C’è un modo per rilevare se un’immagine è sfocata?

Mi stavo chiedendo se c’è un modo per determinare se un’immagine è sfocata o meno analizzando i dati dell’immagine.

sì. calcolare il fft e analizzare il risultato. La trasformata di Fourier ti dice quali frequenze sono presenti nell’immagine. Se c’è una bassa quantità di alte frequenze, l’immagine è sfocata.

La definizione dei termini “basso” e “alto” dipende da te.

modifica : come indicato nei commenti, se vuoi un singolo float che rappresenti la sfocatura di una determinata immagine, devi calcolare una metrica adatta.

la risposta di nikie fornisce una tale metrica. Trasforma l’immagine con un kernel laplaciano:

1 1 -4 1 1 

E utilizzare una metrica massima robusta sull’output per ottenere un numero che è ansible utilizzare per la soglia. Cerca di evitare di lisciare troppo le immagini prima di calcolare il laplaciano, perché scoprirai solo che un’immagine liscia è davvero sfocata :-).

Un altro modo molto semplice per stimare la nitidezza di un’immagine consiste nell’utilizzare un filtro Laplace (o LoG) e selezionare semplicemente il valore massimo. Usare una misura robusta come un quantile del 99,9% è probabilmente meglio se ci si aspetta un rumore (es. Scegliere il N-massimo contrasto invece del più alto contrasto). Se ti aspetti una variazione della luminosità dell’immagine, dovresti anche includere un passo di pre-elaborazione per normalizzare la luminosità dell’immagine / contrasto (eg equalizzazione dell’istogramma).

Ho implementato il suggerimento di Simon e questo in Mathematica, e l’ho provato con alcune immagini di prova:

prova le immagini

Il primo test sfoca le immagini di prova usando un filtro gaussiano con una dimensione del kernel variabile, quindi calcola la FFT dell’immagine sfocata e prende la media delle frequenze più alte del 90%:

 testFft[img_] := Table[ ( blurred = GaussianFilter[img, r]; fft = Fourier[ImageData[blurred]]; {w, h} = Dimensions[fft]; windowSize = Round[w/2.1]; Mean[Flatten[(Abs[ fft[[w/2 - windowSize ;; w/2 + windowSize, h/2 - windowSize ;; h/2 + windowSize]]])]] ), {r, 0, 10, 0.5}] 

Risultato in una trama logaritmica:

risultato fft

Le 5 linee rappresentano le 5 immagini di prova, l’asse X rappresenta il raggio del filtro gaussiano. I grafici stanno diminuendo, quindi la FFT è una buona misura per la nitidezza.

Questo è il codice per lo stimatore di sfocatura “high LoG”: applica semplicemente un filtro LoG e restituisce il pixel più luminoso nel risultato del filtro:

 testLaplacian[img_] := Table[ ( blurred = GaussianFilter[img, r]; Max[Flatten[ImageData[LaplacianGaussianFilter[blurred, 1]]]]; ), {r, 0, 10, 0.5}] 

Risultato in una trama logaritmica:

risultato di laplace

Lo spread per le immagini non sfocate è un po ‘migliore qui (2.5 vs 3.3), principalmente perché questo metodo utilizza solo il contrasto più forte nell’immagine, mentre la FFT è essenzialmente una media sull’intera immagine. Anche le funzioni diminuiscono più velocemente, quindi potrebbe essere più semplice impostare una soglia “sfocata”.

Durante alcuni lavori con un objective con messa a fuoco automatica, mi sono imbattuto in questo utile insieme di algoritmi per il rilevamento della messa a fuoco dell’immagine . È implementato in MATLAB, ma la maggior parte delle funzioni è abbastanza facile da portare a OpenCV con filter2D .

Si tratta fondamentalmente di un’implementazione di sondaggi di molti algoritmi di misurazione della messa a fuoco. Se si desidera leggere i documenti originali, i riferimenti agli autori degli algoritmi sono forniti nel codice. Il documento del 2012 di Pertuz, et al. L’analisi degli operatori di misura della messa a fuoco per la forma a fuoco (SFF) fornisce una buona rassegna di tutte queste misure e delle loro prestazioni (sia in termini di velocità che di accuratezza applicate all’SFF).

EDIT: aggiunto il codice MATLAB nel caso in cui il link muoia.

 function FM = fmeasure(Image, Measure, ROI) %This function measures the relative degree of focus of %an image. It may be invoked as: % % FM = fmeasure(Image, Method, ROI) % %Where % Image, is a grayscale image and FM is the computed % focus value. % Method, is the focus measure algorithm as a string. % see 'operators.txt' for a list of focus % measure methods. % ROI, Image ROI as a rectangle [xo yo width heigth]. % if an empty argument is passed, the whole % image is processed. % % Said Pertuz % Abr/2010 if ~isempty(ROI) Image = imcrop(Image, ROI); end WSize = 15; % Size of local window (only some operators) switch upper(Measure) case 'ACMO' % Absolute Central Moment (Shirvaikar2004) if ~isinteger(Image), Image = im2uint8(Image); end FM = AcMomentum(Image); case 'BREN' % Brenner's (Santos97) [MN] = size(Image); DH = Image; DV = Image; DH(1:M-2,:) = diff(Image,2,1); DV(:,1:N-2) = diff(Image,2,2); FM = max(DH, DV); FM = FM.^2; FM = mean2(FM); case 'CONT' % Image contrast (Nanda2001) ImContrast = inline('sum(abs(x(:)-x(5)))'); FM = nlfilter(Image, [3 3], ImContrast); FM = mean2(FM); case 'CURV' % Image Curvature (Helmli2001) if ~isinteger(Image), Image = im2uint8(Image); end M1 = [-1 0 1;-1 0 1;-1 0 1]; M2 = [1 0 1;1 0 1;1 0 1]; P0 = imfilter(Image, M1, 'replicate', 'conv')/6; P1 = imfilter(Image, M1', 'replicate', 'conv')/6; P2 = 3*imfilter(Image, M2, 'replicate', 'conv')/10 ... -imfilter(Image, M2', 'replicate', 'conv')/5; P3 = -imfilter(Image, M2, 'replicate', 'conv')/5 ... +3*imfilter(Image, M2, 'replicate', 'conv')/10; FM = abs(P0) + abs(P1) + abs(P2) + abs(P3); FM = mean2(FM); case 'DCTE' % DCT energy ratio (Shen2006) FM = nlfilter(Image, [8 8], @DctRatio); FM = mean2(FM); case 'DCTR' % DCT reduced energy ratio (Lee2009) FM = nlfilter(Image, [8 8], @ReRatio); FM = mean2(FM); case 'GDER' % Gaussian derivative (Geusebroek2000) N = floor(WSize/2); sig = N/2.5; [x,y] = meshgrid(-N:N, -N:N); G = exp(-(x.^2+y.^2)/(2*sig^2))/(2*pi*sig); Gx = -x.*G/(sig^2);Gx = Gx/sum(Gx(:)); Gy = -y.*G/(sig^2);Gy = Gy/sum(Gy(:)); Rx = imfilter(double(Image), Gx, 'conv', 'replicate'); Ry = imfilter(double(Image), Gy, 'conv', 'replicate'); FM = Rx.^2+Ry.^2; FM = mean2(FM); case 'GLVA' % Graylevel variance (Krotkov86) FM = std2(Image); case 'GLLV' %Graylevel local variance (Pech2000) LVar = stdfilt(Image, ones(WSize,WSize)).^2; FM = std2(LVar)^2; case 'GLVN' % Normalized GLV (Santos97) FM = std2(Image)^2/mean2(Image); case 'GRAE' % Energy of gradient (Subbarao92a) Ix = Image; Iy = Image; Iy(1:end-1,:) = diff(Image, 1, 1); Ix(:,1:end-1) = diff(Image, 1, 2); FM = Ix.^2 + Iy.^2; FM = mean2(FM); case 'GRAT' % Thresholded gradient (Snatos97) Th = 0; %Threshold Ix = Image; Iy = Image; Iy(1:end-1,:) = diff(Image, 1, 1); Ix(:,1:end-1) = diff(Image, 1, 2); FM = max(abs(Ix), abs(Iy)); FM(FMImage); FM = 1./R1; FM(index) = R1(index); FM = mean2(FM); case 'HISE' % Histogram entropy (Krotkov86) FM = entropy(Image); case 'HISR' % Histogram range (Firestone91) FM = max(Image(:))-min(Image(:)); case 'LAPE' % Energy of laplacian (Subbarao92a) LAP = fspecial('laplacian'); FM = imfilter(Image, LAP, 'replicate', 'conv'); FM = mean2(FM.^2); case 'LAPM' % Modified Laplacian (Nayar89) M = [-1 2 -1]; Lx = imfilter(Image, M, 'replicate', 'conv'); Ly = imfilter(Image, M', 'replicate', 'conv'); FM = abs(Lx) + abs(Ly); FM = mean2(FM); case 'LAPV' % Variance of laplacian (Pech2000) LAP = fspecial('laplacian'); ILAP = imfilter(Image, LAP, 'replicate', 'conv'); FM = std2(ILAP)^2; case 'LAPD' % Diagonal laplacian (Thelen2009) M1 = [-1 2 -1]; M2 = [0 0 -1;0 2 0;-1 0 0]/sqrt(2); M3 = [-1 0 0;0 2 0;0 0 -1]/sqrt(2); F1 = imfilter(Image, M1, 'replicate', 'conv'); F2 = imfilter(Image, M2, 'replicate', 'conv'); F3 = imfilter(Image, M3, 'replicate', 'conv'); F4 = imfilter(Image, M1', 'replicate', 'conv'); FM = abs(F1) + abs(F2) + abs(F3) + abs(F4); FM = mean2(FM); case 'SFIL' %Steerable filters (Minhas2009) % Angles = [0 45 90 135 180 225 270 315]; N = floor(WSize/2); sig = N/2.5; [x,y] = meshgrid(-N:N, -N:N); G = exp(-(x.^2+y.^2)/(2*sig^2))/(2*pi*sig); Gx = -x.*G/(sig^2);Gx = Gx/sum(Gx(:)); Gy = -y.*G/(sig^2);Gy = Gy/sum(Gy(:)); R(:,:,1) = imfilter(double(Image), Gx, 'conv', 'replicate'); R(:,:,2) = imfilter(double(Image), Gy, 'conv', 'replicate'); R(:,:,3) = cosd(45)*R(:,:,1)+sind(45)*R(:,:,2); R(:,:,4) = cosd(135)*R(:,:,1)+sind(135)*R(:,:,2); R(:,:,5) = cosd(180)*R(:,:,1)+sind(180)*R(:,:,2); R(:,:,6) = cosd(225)*R(:,:,1)+sind(225)*R(:,:,2); R(:,:,7) = cosd(270)*R(:,:,1)+sind(270)*R(:,:,2); R(:,:,7) = cosd(315)*R(:,:,1)+sind(315)*R(:,:,2); FM = max(R,[],3); FM = mean2(FM); case 'SFRQ' % Spatial frequency (Eskicioglu95) Ix = Image; Iy = Image; Ix(:,1:end-1) = diff(Image, 1, 2); Iy(1:end-1,:) = diff(Image, 1, 1); FM = mean2(sqrt(double(Iy.^2+Ix.^2))); case 'TENG'% Tenengrad (Krotkov86) Sx = fspecial('sobel'); Gx = imfilter(double(Image), Sx, 'replicate', 'conv'); Gy = imfilter(double(Image), Sx', 'replicate', 'conv'); FM = Gx.^2 + Gy.^2; FM = mean2(FM); case 'TENV' % Tenengrad variance (Pech2000) Sx = fspecial('sobel'); Gx = imfilter(double(Image), Sx, 'replicate', 'conv'); Gy = imfilter(double(Image), Sx', 'replicate', 'conv'); G = Gx.^2 + Gy.^2; FM = std2(G)^2; case 'VOLA' % Vollath's correlation (Santos97) Image = double(Image); I1 = Image; I1(1:end-1,:) = Image(2:end,:); I2 = Image; I2(1:end-2,:) = Image(3:end,:); Image = Image.*(I1-I2); FM = mean2(Image); case 'WAVS' %Sum of Wavelet coeffs (Yang2003) [C,S] = wavedec2(Image, 1, 'db6'); H = wrcoef2('h', C, S, 'db6', 1); V = wrcoef2('v', C, S, 'db6', 1); D = wrcoef2('d', C, S, 'db6', 1); FM = abs(H) + abs(V) + abs(D); FM = mean2(FM); case 'WAVV' %Variance of Wav...(Yang2003) [C,S] = wavedec2(Image, 1, 'db6'); H = abs(wrcoef2('h', C, S, 'db6', 1)); V = abs(wrcoef2('v', C, S, 'db6', 1)); D = abs(wrcoef2('d', C, S, 'db6', 1)); FM = std2(H)^2+std2(V)+std2(D); case 'WAVR' [C,S] = wavedec2(Image, 3, 'db6'); H = abs(wrcoef2('h', C, S, 'db6', 1)); V = abs(wrcoef2('v', C, S, 'db6', 1)); D = abs(wrcoef2('d', C, S, 'db6', 1)); A1 = abs(wrcoef2('a', C, S, 'db6', 1)); A2 = abs(wrcoef2('a', C, S, 'db6', 2)); A3 = abs(wrcoef2('a', C, S, 'db6', 3)); A = A1 + A2 + A3; WH = H.^2 + V.^2 + D.^2; WH = mean2(WH); WL = mean2(A); FM = WH/WL; otherwise error('Unknown measure %s',upper(Measure)) end end %************************************************************************ function fm = AcMomentum(Image) [MN] = size(Image); Hist = imhist(Image)/(M*N); Hist = abs((0:255)-255*mean2(Image))'.*Hist; fm = sum(Hist); end %****************************************************************** function fm = DctRatio(M) MT = dct2(M).^2; fm = (sum(MT(:))-MT(1,1))/MT(1,1); end %************************************************************************ function fm = ReRatio(M) M = dct2(M); fm = (M(1,2)^2+M(1,3)^2+M(2,1)^2+M(2,2)^2+M(3,1)^2)/(M(1,1)^2); end %****************************************************************** 

Alcuni esempi di versioni OpenCV:

 // OpenCV port of 'LAPM' algorithm (Nayar89) double modifiedLaplacian(const cv::Mat& src) { cv::Mat M = (Mat_(3, 1) < < -1, 2, -1); cv::Mat G = cv::getGaussianKernel(3, -1, CV_64F); cv::Mat Lx; cv::sepFilter2D(src, Lx, CV_64F, M, G); cv::Mat Ly; cv::sepFilter2D(src, Ly, CV_64F, G, M); cv::Mat FM = cv::abs(Lx) + cv::abs(Ly); double focusMeasure = cv::mean(FM).val[0]; return focusMeasure; } // OpenCV port of 'LAPV' algorithm (Pech2000) double varianceOfLaplacian(const cv::Mat& src) { cv::Mat lap; cv::Laplacian(src, lap, CV_64F); cv::Scalar mu, sigma; cv::meanStdDev(lap, mu, sigma); double focusMeasure = sigma.val[0]*sigma.val[0]; return focusMeasure; } // OpenCV port of 'TENG' algorithm (Krotkov86) double tenengrad(const cv::Mat& src, int ksize) { cv::Mat Gx, Gy; cv::Sobel(src, Gx, CV_64F, 1, 0, ksize); cv::Sobel(src, Gy, CV_64F, 0, 1, ksize); cv::Mat FM = Gx.mul(Gx) + Gy.mul(Gy); double focusMeasure = cv::mean(FM).val[0]; return focusMeasure; } // OpenCV port of 'GLVN' algorithm (Santos97) double normalizedGraylevelVariance(const cv::Mat& src) { cv::Scalar mu, sigma; cv::meanStdDev(src, mu, sigma); double focusMeasure = (sigma.val[0]*sigma.val[0]) / mu.val[0]; return focusMeasure; } 

Nessuna garanzia sul fatto che queste misure siano la scelta migliore per il tuo problema, ma se rintracci i documenti associati a queste misure, potrebbero darti più informazioni. Spero che trovi utile il codice! So che l'ho fatto.

Costruire la risposta di Nike. È semplice implementare il metodo basato su laplacian con opencv:

 short GetSharpness(char* data, unsigned int width, unsigned int height) { // assumes that your image is already in planner yuv or 8 bit greyscale IplImage* in = cvCreateImage(cvSize(width,height),IPL_DEPTH_8U,1); IplImage* out = cvCreateImage(cvSize(width,height),IPL_DEPTH_16S,1); memcpy(in->imageData,data,width*height); // aperture size of 1 corresponds to the correct matrix cvLaplace(in, out, 1); short maxLap = -32767; short* imgData = (short*)out->imageData; for(int i =0;i< (out->imageSize/2);i++) { if(imgData[i] > maxLap) maxLap = imgData[i]; } cvReleaseImage(&in); cvReleaseImage(&out); return maxLap; } 

Restituirà un breve che indica la massima nitidezza rilevata, che in base ai miei test su campioni del mondo reale, è un buon indicatore se una fotocamera è a fuoco o meno. Non sorprendentemente, i valori normali dipendono dalla scena, ma molto meno dal metodo FFT, che deve essere elevato a una percentuale di falsi positivi per essere utile nella mia applicazione.

Ho trovato una soluzione completamente diversa. Avevo bisogno di analizzare i fotogrammi video per trovare il più nitido in ogni fotogramma (X). In questo modo, avrei rilevato il motion blur e / o le immagini fuori fuoco.

Ho finito per utilizzare il rilevamento Canny Edge e ho ottenuto risultati MOLTO MOLTO buoni con quasi ogni tipo di video (con il metodo di Nikie, ho avuto problemi con video VHS digitalizzati e video interlacciati pesanti).

Ho ottimizzato le prestazioni impostando una regione di interesse (ROI) sull’immagine originale.

Utilizzando EmguCV:

 //Convert image using Canny using (Image imgCanny = imgOrig.Canny(225, 175)) { //Count the number of pixel representing an edge int nCountCanny = imgCanny.CountNonzero()[0]; //Compute a sharpness grade: //< 1.5 = blurred, in movement //de 1.5 à 6 = acceptable //> 6 =stable, sharp double dSharpness = (nCountCanny * 1000.0 / (imgCanny.Cols * imgCanny.Rows)); } 

Grazie nikie per il grande suggerimento di Laplace. I documenti OpenCV mi hanno indirizzato nella stessa direzione: usando python, cv2 (opencv 2.4.10) e numpy …

gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY) numpy.max(cv2.convertScaleAbs(cv2.Laplacian(gray_image,3)))

il risultato è tra 0-255. Ho trovato qualcosa di più di 200 è molto a fuoco, e per 100, è notevolmente sfocato. il massimo non ottiene mai molto sotto i 20 anche se è completamente sfocato.

Un modo che sto attualmente usando misura la diffusione dei bordi nell’immagine. Cerca questo documento:

 @ARTICLE{Marziliano04perceptualblur, author = {Pina Marziliano and Frederic Dufaux and Stefan Winkler and Touradj Ebrahimi}, title = {Perceptual blur and ringing metrics: Application to JPEG2000,” Signal Process}, journal = {Image Commun}, year = {2004}, pages = {163--172} } 

Di solito è dietro un paywall ma ho visto alcune copie gratuite in giro. Fondamentalmente, localizzano i bordi verticali in un’immagine e poi misurano la larghezza di questi bordi. La media della larghezza fornisce il risultato finale della sfocatura dell’immagine. I bordi più larghi corrispondono a immagini sfocate e viceversa.

Questo problema appartiene al campo della stima della qualità dell’immagine senza riferimento . Se lo cerchi su Google Scholar, avrai molti riferimenti utili.

MODIFICARE

Ecco una trama delle stime di sfocatura ottenute per le 5 immagini nel post di nikie. Valori più alti corrispondono a maggiore sfocatura. Ho usato un filtro gaussiano 11×11 a dimensione fissa e ho variato la deviazione standard (usando il comando convert di imagemagick per ottenere le immagini sfocate).

inserisci la descrizione dell'immagine qui

Se confronti immagini di dimensioni diverse, non dimenticare di normalizzare la larghezza dell’immagine, poiché le immagini più grandi avranno bordi più larghi.

Infine, un problema significativo è la distinzione tra sfocatura artistica e sfocatura indesiderata (causata da mancanza di messa a fuoco, compressione, movimento relativo del sobject alla fotocamera), ma ciò va oltre i semplici approcci come questo. Per un esempio di sfocatura artistica, dai un’occhiata all’immagine di Lenna: il riflesso di Lenna nello specchio è sfocato, ma il suo volto è perfettamente a fuoco. Ciò contribuisce a una stima di sfocatura più elevata per l’immagine di Lenna.

Le risposte sopra hanno chiarito molte cose, ma penso che sia utile fare una distinzione concettuale.

Cosa fare se si riprende un’immagine perfettamente a fuoco di un’immagine sfocata?

Il problema di rilevamento della sfocatura è ben posto quando si ha un riferimento . Se è necessario progettare, ad esempio, un sistema di messa a fuoco automatica, si confronta una sequenza di immagini riprese con diversi gradi di sfocatura o levigatura e si tenta di trovare il punto di sfocatura minima all’interno di questo set. In altre parole è necessario incrociare le varie immagini usando una delle tecniche illustrate sopra (fondamentalmente – con vari livelli possibili di raffinatezza nell’approccio – cercando l’unica immagine con il più alto contenuto di alta frequenza).

Ho provato una soluzione basata sul filtro laplaciano da questo post. Non mi ha aiutato. Quindi, ho provato la soluzione da questo post ed è andata bene per il mio caso (ma è lento):

 import cv2 image = cv2.imread("test.jpeg") height, width = image.shape[:2] gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY) def px(x, y): return int(gray[y, x]) sum = 0 for x in range(width-1): for y in range(height): sum += abs(px(x, y) - px(x+1, y)) 

L’immagine meno sfocata ha il valore di sum massimo!

È inoltre ansible sintonizzare la velocità e la precisione modificando il passo, ad es

questa parte

 for x in range(width - 1): 

puoi sostituire con questo

 for x in range(0, width - 1, 10): 

Il codice Matlab di due metodi che sono stati pubblicati su riviste molto apprezzate (IEEE Transactions on Image Processing) sono disponibili qui: https://ivulab.asu.edu/software

controllare gli algoritmi CPBDM e JNBM. Se si controlla il codice non è molto difficile essere trasferito e, incidentalmente, si basa sul metodo di Marzialiano come caratteristica di base.

l’ho implementato utilizzare fft in MATLAB e controllare l’istogramma della media fft calcolo e std ma anche la funzione adatta può essere fatta

 fa = abs(fftshift(fft(sharp_img))); fb = abs(fftshift(fft(blured_img))); f1=20*log10(0.001+fa); f2=20*log10(0.001+fb); figure,imagesc(f1);title('org') figure,imagesc(f2);title('blur') figure,hist(f1(:),100);title('org') figure,hist(f2(:),100);title('blur') mf1=mean(f1(:)); mf2=mean(f2(:)); mfd1=median(f1(:)); mfd2=median(f2(:)); sf1=std(f1(:)); sf2=std(f2(:)); 

Questo è quello che faccio in Opencv per rilevare la qualità della messa a fuoco in una regione:

 Mat grad; int scale = 1; int delta = 0; int ddepth = CV_8U; Mat grad_x, grad_y; Mat abs_grad_x, abs_grad_y; /// Gradient X Sobel(matFromSensor, grad_x, ddepth, 1, 0, 3, scale, delta, BORDER_DEFAULT); /// Gradient Y Sobel(matFromSensor, grad_y, ddepth, 0, 1, 3, scale, delta, BORDER_DEFAULT); convertScaleAbs(grad_x, abs_grad_x); convertScaleAbs(grad_y, abs_grad_y); addWeighted(abs_grad_x, 0.5, abs_grad_y, 0.5, 0, grad); cv::Scalar mu, sigma; cv::meanStdDev(grad, /* mean */ mu, /*stdev*/ sigma); focusMeasure = mu.val[0] * mu.val[0];