3 #include <Geant4/G4SystemOfUnits.hh>
17 , m_MagFieldScale(magfield_rescale)
19 ifstream
file(filename);
22 cout <<
"Could not open " << filename <<
" exiting now" << endl;
34 for (ix = 0; ix <
nx; ix++)
39 for (iy = 0; iy <
ny; iy++)
52 file.getline(buffer, 256);
53 }
while (buffer[1] !=
'0');
56 double xval, yval, zval, bx, by, bz;
59 for (ix = 0; ix <
nx; ix++)
61 for (iy = 0; iy <
ny; iy++)
63 for (iz = 0; iz <
nz; iz++)
65 file >> xval >> yval >> zval >> bx >> by >> bz;
74 cout <<
"mismatch expected x coordinate: " << save_x
75 <<
" read " << xval << endl;
87 cout <<
"mismatch expected y coordinate: " << save_y
88 <<
" read " << yval << endl;
92 if (ix == 0 && iy == 0 && iz == 0)
117 double x = fabs(point[0]);
118 double y = fabs(point[1]);
122 if (x >=
minx && x <= maxx && y >=
miny && y <= maxy && z >=
minz && z <=
maxz)
126 double xfraction = (x -
minx) /
dx;
127 double yfraction = (y -
miny) /
dy;
128 double zfraction = (z -
minz) /
dz;
132 double xdindex, ydindex, zdindex;
136 double xlocal = (std::modf(xfraction * (
nx - 1), &xdindex));
137 double ylocal = (std::modf(yfraction * (
ny - 1), &ydindex));
138 double zlocal = (std::modf(zfraction * (
nz - 1), &zdindex));
142 int xindex =
static_cast<int>(xdindex);
143 int yindex =
static_cast<int>(ydindex);
144 int zindex =
static_cast<int>(zdindex);
147 Bfield[0] =
xField[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) * (1 - zlocal) +
xField[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) * zlocal +
xField[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal * (1 - zlocal) +
xField[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal * zlocal +
xField[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) * (1 - zlocal) +
xField[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) * zlocal +
xField[xindex + 1][yindex + 1][zindex] * xlocal * ylocal * (1 - zlocal) +
xField[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;
148 Bfield[1] =
yField[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) * (1 - zlocal) +
yField[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) * zlocal +
yField[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal * (1 - zlocal) +
yField[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal * zlocal +
yField[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) * (1 - zlocal) +
yField[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) * zlocal +
yField[xindex + 1][yindex + 1][zindex] * xlocal * ylocal * (1 - zlocal) +
yField[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;
149 Bfield[2] =
zField[xindex][yindex][zindex] * (1 - xlocal) * (1 - ylocal) * (1 - zlocal) +
zField[xindex][yindex][zindex + 1] * (1 - xlocal) * (1 - ylocal) * zlocal +
zField[xindex][yindex + 1][zindex] * (1 - xlocal) * ylocal * (1 - zlocal) +
zField[xindex][yindex + 1][zindex + 1] * (1 - xlocal) * ylocal * zlocal +
zField[xindex + 1][yindex][zindex] * xlocal * (1 - ylocal) * (1 - zlocal) +
zField[xindex + 1][yindex][zindex + 1] * xlocal * (1 - ylocal) * zlocal +
zField[xindex + 1][yindex + 1][zindex] * xlocal * ylocal * (1 - zlocal) +
zField[xindex + 1][yindex + 1][zindex + 1] * xlocal * ylocal * zlocal;