Belle II Software light-2406-ragdoll
BFieldComponent3d Class Reference

The BFieldComponent3d class. More...

#include <BFieldComponent3d.h>

Inheritance diagram for BFieldComponent3d:
Collaboration diagram for BFieldComponent3d:

Public Member Functions

 BFieldComponent3d ()=default
 The BFieldComponent3d constructor.
 
virtual ~BFieldComponent3d ()=default
 The BFieldComponent3d destructor.
 
virtual void initialize () override
 Initializes the magnetic field component.
 
virtual ROOT::Math::XYZVector calculate (const ROOT::Math::XYZVector &point) const override
 Calculates the magnetic field vector at the specified space point.
 
virtual void terminate () override
 Terminates the magnetic field component.
 
void setMapFilename (const std::string &filename)
 Sets the filename of the magnetic field map.
 
void setMapSize (int sizeR, int sizePhi, int sizeZ)
 Sets the size of the magnetic field map.
 
void setMapRegionZ (double minZ, double maxZ, double offset)
 Sets the size of the magnetic field map.
 
void setMapRegionR (double minR, double maxR)
 Sets the size of the magnetic field map.
 
void setGridPitch (double pitchR, double pitchPhi, double pitchZ)
 Sets the grid pitch of the magnetic field map.
 
void setExcludeRegionZ (double minZ, double maxZ)
 Sets the size of the magnetic field map to exclude.
 
void setExcludeRegionR (double minR, double maxR)
 Sets the size of the magnetic field map to exclude.
 
void setErrorRegionR (double minR=-1., double maxR=-1., double errBr=0.0, double errBphi=0.0, double errBz=0.0)
 Sets the size of the magnetic field map to apply error on B.
 
void mirrorPhi (bool mirrorPhi=0.)
 
void doInterpolation (bool interpolate=true)
 
void enableCoordinate (const std::string &mapEnable="rphiz")
 Option to reduce 3D to 2D map (in coordinates, NOT Br, Bphi, Bz components)
 

Private Member Functions

ROOT::Math::XYZVector interpolate (unsigned int ir, unsigned int iphi, unsigned int iz, double wr, double wphi, double wz) const
 Interpolate the value of B-field between (ir, iphi, iz) and (ir+1, iphi+1, iz+1) using weights (wr, wphi, wz)
 

Private Attributes

std::string m_mapFilename {""}
 The filename of the magnetic field map.
 
std::vector< ROOT::Math::XYZVector > m_bmap
 The memory buffer for the magnetic field map.
 
std::string m_mapEnable {"rphiz"}
 Enable different dimension, "rphiz", "rphi", "phiz" or "rz" >
 
bool m_interpolate {true}
 Flag to switch on/off interpolation >
 
int m_mapSize [3]
 The size of the map in r, phi and z.
 
double m_mapRegionZ [2] {0}
 The min and max boundaries of the map region in z.
 
double m_mapOffset {0}
 Offset required because the accelerator group defines the Belle center as zero.
 
double m_mapRegionR [2] {0}
 The min and max boundaries of the map region in r.
 
double m_gridPitch [3] {0}
 The grid pitch in r,phi,z.
 
double m_igridPitch [3] {0}
 The inverted grid pitch in r,phi,z.
 
double m_exRegionZ [2] {0}
 The min and max boundaries of the excluded region in z.
 
double m_exRegionR [2] {0}
 The min and max boundaries of the excluded region in r.
 
bool m_exRegion {true}
 Flag to indicate whether there is a region to exclude.
 
bool m_mirrorPhi {true}
 Flag to indicate whether there is a region to exclude.
 
double m_errRegionR [2] {0}
 The min and max boundaries of the region in r to apply error.
 
double m_errB [3] {0}
 The error Br, Bphi, Bz as a scale factor (B_new = m_errB * B_old).
 

Detailed Description

The BFieldComponent3d class.

This class represents a 3d magnetic field map. The magnetic field map is stored as a grid in cylindrical coordinates. It is defined by a minimum radius and a maximum radius, a minimum z and a maximum z value (phi is assumed to be in 360 degrees range, option to mirror about the x-y plane exist), a pitch size in both, r, z and phi and the number of grid points. The ZOffset is used to account for the fact that the acceleration group defines 0 to be in the center of the detector, while the detector group defines the IP to be the center. The filename points to the zip file containing the magnetic field map.

Definition at line 36 of file BFieldComponent3d.h.

Member Function Documentation

◆ calculate()

ROOT::Math::XYZVector calculate ( const ROOT::Math::XYZVector &  point) const
overridevirtual

Calculates the magnetic field vector at the specified space point.

Parameters
pointThe space point in [cm] at which the magnetic field vector should be calculated.
Returns
The magnetic field vector at the given space point in [T]. Returns a zero vector XYZVector(0,0,0) if the space point lies outside the region described by the component.

Implements BFieldComponentAbs.

Definition at line 187 of file BFieldComponent3d.cc.

188{
189 auto getPhiIndexWeight = [this](double y, double x, double & wphi) -> unsigned int {
190 double phi = fast_atan2_minimax<4>(y, x);
191 wphi = phi * m_igridPitch[1];
192 auto iphi = static_cast<unsigned int>(wphi);
193 iphi = min(iphi, static_cast<unsigned int>(m_mapSize[1] - 2));
194 wphi -= iphi;
195 return iphi;
196 };
197
198 ROOT::Math::XYZVector B(0, 0, 0);
199
200 // If both '3d' and 'Beamline' components are defined in xml file,
201 // '3d' component returns zero field where 'Beamline' component is defined.
202 // If no 'Beamline' component is defined in xml file, the following function will never be called.
203 if (BFieldComponentBeamline::Instance().isInRange(point)) {
204 return B;
205 }
206
207 double z = point.Z();
208 // Check if the point lies inside the magnetic field boundaries
209 if (z < m_mapRegionZ[0] || z > m_mapRegionZ[1]) return B;
210
211 double r2 = point.Perp2();
212 // Check if the point lies inside the magnetic field boundaries
213 if (r2 < m_mapRegionR[0]*m_mapRegionR[0] || r2 >= m_mapRegionR[1]*m_mapRegionR[1]) return B;
214 // Check if the point lies in the exclude region
215 if (m_exRegion && (z >= m_exRegionZ[0]) && (z < m_exRegionZ[1]) &&
216 (r2 >= m_exRegionR[0]*m_exRegionR[0]) && (r2 < m_exRegionR[1]*m_exRegionR[1])) return B;
217
218 // Calculate the lower index of the point in the Z grid
219 // Note that z coordinate is inverted to match ANSYS frame
220 double wz = (m_mapRegionZ[1] - z) * m_igridPitch[2];
221 unsigned int iz = static_cast<int>(wz);
222 iz = min(iz, static_cast<unsigned int>(m_mapSize[2] - 2));
223 wz -= iz;
224
225 if (r2 > 0) {
226 double r = sqrt(r2);
227
228 // Calculate the lower index of the point in the R grid
229 double wr = (r - m_mapRegionR[0]) * m_igridPitch[0];
230 auto ir = static_cast<unsigned int>(wr);
231 ir = min(ir, static_cast<unsigned int>(m_mapSize[0] - 2));
232 wr -= ir;
233
234 // Calculate the lower index of the point in the Phi grid
235 double ay = std::abs(point.Y());
236 double wphi;
237 unsigned int iphi = getPhiIndexWeight(ay, point.X(), wphi);
238
239 // Get B-field values from map
240 ROOT::Math::XYZVector b = interpolate(ir, iphi, iz, wr, wphi, wz); // in cylindrical system
241 double norm = 1 / r;
242 double s = ay * norm, c = point.X() * norm;
243 // Flip sign of By if y<0
244 const double sgny = (point.Y() >= 0) - (point.Y() < 0);
245 // in cartesian system
246 B.SetXYZ(-(b.X() * c - b.Y() * s), -sgny * (b.X() * s + b.Y() * c), b.Z());
247 } else {
248 // Get B-field values from map in cartesian system assuming phi=0, so Bx = Br and By = Bphi
249 B = interpolate(0, 0, iz, 0, 0, wz);
250 }
251
252 return B;
253}
double m_exRegionR[2]
The min and max boundaries of the excluded region in r.
ROOT::Math::XYZVector interpolate(unsigned int ir, unsigned int iphi, unsigned int iz, double wr, double wphi, double wz) const
Interpolate the value of B-field between (ir, iphi, iz) and (ir+1, iphi+1, iz+1) using weights (wr,...
double m_mapRegionZ[2]
The min and max boundaries of the map region in z.
int m_mapSize[3]
The size of the map in r, phi and z.
double m_igridPitch[3]
The inverted grid pitch in r,phi,z.
bool m_exRegion
Flag to indicate whether there is a region to exclude.
double m_mapRegionR[2]
The min and max boundaries of the map region in r.
double m_exRegionZ[2]
The min and max boundaries of the excluded region in z.
static BFieldComponentBeamline & Instance()
BFieldComponentBeamline instance.

◆ doInterpolation()

void doInterpolation ( bool  interpolate = true)
inline
Parameters
interpolateFlag to switch on/off interpolation

Definition at line 137 of file BFieldComponent3d.h.

bool m_interpolate
Flag to switch on/off interpolation >

◆ enableCoordinate()

void enableCoordinate ( const std::string &  mapEnable = "rphiz")
inline

Option to reduce 3D to 2D map (in coordinates, NOT Br, Bphi, Bz components)

Parameters
mapEnableList of dimensions to enable: "rphiz", "rphi", "phiz" or "rz"

Definition at line 143 of file BFieldComponent3d.h.

143{ m_mapEnable = mapEnable; }
std::string m_mapEnable
Enable different dimension, "rphiz", "rphi", "phiz" or "rz" >

◆ initialize()

void initialize ( )
overridevirtual

Initializes the magnetic field component.

This method opens the magnetic field map file.

Reimplemented from BFieldComponentAbs.

Definition at line 27 of file BFieldComponent3d.cc.

28{
29
30 // Input field map
31 if (m_mapFilename.empty()) {
32 B2ERROR("The filename for the 3d magnetic field component is empty !");
33 return;
34 }
35 string fullPath = FileSystem::findFile("/data/" + m_mapFilename);
36 if (!FileSystem::fileExists(fullPath)) {
37 B2ERROR("The 3d magnetic field map file '" << m_mapFilename << "' could not be found !");
38 return;
39 }
40
41 // Options to reduce 3D to 2D
42 if (m_mapEnable != "rphiz" && m_mapEnable != "rphi" && m_mapEnable != "phiz" && m_mapEnable != "rz") {
43 B2ERROR("BField3d:: enabled coordinates must be \"rphiz\", \"rphi\", \"phiz\" or \"rz\"");
44 return;
45 }
46
47 // Excluded region
48 m_exRegion = true;
49 if ((m_exRegionR[0] == m_exRegionR[1]) || (m_exRegionZ[0] == m_exRegionZ[1])) m_exRegion = false;
50
51 // Print initial input parameters
52 B2DEBUG(100, "BField3d:: initial input parameters");
53 B2DEBUG(100, Form(" map filename: %s", m_mapFilename.c_str()));
54 B2DEBUG(100, Form(" map dimension: %s", m_mapEnable.c_str()));
55 if (m_interpolate) { B2DEBUG(100, Form(" map interpolation: on")); }
56 else { B2DEBUG(100, Form(" map interpolation: off")); }
57 B2DEBUG(100, Form(" map r pitch & range: %.2e [%.2e, %.2e] cm", m_gridPitch[0], m_mapRegionR[0], m_mapRegionR[1]));
58 B2DEBUG(100, Form(" map phi pitch: %.2e deg", m_gridPitch[1] * 180 / M_PI));
59 B2DEBUG(100, Form(" map z pitch & range: %.2e [%.2e, %.2e] cm", m_gridPitch[2], m_mapRegionZ[0], m_mapRegionZ[1]));
60 if (m_exRegion) {
61 B2DEBUG(100, Form(" map r excluded region: [%.2e, %.2e] cm", m_exRegionR[0], m_exRegionR[1]));
62 B2DEBUG(100, Form(" map z excluded region: [%.2e, %.2e] cm", m_exRegionZ[0], m_exRegionZ[1]));
63 }
64
65 m_bmap.reserve(m_mapSize[0]*m_mapSize[1]*m_mapSize[2]);
66 // Load B-field map file
67 io::filtering_istream fieldMapFile;
68 fieldMapFile.push(io::gzip_decompressor());
69 fieldMapFile.push(io::file_source(fullPath));
70
71 char tmp[256];
72 for (int k = 0; k < m_mapSize[2]; k++) { // z
73 for (int i = 0; i < m_mapSize[0]; i++) { // r
74 for (int j = 0; j < m_mapSize[1]; j++) { // phi
75 double Br, Bz, Bphi;
76 //r[m] phi[deg] z[m] Br[T] Bphi[T] Bz[T]
77 fieldMapFile.getline(tmp, 256);
78 // sscanf(tmp+33,"%lf %lf %lf",&Br,&Bphi,&Bz);
79 char* next;
80 Br = strtod(tmp + 33, &next);
81 Bphi = strtod(next, &next);
82 Bz = strtod(next, nullptr);
83 m_bmap.emplace_back(-Br, -Bphi, -Bz);
84 }
85 }
86 }
87
88 // Introduce error on B field
89 if (m_errRegionR[0] != m_errRegionR[1]) {
90 auto it = m_bmap.begin();
91 for (int k = 0; k < m_mapSize[2]; k++) { // z
92 double r = m_mapRegionR[0];
93 for (int i = 0; i < m_mapSize[0]; i++, r += m_gridPitch[0]) { // r
94 if (!(r >= m_errRegionR[0] && r < m_errRegionR[1])) { it += m_mapSize[1]; continue;}
95 for (int j = 0; j < m_mapSize[1]; j++) { // phi
96 ROOT::Math::XYZVector& B = *it;
97 B.SetX(B.X() * m_errB[0]);
98 B.SetY(B.Y() * m_errB[1]);
99 B.SetZ(B.Z() * m_errB[2]);
100 }
101 }
102 }
103 }
104
105 m_igridPitch[0] = 1 / m_gridPitch[0];
106 m_igridPitch[1] = 1 / m_gridPitch[1];
107 m_igridPitch[2] = 1 / m_gridPitch[2];
108
109 B2DEBUG(100, Form("BField3d:: final map region & pitch: r [%.2e,%.2e] %.2e, phi %.2e, z [%.2e,%.2e] %.2e",
112 B2DEBUG(100, "Memory consumption: " << m_bmap.size()*sizeof(ROOT::Math::XYZVector) / (1024 * 1024.) << " Mb");
113}
double m_errRegionR[2]
The min and max boundaries of the region in r to apply error.
double m_gridPitch[3]
The grid pitch in r,phi,z.
std::string m_mapFilename
The filename of the magnetic field map.
std::vector< ROOT::Math::XYZVector > m_bmap
The memory buffer for the magnetic field map.
double m_errB[3]
The error Br, Bphi, Bz as a scale factor (B_new = m_errB * B_old).
static std::string findFile(const std::string &path, bool silent=false)
Search for given file or directory in local or central release directory, and return absolute path if...
Definition: FileSystem.cc:151
static bool fileExists(const std::string &filename)
Check if the file with given filename exists.
Definition: FileSystem.cc:32

◆ interpolate()

ROOT::Math::XYZVector interpolate ( unsigned int  ir,
unsigned int  iphi,
unsigned int  iz,
double  wr,
double  wphi,
double  wz 
) const
private

Interpolate the value of B-field between (ir, iphi, iz) and (ir+1, iphi+1, iz+1) using weights (wr, wphi, wz)

Definition at line 259 of file BFieldComponent3d.cc.

261{
262 const unsigned int strideZ = m_mapSize[0] * m_mapSize[1];
263 const unsigned int strideR = m_mapSize[1];
264
265 const double wz0 = 1 - wz1, wr0 = 1 - wr1, wphi0 = 1 - wphi1;
266 const unsigned int j000 = iz * strideZ + ir * strideR + iphi;
267 const unsigned int j001 = j000 + 1;
268 const unsigned int j010 = j000 + strideR;
269 const unsigned int j011 = j001 + strideR;
270 const unsigned int j100 = j000 + strideZ;
271 const unsigned int j101 = j001 + strideZ;
272 const unsigned int j110 = j010 + strideZ;
273 const unsigned int j111 = j011 + strideZ;
274 const double w00 = wphi0 * wr0;
275 const double w10 = wphi0 * wr1;
276 const double w01 = wphi1 * wr0;
277 const double w11 = wphi1 * wr1;
278 const vector<ROOT::Math::XYZVector>& B = m_bmap;
279 return
280 (B[j000] * w00 + B[j001] * w01 + B[j010] * w10 + B[j011] * w11) * wz0 +
281 (B[j100] * w00 + B[j101] * w01 + B[j110] * w10 + B[j111] * w11) * wz1;
282}

◆ mirrorPhi()

void mirrorPhi ( bool  mirrorPhi = 0.)
inline
Parameters
mirrorPhiFlag to enable mirroring in phi about x-z plane

Definition at line 132 of file BFieldComponent3d.h.

bool m_mirrorPhi
Flag to indicate whether there is a region to exclude.
void mirrorPhi(bool mirrorPhi=0.)

◆ setErrorRegionR()

void setErrorRegionR ( double  minR = -1.,
double  maxR = -1.,
double  errBr = 0.0,
double  errBphi = 0.0,
double  errBz = 0.0 
)
inline

Sets the size of the magnetic field map to apply error on B.

Parameters
minRThe left (min) border of the magnetic field map region in r [m] -> [cm].
maxRThe right (max) border of the magnetic field map region in r [m] -> [cm].
errBrThe size of the error on Br, Br_new = errBr*Br [fraction].
errBphiThe size of the error on Bphi, Bphi_new = errBphi*Bphi [fraction].
errBzThe size of the error on Bz, Bz_new = errBz*Bz [fraction].

Definition at line 126 of file BFieldComponent3d.h.

127 { m_errRegionR[0] = minR; m_errRegionR[1] = maxR; m_errB[0] = errBr; m_errB[1] = errBphi; m_errB[2] = errBz; }

◆ setExcludeRegionR()

void setExcludeRegionR ( double  minR,
double  maxR 
)
inline

Sets the size of the magnetic field map to exclude.

Parameters
minRThe left (min) border of the magnetic field map region in r [m] -> [cm].
maxRThe right (max) border of the magnetic field map region in r [m] -> [cm].

Definition at line 116 of file BFieldComponent3d.h.

116{ m_exRegionR[0] = minR; m_exRegionR[1] = maxR; }

◆ setExcludeRegionZ()

void setExcludeRegionZ ( double  minZ,
double  maxZ 
)
inline

Sets the size of the magnetic field map to exclude.

Parameters
minZThe left (min) border of the magnetic field map region in z [m] -> [cm].
maxZThe right (max) border of the magnetic field map region in z [m] -> [cm].

Definition at line 109 of file BFieldComponent3d.h.

109{ m_exRegionZ[0] = minZ + m_mapOffset; m_exRegionZ[1] = maxZ + m_mapOffset; }
double m_mapOffset
Offset required because the accelerator group defines the Belle center as zero.

◆ setGridPitch()

void setGridPitch ( double  pitchR,
double  pitchPhi,
double  pitchZ 
)
inline

Sets the grid pitch of the magnetic field map.

Parameters
pitchRThe pitch of the grid in r [m] -> [cm].
pitchPhiThe pitch of the grid in phi [degrees] -> [rad] converted back to degrees in the code.
pitchZThe pitch of the grid in z [m] -> [cm].

Definition at line 102 of file BFieldComponent3d.h.

102{ m_gridPitch[0] = pitchR; m_gridPitch[1] = pitchPhi; m_gridPitch[2] = pitchZ; }

◆ setMapFilename()

void setMapFilename ( const std::string &  filename)
inline

Sets the filename of the magnetic field map.

Parameters
filenameThe filename of the magnetic field map.

Definition at line 71 of file BFieldComponent3d.h.

71{ m_mapFilename = filename; };

◆ setMapRegionR()

void setMapRegionR ( double  minR,
double  maxR 
)
inline

Sets the size of the magnetic field map.

Parameters
minRThe left (min) border of the magnetic field map region in r [m] -> [cm].
maxRThe right (max) border of the magnetic field map region in r [m] -> [cm].

Definition at line 94 of file BFieldComponent3d.h.

94{ m_mapRegionR[0] = minR; m_mapRegionR[1] = maxR; }

◆ setMapRegionZ()

void setMapRegionZ ( double  minZ,
double  maxZ,
double  offset 
)
inline

Sets the size of the magnetic field map.

Parameters
minZThe left (min) border of the magnetic field map region in z [m] -> [cm].
maxZThe right (max) border of the magnetic field map region in z [m] -> [cm].
offsetThe offset in z [m] -> [cm] which is required because the accelerator group defines the Belle center as zero.

Definition at line 87 of file BFieldComponent3d.h.

87{ m_mapRegionZ[0] = minZ; m_mapRegionZ[1] = maxZ; m_mapOffset = offset; }

◆ setMapSize()

void setMapSize ( int  sizeR,
int  sizePhi,
int  sizeZ 
)
inline

Sets the size of the magnetic field map.

Parameters
sizeRThe number of points in the r direction.
sizePhiThe number of points in the phi direction.
sizeZThe number of points in the z direction.

Definition at line 79 of file BFieldComponent3d.h.

79{ m_mapSize[0] = sizeR; m_mapSize[1] = sizePhi; m_mapSize[2] = sizeZ; }

◆ terminate()

void terminate ( )
overridevirtual

Terminates the magnetic field component.

This method closes the magnetic field map file.

Reimplemented from BFieldComponentAbs.

Definition at line 255 of file BFieldComponent3d.cc.

256{
257}

Member Data Documentation

◆ m_bmap

std::vector<ROOT::Math::XYZVector> m_bmap
private

The memory buffer for the magnetic field map.

Definition at line 156 of file BFieldComponent3d.h.

◆ m_errB

double m_errB[3] {0}
private

The error Br, Bphi, Bz as a scale factor (B_new = m_errB * B_old).

Definition at line 184 of file BFieldComponent3d.h.

◆ m_errRegionR

double m_errRegionR[2] {0}
private

The min and max boundaries of the region in r to apply error.

Definition at line 182 of file BFieldComponent3d.h.

◆ m_exRegion

bool m_exRegion {true}
private

Flag to indicate whether there is a region to exclude.

>

Definition at line 178 of file BFieldComponent3d.h.

◆ m_exRegionR

double m_exRegionR[2] {0}
private

The min and max boundaries of the excluded region in r.

Definition at line 176 of file BFieldComponent3d.h.

◆ m_exRegionZ

double m_exRegionZ[2] {0}
private

The min and max boundaries of the excluded region in z.

Definition at line 174 of file BFieldComponent3d.h.

◆ m_gridPitch

double m_gridPitch[3] {0}
private

The grid pitch in r,phi,z.

Definition at line 170 of file BFieldComponent3d.h.

◆ m_igridPitch

double m_igridPitch[3] {0}
private

The inverted grid pitch in r,phi,z.

Definition at line 172 of file BFieldComponent3d.h.

◆ m_interpolate

bool m_interpolate {true}
private

Flag to switch on/off interpolation >

Definition at line 160 of file BFieldComponent3d.h.

◆ m_mapEnable

std::string m_mapEnable {"rphiz"}
private

Enable different dimension, "rphiz", "rphi", "phiz" or "rz" >

Definition at line 158 of file BFieldComponent3d.h.

◆ m_mapFilename

std::string m_mapFilename {""}
private

The filename of the magnetic field map.

Definition at line 154 of file BFieldComponent3d.h.

◆ m_mapOffset

double m_mapOffset {0}
private

Offset required because the accelerator group defines the Belle center as zero.

Definition at line 166 of file BFieldComponent3d.h.

◆ m_mapRegionR

double m_mapRegionR[2] {0}
private

The min and max boundaries of the map region in r.

Definition at line 168 of file BFieldComponent3d.h.

◆ m_mapRegionZ

double m_mapRegionZ[2] {0}
private

The min and max boundaries of the map region in z.

Definition at line 164 of file BFieldComponent3d.h.

◆ m_mapSize

int m_mapSize[3]
private

The size of the map in r, phi and z.

Definition at line 162 of file BFieldComponent3d.h.

◆ m_mirrorPhi

bool m_mirrorPhi {true}
private

Flag to indicate whether there is a region to exclude.

>

Definition at line 180 of file BFieldComponent3d.h.


The documentation for this class was generated from the following files: