www.pudn.com > NoiseEstimate.rar > NoiseEstimateDlg.cpp
// NoiseEstimateDlg.cpp : 实现文件 // #include "stdafx.h" #include "NoiseEstimate.h" #include "NoiseEstimateDlg.h" #include "HistShow.h" #include "math.h" #include#include #include using namespace std; #ifdef _DEBUG #define new DEBUG_NEW #endif #define blocksize 4 //局域标准差法估计噪声时的分块大小 // 用于应用程序“关于”菜单项的 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() // CNoiseEstimateDlg 对话框 CNoiseEstimateDlg::CNoiseEstimateDlg(CWnd* pParent /*=NULL*/) : CDialog(CNoiseEstimateDlg::IDD, pParent) , m_delta(0) , m_NoiseType(_T("高斯")) , m_HistType(_T("原图")) { m_hIcon = AfxGetApp()->LoadIcon(IDR_MAINFRAME); } void CNoiseEstimateDlg::DoDataExchange(CDataExchange* pDX) { CDialog::DoDataExchange(pDX); DDX_Text(pDX, IDC_EDIT1, m_delta); DDX_CBString(pDX, IDC_COMBO1, m_NoiseType); DDX_CBString(pDX, IDC_COMBO2, m_HistType); } BEGIN_MESSAGE_MAP(CNoiseEstimateDlg, CDialog) ON_WM_SYSCOMMAND() ON_WM_PAINT() ON_WM_QUERYDRAGICON() //}}AFX_MSG_MAP ON_BN_CLICKED(IDC_BUTTONOPEN, &CNoiseEstimateDlg::OnBnClickedButtonopen) ON_BN_CLICKED(IDC_BUTTONSAVE, &CNoiseEstimateDlg::OnBnClickedButtonsave) ON_BN_CLICKED(IDC_BUTTONNOISE, &CNoiseEstimateDlg::OnBnClickedButtonnoise) ON_BN_CLICKED(ID_NOISE_ESTIMATE, &CNoiseEstimateDlg::OnBnClickedNoiseEstimate) ON_BN_CLICKED(IDC_HISTSHOW, &CNoiseEstimateDlg::OnBnClickedHistshow) ON_BN_CLICKED(IDC_NOISE_ESTIMATE2, &CNoiseEstimateDlg::OnBnClickedNoiseEstimate2) END_MESSAGE_MAP() // CNoiseEstimateDlg 消息处理程序 BOOL CNoiseEstimateDlg::OnInitDialog() { CDialog::OnInitDialog(); // 将“关于...”菜单项添加到系统菜单中。 // 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: 在此添加额外的初始化代码 m_picR=m_picB=m_picG=NULL; m_picO=new CFrameDisplay(this,IDC_STATIC); m_picA=new CFrameDisplay(this,IDC_STATIC1); m_sSaveFile=""; return TRUE; // 除非将焦点设置到控件,否则返回 TRUE } void CNoiseEstimateDlg::OnSysCommand(UINT nID, LPARAM lParam) { if ((nID & 0xFFF0) == IDM_ABOUTBOX) { CAboutDlg dlgAbout; dlgAbout.DoModal(); } else { CDialog::OnSysCommand(nID, lParam); } } // 如果向对话框添加最小化按钮,则需要下面的代码 // 来绘制该图标。对于使用文档/视图模型的 MFC 应用程序, // 这将由框架自动完成。 void CNoiseEstimateDlg::OnPaint() { CPaintDC dc(this); m_picO->ShowImage(); m_picA->ShowImage(); } //当用户拖动最小化窗口时系统调用此函数取得光标显示。 // HCURSOR CNoiseEstimateDlg::OnQueryDragIcon() { return static_cast (m_hIcon); } void CNoiseEstimateDlg::OnBnClickedButtonopen() { // TODO: 在此添加控件通知处理程序代码 CString sFileName; CFile fBMP; if (m_picR) free(m_picR); if (m_picG) free(m_picG); if (m_picB) free(m_picB); m_picR=m_picG=m_picB=NULL; //选择并打开一个原始bmp文件 CFileDialog FileDlg(TRUE,NULL,NULL,OFN_HIDEREADONLY,"Picture File(*.bmp)|*.bmp"); if (FileDlg.DoModal()!=IDOK) return; sFileName=FileDlg.GetPathName(); if (!fBMP.Open(sFileName,CFile::modeRead|CFile::typeBinary)) { AfxMessageBox("文件无法打开"); return; } //读入文件,返回pos为图像文件头的长度 if (!Read_BMP(&fBMP)) { AfxMessageBox("文件格式不是无压缩24位色Windows 3格式"); fBMP.Close(); return; } fBMP.Close(); m_picO->Copy_RGB24(m_picB,m_picG,m_picR,m_dWidth,m_dHeight); m_picO->ShowImage(); RedrawWindow(); } void CNoiseEstimateDlg::OnBnClickedButtonsave() { // TODO: 在此添加控件通知处理程序代码 CFileDialog dlg(NULL,"",""); dlg.DoModal();//显示对话框 m_sSaveFile = dlg.GetFileName(); if (m_sSaveFile!="") m_picA->SaveImage24(m_sSaveFile); m_sSaveFile=""; } void CNoiseEstimateDlg::OnBnClickedButtonnoise() { // TODO: 在此添加控件通知处理程序代码 CWaitCursor waitcursor; //UpdateData(TRUE); //CString str; //str.Format("%f",m_delta); //MessageBox(str); //添加噪声 int i,j; //unsigned char *B,*G,*R; UpdateData(TRUE); //ofstream fout; //fout.open( "C:\\Documents and Settings\\Administrator\\noise_50.txt" ); //if ( !fout.is_open() ) //{ // MessageBox( "Can't open file" ); // return; //} B = (unsigned char *)malloc(m_dWidth * m_dHeight); G = (unsigned char *)malloc(m_dWidth * m_dHeight); R = (unsigned char *)malloc(m_dWidth * m_dHeight); if (m_NoiseType == "高斯") { //MessageBox( "程序正在运行中,请等待。。。" ); for(i = 0;i < m_dHeight;i++) { for(j = 0;j < m_dWidth;j++) { /**(B + i * m_dWidth + j) = *(m_picB + i*m_dWidth + j ) + Discrete_Gauss_Noise(0, m_delta);*/ float result_G = (float)(*(m_picG + i * m_dWidth + j)) + Discrete_Gauss_Noise1(0, m_delta); if (result_G > 255) { *(G + i * m_dWidth + j) = 255; } else if (result_G < 0) { *(G + i * m_dWidth + j) = 0; } else { *(G + i * m_dWidth + j) = (unsigned char) result_G; } /**(G + i * m_dWidth + j) = *(m_picG + i * m_dWidth + j) + Discrete_Gauss_Noise(0, m_delta);*/ /**(R + i * m_dWidth + j) = *(m_picR + i * m_dWidth + j) + Discrete_Gauss_Noise(0, m_delta);*/ //fout << setw(5); //fout << Discrete_Gauss_Noise( 0, m_delta ); //下面三条语句用来生成一幅均匀灰度图像,灰度值由用户给出,大小和先打开的图像相同。 /**(B + i * m_dWidth + j) = (unsigned char)m_delta; *(G + i * m_dWidth + j) = (unsigned char)m_delta; *(R + i * m_dWidth + j) = (unsigned char)m_delta;*/ }//end j //fout << setw(5); //fout << endl; }//end i //fout << flush; //fout.close(); //MessageBox( "添加噪声已经结束,可以继续其他操作!" ); ////以下语句用来生成一幅灰度棋盘图像,棋盘格大小为8x8 //for (i = 0;i != m_dHeight;i += 8) //{ // for (j = 0;j != m_dWidth;j += 8) // { // //每格内像素相同,格之间均匀(100-200) // int k = 100 + rand() % 100; // //黑白相间(0或者255) // //int k = ((i + j) % 16 == 0) ? 0 : 255; // for (int s = 0;s != 8;s++) // { // for (int t = 0;t != 8;t++) // { // *(B + (i + s) * m_dWidth + j + t) = k; // *(R + (i + s) * m_dWidth + j + t) = k; // *(G + (i + s) * m_dWidth + j + t) = k; // }//end t // }//end s // }//end j //}//end i } else if(m_NoiseType == "椒盐") { for(i = 0;i < m_dHeight;i++) { for(j = 0;j < m_dWidth;j++) { *(B + i * m_dWidth + j) = *(m_picB + i * m_dWidth + j); *(G + i * m_dWidth + j) = *(m_picG + i * m_dWidth + j); *(R + i * m_dWidth + j) = *(m_picR + i * m_dWidth + j); int noise = rand() % 100; if (noise < 10) { *(B + i * m_dWidth + j) = 0; *(G + i * m_dWidth + j) = 0; *(R + i * m_dWidth + j) = 0; } else if (noise > 40) { *(B + i * m_dWidth + j) = 255; *(G + i * m_dWidth + j) = 255; *(R + i * m_dWidth + j) = 255; } } } } m_picA->Copy_RGB24(B,G,R,m_dWidth,m_dHeight); m_picA->ShowImage(); RedrawWindow(); } void CNoiseEstimateDlg::OnBnClickedNoiseEstimate() { //M2M4方法估计噪声 double kx,ky; long double M2; double N,sigma; M2 = Two_Stage_Squares(G,m_dWidth * m_dHeight); kx = Kurtosis(m_picG,m_dWidth * m_dHeight); ky = Kurtosis(G,m_dWidth * m_dHeight); if (fabs((kx - 3)) < 1) { MessageBox("原图为纯噪声图像,该算法不适用!"); return; } else if ((((ky - 3) / (kx - 3)) < 0)) { MessageBox("((ky - 3) / (kx - 3)) < 0,该算法不适用!"); return; } else if((((ky - 3) / (kx - 3)) > 1)) { MessageBox("((ky - 3) / (kx - 3)) > 1),该算法不适用!"); return; } else { N = M2 * (1 - sqrt((ky - 3) / (kx - 3))); sigma = sqrt(N); CString str; str.Format("%f",sigma ); //str.Format("%f",ky); MessageBox(str); } } void CNoiseEstimateDlg::OnBnClickedNoiseEstimate2() { // 局域标准差法估计噪声 //首先计算出各图像块的局部标准差,然后找出频数最大的那个作为噪声估计 int lsdmax,lsdhistmax; double lsdesitimate; int M = int(m_dHeight * m_dWidth / blocksize / blocksize); int N = blocksize * blocksize; double* variance = new double [m_dHeight * m_dWidth / blocksize / blocksize];//存放计算出的各个局域方差 int* lsd = new int [m_dHeight * m_dWidth / blocksize / blocksize];//存放计算出的各个局域标准差的近似值(取整) unsigned char* buf = new unsigned char [blocksize * blocksize];//临时数组,存放每一个图像块的数据 for (int i = 0; i <= m_dHeight;i += blocksize) { for (int j = 0;j <= m_dWidth;j += blocksize) { for (int s = 0; s != blocksize;s++) { for (int t = 0;t != blocksize;t++) { *(buf + s * blocksize + t) = *(G + (i + s) * blocksize + j + t); }//end t }//end s *(variance + (i / blocksize) * (m_dWidth / blocksize) + (j / blocksize) ) = var(buf);//计算出了局域方差 }//end j }//end i for (int i = 0; i != M;i++) { *(lsd + i) = int(sqrt(*(variance + i)) * 100);//计算出了局域标准差 } //求lsd的最大值 lsdmax = 0;lsdesitimate = 0; for (int i = 0;i != M;i++) { if (lsdmax < *(lsd + i)) { lsdmax = *(lsd + i); } } //lsd的直方图数组lsdhist初始化 for (int i = 0;i != lsdmax;i++) { lsdhist[i] = 0; } //求直方图lsdhist(其实不是真正的直方图,因为不是概率,而是出现的次数) for (int i = 0;i != M;i++) { lsdhist[*(lsd + i)]++; } //求直方图lsdhist的最大值和lsd的估计值 lsdhistmax = 0; for (int i = 0;i != lsdmax;i++) { if (lsdhistmax < lsdhist[i]) { lsdhistmax = lsdhist[i]; lsdesitimate = i; } } CString str; str.Format("%f",lsdesitimate / 100 * (double (sqrt(double((N - 1)/(N - 3)))))); AfxMessageBox(str); } double CNoiseEstimateDlg::var(unsigned char *buf) { //计算数组的方差值,数组大小为图像块的大小:blocksize * blocksize double ave,var; ave = 0;var = 0; for (int i = 0; i != blocksize * blocksize;i++) { ave += *(buf + i); } ave /= (blocksize * blocksize); for (int i = 0;i != blocksize * blocksize;i++) { var += pow((*(buf + i) - ave),2); } var /= (blocksize * blocksize - 1); return var; } void CNoiseEstimateDlg::OnBnClickedHistshow() { //显示直方图 UpdateData(); DWORD height,width,max; DWORD i; unsigned char *p; //求长度 height = m_dHeight; //求宽度 width = m_dWidth; //指向图像数据 if (m_HistType == "原图") { p = m_picG; } else if (m_HistType == "噪声图") { p = G; } if (p == NULL) { MessageBox("请先打开一幅图像。"); return; } //直方图数组初始化 for (i = 0;i < 256;i++) { hist[i] = 0; } //求直方图(其实不是真正的直方图,因为不是概率,而是出现的次数) for (i = 0;i < height * width;i++) { hist[*(p + i)]++; } //求直方图最大值 max = 0; for (i = 0;i < 256;i++) { if (max < hist[i]) { max = hist[i]; } } //将直方图对最大值归一化 for (i = 0;i < 256;i++) { hist[i] = hist[i] * 280 / max; } CHistShow dlg(NULL,hist); dlg.DoModal(); } double CNoiseEstimateDlg::Two_Stage_Squares(unsigned char *Y,long n) { //求一组数的二阶矩(统计量) double M2 = 0; /*n = m_dWidth * m_dHeight;*/ //double EY = 0; //for (int i = 0; i < n; i++) //EY += *(Y + i); //EY = EY / n;//均值 //M2 = 0; for (int i = 0;i < n;i++) //M2 += (*(Y + i) - EY) * (*(Y + i) - EY); M2 += (*(Y + i)) * (*(Y + i)); //由于图像的像素值均值肯定大于0,不满足论文中均值为0的条件,所以这里计算的二阶矩是原点二阶矩,不是中心二阶矩。 M2 = M2 / n; return (M2); } double CNoiseEstimateDlg::Kurtosis(unsigned char *Y,long n) { //求某一块(大小为n)的峰度系数 double k; double M2,M4; int i; //double EY = 0; //for (i = 0; i < n; i++) //EY += *(Y + i); //EY = EY / n;//均值 M2 = 0;M4 = 0; for (i = 0;i < n;i++) { /*M2 += (*(Y + i) - EY) * (*(Y + i) - EY); M4 += (*(Y + i) - EY) * (*(Y + i) - EY) * (*(Y + i) - EY) * (*(Y + i) - EY);*/ M2 += (*(Y + i)) * (*(Y + i));//同样,这里都是原点矩,不是中心矩。 M4 += (*(Y + i)) * (*(Y + i)) * (*(Y + i)) * (*(Y + i)); } M2 = M2 / n; M4 = M4 / n; k = M4 / M2; k = k / M2; return (k); } int CNoiseEstimateDlg::Discrete_Gauss_Noise(double mu, double sigma) { //产生符合所给参数的高斯分布的离散值 int k; double MAX=256; double r,y,t; t = 0; unsigned int number; for (int i = 0;i < MAX * 12;i++) { r = 1.0 * rand_s(&number) / RAND_MAX; t += r - 0.5; } //根据大数定理,上述操作产生了服从标准正态分布的值,范围是:[-MAX*12/2,MAX*12/2] y = mu + sigma * t / sqrt(MAX); if (y >= 0) t = 1; else t = -1; //取绝对值再加上0.5,使-0.5和0.5区间没有值存在(左右边分别移动0.5),避免了取整后0值加倍的情况 k = (int)((fabs(y) + 0.5) * t); return(k); } int CNoiseEstimateDlg::Discrete_Gauss_Noise1(double mu, double sigma) { //产生符合所给参数的高斯分布的离散值 int noise; //noise = floor(gen_normal(mu,sigma)); noise = Poisson(sigma); //noise = floor(Normal(mu,sigma)); //noise = floor(gen_uniform(100,200)); return noise; } BOOL CNoiseEstimateDlg::Read_BMP(CFile *fBMP) { // 读入BMP图像,将各信息赋给响应的变量或指针 unsigned char temp[4]; int i,j,m; long t; fBMP->Read(temp,2);//首先读入两个字节 if ((temp[0] != 'B')&&(temp[1] != 'M')) return(FALSE); //位图格式中文件头结构体的第一个变量是bfType,两个字节,代表文件类型,必须是0x424, //即字符串"MB" fBMP->Read(temp,4); fBMP->Read(temp,4); fBMP->Read(temp,4); fBMP->Read(temp,4); if (Bytes_to_Long(temp,4) != 40) return(FALSE); //读到了位图信息头结构体,第一个变量是biSize,WORD型,表示结构体的大小,必须是40个字节 fBMP->Read(temp,4); m_dWidth = Bytes_to_Long(temp,4); //图像宽度 fBMP->Read(temp,4); m_dHeight = Bytes_to_Long(temp,4); //图像高度 fBMP->Read(temp,2); if (Bytes_to_Long(temp,2) != 1) return(FALSE); //biPlanes必须是1 fBMP->Read(temp,2); if (Bytes_to_Long(temp,2) != 24) return(FALSE); //biBitCount是图像数据位数,这里要求是24,是真彩色图像 fBMP->Read(temp,4); if (Bytes_to_Long(temp,4) != 0) return(FALSE); //biCompression指示位图是否被压缩,这里要求非压缩 fBMP->Read(temp,4); fBMP->Read(temp,4); fBMP->Read(temp,4); fBMP->Read(temp,4); if (Bytes_to_Long(temp,4) != 0) return(FALSE); //biClrUsed指示图像用到的颜色数,为0表示用到的颜色数是2的biBitCount次方 fBMP->Read(temp,4); if (Bytes_to_Long(temp,4) != 0) return(FALSE); //biClrImpression指示图像中重要的颜色数,为0则认为所有的颜色都重要 //分配内存 m_picB = (unsigned char *)malloc(m_dWidth * m_dHeight); m_picG = (unsigned char *)malloc(m_dWidth * m_dHeight); m_picR = (unsigned char *)malloc(m_dWidth * m_dHeight); //读入图像数据 m = (4 - (m_dWidth * 3) % 4) % 4; //避免图像宽度不是3的倍数,造成读入混乱,但是这样会忽略边缘的这些零头 for (i = m_dHeight;i > 0;i--) { //BMP文件从下到上,从左到右排列,最先读到的数据是最下面一行的左边的第一个像素 //这里i是递减的,则在数组中改变了这种排列方式,即变为从上到下,从左到右了 t=(i - 1) * m_dWidth; for (j = 0;j < m_dWidth;j++) { fBMP->Read(temp,3); //注意真彩色图像的3字节图像数据依次代表B,G,R分量 *(m_picB + t + j) = temp[0]; *(m_picG + t + j) = temp[1]; *(m_picR + t + j) = temp[2]; } if (m != 0) fBMP->Read(temp,m); } return(TRUE); } long CNoiseEstimateDlg::Bytes_to_Long(unsigned char *temp, int n) { //把字节型数据-->长整型数据 long k=0; for (;n>0;n--) k=k*256+temp[n-1]; //一个字节一个字节地读,故乘于256 return(k); } void CNoiseEstimateDlg::Histogram_Equalizor(unsigned char *pointer) { //直方图均衡 int i,j,p[256]; long n,q[256]; n = m_dWidth * m_dHeight; for (i = 0;i < 256;i++) p[i] = 0; for (i = 0;i < n;i++) p[*(pointer + i)]++; for (i = 0;i < 256;i++) { q[i]=0; for (j = 0;j <= i;j++) q[i] += p[j]; } for (i = 0;i < n;i++) *(pointer + i) = (unsigned char)(q[*(pointer + i)] * 255 / n); }