source: trunk/source/physics_lists/lists/src/G4HadronInelasticQLHEP.cc @ 1315

Last change on this file since 1315 was 1315, checked in by garnier, 14 years ago

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08
File size: 15.6 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: G4HadronInelasticQLHEP.cc,v 1.4 2010/06/08 08:58:03 gunter Exp $
27// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
28//
29//---------------------------------------------------------------------------
30//
31// ClassName:   G4HadronInelasticQLHEP
32//
33// Author: 11 April 2006 V. Ivanchenko
34//
35// Modified:
36// 05.07.2006 V.Ivanchenko fix problem of initialisation of HP
37//
38//----------------------------------------------------------------------------
39//
40
41#include "G4HadronInelasticQLHEP.hh"
42
43#include "G4HadronInelasticProcess.hh"
44
45#include "G4ParticleDefinition.hh"
46#include "G4ProcessManager.hh"
47
48#include "G4MesonConstructor.hh"
49#include "G4BaryonConstructor.hh"
50
51#include "G4PiNuclearCrossSection.hh"
52
53#include "G4TheoFSGenerator.hh"
54#include "G4GeneratorPrecompoundInterface.hh"
55#include "G4PreCompoundModel.hh"
56#include "G4ExcitationHandler.hh"
57#include "G4QGSMFragmentation.hh"
58#include "G4ExcitedStringDecay.hh"
59
60#include "G4CascadeInterface.hh"
61#include "G4BinaryCascade.hh"
62#include "G4HadronFissionProcess.hh"
63#include "G4HadronCaptureProcess.hh"
64#include "G4LFission.hh"
65#include "G4LCapture.hh"
66
67#include "G4NeutronHPInelastic.hh"
68#include "G4NeutronHPFission.hh"
69#include "G4NeutronHPCapture.hh"
70
71#include "G4LEAntiLambdaInelastic.hh"
72#include "G4LEAntiNeutronInelastic.hh"
73#include "G4LEAntiOmegaMinusInelastic.hh"
74#include "G4LEAntiProtonInelastic.hh"
75#include "G4LEAntiSigmaMinusInelastic.hh"
76#include "G4LEAntiSigmaPlusInelastic.hh"
77#include "G4LEAntiXiMinusInelastic.hh"
78#include "G4LEAntiXiZeroInelastic.hh"
79#include "G4LEKaonMinusInelastic.hh"
80#include "G4LEKaonPlusInelastic.hh"
81#include "G4LEKaonZeroSInelastic.hh"
82#include "G4LEKaonZeroLInelastic.hh"
83#include "G4LENeutronInelastic.hh"
84#include "G4LELambdaInelastic.hh"
85#include "G4LEProtonInelastic.hh"
86#include "G4LEPionPlusInelastic.hh"
87#include "G4LEPionMinusInelastic.hh"
88#include "G4LEOmegaMinusInelastic.hh"
89#include "G4LESigmaMinusInelastic.hh"
90#include "G4LESigmaPlusInelastic.hh"
91#include "G4LEXiMinusInelastic.hh"
92#include "G4LEXiZeroInelastic.hh"
93
94#include "G4HEAntiLambdaInelastic.hh"
95#include "G4HEAntiNeutronInelastic.hh"
96#include "G4HEAntiOmegaMinusInelastic.hh"
97#include "G4HEAntiProtonInelastic.hh"
98#include "G4HEAntiSigmaMinusInelastic.hh"
99#include "G4HEAntiSigmaPlusInelastic.hh"
100#include "G4HEAntiXiMinusInelastic.hh"
101#include "G4HEAntiXiZeroInelastic.hh"
102#include "G4HEKaonMinusInelastic.hh"
103#include "G4HEKaonPlusInelastic.hh"
104#include "G4HEKaonZeroInelastic.hh"
105#include "G4HENeutronInelastic.hh"
106#include "G4HELambdaInelastic.hh"
107#include "G4HEProtonInelastic.hh"
108#include "G4HEPionPlusInelastic.hh"
109#include "G4HEPionMinusInelastic.hh"
110#include "G4HEOmegaMinusInelastic.hh"
111#include "G4HESigmaMinusInelastic.hh"
112#include "G4HESigmaPlusInelastic.hh"
113#include "G4HEXiMinusInelastic.hh"
114#include "G4HEXiZeroInelastic.hh"
115
116G4HadronInelasticQLHEP::G4HadronInelasticQLHEP(G4int ver)
117  : G4VPhysicsConstructor("hInelasticQLHEP"), verbose(ver), qgsFlag(false), 
118    bertFlag(false), bicFlag(false), hpFlag(false), wasActivated(false)
119{
120  if(verbose > 1) G4cout << "### HadronInelasticQLHEP" << G4endl;
121  theCascade = 0;
122  theQGStringDecay = 0;
123  theQGStringModel = 0;
124  thePreEquilib = 0;
125  theHPXSecI = 0;
126  theHPXSecC = 0;
127  theHPXSecF = 0;
128}
129
130G4HadronInelasticQLHEP::G4HadronInelasticQLHEP(const G4String&, 
131    G4int ver, G4bool qgs, G4bool bert, G4bool bic, G4bool hp)
132  : G4VPhysicsConstructor("hInelasticQLHEP"), verbose(ver), qgsFlag(qgs), 
133    bertFlag(bert), bicFlag(bic), hpFlag(hp), wasActivated(false)
134{
135  if(verbose > 1) G4cout << "### HadronInelasticQLHEP" << G4endl;
136  theCascade = 0;
137  theQGStringDecay = 0;
138  theQGStringModel = 0;
139  thePreEquilib = 0;
140  theHPXSecI = 0;
141  theHPXSecC = 0;
142  theHPXSecF = 0;
143}
144
145G4HadronInelasticQLHEP::~G4HadronInelasticQLHEP()
146{
147  if(wasActivated) {
148    delete theCascade;
149    delete theQGStringDecay;
150    delete theQGStringModel;
151    delete thePreEquilib;
152    delete theHPXSecI;
153    delete theHPXSecC;
154    delete theHPXSecF;
155  }
156}
157
158void G4HadronInelasticQLHEP::ConstructParticle()
159{
160  G4MesonConstructor pMesonConstructor;
161  pMesonConstructor.ConstructParticle();
162
163  G4BaryonConstructor pBaryonConstructor;
164  pBaryonConstructor.ConstructParticle();
165}
166
167void G4HadronInelasticQLHEP::ConstructProcess()
168{
169  if(wasActivated) return;
170  wasActivated = true;
171
172  if(verbose > 1) G4cout << "### HadronInelasticQLHEP Construct Process" << G4endl;
173
174  G4double minEneutron   = 0.0*GeV;
175  G4double minELEP       = 0.0*GeV;
176  G4double maxELEP       = 55.*GeV;
177  G4double minEstring    = 12.*GeV;
178  if(hpFlag) minEneutron = 19.5*MeV;
179
180  //Bertini
181  G4HadronicInteraction* theBERT = 0;
182  if(bertFlag) {
183    minELEP = 9.5*GeV;
184    theBERT = new G4CascadeInterface();
185    theBERT->SetMinEnergy(0.0);
186    theBERT->SetMaxEnergy(9.9*GeV);
187  }
188
189  //Binari
190  G4HadronicInteraction* theBIC = 0;
191  if(bicFlag) {
192    minELEP = 9.5*GeV;
193    theBIC = new G4BinaryCascade();
194    theBIC->SetMinEnergy(0.0);
195    theBIC->SetMaxEnergy(9.9*GeV);
196  }
197
198  //QGSP
199  G4TheoFSGenerator* theQGSModel = 0;
200  if(qgsFlag) {
201    maxELEP     = 25.*GeV;
202    theQGSModel = new G4TheoFSGenerator();
203    theQGStringModel  = new G4QGSModel< G4QGSParticipants >;
204    theQGStringDecay  = new G4ExcitedStringDecay(new G4QGSMFragmentation());
205    theQGStringModel->SetFragmentationModel(theQGStringDecay);
206    theCascade = new G4GeneratorPrecompoundInterface;
207    thePreEquilib = new G4PreCompoundModel(new G4ExcitationHandler);
208    theCascade->SetDeExcitation(thePreEquilib);
209    theQGSModel->SetTransport(theCascade);
210    theQGSModel->SetHighEnergyGenerator(theQGStringModel);
211    theQGSModel->SetMinEnergy(minEstring);
212    theQGSModel->SetMaxEnergy(100*TeV);
213  }
214
215  theParticleIterator->reset();
216  while( (*theParticleIterator)() ) {
217    G4ParticleDefinition* particle = theParticleIterator->value();
218    G4String pname = particle->GetParticleName();
219    if(verbose > 1) G4cout << "### HadronInelasticQLHEP:  " << pname << G4endl;
220    if(pname == "anti_lambda"  ||
221       pname == "anti_neutron" ||
222       pname == "anti_omega-"  || 
223       pname == "anti_proton"  || 
224       pname == "anti_sigma-"  || 
225       pname == "anti_sigma+"  || 
226       pname == "anti_xi-"  || 
227       pname == "anti_xi0"  || 
228       pname == "kaon-"     || 
229       pname == "kaon+"     || 
230       pname == "kaon0S"    || 
231       pname == "kaon0L"    || 
232       pname == "lambda"    || 
233       pname == "neutron"   || 
234       pname == "omega-"    || 
235       pname == "pi-"       || 
236       pname == "pi+"       || 
237       pname == "proton"    || 
238       pname == "sigma-"    || 
239       pname == "sigma+"    || 
240       pname == "xi-"       || 
241       pname == "xi0") {
242
243      G4ProcessManager* pmanager = particle->GetProcessManager();
244      G4HadronInelasticProcess* hp = 
245        new G4HadronInelasticProcess("hInelastic", particle);
246      pmanager->AddDiscreteProcess(hp);
247
248      if(pname == "proton") {
249        if(qgsFlag) {
250          Register(particle,hp,theQGSModel,"QGS");
251          hp->AddDataSet(&theXSecP);
252        } else {
253          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
254        }
255        AddLEP(particle, hp, minELEP, maxELEP);
256       
257        if(bicFlag)       Register(particle,hp,theBIC,"Binary");
258        else if(bertFlag) Register(particle,hp,theBERT,"Bertini");
259
260      } else if(pname == "neutron") {
261        if(qgsFlag) {
262          Register(particle,hp,theQGSModel,"QGS");
263          hp->AddDataSet(&theXSecN);
264        } else {
265          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
266        }
267        AddLEP(particle, hp, minELEP, maxELEP);
268
269        G4HadronCaptureProcess* theNeutronCapture = 
270          new G4HadronCaptureProcess("nCapture");
271        G4HadronFissionProcess* theNeutronFission = 
272          new G4HadronFissionProcess("nFission");
273        pmanager->AddDiscreteProcess(theNeutronCapture);
274        pmanager->AddDiscreteProcess(theNeutronFission);
275
276        if(hpFlag) {
277          theHPXSecI = new G4NeutronHPInelasticData;
278          theHPXSecC = new G4NeutronHPCaptureData;
279          theHPXSecF = new G4NeutronHPFissionData;
280          hp->AddDataSet(theHPXSecI);
281          theNeutronCapture->AddDataSet(theHPXSecC);
282          theNeutronFission->AddDataSet(theHPXSecF);
283          G4NeutronHPInelastic* hpi = new G4NeutronHPInelastic();
284          G4NeutronHPCapture* hpc = new G4NeutronHPCapture();
285          G4NeutronHPFission* hpf = new G4NeutronHPFission();
286          Register(particle,hp,hpi,"HP");
287          Register(particle,theNeutronCapture,hpc,"HP");
288          Register(particle,theNeutronFission,hpf,"HP");
289        }
290
291        G4HadronicInteraction* theC = new G4LCapture();
292        theC->SetMinEnergy(minEneutron);
293        Register(particle,theNeutronCapture,theC,"LCapture");
294
295        G4HadronicInteraction* theF = new G4LFission();
296        theF->SetMinEnergy(minEneutron);
297        Register(particle,theNeutronFission,theF,"LFission");
298
299        if(bicFlag) {
300          G4BinaryCascade* theB = new G4BinaryCascade();
301          theB->SetMinEnergy(minEneutron);
302          theB->SetMaxEnergy(9.9*GeV);
303          Register(particle,hp,theB,"Binary");
304        } else if(bertFlag) {
305          G4CascadeInterface* theB = new G4CascadeInterface();
306          theB->SetMinEnergy(minEneutron);
307          theB->SetMaxEnergy(9.9*GeV);
308          Register(particle,hp,theB,"Bertini");
309        }
310
311      } else if(pname == "pi-" || pname == "pi+") {
312        if(qgsFlag) {
313          Register(particle,hp,theQGSModel,"QGS");
314          hp->AddDataSet(&thePiCross);
315        } else {
316          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
317        }
318        AddLEP(particle, hp, minELEP, maxELEP);
319        if(bertFlag) Register(particle,hp,theBERT,"Bertini");
320
321      } else if(pname == "kaon-"     || 
322                pname == "kaon+"     || 
323                pname == "kaon0S"    || 
324                pname == "kaon0L") {
325
326        if(qgsFlag) Register(particle,hp,theQGSModel,"QGS");
327        else        AddHEP(particle, hp, 25.*GeV, 100.*TeV);
328       
329        AddLEP(particle, hp, minELEP, maxELEP);
330        if(bertFlag) Register(particle,hp,theBERT,"Bertini");
331
332      } else {
333
334        AddHEP(particle, hp, 25.*GeV, 100.*TeV);
335        AddLEP(particle, hp, 0.0, 55.*GeV);
336      }
337
338      if(verbose > 1)
339        G4cout << "### HadronInelasticQLHEP: " << hp->GetProcessName()
340               << " added for " << pname << G4endl;
341    }
342  }
343}
344
345void G4HadronInelasticQLHEP::AddLEP(G4ParticleDefinition* particle,
346                                    G4HadronicProcess* hp,
347                                    G4double emin,
348                                    G4double emax)
349{
350  G4HadronicInteraction* hi = 0;
351  G4String pname = particle->GetParticleName();
352
353  if(pname == "anti_lambda"      ) hi = new G4LEAntiLambdaInelastic();
354  else if(pname == "anti_neutron") hi = new G4LEAntiNeutronInelastic();
355  else if(pname == "anti_omega-" ) hi = new G4LEAntiOmegaMinusInelastic();
356  else if(pname == "anti_proton" ) hi = new G4LEAntiProtonInelastic();
357  else if(pname == "anti_sigma-" ) hi = new G4LEAntiSigmaMinusInelastic();
358  else if(pname == "anti_sigma+" ) hi = new G4LEAntiSigmaPlusInelastic();
359  else if(pname == "anti_xi-"    ) hi = new G4LEAntiXiMinusInelastic();
360  else if(pname == "anti_xi0"    ) hi = new G4LEAntiXiZeroInelastic();
361  else if(pname == "kaon-"       ) hi = new G4LEKaonMinusInelastic();
362  else if(pname == "kaon+"       ) hi = new G4LEKaonPlusInelastic();
363  else if(pname == "kaon0S"      ) hi = new G4LEKaonZeroSInelastic();
364  else if(pname == "kaon0L"      ) hi = new G4LEKaonZeroLInelastic();
365  else if(pname == "neutron"     ) hi = new G4LENeutronInelastic();
366  else if(pname == "lambda"      ) hi = new G4LELambdaInelastic();
367  else if(pname == "omega-"      ) hi = new G4LEOmegaMinusInelastic();
368  else if(pname == "proton"      ) hi = new G4LEProtonInelastic();
369  else if(pname == "pi+"         ) hi = new G4LEPionPlusInelastic();
370  else if(pname == "pi-"         ) hi = new G4LEPionMinusInelastic();
371  else if(pname == "sigma-"      ) hi = new G4LESigmaMinusInelastic();
372  else if(pname == "sigma+"      ) hi = new G4LESigmaPlusInelastic();
373  else if(pname == "xi-"         ) hi = new G4LEXiMinusInelastic();
374  else if(pname == "xi0"         ) hi = new G4LEXiZeroInelastic();
375
376  if(hi) {
377    hi->SetMinEnergy(emin);
378    hi->SetMaxEnergy(emax);
379    Register(particle,hp,hi,"LHEP");
380  } else {
381    G4cout << "### G4HadronInelasticTHEO: ERROR - no LHEP model for "
382           << pname << G4endl;
383  }
384}
385
386
387void G4HadronInelasticQLHEP::AddHEP(G4ParticleDefinition* particle,
388                                    G4HadronicProcess* hp,
389                                    G4double emin,
390                                    G4double emax)
391{
392  G4HadronicInteraction* hi = 0;
393  G4String pname = particle->GetParticleName();
394
395  if(pname == "anti_lambda"      ) hi = new G4HEAntiLambdaInelastic();
396  else if(pname == "anti_neutron") hi = new G4HEAntiNeutronInelastic();
397  else if(pname == "anti_omega-" ) hi = new G4HEAntiOmegaMinusInelastic();
398  else if(pname == "anti_proton" ) hi = new G4HEAntiProtonInelastic();
399  else if(pname == "anti_sigma-" ) hi = new G4HEAntiSigmaMinusInelastic();
400  else if(pname == "anti_sigma+" ) hi = new G4HEAntiSigmaPlusInelastic();
401  else if(pname == "anti_xi-"    ) hi = new G4HEAntiXiMinusInelastic();
402  else if(pname == "anti_xi0"    ) hi = new G4HEAntiXiZeroInelastic();
403  else if(pname == "kaon-"       ) hi = new G4HEKaonMinusInelastic();
404  else if(pname == "kaon+"       ) hi = new G4HEKaonPlusInelastic();
405  else if(pname == "kaon0S"      ) hi = new G4HEKaonZeroInelastic();
406  else if(pname == "kaon0L"      ) hi = new G4HEKaonZeroInelastic();
407  else if(pname == "neutron"     ) hi = new G4HENeutronInelastic();
408  else if(pname == "lambda"      ) hi = new G4HELambdaInelastic();
409  else if(pname == "omega-"      ) hi = new G4HEOmegaMinusInelastic();
410  else if(pname == "proton"      ) hi = new G4HEProtonInelastic();
411  else if(pname == "pi+"         ) hi = new G4HEPionPlusInelastic();
412  else if(pname == "pi-"         ) hi = new G4HEPionMinusInelastic();
413  else if(pname == "sigma-"      ) hi = new G4HESigmaMinusInelastic();
414  else if(pname == "sigma+"      ) hi = new G4HESigmaPlusInelastic();
415  else if(pname == "xi-"         ) hi = new G4HEXiMinusInelastic();
416  else if(pname == "xi0"         ) hi = new G4HEXiZeroInelastic();
417
418  if(hi) {
419    hi->SetMinEnergy(emin);
420    hi->SetMaxEnergy(emax);
421    Register(particle,hp,hi,"HEP");
422  } else {
423    G4cout << "### G4HadronInelasticQLHEP: ERROR - no HEP model for "
424           << pname << G4endl;
425  }
426}
427
428void G4HadronInelasticQLHEP::Register(G4ParticleDefinition* p, 
429                                      G4HadronicProcess* hp, 
430                                      G4HadronicInteraction* hi, 
431                                      const G4String& m)
432{
433  hp->RegisterMe(hi);
434  if(verbose > 1)
435    G4cout << "### QLHEP: Register new model " << m
436           << " for " << p->GetParticleName() << " and " 
437           << hp->GetProcessName()
438           << " E(GeV) " << hi->GetMinEnergy()/GeV
439           << " - " << hi->GetMaxEnergy()/GeV << G4endl;
440}
Note: See TracBrowser for help on using the repository browser.