数字图像处理领域的二十四个典型算法

转自:http://blog.csdn.net/v_JULY_v/article/details/6210124

作者:July   二零一一年二月二十六日。
参考:百度百科、维基百科、vc数字图像处理。
--------------------------------------------------
数字图像处理领域的二十四个典型算法及vc实现、第一章
一、256色转灰度图
二、Walsh变换
三、二值化变换
四、阈值变换
五、傅立叶变换
六、离散余弦变换

数字图像处理领域的二十四个典型算法及vc实现、第二章
七、高斯平滑
八、图像平移
九、图像缩放
十、图像旋转
数字图像处理领域的二十四个典型算法及vc实现、第三章

 

      像处理,是对图像进行分析、加工、和处理,使其满足视觉、心理以及其他要求的技术。图像处理是信号处理在图像域上的一个应用。目前大多数的图像是以数字形式存储,因而图像处理很多情况下指数字图像处理。

      本文接下来,简单粗略介绍下数字图像处理领域中的24个经典算法,然后全部算法用vc实现。由于篇幅所限,只给出某一算法的主体代码。

      ok,请细看。

一、256色转灰度图
    算法介绍(
百度百科)
    什么叫灰度图?任何颜色都有红、绿、蓝三原色组成,假如原来某点的颜色为RGB(R,G,B),那么,我们可以通过下面几种方法,将其转换为灰度:   
   1.浮点算法:Gray=R*0.3+G*0.59+B*0.11   
   2.整数方法:Gray=(R*30+G*59+B*11)/100   
  3.移位方法:Gray =(R*28+G*151+B*77)>>8;   
   4.平均值法:Gray=(R+G+B)/3;   
   5.仅取绿色:Gray=G;   
    通过上述任一种方法求得Gray后,将原来的RGB(R,G,B)中的R,G,B统一用Gray替换,形成新的颜色RGB(Gray,Gray,Gray),用它替换原来的RGB(R,G,B)就是灰度图了。

灰度分为256阶。所以,用灰度表示的图像称作灰度图。

    程序实现:
    ok,知道了什么叫灰度图,下面,咱们就来实现此256色灰度图。
这个Convert256toGray(),即是将256色位图转化为灰度图:

void Convert256toGray(HDIB hDIB)
{
 LPSTR lpDIB; 
 // 由DIB句柄得到DIB指针并锁定DIB
 lpDIB = (LPSTR) ::GlobalLock((HGLOBAL)hDIB); 
 // 指向DIB象素数据区的指针
 LPSTR   lpDIBBits;  
 // 指向DIB象素的指针
 BYTE * lpSrc; 
 // 图像宽度
 LONG lWidth; 
 // 图像高度
 LONG   lHeight;  
 // 图像每行的字节数
 LONG lLineBytes; 
 // 指向BITMAPINFO结构的指针(Win3.0)
 LPBITMAPINFO lpbmi; 
 // 指向BITMAPCOREINFO结构的指针
 LPBITMAPCOREINFO lpbmc;
 // 获取指向BITMAPINFO结构的指针(Win3.0)
 lpbmi = (LPBITMAPINFO)lpDIB;  
 // 获取指向BITMAPCOREINFO结构的指针
 lpbmc = (LPBITMAPCOREINFO)lpDIB; 
 // 灰度映射表
 BYTE bMap[256];
 
 // 计算灰度映射表(保存各个颜色的灰度值),并更新DIB调色板
 int i,j;
 for (i = 0; i < 256; i ++)
 {
  // 计算该颜色对应的灰度值
  bMap[i] = (BYTE)(0.299 * lpbmi->bmiColors[i].rgbRed +
   0.587 * lpbmi->bmiColors[i].rgbGreen +
   0.114 * lpbmi->bmiColors[i].rgbBlue + 0.5);   
  // 更新DIB调色板红色分量
  lpbmi->bmiColors[i].rgbRed = i; 
  // 更新DIB调色板绿色分量
  lpbmi->bmiColors[i].rgbGreen = i; 
   // 更新DIB调色板蓝色分量
  lpbmi->bmiColors[i].rgbBlue = i;
  // 更新DIB调色板保留位
  lpbmi->bmiColors[i].rgbReserved = 0;
 }
 // 找到DIB图像象素起始位置
 lpDIBBits = ::FindDIBBits(lpDIB);
 // 获取图像宽度
 lWidth = ::DIBWidth(lpDIB); 
 // 获取图像高度
 lHeight = ::DIBHeight(lpDIB); 
 // 计算图像每行的字节数
 lLineBytes = WIDTHBYTES(lWidth * 8); 
 // 更换每个象素的颜色索引(即按照灰度映射表换成灰度值)

 //逐行扫描
 for(i = 0; i < lHeight; i++)
 {
    //逐列扫描
  for(j = 0; j < lWidth; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   // 变换
   *lpSrc = bMap[*lpSrc];
  }
 } 
 //解除锁定
 ::GlobalUnlock ((HGLOBAL)hDIB);

变换效果(以下若无特别说明,图示的右边部分都是为某一算法变换之后的效果):

 
二、Walsh变换
    算法介绍:
    有关Walsh变换的深入介绍,请看此论文:http://www.informatics.org.cn/doc/ucit200510/ucit20051005.pdf

    程序实现:

函数名称:WALSH()
参数:
double * f    - 指向时域值的指针
double * F    - 指向频域值的指针
r      -2的幂数
返回值:无。
说明:该函数用来实现快速沃尔什-哈达玛变换。
VOID WINAPI WALSH(double *f, double *F, int r)
{
 // 沃尔什-哈达玛变换点数
 LONG count;
 // 循环变量
 int  i,j,k;
 // 中间变量
 int  bfsize,p;
 double *X1,*X2,*X;
 // 计算快速沃尔什变换点数
 count = 1 << r;
 // 分配运算所需的数组
 X1 = new double[count];
 X2 = new double[count];
 // 将时域点写入数组X1
 memcpy(X1, f, sizeof(double) * count);
 
 // 蝶形运算
 for(k = 0; k < r; k++)
 {
  for(j = 0; j < 1<<k; j++)
  {
   bfsize = 1 << (r-k);
   for(i = 0; i < bfsize / 2; i++)
   {
    p = j * bfsize;
    X2[i + p] = X1[i + p] + X1[i + p + bfsize / 2];
    X2[i + p + bfsize / 2] = X1[i + p] - X1[i + p + bfsize / 2];
   }
  }
  // 互换X1和X2  
  X = X1;
  X1 = X2;
  X2 = X;
 }

 // 调整系数
 for(j = 0; j < count; j++)
 {
  p = 0;
  for(i = 0; i < r; i++)
  {
   if (j & (1<<i))
   {
    p += 1 << (r-i-1);
   }
  }

  F[j] = X1[p] / count;
 }
 
 // 释放内存
 delete X1;
 delete X2;
}

函数名称:DIBWalsh1()
参数:
LPSTR lpDIBBits    - 指向源DIB图像指针
LONG  lWidth       - 源图像宽度(象素数)
LONG  lHeight      - 源图像高度(象素数)
返回值:BOOL               - 成功返回TRUE,否则返回FALSE。
说明:该函数用来对图像进行沃尔什-哈达玛变换。于上面不同的是,此处是将二维
矩阵转换成一个列向量,然后对该列向量进行一次一维沃尔什-哈达玛变换。

BOOL WINAPI DIBWalsh1(LPSTR lpDIBBits, LONG lWidth, LONG lHeight)
{
 // 指向源图像的指针
 unsigned char* lpSrc;
 // 循环变量
 LONG i;
 LONG j;
 // 进行付立叶变换的宽度和高度(2的整数次方)
 LONG w;
 LONG h;
 // 中间变量
 double dTemp;
 int  wp;
 int  hp;
  // 图像每行的字节数
 LONG lLineBytes;
 // 计算图像每行的字节数
 lLineBytes = WIDTHBYTES(lWidth * 8);
 
 // 赋初值
 w = 1;
 h = 1;
 wp = 0;
 hp = 0;
 
 // 计算进行离散余弦变换的宽度和高度(2的整数次方)
 while(w * 2 <= lWidth)
 {
  w *= 2;
  wp++;
 }
 
 while(h * 2 <= lHeight)
 {
  h *= 2;
  hp++;
 }
 
 // 分配内存
 double *f = new double[w * h];
 double *F = new double[w * h];
 
 // 列
 for(i = 0; i < w; i++)
 {
  // 行
  for(j = 0; j < h; j++)
  {
   // 指向DIB第j行,第i个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - j) + i;
   
   // 给时域赋值
   f[j + i * w] = *(lpSrc);
  }
 }
 
 // 调用快速沃尔什-哈达玛变换
 WALSH(f, F, wp + hp);
 // 列
 for(i = 0; i < w; i++)
 {
  // 行
  for(j = 0; j < h; j++)
  {
   // 计算频谱
   dTemp = fabs(F[i * w + j] * 1000);
   
   // 判断是否超过255
   if (dTemp > 255)
   {
    // 对于超过的,直接设置为255
    dTemp = 255;
   }
   // 指向DIB第j行,第i个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - j) + i;
   
   // 更新源图像
   * (lpSrc) = (BYTE)(dTemp);
  }
 }
 
 //释放内存
 delete f;
 delete F;

 // 返回
 return TRUE;
}

变换效果:


三、二值化变换
    算法描述:
    二值化是图像分割的一种方法。在二值化图象的时候把大于某个临界灰度值的像素灰度设为灰度極大值,把小于这个值的像素灰度设为灰度極小值,从而实现二值化。
    根据阈值选取的不同,二值化的算法分为固定阈值和自适应阈值。 比较常用的二值化方法则有:双峰法、P参数法、迭代法和OTSU法等。

    程序实现:

void CMyDIPView::OnDraw(CDC* pDC)
{   
 CMyDIPDoc* pDoc = GetDocument();
 ASSERT_VALID(pDoc);
 if(pDoc->m_hDIB == NULL)
  return ;
 // TODO: add draw code for native data here
 int i,j;
    unsigned char *lpSrc;
 LPSTR lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) pDoc->m_hDIB);
 int cxDIB = (int) ::DIBWidth(lpDIB);         // Size of DIB - x
 int cyDIB = (int) ::DIBHeight(lpDIB);        // Size of DIB - y
 LPSTR lpDIBBits=::FindDIBBits (lpDIB);
 // 计算图像每行的字节数
 long lLineBytes = WIDTHBYTES(cxDIB * 8);
 // 每行
 for(i = 0; i < cyDIB; i++)
 {
  // 每列
  for(j = 0; j < cxDIB; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (cyDIB - 1 - i) + j;
   // 计算新的灰度值
   //*(lpSrc) = BYTE(255-*lpSrc);
  }
 }
 ::GlobalUnlock((HGLOBAL) pDoc->m_hDIB);
 CRect rect(0,0,cxDIB,cyDIB), rcDIB(0,0,cxDIB,cyDIB);
 ::PaintDIB(pDC->m_hDC, &rect, pDoc->m_hDIB, &rcDIB, pDoc->m_palDIB);
}

void CMyDIPView::OnMenuitem32778() 
{
 // TODO: Add your command handler code here
 int i,j;
    unsigned char *lpSrc;
 CMyDIPDoc* pDoc = GetDocument();
 ASSERT_VALID(pDoc);
 if(pDoc->m_hDIB == NULL)
  return ;
 LPSTR lpDIB = (LPSTR) ::GlobalLock((HGLOBAL) pDoc->m_hDIB);
 LPSTR lpDIBBits=::FindDIBBits (lpDIB);
 int cxDIB = (int) ::DIBWidth(lpDIB);         // Size of DIB - x
 int cyDIB = (int) ::DIBHeight(lpDIB);        // Size of DIB - y
 long lLineBytes = WIDTHBYTES(cxDIB * 8);     // 计算图像每行的字节数
 const float c1=150,c2=2.5;
 // 每行
 for(i = 0; i < cyDIB; i++)
 {
  // 每列
  for(j = 0; j < cxDIB; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (cyDIB - 1 - i) + j;
   
   // 计算新的灰度值
   if(*lpSrc<122) *lpSrc=BYTE(0);
   else *lpSrc = BYTE(255);
  }
 }
 ::GlobalUnlock((HGLOBAL) pDoc->m_hDIB);
    Invalidate(TRUE);
}

变换效果:


四、阈值变换
 算法描述:
  输入图像像元密度值(灰度、亮度值)按对数函数关系变换为输出图像。 

 程序实现:

//参数说明:
//LPSTR lpDIBBits:指向源DIB图像指针
//LONG  lWidth:源图像宽度(象素数)
//LONG  lHeight:源图像高度(象素数)
//BYTE  bThre:阈值
//程序说明:
//该函数用来对图像进行阈值变换。对于灰度值小于阈值的象素直接设置
灰度值为0;灰度值大于阈值的象素直接设置为255。
BOOL WINAPI ThresholdTrans(LPSTR lpDIBBits, LONG lWidth, LONG lHeight, BYTE bThre)
{
 // 指向源图像的指针
 unsigned char* lpSrc;
 // 循环变量
 LONG i;
 LONG j;
 // 图像每行的字节数
 LONG lLineBytes;
 // 计算图像每行的字节数
 lLineBytes = WIDTHBYTES(lWidth * 8);
 // 每行
 for(i = 0; i < lHeight; i++)
 {
  // 每列
  for(j = 0; j < lWidth; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   // 判断是否小于阈值
   if ((*lpSrc) < bThre)
   {
    // 直接赋值为0
    *lpSrc = 0;
   }
   else
   {
    // 直接赋值为255
    *lpSrc = 255;
   }
  }
 }
 // 返回
 return TRUE;
}


五、傅立叶变换
    算法描述:
    关于此傅里叶变换算法的具体介绍,请参考本BLOG文章:十、从头到尾彻底理解傅里叶变换算法、上

    程序实现:

函数名称:FFT()
参数:
complex<double> * TD - 指向时域数组的指针
complex<double> * FD - 指向频域数组的指针
r      -2的幂数,即迭代次数
返回值:无。
说明:该函数用来实现快速付立叶变换。
VOID WINAPI FFT(complex<double> * TD, complex<double> * FD, int r)
{
 // 付立叶变换点数
 LONG count;
 // 循环变量
 int  i,j,k;
 // 中间变量
 int  bfsize,p;
  // 角度
 double angle;
  complex<double> *W,*X1,*X2,*X;
 // 计算付立叶变换点数
 count = 1 << r;
 
 // 分配运算所需存储器
 W  = new complex<double>[count / 2];
 X1 = new complex<double>[count];
 X2 = new complex<double>[count];
 
 // 计算加权系数
 for(i = 0; i < count / 2; i++)
 {
  angle = -i * PI * 2 / count;
  W[i] = complex<double> (cos(angle), sin(angle));
 }
 
 // 将时域点写入X1
 memcpy(X1, TD, sizeof(complex<double>) * count);
 
 // 采用蝶形算法进行快速付立叶变换
 for(k = 0; k < r; k++)
 {
  for(j = 0; j < 1 << k; j++)
  {
   bfsize = 1 << (r-k);
   for(i = 0; i < bfsize / 2; i++)
   {
    p = j * bfsize;
    X2[i + p] = X1[i + p] + X1[i + p + bfsize / 2];
    X2[i + p + bfsize / 2] = (X1[i + p] - X1[i + p + bfsize / 2]) * W[i * (1<<k)];
   }
  }
  X  = X1;
  X1 = X2;
  X2 = X;
 }
 
 // 重新排序
 for(j = 0; j < count; j++)
 {
  p = 0;
  for(i = 0; i < r; i++)
  {
   if (j&(1<<i))
   {
    p+=1<<(r-i-1);
   }
  }
  FD[j]=X1[p];
 }
 
 // 释放内存
 delete W;
 delete X1;
 delete X2;
}

函数名称:Fourier()
参数:
LPSTR lpDIBBits    - 指向源DIB图像指针
LONG  lWidth       - 源图像宽度(象素数)
LONG  lHeight      - 源图像高度(象素数)
返回值:BOOL               - 成功返回TRUE,否则返回FALSE。
说明:该函数用来对图像进行付立叶变换。
BOOL WINAPI Fourier(LPSTR lpDIBBits, LONG lWidth, LONG lHeight)
{
 // 指向源图像的指针
 unsigned char* lpSrc;
  // 中间变量
 double dTemp;
  // 循环变量
 LONG i;
 LONG j;
 
 // 进行付立叶变换的宽度和高度(2的整数次方)
 LONG w;
 LONG h;
  int  wp;
 int  hp;
 
 // 图像每行的字节数
 LONG lLineBytes;
 
 // 计算图像每行的字节数
 lLineBytes = WIDTHBYTES(lWidth * 8);
 
 // 赋初值
 w = 1;
 h = 1;
 wp = 0;
 hp = 0;
 
 // 计算进行付立叶变换的宽度和高度(2的整数次方)
 while(w * 2 <= lWidth)
 {
  w *= 2;
  wp++;
 }
 
 while(h * 2 <= lHeight)
 {
  h *= 2;
  hp++;
 }
 
 // 分配内存
 complex<double> *TD = new complex<double>[w * h];
 complex<double> *FD = new complex<double>[w * h];
 
 // 行
 for(i = 0; i < h; i++)
 {
  // 列
  for(j = 0; j < w; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   
   // 给时域赋值
   TD[j + w * i] = complex<double>(*(lpSrc), 0);
  }
 }
 
 for(i = 0; i < h; i++)
 {
  // 对y方向进行快速付立叶变换
  FFT(&TD[w * i], &FD[w * i], wp);
 }
 
 // 保存变换结果
 for(i = 0; i < h; i++)
 {
  for(j = 0; j < w; j++)
  {
   TD[i + h * j] = FD[j + w * i];
  }
 }
 
 for(i = 0; i < w; i++)
 {
  // 对x方向进行快速付立叶变换
  FFT(&TD[i * h], &FD[i * h], hp);
 }
 
 // 行
 for(i = 0; i < h; i++)
 {
  // 列
  for(j = 0; j < w; j++)
  {
   // 计算频谱
   dTemp = sqrt(FD[j * h + i].real() * FD[j * h + i].real() + 
             FD[j * h + i].imag() * FD[j * h + i].imag()) / 100;
   // 判断是否超过255
   if (dTemp > 255)
   {
    // 对于超过的,直接设置为255
    dTemp = 255;
   }
    // 指向DIB第(i<h/2 ? i+h/2 : i-h/2)行,第(j<w/2 ? j+w/2 : j-w/2)个象素的指针
   // 此处不直接取i和j,是为了将变换后的原点移到中心
   //lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * 
    (lHeight - 1 - (i<h/2 ? i+h/2 : i-h/2)) + (j<w/2 ? j+w/2 : j-w/2);
   
   // 更新源图像
   * (lpSrc) = (BYTE)(dTemp);
  }
 }
 
 // 删除临时变量
 delete TD;
 delete FD;
 
 // 返回
 return TRUE;
}

变换效果:

    July附注:此傅里叶变换算法,在本BLOG内有深入具体的介绍,请参考本BLOG内其它文章。

 

六、离散余弦变换
    算法描述:
    离散余弦变换(DCT for Discrete Cosine Transform)是与傅里叶变换相关的一种变换,它类似于离散傅里叶变换(DFT for Discrete Fourier Transform),但是只使用实数。
    离散余弦变换相当于一个长度大概是它两倍的离散傅里叶变换,这个离散傅里叶变换是对一个实偶函数进行的(因为一个实偶函数的傅里叶变换仍然是一个实偶函数),在有些变形里面需要将输入或者输出的位置移动半个单位(DCT有8种标准类型,其中4种是常见的)。

    程序实现:

函数名称:FFT()
参数:
complex<double> * TD - 指向时域数组的指针
complex<double> * FD - 指向频域数组的指针
r         -2的幂数,即迭代次数
返回值:无。
说明:该函数用来实现快速付立叶变换。
VOID WINAPI FFT(complex<double> * TD, complex<double> * FD, int r)
{
 // 付立叶变换点数
 LONG count;
 // 循环变量
 int  i,j,k;
  // 中间变量
 int  bfsize,p;
  // 角度
 double angle;
 
 complex<double> *W,*X1,*X2,*X;
  // 计算付立叶变换点数
 count = 1 << r;
 
 // 分配运算所需存储器
 W  = new complex<double>[count / 2];
 X1 = new complex<double>[count];
 X2 = new complex<double>[count];
 
 // 计算加权系数
 for(i = 0; i < count / 2; i++)
 {
  angle = -i * PI * 2 / count;
  W[i] = complex<double> (cos(angle), sin(angle));
 }
 
 // 将时域点写入X1
 memcpy(X1, TD, sizeof(complex<double>) * count);
 
 // 采用蝶形算法进行快速付立叶变换
 for(k = 0; k < r; k++)
 {
  for(j = 0; j < 1 << k; j++)
  {
   bfsize = 1 << (r-k);
   for(i = 0; i < bfsize / 2; i++)
   {
    p = j * bfsize;
    X2[i + p] = X1[i + p] + X1[i + p + bfsize / 2];
    X2[i + p + bfsize / 2] = (X1[i + p] - X1[i + p + bfsize / 2]) * W[i * (1<<k)];
   }
  }
  X  = X1;
  X1 = X2;
  X2 = X;
 }
 
 // 重新排序
 for(j = 0; j < count; j++)
 {
  p = 0;
  for(i = 0; i < r; i++)
  {
   if (j&(1<<i))
   {
    p+=1<<(r-i-1);
   }
  }
  FD[j]=X1[p];
 }
 
 // 释放内存
 delete W;
 delete X1;
 delete X2;
}

函数名称:DCT()
参数:
double * f    - 指向时域值的指针
double * F    - 指向频域值的指针
r      -2的幂数
返回值:无。
说明:该函数用来实现快速离散余弦变换,利用2N点的快速付立叶变换来实现离散余弦变换。
VOID WINAPI DCT(double *f, double *F, int r)
{
 // 离散余弦变换点数
 LONG count;
 // 循环变量
 int  i;
  // 中间变量
 double dTemp;
 
 complex<double> *X;
  // 计算离散余弦变换点数
 count = 1<<r;
 
 // 分配内存
 X = new complex<double>[count*2];
  // 赋初值为0
 memset(X, 0, sizeof(complex<double>) * count * 2);
  // 将时域点写入数组X
 for(i=0;i<count;i++)
 {
  X[i] = complex<double> (f[i], 0);
 }
 
 // 调用快速付立叶变换
 FFT(X,X,r+1);
 // 调整系数
 dTemp = 1/sqrt(count);
  // 求F[0]
 F[0] = X[0].real() * dTemp;
 dTemp *= sqrt(2);
 // 求F[u] 
 for(i = 1; i < count; i++)
 {
  F[i]=(X[i].real() * cos(i*PI/(count*2)) + X[i].imag() * sin(i*PI/(count*2))) * dTemp;
 }
 
 // 释放内存
 delete X;
}

函数名称:DIBDct()
参数:
LPSTR lpDIBBits    - 指向源DIB图像指针
LONG  lWidth       - 源图像宽度(象素数)
LONG  lHeight      - 源图像高度(象素数)
返回值:BOOL               - 成功返回TRUE,否则返回FALSE。
说明:该函数用来对图像进行离散余弦变换。
BOOL WINAPI DIBDct(LPSTR lpDIBBits, LONG lWidth, LONG lHeight)
{
  // 指向源图像的指针
 unsigned char* lpSrc;
  // 循环变量
 LONG i;
 LONG j;
  // 进行付立叶变换的宽度和高度(2的整数次方)
 LONG w;
 LONG h;
  // 中间变量
 double dTemp;
 int  wp;
 int  hp;
 
 // 图像每行的字节数
 LONG lLineBytes;
  // 计算图像每行的字节数
 lLineBytes = WIDTHBYTES(lWidth * 8);
 
 // 赋初值
 w = 1;
 h = 1;
 wp = 0;
 hp = 0;
 
 // 计算进行离散余弦变换的宽度和高度(2的整数次方)
 while(w * 2 <= lWidth)
 {
  w *= 2;
  wp++;
 }
 
 while(h * 2 <= lHeight)
 {
  h *= 2;
  hp++;
 }
  // 分配内存
 double *f = new double[w * h];
 double *F = new double[w * h];
 
 // 行
 for(i = 0; i < h; i++)
 {
  // 列
  for(j = 0; j < w; j++)
  {
   // 指向DIB第i行,第j个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   
   // 给时域赋值
   f[j + i * w] = *(lpSrc);
  }
 }
 
 for(i = 0; i < h; i++)
 {
  // 对y方向进行离散余弦变换
  DCT(&f[w * i], &F[w * i], wp);
 }
 
 // 保存计算结果
 for(i = 0; i < h; i++)
 {
  for(j = 0; j < w; j++)
  {
   f[j * h + i] = F[j + w * i];
  }
 }
 
 for(j = 0; j < w; j++)
 {
  // 对x方向进行离散余弦变换
  DCT(&f[j * h], &F[j * h], hp);
 }
  // 行
 for(i = 0; i < h; i++)
 {
  // 列
  for(j = 0; j < w; j++)
  {
   // 计算频谱
   dTemp = fabs(F[j*h+i]);
   
   // 判断是否超过255
   if (dTemp > 255)
   {
    // 对于超过的,直接设置为255
    dTemp = 255;
   }
   
   // 指向DIB第y行,第x个象素的指针
   lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
   
   // 更新源图像
   * (lpSrc) = (BYTE)(dTemp);
  }
 }
 
 // 释放内存
 delete f;
 delete F;

 // 返回
 return TRUE;
}


    变化效果:

更多见下一章: 数字图像处理领域的二十四个典型算法及vc实现、第二章。本文完。

 版权所有,侵权必究。若需转载,请注明出处。谢谢。


using System; using System.Collections.Generic; using System.ComponentModel; using System.Data; using System.Drawing; using System.Text; using System.Windows.Forms; namespace edge { public partial class Form1 : Form { public Form1() { InitializeComponent(); } private void open_Click(object sender, EventArgs e) { OpenFileDialog opnDlg = new OpenFileDialog(); opnDlg.Filter = "所有图像文件 | *.bmp; *.pcx; *.png; *.jpg; *.gif;" + "*.tif; *.ico; *.dxf; *.cgm; *.cdr; *.wmf; *.eps; *.emf|" + "位图( *.bmp; *.jpg; *.png;...) | *.bmp; *.pcx; *.png; *.jpg; *.gif; *.tif; *.ico|" + "矢量图( *.wmf; *.eps; *.emf;...) | *.dxf; *.cgm; *.cdr; *.wmf; *.eps; *.emf"; opnDlg.Title = "打开图像文件"; opnDlg.ShowHelp = true; if (opnDlg.ShowDialog() == DialogResult.OK) { curFileName = opnDlg.FileName; try { curBitmap = (Bitmap)Image.FromFile(curFileName); } catch (Exception exp) { MessageBox.Show(exp.Message); } } Invalidate(); } private void close_Click(object sender, EventArgs e) { this.Close(); } private void Form1_Paint(object sender, PaintEventArgs e) { Graphics g = e.Graphics; if (curBitmap != null) { g.DrawImage(curBitmap, 160, 20, curBitmap.Width, curBitmap.Height); } } private void mask_Click(object sender, EventArgs e) { if (curBitmap != null) { mask operatorMask = new mask(); if (operatorMask.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); int thresholding = operatorMask.GetThresholding; byte flagMask = operatorMask.GetMask; double[] tempArray = new double[bytes]; double gradX, gradY, grad; switch (flagMask) { case 0://Roberts for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { gradX = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[i * curBitmap.Width + j]; gradY = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)]; grad = Math.Sqrt(gradX * gradX + gradY * gradY); tempArray[i * curBitmap.Width + j] = grad; } } break; case 1://Prewitt for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { gradX = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] - grayValues[i * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)]; gradY = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] + grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; grad = Math.Sqrt(gradX * gradX + gradY * gradY); tempArray[i * curBitmap.Width + j] = grad; } } break; case 2://Sobel for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { gradX = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 2 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 2 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)]; gradY = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 2 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 2 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; grad = Math.Sqrt(gradX * gradX + gradY * gradY); tempArray[i * curBitmap.Width + j] = grad; } } break; case 3://Laplacian1公式(8.4) for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { grad = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - 4 * grayValues[i * curBitmap.Width + j]; tempArray[i * curBitmap.Width + j] = grad; } } break; case 4://Laplacian2公式(8.5) for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { grad = grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 8 * grayValues[i * curBitmap.Width + j]; tempArray[i * curBitmap.Width + j] = grad; } } break; case 5://Laplacian3公式(8.6) for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { grad = -1 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 2 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 2 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 2 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 2 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 4 * grayValues[i * curBitmap.Width + j]; tempArray[i * curBitmap.Width + j] = grad; } } break; case 6://Kirsch for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { grad = 0; gradX = -5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = -5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; gradX = -5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 5 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - 5 * grayValues[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + j] + 3 * grayValues[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; if (gradX > grad) grad = gradX; tempArray[i * curBitmap.Width + j] = grad; } } break; default: MessageBox.Show("无效!"); break; } if (thresholding == 0)//不进行阈值处理 { for (int i = 0; i < bytes; i++) { if (tempArray[i] < 0) grayValues[i] = 0; else { if (tempArray[i] > 255) grayValues[i] = 255; else grayValues[i] = Convert.ToByte(tempArray[i]); } } } else//阈值处理,生成二值边缘图像 { if (flagMask == 3 || flagMask == 4 || flagMask == 5) { zerocross(ref tempArray, out grayValues, thresholding); } else { for (int i = 0; i < bytes; i++) { if (tempArray[i] > thresholding) grayValues[i] = 255; else grayValues[i] = 0; } } } System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } private void gaussian_Click(object sender, EventArgs e) { if (curBitmap != null) { gaussian gaussFilter = new gaussian(); if (gaussFilter.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); double thresholding = gaussFilter.GetThresholding; double sigma = gaussFilter.GetSigma; bool flag = gaussFilter.GetFlag; double[] filt, tempArray; createFilter(out filt, sigma, flag); conv2(ref grayValues, ref filt, out tempArray); zerocross(ref tempArray, out grayValues, thresholding); System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } private void zerocross(ref double[] inputImage, out byte[] outImage, double thresh) { outImage = new byte[curBitmap.Width * curBitmap.Height]; for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { if (inputImage[i * curBitmap.Width + j] < 0 && inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j] > 0 && Math.Abs(inputImage[i * curBitmap.Width + j] - inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j]) > thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + j] < 0 && inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] > 0 && Math.Abs(inputImage[i * curBitmap.Width + j] - inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j]) > thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + j] < 0 && inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] > 0 && Math.Abs(inputImage[i * curBitmap.Width + j] - inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)]) > thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + j] < 0 && inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] > 0 && Math.Abs(inputImage[i * curBitmap.Width + j] - inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)]) > thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + j] == 0) { if (inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j] > 0 && inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] < 0 && Math.Abs(inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j]) > 2 * thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j] < 0 && inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] > 0 && Math.Abs(inputImage[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] - inputImage[((i + 1) % curBitmap.Height) * curBitmap.Width + j]) > 2 * thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] > 0 && inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] < 0 && Math.Abs(inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)]) > 2 * thresh) { outImage[i * curBitmap.Width + j] = 255; } else if (inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] < 0 && inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] > 0 && Math.Abs(inputImage[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] - inputImage[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)]) > 2 * thresh) { outImage[i * curBitmap.Width + j] = 255; } else { outImage[i * curBitmap.Width + j] = 0; } } else { outImage[i * curBitmap.Width + j] = 0; } } } } private void createFilter(out double[] filter, double sigma, bool lod) { double std2 = 2 * sigma * sigma; int radius = Convert.ToInt16(Math.Ceiling(3 * sigma)); int filterWidth = 2 * radius + 1; filter = new double[filterWidth * filterWidth]; double sum = 0, average = 0; if (lod == false) { for (int i = 0; i < radius; i++) { for (int j = 0; j < radius; j++) { int xx = (j - radius) * (j - radius); int yy = (i - radius) * (i - radius); filter[i * filterWidth + j] = (xx + yy - std2) * Math.Exp(-(xx + yy) / std2); sum += 4 * filter[i * filterWidth + j]; } } for (int i = 0; i < radius; i++) { int xx = (i - radius) * (i - radius); filter[i * filterWidth + radius] = (xx - std2) * Math.Exp(-xx / std2); sum += 2 * filter[i * filterWidth + radius]; } for (int j = 0; j < radius; j++) { int yy = (j - radius) * (j - radius); filter[radius * filterWidth + j] = (yy - std2) * Math.Exp(-yy / std2); sum += 2 * filter[radius * filterWidth + j]; } filter[radius * filterWidth + radius] = -std2; sum += filter[radius * filterWidth + radius]; average = sum / filter.Length; for (int i = 0; i < radius; i++) { for (int j = 0; j < radius; j++) { filter[i * filterWidth + j] = filter[i * filterWidth + j] - average; filter[filterWidth - 1 - j + i * filterWidth] = filter[i * filterWidth + j]; filter[filterWidth - 1 - j + (filterWidth - 1 - i) * filterWidth] = filter[i * filterWidth + j]; filter[j + (filterWidth - 1 - i) * filterWidth] = filter[i * filterWidth + j]; } } for (int i = 0; i < radius; i++) { filter[i * filterWidth + radius] = filter[i * filterWidth + radius] - average; filter[(filterWidth - 1 - i) * filterWidth + radius] = filter[i * filterWidth + radius]; } for (int j = 0; j < radius; j++) { filter[radius * filterWidth + j] = filter[radius * filterWidth + j] - average; filter[radius * filterWidth + filterWidth - 1 - j] = filter[radius * filterWidth + j]; } filter[radius * filterWidth + radius] = filter[radius * filterWidth + radius] - average; } else { for (int i = 0; i < radius; i++) { for (int j = 0; j < radius; j++) { int xx = (j - radius) * (j - radius); int yy = (i - radius) * (i - radius); filter[i * filterWidth + j] = 1.6 * Math.Exp(-(xx + yy) * 1.6 * 1.6 / std2) / sigma - Math.Exp(-(xx + yy) / std2) / sigma; sum += 4 * filter[i * filterWidth + j]; } } for (int i = 0; i < radius; i++) { int xx = (i - radius) * (i - radius); filter[i * filterWidth + radius] = 1.6 * Math.Exp(-xx * 1.6 * 1.6 / std2) / sigma - Math.Exp(-xx / std2) / sigma; sum += 2 * filter[i * filterWidth + radius]; } for (int j = 0; j < radius; j++) { int yy = (j - radius) * (j - radius); filter[radius * filterWidth + j] = 1.6 * Math.Exp(-yy * 1.6 * 1.6 / std2) / sigma - Math.Exp(-yy / std2) / sigma; sum += 2 * filter[radius * filterWidth + j]; } filter[radius * filterWidth + radius] = 1.6 / sigma - 1 / sigma; sum += filter[radius * filterWidth + radius]; average = sum / filter.Length; for (int i = 0; i < radius; i++) { for (int j = 0; j < radius; j++) { filter[i * filterWidth + j] = filter[i * filterWidth + j] - average; filter[filterWidth - 1 - j + i * filterWidth] = filter[i * filterWidth + j]; filter[filterWidth - 1 - j + (filterWidth - 1 - i) * filterWidth] = filter[i * filterWidth + j]; filter[j + (filterWidth - 1 - i) * filterWidth] = filter[i * filterWidth + j]; } } for (int i = 0; i < radius; i++) { filter[i * filterWidth + radius] = filter[i * filterWidth + radius] - average; filter[(filterWidth - 1 - i) * filterWidth + radius] = filter[i * filterWidth + radius]; } for (int j = 0; j < radius; j++) { filter[radius * filterWidth + j] = filter[radius * filterWidth + j] - average; filter[radius * filterWidth + filterWidth - 1 - j] = filter[radius * filterWidth + j]; } filter[radius * filterWidth + radius] = filter[radius * filterWidth + radius] - average; } } private void conv2(ref byte[] inputImage, ref double[] mask, out double[] outImage) { int windWidth = Convert.ToInt16(Math.Sqrt(mask.Length)); int radius = windWidth / 2; double temp; outImage = new double[curBitmap.Width * curBitmap.Height]; for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { temp = 0; for (int x = -radius; x <= radius; x++) { for (int y = -radius; y <= radius; y++) { temp += inputImage[((Math.Abs(i + x)) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j + y)) % curBitmap.Width] * mask[(x + radius) * windWidth + y + radius]; } } outImage[i * curBitmap.Width + j] = temp; } } } private void canny_Click(object sender, EventArgs e) { if (curBitmap != null) { canny cannyOp = new canny(); if (cannyOp.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); byte[] thresholding = new byte[2]; thresholding = cannyOp.GetThresh; double sigma = cannyOp.GetSigma; double[] tempArray;// = new double[bytes]; double[] tempImage = new double[bytes]; double[] grad = new double[bytes]; byte[] aLabel = new byte[bytes]; double[] edgeMap = new double[bytes]; Array.Clear(edgeMap, 0, bytes); double gradX, gradY, angle; int rad = Convert.ToInt16(Math.Ceiling(3 * sigma)); for (int i = 0; i < bytes; i++) tempImage[i] = Convert.ToDouble(grayValues[i]); gaussSmooth(tempImage, out tempArray, sigma); for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { gradX = tempArray[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] + 2 * tempArray[i * curBitmap.Width + ((j + 1) % curBitmap.Width)] + tempArray[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - tempArray[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 2 * tempArray[i * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - tempArray[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)]; gradY = tempArray[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] + 2 * tempArray[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] + tempArray[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] - tempArray[((i + 1) % curBitmap.Height) * curBitmap.Width + ((Math.Abs(j - 1)) % curBitmap.Width)] - 2 * tempArray[((i + 1) % curBitmap.Height) * curBitmap.Width + j] - tempArray[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]; grad[i * curBitmap.Width + j] = Math.Sqrt(gradX * gradX + gradY * gradY); angle = Math.Atan2(gradY, gradX); if ((angle >= -1.178097 && angle < 1.178097) || angle >= 2.748894 || angle < -2.748894) aLabel[i * curBitmap.Width + j] = 0; else if ((angle >= 0.392699 && angle < 1.178097) || (angle >= -2.748894 && angle < -1.963495)) aLabel[i * curBitmap.Width + j] = 1; else if ((angle >= -1.178097 && angle < -0.392699) || (angle >= 1.963495 && angle < 2.748894)) aLabel[i * curBitmap.Width + j] = 2; else aLabel[i * curBitmap.Width + j] = 3; } } for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { switch (aLabel[i * curBitmap.Width + j]) { case 3://水平方向 if (grad[i * curBitmap.Width + j] > grad[((Math.Abs(i - 1))%curBitmap.Height) * curBitmap.Width + j] && grad[i * curBitmap.Width + j] > grad[((i + 1)%curBitmap.Height) * curBitmap.Width + j]) edgeMap[i * curBitmap.Width + j] = grad[i * curBitmap.Width + j]; break; case 2://正45度方向 if (grad[i * curBitmap.Width + j] > grad[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] && grad[i * curBitmap.Width + j] > grad[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]) edgeMap[i * curBitmap.Width + j] = grad[i * curBitmap.Width + j]; break; case 1://负45度方向 if (grad[i * curBitmap.Width + j] > grad[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] && grad[i * curBitmap.Width + j] > grad[((i + 1) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)]) edgeMap[i * curBitmap.Width + j] = grad[i * curBitmap.Width + j]; break; case 0://垂直方向 if (grad[i * curBitmap.Width + j] > grad[i * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] && grad[i * curBitmap.Width + j] > grad[i * curBitmap.Width + ((j + 1) % curBitmap.Width)]) edgeMap[i * curBitmap.Width + j] = grad[i * curBitmap.Width + j]; break; default: return; } } } Array.Clear(grayValues, 0, bytes); for (int i = 0; i < curBitmap.Height; i++) { for (int j =0; j < curBitmap.Width; j++) { if (edgeMap[i * curBitmap.Width + j] > thresholding[0]) { grayValues[i * curBitmap.Width + j] = 255; traceEdge(i, j, edgeMap, ref grayValues, thresholding[1]); } } } System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } private void gaussSmooth(double[] inputImage, out double[] outputImage, double sigma) { double std2 = 2 * sigma * sigma; int radius = Convert.ToInt16(Math.Ceiling(3 * sigma)); int filterWidth = 2 * radius + 1; double[] filter = new double[filterWidth]; outputImage = new double[inputImage.Length]; int length = Convert.ToInt16(Math.Sqrt(inputImage.Length)); double[] tempImage = new double[inputImage.Length]; double sum = 0; for (int i = 0; i < filterWidth; i++) { int xx = (i - radius) * (i - radius); filter[i] = Math.Exp(-xx / std2); sum += filter[i]; } for (int i = 0; i < filterWidth; i++) { filter[i] = filter[i] / sum; } for (int i = 0; i < length; i++) { for (int j = 0; j < length; j++) { double temp = 0; for (int k = -radius; k <= radius; k++) { int rem = (Math.Abs(j + k)) % length; temp += inputImage[i * length + rem] * filter[k + radius]; } tempImage[i * length + j] = temp; } } for (int j = 0; j < length; j++) { for (int i = 0; i < length; i++) { double temp = 0; for (int k = -radius; k <= radius; k++) { int rem = (Math.Abs(i + k)) % length; temp += tempImage[rem * length + j] * filter[k + radius]; } outputImage[i * length + j] = temp; } } } private void traceEdge(int k, int l, double[] inputImage, ref byte[] outputImage, byte thrLow) { int[] kOffset = new int[] { 1, 1, 0, -1, -1, -1, 0, 1 }; int[] lOffset = new int[] { 0, 1, 1, 1, 0, -1, -1, -1 }; int kk, ll; for (int p = 0; p < 8; p++) { kk = k + kOffset[p]; kk = Math.Abs(kk) % curBitmap.Height; ll = l + lOffset[p]; ll = Math.Abs(ll) % curBitmap.Width; if (outputImage[ll * curBitmap.Width + kk] != 255) { if (inputImage[ll * curBitmap.Width + kk] > thrLow) { outputImage[ll * curBitmap.Width + kk] = 255; traceEdge(ll, kk, inputImage, ref outputImage, thrLow); } } } } private void morph_Click(object sender, EventArgs e) { if (curBitmap != null) { morphologic grayMor = new morphologic(); if (grayMor.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); byte[] tempArray1 = new byte[bytes]; byte[] tempArray2 = new byte[bytes]; bool flag = grayMor.GetMethod; double thresh = grayMor.GetThresh; byte[] struEle = new byte[25]; struEle = grayMor.GetStruction; int temp; tempArray1 = grayDelation(grayValues, struEle, curBitmap.Height, curBitmap.Width); tempArray2 = grayErode(grayValues, struEle, curBitmap.Height, curBitmap.Width); for (int i = 0; i < bytes; i++) { if (flag == false) temp = (tempArray1[i] - tempArray2[i]) / 2; else temp = (tempArray1[i] + tempArray2[i] - 2 * grayValues[i]) / 2; if (temp > thresh) grayValues[i] = 255; else grayValues[i] = 0; } System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } private byte[] grayDelation(byte[] grayImage, byte[] se, int tHeight, int tWidth) { byte[] tempImage = new byte[grayImage.Length]; for (int i = 0; i < tHeight; i++) { for (int j = 0; j < tWidth; j++) { int[] cou = new int[]{grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + (Math.Abs(j - 2)) % tWidth] * se[0], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + (Math.Abs(j - 1)) % tWidth] * se[1], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + j] * se[2], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + ((j + 1) % tWidth)] * se[3], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + ((j + 2) % tWidth)] * se[4], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * se[5], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * se[6], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + j] * se[7], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + ((j + 1) % tWidth)] * se[8], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + ((j + 2) % tWidth)] * se[9], grayImage[i * tWidth + (Math.Abs(j - 2) % tWidth)] * se[10], grayImage[i * tWidth + (Math.Abs(j - 1) % tWidth)] * se[11], grayImage[i * tWidth + j] * se[12], grayImage[i * tWidth + ((j + 1) % tWidth)] * se[13], grayImage[i * tWidth + ((j + 2) % tWidth)] * se[14], grayImage[((i + 1) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * se[15], grayImage[((i + 1) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * se[16], grayImage[((i + 1) % tHeight) * tWidth + j] * se[17], grayImage[((i + 1) % tHeight) * tWidth + ((j + 1) % tWidth)] * se[18], grayImage[((i + 1) % tHeight) * tWidth + ((j + 2) % tWidth)] * se[19], grayImage[((i + 2) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * se[20], grayImage[((i + 2) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * se[21], grayImage[((i + 2) % tHeight) * tWidth + j] * se[22], grayImage[((i + 2) % tHeight) * tWidth + ((j + 1) % tWidth)] * se[23], grayImage[((i + 2) % tHeight) * tWidth + ((j + 2) % tWidth)] * se[24]}; int maxim = cou[0]; for (int k = 1; k < 25; k++) { if (cou[k] > maxim) { maxim = cou[k]; } } tempImage[i * tWidth + j] = (byte)maxim; } } return tempImage; } private byte[] grayErode(byte[] grayImage, byte[] se, int tHeight, int tWidth) { byte[] tempImage = new byte[grayImage.Length]; byte[] tempSe = new byte[25]; tempSe = (byte[])se.Clone(); for (int k = 0; k < 25; k++) { if (tempSe[k] == 0) tempSe[k] = 255; } for (int i = 0; i < tHeight; i++) { for (int j = 0; j < tWidth; j++) { int[] cou = new int[]{grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + (Math.Abs(j - 2)) % tWidth] * tempSe[0], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + (Math.Abs(j - 1)) % tWidth] * tempSe[1], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + j] * tempSe[2], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + ((j + 1) % tWidth)] * tempSe[3], grayImage[((Math.Abs(i - 2)) % tHeight) * tWidth + ((j + 2) % tWidth)] * tempSe[4], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * tempSe[5], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * tempSe[6], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + j] * tempSe[7], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + ((j + 1) % tWidth)] * tempSe[8], grayImage[((Math.Abs(i - 1)) % tHeight) * tWidth + ((j + 2) % tWidth)] * tempSe[9], grayImage[i * tWidth + (Math.Abs(j - 2) % tWidth)] * tempSe[10], grayImage[i * tWidth + (Math.Abs(j - 1) % tWidth)] * tempSe[11], grayImage[i * tWidth + j] * tempSe[12], grayImage[i * tWidth + ((j + 1) % tWidth)] * tempSe[13], grayImage[i * tWidth + ((j + 2) % tWidth)] * tempSe[14], grayImage[((i + 1) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * tempSe[15], grayImage[((i + 1) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * tempSe[16], grayImage[((i + 1) % tHeight) * tWidth + j] * tempSe[17], grayImage[((i + 1) % tHeight) * tWidth + ((j + 1) % tWidth)] * tempSe[18], grayImage[((i + 1) % tHeight) * tWidth + ((j + 2) % tWidth)] * tempSe[19], grayImage[((i + 2) % tHeight) * tWidth + (Math.Abs(j - 2) % tWidth)] * tempSe[20], grayImage[((i + 2) % tHeight) * tWidth + (Math.Abs(j - 1) % tWidth)] * tempSe[21], grayImage[((i + 2) % tHeight) * tWidth + j] * tempSe[22], grayImage[((i + 2) % tHeight) * tWidth + ((j + 1) % tWidth)] * tempSe[23], grayImage[((i + 2) % tHeight) * tWidth + ((j + 2) % tWidth)] * tempSe[24]}; int minimum = cou[0]; for (int k = 1; k < 25; k++) { if (cou[k] < minimum) { minimum = cou[k]; } } tempImage[i * tWidth + j] = (byte)minimum; } } return tempImage; } private void wavelet_Click(object sender, EventArgs e) { if (curBitmap != null) { wvl wavelet = new wvl(); if (wavelet.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); double[] tempArray1 = new double[bytes]; double[] tempArray2 = new double[bytes]; double[] tempArray3 = new double[bytes]; double[] gradX = new double[bytes]; double[] gradY = new double[bytes]; byte multiscale = wavelet.GetScale; byte[] thresholding = new byte[2]; thresholding = wavelet.GetThresh; for (int i = 0; i < bytes; i++) { tempArray1[i] = Convert.ToDouble(grayValues[i]); } for (int z = 0; z <= multiscale; z++) { double[] p = null; double[] q = null; switch (z) { case 0: p = new double[] { 0.125, 0.375, 0.375, 0.125 }; q = new double[] { -2, 2 }; break; case 1: p = new double[] { 0.125, 0, 0.375, 0, 0.375, 0, 0.125 }; q = new double[] { -2, 0, 2 }; break; case 2: p = new double[] { 0.125, 0, 0, 0, 0.375, 0, 0, 0, 0.375, 0, 0, 0, 0.125 }; q = new double[] { -2, 0, 0, 0, 2 }; break; case 3: p = new double[] { 0.125, 0, 0, 0, 0, 0, 0, 0, 0.375, 0, 0, 0, 0, 0, 0, 0, 0.375, 0, 0, 0, 0, 0, 0, 0, 0.125 }; q = new double[] { -2, 0, 0, 0, 0, 0, 0, 0, 2 }; break; default: return; } int coff = Convert.ToInt16(Math.Pow(2, z) - 1); for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { double[] scl = new double[curBitmap.Width]; double[] wvl = new double[curBitmap.Width]; int temp; scl[j] = 0.0; wvl[j] = 0.0; for (int x = -2 - 2 * coff; x < p.Length - 2 - 2 * coff; x++) { temp = (Math.Abs(j + x)) % curBitmap.Width; scl[j] += p[1 + coff - x] * tempArray1[i * curBitmap.Width + temp]; } for (int x = -1 - coff; x < q.Length - 1 - coff; x++) { temp = (Math.Abs(j + x)) % curBitmap.Width; wvl[j] += q[-x] * tempArray1[i * curBitmap.Width + temp]; } tempArray2[i * curBitmap.Width + j] = scl[j]; gradX[i * curBitmap.Width + j] = wvl[j]; } } for (int i = 0; i < curBitmap.Width; i++) { for (int j = 0; j < curBitmap.Height; j++) { double[] scl = new double[curBitmap.Height]; double[] wvl = new double[curBitmap.Height]; int temp; scl[j] = 0.0; wvl[j] = 0.0; for (int x = -2 - 2 * coff; x < p.Length - 2 - 2 * coff; x++) { temp = (Math.Abs(j + x)) % curBitmap.Height; scl[j] += p[1 + coff - x] * tempArray2[temp * curBitmap.Width + i]; } for (int x = -1 - coff; x < q.Length - 1 - coff; x++) { temp = (Math.Abs(j + x)) % curBitmap.Height; wvl[j] += q[-x] * tempArray1[temp * curBitmap.Width + i]; } tempArray3[j * curBitmap.Width + i] = scl[j]; gradY[j * curBitmap.Width + i] = wvl[j]; } } tempArray1 = (double[])tempArray3.Clone(); } double angle; for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { tempArray1[i * curBitmap.Width + j] = Math.Sqrt(gradX[i * curBitmap.Width + j] * gradX[i * curBitmap.Width + j] + gradY[i * curBitmap.Width + j] * gradY[i * curBitmap.Width + j]); angle = Math.Atan2(gradY[i * curBitmap.Width + j], gradX[i * curBitmap.Width + j]); if ((angle >= -1.178097 && angle < 1.178097) || angle >= 2.748894 || angle < -2.748894) tempArray2[i * curBitmap.Width + j] = 0; else if ((angle >= 0.392699 && angle < 1.178097) || (angle >= -2.748894 && angle < -1.963495)) tempArray2[i * curBitmap.Width + j] = 1; else if ((angle >= -1.178097 && angle < -0.392699) || (angle >= 1.963495 && angle < 2.748894)) tempArray2[i * curBitmap.Width + j] = 2; else tempArray2[i * curBitmap.Width + j] = 3; } } Array.Clear(tempArray3, 0, bytes); for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { switch (Convert.ToInt16(tempArray2[i * curBitmap.Width + j])) { case 3://水平方向 if (tempArray1[i * curBitmap.Width + j] > tempArray1[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + j] && tempArray1[i * curBitmap.Width + j] > tempArray1[((i + 1) % curBitmap.Height) * curBitmap.Width + j]) tempArray3[i * curBitmap.Width + j] = tempArray1[i * curBitmap.Width + j]; break; case 1://正45度方向 if (tempArray1[i * curBitmap.Width + j] > tempArray1[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] && tempArray1[i * curBitmap.Width + j] > tempArray1[((i + 1) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)]) tempArray3[i * curBitmap.Width + j] = tempArray1[i * curBitmap.Width + j]; break; case 2://负45度方向 if (tempArray1[i * curBitmap.Width + j] > tempArray1[((Math.Abs(i - 1)) % curBitmap.Height) * curBitmap.Width + ((j + 1) % curBitmap.Width)] && tempArray1[i * curBitmap.Width + j] > tempArray1[((i + 1) % curBitmap.Height) * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)]) tempArray3[i * curBitmap.Width + j] = tempArray1[i * curBitmap.Width + j]; break; case 0://垂直方向 if (tempArray1[i * curBitmap.Width + j] > tempArray1[i * curBitmap.Width + (Math.Abs(j - 1) % curBitmap.Width)] && tempArray1[i * curBitmap.Width + j] > tempArray1[i * curBitmap.Width + ((j + 1) % curBitmap.Width)]) tempArray3[i * curBitmap.Width + j] = tempArray1[i * curBitmap.Width + j]; break; default: return; } } } Array.Clear(grayValues, 0, bytes); for (int i = 0; i < curBitmap.Height; i++) { for (int j = 0; j < curBitmap.Width; j++) { if (tempArray3[i * curBitmap.Width + j] > thresholding[0]) { grayValues[i * curBitmap.Width + j] = 255; traceEdge(i, j, tempArray3, ref grayValues, thresholding[1]); } } } System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } private void pyramid_Click(object sender, EventArgs e) { if (curBitmap != null) { int series = Convert.ToInt16(Math.Log(curBitmap.Width, 2)); glp pyramid = new glp(series); if (pyramid.ShowDialog() == DialogResult.OK) { Rectangle rect = new Rectangle(0, 0, curBitmap.Width, curBitmap.Height); System.Drawing.Imaging.BitmapData bmpData = curBitmap.LockBits(rect, System.Drawing.Imaging.ImageLockMode.ReadWrite, curBitmap.PixelFormat); IntPtr ptr = bmpData.Scan0; int bytes = curBitmap.Width * curBitmap.Height; byte[] grayValues = new byte[bytes]; System.Runtime.InteropServices.Marshal.Copy(ptr, grayValues, 0, bytes); double thresh = pyramid.GetThresh; byte level = pyramid.GetLevel; double sigma = pyramid.GetSigma; double[][] pyramidImage = new double[level + 1][]; double[][] passImage = new double[level + 1][]; int levelBytes = bytes; for (int k = 0; k < level + 1; k++) { passImage[k] = new double[levelBytes]; pyramidImage[k] = new double[levelBytes]; levelBytes = levelBytes / 4; } for (int i = 0; i < bytes; i++) pyramidImage[0][i] = Convert.ToDouble(grayValues[i]); for (int k = 0; k < level; k++) { double[] tempImage = null; gaussSmooth(pyramidImage[k], out tempImage, sigma); int coff = pyramidImage[k].Length; for (int i = 0; i < coff; i++) { passImage[k][i] = pyramidImage[k][i] - tempImage[i]; int div = i / Convert.ToInt16(Math.Sqrt(coff)); int rem = i % Convert.ToInt16(Math.Sqrt(coff)); if (div % 2 == 0 && rem % 2 == 0) { int j = (int)((div / 2) * Math.Sqrt(pyramidImage[k + 1].Length) + rem / 2); pyramidImage[k + 1][j] = tempImage[i]; } } } for (int k = level - 1; k >= 0; k--) { int coff = pyramidImage[k].Length; for (int i = 0; i < coff; i++) { int div = i / Convert.ToInt16(Math.Sqrt(coff)); int rem = i % Convert.ToInt16(Math.Sqrt(coff)); int j = (int)((div / 2) * Math.Sqrt(pyramidImage[k + 1].Length) + rem / 2); if (div % 2 == 0 && rem % 2 == 0) pyramidImage[k][i] = pyramidImage[k + 1][j]; else pyramidImage[k][i] = 0; } double[] tempImage = null; gaussSmooth(pyramidImage[k], out tempImage, 1); for (int i = 0; i < coff; i++) pyramidImage[k][i] = tempImage[i] + passImage[k][i]; } zerocross(ref pyramidImage[0], out grayValues,thresh); System.Runtime.InteropServices.Marshal.Copy(grayValues, 0, ptr, bytes); curBitmap.UnlockBits(bmpData); } Invalidate(); } } } }
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值