how to enter GLCM and pairs

function [out] = GLCM_Features1(glcmin,pairs) if ((nargin > 2) (nargin == 0)) error('Too many or too few input arguments. Enter GLCM and pairs.'); elseif ( (nargin == 2) ) if ((size(glcmin,1) <= 1) (size(glcmin,2) <= 1)) error('The GLCM should be a 2-D or 3-D matrix.'); elseif ( size(glcmin,1) ~= size(glcmin,2) ) error('Each GLCM should be square with NumLevels rows and NumLevels cols'); end elseif (nargin == 1) % only GLCM is entered pairs = 0; % default is numbers and input 1 for percentage if ((size(glcmin,1) <= 1) (size(glcmin,2) <= 1)) error('The GLCM should be a 2-D or 3-D matrix.'); elseif ( size(glcmin,1) ~= size(glcmin,2) ) error('Each GLCM should be square with NumLevels rows and NumLevels cols'); end end
format long e if (pairs == 1) newn = 1; for nglcm = 1:2:size(glcmin,3) glcm(:,:,newn) = glcmin(:,:,nglcm) + glcmin(:,:,nglcm+1); newn = newn + 1; end elseif (pairs == 0) glcm = glcmin; end
size_glcm_1 = size(glcm,1); size_glcm_2 = size(glcm,2); size_glcm_3 = size(glcm,3); out.autoc = zeros(1,size_glcm_3); % Autocorrelation: [2] out.contr = zeros(1,size_glcm_3); % Contrast: matlab/[1,2] out.corrm = zeros(1,size_glcm_3); % Correlation: matlab out.corrp = zeros(1,size_glcm_3); % Correlation: [1,2] out.cprom = zeros(1,size_glcm_3); % Cluster Prominence: [2] out.cshad = zeros(1,size_glcm_3); % Cluster Shade: [2] out.dissi = zeros(1,size_glcm_3); % Dissimilarity: [2] out.energ = zeros(1,size_glcm_3); % Energy: matlab / [1,2] out.entro = zeros(1,size_glcm_3); % Entropy: [2] out.homom = zeros(1,size_glcm_3); % Homogeneity: matlab out.homop = zeros(1,size_glcm_3); % Homogeneity: [2] out.maxpr = zeros(1,size_glcm_3); % Maximum probability: [2] out.sosvh = zeros(1,size_glcm_3); % Sum of sqaures: Variance [1] out.savgh = zeros(1,size_glcm_3); % Sum average [1] out.svarh = zeros(1,size_glcm_3); % Sum variance [1] out.senth = zeros(1,size_glcm_3); % Sum entropy [1] out.dvarh = zeros(1,size_glcm_3); % Difference variance [4] out.denth = zeros(1,size_glcm_3); % Difference entropy [1] out.inf1h = zeros(1,size_glcm_3); % Information measure of correlation1 [1] out.inf2h = zeros(1,size_glcm_3); % Informaiton measure of correlation2 [1] out.indnc = zeros(1,size_glcm_3); % Inverse difference normalized (INN) [3] out.idmnc = zeros(1,size_glcm_3); % Inverse difference moment normalized [3] glcm_sum = zeros(size_glcm_3,1); glcm_mean = zeros(size_glcm_3,1); glcm_var = zeros(size_glcm_3,1); u_x = zeros(size_glcm_3,1); u_y = zeros(size_glcm_3,1); s_x = zeros(size_glcm_3,1); s_y = zeros(size_glcm_3,1); p_x = zeros(size_glcm_1,size_glcm_3); % Ng x #glcms[1] p_y = zeros(size_glcm_2,size_glcm_3); % Ng x #glcms[1] p_xplusy = zeros((size_glcm_1*2 - 1),size_glcm_3); %[1] p_xminusy = zeros((size_glcm_1),size_glcm_3); %[1] hxy = zeros(size_glcm_3,1); hxy1 = zeros(size_glcm_3,1); hx = zeros(size_glcm_3,1); hy = zeros(size_glcm_3,1); hxy2 = zeros(size_glcm_3,1); for k = 1:size_glcm_3 % number glcms glcm_sum(k) = sum(sum(glcm(:,:,k))); glcm(:,:,k) = glcm(:,:,k)./glcm_sum(k); % Normalize each glcm glcm_mean(k) = mean2(glcm(:,:,k)); % compute mean after norm glcm_var(k) = (std2(glcm(:,:,k)))^2; for i = 1:size_glcm_1
for j = 1:size_glcm_2
out.contr(k) = out.contr(k) + (abs(i - j))^2.*glcm(i,j,k);
out.dissi(k) = out.dissi(k) + (abs(i - j)*glcm(i,j,k));
out.energ(k) = out.energ(k) + (glcm(i,j,k).^2);
out.entro(k) = out.entro(k) - (glcm(i,j,k)*log(glcm(i,j,k) + eps));
out.homom(k) = out.homom(k) + (glcm(i,j,k)/( 1 + abs(i-j) ));
out.homop(k) = out.homop(k) + (glcm(i,j,k)/( 1 + (i - j)^2));
out.sosvh(k) = out.sosvh(k) + glcm(i,j,k)*((i - glcm_mean(k))^2);
out.indnc(k) = out.indnc(k) + (glcm(i,j,k)/( 1 + (abs(i-j)/size_glcm_1) ));
out.idmnc(k) = out.idmnc(k) + (glcm(i,j,k)/( 1 + ((i - j)/size_glcm_1)^2));
u_x(k) = u_x(k) + (i)*glcm(i,j,k); % changed 10/26/08
u_y(k) = u_y(k) + (j)*glcm(i,j,k); % changed 10/26/08
end
end
out.maxpr(k) = max(max(glcm(:,:,k)));
end
for k = 1:size_glcm_3
for i = 1:size_glcm_1
for j = 1:size_glcm_2
p_x(i,k) = p_x(i,k) + glcm(i,j,k);
p_y(i,k) = p_y(i,k) + glcm(j,i,k); % taking i for j and j for i
if (ismember((i + j),[2:2*size_glcm_1]))
p_xplusy((i+j)-1,k) = p_xplusy((i+j)-1,k) + glcm(i,j,k);
end
if (ismember(abs(i-j),[0:(size_glcm_1-1)]))
p_xminusy((abs(i-j))+1,k) = p_xminusy((abs(i-j))+1,k) +...
glcm(i,j,k);
end
end
end
end
for k = 1:(size_glcm_3)
for i = 1:(2*(size_glcm_1)-1)
out.savgh(k) = out.savgh(k) + (i+1)*p_xplusy(i,k);
out.senth(k) = out.senth(k) - (p_xplusy(i,k)*log(p_xplusy(i,k) + eps));
end
end
for k = 1:(size_glcm_3)
for i = 1:(2*(size_glcm_1)-1)
out.svarh(k) = out.svarh(k) + (((i+1) - out.senth(k))^2)*p_xplusy(i,k);
end
end for k = 1:size_glcm_3 for i = 0:(size_glcm_1-1) out.denth(k) = out.denth(k) - (p_xminusy(i+1,k)*log(p_xminusy(i+1,k) + eps)); out.dvarh(k) = out.dvarh(k) + (i^2)*p_xminusy(i+1,k); end end for k = 1:size_glcm_3 hxy(k) = out.entro(k); for i = 1:size_glcm_1
for j = 1:size_glcm_2
hxy1(k) = hxy1(k) - (glcm(i,j,k)*log(p_x(i,k)*p_y(j,k) + eps));
hxy2(k) = hxy2(k) - (p_x(i,k)*p_y(j,k)*log(p_x(i,k)*p_y(j,k) + eps));
end
hx(k) = hx(k) - (p_x(i,k)*log(p_x(i,k) + eps));
hy(k) = hy(k) - (p_y(i,k)*log(p_y(i,k) + eps));
end
out.inf1h(k) = ( hxy(k) - hxy1(k) ) / ( max([hx(k),hy(k)]) );
out.inf2h(k) = ( 1 - exp( -2*( hxy2(k) - hxy(k) ) ) )^0.5;
end
corm = zeros(size_glcm_3,1); corp = zeros(size_glcm_3,1); for k = 1:size_glcm_3 for i = 1:size_glcm_1 for j = 1:size_glcm_2 s_x(k) = s_x(k) + (((i) - u_x(k))^2)*glcm(i,j,k); s_y(k) = s_y(k) + (((j) - u_y(k))^2)*glcm(i,j,k); corp(k) = corp(k) + ((i)*(j)*glcm(i,j,k)); corm(k) = corm(k) + (((i) - u_x(k))*((j) - u_y(k))*glcm(i,j,k)); out.cprom(k) = out.cprom(k) + (((i + j - u_x(k) - u_y(k))^4)*... glcm(i,j,k)); out.cshad(k) = out.cshad(k) + (((i + j - u_x(k) - u_y(k))^3)*... glcm(i,j,k)); end end
s_x(k) = s_x(k) ^ 0.5;
s_y(k) = s_y(k) ^ 0.5;
out.autoc(k) = corp(k);
out.corrp(k) = (corp(k) - u_x(k)*u_y(k))/(s_x(k)*s_y(k));
out.corrm(k) = corm(k) / (s_x(k)*s_y(k));
end

1 Comment

Yep, that's a program. What about it?
(For those who don't know GLCM is an image processing concept called Gray Level Co-ocurrence Matrix and is calculated by a MATLAB function called graycomatrix() - like the FFT, it's normally not a matrix that you'd be expected to type in by hand.)

Sign in to comment.

Answers (0)

Categories

Asked:

on 20 Jul 2012

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!