The README describes the final method and gives the benchmark result. These notes cover the decisions that shaped the implementation and why I kept the direct-projection route.
I began with the obvious pipeline: read a CDM, rotate the RTN covariance into inertial coordinates, move the states to linear closest approach, build the encounter plane, and integrate the Gaussian over the hard-body circle. I expected most of the work to be in the integral.
The integral mattered, but the more useful work was establishing which operations were necessary, which were redundant, and which agreement claims the code could actually support.
Why direct projection stayed
The observation that unlocked the final geometry was not a large error. It was no change at all.
I had written a linear closest-approach step because CARA’s workbook distinguishes an adjusted Pc2D result. The function produced a nonzero time offset and visibly different positions. Pc stayed the same.
My first suspicion was that the shift was incomplete. I checked whether I had moved the wrong object, used the wrong sign, or failed to update something downstream. Then I printed the two quantities that actually enter the Foster calculation: the encounter-plane basis and the projected miss vector. They were unchanged before and after the shift.
That forced me to stop treating CARA’s processing order as the definition of the geometry. The linear correction adds a multiple of relative velocity to relative position. The encounter-plane projection removes that direction. The cross product that orients the plane is unchanged for the same reason.
The full derivation is in the README. The direct route was not chosen first and justified later. It came from trying to explain why an apparently meaningful preprocessing step had no effect on the quantity I was calculating.
I removed the shift from the production calculation. I kept the helper and added a test that calculates the encounter both ways. If those routes ever stop agreeing under the assumptions used here, the test fails. The resulting implementation is smaller than a procedural translation of CARA, but it remains tied to the reference through an explicit geometric equivalence.
Making the tail calculation hold up
Once ordinary cases were close, the smallest CARA cases were still troublesome. Two probabilities are around 1e-168. The straightforward implementation could lose the inner interval probability and return zero.
The tempting response was to keep tightening the numerical quadrature. That would not address the failure. The loss occurred inside the integrand when subtracting two nearly equal normal CDF values. Asking the outer integrator for more samples cannot recover digits that have already disappeared.
I separated that operation from the rest of the integral and handled same-tail intervals in log-CDF space. That keeps both CARA tail cases nonzero without introducing a separate high-precision runtime path. It also preserves the same one-dimensional integration structure across ordinary and extreme cases.
This also changed how I used the corpus. A tiny absolute error at 1e-168 is useful as a numerical stress case. It is not useful as an enormous performance ratio, and I did not want the repository to present it that way.
Turning a comparison into a reproducible benchmark
The numerical table came together before the benchmark command fully represented the claim I wanted to make.
The first benchmark runner took HBR directly from the workbook. The CDM parser also read HBR, but its result never reached the calculation. That meant the headline command could pass even if the parser extracted the wrong radius. I changed the runner to require the CDM value, compare it with the workbook, and then use the parsed value.
The command also printed a result and returned zero regardless of whether the acceptance tolerance was met. That is convenient for an exploratory script and wrong for a validation command used by CI. It now fails its process exit status when any case misses 1e-12.
The external reference had the same problem at a larger scale. Cloning the current CARA repository does not reproduce a fixed experiment. I pinned the commit, workbook digest, case count, exact CDM membership, and cleanliness of the benchmark directory. I added those checks after the calculation was working, but they are part of the claim now. Without them, 53/53 has no stable referent.
None of those changes altered the Pc algorithm. They turned the result into an executable claim: the public command exercises the CDM path, checks a fixed external reference, and fails when the stated tolerance is not met.
A covariance rule I could defend
The CDM test cases hid the first covariance problem because their lower-triangular entries are mirrored into an exactly symmetric matrix. A covariance supplied through the other public path may already have passed through several matrix transformations. It can be symmetric mathematically and still carry a small floating-point skew.
My first symmetry check used a fixed absolute tolerance. It rejected legitimate large covariances after an ordinary frame transformation.
I replaced that with a relative test against the scale of the whole matrix. That fixed the false rejection and introduced a different mistake. One very large variance could make enough room for a material mismatch in an unrelated off-diagonal pair. A concrete matrix with C00 = 1e12, C01 = 0.5, and C10 = 0.0 was accepted and averaged into an apparently clean covariance.
The current check works pair by pair. Each asymmetry is judged against the entries involved and the corresponding diagonal variances. Small transformation roundoff is symmetrized; an unrelated large value elsewhere cannot hide a bad pair.
Even the regression test needed correction. I first used a seeded chain of five orthogonal transformations and asserted that it would exceed the old tolerance. It did on one NumPy and BLAS combination and not on another. The implementation decision was still sound, but the test depended on a platform-specific rounding event. I replaced it with a deterministic covariance placed deliberately between the former and current limits.
That sequence is a better account of the covariance code than simply calling it “scale aware.” The current rule accepts transformation-level roundoff, rejects locally material disagreement, and has a deterministic boundary test that behaves the same way across numerical environments.
Keeping CDM meaning at the parser boundary
The initial CDM parser used OBJECT1 and OBJECT2 as delimiters, then discarded the labels. Two sections both named OBJECT1 therefore became two valid-looking objects. Every array had the expected shape, so the error could continue all the way to a plausible Pc.
The fix was small: preserve the section sequence and require OBJECT1 followed by OBJECT2. Dimensional validation would never catch this class of error. The parser now enforces the meaning while the label is still present.
HBR required the opposite change. I originally made it mandatory when constructing a conjunction because Pc needs a collision radius. Real CDMs do not necessarily carry one. Parsing the state should still succeed; calculating Pc should not proceed until a radius is supplied. Splitting those two decisions made the reader less dependent on the CARA fixture format without inventing a default. Together, the section-order and HBR rules give the reader a useful boundary: preserve valid source information, reject semantic contradictions, and do not manufacture missing physical inputs.
Covering the other end of the probability interval
The 53 cases exercise extremely small probabilities, but they did not expose the other end of the interval. A centered isotropic case with a radius of 8.7 sigma returned 1.0000000000000002 from quadrature. My guard raised because the value was greater than one.
That guard looked rigorous and was wrong by one floating-point step. I added a narrow tolerance at the probability boundary, clip only inside that tolerance, and continue to raise on a material violation. The regression uses the centered isotropic closed form rather than another numerical integration.
This was the final implementation change before preparing the repository. It did not alter any CARA result. It closed a public-API case that the benchmark did not contain and left the same calculation well behaved near both zero and one.
What the finished component establishes
If I started again, I would define the benchmark identity before writing the runner: source commit, corpus membership, workbook digest, input path, comparison column, and failure exit. I would also separate corpus reproduction tests from public-API boundary tests from the beginning. They answer different questions.
I would still keep the direct encounter-plane formulation. It survived the part of the work that mattered: I tried to disprove the equivalence, traced why the explicit shift vanished, and left a test that compares the two routes.
The finished component is not a wrapper around CARA and it is not a line-by-line translation. It independently carries CDM state and covariance through frame conversion, encounter-plane construction, and stable circular Gaussian integration, then reproduces all 53 pinned CARA cases within the published tolerance. The tests also cover geometry invariance, probabilities deep in the Gaussian tail, the near-certainty boundary, covariance roundoff, and CDM semantic errors that can otherwise survive as plausible numbers.
That is a meaningful amount of orbital-mechanics and numerical work in a small repository. The path was not obvious, but the result is now compact, inspectable, and reproducible.