]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/ESD/AliTrackerBase.cxx
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / STEER / ESD / AliTrackerBase.cxx
index 72d0db4aafdf082c9f4cb2e131068841ff69744a..b19efbdf31362d62a3939d94675c5b4ebb7f7e33 100644 (file)
@@ -71,7 +71,11 @@ AliTrackerBase::AliTrackerBase(const AliTrackerBase &atr):
 Double_t AliTrackerBase::GetBz()
 {
   AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
-  if (!fld) return 0.5*kAlmost0Field;
+  if (!fld) {
+    AliFatalClass("Field is not loaded");
+    //if (!fld) 
+    return  0.5*kAlmost0Field;
+  }
   Double_t bz = fld->SolenoidField();
   return TMath::Sign(0.5*kAlmost0Field,bz) + bz;
 }
@@ -82,7 +86,11 @@ Double_t AliTrackerBase::GetBz(const Double_t *r) {
   // Returns Bz (kG) at the point "r" .
   //------------------------------------------------------------------
   AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
-  if (!fld) return  0.5*kAlmost0Field;
+  if (!fld) {
+    AliFatalClass("Field is not loaded");
+    //  if (!fld) 
+    return  0.5*kAlmost0Field;
+  }
   Double_t bz = fld->GetBz(r);
   return  TMath::Sign(0.5*kAlmost0Field,bz) + bz;
 }
@@ -94,9 +102,10 @@ void AliTrackerBase::GetBxByBz(const Double_t r[3], Double_t b[3]) {
   //------------------------------------------------------------------
   AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
   if (!fld) {
-     b[0] = b[1] = 0.;
-     b[2] = 0.5*kAlmost0Field;
-     return;
+    AliFatalClass("Field is not loaded");
+    // b[0] = b[1] = 0.;
+    // b[2] = 0.5*kAlmost0Field;
+    return;
   }
 
   if (fld->IsUniform()) {
@@ -142,7 +151,7 @@ Double_t AliTrackerBase::MeanMaterialBudget(const Double_t *start, const Double_
   for (Int_t i=0;i<6;i++) bparam[i]=0;
 
   if (!gGeoManager) {
-    AliErrorClass("No TGeo\n");
+    AliFatalClass("No TGeo\n");
     return 0.;
   }
   //
@@ -259,13 +268,13 @@ Double_t AliTrackerBase::MeanMaterialBudget(const Double_t *start, const Double_
 
 Bool_t 
 AliTrackerBase::PropagateTrackTo(AliExternalTrackParam *track, Double_t xToGo, 
-                                Double_t mass, Double_t maxStep, Bool_t rotateTo, Double_t maxSnp, Int_t sign, Bool_t addTimeStep){
+                                Double_t mass, Double_t maxStep, Bool_t rotateTo, Double_t maxSnp, Int_t sign, Bool_t addTimeStep, Bool_t correctMaterialBudget){
   //----------------------------------------------------------------
   //
   // Propagates the track to the plane X=xk (cm) using the magnetic field map 
   // and correcting for the crossed material.
   //
-  // mass     - mass used in propagation - used for energy loss correction
+  // mass     - mass used in propagation - used for energy loss correction (if <0 then q = 2)
   // maxStep  - maximal step for propagation
   //
   //  Origin: Marian Ivanov,  Marian.Ivanov@cern.ch
@@ -282,36 +291,42 @@ AliTrackerBase::PropagateTrackTo(AliExternalTrackParam *track, Double_t xToGo,
     track->GetXYZ(xyz0);   //starting global position
 
     Double_t bz=GetBz(xyz0); // getting the local Bz
-
     if (!track->GetXYZAt(x,bz,xyz1)) return kFALSE;   // no prolongation
     xyz1[2]+=kEpsilon; // waiting for bug correction in geo
-
-    if (TMath::Abs(track->GetSnpAt(x,bz)) >= maxSnp) return kFALSE;
+    
+    if (maxSnp>0 && TMath::Abs(track->GetSnpAt(x,bz)) >= maxSnp) return kFALSE;
     if (!track->PropagateTo(x,bz))  return kFALSE;
 
-    MeanMaterialBudget(xyz0,xyz1,param);       
-    Double_t xrho=param[0]*param[4], xx0=param[1];
-    if (sign) {if (sign<0) xrho = -xrho;}  // sign is imposed
-    else { // determine automatically the sign from direction
-      if (dir>0) xrho = -xrho; // outward should be negative
+    if (correctMaterialBudget){
+      MeanMaterialBudget(xyz0,xyz1,param);
+      Double_t xrho=param[0]*param[4], xx0=param[1];
+      if (sign) {if (sign<0) xrho = -xrho;}  // sign is imposed
+      else { // determine automatically the sign from direction
+        if (dir>0) xrho = -xrho; // outward should be negative
+      }
+      //
+      if (!track->CorrectForMeanMaterial(xx0,xrho,mass)) return kFALSE;
     }
-    //
-    if (!track->CorrectForMeanMaterial(xx0,xrho,mass)) return kFALSE;
+    
     if (rotateTo){
-      if (TMath::Abs(track->GetSnp()) >= maxSnp) return kFALSE;
-      track->GetXYZ(xyz0);   // global position
-      Double_t alphan = TMath::ATan2(xyz0[1], xyz0[0]); 
-      //
-      Double_t ca=TMath::Cos(alphan-track->GetAlpha()), 
-               sa=TMath::Sin(alphan-track->GetAlpha());
-      Double_t sf=track->GetSnp(), cf=TMath::Sqrt((1.-sf)*(1.+sf));
-      Double_t sinNew =  sf*ca - cf*sa;
-      if (TMath::Abs(sinNew) >= maxSnp) return kFALSE;
-      if (!track->Rotate(alphan)) return kFALSE;
+      track->GetXYZ(xyz1);   // global position
+      Double_t alphan = TMath::ATan2(xyz1[1], xyz1[0]); 
+      if (maxSnp>0) {
+       if (TMath::Abs(track->GetSnp()) >= maxSnp) return kFALSE;
+        
+       //
+       Double_t ca=TMath::Cos(alphan-track->GetAlpha()), sa=TMath::Sin(alphan-track->GetAlpha());
+       Double_t sf=track->GetSnp(), cf=TMath::Sqrt((1.-sf)*(1.+sf));
+       Double_t sinNew =  sf*ca - cf*sa;
+        if (TMath::Abs(sinNew) >= maxSnp) return kFALSE;
+        
+      }
+      if (!track->AliExternalTrackParam::Rotate(alphan)) return kFALSE;
+      
     }
     xpos = track->GetX();
     if (addTimeStep && track->IsStartedTimeIntegral()) {
-      track->GetXYZ(xyz1);
+      if (!rotateTo) track->GetXYZ(xyz1); // if rotateTo==kTRUE, then xyz1 is already extracted
       Double_t dX=xyz0[0]-xyz1[0],dY=xyz0[1]-xyz1[1],dZ=xyz0[2]-xyz1[2]; 
       Double_t d=TMath::Sqrt(dX*dX + dY*dY + dZ*dZ);
       if (sign) {if (sign>0) d = -d;}  // step sign is imposed, positive means inward direction
@@ -324,6 +339,78 @@ AliTrackerBase::PropagateTrackTo(AliExternalTrackParam *track, Double_t xToGo,
   return kTRUE;
 }
 
+Int_t AliTrackerBase::PropagateTrackTo2(AliExternalTrackParam *track, Double_t xToGo,
+                                        Double_t mass, Double_t maxStep, Bool_t rotateTo, Double_t maxSnp, Int_t sign, Bool_t addTimeStep, Bool_t correctMaterialBudget){
+  //----------------------------------------------------------------
+  //
+  // Propagates the track to the plane X=xk (cm) using the magnetic field map
+  // and correcting for the crossed material.
+  //
+  // mass     - mass used in propagation - used for energy loss correction
+  // maxStep  - maximal step for propagation
+  //
+  //  Origin: Marian Ivanov,  Marian.Ivanov@cern.ch
+  //
+  //----------------------------------------------------------------
+  const Double_t kEpsilon = 0.00001;
+  Double_t xpos     = track->GetX();
+  Int_t dir         = (xpos<xToGo) ? 1:-1;
+  //
+  while ( (xToGo-xpos)*dir > kEpsilon){
+    Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
+    Double_t x    = xpos+step;
+    Double_t xyz0[3],xyz1[3],param[7];
+    track->GetXYZ(xyz0);   //starting global position
+    
+    Double_t bz=GetBz(xyz0); // getting the local Bz
+    if (!track->GetXYZAt(x,bz,xyz1)) return -1;   // no prolongation
+    xyz1[2]+=kEpsilon; // waiting for bug correction in geo
+    
+    if (maxSnp>0 && TMath::Abs(track->GetSnpAt(x,bz)) >= maxSnp) return -2;
+    if (!track->PropagateTo(x,bz))  return -3;
+    
+    if (correctMaterialBudget){
+      MeanMaterialBudget(xyz0,xyz1,param);
+      Double_t xrho=param[0]*param[4], xx0=param[1];
+      if (sign) {if (sign<0) xrho = -xrho;}  // sign is imposed
+      else { // determine automatically the sign from direction
+        if (dir>0) xrho = -xrho; // outward should be negative
+      }
+      //
+      if (!track->CorrectForMeanMaterial(xx0,xrho,mass)) return -4;
+    }
+    
+    if (rotateTo){
+      track->GetXYZ(xyz1);   // global position
+      Double_t alphan = TMath::ATan2(xyz1[1], xyz1[0]);
+      if (maxSnp>0) {
+        if (TMath::Abs(track->GetSnp()) >= maxSnp) return -5;
+        
+        //
+        Double_t ca=TMath::Cos(alphan-track->GetAlpha()), sa=TMath::Sin(alphan-track->GetAlpha());
+        Double_t sf=track->GetSnp(), cf=TMath::Sqrt((1.-sf)*(1.+sf));
+        Double_t sinNew =  sf*ca - cf*sa;
+        if (TMath::Abs(sinNew) >= maxSnp) return -6;
+        
+      }
+      if (!track->AliExternalTrackParam::Rotate(alphan)) return -7;
+      
+    }
+    xpos = track->GetX();
+    if (addTimeStep && track->IsStartedTimeIntegral()) {
+      if (!rotateTo) track->GetXYZ(xyz1); // if rotateTo==kTRUE, then xyz1 is already extracted
+      Double_t dX=xyz0[0]-xyz1[0],dY=xyz0[1]-xyz1[1],dZ=xyz0[2]-xyz1[2];
+      Double_t d=TMath::Sqrt(dX*dX + dY*dY + dZ*dZ);
+      if (sign) {if (sign>0) d = -d;}  // step sign is imposed, positive means inward direction
+      else { // determine automatically the sign from direction
+        if (dir<0) d = -d;
+      }
+      track->AddTimeStep(d);
+    }
+  }
+  return 1;
+}
+
 Bool_t 
 AliTrackerBase::PropagateTrackToBxByBz(AliExternalTrackParam *track,
                                       Double_t xToGo,Double_t mass, Double_t maxStep, Bool_t rotateTo, Double_t maxSnp,Int_t sign, Bool_t addTimeStep){
@@ -333,7 +420,7 @@ AliTrackerBase::PropagateTrackToBxByBz(AliExternalTrackParam *track,
   // taking into account all the three components of the magnetic field 
   // and correcting for the crossed material.
   //
-  // mass     - mass used in propagation - used for energy loss correction
+  // mass     - mass used in propagation - used for energy loss correction (if <0 then q=2)
   // maxStep  - maximal step for propagation
   //
   //  Origin: Marian Ivanov,  Marian.Ivanov@cern.ch
@@ -354,7 +441,7 @@ AliTrackerBase::PropagateTrackToBxByBz(AliExternalTrackParam *track,
     if (!track->GetXYZAt(x,b[2],xyz1)) return kFALSE;   // no prolongation
     xyz1[2]+=kEpsilon; // waiting for bug correction in geo
 
-    if (TMath::Abs(track->GetSnpAt(x,b[2])) >= maxSnp) return kFALSE;
+    if (maxSnp>0 && TMath::Abs(track->GetSnpAt(x,b[2])) >= maxSnp) return kFALSE;
     if (!track->PropagateToBxByBz(x,b))  return kFALSE;
 
     MeanMaterialBudget(xyz0,xyz1,param);    
@@ -366,20 +453,20 @@ AliTrackerBase::PropagateTrackToBxByBz(AliExternalTrackParam *track,
     //
     if (!track->CorrectForMeanMaterial(xx0,xrho,mass)) return kFALSE;
     if (rotateTo){
-      if (TMath::Abs(track->GetSnp()) >= maxSnp) return kFALSE;
-      track->GetXYZ(xyz0);   // global position
-      Double_t alphan = TMath::ATan2(xyz0[1], xyz0[0]); 
-      //
-      Double_t ca=TMath::Cos(alphan-track->GetAlpha()), 
-               sa=TMath::Sin(alphan-track->GetAlpha());
-      Double_t sf=track->GetSnp(), cf=TMath::Sqrt((1.-sf)*(1.+sf));
-      Double_t sinNew =  sf*ca - cf*sa;
-      if (TMath::Abs(sinNew) >= maxSnp) return kFALSE;
-      if (!track->Rotate(alphan)) return kFALSE;
+      track->GetXYZ(xyz1);   // global position
+      Double_t alphan = TMath::ATan2(xyz1[1], xyz1[0]); 
+      if (maxSnp>0) {
+       if (TMath::Abs(track->GetSnp()) >= maxSnp) return kFALSE;
+       Double_t ca=TMath::Cos(alphan-track->GetAlpha()), sa=TMath::Sin(alphan-track->GetAlpha());
+       Double_t sf=track->GetSnp(), cf=TMath::Sqrt((1.-sf)*(1.+sf));
+       Double_t sinNew =  sf*ca - cf*sa;
+       if (TMath::Abs(sinNew) >= maxSnp) return kFALSE;
+      }
+      if (!track->AliExternalTrackParam::Rotate(alphan)) return kFALSE;
     }
     xpos = track->GetX();    
     if (addTimeStep && track->IsStartedTimeIntegral()) {
-      track->GetXYZ(xyz1);
+      if (!rotateTo) track->GetXYZ(xyz1); // if rotateTo==kTRUE, then xyz1 is already extracted
       Double_t dX=xyz0[0]-xyz1[0],dY=xyz0[1]-xyz1[1],dZ=xyz0[2]-xyz1[2]; 
       Double_t d=TMath::Sqrt(dX*dX + dY*dY + dZ*dZ);
       if (sign) {if (sign>0) d = -d;}  // step sign is imposed, positive means inward direction
@@ -477,13 +564,18 @@ Double_t AliTrackerBase::MakeTgl(Double_t x1,Double_t y1,
   // Initial approximation of the tangent of the track dip angle
   //-----------------------------------------------------------------
   //
+  const Double_t kEpsilon =0.00001;
   x2-=x1;
   y2-=y1;
   z2-=z1;
   Double_t d  =  TMath::Sqrt(x2*x2+y2*y2);  // distance  straight line
   if (TMath::Abs(d*c*0.5)>1) return 0;
   Double_t   angle2    = TMath::ASin(d*c*0.5); 
-  angle2  = z2*TMath::Abs(c/(angle2*2.));
+  if (TMath::Abs(angle2)>kEpsilon)  {
+    angle2  = z2*TMath::Abs(c/(angle2*2.));
+  }else{
+    angle2=z2/d;
+  }
   return angle2;
 }