mrgloom 242 Жалоба Опубликовано January 23, 2012 да даже если полным перебором пойду, увеличивая тэмплейт попиксельно по x,y получу скажем данные,сравнимы ли будут найденные максимумы или нет? и еще можно ли основное изображение расширять по бокам черными полосками? п.с. при чем тут ранзак? или имеется ввиду сифт и т.д.? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано February 17, 2012 не могу найти этот документ "Kuglin,C.D. and Hines,D.C. (1975) The phase correlation image alignment method. In Proceedings of the IEEE, International Conference on Cybernetics and Society, pp. 163–165." вот еще по поводу log polar координат, интересно падает ли точность при таком методе или нет. http://stackoverflow.com/questions/4814215/template-matching-with-rotation/4814574#4814574 Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано February 17, 2012 Kuglin,C.D. and Hines,D.C. (1975) - труд древний, найти сложно, но можно добыть фрагменты, например здесь: http://www.rec.ri.cmu.edu/fsr09/papers/FSR2009_0007_e5691a8c07c4146597ec6f4b92004434.pdf Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано March 28, 2012 не пойму как оценить точность. например если я беру изображения целиком 100% размера получаю координаты одно относительно другого х=-2; y= 2512; если размера 10%, то х=1; y=251; -> *10 -> х=10; y=2510; т.е. если числа большие, то погрешность даже при скейле в 10 раз всего 2 пикселя, но по х довольно большая получается относительная. я вообще изначально думал, что уменьшаем в n раз и точность падает в n раз. т.е. если например на увеличении 100% - получил координату 1000 то на увеличении 50% могу получить в диапазоне [499-501] на 25% 250+-4 и т.д. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано April 3, 2012 оказывается еще можно с помощью cvCalcBackProjectPatch т.е. с помощью сравнения гистограмм везде примерно одно и то же и на японском http://blog.csdn.net/fdl19881/article/details/6726438 http://opencv.jp/opencv2-x-samples/backprojectionpatch http://opencv.jp/sample/histogram.html наверно работает быстрее, но менее точно. +для быстрого подсчета гистограмм можно использовать интегральное представление изображений. хз как сделано в opencv 1 Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 14, 2012 нашел метод инвариантный к повороту http://opencv-users.1802565.n2.nabble.com/Template-Matching-Issue-td2862903.html как приспособить для масштаба пока не понял. но там пишут I only care about rotation at this point, which shows up in my y axis of the inverse fft'ed phase correlation, but finding the scale isn't hard, its on the x axis, and finding translation would mean going back out of log-polar and getting the phase correlation of the fft'ed images, which wouldn't be hard given this code. еще там есть параметр M как его подбирать на автомате тоже неясно.(возможно зависит от размера картинки?) ошибается где то на 0.2-0.4 градуса, возможно если правильно подбирать параметры, то можно как то эту ситуацию улучшить. double logpolar_fft(IplImage* img,IplImage* temp) { CvSize imgSize = cvSize(img->width, img->height); IplImage* src_Polar = cvCreateImage( imgSize, img->depth, 1 ); IplImage* temp_Polar = cvCreateImage( imgSize, img->depth, 1 ); // The log polar transform in OpenCV CvPoint2D32f pntCenter = cvPoint2D32f(imgSize.width/2, imgSize.height);//почему такой центр? double M = 120.0;//возможно зависит от размера изображения? // Not really sure what the best value of M should be, but this works well cvZero(src_Polar); cvZero(temp_Polar); cvLogPolar(img,src_Polar,pntCenter,M,CV_INTER_LINEAR+CV_WARP_FILL_OUTLIERS); cvLogPolar(temp,temp_Polar,pntCenter,M,CV_INTER_LINEAR+CV_WARP_FILL_OUTLIERS); //test cvSaveImage("src_Polar.png",src_Polar); cvSaveImage("temp_Polar.png",temp_Polar); // Allocate floating point frames used for DFT (real, imaginary, and complex) IplImage* realInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 1 ); IplImage* imaginaryInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 1 ); IplImage* complexInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 2 ); int nDFTHeight= cvGetOptimalDFTSize( imgSize.height ); int nDFTWidth= cvGetOptimalDFTSize( imgSize.width ); CvMat* src_DFT = cvCreateMat( nDFTHeight, nDFTWidth, CV_64FC2 ); CvMat* temp_DFT = cvCreateMat( nDFTHeight, nDFTWidth, CV_64FC2 ); CvSize dftSize = cvSize(nDFTWidth, nDFTHeight); IplImage* imageRe = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageIm = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageImMag = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageMag = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); CvMat tmp; // Processing of src cvScale(src_Polar,realInput,1.0,0);//просто копирование? cvZero(imaginaryInput); cvMerge(realInput,imaginaryInput,NULL,NULL,complexInput); cvGetSubRect(src_DFT,&tmp,cvRect(0,0,src_Polar->width,src_Polar->height)); cvCopy(complexInput,&tmp,NULL); if (src_DFT->cols>src_Polar->width) { cvGetSubRect(src_DFT,&tmp,cvRect(src_Polar->width,0,src_DFT->cols-src_Polar->width,src_Polar->height)); cvZero(&tmp); } cvDFT(src_DFT,src_DFT,CV_DXT_FORWARD,complexInput->height); cvSplit(src_DFT,imageRe,imageIm,0,0); // Processing of temp cvScale(temp_Polar,realInput,1.0,0); cvMerge(realInput,imaginaryInput,NULL,NULL,complexInput); cvGetSubRect(temp_DFT,&tmp,cvRect(0,0,temp_Polar->width,temp_Polar->height)); cvCopy(complexInput,&tmp,NULL); if (temp_DFT->cols>temp_Polar->width) { cvGetSubRect(temp_DFT,&tmp,cvRect(temp_Polar->width,0,temp_DFT->cols-temp_Polar->width,temp_Polar->height)); cvZero( &tmp ); } cvDFT(temp_DFT,temp_DFT,CV_DXT_FORWARD,complexInput->height); // Multiply spectrums of the scene and the model (use CV_DXT_MUL_CONJ to get correlation instead of convolution) cvMulSpectrums(src_DFT,temp_DFT,src_DFT,CV_DXT_MUL_CONJ); // Split Fourier in real and imaginary parts cvSplit(src_DFT,imageRe,imageIm,0,0); // Compute the magnitude of the spectrum components: Mag = sqrt(Re^2 + Im^2) cvPow( imageRe, imageMag, 2.0 ); cvPow( imageIm, imageImMag, 2.0 ); cvAdd( imageMag, imageImMag, imageMag, NULL ); cvPow( imageMag, imageMag, 0.5 ); // Normalize correlation (Divide real and imaginary components by magnitude) cvDiv(imageRe,imageMag,imageRe,1.0); cvDiv(imageIm,imageMag,imageIm,1.0); cvMerge(imageRe,imageIm,NULL,NULL,src_DFT); // inverse dft cvDFT( src_DFT, src_DFT, CV_DXT_INVERSE_SCALE, complexInput->height ); cvSplit( src_DFT, imageRe, imageIm, 0, 0 ); double min = 0.0; double max = 0.0; CvPoint minloc; CvPoint maxloc; cvMinMaxLoc(imageRe,&min,&max,&minloc,&maxloc,NULL); // check Max location and increase or decrease angle from it int x= maxloc.x; // log range if (x>(imageRe->width/2)) x = x-imageRe->width; // positive or negative values int y=maxloc.y; // angle if (y>(imageRe->height/2)) y = y-imageRe->height; // positive or negative values // The correct correlation angle for log-polar depends on the ratio of the length of the y axis // with 360 degrees in 1 complete revolution double angle= (double)y / ((double)imageRe->height/360.0); return angle; } 1 Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 14, 2012 только непонятно как они это делают вот еще док. http://www-cs.engr.ccny.cuny.edu/~wolberg/pub/icip00.pdf исходные данные вот как должно быть вот что получается в качестве центра брал img->width/2, img->height , проецировал в изображением размером с меньшее изображение, M=20 для обоих. double M = 20.0; CvPoint2D32f pntCenter1 = cvPoint2D32f(img->width/2, img->height/2); CvPoint2D32f pntCenter2 = cvPoint2D32f(temp->width/2, temp->height/2); cvLogPolar(img,src_Polar,pntCenter1,M,CV_INTER_LINEAR+CV_WARP_FILL_OUTLIERS); cvLogPolar(temp,temp_Polar,pntCenter2,M,CV_INTER_LINEAR+CV_WARP_FILL_OUTLIERS); cvSaveImage("src_Polar.png",src_Polar); cvSaveImage("temp_Polar.png",temp_Polar); Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 14, 2012 т.е. проблема мне кажется в том, что изображения должны переводится в логполар из одной точки, т.е. например от уха мужика. и действительно такое предположение очень похоже на правду, только оконтовка не совпадает. это получается, если есть сдвиг то надо считать это всё для каждой точки. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 14, 2012 Мне почему-то кажется, что центр в углу изображения, тогда угловая координата меняется 0-90 градусов. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 14, 2012 и почему в углу? попробовал (0,0) ничего хорошего. получается мы получаем инвариантность к повороту и масштабу, но пробелема с выбором точки отсчёта, т.е. проблема со сдвигом. вроде где то предлагалось брать сначало фурье, а потом его только уже в лог-полар переводить возможно ли такое? и поможет ли? апдейт. нашел вот такую схему.надо попробовать. картинка из An FFT-Based Technique for Translation, Rotation, and Scale-Invariant Image Registration Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 14, 2012 Книжка по поиску шаблонов: набираем в гугле Template Matching Techniques in Computer Vision: Theory and Practice pdf на первой странице ссылка на китайский сайт, с него и берем. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 15, 2012 а что там конкретно? я там связки fft+logpolar не нашел есть только 2.4.3. LOG-POLAR MAPPING Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 15, 2012 Да нет, это не по fft+logpolar, просто к теме. По logpolar+fft вроде есть работа "Rotation-discriminating template matching based on Fourier coefficients of radial projections with robustness to scaling and partial occlusion" но в открытом доступе не нашел. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 15, 2012 что то cvLogPolar не хочет работать с IPL_DEPTH_64F и такой вопрос допустим мы сначала захотели сделать фурье преобразование, сделали, получили complex=Re+Im, если мы хотим потом применить к этому logpolar то надо брать Re часть? или это неправильно? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 15, 2012 попробовал FFT->LOGPOLAR пишут, что This method applies a Fourier transform to images to recover translation. Then a log-polar transformation is applied to the magnitude spectrum and the rotation and scale is recovered by using phase correlation in the log-polar space. This method exploits the fact that by operating on the magnitude spectrum of an image, the translational differences are avoided since the magnitude spectrum of an image and its translated counterpart are identical; only their phase spectrum are different. Furthermore, the log-polar transformation causes rotation and scale to be manifest as translation, whereby phase correlation can be applied to recover the rotation angle and scale factor between the pair of input images. The problem here, though, is that limited scale factors can be determined because large scale factors would alter the frequency content beyond recognition. It should be noted that the maximum scale factor recovered in [6] and [7] is 2.0 and 1.8, respectively. только что такое magnitude spectrum? вроде бы даже получилось, единственное смущает, что похоже лучше всего накладываются если просто совместить изображения без сдвига. альтернативная версия(другой центр) тут уже вроде небольшой сдвиг есть. возможно нужно еще этот хайпасс фильтр добавить(я так понимаю это эквивалентно выделению контуров в пространсвенной области) а как это делается в частотной? т.е. когда у нас есть Re+Im после фурье преобразования? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 15, 2012 попытался сравнить код через старый opencv и код через FFTW+opencv похоже что оба работают правильно, только постоянно какие то непонятные доп. смещения по краям(в коде они закомментированы прибавляются или вычитаются от конечной точки), т.е. я так понимаю в зависимости от того в какой квадрант попадает пик корреляции зависит какая будет формула расчета конечного результата, причем для 2-х методов формулы вроде бы будут разные. эти коды я использую не для поиска темплейта, а для сшивки пары изображений.изображения одного размера. вообщем вопрос стоит так как по точке выдаваемой из функции корреляции получить точку смещения? class Peak { public: CvPoint pt; double maxval; }; Peak old_opencv_FFT(IplImage* src,IplImage* temp) { CvSize imgSize = cvSize(src->width, src->height); // Allocate floating point frames used for DFT (real, imaginary, and complex) IplImage* realInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 1 ); IplImage* imaginaryInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 1 ); IplImage* complexInput = cvCreateImage( imgSize, IPL_DEPTH_64F, 2 ); int nDFTHeight= cvGetOptimalDFTSize( imgSize.height ); int nDFTWidth= cvGetOptimalDFTSize( imgSize.width ); CvMat* src_DFT = cvCreateMat( nDFTHeight, nDFTWidth, CV_64FC2 ); CvMat* temp_DFT = cvCreateMat( nDFTHeight, nDFTWidth, CV_64FC2 ); CvSize dftSize = cvSize(nDFTWidth, nDFTHeight); IplImage* imageRe = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageIm = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageImMag = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); IplImage* imageMag = cvCreateImage( dftSize, IPL_DEPTH_64F, 1 ); CvMat tmp; // Processing of src cvScale(src,realInput,1.0,0); cvZero(imaginaryInput); cvMerge(realInput,imaginaryInput,NULL,NULL,complexInput); cvGetSubRect(src_DFT,&tmp,cvRect(0,0,src->width,src->height)); cvCopy(complexInput,&tmp,NULL); if (src_DFT->cols>src->width) { cvGetSubRect(src_DFT,&tmp,cvRect(src->width,0,src_DFT->cols-src->width,src->height)); cvZero(&tmp); } cvDFT(src_DFT,src_DFT,CV_DXT_FORWARD,complexInput->height); cvSplit(src_DFT,imageRe,imageIm,0,0); // Processing of temp cvScale(temp,realInput,1.0,0); cvMerge(realInput,imaginaryInput,NULL,NULL,complexInput); cvGetSubRect(temp_DFT,&tmp,cvRect(0,0,temp->width,temp->height)); cvCopy(complexInput,&tmp,NULL); if (temp_DFT->cols>temp->width) { cvGetSubRect(temp_DFT,&tmp,cvRect(temp->width,0,temp_DFT->cols-temp->width,temp->height)); cvZero( &tmp ); } cvDFT(temp_DFT,temp_DFT,CV_DXT_FORWARD,complexInput->height); // Multiply spectrums of the scene and the model (use CV_DXT_MUL_CONJ to get correlation instead of convolution) cvMulSpectrums(src_DFT,temp_DFT,src_DFT,CV_DXT_MUL_CONJ); // Split Fourier in real and imaginary parts cvSplit(src_DFT,imageRe,imageIm,0,0); // Compute the magnitude of the spectrum components: Mag = sqrt(Re^2 + Im^2) cvPow( imageRe, imageMag, 2.0 ); cvPow( imageIm, imageImMag, 2.0 ); cvAdd( imageMag, imageImMag, imageMag, NULL ); cvPow( imageMag, imageMag, 0.5 ); // Normalize correlation (Divide real and imaginary components by magnitude) cvDiv(imageRe,imageMag,imageRe,1.0); cvDiv(imageIm,imageMag,imageIm,1.0); cvMerge(imageRe,imageIm,NULL,NULL,src_DFT); // inverse dft cvDFT( src_DFT, src_DFT, CV_DXT_INVERSE_SCALE, complexInput->height ); cvSplit( src_DFT, imageRe, imageIm, 0, 0 ); double minval = 0.0; double maxval = 0.0; CvPoint minloc; CvPoint maxloc; cvMinMaxLoc(imageRe,&minval,&maxval,&minloc,&maxloc,NULL); int x=maxloc.x; // log range //if (x>(imageRe->width/2)) // x = x-imageRe->width; // positive or negative values int y=maxloc.y; // angle //if (y>(imageRe->height/2)) // y = y-imageRe->height; // positive or negative values Peak pk; pk.maxval= maxval; pk.pt=cvPoint(x,y); return pk; } void phase_correlation2D( IplImage* src, IplImage *tpl, IplImage *poc ) { int i, j, k; double tmp; /* get image properties */ int width = src->width; int height = src->height; int step = src->widthStep; int fft_size = width * height; /* setup pointers to images */ uchar *src_data = ( uchar* ) src->imageData; uchar *tpl_data = ( uchar* ) tpl->imageData; double *poc_data = ( double* )poc->imageData; /* allocate FFTW input and output arrays */ fftw_complex *img1 = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height ); fftw_complex *img2 = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height ); fftw_complex *res = ( fftw_complex* )fftw_malloc( sizeof( fftw_complex ) * width * height ); /* setup FFTW plans */ fftw_plan fft_img1 = fftw_plan_dft_2d( height ,width, img1, img1, FFTW_FORWARD, FFTW_ESTIMATE ); fftw_plan fft_img2 = fftw_plan_dft_2d( height ,width, img2, img2, FFTW_FORWARD, FFTW_ESTIMATE ); fftw_plan ifft_res = fftw_plan_dft_2d( height ,width, res, res, FFTW_BACKWARD, FFTW_ESTIMATE ); /* load images' data to FFTW input */ for( i = 0, k = 0 ; i < height ; i++ ) { for( j = 0 ; j < width ; j++, k++ ) { img1[k][0] = ( double )src_data[i * step + j]; img1[k][1] = 0.0; img2[k][0] = ( double )tpl_data[i * step + j]; img2[k][1] = 0.0; } } ///* Hamming window */ //double omega = 2.0*M_PI/(fft_size-1); //double A= 0.54; //double B= 0.46; //for(i=0,k=0;i<height;i++) //{ // for(j=0;j<width;j++,k++) // { // img1[k][0]= (img1[k][0])*(A-B*cos(omega*k)); // img2[k][0]= (img2[k][0])*(A-B*cos(omega*k)); // } //} /* obtain the FFT of img1 */ fftw_execute( fft_img1 ); /* obtain the FFT of img2 */ fftw_execute( fft_img2 ); /* obtain the cross power spectrum */ for( i = 0; i < fft_size ; i++ ) { res[i][0] = ( img2[i][0] * img1[i][0] ) - ( img2[i][1] * ( -img1[i][1] ) ); res[i][1] = ( img2[i][0] * ( -img1[i][1] ) ) + ( img2[i][1] * img1[i][0] ); tmp = sqrt( pow( res[i][0], 2.0 ) + pow( res[i][1], 2.0 ) ); res[i][0] /= tmp; res[i][1] /= tmp; } /* obtain the phase correlation array */ fftw_execute(ifft_res); //normalize and copy to result image for( i = 0 ; i < fft_size ; i++ ) { poc_data[i] = res[i][0] / ( double )fft_size; } /* deallocate FFTW arrays and plans */ fftw_destroy_plan( fft_img1 ); fftw_destroy_plan( fft_img2 ); fftw_destroy_plan( ifft_res ); fftw_free( img1 ); fftw_free( img2 ); fftw_free( res ); } Peak FFTW_test(IplImage* src,IplImage* temp) { clock_t start=clock(); int t_w=temp->width; int t_h=temp->height; /* create a new image, to store phase correlation result */ IplImage* poc = cvCreateImage( cvSize(temp->width,temp->height ), IPL_DEPTH_64F, 1 ); /* get phase correlation of input images */ phase_correlation2D( src, temp, poc ); /* find the maximum value and its location */ CvPoint minloc, maxloc; double minval, maxval; cvMinMaxLoc( poc, &minval, &maxval, &minloc, &maxloc, 0 ); /* IplImage* poc_8= cvCreateImage( cvSize(temp->width, temp->height ), 8, 1 ); cvConvertScale(poc,poc_8,(double)255/(maxval-minval),(double)(-minval)*255/(maxval-minval)); cvSaveImage("poc.png",poc_8); */ cvReleaseImage( &poc ); clock_t end=clock(); int time= end-start; //fprintf( stdout, "Time= %d using clock() \n" ,time/*dt*/ ); //fprintf( stdout, "Maxval at (%d, %d) = %2.4f\n", maxloc.x, maxloc.y, maxval ); CvPoint pt; pt.x= maxloc.x; pt.y= maxloc.y; //4 случая локации точки максимальной кореляции //if(maxloc.x>=0&&maxloc.x<=t_w/2&&maxloc.y>=0&&maxloc.y<=t_h/2) //{ // pt.x= src->width-maxloc.x; // pt.y= -maxloc.y; //} //if(maxloc.x>=t_w/2&&maxloc.x<=t_w&&maxloc.y>=0&&maxloc.y<=t_h/2) //{ // pt.x= src->width-maxloc.x; // pt.y= src->height-maxloc.y; //} //if(maxloc.x>=0&&maxloc.x<=t_w/2&&maxloc.y>=t_h/2&&maxloc.y<=t_h) //{ // /*pt.x= -maxloc.x; // pt.y= -maxloc.y;*/ // pt.x= src->width-maxloc.x; //тут непонятно // pt.y= src->height-maxloc.y; //} //if(maxloc.x>=t_w/2&&maxloc.x<=t_w&&maxloc.y>=t_h/2&&maxloc.y<=t_h) //{ // pt.x= -maxloc.x; // pt.y= src->height-maxloc.y; //} Peak pk; pk.maxval= maxval; pk.pt=pt; return pk; } Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 15, 2012 понял какой M надо брать чтобы всё изображение покрыло M=w/log(sqrt((w/2)^2+(h/2)^2)) Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 15, 2012 Фильтрация изображения в частотной области проводится обнулением части спектра (кольцо, например, даст фильтр с заданной полосой пропускания): http://www.compvision.ru/forum/index.php?showtopic=110 (пост 12) Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 16, 2012 кстати я так и не понял как считать scale он зависит от координаты x в logpolar, но вот какая формула? для angle в коде была, но я не понял как она выводится. и да еще почему то когда я пытаюсь сохранить результат фурье преобразования получается что то не то, пробовал сохранять мнимую и реальную часть, может его надо как то объединять перед сохранением? а наверно всё правильно, только обычно просматривают фурье в логарифмическом масштабе. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 16, 2012 вот это я через внешнюю программу получил фурье в лог виде вот это оно же в полярных координатах This method applies a Fourier transform to images to recover translation. Then a log-polar transformation is applied to the magnitude spectrum and the rotation and scale is recovered by using phase correlation in the log-polar space. This method exploits the fact that by operating on the magnitude spectrum of an image, the translational differences are avoided since the magnitude spectrum of an image and its translated counterpart are identical; only their phase spectrum are different. Furthermore, the log-polar transformation causes rotation and scale to be manifest as translation, whereby phase correlation can be applied to recover the rotation angle and scale factor between the pair of input images. The problem here, though, is that limited scale factors can be determined because large scale factors would alter the frequency content beyond recognition. It should be noted that the maximum scale factor recovered in [6] and [7] is 2.0 and 1.8, respectively. что такое magnitude spectrum? нормированная реальная часть фурье преобразования? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 16, 2012 Корень квадратный из суммы квадратов действительной и мнимой части: mag=sqrt(Im^2+Re^2); Вот мой кусок кода по теме: #include "stdafx.h" #pragma once #include "opencv2/core/gpumat.hpp" #include "opencv2/core/opengl_interop.hpp" #include "opencv2/gpu/gpu.hpp" #include "opencv2/ml/ml.hpp" #include "opencv2/highgui/highgui.hpp" #include "opencv2/contrib/contrib.hpp" #include "opencv2/video/tracking.hpp" #include "opencv2/imgproc/imgproc.hpp" #include "opencv2/features2d/features2d.hpp" #include "opencv2/nonfree/nonfree.hpp" #include "opencv2/calib3d/calib3d.hpp" #include "opencv2/legacy/legacy.hpp" #include <string> #include <iostream> using namespace std; using namespace cv; //---------------------------------------------------------- // Функция для перестановки четвертей изображения местамии // так, чтобы ноль спектра находился в центре изображения. //---------------------------------------------------------- void Recomb(Mat &src,Mat &dst) { int cx = src.cols>>1; int cy = src.rows>>1; Mat tmp; tmp.create(src.size(),src.type()); src(Rect(0, 0, cx, cy)).copyTo(tmp(Rect(cx, cy, cx, cy))); src(Rect(cx, cy, cx, cy)).copyTo(tmp(Rect(0, 0, cx, cy))); src(Rect(cx, 0, cx, cy)).copyTo(tmp(Rect(0, cy, cx, cy))); src(Rect(0, cy, cx, cy)).copyTo(tmp(Rect(cx, 0, cx, cy))); dst=tmp; } //---------------------------------------------------------- // По заданному изображению рассчитывает // действительную и мнимую части спектра Фурье //---------------------------------------------------------- void ForwardFFT(Mat &Src, Mat *FImg) { int M = getOptimalDFTSize( Src.rows ); int N = getOptimalDFTSize( Src.cols ); Mat padded; copyMakeBorder(Src, padded, 0, M - Src.rows, 0, N - Src.cols, BORDER_CONSTANT, Scalar::all(0)); // Создаем комплексное представление изображения // planes[0] содержит само изображение, planes[1] его мнимую часть (заполнено нулями) Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)}; Mat complexImg; merge(planes, 2, complexImg); dft(complexImg, complexImg); // После преобразования результат так-же состоит из действительной и мнимой части split(complexImg, planes); // обрежем спектр, если у него нечетное количество строк или столбцов planes[0] = planes[0](Rect(0, 0, planes[0].cols & -2, planes[0].rows & -2)); planes[1] = planes[1](Rect(0, 0, planes[1].cols & -2, planes[1].rows & -2)); Recomb(planes[0],planes[0]); Recomb(planes[1],planes[1]); // Нормализуем спектр planes[0]/=float(M*N); planes[1]/=float(M*N); FImg[0]=planes[0].clone(); FImg[1]=planes[1].clone(); } //---------------------------------------------------------- // По заданным действительной и мнимой части // спектра Фурье восстанавливает изображение //---------------------------------------------------------- void InverseFFT(Mat *FImg,Mat &Dst) { Recomb(FImg[0],FImg[0]); Recomb(FImg[1],FImg[1]); Mat complexImg; merge(FImg, 2, complexImg); // Производим обратное преобразование Фурье idft(complexImg, complexImg); split(complexImg, FImg); normalize(FImg[0], Dst, 0, 1, CV_MINMAX); } //---------------------------------------------------------- // Раскладывает изображение на амплитуду и фазу спектра Фурье //---------------------------------------------------------- void ForwardFFT_Mag_Phase(Mat &src, Mat &Mag,Mat &Phase) { Mat planes[2]; ForwardFFT(src,planes); Mag.zeros(planes[0].rows,planes[0].cols,CV_32F); Phase.zeros(planes[0].rows,planes[0].cols,CV_32F); cv::cartToPolar(planes[0],planes[1],Mag,Phase); } //---------------------------------------------------------- // По заданным амплитуде и фазе // спектра Фурье восстанавливает изображение //---------------------------------------------------------- void InverseFFT_Mag_Phase(Mat &Mag,Mat &Phase, Mat &dst) { Mat planes[2]; planes[0].create(Mag.rows,Mag.cols,CV_32F); planes[1].create(Mag.rows,Mag.cols,CV_32F); cv::polarToCart(Mag,Phase,planes[0],planes[1]); InverseFFT(planes,dst); } //---------------------------------------------------------- // Точка входа //---------------------------------------------------------- int _tmain(int argc, _TCHAR* argv[]) { // Матрица изображения Mat img; // Ампльтуда спектра Mat Mag; // Фаза спектра Mat Phase; // Грузим изображение img=imread("data/lena.bmp",0); // Покажем размеры изображения cout<<img.size().width<<endl; cout<<img.size().height<<endl; // Раскладываем изображение в спектр ForwardFFT_Mag_Phase(img,Mag,Phase); //---------------------------------------------------------- // Частотный фильтр //---------------------------------------------------------- // нарисуем кольцо посередине int R=100; // Внешний радиус кольца пропускания int r=30; // Внутренний радиус кольца пропускания Mat mask; mask.create(Mag.cols,Mag.rows,CV_32F); int cx = Mag.cols>>1; int cy = Mag.rows>>1; mask=1; cv::circle(mask,cv::Point(cx,cy),R,CV_RGB(0,0,0),-1); cv::circle(mask,cv::Point(cx,cy),r,CV_RGB(1,1,1),-1); //mask=1-mask; // Раскомментировать для инверсии cv::multiply(Mag,mask,Mag); cv::multiply(Phase,mask,Phase); //---------------------------------------------------------- // Обратное преобразование //---------------------------------------------------------- InverseFFT_Mag_Phase(Mag,Phase,img); //---------------------------------------------------------- // Вывод результатов //---------------------------------------------------------- // Преобразуем к виду удобному для отображения в окне Mat LogMag; LogMag.zeros(Mag.rows,Mag.cols,CV_32F); LogMag=(Mag+1); cv::log(LogMag,LogMag); //--------------------------------------------------- imshow("Логарифм амплитуды", LogMag); imshow("Фаза", Phase); imshow("Результат фильтрации", img); //---------------------------------------------------------- // Ждем нажатия кнопки //---------------------------------------------------------- cvWaitKey(0); return 0; } [/code] Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 16, 2012 то что у вас по коду Mag это не mag=sqrt(Im^2+Re^2)? судя по картинкам отсюда это то что обычно у фурье и показывают http://www.students.uwosh.edu/~piehld88/a2dweb/a2d_ex2.htm т.е. то что у вас называется Логарифм амплитуды. хотя непонятно как действует то, что взяли логарифм может и без него надо. кстати нашел ошибку /*mask.create(Mag.cols,Mag.rows,CV_32F);*/ mask.create(Mag.rows,Mag.cols,CV_32F); Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
Smorodov 579 Жалоба Опубликовано August 16, 2012 Да, ошибочка наблюдается, спасибо. Просто не тестировал на прямоугольных картинках. Я использовал CartToPolar и PolarToCart, т.к. они сразу делают то что нужно - мнимую и действительную части преобразуют в фазу и амплитуду и обратно. http://docs.opencv.org/modules/core/doc/operations_on_arrays.html?highlight=carttopolar#cv.CartToPolar Логарифм берется только для того, чтобы была видна картинка, на спектре очень резкий спад, и при приведении к диапазону отображаемому на экране, вы увидите просто одну точку посередине. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 16, 2012 Логарифм берется только для того, чтобы была видна картинка, на спектре очень резкий спад, и при приведении к диапазону отображаемому на экране, вы увидите просто одну точку посередине. ну да получилась лишь небольшая точка, так я не понял всё таки надо переводить в лог(т.е. впринципе это никак не повлияет если мы потом будем смотреть корреляцию)? или там где всё черное тоже цифры только маленькие? Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах
mrgloom 242 Жалоба Опубликовано August 16, 2012 что то я получаю не такую картинку файл как выводится через imshow, причем и через imshow не совпадает с тем что показывает сторонняя программа. Mat t; IplImage iplimg = LogMag; CvPoint minloc, maxloc; double minval, maxval; cvMinMaxLoc(&iplimg, &minval, &maxval, &minloc, &maxloc, 0); LogMag.convertTo(t,t.type(),(double)255/(maxval-minval),(double)(-minval)*255/(maxval-minval)); imwrite("polarToCart_Mag.png",t); т.е. по форме то всё ок тут дело именно с переводом похоже, т.е. большая часть получается черной. Поделиться сообщением Ссылка на сообщение Поделиться на других сайтах