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()
251 for (Int_t i = 0; i<32; i++) { fDig[i] = 0; }
254 ////////////////////////////////////////////////////////////////////////
255 AliTPCdigitRow & AliTPCdigitRow::operator=(const AliTPCdigitRow &digOld)
257 for (Int_t i = 0; i<32; i++) fDig[i] = digOld.fDig[i];
260 ////////////////////////////////////////////////////////////////////////
261 void AliTPCdigitRow::SetRow(Int_t row)
264 // set bit mask for given row
267 // cerr<<"AliTPCdigitRow::SetRow: index "<<row<<" out of bounds."<<endl;
275 ////////////////////////////////////////////////////////////////////////
276 Bool_t AliTPCdigitRow::TestRow(Int_t row) const
279 // return kTRUE if row is on
283 return TESTBIT(fDig[iC],iB);
285 ////////////////////////////////////////////////////////////////////////
286 Int_t AliTPCdigitRow::RowsOn(Int_t upto) const
289 // returns number of rows with a digit
290 // count only rows less equal row number upto
293 for (Int_t i = 0; i<32; i++) {
294 for (Int_t j = 0; j < 8; j++) {
295 if (i*8+j > upto) return total;
296 if (TESTBIT(fDig[i],j)) total++;
301 ////////////////////////////////////////////////////////////////////////
302 void AliTPCdigitRow::Reset()
305 // resets all rows to zero
307 for (Int_t i = 0; i<32; i++) {
312 ////////////////////////////////////////////////////////////////////////
313 Int_t AliTPCdigitRow::Last() const
316 // returns the last row number with a digit
317 // returns -1 if now digits
319 for (Int_t i = 32-1; i>=0; i--) {
320 for (Int_t j = 7; j >= 0; j--) {
321 if TESTBIT(fDig[i],j) return i*8+j;
326 ////////////////////////////////////////////////////////////////////////
327 Int_t AliTPCdigitRow::First() const
330 // returns the first row number with a digit
331 // returns -1 if now digits
333 for (Int_t i = 0; i<32; i++) {
334 for (Int_t j = 0; j < 8; j++) {
335 if (TESTBIT(fDig[i],j)) return i*8+j;
342 void AliMCInfo::Clear(Option_t* ){
346 if (fTPCReferences) fTPCReferences->Clear();
347 if (fITSReferences) fITSReferences->Clear();
348 if (fTRDReferences) fTRDReferences->Clear();
349 if (fTOFReferences) fTOFReferences->Clear();
350 for (Int_t i=0;i<4;i++) fVDist[i]=0;
351 for (Int_t i=0;i<3;i++) fDecayCoord[i]=0;
352 fTRdecay.SetTrack(-1);
355 fPrimPart=-1; // index of primary particle in TreeH
356 fMass=0; // mass of the particle
357 fCharge=0; // charge of the particle
358 fLabel=0; // track label
359 fEventNr=0; // event number
360 fMCtracks=0; // indication of how many times the track is retuturned back
362 fTPCdecay=0; //indicates decay in TPC
363 fRowsWithDigitsInn=0; // number of rows with digits in the inner sectors
364 fRowsWithDigits=0; // number of rows with digits in the outer sectors
365 fRowsTrackLength=0; // last - first row with digit
366 fTPCtrackLength=0; // distance between first and last track reference
367 fPrim=0; // theoretical dedx in tpc according particle momenta and mass
368 fNTPCRef=0; // tpc references counter
369 fNITSRef=0; // ITS references counter
370 fNTRDRef=0; // TRD references counter
371 fNTOFRef=0; // TOF references counter
373 fNTPCRefOut=0; // tpc references counter - out
374 fNITSRefOut=0; // ITS references counter - out
375 fNTRDRefOut=0; // TRD references counter - out
376 fNTOFRefOut=0; // TOF references counter - out
381 void AliMCInfo::CalcTPCrows(TClonesArray * runArrayTR){
383 // Calculates the numebr of the track references for detectors
384 // In case of changing direction - curling tracks - the counter is not increasing
387 // of the first and last point in the TPC
399 Float_t tpcminRadius=TMath::Max(90., fParticle.R());
400 Float_t tpcmaxRadius=TMath::Min(250.,fParticle.R());
403 if (runArrayTR->GetEntriesFast()>0){
405 AliTrackReference *ref0 = (AliTrackReference*)runArrayTR->At(0);
406 Float_t dir = ((ref0->X()*ref0->Px()+ref0->Y()*ref0->Py())>0)? 1.:-1.; //inside or outside
409 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
411 AliTrackReference *ref = (AliTrackReference*)runArrayTR->At(iTrackRef);
412 Float_t newdirection = ((ref->X()*ref->Px()+ref->Y()*ref->Py())>0)? 1.:-1.; //inside or outside
415 if (dir*newdirection<0.0) {
420 if (ref->DetectorId()== AliTrackReference::kTRD){
427 if (ref->DetectorId()== AliTrackReference::kITS){
434 if (ref->DetectorId()== AliTrackReference::kTOF){
442 if (ref->DetectorId()== AliTrackReference::kTPC){
446 fTrackRefOut = (*ref);
447 if (ref->R()>tpcmaxRadius) tpcmaxRadius = ref->R();
448 if (ref->R()<tpcminRadius) tpcminRadius = ref->R();
451 if (ref->DetectorId()== AliTrackReference::kDisappeared){
452 if (nover==0 &&TMath::Abs(ref->Z())<250){
453 tpcmaxRadius = TMath::Min(Float_t(250.),ref->R());
457 } // if we have track refs
458 fTPCtrackLength = tpcmaxRadius-tpcminRadius;
468 void AliMCInfo::Update(TParticle * part, TClonesArray * runArrayTR, Double_t pvertex[4], Int_t label){
475 AliMCInfo *info = this;
476 fVDist[0]=part->Vx()-pvertex[0];
477 fVDist[1]=part->Vy()-pvertex[1];
478 fVDist[2]=part->Vz()-pvertex[2];
479 fVDist[3]=part->R()-pvertex[3];
483 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
484 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
485 if (trackRef->DetectorId()== AliTrackReference::kTPC){
486 TClonesArray & arr = *(info->fTPCReferences);
487 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
489 if (trackRef->DetectorId()== AliTrackReference::kITS){
490 TClonesArray & arr = *(info->fITSReferences);
491 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
493 if (trackRef->DetectorId()== AliTrackReference::kTRD){
494 TClonesArray & arr = *(info->fTRDReferences);
495 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
497 if (trackRef->DetectorId()== AliTrackReference::kTOF){
498 TClonesArray & arr = *(info->fTOFReferences);
499 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
504 if (trackRef->DetectorId()== AliTrackReference::kDisappeared){
505 info->fTRdecay = *trackRef;
510 CalcTPCrows(runArrayTR);