Geant4  10.00.p01
G4VDecayChannel.hh
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 //
27 // $Id: G4VDecayChannel.hh 73762 2013-09-10 12:53:56Z gcosmo $
28 //
29 //
30 // ------------------------------------------------------------
31 // GEANT 4 class header file
32 //
33 // History: first implementation, based on object model of
34 // 27 July 1996 H.Kurashige
35 // 30 May 1997 H.Kurashige
36 // 23 Mar. 2000 H.Weber : add GetAngularMomentum()
37 // ------------------------------------------------------------
38 #ifndef G4VDecayChannel_h
39 #define G4VDecayChannel_h 1
40 
41 #include <cmath>
42 
43 #include "G4ios.hh"
44 #include "globals.hh"
46 
48 class G4DecayProducts;
49 class G4ParticleTable;
50 
102 {
103  // Class Description:
104  //
105  // Abstract class to describe decay kinematics
106 
107  public: // with description
108 
109  // Constructors
110  G4VDecayChannel(const G4String &aName, G4int Verbose = 1);
111  G4VDecayChannel(const G4String &aName,
112  const G4String& theParentName,
113  G4double theBR,
114  G4int theNumberOfDaughters,
115  const G4String& theDaughterName1,
116  const G4String& theDaughterName2 = "",
117  const G4String& theDaughterName3 = "",
118  const G4String& theDaughterName4 = "" );
119 
120  // Destructor
121  virtual ~G4VDecayChannel();
122 
123  // equality operators
124  G4int operator==(const G4VDecayChannel &right) const {return (this == &right);}
125  G4int operator!=(const G4VDecayChannel &right) const {return (this != &right);}
126 
127  // less-than operator is defined for G4DecayTable
128  G4int operator<(const G4VDecayChannel &right) const;
129 
130  virtual G4DecayProducts* DecayIt(G4double parentMass = -1.0) = 0;
131 
132  // get kinematics name
133  const G4String& GetKinematicsName() const;
134  // get branching ratio
135  G4double GetBR() const;
136  // get number of daughter particles
137  G4int GetNumberOfDaughters() const;
138 
139  // get the pointer to the parent particle
141  // get the pointer to a daughter particle
143 
144  // get the angular momentum of the decay
146  // get the name of the parent particle
147  const G4String& GetParentName() const;
148  //get the name of a daughter particle
149  const G4String& GetDaughterName(G4int anIndex) const;
150 
151  // get mass of parent
152  G4double GetParentMass() const;
153  G4double GetDaughterMass(G4int anIndex) const;
154 
155  // set the parent particle (by name or by pointer)
156  void SetParent(const G4ParticleDefinition * particle_type);
157  void SetParent(const G4String &particle_name);
158  // set branching ratio
159  void SetBR(G4double value);
160  // set number of daughter particles
161  void SetNumberOfDaughters(G4int value);
162  // set a daughter particle (by name or by pointer)
163  void SetDaughter(G4int anIndex,
164  const G4ParticleDefinition * particle_type);
165  void SetDaughter(G4int anIndex,
166  const G4String &particle_name);
167 
168  void SetVerboseLevel(G4int value);
169  G4int GetVerboseLevel() const;
170  void DumpInfo();
171 
178 
179  protected: // with description
180 
181  // clear daughters array
182  void ClearDaughtersName();
183 
184  // fill daughters array
185  void FillDaughters();
186  // fill parent
187  void FillParent();
188 
189  protected: // without description
190 
191  // default constructor
192  G4VDecayChannel();
193 
194  // copy constructor and assignment operator
197 
198  private:
199 
200  const G4String& GetNoName() const;
201 
202  protected:
203 
204  // kinematics name
206  // branching ratio [0.0 - 1.0]
208  // number of daughters
210  // parent particle
212  // daughter particles
214 
215  // pointer to particle table
217 
218  // control flag for output message
220  // 0: Silent
221  // 1: Warning message
222  // 2: More
223 
224  static const G4String noName;
225 
234 };
235 
236 // ------------------------------------------------------------
237 // Inline methods
238 // ------------------------------------------------------------
239 
240 inline
242 {
243  return (this->rbranch < right.rbranch);
244 }
245 
246 inline
248  {
249  // pointers to daughter particles are filled, if they are not set yet
250  if (G4MT_daughters == 0) FillDaughters();
251 
252  // get the pointer to a daughter particle
253  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
254  return G4MT_daughters[anIndex];
255  } else {
256  if (verboseLevel>0)
257  G4cout << "G4VDecayChannel::GetDaughter index out of range "<<anIndex<<G4endl;
258  return 0;
259  }
260 }
261 
262 inline
264 {
265  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
266  return *daughters_name[anIndex];
267  } else {
268  if (verboseLevel>0){
269  G4cout << "G4VDecayChannel::GetDaughterName ";
270  G4cout << "index out of range " << anIndex << G4endl;
271  }
272  return GetNoName();
273  }
274 }
275 
276 inline
278 {
279  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
280  return G4MT_daughters_mass[anIndex];
281  } else {
282  if (verboseLevel>0){
283  G4cout << "G4VDecayChannel::GetDaughterMass ";
284  G4cout << "index out of range " << anIndex << G4endl;
285  }
286  return 0.0;
287  }
288 }
289 
290 inline
292 {
293  // the pointer to the parent particle is filled, if it is not set yet
294  if (G4MT_parent == 0) FillParent();
295  // get the pointer to the parent particle
296  return G4MT_parent;
297 }
298 
299 inline
301 {
302  return *parent_name;
303 }
304 
305 inline
307 {
308  return G4MT_parent_mass;
309 }
310 
311 inline
312  void G4VDecayChannel::SetParent(const G4String &particle_name)
313 {
314  if (parent_name != 0) delete parent_name;
315  parent_name = new G4String(particle_name);
316  G4MT_parent = 0;
317 }
318 
319 inline
321 {
322  return numberOfDaughters;
323 }
324 
325 inline
327 
328 inline
329  void G4VDecayChannel::SetBR(G4double value){ rbranch = value; }
330 
331 inline
333 
334 inline
336 
337 inline
339 
342 
343 #endif
G4double GetBR() const
const G4String & GetNoName() const
static const G4String noName
void SetBR(G4double value)
G4double * G4MT_daughters_mass
const G4String & GetParentName() const
G4int GetNumberOfDaughters() const
G4ParticleDefinition * G4MT_parent
G4ParticleDefinition * GetDaughter(G4int anIndex)
G4ParticleDefinition ** G4MT_daughters
const G4String & GetKinematicsName() const
G4double G4MT_parent_mass
G4int operator!=(const G4VDecayChannel &right) const
G4int operator<(const G4VDecayChannel &right) const
G4ParticleTable * particletable
int G4int
Definition: G4Types.hh:78
G4int GetAngularMomentum()
G4int operator==(const G4VDecayChannel &right) const
G4GLOB_DLL std::ostream G4cout
void SetVerboseLevel(G4int value)
G4String kinematics_name
void SetNumberOfDaughters(G4int value)
G4double GetParentMass() const
G4String * parent_name
G4VDecayChannel & operator=(const G4VDecayChannel &)
G4double GetDaughterMass(G4int anIndex) const
const G4String & GetDaughterName(G4int anIndex) const
virtual ~G4VDecayChannel()
G4int GetVerboseLevel() const
void SetParent(const G4ParticleDefinition *particle_type)
void SetDaughter(G4int anIndex, const G4ParticleDefinition *particle_type)
#define G4endl
Definition: G4ios.hh:61
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
double G4double
Definition: G4Types.hh:76
G4ParticleDefinition * GetParent()
G4String ** daughters_name