GNSS data processing script

GNSS data processing scriptΒΆ

Note: This notebook is not reproducible because the source data paths direct to a large database. It is meant to give readers a sense about how the provided GNSS records (Kaskawulsh_v2_2018-mm-dd_to_2018-mm-dd_GPS.csv) were made.

import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import pandas as pd
import geopandas as gpd
import os
from datetime import datetime
from scipy import interpolate
import math
import glob
from pyproj import transform, CRS
from math import degrees, atan2
# from scipy.stats import linregress
import scipy.stats
from netCDF4 import Dataset, num2date
import calendar
from scipy.io import loadmat

pd.set_option('display.max_columns', None)

from IPython.core.display import display, HTML
display(HTML("<style>.container { width:100% !important; }</style>"))

def convert_time(seconds): 
    seconds = seconds % (24 * 3600) 
    hour = seconds // 3600
    seconds %= 3600
    minutes = seconds // 60
    seconds %= 60
      
    return "%d:%02d:%02d" % (hour, minutes, seconds)

def custom_mean(df):
    return df.mean(skipna=False)

def rsquared(x, y):
    """ Return R^2 where x and y are array-like."""
    slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(x, y)
    return r_value**2

def truncate(number, decimals=0):
    """
    Returns a value truncated to a specific number of decimal places.
    """
    if not isinstance(decimals, int):
        raise TypeError("decimal places must be an integer.")
    elif decimals < 0:
        raise ValueError("decimal places has to be 0 or more.")
    elif decimals == 0:
        return math.trunc(number)

    factor = 10.0 ** decimals
    return math.trunc(number * factor) / factor

horizontal_error = 0.12
vertical_error = 0.2
# #read Lower data
lower_path='/Users/willkochtitzky/Projects/Kaskawulsh/GPS_data/lower_cleaned.csv'
#read in PPP corrected data (from NR CAN)
#read lower data
lower_pos_data = pd.read_csv(lower_path)#,usecols=['DIR','HGT(m)','YEAR-MM-DD','HR:MN:SS.SS','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)'])

lower_pos_data['date'] = pd.to_datetime(lower_pos_data['YEAR-MM-DD']+'/'+lower_pos_data['HR:MN:SS.SS'],format='%Y-%m-%d/%H:%M:%S.%f')
lower_pos_data['day_date'] = pd.to_datetime(lower_pos_data['YEAR-MM-DD'],format='%Y-%m-%d')
lower_pos_data = lower_pos_data.drop(columns=['DIR','YEAR-MM-DD','HR:MN:SS.SS'])

# lower_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']] = lower_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']].apply(pd.to_numeric)#change the data types
lower_pos_data = lower_pos_data.sort_values(by=['date'])#sort the data

#filter out good data to inspect the bad data later
lower_pos_data_bad_data = lower_pos_data[(lower_pos_data['SDLAT(95%)']>horizontal_error) & (lower_pos_data['SDLON(95%)']>horizontal_error) & (lower_pos_data['SDHGT(95%)']>vertical_error)]  

lower_pos_data = lower_pos_data[lower_pos_data['SDLAT(95%)']<horizontal_error]
lower_pos_data = lower_pos_data[lower_pos_data['SDLON(95%)']<horizontal_error]
lower_pos_data = lower_pos_data[lower_pos_data['SDHGT(95%)']<vertical_error]

lower_pos_data['year'] = lower_pos_data['date'].dt.year
lower_pos_data['day'] = lower_pos_data['date'].dt.dayofyear
lower_pos_data = lower_pos_data.reset_index() #reset the index for the calculations ahead, this is crucial to make sure your counter is right!

#calculate how much time has occured since 1/1/00 so that you can figure out how much time is between observations
pos_sec_since = []
jan1_2007 = datetime(2007,1,1)
for i in range(0,len(lower_pos_data)):
    pos_sec_since.append((lower_pos_data['date'][i]-jan1_2007).total_seconds())
lower_pos_data['sec_since']=pos_sec_since

lower_pos_data_daily = lower_pos_data.set_index('date').groupby(pd.Grouper(freq='d')).mean() #calculate daily data
lower_pos_data_daily = lower_pos_data_daily.reset_index() #reset the index after you claculate daily data

distance=[float('nan')]
obs_duration_days=[float('nan')]
for i in range(1,len(lower_pos_data_daily)):
    distance.append(np.sqrt((lower_pos_data_daily['UTM_EASTING'][i]-lower_pos_data_daily['UTM_EASTING'][i-1])**2+(lower_pos_data_daily['UTM_NORTHING'][i]-lower_pos_data_daily['UTM_NORTHING'][i-1])**2))   
    obs_duration_days.append((lower_pos_data_daily['sec_since'][i]-lower_pos_data_daily['sec_since'][i-1])/3600/24)   
lower_pos_data_daily['Distance']=distance
lower_pos_data_daily['obs_duration_days']=obs_duration_days
lower_pos_data_daily[['Distance','obs_duration_days']] = lower_pos_data_daily[['Distance','obs_duration_days']].apply(pd.to_numeric)#change the data types
lower_pos_data_daily['Velocity_m_per_d']=lower_pos_data_daily['Distance']/lower_pos_data_daily['obs_duration_days']
#read middle data
#middle_path='/Users/willkochtitzky/Projects/Kaskawulsh/GPS_data/Middle_all_data.csv'
middle_path='/Users/willkochtitzky/Projects/Kaskawulsh/GPS_data/middle_cleaned.csv'
#read in PPP corrected data (from NR CAN)
middle_pos_data = pd.read_csv(middle_path)#,usecols=['DIR','HGT(m)','YEAR-MM-DD','HR:MN:SS.SS','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)'])

middle_pos_data['date'] = pd.to_datetime(middle_pos_data['YEAR-MM-DD']+'/'+middle_pos_data['HR:MN:SS.SS'],format='%Y-%m-%d/%H:%M:%S.%f')
middle_pos_data['day_date'] = pd.to_datetime(middle_pos_data['YEAR-MM-DD'],format='%Y-%m-%d')
middle_pos_data = middle_pos_data.drop(columns=['DIR','YEAR-MM-DD','HR:MN:SS.SS'])

# middle_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']] = middle_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']].apply(pd.to_numeric)#change the data types
middle_pos_data = middle_pos_data.sort_values(by=['date'])#sort the data

#filter out good data to inspect the bad data later
middle_pos_data_bad_data = middle_pos_data[(middle_pos_data['SDLAT(95%)']>horizontal_error) & (middle_pos_data['SDLON(95%)']>horizontal_error) & (middle_pos_data['SDHGT(95%)']>vertical_error)]  

middle_pos_data = middle_pos_data[middle_pos_data['SDLAT(95%)']<horizontal_error]
middle_pos_data = middle_pos_data[middle_pos_data['SDLON(95%)']<horizontal_error]
middle_pos_data = middle_pos_data[middle_pos_data['SDHGT(95%)']<vertical_error]

middle_pos_data['year'] = middle_pos_data['date'].dt.year
middle_pos_data['day'] = middle_pos_data['date'].dt.dayofyear
middle_pos_data = middle_pos_data.reset_index() #reset the index for the calculations ahead, this is crucial to make sure your counter is right!

#calculate how much time has occured since 1/1/00 so that you can figure out how much time is between observations
pos_sec_since = []
jan1_2007 = datetime(2007,1,1)
for i in range(0,len(middle_pos_data)):
    pos_sec_since.append((middle_pos_data['date'][i]-jan1_2007).total_seconds())
middle_pos_data['sec_since']=pos_sec_since

middle_pos_data_daily = middle_pos_data.set_index('date').groupby(pd.Grouper(freq='d')).mean() #calculate daily data
middle_pos_data_daily = middle_pos_data_daily.reset_index() #reset the index after you claculate daily data

distance=[float('nan')]
obs_duration_days=[float('nan')]
for i in range(1,len(middle_pos_data_daily)):
    distance.append(np.sqrt((middle_pos_data_daily['UTM_EASTING'][i]-middle_pos_data_daily['UTM_EASTING'][i-1])**2+(middle_pos_data_daily['UTM_NORTHING'][i]-middle_pos_data_daily['UTM_NORTHING'][i-1])**2))   
    obs_duration_days.append((middle_pos_data_daily['sec_since'][i]-middle_pos_data_daily['sec_since'][i-1])/3600/24)   
middle_pos_data_daily['Distance']=distance
middle_pos_data_daily['obs_duration_days']=obs_duration_days
middle_pos_data_daily[['Distance','obs_duration_days']] = middle_pos_data_daily[['Distance','obs_duration_days']].apply(pd.to_numeric)#change the data types
middle_pos_data_daily['Velocity_m_per_d']=middle_pos_data_daily['Distance']/middle_pos_data_daily['obs_duration_days']
#read Upper data
upper_path='/Users/willkochtitzky/Projects/Kaskawulsh/GPS_data/upper_cleaned.csv'
#read in PPP corrected data (from NR CAN)
upper_pos_data = pd.read_csv(upper_path)#,usecols=['DIR','HGT(m)','YEAR-MM-DD','HR:MN:SS.SS','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)'])
upper_pos_data['date'] = pd.to_datetime(upper_pos_data['YEAR-MM-DD']+'/'+upper_pos_data['HR:MN:SS.SS'],format='%Y-%m-%d/%H:%M:%S.%f')
upper_pos_data['day_date'] = pd.to_datetime(upper_pos_data['YEAR-MM-DD'],format='%Y-%m-%d')
upper_pos_data = upper_pos_data.drop(columns=['DIR','YEAR-MM-DD','HR:MN:SS.SS'])

# upper_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']] = upper_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']].apply(pd.to_numeric)#change the data types
upper_pos_data = upper_pos_data.sort_values(by=['date'])#sort the data

#filter out good data to inspect the bad data later
upper_pos_data_bad_data = upper_pos_data[(upper_pos_data['SDLAT(95%)']>horizontal_error) & (upper_pos_data['SDLON(95%)']>horizontal_error) & (upper_pos_data['SDHGT(95%)']>vertical_error)]  

upper_pos_data = upper_pos_data[upper_pos_data['SDLAT(95%)']<horizontal_error]
upper_pos_data = upper_pos_data[upper_pos_data['SDLON(95%)']<horizontal_error]
upper_pos_data = upper_pos_data[upper_pos_data['SDHGT(95%)']<vertical_error]

upper_pos_data['year'] = upper_pos_data['date'].dt.year
upper_pos_data['day'] = upper_pos_data['date'].dt.dayofyear
upper_pos_data = upper_pos_data.reset_index() #reset the index for the calculations ahead, this is crucial to make sure your counter is right!

#calculate how much time has occured since 1/1/00 so that you can figure out how much time is between observations
pos_sec_since = []
jan1_2007 = datetime(2007,1,1)
for i in range(0,len(upper_pos_data)):
    pos_sec_since.append((upper_pos_data['date'][i]-jan1_2007).total_seconds())
upper_pos_data['sec_since']=pos_sec_since

upper_pos_data_daily = upper_pos_data.set_index('date').groupby(pd.Grouper(freq='d')).mean() #calculate daily data
upper_pos_data_daily = upper_pos_data_daily.reset_index() #reset the index after you claculate daily data

distance=[float('nan')]
obs_duration_days=[float('nan')]
for i in range(1,len(upper_pos_data_daily)):
    distance.append(np.sqrt((upper_pos_data_daily['UTM_EASTING'][i]-upper_pos_data_daily['UTM_EASTING'][i-1])**2+(upper_pos_data_daily['UTM_NORTHING'][i]-upper_pos_data_daily['UTM_NORTHING'][i-1])**2))   
    obs_duration_days.append((upper_pos_data_daily['sec_since'][i]-upper_pos_data_daily['sec_since'][i-1])/3600/24)   
upper_pos_data_daily['Distance']=distance
upper_pos_data_daily['obs_duration_days']=obs_duration_days
upper_pos_data_daily[['Distance','obs_duration_days']] = upper_pos_data_daily[['Distance','obs_duration_days']].apply(pd.to_numeric)#change the data types
upper_pos_data_daily['Velocity_m_per_d']=upper_pos_data_daily['Distance']/upper_pos_data_daily['obs_duration_days']
#read Arm data
arm_path='/Users/willkochtitzky/Projects/Kaskawulsh/GPS_data/arm_cleaned.csv'
#read in PPP corrected data (from NR CAN)
arm_pos_data = pd.read_csv(arm_path)#,usecols=['DIR','HGT(m)','YEAR-MM-DD','HR:MN:SS.SS','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)'])
arm_pos_data['date'] = pd.to_datetime(arm_pos_data['YEAR-MM-DD']+'/'+arm_pos_data['HR:MN:SS.SS'],format='%Y-%m-%d/%H:%M:%S.%f')
arm_pos_data['day_date'] = pd.to_datetime(arm_pos_data['YEAR-MM-DD'],format='%Y-%m-%d')
arm_pos_data = arm_pos_data.drop(columns=['DIR','YEAR-MM-DD','HR:MN:SS.SS'])

# arm_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']] = arm_pos_data[['HGT(m)','UTM_EASTING','UTM_NORTHING','GDOP','RMSC(m)','RMSP(m)','SDLAT(95%)','SDLON(95%)','SDHGT(95%)','SDCLK(95%)','SDTZD(95%)']].apply(pd.to_numeric)#change the data types
arm_pos_data = arm_pos_data.sort_values(by=['date'])#sort the data

#filter out good data to inspect the bad data later
arm_pos_data_bad_data = arm_pos_data[(arm_pos_data['SDLAT(95%)']>horizontal_error) & (arm_pos_data['SDLON(95%)']>horizontal_error) & (arm_pos_data['SDHGT(95%)']>vertical_error)]

arm_pos_data = arm_pos_data[arm_pos_data['SDLAT(95%)']<horizontal_error]
arm_pos_data = arm_pos_data[arm_pos_data['SDLON(95%)']<horizontal_error]
arm_pos_data = arm_pos_data[arm_pos_data['SDHGT(95%)']<vertical_error]

arm_pos_data['year'] = arm_pos_data['date'].dt.year
arm_pos_data['day'] = arm_pos_data['date'].dt.dayofyear
arm_pos_data = arm_pos_data.reset_index() #reset the index for the calculations ahead, this is crucial to make sure your counter is right!

#calculate how much time has occured since 1/1/00 so that you can figure out how much time is between observations
pos_sec_since = []
jan1_2007 = datetime(2007,1,1)
for i in range(0,len(arm_pos_data)):
    pos_sec_since.append((arm_pos_data['date'][i]-jan1_2007).total_seconds())
arm_pos_data['sec_since']=pos_sec_since

arm_pos_data_daily = arm_pos_data.set_index('date').groupby(pd.Grouper(freq='d')).mean() #calculate daily data
arm_pos_data_daily = arm_pos_data_daily.reset_index() #reset the index after you claculate daily data

distance=[float('nan')]
obs_duration_days=[float('nan')]
for i in range(1,len(arm_pos_data_daily)):
    distance.append(np.sqrt((arm_pos_data_daily['UTM_EASTING'][i]-arm_pos_data_daily['UTM_EASTING'][i-1])**2+(arm_pos_data_daily['UTM_NORTHING'][i]-arm_pos_data_daily['UTM_NORTHING'][i-1])**2))   
    obs_duration_days.append((arm_pos_data_daily['sec_since'][i]-arm_pos_data_daily['sec_since'][i-1])/3600/24)   
arm_pos_data_daily['Distance']=distance
arm_pos_data_daily['obs_duration_days']=obs_duration_days
arm_pos_data_daily[['Distance','obs_duration_days']] = arm_pos_data_daily[['Distance','obs_duration_days']].apply(pd.to_numeric)#change the data types
arm_pos_data_daily['Velocity_m_per_d']=arm_pos_data_daily['Distance']/arm_pos_data_daily['obs_duration_days']
#calculate velocity between two dates:

# Whyjay project
start_dates = ['2018-05-08', '2018-09-03', '2018-08-18', '2018-08-02', 
               '2018-04-28', '2018-04-12', '2018-05-23', '2018-04-21', 
               '2018-04-05', '2018-03-04', '2018-09-30', '2018-09-20', 
               '2018-09-10', '2018-08-31', '2018-08-11', '2018-08-01', 
               '2018-07-27', '2018-07-22', '2018-06-27', '2018-06-12', 
               '2018-05-23', '2018-05-18', '2018-05-08', '2018-03-29', 
               '2018-03-14', '2018-03-04', '2018-09-17', '2018-08-18', 
               '2018-07-29', '2018-07-24', '2018-07-04', '2018-06-19', 
               '2018-05-15', '2018-03-16', '2018-03-06']
end_dates = ['2018-06-27', '2018-10-05', '2018-09-03', '2018-08-18', 
             '2018-08-02', '2018-04-28', '2018-06-08', '2018-05-23', 
             '2018-04-21', '2018-04-05', '2018-10-05', '2018-09-30', 
             '2018-09-20', '2018-09-10', '2018-08-31', '2018-08-11', 
             '2018-08-01', '2018-07-27', '2018-07-22', '2018-06-27', 
             '2018-06-12', '2018-05-23', '2018-05-18', '2018-05-08', 
             '2018-03-29', '2018-03-14', '2018-10-02', '2018-09-17', 
             '2018-08-18', '2018-07-29', '2018-07-24', '2018-07-04', 
             '2018-06-19', '2018-05-15', '2018-03-16']


for i in range(0,len(start_dates)):
    vel_date_start=pd.to_datetime(start_dates[i],format='%Y-%m-%d')
    vel_date_end=pd.to_datetime(end_dates[i],format='%Y-%m-%d')

    #Lower
    for i in range(0,len(lower_pos_data_daily)):
        if lower_pos_data_daily['date'][i]==vel_date_start:
            start_pos= lower_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]
        if lower_pos_data_daily['date'][i]==vel_date_end:
            end_pos=lower_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]

    distance_traveled = np.sqrt((start_pos['UTM_EASTING']-end_pos['UTM_EASTING'])**2+(start_pos['UTM_NORTHING']-end_pos['UTM_NORTHING'])**2)
    time_between_obs = (end_pos['sec_since']-start_pos['sec_since'])/3600/24

    velocity_obs_lower=pd.DataFrame({'date1':[vel_date_start],'date2':[vel_date_end],'start_easting':[start_pos['UTM_EASTING']],'start_northing':[start_pos['UTM_NORTHING']],'end_easting':[end_pos['UTM_EASTING']],'end_northing':[end_pos['UTM_NORTHING']],'distance_traveled (m)':[distance_traveled],'velocity (m/d)':[distance_traveled/time_between_obs]},index=[0])
    velocity_obs_lower['label']='lower'
    
    #Middle
    for i in range(0,len(middle_pos_data_daily)):
        if middle_pos_data_daily['date'][i]==vel_date_start:
            start_pos= middle_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]
        if middle_pos_data_daily['date'][i]==vel_date_end:
            end_pos=middle_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]

    distance_traveled = np.sqrt((start_pos['UTM_EASTING']-end_pos['UTM_EASTING'])**2+(start_pos['UTM_NORTHING']-end_pos['UTM_NORTHING'])**2)
    time_between_obs = (end_pos['sec_since']-start_pos['sec_since'])/3600/24

    velocity_obs_middle=pd.DataFrame({'date1':[vel_date_start],'date2':[vel_date_end],'start_easting':[start_pos['UTM_EASTING']],'start_northing':[start_pos['UTM_NORTHING']],'end_easting':[end_pos['UTM_EASTING']],'end_northing':[end_pos['UTM_NORTHING']],'distance_traveled (m)':[distance_traveled],'velocity (m/d)':[distance_traveled/time_between_obs]},index=[1])
    velocity_obs_middle['label']='middle'
    
    #Upper
    for i in range(0,len(upper_pos_data_daily)):
        if upper_pos_data_daily['date'][i]==vel_date_start:
            start_pos= upper_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]
        if upper_pos_data_daily['date'][i]==vel_date_end:
            end_pos=upper_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]

    distance_traveled = np.sqrt((start_pos['UTM_EASTING']-end_pos['UTM_EASTING'])**2+(start_pos['UTM_NORTHING']-end_pos['UTM_NORTHING'])**2)
    time_between_obs = (end_pos['sec_since']-start_pos['sec_since'])/3600/24

    velocity_obs_upper=pd.DataFrame({'date1':[vel_date_start],'date2':[vel_date_end],'start_easting':[start_pos['UTM_EASTING']],'start_northing':[start_pos['UTM_NORTHING']],'end_easting':[end_pos['UTM_EASTING']],'end_northing':[end_pos['UTM_NORTHING']],'distance_traveled (m)':[distance_traveled],'velocity (m/d)':[distance_traveled/time_between_obs]},index=[2])
    velocity_obs_upper['label']='upper'
    
    #Arm
    for i in range(0,len(arm_pos_data_daily)):
        if arm_pos_data_daily['date'][i]==vel_date_start:
            start_pos= arm_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]
        if arm_pos_data_daily['date'][i]==vel_date_end:
            end_pos=arm_pos_data_daily[['date','UTM_EASTING','UTM_NORTHING','HGT(m)','sec_since']].iloc[i]

    distance_traveled = np.sqrt((start_pos['UTM_EASTING']-end_pos['UTM_EASTING'])**2+(start_pos['UTM_NORTHING']-end_pos['UTM_NORTHING'])**2)
    time_between_obs = (end_pos['sec_since']-start_pos['sec_since'])/3600/24

    velocity_obs_arm=pd.DataFrame({'date1':[vel_date_start],'date2':[vel_date_end],'start_easting':[start_pos['UTM_EASTING']],'start_northing':[start_pos['UTM_NORTHING']],'end_easting':[end_pos['UTM_EASTING']],'end_northing':[end_pos['UTM_NORTHING']],'distance_traveled':[distance_traveled],'velocity':[distance_traveled/time_between_obs]},index=[3])
    velocity_obs_arm['label']='arm'
    
    # print(velocity_obs_lower)
    # print(velocity_obs_middle)
    # print(velocity_obs_upper)
    # print(velocity_obs_arm)

#     velocity_data_between_two_dates = pd.concat([velocity_obs_lower,velocity_obs_middle,velocity_obs_upper],axis=0)
    velocity_data_between_two_dates = pd.concat([velocity_obs_lower,velocity_obs_middle,velocity_obs_upper,velocity_obs_arm],axis=0)

    csv_name='Kaskawulsh_v2_'+vel_date_start.strftime('%Y-%m-%d')+'_to_'+vel_date_end.strftime('%Y-%m-%d')+'_GPS.csv'
    velocity_data_between_two_dates.to_csv(csv_name)