本文共 6652 字，大约阅读时间需要 22 分钟。

事实上这篇文章并非原创，代码都是别人写的，可是在我的机子上有些地方不能run，我做了一丁点的修改，所以就把文章设置为原创了。另外，最初的参照博客已经看不到了，我看到的已经是别人转的，所以我有必要贴一下别人的链接以示尊重。

第一个链接：  该文章已被加密，打不开了幸亏有别人已经转载博主的文章，我的这篇文章就是从第二篇文章那里复制加修改得到的，代码整理放在CSDN博客上看起来更加美观，且亲测可直接运行。

第二个链接：

对应原文为：

下面分模块贴代码，代码是我在第二个链接的基础上稍作修改得到的，亲测可直接运行，没有任何问题。

第一个指标 Coarseness，粗糙度

%调用举例：%image=rgb2gray(imread('example.jpg'));%f=coarseness(image,5)function Fcrs = coarseness( graypic,kmax )%graphic为待处理的灰度图像，2^kmax为最大窗口[h,w]=size(graypic); %获取图片大小A=zeros(h,w,2^kmax); %平均灰度值矩阵A%计算有效可计算范围内每个点的2^k邻域内的平均灰度值for i=2^(kmax-1)+1:h-2^(kmax-1)    for j=2^(kmax-1)+1:w-2^(kmax-1)        for k=1:kmax            A(i,j,k)=mean2(graypic(i-2^(k-1):i+2^(k-1)-1,j-2^(k-1):j+2^(k-1)-1));        end    endend%对每个像素点计算在水平和垂直方向上不重叠窗口之间的Ak差for i=1+2^(kmax-1):h-2^(kmax-1)    for j=1+2^(kmax-1):w-2^(kmax-1)        for k=1:kmax            Eh(i,j,k)=abs(A(i+2^(k-1),j,k)-A(i-2^(k-1),j));            Ev(i,j,k)=abs(A(i,j+2^(k-1),k)-A(i,j-2^(k-1)));        end    endend%对每个像素点计算使E达到最大值的kfor i=2^(kmax-1)+1:h-2^(kmax-1)    for j=2^(kmax-1)+1:w-2^(kmax-1)        [maxEh,p]=max(Eh(i,j,:));        [maxEv,q]=max(Ev(i,j,:));        if maxEh>maxEv            maxkk=p;        else            maxkk=q;        end        Sbest(i,j)=2^maxkk; %每个像素点的最优窗口大小为2^maxkk    endend%所有Sbest的均值作为整幅图片的粗糙度Fcrs=mean2(Sbest);end

%调用举例：%注意这个函数因为涉及到方差，要求输入类型为double，因此我这里在源代码上做了适当的修改%image=rgb2gray(imread('example.jpg'));%f=contrast(image)function Fcon=contrast(graypic) %graypic为待处理的灰度图片graypic=double(graypic);%这一句我自己做了修改，否则原博文中的代码不能直接运行x=graypic(:); %二维向量一维化M4=mean((x-mean(x)).^4); %四阶矩delta2=var(x,1); %方差alfa4=M4/(delta2^2); %峰度delta=std(x,1); %标准差Fcon=delta/(alfa4^(1/4)); %对比度end

%调用举例：%image=rgb2gray(imread('example.jpg'));%[Fdir,sita]=directionality(image)%sita为各像素点的角度矩阵，在线性度中会用到，所以这里作为结果返回function [Fdir,sita]=directionality(graypic)[h w]=size(graypic);%两个方向的卷积矩阵GradientH=[-1 0 1;-1 0 1;-1 0 1];GradientV=[ 1 1 1;0 0 0;-1 -1 -1];%卷积，取有效结果矩阵MHconv=conv2(graypic,GradientH);MH=MHconv(3:h,3:w);MVconv=conv2(graypic,GradientV);MV=MVconv(3:h,3:w)%向量模MG=(abs(MH)+abs(MV))./2;%有效矩阵大小validH=h-2;validW=w-2%各像素点的方向for i=1:validH    for j=1:validW        sita(i,j)=atan(MV(i,j)/MH(i,j))+(pi/2);    endendn=16;t=12;Nsita=zeros(1,n);%构造方向的统计直方图for i=1:validH    for j=1:validW        for k=1:n            if sita(i,j)>=(2*(k-1)*pi/2/n) && sita(i,j)<((2*(k-1)+1)*pi/2/n) && MG(i,j)>=t                Nsita(k)=Nsita(k)+1;            end        end    endendfor k=1:n    HD(k)=Nsita(k)/sum(Nsita(:));end%假设每幅图片只有一个方向峰值，为计算方便简化了原著[maxvalue,FIp]=max(HD);Fdir=0;for k=1:n    Fdir=Fdir+(k-FIp)^2*HD(k);%公式与原著有改动endend

%调用举例：%image=rgb2gray(imread('example.jpg'));%Flin=linelikeness(image,sita,4) %sita为directionality.m返回的结果function Flin=linelikeness(graypic,sita,d) %d为共生矩阵计算时的像素间隔距离n=16;[h,w]=size(graypic);%构造方向共生矩阵PDd1=zeros(n,n);PDd2=zeros(n,n);PDd3=zeros(n,n);PDd4=zeros(n,n);PDd5=zeros(n,n);PDd6=zeros(n,n);PDd7=zeros(n,n);PDd8=zeros(n,n);for i=d+1:h-d-2    for j=d+1:w-d-2        for m1=1:n            for m2=1:n                %下方向                 if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i+d,j)>=(2*(m2-1)*pi/2/n) && sita(i+d,j)<((2*(m2-1)+1)*pi/2/n))                    PDd1(m1,m2)=PDd1(m1,m2)+1;                end                %上方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i-d,j)>=(2*(m2-1)*pi/2/n) && sita(i-d,j)<((2*(m2-1)+1)*pi/2/n))                    PDd2(m1,m2)=PDd2(m1,m2)+1;                end                %右方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i,j+d)>=(2*(m2-1)*pi/2/n) && sita(i,j+d)<((2*(m2-1)+1)*pi/2/n))                    PDd3(m1,m2)=PDd3(m1,m2)+1;                end                %左方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i,j-d)>=(2*(m2-1)*pi/2/n) && sita(i,j-d)<((2*(m2-1)+1)*pi/2/n))                    PDd4(m1,m2)=PDd4(m1,m2)+1;                end                %右下方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i+d,j+d)>=(2*(m2-1)*pi/2/n) && sita(i+d,j+d)<((2*(m2-1)+1)*pi/2/n))                    PDd5(m1,m2)=PDd5(m1,m2)+1;                end                %右上方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i-d,j+d)>=(2*(m2-1)*pi/2/n) && sita(i-d,j+d)<((2*(m2-1)+1)*pi/2/n))                    PDd6(m1,m2)=PDd6(m1,m2)+1;                end                %左下方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i+d,j-d)>=(2*(m2-1)*pi/2/n) && sita(i+d,j-d)<((2*(m2-1)+1)*pi/2/n))                    PDd7(m1,m2)=PDd7(m1,m2)+1;                end                %左上方向                if (sita(i,j)>=(2*(m1-1)*pi/2/n) && sita(i,j)<((2*(m1-1)+1)*pi/2/n)) && (sita(i-d,j-d)>=(2*(m2-1)*pi/2/n) && sita(i-d,j-d)<((2*(m2-1)+1)*pi/2/n))                    PDd8(m1,m2)=PDd8(m1,m2)+1;                end            end        end    endendf=zeros(1,8);g=zeros(1,8);for i=1:n    for j=1:n        f(1)=f(1)+PDd1(i,j)*cos((i-j)*2*pi/n);        g(1)=g(1)+PDd1(i,j);        f(2)=f(2)+PDd2(i,j)*cos((i-j)*2*pi/n);        g(2)=g(2)+PDd2(i,j);        f(3)=f(3)+PDd3(i,j)*cos((i-j)*2*pi/n);        g(3)=g(3)+PDd3(i,j);        f(4)=f(4)+PDd4(i,j)*cos((i-j)*2*pi/n);        g(4)=g(4)+PDd4(i,j);        f(5)=f(5)+PDd5(i,j)*cos((i-j)*2*pi/n);        g(5)=g(5)+PDd5(i,j);        f(6)=f(6)+PDd6(i,j)*cos((i-j)*2*pi/n);        g(6)=g(6)+PDd6(i,j);        f(7)=f(7)+PDd7(i,j)*cos((i-j)*2*pi/n);        g(7)=g(7)+PDd7(i,j);        f(8)=f(8)+PDd8(i,j)*cos((i-j)*2*pi/n);        g(8)=g(4)+PDd8(i,j);    endendtempM=f./g;Flin=max(tempM);%取8个方向的线性度最大值作为图片的线性度end

%调用举例：%image=rgb2gray(imread('example.jpg'));%Freg=regularity(image,64)function Freg=regularity(graypic,windowsize) %windowsize为计算规则度的子窗口大小[h,w]=size(graypic);k=0;for i=1:windowsize:h-windowsize    for j=1:windowsize:w-windowsize        k=k+1;        crs(k)=coarseness(graypic(i:i+windowsize-1,j:j+windowsize-1),5); %粗糙度        con(k)=contrast(graypic(i:i+windowsize-1,j:j+windowsize-1)); %对比度        [dire(k),sita]=directionality(graypic(i:i+windowsize-1,j:j+windowsize-1));%方向度        lin=linelikeness(graypic(i:i+windowsize-1,j:j+windowsize-1),sita,4)*10; %线性度，*10与crs、con、dire同量级化    endend%求上述各参数的标准差Dcrs=std(crs,1);Dcon=std(con,1);Ddir=std(dire,1);Dlin=std(lin,1);Freg=1-(Dcrs+Dcon+Ddir+Dlin)/4/100;end

第六个指标 Roughness，粗略度

%粗略度Frgh=Fcrs+Fcon;

有问题可发邮件至 jzwangATbjtuDOTeduDOTcn 讨论交流。注：代码非我所写，我只做丁点修改使之可在我的机器上顺利运行。向原作致谢，如有侵权马上删掉！

转载地址：https://jianzhuwang.blog.csdn.net/article/details/51584535 如侵犯您的版权，请留言回复原文章的地址，我们会给您删除此文章，给您带来不便请您谅解！