1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////
20 Origin: marian.ivanov@cern.ch
21 Container classes with MC infomation
23 The AliMCInfo contains the information about the particles properties
24 during transportation throuch ALICE Detector
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
35 some additional information usable for tree draw - TO SPEED UP tree queries
36 IMPORTANT FOR PROOF FAST PROTOTYPING ANALYSIS
43 #if !defined(__CINT__) || defined(__MAKECINT__)
47 #include "TClonesArray.h"
48 #include "TDatabasePDG.h"
50 #include "AliTrackReference.h"
51 #include "AliMCInfo.h"
52 #include "AliMathBase.h"
58 ClassImp(AliTPCdigitRow)
63 ////////////////////////////////////////////////////////////////////////
64 AliMCInfo::AliMCInfo():
77 fRowsWithDigitsInn(0),
83 fNTPCRef(0), // tpc references counter
84 fNITSRef(0), // ITS references counter
85 fNTRDRef(0), // TRD references counter
86 fNTOFRef(0), // TOF references counter
87 fNTPCRefOut(0), // tpc references counter
88 fNITSRefOut(0), // ITS references counter
89 fNTRDRefOut(0), // TRD references counter
90 fNTOFRefOut(0), // TOF references counter
97 // Default constructor
99 fTPCReferences = new TClonesArray("AliTrackReference",10);
100 fITSReferences = new TClonesArray("AliTrackReference",10);
101 fTRDReferences = new TClonesArray("AliTrackReference",10);
102 fTOFReferences = new TClonesArray("AliTrackReference",10);
103 fTRdecay.SetTrack(-1);
104 for (Int_t i=0;i<4;i++) fVDist[i]=0;
105 for (Int_t i=0;i<3;i++) fDecayCoord[i]=0;
109 AliMCInfo::AliMCInfo(const AliMCInfo& info):
111 fTrackRef(info.fTrackRef),
112 fTrackRefOut(info.fTrackRefOut),
113 fTRdecay(info.GetTRdecay()),
114 fPrimPart(info.fPrimPart),
115 fParticle(info.fParticle),
116 fMass(info.GetMass()),
117 fCharge(info.fCharge),
119 fEventNr(info.fEventNr),
120 fMCtracks(info.fMCtracks),
122 fTPCdecay(info.fTPCdecay),
123 fRowsWithDigitsInn(info.fRowsWithDigitsInn),
124 fRowsWithDigits(info.fRowsWithDigits),
125 fRowsTrackLength(info.fRowsTrackLength),
126 fTPCtrackLength(info.fTPCtrackLength),
128 fTPCRow(info.fTPCRow),
129 fNTPCRef(info.fNTPCRef), // tpc references counter
130 fNITSRef(info.fNITSRef), // ITS references counter
131 fNTRDRef(info.fNTRDRef), // TRD references counter
132 fNTOFRef(info.fNTOFRef), // TOF references counter
133 fNTPCRefOut(info.fNTPCRefOut), // tpc references counter
134 fNITSRefOut(info.fNITSRefOut), // ITS references counter
135 fNTRDRefOut(info.fNTRDRefOut), // TRD references counter
136 fNTOFRefOut(info.fNTOFRefOut), // TOF references counter
145 fTPCReferences = (TClonesArray*)info.fTPCReferences->Clone();
146 fITSReferences = (TClonesArray*)info.fITSReferences->Clone();
147 fTRDReferences = (TClonesArray*)info.fTRDReferences->Clone();
148 fTOFReferences = (TClonesArray*)info.fTOFReferences->Clone();
149 for (Int_t i=0;i<4;i++) fVDist[i]=info.fVDist[i];
150 for (Int_t i=0;i<3;i++) fDecayCoord[i]=info.fDecayCoord[i];
155 AliMCInfo& AliMCInfo::operator=(const AliMCInfo& info) {
157 // Assignment operator
160 new (this) AliMCInfo(info);
165 AliMCInfo::~AliMCInfo()
168 // Destructor of the class
170 if (fTPCReferences) {
171 delete fTPCReferences;
174 delete fITSReferences;
177 delete fTRDReferences;
180 delete fTOFReferences;
188 void AliMCInfo::Update()
192 // Calculates some derived variables
199 if (fTRdecay.GetTrack()>0){
200 fDecayCoord[0] = fTRdecay.X();
201 fDecayCoord[1] = fTRdecay.Y();
202 fDecayCoord[2] = fTRdecay.Z();
203 if ( (fTRdecay.R()<250)&&(fTRdecay.R()>85) && (TMath::Abs(fTRdecay.Z())<250) ){
214 //digits information update
215 fRowsWithDigits = fTPCRow.RowsOn();
216 fRowsWithDigitsInn = fTPCRow.RowsOn(63); // 63 = number of inner rows
217 fRowsTrackLength = fTPCRow.Last() - fTPCRow.First();
220 // calculate primary ionization per cm
221 if (fParticle.GetPDG()){
222 fMass = fParticle.GetMass();
223 fCharge = fParticle.GetPDG()->Charge();
224 if (fTPCReferences->GetEntriesFast()>0){
225 fTrackRef = *((AliTrackReference*)fTPCReferences->At(0));
228 Float_t p = TMath::Sqrt(fTrackRef.Px()*fTrackRef.Px()+
229 fTrackRef.Py()*fTrackRef.Py()+
230 fTrackRef.Pz()*fTrackRef.Pz());
232 Float_t betagama = p /fMass;
233 fPrim = AliMathBase::BetheBlochAleph(betagama);
248 ////////////////////////////////////////////////////////////////////////
249 AliTPCdigitRow::AliTPCdigitRow()
253 ////////////////////////////////////////////////////////////////////////
254 AliTPCdigitRow & AliTPCdigitRow::operator=(const AliTPCdigitRow &digOld)
256 for (Int_t i = 0; i<32; i++) fDig[i] = digOld.fDig[i];
259 ////////////////////////////////////////////////////////////////////////
260 void AliTPCdigitRow::SetRow(Int_t row)
263 // set bit mask for given row
266 // cerr<<"AliTPCdigitRow::SetRow: index "<<row<<" out of bounds."<<endl;
274 ////////////////////////////////////////////////////////////////////////
275 Bool_t AliTPCdigitRow::TestRow(Int_t row) const
278 // return kTRUE if row is on
282 return TESTBIT(fDig[iC],iB);
284 ////////////////////////////////////////////////////////////////////////
285 Int_t AliTPCdigitRow::RowsOn(Int_t upto) const
288 // returns number of rows with a digit
289 // count only rows less equal row number upto
292 for (Int_t i = 0; i<32; i++) {
293 for (Int_t j = 0; j < 8; j++) {
294 if (i*8+j > upto) return total;
295 if (TESTBIT(fDig[i],j)) total++;
300 ////////////////////////////////////////////////////////////////////////
301 void AliTPCdigitRow::Reset()
304 // resets all rows to zero
306 for (Int_t i = 0; i<32; i++) {
310 ////////////////////////////////////////////////////////////////////////
311 Int_t AliTPCdigitRow::Last() const
314 // returns the last row number with a digit
315 // returns -1 if now digits
317 for (Int_t i = 32-1; i>=0; i--) {
318 for (Int_t j = 7; j >= 0; j--) {
319 if TESTBIT(fDig[i],j) return i*8+j;
324 ////////////////////////////////////////////////////////////////////////
325 Int_t AliTPCdigitRow::First() const
328 // returns the first row number with a digit
329 // returns -1 if now digits
331 for (Int_t i = 0; i<32; i++) {
332 for (Int_t j = 0; j < 8; j++) {
333 if (TESTBIT(fDig[i],j)) return i*8+j;
340 void AliMCInfo::Clear(Option_t* ){
344 if (fTPCReferences) fTPCReferences->Clear();
345 if (fITSReferences) fITSReferences->Clear();
346 if (fTRDReferences) fTRDReferences->Clear();
347 if (fTOFReferences) fTOFReferences->Clear();
348 for (Int_t i=0;i<4;i++) fVDist[i]=0;
349 for (Int_t i=0;i<3;i++) fDecayCoord[i]=0;
350 fTRdecay.SetTrack(-1);
353 fPrimPart=-1; // index of primary particle in TreeH
354 fMass=0; // mass of the particle
355 fCharge=0; // charge of the particle
356 fLabel=0; // track label
357 fEventNr=0; // event number
358 fMCtracks=0; // indication of how many times the track is retuturned back
360 fTPCdecay=0; //indicates decay in TPC
361 fRowsWithDigitsInn=0; // number of rows with digits in the inner sectors
362 fRowsWithDigits=0; // number of rows with digits in the outer sectors
363 fRowsTrackLength=0; // last - first row with digit
364 fTPCtrackLength=0; // distance between first and last track reference
365 fPrim=0; // theoretical dedx in tpc according particle momenta and mass
366 fNTPCRef=0; // tpc references counter
367 fNITSRef=0; // ITS references counter
368 fNTRDRef=0; // TRD references counter
369 fNTOFRef=0; // TOF references counter
371 fNTPCRefOut=0; // tpc references counter - out
372 fNITSRefOut=0; // ITS references counter - out
373 fNTRDRefOut=0; // TRD references counter - out
374 fNTOFRefOut=0; // TOF references counter - out
379 void AliMCInfo::CalcTPCrows(TClonesArray * runArrayTR){
381 // Calculates the numebr of the track references for detectors
382 // In case of changing direction - curling tracks - the counter is not increasing
385 // of the first and last point in the TPC
397 Float_t tpcminRadius=TMath::Max(90., fParticle.R());
398 Float_t tpcmaxRadius=TMath::Min(250.,fParticle.R());
401 if (runArrayTR->GetEntriesFast()>0){
403 AliTrackReference *ref0 = (AliTrackReference*)runArrayTR->At(0);
404 Float_t dir = ((ref0->X()*ref0->Px()+ref0->Y()*ref0->Py())>0)? 1.:-1.; //inside or outside
407 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
409 AliTrackReference *ref = (AliTrackReference*)runArrayTR->At(iTrackRef);
410 Float_t newdirection = ((ref->X()*ref->Px()+ref->Y()*ref->Py())>0)? 1.:-1.; //inside or outside
413 if (dir*newdirection<0.0) {
418 if (ref->DetectorId()== AliTrackReference::kTRD){
425 if (ref->DetectorId()== AliTrackReference::kITS){
432 if (ref->DetectorId()== AliTrackReference::kTOF){
440 if (ref->DetectorId()== AliTrackReference::kTPC){
444 fTrackRefOut = (*ref);
445 if (ref->R()>tpcmaxRadius) tpcmaxRadius = ref->R();
446 if (ref->R()<tpcminRadius) tpcminRadius = ref->R();
449 if (ref->DetectorId()== AliTrackReference::kDisappeared){
450 if (nover==0 &&TMath::Abs(ref->Z())<250){
451 tpcmaxRadius = TMath::Min(Float_t(250.),ref->R());
455 } // if we have track refs
456 fTPCtrackLength = tpcmaxRadius-tpcminRadius;
466 void AliMCInfo::Update(TParticle * part, TClonesArray * runArrayTR, Double_t pvertex[4], Int_t label){
473 AliMCInfo *info = this;
474 fVDist[0]=part->Vx()-pvertex[0];
475 fVDist[1]=part->Vy()-pvertex[1];
476 fVDist[2]=part->Vz()-pvertex[2];
477 fVDist[3]=part->R()-pvertex[3];
481 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
482 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
483 if (trackRef->DetectorId()== AliTrackReference::kTPC){
484 TClonesArray & arr = *(info->fTPCReferences);
485 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
487 if (trackRef->DetectorId()== AliTrackReference::kITS){
488 TClonesArray & arr = *(info->fITSReferences);
489 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
491 if (trackRef->DetectorId()== AliTrackReference::kTRD){
492 TClonesArray & arr = *(info->fTRDReferences);
493 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
495 if (trackRef->DetectorId()== AliTrackReference::kTOF){
496 TClonesArray & arr = *(info->fTOFReferences);
497 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
502 if (trackRef->DetectorId()== AliTrackReference::kDisappeared){
503 info->fTRdecay = *trackRef;
508 CalcTPCrows(runArrayTR);