怎么拟合一条闭合曲线

dongsheng01

比如最简单的 一个长方形 只知道每边上的有限个点 怎么把这个闭合曲线拟合出来？

messenger

这个问题以前讨论过，搜索一下以前的贴子

dongsheng01

2# messenger 我想使用插值得到更多的点
但是在使用那个插值函数时 总是出现错误提示 是怎么回事 附上代码 是一个很简单的长方形
The values of X should be distinct.

x=[0;1;2;2;2;1;0;0];
y=[0;0;0;1;2;2;2;1];
points=[x;y];
x0=mean(x);y0=mean(y); % the center point
r=sqrt((x-x0).^2+(y-y0).^2); % the radius
theta=acos((x-x0)./r); %the angle
index=(y-y0)<0;
theta(index)=2*pi-theta(index); %计算各个点的角度
[theta,i1]=sort(theta);
r1=r(i1);
theta1=min(theta):2*pi/100:max(theta);
r2=interp1(theta,[r1],theta1,'spline');  
xi=x0+r2.*cos(theta1);
yi=y0+r2.*sin(theta1);
points=[[xi;xi(1)],[yi;yi(1)]]';
curve=cscvn(points);
fnplt(curve);
hold on;
plot(points(1,:),points(2,:),'o')
hold off

messenger

points在数组合并时，数组的行列搞反了。

1. 
   
2. x=[0;1;2;2;2;1;0;0];
   
3. y=[0;0;0;1;2;2;2;1];
   
4. points=[x;y];
   
5. x0=mean(x);y0=mean(y); % the center point
   
6. r=sqrt((x-x0).^2+(y-y0).^2); % the radius
   
7. theta=acos((x-x0)./r); %the angle
   
8. index=(y-y0)<0;
   
9. theta(index)=2*pi-theta(index); %计算各个点的角度
   
10. [theta,i1]=sort(theta);
    
11. r1=r(i1);
    
12. theta1=min(theta):2*pi/100:max(theta);
    
13. r2=interp1(theta,[r1],theta1,'spline');  
    
14. xi=x0+r2.*cos(theta1);
    
15. yi=y0+r2.*sin(theta1);
    
16. %points=[[xi;xi(1)],[yi;yi(1)]]';
    
17. points=[[xi xi(1)];[yi yi(1)]];
    
18. curve=cscvn(points);
    
19. fnplt(curve);
    
20. hold on;
    
21. plot(points(1,:),points(2,:),'yo')
    
22. hold off
    

复制代码

dongsheng01

4# messenger
哦 谢谢 我给出的点 只想拟合出长方形 但是从拟合效果来看 似乎不是很好 是不是插值函数 不应该选择三次样条

starcool

学习一下啊

dongsheng01

4# messenger
怎么求plot拟合的封闭图形的面积 area 只是把这个面积显示出来

nwcwww

4# messenger
怎么求plot拟合的封闭图形的面积 area 只是把这个面积显示出来
dongsheng01 发表于 2011-5-4 08:33

area = 1/2*abs(sum(x.*y([2:end,1])-y.*x([2:end,1])));

dongsheng01

8# nwcwww
x=[-1;1;1;-1;0;0;1;-1];
y=[1; 1; -1; -1;1;-1;0;0];
points=[x;y];
[m,n]=size(y);
x0=mean(x);y0=mean(y); % the center point
r=sqrt((x-x0).^2+(y-y0).^2); % the radius
theta=acos((x-x0)./r); %the angle
index=(y-y0)<0;
theta(index)=2*pi-theta(index); %计算各个点的角度
[theta,i1]=sort(theta);
r1=r(i1);
i=0;s=0;
area = 1/2*abs(sum(x.*y([2:end,1])-y.*x([2:end,1])));
xi=x0+r1.*cos(theta);
yi=y0+r1.*sin(theta);
points=[[xi;xi(1)],[yi;yi(1)]]';
curve=cscvn(points);
fnplt(curve);
hold on;
plot(points(1,:),points(2,:),'o')
area(points(1,:),points(2,:))
hold off
精度不够啊 这是一个长方形 面积应该是4 算出来的是3.5

nwcwww

9# dongsheng01


点要顺时针或者逆时针排列的。
你之前那个例子就没问题，但现在点的顺序比较乱，直接plot(x,y)能看到很多交叉线。
还是这几个点，顺序改变一下就行
>> x=[-1 -1 -1 0 1 1 1 0];y=[-1 0 1 1 1 0 -1 -1];
>> S=1/2*abs(sum(x.*y([2:end,1])-y.*x([2:end,1])))

S =

     4


补充说明下，如果顺序已经排好了你也可以直接用polyarea算面积：
>> x=[-1 -1 -1 0 1 1 1 0];y=[-1 0 1 1 1 0 -1 -1];S=polyarea(x,y)

S =

     4

如果像你9楼那样顺序是乱的，考虑用convhull:
x=[-1;1;1;-1;0;0;1;-1];y=[1; 1; -1; -1;1;-1;0;0];[~,S]=convhull(x,y)

S =

     4

dongsheng01

10# nwcwww 谢谢解答  已经使用polyarea 函数解决

lujianfeng8484

长见识了，得认真学习

sxgywjj2009

学习一下，真是高手如云啊！

dztaobao

恩 是高手啊 我正找呢写了