Add bulk SAR mocomp operator

docs_input/api/signalimage/radar/sar_bulk_mocomp.rst

@@ -0,0 +1,69 @@
+.. _sar_bulk_mocomp_func:
+
+sar_bulk_mocomp
+################
+
+Apply bulk motion compensation to complex SAR phase-history data in the FX
+domain. sar_bulk_mocomp is currently in the matx::experimental namespace
+because its API is subject to change.
+
+The operator is lazy, allocation-free, and composable with other MatX
+expressions. The final two FX dimensions are interpreted as pulse and frequency
+sample, and any leading dimensions are batches:
+
+.. code-block:: text
+
+   fx:           [batch..., pulses, samples]
+   range_offset: [batch..., pulses]
+   result:       [batch..., pulses, samples]
+
+For sample index :math:`n`, the frequency is
+
+.. math::
+
+   f_n = f_{\mathrm{ref}} +
+         \left(n - \left\lfloor N/2 \right\rfloor\right)\Delta f,
+
+where :math:`f_{\mathrm{ref}}` is
+SarBulkMocompParams::phase_reference_frequency and :math:`\Delta f` is
+SarBulkMocompParams::sample_frequency_spacing. For differential reference
+range :math:`\Delta R_p`, the applied correction is
+
+.. math::
+
+   \exp\left(-j\,\mathrm{sgn}\,\frac{4\pi}{c}\,f_n\,\Delta R_p\right).
+
+The one-range overload assumes an initial reference range of zero. The
+two-range overload changes reference from initial_reference_range to
+target_reference_range and uses
+:math:`\Delta R_p = R_{\mathrm{target},p} - R_{\mathrm{initial},p}`.
+
+Examples
+~~~~~~~~
+
+.. literalinclude:: ../../../../test/00_operators/sar_bulk_mocomp_test.cu
+   :language: cpp
+   :start-after: example-begin sar-bulk-mocomp-1
+   :end-before: example-end sar-bulk-mocomp-1
+   :dedent:
+
+Phase arithmetic follows the range-offset value type. Float ranges use
+single-precision frequency, phase, and trigonometric arithmetic; double ranges
+use double precision for phase construction and argument reduction. With CUDA,
+when double ranges are combined with ``cuda::std::complex<float>`` FX data, the
+bounded, reduced angle is converted to float for a fast single-precision
+``sincospif`` evaluation. Double-complex FX data retains the full double-precision
+trigonometric path. The parameter frequencies are stored as double and narrowed
+for the single-precision range path. Other range types, including fltflt and
+half precision, are not currently supported. FX values must be
+cuda::std::complex<float> or cuda::std::complex<double>.
+
+Host, CUDA, and CUDAJIT executors are supported. Distributed expressions and
+dynamic-rank inputs are not currently supported.
+
+.. versionadded:: head
+
+.. doxygenfunction:: sar_bulk_mocomp(const FxOp &fx, const TargetRangeOp &target_reference_range, const SarBulkMocompParams &params)
+.. doxygenfunction:: sar_bulk_mocomp(const FxOp &fx, const InitialRangeOp &initial_reference_range, const TargetRangeOp &target_reference_range, const SarBulkMocompParams &params)
+.. doxygenstruct:: matx::experimental::SarBulkMocompParams
+   :members:

docs_input/executor_compatibility.rst

@@ -190,6 +190,7 @@ existing operators do not implicitly become distributed operations.
    "round", "|yes|", "|yes|", "|yes|", "|no|", "Element-wise expression."
    "rsqrt", "|yes|", "|yes|", "|yes|", "|no|", "Element-wise expression."
    "sar_bp", "|no|", "|yes|", "|no|", "|no|", "CUDA-only SAR backprojection transform."
+   "sar_bulk_mocomp", "|yes|", "|yes|", "|yes|", "|no|", "Experimental element-wise SAR bulk motion-compensation expression."
    "select", "|yes|", "|yes|", "|yes|", "|no|", "Selection expression."
    "shift", "|yes|", "|yes|", "|yes|", "|no|", "View/reindex expression."
    "sign", "|yes|", "|yes|", "|yes|", "|no|", "Element-wise expression."

examples/sarbp/README.md

@@ -115,6 +115,7 @@ real/imag, row-major), written to `output_image.raw` in this example.
 | `-o OUTPUT` | Output file path (default: input with `.raw` extension) |
 | `-u N` | Range upsample factor via zero-padding (default: 1) |
 | `-w {hamming,none}` | Window for range compression (default: hamming) |
+| `--bulk-mocomp` | Apply bulk motion compensation to FX-domain input using the per-pulse `range_to_mcp` values stored by the CPHD converter |
 | `-b {auto,all,0,N}` | Pulses per processing block. `auto` uses the GPU L2 cache size to choose a block size; `all` and `0` use all pulses (default: auto) |
 | `--image-tiles N` | Process the image as N x N tiles during backprojection (default: 1) |
 | `--taylor-fast-third-order` | Add the third-order range term when using `--precision taylor_fast` |
@@ -124,6 +125,11 @@ real/imag, row-major), written to `output_image.raw` in this example.
 | `--gold FILE` | Validate the output against a golden image (raw `complex<float>`, same format the example writes) and report accuracy metrics (3x3-window correlation and signal-to-error ratio) |
 | `--cmap FILE` | Write the correlation map (raw `float32`) to a file; requires `--gold` |
 
+`--bulk-mocomp` assumes the FX-domain input has not already been bulk motion
+compensated to the stored reference ranges. Applying bulk motion compensation
+to data that has already been motion compensated is not idempotent and results
+in an incorrect phase reference.
+
 The `--precision` flag controls the arithmetic used by the `sar_bp` operator. For spaceborne SAR, `float` does not provide enough precision to store fractional wavelengths at the range-to-MCP magnitudes (hundreds of km), so pure `float` is not sufficient to produce focused images. The available modes are:
 
 - `double` -- full double-precision arithmetic. Most accurate.