convert_waveform_coordinates: apply source_redshift to input masses (mu1/mu2 frame)

.travis/test-core-units.sh

@@ -111,6 +111,8 @@ FILES=(
   "$C/test/test_lisa_ini_contract.py"
   "$C/test/test_tracer_placement_gp.py"
   "$C/test/waveforms/test_uv_symmetry.py"
+  # -- coordinate conversion (source-frame CIP input -> detector-frame fit coordinates)
+  "$C/test/test_convert_coordinates_source_redshift.py"
 )
 
 # A manifest entry that stops existing is a SILENT no-op: the gate keeps passing while
@@ -163,6 +165,8 @@ done
 #            fed odeint's float probe to len(x) pdfs; the t_ref wiring in all three ILE drivers)
 #   378/366  + test_response_order.py (8 tests: SNR tightening, Halton independence,
 #            compound-axis semantics, reference resolution/tail charging, and bank preflight)
+#   393/381  + test_convert_coordinates_source_redshift.py (13 tests: vectorized
+#            convert_waveform_coordinates vs per-row extract_param at source_redshift 0 and 0.3)
 #
 # RAISE these when files are added: a floor left at the old value passes while covering less,
 # which is the failure this gate exists to catch.
@@ -178,12 +182,12 @@ done
 # runner's closure the count falls back to 347 and still passes.  Pinning 350 would turn an
 # unrelated dependency change into a red gate.
 # Two XML/grid template-finalization regressions, with no added skips.
-EXPECTED_TESTS=380
+EXPECTED_TESTS=393
 # Outcomes, not just exit status: a collection floor cannot see a test that collects, runs and
 # asserts nothing, and a pytest.skip can quietly absorb a lost gate.  The 12 skips are
 # environment legs -- cupy in test_seeding_reproducibility, device legs in
 # test_dslice_device_native, and the xfail in test_uv_symmetry.
-EXPECTED_PASSED=368
+EXPECTED_PASSED=381
 MAX_SKIPPED=12
 
 # The floors must be INTEGERS, and this is checked rather than assumed.  `[ 347 -lt FOO ]` does

MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py

@@ -5919,6 +5919,17 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_
       - source_redshift: if nonzero, convert m1 -> m1 (1+z)=m_z, as fit is done in the detector frame.  We are **assuming source-frame sampling**
     """
     x_out = np.zeros( (len(x_in), len(coord_names) ) )
+    # Source-frame input: move the mass coordinates to the detector frame ONCE, here, so every vectorized branch below
+    # and the per-row fallback see detector-frame masses.  Other mass-dimensional inputs (e.g. mu1, mu2) do not scale
+    # by (1+z); for those, the per-row fallback applies the redshift to P.m1, P.m2 instead.
+    redshift_applied_to_input = False
+    if source_redshift:
+        mass_names_in = [p for p in ['mc', 'mc_ecc', 'm1', 'm2', 'mtot'] if p in low_level_coord_names]
+        if mass_names_in:
+            x_in = np.array(x_in, dtype=float)   # copy: never rescale the caller's array
+            for p in mass_names_in:
+                x_in[:, low_level_coord_names.index(p)] *= (1+source_redshift)
+            redshift_applied_to_input = True
     # Check for trivial identity transformations and do those by direct copy, then remove those from the list of output coord names
     coord_names_reduced = coord_names.copy() 
     for p in low_level_coord_names:
@@ -6341,13 +6352,6 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_
                 coord_names_reduced.remove('DeltaLambdaTilde')
 
 
-    # perform any mass conversions needed, so output in detector frame given input in source frame
-    if source_redshift:
-        for name in ['mc', 'm1', 'm2']:
-            if name in coord_names:
-                indx_name = coord_names.index(name)
-                x_out[:,indx_name] *= (1+source_redshift)
-
     # return if we don't need to do any more conversions (e.g., if we only have --parameter specification)
     if len(coord_names_reduced)<1:
         return x_out
@@ -6361,8 +6365,9 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_
             if low_level_coord_names[indx] != 'chi_pavg':
                 P.assign_param( low_level_coord_names[indx], x_in[indx_out,indx])            
         # Apply redshift: assume input is source-frame mass, convert m1 -> m1(1+z) = m1_z, as fit used detector frame
-        P.m1 = P.m1*(1+source_redshift)
-        P.m2 = P.m2*(1+source_redshift)
+        if source_redshift and not redshift_applied_to_input:
+            P.m1 = P.m1*(1+source_redshift)
+            P.m2 = P.m2*(1+source_redshift)
         for indx in np.arange(len(coord_names_reduced)):
             p = coord_names_reduced[indx]
             indx_p_out= coord_names.index(p)

MonteCarloMarginalizeCode/Code/test/test_convert_coordinates_source_redshift.py

@@ -0,0 +1,84 @@
+"""convert_waveform_coordinates(..., source_redshift=z) must return detector-frame values.
+
+Reference: the per-row ChooseWaveformParams path the function itself falls back to
+(assign source-frame low-level coordinates, scale P.m1, P.m2 by 1+z, extract_param).
+Masses are in Msun, as CIP passes them.
+"""
+import numpy as np
+import pytest
+
+from RIFT import lalsimutils as lsu
+
+N = 50
+
+
+def _draw(low_level_coord_names, rng):
+    ranges = {
+        'mc': (5., 40.), 'delta_mc': (0.05, 0.85), 'eta': (0.1, 0.249),
+        's1z': (-0.9, 0.9), 's2z': (-0.9, 0.9),
+        'chi1': (0., 0.99), 'chi2': (0., 0.99),
+        'cos_theta1': (-1., 1.), 'cos_theta2': (-1., 1.),
+        'phi1': (0., 2*np.pi), 'phi2': (0., 2*np.pi),
+        's1z_bar': (-0.9, 0.9), 's2z_bar': (-0.9, 0.9),
+        'chi1_perp_bar': (0., 0.9), 'chi2_perp_bar': (0., 0.9),
+        'lambda1': (0., 1000.), 'lambda2': (0., 1000.),
+    }
+    return np.column_stack([rng.uniform(*ranges[p], N) for p in low_level_coord_names])
+
+
+def _per_row(x_in, coord_names, low_level_coord_names, z):
+    out = np.zeros((len(x_in), len(coord_names)))
+    for i, row in enumerate(x_in):
+        P = lsu.ChooseWaveformParams()
+        for p, val in zip(low_level_coord_names, row):
+            P.assign_param(p, val)
+        P.m1 *= 1 + z
+        P.m2 *= 1 + z
+        out[i] = [P.extract_param(p) for p in coord_names]
+    return out
+
+
+CASES = [
+    # spherical spins: vectorized mu1, mu2 (the case PR #204's RF features use)
+    (['mc', 'delta_mc', 'mu1', 'mu2', 'xi', 'chiMinus', 's1x', 's1y', 's2x', 's2y'],
+     ['mc', 'delta_mc', 'chi1', 'cos_theta1', 'phi1', 'chi2', 'cos_theta2', 'phi2']),
+    # aligned cartesian: vectorized mu1, mu2
+    (['mu1', 'mu2', 'delta_mc', 'chiMinus', 'mc'], ['mc', 'delta_mc', 's1z', 's2z']),
+    # vectorized m1, m2 from mc, eta
+    (['mc', 'm1', 'm2', 'eta', 'xi'], ['mc', 'eta', 's1z', 's2z']),
+    # pseudo-cylindrical spins: vectorized mu1, mu2.  s*z_bar precedes chi*_perp_bar because
+    # assign_param('chi1_perp_bar') reads the current s1z; the reverse order builds a different spin.
+    (['mu1', 'mu2', 'delta_mc', 'xi', 'chiMinus', 'chi_p'],
+     ['mc', 'delta_mc', 's1z_bar', 'chi1_perp_bar', 'phi1', 's2z_bar', 'chi2_perp_bar', 'phi2']),
+    # tidal: vectorized mu1, mu2 next to LambdaTilde
+    (['mu1', 'mu2', 'delta_mc', 'LambdaTilde', 'DeltaLambdaTilde'],
+     ['mc', 'delta_mc', 's1z', 's2z', 'lambda1', 'lambda2']),
+    # mtot and q go through the per-row fallback: the redshift must be applied once, not twice
+    (['mc', 'mtot', 'q', 'm1'], ['mc', 'delta_mc', 's1z', 's2z']),
+]
+
+
+@pytest.mark.parametrize("z", [0.0, 0.3])
+@pytest.mark.parametrize("coord_names,low_level_coord_names", CASES)
+def test_vectorized_matches_per_row(coord_names, low_level_coord_names, z):
+    rng = np.random.default_rng(20261004)
+    x_in = _draw(low_level_coord_names, rng)
+    x_in_before = x_in.copy()
+    got = lsu.convert_waveform_coordinates(x_in, coord_names=coord_names,
+                                           low_level_coord_names=low_level_coord_names,
+                                           source_redshift=z)
+    ref = _per_row(x_in, coord_names, low_level_coord_names, z)
+    np.testing.assert_array_equal(x_in, x_in_before)  # caller's array is not rescaled in place
+    np.testing.assert_allclose(got, ref, rtol=1e-9, atol=1e-12,
+                               err_msg=f"z={z} coords={coord_names}")
+
+
+def test_mass_columns_scale_on_every_row():
+    rng = np.random.default_rng(1)
+    low = ['mc', 'eta', 's1z', 's2z']
+    x_in = _draw(low, rng)
+    names = ['mc', 'm1', 'm2']
+    z0 = lsu.convert_waveform_coordinates(x_in, coord_names=names, low_level_coord_names=low)
+    z1 = lsu.convert_waveform_coordinates(x_in, coord_names=names, low_level_coord_names=low,
+                                          source_redshift=0.3)
+    np.testing.assert_allclose(z1, 1.3*z0, rtol=1e-12)