Geant4  10.00.p02
G4TwoBodyAngularDist.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 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4TwoBodyAngularDist.cc 71652 2013-06-19 17:20:45Z mkelsey $
27 // Author: Michael Kelsey (SLAC)
28 // Date: 21 February 2013
29 //
30 // Description: Singleton class to evaluate two-body angular distribution
31 // functions based on intial/final state codes.
32 //
33 // 20130307 M. Kelsey -- Add verbosity interface
34 // 20130422 M. Kelsey -- Add three-body distributions, for temporary use
35 // 20130619 Change singleton instance to be thread-local, to avoid collisions.
36 
37 #include "G4TwoBodyAngularDist.hh"
38 #include "G4GamP2NPipAngDst.hh"
39 #include "G4GamP2PPi0AngDst.hh"
40 #include "G4GammaNuclAngDst.hh"
41 #include "G4PP2PPAngDst.hh"
42 #include "G4NP2NPAngDst.hh"
43 #include "G4HadNElastic1AngDst.hh"
44 #include "G4HadNElastic2AngDst.hh"
45 #include "G4InuclParticleNames.hh"
46 #include "G4NuclNuclAngDst.hh"
47 #include "G4HadNucl3BodyAngDst.hh"
48 #include "G4NuclNucl3BodyAngDst.hh"
49 #include "G4PiNInelasticAngDst.hh"
50 #include "G4InuclParticleNames.hh"
51 using namespace G4InuclParticleNames;
52 
53 // Singleton is created at first invocation
54 
56 
58  if (!theInstance) theInstance = new G4TwoBodyAngularDist;
59  return theInstance;
60 }
61 
62 // Constructor and destructor
63 
65  : gp_npip(new G4GamP2NPipAngDst), gp_ppi0(new G4GamP2PPi0AngDst),
66  ppAngDst(new G4PP2PPAngDst), npAngDst(new G4NP2NPAngDst),
67  nnAngDst(new G4NuclNuclAngDst), qxAngDst(new G4PiNInelasticAngDst),
68  hn1AngDst(new G4HadNElastic1AngDst), hn2AngDst(new G4HadNElastic2AngDst),
69  gnAngDst(new G4GammaNuclAngDst), hn3BodyDst(new G4HadNucl3BodyAngDst),
70  nn3BodyDst(new G4NuclNucl3BodyAngDst)
71 {;}
72 
74  delete gp_npip;
75  delete gp_ppi0;
76  delete ppAngDst;
77  delete nnAngDst;
78  delete qxAngDst;
79  delete hn1AngDst;
80  delete hn2AngDst;
81  delete gnAngDst;
82  delete npAngDst;
83  delete hn3BodyDst;
84  delete nn3BodyDst;
85 }
86 
87 
88 // Set verbosity for all generators (const-cast required)
89 
91  const_cast<G4TwoBodyAngularDist*>(GetInstance())->passVerbose(verbose);
92 }
93 
95  if (gp_npip) gp_npip->setVerboseLevel(verbose);
96  if (gp_ppi0) gp_ppi0->setVerboseLevel(verbose);
97  if (ppAngDst) ppAngDst->setVerboseLevel(verbose);
98  if (nnAngDst) nnAngDst->setVerboseLevel(verbose);
99  if (qxAngDst) qxAngDst->setVerboseLevel(verbose);
100  if (hn1AngDst) hn1AngDst->setVerboseLevel(verbose);
101  if (hn2AngDst) hn2AngDst->setVerboseLevel(verbose);
102  if (gnAngDst) gnAngDst->setVerboseLevel(verbose);
103  if (npAngDst) npAngDst->setVerboseLevel(verbose);
104  if (hn3BodyDst) hn3BodyDst->setVerboseLevel(verbose);
105  if (nn3BodyDst) nn3BodyDst->setVerboseLevel(verbose);
106 }
107 
108 
109 // Return appropriate distribution generator for specified interaction
110 
111 const G4VTwoBodyAngDst*
113  // TEMPORARY: Three-body distributions for hN/NN
114  if (fs==0 && kw==0) {
115  if (is == pro*pro || is == pro*neu || is == neu*neu) return nn3BodyDst;
116  else return hn3BodyDst;
117  }
118 
119  // gamma-nucleon -> nucleon pi0
120  if ((is == gam*pro && fs == pro*pi0) ||
121  (is == gam*neu && fs == neu*pi0)) {
122  return gp_ppi0;
123  }
124 
125  // gamma-nucleon charge exchange
126  if ((is == gam*pro && fs == neu*pip) ||
127  (is == gam*neu && fs == pro*pim)) {
128  return gp_npip;
129  }
130 
131  // pp and nn elastic
132  if (is == pro*pro || is == neu*neu) return ppAngDst;
133 
134  // np and pn elastic
135  if (is == pro*neu) return npAngDst;
136 
137  // hyperon-nucleon
138  if (is == pro*lam || is == pro*sp || is == pro*s0 ||
139  is == pro*sm || is == pro*xi0 || is == pro*xim ||
140  is == pro*om ||
141  is == neu*lam || is == neu*sp || is == neu*s0 ||
142  is == neu*sm || is == neu*xi0 || is == neu*xim ||
143  is == neu*om) {
144  return nnAngDst;
145  }
146 
147  // gamma p -> pi+ n, gamma p -> pi0 p, gamma p -> K Y (and isospin variants)
148  if (kw == 2 && (is == pro*gam || is == neu*gam)) {
149  return gnAngDst;
150  }
151 
152  // pion-nucleon charge/strangeness exchange
153  if (kw == 2) {
154  return qxAngDst;
155  }
156 
157  // pi+p, pi0p, gamma p, k+p, k0bp, pi-n, pi0n, gamma n, k-n, or k0n
158  if (is == pro*pip || is == pro*pi0 || is == pro*gam ||
159  is == pro*kpl || is == pro*k0b ||
160  is == neu*pim || is == neu*pi0 || is == neu*gam ||
161  is == neu*kmi || is == neu*k0) {
162  return hn1AngDst;
163  }
164 
165  // pi-p, pi+n, k-p, k0bn, k+n, or k0p
166  if (is == pro*pim || is == pro*kmi || is == pro*k0 ||
167  is == neu*pip || is == neu*kpl || is == neu*k0b) {
168  return hn2AngDst;
169  }
170 
171  // Invalid interaction
172  return 0;
173 }
G4NuclNucl3BodyAngDst * nn3BodyDst
static G4ThreadLocal G4TwoBodyAngularDist * theInstance
G4GamP2PPi0AngDst * gp_ppi0
G4PiNInelasticAngDst * qxAngDst
static const G4TwoBodyAngularDist * GetInstance()
#define G4ThreadLocal
Definition: tls.hh:52
int G4int
Definition: G4Types.hh:78
G4HadNElastic2AngDst * hn2AngDst
G4NuclNuclAngDst * nnAngDst
G4GamP2NPipAngDst * gp_npip
const G4VTwoBodyAngDst * ChooseDist(G4int is, G4int fs, G4int kw) const
virtual void setVerboseLevel(G4int verbose=0)
static void setVerboseLevel(G4int vb=0)
G4HadNElastic1AngDst * hn1AngDst
G4GammaNuclAngDst * gnAngDst
G4HadNucl3BodyAngDst * hn3BodyDst
void passVerbose(G4int verbose)