Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 3852) +++ /trunk/MagicSoft/Mars/Changelog (revision 3853) @@ -36,4 +36,14 @@ 2004/04/27: Markus Gaug + + * mjobs/MJCalibration.[h,cc] + - put the three camera drawing facilities into a separate base + class MGCamDisplays + + * mjobs/MGCamDisplays.[h,cc] + * mjobs/Makefile + * mjobs/JobsLinkDef.h + - camera drawing facilities: CamDraw, DrawProjection and + DrawRadialProfile, used by all jobs * mcalib/MCalibrationChargeCalc.cc Index: /trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h (revision 3852) +++ /trunk/MagicSoft/Mars/mjobs/JobsLinkDef.h (revision 3853) @@ -9,3 +9,5 @@ #pragma link C++ class MJExtractSignal+; +#pragma link C++ class MGCamDisplays+; + #endif Index: /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.cc =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.cc (revision 3853) +++ /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.cc (revision 3853) @@ -0,0 +1,301 @@ +/* ======================================================================== *\ +! +! * +! * This file is part of MARS, the MAGIC Analysis and Reconstruction +! * Software. It is distributed to you in the hope that it can be a useful +! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. +! * It is distributed WITHOUT ANY WARRANTY. +! * +! * Permission to use, copy, modify and distribute this software and its +! * documentation for any purpose is hereby granted without fee, +! * provided that the above copyright notice appear in all copies and +! * that both that copyright notice and this permission notice appear +! * in supporting documentation. It is provided "as is" without express +! * or implied warranty. +! * +! +! +! Author(s): Thomas Bretz, 1/2004 +! Markus Gaug, 02/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ +///////////////////////////////////////////////////////////////////////////// +// +// MGCamDisplays +// +// Graphical interfaces to display the camera with fits and projections +// +///////////////////////////////////////////////////////////////////////////// +#include "MGCamDisplays.h" + +#include +#include + +#include "MHCamera.h" +#include "MGeomCam.h" +#include "TPad.h" +#include "TProfile.h" +#include "TF1.h" + +#include "MStatusDisplay.h" + +ClassImp(MGCamDisplays); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MGCamDisplays::MGCamDisplays() +{ +} + +// -------------------------------------------------------------------------- +// +// Draw the MHCamera into the MStatusDisplay: +// +// 1) Draw it as histogram (MHCamera::DrawCopy("hist") +// 2) Draw it as a camera, with MHCamera::SetPrettyPalette() set. +// 3) If "rad" is not zero, draw its values vs. the radius from the camera center. +// (DrawRadialProfile()) +// 4) Depending on the variable "fit", draw the values projection on the y-axis +// (DrawProjection()): +// 0: don't draw +// 1: Draw fit to Single Gauss (for distributions flat-fielded over the whole camera) +// 2: Draw and fit to Double Gauss (for distributions different for inner and outer pixels) +// 3: Draw and fit to Triple Gauss (for distributions with inner, outer pixels and outliers) +// 4: Draw and fit to Polynomial grade 0: (for the probability distributions) +// >4: Draw and don;t fit. +// +void MGCamDisplays::CamDraw(TCanvas &c, const Int_t x, const Int_t y, const MHCamera &cam1, + const Int_t fit, const Int_t rad) +{ + + c.cd(x); + gPad->SetBorderMode(0); + gPad->SetTicks(); + MHCamera *obj1=(MHCamera*)cam1.DrawCopy("hist"); + obj1->SetDirectory(NULL); + + c.cd(x+y); + gPad->SetBorderMode(0); + obj1->SetPrettyPalette(); + obj1->Draw(); + + if (rad) + { + c.cd(x+2*y); + gPad->SetBorderMode(0); + gPad->SetTicks(); + DrawRadialProfile(obj1); + } + + + if (!fit) + return; + + c.cd(rad ? x+3*y : x+2*y); + gPad->SetBorderMode(0); + gPad->SetTicks(); + DrawProjection(obj1, fit); +} + +// -------------------------------------------------------------------------- +// +// Draw a projection of MHCamera onto the y-axis values. Depending on the +// variable fit, the following fits are performed: +// +// 1: Single Gauss (for distributions flat-fielded over the whole camera) +// 2: Double Gauss (for distributions different for inner and outer pixels) +// 3: Triple Gauss (for distributions with inner, outer pixels and outliers) +// 4: flat (for the probability distributions) +// +// Moreover, sectors 6,1 and 2 of the camera and sectors 3,4 and 5 are +// drawn separately, for inner and outer pixels. +// +void MGCamDisplays::DrawProjection(MHCamera *obj, Int_t fit) const +{ + + TH1D *obj2 = (TH1D*)obj->Projection(obj->GetName()); + obj2->SetDirectory(0); + obj2->Draw(); + obj2->SetBit(kCanDelete); + + if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic")) + { + TArrayI s0(3); + s0[0] = 6; + s0[1] = 1; + s0[2] = 2; + + TArrayI s1(3); + s1[0] = 3; + s1[1] = 4; + s1[2] = 5; + + TArrayI inner(1); + inner[0] = 0; + + TArrayI outer(1); + outer[0] = 1; + + // Just to get the right (maximum) binning + TH1D *half[4]; + half[0] = obj->ProjectionS(s0, inner, "Sector 6-1-2 Inner"); + half[1] = obj->ProjectionS(s1, inner, "Sector 3-4-5 Inner"); + half[2] = obj->ProjectionS(s0, outer, "Sector 6-1-2 Outer"); + half[3] = obj->ProjectionS(s1, outer, "Sector 3-4-5 Outer"); + + for (int i=0; i<4; i++) + { + half[i]->SetLineColor(kRed+i); + half[i]->SetDirectory(0); + half[i]->SetBit(kCanDelete); + half[i]->Draw("same"); + } + } + + const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); + const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); + const Double_t integ = obj2->Integral("width")/2.5; + const Double_t mean = obj2->GetMean(); + const Double_t rms = obj2->GetRMS(); + const Double_t width = max-min; + + const TString dgausformula = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])" + "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; + + const TString tgausformula = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])" + "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])" + "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])"; + TF1 *f=0; + switch (fit) + { + case 1: + f = new TF1("sgaus", "gaus(0)", min, max); + f->SetLineColor(kYellow); + f->SetBit(kCanDelete); + f->SetParNames("Area", "#mu", "#sigma"); + f->SetParameters(integ/rms, mean, rms); + f->SetParLimits(0, 0, integ); + f->SetParLimits(1, min, max); + f->SetParLimits(2, 0, width/1.5); + + obj2->Fit(f, "QLR"); + break; + + case 2: + f = new TF1("dgaus",dgausformula.Data(),min,max); + f->SetLineColor(kYellow); + f->SetBit(kCanDelete); + f->SetParNames("A_{tot}", "#mu1", "#sigma1", "A2", "#mu2", "#sigma2"); + f->SetParameters(integ,(min+mean)/2.,width/4., + integ/width/2.,(max+mean)/2.,width/4.); + // The left-sided Gauss + f->SetParLimits(0,integ-1.5 , integ+1.5); + f->SetParLimits(1,min+(width/10.), mean); + f->SetParLimits(2,0 , width/2.); + // The right-sided Gauss + f->SetParLimits(3,0 , integ); + f->SetParLimits(4,mean, max-(width/10.)); + f->SetParLimits(5,0 , width/2.); + obj2->Fit(f,"QLRM"); + break; + + case 3: + f = new TF1("tgaus",tgausformula.Data(),min,max); + f->SetLineColor(kYellow); + f->SetBit(kCanDelete); + f->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}", + "A_{2}","#mu_{2}","#sigma_{2}", + "A_{3}","#mu_{3}","#sigma_{3}"); + f->SetParameters(integ,(min+mean)/2,width/4., + integ/width/3.,(max+mean)/2.,width/4., + integ/width/3.,mean,width/2.); + // The left-sided Gauss + f->SetParLimits(0,integ-1.5,integ+1.5); + f->SetParLimits(1,min+(width/10.),mean); + f->SetParLimits(2,width/15.,width/2.); + // The right-sided Gauss + f->SetParLimits(3,0.,integ); + f->SetParLimits(4,mean,max-(width/10.)); + f->SetParLimits(5,width/15.,width/2.); + // The Gauss describing the outliers + f->SetParLimits(6,0.,integ); + f->SetParLimits(7,min,max); + f->SetParLimits(8,width/4.,width/1.5); + obj2->Fit(f,"QLRM"); + break; + + case 4: + obj2->Fit("pol0", "Q"); + obj2->GetFunction("pol0")->SetLineColor(kYellow); + break; + + case 9: + break; + + default: + obj2->Fit("gaus", "Q"); + obj2->GetFunction("gaus")->SetLineColor(kYellow); + break; + } +} + +// -------------------------------------------------------------------------- +// +// Draw a projection of MHCamera vs. the radius from the central pixel. +// +// The inner and outer pixels are drawn separately, both fitted by a polynomial +// of grade 1. +// +void MGCamDisplays::DrawRadialProfile(MHCamera *obj) const +{ + + TProfile *obj2 = (TProfile*)obj->RadialProfile(obj->GetName()); + obj2->SetDirectory(0); + obj2->Draw(); + obj2->SetBit(kCanDelete); + + if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic")) + { + + TArrayI s0(6); + s0[0] = 1; + s0[1] = 2; + s0[2] = 3; + s0[3] = 4; + s0[4] = 5; + s0[5] = 6; + + TArrayI inner(1); + inner[0] = 0; + + TArrayI outer(1); + outer[0] = 1; + + // Just to get the right (maximum) binning + TProfile *half[2]; + half[0] = obj->RadialProfileS(s0, inner,Form("%s%s",obj->GetName(),"Inner")); + half[1] = obj->RadialProfileS(s0, outer,Form("%s%s",obj->GetName(),"Outer")); + + for (Int_t i=0; i<2; i++) + { + Double_t min = obj->GetGeomCam().GetMinRadius(i); + Double_t max = obj->GetGeomCam().GetMaxRadius(i); + + half[i]->SetLineColor(kRed+i); + half[i]->SetDirectory(0); + half[i]->SetBit(kCanDelete); + half[i]->Draw("same"); + half[i]->Fit("pol1","Q","",min,max); + half[i]->GetFunction("pol1")->SetLineColor(kRed+i); + half[i]->GetFunction("pol1")->SetLineWidth(1); + } + } +} + Index: /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.h (revision 3853) +++ /trunk/MagicSoft/Mars/mjobs/MGCamDisplays.h (revision 3853) @@ -0,0 +1,27 @@ +#ifndef MARS_MGCamDisplays +#define MARS_MGCamDisplays + +#ifndef MARS_MAGIC +#include "MAGIC.h" +#endif + +class TCanvas; +class MHCamera; +class MCamEvent; +class MGCamDisplays +{ +protected: + + void DrawProjection ( MHCamera *obj, Int_t fit) const; // Draw projection of pixels values + void DrawRadialProfile( MHCamera *obj) const; // Draw projection of pixels values onto camera radius + void CamDraw(TCanvas &c, const Int_t x, const Int_t y, const MHCamera &cam1, + const Int_t fit, const Int_t rad=0); // Draw histograms into the MStatusDisplay + +public: + + MGCamDisplays(); + + ClassDef(MGCamDisplays, 0) // Tools to display the camera in an MStatusDisplay +}; + +#endif Index: /trunk/MagicSoft/Mars/mjobs/MJCalibration.cc =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibration.cc (revision 3852) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibration.cc (revision 3853) @@ -158,54 +158,4 @@ } -// -------------------------------------------------------------------------- -// -// Draw the MHCamera into the MStatusDisplay: -// -// 1) Draw it as histogram (MHCamera::DrawCopy("hist") -// 2) Draw it as a camera, with MHCamera::SetPrettyPalette() set. -// 3) If "rad" is not zero, draw its values vs. the radius from the camera center. -// (DrawRadialProfile()) -// 4) Depending on the variable "fit", draw the values projection on the y-axis -// (DrawProjection()): -// 0: don't draw -// 1: Draw fit to Single Gauss (for distributions flat-fielded over the whole camera) -// 2: Draw and fit to Double Gauss (for distributions different for inner and outer pixels) -// 3: Draw and fit to Triple Gauss (for distributions with inner, outer pixels and outliers) -// 4: Draw and fit to Polynomial grade 0: (for the probability distributions) -// >4: Draw and don;t fit. -// -void MJCalibration::CamDraw(TCanvas &c, const Int_t x, const Int_t y, const MHCamera &cam1, - const Int_t fit, const Int_t rad) -{ - - c.cd(x); - gPad->SetBorderMode(0); - gPad->SetTicks(); - MHCamera *obj1=(MHCamera*)cam1.DrawCopy("hist"); - obj1->SetDirectory(NULL); - - c.cd(x+y); - gPad->SetBorderMode(0); - obj1->SetPrettyPalette(); - obj1->AddNotify(&fCalibrationCam); - obj1->Draw(); - - if (rad) - { - c.cd(x+2*y); - gPad->SetBorderMode(0); - gPad->SetTicks(); - DrawRadialProfile(obj1); - } - - - if (!fit) - return; - - c.cd(rad ? x+3*y : x+2*y); - gPad->SetBorderMode(0); - gPad->SetTicks(); - DrawProjection(obj1, fit); -} // -------------------------------------------------------------------------- @@ -677,198 +627,4 @@ } -// -------------------------------------------------------------------------- -// -// Draw a projection of MHCamera onto the y-axis values. Depending on the -// variable fit, the following fits are performed: -// -// 1: Single Gauss (for distributions flat-fielded over the whole camera) -// 2: Double Gauss (for distributions different for inner and outer pixels) -// 3: Triple Gauss (for distributions with inner, outer pixels and outliers) -// 4: flat (for the probability distributions) -// -// Moreover, sectors 6,1 and 2 of the camera and sectors 3,4 and 5 are -// drawn separately, for inner and outer pixels. -// -void MJCalibration::DrawProjection(MHCamera *obj, Int_t fit) const -{ - - TH1D *obj2 = (TH1D*)obj->Projection(obj->GetName()); - obj2->SetDirectory(0); - obj2->Draw(); - obj2->SetBit(kCanDelete); - - if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic")) - { - TArrayI s0(3); - s0[0] = 6; - s0[1] = 1; - s0[2] = 2; - - TArrayI s1(3); - s1[0] = 3; - s1[1] = 4; - s1[2] = 5; - - TArrayI inner(1); - inner[0] = 0; - - TArrayI outer(1); - outer[0] = 1; - - // Just to get the right (maximum) binning - TH1D *half[4]; - half[0] = obj->ProjectionS(s0, inner, "Sector 6-1-2 Inner"); - half[1] = obj->ProjectionS(s1, inner, "Sector 3-4-5 Inner"); - half[2] = obj->ProjectionS(s0, outer, "Sector 6-1-2 Outer"); - half[3] = obj->ProjectionS(s1, outer, "Sector 3-4-5 Outer"); - - for (int i=0; i<4; i++) - { - half[i]->SetLineColor(kRed+i); - half[i]->SetDirectory(0); - half[i]->SetBit(kCanDelete); - half[i]->Draw("same"); - } - } - - const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); - const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); - const Double_t integ = obj2->Integral("width")/2.5; - const Double_t mean = obj2->GetMean(); - const Double_t rms = obj2->GetRMS(); - const Double_t width = max-min; - - const TString dgausformula = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])" - "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; - - const TString tgausformula = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])" - "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])" - "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])"; - TF1 *f=0; - switch (fit) - { - case 1: - f = new TF1("sgaus", "gaus(0)", min, max); - f->SetLineColor(kYellow); - f->SetBit(kCanDelete); - f->SetParNames("Area", "#mu", "#sigma"); - f->SetParameters(integ/rms, mean, rms); - f->SetParLimits(0, 0, integ); - f->SetParLimits(1, min, max); - f->SetParLimits(2, 0, width/1.5); - - obj2->Fit(f, "QLR"); - break; - - case 2: - f = new TF1("dgaus",dgausformula.Data(),min,max); - f->SetLineColor(kYellow); - f->SetBit(kCanDelete); - f->SetParNames("A_{tot}", "#mu1", "#sigma1", "A2", "#mu2", "#sigma2"); - f->SetParameters(integ,(min+mean)/2.,width/4., - integ/width/2.,(max+mean)/2.,width/4.); - // The left-sided Gauss - f->SetParLimits(0,integ-1.5 , integ+1.5); - f->SetParLimits(1,min+(width/10.), mean); - f->SetParLimits(2,0 , width/2.); - // The right-sided Gauss - f->SetParLimits(3,0 , integ); - f->SetParLimits(4,mean, max-(width/10.)); - f->SetParLimits(5,0 , width/2.); - obj2->Fit(f,"QLRM"); - break; - - case 3: - f = new TF1("tgaus",tgausformula.Data(),min,max); - f->SetLineColor(kYellow); - f->SetBit(kCanDelete); - f->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}", - "A_{2}","#mu_{2}","#sigma_{2}", - "A_{3}","#mu_{3}","#sigma_{3}"); - f->SetParameters(integ,(min+mean)/2,width/4., - integ/width/3.,(max+mean)/2.,width/4., - integ/width/3.,mean,width/2.); - // The left-sided Gauss - f->SetParLimits(0,integ-1.5,integ+1.5); - f->SetParLimits(1,min+(width/10.),mean); - f->SetParLimits(2,width/15.,width/2.); - // The right-sided Gauss - f->SetParLimits(3,0.,integ); - f->SetParLimits(4,mean,max-(width/10.)); - f->SetParLimits(5,width/15.,width/2.); - // The Gauss describing the outliers - f->SetParLimits(6,0.,integ); - f->SetParLimits(7,min,max); - f->SetParLimits(8,width/4.,width/1.5); - obj2->Fit(f,"QLRM"); - break; - - case 4: - obj2->Fit("pol0", "Q"); - obj2->GetFunction("pol0")->SetLineColor(kYellow); - break; - - case 9: - break; - - default: - obj2->Fit("gaus", "Q"); - obj2->GetFunction("gaus")->SetLineColor(kYellow); - break; - } -} - -// -------------------------------------------------------------------------- -// -// Draw a projection of MHCamera vs. the radius from the central pixel. -// -// The inner and outer pixels are drawn separately, both fitted by a polynomial -// of grade 1. -// -void MJCalibration::DrawRadialProfile(MHCamera *obj) const -{ - - TProfile *obj2 = (TProfile*)obj->RadialProfile(obj->GetName()); - obj2->SetDirectory(0); - obj2->Draw(); - obj2->SetBit(kCanDelete); - - if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic")) - { - - TArrayI s0(6); - s0[0] = 1; - s0[1] = 2; - s0[2] = 3; - s0[3] = 4; - s0[4] = 5; - s0[5] = 6; - - TArrayI inner(1); - inner[0] = 0; - - TArrayI outer(1); - outer[0] = 1; - - // Just to get the right (maximum) binning - TProfile *half[2]; - half[0] = obj->RadialProfileS(s0, inner,Form("%s%s",obj->GetName(),"Inner")); - half[1] = obj->RadialProfileS(s0, outer,Form("%s%s",obj->GetName(),"Outer")); - - for (Int_t i=0; i<2; i++) - { - Double_t min = obj->GetGeomCam().GetMinRadius(i); - Double_t max = obj->GetGeomCam().GetMaxRadius(i); - - half[i]->SetLineColor(kRed+i); - half[i]->SetDirectory(0); - half[i]->SetBit(kCanDelete); - half[i]->Draw("same"); - half[i]->Fit("pol1","Q","",min,max); - half[i]->GetFunction("pol1")->SetLineColor(kRed+i); - half[i]->GetFunction("pol1")->SetLineWidth(1); - } - } -} Index: /trunk/MagicSoft/Mars/mjobs/MJCalibration.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibration.h (revision 3852) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibration.h (revision 3853) @@ -14,12 +14,12 @@ #include "MBadPixelsCam.h" #endif +#ifndef MARS_MGCamDisplays +#include "MGCamDisplays.h" +#endif -class TCanvas; -class MHCamera; -class MCamEvent; class MRunIter; class MParList; class MPedestalCam; -class MJCalibration : public MParContainer +class MJCalibration : public MParContainer, public MGCamDisplays { private: @@ -46,9 +46,4 @@ UInt_t fExtractorLevel; // Level signal extractor (e.g. MExtractSignal2) UInt_t fArrivalTimeLevel; // Level arr. time extractor (e.g. MArrivalTimeCalc2) - - void DrawProjection ( MHCamera *obj, Int_t fit) const; // Draw projection of pixels values - void DrawRadialProfile( MHCamera *obj) const; // Draw projection of pixels values onto camera radius - void CamDraw(TCanvas &c, const Int_t x, const Int_t y, const MHCamera &cam1, - const Int_t fit, const Int_t rad=0); // Draw histograms into the MStatusDisplay void DisplayResult(MParList &plist); Index: /trunk/MagicSoft/Mars/mjobs/Makefile =================================================================== --- /trunk/MagicSoft/Mars/mjobs/Makefile (revision 3852) +++ /trunk/MagicSoft/Mars/mjobs/Makefile (revision 3853) @@ -34,5 +34,6 @@ SRCFILES = MJPedestal.cc \ MJCalibration.cc \ - MJExtractSignal.cc + MJExtractSignal.cc \ + MGCamDisplays.cc SRCS = $(SRCFILES)