LCFIPlus  0.6.5
lcfiplus.h
Go to the documentation of this file.
1 #ifndef lcfiplus_h
2 #define lcfiplus_h 1
3 
4 #include <vector>
5 #include <map>
6 #include <iostream>
7 #include <sstream>
8 #include <exception>
9 #include "Rtypes.h"
10 #include "TObject.h"
11 #include "TLorentzVector.h"
12 #include "TVector3.h"
13 
14 #include "EventStore.h"
15 #include "EVENT/Cluster.h"
16 
17 #include <typeinfo>
18 
19 using namespace std;
20 
21 namespace lcfiplus {
22 
23 // used for: track parameters, track covariance matrix, calorimeter subsystems
24 namespace tpar {
25 enum par { d0=0, z0, ph, om, td, parN };
26 enum cov { d0d0=0, d0ph, phph, d0om, phom, omom, d0z0,
28  };
29 enum hit { VTX=0, FTD, SIT, TPC, SET, ETD, hitN };
30 enum calo { ecal=0, hcal, yoke, lcal, lhcal, bcal, caloN };
31 }
32 
33 template <class T> struct DeleteVector {
34  bool operator() (T x) const {
35  delete x;
36  return true;
37  }
38 };
39 
40 template <class T> vector<const T*>* constVector(vector<T*>* ptr) {
41  return reinterpret_cast<vector<const T*> *>(ptr);
42 }
43 template <class T> vector<const T*>& constVector(vector<T*>& ref) {
44  return *reinterpret_cast<vector<const T*> *>(&ref);
45 }
46 template <class T> const vector<const T*>* constVector(const vector<T*>* ptr) {
47  return reinterpret_cast<const vector<const T*> *>(ptr);
48 }
49 template <class T> const vector<const T*>& constVector(const vector<T*>& ref) {
50  return *reinterpret_cast<const vector<const T*> *>(&ref);
51 }
52 
53 // used for error rescaling of the track parameters
54 // in bins of cosTheta and momentum
55 struct ErrorRescale {
56  int var; // 0: d0, 1: z0
57  double cmin; // min cosTheta
58  double cmax; // max cosTheta
59  double pmin; // min momentum
60  double pmax; // max momentum
61  double pull; // pull - scale factor
62  double pullerr; // pull error
63 };
64 
65 // basic types
66 class Event;
67 class Track;
68 class Neutral;
69 class MCParticle;
70 class MCColorSinglet;
71 class Vertex;
72 class Jet;
73 
74 // const vectors
75 typedef const vector<const Track*> TrackVec;
76 typedef const vector<const Neutral*> NeutralVec;
77 typedef const vector<const MCParticle*> MCParticleVec;
78 typedef const vector<const MCColorSinglet*> MCColorSingletVec;
79 typedef const vector<const Vertex*> VertexVec;
80 typedef const vector<const Jet*> JetVec;
81 
82 // const iterators
83 typedef vector<const Track*>::const_iterator TrackVecIte;
84 typedef vector<const Neutral*>::const_iterator NeutralVecIte;
85 typedef vector<const MCParticle*>::const_iterator MCParticleVecIte;
86 typedef vector<const MCColorSinglet*>::const_iterator MCColorSingletVecIte;
87 typedef vector<const Vertex*>::const_iterator VertexVecIte;
88 typedef vector<const Jet*>::const_iterator JetVecIte;
89 
90 class Exception : std::exception {
91  public:
92  Exception(const char* message)
93  : _message(message)
94  {
95  }
96  virtual ~Exception()throw() {}
97 
98  virtual const char* what() const throw() {
99  return _message.c_str();
100  }
101  void Print() {
102  cerr << "lcfiplus_Exception: " << _message << endl;
103  }
104 
105  private:
106  string _message;
107 };
108 
109 // global parameters
110 class Globals {
111  public:
112  ~Globals();
113  static Globals* Instance();
114  void setBField(double bField) {
115  _bField=bField;
116  }
117  double getBField()const {
118  return _bField;
119  }
120  void setBeamSizeX(double beamSizeX) {
121  _beamSizeX = beamSizeX;
122  }
123  double getBeamSizeX()const {
124  return _beamSizeX;
125  }
126  void setBeamSizeY(double beamSizeY) {
127  _beamSizeY = beamSizeY;
128  }
129  double getBeamSizeY()const {
130  return _beamSizeY;
131  }
132  void setBeamSizeZ(double beamSizeZ) {
133  _beamSizeZ = beamSizeZ;
134  }
135  double getBeamSizeZ()const {
136  return _beamSizeZ;
137  }
138 
139  private:
140  // value of b field to be read from GEAR
141  double _bField;
142  double _beamSizeX;
143  double _beamSizeY;
144  double _beamSizeZ;
145  static Globals* _theInstance;
146 
147  Globals();
148 };
149 
150 
151 // generic parameter class
152 class Parameters {
153  public:
154  Parameters(bool as = true)
155  : _map()
156  , _allowstring(as)
157  {}
158 
159  // destructor for destroying objects in the map
161  clear();
162  }
163  // copy operator/constructor
164  Parameters(const Parameters& ref)
165  : _map(ref._map)
166  , _allowstring(ref._allowstring)
167  {
168  *this = ref;
169  }
170  Parameters& operator =(const Parameters& ref); // in lcfiplus.cc
171 
172  // fetch for non-vector
173  template<typename T> void fetch(const char* key, T& ret, T def = T())const {
174  if (_map.find(key) == _map.end()) {
175  cout << "Parameter " << key << " not found." << endl;
176  ret = def;
177  } else if (_map.find(key)->second.first == &typeid(T))ret = *(T*)(_map.find(key)->second.second);
178  else if (_allowstring && _map.find(key)->second.first == &typeid(string)) {
179  istringstream str(*(string*)_map.find(key)->second.second);
180  str >> ret;
181  } else if (_allowstring && _map.find(key)->second.first == &typeid(vector<string>)) {
182  istringstream str((*(const vector<string>*)_map.find(key)->second.second)[0]);
183  str >> ret;
184  } else {
185  cout << "Parameter type invalid: key = " << key << ", type = " << _map.find(key)->second.first->name() << " vs " << typeid(T).name() << endl;
186  throw (Exception("Parameter type invalid."));
187  }
188  return;
189  }
190 
191  // fetch for vector
192  template<typename T> void fetchArray(const char* key, vector<T>& ret)const {
193  if (_map.find(key) == _map.end()) {
194  cout << "Parameter " << key << " not found." << endl;
195  } else if (_map.find(key)->second.first == &typeid(vector<T>))ret = *(vector<T>*)(_map.find(key)->second.second);
196  else if (_map.find(key)->second.first == &typeid(T))ret.push_back(*(T*)(_map.find(key)->second.second));
197  else if (_allowstring && _map.find(key)->second.first == &typeid(string)) {
198  ret.push_back(T());
199 
200  istringstream str(*(string*)_map.find(key)->second.second);
201  str >> ret[0];
202  } else if (_allowstring && _map.find(key)->second.first == &typeid(vector<string>)) {
203  const vector<string>* svec = (const vector<string>*)_map.find(key)->second.second;
204  ret.resize(svec->size());
205  for (unsigned int n=0; n<svec->size(); n++) {
206  istringstream str((*svec)[n]);
207  str >> ret[n];
208  }
209  } else {
210  cout << "Parameter type invalid: key = " << key << ", type = " << _map.find(key)->second.first->name() << " vs " << typeid(T).name() << endl;
211  throw (Exception("Parameter type invalid."));
212  }
213  }
214 
215  public:
216  // for non-vector only (if string parameter)
217  template<typename T> T get(const char* key, T def = T())const {
218  T ret;
219  fetch(key, ret, def);
220  return ret;
221  }
222 
223  template<typename T> std::vector<T> getVec(const char* key)const {
224  std::vector<T> ret;
225  fetchArray(key, ret);
226  return ret;
227  }
228 
229  bool exist(const char* key)const {
230  return _map.find(key) != _map.end();
231  }
232 
233  template<typename T> void add(const char* key, T data) { // change T& data to T data for convenience
234  if (_map.find(key) != _map.end())throw(Exception("Double entry."));
235 
236  _map[key] = pair<const type_info*, void*>(&typeid(T), new T(data));
237 
238  }
239  void remove(const char* key) {
240  remove(key, true);
241  }
242  void clear();
243 
244  template<typename T> void assign(const char* key, T data) {
245  if (_map.find(key) == _map.end())throw(Exception("Parameters::assign(): the key has not been registered."));
246  else if (_map.find(key)->second.first != &typeid(T))throw(Exception("Parameters::assign(): the value type is imcompatible."));
247 
248  // assign
249  *static_cast<T*>(_map.find(key)->second.second) = data;
250  }
251 
252  // direct accessors: avoiding copy to handle names
253  const map<string, pair<const type_info*, void*> >& paramMap()const {
254  return _map;
255  }
256 
257  private:
258  void remove(const string& key, bool delmap);
259 
260  template<typename T> void deleteData(const string& key) {
261  T* p = static_cast<T*>(_map[key].second);
262  delete p;
263  }
264 
265  map<string, pair<const type_info*, void*> > _map;
266  bool _allowstring;
267 };
268 
269 class Algorithm {
270  public:
271  Algorithm() : _param(0) {
272  }
273  virtual ~Algorithm() {}
274 
275  virtual void init(Parameters* param) {
276  _param = param;
277  }
278  virtual void process() = 0;
279  virtual void end() {}
280 
281  protected:
283  return _param;
284  }
285  private:
286  Parameters* _param;
287 
288  ClassDef(Algorithm,1);
289 };
290 
291 class Event : public EventStore {
292  public:
293  ~Event();
294 
295  // singleton
296  static Event* Instance();
297 
298  // standard collections retrievers
299  // for other collections: use EventStore::Get()
300  const vector<const Track*>& getTracks(const char* trackname = 0) const;
301  const vector<const Neutral*>& getNeutrals(const char* neutralname = 0) const;
302  const vector<const MCParticle*>& getMCParticles(const char* mcpname = 0) const;
303  const vector<const MCColorSinglet*>& getMCColorSinglets(const char* mcpname = 0) const;
304  const Vertex* getPrimaryVertex(const char* privtxname = 0) const;
305  const vector<const Vertex*>& getSecondaryVertices(const char* secvtxname = 0) const;
306  const vector<const Jet*>& getJets(const char* jetname = 0) const;
307 
308  // standard collection name setter/getter
309  void setDefaultTracks(const char* name) {
310  _defTrackName = name;
311  }
312  void setDefaultNeutrals(const char* name) {
313  _defNeutralName = name;
314  }
315  void setDefaultMCParticles(const char* name) {
316  _defMcpName = name;
317  }
318  void setDefaultPrimaryVertex(const char* name) {
319  _defPriVtxName = name;
320  }
321  void setDefaultSecondaryVertices(const char* name) {
322  _defSecVtxName = name;
323  }
324  void setDefaultJets(const char* name) {
325  _defJetName = name;
326  }
327 
328  const char* getDefaultTracks()const {
329  return _defTrackName.c_str();
330  }
331  const char* getDefaultNeutrals()const {
332  return _defNeutralName.c_str();
333  }
334  const char* getDefaultMCParticles()const {
335  return _defMcpName.c_str();
336  }
337  const char* getDefaultPrimaryVertex()const {
338  return _defPriVtxName.c_str();
339  }
340  const char* getDefaultSecondaryVertices()const {
341  return _defSecVtxName.c_str();
342  }
343  const char* getDefaultJets()const {
344  return _defJetName.c_str();
345  }
346 
347  // utility functions for MCParticles
348  const MCParticle* getMCParticle(int id) const;
349  const MCParticle* getMCParticle(const Track* trk) const;
350 
351  vector<const MCParticle*> mcGetColorStrings() const;
352  int mcNumberOfB() const;
353  int mcNumberOfC() const;
354  vector<const MCParticle*> mcGetSemiStableBs() const;
355  vector<const MCParticle*> mcGetSemiStableCs() const;
356  vector<const MCParticle*> mcGetSemiStableBCs(bool separatebc) const;
357  int mcFindParent(MCParticleVec& vec, const MCParticle* p) const;
358 
359 // void rescaleErrors();
360 
361  private:
362 
363  // singleton
364  static Event* _theInstance;
365 
366  // private constructor
367  Event();
368 
369  // members
370 // int id;
371 // vector<ErrorRescale> _errRescaleTable;
372 // void readErrorRescaleTable();
373 // bool _rescaleTableRead;
374 
375  string _defTrackName;
376  string _defNeutralName;
377  string _defMcpName;
378  string _defPriVtxName;
379  string _defSecVtxName;
380  string _defJetName;
381 
382 };
383 
384 class Track : public TLorentzVector {//, protected TrackData {//, public EventPointer {
385  public:
386  // constructor
387  // Track() {}
388 // Track(const TrackData& d);
389  ~Track() {}
390 
391  int getId() const {
392  return _id;
393  }
394  void setId(int id) {
395  _id = id;
396  }
397 
398  const lcfiplus::MCParticle* getMcp() const {
399  return _mcp;
400  }
401  void setMcp(const lcfiplus::MCParticle* mcp) {
402  _mcp = mcp;
403  }
404 
405  int getPDG() const {
406  return _pdg;
407  }
408  void setPDG(int pdg) {
409  _pdg = pdg;
410  }
411 
412  double getCharge() const {
413  return _charge;
414  }
415  void setCharge(double charge) {
416  _charge = charge;
417  }
418 
419  double getD0() const {
420  return _par[tpar::d0];
421  }
422  double getZ0() const {
423  return _par[tpar::z0];
424  }
425  double getPhi() const {
426  return _par[tpar::ph];
427  }
428  double getOmega() const {
429  return _par[tpar::om];
430  }
431  double getTanLambda() const {
432  return _par[tpar::td];
433  }
434  void setHelix(double d0, double z0, double phi, double omega, double tanLambda) {
435  _par[tpar::d0] = d0;
436  _par[tpar::z0] = z0;
437  _par[tpar::ph] = phi;
438  _par[tpar::om] = omega;
439  _par[tpar::td] = tanLambda;
440  }
441  void setHelix(double* par) {
442  memcpy(_par, par, sizeof(_par));
443  }
444 
445  const double* getCovMatrix() const {
446  return _cov;
447  }
448  void setCovMatrix(double* cov);
449 
450  int getVtxHits() const {
451  return _nhits[tpar::VTX];
452  }
453  int getFtdHits() const {
454  return _nhits[tpar::FTD];
455  }
456  int getSitHits() const {
457  return _nhits[tpar::SIT];
458  }
459  int getTpcHits() const {
460  return _nhits[tpar::TPC];
461  }
462  int getSetHits() const {
463  return _nhits[tpar::SET];
464  }
465  int getEtdHits() const {
466  return _nhits[tpar::ETD];
467  }
468  void setTrackHits(int* hits) {
469  memcpy(_nhits, hits, sizeof(_nhits));
470  }
471 
472  const double* getCaloEdep() const {
473  return _calo;
474  }
475  void setCaloEdep(double* calo) {
476  memcpy(_calo, calo, sizeof(_calo));
477  }
478 
479  double getRadiusOfInnermostHit() const {
480  return _rimh;
481  }
482  void setRadiusOfInnermostHit(double rimh) {
483  _rimh = rimh;
484  }
485 
486  double getChi2() const {
487  return _chi2;
488  }
489  void setChi2(double chi2) {
490  _chi2 = chi2;
491  }
492 
493  int getNdf() const {
494  return _ndf;
495  }
496  void setNdf(int ndf) {
497  _ndf = ndf;
498  }
499 
500  void setFlightLength(double flt) const {
501  _flt = flt;
502  }
503  double getFlightLength() const {
504  return _flt;
505  }
506 
507  void setParticleIDProbability(string parName, double pidProbability) {
508  _pidProbability[parName] = pidProbability;
509  }
510 
511  double getParticleIDProbability(const char* parName) const {
512  string partName = parName;
513  double prob=0.0;
514  if (_pidProbability.find(partName) != _pidProbability.end()){
515  prob = _pidProbability.at(partName); // [] is not a const function
516  }
517  return prob;
518  }
519 
520  //for track energy correction
521  void setCorrEnergy(double mass) {
522  _correnergy = sqrt(this->P()*this->P() + mass * mass);
523  }
524 
525  double getCorrEnergy() const {
526  return _correnergy;
527  }
528 
529  void swapEnergy() {
530  double tmpe = this->E();
531  this->SetE(_correnergy);
532  _correnergy = tmpe;
533  }
534 
535  //for storing BNess
536  void setBNess(double bness) {
537  _bness = bness;
538  }
539 
540  double getBNess() const {
541  return _bness;
542  }
543 
544  double getCNess() const {
545  return _cness;
546  }
547 
548  //for storing CNess
549  void setCNess(double cness) {
550  _cness = cness;
551  }
552  double getX() const;
553  double getY() const;
554  double getZ() const;
555  TVector3 getPos()const {
556  return TVector3(getX(), getY(), getZ());
557  }
558  TVector3 momentumAtVertex( const Vertex* vtx ) const;
559 
560  private:
561  mutable double _flt;
562 
563  int _id;
564  const lcfiplus::MCParticle* _mcp;
565  int _pdg;
566  double _charge;
567 
568  // track parameter
569  double _par[tpar::parN];
570  // covariance matrix
571  double _cov[tpar::covN];
572  int _nhits[tpar::hitN];
573  double _calo[tpar::caloN];
574 
575  double _rimh; // radius of innermost hit
576  double _chi2; // from track fit
577  int _ndf;
578 
579  //ParticleID posterior probability
580  map<string, double> _pidProbability;
581  double _correnergy;
582 
583  //BNess
584  double _bness, _cness;
585 
586  ClassDef(lcfiplus::Track, 2);
587 };
588 
589 class Neutral : public TLorentzVector { // : public LorentzVector, protected NeutralData{//, public EventPointer {
590  public:
591  // constructor
592 // Neutral(const NeutralData& d);
593  Neutral() : _id(0), _pdg(0), _isV0(false), _mcp(0), _clstr() {}
594  ~Neutral() {}
595 
596  int getId() const {
597  return _id;
598  }
599  void setId(int id) {
600  _id = id;
601  }
602  const lcfiplus::MCParticle* getMcp() const {
603  return _mcp;
604  }
605  void setMcp(const lcfiplus::MCParticle* mcp) {
606  _mcp = mcp;
607  }
608  int getPDG() const {
609  return _pdg;
610  }
611  void setPDG(int pdg) {
612  _pdg = pdg;
613  }
614 
615  const double* getCaloEdep() const {
616  return _calo;
617  }
618  //TODO: range check
619  void setCaloEdep(double* calo) {
620  memcpy(_calo, calo, sizeof(_calo));
621  }
622 
623  bool isV0() const {
624  return _isV0;
625  }
626  void setV0(bool v0=true) {
627  _isV0 = v0;
628  }
629 
630  void setClusters(EVENT::ClusterVec clu) { _clstr = clu; }
631  EVENT::ClusterVec getClusters() const { return _clstr; }
632 
633  private:
634  int _id;
635  int _pdg;
636  int _isV0;
637  const lcfiplus::MCParticle* _mcp;
638 
639  double _calo[tpar::caloN];
640 
641  EVENT::ClusterVec _clstr;
642 
643  ClassDef(lcfiplus::Neutral, 3);
644 };
645 
646 class MCParticle : public TLorentzVector {
647 
648  public:
649 
650  //ctor/dtors ////////////////////////////////////////////////////////
651  MCParticle(int id, int pdg, MCParticle* parent, double charge, const TLorentzVector& p, const TVector3& v) {
652  Init(id, pdg, parent, charge, p, v);
653  }
656 
657  // initialization: non-const parent is needed because of adding me to the daughter list
658  void Init(int id, int pdg, MCParticle* parent, double charge, const TLorentzVector& p, const TVector3& v);
659 
660  // simple accessors ///////////////////////////////////////////////////
661  int getId() const {
662  return _id;
663  }
664  void setId(int id) {
665  _id = id;
666  }
667 
668  int getPDG() const {
669  return _pdg;
670  }
671  void setPDG(int pdg) {
672  _pdg = pdg;
673  }
674 
675  double getCharge() const {
676  return _charge;
677  }
678  void setCharge(double charge) {
679  _charge = charge;
680  }
681 
682  const TVector3& getVertex() const {
683  return _vertex;
684  }
685  void setVertex(const TVector3& v) {
686  _vertex = v;
687  }
688  const TVector3& getEndVertex()const;
689 
690  const MCParticle* getParent() const {
691  return _parent;
692  }
693  void setParent(const MCParticle* parent) {
694  _parent = parent;
695  }
696  //vector<MCParticle*> getParents() const;
697 
698  const vector<const MCParticle*>& getDaughters() const {
699  return _dau;
700  }
701  void addDaughter(const MCParticle* mcp);
702 
703  bool isParent(const MCParticle* mcp)const;
704 
705  // more intelligent accessors
706  int getFlavor() const;
707  const MCParticle* getColorString()const;
708 
709  const MCColorSinglet* getColorSinglet(const vector<const MCColorSinglet*>* pcs)const;
710 
711  bool isStableTrack() const;
712  bool isStable() const;
713  bool isSemiStableB() const;
714  bool isSemiStableC() const;
715  bool isSemiStableS() const;
716  bool isSemiStable() const;
717  const MCParticle* semileptonicDecay() const; // return MCParticle of lepton
718  int getFlavorTagCategory() const;
719  const MCParticle* getSemiStableParent() const;
720  const MCParticle* getSemiStableBParent() const;
721  const MCParticle* getSemiStableCParent() const;
722 
723  // end points
724  double decayDistance()const;
725  const MCParticle* findDauForDecay()const;
726  double getEx()const;
727  double getEy()const;
728  double getEz()const;
729 
730  vector<const lcfiplus::MCParticle*> promptTracks()const;
731 
732  // helix stuff
733  double getD0()const;
734  double getZ0()const;
735  double getPhi()const;
736  double getOmega()const;
737  double getTanLambda()const;
738 
739  private:
740  // basic properties
741  int _id;
742  int _pdg;
743  double _charge;
744  const MCParticle* _parent;
745  TVector3 _vertex;
746  // list of daughters
747  vector<const lcfiplus::MCParticle*> _dau;
748 
749  mutable const MCParticle* _dauForDecay;
750 
751  ClassDef(lcfiplus::MCParticle, 2);
752 };
753 
755  public:
756  const MCParticle* getMcp()const {
757  return _mcp;
758  }
759 
760  const MCParticle* _mcp;
761 
762  vector<const lcfiplus::MCParticle*> _initials;
763  vector<const lcfiplus::MCParticle*> _finalstrings;
764  vector<const lcfiplus::MCParticle*> _qqgluons;
765  vector<const lcfiplus::MCParticle*> _finalcolorsinglets;
766  vector<const lcfiplus::MCParticle*> _realparticles;
767 
768  ClassDefNV(lcfiplus::MCColorSinglet, 1);
769 };
770 
771 class Vertex {
772  friend class LcfiInterface;
773 
774  public:
775  enum vtx { xx=0, xy, yy, xz, yz, zz };
776 
777  Vertex() : _id(-1), _chi2(0), _prob(0), _x(0), _y(0), _z(0), _isPrimary(false), _tracks(), _chi2Tracks(), _vertexingname(), _rvtxmom(), _pi0mom(), _rvtxmass(-1), _npi0(-1) {}
778  Vertex(const double chi2, const double prob, const double x, const double y, const double z, const double cov[6], bool isPrim)
779  : _id(-1), _chi2(chi2), _prob(prob), _x(x), _y(y), _z(z), _isPrimary(isPrim), _tracks(), _chi2Tracks(), _vertexingname(), _rvtxmom(), _pi0mom(), _rvtxmass(-1), _npi0(-1) {
780  if (cov == 0) {
781  memset(_cov,0, sizeof(_cov));
782  } else {
783  memcpy(_cov, cov, sizeof(_cov));
784  }
785  }
786 
787  // id is not copied
788  Vertex(const Vertex& from) : _id(-1), _chi2(from._chi2), _prob(from._prob), _x(from._x), _y(from._y), _z(from._z), _isPrimary(from._isPrimary), _tracks(from._tracks), _chi2Tracks(from._chi2Tracks), _vertexingname(from._vertexingname), _rvtxmom(from._rvtxmom), _pi0mom(from._pi0mom), _rvtxmass(from._rvtxmass), _npi0(from._npi0) {
789  memcpy(_cov, from._cov, sizeof(_cov));
790  }
791 
792  ~Vertex() {};
793  void add(const Track* trk);
794  void add(const Track* trk,double chi2);
795  void setId(int id)const {
796  _id = id;
797  }
798 
799  int getId() const {
800  return _id;
801  }
802  double getChi2() const {
803  return _chi2;
804  }
805  double getProb() const {
806  return _prob;
807  }
808  double getX() const {
809  return _x;
810  }
811  double getY() const {
812  return _y;
813  }
814  double getZ() const {
815  return _z;
816  }
817  TVector3 getPos() const {
818  return TVector3(_x, _y,_z);
819  }
820  const double* getCov() const {
821  return _cov;
822  }
823  const vector<const Track*>& getTracks() const {
824  return _tracks;
825  }
826  const map<const lcfiplus::Track*, double>& getTracksChi2Map() const {
827  return _chi2Tracks;
828  }
829  double getChi2Track(const Track* tr)const {
830  map<const Track*,double>::const_iterator it = _chi2Tracks.find(tr);
831  if (it != _chi2Tracks.end())return it->second;
832  else return -1;
833  }
834  const Track* getWorstTrack() const;
835  double getChi2TrackFit(const Track* tr, int mode=1)const;
836 
837  double length(const Vertex* primary=0) const;
838  double significance(const Vertex* primary) const;
839 
840  bool isPrimary()const {
841  return _isPrimary;
842  }
843  void setPrimary(bool isPrim) {
844  _isPrimary = isPrim;
845  }
846 
847  double getPparallel(const TVector3& axis)const;
848  double getVertexMass(const Vertex* daughter = 0, const TVector3* paxis = 0, const double dmass = 1.87, double* ppt = 0, double* pp = 0)const;
849  double getVertexAngle(const Vertex* daughter, const Vertex* primary = 0)const;
850 
851  const MCParticle* getMcp()const;
852 
853  double dirdot(const Vertex* primary=0) const;
854  bool passesV0selection(const Vertex* primary=0) const;
855  TLorentzVector getFourMomentum() const;
856 
857  //for AVF
858  void setVertexingName(string vtxname){ _vertexingname = vtxname;}
859  string getVertexingName() const { return _vertexingname;}
861 
862  //for vertex mass recovery
863  double getRecoveredVertexMass() const {return _rvtxmass;}
864  TLorentzVector getRecoveredFourMomentum() const {return _rvtxmom;}
865  TLorentzVector getPi0sFourMomentum() const {return _pi0mom;}
866  int getNPi0() const {return _npi0;}
867 
868  void setRecoveredVertexMass(double rvtxmass){_rvtxmass = rvtxmass;}
869  void setRecoveredVertexMass() {_rvtxmass = _rvtxmom.M();}
870  void setRecoveredFourMomentum(TLorentzVector rvtxmom) {_rvtxmom = rvtxmom;}
871  void setPi0sFourMomentum(TLorentzVector pi0mom) {_pi0mom = pi0mom;}
872  void setNPi0(int npi0) {_npi0 = npi0;}
874 
875  void Print()const;
876  static int dist_sort(const Vertex* a, const Vertex* b);
877  bool covIsGood() const;
878 
879  private:
880  mutable int _id; // necessary for LCIO relation; -1 if not used, mutable for modification in ConvertVertex()
881 
882  double _chi2;
883  double _prob;
884  double _x;
885  double _y;
886  double _z;
887  double _cov[6];
888 
889  bool _isPrimary;
890 
891  vector<const lcfiplus::Track*> _tracks;
892  map<const Track*, double> _chi2Tracks;
893 
894 
895  //for AVF
896  string _vertexingname;
897 
898  //for vertex mass recovery
899  TLorentzVector _rvtxmom;
900  TLorentzVector _pi0mom;
901  double _rvtxmass;
902  int _npi0;
904 
905  ClassDefNV(Vertex, 1);
906 };
907 class MCVertex {
908  public:
909  MCVertex() : _parent(0), _recoVtx(0) {}
910  void setParent(const MCParticle* parent) {
911  _parent = parent;
912  }
913  const MCParticle* getParent() const {
914  return _parent;
915  }
916  void add(const MCParticle* mcp) {
917  _dauVec.push_back(mcp);
918  }
919  void add(const Track* trk) {
920  _recoTrks.push_back(trk);
921  }
922  void setRecoVertex(const Vertex* vtx) {
923  _recoVtx = vtx;
924  }
925  const vector<const MCParticle*>& getDaughters() const {
926  return _dauVec;
927  }
928  const vector<const Track*>& getRecoTracks() const {
929  return _recoTrks;
930  }
931  const Vertex* getRecoVertex() const {
932  return _recoVtx;
933  }
934 
935  const TVector3& getPos()const {
936  return _pos;
937  }
938  void setPos(const TVector3& pos) {
939  _pos = pos;
940  }
941 
942  private:
943  TVector3 _pos;
944  const MCParticle* _parent;
945  vector<const MCParticle*> _dauVec;
946 
947  vector<const Track*> _recoTrks;
948  const Vertex* _recoVtx;
949 
950  ClassDefNV(MCVertex, 1);
951 };
952 
953 class TrackPocaXY {
954 
955  public:
956  TrackPocaXY(const Track* trk, const Vertex* vtx);
957  TrackPocaXY(const Track* trk, const TVector3& v);
958  bool success() {
959  return _success;
960  }
961  double getFlightLength() {
962  return _flt;
963  }
964  double getPoca() {
965  return _poca;
966  }
967 
968  double operator()(double* flt) const;
969 
970  private:
971  void minimize();
972  bool _success;
973  const Track* _trk;
974  double _flt;
975  double _poca;
976  TVector3 _v;
977 
978 };
979 
980 class Jet : public TLorentzVector {
981  public:
982  // constructors
983  Jet() : _id(-1) {};
984  Jet(const Track* trk);
985  Jet(const Neutral* neutral);
986  Jet(const Vertex* vtx) : _id(-1) {
987  add(vtx);
988  }
989  Jet(const Jet& from, bool extractVertex = false);
990  ~Jet() {};
991 
992  void setId(int id)const {
993  _id = id;
994  }
995  int getId()const {
996  return _id;
997  }
998 
999  const vector<const Track*>& getTracks() const {
1000  return _tracks;
1001  }
1002  const vector<const Neutral*>& getNeutrals() const {
1003  return _neutrals;
1004  }
1005  const vector<const Vertex*>& getVertices() const {
1006  return _vertices;
1007  }
1008 
1010  vector<const Vertex*> getVerticesForFT() const;
1011 
1015  vector<const Track*> getAllTracks(bool withoutV0=false) const;
1016 
1017  // methods
1018  void add(const Jet& jet);
1019  void add(const Track* trk) {
1020  _tracks.push_back(trk);
1021  *(TLorentzVector*)this += *(TLorentzVector*)trk;
1022  }
1023  void add(const Neutral* neut) {
1024  _neutrals.push_back(neut);
1025  *(TLorentzVector*)this += *(TLorentzVector*)neut;
1026  }
1027  void add(const Vertex* vtx, bool removeTracks = true) {
1028  _vertices.push_back(vtx);
1029  for (unsigned int i=0; i<vtx->getTracks().size(); i++) {
1030  *(TLorentzVector*)this += *(TLorentzVector*)(vtx->getTracks()[i]);
1031 
1032  if (removeTracks) {
1033  vector<const Track*>::iterator it;
1034  if ((it = find(_tracks.begin(), _tracks.end(), vtx->getTracks()[i])) != _tracks.end()) {
1035  *(TLorentzVector*)this -= *(TLorentzVector*)(*it);
1036  _tracks.erase(it);
1037  }
1038  }
1039  }
1040  }
1041  //void recalculate();
1042 
1043  static int sort_by_energy(const Jet* a, const Jet* b) {
1044  return (a->E() > b->E());
1045  }
1046 
1047  double sphericity() const;
1048 
1049  // parameter contrl
1050 // void addParam(const char *paramname, Parameters &param, bool forcereset = false){
1051  void addParam(const char* paramname, Parameters& param)const {
1052  if (_params.count(paramname) == 1)throw(Exception("Jet::addParam: parameter of the specified name has been already registered."));
1053  _params[paramname] = param;
1054  }
1055 
1056  const Parameters* getParam(const char* paramname)const {
1057  if (_params.count(paramname) == 0)throw(Exception("Jet::getParam: parameter of the specified name is not registered."));
1058  return &(_params.find(paramname)->second);
1059  }
1060 
1061  const map<string, Parameters>& params()const {
1062  return _params;
1063  }
1064 
1066  TLorentzVector v;
1067  vector<const Track*> tr = getAllTracks();
1068  for (unsigned int i=0; i<tr.size(); i++)
1069  v += (*tr[i]);
1070  for (unsigned int i=0; i<_neutrals.size(); i++)
1071  v += (*_neutrals[i]);
1072 
1073  *(TLorentzVector*)this = v;
1074  }
1075 
1076  private:
1077  vector<const Track*> _tracks;
1078  vector<const Neutral*> _neutrals;
1079  vector<const Vertex*> _vertices;
1080 
1081  mutable map<string, Parameters> _params;
1082 
1083  mutable int _id;
1084 
1085  ClassDefNV(Jet, 1);
1086 };
1087 
1088 
1089 }
1090 
1091 #endif
double getProb() const
Definition: lcfiplus.h:805
void addParam(const char *paramname, Parameters &param) const
Definition: lcfiplus.h:1051
void setDefaultMCParticles(const char *name)
Definition: lcfiplus.h:315
void setBeamSizeY(double beamSizeY)
Definition: lcfiplus.h:126
Definition: lcfiplus.h:29
Definition: lcfiplus.h:25
double getRadiusOfInnermostHit() const
Definition: lcfiplus.h:479
vector< const Track * >::const_iterator TrackVecIte
Definition: lcfiplus.h:83
void assign(const char *key, T data)
Definition: lcfiplus.h:244
Definition: lcfiplus.h:27
MCParticle()
Definition: lcfiplus.h:654
Definition: lcfiplus.h:907
const vector< const Track * > & getTracks() const
Definition: lcfiplus.h:999
const MCParticle * _mcp
Definition: lcfiplus.h:760
void add(const Neutral *neut)
Definition: lcfiplus.h:1023
double getZ0() const
Definition: lcfiplus.h:422
int getVtxHits() const
Definition: lcfiplus.h:450
void setId(int id) const
Definition: lcfiplus.h:992
void setRadiusOfInnermostHit(double rimh)
Definition: lcfiplus.h:482
vector< const Jet * >::const_iterator JetVecIte
Definition: lcfiplus.h:88
void setNPi0(int npi0)
Definition: lcfiplus.h:872
Definition: lcfiplus.h:25
int getTpcHits() const
Definition: lcfiplus.h:459
const vector< const Neutral * > NeutralVec
Definition: lcfiplus.h:76
Definition: lcfiplus.h:771
void setMcp(const lcfiplus::MCParticle *mcp)
Definition: lcfiplus.h:605
~Jet()
Definition: lcfiplus.h:990
int getId() const
Definition: lcfiplus.h:799
const vector< const MCParticle * > & getDaughters() const
Definition: lcfiplus.h:925
void fetch(const char *key, T &ret, T def=T()) const
Definition: lcfiplus.h:173
TLorentzVector getRecoveredFourMomentum() const
Definition: lcfiplus.h:864
int var
Definition: lcfiplus.h:56
Definition: lcfiplus.h:29
int getEtdHits() const
Definition: lcfiplus.h:465
void setFlightLength(double flt) const
Definition: lcfiplus.h:500
void add(const Track *trk)
Definition: lcfiplus.h:1019
void setDefaultTracks(const char *name)
Definition: lcfiplus.h:309
void setV0(bool v0=true)
Definition: lcfiplus.h:626
~Vertex()
Definition: lcfiplus.h:792
const vector< const T * > & constVector(const vector< T * > &ref)
Definition: lcfiplus.h:49
void setCNess(double cness)
Definition: lcfiplus.h:549
void fetchArray(const char *key, vector< T > &ret) const
Definition: lcfiplus.h:192
vector< const Vertex * >::const_iterator VertexVecIte
Definition: lcfiplus.h:87
void setBNess(double bness)
Definition: lcfiplus.h:536
EVENT::ClusterVec getClusters() const
Definition: lcfiplus.h:631
double pull
Definition: lcfiplus.h:61
Definition: lcfiplus.h:90
Jet()
Definition: lcfiplus.h:983
const MCParticle * getParent() const
Definition: lcfiplus.h:913
const double * getCovMatrix() const
Definition: lcfiplus.h:445
Definition: lcfiplus.h:27
Parameters * GetParam() const
Definition: lcfiplus.h:282
Definition: lcfiplus.h:384
Definition: lcfiplus.h:646
int getId() const
Definition: lcfiplus.h:661
~Neutral()
Definition: lcfiplus.h:594
void setPos(const TVector3 &pos)
Definition: lcfiplus.h:938
Definition: lcfiplus.h:953
Definition: lcfiplus.h:26
Parameters(const Parameters &ref)
Definition: lcfiplus.h:164
Definition: lcfiplus.h:26
double getZ() const
Definition: lcfiplus.h:814
Exception(const char *message)
Definition: lcfiplus.h:92
Definition: lcfiplus.h:27
Definition: lcfiplus.h:25
double getCharge() const
Definition: lcfiplus.h:412
~Track()
Definition: lcfiplus.h:389
void setRecoveredFourMomentum(TLorentzVector rvtxmom)
Definition: lcfiplus.h:870
void setPDG(int pdg)
Definition: lcfiplus.h:671
int getSetHits() const
Definition: lcfiplus.h:462
void setCharge(double charge)
Definition: lcfiplus.h:678
void Print()
Definition: lcfiplus.h:101
int getSitHits() const
Definition: lcfiplus.h:456
void add(const Track *trk)
Definition: lcfiplus.h:919
Definition: lcfiplus.h:26
double getPoca()
Definition: lcfiplus.h:964
Definition: lcfiplus.h:26
Vertex()
Definition: lcfiplus.h:777
const char * getDefaultPrimaryVertex() const
Definition: lcfiplus.h:337
vector< const Neutral * >::const_iterator NeutralVecIte
Definition: lcfiplus.h:84
void setHelix(double d0, double z0, double phi, double omega, double tanLambda)
Definition: lcfiplus.h:434
virtual const char * what() const
Definition: lcfiplus.h:98
const char * getDefaultSecondaryVertices() const
Definition: lcfiplus.h:340
const lcfiplus::MCParticle * getMcp() const
Definition: lcfiplus.h:602
Definition: lcfiplus.h:29
void setPrimary(bool isPrim)
Definition: lcfiplus.h:843
Definition: lcfiplus.h:29
void add(const MCParticle *mcp)
Definition: lcfiplus.h:916
void setHelix(double *par)
Definition: lcfiplus.h:441
const map< const lcfiplus::Track *, double > & getTracksChi2Map() const
Definition: lcfiplus.h:826
Vertex(const Vertex &from)
Definition: lcfiplus.h:788
MCVertex()
Definition: lcfiplus.h:909
TVector3 getPos() const
Definition: lcfiplus.h:817
Definition: lcfiplus.h:27
Definition: lcfiplus.h:26
void setPDG(int pdg)
Definition: lcfiplus.h:611
void add(const Vertex *vtx, bool removeTracks=true)
Definition: lcfiplus.h:1027
Definition: LcfiInterface.h:49
Definition: lcfiplus.h:269
void setClusters(EVENT::ClusterVec clu)
Definition: lcfiplus.h:630
void setId(int id) const
Definition: lcfiplus.h:795
Definition: lcfiplus.h:27
vector< const lcfiplus::MCParticle * > _finalcolorsinglets
Definition: lcfiplus.h:765
void setVertexingName(string vtxname)
Definition: lcfiplus.h:858
double cmax
Definition: lcfiplus.h:58
int getId() const
Definition: lcfiplus.h:995
bool isV0() const
Definition: lcfiplus.h:623
const char * getDefaultNeutrals() const
Definition: lcfiplus.h:331
void setMcp(const lcfiplus::MCParticle *mcp)
Definition: lcfiplus.h:401
Jet(const Vertex *vtx)
Definition: lcfiplus.h:986
Definition: lcfiplus.h:152
const map< string, pair< const type_info *, void * > > & paramMap() const
Definition: lcfiplus.h:253
vector< const lcfiplus::MCParticle * > _finalstrings
Definition: lcfiplus.h:763
double getBeamSizeX() const
Definition: lcfiplus.h:123
void setPi0sFourMomentum(TLorentzVector pi0mom)
Definition: lcfiplus.h:871
void recalcFourMomentum()
Definition: lcfiplus.h:1065
void setBeamSizeX(double beamSizeX)
Definition: lcfiplus.h:120
const vector< const MCParticle * > & getDaughters() const
Definition: lcfiplus.h:698
Definition: lcfiplus.h:33
const TVector3 & getPos() const
Definition: lcfiplus.h:935
Definition: lcfiplus.h:29
Definition: lcfiplus.h:27
const lcfiplus::MCParticle * getMcp() const
Definition: lcfiplus.h:398
Definition: lcfiplus.h:29
Definition: lcfiplus.h:26
virtual ~Algorithm()
Definition: lcfiplus.h:273
Vertex(const double chi2, const double prob, const double x, const double y, const double z, const double cov[6], bool isPrim)
Definition: lcfiplus.h:778
TLorentzVector getPi0sFourMomentum() const
Definition: lcfiplus.h:865
vector< const lcfiplus::MCParticle * > _initials
Definition: lcfiplus.h:762
Definition: lcfiplus.h:29
Definition: lcfiplus.h:30
double getD0() const
Definition: lcfiplus.h:419
std::vector< T > getVec(const char *key) const
Definition: lcfiplus.h:223
const map< string, Parameters > & params() const
Definition: lcfiplus.h:1061
double getTanLambda() const
Definition: lcfiplus.h:431
A simple named storage for event data.
Definition: EventStore.h:42
int getPDG() const
Definition: lcfiplus.h:405
const vector< const Neutral * > & getNeutrals() const
Definition: lcfiplus.h:1002
double getPhi() const
Definition: lcfiplus.h:425
double getChi2Track(const Track *tr) const
Definition: lcfiplus.h:829
Algorithm()
Definition: lcfiplus.h:271
void setNdf(int ndf)
Definition: lcfiplus.h:496
double cmin
Definition: lcfiplus.h:57
void setRecoVertex(const Vertex *vtx)
Definition: lcfiplus.h:922
void setCharge(double charge)
Definition: lcfiplus.h:415
void setRecoveredVertexMass(double rvtxmass)
Definition: lcfiplus.h:868
const double * getCaloEdep() const
Definition: lcfiplus.h:472
MCParticle(int id, int pdg, MCParticle *parent, double charge, const TLorentzVector &p, const TVector3 &v)
Definition: lcfiplus.h:651
~Parameters()
Definition: lcfiplus.h:160
void setRecoveredVertexMass()
Definition: lcfiplus.h:869
const vector< const Track * > & getRecoTracks() const
Definition: lcfiplus.h:928
int getPDG() const
Definition: lcfiplus.h:608
void setPDG(int pdg)
Definition: lcfiplus.h:408
void setParent(const MCParticle *parent)
Definition: lcfiplus.h:693
vector< const lcfiplus::MCParticle * > _realparticles
Definition: lcfiplus.h:766
Definition: lcfiplus.h:30
const double * getCaloEdep() const
Definition: lcfiplus.h:615
Definition: lcfiplus.h:754
Definition: lcfiplus.h:25
const vector< const Vertex * > & getVertices() const
Definition: lcfiplus.h:1005
void setId(int id)
Definition: lcfiplus.h:394
Definition: lcfiplus.h:30
double getRecoveredVertexMass() const
Definition: lcfiplus.h:863
double getCorrEnergy() const
Definition: lcfiplus.h:525
int getNPi0() const
Definition: lcfiplus.h:866
int getFtdHits() const
Definition: lcfiplus.h:453
const Parameters * getParam(const char *paramname) const
Definition: lcfiplus.h:1056
Definition: lcfiplus.h:980
vector< const lcfiplus::MCParticle * > _qqgluons
Definition: lcfiplus.h:764
double getFlightLength()
Definition: lcfiplus.h:961
void setCaloEdep(double *calo)
Definition: lcfiplus.h:619
double getBeamSizeY() const
Definition: lcfiplus.h:129
Parameters(bool as=true)
Definition: lcfiplus.h:154
calo
Definition: lcfiplus.h:30
const vector< const Jet * > JetVec
Definition: lcfiplus.h:80
int getPDG() const
Definition: lcfiplus.h:668
void setVertex(const TVector3 &v)
Definition: lcfiplus.h:685
const vector< const Track * > TrackVec
Definition: lcfiplus.h:72
void swapEnergy()
Definition: lcfiplus.h:529
const TVector3 & getVertex() const
Definition: lcfiplus.h:682
bool success()
Definition: lcfiplus.h:958
void setParent(const MCParticle *parent)
Definition: lcfiplus.h:910
Definition: lcfiplus.h:291
double getCNess() const
Definition: lcfiplus.h:544
double pullerr
Definition: lcfiplus.h:62
const char * getDefaultJets() const
Definition: lcfiplus.h:343
Definition: lcfiplus.h:55
vector< const MCParticle * >::const_iterator MCParticleVecIte
Definition: lcfiplus.h:85
Definition: lcfiplus.h:110
void setBField(double bField)
Definition: lcfiplus.h:114
virtual void end()
Definition: lcfiplus.h:279
int getNdf() const
Definition: lcfiplus.h:493
const char * getDefaultMCParticles() const
Definition: lcfiplus.h:334
void setParticleIDProbability(string parName, double pidProbability)
Definition: lcfiplus.h:507
virtual void init(Parameters *param)
Definition: lcfiplus.h:275
const char * getDefaultTracks() const
Definition: lcfiplus.h:328
vector< const MCColorSinglet * >::const_iterator MCColorSingletVecIte
Definition: lcfiplus.h:86
string getVertexingName() const
Definition: lcfiplus.h:859
double getChi2() const
Definition: lcfiplus.h:486
void setId(int id)
Definition: lcfiplus.h:599
void setDefaultPrimaryVertex(const char *name)
Definition: lcfiplus.h:318
void setDefaultNeutrals(const char *name)
Definition: lcfiplus.h:312
bool exist(const char *key) const
Definition: lcfiplus.h:229
cov
Definition: lcfiplus.h:26
double getFlightLength() const
Definition: lcfiplus.h:503
const double * getCov() const
Definition: lcfiplus.h:820
double getParticleIDProbability(const char *parName) const
Definition: lcfiplus.h:511
Neutral()
Definition: lcfiplus.h:593
double getCharge() const
Definition: lcfiplus.h:675
void setDefaultSecondaryVertices(const char *name)
Definition: lcfiplus.h:321
void setCorrEnergy(double mass)
Definition: lcfiplus.h:521
double getY() const
Definition: lcfiplus.h:811
virtual ~Exception()
Definition: lcfiplus.h:96
static int sort_by_energy(const Jet *a, const Jet *b)
Definition: lcfiplus.h:1043
void setTrackHits(int *hits)
Definition: lcfiplus.h:468
int getId() const
Definition: lcfiplus.h:596
vtx
Definition: lcfiplus.h:775
~MCParticle()
Definition: lcfiplus.h:655
double getOmega() const
Definition: lcfiplus.h:428
const vector< const MCParticle * > MCParticleVec
Definition: lcfiplus.h:77
Definition: lcfiplus.h:25
double getBField() const
Definition: lcfiplus.h:117
double getBeamSizeZ() const
Definition: lcfiplus.h:135
double pmin
Definition: lcfiplus.h:59
Definition: lcfiplus.h:589
Definition: lcfiplus.h:30
Definition: lcfiplus.h:30
Definition: lcfiplus.h:27
int getId() const
Definition: lcfiplus.h:391
void setChi2(double chi2)
Definition: lcfiplus.h:489
double getX() const
Definition: lcfiplus.h:808
const MCParticle * getMcp() const
Definition: lcfiplus.h:756
const vector< const Track * > & getTracks() const
Definition: lcfiplus.h:823
Definition: lcfiplus.h:27
Definition: lcfiplus.h:27
TVector3 getPos() const
Definition: lcfiplus.h:555
double pmax
Definition: lcfiplus.h:60
double getBNess() const
Definition: lcfiplus.h:540
Definition: lcfiplus.h:26
void add(const char *key, T data)
Definition: lcfiplus.h:233
hit
Definition: lcfiplus.h:29
Definition: lcfiplus.h:30
par
Definition: lcfiplus.h:25
const vector< const Vertex * > VertexVec
Definition: lcfiplus.h:79
const MCParticle * getParent() const
Definition: lcfiplus.h:690
Definition: lcfiplus.h:25
double getChi2() const
Definition: lcfiplus.h:802
Definition: lcfiplus.h:30
bool isPrimary() const
Definition: lcfiplus.h:840
void setDefaultJets(const char *name)
Definition: lcfiplus.h:324
const vector< const MCColorSinglet * > MCColorSingletVec
Definition: lcfiplus.h:78
const Vertex * getRecoVertex() const
Definition: lcfiplus.h:931
void setId(int id)
Definition: lcfiplus.h:664
void setBeamSizeZ(double beamSizeZ)
Definition: lcfiplus.h:132
void setCaloEdep(double *calo)
Definition: lcfiplus.h:475