Merge branch 'isce-framework:develop' into develop

Author: gshiroma

python/packages/nisar/pointing/doppler_lut_from_raw.py

@@ -21,7 +21,8 @@
 def doppler_lut_from_raw(raw, *, freq_band='A', txrx_pol=None,
                          orbit=None, attitude=None, ant=None,
                          dem=None, num_rgb_avg=8, az_block_dur=4.0,
-                         time_interval=2.0, dop_method='CDE', subband=False,
+                         rangeline_limit=None, time_interval=2.0,
+                         dop_method='CDE', subband=False,
                          polyfit_deg=3, polyfit=False, exclude_beams=None,
                          out_path='.', plot=False, logger=None):
     """Generates 2-D Doppler LUT as a function of slant range and azimuth time.
@@ -77,6 +78,10 @@ def doppler_lut_from_raw(raw, *, freq_band='A', txrx_pol=None,
     az_block_dur : float, default=4.0
         Azimuth block duration in seconds defining time-domain correlator
         length used in Doppler estimator.
+    rangeline_limit : tuple[int | None, int | None], optional
+        0-based range line [start, stop) indices to limit echo range lines to
+        be processed. Default is all range lines. The start/stop will be
+        limited to within [0, total rangelines)!
     time_interval : float, default=2.0
         Time stamp interval between azimuth blocks in seconds.
         It should not be larger than "az_block_dur".
@@ -269,10 +274,10 @@ def doppler_lut_from_raw(raw, *, freq_band='A', txrx_pol=None,
                 logger.warning(
                     f'Excluded beams {exclude_beams} is out of range of '
                     f'active RX channels {list_rx_active}!'
-                    )
+                )
                 exclude_beams.intersection_update(list_rx_active)
                 logger.warning(
-                        f'Updated list of excluded beams -> {exclude_beams}')
+                    f'Updated list of excluded beams -> {exclude_beams}')
         # initialize DEM object if None
         if dem is None:
             dem = DEMInterpolator()
@@ -417,9 +422,26 @@ def time_dop_est(echo, prf, lag=1, axis=None):
     sr_stop = sr_start + (num_blk_rg - 1) * sr_spacing
     slrg_per_blk = np.linspace(sr_start, sr_stop, num=num_blk_rg)
 
+    # set range line slice of echo
+    rgl_start_stop = [0, tot_pulses]
+    if rangeline_limit is not None:
+        if rangeline_limit[0] is not None:
+            rgl_start_stop[0] = min(
+                max(rangeline_limit[0], 0),
+                tot_pulses - 1)
+        if rangeline_limit[1] is not None:
+            rgl_start_stop[1] = max(
+                min(rangeline_limit[1], tot_pulses),
+                rgl_start_stop[0] + 1)
+    logger.info('[start, stop) range line index to be processed -> '
+                f'{rgl_start_stop}')
+    # update number of range lines to be processed
+    num_rgls = rgl_start_stop[1] - rgl_start_stop[0]
+    logger.info(f'Number of range lines to be processed -> {num_rgls}')
+
     # form the blocks of range lines / azimuth bins
     len_az_blk_dur, len_tm_int, num_blk_az = _get_az_block_interval_len(
-        tot_pulses, az_block_dur, prf, time_interval)
+        num_rgls, az_block_dur, prf, time_interval)
 
     logger.info(
         f'Final full azimuth block duration -> {len_az_blk_dur/prf:.3f} (sec)')
@@ -434,7 +456,9 @@ def time_dop_est(echo, prf, lag=1, axis=None):
     logger.info(f'Total number of azimuth blocks -> {num_blk_az}')
 
     slice_lines = _azblk_slice_gen(
-        tot_pulses, len_az_blk_dur, len_tm_int, num_blk_az)
+        num_rgls, len_az_blk_dur, len_tm_int,
+        num_blk_az, idx_start=rgl_start_stop[0]
+    )
 
     # parse valid subswath index for all range lines used later
     valid_sbsw_all = raw.getSubSwaths(freq_band, txrx_pol[0])
@@ -447,8 +471,10 @@ def time_dop_est(echo, prf, lag=1, axis=None):
 
     # initialize the azimuth time block and set an intermediate var
     half_az_blk_dur = (len_az_blk_dur - 1) / 2
-    az_time_blk = np.full(num_blk_az, az_time[0] + half_az_blk_dur * pri,
-                          dtype=float)
+    az_time_blk = np.full(
+        num_blk_az,
+        az_time[rgl_start_stop[0]] + half_az_blk_dur * pri,
+        dtype=float)
     tm_int_pri_prod = len_tm_int * pri
 
     # doppler centroid map is azimuth block by slant-range block
@@ -479,7 +505,7 @@ def time_dop_est(echo, prf, lag=1, axis=None):
             echo, rgb_limits = form_single_tap_dbf_echo(
                 raw_dset, slice_line, el_trans, az_trans,
                 pos_mid, vel_mid, quat_mid, sr_lsp, dem
-                )
+            )
             if exclude_beams is not None:
                 mask_bad = np.zeros(rgb_limits[-1], dtype='bool')
                 for n_b in exclude_beams:
@@ -489,7 +515,7 @@ def time_dop_est(echo, prf, lag=1, axis=None):
                         f'For AZ block # {n_azblk + 1}, exclude slant ranges '
                         f'(m, m) within -> ({sr_lsp[rgb_slice.start] :.3f}, '
                         f'{sr_lsp[rgb_slice.stop - 1]:.3f}) '
-                        )
+                    )
         else:  # single channel
             echo = raw_dset[slice_line]
 
@@ -844,7 +870,8 @@ def _get_az_block_interval_len(num_pls: int, az_block_dur: float, prf: float,
     # if the block_dur + interval is too large then raise an exception!
     if (len_tm_int + len_az_blk_dur) > num_pls:
         raise ValueError(
-            'Sum of azimuth block duration and time interval is large than '
+            f'Sum of azimuth block duration {len_az_blk_dur / prf} (sec) '
+            f'and time interval {len_tm_int / prf} (sec) is longer than '
             f'echo duration {(num_pls - 1) / prf} (sec)!'
         )
     # get number of blocks which must be at least 2!
@@ -858,15 +885,14 @@ def _get_az_block_interval_len(num_pls: int, az_block_dur: float, prf: float,
 
 
 def _azblk_slice_gen(num_pls: int, len_az_blk_dur: int, len_tm_int: int,
-                     num_blk_az: int):
+                     num_blk_az: int, idx_start: int = 0):
     """Slice index generator for azimuth blocks/range lines"""
     # generate slice for azimuth block indexing
-    i_str = 0
-    i_stp = len_az_blk_dur
-    for bb in range(num_blk_az):
+    i_max = idx_start + num_pls
+    for i_str in range(
+            idx_start, idx_start + num_blk_az * len_tm_int, len_tm_int):
+        i_stp = min(i_str + len_az_blk_dur, i_max)
         yield slice(i_str, i_stp)
-        i_str += len_tm_int
-        i_stp = min(i_str + len_az_blk_dur, num_pls)
 
 
 def _form_mask_valid_range(tot_rgbs, rgb_valid_sbsw):

python/packages/nisar/pointing/el_null_range_from_raw_ant.py

@@ -55,6 +55,7 @@ def el_null_range_from_raw_ant(
         orbit=None,
         attitude=None,
         az_block_dur=3.0,
+        rangeline_limit=None,
         apply_caltone=False,
         imbalances_right2left=None,
         sample_delays_wrt_left=None,
@@ -98,6 +99,10 @@ def el_null_range_from_raw_ant(
         duration of echo if it is too large.
         The min block duration must be equal or larger than nominal mean
         PRI (pulse repetition interval).
+    rangeline_limit : tuple[int | None, int | None], optional
+        0-based range line [start, stop) indices to limit echo range lines to
+        be processed. Default is all range lines. The start/stop will be
+        limited to within [0, total rangelines)!
     apply_caltone : bool, default=False
         Apply caltone coefficients to RX channels prior to Null formation.
     imbalances_right2left : sequence or array of complex float, optional
@@ -226,9 +231,25 @@ def el_null_range_from_raw_ant(
     if raw_dset.ndim != 3:
         raise RuntimeError(
             'Expected Multi-channel Raw echo aka Diagnostic Mode #2 (DM2)!')
-    num_channels, num_rgls, num_rgbs = raw_dset.shape
+    num_channels, num_rgls_tot, num_rgbs = raw_dset.shape
     logger.info('Shape of the echo data (channels, pulses, ranges) -> '
-                f'({num_channels, num_rgls, num_rgbs})')
+                f'({num_channels, num_rgls_tot, num_rgbs})')
+    # set range line slice of echo
+    rgl_start_stop = [0, num_rgls_tot]
+    if rangeline_limit is not None:
+        if rangeline_limit[0] is not None:
+            rgl_start_stop[0] = min(
+                max(rangeline_limit[0], 0),
+                num_rgls_tot - 1)
+        if rangeline_limit[1] is not None:
+            rgl_start_stop[1] = max(
+                min(rangeline_limit[1], num_rgls_tot),
+                rgl_start_stop[0] + 1)
+    logger.info('[start, stop) range line index to be processed -> '
+                f'{rgl_start_stop}')
+    # update number of range lines to be processed
+    num_rgls = rgl_start_stop[1] - rgl_start_stop[0]
+    logger.info(f'Number of range lines to be processed -> {num_rgls}')
 
     # number of nulls
     num_nulls = num_channels - 1
@@ -277,8 +298,11 @@ def el_null_range_from_raw_ant(
     # [1, num_rgls]
     num_rgl_block = round(az_block_dur * prf)
     if num_rgl_block > num_rgls:
-        logger.warning('Azimuth block duration exceeds max duration. '
-                       'It will be limited to total azimuth duration!')
+        logger.warning(
+            f'Azimuth block duration w/ lines # {num_rgl_block} exceeds max '
+            f'duration w/ lines # {num_rgls}. It will be limited to total '
+            'azimuth duration!'
+        )
         num_rgl_block = num_rgls
     if num_rgl_block < 1:
         raise ValueError('Azimuth block duration is smaller than mean PRI!')
@@ -393,7 +417,9 @@ def el_null_range_from_raw_ant(
 
     # get pair number and range line slice generator
     pairnum_rglslice = _pair_num_rgl_slice_gen(
-        num_rgls, num_nulls, num_azimuth_block, num_rgl_block)
+        num_rgls, num_nulls, num_azimuth_block, num_rgl_block,
+        idx_start=rgl_start_stop[0]
+    )
 
     # check if there is matplotlib package needed for plotting if requested
     if plot and plt is None:
@@ -583,7 +609,7 @@ def form_overlap_antenna_pairs(ant, list_rx, rx_pol, is_half_p2p=False):
 
 
 def _pair_num_rgl_slice_gen(num_rgls, num_nulls, num_azimuth_block,
-                            num_rgl_block):
+                            num_rgl_block, idx_start=0):
     """Generates pair number and respective range line slice.
 
     Parameters
@@ -596,6 +622,8 @@ def _pair_num_rgl_slice_gen(num_rgls, num_nulls, num_azimuth_block,
         Number of azimuth blocks.
     num_rgl_block : int
         Number of range lines per azimuth block
+    idx_start : int, default=0
+        Start range line index.
 
     Yields
     ------
@@ -607,16 +635,17 @@ def _pair_num_rgl_slice_gen(num_rgls, num_nulls, num_azimuth_block,
         AZ block number starting from 1
 
     """
-    i_start = 0
-    i_stop = 0
-    for cc in range(1, num_azimuth_block + 1):
+    for cc, i_start in enumerate(
+            range(idx_start,
+                  idx_start + num_azimuth_block * num_rgl_block, num_rgl_block
+                  ),
+            start=1):
         if cc == num_azimuth_block:
-            i_stop = num_rgls
+            i_stop = idx_start + num_rgls
         else:
-            i_stop += num_rgl_block
+            i_stop = i_start + num_rgl_block
         for nn in range(num_nulls):
             yield nn, slice(i_start, i_stop), cc
-        i_start += num_rgl_block
 
 
 def _is_null_valid(rgb_p2p, rgb_valid_sbsw):

python/packages/nisar/pointing/el_rising_edge_from_raw_ant.py

@@ -23,8 +23,8 @@
 
 def el_rising_edge_from_raw_ant(raw, ant, *, dem_interp=None,
                                 freq_band='A', txrx_pol=None,
-                                orbit=None, attitude=None,
-                                az_block_dur=3.0, beam_num=None,
+                                orbit=None, attitude=None, az_block_dur=3.0,
+                                rangeline_limit=None, beam_num=None,
                                 dbf_pow_norm=True, apply_weight=True,
                                 plot=False, out_path='.', logger=None):
     """
@@ -63,6 +63,10 @@ def el_rising_edge_from_raw_ant(raw, ant, *, dem_interp=None,
         duration of echo if it is too large.
         The min block duration must be equal or larger than nominal mean
         PRI (pulse repetition interval).
+    rangeline_limit : tuple[int | None, int | None], optional
+        0-based range line [start, stop) indices to limit echo range lines to
+        be processed. Default is all range lines. The start/stop will be
+        limited to within [0, total rangelines)!
     beam_num : int, optional
         Beam number used for fetching a desired beam from antenna object simply
         for single-channel raw echo. It will be ignored for SweepSAR case!
@@ -223,9 +227,25 @@ def el_rising_edge_from_raw_ant(raw, ant, *, dem_interp=None,
 
     # Get raw dataset
     raw_dset = raw.getRawDataset(freq_band, txrx_pol)
-    num_rgls, num_rgbs = raw_dset.shape
+    num_rgls_tot, num_rgbs = raw_dset.shape
     logger.info('Shape of the Raw echo data (pulses, ranges) -> '
-                f'({num_rgls, num_rgbs})')
+                f'({num_rgls_tot, num_rgbs})')
+    # set range line slice of echo
+    rgl_start_stop = [0, num_rgls_tot]
+    if rangeline_limit is not None:
+        if rangeline_limit[0] is not None:
+            rgl_start_stop[0] = min(
+                max(rangeline_limit[0], 0),
+                num_rgls_tot - 1)
+        if rangeline_limit[1] is not None:
+            rgl_start_stop[1] = max(
+                min(rangeline_limit[1], num_rgls_tot),
+                rgl_start_stop[0] + 1)
+    logger.info('(start, stop) range line index to be processed -> '
+                f'{rgl_start_stop}')
+    # update number of range lines to be processed
+    num_rgls = rgl_start_stop[1] - rgl_start_stop[0]
+    logger.info(f'Number of range lines to be processed -> {num_rgls}')
 
     # EL angle margin  and resolution in radians
     el_margin = np.deg2rad(el_margin_deg)
@@ -509,7 +529,9 @@ def _rgb_dbf2echo(nrgb_dbf):
             tx_trm_info, norm=True)
 
     # generate rangeline slices
-    rgl_slices = _rgl_slice_gen(num_rgls, num_azimuth_block, num_rgl_block)
+    rgl_slices = _rgl_slice_gen(
+        num_rgls, num_azimuth_block, num_rgl_block, idx_start=rgl_start_stop[0]
+    )
 
     # build ElPatternEst object used for 2-way power pattern est from echo
     # within only rising edge with 3rd-order to be used in roll cost function
@@ -594,10 +616,19 @@ def _rgb_dbf2echo(nrgb_dbf):
         idx_el_ant_str = bisect.bisect_left(rx_beams_el.angle, el_ant_min)
         idx_el_ant_stp = bisect.bisect_right(rx_beams_el.angle, el_ant_max)
         el_ant_slice = slice(idx_el_ant_str, idx_el_ant_stp)
+        logger.info(
+            '(min, max) EL angle coverage used from antenna file (deg, deg) '
+            f'-> ({np.rad2deg(el_ant_min):.3f}, {np.rad2deg(el_ant_max):.3f})'
+        )
         ant_el = rx_beams_el.angle[el_ant_slice]
+        if ant_el.size < 2:
+            raise RuntimeError(
+                f'Not enough antenna EL angles w/ size {ant_el.size} to '
+                'polyfit EL vs slant range'
+            )
         logger.info(
-            f'(min, max) antenna EL coverage -> ({np.rad2deg(ant_el[0]):.3f}'
-            f', {np.rad2deg(ant_el[-1]):.3f}) (deg, deg)'
+            '(min, max) antenna EL coverage used in edge method (deg, deg) '
+            f'-> ({np.rad2deg(ant_el[0]):.3f}, {np.rad2deg(ant_el[-1]):.3f})'
         )
         # convert antenna EL angles to slant range and then compute weighting
         # factor based on relative SNR used in weighing cost function of
@@ -613,12 +644,9 @@ def _rgb_dbf2echo(nrgb_dbf):
         # build polyfit coeffs to represent EL (deg) as a function
         # of slant range (km) with either first or second order within
         # a relatively small rising edge region!
-        if ant_el.size < 2:
-            raise RuntimeError(
-                'Not enough antenna EL angles to polyfit EL vs slant range')
         pf_sr2el = np.polyfit(
             ant_sr * 1e-3, np.rad2deg(ant_el), deg=min(ant_el.size - 1, 2)
-            )
+        )
         logger.info('Slantrange-to-elevation polyfit coeffs for AZ block # '
                     f'{n_azblk} -> {pf_sr2el}')
         def _slantrange_to_elevation(
@@ -647,7 +675,7 @@ def _slantrange_to_elevation(
             slice_el_all = _dbf_windows_to_ant_el_coverge_slice(
                 dwp_dbf[s_rgl], wl_dbf[s_rgl], rgb_first, sr_start,
                 sr_spacing, ant_el, _rgb_dbf2echo, _slantrange_to_elevation
-                )
+            )
             # exclude the portion outside the EL coverage by zeroing out
             # its values.
             for cc, slice_el in enumerate(slice_el_all):
@@ -986,7 +1014,7 @@ def _form_ant2way_sweepsar(ant_pat_tx, ant_pat_rx, ant_el, tx_wgt,
     return pow2db(ant_powpat_2way)
 
 
-def _rgl_slice_gen(num_rgls, num_azimuth_block, num_rgl_block):
+def _rgl_slice_gen(num_rgls, num_azimuth_block, num_rgl_block, idx_start=0):
     """Generates pair number and respective range line slice.
 
     Parameters
@@ -997,22 +1025,25 @@ def _rgl_slice_gen(num_rgls, num_azimuth_block, num_rgl_block):
         Number of azimuth blocks.
     num_rgl_block : int
         Number of range lines per azimuth block
+    idx_start : int, default=0
+        Start index.
 
     Yields
     ------
     slice
         Range line slices
 
     """
-    i_start = 0
-    i_stop = 0
-    for cc in range(1, num_azimuth_block + 1):
+    for cc, i_start in enumerate(
+            range(idx_start,
+                  idx_start + num_azimuth_block * num_rgl_block,
+                  num_rgl_block),
+            start=1):
         if cc == num_azimuth_block:
-            i_stop = num_rgls
+            i_stop = idx_start + num_rgls
         else:
-            i_stop += num_rgl_block
+            i_stop = i_start + num_rgl_block
         yield slice(i_start, i_stop)
-        i_start += num_rgl_block
 
 
 def _is_rising_edge_valid(rgb_fl, rgb_valid_sbsw):