]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Moved old AliMagFCheb to AliMagF, small fixes/optimizations in field classes
authorfca <fca@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sun, 1 Feb 2009 11:38:48 +0000 (11:38 +0000)
committerfca <fca@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sun, 1 Feb 2009 11:38:48 +0000 (11:38 +0000)
STEER/AliCheb3DCalc.h
STEER/AliMagF.cxx
STEER/AliMagF.h
STEER/AliMagWrapCheb.cxx
STEER/AliMagWrapCheb.h
STEER/CMake_libSTEER.txt
STEER/CMake_libSTEERBase.txt
STEER/STEERBaseLinkDef.h
STEER/STEERLinkDef.h
STEER/libSTEER.pkg
STEER/libSTEERBase.pkg

index 8cb0293d609322812e7336baa2c6c7577ab26a09..c6f54b284bae56584283d232bdd6284fc43977a5 100644 (file)
@@ -94,6 +94,7 @@ class AliCheb3DCalc: public TNamed
 inline Float_t AliCheb3DCalc::ChebEval1D(Float_t  x, const Float_t * array, int ncf ) 
 {
   // evaluate 1D Chebyshev parameterization. x is the argument mapped to [-1:1] interval
+  if (ncf<=0) return 0;
   Float_t b0, b1, b2, x2 = x+x;
   b0 = array[--ncf]; 
   b1 = b2 = 0;
index 14a73562045b6d982669f5ae48fabbb535337a0e..0776ff5c7a236f39012b7b0b5bbab37046944870 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-/* $Id$ */
 
-//----------------------------------------------------------------------
-// Basic magnetic field class
-// Used in all the detectors, and also in the traking classes
-// Author:
-//----------------------------------------------------------------------
+#include <TClass.h>
+#include <TFile.h>
+#include <TSystem.h>
 
-#include "AliLog.h"
 #include "AliMagF.h"
+#include "AliMagWrapCheb.h"
+#include "AliLog.h"
 
 ClassImp(AliMagF)
 
+const Double_t AliMagF::fgkSol2DipZ   =  -700.;
+const Double_t AliMagF::fgkBMachineZ1 =   919.;
+const Double_t AliMagF::fgkBMachineZ2 = -1972.;
+
 //_______________________________________________________________________
 AliMagF::AliMagF():
-  fMap(0),
-  fType(0),
+  TVirtualMagField(),
+  fMeasuredMap(0),
+  fMapType(k5kG),
+  fSolenoid(0),
+  fBeamType(kNoBeamField),
+  fBeamEnergy(0),
+  fCompensator(kFALSE),
+  //
   fInteg(0),
-  fPrecInteg(1),
-  fFactor(0),
-  fMax(0),
-  fReadField(1)
-{
+  fPrecInteg(0),
+  fFactorSol(1.),
+  fFactorDip(1.),
+  fMax(15),
+  fDipoleOFF(kFALSE),
   //
+  fQuadGradient(0),
+  fDipoleField(0),
+  fCCorrField(0), 
+  fACorr1Field(0),
+  fACorr2Field(0),
+  fParNames("","")
+{
   // Default constructor
   //
 }
 
 //_______________________________________________________________________
-AliMagF::AliMagF(const char *name, const char *title, Int_t integ, 
-                 Float_t factor, Float_t fmax):
-  TNamed(name,title),
-  fMap(0),
-  fType(0),
-  fInteg(0),
+AliMagF::AliMagF(const char *name, const char* title, Int_t integ, 
+                Double_t factorSol, Double_t factorDip, 
+                Double_t fmax, BMap_t maptype, const char* path,
+                BeamType_t bt, Double_t be, Bool_t compensator):
+  TVirtualMagField(name),
+  fMeasuredMap(0),
+  fMapType(maptype),
+  fSolenoid(0),
+  fBeamType(bt),
+  fBeamEnergy(be),
+  fCompensator(compensator),
+  //
+  fInteg(integ),
   fPrecInteg(1),
-  fFactor(factor),
+  fFactorSol(1.),
+  fFactorDip(1.),
   fMax(fmax),
-  fReadField(1)
+  fDipoleOFF(factorDip==0.),
+  //
+  fQuadGradient(0),
+  fDipoleField(0),
+  fCCorrField(0), 
+  fACorr1Field(0),
+  fACorr2Field(0),
+  fParNames("","")
 {
   //
-  // Standard constructor
+  SetTitle(title);
+  if(integ<0 || integ > 2) {
+    AliWarning(Form("Invalid magnetic field flag: %5d; Helix tracking chosen instead",integ));
+    fInteg = 2;
+  }
+  if (fInteg == 0) fPrecInteg = 0;
   //
-    if(integ<0 || integ > 2) {
-      AliWarning(Form(
-              "Invalid magnetic field flag: %5d; Helix tracking chosen instead"
-              ,integ));
-      fInteg = 2;
-    } else {
-      fInteg = integ;
-    }
-   
-    if (fInteg == 0) fPrecInteg = 0;
-    
-    fType = kUndef;
-    //
+  const char* parname = 0;
+  //  
+  if (fMapType == k2kG) {
+    fSolenoid = 2.;
+    parname = fDipoleOFF ? "Sol12_Dip0_Hole":"Sol12_Dip6_Hole";
+  } else if (fMapType == k5kG) {
+    fSolenoid = 5.;
+    parname = fDipoleOFF ? "Sol30_Dip0_Hole":"Sol30_Dip6_Hole";
+  } else if (fMapType == k5kGUniform) {
+    fSolenoid = 5.;
+    parname = "Sol30_Dip6_Uniform";
+  } else {
+    AliFatal(Form("Unknown field identifier %d is requested\n",fMapType)); 
+  }
+  //
+  SetDataFileName(path);
+  SetParamName(parname);
+  //
+  SetFactorSol(factorSol);
+  SetFactorDip(factorDip);
+  LoadParameterization();
+  InitMachineField(fBeamType,fBeamEnergy);
 }
 
 //_______________________________________________________________________
 AliMagF::AliMagF(const AliMagF &src):
-  TNamed(src),
-  fMap(src.fMap),
-  fType(src.fType),
+  TVirtualMagField(src),
+  fMeasuredMap(0),
+  fMapType(src.fMapType),
+  fSolenoid(src.fSolenoid),
+  fBeamType(src.fBeamType),
+  fBeamEnergy(src.fBeamEnergy),
+  fCompensator(src.fCompensator),
   fInteg(src.fInteg),
   fPrecInteg(src.fPrecInteg),
-  fFactor(src.fFactor),
+  fFactorSol(src.fFactorSol),
+  fFactorDip(src.fFactorDip),
   fMax(src.fMax),
-  fReadField(src.fReadField)
+  fDipoleOFF(src.fDipoleOFF),
+  fQuadGradient(src.fQuadGradient),
+  fDipoleField(src.fDipoleField),
+  fCCorrField(src.fCCorrField), 
+  fACorr1Field(src.fACorr1Field),
+  fACorr2Field(src.fACorr2Field),
+  fParNames(src.fParNames)
 {
-    // Copy constructor
+  if (src.fMeasuredMap) fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
 }
 
 //_______________________________________________________________________
-void AliMagF::Field(const Float_t*, Float_t *b) const
+AliMagF::~AliMagF()
 {
+  delete fMeasuredMap;
+}
+
+//_______________________________________________________________________
+Bool_t AliMagF::LoadParameterization()
+{
+  if (fMeasuredMap) {
+    AliError(Form("Field data %s are already loaded from %s\n",GetParamName(),GetDataFileName()));
+    return kTRUE;
+  }
   //
-  // Method to return the field in one point -- dummy in this case
+  char* fname = gSystem->ExpandPathName(GetDataFileName());
+  TFile* file = TFile::Open(fname);
+  if (!file) {
+    AliError(Form("Failed to open magnetic field data file %s\n",fname)); 
+    return kFALSE;
+  }
   //
-  AliWarning("Undefined MagF Field called, returning 0");
-  b[0]=b[1]=b[2]=0;
+  fMeasuredMap = dynamic_cast<AliMagWrapCheb*>(file->Get(GetParamName()));
+  if (!fMeasuredMap) {
+    AliError(Form("Did not find field %s in %s\n",GetParamName(),fname)); 
+    return kFALSE;
+  }
+  file->Close();
+  delete file;
+  return kTRUE;
 }
 
+
 //_______________________________________________________________________
-void AliMagF::Field(const double*, double *b) const
+void AliMagF::Field(const Double_t *xyz, Double_t *b)
 {
+  // Method to calculate the field at point  xyz
   //
-  // Method to return the field in one point -- dummy in this case
+  b[0]=b[1]=b[2]=0.0;
+  if (xyz[2] > fgkBMachineZ1 || xyz[2] < fgkBMachineZ2) MachineField(xyz, b);
+  else if (fMeasuredMap) {
+    fMeasuredMap->Field(xyz,b);
+    if (xyz[2]>fgkSol2DipZ || fDipoleOFF) for (int i=3;i--;) b[i] *= fFactorSol;
+    else                                  for (int i=3;i--;) b[i] *= fFactorDip;
+  }
   //
-  AliWarning("Undefined MagF Field called, returning 0");
-  b[0]=b[1]=b[2]=0;
 }
 
 //_______________________________________________________________________
-void AliMagF::GetTPCInt(const Float_t *, Float_t *b) const
+Double_t AliMagF::GetBz(const Double_t *xyz) const
 {
-//
-// Obtain the integral of the field components in the TPC from given point
-// to the closest cathod plane
-//
-  AliWarning("Undefined MagF TPCIntegral called, returning 0");
-  b[0]=b[1]=b[2]=0;
+  // Method to calculate the field at point  xyz
+  //
+  if (xyz[2] <= fgkBMachineZ1 && xyz[2] >= fgkBMachineZ2) return fMeasuredMap->GetBz(xyz);
+  else {
+    double b[3] = {0,0,0};
+    MachineField(xyz, b);
+    return b[2];
+  }
+  //
 }
 
 //_______________________________________________________________________
-void AliMagF::GetTPCIntCyl(const Float_t *, Float_t *b) const
+AliMagF& AliMagF::operator=(const AliMagF& src)
 {
-//    
-// Obtain the integral of the field components in the TPC from given point
-// to the closest cathod plane
-//
-  AliWarning("Undefined MagF TPCIntegral called, returning 0");
-  b[0]=b[1]=b[2]=0;
+  if (this != &src && src.fMeasuredMap) { 
+    if (fMeasuredMap) delete fMeasuredMap;
+    fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
+    SetName(src.GetName());
+    fSolenoid    = src.fSolenoid;
+    fBeamType    = src.fBeamType;
+    fBeamEnergy  = src.fBeamEnergy;
+    fCompensator = src.fCompensator;
+    fInteg       = src.fInteg;
+    fPrecInteg   = src.fPrecInteg;
+    fFactorSol   = src.fFactorSol;
+    fFactorDip   = src.fFactorDip;
+    fMax         = src.fMax;
+    fDipoleOFF   = src.fDipoleOFF;
+    fParNames    = src.fParNames;
+  }
+  return *this;
 }
 
 //_______________________________________________________________________
-AliMagF& AliMagF::operator=(const AliMagF& rhs)
+void AliMagF::InitMachineField(BeamType_t btype, Double_t benergy)
 {
-    // Asignment operator
-    fMap       = rhs.fMap;
-    fType      = rhs.fType;
-    fInteg     = rhs.fInteg;
-    fPrecInteg = rhs.fPrecInteg;
-    fFactor    = rhs.fFactor;
-    fMax       = rhs.fMax;
-    fReadField = rhs.fReadField;
-    return *this;
+  const double kToler = 0.1;
+  if (btype==kNoBeamField) {
+    fQuadGradient = fDipoleField = fCCorrField = fACorr1Field = fACorr2Field = 0.;
+  }
+  //
+  else if (btype==kBeamTypepp && TMath::Abs(1.-benergy/5000.)<kToler ){
+    // p-p @ 5+5 TeV
+    fQuadGradient = 15.7145;
+    fDipoleField  = 27.0558;
+    // SIDE C
+    fCCorrField   = 9.7017;
+    // SIDE A
+    fACorr1Field  = -13.2143;
+    fACorr2Field  = -11.9909;
+  } else if (btype == kBeamTypepp && TMath::Abs(1.-benergy/450.)<kToler) {
+    // p-p 0.45+0.45 TeV
+    Double_t const kEnergyRatio = benergy / 7000.;
+    fQuadGradient = 22.0002 * kEnergyRatio;
+    fDipoleField  = 37.8781 * kEnergyRatio;
+    // SIDE C
+    fCCorrField   =  9.6908;
+    // SIDE A
+    fACorr1Field  = -13.2014;
+    fACorr2Field  = -9.6908;
+  } else if ( (btype == kBeamTypepp && TMath::Abs(1.-benergy/7000.)<kToler) || 
+             (fBeamType == kBeamTypeAA) ) {
+    // Pb-Pb @ 2.7+2.7 TeV or p-p @ 7+7 TeV
+    fQuadGradient = 22.0002;
+    fDipoleField  = 37.8781;
+    // SIDE C
+    fCCorrField   = 9.6908;
+    // SIDE A
+    fACorr1Field  = -13.2014;
+    fACorr2Field  = -9.6908;
+  }
+  //
 }
 
-void AliMagF::SetPrecInteg(Int_t integ)
+//_______________________________________________________________________
+void AliMagF::MachineField(const Double_t *x, Double_t *b) const
 {
-    if (fInteg > 0) {
-       fPrecInteg = integ;
+  // ---- This is the ZDC part
+  const Double_t kCCorrBegin = fgkBMachineZ2-0.5,kCCorrEnd = kCCorrBegin - 153., kCCorrSqRadius = 4.5*4.5;
+  //
+  const Double_t kCTripletBegin  = -2296.5;
+  const Double_t kCQ1Begin = kCTripletBegin,        kCQ1End = kCQ1Begin-637., kCQ1SqRadius = 3.5*3.5;
+  const Double_t kCQ2Begin = kCTripletBegin-908.5,  kCQ2End = kCQ2Begin-550., kCQ2SqRadius = 3.5*3.5;
+  const Double_t kCQ3Begin = kCTripletBegin-1558.5, kCQ3End = kCQ3Begin-550., kCQ3SqRadius = 3.5*3.5;
+  const Double_t kCQ4Begin = kCTripletBegin-2400.,  kCQ4End = kCQ4Begin-637., kCQ4SqRadius = 3.5*3.5;
+  //
+  const Double_t kCD1Begin = -5838.3,  kCD1End = kCD1Begin-945., kCD1SqRadius = 4.5*4.5;
+  const Double_t kCD2Begin = -12167.8, kCD2End = kCD2Begin-945., kCD2SqRadius = 4.5*4.5;
+  const Double_t kCD2XCentre1 = -9.7;
+  const Double_t kCD2XCentre2 =  9.7;
+  //
+  // -> SIDE A
+  // NB -> kACorr1Begin = 919. to be checked
+  const Double_t kACorr1Begin = fgkBMachineZ1, kACorr1End = kACorr1Begin+260., kCCorr1SqRadius = 4.*4.;
+  const Double_t kACorr2Begin = -fgkBMachineZ2 + 0.5, kACorr2End = kACorr2Begin+153., kCCorr2SqRadius = 4.5*4.5;
+  const Double_t kATripletBegin  = 2296.5;
+  const Double_t kAQ1Begin = kATripletBegin,   kAQ1End = kAQ1Begin+637., kAQ1SqRadius = 3.5*3.5;
+  const Double_t kAQ2Begin = kATripletBegin+908.5,  kAQ2End = kAQ2Begin+550., kAQ2SqRadius = 3.5*3.5;
+  const Double_t kAQ3Begin = kATripletBegin+1558.5, kAQ3End = kAQ3Begin+550., kAQ3SqRadius = 3.5*3.5;
+  const Double_t kAQ4Begin = kATripletBegin+2400.,  kAQ4End = kAQ4Begin+637., kAQ4SqRadius = 3.5*3.5;
+  //
+  const Double_t kAD1Begin = 5838.3,  kAD1End = kAD1Begin+945., kAD1SqRadius = 3.375*3.375;
+  const Double_t kAD2Begin = 12167.8, kAD2End = kAD2Begin+945., kAD2SqRadius = 3.75*3.75;
+  const Double_t kAD2XCentre1 = -9.4;
+  const Double_t kAD2XCentre2 =  9.4;
+  //
+  double rad2 = x[0] * x[0] + x[1] * x[1];
+  //
+  // SIDE C **************************************************
+  if(x[2]<0.){  
+    if(x[2] < kCCorrBegin && x[2] > kCCorrEnd && rad2 < kCCorrSqRadius){
+      b[0] = fCCorrField;
+      b[1] = 0.;
+      b[2] = 0.;
+    } 
+    else if(x[2] < kCQ1Begin && x[2] > kCQ1End && rad2 < kCQ1SqRadius){
+      b[0] = fQuadGradient*x[1];
+      b[1] = fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] < kCQ2Begin && x[2] > kCQ2End && rad2 < kCQ2SqRadius){
+      b[0] = -fQuadGradient*x[1];
+      b[1] = -fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] < kCQ3Begin && x[2] > kCQ3End && rad2 < kCQ3SqRadius){
+      b[0] = -fQuadGradient*x[1];
+      b[1] = -fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] < kCQ4Begin && x[2] > kCQ4End && rad2 < kCQ4SqRadius){
+      b[0] = fQuadGradient*x[1];
+      b[1] = fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] < kCD1Begin && x[2] > kCD1End && rad2 < kCD1SqRadius){
+      b[1] = fDipoleField;
+      b[2] = 0.;
+      b[2] = 0.;
+    }
+    else if(x[2] < kCD2Begin && x[2] > kCD2End){
+      if(((x[0]-kCD2XCentre1)*(x[0]-kCD2XCentre1)+(x[1]*x[1]))<kCD2SqRadius
+        || ((x[0]-kCD2XCentre2)*(x[0]-kCD2XCentre2)+(x[1]*x[1]))<kCD2SqRadius){
+       b[1] = -fDipoleField;
+       b[2] = 0.;
+       b[2] = 0.;
+      }
+    }
+  }
+  //
+  // SIDE A **************************************************
+  else{        
+    if(fCompensator && (x[2] > kACorr1Begin && x[2] < kACorr1End) && rad2 < kCCorr1SqRadius) {
+      // Compensator magnet at z = 1075 m 
+      b[0] = fACorr1Field;
+      b[1] = 0.;
+      b[2] = 0.;
+      return;
+    }
+    //
+    if(x[2] > kACorr2Begin && x[2] < kACorr2End && rad2 < kCCorr2SqRadius){
+      b[0] = fACorr2Field;
+      b[1] = 0.;
+      b[2] = 0.;
+    }          
+    else if(x[2] > kAQ1Begin && x[2] < kAQ1End && rad2 < kAQ1SqRadius){
+      // First quadrupole of inner triplet de-focussing in x-direction
+      b[0] = -fQuadGradient*x[1];
+      b[1] = -fQuadGradient*x[0];
+      b[2] = 0.;
     }
-    else if (integ != 0)
-    {
-       AliWarning("Precision integration flag set to 0 \n");
-       fPrecInteg = 0;
+    else if(x[2] > kAQ2Begin && x[2] < kAQ2End && rad2 < kAQ2SqRadius){
+      b[0] = fQuadGradient*x[1];
+      b[1] = fQuadGradient*x[0];
+      b[2] = 0.;
     }
+    else if(x[2] > kAQ3Begin && x[2] < kAQ3End && rad2 < kAQ3SqRadius){
+      b[0] = fQuadGradient*x[1];
+      b[1] = fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] > kAQ4Begin && x[2] < kAQ4End && rad2 < kAQ4SqRadius){
+      b[0] = -fQuadGradient*x[1];
+      b[1] = -fQuadGradient*x[0];
+      b[2] = 0.;
+    }
+    else if(x[2] > kAD1Begin && x[2] < kAD1End && rad2 < kAD1SqRadius){
+      b[0] = 0.;
+      b[1] = -fDipoleField;
+      b[2] = 0.;
+    }
+    else if(x[2] > kAD2Begin && x[2] < kAD2End){
+      if(((x[0]-kAD2XCentre1)*(x[0]-kAD2XCentre1)+(x[1]*x[1])) < kAD2SqRadius
+        || ((x[0]-kAD2XCentre2)*(x[0]-kAD2XCentre2)+(x[1]*x[1])) < kAD2SqRadius){
+       b[1] = fDipoleField;
+      }
+    }
+  }
+}
+
+//_______________________________________________________________________
+void AliMagF::GetTPCInt(const Double_t *xyz, Double_t *b) const
+{
+  // Method to calculate the integral of magnetic integral from xyz to nearest cathode plane
+  b[0]=b[1]=b[2]=0.0;
+  if (fMeasuredMap) {
+    fMeasuredMap->GetTPCInt(xyz,b);
+    for (int i=3;i--;) b[i] *= fFactorSol;
+  }
+}
+
+//_______________________________________________________________________
+void AliMagF::GetTPCIntCyl(const Double_t *rphiz, Double_t *b) const
+{
+  // Method to calculate the integral of magnetic integral from point to nearest cathode plane
+  // in cylindrical coordiates ( -pi<phi<pi convention )
+  b[0]=b[1]=b[2]=0.0;
+  if (fMeasuredMap) {
+    fMeasuredMap->GetTPCIntCyl(rphiz,b);
+    for (int i=3;i--;) b[i] *= fFactorSol;
+  }
 }
index 3822d477e2aa8b8b64ffdb023bfd0f1fc16874d9..dd5ecc4032ac050e6175a84e0297fda62932e80f 100644 (file)
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
-/* $Id$ */
+//
+// Interface between the TVirtualMagField and AliMagWrapCheb: wrapper to
+// the set of magnetic field data + Tosca parameterization by 
+// Chebyshev polynomials
+// 
+// Author: ruben.shahoyan@cern.ch
+//
 
-//----------------------------------------------------------------------
-// Basic magnetic field class
-// Used in all the detectors, and also in the traking classes
-// Author:
-//----------------------------------------------------------------------
+//#include <TGeoGlobalMagField.h>
+#include <TVirtualMagField.h>
+class AliMagWrapCheb;
 
-#include "TNamed.h"
-
-enum Field_t {kUndef=1, kConst=1, kConMesh=2, kDipoMap=3};
-
-class AliMagF : public TNamed {
-
-public:
+class AliMagF : public TVirtualMagField
+{
+ public:
+  enum BMap_t      {k2kG, k5kG, k5kGUniform};
+  enum BeamType_t  {kBeamTypeAA, kBeamTypepp, kNoBeamField};
+  //
   AliMagF();
-  AliMagF(const char *name, const char *title, Int_t integ, 
-         Float_t factor = 1., Float_t fmax = 10.);
-  AliMagF(const AliMagF& maps);
-  virtual ~AliMagF() {}
-  AliMagF& operator=(const AliMagF& rhs);
-  virtual void    Field(const Float_t *x, Float_t *b)                  const;
-  virtual void    Field(const Double_t *x, Double_t *b)                const;
-  virtual void    GetTPCInt(const Float_t *xyz, Float_t *b)        const;
-  virtual void    GetTPCIntCyl(const Float_t *rphiz, Float_t *b)   const;
-  virtual Int_t   Type() const {return fType;}
-  virtual Float_t Max() const {return fMax;}
-  virtual Int_t   Map() const {return fMap;}
-  virtual Int_t   Integ() const {return fInteg;}
-  virtual Int_t   PrecInteg() const {return fPrecInteg;}  
-  virtual Float_t Factor() const {return fFactor;}
-  virtual void    ReadField() {}
-  virtual Float_t SolenoidField() const {return 2.;}
-  virtual void    SetPrecInteg(Int_t integ);
-  virtual void    SetReadField(Bool_t flag = kTRUE) {fReadField = flag;}
+  AliMagF(const char *name, const char* title, Int_t integ, 
+         Double_t factorSol=1., Double_t factorDip=1., 
+         Double_t fmax=15, BMap_t maptype = k5kG,
+         const char* path="$(ALICE_ROOT)/data/maps/mfchebKGI_sym.root",
+         BeamType_t btype=kBeamTypepp, Double_t benergy=7000., Bool_t compensator=kFALSE);
+  AliMagF(const AliMagF& src);             
+  AliMagF& operator=(const AliMagF& src);
+  virtual ~AliMagF();
+  //
+  virtual void Field(const Double_t *x, Double_t *b);
+  void       GetTPCInt(const Double_t *xyz, Double_t *b)        const;
+  void       GetTPCIntCyl(const Double_t *rphiz, Double_t *b)   const;
+  Double_t   GetBz(const Double_t *xyz)                         const;
+  //
+  AliMagWrapCheb* GetMeasuredMap()                              const {return fMeasuredMap;}
+  //
+  // former AliMagF methods or their aliases
+  void         SetFactorSol(Float_t fc=1.)                            {fFactorSol = fc;}
+  void         SetFactorDip(Float_t fc=1.)                            {fFactorDip = fc;}
+  Double_t     GetFactorSol()                                   const {return fFactorSol;}
+  Double_t     GetFactorDip()                                   const {return fFactorSol;}
+  Double_t     Factor()                                         const {return GetFactorSol();}
+  Bool_t       IsUniform()                                      const {return fMapType == k5kGUniform;}
+  //
+  void         MachineField(const Double_t  *x, Double_t *b)    const;
+  BMap_t       GetMapType()                                     const {return fMapType;}
+  Bool_t       GetCompensator()                                 const {return fCompensator;}
+  BeamType_t   GetBeamType()                                    const {return fBeamType;}
+  Double_t     GetBeamEnergy()                                  const {return fBeamEnergy;}
+  Double_t     Max()                                            const {return fMax;}
+  Int_t        Integ()                                          const {return fInteg;}
+  Int_t        PrecInteg()                                      const {return fPrecInteg;}  
+  Double_t     SolenoidField()                                  const {return -fFactorSol*fSolenoid;}
+  //
+  Char_t*      GetDataFileName()                                const {return (Char_t*)fParNames.GetName();}
+  Char_t*      GetParamName()                                   const {return (Char_t*)fParNames.GetTitle();}
+  void         SetDataFileName(const Char_t* nm)                      {fParNames.SetName(nm);}
+  void         SetParamName(const Char_t* nm)                         {fParNames.SetTitle(nm);}
+  //
+  Bool_t       LoadParameterization();
+  //
  protected:
-  Int_t     fMap;       // Field Map identifier
-  Int_t     fType;      // Mag Field type
-  Int_t     fInteg;     // Default integration method as indicated in Geant
-  Int_t     fPrecInteg; // Alternative integration method, e.g. for higher precision
-  Float_t   fFactor;    // Multiplicative factor
-  Float_t   fMax;       // Max Field as indicated in Geant
-  Bool_t    fReadField; // Flag for reading the field from file (if available) 
-  ClassDef(AliMagF,5)   //Base class for all Alice MagField
+  // not supposed to be changed during the run, set only at the initialization via constructor
+  void         InitMachineField(BeamType_t btype, Double_t benergy);
+  void         SetBeamType(BeamType_t type)                           {fBeamType = type;}
+  void         SetBeamEnergy(Float_t energy)                          {fBeamEnergy = energy;}
+  void         SetCompensatorMagnet(Bool_t flag)                      {fCompensator = flag;}
+  //
+ protected:
+  AliMagWrapCheb*  fMeasuredMap;     //! Measured part of the field map
+  BMap_t           fMapType;         // field map type
+  Double_t         fSolenoid;        // Solenoid field setting
+  BeamType_t       fBeamType;        // Beam type: A-A (fBeamType=0) or p-p (fBeamType=1)
+  Double_t         fBeamEnergy;      // Beam energy in GeV
+  Bool_t           fCompensator;     // Flag for compensator magnetic field (kTrue -> ON)
+  // 
+  Int_t            fInteg;           // Default integration method as indicated in Geant
+  Int_t            fPrecInteg;       // Alternative integration method, e.g. for higher precision
+  Double_t         fFactorSol;       // Multiplicative factor for solenoid
+  Double_t         fFactorDip;       // Multiplicative factor for dipole
+  Double_t         fMax;             // Max Field as indicated in Geant
+  Bool_t           fDipoleOFF;       // Dipole ON/OFF flag
+  //
+  Double_t         fQuadGradient;    // Gradient field for inner triplet quadrupoles
+  Double_t         fDipoleField;     // Field value for D1 and D2 dipoles
+  Double_t         fCCorrField;      // Side C 2nd compensator field
+  Double_t         fACorr1Field;     // Side A 1st compensator field 
+  Double_t         fACorr2Field;     // Side A 2nd compensator field
+  //
+  TNamed           fParNames;        // file and parameterization loadad
+  //
+  static const Double_t  fgkSol2DipZ;    // conventional Z of transition from L3 to Dipole field 
+  static const Double_t  fgkBMachineZ1;  // Min Z of the LHC mag field range (to be checked?)
+  static const Double_t  fgkBMachineZ2;  // Max Z of the LHC mag field range
+  //   
+  ClassDef(AliMagF, 1)           // Class for all Alice MagField wrapper for measured data + Tosca parameterization
 };
 
+
 #endif
index 05aafcf63b89aa0d4bb95c1a6684117ce4fae095..05e4b04d33aebc91992b85e02f83bf9015557093 100644 (file)
@@ -239,11 +239,11 @@ void AliMagWrapCheb::Clear(const Option_t *)
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::Field(const float *xyz, float *b) const
+void AliMagWrapCheb::Field(const Double_t *xyz, Double_t *b) const
 {
   // compute field in cartesian coordinates. If point is outside of the parameterized region
   // get it at closest valid point
-  static float rphiz[3];
+  static Double_t rphiz[3];
   //
 #ifndef _BRING_TO_BOUNDARY_  // exact matching to fitted volume is requested
   if ( !(xyz[2]>=GetMinZSol()&&xyz[2]<=GetMaxZSol()) && 
@@ -274,46 +274,100 @@ void AliMagWrapCheb::Field(const float *xyz, float *b) const
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::Field(const Double_t *xyz, Double_t *b) const
+Double_t AliMagWrapCheb::GetBz(const Double_t *xyz) const
 {
-  // compute field in cartesian coordinates. If point is outside of the parameterized region
+  // compute Bz for the point in cartesian coordinates. If point is outside of the parameterized region
   // get it at closest valid point
   static Double_t rphiz[3];
   //
 #ifndef _BRING_TO_BOUNDARY_  // exact matching to fitted volume is requested
   if ( !(xyz[2]>=GetMinZSol()&&xyz[2]<=GetMaxZSol()) && 
-       !(xyz[2]>=GetMinZDip()&&xyz[2]<=GetMaxZDip())  ) {for (int i=3;i--;) b[i]=0; return;}
+       !(xyz[2]>=GetMinZDip()&&xyz[2]<=GetMaxZDip())  ) return 0;
 #endif
   //
   if (xyz[2]<fMaxZDip) {    // dipole part?
 #ifndef _BRING_TO_BOUNDARY_
     AliCheb3D* par = GetParamDip(FindDipSegment(xyz));
-    if (par->IsInside(xyz)) {par->Eval(xyz,b); return;}
-    for (int i=3;i--;) b[i]=0; return;
+    if (par->IsInside(xyz)) return par->Eval(xyz,2);
+    else return 0;
 #else
-    GetParamDip(FindDipSegment(xyz))->Eval(xyz,b); return;  
+    return GetParamDip(FindDipSegment(xyz))->Eval(xyz,2);
 #endif
   }
-  //
   // Sol region: convert coordinates to cyl system
   CartToCyl(xyz,rphiz);
-#ifndef _BRING_TO_BOUNDARY_
-  if (rphiz[0]>GetMaxRSol()) {for (int i=3;i--;) b[i]=0; return;}
-#endif
+  return FieldCylSolBz(rphiz);
+}
+
+
+//__________________________________________________________________________________________
+void AliMagWrapCheb::Print(Option_t *) const
+{
+  // print info
+  printf("Alice magnetic field parameterized by Chebyshev polynomials\n");
+  printf("Segmentation for Solenoid (%+.2f<Z<%+.2f cm | R<%.2f cm)\n",fMinZSol,fMaxZSol,fMaxRSol);
   //
-  FieldCylSol(rphiz,b);
+  if (fParamsSol) fParamsSol->Print();
+  /*
+  for (int iz=0;iz<fNSegZSol;iz++) {
+    AliCheb3D* param = GetParamSol( fSegZIdSol[iz] );
+    printf("*** Z Segment %2d (%+7.2f<Z<%+7.2f)\t***\n",iz,param->GetBoundMin(2),param->GetBoundMax(2));
+    for (int ir=0;ir<fNSegRSol[iz];ir++) {
+      param = GetParamSol( fSegZIdSol[iz]+ir );
+      printf("    R Segment %2d (%+7.2f<R<%+7.2f, Precision: %.1e) (ID=%2d)\n",ir, param->GetBoundMin(0),
+            param->GetBoundMax(0),param->GetPrecision(),fSegZIdSol[iz]+ir);
+    }
+  }
+  */
   //
-  // convert field to cartesian system
-  CylToCartCylB(rphiz, b,b);
+  printf("Segmentation for TPC field integral (%+.2f<Z<%+.2f cm | R<%.2f cm)\n",fMinZTPCInt,fMaxZTPCInt,fMaxRTPCInt);
   //
+  if (fParamsTPCInt) fParamsTPCInt->Print();
+  /*
+  for (int iz=0;iz<fNSegZTPCInt;iz++) {
+    AliCheb3D* param = GetParamTPCInt( fSegZIdTPCInt[iz] );
+    printf("*** Z Segment %2d (%+7.2f<Z<%+7.2f)\t***\n",iz,param->GetBoundMin(2),param->GetBoundMax(2));
+    for (int ir=0;ir<fNSegRTPCInt[iz];ir++) {
+      param = GetParamTPCInt( fSegZIdTPCInt[iz]+ir );
+      printf("    R Segment %2d (%+7.2f<R<%+7.2f, Precision: %.1e) (ID=%2d)\n",ir, param->GetBoundMin(0),
+            param->GetBoundMax(0),param->GetPrecision(),fSegZIdTPCInt[iz]+ir);
+    }
+  }
+  */
+  //
+  printf("Segmentation for Dipole (%+.2f<Z<%+.2f cm)\n",fMinZDip,fMaxZDip);
+  if (fParamsDip) fParamsDip->Print();
+  //
+  
+
+}
+
+
+//__________________________________________________________________________________________________
+Int_t    AliMagWrapCheb::FindDipSegment(const Double_t *xyz) const 
+{
+  // find the segment containing point xyz. If it is outside find the closest segment 
+  int xid,yid,zid = TMath::BinarySearch(fNZSegDip,fSegZDip,(Float_t)xyz[2]); // find zsegment
+  int ysegBeg = fBegSegYDip[zid];
+  //
+  for (yid=0;yid<fNSegYDip[zid];yid++) if (xyz[1]<fSegYDip[ysegBeg+yid]) break;
+  if ( --yid < 0 ) yid = 0;
+  yid +=  ysegBeg;
+  //
+  int xsegBeg = fBegSegXDip[yid];
+  for (xid=0;xid<fNSegXDip[yid];xid++) if (xyz[0]<fSegXDip[xsegBeg+xid]) break;
+  if ( --xid < 0) xid = 0;
+  xid +=  xsegBeg;
+  //
+  return fSegIDDip[xid];
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::GetTPCInt(const Float_t *xyz, Float_t *b) const
+void AliMagWrapCheb::GetTPCInt(const Double_t *xyz, Double_t *b) const
 {
   // compute TPC region field integral in cartesian coordinates.
   // If point is outside of the parameterized region get it at closeset valid point
-  static float rphiz[3];
+  static Double_t rphiz[3];
   //
   // TPCInt region
   // convert coordinates to cyl system
@@ -331,7 +385,7 @@ void AliMagWrapCheb::GetTPCInt(const Float_t *xyz, Float_t *b) const
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::FieldCylSol(const float *rphiz, float *b) const
+void AliMagWrapCheb::FieldCylSol(const Double_t *rphiz, Double_t *b) const
 {
   // compute Solenoid field in Cylindircal coordinates
   // note: if the point is outside the volume get the field in closest parameterized point
@@ -345,7 +399,7 @@ void AliMagWrapCheb::FieldCylSol(const float *rphiz, float *b) const
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::FieldCylSol(const Double_t *rphiz, Double_t *b) const
+Double_t AliMagWrapCheb::FieldCylSolBz(const Double_t *rphiz) const
 {
   // compute Solenoid field in Cylindircal coordinates
   // note: if the point is outside the volume get the field in closest parameterized point
@@ -354,12 +408,12 @@ void AliMagWrapCheb::FieldCylSol(const Double_t *rphiz, Double_t *b) const
   int SolRId = fSegZIdSol[SolZId];        // first R segment for this Z
   int SolRMax = SolRId + fNSegRSol[SolZId];
   while (rphiz[0]>fSegRSol[SolRId] && SolRId<SolRMax-1) ++SolRId;    // find R segment
-  GetParamSol( SolRId )->Eval(rphiz,b);
+  return GetParamSol( SolRId )->Eval(rphiz,2);
   //
 }
 
 //__________________________________________________________________________________________
-void AliMagWrapCheb::GetTPCIntCyl(const Float_t *rphiz, Float_t *b) const
+void AliMagWrapCheb::GetTPCIntCyl(const Double_t *rphiz, Double_t *b) const
 {
   // compute field integral in TPC region in Cylindircal coordinates
   // note: the check for the point being inside the parameterized region is done outside
@@ -373,47 +427,6 @@ void AliMagWrapCheb::GetTPCIntCyl(const Float_t *rphiz, Float_t *b) const
 }
 
 
-//__________________________________________________________________________________________
-void AliMagWrapCheb::Print(Option_t *) const
-{
-  // print info
-  printf("Alice magnetic field parameterized by Chebyshev polynomials\n");
-  printf("Segmentation for Solenoid (%+.2f<Z<%+.2f cm | R<%.2f cm)\n",fMinZSol,fMaxZSol,fMaxRSol);
-  //
-  if (fParamsSol) fParamsSol->Print();
-  /*
-  for (int iz=0;iz<fNSegZSol;iz++) {
-    AliCheb3D* param = GetParamSol( fSegZIdSol[iz] );
-    printf("*** Z Segment %2d (%+7.2f<Z<%+7.2f)\t***\n",iz,param->GetBoundMin(2),param->GetBoundMax(2));
-    for (int ir=0;ir<fNSegRSol[iz];ir++) {
-      param = GetParamSol( fSegZIdSol[iz]+ir );
-      printf("    R Segment %2d (%+7.2f<R<%+7.2f, Precision: %.1e) (ID=%2d)\n",ir, param->GetBoundMin(0),
-            param->GetBoundMax(0),param->GetPrecision(),fSegZIdSol[iz]+ir);
-    }
-  }
-  */
-  //
-  printf("Segmentation for TPC field integral (%+.2f<Z<%+.2f cm | R<%.2f cm)\n",fMinZTPCInt,fMaxZTPCInt,fMaxRTPCInt);
-  //
-  if (fParamsTPCInt) fParamsTPCInt->Print();
-  /*
-  for (int iz=0;iz<fNSegZTPCInt;iz++) {
-    AliCheb3D* param = GetParamTPCInt( fSegZIdTPCInt[iz] );
-    printf("*** Z Segment %2d (%+7.2f<Z<%+7.2f)\t***\n",iz,param->GetBoundMin(2),param->GetBoundMax(2));
-    for (int ir=0;ir<fNSegRTPCInt[iz];ir++) {
-      param = GetParamTPCInt( fSegZIdTPCInt[iz]+ir );
-      printf("    R Segment %2d (%+7.2f<R<%+7.2f, Precision: %.1e) (ID=%2d)\n",ir, param->GetBoundMin(0),
-            param->GetBoundMax(0),param->GetPrecision(),fSegZIdTPCInt[iz]+ir);
-    }
-  }
-  */
-  //
-  printf("Segmentation for Dipole (%+.2f<Z<%+.2f cm)\n",fMinZDip,fMaxZDip);
-  if (fParamsDip) fParamsDip->Print();
-  //
-  
-
-}
 #ifdef  _INC_CREATION_ALICHEB3D_
 //_______________________________________________
 void AliMagWrapCheb::LoadData(const char* inpfile)
@@ -778,7 +791,7 @@ void AliMagWrapCheb::BuildTableTPCInt()
   //
   if (fNParamsTPCInt<1) return;
   fNSegZTPCInt = 0;
-  fSegRTPCInt   = new Float_t[fNParamsTPCInt];
+  fSegRTPCInt     = new Float_t[fNParamsTPCInt];
   float *tmpbufF  = new float[fNParamsTPCInt+1];
   int   *tmpbufI  = new int[fNParamsTPCInt+1];
   int   *tmpbufI1 = new int[fNParamsTPCInt+1];
index 195dc50d71b8899ca24bddf2ff3ce9ce699c964e..942f5abdc49e113c249b00af9e1b29f20dd5c58e 100644 (file)
@@ -5,9 +5,9 @@
 //                                                                               //
 //  Wrapper for the set of mag.field parameterizations by Chebyshev polinomials  //
 //  To obtain the field in cartesian coordinates/components use                  //
-//    Field(float* xyz, float* bxyz);                                            //
+//    Field(double* xyz, double* bxyz);                                          //
 //  For cylindrical coordinates/components:                                      //
-//    FieldCyl(float* rphiz, float* brphiz)                                      //
+//    FieldCyl(double* rphiz, double* brphiz)                                    //
 //                                                                               //
 //  The solenoid part is parameterized in the volume  R<500, -550<Z<550 cm       //
 //                                                                               //
@@ -25,9 +25,9 @@
 //                                                                               //
 //  To obtain the field integral in the TPC region from given point to nearest   //
 //  cathod plane (+- 250 cm) use:                                                //
-//  GetTPCInt(float* xyz, float* bxyz);  for Cartesian frame                     //
+//  GetTPCInt(double* xyz, double* bxyz);  for Cartesian frame                   //
 //  or                                                                           //
-//  GetTPCIntCyl(Float_t *rphiz, Float_t *b); for Cylindrical frame              //
+//  GetTPCIntCyl(Double_t *rphiz, Double_t *b); for Cylindrical frame            //
 //                                                                               //
 //                                                                               //
 //  The units are kiloGauss and cm.                                              //
@@ -58,7 +58,7 @@ class AliMagWrapCheb: public TNamed
   //
   Int_t      GetNParamsSol()                              const {return fNParamsSol;}
   Int_t      GetNSegZSol()                                const {return fNSegZSol;}
-  Float_t*     GetSegZSol() const {return fSegZSol;}
+  Float_t*   GetSegZSol() const {return fSegZSol;}
   //
   Int_t      GetNParamsTPCInt()                           const {return fNParamsTPCInt;}
   Int_t      GetNSegZTPCInt()                             const {return fNSegZTPCInt;}
@@ -84,36 +84,26 @@ class AliMagWrapCheb: public TNamed
   //
   virtual void Print(Option_t * = "")                     const;
   //
-  virtual void Field(const Float_t *xyz, Float_t *b)                const;
-  virtual void Field(const Double_t *xyz, Double_t *b)              const;
+  virtual void Field(const Double_t *xyz, Double_t *b)    const;
+  Double_t     GetBz(const Double_t *xyz)                 const;
   //
-  virtual void FieldCyl(const Float_t *rphiz, Float_t *b)     const;
-  virtual void FieldCyl(const Double_t *rphiz, Double_t *b)   const;
+  void FieldCyl(const Double_t *rphiz, Double_t  *b)      const;
+  void GetTPCInt(const Double_t *xyz, Double_t *b)        const;
+  void GetTPCIntCyl(const Double_t *rphiz, Double_t *b)   const;
   //
-  virtual void GetTPCInt(const Float_t *xyz, Float_t *b)        const;
-  virtual void GetTPCIntCyl(const Float_t *rphiz, Float_t *b)   const;
-  //
-  template <class T>
-    Int_t      FindDipSegment(const T *xyz)               const; 
-  //
-  template <class T>
-    static void CylToCartCylB(const T *rphiz, const T *brphiz,T *bxyz);
-  template <class T>
-    static void CylToCartCartB(const T *xyz,  const T *brphiz,T *bxyz);
-  template <class T>
-    static void CartToCylCartB(const T *xyz,  const T *bxyz,  T *brphiz);
-  template <class T>  
-    static void CartToCylCylB(const T *rphiz, const T *bxyz,  T *brphiz);
-  template <class T>
-    static void CartToCyl(const T *xyz,  T *rphiz);
-  template <class T>
-    static void CylToCart(const T *rphiz,T *xyz);
+  Int_t       FindDipSegment(const Double_t *xyz)         const; 
+  static void CylToCartCylB(const Double_t *rphiz, const Double_t *brphiz,Double_t *bxyz);
+  static void CylToCartCartB(const Double_t *xyz,  const Double_t *brphiz,Double_t *bxyz);
+  static void CartToCylCartB(const Double_t *xyz,  const Double_t *bxyz,  Double_t *brphiz);
+  static void CartToCylCylB(const Double_t *rphiz, const Double_t *bxyz,  Double_t *brphiz);
+  static void CartToCyl(const Double_t *xyz,  Double_t *rphiz);
+  static void CylToCart(const Double_t *rphiz,Double_t *xyz);
   //
 #ifdef  _INC_CREATION_ALICHEB3D_                          // see AliCheb3D.h for explanation
   void         LoadData(const char* inpfile);
   //
   AliMagWrapCheb(const char* inputFile);
-  void       SaveData(const char* outfile)              const;
+  void       SaveData(const char* outfile)                const;
   Int_t      SegmentDipDimension(Float_t** seg,const TObjArray* par,int npar, int dim, 
                                 Float_t xmn,Float_t xmx,Float_t ymn,Float_t ymx,Float_t zmn,Float_t zmx);
   //
@@ -129,8 +119,8 @@ class AliMagWrapCheb: public TNamed
 #endif
   //
  protected:
-    virtual void FieldCylSol(const Float_t *rphiz, Float_t *b)      const;
-    virtual void FieldCylSol(const Double_t *rphiz, Double_t *b)    const;
+  void     FieldCylSol(const Double_t *rphiz, Double_t *b)    const;
+  Double_t FieldCylSolBz(const Double_t *rphiz)               const;
   //
  protected:
   //
@@ -182,20 +172,11 @@ class AliMagWrapCheb: public TNamed
   TObjArray* fParamsDip;             // Parameterization pieces for Dipole field
   TObjArray* fParamsTPCInt;          // Parameterization pieces for Solenoid field integrals in TPC region
   //
-  ClassDef(AliMagWrapCheb,3)            // Wrapper class for the set of Chebishev parameterizations of Alice mag.field
+  ClassDef(AliMagWrapCheb,4)         // Wrapper class for the set of Chebishev parameterizations of Alice mag.field
   //
  };
 
 
-//__________________________________________________________________________________________
-inline void AliMagWrapCheb::FieldCyl(const Float_t *rphiz, Float_t *b) const
-{
-  // compute field in Cylindircal coordinates
-  //  if (rphiz[2]<GetMinZSol() || rphiz[2]>GetMaxZSol() || rphiz[0]>GetMaxRSol()) {for (int i=3;i--;) b[i]=0; return;}
-  FieldCylSol(rphiz,b);
-}
-
-
 //__________________________________________________________________________________________
 inline void AliMagWrapCheb::FieldCyl(const Double_t *rphiz, Double_t *b) const
 {
@@ -205,12 +186,11 @@ inline void AliMagWrapCheb::FieldCyl(const Double_t *rphiz, Double_t *b) const
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCartCylB(const T *rphiz, const T *brphiz,T *bxyz)
+inline void AliMagWrapCheb::CylToCartCylB(const Double_t *rphiz, const Double_t *brphiz,Double_t *bxyz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cyl.system
-  T btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
-  T psiPLUSphi = TMath::ATan2(brphiz[1],brphiz[0]) + rphiz[1];
+  Double_t btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
+  Double_t psiPLUSphi = TMath::ATan2(brphiz[1],brphiz[0]) + rphiz[1];
   bxyz[0] = btr*TMath::Cos(psiPLUSphi);
   bxyz[1] = btr*TMath::Sin(psiPLUSphi);
   bxyz[2] = brphiz[2];
@@ -218,12 +198,11 @@ inline void AliMagWrapCheb::CylToCartCylB(const T *rphiz, const T *brphiz,T *bxy
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCartCartB(const T *xyz, const T *brphiz, T *bxyz)
+inline void AliMagWrapCheb::CylToCartCartB(const Double_t* xyz, const Double_t *brphiz, Double_t *bxyz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cart.system
-  T btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
-  T phiPLUSpsi = TMath::ATan2(xyz[1],xyz[0]) +  TMath::ATan2(brphiz[1],brphiz[0]);
+  Double_t btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
+  Double_t phiPLUSpsi = TMath::ATan2(xyz[1],xyz[0]) +  TMath::ATan2(brphiz[1],brphiz[0]);
   bxyz[0] = btr*TMath::Cos(phiPLUSpsi);
   bxyz[1] = btr*TMath::Sin(phiPLUSpsi);
   bxyz[2] = brphiz[2];
@@ -231,12 +210,11 @@ inline void AliMagWrapCheb::CylToCartCartB(const T *xyz, const T *brphiz, T *bxy
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCylCartB(const T *xyz, const T *bxyz, T *brphiz)
+inline void AliMagWrapCheb::CartToCylCartB(const Double_t *xyz, const Double_t *bxyz, Double_t *brphiz)
 {
   // convert field in cylindrical coordinates to cartesian system, poin is in cart.system
-  T btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
-  T psiMINphi = TMath::ATan2(bxyz[1],bxyz[0]) - TMath::ATan2(xyz[1],xyz[0]);
+  Double_t btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
+  Double_t psiMINphi = TMath::ATan2(bxyz[1],bxyz[0]) - TMath::ATan2(xyz[1],xyz[0]);
   //
   brphiz[0] = btr*TMath::Cos(psiMINphi);
   brphiz[1] = btr*TMath::Sin(psiMINphi);
@@ -245,12 +223,11 @@ inline void AliMagWrapCheb::CartToCylCartB(const T *xyz, const T *bxyz, T *brphi
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCylCylB(const T *rphiz, const T *bxyz, T *brphiz)
+inline void AliMagWrapCheb::CartToCylCylB(const Double_t *rphiz, const Double_t *bxyz, Double_t *brphiz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cyl.system
-  T btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
-  T psiMINphi =  TMath::ATan2(bxyz[1],bxyz[0]) - rphiz[1];
+  Double_t btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
+  Double_t psiMINphi =  TMath::ATan2(bxyz[1],bxyz[0]) - rphiz[1];
   brphiz[0] = btr*TMath::Cos(psiMINphi);
   brphiz[1] = btr*TMath::Sin(psiMINphi);
   brphiz[2] = bxyz[2];
@@ -258,8 +235,7 @@ inline void AliMagWrapCheb::CartToCylCylB(const T *rphiz, const T *bxyz, T *brph
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCyl(const T *xyz,T *rphiz)
+inline void AliMagWrapCheb::CartToCyl(const Double_t *xyz, Double_t *rphiz)
 {
   rphiz[0] = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
   rphiz[1] = TMath::ATan2(xyz[1],xyz[0]);
@@ -267,32 +243,11 @@ inline void AliMagWrapCheb::CartToCyl(const T *xyz,T *rphiz)
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCart(const T *rphiz, T *xyz)
+inline void AliMagWrapCheb::CylToCart(const Double_t *rphiz, Double_t *xyz)
 {
   xyz[0] = rphiz[0]*TMath::Cos(rphiz[1]);
   xyz[1] = rphiz[0]*TMath::Sin(rphiz[1]);
   xyz[2] = rphiz[2];
 }
 
-//__________________________________________________________________________________________________
-template <class T>
-Int_t    AliMagWrapCheb::FindDipSegment(const T *xyz) const 
-{
-  // find the segment containing point xyz. If it is outside find the closest segment 
-  int xid,yid,zid = TMath::BinarySearch(fNZSegDip,fSegZDip,(Float_t)xyz[2]); // find zsegment
-  int ysegBeg = fBegSegYDip[zid];
-  //
-  for (yid=0;yid<fNSegYDip[zid];yid++) if (xyz[1]<fSegYDip[ysegBeg+yid]) break;
-  if ( --yid < 0 ) yid = 0;
-  yid +=  ysegBeg;
-  //
-  int xsegBeg = fBegSegXDip[yid];
-  for (xid=0;xid<fNSegXDip[yid];xid++) if (xyz[0]<fSegXDip[xsegBeg+xid]) break;
-  if ( --xid < 0) xid = 0;
-  xid +=  xsegBeg;
-  //
-  return fSegIDDip[xid];
-}
-
 #endif
index a34addb6b8d8ab8db85a1b4f7b5acefca424b044..d9946551ec53ca427f4f738ad6e2234a351a1335 100644 (file)
@@ -8,15 +8,13 @@ AliDigit.cxx  AliHit.cxx
 AliRun.cxx AliGenerator.cxx AliVertexGenerator.cxx 
 AliLego.cxx    AliModule.cxx   AliDigitNew.cxx 
 AliGeometry.cxx        AliRecPoint.cxx 
-AliHitMap.cxx  AliMagFC.cxx    AliMagFCM.cxx 
-AliMagFDM.cxx  AliLegoGenerator.cxx AliLegoGeneratorXYZ.cxx
+AliHitMap.cxx  AliLegoGenerator.cxx AliLegoGeneratorXYZ.cxx
 AliLegoGeneratorPhiZ.cxx AliLegoGeneratorEta.cxx AliLegoGeneratorEtaR.cxx 
 AliRndm.cxx 
 AliDebugVolume.cxx 
 AliConfig.cxx 
 AliRunDigitizer.cxx AliDigitizer.cxx AliStream.cxx 
-AliMergeCombi.cxx AliMagFMaps.cxx AliFieldMap.cxx 
-AliMagFCheb.cxx AliCheb3D.cxx AliCheb3DCalc.cxx
+AliMergeCombi.cxx 
 AliGausCorr.cxx 
 AliTrackMap.cxx AliTrackMapper.cxx 
 AliMemoryWatcher.cxx 
@@ -46,8 +44,6 @@ AliMathBase.cxx AliSignalProcesor.cxx
 AliTracker.cxx AliCluster.cxx AliCluster3D.cxx 
 AliHelix.cxx AliV0.cxx AliKink.cxx 
 AliSelectorRL.cxx 
-AliMagFMapsV1.cxx 
-AliMagWrapCheb.cxx 
 AliSurveyObj.cxx AliSurveyPoint.cxx 
 AliSurveyToAlignObjs.cxx 
 AliFstream.cxx 
@@ -69,6 +65,8 @@ AliGRPRecoParam.cxx
 AliRelAlignerKalman.cxx 
 AliESDTagCreator.cxx 
 AliGRPObject.cxx
+AliMillePede2.cxx AliMatrixSq.cxx AliMatrixSparse.cxx 
+AliSymMatrix.cxx AliMinResSolve.cxx
 )
 
 
index 0eff6689d5e15870c3f64b314b6448f28ce3214d..143c95b2ddd561d96ba1f5960639bbcb98da599d 100644 (file)
@@ -20,7 +20,7 @@ set(SRCS
   AliStack.cxx AliMCEventHandler.cxx AliInputEventHandler.cxx
   AliTrackReference.cxx AliSysInfo.cxx
   AliMCEvent.cxx AliMCParticle.cxx
-  AliMagF.cxx
+  AliMagF.cxx AliMagF.cxx AliMagWrapCheb.cxx AliCheb3D.cxx AliCheb3DCalc.cxx
   AliCodeTimer.cxx
   AliPDG.cxx
   AliQA.cxx
index 350c3c5ff29114cbfea6614a71797b81f4a3f0bb..b918aacde36e105dfc98f4b68b087fd4c4b44dcd 100644 (file)
 #pragma link C++ class AliMCParticle+;
 
 #pragma link C++ class  AliMagF+;
+#pragma link C++ class  AliMagWrapCheb+;
+#pragma link C++ class  AliCheb3DCalc+;
+#pragma link C++ class  AliCheb3D+;
+
 
 #pragma link C++ class AliCodeTimer+;
 #pragma link C++ class AliCodeTimer::AliPair+;
index 0fe17af0540ec7f9552d419ae8f856edb75689ef..e0bf1b9b1954947426c662a4703b967e89aa1600 100644 (file)
 #pragma link C++ class  AliDetector+;
 #pragma link C++ class  AliDigit+;
 #pragma link C++ class  AliHit+;
-#pragma link C++ class  AliMagFC+;
-#pragma link C++ class  AliMagFCM+;
-#pragma link C++ class  AliMagFMaps-;
-#pragma link C++ class  AliMagFMapsV1+;
-#pragma link C++ class  AliMagFDM+;
-#pragma link C++ class  AliMagFCheb+;
-#pragma link C++ class  AliCheb3DCalc+;
-#pragma link C++ class  AliCheb3D+;
-#pragma link C++ class  AliMagWrapCheb+;
 #pragma link C++ class  AliLego+;
 #pragma link C++ class  AliLegoGenerator+;
 #pragma link C++ class  AliLegoGeneratorXYZ+;
@@ -47,7 +38,6 @@
 #pragma link C++ class  AliRunDigitizer+;
 #pragma link C++ class  AliStream+;
 #pragma link C++ class  AliMergeCombi+;
-#pragma link C++ class  AliFieldMap-;
 #pragma link C++ class  AliGausCorr+;
 #pragma link C++ class  AliLoader+;
 #pragma link C++ class  AliDataLoader+;
 
 #pragma link C++ class AliGRPObject+;
 
+#pragma link C++ class AliMillePede2+;
+#pragma link C++ class AliMillePedeRecord+;
+#pragma link C++ class AliMinResSolve+;
+#pragma link C++ class AliMatrixSparse+;
+#pragma link C++ class AliMatrixSq+;
+#pragma link C++ class AliSymMatrix+;
+
+
 #endif
index 2a367dbe012a0479f4a40e2232edf8b7e3f20943..3fc0d8424b0e11efa9911ae039608164143be9ca 100644 (file)
@@ -7,15 +7,13 @@ AliDigit.cxx  AliHit.cxx      \
 AliRun.cxx AliGenerator.cxx AliVertexGenerator.cxx \
 AliLego.cxx    AliModule.cxx   AliDigitNew.cxx \
 AliGeometry.cxx        AliRecPoint.cxx \
-AliHitMap.cxx  AliMagFC.cxx    AliMagFCM.cxx \
-AliMagFDM.cxx  AliLegoGenerator.cxx AliLegoGeneratorXYZ.cxx\
+AliHitMap.cxx AliLegoGenerator.cxx AliLegoGeneratorXYZ.cxx\
 AliLegoGeneratorPhiZ.cxx AliLegoGeneratorEta.cxx AliLegoGeneratorEtaR.cxx \
 AliRndm.cxx \
 AliDebugVolume.cxx \
 AliConfig.cxx \
 AliRunDigitizer.cxx AliDigitizer.cxx AliStream.cxx \
-AliMergeCombi.cxx AliMagFMaps.cxx AliFieldMap.cxx \
-AliMagFCheb.cxx AliCheb3D.cxx AliCheb3DCalc.cxx\
+AliMergeCombi.cxx \
 AliGausCorr.cxx \
 AliTrackMap.cxx AliTrackMapper.cxx \
 AliMemoryWatcher.cxx \
@@ -45,8 +43,6 @@ AliMathBase.cxx AliSignalProcesor.cxx \
 AliTracker.cxx AliCluster.cxx AliCluster3D.cxx \
 AliHelix.cxx AliV0.cxx AliKink.cxx \
 AliSelectorRL.cxx \
-AliMagFMapsV1.cxx \
-AliMagWrapCheb.cxx \
 AliSurveyObj.cxx AliSurveyPoint.cxx \
 AliSurveyToAlignObjs.cxx \
 AliFstream.cxx \
@@ -67,8 +63,8 @@ AliTriggerRunScalers.cxx AliGRPPreprocessor.cxx \
 AliGRPRecoParam.cxx \
 AliRelAlignerKalman.cxx \
 AliESDTagCreator.cxx \
-AliGRPObject.cxx
-
+AliGRPObject.cxx \
+AliMillePede2.cxx AliMatrixSq.cxx AliMatrixSparse.cxx AliSymMatrix.cxx AliMinResSolve.cxx
 
 HDRS:= $(SRCS:.cxx=.h) 
 
index 8dbfd1488439508ea7689577977af7b0dd8fa4a4..43c56b261c92f6fe2519439fada3744b63b22322 100644 (file)
@@ -19,7 +19,7 @@ SRCS = AliVParticle.cxx \
        AliStack.cxx AliMCEventHandler.cxx AliInputEventHandler.cxx \
        AliTrackReference.cxx AliSysInfo.cxx \
        AliMCEvent.cxx AliMCParticle.cxx \
-       AliMagF.cxx \
+       AliMagF.cxx AliMagWrapCheb.cxx AliCheb3D.cxx AliCheb3DCalc.cxx \
        AliCodeTimer.cxx \
        AliPDG.cxx \
        AliQA.cxx \