# Creating an irregular point grid with provided row and column numbers in PostGIS

I am looking to create an irregular grid of points via PostGIS, using the four provided corner (point) coordinates, and number of point-rows, and point-columns.

Four corners may or may not be perpendicular to each other (see image, 7 rows, 5 columns). Terrain is considered flat. What might be the easiest approach to solve this, within PostGIS ?

• Can you provide an example image/sketch of what you're trying to achieve?
– Erik
Aug 8, 2022 at 11:38
• Automation of regular grids is pretty basic, bur, by definition, an irregular grid would seem automation-resistant. Aug 8, 2022 at 11:42
– Erik
Aug 8, 2022 at 13:43
• My approach was first creating four sides (ST_InterpolatePoints), then filling the inside, row by row. But not sure if this can be applicable in SQL/PostGIS limits. Aug 9, 2022 at 6:34
• It is still not clear from your question what you would like to have in the output, for example: 1) only geometry from points as in the picture, or 2) geometry from points and numbered columns and rows calculated from one corner (1 or 2 or 3 or 4)...This is important... Aug 10, 2022 at 19:31

So, there are many ways to solve your question and this approach is one of them.

Create a fun function called `ST_RegularPointsGridOfCornerPoints`

``````DROP FUNCTION ST_RegularPointsGridOfCornerPoints

CREATE OR REPLACE FUNCTION ST_RegularPointsGridOfCornerPoints(
geom GEOMETRY,
r bigint,
c bigint)
RETURNS GEOMETRY AS
\$BODY\$
WITH
tbla AS (SELECT ST_Boundary(ST_Union(geom)) geom FROM (SELECT ((ST_DelaunayTriangles(ST_Collect(geom)))) geom) foo),
tblb AS (SELECT row_number() over() AS id,
ST_MakeLine(pt1, pt2) geom FROM (SELECT ST_PointN(geom, generate_series(1, ST_NPoints(geom)-1)) pt1,
ST_PointN(geom, generate_series(2, ST_NPoints(geom))) pt2 FROM tbla) AS geom),
tblc AS (SELECT generate_series (0,r-1) as steps),
tbld AS (SELECT steps AS stp1, ST_LineInterpolatePoint(geom, steps/(SELECT count(steps)::float-1 FROM tblc)) geom1 FROM tblc, tblb WHERE tblb.id
IN (2) GROUP BY tblc.steps, geom),
tble AS (SELECT steps AS stp2, ST_LineInterpolatePoint(ST_Reverse(geom), steps/(SELECT count(steps)::float-1 FROM tblc)) geom2 FROM tblc, tblb
WHERE tblb.id IN (4) GROUP BY tblc.steps, geom),
tblf AS (SELECT row_number() over() AS id, ST_MakeLine(geom1, geom2) geom FROM tbld JOIN tble ON true AND stp1=stp2),
tblg AS (SELECT generate_series (0,c-1) as steps)
(SELECT ST_LineInterpolatePoint(geom, steps/(SELECT count(steps-1)::float-1 FROM tblg)) geom FROM tblg, tblf geom);
\$BODY\$
LANGUAGE SQL
``````

Run

`SELECT ST_RegularPointsGridOfCornerPoints(ST_Union(geom), 7, 5) geom FROM <name_table>`

See the result - Unfortunately something went wrong and it works not on all versions of PostgreSQL builds (For example, for PostgreSQL 14.0, compiled by Visual C++ build 1914, 64-bit and higher should work :-))... Remember my comment, its future hasn't come yet :-(...

As a consequence, for now, run the body of the function as a CTE and set the required values of columns and rows, for example, as specified in your question for your example. The architecture of the SQL-code is shown below:

``````create table <name_table> AS
WITH
tbla AS (SELECT ST_Boundary(ST_Union(geom)) geom FROM (SELECT ((ST_DelaunayTriangles(ST_Collect(geom)))) geom FROM layer_1) foo),
tblb AS (SELECT row_number() over() AS id,
ST_MakeLine(pt1, pt2) geom FROM (SELECT ST_PointN(geom, generate_series(1, ST_NPoints(geom)-1)) pt1,
ST_PointN(geom, generate_series(2, ST_NPoints(geom))) pt2 FROM tbla) AS geom),
tblc AS (SELECT generate_series (0,4) as steps),
tbld AS (SELECT steps AS stp1, ST_LineInterpolatePoint(geom, steps/(SELECT count(steps)::float-1 FROM tblc)) geom1 FROM tblc, tblb WHERE tblb.id
IN (2) GROUP BY tblc.steps, geom),
tble AS (SELECT steps AS stp2, ST_LineInterpolatePoint(ST_Reverse(geom), steps/(SELECT count(steps)::float-1 FROM tblc)) geom2 FROM tblc, tblb
WHERE tblb.id IN (4) GROUP BY tblc.steps, geom),
tblf AS (SELECT row_number() over() AS id, ST_MakeLine(geom1, geom2) geom FROM tbld JOIN tble ON true AND stp1=stp2),
tblg AS (SELECT generate_series (0,6) as steps)
(SELECT ST_LineInterpolatePoint(geom, steps/(SELECT count(steps-1)::float-1 FROM tblg)) geom FROM tblg, tblf);
``````

The figure below shows the result, you should get the same one for yourself...

The figure

Unfortunately, I only fancy fun and customizable functions and they are not always simple :-(...

P.S. In the following questions, try to present the SQL-code and an explanation of what prevented you from getting the expected result...

P.S.

Everything is fine now, the geospatial function works as expected :-)! Create Geo-SQL function:

``````CREATE OR REPLACE FUNCTION ST_RegularPointsGridOfCornerPoints(
geom GEOMETRY,
r bigint,
c bigint)
RETURNS TABLE (geom GEOMETRY) AS
\$BODY\$
WITH
tbla AS (SELECT ST_Boundary(ST_Union(geom)) geom FROM (SELECT ((ST_DelaunayTriangles(ST_Collect(geom)))) geom) foo),
tblb AS (SELECT row_number() over() AS id,
ST_MakeLine(pt1, pt2) geom FROM (SELECT ST_PointN(geom, generate_series(1, ST_NPoints(geom)-1)) pt1,
ST_PointN(geom, generate_series(2, ST_NPoints(geom))) pt2 FROM tbla) AS geom),
tblc AS (SELECT generate_series (0,r-1) as steps),
tbld AS (SELECT steps AS stp1, ST_LineInterpolatePoint(geom, steps/(SELECT count(steps)::float-1 FROM tblc)) geom1 FROM tblc, tblb WHERE tblb.id
IN (2) GROUP BY tblc.steps, geom),
tble AS (SELECT steps AS stp2, ST_LineInterpolatePoint(ST_Reverse(geom), steps/(SELECT count(steps)::float-1 FROM tblc)) geom2 FROM tblc, tblb
WHERE tblb.id IN (4) GROUP BY tblc.steps, geom),
tblf AS (SELECT row_number() over() AS id, ST_MakeLine(geom1, geom2) geom FROM tbld JOIN tble ON true AND stp1=stp2),
tblg AS (SELECT generate_series (0,c-1) as steps)
(SELECT ST_LineInterpolatePoint(geom, steps/(SELECT count(steps-1)::float-1 FROM tblg)) geom FROM tblg, tblf geom);
\$BODY\$
LANGUAGE SQL
``````

Run:

`SELECT ST_RegularPointsGridOfCornerPoints(geom, 7, 5) geom FROM <name_poly_table>`

(-: FOGS :-)...

Translated with www.DeepL.com/Translator (free version)

• This is working exceptionally well, thank you very much ! Aug 13, 2022 at 11:17
• I'm glad I helped you and others and thank you for the question... 8-)... Aug 15, 2022 at 17:24

An elegant way to do this is to transform a square grid of points into the required quadrilateral. Since the transformation does not preserve parallel lines, a projective transformation is normally required. This is complex to derive. But this answer describes a clever and simple transformation of the unit square to an arbitrary quadrilateral. The diagram below shows how it uses a combination of three vectors derived from the vertices of the quadrilateral.

Here's SQL implementing the transformation, with an example quadrilateral:

``````WITH quad AS (SELECT
5 AS LLx, 5 AS LLy,
10 AS ULx, 30 AS ULy,
25 AS URx, 25 AS URy,
30 AS LRx, 0 AS LRy
),
vec AS (
SELECT  LLx AS ox, LLy AS oy,
ULx - LLx AS ux, ULy - LLy AS uy,
LRx - LLx AS vx, LRy - LLy As vy,
URx - LLx - ((ULX - LLx) + (LRx - LLx)) AS wx,
URy - LLy - ((ULy - LLy) + (LRy - LLy)) AS wy