Skip to content

Commit bae360e

Browse files
maureeungaromaureeungaro
authored
Merge pull request #145 from gemc-eic/master
add cartesian_3d field support
2 parents 0f1d021 + 232e52d commit bae360e

File tree

7 files changed

+339
-6
lines changed

7 files changed

+339
-6
lines changed

‎SConstruct

+2-1
Original file line numberDiff line numberDiff line change
@@ -93,7 +93,8 @@ field_sources = Split("""
9393
fields/multipoleField.cc
9494
fields/symmetries/dipole.cc
9595
fields/symmetries/cylindrical.cc
96-
fields/symmetries/phi-segmented.cc""")
96+
fields/symmetries/phi-segmented.cc
97+
fields/symmetries/cartesian_3d.cc""")
9798
env.Library(source = field_sources, target = "lib/gfields")
9899

99100

‎fields/asciiField.cc

+8-3
Original file line numberDiff line numberDiff line change
@@ -192,15 +192,20 @@ void asciiField::loadFieldMap(gMappedField* map, double v)
192192
// dipole field
193193
if(map->symmetry == "dipole-x" || map->symmetry == "dipole-y" || map->symmetry == "dipole-z")
194194
loadFieldMap_Dipole(map, v);
195-
195+
196196
// cylindrical field
197-
if(map->symmetry == "cylindrical-x" || map->symmetry == "cylindrical-y" || map->symmetry == "cylindrical-z")
197+
else if(map->symmetry == "cylindrical-x" || map->symmetry == "cylindrical-y" || map->symmetry == "cylindrical-z")
198198
loadFieldMap_Cylindrical(map, v);
199199

200200
// phi-segmented field
201201

202-
if(map->symmetry == "phi-segmented")
202+
else if(map->symmetry == "phi-segmented")
203203
loadFieldMap_phiSegmented(map, v);
204+
205+
else if(map->symmetry == "cartesian_3D" || map->symmetry == "cartesian_3D_quadrant")
206+
loadFieldMap_cartesian3d(map, v);
207+
208+
else {cout << "can't recognize the field symmetry "<< map->symmetry << endl; exit(0);}
204209
}
205210

206211

‎fields/asciiField.h

+2-1
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,8 @@ class asciiField : public fieldFactory
1919

2020
void loadFieldMap_Dipole(gMappedField*, double); // load dipole field map
2121
void loadFieldMap_Cylindrical(gMappedField*, double); // load cylindrical field map
22-
void loadFieldMap_phiSegmented(gMappedField*, double); // load cylindrical field map
22+
void loadFieldMap_phiSegmented(gMappedField*, double); // load phiSegmented field map
23+
void loadFieldMap_cartesian3d(gMappedField*, double); // load cartesian3d field map
2324

2425
static fieldFactory *createFieldFactory()
2526
{

‎fields/fieldFactory.h

+3-1
Original file line numberDiff line numberDiff line change
@@ -20,7 +20,9 @@ class fieldFactory
2020
virtual void loadFieldMap(gMappedField*, double) = 0; // load field map. This is a dispatcher function for the various types of fields below
2121

2222
virtual void loadFieldMap_Dipole(gMappedField*, double) = 0; // load 1D dipole field depending on 2 coordinates (transverse and longitudinal)
23-
virtual void loadFieldMap_Cylindrical(gMappedField*, double) = 0; // load phi-symmetric field depending on 2 coordinates (transverse and longitudinal)
23+
virtual void loadFieldMap_Cylindrical(gMappedField*, double) = 0; // load cylindrical field depending on 2 coordinates (transverse and longitudinal)
24+
virtual void loadFieldMap_phiSegmented(gMappedField*, double) = 0; // load phi-symmetric field depending on 3 coordinates (azimuthal,transverse and longitudinal)
25+
virtual void loadFieldMap_cartesian3d(gMappedField*, double) = 0; // load cartesian_3D field depending on 3 coordinates (XX,YY,ZZ)
2426

2527
string factoryType;
2628

‎fields/mappedField.cc

+37
Original file line numberDiff line numberDiff line change
@@ -46,7 +46,12 @@ void gMappedField::GetFieldValue(const double point[3], double *bField) const
4646
// phi-segmented
4747
} else if(symmetry == "phi-segmented") {
4848
GetFieldValue_phiSegmented(rpoint, bField, FIRST_ONLY);
49+
// cartesian_3D
50+
} else if(symmetry == "cartesian_3D" || symmetry == "cartesian_3D_quadrant"){
51+
GetFieldValue_cartesian3d(rpoint, bField, FIRST_ONLY);
4952
}
53+
54+
5055

5156
if(verbosity == 99) FIRST_ONLY = 99;
5257

@@ -89,13 +94,16 @@ void gMappedField::initializeMap()
8994
if(symmetry == "dipole-x" || symmetry == "dipole-y" || symmetry == "dipole-z")
9095
{
9196
startMap = new double[2];
97+
endMap = new double[2];
9298
cellSize = new double[2];
9399
np = new unsigned int[2];
94100

95101
np[0] = getCoordinateWithName("longitudinal").np;
96102
np[1] = getCoordinateWithName("transverse").np;
97103
startMap[0] = getCoordinateWithName("longitudinal").min;
98104
startMap[1] = getCoordinateWithName("transverse").min;
105+
endMap[0] = getCoordinateWithName("longitudinal").max;
106+
endMap[1] = getCoordinateWithName("transverse").max;
99107
cellSize[0] = (getCoordinateWithName("longitudinal").max - startMap[0]) / (np[0] - 1);
100108
cellSize[1] = (getCoordinateWithName("transverse").max - startMap[1]) / (np[1] - 1);
101109
}
@@ -105,13 +113,16 @@ void gMappedField::initializeMap()
105113
if(symmetry == "cylindrical-x" || symmetry == "cylindrical-y" || symmetry == "cylindrical-z")
106114
{
107115
startMap = new double[2];
116+
endMap = new double[2];
108117
cellSize = new double[2];
109118
np = new unsigned int[2];
110119

111120
np[0] = getCoordinateWithName("transverse").np;
112121
np[1] = getCoordinateWithName("longitudinal").np;
113122
startMap[0] = getCoordinateWithName("transverse").min;
114123
startMap[1] = getCoordinateWithName("longitudinal").min;
124+
endMap[0] = getCoordinateWithName("longitudinal").max;
125+
endMap[1] = getCoordinateWithName("transverse").max;
115126
cellSize[0] = (getCoordinateWithName("transverse").max - startMap[0]) / (np[0] - 1);
116127
cellSize[1] = (getCoordinateWithName("longitudinal").max - startMap[1]) / (np[1] - 1);
117128
}
@@ -121,6 +132,7 @@ void gMappedField::initializeMap()
121132
if(symmetry == "phi-segmented")
122133
{
123134
startMap = new double[3];
135+
endMap = new double[3];
124136
cellSize = new double[3];
125137
np = new unsigned int[3];
126138

@@ -130,11 +142,36 @@ void gMappedField::initializeMap()
130142
startMap[0] = getCoordinateWithName("azimuthal").min;
131143
startMap[1] = getCoordinateWithName("transverse").min;
132144
startMap[2] = getCoordinateWithName("longitudinal").min;
145+
endMap[0] = getCoordinateWithName("azimuthal").max;
146+
endMap[1] = getCoordinateWithName("transverse").max;
147+
endMap[2] = getCoordinateWithName("longitudinal").max;
133148
cellSize[0] = (getCoordinateWithName("azimuthal").max - startMap[0]) / (np[0] - 1);
134149
cellSize[1] = (getCoordinateWithName("transverse").max - startMap[1]) / (np[1] - 1);
135150
cellSize[2] = (getCoordinateWithName("longitudinal").max - startMap[2]) / (np[2] - 1);
136151
}
137152

153+
//cartesian_3D
154+
if(symmetry == "cartesian_3D" || symmetry == "cartesian_3D_quadrant")
155+
{
156+
startMap = new double[3];
157+
endMap = new double[3];
158+
cellSize = new double[3];
159+
np = new unsigned int[3];
160+
161+
np[0] = getCoordinateWithName("X").np;
162+
np[1] = getCoordinateWithName("Y").np;
163+
np[2] = getCoordinateWithName("Z").np;
164+
startMap[0] = getCoordinateWithName("X").min;
165+
startMap[1] = getCoordinateWithName("Y").min;
166+
startMap[2] = getCoordinateWithName("Z").min;
167+
endMap[0] = getCoordinateWithName("X").max;
168+
endMap[1] = getCoordinateWithName("Y").max;
169+
endMap[2] = getCoordinateWithName("Z").max;
170+
cellSize[0] = (getCoordinateWithName("X").max - startMap[0]) / (np[0] - 1);
171+
cellSize[1] = (getCoordinateWithName("Y").max - startMap[1]) / (np[1] - 1);
172+
cellSize[2] = (getCoordinateWithName("Z").max - startMap[2]) / (np[2] - 1);
173+
}
174+
138175
// setting rotation sin and cosines
139176
sinAlpha = sin(mapRotation[0]);
140177
cosAlhpa = cos(mapRotation[0]);

‎fields/mappedField.h

+2
Original file line numberDiff line numberDiff line change
@@ -119,6 +119,7 @@ class gMappedField : public G4MagneticField
119119
// it will avoid the time spend in retrieving
120120
// the infos in GetFieldValue
121121
double *startMap;
122+
double *endMap;
122123
double *cellSize;
123124
unsigned int *np;
124125
void initializeMap();
@@ -132,6 +133,7 @@ class gMappedField : public G4MagneticField
132133
void GetFieldValue_Dipole( const double x[3], double *Bfield, int FIRST_ONLY) const;
133134
void GetFieldValue_Cylindrical( const double x[3], double *Bfield, int FIRST_ONLY) const;
134135
void GetFieldValue_phiSegmented( const double x[3], double *Bfield, int FIRST_ONLY) const;
136+
void GetFieldValue_cartesian3d( const double x[3], double *Bfield, int FIRST_ONLY) const;
135137

136138

137139
// we want to rotate the field (axes), not the point

0 commit comments

Comments
 (0)