Home arrow Site Info arrow Latest Portal Content arrow PYTHON! derive speed from GPS
PYTHON! derive speed from GPS PDF Print E-mail

Written by Kevin Bell,

Salt Lake City uses GPS to help gather travel time data to help with traffic issues. After converting from .gpx to a feature class, I have the X, Y, and timestamp on each point from the track log. The code below will add fields that are populated with the number of seconds between points, as well as the speed between a point and it's preceeding point.

This code is kind of fun. Calling the addGPSdata function on the feature class first checks for the needed additional data fields and adds them if needed. Next a dictionary of FID, X, Y, Time is built up in memory for accessing later. A cursor is created to update all of the rows in the feature class and for each row, the time and distance is compared and written into the appropriate field.

NOTE:  all of the time field data's parsing is hardcoded and would need to be altered for a specific time format, and the cursors are hardcoded as well

I duct taped a GPS to my dog and he can run 42 MPH! 

######################################################### 

#Name:     deriveSpeedFromGPS.py
#Author:    Kevin Bell
#Date:       20080108
#Purpose:  add speed/seconds to gps data

import time
import arcgisscripting

#--------------------------------------
def _getTimeCode(inFID,d):  # TIP: google for python string slicing!
    '''Parse a time & return seconds from epoch'''
    t = d[inFID][2]
    y = int(t[:4])
    m = int(t[5:7])
    day = int(t[8:10])
    h = int(t[11:13])
    minute = int(t[14:16])
    sec = int(t[17:19])
    timestamp = y,m,day,h,minute,sec,-1, -1, -1
    timeCode = time.mktime(timestamp)
    return timeCode

def _getDistanceFromPreceeding(inFID,d):  # TIP: google for python dictionaries
    '''return the distance from the point that preceeds the current point'''
    x1 = d[inFID][0]
    y1 = d[inFID][1]
    inFID2 = inFID - 1
    x2 = d[inFID2][0]
    y2 = d[inFID2][1]
    xDelta = abs(x2 - x1)
    yDelta = abs(y2 - y1)
    return round(xDelta + yDelta, 2)

def _buildXYTdict(fc):  # TIP: see the GP cursor help topics
    '''Build a dict of XYT's needed later to compare point against point'''
    d = {}
    cur = gp.SearchCursor(fc)
    row = cur.Next()
    while row:
        d[row.FID] = [row.Y_PROJ, row.X_PROJ, row.TIME]
        row = cur.Next()
    del cur
    return d

def _checkFields(fc, field):
    '''check for fields needed to hold new data, add them if necessary '''
    fldList = gp.ListFields(fc, field)
    if fldList.Next() == None:
        gp.AddField_management(fc, field, "FLOAT", 9, 2)
        print "added field: " + field


def addGPSdata(fc):
    '''add SpeedMPH and seconds to FC... this is relative to its preceeding point'''
    _checkFields(fc, "seconds")
    _checkFields(fc, "SpeedMPH")
    d = _buildXYTdict(fc)
   
    updcur = gp.UpdateCursor(fc)
    updrow = updcur.Next()
    while updrow:
        theFID = updrow.FID
        if theFID != 0: # if it's the 1st point, there's nothing to compare to...
            y = theFID - 1
            curTC = _getTimeCode(theFID, d)
            prevTC = _getTimeCode(y,d)
            secs = curTC - prevTC;                                                   print str(secs) + " seconds"
            meters = _getDistanceFromPreceeding(theFID, d);      print str(meters) + " meters"
            speed = (meters/secs) * 2.23693629;                           print str(round(speed,2)) + " mph" + "\n"

            updrow.seconds = secs
            updrow.SpeedMPH = round(speed,2)
            updcur.UpdateRow(updrow)
            updrow = updcur.Next()
        else:
            updrow = updcur.Next()
    del updcur

#----------------------------------------------------------


gp = arcgisscripting.create()
gp.workspace = r"C:\dog"
fc = "dogNAD83_UTM12N.shp"

addGPSdata(fc);                     print "...added GPS data"

del gp
print "...all done."


Users' Comments  
 

No comment posted

Add your comment

10, Jan. 2008
Last Updated ( 10, Jan. 2008 )
 
< Prev   Next >

AGRC Contacts | UGIC Contacts

feed image feed image

Utah GIS Portal © 2009 AGRC

Optimized for