灰度共生矩阵生成
灰度共生矩阵生成分析及C++实现
发布于 2021610|更新于 202192|遵循 CC BY-NC-SA 4.0 许可

算法目的

通过创建灰度共生矩阵GLCM来表示纹理的空间分布情况。

算法原理

  1. 选择不同方向(0° , 45°, 90°, 135°)

  2. 对原图灰度值进行量化,如8级、16级等

  3. 遍历所有像元,统计特定步长下不同灰度级出现的次数。

    灰度共生矩阵
  4. 输出获得的灰度共生矩阵。

分析过程

在MFC工程中添加GDAL库,并在程序中引用并注册。

cpp
复制代码
1#include "gdal/include/gdal_priv.h" 2#include "gdal/include/gdalwarper.h" 3#include "gdal/include/ogrsf_frmts.h" 4#pragma comment(lib, "gdal/lib/gdal.lib")

初始化对话框时,注册驱动器并禁用UTF8编码。

cpp
复制代码
1GDALAllRegister(); 2CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");

获取图像的灰度信息,这里创建了一维数组pData用于存储数据。

cpp
复制代码
1GDALDriver* poDriver=GetGDALDriverManager()->GetDriverByName("GTIFF"); 2GDALDataset* poDataset = (GDALDataset*)GDALOpen(inPath, GA_ReadOnly); 3GDALRasterBand* poBand = poDataset->GetRasterBand(band+1); 4int xSize = poDataset->GetRasterXSize(); 5int ySize = poDataset->GetRasterYSize(); 6int* pData = new int[xSize * ySize]; 7poBand->RasterIO(GF_Read, 0, 0, xSize, ySize, pData, xSize, ySize, GDT_Int32, 0, 0);

对原图进行量化,灰度级数以gray表示;同时根据灰度级数创建最终输出的灰度共生矩阵。

cpp
复制代码
1int max = pData[0]; 2for (int i = 0; i < xSize * ySize; i++) { 3if (pData[i] > max) max = pData[i]; 4} 5int* pDataClass = new int[xSize * ySize]; 6for (int i = 0; i < xSize * ySize; i++) { 7 pDataClass[i] = (int)(pData[i] / (double)(max) * gray); 8 if (max == pData[i]) pDataClass[i] = gray - 1; 9} 10 11int* glcm = new int[gray * gray]; 12for (int i = 0; i < gray * gray; i++) glcm[i] = 0;

根据不同的方向,进行灰度共生矩阵的统计。

cpp
复制代码
1int p1 = 0, p2 = 0; 2if (strcmp("0°", ori) == 0) { 3 for (int i = 0; i < ySize; i++) { 4 for (int j = 0; j < xSize - step; j++) { 5 p1 = pDataClass[i * xSize + j]; 6 p2 = pDataClass[i * xSize + j + step]; 7 glcm[p1 * gray + p2]++; 8 } 9 } 10} 11else if (strcmp("45°", ori) == 0) { 12 for (int i = step; i < ySize; i++) { 13 for (int j = 0; j < xSize - step; j++) { 14 p1 = pDataClass[i * xSize + j]; 15 p2 = pDataClass[(i - step) * xSize + j + step]; 16 glcm[p1 * gray + p2]++; 17 } 18 } 19} 20else if (strcmp("90°", ori) == 0) { 21 for (int i = step; i < ySize; i++) { 22 for (int j = 0; j < xSize; j++) { 23 p1 = pDataClass[i * xSize + j]; 24 p2 = pDataClass[(i - step) * xSize + j]; 25 glcm[p1 * gray + p2]++; 26 } 27 } 28} 29else if (strcmp("135°", ori) == 0) { 30 for (int i = step; i < ySize; i++) { 31 for (int j = step; j < xSize; j++) { 32 p1 = pDataClass[i * xSize + j]; 33 p2 = pDataClass[(i - step) * xSize + j - step]; 34 glcm[p1 * gray + p2]++; 35 } 36 } 37}

处理完毕,写进文本文档里输出,完成任务。

cpp
复制代码
1std::ofstream outfile; 2outfile.open(outPath); 3for (int i = 0; i < gray; i++) { 4 for (int j = 0; j < gray; j++) { 5 outfile << glcm[i * gray + j] << '\t'; 6 } 7 outfile << '\n'; 8}
Comments