diff --git a/CosmicReco/src/LineFinder_module.cc b/CosmicReco/src/LineFinder_module.cc index 9418aa4382..9a7a865246 100644 --- a/CosmicReco/src/LineFinder_module.cc +++ b/CosmicReco/src/LineFinder_module.cc @@ -190,8 +190,8 @@ int LineFinder::findLine(const ComboHitCollection& shC, std::vectorgetStraw(shC[kloc].strawId()); TwoLinePCA pca( strawk.getMidPoint(), strawk.getDirection(), ipos, newdir); - double dist = (pca.point1()-strawk.getMidPoint()).mag(); - if (pca.dca() < _maxDOCA && dist < strawk.halfLength()){ + double dist = (pca.point1()-strawk.getMidPoint()).dot(strawk.getDirection()); + if (pca.dca() < _maxDOCA && fabs(dist) < strawk.halfLength()){ count += 1; ll += pow(dist-shC[kloc].wireDist(),2)/shC[kloc].wireVar(); } diff --git a/Mu2eKinKal/inc/KKFit.hh b/Mu2eKinKal/inc/KKFit.hh index 36f5d2bb87..3b4dee0c0e 100644 --- a/Mu2eKinKal/inc/KKFit.hh +++ b/Mu2eKinKal/inc/KKFit.hh @@ -82,6 +82,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; // using KKPANELHIT = KKPanelHit; using MEAS = KinKal::Hit; using MEASPTR = std::shared_ptr; @@ -99,10 +102,12 @@ namespace mu2e { bool makeStrawHits(Tracker const& tracker,StrawResponse const& strawresponse, BFieldMap const& kkbf, PTRAJ const& ptraj, ComboHitCollection const& chcol, StrawHitIndexCollection const& strawHitIdxs, KKSTRAWHITCOL& hits, KKSTRAWXINGCOL& exings) const; + void makeSeedParamHit(KTRAJ const& seedtraj, std::vector const& paramconstraints, PARAMHITCOL& paramhits) const; // regrow KKTrack components from a KalSeed bool regrowComponents(KalSeed const& kseed, ComboHitCollection const& chcol, mu2e::IndexMap const& strawindexmap, Tracker const& tracker,Calorimeter const& calo, StrawResponse const& strawresponse,BFieldMap const& kkbf, - PTRAJPTR& ptraj, KKSTRAWHITCOL& strawhits, KKCALOHITCOL& calohits, KKSTRAWXINGCOL& exings, DOMAINCOL& domains) const; + PTRAJPTR& ptraj, KKSTRAWHITCOL& strawhits, KKCALOHITCOL& calohits, PARAMHITCOL& paramhits, + KKSTRAWXINGCOL& exings, DOMAINCOL& domains) const; std::shared_ptr caloAxis(CaloCluster const& cluster, Calorimeter const& calo) const; // should come from CaloCluster TODO bool makeCaloHit(CCPtr const& cluster, Calorimeter const& calo, PTRAJ const& pktraj, KKCALOHITCOL& hits) const; // extend a track with a new configuration, optionally searching for and adding hits and straw material @@ -251,11 +256,29 @@ namespace mu2e { return ngood >= minNStrawHits_; } + template void KKFit::makeSeedParamHit(KTRAJ const& seedtraj, std::vector const& paramconstraints, PARAMHITCOL& paramhits) const { + std::array mask = {false}; + KinKal::Parameters cparams = seedtraj.params(); + for (size_t ipar=0;ipar 0){ + mask[ipar] = true; + cparams.covariance()[ipar][ipar] = paramconstraints[ipar]*paramconstraints[ipar]; + }else{ + cparams.covariance()[ipar][ipar] = 1.0; // otherwise inversion fails + } + } + paramhits.push_back(std::make_shared(seedtraj.range().mid(),seedtraj,cparams,mask)); + } + template bool KKFit::regrowComponents(KalSeed const& kseed, // primary event input ComboHitCollection const& chcol, mu2e::IndexMap const& strawindexmap, // ancillary event input Tracker const& tracker,Calorimeter const& calo, // geometries StrawResponse const& strawresponse,BFieldMap const& kkbf, // other conditions - PTRAJPTR& ptraj, KKSTRAWHITCOL& strawhits, KKCALOHITCOL& calohits, KKSTRAWXINGCOL& exings, DOMAINCOL& domains) const { // return values + PTRAJPTR& ptraj, KKSTRAWHITCOL& strawhits, KKCALOHITCOL& calohits, PARAMHITCOL& paramhits, + KKSTRAWXINGCOL& exings, DOMAINCOL& domains) const { // return values GeomHandle kkmat_h; auto const& smat = kkmat_h->strawMaterial(); unsigned ngood(0), nactive(0), nsactive(0); @@ -282,6 +305,9 @@ namespace mu2e { if(kseed.caloCluster()){ makeCaloHit(kseed.caloCluster(),calo,*ptraj,calohits); } + for (auto const& paramhit : kseed.paramHits()){ + paramhits.push_back(std::make_shared(paramhit.time(),*ptraj,paramhit.params(),paramhit.pmask())); + } if(matcorr_){ // add Straw Xings for straws without hits for(auto const& sx : kseed.straws()){ @@ -724,6 +750,9 @@ namespace mu2e { calohit->unbiasedClosestApproach().tpData(), ctres,ca.particleTraj().momentum3(ca.particleToca())); } + for (auto const& paramhit : kktrk.paramHits()){ + kseed._paramhits.emplace_back(paramhit->time(),paramhit->constraintParameters(),paramhit->constraintMask()); + } kseed._straws.reserve(kktrk.strawXings().size()); for(auto const& sxing : kktrk.strawXings()) { // create and fill the flag diff --git a/Mu2eKinKal/inc/KKTrack.hh b/Mu2eKinKal/inc/KKTrack.hh index 58f495cc8d..20f772ba59 100644 --- a/Mu2eKinKal/inc/KKTrack.hh +++ b/Mu2eKinKal/inc/KKTrack.hh @@ -36,6 +36,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using MEAS = KinKal::Hit; using MEASPTR = std::shared_ptr; using MEASCOL = std::vector; @@ -60,10 +63,10 @@ namespace mu2e { using DOMAINCOL = std::set; // construct from configuration, fit environment, and hits and materials KKTrack(Config const& config, BFieldMap const& bfield, KTRAJ const& seedtraj, PDGCode::type tpart, KKSTRAWHITCLUSTERER const& shclusterer, - KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits, std::array constraints = {0}); + KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits, PARAMHITCOL const& paramhits); // construct from regrown constituants KKTrack(Config const& config, BFieldMap const& bfield, PDGCode::type tpart, - PKTRAJPTR& fittraj, KKSTRAWHITCOL& strawhits, KKSTRAWXINGCOL& strawxings, KKCALOHITCOL& calohits, DOMAINCOL& domains); + PKTRAJPTR& fittraj, KKSTRAWHITCOL& strawhits, KKSTRAWXINGCOL& strawxings, KKCALOHITCOL& calohits, PARAMHITCOL const& paramhits, DOMAINCOL& domains); // copy constructor KKTrack(KKTrack const& rhs, CloneContext& context): KinKal::Track(rhs, context), tpart_(rhs.fitParticle()), @@ -118,7 +121,7 @@ namespace mu2e { // extend the track according to new configuration, hits, and/or exings void extendTrack(Config const& config, - KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits ); + KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits); // extend the track to cover a new set of material Xings. This will reuse the existing config object void extendTrack(EXINGCOL const& xings); // add IPA Xing @@ -135,6 +138,7 @@ namespace mu2e { KKSTRAWHITCOL const& strawHits() const { return strawhits_; } KKSTRAWHITCLUSTERCOL const& strawHitClusters() const { return strawhitclusters_; } KKSTRAWXINGCOL const& strawXings() const { return strawxings_; } + PARAMHITCOL const& paramHits() const { return paramhits_; } KKIPAXINGCOL const& IPAXings() const { return ipaxings_; } KKSTXINGCOL const& STXings() const { return stxings_; } KKCRVXINGCOL const& CRVXings() const { return crvxings_; } @@ -153,13 +157,14 @@ namespace mu2e { KKSTRAWHITCOL strawhits_; // straw hits used in this fit KKSTRAWXINGCOL strawxings_; // straw material crossings used in this fit KKCALOHITCOL calohits_; // calo hits used in this fit + PARAMHITCOL paramhits_; KKIPAXINGCOL ipaxings_; // ipa material crossings used in extrapolation KKSTXINGCOL stxings_; // stopping target material crossings used in extrapolation KKCRVXINGCOL crvxings_; // crv crossings using in extrapolation KKINTERCOL inters_; // other recorded intersections KKSTRAWHITCLUSTERCOL strawhitclusters_; // straw hit clusters used in this fit // utility function to convert to generic types - void convertTypes( KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings,KKCALOHITCOL const& calohits, + void convertTypes( KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings,KKCALOHITCOL const& calohits, PARAMHITCOL const& paramhits, MEASCOL& hits, EXINGCOL& exings); // add hits to clusters void addHitClusters(KKSTRAWHITCOL const& strawhits,KKSTRAWXINGCOL const& strawxings, MEASCOL& hits); @@ -170,9 +175,9 @@ namespace mu2e { KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits, - std::array constraints) : + PARAMHITCOL const& paramhits) : KinKal::Track(config,bfield), tpart_(tpart), shclusterer_(shclusterer), - strawhits_(strawhits), strawxings_(strawxings), calohits_(calohits) { + strawhits_(strawhits), strawxings_(strawxings), calohits_(calohits), paramhits_(paramhits) { MEASCOL hits; // polymorphic container of hits EXINGCOL exings; // polymorphic container of detector element crossings // add the hits to clusters, as required @@ -183,44 +188,22 @@ namespace mu2e { shcluster->print(std::cout,this->config().plevel_); } } - convertTypes(strawhits_, strawxings_, calohits_, hits, exings); + convertTypes(strawhits_, strawxings_, calohits_, paramhits_, hits, exings); - std::array mask = {false}; - bool constraining = false; - for (size_t i=0;i 0){ - mask[i] = true; - constraining = true; - } - } - if (constraining){ - KinKal::Parameters cparams = seedtraj.params(); - for (size_t ipar=0;ipar>(seedtraj.range().mid(),seedtraj,cparams,mask)); - } // now fit these this->fit(hits, exings, seedtraj); } template KKTrack::KKTrack(Config const& config, BFieldMap const& bfield, PDGCode::type tpart, - PKTRAJPTR& fittraj, KKSTRAWHITCOL& strawhits, KKSTRAWXINGCOL& strawxings, KKCALOHITCOL& calohits, DOMAINCOL& domains) : + PKTRAJPTR& fittraj, KKSTRAWHITCOL& strawhits, KKSTRAWXINGCOL& strawxings, KKCALOHITCOL& calohits, + PARAMHITCOL const& paramhits, DOMAINCOL& domains) : KinKal::Track(config,bfield), tpart_(tpart),shclusterer_(StrawIdMask::none,0,0.0), - strawhits_(strawhits), strawxings_(strawxings), calohits_(calohits) { + strawhits_(strawhits), strawxings_(strawxings), calohits_(calohits), paramhits_(paramhits) { // convert types MEASCOL hits; // polymorphic container of hits EXINGCOL exings; // polymorphic container of detector element crossings - convertTypes(strawhits_, strawxings_, calohits_, hits, exings); + convertTypes(strawhits_, strawxings_, calohits_, paramhits_, hits, exings); this->fit(hits,exings,domains,fittraj); // add ParameterHit, intersect, extrapolate, ... TODO } @@ -265,11 +248,13 @@ namespace mu2e { KKSTRAWHITCOL const& strawhits, KKSTRAWXINGCOL const& strawxings, KKCALOHITCOL const& calohits, + PARAMHITCOL const& paramhits, MEASCOL& hits, EXINGCOL& exings) { - hits.reserve(strawhits_.size() + calohits_.size()); + hits.reserve(strawhits_.size() + calohits_.size()+paramhits_.size()); exings.reserve(strawxings_.size()); for(auto const& strawhit : strawhits)hits.emplace_back(std::static_pointer_cast(strawhit)); for(auto const& calohit : calohits)hits.emplace_back(std::static_pointer_cast(calohit)); + for(auto const& paramhit : paramhits)hits.emplace_back(std::static_pointer_cast(paramhit)); for(auto const& strawxing : strawxings)exings.emplace_back(std::static_pointer_cast(strawxing)); } @@ -294,7 +279,8 @@ namespace mu2e { } if(nhit != strawhits_.size()+strawhits.size()) std::cout << "cluster hit sum doesn't match " << nhit << " " << strawhits_.size() << std::endl; } - convertTypes(strawhits,strawxings,calohits,hits,exings); + PARAMHITCOL paramhits; + convertTypes(strawhits,strawxings,calohits,paramhits,hits,exings); this->extend(config,hits,exings); // store the new hits strawhits_.reserve(strawhits_.size()+strawhits.size()); diff --git a/Mu2eKinKal/src/CentralHelixFit_module.cc b/Mu2eKinKal/src/CentralHelixFit_module.cc index ed57b279c1..ba96bd6135 100644 --- a/Mu2eKinKal/src/CentralHelixFit_module.cc +++ b/Mu2eKinKal/src/CentralHelixFit_module.cc @@ -91,6 +91,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using KKFIT = KKFit; using KinKal::VEC3; using KinKal::DMAT; @@ -113,7 +116,7 @@ namespace mu2e { fhicl::OptionalAtom fixedBField { Name("ConstantBField"), Comment("Constant BField value") }; fhicl::Atom seedMom { Name("SeedMomentum"), Comment("Seed momentum") }; fhicl::Atom seedCharge { Name("SeedCharge"), Comment("Seed charge MAGNITUDE, in electron charge units") }; - fhicl::Sequence parconst { Name("ParameterConstraints"), Comment("External constraint on parameters to seed values (rms, various units)") }; + fhicl::Sequence parconst { Name("ParameterConstraints"), Comment("External constraint on parameters to seed values (rms, various units)") }; fhicl::Sequence sampleSurfaces { Name("SampleSurfaces"), Comment("When creating the KalSeed, sample the fit at these surfaces") }; fhicl::Atom sampleInRange { Name("SampleInRange"), Comment("Require sample times to be inside the fit trajectory time range") }; fhicl::Atom sampleInBounds { Name("SampleInBounds"), Comment("Require sample intersection point be inside surface bounds (within tolerance)") }; @@ -165,6 +168,7 @@ namespace mu2e { bool fixedfield_; // double seedMom_; int seedCharge_; + std::vector paramconstraints_; double intertol_; // surface intersection tolerance (mm) double sampletbuff_; // simple time buffer; replace this with extrapolation TODO bool useFitCharge_; // Set the PDG particle to agree with the fit charge @@ -172,7 +176,7 @@ namespace mu2e { bool sampleinrange_, sampleinbounds_; // require samples to be in range or on surface SurfaceIdCollection ssids_; KinKalGeom::SurfacePairCollection surfacess_to_sample_; // surfaces to sample the fit - std::array paramconstraints_; + bool constraining_; }; CentralHelixFit::CentralHelixFit(const Parameters& settings) : art::EDProducer{settings}, @@ -189,6 +193,7 @@ namespace mu2e { fixedfield_(false), seedMom_(settings().modSettings().seedMom()), seedCharge_(settings().modSettings().seedCharge()), + paramconstraints_(settings().modSettings().parconst()), intertol_(settings().modSettings().interTol()), sampletbuff_(settings().modSettings().sampleTBuff()), useFitCharge_(settings().modSettings().useFitCharge()), @@ -208,9 +213,15 @@ namespace mu2e { for(size_t ipar=0;ipar < seederrors.size(); ++ipar){ seedcov_[ipar][ipar] = seederrors[ipar]*seederrors[ipar]; } - auto const& parerrors = settings().modSettings().parconst(); - for (size_t ipar=0;ipar 0) + constraining_ = true; + } + }else if (paramconstraints_.size() > 0){ + throw cet::exception("RECO")<<"mu2e::CentralHelixFit: Parameter constraint configuration error"<< endl; + } if(print_ > 0) std::cout << "Fit " << config_ << "Extension " << exconfig_; double bz(0.0); @@ -330,8 +341,14 @@ namespace mu2e { try { seedtraj.range() = kkfit_.range(strawhits,calohits,strawxings); + // create parameter constraint, use seedtrajectory parameters + PARAMHITCOL paramhits; + if (constraining_){ + kkfit_.makeSeedParamHit(seedtraj,paramconstraints_,paramhits); + } + // create and fit the track - auto ktrk = make_unique(config_,*kkbf_,seedtraj,fitpart,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits,paramconstraints_); + auto ktrk = make_unique(config_,*kkbf_,seedtraj,fitpart,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits,paramhits); // Check the fit auto goodfit = goodFit(*ktrk); // if we have an extension schedule, extend. diff --git a/Mu2eKinKal/src/KinematicLineFit_module.cc b/Mu2eKinKal/src/KinematicLineFit_module.cc index 76daf75111..293240b649 100644 --- a/Mu2eKinKal/src/KinematicLineFit_module.cc +++ b/Mu2eKinKal/src/KinematicLineFit_module.cc @@ -79,6 +79,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using MEAS = KinKal::Hit; using MEASPTR = std::shared_ptr; using MEASCOL = std::vector; @@ -111,7 +114,7 @@ namespace mu2e { struct KKLineModuleConfig : public KKModuleConfig { fhicl::Sequence seedCollections {Name("CosmicTrackSeedCollections"), Comment("Seed fit collections to be processed ") }; fhicl::Atom seedmom { Name("SeedMomentum"), Comment("Initial momentum value")}; - fhicl::Sequence paramconstraints { Name("ParameterConstraints"), Comment("Sigma of direct gaussian constraints on each parameter (0=no constraint)")}; + fhicl::Sequence paramconstraints { Name("ParameterConstraints"), Comment("Sigma of direct gaussian constraints on each parameter (0=no constraint)")}; fhicl::Sequence sampleSurfaces { Name("SampleSurfaces"), Comment("When creating the KalSeed, sample the fit at these surfaces") }; fhicl::Atom sampleInRange { Name("SampleInRange"), Comment("Require sample times to be inside the fit trajectory time range") }; fhicl::Atom sampleInBounds { Name("SampleInBounds"), Comment("Require sample intersection point be inside surface bounds (within tolerance)") }; @@ -150,11 +153,12 @@ namespace mu2e { ProditionsHandle alignedTracker_h_; int print_; float seedmom_; + std::vector paramconstraints_; PDGCode::type fpart_; TrkFitDirection fdir_; KKFIT kkfit_; // fit helper DMAT seedcov_; // seed covariance matrix - std::array paramconstraints_; + bool constraining_; double mass_; // particle mass int charge_; // particle charge std::unique_ptr kkbf_; @@ -175,6 +179,7 @@ namespace mu2e { saveall_(settings().modSettings().saveAll()), print_(settings().modSettings().printLevel()), seedmom_(settings().modSettings().seedmom()), + paramconstraints_(settings().modSettings().paramconstraints()), fpart_(static_cast(settings().modSettings().fitParticle())), kkfit_(settings().mu2eSettings()), intertol_(settings().modSettings().interTol()), @@ -198,14 +203,13 @@ namespace mu2e { for(size_t ipar=0;ipar < seederrors.size(); ++ipar){ seedcov_[ipar][ipar] = seederrors[ipar]*seederrors[ipar]; } - auto const& tempconstraints = settings().modSettings().paramconstraints(); - if (tempconstraints.size() == 0){ - for (size_t ipar=0;ipar 0) + constraining_ = true; + } + }else if (paramconstraints_.size() > 0){ throw cet::exception("RECO")<<"mu2e::KinematicLineFit: Parameter constraint configuration error"<< endl; } for(auto const& sidname : settings().modSettings().sampleSurfaces()) { @@ -288,12 +292,19 @@ namespace mu2e { } // set the seed range given the hit TPOCA values seedtraj.range() = kkfit_.range(strawhits,calohits, strawxings); + + // create parameter constraint, use seedtrajectory parameters + PARAMHITCOL paramhits; + if (constraining_){ + kkfit_.makeSeedParamHit(seedtraj,paramconstraints_,paramhits); + } + if(print_ > 0){ //std::cout << "Seed line parameters " << hseed.track() << std::endl; seedtraj.print(std::cout,print_); } // create and fit the track - auto kktrk = make_unique(config_,*kkbf_,seedtraj,fpart_,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits,paramconstraints_); + auto kktrk = make_unique(config_,*kkbf_,seedtraj,fpart_,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits,paramhits); auto goodfit = goodFit(*kktrk); if(goodfit && exconfig_.schedule().size() > 0){ kkfit_.extendTrack(exconfig_,*kkbf_, *tracker,*strawresponse, chcol, *calo_h, cc_H, *kktrk ); diff --git a/Mu2eKinKal/src/LoopHelixFit_module.cc b/Mu2eKinKal/src/LoopHelixFit_module.cc index 0bcd575bcd..5f9056cf7f 100644 --- a/Mu2eKinKal/src/LoopHelixFit_module.cc +++ b/Mu2eKinKal/src/LoopHelixFit_module.cc @@ -84,6 +84,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using KKFIT = KKFit; using KinKal::VEC3; using KinKal::DMAT; @@ -338,10 +341,11 @@ namespace mu2e { if (kkfit_.useCalo() && hseed.caloCluster().isNonnull()) { kkfit_.makeCaloHit(hseed.caloCluster(),*calo_h, pseedtraj, calohits); } + PARAMHITCOL paramhits; // set the seed range given the hits and xings seedtraj.range() = kkfit_.range(strawhits,calohits,strawxings); // create and fit the track - auto ktrk = make_unique(config_,*kkbf_,seedtraj,fitpart,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits); + auto ktrk = make_unique(config_,*kkbf_,seedtraj,fitpart,kkfit_.strawHitClusterer(),strawhits,strawxings,calohits,paramhits); if(!ktrk) // check that the track exists throw cet::exception("RECO")<<"mu2e::LoopHelixFit: Track fit was performed but no track is found\n"; diff --git a/Mu2eKinKal/src/RegrowKinematicLine_module.cc b/Mu2eKinKal/src/RegrowKinematicLine_module.cc index ff02d69854..c4438e581a 100644 --- a/Mu2eKinKal/src/RegrowKinematicLine_module.cc +++ b/Mu2eKinKal/src/RegrowKinematicLine_module.cc @@ -116,6 +116,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using KKFIT = KKFit; using MEAS = KinKal::Hit; @@ -209,13 +212,14 @@ namespace mu2e { KKSTRAWXINGCOL strawxings; strawxings.reserve(kseed.straws().size()); KKCALOHITCOL calohits; + PARAMHITCOL paramhits; DOMAINCOL domains; // create the trajectory. This is done here to be strongly typed auto goodhits = kkfit_.regrowComponents(kseed, chcol, indexmap, *tracker,*calo_h,*strawresponse,*kkbf_, - trajptr, strawhits, calohits, strawxings, domains); + trajptr, strawhits, calohits, paramhits, strawxings, domains); if(debug_ > 0){ - std::cout << "Regrew " << strawhits.size() << " straw hits, " << strawxings.size() << " straw xings, " << calohits.size() << " CaloHits and " << domains.size() << " domains, status = " << goodhits << std::endl; + std::cout << "Regrew " << strawhits.size() << " straw hits, " << strawxings.size() << " straw xings, " << calohits.size() << " CaloHits, " << paramhits.size() << " ParameterHits, and " << domains.size() << " domains, status = " << goodhits << std::endl; } if(debug_ > 2){ unsigned nhactive(0); @@ -228,7 +232,7 @@ namespace mu2e { if(debug_ > 5)static_cast*>(trajptr.get())->print(std::cout,2); if(goodhits){ // create the KKTrack from these components - auto ktrk = std::make_unique(config_,*kkbf_,kseed.particle(),trajptr,strawhits,strawxings,calohits,domains); + auto ktrk = std::make_unique(config_,*kkbf_,kseed.particle(),trajptr,strawhits,strawxings,calohits,paramhits,domains); if(ktrk && ktrk->fitStatus().usable()){ if(debug_ > 0) std::cout << "RegrowKinematicLine: successful track refit" << std::endl; // extrapolate as requested diff --git a/Mu2eKinKal/src/RegrowLoopHelix_module.cc b/Mu2eKinKal/src/RegrowLoopHelix_module.cc index 085ed2166d..e4a3f93a32 100644 --- a/Mu2eKinKal/src/RegrowLoopHelix_module.cc +++ b/Mu2eKinKal/src/RegrowLoopHelix_module.cc @@ -119,6 +119,9 @@ namespace mu2e { using KKCALOHIT = KKCaloHit; using KKCALOHITPTR = std::shared_ptr; using KKCALOHITCOL = std::vector; + using PARAMHIT = KinKal::ParameterHit; + using PARAMHITPTR = std::shared_ptr; + using PARAMHITCOL = std::vector; using KKFIT = KKFit; using MEAS = KinKal::Hit; @@ -216,13 +219,14 @@ namespace mu2e { KKSTRAWXINGCOL strawxings; strawxings.reserve(kseed.straws().size()); KKCALOHITCOL calohits; + PARAMHITCOL paramhits; DOMAINCOL domains; // create the trajectory. This is done here to be strongly typed auto goodhits = kkfit_.regrowComponents(kseed, chcol, indexmap, *tracker,*calo_h,*strawresponse,*kkbf_, - trajptr, strawhits, calohits, strawxings, domains); + trajptr, strawhits, calohits, paramhits, strawxings, domains); if(debug_ > 1){ - std::cout << "Regrew " << strawhits.size() << " straw hits, " << strawxings.size() << " straw xings, " << calohits.size() << " CaloHits and " << domains.size() << " domains, status = " << goodhits << std::endl; + std::cout << "Regrew " << strawhits.size() << " straw hits, " << strawxings.size() << " straw xings, " << calohits.size() << " CaloHits, " << paramhits.size() << " ParameterHits, and " << domains.size() << " domains, status = " << goodhits << std::endl; } if(debug_ > 2){ unsigned nhactive(0); @@ -235,7 +239,7 @@ namespace mu2e { // require hits and consistent BField domains if(goodhits && (domains.size() > 0 || !config_.bfcorr_)){ // create the KKTrack from these components - auto ktrk = std::make_unique(config_,*kkbf_,kseed.particle(),trajptr,strawhits,strawxings,calohits,domains); + auto ktrk = std::make_unique(config_,*kkbf_,kseed.particle(),trajptr,strawhits,strawxings,calohits,paramhits,domains); if(ktrk && ktrk->fitStatus().usable()){ if(debug_ > 0) std::cout << "RegrowLoopHelix: successful track refit" << std::endl; if(extend_)kkfit_.extendTrack(config_,*kkbf_, *tracker,*strawresponse, chcol, *calo_h, cc_H , *ktrk ); diff --git a/RecoDataProducts/inc/KalSeed.hh b/RecoDataProducts/inc/KalSeed.hh index 56f72e9ed9..ab01e5db8f 100644 --- a/RecoDataProducts/inc/KalSeed.hh +++ b/RecoDataProducts/inc/KalSeed.hh @@ -11,6 +11,7 @@ #include "Offline/RecoDataProducts/inc/TrkStrawHitSeed.hh" #include "Offline/RecoDataProducts/inc/TrkStrawHitCalib.hh" #include "Offline/RecoDataProducts/inc/TrkCaloHitSeed.hh" +#include "Offline/RecoDataProducts/inc/TrkParamHitSeed.hh" #include "Offline/RecoDataProducts/inc/TrkStraw.hh" #include "Offline/RecoDataProducts/inc/KalSegment.hh" #include "Offline/RecoDataProducts/inc/KalIntersection.hh" @@ -45,6 +46,7 @@ namespace mu2e { auto const& hitCalibInfos() const { return _hitcalibs;} auto const& caloHit() const { return _chit; } auto const& straws() const { return _straws;} + auto const& paramHits() const { return _paramhits;} auto const& segments() const { return _segments; } auto const& intersections() const { return _inters; } auto const& status() const { return _status; } @@ -90,6 +92,7 @@ namespace mu2e { std::vector _hits; // hit seeds for all the hits used in this fit std::vector _hitcalibs; // extra calibration/alignment info std::vector _straws; // straws interesected by this fit + std::vector _paramhits; // parameter constraint hits std::vector _domainbounds; // domain time boundaries TrkCaloHitSeed _chit; // CaloCluster-based hit. If it has no CaloCluster, this has no content // static value used in regrowing diff --git a/RecoDataProducts/inc/TrkParamHitSeed.hh b/RecoDataProducts/inc/TrkParamHitSeed.hh new file mode 100644 index 0000000000..b3d2ed899c --- /dev/null +++ b/RecoDataProducts/inc/TrkParamHitSeed.hh @@ -0,0 +1,29 @@ +// +// Persistent representation of a KinKal ParameterHit, used in the +// persistent representation of Kalman Fit +// +#ifndef RecoDataProducts_TrkParamHitSeed_HH +#define RecoDataProducts_TrkParamHitSeed_HH +#include "KinKal/Detector/ParameterHit.hh" +namespace mu2e { + struct TrkParamHitSeed { + using PMASK = std::array; // parameter mask + // default constructor: initialization on declaration + TrkParamHitSeed() {} + //KinKal constructor + TrkParamHitSeed(double time, KinKal::Parameters const& params, PMASK const& pmask) : + time_(time), params_(params), pmask_(pmask) {} + + double time() const { return time_; } + KinKal::Parameters const& params() const { return params_; } + PMASK const& pmask() const { return pmask_; } + + // + // Payload + // + double time_; // time of this constraint: must be supplied on construction and does not change + KinKal::Parameters params_; // constraint parameters with covariance + PMASK pmask_; // subset of parmeters to constrain + }; +} +#endif diff --git a/RecoDataProducts/src/classes_def.xml b/RecoDataProducts/src/classes_def.xml index 41fda60dd4..51b6bf092e 100644 --- a/RecoDataProducts/src/classes_def.xml +++ b/RecoDataProducts/src/classes_def.xml @@ -195,6 +195,7 @@ + @@ -202,6 +203,7 @@ +