gamatchdoc.cpp

数字摄影测量相关算法编程,在VC环境下实现影像特征点提取.

// GaMatchDoc.cpp : implementation of the CGaMatchDoc class
//

#include "stdafx.h"
#include "GaMatch.h"
#include "NCMatch.h"
#include "GaMatchDoc.h"
#include "Dibapi.h"
#include "MainFrm.h"
#include "ChildFrm.h"
#include "GaMatchView.h"
#include  "math.h"
#include "windowsx.h"
#include "TwotoOneDlg.h"
#include "FangsheDlg.h"
#include "basemath.h"
#include "ProcessDlg.h"
#include "FailureDlg.h"
#include <stdlib.h>
#include <time.h>
//#include "math.h"
#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif
short SizeMatch;
short nTask;//任务
short SizeSearch;
CPoint CurPt;//预测点坐标
/////////////////////////////////////////////////////////////////////////////
// CGaMatchDoc

IMPLEMENT_DYNCREATE(CGaMatchDoc, CDocument)

BEGIN_MESSAGE_MAP(CGaMatchDoc, CDocument)
	//{{AFX_MSG_MAP(CGaMatchDoc)
	ON_COMMAND(ID_GOOD_CHOOSE_GA, OnGoodChooseGa)
	//}}AFX_MSG_MAP
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CGaMatchDoc construction/destruction

CGaMatchDoc::CGaMatchDoc()
{
	long xx,yy,xx1,xx2,yy1,xx3;	
	char *zzss,*zzss1,*zzss2;
	char a;
	CString xxs;
    SizeCollect=128;
	for(int i=0;i<6;i++)
		rightmax[i]=0.0;
	nTask=0;
    ZoomFractor=1;
	transL=transR=NULL;
    left_right=0; 
	againmatch=0;
	againding=0;
	SizeMatch=64;
    overlap=65;
	SizeSearch=256;
//	pen=new CPen(PS_SOLID,1,RGB(255,0,0));
	HANDLE 	hMem1  = GlobalAlloc(GMEM_DISCARDABLE | GMEM_MOVEABLE,SizeMatch*(SizeMatch+1)+1);
	llii = (BYTE*)GlobalLock(hMem1);
	GlobalUnlock(hMem1);
//	llii=new BYTE[SizeMatch*SizeMatch];
//    rrii=new BYTE[SizeSearch*SizeSearch+SizeSearch/2];
	HANDLE 	hMem  = GlobalAlloc(GMEM_DISCARDABLE | GMEM_MOVEABLE,(SizeSearch+SizeMatch)*(SizeSearch+1));
	rrii = (BYTE*)GlobalLock(hMem);
	GlobalUnlock(hMem);
}

CGaMatchDoc::~CGaMatchDoc()
{ 

//	delete pen;
/*		if (rrii != NULL)
	{
		 GlobalFreePtr(rrii);
	}
	if (llii != NULL)
	{
       delete llii;
	}
*/


}

BOOL CGaMatchDoc::OnNewDocument()
{
	if (!CDocument::OnNewDocument())
		return FALSE;

	// TODO: add reinitialization code here
	// (SDI documents will reuse this document)

	return TRUE;
}



/////////////////////////////////////////////////////////////////////////////
// CGaMatchDoc serialization

void CGaMatchDoc::Serialize(CArchive& ar)
{
	if (ar.IsStoring())
	{
		// TODO: add storing code here
	}
	else
	{
		// TODO: add loading code here
	}
}

/////////////////////////////////////////////////////////////////////////////
// CGaMatchDoc diagnostics

#ifdef _DEBUG
void CGaMatchDoc::AssertValid() const
{
	CDocument::AssertValid();
}

void CGaMatchDoc::Dump(CDumpContext& dc) const
{
	CDocument::Dump(dc);
}
#endif //_DEBUG

/////////////////////////////////////////////////////////////////////////////
// CGaMatchDoc commands
void  CGaMatchDoc::InitDIBData(CPalette **m_palDIB,HDIB *m_hDIB,CSize *m_sizeDoc)
{
	if (*m_palDIB != NULL)
	{
//		delete *m_palDIB;
		*m_palDIB = NULL;
	}
	if (*m_hDIB == NULL)
	{
		return;
	}
	LPSTR lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) *m_hDIB);
	if (::DIBWidth(lpDIB) > INT_MAX ||::DIBHeight(lpDIB) > INT_MAX)
	{
		::GlobalUnlock((HGLOBAL) *m_hDIB);
		::GlobalFree((HGLOBAL) *m_hDIB);
		*m_hDIB = NULL;
		AfxMessageBox("DIB_TOO_BIG");
		return;
	}
	*m_sizeDoc = CSize((int) ::DIBWidth(lpDIB), (int) ::DIBHeight(lpDIB));
	::GlobalUnlock((HGLOBAL) *m_hDIB);
	// Create copy of palette
	*m_palDIB = new CPalette;
	if (*m_palDIB == NULL)
	{
		// we must be really low on memory
		::GlobalFree((HGLOBAL) *m_hDIB);
		*m_hDIB = NULL;
		return;
	}
	if (::CreateDIBPalette(*m_hDIB, *m_palDIB) == NULL)
	{
		// DIB may not have a palette
		delete *m_palDIB;
		*m_palDIB = NULL;
		return;
	}

}
BOOL CGaMatchDoc::OpenModel(LPCTSTR LImgFile,LPCTSTR RImgFile)
{   
	/////////////////////////////////////////////////////
	//open  files
    CGaMatchApp *CRPapp=(CGaMatchApp *)AfxGetApp();
    
	if(!MyReadImage(LImgFile,&LImg_palette,&LImage,&LImg_size))
		return FALSE;
	
/*	if(!MyReadImage(LImgFile,&LFImg_palette,&LFImage,&LFImg_size))
		return FALSE;
*/	
	if(!MyReadImage(RImgFile,&RImg_palette,&RImage,&RImg_size))
		return FALSE;
/*	if(!MyReadImage(RImgFile,&RFImg_palette,&RFImage,&RFImg_size))
		return FALSE;*/
/*	SetImageInfo((BYTE*)GetDIBData(LFImage),LImg_size.cx,LImg_size.cy,
			     (BYTE*)GetDIBData(RFImage),RImg_size.cx,RImg_size.cy);*/
	UpdateAllViews(NULL,0L,NULL);	
	d_choice=0;
	UpdateAllViews(NULL,0L,NULL);
	return TRUE;
}
BOOL CGaMatchDoc::MyReadImage(LPCTSTR lpszPathName,CPalette **m_palDIB,HDIB *m_hDIB,CSize *m_sizeDoc)
{
	CFile file;
	CFileException fe;
	CString Pathname,namel,namer;
	
	char TPathName[100];
	strcpy(TPathName,lpszPathName);
	strupr(TPathName);
	
	char *ext = strrchr(TPathName, '.');
	if(ext!=NULL)
	{
    if(strcmp(ext, ".BMP") == 0)
		;
	else
        return FALSE;
	}
	else
        return FALSE;
	
	
	if (!file.Open(lpszPathName, CFile::modeRead | CFile::shareDenyWrite, &fe))
	{
		ReportSaveLoadException(lpszPathName, &fe,
			FALSE, AFX_IDP_FAILED_TO_OPEN_DOC);
		return FALSE;
	}
	
	
	// replace calls to Serialize with ReadDIBFile function
	::GlobalUnlock((HGLOBAL) *m_hDIB);
	::GlobalFree((HGLOBAL) *m_hDIB);	
	
	BeginWaitCursor();
	TRY
	{
	    *m_hDIB =::ReadDIBFile(file);
	}

	CATCH (CFileException, eLoad)
	{
		file.Abort(); // will not throw an exception
		EndWaitCursor();
		ReportSaveLoadException(lpszPathName, eLoad,
			FALSE, AFX_IDP_FAILED_TO_OPEN_DOC);
		*m_hDIB = NULL;
		file.Close();
		return FALSE;
	}
	END_CATCH

	InitDIBData(m_palDIB,m_hDIB,m_sizeDoc);
	    
	if (*m_hDIB == NULL)
	{
		// may not be DIB format
		EndWaitCursor();
		AfxMessageBox("CANNOT_LOAD_DIB");
		file.Close();
		return FALSE;
	}
	file.Close();
	EndWaitCursor();
	return TRUE;

}
HDIB CGaMatchDoc::GetLImage(int choice) const
{
	switch(choice)
	{
	  case 0:	return LImage;	break;
	  case 1:	return LFImage;	break;
	  default:  return 0;
	}
}

HDIB CGaMatchDoc::GetRImage(int choice) const
{
	switch(choice)
	{
	  case 0:	return RImage;	break;
	  case 1:	return RFImage;	break;
	  default:  return 0;
	}
}


//Get Image Color Palette Memory point
CPalette* CGaMatchDoc::GetLImgPalette(int choice) const
{
	switch(choice)
	{
	  case 0:	return LImg_palette;	break;
	  case 1:	return LFImg_palette;	break;
	  default:  return 0;
	}
}

CPalette* CGaMatchDoc::GetRImgPalette(int choice) const
{
	switch(choice)
	{
	  case 0:	return RImg_palette;	break;
	  case 1:	return RFImg_palette;	break;
	  default:  return 0;
	};
}

//Get Size of Images	
CSize CGaMatchDoc::GetLImgSize(int choice) const
{
	CSize ZeroSize(0,0);
	switch(choice)
	{
	  case 0:	return LImg_size;	break;
	  case 1:	return LFImg_size;	break;
	  default:  return ZeroSize;
	};
}

CSize CGaMatchDoc::GetRImgSize(int choice) const
{
	CSize ZeroSize(0,0);
	switch(choice)
	{
	  case 0:	return RImg_size;	break;
	  case 1:	return RFImg_size;	break;
	  default:  return ZeroSize;
	};
}
void CGaMatchDoc::GetmatchWindow(int x,int y)
{
	int cxDIB,cyDIB;
	int i,j;
	LPSTR lpDIB;

//    HDIB hDIB =GetLImage(0);
	lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) LImage);
//	cxDIB = (int) ::DIBWidth(lpDIB);         // Size of DIB - x
//	cyDIB = (int) ::DIBHeight(lpDIB);     // Size of DIB - y
	::GlobalUnlock((HGLOBAL) LImage); 
	
	if(x-(1-1.0*overlap/100)*LImg_size.cx<SizeMatch/2)
		            x=SizeMatch/2+(1-1.0*overlap/100)*LImg_size.cx;
	if(y<SizeMatch/2) y=SizeMatch/2;
	if(x>LImg_size.cx-SizeMatch/2-1) x=LImg_size.cx-SizeMatch/2-1;
	if(y>LImg_size.cy-SizeMatch/2-1) y=LImg_size.cy-SizeMatch/2-1;
	for(i=-SizeMatch/2;i<=SizeMatch/2;i++)
		for(j=-SizeMatch/2;j<=SizeMatch/2;j++)
			llii[(i+SizeMatch/2)*SizeMatch+j+SizeMatch/2]
							=lpDIB[(LImg_size.cy-y-i-1)*LImg_size.cx+x+j+1064]; 
}
void CGaMatchDoc::GetSearchWindow(int x,int y)
{
	int cxDIB,cyDIB;
	int i,j;
	LPSTR lpDIB;
	HDIB hDIB =GetRImage(0);
	lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) hDIB);
	cxDIB = (int) ::DIBWidth(lpDIB);         // Size of DIB - x
	cyDIB = (int) ::DIBHeight(lpDIB);        // Size of DIB - y
	::GlobalUnlock((HGLOBAL) hDIB);
	if(x<SizeSearch/2) x=SizeSearch/2;
	if(y<SizeSearch/2) y=SizeSearch/2;
	if(x>cxDIB-SizeSearch/2-1) x=cxDIB-SizeSearch/2-1;
	if(y>cyDIB-SizeSearch/2-1) y=cyDIB-SizeSearch/2-1;
	CurPt.x=x;
	CurPt.y=y;
	for(i=-SizeSearch/2;i<=SizeSearch/2;i++)
		for(j=-SizeSearch/2;j<=SizeSearch/2;j++)
			rrii[(i+SizeSearch/2)*SizeSearch+j+SizeSearch/2]
							=lpDIB[(cyDIB-y-i-1)*cxDIB+x+j+1064];

}
void CGaMatchDoc::GetSearchWind(int x,int y)//考虑航向60%的重叠度
{
	int cxDIB,cyDIB;
	int i,j;
	LPSTR lpDIB;
	HDIB hDIB =GetRImage(0);
	lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) hDIB);
	cxDIB = (int) ::DIBWidth(lpDIB);         // Size of DIB - x
	cyDIB = (int) ::DIBHeight(lpDIB);        // Size of DIB - y
	::GlobalUnlock((HGLOBAL) hDIB);
	x=x-(1-1.0*overlap/100)*cxDIB;
	if(x<SizeSearch/2) x=SizeSearch/2;
	if(y<SizeSearch/2) y=SizeSearch/2;
	if(x>cxDIB-SizeSearch/2-1) x=cxDIB-SizeSearch/2-1;
	if(y>cyDIB-SizeSearch/2-1) y=cyDIB-SizeSearch/2-1;
	CurPt.x=x;
	CurPt.y=y;
	for(i=-SizeSearch/2;i<=SizeSearch/2;i++)
		for(j=-SizeSearch/2;j<=SizeSearch/2;j++)
			rrii[(i+SizeSearch/2)*SizeSearch+j+SizeSearch/2]
							=lpDIB[(cyDIB-y-i-1)*cxDIB+x+j+1064];

}
void CGaMatchDoc::OnGoodChooseGa() 
{
   nTask=10;	
}
void CGaMatchDoc::GoodChooseGa(int x,int y) //择优选取良种
{
	BeginWaitCursor();
	int i,j,step,*matchii,*matchjj;
	int random1,random2,random[16];
	double adnn,maxr=0.0,maxl=0.0,adsum,*matchadnn;
	int lx,ly,rx,ry;
    clock_t start, finish;
	lx=x;ly=y;
    GetmatchWindow(lx,ly);
	rx=x;ry=y;
    GetSearchWind(rx,ry);
	int kk=0,kkk,nn,kk0=0;
	step=SizeMatch/8;//4,16
	nn=((SizeSearch-SizeMatch+1)/step+1)*((SizeSearch-SizeMatch+1)/step+1);
	matchii=new int [nn+2];
	matchjj=new int [nn+2];
	matchadnn=new double [nn+2];
	start = clock();//获得开始时间
    //初始种群的选取
	float sumgray=0.0;
	float avergray;
	float aveA=average(llii,SizeMatch,0,0,SizeMatch);
	for(i=0;i<SizeSearch-SizeMatch+1;i+=1)
	{
		for(j=0;j<SizeSearch-SizeMatch+1;j+=1)
		{
			float aveB=average(rrii,SizeSearch,i,j,SizeMatch);
			float sigmaA=sigma(llii,SizeMatch,0,0,aveA,SizeMatch);		 //方差
			float sigmaB=sigma(rrii,SizeSearch,i,j,aveB,SizeMatch);
			float C_AB=sigma_AB(rrii,i,j,aveB,SizeSearch		 // 协方差
								,llii,aveA,SizeMatch);
			float rel_k=C_AB/sqrt(sigmaA*sigmaB);		//相关系数
			rel_k=fabs(rel_k);
		    if(rel_k>=maxr)//寻找最大值
			{
			maxr=rel_k;
			matchi=i+CurPt.y-SizeSearch/2+SizeMatch/2;
			matchj=j+CurPt.x-SizeSearch/2+SizeMatch/2;
			}
		}
	}
	finish = clock();  //结束时间
	EndWaitCursor();
	if(nTask==10)
	{
	POSITION pos=GetFirstViewPosition();
	CGaMatchView* pView=(CGaMatchView* )GetNextView(pos);
	CPoint pt0=pView->GetDeviceScrollPosition();
	CDC* dc=pView->GetDC();
	CPen pen1(PS_SOLID,1,RGB(255,0,0));
	dc->SelectObject(&pen1);
	dc->MoveTo(matchj-pt0.x-5,matchi-pt0.y);
	dc->LineTo(matchj-pt0.x+5,matchi-pt0.y);
	dc->MoveTo(matchj-pt0.x,matchi-pt0.y-5);
	dc->LineTo(matchj-pt0.x,matchi-pt0.y+5);   
	static int num=1;
	char string[30];
	sprintf(string,"%d",num);
    dc->TextOut(matchj-pt0.x+5,matchi-pt0.y+5,string);
	num++;
    double duration = (double)(finish - start); // CLOCKS_PER_SEC;计算经历时间段  
	char ss[80];
	wsprintf(ss,"i=%d,j=%d,maxr=%d,time_span=%d,nn=%d",matchi,matchj,(int)(maxr*100),(int)duration,(int)nn*30);
	AfxMessageBox(ss);
	}
	}