4 #include <phparameter/PHParameters.h>
14 #include <Geant4/G4ChordFinder.hh>
15 #include <Geant4/G4ClassicalRK4.hh>
16 #include <Geant4/G4FieldManager.hh>
17 #include <Geant4/G4LogicalVolume.hh>
18 #include <Geant4/G4Mag_UsualEqRhs.hh>
19 #include <Geant4/G4MagneticField.hh>
20 #include <Geant4/G4PVPlacement.hh>
21 #include <Geant4/G4PhysicalConstants.hh>
22 #include <Geant4/G4QuadrupoleMagField.hh>
23 #include <Geant4/G4RotationMatrix.hh>
24 #include <Geant4/G4SystemOfUnits.hh>
25 #include <Geant4/G4ThreeVector.hh>
26 #include <Geant4/G4Transform3D.hh>
27 #include <Geant4/G4Tubs.hh>
28 #include <Geant4/G4Types.hh>
29 #include <Geant4/G4UniformMagField.hh>
31 #include <CLHEP/Units/SystemOfUnits.h>
43 , m_Params(parameters)
45 , m_MagnetId(magnet_id)
76 G4RotationMatrix *rotm =
new G4RotationMatrix();
82 G4VSolid *magnet_mother_solid =
new G4Tubs(
GetName().append(
"_Mother"),
86 G4LogicalVolume *magnet_mother_logic =
new G4LogicalVolume(magnet_mother_solid,
89 nullptr,
nullptr,
nullptr);
92 new G4PVPlacement(G4Transform3D(*rotm,
origin),
97 G4ThreeVector field_origin(
origin);
98 G4RotationMatrix *field_rotm = rotm;
105 std::cout << __PRETTY_FUNCTION__ <<
": set field using the global coordinate system, as the magnet is a daughter vol. of "
116 field_rotm =
new G4RotationMatrix();
123 if (magnettype ==
"DIPOLE")
131 field.transform(*field_rotm);
135 std::cout <<
"Creating DIPOLE with field x: " <<
field.x() / tesla
136 <<
", y: " <<
field.y() / tesla
137 <<
", z: " <<
field.z() / tesla
138 <<
" and name " <<
GetName() << std::endl;
141 else if (magnettype ==
"QUADRUPOLE")
148 m_magField =
new G4QuadrupoleMagField(fieldGradient, field_origin, field_rotm);
152 std::cout <<
"Creating QUADRUPOLE with gradient " << fieldGradient <<
" and name " <<
GetName() << std::endl;
153 std::cout <<
"at x, y, z = " << field_origin.x() <<
" , " << field_origin.y() <<
" , " << field_origin.z() << std::endl;
154 std::cout <<
"with rotation around x, y, z axis by: " << field_rotm->phiX() <<
", " << field_rotm->phiY() <<
", " << field_rotm->phiZ() << std::endl;
159 std::cout << __PRETTY_FUNCTION__ <<
" : Fatal error, magnet type of " << magnettype <<
" is not yet supported" << std::endl;
165 std::cout <<
PHWHERE <<
" No magnetic field specified for " <<
GetName()
166 <<
" of type " << magnettype << std::endl;
170 G4VSolid *magnet_iron_solid =
new G4Tubs(
GetName().append(
"_Solid"),
174 G4LogicalVolume *magnet_iron_logic =
new G4LogicalVolume(magnet_iron_solid,
177 nullptr,
nullptr,
nullptr);
185 std::cout <<
"Magnet skin thickness " <<
m_Params->
get_double_param(
"skin_thickness") <<
" too large, exceeds 2xmagnet thickness "
191 G4VSolid *magnet_core_solid =
new G4Tubs(
GetName().append(
"_Core"),
195 G4LogicalVolume *magnet_core_logic =
new G4LogicalVolume(magnet_core_solid,
198 nullptr,
nullptr,
nullptr);
211 G4VSolid *magnet_field_solid =
new G4Tubs(
GetName().append(
"_Field_Solid"),
219 nullptr,
nullptr,
nullptr);
229 magnet_physi =
new G4PVPlacement(
nullptr, G4ThreeVector(0, 0, 0),
246 std::cout << __PRETTY_FUNCTION__ <<
": set field to vol " <<
m_magnetFieldLogic->GetName() <<
" that include "
251 G4Mag_UsualEqRhs *localEquation =
new G4Mag_UsualEqRhs(
m_magField);
252 G4ClassicalRK4 *localStepper =
new G4ClassicalRK4(localEquation);
253 G4double minStep = 0.25 *
mm;
254 G4ChordFinder *localChordFinder =
new G4ChordFinder(
m_magField, minStep, localStepper);
256 G4FieldManager *fieldMgr =
new G4FieldManager();
258 fieldMgr->SetChordFinder(localChordFinder);