Correcting GPS Coordinates

Today I was given an excel sheet containing an X and Y column of GPS coordinates that someone in another department had attempted to display in ArcMap. They were successful in that they got some points to display, but these points were located somewhere in the Atlantic, I was asked to 'fix' it.

In these situations more information is always useful, but none was forthcoming, so I made an educated guess. Often when GPS units capture data in British National Grid, they miss out the first number of each X and Y value as these correspond to the map sheet identifying letter and rarely change for the average user. These letters have numerical equivalents, so I tacked on the correct numbers to each coordinate pair and displayed them in Arc. Lo and behold, they were now in the right place.

I decided to knock together a quick Python script to automate this process, which gave me the opportunity to learn how to interface with excel sheets directly (I would usually export to *.csv or *.dbf) using the xlrd module and refresh my memory about the idiosyncrasies of the wonderful pyshp module. Below is the code I wrote in an hour or so. Originally I was going to write a small interface for it and pass it on to the department, but that is a task for another day. It's a bit ugly, but it does what is needed.

Square Buffers With FME

Using FME has a lot of advantages, however every now and then a function which comes as standard in Quantum or Arc is not available. Today I was trying to clip a raster image of an OS map tile to a buffer of the minimum bounding rectangle of my data, so that I have some context map surrounding my areas of interest.

I knocked together what I assumed would do the trick at the end of my workflow, using the Bounding Box Accumulator to get a bounding box, the Bufferer to this polygon buffer by 500 metres and a clipper to clip the map tile to the buffered polygon.

This resulted in a rounded map output, which although not incorrect, looked a bit odd to me. Tinkering with the options of the Bufferer did not yield any different buffer shapes, so I began to think creatively. Using square end caps, it is possible to create rectangular buffers around lines, so the objective became to explode the bounding box out into its individual lines.

I achieved this through the TopologyBuilder transformer which allows an area to be input and a line output, and by setting the Maximum Coords Per Line parameter to 2, four individual lines are created. These lines could then be buffered to 500 metres each using square end caps to create 4 overlapping rectangles which covered my desired buffer extent.

Blue lines are the individual buffers and the red polygon is the original bounding box.

The final stage was to aggregate these individual polygons into one larger shape, which can then be used as the clipper to cut out the section of OS map I required. The final workflow is more convoluted than the original plan, but it works without any problems to give me the square buffer I required.

Update: I noticed an error with this workbench where on larger areas the centre of the clipping polygon is not filled in, creating a doughnut. This is easily fixed by connecting the TopologyBuilder's area output into the aggregator, as the orginal bounding box will always cover the centre of the tile.

Python's Null Values

A problem which I and many others have encountered when using python for GIS tasks is the identification of python's different ways of representing an empty value, generally described in other languages as NULL.

For example, in SQL, to select a series of rows from a table where an attribute is NULL is performed using the command:

SELECT * FROM [TABLE_NAME] WHERE [FIELD] IS NULL

This query will return any row in which FIELD is empty, regardless of why it is empty.

In python things are a bit more complex, there are three possible types of null which python can deal with and knowing which null type you are using is important in avoiding seemingly inexplicable errors from cropping up.

The three types are:

1. None - the python equivalent of SQL's NULL
2. "" - an empty string
3. An unassigned variable

These three types can be demonstrated at the python shell:

Clearly an unassigned variable is useless as it generates an error, but the x and y values can both be used, but are they equal?

>>> x == y
False
>>> 

They may conceptually be the same thing, but the interpreter considers them to be distinct, which can present problems when trying to set loop conditions, or run models until a null value is filled or encountered. The key is to be clear at the outset which kind of null values you will be encountering.

Welcome

I have set this blog up as a means for me to keep track of code I write, write small how to guides for my own reference, list useful tools I come across and generally archive interesting things. If its of any use to other people then that's a bonus.

I expect to cover topics including:

• Python
• General GIS
• FME
• arcpy
• LiDAR
• Geomorphology
• ArcGIS and ArcPad
• MapGuide