Thread Subject:

how to find pixels inside ellipse ?

hi all

I have an ellipse region in an image and i know the size of
the ellipse axes and its orientation

how can i extract all the pixels inside this ellipse region .

thanks

"elflord woods" <elflord.woods@gmail.com> wrote in message <g1bsct$ggh
$1@fred.mathworks.com>...
> hi all
>
> I have an ellipse region in an image and i know the size of
> the ellipse axes and its orientation
>
> how can i extract all the pixels inside this ellipse region .
>
> thanks
------------------
  I assume you also know the location of the ellipse's center. Call that (x0,y0).
Let t be the counterclockwise angle the major axis makes with respect to the
x-axis. Let a and b be the semi-major and semi-minor axes, respectively. If
P = (x,y) is an arbitrary point then do this:

 X = (x-x0)*cos(t)+(y-y0)*sin(t); % Translate and rotate coords.
 Y = -(x-x0)*sin(t)+(y-y0)*cos(t); % to align with ellipse

If

 X^2/a^2+Y^2/b^2

is less than 1, the point P lies inside the ellipse. If it equals 1, it is right on
the ellipse. If it is greater than 1, P is outside.

Roger Stafford

"Roger Stafford" wrote in message <g1c6l2$b04$1@fred.mathworks.com>...
> "elflord woods" <elflord.woods@gmail.com> wrote in message <g1bsct$ggh
> $1@fred.mathworks.com>...
> > hi all
> >
> > I have an ellipse region in an image and i know the size of
> > the ellipse axes and its orientation
> >
> > how can i extract all the pixels inside this ellipse region .
> >
> > thanks
> ------------------
> I assume you also know the location of the ellipse's center. Call that (x0,y0).
> Let t be the counterclockwise angle the major axis makes with respect to the
> x-axis. Let a and b be the semi-major and semi-minor axes, respectively. If
> P = (x,y) is an arbitrary point then do this:
>
> X = (x-x0)*cos(t)+(y-y0)*sin(t); % Translate and rotate coords.
> Y = -(x-x0)*sin(t)+(y-y0)*cos(t); % to align with ellipse
>
> If
>
> X^2/a^2+Y^2/b^2
>
> is less than 1, the point P lies inside the ellipse. If it equals 1, it is right on
> the ellipse. If it is greater than 1, P is outside.
>
> Roger Stafford
>

I was curious how to do this with lat lon points. I have data describing a rotated ellipse (the center of the ellipse in lat lon coordinates, the lengths of the major and minor axes in kilometers, and the angle that the ellipse is oriented). I do not know the location of the foci, but assume there is a way to figure them out somehow. How can I determine if a specific lat lon point is within this ellipse?

Thanks,
Cody

"Cody Oppermann" wrote in message <j2ellm$gt5$1@newscl01ah.mathworks.com>...
> "Roger Stafford" wrote in message <g1c6l2$b04$1@fred.mathworks.com>...
> > I assume you also know the location of the ellipse's center. Call that (x0,y0).
> > Let t be the counterclockwise angle the major axis makes with respect to the
> > x-axis. Let a and b be the semi-major and semi-minor axes, respectively. If
> > P = (x,y) is an arbitrary point then do this:
> >
> > X = (x-x0)*cos(t)+(y-y0)*sin(t); % Translate and rotate coords.
> > Y = -(x-x0)*sin(t)+(y-y0)*cos(t); % to align with ellipse
> >
> > If
> >
> > X^2/a^2+Y^2/b^2
> >
> > is less than 1, the point P lies inside the ellipse. If it equals 1, it is right on
> > the ellipse. If it is greater than 1, P is outside.
> >
> > Roger Stafford
>
> I was curious how to do this with lat lon points. I have data describing a rotated ellipse (the center of the ellipse in lat lon coordinates, the lengths of the major and minor axes in kilometers, and the angle that the ellipse is oriented). I do not know the location of the foci, but assume there is a way to figure them out somehow. How can I determine if a specific lat lon point is within this ellipse?
>
> Thanks,
> Cody
- - - - - - - - - -
  It depends on the size of the ellipse and your accuracy requirements. The easy way is to approximate the region in the vicinity of the ellipse with a flat plane (as in a map quadrangle) and convert longitude/latitude values of all points to be tested over to east/west and north/south kilometer coordinates relative to the position of the ellipse's center. At the center's latitude, longitude differences correspond approximately to a number of kilometers in the east/west direction as determined by the cosine of the latitude. Latitude differences always correspond to a fixed difference in north/south kilometers (assuming a spherical earth.) Then solve the problem as given earlier in this thread.

  On the other hand if the above is not sufficiently accurate, you are faced with a problem in spherical geometry. To begin with you need to define what an "ellipse" is to be on an (approximately) spherical earth. It can no longer be a planar curve - all planar intersections with a sphere are circles. I won't proceed any further along this line at the moment, assuming hopefully that the procedure in the first paragraph is sufficiently accurate for your purposes.

Roger Stafford

"Roger Stafford" wrote in message <j2eomo$q1e$1@newscl01ah.mathworks.com>...
> "Cody Oppermann" wrote in message <j2ellm$gt5$1@newscl01ah.mathworks.com>...
> > "Roger Stafford" wrote in message <g1c6l2$b04$1@fred.mathworks.com>...
> > > I assume you also know the location of the ellipse's center. Call that (x0,y0).
> > > Let t be the counterclockwise angle the major axis makes with respect to the
> > > x-axis. Let a and b be the semi-major and semi-minor axes, respectively. If
> > > P = (x,y) is an arbitrary point then do this:
> > >
> > > X = (x-x0)*cos(t)+(y-y0)*sin(t); % Translate and rotate coords.
> > > Y = -(x-x0)*sin(t)+(y-y0)*cos(t); % to align with ellipse
> > >
> > > If
> > >
> > > X^2/a^2+Y^2/b^2
> > >
> > > is less than 1, the point P lies inside the ellipse. If it equals 1, it is right on
> > > the ellipse. If it is greater than 1, P is outside.
> > >
> > > Roger Stafford
> >
> > I was curious how to do this with lat lon points. I have data describing a rotated ellipse (the center of the ellipse in lat lon coordinates, the lengths of the major and minor axes in kilometers, and the angle that the ellipse is oriented). I do not know the location of the foci, but assume there is a way to figure them out somehow. How can I determine if a specific lat lon point is within this ellipse?
> >
> > Thanks,
> > Cody
> - - - - - - - - - -
> It depends on the size of the ellipse and your accuracy requirements. The easy way is to approximate the region in the vicinity of the ellipse with a flat plane (as in a map quadrangle) and convert longitude/latitude values of all points to be tested over to east/west and north/south kilometer coordinates relative to the position of the ellipse's center. At the center's latitude, longitude differences correspond approximately to a number of kilometers in the east/west direction as determined by the cosine of the latitude. Latitude differences always correspond to a fixed difference in north/south kilometers (assuming a spherical earth.) Then solve the problem as given earlier in this thread.
>
> On the other hand if the above is not sufficiently accurate, you are faced with a problem in spherical geometry. To begin with you need to define what an "ellipse" is to be on an (approximately) spherical earth. It can no longer be a planar curve - all planar intersections with a sphere are circles. I won't proceed any further along this line at the moment, assuming hopefully that the procedure in the first paragraph is sufficiently accurate for your purposes.
>
> Roger Stafford

Thanks Roger, I was planning on using the flat plane assumption as my ellipses are relatively small. Just to clarify though, say I have a point, for example, at 30°,-30° and the center of my ellipse is at, say, 29°,-29°. Would it be most accurate if x-x0 is the distance from 29°,-29° to 29°,-30° (as opposed to 30°,-30° to 30°,-29°) and y-y0 is the distance from 29°,-29° to 30°,-29° (as opposed to 30°,-30° to 29°,-30°)?
Thanks again.
30,-30 • • 30,-29
29,-30 • • 29,-29

"Cody Oppermann" wrote in message <j2gopq$78p$1@newscl01ah.mathworks.com>...
> Thanks Roger, I was planning on using the flat plane assumption as my ellipses are relatively small. Just to clarify though, say I have a point, for example, at 30°,-30° and the center of my ellipse is at, say, 29°,-29°. Would it be most accurate if x-x0 is the distance from 29°,-29° to 29°,-30° (as opposed to 30°,-30° to 30°,-29°) and y-y0 is the distance from 29°,-29° to 30°,-29° (as opposed to 30°,-30° to 29°,-30°)?
> Thanks again.
> 30,-30 • • 30,-29
> 29,-30 • • 29,-29
- - - - - - - - -
  That's a good question, Cody, but I'm afraid I don't know the answer. However, in any case I am not satisfied with the method I described to you and am now suggesting an approach less prone to error.

  Suppose we define an "ellipse" on the earth's (supposedly) spherical surface as the locus of all points on that surface for which the sum of their great circle arc length distances from two fixed surface focal points is a given constant. If the location of such foci for your ellipse can be found, it is easy to determine whether points lie within it by comparing such a sum of arc lengths with the given constant.

  Transform all latitude, longitude values for points to earth-centered cartesian coordinates using the equations:

 x = R*cosd(theta)*cosd(phi)
 y = R*cosd(theta)*sind(phi)
 z = R*sind(theta)

where R is the earth's radius, theta is the (north-positive) latitude, and phi is the (east-positive) longitude and where theta and phi are in degrees.

  Let v1 = [x1,y1,z1] and v2 = [x2,y2,z2] be vectors for the two given foci and let v = [x,y,z] be the vector for an arbitrary point to be tested. (These are vectors pointing from the earth's center to the surface points.) Then the sum of the two great circle arc lengths referred to above is:

 t = R*(atan2(norm(cross(v,v1)),dot(v,v1)) + ...
        atan2(norm(cross(v,v2)),dot(v,v2)));

The point (x,y,z) will lie within the ellipse whenever t is less than or equal to the above given constant.

  (Note: You can actually assume that R is equal to one in all the above provided you make the corresponding correction to the assumed arc length constant sum for the ellipse.)

  I would venture to guess that, provided you use double precision, the major source of inaccuracy in the above would be the assumption that the earth is a sphere, as opposed to an oblate spheroid. However the computations for the latter would be much more complicated.

Roger Stafford

"Roger Stafford" wrote in message <j2hcup$fh9$1@newscl01ah.mathworks.com>...
> "Cody Oppermann" wrote in message <j2gopq$78p$1@newscl01ah.mathworks.com>...
> > Thanks Roger, I was planning on using the flat plane assumption as my ellipses are relatively small. Just to clarify though, say I have a point, for example, at 30°,-30° and the center of my ellipse is at, say, 29°,-29°. Would it be most accurate if x-x0 is the distance from 29°,-29° to 29°,-30° (as opposed to 30°,-30° to 30°,-29°) and y-y0 is the distance from 29°,-29° to 30°,-29° (as opposed to 30°,-30° to 29°,-30°)?
> > Thanks again.
> > 30,-30 • • 30,-29
> > 29,-30 • • 29,-29
> - - - - - - - - -
> That's a good question, Cody, but I'm afraid I don't know the answer. However, in any case I am not satisfied with the method I described to you and am now suggesting an approach less prone to error.
>
> Suppose we define an "ellipse" on the earth's (supposedly) spherical surface as the locus of all points on that surface for which the sum of their great circle arc length distances from two fixed surface focal points is a given constant. If the location of such foci for your ellipse can be found, it is easy to determine whether points lie within it by comparing such a sum of arc lengths with the given constant.
>
> Transform all latitude, longitude values for points to earth-centered cartesian coordinates using the equations:
>
> x = R*cosd(theta)*cosd(phi)
> y = R*cosd(theta)*sind(phi)
> z = R*sind(theta)
>
> where R is the earth's radius, theta is the (north-positive) latitude, and phi is the (east-positive) longitude and where theta and phi are in degrees.
>
> Let v1 = [x1,y1,z1] and v2 = [x2,y2,z2] be vectors for the two given foci and let v = [x,y,z] be the vector for an arbitrary point to be tested. (These are vectors pointing from the earth's center to the surface points.) Then the sum of the two great circle arc lengths referred to above is:
>
> t = R*(atan2(norm(cross(v,v1)),dot(v,v1)) + ...
> atan2(norm(cross(v,v2)),dot(v,v2)));
>
> The point (x,y,z) will lie within the ellipse whenever t is less than or equal to the above given constant.
>
> (Note: You can actually assume that R is equal to one in all the above provided you make the corresponding correction to the assumed arc length constant sum for the ellipse.)
>
> I would venture to guess that, provided you use double precision, the major source of inaccuracy in the above would be the assumption that the earth is a sphere, as opposed to an oblate spheroid. However the computations for the latter would be much more complicated.
>
> Roger Stafford

This sounds reasonable, Roger, but I think Cody's problem can be solved more easily (and just as easily on an oblate spheroid as on a sphere) by using the distance function in Mapping Toolbox.

The Mapping Toolbox function ellipse1 defines an ellipse on the surface of a sphere or spheroid by generalizing the parametric form of a planar ellipse:

x = a*cos(theta)
y = b*sin(theta)

where the center is a (0,0), the semimajor and semiminor axes are a and b, and the major axis is aligned with the X-axis in the plane. Taking the same approach to this problem, consider an ellipse on a sphere or spheroid defined by its center -- at (lat0,lon0), semimajor and semiminor axes a and b, and a rotation angle between the axes and the local meridian on which the center falls (the meridian for which lon = lon0).

1. Use the distance function to compute the great circle (geodesic) distance d from the center (lat0, lon0) to a point in question, along with the azimuth angle (the second output of distance -- no need for a separate call to the azimuth function).

2. Determine theta from the computed azimuth and the rotation angle for the ellipse.

3. The distance to the elliptical curve, following a geodesic that departs from the center at the computed azimuth, is simply defined as e = hypot(a*cos(theta), b*sin(theta)).

4. If d < e, the point is inside the ellipse. If d == e the point is on the ellipse. If d > e the point is outside the ellipse.

Hope this helps,

Rob Comer
Mapping Toolbox Development
MathWorks

Interesting we have two definitions of ellipse on spherical surface:

1. sum geodesic distances to two focii is constant (Roger's)
2. Generalize of (x/a)^2 + (y/b)^2 = 1 in Mapping tool box (Rob's)

I would guess for many applications using satellite data, the intersection between a 3D cone and spherical surface would be also considered as "elliptical". This is a third definition.

I would guess the three are really distinct.

Which one is the definition is most appropriate to Cody's question?

Bruno

"Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2ja2d$ebk$1@newscl01ah.mathworks.com>...
> This sounds reasonable, Roger, but I think Cody's problem can be solved more easily (and just as easily on an oblate spheroid as on a sphere) by using the distance function in Mapping Toolbox.
>
> The Mapping Toolbox function ellipse1 defines an ellipse on the surface of a sphere or spheroid by generalizing the parametric form of a planar ellipse:
>
> x = a*cos(theta)
> y = b*sin(theta)
>
> where the center is a (0,0), the semimajor and semiminor axes are a and b, and the major axis is aligned with the X-axis in the plane. Taking the same approach to this problem, consider an ellipse on a sphere or spheroid defined by its center -- at (lat0,lon0), semimajor and semiminor axes a and b, and a rotation angle between the axes and the local meridian on which the center falls (the meridian for which lon = lon0).
>
> 1. Use the distance function to compute the great circle (geodesic) distance d from the center (lat0, lon0) to a point in question, along with the azimuth angle (the second output of distance -- no need for a separate call to the azimuth function).
>
> 2. Determine theta from the computed azimuth and the rotation angle for the ellipse.
>
> 3. The distance to the elliptical curve, following a geodesic that departs from the center at the computed azimuth, is simply defined as e = hypot(a*cos(theta), b*sin(theta)).
>
> 4. If d < e, the point is inside the ellipse. If d == e the point is on the ellipse. If d > e the point is outside the ellipse.
>
> Hope this helps,
>
> Rob Comer
> Mapping Toolbox Development
> MathWorks
- - - - - - - - -
Hello Rob,

  If I understand your parametric method correctly, it is equivalent to generalizing the polar coordinate form of a planar ellipse. The cartesian coordinate equation of an ellipse

 x^2/a^2 + y^2/b^2 = 1

has the polar coordinate form

 r = 1/sqrt((cos(phi)/a)^2+(sin(phi)/b)^2)

where x = r*cos(phi), y = r*sin(phi).

  As applied to a spherical or ellipsoidal surface this latter equation can be used to generalize r to a geodesic distance from the ellipse center to a point on the ellipse and phi to the angle between that geodesic and the major (or minor) axis at the center of the ellipse. As you say, the Mapping Toolbox can readily be used for this purpose.

  I agree with Bruno that this definition of an ellipse may well differ from one involving a constant sum of distances from two fixed foci as applied to a spherical or ellipsoidal surface.

Roger Stafford

"Roger Stafford" wrote in message <j2l4p3$cus$1@newscl01ah.mathworks.com>...
> "Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2ja2d$ebk$1@newscl01ah.mathworks.com>...
> > This sounds reasonable, Roger, but I think Cody's problem can be solved more easily (and just as easily on an oblate spheroid as on a sphere) by using the distance function in Mapping Toolbox.
> >
> > The Mapping Toolbox function ellipse1 defines an ellipse on the surface of a sphere or spheroid by generalizing the parametric form of a planar ellipse:
> >
> > x = a*cos(theta)
> > y = b*sin(theta)
> >
> > where the center is a (0,0), the semimajor and semiminor axes are a and b, and the major axis is aligned with the X-axis in the plane. Taking the same approach to this problem, consider an ellipse on a sphere or spheroid defined by its center -- at (lat0,lon0), semimajor and semiminor axes a and b, and a rotation angle between the axes and the local meridian on which the center falls (the meridian for which lon = lon0).
> >
> > 1. Use the distance function to compute the great circle (geodesic) distance d from the center (lat0, lon0) to a point in question, along with the azimuth angle (the second output of distance -- no need for a separate call to the azimuth function).
> >
> > 2. Determine theta from the computed azimuth and the rotation angle for the ellipse.
> >
> > 3. The distance to the elliptical curve, following a geodesic that departs from the center at the computed azimuth, is simply defined as e = hypot(a*cos(theta), b*sin(theta)).
> >
> > 4. If d < e, the point is inside the ellipse. If d == e the point is on the ellipse. If d > e the point is outside the ellipse.
> >
> > Hope this helps,
> >
> > Rob Comer
> > Mapping Toolbox Development
> > MathWorks
> - - - - - - - - -
> Hello Rob,
>
> If I understand your parametric method correctly, it is equivalent to generalizing the polar coordinate form of a planar ellipse. The cartesian coordinate equation of an ellipse
>
> x^2/a^2 + y^2/b^2 = 1
>
> has the polar coordinate form
>
> r = 1/sqrt((cos(phi)/a)^2+(sin(phi)/b)^2)
>
> where x = r*cos(phi), y = r*sin(phi).
>
> As applied to a spherical or ellipsoidal surface this latter equation can be used to generalize r to a geodesic distance from the ellipse center to a point on the ellipse and phi to the angle between that geodesic and the major (or minor) axis at the center of the ellipse. As you say, the Mapping Toolbox can readily be used for this purpose.

Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob

> I agree with Bruno that this definition of an ellipse may well differ from one involving a constant sum of distances from two fixed foci as applied to a spherical or ellipsoidal surface.

I also agree with Bruno on that. Of course, for small ellipses (relative to the size of the planet) they are similar. And the definition used in ellipse1 is pragmatic and relatively easy to compute. The main application of which I'm aware is to indicate the degree of uncertainty with respect to the location of a point on the surface -- an earthquake epicenter, for example. Here a and b relate to the eigenvalues of a covariance matrix, and everything is somewhat approximate to begin with. (Cody, what's your application?)

Intersecting a cone with the surface, as Bruno suggested, is a legitimate and interesting problem also. I'd have to study it before declaring the answer to be an ellipse. Certainly the cone can miss the planet entirely, resulting in the empty set, or partially, resulting in a curve that has sharp cusps where it needs to "shortcut" back along a great circle.
-- Rob

"Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2lfef$coo$1@newscl01ah.mathworks.com>...
> "Roger Stafford" wrote in message <j2l4p3$cus$1@newscl01ah.mathworks.com>...
> > "Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2ja2d$ebk$1@newscl01ah.mathworks.com>...
> > > This sounds reasonable, Roger, but I think Cody's problem can be solved more easily (and just as easily on an oblate spheroid as on a sphere) by using the distance function in Mapping Toolbox.
> > >
> > > The Mapping Toolbox function ellipse1 defines an ellipse on the surface of a sphere or spheroid by generalizing the parametric form of a planar ellipse:
> > >
> > > x = a*cos(theta)
> > > y = b*sin(theta)
> > >
> > > where the center is a (0,0), the semimajor and semiminor axes are a and b, and the major axis is aligned with the X-axis in the plane. Taking the same approach to this problem, consider an ellipse on a sphere or spheroid defined by its center -- at (lat0,lon0), semimajor and semiminor axes a and b, and a rotation angle between the axes and the local meridian on which the center falls (the meridian for which lon = lon0).
> > >
> > > 1. Use the distance function to compute the great circle (geodesic) distance d from the center (lat0, lon0) to a point in question, along with the azimuth angle (the second output of distance -- no need for a separate call to the azimuth function).
> > >
> > > 2. Determine theta from the computed azimuth and the rotation angle for the ellipse.
> > >
> > > 3. The distance to the elliptical curve, following a geodesic that departs from the center at the computed azimuth, is simply defined as e = hypot(a*cos(theta), b*sin(theta)).
> > >
> > > 4. If d < e, the point is inside the ellipse. If d == e the point is on the ellipse. If d > e the point is outside the ellipse.
> > >
> > > Hope this helps,
> > >
> > > Rob Comer
> > > Mapping Toolbox Development
> > > MathWorks
> > - - - - - - - - -
> > Hello Rob,
> >
> > If I understand your parametric method correctly, it is equivalent to generalizing the polar coordinate form of a planar ellipse. The cartesian coordinate equation of an ellipse
> >
> > x^2/a^2 + y^2/b^2 = 1
> >
> > has the polar coordinate form
> >
> > r = 1/sqrt((cos(phi)/a)^2+(sin(phi)/b)^2)
> >
> > where x = r*cos(phi), y = r*sin(phi).
> >
> > As applied to a spherical or ellipsoidal surface this latter equation can be used to generalize r to a geodesic distance from the ellipse center to a point on the ellipse and phi to the angle between that geodesic and the major (or minor) axis at the center of the ellipse. As you say, the Mapping Toolbox can readily be used for this purpose.
>
> Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob
>
> > I agree with Bruno that this definition of an ellipse may well differ from one involving a constant sum of distances from two fixed foci as applied to a spherical or ellipsoidal surface.
>
> I also agree with Bruno on that. Of course, for small ellipses (relative to the size of the planet) they are similar. And the definition used in ellipse1 is pragmatic and relatively easy to compute. The main application of which I'm aware is to indicate the degree of uncertainty with respect to the location of a point on the surface -- an earthquake epicenter, for example. Here a and b relate to the eigenvalues of a covariance matrix, and everything is somewhat approximate to begin with. (Cody, what's your application?)
>
> Intersecting a cone with the surface, as Bruno suggested, is a legitimate and interesting problem also. I'd have to study it before declaring the answer to be an ellipse. Certainly the cone can miss the planet entirely, resulting in the empty set, or partially, resulting in a curve that has sharp cusps where it needs to "shortcut" back along a great circle.
> -- Rob

After bringing my head around geometry again, I was able to figure something out that would work for me. I don't have the latest version of Matlab so I don't have the ellipse1 function, but I was able to locate the lat lon of the foci and then as long as the sum of the distances from the point of interest to each foci is less than 2a (the major axis radii) then it was within the ellipse. Thanks for all the help!
-Cody

"Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2lfef$coo$1@newscl01ah.mathworks.com>...
> Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob
> .........
- - - - - - - - - -
  I am not sure what you mean by "figuratively", Rob. The equivalence is mathematically exact if I understand your method.

  In your step 2. you say "Determine theta from the computed azimuth and the rotation angle for the ellipse" where theta is the parameter in the equations

 x = a*cos(theta)
 y = b*sin(theta)

Combining your azimuth angle with the ellipse's rotation angle should give an angle phi between the (geodesic) path and the major axis. To get from angle phi to the parameter theta it is necessary to make a transformation since these are very different quantities:

 tan(theta) = a/b*y/x = a/b*tan(phi)
 theta = atan(a/b*tan(phi))

This what I understood you meant to be done in step 2. (If that is not the case then I simply don't understand your step 2.) Then in step 3. you do

 e = sqrt((a*cos(theta))^2+(b*sin(theta))^2)

  This is the precise equivalent of the polar coordinate procedure I outlined where r is e. It is not necessary to bring the parameter theta into the computation since the distance can be computed directly from a, b, and phi by the simple formula I gave. The results will be exactly the same.

Roger Stafford

"Roger Stafford" wrote in message <j2m9gp$c3s$1@newscl01ah.mathworks.com>...
> "Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2lfef$coo$1@newscl01ah.mathworks.com>...
> > Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob
> > .........
> - - - - - - - - - -
> I am not sure what you mean by "figuratively", Rob. The equivalence is mathematically exact if I understand your method.
>
> In your step 2. you say "Determine theta from the computed azimuth and the rotation angle for the ellipse" where theta is the parameter in the equations
>
> x = a*cos(theta)
> y = b*sin(theta)
>
> Combining your azimuth angle with the ellipse's rotation angle should give an angle phi between the (geodesic) path and the major axis. To get from angle phi to the parameter theta it is necessary to make a transformation since these are very different quantities:
>
> tan(theta) = a/b*y/x = a/b*tan(phi)
> theta = atan(a/b*tan(phi))
>
> This what I understood you meant to be done in step 2. (If that is not the case then I simply don't understand your step 2.) Then in step 3. you do
>
> e = sqrt((a*cos(theta))^2+(b*sin(theta))^2)
>
> This is the precise equivalent of the polar coordinate procedure I outlined where r is e. It is not necessary to bring the parameter theta into the computation since the distance can be computed directly from a, b, and phi by the simple formula I gave. The results will be exactly the same.
>
> Roger Stafford

By figuratively I just meant that I hadn't taken time to check your expression for r and be certain that phi was the polar angle. Now I have. It looks correct, as does your expression for theta in terms of phi.

Your understanding of my step 2 is correct -- except that it's probably most natural to use an azimuth angle rather than a polar angle. If you have Mapping Toolbox, look at the az2parametric subfunction in ellipse1.m.

Cheers,
Rob

"Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2qsqg$583$1@newscl01ah.mathworks.com>...
> By figuratively I just meant that I hadn't taken time to check your expression for r and be certain that phi was the polar angle. Now I have. It looks correct, as does your expression for theta in terms of phi.
>
> Your understanding of my step 2 is correct -- except that it's probably most natural to use an azimuth angle rather than a polar angle. If you have Mapping Toolbox, look at the az2parametric subfunction in ellipse1.m.
>
> Cheers,
> Rob
- - - - - - - - - -
Hello again Rob,

  Now that I am sure of the method used in ellipse1 (I don't have it on my system) I have compared the curve presumably generated by it with the one defined in terms of two foci. They do generate different curves for a spherical surface, though these are rather close to one another for reasonably sized major and minor axes relative to the radius of the sphere. The following is an example.

  Let the major axis length be 1 radian and the minor axis .5 radians. For the curve given by ellipse1 choose a point P making an angle at the center of .3 radians with respect to the major axis.

 a = 1;
 b = .5;
 phi = .3;

Then we have for the arc distance from the center to P from ellipse1

 r = 1/sqrt(cos(phi)^2/a^2+sin(phi)^2/b^2);

If two points along the major axis are to be valid foci, the sum of the arc distances from them to the end of the minor axis must equal the sum of the arc distances from them to the end of the major axis, so the foci must be located on the major axis a distance of c on either side of the center where c is:

 c = acos(cos(a)/cos(b));

The arc distances from P to these two foci would then be:

 d1 = acos(cos(r)*cos(c)+sin(r)*sin(c)*cos(phi));
 d2 = acos(cos(r)*cos(c)-sin(r)*sin(c)*cos(phi));

These calculations result in

 d1+d2 = 2.00481761

whereas they should have added up to 2*a = 2.0 if this were an "ellipse" according to the foci definition. The curves are therefore very close but definitely not the same.

Roger Stafford

"Roger Stafford" wrote in message <j2m9gp$c3s$1@newscl01ah.mathworks.com>...
> "Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2lfef$coo$1@newscl01ah.mathworks.com>...
> > Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob
> > .........
> - - - - - - - - - -
> I am not sure what you mean by "figuratively", Rob. The equivalence is mathematically exact if I understand your method.
>
> In your step 2. you say "Determine theta from the computed azimuth and the rotation angle for the ellipse" where theta is the parameter in the equations
>
> x = a*cos(theta)
> y = b*sin(theta)
>
> Combining your azimuth angle with the ellipse's rotation angle should give an angle phi between the (geodesic) path and the major axis. To get from angle phi to the parameter theta it is necessary to make a transformation since these are very different quantities:
>
> tan(theta) = a/b*y/x = a/b*tan(phi)
> theta = atan(a/b*tan(phi))
>
> This what I understood you meant to be done in step 2. (If that is not the case then I simply don't understand your step 2.) Then in step 3. you do
>
> e = sqrt((a*cos(theta))^2+(b*sin(theta))^2)
>
> This is the precise equivalent of the polar coordinate procedure I outlined where r is e. It is not necessary to bring the parameter theta into the computation since the distance can be computed directly from a, b, and phi by the simple formula I gave. The results will be exactly the same.
>
> Roger Stafford

Hi Roger,
How do I calculate PHI ? I know the ( lat_center, lon_center) of center of ellipse, lat lon of point to be tested ( lat_point,lon_ppoint) and also know another point through which the major axis passes ( lat_another,lon_another).
Thanks in advance
GlanPlon

Until Roger answers you, you might be interested in the Pick of the
Week this week that talks about ellipses:
http://blogs.mathworks.com/pick/2011/10/21/drawing-ellipses/

"glanplon" wrote in message <j878ve$770$1@newscl01ah.mathworks.com>...
> "Roger Stafford" wrote in message <j2m9gp$c3s$1@newscl01ah.mathworks.com>...
> > "Rob Comer" <rob.comer.nospam@mathworks.com> wrote in message <j2lfef$coo$1@newscl01ah.mathworks.com>...
> > > Roger, I agree that this is the case at least figuratively. One does need take into account that the polar angle and parametric angle are not equivalent except when the ellipse is a circle (just as the geocentric and parametric latitudes on an oblate spheroid are not the same). -- Rob
> > > .........
> > - - - - - - - - - -
> > I am not sure what you mean by "figuratively", Rob. The equivalence is mathematically exact if I understand your method.
> >
> > In your step 2. you say "Determine theta from the computed azimuth and the rotation angle for the ellipse" where theta is the parameter in the equations
> >
> > x = a*cos(theta)
> > y = b*sin(theta)
> >
> > Combining your azimuth angle with the ellipse's rotation angle should give an angle phi between the (geodesic) path and the major axis. To get from angle phi to the parameter theta it is necessary to make a transformation since these are very different quantities:
> >
> > tan(theta) = a/b*y/x = a/b*tan(phi)
> > theta = atan(a/b*tan(phi))
> >
> > This what I understood you meant to be done in step 2. (If that is not the case then I simply don't understand your step 2.) Then in step 3. you do
> >
> > e = sqrt((a*cos(theta))^2+(b*sin(theta))^2)
> >
> > This is the precise equivalent of the polar coordinate procedure I outlined where r is e. It is not necessary to bring the parameter theta into the computation since the distance can be computed directly from a, b, and phi by the simple formula I gave. The results will be exactly the same.
> >
> > Roger Stafford
>
> Hi Roger,
> How do I calculate PHI ? I know the ( lat_center, lon_center) of center of ellipse, lat lon of point to be tested ( lat_point,lon_ppoint) and also know another point through which the major axis passes ( lat_another,lon_another).
> Thanks in advance
> GlanPlon

Hi ,
Is PHI the dihedral angle between planes passing through center-point and center-anotherpoint. Can this dihedral be calculated by using method of computation at
http://en.wikipedia.org/wiki/Dihedral_angle
GlanPlon

ImageAnalyst <imageanalyst@mailinator.com> wrote in message <881adacb-23e9-4ded-b9df-0d49c2d924c4@u10g2000prl.googlegroups.com>...
> Until Roger answers you, you might be interested in the Pick of the
> Week this week that talks about ellipses:
> http://blogs.mathworks.com/pick/2011/10/21/drawing-ellipses/

Anybody there....