More with the GDAL/OGR perl bindings

So my last post talked about my struggles to finally get something saved in the database using the native perl bindings into the GDAL/OGR library. Once I got that working and pushed out the post, I immediately started loading up multiple files and playing around with the data. One thing I noticed was that it was impossible to separate different “trips” within the data without playing around with space and time. What I wanted was an easy way to flag each batch of points with a field identifying the run.

The auto-generated schema for the GPX data looks like this:

d testogr.track_points
                                              Table "testogr.track_points"
       Column       |           Type           |                               Modifiers                                
--------------------+--------------------------+------------------------------------------------------------------------
 ogc_fid            | integer                  | not null default nextval('testogr.track_points_ogc_fid_seq'::regclass)
 wkb_geometry       | geometry(Point,4326)     | 
 track_fid          | integer                  | 
 track_seg_id       | integer                  | 
 track_seg_point_id | integer                  | 
 ele                | double precision         | 
 time               | timestamp with time zone | 
 magvar             | double precision         | 
 geoidheight        | double precision         | 
 name               | character varying        | 
 cmt                | character varying        | 
 desc               | character varying        | 
 src                | character varying        | 
 link1_href         | character varying        | 
 link1_text         | character varying        | 
 link1_type         | character varying        | 
 link2_href         | character varying        | 
 link2_text         | character varying        | 
 link2_type         | character varying        | 
 sym                | character varying        | 
 type               | character varying        | 
 fix                | character varying        | 
 sat                | integer                  | 
 hdop               | double precision         | 
 vdop               | double precision         | 
 pdop               | double precision         | 
 ageofdgpsdata      | double precision         | 
 dgpsid             | integer                  | 
 speed              | double precision         | 
Indexes:
    "track_points_pkey" PRIMARY KEY, btree (ogc_fid)
    "track_points_wkb_geometry_geom_idx" gist (wkb_geometry)

There are three fields that are completely blank: src, desc, and name. I decided to use src to identify the source of the data as the file name it came from.

First I modified my previous program to parse the command line options using Getopt::Long. I don’t use all of its power in this example, but in the past I’ve been well served by starting with that in case the script grows and mutates.

With Getopt::Long, I understand there are ways to input a list of things into the arguments. You can have multiple invocations of the same option, for example, --file mydata.gpx --file moredata.gpx, or you can input them as a comma separated list and follow the recipe in the perldoc for the module. However, I wanted to use a glob, like –file data/*.gpx, so I instead decided to just stick all the files after a double dash on the command line. So really, in the following code, I’m only using Getopt::Long to parse out a –help command! However, it’s there if I need to expand functionality in the future.

use strict;
use warnings;
use Carp;

use Geo::GDAL;
use Data::Dumper;

use Getopt::Long;
use Pod::Usage;

my $man = 0;
my $help = 0;

my @files;

my $result = GetOptions(
    'help|?' => $help,
    ) or pod2usage(2);

pod2usage(-exitval => 0, -verbose => 2) if $help;

@files = @ARGV;
...

With that, I have all of my input files in an array, and I can loop over them and store the filename in the source field in the db by using $new_feature->SetField('src',$_);, as follows:

foreach (@files){

    my $ds = Geo::OGR::Open($_);

    my $layer         = $ds->Layer($layer_name);
    my $feature_count = $layer->GetFeatureCount();
    carp "$layer_name, $feature_count";
    if ( $feature_count < 10 ) {
        next;
    }

    carp "saving $_ to pg";

    # now append each feature
    my $x = 0;
    $pg_layer->StartTransaction();
    while ( my $feature = $layer->GetNextFeature() ) {

        my $new_feature = Geo::OGR::Feature->new($defn);
        $new_feature->SetFrom($feature);

        # write the filename as the src field, for making lines later
        $new_feature->SetField('src',$_);

        my $pgf = $pg_layer->CreateFeature($new_feature);

        $x += 1;
        if ( $x % 128 == 0 ) {
            carp $x;
            # uncomment the following to crash your program
            # $pg_layer->CommitTransaction();
            # StartTransaction() seems to auto commit prior transaction?
            $pg_layer->StartTransaction(); 
            $x = 0;
        }

    }
    if ($x) {
        carp "all done, $x remaining";
        $pg_layer->CommitTransaction(); # this one doesn't crash for some reason
        carp "last transaction committed";
    }
}

That does its magic, and the database now has distinct groups of points. Now if you want to make “lines” out of those points, you can do this in PostGIS:

SELECT ST_MakeLine(wkb_geometry ORDER BY track_seg_point_id ASC) AS linegeom, src
INTO table testogr.lines
FROM testogr.track_points
GROUP BY src;

Et voila

QGIS rendering the new lines table, on top of OSM lines data

QGIS rendering the new lines table, on top of OSM lines data

Of course, that isn’t at all helpful, as I want to see speeds, not just the lines. Next step is to try to figure out how to add a measure to each point, and then collect those (X,Y,M) type points into a line with a measure dimension. I guess that will be my next post.

Advertisements

2 thoughts on “More with the GDAL/OGR perl bindings

  1. Hi thanks for your excellent articles. Finally i have found someone interested in geospatial information using Perl. A few years ago i developed a database for hiking the Bruce Trail in Ontario Canada. It took a lot of manual work since the trail about 900 km long is always changing. The cartographer uses GPX files that he imports into ESRI and creates paper maps for sale by the Bruce Trail Conservancy. So I had to extract these into text files to load my MySQL database. The guts of my code matched up the intersections of the side trails and main trail to allow planning of interesting hikes. The trail has specific parking points so the database allowed for the planning of hikes with doable distances and describable elevation changes. If you are interested i’d like to discuss your efforts further to see if there is any synergy in our objectives.

    Regards
    John

  2. John,

    Your project sounds pretty cool, and certainly more fun collecting the data that my work!

    My application area is street traffic. I’ve done similar operations to those you describe in the past, but usually my geospatial library of choice is PostgreSQL/PostGIS, with more sophisticated geostatistics operations being done in R (kriging, etc). For example, years ago I wrote some code to identify repeated routes and destinations for my travel survey respondents, with the idea of tuning recommended routes to skew towards the types of routes the person has demonstrated doing in the past (if you never get off the freeway, don’t recommend that; if you always get gasoline in wealthy suburbs, don’t offer guidance to a gas station in Watts; etc.)

    I’ve actually never really worked with altitude, if you can believe it.

    This project I’m working on now is a lot more prosaic. We’re evaluating a new device/technology for estimating speeds on surface streets. To get a ground truth, I drove up and down the test road 20 times with my cellphone using OSMTracker. This produced my GPX files that I then loaded into PostgreSQL.

    More than happy to share my code and techniques, but I haven’t pushed this stuff up to github yet. Need to write some tests first.

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s