Trouve le point le plus proche sur une polyligne d'un point donné

Closest Point on a Polyline

Paul Crisp (Paul.CRISP@syntegra.com) posted ro mapinfo-l on July 3, 2002 some VB code to find the closest point from a given point to a polyline (it could be the “foot” of the shortest perpendicular line from the point to every segment of the polyline, or an end point of the polyline if no perpendicular can be traced that would intersects a segment between its 2 nodes). As the translation to pure MapBasic code is not obvious, I have decided, with Paul's entire support, to offer my solution for inspection and download.

The essence of Paul's code is contained in 3 functions: returnPointOfShortestPerpendicular(), FindClosestPointOnLineSegment() and PointIsBetween(). These functions need also a sub, ErrorHandler, that is not despite its name a standard MB procedure. It is useful to locate the source of an error in one of the functions.

I have added a small interface that shows a way to use these functions and that allows for some experimenting but is not a “production” tool. It is mainly made of “traps” to make sure that valid data is entered. It does not exclude working with “degrees” tables but the results are probably incorrect then.

The way I have used the SelChangedHandler is most primitive but it suffices for that kind of application. I would be most interested to learn of another way to code for successive selections of different objects.

You can download the following code as a MB listing in the ClstPntOnPline.zip file

'

'Demo program showing a way to find the closest point on a polyline to a given point
'
'using the VB code offered by Paul Crisp (Paul.CRISP@syntegra.com) July 3, 2002
'adapted to MapBasic and presented with a simple interface by
'Jacques Paris (jacques@paris-pc-gis.com) July 3, 2002
'
'

Include “mapbasic.def”
Type tPoint
X as float
Y as float
end type

Declare Function ReturnPointOfShortestPerpendicular(clsPoint As tPoint,clsNodes() As tPoint, closestPoint As tPoint) As Logical
Declare Function FindClosestPointOnLineSegment(ptNearbyPoint As tPoint, ptLineSegmentEnd1 As tPoint, ptLineSegmentEnd2 As tPoint, ptClosestPoint As tPoint) As Logical
Declare Function PointIsBetween(pt1 As tPoint, pt2 As tPoint, ptInBetween As tPoint) As Logical
Declare sub errorhandler(byval a as string)
Declare sub SelChangedHandler
Declare sub main

dim opline, opoint as object
dim clsPoint as tpoint
dim clsNodes() as tpoint
dim closestPoint as tpoint
dim my_table as string
dim tcall as smallint

'

sub main
'

if numwindows()=0 then
note“A window must be open”
exit sub
end if

if int(windowinfo(frontwindow(),3))<>1 then
note “Fontwindow must be a mapper”
exit sub
end if

note “Warning: this demo program will not return correct results if the ”+chr$(10)+“base table uses an 'earth, not projected' (in degrees) coordsys”

tcall=1

run menu command 304

tcall=2

note “Select a PolyLine”

end sub

'

sub SelChangedHandler
'

dim nnodes, i as integer

do case tcall

case 1

exit sub

case 2

'A polyline must be selected

opline=selection.obj

if int(objectinfo(opline,1))<>4 then
note “Selected object is not a Polyline”
exit sub
end if

if int(objectinfo(opline,21))>1 then
Note “Works only on single section polylines”
exit sub
end if

my_table=selectioninfo(1)

set coordsys table my_table

nnodes=objectinfo(opline,20)

redim clsNodes(nnodes)

for i=1 to nnodes
clsNodes(i).X=objectnodeX(opline,1,i)
clsNodes(i).Y=objectnodeY(opline,1,i)
next

tcall=3

run menu command 304
exit sub

case 3

tcall=4

note “Select a Point”

exit sub

case 4

'A point must be selected

opoint=selection.obj

if int(objectinfo(opoint,1))<>5 then

note “Selected object is not a Point”

exit sub

end if

if selectioninfo(1)<>my_table then
note “PolyLine and Point must be in the same table”
exit sub
end if

print selectioninfo(1)

clsPoint.X=centroidX(opoint)
clsPoint.Y=centroidY(opoint)

'Calling the minimum distance function

if not ReturnPointOfShortestPerpendicular(clsPoint, clsNodes(), closestPoint) then
note “Function failed to find closest point”
exit sub
end if

insert into my_table (obj) Values ( createpoint(closestPoint.X, closestPoint.Y) )

end case

end program

end sub

'

sub errorhandler (byval a as string)
'

note “Error encountered in ”+a+chr$(10)+chr$(10)+err()+“ ”+error$()

end sub

'

Function ReturnPointOfShortestPerpendicular(clsPoint As tPoint,
clsNodes() As tPoint, closestPoint As tPoint) As Logical
'

' Purpose: Routine to return the end point of the perpendicular
' on to nearest line segment from a point
' Inputs: Point, points collection from line
'
' Returns: point
'*

OnError GoTo labelError

Dim i As Integer
Dim ptTemp As tPoint
Dim intCount As Integer
Dim iClosestPoint As Integer
Dim dblShortestDist As Float
Dim dblTempDist As Float
Dim ptLineSegmentEnd1, ptLineSegmentEnd2 as tPoint

intCount = ubound(clsNodes)

'Treat line feature line as a series of lines and get the perpendicular onto each
'take the shortest as the return value
'Initialise the values from the first line segment then compare each of the rest
'Give dblShortestDist a value just in case

closestPoint.X=clsNodes(1).X
closestPoint.Y=clsNodes(1).Y

dblShortestDist = Sqr¹⁾
dblReturnY = (dblM * dblReturnX) + dblC

End If

ptClosestPoint.X = dblReturnX
ptClosestPoint.y = dblReturnY

'Check this closest point is actually between the other two

If Not PointIsBetween(ptLineSegmentEnd1, ptLineSegmentEnd2, ptClosestPoint) Then

'The nearest point on that segment is actually one of the vertices

dblTempDist1 = sqr ((ptClosestPoint.X - ptLineSegmentEnd1.X) ^ 2 +(ptClosestPoint.Y - ptLineSegmentEnd1.Y) ^ 2)

dblTempDist2 = sqr ((ptClosestPoint.X - ptLineSegmentEnd2.X) ^ 2 +(ptClosestPoint.Y - ptLineSegmentEnd2.Y) ^ 2)

If dblTempDist1 < dblTempDist2 Then
ptClosestPoint.X = ptLineSegmentEnd1.X
ptClosestPoint.Y = ptLineSegmentEnd1.Y

Else
ptClosestPoint.X = ptLineSegmentEnd2.X
ptClosestPoint.Y = ptLineSegmentEnd2.Y

End If 'Vertex 1 is close

End If 'Returned point is not on line segment

FindClosestPointOnLineSegment = True

Exit Function

labelError:

FindClosestPointOnLineSegment = False

Call ErrorHandler(“FindClosestPointOnLineSegment”)

End Function

'

Function PointIsBetween(pt1 As tPoint, pt2 As tPoint, ptInBetween As tPoint) As Logical
'

' Purpose: Routine to check whether a point falls between two others
' Inputs: Three Points
'
' Returns: Logical
'*

OnError GoTo labelError

Dim dblMapX As Float
Dim dblMapX2 As Float
Dim dblMapY As Float
Dim dblMapY2 As Float

PointIsBetween = False

'Sort the values so x2, y2 are the larger

If pt1.X > pt2.X Then
dblMapX = pt2.X
dblMapX2 = pt1.X

Else
dblMapX = pt1.X
dblMapX2 = pt2.X

End If

If pt1.Y > pt2.Y Then
dblMapY = pt2.Y
dblMapY2 = pt1.Y

Else
dblMapY = pt1.Y
dblMapY2 = pt2.Y

End If

'Use a 10cm tolerance and see if the point falls between the other two

If (ptInBetween.X + 0.1 >= dblMapX) And (ptInBetween.X - 0.1 ⇐ dblMapX2) And
(ptInBetween.Y + 0.1 >= dblMapY) And (ptInBetween.Y - 0.1 ⇐ dblMapY2) Then
PointIsBetween = True
End If

Exit Function

labelError:

PointIsBetween = False

Call ErrorHandler(“GenericUtils: PointIsBetween”)

End Function

'================================================================================= [/code]