source: trunk/source/processes/hadronic/management/src/G4HadronicProcessStore.cc@ 1036

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

update to geant4.9.2
File size: 18.3 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// $Id: G4HadronicProcessStore.cc,v 1.7 2008/10/22 07:58:20 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-02 $
28//
29// -------------------------------------------------------------------
30//
31// GEANT4 Class file
32//
33//
34// File name: G4HadronicProcessStore
35//
36// Author: Vladimir Ivanchenko
37//
38// Creation date: 09.05.2008
39//
40// Modifications:
41//
42//
43// Class Description:
44//
45// -------------------------------------------------------------------
46//
47//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
48//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
49
50#include "G4HadronicProcessStore.hh"
51#include "G4Element.hh"
52#include "G4ProcessManager.hh"
53#include "G4Electron.hh"
54#include "G4Proton.hh"
55
56G4HadronicProcessStore* G4HadronicProcessStore::theInstance = 0;
57
58//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
59
60G4HadronicProcessStore* G4HadronicProcessStore::Instance()
61{
62 if(0 == theInstance) {
63 static G4HadronicProcessStore manager;
64 theInstance = &manager;
65 }
66 return theInstance;
67}
68
69//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
70
71G4HadronicProcessStore::~G4HadronicProcessStore()
72{
73 /*
74 for (G4int i=0; i<n_proc; i++) {
75 if( process[i] ) delete process[i];
76 }
77 */
78}
79
80//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
81
82G4HadronicProcessStore::G4HadronicProcessStore()
83{
84 n_proc = 0;
85 n_part = 0;
86 n_model= 0;
87 n_extra= 0;
88 currentProcess = 0;
89 currentParticle = 0;
90 verbose = 1;
91 buildTableStart = true;
92}
93
94//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
95
96G4double G4HadronicProcessStore::GetElasticCrossSectionPerVolume(
97 const G4ParticleDefinition *aParticle,
98 G4double kineticEnergy,
99 const G4Material *material)
100{
101 G4double cross = 0.0;
102 const G4ElementVector* theElementVector = material->GetElementVector();
103 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
104 size_t nelm = material->GetNumberOfElements();
105 for (size_t i=0; i<nelm; i++) {
106 const G4Element* elm = (*theElementVector)[i];
107 cross += theAtomNumDensityVector[i]*
108 GetElasticCrossSectionPerAtom(aParticle,kineticEnergy,elm);
109 }
110 return cross;
111}
112
113//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
114
115G4double G4HadronicProcessStore::GetElasticCrossSectionPerAtom(
116 const G4ParticleDefinition *aParticle,
117 G4double kineticEnergy,
118 const G4Element *anElement)
119{
120 G4HadronicProcess* hp = FindProcess(aParticle, fHadronElastic);
121 localDP.SetKineticEnergy(kineticEnergy);
122 G4double cross = 0.0;
123 if(hp) cross = hp->GetMicroscopicCrossSection(&localDP,
124 anElement,
125 STP_Temperature);
126 return cross;
127}
128
129//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
130
131G4double G4HadronicProcessStore::GetElasticCrossSectionPerIsotope(
132 const G4ParticleDefinition*,
133 G4double,
134 G4int, G4int)
135{
136 return 0.0;
137}
138
139//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
140
141G4double G4HadronicProcessStore::GetInelasticCrossSectionPerVolume(
142 const G4ParticleDefinition *aParticle,
143 G4double kineticEnergy,
144 const G4Material *material)
145{
146 G4double cross = 0.0;
147 const G4ElementVector* theElementVector = material->GetElementVector();
148 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
149 size_t nelm = material->GetNumberOfElements();
150 for (size_t i=0; i<nelm; i++) {
151 const G4Element* elm = (*theElementVector)[i];
152 cross += theAtomNumDensityVector[i]*
153 GetInelasticCrossSectionPerAtom(aParticle,kineticEnergy,elm);
154 }
155 return cross;
156}
157
158//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
159
160G4double G4HadronicProcessStore::GetInelasticCrossSectionPerAtom(
161 const G4ParticleDefinition *aParticle,
162 G4double kineticEnergy,
163 const G4Element *anElement)
164{
165 G4HadronicProcess* hp = FindProcess(aParticle, fHadronInelastic);
166 localDP.SetDefinition(const_cast<G4ParticleDefinition*>(aParticle));
167 localDP.SetKineticEnergy(kineticEnergy);
168 G4double cross = 0.0;
169 if(hp) cross = hp->GetMicroscopicCrossSection(&localDP,
170 anElement,
171 STP_Temperature);
172 return cross;
173}
174
175//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
176
177G4double G4HadronicProcessStore::GetInelasticCrossSectionPerIsotope(
178 const G4ParticleDefinition *,
179 G4double,
180 G4int, G4int)
181{
182 return 0.0;
183}
184
185//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
186
187G4double G4HadronicProcessStore::GetCaptureCrossSectionPerVolume(
188 const G4ParticleDefinition *aParticle,
189 G4double kineticEnergy,
190 const G4Material *material)
191{
192 G4double cross = 0.0;
193 const G4ElementVector* theElementVector = material->GetElementVector();
194 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
195 size_t nelm = material->GetNumberOfElements();
196 for (size_t i=0; i<nelm; i++) {
197 const G4Element* elm = (*theElementVector)[i];
198 cross += theAtomNumDensityVector[i]*
199 GetCaptureCrossSectionPerAtom(aParticle,kineticEnergy,elm);
200 }
201 return cross;
202}
203
204//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
205
206G4double G4HadronicProcessStore::GetCaptureCrossSectionPerAtom(
207 const G4ParticleDefinition *aParticle,
208 G4double kineticEnergy,
209 const G4Element *anElement)
210{
211 G4HadronicProcess* hp = FindProcess(aParticle, fCapture);
212 localDP.SetDefinition(const_cast<G4ParticleDefinition*>(aParticle));
213 localDP.SetKineticEnergy(kineticEnergy);
214 G4double cross = 0.0;
215 if(hp) cross = hp->GetMicroscopicCrossSection(&localDP,
216 anElement,
217 STP_Temperature);
218 return cross;
219}
220
221//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
222
223G4double G4HadronicProcessStore::GetCaptureCrossSectionPerIsotope(
224 const G4ParticleDefinition *,
225 G4double,
226 G4int, G4int)
227{
228 return 0.0;
229}
230
231//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
232
233G4double G4HadronicProcessStore::GetFissionCrossSectionPerVolume(
234 const G4ParticleDefinition *aParticle,
235 G4double kineticEnergy,
236 const G4Material *material)
237{
238 G4double cross = 0.0;
239 const G4ElementVector* theElementVector = material->GetElementVector();
240 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
241 size_t nelm = material->GetNumberOfElements();
242 for (size_t i=0; i<nelm; i++) {
243 const G4Element* elm = (*theElementVector)[i];
244 cross += theAtomNumDensityVector[i]*
245 GetFissionCrossSectionPerAtom(aParticle,kineticEnergy,elm);
246 }
247 return cross;
248}
249
250//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
251
252G4double G4HadronicProcessStore::GetFissionCrossSectionPerAtom(
253 const G4ParticleDefinition *aParticle,
254 G4double kineticEnergy,
255 const G4Element *anElement)
256{
257 G4HadronicProcess* hp = FindProcess(aParticle, fFission);
258 localDP.SetDefinition(const_cast<G4ParticleDefinition*>(aParticle));
259 localDP.SetKineticEnergy(kineticEnergy);
260 G4double cross = 0.0;
261 if(hp) cross = hp->GetMicroscopicCrossSection(&localDP,
262 anElement,
263 STP_Temperature);
264 return cross;
265}
266
267//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
268
269G4double G4HadronicProcessStore::GetFissionCrossSectionPerIsotope(
270 const G4ParticleDefinition *,
271 G4double,
272 G4int, G4int)
273{
274 return 0.0;
275}
276
277//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
278
279G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerVolume(
280 const G4ParticleDefinition *aParticle,
281 G4double kineticEnergy,
282 const G4Material *material)
283{
284 G4double cross = 0.0;
285 const G4ElementVector* theElementVector = material->GetElementVector();
286 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
287 size_t nelm = material->GetNumberOfElements();
288 for (size_t i=0; i<nelm; i++) {
289 const G4Element* elm = (*theElementVector)[i];
290 cross += theAtomNumDensityVector[i]*
291 GetChargeExchangeCrossSectionPerAtom(aParticle,kineticEnergy,elm);
292 }
293 return cross;
294}
295
296//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
297
298G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerAtom(
299 const G4ParticleDefinition *aParticle,
300 G4double kineticEnergy,
301 const G4Element *anElement)
302{
303 G4HadronicProcess* hp = FindProcess(aParticle, fChargeExchange);
304 localDP.SetDefinition(const_cast<G4ParticleDefinition*>(aParticle));
305 localDP.SetKineticEnergy(kineticEnergy);
306 G4double cross = 0.0;
307 if(hp) cross = hp->GetMicroscopicCrossSection(&localDP,
308 anElement,
309 STP_Temperature);
310 return cross;
311}
312
313//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
314
315G4double G4HadronicProcessStore::GetChargeExchangeCrossSectionPerIsotope(
316 const G4ParticleDefinition *,
317 G4double,
318 G4int, G4int)
319{
320 return 0.0;
321}
322
323//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
324
325void G4HadronicProcessStore::Register(G4HadronicProcess* proc)
326{
327 for(G4int i=0; i<n_proc; i++) {if(process[i] == proc) return;}
328
329 n_proc++;
330 process.push_back(proc);
331}
332
333//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
334
335void G4HadronicProcessStore::RegisterParticle(G4HadronicProcess* proc,
336 const G4ParticleDefinition* part)
337{
338 G4int i=0;
339 for(; i<n_proc; i++) {if(process[i] == proc) break;}
340 G4int j=0;
341 for(; j<n_part; j++) {if(particle[j] == part) break;}
342
343 if(j == n_part) {
344 n_part++;
345 particle.push_back(part);
346 wasPrinted.push_back(0);
347 }
348
349 // the pair should be added?
350 if(i < n_proc) {
351 std::multimap<PD,HP,std::less<PD> >::iterator it;
352 for(it=p_map.lower_bound(part); it!=p_map.upper_bound(part); ++it) {
353 if(it->first == part) {
354 HP process = (it->second);
355 if(proc == process) return;
356 }
357 }
358 }
359
360 p_map.insert(std::multimap<PD,HP>::value_type(part,proc));
361}
362
363//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
364
365void G4HadronicProcessStore::RegisterInteraction(G4HadronicProcess* proc,
366 G4HadronicInteraction* mod)
367{
368 G4int i=0;
369 for(; i<n_proc; i++) {if(process[i] == proc) break;}
370 G4int k=0;
371 for(; k<n_model; k++) {if(model[k] == mod) break;}
372
373 m_map.insert(std::multimap<HP,HI>::value_type(proc,mod));
374
375 if(k == n_model) {
376 n_model++;
377 model.push_back(mod);
378 modelName.push_back(mod->GetModelName());
379 }
380}
381
382//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
383
384void G4HadronicProcessStore::DeRegister(G4HadronicProcess* proc)
385{
386 for(G4int i=0; i<n_proc; i++) {
387 if(process[i] == proc) {
388 process[i] = 0;
389 break;
390 }
391 }
392}
393
394//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
395
396void G4HadronicProcessStore::RegisterExtraProcess(G4VProcess* proc)
397{
398 for(G4int i=0; i<n_extra; i++) {if(extraProcess[i] == proc) return;}
399
400 n_extra++;
401 extraProcess.push_back(proc);
402}
403
404//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
405
406void G4HadronicProcessStore::RegisterParticleForExtraProcess(
407 G4VProcess* proc,
408 const G4ParticleDefinition* part)
409{
410 G4int i=0;
411 for(; i<n_extra; i++) {if(extraProcess[i] == proc) break;}
412 G4int j=0;
413 for(; j<n_part; j++) {if(particle[j] == part) break;}
414
415 if(j == n_part) {
416 n_part++;
417 particle.push_back(part);
418 wasPrinted.push_back(0);
419 }
420
421 // the pair should be added?
422 if(i < n_extra) {
423 std::multimap<PD,G4VProcess*,std::less<PD> >::iterator it;
424 for(it=ep_map.lower_bound(part); it!=ep_map.upper_bound(part); ++it) {
425 if(it->first == part) {
426 G4VProcess* process = (it->second);
427 if(proc == process) return;
428 }
429 }
430 }
431
432 ep_map.insert(std::multimap<PD,G4VProcess*>::value_type(part,proc));
433}
434
435//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
436
437void G4HadronicProcessStore::DeRegisterExtraProcess(G4VProcess* proc)
438{
439 for(G4int i=0; i<n_extra; i++) {
440 if(extraProcess[i] == proc) {
441 extraProcess[i] = 0;
442 break;
443 }
444 }
445}
446
447//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
448
449void G4HadronicProcessStore::PrintInfo(const G4ParticleDefinition* part)
450{
451 if(buildTableStart && part == particle[n_part - 1]) {
452 buildTableStart = false;
453 Dump(verbose);
454 }
455}
456
457//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
458
459void G4HadronicProcessStore::Dump(G4int level)
460{
461 if(level > 0) {
462 G4cout << "=============================================================="
463 << "=============================="
464 << G4endl;
465 G4cout << " HADRONIC PROCESSES SUMMARY (verbose level " << level
466 << ")" << G4endl;
467 }
468 for(G4int i=0; i<n_part; i++) {
469 PD part = particle[i];
470 G4String pname = part->GetParticleName();
471 G4bool yes = false;
472 if(level >= 2) yes = true;
473 else if(level == 1 && (pname == "proton" ||
474 pname == "neutron" ||
475 pname == "pi+" ||
476 pname == "pi-" ||
477 pname == "gamma" ||
478 pname == "e-" ||
479 pname == "mu-" ||
480 pname == "kaon+" ||
481 pname == "kaon-" ||
482 pname == "lambda" ||
483 pname == "anti_neutron" ||
484 pname == "anti_proton")) yes = true;
485 if(yes) {
486 // main processes
487 std::multimap<PD,HP,std::less<PD> >::iterator it;
488 for(it=p_map.lower_bound(part); it!=p_map.upper_bound(part); ++it) {
489 if(it->first == part) {
490 HP proc = (it->second);
491 G4int j=0;
492 for(; j<n_proc; j++) {
493 if(process[j] == proc) {
494 Print(j, i);
495 }
496 }
497 }
498 }
499 // extra processes
500 std::multimap<PD,G4VProcess*,std::less<PD> >::iterator itp;
501 for(itp=ep_map.lower_bound(part); itp!=ep_map.upper_bound(part); ++itp) {
502 if(itp->first == part) {
503 G4VProcess* proc = (itp->second);
504 if(wasPrinted[i] == 0) {
505 wasPrinted[i] = 1;
506 G4cout<<G4endl;
507 G4cout << " Hadronic Processes for <"
508 <<part->GetParticleName() << ">" << G4endl;
509 }
510 G4cout << " " << proc->GetProcessName() << G4endl;
511 }
512 }
513 }
514 }
515 if(level > 0) {
516 G4cout << "=============================================================="
517 << "=============================="
518 << G4endl;
519 }
520}
521
522//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
523
524void G4HadronicProcessStore::Print(G4int idxProc, G4int idxPart)
525{
526 G4HadronicProcess* proc = process[idxProc];
527 const G4ParticleDefinition* part = particle[idxPart];
528 if(wasPrinted[idxPart] == 0) {
529 wasPrinted[idxPart] = 1;
530 G4cout<<G4endl;
531 G4cout << " Hadronic Processes for <"
532 <<part->GetParticleName() << ">" << G4endl;
533 }
534 HI hi = 0;
535 G4bool first;
536 std::multimap<HP,HI,std::less<HP> >::iterator ih;
537 G4cout << std::setw(20) << proc->GetProcessName()
538 << " Models: ";
539 first = true;
540 for(ih=m_map.lower_bound(proc); ih!=m_map.upper_bound(proc); ++ih) {
541 if(ih->first == proc) {
542 hi = ih->second;
543 G4int i=0;
544 for(; i<n_model; i++) {
545 if(model[i] == hi) break;
546 }
547 if(!first) G4cout << " ";
548 first = false;
549 G4cout << std::setw(25) << modelName[i]
550 << ": Emin(GeV)= "
551 << std::setw(5) << hi->GetMinEnergy()/GeV
552 << " Emax(GeV)= "
553 << hi->GetMaxEnergy()/GeV
554 << G4endl;
555 }
556 }
557}
558
559//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
560
561void G4HadronicProcessStore::SetVerbose(G4int val)
562{
563 verbose = val;
564 G4int i;
565 for(i=0; i<n_proc; i++) {
566 if(process[i]) process[i]->SetVerboseLevel(val);
567 }
568 for(i=0; i<n_model; i++) {
569 if(model[i]) model[i]->SetVerboseLevel(val);
570 }
571}
572
573//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
574
575G4int G4HadronicProcessStore::GetVerbose()
576{
577 return verbose;
578}
579
580//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
581
582G4HadronicProcess* G4HadronicProcessStore::FindProcess(
583 const G4ParticleDefinition* part, G4HadronicProcessType subType)
584{
585 bool isNew = false;
586 G4HadronicProcess* hp = 0;
587
588 if(part != currentParticle) {
589 isNew = true;
590 currentParticle = part;
591 localDP.SetDefinition(const_cast<G4ParticleDefinition*>(part));
592 } else if(!currentProcess) {
593 isNew = true;
594 } else if(subType == currentProcess->GetProcessSubType()) {
595 hp = currentProcess;
596 } else {
597 isNew = true;
598 }
599
600 if(isNew) {
601 std::multimap<PD,HP,std::less<PD> >::iterator it;
602 for(it=p_map.lower_bound(part); it!=p_map.upper_bound(part); ++it) {
603 if(it->first == part && subType == (it->second)->GetProcessSubType()) {
604 hp = it->second;
605 break;
606 }
607 }
608 currentProcess = hp;
609 }
610
611 return hp;
612}
613
614//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
Note: See TracBrowser for help on using the repository browser.