2012-04-01 5 views
1

Я должен определить, какой фильтр использовался на случайном изображении - существует ли общий способ обнаружения правильного (гауссовый, prewitt, sobel, average, ...) или было бы разумно кодировать некоторые что-то вроде «грубой силы»?обнаружение фильтра

Я попытался найти его с помощью Matlab, но у меня нет клея, как искать больше эффективно. В настоящий момент это похоже на поиск иглы в стоге. Я также подумал о том, чтобы использовать некоторые скрипты bash и imagemagick, но это было бы для ресурса голодным.

я учил это не было бы проблемой, но это очень много времени, чтобы угадать фильтр и попробовать, как этот

f = fspecial('gaussian', [3 3], 1); 
res = imfilter(orginal, f); 
corr2(res, pic); 
+1

Если у вас есть обработанное изображение и исходное изображение, то это очень просто определить коэффициенты фильтра, и оттуда можно классифицировать ядро ​​фильтра. –

+0

да, у меня есть оба, но как я могу определить коэффициенты? с матлабом? – NaN

+1

Вы можете deconvolve, взяв БПФ обоих изображений, разделите, а затем возьмите обратный БПФ результата. См. http://www.mathworks.com/matlabcentral/fileexchange/5465-fast-deconvolution –

ответ

4

Пусть е быть исходное изображение, г быть фильтруют один, а ч фильтр применяется к е, так что:

f * h = g 

Передача тха t в частотную область:

F.H = G, so H = G/F 

Проблема в том, что инвертирование F ОЧЕНЬ чувствительно к шуму.

Как осуществить это в MATLAB:

close all; 
f = imread('cameraman.tif'); 
[x,y] = size(f); 
figure,imshow(f); 
h = fspecial('motion', 20, 40); % abitrary filter just for testing the algorithm 
F = fft2(f); 
H = fft2(h,x,y); 
G = F.*H; 
g = ifft2(G); % the filtered image 
figure, imshow(g/max(g(:))); 
% Inverting the original image 
epsilon = 10^(-10); 
small_values = find(abs(F)<epsilon); 
F(small_values) = epsilon; 
F_i = ones(x,y)./F; 
H_calculated = G.*F_i; 

h_calculated = ifft2(H_calculated); 

% remove really small values to try to infer the original size of h 
r = sum(h_calculated,1)<epsilon; 
c = sum(h_calculated,2)<epsilon; 
h_real = h_calculated(~r,~c); 

% Calculate error 
% redo the filtering with the found filter 
figure,g_comp = ifft2(fft2(f).*fft2(h_real,x,y)); 
imshow(g_comp/max(g_comp(:))); 
rmse = sqrt(mean(mean((double(g_comp) - double(g)).^2,2),1)) 

редактировать: Просто объяснить эпсилон часть:

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

Математически, что было сделано это:

For all F < epsilon, F = epsilon 
+3

Он также очень чувствителен к спектральным нулям. –

+2

Да, я знаю, вот почему я усекаю его. Для меня реальной проблемой с этим подходом является шум. Но это определенно лучше, чем ничего ... – Castilho