| 查看: 1102 | 回复: 3 | ||
梦魇拾柒木虫 (著名写手)
|
[求助]
Matlab程序求助 已有1人参与
|
|
function multifractal(A) format long g L=length(A); i=1; modify=1; tmin=2; % 边框间距，“※” tmax=10; ttmin=-10; ttmax=10; % 自定义 q 的范围 for r=tmin:1:tmax c(i,1)=mod(L,r); i=i+1; end c'; % 计算不能被边长r整除的余数 a=L-c'; % 计算并剔除掉不能被边长r整除的原始数据 n=length(a); % 求解格网化边长的个数，即为 n TT=[]; j=1; r=tmin; % 自定义项，“※-2” for i=1:1:n % 即n=25-10+1，自定义的结果 B=A(1:a(i),1); U=reshape(B,r,length(B)/r); T=mean(U); T=T'.*r^3; TT(1: length(T),i)=[T]; modifying(modify,1)=length(T); r=r+1; if r>= tmin+n % 自定义项，“※-3” break; % 边长超过10+n，超过初始限制，则程序自动终止 end modify=modify+1; end modifying; TT= nthroot(TT,1); % 或者不缩小 % TT; TT=nonzeros(TT); for cugb=1:1:n modifying_modifying(cugb,1)=sum(modifying(1:cugb)); % 有影响的新加卷 % & * % ￥ # @ ！） …… end modifying_modifying; j=1; % q 为任意数，这里取1到n，为n，与 k取值保持一致，q过大，计算机无法 %识别，默认为无穷大，q过小，结果接近0，则意义不明确 for q=ttmin:ttmax %这里取 q=-10:1:10 for k=1:1:n X=TT(1:modifying_modifying(k,1),1).^q; if k>1 X=TT(modifying_modifying(k-1,1)+1:modifying_modifying(k,1),1).^q; end t=sum(X); XX(k,j)=[t]; % 这里用到两个循环，即考虑到了幂函数，又需考虑求和 end j=j+1; end XX; % 得到质量分配函数，Xq(ξ)， X=log(tmin:1:tmax); % X=log(tmax:-1:tmin); % 此系以前的自定义输入结果，“※—4” Y=log(XX); % figure(1) % plot(X',Y,'o-k') % 至此，计算多重分形谱的第一步，分配函数构建完毕 side_length= tmin:1:tmax; % 自定义网格边长，“※—5” side_length=side_length'; q=ttmin:ttmax; % q=-5:1:n-10，q=-10:1:10 m=1; [ha,hb]=size(XX); for i=1:1:hb % XX=XX’; % or not s=XX(:,i); % XX b=polyfit(log(side_length),log(s),1); % 在对数尺度下计算斜率 slope(m,1)=b(1,1); m=m+1; end slope; % 这里的Slope即为质量指数，τ(q) N=polyfit(q', slope,1); plot(q', slope) % 此步是考察τ(q)-q 之间的关系， a=diff(slope)./diff(q'); % 第三步计算，diff函数求偏导确实少一列 q=q'; f_a=a.*q(1:end-1)-slope(1:end-1); % figure(2) % a=sort(a,'ascend'); plot(a,f_a,'o-k') xlabel('α','FontSize',12); ylabel('f(α)','FontSize',12); % polyfit(a,f_a,3) a=sort(a,'ascend'); % a+1 % f_a+1 运行后出现??? Attempted to access modifying_modifying(0,1); index must be a positive integer or logical. Error in ==> multifractal at 47 X=TT(1:modifying_modifying(k-1,1),1).^q;哪位童鞋能帮我看看嘛 |
回复此楼» 猜你喜欢
售SCI一区文章，我:8.O.551.O.5.4,科目全,可伽急
已经有6人回复
-
售SCI文章，我:8O.5.5.1O.54,科目全,可十急
已经有8人回复
-
售SCI文章，我:8O.5.5.1O.54,科目全,可十急
已经有5人回复
-
售SCI一区T0P文章，我:8O.55.1.O.5.4,科目齐全,可+急
已经有10人回复
-
售一区SCI文章T0P，我:8O.551.O54,科目全,可十急
已经有6人回复
-
售SCI一区T0P文章，我:8O.55.1.O.54,科目全,可伽急
已经有9人回复
-
售SCI一区T0P文章，我:8.O55.1.O.54,科目全,可十急
已经有8人回复
-
国社科又开始会评了，不知道这次命运如何
已经有14人回复
-
售SCI一区文章，我:8.O.55.1.O.54,科目齐全,可伽急
已经有8人回复
-
麻烦专家们看看评委们的意见（F口面上）
已经有13人回复
1107040177
新虫 (小有名气)
- 应助: 0 (幼儿园)
- 金币: 386.4
- 散金: 10
- 帖子: 110
- 在线: 7.5小时
- 虫号: 1888288
- 注册: 2012-07-11
- 专业: 食品科学基础
2楼2019-05-29 17:37:25
3楼2019-05-29 20:59:31
梦魇拾柒
木虫 (著名写手)
- 应助: 13 (小学生)
- 金币: 2240.9
- 散金: 6355
- 红花: 4
- 帖子: 1033
- 在线: 240.4小时
- 虫号: 987900
- 注册: 2010-04-02
- 性别: GG
- 专业: 地质学
送红花一朵 |
k=1:1:n X=TT(1:modifying_modifying(k,1),1).^q; if k>1 X=TT(modifying_modifying(k，1)+1:modifying_modifying(k,1),1).^q; 有点小错误。应为（k,1）,是一维数组。。兄台能帮我调试一下，A为一维列向量即可。。 发自小木虫Android客户端 |
4楼2019-05-29 23:32:26
30