Find which side of a line a point is on

[ChrisMathers]

How can I find which side of a line a point which isnt on the line falls. I have a road and need to see which side of the line a point is so I know which field to pull a value from between R_ESN and L_ESN.



Also, long time no post. The forums are looking snazzy these days!

[KimOllivier · Accepted solution]

This spatial query is a glaring gap in ArcGIS tools. The function is available in Avenue (for ArcView 3.x)

But we have Python integrated in ArcGIS and even NUMPY, so anything is possible now!
What you have to do is use the NEAR tool with the option of adding the X,Y coordinates of the intersection with the line.

With this intersection point and the original point you can calculate the side of the polyline using a vector cross product 🙂 You use the sign of the vector product for an indication of left/right. To generate a simplified vector of the polyline for this you will need the lastPoint of the polyline as well, which is easy enough to get from the shape properties. Having the lastPoint gives you the direction of the polyline.

To see how this works, just remember the Right Hand Rule. Hold out your right hand and spread your thumb and first two fingers perpendicular as a 3D axis. If the first finger is the polyline, the second finger the point, then the thumb pointing up (+ve) indicates the point is on the left. Turn your hand upside down so the second finger is on the right, then the thumb pointing down is now negative (-ve).

Numpy has a cross product function (to save us the headache of matrix arithmetic.)

import numpy as np a = np.array([1,0,0])   b = np.array([0,1,0])   print np.cross(a,b)

There are a couple of assumptions that are critical for this vector algebra to work as expected.
No multipart shapes
No curves

I agree with another poster that these are evil additions to the data model that break most analysis tools.

# LeftRight.py # Vector Cross Product # to find which side a point is of a polyline # Kim Ollivier 18 Feb 2012  import numpy as np  # Get the points from original point, Near and shape.lastPoint # eg xPnt = (1,20) intPnt = (11,9) lastPnt = (5,20) # make two 3D vectors relative to the intersection point # for a local origin at the intersection point a = np.array([xPnt[0]-intPnt[0],xPnt[1]-intPnt[1],0])   b = np.array([lastPnt[0]-intPnt[0],lastPnt[1]-intPnt[1],0])   c = np.cross(a,b) if c[2] > 0:     print "Side is left" else:     print "Side is right" print c

I haven't written a general tool for this, but it would be a good candidate for the Resource Centre (now deprecated). Maybe ArcScripts could be reinstated?
Competition anyone?

[RichardFairhurst]

I know this is the Python forum, but I wanted to mention that in ArcObjects that the QueryPointAndDistance method of the ICurve interface reports not only the Near Point, but also reports the distance along the line where the near point fell (as either a percentage of the line or a linear distance), the distance of offset from the line, and the side of the line the point falls on (left/right) in one operation.  It is very efficient and cool, and potentially worth the trouble of learning how to use ArcObjects within Python if you need to do a lot of linear analysis.  It works with the evil true curves also and polylines of any shape.  I have used it with Address Range validation and assignment and other Linear Referencing operations.

[ChrisMathers]

Thanks for the help Kimo. Im going to have to look into ArcObjects eventually. I dont even have time to get on here anymore much less learn a new language!

[DuncanHornby]

I just want to compliment Richards answer. When a python script is run from within ArcGIS as a scripting tool it's possible to call arcobjects within the calculate field tool. I've used this approach sometimes when I wanted to access the fine grade geometry tools that ArcObjects has. The QueryPointAndDistance method on IPolyline is certainly a very useful bit of wizardry I use it regularly!

[KimOllivier]

I just want to compliment Richards answer. When a python script is run from within ArcGIS as a scripting tool it's possible to call arcobjects within the calculate field tool. I've used this approach sometimes when I wanted to access the fine grade geometry tools that ArcObjects has. The QueryPointAndDistance method on IPolyline is certainly a very useful bit of wizardry I use it regularly!

I did manage to get an example ArcObjects call running for 9.4 presented by Mark Cederholm, but everything has changed at 10.
Could you post a snippet to give us a hint? It would be a lot easier to use an existing function.
What a pity that the Near tool (and many others) is only a subset of ArcObjects functionality.

Update
Mark has updated his demo for version 10 at http://www.pierssen.com/arcgis10/python.htm

[DuncanHornby]

Kimo,

It looks like I am going to have to apologies here. What I said is correct for ArcGIS 9.3 but ESRI has disappointingly removed this functionality from ArcGIS 10. I was not aware that they had killed this off until I just tried it. It's discussed here and it seems to be something to do with a limited version of VBA.

So unless you can figure out how to access arcobjects from python (here is a presentation that discusses this, you may have seen it already) you'll need to get to grips with ArcObjects in VB or c#.

It's a real shame that ESRI have removed this ability. From my developers point of view about the only thing ESRI have done good with the release of version 10 is the ESRIAddIn, this makes deploying buttons and toolbars developed in .Net a whole lot easier.

So sorry for getting your hopes up, thats "progress" for you...

Duncan

[KimOllivier]

It looks like Python wins over ArcObjects for development speed. 🙂 Which is exactly why I use it.
I always mean to get back into C++ some day, but the huge 8.02 Object Model (still up on the wall) puts me off.

Since the cursor speed has been improved at 10.1 a script might be fast enough to be practical. The limitations I mentioned are only because of the shortcut to get the lastPoint for the line direction. You can eliminate these cases (if present) by exploding to a temporary layer.
If the next vertex (instead of the last) was retrieved then it would handle multi-parts and another tweak is possible for the unusual case of imbedded parametric curves. Or just explode them too. The answer could be put back on the original featureclass.

I will have to finish the snippet above to generalise it.

[KimOllivier]

I have found a geoprocessing tool which returns the side : LocateFeaturesAlongRoutes has a postive or negative distance to indicate the left or right side of a route.

There is an optional flag to measure the side from the route direction or arc direction.

So all you have to do is build a route system with an id of each arc for the route_id and then use it to find the side.