Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DPMJET2_5Model.cc
Go to the documentation of this file.
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 // * *
21 // * Parts of this code which have been developed by QinetiQ Ltd *
22 // * under contract to the European Space Agency (ESA) are the *
23 // * intellectual property of ESA. Rights to use, copy, modify and *
24 // * redistribute this software for general public use are granted *
25 // * in compliance with any licensing, distribution and development *
26 // * policy adopted by the Geant4 Collaboration. This code has been *
27 // * written by QinetiQ Ltd for the European Space Agency, under ESA *
28 // * contract 19770/06/NL/JD (Technology Research Programme). *
29 // * *
30 // * By using, copying, modifying or distributing the software (or *
31 // * any work based on the software) you agree to acknowledge its *
32 // * use in resulting scientific publications, and indicate your *
33 // * acceptance of all terms of the Geant4 Software license. *
34 // ********************************************************************
35 //
36 /// \file hadronic/Hadr02/src/G4DPMJET2_5Model.cc
37 /// \brief Implementation of the G4DPMJET2_5Model class
38 //
39 // $Id: G4DPMJET2_5Model.cc 77519 2013-11-25 10:54:57Z gcosmo $
40 //
41 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
42 //
43 // MODULE: G4DPMJET2_5Model.cc
44 //
45 // Version: 0.B
46 // Date: 02/04/08
47 // Author: P R Truscott
48 // Organisation: QinetiQ Ltd, UK
49 // Customer: ESA/ESTEC, NOORDWIJK
50 // Contract: 19770/06/NL/JD
51 //
52 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
53 ///////////////////////////////////////////////////////////////////////////////
54 //
55 #ifdef G4_USE_DPMJET
56 
57 
58 #include "G4DPMJET2_5Model.hh"
59 #include "G4GlaubAADataSetHandler.hh"
60 
61 #include "G4ExcitationHandler.hh"
62 #include "G4Evaporation.hh"
63 #include "G4FermiBreakUp.hh"
64 #include "G4PhotonEvaporation.hh"
65 #include "G4PreCompoundModel.hh"
66 #include "G4ParticleDefinition.hh"
67 #include "G4ParticleTable.hh"
68 #include "G4DynamicParticle.hh"
69 #include "Randomize.hh"
70 #include "G4Fragment.hh"
71 #include "G4VNuclearDensity.hh"
72 #include "G4NuclearShellModelDensity.hh"
73 #include "G4NuclearFermiDensity.hh"
74 #include "G4FermiMomentum.hh"
75 #include "G4ReactionProductVector.hh"
76 #include "G4LorentzVector.hh"
77 #include "G4ParticleMomentum.hh"
78 #include "G4Poisson.hh"
79 #include "G4ParticleTable.hh"
80 #include "G4IonTable.hh"
81 #include "G4LorentzVector.hh"
82 #include "G4HadTmpUtil.hh"
83 #include "globals.hh"
84 #include "G4PhysicalConstants.hh"
85 #include "G4SystemOfUnits.hh"
86 
87 #include <fstream>
88 
89 #include "G4DPMJET2_5Interface.hh"
90 
91 /////////////////////////////////////////////////////////////////////////////////
92 //
93 // Constructor without arguments
94 //
95 // This constructor uses a default pre-compound. It initialises the
96 // variables (including the de-excitation), but note that much of the work is
97 // done in the member function Initialise(), which is dedicated to
98 // initialising variables in DPMJET-II.5 to class-default values.
99 //
100 
101 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
102 
103 G4DPMJET2_5Model::G4DPMJET2_5Model () : G4HadronicInteraction("DPMJET2_5")
104 {
105 //
106 // Set the default verbose level to 0 - no output.
107 //
108  SetVerboseLevel(1);
109 //
110 #ifdef G4VERBOSE
111  if (GetVerboseLevel()>0) {
112  G4cout <<"G4DPMJET2_5Model default constructor" <<G4endl;
113  }
114 #endif
115 //
116 //
117 // Send message to stdout to advise that the G4DPMJET2_5Model model
118 // is being used.
119 //
120  theInitType = DEFAULT;
121  PrintWelcomeMessage();
122 //
123 // Use the precompound model for nuclear de-excitation.
124 //
125  theExcitationHandler = 0;
126  SetDefaultPreCompoundModel();
127 //
128 //
129 // Set the minimum and maximum range for the model (despite nomanclature, this
130 // is in energy per nucleon number).
131 //
132  SetMinEnergy(5.0*GeV);
133  SetMaxEnergy(1000.0*TeV);
134 
135  SetEnergyMomentumCheckLevels(100*perCent, 1*TeV);
136 
137 //
138 //
139 // Initialise the DPMJET model - this effectively does what the DPMJET
140 // subroutine DMINIT does without reading an input file.
141 //
142  debug = false;
143  debug_level = 0;
144  lunber = 14;
145  dpmver = 2.5;
146 
147  LFALSE = 0;
148  LTRUE = 1;
149 
150  Initialise ();
151 //
152 //
153 // Next bit directs how and how many Glauber data sets are loaded
154 // or created.
155 //
156  theGlauberDataSetHandler = G4GlaubAADataSetHandler::getInstance();
157  theGlauberDataSetHandler->SetMaxGlauberDataSets(-1);
158 
159  theParticleTable = G4ParticleTable::GetParticleTable();
160  theIonTable = const_cast <G4IonTable *> (theParticleTable->GetIonTable());
161 //
162 //
163 }
164 
165 const std::pair<G4double, G4double>
166 G4DPMJET2_5Model::GetFatalEnergyCheckLevels() const
167 {
168  // default level of Check
169  return std::pair<G4double, G4double>(100.*perCent, 500 * GeV);
170 }
171 
172 /////////////////////////////////////////////////////////////////////////////////
173 //
174 // Constructor with DPMJET-II.5 initialisation type
175 //
176 // This constructor uses a default pre-compound. It initialises the
177 // variables (including the de-excitation), but note that much of the work is
178 // done in the member function Initialise(), which is dedicated to
179 // initialising variables in DPMJET-II.5. The user is able to define whether to
180 // use the default values, DPMJET-II.5 settings or DPMJET-III settings.
181 //
182 G4DPMJET2_5Model::G4DPMJET2_5Model (const G4DPMJET2_5InitialisationType initType)
183 {
184 //
185 // Set the default verbose level to 0 - no output.
186 //
187  SetVerboseLevel(1);
188 //
189 #ifdef G4VERBOSE
190  if (GetVerboseLevel()>0) {
191  G4cout <<"G4DPMJET2_5Model constructor #1 " <<G4endl;
192  }
193 #endif
194 //
195 //
196 // Send message to stdout to advise that the G4DPMJET2_5Model model
197 // is being used.
198 //
199  theInitType = initType;
200  PrintWelcomeMessage();
201 //
202 // Use the precompound model for nuclear de-excitation.
203 //
204  theExcitationHandler = 0;
205  SetDefaultPreCompoundModel();
206 //
207 //
208 // Set the minimum and maximum range for the model (despite nomanclature, this
209 // is in energy per nucleon number).
210 //
211  SetMinEnergy(5.0*GeV);
212  SetMaxEnergy(1000.0*TeV);
213 //
214 //
215 // Initialise the DPMJET model - this effectively does what the DPMJET
216 // subroutine DMINIT does without reading an input file.
217 //
218  debug = false;
219  debug_level = 0;
220  lunber = 14;
221  dpmver = 2.5;
222 
223  LFALSE = 0;
224  LTRUE = 1;
225 
226  Initialise ();
227 //
228 //
229 // Next bit directs how and how many Glauber data sets are loaded
230 // or created.
231 //
232  theGlauberDataSetHandler = G4GlaubAADataSetHandler::getInstance();
233  theGlauberDataSetHandler->SetMaxGlauberDataSets (-1);
234 
235  theParticleTable = G4ParticleTable::GetParticleTable();
236  theIonTable = theParticleTable->GetIonTable();
237 //
238 //
239 }
240 ////////////////////////////////////////////////////////////////////////////////
241 //
242 // Constructor with de-excitation handler.
243 //
244 // This constructor uses the user-provided de-excitation handler. However, it
245 // is possible for the use to provide a NULL pointer, in which case, it is
246 // assumed that the user doesn't want to simulate de-excitation - USER, BEWARE!
247 //
248 // The member function initialises the variables (including the de-excitation),
249 // but note that much of the work is done in the member function Initialise(),
250 // which is dedicated to initialising variables in DPMJET-II.5.
251 //
252 G4DPMJET2_5Model::G4DPMJET2_5Model (G4ExcitationHandler *aExcitationHandler,
253  const G4DPMJET2_5InitialisationType initType)
254 {
255 //
256 // Set the default verbose level to 0 - no output.
257 //
258  SetVerboseLevel(1);
259 //
260 // Send message to stdout to advise that the G4DPMJET2_5Model model
261 // is being used.
262 //
263 #ifdef G4VERBOSE
264  if (GetVerboseLevel()>0) {
265  G4cout <<"G4DPMJET2_5Model constructor #2 " <<G4endl;
266  }
267 #endif
268  theInitType = initType;
269  PrintWelcomeMessage();
270 //
271 // The user is able to provide the excitation handler.
272 //
273  theExcitationHandler = aExcitationHandler;
274  thePreComp = 0;
275 //
276 //
277 // Set the minimum and maximum range for the model (despite nomanclature, this
278 // is in energy per nucleon number).
279 //
280  SetMinEnergy(5.0*GeV);
281  SetMaxEnergy(1000.0*TeV);
282 //
283 //
284 // Initialise the DPMJET model - this effectively does what the DPMJET
285 // subroutine DMINIT does without reading an input file.
286 //
287  debug = false;
288  debug_level = 0;
289  lunber = 14;
290  dpmver = 2.5;
291 
292  LFALSE = 0;
293  LTRUE = 1;
294 
295  Initialise ();
296 //
297 //
298 // Next bit directs how and how many Glauber data sets are loaded
299 // or created.
300 //
301  theGlauberDataSetHandler = G4GlaubAADataSetHandler::getInstance();
302  theGlauberDataSetHandler->SetMaxGlauberDataSets (-1);
303 
304  theParticleTable = G4ParticleTable::GetParticleTable();
305  theIonTable = theParticleTable->GetIonTable();
306 //
307 //
308 }
309 ////////////////////////////////////////////////////////////////////////////////
310 //
311 // Constructor with pre-compound model.
312 //
313 // This constructor uses the user-provided pre-equilibrium model. However, it
314 // is possible for the use to provide a NULL pointer, in which case, it is
315 // assumed that the user doesn't want to simulate pre-equilibrium. - USER,
316 // BEWARE!
317 //
318 // The member function initialises the variables (including the de-excitation),
319 // but note that much of the work is done in the member function Initialise(),
320 // which is dedicated to initialising variables in DPMJET-II.5.
321 //
322 G4DPMJET2_5Model::G4DPMJET2_5Model (G4VPreCompoundModel *aPreComp,
323  const G4DPMJET2_5InitialisationType initType)
324 {
325 //
326 // Set the default verbose level to 0 - no output.
327 //
328  SetVerboseLevel(1);
329 //
330 // Send message to stdout to advise that the G4DPMJET2_5Model model
331 // is being used.
332 //
333 #ifdef G4VERBOSE
334  if (GetVerboseLevel()>0) {
335  G4cout <<"G4DPMJET2_5Model constructor #3 " <<G4endl;
336  }
337 #endif
338  theInitType = initType;
339  PrintWelcomeMessage();
340 //
341 // The user is able to provide the pre-compound model.
342 //
343  theExcitationHandler = 0;
344  thePreComp = aPreComp;
345 //
346 //
347 // Set the minimum and maximum range for the model (despite nomanclature, this
348 // is in energy per nucleon number).
349 //
350  SetMinEnergy(5.0*GeV);
351  SetMaxEnergy(1000.0*TeV);
352 //
353 //
354 // Initialise the DPMJET model - this effectively does what the DPMJET
355 // subroutine DMINIT does without reading an input file.
356 //
357  debug = false;
358  debug_level = 0;
359  lunber = 14;
360  dpmver = 2.5;
361 
362  LFALSE = 0;
363  LTRUE = 1;
364 
365  Initialise ();
366 //
367 //
368 // Next bit directs how and how many Glauber data sets are loaded
369 // or created.
370 //
371  theGlauberDataSetHandler = G4GlaubAADataSetHandler::getInstance();
372  theGlauberDataSetHandler->SetMaxGlauberDataSets (-1);
373 
374  theParticleTable = G4ParticleTable::GetParticleTable();
375  theIonTable = theParticleTable->GetIonTable();
376 //
377 //
378 }
379 ////////////////////////////////////////////////////////////////////////////////
380 //
381 // Destructor
382 //
383 G4DPMJET2_5Model::~G4DPMJET2_5Model ()
384 {
385 //
386 //
387 // The destructor doesn't have to do a great deal!
388 //
389  if (theExcitationHandler) delete theExcitationHandler;
390  if (thePreComp) delete thePreComp;
391 // delete theGlauberDataSetHandler;
392  theGlauberDataSetHandler->UnloadAllGlauberData();
393 }
394 ////////////////////////////////////////////////////////////////////////////////
395 //
396 // IsApplicable
397 //
398 // This member function simply determines whether there is relevant information
399 // in Glauber data for this projectile and target, and if the nucleus is
400 // sensible.
401 //
402 G4bool G4DPMJET2_5Model::IsApplicable (
403  const G4HadProjectile &theTrack, G4Nucleus &theTarget)
404 {
405 //
406 //
407 // Get relevant information about the projectile and target (A, Z)
408 //
409  const G4ParticleDefinition *definitionP = theTrack.GetDefinition();
410  G4int AP = definitionP->GetBaryonNumber();
411  G4int ZP = G4int(definitionP->GetPDGCharge()/eplus + 0.5);
412  G4int AT = theTarget.GetA_asInt();
413  G4int ZT = theTarget.GetZ_asInt();
414 
415  if (AP >= 2 && ZP >= 1 && AT >= 2 && ZT >=1) {
416  return theGlauberDataSetHandler->IsGlauberDataSetAvailable(AP,AT);
417  }
418  else {
419  return false;
420  }
421 }
422 ////////////////////////////////////////////////////////////////////////////////
423 //
424 // ApplyYourself
425 //
426 // Member function to process an event, and get information about the products.
427 // The phases are:
428 //
429 // (1) Determine the information about the projectile and target.
430 //
431 // (2) Identify to the Glauber data set handler which data need to be used. If
432 // the GDSH finds there are no Glauber profile data for the collision,the
433 // product is the unchanged projectile.
434 //
435 // (3) Initialise further common-block variables in DPMJET-II.5, and perform
436 // FORTRAN calls (note this is taken largely from an interpretation of CORSIKA
437 // and DPMJET-II.5 FORTRAN ... I think there's duplication in some of the
438 // initialisation on top of what Initialise() does, but CORSIKA has similar
439 // duplication. I'm not confident at removing this out at the moment.
440 //
441 // (4) Call the DPMJET-II.5 FORTRAN subroutine DPMEVT.
442 //
443 // (5) Pick out the final state particles and nuclei. In the case of nuclei
444 // use the de-excitation handler if one has been defined. Transfer all these
445 // particles to the final-state vector.
446 //
447 // (6) if very-verbose output is demanded by the user, there is a print-out
448 // of the total energy and total momentum before and after collision.
449 //
450 G4HadFinalState *G4DPMJET2_5Model::ApplyYourself (
451  const G4HadProjectile &theTrack, G4Nucleus &theTarget)
452 {
453 //
454 //
455 // The secondaries will be returned in G4HadFinalState &theParticleChange -
456 // initialise this. The original track will always be discontinued and
457 // secondaries followed.
458 //
459  theParticleChange.Clear();
460  theParticleChange.SetStatusChange(stopAndKill);
461 //
462 //
463 // Get relevant information about the projectile and target (A, Z, energy/nuc,
464 // momentum, etc).
465 //
466  const G4ParticleDefinition *definitionP = theTrack.GetDefinition();
467  G4int AP = definitionP->GetBaryonNumber();
468  G4int ZP = G4int(definitionP->GetPDGCharge()/eplus+0.5);
469  G4double M = definitionP->GetPDGMass();
470  G4ThreeVector pP = theTrack.Get4Momentum().vect();
471  G4double T = theTrack.GetKineticEnergy()/G4double(AP);
472  // Units are MeV/nuc
473  G4double E = theTrack.GetTotalEnergy()/G4double(AP);
474  // Units are MeV/nuc
475  G4int AT = theTarget.GetA_asInt();
476  G4int ZT = theTarget.GetZ_asInt();
477  G4double mpnt = theTarget.AtomicMass(AT, ZT);
478  G4double TotalEPre = theTrack.GetTotalEnergy() + mpnt;
479  // theTarget.AtomicMass(AT, ZT) + theTarget.GetEnergyDeposit();
480 // G4LorentzRotation transformToLab =
481 // (const_cast <G4HadProjectile*> (&theTrack))->GetTrafoToLab();
482 //
483 //
484 // Output relevant information on initial conditions if verbose. Note that
485 // most of the verbsse output is dealt with through private member function
486 // calls.
487 //
488  if (verboseLevel >= 2)
489  {
490  G4cout <<"########################################"
491  <<"########################################"
492  <<G4endl;
493  G4cout.precision(6);
494  G4cout <<"IN G4DPMJET2_5Model::ApplyYourself" <<G4endl;
495  G4cout <<"START OF EVENT" <<G4endl;
496  G4cout <<"Initial projectile (A,Z) = (" <<AP <<", " <<ZP <<")" <<G4endl;
497  G4cout <<"Initial target (A,Z) = (" <<AT <<", " <<ZT <<")" <<G4endl;
498  G4cout <<"Projectile momentum = " <<pP/MeV <<" MeV/c" <<G4endl;
499  G4cout <<"Total energy = " <<TotalEPre/MeV <<" MeV" <<G4endl;
500  G4cout <<"Kinetic energy/nuc = " <<T/MeV <<" MeV" <<G4endl;
501  }
502 //
503 //
504 // Setup variables and call the DPMJET model. There is a significant amount of
505 // initialisation which still needs to be done, based on DPMJET-II.5 main
506 // program and card reader subroutine, and the CORSIKA implementation.
507 //
508  G4int AP1 = AP;
509  G4int ZP1 = ZP;
510  G4int AT1 = AT;
511  G4int ZT1 = ZT;
512 //
513 //
514 // vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
515 // ***** WARNING *****
516 // The following is a provisional "catch" for ions with A==Z. The de-excitation
517 // and precompound model can produce such nuclei, although they should decay
518 // into free protons. For the moment, DPMJET-II.5 doesn't treat them ... i.e.
519 // the FORTRAN code would crash. Therefore, return such ions without
520 // nuclear interactions.
521 //
522  if (AP1 > 1 && AP1 == ZP1) {
523  theParticleChange.SetStatusChange(isAlive);
524  theParticleChange.SetEnergyChange(theTrack.GetKineticEnergy());
525  theParticleChange.SetMomentumChange(theTrack.Get4Momentum().vect().unit());
526  if (verboseLevel >= 2) {
527  G4cout <<"PROJECTILE WITH AP = " <<AP1 <<" == ZP = " <<ZP1
528  <<" REJECTED" <<G4endl;
529  G4cout <<"########################################"
530  <<"########################################"
531  <<G4endl;
532  }
533  return &theParticleChange;
534  }
535 // ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
536 
537  nucc_.ibproj = 1; // IBPROJ = 1
538  nucc_.ijproj = 1; // IJPROJ = 1
539  collis_.ijprox = 1; // IJPROX = 1
540  nucc_.ip = AP1; // IP = IP_P
541  nucc_.ipz = ZP1; // IPZ = IPZ_P
542  nucc_.it = AT1; // IT = IT_P
543  nucc_.itz = ZT1; // ITZ = ITZ_P
544  collis_.ijtar = 1; // IJTAR=1
545 //
546 //
547 // Note that epn is deliberately given units of GeV/nuc, because of the units
548 // used in DPMJET-II.5.
549 //
550  G4double epn = E / GeV;
551  G4double mpn = M / (GeV*AP); // Projectile mass per nucleon in GeV/(c2*nuc)
552 // G4double gamma = epn / mpn;
553 // G4double elab = epn * M / GeV; // ELAB = EPN * AMPRO_P
554 
555  diffra_.isingd = ISINGD;
556  user2_.isingx = ISINGX;
557  user2_.idubld = IDUBLD;
558  user2_.sdfrac = SDFRAC;
559 
560 // G4double amu_c2GeV = amu_c2 / GeV;
561  G4double ppn = std::sqrt((epn-mpn)*(epn+mpn));
562 // Units of GeV/(c*nuc)
563 // PPN = SQRT( (EPN-AMPROJ)*(EPN+AMPROJ) )
564 //
565 //
566 // Before setting the remainder of the variables for DPMJET-II.5, check for
567 // appropriate Glauber data. If the value returned is false, then set the
568 // change object to the source particle.
569 //
570  if (!(theGlauberDataSetHandler->SetCurrentGlauberDataSet(AP1,AT1,ppn))) {
571  theParticleChange.SetStatusChange(isAlive);
572  theParticleChange.SetEnergyChange(theTrack.GetKineticEnergy());
573  theParticleChange.SetMomentumChange(theTrack.Get4Momentum().vect().unit());
574  return &theParticleChange;
575  }
576  dtumat_.ntaxx[0] = AT1;
577  dtumat_.nztaxx[0] = ZT1;
578  dtumat_.nprxx[0] = AP1;
579  dtumat_.nzprxx[0] = ZP1;
580 //
581 //
582 // Set the remainder of the variables for DPMJET-II.5 FORTRAN.
583 //
584  nncms_.pproj = ppn; // PPROJ = PPN
585  nncms_.eproj = epn; // EPROJ = EPN
586  mpnt /= (AT * GeV); // Mass per nuclon of target in GeV/(c2*nuc)
587  nncms_.umo = std::sqrt(mpn*mpn + mpnt*mpnt + 2.0*mpnt*epn);
588 // UMO = SQRT( AMPROJ**2 + AMTAR**2 +2.D0*AMTAR*EPROJ )
589 // Note I believe this equation is only correct
590 // if the subsequent equations (for pTthr)
591 // needs the Ecm for the NUCLEON-NUCLEON system
592  user2_.cmener = nncms_.umo; // CMENER = UMO
593  collis_.s = nncms_.umo * nncms_.umo;
594  // SS = UMO**2
595  collis_.ptthr = 3.0;
596  if (strufu_.istrut == 1)
597  {
598  collis_.ptthr = 2.1 + 0.15*std::pow(std::log10(user2_.cmener/50.),3.0);
599  // PTTHR = 2.1D0+0.15D0*(LOG10(CMENER/50.))**3
600  }
601  else if (strufu_.istrut == 2)
602  {
603  collis_.ptthr = 2.5 + 0.12*std::pow(std::log10(user2_.cmener/50.),3.0);
604  // PTTHR = 2.5D0+0.12D0*(LOG10(CMENER/50.))**3
605  }
606  collis_.ptthr2 = collis_.ptthr; // PTTHR2 = PTTHR
607  nncms_.gamcm = (epn + mpnt) / nncms_.umo;
608  // GAMCM = (EPROJ+AMTAR)/UMO
609  // Note I believe this equation is only correct
610  // if the subsequent equations (for pTthr)
611  // need the Ecm for the NUCLEON-NUCLEON system
612  nncms_.bgcm = ppn / nncms_.umo; // PPROJ/UMO
613  nncms_.pcm = nncms_.gamcm*ppn - nncms_.bgcm*epn;
614  // PCM = GAMCM*PPROJ - BGCM*EPROJ
615  sigma_.sigsof = 37.8 * std::pow(collis_.s,0.076);
616  // ALFA = 1.076D0
617  // A = 37.8D0
618  // SIGSOF = A * SS**(ALFA-1.D0)
619  seasu3_.seasq = SEASQ; // SEASQ = 0.50D0
620  xseadi_.ssmima = SSMIMA; // SSMIMA = 1.201D0
621  xseadi_.ssmimq = xseadi_.ssmima * xseadi_.ssmima;
622  // SSMIMQ = SSMIMA**2
623  taufo_.taufor = TAUFOR;
624  taufo_.ktauge = KTAUGE;
625 
626  if ( theInitType == DEFAULT ) {
627  final_.ifinal = 0; // IFINAL = 1
628  evappp_.ievap = 0; // IEVAP = 0
629  parevt_.levprt = LTRUE; // LEVPRT = .FALSE.
630  parevt_.ilvmod = 1; // ILVMOD = 1
631  parevt_.ldeexg = LFALSE; // LDEEXG = .FALSE.
632  parevt_.lheavy = LFALSE; // LHEAVY = .FALSE.
633  frbkcm_.lfrmbk = LFALSE; // LFRMBK = .FALSE.
634  inpflg_.ifiss = 0; // IFISS = 0
635  } else if ( theInitType == CORSIKA ) {
636  final_.ifinal = 0; // IFINAL = 1
637  evappp_.ievap = 0; // IEVAP = 0
638  parevt_.levprt = LTRUE; // LEVPRT = .FALSE.
639  parevt_.ilvmod = 1; // ILVMOD = 1
640  parevt_.ldeexg = LFALSE; // LDEEXG = .FALSE.
641  parevt_.lheavy = LTRUE; // LHEAVY = .FALSE.
642  frbkcm_.lfrmbk = LFALSE; // LFRMBK = .FALSE.
643  inpflg_.ifiss = 0; // IFISS = 0
644  }
645  else if ( theInitType == DPM2_5 ) {
646  final_.ifinal = 0; // IFINAL = 0
647  evappp_.ievap = 0; // IEVAP = 0
648  parevt_.levprt = LTRUE;
649 // LEVPRT = .TRUE. NOTE: THIS IS AT ODDS WITH WHAT'S IN DPMJET-II.5,
650 // BUT IF NOT SET, ALL EVENTS GET REJECTED.
651  parevt_.ilvmod = 1; // ILVMOD = 1
652  parevt_.ldeexg = LFALSE; // LDEEXG = .FALSE.
653  parevt_.lheavy = LFALSE; // LHEAVY = .FALSE.
654  frbkcm_.lfrmbk = LFALSE; // LFRMBK = .FALSE.
655  inpflg_.ifiss = 0; // IFISS = 0
656  }
657  else if ( theInitType == DPM3 ) {
658  final_.ifinal = 0; // IFINAL = 0
659  evappp_.ievap = 0; // IEVAP = 0
660  parevt_.levprt = LTRUE;
661 // LEVPRT = .TRUE. NOTE: THIS IS AT ODDS WITH WHAT'S IN DPMJET-II.5,
662 // BUT IF NOT SET, ALL EVENTS GET REJECTED.
663  parevt_.ilvmod = 1; // ILVMOD = 1
664  parevt_.ldeexg = LFALSE; // LDEEXG = .FALSE.
665  parevt_.lheavy = LFALSE; // LHEAVY = .FALSE.
666  frbkcm_.lfrmbk = LFALSE; // LFRMBK = .FALSE.
667  inpflg_.ifiss = 0; // IFISS = 0
668  }
669  xsecpt_.ptcut = collis_.ptthr; // PTCUT = PTTHR
670 
671  G4double dsig1[maxpro+1];
672  csj1mi_ (&xsecpt_.ptcut, &dsig1[0]); // CALL CSJ1MI(PTCUT,DSIG1)
673  xsecpt_.dsigh = dsig1[0]; // SIG1 = DSIG1(0)
674  // DSIGH = SIG1
675  G4int i = 0;
676  G4double pt = 0.0;
677  samppt_ (&i,&pt); // SAMPPT(0,PT)
678  collap_.s3 = collis_.s; // S3 = SS
679  collap_.ijproj1 = collis_.ijprox; // IJPROJ1 = IJPROX
680  collap_.ijtar1 = collis_.ijtar; // IJTAR1 = IJTAR
681  collap_.ptthr1 = collis_.ptthr; // PTTHR1 = PTTHR
682  collap_.iophrd1 = collis_.iophrd; // IOPHRD1 = IOPHRD
683  collap_.ijprlu1 = collis_.ijprlu; // IJPRLU1 = IJPRLU
684  collap_.ijtalu1 = collis_.ijtalu; // IJTALU1 = IJTALU
685  collap_.ptthr3 = collis_.ptthr2; // PTTHR3 = PTTHR2
686 
687  G4int iiipro = nucc_.ijproj;
688  // IIPROJ = IJPROJ & IIIPRO = IIPROJ
689  G4int iiitar = nucc_.ijtarg; // IITARG = IJTARG
690  G4int kkmat = 1;
691  G4int nhkkh1 = 1; // NHKKH1 = 1
692 
693  for (i=0; i<8; i++) user1_.projty[i] = paname_.btype[iiipro][i];
694  // PROJTY=BTYPE(IPROJ)
695  G4int irej = 1;
696  G4int evtcnt = 0;
697 
698  do {
699 /* bufueh_.annvv = 0.001; // ANNVV = 0.001
700  bufueh_.annss = 0.001; // ANNSS = 0.001
701  bufueh_.annsv = 0.001; // ANNSV = 0.001
702  bufueh_.annvs = 0.001; // ANNVS = 0.001
703  bufueh_.anncc = 0.001; // ANNCC = 0.001
704  bufueh_.anndv = 0.001; // ANNDV = 0.001
705  bufueh_.annvd = 0.001; // ANNVD = 0.001
706  bufueh_.annds = 0.001; // ANNDS = 0.001
707  bufueh_.annsd = 0.001; // ANNSD = 0.001
708  bufueh_.annhh = 0.001; // ANNHH = 0.001
709  bufueh_.annzz = 0.001; // ANNZZ = 0.001
710  bufueh_.anndi = 0.001; // ANNDI = 0.001
711  bufueh_.annzd = 0.001; // ANNZD = 0.001
712  bufueh_.anndz = 0.001; // ANNDZ = 0.001
713  bufueh_.ptvv = 0.0; // PTVV = 0.
714  bufueh_.ptss = 0.0; // PTSS = 0.
715  bufueh_.ptsv = 0.0; // PTSV = 0.
716  bufueh_.ptvs = 0.0; // PTVS = 0.
717  bufueh_.ptcc = 0.0; // PTCC = 0.
718  bufueh_.ptdv = 0.0; // PTDV = 0.
719  bufueh_.ptvd = 0.0; // PTVD = 0.
720  bufueh_.ptds = 0.0; // PTDS = 0.
721  bufueh_.ptsd = 0.0; // PTSD = 0.
722  bufueh_.pthh = 0.0; // PTHH = 0.
723  bufueh_.ptzz = 0.0; // PTZZ = 0.
724  bufueh_.ptdi = 0.0; // PTDI = 0.
725  bufueh_.ptzd = 0.0; // PTZD = 0.
726  bufueh_.ptdz = 0.0; // PTDZ = 0.
727  bufueh_.eevv = 0.0; // EEVV = 0.
728  bufueh_.eess = 0.0; // EESS = 0.
729  bufueh_.eesv = 0.0; // EESV = 0.
730  bufueh_.eevs = 0.0; // EEVS = 0.
731  bufueh_.eecc = 0.0; // EECC = 0.
732  bufueh_.eedv = 0.0; // EEDV = 0.
733  bufueh_.eevd = 0.0; // EEVD = 0.
734  bufueh_.eeds = 0.0; // EEDS = 0.
735  bufueh_.eesd = 0.0; // EESD = 0.
736  bufueh_.eehh = 0.0; // EEHH = 0.
737  bufueh_.eezz = 0.0; // EEZZ = 0.
738  bufueh_.eedi = 0.0; // EEDI = 0.
739  bufueh_.eezd = 0.0; // EEZD = 0.
740  bufueh_.eedz = 0.0; // EEDZ = 0.
741  ncouch_.acouvv = 0.0; // ACOUVV = 0.
742  ncouch_.acouss = 0.0; // ACOUSS = 0.
743  ncouch_.acousv = 0.0; // ACOUSV = 0.
744  ncouch_.acouvs = 0.0; // ACOUVS = 0.
745  ncouch_.acouzz = 0.0; // ACOUZZ = 0.
746  ncouch_.acouhh = 0.0; // ACOUHH = 0.
747  ncouch_.acouds = 0.0; // ACOUDS = 0.
748  ncouch_.acousd = 0.0; // ACOUSD = 0.
749  ncouch_.acoudz = 0.0; // ACOUDZ = 0.
750  ncouch_.acouzd = 0.0; // ACOUZD = 0.
751  ncouch_.acoudi = 0.0; // ACOUDI = 0.
752  ncouch_.acoudv = 0.0; // ACOUDV = 0.
753  ncouch_.acouvd = 0.0; // ACOUVD = 0.
754  ncouch_.acoucc = 0.0; // ACOUCC = 0.*/
755  if (evtcnt > 0)
756  G4cout <<"REJECTED KKINC EVENT. RETRY # = " <<evtcnt <<G4endl;
757 //
758 //
759 // Generate an event using DPMJET II-5. NOTE this is a call to dpmevt, and
760 // could possibly be replaced by call to kkinc_. After the call,
761 // ResetCurrentGlauberDataSet is used to delete any temporary Glauber data
762 // set if one had to be set up.
763 //
764  //G4cout << "Call to kkinc_" << G4endl;
765  kkinc_ (&epn, &AT1, &ZT1, &AP1, &ZP1, &iiipro, &kkmat, &iiitar, &nhkkh1,
766  &irej);
767  // CALL KKINC(EPN,IIT,IITZ,IIP,IIPZ,IIPROJ,KKMAT,
768  // * IITARG,NHKKH1,IREJ)
769 // dpmevt_ (&elabt, &iiipro, &AP1, &ZP1, &AT1, &ZT1, &kkmat, &nhkkh1);
770  } while (irej == 1 && ++evtcnt <100);
771  //G4cout << "Call to reset G-data" << G4endl;
772 
773  theGlauberDataSetHandler->ResetCurrentGlauberDataSet();
774 //
775 //
776 // If the event has been rejected more than 100 times, then set the track
777 // as still active and return to the calling routine.
778 //
779  if (irej == 1) {
780  theParticleChange.SetStatusChange(isAlive);
781  theParticleChange.SetEnergyChange(theTrack.GetKineticEnergy());
782  theParticleChange.SetMomentumChange(theTrack.Get4Momentum().vect().unit());
783  if (verboseLevel >= 2) {
784  G4cout <<"Event rejected and original track maintained" <<G4endl;
785  G4cout <<"########################################"
786  <<"########################################"
787  <<G4endl;
788  }
789  return &theParticleChange;
790  }
791 //
792 //
793 // Determine number of final state particles (including nuclear fragments,
794 // and load into the particle change if you can identify the particles.
795 //
796  G4int n = hkkevt_.nhkk;
797  G4int M1 = 0;
798  G4Fragment *fragment = 0;
799  if (verboseLevel >= 2) DumpVerboseInformation1 (n);
800 //
801 //
802 // Now go through each of the secondaries and add to theParticleChange.
803 //
804  for (G4int ii=0; ii<n; ii++)
805  {
806  if (hkkevt_.isthkk[ii]==1 || hkkevt_.isthkk[ii]==-1)
807  {
808 //
809 // Particle is a final state secondary and not a nucleus.
810 // Determine what this secondary particle is, and if valid, load dynamic
811 // parameters.
812 //
813  G4ParticleDefinition* theParticle =
814  theParticleTable->FindParticle(hkkevt_.idhkk[i]);
815  if (theParticle)
816  {
817  G4double px = hkkevt_.phkk[i][0] * GeV;
818  G4double py = hkkevt_.phkk[i][1] * GeV;
819  G4double pz = hkkevt_.phkk[i][2] * GeV;
820  G4double et = hkkevt_.phkk[i][3] * GeV;
821 // G4LorentzVector lv = transformToLab * G4LorentzVector(px,py,pz,et);
822  G4LorentzVector lv = G4LorentzVector(px,py,pz,et);
823 
824  G4DynamicParticle *theDynamicParticle =
825  new G4DynamicParticle(theParticle,lv);
826  theParticleChange.AddSecondary (theDynamicParticle);
827 
828  if (verboseLevel >= 2)
829  DumpVerboseInformation2 (theParticle->GetParticleName(),
830  lv.vect(), et, theDynamicParticle->GetKineticEnergy(),pP);
831  }
832  }
833  else if (hkkevt_.idhkk[i]==80000 && hkkevt_.isthkk[i]==1001)
834  {
835 //
836 //
837 // Particle is a secondary nucleus. Determine the details of the nuclear
838 // fragment prior to de-excitation. (Note that the 1 eV in the total energy
839 // is a safety factor to avoid any possibility of negative rest mass energy.)
840 // Note also that we don't full trust the energy provided by the DPMJET-II.5,
841 // and there it's based on the Geant4-determined ion rest-mass.
842 //
843  G4int nucA = extevt_.idres[i];
844  G4int nucZ = extevt_.idxres[i];
845  if (nucA>0 && nucZ>0) {
846  M1++;
847  fragment = 0;
848  G4double px = hkkevt_.phkk[i][0] * GeV;
849  G4double py = hkkevt_.phkk[i][1] * GeV;
850  G4double pz = hkkevt_.phkk[i][2] * GeV;
851  G4double ionMass = theIonTable->GetIonMass(nucZ,nucA);
852  //GetIonMass(nucZ,nucA) + nucex[i]; // check how to get this energy
853  G4double dpmMass = hkkevt_.phkk[i][4] * GeV;
854  if (dpmMass > ionMass) ionMass = dpmMass;
855  G4double et = std::sqrt(px*px + py*py + pz*pz + ionMass*ionMass);
856  G4LorentzVector lv = G4LorentzVector(px,py,pz,et+1.0*eV);
857  fragment = new G4Fragment(nucA, nucZ, lv);
858  if (verboseLevel >= 2)
859  DumpVerboseInformation3 (M1, nucA, nucZ, lv.vect(),
860  et, et-ionMass, pP);
861 //
862 //
863 // Now we can decay the nuclear fragment if present. The secondaries are
864 // collected and boosted as well. The priority is to use a pre-compound
865 // de-excitation, otherwise the standard excitation-handler is used.
866 //
867  if (fragment && (thePreComp || theExcitationHandler))
868  {
869  G4ReactionProductVector *products = 0;
870  if (thePreComp && fragment->GetA() > 1.5)
871  products = thePreComp->DeExcite(*fragment);
872  else
873  products = theExcitationHandler->BreakItUp(*fragment);
874  delete fragment;
875  fragment = 0;
876  G4ReactionProductVector::iterator iter;
877  for (iter = products->begin(); iter != products->end(); ++iter)
878  {
879  G4DynamicParticle *secondary =
880  new G4DynamicParticle((*iter)->GetDefinition(),
881  (*iter)->GetTotalEnergy(), (*iter)->GetMomentum());
882  theParticleChange.AddSecondary (secondary);
883  G4String particleName =
884  (*iter)->GetDefinition()->GetParticleName();
885  delete (*iter);
886  if (verboseLevel >= 2) {
887  if (particleName.find("[",0) < particleName.size())
888  DumpVerboseInformation4 (m, particleName,
889  secondary->GetMomentum(),
890  secondary->GetTotalEnergy(),
891  secondary->GetKineticEnergy(), pP);
892  else DumpVerboseInformation2 (particleName,
893  secondary->GetMomentum(),
894  secondary->GetTotalEnergy(),
895  secondary->GetKineticEnergy(),
896  pP);
897  }
898  }
899  delete products;
900  }
901  if (fragment != 0)
902  {
903 //
904 //
905 // Add the excited fragment to the product vector. Note that this is temporary
906 // since we should at the atomic excitation to strip all electrons ... i.e. it's
907 // actually a bare nucleus not an atom in the ground state.
908 //
909  G4ParticleDefinition *theParticleDefinition = theIonTable->
910  GetIon(nucZ,nucA);
911  G4DynamicParticle *theDynamicParticle =
912  new G4DynamicParticle(theParticleDefinition,lv);
913  theParticleChange.AddSecondary (theDynamicParticle);
914  delete fragment;
915  fragment = 0;
916  }
917  }
918  }
919  }
920 
921  if (verboseLevel >= 3) {
922 //
923 //
924 // Calculate and display the energy and momenta before and after the collision.
925 // Everything is calculated for the lab frame.
926 //
927  G4double TotalEPost = 0.0;
928  G4ThreeVector TotalPPost;
929  G4double charge = 0.0;
930  G4int baryon = 0;
931  G4int lepton = 0;
932 // G4int parity = 0;
933  G4int nSecondaries = theParticleChange.GetNumberOfSecondaries();
934  for (G4int j=0; j<nSecondaries; j++) {
935  TotalEPost += theParticleChange.GetSecondary(j)->
936  GetParticle()->GetTotalEnergy();
937  TotalPPost += theParticleChange.GetSecondary(j)->
938  GetParticle()->GetMomentum();
939  G4ParticleDefinition *theParticle = theParticleChange.GetSecondary(j)->
940  GetParticle()->GetDefinition();
941  charge += theParticle->GetPDGCharge();
942  baryon += theParticle->GetBaryonNumber();
943  lepton += theParticle->GetLeptonNumber();
944 // parity += theParticle->GetPDGiParity();
945  }
946  G4cout <<"----------------------------------------"
947  <<"----------------------------------------"
948  <<G4endl;
949  G4cout <<"Total energy before collision = " <<TotalEPre/MeV
950  <<" MeV" <<G4endl;
951  G4cout <<"Total energy after collision = " <<TotalEPost/MeV
952  <<" MeV" <<G4endl;
953  G4cout <<"Total momentum before collision = " <<pP/MeV
954  <<" MeV/c" <<G4endl;
955  G4cout <<"Total momentum after collision = " <<TotalPPost/MeV
956  <<" MeV/c" <<G4endl;
957  if (verboseLevel >= 4) {
958  G4cout <<"Total charge before collision = " <<(ZP+ZT)*eplus
959  <<G4endl;
960  G4cout <<"Total charge after collision = " <<charge
961  <<G4endl;
962  G4cout <<"Total baryon number before collision = "<<AP+AT
963  <<G4endl;
964  G4cout <<"Total baryon number after collision = "<<baryon
965  <<G4endl;
966  G4cout <<"Total lepton number before collision = 0"
967  <<G4endl;
968  G4cout <<"Total lepton number after collision = "<<lepton
969  <<G4endl;
970  }
971  G4cout <<"----------------------------------------"
972  <<"----------------------------------------"
973  <<G4endl;
974  }
975 
976  if (verboseLevel >= 2)
977  G4cout <<"########################################"
978  <<"########################################"
979  <<G4endl;
980 
981  return &theParticleChange;
982 }
983 ////////////////////////////////////////////////////////////////////////////////
984 //
985 // SetNoDeexcitation
986 //
987 // Deletes an exiting de-excitation handlers and zeros the pointer. This
988 // allows the simulation to run without any de-excitation, or just pre-compound
989 // model only.
990 //
991 // Note that you need to separately SetNoPreCompoundModel, if you REALLY don't
992 // want any nuclear de-excitation. But please only do this to understand the
993 // contribution of de-excitation/pre-equilibrium to the full simulation.
994 // Running without de-excitation or pre-compound is physically unrealistic!
995 //
996 //
997 void G4DPMJET2_5Model::SetNoDeexcitation ()
998 {
999  if (theExcitationHandler)
1000  {
1001  delete theExcitationHandler;
1002  theExcitationHandler = 0;
1003  }
1004 }
1005 ////////////////////////////////////////////////////////////////////////////////
1006 //
1007 // SetNoPreCompoundModel
1008 //
1009 // Deletes an exiting pre-equilibrium model and zeros the pointer. This
1010 // allows the simulation to run without any pre-compound, or just -de-excitation
1011 // model only.
1012 //
1013 // Note that you need to separately SetNoDeexcitation, if you REALLY don't want
1014 // any nuclear de-excitation. But please only do this to understand the
1015 // contribution of de-excitation/pre-equilibrium to the full simulation.
1016 // Running without de-excitation or pre-compound is physically unrealistic!
1017 //
1018 //
1019 void G4DPMJET2_5Model::SetNoPreCompoundModel ()
1020 {
1021  if (thePreComp)
1022  {
1023  delete thePreComp;
1024  thePreComp = 0;
1025  }
1026 }
1027 ////////////////////////////////////////////////////////////////////////////////
1028 //
1029 // SetDefaultDeexcitation
1030 //
1031 // Note that this is used by the default constructor, as well as can be called
1032 // directly by the user.
1033 //
1034 void G4DPMJET2_5Model::SetDefaultDeexcitation ()
1035 {
1036 // SetNoDeexcitation();
1037 
1038  theExcitationHandler = new G4ExcitationHandler;
1039  G4Evaporation * theEvaporation = new G4Evaporation;
1040  G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
1041  G4PhotonEvaporation* thePhotonEvap = new G4PhotonEvaporation;
1042  theExcitationHandler->SetEvaporation(theEvaporation);
1043  theExcitationHandler->SetFermiModel(theFermiBreakUp);
1044  theExcitationHandler->SetMaxAandZForFermiBreakUp(17, 9);
1045  theExcitationHandler->SetFermiModel(theFermiBreakUp);
1046  theExcitationHandler->SetPhotonEvaporation(thePhotonEvap);
1047 }
1048 ////////////////////////////////////////////////////////////////////////////////
1049 //
1050 // SetDefaultPreCompoundModel
1051 //
1052 // Note that this is used by the default constructor, as well as can be called
1053 // directly by the user.
1054 //
1055 void G4DPMJET2_5Model::SetDefaultPreCompoundModel ()
1056 {
1057 // SetNoPreCompoundModel();
1058 
1059  G4ExcitationHandler *anExcitationHandler = new G4ExcitationHandler;
1060  G4Evaporation * theEvaporation = new G4Evaporation;
1061  G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
1062  G4PhotonEvaporation* thePhotonEvap = new G4PhotonEvaporation;
1063  anExcitationHandler->SetEvaporation(theEvaporation);
1064  anExcitationHandler->SetFermiModel(theFermiBreakUp);
1065  anExcitationHandler->SetMaxAandZForFermiBreakUp(17, 9);
1066  anExcitationHandler->SetPhotonEvaporation(thePhotonEvap);
1067 
1068  thePreComp = new G4PreCompoundModel(anExcitationHandler);
1069 }
1070 ////////////////////////////////////////////////////////////////////////////////
1071 //
1072 // SetVerboseFortranOutput
1073 //
1074 G4bool G4DPMJET2_5Model::SetVerboseFortranOutput (const G4String filename)
1075 {
1076  g4dpmjet_close_fort6_ ();
1077  if (filename == "" || filename == "stdo" ||
1078  filename == "stdout" || filename == "std::out" )
1079  {
1080  verboseFortranFile = "std::out";
1081  return true;
1082  }
1083  else
1084  {
1085  G4int namelen = filename.length();
1086  char *ptr = new char[namelen+1];
1087  filename.copy(ptr,namelen,0);
1088  ptr[namelen] = '\0';
1089 // char *ptr = 0;
1090 // ptr = const_cast<char*> (filename.c_str());
1091  ftnlogical opened = LFALSE;
1092  g4dpmjet_open_fort6_ (&namelen, &opened, ptr);
1093  delete [] ptr;
1094  if (opened == LTRUE) {
1095  verboseFortranFile = filename;
1096  return true;
1097  } else {
1098  verboseFortranFile = "std::out";
1099  return false;
1100  }
1101  }
1102 }
1103 ////////////////////////////////////////////////////////////////////////////////
1104 //
1105 // PrintWelcomeMessage
1106 //
1107 void G4DPMJET2_5Model::PrintWelcomeMessage () const
1108 {
1109  G4cout <<G4endl;
1110  G4cout <<" *****************************************************************"
1111  <<G4endl;
1112  G4cout <<" Interface to DPMJET2.5 for nuclear-nuclear interactions activated"
1113  <<G4endl;
1114  G4cout <<" Version number : 00.00.0B File date : 23/05/08" <<G4endl;
1115  G4cout <<" (Interface written by QinetiQ Ltd for the European Space Agency)"
1116  <<G4endl;
1117  G4cout <<G4endl;
1118  G4cout <<" Initialisation of DPMJET-II.5 variables will be according to "
1119  <<theInitType <<G4endl;
1120  G4cout <<" *****************************************************************"
1121  <<G4endl;
1122  G4cout << G4endl;
1123 
1124  return;
1125 }
1126 ////////////////////////////////////////////////////////////////////////////////
1127 //
1128 // DumpVerboseInformation1
1129 //
1130 // Dumps raw information about the DPMJET-II.5 simulation if verbosity set
1131 // to 4 or more.
1132 //
1133 void G4DPMJET2_5Model::DumpVerboseInformation1 (const G4int n) const
1134 {
1135  G4cout <<"----------------------------------------"
1136  <<"----------------------------------------" <<G4endl;
1137  G4cout <<n <<" INTERMEDIATE AND FINAL-STATE SECONDARIES PRODUCED" <<G4endl;
1138  if (verboseLevel >= 4)
1139  {
1140  G4cout <<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"
1141  <<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" <<G4endl;
1142  G4cout <<"ORIGINAL DPMJET-II.5 OUTPUT FOR EVENT:" <<G4endl;
1143  G4cout <<"Note that (1) the particles are yet to be transformed according"
1144  <<G4endl;
1145  G4cout <<" to incident particle direction" <<G4endl;
1146  G4cout <<" (2) the units of energy, momentum and mass are GeV,"
1147  <<G4endl;
1148  G4cout <<" GeV/c and GeV/c^2 respectively" <<G4endl;
1149  G4cout <<" I"
1150  <<" ISTHKK"
1151  <<" IDHKK"
1152  <<" IDRES"
1153  <<" IDXRES"
1154  <<" PX"
1155  <<" PY"
1156  <<" PZ"
1157  <<" TOTAL ENERGY"
1158  <<" MASS"
1159  <<G4endl;
1160  for (G4int i=0; i<n; i++)
1161  {
1162  G4cout.unsetf(std::ios::scientific);
1163  G4cout.setf(std::ios::fixed|std::ios::right|std::ios::adjustfield);
1164  G4cout.precision(0);
1165  G4cout <<std::setw(5) <<i
1166  <<std::setw(10) <<hkkevt_.isthkk[i]
1167  <<std::setw(10) <<hkkevt_.idhkk[i]
1168  <<std::setw(10) <<extevt_.idres[i]
1169  <<std::setw(10) <<extevt_.idxres[i];
1170  G4cout.unsetf(std::ios::fixed);
1171  G4cout.setf(std::ios::scientific|std::ios::right|std::ios::adjustfield);
1172  G4cout.precision(7);
1173  G4cout <<std::setw(15) <<hkkevt_.phkk[i][0]
1174  <<std::setw(15) <<hkkevt_.phkk[i][1]
1175  <<std::setw(15) <<hkkevt_.phkk[i][2]
1176  <<std::setw(15) <<hkkevt_.phkk[i][3]
1177  <<std::setw(15) <<hkkevt_.phkk[i][4]
1178  <<G4endl;
1179  }
1180  G4cout <<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"
1181  <<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" <<G4endl;
1182  }
1183  G4cout.setf(std::ios::fixed);
1184  G4cout <<" THE FOLLOWING LISTS ONLY THE FINAL-STATE SECONDARIES" <<G4endl;
1185 }
1186 
1187 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1188 
1189 void G4DPMJET2_5Model::DumpVerboseInformation2
1190  (const G4String particleName, const G4ThreeVector p,
1191  const G4double E, const G4double T, const G4ThreeVector pinit) const
1192 {
1193  G4cout <<"Name = " <<particleName <<G4endl;
1194  G4cout <<" Momentum = " <<p/MeV <<" MeV/c" <<G4endl;
1195  G4cout <<" T. Energy = " <<E/MeV <<" MeV" <<G4endl;
1196  G4cout <<" K. Energy = " <<T/MeV <<" MeV" <<G4endl;
1197  if (verboseLevel >= 3)
1198  {
1199  G4ThreeVector axis = pinit.unit();
1200  G4double pz = p.dot(axis);
1201  G4cout <<" Transverse mass = " <<std::sqrt(E*E-pz*pz)/MeV
1202  <<" MeV"
1203  <<G4endl;
1204  G4cout <<" Rapidity = "
1205  <<0.5*std::log((E+pz)/(E-pz)) <<G4endl;
1206  }
1207 }
1208 
1209 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1210 
1211 void G4DPMJET2_5Model::DumpVerboseInformation3 (const G4int i,
1212  const G4int A, const G4int Z, const G4ThreeVector p,
1213  const G4double E, const G4double T, const G4ThreeVector pinit) const
1214 {
1215  G4cout <<"----------------------------------------"
1216  <<"----------------------------------------" <<G4endl;
1217  G4cout <<"The nuclear fragment #" <<i <<" before" <<G4endl;
1218  G4cout <<"----------------------------------------"
1219  <<"----------------------------------------" <<G4endl;
1220 
1221  std::ostringstream tmpStream;
1222  tmpStream <<"(A = " <<A <<", Z = " <<Z <<")";
1223  G4String AZ = tmpStream.str();
1224 
1225  DumpVerboseInformation2(AZ, p, E, T, pinit);
1226 }
1227 
1228 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1229 
1230 void G4DPMJET2_5Model::DumpVerboseInformation4 (const G4int i,
1231  const G4String particleName, const G4ThreeVector p,
1232  const G4double E, const G4double T, const G4ThreeVector pinit) const
1233 {
1234  G4cout <<"----------------------------------------"
1235  <<"----------------------------------------" <<G4endl;
1236  G4cout <<"The nuclear fragment #" <<i <<" after" <<G4endl;
1237  G4cout <<"----------------------------------------"
1238  <<"----------------------------------------" <<G4endl;
1239 
1240  DumpVerboseInformation2(particleName, p, E, T, pinit);
1241 }
1242 
1243 ////////////////////////////////////////////////////////////////////////////////
1244 //
1245 // Initialise
1246 // This is intended to do exactly what the main program and subroutine DMINIT
1247 // attempt to do with variable initialisation. i've included most of the
1248 // FORTRAN source code for reference. Note that I'm certain that many of these
1249 // lines relate to tallying of events for output by the standalone version
1250 // of DPMJET-II.5, but for the moment, I'd rather initialise those variables,
1251 // just in case not doing so has an adverse effect.
1252 //
1253 void G4DPMJET2_5Model::Initialise ()
1254 {
1255 //
1256 //
1257 // This first line is intended to make sure the block data statements are
1258 // executed, since we're not running from a FORTRAN main program.
1259 //
1260  g4dpmjet_initialise_block_data_ ();
1261 //
1262 //
1263  dpar_.aam[4] = 0.001; // AAM(5)=0.001D0
1264  dpar_.aam[5] = 0.001; // AAM(6)=0.001D0
1265  dpar_.aam[132] = 0.001; // AAM(133)=0.001D0
1266  dpar_.aam[133] = 0.001; // AAM(134)=0.001D0
1267  dpar_.aam[134] = 0.001; // AAM(135)=0.001D0
1268  dpar_.aam[135] = 0.001; // AAM(136)=0.001D0
1269 
1270  vxsvd_.nxsp = 0; // NXSP=0
1271  vxsvd_.nxst = 0; // NXST=0
1272  vxsvd_.nxsap = 0; // NXSAP=0
1273  vxsvd_.nxsat = 0; // NXSAT=0
1274  vxsvd_.nxvp = 0; // NXVP=0
1275  vxsvd_.nxvt = 0; // NXVT=0
1276  vxsvd_.nxdp = 0; // NXDP=0
1277  vxsvd_.nxdt = 0; // NXDT=0
1278 
1279  for (G4int i=0; i<50; i++)
1280  {
1281  vxsvd_.vxsp[i] = 1.0E-08; // VXSP(II)=1.D-8
1282  vxsvd_.vxst[i] = 1.0E-08; // VXST(II)=1.D-8
1283  vxsvd_.vxsap[i] = 1.0E-08; // VXSAP(II)=1.D-8
1284  vxsvd_.vxsat[i] = 1.0E-08; // VXSAT(II)=1.D-8
1285  vxsvd_.vxvp[i] = 1.0E-08; // VXVP(II)=1.D-8
1286  vxsvd_.vxvt[i] = 1.0E-08; // VXVT(II)=1.D-8
1287  vxsvd_.vxdp[i] = 1.0E-08; // VXDP(II)=1.D-8
1288  vxsvd_.vxst[i] = 1.0E-08; // VXST(II)=1.D-8
1289  }
1290 
1291  if (debug)
1292  {
1293 #ifdef G4VERBOSE
1294  if (GetVerboseLevel()>0) {
1295  G4cout <<"AT G4DPMJET2_5Model::Initialise:" <<G4endl;
1296  }
1297 #endif
1298  dprin_.ipri = debug_level; // IPRI = LEVLDB
1299  dprin_.ipev = debug_level; // IPEV = LEVLDB
1300  dprin_.ippa = debug_level; // IPPA = LEVLDB
1301  dprin_.ipco = debug_level; // IPCO = LEVLDB
1302  dprin_.init = debug_level; // INIT = LEVLDB
1303  dprin_.iphkk = debug_level; // IPHKK = LEVLDB
1304  pydat1_.mstu[25] = 10; // MSTU(26) = 10
1305  }
1306  else
1307  {
1308  dprin_.ipri = 0; // IPRI = 0
1309  dprin_.ipev = 0; // IPEV = 0
1310  dprin_.ippa = 0; // IPPA = 0
1311  dprin_.ipco =-2; // IPCO = -2
1312  dprin_.init = 0; // INIT = 0
1313  dprin_.iphkk = 0; // IPHKK = 0
1314  pydat1_.mstu[25] = 0; // MSTU(26) = 10
1315  }
1316 
1317  for (G4int i=0; i<7; i++)
1318  {
1319  diqrej_.idiqre[i] = 0;
1320  diqrej_.idiqrz[i] = 0;
1321  }
1322 
1323  for (G4int i=0; i<3; i++)
1324  {
1325  diqrej_.idvre[i] = 0;
1326  diqrej_.ivdre[i] = 0;
1327  diqrej_.idsre[i] = 0;
1328  diqrej_.isdre[i] = 0;
1329  diqrej_.idzre[i] = 0;
1330  diqrej_.izdre[i] = 0;
1331  }
1332 
1333  diqsum_.ndvuu = 0; // NDVUU = 0
1334  diqsum_.ndvus = 0; // NDVUS = 0
1335  diqsum_.ndvss = 0; // NDVSS = 0
1336  diqsum_.nvduu = 0; // NVDUU = 0
1337  diqsum_.nvdus = 0; // NVDUS = 0
1338  diqsum_.nvdss = 0; // NVDSS = 0
1339  diqsum_.ndsuu = 0; // NDSUU = 0
1340  diqsum_.ndsus = 0; // NDSUS = 0
1341  diqsum_.ndsss = 0; // NDSSS = 0
1342  diqsum_.nsduu = 0; // NSDUU = 0
1343  diqsum_.nsdus = 0; // NSDUS = 0
1344  diqsum_.nsdss = 0; // NSDSS = 0
1345  diqsum_.ndzuu = 0; // NDZUU = 0
1346  diqsum_.ndzus = 0; // NDZUS = 0
1347  diqsum_.ndzss = 0; // NDZSS = 0
1348  diqsum_.nzduu = 0; // NZDUU = 0
1349  diqsum_.nzdus = 0; // NZDUS = 0
1350  diqsum_.nzdss = 0; // NZDSS = 0
1351  diqsum_.nadvuu = 0; // NADVUU = 0
1352  diqsum_.nadvus = 0; // NADVUS = 0
1353  diqsum_.nadvss = 0; // NADVSS = 0
1354  diqsum_.navduu = 0; // NAVDUU = 0
1355  diqsum_.navdus = 0; // NAVDUS = 0
1356  diqsum_.navdss = 0; // NAVDSS = 0
1357  diqsum_.nadsuu = 0; // NADSUU = 0
1358  diqsum_.nadsus = 0; // NADSUS = 0
1359  diqsum_.nadsss = 0; // NADSSS = 0
1360  diqsum_.nasduu = 0; // NASDUU = 0
1361  diqsum_.nasdus = 0; // NASDUS = 0
1362  diqsum_.nasdss = 0; // NASDSS = 0
1363  diqsum_.nadzuu = 0; // NADZUU = 0
1364  diqsum_.nadzus = 0; // NADZUS = 0
1365  diqsum_.nadzss = 0; // NADZSS = 0
1366  diqsum_.nazduu = 0; // NAZDUU = 0
1367  diqsum_.nazdus = 0; // NAZDUS = 0
1368  diqsum_.nazdss = 0; // NAZDSS = 0
1369  hdjase_.nhse1 = 0; // NHSE1 = 0
1370  hdjase_.nhse2 = 0; // NHSE2 = 0
1371  hdjase_.nhse3 = 0; // NHSE3 = 0
1372  hdjase_.nhase1 = 0; // NHASE1 = 0
1373  hdjase_.nhase2 = 0; // NHASE2 = 0
1374  hdjase_.nhase3 = 0; // NHASE3 = 0
1375 //
1376 //
1377 // Parton pt distribution.
1378 //
1379  G4int i = 1;
1380  G4double pt1 = 0.0;
1381  G4double pt2 = 0.0;
1382  G4int ipt = 0;
1383  G4int nevt = 0;
1384  parpt_ (&i,&pt1,&pt2,&ipt,&nevt);
1385  // CALL PARPT(1,PT1,PT2,IPT,NEVT)
1386 //
1387 //
1388 // Initialise BAMJET, DECAY and HADRIN.
1389 //
1390  ddatar_ (); // CALL DDATAR
1391  dhadde_ (); // CALL DHADDE
1392  dchant_ (); // CALL DCHANT
1393  dchanh_ (); // CALL DCHANH
1394 
1395  G4double epn = 0.0;
1396  G4double ppn = 0.0;
1397  defaul_ (&epn,&ppn); // CALL DEFAUL(EPN,PPN)
1398  defaux_ (&epn,&ppn); // CALL DEFAUX(EPN,PPN)
1399 
1400  coulo_.icoul = 1; // ICOUL = 1
1401  nuclea_.icoull = 1; // ICOULL = 1
1402  edens_.ieden = 0; // IEDEN = 0
1403  dprin_.itopd = 0; // ITOPD = 0
1404 
1405  if ( theInitType == DEFAULT ) {
1406  TAUFOR = 5.0E+00;
1407  KTAUGE = 25;
1408  } else if ( theInitType == CORSIKA ) {
1409  TAUFOR = 5.0E+00; // TAUFOR = 5.D0
1410  KTAUGE = 25; // KTAUGE = 25
1411  } else if ( theInitType == DPM2_5 ) {
1412  TAUFOR = 105.0E+00; // TAUFOR = 105.D0
1413  KTAUGE = 10; // KTAUGE = 10
1414  } else if ( theInitType == DPM3 ) {
1415  TAUFOR = 3.5E+00;
1416  KTAUGE = 10;
1417  }
1418  taufo_.taufor = TAUFOR;
1419  taufo_.ktauge = KTAUGE;
1420 
1421  ITAUVE = 1; // ITAUVE = 1
1422  taufo_.itauve = ITAUVE;
1423  taufo_.incmod = 1; // INCMOD = 1
1424 //
1425 //
1426 // Definition of soft quark distributions, Fermi, Pauli
1427 //
1428  G4double xseaco = 1.0; // XSEACO = 1.00D0
1429  xseadi_.xseacu = 1.05 - xseaco; // XSEACU = 1.05D0-XSEACO
1430 
1431  if ( theInitType == DEFAULT ) {
1432  UNON = 3.50;
1433  UNOM = 1.11;
1434  UNOSEA = 5.0;
1435  droppt_.fermp = LTRUE;
1436  nucimp_.fermod = 0.6;
1437  } else if ( theInitType == CORSIKA ) {
1438  UNON = 3.50; // UNON = 3.50D0
1439  UNOM = 1.11; // UNOM = 1.11D0
1440  UNOSEA = 5.0; // UNOSEA = 5.0D0
1441  droppt_.fermp = LTRUE; // FERMP = .TRUE.
1442  nucimp_.fermod = 0.6; // FERMOD = 0.6D0
1443  } else if ( theInitType == DPM2_5 ) {
1444  UNON = 3.50; // UNON = 3.50D0
1445  UNOM = 1.11; // UNOM = 1.11D0
1446  UNOSEA = 5.0; // UNOSEA = 5.0D0
1447  droppt_.fermp = LTRUE; // FERMP = .TRUE.
1448  nucimp_.fermod = 0.6; // FERMOD = 0.6D0
1449  } else if ( theInitType == DPM3 ) {
1450  UNON = 2.00;
1451  UNOM = 1.5;
1452  UNOSEA = 5.0;
1453  droppt_.fermp = LTRUE;
1454  nucimp_.fermod = 0.55;
1455  }
1456  xseadi_.unon = UNON;
1457  xseadi_.unom = UNOM;
1458  xseadi_.unosea = UNOSEA;
1459 
1460  nuclea_.fermdd = 0.6; // FERMDD = 0.6D0
1461  ferfor_.iferfo = 1; // IFERFO = 1
1462  dprin_.ipaupr = 0; // IPAUPR = 0
1463  droppt_.lpauli = LTRUE; // LPAULI = .TRUE.
1464 //
1465 //
1466 // Definition of cuts for x-sampling.
1467 //
1468  if ( theInitType == DEFAULT ) {
1469  CVQ = 1.8;
1470  CDQ = 2.0;
1471  CSEA = 0.5;
1472  SSMIMA = 0.9;
1473  } else if ( theInitType == CORSIKA ) {
1474  CVQ = 1.8; // CVQ = 1.8D0
1475  CDQ = 2.0; // CDQ = 2.0D0
1476  CSEA = 0.5; // CSEA = 0.5D0
1477  SSMIMA = 0.901; // SSMIMA = 0.901D0
1478  } else if ( theInitType == DPM2_5 ) {
1479  CVQ = 1.8; // CVQ = 1.8D0
1480  CDQ = 2.0; // CDQ = 2.0D0
1481  CSEA = 0.5; // CSEA = 0.5D0
1482  SSMIMA = 1.201; // SSMIMA = 1.201D0
1483  } else if ( theInitType == DPM3 ) {
1484  CVQ = 1.0;
1485  CDQ = 2.0;
1486  CSEA = 0.1;
1487  SSMIMA = 0.14;
1488  }
1489  xseadi_.cvq = CVQ;
1490  xseadi_.cdq = CDQ;
1491  xseadi_.csea = CSEA;
1492  xseadi_.ssmima = SSMIMA;
1493 
1494  xseadi_.ssmimq = xseadi_.ssmima * xseadi_.ssmima;
1495  // SSMIMQ = SSMIMA**2
1496  if ( theInitType == DEFAULT ) {
1497  VVMTHR = 0.0;
1498  } else if ( theInitType == CORSIKA ) {
1499  VVMTHR = 0.0; // VVMTHR = 0.D0
1500  } else if ( theInitType == DPM2_5 ) {
1501  VVMTHR = 0.0; // VVMTHR = 0.D0
1502  } else if ( theInitType == DPM3 ) {
1503  VVMTHR = 2.0;
1504  }
1505  xseadi_.vvmthr = VVMTHR;
1506 //
1507 //
1508 // There is a final call. Set recombin, seasu3, coninpt, allpart and interdpm
1509 //
1510  final_.ifinal = 0; // IFINAL = 0
1511  recom_.irecom = 0; // IRECOM = 0
1512  seadiq_.lseadi = LTRUE; // LSEADI = .TRUE.
1513 
1514  if ( theInitType == DEFAULT ) {
1515  SEASQ = 0.5;
1516  MKCRON = 1;
1517  CRONCO = 0.64;
1518  } else if ( theInitType == CORSIKA ) {
1519  SEASQ = 0.5; // SEASQ = 0.50D0
1520  MKCRON = 0; // MKCRON = 0
1521  CRONCO = 0.0; // CRONCO = 0.00D0
1522  } else if ( theInitType == DPM2_5 ) {
1523  SEASQ = 0.5; // SEASQ = 0.50D0
1524  MKCRON = 1; // MKCRON = 1
1525  CRONCO = 0.64; // CRONCO = 0.64D0
1526  } else if ( theInitType == DPM3 ) {
1527  SEASQ = 1.0;
1528  MKCRON = 1;
1529  CRONCO = 0.64;
1530  }
1531  seasu3_.seasq = SEASQ;
1532  cronin_.mkcron = MKCRON;
1533  cronin_.cronco = CRONCO;
1534 
1535  droppt_.ihada = LTRUE; // IHADA = .TRUE.
1536 
1537 // inxdpm_.intdpm = 0; // INTDPM = 0
1538 // Note FORTRAN initialises the variable IROEH = 0, but this isn't used
1539 // nor is it contained within a common block.
1540 //
1541 //
1542 // Definition for popcork, casadiqu, popcorse
1543 //
1544  popcck_.pdbck = 0.0; // PDBCK = 0.D0
1545  popcck_.ijpock = 0; // IJPOCK = 0
1546  casadi_.icasad = 1; // ICASAD = 1
1547  casadi_.casaxx = 0.05;
1548  // CASAXX = 0.05D0 ! corrected Nov. 2001
1549  popcck_.pdbse = 0.45;
1550  // PDBSE = 0.45D0 ! with baryon stopping
1551  popcck_.pdbseu = 0.45;
1552  // PDBSEU = 0.45D0 ! with baryon stopping
1553  // PDBSE = 0.D0 ! without baryon stopping
1554  // PDBSEU = 0.D0 ! without baryon stopping
1555  popcck_.irejck = 0; // IREJCK = 0
1556  popcck_.irejse = 0; // IREJSE = 0
1557  popcck_.irejs3 = 0; // IREJS3 = 0
1558  popcck_.irejs0 = 0; // IREJS0 = 0
1559  popcck_.ick4 = 0; // ICK4 = 0
1560  popcck_.ise4 = 0; // ISE4 = 0
1561  popcck_.ise43 = 0; // ISE43 = 0
1562  popcck_.ihad4 = 0; // IHAD4 = 0
1563  popcck_.ick6 = 0; // ICK6 = 0
1564  popcck_.ise6 = 0; // ISE6 = 0
1565  popcck_.ise63 = 0; // ISE63 = 0
1566  popcck_.ihad6 = 0; // IHAD6 = 0
1567  popcck_.irejsa = 0; // IREJSA = 0
1568  popcck_.ireja3 = 0; // IREJA3 = 0
1569  popcck_.ireja0 = 0; // IREJA0 = 0
1570  popcck_.isea4 = 0; // ISEA4 = 0
1571  popcck_.isea43 = 0; // ISEA43 = 0
1572  popcck_.ihada4 = 0; // IHADA4 = 0
1573  popcck_.isea6 = 0; // ISEA6 = 0
1574  popcck_.isea63 = 0; // ISEA63 = 0
1575  popcck_.ihada6 = 0; // IHADA6 = 0
1576 
1577  if ( theInitType == DEFAULT || theInitType == CORSIKA ) {
1578  popcor_.pdb = 0.1; // PDB = 0.10D0
1579  } else if ( theInitType == DPM2_5 ) {
1580  popcor_.pdb = 0.1; // PDB = 0.10D0
1581  } else if ( theInitType == DPM3 ) {
1582  popcor_.pdb = 0.15;
1583  }
1584  popcor_.ajsdef = 0.0; // AJSDEF = 0.D0
1585 //
1586 //
1587 // Definition of fluctuat, intpt, hadroniz, diquarks, singlech, evapor
1588 // (Charmed particles set to decay : IHADRINZ>=2
1589 //
1590  fluctu_.ifluct = 0; // IFLUCT = 0
1591  droppt_.intpt = LTRUE; // INTPT = .TRUE.
1592  colle_.ihadrz = 2; // IHADRZ = 2
1593  ifragm_.ifrag = 1; // IFRAG = 1
1594  promu_.ipromu = 1; // IPROMU = 1
1595  if (colle_.ihadrz >= 2)
1596  {
1597  ifragm_.ifrag = colle_.ihadrz - 1;
1598  // IFRAG = IHADRZ-1
1599  lundin_ (); // CALL LUNDIN
1600  }
1601  diquax_.idiqua = 1; // IDIQUA = 1
1602  diquax_.idiquu = 1; // IDIQUU = 1
1603  diquax_.amedd = 0.9; // AMEDD = 0.9D0
1604  sincha_.isicha = 0; // ISICHA = 0
1605  evappp_.ievap = 0; // IEVAP = 0
1606  seaqxx_.seaqx = 0.5; // SEAQX = 0.5D0
1607  seaqxx_.seaqxn = 0.5; // SEAQXN = 0.5D0
1608  kglaub_.jglaub = 2; // JGLAUB = 2
1609  hadthr_.ehadth = 5.0; // EHADTH = 5.D0
1610 //
1611 //
1612 // Definitions for hbook, pomtable, cmhisto, central, strucfun
1613 //
1614  hboo_.ihbook = 1; // IHBOOK = 1
1615  pomtab_.ipomta = 1; // IPOMTA = 1
1616  cmhico_. cmhis = 0.0;
1617  // CMHIS = 0.0D+00 ! Lab System
1618  zentra_.icentr = 0; // ICENTR = 0
1619  user2_.istruf = 222; // ISTRUF = 222
1620  strufu_.istrum = 0; // ISTRUM = 0
1621  strufu_.istrut = user2_.istruf / 100;
1622  // ISTRUT = ISTRUF/100
1623  user2_.istruf = user2_.istruf - strufu_.istrut*100;
1624  // ISTRUF = ISTRUF-ISTRUT*100
1625  strufu_.istrum = user2_.istruf; // ISTRUM = ISTRUF
1626 
1627  if ( theInitType == DEFAULT ) {
1628  ISINGD = 1;
1629  ISINGX = 1;
1630  IDUBLD = 0;
1631  SDFRAC = 1.0;
1632  } else if ( theInitType == CORSIKA ) {
1633  ISINGD = 0; // ISINGD = 0
1634  ISINGX = 0; // ISINGX = 0
1635  IDUBLD = 0; // IDUBLD = 0
1636  SDFRAC = 0.0; // SDFRAC = 0.
1637  } else if ( theInitType == DPM2_5 ) {
1638  ISINGD = 1; // ISINGD = 1
1639  ISINGX = 1; // ISINGX = 1
1640  IDUBLD = 0; // IDUBLD = 0
1641  SDFRAC = 1.0; // SDFRAC = 1.
1642  } else if ( theInitType == DPM3 ) {
1643  ISINGD = 0;
1644  ISINGX = 1;
1645  IDUBLD = 0;
1646  SDFRAC = 1.0;
1647  }
1648  diffra_.isingd = ISINGD;
1649  user2_.isingx = ISINGX;
1650  user2_.idubld = IDUBLD;
1651  user2_.sdfrac = SDFRAC;
1652 //
1653 //
1654 // Definitions for start, inforeje, gluxplit, partev, sampt
1655 // NEVNTS passed to DTMAI as argument, NEVHAD in COMMON
1656 //
1657  colle_.nfile = 0; // NFILE = 0
1658  nstari_.nstart = 1; // NSTART = 1
1659  stars_.istar2 = 0; // ISTAR2 = 0
1660  stars_.istar3 = 0; // ISTAR3 = 0
1661  user2_.ptlar = 2.0; // PTLAR = 2.D0
1662 
1663 // G4int iglaub = 0;
1664 // IGLAUB = 0 !! Note that this is just a local variable
1665 
1666 // infore_.ifrej = 0;
1667 // IFREJ = 0 Note rejection diagnostics not required
1668  gluspl_.nugluu = 1; // NUGLUU = 1
1669  gluspl_.nsgluu = 0; // NSGLUU = 0
1670  colle_.nvers = 1; // NVERS = 1
1671  ptsamp_.isampt = 4; // ISAMPT = 4
1672 //
1673 //
1674 // Definitions for selhard, sigmapom, pshow, secinter.
1675 //
1676  dropjj_.dropjt = 0.0; // DROPJT = 0.D0
1677  collis_.iophrd = 2; // IOPHRD = 2
1678  collis_.ptthr = 3.0; // PTTHR = 3.D0
1679 // collis_.ptthr2 = collis_.ptthr; // PTTHR2 = PTTHR
1680  user2_.cmener = 100.0; // CMENER = 100.D0
1681  if (strufu_.istrut == 1)
1682  {
1683  collis_.ptthr = 2.1 + 0.15*std::pow(std::log10(user2_.cmener/50.),3.0);
1684  // PTTHR = 2.1D0+0.15D0*(LOG10(CMENER/50.))**3
1685  }
1686  else if (strufu_.istrut == 2)
1687  {
1688  collis_.ptthr = 2.5 + 0.12*std::pow(std::log10(user2_.cmener/50.),3.0);
1689  // PTTHR = 2.5D0+0.12D0*(LOG10(CMENER/50.))**3
1690  }
1691  collis_.ptthr2 = collis_.ptthr; // PTTHR2 = PTTHR
1692  pomtyp_.ipim = 2; // IPIM = 2
1693  pomtyp_.icon = 48; // ICON = 48
1694  pomtyp_.isig = 10; // ISIG = 10
1695  pomtyp_.lmax = 30; // LMAX = 30
1696  pomtyp_.mmax = 100; // MMAX = 100
1697  pomtyp_.nmax = 2; // NMAX = 2
1698  pomtyp_.difel = 0.0; // DIFEL = 0.D0
1699  pomtyp_.difnu = 1.0; // DIFNU = 1.D0
1700  pshow_.ipshow = 1; // IPSHOW = 1
1701  secint_.isecin = 0; // ISECIN = 0
1702 //
1703 //
1704 // This next bit is associated with evaporation. I'm not sure if it's needed
1705 // as IEVAP = 0, but will initialise in any case.
1706 //
1707  if ( theInitType == DEFAULT ) {
1708  parevt_.levprt = LTRUE;
1709  parevt_.ilvmod = 1;
1710  parevt_.ldeexg = LFALSE;
1711  parevt_.lheavy = LFALSE;
1712  frbkcm_.lfrmbk = LFALSE;
1713  inpflg_.ifiss = 0;
1714  } else if ( theInitType == CORSIKA ) {
1715  parevt_.levprt = LTRUE; // LEVPRT = .TRUE.
1716  parevt_.ilvmod = 1; // ILVMOD = 1
1717  parevt_.ldeexg = LTRUE; // LDEEXG = .TRUE.
1718  parevt_.lheavy = LTRUE; // LHEAVY = .TRUE.
1719  frbkcm_.lfrmbk = LTRUE; // LFRMBK = .TRUE.
1720  inpflg_.ifiss = 0; // IFISS = 0
1721  } else if ( theInitType == DPM2_5 ) {
1722  parevt_.levprt = LFALSE; // LEVPRT = .FALSE.
1723  parevt_.ilvmod = 1; // ILVMOD = 1
1724  parevt_.ldeexg = LFALSE; // LDEEXG = .FALSE.
1725  parevt_.lheavy = LFALSE; // LHEAVY = .FALSE.
1726  frbkcm_.lfrmbk = LFALSE; // LFRMBK = .FALSE.
1727  inpflg_.ifiss = 0; // IFISS = 0
1728  } else if ( theInitType == DPM3 ) {
1729  parevt_.levprt = LFALSE;
1730  parevt_.ilvmod = 1;
1731  parevt_.ldeexg = LFALSE;
1732  parevt_.lheavy = LFALSE;
1733  frbkcm_.lfrmbk = LFALSE;
1734  inpflg_.ifiss = 0;
1735  }
1736 
1737  hettp_.nbertp = lunber; // NBERTP = LUNBER
1738 
1739  verboseFortranFile = "fort.6";
1740  G4int namelen = verboseFortranFile.length();
1741  char *ptr1 = new char[namelen+1];
1742  verboseFortranFile.copy(ptr1,namelen,0);
1743  ptr1[namelen] = '\0';
1744 // ptr1 = const_cast<char*> (verboseFortranFile.c_str());
1745  ftnlogical opened = LFALSE;
1746  g4dpmjet_open_fort6_ (&namelen, &opened, ptr1);
1747  delete [] ptr1;
1748  if (opened == LFALSE)
1749  {
1750  G4cout <<"ATTEMPTED TO OPEN fort.6 TO OUTPUT VERBOSE FORTRAN TEXT"
1751  <<G4endl;
1752  G4cout <<"HOWEVER THIS WAS NOT POSSIBLE" <<G4endl;
1753  }
1754 
1755 #ifdef G4VERBOSE
1756  if (GetVerboseLevel()>0) {
1757  G4cout <<"AT G4DPMJET2_5Model::Initialise: before NUCLEAR.BIN" <<G4endl;
1758  G4cout <<"OPENING NUCLEAR.BIN ON FILE UNIT " <<lunber <<G4endl;
1759  }
1760 #endif
1761  if ( !getenv("G4DPMJET2_5DATA") )
1762  {
1763  G4cout <<"ENVIRONMENT VARIABLE G4DPMJET2_5DATA NOT SET " <<G4endl;
1764  throw G4HadronicException(__FILE__, __LINE__,
1765  "Please setenv G4DPMJET2_5DATA to point to the dpmjet2.5 data files.");
1766  }
1767  defaultDirName = G4String(getenv("G4DPMJET2_5DATA")) + "/NUCLEAR.BIN";
1768  namelen = defaultDirName.length();
1769  ptr1 = new char[namelen+1];
1770  defaultDirName.copy(ptr1,namelen,0);
1771  ptr1[namelen] = '\0';
1772 // ptr1 = const_cast<char*> (defaultDirName.c_str());
1773  opened = LFALSE;
1774  g4dpmjet_open_nuclear_bin_ (&namelen, &lunber, &opened, ptr1);
1775  delete [] ptr1;
1776  if (opened == LFALSE)
1777  {
1778 //
1779 //
1780 // Problems with locating NUCLEAR.BIN file.
1781 //
1782  G4cout <<"NUCLEAR.BIN FILE NOT FOUND IN DIRECTORY " <<defaultDirName
1783  <<G4endl;
1784  throw G4HadronicException(__FILE__, __LINE__,
1785  "NUCLEAR.BIN file not present.");
1786  }
1787 #ifdef G4VERBOSE
1788  if (GetVerboseLevel()>0)
1789  G4cout <<"AT G4DPMJET2_5Model::Initialise: after NUCLEAR.BIN" <<G4endl;
1790 #endif
1791  G4cout << "CALL BERTTP" << G4endl;
1792  berttp_ (); // CALL BERTTP
1793  G4cout << "CALL BERTTP done" << G4endl;
1794  if (evappp_.ievap == 1)
1795  {
1796 #ifdef G4VERBOSE
1797  if (GetVerboseLevel()>0)
1798  G4cout <<"AT G4DPMJET2_5Model::Initialise: before INCINI" <<G4endl;
1799 #endif
1800  incini_ (); // CALL INCINI
1801 #ifdef G4VERBOSE
1802  if (GetVerboseLevel()>0)
1803  G4cout <<"AT G4DPMJET2_5Model::Initialise: after INCINI" <<G4endl;
1804 #endif
1805  }
1806  g4dpmjet_close_nuclear_bin_ (&lunber);
1807 #ifdef G4VERBOSE
1808  if (GetVerboseLevel()>0)
1809  G4cout <<"AT G4DPMJET2_5Model::Initialise: NUCLEAR.BIN closed" <<G4endl;
1810 #endif
1811 
1812  ptlarg_.xsmax = 0.8; // XSMAX = 0.8D0
1813  dprin_.itopd = 0; // ITOPD = 0
1814 
1815  G4int iseed1 = 0;
1816  G4int iseed2 = 0;
1817  rd2out_ (&iseed1,&iseed2);
1818 //
1819 //
1820 // This next bit outputs important variables if debug is switched on.
1821 //
1822 #ifdef G4VERBOSE
1823  if (GetVerboseLevel()>0) {
1824  G4cout <<"AT G4DPMJET2_5Model::Initialise:" <<G4endl;
1825  G4cout <<"Printout of important Parameters before DPMJET2.5" <<G4endl;
1826  G4cout <<"Please note for DPMJET input all numbers are floating point!"
1827  <<G4endl;
1828  G4cout <<"PROJPAR " <<nucc_.ip <<" "
1829  <<nucc_.ipz <<G4endl;
1830  G4cout <<"TARPAR " <<nucc_.it <<" "
1831  <<nucc_.itz <<G4endl;
1832  G4cout <<"MOMENTUM " <<ppn <<G4endl;
1833  G4cout <<"ENERGY " <<epn <<G4endl;
1834  G4cout <<"CMENERGY " <<nncms_.umo <<G4endl;
1835  G4cout <<"NOFINALE " <<final_.ifinal <<G4endl;
1836  G4cout <<"EVAPORAT " <<evappp_.ievap <<G4endl;
1837  G4cout <<"OUTLEVEL " <<dprin_.ipri <<" "
1838  <<dprin_.ipev <<" "
1839  <<dprin_.ippa <<" "
1840  <<dprin_.ipco <<" "
1841  <<dprin_.init <<" "
1842  <<dprin_.iphkk <<G4endl;
1843  G4cout <<"RANDOMIZ " <<iseed1 <<" "
1844  <<iseed2 <<G4endl;
1845  G4cout <<"STRUCFUN " <<user2_.istruf+100*strufu_.istrut <<G4endl;
1846  G4cout <<"SAMPT " <<ptsamp_.isampt <<G4endl;
1847  G4cout <<"SELHARD " <<0 <<" "
1848  <<collis_.iophrd <<" "
1849  <<0 <<" "
1850  <<dropjj_.dropjt <<" "
1851  <<collis_.ptthr <<" "
1852  <<collis_.ptthr2 << G4endl;
1853  G4cout <<"SIGMAPOM " <<0 <<" "
1854  <<pomtyp_.isig <<" "
1855  <<pomtyp_.ipim + 10*pomtyp_.icon <<" "
1856  <<pomtyp_.lmax <<" "
1857  <<pomtyp_.mmax <<" "
1858  <<pomtyp_.nmax <<G4endl;
1859  G4cout <<"PSHOWER " <<pshow_.ipshow <<G4endl;
1860  G4cout <<"CENTRAL " <<zentra_.icentr <<G4endl;
1861  G4cout <<"CMHISTO " <<cmhico_.cmhis <<G4endl;
1862  G4cout <<"SEASU3 " <<seasu3_.seasq <<G4endl;
1863  G4cout <<"RECOMBIN " <<recom_.irecom <<G4endl;
1864  G4cout <<"SINGDIFF " <<diffra_.isingd <<G4endl;
1865  G4cout <<"TAUFOR " <<taufo_.taufor <<" "
1866  <<taufo_.ktauge <<" "
1867  <<taufo_.itauve <<G4endl;
1868  G4cout <<"POPCORN " <<popcor_.pdb <<G4endl;
1869  G4cout <<"POPCORCK " <<popcck_.ijpock <<" "
1870  <<popcck_.pdbck <<G4endl;
1871  G4cout <<"POPCORSE " <<popcck_.pdbse <<" "
1872  <<popcck_.pdbseu <<G4endl;
1873  G4cout <<"CASADIQU " <<casadi_.icasad <<" "
1874  <<casadi_.casaxx <<G4endl;
1875  G4cout <<"DIQUARKS " <<diquax_.idiqua <<" "
1876  <<diquax_.idiquu <<" "
1877  <<diquax_.amedd <<G4endl;
1878  G4cout <<"HADRONIZ " <<colle_.ihadrz <<G4endl;
1879  G4cout <<"INTPT " <<droppt_.intpt <<G4endl;
1880  G4cout <<"PAULI " <<droppt_.lpauli <<G4endl;
1881  G4cout <<"FERMI " <<droppt_.fermp <<" "
1882  <<nucimp_.fermod <<G4endl;
1883  G4cout <<"CRONINPT " <<cronin_.mkcron <<" "
1884  <<cronin_.cronco <<G4endl;
1885  G4cout <<"SEADISTR " <<xseadi_.xseacu+0.95 <<" "
1886  <<xseadi_.unon <<" "
1887  <<xseadi_.unom <<" "
1888  <<xseadi_.unosea <<G4endl;
1889  G4cout <<"SEAQUARK " <<seaqxx_.seaqx <<" "
1890  <<seaqxx_.seaqxn <<G4endl;
1891  G4cout <<"SECINTER " <<secint_.isecin <<G4endl;
1892  G4cout <<"XCUTS " <<xseadi_.cvq <<" "
1893  <<xseadi_.cdq <<" "
1894  <<xseadi_.csea <<" "
1895  <<xseadi_.ssmima <<G4endl;
1896  }
1897 #endif
1898 
1899  bufues_.bnnvv = 0.001; // BNNVV=0.001
1900  bufues_.bnnss = 0.001; // BNNSS=0.001
1901  bufues_.bnnsv = 0.001; // BNNSV=0.001
1902  bufues_.bnnvs = 0.001; // BNNVS=0.001
1903  bufues_.bnncc = 0.001; // BNNCC=0.001
1904  bufues_.bnndv = 0.001; // BNNDV=0.001
1905  bufues_.bnnvd = 0.001; // BNNVD=0.001
1906  bufues_.bnnds = 0.001; // BNNDS=0.001
1907  bufues_.bnnsd = 0.001; // BNNSD=0.001
1908  bufues_.bnnhh = 0.001; // BNNHH=0.001
1909  bufues_.bnnzz = 0.001; // BNNZZ=0.001
1910  bufues_.bnndi = 0.001; // BNNDI=0.001
1911  bufues_.bnnzd = 0.001; // BNNZD=0.001
1912  bufues_.bnndz = 0.001; // BNNDZ=0.001
1913  bufues_.bptvv = 0.0; // BPTVV=0.
1914  bufues_.bptss = 0.0; // BPTSS=0.
1915  bufues_.bptsv = 0.0; // BPTSV=0.
1916  bufues_.bptvs = 0.0; // BPTVS=0.
1917  bufues_.bptcc = 0.0; // BPTCC=0.
1918  bufues_.bptdv = 0.0; // BPTDV=0.
1919  bufues_.bptvd = 0.0; // BPTVD=0.
1920  bufues_.bptds = 0.0; // BPTDS=0.
1921  bufues_.bptsd = 0.0; // BPTSD=0.
1922  bufues_.bpthh = 0.0; // BPTHH=0.
1923  bufues_.bptzz = 0.0; // BPTZZ=0.
1924  bufues_.bptdi = 0.0; // BPTDI=0.
1925  bufues_.bptzd = 0.0; // BPTZD=0.
1926  bufues_.bptdz = 0.0; // BPTDZ=0.
1927  bufues_.beevv = 0.0; // BEEVV=0.
1928  bufues_.beess = 0.0; // BEESS=0.
1929  bufues_.beesv = 0.0; // BEESV=0.
1930  bufues_.beevs = 0.0; // BEEVS=0.
1931  bufues_.beecc = 0.0; // BEECC=0.
1932  bufues_.beedv = 0.0; // BEEDV=0.
1933  bufues_.beevd = 0.0; // BEEVD=0.
1934  bufues_.beeds = 0.0; // BEEDS=0.
1935  bufues_.beesd = 0.0; // BEESD=0.
1936  bufues_.beehh = 0.0; // BEEHH=0.
1937  bufues_.beezz = 0.0; // BEEZZ=0.
1938  bufues_.beedi = 0.0; // BEEDI=0.
1939  bufues_.beezd = 0.0; // BEEZD=0.
1940  bufues_.beedz = 0.0; // BEEDZ=0.
1941  ncoucs_.bcouvv = 0.0; // BCOUVV=0.
1942  ncoucs_.bcouss = 0.0; // BCOUSS=0.
1943  ncoucs_.bcousv = 0.0; // BCOUSV=0.
1944  ncoucs_.bcouvs = 0.0; // BCOUVS=0.
1945  ncoucs_.bcouzz = 0.0; // BCOUZZ=0.
1946  ncoucs_.bcouhh = 0.0; // BCOUHH=0.
1947  ncoucs_.bcouds = 0.0; // BCOUDS=0.
1948  ncoucs_.bcousd = 0.0; // BCOUSD=0.
1949  ncoucs_.bcoudz = 0.0; // BCOUDZ=0.
1950  ncoucs_.bcouzd = 0.0; // BCOUZD=0.
1951  ncoucs_.bcoudi = 0.0; // BCOUDI=0.
1952  ncoucs_.bcoudv = 0.0; // BCOUDV=0.
1953  ncoucs_.bcouvd = 0.0; // BCOUVD=0.
1954  ncoucs_.bcoucc = 0.0; // BCOUCC=0.
1955 //
1956 //
1957 // Initialisation of
1958 // ANNVV, ANNSS ... ANNDZ
1959 // PTVV, PTSS ... PTDZ
1960 // EEVV, EESS ... EEDZ
1961 // ACOUVV, ACOUSS, ... ACOUCC
1962 // now all moved to ApplyYourself member function.
1963 //
1964 // droppt_.ipadis = LFALSE; // IPADIS = .FALSE.
1965 // droppt_.ihadvv = LFALSE; // IHADVV = .FALSE.
1966 // droppt_.ihadsv = LFALSE; // IHADSV = .FALSE.
1967 // droppt_.ihadvs = LFALSE; // IHADVS = .FALSE.
1968 
1969  nucc_.ijtarg = 1; // IJTARG=1
1970 
1971 //
1972 //
1973 // The following commented out since it seems to have more to do with histogram
1974 // generation.
1975 //
1976 // i = 1;
1977 // G4int idummy;
1978 // distr_ (&i,&nucc_.ijproj,&ppn,&idummy);
1979 // CALL DISTR( 1,IJPROJ,PPN,IDUMMY )
1980  if (pomtyp_.ipim == 2) {prblm2_ (&user2_.cmener);}
1981 //
1982 //
1983 // Initialise hard scattering & transverse momentum for soft scattering
1984 //
1985  i = 0;
1986  jtdtu_ (&i); // CALL JTDTU( 0 )
1987  i = 0;
1988  G4double pt;
1989  samppt_ (&i,&pt); // CALL SAMPPT(0,PT)
1990 }
1991 ////////////////////////////////////////////////////////////////////////////////
1992 //
1993 #endif
cmhico_
struct ccdpm25cmhico cmhico_
G4LorentzVector.hh
G4GlaubAADataSetHandler::IsGlauberDataSetAvailable
G4bool IsGlauberDataSetAvailable(const G4int AP, const G4int AT) const
G4Nucleus::GetA_asInt
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
sigma_
struct ccdpm25sigma sigma_
G4Nucleus::AtomicMass
G4double AtomicMass(const G4double A, const G4double Z) const
Definition: G4Nucleus.cc:240
G4GlaubAADataSetHandler::SetMaxGlauberDataSets
void SetMaxGlauberDataSets(const G4int n)
ccdpm25hkkevt::phkk
G4double phkk[nmxhkk][5]
Definition: G4DPMJET2_5Interface.hh:468
colle_
struct ccdpm25colle colle_
ncoucs_
struct ccdpm25ncoucs ncoucs_
stopAndKill
Definition: G4HadFinalState.hh:42
right
Definition: F04UserTrackInformation.hh:37
evappp_
struct ccdpm25evappp evappp_
pydat1_
struct ccdpm25pydat1 pydat1_
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4DynamicParticle::GetTotalEnergy
G4double GetTotalEnergy() const
CLHEP::Hep3Vector::dot
double dot(const Hep3Vector &) const
pomtyp_
struct ccdpm25pomtyp pomtyp_
nncms_
struct ccdpm25nncms nncms_
G4Fragment::GetA
G4double GetA() const
Definition: G4Fragment.hh:303
diqsum_
struct ccdpm25diqsum diqsum_
G4ExcitationHandler.hh
G4VNuclearDensity.hh
CORSIKA
Definition: G4DPMJET2_5Model.hh:85
parevt_
struct ccdpm25parevt parevt_
p
const char * p
Definition: xmltok.h:285
G4DPMJET2_5Model::SetVerboseFortranOutput
G4bool SetVerboseFortranOutput(const G4String filename)
dpar_
struct ccdpm25dpar dpar_
final_
struct ccdpm25final final_
paname_
struct ccdpm25paname paname_
G4DynamicParticle.hh
strufu_
struct ccdpm25strufu strufu_
G4DPMJET2_5Model::~G4DPMJET2_5Model
virtual ~G4DPMJET2_5Model()
g4dpmjet_open_nuclear_bin_
void g4dpmjet_open_nuclear_bin_(int *, int *, int *, char *)
dchant_
void dchant_()
bufues_
struct ccdpm25bufues bufues_
nstari_
struct ccdpm25nstari nstari_
secint_
struct ccdpm25secint secint_
dhadde_
void dhadde_()
G4int
int G4int
Definition: G4Types.hh:78
taufo_
struct ccdpm25taufo taufo_
rd2out_
double rd2out_(int *, int *)
user1_
struct ccdpm25user1 user1_
coulo_
struct ccdpm25coulo coulo_
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:159
g4dpmjet_close_nuclear_bin_
void g4dpmjet_close_nuclear_bin_(int *)
G4PhotonEvaporation.hh
G4IonTable.hh
ferfor_
struct ccdpm25ferfor ferfor_
xseadi_
struct ccdpm25xseadi xseadi_
cronin_
struct ccdpm25cronin cronin_
G4FermiMomentum.hh
G4HadFinalState::SetStatusChange
void SetStatusChange(G4HadFinalStateStatus aS)
Definition: G4HadFinalState.cc:81
pomtab_
struct ccdpm25pomtab pomtab_
G4ReactionProductVector
std::vector< G4ReactionProduct * > G4ReactionProductVector
Definition: G4ReactionProductVector.hh:43
dtumat_
struct ccdpm25dtumat dtumat_
G4NuclearShellModelDensity.hh
g4dpmjet_open_fort6_
void g4dpmjet_open_fort6_(int *, int *, char *)
seadiq_
struct ccdpm25seadiq seadiq_
G4ReactionProductVector.hh
G4HadronicInteraction::SetMinEnergy
void SetMinEnergy(G4double anEnergy)
Definition: G4HadronicInteraction.hh:90
DPM2_5
Definition: G4DPMJET2_5Model.hh:85
hboo_
struct ccdpm25hboo hboo_
g4dpmjet_close_fort6_
void g4dpmjet_close_fort6_()
vxsvd_
struct ccdpm25vxsvd vxsvd_
ptsamp_
struct ccdpm25ptsamp ptsamp_
G4PreCompoundModel.hh
G4ParticleTable::GetIonTable
G4IonTable * GetIonTable() const
Definition: G4ParticleTable.cc:712
CLHEP::HepLorentzVector::vect
Hep3Vector vect() const
G4DPMJET2_5Model::SetDefaultPreCompoundModel
void SetDefaultPreCompoundModel()
prblm2_
void prblm2_(double *)
G4cout
G4GLOB_DLL std::ostream G4cout
jtdtu_
void jtdtu_(int *)
G4DPMJET2_5Model::G4DPMJET2_5Model
G4DPMJET2_5Model()
G4HadronicInteraction::GetVerboseLevel
virtual G4int GetVerboseLevel() const
Definition: G4HadronicInteraction.hh:113
G4HadProjectile::GetDefinition
const G4ParticleDefinition * GetDefinition() const
Definition: G4HadProjectile.hh:81
G4DPMJET2_5Model::IsApplicable
virtual G4bool IsApplicable(const G4HadProjectile &theTrack, G4Nucleus &theTarget)
zentra_
struct ccdpm25zentra zentra_
ftnlogical
G4int ftnlogical
Definition: G4DPMJET2_5Interface.hh:85
G4ParticleTable.hh
casadi_
struct ccdpm25casadi casadi_
G4Fragment.hh
G4DPMJET2_5Model::SetVerboseLevel
void SetVerboseLevel(G4int)
Definition: G4DPMJET2_5Model.hh:238
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Poisson.hh
incini_
void incini_()
G4HadProjectile::GetKineticEnergy
G4double GetKineticEnergy() const
Definition: G4HadProjectile.hh:106
G4ExcitationHandler::SetFermiModel
void SetFermiModel(G4VFermiBreakUp *ptr)
Definition: G4ExcitationHandler.cc:476
xsecpt_
struct ccdpm25xsecpt xsecpt_
inpflg_
struct ccdpm25inpflg inpflg_
nucimp_
struct ccdpm25nucimp nucimp_
theTarget
G4ErrorTarget * theTarget
Definition: errprop.cc:59
G4DPMJET2_5Model::ApplyYourself
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &, G4Nucleus &)
G4ParticleDefinition.hh
csj1mi_
void csj1mi_(double *, double *)
Randomize.hh
sincha_
struct ccdpm25sincha sincha_
extevt_
struct ccdpm25extevt extevt_
hdjase_
struct ccdpm25hdjase hdjase_
dropjj_
struct ccdpm25dropjj dropjj_
n
const G4int n
Definition: G4UrQMD1_3Interface.hh:144
seasu3_
struct ccdpm25seasu3 seasu3_
nucc_
struct ccdpm25nucc nucc_
maxpro
const int maxpro
Definition: G4DPMJET2_5Interface.hh:67
G4DPMJET2_5Model::SetNoPreCompoundModel
void SetNoPreCompoundModel()
ptlarg_
struct ccdpm25ptlarg ptlarg_
G4HadProjectile::Get4Momentum
const G4LorentzVector & Get4Momentum() const
Definition: G4HadProjectile.hh:86
DPM3
Definition: G4DPMJET2_5Model.hh:85
G4ParticleMomentum.hh
ccdpm25extevt::idxres
G4int idxres[nmxhkk]
Definition: G4DPMJET2_5Interface.hh:461
droppt_
struct ccdpm25droppt droppt_
kglaub_
struct ccdpm25kglaub kglaub_
ccdpm25dpar::aam
G4double aam[210]
Definition: G4DPMJET2_5Interface.hh:455
g4dpmjet_initialise_block_data_
void g4dpmjet_initialise_block_data_()
parpt_
void parpt_(int *, double *, double *, int *, int *)
G4FermiBreakUp.hh
G4PhysicalConstants.hh
dchanh_
void dchanh_()
kkinc_
void kkinc_(double *, int *, int *, int *, int *, int *, int *, int *, int *, int *)
popcck_
struct ccdpm25popcck popcck_
G4HadFinalState::SetEnergyChange
void SetEnergyChange(G4double anEnergy)
Definition: G4HadFinalState.cc:42
stars_
struct ccdpm25stars stars_
defaul_
void defaul_(double *, double *)
isAlive
Definition: G4HadFinalState.hh:42
frbkcm_
struct ccdpm25frbkcm frbkcm_
collap_
struct ccdpm25collap collap_
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:161
G4ParticleTable::GetParticleTable
static G4ParticleTable * GetParticleTable()
Definition: G4ParticleTable.cc:95
lundin_
void lundin_()
diffra_
struct ccdpm25diffra diffra_
debug
#define debug
Definition: CCalPrimaryGeneratorAction.cc:46
ifragm_
struct ccdpm25ifragm ifragm_
G4ExcitationHandler::SetMaxAandZForFermiBreakUp
void SetMaxAandZForFermiBreakUp(G4int anA, G4int aZ)
Definition: G4ExcitationHandler.cc:503
CLHEP::Hep3Vector::unit
Hep3Vector unit() const
user2_
struct ccdpm25user2 user2_
G4GlaubAADataSetHandler::SetCurrentGlauberDataSet
G4bool SetCurrentGlauberDataSet(const G4int AP, const G4int AT, const G4double ppn=0.0)
G4Nucleus::GetZ_asInt
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
DEFAULT
Definition: G4DPMJET2_5Model.hh:85
collis_
struct ccdpm25collis collis_
G4DPMJET2_5Interface.hh
Definition of the G4DPMJET2_5Interface class.
G4ExcitationHandler::SetEvaporation
void SetEvaporation(G4VEvaporation *ptr)
Definition: G4ExcitationHandler.cc:456
diqrej_
struct ccdpm25diqrej diqrej_
G4HadronicInteraction::SetMaxEnergy
void SetMaxEnergy(const G4double anEnergy)
Definition: G4HadronicInteraction.hh:103
python.hepunit.perCent
float perCent
Definition: hepunit.py:239
G4endl
#define G4endl
Definition: G4ios.hh:61
hkkevt_
struct ccdpm25hkkevt hkkevt_
seaqxx_
struct ccdpm25seaqxx seaqxx_
G4DPMJET2_5Model::GetFatalEnergyCheckLevels
virtual const std::pair< G4double, G4double > GetFatalEnergyCheckLevels() const
nuclea_
struct ccdpm25nuclea nuclea_
ccdpm25hkkevt::isthkk
G4int isthkk[nmxhkk]
Definition: G4DPMJET2_5Interface.hh:466
G4GlaubAADataSetHandler.hh
Definition of the G4GlaubAADataSetHandler class.
pshow_
struct ccdpm25pshow pshow_
ddatar_
void ddatar_()
G4DPMJET2_5InitialisationType
G4DPMJET2_5InitialisationType
Definition: G4DPMJET2_5Model.hh:85
G4DPMJET2_5Model::SetDefaultDeexcitation
void SetDefaultDeexcitation()
G4GlaubAADataSetHandler::getInstance
static G4GlaubAADataSetHandler * getInstance()
G4double
double G4double
Definition: G4Types.hh:76
G4SystemOfUnits.hh
dprin_
struct ccdpm25dprin dprin_
edens_
struct ccdpm25edens edens_
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:163
G4ExcitationHandler::SetPhotonEvaporation
void SetPhotonEvaporation(G4VEvaporationChannel *ptr)
Definition: G4ExcitationHandler.cc:485
promu_
struct ccdpm25promu promu_
defaux_
void defaux_(double *, double *)
berttp_
void berttp_()
G4HadFinalState::SetMomentumChange
void SetMomentumChange(const G4ThreeVector &aV)
Definition: G4HadFinalState.cc:54
G4DPMJET2_5Model::SetNoDeexcitation
void SetNoDeexcitation()
gluspl_
struct ccdpm25gluspl gluspl_
diquax_
struct ccdpm25diquax diquax_
G4NuclearFermiDensity.hh
popcor_
struct ccdpm25popcor popcor_
fluctu_
struct ccdpm25fluctu fluctu_
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
hettp_
struct ccdpm25hettp hettp_
samppt_
void samppt_(int *, double *)
G4DPMJET2_5Model.hh
Definition of the G4DPMJET2_5Model class.
G4Evaporation.hh
hadthr_
struct ccdpm25hadthr hadthr_
G4HadProjectile::GetTotalEnergy
G4double GetTotalEnergy() const
Definition: G4HadProjectile.hh:96
G4LorentzVector
CLHEP::HepLorentzVector G4LorentzVector
Definition: G4LorentzVector.hh:35
recom_
struct ccdpm25recom recom_
G4HadTmpUtil.hh
Generated on Wed Apr 30 2014 15:54:42 for Geant4.10 by   doxygen 1.8.7