Fix issues with times near 00:00Z. Add log parser.

chasemapper/geometry.py

@@ -30,7 +30,7 @@ class GenericTrack(object):
         self.ASCENT_AVERAGING = ascent_averaging
         # Payload state.
         self.landing_rate = landing_rate
-        self.ascent_rate = 5.0
+        self.ascent_rate = 0.0
         self.heading = 0.0
         self.speed = 0.0
         self.is_descending = False
@@ -87,7 +87,7 @@ class GenericTrack(object):
     def calculate_ascent_rate(self):
         ''' Calculate the ascent/descent rate of the payload based on the available data '''
         if len(self.track_history) <= 1:
-            return 5.0
+            return 0.0
         elif len(self.track_history) == 2:
             # Basic ascent rate case - only 2 samples.
             _time_delta = (self.track_history[-1][0] - self.track_history[-2][0]).total_seconds()

chasemapper/listeners.py

@@ -12,10 +12,54 @@
 import socket, json, sys, traceback
 from threading import Thread
 from dateutil.parser import parse
-from datetime import datetime
+from datetime import datetime, timedelta
 
 MAX_JSON_LEN = 2048
 
+
+def fix_datetime(datetime_str, local_dt_str = None):
+    '''
+    Given a HH:MM:SS string from an telemetry sentence, produce a complete timestamp, using the current system time as a guide for the date.
+    '''
+
+    if local_dt_str is None:
+        _now = datetime.utcnow()
+    else:
+        _now = parse(local_dt_str)
+
+    # Are we in the rollover window?
+    if _now.hour == 23 or _now.hour == 0:
+        _outside_window = False
+    else:
+        _outside_window = True
+    
+    # Append on a timezone indicator if the time doesn't have one.
+    if datetime_str.endswith('Z') or datetime_str.endswith('+00:00'):
+        pass
+    else:
+        datetime_str += "Z"
+
+
+    # Parsing just a HH:MM:SS will return a datetime object with the year, month and day replaced by values in the 'default'
+    # argument. 
+    _telem_dt = parse(datetime_str, default=_now)
+
+    if _outside_window:
+        # We are outside the day-rollover window, and can safely use the current zulu date.
+        return _telem_dt
+    else:
+        # We are within the window, and need to adjust the day backwards or forwards based on the sonde time.
+        if _telem_dt.hour == 23 and _now.hour == 0:
+            # Assume system clock running slightly fast, and subtract a day from the telemetry date.
+            _telem_dt = _telem_dt - timedelta(days=1)
+
+        elif _telem_dt.hour == 00 and _now.hour == 23:
+            # System clock running slow. Add a day.
+            _telem_dt = _telem_dt + timedelta(days=1)
+
+        return _telem_dt
+
+
 class UDPListener(object):
     ''' UDP Broadcast Packet Listener 
     Listens for Horuslib UDP broadcast packets, and passes them onto a callback function
@@ -191,6 +235,8 @@ class OziListener(object):
         _full_time = datetime.utcnow().strftime("%Y-%m-%dT") + _short_time + "Z"
         _time_dt = parse(_full_time)
 
+        _time_dt = fix_datetime(_short_time)
+
         _output = {
             'time'  : _time_dt,
             'lat'   : _lat,

chasemapper/logger.py

@@ -9,6 +9,7 @@ import datetime
 import json
 import logging
 import os
+import pytz
 import time
 from threading import Thread
 try:
@@ -59,7 +60,7 @@ class ChaseLogger(object):
         """
 
         data['log_type'] = 'CAR POSITION'
-        data['log_time'] = datetime.datetime.utcnow().isoformat()
+        data['log_time'] = pytz.utc.localize(datetime.datetime.utcnow()).isoformat()
 
         # Convert the input datetime object into a string.
         data['time'] = data['time'].isoformat()
@@ -75,7 +76,7 @@ class ChaseLogger(object):
         """
 
         data['log_type'] = 'BALLOON TELEMETRY'
-        data['log_time'] = datetime.datetime.utcnow().isoformat()
+        data['log_time'] = pytz.utc.localize(datetime.datetime.utcnow()).isoformat()
 
         # Convert the input datetime object into a string.
         data['time'] = data['time_dt'].isoformat()
@@ -93,7 +94,7 @@ class ChaseLogger(object):
         """ Log a prediction run """
 
         data['log_type'] = 'PREDICTION'
-        data['log_time'] = datetime.datetime.utcnow().isoformat()
+        data['log_time'] = pytz.utc.localize(datetime.datetime.utcnow()).isoformat()
 
 
         # Add it to the queue if we are running.

horusmapper.py

@@ -10,6 +10,7 @@ import logging
 import flask
 from flask_socketio import SocketIO
 import os.path
+import pytz
 import sys
 import time
 import traceback
@@ -23,7 +24,7 @@ from chasemapper.geometry import *
 from chasemapper.gps import SerialGPS
 from chasemapper.gpsd import GPSDAdaptor
 from chasemapper.atmosphere import time_to_landing
-from chasemapper.listeners import OziListener, UDPListener
+from chasemapper.listeners import OziListener, UDPListener, fix_datetime
 from chasemapper.predictor import predictor_spawn_download, model_download_running
 from chasemapper.habitat import HabitatChaseUploader
 from chasemapper.logger import ChaseLogger
@@ -550,11 +551,13 @@ def udp_listener_summary_callback(data):
 
     # Process the 'short time' value if we have been provided it.
     if 'time' in data.keys():
-        _full_time = datetime.utcnow().strftime("%Y-%m-%dT") + data['time'] + "Z"
+        output['time_dt'] = fix_datetime(data['time'])
-        output['time_dt'] = parse(_full_time)
+        #_full_time = datetime.utcnow().strftime("%Y-%m-%dT") + data['time'] + "Z"
+        #output['time_dt'] = parse(_full_time)
     else:
         # Otherwise use the current UTC time.
-        output['time_dt'] = datetime.utcnow()
+        
+        output['time_dt'] = pytz.utc.localize(datetime.utcnow())
 
     # Copy out any extra fields that we want to pass on to the GUI.
     for _field in EXTRA_FIELDS:
@@ -576,7 +579,7 @@ def udp_listener_car_callback(data):
     _lon = data['longitude']
     _alt = data['altitude']
     _comment = "CAR"
-    _time_dt = datetime.utcnow()
+    _time_dt = pytz.utc.localize(datetime.utcnow())#datetime.utcnow()
 
     logging.debug("Car Position: %.5f, %.5f" % (_lat, _lon))
 

log_parse.py

@@ -0,0 +1,283 @@
+#!/usr/bin/env python2.7
+#
+#   Project Horus - Browser-Based Chase Mapper
+#   Log File Parsing
+#
+#   Copyright (C) 2018  Mark Jessop <vk5qi@rfhead.net>
+#   Released under GNU GPL v3 or later
+#
+import argparse
+import datetime
+import json
+import logging
+import sys
+import numpy as np
+import matplotlib.pyplot as plt
+from chasemapper.earthmaths import *
+from chasemapper.geometry import *
+from dateutil.parser import parse
+
+
+def read_file(filename):
+    """ Read log file, and output an array of dicts. """
+    _output = []
+
+    _f = open(filename, 'r')
+    for _line in _f:
+        try:
+            _data = json.loads(_line)
+            _output.append(_data)
+        except Exception as e:
+            logging.debug("Error reading line: %s" % str(e))
+
+    logging.info("Read %d log entries." % len(_output))
+
+    return _output
+
+
+
+def extract_data(log_entries):
+    """ Step through the log entries, and extract:
+    - Car position telemetry
+    - Balloon positions
+    - Predictions 
+    """
+
+    # There's only ever going to be one car showing up on the map, so we just output a list.
+    _car = []
+    # We might have more than one balloon though, so we use a dictionary, with one entry per callsign.
+    _telemetry = {}
+
+    for _entry in log_entries:
+
+        if _entry['log_type'] == "CAR POSITION":
+            _car.append(_entry)
+
+        elif _entry['log_type'] == "BALLOON TELEMETRY":
+            # Extract the callsign.
+            _call = _entry['callsign']
+
+            if _call not in _telemetry:
+                _telemetry[_call] = {'telemetry': [], 'predictions': []}
+
+            _telemetry[_call]['telemetry'].append(_entry)
+
+        elif _entry['log_type'] == "PREDICTION":
+            # Extract the callsign.
+            _call = _entry['callsign']
+
+            if _call not in _telemetry:
+                _telemetry[_call] = {'telemetry': [], 'predictions': []}
+
+            _telemetry[_call]['predictions'].append(_entry)
+
+    logging.info("Extracted %d Car Positions" % len(_car))
+    for _call in _telemetry:
+        logging.info("Callsign %s: Extracted %d telemetry positions, %d predictions." % (_call, len(_telemetry[_call]['telemetry']), len(_telemetry[_call]['predictions'])))
+
+    return (_car, _telemetry)
+
+
+
+def flight_stats(telemetry, ascent_threshold = 3.0, descent_threshold=-5.0,  landing_threshold = 0.5):
+    """ Process a set of balloon telemetry, and calculate some statistics about the flight """
+
+    asc_rate_avg_length = 5
+
+    _stats = {
+        'ascent_rates': np.array([]),
+        'positions': []
+    }
+
+    _flight_segment = "UNKNOWN"
+
+    _track = GenericTrack() 
+
+    for _entry in telemetry:
+        # Produce a dict which we can pass into the GenericTrack object.
+        _position = {
+            'time': parse(_entry['time']),
+            'lat': _entry['lat'],
+            'lon': _entry['lon'],
+            'alt': _entry['alt']
+        }
+
+        _stats['positions'].append([_position['time'], _position['lat'], _position['lon'], _position['alt']])
+        _state = _track.add_telemetry(_position)
+
+        if len(_stats['ascent_rates']) < asc_rate_avg_length:
+            # Not enough data to make any judgements about the state of the flight yet.
+            _stats['ascent_rates'] = np.append(_stats['ascent_rates'], _state['ascent_rate'])
+        else:
+            # Roll the array, and add the new value on the end.
+            _stats['ascent_rates'] = np.roll(_stats['ascent_rates'], -1)
+            _stats['ascent_rates'][-1] = _state['ascent_rate']
+
+            _mean_asc_rate = np.mean(_stats['ascent_rates'])
+
+            if _flight_segment == "UNKNOWN":
+                # Make a determination on flight state based on what we know now.
+                if _mean_asc_rate > ascent_threshold:
+                    _flight_segment = "ASCENT"
+                    _stats['launch'] = [_state['time'], _state['lat'], _state['lon'], _state['alt']]
+                    logging.info("Detected Launch: %s, %.5f, %.5f, %dm" % 
+                        (_state['time'].isoformat(), _state['lat'], _state['lon'], _state['alt']))
+                elif _mean_asc_rate < descent_threshold:
+                    _flight_segment = "DESCENT"
+                else:
+                    pass
+
+            if _flight_segment == "ASCENT":
+                if _track.track_history[-1][3] < _track.track_history[-2][3]:
+                    # Possible detection of burst.
+                    if 'burst_position' not in _stats:
+                        _stats['burst_position'] = _track.track_history[-2]
+                        logging.info("Detected Burst: %s, %.5f, %.5f, %dm" % (
+                            _stats['burst_position'][0].isoformat(),
+                            _stats['burst_position'][1],
+                            _stats['burst_position'][2],
+                            _stats['burst_position'][3]
+                            ))
+
+                if _mean_asc_rate < descent_threshold:
+                    _flight_segment = "DESCENT"
+                    continue
+
+            if _flight_segment == "DESCENT":
+                if abs(_mean_asc_rate) < landing_threshold:
+                    _stats['landing'] = [parse(_entry['log_time']), _state['lat'], _state['lon'], _state['alt']]
+                    logging.info("Detected Landing: %s, %.5f, %.5f, %dm" % 
+                        (_entry['log_time'], _state['lat'], _state['lon'], _state['alt']))
+
+
+                    _flight_segment = "LANDED"
+                    return _stats
+
+    return _stats
+
+
+def calculate_predictor_error(predictions, landing_time, lat, lon, alt):
+    """ Process a list of predictions, and determine the landing position error for each one """
+
+
+    _output = []
+
+    _landing = (lat, lon, alt)
+
+
+    for _predict in predictions:
+        _predict_time = _predict['log_time']
+
+        if parse(_predict_time) > (landing_time-datetime.timedelta(0,30)):
+            break
+
+        _predict_altitude = _predict['pred_path'][0][2]
+        _predict_landing = (_predict['pred_landing'][0], _predict['pred_landing'][1], _predict['pred_landing'][2])
+
+
+        _pos_info = position_info(_landing, _predict_landing)
+
+        logging.info("Prediction %s: Altitude %d, Predicted Landing: %.4f, %.4f Prediction Error: %.1f km, %s" % (
+            _predict_time,
+            int(_predict_altitude),
+            _predict['pred_landing'][0],
+            _predict['pred_landing'][1],
+            (_pos_info['great_circle_distance']/1000.0),
+            bearing_to_cardinal(_pos_info['bearing'])
+            ))
+
+        _output.append([
+            parse(_predict_time+"+00:00"),
+            _pos_info['great_circle_distance']/1000.0,
+            _pos_info['bearing'],
+            ])
+
+    return _output
+
+
+
+def plot_predictor_error(flight_stats, predictor_errors, callsign = ""):
+
+    # Get launch time.
+    _launch_time = flight_stats['launch'][0]
+
+    # Generate datasets of time-since-launch and altitude.
+    _flight_time = []
+    _flight_alt = []
+    for _entry in flight_stats['positions']:
+        _ft = (_entry[0]-_launch_time).total_seconds()/60.0
+        if _ft > 0:
+            _flight_time.append(_ft)
+            _flight_alt.append(_entry[3])
+
+    # Generate datasets of time-since-launch and altitude.
+    _predict_time = []
+    _predict_error = []
+    for _entry in predictor_errors:
+        _ft = (_entry[0]-_launch_time).total_seconds()/60.0
+        if _ft > 0:
+            _predict_time.append(_ft)
+            _predict_error.append(_entry[1])
+
+
+
+    # Altitude vs Time
+    plt.figure()
+    plt.plot(_flight_time, _flight_alt)
+    plt.grid()
+    plt.xlabel("Time (minutes)")
+    plt.ylabel("Altitude (metres)")
+    plt.title("Flight Profile - %s" % callsign)
+
+    # Prediction error vs time.
+    plt.figure()
+    plt.plot(_predict_time, _predict_error)
+    plt.xlabel("Time (minutes)")
+    plt.ylabel("Landing Prediction Error (km)")
+    plt.title("Landing Prediction Error - %s" % callsign)
+    plt.grid()
+    plt.show()
+
+
+
+if __name__ == "__main__":
+    parser = argparse.ArgumentParser()
+    parser.add_argument("filename", type=str, help="Input log file.")
+    parser.add_argument("-c", "--config", type=str, default="horusmapper.cfg", help="Configuration file.")
+    parser.add_argument("-v", "--verbose", action="store_true", default=False, help="Verbose output.")
+    parser.add_argument("--predict-error", action="store_true", default=False, help="Calculate Prediction Error.")
+    args = parser.parse_args()
+
+    # Configure logging
+    if args.verbose:
+        _log_level = logging.DEBUG
+    else:
+        _log_level = logging.INFO
+
+    logging.basicConfig(format='%(asctime)s %(levelname)s:%(message)s', stream=sys.stdout, level=_log_level)
+
+
+    _log_entries = read_file(args.filename)
+
+    _car, _telemetry = extract_data(_log_entries)
+
+    for _call in _telemetry:
+        logging.info("Processing Callsign: %s" % _call)
+        _stats = flight_stats(_telemetry[_call]['telemetry'])
+
+        if ('landing' in _stats) and ('launch' in _stats):
+            _total_flight = position_info((_stats['launch'][1],_stats['launch'][2],_stats['launch'][3]),(_stats['landing'][1], _stats['landing'][2], _stats['landing'][3]))
+            logging.info("%s Flight Distance: %.2f km" % (_call, _total_flight['great_circle_distance']/1000.0))
+
+        if args.predict_error:
+            if 'landing' in _stats:
+                _time = _stats['landing'][0]
+                _lat = _stats['landing'][1]
+                _lon = _stats['landing'][2]
+                _alt = _stats['landing'][3]
+                _predict_errors = calculate_predictor_error(_telemetry[_call]['predictions'], _time, _lat, _lon, _alt)
+
+                plot_predictor_error(_stats, _predict_errors, _call)
+
+
+