remove obselete clean QA macros command
[u/mrichter/AliRoot.git] / PWG1 / AliMCInfo.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16
17 ///////////////////////////////////////////////////////////////////////////
18 /*
19
20 Origin: marian.ivanov@cern.ch
21 Container classes with MC infomation
22
23 The AliMCInfo contains the information about the particles properties 
24 during transportation throuch ALICE Detector
25
26 The base Information :
27 TParticle - fParticle             -  properties of the particle at creation point
28 AliTrackReference - fXXXRefernces -  TClonesArray of refernces in differnt detectors
29 fNXXXRef                          -  number of the track refernces in differnt detectors
30 AliTPCdigitRow    - fTPCRow       -  the map of the hitted rows - (will be repalced by TBits)
31 fRowsWith*                        - number of rows hitted by particle
32 fMCtracks         -               - number of turn over of the track inside of the TPC
33
34 ++++
35 some additional information usable for tree draw - TO SPEED UP tree queries
36 IMPORTANT FOR PROOF FAST PROTOTYPING ANALYSIS 
37                                       
38
39
40
41 */
42
43 #if !defined(__CINT__) || defined(__MAKECINT__)
44 #include <stdio.h>
45 //ROOT includes
46 #include "Rtypes.h"
47 #include "TClonesArray.h"
48 //ALIROOT includes
49 #include "AliTrackReference.h"
50 #include "AliMCInfo.h" 
51 #include "AliMathBase.h" 
52 #endif
53
54 //
55 // 
56
57 ClassImp(AliTPCdigitRow)
58 ClassImp(AliMCInfo)
59
60
61
62 ////////////////////////////////////////////////////////////////////////
63 AliMCInfo::AliMCInfo():
64   fTrackRef(),
65   fTrackRefOut(),
66   fTRdecay(),
67   fPrimPart(0),
68   fParticle(),
69   fMass(0),
70   fCharge(0),
71   fLabel(0),
72   fEventNr(),
73   fMCtracks(),
74   fPdg(0),
75   fTPCdecay(0),
76   fRowsWithDigitsInn(0),
77   fRowsWithDigits(0),
78   fRowsTrackLength(0),
79   fTPCtrackLength(-1),
80   fPrim(0),
81   fTPCRow(), 
82   fNTPCRef(0),                    // tpc references counter
83   fNITSRef(0),                    // ITS references counter
84   fNTRDRef(0),                    // TRD references counter
85   fNTOFRef(0),                    // TOF references counter
86   fNTPCRefOut(0),                    // tpc references counter
87   fNITSRefOut(0),                    // ITS references counter
88   fNTRDRefOut(0),                    // TRD references counter
89   fNTOFRefOut(0),                    // TOF references counter
90   fTPCReferences(0),
91   fITSReferences(0),
92   fTRDReferences(0),
93   fTOFReferences(0)
94 {
95   //
96   // Default constructor
97   //
98   fTPCReferences  = new TClonesArray("AliTrackReference",10);
99   fITSReferences  = new TClonesArray("AliTrackReference",10);
100   fTRDReferences  = new TClonesArray("AliTrackReference",10);
101   fTOFReferences  = new TClonesArray("AliTrackReference",10);
102   fTRdecay.SetTrack(-1);
103   for (Int_t i=0;i<4;i++) fVDist[i]=0;
104   for (Int_t i=0;i<3;i++) fDecayCoord[i]=0;
105   fCharge = 0;
106 }
107
108 AliMCInfo::AliMCInfo(const AliMCInfo& info):
109   TObject(),
110   fTrackRef(info.fTrackRef),
111   fTrackRefOut(info.fTrackRefOut),
112   fTRdecay(info.GetTRdecay()),
113   fPrimPart(info.fPrimPart),
114   fParticle(info.fParticle),
115   fMass(info.GetMass()),
116   fCharge(info.fCharge),
117   fLabel(info.fLabel),
118   fEventNr(info.fEventNr),
119   fMCtracks(info.fMCtracks),
120   fPdg(info.fPdg),
121   fTPCdecay(info.fTPCdecay),
122   fRowsWithDigitsInn(info.fRowsWithDigitsInn),
123   fRowsWithDigits(info.fRowsWithDigits),
124   fRowsTrackLength(info.fRowsTrackLength),
125   fTPCtrackLength(info.fTPCtrackLength),
126   fPrim(info.fPrim),
127   fTPCRow(info.fTPCRow), 
128   fNTPCRef(info.fNTPCRef),                    // tpc references counter
129   fNITSRef(info.fNITSRef),                    // ITS references counter
130   fNTRDRef(info.fNTRDRef),                    // TRD references counter
131   fNTOFRef(info.fNTOFRef),                    // TOF references counter
132   fNTPCRefOut(info.fNTPCRefOut),                    // tpc references counter
133   fNITSRefOut(info.fNITSRefOut),                    // ITS references counter
134   fNTRDRefOut(info.fNTRDRefOut),                    // TRD references counter
135   fNTOFRefOut(info.fNTOFRefOut),                    // TOF references counter
136   fTPCReferences(0),
137   fITSReferences(0),
138   fTRDReferences(0),
139   fTOFReferences(0)
140 {
141   //
142   // copy constructor
143   //
144   fTPCReferences = (TClonesArray*)info.fTPCReferences->Clone();
145   fITSReferences = (TClonesArray*)info.fITSReferences->Clone();
146   fTRDReferences = (TClonesArray*)info.fTRDReferences->Clone();
147   fTOFReferences = (TClonesArray*)info.fTOFReferences->Clone();
148   for (Int_t i=0;i<4;i++) fVDist[i]=info.fVDist[i];
149   for (Int_t i=0;i<3;i++) fDecayCoord[i]=info.fDecayCoord[i];
150
151 }
152
153
154 AliMCInfo& AliMCInfo::operator=(const AliMCInfo& info) { 
155   //
156   // Assignment operator
157   //
158   this->~AliMCInfo();
159   new (this) AliMCInfo(info);
160   return *this;
161 }
162
163
164 AliMCInfo::~AliMCInfo()
165 {
166   //
167   // Destructor of the class
168   //
169   if (fTPCReferences) {
170     delete fTPCReferences;
171   }
172   if (fITSReferences){
173     delete fITSReferences;
174   }
175   if (fTRDReferences){    
176     delete fTRDReferences;  
177   }
178   if (fTOFReferences){    
179     delete fTOFReferences;  
180   }
181
182 }
183
184
185
186
187 void AliMCInfo::Update()
188 {
189   //
190   // Update MC info
191   // Calculates some derived variables
192   //
193   //
194   fMCtracks =1;
195   //
196   // decay info
197   fTPCdecay=kFALSE;
198   if (fTRdecay.GetTrack()>0){
199     fDecayCoord[0] = fTRdecay.X();
200     fDecayCoord[1] = fTRdecay.Y();
201     fDecayCoord[2] = fTRdecay.Z();
202     if ( (fTRdecay.R()<250)&&(fTRdecay.R()>85) && (TMath::Abs(fTRdecay.Z())<250) ){
203       fTPCdecay=kTRUE;     
204     }
205     else{
206       fDecayCoord[0] = 0;
207       fDecayCoord[1] = 0;
208       fDecayCoord[2] = 0;
209     }
210   }
211   //
212   //
213   //digits information update
214   fRowsWithDigits    = fTPCRow.RowsOn();    
215   fRowsWithDigitsInn = fTPCRow.RowsOn(63); // 63 = number of inner rows
216   fRowsTrackLength   = fTPCRow.Last() - fTPCRow.First();
217   //
218   //
219   // calculate primary ionization per cm  
220   if (fParticle.GetPDG()){
221     fMass = fParticle.GetMass();  
222     fCharge = fParticle.GetPDG()->Charge();
223     if (fTPCReferences->GetEntriesFast()>0){
224       fTrackRef = *((AliTrackReference*)fTPCReferences->At(0));
225     }
226     if (fMass>0){
227       Float_t p = TMath::Sqrt(fTrackRef.Px()*fTrackRef.Px()+
228                               fTrackRef.Py()*fTrackRef.Py()+
229                               fTrackRef.Pz()*fTrackRef.Pz());    
230       if (p>0.001){
231         Float_t betagama = p /fMass;
232         fPrim = AliMathBase::BetheBlochAleph(betagama);
233       }else fPrim=0;
234     }
235   }else{
236     fMass =0;
237     fPrim =0;
238   }  
239   
240
241
242 }
243
244
245
246
247 ////////////////////////////////////////////////////////////////////////
248 AliTPCdigitRow::AliTPCdigitRow()
249 {
250   Reset();
251 }
252 ////////////////////////////////////////////////////////////////////////
253 AliTPCdigitRow & AliTPCdigitRow::operator=(const AliTPCdigitRow &digOld)
254 {
255   for (Int_t i = 0; i<32; i++) fDig[i] = digOld.fDig[i];
256   return (*this);
257 }
258 ////////////////////////////////////////////////////////////////////////
259 void AliTPCdigitRow::SetRow(Int_t row) 
260 {
261   //
262   // set bit mask for given row
263   //
264   if (row >= 8*32) {
265     //    cerr<<"AliTPCdigitRow::SetRow: index "<<row<<" out of bounds."<<endl;
266     return;
267   }
268   Int_t iC = row/8;
269   Int_t iB = row%8;
270   SETBIT(fDig[iC],iB);
271 }
272
273 ////////////////////////////////////////////////////////////////////////
274 Bool_t AliTPCdigitRow::TestRow(Int_t row) const
275 {
276 //
277 // return kTRUE if row is on
278 //
279   Int_t iC = row/8;
280   Int_t iB = row%8;
281   return TESTBIT(fDig[iC],iB);
282 }
283 ////////////////////////////////////////////////////////////////////////
284 Int_t AliTPCdigitRow::RowsOn(Int_t upto) const
285 {
286 //
287 // returns number of rows with a digit  
288 // count only rows less equal row number upto
289 //
290   Int_t total = 0;
291   for (Int_t i = 0; i<32; i++) {
292     for (Int_t j = 0; j < 8; j++) {
293       if (i*8+j > upto) return total;
294       if (TESTBIT(fDig[i],j))  total++;
295     }
296   }
297   return total;
298 }
299 ////////////////////////////////////////////////////////////////////////
300 void AliTPCdigitRow::Reset()
301 {
302 //
303 // resets all rows to zero
304 //
305   for (Int_t i = 0; i<32; i++) {
306     fDig[i] <<= 8;
307   }
308 }
309 ////////////////////////////////////////////////////////////////////////
310 Int_t AliTPCdigitRow::Last() const
311 {
312 //
313 // returns the last row number with a digit
314 // returns -1 if now digits 
315 //
316   for (Int_t i = 32-1; i>=0; i--) {
317     for (Int_t j = 7; j >= 0; j--) {
318       if TESTBIT(fDig[i],j) return i*8+j;
319     }
320   }
321   return -1;
322 }
323 ////////////////////////////////////////////////////////////////////////
324 Int_t AliTPCdigitRow::First() const
325 {
326 //
327 // returns the first row number with a digit
328 // returns -1 if now digits 
329 //
330   for (Int_t i = 0; i<32; i++) {
331     for (Int_t j = 0; j < 8; j++) {
332       if (TESTBIT(fDig[i],j)) return i*8+j;
333     }
334   }
335   return -1;
336 }
337
338
339 void AliMCInfo::Clear(Option_t* ){
340   //
341   // clear the content
342   //
343   if (fTPCReferences) fTPCReferences->Clear();
344   if (fITSReferences) fITSReferences->Clear();
345   if (fTRDReferences) fTRDReferences->Clear();
346   if (fTOFReferences) fTOFReferences->Clear();
347   for (Int_t i=0;i<4;i++) fVDist[i]=0;
348   for (Int_t i=0;i<3;i++) fDecayCoord[i]=0;
349   fTRdecay.SetTrack(-1);
350   fCharge = 0;
351   fEventNr=-1;
352   fPrimPart=-1;               // index of primary particle in TreeH
353   fMass=0;                   // mass of the particle
354   fCharge=0;                 // charge of the particle
355   fLabel=0;                  // track label
356   fEventNr=0;                // event number
357   fMCtracks=0;               // indication of how many times the track is retuturned back
358   fPdg=0;                        //pdg code
359   fTPCdecay=0;                  //indicates decay in TPC
360   fRowsWithDigitsInn=0;          // number of rows with digits in the inner sectors
361   fRowsWithDigits=0;             // number of rows with digits in the outer sectors
362   fRowsTrackLength=0;            // last - first row with digit
363   fTPCtrackLength=0;           // distance between first and last track reference
364   fPrim=0;                     // theoretical dedx in tpc according particle momenta and mass
365   fNTPCRef=0;                    // tpc references counter
366   fNITSRef=0;                    // ITS references counter
367   fNTRDRef=0;                    // TRD references counter
368   fNTOFRef=0;                    // TOF references counter
369   //
370   fNTPCRefOut=0;                    // tpc references counter - out
371   fNITSRefOut=0;                    // ITS references counter - out
372   fNTRDRefOut=0;                    // TRD references counter - out
373   fNTOFRefOut=0;                    // TOF references counter - out
374 }
375
376
377
378 void AliMCInfo::CalcTPCrows(TClonesArray * runArrayTR){
379   //
380   // Calculates the numebr of the track references for detectors
381   // In case of changing direction - curling tracks - the counter is not increasing
382   //
383   // Rough calculation 
384   // of the first and last point in the TPC  
385   //
386   fNTPCRef = 0;
387   fNITSRef = 0;
388   fNTRDRef = 0;
389   fNTOFRef = 0;  
390   fNTPCRefOut = 0;
391   fNITSRefOut = 0;
392   fNTRDRefOut = 0;
393   fNTOFRefOut = 0;  
394   //
395   //
396   Float_t tpcminRadius=TMath::Max(90., fParticle.R());
397   Float_t tpcmaxRadius=TMath::Min(250.,fParticle.R());  
398   Int_t nover=-1;
399
400   if (runArrayTR->GetEntriesFast()>0){  
401     nover=0;
402     AliTrackReference *ref0 = (AliTrackReference*)runArrayTR->At(0);    
403     Float_t dir = ((ref0->X()*ref0->Px()+ref0->Y()*ref0->Py())>0)? 1.:-1.; //inside or outside
404     //
405     //    
406     for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
407       //
408       AliTrackReference *ref = (AliTrackReference*)runArrayTR->At(iTrackRef);  
409       Float_t newdirection = ((ref->X()*ref->Px()+ref->Y()*ref->Py())>0)? 1.:-1.; //inside or outside
410       //
411       //
412       if (dir*newdirection<0.0) {
413         nover++;
414         dir = newdirection;
415       }
416       //
417       if (ref->DetectorId()== AliTrackReference::kTRD){
418         if (nover==0) {
419           tpcmaxRadius =250;
420           fNTRDRefOut++;
421         }
422         fNTRDRef++;
423       }
424       if (ref->DetectorId()== AliTrackReference::kITS){ 
425         if (nover==0) {
426           fNITSRefOut++;
427           tpcminRadius =90;
428         }
429         fNITSRef++;
430       }
431       if (ref->DetectorId()== AliTrackReference::kTOF){
432         fNTOFRef++;
433         if (nover==0) {
434           fNTOFRefOut++;
435           tpcmaxRadius =250;
436         }
437       }
438       //
439       if (ref->DetectorId()== AliTrackReference::kTPC){
440         fNTPCRef++;
441         if (nover==0) {
442           fNTPCRefOut++;
443           fTrackRefOut = (*ref);
444           if (ref->R()>tpcmaxRadius) tpcmaxRadius = ref->R();
445           if (ref->R()<tpcminRadius) tpcminRadius = ref->R();
446         }
447       }
448       if (ref->DetectorId()== AliTrackReference::kDisappeared){
449         if (nover==0 &&TMath::Abs(ref->Z())<250){
450           tpcmaxRadius = TMath::Min(Float_t(250.),ref->R());
451         }
452       }
453     }
454   }  // if we have track refs
455   fTPCtrackLength = tpcmaxRadius-tpcminRadius;
456   fMCtracks=nover+1;
457 }
458
459 //
460 //
461 //
462
463
464
465 void AliMCInfo::Update(TParticle * part, TClonesArray * runArrayTR, Double_t pvertex[4], Int_t label){
466   //
467   //
468   //
469   Clear();
470   fLabel=label;
471   fParticle=(*part);
472   AliMCInfo *info = this;
473   fVDist[0]=part->Vx()-pvertex[0];
474   fVDist[1]=part->Vy()-pvertex[1];
475   fVDist[2]=part->Vz()-pvertex[2];
476   fVDist[3]=part->R()-pvertex[3];
477   //
478   //
479   //
480   for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
481     AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);         
482     if (trackRef->DetectorId()== AliTrackReference::kTPC){     
483       TClonesArray & arr = *(info->fTPCReferences);
484       new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);     
485     }
486     if (trackRef->DetectorId()== AliTrackReference::kITS){
487       TClonesArray & arr = *(info->fITSReferences);
488       new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);     
489     }
490     if (trackRef->DetectorId()== AliTrackReference::kTRD){
491       TClonesArray & arr = *(info->fTRDReferences);
492       new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);     
493     }
494     if (trackRef->DetectorId()== AliTrackReference::kTOF){
495       TClonesArray & arr = *(info->fTOFReferences);
496       new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);     
497     }
498     //
499     // decay case
500     //
501     if (trackRef->DetectorId()== AliTrackReference::kDisappeared){
502       info->fTRdecay = *trackRef;       
503     }
504   }
505   //
506   Update();
507   CalcTPCrows(runArrayTR);
508 }
509