on Code

Creating a GIS Service Oriented Architecture

Last month we talked about what a service oriented architecture is. This month I’m going to talk about making a service oriented architecture for GIS. While I’m using ASP.NET for the actual services for convenience (you can make web services in a variety of langauges), the GIS back-end of our SOA is free, open source software. Specifically, PostgreSQL with PostGIS.

Let’s take a look at some of the kinds of things we might want to ask or have a GIS service do:

  • What spatial layers are available?
  • What fields are available for a given layer?
  • Query a layer based on it’s attributes.
  • What is the extent of a selected feature or layer?
  • Project some coordinates from one projection to another.
  • Intersect a point with a feature layer and return the results.
  • Intersect a feature or features of one layer with another layer and return the results.
  • Buffer a point a given distance and return intersected features.
  • Buffer a feature
    or features a given distance and return intersected features from that or another layer.
  • Etc.
If you’ve done GIS programming, you know automating tasks like these can be mildly unpleasant. With PostgreSQL and PostGIS, however, it’s a piece of cake. Each of those tasks is actually a single SQL call.

First, let’s load a couple of GIS layers in to PostgreSQL - tax parcels (parcel.shp), voting precincts (baseprec.shp), and polling locations (polling.shp). We’ll look at tax parcels as an example, as the code is basically the same. shp2pgsql and psql are executables included with PostgreSQL and PostGIS.

First, we turn the shape file into a sql file with this line of code. We’re making the name of the postgresql table tax_parcels, the -d options tells PostgreSQL to delete the table if it already exists, and the -s indicates the projection, and the -i option tells it to use standard 32-bit integers.
shp2pgsql d:\pgrefresh\parcel.shp tax_parcels > d:\pgrefresh\parcel.sql -d -s 2264 -i

Then we’ll tell PostgreSQL to load the table and create a spatial index on the geometry field.
psql -d GISData -U username -f d:\pgrefresh\parcel.sql
psql -d GISData -U username -c “create index idx_tax_parcels_geo on tax_parcels using GIST (the_geom GIST_GEOMETRY_OPS)”

We can also add regular indexes on other fields:
psql -d GISData -U username -c “create index idx_tax_parcels_pid on tax_parcels using btree(pid)”

Finally, we tell PostgreSQL to do some maintenance to speed things up.
psql -d GISData -U username -c “vacuum analyze tax_parcels”</
span>

That’s it. Simply do that for each layer to load. It’s relatively to put in a batch file and automate.

Now let’s look at some SQL that will answer some of our questions. First, let’s list all of our layers:
select f_table_name as LAYER_NAME from geometry_columns

Not so bad. How about something more challenging - getting the extent of some spatial features.
select extent(the_geom) as extent from tax_parcels where pid=’11111111’

Not bad either. Extent is a database function from PostGIS. Now let’s look at something obstensibly harder, intersecting the features of one layer with another layer.
select s.precno from voting_precincts as s, tax_parcels as o where s.the_geom && o.the_geom and intersects(s.the_geom, o.the_geom) and o.pid=’11111111’

Here we’re returning the voting precinct for parcel 11111111. The intersects function does what it sounds like it does. The && function makes the query run faster by using
the minimum bounding rectangles of each layer to reduce the set the intersect function has to work on.

These are simple SQL calls are return a standard recordset. All of the queries we talked about can be performed with a single line of SQL.

Now all we need to do is put that SQL in to a web service. First we’ll make a function to connect to PostgreSQL and return a recordset:
Private Function ReturnDataset(ByVal geoSQL As String) As DataSet
Dim conn As NpgsqlConnection = New NpgsqlConnection(“Server=ntgisnet1;Port=5432;User Id=username;Password=password;Database=database;”)
conn.Open()
Dim da As New NpgsqlDataAdapter(geoSQL, conn)
Dim ds As New DataSet(“GEOPROCESSING”)
da.Fill(ds, “QUERY_RESULTS”)
Return ds
End Function

Then we’ll make our feature-on-feature overlay web function:

<WebMethod>_
Public Function FeatureOverlay(ByVal strFields As String, ByVal geoSourceTable As String, ByVal geoOverlayTable As String, ByVal strWhereClause As String) As String
Try
Return ReturnDataset(“select “ & strFields & “ from “ & geoSourceTable & “ as s, “ & geoOverlayTable & “ as o where s.the_geom && o.the_geom and intersects(s.the_geom, o.the_geom) and “ & strWhereClause).GetXml()
Catch ex As Exception
Return ReturnError(Err.Description).GetXml()
End Try
End Function

Of course, you’ll want to add some code to trap SQL injections and the like (a good coding practice even for internal services), but that’s basically it. Clean XML is
returned with an answer to a feature-on-feature overlay question. What fire district am I in? What voting precinct? Am I in a floodplain? What is my zoning designation? All kinds of questions can be answered with this simple service.

Think about what you’ve just done for a minute. How would your customer normally answer this GIS question? It would involve some desktop software ($), a PC with the muscle to run said software ($), training to use GIS desktop software ($), an enterprise database server like SDE and Oracle/SQL Server/etc. ($), and probably someone to write some code ($), for every single instance you want that function (majorly inefficient - $$$$$).

Now look what we have. Any programmer (including non-GIS programmers) can write a simple interface to a web service and answer the customer’s GIS question. No expensive desktop software, no training, no expensive database back end, and the solution written one time and is available enterprise-wide. You can even access
the web service in programs like Word or Excel, making a couple of text boxes and a submit button and putting the answer in a printable form.

We’re currently using this SOA approach throughout the county to extend GIS functionality to traditional IT products. Previously we might write an entire ArcIMS application so a desktop application could open a browser window to answer some GIS questions, which that user may have to transcribe back to the original application manually. Now we can embed that functionality right into the application itself. Users may not even know they are using GIS.

Which, when you come down to it, is how it should be.