[mock_uss, monitorlib] Add EGM96 datum conversion (interuss#342)

* Add EGM96 datum conversion

* make format

* Address comments and `make format`

monitoring/mock_uss/riddp/routes_observation.py

@@ -24,6 +24,7 @@
 from .config import KEY_RID_VERSION
 from .database import db
 from monitoring.monitorlib.formatting import limit_resolution
+from monitoring.monitorlib.geo import egm96_geoid_offset
 
 
 def _make_flight_observation(
@@ -61,6 +62,8 @@ def _make_flight_observation(
         paths.append(current_path)
 
     p = flight.most_recent_position
+    msl_alt = p.alt - egm96_geoid_offset(s2sphere.LatLng.from_degrees(p.lat, p.lng))
+    # TODO: Return msl_alt in observation information
     current_state = observation_api.CurrentState(
         timestamp=p.time.isoformat(),
         operational_status=flight.operational_status,

monitoring/monitorlib/assets/README.md

@@ -0,0 +1,9 @@
+# monitorlib assets
+
+This folder contains non-code artifacts used by monitorlib.
+
+## WW15MGH.DAC
+
+EGM96 geoid offsets from WGS84 ellipsoid in binary 15-minute grid format.
+
+Previously hosted at: http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/binary/WW15MGH.DAC

monitoring/monitorlib/geo.py

@@ -1,10 +1,13 @@
 from __future__ import annotations
 import math
 from enum import Enum
+import os
 from typing import List, Tuple, Union, Optional
 
 from implicitdict import ImplicitDict
+import numpy as np
 import s2sphere
+from scipy.interpolate import RectBivariateSpline as Spline
 import shapely.geometry
 from uas_standards.astm.f3548.v21 import api as f3548v21
 from uas_standards.astm.f3411.v19 import api as f3411v19
@@ -462,3 +465,39 @@ def to_vertices(self) -> List[s2sphere.LatLng]:
 
 def latitude_degrees(distance_meters: float) -> float:
     return 360 * distance_meters / EARTH_CIRCUMFERENCE_M
+
+
+_egm96: Optional[Spline] = None
+"""Cached EGM96 geoid interpolation function with inverted latitude"""
+
+
+def egm96_geoid_offset(p: s2sphere.LatLng) -> float:
+    """Estimate the EGM96 geoid height above the WGS84 ellipsoid.
+
+    Args:
+        p: Point where offset should be estimated.
+
+    Returns: Meters above WGS84 ellipsoid of the EGM96 geoid at p.
+    """
+    global _egm96
+    if _egm96 is None:
+        grid_size = 0.25  # degrees
+        # Latitude data is [90, -90] degrees
+        lats = np.arange(-90, 90 + grid_size / 2, grid_size)
+        # Longitude data is [0, 360) degrees
+        lngs = np.arange(0, 360, grid_size)
+        grid_path = os.path.join(os.path.dirname(__file__), "assets/WW15MGH.DAC")
+        grid = np.fromfile(grid_path, ">i2").reshape(lats.size, lngs.size) / 100
+        _egm96 = Spline(lats, lngs, grid)
+    lng = math.fmod(p.lng().degrees, 360)
+    while lng < 0:
+        lng += 360
+    lat = p.lat().degrees
+    if lat < -90 or lat > 90:
+        raise ValueError(f"Cannot compute EGM96 geoid offset at latitude {lat} degrees")
+
+    # Negative latitude because the grid file lists offsets from 90 to -90
+    # degrees latitude, but Splines must have increasing X so latitudes must be
+    # listed -90 to 90.  Since latitude data are symmetric, we can simply
+    # convert "-90 to 90" to "90 to -90" by inverting the requested latitude.
+    return _egm96.ev(-lat, lng)

requirements.txt

@@ -25,6 +25,7 @@ locust==1.2.2  # loadtest
 loguru==0.6.0
 lxml==4.9.1
 marko==1.2.2  # uss_qualifier
+numpy==1.24.4
 pem==21.2.0  # deployment_manager
 pvlib==0.10.1
 pyjwt==2.4.0
@@ -36,6 +37,7 @@ pytest-mock==3.6.1
 pyyaml==6.0.1
 requests==2.31.0
 s2sphere==0.2.5
+scipy==1.10.1
 shapely==1.7.1
 structlog==21.5.0  # deployment_manager
 termcolor==1.1.0

schemas/monitoring/monitorlib/clients/flight_planning/flight_info_template/BasicFlightPlanInformationTemplate.json

@@ -7,6 +7,13 @@
       "description": "Path to content that replaces the $ref",
       "type": "string"
     },
+    "area": {
+      "description": "User intends to or may fly anywhere in this entire area.",
+      "items": {
+        "$ref": "../../../geotemporal/Volume4DTemplate.json"
+      },
+      "type": "array"
+    },
     "uas_state": {
       "description": "State of the user's UAS associated with this flight plan.",
       "enum": [
@@ -26,6 +33,7 @@
     }
   },
   "required": [
+    "area",
     "uas_state",
     "usage_state"
   ],

schemas/monitoring/monitorlib/geo/Altitude.json

@@ -0,0 +1,34 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/geo/Altitude.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.geo.Altitude, as defined in monitoring/monitorlib/geo.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "reference": {
+      "enum": [
+        "W84",
+        "SFC"
+      ],
+      "type": "string"
+    },
+    "units": {
+      "enum": [
+        "M",
+        "FT"
+      ],
+      "type": "string"
+    },
+    "value": {
+      "type": "number"
+    }
+  },
+  "required": [
+    "reference",
+    "units",
+    "value"
+  ],
+  "type": "object"
+}

schemas/monitoring/monitorlib/geo/Circle.json

@@ -0,0 +1,22 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/geo/Circle.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.geo.Circle, as defined in monitoring/monitorlib/geo.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "center": {
+      "$ref": "LatLngPoint.json"
+    },
+    "radius": {
+      "$ref": "Radius.json"
+    }
+  },
+  "required": [
+    "center",
+    "radius"
+  ],
+  "type": "object"
+}

schemas/monitoring/monitorlib/geo/Polygon.json

@@ -0,0 +1,21 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/geo/Polygon.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.geo.Polygon, as defined in monitoring/monitorlib/geo.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "vertices": {
+      "items": {
+        "$ref": "LatLngPoint.json"
+      },
+      "type": "array"
+    }
+  },
+  "required": [
+    "vertices"
+  ],
+  "type": "object"
+}

schemas/monitoring/monitorlib/geo/Radius.json

@@ -0,0 +1,26 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/geo/Radius.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.geo.Radius, as defined in monitoring/monitorlib/geo.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "units": {
+      "enum": [
+        "M",
+        "FT"
+      ],
+      "type": "string"
+    },
+    "value": {
+      "type": "number"
+    }
+  },
+  "required": [
+    "units",
+    "value"
+  ],
+  "type": "object"
+}

schemas/monitoring/monitorlib/geotemporal/Volume4DTemplate.json

@@ -0,0 +1,86 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/geotemporal/Volume4DTemplate.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.geotemporal.Volume4DTemplate, as defined in monitoring/monitorlib/geotemporal.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "altitude_lower": {
+      "description": "The minimum altitude at which the virtual user will fly while using this volume for their flight.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../geo/Altitude.json"
+        }
+      ]
+    },
+    "altitude_upper": {
+      "description": "The maximum altitude at which the virtual user will fly while using this volume for their flight.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../geo/Altitude.json"
+        }
+      ]
+    },
+    "duration": {
+      "description": "If only one of start_time and end_time is specified, then the other time should be separated from the specified time by this amount.  May not be defined in both start_time and end_time are defined.",
+      "format": "duration",
+      "type": [
+        "string",
+        "null"
+      ]
+    },
+    "end_time": {
+      "description": "The time at which the virtual user will be finished using the specified geospatial area for their flight.  May not be defined if duration and start_time are defined.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../temporal/TestTime.json"
+        }
+      ]
+    },
+    "outline_circle": {
+      "description": "Circular outline/footprint of the specified area.  May not be defined if outline_polygon is defined.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../geo/Circle.json"
+        }
+      ]
+    },
+    "outline_polygon": {
+      "description": "Polygonal 2D outline/footprint of the specified area.  May not be defined if outline_circle is defined.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../geo/Polygon.json"
+        }
+      ]
+    },
+    "start_time": {
+      "description": "The time at which the virtual user may start using the specified geospatial area for their flight.  May not be defined if duration and end_time are defined.",
+      "oneOf": [
+        {
+          "type": "null"
+        },
+        {
+          "$ref": "../temporal/TestTime.json"
+        }
+      ]
+    }
+  },
+  "type": "object"
+}

schemas/monitoring/monitorlib/temporal/NextDay.json

@@ -0,0 +1,43 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/temporal/NextDay.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.temporal.NextDay, as defined in monitoring/monitorlib/temporal.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "days_of_the_week": {
+      "description": "Acceptable days of the week.  Omit to indicate that any day of the week is acceptable.",
+      "items": {
+        "enum": [
+          "Mo",
+          "Tu",
+          "We",
+          "Th",
+          "Fr",
+          "Sa",
+          "Su"
+        ],
+        "type": "string"
+      },
+      "type": [
+        "array",
+        "null"
+      ]
+    },
+    "starting_from": {
+      "$ref": "TestTime.json",
+      "description": "The time after which the first instance of one of the days should be found."
+    },
+    "time_zone": {
+      "description": "Time zone in which \"day\" is understood.  Examples:\n  * \"local\" (local time of machine running this code)\n  * \"Z\" (Zulu time)\n  * \"-08:00\" (ISO time zone)\n  * \"US/Pacific\" (IANA time zone)",
+      "type": "string"
+    }
+  },
+  "required": [
+    "starting_from",
+    "time_zone"
+  ],
+  "type": "object"
+}

schemas/monitoring/monitorlib/temporal/NextSunPosition.json

@@ -0,0 +1,29 @@
+{
+  "$id": "https://github.com/interuss/monitoring/blob/main/schemas/monitoring/monitorlib/temporal/NextSunPosition.json",
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "description": "monitoring.monitorlib.temporal.NextSunPosition, as defined in monitoring/monitorlib/temporal.py",
+  "properties": {
+    "$ref": {
+      "description": "Path to content that replaces the $ref",
+      "type": "string"
+    },
+    "elevation_deg": {
+      "description": "Elevation of the center of the sun above horizontal, in degrees.",
+      "type": "number"
+    },
+    "observed_from": {
+      "$ref": "../geo/LatLngPoint.json",
+      "description": "The location on earth observing the sun."
+    },
+    "starting_from": {
+      "$ref": "TestTime.json",
+      "description": "The time after which the first time the sun is at the specified position should be found."
+    }
+  },
+  "required": [
+    "elevation_deg",
+    "observed_from",
+    "starting_from"
+  ],
+  "type": "object"
+}