OpenCV: homomorphic filter

Question

I want to use a homomorphic filter to work on an underwater image. I tried to encode it with codes found on the Internet, but I always have a black image ... I tried to normalize my result, but did not work.

Here are my functions:

void HomomorphicFilter::butterworth_homomorphic_filter(Mat &dft_Filter, int D, int n, float high_h_v_TB, float low_h_v_TB) { Mat single(dft_Filter.rows, dft_Filter.cols, CV_32F); Point centre = Point(dft_Filter.rows/2, dft_Filter.cols/2); double radius; float upper = (high_h_v_TB * 0.01); float lower = (low_h_v_TB * 0.01); //create essentially create a butterworth highpass filter //with additional scaling and offset for(int i = 0; i < dft_Filter.rows; i++) { for(int j = 0; j < dft_Filter.cols; j++) { radius = (double) sqrt(pow((i - centre.x), 2.0) + pow((double) (j - centre.y), 2.0)); single.at<float>(i,j) =((upper - lower) * (1/(1 + pow((double) (D/radius), (double) (2*n))))) + lower; } } //normalize(single, single, 0, 1, CV_MINMAX); //Apply filter mulSpectrums( dft_Filter, single, dft_Filter, 0); } void HomomorphicFilter::Shifting_DFT(Mat &fImage) { //For visualization purposes we may also rearrange the quadrants of the result, so that the origin (0,0), corresponds to the image center. Mat tmp, q0, q1, q2, q3; /*First crop the image, if it has an odd number of rows or columns. Operator & bit to bit by -2 (two complement : -2 = 111111111....10) to eliminate the first bit 2^0 (In case of odd number on row or col, we take the even number in below)*/ fImage = fImage(Rect(0, 0, fImage.cols & -2, fImage.rows & -2)); int cx = fImage.cols/2; int cy = fImage.rows/2; /*Rearrange the quadrants of Fourier image so that the origin is at the image center*/ q0 = fImage(Rect(0, 0, cx, cy)); q1 = fImage(Rect(cx, 0, cx, cy)); q2 = fImage(Rect(0, cy, cx, cy)); q3 = fImage(Rect(cx, cy, cx, cy)); /*We reverse each quadrant of the frame with its other quadrant diagonally opposite*/ /*We reverse q0 and q3*/ q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); /*We reverse q1 and q2*/ q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); } void HomomorphicFilter::Fourier_Transform(Mat frame_bw, Mat &image_phase, Mat &image_mag) { Mat frame_log; frame_bw.convertTo(frame_log, CV_32F); /*Take the natural log of the input (compute log(1 + Mag)*/ frame_log += 1; log( frame_log, frame_log); // log(1 + Mag) /*2. Expand the image to an optimal size The performance of the DFT depends of the image size. It tends to be the fastest for image sizes that are multiple of 2, 3 or 5. We can use the copyMakeBorder() function to expand the borders of an image.*/ Mat padded; int M = getOptimalDFTSize(frame_log.rows); int N = getOptimalDFTSize(frame_log.cols); copyMakeBorder(frame_log, padded, 0, M - frame_log.rows, 0, N - frame_log.cols, BORDER_CONSTANT, Scalar::all(0)); /*Make place for both the complex and real values The result of the DFT is a complex. Then the result is 2 images (Imaginary + Real), and the frequency domains range is much larger than the spatial one. Therefore we need to store in float ! That why we will convert our input image "padded" to float and expand it to another channel to hold the complex values. Planes is an arrow of 2 matrix (planes[0] = Real part, planes[1] = Imaginary part)*/ Mat image_planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)}; Mat image_complex; /*Creates one multichannel array out of several single-channel ones.*/ merge(image_planes, 2, image_complex); /*Make the DFT The result of thee DFT is a complex image : "image_complex"*/ dft(image_complex, image_complex); /***************************/ //Create spectrum magnitude// /***************************/ /*Transform the real and complex values to magnitude NB: We separe Real part to Imaginary part*/ split(image_complex, image_planes); //Starting with this part we have the real part of the image in planes[0] and the imaginary in planes[1] phase(image_planes[0], image_planes[1], image_phase); magnitude(image_planes[0], image_planes[1], image_mag); } void HomomorphicFilter::Inv_Fourier_Transform(Mat image_phase, Mat image_mag, Mat &inverseTransform) { /*Calculates x and y coordinates of 2D vectors from their magnitude and angle.*/ Mat result_planes[2]; polarToCart(image_mag, image_phase, result_planes[0], result_planes[1]); /*Creates one multichannel array out of several single-channel ones.*/ Mat result_complex; merge(result_planes, 2, result_complex); /*Make the IDFT*/ dft(result_complex, inverseTransform, DFT_INVERSE|DFT_REAL_OUTPUT); /*Take the exponential*/ exp(inverseTransform, inverseTransform); } 

and here is my main code:

  /**************************/ /****Homomorphic filter****/ /**************************/ /**********************************************/ //Getting the frequency and magnitude of image// /**********************************************/ Mat image_phase, image_mag; HomomorphicFilter().Fourier_Transform(frame_bw, image_phase, image_mag); /******************/ //Shifting the DFT// /******************/ HomomorphicFilter().Shifting_DFT(image_mag); /********************************/ //Butterworth homomorphic filter// /********************************/ int high_h_v_TB = 101; int low_h_v_TB = 99; int D = 10;// radius of band pass filter parameter int order = 2;// order of band pass filter parameter HomomorphicFilter().butterworth_homomorphic_filter(image_mag, D, order, high_h_v_TB, low_h_v_TB); /******************/ //Shifting the DFT// /******************/ HomomorphicFilter().Shifting_DFT(image_mag); /*******************************/ //Inv Discret Fourier Transform// /*******************************/ Mat inverseTransform; HomomorphicFilter().Inv_Fourier_Transform(image_phase, image_mag, inverseTransform); imshow("Result", inverseTransform); 

If someone can explain my mistakes to me, I would really appreciate :). Thanks and sorry for my bad english.

EDIT: Now I have something, but it's not perfect ... I changed 2 things in my code. I applied log (mag + 1) after dft and not on the input image. I removed exp () after idft.

here are the results (I can publish only 2 links ...):

my input image: final result:

having seen several topics, I find similar results for the butterworth filter and for my size after dft / shifting. Unfortunately, my end result is not very good. Why do I have so much "noise"?

Answer 1

I did this method to balance the backlight when changing the camera, thereby causing the image to darken.

I tried FFT on frequency to filter the image! it works. but use too much time (image 2750 * 3680RGB). I do this in the spatial domain.

here is my code!

 //IplImage *imgSrcI=cvLoadImage("E:\\lean.jpg",-1); Mat imgSrcM(imgSrc,true); Mat imgDstM; Mat imgGray; Mat imgHls; vector<Mat> vHls; Mat imgTemp1=Mat::zeros(imgSrcM.size(),CV_64FC1); Mat imgTemp2=Mat::zeros(imgSrcM.size(),CV_64FC1); if(imgSrcM.channels()==1) { imgGray=imgSrcM.clone(); } else if (imgSrcM.channels()==3) { cvtColor(imgSrcM, imgHls, CV_BGR2HLS); split(imgHls, vHls); imgGray=vHls.at(1); } else { return -1; } imgGray.convertTo(imgTemp1,CV_64FC1); imgTemp1=imgTemp1+0.0001; log(imgTemp1,imgTemp1); GaussianBlur(imgTemp1, imgTemp2, Size(21, 21), 0.1, 0.1, BORDER_DEFAULT);//imgTemp2是低通滤波的结果imgTemp1 = (imgTemp1 - imgTemp2);//imgTemp1是对数减低通的高通addWeighted(imgTemp2, 0.7, imgTemp1, 1.4, 1, imgTemp1, -1);//imgTemp1是压制低频增强高频的结构exp(imgTemp1,imgTemp1); normalize(imgTemp1,imgTemp1,0,1,NORM_MINMAX); imgTemp1=imgTemp1*255; imgTemp1.convertTo(imgGray, CV_8UC1); //imwrite("E:\\leanImgGray.jpg",imgGray); if (imgSrcM.channels()==3) { vHls.at(1)=imgGray; merge(vHls,imgHls); cvtColor(imgHls, imgDstM, CV_HLS2BGR); } else if (imgSrcM.channels()==1) { imgDstM=imgGray.clone(); } cvCopy(&(IplImage)imgDstM,imgDst); //cvShowImage("jpg",imgDst); return 0; 

Answer 2

I took your code, adjusted it in several places and got decent results as an output for a homographic filter.

Here are the corrections I made.

1) Instead of working only with image_mag, work on the full output of the FFT.

2) your filter values ​​high_h_v_TB = 101 and low_h_v_TB = 99 had practically no effect on filtering.

Here are the values ​​I used.

 int high_h_v_TB = 100; int low_h_v_TB = 20; int D = 10;// radius of band pass filter parameter int order = 4; 

Here is my main code

 //float_img == grayscale image in 0-1 scale Mat log_img; log(float_img, log_img); Mat fft_phase, fft_mag; Mat fft_complex; HomomorphicFilter::Fourier_Transform(log_img, fft_complex); HomomorphicFilter::ShiftFFT(fft_complex); int high_h_v_TB = 100; int low_h_v_TB = 30; int D = 10;// radius of band pass filter parameter int order = 4; //get a butterworth filter of same image size as the input image //dont call mulSpectrums yet, just get the filter of correct size Mat butterWorthFreqDomain; HomomorphicFilter::ButterworthFilter(fft_complex.size(), butterWorthFreqDomain, D, order, high_h_v_TB, low_h_v_TB); //this should match fft_complex in size and type //and is what we will be using for 'mulSpectrums' call Mat butterworth_complex; //make two channels to match fft_complex Mat butterworth_channels[] = {Mat_<float>(butterWorthFreqDomain.size()), Mat::zeros(butterWorthFreqDomain.size(), CV_32F)}; merge(butterworth_channels, 2, butterworth_complex); //do mulSpectrums on the full fft mulSpectrums(fft_complex, butterworth_complex, fft_complex, 0); //shift back the output HomomorphicFilter::ShiftFFT(fft_complex); Mat log_img_out; HomomorphicFilter::Inv_Fourier_Transform(fft_complex, log_img_out); Mat float_img_out; exp(log_img_out, float_img_out); //float_img_out is gray in 0-1 range 

Here is my conclusion.