Offset points in Mapserver and PostGIS

Recently, I had to find a way to offset point features that have been created using dynamic segmentation in PostGIS and display them through Mapserver.

Obviously, the points would have to be offset perpendicularly to the line. Mapserver includes an offset function through the STYLE object which works well for lines. By defining the offset as OFFSET [x] –99 the line geometry will be rendered as shifted by x SIZEUNITS perpendicular to the original line geometry. But this technique wouldn’t work for points.

So the idea was to do the offset server-side but WITHOUT changing the actual geometry as the offset was for display purposes only.

To do this I used the st_line_interpolate_point_with_offset function found here and then created a simple wrapper around it:

CREATE OR REPLACE FUNCTION offset_point(
	pnt_geometry geometry,
	line_geometry geometry,
	offset_val double precision)
    RETURNS geometry
    LANGUAGE 'plpgsql'

    COST 100
    VOLATILE 
AS $BODY$

 DECLARE
 loc double precision;
 offset_geom geometry;
  BEGIN
  loc:=ST_LineLocatePoint(line_geometry, pnt_geometry);
	--st_line_interpolate_point taken from: https://trac.osgeo.org/postgis/wiki/UsersWikiExamplesInterpolateWithOffset	
	offset_geom:=st_line_interpolate_point(line_geometry, loc, offset_val);
	--raise notice 'offset geom: %1', ST_AsText(offset_geom);
	return offset_geom;
  END;

$BODY$;

Then, in the mapfile, I set the DATA statement to be:

CONNECTIONTYPE POSTGIS
    DATA "mh_geom from (select mh_gid, offset_point(mh_geom, 'water_pipe','wp_geom',%offsetval%) as
      mh_geom from manhole g) as
      subquery using unique mh_gid using srid=2100"

I also added a VALIDATION block to ensure the offsetval substitution variable was a number:

 
VALIDATION 
    "offsetval"	"-?[0-9]{0,10}"
END

I could then fire a mapserver/WMS request appending the offsetval as a parameter:


http://localhost/mapserver/mapserv?map=test_offsets.map... &offsetval=10

With an offset of 10:

image

Or, with an offset of –10:

image

Hope this helps someone. Happy coding!

Advertisements

Look Ma! No GIS!

Before you light’ em torches, let me tell you: Not really, no. But the following are a couple of posts I wrote on LinkedIn explaining why a “GIS-only” approach is not ideal when it comes to Linear Referencing and Transportation networks. So, for your reading pleasure:

When GIS is just not enough – Part 1

When GIS is just not enough – Part 2

And the Tech Talk video

 

More linear referencing functions for SQL Server 2008

Following Geographika’s LRS add-ons for  SQL Spatial Tools, I have added a new function to locate the closest point along a geometry LineString from a given point. The function called  LineLocatePoint, is the equivalent of PostGIS ST_line_locate_point function. It will return a number between 0 and 1, representing the location of the closest point along a geometry LineString to the input point. In other words, “projecting” the input point to the line. You can use this function in conjunction with the existing LocateAlongGeom function to retrieve the actual point geometry. This is a schematic of how this function works:

image

Thanks to Seth, a.k.a. Geographika, I have uploaded the updated Spatial Tools code to the BitBucket repository (https://bitbucket.org/geographika/sql-server-spatial-tools) along with the register and unregister scripts.

To test the function use the following script:

select dbo.LineLocatePoint(geometry::STGeomFromText(
                  'LINESTRING (411714.523521 4493882.696417, 
                 411737.898344 4493862.196547,
                411753.991966 4493838.696703, 411759.616915 4493820.696827,
                411753.179452 4493807.696923, 411740.179536 4493794.697023)',0) 
             ,geometry::STGeomFromText('POINT (411740.179536 4493794.697023)',0))

One important caveat: If for some reason the point cannot be projected onto the line, the function will return a nice, clean 0.5 value i.e.. the mid-section. If you get 0.5, you are either very lucky, OR (more likely) there was some kind of error or division by 0 due to the way the line and point are located in relation to each other. I did this on purpose since I use this function for some address geocoding routines where I always want a value returned.