24小时热门版块排行榜    

查看: 1781  |  回复: 5
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

kobybry_cn

铜虫 (初入文坛)

[求助] 这个程序求多重分形维数对吗 已有2人参与

clc;
clear all;
close all;
%%
%indata=inputdlg({'dla'});
a = imread('1.jpg');%%%%此处的1.jpg就是一条直线
[rows, cols] = size(a);
figure;imshow(a);
npix = sum(sum(a));
%% calculates niL which is the number of pixels in the ith box of size L
% ideas from boxcount.m by F. Moisy have been borrowed here
width = rows;
p =log(width)/log(2);
max_boxes = power(rows,2)/power(2,2);
nL = double(zeros(max_boxes,p));
for g=(p-1):-1:0
            siz = 2^(p-g);
            sizm1 = siz - 1;
            index = log2(siz);
            count = 0;
            for i=1:sizwidth-siz+1)
                for j=1:sizwidth-siz+1)
                    count = count + 1;
                    sums = sum(sum(a(i:i+sizm1,j:j+sizm1)));
                    nL(count,index) = sums;
                end
            end
end
%%
qran = 1;
logl = zeros(p,1);

for l=1:p
    logl(l) = log(power(2,l));
end
%% normalized masses
pL = double(zeros(max_boxes,p));
for l=1:p
    nboxes = power(rows,2)/power(power(2,l),2);
    norm = sum(nL(1:nboxes,l));
    if(norm ~= npix)
        display('error');
    end
    for i=1:nboxes
        pL(i,l) = nL(i,l)/norm;
    end
end
%%
%falpha, alpha
for l=1:p
   
    count = 0;
    nboxes = power(rows,2)/power(power(2,l),2);
   
    for q =-qran:+0.1:qran      
        
        %denominator of muiql
        qsum = 0.0;
        for i=1:nboxes
            if(pL(i,l) ~= 0)
                qsum = qsum + power(pL(i,l),q);
            end
        end

        fqnum = 0.0;
        aqnum = 0.0;
        smuiqL = 0.0;
        for i=1:nboxes
            if(pL(i,l) ~= 0)
                  muiqL = power(pL(i,l),q)/qsum;
                  fqnum = fqnum + (muiqL * log(muiqL));
                  aqnum = aqnum + (muiqL * log(pL(i,l)));
                  smuiqL = smuiqL + muiqL;
            end
        end
        if(uint8(smuiqL)~=1)
            display('error');
        end
        
        count = count + 1;
        fql(l,count) = fqnum;
        aql(l,count) = aqnum;
        qval(count) = q;
    end

end
%%
% alpha_q
for i=1:count
     line = polyfit(logl,aql(:,i),1);
     aq(i) = line(1);
     yfit = polyval(line,logl);
     sse = sum(power(aql(:,i)-yfit,2));
     sst = sum(power(aql(:,i)-mean(aql(:,i)),2));
     ar2(i) = 1+(sse/sst);
end
% f_q
for i=1:count
     line = polyfit(logl,fql(:,i),1);
     fq(i) = line(1);
     yfit = polyval(line,logl);
     sse = sum(power(fql(:,i)-yfit,2));
     sst = sum(power(fql(:,i)-mean(fql(:,i)),2));
     fr2(i) = 1+(sse/sst);
end

figure;plot(aq,fq,'r');
xlabel('\alpha(q)','FontSize',14);
ylabel('f(q)','FontSize',14);

[ Last edited by kobybry_cn on 2014-1-13 at 22:20 ]
回复此楼
climbing ing ing……
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
6楼2015-11-09 12:33:21
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 6 个回答

戴钢盔的猪头

木虫 (知名作家)

【答案】应助回帖

兄弟,你这么发帖不要怪没人回你,DEBUG比直接写新代码还要让人头痛!你自己运行下,如果有错,把错误报告贴上来,别人才好帮你分析。

[ 发自手机版 http://muchong.com/3g ]
3楼2014-01-24 18:59:31
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

kobybry_cn

铜虫 (初入文坛)

呃,这个程序运行是正常的,但我觉得结果是不对的,因为结果不应该超过2……
climbing ing ing……
4楼2014-01-25 23:00:59
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

冰凝香

新虫 (初入文坛)

【答案】应助回帖

楼主,请问这个程序我为什么运行不了呢 sums = sum(sum(a(i:i+sizm1,j:j+sizm1)));此处超过了指数范围
5楼2014-03-10 09:44:43
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考博] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +4 k0dTPqJtl0jt 2026-08-14 4/200 2026-08-15 01:40 by 4wMiSEwB6436
[考研] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 7lpolszZVXgi 2026-08-14 5/250 2026-08-15 00:52 by 4wMiSEwB6436
[考研] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +3 HFw0lei2R37i 2026-08-14 7/350 2026-08-15 00:40 by 4wMiSEwB6436
[论文投稿] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +3 HFw0lei2R37i 2026-08-14 5/250 2026-08-15 00:16 by 4wMiSEwB6436
[基金申请] 是这周出结果还是下周出结果? +4 yuleib84 2026-08-11 4/200 2026-08-14 23:05 by lfy8008
[基金申请] 欢迎发来filecode的Mz6后的代码验证其规律 +23 医学老男孩 2026-08-13 49/2450 2026-08-14 20:02 by zhaifei
[基金申请] 哪位老哥知道今年的国自然具体哪一天放榜? +5 Ldrop2023 2026-08-13 5/250 2026-08-14 18:48 by ssxclkj
[基金申请] 好奇怪的filecode +6 布布和一二 2026-08-08 7/350 2026-08-14 13:38 by Equinoxhua
[论文投稿] 职称评审,求友友推荐见刊最快的期刊 +5 工厂打螺丝 2026-08-08 5/250 2026-08-14 11:03 by 玖戈弋
[基金申请] filecode +15 documentary 2026-08-10 17/850 2026-08-14 10:08 by kissu88
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +10 家与远方 2026-08-10 15/750 2026-08-14 02:08 by 绵羊哥哥
[文学芳草园] 阿姨 +4 汪汪锅 2026-08-09 4/200 2026-08-13 19:43 by arzu_hma
[基金申请] 静等基金结果 +8 gjjjzhong 2026-08-10 21/1050 2026-08-13 17:56 by 且听虎啸
[基金申请] 重要来源:本周末出结果 +10 瞬息宇宙 2026-08-12 10/500 2026-08-13 15:46 by likettle
[基金申请] 分享一下我之前已中青C的计划书的filecode +4 布布和一二 2026-08-11 5/250 2026-08-13 12:56 by cratir
[基金申请] 结合人工智能,周易传统文化,filecode打分制来了,3分以上希望很大。 +3 Tide man 2026-08-12 4/200 2026-08-13 08:35 by ZJTJZ
[基金申请] 帮忙看看fileCode +7 wwncly 2026-08-10 13/650 2026-08-11 19:36 by 冰心玉壶晴
[基金申请] 什么时候出结果,有咨询渠道??? +3 Tide man 2026-08-11 3/150 2026-08-11 17:54 by kudofaye
[基金申请] 国自然结果 +4 Vierhys 2026-08-10 8/400 2026-08-10 15:06 by Vierhys
[基金申请] 这样的filecode谁见过 +11 布布和一二 2026-08-08 22/1100 2026-08-10 11:10 by wmfsnow
信息提示
请填处理意见