Geant4  10.02.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 93024 2015-09-30 15:59:32Z 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"
45 #include "G4ThreeVector.hh"
47 
49 class G4DecayProducts;
50 class G4ParticleTable;
51 
106 {
107  // Class Description:
108  //
109  // Abstract class to describe decay kinematics
110 
111  public: // with description
112 
113  // Constructors
114  G4VDecayChannel(const G4String &aName, G4int Verbose = 1);
115  G4VDecayChannel(const G4String &aName,
116  const G4String& theParentName,
117  G4double theBR,
118  G4int theNumberOfDaughters,
119  const G4String& theDaughterName1,
120  const G4String& theDaughterName2 = "",
121  const G4String& theDaughterName3 = "",
122  const G4String& theDaughterName4 = "" );
123 
124  // Destructor
125  virtual ~G4VDecayChannel();
126 
127  // equality operators
128  G4int operator==(const G4VDecayChannel &right) const {return (this == &right);}
129  G4int operator!=(const G4VDecayChannel &right) const {return (this != &right);}
130 
131  // less-than operator is defined for G4DecayTable
132  G4int operator<(const G4VDecayChannel &right) const;
133 
134  virtual G4DecayProducts* DecayIt(G4double parentMass = -1.0) = 0;
135 
136  // get kinematics name
137  const G4String& GetKinematicsName() const;
138  // get branching ratio
139  G4double GetBR() const;
140  // get number of daughter particles
141  G4int GetNumberOfDaughters() const;
142 
143  // get the pointer to the parent particle
145  // get the pointer to a daughter particle
147 
148  // get the angular momentum of the decay
150  // get the name of the parent particle
151  const G4String& GetParentName() const;
152  //get the name of a daughter particle
153  const G4String& GetDaughterName(G4int anIndex) const;
154 
155  // get mass of parent
156  G4double GetParentMass() const;
157  G4double GetDaughterMass(G4int anIndex) const;
158 
159  // set the parent particle (by name or by pointer)
160  void SetParent(const G4ParticleDefinition * particle_type);
161  void SetParent(const G4String &particle_name);
162  // set branching ratio
163  void SetBR(G4double value);
164  // set number of daughter particles
165  void SetNumberOfDaughters(G4int value);
166  // set a daughter particle (by name or by pointer)
167  void SetDaughter(G4int anIndex,
168  const G4ParticleDefinition * particle_type);
169  void SetDaughter(G4int anIndex,
170  const G4String &particle_name);
171 
172  void SetVerboseLevel(G4int value);
173  G4int GetVerboseLevel() const;
174  void DumpInfo();
175 
176  G4double GetRangeMass() const;
177  void SetRangeMass(G4double val);
178  virtual G4bool IsOKWithParentMass(G4double parentMass);
179 
180  void SetPolarization(const G4ThreeVector&);
181  const G4ThreeVector& GetPolarization() const;
182 
188 
189  protected: // with description
190 
191  // clear daughters array
192  void ClearDaughtersName();
193 
194  // fill daughters array
195  void FillDaughters();
196  // fill parent
197  void FillParent();
198 
199  protected: // without description
200 
201  // default constructor
202  G4VDecayChannel();
203 
204  // copy constructor and assignment operator
207 
208  private:
209 
210  const G4String& GetNoName() const;
211 
212  protected:
213  G4double DynamicalMass(G4double massPDG, G4double width, G4double maxDev = +1.) const;
214 
215  protected:
216 
217  // kinematics name
219  // branching ratio [0.0 - 1.0]
221  // number of daughters
223  // parent particle
225  // daughter particles
227 
228  // range of mass allowed in decay
230 
231  // polarization of the parent particle
233 
234  // pointer to particle table
236 
237  // control flag for output message
239  // 0: Silent
240  // 1: Warning message
241  // 2: More
242 
243  static const G4String noName;
244 
254 };
255 
256 // ------------------------------------------------------------
257 // Inline methods
258 // ------------------------------------------------------------
259 
260 inline
262 {
263  return (this->rbranch < right.rbranch);
264 }
265 
266 inline
268  {
269  // pointers to daughter particles are filled, if they are not set yet
270  if (G4MT_daughters == 0) FillDaughters();
271 
272  // get the pointer to a daughter particle
273  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
274  return G4MT_daughters[anIndex];
275  } else {
276  if (verboseLevel>0)
277  G4cout << "G4VDecayChannel::GetDaughter index out of range "<<anIndex<<G4endl;
278  return 0;
279  }
280 }
281 
282 inline
284 {
285  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
286  return *daughters_name[anIndex];
287  } else {
288  if (verboseLevel>0){
289  G4cout << "G4VDecayChannel::GetDaughterName ";
290  G4cout << "index out of range " << anIndex << G4endl;
291  }
292  return GetNoName();
293  }
294 }
295 
296 inline
298 {
299  if ( (anIndex>=0) && (anIndex<numberOfDaughters) ) {
300  return G4MT_daughters_mass[anIndex];
301  } else {
302  if (verboseLevel>0){
303  G4cout << "G4VDecayChannel::GetDaughterMass ";
304  G4cout << "index out of range " << anIndex << G4endl;
305  }
306  return 0.0;
307  }
308 }
309 
310 inline
312 {
313  // the pointer to the parent particle is filled, if it is not set yet
314  if (G4MT_parent == 0) FillParent();
315  // get the pointer to the parent particle
316  return G4MT_parent;
317 }
318 
319 inline
321 {
322  return *parent_name;
323 }
324 
325 inline
327 {
328  return G4MT_parent_mass;
329 }
330 
331 inline
332  void G4VDecayChannel::SetParent(const G4String &particle_name)
333 {
334  if (parent_name != 0) delete parent_name;
335  parent_name = new G4String(particle_name);
336  G4MT_parent = 0;
337 }
338 
339 inline
341 {
342  return numberOfDaughters;
343 }
344 
345 inline
347 
348 inline
349  void G4VDecayChannel::SetBR(G4double value){ rbranch = value; }
350 
351 inline
353 
354 inline
356 
357 inline
359 
360 inline
362 
363 inline
364  void G4VDecayChannel::SetRangeMass(G4double val){ if(val>=0.) rangeMass=val; }
365 
366 inline
368 {
369  parent_polarization = polar;
370 }
371 
372 inline
374 {
375  return parent_polarization;
376 }
377 
378 
381 
382 #endif
G4double GetBR() const
const G4String & GetNoName() const
static const G4String noName
void SetBR(G4double value)
G4double * G4MT_daughters_mass
CLHEP::Hep3Vector G4ThreeVector
const G4String & GetParentName() const
G4int GetNumberOfDaughters() const
G4ParticleDefinition * G4MT_parent
G4ParticleDefinition * GetDaughter(G4int anIndex)
G4ParticleDefinition ** G4MT_daughters
const G4String & GetKinematicsName() const
#define width
G4double G4MT_parent_mass
virtual G4bool IsOKWithParentMass(G4double parentMass)
G4int operator!=(const G4VDecayChannel &right) const
G4int operator<(const G4VDecayChannel &right) const
G4ParticleTable * particletable
void SetRangeMass(G4double val)
int G4int
Definition: G4Types.hh:78
G4int GetAngularMomentum()
G4int operator==(const G4VDecayChannel &right) const
G4GLOB_DLL std::ostream G4cout
void SetVerboseLevel(G4int value)
bool G4bool
Definition: G4Types.hh:79
G4String kinematics_name
void SetNumberOfDaughters(G4int value)
G4double GetParentMass() const
G4String * parent_name
G4VDecayChannel & operator=(const G4VDecayChannel &)
G4ThreeVector parent_polarization
G4double GetDaughterMass(G4int anIndex) const
const G4String & GetDaughterName(G4int anIndex) const
virtual ~G4VDecayChannel()
G4int GetVerboseLevel() const
G4double * G4MT_daughters_width
void SetParent(const G4ParticleDefinition *particle_type)
void SetDaughter(G4int anIndex, const G4ParticleDefinition *particle_type)
#define G4endl
Definition: G4ios.hh:61
G4double GetRangeMass() const
virtual G4DecayProducts * DecayIt(G4double parentMass=-1.0)=0
G4double DynamicalMass(G4double massPDG, G4double width, G4double maxDev=+1.) const
double G4double
Definition: G4Types.hh:76
void SetPolarization(const G4ThreeVector &)
const G4ThreeVector & GetPolarization() const
G4ParticleDefinition * GetParent()
G4String ** daughters_name