dlgreg.cpp

《Visual C++数字图像获取 处理及实践应用》一书的源代码

	// 循环变量
	int i,j;

	// 如果特征选取标志位为TRUE，则进行特征点的选取和配准，否则退出
	if(!m_bChoseFeature){
		return;
	}

	// 待配准图象的特征选取区域，在这里选择特征点的选择区域要比图象的区
	// 域小一圈
	CRect rectChoose;
	rectChoose.bottom = m_rectResltImage.bottom - 5;
	rectChoose.top    = m_rectResltImage.top + 5;
	rectChoose.left   = m_rectResltImage.left + 5;
	rectChoose.right  = m_rectResltImage.right - 5;

	// 特征点的区域
	CRect rectFeature;

	// 标志位，表示此点是否是已经选择的特征点
	BOOL bFlag = FALSE;

	// 判断此点是否合法，并判断此点是否已经选择，如果是，则去掉此点
	if(rectChoose.PtInRect(point))
	{
		// 如果所选择的特征点是以前的特征点，则去掉此点
		for( i = 0; i<m_nChsFeatureNum; i++){

			// 选择特征点的显示区域，以便对特征点进行取舍
			rectFeature.bottom = m_pPointSampl[i].y + m_rectResltImage.top + 5;
			rectFeature.top    = m_pPointSampl[i].y + m_rectResltImage.top - 5;
			rectFeature.left   = m_pPointSampl[i].x + m_rectResltImage.left- 5;
			rectFeature.right  = m_pPointSampl[i].x + m_rectResltImage.left+ 5;

			// 判断所选择的特征点是否为原来选择的特征点
			// 如果是，则去掉此特征点
			if(rectFeature.PtInRect(point)){
				// 将后面的特征点向前移动一位，去掉所选择的此特征点
				for(j=i; j<m_nChsFeatureNum-1; j++){
					m_pPointSampl[j] = m_pPointSampl[j+1];
					m_pPointBase[j]  = m_pPointBase[j+1];						
				}

				// 将特征点的计数减一
				m_nChsFeatureNum--;
				
				// 更新显示
				Invalidate();
				
				// 设置标志位
				bFlag = TRUE;

				// 退出
				return;
			}
		}
		
		// 在判断特征点是否已经选取够了
		if(m_nChsFeatureNum == 3){
			AfxMessageBox("你已经选取了3个特征点，如果要继续选取，你可以去掉配准不正确的特征点在进行选取");
			return;
		}

		// 如果此点是需要选取的，则进行相关操作
		if(!bFlag){
			// 将此待配准特征点选取，注意特征点的坐标是以图象的左上角为原点确定的
			m_pPointSampl[m_nChsFeatureNum].x = point.x - m_rectResltImage.left;
			m_pPointSampl[m_nChsFeatureNum].y = point.y - m_rectResltImage.top;
			
			// 配准此特征点
			m_pPointBase[m_nChsFeatureNum] = FindMatchPoint(m_pDibInit, m_pDibSamp, 
				m_pPointSampl[m_nChsFeatureNum]);

			// 将特征点计数加一
			m_nChsFeatureNum++;
		}
			
		
	}	

	// 更新显示
	Invalidate();

	CDialog::OnLButtonUp(nFlags, point);
}

/*************************************************************************
 *
 * \函数名称：
 *   OnLButtonMove()
 *
 * \输入参数:
 *   无
 *
 * \返回值:
 *   无
 *
 * \说明:
 *   如果特征点选取标志位为TRUE，则该函数将待配准区域的鼠标设置为十字形状，
 *以便能更精确的定位特征点。
 *
 *************************************************************************
 */
void CDlgReg::OnMouseMove(UINT nFlags, CPoint point) 
{
	// 如果不是在特征点选取状态，则不进行相关的操作
	if(m_bChoseFeature){
		// 如果鼠标在待配准图象区域中，则更改鼠标形状
		if(m_rectResltImage.PtInRect(point))
		{
			::SetCursor(LoadCursor(NULL,IDC_CROSS));
		}
		else
		{
			::SetCursor(LoadCursor(NULL,IDC_ARROW));
		}
	}
	CDialog::OnMouseMove(nFlags, point);
}

/*************************************************************************
 *
 * \函数名称：
 *   FindMatchPoint()
 *
 * \输入参数:
 *   CDib*	pDibBase	- 基准图象指针
 *   CDib*	pDibSamp	- 待配准图象指针
 *   CPoint	pointSamp	- 待配准特征点的位置
 *
 * \返回值:
 *   CPoint			- 返回在基准图象中配准的特征点的位置
 *
 * \说明:
 *   该函数根据待待配准图象中的特征点位置在基准图象中寻找配准特征点，并将
 *配准的特征点位置返回。在配准的过程中，采取的是块配准的方法进行配准，在这里
 *选取块的大小为7*7。搜索的方法为全局搜索。
 *
 *************************************************************************
 */
CPoint CDlgReg::FindMatchPoint(CDib* pDibBase, CDib* pDibSamp, 
			       CPoint pointSamp)
{
	// 循环变量
	int i,j,m,n;

	// 临时变量
	int nX, nY;

	// 配准特征点位置
	CPoint pointBase;

	// 配准数据块的尺寸
	int nBlockLen = 7;
	int nBlockHalfLen =3;

	// 基准图象数据指针	
	unsigned char* pBase;
	pBase = (unsigned char *)pDibBase->m_lpImage;
	
	// 待配准图象数据指针
	unsigned char* pSamp;
	pSamp = (unsigned char *)pDibSamp->m_lpImage;		

	// 特征点位置的数据配准块
	unsigned char* pUnchSampBlock;
	pUnchSampBlock = new unsigned char[nBlockLen*nBlockLen];

	// 临时分配内存，用于存放配准数据块
	unsigned char* pUnchBaseBlock;
	pUnchBaseBlock = new unsigned char[nBlockLen*nBlockLen];

	// 相似度
	double dbCor;

	// 最大相似度
	double dbMaxCor = 0;

	// 基准图象的存储大小
	CSize sizeBaseImg;
	sizeBaseImg = pDibBase->GetDibSaveDim();

	// 待配准图象的存储大小
	CSize sizeSampImg;
	sizeSampImg = pDibSamp->GetDibSaveDim();

	// 从待配准图象中提取以特征点为中心的nBlockLen*nBlockLen0的数据块
	for(i=-nBlockHalfLen; i<=nBlockHalfLen; i++){
		for(j=-nBlockHalfLen; j<=nBlockHalfLen; j++){
			// 计算此点在图象中的位置
			nX = pointSamp.x + i;
			nY = sizeSampImg.cy - pointSamp.y + j +1;

			// 提取图象数据
			pUnchSampBlock[(j+nBlockHalfLen)*nBlockLen + (i+nBlockHalfLen)] = 
				pSamp[nY*sizeSampImg.cx + nX];				
		}
	}

	// 基准图象的高度和宽度
	int nBaseImgHeight, nBaseImgWidth;
	nBaseImgHeight = pDibBase->m_lpBMIH->biHeight;
	nBaseImgWidth  = pDibBase->m_lpBMIH->biWidth;

	// 在基准图象中寻找配准特征点，采取的搜索方法为全局搜索
	for(m = nBlockHalfLen; m< nBaseImgHeight-nBlockHalfLen; m++){
		for(n=nBlockHalfLen; n<nBaseImgWidth-nBlockHalfLen; n++){
			// 提取以此点为中心，大小为nBlockLen*nBlockLen的数据块
			for(i=-nBlockHalfLen; i<=nBlockHalfLen; i++){
				for(j=-nBlockHalfLen; j<=nBlockHalfLen; j++){
					// 计算此点在图象中存储的位置
					nX = n + i;
					nY = sizeBaseImg.cy - m + j - 1;

					// 提取图象数据
					pUnchBaseBlock[(j+nBlockHalfLen)*nBlockLen + (i+nBlockHalfLen)] = 
						pBase[nY*sizeBaseImg.cx + nX];				
				}
			}

			// 对这两个数据块进行配准，计算相似度
			dbCor = CalCorrelation(pUnchBaseBlock, pUnchSampBlock, nBlockLen);

			// 判断是否为最大相似度，如果是，则记录此相似度和配准特征点位置
			if(dbCor > dbMaxCor){
				dbMaxCor  = dbCor;
				pointBase.x = n;
				pointBase.y = m;
			}

		}
	}
	
	return pointBase;
}

/*************************************************************************
 *
 * \函数名称：
 *   CalCorrelation()
 *
 * \输入参数:
 *   unsigned char*	pBase		- 基准图象数据指针
 *   unsigned char*	pSamp		- 待配准图象数据指针
 *   int		nBlockLen	- 配准数据块的尺度大小
 *
 * \返回值:
 *   double			- 返回两个数据块配准的相似度
 *
 * \说明:
 *   该函数对给定的两个大小为nBlockLen*nBlockLen的数据块，计算两者之间的
 *的配准相似度。其中，去掉均值以消除亮度变换的影响。
 *
 *************************************************************************
 */
double CDlgReg::CalCorrelation(unsigned char* pBase, unsigned char* pSamp, int nBlockLen)
{
	// 临时变量
	double dbSelfBase=0,dbSelfSamp=0;

	// 相似度
	double dbCor=0;

	// 块均值
	double dbMeanBase=0,dbMeanSamp=0;

	// 计算两个块的平均值
	for(int i=0;i<nBlockLen;i++)
		for(int j=0;j<nBlockLen;j++)
		{
			dbMeanBase += pBase[j*nBlockLen + i];
			dbMeanSamp += pSamp[j*nBlockLen + i];
		}
	dbMeanBase = dbMeanBase/(nBlockLen*nBlockLen);
	dbMeanSamp = dbMeanSamp/(nBlockLen*nBlockLen);

	// 求取相似度
	for(i=0;i<nBlockLen;i++)
		for(int j=0;j<nBlockLen;j++)
		{
			dbSelfBase += (pBase[j*nBlockLen + i] - dbMeanBase)*
				(pBase[j*nBlockLen + i]  - dbMeanBase);
			dbSelfSamp += (pSamp[j*nBlockLen + i] - dbMeanSamp)*
				(pSamp[j*nBlockLen + i]  - dbMeanSamp);
			dbCor += (pBase[j*nBlockLen + i] - dbMeanBase)*
				(pSamp[j*nBlockLen + i]  - dbMeanSamp);
		}
	dbCor = dbCor / sqrt(dbSelfBase * dbSelfSamp);

	// 返回相似度
	return dbCor;
}

/*************************************************************************
 *
 * \函数名称：
 *   OnRegReg()
 *
 * \输入参数:
 *   无
 *
 * \返回值:
 *   无
 *
 * \说明:
 *   该函数对两幅图象m_pDibInit和m_pDibSamp进行配准。其中，需要先确定特征点，
 *特征点至少需要选取3个以上，这一步骤需要在“选取特征点”步骤中进行。根据选取
 *的特征点，计算得到仿射参数。最后，将配准的两幅图象拼接起来显示。
 *
 *************************************************************************
 */
void CDlgReg::OnRegReg() 
{
	// 进行配准之前，判断是否已经选取特征点
	if(!m_bChoseFeature){
		AfxMessageBox("还没有选取特征点，请先选取特征点");
		return;
	}

	// 如果选取的特征点数目不够，也不能进行处理
	if(m_nChsFeatureNum < 3){
		AfxMessageBox("请选择三个特征点，再进行图象配准");
		return;
	}

	// 仿射变换系数，仿射变换的系数为6个。此系数为从基准图象到待配准图象的变换系数
	double *pDbBs2SpAffPara;
	pDbBs2SpAffPara = new double[2*3];

	// 仿射变换系数，仿射变换的系数为6个。此系数为从待配准图象到基准图象的变换系数
	double *pDbSp2BsAffPara;
	pDbSp2BsAffPara = new double[2*3];

	// 如果已经选取了特征点，并得到了配准的特征点，则可以计算仿射变换系数了
	GetAffinePara(m_pPointBase, m_pPointSampl, pDbSp2BsAffPara);
	
	// 计算从待配准图象到基准图象的仿射变换系数
	GetAffinePara(m_pPointSampl, m_pPointBase, pDbBs2SpAffPara);

	// 计算图象经过仿射变换后的尺寸大小
	CRect rectAftAff;
	rectAftAff = GetAftAffDim(pDbSp2BsAffPara);
	
	// 根据图象的尺寸大小创建新的图象
	//CreateDIB();
	int nNewImgWidth, nNewImgHeight;
	nNewImgWidth  = rectAftAff.right - rectAftAff.left;
	nNewImgHeight = rectAftAff.bottom- rectAftAff.top;

	
	// 将基准图象放入新的图象中
	LPBYTE lpBaseImg;
	lpBaseImg = SetBaseImgAftAff(rectAftAff);

	// 将待配准图象放入新的图象中
	LPBYTE lpSampImg;
	lpSampImg = SetSampImgAftAff(pDbBs2SpAffPara, rectAftAff);

	// 将此图象用CDib类封装
	m_pDibResult = new CDib(CSize(nNewImgWidth,nNewImgHeight), 8);

	// 计算结果图象的存储大小尺寸
	CSize sizeSaveResult;
	sizeSaveResult = m_pDibResult->GetDibSaveDim();

	// 拷贝调色板
	memcpy(m_pDibResult->m_lpvColorTable, m_pDibInit->m_lpvColorTable, m_pDibResult->m_nColorTableEntries*sizeof(RGBQUAD));

	// 应用调色板
	m_pDibResult->MakePalette();
	
	// 分配内存给合并后的图象
	LPBYTE lpImgResult;
	lpImgResult = (LPBYTE)new unsigned char[sizeSaveResult.cx * sizeSaveResult.cy];

	// 对图象进行赋值
	for( int i=0; i<nNewImgWidth; i++)
		for( int j=0; j<nNewImgHeight; j++){
			int nX = i;
			int nY = sizeSaveResult.cy - j - 1;
			lpImgResult[nY*sizeSaveResult.cx + nX] = (lpBaseImg[j*nNewImgWidth + i] +lpSampImg[j*nNewImgWidth + i])/2;
		}

	// 将指针赋值给CDib类的数据
	m_pDibResult->m_lpImage = lpImgResult;

	// 显示合并后的图象
	CDlgAftReg* pDlg;
	pDlg = new CDlgAftReg(NULL, m_pDibResult);
	pDlg->DoModal();
	
	// 删除对象
	delete pDlg;

	// 释放已分配内存
	delete[]lpBaseImg;
	delete[]lpSampImg;
	delete[]pDbBs2SpAffPara;
	delete[]pDbSp2BsAffPara;
}

/*************************************************************************
 *
 * \函数名称：
 *   GetAffinePara()
 *
 * \输入参数:
 *   double  *pDbBs2SpAffPara	- 用于存放基准图象到待配准图象的仿射变换系数
 *   double  *pDbSp2BsAffPara	- 用于存放待配准图象到基准图象的仿射变换系数
 *
 * \返回值:
 *   无
 *
 * \说明:
 *   该函数根据得到的三对配准的特征点，计算仿射变换系数。得到的仿射变换系数
 *存放在两个输入参数所指向的内存中。
 *
 *************************************************************************
 */
void CDlgReg::GetAffinePara(CPoint* pPointBase, CPoint* pPointSampl, double* pDbAffPara)
{
	// pDbBMatrix中存放的是基准图象中特征点的坐标，
	// 大小为2*m_nChsFeatureNum，前m_nChsFeatureNum为X坐标
	double *pDbBMatrix;
	pDbBMatrix = new double[2*m_nChsFeatureNum];

	// pDbSMatrix中存放的是待配准图象中特征点的扩展坐标
	// 大小为3*m_nChsFeatureNum，其中前m_nChsFeatureNum为X坐标
	// 中间m_nChsFeatureNum个为Y坐标，后面m_nChsFeatureNum为1
	double *pDbSMatrix;
	pDbSMatrix = new double[3*m_nChsFeatureNum];

	// pDbSMatrixT中存放的pDbSMatrix的转置矩阵，
	// 大小为m_nChsFeatureNum*3
	double *pDbSMatrixT;
	pDbSMatrixT = new double[m_nChsFeatureNum*3];

	// pDbInvMatrix为临时变量，存放的是pDbSMatrix*pDbSMatrixT的逆
	// 大小为3*3
	double *pDbInvMatrix;
	pDbInvMatrix = new double[3*3];

	// 临时内存
	double *pDbTemp;
	pDbTemp = new double[2*3];

	// 循环变量
	int count;

	// 给矩阵赋值
	for(count = 0; count<m_nChsFeatureNum; count++)
	{
		pDbBMatrix[0*m_nChsFeatureNum + count] = pPointBase[count].x;
		pDbBMatrix[1*m_nChsFeatureNum + count] = pPointBase[count].y;
		pDbSMatrix[0*m_nChsFeatureNum + count] = pPointSampl[count].x;
		pDbSMatrix[1*m_nChsFeatureNum + count] = pPointSampl[count].y;
		pDbSMatrix[2*m_nChsFeatureNum + count] = 1;
		pDbSMatrixT[count*m_nChsFeatureNum + 0] = pPointSampl[count].x;
		pDbSMatrixT[count*m_nChsFeatureNum + 1] = pPointSampl[count].y;