testmcazerodlg.cpp

最小二乘拟和工具

// TestMCAZeroDlg.cpp : 实现文件
//

#include "stdafx.h"
#include "TestMCAZero.h"
#include "TestMCAZeroDlg.h"
#include ".\testmcazerodlg.h"
#include <math.h>
#include "mkl.h"

#ifdef _DEBUG
#define new DEBUG_NEW
#endif


// 用于应用程序“关于”菜单项的 CAboutDlg 对话框

class CAboutDlg : public CDialog
{
public:
	CAboutDlg();

// 对话框数据
	enum { IDD = IDD_ABOUTBOX };

	protected:
	virtual void DoDataExchange(CDataExchange* pDX);    // DDX/DDV 支持

// 实现
protected:
	DECLARE_MESSAGE_MAP()
};

CAboutDlg::CAboutDlg() : CDialog(CAboutDlg::IDD)
{
}

void CAboutDlg::DoDataExchange(CDataExchange* pDX)
{
	CDialog::DoDataExchange(pDX);
}

BEGIN_MESSAGE_MAP(CAboutDlg, CDialog)
END_MESSAGE_MAP()


// CTestMCAZeroDlg 对话框


void Analysis(double* A,int n,int m,double* W,double* C,double* result)
{
	ASSERT(n!=0);
	ASSERT(m!=0);

	double* U=new double[m*m];
	double* V=new double[m];
	double* T=new double[m*n];
	int lwork=2*m;
	double* work=new double[lwork];
	int* ipiv=new int[m];
	int info,i,j;

	memset(U,0,sizeof(double)*m*m);
	memset(T,0,sizeof(double)*m*n);
	memset(V,0,sizeof(double)*m);
	memset(ipiv,0,sizeof(int)*m);
	memset(result,0,sizeof(double)*m);

	for(i=0;i<m;i++)
		for(j=0;j<n;j++){
			T[i*n+j]=A[j*m+i]*W[j];			
		}
		cblas_dgemm(CblasRowMajor,CblasNoTrans,CblasNoTrans,m,m,n,1.0,T,n,A,m,0,U,m);

		dgetrf(&m,&m,U,&m,ipiv,&info);
		dgetri(&m,U,&m,ipiv,work,&lwork,&info);

		cblas_dgemv(CblasRowMajor,CblasNoTrans,m,n,1.0,T,n,C,1,0,V,1);
		cblas_dgemv(CblasRowMajor,CblasNoTrans,m,m,1.0,U,m,V,1,0,result,1);

		delete[] U;//U=NULL;
		delete[] V;//V=NULL;
		delete[] T;//T=NULL;
		delete[] work;//work=NULL;
		delete[] ipiv;//ipiv=NULL;
}

CTestMCAZeroDlg::CTestMCAZeroDlg(CWnd* pParent /*=NULL*/)
	: CDialog(CTestMCAZeroDlg::IDD, pParent)
	, m_channel(0)
	, m_vol(0)
	, bTwo(FALSE)
{
	m_hIcon = AfxGetApp()->LoadIcon(IDR_MAINFRAME);
	a0=a1=a2=0;
}

CTestMCAZeroDlg::~CTestMCAZeroDlg()
{
	for(int i=m_cv_list.GetItemCount()-1;i>=0;i--)
	{
		DeleteCV(i);
	} 
}

void CTestMCAZeroDlg::DoDataExchange(CDataExchange* pDX)
{
	CDialog::DoDataExchange(pDX);
	DDX_Text(pDX, IDC_EDIT_C, m_channel);
	DDX_Text(pDX, IDC_EDIT_V, m_vol);
	DDX_Control(pDX, IDC_LIST, m_cv_list);
	DDX_Check(pDX, IDC_CHECK1, bTwo);
}

BEGIN_MESSAGE_MAP(CTestMCAZeroDlg, CDialog)
	ON_WM_SYSCOMMAND()
	ON_WM_PAINT()
	ON_WM_QUERYDRAGICON()
	//}}AFX_MSG_MAP
//	ON_STN_CLICKED(IDC_STATIC_A1, OnStnClickedStaticA1)
ON_BN_CLICKED(IDC_BT_ADD, OnBnClickedBtAdd)
ON_BN_CLICKED(IDC_BT_CLEAR, OnBnClickedBtClear)
ON_BN_CLICKED(IDC_BT_DEL, OnBnClickedBtDel)
//ON_BN_CLICKED(IDC_CHECK1, OnBnClickedCheck1)
ON_BN_CLICKED(IDC_BT_CALC, OnBnClickedBtCalc)
END_MESSAGE_MAP()


// CTestMCAZeroDlg 消息处理程序

BOOL CTestMCAZeroDlg::OnInitDialog()
{
	CDialog::OnInitDialog();

	m_cv_list.SetExtendedStyle(m_cv_list.GetExtendedStyle()|LVS_EX_FULLROWSELECT);
	m_cv_list.InsertColumn(0, "道址", LVCFMT_LEFT,70);
	m_cv_list.InsertColumn(1, "电压", LVCFMT_LEFT,70);

	// 将\“关于...\”菜单项添加到系统菜单中。

	// IDM_ABOUTBOX 必须在系统命令范围内。
	ASSERT((IDM_ABOUTBOX & 0xFFF0) == IDM_ABOUTBOX);
	ASSERT(IDM_ABOUTBOX < 0xF000);

	CMenu* pSysMenu = GetSystemMenu(FALSE);
	if (pSysMenu != NULL)
	{
		CString strAboutMenu;
		strAboutMenu.LoadString(IDS_ABOUTBOX);
		if (!strAboutMenu.IsEmpty())
		{
			pSysMenu->AppendMenu(MF_SEPARATOR);
			pSysMenu->AppendMenu(MF_STRING, IDM_ABOUTBOX, strAboutMenu);
		}
	}

	// 设置此对话框的图标。当应用程序主窗口不是对话框时，框架将自动
	//  执行此操作
	SetIcon(m_hIcon, TRUE);			// 设置大图标
	SetIcon(m_hIcon, FALSE);		// 设置小图标

	// TODO: 在此添加额外的初始化代码
	
	return TRUE;  // 除非设置了控件的焦点，否则返回 TRUE
}

void CTestMCAZeroDlg::OnSysCommand(UINT nID, LPARAM lParam)
{
	if ((nID & 0xFFF0) == IDM_ABOUTBOX)
	{
		CAboutDlg dlgAbout;
		dlgAbout.DoModal();
	}
	else
	{
		CDialog::OnSysCommand(nID, lParam);
	}
}

// 如果向对话框添加最小化按钮，则需要下面的代码
//  来绘制该图标。对于使用文档/视图模型的 MFC 应用程序，
//  这将由框架自动完成。

void CTestMCAZeroDlg::OnPaint() 
{
	if (IsIconic())
	{
		CPaintDC dc(this); // 用于绘制的设备上下文

		SendMessage(WM_ICONERASEBKGND, reinterpret_cast<WPARAM>(dc.GetSafeHdc()), 0);

		// 使图标在工作矩形中居中
		int cxIcon = GetSystemMetrics(SM_CXICON);
		int cyIcon = GetSystemMetrics(SM_CYICON);
		CRect rect;
		GetClientRect(&rect);
		int x = (rect.Width() - cxIcon + 1) / 2;
		int y = (rect.Height() - cyIcon + 1) / 2;

		// 绘制图标
		dc.DrawIcon(x, y, m_hIcon);
	}
	else
	{
		OnBnClickedBtCalc();
		CDialog::OnPaint();
	}
}

//当用户拖动最小化窗口时系统调用此函数取得光标显示。
HCURSOR CTestMCAZeroDlg::OnQueryDragIcon()
{
	return static_cast<HCURSOR>(m_hIcon);
}

void CTestMCAZeroDlg::AddCV(ChannelVoltage* pCV)
{
	CString str;
	int index=m_cv_list.GetItemCount();
	str.Format("%.2f",pCV->m_channel);
	m_cv_list.InsertItem(index,str);//name
	str.Format("%.2f",pCV->m_vol);
	m_cv_list.SetItemText(index,1,str);//z
	m_cv_list.SetItemData(index,(DWORD)pCV);
}

void CTestMCAZeroDlg::DeleteCV(int index)
{
	ChannelVoltage* pCV=(ChannelVoltage*)m_cv_list.GetItemData(index);
	m_cv_list.DeleteItem(index);
	delete pCV;
}

void CTestMCAZeroDlg::OnBnClickedBtAdd()
{
	UpdateData();
	ChannelVoltage* pCV=new ChannelVoltage;
	pCV->m_channel=m_channel;
	pCV->m_vol=m_vol;
	AddCV(pCV);
	OnBnClickedBtCalc();
}

void CTestMCAZeroDlg::OnBnClickedBtClear()
{
	for(int i=m_cv_list.GetItemCount()-1;i>=0;i--)
	{
		DeleteCV(i);
	}
	OnBnClickedBtCalc();
}

void CTestMCAZeroDlg::OnBnClickedBtDel()
{
	POSITION pos = m_cv_list.GetFirstSelectedItemPosition();
	if(pos)
	{
		int nItem = m_cv_list.GetNextSelectedItem(pos);
		DeleteCV(nItem);
		OnBnClickedBtCalc();
	}	
}
void CTestMCAZeroDlg::OnBnClickedBtCalc()
{
	UpdateData();

	int count=m_cv_list.GetItemCount();
	if(count<2)return;
	
	int i,j;
	double* voltage=new double[count];
	double* address=new double[count];	
	ChannelVoltage* pCV;

	for(i=0;i<count;i++)
	{
		pCV=(ChannelVoltage*)m_cv_list.GetItemData(i);
		address[i]=pCV->m_channel;
		voltage[i]=pCV->m_vol;
	}

	int n=count;
	int m=bTwo?3:2;
	double* matrixA=new double[n*m];//set A;
	for(i=0;i<m;i++)
		for(j=0;j<n;j++)
			matrixA[j*m+i]=pow(voltage[j],m-i-1);
	double* matrixW=new double[n];//set matrixW every unit to 1
	for(i=0;i<n;i++)matrixW[i]=1;

	double* result=new double[m];
	memset(result,0,sizeof(double)*m);

	Analysis(matrixA,n,m,matrixW,address,result);

	if(bTwo)
	{
		a2=result[0];
		a1=result[1];
		a0=result[2];
	}
	else
	{
		a2=0;
		a1=result[0];
		a0=result[1];
	}

	CString str;
	str.Format("A2=%.7f",a2);
	SetDlgItemText(IDC_STATIC_A2,str);
	str.Format("A1=%.7f",a1);
	SetDlgItemText(IDC_STATIC_A1,str);
	str.Format("A0=%.7f",a0);
	SetDlgItemText(IDC_STATIC_A0,str);

	delete[] voltage;
	delete[] address;
	delete[] matrixA;
	delete[] matrixW;
	delete[] result;

	Reflash();
}

void CTestMCAZeroDlg::Reflash(void)
{
	UpdateData();

	CWnd* pWnd;
	CDC* pDC;
	CRect rect;

	pWnd =GetDlgItem(IDC_PICT);
	pDC=pWnd->GetDC();
	pWnd->GetClientRect(&rect);
	rect.DeflateRect(1,1);

	CDC memDC;
	memDC.CreateCompatibleDC(pDC);
	CBitmap m_bitmap;
	m_bitmap.CreateCompatibleBitmap(pDC, rect.Width(), rect.Height());
	CBitmap* m_oldBitmap = memDC.SelectObject(&m_bitmap);

	memDC.FillSolidRect(rect,RGB(255,255,255));
	int count=m_cv_list.GetItemCount();
	int i;
	ChannelVoltage* pCV;
	double minchannel=2048;
	double maxchannel=0;

	double minvol=30;
	double maxvol=0;
	for(i=0;i<count;i++)
	{
		pCV=(ChannelVoltage*)m_cv_list.GetItemData(i);
		minchannel=min(minchannel,pCV->m_channel);
		maxchannel=max(maxchannel,pCV->m_channel);
		minvol=min(minvol,pCV->m_vol);
		maxvol=max(maxvol,pCV->m_vol);
	}

	if(count>1)
	{
		double dx=(rect.Width()-10)/(maxvol-minvol);
		double dy=(rect.Height()-10)/(maxchannel-minchannel);

		pCV=(ChannelVoltage*)m_cv_list.GetItemData(0);
		double vol=pCV->m_vol;
		double channel=a2*vol*vol+a1*vol+a0;
		memDC.MoveTo(5+dx*(vol-minvol),rect.bottom-5-dy*(channel-minchannel));
		memDC.Ellipse(5+dx*(vol-minvol)-2,rect.bottom-5-dy*(pCV->m_channel-minchannel)-2,5+dx*(vol-minvol)+2,rect.bottom-5-dy*(pCV->m_channel-minchannel)+2);
		for(i=1;i<count;i++)
		{
			pCV=(ChannelVoltage*)m_cv_list.GetItemData(i);
			vol=pCV->m_vol;
			channel=a2*vol*vol+a1*vol+a0;
			memDC.LineTo(5+dx*(vol-minvol),rect.bottom-5-dy*(channel-minchannel));
			memDC.Ellipse(5+dx*(vol-minvol)-2,rect.bottom-5-dy*(pCV->m_channel-minchannel)-2,5+dx*(vol-minvol)+2,rect.bottom-5-dy*(pCV->m_channel-minchannel)+2);
		}
	}
	
	pDC->BitBlt(rect.left, rect.top, rect.Width(), rect.Height(),&memDC, rect.left, rect.top, SRCCOPY);
	memDC.SelectObject(m_oldBitmap);
	pWnd->ReleaseDC(pDC);
}