Changeset 17148 for trunk/Mars/msimcamera/MSimCamera.cc

Timestamp:

09/09/13 12:39:18 (11 years ago)

Author:

ftemme

Message:

Merging changes from the MC branch in the trunk

Location:

trunk/Mars

Files:

• . (modified) (1 prop)
• msimcamera/MSimCamera.cc (modified) (8 diffs)

trunk/Mars/msimcamera/MSimCamera.cc

@@ -43 +43 @@
 #include <TF1.h>
 #include <TRandom.h>
+#include <TRandom3.h>
 
 #include "MLog.h"
@@ -63 +64 @@
 #include "MMcEvt.hxx"            // To be replaced by a CheObs class
 #include "MRawRunHeader.h"
+
+#include "math.h"
 
 ClassImp(MSimCamera);
@@ -75 +78 @@
     : fEvt(0), fStat(0), fRunHeader(0), fElectronicNoise(0), fGain(0),
     fCamera(0), fMcEvt(0), fSpline(0), fBaselineGain(kFALSE),
-    fDefaultOffset(-1), fDefaultNoise(-1), fDefaultGain(-1)
+      fDefaultOffset(-1), fDefaultNoise(-1), fDefaultGain(-1), fACFudgeFactor(0),
+      fACTimeConstant(0)
 
 {
@@ -135 +139 @@
         return kFALSE;
 
+    fAccidentalPhotons = (MPedestalCam*)pList->FindObject("AccidentalPhotonRates","MPedestalCam");
+    if(!fAccidentalPhotons)
+    {
+        *fLog << err << "AccidentalPhotonRates [MPedestalCam] not found... aborting." << endl;
+        return kFALSE;
+    }
+
     MParSpline *pulse = (MParSpline*)pList->FindObject("PulseShape", "MParSpline");
     if (!pulse)
@@ -142 +153 @@
     }
 
-    if (fRunHeader->GetFreqSampling()!=1000)
-    {
-        *fLog << err  << "ERROR - Sampling frequencies others than 1GHz are not yet supported." << endl;
-        *fLog << warn << "FIXME - SCALE MPulsShape WITH THE SAMPLING FREQUENCY." << endl;
-        return kFALSE;
-    }
+//    if (fRunHeader->GetFreqSampling()!=1000)
+//    {
+//        *fLog << err  << "ERROR - Sampling frequencies others than 1GHz are not yet supported." << endl;
+//        *fLog << warn << "FIXME - SCALE MPulsShape WITH THE SAMPLING FREQUENCY." << endl;
+//        return kFALSE;
+//    }
 
     fSpline = pulse->GetSpline();
@@ -238 +249 @@
         const Double_t rms = pix.GetPedestalRms();
 
+        // FTemme: Implementation of AC-coupling:
+        // to calculate the value of the accoupling per slice I use the
+        // following equation:
+        // accouplingPerSlice = accidentalPhotonRate * (1 + crossTalkProb)
+        //       * areaOfOnePulse / samplingRate;
+        // Therefore I need the following variables
+        Double_t accidentalPhotonRate   = 0; // [MHz]
+        Float_t crossTalkProb           = 0; // [1]
+        Double_t areaOfOnePulse         = 0; // [ADC-Counts * s]
+        Double_t samplingRate           = 0; // [slices * MHz]
+
+        // The accidental photon rate is stored in GHz, so we have to multiply
+        // with 1E3 to get MHz:
+        MPedestalPix &accPhoPix         = (*fAccidentalPhotons)[i];
+        accidentalPhotonRate            = accPhoPix.GetPedestal() * 1E3;
+        Double_t currentAccidentalPhotonRate = accidentalPhotonRate;
+
+        if(fACTimeConstant!=0)
+        {
+            Double_t accidentalPhotons      = fACTimeConstant * accidentalPhotonRate;
+            Double_t sigmaAccidentalPhotons = sqrt(accidentalPhotons);
+
+            TRandom3 *random = new TRandom3(0);
+            Double_t gaus = random->Gaus(accidentalPhotons,sigmaAccidentalPhotons);
+            currentAccidentalPhotonRate = gaus / fACTimeConstant;
+        }
+
+        // I don't know how to get the variable fCrosstalkProb from
+        // the class APD (see MAvalanchePhotoDiode.h), because there is no
+        // getter for the APD array(fAPDs) in MSimAPD.
+        // So I set the crossTalkProb hardcoded to the value 0.15, which is
+        // equal to the value of the apd of type 4
+        crossTalkProb                   = 0.15;
+
+        // To get the area of one Pulse, I only need to calculate the Integral
+        // of the Pulse Shape, which is stored in fSpline. Because the spline is
+        // normalized to a maximal amplitude of 1.0, I had to multiply it with
+        // the Default gain:
+        areaOfOnePulse                  = fSpline->Integral() * fDefaultGain;
+
+        // The sampling rate I get from the RunHeader:
+        samplingRate                    = fRunHeader->GetFreqSampling();
+
+        Double_t accouplingPerSlice = currentAccidentalPhotonRate
+                * (1 + crossTalkProb + fACFudgeFactor)
+                * areaOfOnePulse / samplingRate;
+
+        // The accoupling is substracted from the timeline by decreasing the
+        // mean of the gaussian noise which is added
+
         if (!fBaselineGain)
         {
-            (*fCamera)[i].AddGaussianNoise(rms, val);
+            (*fCamera)[i].AddGaussianNoise(rms, val - accouplingPerSlice);
             continue;
         }
@@ -250 +311 @@
         // FIXME: How stable is the offset?
         // FIXME: Should we write a container AppliedGain for MSImTrigger?
-        (*fCamera)[i].AddGaussianNoise(rms*gain, val*gain);
+
+        (*fCamera)[i].AddGaussianNoise(rms*gain, (val - accouplingPerSlice)*gain);
     }
 
@@ -323 +385 @@
         fDefaultGain = GetEnvValue(env, prefix, "DefaultGain", fDefaultGain);
     }
+    if (IsEnvDefined(env, prefix, "fACFudgeFactor", print))
+    {
+        rc = kTRUE;
+        fACFudgeFactor = GetEnvValue(env, prefix, "fACFudgeFactor", fACFudgeFactor);
+    }
+    if (IsEnvDefined(env, prefix, "fACTimeConstant", print))
+    {
+        rc = kTRUE;
+        fACTimeConstant = GetEnvValue(env, prefix, "fACTimeConstant", fACTimeConstant);
+    }
 
     return rc;