Perpendicular distance from Point to Line

704 views
Skip to first unread message

cathal coffey

unread,
Jun 21, 2011, 9:49:14 AM6/21/11
to sympy
Hello,

I have been using PostGis and Django Geos for a long time now, however
today they both failed me. This prompted me to search and find sympy,
so here I am.

I am trying to find the perpendicular distance from a Line to a Point.
PostGis and Django Geos will only let me calculate the distance from a
Point to a Line Segment i.e. the distance from the closest point on a
Line Segment to a Point.

I have tried the following in sympy but it is not returning what I
expect.

-----------------------------------------------------------------------------------------------------
from sympy.geometry import *

l = Line(p1, p2) # My line is defined using two points.
p = Point(x1, y1) # The point I wish to find the perpendicular
distance to.

# I expected the below to return a Line Segment connecting l and p,
the length of this Line Segment being the distance I want.
s = l.perpendicular_segment(p)
distance = s.length

However I am getting the below output.
In [9]: l.perpendicular_segment(p)
Out[9]: Point(Integer(123), Integer(123))
-----------------------------------------------------------------------------------------------------

What am I doing wrong?

Kind regards,
Cathal

Chris Smith

unread,
Jun 21, 2011, 11:25:40 AM6/21/11
to sy...@googlegroups.com
    >>> help(Line.perpendicular_segment)
    
    Help on method perpendicular_segment in module     sympy.geometry.line:

    perpendicular_segment(self, p) unbound     sympy.geometry.line.Line method
        Returns a new Segment which connects p to a point on this     linear
        entity and is also perpendicular to this line. Returns p     itself
        if p is on this linear entity.
    >>> l=Line(Point(1,2), Point(2,3))
    >>> l.perpendicular_segment(Point(0,0))
    Segment(Point(Rational(-1, 2), Half(1, 2)), Point(Zero,     Zero))
    >>> p=Point(3,4)
    
    >>> p in l
    
    True
    >>> l.perpendicular_segment(p)
    Point(Integer(3), Integer(4))

Looks like your point in on the line.

cathal coffey

unread,
Jun 21, 2011, 11:57:26 AM6/21/11
to sy...@googlegroups.com
Chris thank you for your quick reply.
The original example I posted used test data, below is as far as I
have gotten with real numbers.

from sympy.geometry import *

# Define a Line given two points p1 and p2.
p1 = Point(-694310.42867240659, 7051910.5392115889)
p2 = Point(-694286.161023414, 7051916.4164796043)
l = Line(p1, p2)

# Get the Segment perpendicular to l given p.
p = Point(-694210.80222825997, 7051914.7562929001)


s = l.perpendicular_segment(p)
distance = s.length

This throws the following error
In [61]: s = l.perpendicular_segment(p)
ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (102, 0))

In my attempt to figure out why this error is happening I devised the
below test.

l2 = l.perpendicular_line(p)
l.is_perpendicular(l2)
Out[56]: False

You might ask why my points have such large co-ordinates? They are
Points using the Google Spherical Merchant projection. As I said
originally I have been using PostGis and Django Geos for a long time.
This is the format of data that I work with.

Any ideas?

Kind Regards,
Cathal

Chris Smith

unread,
Jun 21, 2011, 3:23:05 PM6/21/11
to sy...@googlegroups.com
I think we should not warn about floats; we should just convert them in routines such as this. Until then, try change your floats to Rational, e.g.

    >>> def rat(d):
    ...  return S(str(d), rational=True)
    ...
    >>> p1 = Point(rat(-694310.42867240659), rat(7051910.5392115889))
    >>> p2 = Point(rat(-694286.161023414), rat(7051916.4164796043))
    >>> l=Line(p1,p2)
    >>> p = Point(rat(-694210.80222825997), rat(7051914.7562929001))
    >>> s = l.perpendicular_segment(p)
    >>> s.length
    862844452291*12723695727349**(1/2)/159046196591862500
    >>> _.n()
    19.3515546455920

Both S and nsimplify have a rational flag that can be used to change numbers to Rational form. S requires you to convert the number to a string first while nsimplify does not (so I could have used nsimplify in rat, too).

    >>> nsimplify(2.3, rational=1)
    23/10

Aaron Meurer

unread,
Jun 21, 2011, 10:39:03 PM6/21/11
to sy...@googlegroups.com

What is the Google Spherical Merchant projection? I searched for it
on Google, but the only exact match I found was this thread.

Aaron Meurer

>
> Any ideas?
>
> Kind Regards,
> Cathal
>

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to sympy+un...@googlegroups.com.
> For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
>
>

cathal coffey

unread,
Jun 22, 2011, 5:06:10 AM6/22/11
to sy...@googlegroups.com
Chris,

when I execute your code I get the error: TypeError: sympify() got an
unexpected keyword argument 'rational'.

I installed sympy using apt-get install python-sympy and synaptic
reports 0.6.7-1.1 as the installed version.

When I execute ?S I don't see any options for rational.
Call def: S(a, locals=None, convert_xor=True)

The same is true for nsimplify(2.3, rational=1)

from sympy import nsimplify
nsimplify(2.3, rational=1)

TypeError: nsimplify() got an unexpected keyword argument 'rational'

Again executing ?nsimplify mentions no rational param.
Definition: nsimplify(expr, constants=[], tolerance=None, full=False)

Any ideas?
Cathal

from sympy.geometry import *
from sympy import S

def rat(d):
   return S(str(d), rational=True)

p1 = Point(rat(-694310.42867240659), rat(7051910.5392115889))
TypeError: sympify() got an unexpected keyword argument 'rational'

cathal coffey

unread,
Jun 22, 2011, 6:04:00 AM6/22/11
to sy...@googlegroups.com
Chris,

I got it to work, solution is below. Thanks for all your help... sympy rocks!!!.

Cathal

import sympy
from sympy import Rational

# Convert to Rational
def rat(d):
return Rational(str(d))

# Calculate the perpendicular distance from the point (lat, lon) to
the line defined by the points (x1, y1) and (x2, y2).
def perpendicular_distance(x1, y1, x2, y2, lat, lon):
p1 = sympy.geometry.Point(self.rat(x1), self.rat(y1))
p2 = sympy.geometry.Point(self.rat(x2), self.rat(y2))
l = sympy.geometry.Line(p1,p2)
p = sympy.geometry.Point(self.rat(lat), self.rat(lon))

s = l.perpendicular_segment(p)
return s.length.n()

cathal coffey

unread,
Jun 22, 2011, 5:08:51 AM6/22/11
to sy...@googlegroups.com
Aaron,

sorry its not actually a Google thing.

http://en.wikipedia.org/wiki/Mercator_projection

Cathal

Chris Smith

unread,
Jun 22, 2011, 7:57:52 AM6/22/11
to sy...@googlegroups.com
Sometimes you have two sympy installations, e.g. one in site-packages and one somewhere else. Figuring out which one is getting used is the trick. I haven't really spent the time to figure out how to get my IDE to use the one that is in my git directory. I just know what I have to do to use the git one to work. Others may be of more help.

I see that you found Rational(str(...)) to work. That's good, too. If you ever get nsimplify or S to work they are good for situations where the Float is buried in an expression. But your application can easily target the Float so Rational works more directly.

Glad it's working. One thing I like about python is the way you can write an adapting function so easily (as you have done) to deal with a new situation.

Best regards,
 Chris

Aaron S. Meurer

unread,
Jun 22, 2011, 5:22:09 PM6/22/11
to sy...@googlegroups.com
I think he was not using the git version at all. To get those features Chris mentioned, you need to use the development version, or else use the 0.7.0 version that will be coming out in a few days.

Aaron Meurer

Reply all
Reply to author
Forward
0 new messages