#include "AliMUONTrackParam.h"
#include "AliMUONConstants.h"
#include "AliMUONReconstructor.h"
-#include "AliMUONRecoParam.h"
#include "AliMagF.h"
-#include <TMath.h>
+#include <TGeoGlobalMagField.h>
#include <TGeoManager.h>
+#include <TMath.h>
#include <Riostream.h>
ClassImp(AliMUONTrackExtrap) // Class implementation in ROOT context
/// \endcond
-const AliMagF* AliMUONTrackExtrap::fgkField = 0x0;
+const Double_t AliMUONTrackExtrap::fgkSimpleBPosition = 0.5 * (AliMUONConstants::CoilZ() + AliMUONConstants::YokeZ());
+const Double_t AliMUONTrackExtrap::fgkSimpleBLength = 0.5 * (AliMUONConstants::CoilL() + AliMUONConstants::YokeL());
+ Double_t AliMUONTrackExtrap::fgSimpleBValue = 0.;
+ Bool_t AliMUONTrackExtrap::fgFieldON = kFALSE;
const Bool_t AliMUONTrackExtrap::fgkUseHelix = kFALSE;
const Int_t AliMUONTrackExtrap::fgkMaxStepNumber = 5000;
const Double_t AliMUONTrackExtrap::fgkHelixStepLength = 6.;
const Double_t AliMUONTrackExtrap::fgkRungeKuttaMaxResidue = 0.002;
+//__________________________________________________________________________
+void AliMUONTrackExtrap::SetField()
+{
+ // set field on/off flag
+ // set field at the centre of the dipole
+ const Double_t x[3] = {50.,50.,fgkSimpleBPosition};
+ Double_t b[3] = {0.,0.,0.};
+ TGeoGlobalMagField::Instance()->Field(x,b);
+ fgSimpleBValue = b[0];
+ fgFieldON = fgSimpleBValue ? kTRUE : kFALSE;
+
+}
+
//__________________________________________________________________________
Double_t AliMUONTrackExtrap::GetImpactParamFromBendingMomentum(Double_t bendingMomentum)
{
if (bendingMomentum == 0.) return 1.e10;
const Double_t kCorrectionFactor = 0.9; // impact parameter is 10% overestimated
- Double_t simpleBPosition = 0.5 * (AliMUONConstants::CoilZ() + AliMUONConstants::YokeZ());
- Double_t simpleBLength = 0.5 * (AliMUONConstants::CoilL() + AliMUONConstants::YokeL());
- Float_t b[3], x[3] = {50.,50.,(Float_t) simpleBPosition};
- if (fgkField) fgkField->Field(x,b);
- else {
- cout<<"F-AliMUONTrackExtrap::GetField: fgkField = 0x0"<<endl;
- exit(-1);
- }
- Double_t simpleBValue = (Double_t) b[0];
- return kCorrectionFactor * (-0.0003 * simpleBValue * simpleBLength * simpleBPosition / bendingMomentum);
+ return kCorrectionFactor * (-0.0003 * fgSimpleBValue * fgkSimpleBLength * fgkSimpleBPosition / bendingMomentum);
}
//__________________________________________________________________________
-Double_t AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(Double_t impactParam)
+Double_t
+AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(Double_t impactParam)
{
/// Returns signed bending momentum in bending plane (GeV/c),
/// the sign being the sign of the charge for particles moving forward in Z,
if (impactParam == 0.) return 1.e10;
const Double_t kCorrectionFactor = 1.1; // bending momentum is 10% underestimated
- Double_t simpleBPosition = 0.5 * (AliMUONConstants::CoilZ() + AliMUONConstants::YokeZ());
- Double_t simpleBLength = 0.5 * (AliMUONConstants::CoilL() + AliMUONConstants::YokeL());
- Float_t b[3], x[3] = {50.,50.,(Float_t) simpleBPosition};
- if (fgkField) fgkField->Field(x,b);
- else {
- cout<<"F-AliMUONTrackExtrap::GetField: fgkField = 0x0"<<endl;
- exit(-1);
- }
- Double_t simpleBValue = (Double_t) b[0];
- if (TMath::Abs(simpleBValue) > 0.01) return kCorrectionFactor * (-0.0003 * simpleBValue * simpleBLength * simpleBPosition / impactParam);
- else return AliMUONReconstructor::GetRecoParam()->GetMostProbBendingMomentum();
+ if (fgFieldON)
+ {
+ return kCorrectionFactor * (-0.0003 * fgSimpleBValue * fgkSimpleBLength * fgkSimpleBPosition / impactParam);
+ }
+ else
+ {
+ return AliMUONConstants::GetMostProbBendingMomentum();
+ }
}
//__________________________________________________________________________
-void AliMUONTrackExtrap::LinearExtrapToZ(AliMUONTrackParam* trackParam, Double_t zEnd)
+void AliMUONTrackExtrap::LinearExtrapToZ(AliMUONTrackParam* trackParam, Double_t zEnd, Bool_t updatePropagator)
{
/// Track parameters (and their covariances if any) linearly extrapolated to the plane at "zEnd".
/// On return, results from the extrapolation are updated in trackParam.
paramCov(2,3) += dZ * paramCov(3,3);
paramCov(3,2) = paramCov(2,3);
trackParam->SetCovariances(paramCov);
+
+ // Update the propagator if required
+ if (updatePropagator) {
+ TMatrixD jacob(5,5);
+ jacob.UnitMatrix();
+ jacob(0,1) = dZ;
+ jacob(2,3) = dZ;
+ trackParam->UpdatePropagator(jacob);
+ }
+
}
}
{
/// Interface to track parameter extrapolation to the plane at "Z" using Helix or Rungekutta algorithm.
/// On return, the track parameters resulting from the extrapolation are updated in trackParam.
- if (fgkUseHelix) AliMUONTrackExtrap::ExtrapToZHelix(trackParam,zEnd);
+ if (!fgFieldON) AliMUONTrackExtrap::LinearExtrapToZ(trackParam,zEnd);
+ else if (fgkUseHelix) AliMUONTrackExtrap::ExtrapToZHelix(trackParam,zEnd);
else AliMUONTrackExtrap::ExtrapToZRungekutta(trackParam,zEnd);
}
trackParam->SetNonBendingCoor(v3[0]); // X
trackParam->SetBendingCoor(v3[1]); // Y
trackParam->SetZ(v3[2]); // Z
- Double_t pYZ = v3[6] * TMath::Sqrt(1.0 - v3[3] * v3[3]);
+ Double_t pYZ = v3[6] * TMath::Sqrt((1.-v3[3])*(1.+v3[3]));
trackParam->SetInverseBendingMomentum(charge/pYZ);
trackParam->SetBendingSlope(v3[4]/v3[5]);
trackParam->SetNonBendingSlope(v3[3]/v3[5]);
if (trackParam->GetZ() == zEnd) return; // nothing to be done if same z
+ if (!fgFieldON) { // linear extrapolation if no magnetic field
+ AliMUONTrackExtrap::LinearExtrapToZ(trackParam,zEnd,updatePropagator);
+ return;
+ }
+
// No need to propagate the covariance matrix if it does not exist
if (!trackParam->CovariancesExist()) {
cout<<"W-AliMUONTrackExtrap::ExtrapToZCov: Covariance matrix does not exist"<<endl;
xyz[2] = vect[kiz] + 0.5 * step * vect[kipz];
//cmodif: call gufld (xyz, h) changed into:
- GetField (xyz, h);
+ TGeoGlobalMagField::Instance()->Field(xyz,h);
h2xy = h[0]*h[0] + h[1]*h[1];
h[3] = h[2]*h[2]+ h2xy;
rest = step - tl;
if (TMath::Abs(h) > TMath::Abs(rest)) h = rest;
//cmodif: call gufld(vout,f) changed into:
-
- GetField(vout,f);
+ TGeoGlobalMagField::Instance()->Field(vout,f);
// *
// * start of integration
xyzt[2] = zt;
//cmodif: call gufld(xyzt,f) changed into:
- GetField(xyzt,f);
+ TGeoGlobalMagField::Instance()->Field(xyzt,f);
at = a + secxs[0];
bt = b + secys[0];
xyzt[2] = zt;
//cmodif: call gufld(xyzt,f) changed into:
- GetField(xyzt,f);
+ TGeoGlobalMagField::Instance()->Field(xyzt,f);
z = z + (c + (seczs[0] + seczs[1] + seczs[2]) * kthird) * h;
y = y + (b + (secys[0] + secys[1] + secys[2]) * kthird) * h;
return;
}
-//___________________________________________________________
-void AliMUONTrackExtrap::GetField(Double_t *Position, Double_t *Field)
-{
- /// interface for arguments in double precision (Why ? ChF)
- Float_t x[3], b[3];
-
- x[0] = Position[0]; x[1] = Position[1]; x[2] = Position[2];
-
- if (fgkField) fgkField->Field(x,b);
- else {
- cout<<"F-AliMUONTrackExtrap::GetField: fgkField = 0x0"<<endl;
- exit(-1);
- }
-
- Field[0] = b[0]; Field[1] = b[1]; Field[2] = b[2];
-
- return;
-}