The principle is like that...
Firstly, to obtain 2 adjacent images
extract good feature to track.
Match each features between 2 images
Get Homography matrix
Warp current image to old image by using H.
However, H is multiplied by cumulatively.
Optical flow source code is based on http://feelmare.blogspot.kr/2012/10/optical-flow-sample-source-code-using.html
This is the result of video stabilizing.
This is example source code. ---
#define MAX_COUNT 250 #define DELAY_T 3 #define PI 3.1415 void main() { ////////////////////////////////////////////////////////////////////////// //image class IplImage* image = 0; //T, T-1 image IplImage* current_Img = 0; IplImage* Old_Img = 0; //Optical Image IplImage * imgA=0; IplImage * imgB=0; //Video Load CvCapture * capture = cvCreateFileCapture("cam1.wmv"); //cvCaptureFromCAM(0); //cvCreateFileCapture("1.avi"); //Window cvNamedWindow( "Origin OpticalFlow" , WINDOW_NORMAL); //cvNamedWindow( "RealOrigin" , WINDOW_NORMAL); ////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////// //Optical Flow Variables IplImage * eig_image=0; IplImage * tmp_image=0; int corner_count = MAX_COUNT; CvPoint2D32f* cornersA = new CvPoint2D32f[ MAX_COUNT ]; CvPoint2D32f * cornersB = new CvPoint2D32f[ MAX_COUNT ]; CvSize img_sz; int win_size=20; IplImage* pyrA=0; IplImage* pyrB=0; char features_found[ MAX_COUNT ]; float feature_errors[ MAX_COUNT ]; ////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////// //Variables for time different video int one_zero=0; //int t_delay=0; double gH[9]={1,0,0, 0,1,0, 0,0,1}; CvMat gmxH = cvMat(3, 3, CV_64F, gH); //Routine Start while(1) { //capture a frame form cam if( cvGrabFrame( capture ) == 0 ) break; //image = cvRetrieveFrame( capture ); //cvShowImage("RealOrigin", image ); //Image Create if(Old_Img == 0) { image = cvRetrieveFrame( capture ); current_Img = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels); memcpy(current_Img->imageData, image->imageData, sizeof(char)*image->imageSize ); Old_Img = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels); one_zero=1; } if(one_zero == 0 ) { if(eig_image == 0) { eig_image = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels); tmp_image = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels); } //copy to image class memcpy(Old_Img->imageData, current_Img->imageData, sizeof(char)*image->imageSize ); image = cvRetrieveFrame( capture ); memcpy(current_Img->imageData, image->imageData, sizeof(char)*image->imageSize ); ////////////////////////////////////////////////////////////////////////// //Create image for Optical flow if(imgA == 0) { imgA = cvCreateImage( cvSize(image->width, image->height), IPL_DEPTH_8U, 1); imgB = cvCreateImage( cvSize(image->width, image->height), IPL_DEPTH_8U, 1); } //RGB to Gray for Optical Flow cvCvtColor(current_Img, imgA, CV_BGR2GRAY); cvCvtColor(Old_Img, imgB, CV_BGR2GRAY); //extract features cvGoodFeaturesToTrack(imgA, eig_image, tmp_image, cornersA, &corner_count, 0.01, 5.0, 0, 3, 0, 0.04); cvFindCornerSubPix(imgA, cornersA, corner_count, cvSize(win_size, win_size), cvSize(-1, -1), cvTermCriteria(CV_TERMCRIT_ITER|CV_TERMCRIT_EPS, 20, 0.03)); CvSize pyr_sz = cvSize( imgA->width+8, imgB->height/3 ); if( pyrA == 0) { pyrA = cvCreateImage( pyr_sz, IPL_DEPTH_32F, 1); pyrB = cvCreateImage( pyr_sz, IPL_DEPTH_32F, 1); } //Optical flow cvCalcOpticalFlowPyrLK( imgA, imgB, pyrA, pyrB, cornersA, cornersB, corner_count, cvSize(win_size, win_size), 5, features_found, feature_errors, cvTermCriteria( CV_TERMCRIT_ITER | CV_TERMCRIT_EPS, 20, 0.3), 0); ///////////////////////////////////////////////////////////////////////// int fCount=0; for(int i=0; i < corner_count; ++i) { if( features_found[i] == 0 || feature_errors[i] > MAX_COUNT ) continue; fCount++; ////////////////////////////////////////////////////////////////////////// //Vector Length //float fVecLength = sqrt((float)((cornersA[i].x-cornersB[i].x)*(cornersA[i].x-cornersB[i].x)+(cornersA[i].y-cornersB[i].y)*(cornersA[i].y-cornersB[i].y))); //Vector Angle //float fVecSetha = fabs( atan2((float)(cornersB[i].y-cornersA[i].y), (float)(cornersB[i].x-cornersA[i].x)) * 180/PI ); //cvLine( image, cvPoint(cornersA[i].x, cornersA[i].y), cvPoint(cornersB[i].x, cornersA[i].y), CV_RGB(0, 255, 0), 2); } printf("%d \n", fCount); int inI=0; CvPoint2D32f* pt1 = new CvPoint2D32f[ fCount ]; CvPoint2D32f * pt2 = new CvPoint2D32f[ fCount ]; for(int i=0; i < corner_count; ++i) { if( features_found[i] == 0 || feature_errors[i] > MAX_COUNT ) continue; pt1[inI] = cornersA[i]; pt2[inI] = cornersB[i]; cvLine( image, cvPoint(pt1[inI].x, pt1[inI].y), cvPoint(pt2[inI].x, pt2[inI].y), CV_RGB(0, 255, 0), 2); inI++; } //FindHomography CvMat M1, M2; double H[9]; CvMat mxH = cvMat(3, 3, CV_64F, H); M1 = cvMat(1, fCount, CV_32FC2, pt1); M2 = cvMat(1, fCount, CV_32FC2, pt2); //M2 = H*M1 , old = H*current if( !cvFindHomography(&M1, &M2, &mxH, CV_RANSAC, 2)) //if( !cvFindHomography(&M1, &M2, &mxH, CV_RANSAC, 2)) { printf("Find Homography Fail!\n"); }else{ //printf(" %lf %lf %lf \n %lf %lf %lf \n %lf %lf %lf\n", H[0], H[1], H[2], H[3], H[4], H[5], H[6], H[7], H[8] ); } delete pt1; delete pt2; //warping by H //warpAffine(warped_2,warped_3,Transform_avg,Size( reSizeMat.cols, reSizeMat.rows)); //warpPerspective(cameraFrame2, WarpImg, H, Size(WarpImg.cols, WarpImg.rows)); IplImage* WarpImg = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels);; //cvCreateImage(cvSize(T1Img->width*2, T1Img->height*2), T1Img->depth, T1Img->nChannels); cvMatMul( &gmxH, &mxH, &gmxH); // Ma*Mb -> Mc printf(" %lf %lf %lf \n %lf %lf %lf \n %lf %lf %lf\n", H[0], H[1], H[2], H[3], H[4], H[5], H[6], H[7], H[8] ); printf(" -----\n"); printf(" %lf %lf %lf \n %lf %lf %lf \n %lf %lf %lf\n\n\n", gH[0], gH[1], gH[2], gH[3], gH[4], gH[5], gH[6], gH[7], gH[8] ); //cvWarpAffine(current_Img, WarpImg, &gmxH); cvWarpPerspective(current_Img, WarpImg, &gmxH); //cvWarpPerspective(Old_Img, WarpImg, &mxH); //display cvNamedWindow("Stabilizing",WINDOW_NORMAL ); cvShowImage("Stabilizing", WarpImg); cvReleaseImage(&WarpImg); // //printf("[%d] - Sheta:%lf, Length:%lf\n",i , fVecSetha, fVecLength); //cvWaitKey(0); ////////////////////////////////////////////////////////////////////////// } cvShowImage( "Origin OpticalFlow", image); ////////////////////////////////////////////////////////////////////////// //time delay one_zero++; if( (one_zero % DELAY_T ) == 0) { one_zero=0; } //break if( cvWaitKey(10) >= 0 ) break; } //release capture point cvReleaseCapture(&capture); //close the window cvDestroyWindow( "Origin" ); cvReleaseImage(&Old_Img); ////////////////////////////////////////////////////////////////////////// cvReleaseImage(&imgA); cvReleaseImage(&imgB); cvReleaseImage(&eig_image); cvReleaseImage(&tmp_image); delete cornersA; delete cornersB; cvReleaseImage(&pyrA); cvReleaseImage(&pyrB); ////////////////////////////////////////////////////////////////////////// }---