From 0f933e7a7a2a1f772eb4feba21e03303d29a80e1 Mon Sep 17 00:00:00 2001 From: Richard O'Shaughnessy Date: Wed, 30 Sep 2026 14:01:38 -0700 Subject: [PATCH 1/4] lalsimutils: phi12 and chi_p_vec, and vectorized in-plane spin coordinates rift_O4d port of oshaughn/research-projects-RIT#203. extract_param gains phi12 (listed in valid_params but not implemented) and chi_p_vec, the vector-sum analogue of chi_p (same A1, A2 weights). The spherical branch of convert_waveform_coordinates builds chi1_perp, chi2_perp, phi12, SOverM2_perp, DeltaOverM2_perp and chi_p_vec vectorized, when spin_convention is "L"; otherwise they fall through to extract_param as before. test_ring_coordinates.py joins the core-units gate; its floors are re-measured at 600/587. Co-Authored-By: Claude Opus 5.5 --- .travis/test-core-units.sh | 10 +++- .../Code/RIFT/lalsimutils.py | 35 ++++++++++- .../Code/test/test_ring_coordinates.py | 60 +++++++++++++++++++ 3 files changed, 102 insertions(+), 3 deletions(-) create mode 100644 MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py diff --git a/.travis/test-core-units.sh b/.travis/test-core-units.sh index 6f4733768..956b3294d 100755 --- a/.travis/test-core-units.sh +++ b/.travis/test-core-units.sh @@ -159,6 +159,8 @@ FILES=( "$C/test/test_lisa_ini_contract.py" "$C/test/test_tracer_placement_gp.py" "$C/test/waveforms/test_uv_symmetry.py" + # -- coordinates: vectorized in-plane spin / ring coordinates agree with extract_param + "$C/test/test_ring_coordinates.py" ) # A manifest entry that stops existing is a SILENT no-op: the gate keeps passing while @@ -320,7 +322,11 @@ 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. # The detector-network sky mapping adds two passing tests and no skips. -EXPECTED_TESTS=592 +# test_ring_coordinates.py adds two more, no skips. RE-MEASURED with it on CIT (ldas-grid; +# `import cupy` FAILS there, numpy backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0: +# per-file 600, junit 603 collected / 590 passed / 13 skipped / 0 failed, of which 3 are +# subtests. So the plugin-free floors are 600/587; the old 592/579 sat 6 below the tree. +EXPECTED_TESTS=600 # 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 13 skips are # environment legs -- cupy in test_seeding_reproducibility, device legs in @@ -329,7 +335,7 @@ EXPECTED_TESTS=592 # (mcsamplerNFlow is an optional dependency and is absent from the IGWN environment), and # the xfail in test_uv_symmetry. test_eos_portfolio_sampler.py adds 12 tests and # test_cip_portfolio_members.py 4, none of them skips. -EXPECTED_PASSED=579 +EXPECTED_PASSED=587 MAX_SKIPPED=13 # The floors must be INTEGERS, and this is checked rather than assumed. `[ 347 -lt FOO ]` does diff --git a/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py b/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py index 0c2f48111..7eaf427a6 100644 --- a/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py +++ b/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py @@ -414,7 +414,7 @@ def lsu_StringFromPNOrder(order): # Class to hold arguments of ChooseWaveform functions # -valid_params = ['m1', 'm2', 's1x', 's1y', 's1z', 's2x', 's2y', 's2z', 'chi1_perp', 'chi2_perp', 'chi1_perp_bar', 'chi2_perp_bar','chi1_perp_u', 'chi2_perp_u', 's1z_bar', 's2z_bar', 'lambda1', 'lambda2', 'theta','phi', 'phiref', 'psi', 'incl', 'tref', 'dist', 'mc', 'mc_ecc', 'eta', 'delta_mc', 'chi1', 'chi2', 'thetaJN', 'phiJL', 'theta1', 'theta2', 'cos_theta1', 'cos_theta2', 'theta1_Jfix', 'theta2_Jfix', 'psiJ', 'beta', 'cos_beta', 'sin_phiJL', 'cos_phiJL', 'phi12', 'phi1', 'phi2', 'LambdaTilde', 'DeltaLambdaTilde', 'lambda_plus', 'lambda_minus', 'q', 'mtot','xi','chiz_plus', 'chiz_minus', 'chieff_aligned','fmin','fref', "SOverM2_perp", "SOverM2_L", "DeltaOverM2_perp", "DeltaOverM2_L", "shu","ampO", "phaseO",'eccentricity','eccentricity_squared','eccentricity_ln', 'chi_pavg','mu1','mu2','eos_table_index','meanPerAno','a6c','E0','p_phi0','hypclass'] +valid_params = ['m1', 'm2', 's1x', 's1y', 's1z', 's2x', 's2y', 's2z', 'chi1_perp', 'chi2_perp', 'chi1_perp_bar', 'chi2_perp_bar','chi1_perp_u', 'chi2_perp_u', 's1z_bar', 's2z_bar', 'lambda1', 'lambda2', 'theta','phi', 'phiref', 'psi', 'incl', 'tref', 'dist', 'mc', 'mc_ecc', 'eta', 'delta_mc', 'chi1', 'chi2', 'thetaJN', 'phiJL', 'theta1', 'theta2', 'cos_theta1', 'cos_theta2', 'theta1_Jfix', 'theta2_Jfix', 'psiJ', 'beta', 'cos_beta', 'sin_phiJL', 'cos_phiJL', 'phi12', 'phi1', 'phi2', 'LambdaTilde', 'DeltaLambdaTilde', 'lambda_plus', 'lambda_minus', 'q', 'mtot','xi','chiz_plus', 'chiz_minus', 'chieff_aligned','fmin','fref', "SOverM2_perp", "SOverM2_L", "DeltaOverM2_perp", "DeltaOverM2_L", "shu","ampO", "phaseO",'eccentricity','eccentricity_squared','eccentricity_ln', 'chi_pavg','chi_p_vec','mu1','mu2','eos_table_index','meanPerAno','a6c','E0','p_phi0','hypclass'] # so far, used for puffball, to prevent insanity (infinite growth) and/or death to downselect # - note we also provide for extrinsic: RA (phi), phiref, psi, just in case we need it in the future @@ -441,6 +441,7 @@ def lsu_StringFromPNOrder(order): "DeltaOverM2_perp" : r"$\Delta_\perp$", "DeltaOverM2_L" : r"$\Delta_{||}$", "SOverM2_perp" : r"$S_\perp$", + "chi_p_vec" : r"$\chi_{p,{\rm vec}}$", "SOverM2_L" : r"$S_{||}$", "eta": r"$\eta$", "chi_eff": r"$\chi_{eff}$", @@ -1390,6 +1391,16 @@ def extract_param(self,p): S2p = (m2**2 * chi2)[:2] Sp = np.max([np.linalg.norm( A1*S1p), np.linalg.norm(A2*S2p)]) return Sp/(A1*m1**2) # divide by term for *larger* BH + if p == 'phi12': + # azimuth of spin 2's in-plane component relative to spin 1's, in [0, 2 pi), L frame + return np.mod(np.arctan2(self.s2y, self.s2x) - np.arctan2(self.s1y, self.s1x), 2*np.pi) + if p == 'chi_p_vec': + # vector-sum (ring) analogue of chi_p: |A1 S1perp + A2 S2perp| / (A1 m1^2), same A1, A2 as chi_p. + # chi_p keeps the larger of the two terms; this keeps their vector sum, so it depends on phi12 + # (in-plane spins that cancel give a small value). L frame. + q = self.m2/self.m1 + A1 = (2+ 3.*q/2); A2 = (2+3./(2*q)) + return np.abs( (self.s1x + 1j*self.s1y) + (A2/A1)*q**2*(self.s2x + 1j*self.s2y) ) if p == 'chi_pavg': if (abs(self.s1x) < 1e-4 and abs(self.s1y) < 1e-4 and abs(self.s2x) < 1e-4 and abs(self.s2y) < 1e-4): chipavg = 0.0 @@ -6132,6 +6143,28 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_ x_out[:,indx_pout_s2y] = x_in[:,indx_chi2]*sintheta2*sinphi2 coord_names_reduced.remove('s2x') coord_names_reduced.remove('s2y') + # in-plane magnitudes, relative azimuth, and ring coordinates, vectorized. L frame only: + # for any other spin_convention these fall through to extract_param, as before. + ring_names = ['chi1_perp', 'chi2_perp', 'phi12', 'SOverM2_perp', 'DeltaOverM2_perp', 'chi_p_vec'] + if spin_convention == "L" and any(p in coord_names_reduced for p in ring_names): + indx_phi1 = low_level_coord_names.index('phi1') + indx_phi2 = low_level_coord_names.index('phi2') + chi1_perp = x_in[:,indx_chi1]*np.sqrt(1-x_in[:,indx_ct1]**2) + chi2_perp = x_in[:,indx_chi2]*np.sqrt(1-x_in[:,indx_ct2]**2) + v1 = chi1_perp*np.exp(1j*x_in[:,indx_phi1]) + v2 = chi2_perp*np.exp(1j*x_in[:,indx_phi2]) + mtot_vals = m1_vals + m2_vals + q_vals = m2_vals/m1_vals + A1 = 2 + 1.5*q_vals; A2 = 2 + 1.5/q_vals + ring_vals = {'chi1_perp': chi1_perp, 'chi2_perp': chi2_perp, + 'phi12': np.mod(x_in[:,indx_phi2] - x_in[:,indx_phi1], 2*np.pi), + 'SOverM2_perp': np.abs(v1*m1_vals**2 + v2*m2_vals**2)/mtot_vals**2, + 'DeltaOverM2_perp': np.abs(v1*m1_vals - v2*m2_vals)/mtot_vals, + 'chi_p_vec': np.abs(v1 + (A2/A1)*q_vals**2*v2)} + for p in ring_names: + if p in coord_names_reduced: + x_out[:,coord_names.index(p)] = ring_vals[p] + coord_names_reduced.remove(p) # Spin pseudo-cylindrical coordinate names, standard framing if ('s1z_bar' in low_level_coord_names) and ('phi1' in low_level_coord_names) and ('s2z_bar' in low_level_coord_names) and ('phi2' in low_level_coord_names) and ('mc' in low_level_coord_names) and ('eta' in low_level_coord_names or 'delta_mc' in low_level_coord_names): diff --git a/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py b/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py new file mode 100644 index 000000000..9580302b7 --- /dev/null +++ b/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py @@ -0,0 +1,60 @@ +"""phi12 and the ring coordinates (chi1_perp, chi2_perp, SOverM2_perp, DeltaOverM2_perp, chi_p_vec). + +The vectorized spherical-spin path of convert_waveform_coordinates must agree with +ChooseWaveformParams.extract_param, and must not fall through to the per-row loop. +""" + +import numpy as np +import lal + +from RIFT import lalsimutils + +RING = ['chi1_perp', 'chi2_perp', 'phi12', 'SOverM2_perp', 'DeltaOverM2_perp', 'chi_p_vec'] +LOW = ['mc', 'delta_mc', 'chi1', 'cos_theta1', 'phi1', 'chi2', 'cos_theta2', 'phi2'] + + +def _draws(n=300, seed=4): + rng = np.random.default_rng(seed) + return np.column_stack([rng.uniform(5, 30, n), rng.uniform(0.01, 0.8, n), rng.uniform(0, 0.99, n), + rng.uniform(-1, 1, n), rng.uniform(0, 2 * np.pi, n), rng.uniform(0, 0.99, n), + rng.uniform(-1, 1, n), rng.uniform(0, 2 * np.pi, n)]) + + +def _params(row): + mc, dmc, c1, ct1, p1, c2, ct2, p2 = row + m1, m2 = lalsimutils.m1m2(mc, 0.25 * (1 - dmc ** 2)) + P = lalsimutils.ChooseWaveformParams() + P.m1, P.m2 = m1 * lal.MSUN_SI, m2 * lal.MSUN_SI + s1 = c1 * np.sqrt(1 - ct1 ** 2) + s2 = c2 * np.sqrt(1 - ct2 ** 2) + P.s1x, P.s1y, P.s1z = s1 * np.cos(p1), s1 * np.sin(p1), c1 * ct1 + P.s2x, P.s2y, P.s2z = s2 * np.cos(p2), s2 * np.sin(p2), c2 * ct2 + return P + + +def test_vectorized_matches_extract_param(capsys): + x = _draws() + y = lalsimutils.convert_waveform_coordinates(x, coord_names=RING, low_level_coord_names=LOW) + assert "Fallthrough" not in capsys.readouterr().out + for i, row in enumerate(x): + P = _params(row) + for j, name in enumerate(RING): + ref = P.extract_param(name) + d = abs(y[i, j] - ref) + if name == 'phi12': + d = min(d, 2 * np.pi - d) + assert d < 1e-9, (name, y[i, j], ref) + + +def test_chi_p_vec_limits(): + # one in-plane spin: chi_p_vec equals chi1_perp. Antiparallel in-plane spins of the weighted size cancel. + P = lalsimutils.ChooseWaveformParams() + P.m1, P.m2 = 10 * lal.MSUN_SI, 5 * lal.MSUN_SI + P.s1x, P.s1y, P.s2x, P.s2y = 0.4, 0.0, 0.0, 0.0 + assert abs(P.extract_param('chi_p_vec') - 0.4) < 1e-12 + q = 0.5 + A1 = 2 + 1.5 * q + A2 = 2 + 1.5 / q + P.s2x = -0.4 / ((A2 / A1) * q ** 2) + assert P.extract_param('chi_p_vec') < 1e-12 + assert abs(P.extract_param('phi12') - np.pi) < 1e-12 From 9e551386bb8123da593fdffeb38d2513789d1869 Mon Sep 17 00:00:00 2001 From: Richard O'Shaughnessy Date: Wed, 30 Sep 2026 15:59:23 -0700 Subject: [PATCH 2/4] lalsimutils: review fixes for the ring coordinates (object arrays, Kerr, phi12) From adversarial review of this PR and of the O4c version (#203 upstream): - CIP's default sampler passes an object array of python floats; the ring block's ufuncs raised TypeError, so CIP exited 1 with these fit coordinates. The block now casts to float. - The ring block can end convert_waveform_coordinates before the per-row fallthrough, whose enforce_kerr rule then never ran. The block applies it. - phi12 is 0 in both paths when an in-plane spin vanishes (they disagreed), cannot return exactly 2 pi, and joins periodic_params. - chi_p_vec leaves valid_params: grid readers assign_param every listed column. - Tests: phi12 at pi/3 in both paths, object arrays, zero in-plane spin, Kerr, valid_params. Core-units floors re-measured: 605/592. Co-Authored-By: Claude Opus 5.5 --- .travis/test-core-units.sh | 6 ++- .../Code/RIFT/lalsimutils.py | 40 ++++++++++++------ .../Code/test/test_ring_coordinates.py | 41 +++++++++++++++++++ 3 files changed, 72 insertions(+), 15 deletions(-) diff --git a/.travis/test-core-units.sh b/.travis/test-core-units.sh index 956b3294d..5611b97a0 100755 --- a/.travis/test-core-units.sh +++ b/.travis/test-core-units.sh @@ -326,7 +326,9 @@ done # `import cupy` FAILS there, numpy backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0: # per-file 600, junit 603 collected / 590 passed / 13 skipped / 0 failed, of which 3 are # subtests. So the plugin-free floors are 600/587; the old 592/579 sat 6 below the tree. -EXPECTED_TESTS=600 +# Review of #377 added five more ring-coordinate tests, no skips. RE-MEASURED the same way: +# junit 608 collected / 595 passed / 13 skipped / 0 failed (3 subtests), so the floors are 605/592. +EXPECTED_TESTS=605 # 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 13 skips are # environment legs -- cupy in test_seeding_reproducibility, device legs in @@ -335,7 +337,7 @@ EXPECTED_TESTS=600 # (mcsamplerNFlow is an optional dependency and is absent from the IGWN environment), and # the xfail in test_uv_symmetry. test_eos_portfolio_sampler.py adds 12 tests and # test_cip_portfolio_members.py 4, none of them skips. -EXPECTED_PASSED=587 +EXPECTED_PASSED=592 MAX_SKIPPED=13 # The floors must be INTEGERS, and this is checked rather than assumed. `[ 347 -lt FOO ]` does diff --git a/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py b/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py index 7eaf427a6..e6a0790b8 100644 --- a/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py +++ b/MonteCarloMarginalizeCode/Code/RIFT/lalsimutils.py @@ -414,11 +414,11 @@ def lsu_StringFromPNOrder(order): # Class to hold arguments of ChooseWaveform functions # -valid_params = ['m1', 'm2', 's1x', 's1y', 's1z', 's2x', 's2y', 's2z', 'chi1_perp', 'chi2_perp', 'chi1_perp_bar', 'chi2_perp_bar','chi1_perp_u', 'chi2_perp_u', 's1z_bar', 's2z_bar', 'lambda1', 'lambda2', 'theta','phi', 'phiref', 'psi', 'incl', 'tref', 'dist', 'mc', 'mc_ecc', 'eta', 'delta_mc', 'chi1', 'chi2', 'thetaJN', 'phiJL', 'theta1', 'theta2', 'cos_theta1', 'cos_theta2', 'theta1_Jfix', 'theta2_Jfix', 'psiJ', 'beta', 'cos_beta', 'sin_phiJL', 'cos_phiJL', 'phi12', 'phi1', 'phi2', 'LambdaTilde', 'DeltaLambdaTilde', 'lambda_plus', 'lambda_minus', 'q', 'mtot','xi','chiz_plus', 'chiz_minus', 'chieff_aligned','fmin','fref', "SOverM2_perp", "SOverM2_L", "DeltaOverM2_perp", "DeltaOverM2_L", "shu","ampO", "phaseO",'eccentricity','eccentricity_squared','eccentricity_ln', 'chi_pavg','chi_p_vec','mu1','mu2','eos_table_index','meanPerAno','a6c','E0','p_phi0','hypclass'] +valid_params = ['m1', 'm2', 's1x', 's1y', 's1z', 's2x', 's2y', 's2z', 'chi1_perp', 'chi2_perp', 'chi1_perp_bar', 'chi2_perp_bar','chi1_perp_u', 'chi2_perp_u', 's1z_bar', 's2z_bar', 'lambda1', 'lambda2', 'theta','phi', 'phiref', 'psi', 'incl', 'tref', 'dist', 'mc', 'mc_ecc', 'eta', 'delta_mc', 'chi1', 'chi2', 'thetaJN', 'phiJL', 'theta1', 'theta2', 'cos_theta1', 'cos_theta2', 'theta1_Jfix', 'theta2_Jfix', 'psiJ', 'beta', 'cos_beta', 'sin_phiJL', 'cos_phiJL', 'phi12', 'phi1', 'phi2', 'LambdaTilde', 'DeltaLambdaTilde', 'lambda_plus', 'lambda_minus', 'q', 'mtot','xi','chiz_plus', 'chiz_minus', 'chieff_aligned','fmin','fref', "SOverM2_perp", "SOverM2_L", "DeltaOverM2_perp", "DeltaOverM2_L", "shu","ampO", "phaseO",'eccentricity','eccentricity_squared','eccentricity_ln', 'chi_pavg','mu1','mu2','eos_table_index','meanPerAno','a6c','E0','p_phi0','hypclass'] # so far, used for puffball, to prevent insanity (infinite growth) and/or death to downselect # - note we also provide for extrinsic: RA (phi), phiref, psi, just in case we need it in the future -periodic_params = {'phi1':2*np.pi, 'phi2':2*np.pi, 'phiref':2*np.pi, 'psi':np.pi, 'meanPerAno':2*np.pi, 'phi':2*np.pi, 'phiJL':2*np.pi, 'psiJ':2*np.pi} +periodic_params = {'phi1':2*np.pi, 'phi2':2*np.pi, 'phi12':2*np.pi, 'phiref':2*np.pi, 'psi':np.pi, 'meanPerAno':2*np.pi, 'phi':2*np.pi, 'phiJL':2*np.pi, 'psiJ':2*np.pi} tex_dictionary = { "mtot": r'$M$', @@ -1392,8 +1392,12 @@ def extract_param(self,p): Sp = np.max([np.linalg.norm( A1*S1p), np.linalg.norm(A2*S2p)]) return Sp/(A1*m1**2) # divide by term for *larger* BH if p == 'phi12': - # azimuth of spin 2's in-plane component relative to spin 1's, in [0, 2 pi), L frame - return np.mod(np.arctan2(self.s2y, self.s2x) - np.arctan2(self.s1y, self.s1x), 2*np.pi) + # azimuth of spin 2's in-plane component relative to spin 1's, in [0, 2 pi), L frame. + # Undefined if either in-plane component vanishes; 0 is returned then (same as the vectorized path). + if np.hypot(self.s1x, self.s1y) == 0 or np.hypot(self.s2x, self.s2y) == 0: + return 0. + val = np.mod(np.arctan2(self.s2y, self.s2x) - np.arctan2(self.s1y, self.s1x), 2*np.pi) + return 0. if val >= 2*np.pi else val # np.mod of a tiny negative number rounds to 2 pi if p == 'chi_p_vec': # vector-sum (ring) analogue of chi_p: |A1 S1perp + A2 S2perp| / (A1 m1^2), same A1, A2 as chi_p. # chi_p keeps the larger of the two terms; this keeps their vector sum, so it depends on phi12 @@ -5930,6 +5934,7 @@ 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) ) ) + kerr_violation_ring = None # set by the vectorized in-plane block, which can end the conversion early # 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: @@ -6147,24 +6152,31 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_ # for any other spin_convention these fall through to extract_param, as before. ring_names = ['chi1_perp', 'chi2_perp', 'phi12', 'SOverM2_perp', 'DeltaOverM2_perp', 'chi_p_vec'] if spin_convention == "L" and any(p in coord_names_reduced for p in ring_names): + # CIP's default sampler passes an object array of python floats; ufuncs need a float array + xf = np.asarray(x_in, dtype=float) + m1f = np.asarray(m1_vals, dtype=float); m2f = np.asarray(m2_vals, dtype=float) indx_phi1 = low_level_coord_names.index('phi1') indx_phi2 = low_level_coord_names.index('phi2') - chi1_perp = x_in[:,indx_chi1]*np.sqrt(1-x_in[:,indx_ct1]**2) - chi2_perp = x_in[:,indx_chi2]*np.sqrt(1-x_in[:,indx_ct2]**2) - v1 = chi1_perp*np.exp(1j*x_in[:,indx_phi1]) - v2 = chi2_perp*np.exp(1j*x_in[:,indx_phi2]) - mtot_vals = m1_vals + m2_vals - q_vals = m2_vals/m1_vals + chi1_perp = xf[:,indx_chi1]*np.sqrt(1-xf[:,indx_ct1]**2) + chi2_perp = xf[:,indx_chi2]*np.sqrt(1-xf[:,indx_ct2]**2) + v1 = chi1_perp*np.exp(1j*xf[:,indx_phi1]) + v2 = chi2_perp*np.exp(1j*xf[:,indx_phi2]) + mtot_vals = m1f + m2f + q_vals = m2f/m1f A1 = 2 + 1.5*q_vals; A2 = 2 + 1.5/q_vals + phi12 = np.where((chi1_perp > 0) & (chi2_perp > 0), np.mod(xf[:,indx_phi2] - xf[:,indx_phi1], 2*np.pi), 0.) ring_vals = {'chi1_perp': chi1_perp, 'chi2_perp': chi2_perp, - 'phi12': np.mod(x_in[:,indx_phi2] - x_in[:,indx_phi1], 2*np.pi), - 'SOverM2_perp': np.abs(v1*m1_vals**2 + v2*m2_vals**2)/mtot_vals**2, - 'DeltaOverM2_perp': np.abs(v1*m1_vals - v2*m2_vals)/mtot_vals, + 'phi12': np.where(phi12 >= 2*np.pi, 0., phi12), + 'SOverM2_perp': np.abs(v1*m1f**2 + v2*m2f**2)/mtot_vals**2, + 'DeltaOverM2_perp': np.abs(v1*m1f - v2*m2f)/mtot_vals, 'chi_p_vec': np.abs(v1 + (A2/A1)*q_vals**2*v2)} for p in ring_names: if p in coord_names_reduced: x_out[:,coord_names.index(p)] = ring_vals[p] coord_names_reduced.remove(p) + if enforce_kerr: + # same rule as the per-row fallthrough below, which this block can bypass + kerr_violation_ring = (xf[:,indx_chi1] > 1) | (xf[:,indx_chi2] > 1) # Spin pseudo-cylindrical coordinate names, standard framing if ('s1z_bar' in low_level_coord_names) and ('phi1' in low_level_coord_names) and ('s2z_bar' in low_level_coord_names) and ('phi2' in low_level_coord_names) and ('mc' in low_level_coord_names) and ('eta' in low_level_coord_names or 'delta_mc' in low_level_coord_names): @@ -6382,6 +6394,8 @@ def convert_waveform_coordinates(x_in,coord_names=['mc', 'eta'],low_level_coord_ 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 kerr_violation_ring is not None: + x_out[kerr_violation_ring] = -np.inf if len(coord_names_reduced)<1: return x_out diff --git a/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py b/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py index 9580302b7..c607c57e2 100644 --- a/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py +++ b/MonteCarloMarginalizeCode/Code/test/test_ring_coordinates.py @@ -58,3 +58,44 @@ def test_chi_p_vec_limits(): P.s2x = -0.4 / ((A2 / A1) * q ** 2) assert P.extract_param('chi_p_vec') < 1e-12 assert abs(P.extract_param('phi12') - np.pi) < 1e-12 + + +def _one(*row): + return np.array([row], dtype=float) + + +def test_phi12_sign_both_paths(): + # spin 2 sixty degrees ahead of spin 1: phi12 = pi/3, not 5 pi/3 (a sign flip in both paths fails here) + row = (10., 0.3, 0.5, 0., 0.2, 0.5, 0., 0.2 + np.pi / 3) + y = lalsimutils.convert_waveform_coordinates(_one(*row), coord_names=['phi12'], low_level_coord_names=LOW) + assert abs(y[0, 0] - np.pi / 3) < 1e-12 + assert abs(_params(row).extract_param('phi12') - np.pi / 3) < 1e-12 + + +def test_object_array_input(): + # CIP's default sampler (adaptive_cartesian) passes an object array of python floats + x = _draws(20) + y = lalsimutils.convert_waveform_coordinates(x.astype(object), coord_names=RING, low_level_coord_names=LOW) + y_ref = lalsimutils.convert_waveform_coordinates(x, coord_names=RING, low_level_coord_names=LOW) + assert np.allclose(np.asarray(y, dtype=float), y_ref, atol=1e-12) + + +def test_phi12_zero_inplane_spin(): + # phi12 is undefined without an in-plane component; both paths return 0 + for row in [(10., 0.3, 0.0, 0.5, 1.5, 0.5, 0.2, 2.5), (10., 0.3, 0.5, 0.2, 1.5, 0.5, -1.0, 2.5)]: + y = lalsimutils.convert_waveform_coordinates(_one(*row), coord_names=['phi12'], low_level_coord_names=LOW) + assert y[0, 0] == 0. + assert _params(row).extract_param('phi12') == 0. + + +def test_enforce_kerr_in_ring_block(): + # the vectorized block can end the conversion early; it applies the fallthrough's Kerr rule itself + x = np.vstack([_one(10., 0.3, 1.2, 0.2, 1.5, 0.5, 0.2, 2.5), _one(10., 0.3, 0.5, 0.2, 1.5, 0.5, 0.2, 2.5)]) + y = lalsimutils.convert_waveform_coordinates(x, coord_names=['mc', 'delta_mc', 'chi1_perp', 'chi_p_vec'], + low_level_coord_names=LOW, enforce_kerr=True) + assert np.all(y[0] == -np.inf) and np.all(np.isfinite(y[1])) + + +def test_chi_p_vec_not_assignable(): + # grid readers call assign_param on every column named in valid_params; chi_p_vec is derived only + assert 'chi_p_vec' not in lalsimutils.valid_params From 54eadb1cce526cbac02df367502faf941079a6e6 Mon Sep 17 00:00:00 2001 From: Session Router Gate Date: Thu, 1 Oct 2026 00:33:02 +0000 Subject: [PATCH 3/4] Address automated review findings for PR #377 --- .travis/test-core-units.sh | 34 +++++++++++++++++++--------------- 1 file changed, 19 insertions(+), 15 deletions(-) diff --git a/.travis/test-core-units.sh b/.travis/test-core-units.sh index d2e407ada..25324bfe7 100755 --- a/.travis/test-core-units.sh +++ b/.travis/test-core-units.sh @@ -315,12 +315,17 @@ done # RIFT_COREUNIT_PYTHON pointed at the IGWN interpreter: per-file 600 over 52 files, # junit 603 collected / 590 passed / 13 skipped / 0 failed. # -# 613/600 test_ring_coordinates.py added (7 tests: vectorized in-plane spin and ring -# coordinates, phi12, chi_p_vec, object-array input, Kerr rule), merged with -# test_complex_overlap_interpolate_max.py above. MEASURED on CIT (ldas-grid; numpy -# backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0, with RIFT_COREUNIT_PYTHON -# pointed at the IGWN interpreter: junit 616 collected / 603 passed / 13 skipped / -# 0 failed, of which 3 are subtests. +# 607/594 test_ring_coordinates.py added (7 tests: vectorized in-plane spin and ring +# coordinates, phi12, chi_p_vec, object-array input, Kerr rule). The 600/587 row +# above is the base manifest INCLUDING test_complex_overlap_interpolate_max.py, +# measured over all 52 of its files; this branch adds exactly ONE file on top of +# it, and its seven entry points are plain test_* functions -- no parametrization, +# no classes, no skips. So 607/594. +# A merge resolution first wrote 613/600, carrying #375's eight tests a second +# time on top of a 600/587 that already contained them -- six above the tree. A +# floor set too high fails a run in which every test passes, which is not the +# failure this gate is for; recount against the base row, never re-add a file the +# base already counts. # # 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. @@ -336,14 +341,13 @@ 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. # The detector-network sky mapping adds two passing tests and no skips. -# test_ring_coordinates.py adds two more, no skips. RE-MEASURED with it on CIT (ldas-grid; -# `import cupy` FAILS there, numpy backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0: -# per-file 600, junit 603 collected / 590 passed / 13 skipped / 0 failed, of which 3 are -# subtests. So the plugin-free floors are 600/587; the old 592/579 sat 6 below the tree. -# Review of #377 added five more ring-coordinate tests, no skips (605/592 measured). -# test_complex_overlap_interpolate_max.py (#375) adds 8 passing tests and no skips. -# Merged with #375: 613/600 (see the table above). -EXPECTED_TESTS=613 +# test_complex_overlap_interpolate_max.py (#375) is part of the BASE manifest and is already +# inside the 600/587 row of the table above (per-file 600 over 52 files, junit 603 collected / +# 590 passed / 13 skipped / 0 failed, of which 3 are subtests). This branch adds one file, +# test_ring_coordinates.py: seven plain test_* functions, no parametrization, no skips. So the +# plugin-free floors are 600+7 / 587+7 = 607/594, and MAX_SKIPPED stays 13. #375's eight tests +# are NOT added again here -- doing that is what put the floors six above the tree. +EXPECTED_TESTS=607 # 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 13 skips are # environment legs -- cupy in test_seeding_reproducibility, device legs in @@ -352,7 +356,7 @@ EXPECTED_TESTS=613 # (mcsamplerNFlow is an optional dependency and is absent from the IGWN environment), and # the xfail in test_uv_symmetry. test_eos_portfolio_sampler.py adds 12 tests and # test_cip_portfolio_members.py 4, none of them skips. -EXPECTED_PASSED=600 +EXPECTED_PASSED=594 MAX_SKIPPED=13 # The floors must be INTEGERS, and this is checked rather than assumed. `[ 347 -lt FOO ]` does From aae57ce2576555a8ea78ff09d2774d4bb65ece12 Mon Sep 17 00:00:00 2001 From: Richard O'Shaughnessy Date: Thu, 1 Oct 2026 13:56:13 -0700 Subject: [PATCH 4/4] Restore core-unit floors 613/600 (revert 54eadb1c) 54eadb1c lowered the floors to 607/594 on the premise that the merge counted #375's tests twice. CI on that very commit counted 616 collected / 603 passed / 13 skipped (3 subtests), i.e. 613/600, so 607/594 sat six below the tree. Co-Authored-By: Claude Opus 5.5 --- .travis/test-core-units.sh | 34 +++++++++++++++------------------- 1 file changed, 15 insertions(+), 19 deletions(-) diff --git a/.travis/test-core-units.sh b/.travis/test-core-units.sh index 25324bfe7..d2e407ada 100755 --- a/.travis/test-core-units.sh +++ b/.travis/test-core-units.sh @@ -315,17 +315,12 @@ done # RIFT_COREUNIT_PYTHON pointed at the IGWN interpreter: per-file 600 over 52 files, # junit 603 collected / 590 passed / 13 skipped / 0 failed. # -# 607/594 test_ring_coordinates.py added (7 tests: vectorized in-plane spin and ring -# coordinates, phi12, chi_p_vec, object-array input, Kerr rule). The 600/587 row -# above is the base manifest INCLUDING test_complex_overlap_interpolate_max.py, -# measured over all 52 of its files; this branch adds exactly ONE file on top of -# it, and its seven entry points are plain test_* functions -- no parametrization, -# no classes, no skips. So 607/594. -# A merge resolution first wrote 613/600, carrying #375's eight tests a second -# time on top of a 600/587 that already contained them -- six above the tree. A -# floor set too high fails a run in which every test passes, which is not the -# failure this gate is for; recount against the base row, never re-add a file the -# base already counts. +# 613/600 test_ring_coordinates.py added (7 tests: vectorized in-plane spin and ring +# coordinates, phi12, chi_p_vec, object-array input, Kerr rule), merged with +# test_complex_overlap_interpolate_max.py above. MEASURED on CIT (ldas-grid; numpy +# backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0, with RIFT_COREUNIT_PYTHON +# pointed at the IGWN interpreter: junit 616 collected / 603 passed / 13 skipped / +# 0 failed, of which 3 are subtests. # # 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. @@ -341,13 +336,14 @@ 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. # The detector-network sky mapping adds two passing tests and no skips. -# test_complex_overlap_interpolate_max.py (#375) is part of the BASE manifest and is already -# inside the 600/587 row of the table above (per-file 600 over 52 files, junit 603 collected / -# 590 passed / 13 skipped / 0 failed, of which 3 are subtests). This branch adds one file, -# test_ring_coordinates.py: seven plain test_* functions, no parametrization, no skips. So the -# plugin-free floors are 600+7 / 587+7 = 607/594, and MAX_SKIPPED stays 13. #375's eight tests -# are NOT added again here -- doing that is what put the floors six above the tree. -EXPECTED_TESTS=607 +# test_ring_coordinates.py adds two more, no skips. RE-MEASURED with it on CIT (ldas-grid; +# `import cupy` FAILS there, numpy backend) 2026-09-30, IGWN conda python 3.11 / lal 7.7.0: +# per-file 600, junit 603 collected / 590 passed / 13 skipped / 0 failed, of which 3 are +# subtests. So the plugin-free floors are 600/587; the old 592/579 sat 6 below the tree. +# Review of #377 added five more ring-coordinate tests, no skips (605/592 measured). +# test_complex_overlap_interpolate_max.py (#375) adds 8 passing tests and no skips. +# Merged with #375: 613/600 (see the table above). +EXPECTED_TESTS=613 # 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 13 skips are # environment legs -- cupy in test_seeding_reproducibility, device legs in @@ -356,7 +352,7 @@ EXPECTED_TESTS=607 # (mcsamplerNFlow is an optional dependency and is absent from the IGWN environment), and # the xfail in test_uv_symmetry. test_eos_portfolio_sampler.py adds 12 tests and # test_cip_portfolio_members.py 4, none of them skips. -EXPECTED_PASSED=594 +EXPECTED_PASSED=600 MAX_SKIPPED=13 # The floors must be INTEGERS, and this is checked rather than assumed. `[ 347 -lt FOO ]` does