STEER/AliMagF.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16
  17 #include <TClass.h>
  18 #include <TFile.h>
  19 #include <TSystem.h>
  20 #include <TPRegexp.h>
  21
  22 #include "AliMagF.h"
  23 #include "AliMagWrapCheb.h"
  24 #include "AliLog.h"
  25
  26 ClassImp(AliMagF)
  27
  28 const Double_t AliMagF::fgkSol2DipZ    =  -700.;
  29 const UShort_t AliMagF::fgkPolarityConvention = AliMagF::kConvLHC;
  30 /*
  31  Explanation for polarity conventions: these are the mapping between the
  32  current signs and main field components in L3 (Bz) and Dipole (Bx) (in Alice frame)
  33  1) kConvMap2005: used for the field mapping in 2005
  34  positive L3  current -> negative Bz
  35  positive Dip current -> positive Bx
  36  2) kConvMapDCS2008: defined by the microswitches/cabling of power converters as of 2008 - 1st half 2009
  37  positive L3  current -> positive Bz
  38  positive Dip current -> positive Bx
  39  3) kConvLHC : defined by LHC
  40  positive L3  current -> positive Bz
  41  positive Dip current -> negative Bx
  42
  43  Note: only "negative Bz(L3) with postive Bx(Dipole)" and its inverse was mapped in 2005. Hence
  44  the GRP Manager will reject the runs with the current combinations (in the convention defined by the
  45  static Int_t AliMagF::GetPolarityConvention()) which do not lead to such field polarities.
  46
  47  -----------------------------------------------
  48
  49  Explanation on integrals in the TPC region
  50  GetTPCInt(xyz,b) and GetTPCRatInt(xyz,b) give integrals from point (x,y,z) to point (x,y,0)
  51  (irrespectively of the z sign) of the following:
  52  TPCInt:    b contains int{bx}, int{by}, int{bz}
  53  TPCRatInt: b contains int{bx/bz}, int{by/bz}, int{(bx/bz)^2+(by/bz)^2}
  54
  55  The same applies to integral in cylindrical coordinates:
  56  GetTPCIntCyl(rphiz,b)
  57  GetTPCIntRatCyl(rphiz,b)
  58  They accept the R,Phi,Z coordinate (-pi<phi<pi) and return the field
  59  integrals in cyl. coordinates.
  60
  61  Thus, to compute the integral from arbitrary xy_z1 to xy_z2, one should take
  62  b = b1-b2 with b1 and b2 coming from GetTPCInt(xy_z1,b1) and GetTPCInt(xy_z2,b2)
  63
  64  Note: the integrals are defined for the range -300<Z<300 and 0<R<300
  65 */
  66 //_______________________________________________________________________
  67 AliMagF::AliMagF():
  68   TVirtualMagField(),
  69   fMeasuredMap(0),
  70   fMapType(k5kG),
  71   fSolenoid(0),
  72   fBeamType(kNoBeamField),
  73   fBeamEnergy(0),
  74   //
  75   fInteg(0),
  76   fPrecInteg(0),
  77   fFactorSol(1.),
  78   fFactorDip(1.),
  79   fMax(15),
  80   fDipoleOFF(kFALSE),
  81   //
  82   fQuadGradient(0),
  83   fDipoleField(0),
  84   fCCorrField(0),
  85   fACorr1Field(0),
  86   fACorr2Field(0),
  87   fParNames("","")
  88 {
  89   // Default constructor
  90   //
  91 }
  92
  93 //_______________________________________________________________________
  94 AliMagF::AliMagF(const char *name, const char* title, Double_t factorSol, Double_t factorDip,
  95                  BMap_t maptype, BeamType_t bt, Double_t be,Int_t integ, Double_t fmax, const char* path):
  96   TVirtualMagField(name),
  97   fMeasuredMap(0),
  98   fMapType(maptype),
  99   fSolenoid(0),
 100   fBeamType(bt),
 101   fBeamEnergy(be),
 102   //
 103   fInteg(integ),
 104   fPrecInteg(1),
 105   fFactorSol(1.),
 106   fFactorDip(1.),
 107   fMax(fmax),
 108   fDipoleOFF(factorDip==0.),
 109   //
 110   fQuadGradient(0),
 111   fDipoleField(0),
 112   fCCorrField(0),
 113   fACorr1Field(0),
 114   fACorr2Field(0),
 115   fParNames("","")
 116 {
 117   // Initialize the field with Geant integration option "integ" and max field "fmax,
 118   // Impose scaling of parameterized L3 field by factorSol and of dipole by factorDip.
 119   // The "be" is the energy of the beam in GeV/nucleon
 120   //
 121   SetTitle(title);
 122   if(integ<0 || integ > 2) {
 123     AliWarning(Form("Invalid magnetic field flag: %5d; Helix tracking chosen instead",integ));
 124     fInteg = 2;
 125   }
 126   if (fInteg == 0) fPrecInteg = 0;
 127   //
 128   if (fBeamEnergy<=0 && fBeamType!=kNoBeamField) {
 129     if      (fBeamType == kBeamTypepp) fBeamEnergy = 7000.; // max proton energy
 130     else if (fBeamType == kBeamTypeAA) fBeamEnergy = 2760;  // max PbPb energy
 131     AliInfo("Maximim possible beam energy for requested beam is assumed");
 132   }
 133   const char* parname = 0;
 134   //
 135   if      (fMapType == k2kG) parname = fDipoleOFF ? "Sol12_Dip0_Hole":"Sol12_Dip6_Hole";
 136   else if (fMapType == k5kG) parname = fDipoleOFF ? "Sol30_Dip0_Hole":"Sol30_Dip6_Hole";
 137   else if (fMapType == k5kGUniform) parname = "Sol30_Dip6_Uniform";
 138   else AliFatal(Form("Unknown field identifier %d is requested\n",fMapType));
 139   //
 140   SetDataFileName(path);
 141   SetParamName(parname);
 142   //
 143   LoadParameterization();
 144   InitMachineField(fBeamType,fBeamEnergy);
 145   double xyz[3]={0.,0.,0.};
 146   fSolenoid = GetBz(xyz);
 147   SetFactorSol(factorSol);
 148   SetFactorDip(factorDip);
 149   Print("a");
 150 }
 151
 152 //_______________________________________________________________________
 153 AliMagF::AliMagF(const AliMagF &src):
 154   TVirtualMagField(src),
 155   fMeasuredMap(0),
 156   fMapType(src.fMapType),
 157   fSolenoid(src.fSolenoid),
 158   fBeamType(src.fBeamType),
 159   fBeamEnergy(src.fBeamEnergy),
 160   fInteg(src.fInteg),
 161   fPrecInteg(src.fPrecInteg),
 162   fFactorSol(src.fFactorSol),
 163   fFactorDip(src.fFactorDip),
 164   fMax(src.fMax),
 165   fDipoleOFF(src.fDipoleOFF),
 166   fQuadGradient(src.fQuadGradient),
 167   fDipoleField(src.fDipoleField),
 168   fCCorrField(src.fCCorrField),
 169   fACorr1Field(src.fACorr1Field),
 170   fACorr2Field(src.fACorr2Field),
 171   fParNames(src.fParNames)
 172 {
 173   if (src.fMeasuredMap) fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
 174 }
 175
 176 //_______________________________________________________________________
 177 AliMagF::~AliMagF()
 178 {
 179   delete fMeasuredMap;
 180 }
 181
 182 //_______________________________________________________________________
 183 Bool_t AliMagF::LoadParameterization()
 184 {
 185   if (fMeasuredMap) {
 186     AliFatal(Form("Field data %s are already loaded from %s\n",GetParamName(),GetDataFileName()));
 187   }
 188   //
 189   char* fname = gSystem->ExpandPathName(GetDataFileName());
 190   TFile* file = TFile::Open(fname);
 191   if (!file) {
 192     AliFatal(Form("Failed to open magnetic field data file %s\n",fname));
 193   }
 194   //
 195   fMeasuredMap = dynamic_cast<AliMagWrapCheb*>(file->Get(GetParamName()));
 196   if (!fMeasuredMap) {
 197     AliFatal(Form("Did not find field %s in %s\n",GetParamName(),fname));
 198   }
 199   file->Close();
 200   delete file;
 201   return kTRUE;
 202 }
 203
 204
 205 //_______________________________________________________________________
 206 void AliMagF::Field(const Double_t *xyz, Double_t *b)
 207 {
 208   // Method to calculate the field at point  xyz
 209   //
 210   //  b[0]=b[1]=b[2]=0.0;
 211   if (fMeasuredMap && xyz[2]>fMeasuredMap->GetMinZ() && xyz[2]<fMeasuredMap->GetMaxZ()) {
 212     fMeasuredMap->Field(xyz,b);
 213     if (xyz[2]>fgkSol2DipZ || fDipoleOFF) for (int i=3;i--;) b[i] *= fFactorSol;
 214     else                                  for (int i=3;i--;) b[i] *= fFactorDip;
 215   }
 216   else MachineField(xyz, b);
 217   //
 218 }
 219
 220 //_______________________________________________________________________
 221 Double_t AliMagF::GetBz(const Double_t *xyz) const
 222 {
 223   // Method to calculate the field at point  xyz
 224   //
 225   if (fMeasuredMap && xyz[2]>fMeasuredMap->GetMinZ() && xyz[2]<fMeasuredMap->GetMaxZ()) {
 226     double bz = fMeasuredMap->GetBz(xyz);
 227     return (xyz[2]>fgkSol2DipZ || fDipoleOFF) ? bz*fFactorSol : bz*fFactorDip;
 228   }
 229   else return 0.;
 230 }
 231
 232 //_______________________________________________________________________
 233 AliMagF& AliMagF::operator=(const AliMagF& src)
 234 {
 235   if (this != &src && src.fMeasuredMap) {
 236     if (fMeasuredMap) delete fMeasuredMap;
 237     fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
 238     SetName(src.GetName());
 239     fSolenoid    = src.fSolenoid;
 240     fBeamType    = src.fBeamType;
 241     fBeamEnergy  = src.fBeamEnergy;
 242     fInteg       = src.fInteg;
 243     fPrecInteg   = src.fPrecInteg;
 244     fFactorSol   = src.fFactorSol;
 245     fFactorDip   = src.fFactorDip;
 246     fMax         = src.fMax;
 247     fDipoleOFF   = src.fDipoleOFF;
 248     fParNames    = src.fParNames;
 249   }
 250   return *this;
 251 }
 252
 253 //_______________________________________________________________________
 254 void AliMagF::InitMachineField(BeamType_t btype, Double_t benergy)
 255 {
 256   if (btype==kNoBeamField) {
 257     fQuadGradient = fDipoleField = fCCorrField = fACorr1Field = fACorr2Field = 0.;
 258     return;
 259   }
 260   //
 261   double rigScale = benergy/7000.;   // scale according to ratio of E/Enominal
 262   // for ions assume PbPb (with energy provided per nucleon) and account for A/Z
 263   if (btype == kBeamTypeAA) rigScale *= 208./82.;
 264   //
 265   fQuadGradient = 22.0002*rigScale;
 266   fDipoleField  = 37.8781*rigScale;
 267   //
 268   // SIDE C
 269   fCCorrField   = -9.6980;
 270   // SIDE A
 271   fACorr1Field  = -13.2247;
 272   fACorr2Field  =  11.7905;
 273   //
 274 }
 275
 276 //_______________________________________________________________________
 277 void AliMagF::MachineField(const Double_t *x, Double_t *b) const
 278 {
 279   // ---- This is the ZDC part
 280   // Compansators for Alice Muon Arm Dipole
 281   const Double_t kBComp1CZ = 1075., kBComp1hDZ = 260./2., kBComp1SqR = 4.0*4.0;
 282   const Double_t kBComp2CZ = 2049., kBComp2hDZ = 153./2., kBComp2SqR = 4.5*4.5;
 283   //
 284   const Double_t kTripQ1CZ = 2615., kTripQ1hDZ = 637./2., kTripQ1SqR = 3.5*3.5;
 285   const Double_t kTripQ2CZ = 3480., kTripQ2hDZ = 550./2., kTripQ2SqR = 3.5*3.5;
 286   const Double_t kTripQ3CZ = 4130., kTripQ3hDZ = 550./2., kTripQ3SqR = 3.5*3.5;
 287   const Double_t kTripQ4CZ = 5015., kTripQ4hDZ = 637./2., kTripQ4SqR = 3.5*3.5;
 288   //
 289   const Double_t kDip1CZ = 6310.8,  kDip1hDZ = 945./2., kDip1SqRC = 4.5*4.5, kDip1SqRA = 3.375*3.375;
 290   const Double_t kDip2CZ = 12640.3, kDip2hDZ = 945./2., kDip2SqRC = 4.5*4.5, kDip2SqRA = 3.75*3.75;
 291   const Double_t kDip2DXC = 9.7, kDip2DXA = 9.4;
 292   //
 293   double rad2 = x[0] * x[0] + x[1] * x[1];
 294   //
 295   b[0] = b[1] = b[2] = 0;
 296   //
 297   // SIDE C **************************************************
 298   if(x[2]<0.){
 299     if(TMath::Abs(x[2]+kBComp2CZ)<kBComp2hDZ && rad2 < kBComp2SqR){
 300       b[0] = fCCorrField*fFactorDip;
 301     }
 302     else if(TMath::Abs(x[2]+kTripQ1CZ)<kTripQ1hDZ && rad2 < kTripQ1SqR){
 303       b[0] = fQuadGradient*x[1];
 304       b[1] = fQuadGradient*x[0];
 305     }
 306     else if(TMath::Abs(x[2]+kTripQ2CZ)<kTripQ2hDZ && rad2 < kTripQ2SqR){
 307       b[0] = -fQuadGradient*x[1];
 308       b[1] = -fQuadGradient*x[0];
 309     }
 310     else if(TMath::Abs(x[2]+kTripQ3CZ)<kTripQ3hDZ && rad2 < kTripQ3SqR){
 311       b[0] = -fQuadGradient*x[1];
 312       b[1] = -fQuadGradient*x[0];
 313     }
 314     else if(TMath::Abs(x[2]+kTripQ4CZ)<kTripQ4hDZ && rad2 < kTripQ4SqR){
 315       b[0] = fQuadGradient*x[1];
 316       b[1] = fQuadGradient*x[0];
 317     }
 318     else if(TMath::Abs(x[2]+kDip1CZ)<kDip1hDZ && rad2 < kDip1SqRC){
 319       b[1] = fDipoleField;
 320     }
 321     else if(TMath::Abs(x[2]+kDip2CZ)<kDip2hDZ && rad2 < kDip2SqRC) {
 322       double dxabs = TMath::Abs(x[0])-kDip2DXC;
 323       if ( (dxabs*dxabs + x[1]*x[1])<kDip2SqRC) {
 324         b[1] = -fDipoleField;
 325       }
 326     }
 327   }
 328   //
 329   // SIDE A **************************************************
 330   else{
 331     if(TMath::Abs(x[2]-kBComp1CZ)<kBComp1hDZ && rad2 < kBComp1SqR) {
 332       // Compensator magnet at z = 1075 m
 333       b[0] = fACorr1Field*fFactorDip;
 334     }
 335     //
 336     if(TMath::Abs(x[2]-kBComp2CZ)<kBComp2hDZ && rad2 < kBComp2SqR){
 337       b[0] = fACorr2Field*fFactorDip;
 338     }
 339     else if(TMath::Abs(x[2]-kTripQ1CZ)<kTripQ1hDZ && rad2 < kTripQ1SqR){
 340       b[0] = -fQuadGradient*x[1];
 341       b[1] = -fQuadGradient*x[0];
 342     }
 343     else if(TMath::Abs(x[2]-kTripQ2CZ)<kTripQ2hDZ && rad2 < kTripQ2SqR){
 344       b[0] =  fQuadGradient*x[1];
 345       b[1] =  fQuadGradient*x[0];
 346     }
 347     else if(TMath::Abs(x[2]-kTripQ3CZ)<kTripQ3hDZ && rad2 < kTripQ3SqR){
 348       b[0] =  fQuadGradient*x[1];
 349       b[1] =  fQuadGradient*x[0];
 350     }
 351     else if(TMath::Abs(x[2]-kTripQ4CZ)<kTripQ4hDZ && rad2 < kTripQ4SqR){
 352       b[0] = -fQuadGradient*x[1];
 353       b[1] = -fQuadGradient*x[0];
 354     }
 355     else if(TMath::Abs(x[2]-kDip1CZ)<kDip1hDZ && rad2 < kDip1SqRA){
 356       b[1] = -fDipoleField;
 357     }
 358     else if(TMath::Abs(x[2]-kDip2CZ)<kDip2hDZ && rad2 < kDip2SqRA) {
 359       double dxabs = TMath::Abs(x[0])-kDip2DXA;
 360       if ( (dxabs*dxabs + x[1]*x[1])<kDip2SqRA) {
 361         b[1] = fDipoleField;
 362       }
 363     }
 364   }
 365   //
 366 }
 367
 368 //_______________________________________________________________________
 369 void AliMagF::GetTPCInt(const Double_t *xyz, Double_t *b) const
 370 {
 371   // Method to calculate the integral_0^z of br,bt,bz
 372   b[0]=b[1]=b[2]=0.0;
 373   if (fMeasuredMap) {
 374     fMeasuredMap->GetTPCInt(xyz,b);
 375     for (int i=3;i--;) b[i] *= fFactorSol;
 376   }
 377 }
 378
 379 //_______________________________________________________________________
 380 void AliMagF::GetTPCRatInt(const Double_t *xyz, Double_t *b) const
 381 {
 382   // Method to calculate the integral_0^z of bx/bz,by/bz and (bx/bz)^2+(by/bz)^2
 383   b[0]=b[1]=b[2]=0.0;
 384   if (fMeasuredMap) {
 385     fMeasuredMap->GetTPCRatInt(xyz,b);
 386     b[2] /= 100;
 387   }
 388 }
 389
 390 //_______________________________________________________________________
 391 void AliMagF::GetTPCIntCyl(const Double_t *rphiz, Double_t *b) const
 392 {
 393   // Method to calculate the integral_0^z of br,bt,bz
 394   // in cylindrical coordiates ( -pi<phi<pi convention )
 395   b[0]=b[1]=b[2]=0.0;
 396   if (fMeasuredMap) {
 397     fMeasuredMap->GetTPCIntCyl(rphiz,b);
 398     for (int i=3;i--;) b[i] *= fFactorSol;
 399   }
 400 }
 401
 402 //_______________________________________________________________________
 403 void AliMagF::GetTPCRatIntCyl(const Double_t *rphiz, Double_t *b) const
 404 {
 405   // Method to calculate the integral_0^z of bx/bz,by/bz and (bx/bz)^2+(by/bz)^2
 406   // in cylindrical coordiates ( -pi<phi<pi convention )
 407   b[0]=b[1]=b[2]=0.0;
 408   if (fMeasuredMap) {
 409     fMeasuredMap->GetTPCRatIntCyl(rphiz,b);
 410     b[2] /= 100;
 411   }
 412 }
 413
 414 //_______________________________________________________________________
 415 void AliMagF::SetFactorSol(Float_t fc)
 416 {
 417   // set the sign/scale of the current in the L3 according to fgkPolarityConvention
 418   switch (fgkPolarityConvention) {
 419   case kConvDCS2008: fFactorSol = -fc; break;
 420   case kConvLHC    : fFactorSol = -fc; break;
 421   default          : fFactorSol =  fc; break;  // case kConvMap2005: fFactorSol =  fc; break;
 422   }
 423 }
 424
 425 //_______________________________________________________________________
 426 void AliMagF::SetFactorDip(Float_t fc)
 427 {
 428   // set the sign*scale of the current in the Dipole according to fgkPolarityConvention
 429   switch (fgkPolarityConvention) {
 430   case kConvDCS2008: fFactorDip =  fc; break;
 431   case kConvLHC    : fFactorDip = -fc; break;
 432   default          : fFactorDip =  fc; break;  // case kConvMap2005: fFactorDip =  fc; break;
 433   }
 434 }
 435
 436 //_______________________________________________________________________
 437 Double_t AliMagF::GetFactorSol() const
 438 {
 439   // return the sign*scale of the current in the Dipole according to fgkPolarityConventionthe
 440   switch (fgkPolarityConvention) {
 441   case kConvDCS2008: return -fFactorSol;
 442   case kConvLHC    : return -fFactorSol;
 443   default          : return  fFactorSol;       //  case kConvMap2005: return  fFactorSol;
 444   }
 445 }
 446
 447 //_______________________________________________________________________
 448 Double_t AliMagF::GetFactorDip() const
 449 {
 450   // return the sign*scale of the current in the Dipole according to fgkPolarityConventionthe
 451   switch (fgkPolarityConvention) {
 452   case kConvDCS2008: return  fFactorDip;
 453   case kConvLHC    : return -fFactorDip;
 454   default          : return  fFactorDip;       //  case kConvMap2005: return  fFactorDip;
 455   }
 456 }
 457
 458 //_____________________________________________________________________________
 459 AliMagF* AliMagF::CreateFieldMap(Float_t l3Cur, Float_t diCur, Int_t convention, Bool_t uniform,
 460                                  Float_t beamenergy, const Char_t *beamtype, const Char_t *path)
 461 {
 462   //------------------------------------------------
 463   // The magnetic field map, defined externally...
 464   // L3 current 30000 A  -> 0.5 T
 465   // L3 current 12000 A  -> 0.2 T
 466   // dipole current 6000 A
 467   // The polarities must match the convention (LHC or DCS2008)
 468   // unless the special uniform map was used for MC
 469   //------------------------------------------------
 470   const Float_t l3NominalCurrent1=30000.; // (A)
 471   const Float_t l3NominalCurrent2=12000.; // (A)
 472   const Float_t diNominalCurrent =6000. ; // (A)
 473
 474   const Float_t tolerance=0.03; // relative current tolerance
 475   const Float_t zero=77.;       // "zero" current (A)
 476   //
 477   BMap_t map = k5kG;
 478   double sclL3,sclDip;
 479   //
 480   Float_t l3Pol = l3Cur > 0 ? 1:-1;
 481   Float_t diPol = diCur > 0 ? 1:-1;
 482
 483   l3Cur = TMath::Abs(l3Cur);
 484   diCur = TMath::Abs(diCur);
 485   //
 486   if (TMath::Abs((sclDip=diCur/diNominalCurrent)-1.) > tolerance && !uniform) {
 487     if (diCur <= zero) sclDip = 0.; // some small current.. -> Dipole OFF
 488     else {
 489       AliFatalGeneral("AliMagF",Form("Wrong dipole current (%f A)!",diCur));
 490     }
 491   }
 492   //
 493   if (uniform) {
 494     // special treatment of special MC with uniform mag field (normalized to 0.5 T)
 495     // no check for scaling/polarities are done
 496     map   = k5kGUniform;
 497     sclL3 = l3Cur/l3NominalCurrent1;
 498   }
 499   else {
 500     if      (TMath::Abs((sclL3=l3Cur/l3NominalCurrent1)-1.) < tolerance) map  = k5kG;
 501     else if (TMath::Abs((sclL3=l3Cur/l3NominalCurrent2)-1.) < tolerance) map  = k2kG;
 502     else if (l3Cur <= zero && diCur<=zero)   { sclL3=0; sclDip=0; map  = k5kGUniform;}
 503     else {
 504       AliFatalGeneral("AliMagF",Form("Wrong L3 current (%f A)!",l3Cur));
 505     }
 506   }
 507   //
 508   if (sclDip!=0 && map!=k5kGUniform) {
 509     if ( (l3Cur<=zero) || ((convention==kConvLHC && l3Pol!=diPol) || (convention==kConvDCS2008 && l3Pol==diPol)) ) {
 510       AliFatalGeneral("AliMagF",Form("Wrong combination for L3/Dipole polarities (%c/%c) for convention %d",
 511                                      l3Pol>0?'+':'-',diPol>0?'+':'-',GetPolarityConvention()));
 512     }
 513   }
 514   //
 515   if (l3Pol<0) sclL3  = -sclL3;
 516   if (diPol<0) sclDip = -sclDip;
 517   //
 518   BeamType_t btype = kNoBeamField;
 519   TString btypestr = beamtype;
 520   btypestr.ToLower();
 521   TPRegexp protonBeam("(proton|p)\\s*-?\\s*\\1");
 522   TPRegexp ionBeam("(lead|pb|ion|a)\\s*-?\\s*\\1");
 523   if (btypestr.Contains(ionBeam)) btype = kBeamTypeAA;
 524   else if (btypestr.Contains(protonBeam)) btype = kBeamTypepp;
 525   else AliInfoGeneral("AliMagF",Form("Assume no LHC magnet field for the beam type %s, ",beamtype));
 526   char ttl[80];
 527   snprintf(ttl,79,"L3: %+5d Dip: %+4d kA; %s | Polarities in %s convention",(int)TMath::Sign(l3Cur,float(sclL3)),
 528           (int)TMath::Sign(diCur,float(sclDip)),uniform ? " Constant":"",
 529           convention==kConvLHC ? "LHC":"DCS2008");
 530   // LHC and DCS08 conventions have opposite dipole polarities
 531   if ( GetPolarityConvention() != convention) sclDip = -sclDip;
 532   //
 533   return new AliMagF("MagneticFieldMap", ttl,sclL3,sclDip,map,btype,beamenergy,2,10.,path);
 534   //
 535 }
 536
 537 //_____________________________________________________________________________
 538 const char*  AliMagF::GetBeamTypeText() const
 539 {
 540   const char *beamNA  = "No Beam";
 541   const char *beamPP  = "p-p";
 542   const char *beamPbPb= "Pb-Pb";
 543   switch ( fBeamType ) {
 544   case kBeamTypepp : return beamPP;
 545   case kBeamTypeAA : return beamPbPb;
 546   case kNoBeamField:
 547   default:           return beamNA;
 548   }
 549 }
 550
 551 //_____________________________________________________________________________
 552 void AliMagF::Print(Option_t *opt) const
 553 {
 554   // print short or long info
 555   TString opts = opt; opts.ToLower();
 556   AliInfo(Form("%s:%s",GetName(),GetTitle()));
 557   AliInfo(Form("Solenoid (%+.2f*)%.0f kG, Dipole %s (%+.2f) %s",
 558                GetFactorSol(),(fMapType==k5kG||fMapType==k5kGUniform)?5.:2.,
 559                fDipoleOFF ? "OFF":"ON",GetFactorDip(),fMapType==k5kGUniform?" |Constant Field!":""));
 560   if (opts.Contains("a")) {
 561     AliInfo(Form("Machine B fields for %s beam (%.0f GeV): QGrad: %.4f Dipole: %.4f",
 562                  fBeamType==kBeamTypeAA ? "A-A":(fBeamType==kBeamTypepp ? "p-p":"OFF"),
 563                  fBeamEnergy,fQuadGradient,fDipoleField));
 564     AliInfo(Form("Uses %s of %s",GetParamName(),GetDataFileName()));
 565   }
 566 }