Странное поведение при гауссовом размытии изображения с использованием большого радиуса / стандартного отклонения

Question

См редактировать

Я попытался реализовать алгоритм размытия по Гауссу самостоятельно в MATLAB вместо использования встроенного алгоритма из соображений его детального понимания.

Я нашел интересную реализацию, и кто-то уже спрашивал, как кодировать такой алгоритм. Так что это не бизнес.

Кроме того, я использую следующую формулу для вычисления стандартного отклонения заданного радиуса, как это делает GIMP:

stdDeviation = sqrt(-(radius^2) / (2*log10(1 / 255)));

Мой алгоритм работает для небольших значений радиуса (например, 3,5,7) без каких-либо проблем (по крайней мере, вы не видите разницу). Если я попытаюсь размыть изображение с радиусом 21, получится следующее:

По сравнению с GIMP / MATLAB imgaussfilt(A,sigma) выход:

Очевидно, что алгоритмы вычисляют не одно и то же (или, по крайней мере, похожее) выходное изображение. Что такое GIMP / MATLAB? imgaussfilt(A,sigma) кроме этого?

Границу изображения можно пренебречь. Я знаю об этой проблеме. Но я не понимаю происхождение "странных полос пикселей" в моем выходном изображении.

По причинам полноты, вот исходный код:

function y = gaussianBlurSepTest(inputImage)
% radius in pixel; RADIUS MUST BE ODD! 
radius = 21;
% approximate value for standard deviation
stdDeviation = sqrt(-(radius^2) / (2*log10(1 / 255)));

ind = -floor(radius/2):floor(radius/2);
[X, Y] = meshgrid(ind, ind);
h = exp(-(X.^2 + Y.^2) / (2*stdDeviation*stdDeviation));
h = h / sum(h(:));

redChannel = inputImage(:,:,1);
greenChannel = inputImage(:,:,2);
blueChannel = inputImage(:,:,3);

redBlurred = conv2(redChannel, h);
greenBlurred = conv2(greenChannel, h);
blueBlurred = conv2(blueChannel, h);

y = cat(3, uint8(redBlurred), uint8(greenBlurred), uint8(blueBlurred));

РЕДАКТИРОВАТЬ:

Для полноты картины и для помощи другим: я применил модификации Эрфана. Результат теперь намного лучше, но все еще можно увидеть очевидную разницу в вычислениях GIMP. Результат GIMP выглядит более "гладким".

Реализованный алгоритм:

function y = gaussianBlurSepTest(inputImage)
radius = 21;
stdDeviation = sqrt(-(radius^2) / (2*log10(1 / 255)));
ind = -floor(radius/2):floor(radius/2);
[X, Y] = meshgrid(ind, ind);
[~, R] = cart2pol(X, Y);  % Here R is defined
h = exp(-R.^2 / (2*stdDeviation*stdDeviation));
h = h / sum(h(:));
h(R > radius/2) = 0;
h = h / sum(h(:));
redChannel = inputImage(:,:,1);
greenChannel = inputImage(:,:,2);
blueChannel = inputImage(:,:,3);
redBlurred = conv2(redChannel, h);
greenBlurred = conv2(greenChannel, h);
blueBlurred = conv2(blueChannel, h);
y = cat(3, uint8(redBlurred), uint8(greenBlurred), uint8(blueBlurred));

Результат:

Результат GIMP:

С точки зрения полного ответа на вопрос и помощи другим в том же вопросе, я думаю, было бы полезно спросить о причинах разногласий.

Спасибо!

Answer 1

h отвечает за появившиеся горизонтальные и вертикальные полосы. Хотя вы определили h с угловой симметрией, как видно на следующем графике, его границы нарушают эту симметрию:

Чтобы сделать все правильно, вы можете обрезать h на правильном радиусе:

Я применил модификацию к вашей функции, и теперь она должна дать гораздо лучшие результаты:

function y = gaussianBlurSepTest(inputImage, radius)
stdDeviation = sqrt(-(radius^2) / (2*log10(1 / 255)));
ind = -floor(radius/2):floor(radius/2);
[X, Y] = meshgrid(ind, ind);
[~, R] = cart2pol(X, Y);  % Here R is defined
h = exp(-R.^2 / (2*stdDeviation*stdDeviation));
h = h / sum(h(:));
h(R > radius/2) = 0;  % And here h is truncated. The rest is the same.

Вот мой тест. Мое изображение:

После вашего gaussianBlurSepTest (радиус = 35):

После модифицированной функции:

Примечание: вывод немного темнее. Если это проблема, вы можете повторно нормализовать stdDeviation или расширить свою сетку.

Answer 2

Я думаю, что виновник conv2 инструкция, и это потому, что ваш тип изображения uint8 или же uint16, когда вы фильтруете свое изображение, вы должны использовать float или же double тип. Размер вашего фильтра настолько велик, что результат применения окна коэффициентов к пикселям превышает 8-битный или 16-битный целочисленный диапазон (255 или 65535). Я предлагаю попробовать привести ваше изображение к double перед применением вашего фильтра и после применения фильтра верните его uint8 как показано ниже:

   Output=gaussianBlurSepTest(double(inputImage));
   Output=uint8(Output);

Удачи.