- You will see updates in your activity feed.
- You may receive emails, depending on your notification preferences.

35 views (last 30 days)

Show older comments

Please can anyone tell me how to do it?Any way will be fine

Walter Roberson
on 23 Jan 2017

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
on 23 Jan 2017

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
on 28 Jan 2017

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
on 28 Jan 2017

Walter Roberson
on 28 Jan 2017

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
on 29 Jan 2017

Image Analyst
on 28 Jan 2017

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-');

L K
on 29 Jan 2017

L K
on 29 Jan 2017

Image Analyst
on 29 Jan 2017

L K
on 30 Jan 2017

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
on 30 Jan 2017

Walter Roberson
on 31 Jan 2017

Requires R2014b or later.

For releases before that look in the File Exchange for "alpha shapes"

L K
on 31 Jan 2017

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
on 31 Jan 2017

L K
on 31 Jan 2017

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
on 31 Jan 2017

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
on 1 Feb 2017

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
on 1 Feb 2017

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
on 2 Feb 2017

L K
on 10 Feb 2017

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
on 10 Feb 2017

Walter Roberson
on 10 Feb 2017

L K
on 22 Mar 2017

Walter Roberson
on 22 Mar 2017

However, it is my opinion that the result you are coming up with has no useful meaning.

Walter Roberson
on 22 Mar 2017

L K
on 22 Mar 2017

Image Analyst
on 22 Mar 2017

L K
on 23 Mar 2017

L K
on 23 Mar 2017

Walter Roberson
on 23 Jan 2017

Walter Roberson
on 28 Jan 2017

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
on 29 Jan 2017

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
on 29 Jan 2017

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
on 29 Jan 2017

Walter Roberson
on 30 Jan 2017

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
on 31 Jan 2017

Walter Roberson
on 31 Jan 2017

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
on 31 Jan 2017

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);

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!Unable to complete the action because of changes made to the page. Reload the page to see its updated state.

Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .

Select web siteYou can also select a web site from the following list:

Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.

- 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)