<DIV>Thanks for posting this. Let me make sure that I have this correct - as we may need something similar to this in an upcoming project. It sounds like this function essentially "snaps" a point to the nearest line. Correct? If so, does it snap the point to the nearest VERTEX on the line? Or will it place the point at the appropriate location between two vertex's on the line? I just ask because an Avenue script ArcView 3.2 would place the point at the nearest vertex and therefore not do what we need. </DIV>
<DIV> </DIV>
<DIV>Thanks again for sharing!</DIV>
<DIV> </DIV>
<DIV>AS<BR><BR><B><I>Stephen Woodbridge <woodbri@swoodbridge.com></I></B> wrote:</DIV>
<BLOCKQUOTE class=replbq style="PADDING-LEFT: 5px; MARGIN-LEFT: 5px; BORDER-LEFT: #1010ff 2px solid">Hi all,<BR><BR>Someone on the list asked if they could get the point location of where <BR>on a line object the shortest distance to that line was from a point. <BR>The distance(line, point) will tell you the shortest distance, but it <BR>does not project the point onto the line at the point where the shortest <BR>distance is located.<BR><BR>Then I needed the same thing for a project at Where2GetIt.com and it <BR>seemed like this functionality should be a part of PostGIS. So attached <BR>is the code I wrote for our needs. There are some examples of its use. <BR>It has an assumption that you are working in decimal degrees, but it <BR>could be modified easily to work in any space. I will leave that up to <BR>others.<BR><BR>The algorithm is very straight forward. Given a linestring and a point, <BR>project the point onto each segment of the linestring and save the <BR>projection of the
closest point. Then compute the distance from the <BR>start of the linestring to the projection point and return an array of <BR>pointx, pointy, distance.<BR><BR>Paul, please feel free to add this or something like it to PostGIS.<BR><BR>Enjoy, I would be open to any feedback or improvements any might suggest.<BR><BR>-Steve<BR>CREATE OR REPLACE FUNCTION x_y_offset_from_line_pnt (geometry, geometry) RETURNS float[] AS '<BR>/*<BR>* {x, y, offset} = x_y_offset_from_line_pnt(LINESTRING, POINT($lon $lat))<BR>*/<BR>DECLARE<BR>line ALIAS FOR $1;<BR>pnt ALIAS FOR $2;<BR>nline geometry;<BR>ipnt geometry;<BR>bestp geometry;<BR>pa geometry;<BR>pb geometry;<BR>npts integer;<BR>besti integer;<BR>seglen float;<BR>dist float;<BR>best float;<BR>offset float;<BR>U float;<BR>--<BR>BEGIN<BR>IF NOT (GeometryType(line) = ''LINESTRING'' ||<BR>GeometryType(line) = ''MULTILINESTRING'' ) AND<BR>GeometryType(pnt) != ''POINT'' THEN<BR>RETURN NULL;<BR>END IF;<BR>/*<BR>* iterate over the linestring segments and
find the segment and <BR>* point on the linestring that are closest to the pnt<BR>*/<BR>npts := npoints(line);<BR>best := 360.0;<BR>FOR i IN 1 .. npts-1 LOOP<BR>pa := PointN(line,i);<BR>pb := PointN(line,i+1);<BR>seglen := distance(pa, pb);<BR>IF (seglen > 0.0) THEN<BR>U := ( (X(pnt)-X(pa)) * (X(pb)-X(pa)) +<BR>(Y(pnt)-Y(pa)) * (Y(pb)-Y(pb)) ) /<BR>(seglen*seglen);<BR>IF U < 0.0 THEN<BR>ipnt := pa;<BR>ELSIF U > 1.0 THEN<BR>ipnt := pb;<BR>ELSE<BR>ipnt := MakePoint(<BR>X(pa) + U*(X(pb) - X(pa)),<BR>Y(pa) + U*(Y(pb) - Y(pa)) );<BR>END IF;<BR>dist := distance(pnt, ipnt);<BR>IF dist < best THEN<BR>best := dist;<BR>besti := i;<BR>bestp := ipnt;<BR>END IF;<BR>END IF;<BR>END LOOP;<BR>/*<BR>* we should now have the best point so calcuate the offset<BR>* and return the results.<BR>*/<BR>IF besti = 1 THEN<BR>nline = MakeLine(PointN(line,1), bestp);<BR>ELSE<BR>nline = MakeLine(PointN(line,1), PointN(line,2));<BR>IF besti > 2 THEN<BR>FOR i IN 3 .. besti LOOP<BR>pa :=
PointN(line,i);<BR>nline := AddPoint(nline, pa, -1);<BR>END LOOP;<BR>END IF;<BR>nline := AddPoint(nline, bestp, -1);<BR>END IF;<BR>offset := length_spheroid(nline,''SPHEROID["GRS_1980", 6378137, 298.257222101]'');<BR>RETURN ARRAY[X(bestp), Y(bestp), offset];<BR>END;<BR>' LANGUAGE plpgsql;<BR><BR>-- here are some test calls to the function<BR><BR>select x_y_offset_from_line_pnt(GeomFromText('LINESTRING(0 0,1 0)',-1), GeomFromText('POINT(0.5 0.1)', -1));<BR><BR>select x_y_offset_from_line_pnt(GeomFromText('LINESTRING(0 0,0.1 0,0.3 .1,.4 .15,.6 .2)',-1), GeomFromText('POINT(0.5 0.1)', -1));<BR><BR>select link_id, the_geom, min(distance(the_geom, GeomFromText('POINT(-71.38900904888 42.619402684661)', -1))) as dist<BR>from <BR>rgeo <BR>where<BR>expand(GeomFromText('POINT(-71.38900904888 42.619402684661)', -1), 0.0013) &&<BR>the_geom<BR>group by the_geom order by dist limit 1;<BR><BR>select (x_y_offset_from_line_pnt(<BR>(select the_geom from rgeo where link_id=$segment_id),
<BR>GeomFromText('POINT(-71.38900904888 42.619402684661)', -1))<BR>)[3] as offset;<BR>_______________________________________________<BR>postgis-users mailing list<BR>postgis-users@postgis.refractions.net<BR>http://postgis.refractions.net/mailman/listinfo/postgis-users<BR></BLOCKQUOTE><p>
<hr size=1>Yahoo! Sports<br>
<a href="http://pa.yahoo.com/*http://us.rd.yahoo.com/evt=33539/*http://football.fantasysports.yahoo.com?ovchn=YAH&ovcpn=Integration&ovcrn=Mail+footer&ovrfd=YAH&ovtac=AD ">Rekindle the Rivalries. Sign up for Fantasy Football</a>