Fixes ISIS CSM serial issue and adds geoparquet support.

CHANGELOG.md

@@ -37,6 +37,7 @@ release.
 ### Added
 - Hard coded support for `csminit`ed cubes that have serial numbers that differ from `spiceinit`ed cubes.
 - Added support for scale and offset for GeoDataset objects.
+- The ability to write ISIS formated control networks to geoparquet. Geoparquet files are smaller, read/write faster, and interoperate in spatial tools like QGIS.
 
 ### Fixed
 - Fixed a bug where scale and offset were being applied backwards to GeoDataset objects.

plio/io/io_controlnetwork.py

@@ -2,15 +2,18 @@
 from time import gmtime, strftime
 import warnings
 
+import geopandas as gpd
 import pandas as pd
 import numpy as np
+import pyproj
 import pvl
 import struct
 
 from plio.io import ControlNetFileV0002_pb2 as cnf
 from plio.io import ControlNetFileHeaderV0005_pb2 as cnh5
 from plio.io import ControlPointFileEntryV0005_pb2 as cnp5
 from plio.utils.utils import xstr, find_in_dict
+from plio.utils.df_manipulation import add_latlon_to_network
 
 HEADERSTARTBYTE = 65536
 DEFAULTUSERNAME = 'None'
@@ -98,7 +101,6 @@ def to_protobuf(self, version=2):
     def from_protobuf(cls, protobuf):
         return cls(protobuf.doubleDataType, protobuf.doubleDataValue)
 
-
 class IsisControlNetwork(pd.DataFrame):
 
     # normal properties
@@ -500,3 +502,27 @@ def create_pvl_header(self, version, headerstartbyte,
         ])
 
         return pvl.dumps(header, encoder=encoder) + "\n"
+
+def to_parquet(obj, path, crs='IAU_2015:49900'):
+    def prepare_for_parquet(obj):
+        obj.loc[:, 'aprioriCovar'] = obj['aprioriCovar'].apply(lambda x: list(x))
+        obj.loc[:, 'adjustedCovar'] = obj['adjustedCovar'].apply(lambda x: list(x))
+        obj.loc[:, 'pointLog'] = obj['pointLog'].apply(lambda x: list(x))
+
+        if not ('lon' in obj.columns) or not ('lat' in obj.columns):
+            add_latlon_to_network(obj, crs)
+
+    prepare_for_parquet(obj)
+
+    network = gpd.GeoDataFrame(obj, geometry=gpd.points_from_xy(obj.lon, obj.lat))
+    network = network.set_crs(crs=crs, inplace=True)
+    network.to_parquet(path)
+
+def from_parquet(path):
+    return gpd.read_parquet(path)
+
+def isis2parquet(inpath, outpath, crs='IAU_2015:49900'):
+    to_parquet(from_isis(inpath), outpath, crs=crs)
+
+def parquet2isis(inpath, outpath, targetname, **kwargs):
+    to_isis(from_parquet(inpath), outpath, targetname, **kwargs)

plio/io/isis_serial_number.py

@@ -106,7 +106,7 @@ def generate_serial_number(label):
                 serial_entry = search_translation['*']
             if isinstance(serial_entry, datetime.datetime):
                 # PVL returns datetime objects now. Convert these to string and strip trailing zeros on microseconds.
-                serial_entry = serial_entry.strftime('%Y-%m-%dT%H:%M:%S.%f').rstrip('0')                
+                serial_entry = serial_entry.strftime('%Y-%m-%dT%H:%M:%SZ%f').rstrip('0')                
             serial_number.append(serial_entry)
         except:
             pass

plio/utils/df_manipulation.py

@@ -0,0 +1,31 @@
+import pyproj
+from plio.utils.proj import parse_radii_from_proj
+
+def add_latlon_to_network(network, crs):
+    """
+    Adds a 'lat' and 'lon' columns to a control network data frame. This
+    function modifies the dataframe in place.
+
+    Parameters
+    ----------
+    network : pd.DataFrame
+            A pandas dataframe with 'adjustedX', 'adjustedY', and 'adjustedZ' columns
+
+    semimajor_axis : float
+                    The semimajor axis of the ellipsoid for the body. Defaults to 3396190 (Mars).
+
+    semiminor_axis : float
+                    The semiminor axis of the ellipsoid for the body. Defaults to 3396190 (Mars).
+    """
+    semimajor_axis, semiminor_axis = parse_radii_from_proj(crs)
+    proj_str = f"""
+        +proj=pipeline
+        +step +proj=cart +a={semimajor_axis} +b={semiminor_axis} +inv
+    """
+    geocent2latlon = pyproj.transformer.Transformer.from_pipeline(proj_str)
+    lon, lat, _ = geocent2latlon.transform(network['adjustedX'],
+                                        network['adjustedY'],
+                                        network['adjustedZ'],
+                                        errcheck=True)
+    network.loc[:, 'lon'] = lon
+    network.loc[:, 'lat'] = lat

plio/utils/proj.py

@@ -0,0 +1,7 @@
+import pyproj
+
+def parse_radii_from_proj(crs):
+    prj = pyproj.CRS(crs)
+    semimajor_axis = prj.ellipsoid.semi_major_metre
+    semiminor_axis = prj.ellipsoid.semi_minor_metre
+    return semimajor_axis, semiminor_axis