In reconstruction:

 - Correct calculations of energy loss and multiple scattering in absorber.
 - Avoid re-extrapolating tracks to vertex in MUONRecoCheck.C
(Diego and Philippe P.)

diff --git a/MUON/AliMUONTrackExtrap.cxx b/MUON/AliMUONTrackExtrap.cxx
index 3ab3ebee54f560ed7d687d10b78107c4d9c46575..09bef71649ee2bc5019bcd51a75aa6f673e945f5 100644 (file)
--- a/MUON/AliMUONTrackExtrap.cxx
+++ b/MUON/AliMUONTrackExtrap.cxx
@@ -561,6 +561,7 @@ Bool_t AliMUONTrackExtrap::GetAbsorberCorrectionParam(Double_t trackXYZIn[3], Do
   Double_t x0 = 0.;  // radiation-length (cm-1)
   Double_t atomicA = 0.; // A of material
   Double_t atomicZ = 0.; // Z of material
+  Double_t atomicZoverA = 0.; // Z/A of material
   Double_t localPathLength = 0;
   Double_t remainingPathLength = pathLength;
   Double_t zB = trackXYZIn[2];
@@ -570,9 +571,19 @@ Bool_t AliMUONTrackExtrap::GetAbsorberCorrectionParam(Double_t trackXYZIn[3], Do
     TGeoMaterial *material = currentnode->GetVolume()->GetMedium()->GetMaterial();
     rho = material->GetDensity();
     x0 = material->GetRadLen();
-    if (!material->IsMixture()) x0 /= rho; // different normalization in the modeler for mixture
     atomicA = material->GetA();
     atomicZ = material->GetZ();
+    if(material->IsMixture()){
+      TGeoMixture * mixture = (TGeoMixture*)material;
+      atomicZoverA = 0.;
+      Double_t sum = 0.;
+      for (Int_t iel=0;iel<mixture->GetNelements();iel++){
+       sum  += mixture->GetWmixt()[iel];
+       atomicZoverA += mixture->GetZmixt()[iel]*mixture->GetWmixt()[iel]/mixture->GetAmixt()[iel];
+      }
+      atomicZoverA/=sum;
+    }
+    else atomicZoverA = atomicZ/atomicA;
     
     // Get path length within this material
     gGeoManager->FindNextBoundary(remainingPathLength);
@@ -607,8 +618,8 @@ Bool_t AliMUONTrackExtrap::GetAbsorberCorrectionParam(Double_t trackXYZIn[3], Do
     f1 += (dzE*dzE - dzB*dzB) / b[2] / b[2] / x0 / 2.;
     f2 += (dzE*dzE*dzE - dzB*dzB*dzB) / b[2] / b[2] / b[2] / x0 / 3.;
     meanRho += localPathLength * rho;
-    totalELoss += BetheBloch(pTotal, localPathLength, rho, atomicA, atomicZ);
-    sigmaELoss2 += EnergyLossFluctuation2(pTotal, localPathLength, rho, atomicA, atomicZ);
+    totalELoss += BetheBloch(pTotal, localPathLength, rho, atomicZ, atomicZoverA);
+    sigmaELoss2 += EnergyLossFluctuation2(pTotal, localPathLength, rho, atomicZoverA);
     
     // prepare next step
     zB = zE;
@@ -883,7 +894,7 @@ Double_t AliMUONTrackExtrap::TotalMomentumEnergyLoss(AliMUONTrackParam* trackPar
 }
 
 //__________________________________________________________________________
-Double_t AliMUONTrackExtrap::BetheBloch(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ)
+Double_t AliMUONTrackExtrap::BetheBloch(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicZ, Double_t atomicZoverA)
 {
   /// Returns the mean total momentum energy loss of muon with total momentum='pTotal'
   /// in the absorber layer of lenght='pathLength', density='rho', A='atomicA' and Z='atomicZ'
@@ -894,11 +905,11 @@ Double_t AliMUONTrackExtrap::BetheBloch(Double_t pTotal, Double_t pathLength, Do
   if (atomicZ < 13) i = (12. * atomicZ + 7.) * 1.e-9;
   else i = (9.76 * atomicZ + 58.8 * TMath::Power(atomicZ,-0.19)) * 1.e-9;
   
-  return pathLength * rho * AliExternalTrackParam::BetheBlochGeant(pTotal/muMass, rho, 0.20, 3.00, i, atomicZ/atomicA);
+  return pathLength * rho * AliExternalTrackParam::BetheBlochGeant(pTotal/muMass, rho, 0.20, 3.00, i, atomicZoverA);
 }
 
 //__________________________________________________________________________
-Double_t AliMUONTrackExtrap::EnergyLossFluctuation2(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ)
+Double_t AliMUONTrackExtrap::EnergyLossFluctuation2(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicZoverA)
 {
   /// Returns the total momentum energy loss fluctuation of muon with total momentum='pTotal'
   /// in the absorber layer of lenght='pathLength', density='rho', A='atomicA' and Z='atomicZ'
@@ -908,7 +919,7 @@ Double_t AliMUONTrackExtrap::EnergyLossFluctuation2(Double_t pTotal, Double_t pa
   Double_t p2=pTotal*pTotal;
   Double_t beta2=p2/(p2 + muMass*muMass);
   
-  Double_t fwhm = 2. * k * rho * pathLength * atomicZ / atomicA / beta2; // FWHM of the energy loss Landau distribution
+  Double_t fwhm = 2. * k * rho * pathLength * atomicZoverA / beta2; // FWHM of the energy loss Landau distribution
   Double_t sigma2 = fwhm * fwhm / (8.*log(2.)); // gaussian: fwmh = 2 * srqt(2*ln(2)) * sigma (i.e. fwmh = 2.35 * sigma)
   
   //sigma2 = k * rho * pathLength * atomicZ / atomicA * eMass; // sigma2 of the energy loss gaussian distribution

diff --git a/MUON/AliMUONTrackExtrap.h b/MUON/AliMUONTrackExtrap.h
index 85a7042a223b355088afab776e3ad3911ba76d6a..9a9809ba93a1e7649f5d4d4ae8388f9404fde280 100644 (file)
--- a/MUON/AliMUONTrackExtrap.h
+++ b/MUON/AliMUONTrackExtrap.h
@@ -114,8 +114,8 @@ class AliMUONTrackExtrap : public TObject
                                            Double_t &pathLength, Double_t &f0, Double_t &f1, Double_t &f2,
                                            Double_t &meanRho, Double_t &totalELoss, Double_t &sigmaELoss2);
   
-  static Double_t BetheBloch(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ);
-  static Double_t EnergyLossFluctuation2(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ);
+  static Double_t BetheBloch(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicZ, Double_t atomicZoverA);
+  static Double_t EnergyLossFluctuation2(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicZoverA);
   
   static void Cov2CovP(const TMatrixD &param, TMatrixD &cov);
   static void CovP2Cov(const TMatrixD &param, TMatrixD &cov);

diff --git a/MUON/MUONRecoCheck.C b/MUON/MUONRecoCheck.C
index a149a188a32fb75a0f7a99eecf787a2d6a964afe..76ed952f317332ec97c442f1169d275d2664f563 100644 (file)
--- a/MUON/MUONRecoCheck.C
+++ b/MUON/MUONRecoCheck.C
@@ -197,8 +197,7 @@ void MUONRecoCheck (Int_t nEvent = -1, char* geoFilename = "geometry.root",
         
         //     printf(" Ref. track at vertex: x,y,z: %f %f %f px,py,pz,p: %f %f %f %f \n",x1,y1,z1,pX1,pY1,pZ1,p1);
         trackReco = trackOK;
-        trackParam = new AliMUONTrackParam(*((AliMUONTrackParam*)(trackReco->GetTrackParamAtCluster()->First())));
-        AliMUONTrackExtrap::ExtrapToVertex(trackParam,x1,y1,z1,0.,0.);
+       trackParam = new AliMUONTrackParam(*((AliMUONTrackParam*)(trackReco->GetTrackParamAtVertex())));
         x2 = trackParam->GetNonBendingCoor();
         y2 = trackParam->GetBendingCoor();
         z2 = trackParam->GetZ();