STEER/STEERBase/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     else if (fBeamType == kBeamTypepA || fBeamType == kBeamTypeAp) fBeamEnergy = 2760;  // same rigitiy max PbPb energy
 132     AliInfo("Maximim possible beam energy for requested beam is assumed");
 133   }
 134   const char* parname = 0;
 135   //
 136   if      (fMapType == k2kG) parname = fDipoleOFF ? "Sol12_Dip0_Hole":"Sol12_Dip6_Hole";
 137   else if (fMapType == k5kG) parname = fDipoleOFF ? "Sol30_Dip0_Hole":"Sol30_Dip6_Hole";
 138   else if (fMapType == k5kGUniform) parname = "Sol30_Dip6_Uniform";
 139   else AliFatal(Form("Unknown field identifier %d is requested\n",fMapType));
 140   //
 141   SetDataFileName(path);
 142   SetParamName(parname);
 143   //
 144   LoadParameterization();
 145   InitMachineField(fBeamType,fBeamEnergy);
 146   double xyz[3]={0.,0.,0.};
 147   fSolenoid = GetBz(xyz);
 148   SetFactorSol(factorSol);
 149   SetFactorDip(factorDip);
 150   Print("a");
 151 }
 152
 153 //_______________________________________________________________________
 154 AliMagF::AliMagF(const AliMagF &src):
 155   TVirtualMagField(src),
 156   fMeasuredMap(0),
 157   fMapType(src.fMapType),
 158   fSolenoid(src.fSolenoid),
 159   fBeamType(src.fBeamType),
 160   fBeamEnergy(src.fBeamEnergy),
 161   fInteg(src.fInteg),
 162   fPrecInteg(src.fPrecInteg),
 163   fFactorSol(src.fFactorSol),
 164   fFactorDip(src.fFactorDip),
 165   fMax(src.fMax),
 166   fDipoleOFF(src.fDipoleOFF),
 167   fQuadGradient(src.fQuadGradient),
 168   fDipoleField(src.fDipoleField),
 169   fCCorrField(src.fCCorrField),
 170   fACorr1Field(src.fACorr1Field),
 171   fACorr2Field(src.fACorr2Field),
 172   fParNames(src.fParNames)
 173 {
 174   if (src.fMeasuredMap) fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
 175 }
 176
 177 //_______________________________________________________________________
 178 AliMagF::~AliMagF()
 179 {
 180   delete fMeasuredMap;
 181 }
 182
 183 //_______________________________________________________________________
 184 Bool_t AliMagF::LoadParameterization()
 185 {
 186   if (fMeasuredMap) {
 187     AliFatal(Form("Field data %s are already loaded from %s\n",GetParamName(),GetDataFileName()));
 188   }
 189   //
 190   char* fname = gSystem->ExpandPathName(GetDataFileName());
 191   TFile* file = TFile::Open(fname);
 192   if (!file) {
 193     AliFatal(Form("Failed to open magnetic field data file %s\n",fname));
 194   }
 195   //
 196   fMeasuredMap = dynamic_cast<AliMagWrapCheb*>(file->Get(GetParamName()));
 197   if (!fMeasuredMap) {
 198     AliFatal(Form("Did not find field %s in %s\n",GetParamName(),fname));
 199   }
 200   file->Close();
 201   delete file;
 202   return kTRUE;
 203 }
 204
 205
 206 //_______________________________________________________________________
 207 void AliMagF::Field(const Double_t *xyz, Double_t *b)
 208 {
 209   // Method to calculate the field at point  xyz
 210   //
 211   //  b[0]=b[1]=b[2]=0.0;
 212   if (fMeasuredMap && xyz[2]>fMeasuredMap->GetMinZ() && xyz[2]<fMeasuredMap->GetMaxZ()) {
 213     fMeasuredMap->Field(xyz,b);
 214     if (xyz[2]>fgkSol2DipZ || fDipoleOFF) for (int i=3;i--;) b[i] *= fFactorSol;
 215     else                                  for (int i=3;i--;) b[i] *= fFactorDip;
 216   }
 217   else MachineField(xyz, b);
 218   //
 219 }
 220
 221 //_______________________________________________________________________
 222 Double_t AliMagF::GetBz(const Double_t *xyz) const
 223 {
 224   // Method to calculate the field at point  xyz
 225   //
 226   if (fMeasuredMap && xyz[2]>fMeasuredMap->GetMinZ() && xyz[2]<fMeasuredMap->GetMaxZ()) {
 227     double bz = fMeasuredMap->GetBz(xyz);
 228     return (xyz[2]>fgkSol2DipZ || fDipoleOFF) ? bz*fFactorSol : bz*fFactorDip;
 229   }
 230   else return 0.;
 231 }
 232
 233 //_______________________________________________________________________
 234 AliMagF& AliMagF::operator=(const AliMagF& src)
 235 {
 236   if (this != &src) {
 237     if (src.fMeasuredMap) {
 238       if (fMeasuredMap) delete fMeasuredMap;
 239       fMeasuredMap = new AliMagWrapCheb(*src.fMeasuredMap);
 240     }
 241     SetName(src.GetName());
 242     fSolenoid    = src.fSolenoid;
 243     fBeamType    = src.fBeamType;
 244     fBeamEnergy  = src.fBeamEnergy;
 245     fInteg       = src.fInteg;
 246     fPrecInteg   = src.fPrecInteg;
 247     fFactorSol   = src.fFactorSol;
 248     fFactorDip   = src.fFactorDip;
 249     fMax         = src.fMax;
 250     fDipoleOFF   = src.fDipoleOFF;
 251     fParNames    = src.fParNames;
 252   }
 253   return *this;
 254 }
 255
 256 //_______________________________________________________________________
 257 void AliMagF::InitMachineField(BeamType_t btype, Double_t benergy)
 258 {
 259   if (btype==kNoBeamField) {
 260     fQuadGradient = fDipoleField = fCCorrField = fACorr1Field = fACorr2Field = 0.;
 261     return;
 262   }
 263   //
 264   double rigScale = benergy/7000.;   // scale according to ratio of E/Enominal
 265   // for ions assume PbPb (with energy provided per nucleon) and account for A/Z
 266   if (btype==kBeamTypeAA/* || btype==kBeamTypepA || btype==kBeamTypeAp */) rigScale *= 208./82.;
 267   // Attention: in p-Pb the energy recorded in the GRP is the PROTON energy, no rigidity
 268   // rescaling is needed
 269   //
 270   fQuadGradient = 22.0002*rigScale;
 271   fDipoleField  = 37.8781*rigScale;
 272   //
 273   // SIDE C
 274   fCCorrField   = -9.6980;
 275   // SIDE A
 276   fACorr1Field  = -13.2247;
 277   fACorr2Field  =  11.7905;
 278   //
 279 }
 280
 281 //_______________________________________________________________________
 282 void AliMagF::MachineField(const Double_t *x, Double_t *b) const
 283 {
 284   // ---- This is the ZDC part
 285   // Compansators for Alice Muon Arm Dipole
 286   const Double_t kBComp1CZ = 1075., kBComp1hDZ = 260./2., kBComp1SqR = 4.0*4.0;
 287   const Double_t kBComp2CZ = 2049., kBComp2hDZ = 153./2., kBComp2SqR = 4.5*4.5;
 288   //
 289   const Double_t kTripQ1CZ = 2615., kTripQ1hDZ = 637./2., kTripQ1SqR = 3.5*3.5;
 290   const Double_t kTripQ2CZ = 3480., kTripQ2hDZ = 550./2., kTripQ2SqR = 3.5*3.5;
 291   const Double_t kTripQ3CZ = 4130., kTripQ3hDZ = 550./2., kTripQ3SqR = 3.5*3.5;
 292   const Double_t kTripQ4CZ = 5015., kTripQ4hDZ = 637./2., kTripQ4SqR = 3.5*3.5;
 293   //
 294   const Double_t kDip1CZ = 6310.8,  kDip1hDZ = 945./2., kDip1SqRC = 4.5*4.5, kDip1SqRA = 3.375*3.375;
 295   const Double_t kDip2CZ = 12640.3, kDip2hDZ = 945./2., kDip2SqRC = 4.5*4.5, kDip2SqRA = 3.75*3.75;
 296   const Double_t kDip2DXC = 9.7, kDip2DXA = 9.4;
 297   //
 298   double rad2 = x[0] * x[0] + x[1] * x[1];
 299   //
 300   b[0] = b[1] = b[2] = 0;
 301   //
 302   // SIDE C **************************************************
 303   if(x[2]<0.){
 304     if(TMath::Abs(x[2]+kBComp2CZ)<kBComp2hDZ && rad2 < kBComp2SqR){
 305       b[0] = fCCorrField*fFactorDip;
 306     }
 307     else if(TMath::Abs(x[2]+kTripQ1CZ)<kTripQ1hDZ && rad2 < kTripQ1SqR){
 308       b[0] = fQuadGradient*x[1];
 309       b[1] = fQuadGradient*x[0];
 310     }
 311     else if(TMath::Abs(x[2]+kTripQ2CZ)<kTripQ2hDZ && rad2 < kTripQ2SqR){
 312       b[0] = -fQuadGradient*x[1];
 313       b[1] = -fQuadGradient*x[0];
 314     }
 315     else if(TMath::Abs(x[2]+kTripQ3CZ)<kTripQ3hDZ && rad2 < kTripQ3SqR){
 316       b[0] = -fQuadGradient*x[1];
 317       b[1] = -fQuadGradient*x[0];
 318     }
 319     else if(TMath::Abs(x[2]+kTripQ4CZ)<kTripQ4hDZ && rad2 < kTripQ4SqR){
 320       b[0] = fQuadGradient*x[1];
 321       b[1] = fQuadGradient*x[0];
 322     }
 323     else if(TMath::Abs(x[2]+kDip1CZ)<kDip1hDZ && rad2 < kDip1SqRC){
 324       b[1] = fDipoleField;
 325     }
 326     else if(TMath::Abs(x[2]+kDip2CZ)<kDip2hDZ && rad2 < kDip2SqRC) {
 327       double dxabs = TMath::Abs(x[0])-kDip2DXC;
 328       if ( (dxabs*dxabs + x[1]*x[1])<kDip2SqRC) {
 329         b[1] = -fDipoleField;
 330       }
 331     }
 332   }
 333   //
 334   // SIDE A **************************************************
 335   else{
 336     if(TMath::Abs(x[2]-kBComp1CZ)<kBComp1hDZ && rad2 < kBComp1SqR) {
 337       // Compensator magnet at z = 1075 m
 338       b[0] = fACorr1Field*fFactorDip;
 339     }
 340     //
 341     if(TMath::Abs(x[2]-kBComp2CZ)<kBComp2hDZ && rad2 < kBComp2SqR){
 342       b[0] = fACorr2Field*fFactorDip;
 343     }
 344     else if(TMath::Abs(x[2]-kTripQ1CZ)<kTripQ1hDZ && rad2 < kTripQ1SqR){
 345       b[0] = -fQuadGradient*x[1];
 346       b[1] = -fQuadGradient*x[0];
 347     }
 348     else if(TMath::Abs(x[2]-kTripQ2CZ)<kTripQ2hDZ && rad2 < kTripQ2SqR){
 349       b[0] =  fQuadGradient*x[1];
 350       b[1] =  fQuadGradient*x[0];
 351     }
 352     else if(TMath::Abs(x[2]-kTripQ3CZ)<kTripQ3hDZ && rad2 < kTripQ3SqR){
 353       b[0] =  fQuadGradient*x[1];
 354       b[1] =  fQuadGradient*x[0];
 355     }
 356     else if(TMath::Abs(x[2]-kTripQ4CZ)<kTripQ4hDZ && rad2 < kTripQ4SqR){
 357       b[0] = -fQuadGradient*x[1];
 358       b[1] = -fQuadGradient*x[0];
 359     }
 360     else if(TMath::Abs(x[2]-kDip1CZ)<kDip1hDZ && rad2 < kDip1SqRA){
 361       b[1] = -fDipoleField;
 362     }
 363     else if(TMath::Abs(x[2]-kDip2CZ)<kDip2hDZ && rad2 < kDip2SqRA) {
 364       double dxabs = TMath::Abs(x[0])-kDip2DXA;
 365       if ( (dxabs*dxabs + x[1]*x[1])<kDip2SqRA) {
 366         b[1] = fDipoleField;
 367       }
 368     }
 369   }
 370   //
 371 }
 372
 373 //_______________________________________________________________________
 374 void AliMagF::GetTPCInt(const Double_t *xyz, Double_t *b) const
 375 {
 376   // Method to calculate the integral_0^z of br,bt,bz
 377   b[0]=b[1]=b[2]=0.0;
 378   if (fMeasuredMap) {
 379     fMeasuredMap->GetTPCInt(xyz,b);
 380     for (int i=3;i--;) b[i] *= fFactorSol;
 381   }
 382 }
 383
 384 //_______________________________________________________________________
 385 void AliMagF::GetTPCRatInt(const Double_t *xyz, Double_t *b) const
 386 {
 387   // Method to calculate the integral_0^z of bx/bz,by/bz and (bx/bz)^2+(by/bz)^2
 388   b[0]=b[1]=b[2]=0.0;
 389   if (fMeasuredMap) {
 390     fMeasuredMap->GetTPCRatInt(xyz,b);
 391     b[2] /= 100;
 392   }
 393 }
 394
 395 //_______________________________________________________________________
 396 void AliMagF::GetTPCIntCyl(const Double_t *rphiz, Double_t *b) const
 397 {
 398   // Method to calculate the integral_0^z of br,bt,bz
 399   // in cylindrical coordiates ( -pi<phi<pi convention )
 400   b[0]=b[1]=b[2]=0.0;
 401   if (fMeasuredMap) {
 402     fMeasuredMap->GetTPCIntCyl(rphiz,b);
 403     for (int i=3;i--;) b[i] *= fFactorSol;
 404   }
 405 }
 406
 407 //_______________________________________________________________________
 408 void AliMagF::GetTPCRatIntCyl(const Double_t *rphiz, Double_t *b) const
 409 {
 410   // Method to calculate the integral_0^z of bx/bz,by/bz and (bx/bz)^2+(by/bz)^2
 411   // in cylindrical coordiates ( -pi<phi<pi convention )
 412   b[0]=b[1]=b[2]=0.0;
 413   if (fMeasuredMap) {
 414     fMeasuredMap->GetTPCRatIntCyl(rphiz,b);
 415     b[2] /= 100;
 416   }
 417 }
 418
 419 //_______________________________________________________________________
 420 void AliMagF::SetFactorSol(Float_t fc)
 421 {
 422   // set the sign/scale of the current in the L3 according to fgkPolarityConvention
 423   switch (fgkPolarityConvention) {
 424   case kConvDCS2008: fFactorSol = -fc; break;
 425   case kConvLHC    : fFactorSol = -fc; break;
 426   default          : fFactorSol =  fc; break;  // case kConvMap2005: fFactorSol =  fc; break;
 427   }
 428 }
 429
 430 //_______________________________________________________________________
 431 void AliMagF::SetFactorDip(Float_t fc)
 432 {
 433   // set the sign*scale of the current in the Dipole according to fgkPolarityConvention
 434   switch (fgkPolarityConvention) {
 435   case kConvDCS2008: fFactorDip =  fc; break;
 436   case kConvLHC    : fFactorDip = -fc; break;
 437   default          : fFactorDip =  fc; break;  // case kConvMap2005: fFactorDip =  fc; break;
 438   }
 439 }
 440
 441 //_______________________________________________________________________
 442 Double_t AliMagF::GetFactorSol() const
 443 {
 444   // return the sign*scale of the current in the Dipole according to fgkPolarityConventionthe
 445   switch (fgkPolarityConvention) {
 446   case kConvDCS2008: return -fFactorSol;
 447   case kConvLHC    : return -fFactorSol;
 448   default          : return  fFactorSol;       //  case kConvMap2005: return  fFactorSol;
 449   }
 450 }
 451
 452 //_______________________________________________________________________
 453 Double_t AliMagF::GetFactorDip() const
 454 {
 455   // return the sign*scale of the current in the Dipole according to fgkPolarityConventionthe
 456   switch (fgkPolarityConvention) {
 457   case kConvDCS2008: return  fFactorDip;
 458   case kConvLHC    : return -fFactorDip;
 459   default          : return  fFactorDip;       //  case kConvMap2005: return  fFactorDip;
 460   }
 461 }
 462
 463 //_____________________________________________________________________________
 464 AliMagF* AliMagF::CreateFieldMap(Float_t l3Cur, Float_t diCur, Int_t convention, Bool_t uniform,
 465                                  Float_t beamenergy, const Char_t *beamtype, const Char_t *path)
 466 {
 467   //------------------------------------------------
 468   // The magnetic field map, defined externally...
 469   // L3 current 30000 A  -> 0.5 T
 470   // L3 current 12000 A  -> 0.2 T
 471   // dipole current 6000 A
 472   // The polarities must match the convention (LHC or DCS2008)
 473   // unless the special uniform map was used for MC
 474   //------------------------------------------------
 475   const Float_t l3NominalCurrent1=30000.; // (A)
 476   const Float_t l3NominalCurrent2=12000.; // (A)
 477   const Float_t diNominalCurrent =6000. ; // (A)
 478
 479   const Float_t tolerance=0.03; // relative current tolerance
 480   const Float_t zero=77.;       // "zero" current (A)
 481   //
 482   BMap_t map = k5kG;
 483   double sclL3,sclDip;
 484   //
 485   Float_t l3Pol = l3Cur > 0 ? 1:-1;
 486   Float_t diPol = diCur > 0 ? 1:-1;
 487
 488   l3Cur = TMath::Abs(l3Cur);
 489   diCur = TMath::Abs(diCur);
 490   //
 491   if (TMath::Abs((sclDip=diCur/diNominalCurrent)-1.) > tolerance && !uniform) {
 492     if (diCur <= zero) sclDip = 0.; // some small current.. -> Dipole OFF
 493     else {
 494       AliFatalGeneral("AliMagF",Form("Wrong dipole current (%f A)!",diCur));
 495     }
 496   }
 497   //
 498   if (uniform) {
 499     // special treatment of special MC with uniform mag field (normalized to 0.5 T)
 500     // no check for scaling/polarities are done
 501     map   = k5kGUniform;
 502     sclL3 = l3Cur/l3NominalCurrent1;
 503   }
 504   else {
 505     if      (TMath::Abs((sclL3=l3Cur/l3NominalCurrent1)-1.) < tolerance) map  = k5kG;
 506     else if (TMath::Abs((sclL3=l3Cur/l3NominalCurrent2)-1.) < tolerance) map  = k2kG;
 507     else if (l3Cur <= zero && diCur<=zero)   { sclL3=0; sclDip=0; map  = k5kGUniform;}
 508     else {
 509       AliFatalGeneral("AliMagF",Form("Wrong L3 current (%f A)!",l3Cur));
 510     }
 511   }
 512   //
 513   if (sclDip!=0 && map!=k5kGUniform) {
 514     if ( (l3Cur<=zero) || ((convention==kConvLHC && l3Pol!=diPol) || (convention==kConvDCS2008 && l3Pol==diPol)) ) {
 515       AliFatalGeneral("AliMagF",Form("Wrong combination for L3/Dipole polarities (%c/%c) for convention %d",
 516                                      l3Pol>0?'+':'-',diPol>0?'+':'-',GetPolarityConvention()));
 517     }
 518   }
 519   //
 520   if (l3Pol<0) sclL3  = -sclL3;
 521   if (diPol<0) sclDip = -sclDip;
 522   //
 523   BeamType_t btype = kNoBeamField;
 524   TString btypestr = beamtype;
 525   btypestr.ToLower();
 526   TPRegexp protonBeam("(proton|p)\\s*-?\\s*\\1");
 527   TPRegexp ionBeam("(lead|pb|ion|a|A)\\s*-?\\s*\\1");
 528   TPRegexp protonionBeam("(proton|p)\\s*-?\\s*(lead|pb|ion|a|A)");
 529   TPRegexp ionprotonBeam("(lead|pb|ion|a|A)\\s*-?\\s*(proton|p)");
 530   if (btypestr.Contains(ionBeam)) btype = kBeamTypeAA;
 531   else if (btypestr.Contains(protonBeam)) btype = kBeamTypepp;
 532   else if (btypestr.Contains(protonionBeam)) btype = kBeamTypepA;
 533   else if (btypestr.Contains(ionprotonBeam)) btype = kBeamTypeAp;
 534   else AliInfoGeneral("AliMagF",Form("Assume no LHC magnet field for the beam type %s, ",beamtype));
 535   char ttl[80];
 536   snprintf(ttl,79,"L3: %+5d Dip: %+4d kA; %s | Polarities in %s convention",(int)TMath::Sign(l3Cur,float(sclL3)),
 537           (int)TMath::Sign(diCur,float(sclDip)),uniform ? " Constant":"",
 538           convention==kConvLHC ? "LHC":"DCS2008");
 539   // LHC and DCS08 conventions have opposite dipole polarities
 540   if ( GetPolarityConvention() != convention) sclDip = -sclDip;
 541   //
 542   return new AliMagF("MagneticFieldMap", ttl,sclL3,sclDip,map,btype,beamenergy,2,10.,path);
 543   //
 544 }
 545
 546 //_____________________________________________________________________________
 547 const char*  AliMagF::GetBeamTypeText() const
 548 {
 549   // beam type in text form
 550   const char *beamNA  = "No Beam";
 551   const char *beamPP  = "p-p";
 552   const char *beamPbPb= "A-A";
 553   const char *beamPPb = "p-A";
 554   const char *beamPbP = "A-p";
 555   switch ( fBeamType ) {
 556   case kBeamTypepp : return beamPP;
 557   case kBeamTypeAA : return beamPbPb;
 558   case kBeamTypepA : return beamPPb;
 559   case kBeamTypeAp : return beamPbP;
 560   case kNoBeamField:
 561   default:           return beamNA;
 562   }
 563 }
 564
 565 //_____________________________________________________________________________
 566 void AliMagF::Print(Option_t *opt) const
 567 {
 568   // print short or long info
 569   TString opts = opt; opts.ToLower();
 570   AliInfo(Form("%s:%s",GetName(),GetTitle()));
 571   AliInfo(Form("Solenoid (%+.2f*)%.0f kG, Dipole %s (%+.2f) %s",
 572                GetFactorSol(),(fMapType==k5kG||fMapType==k5kGUniform)?5.:2.,
 573                fDipoleOFF ? "OFF":"ON",GetFactorDip(),fMapType==k5kGUniform?" |Constant Field!":""));
 574   if (opts.Contains("a")) {
 575     AliInfo(Form("Machine B fields for %s beam (%.0f GeV): QGrad: %.4f Dipole: %.4f",
 576                  GetBeamTypeText(),
 577                  fBeamEnergy,fQuadGradient,fDipoleField));
 578     AliInfo(Form("Uses %s of %s",GetParamName(),GetDataFileName()));
 579   }
 580 }