[postgis-users] Degrees of a linestring segment

[Bruce Rindahl]

Flavio
This will not work because line_interpolate_point and line_locate_point do
not give the vertex locations, only points and locations along the line.
Say you have a linestring with three segments (and four points).  If you use
line_locate_point(the_line,my_point) where my_line is your line and my_point
is your point that is close to the middle of the line, the result will be
0.5 or half way along the line.  It will not tell you which vertex it is
closest to.  It you then use line_interpolate_point(my_line,0.5) you will
get a new point halfway along your line or in the middle of the middle
segment.
The good news is it can be done with a new function called
line_locate_vertex(linestring,point) that will give you the closest vertex
to the point.  I wrote this to split polygons into four sides for section
lines recently.  The function is:
--------------------------------------------------------------------------

-- Function: line_locate_vertex(geometry, geometry)

-- DROP FUNCTION line_locate_vertex(geometry, geometry);

CREATE OR REPLACE FUNCTION line_locate_vertex(geometry, geometry)
  RETURNS integer AS
$BODY$
  DECLARE     
	dist double precision;
	n_dist double precision;
	i integer;
	p integer;
  BEGIN
	dist = 9999999999999;
     
   IF (GeometryType($1) != 'LINESTRING') THEN
	raise notice 'First Geometry type must be a LINESTRING';
	RETURN -1;
   END IF;
   IF (GeometryType($2) != 'POINT') THEN
	raise notice 'Second Geometry type must be a POINT';
	RETURN -1;
   END IF;

	FOR i IN 1 .. numPoints($1)
	LOOP
		n_dist = distance($2,pointN($1,i));
		IF (n_dist < dist) THEN
			p = i;
			dist = n_dist;
		END IF;
	END LOOP;

    RETURN p;
  
  END;
$BODY$
  LANGUAGE 'plpgsql' VOLATILE;
ALTER FUNCTION line_locate_vertex(geometry, geometry) OWNER TO postgres;

-------------------------------------------------------------------------

Using this function, the result of line_locate_vertex(my_line,my_point) will
give the index of the vertex closest to your point.  You then have to figure
out if you need the next vertex or the previous vertex for your line
segment.  This can be done using line_locate_point.  First get the location
your closest vertex along the line:

	index = line_locate_vertex(my_line,my_point)

Then find the location of that vertex along the line:

	indx_loc = line_locate_point(pointN(my_line,index))

Then find the location of your point along the line:

	pnt_loc = line_locate_point(my_line,my_point)

If pn_loc is less than indx_loc your line segment goes from
pointN(my_line,index-1) to pointN(my_line,index).

If pn_loc is greater than indx_loc your line segment goes from
pointN(my_line,index) to pointN(my_line,index+1).

Hope this helps!

Bruce Rindahl

 

-----Original Message-----
From: postgis-users-bounces at postgis.refractions.net
[mailto:postgis-users-bounces at postgis.refractions.net] On Behalf Of Flavio
Perri
Sent: Friday, April 27, 2007 3:26 AM
To: postgis-users at postgis.refractions.net
Subject: [postgis-users] Degrees of a linestring segment


Hi all, i'm trying to calculate the degrees of a linestring, but i need the
degree of the closest segment that contains a point. I explain better.

I have a point near a line. I use the line_interpolate_point() and
line_locate_point() to locate the point on the linestring. Then I make two
segments with Makeline from the linestring:
1) StartPoint - interpolated point
2) (Interpolated point+1) - EndPoint

With PointN() i get the (n-1)th point of (1) and the first point of (2) and
finally i can call degrees(azimuth(P1,P2)) to retrieve the degrees.

Is this way ok? The fact is that sometime the return value is NULL because
the two points are the same.

This is the query i used:

      SELECT
degrees(azimuth(PointN(MakeLine(StartPoint(l.the_geom),line_interpolate_poin
t(l.the_geom,
 line_locate_point(l.the_geom,
a.point_geom))),(numPoints(MakeLine(StartPoint(l.the_geom),line_interpolate_
point(l.the_geom,
 line_locate_point(l.the_geom,
a.point_geom))))-1)),PointN(MakeLine(line_interpolate_point(l.the_geom,
line_locate_point(l.the_geom, a.point_geom)), EndPoint(l.the_geom)),1))) as
degrees from pointTable a, LineTable l

Bye,
Flavio

_______________________________________________
postgis-users mailing list
postgis-users at postgis.refractions.net
http://postgis.refractions.net/mailman/listinfo/postgis-users