纹理特征：灰度共生矩阵

(1).CalTexture.m文件

function T = CalTexture(Image)
[M,N,O] = size(Image);%% 得到图像的宽，高和通道
%--------------------------------------------------------------------------
%1.原始图像灰度化
%--------------------------------------------------------------------------
if O==3
    Gray = double(0.3*Image(:,:,1)+0.59*Image(:,:,2)+0.11*Image(:,:,3));
else
    if O==1
        Gray =Image;
    end;
end;
% imshow(Gray);
%--------------------------------------------------------------------------
%2.Î灰度共生矩阵分块，减少计算量
%--------------------------------------------------------------------------
Layer =8;
for i = 1:M
    for j = 1:N
        for n = 1:256/Layer
            if (n-1)*Layer <= Gray(i,j) && Gray(i,j)<=((n-1)*Layer +Layer-1)
                Gray(i,j) = n-1;
            end
        end
    end
end
%--------------------------------------------------------------------------
%3.计算四个方向的灰度矩阵0,45,90,135
%--------------------------------------------------------------------------
P = zeros(Layer,Layer,4);
dis = 5;
for m = 1:Layer
    for n = 1:Layer
        for i = 1:M
            for j = 1:N
                if j<N-dis&&Gray(i,j)==m-1&&Gray(i,j+dis)==n-1
                    P(m,n,1) = P(m,n,1)+1;
                    %P(n,m,1) = P(m,n,1)+1;
                    P(n,m,1) = P(n,m,1)+1;
                end
                if i<M-dis&&j<N-dis&&Gray(i,j)==m-1&&Gray(i+dis,j+dis)==n-1
                    P(m,n,2) = P(m,n,2)+1;
                    P(n,m,2) = P(n,m,2)+1;
                end
                if i<M-dis&&Gray(i,j)==m-1&&Gray(i+dis,j)==n-1
                    P(m,n,3) = P(m,n,3)+1;
                    P(n,m,3) = P(n,m,3)+1;
                end
                if i<M-dis&&j>dis&&Gray(i,j)==m-1&&Gray(i+dis,j-dis)==n-1
                    P(m,n,4) = P(m,n,4)+1;
                    P(n,m,4) = P(n,m,4)+1;
                end
            end
        end
        if m==n
            P(m,n,:) = P(m,n,:)*2;
        end
    end
end

%
% for n=1:4
%     imshow(p(n));
% end

%%---------------------------------------------------------
% ¶灰度矩阵归一化
%%---------------------------------------------------------
for n = 1:4
   P(:,:,n) = P(:,:,n)/sum(sum(P(:,:,n)));
end
%--------------------------------------------------------------------------
%4.¶计算灰度矩阵的特征
%--------------------------------------------------------------------------
H = zeros(1,4);
I = H;
Ux = H;      Uy = H;
deltaX= H;  deltaY = H;
C =H;
L =H;
for n = 1:4
    E(n) = sum(sum(P(:,:,n).^2)); %%计算能量
    for i = 1:Layer
        for j = 1:Layer
            if P(i,j,n)~=0
                H(n) = -P(i,j,n)*log(P(i,j,n))+H(n); %%计算熵
            end
            I(n) = (i-j)^2*P(i,j,n)+I(n);  %%计算中心矩
            L(n)  = P(i,j,n)/(1+(i-j).^2)+L(n);%计算惯性
            
            Ux(n) = i*P(i,j,n)+Ux(n); %计算X均值
            Uy(n) = j*P(i,j,n)+Uy(n); %计算Y均值
     
        end
    end
end
for n = 1:4
    for i = 1:Layer
        for j = 1:Layer
            deltaX(n) = (i-Ux(n))^2*P(i,j,n)+deltaX(n); %
            deltaY(n) = (j-Uy(n))^2*P(i,j,n)+deltaY(n); %
           % C(n) = i*j*P(i,j,n)+C(n);
           C(n) = (i-Ux(n))*(j-Uy(n))*P(i,j,n)+C(n);%计算相关性
        end
    end
    %C(n) = (C(n)-Ux(n)*Uy(n))/deltaX(n)/deltaY(n); %Ïà¹ØÐÔ
    C(n)= C(n)/deltaX(n)/deltaY(n);
end

%--------------------------------------------------------------------------
T(1,1) = mean(E);   T(2,1) = sqrt(cov(E));
T(1,2) = mean(H);   T(2,2) = sqrt(cov(H));
T(1,3) = mean(I);   T(2,3) = sqrt(cov(I));
T(1,4) = mean(C);   T(2,4) = sqrt(cov(C));

//输出结果到数组，有四组结果

for n =1:4;
    T(3,n)=C(n);
    T(4,n)=E(n);
    T(5,n)=H(n);
    T(6,n)=I(n);
    T(7,n)=L(n);
end;