24小时热门版块排行榜    

查看: 614  |  回复: 2
本帖产生 1 个 仿真EPI ,点击这里进行查看

飞扬2282

荣誉版主 (著名写手)

[交流] 【求助】再次求助 已有1人参与

http://muchong.com/bbs/viewthread.php?tid=1952310

6楼的程序。

现在要求略有变化,数据存储在C:\TDdownload\data.dat 文件中,数据为512x512x
90,即90幅数据,请问如何读取?
应该在6楼的程序上修改即可?
回复此楼

» 猜你喜欢

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

anyuezhiji

银虫 (正式写手)

星空行者

★ ★ ★
飞扬2282(金币+30):高手啊,再次感谢,牛人! 2010-04-23 22:29
adu886886(金币+3, 仿真EPI+1):AN版造诣很好啊 2010-04-24 11:11
代码如下:

可以和http://muchong.com/bbs/viewthread.php?tid=1952310
6楼做比较
引用回帖:
function Find_Center_of_circles2
%%Find Center of circles from picture
%%PSL@CSU
%%QQ:547423688
%%Email:anyuezhiji@qq.com
%%Edit @ 2010.4.23
path='F:\Program Files\MATLAB\work\CurrentWork\test\';
filename='zhanli 512 512 90-float.dat';
outfile='CResult.txt';
clear CResult;
fid = fopen([path filename],'rb','ieee-le');
n=0;
maxcm=0;
flag=1;
for i=1:90
   fname=['Rusult',num2str(i)];
   disp(['正在处理第',num2str(i),'张图。。。']);
   n=n+1;
   img = fread(fid,[512, 512],'float32');
   %figure
   %image(img);
   %imshow(img, []);
       minp=max(max(img));
       mask=img>minp/3;
       mask = bwmorph(mask,'clean');
       mask = bwmorph(mask,'majority',inf);
   %figure;imshow(mask);
   CResult(n).fname=fname;
   CResult(n).R=[];
   CResult(n).C=[];
   [L,m]=bwlabel(mask);
   %m
   if flag
   for i=1:m
   mask1=mask;
   mask1(find(L~=i))=0;
   BW = edge(uint8(mask1),'sobel');
   %figure,imshow(BW);
   [r,c]=find(BW~=0);
   %eval(['BW',num2str(i),'=BW;']);
   %eval(['edge',num2str(i),'=[r c];']);
   CircleEdge=[r,c];
   [r,center]=findcenter(CircleEdge,2);
   CResult(n).R=[CResult(n).R;round(r)];
   CResult(n).C=[CResult(n).C;round(center([2 1]))];
   end
   if m>maxcm
   maxcm=m;
   end
   CResult(n).num=m;
   [tc,p]=sort(CResult(n).C(:,2));
   CResult(n).C=CResult(n).C(p,;
   fidout=fopen([path,outfile],'a+');
   fprintf(fidout,'%-8s',CResult(n).fname);
   cc=CResult(n).C;
   for j=1:m
   fprintf(fidout,'\t\t圆心%d:(%5d,%5d)',j,cc(j,1),cc(j,2));
   end
   fprintf(fidout,'\n');
   fclose(fidout);
   end%end if flag
end
status = fclose(fid);
if flag
fidout=fopen([path,outfile],'a+');
fprintf(fidout,'\n');
fprintf(fidout,'\n');
for j=1:maxcm
fprintf(fidout,'%%圆心%d坐标\n',j);
fprintf(fidout,'x%d=[',j);
for i=1:n
if j'%5d',CResult(i).C(j,1));
else
fprintf(fidout,'  nan');
end
end
fprintf(fidout,']'';\n');
fprintf(fidout,'y%d=[',j);
for i=1:n
if j'%5d',CResult(i).C(j,2));
else
fprintf(fidout,'  nan');
end
end
fprintf(fidout,']'';\n');
end
fprintf(fidout,'\n');
fprintf(fidout,'\n');
fclose(fidout);
end%end if flag
function [R,C]=findcenter(edge,mode)
%mode1 相加取中点
%mode2 最小二乘拟合圆x^2+y^2+a(1)*x+a(2)*y+a(3)=0
%mode3 用遗传算法,找到一点和一个半径R,各点到这点的距离与R之差的平方和最小
[m,n]= size(edge);
switch mode
    case 1
    C=sum(edge)/m;
    R=sum(sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2))/m;
    case 2
    A=[edge(:,1) edge(:,2) ones(m,1)];
    B=-[edge(:,1).*edge(:,1)+edge(:,2).*edge(:,2)];
    a=A\B;
    C=-.5*a([1 2])';
    R = sqrt((a(1)^2+a(2)^2)/4-a(3));
    case 3
    save('edge_m','edge','m');
    %[x, fval, reason]=ga(fitnessfcn,nvars,A,b,Aeq,beq,LB,UB,nonlcon,options)
    C=sum(edge)/m;
    R=sum(sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2))/m;
    options = gaoptimset('Generations',200,'InitialPopulation',[R,C],...
        'PopulationSize',50,'TolFun',1e-7,'TolCon',1e-7,'MutationFcn',@mutationadaptfeasible );
    [sol, fval, exitflag] = ga(@f,3,[],[],[],[],[0,min(edge)],[2*R,max(edge)],[],options);
    R=sol(1);
    C=sol([2,3]);
    fval
end
%适应度函数的matlab代码
function [eval]=f(sol)
load edge_m
R=sol(1);
C=sol([2,3]);
eval=sum((sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2)-R).^2)/m;

暗月下没有留下风的痕迹,但它已经寂然飘逝。。By<暗月之寂>:tiger38:
2楼2010-04-23 21:53:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

anyuezhiji

银虫 (正式写手)

星空行者

★ ★
小木虫(金币+0.5):给个红包,谢谢回帖交流
飞扬2282(金币+20):追加金币,感谢帮助 2010-04-24 07:56
adu886886(金币+1):尽心尽力。赞一个 2010-04-24 11:11
修改版 读数据 计算坐标
包含生成jpg gif avi 等附加功能
引用回帖:
function Find_Center_of_circles2
%%Find Center of circles from picture
%%PSL@CSU
%%QQ:547423688
%%Email:anyuezhiji@qq.com
%%Edit @ 2010.4.23
%%Modify @ 2010.4.24
clear all
path='F:\Program Files\MATLAB\work\CurrentWork\test\';
filename='zhanli 512 512 90-float.dat';
outfile='CResult.txt';
clear CResult;
fid = fopen([path filename],'rb','ieee-le');
n=0;
maxcm=0;
flag=1;%是否执行计算
savemov=0;%是否生成avi文件,1为是,0为否
savejpg=0;%是否生成jpg文件,1为是,0为否
savegif=1;%是否生成gif文件,1为是,0为否
if(savemov==1)
mov=avifile([path,'CResult.avi']);
end
%[header,count]=fread(fid,4100,'uint8');
for i=1:90
   fname=['Rusult',num2str(i)];
   disp(['正在处理第',num2str(i),'张图。。。']);
   n=n+1;
      img = fread(fid,[512, 512],'float32');
      img=img';
      img=img(512:-1:1,;
      %img=img(:,512:-1:1);
   fig=figure(1);
   %image(img);
   imshow(img, []);
   %movie
   pause(0.01);
   %drawnow;%刷新
   if(savegif==1)
   %f=getframe(fig);
   %imind=frame2im(f);
   %[imind,cm] = rgb2ind(imind,256);
   if n==1
       %imwrite(imind,cm,[path 'CResult.gif'],'gif', 'Loopcount',inf,'DelayTime',0.1);
       imwrite(im2uint8(img),[path 'CResult.gif'],'gif', 'Loopcount',inf,'DelayTime',0.1);
   else
       %imwrite(imind,cm,[path 'CResult.gif'],'gif','WriteMode','append','DelayTime',0.1);
       imwrite(im2uint8(img),[path 'CResult.gif'],'gif','WriteMode','append','DelayTime',0.1);
   end
   end
   if(savejpg==1)
   print(1,'-djpeg',[path 'Result_' fname '.jpg']);
   end
   if(savemov==1)
   f=getframe(fig);
   mov=addframe(mov,f);
   end
       minp=max(max(img));
       mask=img>minp/3;
       mask = bwmorph(mask,'clean');
       mask = bwmorph(mask,'majority',inf);
   %figure;imshow(mask);
   CResult(n).fname=fname;
   CResult(n).R=[];
   CResult(n).C=[];
   [L,m]=bwlabel(mask);
   %m
   if flag
   for i=1:m
   mask1=mask;
   mask1(find(L~=i))=0;
   BW = edge(uint8(mask1),'sobel');
   %figure,imshow(BW);
   [r,c]=find(BW~=0);
   %eval(['BW',num2str(i),'=BW;']);
   %eval(['edge',num2str(i),'=[r c];']);
   CircleEdge=[r,c];
   [r,center]=findcenter(CircleEdge,2);
   CResult(n).R=[CResult(n).R;round(r)];
   CResult(n).C=[CResult(n).C;round(center([2 1]))];
   end
   if m>maxcm
   maxcm=m;
   end
   CResult(n).num=m;
   [tc,p]=sort(CResult(n).C(:,2));
   CResult(n).C=CResult(n).C(p,;
   fidout=fopen([path,outfile],'a+');
   fprintf(fidout,'%-8s',CResult(n).fname);
   cc=CResult(n).C;
   for j=1:m
   fprintf(fidout,'\t\t圆心%d:(%5d,%5d)',j,cc(j,1),cc(j,2));
   end
   fprintf(fidout,'\n');
   fclose(fidout);
   end%end if flag
end
if(savemov==1)
mov=close(mov);
end
status = fclose(fid);
if flag
fidout=fopen([path,outfile],'a+');
fprintf(fidout,'\n');
fprintf(fidout,'\n');
for j=1:maxcm
fprintf(fidout,'%%圆心%d坐标\n',j);
fprintf(fidout,'x%d=[',j);
for i=1:n
if j'%5d',CResult(i).C(j,1));
else
fprintf(fidout,'  nan');
end
end
fprintf(fidout,']'';\n');
fprintf(fidout,'y%d=[',j);
for i=1:n
if j'%5d',CResult(i).C(j,2));
else
fprintf(fidout,'  nan');
end
end
fprintf(fidout,']'';\n');
end
fprintf(fidout,'\n');
fprintf(fidout,'\n');
fclose(fidout);
end%end if flag
function [R,C]=findcenter(edge,mode)
%mode1 相加取中点
%mode2 最小二乘拟合圆x^2+y^2+a(1)*x+a(2)*y+a(3)=0
%mode3 用遗传算法,找到一点和一个半径R,各点到这点的距离与R之差的平方和最小
[m,n]= size(edge);
switch mode
    case 1
    C=sum(edge)/m;
    R=sum(sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2))/m;
    case 2
    A=[edge(:,1) edge(:,2) ones(m,1)];
    B=-[edge(:,1).*edge(:,1)+edge(:,2).*edge(:,2)];
    a=A\B;
    C=-.5*a([1 2])';
    R = sqrt((a(1)^2+a(2)^2)/4-a(3));
    case 3
    save('edge_m','edge','m');
    %[x, fval, reason]=ga(fitnessfcn,nvars,A,b,Aeq,beq,LB,UB,nonlcon,options)
    C=sum(edge)/m;
    R=sum(sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2))/m;
    options = gaoptimset('Generations',200,'InitialPopulation',[R,C],...
        'PopulationSize',50,'TolFun',1e-7,'TolCon',1e-7,'MutationFcn',@mutationadaptfeasible );
    [sol, fval, exitflag] = ga(@f,3,[],[],[],[],[0,min(edge)],[2*R,max(edge)],[],options);
    R=sol(1);
    C=sol([2,3]);
    fval
end
%适应度函数的matlab代码
function [eval]=f(sol)
load edge_m
R=sol(1);
C=sol([2,3]);
eval=sum((sqrt((edge(:,1)-C(1)).^2+(edge(:,2)-C(2)).^2)-R).^2)/m;

运行结果
引用回帖:
正在处理第1张图。。。
正在处理第2张图。。。
正在处理第3张图。。。
。。。。。。

结果写入CResult.txt
引用回帖:
Rusult1                 圆心1:(  257,   99)                圆心2:(  257,  457)
Rusult1                 圆心1:(  257,   99)                圆心2:(  257,  457)
Rusult2                 圆心1:(  268,   99)                圆心2:(  243,  457)
。。。。。。

%圆心1坐标               
x1=[  257  268  279  290 。。。。。。]';
y1=[   99   99   99   99 。。。。。。]';
%圆心2坐标               
x2=[  257  243  229  215 。。。。。。]';
y2=[  457  457  457  457 。。。。。。]';

静态图
引用回帖:

。。。。。。

gif动画
引用回帖:

运动轨迹
引用回帖:

avi。。 略。。。。
暗月下没有留下风的痕迹,但它已经寂然飘逝。。By<暗月之寂>:tiger38:
3楼2010-04-24 01:48:55
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 飞扬2282 的主题更新
普通表情 高级回复 (可上传附件)
信息提示
请填处理意见