more comments/docstrings

add environment yaml to help create clean/minimal environment
add bounds check to margin application

.gitignore

@@ -151,4 +151,12 @@ cython_debug/
 #  option (not recommended) you can uncomment the following to ignore the entire idea folder.
 #.idea/
 
-
+# unit test files/directories
+tests/data/geo_cslc_s1.yaml
+tests/data/2022-10-16_0000_Rosamond-corner-reflectors.csv
+tests/data/S1A_IW_SLC__1SDV_20221016T015043_20221016T015111_045461_056FC0_6681.zip
+tests/data/orbits
+tests/data/test_dem.tiff
+tests/data/test_burst_map.sqlite3
+tests/product
+tests/scratch

compass_env.yml

@@ -0,0 +1,19 @@
+name: compass
+channels:
+  - conda-forge
+  - nodefaults
+dependencies:
+  # isce3 can be added if forge version is sufficient
+  - gdal>=3.0
+  - h5py
+  - lxml
+  - numpy
+  - pandas
+  - pyproj
+  - pytest
+  - python>=3.9
+  - requests
+  - ruamel.yaml
+  - scipy
+  - shapely
+  - yamale

src/compass/utils/runconfig.py

@@ -206,16 +206,14 @@ def runconfig_to_bursts(cfg: SimpleNamespace) -> list[Sentinel1BurstSlc]:
         for pol, i_subswath in pol_subswath_index_pairs:
 
             # loop over burst objs extracted from SAFE zip
-            loaded_bursts =  load_bursts(safe_file, orbit_path, i_subswath, pol)
-            for burst in loaded_bursts:
+            for burst in load_bursts(safe_file, orbit_path, i_subswath, pol):
                 # get burst ID
                 burst_id = str(burst.burst_id)
 
                 # include ALL bursts if no burst IDs given
                 # is burst_id wanted? skip if not given in config
                 if (cfg.input_file_group.burst_id is not None and
                         burst_id not in cfg.input_file_group.burst_id):
-                    print(burst_id, cfg.input_file_group.burst_id, burst_id in cfg.input_file_group.burst_id)
                     continue
 
                 # get polarization and save as tuple with burst ID
@@ -239,7 +237,6 @@ def runconfig_to_bursts(cfg: SimpleNamespace) -> list[Sentinel1BurstSlc]:
                 # if reference burst, ok to add if id+pol combo does not exist
                 # no duplicate id+pol combos for reference bursts
                 if not_ref or not burst_id_pol_exist:
-                    print(f'{burst_id} found')
                     burst_ids_found.append(burst_id)
                     bursts.append(burst)
 

tests/conftest.py

@@ -9,25 +9,48 @@
 
 from compass.utils import iono
 
+
 def download_if_needed(local_path):
-    # check if test inputs and reference files exists; download if not found.
+    '''
+    Check if given path to file exists. Download if it from zenodo does not.
+
+    Parameters
+    ----------
+    local_path: str
+        Path to file
+    '''
+    # return if file is found
     if os.path.isfile(local_path):
         return
 
     check_internet_connection()
 
-    dataset_url = 'https://zenodo.org/record/7668411/files/'
     dst_dir, file_name = os.path.split(local_path)
+
+    # create destination directory if it does not exist
     if dst_dir:
         os.makedirs(dst_dir, exist_ok=True)
+
+    # download data
+    dataset_url = 'https://zenodo.org/record/7668411/files/'
     target_url = f'{dataset_url}/{file_name}'
     with open(local_path, 'wb') as f:
         f.write(requests.get(target_url).content)
 
+
 @pytest.fixture(scope="session")
-def unit_test_paths():
-    test_paths = types.SimpleNamespace()
+def geocode_slc_params():
+    '''
+    Parameters to be used by geocode SLC unit test
+
+    Returns
+    -------
+    test_params: SimpleNamespace
+        SimpleNamespace containing geocode SLC unit test parameters
+    '''
+    test_params = types.SimpleNamespace()
 
+    # burst ID and date of burst
     burst_id = 't064_135523_iw2'
     b_date = '20221016'
 
@@ -39,34 +62,40 @@ def unit_test_paths():
 
     # paths for template and actual runconfig
     gslc_template_path = f'{test_data_path}/geo_cslc_s1_template.yaml'
-    test_paths.gslc_cfg_path = f'{test_data_path}/geo_cslc_s1.yaml'
+    test_params.gslc_cfg_path = f'{test_data_path}/geo_cslc_s1.yaml'
 
     # read runconfig template, replace pieces, write to runconfig
     with open(gslc_template_path, 'r') as f_template, \
-            open(test_paths.gslc_cfg_path, 'w') as f_cfg:
+            open(test_params.gslc_cfg_path, 'w') as f_cfg:
         cfg = f_template.read().replace('@TEST_PATH@', str(test_path)).\
             replace('@DATA_PATH@', test_data_path).\
             replace('@BURST_ID@', burst_id)
         f_cfg.write(cfg)
 
-    # check for files and download as needed
+    # files needed for geocode SLC unit test
     test_files = ['S1A_IW_SLC__1SDV_20221016T015043_20221016T015111_045461_056FC0_6681.zip',
                   'orbits/S1A_OPER_AUX_POEORB_OPOD_20221105T083813_V20221015T225942_20221017T005942.EOF',
                   'test_dem.tiff', 'test_burst_map.sqlite3',
                   '2022-10-16_0000_Rosamond-corner-reflectors.csv']
     test_files = [f'{test_data_path}/{test_file}' for test_file in test_files]
 
-    # parallel download of test files
+    # parallel download of test files (if necessary)
     pool = mp.Pool(len(test_files))
     _ = pool.map(download_if_needed, test_files)
     pool.close()
     pool.join()
 
-    test_paths.corner_coord_csv_path = test_files[-1]
-    test_paths.output_hdf5 = f'{test_path}/product/{burst_id}/{b_date}/{burst_id}_{b_date}.h5'
-    test_paths.grid_group_path = '/science/SENTINEL1/CSLC/grids'
+    # path to file containing corner reflectors
+    test_params.corner_coord_csv_path = test_files[-1]
+
+    # path the output HDF5
+    test_params.output_hdf5 = f'{test_path}/product/{burst_id}/{b_date}/{burst_id}_{b_date}.h5'
 
-    return test_paths
+    # path to groups and datasets in output HDF5
+    test_params.grid_group_path = '/science/SENTINEL1/CSLC/grids'
+    test_params.raster_path = f'{test_params.grid_group_path}/VV'
+
+    return test_params
 
 @pytest.fixture(scope='session')
 def ionex_params(download_data=True):
@@ -81,9 +110,8 @@ def ionex_params(download_data=True):
 
     Returns
     -------
-    tec_file: str
-        Path to local or downloaded TEC file to
-        use in the unit test
+    tec_file: SimpleNamespace
+        SimpleNamespace containing parameters needed for ionex unit test
     '''
     test_params = types.SimpleNamespace()
 

tests/test_s1_geocode_slc.py

@@ -8,27 +8,55 @@
 from compass.utils.geo_runconfig import GeoRunConfig
 
 
-def test_geocode_slc_run(unit_test_paths):
+def test_geocode_slc_run(geocode_slc_params):
     '''
-    run s1_geocode_slc to ensure it does not crash
+    Run s1_geocode_slc to ensure it does not crash
+
+    Parameters
+    ----------
+    geocode_slc_params: SimpleNamespace
+        SimpleNamespace containing geocode SLC unit test parameters
     '''
     # load yaml to cfg
-    cfg = GeoRunConfig.load_from_yaml(unit_test_paths.gslc_cfg_path,
+    cfg = GeoRunConfig.load_from_yaml(geocode_slc_params.gslc_cfg_path,
                                       workflow_name='s1_cslc_geo')
 
     # pass cfg to s1_geocode_slc
     s1_geocode_slc.run(cfg)
 
-def get_nearest_index(arr, val):
+
+def _get_nearest_index(arr, val):
+    '''
+    Find index of element in given array closest to given value
+
+    Parameters
+    ----------
+    arr: np.ndarray
+        1D array to be searched
+    val: float
+        Number to be searched for
+
+    Returns
+    -------
+    _: int
+        Index of element in arr where val is closest
+    '''
     return np.abs(arr - val).argmin()
 
-def get_reflectors_extents_slice(unit_test_paths, margin=50):
+
+def _get_reflectors_bounding_slice(geocode_slc_params):
     '''
-    get max and min lat, lon
+    Get latitude, longitude slice that contains all the corner reflectors in
+    CSV list of corner reflectors
+
+    Parameters
+    ----------
+    geocode_slc_params: SimpleNamespace
+        SimpleNamespace containing geocode SLC unit test parameters
     '''
     # extract from HDF5
-    with h5py.File(unit_test_paths.output_hdf5, 'r') as h5_obj:
-        grid_group = h5_obj[unit_test_paths.grid_group_path]
+    with h5py.File(geocode_slc_params.output_hdf5, 'r') as h5_obj:
+        grid_group = h5_obj[geocode_slc_params.grid_group_path]
 
         # create projection to covert from UTM to LLH
         epsg = int(grid_group['projection'][()])
@@ -43,33 +71,50 @@ def get_reflectors_extents_slice(unit_test_paths, margin=50):
         lats = np.array([np.degrees(proj.inverse([x_coords_utm[0], y, 0])[1])
                          for y in y_coords_utm])
 
+        # get array shape for later check of slice with margins applied
+        height, width = h5_obj[geocode_slc_params.raster_path].shape
+
     # extract all lat/lon corner reflector coordinates
     corner_lats = []
     corner_lons = []
-    with open(unit_test_paths.corner_coord_csv_path, 'r') as csvfile:
+    with open(geocode_slc_params.corner_coord_csv_path, 'r') as csvfile:
         corner_reader = csv.DictReader(csvfile)
         for row in corner_reader:
             corner_lats.append(float(row['Latitude (deg)']))
             corner_lons.append(float(row['Longitude (deg)']))
 
-    i_max_lat = get_nearest_index(lats, np.max(corner_lats))
-    i_min_lat = get_nearest_index(lats, np.min(corner_lats))
-    i_max_lon = get_nearest_index(lons, np.max(corner_lons))
-    i_min_lon = get_nearest_index(lons, np.min(corner_lons))
+    # find nearest index for min/max of lats/lons and apply margin
+    # apply margin to bounding box and ensure raster bounds are not exceeded
+    # application of margin y indices reversed due descending order lats vector
+    margin = 50
+    i_max_y = max(_get_nearest_index(lats, np.max(corner_lats)) - margin, 0)
+    i_min_y = min(_get_nearest_index(lats, np.min(corner_lats)) + margin,
+                  height - 1)
+    i_max_x = min(_get_nearest_index(lons, np.max(corner_lons)) + margin,
+                  width - 1)
+    i_min_x = max(_get_nearest_index(lons, np.min(corner_lons)) - margin, 0)
 
-    return np.s_[i_max_lat - margin:i_min_lat + margin,
-                 i_min_lon - margin:i_max_lon + margin]
+    # return as slice
+    # y indices reversed to account for descending order lats vector
+    return np.s_[i_max_y:i_min_y, i_min_x:i_max_x]
 
-def test_geocode_slc_validate(unit_test_paths):
+
+def test_geocode_slc_validate(geocode_slc_params):
     '''
-    check for reflectors in geocoded output
+    Check for reflectors in geocoded output
+
+    Parameters
+    ----------
+    geocode_slc_params: SimpleNamespace
+        SimpleNamespace containing geocode SLC unit test parameters
     '''
-    s_ = get_reflectors_extents_slice(unit_test_paths)
+    # get slice where corner reflectors should be
+    s_ = _get_reflectors_bounding_slice(geocode_slc_params)
 
-    with h5py.File(unit_test_paths.output_hdf5, 'r') as h5_obj:
-        src_path = f'{unit_test_paths.grid_group_path}/VV'
-        arr = h5_obj[src_path][()][s_]
-        print(arr.shape, s_)
+    # slice raster array
+    with h5py.File(geocode_slc_params.output_hdf5, 'r') as h5_obj:
+        arr = h5_obj[geocode_slc_params.raster_path][()][s_]
 
+    # check for bright spots in sliced array
     corner_reflector_threshold = 3e3
     assert np.any(arr > corner_reflector_threshold)