Stability controls for the bfcorr BField domain walk - #224
Conversation
…cosmics Robustness fixes to the bfcorr=true CentralHelix domain machinery for near-vertical cosmic tracks that exit the DS bore (|B| -> 0): - Config: add mindtstep_ (min domain step), domainmargin_ (clamp domains to the active range +/- margin), and minfield_ (stop extrapolation before the dPardB pole). Defaults preserve legacy behaviour (margin = max(), minfield = 0). - createDomains/extendDomains: floor the domain step at mindtstep_ so a vanishing rangeInTolerance in a high-gradient/low-momentum region can no longer spawn an unbounded number of micro-domains (CPU hang / OOM). Clamp domain bounds to the active range +/- domainmargin_ so a uniform-field track just past the tracker doesn't sample Bz->0 at a domain midpoint and blow the CentralHelix reference off by metres via the singular bfrac/(1+bfrac) correction. - createEffects: connect only domains that actually abut, instead of throwing "Invalid domains" across a gap between separate contiguous domain blocks. - extrapolate: break before |B| < minfield_, tested at the new domain midpoint (not the frontier), handing the field-free region to a straight line tail. - convertSeed: guard against an empty fittraj_ when a degenerate active range yields zero domains -- createEffects would otherwise build a Measurement against an empty PiecewiseTrajectory and throw std::length_error, aborting the whole art event. Treat it as an unfittable track (outsidemap) so the module drops it cleanly. Includes temporary domain-count / empty-seed diagnostics (marked "remove before PR"). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
…gnostics - Config: new maxdomains_ (default unlimited, preserving legacy behaviour) printed in operator<<. A hard cap on the number of BField domains a single fit may accumulate. - Track.hh createDomains/extendDomains: when domains exceed maxdomains_, throw -- caught by the existing fit()/processEnds() handlers, recording the fit as failed so the unusable track is dropped. A diverging low-momentum track can otherwise build ~1e5 domains (each a KKDW effect + traj piece) -> wasted CPU + ~GB memory before being dropped anyway. Validated over 112 cosmic files: usable tracks peak at 368 domains, runaways at 1e4-1e5; a cap of 1000 (set in the reco config) cleanly separates them. - Remove the temporary domain-count / empty-seed cout diagnostics added during the investigation (the empty-seed outsidemap guard itself is retained). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
brownd1978
left a comment
brownd1978
left a comment
There was a problem hiding this comment.
Thanks for this PR. You've clearly understood a lot about how the fit works and how it can fail, and made important improvements.
These fixes make CreateDomains quite long, which suggests it should be refactored. I see big blocks that look very similar: maybe those need to be encapsulated in new functions or lambdas. I wonder if some protection could be moved to BFieldMap. Perhaps BFieldMap could define 'boundary' and '0-field' regions on construction, so that those wouldn't have to be test for by each track. It might help to add a createDomain function that returns a guaranteed usable domain, which could deal with things like minimum field.
Protections connected with sign flipping could maybe be moved to CentralHelix; that should allow omega to change sign while keeping the same particle charge, at least within the same piecetraj.
I suspect the config can be simplified; perhaps a min and max domain range are sufficient. Limits on physical ranges (time) have clearer interpretation than limits on container size.
I believe the minimum field value could go into the BFieldMap constructor: LH and KL won't ever care, so it's fine to set it to a single, global value.
I'm also curious if unusable inputs could be identified and failed before any field calculations, to avoid burdening downstream code with testing.
It's always good to catch failure conditions, but some of the throw/catch blocks seem designed to handle common/routine conditions, which might be clearer to implement as normal functions.
The term 'legacy' in comments won't mean much in a few years: referring to conditions that make the block relevant will be more useful.
I'm happy to discuss if any points above aren't clear.
| unsigned maxdomains_ = std::numeric_limits<unsigned>::max(); // hard cap on BField domains per fit; over-cap fits fail cleanly (max = unlimited/legacy) | ||
| unsigned minndof_ = 5; // minimum number of DOFs to continue fit | ||
| bool bfcorr_ = true; // whether to make BFieldMap corrections in the fit | ||
| bool zerofield_extrap_ = false; // if true: (1) bfcorr extrapolate() hands off to geometric free-particle continuation outside the map / below minfield_; (2) createDomains stops DomainWalls at that edge instead of failing Extension; (3) replaceDomains charge/mass mismatch soft-keeps the prior usable fit (CH cosmic CRV). Default false preserves legacy LH/CH behaviour. |
There was a problem hiding this comment.
Would zerofield_extrap_ == true result in different behavior for minfield_ = 0? If not, this config could be replaced with (minfield_ > 0)
| if(olditer != last)++olditer; | ||
| } while(olditer != last); | ||
| } | ||
| } catch (...) { |
There was a problem hiding this comment.
what function in the try block is throwing?
| // Degenerate active range (or DomainMargin confinement dropping every domain) can leave | ||
| // createDomains with zero domains → empty fittraj_. createEffects would then throw | ||
| // std::length_error("Empty PiecewiseTrajectory!") and abort the art event. Soft-fail instead. | ||
| if(fittraj_->pieces().empty()){ |
There was a problem hiding this comment.
Can you test for this failure back at line 337?
2 participants
Summary
Bringing up
bfcorr=trueCentralHelix→CRV extrapolation for Mu2e cosmics exposed four ways the BField"domain walk" can fail: an OOM, a multi-minute CPU spin, a lost fit, and an aborted event. They share
one root cause — the domain step has no lower bound and the domain set no upper bound — compounded by
the fact that the CentralHelix parameterization is singular as |B|→0 (
p_T ∝ B/ω; thedPardB ~ bfrac/(1+bfrac)correction has a pole atbfrac = −1). Near-vertical cosmics leaving the DS bore justpast the tracker drive the walk straight into that region.
The failure modes
rangeInTolerancereturns a vanishingdt, andaddDomainappends an effect + traj piece per step without bound.BCorrTolerancedrivescreateDomains'trange→0, sotstartstops advancing and the loop spins emitting deduped domains into astd::set.throw std::invalid_argument("Invalid domains")extendbuilds extension domains in independentcreateDomainscalls, each with its own ±0.5·trangemargin; the blocks don't tile contiguously andcreateEffectsrequires abutting domains.throw std::length_error("Empty PiecewiseTrajectory!")New
Configknobsmindtstep_0.0nsmax(rangeInTolerance, 0) == rangeInTolerancerangeInTolerance→0minfield_0.0Tminfield_ > 0|B|drops below it, before thedPardBpoledomainmargin_DBL_MAXnsmax(x, x−DBL_MAX) == xmaxdomains_UINT_MAXsize() > UINT_MAXnever truezerofield_extrap_falseminfield_Always-on changes
Four changes are not behind a knob. Each only fires where the previous code had undefined behavior or an
uncaught throw, so no configuration reaches them on a well-behaved track:
detectorRangereturns a nullTimeRangewhen there are no active hits/xings, instead of constructingan inverted range that throws.
container.
length_errorof Mat #4.replaceDomainsbuilds the replacement trajectory before clearingdomains_/effects_, so a failedappend leaves the track intact rather than half-mutated. Behavior-preserving on the success path.
Verification
Opt-in / no-effect-when-unconfigured. Mu2e Offline built twice against the same source base, differing
only by this branch versus its base
2011f34, run on identical input with no knob configured:reco/OnSpill.fcldefaults)rtol = 0)KalSeedBField domain boundsbfcorr=false, and every change here is inside abfcorr_branchStability. Event 5160 under a 4 GB address-space cap, knobs off versus on, with nothing else changed:
FatalRootErrorfinalizing a 367,596-domain runaway