neichaview.cpp

根据已知点绘拟合出附近点的高程

// NeiChaView.cpp : implementation of the CNeiChaView class
//

#include "stdafx.h"
#include "NeiCha.h"

#include "NeiChaDoc.h"
#include "NeiChaView.h"

#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView

IMPLEMENT_DYNCREATE(CNeiChaView, CView)

BEGIN_MESSAGE_MAP(CNeiChaView, CView)
	//{{AFX_MSG_MAP(CNeiChaView)
	ON_COMMAND(ID_SETCENTER, OnSetcenter)
	ON_WM_MOUSEMOVE()
	ON_COMMAND(ID_COMPUTE, OnCompute)
	//}}AFX_MSG_MAP
	// Standard printing commands
	ON_COMMAND(ID_FILE_PRINT, CView::OnFilePrint)
	ON_COMMAND(ID_FILE_PRINT_DIRECT, CView::OnFilePrint)
	ON_COMMAND(ID_FILE_PRINT_PREVIEW, CView::OnFilePrintPreview)
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView construction/destruction

CNeiChaView::CNeiChaView()
{
	// TODO: add construction code here
	m_knowPX[0]=102;
	m_knowPX[1]=109;
	m_knowPX[2]=105;
	m_knowPX[3]=103;
	m_knowPX[4]=108;
	m_knowPX[5]=105;
	m_knowPX[6]=115;
	m_knowPX[7]=118;
	m_knowPX[8]=116;
    m_knowPX[9]=113;
	m_knowPY[0]=110;
	m_knowPY[1]=113;
	m_knowPY[2]=115;
	m_knowPY[3]=103;
	m_knowPY[4]=105;
	m_knowPY[5]=108;
	m_knowPY[6]=104;
	m_knowPY[7]=108;
	m_knowPY[8]=113;
	m_knowPY[9]=118;
	m_knowPZ[0]=15;
    m_knowPZ[1]=18;
	m_knowPZ[2]=19;
	m_knowPZ[3]=17;
	m_knowPZ[4]=21;
	m_knowPZ[5]=15;
	m_knowPZ[6]=20;
	m_knowPZ[7]=15;
	m_knowPZ[8]=17;
	m_knowPZ[9]=22;
	/*CDC *pDC=GetDC();
    CPen newpen1(PS_SOLID,1,RGB(0,255,0));
	CPen *mypen=pDC->SelectObject(&newpen1);
	for(int i=0;i<10;i++)
	{
	  pDC->MoveTo(m_knowPX[i],m_knowPY[i]);
	  pDC->LineTo(m_knowPX[i],m_knowPY[i]);
	}*/
}
/********************************************************************
函数名：make_M
传入数据：以P点为中心的已知点数据X,Y,Z
得到数据：曲面拟合的系数X*X,XY,Y*Y,X,Y,-1
作用：获得数组M
*********************************************************************/
void CNeiChaView::make_M(double arrX[10],double arrY[10],double arrZ[10])
{
	for(int i=0;i<10;i++)
	{
		m_M[i][0]=arrX[i]*arrX[i];
		m_M[i][1]=arrX[i]*arrY[i];
		m_M[i][2]=arrY[i]*arrY[i];
		m_M[i][3]=arrX[i];
		m_M[i][4]=arrY[i];
		m_M[i][5]=1;
	}
}
void CNeiChaView::Get_TM(double arr[10][6])
{
	for(int i=0;i<10;i++)
		for(int j=0;j<6;j++)
			m_TM[j][i]=m_M[i][j];
}
void CNeiChaView::computeP(double arr1[10],double arr2[10])
{
	for(int i=0;i<10;i++)
		for(int j=0;j<10;j++)
		{
			if(i==j)
            {
				m_P[i][j]=(arr1[j]*arr1[j]+arr2[j]*arr2[j]);
				m_P[i][j]=1.0/m_P[i][j];
			}
			else 
			m_P[i][j]=0.0;
		}
}
void CNeiChaView::MultMTxP(double MT[6][10],double P[10][10])
{
	for(int i=0;i<6;i++)
	{   
		
		for(int j=0;j<10;j++)
		{    
			m_MTP[i][j]=0.0;
			for (int k=0;k<10;k++)
				m_MTP[i][j]+=MT[i][k]*P[k][j];
		}
	}
}
void CNeiChaView::MultMTPxM(double arr1[6][10],double arr2[10][6])
{
	for(int i=0;i<6;i++)
	{
		for(int j=0;j<6;j++)
		{
			m_MTPM[i][j]=0.0;
			for(int k=0;k<10;k++)
				m_MTPM[i][j]+=arr1[i][k]*arr2[k][j];
		}
	}
}
/**********************************************************
以下程序实现矩阵求逆,数组nis[]标记交换过程中的行数,njs[]标
记交换过程的列数,用高斯全选主元法,用d做中间变量,p做最大主元素
数组nis[]标记行位置,njs[]标记列位置
	***********************************************************/
void CNeiChaView::make_Re(int demsion,double matx[6][6])
{	
	int nis[10],njs[10];
 int i,j,k;
 double d,p=0;
 /*选出最大的元素*/
 for(k=0;k<demsion;k++)
 { 
  d=0.0;
  for(i=k;i<demsion;i++)
  {
   for(j=k;j<demsion;j++)
   {
    p=fabs(matx[i][j]);
   if(p>d)
   {d=p;nis[k]=i;njs[k]=j;}
   }
  // if (d+1.0==1.0)
           //return(0);
  }
  /*实现行列交换,使最大元素到主元位置*/
   
   if(nis[k]!=k)
   {
    for(j=0;j<demsion;j++)
    {p=matx[k][j];matx[k][j]=matx[nis[k]][j];matx[nis[k]][j]=p;}
   }
   if(njs[k]!=k)
   {
    for(i=0;i<demsion;i++)
    {p=matx[i][k];matx[i][k]=matx[i][njs[k]];matx[i][njs[k]]=p;}
   }    
   matx[k][k]=1.0/matx[k][k];
   for(j=0;j<demsion;j++)
    if(j!=k) 
     matx[k][j]=matx[k][j]*matx[k][k];
    for(i=0;i<demsion;i++)
  if(i!=k){
     for(j=0;j<demsion;j++)
      if(j!=k)
   {matx[i][j]=matx[i][j]-matx[i][k]*matx[k][j];}
  }
    
      for(i=0;i<demsion;i++)
       if(i!=k)
        matx[i][k]=(-1)*matx[i][k]*matx[k][k];
       
 }
 /*调整恢复行列次序*/
 for(k=demsion-1;k>=0;k--)
 {
  if(njs[k]!=k)
  
   for(j=0;j<demsion;j++)
   {p=matx[k][j];matx[k][j]=matx[njs[k]][j];matx[njs[k]][j]=p;}
  
   if(nis[k]!=k)
   
    for(i=0;i<demsion;i++){p=matx[i][k];matx[i][k]=matx[i][nis[k]];matx[i][nis[k]]=p;}
   
 }
	for(i=0;i<demsion;i++)
		for(j=0;j<demsion;j++)
			m_ReM[i][j]=matx[i][j];
}
void CNeiChaView::GetX()
{
	
	double arr[6][10];
	int i,j,k;
	/*input.DoModal();
	UpdateData(true);
	m_PX=input.m_Px;
	m_PY=input.m_Py;
	UpdateData(false);*/
	for(i=0;i<10;i++)
    {
		m_ppx[i]=m_knowPX[i]-m_PX;
		m_ppy[i]=m_knowPY[i]-m_PY;
		
	}
	make_M(m_ppx,m_ppy,m_ppz);
	computeP(m_ppx,m_ppy);
	Get_TM(m_M);
	MultMTxP(m_TM,m_P);
	MultMTPxM(m_MTP,m_M);
	make_Re(6,m_MTPM);
	for (i=0;i<6;i++)
	{
		for (j=0;j<10;j++)
		{
			arr[i][j]=0.0;
			for(k=0;k<6;k++)
				arr[i][j]+=m_ReM[i][k]*m_MTP[k][j];
		}
	}
	for(i=0;i<6;i++)
	{
		m_X[i]=0.0;
		
			for(k=0;k<10;k++)
				m_X[i]+=arr[i][k]*m_knowPZ[k];
	
	}
	m_PZ=m_X[5];
}
CNeiChaView::~CNeiChaView()
{
}

BOOL CNeiChaView::PreCreateWindow(CREATESTRUCT& cs)
{
	// TODO: Modify the Window class or styles here by modifying
	//  the CREATESTRUCT cs

	return CView::PreCreateWindow(cs);
}

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView drawing

void CNeiChaView::OnDraw(CDC* pDC)
{
	CNeiChaDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	// TODO: add draw code for native data here
}

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView printing

BOOL CNeiChaView::OnPreparePrinting(CPrintInfo* pInfo)
{
	// default preparation
	return DoPreparePrinting(pInfo);
}

void CNeiChaView::OnBeginPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)
{
	// TODO: add extra initialization before printing
}

void CNeiChaView::OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)
{
	// TODO: add cleanup after printing
}

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView diagnostics

#ifdef _DEBUG
void CNeiChaView::AssertValid() const
{
	CView::AssertValid();
}

void CNeiChaView::Dump(CDumpContext& dc) const
{
	CView::Dump(dc);
}

CNeiChaDoc* CNeiChaView::GetDocument() // non-debug version is inline
{
	ASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CNeiChaDoc)));
	return (CNeiChaDoc*)m_pDocument;
}
#endif //_DEBUG

/////////////////////////////////////////////////////////////////////////////
// CNeiChaView message handlers

void CNeiChaView::OnSetcenter() 
{
	// TODO: Add your command handler code here
	
	CDC *pDC=GetDC();
    CPen newpen1(PS_SOLID,1,RGB(0,255,0));
	CPen newpen2(PS_SOLID,10,RGB(255,0,0));
	CPen *mypen1=pDC->SelectObject(&newpen1);
	for(int j=0;j<10;j++)
	{
		//pDC->MoveTo(m_knowPX[j],m_knowPY[j]);
		//pDC->LineTo(m_knowPX[j],m_knowPY[j]);
		pDC->SetPixel(CPoint(m_knowPX[j],m_knowPY[j]),RGB(255,0,0));
	}
	CPen *mypen=pDC->SelectObject(&newpen2);
	input.DoModal();
    UpdateData(true);
	m_PX=input.m_Px;
	m_PY=input.m_Py;
	UpdateData(false);
	for(int i=0;i<10;i++)
    {
		m_ppx[i]=m_knowPX[i]-m_PX;
		m_ppy[i]=m_knowPY[i]-m_PY;
		//pDC->MoveTo(m_ppx[i],m_ppy[i]);
		//pDC->LineTo(m_ppx[i],m_ppy[i]);
		pDC->SetPixel(CPoint(m_ppx[i]+200,m_ppy[i]+200),RGB(0,255,0));
	}
}

void CNeiChaView::OnMouseMove(UINT nFlags, CPoint point) 
{
	// TODO: Add your message handler code here and/or call default
	CString str;
	CMainFrame *pFrame=(CMainFrame*)AfxGetApp()->m_pMainWnd;
	//获取状态栏的指针
	CStatusBar *pStatus=&pFrame->m_wndStatusBar;
	//CMoPoint mapPt(m_map.ToMapPoint(X,Y));
	if(pStatus)
	{
		str.Format("坐标x=%f",point.x);//显示x坐标
		pStatus->SetPaneText(1,str);
		str.Format("坐标y=%f",point.y);//显示y坐标
		pStatus->SetPaneText(2,str);
	}
	
	CView::OnMouseMove(nFlags, point);
}

void CNeiChaView::OnCompute() 
{
	// TODO: Add your command handler code here
	CDC *pDC=GetDC();
	CString expr,expr1,expr2,expr3,expr4,expr5,expr6;
	GetX();
	expr.Format("'P='%f",m_PZ);
	expr1.Format("'A='%f",m_X[0]);
	expr2.Format("'B='%f",m_X[1]);
	expr3.Format("'C='%f",m_X[2]);
	expr4.Format("'D='%f",m_X[3]);
	expr5.Format("'E='%f",m_X[4]);
	expr6.Format("'F='%f",m_X[5]);
	pDC->TextOut(200,200,expr);
	pDC->TextOut(200,220,expr1);
	pDC->TextOut(200,240,expr2);
	pDC->TextOut(200,260,expr3);
	pDC->TextOut(200,280,expr4);
	pDC->TextOut(200,300,expr5);
	pDC->TextOut(200,320,expr6);
}