Changeset 18450 for trunk/Mars/mjobs


Ignore:
Timestamp:
03/08/16 15:57:17 (9 years ago)
Author:
ftemme
Message:
Restructured the whole MJSimulation::Process method and added a lot of comments for a better overview what is done in ceres
Location:
trunk/Mars/mjobs
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/Mars/mjobs/MJSimulation.cc

    r18333 r18450  
    5656#include <TEnv.h>
    5757#include <TCanvas.h>
     58#include <iostream>
    5859
    5960// Core
    … …  
    144145    return kTRUE;
    145146}
    146 /*
    147 TString MJSimulation::GetOutFile(const MSequence &seq) const
    148 {
    149     return seq.IsValid() ? Form("ceres%08d.root", seq.GetSequence()) : "ceres.root";
    150 }
    151 */
    152147
    153148Bool_t MJSimulation::WriteResult(const MParList &plist, const MSequence &seq, Int_t num)
    … …  
    337332}
    338333
     334// Setup the header accordingly to the used operation mode
     335void MJSimulation::SetupHeaderOperationMode(MRawRunHeader &header) const
     336{
     337    switch (fOperationMode)
     338    {
     339    case kModeData:
     340        header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTData);
     341        header.SetRunInfo(0, fRunNumber<0 ? 3 : fRunNumber);
     342        break;
     343
     344    case kModePed:
     345        header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPedestal);
     346        header.SetSourceInfo("Pedestal");
     347        header.SetRunInfo(0, fRunNumber<0 ? 1 : fRunNumber);
     348        break;
     349
     350    case kModeCal:
     351        header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTCalibration);
     352        header.SetSourceInfo("Calibration");
     353        header.SetRunInfo(0, fRunNumber<0 ? 2 : fRunNumber);
     354        break;
     355
     356    case kModePointRun:
     357        header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPointRun);
     358        header.SetRunInfo(0, fRunNumber<0 ? 0 : fRunNumber);
     359        break;
     360    }
     361}
     362
     363void MJSimulation::CreateBinningObjects(MParList &plist) const
     364{
     365    MBinning *binse = (MBinning*) plist.FindCreateObj("MBinning","BinningEnergy");
     366    binse->SetEdges(120,      1,  1000000,"log");
     367    MBinning *binsth = (MBinning*) plist.FindCreateObj("MBinning","BinningThreshold");
     368    binsth->SetEdges( 60,   0.9,   900000, "log");
     369    MBinning *binsee = (MBinning*) plist.FindCreateObj("MBinning","BinningEnergyEst");
     370    binsee->SetEdges( 36,   0.9,   900000, "log");
     371    MBinning *binss = (MBinning*) plist.FindCreateObj("MBinning","BinningSize");
     372    binss->SetEdges( 100,     1, 10000000,      "log");
     373    MBinning *binsi = (MBinning*) plist.FindCreateObj("MBinning","BinningImpact");
     374    binsi->SetEdges(  32,     0,      800);
     375    MBinning *binsh = (MBinning*) plist.FindCreateObj("MBinning","BinningHeight");
     376    binsh->SetEdges( 150,     0,       50);
     377    MBinning *binsaz = (MBinning*) plist.FindCreateObj("MBinning","BinningAz");
     378    binsaz->SetEdges(720,  -360,      360);
     379    MBinning *binszd = (MBinning*) plist.FindCreateObj("MBinning","BinningZd");
     380    binszd->SetEdges( 70,     0,       70);
     381    MBinning *binsvc = (MBinning*) plist.FindCreateObj("MBinning","BinningViewCone");
     382    binsvc->SetEdges( 35,     0,        7);
     383    MBinning *binsel = (MBinning*) plist.FindCreateObj("MBinning","BinningTotLength");
     384    binsel->SetEdges(150,     0,       50);
     385    MBinning *binsew = (MBinning*) plist.FindCreateObj("MBinning","BinningMedLength");
     386    binsew->SetEdges(150,     0,       15);
     387    plist.FindCreateObj("MBinning","BinningDistC");
     388    plist.FindCreateObj("MBinning","BinningDist");
     389    plist.FindCreateObj("MBinning","BinningTrigPos");
     390}
     391
    339392Bool_t MJSimulation::Process(const MArgs &args, const MSequence &seq)
    340393{
     394    // --------------------------------------------------------------------------------
     395    // --------------------------------------------------------------------------------
     396    // Initialization
     397    // --------------------------------------------------------------------------------
     398    // --------------------------------------------------------------------------------
     399
     400    // --------------------------------------------------------------------------------
     401    // Logging output:
     402    // --------------------------------------------------------------------------------
     403    // - description of the job itself
     404    // - list of the to processed sequence
    341405    *fLog << inf;
    342406    fLog->Separator(GetDescriptor());
    … …  
    357421        return kFALSE;
    358422    }
    359 
    360     // --------------------------------------------------------------------------------
    361     // Setup Parlist
     423    // --------------------------------------------------------------------------------
     424    // Setup MParList and MTasklist
     425    // --------------------------------------------------------------------------------
    362426    MParList plist;
    363427    plist.AddToList(this); // take care of fDisplay!
    … …  
    365429    MTaskList tasks;
    366430    plist.AddToList(&tasks);
    367     // --------------------------------------------------------------------------------
    368 
     431
     432    // --------------------------------------------------------------------------------
     433    // --------------------------------------------------------------------------------
     434    // Parameter Container Setup
     435    // --------------------------------------------------------------------------------
     436    // --------------------------------------------------------------------------------
     437
     438    // --------------------------------------------------------------------------------
     439    // Setup container for the reflector, the cones and the camera geometry
     440    // --------------------------------------------------------------------------------
    369441    // FIXME: Allow empty GeomCones!
    370442    MParEnv env0("Reflector");
    … …  
    377449    plist.AddToList(&env1);
    378450    plist.AddToList(&env2);
    379 
     451    // --------------------------------------------------------------------------------
     452    // Setup container for the "camera array" of the extracted number of photons
     453    // --------------------------------------------------------------------------------
     454    // from ExtractorRndm
    380455    plist.FindCreateObj("MPedPhotCam", "MPedPhotFromExtractorRndm");
    381 
     456    // --------------------------------------------------------------------------------
     457    // Setup spline containers for:
     458    // --------------------------------------------------------------------------------
     459    // - the pulse shape
     460    // - the different photon acceptance splines
    382461    MParSpline shape("PulseShape");
    383     plist.AddToList(&shape);
    384 
    385     // *** FIXME *** FIXME *** FIXME ***
    386     plist.FindCreateObj("MMcRunHeader");
    387 
    388     MRawRunHeader header;
    389     header.SetValidMagicNumber();
    390     //header.InitFadcType(3);
    391 
    392     switch (fOperationMode)
    393     {
    394     case kModeData:
    395         header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTData);
    396         header.SetRunInfo(0, fRunNumber<0 ? 3 : fRunNumber);
    397         break;
    398 
    399     case kModePed:
    400         header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPedestal);
    401         header.SetSourceInfo("Pedestal");
    402         header.SetRunInfo(0, fRunNumber<0 ? 1 : fRunNumber);
    403         break;
    404 
    405     case kModeCal:
    406         header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTCalibration);
    407         header.SetSourceInfo("Calibration");
    408         header.SetRunInfo(0, fRunNumber<0 ? 2 : fRunNumber);
    409         break;
    410 
    411     case kModePointRun:
    412         header.SetRunType(MRawRunHeader::kRTMonteCarlo|MRawRunHeader::kRTPointRun);
    413         header.SetRunInfo(0, fRunNumber<0 ? 0 : fRunNumber);
    414         break;
    415     }
    416 
    417     // FIXME: Move to MSimPointingPos, MSimCalibrationSignal
    418     //        Can we use this as input for MSimPointingPos?
    419     header.SetObservation("On", "MonteCarlo");
    420     plist.AddToList(&header);
    421 
    422     MCorsikaRead read;
    423     read.SetForceMode(fForceMode);
    424 
    425     if (!seq.IsValid())
    426     {
    427         for (int i=0; i<args.GetNumArguments(); i++)
    428             read.AddFile(args.GetArgumentStr(i));
    429     }
    430     else
    431         read.AddFiles(iter);
    432 
    433     MContinue precut("", "PreCut");
    434     precut.IsInverted();
    435     precut.SetAllowEmpty();
    436 
    437     MSimMMCS simmmcs;
    438 
    439462    MParSpline splinepde("PhotonDetectionEfficiency");
    440463    MParSpline splinemirror("MirrorReflectivity");
    … …  
    442465    MParSpline splinecones2("ConesTransmission");
    443466    MParSpline splineAdditionalPhotonAcceptance("AdditionalPhotonAcceptance");
     467    plist.AddToList(&shape);
    444468    plist.AddToList(&splinepde);
    445469    plist.AddToList(&splinemirror);
    … …  
    448472    plist.AddToList(&splineAdditionalPhotonAcceptance);
    449473
     474    // --------------------------------------------------------------------------------
     475    // Setup container for the MC run header and the Raw run header
     476    // --------------------------------------------------------------------------------
     477    plist.FindCreateObj("MMcRunHeader");
     478
     479    MRawRunHeader header;
     480    header.SetValidMagicNumber();
     481    SetupHeaderOperationMode(header);
     482    // FIXME: Move to MSimPointingPos, MSimCalibrationSignal
     483    //        Can we use this as input for MSimPointingPos?
     484    header.SetObservation("On", "MonteCarlo");
     485    plist.AddToList(&header);
     486
     487    // --------------------------------------------------------------------------------
     488    // Setup container for the intended pulse position and for the trigger position
     489    // --------------------------------------------------------------------------------
     490    MParameterD pulpos("IntendedPulsePos");
     491    // FIXME: Set a default which could be 1/3 of the digitization window
     492    //    pulpos.SetVal(40);  // [ns]
     493    plist.AddToList(&pulpos);
     494
     495    MParameterD trigger("TriggerPos");
     496    trigger.SetVal(0);
     497    plist.AddToList(&trigger);
     498
     499
     500    // --------------------------------------------------------------------------------
     501    // Setup container for the fix time offset, for residual time spread and for gapd time jitter
     502    // --------------------------------------------------------------------------------
     503    // Dominik and Sebastian on: fix time offsets
     504    MMatrix fix_time_offsets_between_pixels_in_ns(
     505        "MFixTimeOffset","titel"
     506    );
     507    plist.AddToList(&fix_time_offsets_between_pixels_in_ns);
     508
     509    // Jens Buss on: residual time spread
     510    MParameterD resTimeSpread("ResidualTimeSpread");
     511    resTimeSpread.SetVal(0.0);
     512    plist.AddToList(&resTimeSpread);
     513    // Jens Buss on: residual time spread
     514    MParameterD gapdTimeJitter("GapdTimeJitter");
     515    gapdTimeJitter.SetVal(0.0);
     516    plist.AddToList(&gapdTimeJitter);
     517
     518    // --------------------------------------------------------------------------------
     519    // Creation of binning objects and apply default settings
     520    // --------------------------------------------------------------------------------
     521    CreateBinningObjects(plist);
     522
     523    // --------------------------------------------------------------------------------
     524    // --------------------------------------------------------------------------------
     525    // Tasks Setup (Reading, Absorption, Reflector)
     526    // --------------------------------------------------------------------------------
     527    // --------------------------------------------------------------------------------
     528
     529    // --------------------------------------------------------------------------------
     530    // Corsika read setup
     531    // --------------------------------------------------------------------------------
     532    MCorsikaRead read;
     533    read.SetForceMode(fForceMode);
     534
     535    if (!seq.IsValid())
     536    {
     537        for (int i=0; i<args.GetNumArguments(); i++)
     538            read.AddFile(args.GetArgumentStr(i));
     539    }
     540    else
     541        read.AddFiles(iter);
     542
     543    // --------------------------------------------------------------------------------
     544    // Precut
     545    // --------------------------------------------------------------------------------
     546    MContinue precut("", "PreCut");
     547    precut.IsInverted();
     548    precut.SetAllowEmpty();
     549
     550    // --------------------------------------------------------------------------------
     551    // MSimMMCS (don't know what this is doing) and calculation of the incident angle
     552    // --------------------------------------------------------------------------------
     553    MSimMMCS simmmcs;
    450554    const TString sin2 = "sin(MCorsikaEvtHeader.fZd)*sin(MCorsikaRunHeader.fZdMin*TMath::DegToRad())";
    451555    const TString cos2 = "cos(MCorsikaEvtHeader.fZd)*cos(MCorsikaRunHeader.fZdMin*TMath::DegToRad())";
    … …  
    457561    calcangle.SetNameParameter("IncidentAngle");
    458562
     563    // --------------------------------------------------------------------------------
     564    // Absorption tasks
     565    // --------------------------------------------------------------------------------
     566    // - Atmosphere (simatm)
     567    // - photon detection efficiency (absapd)
     568    // - mirror reflectivity (absmir)
     569    // - angular acceptance of winston cones (cones)
     570    // - transmission of winston cones (cones2)
     571    // - additional photon acceptance (only motivated, not measured) (additionalPhotonAcceptance)
    459572    MSimAtmosphere simatm;
    460573    MSimAbsorption absapd("SimPhotonDetectionEfficiency");
    … …  
    470583    additionalPhotonAcceptance.SetParName("AdditionalPhotonAcceptance");
    471584 
     585    // --------------------------------------------------------------------------------
     586    // Simulating pointing position and simulating reflector
     587    // --------------------------------------------------------------------------------
    472588    MSimPointingPos pointing;
    473589
    … …  
    477593//  MSimStarField stars;
    478594
     595    // --------------------------------------------------------------------------------
     596    // Continue tasks
     597    // --------------------------------------------------------------------------------
     598    // - Processing of the current events stops, if there aren't at least two photons (see cont3)
    479599    MContinue cont1("MPhotonEvent.GetNumPhotons<1", "ContEmpty1");
    480600    MContinue cont2("MPhotonEvent.GetNumPhotons<1", "ContEmpty2");
    … …  
    484604    cont3.SetAllowEmpty();
    485605
    486     // -------------------------------------------------------------------
    487 
    488     MBinning binse( 120,     1,  1000000, "BinningEnergy",    "log");
    489     MBinning binsth( 60,   0.9,   900000, "BinningThreshold", "log");
    490     MBinning binsee( 36,   0.9,   900000, "BinningEnergyEst", "log");
    491     MBinning binss( 100,     1, 10000000, "BinningSize",      "log");
    492 //    MBinning binsi( 100,  -500,      500, "BinningImpact");
    493     MBinning binsi(  32,     0,      800, "BinningImpact");
    494     MBinning binsh( 150,     0,       50, "BinningHeight");
    495     MBinning binsaz(720,  -360,      360, "BinningAz");
    496     MBinning binszd( 70,     0,       70, "BinningZd");
    497     MBinning binsvc( 35,     0,        7, "BinningViewCone");
    498     MBinning binsel(150,     0,       50, "BinningTotLength");
    499     MBinning binsew(150,     0,       15, "BinningMedLength");
    500     MBinning binsd("BinningDistC");
    501     MBinning binsd0("BinningDist");
    502     MBinning binstr("BinningTrigPos");
    503 
    504     plist.AddToList(&binsee);
    505     plist.AddToList(&binse);
    506     plist.AddToList(&binsth);
    507     plist.AddToList(&binss);
    508     plist.AddToList(&binsi);
    509     plist.AddToList(&binsh);
    510     plist.AddToList(&binszd);
    511     plist.AddToList(&binsaz);
    512     plist.AddToList(&binsvc);
    513     plist.AddToList(&binsel);
    514     plist.AddToList(&binsew);
    515     plist.AddToList(&binstr);
    516     plist.AddToList(&binsd);
    517     plist.AddToList(&binsd0);
    518 
     606    // --------------------------------------------------------------------------------
     607    // --------------------------------------------------------------------------------
     608    // Tasks Setup (Cones, SiPM-Simulation, RandomPhotons, Camera, Trigger, DAQ)
     609    // --------------------------------------------------------------------------------
     610    // --------------------------------------------------------------------------------
     611
     612    // --------------------------------------------------------------------------------
     613    // GeomApply, SimPSF, SimGeomCam, SimCalibrationSignal
     614    // --------------------------------------------------------------------------------
     615    // - MGeomApply only resizes all MGeomCam parameter container to the actual
     616    //   number of pixels
     617    // - MSimPSF applies an additional smearing of the photons on the camera window
     618    //   on default it is switched of (set to 0), but can be applied by setting:
     619    //   MSimPSF.fSigma: <value>
     620    //   in the config file.
     621    //   Be aware, that there is also a smearing of the photons implemented in
     622    //   MSimReflector by setting:
     623    //   MReflector.SetSigmaPSF: <value>
     624    //   in the config file. This is at the moment done in the default config file.
     625    // - MSimGeomCam identifies which pixel was hit by which photon
     626    // - MSimCalibrationSignal is only used for simulated calibration runs
     627    MGeomApply apply;
     628
     629    MSimPSF simpsf;
     630
     631    MSimGeomCam simgeom;
     632    simgeom.SetNameGeomCam("GeomCones");
     633
     634    MSimCalibrationSignal simcal;
     635    simcal.SetNameGeomCam("GeomCones");
     636
     637    // --------------------------------------------------------------------------------
     638    // Simulation of random photons
     639    // --------------------------------------------------------------------------------
     640    // - be aware that here default values for the nsb rate and for the dark
     641    //   count rate are set. They will be overwritten if the values are defined
     642    //   in the config file
     643    // - if you want to simulate a specific nsb rate you have to set the filename
     644    //   in the config file to an empty string and then you can specify the NSB rate:
     645    //   MSimRandomPhotons.FileNameNSB:
     646    //   MSimRandomPhotons.FrequencyNSB: 0.0
     647
     648    // Sky Quality Meter: 21.82M = 2.02e-4cd/m^2
     649    // 1cd = 12.566 lm / 4pi sr
     650    // FIXME: Simulate photons before cones and QE!
     651    MSimRandomPhotons  simnsb;
     652    simnsb.SetFreq(5.8, 0.004);  // ph/m^2/nm/sr/ns NSB, 4MHz dark count rate
     653    simnsb.SetNameGeomCam("GeomCones");
     654    // FIXME: How do we handle star-light? We need the rate but we also
     655    //        need to process it through the mirror!
     656
     657    // --------------------------------------------------------------------------------
     658    // Simulation from the SiPM to the DAQ
     659    // --------------------------------------------------------------------------------
     660    // - MSimAPD simulates the whole behaviour of the SiPMs:
     661    //   (dead time of cells, crosstalk, afterpulses)
     662    // - MSimExcessNoise adds a spread on the weight of each signal in the SiPMs
     663    // - MSimBundlePhotons restructured the photon list of MPhotonEvent via a look
     664    //   up table. At the moment the tasks is skipped, due to the fact that there
     665    //   is no file name specified in the config file
     666    // - MSimSignalCam only fills a SignalCam container for the cherenkov photons
     667    // - MSimCamera simulates the behaviour of the camera:
     668    //   (electronic noise, accoupling, time smearing and offset for the photons,
     669    //   pulses injected by all photons)
     670    // - MSimTrigger simulates the behaviour of the Trigger:
     671    //   (Adding patch voltages, applying clipping, applying discriminator, applying
     672    //   time over threshold and a possible trigger logic)
     673    // - MContinue conttrig stops the processing of the current event if the trigger
     674    //   position is negativ (and therefore not valid)
     675    // - MSimReadout simulates the behaviour of the readout:
     676    //   (Digitization and saturation of the ADC)
     677    MSimAPD simapd;
     678    simapd.SetNameGeomCam("GeomCones");
     679
     680    MSimExcessNoise   simexcnoise;
     681    MSimBundlePhotons simsum;
     682    MSimSignalCam     simsignal;
     683    MSimCamera        simcam;
     684    MSimTrigger       simtrig;
     685    MContinue conttrig("TriggerPos.fVal<0", "ContTrigger");
     686    MSimReadout       simdaq;
     687
     688    // --------------------------------------------------------------------------------
     689    // Standard analysis with hillas parameters for the true cherenkov photons
     690    // --------------------------------------------------------------------------------
     691    // Remove isolated pixels
     692    MImgCleanStd clean(0, 0);
     693    //clean.SetCleanLvl0(0); // The level above which isolated pixels are kept
     694    clean.SetCleanRings(0);
     695    clean.SetMethod(MImgCleanStd::kAbsolute);
     696
     697    //clean.SetNamePedPhotCam("MPedPhotFromExtractorRndm");
     698
     699    MHillasCalc hcalc;
     700    hcalc.Disable(MHillasCalc::kCalcConc);
     701
     702
     703    // --------------------------------------------------------------------------------
     704    // Setup of histogram containers and the corresponding fill tasks
     705    // --------------------------------------------------------------------------------
     706    // Here is no functionality for the simulation, it is only visualization via
     707    // showplot
    519708    MHn mhn1, mhn2, mhn3, mhn4;
    520709    SetupHist(mhn1);
    … …  
    565754    plist.AddToList(&plane4);
    566755    plist.AddToList(&planeF1);
    567     plist.AddToList(&planeF2);;
     756    plist.AddToList(&planeF2);
    568757
    569758    //MHPSF psf;
    … …  
    589778    //fillP.SetNameTab("PSF",        "Photon distribution after cones for all mirrors");
    590779
    591     // -------------------------------------------------------------------
    592 
     780
     781    MHCamEvent evt0a(/*10*/3, "Signal",    "Average signal (absolute);;S [ph]");
     782    MHCamEvent evt0c(/*10*/3, "MaxSignal", "Maximum signal (absolute);;S [ph]");
     783    MHCamEvent evt0d(/*11*/8, "ArrTm",     "Time after first photon;;T [ns]");
     784    evt0a.SetErrorSpread(kFALSE);
     785    evt0c.SetCollectMax();
     786
     787    MContinue cut("", "Cut");
     788
     789    MFillH fillx0a(&evt0a,             "MSignalCam",      "FillAvgSignal");
     790    MFillH fillx0c(&evt0c,             "MSignalCam",      "FillMaxSignal");
     791    MFillH fillx0d(&evt0d,             "MSignalCam",      "FillArrTm");
     792    MFillH fillx1("MHHillas",          "MHillas",         "FillHillas");
     793    MFillH fillx3("MHHillasSrc",       "MHillasSrc",      "FillHillasSrc");
     794    MFillH fillx4("MHNewImagePar",     "MNewImagePar",    "FillNewImagePar");
     795    MFillH fillth("MHThreshold",       "",                "FillThreshold");
     796    MFillH fillca("MHCollectionArea",  "",                "FillTrigArea");
     797
     798    fillth.SetNameTab("Threshold");
     799    fillca.SetNameTab("TrigArea");
     800
     801    // --------------------------------------------------------------------------------
     802    // Formating of the outputfilepath
     803    // --------------------------------------------------------------------------------
    593804    if (!fFileOut.IsNull())
    594805    {
    … …  
    598809            fFileOut = fFileOut.Remove(dot);
    599810    }
    600 
    601     // -------------------------------------------------------------------
    602 
    603811    const char *fmt = fFileOut.IsNull() ?
    604812        Form("s/cer([0-9]+)([0-9][0-9][0-9])/%s\\/%08d.$2%%s_MonteCarlo$1.root/", Esc(fPathOut).Data(), header.GetRunNumber()) :
    … …  
    610818    TString rule3(Form(fmt, Form("_%c", header.GetRunTypeChar())));
    611819
     820    // --------------------------------------------------------------------------------
     821    // Setup of the outputfile tasks
     822    // --------------------------------------------------------------------------------
    612823    MWriteRootFile write4a( 2, rule4, fOverwrite?"RECREATE":"NEW", "Star file");
    613824    MWriteRootFile write4b( 2, rule4, fOverwrite?"RECREATE":"NEW", "Star file");
    … …  
    682893    write4b.AddContainer("MMcEvtBasic", "OriginalMC");
    683894
    684     // -------------------------------------------------------------------
    685 
    686     MGeomApply apply;
    687 
    688     MSimPSF simpsf;
    689 
    690     MSimGeomCam simgeom;
    691     simgeom.SetNameGeomCam("GeomCones");
    692 
    693     MSimCalibrationSignal simcal;
    694     simcal.SetNameGeomCam("GeomCones");
    695 
    696     // Sky Quality Meter: 21.82M = 2.02e-4cd/m^2
    697     // 1cd = 12.566 lm / 4pi sr
    698 
    699     // FIXME: Simulate photons before cones and QE!
    700 
    701     MSimRandomPhotons  simnsb;
    702     simnsb.SetFreq(5.8, 0.004);  // ph/m^2/nm/sr/ns NSB, 4MHz dark count rate
    703     simnsb.SetNameGeomCam("GeomCones");
    704 
    705     // FIXME: How do we handle star-light? We need the rate but we also
    706     //        need to process it through the mirror!
    707 
    708     MSimAPD simapd;
    709     simapd.SetNameGeomCam("GeomCones");
    710 
    711     MSimExcessNoise   simexcnoise;
    712     MSimBundlePhotons simsum;
    713     MSimSignalCam     simsignal;
    714     MSimCamera        simcam;
    715     MSimTrigger       simtrig;
    716     MSimReadout       simdaq;
    717 
    718     MContinue conttrig("TriggerPos.fVal<0", "ContTrigger");
    719 
    720     MParameterD pulpos("IntendedPulsePos");
    721     // FIXME: Set a default which could be 1/3 of the digitization window
    722     //    pulpos.SetVal(40);  // [ns]
    723     plist.AddToList(&pulpos);
    724 
    725     MParameterD trigger("TriggerPos");
    726     trigger.SetVal(0);
    727     plist.AddToList(&trigger);
    728    
    729     // -------------------------------------------------------------------
    730     // Dominik and Sebastian on: fix time offsets
    731     MMatrix fix_time_offsets_between_pixels_in_ns(
    732         "MFixTimeOffset","titel"
    733     );
    734     plist.AddToList(&fix_time_offsets_between_pixels_in_ns);
    735 
    736     // -------------------------------------------------------------------
    737 
    738     // Jens Buss on: residual time spread
    739     MParameterD resTimeSpread("ResidualTimeSpread");
    740     resTimeSpread.SetVal(0.0);
    741     plist.AddToList(&resTimeSpread);
    742 
    743     // -------------------------------------------------------------------
    744 
    745     // Jens Buss on: residual time spread
    746     MParameterD gapdTimeJitter("GapdTimeJitter");
    747     gapdTimeJitter.SetVal(0.0);
    748     plist.AddToList(&gapdTimeJitter);
    749 
    750     // -------------------------------------------------------------------
    751 
    752     // Remove isolated pixels
    753     MImgCleanStd clean(0, 0);
    754     //clean.SetCleanLvl0(0); // The level above which isolated pixels are kept
    755     clean.SetCleanRings(0);
    756     clean.SetMethod(MImgCleanStd::kAbsolute);
    757 
    758     //clean.SetNamePedPhotCam("MPedPhotFromExtractorRndm");
    759 
    760     MHillasCalc hcalc;
    761     hcalc.Disable(MHillasCalc::kCalcConc);
    762 
    763     MHCamEvent evt0a(/*10*/3, "Signal",    "Average signal (absolute);;S [ph]");
    764     MHCamEvent evt0c(/*10*/3, "MaxSignal", "Maximum signal (absolute);;S [ph]");
    765     MHCamEvent evt0d(/*11*/8, "ArrTm",     "Time after first photon;;T [ns]");
    766     evt0a.SetErrorSpread(kFALSE);
    767     evt0c.SetCollectMax();
    768 
    769     MContinue cut("", "Cut");
    770 
    771     MFillH fillx0a(&evt0a,             "MSignalCam",      "FillAvgSignal");
    772     MFillH fillx0c(&evt0c,             "MSignalCam",      "FillMaxSignal");
    773     MFillH fillx0d(&evt0d,             "MSignalCam",      "FillArrTm");
    774     MFillH fillx1("MHHillas",          "MHillas",         "FillHillas");
    775     MFillH fillx3("MHHillasSrc",       "MHillasSrc",      "FillHillasSrc");
    776     MFillH fillx4("MHNewImagePar",     "MNewImagePar",    "FillNewImagePar");
    777     MFillH fillth("MHThreshold",       "",                "FillThreshold");
    778     MFillH fillca("MHCollectionArea",  "",                "FillTrigArea");
    779 
    780     fillth.SetNameTab("Threshold");
    781     fillca.SetNameTab("TrigArea");
    782 
    783     // -------------------------------------------------------------------
     895
     896    // --------------------------------------------------------------------------------
     897    // --------------------------------------------------------------------------------
     898    // Filling of tasks in tasklist
     899    // --------------------------------------------------------------------------------
     900    // --------------------------------------------------------------------------------
    784901
    785902    if (header.IsDataRun())
    786903    {
    787         tasks.AddToList(&read);
    788         tasks.AddToList(&precut);
    789         tasks.AddToList(&pointing);
    790         tasks.AddToList(&simmmcs);
    791         if (!fPathOut.IsNull() && !HasNullOut())
     904        tasks.AddToList(&read);  // Reading Corsika
     905        tasks.AddToList(&precut);  // Precut
     906        tasks.AddToList(&pointing);  // Simulating pointing
     907        tasks.AddToList(&simmmcs);  // Simulating MMCS
     908        if (!fPathOut.IsNull() && !HasNullOut())  // Write Tasks for corsika infos
    792909        {
    793910            //tasks.AddToList(&write1b);
    … …  
    801918        //        tasks.AddToList(&stars);
    802919        if (1)
    803             tasks.AddToList(&simatm); // Here because before fillh1
    804         tasks.AddToList(&calcangle);
    805         tasks.AddToList(&fillh1);
    806         tasks.AddToList(&fillG);
     920            tasks.AddToList(&simatm); // Here because before fillh1  // atmosphere absorption
     921        tasks.AddToList(&calcangle);  // calculation incident angle
     922        tasks.AddToList(&fillh1);  // fill histogram task
     923        tasks.AddToList(&fillG);  // fill histogram task
    807924        if (!header.IsPointRun())
    808925        {
    809             tasks.AddToList(&absapd);
    810             tasks.AddToList(&absmir);
    811             tasks.AddToList(&additionalPhotonAcceptance);
     926            tasks.AddToList(&absapd);  // photon detection efficiency of the apds
     927            tasks.AddToList(&absmir);  // mirror reflectivity
     928            tasks.AddToList(&additionalPhotonAcceptance);  // addition photon acceptance
    812929            if (0)
    813                 tasks.AddToList(&simatm); // FASTER?
     930                tasks.AddToList(&simatm); // FASTER?  // be aware this is 'commented'
    814931        }
    815         tasks.AddToList(&reflect);
     932        tasks.AddToList(&reflect);  // Simulation of the reflector
    816933        if (!header.IsPointRun())
    817934        {
    818             tasks.AddToList(&fill0);
     935            tasks.AddToList(&fill0);  // fill histogram task
    819936            //tasks.AddToList(&fill1);
    820             tasks.AddToList(&fill2);
    821             tasks.AddToList(&fill3);
    822             tasks.AddToList(&fill4);
    823             tasks.AddToList(&fillF1);
     937            tasks.AddToList(&fill2);  // fill histogram task
     938            tasks.AddToList(&fill3);  // fill histogram task
     939            tasks.AddToList(&fill4);  // fill histogram task
     940            tasks.AddToList(&fillF1);  // fill histogram task
    824941        }
    825         tasks.AddToList(&cones);
    826         tasks.AddToList(&cones2);
     942        tasks.AddToList(&cones);  // angular acceptance of winston cones
     943        tasks.AddToList(&cones2);  // transmission of winston cones
    827944        //if (header.IsPointRun())
    828945        //    tasks.AddToList(&fillP);
    829         tasks.AddToList(&cont1);
     946        tasks.AddToList(&cont1);  // continue if at least 2 photons
    830947    }
    831948    // -------------------------------
    832     tasks.AddToList(&apply);
     949    tasks.AddToList(&apply);  // apply geometry to MGeomCam containers
    833950    if (header.IsDataRun())
    834951    {
    835         tasks.AddToList(&simpsf);
    836         tasks.AddToList(&simgeom);
    837         tasks.AddToList(&cont2);
     952        tasks.AddToList(&simpsf);  // simulates additional psf (NOT used in default mode)
     953        tasks.AddToList(&simgeom);  // tag which pixel is hit by which photon
     954        tasks.AddToList(&cont2);  // continue if at least 2 photons
    838955        if (!header.IsPointRun())
    839956        {
    840             tasks.AddToList(&fillF2);
    841             tasks.AddToList(&fillh2);
     957            tasks.AddToList(&fillF2);  // fill histogram task
     958            tasks.AddToList(&fillh2);  // fill histogram task
    842959        }
    843         tasks.AddToList(&cont3);
     960        tasks.AddToList(&cont3);  // continue if at least 3 photons
    844961    }
    845962    if (fCamera)
    846963    {
    847964        if (header.IsPedestalRun() || header.IsCalibrationRun())
    848             tasks.AddToList(&simcal);
    849         tasks.AddToList(&simnsb);
    850         tasks.AddToList(&simapd);
    851         tasks.AddToList(&simexcnoise);
    852     }
    853     tasks.AddToList(&simsum);
     965            tasks.AddToList(&simcal);  // add calibration signal for calibration runs
     966        tasks.AddToList(&simnsb);  // simulate nsb and dark counts
     967        tasks.AddToList(&simapd);  // simulate SiPM behaviour (dead time, crosstalk ...)
     968        tasks.AddToList(&simexcnoise);  // add excess noise
     969    }
     970    tasks.AddToList(&simsum);  // bundle photons (NOT used in default mode)
    854971    if (fCamera)
    855972    {
    856         tasks.AddToList(&simcam);
     973        tasks.AddToList(&simcam);  // simulate camera behaviour (creates analog signal)
    857974        if (header.IsDataRun() || fForceTrigger)
    858             tasks.AddToList(&simtrig);
    859         tasks.AddToList(&conttrig);
    860         tasks.AddToList(&simdaq);
    861     }
    862     tasks.AddToList(&simsignal);  // What do we do if signal<0?
     975            tasks.AddToList(&simtrig);  // simulate trigger
     976        tasks.AddToList(&conttrig);  // continue if trigger pos is valid
     977        tasks.AddToList(&simdaq);  // simulate data aquisition
     978    }
     979    tasks.AddToList(&simsignal);  // What do we do if signal<0?  // fill MSimSignal container
    863980    if (!fPathOut.IsNull() && !HasNullOut())
    864981    {
    865982        //tasks.AddToList(&write1a);
    866983        if (!header.IsPedestalRun())
    867             tasks.AddToList(&write2a);
     984            tasks.AddToList(&write2a);  // outputtask
    868985        if (fCamera)
    869             tasks.AddToList(&write3a);
     986            tasks.AddToList(&write3a);  // outputtask (this is the raw data writing tasks)
    870987    }
    871988    // -------------------------------
    … …  
    874991        // FIXME: MHCollectionArea Trigger Area!
    875992        if (header.IsDataRun())
    876             tasks.AddToList(&fillh3);
    877         tasks.AddToList(&filltp);
     993            tasks.AddToList(&fillh3);  // fill histogram task
     994        tasks.AddToList(&filltp);  // fill histogram task
    878995    }
    879996    if (header.IsDataRun())
    880997    {
    881         tasks.AddToList(&fillew);
    882         tasks.AddToList(&filled);
     998        tasks.AddToList(&fillew);  // fill histogram task
     999        tasks.AddToList(&filled);  // fill histogram task
    8831000    }
    8841001    if (!header.IsPedestalRun())
    8851002    {
    886         tasks.AddToList(&fillx0a);
    887         tasks.AddToList(&fillx0c);
     1003        tasks.AddToList(&fillx0a);  // fill histogram task
     1004        tasks.AddToList(&fillx0c);  // fill histogram task
    8881005    }
    8891006    if (header.IsDataRun())
    8901007    {
    891         tasks.AddToList(&clean);
    892         tasks.AddToList(&hcalc);
     1008        tasks.AddToList(&clean);  // Cleaning for standard analysis
     1009        tasks.AddToList(&hcalc);  // Hillas parameter calculation
    8931010        tasks.AddToList(&cut);
    894         tasks.AddToList(&fillx0d);
    895         tasks.AddToList(&fillx1);
     1011        tasks.AddToList(&fillx0d);  // fill histogram task
     1012        tasks.AddToList(&fillx1);  // fill histogram task
    8961013        //tasks.AddToList(&fillx2);
    897         tasks.AddToList(&fillx3);
    898         tasks.AddToList(&fillx4);
     1014        tasks.AddToList(&fillx3);  // fill histogram task
     1015        tasks.AddToList(&fillx4);  // fill histogram task
    8991016        if (!HasNullOut())
    9001017            tasks.AddToList(&write4a);
    9011018        //tasks.AddToList(&fillx5);
    9021019
    903         tasks.AddToList(&fillh4);
    904         tasks.AddToList(&fillth);
    905         tasks.AddToList(&fillca);
    906     }
    907     //-------------------------------------------
    908 
     1020        tasks.AddToList(&fillh4);  // fill histogram task
     1021        tasks.AddToList(&fillth);  // fill histogram task
     1022        tasks.AddToList(&fillca);  // fill histogram task
     1023    }
     1024
     1025
     1026    // --------------------------------------------------------------------------------
     1027    // --------------------------------------------------------------------------------
     1028    // Event loop and processing display
     1029    // --------------------------------------------------------------------------------
     1030    // --------------------------------------------------------------------------------
     1031
     1032    // --------------------------------------------------------------------------------
     1033    // Creation of event loop
     1034    // --------------------------------------------------------------------------------
    9091035    tasks.SetAccelerator(MTask::kAccDontReset|MTask::kAccDontTime);
    9101036
    … …  
    9181044
    9191045    // FIXME: If pedestal file given use the values from this file
    920     //-------------------------------------------
     1046
     1047
     1048    // --------------------------------------------------------------------------------
     1049    // Preparation of the graphicalk display which is shown while processing
     1050    // --------------------------------------------------------------------------------
    9211051
    9221052    MGeomCam *cam = static_cast<MGeomCam*>(env2.GetCont());
    9231053
    924     if (binstr.IsDefault())
    925         binstr.SetEdgesLin(150, -shape.GetWidth(),
     1054
     1055    MBinning *binsd  = (MBinning*) plist.FindObject("BinningDistC");
     1056    MBinning *binsd0 = (MBinning*) plist.FindObject("BinningDist");
     1057    MBinning *binstr = (MBinning*) plist.FindObject("BinningTrigPos");
     1058    if (binstr->IsDefault())
     1059        binstr->SetEdgesLin(150, -shape.GetWidth(),
    9261060                           header.GetFreqSampling()/1000.*header.GetNumSamples()+shape.GetWidth());
    9271061
    928     if (binsd.IsDefault() && cam)
    929         binsd.SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());
    930 
    931     if (binsd0.IsDefault() && cam)
    932         binsd0.SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());
     1062    if (binsd->IsDefault() && cam)
     1063        binsd->SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());
     1064
     1065    if (binsd0->IsDefault() && cam)
     1066        binsd0->SetEdgesLin(100, 0, cam->GetMaxRadius()*cam->GetConvMm2Deg());
    9331067
    9341068
    … …  
    10451179    }
    10461180
    1047     //-------------------------------------------
     1181    // --------------------------------------------------------------------------------
     1182    // Perform event loop
     1183    // --------------------------------------------------------------------------------
    10481184
    10491185    // Execute first analysis

  • trunk/Mars/mjobs/MJSimulation.h

    r17866 r18450  
    4040    void SetupHeaderKeys(MWriteFitsFile& write, MRawRunHeader &header) const;
    4141    void SetupVetoColumns(MWriteFitsFile& write) const;
     42    void SetupHeaderOperationMode(MRawRunHeader &header) const;
     43    void CreateBinningObjects(MParList &plist) const;
    4244
    4345public:
Note: See TracChangeset for help on using the changeset viewer.