01年數(shù)學(xué)建模A題程序

1、fork=1:length(xx)fork=1:length(xx)程序1：計(jì)算最大內(nèi)切圓半徑的粗略范圍。運(yùn)行前請先將程序中“F:01A圖片”改為圖片在你電腦中的實(shí)際儲存路徑。運(yùn)行時(shí)用鼠標(biāo)在切片圖中點(diǎn)出內(nèi)切圓的大致位置；運(yùn)行后顯示的r是半徑估計(jì)區(qū)間的下限，rl是上限。clear,clcima=imread(F:01A圖片);%讀入一張圖片figureimshow(ima)%顯示這張圖像j,i=ginput(l);%在坐標(biāo)上取圓心i=round(i);j=round(j);r=l;s=0;xxyy=yuan(j,i,r);s=s+ima(yy(k),xx(k);endwhiles=0r=r+1;

2、xxyy=yuan(j,i,r);s=0;fork=1:length(xx)s=s+ima(yy(k),xx(k);endendr=r-1%rminr1=sqrt(sum(sum(ones(size(ima)-sum(sum(ima)/pi)程序2：計(jì)算最大內(nèi)切圓半徑、圓心坐標(biāo)，作投影圖與三維圖。運(yùn)行前同樣要修改路徑。forj=3:3:510forj=3:3:510運(yùn)行時(shí)要調(diào)用函數(shù)，在命令窗口當(dāng)tt=100時(shí)程序運(yùn)行完畢，運(yùn)行時(shí)間一般不超過2分鐘。運(yùn)行結(jié)束后，工作空間中，aax,aay為利用計(jì)算出的數(shù)據(jù)擬合后的球心橫、縱坐標(biāo)。ar為從每張圖計(jì)算出的最大內(nèi)切圓半徑。clear,clcr=28;t

3、t=0;forii=0:99tt=tt+1ima二imread(F:01A圖片int2str(ttT).bmp);ima1=ima;ima1(1512,:)=0;ima1(:,1512)=0;fori=3:3:510ifsum(sum(ima(i-1:i+1,j-2j+2)+sum(sum(ima(i-2fork1=1:length(xx)fork1=1:length(xx)i+2,j-1:j+1)=12ima1(i-1:i+1,j-1:j+1)=0;endendend%eg=edge(ima,sobel);求出切片的邊緣點(diǎn)的位置%yeg,xeg=find(eg);邊緣點(diǎn)的橫縱坐標(biāo)yeg,xe

4、g=find(ima1);fork=1:length(xeg)xxyy=yuan(xeg(k),yeg(k),r);ifrem(k,5)=1fl=floor(r/sqrt(2);ima(yeg(k)-fl:yeg(k)+fl,xeg(k)-fl:xeg(k)+fl)=1;endima(yy(kl),xx(kl)=l;%剔除區(qū)域C以外的點(diǎn)endendyc,xc=find(1-ima);len2=500000;len22=0;fork=1:length(xc)len2=500000;fork1=1:length(xeg)t=(xc(k)-xeg(k1)*(xc(k)-xeg(k1)+(yc(k)-

5、yeg(k1)*(yc(k)-yeg(k1);iftlen22len22=len2;z=k;endendlen2=500000;len22=0;fork=:fork0=:len2=500000;fork1=1:length(xeg)xz=xc(z)+k0;yz=yc(z)+k;t=(xz-xeg(k1)*(xz-xeg(k1)+(yz-yeg(k1)*(yz-yeg(k1);iftlen22len22=len2;kx=k0;ky=k;endendendax(tt)=yc(z)+ky-256;ay(tt)=xc(z)+kx-256;ar(tt)=sqrt(len22);endz=0:99;A=p

6、olyfit(z,ax,7);B=polyfit(z,ay,7);aax=polyval(A,z);aay=polyval(B,z);figure(1)plot(aax,aay)gridonfigure(2)plot(aax,z)gridonfigure(3)plot(aay,z)gridonfigure(4)plot3(aax,aay,z,r)gridon函數(shù)functionxx,yy=yuan(x,y,r)k=1;forj=x-r+1:x-1fori=y-r+1:y-1if(x-j)*(x-j)+(y-i)*(y-i)r*rxx(k)=j;yy(k)=i;k=k+1;break;endendendfori=y-r+1:y-1forj=x-r+1:x-1if(x-j)*(x-j)+(y-i)*(y-i)r*rxx(k)=j;yy(k)=i;k=k+1;break;endendendxx(k:k+1)=x-r+1x;yy(k:k+