Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:59:20

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // G4VDecayChannel
0027 //
0028 // Class description:
0029 //
0030 // Abstract class to describe decay kinematics
0031 
0032 // Author: H.Kurashige, 27 July 1996
0033 // --------------------------------------------------------------------
0034 #ifndef G4VDecayChannel_hh
0035 #define G4VDecayChannel_hh 1
0036 
0037 #include "G4AutoLock.hh"
0038 #include "G4Threading.hh"
0039 #include "G4ThreeVector.hh"
0040 #include "G4ios.hh"
0041 #include "globals.hh"
0042 
0043 #include <cmath>
0044 
0045 class G4ParticleDefinition;
0046 class G4DecayProducts;
0047 class G4ParticleTable;
0048 
0049 class G4VDecayChannel
0050 {
0051   public:
0052     // Constructors
0053     G4VDecayChannel(const G4String& aName, G4int Verbose = 1);
0054     G4VDecayChannel(const G4String& aName, const G4String& theParentName, G4double theBR,
0055                     G4int theNumberOfDaughters, const G4String& theDaughterName1,
0056                     const G4String& theDaughterName2 = "", const G4String& theDaughterName3 = "",
0057                     const G4String& theDaughterName4 = "", const G4String& theDaughterName5 = "");
0058 
0059     // Destructor
0060     virtual ~G4VDecayChannel();
0061 
0062     // Equality operators
0063     G4bool operator==(const G4VDecayChannel& r) const { return (this == &r); }
0064     G4bool operator!=(const G4VDecayChannel& r) const { return (this != &r); }
0065 
0066     // Less-than operator is defined for G4DecayTable
0067     inline G4bool operator<(const G4VDecayChannel& right) const;
0068 
0069     virtual G4DecayProducts* DecayIt(G4double parentMass = -1.0) = 0;
0070 
0071     // Get kinematics name
0072     inline const G4String& GetKinematicsName() const;
0073 
0074     // Get branching ratio
0075     inline G4double GetBR() const;
0076 
0077     // Get number of daughter particles
0078     inline G4int GetNumberOfDaughters() const;
0079 
0080     // Get the pointer to the parent particle
0081     inline G4ParticleDefinition* GetParent();
0082 
0083     // Get the pointer to a daughter particle
0084     inline G4ParticleDefinition* GetDaughter(G4int anIndex);
0085 
0086     // Get the angular momentum of the decay
0087     G4int GetAngularMomentum();
0088 
0089     // Get the name of the parent particle
0090     inline const G4String& GetParentName() const;
0091 
0092     // Get the name of a daughter particle
0093     inline const G4String& GetDaughterName(G4int anIndex) const;
0094 
0095     // Get mass of parent
0096     inline G4double GetParentMass() const;
0097 
0098     // Get mass of daughter
0099     inline G4double GetDaughterMass(G4int anIndex) const;
0100 
0101     // Set the parent particle (by name or by pointer)
0102     void SetParent(const G4ParticleDefinition* particle_type);
0103     inline void SetParent(const G4String& particle_name);
0104 
0105     // Set branching ratio
0106     void SetBR(G4double value);
0107 
0108     // Set number of daughter particles
0109     void SetNumberOfDaughters(G4int value);
0110 
0111     // Set a daughter particle (by name or by pointer)
0112     void SetDaughter(G4int anIndex, const G4ParticleDefinition* particle_type);
0113     void SetDaughter(G4int anIndex, const G4String& particle_name);
0114 
0115     inline void SetVerboseLevel(G4int value);
0116     inline G4int GetVerboseLevel() const;
0117     void DumpInfo();
0118 
0119     inline G4double GetRangeMass() const;
0120     inline void SetRangeMass(G4double val);
0121     virtual G4bool IsOKWithParentMass(G4double parentMass);
0122 
0123     void SetPolarization(const G4ThreeVector&);
0124     inline const G4ThreeVector& GetPolarization() const;
0125 
0126   protected:
0127     // Default constructor
0128     G4VDecayChannel();
0129 
0130     // Copy constructor and assignment operator
0131     G4VDecayChannel(const G4VDecayChannel&);
0132     G4VDecayChannel& operator=(const G4VDecayChannel&);
0133 
0134     // Clear daughters array
0135     void ClearDaughtersName();
0136 
0137     inline void CheckAndFillDaughters();
0138     inline void CheckAndFillParent();
0139 
0140     G4double DynamicalMass(G4double massPDG, G4double width, G4double maxDev = 1.0) const;
0141 
0142   protected:
0143     // Kinematics name
0144     G4String kinematics_name = "";
0145     // Branching ratio  [0.0 - 1.0]
0146     G4double rbranch = 0.0;
0147     // Parent particle
0148     G4String* parent_name = nullptr;
0149     // Daughter particles
0150     G4String** daughters_name = nullptr;
0151 
0152     // Range of mass allowed in decay
0153     G4double rangeMass = 2.5;
0154 
0155     // Polarization of the parent particle
0156     G4ThreeVector parent_polarization;
0157 
0158     // Pointer to particle table
0159     G4ParticleTable* particletable = nullptr;
0160 
0161     static const G4String noName;
0162 
0163     G4ParticleDefinition* G4MT_parent = nullptr;
0164     G4ParticleDefinition** G4MT_daughters = nullptr;
0165     G4double G4MT_parent_mass = 0.0;
0166     G4double* G4MT_daughters_mass = nullptr;
0167     G4double* G4MT_daughters_width = nullptr;
0168     G4Mutex daughtersMutex;
0169     G4Mutex parentMutex;
0170 
0171     // Number of daughters
0172     G4int numberOfDaughters = 0;
0173 
0174     // Control flag for output message
0175     //  0: Silent
0176     //  1: Warning message
0177     //  2: More
0178     G4int verboseLevel = 1;
0179 
0180   private:
0181     // Fill daughters array
0182     void FillDaughters();
0183 
0184     // Fill parent
0185     void FillParent();
0186 
0187     const G4String& GetNoName() const;
0188 };
0189 
0190 // ------------------------------------------------------------
0191 // Inline methods
0192 // ------------------------------------------------------------
0193 
0194 inline G4bool G4VDecayChannel::operator<(const G4VDecayChannel& right) const
0195 {
0196   return (this->rbranch < right.rbranch);
0197 }
0198 
0199 inline G4ParticleDefinition* G4VDecayChannel::GetDaughter(G4int anIndex)
0200 {
0201   // pointers to daughter particles are filled, if they are not set yet
0202   CheckAndFillDaughters();
0203 
0204   // get the pointer to a daughter particle
0205   if ((anIndex >= 0) && (anIndex < numberOfDaughters)) {
0206     return G4MT_daughters[anIndex];
0207   }
0208 
0209   if (verboseLevel > 0)
0210     G4cout << "G4VDecayChannel::GetDaughter  index out of range " << anIndex << G4endl;
0211   return nullptr;
0212 }
0213 
0214 inline const G4String& G4VDecayChannel::GetDaughterName(G4int anIndex) const
0215 {
0216   if ((anIndex >= 0) && (anIndex < numberOfDaughters)) {
0217     return *daughters_name[anIndex];
0218   }
0219 
0220   if (verboseLevel > 0) {
0221     G4cout << "G4VDecayChannel::GetDaughterName ";
0222     G4cout << "index out of range " << anIndex << G4endl;
0223   }
0224   return GetNoName();
0225 }
0226 
0227 inline G4double G4VDecayChannel::GetDaughterMass(G4int anIndex) const
0228 {
0229   if ((anIndex >= 0) && (anIndex < numberOfDaughters)) {
0230     return G4MT_daughters_mass[anIndex];
0231   }
0232 
0233   if (verboseLevel > 0) {
0234     G4cout << "G4VDecayChannel::GetDaughterMass ";
0235     G4cout << "index out of range " << anIndex << G4endl;
0236   }
0237   return 0.0;
0238 }
0239 
0240 inline G4ParticleDefinition* G4VDecayChannel::GetParent()
0241 {
0242   // the pointer to the parent particle is filled, if it is not set yet
0243   CheckAndFillParent();
0244   // get the pointer to the parent particle
0245   return G4MT_parent;
0246 }
0247 
0248 inline const G4String& G4VDecayChannel::GetParentName() const
0249 {
0250   return *parent_name;
0251 }
0252 
0253 inline G4double G4VDecayChannel::GetParentMass() const
0254 {
0255   return G4MT_parent_mass;
0256 }
0257 
0258 inline void G4VDecayChannel::SetParent(const G4String& particle_name)
0259 {
0260   delete parent_name;
0261   parent_name = new G4String(particle_name);
0262   G4MT_parent = nullptr;
0263 }
0264 
0265 inline G4int G4VDecayChannel::GetNumberOfDaughters() const
0266 {
0267   return numberOfDaughters;
0268 }
0269 
0270 inline const G4String& G4VDecayChannel::GetKinematicsName() const
0271 {
0272   return kinematics_name;
0273 }
0274 
0275 inline G4double G4VDecayChannel::GetBR() const
0276 {
0277   return rbranch;
0278 }
0279 
0280 inline void G4VDecayChannel::SetVerboseLevel(G4int value)
0281 {
0282   verboseLevel = value;
0283 }
0284 
0285 inline G4int G4VDecayChannel::GetVerboseLevel() const
0286 {
0287   return verboseLevel;
0288 }
0289 
0290 inline G4double G4VDecayChannel::GetRangeMass() const
0291 {
0292   return rangeMass;
0293 }
0294 
0295 inline void G4VDecayChannel::SetRangeMass(G4double val)
0296 {
0297   if (val >= 0.) rangeMass = val;
0298 }
0299 
0300 inline void G4VDecayChannel::SetPolarization(const G4ThreeVector& polar)
0301 {
0302   parent_polarization = polar;
0303 }
0304 
0305 inline const G4ThreeVector& G4VDecayChannel::GetPolarization() const
0306 {
0307   return parent_polarization;
0308 }
0309 
0310 inline void G4VDecayChannel::CheckAndFillDaughters()
0311 {
0312   G4AutoLock l(&daughtersMutex);
0313   if (G4MT_daughters == nullptr) {
0314     l.unlock();
0315     FillDaughters();
0316   }
0317 }
0318 
0319 inline void G4VDecayChannel::CheckAndFillParent()
0320 {
0321   G4AutoLock l(&parentMutex);
0322   if (G4MT_parent == nullptr) {
0323     l.unlock();
0324     FillParent();
0325   }
0326 }
0327 
0328 #endif