Merge pull request DOI-USGS#130 from amystamile-usgs/sensor-util-func…

…tions

Adds Sensor Util Functions

examples/sensor_utils.ipynb

@@ -0,0 +1,312 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Sensor Utils\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import os\n",
+    "\n",
+    "from csmapi import csmapi\n",
+    "from knoten import csm, sensor_utils\n",
+    "\n",
+    "from knoten.shape import Ellipsoid\n",
+    "from knoten.illuminator import Illuminator\n",
+    "\n",
+    "import ale\n",
+    "import json"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Create a usgscsm sensor model"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "fileName = \"data/N1573082850_1.cub\"\n",
+    "\n",
+    "kernels = ale.util.generate_kernels_from_cube(fileName, expand=True)\n",
+    "isd_string = ale.loads(fileName, props={'kernels': kernels})\n",
+    "csm_isd = os.path.splitext(fileName)[0] + '.json'\n",
+    "\n",
+    "with open(csm_isd, 'w') as isd_file:\n",
+    "    isd_file.write(isd_string)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Run Sensor Utils with usgscsm sensor model and image point"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "camera = csm.create_csm(csm_isd)\n",
+    "image_pt = csmapi.ImageCoord(511.5, 511.5)\n",
+    "shape = Ellipsoid.from_csm_sensor(camera)\n",
+    "illuminator = Illuminator()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "38.87212509629895"
+      ]
+     },
+     "execution_count": 4,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "phaseAngle = sensor_utils.phase_angle(image_pt, camera, shape, illuminator)\n",
+    "\n",
+    "phaseAngle"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "49.60309924893989"
+      ]
+     },
+     "execution_count": 5,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "emissionAngle = sensor_utils.emission_angle(image_pt, camera, shape)\n",
+    "\n",
+    "emissionAngle"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "2903512972.146115"
+      ]
+     },
+     "execution_count": 6,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "slantDistance = sensor_utils.slant_distance(image_pt, camera, shape)\n",
+    "\n",
+    "slantDistance"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "2943536048.858226"
+      ]
+     },
+     "execution_count": 7,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "targetCenterDistance = sensor_utils.target_center_distance(image_pt, camera)\n",
+    "\n",
+    "targetCenterDistance"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "LatLon(lat=3.2229625890973583, lon=258.6197326526089)"
+      ]
+     },
+     "execution_count": 8,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "subSpacecraftPoint = sensor_utils.sub_spacecraft_point(image_pt, camera)\n",
+    "\n",
+    "subSpacecraftPoint"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 9,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "59096282.024265066"
+      ]
+     },
+     "execution_count": 9,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "localRadius = sensor_utils.local_radius(image_pt, camera, shape)\n",
+    "\n",
+    "localRadius"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 10,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "(79.34815579474038, -2.7790780986459485)"
+      ]
+     },
+     "execution_count": 10,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "rightAscDec = sensor_utils.right_ascension_declination(image_pt, camera)\n",
+    "\n",
+    "rightAscDec"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 11,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "17397.96094194587"
+      ]
+     },
+     "execution_count": 11,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "lineResolution = sensor_utils.line_resolution(image_pt, camera, shape)\n",
+    "\n",
+    "lineResolution"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 12,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "17397.93370038153"
+      ]
+     },
+     "execution_count": 12,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "sampleResolution = sensor_utils.sample_resolution(image_pt, camera, shape)\n",
+    "\n",
+    "sampleResolution"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 13,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "17397.9473211637"
+      ]
+     },
+     "execution_count": 13,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "pixelResolution = sensor_utils.pixel_resolution(image_pt, camera, shape)\n",
+    "\n",
+    "pixelResolution"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.9.18"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}

knoten/bundle.py

@@ -212,6 +212,35 @@ def compute_ground_partials(sensor, ground_pt):
     partials = np.array(sensor.computeGroundPartials(csm_ground))
     return np.reshape(partials, (2, 3))
 
+def compute_image_partials(sensor, ground_pt):
+    """
+    Compute the partial derivatives of the ground point with respect to
+    the line and sample at a ground point.
+
+    These are not normally available from the CSM model, so we use
+    csm::RasterGM::computeGroundPartials to get the Jacobian of the ground to
+    image transformation. Then we use the pseudoinverse of that to get the
+    Jacobian of the image to ground transformation.
+
+    Parameters
+    ----------
+    sensor : CSM sensor
+             The CSM sensor model
+    ground_pt : array
+                The (x, y, z) ground point to compute the partial derivatives W.R.T.
+
+    Returns
+    -------
+     : array
+       The partial derivatives of the image to ground transformation
+    """
+    if isinstance(ground_pt, csmapi.EcefCoord):
+        ground_pt = [ground_pt.x, ground_pt.y, ground_pt.z]
+    ground_matrix = compute_ground_partials(sensor, ground_pt)
+    image_matrix = np.linalg.pinv(ground_matrix)
+
+    return image_matrix.flatten()
+
 def compute_coefficient_columns(network, sensors, parameters):
     """
     Compute the columns for different coefficients

knoten/csm.py

@@ -73,6 +73,38 @@ def create_camera(label, url='http://pfeffer.wr.usgs.gov/api/1.0/pds/'):
     if plugin.canModelBeConstructedFromISD(isd, model_name):
         model = plugin.constructModelFromISD(isd, model_name)
         return model
+
+def get_state(sensor, image_pt):
+    """
+    Get the state of the sensor model at a given image point.
+
+    Parameters
+    ----------
+    sensor : object
+             A CSM compliant sensor model object
+
+    image_pt : tuple
+               Pair of x, y (sample, line) coordinates in pixel space
+
+    Returns
+    -------
+    : dict
+        Dictionary containing lookVec, sensorPos, sensorTime, and imagePoint
+    """
+    if not isinstance(sensor, csmapi.RasterGM):
+        raise TypeError("inputted sensor not a csm.RasterGM object")
+
+    sensor_time = sensor.getImageTime(image_pt)
+    locus = sensor.imageToRemoteImagingLocus(image_pt)
+    sensor_position = sensor.getSensorPosition(image_pt)
+
+    sensor_state = {
+        "lookVec": locus.direction,
+        "sensorPos": sensor_position,
+        "sensorTime": sensor_time,
+        "imagePoint": image_pt
+    }
+    return sensor_state
 
 def _from_state(state, verbose):
     with open(state, 'r') as stream: