灰度共生矩阵(GLCM)并计算能量、熵、惯性矩、相关性(matlab)(待总结)

灰度共生矩阵(GLCM)并计算能量、熵、惯性矩、相关性(matlab)(待总结)

关于灰度共生矩阵的介绍可参考

http://blog.csdn.net/chuminnan2010/article/details/22035751

http://blog.csdn.net/xuezhisd/article/details/8908824

http://blog.csdn.net/xuexiang0704/article/details/8713204

http://cn.mathworks.com/help/images/ref/imlincomb.html?refresh=true

理论介绍

http://blog.csdn.net/lskyne/article/details/8659225 
http://blog.csdn.net/light_lj/article/details/26098815 
http://blog.csdn.net/kezunhai/article/details/42001477

共生矩阵的物理意义 
http://blog.csdn.net/light_lj/article/details/26098815

下面给出不同的NumLevels的例子

  1. gray = [1 1 5 6 8
  2. 2 3 5 7 1
  3. 4 5 7 1 2
  4. 8 5 1 2 5];
  5. GLCM = graycomatrix(gray,'GrayLimits',[])
  6. GLCM =
  7. 1 2 0 0 1 0 0 0
  8. 0 0 1 0 1 0 0 0
  9. 0 0 0 0 1 0 0 0
  10. 0 0 0 0 1 0 0 0
  11. 1 0 0 0 0 1 2 0
  12. 0 0 0 0 0 0 0 1
  13. 2 0 0 0 0 0 0 0
  14. 0 0 0 0 1 0 0 0
  15. GLCM = graycomatrix(gray, 'GrayLimits',[],'offset', [0 1])
  16. GLCM =
  17. 1 2 0 0 1 0 0 0
  18. 0 0 1 0 1 0 0 0
  19. 0 0 0 0 1 0 0 0
  20. 0 0 0 0 1 0 0 0
  21. 1 0 0 0 0 1 2 0
  22. 0 0 0 0 0 0 0 1
  23. 2 0 0 0 0 0 0 0
  24. 0 0 0 0 1 0 0 0
  25. GLCM = graycomatrix(gray, 'GrayLimits',[],'offset', [0 1],'NumLevels',8)
  26. GLCM =
  27. 1 2 0 0 1 0 0 0
  28. 0 0 1 0 1 0 0 0
  29. 0 0 0 0 1 0 0 0
  30. 0 0 0 0 1 0 0 0
  31. 1 0 0 0 0 1 2 0
  32. 0 0 0 0 0 0 0 1
  33. 2 0 0 0 0 0 0 0
  34. 0 0 0 0 1 0 0 0
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43

灰度共生矩阵(GLCM)并计算能量、熵、惯性矩、相关性(matlab)(待总结)

上图显示了如何求解灰度共生矩阵,以(1,1)点为例,GLCM(1,1)值为1说明只有一对灰度为1的像素水平相邻。GLCM(1,2)值为2,是因为有两对灰度为1和2的像素水平相邻。(1,5)出现一次,所以在(1.5)位置上标记1,没出现(1,6)所以为0;

上面所有的参数都是默认设置,NumLevels=8, 下面考虑NumLevels=3 的情况

  1. gray = [1 1 5 6 8
  2. 2 3 5 7 1
  3. 4 5 7 1 2
  4. 8 5 1 2 5];
  5. GL(2) = max(max(gray));
  6. GL(1) = min(min(gray));
  7. if GL(2) == GL(1)
  8. SI = ones(size(gray));
  9. else
  10. slope = NumLevels/(GL(2) - GL(1));
  11. intercept = 1 - (slope*(GL(1)));
  12. SI = floor(imlincomb(slope,gray,intercept,'double'));
  13. end
  14. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  15. SI =
  16. 1 1 2 3 4
  17. 1 1 2 3 1
  18. 2 2 3 1 1
  19. 4 2 1 1 2
  20. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  21. SI(SI > NumLevels) = NumLevels;
  22. SI(SI < 1) = 1;
  23. SI =
  24. 1 1 2 3 3
  25. 1 1 2 3 1
  26. 2 2 3 1 1
  27. 3 2 1 1 2
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33

上面给出了如何将初始矩阵gray变成3阶的灰度级,SI就是gray的3阶灰度级矩阵。

灰度共生矩阵(GLCM)并计算能量、熵、惯性矩、相关性(matlab)(待总结)

上图显示了如何求解3级灰度共生矩阵,以(1,1)点为例,GLCM(1,1)值为4说明只有4对灰度为1的像素水平相邻。GLCM(3,1)值为2,是因为有两对灰度为3和1的像素水平相邻。(2,1)出现一次,所以在(2.1)位置上标记1,没出现(1,3)所以为0;

  1. gray = [1 1 5 6 8
  2. 2 3 5 7 1
  3. 4 5 7 1 2
  4. 8 5 1 2 5];
  5. [GLCM,SI] = graycomatrix(gray,'NumLevels',3,'G',[])
  6. GLCM =
  7. 4 3 0
  8. 1 1 3
  9. 2 1 1
  10. SI =
  11. 1 1 2 3 3
  12. 1 1 2 3 1
  13. 2 2 3 1 1
  14. 3 2 1 1 2
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19

用3阶灰度级去计算四个共生矩阵P,取距离为1,角度分别为0,45,90,135

  1. gray = [1 1 5 6 8
  2. 2 3 5 7 1
  3. 4 5 7 1 2
  4. 8 5 1 2 5];
  5. offsets = [0 1;-1 1;-1 0;-1 -1];
  6. m = 3; % 3阶灰度级
  7. [GLCMS,SI] = graycomatrix(gray,'GrayLimits',[],'Of',offsets,'NumLevels',m);
  8. P = GLCMS;
  9. [kk,ll,mm] = size(P);
  10. % 对共生矩阵归一化
  11. %---------------------------------------------------------
  12. for n = 1:mm
  13. P(:,:,n) = P(:,:,n)/sum(sum(P(:,:,n)));
  14. end
  15. %-----------------------------------------------------------
  16. %对共生矩阵计算能量、熵、惯性矩、相关4个纹理参数
  17. %-----------------------------------------------------------
  18. H = zeros(1,mm);
  19. I = H;
  20. Ux = H; Uy = H;
  21. deltaX= H; deltaY = H;
  22. C =H;
  23. for n = 1:mm
  24. E(n) = sum(sum(P(:,:,n).^2)); %能量
  25. for i = 1:kk
  26. for j = 1:ll
  27. if P(i,j,n)~=0
  28. H(n) = -P(i,j,n)*log(P(i,j,n))+H(n); %熵
  29. end
  30. I(n) = (i-j)^2*P(i,j,n)+I(n); %惯性矩
  31. Ux(n) = i*P(i,j,n)+Ux(n); %相关性中μx
  32. Uy(n) = j*P(i,j,n)+Uy(n); %相关性中μy
  33. end
  34. end
  35. end
  36. for n = 1:mm
  37. for i = 1:kk
  38. for j = 1:ll
  39. deltaX(n) = (i-Ux(n))^2*P(i,j,n)+deltaX(n); %相关性中σx
  40. deltaY(n) = (j-Uy(n))^2*P(i,j,n)+deltaY(n); %相关性中σy
  41. C(n) = i*j*P(i,j,n)+C(n);
  42. end
  43. end
  44. C(n) = (C(n)-Ux(n)*Uy(n))/deltaX(n)/deltaY(n); %相关性
  45. end
  46. %--------------------------------------------------------------------------
  47. %求能量、熵、惯性矩、相关的均值和标准差作为最终8维纹理特征
  48. %--------------------------------------------------------------------------
  49. a1 = mean(E)
  50. b1 = sqrt(cov(E))
  51. a2 = mean(H)
  52. b2 = sqrt(cov(H))
  53. a3 = mean(I)
  54. b3 = sqrt(cov(I))
  55. a4 = mean(C)
  56. b4 = sqrt(cov(C))
  57. sprintf('0,45,90,135方向上的能量依次为: %f, %f, %f, %f',E(1),E(2),E(3),E(4)) % 输出数据;
  58. sprintf('0,45,90,135方向上的熵依次为: %f, %f, %f, %f',H(1),H(2),H(3),H(4)) % 输出数据;
  59. sprintf('0,45,90,135方向上的惯性矩依次为: %f, %f, %f, %f',I(1),I(2),I(3),I(4)) % 输出数据;
  60. sprintf('0,45,90,135方向上的相关性依次为: %f, %f, %f, %f',C(1),C(2),C(3),C(4)) % 输出数据;
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43
  • 44
  • 45
  • 46
  • 47
  • 48
  • 49
  • 50
  • 51
  • 52
  • 53
  • 54
  • 55
  • 56
  • 57
  • 58
  • 59
  • 60
  • 61
  • 62
  • 63
  • 64
  • 65
  • 66
  • 67
  • 68
  • 69
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43
  • 44
  • 45
  • 46
  • 47
  • 48
  • 49
  • 50
  • 51
  • 52
  • 53
  • 54
  • 55
  • 56
  • 57
  • 58
  • 59
  • 60
  • 61
  • 62
  • 63
  • 64
  • 65
  • 66
  • 67
  • 68
  • 69

下面给出的是默认的设置下求四个方向的灰度共生矩阵,NumLevels=8.

  1. gray = [1 1 5 6 8
  2. 2 3 5 7 1
  3. 4 5 7 1 2
  4. 8 5 1 2 5];
  5. offsets = [0 1;-1 1;-1 0;-1 -1];
  6. [GLCMS,SI] = graycomatrix(gray,'GrayLimits',[],'Of',offsets)
  7. GLCMS(:,:,1) =
  8. 1 2 0 0 1 0 0 0
  9. 0 0 1 0 1 0 0 0
  10. 0 0 0 0 1 0 0 0
  11. 0 0 0 0 1 0 0 0
  12. 1 0 0 0 0 1 2 0
  13. 0 0 0 0 0 0 0 1
  14. 2 0 0 0 0 0 0 0
  15. 0 0 0 0 1 0 0 0
  16. GLCMS(:,:,2) =
  17. 2 0 0 0 0 0 0 0
  18. 1 1 0 0 0 0 0 0
  19. 0 0 0 0 1 0 0 0
  20. 0 0 1 0 0 0 0 0
  21. 0 0 0 0 1 1 1 0
  22. 0 0 0 0 0 0 0 0
  23. 0 0 0 0 0 0 1 1
  24. 0 0 0 0 1 0 0 0
  25. GLCMS(:,:,3) =
  26. 0 0 0 0 0 0 2 1
  27. 3 0 0 0 0 0 0 0
  28. 1 0 0 0 0 0 0 0
  29. 0 1 0 0 0 0 0 0
  30. 0 1 1 0 2 0 0 0
  31. 0 0 0 0 0 0 0 0
  32. 0 0 0 0 1 1 0 0
  33. 0 0 0 1 0 0 0 0
  34. GLCMS(:,:,4) =
  35. 0 0 0 0 2 1 0 0
  36. 0 0 0 0 0 0 2 0
  37. 1 0 0 0 0 0 0 0
  38. 0 0 0 0 0 0 0 0
  39. 2 1 0 1 0 0 0 0
  40. 0 0 0 0 0 0 0 0
  41. 0 0 1 0 1 0 0 0
  42. 0 0 0 0 0 0 0 0
  43. SI =
  44. 1 1 5 6 8
  45. 2 3 5 7 1
  46. 4 5 7 1 2
  47. 8 5 1 2 5
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43
  • 44
  • 45
  • 46
  • 47
  • 48
  • 49
  • 50
  • 51
  • 52
  • 53
  • 54
  • 55
  • 56
  • 57
  • 58
  • 59
  • 60
  • 61
  • 62
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 19
  • 20
  • 21
  • 22
  • 23
  • 24
  • 25
  • 26
  • 27
  • 28
  • 29
  • 30
  • 31
  • 32
  • 33
  • 34
  • 35
  • 36
  • 37
  • 38
  • 39
  • 40
  • 41
  • 42
  • 43
  • 44
  • 45
  • 46
  • 47
  • 48
  • 49
  • 50
  • 51
  • 52
  • 53
  • 54
  • 55
  • 56
  • 57
  • 58
  • 59
  • 60
  • 61
  • 62

下面给出更多的关于灰度共生矩阵的特征

  1. I = imread('circuit.tif');
  2. GLCM2 = graycomatrix(I,'GrayLimits',[],'Offset',[0 1;-1 1;-1 0;-1 -1]);
  3. stats = GLCM_Features1(GLCM2,0)
  • 1
  • 2
  • 3
  • 1
  • 2
  • 3
  1. I = imread('circuit.tif');
  2. AllGLCMFeatureb= caluateglcmfeature(I);
  • 1
  • 2
  • 1
  • 2
  1. function [newglcm] = caluateglcmfeature(I);
  2. GLCM2 = graycomatrix(I,'GrayLimits',[],'Offset',[0 1;-1 1;-1 0;-1 -1]);
  3. stats = GLCM_Features1(GLCM2,0)
  4. newstats = struct2cell(stats);
  5. [statsx,statsy] = size(newstats);
  6. newglcm = [];
  7. for i = 1:statsx
  8. element = [];
  9. newelement = [];
  10. glcm = [];
  11. element = newstats(i,1);
  12. newelement = element{1,1};
  13. average = mean(newelement);
  14. variance = sqrt(cov(newelement));
  15. glcm = [newelement average variance];
  16. newglcm = [newglcm glcm];
  17. end
  18. end
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  • 1
  • 2
  • 3
  • 4
  • 5
  • 6
  • 7
  • 8
  • 9
  • 10
  • 11
  • 12
  • 13
  • 14
  • 15
  • 16
  • 17
  • 18
  1. % http://www.mathworks.com/matlabcentral/fileexchange/22187-glcm-texture-features
  2. % This code from the upper website but we have changed some.
  3. function [out] = GLCM_Features1(glcmin,pairs)
  4. % GLCM_Features1 helps to calculate the features from the different GLCMs
  5. % that are input to the function. The GLCMs are stored in a i x j x n
  6. % matrix, where n is the number of GLCMs calculated usually due to the
  7. % different orientation and displacements used in the algorithm. Usually
  8. % the values i and j are equal to 'NumLevels' parameter of the GLCM
  9. % computing function graycomatrix(). Note that matlab quantization values
  10. % belong to the set {1,..., NumLevels} and not from {0,...,(NumLevels-1)}
  11. % as provided in some references
  12. % http://www.mathworks.com/access/helpdesk/help/toolbox/images/graycomatrix
  13. % .html
  14. %
  15. % Although there is a function graycoprops() in Matlab Image Processing
  16. % Toolbox that computes four parameters Contrast, Correlation, Energy,
  17. % and Homogeneity. The paper by Haralick suggests a few more parameters
  18. % that are also computed here. The code is not fully vectorized and hence
  19. % is not an efficient implementation but it is easy to add new features
  20. % based on the GLCM using this code. Takes care of 3 dimensional glcms
  21. % (multiple glcms in a single 3D array)
  22. %
  23. % If you find that the values obtained are different from what you expect
  24. % or if you think there is a different formula that needs to be used
  25. % from the ones used in this code please let me know.
  26. % A few questions which I have are listed in the link
  27. % http://www.mathworks.com/matlabcentral/newsreader/view_thread/239608
  28. %
  29. % I plan to submit a vectorized version of the code later and provide
  30. % updates based on replies to the above link and this initial code.
  31. %
  32. % Features computed
  33. % Autocorrelation: [2] (out.autoc)
  34. % Contrast: matlab/[1,2] (out.contr)
  35. % Correlation: matlab (out.corrm)
  36. % Correlation: [1,2] (out.corrp)
  37. % Cluster Prominence: [2] (out.cprom)
  38. % Cluster Shade: [2] (out.cshad)
  39. % Dissimilarity: [2] (out.dissi)
  40. % Energy: matlab / [1,2] (out.energ)
  41. % Entropy: [2] (out.entro)
  42. % Homogeneity: matlab (out.homom)
  43. % Homogeneity: [2] (out.homop)
  44. % Maximum probability: [2] (out.maxpr)
  45. % Sum of sqaures: Variance [1] (out.sosvh)
  46. % Sum average [1] (out.savgh)
  47. % Sum variance [1] (out.svarh)
  48. % Sum entropy [1] (out.senth)
  49. % Difference variance [1] (out.dvarh)
  50. % Difference entropy [1] (out.denth)
  51. % Information measure of correlation1 [1] (out.inf1h)
  52. % Informaiton measure of correlation2 [1] (out.inf2h)
  53. % Inverse difference (INV) is homom [3] (out.homom)
  54. % Inverse difference normalized (INN) [3] (out.indnc)
  55. % Inverse difference moment normalized [3] (out.idmnc)
  56. %
  57. % The maximal correlation coefficient was not calculated due to
  58. % computational instability
  59. % http://murphylab.web.cmu.edu/publications/boland/boland_node26.html
  60. %
  61. % Formulae from MATLAB site (some look different from
  62. % the paper by Haralick but are equivalent and give same results)
  63. % Example formulae:
  64. % Contrast = sum_i(sum_j( (i-j)^2 * p(i,j) ) ) (same in matlab/paper)
  65. % Correlation = sum_i( sum_j( (i - u_i)(j - u_j)p(i,j)/(s_i.s_j) ) ) (m)
  66. % Correlation = sum_i( sum_j( ((ij)p(i,j) - u_x.u_y) / (s_x.s_y) ) ) (p[2])
  67. % Energy = sum_i( sum_j( p(i,j)^2 ) ) (same in matlab/paper)
  68. % Homogeneity = sum_i( sum_j( p(i,j) / (1 + |i-j|) ) ) (as in matlab)
  69. % Homogeneity = sum_i( sum_j( p(i,j) / (1 + (i-j)^2) ) ) (as in paper)
  70. %
  71. % Where:
  72. % u_i = u_x = sum_i( sum_j( i.p(i,j) ) ) (in paper [2])
  73. % u_j = u_y = sum_i( sum_j( j.p(i,j) ) ) (in paper [2])
  74. % s_i = s_x = sum_i( sum_j( (i - u_x)^2.p(i,j) ) ) (in paper [2])
  75. % s_j = s_y = sum_i( sum_j( (j - u_y)^2.p(i,j) ) ) (in paper [2])
  76. %
  77. %
  78. % Normalize the glcm:
  79. % Compute the sum of all the values in each glcm in the array and divide
  80. % each element by it sum
  81. %
  82. % Haralick uses 'Symmetric' = true in computing the glcm
  83. % There is no Symmetric flag in the Matlab version I use hence
  84. % I add the diagonally opposite pairs to obtain the Haralick glcm
  85. % Here it is assumed that the diagonally opposite orientations are paired
  86. % one after the other in the matrix
  87. % If the above assumption is true with respect to the input glcm then
  88. % setting the flag 'pairs' to 1 will compute the final glcms that would result
  89. % by setting 'Symmetric' to true. If your glcm is computed using the
  90. % Matlab version with 'Symmetric' flag you can set the flag 'pairs' to 0
  91. %
  92. % References:
  93. % 1. R. M. Haralick, K. Shanmugam, and I. Dinstein, Textural Features of
  94. % Image Classification, IEEE Transactions on Systems, Man and Cybernetics,
  95. % vol. SMC-3, no. 6, Nov. 1973
  96. % 2. L. Soh and C. Tsatsoulis, Texture Analysis of SAR Sea Ice Imagery
  97. % Using Gray Level Co-Occurrence Matrices, IEEE Transactions on Geoscience
  98. % and Remote Sensing, vol. 37, no. 2, March 1999.
  99. % 3. D A. Clausi, An analysis of co-occurrence texture statistics as a
  100. % function of grey level quantization, Can. J. Remote Sensing, vol. 28, no.
  101. % 1, pp. 45-62, 2002
  102. % 4. http://murphylab.web.cmu.edu/publications/boland/boland_node26.html
  103. %
  104. %
  105. % Example:
  106. %
  107. % Usage is similar to graycoprops() but needs extra parameter 'pairs' apart
  108. % from the GLCM as input
  109. % I = imread('circuit.tif');
  110. % GLCM2 = graycomatrix(I,'GrayLimits',[],'Offset',[0 1;-1 1;-1 0;-1 -1]);
  111. % stats = GLCM_Features1(GLCM2,0)
  112. % The output is a structure containing all the parameters for the different
  113. % GLCMs
  114. %
  115. % [Avinash Uppuluri: [email protected]: Last modified: 11/20/08]
  116. % If 'pairs' not entered: set pairs to 0
  117. if ((nargin > 2) || (nargin == 0))
  118. error('Too many or too few input arguments. Enter GLCM and pairs.');
  119. elseif ( (nargin == 2) )
  120. if ((size(glcmin,1) <= 1) || (size(glcmin,2) <= 1))
  121. error('The GLCM should be a 2-D or 3-D matrix.');
  122. elseif ( size(glcmin,1) ~= size(glcmin,2) )
  123. error('Each GLCM should be square with NumLevels rows and NumLevels cols');
  124. end
  125. elseif (nargin == 1) % only GLCM is entered
  126. pairs = 0; % default is numbers and input 1 for percentage
  127. if ((size(glcmin,1) <= 1) || (size(glcmin,2) <= 1))
  128. error('The GLCM should be a 2-D or 3-D matrix.');
  129. elseif ( size(glcmin,1) ~= size(glcmin,2) )
  130. error('Each GLCM should be square with NumLevels rows and NumLevels cols');
  131. end
  132. end
  133. format long e
  134. if (pairs == 1)
  135. newn = 1;
  136. for nglcm = 1:2:size(glcmin,3)
  137. glcm(:,:,newn) = glcmin(:,:,nglcm) + glcmin(:,:,nglcm+1);
  138. newn = newn + 1;
  139. end
  140. elseif (pairs == 0)
  141. glcm = glcmin;
  142. end
  143. size_glcm_1 = size(glcm,1);
  144. size_glcm_2 = size(glcm,2);
  145. size_glcm_3 = size(glcm,3);
  146. % checked
  147. out.autoc = zeros(1,size_glcm_3); % Autocorrelation: [2]
  148. out.contr = zeros(1,size_glcm_3); % Contrast: matlab/[1,2]
  149. out.corrm = zeros(1,size_glcm_3); % Correlation: matlab
  150. out.corrp = zeros(1,size_glcm_3); % Correlation: [1,2]
  151. out.cprom = zeros(1,size_glcm_3); % Cluster Prominence: [2]
  152. out.cshad = zeros(1,size_glcm_3); % Cluster Shade: [2]
  153. out.dissi = zeros(1,size_glcm_3); % Dissimilarity: [2]
  154. out.energ = zeros(1,size_glcm_3); % Energy: matlab / [1,2]
  155. out.entro = zeros(1,size_glcm_3); % Entropy: [2]
  156. out.homom = zeros(1,size_glcm_3); % Homogeneity: matlab
  157. out.homop = zeros(1,size_glcm_3); % Homogeneity: [2]
  158. out.maxpr = zeros(1,size_glcm_3); % Maximum probability: [2]
  159. out.sosvh = zeros(1,size_glcm_3); % Sum of sqaures: Variance [1]
  160. out.savgh = zeros(1,size_glcm_3); % Sum average [1]
  161. out.svarh = zeros(1,size_glcm_3); % Sum variance [1]
  162. out.senth = zeros(1,size_glcm_3); % Sum entropy [1]
  163. out.dvarh = zeros(1,size_glcm_3); % Difference variance [4]
  164. %out.dvarh2 = zeros(1,size_glcm_3); % Difference variance [1]
  165. out.denth = zeros(1,size_glcm_3); % Difference entropy [1]
  166. out.inf1h = zeros(1,size_glcm_3); % Information measure of correlation1 [1]
  167. out.inf2h = zeros(1,size_glcm_3); % Informaiton measure of correlation2 [1]
  168. %out.mxcch = zeros(1,size_glcm_3);% maximal correlation coefficient [1]
  169. %out.invdc = zeros(1,size_glcm_3);% Inverse difference (INV) is homom [3]
  170. out.indnc = zeros(1,size_glcm_3); % Inverse difference normalized (INN) [3]
  171. out.idmnc = zeros(1,size_glcm_3); % Inverse difference moment normalized [3]
  172. % correlation with alternate definition of u and s
  173. %out.corrm2 = zeros(1,size_glcm_3); % Correlation: matlab
  174. %out.corrp2 = zeros(1,size_glcm_3); % Correlation: [1,2]
  175. glcm_sum = zeros(size_glcm_3,1);
  176. glcm_mean = zeros(size_glcm_3,1);
  177. glcm_var = zeros(size_glcm_3,1);
  178. % http://www.fp.ucalgary.ca/mhallbey/glcm_mean.htm confuses the range of
  179. % i and j used in calculating the means and standard deviations.
  180. % As of now I am not sure if the range of i and j should be [1:Ng] or
  181. % [0:Ng-1]. I am working on obtaining the values of mean and std that get
  182. % the values of correlation that are provided by matlab.
  183. u_x = zeros(size_glcm_3,1);
  184. u_y = zeros(size_glcm_3,1);
  185. s_x = zeros(size_glcm_3,1);
  186. s_y = zeros(size_glcm_3,1);
  187. % % alternate values of u and s
  188. % u_x2 = zeros(size_glcm_3,1);
  189. % u_y2 = zeros(size_glcm_3,1);
  190. % s_x2 = zeros(size_glcm_3,1);
  191. % s_y2 = zeros(size_glcm_3,1);
  192. % checked p_x p_y p_xplusy p_xminusy
  193. p_x = zeros(size_glcm_1,size_glcm_3); % Ng x #glcms[1]
  194. p_y = zeros(size_glcm_2,size_glcm_3); % Ng x #glcms[1]
  195. p_xplusy = zeros((size_glcm_1*2 - 1),size_glcm_3); %[1]
  196. p_xminusy = zeros((size_glcm_1),size_glcm_3); %[1]
  197. % checked hxy hxy1 hxy2 hx hy
  198. hxy = zeros(size_glcm_3,1);
  199. hxy1 = zeros(size_glcm_3,1);
  200. hx = zeros(size_glcm_3,1);
  201. hy = zeros(size_glcm_3,1);
  202. hxy2 = zeros(size_glcm_3,1);
  203. %Q = zeros(size(glcm));