added projection code

ale/base/base.py

@@ -1,6 +1,9 @@
 import pvl
 import json
 
+import tempfile
+import os 
+
 class Driver():
     """
     Base class for all Drivers.
@@ -323,3 +326,60 @@ def center_ephemeris_time(self):
     @property
     def short_mission_name(self):
         return self.__module__.split('.')[-1].split('_')[0]
+
+    @property 
+    def projection(self):
+        if not hasattr(self, "_projection"): 
+            try: 
+              from osgeo import gdal 
+            except: 
+                self._projection = None
+                return self._projection
+
+            if isinstance(self._file, pvl.PVLModule):
+                # save it to a temp folder
+                with tempfile.NamedTemporaryFile() as tmp:
+                    tmp.write(pvl.dumps(self._file)) 
+
+                    geodata = gdal.Open(tempfile.name)
+                    self._projection = geodata.GetSpatialRef()
+            else: 
+                # should be a path
+                if not os.path.exists(self._file): 
+                    self._projection = None 
+                else: 
+                    geodata = gdal.Open(self._file)
+                    self._projection = geodata.GetSpatialRef()
+
+            # is None if not projected
+            if self._projection: 
+                self._projection = self._projection.ExportToProj4()
+
+        return self._projection
+
+
+    @property 
+    def geotransform(self):
+        if not hasattr(self, "_geotransform"): 
+            try: 
+              from osgeo import gdal 
+            except: 
+                self._geotransform = (0.0, 1.0, 0.0, 0.0, 0.0, 1.0)
+                return self._geotransform
+
+            if isinstance(self._file, pvl.PVLModule):
+                # save it to a temp folder
+                with tempfile.NamedTemporaryFile() as tmp:
+                    tmp.write(pvl.dumps(self._file)) 
+
+                    geodata = gdal.Open(tempfile.name)
+                    self._geotransform = geodata.GetGeoTransform()
+            else: 
+                # should be a path
+                if not os.path.exists(self._file): 
+                    self._geotransform = (0.0, 1.0, 0.0, 0.0, 0.0, 1.0) 
+                else: 
+                    geodata = gdal.Open(self._file)
+                    self._geotransform = geodata.GetGeoTransform()
+
+        return self._geotransform

ale/formatters/usgscsm_formatter.py

@@ -51,6 +51,9 @@ def to_usgscsm(driver):
         'unit' : 'm'
     }
 
+    isd_data["projection"] = driver.projection
+    isd_data["geotransform"] = driver.geotransform
+
     # shared isd keywords for Framer and Linescanner
     if isinstance(driver, LineScanner) or isinstance(driver, Framer):
         # exterior orientation for just Framer and LineScanner

environment.yml

@@ -7,6 +7,7 @@ dependencies:
   - cmake>=3.15
   - pytest
   - eigen
+  - gdal 
   - jupyter
   - nlohmann_json
   - numpy

recipe/meta.yaml

@@ -31,7 +31,9 @@ requirements:
     - scipy >=1.4.0
     - spiceypy >=4.0.1
     - pyyaml
-
+  run_contrained: 
+    - gdal 
+
 test:
   imports:
     - ale

tests/pytests/data/B10_013341_1010_XN_79S172W/B10_013341_1010_XN_79S172W_isis3.lbl

@@ -69,6 +69,25 @@ Object = IsisCube
     Source                    = isis
   End_Group
 
+  Group = Mapping
+    ProjectionName     = Sinusoidal
+    CenterLongitude    = 148.36859083039
+    TargetName         = MARS
+    EquatorialRadius   = 3396190.0 <meters>
+    PolarRadius        = 3376200.0 <meters>
+    LatitudeType       = Planetocentric
+    LongitudeDirection = PositiveEast
+    LongitudeDomain    = 360
+    MinimumLatitude    = 63.636322793577
+    MaximumLatitude    = 87.296295823424
+    MinimumLongitude   = 139.6658284858
+    MaximumLongitude   = 157.07135317498
+    UpperLeftCornerX   = -219771.1526456 <meters>
+    UpperLeftCornerY   = 5175537.8728989 <meters>
+    PixelResolution    = 1455.4380969907 <meters/pixel>
+    Scale              = 40.726361118253 <pixels/degree>
+  End_Group
+
   Group = AlphaCube
     AlphaSamples        = 5000
     AlphaLines          = 24576

tests/pytests/test_load.py

@@ -73,12 +73,10 @@ def test_load_mes_from_metakernels(tmpdir, monkeypatch, mess_kernels):
         mk_file.write(mk_str)
 
     usgscsm_isd_obj = ale.load(label_file, verbose=True)
-
     assert usgscsm_isd_obj['name_platform'] == 'MESSENGER'
     assert usgscsm_isd_obj['name_sensor'] == 'MERCURY DUAL IMAGING SYSTEM NARROW ANGLE CAMERA'
     assert usgscsm_isd_obj['name_model'] == 'USGS_ASTRO_FRAME_SENSOR_MODEL'
 
-
 def test_load_mes_with_no_metakernels(tmpdir, monkeypatch, mess_kernels):
     monkeypatch.setenv('ALESPICEROOT', str(tmpdir))
 

tests/pytests/test_mex_drivers.py

@@ -202,7 +202,9 @@ def usgscsm_compare_dict():
               "t0_ephemeris": -101.83713859319687,
               "dt_ephemeris": 40.734855437278746,
               "t0_quaternion": -101.83713859319687,
-              "dt_quaternion": 40.734855437278746
+              "dt_quaternion": 40.734855437278746,
+              "projection" : None,
+              "geotransform" : (0.0, 1.0, 0.0, 0.0, 0.0, 1.0)
               },
 
         "isis" :

tests/pytests/test_usgscsm_formatter.py

@@ -9,6 +9,8 @@
 from ale.base.data_naif import NaifSpice
 from ale.rotation import ConstantRotation, TimeDependentRotation
 
+from conftest import get_image_label
+
 class TestDriver(Driver, NaifSpice):
     """
     Test Driver implementation with dummy values
@@ -342,3 +344,20 @@ def test_line_scan_sun_position(test_line_scan_driver):
     assert sun_position_obj['positions'] == [[0, 1, 2], [3, 4, 5]]
     assert sun_position_obj['velocities'] == [[0, -1, -2], [-3, -4, -5]]
     assert sun_position_obj['unit'] == 'm'
+
+def test_no_projection(test_frame_driver):
+    isd = usgscsm_formatter.to_usgscsm(test_frame_driver)
+    # isn't using real projection so it should be None
+    assert isd['projection'] == None
+
+def test_isis_projection():
+    isd = usgscsm_formatter.to_usgscsm(TestLineScanner(get_image_label('B10_013341_1010_XN_79S172W', "isis3")))
+    assert isd["projection"] == "+proj=sinu +lon_0=148.36859083039 +x_0=0 +y_0=0 +R=3396190 +units=m +no_defs"
+
+
+def test_isis_geotransform():
+    isd = usgscsm_formatter.to_usgscsm(TestLineScanner(get_image_label('B10_013341_1010_XN_79S172W', "isis3")))
+    expected = (-219771.1526456, 1455.4380969907, 0.0, 5175537.8728989, 0.0, -1455.4380969907)
+    for value, truth in zip(isd["geotransform"], expected):
+        pytest.approx(value, truth)
+