%%%%%%%%%%%%%%%%%%%%%%%%
%由已知分形参数模拟分形轮廓
L=0.0005;
l=0.000000204;
N=floor(L/l);
d=1.597;
g=2.86e-10;
gama=1.5;
nl=floor(log(1/L)/log(gama));
nmax=floor(log(1/l)/log(gama));
x=0:l:L; %采样区间为[0,500],采样间隔为l=0.204 *注*在此,x的取值对轮廓有影响,对频谱也有影响
z=0;
for n=nl:nmax
g=2.86*10^(-10);gama=1.5;d=1.597;
h=g^(d-1)*(cos(2*pi*gama.^n*x)/gama^((2-d)*n));
z=h+z;
end
plot(x,z,'r')
%%%%%%%%%%%%%%%%%
%结构函数法求分形维数与分形粗糙度
N=length(z);
S=[];
for n=0:1:N-1
sum=0;
for i=1:1:N-n
sum=sum+(z(i+n)-z(i))^2;
end
S=[S sum/(N-n)];
end
k=2:1:150;
t=log10(k.*l);
S=log10(S(k));
plot(t,S,'bo')
hold on
p=polyfit(t,S,1);
Sa=polyval(p,t);
plot(t,Sa)
K=p(1);
B=p(2);
D=(4-K)/2;
C=gamma(2*D-3)*sin((D-1.5)/pi)/((4-2*D)*log(gama));
G=10^((B-log10(C))/(2*D-2));
%%%%%%%%%%%%%%%%%%%%%%%%
%结果为D =1.585004727033449 (D-1.597)/1.597=-0.007511128970915 C =0.438205817882719 G=5.300056019492621e-09
%k=2:1:200时,D =1.598924116746748 (D-1.597)/1.597=0.001204832026768 C =0.449502607658423 G =6.102802884236269e-09
%k=2:1:250时,D =1.605955456309375 (D-1.597)/1.597=0.005607674583203 C =0.455695985622856 G =6.554769875068262e-09
%k=2:1:300时,D =1.610228510638337 (D-1.597)/1.597=0.008283350431019 C =0.459627811261989 G =6.843935252885409e-09
%k=2:1:1000时,D =1.649746790392541 (D-1.597)/1.597=0.033028672756757 C =0.502744824303370 G =1.027739980890304e-08
评论2