splinedemoview.cpp

三次样条插值程序演示。在SplineDemoView.cpp文件中的函数BOOL CSplineDemoView::ZSpline3()和函数double Spline3(double)是实现三次样条

// SplineDemoView.cpp : implementation of the CSplineDemoView class
//

#include "stdafx.h"
#include "SplineDemo.h"

#include "SplineDemoDoc.h"
#include "SplineDemoView.h"

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

/////////////////////////////////////////////////////////////////////////////
// CSplineDemoView

IMPLEMENT_DYNCREATE(CSplineDemoView, CView)

BEGIN_MESSAGE_MAP(CSplineDemoView, CView)
	//{{AFX_MSG_MAP(CSplineDemoView)
	ON_WM_LBUTTONDOWN()
	ON_WM_RBUTTONDOWN()
	ON_WM_LBUTTONDBLCLK()
	//}}AFX_MSG_MAP
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CSplineDemoView construction/destruction

CSplineDemoView::CSplineDemoView()
{
	// TODO: add construction code here

	N=0; NN=0;
	x[0]=0.0;
	h=NULL;
	u=NULL;
	v=NULL;
	z=NULL;
	spx=NULL;
	spy=NULL;

}

CSplineDemoView::~CSplineDemoView()
{
}

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

	return CView::PreCreateWindow(cs);
}

/////////////////////////////////////////////////////////////////////////////
// CSplineDemoView drawing

void CSplineDemoView::OnDraw(CDC* pDC)
{
	CSplineDemoDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	// TODO: add draw code for native data here
	pDC->TextOut(0,0,"3 Spline Demo: 三次样条插值方法演示");
	pDC->TextOut(0,20,"click mouse,(from left to right): LBtn:set point, RBtn:done, Double click: clear rectangle.");
	GetClientRect(&DRect);
	DRect.top+=40;
	DRect.bottom-=20;
	DRect.left+=20;
	DRect.right-=20;
	pDC->Rectangle(DRect);

	CPoint point;
	if(N>0)
	{
		for(int i=1;i<=N;i++)
		{
			point.x=DRect.left+(int)(DRect.Width()*x[i]/100.0);
			point.y=DRect.bottom-(int)(DRect.Height()*y[i]/100.0);
			DrawPoint(point.x,point.y,4,pDC);
		}
		if(spx!=NULL&&spy!=NULL)
		{
			point.x=DRect.left+(int)(DRect.Width()*spx[0]/100.0);
			point.y=DRect.bottom-(int)(DRect.Height()*spy[0]/100.0);
			pDC->MoveTo(point);
			for(int j=1;j<NN;j++)
			{
				point.x=DRect.left+(int)(DRect.Width()*spx[j]/100.0);
				point.y=DRect.bottom-(int)(DRect.Height()*spy[j]/100.0);
				pDC->LineTo(point);
			}
		}
	}	
}

/////////////////////////////////////////////////////////////////////////////
// CSplineDemoView diagnostics

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

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

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

/////////////////////////////////////////////////////////////////////////////
// CSplineDemoView message handlers

void CSplineDemoView::DrawPoint(int x,int y,int size,CDC* pdc)
{
	if(size<1)size=1;
	if(size>50)size=50;
	CRect rect;
	rect.top=y-size/2;
	rect.bottom=y+size/2;
	rect.left=x-size/2;
	rect.right=x+size/2;
	CBrush brush;
	brush.CreateSolidBrush(RGB(0,0,0));
	pdc->FillRect(rect,&brush);

}

void CSplineDemoView::OnLButtonDown(UINT nFlags, CPoint point) 
{//set x and y value ,N<100,N=0,...,99.  N is the sum of datas
	double x1;
	if(point.x>DRect.left&&point.x<DRect.right&&point.y>DRect.top&&point.y<DRect.bottom&&N<99)
	{
		x1=(point.x-DRect.left)*100.0/DRect.Width();//x=0.0--100.0
		if(x1>x[N])
		{	N++;			
			x[N]=x1;
			y[N]=(DRect.bottom-point.y)*100.0/DRect.Height();//y=0.0--100.0
		}
	}		
	CView::OnLButtonDown(nFlags, point);
	this->Invalidate();
}

void CSplineDemoView::OnRButtonDown(UINT nFlags, CPoint point) 
{
	CView::OnRButtonDown(nFlags, point);
	//do spline ,if N==0,no data,if N<3,no enough data
	if(ZSpline3()==TRUE)
	{
		NN=N*10+1;
		spx=(double*)new double[NN];
		spy=(double*)new double[NN];
		int i;
		for(i=0;i<NN;i++)
		{
			spx[i]=(x[N]-x[1])*i/(N*10)+x[1];
			spy[i]=Spline3(spx[i]);
		}
	}

	this->Invalidate();
}

void CSplineDemoView::OnLButtonDblClk(UINT nFlags, CPoint point) 
{

	N=0;
	spx=NULL;
	spy=NULL;
	
	this->Invalidate();
	CView::OnLButtonDblClk(nFlags, point);
}

BOOL CSplineDemoView::ZSpline3()
{
	if(N<=3)
	{
		MessageBox("The sum of data is less than 4");
		return FALSE;
	}
	double *b=NULL;
	h=(double*)new double[N+1];//下标从1 到N
	u=(double*)new double[N+1];
	v=(double*)new double[N+1];
	z=(double*)new double[N+1];
	b=(double*)new double[N+1];
	if(!(h&&u&&v&&z&&b))
	{
		MessageBox("内存分配失败");
		return FALSE;
	}
	int i;
	for(i=1;i<N;i++)
	{
		h[i]=x[i+1]-x[i];
		b[i]=(y[i+1]-y[i])/h[i];
	}
	u[2]=2*(h[1]+h[2]);
	v[2]=6*(b[2]-b[1]);
	for(i=3;i<N;i++)
	{
		u[i]=2*(h[i-1]+h[i])-h[i-1]*h[i-1]/u[i-1];
		v[i]=6*(b[i]-b[i-1])-h[i-1]*v[i-1]/u[i-1];
	}
	z[N]=0;
	z[1]=0;
	for(i=N-1;i>=2;i--)
		z[i]=(v[i]-h[i]*z[i+1])/u[i];


	return TRUE;
}

double CSplineDemoView::Spline3(double xx)
{
	int i;
	double yy,x_ti,ti1_x;
	if(N>3&&h!=NULL&&z!=NULL)
	{
		i=1;
		if(xx<x[1])return 0.0;
		while(xx>=x[i+1]&&i+1<N)
			i++;
		if(xx>x[N])return 0.0;
		x_ti=xx-x[i];
		ti1_x=x[i+1]-xx;
		yy=z[i+1]*x_ti*x_ti*x_ti/(6*h[i]) + z[i]*ti1_x*ti1_x*ti1_x/(6*h[i])
			+ (y[i+1]/h[i]-h[i]*z[i+1]/6.0)*x_ti + (y[i]/h[i]-h[i]*z[i]/6.0)*ti1_x;	
	}
	else
		yy=0.0;
	return yy;
}