| Geant4 Cross Reference |
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 // ClassName: G4HadronPhysicsQGSP_BERT
30 //
31 // Author: 2002 J.P. Wellisch
32 //
33 // Modified:
34 // 23.11.2005 G.Folger: migration to non static particles
35 // 08.06.2006 V.Ivanchenko: remove stopping
36 // 20.06.2006 G.Folger: Bertini applies to Kaons, i.e. use SetMinEnergy instead of SetMinPionEnergy
37 // 25.04.2007 G.Folger: Add code for quasielastic
38 // 10.12.2007 G.Folger: Add projectilediffrative option for proton/neutron, off by default
39 // 31.10.2012 A.Ribon: Use G4MiscBuilder
40 // 19.03.2013 A.Ribon: Replace LEP with FTFP
41 // 05.05.2020 A.Ribon: Use QGSP for antibaryons at high energies
42 // 07.05.2020 A.Ribon: Use QGSP for hyperons (and anti-hyperons) at high energies
43 //
44 //----------------------------------------------------------------------------
45 //
46 #include <iomanip>
47
48 #include "G4HadronPhysicsQGSP_BERT.hh"
49
50 #include "globals.hh"
51 #include "G4ios.hh"
52 #include "G4SystemOfUnits.hh"
53 #include "G4ParticleDefinition.hh"
54 #include "G4ParticleTable.hh"
55
56 #include "G4PiKBuilder.hh"
57 #include "G4FTFPPiKBuilder.hh"
58 #include "G4QGSPPiKBuilder.hh"
59 #include "G4BertiniPiKBuilder.hh"
60
61 #include "G4ProtonBuilder.hh"
62 #include "G4FTFPProtonBuilder.hh"
63 #include "G4QGSPProtonBuilder.hh"
64 #include "G4BertiniProtonBuilder.hh"
65
66 #include "G4NeutronBuilder.hh"
67 #include "G4FTFPNeutronBuilder.hh"
68 #include "G4QGSPNeutronBuilder.hh"
69 #include "G4BertiniNeutronBuilder.hh"
70
71 #include "G4MesonConstructor.hh"
72 #include "G4BaryonConstructor.hh"
73 #include "G4ShortLivedConstructor.hh"
74 #include "G4IonConstructor.hh"
75
76 #include "G4NeutronRadCapture.hh"
77 #include "G4NeutronInelasticXS.hh"
78 #include "G4NeutronCaptureXS.hh"
79
80 #include "G4PhysListUtil.hh"
81 #include "G4HadParticles.hh"
82 #include "G4HadronicParameters.hh"
83 #include "G4HadronicBuilder.hh"
84 #include "G4BuilderType.hh"
85
86 #include "G4PhysicsConstructorFactory.hh"
87 //
88 G4_DECLARE_PHYSCONSTR_FACTORY(G4HadronPhysicsQGSP_BERT);
89
90 G4HadronPhysicsQGSP_BERT::G4HadronPhysicsQGSP_BERT(G4int verb)
91 : G4HadronPhysicsQGSP_BERT("hInelastic QGSP_BERT",true)
92 {
93 G4HadronicParameters::Instance()->SetVerboseLevel(verb);
94 }
95
96 G4HadronPhysicsQGSP_BERT::G4HadronPhysicsQGSP_BERT(const G4String& name, G4bool)
97 : G4VPhysicsConstructor(name)
98 {
99 SetPhysicsType(bHadronInelastic);
100 QuasiElasticFTF= false; // Use built-in quasi-elastic (not add-on)
101 QuasiElasticQGS= true; // For QGS, it must use it.
102 G4HadronicParameters* param = G4HadronicParameters::Instance();
103 minQGSP_proton = minQGSP_neutron = minQGSP_pik =
104 param->GetMinEnergyTransitionQGS_FTF();
105 maxFTFP_proton = maxFTFP_neutron = maxFTFP_pik =
106 param->GetMaxEnergyTransitionQGS_FTF();
107 minFTFP_proton = minFTFP_neutron = minFTFP_pik =
108 param->GetMinEnergyTransitionFTF_Cascade();
109 maxBERT_proton = maxBERT_neutron = maxBERT_pik =
110 param->GetMaxEnergyTransitionFTF_Cascade();
111 minBERT_proton = minBERT_neutron = minBERT_pik = 0.0;
112 param->SetEnableBCParticles(true);
113 }
114
115 void G4HadronPhysicsQGSP_BERT::CreateModels()
116 {
117 Neutron();
118 Proton();
119 Pion();
120 Others();
121 }
122
123 void G4HadronPhysicsQGSP_BERT::Neutron()
124 {
125 G4HadronicParameters* param = G4HadronicParameters::Instance();
126 G4bool useFactorXS = param->ApplyFactorXS();
127 //General schema:
128 // 1) Create a builder
129 // 2) Call AddBuilder
130 // 3) Configure the builder, possibly with sub-builders
131 // 4) Call builder->Build()
132 auto neu = new G4NeutronBuilder;
133 AddBuilder(neu);
134 auto qgs = new G4QGSPNeutronBuilder(QuasiElasticQGS);
135 AddBuilder(qgs);
136 qgs->SetMinEnergy(minQGSP_neutron);
137 neu->RegisterMe(qgs);
138 auto ftf = new G4FTFPNeutronBuilder(QuasiElasticFTF);
139 AddBuilder(ftf);
140 ftf->SetMinEnergy(minFTFP_neutron);
141 ftf->SetMaxEnergy(maxFTFP_neutron);
142 neu->RegisterMe(ftf);
143 auto bert = new G4BertiniNeutronBuilder;
144 AddBuilder(bert);
145 bert->SetMinEnergy(minBERT_neutron);
146 bert->SetMaxEnergy(maxBERT_neutron);
147 neu->RegisterMe(bert);
148 neu->Build();
149
150 const G4ParticleDefinition* neutron = G4Neutron::Neutron();
151 G4HadronicProcess* inel = G4PhysListUtil::FindInelasticProcess(neutron);
152 if(nullptr != inel) {
153 inel->AddDataSet(new G4NeutronInelasticXS());
154 if( useFactorXS ) inel->MultiplyCrossSectionBy( param->XSFactorNucleonInelastic() );
155 }
156 G4HadronicProcess* capture = G4PhysListUtil::FindCaptureProcess(neutron);
157 if (nullptr != capture) {
158 capture->RegisterMe(new G4NeutronRadCapture());
159 }
160 }
161
162 void G4HadronPhysicsQGSP_BERT::Proton()
163 {
164 G4HadronicParameters* param = G4HadronicParameters::Instance();
165 G4bool useFactorXS = param->ApplyFactorXS();
166
167 auto pro = new G4ProtonBuilder;
168 AddBuilder(pro);
169 auto qgs = new G4QGSPProtonBuilder(QuasiElasticQGS);
170 AddBuilder(qgs);
171 qgs->SetMinEnergy(minQGSP_proton);
172 pro->RegisterMe(qgs);
173 auto ftf = new G4FTFPProtonBuilder(QuasiElasticFTF);
174 AddBuilder(ftf);
175 ftf->SetMinEnergy(minFTFP_proton);
176 ftf->SetMaxEnergy(maxFTFP_proton);
177 pro->RegisterMe(ftf);
178 auto bert = new G4BertiniProtonBuilder;
179 AddBuilder(bert);
180 bert->SetMinEnergy(minBERT_proton);
181 bert->SetMaxEnergy(maxBERT_proton);
182 pro->RegisterMe(bert);
183 pro->Build();
184
185 const G4ParticleDefinition* proton = G4Proton::Proton();
186 G4HadronicProcess* inel = G4PhysListUtil::FindInelasticProcess(proton);
187 if(nullptr != inel) {
188 if( useFactorXS ) inel->MultiplyCrossSectionBy( param->XSFactorNucleonInelastic() );
189 }
190 }
191
192 void G4HadronPhysicsQGSP_BERT::Pion()
193 {
194 G4HadronicParameters* param = G4HadronicParameters::Instance();
195 G4bool useFactorXS = param->ApplyFactorXS();
196
197 auto pik = new G4PiKBuilder;
198 AddBuilder(pik);
199 auto qgs = new G4QGSPPiKBuilder(QuasiElasticQGS);
200 AddBuilder(qgs);
201 qgs->SetMinEnergy(minQGSP_pik);
202 pik->RegisterMe(qgs);
203 auto ftf = new G4FTFPPiKBuilder(QuasiElasticFTF);
204 AddBuilder(ftf);
205 ftf->SetMinEnergy(minFTFP_pik);
206 ftf->SetMaxEnergy(maxFTFP_pik);
207 pik->RegisterMe(ftf);
208 auto bert = new G4BertiniPiKBuilder;
209 AddBuilder(bert);
210 bert->SetMinEnergy(minBERT_pik);
211 bert->SetMaxEnergy(maxBERT_pik);
212 pik->RegisterMe(bert);
213 pik->Build();
214
215 // add cross section factor
216 if( useFactorXS ) {
217 const G4ParticleDefinition* pion = G4PionPlus::PionPlus();
218 G4HadronicProcess* inel = G4PhysListUtil::FindInelasticProcess(pion);
219 if(nullptr != inel) {
220 inel->MultiplyCrossSectionBy( param->XSFactorPionInelastic() );
221 }
222 pion = G4PionMinus::PionMinus();
223 inel = G4PhysListUtil::FindInelasticProcess(pion);
224 if(inel) {
225 inel->MultiplyCrossSectionBy( param->XSFactorPionInelastic() );
226 }
227 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
228 for( auto & pdg : G4HadParticles::GetKaons() ) {
229 auto part = table->FindParticle( pdg );
230 if ( part == nullptr ) { continue; }
231 inel = G4PhysListUtil::FindInelasticProcess(part);
232 if(nullptr != inel) {
233 inel->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
234 }
235 }
236 }
237 }
238
239 void G4HadronPhysicsQGSP_BERT::Others()
240 {
241 G4HadronicParameters* param = G4HadronicParameters::Instance();
242
243 // high energy particles
244 if( param->GetMaxEnergy() > param->EnergyThresholdForHeavyHadrons() ) {
245
246 // anti light ions
247 G4HadronicBuilder::BuildAntiLightIonsFTFP();
248
249 // hyperons
250 G4HadronicBuilder::BuildHyperonsQGSP_FTFP_BERT(true);
251
252 // b-, c- baryons and mesons
253 if( param->EnableBCParticles() ) {
254 G4HadronicBuilder::BuildBCHadronsQGSP_FTFP_BERT(true);
255 }
256 }
257 }
258
259 G4HadronPhysicsQGSP_BERT::~G4HadronPhysicsQGSP_BERT()
260 {}
261
262 void G4HadronPhysicsQGSP_BERT::ConstructParticle()
263 {
264 G4MesonConstructor pMesonConstructor;
265 pMesonConstructor.ConstructParticle();
266
267 G4BaryonConstructor pBaryonConstructor;
268 pBaryonConstructor.ConstructParticle();
269
270 G4ShortLivedConstructor pShortLivedConstructor;
271 pShortLivedConstructor.ConstructParticle();
272
273 G4IonConstructor pIonConstructor;
274 pIonConstructor.ConstructParticle();
275 }
276
277 void G4HadronPhysicsQGSP_BERT::ConstructProcess()
278 {
279 // allow changing of parameters at PreInit
280 G4HadronicParameters* param = G4HadronicParameters::Instance();
281 minQGSP_proton = minQGSP_neutron = minQGSP_pik =
282 param->GetMinEnergyTransitionQGS_FTF();
283 maxFTFP_proton = maxFTFP_neutron = maxFTFP_pik =
284 param->GetMaxEnergyTransitionQGS_FTF();
285 minFTFP_proton = minFTFP_neutron = minFTFP_pik =
286 param->GetMinEnergyTransitionFTF_Cascade();
287 maxBERT_pik = param->GetMaxEnergyTransitionFTF_Cascade();
288
289 if(G4Threading::IsMasterThread() && param->GetVerboseLevel() > 0) {
290 DumpBanner();
291 }
292 CreateModels();
293 }
294
295 void G4HadronPhysicsQGSP_BERT::DumpBanner()
296 {
297 G4cout << G4endl;
298 G4cout << " " << GetPhysicsName() << " Thresholds: " << G4endl;
299 G4cout << " 1) between BERT and FTF/P over the interval "
300 << minFTFP_proton/GeV << " to " << maxBERT_proton/GeV << " GeV. " << G4endl;
301 G4cout << " 2) between FTF/P and QGS/P over the interval "
302 << minQGSP_proton/GeV << " to " << maxFTFP_proton/GeV << " GeV. " << G4endl;
303 G4cout << " -- quasiElastic: " << QuasiElasticQGS << " for QGS "
304 << " and " << QuasiElasticFTF << " for FTF" << G4endl;
305 }
306
307
308