Changeset 1315 for trunk/source/processes/electromagnetic/xrays/src/G4Scintillation.cc

Timestamp:

Jun 18, 2010, 11:42:07 AM (14 years ago)

Author:

garnier

Message:

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08

File:

• trunk/source/processes/electromagnetic/xrays/src/G4Scintillation.cc (modified) (10 diffs)

trunk/source/processes/electromagnetic/xrays/src/G4Scintillation.cc

@@ -25 +25 @@
 //
 //
-// $Id: G4Scintillation.cc,v 1.30 2008/10/22 01:19:11 gum Exp $
-// GEANT4 tag $Name: geant4-09-03 $
+// $Id: G4Scintillation.cc,v 1.31 2010/05/27 20:49:40 gum Exp $
+// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
 //
 ////////////////////////////////////////////////////////////////////////
@@ -37 +37 @@
 // Created:     1998-11-07  
 // Author:      Peter Gumplinger
-// Updated:     2005-08-17 by Peter Gumplinger
+// Updated:     2010-92-22 by Peter Gumplinger
+//              > scintillation rise time included, thanks to
+//              > Martin Goettlich/DESY
+//              2005-08-17 by Peter Gumplinger
 //              > change variable name MeanNumPhotons -> MeanNumberOfPhotons
 //              2005-07-28 by Peter Gumplinger
@@ -93 +96 @@
 
         fTrackSecondariesFirst = false;
+        fFiniteRiseTime = false;
 
         YieldFactor = 1.0;
@@ -248 +252 @@
 
             G4double ScintillationTime = 0.*ns;
+            G4double ScintillationRiseTime = 0.*ns;
             G4PhysicsOrderedFreeVector* ScintillationIntegral = NULL;
 
@@ -255 +260 @@
                    ScintillationTime   = aMaterialPropertiesTable->
                                            GetConstProperty("FASTTIMECONSTANT");
+                   if (fFiniteRiseTime) {
+                      ScintillationRiseTime = aMaterialPropertiesTable->
+                                  GetConstProperty("FASTSCINTILLATIONRISETIME");
+                   }
                    ScintillationIntegral =
                    (G4PhysicsOrderedFreeVector*)((*theFastIntegralTable)(materialIndex));
@@ -261 +270 @@
                    ScintillationTime   = aMaterialPropertiesTable->
                                            GetConstProperty("SLOWTIMECONSTANT");
+                   if (fFiniteRiseTime) {
+                      ScintillationRiseTime = aMaterialPropertiesTable->
+                                  GetConstProperty("SLOWSCINTILLATIONRISETIME");                   }
                    ScintillationIntegral =
                    (G4PhysicsOrderedFreeVector*)((*theSlowIntegralTable)(materialIndex));
@@ -276 +288 @@
                  ScintillationTime   = aMaterialPropertiesTable->
                                           GetConstProperty("FASTTIMECONSTANT");
+                 if (fFiniteRiseTime) {
+                      ScintillationRiseTime = aMaterialPropertiesTable->
+                                 GetConstProperty("FASTSCINTILLATIONRISETIME");
+                 }
                  ScintillationIntegral =
                   (G4PhysicsOrderedFreeVector*)((*theFastIntegralTable)(materialIndex));
@@ -284 +300 @@
                ScintillationTime   =   aMaterialPropertiesTable->
                                           GetConstProperty("SLOWTIMECONSTANT");
+               if (fFiniteRiseTime) {
+                    ScintillationRiseTime = aMaterialPropertiesTable->
+                                 GetConstProperty("SLOWSCINTILLATIONRISETIME");
+               }
                ScintillationIntegral =
                   (G4PhysicsOrderedFreeVector*)((*theSlowIntegralTable)(materialIndex));
@@ -369 +389 @@
                          pPostStepPoint->GetVelocity())/2.);
 
-                deltaTime = deltaTime - 
-                            ScintillationTime * log( G4UniformRand() );
+                // emission time distribution
+                if (ScintillationRiseTime==0.0) {
+                   deltaTime = deltaTime - 
+                          ScintillationTime * log( G4UniformRand() );
+                } else {
+                   deltaTime = deltaTime +
+                          sample_time(ScintillationRiseTime, ScintillationTime);
+                }
 
                 G4double aSecondaryTime = t0 + deltaTime;
@@ -592 +618 @@
 
 }
+
+G4double G4Scintillation::sample_time(G4double tau1, G4double tau2)
+{
+// tau1: rise time and tau2: decay time
+
+        while(1) {
+          // two random numbers
+          G4double ran1 = G4UniformRand();
+          G4double ran2 = G4UniformRand();
+          //
+          // exponential distribution as envelope function: very efficient
+          //
+          G4double d = (tau1+tau2)/tau2;
+          // make sure the envelope function is 
+          // always larger than the bi-exponential
+          G4double t = -1.0*tau2*log(1-ran1);
+          G4double g = d*single_exp(t,tau2);
+          if (ran2 <= bi_exp(t,tau1,tau2)/g) return t;
+        }
+        return -1.0;
+}