Use const for mag field input

[u/mrichter/AliRoot.git] / STEER / AliMagFMaps.cxx

diff --git a/STEER/AliMagFMaps.cxx b/STEER/AliMagFMaps.cxx
index 1c78a65184dc4e136d9d02b5134d7d8dcb78d4ff..7db65add35693b42ac9607d3d04a120f07433351 100644 (file)
--- a/STEER/AliMagFMaps.cxx
+++ b/STEER/AliMagFMaps.cxx
@@ -22,7 +22,9 @@
 // Author: Andreas Morsch <andreas.morsch@cern.ch>
 //------------------------------------------------------------------------
 
+#include <TClass.h>
 #include <TFile.h>
+#include <TMath.h>
 #include <TSystem.h>
 
 #include "AliLog.h"
@@ -73,6 +75,7 @@ AliMagFMaps::AliMagFMaps(const char *name, const char *title, Int_t integ,
 AliMagFMaps::AliMagFMaps(const AliMagFMaps &magf):
   AliMagFC(magf),
   fSolenoid(0),
+  fSolenoidUser(0),
   fL3Option(0)
 {
   //
@@ -99,10 +102,12 @@ void AliMagFMaps::ReadField()
   //
   //  don't read twice
   //
-    if (!fgReadField) return;
-    fgReadField = 0;
+    if (!fReadField) return;
+    fReadField = kFALSE;
   //    
-  char* fname;
+   
+
+  char* fname = 0;
   TFile* file = 0;
   if (fMap == k2kG) {
       fSolenoid = 2.;
@@ -111,17 +116,19 @@ void AliMagFMaps::ReadField()
       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B02"));
       file->Close();
       delete file;
+      delete [] fname; // see documentation to expand filenam
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB02.root");
       file = new TFile(fname);
       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB02"));
       file->Close();
       delete file;;
-      
+      delete [] fname;
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB02.root");
       file = new TFile(fname);
       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB02"));
       file->Close();
       delete file;
+      delete [] fname;
   } else if (fMap == k4kG) {
       fSolenoid = 4.;
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/L3B04.root");
@@ -129,16 +136,19 @@ void AliMagFMaps::ReadField()
       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B04"));
       file->Close();
       delete file;
+      delete [] fname;
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB04.root");
       file = new TFile(fname);
       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB04"));
       file->Close();
+      delete [] fname;
       delete file;
       
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB04.root");
       file = new TFile(fname);
       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB04"));
       file->Close();
+      delete [] fname;
       delete file;
   } else if (fMap == k5kG) {
       fSolenoid = 5.;
@@ -146,18 +156,20 @@ void AliMagFMaps::ReadField()
       file = new TFile(fname);
       fFieldMap[0] = dynamic_cast<AliFieldMap*>(file->Get("L3B05"));
       file->Close();
-      delete file;
+      delete file;      
+      delete [] fname;
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/DipB05.root");
       file = new TFile(fname);
       fFieldMap[1] = dynamic_cast<AliFieldMap*>(file->Get("DipB05"));
       file->Close();
       delete file;
-      
+      delete [] fname;
       fname = gSystem->ExpandPathName("$(ALICE_ROOT)/data/maps/ExtB05.root");
       file = new TFile(fname);
       fFieldMap[2] = dynamic_cast<AliFieldMap*>(file->Get("ExtB05"));
       file->Close();
       delete file;
+      delete [] fname;
   }
 
   if (!fL3Option) fSolenoidUser = fSolenoid;
@@ -171,11 +183,11 @@ Float_t AliMagFMaps::SolenoidField() const
   // Returns max. L3 (solenoid) field strength 
   // according to field map setting 
   //
-  return fSolenoid;
+  return -Factor()*fSolenoid;
 }
 
 //_______________________________________________________________________
-void AliMagFMaps::Field(Float_t *x, Float_t *b) const
+void AliMagFMaps::Field(const float *x, float *b) const
 {
   //
   // Method to calculate the magnetic field
@@ -207,7 +219,7 @@ void AliMagFMaps::Field(Float_t *x, Float_t *b) const
       //
       //     Constant L3 field , if this option was selected
       //
-         b[2] = - fSolenoid;
+       b[2] = (- fSolenoid)*fFactor;
          return;
       } 
   } else if (fFieldMap[1]->Inside(xm[0], xm[1], xm[2])) {
@@ -217,20 +229,74 @@ void AliMagFMaps::Field(Float_t *x, Float_t *b) const
   }
   
   if(map){
-    map->Field(xm,b);
-    b[0] = - b[0];
-    b[2] = - b[2];
-
+      map->Field(xm,b);
+      b[0] = - b[0];
+      b[2] = - b[2];
+      
+      if(fFactor!=1) {
+         b[0]*=fFactor;
+         b[1]*=fFactor;
+         b[2]*=fFactor;
+      }
   } else {
       //This is the ZDC part
       ZDCField(x, b);
   }
+}
 
+//_______________________________________________________________________
+void AliMagFMaps::Field(const double *x, double *b) const
+{
+  //
+  // Method to calculate the magnetic field
+  //
+  // --- find the position in the grid ---
+  
+ //    if (!fFieldRead) ReadField();
+
+  //
+  // Field Maps have been calculated for the coordinate system in which 
+  // the Muon Spectrometer is placed at z > 0
+  // Transform coordinates corresponingly 
+  //
   
-  if(fFactor!=1) {
-      b[0]*=fFactor;
-      b[1]*=fFactor;
-      b[2]*=fFactor;
+  b[0]=b[1]=b[2]=0;
+  double xm[3];
+  xm[0] = - x[0];
+  xm[1] =   x[1];
+  xm[2] = - x[2];
+  
+  AliFieldMap* map = 0;
+  if (fFieldMap[0]->Inside(xm[0], xm[1], xm[2])) {
+      map = fFieldMap[0];
+      Float_t r = TMath::Sqrt(xm[0] * xm[0] + xm[1] * xm[1]);
+      
+      if (!fL3Option && TMath::Abs(xm[2]) < 370. && r < 550.) {
+      //
+      //     Constant L3 field , if this option was selected
+      //
+       b[2] = (- fSolenoid)*fFactor;
+         return;
+      } 
+  } else if (fFieldMap[1]->Inside(xm[0], xm[1], xm[2])) {
+    map = fFieldMap[1];
+  } else if (fFieldMap[2]->Inside(xm[0], xm[1], xm[2])) {
+    map = fFieldMap[2];
+  }
+  
+  if(map){
+      map->Field(xm,b);
+      b[0] = - b[0];
+      b[2] = - b[2];
+      
+      if(fFactor!=1) {
+         b[0]*=fFactor;
+         b[1]*=fFactor;
+         b[2]*=fFactor;
+      }
+  } else {
+      //This is the ZDC part
+      ZDCField(x, b);
   }
 }
 
@@ -248,9 +314,11 @@ void AliMagFMaps::Streamer(TBuffer &R__b)
 {
   // Stream an object of class AliMagFMaps.
   if (R__b.IsReading()) {
-    AliMagFMaps::Class()->ReadBuffer(R__b, this);
-    ReadField();
+   R__b.ReadClassBuffer(AliMagFMaps::Class(), this);
+   // 
+   fReadField = kTRUE;
+   ReadField();
   } else {
-    AliMagFMaps::Class()->WriteBuffer(R__b, this);
+    R__b.WriteClassBuffer(AliMagFMaps::Class(), this);
   }
 }