Skip to content

Commit

Permalink
add cartesian_3d field support
Browse files Browse the repository at this point in the history
  • Loading branch information
Zhiwen Zhao committed Sep 6, 2018
1 parent 78c4c42 commit 232e52d
Show file tree
Hide file tree
Showing 7 changed files with 339 additions and 6 deletions.
3 changes: 2 additions & 1 deletion SConstruct
View file
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,8 @@ field_sources = Split("""
fields/multipoleField.cc
fields/symmetries/dipole.cc
fields/symmetries/cylindrical.cc
fields/symmetries/phi-segmented.cc""")
fields/symmetries/phi-segmented.cc
fields/symmetries/cartesian_3d.cc""")
env.Library(source = field_sources, target = "lib/gfields")


Expand Down
11 changes: 8 additions & 3 deletions fields/asciiField.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -192,15 +192,20 @@ void asciiField::loadFieldMap(gMappedField* map, double v)
// dipole field
if(map->symmetry == "dipole-x" || map->symmetry == "dipole-y" || map->symmetry == "dipole-z")
loadFieldMap_Dipole(map, v);

// cylindrical field
if(map->symmetry == "cylindrical-x" || map->symmetry == "cylindrical-y" || map->symmetry == "cylindrical-z")
else if(map->symmetry == "cylindrical-x" || map->symmetry == "cylindrical-y" || map->symmetry == "cylindrical-z")
loadFieldMap_Cylindrical(map, v);

// phi-segmented field

if(map->symmetry == "phi-segmented")
else if(map->symmetry == "phi-segmented")
loadFieldMap_phiSegmented(map, v);

else if(map->symmetry == "cartesian_3D" || map->symmetry == "cartesian_3D_quadrant")
loadFieldMap_cartesian3d(map, v);

else {cout << "can't recognize the field symmetry "<< map->symmetry << endl; exit(0);}
}


Expand Down
3 changes: 2 additions & 1 deletion fields/asciiField.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@ class asciiField : public fieldFactory

void loadFieldMap_Dipole(gMappedField*, double); // load dipole field map
void loadFieldMap_Cylindrical(gMappedField*, double); // load cylindrical field map
void loadFieldMap_phiSegmented(gMappedField*, double); // load cylindrical field map
void loadFieldMap_phiSegmented(gMappedField*, double); // load phiSegmented field map
void loadFieldMap_cartesian3d(gMappedField*, double); // load cartesian3d field map

static fieldFactory *createFieldFactory()
{
Expand Down
4 changes: 3 additions & 1 deletion fields/fieldFactory.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,9 @@ class fieldFactory
virtual void loadFieldMap(gMappedField*, double) = 0; // load field map. This is a dispatcher function for the various types of fields below

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

string factoryType;

Expand Down
37 changes: 37 additions & 0 deletions fields/mappedField.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,12 @@ void gMappedField::GetFieldValue(const double point[3], double *bField) const
// phi-segmented
} else if(symmetry == "phi-segmented") {
GetFieldValue_phiSegmented(rpoint, bField, FIRST_ONLY);
// cartesian_3D
} else if(symmetry == "cartesian_3D" || symmetry == "cartesian_3D_quadrant"){
GetFieldValue_cartesian3d(rpoint, bField, FIRST_ONLY);
}



if(verbosity == 99) FIRST_ONLY = 99;

Expand Down Expand Up @@ -89,13 +94,16 @@ void gMappedField::initializeMap()
if(symmetry == "dipole-x" || symmetry == "dipole-y" || symmetry == "dipole-z")
{
startMap = new double[2];
endMap = new double[2];
cellSize = new double[2];
np = new unsigned int[2];

np[0] = getCoordinateWithName("longitudinal").np;
np[1] = getCoordinateWithName("transverse").np;
startMap[0] = getCoordinateWithName("longitudinal").min;
startMap[1] = getCoordinateWithName("transverse").min;
endMap[0] = getCoordinateWithName("longitudinal").max;
endMap[1] = getCoordinateWithName("transverse").max;
cellSize[0] = (getCoordinateWithName("longitudinal").max - startMap[0]) / (np[0] - 1);
cellSize[1] = (getCoordinateWithName("transverse").max - startMap[1]) / (np[1] - 1);
}
Expand All @@ -105,13 +113,16 @@ void gMappedField::initializeMap()
if(symmetry == "cylindrical-x" || symmetry == "cylindrical-y" || symmetry == "cylindrical-z")
{
startMap = new double[2];
endMap = new double[2];
cellSize = new double[2];
np = new unsigned int[2];

np[0] = getCoordinateWithName("transverse").np;
np[1] = getCoordinateWithName("longitudinal").np;
startMap[0] = getCoordinateWithName("transverse").min;
startMap[1] = getCoordinateWithName("longitudinal").min;
endMap[0] = getCoordinateWithName("longitudinal").max;
endMap[1] = getCoordinateWithName("transverse").max;
cellSize[0] = (getCoordinateWithName("transverse").max - startMap[0]) / (np[0] - 1);
cellSize[1] = (getCoordinateWithName("longitudinal").max - startMap[1]) / (np[1] - 1);
}
Expand All @@ -121,6 +132,7 @@ void gMappedField::initializeMap()
if(symmetry == "phi-segmented")
{
startMap = new double[3];
endMap = new double[3];
cellSize = new double[3];
np = new unsigned int[3];

Expand All @@ -130,11 +142,36 @@ void gMappedField::initializeMap()
startMap[0] = getCoordinateWithName("azimuthal").min;
startMap[1] = getCoordinateWithName("transverse").min;
startMap[2] = getCoordinateWithName("longitudinal").min;
endMap[0] = getCoordinateWithName("azimuthal").max;
endMap[1] = getCoordinateWithName("transverse").max;
endMap[2] = getCoordinateWithName("longitudinal").max;
cellSize[0] = (getCoordinateWithName("azimuthal").max - startMap[0]) / (np[0] - 1);
cellSize[1] = (getCoordinateWithName("transverse").max - startMap[1]) / (np[1] - 1);
cellSize[2] = (getCoordinateWithName("longitudinal").max - startMap[2]) / (np[2] - 1);
}

//cartesian_3D
if(symmetry == "cartesian_3D" || symmetry == "cartesian_3D_quadrant")
{
startMap = new double[3];
endMap = new double[3];
cellSize = new double[3];
np = new unsigned int[3];

np[0] = getCoordinateWithName("X").np;
np[1] = getCoordinateWithName("Y").np;
np[2] = getCoordinateWithName("Z").np;
startMap[0] = getCoordinateWithName("X").min;
startMap[1] = getCoordinateWithName("Y").min;
startMap[2] = getCoordinateWithName("Z").min;
endMap[0] = getCoordinateWithName("X").max;
endMap[1] = getCoordinateWithName("Y").max;
endMap[2] = getCoordinateWithName("Z").max;
cellSize[0] = (getCoordinateWithName("X").max - startMap[0]) / (np[0] - 1);
cellSize[1] = (getCoordinateWithName("Y").max - startMap[1]) / (np[1] - 1);
cellSize[2] = (getCoordinateWithName("Z").max - startMap[2]) / (np[2] - 1);
}

// setting rotation sin and cosines
sinAlpha = sin(mapRotation[0]);
cosAlhpa = cos(mapRotation[0]);
Expand Down
2 changes: 2 additions & 0 deletions fields/mappedField.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -119,6 +119,7 @@ class gMappedField : public G4MagneticField
// it will avoid the time spend in retrieving
// the infos in GetFieldValue
double *startMap;
double *endMap;
double *cellSize;
unsigned int *np;
void initializeMap();
Expand All @@ -132,6 +133,7 @@ class gMappedField : public G4MagneticField
void GetFieldValue_Dipole( const double x[3], double *Bfield, int FIRST_ONLY) const;
void GetFieldValue_Cylindrical( const double x[3], double *Bfield, int FIRST_ONLY) const;
void GetFieldValue_phiSegmented( const double x[3], double *Bfield, int FIRST_ONLY) const;
void GetFieldValue_cartesian3d( const double x[3], double *Bfield, int FIRST_ONLY) const;


// we want to rotate the field (axes), not the point
Expand Down
Loading

0 comments on commit 232e52d

Please sign in to comment.