| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 /// \file field/field04/src/F04PhysicsList.cc <<
28 /// \brief Implementation of the F04PhysicsLis <<
29 // 27 //
30 28
31 #include "F04PhysicsList.hh" 29 #include "F04PhysicsList.hh"
32 <<
33 #include "F04PhysicsListMessenger.hh" 30 #include "F04PhysicsListMessenger.hh"
34 31
35 #include "G4LossTableManager.hh" << 32 #include "F04ExtraPhysics.hh"
36 #include "G4OpticalPhysics.hh" 33 #include "G4OpticalPhysics.hh"
37 #include "G4ParticleTable.hh" << 34
38 #include "G4ParticleTypes.hh" << 35 #include "G4LossTableManager.hh"
>> 36
39 #include "G4ProcessManager.hh" 37 #include "G4ProcessManager.hh"
40 #include "G4StepLimiterPhysics.hh" << 38 #include "G4ParticleTypes.hh"
>> 39 #include "G4ParticleTable.hh"
41 40
42 #include "F04StepMax.hh" << 41 #include "G4PhysListFactory.hh"
43 #include "FTFP_BERT.hh" <<
44 #include "QGSP_BERT.hh" <<
45 42
46 #include "G4DecayTable.hh" <<
47 #include "G4DecayWithSpin.hh" <<
48 #include "G4Electron.hh" <<
49 #include "G4Gamma.hh" 43 #include "G4Gamma.hh"
50 #include "G4MuMinusCapturePrecompound.hh" << 44 #include "G4Electron.hh"
51 #include "G4MuonDecayChannelWithSpin.hh" <<
52 #include "G4MuonMinusCapture.hh" <<
53 #include "G4MuonRadiativeDecayChannelWithSpin. <<
54 #include "G4PionDecayMakeSpin.hh" <<
55 #include "G4Positron.hh" 45 #include "G4Positron.hh"
>> 46
>> 47 #include "F04StepMax.hh"
>> 48
56 #include "G4ProcessTable.hh" 49 #include "G4ProcessTable.hh"
57 #include "G4SystemOfUnits.hh" <<
58 50
59 G4ThreadLocal F04StepMax* F04PhysicsList::fSte << 51 #include "G4PionDecayMakeSpin.hh"
>> 52 #include "G4DecayWithSpin.hh"
60 53
61 //....oooOO0OOooo........oooOO0OOooo........oo << 54 #include "G4DecayTable.hh"
>> 55 #include "G4MuonDecayChannelWithSpin.hh"
>> 56 #include "G4MuonRadiativeDecayChannelWithSpin.hh"
62 57
63 F04PhysicsList::F04PhysicsList(const G4String& << 58 F04PhysicsList::F04PhysicsList(G4String physName) : G4VModularPhysicsList()
64 { 59 {
65 G4LossTableManager::Instance(); << 60 G4LossTableManager::Instance();
66 61
67 defaultCutValue = 1. * mm; << 62 defaultCutValue = 1.*mm;
>> 63 fCutForGamma = defaultCutValue;
>> 64 fCutForElectron = defaultCutValue;
>> 65 fCutForPositron = defaultCutValue;
68 66
69 fMessenger = new F04PhysicsListMessenger(thi << 67 fMessenger = new F04PhysicsListMessenger(this);
70 68
71 SetVerboseLevel(1); << 69 SetVerboseLevel(1);
72 70
73 // G4PhysListFactory factory; << 71 G4PhysListFactory factory;
74 G4VModularPhysicsList* phys = nullptr; << 72 G4VModularPhysicsList* phys = NULL;
75 if (physName == "QGSP_BERT") { <<
76 phys = new QGSP_BERT; <<
77 } <<
78 else { <<
79 phys = new FTFP_BERT; <<
80 } <<
81 73
82 // if (factory.IsReferencePhysList(physNa << 74 if (factory.IsReferencePhysList(physName))
83 // phys =factory.GetReferencePhysList( << 75 phys =factory.GetReferencePhysList(physName);
84 76
85 // Physics List is defined via environment v << 77 // Physics List is defined via environment variable PHYSLIST
86 // if (!phys) phys = factory.ReferencePhy << 78 if (!phys) phys = factory.ReferencePhysList();
87 79
88 if (!phys) << 80 if (!phys) G4Exception("WLSPhysicsList::WLSPhysicsList","InvalidSetup",
89 G4Exception("F04PhysicsList::F04PhysicsLis << 81 FatalException,"PhysicsList does not exist");
90 "PhysicsList does not exist"); <<
91 82
92 for (G4int i = 0;; ++i) { << 83 for (G4int i = 0; ; ++i) {
93 auto elem = const_cast<G4VPhysicsConstruct << 84 G4VPhysicsConstructor* elem =
94 if (elem == nullptr) break; << 85 const_cast<G4VPhysicsConstructor*> (phys->GetPhysics(i));
95 G4cout << "RegisterPhysics: " << elem->Get << 86 if (elem == NULL) break;
96 RegisterPhysics(elem); << 87 G4cout << "RegisterPhysics: " << elem->GetPhysicsName() << G4endl;
97 } << 88 RegisterPhysics(elem);
>> 89 }
98 90
99 RegisterPhysics(new G4StepLimiterPhysics()); << 91 RegisterPhysics(new F04ExtraPhysics());
100 RegisterPhysics(new G4OpticalPhysics()); << 92 RegisterPhysics(new G4OpticalPhysics());
101 } <<
102 93
103 //....oooOO0OOooo........oooOO0OOooo........oo << 94 stepMaxProcess = new F04StepMax();
>> 95 }
104 96
105 F04PhysicsList::~F04PhysicsList() 97 F04PhysicsList::~F04PhysicsList()
106 { 98 {
107 delete fMessenger; << 99 delete fMessenger;
108 } <<
109 100
110 //....oooOO0OOooo........oooOO0OOooo........oo << 101 delete stepMaxProcess;
>> 102 }
111 103
112 void F04PhysicsList::ConstructParticle() 104 void F04PhysicsList::ConstructParticle()
113 { 105 {
114 G4VModularPhysicsList::ConstructParticle(); << 106 G4VModularPhysicsList::ConstructParticle();
115 <<
116 G4GenericIon::GenericIonDefinition(); <<
117 107
118 auto muonPlusDecayTable = new G4DecayTable() << 108 G4DecayTable* MuonPlusDecayTable = new G4DecayTable();
119 muonPlusDecayTable->Insert(new G4MuonDecayCh << 109 MuonPlusDecayTable -> Insert(new
120 muonPlusDecayTable->Insert(new G4MuonRadiati << 110 G4MuonDecayChannelWithSpin("mu+",0.986));
121 G4MuonPlus::MuonPlusDefinition()->SetDecayTa << 111 MuonPlusDecayTable -> Insert(new
122 << 112 G4MuonRadiativeDecayChannelWithSpin("mu+",0.014));
123 auto muonMinusDecayTable = new G4DecayTable( << 113 G4MuonPlus::MuonPlusDefinition() -> SetDecayTable(MuonPlusDecayTable);
124 muonMinusDecayTable->Insert(new G4MuonDecayC << 114
125 muonMinusDecayTable->Insert(new G4MuonRadiat << 115 G4DecayTable* MuonMinusDecayTable = new G4DecayTable();
126 G4MuonMinus::MuonMinusDefinition()->SetDecay << 116 MuonMinusDecayTable -> Insert(new
>> 117 G4MuonDecayChannelWithSpin("mu-",0.986));
>> 118 MuonMinusDecayTable -> Insert(new
>> 119 G4MuonRadiativeDecayChannelWithSpin("mu-",0.014));
>> 120 G4MuonMinus::MuonMinusDefinition() -> SetDecayTable(MuonMinusDecayTable);
127 } 121 }
128 122
129 //....oooOO0OOooo........oooOO0OOooo........oo <<
130 <<
131 void F04PhysicsList::ConstructProcess() 123 void F04PhysicsList::ConstructProcess()
132 { 124 {
133 G4VModularPhysicsList::ConstructProcess(); << 125 G4VModularPhysicsList::ConstructProcess();
134 126
135 fStepMaxProcess = new F04StepMax(); << 127 G4DecayWithSpin* decayWithSpin = new G4DecayWithSpin();
136 128
137 auto decayWithSpin = new G4DecayWithSpin(); << 129 G4ProcessTable* processTable = G4ProcessTable::GetProcessTable();
138 130
139 G4ProcessTable* processTable = G4ProcessTabl << 131 G4VProcess* decay;
>> 132 decay = processTable->FindProcess("Decay",G4MuonPlus::MuonPlus());
140 133
141 G4VProcess* decay; << 134 G4ProcessManager* fManager;
142 decay = processTable->FindProcess("Decay", G << 135 fManager = G4MuonPlus::MuonPlus()->GetProcessManager();
143 136
144 G4ProcessManager* pmanager; << 137 if (fManager) {
145 pmanager = G4MuonPlus::MuonPlus()->GetProces << 138 if (decay) fManager->RemoveProcess(decay);
146 << 139 fManager->AddProcess(decayWithSpin);
147 if (pmanager) { << 140 // set ordering for PostStepDoIt and AtRestDoIt
148 if (decay) pmanager->RemoveProcess(decay); << 141 fManager ->SetProcessOrdering(decayWithSpin, idxPostStep);
149 pmanager->AddProcess(decayWithSpin); << 142 fManager ->SetProcessOrdering(decayWithSpin, idxAtRest);
150 // set ordering for PostStepDoIt and AtRes << 143 }
151 pmanager->SetProcessOrdering(decayWithSpin <<
152 pmanager->SetProcessOrdering(decayWithSpin <<
153 } <<
154 <<
155 decay = processTable->FindProcess("Decay", G <<
156 <<
157 pmanager = G4MuonMinus::MuonMinus()->GetProc <<
158 144
159 if (pmanager) { << 145 decay = processTable->FindProcess("Decay",G4MuonMinus::MuonMinus());
160 if (decay) pmanager->RemoveProcess(decay); <<
161 pmanager->AddProcess(decayWithSpin); <<
162 // set ordering for PostStepDoIt and AtRes <<
163 pmanager->SetProcessOrdering(decayWithSpin <<
164 pmanager->SetProcessOrdering(decayWithSpin <<
165 } <<
166 146
167 G4VProcess* process = processTable->FindProc << 147 fManager = G4MuonMinus::MuonMinus()->GetProcessManager();
168 148
169 if (pmanager) { << 149 if (fManager) {
170 if (process) pmanager->RemoveProcess(proce << 150 if (decay) fManager->RemoveProcess(decay);
171 process = new G4MuonMinusCapture(new G4MuM << 151 fManager->AddProcess(decayWithSpin);
172 pmanager->AddRestProcess(process); << 152 // set ordering for PostStepDoIt and AtRestDoIt
173 } << 153 fManager ->SetProcessOrdering(decayWithSpin, idxPostStep);
>> 154 fManager ->SetProcessOrdering(decayWithSpin, idxAtRest);
>> 155 }
174 156
175 auto poldecay = new G4PionDecayMakeSpin(); << 157 G4PionDecayMakeSpin* poldecay = new G4PionDecayMakeSpin();
176 158
177 decay = processTable->FindProcess("Decay", G << 159 decay = processTable->FindProcess("Decay",G4PionPlus::PionPlus());
178 160
179 pmanager = G4PionPlus::PionPlus()->GetProces << 161 fManager = G4PionPlus::PionPlus()->GetProcessManager();
180 162
181 if (pmanager) { << 163 if (fManager) {
182 if (decay) pmanager->RemoveProcess(decay); << 164 if (decay) fManager->RemoveProcess(decay);
183 pmanager->AddProcess(poldecay); << 165 fManager->AddProcess(poldecay);
184 // set ordering for PostStepDoIt and AtRes << 166 // set ordering for PostStepDoIt and AtRestDoIt
185 pmanager->SetProcessOrdering(poldecay, idx << 167 fManager ->SetProcessOrdering(poldecay, idxPostStep);
186 pmanager->SetProcessOrdering(poldecay, idx << 168 fManager ->SetProcessOrdering(poldecay, idxAtRest);
187 } << 169 }
188 170
189 decay = processTable->FindProcess("Decay", G << 171 decay = processTable->FindProcess("Decay",G4PionMinus::PionMinus());
190 172
191 pmanager = G4PionMinus::PionMinus()->GetProc << 173 fManager = G4PionMinus::PionMinus()->GetProcessManager();
192 174
193 if (pmanager) { << 175 if (fManager) {
194 if (decay) pmanager->RemoveProcess(decay); << 176 if (decay) fManager->RemoveProcess(decay);
195 pmanager->AddProcess(poldecay); << 177 fManager->AddProcess(poldecay);
196 // set ordering for PostStepDoIt and AtRes << 178 // set ordering for PostStepDoIt and AtRestDoIt
197 pmanager->SetProcessOrdering(poldecay, idx << 179 fManager ->SetProcessOrdering(poldecay, idxPostStep);
198 pmanager->SetProcessOrdering(poldecay, idx << 180 fManager ->SetProcessOrdering(poldecay, idxAtRest);
199 } << 181 }
200 182
201 AddStepMax(); << 183 AddStepMax();
202 } 184 }
203 185
204 /* 186 /*
205 //....oooOO0OOooo........oooOO0OOooo........oo <<
206 <<
207 void F04PhysicsList::RemoveFromPhysicsList(con 187 void F04PhysicsList::RemoveFromPhysicsList(const G4String& name)
208 { 188 {
209 G4bool success = false; 189 G4bool success = false;
210 for (G4PhysConstVector::iterator p = phys 190 for (G4PhysConstVector::iterator p = physicsVector->begin();
211 p != phys 191 p != physicsVector->end(); ++p) {
212 G4VPhysicsConstructor* e = (*p); 192 G4VPhysicsConstructor* e = (*p);
213 if (e->GetPhysicsName() == name) { 193 if (e->GetPhysicsName() == name) {
214 physicsVector->erase(p); 194 physicsVector->erase(p);
215 success = true; 195 success = true;
216 break; 196 break;
217 } 197 }
218 } 198 }
219 if (!success) { 199 if (!success) {
220 std::ostringstream message; 200 std::ostringstream message;
221 message << "PhysicsList::RemoveFromPhys 201 message << "PhysicsList::RemoveFromPhysicsList "<< name << "not found";
222 G4Exception(message.str().c_str()); 202 G4Exception(message.str().c_str());
223 } 203 }
224 } 204 }
225 205
226 //....oooOO0OOooo........oooOO0OOooo........oo <<
227 <<
228 void F04PhysicsList::ClearPhysics() 206 void F04PhysicsList::ClearPhysics()
229 { 207 {
230 for (G4PhysConstVector::iterator p = phys 208 for (G4PhysConstVector::iterator p = physicsVector->begin();
231 p != phys 209 p != physicsVector->end(); ++p) {
232 delete (*p); 210 delete (*p);
233 } 211 }
234 physicsVector->clear(); 212 physicsVector->clear();
235 } 213 }
236 */ 214 */
237 215
238 //....oooOO0OOooo........oooOO0OOooo........oo << 216 void F04PhysicsList::SetCuts()
>> 217 {
>> 218 if (verboseLevel >0) {
>> 219 G4cout << "F04PhysicsList::SetCuts:";
>> 220 G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length")
>> 221 << G4endl;
>> 222 }
239 223
240 void F04PhysicsList::SetStepMax(G4double step) << 224 // set cut values for gamma at first and for e- second and next for e+,
>> 225 // because some processes for e+/e- need cut values for gamma
>> 226 SetCutValue(fCutForGamma, "gamma");
>> 227 SetCutValue(fCutForElectron, "e-");
>> 228 SetCutValue(fCutForPositron, "e+");
>> 229
>> 230 if (verboseLevel>0) DumpCutValuesTable();
>> 231 }
>> 232
>> 233 void F04PhysicsList::SetCutForGamma(G4double cut)
>> 234 {
>> 235 fCutForGamma = cut;
>> 236 SetParticleCuts(fCutForGamma, G4Gamma::Gamma());
>> 237 }
>> 238
>> 239 void F04PhysicsList::SetCutForElectron(G4double cut)
241 { 240 {
242 fMaxChargedStep = step; << 241 fCutForElectron = cut;
243 fStepMaxProcess->SetStepMax(fMaxChargedStep) << 242 SetParticleCuts(fCutForElectron, G4Electron::Electron());
244 } 243 }
245 244
246 //....oooOO0OOooo........oooOO0OOooo........oo << 245 void F04PhysicsList::SetCutForPositron(G4double cut)
>> 246 {
>> 247 fCutForPositron = cut;
>> 248 SetParticleCuts(fCutForPositron, G4Positron::Positron());
>> 249 }
247 250
248 F04StepMax* F04PhysicsList::GetStepMaxProcess( << 251 void F04PhysicsList::SetStepMax(G4double step)
249 { 252 {
250 return fStepMaxProcess; << 253 MaxChargedStep = step ;
>> 254 stepMaxProcess->SetStepMax(MaxChargedStep);
251 } 255 }
252 256
253 //....oooOO0OOooo........oooOO0OOooo........oo << 257 F04StepMax* F04PhysicsList::GetStepMaxProcess()
>> 258 {
>> 259 return stepMaxProcess;
>> 260 }
254 261
255 void F04PhysicsList::AddStepMax() 262 void F04PhysicsList::AddStepMax()
256 { 263 {
257 // Step limitation seen as a process 264 // Step limitation seen as a process
258 265
259 auto particleIterator = GetParticleIterator( << 266 theParticleIterator->reset();
260 particleIterator->reset(); << 267 while ((*theParticleIterator)()){
261 while ((*particleIterator)()) { << 268 G4ParticleDefinition* particle = theParticleIterator->value();
262 G4ParticleDefinition* particle = particleI << 269 G4ProcessManager* pmanager = particle->GetProcessManager();
263 G4ProcessManager* pmanager = particle->Get << 270
264 << 271 if (stepMaxProcess->IsApplicable(*particle) && !particle->IsShortLived())
265 if (fStepMaxProcess->IsApplicable(*particl << 272 {
266 if (pmanager) pmanager->AddDiscreteProce << 273 if (pmanager) pmanager ->AddDiscreteProcess(stepMaxProcess);
267 } << 274 }
268 } 275 }
269 } 276 }
270 277