/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

DsG4Scintillation.h
Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 //  * DISCLAIMER                                                       *
00004 //  *                                                                  *
00005 //  * The following disclaimer summarizes all the specific disclaimers *
00006 //  * of contributors to this software. The specific disclaimers,which *
00007 //  * govern, are listed with their locations in:                      *
00008 //  *   http://cern.ch/geant4/license                                  *
00009 //  *                                                                  *
00010 //  * Neither the authors of this software system, nor their employing *
00011 //  * institutes,nor the agencies providing financial support for this *
00012 //  * work  make  any representation or  warranty, express or implied, *
00013 //  * regarding  this  software system or assume any liability for its *
00014 //  * use.                                                             *
00015 //  *                                                                  *
00016 //  * This  code  implementation is the  intellectual property  of the *
00017 //  * GEANT4 collaboration.                                            *
00018 //  * By copying,  distributing  or modifying the Program (or any work *
00019 //  * based  on  the Program)  you indicate  your  acceptance of  this *
00020 //  * statement, and all its terms.                                    *
00021 //  ********************************************************************
00022 // 
00023 // 
00024 // 
00025 // //////////////////////////////////////////////////////////////////////
00026 //  Scintillation Light Class Implementation
00027 // //////////////////////////////////////////////////////////////////////
00028 // 
00029 //  File:        G4Scintillation.cc 
00030 //  Description: RestDiscrete Process - Generation of Scintillation Photons
00031 //  Version:     1.0
00032 //  Created:     1998-11-07  
00033 //  Author:      Peter Gumplinger
00034 //  Updated:     2005-08-17 by Peter Gumplinger
00035 //               > change variable name MeanNumPhotons -> MeanNumberOfPhotons
00036 //               2005-07-28 by Peter Gumplinger
00037 //               > add G4ProcessType to constructor
00038 //               2004-08-05 by Peter Gumplinger
00039 //               > changed StronglyForced back to Forced in GetMeanLifeTime
00040 //               2002-11-21 by Peter Gumplinger
00041 //               > change to use G4Poisson for small MeanNumberOfPhotons
00042 //               2002-11-07 by Peter Gumplinger
00043 //               > now allow for fast and slow scintillation component
00044 //               2002-11-05 by Peter Gumplinger
00045 //               > now use scintillation constants from G4Material
00046 //               2002-05-09 by Peter Gumplinger
00047 //               > use only the PostStepPoint location for the origin of
00048 //                scintillation photons when energy is lost to the medium
00049 //                by a neutral particle
00050 //                2000-09-18 by Peter Gumplinger
00051 //               > change: aSecondaryPosition=x0+rand*aStep.GetDeltaPosition();
00052 //                aSecondaryTrack->SetTouchable(0);
00053 //                2001-09-17, migration of Materials to pure STL (mma) 
00054 //                2003-06-03, V.Ivanchenko fix compilation warnings
00055 //    
00056 //mail:        gum@triumf.ca
00057 //               
00059 
00060 //-----------------------------------------------------------------
00061 // DsG4Scintillation
00062 // 
00063 // A class modified from G4Scintillation.
00064 // Birks' law is used to calculate the scintillation photon number 
00065 // Author: Liang Zhan, 2006/01/27
00066 // Modified: bv@bnl.gov, 2008/4/16 for DetSim
00067 //--------------------------------------------------------------
00068 
00069 #ifndef DsG4Scintillation_h
00070 #define DsG4Scintillation_h 1
00071 
00072 #include "globals.hh"
00073 #include "templates.hh"
00074 #include "Randomize.hh"
00075 #include "G4Poisson.hh"
00076 #include "G4ThreeVector.hh"
00077 #include "G4ParticleMomentum.hh"
00078 #include "G4Step.hh"
00079 #include "G4VRestDiscreteProcess.hh"
00080 #include "G4OpticalPhoton.hh"
00081 #include "G4DynamicParticle.hh"
00082 #include "G4Material.hh" 
00083 #include "G4PhysicsTable.hh"
00084 #include "G4MaterialPropertiesTable.hh"
00085 #include "G4PhysicsOrderedFreeVector.hh"
00086 #include "G4UImessenger.hh"
00087 #include "DsPhysConsOptical.h"
00088 #include <fstream>
00089 #include <iostream>
00090 class G4UIcommand;
00091 class G4UIdirectory;
00092 
00094 // Class Description:
00095 // RestDiscrete Process - Generation of Scintillation Photons.
00096 // Class inherits publicly from G4VRestDiscreteProcess.
00097 // Class Description - End:
00098 
00100 // Class Definition
00102 
00103 class DsG4Scintillation : public G4VRestDiscreteProcess, public G4UImessenger
00104 { //too lazy to create another UImessenger class
00105 
00106 private:
00107 
00109         // Operators
00111 
00112         // DsG4Scintillation& operator=(const DsG4Scintillation &right);
00113 
00114 public: // Without description
00115 
00117         // Constructors and Destructor
00119 
00120         DsG4Scintillation(const G4String& processName = "Scintillation",
00121                           G4ProcessType type = fElectromagnetic);
00122 
00123         // DsG4Scintillation(const DsG4Scintillation &right);
00124 
00125         ~DsG4Scintillation();   
00126 
00128         // Methods
00130 
00131 public: // With description
00132 
00133         // DsG4Scintillation Process has both PostStepDoIt (for energy 
00134         // deposition of particles in flight) and AtRestDoIt (for energy
00135         // given to the medium by particles at rest)
00136 
00137         G4bool IsApplicable(const G4ParticleDefinition& aParticleType);
00138         // Returns true -> 'is applicable', for any particle type except
00139         // for an 'opticalphoton' 
00140 
00141         G4double GetMeanFreePath(const G4Track& aTrack,
00142                                        G4double ,
00143                                        G4ForceCondition* );
00144         // Returns infinity; i. e. the process does not limit the step,
00145         // but sets the 'StronglyForced' condition for the DoIt to be 
00146         // invoked at every step.
00147 
00148         G4double GetMeanLifeTime(const G4Track& aTrack,
00149                                  G4ForceCondition* );
00150         // Returns infinity; i. e. the process does not limit the time,
00151         // but sets the 'StronglyForced' condition for the DoIt to be
00152         // invoked at every step.
00153 
00154         G4VParticleChange* PostStepDoIt(const G4Track& aTrack, 
00155                                         const G4Step&  aStep);
00156         G4VParticleChange* AtRestDoIt (const G4Track& aTrack,
00157                                        const G4Step& aStep);
00158 
00159         // These are the methods implementing the scintillation process.
00160 
00161         void SetTrackSecondariesFirst(const G4bool state);
00162         // If set, the primary particle tracking is interrupted and any
00163         // produced scintillation photons are tracked next. When all 
00164         // have been tracked, the tracking of the primary resumes.
00165 
00166         G4bool GetTrackSecondariesFirst() const;
00167         // Returns the boolean flag for tracking secondaries first.
00168         
00169         void SetScintillationYieldFactor(const G4double yieldfactor);
00170         // Called to set the scintillation photon yield factor, needed when
00171         // the yield is different for different types of particles. This
00172         // scales the yield obtained from the G4MaterialPropertiesTable.
00173         G4double GetScintillationYieldFactor() const;
00174         // Returns the photon yield factor.
00175 
00176        void SetUseFastMu300nsTrick(const G4bool fastMu300nsTrick);
00177        G4bool GetUseFastMu300nsTrick() const;
00178 
00179 
00180         void SetScintillationExcitationRatio(const G4double excitationratio);
00181         // Called to set the scintillation exciation ratio, needed when
00182         // the scintillation level excitation is different for different
00183         // types of particles. This overwrites the YieldRatio obtained
00184         // from the G4MaterialPropertiesTable.
00185 
00186         G4double GetScintillationExcitationRatio() const;
00187         // Returns the scintillation level excitation ratio.
00188 
00189         G4PhysicsTable* GetFastIntegralTable() const;
00190         // Returns the address of the fast scintillation integral table.
00191 
00192         G4PhysicsTable* GetSlowIntegralTable() const;
00193         // Returns the address of the slow scintillation integral table.
00194 
00195         G4PhysicsTable* GetReemissionIntegralTable() const;
00196         // Returns the address of the reemission integral table.
00197 
00198         void DumpPhysicsTable() const;
00199         // Prints the fast and slow scintillation integral tables.
00200 
00201         // Configuration
00202         G4double GetPhotonWeight() const { return fPhotonWeight; }
00203         void SetPhotonWeight(G4double weight) { fPhotonWeight = weight; }
00204         void SetDoReemission(bool tf = true) { doReemission = tf; }
00205         bool GetDoReemission() { return doReemission; }
00206         void SetDoBothProcess(bool tf = true) { doBothProcess = tf; }
00207         bool GetDoBothProcess() { return doBothProcess; }
00208 
00209         G4bool GetApplyPreQE() const { return fApplyPreQE; }
00210         void SetApplyPreQE(G4bool a) { fApplyPreQE = a; }
00211 
00212         G4double GetPreQE() const { return fPreQE; }
00213         void SetPreQE(G4double a) { fPreQE = a; }
00214 
00215         void SetBirksConstant1(double c1) {birksConstant1 = c1;}
00216         double GetBirksConstant1() {return birksConstant1;}
00217         void SetBirksConstant2(double c2) {birksConstant2 = c2;}
00218         double GetBirksConstant2() {return birksConstant2;}
00219 
00220         void SetGammaSlowerTimeConstant(double st) { gammaSlowerTime = st;}
00221         double GetGammaSlowerTimeConstant() {return gammaSlowerTime;}
00222 
00223         void SetGammaSlowerRatio(double sr) { gammaSlowerRatio = sr;}
00224         double GetGammaSlowerRatio() {return gammaSlowerRatio;}
00225 
00226         void SetNeutronSlowerTimeConstant(double st) { neutronSlowerTime = st;}
00227         double GetNeutronSlowerTimeConstant() {return neutronSlowerTime;}
00228 
00229         void SetNeutronSlowerRatio(double sr) { neutronSlowerRatio = sr;}
00230         double GetNeutronSlowerRatio() {return neutronSlowerRatio;}
00231         
00232         void SetAlphaSlowerTimeConstant(double st) { alphaSlowerTime = st;}
00233         double GetAlphaSlowerTimeConstant() {return alphaSlowerTime;}
00234 
00235         void SetAlphaSlowerRatio(double sr) { alphaSlowerRatio = sr;}
00236         double GetAlphaSlowerRatio() {return alphaSlowerRatio;}
00237         
00238        
00239         // Don't actually do anything.  This is needed, as apposed to
00240         // simply not including the scintilation process, because
00241         // w/out this process no optical photons make it into the
00242         // photocathode (???) 
00243         void SetNoOp(bool tf = true) { m_noop = tf; }
00244 private:
00245 
00246         void BuildThePhysicsTable();
00247         // It builds either the fast or slow scintillation integral table; 
00248         // or both. 
00249 
00251         // Class Data Members
00253 
00254 protected:
00255 
00256         G4PhysicsTable* theSlowIntegralTable;
00257         G4PhysicsTable* theFastIntegralTable;
00258         G4PhysicsTable* theReemissionIntegralTable;
00259 
00260         // on/off flag for absorbed opticalphoton reemission
00261         G4bool doReemission;
00262         // choose only reemission of Cerenkov or both of Cerenkov and Scintillaton;
00263         G4bool doBothProcess;
00264 
00265         // Birks constant C1 and C2
00266         double  birksConstant1;
00267         double  birksConstant2;
00268 
00269         double  slowerTimeConstant;
00270         double  slowerRatio;
00271         
00272         double gammaSlowerTime;
00273         double gammaSlowerRatio;
00274         double neutronSlowerTime;
00275         double neutronSlowerRatio;
00276         double alphaSlowerTime;
00277         double alphaSlowerRatio;
00278 
00279 
00280 private:
00281 
00282         G4bool fTrackSecondariesFirst;
00283 
00284         G4double YieldFactor;
00285         G4bool FastMu300nsTrick;
00286         G4double ExcitationRatio;
00287   
00288         //mean number of true photons per secondary track in GLG4Scint
00289         G4double fPhotonWeight;
00290         G4bool   fApplyPreQE;
00291         G4double fPreQE;
00292         bool m_noop;
00293 
00294 };
00295 
00297 // Inline methods
00299 
00300 inline 
00301 G4bool DsG4Scintillation::IsApplicable(const G4ParticleDefinition& aParticleType)
00302 {
00303         if (aParticleType.GetParticleName() == "opticalphoton"){
00304            return true;
00305         } else {
00306            return true;
00307         }
00308 }
00309 
00310 inline 
00311 void DsG4Scintillation::SetTrackSecondariesFirst(const G4bool state) 
00312 { 
00313         fTrackSecondariesFirst = state;
00314 }
00315 
00316 inline
00317 G4bool DsG4Scintillation::GetTrackSecondariesFirst() const
00318 {
00319         return fTrackSecondariesFirst;
00320 }
00321 
00322 inline
00323 void DsG4Scintillation::SetScintillationYieldFactor(const G4double yieldfactor)
00324 {
00325         YieldFactor = yieldfactor;
00326 }
00327 
00328 
00329 inline
00330 G4double DsG4Scintillation::GetScintillationYieldFactor() const
00331 {
00332         return YieldFactor;
00333 }
00334 
00335 inline 
00336 void DsG4Scintillation::SetUseFastMu300nsTrick(const G4bool fastMu300nsTrick)
00337 {
00338         FastMu300nsTrick = fastMu300nsTrick;
00339 }
00340 
00341 inline 
00342 G4bool DsG4Scintillation::GetUseFastMu300nsTrick() const
00343 {
00344       return  FastMu300nsTrick;
00345 }
00346 
00347 
00348 
00349 
00350 inline
00351 void DsG4Scintillation::SetScintillationExcitationRatio(const G4double excitationratio)
00352 {
00353         ExcitationRatio = excitationratio;
00354 }
00355 
00356 inline
00357 G4double DsG4Scintillation::GetScintillationExcitationRatio() const
00358 {
00359         return ExcitationRatio;
00360 }
00361 
00362 inline
00363 G4PhysicsTable* DsG4Scintillation::GetSlowIntegralTable() const
00364 {
00365         return theSlowIntegralTable;
00366 }
00367 
00368 inline
00369 G4PhysicsTable* DsG4Scintillation::GetFastIntegralTable() const
00370 {
00371         return theFastIntegralTable;
00372 }
00373 
00374 inline
00375 G4PhysicsTable* DsG4Scintillation::GetReemissionIntegralTable() const
00376 {
00377         return theReemissionIntegralTable;
00378 }
00379 
00380 inline
00381 void DsG4Scintillation::DumpPhysicsTable() const
00382 {
00383         if (theFastIntegralTable) {
00384            G4int PhysicsTableSize = theFastIntegralTable->entries();
00385            G4PhysicsOrderedFreeVector *v;
00386 
00387            for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00388            {
00389                 v = (G4PhysicsOrderedFreeVector*)(*theFastIntegralTable)[i];
00390                 v->DumpValues();
00391            }
00392          }
00393 
00394         if (theSlowIntegralTable) {
00395            G4int PhysicsTableSize = theSlowIntegralTable->entries();
00396            G4PhysicsOrderedFreeVector *v;
00397 
00398            for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00399            {
00400                 v = (G4PhysicsOrderedFreeVector*)(*theSlowIntegralTable)[i];
00401                 v->DumpValues();
00402            }
00403          }
00404 
00405         if (theReemissionIntegralTable) {
00406            G4int PhysicsTableSize = theReemissionIntegralTable->entries();
00407            G4PhysicsOrderedFreeVector *v;
00408 
00409            for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00410            {
00411                 v = (G4PhysicsOrderedFreeVector*)(*theReemissionIntegralTable)[i];
00412                 v->DumpValues();
00413            }
00414          }
00415 }
00416 
00417 #endif /* DsG4Scintillation_h */
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:17:58 for DetSim by doxygen 1.7.4