From 200e4de2ed7984052b846627d7a90efa58be1eff Mon Sep 17 00:00:00 2001 From: jose-luis-rs Date: Thu, 23 Feb 2023 21:37:39 +0100 Subject: [PATCH] GLAD: Added new setter SetFieldfromCurrent to scale the field Minor update Finishing the update Revert clang-format Revert bools clang-format --- field/R3BGladFieldMap.cxx | 209 +++++++++------ field/R3BGladFieldMap.h | 321 ++++++++++++----------- glad/R3BGladMagnet.cxx | 2 +- r3bbase/R3BLogger.h | 37 ++- r3bgen/R3BINCLRootGenerator.h | 8 +- twim/ana/R3BTwimVertexReconstruction.cxx | 14 +- 6 files changed, 320 insertions(+), 271 deletions(-) diff --git a/field/R3BGladFieldMap.cxx b/field/R3BGladFieldMap.cxx index 1c9d67620..b84d863bf 100644 --- a/field/R3BGladFieldMap.cxx +++ b/field/R3BGladFieldMap.cxx @@ -24,6 +24,7 @@ #include "TMath.h" #include "R3BGladFieldMap.h" +#include "R3BLogger.h" using std::cerr; using std::cout; @@ -55,14 +56,14 @@ R3BGladFieldMap::R3BGladFieldMap() // ------------- Standard constructor --------------------------------- // -R3BGladFieldMap::R3BGladFieldMap(const char* mapName, const char* fileType) - : FairField(mapName) +R3BGladFieldMap::R3BGladFieldMap(const TString& mapName, const TString& fileType) + : FairField(mapName.Data()) { // Default field positions (in cm) in lab: // between target position (0,0,0) and GLAD rotation point (field origin) // Override these values by calling SetPosition(x,y,z) before Init() fPosX = 0.; - fPosY = 0.; + fPosY = 1.75; fPosZ = 163.4; // Default translation vector of the local filed coordiantes @@ -74,39 +75,73 @@ R3BGladFieldMap::R3BGladFieldMap(const char* mapName, const char* fileType) fYAngle = -14.; fZAngle = 0.; - fXmin = fYmin = fZmin = 0.; - fXmax = fYmax = fZmax = 0.; - fXstep = fYstep = fZstep = 0.; - fNx = fNy = fNz = 0; + fXmin = 0; + fYmin = 0; + fZmin = 0.; + fXmax = 0; + fYmax = 0; + fZmax = 0.; + fXstep = 0; + fYstep = 0; + fZstep = 0.; + fNx = 0; + fNy = 0; + fNz = 0; fScale = 1.; - fBx = fBy = fBz = NULL; + fBx = nullptr; + fBy = nullptr; + fBz = nullptr; fName = mapName; TString dir = getenv("VMCWORKDIR"); fFileName = dir + "/field/magField/R3B/" + mapName; if (fileType[0] == 'R') + { fFileName += ".root"; + } else + { fFileName += ".dat"; + } fType = 2; } +void R3BGladFieldMap::SetFieldfromCurrent(double current) +{ + R3BLOG_IF(fatal, TMath::Abs(current) > 3583.81, "GLAD current cannot be larger than 3583.81 A."); + fScale = current / 3583.81; + R3BLOG(info, "GLAD current set to " << current << " A, which corresponds to a scaling factor of " << fScale); + return; +} + // ------------ Constructor from R3BGladFieldPar ----------------------- // R3BGladFieldMap::R3BGladFieldMap(R3BFieldPar* fieldPar) { fType = 2; - fPosX = fPosY = fPosZ = 0.; - fXAngle = fYAngle = fZAngle = 0.; - fXmin = fYmin = fZmin = 0.; - fXmax = fYmax = fZmax = 0.; - fXstep = fYstep = fZstep = 0.; - fNx = fNy = fNz = 0; + fPosX = 0; + fPosY = 0; + fPosZ = 0.; + fXAngle = 0; + fYAngle = 0; + fZAngle = 0.; + fXmin = 0; + fYmin = 0; + fZmin = 0.; + fXmax = 0; + fYmax = 0; + fZmax = 0.; + fXstep = 0; + fYstep = 0; + fZstep = 0.; + fNx = 0; + fNy = 0; + fNz = 0; fScale = 1.; fBx = fBy = fBz = NULL; if (!fieldPar) { - cerr << "-W- R3BGladFieldConst::R3BGladFieldMap: empty parameter container!" << endl; + R3BLOG(warn, "empty parameter container!"); fName = ""; fFileName = ""; fType = 1; @@ -150,8 +185,7 @@ void R3BGladFieldMap::Init() ReadRootFile(fFileName); else { - cerr << "-E- R3BGladFieldMap::Init: No proper file name defined! (" << fFileName << ")" << endl; - LOG(fatal) << "Init: No proper file name"; + R3BLOG(fatal, "No proper file name defined! (" << fFileName.Data() << ")"); } Print(); } @@ -160,7 +194,7 @@ void R3BGladFieldMap::Init() Double_t R3BGladFieldMap::GetBx(Double_t x, Double_t y, Double_t z) { TVector3 B = GetBtrans(x, y, z); - Double_t val = B.X() * 10.000000; + Double_t val = B.X() * 10.0; return (val); // should be in kGaus units } // ----------- Get y component of the field --------------------------- @@ -168,7 +202,7 @@ Double_t R3BGladFieldMap::GetBx(Double_t x, Double_t y, Double_t z) Double_t R3BGladFieldMap::GetBy(Double_t x, Double_t y, Double_t z) { TVector3 B = GetBtrans(x, y, z); - Double_t val = B.Y() * 10.000000; + Double_t val = B.Y() * 10.0; return (val); // should be in kGaus units } @@ -177,7 +211,7 @@ Double_t R3BGladFieldMap::GetBy(Double_t x, Double_t y, Double_t z) Double_t R3BGladFieldMap::GetBz(Double_t x, Double_t y, Double_t z) { TVector3 B = GetBtrans(x, y, z); - Double_t val = B.Z() * 10.000000; + Double_t val = B.Z() * 10.0; return (val); // should be in kGaus units } @@ -317,14 +351,14 @@ Bool_t R3BGladFieldMap::IsInside(Double_t x, // ---------- Write the map to an ASCII file -------------------------- // -void R3BGladFieldMap::WriteAsciiFile(const char* fileName) +void R3BGladFieldMap::WriteAsciiFile(const TString& fileName) { // Open file - cout << "-I- R3BGladFieldMap: Writing field map to ASCII file " << fileName << endl; + R3BLOG(info, "Writing field map to ASCII file " << fileName.Data()); ofstream mapFile(fileName); if (!mapFile.is_open()) { - cerr << "-E- R3BGladFieldMap:WriteAsciiFile: Could not open file! " << endl; + R3BLOG(error, "Could not open file!"); return; } @@ -345,7 +379,7 @@ void R3BGladFieldMap::WriteAsciiFile(const char* fileName) Double_t factor = 10. * fScale; // Takes out scaling and converts kG->T cout << right; Int_t nTot = fNx * fNy * fNz; - cout << "-I- R3BGladFieldMap: " << fNx * fNy * fNz << " entries to write... " << setw(3) << 0 << " % "; + R3BLOG(info, fNx * fNy * fNz << " entries to write... " << setw(3) << 0 << " % "); Int_t index = 0; div_t modul; Int_t iDiv = TMath::Nint(nTot / 100.); @@ -367,7 +401,7 @@ void R3BGladFieldMap::WriteAsciiFile(const char* fileName) } // z-Loop } // y-Loop } // x-Loop - cout << " " << index + 1 << " written" << endl; + R3BLOG(info, " " << index + 1 << " written"); mapFile.close(); } @@ -383,46 +417,65 @@ void R3BGladFieldMap::SetPosition(Double_t x, Double_t y, Double_t z) // --------- Screen output -------------------------------------------- // -void R3BGladFieldMap::Print(Option_t* option) const +void R3BGladFieldMap::Print(Option_t*) const { TString type = "Map"; if (fType == 2) + { type = "Map sym2"; - if (fType == 3) + } + else if (fType == 3) + { type = "Map sym3"; - cout << "======================================================" << endl; - cout.precision(4); - cout << showpoint; - cout << "---- " << fTitle << " : " << fName << endl; - cout << "----" << endl; - cout << "---- Field type : " << type << endl; - cout << "----" << endl; - cout << "---- Field map grid : " << endl; - cout << "---- x = " << setw(4) << fXmin << " to " << setw(4) << fXmax << " cm, " << fNx - << " grid points, dx = " << fXstep << " cm" << endl; - cout << "---- y = " << setw(4) << fYmin << " to " << setw(4) << fYmax << " cm, " << fNy - << " grid points, dy = " << fYstep << " cm" << endl; - cout << "---- z = " << setw(4) << fZmin << " to " << setw(4) << fZmax << " cm, " << fNz - << " grid points, dz = " << fZstep << " cm" << endl; - cout << endl; - cout << "---- Field centre position: ( " << setw(6) << fPosX << ", " << setw(6) << fPosY << ", " << setw(6) - << fPosZ << ") cm" << endl; - cout << "---- Field rotation X: " << setw(6) << fXAngle << " deg" << endl; - cout << "---- Field rotation Y: " << setw(6) << fYAngle << " deg" << endl; - cout << "---- Field rotation Z: " << setw(6) << fZAngle << " deg" << endl; - cout << "---- Field scaling factor: " << fScale << endl; - cout << "======================================================" << endl; + } + + std::stringstream sprint; + sprint << endl; + sprint << "======================================================" << endl; + sprint.precision(4); + sprint << showpoint; + sprint << "---- " << fTitle << " : " << fName << endl; + sprint << "----" << endl; + sprint << "---- Field type : " << type << endl; + sprint << "----" << endl; + sprint << "---- Field map grid : " << endl; + sprint << "---- x = " << setw(4) << fXmin << " to " << setw(4) << fXmax << " cm, " << fNx + << " grid points, dx = " << fXstep << " cm" << endl; + sprint << "---- y = " << setw(4) << fYmin << " to " << setw(4) << fYmax << " cm, " << fNy + << " grid points, dy = " << fYstep << " cm" << endl; + sprint << "---- z = " << setw(4) << fZmin << " to " << setw(4) << fZmax << " cm, " << fNz + << " grid points, dz = " << fZstep << " cm" << endl; + sprint << endl; + sprint << "---- Field centre position: ( " << setw(6) << fPosX << ", " << setw(6) << fPosY << ", " << setw(6) + << fPosZ << ") cm" << endl; + sprint << "---- Field rotation X: " << setw(6) << fXAngle << " deg" << endl; + sprint << "---- Field rotation Y: " << setw(6) << fYAngle << " deg" << endl; + sprint << "---- Field rotation Z: " << setw(6) << fZAngle << " deg" << endl; + sprint << "---- Field scaling factor: " << fScale << endl; + sprint << "======================================================" << endl; + + R3BLOG(info, sprint.str()); } // --------- Reset parameters and data (private) ---------------------- // void R3BGladFieldMap::Reset() { - fPosX = fPosY = fPosZ = 0.; - fXmin = fYmin = fZmin = 0.; - fXmax = fYmax = fZmax = 0.; - fXstep = fYstep = fZstep = 0.; - fNx = fNy = fNz = 0; + fPosX = 0.; + fPosY = 0.; + fPosZ = 0.; + fXmin = 0.; + fYmin = 0.; + fZmin = 0.; + fXmax = 0.; + fYmax = 0.; + fZmax = 0.; + fXstep = 0.; + fYstep = 0.; + fZstep = 0.; + fNx = 0; + fNy = 0; + fNz = 0; fScale = 1.; if (fBx) { @@ -443,29 +496,34 @@ void R3BGladFieldMap::Reset() // ----- Read field map from ASCII file (private) --------------------- // -void R3BGladFieldMap::ReadAsciiFile(const char* fileName) +void R3BGladFieldMap::ReadAsciiFile(const TString& fileName) { Double_t bx = 0., by = 0., bz = 0.; Double_t ax = 0., ay = 0., az = 0.; // Open file - cout << "-I- R3BGladFieldMap: Reading field map from ASCII file " << fileName << endl; + R3BLOG(info, "Reading field map from ASCII file " << fileName.Data()); ifstream mapFile(fileName); if (!mapFile.is_open()) { - cerr << "-E- R3BGladFieldMap:ReadAsciiFile: Could not open file! " << endl; - LOG(fatal) << "ReadAsciiFile: Could not open file"; + R3BLOG(fatal, "Could not open file! "); } // Read map type TString type; mapFile >> type; Int_t iType = 0; if (type == "nosym") + { iType = 1; - if (type == "sym2") + } + else if (type == "sym2") + { iType = 2; - if (type == "sym3") + } + else if (type == "sym3") + { iType = 3; + } // Read grid parameters mapFile >> fXmin >> fXmax >> fNx; @@ -487,7 +545,7 @@ void R3BGladFieldMap::ReadAsciiFile(const char* fileName) Double_t factor = fScale * 10.; // Factor 10 for T -> kG cout << right; Int_t nTot = fNx * fNy * fNz; - cout << "-I- R3BGladFieldMap: " << nTot << " entries to read... " << setw(3) << 0 << " % "; + R3BLOG(info, nTot << " entries to read... " << setw(3) << 0 << " % "); Int_t index = 0; div_t modul; Int_t iDiv = TMath::Nint(nTot / 100.); @@ -498,8 +556,9 @@ void R3BGladFieldMap::ReadAsciiFile(const char* fileName) for (Int_t iz = 0; iz < fNz; iz++) { if (!mapFile.good()) - cerr << "-E- R3BGladFieldMap::ReadAsciiFile: " - << "I/O Error at " << ix << " " << iy << " " << iz << endl; + { + R3BLOG(error, "I/O Error at " << ix << " " << iy << " " << iz); + } index = ix * fNy * fNz + iy * fNz + iz; modul = div(index, iDiv); if (modul.rem == 0) @@ -527,9 +586,9 @@ void R3BGladFieldMap::ReadAsciiFile(const char* fileName) // cout << "-I- " << bx << " : " << by << " : " << bz << " : " << endl; if (mapFile.eof()) { - cerr << endl - << "-E- R3BGladFieldMap::ReadAsciiFile: EOF" - << " reached at " << ix << " " << iy << " " << iz << endl; + R3BLOG(error, + " EOF" + << " reached at " << ix << " " << iy << " " << iz); mapFile.close(); break; } @@ -537,22 +596,20 @@ void R3BGladFieldMap::ReadAsciiFile(const char* fileName) } // y-Loop } // x-Loop - cout << " " << index + 1 << " read" << endl; - + R3BLOG(info, " " << index + 1 << " read"); mapFile.close(); } // ----- Read field map from ROOT file (private) --------------------- // -void R3BGladFieldMap::ReadRootFile(const char* fileName) +void R3BGladFieldMap::ReadRootFile(const TString& fileName) { // Opening root file - cout << "-I- R3BGladFieldMap: Reading field map from ROOT file " << fileName; + R3BLOG(info, "Reading field map from ROOT file " << fileName.Data()); fFile = new TFile(fileName, "READ"); if (!(fFile->IsOpen())) { - cerr << "-E- R3BGladFieldMap::ReadRootFile: Cannot read from file! " << endl; - LOG(fatal) << "ReadRootFile: Cannot read from file"; + R3BLOG(fatal, "Cannot read from file!"); } TTree* fTreeMap = NULL; @@ -560,8 +617,8 @@ void R3BGladFieldMap::ReadRootFile(const char* fileName) if (!fTreeMap) { - cerr << "-E- R3BGladFieldMap::ReadRootFile: no TTree named 'tree' found in the file" << fileName; - LOG(fatal) << "No field map data"; + R3BLOG(error, "No TTree named 'tree' found in the file " << fileName.Data()); + R3BLOG(fatal, "No field map data"); } Double_t tBx, tBy, tBz; // branches @@ -596,7 +653,7 @@ void R3BGladFieldMap::ReadRootFile(const char* fileName) Long64_t Nentries = fTreeMap->GetEntries(); - cout << "\n-I- Reading GLAD field data from root tree" << endl; + R3BLOG(info, "Reading GLAD field data from root tree"); TVector3 fBvec; @@ -613,7 +670,7 @@ void R3BGladFieldMap::ReadRootFile(const char* fileName) fBy->AddAt(fBvec.Y(), ev - 3); fBz->AddAt(fBvec.Z(), ev - 3); } - cout << "\n-I- Finished reading root tree" << endl; + R3BLOG(info, "Finished reading root tree"); return; } diff --git a/field/R3BGladFieldMap.h b/field/R3BGladFieldMap.h index 899322419..75ef20d87 100644 --- a/field/R3BGladFieldMap.h +++ b/field/R3BGladFieldMap.h @@ -17,172 +17,173 @@ #include "FairField.h" #include "R3BFieldPar.h" #include "TRotation.h" -#include "TVector3.h" +#include "TString.h" #include "TTree.h" +#include "TVector3.h" class TArrayD; class R3BGladFieldMap : public FairField { - - public: - /** Default constructor **/ - R3BGladFieldMap(); - - /** Standard constructor - ** @param name Name of field map file without extension - ** @param fileType R = ROOT file, A = ASCII - **/ - R3BGladFieldMap(const char* mapName, const char* fileType = "A"); - - /** Constructor from R3BGladFieldPar **/ - R3BGladFieldMap(R3BFieldPar* fieldPar); - - /** Constructor from R3BGladFieldMapCreator **/ - // R3BGladFieldMap(R3BGladFieldMapCreator* creator); - - /** Destructor **/ - virtual ~R3BGladFieldMap(); - - /** Initialisation (read map from file) **/ - virtual void Init(); - - /** Get the field components at a certain point - ** @param x,y,z Point coordinates (global) [cm] - ** @value Bx,By,Bz Field components [kG] - **/ - virtual Double_t GetBx(Double_t x, Double_t y, Double_t z); - virtual Double_t GetBy(Double_t x, Double_t y, Double_t z); - virtual Double_t GetBz(Double_t x, Double_t y, Double_t z); - - //special function for the field transformation - TVector3 GetBtrans(Double_t x, Double_t y, Double_t z); - - /** Determine whether a point is inside the field map - ** @param x,y,z Point coordinates (global) [cm] - ** @param ix,iy,iz (return) Grid cell - ** @param dx,dy,dz (return) Distance from grid point [cm] if inside - ** @value kTRUE if inside map, else kFALSE - **/ - virtual Bool_t IsInside(Double_t x, - Double_t y, - Double_t z, - Int_t& ix, - Int_t& iy, - Int_t& iz, - Double_t& dx, - Double_t& dy, - Double_t& dz); - - /** Write the field map to an ASCII file **/ - void WriteAsciiFile(const char* fileName); - - /** Set the position (in cm) of the field origin in lab **/ - virtual void SetPosition(Double_t x, Double_t y, Double_t z); - - /** Set a global field scaling factor **/ - virtual void SetScale(Double_t factor) { fScale = factor; } - - /* Set Euler rotation angles of the field (in degrees) - * default fYAngle = -14 deg, fXAngle=0, fZAngle=0 */ - virtual void SetXAngle(Double_t a) { fXAngle = a; }; - virtual void SetYAngle(Double_t a) { fYAngle = a; }; - virtual void SetZAngle(Double_t a) { fZAngle = a; }; - - /** Accessors to field parameters in local coordinate system **/ - Double_t GetXmin() const { return fXmin; } - Double_t GetYmin() const { return fYmin; } - Double_t GetZmin() const { return fZmin; } - Double_t GetXmax() const { return fXmax; } - Double_t GetYmax() const { return fYmax; } - Double_t GetZmax() const { return fZmax; } - Double_t GetXstep() const { return fXstep; } - Double_t GetYstep() const { return fYstep; } - Double_t GetZstep() const { return fZstep; } - Int_t GetNx() const { return fNx; } - Int_t GetNy() const { return fNy; } - Int_t GetNz() const { return fNz; } - - /** Accessor to field centre position in lab system **/ - Double_t GetPositionX() const { return fPosX; } - Double_t GetPositionY() const { return fPosY; } - Double_t GetPositionZ() const { return fPosZ; } - - /** Accessor to field rotation **/ - Double_t GetXAngle() const { return fXAngle; } - Double_t GetYAngle() const { return fYAngle; } - Double_t GetZAngle() const { return fZAngle; } - - /** Accessor to global scaling factor **/ - Double_t GetScale() const { return fScale; } - - /** Accessors to the field value arrays **/ - TArrayD* GetBx() const { return fBx; } - TArrayD* GetBy() const { return fBy; } - TArrayD* GetBz() const { return fBz; } - - /** Accessor to field map file **/ - TString GetFileName() { return fFileName; } - - /** Screen output **/ - virtual void Print(Option_t* option = "") const; - - protected: - /** Reset the field parameters and data **/ - void Reset(); - - /** Read the field map from an ASCII file **/ - void ReadAsciiFile(const char* fileName); - - /** Read field map from a ROOT file **/ - void ReadRootFile(const char* fileName); - - /** Get field values by interpolation of the grid. - ** @param dx,dy,dz Relative distance from grid point [cell units] - **/ - Double_t Interpolate(Double_t dx, Double_t dy, Double_t dz); - - /** Map file name **/ - TString fFileName; - - /** Global scaling factor (w.r.t. map on file) **/ - Double_t fScale; - - /** Field centre position in global coordinates **/ - Double_t fPosX, fPosY, fPosZ; - - /** Euler rotations of the field around field local X, Y, Z axis **/ - Double_t fXAngle; - Double_t fYAngle; - Double_t fZAngle; - - /** Field limits in local coordinate system **/ - Double_t fXmin, fXmax, fXstep; - Double_t fYmin, fYmax, fYstep; - Double_t fZmin, fZmax, fZstep; - - /** Number of grid points **/ - Int_t fNx, fNy, fNz; // - - /** Arrays of doubles with the field values **/ - TArrayD* fBx; //! - TArrayD* fBy; //! - TArrayD* fBz; //! - - /** Variables for temporary storage - ** Used in the very frequently called method GetFieldValue **/ - Double_t fHa[2][2][2]; //! Field at corners of a grid cell - Double_t fHb[2][2]; //! Interpolated field (2-dim) - Double_t fHc[2]; //! Interpolated field (1-dim) - - // local transformation - TRotation* gRot; //! - TVector3* gTrans; //! - - //TTree with the map data when reading a ROOT file - TFile* fFile; // root file with the map data - - ClassDef(R3BGladFieldMap, 4) + public: + /** Default constructor **/ + R3BGladFieldMap(); + + /** Standard constructor + ** @param name Name of field map file without extension + ** @param fileType R = ROOT file, A = ASCII + **/ + R3BGladFieldMap(const TString& mapName, const TString& fileType = "A"); + + /** Constructor from R3BGladFieldPar **/ + R3BGladFieldMap(R3BFieldPar* fieldPar); + + /** Constructor from R3BGladFieldMapCreator **/ + // R3BGladFieldMap(R3BGladFieldMapCreator* creator); + + /** Destructor **/ + virtual ~R3BGladFieldMap(); + + /** Initialisation (read map from file) **/ + virtual void Init(); + + /** Get the field components at a certain point + ** @param x,y,z Point coordinates (global) [cm] + ** @value Bx,By,Bz Field components [kG] + **/ + virtual Double_t GetBx(Double_t x, Double_t y, Double_t z); + virtual Double_t GetBy(Double_t x, Double_t y, Double_t z); + virtual Double_t GetBz(Double_t x, Double_t y, Double_t z); + + // special function for the field transformation + TVector3 GetBtrans(Double_t x, Double_t y, Double_t z); + + /** Determine whether a point is inside the field map + ** @param x,y,z Point coordinates (global) [cm] + ** @param ix,iy,iz (return) Grid cell + ** @param dx,dy,dz (return) Distance from grid point [cm] if inside + ** @value kTRUE if inside map, else kFALSE + **/ + virtual Bool_t IsInside(Double_t x, + Double_t y, + Double_t z, + Int_t& ix, + Int_t& iy, + Int_t& iz, + Double_t& dx, + Double_t& dy, + Double_t& dz); + + /** Write the field map to an ASCII file **/ + void WriteAsciiFile(const TString& fileName); + + /** Set the position (in cm) of the field origin in lab **/ + virtual void SetPosition(Double_t x, Double_t y, Double_t z); + + /** Set a global field scaling factor **/ + virtual void SetScale(Double_t factor) { fScale = factor; } + virtual void SetFieldfromCurrent(double current); + + /* Set Euler rotation angles of the field (in degrees) + * default fYAngle = -14 deg, fXAngle=0, fZAngle=0 */ + virtual void SetXAngle(Double_t a) { fXAngle = a; }; + virtual void SetYAngle(Double_t a) { fYAngle = a; }; + virtual void SetZAngle(Double_t a) { fZAngle = a; }; + + /** Accessors to field parameters in local coordinate system **/ + Double_t GetXmin() const { return fXmin; } + Double_t GetYmin() const { return fYmin; } + Double_t GetZmin() const { return fZmin; } + Double_t GetXmax() const { return fXmax; } + Double_t GetYmax() const { return fYmax; } + Double_t GetZmax() const { return fZmax; } + Double_t GetXstep() const { return fXstep; } + Double_t GetYstep() const { return fYstep; } + Double_t GetZstep() const { return fZstep; } + Int_t GetNx() const { return fNx; } + Int_t GetNy() const { return fNy; } + Int_t GetNz() const { return fNz; } + + /** Accessor to field centre position in lab system **/ + Double_t GetPositionX() const { return fPosX; } + Double_t GetPositionY() const { return fPosY; } + Double_t GetPositionZ() const { return fPosZ; } + + /** Accessor to field rotation **/ + Double_t GetXAngle() const { return fXAngle; } + Double_t GetYAngle() const { return fYAngle; } + Double_t GetZAngle() const { return fZAngle; } + + /** Accessor to global scaling factor **/ + Double_t GetScale() const { return fScale; } + + /** Accessors to the field value arrays **/ + TArrayD* GetBx() const { return fBx; } + TArrayD* GetBy() const { return fBy; } + TArrayD* GetBz() const { return fBz; } + + /** Accessor to field map file **/ + TString GetFileName() { return fFileName; } + + /** Screen output **/ + virtual void Print(Option_t* option = "") const; + + protected: + /** Reset the field parameters and data **/ + void Reset(); + + /** Read the field map from an ASCII file **/ + void ReadAsciiFile(const TString& fileName); + + /** Read field map from a ROOT file **/ + void ReadRootFile(const TString& fileName); + + /** Get field values by interpolation of the grid. + ** @param dx,dy,dz Relative distance from grid point [cell units] + **/ + Double_t Interpolate(Double_t dx, Double_t dy, Double_t dz); + + /** Map file name **/ + TString fFileName; + + /** Global scaling factor (w.r.t. map on file) **/ + Double_t fScale; + + /** Field centre position in global coordinates **/ + Double_t fPosX, fPosY, fPosZ; + + /** Euler rotations of the field around field local X, Y, Z axis **/ + Double_t fXAngle; + Double_t fYAngle; + Double_t fZAngle; + + /** Field limits in local coordinate system **/ + Double_t fXmin, fXmax, fXstep; + Double_t fYmin, fYmax, fYstep; + Double_t fZmin, fZmax, fZstep; + + /** Number of grid points **/ + Int_t fNx, fNy, fNz; // + + /** Arrays of doubles with the field values **/ + TArrayD* fBx; //! + TArrayD* fBy; //! + TArrayD* fBz; //! + + /** Variables for temporary storage + ** Used in the very frequently called method GetFieldValue **/ + Double_t fHa[2][2][2]; //! Field at corners of a grid cell + Double_t fHb[2][2]; //! Interpolated field (2-dim) + Double_t fHc[2]; //! Interpolated field (1-dim) + + // local transformation + TRotation* gRot; //! + TVector3* gTrans; //! + + // TTree with the map data when reading a ROOT file + TFile* fFile; // root file with the map data + + ClassDef(R3BGladFieldMap, 5) }; #endif diff --git a/glad/R3BGladMagnet.cxx b/glad/R3BGladMagnet.cxx index d8aa986d7..69d13b559 100644 --- a/glad/R3BGladMagnet.cxx +++ b/glad/R3BGladMagnet.cxx @@ -18,7 +18,7 @@ // for the geometry creation (v17) and (v2023.1). // These will move also old files. const Double_t __GLAD_POS_DX = -42.0; // offset on the Z axis -const Double_t __GLAD_POS_DY = 2.0; // offset on the Y axis (2cm with respect to the beam line) +const Double_t __GLAD_POS_DY = 1.75; // offset on the Y axis (2cm with respect to the beam line) const Double_t __GLAD_POS_DZ = 308.8; // offset on the Z axis (distance from target) const Double_t __GLAD_ROT = 14; // rotation on the -Y axis const TString __GLAD_NAME = "Glad Magnet"; diff --git a/r3bbase/R3BLogger.h b/r3bbase/R3BLogger.h index faa7ee0ad..214d355e2 100644 --- a/r3bbase/R3BLogger.h +++ b/r3bbase/R3BLogger.h @@ -32,27 +32,26 @@ class R3BLogger; class R3BLogger : public FairLogger { public: -#define R3BLOG(severity, x) \ - if (true) \ - { \ - std::string fileName(__FILE__); \ - std::stringstream ss; \ - ss << fileName.substr(fileName.find_last_of("/") + 1) << ":" << __LINE__ << ":" << __FUNCTION__ << "(): "; \ - LOG(severity) << ss.str() << x; \ - } \ - else \ +#define R3BLOG(severity, x) \ + if (true) \ + { \ + std::string fN(__FILE__); \ + std::stringstream ss; \ + ss << fN.substr(fN.find_last_of("/") + 1) << ":" << __LINE__ << ":" << __FUNCTION__ << "(): "; \ + LOG(severity) << ss.str() << x; \ + } \ + else \ (void)0 -#define R3BLOG_IF(severity, condition, x) \ - if (true) \ - { \ - std::string fileNameif(__FILE__); \ - std::stringstream ssif; \ - ssif << fileNameif.substr(fileNameif.find_last_of("/") + 1) << ":" << __LINE__ << ":" << __FUNCTION__ \ - << "(): "; \ - LOG_IF(severity, condition) << ssif.str() << x; \ - } \ - else \ +#define R3BLOG_IF(severity, condition, x) \ + if (true) \ + { \ + std::string fNif(__FILE__); \ + std::stringstream ssif; \ + ssif << fNif.substr(fNif.find_last_of("/") + 1) << ":" << __LINE__ << ":" << __FUNCTION__ << "(): "; \ + LOG_IF(severity, condition) << ssif.str() << x; \ + } \ + else \ (void)0 private: diff --git a/r3bgen/R3BINCLRootGenerator.h b/r3bgen/R3BINCLRootGenerator.h index 566aef83c..e505dd643 100644 --- a/r3bgen/R3BINCLRootGenerator.h +++ b/r3bgen/R3BINCLRootGenerator.h @@ -63,17 +63,17 @@ class R3BINCLRootGenerator : public FairGenerator /** ** Method to simulate only fission events **/ - void SetOnlyFission(Bool_t Opt) { fOnlyFission = Opt; } + void SetOnlyFission(Bool_t Opt = true) { fOnlyFission = Opt; } /** ** Method to simulate only spallation events **/ - void SetOnlySpallation(Bool_t Opt) { fOnlySpallation = Opt; } + void SetOnlySpallation(Bool_t Opt = true) { fOnlySpallation = Opt; } /** ** Method to simulate only p2p-fission events **/ - void SetOnlyP2pFission(Bool_t Opt) { fOnlyP2pFission = Opt; } + void SetOnlyP2pFission(Bool_t Opt = true) { fOnlyP2pFission = Opt; } private: TString fFileName; // Input file name @@ -86,7 +86,7 @@ class R3BINCLRootGenerator : public FairGenerator ** any ion needed. TODO: Should not be needed by FairRoot. **/ void RegisterIons(); - inline Int_t GetIonPdgId(int z, int a) { return 1000000000 + 10 * 1000 * z + 10 * a; } + inline Int_t GetIonPdgId(int z, int a) { return 1000000000 + 10000 * z + 10 * a; } Int_t fEvt; TTree* Tree; diff --git a/twim/ana/R3BTwimVertexReconstruction.cxx b/twim/ana/R3BTwimVertexReconstruction.cxx index 99c6695df..978e066c0 100644 --- a/twim/ana/R3BTwimVertexReconstruction.cxx +++ b/twim/ana/R3BTwimVertexReconstruction.cxx @@ -59,23 +59,15 @@ InitStatus R3BTwimVertexReconstruction::Init() { R3BLOG(info, ""); FairRootManager* rootManager = FairRootManager::Instance(); - if (!rootManager) - { - R3BLOG(fatal, "FairRootManager not found"); - return kFATAL; - } + R3BLOG_IF(fatal, !rootManager, "FairRootManager not found"); header = dynamic_cast(rootManager->GetObject("EventHeader.")); - R3BLOG_IF(error, !header, "EventHeadder. not found"); + R3BLOG_IF(warn, !header, "EventHeadder. not found"); // INPUT DATA // get access to twim hit data fTwimHitDataCA = dynamic_cast(rootManager->GetObject("TwimHitData")); - if (!fTwimHitDataCA) - { - R3BLOG(fatal, "TwimHitData not found"); - return kFATAL; - } + R3BLOG_IF(fatal, !fTwimHitDataCA, "TwimHitData not found"); return kSUCCESS; }