[postgis-users] Getting dimensions of polygon sides

Obe, Regina robe.dnd at cityofboston.gov
Mon Jun 16 06:11:58 PDT 2008


Bruce,

Thanks for the suggestion.  After some struggling with my highschool
trig and some help,  came up with this function.

CREATE OR REPLACE FUNCTION upgis_lineshift(centerline geometry, dist
double precision)
RETURNS geometry AS
$$
DECLARE 
	delx float;
	dely float;
	x0 float;
	y0 float;
	x1 float;
	y1 float;
	az float;
	dir integer;
	line geometry;
BEGIN
	az := ST_Azimuth(ST_StartPoint(centerline),
ST_EndPoint(centerline));
	dir := CASE WHEN az < pi() THEN -1 ELSE 1 END;
	delx := ABS(COS(az)) * dist * dir;
	dely := ABS(SIN(az)) * dist * dir;

	x0 := ST_X(ST_StartPoint(centerline));
	y0 := ST_Y(ST_StartPoint(centerline));
	x1 := ST_X(ST_EndPoint(centerline));
	y1 := ST_Y(ST_EndPoint(centerline));

	IF az > pi()/2 AND az < pi() OR az > 3 * pi()/2 THEN
		line := ST_Translate(centerline, delx, dely) ;
	ELSE
		line := ST_Translate(centerline, -delx, dely);
	END IF;

	RETURN line;
END;
$$

LANGUAGE 'plpgsql' IMMUTABLE;
COMMENT ON FUNCTION upgis_lineshift(geometry, double precision) IS
'Takes a 2D line string and shifts it dist units along 
the perpendicular defined by the straight line between the start and end
point 
Convention: (right is positive and left is negative.  right being
defined as to right of observer
standing at start point and looking down the end point)';

Using it like so
  INSERT INTO assessing.parcdimstime(pid, pid_year, the_geom)
	SELECT b.pid, b.pid_year, ST_MakeLine(ST_PointN(the_boundary,n),
ST_PointN(the_boundary, n + 1)) As the_side
	FROM (SELECT pid, pid_year, the_boundary,
ST_NumPoints(the_boundary) As nump
		FROM (SELECT ap.parcelid as pid, 2009 as pid_year,
ST_Simplify(ST_Boundary(ST_ForceRHR(ap.the_geom)), 1) As
the_boundary
			FROM dnd.rems_survey ap  ) p) As b
				CROSS JOIN generate_series(1, 500) n
	WHERE n < b.nump;

UPDATE assessing.parcdimstime SET sd_length = CAST(ST_Length(the_geom)
As numeric(8,2)),the_linepos = upgis_lineshift(the_geom,3);


That seems to fall in the polygon as you said.  Its strange I never
thought about lines having direction.  Then I snapped the deed
annotations to my lines with another piece of sql code.  See snapshot.

The dashed lines is what I get from applying the shift function to
create a line 3 feet into the parcel.  The blank text is what I
calculated the edge to be.  The red is what the legal deeds say it
should be.

Which brings me to a bit of a sad and happy note.  Thanks Frank
Warmerdam by the way.

I had originally thought that the annotations were generated by some
sophisticated ESRI software.  Turns out they were hand input from paper
legal deeds which were stored in intergraph and then migrated to ArcGIS
and they don't match with my calculated measures.  I asked why these
numbers are off.  Thats a long story involving paper maps scanned in
and warped because paper warps, alignment with aerials etc.  So short
answer the lengths generated from the geometries can not be fully
trusted.  Though they have been trying for years to have them agree.
Seems like an interesting operations research problem.

The reason I went thru this exercercise was that I asked 2 people using
ArcMap to export it into a regular point/line format that the ArcPad
version we had could deal with.  Neither of them was able to and they
even asked ESRI support for help.  One said he was able to in 9.1 but
for some reason couldn't in 9.2  After discovering that ArcPad actually
uses that ESRI MDB geodatabase format but can't read the annotation
layer for some reason and then remembering that Ogr2Ogr can read ESRI
geodatabase format,  I gave Ogr2Ogr  a try and it with one line - dumped
out the annotations as a nice shape and also a nice PostgreSQL format
that seemed pretty trivial to work with.  It even had an angle field in
it for rotation of the text.

The lines I used were pretty simple
ogr2ogr -f "ESRI Shapefile" parcdims2008.shp myGeodatabase.mdb
dimensions1
ogr2ogr -f "PostgreSQL" PG:"host=localhost user=someuser dbname=somedb"
mygeodatabase.mdb dimensions1


Can someone who has ArcMap 9.2 verify for me if it can or can not export
annotations in ESRI Shapefile format.  The idea that it can't is totally
baffling to me.  So I'm guessing its some misunderstanding. I still have
to install my ArcMap 9.2 and become an addict, but I think it would take
me a while to get up to speed since all the menu options are too much
for my little brain to handle.

Thanks,
Regina


 

-----Original Message-----
From: postgis-users-bounces at postgis.refractions.net
[mailto:postgis-users-bounces at postgis.refractions.net] On Behalf Of
Bruce Rindahl
Sent: Thursday, June 12, 2008 1:04 PM
To: PostGIS Users Discussion
Subject: Re: [postgis-users] Getting dimensions of polygon sides

Obe, Regina wrote:
> Anyrate my new problem:  This whole exercise just made me realize they
> must have put these polygons together with the constituent linework of
> the edges.  Because if I look at the abutting parcels those extra dims
> represent sides of those smaller parcels.  My simplify solution works
> great if you are just looking at one parcel, but when you look at
> abutters - it shall we say, looks crowded.  So my thought is to put
the
> position of my length annotation inside each parcel it represents.
>
> I'm lost how to do that.  So basically jimmy the ST_Centroid(the_geom)
> of my above such that it always sits inside the parcel boundary.
>
> Thanks,
> Regina
>   
I am not sure how you are placing the annotation, but this could be 
interesting....

Looking at the aftersimplify.png, you have the length of each line 
segment (and the segment itself).  The center point of each segment is 
easy.  Now you need to define a line segment offset by the text height 
of your annotation and translate it INSIDE the polygon.  The you draw 
(label?, annotate? ....) the text along that line.  If you first update 
the polygon with forceRHR() then INSIDE is always to the right of the 
line segment.  You can easily get the azmith of the linesegment from the

start and end points, add 90 degrees and compute the x and y values for 
the translate() function.  This should be basic high school trig - ;-) .

Is this what you are thinking of?

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


-----------------------------------------
The substance of this message, including any attachments, may be
confidential, legally privileged and/or exempt from disclosure
pursuant to Massachusetts law. It is intended
solely for the addressee. If you received this in error, please
contact the sender and delete the material from any computer.
-------------- next part --------------
A non-text attachment was scrubbed...
Name: labelinside.png
Type: image/png
Size: 18220 bytes
Desc: labelinside.png
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20080616/bda5bdcf/attachment.png>


More information about the postgis-users mailing list