人体骨骼核扫描图的混合空间增强

1、数字图像处理混合空间增强法学习报告学生姓名： 邹晓敏学 号：6103313025 专业班级：数媒131一、 学习目的与内容：对混合模糊图像进行处理， 为得到一个满意的结果，对其应用多种互补的图像增强技术， 最终得到良好效果。本次的数字图像处理实验是对人体骨骼核扫描图片进行混合增强处理， 通过混合增强突出骨骼的更多细节。原图中骨骼比较模糊，边缘不够清晰，对比度不够，对此应该采取的策略是首先用拉普拉斯变换突出图像中的小细节，然后用梯度法突出其边缘。平滑后的拉普拉斯变换将用于掩蔽拉普拉斯图像，最后用灰度变换来扩展图像的灰度动态范围。二、算法介绍：软件开发环境：基于vc+ MFC实现软件编写。具体步骤

2、：(1)拉普拉斯算子拉普拉斯算子是最简单的各向同性微分算子，具有旋转不变性。一个二维图像函数的拉普拉斯变换是各向同性的二阶导数，定义为：为了更适合于数字图像处理，将该方程表示为离散形式+共工一1,¥)+/(工,¥+1)十六工/一1)-4/(工,中另外，拉普拉斯算子还可以表示成模板的形式，一般增强技术对于陡峭的边缘和缓慢变化的边缘很难确定其边缘线的位置。 定，对孤立点或端点更为敏感， 目的的场合。1a1.。A 1orlfo3Oi10-lj4r-I0,工 0但此算子却可用二次微分正峰和负峰之间的过零点来确因此特别适用于以突出图像中的孤立点、孤立线或线端点为21, lr 1lx电

3、 1L 14 Tf -L-L 瓦-1-1 r ' lr -1用以上算子与原图进行卷积，得到中心像素点与邻域像素点的差值。得到的拉普拉斯变换图M1 , M1与原图相加得到 M2,锐化原图。(2)边缘检测一幅图像的边缘是通过一阶和二阶数字导数来检测得到的。边缘的宽度取决于从初始灰度级跃变到最终灰度级的斜坡的长度。 这个长度又取决于斜度，而斜度又取决于模糊程度。 所以，我们可知，模糊的边缘使其变粗，而清晰的边缘使其变得较细。一幅数字图像的一阶导数是基于各种二维梯度的近似值。图像f (x, y)在位置(x, y)的梯度定义为下列向量：13梯度向量指向在坐标 (x，y)的f的最大变化率方向。在边

4、缘检测中，一个重要的量是向量的大小，用 Vf表示，。一般来讲称' f为梯度。考虑到图片中的噪声，所以在此处+1+ 2+ 1000-1-2-1Gy该算子包含两组3x3的矩阵，分别为横向及纵向， 将之与图像作平面卷积，即可分别得出横向及纵向的亮度差分近似值。如果以A代表原始图像，Gx及Gy分别代表经横向及纵向边缘检测的图像灰度值，其公式如下：-1 0 +1G*= 2 0 +2 *A-1 0 +1+1 +2 +1and Gv =000 * A-1 -2 -1sobel算子运算结果得到 M3,保留边缘去除噪声根据最后运行的结果显示，sobel图像的边缘要比拉普拉斯图像的边缘突出很多。拉普拉斯变

5、换作为一种二阶微分算子，在图像细节的增强处理方面有明显的优点，但拉普拉斯变换与梯度变换相比较会产生更多的噪声。其中位于平滑区域的噪声非常显眼，梯度变换在灰度变换的区域的相应要比拉普拉斯更强烈，而梯度变换对噪声和小细节的响应要比拉普拉斯变换弱，而且可以通过均值滤波器对其进行平滑处理而进一步的降低，对梯度图像进行平滑处理并用拉普拉斯图像与他相乘。处理后的结果在灰度变化强的区域仍然保留细节，而在灰 度变化相对平坦的区域则减少噪声， 这种处理可以看做是将拉普拉斯变换与梯度变换的有点结合。（ 3）平滑处理均值滤波是典型的线性滤波算法， 它是指在图像上对目标像素给一个模板， 该模板包括了其周围的临近像素

6、（以目标象素为中心的周围 8 个像素， 构成一个滤波模板， 即去掉目标 像素本身） ，再用模板中的全体像素的平均值来代替原来像素值。一个 5*5 的均值滤波器对图像进行平滑处理， 5*5 的均值滤波器对图像进行平滑处理得到图像M4再由M2和M4相乘形成掩蔽图像M5（ 4）锐化图像拉普拉斯与平滑后的梯度图像乘积的结果可以看出，强边缘的优势和可见噪声的减少。 将乘积图像与原始图像相加就得到了锐化图像。 与原始图像相比， 锐化处理后的图像细节增加很明显。 单独使用拉普拉斯或者梯度变换不能达到这种效果。 而且以上的锐化过程没有影响图像的灰度动态变化范围。 进行增强处理的最后一步就是扩大锐化图像的灰度动

7、态范围。 将 M5 和原图 M 求和锐化得到 M6 ， 再对 M6 做幂律变换， 扩散灰度范围， 提升对比度，但幂次必须小于1。得到最后的混合空间加强图像 M7 。三、 源程序核心代码及运行结果/ 拉普拉斯 (Laplacian )算子void CImageProcessingDoc:LaplaceOperator(CImage * pSrcImage, CImage * pDstImage, int laplacetype,bool bCalibration) int x;int y;int i;int j;int mWidth = pSrcImage->GetWidth();int

8、mHeight = pSrcImage->GetHeight();int mRowBytes = pSrcImage->GetPitch();BYTE *rImageData = (BYTE *)pDstImage->GetBits();int *pBitData = NULL;pBitData = new int3 * mHeight*mWidth;BYTE * pData = NULL;pData = new BYTE9;int Laplace033 = 0,1,0,1,-4,1 , 0,1,0 ;int Laplace133 = 1,1,1,1,-8, 1, 1,1,1

9、 ;int Laplace233 = 0,-1,0,-1,4,-1, 0,-1,0 ;int Laplace333 = -1,-1,-1,-1,8,-1,-1,-1,-1 ;int *p;switch (laplacetype)case 1: p = (int *)Laplace0;break;case 2: p = (int *)Laplace1;break;case 3: p = (int *)Laplace2;break;case 4: p = (int *)Laplace3;break;default:return;int laplace;for (y = 0; y < mHei

10、ght; y+)for (x = 0; x < mWidth; x+)CopyAreaBlueFromImgBits(pSrcImage, pData, x, y, 1, 1);laplace = 0;for (j = 0; j < 3; j+)for (i = 0; i < 3; i+) laplace += pDataj * 3 + i * pj*3 + i;pBitData3 * y * mWidth + x * 3 = laplace;CopyAreaGreenFromImgBits(pSrcImage, pData, x, y, 1, 1); laplace = 0

11、;for (j = 0; j < 3; j+)for (i = 0; i < 3; i+)laplace += pDataj * 3 + i * pj*3 + i;pBitData3 * y * mWidth + x * 3 + 1 = laplace;CopyAreaRedFromImgBits(pSrcImage, pData, x, y, 1, 1);laplace = 0;for (j = 0; j < 3; j+)for (i = 0; i < 3; i+)laplace += pDataj * 3 + i * pj*3 + i;pBitData3 * y *

12、 mWidth + x * 3 + 2 = laplace;int minBlue = pBitData0;int minGreen = pBitData1;int minRed = pBitData2;for (y = 0; y < mHeight; y+)for (x = 0; x < mWidth; x+)if (pBitData3 * y * mWidth + x * 3 < minBlue)minBlue = pBitData3 * y * mWidth + x * 3;if (pBitData3 * y * mWidth + x * 3 + 1 < minG

13、reen)minGreen = pBitData3 * y * mWidth + x * 3 + 1;if (pBitData3 * y * mWidth + x * 3 + 2 < minRed) minRed = pBitData3 * y * mWidth + x * 3 + 2;for (y = 0; y < mHeight; y+)for (x = 0; x < mWidth; x+)laplace = pBitData3 * y * mWidth + x * 3;laplace = bCalibration ?(laplace + minBlue>255?2

14、55: laplace + minBlue): (laplace > 255) ? 255 : (laplace < 0 ? 0 : laplace);rImageDatay*mRowBytes + x * 3 = laplace;laplace = pBitData3 * y * mWidth + x * 3 + 1;laplace = bCalibration ? (laplace + minBlue>255 ? 255 : laplace + minBlue) (laplace > 255) ? 255 : (laplace < 0 ? 0 : laplac

15、e);rImageDatay*mRowBytes + x * 3 + 1 = laplace;laplace = pBitData3 * y * mWidth + x * 3 + 2;laplace = bCalibration ? (laplace + minBlue>255 ? 255 : laplace + minBlue) (laplace > 255) ? 255 : (laplace < 0 ? 0 : laplace);rImageDatay*mRowBytes + x * 3 + 2 = laplace;if (pData)delete pData;pData

16、 = NULL;if (pBitData)delete pBitData;pBitData = NULL;/Sobel 算子void CImageProcessingDoc:SobelOperator(CImage * pSrcImage, CImage * pDstImage)int x;int y;int dx;int dy;int sobel;int mWidth = pSrcImage->GetWidth();int mHeight = pSrcImage->GetHeight();int mRowBytes = pSrcImage->GetPitch();BYTE

17、*sImageData = (BYTE *)pSrcImage->GetBits();BYTE *rImageData = (BYTE *)pDstImage->GetBits();for (y = 1; y < mHeight - 1; y+)for (x = 1; x < mWidth - 1; x+)dx = sImageData(y - 1)*mRowBytes + (x - 1) * 3 + 0 - sImageData(y + 1)*mRowBytes + (x - 1) * 3 + 0+ 2 * (sImageData(y - 1)*mRowBytes +

18、 x * 3 + 0 - sImageData(y + 1)*mRowBytes + x * 3 + 0)+ sImageData(y - 1)*mRowBytes + (x + 1) * 3 + 0 - sImageData(y + 1)*mRowBytes + (x + 1) * 3 + 0;dy = sImageData(y - 1)*mRowBytes + (x - 1) * 3 + 0 - sImageData(y - 1)*mRowBytes + (x + 1) * 3 + 0+ 2 * (sImageDatay*mRowBytes + (x - 1) * 3 + 0 - sIma

19、geDatay*mRowBytes +(x - 1) * 3 + 0)+ sImageData(y + 1)*mRowBytes + (x - 1) * 3 + 0 - sImageData(y + 1)*mRowBytes + (x + 1) * 3 + 0;sobel = (int)(sqrt(dx*dx*1.0 + dy*dy);sobel = (sobel > 255) ? 255 : sobel;rImageDatay*mRowBytes + x * 3 + 0 = sobel;dx = sImageData(y - 1)*mRowBytes + (x - 1) * 3 + 1

20、 - sImageData(y + 1)*mRowBytes + (x - 1) * 3 + 1+ 2 * (sImageData(y - 1)*mRowBytes + x * 3 + 1 - sImageData(y +1)*mRowBytes + x * 3 + 1)+ sImageData(y - 1)*mRowBytes + (x + 1) * 3 + 1 - sImageData(y + 1)*mRowBytes + (x + 1) * 3 + 1;dy = sImageData(y - 1)*mRowBytes + (x - 1) * 3 + 1 - sImageData(y -

21、1)*mRowBytes + (x + 1) * 3 + 1+ 2 * (sImageDatay*mRowBytes + (x - 1) * 3 + 1 - sImageDatay*mRowBytes +(x - 1) * 3 + 1)+ sImageData(y + 1)*mRowBytes + (x - 1) * 3 + 1 - sImageData(y + 1)*mRowBytes + (x + 1) * 3 + 1;sobel = (int)(sqrt(dx*dx*1.0 + dy*dy);sobel = (sobel > 255) ? 255 : sobel;1) * 3 +

22、2 - sImageData(y +rImageDatay*mRowBytes + x * 3 + 1 = sobel;dx = sImageData(y - 1)*mRowBytes + (x 1)*mRowBytes + (x - 1) * 3 + 2+ 2 * (sImageData(y1)*mRowBytes + x * 3 + 2)1)*mRowBytes+ x * 3 + 2 - sImageData(y +sImageData(y +sImageData(y -+ sImageData(y - 1)*mRowBytes + (x + 1) * 3 + 2 1)*mRowBytes

23、 + (x + 1) * 3 + 2;dy = sImageData(y - 1)*mRowBytes + (x - 1) * 3 + 21)*mRowBytes + (x + 1) * 3 + 2+ 2 * (sImageDatay*mRowBytes + (x - 1) * 3 + 2 - sImageDatay*mRowBytes +(x - 1) * 3 + 2)+ sImageData(y + 1)*mRowBytes + (x - 1) * 3 + 2 - sImageData(y +1)*mRowBytes + (x + 1) * 3 + 2;sobel = (int)(sqrt

24、(dx*dx*1.0 + dy*dy);sobel = (sobel > 255) ? 255 : sobel;rImageDatay*mRowBytes + x * 3 + 2 = sobel;/均值滤波void CImageProcessingDoc:MeanOperator(CImage * pSrcImage, CImage * pDstImage, int m, int n)int x;int y;int mWidth = pSrcImage->GetWidth();int mHeight = pSrcImage->GetHeight();int mRowBytes

25、 = pSrcImage->GetPitch();BYTE *rImageData = (BYTE *)pDstImage->GetBits();double sumB;double sumG;double sumR;int nums = (2 * m + 1)*(2 * n + 1);BYTE * pData = NULL;pData = new BYTE3 * nums;int lineBytes = 3 * (2 * n + 1);for (y = 0; y < mHeight; y+) for (x = 0; x < mWidth; x+)CopyAreaRGB

26、FromImgBits(pSrcImage, pData, x, y, m, n); sumB = sumG = sumR = 0;for (int j = -m; j <= m; j+)for (int i = -n; i <= n; i+)sumB += pData(j + m)*lineBytes + (i + n) * 3 + 0;sumG += pData(j + m)*lineBytes + (i + n) * 3 + 1;sumR += pData(j + m)*lineBytes + (i + n) * 3 + 2;rImageDatay*mRowBytes + x

27、 * 3 = (int)(sumB/nums);rImageDatay*mRowBytes + x * 3 + 1 = (int)(sumG / nums);rImageDatay*mRowBytes + x * 3 + 2 = (int)(sumR / nums);if (pData)delete pData; pData = NULL;/幂律变换void CImageProcessingDoc:OnPowerTransform()/ TODO: 在此添加命令处理程序代码if (m_sImage.IsNull()AfxMessageBox(_T(" 请打开要处理的图像！ "

28、;); return;CDialogPowerParameters mdlg;if (mdlg.DoModal() = IDCANCEL)return;int x;int y;int nRows = m_sImage.GetHeight();int nCols = m_sImage.GetWidth();int nBPP = m_sImage.GetBPP();BYTE nRed, nGreen, nBlue;double dGrey;m_rImage.Destroy();if (!m_rImage.Create(nCols, nRows, nBPP) return;BYTE *sImageD

29、ata = (BYTE *)m_sImage.GetBits();BYTE *rImageData = (BYTE *)m_rImage.GetBits();int nRowBytes = m_sImage.GetPitch();for (y = 0; y < nRows; y+)for (x = 0; x < nCols; x+)nBlue = sImageDatay * nRowBytes + x * 3 + 0;nGreen = sImageDatay * nRowBytes + x * 3 + 1;nRed = sImageDatay * nRowBytes + x * 3

30、 + 2;dGrey = 0.299*nRed + 0.587*nGreen + 0.114*nBlue;dGrey = 255 * mdlg.m_power_c * exp(mdlg.m_power_gama * log(dGrey / 255 + 0.00001);rImageDatay * nRowBytes + x * 3 = (BYTE)dGrey;rImageDatay * nRowBytes + x * 3 + 1 = (BYTE)dGrey;rImageDatay * nRowBytes + x * 3 + 2 = (BYTE)dGrey;UpdateAllViews(NULL

31、);/图像相加void CImageProcessingDoc:OnImageAdd() CImage lsImage;CString strFilter;CSimpleArray<GUID> aguidFileTypes;HRESULT hResult;/ 获取 CImage 支持的图像文件的过滤字符串hResult = m_sImage.GetExporterFilterString(strFilter, aguidFileTypes, _T("All Image Files");if (FAILED(hResult)AfxMessageBox(_T(&qu

32、ot;GetExporterFilter 函数失败 !"); return;CFileDialog dlg(TRUE, NULL, NULL, OFN_FILEMUSTEXIST, strFilter);if (IDOK != dlg.DoModal()return;lsImage.Destroy();/ 将外部图像文件装载到 CImage 对象中hResult = lsImage.Load(dlg.GetPathName();if( ( lsImage.GetHeight() != m_sImage.GetHeight() ) | ( lsImage.GetWidth() != m

33、_sImage.GetWidth() ) )lsImage.Destroy();AfxMessageBox(_T(" 本软件只处理相同大小图像的逻辑与代数运算！ ");return;CDialogImageAddRatio mdlg;if (mdlg.DoModal() != IDOK)lsImage.Destroy();return;double lambuda = mdlg.addRatio / 100.0;int x; int y;int nRows = m_sImage.GetHeight();int nCols = m_sImage.GetWidth();int

34、nBPP = m_sImage.GetBPP();m_rImage.Destroy();if (!m_rImage.Create(nCols, nRows, nBPP)lsImage.Destroy();return;BYTE *sImageData = (BYTE *)m_sImage.GetBits();BYTE *lsImageData = (BYTE *)lsImage.GetBits();BYTE *rImageData = (BYTE *)m_rImage.GetBits();int nRowBytes = m_sImage.GetPitch();int sRowBytes = m

35、_sImage.GetPitch();for (y = 0; y < nRows; y+)for (x = 0; x < nCols ; x+)rImageDatay * nRowBytes + x * 3 + 0 = (int)(lambuda * sImageDatay * nRowBytes + x * 3 + 0 + (1- lambuda) * lsImageDatay * nRowBytes + x * 3 + 0 );rImageDatay * nRowBytes + x * 3 + 1 = (int)(lambuda * sImageDatay * nRowByte

36、s + x * 3 + 1 + (1 - lambuda) * lsImageDatay * nRowBytes + x * 3 + 1);rImageDatay * nRowBytes + x * 3 + 2 = (int)(lambuda * sImageDatay * nRowBytes + x * 3 + 2 + (1 - lambuda) * lsImageDatay * nRowBytes + x * 3 + 2);lsImage.Destroy();UpdateAllViews(NULL);/图像相乘void CImageProcessingDoc:OnImageMul() CI

37、mage lsImage;CString strFilter;CSimpleArray<GUID> aguidFileTypes;HRESULT hResult;/ 获取 CImage 支持的图像文件的过滤字符串hResult = m_sImage.GetExporterFilterString(strFilter, aguidFileTypes, _T("All Image Files");if (FAILED(hResult)AfxMessageBox(_T("GetExporterFilter 函数失败 !"); return; CFi

38、leDialog dlg(TRUE, NULL, NULL, OFN_FILEMUSTEXIST, strFilter);if (IDOK != dlg.DoModal()return;lsImage.Destroy();/ 将外部图像文件装载到 CImage 对象中hResult = lsImage.Load(dlg.GetPathName();int x;int y;int nRows = m_sImage.GetHeight();int nCols = m_sImage.GetWidth();int nBPP = m_sImage.GetBPP();m_rImage.Destroy();if (!m_rImage.Create(nCols, nRows, nBPP)lsIm