Changeset 8935

Timestamp:

06/11/08 12:32:36 (17 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars

Files:

• Changelog (modified) (2 diffs)
• NEWS (modified) (2 diffs)
• mhflux/MHEnergyEst.cc (modified) (16 diffs)
• mhflux/MHEnergyEst.h (modified) (3 diffs)

trunk/MagicSoft/Mars/Changelog

@@ -20 +20 @@
 
 
+ 2008/06/11 Thomas Bretz
+
+   * mhflux/MHEnergyEst.[h,cc]:
+     - finally replaced fResolution by more correct histograms
+     - some code cleanup in projecting, profiling and drawing
+     - increased class version number by one
+
+   * mbase/MStatusDisplay.cc:
+     - remove the embedded canvas from the global list to prevent
+       global access to it
+
+
+
  2008/06/10 Thomas Bretz
 
@@ -41 +54 @@
      - replaced some AddContainer by the new AddTree
      - added Pyrometer information to in- and output, respectively
+
+   * datacenter/macros/fillstar.C:
+     - added new columns fAvgHumidity, fAvgCloudiness, fRmsCloudiness
+       and fAvgTempSky
+
+   * mhist/MHWeather.[h,cc]:
+     - removed the display of the solar radiation which was 
+       never working
+     - added display of the pyrometer data to the display
+     - reorganized display
+
+   * mjobs/MJStar.cc:
+     - added filling of the weather data also from the pyrometer branch
 
 

trunk/MagicSoft/Mars/NEWS

@@ -61 +61 @@
    * The code has been prepared for compilation with root 5.18/00d
 
+   * The MHEnergyEst histogram now shows the distribution of 
+     (Eest-Emc)/Est and the distributions (Eest-Emc)/Eest vs. Eest
+     and (Eest-Emc)/Emc vs Emc.
+
  ;merpp
 
@@ -119 +123 @@
     * The effective on-time calculation doesn't use events with only
       sum-trigger anymore
+
+    * The data in the MHWeather tab has been reorganized. The never 
+      working solar radiation has been removed and the data from
+      the pyrometer (cloudiness, air and sky temperature) is 
+      displayed in addition.
 
  ;ganymed/sponde

trunk/MagicSoft/Mars/mhflux/MHEnergyEst.cc

@@ -20 +20 @@
 !   Author(s): Thomas Bretz 1/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
 !
-!   Copyright: MAGIC Software Development, 2000-2005
+!   Copyright: MAGIC Software Development, 2000-2008
 !
 !
@@ -31 +31 @@
 //  calculates the migration matrix E-est vs. E-true
 //  for different bins in Theta
+//
+//  Class Version 2:
+//  - fHResolution
+//  + fHResolutionEst
+//  + fHResolutionMC
 //
 //////////////////////////////////////////////////////////////////////////////
@@ -82 +87 @@
     fHEnergy.SetZTitle("\\Theta [\\circ]");
 
-    fHResolution.SetDirectory(NULL);
-    fHResolution.SetName("EnergyRes");
-    fHResolution.SetTitle("Histogram in \\Delta E/E vs E_{est} and E_{mc}");
-    fHResolution.SetXTitle("E_{est} [GeV]");
-    fHResolution.SetYTitle("E_{mc} [GeV]");
-    fHResolution.SetZTitle("E_{est}/E_{mc}-1");
+    fHResolutionEst.SetDirectory(NULL);
+    fHResolutionEst.SetName("ResEnergyEst");
+    fHResolutionEst.SetTitle("Histogram in \\Delta E/E_{est} vs E_{est}");
+    fHResolutionEst.SetXTitle("E_{est} [GeV]");
+    fHResolutionEst.SetYTitle("1-E_{mc}/E_{est}");
+
+    fHResolutionMC.SetDirectory(NULL);
+    fHResolutionMC.SetName("ResEnergyMC");
+    fHResolutionMC.SetTitle("Histogram in \\Delta E/E_{mc} vs E_{mc}");
+    fHResolutionMC.SetXTitle("E_{mc} [GeV]");
+    fHResolutionMC.SetYTitle("E_{est}/E_{mc}-1");
 
     fHImpact.SetDirectory(NULL);
@@ -105 +115 @@
     SetBinning(&fHEnergy,     &binse, &binse, &binst);
     SetBinning(&fHImpact,     &binsi, &binsr);
-    SetBinning(&fHResolution, &binse, &binse, &binsr);
+
+    SetBinning(&fHResolutionEst, &binse, &binsr);
+    SetBinning(&fHResolutionMC,  &binse, &binsr);
 
     // For some unknown reasons this must be called after
     // the binning has been initialized at least once
     fHEnergy.Sumw2();
-    fHResolution.Sumw2();
     fHImpact.Sumw2();
+    fHResolutionEst.Sumw2();
+    fHResolutionMC.Sumw2();
 }
 
@@ -154 +167 @@
     SetBinning(&fHEnergy,     &binse, &binse, &binst);
     SetBinning(&fHImpact,     &binsi, &binsr);
-    SetBinning(&fHResolution, &binse, &binse, &binsr);
+
+    SetBinning(&fHResolutionEst, &binse, &binsr);
+    SetBinning(&fHResolutionMC,  &binse, &binsr);
 
     fChisq = 0;
@@ -161 +176 @@
     fHEnergy.Reset();
     fHImpact.Reset();
-    fHResolution.Reset();
+
+    fHResolutionEst.Reset();
+    fHResolutionMC.Reset();
 
     return kTRUE;
@@ -176 +193 @@
     const Double_t imp   = fMatrix ? GetVal(1) : fMcEvt->GetImpact()/100;
     const Double_t theta = fMatrix ? GetVal(2) : fMcEvt->GetTelescopeTheta()*TMath::RadToDeg();
-    const Double_t resE  = (eest-etru)/etru;
+
+    const Double_t resEst  = (eest-etru)/eest;
+    const Double_t resMC   = (eest-etru)/etru;
 
     fHEnergy.Fill(eest, etru, theta, w);
-    fHResolution.Fill(eest, etru, resE, w);
-    fHImpact.Fill(imp, resE, w);
+    fHImpact.Fill(imp, resEst, w);
+
+    fHResolutionEst.Fill(eest, resEst, w);
+    fHResolutionMC.Fill(etru, resMC, w);
 
     // For the fit we use a different quantity
@@ -223 +244 @@
     }
 
-    TH1D *h = (TH1D*)fHResolution.ProjectionZ("Resolution");
+    TH1D *h = (TH1D*)fHResolutionEst.ProjectionY("Dummy", -1, -1, "s");
     h->Fit("gaus", "Q0", "", -1.0, 0.25);
 
@@ -270 +291 @@
     pad->cd(1);
 
-    TH1D *hx=0;
-    TH1D *hy=0;
+    TH1 *hx=0;
+    TH1 *hy=0;
 
     if (pad->GetPad(1))
@@ -277 +298 @@
         pad->GetPad(1)->cd(1);
 
-        if ((hx=(TH1D*)gPad->FindObject("EnergyEst_ex")))
-        {
-            TH1D *h2 = (TH1D*)fHEnergy.Project3D("dum_ex");
-            hx->Reset();
-            hx->Add(h2);
-            delete h2;
-        }
-
-        if ((hy=(TH1D*)gPad->FindObject("EnergyEst_ey")))
-        {
-            TH1D *h2 = (TH1D*)fHEnergy.Project3D("dum_ey");
-            hy->Reset();
-            hy->Add(h2);
-            delete h2;
-        }
+        if (gPad->FindObject("EnergyEst_ex"))
+            hx = fHEnergy.Project3D("ex");
+
+        if (gPad->FindObject("EnergyEst_ey"))
+            hy = fHEnergy.Project3D("ey");
 
         if (hx && hy)
         {
+            hx->SetLineColor(kBlue);
+            hx->SetMarkerColor(kBlue);
             hy->SetMaximum();
             hy->SetMaximum(TMath::Max(hx->GetMaximum(), hy->GetMaximum())*1.2);
@@ -305 +318 @@
             pad->GetPad(1)->GetPad(2)->cd(1);
             if ((hx=(TH1D*)gPad->FindObject("EnergyEst_ez")))
-            {
-                TH1D *h2 = (TH1D*)fHEnergy.Project3D("dum_ez");
-                hx->Reset();
-                hx->Add(h2);
-                delete h2;
-            }
+                fHEnergy.Project3D("ez");
 
             pad->GetPad(1)->GetPad(2)->cd(2);
-            hx = (TH1D*)fHResolution.ProjectionZ("Resolution");
+            hx = (TH1D*)fHResolutionEst.ProjectionY("Resolution", -1, -1, "e");
             TPaveStats *stats = dynamic_cast<TPaveStats*>(hx->FindObject("stats"));
             if (stats)
@@ -323 +331 @@
 
             hx->Fit("gaus", "Q", "", -1.0, 0.25);
+            hx->GetFunction("gaus")->SetLineColor(kBlue);
+            hx->GetFunction("gaus")->SetLineWidth(2);
             gPad=NULL;
             gStyle->SetOptFit(101);
@@ -331 +341 @@
     {
         pad->GetPad(2)->cd(1);
-        UpdatePlot(fHEnergy, "yx", kTRUE);
+        if (gPad->FindObject("EnergyEst_yx"))
+        {
+            TH2D *hyx = static_cast<TH2D*>(fHEnergy.Project3D("yx"));
+            UpdateProf(*hyx, kTRUE);
+        }
 
         TLine *l = (TLine*)gPad->FindObject("TLine");
@@ -346 +360 @@
 
         pad->GetPad(2)->cd(2);
-        UpdatePlot(fHResolution, "zy");
+        UpdateProf(fHResolutionEst, kFALSE);
 
         pad->GetPad(2)->cd(3);
-        UpdatePlot(fHResolution, "zx");
+        UpdateProf(fHResolutionMC, kFALSE);
     }
 }
 
-void MHEnergyEst::UpdatePlot(TH3 &h, const char *how, Bool_t logy)
-{
-    TH2D *hyx=0;
-    if (!(hyx=(TH2D*)gPad->FindObject(MString::Format("%s_%s", h.GetName(), how))))
+void MHEnergyEst::UpdateProf(TH2 &h, Bool_t logy)
+{
+    const TString pname = MString::Format("Prof%s", h.GetName());
+
+    if (!gPad->FindObject(pname))
         return;
 
-    TH2D *h2 = (TH2D*)h.Project3D(MString::Format("dum_%s", how));
-    hyx->Reset();
-    hyx->Add(h2);
-    delete h2;
-
-    TH1D *hx = 0;
-    if ((hx=(TH1D*)gPad->FindObject(MString::Format("Prof%s", h.GetName()))))
-    {
-        hx = hyx->ProfileX(Form("Prof%s", h.GetName()), -1, -1, "s");
-
-        if (logy && hx->GetMaximum()>0)
-            gPad->SetLogy();
-    }
-}
-
-TH1 *MHEnergyEst::MakePlot(TH3 &h, const char *how)
-{
-    gPad->SetBorderMode(0);
+    TH1D *hx = h.ProfileX(pname, -1, -1, "s");
+    hx->SetLineColor(kBlue);
+    hx->SetMarkerColor(kBlue);
+
+    if (logy && hx->GetMaximum()>0)
+        gPad->SetLogy();
+}
+
+TH1 *MHEnergyEst::MakeProj(const char *how)
+{
+    TH1 *p = fHEnergy.Project3D(how);
+    p->SetDirectory(NULL);
+    p->SetBit(kCanDelete);
+    p->SetBit(TH1::kNoStats);
+    p->SetMarkerStyle(kFullDotMedium);
+    p->SetLineColor(kBlue);
+
+    return p;
+}
+
+TH1 *MHEnergyEst::MakeProf(TH2 &h)
+{
+    TH1 *p = h.ProfileX(Form("Prof%s", h.GetName()), -1, -1, "s");
+    p->SetDirectory(NULL);
+    p->SetBit(kCanDelete);
+    p->SetLineWidth(2);
+    p->SetLineColor(kBlue);
+    p->SetFillStyle(4000);
+    p->SetStats(kFALSE);
+
+    return p;
+}
+
+// --------------------------------------------------------------------------
+//
+// Draw the histogram
+//
+void MHEnergyEst::Draw(Option_t *opt)
+{
+    TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
+
+    // Do the projection before painting the histograms into
+    // the individual pads
+    AppendPad("");
+
+    pad->SetBorderMode(0);
+    pad->Divide(2, 1, 1e-10, 1e-10);
+
+    TH1 *h;
+
+    // ----------------------------------------
+
+    pad->cd(1);
+    gPad->SetBorderMode(0);
+
+    gPad->Divide(1, 2, 1e-10, 1e-10);
+
+    TVirtualPad *pad2 = gPad;
+
+    // ----------------------------------------
+
+    pad2->cd(1);
+    gPad->SetBorderMode(0);
+
+    gPad->SetGridx();
+    gPad->SetGridy();
+    gPad->SetLogx();
+
+    h = MakeProj("ey");
+    h->SetTitle("Energy disribution: Monte Carlo E_{mc} (black), Estimated E_{est} (blue)");
+    h->SetXTitle("E [GeV]"); // E_mc
+    h->SetYTitle("Counts");
+    h->Draw();
+
+    h = MakeProj("ex");
+    h->SetLineColor(kBlue);
+    h->SetMarkerColor(kBlue);
+    h->Draw("same");
+
+    // ----------------------------------------
+
+    pad2->cd(2);
+    gPad->SetBorderMode(0);
+
+    TVirtualPad *pad3 = gPad;
+    pad3->Divide(2, 1, 1e-10, 1e-10);
+    pad3->cd(1);
+    gPad->SetBorderMode(0);
+    gPad->SetGridx();
+    gPad->SetGridy();
+
+    h = MakeProj("ez");
+    h->SetTitle("Zenith Angle Distribution");
+    h->GetXaxis()->SetMoreLogLabels();
+    h->GetXaxis()->SetNoExponent();
+    h->Draw();
+
+    // ----------------------------------------
+
+    pad3->cd(2);
+    gPad->SetBorderMode(0);
+    gPad->SetGridx();
+    gPad->SetGridy();
+
+    h = fHResolutionEst.ProjectionY("_py");
+    h->SetTitle("Distribution of \\Delta E/E_{est}");
+    h->SetDirectory(NULL);
+    h->SetBit(kCanDelete);
+    h->GetXaxis()->SetRangeUser(-1.3, 1.3);
+    h->Draw();
+    // ----------------------------------------
+
+    pad->cd(2);
+    gPad->SetBorderMode(0);
+
+    gPad->Divide(1, 3, 1e-10, 1e-10);
+    pad2 = gPad;
+
+    // ----------------------------------------
+
+    pad2->cd(1);
+    gPad->SetBorderMode(0);
+    gPad->SetLogy();
     gPad->SetLogx();
     gPad->SetGridx();
@@ -384 +504 @@
     //gROOT->GetListOfCleanups()->Add(gPad); // WHY?
 
-    TH2D *h2 = (TH2D*)h.Project3D(how);
+    TH2D *h2 = (TH2D*)fHEnergy.Project3D("yx");
     h2->SetDirectory(NULL);
     h2->SetBit(kCanDelete);
     h2->SetFillColor(kBlue);
-    h2->SetLineColor(kRed);
-
-    TH1D *h1 = h2->ProfileX(Form("Prof%s", h.GetName()), -1, -1, "s");
-    h1->SetDirectory(NULL);
-    h1->SetBit(kCanDelete);
-    h1->SetLineWidth(2);
-    h1->SetLineColor(kBlue);
-    h1->SetFillStyle(4000);
-    h1->SetStats(kFALSE);
+    h2->SetTitle("Estimated energy E_{est} vs Monte Carlo energy E_{mc}");
+
+    TH1 *h1 = MakeProf(*h2);
 
     h2->Draw("");
-    h1->Draw("E0 hist C same");
-//    h1->Draw("Chistsame");
-
-    return h2;
-}
-
-
-// --------------------------------------------------------------------------
-//
-// Draw the histogram
-//
-void MHEnergyEst::Draw(Option_t *opt)
-{
-    TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
-
-    // Do the projection before painting the histograms into
-    // the individual pads
-    AppendPad("");
-
-    pad->SetBorderMode(0);
-    pad->Divide(2, 1, 1e-10, 1e-10);
-
-    TH1 *h;
-
-    pad->cd(1);
-    gPad->SetBorderMode(0);
-
-    gPad->Divide(1, 2, 1e-10, 1e-10);
-
-    TVirtualPad *pad2 = gPad;
-
-    pad2->cd(1);
-    gPad->SetBorderMode(0);
-
-    gPad->SetGridx();
-    gPad->SetGridy();
-    gPad->SetLogx();
-    h = (TH1D*)fHEnergy.Project3D("ey");
-    h->SetBit(TH1::kNoStats);
-    h->SetTitle("Energy disribution: Monte Carlo E_{mc} (black), Estimated E_{est} (blue)");
-    h->SetXTitle("E [GeV]"); // E_mc
-    h->SetYTitle("Counts");
-    h->SetBit(kCanDelete);
-    h->SetDirectory(NULL);
-    h->SetMarkerStyle(kFullDotMedium);
-    h->Draw();
-
-    h = (TH1D*)fHEnergy.Project3D("ex");
-    h->SetBit(TH1::kNoTitle|TH1::kNoStats);
-    h->SetXTitle("E [GeV]"); // E_est
-    h->SetYTitle("Counts");
-    h->SetBit(kCanDelete);
-    h->SetDirectory(NULL);
-    h->SetMarkerStyle(kFullDotMedium);
-    h->SetLineColor(kBlue);
-    h->SetMarkerColor(kBlue);
-    h->Draw("same");
-
-    // FIXME: LEGEND
-
-    pad2->cd(2);
-    gPad->SetBorderMode(0);
-
-    TVirtualPad *pad3 = gPad;
-    pad3->Divide(2, 1, 1e-10, 1e-10);
-    pad3->cd(1);
-    gPad->SetBorderMode(0);
-    gPad->SetGridx();
-    gPad->SetGridy();
-    h = fHEnergy.Project3D("ez");
-    h->SetTitle("Zenith Angle Distribution");
-    h->SetBit(TH1::kNoStats);
-    h->SetDirectory(NULL);
-    h->SetXTitle("\\Theta [\\circ]");
-    h->SetBit(kCanDelete);
-    h->Draw();
-
-    pad3->cd(2);
-    gPad->SetBorderMode(0);
-    gPad->SetGridx();
-    gPad->SetGridy();
-    /*
-    h = fHImpact.ProjectionX("Impact", -1, -1, "e");
-    h->SetBit(TH1::kNoStats);
-    h->SetTitle("Distribution of Impact");
-    h->SetDirectory(NULL);
-    h->SetXTitle("Impact [m]");
-    h->SetBit(kCanDelete);
-    h->Draw();*/
-    // ----------------------------------------
-    h = fHResolution.ProjectionZ("Resolution");
-    //h->SetXTitle("\\frac{lg(E_{est}) - lg(E_{mc})}{lg(E_{mc})}");
-    h->SetYTitle("Counts");
-    h->SetTitle("Distribution of \\Delta E/E");
-    h->SetDirectory(NULL);
-    h->SetBit(kCanDelete);
-    h->GetXaxis()->SetRangeUser(-1.3, 1.3);
-    h->Draw("");
-    //h->Fit("gaus");
-    // ----------------------------------------
-
-    pad->cd(2);
-    gPad->SetBorderMode(0);
-
-    gPad->Divide(1, 3, 1e-10, 1e-10);
-    pad2 = gPad;
-
-    pad2->cd(1);
-    gPad->SetLogy();
-    h = MakePlot(fHEnergy, "xy");
-    h->SetXTitle("E_{mc} [GeV]");
-    h->SetYTitle("E_{est} [GeV]");
-    h->SetTitle("Estimated energy E_{est} vs Monte Carlo energy E_{mc}");
+    h1->Draw("E0 same");
 
     TLine line;
@@ -520 +522 @@
     line.SetLineStyle(kDashed);
 
+    // ----------------------------------------
+
     pad2->cd(2);
-    h = MakePlot(fHResolution, "zy");
-    h->SetXTitle("E_{mc} [GeV]");
-    h->SetYTitle("(E_{est}/E_{mc}-1");
-    h->SetTitle("Energy resolution \\Delta E/E vs Monte Carlo Energy E_{mc}");
-    h->SetMinimum(-1.3);
-    h->SetMaximum(1.3);
-
-    line.DrawLine(h->GetXaxis()->GetXmin(), 0, h->GetXaxis()->GetXmax(), 0);
+    gPad->SetBorderMode(0);
+    gPad->SetLogx();
+    gPad->SetGridx();
+    gPad->SetGridy();
+    fHResolutionEst.Draw();
+    MakeProf(fHResolutionEst)->Draw("E0 same");
+
+    fHResolutionEst.GetXaxis()->SetMoreLogLabels();
+    fHResolutionEst.GetXaxis()->SetNoExponent();
+
+    line.DrawLine(fHResolutionEst.GetXaxis()->GetXmin(), 0, fHResolutionEst.GetXaxis()->GetXmax(), 0);
+
+    // ----------------------------------------
 
     pad2->cd(3);
-    h = MakePlot(fHResolution, "zx");
-    h->SetXTitle("E_{est} [GeV]");
-    h->SetYTitle("(E_{est}/E_{mc}-1");
-    h->SetTitle("Energy resolution \\Delta E/E vs estimated Energy E_{est}");
-    h->SetMinimum(-1.3);
-    h->SetMaximum(1.3);
-
-    line.DrawLine(h->GetXaxis()->GetXmin(), 0, h->GetXaxis()->GetXmax(), 0);
+    gPad->SetBorderMode(0);
+    gPad->SetLogx();
+    gPad->SetGridx();
+    gPad->SetGridy();
+    fHResolutionMC.Draw();
+    MakeProf(fHResolutionMC)->Draw("E0 same");
+
+    fHResolutionMC.GetXaxis()->SetMoreLogLabels();
+    fHResolutionMC.GetXaxis()->SetNoExponent();
+
+    line.DrawLine(fHResolutionMC.GetXaxis()->GetXmin(), 0, fHResolutionMC.GetXaxis()->GetXmax(), 0);
 }
 

trunk/MagicSoft/Mars/mhflux/MHEnergyEst.h

@@ -33 +33 @@
 
     TH3D fHEnergy;
-    TH3D fHResolution;
+    TH2D fHResolutionEst;
+    TH2D fHResolutionMC;
     TH2D fHImpact;
 
@@ -39 +40 @@
     Double_t fBias;
 
-    TH1 *MakePlot(TH3 &h, const char *how);
-    void UpdatePlot(TH3 &h, const char *how, Bool_t logy=kFALSE);
+    TH1 *MakeProj(const char *how);
+    TH1 *MakeProf(TH2 &h);
+    void UpdateProf(TH2 &h, Bool_t logy);
 
     Double_t GetVal(Int_t i) const;
@@ -62 +64 @@
     void Print(Option_t *o="") const;
 
-    //ClassDef(MHEnergyEst, 2) //
-    ClassDef(MHEnergyEst, 1) // Histogram for the result of the energy reconstruction
+    ClassDef(MHEnergyEst, 2) // Histogram for the result of the energy reconstruction
 };