08 Demo: Trajectory Analysis#

UW Geospatial Data Analysis
CEE467/CEWA567
David Shean

Pandas Time Series#

Topics to cover#

  • Ingesting and DateTime Index

  • Resampling

  • Indexing

import os
import pandas as pd
import geopandas as gpd

import hvplot.pandas

datadir = './trajectory_data'
fn1=os.path.join(datadir, 'Shean_OASfM_20130909_PosnPnt.csv')
fn2=os.path.join(datadir, 'Shean_OASfM_20130910_PosnPnt.csv')

df = pd.read_csv(fn1)

df

#This is a string
df['GPS Date'].iloc[0]

pd.to_datetime(df['GPS Date'])

pd.to_datetime(df['GPS Date']).iloc[0]

pd.to_datetime(df['GPS Date'] + ' ' + df['GPS Time'])

#read_csv can actually run `to_datetime` if provided with appropriate column names
#pd.read_csv?

df = pd.read_csv(fn1, parse_dates=[['GPS Date', 'GPS Time']])

df.columns

df

df.drop(columns=['Latitude.1', 'Longitude.1'], inplace=True)

df.rename(columns={'GPS Date_GPS Time':'GPS DateTime'}, inplace=True)

df.set_index('GPS DateTime', inplace=True)

df.head()

df.index

#Define function to do the above wrangling
def parse_trimblexh_export(fn):
    df = pd.read_csv(fn, parse_dates=[['GPS Date', 'GPS Time']])
    df.drop(columns=['Latitude.1', 'Longitude.1'], inplace=True)
    df.rename(columns={'GPS Date_GPS Time':'GPS DateTime'}, inplace=True)
    df.set_index('GPS DateTime', inplace=True)
    return df

df1 = parse_trimblexh_export(fn1)
df2 = parse_trimblexh_export(fn2)

print(df1.shape)
print(df2.shape)

df = pd.concat([df1, df2])

df

df1.columns

#%matplotlib widget

df1['HAE'].plot();

df1['HAE'].hvplot()

df.hvplot()

df['2013-09-10 20:26:57':'2013-09-10 20:26:59']

df.loc['2013-09-09']

gdf1 = gpd.GeoDataFrame(df1, crs='EPSG:4326', geometry=gpd.points_from_xy(df1['Longitude'], df1['Latitude']))

gdf = gpd.GeoDataFrame(df, crs='EPSG:4326', geometry=gpd.points_from_xy(df['Longitude'], df['Latitude']))

gdf1

gdf1.hvplot(aspect='equal')

gdf1_utm = gdf1.to_crs('EPSG:32610')

gdf1_utm.hvplot(aspect='equal', c='HAE')

#https://stackoverflow.com/questions/40452759/pandas-latitude-longitude-to-distance-between-successive-rows

gdf1_utm.head()

Horizontal Velocity (ground speed)#

  • Could also look at 3-D vector displacement, vertical velocity, etc.

gdf1_utm.distance(gdf1_utm.shift(1))

gdf1_utm['dm'] = gdf1_utm.distance(gdf1_utm.shift(1))

gdf1_utm['dm'].mean()

gdf1_utm.index

gdf1_utm.index.to_series().diff()

gdf1_utm['dt'] = gdf1_utm.index.to_series().diff().dt.total_seconds()

gdf1_utm['dt']

gdf1_utm['v'] = gdf1_utm['dm']/gdf1_utm['dt']

gdf1_utm['v'].hvplot()

gdf1_utm['a'] = gdf1_utm['v'].diff()

from geoviews import tile_sources as gvts
#map_tiles = gvts.StamenTerrain
map_tiles = gvts.EsriImagery

map_tiles

#kw = {'width':600, 'height':600, 'data_aspect':1, 'alpha':1.0}
kw = {'width':400, 'height':400, 'data_aspect':1, 'alpha':1.0}
kw['colorbar'] = True
kw['cmap'] = 'inferno'
#kw['hover'] = False
#kw['datashade'] = True

Overlay points on map tiles (subset of points)#

  • To combine layers in a single plot using holoviews, you use the asterisk (*) to overlay the two objects

    • You can use the plus sign (+) to build a layout with two separate subplots

  • Currently need to use the geo=True option and revert back to the glas_gdf GeoDataFrame with lat/lon geometry, not our reprojected points in glas_gdf_aea.

    • In principle, hvplot/geoviews should work with our projected points, but there are some residual issues, and it appears that some underlying code somewhere in the stack is assuming lat/lon

map_tiles * gdf1_utm[::3].to_crs('EPSG:3857').hvplot(c='HAE', size=1, **kw)

gdf1_utm[::3].shape

Using datashader to efficiently render all points#

map_tiles * gdf1_utm.to_crs('EPSG:3857').hvplot(c='HAE', datashade=True, dynamic=True, **kw)

Other useful analysis#

  • Distance from starting point

  • Time spent in different elevation bins

  • Altitude above ground - sample DEM

  • Geotagging photos based on timestamp