Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/physics_lists/builders/src/G4HadronicBuilder.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

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