[postgis-tickets] [SCM] PostGIS branch master updated. 4332fa3330ea651317d3c9dc3f5bd2092ba9241d

git at osgeo.org git at osgeo.org
Thu Jan 16 15:34:49 PST 2020


This is an automated email from the git hooks/post-receive script. It was
generated because a ref change was pushed to the repository containing
the project "PostGIS".

The branch, master has been updated
       via  4332fa3330ea651317d3c9dc3f5bd2092ba9241d (commit)
       via  2feb73f7b6e207e1cdb46f022b24b66b88bcbac2 (commit)
       via  7213ee370d15ebdf56a28c450fc7542c237a38a9 (commit)
       via  27aef0e58aeba2ea49eed4dcffa01f953ffcad12 (commit)
       via  af1d1e36ef2a00e75fae1481b6e974391f7c23c2 (commit)
       via  e54690bc393cea5f49653a4496dcaff3050ccc29 (commit)
       via  ebf2cc2256095c38ee62b21a4fb96bcd3151dd0f (commit)
       via  aa253960ffdc83b031aa7c598c73ed23dbb09dbd (commit)
       via  cba68ed602c4268ea50a36251e5417fa5b5ed6f1 (commit)
       via  c3de9f82a5e80fb8684b7d70f8d7ee977b8ab27f (commit)
       via  3dfb4c7defa41f0c7227172f6ae7a112ccfb3e0e (commit)
       via  6a1e56958834957bcde1fcf3c502af30e2f02389 (commit)
      from  51f870db141d468a7e5f2b3eeb10a7e8bd5fda7f (commit)

Those revisions listed above that are new to this repository have
not appeared on any other notification email; so we list those
revisions in full, below.

- Log -----------------------------------------------------------------
commit 4332fa3330ea651317d3c9dc3f5bd2092ba9241d
Merge: 2feb73f 51f870d
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 15:24:29 2020 -0800

    Merge branch 'generator' into svn-trunk-hexgrid


commit 2feb73f7b6e207e1cdb46f022b24b66b88bcbac2
Merge: 7213ee3 e5ea6e1
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 15:14:45 2020 -0800

    Merge branch 'master' of https://git.osgeo.org/gogs/postgis/postgis into svn-trunk-hexgrid


commit 7213ee370d15ebdf56a28c450fc7542c237a38a9
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 15:10:47 2020 -0800

    docbook dtd issues?

diff --git a/doc/reference_constructor.xml b/doc/reference_constructor.xml
index 6565fda..1b38b5f 100644
--- a/doc/reference_constructor.xml
+++ b/doc/reference_constructor.xml
@@ -838,18 +838,18 @@ SELECT ST_AsText( ST_TileEnvelope(3, 1, 1, ST_MakeEnvelope(-180, -90, 180, 90, 4
 			This function answers the question: what hexagons in a given Tiling(SRS, Size)
 			overlap with a given bounds.</para>
 
-			<inlinemediaobject><imageobject>
+			<para><inlinemediaobject><imageobject>
 				<imagedata fileref='images/st_hexagongrid01.png' />
-			</imageobject></inlinemediaobject>
+			</imageobject></inlinemediaobject></para>
 
 			<para>The SRS for the output hexagons is the SRS provided by the bounds geometry.</para>
 			<para>Doubling or tripling the edge size of the hexagon generates a new parent tiling that
 			fits with the origin tiling. Unfortunately, it is not possible to generate parent
 			hexagon tilings that the child tiles perfectly fit inside.</para>
 
-			<inlinemediaobject><imageobject>
+			<para><inlinemediaobject><imageobject>
 				<imagedata fileref='images/st_hexagongrid02.png' />
-			</imageobject></inlinemediaobject>
+			</imageobject></inlinemediaobject></para>
 
 			<para>Availability: 3.1</para>
 
@@ -874,9 +874,9 @@ SELECT ST_AsText( ST_TileEnvelope(3, 1, 1, ST_MakeEnvelope(-180, -90, 180, 90, 4
 		<para>If we generate a set of hexagons for each polygon boundary and filter
 			out those that do not intersect their hexagons, we end up with a tiling for
 			each polygon.</para>
-		<inlinemediaobject><imageobject>
+		<para><inlinemediaobject><imageobject>
 			<imagedata fileref='images/st_hexagongrid02.png' />
-		</imageobject></inlinemediaobject>
+		</imageobject></inlinemediaobject></para>
 		<para>Tiling states results in a hexagon coverage of each state, and multiple
 			hexagons overlapping at the borders between states.</para>
 		 <programlisting>WITH hexes AS (

commit 27aef0e58aeba2ea49eed4dcffa01f953ffcad12
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 14:47:48 2020 -0800

    Documentation and test cases

diff --git a/doc/html/images/st_hexagongrid01.png b/doc/html/images/st_hexagongrid01.png
new file mode 100644
index 0000000..05e6e32
Binary files /dev/null and b/doc/html/images/st_hexagongrid01.png differ
diff --git a/doc/html/images/st_hexagongrid02.png b/doc/html/images/st_hexagongrid02.png
new file mode 100644
index 0000000..80b7262
Binary files /dev/null and b/doc/html/images/st_hexagongrid02.png differ
diff --git a/doc/html/images/st_hexagongrid03.png b/doc/html/images/st_hexagongrid03.png
new file mode 100644
index 0000000..18f9fff
Binary files /dev/null and b/doc/html/images/st_hexagongrid03.png differ
diff --git a/doc/reference_constructor.xml b/doc/reference_constructor.xml
index 4b8e9a7..6565fda 100644
--- a/doc/reference_constructor.xml
+++ b/doc/reference_constructor.xml
@@ -809,4 +809,150 @@ SELECT ST_AsText( ST_TileEnvelope(3, 1, 1, ST_MakeEnvelope(-180, -90, 180, 90, 4
 		</refsection>
 	</refentry>
 
+
+
+	<refentry id="ST_HexagonGrid">
+		<refnamediv>
+		<refname>ST_HexagonGrid</refname>
+		<refpurpose>Returns a set of hexagons and cell indices that completely cover the bounds of the geometry argument.</refpurpose>
+		</refnamediv>
+
+		<refsynopsisdiv>
+		<funcsynopsis>
+		  <funcprototype>
+			<funcdef>setof record <function>ST_HexagonGrid</function></funcdef>
+			<paramdef><type>float8</type> <parameter>size</parameter></paramdef>
+			<paramdef><type>geometry</type> <parameter>bounds</parameter></paramdef>
+		  </funcprototype>
+		 </funcsynopsis>
+		</refsynopsisdiv>
+
+		<refsection>
+			<title>Description</title>
+
+			<para>Starts with the concept of a hexagon tiling of the plane.
+			(Not a hexagon tiling of the globe, this is not the
+			<ulink url="https://github.com/uber/h3">H3</ulink> tiling scheme.)
+			For a given planar SRS, and a given edge size, starting at the origin of the SRS,
+			there is one unique hexagonal tiling of the plane, Tiling(SRS, Size).
+			This function answers the question: what hexagons in a given Tiling(SRS, Size)
+			overlap with a given bounds.</para>
+
+			<inlinemediaobject><imageobject>
+				<imagedata fileref='images/st_hexagongrid01.png' />
+			</imageobject></inlinemediaobject>
+
+			<para>The SRS for the output hexagons is the SRS provided by the bounds geometry.</para>
+			<para>Doubling or tripling the edge size of the hexagon generates a new parent tiling that
+			fits with the origin tiling. Unfortunately, it is not possible to generate parent
+			hexagon tilings that the child tiles perfectly fit inside.</para>
+
+			<inlinemediaobject><imageobject>
+				<imagedata fileref='images/st_hexagongrid02.png' />
+			</imageobject></inlinemediaobject>
+
+			<para>Availability: 3.1</para>
+
+		</refsection>
+
+		<refsection>
+		<title>Example: Counting points in hexagons</title>
+		<para>To do a point summary against a hexagonal tiling, generate a hexagon grid using the
+		extent of the points as the bounds, then spatially join to that grid.</para>
+		 <programlisting>WITH bounds AS (
+ 	SELECT ST_SetSRID(ST_EstimatedExtent('pointtable', 'geom'), 3857) AS geom
+ )
+ SELECT Count(*), hexes.geom
+ FROM pointtable pts
+ JOIN ST_HexagonGrid(1000, bounds.geom) hexes
+ GROUP BY hexes.geom
+</programlisting>
+		</refsection>
+
+		<refsection>
+		<title>Example: Generating hex coverage of polygons</title>
+		<para>If we generate a set of hexagons for each polygon boundary and filter
+			out those that do not intersect their hexagons, we end up with a tiling for
+			each polygon.</para>
+		<inlinemediaobject><imageobject>
+			<imagedata fileref='images/st_hexagongrid02.png' />
+		</imageobject></inlinemediaobject>
+		<para>Tiling states results in a hexagon coverage of each state, and multiple
+			hexagons overlapping at the borders between states.</para>
+		 <programlisting>WITH hexes AS (
+	SELECT admin1.gid,
+	       (ST_HexagonGrid(100000, admin1.geom)).* AS hex
+	FROM admin1
+	WHERE adm0_a3 = 'USA'
+)
+SELECT hexes.*
+FROM hexes
+JOIN admin1 USING (gid)
+WHERE ST_Intersects(admin1.geom, hexes.geom)
+</programlisting>
+		</refsection>
+		<refsection>
+			<title>See Also</title>
+			<para><xref linkend="ST_TileEnvelope" /><xref linkend="ST_SquareGrid" /></para>
+		</refsection>
+	</refentry>
+
+
+	<refentry id="ST_SquareGrid">
+		<refnamediv>
+		<refname>ST_SquareGrid</refname>
+		<refpurpose>Returns a set of grid squares and cell indices that completely cover the bounds of the geometry argument.</refpurpose>
+		</refnamediv>
+
+		<refsynopsisdiv>
+		<funcsynopsis>
+		  <funcprototype>
+			<funcdef>setof record <function>ST_SquareGrid</function></funcdef>
+			<paramdef><type>float8</type> <parameter>size</parameter></paramdef>
+			<paramdef><type>geometry</type> <parameter>bounds</parameter></paramdef>
+		  </funcprototype>
+		 </funcsynopsis>
+		</refsynopsisdiv>
+
+		<refsection>
+			<title>Description</title>
+
+			<para>Starts with the concept of a square tiling of the plane.
+			For a given planar SRS, and a given edge size, starting at the origin of the SRS,
+			there is one unique square tiling of the plane, Tiling(SRS, Size).
+			This function answers the question: what grids in a given Tiling(SRS, Size)
+			overlap with a given bounds.</para>
+
+			<para>The SRS for the output squares is the SRS provided by the bounds geometry.</para>
+			<para>Doubling or edge size of the square generates a new parent tiling that
+			perfectly fits with the original tiling. Standard web map tilings in mercator
+			are just powers-of-two square grids in the mercator plane.</para>
+
+			<para>Availability: 3.1</para>
+
+		</refsection>
+
+		<refsection>
+		<title>Example: Counting points in squares</title>
+		<para>To do a point summary against a square tiling, generate a square grid using the
+		extent of the points as the bounds, then spatially join to that grid.</para>
+		 <programlisting>WITH bounds AS (
+ 	SELECT ST_SetSRID(ST_EstimatedExtent('pointtable', 'geom'), 3857) AS geom
+ )
+ SELECT Count(*), squares.geom
+ FROM pointtable pts
+ JOIN ST_HexagonGrid(1000, bounds.geom) squares
+ GROUP BY squares.geom
+</programlisting>
+		</refsection>
+
+		<refsection>
+			<title>See Also</title>
+			<para><xref linkend="ST_TileEnvelope" /><xref linkend="ST_HexagonGrid" /></para>
+		</refsection>
+	</refentry>
+
+
+
+
   </sect1>
diff --git a/regress/core/regress.sql b/regress/core/regress.sql
index 6f5eea3..f77b1c5 100644
--- a/regress/core/regress.sql
+++ b/regress/core/regress.sql
@@ -286,6 +286,27 @@ select '252', ST_AsEWKT(ST_TileEnvelope(10,300,387, margin => 0.5));
 select '253', ST_AsEWKT(ST_TileEnvelope(10,300,387, margin => 2));
 select '254', ST_AsEWKT(ST_TileEnvelope(10,300,387, margin => -0.3));
 
+-- ST_Hexagon()
+select '300', ST_AsEWKT(ST_SnapToGrid(ST_Hexagon(10, 0, 0),0.00001));
+select '301', ST_AsEWKT(ST_SnapToGrid(ST_Hexagon(10, 1, 1),0.00001));
+select '302', ST_AsEWKT(ST_SnapToGrid(ST_Hexagon(10, -1, -1),0.00001));
+select '303', ST_AsEWKT(ST_SnapToGrid(ST_Hexagon(10, 100, -100),0.00001));
+-- ST_Square()
+select '304', ST_AsEWKT(ST_SnapToGrid(ST_Square(10, 0, 0),0.00001));
+select '305', ST_AsEWKT(ST_SnapToGrid(ST_Square(10, 1, 1),0.00001));
+select '306', ST_AsEWKT(ST_SnapToGrid(ST_Square(10, -1, -1),0.00001));
+select '307', ST_AsEWKT(ST_SnapToGrid(ST_Square(10, 100, -100),0.00001));
+-- ST_HexagonGrid()
+select '308', Count(*) FROM ST_HexagonGrid(100000, ST_TileEnvelope(4, 7, 7));
+select '309', Count(*) FROM ST_HexagonGrid(100000, ST_TileEnvelope(4, 7, 7)) hex, ST_TileEnvelope(4, 7, 7) tile WHERE NOT ST_Intersects(hex.geom, tile);
+select '310', Count(*) FROM ST_HexagonGrid(100000, ST_TileEnvelope(4, 2, 2));
+select '311', Count(*) FROM ST_HexagonGrid(100000, ST_TileEnvelope(4, 2, 2)) hex, ST_TileEnvelope(4, 2, 2) tile WHERE NOT ST_Intersects(hex.geom, tile);
+-- ST_SquareGrid()
+select '312', Count(*) FROM ST_SquareGrid(100000, ST_TileEnvelope(4, 7, 7));
+select '313', Count(*) FROM ST_SquareGrid(100000, ST_TileEnvelope(4, 7, 7)) hex, ST_TileEnvelope(4, 7, 7) tile WHERE NOT ST_Intersects(hex.geom, tile);
+select '314', Count(*) FROM ST_SquareGrid(100000, ST_TileEnvelope(4, 2, 2));
+select '315', Count(*) FROM ST_SquareGrid(100000, ST_TileEnvelope(4, 2, 2)) hex, ST_TileEnvelope(4, 2, 2) tile WHERE NOT ST_Intersects(hex.geom, tile);
+
 -- Drop test table
 DROP table test;
 
diff --git a/regress/core/regress_expected b/regress/core/regress_expected
index f07b639..e7cef06 100644
--- a/regress/core/regress_expected
+++ b/regress/core/regress_expected
@@ -214,3 +214,19 @@ ERROR:  ST_TileEnvelope: Margin must not be less than -50%, margin=-0.510000
 252|SRID=3857;POLYGON((-8316348.67742708 4833266.17252821,-8316348.67742708 4911537.68949223,-8238077.16046306 4911537.68949223,-8238077.16046306 4833266.17252821,-8316348.67742708 4833266.17252821))
 253|SRID=3857;POLYGON((-8375052.31515009 4774562.53480519,-8375052.31515009 4970241.32721524,-8179373.52274004 4970241.32721524,-8179373.52274004 4774562.53480519,-8375052.31515009 4774562.53480519))
 254|SRID=3857;POLYGON((-8285040.07064147 4864574.77931381,-8285040.07064147 4880229.08270662,-8269385.76724866 4880229.08270662,-8269385.76724866 4864574.77931381,-8285040.07064147 4864574.77931381))
+300|POLYGON((-10 0,-5 -8.66025,5 -8.66025,10 0,5 8.66025,-5 8.66025,-10 0))
+301|POLYGON((5 25.98076,10 17.32051,20 17.32051,25 25.98076,20 34.64102,10 34.64102,5 25.98076))
+302|POLYGON((-25 -8.66025,-20 -17.32051,-10 -17.32051,-5 -8.66025,-10 0,-20 0,-25 -8.66025))
+303|POLYGON((1490 -1732.05081,1495 -1740.71106,1505 -1740.71106,1510 -1732.05081,1505 -1723.39055,1495 -1723.39055,1490 -1732.05081))
+304|POLYGON((0 0,0 10,10 10,10 0,0 0))
+305|POLYGON((10 10,10 20,20 20,20 10,10 10))
+306|POLYGON((-10 -10,-10 0,0 0,0 -10,-10 -10))
+307|POLYGON((1000 -1000,1000 -990,1010 -990,1010 -1000,1000 -1000))
+308|270
+309|0
+310|279
+311|0
+312|702
+313|0
+314|676
+315|0

commit af1d1e36ef2a00e75fae1481b6e974391f7c23c2
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 11:33:25 2020 -0800

    Fix up boundary logic for hex columns

diff --git a/postgis/lwgeom_generate_grid.c b/postgis/lwgeom_generate_grid.c
index c2f50e6..d59c7ec 100644
--- a/postgis/lwgeom_generate_grid.c
+++ b/postgis/lwgeom_generate_grid.c
@@ -125,12 +125,12 @@ hexagon_grid_state(double size, const GBOX *gbox, int32_t srid)
 	/* Column address is just width / column width, with an */
 	/* adjustment to account for partial overlaps */
 	state->column_min = floor(gbox->xmin / col_width);
-	if(gbox->xmin - state->column_min * size > size)
+	if(gbox->xmin - state->column_min * col_width > size)
 		state->column_min++;
 
 	state->column_max = ceil(gbox->xmax / col_width);
-	if(state->column_max * size - gbox->xmax > size)
-		state->column_max++;
+	if(state->column_max * col_width - gbox->xmax > size)
+		state->column_max--;
 
 	/* Row address range depends on the odd/even column we're on */
 	state->row_min_even = floor(gbox->ymin/row_height + 0.5);

commit e54690bc393cea5f49653a4496dcaff3050ccc29
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Jan 16 10:45:06 2020 -0800

    Flip around logic for handling of 'bounds'.
    Bounds is no longer the domain of the grid, instead the grid always
    starts at origin 0,0 in whatever SRS the bounds are. Bounds is the
    rectangle-to-be-filled with shapes, and only those shapes needed to
    fill the bounds are instantiated, along with their cell coordinates.

diff --git a/postgis/lwgeom_generate_grid.c b/postgis/lwgeom_generate_grid.c
index 854c825..c2f50e6 100644
--- a/postgis/lwgeom_generate_grid.c
+++ b/postgis/lwgeom_generate_grid.c
@@ -43,30 +43,49 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS);
 
 /* ********* ********* ********* ********* ********* ********* ********* ********* */
 
-/* XXX TO DO: what about negative sizes? */
+typedef enum
+{
+	SHAPE_SQUARE,
+	SHAPE_HEXAGON,
+	SHAPE_TRIANGLE
+} GeometryShape;
 
-/* Ratio of hexagon edge length to height = 2*cos(pi/6) */
-#define HXR 1.7320508075688774
+typedef struct GeometryGridState
+{
+	GeometryShape cell_shape;
+	bool done;
+	GBOX bounds;
+	int32_t srid;
+	double size;
+	int32_t i, j;
+}
+GeometryGridState;
 
-static const double hex_x[] = {0.0, 1.0, 1.5,     1.0, 0.0, -0.5,     0.0};
-static const double hex_y[] = {0.0, 0.0, 0.5*HXR, HXR, HXR,  0.5*HXR, 0.0};
+/* ********* ********* ********* ********* ********* ********* ********* ********* */
+
+/* Ratio of hexagon edge length to edge->center vertical distance = cos(M_PI/6.0) */
+#define H 0.8660254037844387
+
+/* Build origin hexagon centered around origin point */
+static const double hex_x[] = {-1.0, -0.5,  0.5, 1.0, 0.5, -0.5, -1.0};
+static const double hex_y[] = { 0.0, -1*H, -1*H, 0.0,   H,    H,  0.0};
 
 static LWGEOM *
-hexagon(double origin_x, double origin_y, double edge, int cell_i, int cell_j, uint32_t srid)
+hexagon(double origin_x, double origin_y, double size, int cell_i, int cell_j, uint32_t srid)
 {
-	double height = edge * HXR;
-	POINT4D pt = {0, 0, 0, 0};
+	double height = size * 2 * H;
+	POINT4D pt;
 	POINTARRAY **ppa = lwalloc(sizeof(POINTARRAY*));
 	POINTARRAY *pa = ptarray_construct(0, 0, 7);
 	uint32_t i;
 
-	double offset_x = origin_x + (1.5 * edge * cell_i);
-	double offset_y = origin_y + (height * cell_j) + (0.5 * height * (cell_i % 2));
+	double offset_x = origin_x + (1.5 * size * cell_i);
+	double offset_y = origin_y + (height * cell_j) + (0.5 * height * (abs(cell_i) % 2));
 
 	for (i = 0; i < 7; ++i)
 	{
-		pt.x = edge * hex_x[i] + offset_x;
-		pt.y = edge * hex_y[i] + offset_y;
+		pt.x = size * hex_x[i] + offset_x;
+		pt.y = size * hex_y[i] + offset_y;
 		ptarray_set_point4d(pa, i, &pt);
 	}
 
@@ -74,17 +93,76 @@ hexagon(double origin_x, double origin_y, double edge, int cell_i, int cell_j, u
 	return lwpoly_as_lwgeom(lwpoly_construct(srid, NULL, 1 /* nrings */, ppa));
 }
 
-static int
-hexagon_grid_size(double edge, const GBOX *gbox, int *cell_width, int *cell_height)
+
+typedef struct HexagonGridState
 {
-	double hex_height = edge * HXR;
-	double bounds_width  = gbox->xmax - gbox->xmin;
-	double bounds_height = gbox->ymax - gbox->ymin;
-	int width  = (int)floor(bounds_width/(1.5*edge)) + 1;
-	int height = (int)floor(bounds_height/hex_height) + 1;
-	if (cell_width) *cell_width = width;
-	if (cell_height) *cell_height = height;
-	return width * height;
+	GeometryShape cell_shape;
+	bool done;
+	GBOX bounds;
+	int32_t srid;
+	double size;
+	int32_t i, j;
+	int32_t column_min, column_max;
+	int32_t row_min_odd, row_max_odd;
+	int32_t row_min_even, row_max_even;
+}
+HexagonGridState;
+
+static HexagonGridState *
+hexagon_grid_state(double size, const GBOX *gbox, int32_t srid)
+{
+	HexagonGridState *state = palloc0(sizeof(HexagonGridState));
+	double col_width = 1.5 * size;
+	double row_height = size * 2 * H;
+
+	/* fill in state */
+	state->cell_shape = SHAPE_HEXAGON;
+	state->size = size;
+	state->srid = srid;
+	state->done = false;
+	state->bounds = *gbox;
+
+	/* Column address is just width / column width, with an */
+	/* adjustment to account for partial overlaps */
+	state->column_min = floor(gbox->xmin / col_width);
+	if(gbox->xmin - state->column_min * size > size)
+		state->column_min++;
+
+	state->column_max = ceil(gbox->xmax / col_width);
+	if(state->column_max * size - gbox->xmax > size)
+		state->column_max++;
+
+	/* Row address range depends on the odd/even column we're on */
+	state->row_min_even = floor(gbox->ymin/row_height + 0.5);
+	state->row_max_even = floor(gbox->ymax/row_height + 0.5);
+	state->row_min_odd  = floor(gbox->ymin/row_height);
+	state->row_max_odd  = floor(gbox->ymax/row_height);
+
+	/* Set initial state */
+	state->i = state->column_min;
+	state->j = (state->i % 2) ? state->row_min_odd : state->row_min_even;
+
+	return state;
+}
+
+static bool
+hexagon_state_next(HexagonGridState *state)
+{
+	if (!state || state->done) return false;
+	/* Move up one row */
+	state->j++;
+	/* Off the end, increment column counter, reset row counter back to (appropriate) minimum */
+	if (state->j > ((state->i % 2) ? state->row_max_odd : state->row_max_even))
+	{
+		state->i++;
+		state->j = ((state->i % 2) ? state->row_min_odd : state->row_min_even);
+	}
+	/* Column counter over max means we have used up all combinations */
+	if (state->i > state->column_max)
+	{
+		state->done = true;
+	}
+	return !state->done;
 }
 
 /* ********* ********* ********* ********* ********* ********* ********* ********* */
@@ -99,48 +177,67 @@ square(double origin_x, double origin_y, double size, int cell_i, int cell_j, ui
 	return (LWGEOM*)lwpoly_construct_envelope(srid, ll_x, ll_y, ur_x, ur_y);
 }
 
-static int
-square_grid_size(double edge, const GBOX *gbox, int *cell_width, int *cell_height)
+typedef struct SquareGridState
 {
-	int width  = (int)floor((gbox->xmax - gbox->xmin)/edge) + 1;
-	int height = (int)floor((gbox->ymax - gbox->ymin)/edge) + 1;
-	if (cell_width) *cell_width = width;
-	if (cell_height) *cell_height = height;
-	return width * height;
+	GeometryShape cell_shape;
+	bool done;
+	GBOX bounds;
+	int32_t srid;
+	double size;
+	int32_t i, j;
+	int32_t column_min, column_max;
+	int32_t row_min, row_max;
 }
+SquareGridState;
 
-/* ********* ********* ********* ********* ********* ********* ********* ********* */
-
-typedef enum
+static SquareGridState *
+square_grid_state(double size, const GBOX *gbox, int32_t srid)
 {
-	SHAPE_SQUARE,
-	SHAPE_HEXAGON,
-	SHAPE_TRIANGLE
-} GeometryShape;
+	SquareGridState *state = palloc0(sizeof(SquareGridState));
+
+	/* fill in state */
+	state->cell_shape = SHAPE_SQUARE;
+	state->size = size;
+	state->srid = srid;
+	state->done = false;
+	state->bounds = *gbox;
+
+	state->column_min = floor(gbox->xmin / size);
+	state->column_max = floor(gbox->xmax / size);
+	state->row_min = floor(gbox->ymin / size);
+	state->row_max = floor(gbox->ymax / size);
+	state->i = state->column_min;
+	state->j = state->row_min;
+	return state;
+}
 
-typedef struct GeometryGridState
+static bool
+square_state_next(SquareGridState *state)
 {
-	int64_t cell_current;
-	int64_t cell_count;
-	int32_t cell_width;
-	int32_t cell_height;
-	GBOX bounds;
-	double size;
-	int32_t srid;
-	GeometryShape cell_shape;
+	if (!state || state->done) return false;
+	/* Move up one row */
+	state->j++;
+	/* Off the end, increment column counter, reset row counter back to (appropriate) minimum */
+	if (state->j > state->row_max)
+	{
+		state->i++;
+		state->j = state->row_min;
+	}
+	/* Column counter over max means we have used up all combinations */
+	if (state->i > state->column_max)
+	{
+		state->done = true;
+	}
+	return !state->done;
 }
-GeometryGridState;
 
 /**
-* ST_HexagonGrid(size float8 DEFAULT 1.0,
-*		bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
-* ST_SquareGrid(size float8 DEFAULT 1.0,
-*		bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
+* ST_HexagonGrid(size float8, bounds geometry)
+* ST_SquareGrid(size float8, bounds geometry)
 */
 PG_FUNCTION_INFO_V1(ST_ShapeGrid);
 Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 {
-	int32_t i, j;
 	FuncCallContext *funcctx;
 	MemoryContext oldcontext, newcontext;
 
@@ -157,12 +254,12 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 	{
 		const char *func_name;
 		double bounds_width, bounds_height;
-		int gbounds_is_empty;
+		char gbounds_is_empty;
 		GBOX bounds;
 		double size;
 		funcctx = SRF_FIRSTCALL_INIT();
 
-		/* get a local copy of what we're doing a dump points on */
+		/* Get a local copy of the bounds we are going to fill with shapes */
 		gbounds = PG_GETARG_GSERIALIZED_P(1);
 		size = PG_GETARG_FLOAT8(0);
 
@@ -181,13 +278,6 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 		newcontext = funcctx->multi_call_memory_ctx;
 		oldcontext = MemoryContextSwitchTo(newcontext);
 
-		/* Create function state */
-		state = lwalloc(sizeof(GeometryGridState));
-		state->cell_current = 0;
-		state->bounds = bounds;
-		state->srid = gserialized_get_srid(gbounds);
-		state->size = size;
-
 		/*
 		* Support both hexagon and square grids with one function,
 		* by checking the calling signature up front.
@@ -195,15 +285,11 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 		func_name = get_func_name(fcinfo->flinfo->fn_oid);
 		if (strcmp(func_name, "st_hexagongrid") == 0)
 		{
-			state->cell_shape = SHAPE_HEXAGON;
-			state->cell_count = hexagon_grid_size(size, &state->bounds,
-			                                      &state->cell_width, &state->cell_height);
+			state = (GeometryGridState*)hexagon_grid_state(size, &bounds, gserialized_get_srid(gbounds));
 		}
 		else if (strcmp(func_name, "st_squaregrid") == 0)
 		{
-			state->cell_shape = SHAPE_SQUARE;
-			state->cell_count = square_grid_size(size, &state->bounds,
-			                                     &state->cell_width, &state->cell_height);
+			state = (GeometryGridState*)square_grid_state(size, &bounds, gserialized_get_srid(gbounds));
 		}
 		else
 		{
@@ -232,22 +318,27 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 	state = funcctx->user_fctx;
 
 	/* Stop when we've used up all the grid squares */
-	if (state->cell_current >= state->cell_count)
+	if (state->done)
 	{
 		SRF_RETURN_DONE(funcctx);
 	}
 
-	i = state->cell_current % state->cell_width;
-	j = state->cell_current / state->cell_width;
+	/* Store grid cell coordinates */
+	tuple_arr[1] = Int32GetDatum(state->i);
+	tuple_arr[2] = Int32GetDatum(state->j);
 
 	/* Generate geometry */
 	switch (state->cell_shape)
 	{
 		case SHAPE_HEXAGON:
-			lwgeom = hexagon(state->bounds.xmin, state->bounds.ymin, state->size, i, j, state->srid);
+			lwgeom = hexagon(0.0, 0.0, state->size, state->i, state->j, state->srid);
+			/* Increment to next cell */
+			hexagon_state_next((HexagonGridState*)state);
 			break;
 		case SHAPE_SQUARE:
-			lwgeom = square(state->bounds.xmin, state->bounds.ymin, state->size, i, j, state->srid);
+			lwgeom = square(0.0, 0.0, state->size, state->i, state->j, state->srid);
+			/* Increment to next cell */
+			square_state_next((SquareGridState*)state);
 			break;
 		default:
 			ereport(ERROR, (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
@@ -257,30 +348,24 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 	tuple_arr[0] = PointerGetDatum(geometry_serialize(lwgeom));
 	lwfree(lwgeom);
 
-	/* Store grid cell coordinates */
-	tuple_arr[1] = Int32GetDatum(i);
-	tuple_arr[2] = Int32GetDatum(j);
-
+	/* Form tuple and return */
 	tuple = heap_form_tuple(funcctx->tuple_desc, tuple_arr, isnull);
 	result = HeapTupleGetDatum(tuple);
-	state->cell_current++;
 	SRF_RETURN_NEXT(funcctx, result);
 }
 
 /**
-* ST_Hexagon(size double,
-*            origin geometry default 'POINT(0 0)',
-*            int cellx default 0, int celly default 0)
+* ST_Hexagon(size double, cell_i integer default 0, cell_j integer default 0, origin geometry default 'POINT(0 0)')
 */
 PG_FUNCTION_INFO_V1(ST_Hexagon);
 Datum ST_Hexagon(PG_FUNCTION_ARGS)
 {
 	LWPOINT *lwpt;
 	double size = PG_GETARG_FLOAT8(0);
-	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(1);
 	GSERIALIZED *ghex;
-	int cell_i = PG_GETARG_INT32(2);
-	int cell_j = PG_GETARG_INT32(3);
+	int cell_i = PG_GETARG_INT32(1);
+	int cell_j = PG_GETARG_INT32(2);
+	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(3);
 	LWGEOM *lwhex;
 	LWGEOM *lworigin = lwgeom_from_gserialized(gorigin);
 
@@ -308,19 +393,17 @@ Datum ST_Hexagon(PG_FUNCTION_ARGS)
 
 
 /**
-* ST_Square(size double,
-*           origin geometry default 'POINT(0 0)',
-*           int cellx default 0, int celly default 0)
+* ST_Square(size double, cell_i integer, cell_j integer, origin geometry default 'POINT(0 0)')
 */
 PG_FUNCTION_INFO_V1(ST_Square);
 Datum ST_Square(PG_FUNCTION_ARGS)
 {
 	LWPOINT *lwpt;
 	double size = PG_GETARG_FLOAT8(0);
-	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(1);
 	GSERIALIZED *gsqr;
-	int cell_i = PG_GETARG_INT32(2);
-	int cell_j = PG_GETARG_INT32(3);
+	int cell_i = PG_GETARG_INT32(1);
+	int cell_j = PG_GETARG_INT32(2);
+	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(3);
 	LWGEOM *lwsqr;
 	LWGEOM *lworigin = lwgeom_from_gserialized(gorigin);
 
diff --git a/postgis/postgis.sql.in b/postgis/postgis.sql.in
index 2933c08..62cc5ad 100644
--- a/postgis/postgis.sql.in
+++ b/postgis/postgis.sql.in
@@ -5975,42 +5975,39 @@ CREATE OR REPLACE FUNCTION ST_InterpolatePoint(Line geometry, Point geometry)
 ---
 
 -- Availability: 3.0.0
-CREATE OR REPLACE FUNCTION ST_Hexagon(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
+CREATE OR REPLACE FUNCTION ST_Hexagon(size float8, cell_i integer, cell_j integer, origin geometry DEFAULT 'POINT(0 0)')
 	RETURNS geometry
 	AS 'MODULE_PATHNAME', 'ST_Hexagon'
-	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	LANGUAGE 'c' IMMUTABLE STRICT
+	PARALLEL SAFE
 	_COST_MEDIUM;
 
 -- Availability: 3.0.0
-CREATE OR REPLACE FUNCTION ST_Square(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
+CREATE OR REPLACE FUNCTION ST_Square(size float8, cell_i integer, cell_j integer, origin geometry DEFAULT 'POINT(0 0)')
 	RETURNS geometry
 	AS 'MODULE_PATHNAME', 'ST_Square'
-	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	LANGUAGE 'c' IMMUTABLE STRICT
+	PARALLEL SAFE
 	_COST_MEDIUM;
 
--- Availability: 1.0.0
-CREATE TYPE geometry_grid AS (
-	geom geometry,
-	i integer,
-	j integer
-);
-
 -- Availability: 3.0.0
-CREATE OR REPLACE FUNCTION ST_HexagonGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
-    RETURNS SETOF geometry_grid
+CREATE OR REPLACE FUNCTION ST_HexagonGrid(size float8, bounds geometry, OUT geom geometry, OUT i integer, OUT j integer)
+    RETURNS SETOF record
 	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
-	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	LANGUAGE 'c' IMMUTABLE STRICT
+	PARALLEL SAFE
 	_COST_MEDIUM;
 
 -- Availability: 3.0.0
-CREATE OR REPLACE FUNCTION ST_SquareGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
-    RETURNS SETOF geometry_grid
+CREATE OR REPLACE FUNCTION ST_SquareGrid(size float8, bounds geometry, OUT geom geometry, OUT i integer, OUT j integer)
+    RETURNS SETOF record
 	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
-	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	LANGUAGE 'c' IMMUTABLE STRICT
+	PARALLEL SAFE
 	_COST_MEDIUM;
 
 
--- moved to separate file cause its invovled
+-- moved to separate file cause its involved
 #include "postgis_brin.sql.in"
 
 ---------------------------------------------------------------

commit ebf2cc2256095c38ee62b21a4fb96bcd3151dd0f
Merge: aa25396 2ee8e58
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Fri Jan 10 14:57:11 2020 -0800

    Merge branch 'master' of https://git.osgeo.org/gogs/postgis/postgis into svn-trunk-hexgrid

diff --cc postgis/postgis.sql.in
index 773ecb6,ded515c..2933c08
--- a/postgis/postgis.sql.in
+++ b/postgis/postgis.sql.in
@@@ -5988,49 -5967,9 +5967,49 @@@ CREATE OR REPLACE FUNCTION ST_LocateBet
  CREATE OR REPLACE FUNCTION ST_InterpolatePoint(Line geometry, Point geometry)
  	RETURNS float8
  	AS 'MODULE_PATHNAME', 'ST_InterpolatePoint'
- 	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+ 	LANGUAGE 'c' IMMUTABLE STRICT PARALLEL SAFE
  	_COST_MEDIUM;
  
 +---------------------------------------------------------------
 +-- Grid / Hexagon coverage functions
 +---
 +
 +-- Availability: 3.0.0
 +CREATE OR REPLACE FUNCTION ST_Hexagon(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
 +	RETURNS geometry
 +	AS 'MODULE_PATHNAME', 'ST_Hexagon'
 +	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
 +	_COST_MEDIUM;
 +
 +-- Availability: 3.0.0
 +CREATE OR REPLACE FUNCTION ST_Square(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
 +	RETURNS geometry
 +	AS 'MODULE_PATHNAME', 'ST_Square'
 +	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
 +	_COST_MEDIUM;
 +
 +-- Availability: 1.0.0
 +CREATE TYPE geometry_grid AS (
 +	geom geometry,
 +	i integer,
 +	j integer
 +);
 +
 +-- Availability: 3.0.0
 +CREATE OR REPLACE FUNCTION ST_HexagonGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
 +    RETURNS SETOF geometry_grid
 +	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
 +	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
 +	_COST_MEDIUM;
 +
 +-- Availability: 3.0.0
 +CREATE OR REPLACE FUNCTION ST_SquareGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
 +    RETURNS SETOF geometry_grid
 +	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
 +	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
 +	_COST_MEDIUM;
 +
 +
  -- moved to separate file cause its invovled
  #include "postgis_brin.sql.in"
  

commit aa253960ffdc83b031aa7c598c73ed23dbb09dbd
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Tue Sep 3 15:27:53 2019 -0700

    Deal with some PR comments

diff --git a/postgis/Makefile.in b/postgis/Makefile.in
index d821041..237dd38 100644
--- a/postgis/Makefile.in
+++ b/postgis/Makefile.in
@@ -91,7 +91,7 @@ PG_OBJS= \
 	lwgeom_geos_prepared.o \
 	lwgeom_geos_clean.o \
 	lwgeom_geos_relatematch.o \
-	lwgeom_hexgrid.o \
+	lwgeom_generate_grid.o \
 	lwgeom_export.o \
 	lwgeom_in_gml.o \
 	lwgeom_in_kml.o \
diff --git a/postgis/lwgeom_hexgrid.c b/postgis/lwgeom_generate_grid.c
similarity index 89%
rename from postgis/lwgeom_hexgrid.c
rename to postgis/lwgeom_generate_grid.c
index 61278b4..854c825 100644
--- a/postgis/lwgeom_hexgrid.c
+++ b/postgis/lwgeom_generate_grid.c
@@ -18,9 +18,7 @@
  *
  **********************************************************************
  *
- * Copyright 2009 Mark Cave-Ayland <mark.cave-ayland at siriusit.co.uk>
- * Copyright 2009-2017 Paul Ramsey <pramsey at cleverelephant.ca>
- * Copyright 2018 Darafei Praliaskouski <me at komzpa.net>
+ * Copyright 2018 Paul Ramsey <pramsey at cleverelephant.ca>
  *
  **********************************************************************/
 
@@ -57,14 +55,13 @@ static LWGEOM *
 hexagon(double origin_x, double origin_y, double edge, int cell_i, int cell_j, uint32_t srid)
 {
 	double height = edge * HXR;
-	POINT4D pt;
+	POINT4D pt = {0, 0, 0, 0};
 	POINTARRAY **ppa = lwalloc(sizeof(POINTARRAY*));
 	POINTARRAY *pa = ptarray_construct(0, 0, 7);
 	uint32_t i;
 
 	double offset_x = origin_x + (1.5 * edge * cell_i);
 	double offset_y = origin_y + (height * cell_j) + (0.5 * height * (cell_i % 2));
-	pt.z = pt.m = 0.0;
 
 	for (i = 0; i < 7; ++i)
 	{
@@ -92,30 +89,14 @@ hexagon_grid_size(double edge, const GBOX *gbox, int *cell_width, int *cell_heig
 
 /* ********* ********* ********* ********* ********* ********* ********* ********* */
 
-static const double sqr_x[] = {0.0, 1.0, 1.0, 0.0, 0.0};
-static const double sqr_y[] = {0.0, 0.0, 1.0, 1.0, 0.0};
-
 static LWGEOM *
-square(double origin_x, double origin_y, double E, int cell_i, int cell_j, uint32_t srid)
+square(double origin_x, double origin_y, double size, int cell_i, int cell_j, uint32_t srid)
 {
-	POINT4D pt;
-	POINTARRAY **ppa = lwalloc(sizeof(POINTARRAY*));
-	POINTARRAY *pa = ptarray_construct(0, 0, 5);
-	uint32_t i;
-
-	double offset_x = origin_x + (E * cell_i);
-	double offset_y = origin_y + (E * cell_j);
-	pt.z = pt.m = 0.0;
-
-	for (i = 0; i < 5; ++i)
-	{
-		pt.x = E * sqr_x[i] + offset_x;
-		pt.y = E * sqr_y[i] + offset_y;
-		ptarray_set_point4d(pa, i, &pt);
-	}
-
-	ppa[0] = pa;
-	return lwpoly_as_lwgeom(lwpoly_construct(srid, NULL, 1 /* nrings */, ppa));
+	double ll_x = origin_x + (size * cell_i);
+	double ll_y = origin_y + (size * cell_j);
+	double ur_x = origin_x + (size * (cell_i + 1));
+	double ur_y = origin_y + (size * (cell_j + 1));
+	return (LWGEOM*)lwpoly_construct_envelope(srid, ll_x, ll_y, ur_x, ur_y);
 }
 
 static int
@@ -139,25 +120,27 @@ typedef enum
 
 typedef struct GeometryGridState
 {
-	int cell_current;
-	int cell_count;
-	int cell_width;
-	int cell_height;
+	int64_t cell_current;
+	int64_t cell_count;
+	int32_t cell_width;
+	int32_t cell_height;
 	GBOX bounds;
 	double size;
-	uint32_t srid;
+	int32_t srid;
 	GeometryShape cell_shape;
 }
 GeometryGridState;
 
 /**
-* ST_HexagonGrid(size double default 1.0,
-*            bounds geometry default 'LINESTRING(0 0, 100 100)')
+* ST_HexagonGrid(size float8 DEFAULT 1.0,
+*		bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
+* ST_SquareGrid(size float8 DEFAULT 1.0,
+*		bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
 */
 PG_FUNCTION_INFO_V1(ST_ShapeGrid);
 Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 {
-	int i, j;
+	int32_t i, j;
 	FuncCallContext *funcctx;
 	MemoryContext oldcontext, newcontext;
 
@@ -205,6 +188,10 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 		state->srid = gserialized_get_srid(gbounds);
 		state->size = size;
 
+		/*
+		* Support both hexagon and square grids with one function,
+		* by checking the calling signature up front.
+		*/
 		func_name = get_func_name(fcinfo->flinfo->fn_oid);
 		if (strcmp(func_name, "st_hexagongrid") == 0)
 		{

commit cba68ed602c4268ea50a36251e5417fa5b5ed6f1
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Aug 22 20:29:10 2019 -0700

    push

diff --git a/postgis/lwgeom_hexgrid.c b/postgis/lwgeom_hexgrid.c
index 0ab720d..61278b4 100644
--- a/postgis/lwgeom_hexgrid.c
+++ b/postgis/lwgeom_hexgrid.c
@@ -163,8 +163,8 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 
 	GSERIALIZED *gbounds;
 	GeometryGridState *state;
-	LWGEOM *lwgeom;
 
+	LWGEOM *lwgeom;
 	bool isnull[3] = {0,0,0}; /* needed to say no value is null */
 	Datum tuple_arr[3]; /* used to construct the composite return value */
 	HeapTuple tuple;

commit c3de9f82a5e80fb8684b7d70f8d7ee977b8ab27f
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Aug 22 20:21:42 2019 -0700

    mixed decls and code

diff --git a/postgis/lwgeom_hexgrid.c b/postgis/lwgeom_hexgrid.c
index 8024032..0ab720d 100644
--- a/postgis/lwgeom_hexgrid.c
+++ b/postgis/lwgeom_hexgrid.c
@@ -163,6 +163,7 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 
 	GSERIALIZED *gbounds;
 	GeometryGridState *state;
+	LWGEOM *lwgeom;
 
 	bool isnull[3] = {0,0,0}; /* needed to say no value is null */
 	Datum tuple_arr[3]; /* used to construct the composite return value */
@@ -253,7 +254,6 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 	j = state->cell_current / state->cell_width;
 
 	/* Generate geometry */
-	LWGEOM *lwgeom;
 	switch (state->cell_shape)
 	{
 		case SHAPE_HEXAGON:

commit 3dfb4c7defa41f0c7227172f6ae7a112ccfb3e0e
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Aug 22 20:12:35 2019 -0700

    CI build issues

diff --git a/postgis/lwgeom_hexgrid.c b/postgis/lwgeom_hexgrid.c
index 45b6769..8024032 100644
--- a/postgis/lwgeom_hexgrid.c
+++ b/postgis/lwgeom_hexgrid.c
@@ -27,6 +27,7 @@
 #include "postgres.h"
 #include "fmgr.h"
 #include "funcapi.h"
+#include "access/htup_details.h"
 
 #include "../postgis_config.h"
 #include "lwgeom_pg.h"
@@ -156,11 +157,11 @@ GeometryGridState;
 PG_FUNCTION_INFO_V1(ST_ShapeGrid);
 Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 {
+	int i, j;
 	FuncCallContext *funcctx;
 	MemoryContext oldcontext, newcontext;
 
 	GSERIALIZED *gbounds;
-	GBOX bounds;
 	GeometryGridState *state;
 
 	bool isnull[3] = {0,0,0}; /* needed to say no value is null */
@@ -170,15 +171,20 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 
 	if (SRF_IS_FIRSTCALL())
 	{
+		const char *func_name;
+		double bounds_width, bounds_height;
+		int gbounds_is_empty;
+		GBOX bounds;
+		double size;
 		funcctx = SRF_FIRSTCALL_INIT();
 
 		/* get a local copy of what we're doing a dump points on */
 		gbounds = PG_GETARG_GSERIALIZED_P(1);
-		double size = PG_GETARG_FLOAT8(0);
+		size = PG_GETARG_FLOAT8(0);
 
-		int gbounds_is_empty = (gserialized_get_gbox_p(gbounds, &bounds) == LW_FAILURE);
-		double bounds_width = bounds.xmax - bounds.xmin;
-		double bounds_height = bounds.ymax - bounds.ymin;
+		gbounds_is_empty = (gserialized_get_gbox_p(gbounds, &bounds) == LW_FAILURE);
+		bounds_width = bounds.xmax - bounds.xmin;
+		bounds_height = bounds.ymax - bounds.ymin;
 
 		/* quick opt-out if we get nonsensical inputs  */
 		if (size <= 0.0 || gbounds_is_empty ||
@@ -198,7 +204,7 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 		state->srid = gserialized_get_srid(gbounds);
 		state->size = size;
 
-		const char *func_name = get_func_name(fcinfo->flinfo->fn_oid);
+		func_name = get_func_name(fcinfo->flinfo->fn_oid);
 		if (strcmp(func_name, "st_hexagongrid") == 0)
 		{
 			state->cell_shape = SHAPE_HEXAGON;
@@ -243,8 +249,8 @@ Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
 		SRF_RETURN_DONE(funcctx);
 	}
 
-	int i = state->cell_current % state->cell_width;
-	int j = state->cell_current / state->cell_width;
+	i = state->cell_current % state->cell_width;
+	j = state->cell_current / state->cell_width;
 
 	/* Generate geometry */
 	LWGEOM *lwgeom;

commit 6a1e56958834957bcde1fcf3c502af30e2f02389
Author: Paul Ramsey <pramsey at cleverelephant.ca>
Date:   Thu Aug 22 15:47:57 2019 -0700

    First draft, still needs documentation

diff --git a/postgis/Makefile.in b/postgis/Makefile.in
index 76a2c01..d821041 100644
--- a/postgis/Makefile.in
+++ b/postgis/Makefile.in
@@ -91,6 +91,7 @@ PG_OBJS= \
 	lwgeom_geos_prepared.o \
 	lwgeom_geos_clean.o \
 	lwgeom_geos_relatematch.o \
+	lwgeom_hexgrid.o \
 	lwgeom_export.o \
 	lwgeom_in_gml.o \
 	lwgeom_in_kml.o \
diff --git a/postgis/lwgeom_hexgrid.c b/postgis/lwgeom_hexgrid.c
new file mode 100644
index 0000000..45b6769
--- /dev/null
+++ b/postgis/lwgeom_hexgrid.c
@@ -0,0 +1,354 @@
+/**********************************************************************
+ *
+ * PostGIS - Spatial Types for PostgreSQL
+ * http://postgis.net
+ *
+ * PostGIS is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * PostGIS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with PostGIS.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ **********************************************************************
+ *
+ * Copyright 2009 Mark Cave-Ayland <mark.cave-ayland at siriusit.co.uk>
+ * Copyright 2009-2017 Paul Ramsey <pramsey at cleverelephant.ca>
+ * Copyright 2018 Darafei Praliaskouski <me at komzpa.net>
+ *
+ **********************************************************************/
+
+#include "postgres.h"
+#include "fmgr.h"
+#include "funcapi.h"
+
+#include "../postgis_config.h"
+#include "lwgeom_pg.h"
+#include "liblwgeom.h"
+
+#include <math.h>
+#include <float.h>
+#include <string.h>
+#include <stdio.h>
+
+
+Datum ST_Hexagon(PG_FUNCTION_ARGS);
+Datum ST_Square(PG_FUNCTION_ARGS);
+Datum ST_ShapeGrid(PG_FUNCTION_ARGS);
+
+/* ********* ********* ********* ********* ********* ********* ********* ********* */
+
+/* XXX TO DO: what about negative sizes? */
+
+/* Ratio of hexagon edge length to height = 2*cos(pi/6) */
+#define HXR 1.7320508075688774
+
+static const double hex_x[] = {0.0, 1.0, 1.5,     1.0, 0.0, -0.5,     0.0};
+static const double hex_y[] = {0.0, 0.0, 0.5*HXR, HXR, HXR,  0.5*HXR, 0.0};
+
+static LWGEOM *
+hexagon(double origin_x, double origin_y, double edge, int cell_i, int cell_j, uint32_t srid)
+{
+	double height = edge * HXR;
+	POINT4D pt;
+	POINTARRAY **ppa = lwalloc(sizeof(POINTARRAY*));
+	POINTARRAY *pa = ptarray_construct(0, 0, 7);
+	uint32_t i;
+
+	double offset_x = origin_x + (1.5 * edge * cell_i);
+	double offset_y = origin_y + (height * cell_j) + (0.5 * height * (cell_i % 2));
+	pt.z = pt.m = 0.0;
+
+	for (i = 0; i < 7; ++i)
+	{
+		pt.x = edge * hex_x[i] + offset_x;
+		pt.y = edge * hex_y[i] + offset_y;
+		ptarray_set_point4d(pa, i, &pt);
+	}
+
+	ppa[0] = pa;
+	return lwpoly_as_lwgeom(lwpoly_construct(srid, NULL, 1 /* nrings */, ppa));
+}
+
+static int
+hexagon_grid_size(double edge, const GBOX *gbox, int *cell_width, int *cell_height)
+{
+	double hex_height = edge * HXR;
+	double bounds_width  = gbox->xmax - gbox->xmin;
+	double bounds_height = gbox->ymax - gbox->ymin;
+	int width  = (int)floor(bounds_width/(1.5*edge)) + 1;
+	int height = (int)floor(bounds_height/hex_height) + 1;
+	if (cell_width) *cell_width = width;
+	if (cell_height) *cell_height = height;
+	return width * height;
+}
+
+/* ********* ********* ********* ********* ********* ********* ********* ********* */
+
+static const double sqr_x[] = {0.0, 1.0, 1.0, 0.0, 0.0};
+static const double sqr_y[] = {0.0, 0.0, 1.0, 1.0, 0.0};
+
+static LWGEOM *
+square(double origin_x, double origin_y, double E, int cell_i, int cell_j, uint32_t srid)
+{
+	POINT4D pt;
+	POINTARRAY **ppa = lwalloc(sizeof(POINTARRAY*));
+	POINTARRAY *pa = ptarray_construct(0, 0, 5);
+	uint32_t i;
+
+	double offset_x = origin_x + (E * cell_i);
+	double offset_y = origin_y + (E * cell_j);
+	pt.z = pt.m = 0.0;
+
+	for (i = 0; i < 5; ++i)
+	{
+		pt.x = E * sqr_x[i] + offset_x;
+		pt.y = E * sqr_y[i] + offset_y;
+		ptarray_set_point4d(pa, i, &pt);
+	}
+
+	ppa[0] = pa;
+	return lwpoly_as_lwgeom(lwpoly_construct(srid, NULL, 1 /* nrings */, ppa));
+}
+
+static int
+square_grid_size(double edge, const GBOX *gbox, int *cell_width, int *cell_height)
+{
+	int width  = (int)floor((gbox->xmax - gbox->xmin)/edge) + 1;
+	int height = (int)floor((gbox->ymax - gbox->ymin)/edge) + 1;
+	if (cell_width) *cell_width = width;
+	if (cell_height) *cell_height = height;
+	return width * height;
+}
+
+/* ********* ********* ********* ********* ********* ********* ********* ********* */
+
+typedef enum
+{
+	SHAPE_SQUARE,
+	SHAPE_HEXAGON,
+	SHAPE_TRIANGLE
+} GeometryShape;
+
+typedef struct GeometryGridState
+{
+	int cell_current;
+	int cell_count;
+	int cell_width;
+	int cell_height;
+	GBOX bounds;
+	double size;
+	uint32_t srid;
+	GeometryShape cell_shape;
+}
+GeometryGridState;
+
+/**
+* ST_HexagonGrid(size double default 1.0,
+*            bounds geometry default 'LINESTRING(0 0, 100 100)')
+*/
+PG_FUNCTION_INFO_V1(ST_ShapeGrid);
+Datum ST_ShapeGrid(PG_FUNCTION_ARGS)
+{
+	FuncCallContext *funcctx;
+	MemoryContext oldcontext, newcontext;
+
+	GSERIALIZED *gbounds;
+	GBOX bounds;
+	GeometryGridState *state;
+
+	bool isnull[3] = {0,0,0}; /* needed to say no value is null */
+	Datum tuple_arr[3]; /* used to construct the composite return value */
+	HeapTuple tuple;
+	Datum result; /* the actual composite return value */
+
+	if (SRF_IS_FIRSTCALL())
+	{
+		funcctx = SRF_FIRSTCALL_INIT();
+
+		/* get a local copy of what we're doing a dump points on */
+		gbounds = PG_GETARG_GSERIALIZED_P(1);
+		double size = PG_GETARG_FLOAT8(0);
+
+		int gbounds_is_empty = (gserialized_get_gbox_p(gbounds, &bounds) == LW_FAILURE);
+		double bounds_width = bounds.xmax - bounds.xmin;
+		double bounds_height = bounds.ymax - bounds.ymin;
+
+		/* quick opt-out if we get nonsensical inputs  */
+		if (size <= 0.0 || gbounds_is_empty ||
+		    bounds_width <= 0.0 || bounds_height <= 0.0)
+		{
+			funcctx = SRF_PERCALL_SETUP();
+			SRF_RETURN_DONE(funcctx);
+		}
+
+		newcontext = funcctx->multi_call_memory_ctx;
+		oldcontext = MemoryContextSwitchTo(newcontext);
+
+		/* Create function state */
+		state = lwalloc(sizeof(GeometryGridState));
+		state->cell_current = 0;
+		state->bounds = bounds;
+		state->srid = gserialized_get_srid(gbounds);
+		state->size = size;
+
+		const char *func_name = get_func_name(fcinfo->flinfo->fn_oid);
+		if (strcmp(func_name, "st_hexagongrid") == 0)
+		{
+			state->cell_shape = SHAPE_HEXAGON;
+			state->cell_count = hexagon_grid_size(size, &state->bounds,
+			                                      &state->cell_width, &state->cell_height);
+		}
+		else if (strcmp(func_name, "st_squaregrid") == 0)
+		{
+			state->cell_shape = SHAPE_SQUARE;
+			state->cell_count = square_grid_size(size, &state->bounds,
+			                                     &state->cell_width, &state->cell_height);
+		}
+		else
+		{
+			ereport(ERROR, (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
+				errmsg("%s called from unsupported functional context '%s'", __func__, func_name)));
+		}
+
+		funcctx->user_fctx = state;
+
+		/* get tuple description for return type */
+		if (get_call_result_type(fcinfo, 0, &funcctx->tuple_desc) != TYPEFUNC_COMPOSITE)
+		{
+			ereport(ERROR, (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
+				errmsg("set-valued function called in context that cannot accept a set")));
+		}
+
+		BlessTupleDesc(funcctx->tuple_desc);
+		MemoryContextSwitchTo(oldcontext);
+	}
+
+	/* stuff done on every call of the function */
+	funcctx = SRF_PERCALL_SETUP();
+	newcontext = funcctx->multi_call_memory_ctx;
+
+	/* get state */
+	state = funcctx->user_fctx;
+
+	/* Stop when we've used up all the grid squares */
+	if (state->cell_current >= state->cell_count)
+	{
+		SRF_RETURN_DONE(funcctx);
+	}
+
+	int i = state->cell_current % state->cell_width;
+	int j = state->cell_current / state->cell_width;
+
+	/* Generate geometry */
+	LWGEOM *lwgeom;
+	switch (state->cell_shape)
+	{
+		case SHAPE_HEXAGON:
+			lwgeom = hexagon(state->bounds.xmin, state->bounds.ymin, state->size, i, j, state->srid);
+			break;
+		case SHAPE_SQUARE:
+			lwgeom = square(state->bounds.xmin, state->bounds.ymin, state->size, i, j, state->srid);
+			break;
+		default:
+			ereport(ERROR, (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
+				errmsg("%s called from with unsupported shape '%d'", __func__, state->cell_shape)));
+	}
+
+	tuple_arr[0] = PointerGetDatum(geometry_serialize(lwgeom));
+	lwfree(lwgeom);
+
+	/* Store grid cell coordinates */
+	tuple_arr[1] = Int32GetDatum(i);
+	tuple_arr[2] = Int32GetDatum(j);
+
+	tuple = heap_form_tuple(funcctx->tuple_desc, tuple_arr, isnull);
+	result = HeapTupleGetDatum(tuple);
+	state->cell_current++;
+	SRF_RETURN_NEXT(funcctx, result);
+}
+
+/**
+* ST_Hexagon(size double,
+*            origin geometry default 'POINT(0 0)',
+*            int cellx default 0, int celly default 0)
+*/
+PG_FUNCTION_INFO_V1(ST_Hexagon);
+Datum ST_Hexagon(PG_FUNCTION_ARGS)
+{
+	LWPOINT *lwpt;
+	double size = PG_GETARG_FLOAT8(0);
+	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(1);
+	GSERIALIZED *ghex;
+	int cell_i = PG_GETARG_INT32(2);
+	int cell_j = PG_GETARG_INT32(3);
+	LWGEOM *lwhex;
+	LWGEOM *lworigin = lwgeom_from_gserialized(gorigin);
+
+	if (lwgeom_is_empty(lworigin))
+	{
+		elog(ERROR, "%s: origin point is empty", __func__);
+		PG_RETURN_NULL();
+	}
+
+	lwpt = lwgeom_as_lwpoint(lworigin);
+	if (!lwpt)
+	{
+		elog(ERROR, "%s: origin argument is not a point", __func__);
+		PG_RETURN_NULL();
+	}
+
+	lwhex = hexagon(lwpoint_get_x(lwpt), lwpoint_get_y(lwpt),
+	                size, cell_i, cell_j,
+		            lwgeom_get_srid(lworigin));
+
+	ghex = geometry_serialize(lwhex);
+	PG_FREE_IF_COPY(gorigin, 1);
+	PG_RETURN_POINTER(ghex);
+}
+
+
+/**
+* ST_Square(size double,
+*           origin geometry default 'POINT(0 0)',
+*           int cellx default 0, int celly default 0)
+*/
+PG_FUNCTION_INFO_V1(ST_Square);
+Datum ST_Square(PG_FUNCTION_ARGS)
+{
+	LWPOINT *lwpt;
+	double size = PG_GETARG_FLOAT8(0);
+	GSERIALIZED *gorigin = PG_GETARG_GSERIALIZED_P(1);
+	GSERIALIZED *gsqr;
+	int cell_i = PG_GETARG_INT32(2);
+	int cell_j = PG_GETARG_INT32(3);
+	LWGEOM *lwsqr;
+	LWGEOM *lworigin = lwgeom_from_gserialized(gorigin);
+
+	if (lwgeom_is_empty(lworigin))
+	{
+		elog(ERROR, "%s: origin point is empty", __func__);
+		PG_RETURN_NULL();
+	}
+
+	lwpt = lwgeom_as_lwpoint(lworigin);
+	if (!lwpt)
+	{
+		elog(ERROR, "%s: origin argument is not a point", __func__);
+		PG_RETURN_NULL();
+	}
+
+	lwsqr = square(lwpoint_get_x(lwpt), lwpoint_get_y(lwpt),
+	               size, cell_i, cell_j,
+		           lwgeom_get_srid(lworigin));
+
+	gsqr = geometry_serialize(lwsqr);
+	PG_FREE_IF_COPY(gorigin, 1);
+	PG_RETURN_POINTER(gsqr);
+}
diff --git a/postgis/postgis.sql.in b/postgis/postgis.sql.in
index 1dc1b1e..773ecb6 100644
--- a/postgis/postgis.sql.in
+++ b/postgis/postgis.sql.in
@@ -5991,6 +5991,46 @@ CREATE OR REPLACE FUNCTION ST_InterpolatePoint(Line geometry, Point geometry)
 	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
 	_COST_MEDIUM;
 
+---------------------------------------------------------------
+-- Grid / Hexagon coverage functions
+---
+
+-- Availability: 3.0.0
+CREATE OR REPLACE FUNCTION ST_Hexagon(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
+	RETURNS geometry
+	AS 'MODULE_PATHNAME', 'ST_Hexagon'
+	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	_COST_MEDIUM;
+
+-- Availability: 3.0.0
+CREATE OR REPLACE FUNCTION ST_Square(size float8, origin geometry DEFAULT 'POINT(0 0)', cell_i int4 DEFAULT 0, cell_j int4 DEFAULT 0)
+	RETURNS geometry
+	AS 'MODULE_PATHNAME', 'ST_Square'
+	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	_COST_MEDIUM;
+
+-- Availability: 1.0.0
+CREATE TYPE geometry_grid AS (
+	geom geometry,
+	i integer,
+	j integer
+);
+
+-- Availability: 3.0.0
+CREATE OR REPLACE FUNCTION ST_HexagonGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
+    RETURNS SETOF geometry_grid
+	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
+	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	_COST_MEDIUM;
+
+-- Availability: 3.0.0
+CREATE OR REPLACE FUNCTION ST_SquareGrid(size float8 DEFAULT 1.0, bounds geometry DEFAULT 'LINESTRING(0 0, 100 100)')
+    RETURNS SETOF geometry_grid
+	AS 'MODULE_PATHNAME', 'ST_ShapeGrid'
+	LANGUAGE 'c' IMMUTABLE STRICT _PARALLEL
+	_COST_MEDIUM;
+
+
 -- moved to separate file cause its invovled
 #include "postgis_brin.sql.in"
 

-----------------------------------------------------------------------

Summary of changes:
 doc/html/images/st_hexagongrid01.png | Bin 0 -> 94623 bytes
 doc/html/images/st_hexagongrid02.png | Bin 0 -> 20324 bytes
 doc/html/images/st_hexagongrid03.png | Bin 0 -> 275241 bytes
 doc/reference_constructor.xml        | 146 ++++++++++++
 postgis/Makefile.in                  |   1 +
 postgis/lwgeom_generate_grid.c       | 430 +++++++++++++++++++++++++++++++++++
 postgis/postgis.sql.in               |  39 +++-
 regress/core/regress.sql             |  21 ++
 regress/core/regress_expected        |  16 ++
 9 files changed, 652 insertions(+), 1 deletion(-)
 create mode 100644 doc/html/images/st_hexagongrid01.png
 create mode 100644 doc/html/images/st_hexagongrid02.png
 create mode 100644 doc/html/images/st_hexagongrid03.png
 create mode 100644 postgis/lwgeom_generate_grid.c


hooks/post-receive
-- 
PostGIS


More information about the postgis-tickets mailing list