CbmRoot
CbmAnaConversionGlobalFunctions.h
Go to the documentation of this file.
1 
6 #ifndef CBM_ANA_CONVERSION_GLOBAL_FUNCTIONS
7 #define CBM_ANA_CONVERSION_GLOBAL_FUNCTIONS
8 
9 #include "CbmGlobalTrack.h"
10 #include "CbmMCTrack.h"
11 #include "CbmRichElectronIdAnn.h"
12 #include "CbmRichRing.h"
13 
14 #include "FairRootManager.h"
15 
16 #include "TClonesArray.h"
17 #include "TLorentzVector.h"
18 #include "TMath.h"
19 
20 #define M2E 2.6112004954086e-7
21 
23 public:
24  /*
25  * Calculate ANN response for a given globaltrack index and momentum
26  */
27  static Double_t ElectronANNvalue(Int_t globalTrackIndex, Double_t momentum) {
28  FairRootManager* ioman = FairRootManager::Instance();
29  if (NULL == ioman) return -2;
30  TClonesArray* fGlobalTracks =
31  (TClonesArray*) ioman->GetObject("GlobalTrack");
32  if (NULL == fGlobalTracks) return -2;
33  TClonesArray* fRichRings = (TClonesArray*) ioman->GetObject("RichRing");
34  if (NULL == fRichRings) return -2;
35 
36  //if (NULL == fGlobalTracks || NULL == fRichRings) return -2;
37  //CbmGlobalTrack* globalTrack = (CbmGlobalTrack*) fGlobalTracks->At(globalTrackIndex);
38  const CbmGlobalTrack* globalTrack =
39  static_cast<const CbmGlobalTrack*>(fGlobalTracks->At(globalTrackIndex));
40  Int_t richId = globalTrack->GetRichRingIndex();
41  if (richId < 0) return -2;
42  CbmRichRing* ring = static_cast<CbmRichRing*>(fRichRings->At(richId));
43  if (NULL == ring) return -2;
44 
46  globalTrackIndex, momentum);
47  return ann;
48  }
49 
50  /*
51  * Checks, whether a track is an electron based on ANN output for a given globaltrack index and momentum, and for a given ANNcut value
52  */
53  static Double_t IsRICHElectronANN(Int_t globalTrackIndex,
54  Double_t momentum,
55  Double_t ANNcut) {
56  FairRootManager* ioman = FairRootManager::Instance();
57  if (NULL == ioman) return -2;
58  TClonesArray* fGlobalTracks =
59  (TClonesArray*) ioman->GetObject("GlobalTrack");
60  if (NULL == fGlobalTracks) return -2;
61  TClonesArray* fRichRings = (TClonesArray*) ioman->GetObject("RichRing");
62  if (NULL == fRichRings) return -2;
63 
64 
65  //if (NULL == fGlobalTracks || NULL == fRichRings) return -2;
66  //CbmGlobalTrack* globalTrack = (CbmGlobalTrack*) fGlobalTracks->At(globalTrackIndex);
67  const CbmGlobalTrack* globalTrack =
68  static_cast<const CbmGlobalTrack*>(fGlobalTracks->At(globalTrackIndex));
69  Int_t richId = globalTrack->GetRichRingIndex();
70  if (richId < 0) return -2;
71  CbmRichRing* ring = static_cast<CbmRichRing*>(fRichRings->At(richId));
72  if (NULL == ring) return -2;
73 
75  globalTrackIndex, momentum);
76  if (ann > ANNcut)
77  return true;
78  else
79  return false;
80  }
81 
82 
83  static Double_t OpeningAngleBetweenGamma(const TVector3 part11,
84  const TVector3 part12,
85  const TVector3 part21,
86  const TVector3 part22) {
87  Double_t openingAngle = 0;
88 
89  Double_t energy11 = TMath::Sqrt(part11.Mag2() + M2E);
90  TLorentzVector lorVec11(part11, energy11);
91 
92  Double_t energy12 = TMath::Sqrt(part12.Mag2() + M2E);
93  TLorentzVector lorVec12(part12, energy12);
94 
95  Double_t energy21 = TMath::Sqrt(part21.Mag2() + M2E);
96  TLorentzVector lorVec21(part21, energy21);
97 
98  Double_t energy22 = TMath::Sqrt(part22.Mag2() + M2E);
99  TLorentzVector lorVec22(part22, energy22);
100 
101  TLorentzVector gamma1 = lorVec11 + lorVec12;
102  TLorentzVector gamma2 = lorVec21 + lorVec22;
103 
104  Double_t angle = gamma1.Angle(gamma2.Vect());
105  openingAngle = 180. * angle / TMath::Pi();
106 
107  return openingAngle;
108  }
109 };
110 
111 #endif
CbmAnaConversionGlobalFunctions::IsRICHElectronANN
static Double_t IsRICHElectronANN(Int_t globalTrackIndex, Double_t momentum, Double_t ANNcut)
Definition: CbmAnaConversionGlobalFunctions.h:53
CbmGlobalTrack::GetRichRingIndex
Int_t GetRichRingIndex() const
Definition: CbmGlobalTrack.h:41
CbmAnaConversionGlobalFunctions::OpeningAngleBetweenGamma
static Double_t OpeningAngleBetweenGamma(const TVector3 part11, const TVector3 part12, const TVector3 part21, const TVector3 part22)
Definition: CbmAnaConversionGlobalFunctions.h:83
CbmGlobalTrack.h
CbmRichRing
Definition: CbmRichRing.h:17
CbmRichRing.h
CbmAnaConversionGlobalFunctions
Definition: CbmAnaConversionGlobalFunctions.h:22
xMath::Pi
double Pi()
Definition: xMath.h:5
M2E
#define M2E
Definition: CbmAnaConversionGlobalFunctions.h:20
CbmGlobalTrack
Definition: CbmGlobalTrack.h:26
CbmRichElectronIdAnn::GetInstance
static CbmRichElectronIdAnn & GetInstance()
Definition: CbmRichElectronIdAnn.h:43
CbmMCTrack.h
CbmRichElectronIdAnn::CalculateAnnValue
double CalculateAnnValue(int globalTrackIndex, double momentum)
Calculate output value of the ANN.
Definition: CbmRichElectronIdAnn.cxx:84
CbmAnaConversionGlobalFunctions::ElectronANNvalue
static Double_t ElectronANNvalue(Int_t globalTrackIndex, Double_t momentum)
Definition: CbmAnaConversionGlobalFunctions.h:27
CbmRichElectronIdAnn.h
Implementation of the electron identification algorithm in the RICH detector using Artificial Neural ...