5

I have a PostGIS database with a table called 'sites' that includes a field sitename and a field geometry.

If I select, e.g. 

select sitename, ST_AsGeoJSON(geometry) from sites limit 1;

I get:

Plot 7 | {"type":"Polygon","coordinates": [[[-111.974951718315,33.0745275800293,353], [-111.974935360975,33.0745276003588,353], [-111.974935370263,33.0745546649121,353], [-111.974951727608,33.0745546445826,353],[-111.974951718315,33.0745275800293,353]]]}

These are rectangles oriented very close to N-S. I would like to divide each plot in half along (or close to) its N-S axis. Is there an easy way to do this?

Solutions in PostGIS-SQL preferred. Ultimately I plan to write insert statements like

insert into sites (name, geometry) 
   values ('Site 1 E', (select ??east half of geometry?? from sites where name = 'Site 1'))

The closest I've gotten is ST_Split, which requires the polygon as well as the 'blade' (line used to split), but it is not clear how I can compute the blade from the geometry.

David LeBauer
  • 554
  • 3
  • 5
  • 19

1 Answers1

4

I would like to divide each plot in half along (or close to) its N-S axis. Is there an easy way to do this? [...] can compute the blade from the geometry.

Something like?

xmin+((xmax-xmin)/2)

Here we create a sample data set foo. In it we create a polygon, we calculate the xaxis to bisect the polygon with, and we create the blade.

Then we do the bisection and display the results with ST_AsText

WITH foo AS (
  SELECT geom, blade
  FROM ST_MakeEnvelope(-10,-10,10,10) AS geom
  CROSS JOIN LATERAL ( SELECT ST_xMin(geom) + (ST_xMax(geom) - ST_xMin(geom)) / 2 ) AS axis(x)
  CROSS JOIN LATERAL ST_MakeLine(
    ST_MakePoint(axis.x, ST_yMin(geom)),
    ST_MakePoint(axis.x, ST_yMax(geom))
  ) AS blade
)
SELECT ST_AsText(geom),
  ST_AsText(blade),
  ST_AsText( ST_Split(geom, blade) )
FROM foo;
                   st_astext                    |       st_astext        |                                                st_astext                                                 
------------------------------------------------+------------------------+----------------------------------------------------------------------------------------------------------
 POLYGON((-10 -10,-10 10,10 10,10 -10,-10 -10)) | LINESTRING(0 -10,0 10) | GEOMETRYCOLLECTION(POLYGON((-10 -10,-10 10,0 10,0 -10,-10 -10)),POLYGON((0 10,10 10,10 -10,0 -10,0 10)))
Evan Carroll
  • 7,071
  • 2
  • 32
  • 58
  • 1
    How about ST_Centroid ? – dbaston May 10 '17 at 00:22
  • 1
    @dbaston that sounds more cpuintensive, I showed you how to read xMin/xMax/yMin/yMax but for the purpose of creating a blade to split something, you need only two of the four. – Evan Carroll May 10 '17 at 00:24
  • Maybe it sounds CPU-intensive, but it isn't, especially in the context of ST_Split. Whether it should be used depends on the requirements, I guess. I was suggesting centroid in order to get an equal-area split. – dbaston May 10 '17 at 12:23
  • Just noticed that the inputs are strictly rectangles, which makes my above point irrelevant. – dbaston May 10 '17 at 12:24
  • @EvanCarroll this was great, I learned a lot ... thanks! Though would be awesome if the query generated insert statements. (for future reference) to generate an insert statement I had to 1) ST_Dump(ST_Split to the polygons, but 2) had to use ST_xmax to sort 3) and remove split ({1},, ({2} and ) generated from the dump; 4) ST_SetSRID to avoid error; my final query https://gist.github.com/dlebauer/f586e380fa4cfadc25551cd1df487b7d got me most of the way but still required spreadsheet-foo. – David LeBauer May 17 '17 at 22:46
  • Note that xmin+((xmax-xmin)/2) = xmin+xmax/2-xmin/2 = ((xmax+xmin)/2). – jpmc26 Jan 30 '19 at 17:22