CbmRoot
CbmTrdSetTracksPidWkn.cxx
Go to the documentation of this file.
1 
2 // -------------------------------------------------------------------------
3 // ----- CbmTrdSetTracksPidWkn source file -----
4 // ----- Created 13/02/07 by F.Uhlig -----
5 //------ Last modification 01/07/18 by O.Derenovskaya
6 // -------------------------------------------------------------------------
8 
9 #include "CbmTrdHit.h"
10 #include "CbmTrdTrack.h"
11 
12 #include "FairRootManager.h"
13 
14 #include "TClonesArray.h"
15 #include "TMath.h"
16 
17 #include <iostream>
18 #include <vector>
19 
20 using std::cout;
21 using std::endl;
22 
23 
24 // ----- Default constructor -------------------------------------------
26  : CbmTrdSetTracksPidWkn("TrdSetTracksPidWkn", "") {}
27 // -------------------------------------------------------------------------
28 
29 
30 // ----- Standard constructor ------------------------------------------
31 CbmTrdSetTracksPidWkn::CbmTrdSetTracksPidWkn(const char* name, const char*)
32  : FairTask(name)
33  , fnSet(0)
34  , fdegWkn(0)
35  , fNofTracks(0)
36  , fk1(0)
37  , fwkn0(0)
38  , fEmp(0)
39  , fXi(0)
40  , fWmin(0)
41  , fWmax(0)
42  , fDiff(0)
43  , fSISType("sis100")
44  , fTrackArray(NULL)
45  , fTrdHitArray(NULL)
46 
47 
48 {}
49 // -------------------------------------------------------------------------
50 
51 
52 // ----- Destructor ----------------------------------------------------
54 // -------------------------------------------------------------------------
55 
56 
57 // ----- SetParContainers -------------------------------------------------
59 // -------------------------------------------------------------------------
60 
61 
62 // ----- Public method Init (abstract in base class) --------------------
64 
65  // Get and check FairRootManager
66  FairRootManager* ioman = FairRootManager::Instance();
67  if (!ioman) {
68  cout << "-E- CbmTrdSetTracksPidWkn::Init: "
69  << "RootManager not instantised!" << endl;
70  return kFATAL;
71  }
72 
73  // Get TrdTrack array
74  fTrackArray = (TClonesArray*) ioman->GetObject("TrdTrack");
75  if (!fTrackArray) {
76  cout << "-E- CbmTrdSetTracksPidWkn::Init: No TrdTrack array!" << endl;
77  return kERROR;
78  }
79 
80  // Get TrdHit array
81  fTrdHitArray = (TClonesArray*) ioman->GetObject("TrdHit");
82  if (!fTrdHitArray) {
83  cout << "-E- CbmTrdSetTracksPidWkn::Init: No TrdHit array!" << endl;
84  return kERROR;
85  }
86 
87  SetParameters();
88 
89  return kSUCCESS;
90 }
91 // -------------------------------------------------------------------------
92 
93 
94 // ----- Public method Exec --------------------------------------------
95 void CbmTrdSetTracksPidWkn::Exec(Option_t*) {
96 
97  if (!fTrackArray) return;
98 
99  Int_t nTracks = fTrackArray->GetEntriesFast();
100 
101  Double_t result_wkn;
102  Int_t NHits;
103 
104  fNofTracks = 0;
105 
106 
107  std::vector<float> eLossVectorTmp;
108  std::vector<float> eLossVector;
109 
110  for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
111  eLossVectorTmp.clear();
112  eLossVector.clear();
113  CbmTrdTrack* pTrack = (CbmTrdTrack*) fTrackArray->At(iTrack);
114  NHits = pTrack->GetNofHits();
115 
116  // Set total energy loss member in CbmTrdTack object if not already done
117  assert(pTrack);
118  if (pTrack->GetELoss() < 0.) {
119  double_t sumEloss = 0.;
120  for (Int_t iHit = 0; iHit < NHits; iHit++) {
121  CbmTrdHit* hit =
122  (CbmTrdHit*) fTrdHitArray->At(pTrack->GetHitIndex(iHit));
123  assert(hit);
124  sumEloss += hit->GetELoss();
125  }
126  pTrack->SetELoss(sumEloss);
127  } //? Track energy loss not set
128 
129  // Up to now only for tracks with twelve hits the Wkn can be calculated
130  // This should be changed in the future.
131  if (NHits < fnSet) {
132  fNofTracks++;
133  pTrack->SetPidWkn(-2.);
134  continue;
135  }
136 
137  for (Int_t iTRD = 0; iTRD < NHits; iTRD++) {
138  Int_t TRDindex = pTrack->GetHitIndex(iTRD);
139  CbmTrdHit* trdHit = (CbmTrdHit*) fTrdHitArray->At(TRDindex);
140  eLossVectorTmp.push_back((trdHit->GetELoss()) * 1000000);
141  }
142 
143  // calculate the lambda parameter for each TRD layer
144  for (unsigned int jVec = 0; jVec < eLossVectorTmp.size(); jVec++)
145  eLossVectorTmp[jVec] = (eLossVectorTmp[jVec] - fEmp) / fXi - 0.225;
146 
147  sort(eLossVectorTmp.begin(), eLossVectorTmp.end());
148 
149  for (unsigned int jVec = 0; jVec < eLossVectorTmp.size(); jVec++) {
150  eLossVectorTmp[jVec] = TMath::LandauI(eLossVectorTmp[jVec]);
151  }
152 
153  for (int iHit = 0; iHit < fnSet; iHit++)
154  eLossVector.push_back(eLossVectorTmp[NHits - fnSet + iHit]);
155 
156  // calculate the Wkn and add the information to the TrdTrack
157 
158  Double_t S = 0, ty = 0, ti = 0;
159 
160  for (Int_t i = 0; i < fnSet; i++) {
161  ty = eLossVector[i];
162  ti = i;
163  S += pow((ti - 1) / fnSet - ty, fk1) - pow(ti / fnSet - ty, fk1);
164  }
165  result_wkn = -fwkn0 * S;
166  // cout<<result_wkn<<endl;
167  pTrack->SetPidWkn(result_wkn);
168  }
169 }
170 // -------------------------------------------------------------------------
171 
172 // ----- Public method Finish ------------------------------------------
174 // -------------------------------------------------------------------------
175 
177  if (fSISType == "sis300") {
178  fnSet = 5; // number of the layers with TR
179  fdegWkn = 4; // statistics degree
180  fEmp = 1.06;
181  fXi = 0.57;
182  }
183 
184  if (fSISType == "sis100") {
185  fnSet = 3; // number of the layers with TR
186  fdegWkn = 2; // statistics degree
187  fEmp = 3.5;
188  fXi = 5.0;
189  }
190 
191 
192  fk1 = fdegWkn + 1;
193  fwkn0 = pow(fnSet, 0.5 * fdegWkn) / fk1;
194  fWmin = 1 / (pow(2, fdegWkn) * pow(fnSet, fdegWkn / 2) * (fdegWkn + 1));
195  fWmax = pow(fnSet, fdegWkn / 2) / (fdegWkn + 1);
196  fDiff = fWmax - fWmin;
197 }
198 
CbmTrack::GetNofHits
virtual Int_t GetNofHits() const
Definition: CbmTrack.h:53
CbmTrdSetTracksPidWkn.h
CbmTrdSetTracksPidWkn::fnSet
int fnSet
Definition: CbmTrdSetTracksPidWkn.h:77
CbmTrdSetTracksPidWkn::fWmin
float fWmin
Definition: CbmTrdSetTracksPidWkn.h:83
CbmTrdHit::GetELoss
Double_t GetELoss() const
Definition: CbmTrdHit.h:79
CbmTrdSetTracksPidWkn::Exec
virtual void Exec(Option_t *opt)
Definition: CbmTrdSetTracksPidWkn.cxx:95
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmTrdHit
data class for a reconstructed Energy-4D measurement in the TRD
Definition: CbmTrdHit.h:35
CbmTrdSetTracksPidWkn::fwkn0
float fwkn0
Definition: CbmTrdSetTracksPidWkn.h:80
CbmTrdSetTracksPidWkn::~CbmTrdSetTracksPidWkn
virtual ~CbmTrdSetTracksPidWkn()
Definition: CbmTrdSetTracksPidWkn.cxx:53
CbmTrdSetTracksPidWkn::fNofTracks
int fNofTracks
Definition: CbmTrdSetTracksPidWkn.h:79
CbmTrdSetTracksPidWkn::fEmp
float fEmp
Definition: CbmTrdSetTracksPidWkn.h:81
CbmTrdTrack::GetELoss
Double_t GetELoss() const
Definition: CbmTrdTrack.h:33
CbmTrdSetTracksPidWkn::Finish
virtual void Finish()
Definition: CbmTrdSetTracksPidWkn.cxx:173
CbmTrack::GetHitIndex
Int_t GetHitIndex(Int_t iHit) const
Definition: CbmTrack.h:54
CbmTrdSetTracksPidWkn::fTrdHitArray
TClonesArray * fTrdHitArray
Definition: CbmTrdSetTracksPidWkn.h:93
CbmTrdSetTracksPidWkn::fSISType
std::string fSISType
Definition: CbmTrdSetTracksPidWkn.h:90
CbmTrdHit.h
Class for hits in TRD detector.
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmTrdSetTracksPidWkn::fXi
float fXi
Definition: CbmTrdSetTracksPidWkn.h:82
CbmTrdSetTracksPidWkn::fTrackArray
TClonesArray * fTrackArray
Definition: CbmTrdSetTracksPidWkn.h:92
CbmTrdSetTracksPidWkn::fWmax
float fWmax
Definition: CbmTrdSetTracksPidWkn.h:83
CbmTrdTrack
Definition: CbmTrdTrack.h:22
CbmTrdTrack::SetPidWkn
void SetPidWkn(Double_t pid)
Definition: CbmTrdTrack.h:41
CbmTrdTrack::SetELoss
void SetELoss(Double_t eLoss)
Definition: CbmTrdTrack.h:43
CbmTrdSetTracksPidWkn::SetParContainers
virtual void SetParContainers()
Definition: CbmTrdSetTracksPidWkn.cxx:58
CbmTrdSetTracksPidWkn::fDiff
float fDiff
Definition: CbmTrdSetTracksPidWkn.h:83
CbmTrdTrack.h
CbmTrdSetTracksPidWkn::fk1
float fk1
Definition: CbmTrdSetTracksPidWkn.h:80
CbmTrdSetTracksPidWkn::CbmTrdSetTracksPidWkn
CbmTrdSetTracksPidWkn()
Definition: CbmTrdSetTracksPidWkn.cxx:25
CbmTrdSetTracksPidWkn::Init
virtual InitStatus Init()
Definition: CbmTrdSetTracksPidWkn.cxx:63
CbmTrdSetTracksPidWkn::fdegWkn
int fdegWkn
Definition: CbmTrdSetTracksPidWkn.h:78
CbmTrdSetTracksPidWkn::SetParameters
void SetParameters()
Definition: CbmTrdSetTracksPidWkn.cxx:176
CbmTrdSetTracksPidWkn
Definition: CbmTrdSetTracksPidWkn.h:36