Thread Subject: Deciding which side a numerically generated line you are on

I am looking at a number of functions of two variables (say x and y), and I am restricting my attention to x \in [-pi/2, pi/2] and y \in [-pi, pi]. Call this set S \subset R^2. The problem is such that I really only want to look at some of the data in this subset of the plane. I numerically generate some boundary curves that separate feasible (x,y) coordinates from infeasible ones (by solving a system of two equations in two unknowns), and I need to be able to decide if an arbitrary point in S is on one side or the other side of the numerically generated curve. I have the points on the curve as a 2xn array, each column of which is a point (x,y) on the curve.

Here is a link to what the boundary curves look like:
http://www.dlpeterson.com/feasible_lean_steer.pdf

The reason for why I need to do this is because on one side of my curve I can't compute a number of the quantities and if I try, I run into a slew of other numerical non-convergence issues due to the inherent infeasibility of points on that side of the line.

Any ideas?

Thanks,
Luke

"Dale" <hazelREMOVETHISnusse@gmail.com> wrote in message <gkr5ue$8b5$1@fred.mathworks.com>...
> I am looking at a number of functions of two variables (say x and y), and I am restricting my attention to x \in [-pi/2, pi/2] and y \in [-pi, pi]. Call this set S \subset R^2. The problem is such that I really only want to look at some of the data in this subset of the plane. I numerically generate some boundary curves that separate feasible (x,y) coordinates from infeasible ones (by solving a system of two equations in two unknowns), and I need to be able to decide if an arbitrary point in S is on one side or the other side of the numerically generated curve. I have the points on the curve as a 2xn array, each column of which is a point (x,y) on the curve.
>
> Here is a link to what the boundary curves look like:
> http://www.dlpeterson.com/feasible_lean_steer.pdf
>
> The reason for why I need to do this is because on one side of my curve I can't compute a number of the quantities and if I try, I run into a slew of other numerical non-convergence issues due to the inherent infeasibility of points on that side of the line.
>
> Any ideas?
>
> Thanks,
> Luke

Hello Luke,

Just for brainstorming purposes.

You know that the power of a point inside a circle is negative. It is zero if it is on the boundary and it is positive if it is outside. This power is given by:

F(x,y) = (x-xc)^2 + (y-yc)^2 - r^2

for a circle of radius r centered at (xc,yc).

Looking at the curves you have given, if you could use their equations in a similar manner, you could define the power of a point with respect to a given curve. Heuristically, I would expect to see negative power values in the upper right infeasible region with respect to the again upper right curve because the points are in the "interior" region of that curve.

Just an idea. Hope it helps a bit.

Thanks.

"Sadik " <sadik.hava@gmail.com> wrote in message <gkr7it$4s4$1@fred.mathworks.com>...
> "Dale" <hazelREMOVETHISnusse@gmail.com> wrote in message <gkr5ue$8b5$1@fred.mathworks.com>...
> > I am looking at a number of functions of two variables (say x and y), and I am restricting my attention to x \in [-pi/2, pi/2] and y \in [-pi, pi]. Call this set S \subset R^2. The problem is such that I really only want to look at some of the data in this subset of the plane. I numerically generate some boundary curves that separate feasible (x,y) coordinates from infeasible ones (by solving a system of two equations in two unknowns), and I need to be able to decide if an arbitrary point in S is on one side or the other side of the numerically generated curve. I have the points on the curve as a 2xn array, each column of which is a point (x,y) on the curve.
> >
> > Here is a link to what the boundary curves look like:
> > http://www.dlpeterson.com/feasible_lean_steer.pdf
> >
> > The reason for why I need to do this is because on one side of my curve I can't compute a number of the quantities and if I try, I run into a slew of other numerical non-convergence issues due to the inherent infeasibility of points on that side of the line.
> >
> > Any ideas?
> >
> > Thanks,
> > Luke
>
> Hello Luke,
>
> Just for brainstorming purposes.
>
> You know that the power of a point inside a circle is negative. It is zero if it is on the boundary and it is positive if it is outside. This power is given by:
>
> F(x,y) = (x-xc)^2 + (y-yc)^2 - r^2
>
> for a circle of radius r centered at (xc,yc).
>
> Looking at the curves you have given, if you could use their equations in a similar manner, you could define the power of a point with respect to a given curve. Heuristically, I would expect to see negative power values in the upper right infeasible region with respect to the again upper right curve because the points are in the "interior" region of that curve.
>
> Just an idea. Hope it helps a bit.
>
> Thanks.

By the way, you have said that the curves are generated numerically. So if you wish, you could fit a polynomial to obtain such equations.

If you think of the simplest parabola y = x^2, then the power should be F = x^2-y
so that every point in the interior of the parabola will give negative values. Those on the parabola will give zero and the exterior points will yield positive values.

"Dale" <hazelREMOVETHISnusse@gmail.com> wrote in message <gkr5ue$8b5$1@fred.mathworks.com>...
>..... I need to be able to decide if an arbitrary point in S is on one side or the other side of the numerically generated curve. .......

  Your question only has true significance if the "boundary curves" are closed. Otherwise you would be able to move a test point continuously from one side of a point on the curve to the other side without ever crossing the curve at any point - there could be no inside or outside. If the boundary curve consists of a discrete set of vertex points and the segments connecting them, forming a closed polygon, matlab can give a definite answer to this question.

  Suppose P = (x,y) is a test point and P1 = (x1,y1) and P2 = (x2,y2) are two successive vertices along the polygon. Then the angle at P in triangle P1P2P is given by

 a12 = atan2((x-x1)*(y-y2)-(y-y1)*(x-x2),(x-x1)*(x-x2)+(y-y1)*(y-y2));

where it is understood that a12 will be positive if the path P --> P1 --> P2 --> P is in a counterclockwise direction and negative if clockwise.

  If these angles are summed for all successive line segments of the polygon, and if the point P is interior to the polygon, these angles will have a sum of 2*pi for a counterclockwise direction around the polygon and -2*pi for a clockwise direction. On the other hand the sum will be zero for all exterior points. This does not depend on the polygon being convex, though I have assumed that the boundary curve is simple in the sense that it never crosses over itself.

  Note that 'atan2 ' can accept vectors as inputs so these calculations can be done in a vectorized manner (without loops) for an entire polygon relative to a point (x,y).

Roger Stafford

This makes perfect sense, thank you! I'll let you know once I implement it, hopefully it should be as straightforward as you made it sound :)
Thanks,
~Luke

"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message <gkrcu6$sh6$1@fred.mathworks.com>...
> "Dale" <hazelREMOVETHISnusse@gmail.com> wrote in message <gkr5ue$8b5$1@fred.mathworks.com>...
> >..... I need to be able to decide if an arbitrary point in S is on one side or the other side of the numerically generated curve. .......
>
> Your question only has true significance if the "boundary curves" are closed. Otherwise you would be able to move a test point continuously from one side of a point on the curve to the other side without ever crossing the curve at any point - there could be no inside or outside. If the boundary curve consists of a discrete set of vertex points and the segments connecting them, forming a closed polygon, matlab can give a definite answer to this question.
>
> Suppose P = (x,y) is a test point and P1 = (x1,y1) and P2 = (x2,y2) are two successive vertices along the polygon. Then the angle at P in triangle P1P2P is given by
>
> a12 = atan2((x-x1)*(y-y2)-(y-y1)*(x-x2),(x-x1)*(x-x2)+(y-y1)*(y-y2));
>
> where it is understood that a12 will be positive if the path P --> P1 --> P2 --> P is in a counterclockwise direction and negative if clockwise.
>
> If these angles are summed for all successive line segments of the polygon, and if the point P is interior to the polygon, these angles will have a sum of 2*pi for a counterclockwise direction around the polygon and -2*pi for a clockwise direction. On the other hand the sum will be zero for all exterior points. This does not depend on the polygon being convex, though I have assumed that the boundary curve is simple in the sense that it never crosses over itself.
>
> Note that 'atan2 ' can accept vectors as inputs so these calculations can be done in a vectorized manner (without loops) for an entire polygon relative to a point (x,y).
>
> Roger Stafford

"Dale" <hazelREMOVETHISnusse@gmail.com> wrote in message <gkrkbl$21h$1@fred.mathworks.com>...
> This makes perfect sense, thank you! I'll let you know once I implement it, hopefully it should be as straightforward as you made it sound :)
> Thanks,
> ~Luke

  Dale, I have read your original question more carefully and realize that my remarks really apply to situations without bounds. In your case of being confined to a rectangular region, when and if your boundary curve (polygon) meets a bound of the rectangular region, you would want to continue the polygon on around one side or the other of the rectangle so as to be closed. Then the procedure I described could be applied.

  Second, you said "I numerically generate some boundary curves that separate feasible (x,y) coordinates from infeasible ones (by solving a system of two equations in two unknowns)" and also "on one side of my curve I can't compute a number of the quantities." You probably meant one equation in two unknowns here so as to generate a one-dimensional curve. I take it that your presumably single equation is such that for many points it makes no sense to attempt to evaluate the expressions on the two side of it. Otherwise, a direct evaluation of these expressions would seem to be the most efficient way of determining which side of the curve a given point (x,y) lies on, by the simple expedient of which expression is the greater.

  Finally I will give you my version of implementing the procedure I described earlier. Let X and Y be column vectors of the coordinates of successive vertices of a closed polygon which is to be one of your boundary curves. Let x and y be scalar coordinates of a test point P = (x,y).

 x1 = [X(end);X}; y1 = [Y(end);Y};
 x2 = [X;X(1)]; y2 = [Y;Y(1)];
 a = sum(atan2((x-x1).*(y-y2)-(y-y1).*(x-x2), ...
               (x-x1).*(x-x2)+(y-y1).*(y-y2)));

The value of angle 'a' then shows whether (x,y) is inside or outside.

Roger Stafford

"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message <gkrp96$8rn$1@fred.mathworks.com>...
> ......
> Finally I will give you my version of implementing the procedure I described earlier. Let X and Y be column vectors of the coordinates of successive vertices of a closed polygon which is to be one of your boundary curves. Let x and y be scalar coordinates of a test point P = (x,y).
>
> x1 = [X(end);X}; y1 = [Y(end);Y};
> x2 = [X;X(1)]; y2 = [Y;Y(1)];
> a = sum(atan2((x-x1).*(y-y2)-(y-y1).*(x-x2), ...
> (x-x1).*(x-x2)+(y-y1).*(y-y2)));
>
> The value of angle 'a' then shows whether (x,y) is inside or outside.

  I wish to amend that code I sent you, Dale. One of the segments was repeated and this would give an error. It should read as follows.

  Let X and Y be column vectors of the coordinates of successive vertices of a closed polygon which is to be one of your boundary curves. Let x and y be scalar coordinates of a test point P = (x,y).

 X1 = X([2:end,1]); Y1 = Y([2:end,1]);
 a = sum(atan2((x-X).*(y-Y1)-(y-Y).*(x-X1), ...
               (x-X).*(x-X1)+(y-Y).*(y-Y1)));

Roger Stafford

On 17 Jan, 00:41, "Dale" <hazelREMOVETHISnu...@gmail.com> wrote:
> I am looking at a number of functions of two variables (say x and y), and=
 I am restricting my attention to x \in [-pi/2, pi/2] and y \in [-pi, pi]. =
=A0Call this set S \subset R^2. =A0The problem is such that I really only w=
ant to look at some of the data in this subset of the plane. =A0I numerical=
ly generate some boundary curves that separate feasible (x,y) coordinates f=
rom infeasible ones (by solving a system of two equations in two unknowns),=
 and I need to be able to decide if an arbitrary point in S is on one side =
or the other side of the numerically generated curve. =A0I have the points =
on the curve as a 2xn array, each column of which is a point (x,y) on the c=
urve. =A0

The terminology you use would suggest that this is a constrained
linear programming problem. There are well-known methods to solve
such problems, but they are a bit complicated to learn and
implement. Linear programming tecniques work in any dimension.

The 2D nature of your curves allows a different approach. You
could view this as a geometric problem. In this case you check
the signs of the vector cross products

c_n =3D a_n x b_n

where a_n is the vector from a point P1 on the line to the
point under test P, and b_n is a vector from the point P1
to another point P2 on the line. If all the lines are
directed, the signs of all the cross products c_n will
be the same if P lies inside a convex region bounded
by the lines.

Rune

Rune Allnor <allnor@tele.ntnu.no> wrote in message <6a95faa5-b673-4fa8-804d-dd9ae4da773d@b38g2000prf.googlegroups.com>...

>
> The terminology you use would suggest that this is a constrained
> linear programming problem.

I don't think it is. The feasible set in LP is convex polygonal. OP shows a non-convex feasible set in his example.

Bruno