| 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 // Geant4 class G4HadronicBuilder 27 // Geant4 class G4HadronicBuilder
28 // 28 //
29 // Author V.Ivanchenko 14.05.2020 29 // Author V.Ivanchenko 14.05.2020
30 // 30 //
31 31
32 #include "G4HadronicBuilder.hh" 32 #include "G4HadronicBuilder.hh"
33 #include "G4HadParticles.hh" 33 #include "G4HadParticles.hh"
34 #include "G4HadProcesses.hh" 34 #include "G4HadProcesses.hh"
35 35
36 #include "G4ParticleDefinition.hh" 36 #include "G4ParticleDefinition.hh"
37 #include "G4ParticleTable.hh" 37 #include "G4ParticleTable.hh"
38 #include "G4PhysicsListHelper.hh" 38 #include "G4PhysicsListHelper.hh"
39 #include "G4SystemOfUnits.hh" 39 #include "G4SystemOfUnits.hh"
40 40
41 #include "G4HadronicParameters.hh" 41 #include "G4HadronicParameters.hh"
42 42
43 #include "G4TheoFSGenerator.hh" 43 #include "G4TheoFSGenerator.hh"
44 #include "G4FTFModel.hh" 44 #include "G4FTFModel.hh"
45 #include "G4ExcitedStringDecay.hh" 45 #include "G4ExcitedStringDecay.hh"
46 #include "G4GeneratorPrecompoundInterface.hh" 46 #include "G4GeneratorPrecompoundInterface.hh"
47 47
48 #include "G4QGSModel.hh" 48 #include "G4QGSModel.hh"
49 #include "G4QGSParticipants.hh" 49 #include "G4QGSParticipants.hh"
50 #include "G4QGSMFragmentation.hh" 50 #include "G4QGSMFragmentation.hh"
51 #include "G4QuasiElasticChannel.hh" 51 #include "G4QuasiElasticChannel.hh"
52 52
53 #include "G4CascadeInterface.hh" 53 #include "G4CascadeInterface.hh"
54 #include "G4CrossSectionDataSetRegistry.hh" 54 #include "G4CrossSectionDataSetRegistry.hh"
55 #include "G4CrossSectionInelastic.hh" 55 #include "G4CrossSectionInelastic.hh"
56 #include "G4CrossSectionElastic.hh" 56 #include "G4CrossSectionElastic.hh"
57 #include "G4HadronElastic.hh" 57 #include "G4HadronElastic.hh"
58 #include "G4CrossSectionDataSetRegistry.hh" 58 #include "G4CrossSectionDataSetRegistry.hh"
59 59
60 #include "G4HadronElasticProcess.hh" 60 #include "G4HadronElasticProcess.hh"
61 #include "G4HadronInelasticProcess.hh" 61 #include "G4HadronInelasticProcess.hh"
62 62
63 #include "G4DecayTable.hh" 63 #include "G4DecayTable.hh"
64 #include "G4VDecayChannel.hh" 64 #include "G4VDecayChannel.hh"
65 #include "G4PhaseSpaceDecayChannel.hh" 65 #include "G4PhaseSpaceDecayChannel.hh"
66 66
67 #include "G4PreCompoundModel.hh" 67 #include "G4PreCompoundModel.hh"
68 #include "G4INCLXXInterface.hh" 68 #include "G4INCLXXInterface.hh"
69 #include "G4ComponentAntiNuclNuclearXS.hh" 69 #include "G4ComponentAntiNuclNuclearXS.hh"
70 70
71 71
72 72
73 void G4HadronicBuilder::BuildFTFP_BERT(const s 73 void G4HadronicBuilder::BuildFTFP_BERT(const std::vector<G4int>& partList,
74 G4bool 74 G4bool bert, const G4String& xsName) {
75 75
76 G4HadronicParameters* param = G4HadronicPara 76 G4HadronicParameters* param = G4HadronicParameters::Instance();
77 G4PhysicsListHelper* ph = G4PhysicsListHelpe 77 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
78 78
79 auto theModel = new G4TheoFSGenerator("FTFP" 79 auto theModel = new G4TheoFSGenerator("FTFP");
80 auto theStringModel = new G4FTFModel(); 80 auto theStringModel = new G4FTFModel();
81 theStringModel->SetFragmentationModel(new G4 81 theStringModel->SetFragmentationModel(new G4ExcitedStringDecay());
82 theModel->SetHighEnergyGenerator( theStringM 82 theModel->SetHighEnergyGenerator( theStringModel );
83 theModel->SetTransport( new G4GeneratorPreco 83 theModel->SetTransport( new G4GeneratorPrecompoundInterface() );
84 theModel->SetMaxEnergy( param->GetMaxEnergy( 84 theModel->SetMaxEnergy( param->GetMaxEnergy() );
85 85
86 G4CascadeInterface* theCascade = nullptr; 86 G4CascadeInterface* theCascade = nullptr;
87 if(bert) { 87 if(bert) {
88 theCascade = new G4CascadeInterface(); 88 theCascade = new G4CascadeInterface();
89 theCascade->SetMaxEnergy( param->GetMaxEne 89 theCascade->SetMaxEnergy( param->GetMaxEnergyTransitionFTF_Cascade() );
90 theModel->SetMinEnergy( param->GetMinEnerg 90 theModel->SetMinEnergy( param->GetMinEnergyTransitionFTF_Cascade() );
91 } 91 }
92 92
93 auto xsinel = G4HadProcesses::InelasticXS( x 93 auto xsinel = G4HadProcesses::InelasticXS( xsName );
94 94
95 G4ParticleTable* table = G4ParticleTable::Ge 95 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
96 for( auto & pdg : partList ) { 96 for( auto & pdg : partList ) {
97 97
98 auto part = table->FindParticle( pdg ); 98 auto part = table->FindParticle( pdg );
99 if ( part == nullptr ) { continue; } 99 if ( part == nullptr ) { continue; }
100 100
101 auto hadi = new G4HadronInelasticProcess( 101 auto hadi = new G4HadronInelasticProcess( part->GetParticleName()+"Inelastic", part );
102 hadi->AddDataSet( xsinel ); 102 hadi->AddDataSet( xsinel );
103 hadi->RegisterMe( theModel ); 103 hadi->RegisterMe( theModel );
104 if( theCascade != nullptr ) hadi->Register 104 if( theCascade != nullptr ) hadi->RegisterMe( theCascade );
105 if( param->ApplyFactorXS() ) hadi->Multipl 105 if( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
106 ph->RegisterProcess(hadi, part); 106 ph->RegisterProcess(hadi, part);
107 } 107 }
108 } 108 }
109 109
110 void G4HadronicBuilder::BuildFTFQGSP_BERT(cons 110 void G4HadronicBuilder::BuildFTFQGSP_BERT(const std::vector<G4int>& partList,
111 G4bo 111 G4bool bert, const G4String& xsName) {
112 112
113 G4HadronicParameters* param = G4HadronicPara 113 G4HadronicParameters* param = G4HadronicParameters::Instance();
114 G4PhysicsListHelper* ph = G4PhysicsListHelpe 114 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
115 115
116 auto theModel = new G4TheoFSGenerator("FTFQG 116 auto theModel = new G4TheoFSGenerator("FTFQGSP");
117 auto theStringModel = new G4FTFModel(); 117 auto theStringModel = new G4FTFModel();
118 theStringModel->SetFragmentationModel(new G4 118 theStringModel->SetFragmentationModel(new G4ExcitedStringDecay( new G4QGSMFragmentation() ) );
119 theModel->SetHighEnergyGenerator( theStringM 119 theModel->SetHighEnergyGenerator( theStringModel );
120 theModel->SetTransport( new G4GeneratorPreco 120 theModel->SetTransport( new G4GeneratorPrecompoundInterface() );
121 theModel->SetMaxEnergy( param->GetMaxEnergy( 121 theModel->SetMaxEnergy( param->GetMaxEnergy() );
122 122
123 G4CascadeInterface* theCascade = nullptr; 123 G4CascadeInterface* theCascade = nullptr;
124 if(bert) { 124 if(bert) {
125 theCascade = new G4CascadeInterface(); 125 theCascade = new G4CascadeInterface();
126 theCascade->SetMaxEnergy( param->GetMaxEne 126 theCascade->SetMaxEnergy( param->GetMaxEnergyTransitionFTF_Cascade() );
127 theModel->SetMinEnergy( param->GetMinEnerg 127 theModel->SetMinEnergy( param->GetMinEnergyTransitionFTF_Cascade() );
128 } 128 }
129 129
130 auto xsinel = G4HadProcesses::InelasticXS( x 130 auto xsinel = G4HadProcesses::InelasticXS( xsName );
131 131
132 G4ParticleTable* table = G4ParticleTable::Ge 132 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
133 for( auto & pdg : partList ) { 133 for( auto & pdg : partList ) {
134 134
135 auto part = table->FindParticle( pdg ); 135 auto part = table->FindParticle( pdg );
136 if ( part == nullptr ) { continue; } 136 if ( part == nullptr ) { continue; }
137 137
138 auto hadi = new G4HadronInelasticProcess( 138 auto hadi = new G4HadronInelasticProcess( part->GetParticleName()+"Inelastic", part );
139 hadi->AddDataSet( xsinel ); 139 hadi->AddDataSet( xsinel );
140 hadi->RegisterMe( theModel ); 140 hadi->RegisterMe( theModel );
141 if( theCascade != nullptr ) hadi->Register 141 if( theCascade != nullptr ) hadi->RegisterMe( theCascade );
142 if( param->ApplyFactorXS() ) hadi->Multipl 142 if( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
143 ph->RegisterProcess(hadi, part); 143 ph->RegisterProcess(hadi, part);
144 } 144 }
145 } 145 }
146 146
147 void G4HadronicBuilder::BuildQGSP_FTFP_BERT(co 147 void G4HadronicBuilder::BuildQGSP_FTFP_BERT(const std::vector<G4int>& partList,
148 G4 148 G4bool bert, G4bool quasiElastic,
149 co 149 const G4String& xsName) {
150 150
151 G4HadronicParameters* param = G4HadronicPara 151 G4HadronicParameters* param = G4HadronicParameters::Instance();
152 G4PhysicsListHelper* ph = G4PhysicsListHelpe 152 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
153 153
154 auto theTransport = new G4GeneratorPrecompou 154 auto theTransport = new G4GeneratorPrecompoundInterface();
155 155
156 auto theHEModel = new G4TheoFSGenerator("QGS 156 auto theHEModel = new G4TheoFSGenerator("QGSP");
157 G4QGSModel< G4QGSParticipants >* theQGSModel 157 G4QGSModel< G4QGSParticipants >* theQGSModel = new G4QGSModel< G4QGSParticipants >;
158 theQGSModel->SetFragmentationModel( new G4Ex 158 theQGSModel->SetFragmentationModel( new G4ExcitedStringDecay( new G4QGSMFragmentation() ) );
159 theHEModel->SetTransport( theTransport ); 159 theHEModel->SetTransport( theTransport );
160 theHEModel->SetHighEnergyGenerator( theQGSMo 160 theHEModel->SetHighEnergyGenerator( theQGSModel );
161 if (quasiElastic) { 161 if (quasiElastic) {
162 theHEModel->SetQuasiElasticChannel(new G4Q 162 theHEModel->SetQuasiElasticChannel(new G4QuasiElasticChannel());
163 } 163 }
164 theHEModel->SetMinEnergy( param->GetMinEnerg 164 theHEModel->SetMinEnergy( param->GetMinEnergyTransitionQGS_FTF() );
165 theHEModel->SetMaxEnergy( param->GetMaxEnerg 165 theHEModel->SetMaxEnergy( param->GetMaxEnergy() );
166 166
167 auto theLEModel = new G4TheoFSGenerator("FTF 167 auto theLEModel = new G4TheoFSGenerator("FTFP");
168 auto theFTFModel = new G4FTFModel(); 168 auto theFTFModel = new G4FTFModel();
169 theFTFModel->SetFragmentationModel(new G4Exc 169 theFTFModel->SetFragmentationModel(new G4ExcitedStringDecay());
170 theLEModel->SetHighEnergyGenerator( theFTFMo 170 theLEModel->SetHighEnergyGenerator( theFTFModel );
171 theLEModel->SetTransport( theTransport ); 171 theLEModel->SetTransport( theTransport );
172 theLEModel->SetMaxEnergy( param->GetMaxEnerg 172 theLEModel->SetMaxEnergy( param->GetMaxEnergyTransitionQGS_FTF() );
173 173
174 G4CascadeInterface* theCascade = nullptr; 174 G4CascadeInterface* theCascade = nullptr;
175 if(bert) { 175 if(bert) {
176 theCascade = new G4CascadeInterface(); 176 theCascade = new G4CascadeInterface();
177 theCascade->SetMaxEnergy( param->GetMaxEne 177 theCascade->SetMaxEnergy( param->GetMaxEnergyTransitionFTF_Cascade() );
178 theLEModel->SetMinEnergy( param->GetMinEne 178 theLEModel->SetMinEnergy( param->GetMinEnergyTransitionFTF_Cascade() );
179 } 179 }
180 180
181 auto xsinel = G4HadProcesses::InelasticXS( x 181 auto xsinel = G4HadProcesses::InelasticXS( xsName );
182 182
183 G4ParticleTable* table = G4ParticleTable::Ge 183 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
184 for( auto & pdg : partList ) { 184 for( auto & pdg : partList ) {
185 185
186 auto part = table->FindParticle( pdg ); 186 auto part = table->FindParticle( pdg );
187 if ( part == nullptr ) { continue; } 187 if ( part == nullptr ) { continue; }
188 188
189 auto hadi = new G4HadronInelasticProcess( 189 auto hadi = new G4HadronInelasticProcess( part->GetParticleName()+"Inelastic", part );
190 hadi->AddDataSet( xsinel ); 190 hadi->AddDataSet( xsinel );
191 hadi->RegisterMe( theHEModel ); 191 hadi->RegisterMe( theHEModel );
192 hadi->RegisterMe( theLEModel ); 192 hadi->RegisterMe( theLEModel );
193 if(theCascade != nullptr) hadi->RegisterMe 193 if(theCascade != nullptr) hadi->RegisterMe( theCascade );
194 if( param->ApplyFactorXS() ) hadi->Multipl 194 if( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
195 ph->RegisterProcess(hadi, part); 195 ph->RegisterProcess(hadi, part);
196 } 196 }
197 } 197 }
198 198
199 void G4HadronicBuilder::BuildINCLXX(const std: 199 void G4HadronicBuilder::BuildINCLXX(const std::vector<G4int>& partList,
200 G4bool 200 G4bool bert, const G4String& xsName) {
201 201
202 // FTF 202 // FTF
203 G4HadronicParameters* param = G4HadronicPara 203 G4HadronicParameters* param = G4HadronicParameters::Instance();
204 G4PhysicsListHelper* ph = G4PhysicsListHelpe 204 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
205 205
206 auto theModel = new G4TheoFSGenerator("FTFP" 206 auto theModel = new G4TheoFSGenerator("FTFP");
207 auto theStringModel = new G4FTFModel(); 207 auto theStringModel = new G4FTFModel();
208 theStringModel->SetFragmentationModel(new G4 208 theStringModel->SetFragmentationModel(new G4ExcitedStringDecay());
209 theModel->SetHighEnergyGenerator( theStringM 209 theModel->SetHighEnergyGenerator( theStringModel );
210 theModel->SetTransport( new G4GeneratorPreco 210 theModel->SetTransport( new G4GeneratorPrecompoundInterface() );
211 theModel->SetMaxEnergy( param->GetMaxEnergy( 211 theModel->SetMaxEnergy( param->GetMaxEnergy() );
212 212
213 G4CascadeInterface* theCascade = nullptr; 213 G4CascadeInterface* theCascade = nullptr;
214 if(bert) { 214 if(bert) {
215 theCascade = new G4CascadeInterface(); 215 theCascade = new G4CascadeInterface();
216 theCascade->SetMaxEnergy( param->GetMaxEne 216 theCascade->SetMaxEnergy( param->GetMaxEnergyTransitionFTF_Cascade() );
217 theModel->SetMinEnergy( param->GetMinEnerg 217 theModel->SetMinEnergy( param->GetMinEnergyTransitionFTF_Cascade() );
218 } 218 }
219 219
220 // INCLXX 220 // INCLXX
221 auto theModelINCLXX = new G4INCLXXInterface( 221 auto theModelINCLXX = new G4INCLXXInterface();
222 theModelINCLXX->SetMinEnergy( param->GetMinE 222 theModelINCLXX->SetMinEnergy( param->GetMinEnergyINCLXX_Pbar() );
223 theModelINCLXX->SetMaxEnergy( param->GetMaxE 223 theModelINCLXX->SetMaxEnergy( param->GetMaxEnergyINCLXX_Pbar() );
224 224
225 // 225 //
226 auto xsinel = G4HadProcesses::InelasticXS( x 226 auto xsinel = G4HadProcesses::InelasticXS( xsName );
227 227
228 G4ParticleTable* table = G4ParticleTable::Ge 228 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
229 for( auto & pdg : partList ) { 229 for( auto & pdg : partList ) {
230 230
231 auto part = table->FindParticle( pdg ); 231 auto part = table->FindParticle( pdg );
232 if ( part == nullptr ) { continue; } 232 if ( part == nullptr ) { continue; }
233 233
234 auto hadi = new G4HadronInelasticProcess( 234 auto hadi = new G4HadronInelasticProcess( part->GetParticleName()+"Inelastic", part );
235 if( pdg == -2212 ) { // pbar use INCLXX 235 if( pdg == -2212 ) { // pbar use INCLXX
236 hadi->AddDataSet( xsinel ); 236 hadi->AddDataSet( xsinel );
237 hadi->RegisterMe( theModelINCLXX ); 237 hadi->RegisterMe( theModelINCLXX );
238 if( param->ApplyFactorXS() ) hadi->Multi 238 if( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
239 ph->RegisterProcess(hadi, part); 239 ph->RegisterProcess(hadi, part);
240 } else { // other anti-X use FTF 240 } else { // other anti-X use FTF
241 hadi->AddDataSet( xsinel ); 241 hadi->AddDataSet( xsinel );
242 hadi->RegisterMe( theModel ); 242 hadi->RegisterMe( theModel );
243 if( theCascade != nullptr ) hadi->Register 243 if( theCascade != nullptr ) hadi->RegisterMe( theCascade );
244 if( param->ApplyFactorXS() ) hadi->Multipl 244 if( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
245 ph->RegisterProcess(hadi, part); 245 ph->RegisterProcess(hadi, part);
246 } 246 }
247 } 247 }
248 } 248 }
249 249
250 void G4HadronicBuilder::BuildElastic(const std 250 void G4HadronicBuilder::BuildElastic(const std::vector<G4int>& partList) {
251 251
252 G4HadronicParameters* param = G4HadronicPara 252 G4HadronicParameters* param = G4HadronicParameters::Instance();
253 G4PhysicsListHelper* ph = G4PhysicsListHelpe 253 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
254 254
255 auto xsel = G4HadProcesses::ElasticXS("Glaub 255 auto xsel = G4HadProcesses::ElasticXS("Glauber-Gribov");
256 256
257 auto elModel = new G4HadronElastic(); 257 auto elModel = new G4HadronElastic();
258 elModel->SetMaxEnergy( param->GetMaxEnergy() 258 elModel->SetMaxEnergy( param->GetMaxEnergy() );
259 259
260 G4ParticleTable* table = G4ParticleTable::Ge 260 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
261 for( auto & pdg : partList ) { 261 for( auto & pdg : partList ) {
262 262
263 auto part = table->FindParticle( pdg ); 263 auto part = table->FindParticle( pdg );
264 if ( part == nullptr ) { continue; } 264 if ( part == nullptr ) { continue; }
265 265
266 auto hade = new G4HadronElasticProcess(); 266 auto hade = new G4HadronElasticProcess();
267 hade->AddDataSet( xsel ); 267 hade->AddDataSet( xsel );
268 hade->RegisterMe( elModel ); 268 hade->RegisterMe( elModel );
269 if( param->ApplyFactorXS() ) hade->Multipl 269 if( param->ApplyFactorXS() ) hade->MultiplyCrossSectionBy( param->XSFactorHadronElastic() );
270 ph->RegisterProcess(hade, part); 270 ph->RegisterProcess(hade, part);
271 } 271 }
272 } 272 }
273 273
274 void G4HadronicBuilder::BuildHyperonsFTFP_BERT 274 void G4HadronicBuilder::BuildHyperonsFTFP_BERT() {
275 // For hyperons, Bertini is used at low ener 275 // For hyperons, Bertini is used at low energies;
276 // for anti-hyperons, FTFP can be used down 276 // for anti-hyperons, FTFP can be used down to zero kinetic energy.
277 BuildFTFP_BERT(G4HadParticles::GetHyperons() 277 BuildFTFP_BERT(G4HadParticles::GetHyperons(), true, "Glauber-Gribov");
278 BuildFTFP_BERT(G4HadParticles::GetAntiHypero 278 BuildFTFP_BERT(G4HadParticles::GetAntiHyperons(), false, "Glauber-Gribov");
279 } 279 }
280 280
281 void G4HadronicBuilder::BuildHyperonsFTFQGSP_B 281 void G4HadronicBuilder::BuildHyperonsFTFQGSP_BERT() {
282 // For hyperons, Bertini is used at low ener 282 // For hyperons, Bertini is used at low energies;
283 // for anti-hyperons, FTFP can be used down 283 // for anti-hyperons, FTFP can be used down to zero kinetic energy.
284 BuildFTFQGSP_BERT(G4HadParticles::GetHyperon 284 BuildFTFQGSP_BERT(G4HadParticles::GetHyperons(), true, "Glauber-Gribov");
285 BuildFTFQGSP_BERT(G4HadParticles::GetAntiHyp 285 BuildFTFQGSP_BERT(G4HadParticles::GetAntiHyperons(), false, "Glauber-Gribov");
286 } 286 }
287 287
288 void G4HadronicBuilder::BuildHyperonsQGSP_FTFP 288 void G4HadronicBuilder::BuildHyperonsQGSP_FTFP_BERT(G4bool qElastic) {
289 // For hyperons, Bertini is used at low ener 289 // For hyperons, Bertini is used at low energies;
290 // for anti-hyperons, FTFP can be used down 290 // for anti-hyperons, FTFP can be used down to zero kinetic energy.
291 // QGSP is used at high energies in all case 291 // QGSP is used at high energies in all cases.
292 BuildQGSP_FTFP_BERT(G4HadParticles::GetHyper 292 BuildQGSP_FTFP_BERT(G4HadParticles::GetHyperons(), true, qElastic, "Glauber-Gribov");
293 BuildQGSP_FTFP_BERT(G4HadParticles::GetAntiH 293 BuildQGSP_FTFP_BERT(G4HadParticles::GetAntiHyperons(), false, qElastic, "Glauber-Gribov");
294 } 294 }
295 295
296 void G4HadronicBuilder::BuildKaonsFTFP_BERT() 296 void G4HadronicBuilder::BuildKaonsFTFP_BERT() {
297 BuildFTFP_BERT(G4HadParticles::GetKaons(), t 297 BuildFTFP_BERT(G4HadParticles::GetKaons(), true, "Glauber-Gribov");
298 } 298 }
299 299
300 void G4HadronicBuilder::BuildKaonsFTFQGSP_BERT 300 void G4HadronicBuilder::BuildKaonsFTFQGSP_BERT() {
301 BuildFTFQGSP_BERT(G4HadParticles::GetKaons() 301 BuildFTFQGSP_BERT(G4HadParticles::GetKaons(), true, "Glauber-Gribov");
302 } 302 }
303 303
304 void G4HadronicBuilder::BuildKaonsQGSP_FTFP_BE 304 void G4HadronicBuilder::BuildKaonsQGSP_FTFP_BERT(G4bool qElastic) {
305 BuildQGSP_FTFP_BERT(G4HadParticles::GetKaons 305 BuildQGSP_FTFP_BERT(G4HadParticles::GetKaons(), true, qElastic, "Glauber-Gribov");
306 } 306 }
307 307
308 void G4HadronicBuilder::BuildAntiLightIonsFTFP 308 void G4HadronicBuilder::BuildAntiLightIonsFTFP() {
309 BuildFTFP_BERT(G4HadParticles::GetLightAntiI 309 BuildFTFP_BERT(G4HadParticles::GetLightAntiIons(), false, "AntiAGlauber");
310 } 310 }
311 311
312 //void G4HadronicBuilder::BuildAntiLightIonsQG 312 //void G4HadronicBuilder::BuildAntiLightIonsQGSP_FTFP(G4bool qElastic) {
313 // Note: currently QGSP cannot be applied for 313 // Note: currently QGSP cannot be applied for any ion or anti-ion!
314 // BuildQGSP_FTFP_BERT(G4HadParticles::GetLig 314 // BuildQGSP_FTFP_BERT(G4HadParticles::GetLightAntiIons(), false, qElastic, "AntiAGlauber");
315 //} 315 //}
316 316
317 void G4HadronicBuilder::BuildAntiLightIonsINCL 317 void G4HadronicBuilder::BuildAntiLightIonsINCLXX() {
318 BuildINCLXX(G4HadParticles::GetLightAntiIons 318 BuildINCLXX(G4HadParticles::GetLightAntiIons(), false, "AntiAGlauber");
319 } 319 }
320 320
321 void G4HadronicBuilder::BuildBCHadronsFTFP_BER 321 void G4HadronicBuilder::BuildBCHadronsFTFP_BERT() {
322 if( G4HadronicParameters::Instance()->Enable 322 if( G4HadronicParameters::Instance()->EnableBCParticles() ) {
323 // Bertini is not applicable for charm and 323 // Bertini is not applicable for charm and bottom hadrons, therefore FTFP is used
324 // down to zero kinetic energy (but at ver 324 // down to zero kinetic energy (but at very low energies, a dummy model is used
325 // that returns the projectile heavy hadro 325 // that returns the projectile heavy hadron in the final state).
326 BuildFTFP_BERT(G4HadParticles::GetBCHadron 326 BuildFTFP_BERT(G4HadParticles::GetBCHadrons(), false, "Glauber-Gribov");
327 BuildDecayTableForBCHadrons(); 327 BuildDecayTableForBCHadrons();
328 } 328 }
329 } 329 }
330 330
331 void G4HadronicBuilder::BuildBCHadronsFTFQGSP_ 331 void G4HadronicBuilder::BuildBCHadronsFTFQGSP_BERT() {
332 if( G4HadronicParameters::Instance()->Enable 332 if( G4HadronicParameters::Instance()->EnableBCParticles() ) {
333 // Bertini is not applicable for charm and 333 // Bertini is not applicable for charm and bottom hadrons, therefore FTFP is used
334 // down to zero kinetic energy (but at ver 334 // down to zero kinetic energy (but at very low energies, a dummy model is used
335 // that returns the projectile heavy hadro 335 // that returns the projectile heavy hadron in the final state).
336 BuildFTFQGSP_BERT(G4HadParticles::GetBCHad 336 BuildFTFQGSP_BERT(G4HadParticles::GetBCHadrons(), false, "Glauber-Gribov");
337 BuildDecayTableForBCHadrons(); 337 BuildDecayTableForBCHadrons();
338 } 338 }
339 } 339 }
340 340
341 void G4HadronicBuilder::BuildBCHadronsQGSP_FTF 341 void G4HadronicBuilder::BuildBCHadronsQGSP_FTFP_BERT(G4bool qElastic) {
342 if( G4HadronicParameters::Instance()->Enable 342 if( G4HadronicParameters::Instance()->EnableBCParticles() ) {
343 // Bertini is not applicable for charm and 343 // Bertini is not applicable for charm and bottom hadrons, therefore FTFP is used
344 // down to zero kinetic energy (but at ver 344 // down to zero kinetic energy (but at very low energies, a dummy model is used
345 // that returns the projectile heavy hadro 345 // that returns the projectile heavy hadron in the final state).
346 // QGSP is used at high energies in all ca 346 // QGSP is used at high energies in all cases.
347 BuildQGSP_FTFP_BERT(G4HadParticles::GetBCH 347 BuildQGSP_FTFP_BERT(G4HadParticles::GetBCHadrons(), false, qElastic, "Glauber-Gribov");
348 BuildDecayTableForBCHadrons(); 348 BuildDecayTableForBCHadrons();
349 } 349 }
350 } 350 }
351 351
352 void G4HadronicBuilder::BuildDecayTableForBCHa 352 void G4HadronicBuilder::BuildDecayTableForBCHadrons() {
353 // Geant4 does not define the decay of most 353 // Geant4 does not define the decay of most of charmed and bottom hadrons.
354 // The reason is that most of these heavy ha 354 // The reason is that most of these heavy hadrons have many different
355 // decay channels, with a complex dynamics, 355 // decay channels, with a complex dynamics, quite different from the flat
356 // phase space kinematical treatment used in 356 // phase space kinematical treatment used in Geant4 for most of hadronic decays.
357 // High-energy experiments usually use dedic 357 // High-energy experiments usually use dedicated Monte Carlo Event Generators
358 // for the decays of charmed and bottom hadr 358 // for the decays of charmed and bottom hadrons; therefore, these heavy
359 // hadrons, which are passed to Geant4 as pr 359 // hadrons, which are passed to Geant4 as primary tracks, have pre-assigned
360 // decays. Moreover, no charmed or bottom se 360 // decays. Moreover, no charmed or bottom secondary hadrons were created
361 // in Geant4 hadronic interactions before Ge 361 // in Geant4 hadronic interactions before Geant4 10.7.
362 // With the extension of Geant4 hadronic int 362 // With the extension of Geant4 hadronic interactions to charmed and bottom
363 // hadrons, in version Geant4 10.7, we do ne 363 // hadrons, in version Geant4 10.7, we do need to define decays in Geant4
364 // for these heavy hadrons, for two reasons: 364 // for these heavy hadrons, for two reasons:
365 // 1. For testing purposes, unless we pre-as 365 // 1. For testing purposes, unless we pre-assign decays of heavy hadrons
366 // (as the HEP experiments normally do by 366 // (as the HEP experiments normally do by using MC Event Generators);
367 // 2. To avoid crashes (due to missing decay 367 // 2. To avoid crashes (due to missing decay channels) whenever charmed or
368 // bottom secondary hadrons are produced 368 // bottom secondary hadrons are produced by Geant4 hadronic interactions,
369 // even with ordinary (i.e. not heavy) ha 369 // even with ordinary (i.e. not heavy) hadron projectiles, because in
370 // this case we cannot (easily!) pre-assi 370 // this case we cannot (easily!) pre-assign decays to them.
371 // Given that 1. is just a convenience for t 371 // Given that 1. is just a convenience for testing, and 2. happens rather
372 // rarely in practice - because very few pri 372 // rarely in practice - because very few primary energetic (i.e. boosted)
373 // heavy hadrons fly enough to reach the bea 373 // heavy hadrons fly enough to reach the beam pipe or the tracker and
374 // having an inelastic interaction there, an 374 // having an inelastic interaction there, and the very low probability
375 // to create a heavy hadrons from the string 375 // to create a heavy hadrons from the string fragmentation in ordinary
376 // (i.e. not heavy) hadronic interactions - 376 // (i.e. not heavy) hadronic interactions - there is no need in practice
377 // to define accurately the decays of heavy 377 // to define accurately the decays of heavy hadrons in Geant4.
378 // So, for our practical purposes, it is eno 378 // So, for our practical purposes, it is enough to define very simple,
379 // "dummy" decays of charmed and bottom hadr 379 // "dummy" decays of charmed and bottom hadrons.
380 // Here we use a single, fully hadronic chan 380 // Here we use a single, fully hadronic channel, with 2 or 3 or 4
381 // daughters, for each of these heavy hadron 381 // daughters, for each of these heavy hadrons, assigning to this single
382 // decay channel a 100% branching ratio, alt 382 // decay channel a 100% branching ratio, although in reality such a
383 // channel is one between hundreds of possib 383 // channel is one between hundreds of possible ones (and therefore its
384 // real branching ratio is typical of a few 384 // real branching ratio is typical of a few per-cent); moreover, we treat
385 // the decay without any dynamics, i.e. with 385 // the decay without any dynamics, i.e. with a flat phase space kinematical
386 // treatment. 386 // treatment.
387 // Note that some of the charmed and bottom 387 // Note that some of the charmed and bottom hadrons such as SigmaC++,
388 // SigmaC+, SigmaC0, SigmaB+, SigmaB0 and Si 388 // SigmaC+, SigmaC0, SigmaB+, SigmaB0 and SigmaB- have one dominant
389 // decay channel (to LambdaC/B + Pion) which 389 // decay channel (to LambdaC/B + Pion) which is already defined in Geant4.
390 // This is not the case for EtaC, JPsi and U 390 // This is not the case for EtaC, JPsi and Upsilon, whose decays need to
391 // be defined here (although they decay so q 391 // be defined here (although they decay so quickly that their hadronic
392 // interactions can be neglected, as we do f 392 // interactions can be neglected, as we do for Pi0 and Sigma0).
393 // Note that our definition of the decay tab 393 // Note that our definition of the decay tables for these heavy hadrons
394 // do not interfere with the pre-assign deca 394 // do not interfere with the pre-assign decays of primary charmed and
395 // bottom tracks made by the HEP experiments 395 // bottom tracks made by the HEP experiments. In fact, pre-assign decays
396 // have priority over (i.e. override) decay 396 // have priority over (i.e. override) decay tables.
397 static G4bool isFirstCall = true; 397 static G4bool isFirstCall = true;
398 if ( ! isFirstCall ) return; 398 if ( ! isFirstCall ) return;
399 isFirstCall = false; 399 isFirstCall = false;
400 G4ParticleTable* particleTable = G4ParticleT 400 G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
401 for ( auto & pdg : G4HadParticles::GetBCHadr 401 for ( auto & pdg : G4HadParticles::GetBCHadrons() ) {
402 auto part = particleTable->FindParticle( p 402 auto part = particleTable->FindParticle( pdg );
403 if ( part == nullptr ) { 403 if ( part == nullptr ) {
404 G4cout << "G4HadronicBuilder::BuildDecay 404 G4cout << "G4HadronicBuilder::BuildDecayTableForBCHadrons : ERROR ! particlePDG="
405 << pdg << " is not defined !" << 405 << pdg << " is not defined !" << G4endl;
406 continue; 406 continue;
407 } 407 }
408 if ( part->GetDecayTable() ) { 408 if ( part->GetDecayTable() ) {
409 G4cout << "G4HadronicBuilder::BuildDecay 409 G4cout << "G4HadronicBuilder::BuildDecayTableForBCHadrons : WARNING ! particlePDG="
410 << pdg << " has already a decay t 410 << pdg << " has already a decay table defined !" << G4endl;
411 continue; 411 continue;
412 } 412 }
413 G4DecayTable* decayTable = new G4DecayTabl 413 G4DecayTable* decayTable = new G4DecayTable;
414 const G4int numberDecayChannels = 1; 414 const G4int numberDecayChannels = 1;
415 G4VDecayChannel** mode = new G4VDecayChann 415 G4VDecayChannel** mode = new G4VDecayChannel*[ numberDecayChannels ];
416 for ( G4int i = 0; i < numberDecayChannels 416 for ( G4int i = 0; i < numberDecayChannels; ++i ) mode[i] = nullptr;
417 switch ( pdg ) { 417 switch ( pdg ) {
418 // Charmed mesons 418 // Charmed mesons
419 case 411 : // D+ 419 case 411 : // D+
420 mode[0] = new G4PhaseSpaceDecayChannel 420 mode[0] = new G4PhaseSpaceDecayChannel( "D+", 1.0, 3, "kaon-", "pi+", "pi+" );
421 break; 421 break;
422 case -411 : // D- 422 case -411 : // D-
423 mode[0] = new G4PhaseSpaceDecayChannel 423 mode[0] = new G4PhaseSpaceDecayChannel( "D-", 1.0, 3, "kaon+", "pi-", "pi-" );
424 break; 424 break;
425 case 421 : // D0 425 case 421 : // D0
426 mode[0] = new G4PhaseSpaceDecayChannel 426 mode[0] = new G4PhaseSpaceDecayChannel( "D0", 1.0, 3, "kaon-", "pi+", "pi0" );
427 break; 427 break;
428 case -421 : // anti_D0 428 case -421 : // anti_D0
429 mode[0] = new G4PhaseSpaceDecayChannel 429 mode[0] = new G4PhaseSpaceDecayChannel( "anti_D0", 1.0, 3, "kaon+", "pi-", "pi0" );
430 break; 430 break;
431 case 431 : // Ds+ 431 case 431 : // Ds+
432 mode[0] = new G4PhaseSpaceDecayChannel 432 mode[0] = new G4PhaseSpaceDecayChannel( "Ds+", 1.0, 3, "kaon+", "kaon-", "pi+" );
433 break; 433 break;
434 case -431 : // Ds- 434 case -431 : // Ds-
435 mode[0] = new G4PhaseSpaceDecayChannel 435 mode[0] = new G4PhaseSpaceDecayChannel( "Ds-", 1.0, 3, "kaon-", "kaon+", "pi-" );
436 break; 436 break;
437 // Bottom mesons 437 // Bottom mesons
438 case 521 : // B+ 438 case 521 : // B+
439 mode[0] = new G4PhaseSpaceDecayChannel 439 mode[0] = new G4PhaseSpaceDecayChannel( "B+", 1.0, 3, "anti_D0", "pi+", "pi0" );
440 break; 440 break;
441 case -521 : // B- 441 case -521 : // B-
442 mode[0] = new G4PhaseSpaceDecayChannel 442 mode[0] = new G4PhaseSpaceDecayChannel( "B-", 1.0, 3, "D0", "pi-", "pi0" );
443 break; 443 break;
444 case 511 : // B0 444 case 511 : // B0
445 mode[0] = new G4PhaseSpaceDecayChannel 445 mode[0] = new G4PhaseSpaceDecayChannel( "B0", 1.0, 3, "D-", "pi+", "pi0" );
446 break; 446 break;
447 case -511 : // anti_B0 447 case -511 : // anti_B0
448 mode[0] = new G4PhaseSpaceDecayChannel 448 mode[0] = new G4PhaseSpaceDecayChannel( "anti_B0", 1.0, 3, "D+", "pi-", "pi0" );
449 break; 449 break;
450 case 531 : // Bs0 450 case 531 : // Bs0
451 mode[0] = new G4PhaseSpaceDecayChannel 451 mode[0] = new G4PhaseSpaceDecayChannel( "Bs0", 1.0, 3, "Ds-", "pi+", "pi0" );
452 break; 452 break;
453 case -531 : // anti_Bs0 453 case -531 : // anti_Bs0
454 mode[0] = new G4PhaseSpaceDecayChannel 454 mode[0] = new G4PhaseSpaceDecayChannel( "anti_Bs0", 1.0, 3, "Ds+", "pi-", "pi0" );
455 break; 455 break;
456 case 541 : // Bc+ 456 case 541 : // Bc+
457 mode[0] = new G4PhaseSpaceDecayChannel 457 mode[0] = new G4PhaseSpaceDecayChannel( "Bc+", 1.0, 2, "J/psi", "pi+" );
458 break; 458 break;
459 case -541 : // Bc- 459 case -541 : // Bc-
460 mode[0] = new G4PhaseSpaceDecayChannel 460 mode[0] = new G4PhaseSpaceDecayChannel( "Bc-", 1.0, 2, "J/psi", "pi-" );
461 break; 461 break;
462 // Charmed baryons (and anti-baryons) 462 // Charmed baryons (and anti-baryons)
463 case 4122 : // lambda_c+ 463 case 4122 : // lambda_c+
464 mode[0] = new G4PhaseSpaceDecayChannel 464 mode[0] = new G4PhaseSpaceDecayChannel( "lambda_c+", 1.0, 3, "proton", "kaon-", "pi+" );
465 break; 465 break;
466 case -4122 : // anti_lambda_c+ 466 case -4122 : // anti_lambda_c+
467 mode[0] = new G4PhaseSpaceDecayChannel 467 mode[0] = new G4PhaseSpaceDecayChannel( "anti_lambda_c+", 1.0, 3, "anti_proton", "kaon+", "pi-" );
468 break; 468 break;
469 case 4232 : // xi_c+ 469 case 4232 : // xi_c+
470 mode[0] = new G4PhaseSpaceDecayChannel 470 mode[0] = new G4PhaseSpaceDecayChannel( "xi_c+", 1.0, 3, "sigma+", "kaon-", "pi+" );
471 break; 471 break;
472 case -4232 : // anti_xi_c+ 472 case -4232 : // anti_xi_c+
473 mode[0] = new G4PhaseSpaceDecayChannel 473 mode[0] = new G4PhaseSpaceDecayChannel( "anti_xi_c+", 1.0, 3, "anti_sigma+", "kaon+", "pi-" );
474 break; 474 break;
475 case 4132 : // xi_c0 475 case 4132 : // xi_c0
476 mode[0] = new G4PhaseSpaceDecayChannel 476 mode[0] = new G4PhaseSpaceDecayChannel( "xi_c0", 1.0, 3, "lambda", "kaon-", "pi+" );
477 break; 477 break;
478 case -4132 : // anti_xi_c0 478 case -4132 : // anti_xi_c0
479 mode[0] = new G4PhaseSpaceDecayChannel 479 mode[0] = new G4PhaseSpaceDecayChannel( "anti_xi_c0", 1.0, 3, "anti_lambda", "kaon+", "pi-" );
480 break; 480 break;
481 case 4332 : // omega_c0 481 case 4332 : // omega_c0
482 mode[0] = new G4PhaseSpaceDecayChannel 482 mode[0] = new G4PhaseSpaceDecayChannel( "omega_c0", 1.0, 3, "xi0", "kaon-", "pi+" );
483 break; 483 break;
484 case -4332 : // anti_omega_c0 484 case -4332 : // anti_omega_c0
485 mode[0] = new G4PhaseSpaceDecayChannel 485 mode[0] = new G4PhaseSpaceDecayChannel( "anti_omega_c0", 1.0, 3, "anti_xi0", "kaon+", "pi-" );
486 break; 486 break;
487 // Bottom baryons (and anti-baryons) 487 // Bottom baryons (and anti-baryons)
488 case 5122 : // lambda_b 488 case 5122 : // lambda_b
489 mode[0] = new G4PhaseSpaceDecayChannel 489 mode[0] = new G4PhaseSpaceDecayChannel( "lambda_b", 1.0, 4, "lambda_c+", "pi+", "pi-", "pi-" );
490 break; 490 break;
491 case -5122 : // anti_lambda_b 491 case -5122 : // anti_lambda_b
492 mode[0] = new G4PhaseSpaceDecayChannel 492 mode[0] = new G4PhaseSpaceDecayChannel( "anti_lambda_b", 1.0, 4, "anti_lambda_c+", "pi-", "pi+", "pi+" );
493 break; 493 break;
494 case 5232 : // xi_b0 494 case 5232 : // xi_b0
495 mode[0] = new G4PhaseSpaceDecayChannel 495 mode[0] = new G4PhaseSpaceDecayChannel( "xi_b0", 1.0, 3, "lambda_c+", "kaon-", "pi0" );
496 break; 496 break;
497 case -5232 : // anti_xi_b0 497 case -5232 : // anti_xi_b0
498 mode[0] = new G4PhaseSpaceDecayChannel 498 mode[0] = new G4PhaseSpaceDecayChannel( "anti_xi_b0", 1.0, 3, "anti_lambda_c+", "kaon+", "pi0" );
499 break; 499 break;
500 case 5132 : // xi_b- 500 case 5132 : // xi_b-
501 mode[0] = new G4PhaseSpaceDecayChannel 501 mode[0] = new G4PhaseSpaceDecayChannel( "xi_b-", 1.0, 3, "lambda_c+", "kaon-", "pi-" );
502 break; 502 break;
503 case -5132 : // anti_xi_b- 503 case -5132 : // anti_xi_b-
504 mode[0] = new G4PhaseSpaceDecayChannel 504 mode[0] = new G4PhaseSpaceDecayChannel( "anti_xi_b-", 1.0, 3, "anti_lambda_c+", "kaon+", "pi+" );
505 break; 505 break;
506 case 5332 : // omega_b- 506 case 5332 : // omega_b-
507 mode[0] = new G4PhaseSpaceDecayChannel 507 mode[0] = new G4PhaseSpaceDecayChannel( "omega_b-", 1.0, 3, "xi_c+", "kaon-", "pi-" );
508 break; 508 break;
509 case -5332 : // anti_omega_b- 509 case -5332 : // anti_omega_b-
510 mode[0] = new G4PhaseSpaceDecayChannel 510 mode[0] = new G4PhaseSpaceDecayChannel( "anti_omega_b-", 1.0, 3, "anti_xi_c+", "kaon+", "pi+" );
511 break; 511 break;
512 default : 512 default :
513 G4cout << "G4HadronicBuilder::BuildDec 513 G4cout << "G4HadronicBuilder::BuildDecayTableForBCHadrons : UNKNOWN particlePDG=" << pdg << G4endl;
514 } // End of the switch 514 } // End of the switch
515 515
516 for ( G4int index = 0; index < numberDecay 516 for ( G4int index = 0; index < numberDecayChannels; ++index ) decayTable->Insert( mode[index] );
517 delete [] mode; 517 delete [] mode;
518 part->SetDecayTable( decayTable ); 518 part->SetDecayTable( decayTable );
519 } // End of the for loop over heavy hadrons 519 } // End of the for loop over heavy hadrons
520 // Add now the decay for etac, JPsi and Upsi 520 // Add now the decay for etac, JPsi and Upsilon because these can be produced as
521 // secondaries in hadronic interactions, whi 521 // secondaries in hadronic interactions, while they are not part of the heavy
522 // hadrons included in G4HadParticles::GetBC 522 // hadrons included in G4HadParticles::GetBCHadrons() because they live too shortly
523 // and therefore their hadronic interactions 523 // and therefore their hadronic interactions can be neglected (as we do for pi0 and sigma0).
524 if ( ! G4Etac::Definition()->GetDecayTable() 524 if ( ! G4Etac::Definition()->GetDecayTable() ) {
525 G4DecayTable* decayTable = new G4DecayTabl 525 G4DecayTable* decayTable = new G4DecayTable;
526 const G4int numberDecayChannels = 1; 526 const G4int numberDecayChannels = 1;
527 G4VDecayChannel** mode = new G4VDecayChann 527 G4VDecayChannel** mode = new G4VDecayChannel*[ numberDecayChannels ];
528 for ( G4int i = 0; i < numberDecayChannels 528 for ( G4int i = 0; i < numberDecayChannels; ++i ) mode[i] = nullptr;
529 mode[0] = new G4PhaseSpaceDecayChannel( "e 529 mode[0] = new G4PhaseSpaceDecayChannel( "etac", 1.0, 3, "eta", "pi+", "pi-" );
530 for ( G4int index = 0; index < numberDecay 530 for ( G4int index = 0; index < numberDecayChannels; ++index ) decayTable->Insert( mode[index] );
531 delete [] mode; 531 delete [] mode;
532 G4Etac::Definition()->SetDecayTable( decay 532 G4Etac::Definition()->SetDecayTable( decayTable );
533 } 533 }
534 if ( ! G4JPsi::Definition()->GetDecayTable() 534 if ( ! G4JPsi::Definition()->GetDecayTable() ) {
535 G4DecayTable* decayTable = new G4DecayTabl 535 G4DecayTable* decayTable = new G4DecayTable;
536 const G4int numberDecayChannels = 1; 536 const G4int numberDecayChannels = 1;
537 G4VDecayChannel** mode = new G4VDecayChann 537 G4VDecayChannel** mode = new G4VDecayChannel*[ numberDecayChannels ];
538 for ( G4int i = 0; i < numberDecayChannels 538 for ( G4int i = 0; i < numberDecayChannels; ++i ) mode[i] = nullptr;
539 mode[0] = new G4PhaseSpaceDecayChannel( "J 539 mode[0] = new G4PhaseSpaceDecayChannel( "J/psi", 1.0, 3, "pi0", "pi+", "pi-" );
540 for ( G4int index = 0; index < numberDecay 540 for ( G4int index = 0; index < numberDecayChannels; ++index ) decayTable->Insert( mode[index] );
541 delete [] mode; 541 delete [] mode;
542 G4JPsi::Definition()->SetDecayTable( decay 542 G4JPsi::Definition()->SetDecayTable( decayTable );
543 } 543 }
544 if ( ! G4Upsilon::Definition()->GetDecayTabl 544 if ( ! G4Upsilon::Definition()->GetDecayTable() ) {
545 G4DecayTable* decayTable = new G4DecayTabl 545 G4DecayTable* decayTable = new G4DecayTable;
546 const G4int numberDecayChannels = 1; 546 const G4int numberDecayChannels = 1;
547 G4VDecayChannel** mode = new G4VDecayChann 547 G4VDecayChannel** mode = new G4VDecayChannel*[ numberDecayChannels ];
548 for ( G4int i = 0; i < numberDecayChannels 548 for ( G4int i = 0; i < numberDecayChannels; ++i ) mode[i] = nullptr;
549 mode[0] = new G4PhaseSpaceDecayChannel( "U 549 mode[0] = new G4PhaseSpaceDecayChannel( "Upsilon", 1.0, 3, "eta_prime", "pi+", "pi-" );
550 for ( G4int index = 0; index < numberDecay 550 for ( G4int index = 0; index < numberDecayChannels; ++index ) decayTable->Insert( mode[index] );
551 delete [] mode; 551 delete [] mode;
552 G4Upsilon::Definition()->SetDecayTable( de 552 G4Upsilon::Definition()->SetDecayTable( decayTable );
553 } 553 }
554 } 554 }
555 555
556 556
557 void G4HadronicBuilder::BuildHyperNucleiFTFP_B 557 void G4HadronicBuilder::BuildHyperNucleiFTFP_BERT() {
558 if ( G4HadronicParameters::Instance()->Enabl 558 if ( G4HadronicParameters::Instance()->EnableHyperNuclei() ) {
559 // Bertini intra-nuclear cascade model is 559 // Bertini intra-nuclear cascade model is currently not applicable for light
560 // hypernuclei, therefore FTFP is used dow 560 // hypernuclei, therefore FTFP is used down to zero kinetic energy (but at
561 // very low energies, a dummy model is use 561 // very low energies, a dummy model is used that simply returns the projectile
562 // hypernucleus in the final state). 562 // hypernucleus in the final state).
563 BuildFTFP_BERT( G4HadParticles::GetHyperNu 563 BuildFTFP_BERT( G4HadParticles::GetHyperNuclei(), false, "Glauber-Gribov" );
564 } 564 }
565 } 565 }
566 566
567 567
568 void G4HadronicBuilder::BuildHyperAntiNucleiFT 568 void G4HadronicBuilder::BuildHyperAntiNucleiFTFP_BERT() {
569 if ( G4HadronicParameters::Instance()->Enabl 569 if ( G4HadronicParameters::Instance()->EnableHyperNuclei() ) {
570 // FTFP can be used down to zero kinetic e 570 // FTFP can be used down to zero kinetic energy.
571 BuildFTFP_BERT( G4HadParticles::GetHyperAn 571 BuildFTFP_BERT( G4HadParticles::GetHyperAntiNuclei(), false, "AntiAGlauber" );
572 } 572 }
573 } 573 }
574 574
575 575
576 void G4HadronicBuilder::BuildHyperNucleiFTFP_I 576 void G4HadronicBuilder::BuildHyperNucleiFTFP_INCLXX() {
577 if ( G4HadronicParameters::Instance()->Enabl 577 if ( G4HadronicParameters::Instance()->EnableHyperNuclei() ) {
578 BuildFTFP_INCLXX( G4HadParticles::GetHyper 578 BuildFTFP_INCLXX( G4HadParticles::GetHyperNuclei(), "Glauber-Gribov" );
579 } 579 }
580 } 580 }
581 581
582 582
583 void G4HadronicBuilder::BuildFTFP_INCLXX( cons 583 void G4HadronicBuilder::BuildFTFP_INCLXX( const std::vector< G4int >& partList, const G4String& xsName ) {
584 G4HadronicParameters* param = G4HadronicPara 584 G4HadronicParameters* param = G4HadronicParameters::Instance();
585 G4PhysicsListHelper* ph = G4PhysicsListHelpe 585 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
586 auto theTheoFSModel = new G4TheoFSGenerator( 586 auto theTheoFSModel = new G4TheoFSGenerator( "FTFP" );
587 auto theStringModel = new G4FTFModel; 587 auto theStringModel = new G4FTFModel;
588 theStringModel->SetFragmentationModel( new G 588 theStringModel->SetFragmentationModel( new G4ExcitedStringDecay );
589 theTheoFSModel->SetHighEnergyGenerator( theS 589 theTheoFSModel->SetHighEnergyGenerator( theStringModel );
590 theTheoFSModel->SetTransport( new G4Generato 590 theTheoFSModel->SetTransport( new G4GeneratorPrecompoundInterface );
591 theTheoFSModel->SetMaxEnergy( param->GetMaxE 591 theTheoFSModel->SetMaxEnergy( param->GetMaxEnergy() );
592 theTheoFSModel->SetMinEnergy( 15.0*CLHEP::Ge 592 theTheoFSModel->SetMinEnergy( 15.0*CLHEP::GeV );
593 G4VPreCompoundModel* thePrecoModel = new G4P 593 G4VPreCompoundModel* thePrecoModel = new G4PreCompoundModel;
594 thePrecoModel->SetMinEnergy( 0.0 ); 594 thePrecoModel->SetMinEnergy( 0.0 );
595 thePrecoModel->SetMaxEnergy( 2.0*CLHEP::MeV 595 thePrecoModel->SetMaxEnergy( 2.0*CLHEP::MeV );
596 G4INCLXXInterface* theINCLXXModel = new G4IN 596 G4INCLXXInterface* theINCLXXModel = new G4INCLXXInterface( thePrecoModel );
597 theINCLXXModel->SetMinEnergy( 1.0*CLHEP::MeV 597 theINCLXXModel->SetMinEnergy( 1.0*CLHEP::MeV );
598 theINCLXXModel->SetMaxEnergy( 20.0*CLHEP::Ge 598 theINCLXXModel->SetMaxEnergy( 20.0*CLHEP::GeV );
599 auto xsinel = G4HadProcesses::InelasticXS( x 599 auto xsinel = G4HadProcesses::InelasticXS( xsName );
600 G4ParticleTable* table = G4ParticleTable::Ge 600 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
601 for ( auto & pdg : partList ) { 601 for ( auto & pdg : partList ) {
602 auto part = table->FindParticle( pdg ); 602 auto part = table->FindParticle( pdg );
603 if ( part == nullptr ) continue; 603 if ( part == nullptr ) continue;
604 auto hadi = new G4HadronInelasticProcess( 604 auto hadi = new G4HadronInelasticProcess( part->GetParticleName()+"Inelastic", part );
605 hadi->AddDataSet( xsinel ); 605 hadi->AddDataSet( xsinel );
606 hadi->RegisterMe( theTheoFSModel ); 606 hadi->RegisterMe( theTheoFSModel );
607 hadi->RegisterMe( theINCLXXModel ); 607 hadi->RegisterMe( theINCLXXModel );
608 if ( param->ApplyFactorXS() ) hadi->Multip 608 if ( param->ApplyFactorXS() ) hadi->MultiplyCrossSectionBy( param->XSFactorHadronInelastic() );
609 ph->RegisterProcess( hadi, part ); 609 ph->RegisterProcess( hadi, part );
610 } 610 }
611 } 611 }
612 612