[807] | 1 | |
---|
| 2 | ========================================================= |
---|
| 3 | Geant4 - an Object-Oriented Toolkit for Simulation in HEP |
---|
| 4 | ========================================================= |
---|
| 5 | |
---|
| 6 | |
---|
| 7 | |
---|
| 8 | field04 Example |
---|
| 9 | --------------- |
---|
| 10 | |
---|
| 11 | This example shows how to define/use OVERLAPPING field elements |
---|
| 12 | in Geant4. Fields might be either magnetic, electric or both. |
---|
| 13 | |
---|
| 14 | Credit goes to Tom Roberts and Muons Inc. since much of the code |
---|
| 15 | and ideas were taken at liberty from the (GNU GPL) source of |
---|
| 16 | G4BEAMLINE release 1.12. |
---|
| 17 | |
---|
| 18 | http://g4beamline.muonsinc.com |
---|
| 19 | |
---|
| 20 | ************** |
---|
| 21 | *Classes Used* |
---|
| 22 | ************** |
---|
| 23 | |
---|
| 24 | 1 - main() |
---|
| 25 | |
---|
| 26 | Note: |
---|
| 27 | |
---|
| 28 | // runManager->Initialize(); |
---|
| 29 | |
---|
| 30 | is commented out and must be executed with: /run/initialize |
---|
| 31 | either in a macro or on the command line. Note: One can change |
---|
| 32 | the geometry (see below) only BEFORE the initialization! |
---|
| 33 | |
---|
| 34 | Is is possible to assign the PhysicsList via an argument to |
---|
| 35 | program submission: |
---|
| 36 | |
---|
| 37 | field04 -p QGSP_BERT |
---|
| 38 | |
---|
| 39 | and an initial random number seed with: |
---|
| 40 | |
---|
| 41 | field04 field04.in 12345 |
---|
| 42 | |
---|
| 43 | It is also possible to specify more than one macro file which are |
---|
| 44 | then executed in succession: |
---|
| 45 | |
---|
| 46 | field04 field04.in field04.in 12345 |
---|
| 47 | |
---|
| 48 | (Note, in this case the initial random # must also be provided) |
---|
| 49 | |
---|
| 50 | (1) argc holds the number of arguments on the command line. Since |
---|
| 51 | this MUST include at least the program name, argc value is >= 1! |
---|
| 52 | (2) argv[0] is always the name of the program |
---|
| 53 | argv[1] points to the first argument, and so on ... |
---|
| 54 | (3) argv[argc-1] is assumed to be the intial random # seed |
---|
| 55 | (4) the loop over macro files goes from optind = 1 to (argc-1) |
---|
| 56 | |
---|
| 57 | |
---|
| 58 | 2- GEOMETRY DEFINITION |
---|
| 59 | |
---|
| 60 | The geometry consists of two solenoidal magnets: a "CaptureMgnt" |
---|
| 61 | followed by a (bloe-colored "TransferMgnt". By definition, the |
---|
| 62 | axis and center of the "CaptureMgnt" coincide with the "World". The |
---|
| 63 | position of the "TransferMgnt" relative to the downstream end of the |
---|
| 64 | "CaptureMgnt", as well as its axis angle, both may vary. A cylindrical |
---|
| 65 | "Target" is positioned inside the "CaptureMgnt". Its axis can vary |
---|
| 66 | from 0 to 180 deg, and hence also the direction of the incoming |
---|
| 67 | proton beam wrt the "CaptureMgnt"'s axis. A "Degrader" is located |
---|
| 68 | inside the "TransferMgnt", its default position being at the |
---|
| 69 | upstream end of the "TransferMgnt". Finally, also a "TestPlane" is |
---|
| 70 | located inside the "TransferMgnt", by default at its downstream end. |
---|
| 71 | |
---|
| 72 | |
---|
| 73 | The "World" consists of a solid cylinder made of a given material. |
---|
| 74 | (It is the responsibility of the user to make the world |
---|
| 75 | large enough to contain the rest of the geometry!) |
---|
| 76 | |
---|
| 77 | Three parameters define the world : |
---|
| 78 | - the material of the world, |
---|
| 79 | - the world radius, |
---|
| 80 | - the world length. |
---|
| 81 | |
---|
| 82 | Example (default values): |
---|
| 83 | /field04/SetWorldMat G4_AIR |
---|
| 84 | /field04/SetWorldR 5.0 m |
---|
| 85 | /field04/SetWorldZ 50.0 m |
---|
| 86 | |
---|
| 87 | |
---|
| 88 | The "Target" is a solid cylinder made of a given material. |
---|
| 89 | |
---|
| 90 | Five parameters define the target: |
---|
| 91 | - the material of the target, |
---|
| 92 | - the target radius, |
---|
| 93 | - the target thickness, |
---|
| 94 | - the target position inside the "CaptureMgnt", |
---|
| 95 | - the target axis angle relative to that of the "CaptureMgnt". |
---|
| 96 | |
---|
| 97 | Example (default values): |
---|
| 98 | /field04/SetTgtMat G4_W |
---|
| 99 | /field04/SetTgtRad 0.4 cm |
---|
| 100 | /field04/SetTgtThick 16.0 cm |
---|
| 101 | /field04/SetTgtPos 0.0 cm |
---|
| 102 | /field04/SetTgtAng 170 |
---|
| 103 | |
---|
| 104 | |
---|
| 105 | The "Degrader" is a solid cylinder made of a given material. |
---|
| 106 | |
---|
| 107 | Four parameters define the degrader: |
---|
| 108 | - the material of the degrader, |
---|
| 109 | - the degrader radius, |
---|
| 110 | - the degrader thickness, |
---|
| 111 | - the degrader position relative to the "TransferMgnt" center. |
---|
| 112 | |
---|
| 113 | Example (default values): |
---|
| 114 | /field04/SetDgrMat G4_Pb |
---|
| 115 | /field04/SetDgrRad 30.0 cm |
---|
| 116 | /field04/SetDgrThick 0.1 cm |
---|
| 117 | #/field04/SetDgrPos -7.4 m |
---|
| 118 | |
---|
| 119 | |
---|
| 120 | The "CaptureMgnt" is a solenoid (vacuum cylinder). It is either |
---|
| 121 | a two-sided or a one-sided magnetic bottle with the B field |
---|
| 122 | varying linearly from the center value B1 to the edge value B2. |
---|
| 123 | The one-sided 'FocusSolenoid' has the open end at +z and focuses |
---|
| 124 | on the z < 0 side. |
---|
| 125 | |
---|
| 126 | Four parameters define the "CaptureMgnt": |
---|
| 127 | - the magnet radius, |
---|
| 128 | - the magnet length, |
---|
| 129 | - the weaker magnetic field at the center B1 |
---|
| 130 | - the stronger magnetic field at the edge B2 |
---|
| 131 | |
---|
| 132 | Example (default values): |
---|
| 133 | /field04/SetCaptureR 0.6 m |
---|
| 134 | /field04/SetCaptureZ 4.0 m |
---|
| 135 | /field04/SetCaptureB1 2.5 tesla |
---|
| 136 | /field04/SetCaptureB2 5.0 tesla |
---|
| 137 | |
---|
| 138 | |
---|
| 139 | The "TransferMgnt" is a solenoid (vacuum cylinder) with a |
---|
| 140 | constant B-field. When the "TransferMgnt" follows immediately |
---|
| 141 | the "CaptureMgnt", its relative position is at 0 cm. |
---|
| 142 | |
---|
| 143 | Four parameters define the "TransferMgnt": |
---|
| 144 | - the magnet radius, |
---|
| 145 | - the magnet length, |
---|
| 146 | - the magnet field, |
---|
| 147 | - the magnet relative position |
---|
| 148 | (its upstream face wrt the downstream face of the "CaptureMgnt".) |
---|
| 149 | |
---|
| 150 | Example (default values): |
---|
| 151 | /field04/SetTransferR 0.3 m |
---|
| 152 | /field04/SetTransferZ 15.0 m |
---|
| 153 | /field04/SetTransferB 5.0 tesla |
---|
| 154 | /field04/SetTransferP 0.0 m |
---|
| 155 | |
---|
| 156 | The default geometry is constructed in DetectorConstruction class, |
---|
| 157 | but all the parameters can be changed via the commands defined in |
---|
| 158 | the DetectorMessenger class. |
---|
| 159 | |
---|
| 160 | |
---|
| 161 | 3- MATERIAL DEFINITION |
---|
| 162 | |
---|
| 163 | Material definitions are done through the singleton class Materials |
---|
| 164 | which keeps a pointer to the G4NistManager. It has a method |
---|
| 165 | GetMaterial by name (G4String) which in turn invokes the |
---|
| 166 | G4NistManager::FindOrBuildMaterial, and/or G4Material::GetMaterial |
---|
| 167 | methods. It has also a method CreateMaterials which, for materials |
---|
| 168 | absent from the NIST data base, shows how to create them using the |
---|
| 169 | G4NistManager::ConstructNewMaterial method. |
---|
| 170 | |
---|
| 171 | |
---|
| 172 | 4- AN EVENT: THE PRIMARY GENERATOR |
---|
| 173 | |
---|
| 174 | The primary kinematic consists of a single particle which hits the |
---|
| 175 | target perpendicular to its upstream face. The type of the particle |
---|
| 176 | and its energy are set in the PrimaryGeneratorAction class, and can |
---|
| 177 | be changed via the G4 build-in commands of the ParticleGun class. |
---|
| 178 | In addition, there is a rndmFlag, which once set allows the beam to |
---|
| 179 | explore randomly the whole cross section of the target. The default |
---|
| 180 | beam consists of 500 MeV protons, starting at the upstream face of |
---|
| 181 | the target, directed along dx = dy = 0, dz = 1 wrt the target frame. |
---|
| 182 | The default direction should NOT be changed! The arguments of the |
---|
| 183 | x/y/zvertex commands are relative to the target center. |
---|
| 184 | |
---|
| 185 | Example: |
---|
| 186 | /gun/random on |
---|
| 187 | #/gun/xvertex 0 mm |
---|
| 188 | #/gun/yvertex 0 mm |
---|
| 189 | #/gun/zvertex -100 mm |
---|
| 190 | |
---|
| 191 | |
---|
| 192 | 5- DETECTOR RESPONSE |
---|
| 193 | |
---|
| 194 | Information is extracted from the program via SteppingAction |
---|
| 195 | at the TestPlane. |
---|
| 196 | |
---|
| 197 | |
---|
| 198 | 6- PHYSICS |
---|
| 199 | |
---|
| 200 | The PhysicsList is adopted from examples/extended/hadronic/Hadr01. |
---|
| 201 | It uses components which are distributed with Geant4 in |
---|
| 202 | /geant4/physics_lists subdirectory. So, before compiling field04 |
---|
| 203 | it is necessary to compile physics_lists. |
---|
| 204 | |
---|
| 205 | There are two std::vector<G4VPhysicsConstructor*> PhysicsListVectors; |
---|
| 206 | one for EM- and one for HadronPHysics and a user defined process |
---|
| 207 | "StepMax". |
---|
| 208 | |
---|
| 209 | The choice of the physics is provided by the UI command: |
---|
| 210 | |
---|
| 211 | /exp/phys/addPhysics emstandard_opt1 |
---|
| 212 | /exp/phys/addPhysics QGSP_BERT |
---|
| 213 | |
---|
| 214 | To see the list of available configurations one can use: |
---|
| 215 | |
---|
| 216 | /exp/phys/list |
---|
| 217 | |
---|
| 218 | The command: |
---|
| 219 | |
---|
| 220 | /exp/phys/addPhysics PHYSLIST |
---|
| 221 | |
---|
| 222 | allows to download a physics configuration defined by an |
---|
| 223 | environment variable PHYSLIST. |
---|
| 224 | |
---|
| 225 | The cuts for electromagnetic phsyics can be established via: |
---|
| 226 | |
---|
| 227 | /exp/phys/allCuts 1 mm |
---|
| 228 | /exp/phys/gammaCut 0.1 mm |
---|
| 229 | /exp/phys/electronCut 0.2 mm |
---|
| 230 | /exp/phys/positronCut 0.3 mm |
---|
| 231 | |
---|
| 232 | The cut for the StepMax process via: |
---|
| 233 | |
---|
| 234 | /exp/phys/stepMax |
---|
| 235 | |
---|
| 236 | The two PhysicsListVectors can be cleared via: |
---|
| 237 | |
---|
| 238 | /exp/phys/clearEMPhysics |
---|
| 239 | /exp/phys/clearHadronPhysics |
---|
| 240 | |
---|
| 241 | and individual processes can be removed from the vectors via, |
---|
| 242 | for example: |
---|
| 243 | |
---|
| 244 | /exp/phys/removeEMPhysics msc |
---|
| 245 | /exp/phys/removeHadronPhysics gamma_nuc |
---|
| 246 | |
---|
| 247 | Furthermore, the following commands are also available, but |
---|
| 248 | may only be used AFTER /run/initialize |
---|
| 249 | |
---|
| 250 | /process/inactivate msc |
---|
| 251 | /process/activate msc |
---|
| 252 | |
---|
| 253 | The decay of pions can be assigned via (pi -> e nu, pi -> mu nu): |
---|
| 254 | |
---|
| 255 | /decay/pienu |
---|
| 256 | /decay/pimunu |
---|
| 257 | |
---|
| 258 | The pienu assignment includes a small fraction of radiative decay: |
---|
| 259 | e nu gamma (G4PionRadiativeDecayChannel). |
---|
| 260 | |
---|
| 261 | The standard/default muon decay chain is modified to be 98.6% |
---|
| 262 | G4MuonDecayChannelWithSpin and 1.4% G4MuonRadiativeDecayChannelWithSpin |
---|
| 263 | in ConstructParticle(). |
---|
| 264 | |
---|
| 265 | The pion decay process G4PolDecay inherits from G4Decay and implements |
---|
| 266 | the virtual method - empty in the base class - DaughterPolarization |
---|
| 267 | |
---|
| 268 | The muon decay process is G4DecayWithSpin |
---|
| 269 | |
---|
| 270 | |
---|
| 271 | 7- Overlapping Fields |
---|
| 272 | |
---|
| 273 | The GlobalField (a singleton) is instantiated in DetectorConstruction() |
---|
| 274 | and assigned to the global field manager in updateField(): |
---|
| 275 | |
---|
| 276 | fFieldManager = GetGlobalFieldManager(); |
---|
| 277 | fFieldManager->SetDetectorField(this); |
---|
| 278 | |
---|
| 279 | The GlobalField has a std::vector<ElementField*> FieldList |
---|
| 280 | |
---|
| 281 | The field from each individual beamline element is given by a |
---|
| 282 | ElementField object. Any number of overlapping ElementField |
---|
| 283 | objects can be added to the global field. Any element that |
---|
| 284 | represents an element with an EM field must add the appropriate |
---|
| 285 | ElementField to the global GlobalField object. |
---|
| 286 | |
---|
| 287 | Of course, the GlobalField has the method GetFieldValue implemented. |
---|
| 288 | |
---|
| 289 | Before /run/initialize in the macro file or command, the update |
---|
| 290 | field command must have been issued if any of the other following |
---|
| 291 | field commands was employed: |
---|
| 292 | |
---|
| 293 | /field/update |
---|
| 294 | |
---|
| 295 | Other options are: |
---|
| 296 | |
---|
| 297 | /field/setStepperType 4 |
---|
| 298 | /field/setMinStep 10 mm |
---|
| 299 | /field/setDeltaChord 3.0 mm |
---|
| 300 | /field/setDeltaOneStep 0.01 mm |
---|
| 301 | /field/setDeltaIntersection 0.1 mm |
---|
| 302 | /field/setEpsMin 2.5e-7 mm |
---|
| 303 | /field/setEpsMax 0.05 mm |
---|
| 304 | |
---|
| 305 | Each field element has a rectilinear bounding box in global |
---|
| 306 | coordinate space which is checked before a point is verified to |
---|
| 307 | actually be inside the ElementField (isWithin and isOutside). |
---|
| 308 | setGlobalPoint is called 8 times for the corners of the local |
---|
| 309 | bounding box, after a local->global coordinate transform. |
---|
| 310 | |
---|
| 311 | The ElementField is the interface class used by GlobalField to |
---|
| 312 | compute the field value at a given point[]. |
---|
| 313 | |
---|
| 314 | A beamline element, for example the SimpleSolenoid, will derive |
---|
| 315 | from ElementField and implement the computation for the element. |
---|
| 316 | |
---|
| 317 | simpleSolenoid = new SimpleSolenoid(B, l, |
---|
| 318 | logicTransferMgnt,TransferMgntCenter); |
---|
| 319 | |
---|
| 320 | Besides the magnetic field and the length of the simple solenoid, |
---|
| 321 | the constructor needs the knowledge of the G4LogicalVolume for |
---|
| 322 | the beamline element and where its center is located in the |
---|
| 323 | 'World'. |
---|
| 324 | |
---|
| 325 | The ElementField has a G4AffineTransform "global2local" which |
---|
| 326 | allows the quick computation of coordinate transformations. It can |
---|
| 327 | only be determined by knowing the element's coordinate origin in |
---|
| 328 | the global frame and after all of the geometry has been defined. |
---|
| 329 | For this reason, the object is prepared in two stages, through the |
---|
| 330 | constructor providing it with the coordinate center and a pointer |
---|
| 331 | to the G4LogicalVolume. Later the construct() method is called to |
---|
| 332 | calculate the global2local and the bounding box. This can be done |
---|
| 333 | from the RunAction::BeginOfRunAction method, for only then are we |
---|
| 334 | certain that the geometry has been completely built: |
---|
| 335 | |
---|
| 336 | FieldList* fields = GlobalField::getObject()->getFields(); |
---|
| 337 | |
---|
| 338 | if (fields) { |
---|
| 339 | if (fields->size()>0) { |
---|
| 340 | FieldList::iterator i; |
---|
| 341 | for (i=fields->begin(); i!=fields->end(); ++i)(*i)->construct(); |
---|
| 342 | } |
---|
| 343 | } |
---|
| 344 | |
---|
| 345 | The ElementField constructor will also add the derived object into |
---|
| 346 | GlobalField. Finally, its addFieldValue() will add the field value |
---|
| 347 | for this element to field[]. |
---|
| 348 | |
---|
| 349 | |
---|
| 350 | 8- User Action Classes |
---|
| 351 | |
---|
| 352 | RunActionMessenger: |
---|
| 353 | |
---|
| 354 | /rndm/save freq - to save rndm status in external files |
---|
| 355 | 0 not saved |
---|
| 356 | >0 saved on: beginOfRun.rndm |
---|
| 357 | 1 saved on: endOfRun.rndm |
---|
| 358 | 2 saved on: endOfEvent.rndm |
---|
| 359 | /rndm/read random/run0evt8268.rndm |
---|
| 360 | |
---|
| 361 | RunAction: |
---|
| 362 | BeginOfRunAction: Deal with random number storage, |
---|
| 363 | initialization etc. Call the construct() method of |
---|
| 364 | ElementFields in the FieldList of GlobalField object. |
---|
| 365 | EndOfRunAction: random number storage/status printing. |
---|
| 366 | |
---|
| 367 | EventActionMessenger: |
---|
| 368 | /event/setverbose |
---|
| 369 | /event/drawTracks none/charged/all |
---|
| 370 | /event/printModulo |
---|
| 371 | |
---|
| 372 | EventAction(RunAction* RA): |
---|
| 373 | Customized BeginOfEvent printing (frequency: printModulo) |
---|
| 374 | EndofEvent: condition depending drawing of trajectories, |
---|
| 375 | saveEngingStatus and showEngineStatus according to flag |
---|
| 376 | in RunAction |
---|
| 377 | |
---|
| 378 | TrackingAction: |
---|
| 379 | PreUserTrackingAction: Instantiate UserTrackInformation |
---|
| 380 | and set the application TrackStatus. |
---|
| 381 | PostUserTrackingAction: Retreive UserTrackInformation |
---|
| 382 | and decide to save random number status accordingly. |
---|
| 383 | |
---|
| 384 | SteppingActionMessenger: |
---|
| 385 | |
---|
| 386 | SteppingAction: |
---|
| 387 | UserSteppingAction: Kill primary if/when outside Target |
---|
| 388 | volume. Diagnostic/histogram filling for particles at a |
---|
| 389 | TestPlane. Find decay position and when particle |
---|
| 390 | FIRST reverses z-momentum component via using a |
---|
| 391 | UserTrackInformation object. |
---|
| 392 | |
---|
| 393 | StackingAction: |
---|
| 394 | Track only primaries, pi+ or mu+ |
---|
| 395 | |
---|
| 396 | UserTrackInformation: |
---|
| 397 | Keep a application TrackStatus for the track: |
---|
| 398 | undefined, left, right, reverse |
---|
| 399 | |
---|
| 400 | SteppingVerbose: |
---|
| 401 | Only print track header and step information for |
---|
| 402 | pi+ and mu+ |
---|
| 403 | |
---|
| 404 | Trajectory, TrajectoryPoint: |
---|
| 405 | Example of application specific implementations |
---|
| 406 | |
---|
| 407 | 9- HOW TO START ? |
---|
| 408 | |
---|
| 409 | - compile and link to generate an executable |
---|
| 410 | % cd $G4INSTALL/example/extended/field/field04 |
---|
| 411 | % gmake |
---|
| 412 | |
---|
| 413 | - execute field04 in 'batch' mode from macro files e.g. |
---|
| 414 | % $(G4INSTALL)/bin/$(G4SYSTEM)/field04 field04.in |
---|
| 415 | |
---|
| 416 | - execute field04 in 'interactive' mode with visualization e.g. |
---|
| 417 | % $(G4INSTALL)/bin/$(G4SYSTEM)/field04 |
---|
| 418 | .... |
---|
| 419 | Idle> type your commands, for example: |
---|
| 420 | Idle> control/execute field04.in |
---|
| 421 | .... |
---|