House of the rising Moonglum

Have you ever had an excess of GIS vector data, that you just need to simplify? GRASS too confusing? Random CLI tools not fitting the bill? Cartoweb (makers of many fine web GIS tools) has built just the thing, a vertex simplification tool that does all of the hard work of point correspondence for you. Well it turns out, that it actually has a bunch of other programs do that hard work for you, but hey, reuse is king. Unfortunately the instructions are a bit lacking, and I used many online resources to figure out how to do this, so I thought I’d give a little back with instructions for getting all the tools together.

The instructions for the tool say that you need the following:

• a Python interpreter.
• a Postgres Database and the Postgis extension.
• a Python API for Postgres : psycopg
• a lex/yacc Python module : PLY

Welcome to an exciting adventure, we will look at each of those in order. Fortunately OS X already has a python interpreter, so the first item is complete. There are a number of ways to install postgres. Grass comes with a handy installer for both postgres and postgis unfortunately the versions it installs are 7.4.1 and the tool uses features only found in 8. This led me to search elsewhere. Marc Liyanage has a nice prebuilt version for 10.4, or you can use the instructions provided by apple to install it. Since you will need to download the code anyway, I am going to use that method. (See their web page for more complete instructions, I combined those with those from the GRASS site to get the below.)

Download the tarball. 8.1.4 is the latest version at the time of this writing.
Double click on it to decompress it.
> sudo sh
> cd postgres-8.1.4
> ./configure –with-openssl –with-pam –with-krb5 –with-perl –enable-thread-safety –with-includes=/sw/include/ –with-libraries=/sw/lib
> make
> make install
> chown -R username /usr/local/pgsql
> exit

Now you should have postgresql installed. Next is postgis, but lets take a little breather though and make sure our Postgres installation is working.

> /usr/local/pgsql/bin/initdb -D /usr/local/pgsql/data
> /usr/local/pgsql/bin/pg_ctl -D /usr/local/pgsql/data -l /usr/local/pgsql/postgres.log start -o -i
> /usr/local/pgsql/bin/createdb test
> /usr/local/pgsql/bin/psql test

Heck, lets add a user while we are at it:

test=# create user postgres with password ‘p0stgr3s’;
test=# \q

Before we can install postgis, we need to install geos and proj:

Download geos.
> cd ../geos-2.2.3
> ./configure –man-dir=/usr/local/share/man
> make
> sudo make install-strip

Then you need proj4.

Download the tarball. Don’t forget the extra datum grid either.
> cd ../proj-4.4.9
> cp ../proj-datumgrid-1.3/* nad/
> ./configure –mandir=/usr/local/share/man
> make
> sudo make install
> sudo cp -f nad/nzgd2kgrid0005.gsb /usr/local/share/proj/

Now we are finally ready to install Postgis. This is a little tricky, you may want to check out the instructions again. You need to make two changes to the existing Makefile.config.in. Set USE_PROJ=1 and USE_GEOS=1. Then set PROJ_DIR=/usr/local and GEOS_DIR=/usr/local. Finally set USE_ICONV=1. Then open loader/Makefile and set override CFLAGS += -DUSE_ICONV -liconv.

Download the tarball. 1.1.3 as of this writing.
> cd ../postgis-1.1.3
> export PGSQL_SRC=../postgresql-8.1.4/src/
> make
> sudo make install

Now you need to stick the new functionality into your database:

> /usr/local/pgsql/bin/createdb gisdata
> /usr/local/pgsql/bin/createlang plpgsql gisdata
> /usr/local/pgsql/bin/psql -d gisdata -f /usr/local/pgsql/share/contrib/lwpostgis.sql
> /usr/local/pgsql/bin/psql -d gisdata -f /usr/local/pgsql/share/contrib/spatial_ref_sys.sql

OK, now you have postgres and postgis running on your system. Now you need to get all the different layers of map data into your database. I’m going to assume your maps are in .shp format. If not, download GRASS and use the whatever2ogr command line tools to get it in the correct format. Then do:

> /usr/local/bin/shp2pgsql -w layer.shp layer | /usr/local/pgsql/bin/psql gisdata

for each of your data layers. You should probably use: alter table table owner to postgres; to change the owner to be your postgres owner so you can log in. At this point you can check to see if you did everything correctly by grabbing qgis and connecting to your postgis server and downloading a couple layers.

OK, two steps down, next is psycopg. Unfortunately this is a bit tricky to install. An overview is given here. We can’t just dig in though, because there are some header files that we sill need. Also we need the mxDateTime python extension, so…

Download mxBase
> cd ../egenix-mx-base-2.0.6/
> sudo python setup.py install
> cp -R ../postgresql-8.1.4/src/include/catalog /usr/local/pgsql/include/
Download psycopg.
> cd ../psycopg-1.1.20/
> ./configure –with-postgres-libraries=/usr/local/pgsql/lib –with-postgres-includes=/usr/local/pgsql/include –with-mxdatetime-includes=/sw/lib/python2.4/site-packages/mx/DateTime/mxDateTime
> make MACOSX_DEPLOYMENT_TARGET=10.4
> sudo make install

OK, we are getting there. Now we need to do PLY. Fortunately this is just the usual download and install:

Download PLY.
> cd ../ply-1.8
> python setup.py install

Woot! We have all of the precursors installed! Now we need a smoke.

So we should be good to go. Only one problem: the helpful instructions don’t really give you a good example of input, so simplified: install vertex simplifier, build the vertex list, simplify that list. Note that you need to fully qualify the layer name or won’t work. Also the layer must be the layer with the most points if you want to use the same simplification table for multiple layers.

Download vertex simplification program
> cd ../Vertex-Simplification-Program-2.0
> /usr/local/pgsql/bin/psql -f func.sql gisdata
> python vertex.py -H localhost -u postgres -p p0stgr3s -d gisdata -t ‘public.layer’ -v ‘public’ -i ‘gid’ -g ‘the_geom’ -m 1

M is the minimum distance, make it higher for more reduction in size. Go have another smoke. This is going to take a while. It does about one poly every couple seconds. 3,000 polys… 50 minutes! I’ll be back later.

OK, that was a nice break, now you can actually generate the simplified data for each of your layers:

> python simplify.py -H localhost -u postgres -p p0stgr3s -d gisdata -t ‘public.layer’ -v ‘public.vertex1′ -i ‘gid’ -g ‘the_geom’ -m 1

Do that for each layer, and it will make a new geometry column in your layer. Finally you can export your new data using pgsql2shp if you need it in .shp format. Don’t forget to use the -g tag to give it the name of your new data row!

And now you have a nice, simplified map that holds onto its contiguous borders! How cool is that?