Public Member Functions | Private Member Functions | Private Attributes | List of all members
gar::rec::veefinder1 Class Reference
Inheritance diagram for gar::rec::veefinder1:
art::EDProducer art::detail::Producer art::detail::LegacyModule art::Modifier art::ModuleBase art::ProductRegistryHelper

Public Member Functions

 veefinder1 (fhicl::ParameterSet const &p)
 
 veefinder1 (veefinder1 const &)=delete
 
 veefinder1 (veefinder1 &&)=delete
 
veefinder1operator= (veefinder1 const &)=delete
 
veefinder1operator= (veefinder1 &&)=delete
 
void produce (art::Event &e) override
 
- Public Member Functions inherited from art::EDProducer
 EDProducer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDProducer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Producer
virtual ~Producer () noexcept
 
 Producer (fhicl::ParameterSet const &)
 
 Producer (Producer const &)=delete
 
 Producer (Producer &&)=delete
 
Produceroperator= (Producer const &)=delete
 
Produceroperator= (Producer &&)=delete
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Modifier
 ~Modifier () noexcept
 
 Modifier ()
 
 Modifier (Modifier const &)=delete
 
 Modifier (Modifier &&)=delete
 
Modifieroperator= (Modifier const &)=delete
 
Modifieroperator= (Modifier &&)=delete
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Member Functions

bool goodTrack (Track &trk, TrackEnd usebeg)
 
bool goodOtherTrack (Track &trk, TrackEnd usebeg, TVector3 &veepos, TVector3 &veep)
 
int fitVertex (std::vector< TrackPar > &tracks, std::vector< float > &xyz, float &chisquared, std::vector< std::vector< float > > &covmat, double &time, std::vector< TrackEnd > &usebeg, std::vector< float > &dca)
 

Private Attributes

size_t fNTPCClusPerTrack
 min number of TPC clusters for a track to be identified as part of a vee More...
 
float fPCut
 min momentum for a track to be identified as part of a vee More...
 
float fRCut
 DCA cut from track to vee vertex. More...
 
float fDCut
 distance cut from endpoints of tracks More...
 
size_t fRequireNearTrack
 number of nearby tracks to require (primary vertex) More...
 
size_t fNearTrackNTPCClus
 number of nearby tracks' TPC cluster count requirement More...
 
float fNearTrackDMinCut
 minimum distance of nearby track to the vee vertex More...
 
float fNearTrackDMaxCut
 maximum distance of nearby track to the vee vertex More...
 
float fNearTrackAngCut
 pointing requirement of the vee at the nearby track endpoint More...
 
int fPrintLevel
 module print level More...
 
float fOpenAngMinCut
 minimum opening angle cut in 3D More...
 
float fOpenAngMaxCut
 opening angle cut in 3D More...
 
std::string fTrackLabel
 track module label More...
 
bool fRequireOppositeCharge
 whether or not to require opposite-sign tracks. More...
 

Additional Inherited Members

- Public Types inherited from art::EDProducer
using ModuleType = EDProducer
 
using WorkerType = WorkerT< EDProducer >
 
- Public Types inherited from art::detail::Producer
template<typename UserConfig , typename KeysToIgnore = void>
using Table = Modifier::Table< UserConfig, KeysToIgnore >
 
- Public Types inherited from art::Modifier
template<typename UserConfig , typename UserKeysToIgnore = void>
using Table = ProducerTable< UserConfig, detail::ModuleConfig, UserKeysToIgnore >
 
- Static Public Member Functions inherited from art::EDProducer
static void commitEvent (EventPrincipal &ep, Event &e)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollectorconsumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 36 of file veefinder1_module.cc.

Constructor & Destructor Documentation

gar::rec::veefinder1::veefinder1 ( fhicl::ParameterSet const &  p)
explicit

Definition at line 84 of file veefinder1_module.cc.

84  : EDProducer{p}
85  {
86  fTrackLabel = p.get<std::string>("TrackLabel","track");
87  fPCut = p.get<float>("PCut",0.1); // in GeV
88  fRCut = p.get<float>("RCut", 0.1); // in cm
89  fDCut = p.get<float>("DCut",20.); // in cm
90  fOpenAngMinCut = p.get<float>("OpenAngMinCut", 0.05); // in radians
91  fOpenAngMaxCut = p.get<float>("OpenAngMaxCut", 3.09); // in radians
92  fNTPCClusPerTrack = p.get<size_t>("NTPCClusPerTrack",100);
93  fRequireNearTrack = p.get<size_t>("RequireNearTrack",0);
94  fNearTrackNTPCClus = p.get<size_t>("NearTrackNTPCClus",100);
95  fNearTrackDMinCut = p.get<float>("NearTrackDMinCut",2.0);
96  fNearTrackDMaxCut = p.get<float>("NearTrackDMaxCut",20.0);
97  fNearTrackAngCut = p.get<float>("NearTrackAngCut",0.2);
98  fRequireOppositeCharge = p.get<bool>("ReqireOppositeCharge",true);
99  fPrintLevel = p.get<int>("PrintLevel",0);
100 
101  art::InputTag trackTag(fTrackLabel);
102  consumes< std::vector<rec::Track> >(trackTag);
103 
104  // Produce vertex, associated tracks & ends. TrackEnd is defined in Track.h
105  produces< std::vector<rec::Vee> >();
106  produces< art::Assns<rec::Track, rec::Vee, rec::TrackEnd> >();
107  }
float fOpenAngMaxCut
opening angle cut in 3D
float fOpenAngMinCut
minimum opening angle cut in 3D
float fRCut
DCA cut from track to vee vertex.
size_t fNTPCClusPerTrack
min number of TPC clusters for a track to be identified as part of a vee
float fNearTrackAngCut
pointing requirement of the vee at the nearby track endpoint
std::string string
Definition: nybbler.cc:12
int fPrintLevel
module print level
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.h:20
std::string fTrackLabel
track module label
float fDCut
distance cut from endpoints of tracks
float fPCut
min momentum for a track to be identified as part of a vee
p
Definition: test.py:223
bool fRequireOppositeCharge
whether or not to require opposite-sign tracks.
float fNearTrackDMaxCut
maximum distance of nearby track to the vee vertex
size_t fNearTrackNTPCClus
number of nearby tracks&#39; TPC cluster count requirement
float fNearTrackDMinCut
minimum distance of nearby track to the vee vertex
size_t fRequireNearTrack
number of nearby tracks to require (primary vertex)
gar::rec::veefinder1::veefinder1 ( veefinder1 const &  )
delete
gar::rec::veefinder1::veefinder1 ( veefinder1 &&  )
delete

Member Function Documentation

int gar::rec::veefinder1::fitVertex ( std::vector< TrackPar > &  tracks,
std::vector< float > &  xyz,
float &  chisquared,
std::vector< std::vector< float > > &  covmat,
double &  time,
std::vector< TrackEnd > &  usebeg,
std::vector< float > &  dca 
)
private

Definition at line 315 of file veefinder1_module.cc.

323  {
324  // find the ends of the tracks that are closest to the other tracks' ends
325  // pick the end that minimizes the sums of the min(dist) to the other tracks' endpoints
326 
327  // todo -- iterate on this as the tracks are helices and not lines -- move to close point on track and re-do linear extrapolation
328 
329  if (tracks.size() < 2) return(1); // need at least two tracks to make a vertex
330 
331  // find the average time of the tracks
332 
333  time = 0;
334  for (size_t i=0; i<tracks.size(); ++i)
335  {
336  time += tracks.at(i).getTime();
337  }
338  if (tracks.size() > 0)
339  {
340  time /= tracks.size();
341  }
342 
343  xyz.resize(3);
344  dca.clear();
345 
346  std::vector<TVectorF> dir;
347  std::vector<TVectorF> p;
348 
349  for (size_t itrack = 0; itrack < tracks.size(); ++itrack)
350  {
351  TVector3 trackbeg = tracks.at(itrack).getXYZBeg();
352  TVector3 trackend = tracks.at(itrack).getXYZEnd();
353  float phi=0;
354  float si=0;
355  if (usebeg.at(itrack)==TrackEndBeg)
356  {
357  float tmppos[3] = {(float) trackbeg.X(), (float) trackbeg.Y(), (float) trackbeg.Z()};
358  p.emplace_back(3,tmppos);
359  phi = tracks.at(itrack).getTrackParametersBegin()[3];
360  si = TMath::Tan(tracks.at(itrack).getTrackParametersBegin()[4]);
361  }
362  else
363  {
364  float tmppos[3] = {(float) trackend.X(), (float) trackend.Y(), (float) trackend.Z()};
365  p.emplace_back(3,tmppos);
366  phi = tracks.at(itrack).getTrackParametersEnd()[3];
367  si = TMath::Tan(tracks.at(itrack).getTrackParametersEnd()[4]);
368  }
369  TVectorF dtmp(3);
370  dtmp[0] = si;
371  dtmp[1] = TMath::Sin(phi);
372  dtmp[2] = TMath::Cos(phi);
373  float norminv = 1.0/TMath::Sqrt(dtmp.Norm2Sqr());
374  dtmp *= norminv;
375  dir.push_back(dtmp);
376  }
377 
378  TMatrixF I(3,3);
379  I[0][0] = 1;
380  I[1][1] = 1;
381  I[2][2] = 1;
382 
383  TMatrixF A(3,3);
384  TMatrixF VVT(3,3);
385  TMatrixF IMVVT(3,3);
386  TVectorF vsum(3);
387  vsum.Zero();
388  A *= 0;
389 
390  for (size_t itrack=0; itrack < tracks.size(); ++itrack)
391  {
392  for (size_t i=0; i<3; ++i)
393  {
394  for (size_t j=0; j<3; ++j)
395  {
396  VVT[i][j] = dir.at(itrack)[i]*dir.at(itrack)[j];
397  }
398  }
399  IMVVT = I - VVT;
400  A += IMVVT;
401  vsum += IMVVT*p.at(itrack);
402  }
403  double det;
404  A.Invert(&det);
405  if (det == 0) return(1);
406  TVectorF xyzsol = A*vsum;
407  xyz.at(0) = xyzsol[0];
408  xyz.at(1) = xyzsol[1];
409  xyz.at(2) = xyzsol[2];
410 
411  chisquared = 0;
412  for (size_t itrack=0; itrack < tracks.size(); ++itrack)
413  {
414 
415  TVector3 dir3(dir.at(itrack)[0],dir.at(itrack)[1],dir.at(itrack)[2]);
416  TVectorF diff = xyzsol - p.at(itrack);
417  TVector3 diff3(diff[0],diff[1],diff[2]);
418  float dcatrack = (diff3.Cross(dir3)).Mag();
419  dca.push_back(dcatrack);
420  chisquared += dcatrack*dcatrack; // todo -- include track errors in chisquared calc
421  }
422 
423  // to iterate -- extrapolate helical tracks to the closest point to the first found vertex and
424  // run the fitter again
425 
426  // to do: figure out what the vertex covariance matrix is
427 
428  covmat.clear();
429  for (size_t i=0;i<3;++i)
430  {
431  std::vector<float> cmr;
432  for (size_t j=0; j<3; ++j)
433  {
434  cmr.push_back(0);
435  }
436  covmat.push_back(cmr);
437  }
438  return 0;
439 
440  }
string dir
p
Definition: test.py:223
TrackEnd const TrackEndBeg
Definition: Track.h:33
#define A
Definition: memgrp.cpp:38
bool gar::rec::veefinder1::goodOtherTrack ( Track trk,
TrackEnd  usebeg,
TVector3 &  veepos,
TVector3 &  veep 
)
private

Definition at line 468 of file veefinder1_module.cc.

469  {
470  if (trk.NHits() < fNearTrackNTPCClus) // historically-named NHits gets the number of TPC clusters
471  { return false; }
472  TVector3 tep(0,0,0);
473  if (usebeg == TrackEndBeg)
474  {
475  tep.SetXYZ(trk.Vertex()[0],trk.Vertex()[1],trk.Vertex()[2]);
476  }
477  else
478  {
479  tep.SetXYZ(trk.End()[0],trk.End()[1],trk.End()[2]);
480  }
481  TVector3 dv = tep - veepos;
482  float dve = dv.Mag();
483  if (dve < fNearTrackDMinCut || dve > fNearTrackDMaxCut)
484  { return false; }
485  if (dv.Angle(veep) > fNearTrackAngCut)
486  { return false; }
487  return true;
488  }
float fNearTrackAngCut
pointing requirement of the vee at the nearby track endpoint
float fNearTrackDMaxCut
maximum distance of nearby track to the vee vertex
TrackEnd const TrackEndBeg
Definition: Track.h:33
size_t fNearTrackNTPCClus
number of nearby tracks&#39; TPC cluster count requirement
bool gar::rec::veefinder1::goodTrack ( Track trk,
TrackEnd  usebeg 
)
private

Definition at line 445 of file veefinder1_module.cc.

446  {
447  if (trk.NHits() < fNTPCClusPerTrack) // historically-named NHits gets the number of TPC clusters
448  { return false; }
449  if (usebeg == TrackEndBeg)
450  {
451  if (trk.Momentum_beg() < fPCut)
452  { return false; }
453  }
454  else
455  {
456  if (trk.Momentum_end() < fPCut)
457  { return false; }
458  }
459 
460  return true; // if we got here, then the track passed all the cuts
461  }
size_t fNTPCClusPerTrack
min number of TPC clusters for a track to be identified as part of a vee
float fPCut
min momentum for a track to be identified as part of a vee
TrackEnd const TrackEndBeg
Definition: Track.h:33
veefinder1& gar::rec::veefinder1::operator= ( veefinder1 const &  )
delete
veefinder1& gar::rec::veefinder1::operator= ( veefinder1 &&  )
delete
void gar::rec::veefinder1::produce ( art::Event e)
overridevirtual

Implements art::EDProducer.

Definition at line 109 of file veefinder1_module.cc.

109  {
110  std::unique_ptr< std::vector<rec::Vee> >
111  veeCol(new std::vector<rec::Vee>);
112 
113  auto trkVeeAssns = std::make_unique< art::Assns<rec::Track, rec::Vee, rec::TrackEnd> >();
114  //std::unique_ptr< art::Assns<rec::Track, rec::Vee, rec::TrackEnd> >
115  // trkVeeAssns(new art::Assns<rec::Track, rec::Vee, rec::TrackEnd>);
116 
117  auto trackHandle = e.getValidHandle< std::vector<Track> >(fTrackLabel);
118  auto const& tracks = *trackHandle;
119 
120  auto const veePtrMaker = art::PtrMaker<rec::Vee>(e);
121  auto const trackPtrMaker = art::PtrMaker<rec::Track>(e, trackHandle.id());
122 
123  size_t nTrack = tracks.size();
124 
125  const float m_prot = 0.9382723; // in GeV
126  const float m_pi = 0.13957018; // in GeV
127 
128  // For each end of each track, look for candidate vees. Only two tracks per vee, i and j, or this and that
129 
130  if (nTrack>=2)
131  {
132  std::vector<TrackPar> tplist;
133  std::vector<TrackEnd> usebeg;
134  for (size_t iTrack=0; iTrack<nTrack-1; ++iTrack)
135  {
136  Track thisTrack = tracks.at(iTrack);
137  for (int iEnd=0; iEnd<2; ++iEnd)
138  {
139  int thisCharge = 0;
140  TrackEnd iTrackEnd(iEnd);
141  if (!goodTrack(thisTrack,iEnd)) continue;
142  TVector3 thisEndinSpace;
143  TVector3 thisP;
144  if (iTrackEnd==TrackEndBeg)
145  {
146  TVector3 tmpv(thisTrack.Vertex());
147  thisEndinSpace = tmpv;
148  TVector3 tmpv2(thisTrack.VtxDir());
149  tmpv2 *= thisTrack.Momentum_beg();
150  thisP = tmpv2;
151  thisCharge = thisTrack.ChargeBeg();
152  }
153  else
154  {
155  TVector3 tmpv(thisTrack.End());
156  thisEndinSpace = tmpv;
157  TVector3 tmpv2(thisTrack.EndDir());
158  tmpv2 *= thisTrack.Momentum_end();
159  thisP = tmpv2;
160  thisCharge = thisTrack.ChargeEnd();
161  }
162  tplist.clear();
163  tplist.emplace_back(thisTrack);
164  usebeg.clear();
165  usebeg.push_back(iTrackEnd);
166 
167  for (size_t jTrack=iTrack+1; jTrack<nTrack; ++jTrack)
168  {
169  Track thatTrack = tracks.at(jTrack);
170  for (int jEnd = 0; jEnd<2; ++jEnd)
171  {
172  int thatCharge = 0;
173  TrackEnd jTrackEnd(jEnd);
174  if (!goodTrack(thatTrack,jEnd)) continue;
175  TVector3 thatEndinSpace;
176  TVector3 thatP;
177  if (jTrackEnd==TrackEndBeg)
178  {
179  TVector3 tmpv(thatTrack.Vertex());
180  thatEndinSpace = tmpv;
181  TVector3 tmpv2(thatTrack.VtxDir());
182  tmpv2 *= thatTrack.Momentum_beg();
183  thatP = tmpv2;
184  thatCharge = thatTrack.ChargeBeg();
185  }
186  else
187  {
188  TVector3 tmpv(thatTrack.End());
189  thatEndinSpace = tmpv;
190  TVector3 tmpv2(thatTrack.EndDir());
191  tmpv2 *= thatTrack.Momentum_end();
192  thatP = tmpv2;
193  thatCharge = thatTrack.ChargeEnd();
194  }
195 
196  // distance and opening angle cuts
197 
198  if (fPrintLevel>1) std::cout << "veefinder: got a pair of tracks" << std::endl;
199  if ( (thisEndinSpace-thatEndinSpace).Mag() > fDCut ) continue;
200  if (fPrintLevel>1) std::cout << "veefinder: passed the endpoint distance cut" << std::endl;
201  if (fPrintLevel>1) std::cout << "veefinder: opening angle: " << thisP.Angle(thatP) << std::endl;
202  float openang = thisP.Angle(thatP);
203  if ( openang < fOpenAngMinCut) continue;
204  if ( openang > fOpenAngMaxCut) continue;
205  if (fPrintLevel>1) std::cout << "veefinder: passed the opening angle cut" << std::endl;
207  {
208  if (thisCharge == 0 || thatCharge == 0)
209  {
210  if (fPrintLevel>1)
211  {
212  std::cout << "veefinder: One of the track charges is zero: " << thisCharge << " " << thatCharge << std::endl;
213  }
214  }
215  if (thisCharge == thatCharge) continue;
216  if (fPrintLevel>1) std::cout << "veefinder: passed the opposite-charge cut" << std::endl;
217  }
218 
219  tplist.resize(1);
220  tplist.emplace_back(thatTrack);
221  usebeg.resize(1);
222  usebeg.push_back(jTrackEnd);
223 
224  std::vector<float> dcas;
225  std::vector<float> vpos;
226  std::vector<std::vector<float>> covmat;
227  double time;
228  float chisquared;
229  int iret = fitVertex(tplist,vpos,chisquared,covmat,time,usebeg,dcas);
230  if (fPrintLevel>1) std::cout << "veefinder: vertex fitted" << std::endl;
231  if (iret != 0) continue;
232  if (fPrintLevel>1) std::cout << "veefinder: vertex retcode okay" << std::endl;
233  if (dcas[0] > fRCut || dcas[1] > fRCut) continue;
234  if (fPrintLevel>1) std::cout << "veefinder: tracks passed rcut" << std::endl;
235 
236  TVector3 sumP = thisP + thatP;
237  TVector3 vposv(vpos.data());
238 
239  // see if there are nearby tracks that might be from the primary vertex
240 
241  if (fRequireNearTrack > 0)
242  {
243  size_t otherTrackCount = 0;
244  for (size_t kTrack=0; kTrack<nTrack; ++kTrack)
245  {
246  if (kTrack == iTrack || kTrack == jTrack) continue;
247  Track otherTrack = tracks.at(kTrack);
248  for (int kEnd = 0; kEnd<2; ++kEnd)
249  {
250  TrackEnd kTrackEnd(kEnd);
251  if (!goodOtherTrack(otherTrack,kTrackEnd,vposv,sumP)) continue;
252  ++otherTrackCount;
253  break;
254  }
255  }
256  if (otherTrackCount < fRequireNearTrack) continue;
257  }
258  // calculate invariant masses for the three hypotheses
259 
260  std::vector<TLorentzVector> fourmomentumvec;
261  float energy_kshort_hyp =
262  TMath::Sqrt(thisP.Mag2() + m_pi*m_pi) +
263  TMath::Sqrt(thatP.Mag2() + m_pi*m_pi);
264 
265  // sort the lambda hypotheses so that the first one has a positive proton (lambda)
266  // and the second one has a negative proton (lambdabar) if we are requiring opposite charged tracks
267  // if we have two same-charged tracks, the sorting is arbitrary
268 
269  float energy_lambda1_hyp =
270  TMath::Sqrt(thisP.Mag2() + m_prot*m_prot) +
271  TMath::Sqrt(thatP.Mag2() + m_pi*m_pi);
272  float energy_lambda2_hyp =
273  TMath::Sqrt(thisP.Mag2() + m_pi*m_pi) +
274  TMath::Sqrt(thatP.Mag2() + m_prot*m_prot);
275  if (thisCharge == -1)
276  {
277  float tmp = energy_lambda1_hyp;
278  energy_lambda1_hyp = energy_lambda2_hyp;
279  energy_lambda2_hyp = tmp;
280  }
281  fourmomentumvec.emplace_back(sumP,energy_kshort_hyp);
282  fourmomentumvec.emplace_back(sumP,energy_lambda1_hyp);
283  fourmomentumvec.emplace_back(sumP,energy_lambda2_hyp);
284  float covmatvec[9];
285  size_t icounter=0;
286  for (size_t icm=0; icm<3; ++icm)
287  {
288  for (size_t jcm=0;jcm<3; ++jcm)
289  {
290  covmatvec[icounter++] = covmat.at(icm).at(jcm);
291  }
292  }
293  if (fPrintLevel > 1) std::cout << "veefinder: Creating a vee object" << std::endl;
294  veeCol->emplace_back(vpos.data(),covmatvec,chisquared,time,fourmomentumvec.data());
295  auto const veeptr = veePtrMaker(veeCol->size()-1);
296  auto const itrackptr = trackPtrMaker( iTrack );
297  auto const jtrackptr = trackPtrMaker( jTrack );
298  trkVeeAssns->addSingle(itrackptr,veeptr,iTrackEnd);
299  trkVeeAssns->addSingle(jtrackptr,veeptr,jTrackEnd);
300  } // end loop over jEnd
301  } // end loop over jTrack
302  } // end loop over iEnd
303  } // end loop over iTrack
304 
305  } // End test for nTrack>=2; might put empty vertex list & Assns into event
306  e.put(std::move(veeCol));
307  e.put(std::move(trkVeeAssns));
308  } // end produce method.
float fOpenAngMaxCut
opening angle cut in 3D
double m_pi
bool goodOtherTrack(Track &trk, TrackEnd usebeg, TVector3 &veepos, TVector3 &veep)
float fOpenAngMinCut
minimum opening angle cut in 3D
int TrackEnd
Definition: Track.h:32
float fRCut
DCA cut from track to vee vertex.
int fPrintLevel
module print level
std::string fTrackLabel
track module label
float fDCut
distance cut from endpoints of tracks
int fitVertex(std::vector< TrackPar > &tracks, std::vector< float > &xyz, float &chisquared, std::vector< std::vector< float > > &covmat, double &time, std::vector< TrackEnd > &usebeg, std::vector< float > &dca)
const double e
def move(depos, offset)
Definition: depos.py:107
ValidHandle< PROD > getValidHandle(InputTag const &tag) const
Definition: DataViewImpl.h:441
bool goodTrack(Track &trk, TrackEnd usebeg)
string tmp
Definition: languages.py:63
bool fRequireOppositeCharge
whether or not to require opposite-sign tracks.
Definition: tracks.py:1
TrackEnd const TrackEndBeg
Definition: Track.h:33
TrackCollectionProxyElement< TrackCollProxy > Track
Proxy to an element of a proxy collection of recob::Track objects.
Definition: Track.h:1036
QTextStream & endl(QTextStream &s)
size_t fRequireNearTrack
number of nearby tracks to require (primary vertex)

Member Data Documentation

float gar::rec::veefinder1::fDCut
private

distance cut from endpoints of tracks

Definition at line 56 of file veefinder1_module.cc.

float gar::rec::veefinder1::fNearTrackAngCut
private

pointing requirement of the vee at the nearby track endpoint

Definition at line 61 of file veefinder1_module.cc.

float gar::rec::veefinder1::fNearTrackDMaxCut
private

maximum distance of nearby track to the vee vertex

Definition at line 60 of file veefinder1_module.cc.

float gar::rec::veefinder1::fNearTrackDMinCut
private

minimum distance of nearby track to the vee vertex

Definition at line 59 of file veefinder1_module.cc.

size_t gar::rec::veefinder1::fNearTrackNTPCClus
private

number of nearby tracks' TPC cluster count requirement

Definition at line 58 of file veefinder1_module.cc.

size_t gar::rec::veefinder1::fNTPCClusPerTrack
private

min number of TPC clusters for a track to be identified as part of a vee

Definition at line 53 of file veefinder1_module.cc.

float gar::rec::veefinder1::fOpenAngMaxCut
private

opening angle cut in 3D

Definition at line 64 of file veefinder1_module.cc.

float gar::rec::veefinder1::fOpenAngMinCut
private

minimum opening angle cut in 3D

Definition at line 63 of file veefinder1_module.cc.

float gar::rec::veefinder1::fPCut
private

min momentum for a track to be identified as part of a vee

Definition at line 54 of file veefinder1_module.cc.

int gar::rec::veefinder1::fPrintLevel
private

module print level

Definition at line 62 of file veefinder1_module.cc.

float gar::rec::veefinder1::fRCut
private

DCA cut from track to vee vertex.

Definition at line 55 of file veefinder1_module.cc.

size_t gar::rec::veefinder1::fRequireNearTrack
private

number of nearby tracks to require (primary vertex)

Definition at line 57 of file veefinder1_module.cc.

bool gar::rec::veefinder1::fRequireOppositeCharge
private

whether or not to require opposite-sign tracks.

Definition at line 66 of file veefinder1_module.cc.

std::string gar::rec::veefinder1::fTrackLabel
private

track module label

Definition at line 65 of file veefinder1_module.cc.


The documentation for this class was generated from the following file: