From da8bab611348e0afcb2d34dad9f0a1cf89025015 Mon Sep 17 00:00:00 2001
From: Kien Le <kienle@passcal.nmt.edu>
Date: Fri, 21 Apr 2023 13:48:19 -0600
Subject: [PATCH] Fix problem with GPS plot

---
 .../plotting/gps_plot/extract_gps_data.py     | 59 +++++++++++--------
 .../view/plotting/state_of_health_widget.py   |  4 ++
 2 files changed, 37 insertions(+), 26 deletions(-)

diff --git a/sohstationviewer/view/plotting/gps_plot/extract_gps_data.py b/sohstationviewer/view/plotting/gps_plot/extract_gps_data.py
index 18bc5fe15..35f04bcd3 100644
--- a/sohstationviewer/view/plotting/gps_plot/extract_gps_data.py
+++ b/sohstationviewer/view/plotting/gps_plot/extract_gps_data.py
@@ -195,7 +195,7 @@ def get_gps_channel_prefix(data_obj: MSeed, data_type: str) -> Optional[str]:
 
 
 def get_chan_soh_trace_as_dict(data_obj: MSeed, chan: str
-                               ) -> Dict[float, float]:
+                               ) -> Dict[float, List[float]]:
     """
     Get the data of a channel as a dictionary mapping a data point's time to
     its data. Suppose that gps data's spr <= 1.
@@ -204,13 +204,12 @@ def get_chan_soh_trace_as_dict(data_obj: MSeed, chan: str
     :return: a dict that maps the times of channel chan to its data
     """
     chan_data = data_obj.soh_data[data_obj.selected_key][chan]
-    data = chan_data['tracesInfo']['data']
-    times = chan_data['tracesInfo']['times']
-    # Remove the filled in gaps from the data
-    data = data[~data.mask]
-    times = times[~times.mask]
-    data_dict = {time: data
-                 for time, data in np.column_stack((times, data))}
+    traces = chan_data['tracesInfo']
+    times = np.hstack([trace['times'] for trace in traces])
+    data = np.hstack([trace['data'] for trace in traces])
+    data_dict = {}
+    for times_point, data_point in zip(times, data):
+        data_dict.setdefault(times_point, []).append(data_point)
     return data_dict
 
 
@@ -252,13 +251,15 @@ def extract_gps_data_pegasus_centaur(data_obj: MSeed, data_type: str
         # There is no channel for GPS fix type in Pegasus and Centaur data,
         # so we are giving it a dummy value.
         fix_type = 'N/A'
-        lat = la_dict.get(time, None)
-        long = lo_dict.get(time, None)
-        height = el_dict.get(time, None)
+        current_lat = la_dict.get(time, None)
+        current_long = lo_dict.get(time, None)
+        current_height = el_dict.get(time, None)
         # We are ignoring any point that does not have complete location data.
         # It might be possible to ignore points with missing latitude or
         # longitude, seeing as height is not required to plot a GPS point.
-        if lat is None or long is None or height is None:
+        if (current_lat is None or
+                current_long is None or
+                current_height is None):
             continue
         # Convert the location data to the appropriate unit. Centaur and
         # Pegasus dataloggers both store latitude and longitude in microdegree,
@@ -266,20 +267,26 @@ def extract_gps_data_pegasus_centaur(data_obj: MSeed, data_type: str
         # stores elevation in micrometer while Pegasus stores it in centimeter,
         # so we have to use a difference conversion factor for elevation
         # depending on the data type.
-        lat = lat / 1e6
-        long = long / 1e6
-        # The GPS prefix is unique between Pegasus (V) and Centaur (G), so we
-        # can use it in place of the data type.
-        if gps_prefix == 'V':
-            height_factor = 100
-        else:
-            height_factor = 1e6
-        height = height / height_factor
-        height_unit = 'M'
-        formatted_time = UTCDateTime(time).strftime('%Y-%m-%d %H:%M:%S')
-        gps_point = GPSPoint(formatted_time, fix_type, num_sats_used, lat,
-                             long, height, height_unit)
-        extracted_gps_points.append(gps_point)
+        for i, num_sats in enumerate(num_sats_used):
+            try:
+                lat = current_lat[i] / 1e6
+                long = current_long[i] / 1e6
+                # The GPS prefix is unique between Pegasus (V) and Centaur (G),
+                # so we can use it in place of the data type.
+                if gps_prefix == 'V':
+                    height_factor = 100
+                else:
+                    height_factor = 1e6
+                height = current_height[i] / height_factor
+                height_unit = 'M'
+                formatted_time = UTCDateTime(time).strftime(
+                    '%Y-%m-%d %H:%M:%S'
+                )
+                gps_point = GPSPoint(formatted_time, fix_type, num_sats, lat,
+                                     long, height, height_unit)
+                extracted_gps_points.append(gps_point)
+            except IndexError:
+                break
     # We only need to loop through one dictionary. If a time is not
     # available for a channel, the GPS data point at that time would be
     # invalid (it is missing a piece of data). Once we loop through a
diff --git a/sohstationviewer/view/plotting/state_of_health_widget.py b/sohstationviewer/view/plotting/state_of_health_widget.py
index 583721e37..f836ba1e4 100644
--- a/sohstationviewer/view/plotting/state_of_health_widget.py
+++ b/sohstationviewer/view/plotting/state_of_health_widget.py
@@ -138,6 +138,10 @@ class SOHWidget(plotting_widget.PlottingWidget):
         :param first_time: bool - flag identify when the data set is first
             plotted, not zoom in yet
         """
+        # These channels contain GPS data, which we only show in the GPS
+        # dialog.
+        if chan_id in {'GLO', 'VLO', 'VLA', 'GLA', 'VEL', 'VNS', 'GNS', 'GEL'}:
+            return
         chan_db_info = c_data['chan_db_info']
         plot_type = chan_db_info['plotType']
         if c_data['samplerate'] <= 1:
-- 
GitLab