Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 1556) +++ trunk/MagicSoft/Mars/Changelog (revision 1557) @@ -1,3 +1,21 @@ -*-*- END -*-*- + 2002/10/29: Thomas Bretz + + * manalysis/MHadroness.[h,cc]: + - renamed to MHadronness + + * mhist/MHHadroness.[h,cc]: + - renamed to MHHadronness + - small changes to the graphical and text output + + * manalysis/Makefile, manalysis/AnalysisLinkDef.h, macros/multidimdist.C, + manalysis/MMultiDimDistCalc.[h,cc], macros/comprob.C, + mhist/Makefile, mhist/HistLinkDef.h: + - changed Hadroness to Hadronness + + * manalysis/MCerPhotAnal.cc: + - small changes to the code layout + + 2002/10/24: Thomas Bretz Index: trunk/MagicSoft/Mars/macros/comprob.C =================================================================== --- trunk/MagicSoft/Mars/macros/comprob.C (revision 1556) +++ trunk/MagicSoft/Mars/macros/comprob.C (revision 1557) @@ -113,5 +113,5 @@ // fill probabilities (hadroness) into histogram - MFillH fill2("MHHadroness"); + MFillH fill2("MHHadronness"); tlist2.AddToList(&fill2); @@ -124,6 +124,6 @@ // the console window // - plist.FindObject("MHHadroness")->DrawClone(); - plist.FindObject("MHHadroness")->Print(); + plist.FindObject("MHHadronness")->DrawClone(); + plist.FindObject("MHHadronness")->Print(); } Index: trunk/MagicSoft/Mars/macros/multidimdist.C =================================================================== --- trunk/MagicSoft/Mars/macros/multidimdist.C (revision 1556) +++ trunk/MagicSoft/Mars/macros/multidimdist.C (revision 1557) @@ -48,5 +48,5 @@ // CalEvents: Calibration Events // - MReadMarsFile read("Events", "~/Mars/star2.root"); + MReadMarsFile read("Events", "~/Mars/star.root"); read.DisableAutoScheme(); tlist.AddToList(&read); @@ -107,5 +107,6 @@ plist.Replace(&tlist2); - MReadMarsFile read2("Events", "~/Mars/star2.root"); + MReadMarsFile read2("Events", "~/Mars/star.root"); +// read2.AddFile("~/Mars/star2.root"); read2.DisableAutoScheme(); tlist2.AddToList(&read2); @@ -116,5 +117,5 @@ tlist2.AddToList(&calc); - MFillH fillh("MHHadroness"); + MFillH fillh("MHHadronness"); /* @@ -134,5 +135,5 @@ tlist2.PrintStatistics(); - plist.FindObject("MHHadroness")->DrawClone(); - plist.FindObject("MHHadroness")->Print(); + plist.FindObject("MHHadronness")->DrawClone(); + plist.FindObject("MHHadronness")->Print(); } Index: trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1556) +++ trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1557) @@ -21,5 +21,5 @@ #pragma link C++ class MMcPedestalNSBAdd+; -#pragma link C++ class MHadroness+; +#pragma link C++ class MHadronness+; #pragma link C++ class MCompProbCalc+; Index: trunk/MagicSoft/Mars/manalysis/MCerPhotAnal.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MCerPhotAnal.cc (revision 1556) +++ trunk/MagicSoft/Mars/manalysis/MCerPhotAnal.cc (revision 1557) @@ -123,5 +123,5 @@ while (pixel.Next()) - { + { Byte_t *ptr = pixel.GetHiGainSamples(); const Byte_t *end = ptr + fRawEvt->GetNumHiGainSamples(); @@ -134,43 +134,49 @@ Int_t slicelocal=0; - Float_t nphot; - Float_t pedes, sigmaped; - - do - { - sumlocal = 0; - for(Int_t i = 0; i<5;i++) - sumlocal+=*(ptr+i); - if (sumpeak < sumlocal){ - slice=slicelocal; - sumpeak = sumlocal; - } - slicelocal++; + do + { + sumlocal = 0; + for (Int_t i = 0; i<5;i++) + sumlocal += *(ptr+i); + + if (sumpeak < sumlocal) + { + slice=slicelocal; + sumpeak = sumlocal; + } + + slicelocal++; sum += *ptr; - } - while (++ptr != limit); - - do sum += *ptr; - while (++ptr != end); - - pedes = (Float_t)(sum-sumpeak)/(fRawEvt->GetNumHiGainSamples()-5); - nphot=(Float_t)sumpeak - 5.0*pedes; - sigmaped=0; - sumlocal=0; - slicelocal=0; + } while (++ptr != limit); + + do sum += *ptr; + while (++ptr != end); + + Float_t pedes = (Float_t)(sum-sumpeak)/(fRawEvt->GetNumHiGainSamples()-5); + Float_t nphot = (Float_t)sumpeak - 5.0*pedes; + + Float_t sigmaped=0; + + slicelocal=0; + sumlocal = 0; ptr = pixel.GetHiGainSamples(); - do { - if (slicelocal==slice) - ptr+=4; - else{ - sumlocal=*ptr; - sigmaped+=((Float_t)sumlocal-pedes)*((Float_t)sumlocal-pedes); - } - slicelocal++; + do + { + if (slicelocal==slice) + ptr += 4; + else + { + sumlocal = *ptr; + + const Float_t d = (Float_t)sumlocal-pedes; + sigmaped += d*d; + } + slicelocal++; } - while(++ptr != end); - sigmaped/=(fRawEvt->GetNumHiGainSamples()-5); - sigmaped=sqrt(sumlocal); + while (++ptr != end); + + sigmaped /= (fRawEvt->GetNumHiGainSamples()-5); + sigmaped = sqrt(sumlocal); const UInt_t pixid = pixel.GetPixelId(); @@ -188,5 +194,5 @@ fCerPhotEvt->AddPixel(pixid, nphot, sigmaped/2.236); - ped.SetPedestal(pedes,sigmaped); + ped.SetPedestal(pedes, sigmaped); ped.SetPedestalRms(sigmaped/sqrt(fRawEvt->GetNumHiGainSamples()-5), sigmaped/sqrt(2*(fRawEvt->GetNumHiGainSamples()-5))); Index: trunk/MagicSoft/Mars/manalysis/MHadroness.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHadroness.cc (revision 1556) +++ (revision ) @@ -1,69 +1,0 @@ -/* ======================================================================== *\ -! -! * -! * 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 5/2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHadroness -// -// The Value must be in the range [0,1] -// It should be the value used for gamma/hadron seperation. -// For quality histograms look at MHHadroness -// -///////////////////////////////////////////////////////////////////////////// -#include "MHadroness.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHadroness); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MHadroness::MHadroness(const char *name, const char *title) - : fHadroness(-1) -{ - fName = name ? name : "MHadroness"; - fTitle = title ? title : "Storage container for hadroness value"; -} - -// -------------------------------------------------------------------------- -// -// Reset hadroness, -1 indicates: invalid value -// -void MHadroness::Reset() -{ - fHadroness = -1; -} - -// -------------------------------------------------------------------------- -// -// Print the hillas Parameters to *fLog -// -void MHadroness::Print(Option_t *) const -{ - *fLog << all << GetDescriptor() << ": Hadroness = " << fHadroness << endl; -} - Index: trunk/MagicSoft/Mars/manalysis/MHadroness.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHadroness.h (revision 1556) +++ (revision ) @@ -1,26 +1,0 @@ -#ifndef MARS_MHadroness -#define MARS_MHadroness - -#ifndef MARS_MParContainer -#include "MParContainer.h" -#endif - -class MHadroness : public MParContainer -{ -private: - Double_t fHadroness; // [0,1] Hadroness - -public: - MHadroness(const char *name=NULL, const char *title=NULL); - - void Reset(); - - virtual void Print(Option_t *opt=NULL) const; - - Double_t GetHadroness() const { return fHadroness; } - void SetHadroness(Double_t h) { fHadroness = h; } - - ClassDef(MHadroness, 1) // Storage Container for the hadroness -}; - -#endif Index: trunk/MagicSoft/Mars/manalysis/MHadronness.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHadronness.cc (revision 1557) +++ trunk/MagicSoft/Mars/manalysis/MHadronness.cc (revision 1557) @@ -0,0 +1,69 @@ +/* ======================================================================== *\ +! +! * +! * 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 5/2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHadronness +// +// The Value must be in the range [0,1] +// It should be the value used for gamma/hadron seperation. +// For quality histograms look at MHHadronness +// +///////////////////////////////////////////////////////////////////////////// +#include "MHadronness.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MHadronness); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MHadronness::MHadronness(const char *name, const char *title) + : fHadronness(-1) +{ + fName = name ? name : "MHadronness"; + fTitle = title ? title : "Storage container for hadroness value"; +} + +// -------------------------------------------------------------------------- +// +// Reset hadroness, -1 indicates: invalid value +// +void MHadronness::Reset() +{ + fHadronness = -1; +} + +// -------------------------------------------------------------------------- +// +// Print the hillas Parameters to *fLog +// +void MHadronness::Print(Option_t *) const +{ + *fLog << all << GetDescriptor() << ": Hadronness = " << fHadronness << endl; +} + Index: trunk/MagicSoft/Mars/manalysis/MHadronness.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHadronness.h (revision 1557) +++ trunk/MagicSoft/Mars/manalysis/MHadronness.h (revision 1557) @@ -0,0 +1,26 @@ +#ifndef MARS_MHadronness +#define MARS_MHadronness + +#ifndef MARS_MParContainer +#include "MParContainer.h" +#endif + +class MHadronness : public MParContainer +{ +private: + Double_t fHadronness; // [0,1] Hadronness + +public: + MHadronness(const char *name=NULL, const char *title=NULL); + + void Reset(); + + virtual void Print(Option_t *opt=NULL) const; + + Double_t GetHadronness() const { return fHadronness; } + void SetHadronness(Double_t h) { fHadronness = h; } + + ClassDef(MHadronness, 1) // Storage Container for the hadroness +}; + +#endif Index: trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.cc (revision 1556) +++ trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.cc (revision 1557) @@ -54,5 +54,5 @@ #include "MDataArray.h" -#include "MHadroness.h" +#include "MHadronness.h" ClassImp(MMultiDimDistCalc); @@ -149,6 +149,6 @@ } - fHadroness = (MHadroness*)plist->FindCreateObj("MHadroness"); - if (!fHadroness) + fHadronness = (MHadronness*)plist->FindCreateObj("MHadronness"); + if (!fHadronness) return kFALSE; @@ -200,6 +200,6 @@ Double_t dh = fMHadrons->CalcDist(event, numh); - //fHadroness->SetHadroness(dg/(dg+dh)); - fHadroness->SetHadroness(exp(-dh/dg)); + //fHadronness->SetHadronness(dg/(dg+dh)); + fHadronness->SetHadronness(exp(-dh/dg)); return kTRUE; Index: trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.h (revision 1556) +++ trunk/MagicSoft/Mars/manalysis/MMultiDimDistCalc.h (revision 1557) @@ -8,19 +8,19 @@ class MHMatrix; class MParList; -class MHadroness; class MDataArray; +class MHadronness; class MMultiDimDistCalc : public MTask { private: - Int_t fNum; // number of distances used for an avarage - Bool_t fUseKernel; // Flag whether kernel method should be used + Int_t fNum; // number of distances used for an avarage + Bool_t fUseKernel; // Flag whether kernel method should be used - MHMatrix *fMGammas; //! Gammas describing matrix - MHMatrix *fMHadrons; //! Hadrons (non gammas) describing matrix + MHMatrix *fMGammas; //! Gammas describing matrix + MHMatrix *fMHadrons; //! Hadrons (non gammas) describing matrix - MHadroness *fHadroness; //! Output container for calculated hadroness + MHadronness *fHadronness; //! Output container for calculated hadroness - MDataArray *fData; //! Used to store the MDataChains to get the event values + MDataArray *fData; //! Used to store the MDataChains to get the event values void StreamPrimitive(ofstream &out) const; Index: trunk/MagicSoft/Mars/manalysis/Makefile =================================================================== --- trunk/MagicSoft/Mars/manalysis/Makefile (revision 1556) +++ trunk/MagicSoft/Mars/manalysis/Makefile (revision 1557) @@ -39,5 +39,5 @@ MSrcPosCam.cc \ MCameraSmooth.cc \ - MHadroness.cc \ + MHadronness.cc \ MMatrixLoop.cc \ MCompProbCalc.cc \ Index: trunk/MagicSoft/Mars/mhist/HistLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mhist/HistLinkDef.h (revision 1556) +++ trunk/MagicSoft/Mars/mhist/HistLinkDef.h (revision 1557) @@ -34,5 +34,5 @@ #pragma link C++ class MHCompProb+; -#pragma link C++ class MHHadroness+; +#pragma link C++ class MHHadronness+; #pragma link C++ class MHMcRate+; Index: trunk/MagicSoft/Mars/mhist/MHHadroness.cc =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHadroness.cc (revision 1556) +++ (revision ) @@ -1,426 +1,0 @@ -/* ======================================================================== *\ -! -! * -! * 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, 5/2002 -! Author(s): Abelardo Moralejo -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHHadroness -// -// This is histogram is a way to evaluate the quality of a gamma/hadron -// seperation method. It is filled from a MHadroness container, which -// stores a hadroness for the current event. The Value must be in the -// range [0,1]. To fill the histograms correctly the information -// whether it is a gamma or hadron (not a gamma) must be available from -// a MMcEvt container. -// -// In the constructor you can change the number of used bns for the -// evaluation. -// -// The meaning of the histograms (Draw, DrawClone) are the following: -// * Upper Left Corner: -// - black: histogram of all hadronesses for gammas -// - red: histogram of all hadronesses for non gammas -// * Upper Right Corner: -// - black: acceptance of gammas vs. the hadroness -// - red: acceptance of non gammas vs. the hadroness -// - blue: 2D distance of (acceptance_hadrons, acceptances_gammas) -// to optimum (0, 1) -// 1-sqrt(Ag*Ag + (1-Ah)*(1-Ah)) -// * Bottom Left Corner: -// Naive quality factor: Ag/sqrt(Ah) -// * Bottom Right Corner: -// - black: Acceprtance Gammas vs. Acceptance Hadrons -// - blue cross: minimum of distance to (0, 1) -// -//////////////////////////////////////////////////////////////////////////// -#include "MHHadroness.h" - -#include -#include -#include -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" -#include "MBinning.h" -#include "MHadroness.h" - -#include "MMcEvt.hxx" - -ClassImp(MHHadroness); - -// -------------------------------------------------------------------------- -// -// Setup histograms, nbins is the number of bins used for the evaluation. -// The default is 100 bins. -// -MHHadroness::MHHadroness(Int_t nbins, const char *name, const char *title) -{ - // - // set the name and title of this object - // - fName = name ? name : "MHHadroness"; - fTitle = title ? title : "Gamma/Hadron Separation Quality Histograms"; - - fGraph = new TGraph; - fGraph->SetTitle("Acceptance Gammas vs. Hadrons"); - fGraph->SetMaximum(1); - - fGhness = new TH1D("Ghness", "Hadronness", nbins, 0, 1); - fPhness = new TH1D("Phness", "Hadronness", nbins, 0, 1); - fGhness->SetXTitle("Hadroness"); - fPhness->SetXTitle("Hadroness"); - fGhness->SetYTitle("Counts"); - fPhness->SetYTitle("Counts"); - - fIntGhness = new TH1D("AccGammas", "Acceptance", nbins, 0, 1); - fIntPhness = new TH1D("AccHadrons", "Acceptance", nbins, 0, 1); - fIntGhness->SetXTitle("Hadroness"); - fIntPhness->SetXTitle("Hadroness"); - fIntGhness->SetYTitle("Acceptance"); - fIntPhness->SetYTitle("Acceptance"); - - /* - fQfac = new TH1D("Qfac", "Naive Quality factor", nbins, 0, 1); - fQfac->SetXTitle("Hadroness"); - fQfac->SetYTitle("Quality"); - */ - fQfac = new TGraph; - fQfac->SetTitle(" Naive Quality factor "); - - fMinDist = new TH1D("MinDist", "Minimum Dist to (0, 1)", nbins, 0, 1); - fMinDist->SetXTitle("Hadroness"); - fMinDist->SetYTitle("Distance"); - - //fQfac->SetDirectory(NULL); - fGhness->SetDirectory(NULL); - fPhness->SetDirectory(NULL); - fMinDist->SetDirectory(NULL); - fIntGhness->SetDirectory(NULL); - fIntPhness->SetDirectory(NULL); -} - -// -------------------------------------------------------------------------- -// -// Delete the histograms. -// -MHHadroness::~MHHadroness() -{ - delete fGhness; - delete fIntGhness; - delete fPhness; - delete fIntPhness; - delete fQfac; - delete fMinDist; - delete fGraph; -} - -// -------------------------------------------------------------------------- -// -// Setup Filling of the histograms. It needs: -// MMcEvt and MHadroness -// -Bool_t MHHadroness::SetupFill(const MParList *plist) -{ - fMcEvt = (MMcEvt*)plist->FindObject("MMcEvt"); - if (!fMcEvt) - { - *fLog << err << dbginf << "MMcEvt not found... aborting." << endl; - return kFALSE; - } - - fHadroness = (MHadroness*)plist->FindObject("MHadroness"); - if (!fHadroness) - { - *fLog << err << dbginf << "MHadroness not found... aborting." << endl; - return kFALSE; - } - - /* - MBinning* bins = (MBinning*)plist->FindObject("BinningHadroness"); - if (!bins) - { - *fLog << err << dbginf << "BinningHadroness [MBinning] not found... aborting." << endl; - return kFALSE; - } - - SetBinning(&fHist, binsalpha, binsenergy, binstheta); - - fHist.Sumw2(); - */ - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Fill the Hadroness from a MHadroness container into the corresponding -// histogram dependant on the particle id. -// -Bool_t MHHadroness::Fill(const MParContainer *par) -{ - if (fMcEvt->GetPartId()==kGAMMA) - fGhness->Fill(fHadroness->GetHadroness()); - else - fPhness->Fill(fHadroness->GetHadroness()); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Finalize the histograms: -// - integrate the hadroness histograms --> acceptance -// - fill the Minimum Distance histogram (formular see class description) -// - fill the Quality histogram (formular see class description) -// -Bool_t MHHadroness::Finalize() -{ - Int_t n = fGhness->GetNbinsX(); - - fGraph->Set(n); - fQfac->Set(n); - - const Stat_t sumg = fGhness->Integral(1, n+1); - const Stat_t sump = fPhness->Integral(1, n+1); - - for (Int_t i=1; i<=n; i++) - { - const Stat_t ip = fPhness->Integral(1, i)/sump; - const Stat_t ig = fGhness->Integral(1, i)/sumg; - - fIntPhness->SetBinContent(i, ip); - fIntGhness->SetBinContent(i, ig); - - fGraph->SetPoint(i, ip, ig); - - if (ip<=0) - continue; - - fMinDist->SetBinContent(i, 1-sqrt(ip*ip + (ig-1)*(ig-1))); - fQfac->SetPoint(i, ig, ig/sqrt(ip)); - //fQfac->SetBinContent(fQfac->FindBin(ig), ig/sqrt(ip)); - *fLog << all << i << ": " << ig << " " << ig/sqrt(ip) << endl; - } - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Search the corresponding points for the given hadron acceptance (acchad) -// and interpolate the to pointsd (linear) -// -Double_t MHHadroness::GetGammaAcceptance(Double_t acchad) const -{ - const Int_t n = fGraph->GetN(); - const Double_t *x = fGraph->GetX(); - const Double_t *y = fGraph->GetY(); - - Int_t i = 0; - while (iGetEntries() << " Gammas and " << fPhness->GetEntries() << " Hadrons." << endl; - *fLog << "Acc Gammas at 1% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.01)*100) << "%" << endl; - *fLog << "Acc Gammas at 2% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.02)*100) << "%" << endl; - *fLog << "Acc Gammas at 5% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.05)*100) << "%" << endl; - *fLog << "Acc Gammas at 10% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.1)*100) << "%" << endl; - *fLog << "Acc Gammas at 20% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.2)*100) << "%" << endl; - *fLog << "Acc Gammas at 30% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.3)*100) << "%" << endl; - *fLog << "Acc Gammas at 40% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.4)*100) << "%" << endl; - *fLog << "Acc Gammas at 50% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.5)*100) << "%" << endl; - const Int_t h = fMinDist->GetMaximumBin(); - *fLog << "Minimum Distance to (0, 1): " << Form("%.2f", 1-fMinDist->GetMaximum()) << " @ H=" << fMinDist->GetBinCenter(h) << endl; - *fLog << " Acc Gammas = " << Form("%3.0f", fIntGhness->GetBinContent(h)*100) << "%, "; - *fLog << "Acc Hadrons = " << Form("%3.0f", fIntPhness->GetBinContent(h)*100) << "%" << endl; - /* - const Int_t q = GetQfac()->GetMaximumBin(); - *fLog << "Maximum Q-Factor: " << GetQfac()->GetMaximum() << " @ H="; - *fLog << GetQfac()->GetBinCenter(q) << endl;; - *fLog << " Acc Gammas = " << Form("%3.0f", fIntGhness->GetBinContent(q)*100) << "%, "; - *fLog << "Acc Hadrons = " << Form("%3.0f", fIntPhness->GetBinContent(q)*100) << "%" << endl; - */ - *fLog << endl; -} - -// -------------------------------------------------------------------------- -// -// Draw clone of all histograms. (For the Meaning see class description) -// -TObject *MHHadroness::DrawClone(Option_t *opt) const -{ - TCanvas &c = *MakeDefCanvas("Hadroness", fTitle); - c.Divide(2, 2); - - gROOT->SetSelectedPad(NULL); - - c.cd(1); - gStyle->SetOptStat(10); - Getghness()->DrawCopy(); - Getphness()->SetLineColor(kRed); - Getphness()->DrawCopy("same"); - c.cd(2); - gStyle->SetOptStat(0); - Getighness()->DrawCopy(); - Getiphness()->SetLineColor(kRed); - Getiphness()->DrawCopy("same"); - fMinDist->SetLineColor(kBlue); - fMinDist->DrawCopy("same"); - c.cd(3); - //GetQfac()->DrawCopy(); - TGraph &g2 = (TGraph&)*fQfac->DrawClone("AC"); - g2.SetBit(kCanDelete); - gPad->Modified(); - gPad->Update(); - if (g2.GetHistogram()) - { - g2.GetXaxis()->SetRangeUser(0, 1); - g2.GetXaxis()->SetTitle("Acceptance Gammas"); - g2.GetYaxis()->SetTitle("Quality"); - g2.SetMarkerStyle(kFullDotSmall); - g2.Draw("P"); - - gPad->Modified(); - gPad->Update(); - } - c.cd(4); - gPad->Modified(); - gPad->Update(); - TGraph &g = (TGraph&)*fGraph->DrawClone("AC"); - g.SetBit(kCanDelete); - gPad->Modified(); - gPad->Update(); - if (g.GetHistogram()) - { - g.GetXaxis()->SetRangeUser(0, 1); - g.GetXaxis()->SetTitle("Acceptance Hadrons"); - g.GetYaxis()->SetTitle("Acceptance Gammas"); - g.SetMarkerStyle(kFullDotSmall); - g.Draw("P"); - - gPad->Modified(); - gPad->Update(); - } - const Int_t h = fMinDist->GetMaximumBin(); - TMarker *m = new TMarker(fIntPhness->GetBinContent(h), - fIntGhness->GetBinContent(h), kStar); - m->SetMarkerColor(kBlue); - m->SetBit(kCanDelete); - m->Draw(); - return &c; -} - -// -------------------------------------------------------------------------- -// -// Draw all histograms. (For the Meaning see class description) -// -void MHHadroness::Draw(Option_t *) -{ - if (!gPad) - MakeDefCanvas("Hadroness", fTitle); - - gPad->Divide(2, 2); - - gPad->cd(1); - gStyle->SetOptStat(10); - Getghness()->Draw(); - Getphness()->SetLineColor(kRed); - Getphness()->Draw("same"); - gPad->cd(2); - gStyle->SetOptStat(0); - Getighness()->Draw(); - Getiphness()->SetLineColor(kRed); - Getiphness()->Draw("same"); - fMinDist->SetLineColor(kBlue); - fMinDist->Draw("same"); - gPad->cd(3); - //GetQfac()->Draw(); - fQfac->Draw("AC"); - gPad->Modified(); - gPad->Update(); - if (fQfac->GetHistogram()) - { - fQfac->GetXaxis()->SetRangeUser(0, 1); - fQfac->GetXaxis()->SetTitle("Acceptance Gammas"); - fQfac->GetYaxis()->SetTitle("Quality"); - fQfac->SetMarkerStyle(kFullDotSmall); - fQfac->Draw("P"); - - gPad->Modified(); - gPad->Update(); - } - gPad->cd(4); - gPad->Modified(); - gPad->Update(); - fGraph->Draw("AC"); - gPad->Modified(); - gPad->Update(); - if (fGraph->GetHistogram()) - { - fGraph->GetXaxis()->SetRangeUser(0, 1); - fGraph->GetXaxis()->SetTitle("Acceptance Hadrons"); - fGraph->GetYaxis()->SetTitle("Acceptance Gammas"); - fGraph->SetMarkerStyle(kFullDotSmall); - fGraph->Draw("P"); - gPad->Modified(); - gPad->Update(); - } - const Int_t h = fMinDist->GetMaximumBin(); - TMarker *m = new TMarker(fIntPhness->GetBinContent(h), - fIntGhness->GetBinContent(h), kStar); - m->SetMarkerColor(kBlue); - m->SetBit(kCanDelete); - m->Draw(); -} Index: trunk/MagicSoft/Mars/mhist/MHHadroness.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHadroness.h (revision 1556) +++ (revision ) @@ -1,53 +1,0 @@ -#ifndef MARS_MHHadroness -#define MARS_MHHadroness - -#ifndef MARS_MH -#include "MH.h" -#endif - -class TH1D; -class TGraph; -class MParList; -class MMcEvt; -class MHadroness; - -class MHHadroness : public MH -{ -private: - const MMcEvt *fMcEvt; //! - const MHadroness *fHadroness; //! - - TH1D* fPhness; // Hadrons Hadroness - TH1D* fGhness; // Gammas Hadroness - TH1D* fIntPhness; // Hadrons Acceptance - TH1D* fIntGhness; // Gammas Acceptance - TH1D* fMinDist; // Minimum Distance to optimum acceptance - - TGraph *fQfac; // Quality factor - TGraph *fGraph; // gamma acceptance vs. hadron acceptance - -public: - MHHadroness(Int_t nbins=100, const char *name=NULL, const char *title=NULL); - ~MHHadroness(); - - Double_t GetGammaAcceptance(Double_t acchad) const; - - TH1D *Getghness() const { return fGhness; } - TH1D *Getphness() const { return fPhness; } - TH1D *Getighness() const { return fIntGhness; } - TH1D *Getiphness() const { return fIntPhness; } - //TH2D *GetHist() const { return fHist; } - - Bool_t SetupFill(const MParList *plist); - Bool_t Fill(const MParContainer *par); - Bool_t Finalize(); - - void Print(Option_t *option="") const; - - void Draw(Option_t *opt=""); - TObject *DrawClone(Option_t *opt="") const; - - ClassDef(MHHadroness, 1) // Gamma/Hadron Separation Quality Histograms -}; - -#endif Index: trunk/MagicSoft/Mars/mhist/MHHadronness.cc =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHadronness.cc (revision 1557) +++ trunk/MagicSoft/Mars/mhist/MHHadronness.cc (revision 1557) @@ -0,0 +1,492 @@ +/* ======================================================================== *\ +! +! * +! * 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, 5/2002 +! Author(s): Abelardo Moralejo +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHHadronness +// +// This is histogram is a way to evaluate the quality of a gamma/hadron +// seperation method. It is filled from a MHadronness container, which +// stores a hadroness for the current event. The Value must be in the +// range [0,1]. To fill the histograms correctly the information +// whether it is a gamma or hadron (not a gamma) must be available from +// a MMcEvt container. +// +// In the constructor you can change the number of used bns for the +// evaluation. +// +// The meaning of the histograms (Draw, DrawClone) are the following: +// * Upper Left Corner: +// - black: histogram of all hadronesses for gammas +// - red: histogram of all hadronesses for non gammas +// * Upper Right Corner: +// - black: acceptance of gammas vs. the hadroness +// - red: acceptance of non gammas vs. the hadroness +// - blue: 2D distance of (acceptance_hadrons, acceptances_gammas) +// to optimum (0, 1) +// 1-sqrt(Ag*Ag + (1-Ah)*(1-Ah)) +// * Bottom Left Corner: +// Naive quality factor: Ag/sqrt(Ah) +// * Bottom Right Corner: +// - black: Acceprtance Gammas vs. Acceptance Hadrons +// - blue cross: minimum of distance to (0, 1) +// +//////////////////////////////////////////////////////////////////////////// +#include "MHHadronness.h" + +#include +#include +#include +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" +#include "MBinning.h" +#include "MHadronness.h" + +#include "MMcEvt.hxx" + +ClassImp(MHHadronness); + +// -------------------------------------------------------------------------- +// +// Setup histograms, nbins is the number of bins used for the evaluation. +// The default is 100 bins. +// +MHHadronness::MHHadronness(Int_t nbins, const char *name, const char *title) +{ + // + // set the name and title of this object + // + fName = name ? name : "MHHadronness"; + fTitle = title ? title : "Gamma/Hadron Separation Quality Histograms"; + + fGraph = new TGraph; + fGraph->SetTitle("Acceptance Gammas vs. Hadrons"); + fGraph->SetMaximum(1); + + fGhness = new TH1D("Ghness", "Hadronness", nbins, 0, 1); + fPhness = new TH1D("Phness", "Hadronness", nbins, 0, 1); + fGhness->SetXTitle("Hadronness"); + fPhness->SetXTitle("Hadronness"); + fGhness->SetYTitle("Counts"); + fPhness->SetYTitle("Counts"); + + fIntGhness = new TH1D("AccGammas", "Acceptance", nbins, 0, 1); + fIntPhness = new TH1D("AccHadrons", "Acceptance", nbins, 0, 1); + fIntGhness->SetXTitle("Hadronness"); + fIntPhness->SetXTitle("Hadronness"); + fIntGhness->SetYTitle("Acceptance"); + fIntPhness->SetYTitle("Acceptance"); + + /* + fQfac = new TH1D("Qfac", "Naive Quality factor", nbins, 0, 1); + fQfac->SetXTitle("Hadronness"); + fQfac->SetYTitle("Quality"); + */ + fQfac = new TGraph; + fQfac->SetTitle(" Naive Quality factor "); + + fMinDist = new TH1D("MinDist", "Minimum Dist to (0, 1)", nbins, 0, 1); + fMinDist->SetXTitle("Hadronness"); + fMinDist->SetYTitle("Distance"); + + //fQfac->SetDirectory(NULL); + fGhness->SetDirectory(NULL); + fPhness->SetDirectory(NULL); + fMinDist->SetDirectory(NULL); + fIntGhness->SetDirectory(NULL); + fIntPhness->SetDirectory(NULL); +} + +// -------------------------------------------------------------------------- +// +// Delete the histograms. +// +MHHadronness::~MHHadronness() +{ + delete fGhness; + delete fIntGhness; + delete fPhness; + delete fIntPhness; + delete fQfac; + delete fMinDist; + delete fGraph; +} + +// -------------------------------------------------------------------------- +// +// Setup Filling of the histograms. It needs: +// MMcEvt and MHadronness +// +Bool_t MHHadronness::SetupFill(const MParList *plist) +{ + fMcEvt = (MMcEvt*)plist->FindObject("MMcEvt"); + if (!fMcEvt) + { + *fLog << err << dbginf << "MMcEvt not found... aborting." << endl; + return kFALSE; + } + + fHadronness = (MHadronness*)plist->FindObject("MHadronness"); + if (!fHadronness) + { + *fLog << err << dbginf << "MHadronness not found... aborting." << endl; + return kFALSE; + } + + /* + MBinning* bins = (MBinning*)plist->FindObject("BinningHadronness"); + if (!bins) + { + *fLog << err << dbginf << "BinningHadronness [MBinning] not found... aborting." << endl; + return kFALSE; + } + + SetBinning(&fHist, binsalpha, binsenergy, binstheta); + + fHist.Sumw2(); + */ + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Fill the Hadronness from a MHadronness container into the corresponding +// histogram dependant on the particle id. +// +Bool_t MHHadronness::Fill(const MParContainer *par) +{ + if (fMcEvt->GetPartId()==kGAMMA) + fGhness->Fill(fHadronness->GetHadronness()); + else + fPhness->Fill(fHadronness->GetHadronness()); + + return kTRUE; +} + +Float_t MHHadronness::GetQ05() const +{ + Int_t n = fQfac->GetN(); + + Double_t val1x=0; + Double_t val2x=1; + + Double_t val1y=0; + Double_t val2y=0; + + for (Int_t i=1; i<=n; i++) + { + Double_t x, y; + + fQfac->GetPoint(i, x, y); + + if (x<0.5 && x>val1x) + { + val1x = x; + val1y = y; + } + + if (x>0.5 && x acceptance +// - fill the Minimum Distance histogram (formular see class description) +// - fill the Quality histogram (formular see class description) +// +Bool_t MHHadronness::Finalize() +{ + Int_t n = fGhness->GetNbinsX(); + + fGraph->Set(n); + fQfac->Set(n); + + const Stat_t sumg = fGhness->Integral(1, n+1); + const Stat_t sump = fPhness->Integral(1, n+1); + + Float_t max=0; + + for (Int_t i=1; i<=n; i++) + { + const Stat_t ip = fPhness->Integral(1, i)/sump; + const Stat_t ig = fGhness->Integral(1, i)/sumg; + + fIntPhness->SetBinContent(i, ip); + fIntGhness->SetBinContent(i, ig); + + fGraph->SetPoint(i, ip, ig); + + if (ip<=0) + continue; + + fMinDist->SetBinContent(i, 1-sqrt(ip*ip + (ig-1)*(ig-1))); + + Double_t val = ig/sqrt(ip); + fQfac->SetPoint(i, ig, val); + + if (val>max) + max = val; + } + + fQfac->SetMaximum(max+1); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Search the corresponding points for the given hadron acceptance (acchad) +// and interpolate the tow points (linear) +// +Double_t MHHadronness::GetGammaAcceptance(Double_t acchad) const +{ + const Int_t n = fGraph->GetN(); + const Double_t *x = fGraph->GetX(); + const Double_t *y = fGraph->GetY(); + + Int_t i = 0; + while (iGetN(); + const Double_t *x = fGraph->GetX(); + const Double_t *y = fGraph->GetY(); + + Int_t i = 0; + while (iGetEntries() << " Gammas and " << fPhness->GetEntries() << " Hadrons." << endl; + *fLog << "Acc Gammas @ 1% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.01)*100) << "%" << endl; + *fLog << "Acc Gammas @ 2% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.02)*100) << "%" << endl; + *fLog << "Acc Gammas @ 5% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.05)*100) << "%" << endl; + *fLog << "Acc Gammas @ 10% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.1)*100) << "%" << endl; + *fLog << "Acc Gammas @ 20% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.2)*100) << "%" << endl; + *fLog << "Acc Gammas @ 30% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.3)*100) << "%" << endl; + *fLog << "Acc Gammas @ 40% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.4)*100) << "%" << endl; + *fLog << "Acc Gammas @ 50% Hadron-acc: " << Form("%3.0f", GetGammaAcceptance(0.5)*100) << "%" << endl; + + const Int_t h = fMinDist->GetMaximumBin(); + *fLog << "Minimum Distance to (0, 1): " << Form("%.2f", 1-fMinDist->GetMaximum()) << " @ H=" << fMinDist->GetBinCenter(h) << endl; + *fLog << " Acc Gammas = " << Form("%3.0f", fIntGhness->GetBinContent(h)*100) << "%, "; + *fLog << "Acc Hadrons = " << Form("%3.0f", fIntPhness->GetBinContent(h)*100) << "%" << endl; + + *fLog << "Q-Factor @ Acc Gammas=0.5: Q(0.5)=" << Form("%.1f", GetQ05()) << endl; + *fLog << " Acc Hadrons = " << Form("%5.1f", GetHadronAcceptance(0.5)*100) << "%" << endl; + *fLog << endl; +} + +// -------------------------------------------------------------------------- +// +// Draw clone of all histograms. (For the Meaning see class description) +// +TObject *MHHadronness::DrawClone(Option_t *opt) const +{ + TCanvas &c = *MakeDefCanvas("Hadronness", fTitle); + c.Divide(2, 2); + + gROOT->SetSelectedPad(NULL); + + c.cd(1); + gStyle->SetOptStat(10); + Getghness()->DrawCopy(); + Getphness()->SetLineColor(kRed); + Getphness()->DrawCopy("same"); + c.cd(2); + gStyle->SetOptStat(0); + Getighness()->DrawCopy(); + Getiphness()->SetLineColor(kRed); + Getiphness()->DrawCopy("same"); + fMinDist->SetLineColor(kBlue); + fMinDist->DrawCopy("same"); + c.cd(3); + //GetQfac()->DrawCopy(); + TGraph &g2 = (TGraph&)*fQfac->DrawClone("A*"); + g2.SetBit(kCanDelete); + gPad->Modified(); + gPad->Update(); + if (g2.GetHistogram()) + { + g2.GetXaxis()->SetRangeUser(0, 1); + g2.GetXaxis()->SetTitle("Acceptance Gammas"); + g2.GetYaxis()->SetTitle("Quality"); + g2.SetMarkerStyle(kFullDotSmall); + g2.Draw("P"); + + gPad->Modified(); + gPad->Update(); + } + c.cd(4); + gPad->Modified(); + gPad->Update(); + TGraph &g = (TGraph&)*fGraph->DrawClone("AC"); + g.SetBit(kCanDelete); + gPad->Modified(); + gPad->Update(); + if (g.GetHistogram()) + { + g.GetXaxis()->SetRangeUser(0, 1); + g.GetXaxis()->SetTitle("Acceptance Hadrons"); + g.GetYaxis()->SetTitle("Acceptance Gammas"); + g.SetMarkerStyle(kFullDotSmall); + g.Draw("P"); + + gPad->Modified(); + gPad->Update(); + } + const Int_t h = fMinDist->GetMaximumBin(); + TMarker *m = new TMarker(fIntPhness->GetBinContent(h), + fIntGhness->GetBinContent(h), kStar); + m->SetMarkerColor(kBlue); + m->SetBit(kCanDelete); + m->Draw(); + return &c; +} + +// -------------------------------------------------------------------------- +// +// Draw all histograms. (For the Meaning see class description) +// +void MHHadronness::Draw(Option_t *) +{ + if (!gPad) + MakeDefCanvas("Hadronness", fTitle); + + gPad->Divide(2, 2); + + gPad->cd(1); + gStyle->SetOptStat(10); + Getghness()->Draw(); + Getphness()->SetLineColor(kRed); + Getphness()->Draw("same"); + gPad->cd(2); + gStyle->SetOptStat(0); + Getighness()->Draw(); + Getiphness()->SetLineColor(kRed); + Getiphness()->Draw("same"); + fMinDist->SetLineColor(kBlue); + fMinDist->Draw("same"); + gPad->cd(3); + //GetQfac()->Draw(); + fQfac->Draw("A*"); + gPad->Modified(); + gPad->Update(); + if (fQfac->GetHistogram()) + { + fQfac->GetXaxis()->SetRangeUser(0, 1); + fQfac->GetXaxis()->SetTitle("Acceptance Gammas"); + fQfac->GetYaxis()->SetTitle("Quality"); + fQfac->SetMarkerStyle(kFullDotSmall); + fQfac->Draw("P"); + + gPad->Modified(); + gPad->Update(); + } + gPad->cd(4); + gPad->Modified(); + gPad->Update(); + fGraph->Draw("AC"); + gPad->Modified(); + gPad->Update(); + if (fGraph->GetHistogram()) + { + fGraph->GetXaxis()->SetRangeUser(0, 1); + fGraph->GetXaxis()->SetTitle("Acceptance Hadrons"); + fGraph->GetYaxis()->SetTitle("Acceptance Gammas"); + fGraph->SetMarkerStyle(kFullDotSmall); + fGraph->Draw("P"); + gPad->Modified(); + gPad->Update(); + } + const Int_t h = fMinDist->GetMaximumBin(); + TMarker *m = new TMarker(fIntPhness->GetBinContent(h), + fIntGhness->GetBinContent(h), kStar); + m->SetMarkerColor(kBlue); + m->SetBit(kCanDelete); + m->Draw(); +} Index: trunk/MagicSoft/Mars/mhist/MHHadronness.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHadronness.h (revision 1557) +++ trunk/MagicSoft/Mars/mhist/MHHadronness.h (revision 1557) @@ -0,0 +1,56 @@ +#ifndef MARS_MHHadronness +#define MARS_MHHadronness + +#ifndef MARS_MH +#include "MH.h" +#endif + +class TH1D; +class TGraph; +class MParList; +class MMcEvt; +class MHadronness; + +class MHHadronness : public MH +{ +private: + const MMcEvt *fMcEvt; //! + const MHadronness *fHadronness; //! + + TH1D* fPhness; // Hadrons Hadronness + TH1D* fGhness; // Gammas Hadronness + TH1D* fIntPhness; // Hadrons Acceptance + TH1D* fIntGhness; // Gammas Acceptance + TH1D* fMinDist; // Minimum Distance to optimum acceptance + + TGraph *fQfac; // Quality factor + TGraph *fGraph; // gamma acceptance vs. hadron acceptance + +public: + MHHadronness(Int_t nbins=100, const char *name=NULL, const char *title=NULL); + ~MHHadronness(); + + Double_t GetGammaAcceptance(Double_t acchad) const; + Double_t GetHadronAcceptance(Double_t accgam) const; + + TH1D *Getghness() const { return fGhness; } + TH1D *Getphness() const { return fPhness; } + TH1D *Getighness() const { return fIntGhness; } + TH1D *Getiphness() const { return fIntPhness; } + //TH2D *GetHist() const { return fHist; } + + Float_t GetQ05() const; + + Bool_t SetupFill(const MParList *plist); + Bool_t Fill(const MParContainer *par); + Bool_t Finalize(); + + void Print(Option_t *option="") const; + + void Draw(Option_t *opt=""); + TObject *DrawClone(Option_t *opt="") const; + + ClassDef(MHHadronness, 1) // Gamma/Hadron Separation Quality Histograms +}; + +#endif Index: trunk/MagicSoft/Mars/mhist/Makefile =================================================================== --- trunk/MagicSoft/Mars/mhist/Makefile (revision 1556) +++ trunk/MagicSoft/Mars/mhist/Makefile (revision 1557) @@ -49,5 +49,5 @@ MHTimeDiffTheta.cc \ MHCompProb.cc \ - MHHadroness.cc \ + MHHadronness.cc \ MHMcEnergy.cc \ MHMcEfficiency.cc \