www.gusucode.com > VC++遥感原理与数字摄影测量实习源码程序 > VC++遥感原理与数字摄影测量实习源码程序\code\Moravec\Moravec.cpp
// Moravec.cpp : Defines the entry point for the console application. // Download by http://www.NewXing.com #include "stdafx.h" #include "Moravec.h" #include <vector> #include <math.h> #ifdef _DEBUG #define new DEBUG_NEW #undef THIS_FILE static char THIS_FILE[] = __FILE__; #endif #define FLT_MAX 3.402823466e+38F /* max value */ ///////////////////////////////////////////////////////////////////////////// // The one and only application object CWinApp theApp; using namespace std; struct FEATUREPOINT { int x; int y; //坐标 float IV; //兴趣值 }; vector<FEATUREPOINT> RemoveReplicative(vector<FEATUREPOINT> v) { vector<FEATUREPOINT> ret; vector<FEATUREPOINT>::iterator iter = v.begin(); ret.clear(); ret.push_back(*iter); vector<FEATUREPOINT>::iterator iter2; int i(0); BOOL bEXIST(FALSE); for (iter = v.begin(); iter != v.end(); ++iter) { bEXIST = false; for (iter2 = ret.begin(); iter2 != ret.end(); ++iter2) { int a = iter->x; int aa = iter2->x; int b = iter->y; int bb = iter2->y; if (a == aa && b == bb) { bEXIST = TRUE; //存在 break; } } if (bEXIST == FALSE) { ret.push_back(*iter); ++i; } } cout << "i = " << i << endl; return ret; } void SaveBand(int width, int height, int byteCount, int biBitCount, LPBYTE pBits, CString SavePath) { BITMAPFILEHEADER bmfh; BITMAPINFOHEADER bmih; bmfh.bfType = 0x4d42; // 0x42 = "B" 0x4d = "M" bmfh.bfReserved1 = 0; bmfh.bfReserved2 = 0; if (biBitCount == 8) { bmfh.bfOffBits = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER) + 256*4; bmfh.bfSize = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER) + 256*4 + byteCount; } else if (biBitCount == 24) { bmfh.bfOffBits = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER); bmfh.bfSize = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER) + byteCount; } else { MessageBox(NULL, "8bit or 24bit accepted", "错误", MB_ICONERROR); return; } bmih.biBitCount = biBitCount; bmih.biWidth = width; bmih.biHeight = height; bmih.biSize = 40; bmih.biPlanes = 1; bmih.biCompression = BI_RGB; bmih.biSizeImage = bmfh.bfSize; bmih.biClrImportant = 0; bmih.biClrUsed = 0; bmih.biXPelsPerMeter = 0; bmih.biYPelsPerMeter = 0; CFile f; //保存位图 if (f.Open(SavePath, CFile::modeCreate | CFile::modeWrite)) { f.Write(&bmfh, sizeof(bmfh)); f.Write(&bmih, sizeof(bmih)); if (biBitCount == 8) { RGBQUAD rgb[256]; int i(0); for (; i < 256; ++i) { rgb[i].rgbBlue = i; rgb[i].rgbGreen = i; rgb[i].rgbRed = i; rgb[i].rgbReserved = 0; } //修改颜色表255对应的颜色 rgb[255].rgbBlue = 0; rgb[255].rgbGreen = 0; rgb[255].rgbRed = 255; f.Write(rgb, sizeof(RGBQUAD)*256); } f.Write(pBits, byteCount); f.Close(); } else { MessageBox(NULL, "保存图像失败", "错误", MB_ICONERROR); return; } } int _tmain(int argc, TCHAR* argv[], TCHAR* envp[]) { if (argc != 9) { ::MessageBox(NULL,"启动参数不正确,请从GUI.exe或命令行参数启动","",MB_ICONASTERISK); return -1; } system("Mode con: COLS=50 LINES=15"); system("title Moravec_Operator"); int nRetCode = 0; // initialize MFC and print and error on failure if (!AfxWinInit(::GetModuleHandle(NULL), NULL, ::GetCommandLine(), 0)) { // TODO: change error code to suit your needs cerr << _T("Fatal Error: MFC initialization failed") << endl; nRetCode = 1; return nRetCode; } CString srcFile,dstFile,txtFile,threshold,window1,window2,flag1,flag2; //启动参数 srcFile = argv[1]; dstFile = argv[2]; txtFile = argv[3]; threshold = argv[4]; //阈值 window1 = argv[5]; window2 = argv[6]; flag1 = argv[7]; //判断是否在创建结果文件后显示 flag2 = argv[8]; int m_threshold = atoi(threshold); int m_window1 = atoi(window1); int m_window2 = atoi(window2); int m_flag1 = atoi(flag1); int m_flag2 = atoi(flag2); if (GetFileAttributes(srcFile) == -1) { //the module does not exist CString err; err.Format("文件:%s不存在", srcFile); ::MessageBox(NULL,err, "Error", MB_ICONERROR); return -1; } FILE* pSrcFile = NULL; pSrcFile = fopen(srcFile,"rb"); if(pSrcFile == NULL) { printf("原始图像打开失败\n"); nRetCode = 1; return nRetCode; } BITMAPFILEHEADER bmfh; BITMAPINFOHEADER bmih; fread(&bmfh, sizeof(BITMAPFILEHEADER), 1, pSrcFile); fread(&bmih, sizeof(BITMAPINFOHEADER), 1, pSrcFile); if(bmih.biBitCount != 8) { printf("只支持8位位图\n"); nRetCode = 1; return nRetCode; } DWORD width = bmih.biWidth; DWORD height = bmih.biHeight; DWORD _width = (width*8 + 31)/32*4; //saved width in files DWORD pixelCount = _width*height; fseek(pSrcFile, 256*sizeof(RGBQUAD), SEEK_CUR); BYTE* pSrcBits = new BYTE[pixelCount]; fread(pSrcBits,sizeof(BYTE),pixelCount,pSrcFile); fclose(pSrcFile); //开始计算 int i(0), j(0), m(0), n(0); int k = m_window1 / 2; float v1(0), v2(0), v3(0), v4(0); float IV(0); float* pInterest = new float[width*height]; ZeroMemory(pInterest, sizeof(float)*width*height); int Count1(0); for (i = k; i < height - k; ++i) { for (j = k; j < width - k; ++j) { //计算一个点 v1 = v2 = v3 = v4 = 0.0; for (m = 0; m < m_window1-1; ++m) { v1 += (pSrcBits[(i-k+m)*_width + (j)] - pSrcBits[(i-k+1+m)*_width + (j)])* (pSrcBits[(i-k+m)*_width + (j)] - pSrcBits[(i-k+1+m)*_width + (j)]); v2 += (pSrcBits[(i-k+m)*_width + (j-k+m)] - pSrcBits[(i-k+1+m)*_width + (j-k+m+1)])* (pSrcBits[(i-k+m)*_width + (j-k+m)] - pSrcBits[(i-k+1+m)*_width + (j-k+m+1)]); v3 += (pSrcBits[(i)*_width + (j-k+m)] - pSrcBits[(i)*_width + (j-k+1+m)])* (pSrcBits[(i)*_width + (j-k+m)] - pSrcBits[(i)*_width + (j-k+1+m)]); v4 += (pSrcBits[(i-k+m)*_width + (j+k-m)] - pSrcBits[(i-k+m+1)*_width + (j+k-m-1)])* (pSrcBits[(i-k+m)*_width + (j+k-m)] - pSrcBits[(i-k+m+1)*_width + (j+k-m-1)]); } IV = min(v1, v2); IV = min(IV, v3); IV = min(IV, v4); if (IV > m_threshold) { //成为候选点 pInterest[i*width+j] = IV; ++Count1; } } } // cout << "候选点数: " << Count1 << endl; //候选点列表 /* FILE* tempFile1 = fopen("d:\\temp1.txt", "wt"); if (!tempFile1) return -1; for (i = 0; i < height; ++i) { for (j = 0; j < width; ++j) { if (pInterest[i*width+height] == 0) continue; fprintf(tempFile1, "%d\t%d\t%.20f\n", i, j, pInterest[i*width+j]); } } fclose(tempFile1);*/ //候选点列表 bool* bMatrix = new bool[width*height]; memset(bMatrix, 1, sizeof(bool)*width*height); // ZeroMemory(bMatrix, sizeof(bool)*width*height); vector<FEATUREPOINT> FeaturePoint; FEATUREPOINT tempFP; ZeroMemory(&tempFP, sizeof(FEATUREPOINT)); int Count2(0); int k2 = m_window2 / 2; float max2(0); int tempX(0), tempY(0); bool b(false); for(i = k2; i < height - k2; i += m_window2) { for(j = k2; j < width - k2; j += m_window2) { max2 = 0; b = false; for (int m = i-k2; m < i+k2+1; ++m) { for (int n = j-k2; n<j+k2+1; ++n) { if(pInterest[m*width + n]>max2) { max2 = pInterest[m*width + n]; tempY = m; tempX = n; b = true; } } } if (b) { bMatrix[tempY*width+tempX] = 1; tempFP.IV = pInterest[tempY*width+tempX]; tempFP.x = tempX; tempFP.y = tempY; FeaturePoint.push_back(tempFP); ++Count2; //尚未剔除重复点 } } } cout << "未剔除重复: " << FeaturePoint.size() << endl; vector<FEATUREPOINT> FeaturePoint2; FeaturePoint2 = RemoveReplicative(FeaturePoint); // cout << "特征点数: " << Count2 << endl; cout << "剔除重复: " << FeaturePoint2.size() << endl; //保存结果位图 vector<FEATUREPOINT>::iterator iter; for (iter = FeaturePoint2.begin(); iter != FeaturePoint2.end(); ++iter) { int i = iter->y; int j = iter->x; pSrcBits[i*_width+j] = 255; pSrcBits[i*_width+j+1] = 255; pSrcBits[i*_width+j-1] = 255; pSrcBits[(i+1)*_width+j] = 255; pSrcBits[(i-1)*_width+j] = 255; pSrcBits[i*_width+j+2] = 255; pSrcBits[i*_width+j-2] = 255; pSrcBits[(i+2)*_width+j] = 255; pSrcBits[(i-2)*_width+j] = 255; } /* for (i = 0; i < height; ++i) { for (j = 0; j < width; ++j) { if (bMatrix[i*width+j]) { pSrcBits[i*_width+j] = 255; pSrcBits[i*_width+j+1] = 255; pSrcBits[i*_width+j-1] = 255; pSrcBits[(i+1)*_width+j] = 255; pSrcBits[(i-1)*_width+j] = 255; pSrcBits[i*_width+j+2] = 255; pSrcBits[i*_width+j-2] = 255; pSrcBits[(i+2)*_width+j] = 255; pSrcBits[(i-2)*_width+j] = 255; } } }*/ SaveBand(width, height, _width*height, 8, pSrcBits, dstFile); //保存特征点列表 FILE* pTxtFile= fopen(txtFile, "wt"); if (!pTxtFile) return -1; fprintf(pTxtFile, "#\tMoravec算子\n#\t源图像: %s\n#\t结果图: %s\n", srcFile, dstFile); fprintf(pTxtFile, "#\t阈值: %d\n", m_threshold); fprintf(pTxtFile, "#\t窗口1: %d×%d\n#\t窗口2: %d×%d\n", m_window1, m_window1, m_window2, m_window2); fprintf(pTxtFile, "#\t兴趣点数: %d\n", Count1); fprintf(pTxtFile, "#\t特征点数: %d\n\n", FeaturePoint2.size()); fprintf(pTxtFile, "坐标X\t坐标Y\t兴趣值\n"); // vector<FEATUREPOINT>::iterator iter; for (iter = FeaturePoint2.begin(); iter != FeaturePoint2.end(); ++ iter) fprintf(pTxtFile, "%d\t%d\t%.1f\n", iter->x, iter->y, iter->IV); fclose(pTxtFile); //clean up delete[] pInterest; delete[] bMatrix; delete[] pSrcBits; pInterest = NULL; bMatrix = NULL; pSrcBits = NULL; if (m_flag1 == 1) { ShellExecute(NULL, "open", dstFile, NULL, NULL, SW_SHOWNORMAL); } if (m_flag2 == 1) { ShellExecute(NULL, "open", txtFile, NULL, NULL, SW_SHOWNORMAL); } return nRetCode; }