LDMX Software
simcore::MagneticFieldMap3D Class Reference

A 3D B-field map defined as a grid of points with associated B-field values. More...

#include <MagneticFieldMap3D.h>

Public Member Functions

 MagneticFieldMap3D (const char *filename, double xOffset, double yOffset, double zOffset)
 Class constructor.
 
void GetFieldValue (const double point[4], double *bfield) const
 Implementation of primary virtual method from G4MagneticField interface.
 

Private Attributes

vector< vector< vector< double > > > x_field_
 
vector< vector< vector< double > > > y_field_
 
vector< vector< vector< double > > > z_field_
 
int nx_
 
int ny_
 
int nz_
 
double minx_
 
double maxx_
 
double miny_
 
double maxy_
 
double minz_
 
double maxz_
 
double dx_
 
double dy_
 
double dz_
 
double x_offset_
 
double y_offset_
 
double z_offset_
 
bool invert_x_
 
bool invert_y_
 
bool invert_z_
 

Detailed Description

A 3D B-field map defined as a grid of points with associated B-field values.

Note
Values are interpolated to obtain B-field information at points in between.

The header format of the input data file is as follows:

1) Blank line

2) int int int

Number of X, Y, and Z grid points (N_x, N_y, N_z)

3-8) int String

A description of fields in the lines of the map. In our case it is literally:

1 X 2 Y 3 Z 4 BX 5 BY 6 BZ

9) int String

Header terminator

10+) float float float float float float

N_x*N_y*N_z+9

x_ y_ z_ B_x B_y B_z

Original PurgMagTabulatedField3D code developed by: S.Larsson and J. Generowicz.

Definition at line 63 of file MagneticFieldMap3D.h.

Constructor & Destructor Documentation

◆ MagneticFieldMap3D()

simcore::MagneticFieldMap3D::MagneticFieldMap3D ( const char * filename,
double xOffset,
double yOffset,
double zOffset )

Class constructor.

Parameters
[in]filenameThe name of the file defining the B-field grid.
[in]xOffset,yOffset,zOffsetThe offset of the grid's coordinate system.

Definition at line 18 of file MagneticFieldMap3D.cxx.

20 : nx_(0),
21 ny_(0),
22 nz_(0),
23 x_offset_(xOffset),
24 y_offset_(yOffset),
25 z_offset_(zOffset),
26 invert_x_(false),
27 invert_y_(false),
28 invert_z_(false) {
29 ifstream file(filename); // Open the file for reading.
30
31 // Throw an error if file does not exist.
32 if (!file.good()) {
33 EXCEPTION_RAISE("FileDNE", "The field map file '" + std::string(filename) +
34 "' does not exist!");
35 }
36
37 ldmx_log(trace)
38 << "-----------------------------------------------------------";
39 ldmx_log(trace) << " Magnetic Field Map 3D";
40 ldmx_log(trace)
41 << "-----------------------------------------------------------";
42
43 ldmx_log(trace) << "Reading the field grid from " << filename << " ... ";
44 ldmx_log(trace) << " Offsets: " << xOffset << " " << yOffset << " "
45 << zOffset;
46
47 // Ignore first blank line
48 char buffer[256];
49 file.getline(buffer, 256);
50
51 // Read table dimensions
52 file >> nx_ >> ny_ >> nz_; // Note dodgy order
53
54 ldmx_log(trace) << " Number of values: " << nx_ << " " << ny_ << " " << nz_;
55
56 // Set up storage space for table
57 x_field_.resize(nx_);
58 y_field_.resize(nx_);
59 z_field_.resize(nx_);
60 int ix, iy, iz;
61 for (ix = 0; ix < nx_; ix++) {
62 x_field_[ix].resize(ny_);
63 y_field_[ix].resize(ny_);
64 z_field_[ix].resize(ny_);
65 for (iy = 0; iy < ny_; iy++) {
66 x_field_[ix][iy].resize(nz_);
67 y_field_[ix][iy].resize(nz_);
68 z_field_[ix][iy].resize(nz_);
69 }
70 }
71
72 // Ignore other header information
73 // The first line whose second character is '0' is considered to
74 // be the last line of the header.
75 do {
76 file.getline(buffer, 256);
77 } while (buffer[1] != '0');
78
79 // Read in the data
80 double xval, yval, zval, bx, by, bz;
81 for (ix = 0; ix < nx_; ix++) {
82 for (iy = 0; iy < ny_; iy++) {
83 for (iz = 0; iz < nz_; iz++) {
84 file >> xval >> yval >> zval >> bx >> by >> bz;
85 if (ix == 0 && iy == 0 && iz == 0) {
86 minx_ = xval;
87 miny_ = yval;
88 minz_ = zval;
89 }
90 x_field_[ix][iy][iz] = bx;
91 y_field_[ix][iy][iz] = by;
92 z_field_[ix][iy][iz] = bz;
93 }
94 }
95 }
96 file.close();
97
98 maxx_ = xval;
99 maxy_ = yval;
100 maxz_ = zval;
101
102 ldmx_log(trace) << " ... done reading ";
103 ldmx_log(trace) << "Read values of field from file " << filename;
104 ldmx_log(trace) << " Assumed the order: x_, y_, z_, Bx, By, Bz";
105 ldmx_log(trace) << " Min values: " << minx_ << " " << miny_ << " " << minz_
106 << " mm ";
107 ldmx_log(trace) << " Max values: " << maxx_ << " " << maxy_ << " " << maxz_
108 << " mm ";
109 ldmx_log(trace) << " Field offsets: " << x_offset_ << " " << y_offset_ << " "
110 << z_offset_ << " mm ";
111
112 // Should really check that the limits are not the wrong way around.
113 if (maxx_ < minx_) {
114 swap(maxx_, minx_);
115 invert_x_ = true;
116 }
117 if (maxy_ < miny_) {
118 swap(maxy_, miny_);
119 invert_y_ = true;
120 }
121 if (maxz_ < minz_) {
122 swap(maxz_, minz_);
123 invert_z_ = true;
124 }
125
126 ldmx_log(trace) << "After reordering if necessary";
127 ldmx_log(trace) << " Min values: " << minx_ << " " << miny_ << " " << minz_
128 << " mm ";
129 ldmx_log(trace) << " Max values: " << maxx_ << " " << maxy_ << " " << maxz_
130 << " mm ";
131 ;
132
133 dx_ = maxx_ - minx_;
134 dy_ = maxy_ - miny_;
135 dz_ = maxz_ - minz_;
136
137 ldmx_log(trace) << " Range of values: " << dx_ << " " << dy_ << " " << dz_
138 << " mm";
139 ldmx_log(trace) << "Done loading field map";
140 ldmx_log(trace)
141 << "-----------------------------------------------------------";
142}

Member Function Documentation

◆ GetFieldValue()

void simcore::MagneticFieldMap3D::GetFieldValue ( const double point[4],
double * bfield ) const

Implementation of primary virtual method from G4MagneticField interface.

Parameters
[in]pointThe point in 3D space.
[out]bfieldThe output B-field data at the point.

Definition at line 144 of file MagneticFieldMap3D.cxx.

145 {
146 double x = point[0] - x_offset_;
147 double y = point[1] - y_offset_;
148 double z = point[2] - z_offset_;
149 double eps = 1E-6;
150
151 // Check that the point is within the defined region
152 if ((x >= minx_ && (x < maxx_ - eps)) && (y >= miny_ && (y < maxy_ - eps)) &&
153 (z >= minz_ && (z < maxz_ - eps))) {
154 // Position of given point within region, normalized to the range
155 // [0,1]
156 double xfraction = (x - minx_) / dx_;
157 double yfraction = (y - miny_) / dy_;
158 double zfraction = (z - minz_) / dz_;
159
160 if (invert_x_) {
161 xfraction = 1 - xfraction;
162 }
163 if (invert_y_) {
164 yfraction = 1 - yfraction;
165 }
166 if (invert_z_) {
167 zfraction = 1 - zfraction;
168 }
169
170 // Need addresses of these to pass to modf below.
171 // modf uses its second argument as an OUTPUT argument.
172 double xdindex, ydindex, zdindex;
173
174 // Position of the point within the cuboid defined by the
175 // nearest surrounding tabulated points
176 double xlocal = (std::modf(xfraction * (nx_ - 1), &xdindex));
177 double ylocal = (std::modf(yfraction * (ny_ - 1), &ydindex));
178 double zlocal = (std::modf(zfraction * (nz_ - 1), &zdindex));
179
180 // The indices of the nearest tabulated point whose coordinates
181 // are all less than those of the given point
182 int xindex = static_cast<int>(xdindex);
183 int yindex = static_cast<int>(ydindex);
184 int zindex = static_cast<int>(zdindex);
185
186#ifdef DEBUG_INTERPOLATING_FIELD
187 ldmx_log(trace) << "Local x_,y_,z_: " << xlocal << " " << ylocal << " "
188 << zlocal;
189 ldmx_log(trace) << "Index x_,y_,z_: " << xindex << " " << yindex << " "
190 << zindex;
191 double valx0z0, mulx0z0, valx1z0, mulx1z0;
192 double valx0z1, mulx0z1, valx1z1, mulx1z1;
193 valx0z0 = table[xindex][0][zindex];
194 mulx0z0 = (1 - xlocal) * (1 - zlocal);
195 valx1z0 = table[xindex + 1][0][zindex];
196 mulx1z0 = xlocal * (1 - zlocal);
197 valx0z1 = table[xindex][0][zindex + 1];
198 mulx0z1 = (1 - xlocal) * zlocal;
199 valx1z1 = table[xindex + 1][0][zindex + 1];
200 mulx1z1 = xlocal * zlocal;
201#endif
202
203 // Full 3-dimensional version
204 bfield[0] =
205 x_field_[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) *
206 (1 - zlocal) +
207 x_field_[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) *
208 zlocal +
209 x_field_[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal *
210 (1 - zlocal) +
211 x_field_[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal *
212 zlocal +
213 x_field_[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) *
214 (1 - zlocal) +
215 x_field_[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) *
216 zlocal +
217 x_field_[xindex + 1][yindex + 1][zindex] * xlocal * ylocal *
218 (1 - zlocal) +
219 x_field_[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;
220 bfield[1] =
221 y_field_[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) *
222 (1 - zlocal) +
223 y_field_[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) *
224 zlocal +
225 y_field_[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal *
226 (1 - zlocal) +
227 y_field_[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal *
228 zlocal +
229 y_field_[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) *
230 (1 - zlocal) +
231 y_field_[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) *
232 zlocal +
233 y_field_[xindex + 1][yindex + 1][zindex] * xlocal * ylocal *
234 (1 - zlocal) +
235 y_field_[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;
236 bfield[2] =
237 z_field_[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) *
238 (1 - zlocal) +
239 z_field_[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) *
240 zlocal +
241 z_field_[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal *
242 (1 - zlocal) +
243 z_field_[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal *
244 zlocal +
245 z_field_[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) *
246 (1 - zlocal) +
247 z_field_[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) *
248 zlocal +
249 z_field_[xindex + 1][yindex + 1][zindex] * xlocal * ylocal *
250 (1 - zlocal) +
251 z_field_[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;
252
253 } else {
254 bfield[0] = 0.0;
255 bfield[1] = 0.0;
256 bfield[2] = 0.0;
257 }
258}

Member Data Documentation

◆ dx_

double simcore::MagneticFieldMap3D::dx_
private

Definition at line 102 of file MagneticFieldMap3D.h.

◆ dy_

double simcore::MagneticFieldMap3D::dy_
private

Definition at line 102 of file MagneticFieldMap3D.h.

◆ dz_

double simcore::MagneticFieldMap3D::dz_
private

Definition at line 102 of file MagneticFieldMap3D.h.

◆ invert_x_

bool simcore::MagneticFieldMap3D::invert_x_
private

Definition at line 114 of file MagneticFieldMap3D.h.

◆ invert_y_

bool simcore::MagneticFieldMap3D::invert_y_
private

Definition at line 114 of file MagneticFieldMap3D.h.

◆ invert_z_

bool simcore::MagneticFieldMap3D::invert_z_
private

Definition at line 114 of file MagneticFieldMap3D.h.

◆ maxx_

double simcore::MagneticFieldMap3D::maxx_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ maxy_

double simcore::MagneticFieldMap3D::maxy_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ maxz_

double simcore::MagneticFieldMap3D::maxz_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ minx_

double simcore::MagneticFieldMap3D::minx_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ miny_

double simcore::MagneticFieldMap3D::miny_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ minz_

double simcore::MagneticFieldMap3D::minz_
private

Definition at line 97 of file MagneticFieldMap3D.h.

◆ nx_

int simcore::MagneticFieldMap3D::nx_
private

Definition at line 92 of file MagneticFieldMap3D.h.

◆ ny_

int simcore::MagneticFieldMap3D::ny_
private

Definition at line 92 of file MagneticFieldMap3D.h.

◆ nz_

int simcore::MagneticFieldMap3D::nz_
private

Definition at line 92 of file MagneticFieldMap3D.h.

◆ x_field_

vector<vector<vector<double> > > simcore::MagneticFieldMap3D::x_field_
private

Definition at line 85 of file MagneticFieldMap3D.h.

◆ x_offset_

double simcore::MagneticFieldMap3D::x_offset_
private

Definition at line 107 of file MagneticFieldMap3D.h.

◆ y_field_

vector<vector<vector<double> > > simcore::MagneticFieldMap3D::y_field_
private

Definition at line 86 of file MagneticFieldMap3D.h.

◆ y_offset_

double simcore::MagneticFieldMap3D::y_offset_
private

Definition at line 108 of file MagneticFieldMap3D.h.

◆ z_field_

vector<vector<vector<double> > > simcore::MagneticFieldMap3D::z_field_
private

Definition at line 87 of file MagneticFieldMap3D.h.

◆ z_offset_

double simcore::MagneticFieldMap3D::z_offset_
private

Definition at line 109 of file MagneticFieldMap3D.h.


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