Ignore:
Timestamp:
Apr 6, 2009, 12:30:29 PM (15 years ago)
Author:
garnier
Message:

update processes

Location:
trunk/source/processes/hadronic/management
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • trunk/source/processes/hadronic/management/History

    r819 r962  
    1414     * Please list in reverse chronological order (last date on top)
    1515     ---------------------------------------------------------------
     16
     1725 February 2009 Vladimir Ivanchenko (hadr-man-V09-02-04)
     18---------------------------------------------------------
     19- G4HadronicProcessStore - added protection to the Clean method
     20
     2122 February 2009 Vladimir Ivanchenko (hadr-man-V09-02-03)
     22---------------------------------------------------------
     23
     2414 February 2009 Vladimir Ivanchenko (hadr-man-V09-02-02)
     25---------------------------------------------------------
     26- G4HadronicProcessStore - added Clean method and cleanup of model
     27                           and cross section stores
     28
     2906 February 2009 Vladimir Ivanchenko (hadr-man-V09-02-01)
     30---------------------------------------------------------
     31
     3224 January 2009 Vladimir Ivanchenko (hadr-man-V09-02-00)
     33-------------------------------------------------------
     34- G4HadronicProcessStore - added destruction of processes
     35
     3601 December 2008 Dennis Wright (hadr-man-V09-01-10)
     37---------------------------------------------------
     38- G4HadronicProcess - make MeanFreePath() public again.  In future
     39  make it protected again, but make a public method which calls it.
     40
     4121 November 2008 Dennis Wright (hadr-man-V09-01-09)
     42---------------------------------------------------
     43- G4HadronicProcess - remove again method SetDispatch for major release
     44
     4522 October 2008 Vladimir Ivanchenko (hadr-man-V09-01-08)
     46------------------------------------------------------
     47- G4HadronicProcessStore - use G4HadronicProcessType enumerator
     48- G4HadronicProcess - returned back obsolete method SetDispatch for minor release
     49
     5002 October 2008 Dennis Wright (hadr-man-V09-01-07)
     51--------------------------------------------------
     52- create new hadronic process subtype enum G4HadronicProcessType:
     53    enum G4HadronicProcessType
     54    {
     55      fHadronElastic =   111,
     56      fHadronInelastic = 121,
     57      fCapture =         131,
     58      fFission =         141,
     59      fHadronAtRest =    151,
     60      fChargeExchange =  161
     61    };
     62
     63- G4HadronicProcess.hh - add enum to identify process subtypes
     64
     65- G4HadronInelasticProcess.cc - change process sub-type from 12 to fHadronInelastic
     66
     6704 August 2008 Vladimir Ivanchenko (hadr-man-V09-01-06)
     68------------------------------------------------------
     69- G4HadronicProcessStore - improve cout
     70- G4HadronicProcess - cleanup: use method SampleZandA to select an isotope
     71                               do not use home-made NanCheck
     72                               do not check environment variables run-time
     73                               directly fill G4HadronicWhiteBoard
     74                               implement PostStepDoIt and DumpPhysicsTable
     75                               methods
     76- G4HadronInelasticProcess - use methods of the base class
     77                         
     7804 August 2008 Vladimir Ivanchenko (hadr-man-V09-01-05)
     79------------------------------------------------------
     80- G4HadronicProcessStore - use sub-types to access cross sections,
     81                           add proceses following only G4VProcess interface   
     82
     8308 July 2008 Dennis Wright (hadr-man-V09-01-04)
     84-----------------------------------------------
     85- set process sub-type to 12 for G4HadronInelasticProcess
     86
     8709 June 2008 Dennis Wright (hadr-man-V09-01-03)
     88-----------------------------------------------
     89- G4HadronicProcess.cc - turn off error in case of fStopButAlive, but leave
     90  it in place for fStopAndKill, fKillTrackAndSecondaries and fPostponeToNextEvent
     91
     9205 June 2008 Vladimir Ivanchenko (hadr-man-V09-01-02)
     93----------------------------------------------------
     94- G4HadronicProcessStore - comment out destructor
     95
     9619 May 2008 Vladimir Ivanchenko (hadr-man-V09-01-01)
     97----------------------------------------------------
     98- G4HadronicProcessStore - new singleton to keep pointers to all
     99                      hadronic processes
     100- G4HadronicProcess - added PreparPhysicsTable and BuildPhysicsTable
     101                      methods, added registration in G4HadronicProcessStore
     102
     103
     10419 May 2008 Vladimir Ivanchenko (hadr-man-V09-01-00)
     105----------------------------------------------------
     106- G4HadronicProcess - cleanup of the header (add comments, move
     107                      methods for isotope production to the src, make
     108                      GetMeanFreePath protected, remove duplication
     109                      of PostStepDoIt), add default implementation of
     110                      GetMicroscopicCrossSection, comment out check
     111                      IfApplicable in MeanFreePath
    16112
    1711315 October 2007 Dennis Wright (hadr-man-V09-00-00)

  • trunk/source/processes/hadronic/management/include/G4EnergyRangeManager.hh

    r819 r962  
    2626//
    2727// $Id: G4EnergyRangeManager.hh,v 1.9 2006/06/29 19:58:05 gunter Exp $
    28 // GEANT4 tag $Name: $
     28// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2929//
    3030 // Hadronic Process: Energy Range Manager

  • trunk/source/processes/hadronic/management/include/G4HadLeadBias.hh

    r819 r962  
    2424// ********************************************************************
    2525//
    26 // GEANT4 tag $Name: $
     26// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2727//
    2828// --------------------------------------------------------------------

  • trunk/source/processes/hadronic/management/include/G4HadronInelasticProcess.hh

    r819 r962  
    4343
    4444#include "G4HadronicProcess.hh"
    45 //#include "G4LPhysicsFreeVector.hh"
    46 #include "G4HadronCrossSections.hh"
    47 #include "G4CrossSectionDataStore.hh"
    48 #include "G4HadronInelasticDataSet.hh"
    49 #include "G4ParticleChange.hh"
    50  
    5145
    52  class G4HadronInelasticProcess : public G4HadronicProcess
    53  {
    54  public:
     46class G4ParticleDefinition;
     47
     48class G4HadronInelasticProcess : public G4HadronicProcess
     49{
     50public:
    5551   
    56     G4HadronInelasticProcess(
    57      const G4String &processName,
    58      G4ParticleDefinition *aParticle );
     52  G4HadronInelasticProcess(const G4String &processName,
     53                           G4ParticleDefinition *aParticle );
    5954   
    60     virtual ~G4HadronInelasticProcess();
     55  virtual ~G4HadronInelasticProcess();
    6156       
    62     void BuildThePhysicsTable();
    63    
    64     G4bool IsApplicable(const G4ParticleDefinition& aP);
     57  virtual G4bool IsApplicable(const G4ParticleDefinition& aParticleType);
    6558
    66     G4VParticleChange *PostStepDoIt(const G4Track &aTrack, const G4Step &aStep);
     59private:
    6760
    68   private:   
     61  G4ParticleDefinition* theParticle;
    6962
    70     virtual G4double GetMicroscopicCrossSection( const G4DynamicParticle *aParticle,
    71                                                  const G4Element *anElement,
    72                                                  G4double aTemp );
    73    
    74  protected:
    75 
    76     G4ParticleDefinition *theParticle;
    77     G4ParticleChange theParticleChange;
    78  };
     63};
    7964 
    8065#endif

  • trunk/source/processes/hadronic/management/include/G4HadronicProcess.hh

    r819 r962  
    2626//
    2727//
    28  // This is the top level Hadronic Process class
    29  // The inelastic, elastic, capture, and fission processes
    30  // should derive from this class
    31  // This is an abstract base class, since the pure virtual function
    32  // PostStepDoIt has not been defined yet.
    33  // Note:  there is no .cc file
    34  //
    35  // original by H.P.Wellisch
    36  // J.L. Chuma, TRIUMF, 10-Mar-1997
    37  // Last modified: 04-Apr-1997
     28// This is the top level Hadronic Process class
     29// The inelastic, elastic, capture, and fission processes
     30// should derive from this class
     31//
     32// original by H.P.Wellisch
     33// J.L. Chuma, TRIUMF, 10-Mar-1997
     34// Last modified: 04-Apr-1997
     35// 19-May-2008 V.Ivanchenko cleanup and added comments
    3836 
    3937#ifndef G4HadronicProcess_h
    … …  
    5048#include "G4VCrossSectionDataSet.hh"
    5149#include "G4VLeadingParticleBiasing.hh"
    52 #include "G4Delete.hh"
     50//#include "G4Delete.hh"
    5351#include "G4CrossSectionDataStore.hh"
    54 #include "G4HadronicException.hh"
    55 
     52#include "G4HadronicProcessType.hh"
    5653
    5754class G4Track;
    … …  
    6057class G4ParticleChange;
    6158
    62  class G4HadronicProcess : public G4VDiscreteProcess
    63  {
    64  public:
    65    
    66     G4HadronicProcess( const G4String &processName = "Hadronic",
    67                        G4ProcessType   aType = fHadronic );   
    68     virtual ~G4HadronicProcess();
    69 
    70     void RegisterMe( G4HadronicInteraction *a );
    71 
    72     void AddDataSet(G4VCrossSectionDataSet * aDataSet)
    73     {
    74        theCrossSectionDataStore->AddDataSet(aDataSet);
    75     }
     59class G4HadronicProcess : public G4VDiscreteProcess
     60{
     61public:
     62   
     63  G4HadronicProcess( const G4String &processName = "Hadronic",
     64                     G4ProcessType   aType = fHadronic );   
     65
     66  virtual ~G4HadronicProcess();
     67
     68  // register generator of secondaries
     69  void RegisterMe( G4HadronicInteraction *a );
     70
     71  // get cross section per element
     72  virtual
     73  G4double GetMicroscopicCrossSection(const G4DynamicParticle *aParticle,
     74                                      const G4Element *anElement,
     75                                      G4double aTemp );
     76
     77  // generic PostStepDoIt recommended for all derived classes
     78  virtual G4VParticleChange* PostStepDoIt(const G4Track& aTrack,
     79                                          const G4Step& aStep);
     80
     81  // initialisation of physics tables and G4HadronicProcessStore
     82  virtual void PreparePhysicsTable(const G4ParticleDefinition&);
     83
     84  // build physics tables and print out the configuration of the process
     85  virtual void BuildPhysicsTable(const G4ParticleDefinition&);
     86
     87  // dump physics tables
     88  inline void DumpPhysicsTable(const G4ParticleDefinition& p)
     89  { theCrossSectionDataStore->DumpPhysicsTable(p); }
     90
     91  // add cross section data set
     92  inline void AddDataSet(G4VCrossSectionDataSet * aDataSet)
     93  { theCrossSectionDataStore->AddDataSet(aDataSet);}
     94
     95  // access to the manager
     96  inline G4EnergyRangeManager *GetManagerPointer()
     97  { return &theEnergyRangeManager; }
     98         
     99  // get inverse cross section per volume
     100  G4double GetMeanFreePath(const G4Track &aTrack, G4double,
     101                           G4ForceCondition *);
     102
     103protected:   
     104
     105  // reset number of interaction length and save 
     106  virtual void ResetNumberOfInteractionLengthLeft()
     107  { G4VProcess::ResetNumberOfInteractionLengthLeft();
     108    theInitialNumberOfInteractionLength =
     109      G4VProcess::theNumberOfInteractionLengthLeft;
     110  }
     111
     112  // generic method to choose secondary generator
     113  // recommended for all derived classes
     114  inline G4HadronicInteraction *ChooseHadronicInteraction(
     115      G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement )
     116  { return theEnergyRangeManager.GetHadronicInteraction(kineticEnergy,
     117                                                        aMaterial,anElement);
     118  }
     119
     120public:
     121
     122  // Methods for isotope production   
     123  static void EnableIsotopeProductionGlobally();
     124  static void DisableIsotopeProductionGlobally();
     125   
     126  void EnableIsotopeCounting()  {isoIsOnAnyway = 1;}
     127  void DisableIsotopeCounting() {isoIsOnAnyway = -1;}
     128   
     129  void RegisterIsotopeProductionModel(G4VIsotopeProduction * aModel)
     130  { theProductionModels.push_back(aModel); }
     131
     132  static G4IsoParticleChange * GetIsotopeProductionInfo();
     133
     134  void BiasCrossSectionByFactor(G4double aScale);
     135   
     136protected:
     137           
     138  // obsolete method will be removed
     139  inline const G4EnergyRangeManager &GetEnergyRangeManager() const
     140  { return theEnergyRangeManager; }
     141   
     142  // obsolete method will be removed
     143  inline void SetEnergyRangeManager( const G4EnergyRangeManager &value )
     144  { theEnergyRangeManager = value; }
     145
     146  // access to the chosen generator
     147  inline G4HadronicInteraction *GetHadronicInteraction()
     148  { return theInteraction; }
     149   
     150  // access to the cross section data store
     151  inline G4CrossSectionDataStore* GetCrossSectionDataStore()
     152  { return theCrossSectionDataStore; }
     153   
     154  // access to the cross section data set
     155  inline G4double GetLastCrossSection()
     156  { return theLastCrossSection; }
     157
     158private:
     159   
     160  void FillTotalResult(G4HadFinalState * aR, const G4Track & aT);
     161
     162  G4HadFinalState * DoIsotopeCounting(G4HadFinalState * aResult,
     163                                      const G4Track & aTrack,
     164                                      const G4Nucleus & aNucleus);
     165                                         
     166  G4IsoResult * ExtractResidualNucleus(const G4Track & aTrack,
     167                                       const G4Nucleus & aNucleus,
     168                                       G4HadFinalState * aResult);
     169
     170  inline G4double GetTotalNumberOfInteractionLengthTraversed()
     171  { return theInitialNumberOfInteractionLength
     172      -G4VProcess::theNumberOfInteractionLengthLeft;
     173  }
    76174       
    77     virtual G4VParticleChange *PostStepDoIt( const G4Track &aTrack,
    78                                             const G4Step &aStep ) = 0;
    79        
    80     virtual
    81     G4double GetMicroscopicCrossSection(const G4DynamicParticle *aParticle,
    82                                         const G4Element *anElement,
    83                                         G4double aTemp ) = 0;
    84    
    85     G4double GetMeanFreePath(const G4Track &aTrack, G4double,
    86                              G4ForceCondition *);
    87 
    88     // Set methods for isotope production
    89    
    90     static void EnableIsotopeProductionGlobally();
    91     static void DisableIsotopeProductionGlobally();
    92    
    93     void EnableIsotopeCounting()  {isoIsOnAnyway = 1;}
    94     void DisableIsotopeCounting() {isoIsOnAnyway = -1;}
    95    
    96     void RegisterIsotopeProductionModel(G4VIsotopeProduction * aModel)
    97     { theProductionModels.push_back(aModel); }
    98 
    99     static G4IsoParticleChange * GetIsotopeProductionInfo()
    100     {
    101       G4IsoParticleChange * anIsoResult = theIsoResult;
    102       if(theIsoResult) theOldIsoResult = theIsoResult;
    103       theIsoResult = 0;
    104       return anIsoResult;
    105     }
    106 
    107     void BiasCrossSectionByFactor(G4double aScale)
    108     {
    109       xBiasOn = true;
    110       aScaleFactor = aScale;
    111       G4String it = GetProcessName();
    112       if( (it != "PhotonInelastic") &&
    113           (it != "ElectroNuclear") &&
    114           (it != "PositronNuclear") )
    115       {
    116         G4Exception("G4HadronicProcess", "007", FatalException,
    117                     "Cross-section biasing available only for gamma and electro nuclear reactions.");
    118       }
    119       if(aScale<100)
    120       {
    121         G4Exception("G4HadronicProcess", "001", JustWarning,
    122                     "Cross-section bias readjusted to be above safe limit. New value is 100");       
    123         aScaleFactor = 100.;
    124       }
    125     }
    126    
    127  protected:
    128    
    129     virtual void ResetNumberOfInteractionLengthLeft()
    130     {
    131       G4VProcess::theNumberOfInteractionLengthLeft = 
    132                                          -std::log( G4UniformRand() );
    133       theInitialNumberOfInteractionLength =
    134                                  G4VProcess::theNumberOfInteractionLengthLeft;
    135     }
    136 
    137     G4VParticleChange *GeneralPostStepDoIt( const G4Track &aTrack,
    138                                            const G4Step &aStep );
    139    
    140     void SetDispatch( G4HadronicProcess *value )
    141     { dispatch=value; }
    142    
    143     G4Element* ChooseAandZ(const G4DynamicParticle *aParticle,
    144                            const G4Material *aMaterial);
    145 
    146     inline const G4EnergyRangeManager &GetEnergyRangeManager() const
    147     { return theEnergyRangeManager; }
    148    
    149     inline void SetEnergyRangeManager( const G4EnergyRangeManager &value )
    150     { theEnergyRangeManager = value; }
    151 
    152     inline G4HadronicInteraction *ChooseHadronicInteraction(
    153      G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement )
    154     {
    155       G4EnergyRangeManager* ERMan = GetManagerPointer();
    156       if(!ERMan->GetHadronicInteractionCounter())
    157         G4cout<< "*G4HadronicProcess::ChooseHadronicInteraction: process = "
    158               << GetProcessName() << ", nM="
    159               << ERMan->GetHadronicInteractionCounter() << G4endl;
    160       return ERMan->GetHadronicInteraction(kineticEnergy,aMaterial,anElement);
    161     }
    162 
    163     inline G4HadronicInteraction *GetHadronicInteraction()
    164     { return theInteraction; }
    165    
    166  public:
    167     inline G4EnergyRangeManager *GetManagerPointer()
    168     { return &theEnergyRangeManager; }
    169  protected:
    170 
    171     G4CrossSectionDataStore* GetCrossSectionDataStore()
    172     {
    173        return theCrossSectionDataStore;
    174     }
    175    
    176     G4double GetLastCrossSection() {return theLastCrossSection;}
    177  private:
    178    
    179     G4HadFinalState * DoIsotopeCounting(G4HadFinalState * aResult,
    180                                           const G4Track & aTrack,
    181                                           const G4Nucleus & aNucleus);
    182                                          
    183     G4IsoResult * ExtractResidualNucleus(const G4Track & aTrack,
    184                                          const G4Nucleus & aNucleus,
    185                                          G4HadFinalState * aResult);
    186 
    187     G4double GetTotalNumberOfInteractionLengthTraversed()
    188     {
    189       return theInitialNumberOfInteractionLength
    190             -G4VProcess::theNumberOfInteractionLengthLeft;
    191     }
    192    
    193 
    194     void FillTotalResult(G4HadFinalState * aR, const G4Track & aT);
    195    
    196     void SetCrossSectionDataStore(G4CrossSectionDataStore* aDataStore)
    197     {
    198        theCrossSectionDataStore = aDataStore;
    199     }
    200    
    201     G4double XBiasSurvivalProbability();
    202     G4double XBiasSecondaryWeight();
    203    
    204  private:
    205    
    206     G4EnergyRangeManager theEnergyRangeManager;
    207    
    208     G4HadronicInteraction *theInteraction;
    209 
    210     G4CrossSectionDataStore* theCrossSectionDataStore;
    211  
    212     G4Nucleus targetNucleus;
    213    
    214     G4HadronicProcess *dispatch;
    215 
    216 // swiches for isotope production
    217    
    218     static G4bool isoIsEnabled; // true or false; local swich overrides
    219     G4int isoIsOnAnyway; // true(1), false(-1) or default(0)
    220    
    221     G4IsoParticleChange theIsoPC;
    222     std::vector<G4VIsotopeProduction *> theProductionModels;
    223    
    224     std::vector<G4VLeadingParticleBiasing *> theBias;
    225 
    226     static G4IsoParticleChange* theIsoResult;
    227     static G4IsoParticleChange* theOldIsoResult;
    228    
    229     G4ParticleChange* theTotalResult;
    230    
    231     G4double theInitialNumberOfInteractionLength;   
    232 
    233     G4double aScaleFactor;
    234     G4bool xBiasOn;
    235     G4double theLastCrossSection;
    236 
    237     G4int ModelingState;
    238  };
     175  //  inline void SetCrossSectionDataStore(G4CrossSectionDataStore* aDataStore)
     176  //  { theCrossSectionDataStore = aDataStore; }
     177   
     178  G4double XBiasSurvivalProbability();
     179  G4double XBiasSecondaryWeight();
     180   
     181private:
     182   
     183  G4EnergyRangeManager theEnergyRangeManager;
     184   
     185  G4HadronicInteraction *theInteraction;
     186
     187  G4CrossSectionDataStore* theCrossSectionDataStore;
     188 
     189  G4Nucleus targetNucleus;
     190   
     191  G4HadronicProcess *dispatch;
     192
     193  bool G4HadronicProcess_debug_flag;
     194
     195  // swiches for isotope production   
     196  static G4bool isoIsEnabled; // true or false; local swich overrides
     197  G4int isoIsOnAnyway; // true(1), false(-1) or default(0)
     198   
     199  G4IsoParticleChange theIsoPC;
     200  std::vector<G4VIsotopeProduction *> theProductionModels;
     201   
     202  std::vector<G4VLeadingParticleBiasing *> theBias;
     203
     204  static G4IsoParticleChange* theIsoResult;
     205  static G4IsoParticleChange* theOldIsoResult;
     206   
     207  G4ParticleChange* theTotalResult;
     208   
     209  G4double theInitialNumberOfInteractionLength;   
     210
     211  G4double aScaleFactor;
     212  G4bool xBiasOn;
     213  G4double theLastCrossSection;
     214
     215  G4int ModelingState;
     216};
    239217 
    240218#endif

  • trunk/source/processes/hadronic/management/include/G4VLeadingParticleBiasing.hh

    r819 r962  
    2525//
    2626// $Id: G4VLeadingParticleBiasing.hh,v 1.6 2006/06/29 19:58:19 gunter Exp $
    27 // GEANT4 tag $Name: $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// --------------------------------------------------------------------

  • trunk/source/processes/hadronic/management/src/G4EnergyRangeManager.cc

    r819 r962  
    2626//
    2727// $Id: G4EnergyRangeManager.cc,v 1.15 2006/06/29 19:58:21 gunter Exp $
    28 // GEANT4 tag $Name: $
     28// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2929//
    3030 // Hadronic Process: Energy Range Manager

  • trunk/source/processes/hadronic/management/src/G4HadronInelasticProcess.cc

    r819 r962  
    2525//
    2626//
    27 //
    28  // Hadronic Inelastic Process Class
    29  // J.L. Chuma, TRIUMF, 24-Mar-1997
    30  // Last modified: 27-Mar-1997
    31  // J.P. Wellisch: Bug hunting, 23-Apr-97
    32  // Modified by J.L.Chuma 8-Jul-97 to eliminate possible division by zero for sigma
     27// Hadronic Inelastic Process Class
     28// J.L. Chuma, TRIUMF, 24-Mar-1997
     29// Last modified: 27-Mar-1997
     30// J.P. Wellisch: Bug hunting, 23-Apr-97
     31// Modified by J.L.Chuma 8-Jul-97 to eliminate possible division by zero for sigma
    3332//
    3433// 14-APR-98 F.W.Jones: variant G4HadronInelastic process for
    … …  
    3635//
    3736// 17-JUN-98 F.W.Jones: removed extraneous code causing core dump.
     37// 01-SEP-2008 V.Ivanchenko: use methods from the base class
    3838//
    3939 
    4040#include "G4HadronInelasticProcess.hh"
     41#include "G4HadronInelasticDataSet.hh"
    4142#include "G4GenericIon.hh"
    42 #include "G4ProcessManager.hh"
    43 #include "G4ProcessVector.hh"
    44 #include "G4HadronicException.hh"
     43#include "G4ParticleDefinition.hh"
    4544 
    46  void G4HadronInelasticProcess::BuildThePhysicsTable()
    47   {
    48     if (!G4HadronicProcess::GetCrossSectionDataStore()) {
    49       return;
    50     }
    51     G4HadronicProcess::GetCrossSectionDataStore()->BuildPhysicsTable(*theParticle);
    52   }
    53  
    54  G4HadronInelasticProcess::G4HadronInelasticProcess(
    55   const G4String &processName,
    56   G4ParticleDefinition *aParticle ) :
    57    G4HadronicProcess( processName )
    58  {
    59    G4HadronicProcess::AddDataSet(new G4HadronInelasticDataSet);
    60    theParticle = aParticle;
    61  }
     45G4HadronInelasticProcess::G4HadronInelasticProcess(const G4String& processName,
     46                                                   G4ParticleDefinition* aParticle):
     47  G4HadronicProcess(processName)
     48{
     49  SetProcessSubType(fHadronInelastic);
     50  AddDataSet(new G4HadronInelasticDataSet());
     51  theParticle = aParticle;
     52}
    6253
    63  G4HadronInelasticProcess::~G4HadronInelasticProcess() { }
     54G4HadronInelasticProcess::~G4HadronInelasticProcess()
     55{}
    6456
    65  G4VParticleChange *G4HadronInelasticProcess::
    66  PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
    67  {
    68    if(0==GetLastCrossSection()&&!getenv("DebugNeutronHP"))
    69    {
    70      G4cerr << "G4HadronInelasticProcess: called for final state, while cross-section was zero"<<G4endl;
    71      G4cerr << "                          Returning empty particle change...."<<G4endl;
    72      G4double dummy=0;
    73      G4ForceCondition condition;
    74      G4double it = GetMeanFreePath(aTrack, dummy, &condition);
    75      G4cerr << "                          current MeanFreePath is "<<it<<G4endl;
    76      theParticleChange.Initialize(aTrack);
    77      return &theParticleChange;
    78    }
    79    SetDispatch( this );
    80    return G4HadronicProcess::GeneralPostStepDoIt( aTrack, aStep );
    81  }
    82 
    83  G4bool G4HadronInelasticProcess::
    84  IsApplicable(const G4ParticleDefinition& aP)
    85  {
    86     return  theParticle == &aP || theParticle == G4GenericIon::GenericIon();
    87  }
    88 
    89  G4double G4HadronInelasticProcess::GetMicroscopicCrossSection(
    90   const G4DynamicParticle *aParticle,
    91   const G4Element *anElement,
    92   G4double aTemp)
    93   {
    94     // returns the microscopic cross section in GEANT4 internal units
    95    
    96    if (!G4HadronicProcess::GetCrossSectionDataStore())
    97    {
    98       throw G4HadronicException(__FILE__, __LINE__,
    99       "G4HadronInelasticProcess::GetMicroscopicCrossSection: "
    100                   "no CrossSectionDataStore");
    101       return DBL_MIN;
    102    }
    103    return G4HadronicProcess::GetCrossSectionDataStore()->GetCrossSection(aParticle, anElement, aTemp);
    104 
    105   }
    106  
    107  /* end of file */
     57G4bool G4HadronInelasticProcess::IsApplicable(const G4ParticleDefinition& aP)
     58{
     59  return  theParticle == &aP || theParticle == G4GenericIon::GenericIon();
     60}

  • trunk/source/processes/hadronic/management/src/G4HadronicProcess.cc

    r819 r962  
    2828#include "G4Types.hh"
    2929
    30 #include <fstream>
    31 #include <sstream>
    32 #include <stdlib.h>
     30//#include <fstream>
     31//#include <sstream>
     32//#include <stdlib.h>
    3333#include "G4HadronicProcess.hh"
    34 // #include "G4EffectiveCharge.hh"
     34
    3535#include "G4HadProjectile.hh"
    3636#include "G4ElementVector.hh"
    … …  
    4949#include "G4HadronicException.hh"
    5050#include "G4HadReentrentException.hh"
    51 #include "G4HadronicInteractionWrapper.hh"
     51
     52#include "G4HadronicWhiteBoard.hh"
    5253
    5354#include "G4HadSignalHandler.hh"
     55
     56#include "G4HadronicProcessStore.hh"
    5457
    5558#include <typeinfo>
    … …  
    7275void G4HadronicProcess::
    7376DisableIsotopeProductionGlobally() {isoIsEnabled = false;}
     77
     78//////////////////////////////////////////////////////////////////
    7479
    7580G4HadronicProcess::G4HadronicProcess( const G4String &processName,
    … …  
    8186  theTotalResult = new G4ParticleChange();
    8287  theCrossSectionDataStore = new G4CrossSectionDataStore();
     88  G4HadronicProcessStore::Instance()->Register(this);
    8389  aScaleFactor = 1;
    8490  xBiasOn = false;
     91  G4HadronicProcess_debug_flag = false;
    8592  if(getenv("SwitchLeadBiasOn")) theBias.push_back(new G4HadLeadBias());
    8693}
    … …  
    8895G4HadronicProcess::~G4HadronicProcess()
    8996{
     97  G4HadronicProcessStore::Instance()->DeRegister(this);
    9098  delete theTotalResult;
    9199
    … …  
    94102  std::for_each(theBias.begin(), theBias.end(), G4Delete());
    95103 
    96   delete theOldIsoResult; delete theIsoResult;
     104  delete theOldIsoResult;
     105  delete theIsoResult;
    97106  delete theCrossSectionDataStore;
    98107}
    … …  
    107116    "Could not register G4HadronicInteraction");
    108117  }
     118  G4HadronicProcessStore::Instance()->RegisterInteraction(this, a); 
     119}
     120
     121void G4HadronicProcess::PreparePhysicsTable(const G4ParticleDefinition& p)
     122{
     123  if(getenv("G4HadronicProcess_debug")) G4HadronicProcess_debug_flag = true;
     124  G4HadronicProcessStore::Instance()->RegisterParticle(this, &p);
     125}
     126
     127void G4HadronicProcess::BuildPhysicsTable(const G4ParticleDefinition& p)
     128{
     129  theCrossSectionDataStore->BuildPhysicsTable(p);
     130  G4HadronicProcessStore::Instance()->PrintInfo(&p);
    109131}
    110132
    … …  
    112134GetMeanFreePath(const G4Track &aTrack, G4double, G4ForceCondition *)
    113135{
    114   G4double sigma = 0.0;
    115136  try
    116137  {
    117     const G4DynamicParticle *aParticle = aTrack.GetDynamicParticle();
    118     if( !IsApplicable(*aParticle->GetDefinition()))
    119     {
    120       G4cout << "Unrecoverable error: "<<G4endl;
    121       G4ProcessManager * it = aParticle->GetDefinition()->GetProcessManager();
    122       G4ProcessVector * itv = it->GetProcessList();
    123       G4cout <<aParticle->GetDefinition()->GetParticleName()<<
    124       " has the following processes:"<<G4endl;
    125       for(G4int i=0; i<itv->size(); i++)
    126       {
    127         G4cout <<"  "<<(*itv)[i]->GetProcessName()<<G4endl;             
    128       }
    129       G4cout << "for kinetic energy "<<aParticle->GetKineticEnergy()<<G4endl;
    130       G4cout << "and material "<<aTrack.GetMaterial()->GetName()<<G4endl;
    131       G4Exception("G4HadronicProcess", "007", FatalException,
    132       std::string(this->GetProcessName()+
    133       " was called for "+
    134       aParticle->GetDefinition()->GetParticleName()).c_str() );
    135     }
    136     G4Material *aMaterial = aTrack.GetMaterial();
    137     ModelingState = 1;
    138    
    139     sigma = theCrossSectionDataStore->GetCrossSection(aParticle, aMaterial);
    140 
    141     sigma *= aScaleFactor;
    142     theLastCrossSection = sigma;
     138    ModelingState = 1;   
     139    theLastCrossSection = aScaleFactor*
     140      theCrossSectionDataStore->GetCrossSection(aTrack.GetDynamicParticle(),
     141                                                aTrack.GetMaterial());
    143142  }
    144143  catch(G4HadronicException aR)
    … …  
    148147    "G4HadronicProcess::GetMeanFreePath failed");
    149148  }
    150   if( sigma > 0.0 )
    151     return 1.0/sigma;
    152   else
    153     return DBL_MAX;
    154 }
    155 
    156 
    157 G4Element* G4HadronicProcess::ChooseAandZ(
    158 const G4DynamicParticle *aParticle, const G4Material *aMaterial )
    159 {
    160   std::pair<G4double, G4double> ZA =
    161       theCrossSectionDataStore->SelectRandomIsotope(aParticle, aMaterial);
    162   G4double ZZ = ZA.first;
    163   G4double AA = ZA.second;
    164 
    165   targetNucleus.SetParameters(AA, ZZ);
    166 
    167   const G4int numberOfElements = aMaterial->GetNumberOfElements();
    168   const G4ElementVector* theElementVector = aMaterial->GetElementVector();
    169   G4Element* chosen = 0;
    170   for (G4int i = 0; i < numberOfElements; i++) {
    171     chosen = (*theElementVector)[i];
    172     if (chosen->GetZ() == ZZ) break;
    173   }
    174   return chosen;
    175 }
    176 
    177 
    178 struct G4Nancheck{ bool operator()(G4double aV){return (!(aV<1))&&(!(aV>-1));}};
    179 
    180 G4VParticleChange *G4HadronicProcess::GeneralPostStepDoIt(
    181 const G4Track &aTrack, const G4Step &)
     149  G4double res = DBL_MAX;
     150  if( theLastCrossSection > 0.0 ) res = 1.0/theLastCrossSection;
     151  return res;
     152}
     153
     154G4double G4HadronicProcess::
     155GetMicroscopicCrossSection(const G4DynamicParticle *aParticle,
     156                           const G4Element *anElement,
     157                           G4double aTemp )
     158{
     159  return
     160    theCrossSectionDataStore->GetCrossSection(aParticle, anElement, aTemp);
     161}
     162
     163G4VParticleChange *G4HadronicProcess::PostStepDoIt(
     164                   const G4Track &aTrack, const G4Step &)
    182165{
    183166  // Debugging stuff
    184 
    185   bool G4HadronicProcess_debug_flag = false;
    186   if(getenv("G4HadronicProcess_debug")) G4HadronicProcess_debug_flag = true;
    187167  if(G4HadronicProcess_debug_flag)
    188168         std::cout << "@@@@ hadronic process start "<< std::endl;
    189169  // G4cout << theNumberOfInteractionLengthLeft<<G4endl;
    190   #ifndef G4HadSignalHandler_off
     170#ifndef G4HadSignalHandler_off
    191171  G4HadSignalHandler aHandler(G4HadronicProcess_local::G4HadronicProcessHandler_1);
    192   #endif
    193 
    194   if(aTrack.GetTrackStatus() != fAlive && aTrack.GetTrackStatus() != fSuspend)
    195   {
    196     G4cerr << "G4HadronicProcess: track in unusable state - "
    197     <<aTrack.GetTrackStatus()<<G4endl;
    198     G4cerr << "G4HadronicProcess: returning unchanged track "<<G4endl;
    199     G4Exception("G4HadronicProcess", "001", JustWarning, "bailing out");
     172#endif
     173
     174  if(aTrack.GetTrackStatus() != fAlive && aTrack.GetTrackStatus() != fSuspend) {
     175    if (aTrack.GetTrackStatus() == fStopAndKill ||
     176        aTrack.GetTrackStatus() == fKillTrackAndSecondaries ||
     177        aTrack.GetTrackStatus() == fPostponeToNextEvent) {
     178      G4cerr << "G4HadronicProcess: track in unusable state - "
     179             << aTrack.GetTrackStatus() << G4endl;
     180      G4cerr << "G4HadronicProcess: returning unchanged track " << G4endl;
     181      G4Exception("G4HadronicProcess", "001", JustWarning, "bailing out");
     182    }
     183    // No warning for fStopButAlive which is a legal status here
    200184    theTotalResult->Clear();
    201185    theTotalResult->Initialize(aTrack);
    … …  
    203187  }
    204188
    205   const G4DynamicParticle *aParticle = aTrack.GetDynamicParticle();
    206   G4Material *aMaterial = aTrack.GetMaterial();
     189  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
     190  G4Material* aMaterial = aTrack.GetMaterial();
    207191  G4double originalEnergy = aParticle->GetKineticEnergy();
    208192  G4double kineticEnergy = originalEnergy;
    209193
    210   // More debugging
    211 
     194  /*
     195  // It is not needed with standard NaN check
     196  // More debugging 
    212197  G4Nancheck go_wild;
    213198  if(go_wild(originalEnergy) ||
    … …  
    223208    return theTotalResult;
    224209  }
     210  */
    225211
    226212  // Get kinetic energy per nucleon for ions
    227 
    228213  if(aParticle->GetDefinition()->GetBaryonNumber() > 1.5)
    229214          kineticEnergy/=aParticle->GetDefinition()->GetBaryonNumber();
    … …  
    232217  try
    233218  {
    234     anElement = ChooseAandZ( aParticle, aMaterial );
     219    //    anElement = ChooseAandZ( aParticle, aMaterial );
     220    anElement = theCrossSectionDataStore->SampleZandA(aParticle,
     221                                                      aMaterial,
     222                                                      targetNucleus);
    235223  }
    236224  catch(G4HadronicException & aR)
    … …  
    243231    <<aParticle->GetDefinition()->GetParticleName()<<G4endl;
    244232    G4Exception("G4HadronicProcess", "007", FatalException,
    245     "GeneralPostStepDoIt failed on element selection.");
     233    "PostStepDoIt failed on element selection.");
    246234  }
    247235
    … …  
    275263    {
    276264      // Call the interaction
    277 
    278       G4HadronicInteractionWrapper aW;
    279       result = aW.ApplyInteraction(thePro, targetNucleus, theInteraction,
    280                                    GetProcessName(),
    281                                    theInteraction->GetModelName());
     265      result = theInteraction->ApplyYourself( thePro, targetNucleus);
    282266    }
    283267    catch(G4HadReentrentException aR)
    284268    {
     269      G4HadronicWhiteBoard & theBoard = G4HadronicWhiteBoard::Instance();
     270      theBoard.SetProjectile(thePro);
     271      theBoard.SetTargetNucleus(targetNucleus);
     272      theBoard.SetProcessName(GetProcessName());
     273      theBoard.SetModelName(theInteraction->GetModelName());
     274
    285275      aR.Report(G4cout);
    286276      G4cout << " G4HadronicProcess re-entering the ApplyYourself call for "
    … …  
    294284      {
    295285        G4Exception("G4HadronicProcess", "007", FatalException,
    296         "GetHadronicProcess: Reentering ApplyYourself too often - GeneralPostStepDoIt failed."); 
     286        "GetHadronicProcess: Reentering ApplyYourself too often - PostStepDoIt failed."); 
    297287      }
    298288      G4Exception("G4HadronicProcess", "007", FatalException,
    299       "GetHadronicProcess: GeneralPostStepDoIt failed (Reentering ApplyYourself not yet supported.)"); 
     289      "GetHadronicProcess: PostStepDoIt failed (Reentering ApplyYourself not yet supported.)"); 
    300290    }
    301291    catch(G4HadronicException aR)
    302292    {
     293      G4HadronicWhiteBoard & theBoard = G4HadronicWhiteBoard::Instance();
     294      theBoard.SetProjectile(thePro);
     295      theBoard.SetTargetNucleus(targetNucleus);
     296      theBoard.SetProcessName(GetProcessName());
     297      theBoard.SetModelName(theInteraction->GetModelName());
     298
    303299      aR.Report(G4cout);
    304300      G4cout << " G4HadronicProcess failed in ApplyYourself call for"
    … …  
    309305             << aParticle->GetDefinition()->GetParticleName() << G4endl;
    310306      G4Exception("G4HadronicProcess", "007", FatalException,
    311       "GeneralPostStepDoIt failed.");
     307      "PostStepDoIt failed.");
    312308    }
    313309  }
    … …  
    331327  result->SetTrafoToLab(thePro.GetTrafoToLab());
    332328
    333   /*
    334   // Loop over charged ion secondaries
    335 
    336   for(G4int i=0; i<result->GetNumberOfSecondaries(); i++)
    337   {
    338     G4DynamicParticle* aSecTrack = result->GetSecondary(i)->GetParticle();
    339     if(aSecTrack->GetDefinition()->GetPDGCharge()>1.5)
    340     {
    341       G4EffectiveCharge aCalculator;
    342       G4double charge =
    343            aCalculator.GetCharge(aMaterial, aSecTrack->GetKineticEnergy(),
    344                                  aSecTrack->GetDefinition()->GetPDGMass(),
    345                                  aSecTrack->GetDefinition()->GetPDGCharge());
    346       if(getenv("GHADChargeDebug"))
    347       {
    348         std::cout << "Recoil fractional charge is "
    349         << charge/aSecTrack->GetDefinition()->GetPDGCharge()<<" "
    350         << charge <<" "<<aSecTrack->GetDefinition()->GetPDGCharge()<<std::endl;
    351       }
    352       aSecTrack->SetCharge(charge);
    353     }
    354   }
    355   */
    356 
    357329  if(getenv("HadronicDoitLogging") )
    358330  {
    359331    G4cout << "HadronicDoitLogging "
    360     << GetProcessName() <<" "
    361     << aParticle->GetDefinition()->GetPDGEncoding()<<" "
    362     << originalEnergy<<" "
    363     << aParticle->GetMomentum()<<" "
    364     << targetNucleus.GetN()<<" "
    365     << targetNucleus.GetZ()<<" "
    366     << G4endl;
     332           << GetProcessName() <<" "
     333           << aParticle->GetDefinition()->GetPDGEncoding()<<" "
     334           << originalEnergy<<" "
     335           << aParticle->GetMomentum()<<" "
     336           << targetNucleus.GetN()<<" "
     337           << targetNucleus.GetZ()<<" "
     338           << G4endl;
    367339  }
    368340
    … …  
    509481G4HadronicProcess::FillTotalResult(G4HadFinalState * aR, const G4Track & aT)
    510482{
    511   G4Nancheck go_wild;
     483  //  G4Nancheck go_wild;
    512484  theTotalResult->Clear();
    513485  theTotalResult->ProposeLocalEnergyDeposit(0.);
    … …  
    559531      theTotalResult->ProposeParentWeight( XBiasSurvivalProbability()*aT.GetWeight() );
    560532      G4double newWeight = aR->GetWeightChange()*aT.GetWeight();
     533      /*
    561534      if(go_wild(aR->GetEnergyChange()))
    562535      {
    … …  
    571544        "surviving track received NaN momentum."); 
    572545      }
     546      */
    573547      G4double newM=aT.GetDefinition()->GetPDGMass();
    574548      G4double newE=aR->GetEnergyChange() + newM;
    … …  
    585559      if(aR->GetEnergyChange()>-.5)
    586560      {
     561        /*
    587562        if(go_wild(aR->GetEnergyChange()))
    588563        {
    … …  
    590565          "track received NaN energy."); 
    591566        }
     567        */
    592568        theTotalResult->ProposeEnergy(aR->GetEnergyChange());
    593569      }
    … …  
    606582  if(GetProcessName() != "hElastic" && GetProcessName() != "HadronElastic"
    607583     &&  theTotalResult->GetTrackStatus()==fAlive
    608      && aR->GetStatusChange()==isAlive
    609     )
    610   {
     584     && aR->GetStatusChange()==isAlive)
     585    {
    611586    // Use for debugging:   G4double newWeight = theTotalResult->GetParentWeight();
    612587
    … …  
    643618    theM.rotate(rotation, it);
    644619    theM*=aR->GetTrafoToLab();
    645 
     620    /*
    646621    if(go_wild(theM.e()))
    647622    {
    … …  
    656631      "secondary track received NaN momentum."); 
    657632    }
    658 
     633    */
    659634    aR->GetSecondary(i)->GetParticle()->Set4Momentum(theM);
    660635    G4double time = aR->GetSecondary(i)->GetTime();
    … …  
    662637
    663638    G4Track* track = new G4Track(aR->GetSecondary(i)->GetParticle(),
    664     time,
    665     aT.GetPosition());
     639                                time,
     640                                aT.GetPosition());
    666641
    667642    G4double newWeight = aT.GetWeight()*aR->GetSecondary(i)->GetWeight();
    … …  
    679654    // <<G4endl;
    680655    track->SetWeight(newWeight);
     656    /*
    681657    G4double trackDeb = track->GetKineticEnergy();
    682658    if( (  trackDeb<0
    … …  
    689665        <<G4endl;
    690666      }
    691       track->SetTouchableHandle(aT.GetTouchableHandle());
    692       theTotalResult->AddSecondary(track);
     667    */
     668    track->SetTouchableHandle(aT.GetTouchableHandle());
     669    theTotalResult->AddSecondary(track);
    693670  }
    694671
    … …  
    696673  return;
    697674}
     675
     676G4IsoParticleChange* G4HadronicProcess::GetIsotopeProductionInfo()
     677{
     678  G4IsoParticleChange * anIsoResult = theIsoResult;
     679  if(theIsoResult) theOldIsoResult = theIsoResult;
     680  theIsoResult = 0;
     681  return anIsoResult;
     682}
     683
     684void G4HadronicProcess::BiasCrossSectionByFactor(G4double aScale)
     685{
     686  xBiasOn = true;
     687  aScaleFactor = aScale;
     688  G4String it = GetProcessName();
     689  if( (it != "PhotonInelastic") &&
     690      (it != "ElectroNuclear") &&
     691      (it != "PositronNuclear") )
     692    {
     693      G4Exception("G4HadronicProcess", "007", FatalException,
     694                  "Cross-section biasing available only for gamma and electro nuclear reactions.");
     695    }
     696  if(aScale<100)
     697    {
     698      G4Exception("G4HadronicProcess", "001", JustWarning,
     699                  "Cross-section bias readjusted to be above safe limit. New value is 100");
     700      aScaleFactor = 100.;
     701    }
     702}
    698703/* end of file */
Note: See TracChangeset for help on using the changeset viewer.