[postgis-users] st_distance accuracy?

[Sufficool, Stanley]

Following the KISS principle helps sometimes:

CREATE OR REPLACE FUNCTION st_line_interpolate_point(m_geom geometry, m_percent double precision, m_offset double precision)
  RETURNS geometry AS
$BODY$

declare m_seg geometry;
declare p1 geometry;
declare p2 geometry;

BEGIN

m_seg := st_line_substring(m_geom, 0, m_percent);
p1:=st_pointn(m_seg, st_npoints(m_seg)-1);
p2:=st_pointn(m_seg, st_npoints(m_seg));

/* get the delta of the line segment and offset the interpolated point */
/* Rotate deltas by 90 degrees (invert dx/dy ) */
return st_translate(
	st_line_interpolate_point(m_geom, m_percent), 
	((st_y(p2) - st_y(p1)) / st_distance(p1,p2)) * m_offset,
	((st_x(p2) - st_x(p1)) / st_distance(p1,p2)) * (-m_offset)
);

END;

$BODY$
  LANGUAGE 'plpgsql' VOLATILE
  COST 100;
ALTER FUNCTION st_line_interpolate_point(geometry, double precision, double precision) OWNER TO postgres;



> -----Original Message-----
> From: postgis-users-bounces at postgis.refractions.net 
> [mailto:postgis-users-bounces at postgis.refractions.net] On 
> Behalf Of Nicklas Avén
> Sent: Tuesday, October 20, 2009 2:29 PM
> To: PostGIS Users Discussion
> Subject: Re: [postgis-users] st_distance accuracy?
> 
> 
> Sorry wrong thinking from me in last post. But look at it 
> again. You are not measuring along the line but strait from 
> first vertex to each of the others. /nicklas
> 
> -- originalmedd. --
> Ämne:	Re: [postgis-users] st_distance accuracy?
> Från:	nicklas.aven at jordogskog.no
> Datum:		2009-10-20 23:14
> 
> To find the last two vertexes why not use: 
> p1:=st_pointn(m_geom,ST_NPoints(m_geom)-1);
> p2:=st_pointn(m_geom,ST_NPoints(m_geom));
> /Nicklas
> 
> 
> 
> 2009-10-20 nicklas.aven at jordogskog.no wrote:
> 
> I'm not sure I understand what you are trying to do. But as I 
> read it you are comparing the length of the line with the 
> >straight distance from the first to one and each of the 
> following vertexes. > >Is that really what you want?> 
> >Nothing says that the last vertex is the one most far away 
> from the first and the only case when that distance is the 
> same as the length of the line is if the line is straight.
> >> >Or am I missunderstanding you?> >/Nicklas> > >
> > 2009-10-20 Sufficool Stanley wrote:
> >
> > I am attempting to get the last vertex index of a line for a given
> > >measure using the following:
> > >
> > >idx := max(np) from generate_series(1,st_npoints(m_geom)) 
> as np where 
> > >st_distance(st_pointn(m_geom,1), st_pointn(m_geom,np)) <= 
> > >m_distance;
> > >
> > >However, some of my queries are returning the last point on the 
> > >linestring when "m_distance" is less than the st_length() of the 
> > >given geometry.
> > >
> > >Is this an issue with st_distance? I don't get the same issue when 
> > >using Lat/Lon and st_distance_spheroid().
> > >
> > >select postgis_version()
> > >1.4 USE_GEOS=1 USE_PROJ=1 USE_STATS=1
> > >
> > >select version()
> > >PostgreSQL 8.4.0, compiled by Visual C++ build 1400, 32-bit
> > >
> > >
> > >FULL FUNCTION TEXT
> > >------------------------
> > >CREATE OR REPLACE FUNCTION 
> st_line_interpolate_point(m_geom geometry, 
> > >m_percent double precision, m_offset double precision)  RETURNS 
> > >geometry AS $BODY$
> > >
> > >declare idx int;
> > >declare m_distance double precision;
> > >declare p1 geometry;
> > >declare p2 geometry;
> > >
> > >BEGIN
> > >
> > >IF (m_distance > 1 or m_distance < 0) THEN
> > > raise exception 'Percentage must be between 0 and 1';
> > >END IF;
> > >
> > >/* Convert percentage to actual distance measure */ m_distance := 
> > >st_length(m_geom) * m_percent;
> > >
> > >/* Get the last point number on the line that is less than the 
> > >requested distance */
> > >
> > >/* For Feet (NAD1983) */
> > >idx := max(np) from generate_series(1,st_npoints(m_geom)) 
> as np where 
> > >st_distance(st_pointn(m_geom,1), st_pointn(m_geom,np)) <= 
> > >m_distance;
> > >
> > >/* For Lat/Long (TIGER Lines) */
> > >--idx := max(np) from generate_series(1,st_npoints(m_geom)) as np 
> > >where st_distance_sphere(st_pointn(m_geom,1), 
> st_pointn(m_geom,np)) 
> > >< m_distance;
> > >
> > >/* Past the end of the line segment return exception (this should 
> > >never
> > >happen) */
> > >if idx = st_npoints(m_geom) then 
> > > raise exception 'st_line_interpolate_point(geom, dbl, 
> dbl): Attempt to
> > >seek beyond end of line.';
> > >end if;
> > >
> > >/* Grab the 2 points of the matching line segment */ 
> > >p1:=st_pointn(m_geom,idx); p2:=st_pointn(m_geom,idx + 1);
> > >
> > >/* get the delta of the line segment and offset the interpolated 
> > >point */
> > >/* Rotate deltas by 90 degrees (invert dx/dy ) */
> > >/* spheroid projections will need the offset to vary based 
> on lat/lon 
> > >*/ return st_translate(  line_interpolate_point(m_geom, m_percent),
> > > ((st_y(p2) - st_y(p1)) / st_distance(p1,p2)) * m_offset,
> > > ((st_x(p2) - st_x(p1)) / st_distance(p1,p2)) * (-m_offset)
> > >);
> > >
> > >END;
> > >
> > >$BODY$
> > > LANGUAGE 'plpgsql' VOLATILE
> > > COST 100;
> > >ALTER FUNCTION st_line_interpolate_point(geometry, double 
> precision, 
> > >double precision) OWNER TO postgres;
> > >
> > >_______________________________________________
> > >postgis-users mailing list postgis-users at postgis.refractions.net
> > >postgis.refractions.net/mailman/listinfo/postgis-users
> > >
> > >
> 
> _______________________________________________
> postgis-users mailing list postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users
> 
> 
> _______________________________________________
> postgis-users mailing list postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users
>