]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliMCParticle.cxx
AliESDHeader: AliTriggerConfiguration and more trigger scalers added
[u/mrichter/AliRoot.git] / STEER / AliMCParticle.cxx
index ffb377d2f48e2e930b87468b6e2fb2344a783b9b..c4b66304cdc047502f942a4a226ff966394999d1 100644 (file)
@@ -24,6 +24,7 @@
 #include <TRefArray.h>
 
 #include "AliMCParticle.h"
+#include "AliExternalTrackParam.h" 
 
 
 ClassImp(AliMCParticle)
@@ -32,16 +33,25 @@ AliMCParticle::AliMCParticle():
     AliVParticle(),
     fParticle(0),
     fTrackReferences(0),
-    fNTrackRef(0)
+    fNTrackRef(0),
+    fLabel(-1),
+    fMother(-1),
+    fFirstDaughter(-1),
+    fLastDaughter(-1)
 {
     // Constructor
 }
 
-AliMCParticle::AliMCParticle(TParticle* part, TRefArray* rarray):
+    
+AliMCParticle::AliMCParticle(TParticle* part, TRefArray* rarray, Int_t index):
     AliVParticle(),
     fParticle(part),
     fTrackReferences(rarray),
-    fNTrackRef(0)
+    fNTrackRef(0),
+    fLabel(index),
+    fMother(-1),
+    fFirstDaughter(-1),
+    fLastDaughter(-1)
 {
     // Constructor
     if (rarray != 0) {
@@ -54,14 +64,18 @@ AliMCParticle::AliMCParticle(const AliMCParticle& mcPart) :
     AliVParticle(mcPart),
     fParticle(0),    
     fTrackReferences(0),
-    fNTrackRef(0)
+    fNTrackRef(0),
+    fLabel(-1),
+    fMother(-1),
+    fFirstDaughter(-1),
+    fLastDaughter(-1)
 {
 // Copy constructor
 }
 
 AliMCParticle& AliMCParticle::operator=(const AliMCParticle& mcPart)
 { 
-
+// Copy constructor
   if (this!=&mcPart) { 
     AliVParticle::operator=(mcPart);
   }
@@ -75,10 +89,86 @@ AliMCParticle::~AliMCParticle()
   // fParticle should be handled by the user
   // AliStack in case of AliMCEventHandler
   if(fTrackReferences){
-    delete fTrackReferences;
-    fTrackReferences = 0;
+      fTrackReferences->Delete();
+      delete fTrackReferences;
+      fTrackReferences = 0;
   }
 }
 
 
 
+Float_t  AliMCParticle::GetTPCTrackLength(Float_t bz,  Float_t ptmin, Int_t &counter, Float_t deadWidth){
+  //
+  // return track length in geometrically active volume of TPC.
+  // z nad rphi acceptance is included
+  // doesn't take into account dead channel and ExB  
+  // Intput:
+  // trackRefs
+  // bz - magnetic field
+  // deadWidth - dead zone in r-phi
+  // Additional output:
+  // counter   - number of circles
+
+    if (fNTrackRef == 0) return 0.; 
+
+    const Float_t kRMin = 90;
+    const Float_t kRMax = 245;
+    const Float_t kZMax = 250;
+    const Float_t kMinPt= ptmin; 
+
+    Float_t length =0;
+    Int_t nrefs = fNTrackRef;
+
+    
+    AliExternalTrackParam param;
+    Double_t cv[21];
+    for (Int_t i = 0; i < 21; i++) cv[i]=0;
+    counter=0; 
+    //
+    //
+    
+    AliTrackReference *ref0 = (AliTrackReference*) (fTrackReferences->At(0));
+    Float_t direction = 0;
+  //
+    for (Int_t iref = 1; iref < nrefs; iref++){
+       AliTrackReference *ref = (AliTrackReference*) (fTrackReferences->At(iref));
+       if (!ref) continue;
+       if (!ref0 || ref0->DetectorId()!= AliTrackReference::kTPC){
+           ref0 = ref;
+           direction = ((ref0->X() * ref0->Px() + ref0->Y() * ref0->Py()) > 0)? 1. : -1.;
+           continue;
+       }
+       
+       Float_t newdirection = ((ref->X() * ref->Px() + ref->Y() * ref->Py()) > 0)? 1. : -1.;
+       if (newdirection*direction<0) {
+           counter++;  //circle counter 
+           direction = newdirection;
+           continue;
+       }
+       if (counter>0) continue;
+       if (ref0->Pt() < kMinPt) break;
+       Float_t radius0 = TMath::Max(TMath::Min(ref0->R(),kRMax),kRMin);;
+       Float_t radius1 = TMath::Max(TMath::Min(ref->R(),kRMax),kRMin);
+       Double_t xyz[3] = {ref0->X(), ref0->Y(), ref0->Z()};
+       Double_t pxyz[3]= {ref0->Px(), ref0->Py(), ref0->Pz()};
+       Double_t alpha;
+       param.Set(xyz,pxyz,cv,TMath::Nint(fParticle->GetPDG()->Charge()/3.));
+       
+       for (Float_t radius = radius0; radius < radius1; radius+=1){
+           param.GetXYZAt(radius, bz, xyz);
+           if (TMath::Abs(xyz[2]) > kZMax) continue;
+           Float_t gradius = TMath::Sqrt(xyz[1] * xyz[1] + xyz[0] * xyz[0]);
+           if (gradius > kRMax) continue;
+           alpha = TMath::ATan2(xyz[1],xyz[0]);
+           if (alpha<0) alpha += TMath::TwoPi();
+           //
+           Int_t sector   = Int_t(9 * alpha / TMath::Pi());
+           Float_t lalpha = alpha - ((sector + 0.5) * TMath::Pi() / 9.);
+           Float_t dedge  = (TMath::Tan(TMath::Pi() / 18.) - TMath::Abs(TMath::Tan(lalpha))) * gradius; 
+           if (dedge>deadWidth) length++;
+       }
+       if (ref->DetectorId()!= AliTrackReference::kTPC) break; 
+       ref0 = ref;
+    }
+    return length;
+}