Index: trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.cc =================================================================== --- trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.cc (revision 2097) +++ trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.cc (revision 2107) @@ -17,4 +17,5 @@ ! ! Author(s): Wolfgang Wittek 4/2002 +! Abelardo Moralejo 5/2003 ! ! Copyright: MAGIC Software Development, 2000-2002 @@ -45,4 +46,5 @@ #include "MLog.h" #include "MLogManip.h" +#include ClassImp(MHMcEnergyMigration); @@ -54,5 +56,5 @@ // MHMcEnergyMigration::MHMcEnergyMigration(const char *name, const char *title) - : fHist() + : fHist(), fHist2(), fHistImp() { // @@ -63,9 +65,20 @@ fHist.SetDirectory(NULL); - - fHist.SetTitle("3D-plot E-true E-est Theta"); - fHist.SetXTitle("E_{true} [GeV]"); - fHist.SetYTitle("E_{est} [GeV]"); + fHist.SetTitle("3D-plot E_{TRUE} E_{EST} \\Theta"); + fHist.SetXTitle("E_{TRUE} [GeV]"); + fHist.SetYTitle("E_{EST} [GeV]"); fHist.SetZTitle("\\Theta [\\circ]"); + + fHist2.SetDirectory(NULL); + fHist2.SetTitle("3D-plot \\Delta E / E vs E_{TRUE} E_{EST}"); + fHist2.SetXTitle("E_{TRUE} [GeV]"); + fHist2.SetYTitle("E_{EST} [GeV]"); + fHist2.SetZTitle("\\frac{E_{EST} - E_{TRUE}}{E_{TRUE}}"); + + fHistImp.SetDirectory(NULL); + fHistImp.SetTitle("\\Delta E / E vs Impact parameter"); + fHistImp.SetXTitle("Impact parameter (m)"); + fHistImp.SetYTitle("\\frac{E_{EST} - E_{TRUE}}{E_{TRUE}}"); + } @@ -99,6 +112,14 @@ SetBinning(&fHist, binsenergy, binsenergy, binstheta); - fHist.Sumw2(); + + const MBinning* binsde = (MBinning*)plist->FindObject("BinningDE"); + const MBinning* binsimpact = (MBinning*)plist->FindObject("BinningImpact"); + SetBinning(&fHistImp, binsimpact, binsde); + + fHistImp.Sumw2(); + + SetBinning(&fHist2, binsenergy, binsenergy, binsde); + fHist2.Sumw2(); return kTRUE; @@ -124,8 +145,9 @@ gPad->SetLogx(); h = fHist.Project3D("ex_pro"); - h->SetTitle("Distribution of E_{true}"); - h->SetXTitle("E_{true} [GeV]"); + h->SetTitle("Distribution of E_{TRUE}"); + h->SetXTitle("E_{TRUE} [GeV]"); h->SetBit(kCanDelete); h->Draw(opt); + h->SetDirectory(NULL); pad->cd(2); @@ -133,23 +155,133 @@ gPad->SetLogx(); h = fHist.Project3D("ey_pro"); - h->SetTitle("Distribution of E_{est}"); + h->SetTitle("Distribution of E_{EST}"); h->SetXTitle("E_{est} [GeV]"); h->SetBit(kCanDelete); h->Draw(opt); + Double_t minEest = h->GetXaxis()->GetXmin(); + h->SetDirectory(NULL); pad->cd(3); gPad->SetBorderMode(0); - h = fHist.Project3D("ez_pro"); + h = fHist.Project3D("z_pro"); h->SetTitle("Distribution of \\Theta"); h->SetXTitle("\\Theta [\\circ]"); h->SetBit(kCanDelete); + h->SetLineWidth(2); h->Draw(opt); + h->SetDirectory(NULL); pad->cd(4); gPad->SetBorderMode(0); - fHist.Draw(opt); + TH1D* hpx; + hpx = fHistImp.ProjectionX("_px", 1, fHistImp.GetNbinsY(),"e"); + hpx->SetTitle("Distribution of Impact parameter"); + hpx->SetXTitle("Impact parameter (m)"); + hpx->SetBit(kCanDelete); + hpx->Draw(opt); + hpx->SetDirectory(NULL); + fHistImp.SetDirectory(NULL); pad->Modified(); pad->Update(); + + TVirtualPad *pad2 = MakeDefCanvas((TString)GetName()+"2"); + pad2->SetBorderMode(0); + + AppendPad(""); + + pad2->Divide(2,2); + + TH2D *h2; + + pad2->cd(1); + gPad->SetBorderMode(0); + gPad->SetLogx(); + gPad->SetLogy(); + h2 = (TH2D*) fHist.Project3D("xy"); + + TProfile* h2pfx; + h2pfx = h2->ProfileX("_pfx", 1, h2->GetNbinsY(),"S"); + h2pfx->SetXTitle("E_{TRUE} (GeV)"); + h2pfx->SetYTitle("E_{EST} (GeV)"); + h2pfx->SetTitle("E_{EST} vs E_{TRUE}"); + h2pfx->SetBit(kCanDelete); + h2pfx->SetFillColor(41); + h2pfx->SetMinimum(minEest); + h2pfx->SetStats(kFALSE); + h2pfx->Draw("E3"); + + h2->SetBit(kCanDelete); + h2->SetFillColor(1); + h2->Draw("box,same"); + h2->SetDirectory(NULL); + + h2pfx->SetLineColor(2); + h2pfx->SetLineWidth(2); + h2pfx->Draw("L,histo,same"); + h2pfx->SetDirectory(NULL); + + pad2->cd(2); + gPad->SetBorderMode(0); + gPad->SetLogx(); + gPad->SetLeftMargin(0.15); + h2 = (TH2D*) fHist2.Project3D("zx"); + h2->SetBit(kCanDelete); + h2pfx = h2->ProfileX("_pfx", 1, h2->GetNbinsY(),"S"); + h2pfx->SetTitle("\\Delta E / E vs E_{TRUE}"); + h2pfx->SetXTitle("E_{TRUE} (GeV)"); + h2pfx->SetYTitle("\\frac{E_{EST} - E_{TRUE}}{E_{TRUE}}"); + h2pfx->SetBit(kCanDelete); + h2pfx->SetFillColor(41); + h2pfx->GetYaxis()->SetTitleOffset(1.4); + h2pfx->SetStats(kFALSE); + h2pfx->Draw("E3"); + h2->SetFillColor(1); + h2->Draw("same,box"); + h2->SetDirectory(NULL); + h2pfx->SetLineColor(2); + h2pfx->SetLineWidth(2); + h2pfx->Draw("L,histo,same"); + h2pfx->SetDirectory(NULL); + + pad2->cd(3); + gPad->SetBorderMode(0); + gPad->SetLogx(); + gPad->SetLeftMargin(0.15); + h2 = (TH2D*) fHist2.Project3D("zy"); + h2->SetBit(kCanDelete); + h2pfx = h2->ProfileX("_pfx", 1, h2->GetNbinsY(),"S"); + h2pfx->SetTitle("\\Delta E / E vs E_{EST}"); + h2pfx->SetXTitle("E_{TRUE} (GeV)"); + h2pfx->SetYTitle("\\frac{E_{EST} - E_{TRUE}}{E_{TRUE}}"); + h2pfx->SetBit(kCanDelete); + h2pfx->SetFillColor(41); + h2pfx->GetYaxis()->SetTitleOffset(1.4); + h2pfx->SetStats(kFALSE); + h2pfx->SetMinimum(-1.); + h2pfx->SetMaximum(1.); + h2pfx->Draw("E3"); + + h2->SetFillColor(1); + h2->Draw("same,box"); + h2->SetDirectory(NULL); + h2pfx->SetLineColor(2); + h2pfx->SetLineWidth(2); + h2pfx->Draw("L,histo,same"); + h2pfx->SetDirectory(NULL); + + pad2->cd(4); + gPad->SetBorderMode(0); + h = fHist2.ProjectionZ("_pz",1,fHist2.GetNbinsX(),1,fHist2.GetNbinsY(),"e"); + h->SetBit(kCanDelete); + h->Draw(); + h->SetDirectory(NULL); + + pad2->Modified(); + pad2->Update(); + + fHist.SetDirectory(NULL); + fHist2.SetDirectory(NULL); + } @@ -160,9 +292,11 @@ Bool_t MHMcEnergyMigration::Fill(const MParContainer *par, const Stat_t w) { - // MHillasSrc &hil = *(MHillasSrc*)par; - // fHist.Fill(fMcEvt->GetEnergy(), hil.GetSize()); - // get E-true from fMcEvt and E-est from fEnergy - fHist.Fill(fMcEvt->GetEnergy(), fEnergy->GetEnergy(), fMcEvt->GetTheta()*kRad2Deg); + + fHist.Fill(fMcEvt->GetEnergy(), fEnergy->GetEnergy(), fMcEvt->GetTelescopeTheta()*kRad2Deg, w); + + fHist2.Fill(fMcEvt->GetEnergy(), fEnergy->GetEnergy(), (fEnergy->GetEnergy()-fMcEvt->GetEnergy())/fMcEvt->GetEnergy(), w); + + fHistImp.Fill(fMcEvt->GetImpact()/100., (fEnergy->GetEnergy()-fMcEvt->GetEnergy())/fMcEvt->GetEnergy(), w); return kTRUE; Index: trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.h =================================================================== --- trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.h (revision 2097) +++ trunk/MagicSoft/Mars/mhistmc/MHMcEnergyMigration.h (revision 2107) @@ -8,4 +8,8 @@ #ifndef ROOT_TH3 #include +#endif + +#ifndef ROOT_TH2 +#include #endif @@ -21,4 +25,6 @@ TH3D fHist; + TH3D fHist2; + TH2D fHistImp; public: @@ -30,6 +36,8 @@ const TH3D *GetHist() { return &fHist; } const TH3D *GetHist() const { return &fHist; } + TH1 *GetHistByName(const TString name) { return &fHist; } - TH1 *GetHistByName(const TString name) { return &fHist; } + const TH3D *GetHist2() { return &fHist2; } + const TH2D *GetHistImp() { return &fHistImp; } void Draw(Option_t *option=""); @@ -40,5 +48,2 @@ #endif - - -