Planetary mapping with QGIS Part Three – Creating Shapefiles From Planetary Datasets

In Part One, I covered how to set up the QGIS workspace to correspond to a custom ellipsoid (Mercury in this case), rather than a terrestrial one.

Part Two dealt with creating polygons and using the GDAL/OGR command line to transform the resulting shapefiles into raw coordinate lists.

This time, I’ll cover the process of creating shapefiles from raw data and importing them into QGIS.

Let’s say I have a file of tab-delimited values of   longitude    latitude    elevation
131.1758   -11.6734   1.5
131.1759   -11.6722   1.2
131.1794   -11.4167   0.4
130.3352   -6.20119   -0.3

To get this into meter map-coordinates, use the PROJ.4 utility “proj” to convert them, basically reversing the process from Part Two.

proj +proj=eqc +lat_ts=0 +a=2440000 +b=2440000 +no_defs test_lon_lat_elev.tab | sed 's/\t/,/g' > test_lon_lat_elev_proj.csv

Where test_lon_lat_elev.tab is the input file listed above and test_lon_lat_elev_proj.csv is the comma-delimited output file. Note that I used the sed utility to convert the tabs to commas.

The resulting file should look like this:

5586257.05,-497123.81,1.5
5586261.30,-497072.70,1.2
5586410.36,-486191.97,0.4
5550459.23,-264084.09,-0.3

Now we need to add the header row with field descriptions so that the file becomes this:

Longitude,Latitude,Elevation
5586257.05,-497123.81,1.5
5586261.30,-497072.70,1.2
5586410.36,-486191.97,0.4
5550459.23,-264084.09,-0.3

As far as I can tell, PROJ.4 prefers a tab-delimited input file. If you try to use a comma-delimited one, it fails to properly parse the long/lat/alt values.
OGR, on the other hand, needs comma-delimited files as input, since it uses the concept of “drivers” to parse the file, as will be shown next.

Now that the input is properly formatted, we can use OGR to create the shapefiles.  OGR uses the concept of “virtual formats” to describe the data contained in the input file.  These are basically XML-styles that must be prepared for each file to be converted.  A .csv file will use the following VRT file:

<OGRVRTDataSource>
 <OGRVRTLayer name="test_lon_lat_elev_proj">
 <SrcDataSource>test_lon_lat_elev_proj.csv</SrcDataSource>
 <GeometryType>wkbPoint</GeometryType>
 <LayerSRS>NULL</LayerSRS>
 <GeometryField encoding="PointFromColumns" x="Longitude" y="Latitude"/>
 </OGRVRTLayer>
 </OGRVRTDataSource>

In the above .vrt file, test_lon_lat_elev_proj.csv is the file produced above.
OGRVRTLayer name must apparently match exactly with the SrcDataSource field.

Now, using the OGR command line:
ogr2ogr test test_lon_lat_elev.vrt
Will produce a directory called test containing, among other things a .shp shapefile.

If you still have the Mercury global mosaic loaded in the QGIS workspace, you can click on the Open Shapefile button, navigate to the test directory that was created above, and then select the new .shp file. The result should be something similar to the image below.

Note that it only appears that there are two points.  There are actually four, which can be seen if you zoom in on the cluster far enough to resolve them.