計算灰階共生矩陣
共生矩陣用兩個位置的象素的聯合機率密度來定義,它不僅反映亮度的分布特性,也反映具有同樣亮度或接近亮度的象素之間的位置分布特性,是有關圖象亮度變化的二階統計特徵。它是定義一組紋理特徵的基礎。
灰階共生矩陣能反映出圖象灰階關於方向、相鄰間隔、變化幅度的綜合資訊。設f(x,y)為一幅二維數字圖象,其大小為M×N,灰階層級為Ng,則滿足一定空間關係的灰階共生矩陣為:
P(i,j)=#{(x1,y1),(x2,y2)∈M×N|f(x1,y1)=i,f(x2,y2)=j}
其中#(x)表示集合x中的元素個數,顯然P為Ng×Ng的矩陣,若(x1,y1)與(x2,y2)間距離為d,兩者與座標橫軸的夾角為θ,則可以得到各種間距及角度的灰階共生矩陣P(i,j,d,θ)。
class CGlcm<br />{<br />public:<br />CGlcm(void);<br />virtual ~CGlcm(void);<br />public:<br />int *m_pMat1;<br />int *m_pMat2;<br />int *m_pMat3;<br />int *m_pMat4;<br />protected:<br />int *m_pLine[4][256];<br />int m_nMin;<br />int m_nMax;<br />int m_nSum;<br />public:<br />// 計算共生矩陣<br />// 參數:<br />// 1. pImageData: 映像資料指標,單通道,8位。<br />// 2. nLeft,nTop,nWidth,nHeight: 計算共生矩陣的地區。<br />// 3. nWidthStep: 行位移量。<br />// 4. nScale: 尺度。<br />// 5. nReduction: 灰階級壓縮。<br />bool CalGlcm(unsigned char *pImageData, int nLeft, int nTop, int nWidth, int nHeight, int nWidthStep, int nScale, int nReduction);<br />};<br />// 灰階共生矩陣<br />CGlcm::CGlcm(void)<br />: m_pMat1(NULL)<br />, m_pMat2(NULL)<br />, m_pMat3(NULL)<br />, m_pMat4(NULL)<br />, m_nMin(0)<br />, m_nMax(0)<br />, m_nSum(0)<br />{<br />int i, j;<br />unsigned int nSize = sizeof(int) * 256 * 256;<br />// 建立共生矩陣<br />m_pMat1 = (int*) malloc(nSize);<br />m_pMat2 = (int*) malloc(nSize);<br />m_pMat3 = (int*) malloc(nSize);<br />m_pMat4 = (int*) malloc(nSize);<br />if (m_pMat1 && m_pMat2 && m_pMat3 && m_pMat4)<br />{<br />for (i = 0, j = 0; i < 256; i++, j += 256)<br />{<br />m_pLine[0][i] = m_pMat1 + j;<br />m_pLine[1][i] = m_pMat2 + j;<br />m_pLine[2][i] = m_pMat3 + j;<br />m_pLine[3][i] = m_pMat4 + j;<br />}<br />}<br />}<br />CGlcm::~CGlcm(void)<br />{<br />// 釋放共生矩陣<br />if (m_pMat1) free(m_pMat1);<br />if (m_pMat2) free(m_pMat2);<br />if (m_pMat3) free(m_pMat3);<br />if (m_pMat4) free(m_pMat4);<br />}<br />// 計算共生矩陣<br />bool CGlcm::CalGlcm(unsigned char *pImageData, int nLeft, int nTop, int nWidth, int nHeight, int nWidthStep, int nScale, int nReduction)<br />{<br />bool bResult = false;<br />int x0, x1, x2;<br />int y0, y1, y2;<br />int nGray;<br />unsigned int nSize;<br />unsigned char *pLine[3];<br />if (pImageData)<br />{<br />// 灰階最值<br />m_nMin = 0xFF;<br />m_nMax = 0;<br />pLine[1] = pImageData + nWidthStep * nTop + nLeft;<br />for (y1 = 0; y1 < nHeight; y1++)<br />{<br />for (x1 = 0; x1 < nWidth; x1++)<br />{<br />// 灰階級壓縮<br />if (nReduction > 0)<br />{<br />pLine[1][x1] = pLine[1][x1] >> nReduction;<br />}<br />nGray = pLine[1][x1];<br />if (nGray < m_nMin)<br />{<br />m_nMin = nGray;<br />}<br />else if (nGray > m_nMax)<br />{<br />m_nMax = nGray;<br />}<br />}<br />pLine[1] += nWidthStep;<br />}<br />// 累加和<br />m_nSum = nWidth * nHeight * 2;<br />if (m_nMax >= m_nMin)<br />{<br />// 清空記憶體<br />nSize = sizeof(int) * (m_nMax - m_nMin + 1);<br />for (y1 = m_nMin; y1 <= m_nMax; y1++)<br />{<br />memset(&m_pLine[0][y1][m_nMin], 0, nSize);<br />memset(&m_pLine[1][y1][m_nMin], 0, nSize);<br />memset(&m_pLine[2][y1][m_nMin], 0, nSize);<br />memset(&m_pLine[3][y1][m_nMin], 0, nSize);<br />}<br />// 計算共生矩陣<br />pLine[1] = pImageData + nWidthStep * nTop + nLeft;<br />for (y0 = -nScale, y1 = 0, y2 = nScale; y1 < nHeight; y0++, y1++, y2++)<br />{<br />pLine[0] = pLine[1] - nWidthStep;<br />pLine[2] = pLine[1] + nWidthStep;<br />for (x0 = -nScale, x1 = 0, x2 = nScale; x1 < nWidth; x0++, x1++, x2++)<br />{<br />nGray = pLine[1][x1];<br />// 0 度<br />if (x2 < nWidth)<br />{<br />m_pLine[0][nGray][pLine[1][x2]]++;<br />}<br />// 45 度<br />if (y0 > 0 && x2 < nWidth)<br />{<br />m_pLine[1][nGray][pLine[0][x2]]++;<br />}<br />// 90 度<br />if (y0 > 0)<br />{<br />m_pLine[2][nGray][pLine[0][x1]]++;<br />}<br />// 135 度<br />if (y0 > 0 && x0 > 0)<br />{<br />m_pLine[3][nGray][pLine[0][x0]]++;<br />}<br />// 180 度<br />if (x0 > 0)<br />{<br />m_pLine[0][nGray][pLine[1][x0]]++;<br />}<br />// 225 度<br />if (x0 > 0 && y2 < nHeight)<br />{<br />m_pLine[1][nGray][pLine[2][x0]]++;<br />}<br />// 270 度<br />if (y2 < nHeight)<br />{<br />m_pLine[2][nGray][pLine[2][x1]]++;<br />}<br />// 315 度<br />if (y2 < nHeight && x2 < nWidth)<br />{<br />m_pLine[3][nGray][pLine[2][x2]]++;<br />}<br />}<br />pLine[1] += nWidthStep;<br />}<br />bResult = true;<br />}<br />}<br />return bResult;<br />}<br />
提取共生矩陣特徵
為了能更直觀地以共生矩陣描述紋理狀況,從共生矩陣匯出一些反映矩陣狀況的參數,典型的有以下幾種:
- 角
二階矩(ASM):是灰階共生矩陣元素值的平方和,所以也稱能量。它反映了映像灰階分布均勻程度和紋理粗細度。如果共生矩陣的所有值均相等,則ASM值
小;相反,如果其中一些值大而其它值小,則ASM值大。當共生矩陣中元素集中分布時,此時ASM值大。ASM值大表明一種較均一和規則變化的紋理模式。
- 熵(ENT):是映像所具有的資訊量的度量,紋理資訊也屬於映像的資訊,是一個隨機性的度量,當共生矩陣中所有元素有最大的隨機性、空間共生矩陣中所有值幾乎相等時,共生矩陣中元素分散分布時,熵較大。它表示了映像中紋理的非均勻程度或複雜程度。
- 對比(CON):反映了映像的清晰度和紋理溝紋深淺的程度。紋理溝紋越深,其對比越大,視覺效果越清晰;反之,對比小,則溝紋淺,效果模糊。灰階差即對比大的象素對越多,這個值越大。灰階公生矩陣中遠離對角線的元素值越大,CON越大。
- 逆差距(HOM):反映映像紋理的同質性,度量映像紋理局部變化的多少。其值大則說明映像紋理的不同地區間缺少變化,局部非常均勻。
- 相
關(COR):度量空間灰階共生矩陣元素在行或列方向上的相似程度,因此,相關值大小反映了映像中局部灰階相關性。當矩陣元素值均勻相等時,相關值就大;
相反,如果矩陣像元值相差很大則相關值小。如果映像中有水平方向紋理,則水平方向矩陣的COR大於其餘矩陣的COR值。
#include <math.h><br />// 共生矩陣特徵<br />struct GlcmFeature<br />{<br />double ASM[4]; // 角二階矩/能量<br />double ENT[4]; // 熵<br />double CON[4]; // 對比<br />double HOM[4]; // 逆差矩/同質性<br />double COR[4]; // 相關性<br />};<br />// 灰階共生矩陣<br />class CGlcm<br />{<br />public:<br />CGlcm(void);<br />virtual ~CGlcm(void);<br />public:<br />int *m_pMat1;<br />int *m_pMat2;<br />int *m_pMat3;<br />int *m_pMat4;<br />protected:<br />int *m_pLine[4][256];<br />int m_nMin;<br />int m_nMax;<br />int m_nSum;<br />public:<br />// 計算共生矩陣<br />// 參數:<br />// 1. pImageData: 映像資料指標,單通道,8位。<br />// 2. nLeft,nTop,nWidth,nHeight: 計算共生矩陣的地區。<br />// 3. nWidthStep: 行位移量。<br />// 4. nScale: 尺度。<br />// 5. nReduction: 灰階級壓縮。<br />bool CalGlcm(unsigned char *pImageData, int nLeft, int nTop, int nWidth, int nHeight, int nWidthStep, int nScale, int nReduction);<br />// 共生矩陣特徵<br />// 參數:<br />// 1. Feature: 輸出共生矩陣特徵。<br />void GetFeature(GlcmFeature& Feature);<br />};<br />// 共生矩陣特徵<br />void CGlcm::GetFeature(GlcmFeature& Feature)<br />{<br />int x, y;<br />int nTheta;<br />int nValue;<br />int nTemp;<br />double dValue;<br />double dMean, dStdDev;<br />double dSum[256];<br />// 清空記憶體<br />memset(&Feature, 0, sizeof(Feature));<br />// 方向迴圈<br />for (nTheta = 0; nTheta < 4; nTheta++)<br />{<br />dMean = 0;<br />dStdDev = 0;<br />// 清空記憶體<br />memset(dSum, 0, sizeof(dSum));<br />for (y = m_nMin; y <= m_nMax; y++)<br />{<br />for (x = m_nMin; x <= m_nMax; x++)<br />{<br />nValue = m_pLine[nTheta][y][x];<br />if (nValue != 0)<br />{<br />// 歸一化共生矩陣<br />dValue = (double) nValue / (double) m_nSum;<br />nTemp = (x - y) * (x - y);<br />// 角二階矩/能量<br />Feature.ASM[nTheta] += (dValue * dValue);<br />// 熵<br />Feature.ENT[nTheta] -= (dValue * log(dValue));<br />// 對比<br />Feature.CON[nTheta] += (nTemp * dValue);<br />// 逆差矩/同質性<br />Feature.HOM[nTheta] += (dValue / (1 + nTemp));<br />// 相關性<br />Feature.COR[nTheta] += (x * y * dValue);<br />dSum[y] += dValue;<br />}<br />}<br />}<br />for (y = m_nMin; y <= m_nMax; y++)<br />{<br />dMean += (y * dSum[y]);<br />}<br />for (y = m_nMin; y <= m_nMax; y++)<br />{<br />dStdDev += ((y - dMean) * (y - dMean) * dSum[y]);<br />}<br />// 相關性<br />if (abs(dStdDev) > 1e-15)<br />{<br />Feature.COR[nTheta] = (Feature.COR[nTheta] - dMean * dMean) / dStdDev;<br />}<br />else<br />{<br />Feature.COR[nTheta] = 0;<br />}<br />}<br />}