Changeset 9425 for trunk/MagicSoft/Mars/mjobs/MJSimulation.cc

Timestamp:

04/16/09 12:04:29 (16 years ago)

Author:

tbretz

Message:

*** empty log message ***

File:

• trunk/MagicSoft/Mars/mjobs/MJSimulation.cc (modified) (20 diffs)

trunk/MagicSoft/Mars/mjobs/MJSimulation.cc

@@ -108 +108 @@
 #include "MReflector.h"
 #include "MParEnv.h"
-#include "MPulseShape.h"
+#include "MSpline3.h"
+#include "MParSpline.h"
 #include "MGeomCam.h"
 
@@ -152 +153 @@
     //cont.Add(const_cast<MSequence*>(&fSequence));
 
-    cont.Add(plist.FindObject("MPulseShape"));
+    cont.Add(plist.FindObject("PulseShape"));
     cont.Add(plist.FindObject("Reflector"));
     cont.Add(plist.FindObject("MGeomCam"));
@@ -185 +186 @@
     hist.SetLog(kTRUE, kTRUE, kFALSE);
 
-    hist.AddHist("-MCorsikaEvtHeader.fX/100","-MCorsikaEvtHeader.fY/100");
-    hist.SetDrawOption("colz");
-    hist.InitName("Impact;Impact;Impact");
-    hist.InitTitle("Impact;West <--> East [m];South <--> North [m]");
+    /*
+     hist.AddHist("-MCorsikaEvtHeader.fX/100","-MCorsikaEvtHeader.fY/100");
+     hist.SetDrawOption("colz");
+     hist.InitName("Impact;Impact;Impact");
+     hist.InitTitle("Impact;West <--> East [m];South <--> North [m]");
+     hist.SetAutoRange();
+     */
+
+    hist.AddHist("MCorsikaEvtHeader.GetImpact/100");
+    hist.InitName("Impact");
+    hist.InitTitle("Impact;Impact [m]");
     hist.SetAutoRange();
 
@@ -197 +205 @@
     hist.AddHist("(MCorsikaEvtHeader.fAz+MCorsikaRunHeader.fMagneticFieldAz)*TMath::RadToDeg()", "MCorsikaEvtHeader.fZd*TMath::RadToDeg()");
     hist.InitName("SkyOrigin;Az;Zd");
-    hist.InitTitle("Sky Origin;Az [deg];Zd [deg]");
+    hist.InitTitle("Sky Origin;Az [\\deg];Zd [\\deg]");
     hist.SetDrawOption("colz");
     hist.SetAutoRange();
@@ -279 +287 @@
     plist.FindCreateObj("MPedPhotCam", "MPedPhotFromExtractorRndm");
 
-    MPulseShape shape;
+    MParSpline shape("PulseShape");
     plist.AddToList(&shape);
 
@@ -325 +333 @@
     read.SetForceMode(fForceMode);
 
-    for (int i=0; i<args.GetNumArguments();i ++)
+    for (int i=0; i<args.GetNumArguments(); i++)
         read.AddFile(args.GetArgumentStr(i));
 
     MSimMMCS simmmcs;
 
+    MParSpline splinepde("PhotonDetectionEfficiency");
+    MParSpline splinemirror("MirrorReflectivity");
+    MParSpline splinecones("ConesAngularAcceptance");
+    plist.AddToList(&splinepde);
+    plist.AddToList(&splinemirror);
+    plist.AddToList(&splinecones);
+
     MSimAtmosphere simatm;
-    MSimAbsorption absapd("PhotonDetectionEfficiency");
-    MSimAbsorption absmir("MirrorReflectivity");
-    MSimAbsorption cones("ConesAngularAcceptance");
+    MSimAbsorption absapd("SimPhotonDetectionEfficiency");
+    MSimAbsorption absmir("SimMirrorReflectivity");
+    MSimAbsorption cones("SimConesAngularAcceptance");
+    absapd.SetParName("PhotonDetectionEfficiency");
+    absmir.SetParName("MirrorReflectivity");
+    cones.SetParName("ConesAngularAcceptance");
     cones.SetUseTheta();
 
@@ -344 +362 @@
     MContinue cont1("MPhotonEvent.GetNumPhotons<1", "ContEmpty1");
     MContinue cont2("MPhotonEvent.GetNumPhotons<1", "ContEmpty2");
+    MContinue cont3("MPhotonEvent.GetNumPhotons<2", "ContEmpty3");
 
     // -------------------------------------------------------------------
 
     MBinning binse( 100,     1,   100000, "BinningEnergy",    "log");
-    MBinning binsth( 35,   0.9,   900000, "BinningThreshold", "log");
+    MBinning binsth( 70,   0.9,   900000, "BinningThreshold", "log");
     MBinning binsee( 35,   0.9,   900000, "BinningEnergyEst", "log");
     MBinning binss( 100,     1, 10000000, "BinningSize",      "log");
-    MBinning binsi( 100,  -500,      500, "BinningImpact");
+//    MBinning binsi( 100,  -500,      500, "BinningImpact");
+    MBinning binsi(  55,     0,     1100, "BinningImpact");
     MBinning binsh( 150,     0,       50, "BinningHeight");
     MBinning binsaz(720,  -360,      360, "BinningAz");
@@ -377 +397 @@
     plist.AddToList(&binsd0);
 
-    MHn mhn1, mhn2, mhn3;
+    MHn mhn1, mhn2, mhn3, mhn4;
     SetupHist(mhn1);
     SetupHist(mhn2);
     SetupHist(mhn3);
-
-    MH3 mhtp("TriggerPos.fVal-IntendedPulsePos.fVal-MPulseShape.GetPulseWidth");
+    SetupHist(mhn4);
+
+    MH3 mhtp("TriggerPos.fVal-IntendedPulsePos.fVal-PulseShape.GetWidth");
     mhtp.SetName("TrigPos");
     mhtp.SetTitle("Trigger position w.r.t. the first photon hitting a detector");
@@ -397 +418 @@
     MFillH fillh2(&mhn2, "", "FillH2");
     MFillH fillh3(&mhn3, "", "FillH3");
+    MFillH fillh4(&mhn4, "", "FillH4");
     MFillH filltp(&mhtp, "", "FillTriggerPos");
     MFillH fillew(&mhew, "", "FillEvtWidth");
@@ -403 +425 @@
     fillh2.SetNameTab("Detectable", "Distribution of events hit the detector");
     fillh3.SetNameTab("Triggered",  "Distribution of triggered events");
+    fillh4.SetNameTab("Cleaned",    "Distribution after cleaning and cuts");
     filltp.SetNameTab("TrigPos",    "TriggerPosition w.r.t the first photon");
     fillew.SetNameTab("EvtWidth",   "Time between first and last photon hitting a detector");
@@ -491 +514 @@
     simcal.SetNameGeomCam("GeomCones");
 
-    //  Dark Current: 4MHz = 0.004 GHz
-    //  Light of night sky at La Palma: ~ 0.2 ph / cm^2 / sr / ns
-    //  5deg camera: 6e-3 sr
-    //  NSB = 0.2*6e-3
-
-    // MAGIC: 84MHz       / Pixel
-    // DWARF: 14MHz-25Mhz / APD
-
     // FIXME: Simulate photons before cones and QE!
 
     MSimRandomPhotons  simnsb;
-    simnsb.SetFreq(0.0015, 0.04);  // 1.5MHz dark count rate, 40MHZ/cm^2 NSB
+    simnsb.SetFreq(5.8, 0.004);  // ph/m^2/nm/sr/ns NSB, 4MHz dark count rate
     simnsb.SetNameGeomCam("GeomCones");
 
@@ -509 +524 @@
 
     MSimAPD simapd;
-    //simapd.SetFreq(0.04); // Must be identical to MSimRandomPhotons!!
     simapd.SetNameGeomCam("GeomCones");
 
@@ -555 +569 @@
     MFillH fillx0d(&evt0d,             "MSignalCam",      "FillArrTm");
     MFillH fillx1("MHHillas",          "MHillas",         "FillHillas");
-    //MFillH fillx2("MHHillasExt",       "",                "FillHillasExt");
     MFillH fillx3("MHHillasSrc",       "MHillasSrc",      "FillHillasSrc");
+    MFillH fillx4("MHNewImagePar",     "MNewImagePar",    "FillNewImagePar");
     MFillH fillth("MHThreshold",       "",                "FillThreshold");
     MFillH fillca("MHCollectionArea",  "",                "FillTrigArea");
-    //MFillH fillx4("MHImagePar",        "MImagePar",       "FillImagePar");
-    //MFillH fillx5("MHNewImagePar",     "MNewImagePar",    "FillNewImagePar");
 
     fillth.SetNameTab("Threshold");
@@ -603 +615 @@
         }
         tasks.AddToList(&cones);
-        tasks.AddToList(&fillF2);
         //if (header.IsPointRun())
         //    tasks.AddToList(&fillP);
         tasks.AddToList(&cont1);
-        if (!header.IsPointRun())
-            tasks.AddToList(&fillh2);
     }
     // -------------------------------
@@ -617 +626 @@
         tasks.AddToList(&simgeom);
         tasks.AddToList(&cont2);
+        if (!header.IsPointRun())
+        {
+            tasks.AddToList(&fillF2);
+            tasks.AddToList(&fillh2);
+        }
+        tasks.AddToList(&cont3);
     }
     if (fCamera)
@@ -669 +684 @@
         //tasks.AddToList(&fillx2);
         tasks.AddToList(&fillx3);
-        //tasks.AddToList(&fillx4);
+        tasks.AddToList(&fillx4);
         //tasks.AddToList(&fillx5);
     }
     if (header.IsDataRun())
     {
+        tasks.AddToList(&fillh4);
         tasks.AddToList(&fillth);
         tasks.AddToList(&fillca);
@@ -695 +711 @@
 
     if (binstr.IsDefault())
-        binstr.SetEdgesLin(150, -shape.GetPulseWidth(),
-                           header.GetFreqSampling()/1000.*header.GetNumSamples()+shape.GetPulseWidth());
+        binstr.SetEdgesLin(150, -shape.GetWidth(),
+                           header.GetFreqSampling()/1000.*header.GetNumSamples()+shape.GetWidth());
 
     if (binsd.IsDefault() && cam)
@@ -731 +747 @@
             gROOT->SetSelectedPad(0);
             static_cast<MReflector*>(env0.GetCont())->DrawClone("line")->SetBit(kCanDelete);
-
         }
 
@@ -749 +764 @@
                 c->SetBit(kCanDelete);
                 c->Draw();
-
             }
 
@@ -767 +781 @@
             }
         }
+
+        TCanvas &d = fDisplay->AddTab("Info2");
+        d.Divide(2,2);
+
+        d.cd(1);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridx();
+        gPad->SetGridy();
+        gROOT->SetSelectedPad(0);
+        splinepde.DrawClone()->SetBit(kCanDelete);
+
+        d.cd(2);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridx();
+        gPad->SetGridy();
+        gROOT->SetSelectedPad(0);
+        splinemirror.DrawClone()->SetBit(kCanDelete);
+
+        d.cd(3);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridx();
+        gPad->SetGridy();
+        gROOT->SetSelectedPad(0);
+        splinecones.DrawClone()->SetBit(kCanDelete);
+
+        d.cd(4);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridx();
+        gPad->SetGridy();
+        gROOT->SetSelectedPad(0);
+        MParSpline *all = (MParSpline*)splinepde.DrawClone();
+        //all->SetTitle("Combined acceptance");
+        all->SetBit(kCanDelete);
+        all->Multiply(*splinemirror.GetSpline());
     }