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

Last change on this file since 1202 was 988, checked in by garnier, 15 years ago

fichiers manquants

File size: 15.2 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.2 2008/05/19 10:21:34 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-02-ref-02 $
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(const G4String& name, 
117    G4int ver, G4bool qgs, G4bool bert, G4bool bic, G4bool hp)
118  : G4VPhysicsConstructor(name), verbose(ver), qgsFlag(qgs), 
119    bertFlag(bert), bicFlag(bic), hpFlag(hp), wasActivated(false)
120{
121  if(verbose > 1) G4cout << "### HadronInelasticQLHEP" << G4endl;
122  theCascade = 0;
123  theQGStringDecay = 0;
124  theQGStringModel = 0;
125  thePreEquilib = 0;
126  theHPXSecI = 0;
127  theHPXSecC = 0;
128  theHPXSecF = 0;
129}
130
131G4HadronInelasticQLHEP::~G4HadronInelasticQLHEP()
132{
133  if(wasActivated) {
134    delete theCascade;
135    delete theQGStringDecay;
136    delete theQGStringModel;
137    delete thePreEquilib;
138    delete theHPXSecI;
139    delete theHPXSecC;
140    delete theHPXSecF;
141  }
142}
143
144void G4HadronInelasticQLHEP::ConstructParticle()
145{
146  G4MesonConstructor pMesonConstructor;
147  pMesonConstructor.ConstructParticle();
148
149  G4BaryonConstructor pBaryonConstructor;
150  pBaryonConstructor.ConstructParticle();
151}
152
153void G4HadronInelasticQLHEP::ConstructProcess()
154{
155  if(wasActivated) return;
156  wasActivated = true;
157
158  if(verbose > 1) G4cout << "### HadronInelasticQLHEP Construct Process" << G4endl;
159
160  G4double minEneutron   = 0.0*GeV;
161  G4double minELEP       = 0.0*GeV;
162  G4double maxELEP       = 55.*GeV;
163  G4double minEstring    = 12.*GeV;
164  if(hpFlag) minEneutron = 19.5*MeV;
165
166  //Bertini
167  G4HadronicInteraction* theBERT = 0;
168  if(bertFlag) {
169    minELEP = 9.5*GeV;
170    theBERT = new G4CascadeInterface();
171    theBERT->SetMinEnergy(0.0);
172    theBERT->SetMaxEnergy(9.9*GeV);
173  }
174
175  //Binari
176  G4HadronicInteraction* theBIC = 0;
177  if(bicFlag) {
178    minELEP = 9.5*GeV;
179    theBIC = new G4BinaryCascade();
180    theBIC->SetMinEnergy(0.0);
181    theBIC->SetMaxEnergy(9.9*GeV);
182  }
183
184  //QGSP
185  G4TheoFSGenerator* theQGSModel = 0;
186  if(qgsFlag) {
187    maxELEP     = 25.*GeV;
188    theQGSModel = new G4TheoFSGenerator();
189    theQGStringModel  = new G4QGSModel< G4QGSParticipants >;
190    theQGStringDecay  = new G4ExcitedStringDecay(new G4QGSMFragmentation());
191    theQGStringModel->SetFragmentationModel(theQGStringDecay);
192    theCascade = new G4GeneratorPrecompoundInterface;
193    thePreEquilib = new G4PreCompoundModel(new G4ExcitationHandler);
194    theCascade->SetDeExcitation(thePreEquilib);
195    theQGSModel->SetTransport(theCascade);
196    theQGSModel->SetHighEnergyGenerator(theQGStringModel);
197    theQGSModel->SetMinEnergy(minEstring);
198    theQGSModel->SetMaxEnergy(100*TeV);
199  }
200
201  theParticleIterator->reset();
202  while( (*theParticleIterator)() ) {
203    G4ParticleDefinition* particle = theParticleIterator->value();
204    G4String pname = particle->GetParticleName();
205    if(verbose > 1) G4cout << "### HadronInelasticQLHEP:  " << pname << G4endl;
206    if(pname == "anti_lambda"  ||
207       pname == "anti_neutron" ||
208       pname == "anti_omega-"  || 
209       pname == "anti_proton"  || 
210       pname == "anti_sigma-"  || 
211       pname == "anti_sigma+"  || 
212       pname == "anti_xi-"  || 
213       pname == "anti_xi0"  || 
214       pname == "kaon-"     || 
215       pname == "kaon+"     || 
216       pname == "kaon0S"    || 
217       pname == "kaon0L"    || 
218       pname == "lambda"    || 
219       pname == "neutron"   || 
220       pname == "omega-"    || 
221       pname == "pi-"       || 
222       pname == "pi+"       || 
223       pname == "proton"    || 
224       pname == "sigma-"    || 
225       pname == "sigma+"    || 
226       pname == "xi-"       || 
227       pname == "xi0") {
228
229      G4ProcessManager* pmanager = particle->GetProcessManager();
230      G4HadronInelasticProcess* hp = 
231        new G4HadronInelasticProcess("hInelastic", particle);
232      pmanager->AddDiscreteProcess(hp);
233
234      if(pname == "proton") {
235        if(qgsFlag) {
236          Register(particle,hp,theQGSModel,"QGS");
237          hp->AddDataSet(&theXSecP);
238        } else {
239          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
240        }
241        AddLEP(particle, hp, minELEP, maxELEP);
242       
243        if(bicFlag)       Register(particle,hp,theBIC,"Binary");
244        else if(bertFlag) Register(particle,hp,theBERT,"Bertini");
245
246      } else if(pname == "neutron") {
247        if(qgsFlag) {
248          Register(particle,hp,theQGSModel,"QGS");
249          hp->AddDataSet(&theXSecN);
250        } else {
251          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
252        }
253        AddLEP(particle, hp, minELEP, maxELEP);
254
255        G4HadronCaptureProcess* theNeutronCapture = 
256          new G4HadronCaptureProcess("nCapture");
257        G4HadronFissionProcess* theNeutronFission = 
258          new G4HadronFissionProcess("nFission");
259        pmanager->AddDiscreteProcess(theNeutronCapture);
260        pmanager->AddDiscreteProcess(theNeutronFission);
261
262        if(hpFlag) {
263          theHPXSecI = new G4NeutronHPInelasticData;
264          theHPXSecC = new G4NeutronHPCaptureData;
265          theHPXSecF = new G4NeutronHPFissionData;
266          hp->AddDataSet(theHPXSecI);
267          theNeutronCapture->AddDataSet(theHPXSecC);
268          theNeutronFission->AddDataSet(theHPXSecF);
269          G4NeutronHPInelastic* hpi = new G4NeutronHPInelastic();
270          G4NeutronHPCapture* hpc = new G4NeutronHPCapture();
271          G4NeutronHPFission* hpf = new G4NeutronHPFission();
272          Register(particle,hp,hpi,"HP");
273          Register(particle,theNeutronCapture,hpc,"HP");
274          Register(particle,theNeutronFission,hpf,"HP");
275        }
276
277        G4HadronicInteraction* theC = new G4LCapture();
278        theC->SetMinEnergy(minEneutron);
279        Register(particle,theNeutronCapture,theC,"LCapture");
280
281        G4HadronicInteraction* theF = new G4LFission();
282        theF->SetMinEnergy(minEneutron);
283        Register(particle,theNeutronFission,theF,"LFission");
284
285        if(bicFlag) {
286          G4BinaryCascade* theB = new G4BinaryCascade();
287          theB->SetMinEnergy(minEneutron);
288          theB->SetMaxEnergy(9.9*GeV);
289          Register(particle,hp,theB,"Binary");
290        } else if(bertFlag) {
291          G4CascadeInterface* theB = new G4CascadeInterface();
292          theB->SetMinEnergy(minEneutron);
293          theB->SetMaxEnergy(9.9*GeV);
294          Register(particle,hp,theB,"Bertini");
295        }
296
297      } else if(pname == "pi-" || pname == "pi+") {
298        if(qgsFlag) {
299          Register(particle,hp,theQGSModel,"QGS");
300          hp->AddDataSet(&thePiCross);
301        } else {
302          AddHEP(particle, hp, 25.*GeV, 100.*TeV);
303        }
304        AddLEP(particle, hp, minELEP, maxELEP);
305        if(bertFlag) Register(particle,hp,theBERT,"Bertini");
306
307      } else if(pname == "kaon-"     || 
308                pname == "kaon+"     || 
309                pname == "kaon0S"    || 
310                pname == "kaon0L") {
311
312        if(qgsFlag) Register(particle,hp,theQGSModel,"QGS");
313        else        AddHEP(particle, hp, 25.*GeV, 100.*TeV);
314       
315        AddLEP(particle, hp, minELEP, maxELEP);
316        if(bertFlag) Register(particle,hp,theBERT,"Bertini");
317
318      } else {
319
320        AddHEP(particle, hp, 25.*GeV, 100.*TeV);
321        AddLEP(particle, hp, 0.0, 55.*GeV);
322      }
323
324      if(verbose > 1)
325        G4cout << "### HadronInelasticQLHEP: " << hp->GetProcessName()
326               << " added for " << pname << G4endl;
327    }
328  }
329}
330
331void G4HadronInelasticQLHEP::AddLEP(G4ParticleDefinition* particle,
332                                    G4HadronicProcess* hp,
333                                    G4double emin,
334                                    G4double emax)
335{
336  G4HadronicInteraction* hi = 0;
337  G4String pname = particle->GetParticleName();
338
339  if(pname == "anti_lambda"      ) hi = new G4LEAntiLambdaInelastic();
340  else if(pname == "anti_neutron") hi = new G4LEAntiNeutronInelastic();
341  else if(pname == "anti_omega-" ) hi = new G4LEAntiOmegaMinusInelastic();
342  else if(pname == "anti_proton" ) hi = new G4LEAntiProtonInelastic();
343  else if(pname == "anti_sigma-" ) hi = new G4LEAntiSigmaMinusInelastic();
344  else if(pname == "anti_sigma+" ) hi = new G4LEAntiSigmaPlusInelastic();
345  else if(pname == "anti_xi-"    ) hi = new G4LEAntiXiMinusInelastic();
346  else if(pname == "anti_xi0"    ) hi = new G4LEAntiXiZeroInelastic();
347  else if(pname == "kaon-"       ) hi = new G4LEKaonMinusInelastic();
348  else if(pname == "kaon+"       ) hi = new G4LEKaonPlusInelastic();
349  else if(pname == "kaon0S"      ) hi = new G4LEKaonZeroSInelastic();
350  else if(pname == "kaon0L"      ) hi = new G4LEKaonZeroLInelastic();
351  else if(pname == "neutron"     ) hi = new G4LENeutronInelastic();
352  else if(pname == "lambda"      ) hi = new G4LELambdaInelastic();
353  else if(pname == "omega-"      ) hi = new G4LEOmegaMinusInelastic();
354  else if(pname == "proton"      ) hi = new G4LEProtonInelastic();
355  else if(pname == "pi+"         ) hi = new G4LEPionPlusInelastic();
356  else if(pname == "pi-"         ) hi = new G4LEPionMinusInelastic();
357  else if(pname == "sigma-"      ) hi = new G4LESigmaMinusInelastic();
358  else if(pname == "sigma+"      ) hi = new G4LESigmaPlusInelastic();
359  else if(pname == "xi-"         ) hi = new G4LEXiMinusInelastic();
360  else if(pname == "xi0"         ) hi = new G4LEXiZeroInelastic();
361
362  if(hi) {
363    hi->SetMinEnergy(emin);
364    hi->SetMaxEnergy(emax);
365    Register(particle,hp,hi,"LHEP");
366  } else {
367    G4cout << "### G4HadronInelasticTHEO: ERROR - no LHEP model for "
368           << pname << G4endl;
369  }
370}
371
372
373void G4HadronInelasticQLHEP::AddHEP(G4ParticleDefinition* particle,
374                                    G4HadronicProcess* hp,
375                                    G4double emin,
376                                    G4double emax)
377{
378  G4HadronicInteraction* hi = 0;
379  G4String pname = particle->GetParticleName();
380
381  if(pname == "anti_lambda"      ) hi = new G4HEAntiLambdaInelastic();
382  else if(pname == "anti_neutron") hi = new G4HEAntiNeutronInelastic();
383  else if(pname == "anti_omega-" ) hi = new G4HEAntiOmegaMinusInelastic();
384  else if(pname == "anti_proton" ) hi = new G4HEAntiProtonInelastic();
385  else if(pname == "anti_sigma-" ) hi = new G4HEAntiSigmaMinusInelastic();
386  else if(pname == "anti_sigma+" ) hi = new G4HEAntiSigmaPlusInelastic();
387  else if(pname == "anti_xi-"    ) hi = new G4HEAntiXiMinusInelastic();
388  else if(pname == "anti_xi0"    ) hi = new G4HEAntiXiZeroInelastic();
389  else if(pname == "kaon-"       ) hi = new G4HEKaonMinusInelastic();
390  else if(pname == "kaon+"       ) hi = new G4HEKaonPlusInelastic();
391  else if(pname == "kaon0S"      ) hi = new G4HEKaonZeroInelastic();
392  else if(pname == "kaon0L"      ) hi = new G4HEKaonZeroInelastic();
393  else if(pname == "neutron"     ) hi = new G4HENeutronInelastic();
394  else if(pname == "lambda"      ) hi = new G4HELambdaInelastic();
395  else if(pname == "omega-"      ) hi = new G4HEOmegaMinusInelastic();
396  else if(pname == "proton"      ) hi = new G4HEProtonInelastic();
397  else if(pname == "pi+"         ) hi = new G4HEPionPlusInelastic();
398  else if(pname == "pi-"         ) hi = new G4HEPionMinusInelastic();
399  else if(pname == "sigma-"      ) hi = new G4HESigmaMinusInelastic();
400  else if(pname == "sigma+"      ) hi = new G4HESigmaPlusInelastic();
401  else if(pname == "xi-"         ) hi = new G4HEXiMinusInelastic();
402  else if(pname == "xi0"         ) hi = new G4HEXiZeroInelastic();
403
404  if(hi) {
405    hi->SetMinEnergy(emin);
406    hi->SetMaxEnergy(emax);
407    Register(particle,hp,hi,"HEP");
408  } else {
409    G4cout << "### G4HadronInelasticQLHEP: ERROR - no HEP model for "
410           << pname << G4endl;
411  }
412}
413
414void G4HadronInelasticQLHEP::Register(G4ParticleDefinition* p, 
415                                      G4HadronicProcess* hp, 
416                                      G4HadronicInteraction* hi, 
417                                      const G4String& m)
418{
419  hp->RegisterMe(hi);
420  if(verbose > 1)
421    G4cout << "### QLHEP: Register new model " << m
422           << " for " << p->GetParticleName() << " and " 
423           << hp->GetProcessName()
424           << " E(GeV) " << hi->GetMinEnergy()/GeV
425           << " - " << hi->GetMaxEnergy()/GeV << G4endl;
426}
Note: See TracBrowser for help on using the repository browser.