[postgis-users] st_distance accuracy?
Sufficool, Stanley
ssufficool at sbcounty.gov
Tue Oct 20 13:46:32 PDT 2009
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;
More information about the postgis-users
mailing list