I want to fit curve on this polar plot around these points,As shown in the image,Can anyone help me please
14 次查看(过去 30 天)
显示 更早的评论
L K
2017-1-23
Please can anyone tell me how to do it?Any way will be fine
8 个评论
Walter Roberson
2017-1-23
编辑:Walter Roberson
2017-1-23
Is the group at the bottom negative r or large theta?
The ones near the origin, are those the same theta for a few points or is the scale just not large enough to show that they are different theta?
L K
2017-1-23
The bottom are the negative theta values
and at the origin they are different values, I have attached another image with a larger scale(refer img5.1) and img 5.2 is with larger scale for negative theta
L K
2017-1-28
编辑:L K
2017-1-28
Yes ,here it is...it is like at that particular radius the theta is plotted
theta=88,89,90,91,92,94,96,94,90,89,-100,-102,-104,-105,-104,-102,-101,-100
radius=5,7,11,17,26,39,46,44,32,3,0,18,34,32,33,29,28,20
i want an envelope kind of thing around those points both up(for positive values) and down(for negative).. (by envelope i only mean a curve looking like an envelope)
I am attaching a figure.. the dots represent some angles at some radius.... so i want that envelope around those angles by using curve fitting..
Image Analyst
2017-1-28
编辑:Image Analyst
2017-1-28
Do you mean like a pie-shaped sector containing the "bounding box" of all the points? Or do you mean the convex hull?
Walter Roberson
2017-1-28
Is it correct that you want to find an equation in theta that generates approximately radius, together with some bounding values creating a curved box near the projected line, such that the curved box has to contain some specified fraction of the original points? Sort of like if you had computed a line and drawn it on thickly and the thick line had to cover most of the points?
If so, then what portion of the original points need to be covered? 1 standard deviation for example?
Image Analyst
2017-1-29
What black dots? Do you mean the blue stars? And the envelope would be like if you lasooed the blue stars by hand tracing around them? Or used an active contour or alpha shape?
采纳的回答
Image Analyst
2017-1-28
Try this:
theta=[88,89,90,91,92,94,96,94,90,89,-100,-102,-104,-105,-104,-102,-101,-100];
radius=[5,7,11,17,26,39,46,44,32,3,0,18,34,32,33,29,28,20];
polarplot(theta, radius, 'b*');
hold on;
[x, y] = pol2cart(theta, radius);
k = convhull(x, y);
xch = x(k);
ych = y(k);
[thetaCH, rhoCH] = cart2pol(xch, ych);
polarplot(thetaCH, rhoCH, 'ro-');
27 个评论
L K
2017-1-29
Thank you very much !!!! Indeed a great work...!!! Helped me really :) i tried to figure out this lot many times, you made my job easier
L K
2017-1-29
can this be done automatically??? as in if i am given any set of data points...will it work in a similar way?
Image Analyst
2017-1-29
What's your definition of "automatically"? What about my code is not "automatic". There was no manual interaction needed like asking the user a question or asking him to draw anything. Please give the difference between automatic and manual.
L K
2017-1-30
编辑:Walter Roberson
2017-1-30
I was just trying it on some example data points...
here is the code:
VarName1 is just data points from 1 to 10
x=VarName1;
y = abs(sqrt(x));
[xx,yy] = pol2cart(x,y);
k = convhull(xx, yy);
xch = xx(k);
ych = yy(k);
[xCH, yCH] = cart2pol(xch, ych);
plot(xx(k), yy(k), 'ro-',x,y,'b*');
convex hull is like drawing lines on the edge points going outwards unlike concave hull... from what i understood
the output i got is :
L K
2017-1-30
OK I GOT MY FAULT, i had converted it from poltocart,it was not needed as they were not in degrees, my bad!!
Walter Roberson
2017-1-31
Requires R2014b or later.
For releases before that look in the File Exchange for "alpha shapes"
L K
2017-1-31
if i use k = boundary(x, y);
and keep the rest of the code same,I am not getting the output, it gives me an error saying
Error using alphaShape The input points must be 2D or 3D coordinates in numpoints-by-ndim format.
I think it is because i am plotting the points on polar plot ,and then drawing the boundary.
So what change in the code is necessary ?
Image Analyst
2017-1-31
Is your point set different than when I took them to create my answer? What is your current point set now?
L K
2017-1-31
What you showed me is convex hull...
I was trying to use the same data points,for concave hull,by using k = boundary(x, y)
but i am not getting concave hull...please tell me what changes are needed
Walter Roberson
2017-1-31
Please post your current code.
At the time of the call, are your x and y scalars or are they vectors or are they arrays?
L K
2017-2-1
编辑:Walter Roberson
2017-2-1
theta=[88,89,90,91,92,94,96,94,90,89,-100,-102,-104,-105,-104,-102,-101,-100];
radius=[5,7,11,17,26,39,46,44,32,3,0,18,34,32,33,29,28,20];
polar(theta, radius, 'b*');
hold on;
[x, y] = pol2cart(theta, radius);
k = boundary(x, y);
xch = x(k);
ych = y(k);
[thetaCH, rhoCH] = cart2pol(xch, ych);
polar(thetaCH, rhoCH, 'ro-');
here is the code
Walter Roberson
2017-2-1
Remember that your theta is not in radians.
theta=[88,89,90,91,92,94,96,94,90,89,-100,-102,-104,-105,-104,-102,-101,-100];
radius=[5,7,11,17,26,39,46,44,32,3,0,18,34,32,33,29,28,20];
theta_rad = theta * pi/180;
polar(theta_rad, radius, 'b*');
hold on;
[x, y] = pol2cart(theta_rad, radius);
k = boundary(x(:), y(:)); %needs column vectors
xch = x(k);
ych = y(k);
[thetaCH, rhoCH] = cart2pol(xch, ych);
polar(thetaCH, rhoCH, 'ro-');
L K
2017-2-2
Is there any way to give a tolerance say of about 10% to 30% for example around this plot ? so that i can say the points which lie within this tolerance is of so-n-so and which lie outside it is so-n so?
L K
2017-2-10
Hi Image analyst, Can you please help me in how to give a tolerance of 10 or 20% around this plot?
So that i can tell,values which lie in this tolerance belongs to some specific data and outside it to some other data?
Image Analyst
2017-2-10
How about you compute the centroid of the region, and then the radius of all the vertex points, and then just send them outwards 10 or 20%? That's what I'd do.
Walter Roberson
2017-2-10
The poster asked a new Question about this; I posted a solution there. It is a small nuisance because it is polar work.
L K
2017-3-22
The boundary function which i have used to fit a curve,follows what algorithm? or any formula? i needed it badly in order to understand how the boundary function actually goes about , please if any one of you know,,could you help me out?
Walter Roberson
2017-3-22
However, it is my opinion that the result you are coming up with has no useful meaning.
Walter Roberson
2017-3-22
The enclosing boundary that you are coming up with has no useful meaning that I have been able to come up with. I see no reason to use linear boundaries on a curved surface.
L K
2017-3-22
But is there any way in which i can join the points on the curved surface ?Basically i need to fit a curve around those points... so i thought this was quite closely satisfying the need... but how can i form a curved boundary if not linear boundaries?
Image Analyst
2017-3-22
Why does it matter? What are you going to do with a boundary, either curved or linear, if you had it? Maybe you think you need it for some reason but actually you don't. If you tell us what you're going to do with the information, maybe we can say whether that is the best approach or not.
L K
2017-3-23
HI Image analyst, i just have some data points..which i will plot on polar graph and fit boundary around those points...so the next time i run the same thing and if my data points tend to fall outside i will reject it or if within i accept it...this is just overview..... do u have some opinion on this ?Please feel free to share the same ...or if you have some better idea...
L K
2017-3-23
About mathematical boundaries i am not sure,basically what i know is i need to form an curve like around the points...like i mentioned above ..but after the curve is formed ,,is there a way to know its model ?
更多回答(1 个)
Walter Roberson
2017-1-23
If you have the curve fitting toolbox then use cftool. Import your theta and r and play with some models. Or use pol2cart first and import the resulting x and y and fit those.
9 个评论
Walter Roberson
2017-1-28
编辑:Walter Roberson
2017-1-29
My work so far suggests that a reasonable fit is:
R = -67506.3123723029 - 425.950909132369*Theta + 25140.4082290502*Theta^2 + 86472.2200845025 * sin(Theta + 1.56195303516666) + 6006.21747802037 * sin(Theta^2 + 4.93132363492489)
Here, Theta is theta * Pi / 180 -- that is, theta converted to radians.
The model I fitted against was
a0 + a1 * t + b1 * sin(t + c1) + a2 * t^2 + b2 * sin(t^2 + c2)
Here, without loss of generality, c1 and c2 can be restricted to 0 to 2*Pi .
My fitting was showing two distinct areas of solution, which was really slowing it down. After some thought I realized that because sin(theta+Pi) = -sin(theta), that if b1 and c1 are solutions then so are -b1 and c1+Pi. Therefore without loss of generality you can restrict b1 and b2 to be non-negative as well.
Walter Roberson
2017-1-29
Marginally better:
-67507.8055608757 - 425.954107825018*Theta + 25140.966254751*Theta^2 + 86474.0893419859*sin(Theta + 1.56195316282131) + 6006.33822370626*sin(Theta^2 + 4.93132643353749)
The range over which the residue values are "quite close" to the above are
-67508.4024493893 -425.955230769825 25140.7891946898 86473.496119004 1.56195312294972 6006.29986670352 4.93132555131309
to
-67507.3317487232 -425.953060425301 25141.1894099114 86474.8375860513 1.56195321779578 6006.38638459066 4.93132751208825
so there is some free play on a0, a2, b1, and a slight bit on b2. The free play is largest on b1 and a0.
The plot looks like this, where blue is the projected values and red is the actual values.
Walter Roberson
2017-1-29
Visually you do not do all that badly with
- 8.202817146914601 + 31.00545829304457*Theta + 299.819287460388*sin(Theta + 4.628895715827117)
again Theta = theta * pi/180
The other model has a better fit.
This is the model
a0 + a1 * Theta + b1 * sin(Theta + c1)
If you construct
syms a0 a1 b1 c1 f2(t)
f2(t) = a0 + a1*t + b1*sin(t + c1)
obj2 = sum((f2(theta*Pi/180) - radius).^2)
CCV2 = matlabFunction(obj2, 'Vars', {[a0, a1, b1 c1]})
fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi])
then you will probably get the above result. If not, then try the fmincon again; sometimes the rand() puts it into a worse location
L K
2017-1-29
yes i am getting the same prob, but when i tried it again the error went..but am still not getting the above figure... is there some statement to be put?
Walter Roberson
2017-1-30
syms a0 a1 b1 c1 f2(t) Pi
Pi = sym('pi');
f2(t) = a0 + a1*t + b1*sin(t + c1)
obj2 = sum((f2(theta*Pi/180) - radius).^2)
CCV2 = matlabFunction(obj2, 'Vars', {[a0, a1, b1 c1]})
options = struct('MaxFunEvals', 3000);
fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
You might need to repeat the fmincon a small number of times to get the solution I posted above.
L K
2017-1-31
I repeated the fmincon,but still the same error...i am not understanding why am i not getting that result..
Walter Roberson
2017-1-31
theta = [88,89,90,91,92,94,96,94,90,89,-100,-102,-104,-105,-104,-102,-101,-100];
radius = [5,7,11,17,26,39,46,44,32,3,0,18,34,32,33,29,28,20];
syms a0 a1 b1 c1 f2(t) Pi
Pi = sym('pi');
f2(t) = a0 + a1*t + b1*sin(t + c1)
obj2 = sum((f2(theta*Pi/180) - radius).^2)
CCV2 = matlabFunction(obj2, 'Vars', {[a0, a1, b1 c1]})
options = struct('MaxFunEvals', 3000);
fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
f2(t) =
a0 + a1*t + b1*sin(c1 + t)
obj2 =
(a0 + (pi*a1)/2 + b1*sin(c1 + pi/2) - 11)^2 + (a0 - (5*pi*a1)/9 + b1*sin(c1 - (5*pi)/9))^2 + (a0 + (pi*a1)/2 + b1*sin(c1 + pi/2) - 32)^2 + (a0 - (5*pi*a1)/9 + b1*sin(c1 - (5*pi)/9) - 20)^2 + (a0 - (7*pi*a1)/12 + b1*sin(c1 - (7*pi)/12) - 32)^2 + (a0 + (8*pi*a1)/15 + b1*sin(c1 + (8*pi)/15) - 46)^2 + (a0 - (17*pi*a1)/30 + b1*sin(c1 - (17*pi)/30) - 18)^2 + (a0 - (17*pi*a1)/30 + b1*sin(c1 - (17*pi)/30) - 29)^2 + (a0 + (22*pi*a1)/45 + b1*sin(c1 + (22*pi)/45) - 5)^2 + (a0 + (23*pi*a1)/45 + b1*sin(c1 + (23*pi)/45) - 26)^2 + (a0 - (26*pi*a1)/45 + b1*sin(c1 - (26*pi)/45) - 33)^2 + (a0 - (26*pi*a1)/45 + b1*sin(c1 - (26*pi)/45) - 34)^2 + (a0 + (47*pi*a1)/90 + b1*sin(c1 + (47*pi)/90) - 39)^2 + (a0 + (47*pi*a1)/90 + b1*sin(c1 + (47*pi)/90) - 44)^2 + (a0 + (89*pi*a1)/180 + b1*sin(c1 + (89*pi)/180) - 3)^2 + (a0 + (89*pi*a1)/180 + b1*sin(c1 + (89*pi)/180) - 7)^2 + (a0 + (91*pi*a1)/180 + b1*sin(c1 + (91*pi)/180) - 17)^2 + (a0 - (101*pi*a1)/180 + b1*sin(c1 - (101*pi)/180) - 28)^2
CCV2 =
function_handle with value:
@(in1)(in1(:,1)+in1(:,2).*pi.*(1.0./2.0)+in1(:,3).*sin(in1(:,4)+pi.*(1.0./2.0))-1.1e1).^2+(in1(:,1)-in1(:,2).*pi.*(5.0./9.0)+in1(:,3).*sin(in1(:,4)-pi.*(5.0./9.0))).^2+(in1(:,1)+in1(:,2).*pi.*(1.0./2.0)+in1(:,3).*sin(in1(:,4)+pi.*(1.0./2.0))-3.2e1).^2+(in1(:,1)-in1(:,2).*pi.*(5.0./9.0)+in1(:,3).*sin(in1(:,4)-pi.*(5.0./9.0))-2.0e1).^2+(in1(:,1)-in1(:,2).*pi.*(7.0./1.2e1)+in1(:,3).*sin(in1(:,4)-pi.*(7.0./1.2e1))-3.2e1).^2+(in1(:,1)+in1(:,2).*pi.*(8.0./1.5e1)+in1(:,3).*sin(in1(:,4)+pi.*(8.0./1.5e1))-4.6e1).^2+(in1(:,1)-in1(:,2).*pi.*(1.7e1./3.0e1)+in1(:,3).*sin(in1(:,4)-pi.*(1.7e1./3.0e1))-1.8e1).^2+(in1(:,1)-in1(:,2).*pi.*(1.7e1./3.0e1)+in1(:,3).*sin(in1(:,4)-pi.*(1.7e1./3.0e1))-2.9e1).^2+(in1(:,1)+in1(:,2).*pi.*(2.2e1./4.5e1)+in1(:,3).*sin(in1(:,4)+pi.*(2.2e1./4.5e1))-5.0).^2+(in1(:,1)+in1(:,2).*pi.*(2.3e1./4.5e1)+in1(:,3).*sin(in1(:,4)+pi.*(2.3e1./4.5e1))-2.6e1).^2+(in1(:,1)-in1(:,2).*pi.*(2.6e1./4.5e1)+in1(:,3).*sin(in1(:,4)-pi.*(2.6e1./4.5e1))-3.3e1).^2+(in1(:,1)-in1(:,2).*pi.*(2.6e1./4.5e1)+in1(:,3).*sin(in1(:,4)-pi.*(2.6e1./4.5e1))-3.4e1).^2+(in1(:,1)+in1(:,2).*pi.*(4.7e1./9.0e1)+in1(:,3).*sin(in1(:,4)+pi.*(4.7e1./9.0e1))-3.9e1).^2+(in1(:,1)+in1(:,2).*pi.*(4.7e1./9.0e1)+in1(:,3).*sin(in1(:,4)+pi.*(4.7e1./9.0e1))-4.4e1).^2+(in1(:,1)+in1(:,2).*pi.*(8.9e1./1.8e2)+in1(:,3).*sin(in1(:,4)+pi.*(8.9e1./1.8e2))-3.0).^2+(in1(:,1)+in1(:,2).*pi.*(8.9e1./1.8e2)+in1(:,3).*sin(in1(:,4)+pi.*(8.9e1./1.8e2))-7.0).^2+(in1(:,1)+in1(:,2).*pi.*(9.1e1./1.8e2)+in1(:,3).*sin(in1(:,4)+pi.*(9.1e1./1.8e2))-1.7e1).^2+(in1(:,1)-in1(:,2).*pi.*(1.01e2./1.8e2)+in1(:,3).*sin(in1(:,4)-pi.*(1.01e2./1.8e2))-2.8e1).^2
Solver stopped prematurely.
fmincon stopped because it exceeded the function evaluation limit,
options.MaxFunctionEvaluations = 3000 (the selected value).
ans =
22.5101277310399 1.5446077646475 11.4355301636198 4.51809594253506
>> fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
Local minimum found that satisfies the constraints.
Optimization completed because the objective function is non-decreasing in
feasible directions, to within the default value of the optimality tolerance,
and constraints are satisfied to within the default value of the constraint tolerance.
<stopping criteria details>
ans =
23.5818762479881 -0.285742815508204 1.80149433600663e-07 0.424643259717351
>> fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the default value of the step size tolerance and constraints are
satisfied to within the default value of the constraint tolerance.
<stopping criteria details>
ans =
23.5818765447779 -0.285742717659878 3.17390351376664e-08 2.06965743037983
>> fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
Solver stopped prematurely.
fmincon stopped because it exceeded the function evaluation limit,
options.MaxFunctionEvaluations = 3000 (the selected value).
ans =
22.4306120895266 1.6798711087781 12.2504648986848 4.5366287077669
>> fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options)
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the default value of the step size tolerance and constraints are
satisfied to within the default value of the constraint tolerance.
<stopping criteria details>
ans =
-8.20120347875585 31.0148226862162 299.81610540455 4.62883973782363
If you are seeing an error message, please post it.
Walter Roberson
2017-1-31
coeffs = fmincon(CCV2, rand(1,4)*2, [], [], [], [], [-inf -inf 0 0], [inf inf inf 2*pi], [], options);
F2 = subs(f2, {a0, a1, b0, c1}, {coefs(1), coeffs(2), coeffs(3), coeffs(4)};
sort_theta = sort(theta);
projected_r = double(F2(sort_theta * Pi/180));
polar(sort_theta, projected_r);
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Least Squares 的更多信息
标签
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)