Ignore:
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:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/Mars

  • trunk/Mars/msimcamera/MSimCamera.cc

    r10167 r17148  
    4343#include <TF1.h>
    4444#include <TRandom.h>
     45#include <TRandom3.h>
    4546
    4647#include "MLog.h"
     
    6364#include "MMcEvt.hxx"            // To be replaced by a CheObs class
    6465#include "MRawRunHeader.h"
     66
     67#include "math.h"
    6568
    6669ClassImp(MSimCamera);
     
    7578    : fEvt(0), fStat(0), fRunHeader(0), fElectronicNoise(0), fGain(0),
    7679    fCamera(0), fMcEvt(0), fSpline(0), fBaselineGain(kFALSE),
    77     fDefaultOffset(-1), fDefaultNoise(-1), fDefaultGain(-1)
     80      fDefaultOffset(-1), fDefaultNoise(-1), fDefaultGain(-1), fACFudgeFactor(0),
     81      fACTimeConstant(0)
    7882
    7983{
     
    135139        return kFALSE;
    136140
     141    fAccidentalPhotons = (MPedestalCam*)pList->FindObject("AccidentalPhotonRates","MPedestalCam");
     142    if(!fAccidentalPhotons)
     143    {
     144        *fLog << err << "AccidentalPhotonRates [MPedestalCam] not found... aborting." << endl;
     145        return kFALSE;
     146    }
     147
    137148    MParSpline *pulse = (MParSpline*)pList->FindObject("PulseShape", "MParSpline");
    138149    if (!pulse)
     
    142153    }
    143154
    144     if (fRunHeader->GetFreqSampling()!=1000)
    145     {
    146         *fLog << err  << "ERROR - Sampling frequencies others than 1GHz are not yet supported." << endl;
    147         *fLog << warn << "FIXME - SCALE MPulsShape WITH THE SAMPLING FREQUENCY." << endl;
    148         return kFALSE;
    149     }
     155//    if (fRunHeader->GetFreqSampling()!=1000)
     156//    {
     157//        *fLog << err  << "ERROR - Sampling frequencies others than 1GHz are not yet supported." << endl;
     158//        *fLog << warn << "FIXME - SCALE MPulsShape WITH THE SAMPLING FREQUENCY." << endl;
     159//        return kFALSE;
     160//    }
    150161
    151162    fSpline = pulse->GetSpline();
     
    238249        const Double_t rms = pix.GetPedestalRms();
    239250
     251        // FTemme: Implementation of AC-coupling:
     252        // to calculate the value of the accoupling per slice I use the
     253        // following equation:
     254        // accouplingPerSlice = accidentalPhotonRate * (1 + crossTalkProb)
     255        //       * areaOfOnePulse / samplingRate;
     256        // Therefore I need the following variables
     257        Double_t accidentalPhotonRate   = 0; // [MHz]
     258        Float_t crossTalkProb           = 0; // [1]
     259        Double_t areaOfOnePulse         = 0; // [ADC-Counts * s]
     260        Double_t samplingRate           = 0; // [slices * MHz]
     261
     262        // The accidental photon rate is stored in GHz, so we have to multiply
     263        // with 1E3 to get MHz:
     264        MPedestalPix &accPhoPix         = (*fAccidentalPhotons)[i];
     265        accidentalPhotonRate            = accPhoPix.GetPedestal() * 1E3;
     266        Double_t currentAccidentalPhotonRate = accidentalPhotonRate;
     267
     268        if(fACTimeConstant!=0)
     269        {
     270            Double_t accidentalPhotons      = fACTimeConstant * accidentalPhotonRate;
     271            Double_t sigmaAccidentalPhotons = sqrt(accidentalPhotons);
     272
     273            TRandom3 *random = new TRandom3(0);
     274            Double_t gaus = random->Gaus(accidentalPhotons,sigmaAccidentalPhotons);
     275            currentAccidentalPhotonRate = gaus / fACTimeConstant;
     276        }
     277
     278        // I don't know how to get the variable fCrosstalkProb from
     279        // the class APD (see MAvalanchePhotoDiode.h), because there is no
     280        // getter for the APD array(fAPDs) in MSimAPD.
     281        // So I set the crossTalkProb hardcoded to the value 0.15, which is
     282        // equal to the value of the apd of type 4
     283        crossTalkProb                   = 0.15;
     284
     285        // To get the area of one Pulse, I only need to calculate the Integral
     286        // of the Pulse Shape, which is stored in fSpline. Because the spline is
     287        // normalized to a maximal amplitude of 1.0, I had to multiply it with
     288        // the Default gain:
     289        areaOfOnePulse                  = fSpline->Integral() * fDefaultGain;
     290
     291        // The sampling rate I get from the RunHeader:
     292        samplingRate                    = fRunHeader->GetFreqSampling();
     293
     294        Double_t accouplingPerSlice = currentAccidentalPhotonRate
     295                * (1 + crossTalkProb + fACFudgeFactor)
     296                * areaOfOnePulse / samplingRate;
     297
     298        // The accoupling is substracted from the timeline by decreasing the
     299        // mean of the gaussian noise which is added
     300
    240301        if (!fBaselineGain)
    241302        {
    242             (*fCamera)[i].AddGaussianNoise(rms, val);
     303            (*fCamera)[i].AddGaussianNoise(rms, val - accouplingPerSlice);
    243304            continue;
    244305        }
     
    250311        // FIXME: How stable is the offset?
    251312        // FIXME: Should we write a container AppliedGain for MSImTrigger?
    252         (*fCamera)[i].AddGaussianNoise(rms*gain, val*gain);
     313
     314        (*fCamera)[i].AddGaussianNoise(rms*gain, (val - accouplingPerSlice)*gain);
    253315    }
    254316
     
    323385        fDefaultGain = GetEnvValue(env, prefix, "DefaultGain", fDefaultGain);
    324386    }
     387    if (IsEnvDefined(env, prefix, "fACFudgeFactor", print))
     388    {
     389        rc = kTRUE;
     390        fACFudgeFactor = GetEnvValue(env, prefix, "fACFudgeFactor", fACFudgeFactor);
     391    }
     392    if (IsEnvDefined(env, prefix, "fACTimeConstant", print))
     393    {
     394        rc = kTRUE;
     395        fACTimeConstant = GetEnvValue(env, prefix, "fACTimeConstant", fACTimeConstant);
     396    }
    325397
    326398    return rc;
  • trunk/Mars/msimcamera/MSimCamera.h

    r9625 r17148  
    2424    MPedestalCam      *fElectronicNoise; //! Electronic noise (baseline and rms)
    2525    MPedestalCam      *fGain;            //! Electronic noise (baseline and rms)
     26    MPedestalCam      *fAccidentalPhotons;//! Accidental photon rates
    2627
    2728    MAnalogChannels   *fCamera;          //! Output of the analog signals
    2829    MMcEvt            *fMcEvt;           //! For information stored in MMcEvt
    2930
    30     const MSpline3    *fSpline;          // Pulse Shape
     31    const MSpline3    *fSpline;          // Pulse Shape   
    3132
    3233    Bool_t fBaselineGain;  // Should the gain be applied to baseline and electronic noise?
     
    3536    Double_t fDefaultNoise;   // Default noise (RMS) added to all samples
    3637    Double_t fDefaultGain;    // Default gain (multiplication factor to the given/used pulse shape)
     38
     39    Double_t fACFudgeFactor;
     40
     41    Double_t fACTimeConstant;
    3742
    3843    // MParContainer
  • trunk/Mars/msimcamera/MSimReadout.cc

    r9422 r17148  
    145145    fData->InitRead(fRunHeader);
    146146    fData->ResetPixels();
     147    fData->InitStartCells();
    147148    fData->SetIndices();
    148149
     
    202203
    203204    // Check if the whole requested signal can be digitized
    204     if (trig+nslices>=nsamp)
     205    if (trig+nslices>nsamp)
    205206    {
    206207        *fLog << err << "ERROR - Trigger position beyond valid analog signal range." << endl;
Note: See TracChangeset for help on using the changeset viewer.