Producing GeoJSON with SQL

A couple of weeks ago, I had a need to return GeoJSON from PostGIS. This is something I’ve done many times, but I usually do the final formatting in the application/API layer with Javascript or Python. Basically, my workflow has been to SELECT, using the built-in ST_AsGeoJSON function to convert the geometry to GeoJSON. I then process the recordset to features/feature collections in the application logic.

I like to push such logic as far down as I can, so I finally decided to tackle full FeatureCollection creation in the database. The resulting query makes use of JSON functions in PostgreSQL. I’ll walk through the process with some OSM data I imported for a previous post. Let’s jump in.

It’s best to start with a common table expression (CTE) to do any filter, casts, aliases, or other data manipulation. This is generally good practice anyway, but the process I’ll be describing here doesn’t provide the opportunity to do things such as WHERE or LIMIT farther down in the query. So let’s start with a basic CTE:

WITH geodata AS (
SELECT
	*
FROM
	planet_osm_polygon pop
WHERE
	pop.amenity IS NOT NULL)

I’m starting with polygons that have a value in the “amenity” column. This mainly serves to show, for this post, where to do filtering. It would be a better practice to explicitly name the columns to return but I am choosing not to belabor the CTE for this post. It’s standard SQL that doesn’t need further explanation – but always remember to include the geometry column if you choose to list the columns.

From here forward, we’ll rely on two of PostgreSQL’s JSON functions and the ST_AsGeoJSON function in PostGIS. Here is the full SQL that selects data from the CTE above. I’ll discuss what it’s doing on the other side.

select json_build_object(
'type', 'FeatureCollection',
    'features', json_agg(
        json_build_object(
        	'type', 'Feature',
            'geometry', ST_AsGeoJSON(st_transform(way,4326))::json,
            'properties', json_build_object(
            'osm_id', osm_id,
            'area', way_area,
            'name', "name",
            'amenity', amenity
            )
        )
        )
) AS geojson_featurecollection
FROM geodata;

The workhorse of the SQL above is the json_build_object function. This function creates a JSON object from your data, using a list of key/value pairs you provide. In the example above, I am using string literals for the key names (such as ‘type’), but you could source those dynamically. The function coerces whatever you provide as a key name into text.

The value side is any valid SQL value (a literal, a function return, a column) which will be converted based on the behavior of the “to_json” function for its data type. The query above uses all three options. Where the GeoJSON “type” is required, the literal “FeatureCollection” or “Feature” is provided to conform to the GeoJSON spec. There are three functions embedded to do some data manipulation as well – ST_AsGeoJSON to convert the geometry to GeoJSON. This function, in turn, wraps the ST_Transform function to convert the geometry to WGS84, per the GeoJSON spec. The query also uses the json_agg function which I will discuss next.

If you are familiar with the GeoJSON spec, you can see how the json_build_object is working to create the necessary GeoJSON objects. It is used in three places. In the outermost use, it builds the FeatureCollection object. In the middle use, it builds each feature in the collection. In the innermost use, it creates the properties collection for each feature.

Notice that the json_agg function is used to create the feature. This function simply creates a JSON array of the list of objects it is provided. In the query above, this is where the row set from the CTE is processed into features. It is here that I am choosing which columns to return. If you optimize the CTE, you won’t throw data on the floor during this process. Let’s take a look at the full query, followed by what the output looks like in QGIS:

WITH geodata AS (SELECT * FROM planet_osm_polygon pop WHERE pop.amenity is not null LIMIT 2)
SELECT json_build_object(
'type', 'FeatureCollection',
    'features', json_agg(
        json_build_object(
        	'type', 'Feature',
            'geometry', ST_AsGeoJSON(st_transform(way,4326))::json,
            'properties', json_build_object(
            'osm_id', osm_id,
            'area', way_area,
            'name', "name",
            'amenity', amenity
            )
        )
        )
) AS geojson_featurecollection
FROM geodata;

The data returned from this query is at the end of this post. Note that I added a LIMIT clause to make the output smaller for the post. If you paste that output into a file and load it up in QGIS, you should see something like the map above.

This technique helped me push the GeoJSON formatting logic down into the database, where I prefer that it live. It will help me reduce the surface area of my application code a little, which is always a good thing. One of the main advantages to having a fully-functioning GIS at the database level is that it is a lot easier to mix spatial logic with other data-handling, such as JSON manipulation in this case, down at the tier where data management should occur.

As GIS dissolves into the matrix, there will be more opportunities to leverage proper separation of concerns in system design and move away from needing dedicated resources to analyze and process spatial data.

{
   "type":"FeatureCollection",
   "features":[
      {
         "type":"Feature",
         "geometry":{
            "type":"Polygon",
            "coordinates":[
               [
                  [
                     -76.6386504,
                     38.3217006
                  ],
                  [
                     -76.6378594,
                     38.3209138
                  ],
                  [
                     -76.637446,
                     38.3206796
                  ],
                  [
                     -76.6369433,
                     38.3209459
                  ],
                  [
                     -76.6381707,
                     38.3219188
                  ],
                  [
                     -76.6386504,
                     38.3217006
                  ]
               ]
            ]
         },
         "properties":{
            "osm_id":995969469,
            "area":13493.2,
            "name":"Old Saint Aloysius Cemetery",
            "amenity":"grave_yard"
         }
      },
      {
         "type":"Feature",
         "geometry":{
            "type":"Polygon",
            "coordinates":[
               [
                  [
                     -76.6342374,
                     38.2983132
                  ],
                  [
                     -76.6317095,
                     38.2978559
                  ],
                  [
                     -76.6317474,
                     38.2977412
                  ],
                  [
                     -76.6315868,
                     38.2976973
                  ],
                  [
                     -76.6314641,
                     38.2981277
                  ],
                  [
                     -76.6316157,
                     38.2981561
                  ],
                  [
                     -76.6315936,
                     38.2982684
                  ],
                  [
                     -76.6341184,
                     38.2987296
                  ],
                  [
                     -76.6342374,
                     38.2983132
                  ]
               ]
            ]
         },
         "properties":{
            "osm_id":995969244,
            "area":18525,
            "name":null,
            "amenity":"parking"
         }
      }
   ]
}