source: trunk/examples/extended/field/field04/src/F04PhysicsList.cc@ 1036

Last change on this file since 1036 was 807, checked in by garnier, 17 years ago

update
File size: 21.0 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26//
27//
28
29#include "F04PhysicsList.hh"
30#include "F04PhysicsListMessenger.hh"
31
32#include "F04ExtraPhysics.hh"
33//#include "F04OpticalPhysics.hh"
34
35#include "G4DecayPhysics.hh"
36
37#include "G4EmStandardPhysics_option1.hh"
38#include "G4EmStandardPhysics_option2.hh"
39#include "G4EmStandardPhysics.hh"
40
41#include "G4EmExtraPhysics.hh"
42#include "G4EmProcessOptions.hh"
43
44#include "G4HadronElasticPhysics.hh"
45#include "G4HadronQElasticPhysics.hh"
46#include "G4HadronHElasticPhysics.hh"
47#include "G4NeutronTrackingCut.hh"
48#include "G4QStoppingPhysics.hh"
49#include "G4LHEPStoppingPhysics.hh"
50#include "G4IonBinaryCascadePhysics.hh"
51#include "G4IonPhysics.hh"
52
53#include "HadronPhysicsFTFP.hh"
54#include "HadronPhysicsFTFC.hh"
55
56#include "HadronPhysicsLHEP.hh"
57#include "HadronPhysicsLHEP_BERT.hh"
58#include "HadronPhysicsLHEP_EMV.hh"
59#include "HadronPhysicsLHEP_PRECO_HP.hh"
60
61#include "HadronPhysicsQGSC.hh"
62#include "HadronPhysicsQGSC_EFLOW.hh"
63
64#include "HadronPhysicsQGSP.hh"
65#include "HadronPhysicsQGSP_BERT.hh"
66#include "HadronPhysicsQGSP_BERT_HP.hh"
67#include "HadronPhysicsQGSP_BERT_TRV.hh"
68#include "HadronPhysicsQGSP_BIC.hh"
69#include "HadronPhysicsQGSP_BIC_HP.hh"
70
71#include "G4HadronInelasticQBBC.hh"
72#include "G4HadronInelasticQLHEP.hh"
73
74#include "G4IonPhysics.hh"
75
76#include "G4HadronProcessStore.hh"
77
78#include "G4LossTableManager.hh"
79
80#include "G4ProcessManager.hh"
81#include "G4ParticleTypes.hh"
82#include "G4ParticleTable.hh"
83
84#include "G4Gamma.hh"
85#include "G4Electron.hh"
86#include "G4Positron.hh"
87
88#include "F04StepMax.hh"
89
90#include "G4ProcessTable.hh"
91
92#include "G4PionDecayMakeSpin.hh"
93#include "G4DecayWithSpin.hh"
94
95#include "G4DecayTable.hh"
96#include "G4MuonDecayChannelWithSpin.hh"
97#include "G4MuonRadiativeDecayChannelWithSpin.hh"
98
99F04PhysicsList::F04PhysicsList(G4String physicsList) : G4VModularPhysicsList()
100{
101 G4LossTableManager::Instance();
102
103 defaultCutValue = 1.*mm;
104 fCutForGamma = defaultCutValue;
105 fCutForElectron = defaultCutValue;
106 fCutForPositron = defaultCutValue;
107 fDump = false;
108
109 fMessenger = new F04PhysicsListMessenger(this);
110
111 fEMPhysics = new PhysicsListVector();
112 fHadronPhysics = new PhysicsListVector();
113
114 // Particles
115 fParticleList = new G4DecayPhysics("decays");
116
117 // EM physics
118 AddPhysicsList("emstandard_opt1");
119
120 // Set the default hadronic physics.
121 AddPhysicsList(physicsList);
122
123 stepMaxProcess = new F04StepMax();
124}
125
126F04PhysicsList::~F04PhysicsList()
127{
128 delete fMessenger;
129 delete fParticleList;
130
131 ClearEMPhysics();
132 ClearHadronPhysics();
133
134 delete fEMPhysics;
135 delete fHadronPhysics;
136
137 delete stepMaxProcess;
138}
139
140void F04PhysicsList::ClearEMPhysics()
141{
142 for (PhysicsListVector::iterator p = fEMPhysics->begin();
143 p != fEMPhysics->end(); ++p) {
144 delete (*p);
145 }
146 fEMPhysics->clear();
147}
148
149void F04PhysicsList::ClearHadronPhysics()
150{
151 for (PhysicsListVector::iterator p = fHadronPhysics->begin();
152 p != fHadronPhysics->end(); ++p) {
153 delete (*p);
154 }
155 fHadronPhysics->clear();
156}
157
158void F04PhysicsList::ConstructParticle()
159{
160 fParticleList->ConstructParticle();
161
162 G4DecayTable* MuonPlusDecayTable = new G4DecayTable();
163 MuonPlusDecayTable -> Insert(new
164 G4MuonDecayChannelWithSpin("mu+",0.986));
165 MuonPlusDecayTable -> Insert(new
166 G4MuonRadiativeDecayChannelWithSpin("mu+",0.014));
167 G4MuonPlus::MuonPlusDefinition() -> SetDecayTable(MuonPlusDecayTable);
168
169 G4DecayTable* MuonMinusDecayTable = new G4DecayTable();
170 MuonMinusDecayTable -> Insert(new
171 G4MuonDecayChannelWithSpin("mu-",0.986));
172 MuonMinusDecayTable -> Insert(new
173 G4MuonRadiativeDecayChannelWithSpin("mu-",0.014));
174 G4MuonMinus::MuonMinusDefinition() -> SetDecayTable(MuonMinusDecayTable);
175
176}
177
178void F04PhysicsList::ConstructProcess()
179{
180 AddTransportation();
181
182 for (PhysicsListVector::iterator p = fEMPhysics->begin();
183 p != fEMPhysics->end(); ++p) {
184 (*p)->ConstructProcess();
185 }
186
187 fParticleList->ConstructProcess();
188
189 G4DecayWithSpin* decayWithSpin = new G4DecayWithSpin();
190
191 G4ProcessTable* processTable = G4ProcessTable::GetProcessTable();
192
193 G4VProcess* decay;
194 decay = processTable->FindProcess("Decay",G4MuonPlus::MuonPlus());
195
196 G4ProcessManager* fManager;
197 fManager = G4MuonPlus::MuonPlus()->GetProcessManager();
198
199 if (fManager) {
200 if (decay) fManager->RemoveProcess(decay);
201 fManager->AddProcess(decayWithSpin);
202 // set ordering for PostStepDoIt and AtRestDoIt
203 fManager ->SetProcessOrdering(decayWithSpin, idxPostStep);
204 fManager ->SetProcessOrdering(decayWithSpin, idxAtRest);
205 }
206
207 decay = processTable->FindProcess("Decay",G4MuonMinus::MuonMinus());
208
209 fManager = G4MuonMinus::MuonMinus()->GetProcessManager();
210
211 if (fManager) {
212 if (decay) fManager->RemoveProcess(decay);
213 fManager->AddProcess(decayWithSpin);
214 // set ordering for PostStepDoIt and AtRestDoIt
215 fManager ->SetProcessOrdering(decayWithSpin, idxPostStep);
216 fManager ->SetProcessOrdering(decayWithSpin, idxAtRest);
217 }
218
219 G4PionDecayMakeSpin* poldecay = new G4PionDecayMakeSpin();
220
221 decay = processTable->FindProcess("Decay",G4PionPlus::PionPlus());
222
223 fManager = G4PionPlus::PionPlus()->GetProcessManager();
224
225 if (fManager) {
226 if (decay) fManager->RemoveProcess(decay);
227 fManager->AddProcess(poldecay);
228 // set ordering for PostStepDoIt and AtRestDoIt
229 fManager ->SetProcessOrdering(poldecay, idxPostStep);
230 fManager ->SetProcessOrdering(poldecay, idxAtRest);
231 }
232
233 decay = processTable->FindProcess("Decay",G4PionMinus::PionMinus());
234
235 fManager = G4PionMinus::PionMinus()->GetProcessManager();
236
237 if (fManager) {
238 if (decay) fManager->RemoveProcess(decay);
239 fManager->AddProcess(poldecay);
240 // set ordering for PostStepDoIt and AtRestDoIt
241 fManager ->SetProcessOrdering(poldecay, idxPostStep);
242 fManager ->SetProcessOrdering(poldecay, idxAtRest);
243 }
244
245 for (PhysicsListVector::iterator p = fHadronPhysics->begin();
246 p != fHadronPhysics->end();++p) {
247 (*p)->ConstructProcess();
248 }
249
250 AddStepMax();
251
252 if (fDump) G4HadronProcessStore::Instance()->Dump(1);
253}
254
255void F04PhysicsList::RemoveFromEMPhysicsList(const G4String& name)
256{
257 G4bool success = false;
258 for (PhysicsListVector::iterator p = fEMPhysics->begin();
259 p != fEMPhysics->end(); ++p) {
260 G4VPhysicsConstructor* e = (*p);
261 if (e->GetPhysicsName() == name) {
262 fEMPhysics->erase(p);
263 success = true;
264 break;
265 }
266 }
267 if (!success) {
268 std::ostringstream message;
269 message << "PhysicsList::RemoveFromEMPhysicsList "<< name << "not found";
270 G4Exception(message.str().c_str());
271 }
272}
273
274void F04PhysicsList::RemoveFromHadronPhysicsList(const G4String& name)
275{
276 G4bool success = false;
277 for (PhysicsListVector::iterator p = fHadronPhysics->begin();
278 p != fHadronPhysics->end(); ++p) {
279 G4VPhysicsConstructor* e = (*p);
280 if (e->GetPhysicsName() == name) {
281 fHadronPhysics->erase(p);
282 success = true;
283 break;
284 }
285 }
286 if (!success) {
287 G4String message("F04PhysicsList::RemoveFromHadronPhysicsList \"");
288 message += name;
289 message += " not found \"";
290 G4Exception(message);
291 }
292}
293
294void F04PhysicsList::AddPhysicsList(const G4String& name)
295{
296 if (verboseLevel>0) {
297 G4cout << "F04PhysicsList::AddPhysicsList: <" << name
298 << ">" << G4endl;
299 }
300
301 if (name == "emstandard_opt1") {
302
303 ClearEMPhysics();
304 fEMPhysics->push_back(new G4EmStandardPhysics_option1());
305 fEMPhysics->push_back(new F04ExtraPhysics());
306// fEMPhysics->push_back(new F04OpticalPhysics());
307
308 } else if (name == "emstandard_opt2") {
309
310 ClearEMPhysics();
311 fEMPhysics->push_back(new G4EmStandardPhysics_option2());
312 fEMPhysics->push_back(new F04ExtraPhysics());
313// fEMPhysics->push_back(new F04OpticalPhysics());
314
315 } else if (name == "FTFC") {
316
317 SetStandardList(false, false);
318 fHadronPhysics->push_back(new HadronPhysicsFTFC("hadron",true));
319
320 } else if (name == "FTFP") {
321
322 SetStandardList(false, false);
323 fHadronPhysics->push_back(new HadronPhysicsFTFP("hadron",true));
324
325 } else if (name == "FTFP_EMV") {
326
327 AddPhysicsList("emstandard_opt1");
328 AddPhysicsList("FTFP");
329
330 } else if (name == "LHEP") {
331
332 ClearHadronPhysics();
333 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
334 fHadronPhysics->push_back(new G4HadronElasticPhysics("LElastic",
335 verboseLevel,
336 false));
337 fHadronPhysics->push_back(new HadronPhysicsLHEP("hadron"));
338 fHadronPhysics->push_back(new G4IonPhysics("ion"));
339 fDump = true;
340
341 } else if (name == "LHEP_BERT") {
342
343 ClearHadronPhysics();
344 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
345 fHadronPhysics->push_back(new G4HadronElasticPhysics("LElastic",
346 verboseLevel,
347 false));
348 fHadronPhysics->push_back(new HadronPhysicsLHEP_BERT("hadron"));
349 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
350 verboseLevel));
351 fHadronPhysics->push_back(new G4IonPhysics("ion"));
352 fDump = true;
353
354 } else if (name == "LHEP_EMV") {
355
356 ClearHadronPhysics();
357 AddPhysicsList("emstandard_opt1");
358 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
359 fHadronPhysics->push_back(new G4HadronElasticPhysics("LElastic",
360 verboseLevel,
361 false));
362 fHadronPhysics->push_back(new HadronPhysicsLHEP_EMV("hadron"));
363 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
364 verboseLevel));
365 fHadronPhysics->push_back(new G4IonPhysics("ion"));
366 fDump = true;
367
368 } else if (name == "LHEP_PRECO_HP") {
369
370 ClearHadronPhysics();
371 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
372 fHadronPhysics->push_back(new G4HadronElasticPhysics("LElastic",
373 verboseLevel,
374 false));
375 fHadronPhysics->push_back(new HadronPhysicsLHEP_PRECO_HP("hadron"));
376 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
377 verboseLevel));
378 fHadronPhysics->push_back(new G4IonPhysics("ion"));
379 fDump = true;
380
381 } else if (name == "QBBC") {
382
383 SetStandardList(false,true);
384 fHadronPhysics->push_back(new G4HadronInelasticQBBC("QBBC",
385 verboseLevel,
386 false,false,
387 false,false,true));
388
389 } else if (name == "QBBCG") {
390
391 SetStandardList(false, false);
392 fHadronPhysics->push_back(new G4HadronInelasticQBBC("QBBC",
393 verboseLevel,
394 false,false,
395 false,false,false));
396
397 } else if (name == "QBEC") {
398
399 SetStandardList(false,true);
400 fHadronPhysics->push_back(new G4HadronInelasticQBBC("QBBC",
401 verboseLevel,
402 false,true,
403 false,false,true));
404
405 } else if (name == "QBBC_HP") {
406
407 SetStandardList(true,true);
408 fHadronPhysics->push_back(new G4HadronInelasticQBBC("QBBC",
409 verboseLevel,
410 false,false,
411 false,true,true));
412
413 } else if (name == "QBEC_HP") {
414
415 SetStandardList(true,true);
416 fHadronPhysics->push_back(new G4HadronInelasticQBBC("QBBC",
417 verboseLevel,
418 false,true,
419 false,true,true));
420
421 } else if (name == "QGSC") {
422
423 ClearHadronPhysics();
424 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
425 fHadronPhysics->push_back(new G4HadronQElasticPhysics("QElastic",
426 verboseLevel));
427 fHadronPhysics->push_back(new HadronPhysicsQGSC("hadron",true));
428 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
429 verboseLevel,false));
430 fHadronPhysics->push_back(new G4IonPhysics("ion"));
431 fHadronPhysics->push_back(new G4NeutronTrackingCut());
432 fDump = true;
433
434 } else if (name == "QGSC_EFLOW") {
435
436 ClearHadronPhysics();
437 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
438 fHadronPhysics->push_back(new G4HadronQElasticPhysics("QElastic",
439 verboseLevel));
440 fHadronPhysics->push_back(new HadronPhysicsQGSC_EFLOW("hadron",true));
441 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
442 verboseLevel,false));
443 fHadronPhysics->push_back(new G4IonPhysics("ion"));
444 fHadronPhysics->push_back(new G4NeutronTrackingCut());
445 fDump = true;
446
447 } else if (name == "QGSC_EMV") {
448
449 AddPhysicsList("emstandard_opt1");
450 AddPhysicsList("QGSC");
451
452 } else if (name == "QGSP") {
453
454 SetStandardList(false, false);
455 fHadronPhysics->push_back(new HadronPhysicsQGSP("hadron",true));
456
457 } else if (name == "QGSP_BERT") {
458
459 SetStandardList(false, false);
460 fHadronPhysics->push_back(new HadronPhysicsQGSP_BERT("hadron",true));
461
462 } else if (name == "QGSP_BERT_EMV") {
463
464 AddPhysicsList("emstandard_opt1");
465 AddPhysicsList("QGSP_BERT");
466
467 } else if (name == "QGSP_BERT_HP") {
468
469 SetStandardList(true, false);
470 fHadronPhysics->push_back(new HadronPhysicsQGSP_BERT_HP("hadron",true));
471
472 } else if (name == "QGSP_BERT_NQE") {
473
474 SetStandardList(false, false);
475 fHadronPhysics->push_back(new HadronPhysicsQGSP_BERT("hadron",false));
476
477 } else if (name == "QGSP_BERT_TRV") {
478
479 SetStandardList(false, false);
480 fHadronPhysics->push_back(new HadronPhysicsQGSP_BERT_TRV("hadron",true));
481
482 } else if (name == "QGSP_BIC") {
483
484 SetStandardList(false, false);
485 fHadronPhysics->push_back(new HadronPhysicsQGSP_BIC("hadron",true));
486
487 } else if (name == "QGSP_BIC_HP") {
488
489 SetStandardList(true, false);
490 fHadronPhysics->push_back(new HadronPhysicsQGSP_BIC_HP("hadron",true));
491
492 } else if (name == "QGSP_EMV") {
493
494 AddPhysicsList("emstandard_opt1");
495 AddPhysicsList("QGSP");
496
497 } else if (name == "QGSP_EMV_NQE") {
498
499 AddPhysicsList("emstandard_opt1");
500 AddPhysicsList("QGSP_NQE");
501
502 } else if (name == "QGSP_EMX") {
503
504 AddPhysicsList("emstandard_opt2");
505 AddPhysicsList("QGSP");
506
507 } else if (name == "QGSP_NQE") {
508
509 SetStandardList(false, false);
510 fHadronPhysics->push_back(new HadronPhysicsQGSP("hadron",false));
511
512 } else if (name == "QGSP_QEL") {
513
514 ClearHadronPhysics();
515 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
516 fHadronPhysics->push_back(new G4HadronQElasticPhysics("QElastic",
517 verboseLevel));
518 fHadronPhysics->push_back(new HadronPhysicsQGSP("hadron",true));
519 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
520 verboseLevel));
521 fHadronPhysics->push_back(new G4IonPhysics("ion"));
522 fHadronPhysics->push_back(new G4NeutronTrackingCut());
523 fDump = true;
524
525 } else {
526
527 G4cout << "F04PhysicsList::AddPhysicsList: <" << name << ">"
528 << " is not defined" << G4endl;
529 }
530}
531
532void F04PhysicsList::SetStandardList(G4bool flagHP, G4bool glauber)
533{
534 ClearHadronPhysics();
535
536 fHadronPhysics->push_back(new G4EmExtraPhysics("gamma_nuc"));
537 fHadronPhysics->push_back(new G4HadronElasticPhysics("elastic",
538 verboseLevel,
539 flagHP,glauber));
540 fHadronPhysics->push_back(new G4QStoppingPhysics("stopping",
541 verboseLevel));
542 fHadronPhysics->push_back(new G4IonBinaryCascadePhysics("binary_ion"));
543 fHadronPhysics->push_back(new G4NeutronTrackingCut());
544
545 fDump = true;
546}
547
548void F04PhysicsList::SetCuts()
549{
550 if (verboseLevel >0) {
551 G4cout << "F04PhysicsList::SetCuts:";
552 G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length")
553 << G4endl;
554 }
555
556 // set cut values for gamma at first and for e- second and next for e+,
557 // because some processes for e+/e- need cut values for gamma
558 SetCutValue(fCutForGamma, "gamma");
559 SetCutValue(fCutForElectron, "e-");
560 SetCutValue(fCutForPositron, "e+");
561
562 if (verboseLevel>0) DumpCutValuesTable();
563}
564
565void F04PhysicsList::SetCutForGamma(G4double cut)
566{
567 fCutForGamma = cut;
568 SetParticleCuts(fCutForGamma, G4Gamma::Gamma());
569}
570
571void F04PhysicsList::SetCutForElectron(G4double cut)
572{
573 fCutForElectron = cut;
574 SetParticleCuts(fCutForElectron, G4Electron::Electron());
575}
576
577void F04PhysicsList::SetCutForPositron(G4double cut)
578{
579 fCutForPositron = cut;
580 SetParticleCuts(fCutForPositron, G4Positron::Positron());
581}
582
583void F04PhysicsList::List()
584{
585 G4cout << "### PhysicsLists available: FTFC FTFP FTFP_EMV"
586 << " LHEP LHEP_BERT LHEP_EMV"
587 << G4endl;
588 G4cout << " LHEP_PRECO_HP QBBC QBBCG"
589 << " QBEC QBBC_HP QBEC_HP"
590 << G4endl;
591 G4cout << " QGSC QGSC_EFLOW QGSC_EMV"
592 << " QGSP QGSP_BERT QGSP_BER_EMV"
593 << G4endl;
594 G4cout << " QGSP_BERT_HP QGSP_BERT_NQE"
595 << " QGSP_BERT_TRV"
596 << G4endl;
597 G4cout << " QGSP_BIC QGSP_BIC_HP QGSP_EMV"
598 << " QGSP_EMV_NQE"
599 << G4endl;
600 G4cout << " QGSC_EMX QGSP_NQE QGSP_QEL"
601 << G4endl;
602}
603
604void F04PhysicsList::SetStepMax(G4double step)
605{
606 MaxChargedStep = step ;
607 stepMaxProcess->SetStepMax(MaxChargedStep);
608}
609
610F04StepMax* F04PhysicsList::GetStepMaxProcess()
611{
612 return stepMaxProcess;
613}
614
615void F04PhysicsList::AddStepMax()
616{
617 // Step limitation seen as a process
618
619 theParticleIterator->reset();
620 while ((*theParticleIterator)()){
621 G4ParticleDefinition* particle = theParticleIterator->value();
622 G4ProcessManager* pmanager = particle->GetProcessManager();
623
624 if (stepMaxProcess->IsApplicable(*particle) && !particle->IsShortLived())
625 {
626 if (pmanager) pmanager ->AddDiscreteProcess(stepMaxProcess);
627 }
628 }
629}
Note: See TracBrowser for help on using the repository browser.