MagneticField.h

Go to the documentation of this file.

// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
 
 
#ifndef ALICEO2_FIELD_MAGNETICFIELD_H_
#define ALICEO2_FIELD_MAGNETICFIELD_H_
 
#include "FairField.h" // for FairField
#include "Field/MagFieldParam.h"
#include "Field/MagneticWrapperChebyshev.h" // for MagneticWrapperChebyshev
#include "Field/MagFieldFast.h"
#include "TSystem.h"
#include "Rtypes.h" // for Double_t, Char_t, Int_t, Float_t, etc
#include "TNamed.h" // for TNamed
#include <memory>   // for str::unique_ptr
 
class FairParamList;
 
namespace o2
{
namespace field
{
class MagneticWrapperChebyshev;
}
} // namespace o2
namespace o2
{
namespace field
{
 
class MagneticField : public FairField
{
 
 public:
  enum PolarityConvention_t { kConvLHC,
                              kConvDCS2008,
                              kConvMap2005 };
  enum { kOverrideGRP = BIT(14) }; // don't recreate from GRP if set
 
  MagneticField();
 
  MagneticField(const char* name, const char* title, Double_t factorSol = 1., Double_t factorDip = 1.,
                MagFieldParam::BMap_t maptype = MagFieldParam::k5kG,
                MagFieldParam::BeamType_t btype = MagFieldParam::kBeamTypepp, Double_t benergy = -1, Int_t integ = 2,
                Double_t fmax = 15, const std::string path = "$(O2_ROOT)/share/Common/maps/mfchebKGI_sym.root");
 
  MagneticField(const MagFieldParam& param);
 
  MagneticField& operator=(const MagneticField& src);
 
  ~MagneticField() override = default;
 
  static MagneticField* createNominalField(int fld, bool uniform = false);
 
  void CreateField();
 
  void AllowFastField(bool v = true);
 
  bool fastFieldExists() const
  {
    return !(mMapType == MagFieldParam::k5kGUniform || mDipoleOnOffFlag == true);
  }
 
  void rescaleField(float l3Cur, float diCur, bool uniform, int convention = 0);
 
 
  Double_t GetBx(Double_t x, Double_t y, Double_t z) override
  {
    double xyz[3] = {x, y, z}, b[3];
    MagneticField::Field(xyz, b);
    return b[0];
  }
 
  Double_t GetBy(Double_t x, Double_t y, Double_t z) override
  {
    double xyz[3] = {x, y, z}, b[3];
    MagneticField::Field(xyz, b);
    return b[1];
  }
 
  Double_t GetBz(Double_t x, Double_t y, Double_t z) override
  {
    double xyz[3] = {x, y, z};
    return getBz(xyz);
  }
 
  void Field(const Double_t* __restrict__ point, Double_t* __restrict__ bField) override;
 
  void field(const math_utils::Point3D<float> xyz, float bxyz[3])
  {
    double xyzd[3] = {xyz.X(), xyz.Y(), xyz.Z()}, bxyzd[3] = {0};
    Field(xyzd, bxyzd);
    bxyz[0] = bxyzd[0];
    bxyz[1] = bxyzd[1];
    bxyz[2] = bxyzd[2];
  }
 
  void field(const math_utils::Point3D<double> xyz, double bxyz[3])
  {
    double xyzd[3] = {xyz.X(), xyz.Y(), xyz.Z()};
    Field(xyzd, bxyz);
  }
 
  void field(const double* __restrict__ point, double* __restrict__ bField)
  {
    Field(point, bField);
  }
 
  void field(const float* __restrict__ point, float* __restrict__ bField)
  {
    double xyz[3] = {point[0], point[1], point[2]}, bxyz[3] = {0};
    Field(xyz, bxyz);
    bField[0] = bxyz[0];
    bField[1] = bxyz[1];
    bField[2] = bxyz[2];
  }
 
  void GetBxyz(const Double_t p[3], Double_t* b) override { MagneticField::Field(p, b); }
 
  void FillParContainer() override;
 
  void getTPCIntegral(const Double_t* xyz, Double_t* b) const;
 
  void getTPCRatIntegral(const Double_t* xyz, Double_t* b) const;
 
  void getTPCIntegralCylindrical(const Double_t* rphiz, Double_t* b) const;
 
  void getTPCRatIntegralCylindrical(const Double_t* rphiz, Double_t* b) const;
 
  Double_t getBz(const Double_t* xyz) const;
 
  MagneticWrapperChebyshev* getMeasuredMap() const { return mMeasuredMap.get(); }
 
  const MagFieldFast* getFastField() const { return mFastField.get(); }
 
  // Former MagF methods or their aliases
 
  void setFactorSolenoid(float fc = 1.);
 
  void setFactorDipole(float fc = 1.);
 
  Double_t getFactorSolenoid() const;
 
  Double_t getFactorDipole() const;
 
  Double_t Factor() const { return getFactorSolenoid(); }
 
  Double_t getCurrentSolenoid() const
  {
    return getFactorSolenoid() * (mMapType == MagFieldParam::k2kG ? 12000 : 30000);
  }
 
  Double_t getCurrentDipole() const { return getFactorDipole() * 6000; }
 
  Bool_t IsUniform() const { return mMapType == MagFieldParam::k5kGUniform; }
 
  void MachineField(const Double_t* __restrict__ x, Double_t* __restrict__ b) const;
 
  MagFieldParam::BMap_t getMapType() const { return mMapType; }
 
  MagFieldParam::BeamType_t getBeamType() const { return mBeamType; }
 
  const char* getBeamTypeText() const;
 
  Double_t getBeamEnergy() const { return mBeamEnergy; }
 
  Double_t Max() const { return mMaxField; }
 
  Int_t Integral() const { return mDefaultIntegration; }
 
  Int_t precIntegral() const { return mPrecisionInteg; }
 
  Double_t solenoidField() const { return mSolenoid; }
 
  Char_t* getDataFileName() const { return (Char_t*)mParameterNames.GetName(); }
 
  Char_t* getParameterName() const { return (Char_t*)mParameterNames.GetTitle(); }
 
  void setDataFileName(const Char_t* nm) { mParameterNames.SetName(nm); }
 
  void setParameterName(const Char_t* nm) { mParameterNames.SetTitle(nm); }
 
  void Print(Option_t* opt) const override;
 
  Bool_t loadParameterization();
 
  static Int_t getPolarityConvention() { return Int_t(sPolarityConvention); }
 
  static MagFieldParam::BMap_t getFieldMapScale(float& l3, float& dip, bool uniform, int convention = 0);
 
  static MagneticField* createFieldMap(float l3Current = -30000., float diCurrent = -6000., Int_t convention = 0,
                                       Bool_t uniform = kFALSE, float beamenergy = 7000, const Char_t* btype = "pp",
                                       const std::string path = std::string(gSystem->Getenv("VMCWORKDIR")) +
                                                                std::string("/Common/maps/mfchebKGI_sym.root"));
 
 protected:
  // not supposed to be changed during the run, set only at the initialization via constructor
  void initializeMachineField(MagFieldParam::BeamType_t btype, Double_t benergy);
 
  void setBeamType(MagFieldParam::BeamType_t type) { mBeamType = type; }
 
  void setBeamEnergy(float energy) { mBeamEnergy = energy; }
 
 private:
  std::unique_ptr<MagneticWrapperChebyshev> mMeasuredMap; 
  std::unique_ptr<MagFieldFast> mFastField;               // ! optional fast parametrization
  MagFieldParam::BMap_t mMapType;                         
  Double_t mSolenoid;                                     
  MagFieldParam::BeamType_t mBeamType;                    
  Double_t mBeamEnergy;                                   
 
  Int_t mDefaultIntegration;             
  Int_t mPrecisionInteg;                 
  Double_t mMultipicativeFactorSolenoid; 
  Double_t mMultipicativeFactorDipole;   
  Double_t mMaxField;                    
  Bool_t mDipoleOnOffFlag;               
 
  Double_t mQuadrupoleGradient; 
  Double_t mDipoleField;        
  Double_t mCompensatorField2C; 
  Double_t mCompensatorField1A; 
  Double_t mCompensatorField2A; 
 
  TNamed mParameterNames; 
 
  static const Double_t sSolenoidToDipoleZ;  
  static const UShort_t sPolarityConvention; 
 
  MagneticField(const MagneticField& src);
 
  ClassDefOverride(o2::field::MagneticField,
                   3) // Class for all Alice MagField wrapper for measured data + Tosca parameterization
};
} // namespace field
} // namespace o2
 
#endif