Merge pull request #490 from dbekaert/dev

Release 0.4.2

Author: jlmaurer

.github/workflows/deploy-docs.yml

@@ -14,7 +14,7 @@ jobs:
         with:
           fetch-depth: 0
 
-      - uses: mamba-org/provision-with-micromamba@v14
+      - uses: mamba-org/provision-with-micromamba@v15
         with:
           extra-specs: |
             python=3.10

.gitignore

@@ -20,6 +20,9 @@ test/data
 *.cpp
 test/*/geom/*.dem
 test/*/geom/*.nc
+test/test_geom/*.nc
+test/scenario_2/GMAO_Delay_*.csv
+test/GUNW/
 .eggs/
 dist/
 tools/RAiDER.egg-info/

CHANGELOG.md

@@ -6,6 +6,16 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [PEP 440](https://www.python.org/dev/peps/pep-0440/)
 and uses [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
++ Add assert statement to raise error if the delay cube for each SAR date in a GUNW IFG is not written 
+
+## [0.4.2]
+
+### New/Updated Features
++ `calcDelaysGUNW` allows processing with any supported weather model as listed in [`RAiDER.models.allowed.ALLOWED_MODELS`](https://github.com/dbekaert/RAiDER/blob/dev/tools/RAiDER/models/allowed.py).
++ Removed NCMR removed from supported model list till re-tested 
++ `credentials` looks for weather model API credentials RC_file hidden file, and creates it if it does not exists
++ Isolate ISCE3 imports to only those functions that need it.
++ Small bugfixes and updates to docstrings
 
 ## [0.4.1]
 

test/conftest.py

@@ -0,0 +1,19 @@
+import pytest
+
+
+def pytest_addoption(parser):
+    parser.addoption(
+        "--skip-isce3", action="store_true", default=False, help="skip tests which require ISCE3"
+    )
+
+
+def pytest_configure(config):
+    config.addinivalue_line("markers", "isce3: mark test as requiring ISCE3 to run")
+
+
+def pytest_collection_modifyitems(config, items):
+    if config.getoption("--skip-isce3"):
+        skip_isce3 = pytest.mark.skip(reason="--skip-isce3 option given")
+        for item in items:
+            if "isce3" in item.keywords:
+                item.add_marker(skip_isce3)

test/test_GUNW.py

@@ -6,11 +6,13 @@
 import numpy as np
 import rasterio as rio
 import xarray as xr
+import pytest
 
 from test import TEST_DIR
 
 WM = 'GMAO'
 
+@pytest.mark.isce3
 def test_GUNW():
     ## eventually to be implemented
     # home = os.path.expanduser('~')

test/test_credentials.py

@@ -0,0 +1,79 @@
+import os
+from platform import system
+from RAiDER.models import credentials
+
+# Test checking/creating ECMWF_RC CAPI file
+def test_ecmwfApi_createFile():
+    import ecmwfapi
+    
+    #Check extension for hidden files
+    hidden_ext = '_' if system()=="Windows" else '.'
+
+    # Test creation of ~/.ecmwfapirc file
+    ecmwf_file = os.path.expanduser('./') + hidden_ext + credentials.API_FILENAME['ERAI']
+    credentials.check_api('ERAI', 'dummy', 'dummy', './', update_flag=True)
+    assert os.path.exists(ecmwf_file) == True,f'{ecmwf_file} does not exist' 
+
+    # Get existing ECMWF_API_RC env if exist
+    default_ecmwf_file = os.getenv("ECMWF_API_RC_FILE")
+    if default_ecmwf_file is None:
+        default_ecmwf_file = ecmwfapi.api.DEFAULT_RCFILE_PATH
+
+    #Set it to current dir to avoid overwriting ~/.ecmwfapirc file 
+    os.environ["ECMWF_API_RC_FILE"] = ecmwf_file
+    key, url, uid = ecmwfapi.api.get_apikey_values()
+    
+    # Return to default_ecmwf_file and remove local API file
+    os.environ["ECMWF_API_RC_FILE"] = default_ecmwf_file 
+    os.remove(ecmwf_file)
+
+    #Check if API is written correctly
+    assert uid ==  'dummy', f'{ecmwf_file}: UID was not written correctly'
+    assert key == 'dummy', f'{ecmwf_file}: KEY was not written correctly'
+
+
+# Test checking/creating Copernicus Climate Data Store
+#   CDS_RC CAPI file 
+def test_cdsApi_createFile():
+    import cdsapi
+
+    #Check extension for hidden files
+    hidden_ext = '_' if system()=="Windows" else '.'
+    
+    # Test creation of .cdsapirc file in current dir
+    cds_file = os.path.expanduser('./') + hidden_ext + credentials.API_FILENAME['ERA5']
+    credentials.check_api('ERA5', 'dummy', 'dummy', './', update_flag=True)
+    assert os.path.exists(cds_file) == True,f'{cds_file} does not exist' 
+
+    # Check the content
+    cds_credentials = cdsapi.api.read_config(cds_file)
+    uid, key = cds_credentials['key'].split(':')
+
+    # Remove local API file
+    os.remove(cds_file)
+    
+    assert uid ==  'dummy', f'{cds_file}: UID was not written correctly'
+    assert key == 'dummy', f'{cds_file}: KEY was not written correctly'
+
+# Test checking/creating EARTHDATA_RC API file 
+def test_netrcApi_createFile():
+    import netrc
+    
+    #Check extension for hidden files
+    hidden_ext = '_' if system()=="Windows" else '.'
+
+    # Test creation of ~/.cdsapirc file
+    netrc_file = os.path.expanduser('./') + hidden_ext + credentials.API_FILENAME['GMAO']
+    credentials.check_api('GMAO', 'dummy', 'dummy', './', update_flag=True)
+    assert os.path.exists(netrc_file) == True,f'{netrc_file} does not exist' 
+
+    # Check the content
+    host = 'urs.earthdata.nasa.gov'
+    netrc_credentials = netrc.netrc(netrc_file)
+    uid, _, key = netrc_credentials.authenticators(host)
+
+    # Remove local API file
+    os.remove(netrc_file)
+    
+    assert uid ==  'dummy', f'{netrc_file}: UID was not written correctly'
+    assert key == 'dummy', f'{netrc_file}: KEY was not written correctly'

test/test_raytracing.py

@@ -1,49 +1,180 @@
 import datetime
 import pytest
+
 import numpy as np
 
-from RAiDER.losreader import Raytracing
+from pyproj import CRS
+from scipy.interpolate import RegularGridInterpolator as rgi
+
+from RAiDER.delay import _build_cube_ray
+from RAiDER.models.weatherModel import (
+    WeatherModel,
+)
+
+_LON0 = 0
+_LAT0 = 0
+_OMEGA = 0.1 / (180/np.pi)
+
+class MockWeatherModel(WeatherModel):
+    """Implement abstract methods for testing."""
+
+    def __init__(self):
+        super().__init__()
+
+        self._k1 = 1
+        self._k2 = 1
+        self._k3 = 1
+
+        self._Name = "MOCK"
+        self._valid_range = (datetime.datetime(1970, 1, 1), "Present")
+        self._lag_time = datetime.timedelta(days=15)
+
+    def _fetch(self, ll_bounds, time, out):
+        pass
+
+    def load_weather(self, *args, **kwargs):
+        _N_Z = 32
+        self._ys = np.arange(-2,3) + _LAT0
+        self._xs = np.arange(-3,4) + _LON0
+        self._zs = np.linspace(0, 1e5, _N_Z)
+        self._t = np.ones((len(self._ys), len(self._xs), _N_Z))
+        self._e = self._t.copy()
+        self._e[:,3:,:] = 2
+
+        _p = np.arange(31, -1, -1)
+        self._p = np.broadcast_to(_p, self._t.shape)
+
+        self._true_hydro_refr = np.broadcast_to(_p, (self._t.shape))
+        self._true_wet_ztd = 1e-6 * 2 * np.broadcast_to(np.flip(self._zs), (self._t.shape))
+        self._true_wet_ztd[:,3:] = 2 * self._true_wet_ztd[:,3:]
+
+        self._true_hydro_ztd = np.zeros(self._t.shape)
+        for layer in range(len(self._zs)):
+            self._true_hydro_ztd[:,:,layer] = 1e-6 * 0.5 * (self._zs[-1] - self._zs[layer]) * _p[layer]
+        
+        self._true_wet_refr = 2 * np.ones(self._t.shape)
+        self._true_wet_refr[:,3:] = 4
+    
+    def interpWet(self):
+        _ifWet = rgi((self._ys, self._xs, self._zs), self._true_wet_refr)
+        return _ifWet
+    def interpHydro(self):
+        _ifHydro = rgi((self._ys, self._xs, self._zs), self._true_hydro_refr)
+        return _ifHydro
+
+
+@pytest.fixture
+def model():
+    return MockWeatherModel()
+
+
+@pytest.fixture
+def setup_fake_raytracing():
+    '''This sets up a fake orbit for the weather model'''
+
+    import isce3.ext.isce3 as isce
+    from isce3.core import DateTime, TimeDelta
+
+    lat0 = _LAT0 # degrees
+    lon0 = _LON0
+    hsat = 700000.
+    omega = _OMEGA # degrees
+    Nvec = 30
+
+    t0 = DateTime("2017-02-12T01:12:30.0")
+
+    elp = isce.core.Ellipsoid(6378137.,.0066943799901)
+    look = isce.core.LookSide.Left
+
+    sat_hgt = elp.a + hsat
+    sat_lat = np.sin(np.radians(lat0))
+    clat = np.cos(np.radians(lat0))
+
+    svs = []
+    for k in range(Nvec):
+        dt = TimeDelta(100 * k)
+        lon = lon0 + omega * dt.total_seconds()
+
+        pos = []
+        pos.append(sat_hgt * clat * np.cos(np.radians(lon)))
+        pos.append(sat_hgt * clat * np.sin(np.radians(lon)))
+        pos.append(sat_hgt * sat_lat)
+
+        vel = []
+        vel.append(-omega * pos[1])
+        vel.append(omega * pos[0])
+        vel.append(0.)
+
+        epoch = t0 + dt
+
+        svs.append(
+            isce.core.StateVector(epoch,pos, vel)
+        )
+
+    orb = isce.core.Orbit(svs)
+
+    return orb, look, elp, sat_hgt
 
+def solve(R, hsat, ellipsoid, side='left'):
+    # temp = 1.0 + hsat/ellipsoid.a
+    # temp1 = R/ellipsoid.a
+    # temp2 = R/(ellipsoid.a + hsat)
+    t2 = (np.square(hsat) + np.square(R)) - np.square(ellipsoid.a)
+    # cosang = 0.5 * (temp + (1.0/temp) - temp1 * temp2)
+    cosang = 0.5 * t2 / (R * hsat)
+    angdiff = np.arccos(cosang)
+    if side=='left':
+        x = _LAT0 + angdiff
+    else:
+        x = _LAT0 - angdiff
+    return x
 
-def test_Raytracing():
-    lats = np.array([-90, 0, 0, 90])
-    lons = np.array([-90, 0, 90, 180])
-    hgts = np.array([-10, 0, 10, 1000])
 
-    unit_vecs = np.array([[0,0,-1], [1,0,0], [0,1,0], [0,0,1]])
+@pytest.mark.isce3
+def test_llhs(setup_fake_raytracing, model):
+    orb, look_dir, elp, sat_hgt = setup_fake_raytracing
+    llhs = []
+    for k in range(20):
+        tinp = 5 + k * 2
+        rng = 800000 + 1000 * k
+        expLon = _LON0 + _OMEGA * tinp
+        geocentricLat = solve(rng, sat_hgt, elp)
 
-    z = Raytracing()
-    with pytest.raises(RuntimeError):
-        z.setPoints(lats=None)
+        xyz = [
+            elp.a * np.cos(geocentricLat) * np.cos(expLon),
+            elp.a * np.cos(geocentricLat) * np.sin(expLon),
+            elp.a * np.sin(geocentricLat)
+        ]
+        llh = elp.xyz_to_lon_lat(xyz)
+        llhs.append(llh)
 
-    z.setPoints(lats=lats, lons=lons, heights = hgts)
-    assert z._lats.shape == (4,)
-    assert z._lats.shape == z._lons.shape
-    assert np.allclose(z._heights, hgts)
+    assert len(llhs) == 20
 
+@pytest.mark.skip
+@pytest.mark.isce3
+def test_build_cube_ray(setup_fake_raytracing, model):
+    import isce3.ext.isce3 as isce
+    from RAiDER.losreader import state_to_los
 
-def test_toa():
-    lats = np.array([0, 0, 0, 0])
-    lons = np.array([0, 180, 90, -90])
-    hgts = np.array([0, 0, 0, 0])
+    orb, look_dir, elp, _ = setup_fake_raytracing
+    m = model
+    m.load_weather()
 
-    z = Raytracing()
-    z.setPoints(lats=lats, lons=lons, heights=hgts)
+    ys = np.arange(-1,1) + _LAT0
+    xs = np.arange(-1,1) + _LON0
+    zs = np.arange(0, 1e5, 1e3)
 
-    # Mock xyz
-    z._xyz = np.array([[6378137.0, 0.0, 0.0],
-                       [-6378137.0, 0.0, 0.0],
-                       [0.0, 6378137.0, 0.0],
-                       [0.0, -6378137.0, 0.0]])
-    z._look_vecs = np.array([[1, 0, 0],
-                            [-1, 0, 0],
-                            [0, 1, 0],
-                            [0, -1, 0]])
-    toppts = np.array([[6388137.0, 0.0, 0.0],
-                       [-6388137.0, 0.0, 0.0],
-                       [0.0, 6388137.0, 0.0],
-                       [0.0, -6388137.0, 0.0]])
+    _Y, _X, _Z = np.meshgrid(ys, xs, zs)
 
-    topxyz = z.getIntersectionWithHeight(10000.)
+    out_true = np.zeros(_Y.shape)
+    t0 = orb.start_time
+    tm1 = orb.end_time
+    ts = np.arange(t0, tm1 + orb.time.spacing, orb.time.spacing)
+    tlist = [orb.reference_epoch + isce.core.TimeDelta(dt) for dt in ts]
+    ts = np.broadcast_to(tlist, (1, len(tlist))).T
+    svs = np.hstack([ts, orb.position, orb.velocity])
 
-    assert np.allclose(topxyz, toppts)
+    #TODO: Check that the look vectors are not nans
+    lv, xyz = state_to_los(svs, np.stack([_Y.ravel(), _X.ravel(), _Z.ravel()], axis=-1),out="ecef")
+    out = _build_cube_ray(xs, ys, zs, orb, look_dir, CRS(4326), CRS(4326), [m.interpWet(), m.interpHydro()], elp=elp)
+    assert out.shape == out_true.shape

test/test_scenario_3.py

@@ -11,6 +11,7 @@
 
 
 
+@pytest.mark.isce3
 def test_scenario_3():
     SCENARIO_DIR = os.path.join(TEST_DIR, "scenario_3")