Changeset 10093 for trunk/Mars/msimcamera

Timestamp:

01/07/11 15:48:35 (14 years ago)

Author:

tbretz

Message:

Added afterpulse treatment to MAvalanchePhotoDiode and MSimAPD

File:

• trunk/Mars/msimcamera/MSimAPD.cc (modified) (9 diffs)

trunk/Mars/msimcamera/MSimAPD.cc

@@ -37 +37 @@
 //
 // Remark:
-//   - The photons MUST be sorted increasing in time.
 //   - The photon rate used to initialize the APD must match the one used
 //     to "fill" the random photons. (FIXME: This should be stored somewhere)
@@ -152 +151 @@
     //   * Make the arguments a data member of MSimAPD
 
-    Int_t   ncells    = 0;
-    Float_t crosstalk = 0;
-    Float_t deadtime  = 0;
-    Float_t recovery  = 0;
+    Int_t   ncells     = 0;
+    Float_t crosstalk  = 0;
+    Float_t deadtime   = 0;
+    Float_t recovery   = 0;
+    Float_t afterprob1 = 0;
+    Float_t afterprob2 = 0;
 
     switch (fType)
     {
     case 1:
-        ncells    = 30;
-        crosstalk = 0.2;
-        deadtime  = 3;
-        recovery  = 8.75*4;
+        ncells     = 30;
+        crosstalk  = 0.2;
+        deadtime   = 3;
+        recovery   = 8.75*4;
+        afterprob1 = 0;
+        afterprob2 = 0;
         break;
 
     case 2:
-        ncells    = 60;
-        crosstalk = 0.2;
-        deadtime  = 3;
-        recovery  = 8.75;
+        ncells     = 60;
+        crosstalk  = 0.2;
+        deadtime   = 3;
+        recovery   = 8.75;
+        afterprob1 = 0;
+        afterprob2 = 0;
         break;
 
     case 3:
-        ncells    = 60;
-        crosstalk = 0.15;
-        deadtime  = 3;
-        recovery  = 8.75;
+        ncells     = 60;
+        crosstalk  = 0.15;
+        deadtime   = 3;
+        recovery   = 8.75;
+        afterprob1 = 0;
+        afterprob2 = 0;
+        break;
+
+    case 4:
+        ncells     = 60;
+        crosstalk  = 0.15;
+        deadtime   = 3;
+        recovery   = 8.75;
+        afterprob1 = 0.14;
+        afterprob2 = 0.11;
         break;
 
@@ -186 +202 @@
 
     for (UInt_t i=0; i<fGeom->GetNumPixels(); i++)
-        fAPDs.Add(new APD(ncells, crosstalk, deadtime, recovery));
+    {
+        APD *apd = new APD(ncells, crosstalk, deadtime, recovery);
+        apd->SetAfterpulseProb(afterprob1, afterprob2);
+
+        fAPDs.Add(apd);
+    }
 
     return kTRUE;
@@ -198 +219 @@
 Int_t MSimAPD::Process()
 {
-    //const Double_t rate = 40e9;
-    // FIXME: Where do we get this number from??
-    // const Double_t rate = 0.04;
-
     // Make all APDs look neutral for the first hit by a photon according to the
     // average hit rate
@@ -213 +230 @@
         return kERROR;
     }
-/*
-    for (UInt_t idx=0; idx<npix; idx++)
-    {
-        const Double_t freq = fRates ? (*fRates)[idx].GetPedestal() : fFreq;
-        static_cast<APD*>(fAPDs.UncheckedAt(idx))->FillRandom(freq, fStat->GetTimeFirst());
-    }
-*/
+
+    // To ensure that the photons are always sorted when calling
+    // HitRandomCellRelative the array is sorted here. If it is sorted
+    // already nothing will be done since the status is stored.
+    // FIXME: Check that this is true and check that it is really necessary
+    fEvt->Sort();
 
     // This tries to initialize dead and relaxing cells properly. If
@@ -228 +244 @@
     {
         const Double_t freq = fRates ? (*fRates)[idx].GetPedestal() : fFreq;
+
+        // Init creates an empty G-APD, i.e. without external pulses
+        // but the correct history for afterpulses and relaxation.
+        // If it was already initialized once it evolves the G-APD
+        // for a time until the effect of relaxation and afterpulses
+        // is below 0.1%. The also creates the possible afterpulses
+        // of the future and deletes later afterpulses from the list.
+        // After the the time stamp fTime is set to 0.
         static_cast<APD*>(fAPDs.UncheckedAt(idx))->Init(freq);
     }
@@ -240 +264 @@
         MPhotonData &ph = (*fEvt)[i];
 
-        // Get arrival time of photon and idx
+        // Get arrival time of photon wrt to left edge of window and its index
         const Double_t t = ph.GetTime()-fStat->GetTimeFirst();
         const Int_t  idx = ph.GetTag();
@@ -255 +279 @@
             return kERROR;
         }
-        // Simulate hitting the APD (the signal height in effective
-        // "number of photons" is returned)
+
+        // Simulate hitting the APD at a time t after T0 (APD::fTime).
+        // Crosstalk is taken into account and the resulting signal height
+        // in effective "number of photons" is returned. Afterpulses until
+        // this time "hit" the G-APD and newly created afterpulses
+        // are stored in the list of afterpulses
         const Double_t hits = static_cast<APD*>(fAPDs.UncheckedAt(idx))->HitRandomCellRelative(t);
-
-        // FIXME: Make a proper simulation of the excess noise!!!
-        //const Double_t signal = gRandom->Gaus(hits, 0.2*TMath::Sqrt(hits));
 
         // Set the weight to the input
@@ -269 +294 @@
     // simulated time.
     for (UInt_t idx=0; idx<npix; idx++)
-        static_cast<APD*>(fAPDs.UncheckedAt(idx))->IncreaseTime(fStat->GetTimeLast());
+    {
+        APD *a = static_cast<APD*>(fAPDs.UncheckedAt(idx));
+
+        const Double_t end = fStat->GetTimeLast()-fStat->GetTimeFirst();
+
+        // This moves T0 (APD::fTime) at the right edge of our time-
+        // window. For the next event this means that afterpulses of past
+        // noise and afterpulses will be available already.
+        // FIXME: Note, that this might mean that a cosmics-pulse
+        //        might increase the noise above the normal level.
+        a->IncreaseTime(end);
+
+        // Get the afterpulses and add them to the signal
+        TIter Next(&a->GetListOfAfterpulses());
+        Afterpulse *ap = 0;
+        while ((ap=static_cast<Afterpulse*>(Next())))
+        {
+            // Skip afterpulses later than that which have been
+            // already produced
+            if (ap->GetTime()>=end)
+                continue;
+
+            // Add a new photon
+            // FIXME: SLOW!
+            MPhotonData &ph = fEvt->Add();
+
+            // Set source to Artificial (noise), the amplitude produced by the
+            // afterpulse (includes crosstalk), its arrival time
+            // and its amplitude, as well as the index of the channel it
+            // corresponds to.
+            ph.SetPrimary(MMcEvtBasic::kArtificial);
+            ph.SetWeight(ap->GetAmplitude());
+            ph.SetTime(ap->GetTime()+fStat->GetTimeFirst());
+            ph.SetTag(idx);
+        }
+
+        // It seems to make sense to delete the previous afterpulses now
+        // but this is not necessary. We have to loop over them in any
+        // case. So we omit the additional loop for deletion but instead
+        // do the deletion in the next loop at the end of Init()
+        // If needed this could be used to keep the total memory
+        // consumption slighly lower.
+    }
+
+    // Now the newly added afterpulses have to be sorted into the array correctly
+    fEvt->Sort();
 
     return kTRUE;