47 Name(
"inputPFParticleLabel"),
48 Comment(
"Label of recob::PFParticle Collection to be fit")};
50 Name(
"inputTracksLabel"),
51 Comment(
"Label of recob::Track Collection to be fit")};
53 Name(
"inputShowersLabel"),
54 Comment(
"Label of recob::Shower Collection (associated to PFParticles) to be fit")};
56 Name(
"inputCaloLabel"),
57 Comment(
"Label of anab::Calorimetry Collection, matching inputTracksLabel, to be used for " 58 "initial momentum estimate. Used only if momFromCalo is set to true.")};
61 Comment(
"Label of sim::MCTrack Collection to be used for initial momentum estimate. Used " 62 "only if momFromMC is set to true.")};
64 Name(
"inputPidLabel"),
65 Comment(
"Label of anab::ParticleID Collection, matching inputTracksLabel, to be used for " 73 Comment(
"If true extract tracks from inputPFParticleLabel " 74 "collection, if false from inputTracksLabel.")};
77 Comment(
"If true extract showers from inputPFParticleLabel collection.")};
79 Name(
"momFromMSChi2"),
80 Comment(
"Flag used to get initial momentum estimate from " 81 "trkf::TrackMomentumCalculator::GetMomentumMultiScatterChi2().")};
83 Comment(
"Flag used to get initial momentum estimate from " 84 "trkf::TrackMomentumCalculator::GetTrackMomentum().")};
87 Comment(
"Flag used to get initial momentum estimate from inputCaloLabel collection.")};
90 Comment(
"Flag used to get initial momentum estimate from inputMCLabel collection.")};
92 Name(
"momentumInGeV"),
93 Comment(
"Fixed momentum estimate value, to be used when momFromCalo, momFromMSChi2, " 94 "momFromLength and momFromMC are all false, or if the estimate is not available.")};
96 Comment(
"Flag used to get particle ID estimate from corresponding " 97 "PFParticle. Needs trackFromPF=true.")};
99 Name(
"idFromCollection"),
100 Comment(
"Flag used to get particle ID estimate from inputPidLabel collection.")};
103 Comment(
"Default particle id hypothesis in case no valid id is provided either via " 104 "PFParticle or in the ParticleId collection.")};
106 Name(
"dirFromVtxPF"),
107 Comment(
"Assume track direction from Vertex in PFParticle. Needs trackFromPF=true.")};
110 Comment(
"Assume track direction from as the one giving positive " 111 "dot product with vector specified by dirVec.")};
113 Comment(
"Fhicl sequence defining the vector used when " 114 "dirFromVec=true. It must have 3 elements.")};
116 Name(
"alwaysInvertDir"),
117 Comment(
"If true, fit all tracks from end to vertex assuming inverted direction.")};
119 Name(
"produceTrackFitHitInfo"),
120 Comment(
"Option to produce (or not) the detailed TrackFitHitInfo.")};
122 Name(
"produceSpacePoints"),
123 Comment(
"Option to produce (or not) the associated SpacePoints.")};
125 Name(
"keepInputTrajectoryPoints"),
126 Comment(
"Option to keep positions and directions from input trajectory/track. The fit will " 127 "provide only covariance matrices, chi2, ndof, particle Id and absolute momentum. " 128 "It may also modify the trajectory point flags. In order to avoid inconsistencies, " 129 "it has to be used with the following fitter options all set to false: " 130 "sortHitsByPlane, sortOutputHitsMinLength, skipNegProp.")};
168 std::unique_ptr<art::FindManyP<anab::Calorimetry>>
trackCalo;
169 std::unique_ptr<art::FindManyP<anab::ParticleID>>
trackId;
177 const int pId)
const;
178 int setPId(
const unsigned int iTrack,
180 const int pfPid = 0)
const;
198 ,
inputFromPF{
p_().options().trackFromPF() ||
p_().options().showerFromPF()}
203 if (
p_().options().showerFromPF())
208 if (
p_().options().idFromCollection())
212 if (
p_().options().pFromMC() ||
p_().options().dirFromMC())
215 produces<std::vector<recob::Track>>();
216 produces<art::Assns<recob::Track, recob::Hit, recob::TrackHitMeta>>();
217 produces<art::Assns<recob::Track, recob::Hit>>();
218 if (
inputFromPF) { produces<art::Assns<recob::PFParticle, recob::Track>>(); }
219 if (
p_().options().produceTrackFitHitInfo()) {
220 produces<std::vector<std::vector<recob::TrackFitHitInfo>>>();
222 if (
p_().options().produceSpacePoints()) {
223 produces<std::vector<recob::SpacePoint>>();
224 produces<art::Assns<recob::Hit, recob::SpacePoint>>();
228 if (
p_().options().trackFromPF() == 0 &&
p_().options().idFromPF())
230 <<
"Incompatible configuration parameters: cannot use idFromPF=true with trackFromPF=false." 232 if (
p_().options().trackFromPF() == 0 &&
p_().options().dirFromVtxPF())
234 <<
"Incompatible configuration parameters: cannot use dirFromVtxPF=true with " 238 unsigned int nIds = 0;
239 if (
p_().options().idFromPF()) nIds++;
240 if (
p_().options().idFromCollection()) nIds++;
243 <<
"Incompatible configuration parameters: only at most one can be set to true among " 244 "idFromPF and idFromCollection." 248 unsigned int nDirs = 0;
249 if (
p_().options().dirFromVtxPF()) nDirs++;
250 if (
p_().options().dirFromMC()) nDirs++;
251 if (
p_().options().dirFromVec()) nDirs++;
252 if (
p_().options().alwaysInvertDir()) nDirs++;
255 <<
"Incompatible configuration parameters: only at most one can be set to true among " 256 "dirFromVtxPF, dirFromMC, dirFromVec, and alwaysInvertDir." 260 unsigned int nPFroms = 0;
261 if (
p_().options().pFromCalo()) nPFroms++;
262 if (
p_().options().pFromMSChi2()) nPFroms++;
263 if (
p_().options().pFromLength()) nPFroms++;
264 if (
p_().options().pFromMC()) nPFroms++;
267 <<
"Incompatible configuration parameters: only at most one can be set to true among " 268 "pFromCalo, pFromMSChi2, pFromLength, and pFromMC." 272 if (
p_().options().keepInputTrajectoryPoints()) {
273 if (
p_().fitter().sortHitsByPlane() ||
p_().fitter().sortOutputHitsMinLength() ||
274 p_().fitter().skipNegProp()) {
276 <<
"Incompatible configuration parameters: keepInputTrajectoryPoints needs the following " 277 "fitter options all set to false: sortHitsByPlane, sortOutputHitsMinLength, skipNegProp." 282 if (
p_().options().showerFromPF()) {
283 if (nPFroms > 0 || nIds > 0 || nDirs > 0) {
285 <<
"Incompatible configuration parameters: showerFromPF currently does not support " 286 "optional momentum values, particle hypotheses and directions." 294 auto outputTracks = std::make_unique<std::vector<recob::Track>>();
295 auto outputHitsMeta =
296 std::make_unique<art::Assns<recob::Track, recob::Hit, recob::TrackHitMeta>>();
297 auto outputHits = std::make_unique<art::Assns<recob::Track, recob::Hit>>();
298 auto outputHitInfo = std::make_unique<std::vector<std::vector<recob::TrackFitHitInfo>>>();
300 auto const tid = e.
getProductID<std::vector<recob::Track>>();
303 auto outputSpacePoints = std::make_unique<std::vector<recob::SpacePoint>>();
304 auto outputHitSpacePointAssn = std::make_unique<art::Assns<recob::Hit, recob::SpacePoint>>();
305 auto const spid = e.
getProductID<std::vector<recob::SpacePoint>>();
311 if (
p_().options().pFromMC() ||
p_().options().dirFromMC()) {
314 for (
unsigned int iMC = 0; iMC < simTracks->size(); ++iMC) {
317 if (mctrack.
PdgCode() != 13)
continue;
318 if (mctrack.
Process() !=
"primary")
continue;
320 mcdir = TVector3(mctrack.
Start().
Momentum().X() * 0.001 / pMC,
331 auto outputPFAssn = std::make_unique<art::Assns<recob::PFParticle, recob::Track>>();
334 if (
p_().options().trackFromPF())
337 if (
p_().options().showerFromPF())
343 for (
unsigned int iPF = 0; iPF < inputPFParticle->size(); ++iPF) {
345 if (
p_().options().trackFromPF()) {
346 const std::vector<art::Ptr<recob::Track>>& tracks =
assocTracks->at(iPF);
347 auto const& tkHitsAssn =
349 const std::vector<art::Ptr<recob::Vertex>>& vertices =
assocVertices->at(iPF);
351 if (
p_().options().pFromCalo()) {
355 for (
unsigned int iTrack = 0; iTrack < tracks.size(); ++iTrack) {
359 const int pId =
setPId(iTrack,
trackId, inputPFParticle->at(iPF).PdgCode());
361 const bool flipDir =
setDirFlip(track, mcdir, &vertices);
364 std::vector<art::Ptr<recob::Hit>> inHits;
365 for (
auto it = tkHitsAssn.begin(); it != tkHitsAssn.end(); ++it) {
366 if (it->first == ptrack)
367 inHits.push_back(it->second);
368 else if (inHits.size() > 0)
373 std::vector<art::Ptr<recob::Hit>> outHits;
388 if (!fitok)
continue;
390 if (
p_().options().keepInputTrajectoryPoints()) {
394 outputTracks->emplace_back(std::move(outTrack));
397 for (
auto const& trhit : outHits) {
400 outputHitsMeta->addSingle(aptr, trhit, metadata);
401 outputHits->addSingle(aptr, trhit);
409 if (
p_().options().showerFromPF()) {
411 const std::vector<art::Ptr<recob::Shower>>& showers =
assocShowers->at(iPF);
412 if (showers.size() == 0)
continue;
413 auto const& pfClustersAssn =
415 auto const& clHitsAssn =
417 std::vector<art::Ptr<recob::Hit>> inHits;
418 for (
auto itpf = pfClustersAssn.begin(); itpf != pfClustersAssn.end(); ++itpf) {
419 if (itpf->first == pPF) {
421 for (
auto it = clHitsAssn.
begin(); it != clHitsAssn.
end(); ++it) {
422 if (it->first == clust) inHits.push_back(it->second);
425 else if (inHits.
size() > 0)
428 for (
unsigned int iShower = 0; iShower < showers.size(); ++iShower) {
431 std::vector<art::Ptr<recob::Hit>> outHits;
437 auto pid =
p_().options().pdgId();
438 auto mom =
p_().options().pval();
444 std::vector<recob::TrajectoryPointFlags>(),
451 if (!fitok)
continue;
453 outputTracks->emplace_back(std::move(outTrack));
456 for (
auto const& trhit : outHits) {
459 outputHitsMeta->addSingle(aptr, trhit, metadata);
460 outputHits->addSingle(aptr, trhit);
461 if (
p_().options().produceSpacePoints() && outputTracks->back().HasValidPoint(ip)) {
462 auto& tp = outputTracks->back().Trajectory().LocationAtPoint(ip);
463 double fXYZ[3] = {tp.X(), tp.Y(), tp.Z()};
464 double fErrXYZ[6] = {0};
466 outputSpacePoints->emplace_back(std::move(sp));
468 outputHitSpacePointAssn->addSingle(trhit, apsp);
477 e.put(std::move(outputTracks));
478 e.put(std::move(outputHitsMeta));
479 e.put(std::move(outputHits));
480 e.put(std::move(outputPFAssn));
481 if (
p_().options().produceTrackFitHitInfo()) { e.put(std::move(outputHitInfo)); }
482 if (
p_().options().produceSpacePoints()) {
483 e.put(std::move(outputSpacePoints));
484 e.put(std::move(outputHitSpacePointAssn));
493 if (
p_().options().pFromCalo()) {
497 if (
p_().options().idFromCollection()) {
501 for (
unsigned int iTrack = 0; iTrack < inputTracks->size(); ++iTrack) {
507 const bool flipDir =
setDirFlip(track, mcdir);
510 std::vector<art::Ptr<recob::Hit>> inHits;
511 for (
auto it = tkHitsAssn.begin(); it != tkHitsAssn.end(); ++it) {
512 if (it->first == ptrack)
513 inHits.push_back(it->second);
514 else if (inHits.size() > 0)
519 std::vector<art::Ptr<recob::Hit>> outHits;
534 if (!fitok)
continue;
536 if (
p_().options().keepInputTrajectoryPoints()) {
540 outputTracks->emplace_back(std::move(outTrack));
543 for (
auto const& trhit : outHits) {
546 outputHitsMeta->addSingle(aptr, trhit, metadata);
547 outputHits->addSingle(aptr, trhit);
548 if (
p_().options().produceSpacePoints() && outputTracks->back().HasValidPoint(ip)) {
549 auto& tp = outputTracks->back().Trajectory().LocationAtPoint(ip);
550 double fXYZ[3] = {tp.X(), tp.Y(), tp.Z()};
551 double fErrXYZ[6] = {0};
553 outputSpacePoints->emplace_back(std::move(sp));
555 outputHitSpacePointAssn->addSingle(trhit, apsp);
561 e.
put(std::move(outputTracks));
562 e.
put(std::move(outputHitsMeta));
563 e.
put(std::move(outputHits));
564 if (
p_().options().produceTrackFitHitInfo()) { e.
put(std::move(outputHitInfo)); }
565 if (
p_().options().produceSpacePoints()) {
566 e.
put(std::move(outputSpacePoints));
567 e.
put(std::move(outputHitSpacePointAssn));
579 std::vector<Point_t> positions(np);
580 std::vector<Vector_t> momenta(np);
581 std::vector<recob::TrajectoryPointFlags> outFlags(np);
583 for (
unsigned int p = 0; p < np; ++p) {
586 auto op = flag.fromHit();
589 auto mask = flag.mask();
600 std::move(covs.first),
601 std::move(covs.second),
605 for (
auto h : inHits)
606 outHits.push_back(h);
615 double result =
p_().options().pval();
617 else if (
p_().options().pFromLength()) {
620 else if (
p_().options().pFromCalo()) {
622 const std::vector<art::Ptr<anab::Calorimetry>>&
calo =
trackCalo->at(ptrack.
key());
623 double sumenergy = 0.;
625 for (
auto caloit : calo) {
626 if (caloit->KineticEnergy() > 0.) {
627 sumenergy += caloit->KineticEnergy();
631 if (nviews != 0 && sumenergy != 0.) {
633 result = sumenergy / (nviews * 1000.);
636 else if (
p_().options().pFromMC() && pMC > 0.) {
670 if (
p_().options().alwaysInvertDir()) {
return true; }
671 else if (
p_().options().dirFromMC()) {
673 if ((mcdir.X() * tdir.X() + mcdir.Y() * tdir.Y() + mcdir.Z() * tdir.Z()) < 0.) result =
true;
675 else if (
p_().options().dirFromVec()) {
676 std::array<float, 3>
dir =
p_().options().dirVec();
678 if ((dir[0] * tdir.X() + dir[1] * tdir.Y() + dir[2] * tdir.Z()) < 0.) result =
true;
680 else if (
p_().options().trackFromPF() &&
p_().options().dirFromVtxPF() && vertices->size() > 0) {
683 (*vertices)[0]->XYZ(xyz);
686 if (((xyz[0] - te.X()) * (xyz[0] - te.X()) + (xyz[1] - te.Y()) * (xyz[1] - te.Y()) +
687 (xyz[2] - te.Z()) * (xyz[2] - te.Z())) >
688 ((xyz[0] - tv.X()) * (xyz[0] - tv.X()) + (xyz[1] - tv.Y()) * (xyz[1] - tv.Y()) +
689 (xyz[2] - tv.Z()) * (xyz[2] - tv.Z())))
const TVector3 & ShowerStart() const
double VertexMomentum() const
Trajectory_t const & Trajectory() const
Returns the plain trajectory of this object.
Fit tracks using Kalman Filter fit+smooth.
void initTrackFitInfos()
initialize the output vector of TrackFitHitInfos
static constexpr Flag_t NoPoint
The trajectory point is not defined.
int ParticleId() const
Access to various track properties.
ProductID getProductID(std::string const &instance_name="") const
recob::tracking::Point_t Point_t
recob::tracking::Vector_t Vector_t
const recob::TrackTrajectory & Trajectory() const
Access to the stored recob::TrackTrajectory.
art::InputTag pidInputTag
Declaration of signal hit object.
EDProducer(fhicl::ParameterSet const &pset)
recob::tracking::SMatrixSym55 SMatrixSym55
float Chi2() const
Access to various track properties.
size_t NumberTrajectoryPoints() const
Various functions related to the presence and the number of (valid) points.
TrackStatePropagator prop
Vector_t VertexDirection() const
Access to track direction at different points.
bool fitTrack(detinfo::DetectorPropertiesData const &detProp, const recob::TrackTrajectory &traj, int tkID, const SMatrixSym55 &covVtx, const SMatrixSym55 &covEnd, const std::vector< art::Ptr< recob::Hit >> &hits, const double pval, const int pdgid, const bool flipDirection, recob::Track &outTrack, std::vector< art::Ptr< recob::Hit >> &outHits, trkmkr::OptionalOutputs &optionals) const
Fit track starting from TrackTrajectory.
Class for propagation of a trkf::TrackState to a recob::tracking::Plane.
art::InputTag simTrackInputTag
art::InputTag pfParticleInputTag
std::pair< SMatrixSym55, SMatrixSym55 > Covariances() const
Access to position, momentum or covariance at the start and end of the track.
art::InputTag caloInputTag
trkf::TrackKalmanFitter kalmanFitter
double GetMomentumMultiScatterChi2(art::Ptr< recob::Track > const &trk, const bool checkValidPoints=false, const int maxMomentum_MeV=7500)
Calculate muon momentum (GeV) using multiple coulomb scattering. Chi2 minimization of the Highland fo...
std::unique_ptr< art::FindManyP< recob::Vertex > > assocVertices
const_iterator begin() const
Vector_t DirectionAtPoint(size_t i) const
Computes and returns the direction of the trajectory at a point.
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
double setMomValue(art::Ptr< recob::Track > ptrack, const std::unique_ptr< art::FindManyP< anab::Calorimetry >> &trackCalo, const double pMC, const int pId) const
double Length(size_t p=0) const
Access to various track properties.
EDProductGetter const * productGetter(ProductID const pid) const
#define DEFINE_ART_MODULE(klass)
trkf::TrackMomentumCalculator tmc
A trajectory in space reconstructed from hits.
key_type key() const noexcept
int setPId(const unsigned int iTrack, const std::unique_ptr< art::FindManyP< anab::ParticleID >> &trackId, const int pfPid=0) const
art::InputTag trackInputTag
Point_t const & LocationAtPoint(size_t i) const
Returns the position at the specified trajectory point.
art::InputTag showerInputTag
Provides recob::Track data product.
std::unique_ptr< art::FindManyP< anab::ParticleID > > trackId
const TVector3 & Direction() const
KalmanFilterFinalTrackFitter(Parameters const &p)
Declaration of cluster object.
Class def header for mctrack data container.
void produce(art::Event &e) override
int Ndof() const
Access to various track properties.
Point_t const & Vertex() const
Returns the position of the first valid point of the trajectory [cm].
const TLorentzVector & Momentum() const
const_iterator end() const
A trajectory in space reconstructed from hits.
KalmanFilterFinalTrackFitter & operator=(KalmanFilterFinalTrackFitter const &)=delete
ValidHandle< PROD > getValidHandle(InputTag const &tag) const
const SMatrixSym55 & EndCovarianceLocal5D() const
Accessors to track parameters and covariance matrices in Local5D and Global6D coordinates.
std::vector< recob::TrackFitHitInfo > trackFitHitInfos()
get the output vector of TrackFitHitInfos by releasing and moving
std::unique_ptr< art::FindManyP< anab::Calorimetry > > trackCalo
std::unique_ptr< art::FindManyP< recob::Shower > > assocShowers
Point_t const & End() const
Returns the position of the last valid point of the trajectory [cm].
const std::string & Process() const
PointFlags_t const & FlagsAtPoint(size_t i) const
Access to i-th TrajectoryPoint or its Flags.
const MCStep & Start() const
TrackCollectionProxyElement< TrackCollProxy > Track
Proxy to an element of a proxy collection of recob::Track objects.
const SMatrixSym55 & VertexCovarianceLocal5D() const
Accessors to track parameters and covariance matrices in Local5D and Global6D coordinates.
bool setDirFlip(const recob::Track &track, TVector3 &mcdir, const std::vector< art::Ptr< recob::Vertex >> *vertices=0) const
Struct holding optional TrackMaker outputs.
Set of flags pertaining a point of the track.
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
double GetTrackMomentum(double trkrange, int pdg) const
cet::coded_exception< error, detail::translate > exception
std::unique_ptr< art::FindManyP< recob::Track > > assocTracks
void restoreInputPoints(const recob::Trajectory &track, const std::vector< art::Ptr< recob::Hit >> &inHits, recob::Track &outTrack, std::vector< art::Ptr< recob::Hit >> &outHits) const