SQL Geography and GEOJSON

Recently I’ve done a chunk of work with UK mapping of statistical data, which has involved a lot of working with the Geography data type as well as extracting from GEOJSON files.

This has then been plotted in Qlik using a variation of Brian Munz’s Qlikmap which I tweaked to allow me to plot two data layers:

qlikmap

See the git WIP here: QuickMap-QV11

Here is an example of one of the GEOJSON’s from the ONS which can be downloaded from here: ONS Data download

You’ll need SQL16 to be able to use the following code but it works very well indeed:

Select objectid, cty15cd,cty15nm,coordinates,replace(replace(replace(coordinates,'[[[',''),'],[','|'),']]]','') as poly_coordinates
FROM (

SELECT objectid, cty15cd,cty15nm,coordinates

--into mymappingtable
from OPENJSON(@geojson, '$.features')
WITH (
ObjectId int '$.properties.objectid', --NOTE! Case sensitive!!!
cty15cd nvarchar(200) '$.properties.cty15cd',
cty15nm nvarchar(max) '$.properties.cty15nm',
coordinates nvarchar(max) '$.geometry.coordinates' AS JSON )
) a

This will give you the detail from the JSON as well as giving you a set of poly_coordinates.

Sometimes you’re going to end up with a multiple poly shape (a collection of islands for example) so you might need to expand this out a little. The code below will create new rows of data for every multipolygon so that you end up with a clean set of co-ordinates that can be mapped. Excuse any abnormalities or oversights in the code but I have adapted it on the fly to work with the above example

;with tmp(objectid,cty15cd,cty15nm,originalcoordinates, DataItem, coordinates) as (
select objectid,cty15cd,cty15nm,coordinates, LEFT(coordinates, CHARINDEX(']]],[[[',coordinates+']]],[[[')-1),
STUFF(coordinates, 1, CHARINDEX(']]],[[[',coordinates+']]],[[['), '')
from mymappingtable
union all
select objectid,cty15cd,cty15nm,originalcoordinates, LEFT(coordinates, CHARINDEX(']]],[[[',coordinates+']]],[[[')-1),
STUFF(coordinates, 1, CHARINDEX(']]],[[[',coordinates+']]],[[['), '')
from tmp
where coordinates > ''
)
select
replace(cast(objectid as varchar)+'-'+cast((ROW_NUMBER() OVER (PARTITION BY objectid ORDER BY objectid DESC ))-1 as varchar),'-0','') as objectid
,cty15cd
,replace(cty15nm+'-'+cast((ROW_NUMBER() OVER (PARTITION BY objectid ORDER BY objectid DESC ))-1 as varchar),'-0','') as cmlad11nm
,originalcoordinates
,REPLACE(REPLACE(REPLACE(REPLACE(DataItem,']],[[[',''),'],[','|'),'[',''),']','') as poly_coordinates

from tmp
order by objectid
option (maxrecursion 0)

The following bits should work with any SQL version >= 2008 as this is when the geography data type was implemented. This is really easy to work with when you want to check intersects or distances.

To convert the list of polygons to a geography data type you can use (Where poly_coordinates is our column name):

GEOGRAPHY::STGeomFromText(poly_coordinates.MakeValid().STUnion(poly_coordinates.STStartPoint()).STAsText(), 4326)

Testing the as the crow flies distance between two points:

p.geography.STDistance(c.jobgeography)/1609.344)

Leave a Reply