Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 1904) +++ trunk/MagicSoft/Mars/Changelog (revision 1905) @@ -1,3 +1,40 @@ -*-*- END OF LINE -*-*- + + 2003/04/04: Wolfgang Wittek + + * mfilter/Makefile + /FilterLinkDef.h + - added MFCT1SelBasic, MFCT1SelStandard, MFCT1SelFinal, + MFCT1Supercuts + + * mfilter/MFCT1SelBasic.[h,cc] + /MFCT1SelStandard.[h,cc] + /MFCT1SelFinal.[h,cc] + - new; they replace manalysis/MSelBasic.[h,cc] + /MSelStandard.[h,cc] + /MSelFinal.[h,cc] + + * mfilter/MFCT1Supercuts.[h,cc] + - new (supercuts from Daniel for CT1 Mkn421 2001 data) + + * manalysis/MSelBasic.[h,cc] + /MSelStandard.[h,cc] + /MSelFinal.[h,cc] + - deleted + + * manalysis/Makefile + AnalysisLinkDef.h + - MSelBasic, MSelStandard, MSelFinal removed + + * mfileio/MCT1ReadPreProc.[h,cc] + - smeare theta only for MC events + + // The following is not yet commited !!!!!!!!!!!!! + * Mars/Makefile + - in SUBDIRS put mfilter before manalysis + (otherwise get link error) + + + 2003/04/03: Abelardo Moralejo @@ -75,4 +112,5 @@ can be used in the parameters used in the random forest separation) + Index: trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1904) +++ trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1905) @@ -54,7 +54,4 @@ #pragma link C++ class MPadSchweizer+; -#pragma link C++ class MSelBasic+; -#pragma link C++ class MSelStandard+; -#pragma link C++ class MSelFinal+; #pragma link C++ class MPointingCorr+; @@ -71,2 +68,3 @@ + Index: trunk/MagicSoft/Mars/manalysis/MSelBasic.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelBasic.cc (revision 1904) +++ (revision ) @@ -1,249 +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): Wolfgang Wittek 02/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelBasic // -// // -// This is a task to evaluate basic cuts // -// // -// to be called after the calibration (when the number of photons is // -// available for all pixels) // -// // -///////////////////////////////////////////////////////////////////////////// - -#include "MSelBasic.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MCerPhotEvt.h" -#include "MMcEvt.hxx" -#include "MRawRunHeader.h" -#include "MGeomCam.h" -#include "MPedestalCam.h" -#include "MPedestalPix.h" -#include "MGeomPix.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MSelBasic); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MSelBasic::MSelBasic(const char *name, const char *title) -{ - fName = name ? name : "MSelBasic"; - fTitle = title ? title : "Task to evaluate basic cuts"; - - fThetaMin = 0.0; - fThetaMax = 60.0; -} - -// -------------------------------------------------------------------------- -// -// -// -// -// -Bool_t MSelBasic::PreProcess(MParList *pList) -{ - fRawRun = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); - if (!fRawRun) - { - *fLog << dbginf << "MRawRunHeader not found... aborting." << endl; - return kFALSE; - } - - fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); - if (!fMcEvt) - { - *fLog << dbginf << "MMcEvt not found... aborting." << endl; - return kFALSE; - } - - fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fEvt) - { - *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - - fCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fCam) - { - *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; - return kFALSE; - } - - fPed = (MPedestalCam*)pList->FindObject("MPedestalCam"); - if (!fPed) - { - *fLog << dbginf << "MPedestalCam missing in Parameter List... aborting." << endl; - return kFALSE; - } - - memset(fCut, 0, sizeof(fCut)); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Evaluate basic cuts -// -// if cuts are fulfilled : return 0 -// if they are not fullfilled : skip remaining tasks for this event -// -Bool_t MSelBasic::Process() -{ - Int_t rc = 0; - - //if ( fRawRun->GetRunNumber() < 16279 ) - //{ - // rc = 1; - // return kCONTINUE; - //} - - Double_t theta = kRad2Deg*fMcEvt->GetTelescopeTheta(); - if ( theta < fThetaMin ) - { - *fLog << "MSelBasic::Process; Run, Event, Theta = " - << fRawRun->GetRunNumber()<< ", " - << fMcEvt->GetEvtNumber() << ", " << theta << endl; - rc = 1; - } - - else if ( theta > fThetaMax ) - { - rc = 2; - } - - else if ( !SwTrigger() ) - { - //*fLog << "MSelBasic::Process; SwTrigger = " << SwTrigger << endl; - rc = 3; - } - - fCut[rc]++; - - return rc==0 ? kTRUE : kCONTINUE; -} -// -------------------------------------------------------------------------- -// -// Software trigger -// -// require 2 neighboring pixels (which are not in the outermost ring), -// each having at least 'minphotons' photons -// -// -Bool_t MSelBasic::SwTrigger() -{ - // minimum number of photons required - Double_t minphotons = 13.0; - - const Int_t entries = fEvt->GetNumPixels(); - - for (Int_t i=0; iIsPixelUsed(id1) ) continue; - - MGeomPix &gpix1 = (*fCam)[id1]; - if ( gpix1.IsInOutermostRing() ) continue; - - //MCerPhotPix &pix1 = (*fEvt)[id1]; - MCerPhotPix &pix1 = *(fEvt->GetPixById(id1)); - if ( &pix1 == NULL ) - { - *fLog << "MSelBasic::SwTrigger; &pix1 is NULL" << endl; - continue; - } - - Double_t photons1 = pix1.GetNumPhotons(); - if (photons1 >= minphotons) return kTRUE; - } - } - return kFALSE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the Basic selections. -// -Bool_t MSelBasic::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) - << (int)(fCut[1]*100/GetNumExecutions()) - << "%) Evts skipped due to: Zenith angle < " << fThetaMin << endl; - - *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) - << (int)(fCut[2]*100/GetNumExecutions()) - << "%) Evts skipped due to: Zenith angle > " << fThetaMax << endl; - - *fLog << " " << setw(7) << fCut[3] << " (" << setw(3) - << (int)(fCut[3]*100/GetNumExecutions()) - << "%) Evts skipped due to: Software trigger not fullfilled" << endl; - - *fLog << " " << fCut[0] << " (" << (int)(fCut[0]*100/GetNumExecutions()) - << "%) Evts survived Basic selections!" << endl; - *fLog << endl; - - return kTRUE; -} - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MSelBasic.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelBasic.h (revision 1904) +++ (revision ) @@ -1,62 +1,0 @@ -#ifndef MARS_MSelBasic -#define MARS_MSelBasic - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelBasic // -// // -// Task to evaluate basic cuts // -// // -///////////////////////////////////////////////////////////////////////////// - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MGeomCam; -class MPedestalCam; -class MCerPhotEvt; -class MHillas; -class MMcEvt; -class MRawRunHeader; - -class MSelBasic : public MTask -{ -private: - const MPedestalCam *fPed; // Pedestal information - const MGeomCam *fCam; // Camera Geometry - const MCerPhotEvt *fEvt; // Cerenkov Photon Event - const MMcEvt *fMcEvt; - const MRawRunHeader *fRawRun; - - Float_t fThetaMin; - Float_t fThetaMax; - - Int_t fCut[4]; - -public: - MSelBasic(const char *name=NULL, const char *title=NULL); - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - Bool_t PostProcess(); - - Bool_t SwTrigger(); - void SetCuts(Float_t thetamin, Float_t thetamax) - { fThetaMin = thetamin; fThetaMax = thetamax; } - - ClassDef(MSelBasic, 0) // Task to evaluate basic cuts -}; - -#endif - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MSelFinal.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelFinal.cc (revision 1904) +++ (revision ) @@ -1,185 +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): Wolfgang Wittek 02/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelFinal // -// // -// This is a task to evaluate the Final Cuts // -// (these cuts define the final sample of gammas; // -// relevant for the calculation of the effective collection areas) // -// // -// to be called after the calculation of the hadroness // -// // -///////////////////////////////////////////////////////////////////////////// - -#include "math.h" -#include "MSelFinal.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MHillasExt.h" -#include "MHillasSrc.h" -#include "MCerPhotEvt.h" -#include "MMcEvt.hxx" -#include "MGeomCam.h" -#include "MGeomPix.h" -#include "MHadronness.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MSelFinal); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MSelFinal::MSelFinal(const char *HilName, const char *HilSrcName, - const char *name, const char *title) -{ - fName = name ? name : "MSelFinal"; - fTitle = title ? title : "Task to evaluate the Final Cuts"; - - fHilName = HilName; - fHilSrcName = HilSrcName; - - // default values of cuts - fHadronnessMax = 1.0; - fAlphaMax = 100.0; //degrees -} - -// -------------------------------------------------------------------------- -// -// -// -// -// -Bool_t MSelFinal::PreProcess(MParList *pList) -{ - fHil = (MHillasExt*)pList->FindObject(fHilName, "MHillasExt"); - if (!fHil) - { - *fLog << dbginf << "MHillasExt object " << fHilName << " not found... aborting." << endl; - return kFALSE; - } - - fHilSrc = (MHillasSrc*)pList->FindObject(fHilSrcName, "MHillasSrc"); - if (!fHilSrc) - { - *fLog << dbginf << "MHillasSrc object " << fHilSrcName << " not found... aborting." << endl; - return kFALSE; - } - - - fHadronness = (MHadronness*)pList->FindObject("MHadronness"); - if (!fHadronness) - { - *fLog << dbginf << "MHadronness not found... aborting." << endl; - return kFALSE; - } - - - fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); - if (!fMcEvt) - { - *fLog << dbginf << "MMcEvt not found... aborting." << endl; - return kFALSE; - } - - memset(fCut, 0, sizeof(fCut)); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Evaluate final cuts -// -// if cuts are fulfilled : return 0 -// if they are not fullfilled : skip the remaining tasks for this event -// -Bool_t MSelFinal::Process() -{ - //*fLog << "Entry MSelFinal; fHilSrc = " << fHilSrc << endl; - - - - Int_t rc = 0; - - Double_t modalpha = fabs( fHilSrc->GetAlpha() ); - - Double_t h = fHadronness->GetHadronness(); - - if ( h>fHadronnessMax ) - { - //*fLog << "MSelFinal::Process; h, alpha = " << h << ", " - // << fHilSrc->GetAlpha() << endl; - rc = 1; - } - - else if ( modalpha > fAlphaMax ) - { - //*fLog << "MSelFinal::Process; h, alpha = " << h << ", " - // << fHilSrc->GetAlpha() << endl; - rc = 2; - } - - fCut[rc]++; - - return rc==0 ? kTRUE : kCONTINUE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the Final selections. -// -Bool_t MSelFinal::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) - << (int)(fCut[1]*100/GetNumExecutions()) - << "%) Evts skipped due to: hadroness > "<< fHadronnessMax << endl; - - *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) - << (int)(fCut[2]*100/GetNumExecutions()) - << "%) Evts skipped due to: |ALPHA| > " << fAlphaMax - << " [degrees]" << endl; - - *fLog << " " << fCut[0] << " (" - << (int)(fCut[0]*100/GetNumExecutions()) - << "%) Evts survived Final selections!" << endl; - *fLog << endl; - - return kTRUE; -} - - - Index: trunk/MagicSoft/Mars/manalysis/MSelFinal.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelFinal.h (revision 1904) +++ (revision ) @@ -1,66 +1,0 @@ -#ifndef MARS_MSelFinal -#define MARS_MSelFinal - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelFinal // -// // -// Task to evaluate final cuts // -// // -///////////////////////////////////////////////////////////////////////////// - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MGeomCam; -class MCerPhotEvt; -class MHillas; -class MHillasSrc; -class MMcEvt; -class MHadronness; - -class MSelFinal : public MTask -{ -private: - MGeomCam *fCam; // Camera Geometry - MCerPhotEvt *fEvt; // Cerenkov Photon Event - MMcEvt *fMcEvt; - MHillas *fHil; - MHillasSrc *fHilSrc; - MHadronness *fHadronness; - - Double_t fMm2Deg; // conversion mm to degrees in camera - Int_t fCut[3]; - TString fHilName; - TString fHilSrcName; - - Float_t fHadronnessMax; - Float_t fAlphaMax; - -public: - MSelFinal(const char *HilName="MHillas", const char *HilSrcName="MHillasSrc", - const char *name=NULL, const char *title=NULL); - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - Bool_t PostProcess(); - - void SetCuts(Float_t hadmax, Float_t alphamax) - { fHadronnessMax = hadmax; fAlphaMax = alphamax; } - - ClassDef(MSelFinal, 0) // Task to evaluate final cuts -}; - -#endif - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MSelStandard.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelStandard.cc (revision 1904) +++ (revision ) @@ -1,225 +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): Wolfgang Wittek 02/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelStandard // -// // -// This is a task to evaluate the Standard Cuts // -// // -// to be called after the calculation of the image parameters // -// before the g/h separation // -// // -///////////////////////////////////////////////////////////////////////////// - -#include "MSelStandard.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MHillasExt.h" -#include "MHillasSrc.h" -#include "MCerPhotEvt.h" -#include "MMcEvt.hxx" -#include "MGeomCam.h" -#include "MGeomPix.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MSelStandard); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MSelStandard::MSelStandard(const char *HilName, const char *HilSrcName, - const char *name, const char *title) -{ - fName = name ? name : "MSelStandard"; - fTitle = title ? title : "Task to evaluate the Standard Cuts"; - - fHilName = HilName; - fHilSrcName = HilSrcName; - - // default values of cuts - fUsedPixelsMax = 92; - fCorePixelsMin = 4; - fSizeMin = 60; - fDistMin = 0.4; - fDistMax = 1.05; - fLengthMin = 0.0; - fWidthMin = 0.0; -} - -// -------------------------------------------------------------------------- -// -// -// -// -// -Bool_t MSelStandard::PreProcess(MParList *pList) -{ - fHil = (MHillasExt*)pList->FindObject(fHilName, "MHillasExt"); - if (!fHil) - { - *fLog << dbginf << "MHillasExt object " << fHilName << " not found... aborting." << endl; - return kFALSE; - } - - fHilSrc = (MHillasSrc*)pList->FindObject(fHilSrcName, "MHillasSrc"); - if (!fHilSrc) - { - *fLog << dbginf << "MHillasSrc object " << fHilSrcName << " not found... aborting." << endl; - return kFALSE; - } - - - fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); - if (!fMcEvt) - { - *fLog << dbginf << "MMcEvt not found... aborting." << endl; - return kFALSE; - } - - fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fEvt) - { - *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - - fCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fCam) - { - *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; - return kFALSE; - } - fMm2Deg = fCam->GetConvMm2Deg(); - - //*fLog << "fMm2Deg = " << fMm2Deg << endl; - - memset(fCut, 0, sizeof(fCut)); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Evaluate standard cuts -// -// if cuts are fulfilled : return 0 -// if they are not fullfilled : skip the remaining tasks for this event -// -// -Bool_t MSelStandard::Process() -{ - Int_t rc = 0; - - Double_t length = fHil->GetLength() * fMm2Deg; - Double_t width = fHil->GetWidth() * fMm2Deg; - Double_t dist = fHilSrc->GetDist()* fMm2Deg; - //Double_t delta = fHil->GetDelta() * kRad2Deg; - Double_t size = fHil->GetSize(); - Int_t numusedpixels = fHil->GetNumUsedPixels(); - Int_t numcorepixels = fHil->GetNumCorePixels(); - - if ( numusedpixels >= fUsedPixelsMax || numcorepixels <= fCorePixelsMin ) - { - rc = 1; - } - - else if ( size <= fSizeMin ) - { - rc = 2; - } - - else if ( dist< fDistMin || dist > fDistMax ) - { - rc = 3; - } - - else if ( length <= fLengthMin || width <= fWidthMin ) - { - rc = 4; - } - - - fCut[rc]++; - - return rc==0 ? kTRUE : kCONTINUE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the Standard selections. -// -Bool_t MSelStandard::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) - << (int)(fCut[1]*100/GetNumExecutions()) - << "%) Evts skipped due to: Used pixels >= " << fUsedPixelsMax - << " or Core pixels <= " << fCorePixelsMin << endl; - - *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) - << (int)(fCut[2]*100/GetNumExecutions()) - << "%) Evts skipped due to: SIZE <= " << fSizeMin << endl; - - *fLog << " " << setw(7) << fCut[3] << " (" << setw(3) - << (int)(fCut[3]*100/GetNumExecutions()) - << "%) Evts skipped due to: DIST < " << fDistMin - << " or DIST > " << fDistMax << endl; - - *fLog << " " << setw(7) << fCut[4] << " (" << setw(3) - << (int)(fCut[4]*100/GetNumExecutions()) - << "%) Evts skipped due to: LENGTH <= " << fLengthMin - << " or WIDTH <= " << fWidthMin << endl; - - *fLog << " " << fCut[0] << " (" - << (int)(fCut[0]*100/GetNumExecutions()) - << "%) Evts survived Standard selections!" << endl; - *fLog << endl; - - return kTRUE; -} - - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MSelStandard.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MSelStandard.h (revision 1904) +++ (revision ) @@ -1,73 +1,0 @@ -#ifndef MARS_MSelStandard -#define MARS_MSelStandard - -///////////////////////////////////////////////////////////////////////////// -// // -// MSelStandard // -// // -// Task to evaluate standard cuts // -// // -///////////////////////////////////////////////////////////////////////////// - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MGeomCam; -class MCerPhotEvt; -class MHillas; -class MHillasSrc; -class MMcEvt; - -class MSelStandard : public MTask -{ -private: - MGeomCam *fCam; // Camera Geometry - MCerPhotEvt *fEvt; // Cerenkov Photon Event - MMcEvt *fMcEvt; - MHillas *fHil; - MHillasSrc *fHilSrc; - - Double_t fMm2Deg; // conversion mm to degrees in camera - Int_t fCut[5]; - TString fHilName; - TString fHilSrcName; - - Float_t fUsedPixelsMax; - Float_t fCorePixelsMin; - Float_t fSizeMin; - Float_t fDistMin; - Float_t fDistMax; - Float_t fLengthMin; - Float_t fWidthMin; - -public: - MSelStandard(const char *HilName="MHillas", const char *HilSrcName="MHillasSrc", - const char *name=NULL, const char *title=NULL); - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - Bool_t PostProcess(); - - void SetCuts(Float_t usedpixelsmax, Float_t corepixelsmin, - Float_t sizemin, Float_t distmin, Float_t distmax, - Float_t lengthmin, Float_t widthmin) - { fUsedPixelsMax = usedpixelsmax; fCorePixelsMin = corepixelsmin; - fSizeMin = sizemin; fDistMin = distmin; fDistMax = distmax; - fLengthMin = lengthmin; fWidthMin = widthmin; } - - ClassDef(MSelStandard, 0) // Task to evaluate standard cuts -}; - -#endif - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/Makefile =================================================================== --- trunk/MagicSoft/Mars/manalysis/Makefile (revision 1904) +++ trunk/MagicSoft/Mars/manalysis/Makefile (revision 1905) @@ -65,7 +65,4 @@ MPadding.cc \ MPadSchweizer.cc \ - MSelBasic.cc \ - MSelStandard.cc \ - MSelFinal.cc \ MPointingCorr.cc \ MNewImagePar.cc \ Index: trunk/MagicSoft/Mars/mfileio/MCT1ReadPreProc.cc =================================================================== --- trunk/MagicSoft/Mars/mfileio/MCT1ReadPreProc.cc (revision 1904) +++ trunk/MagicSoft/Mars/mfileio/MCT1ReadPreProc.cc (revision 1905) @@ -912,4 +912,7 @@ fThetaOrig->SetVal(theta); + //*fLog << "Theta, smearedTheta = " << theta << ", " << SmearTheta(theta) + // << endl; + // store hour angle fHourAngle->SetVal(event.fhourangle); @@ -931,5 +934,5 @@ fIsMcFile ? event.imcimpact_m*100 : 0, TMath::Pi()/180*event.iaz_arcs/3600, // azimuth (arcseconds) - SmearTheta(theta), + fIsMcFile ? SmearTheta(theta) : theta, 0, /* fTFirst */ 0, /* fTLast */ Index: trunk/MagicSoft/Mars/mfilter/FilterLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mfilter/FilterLinkDef.h (revision 1904) +++ trunk/MagicSoft/Mars/mfilter/FilterLinkDef.h (revision 1905) @@ -16,3 +16,8 @@ #pragma link C++ class MFDataMember+; +#pragma link C++ class MFCT1SelBasic+; +#pragma link C++ class MFCT1SelStandard+; +#pragma link C++ class MFCT1SelFinal+; +#pragma link C++ class MFCT1Supercuts+; + #endif Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.cc =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.cc (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.cc (revision 1905) @@ -0,0 +1,263 @@ +/* ======================================================================== *\ +! +! * +! * 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): Wolfgang Wittek 04/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1SelBasic // +// // +// This is a class to evaluate basic cuts // +// // +// to be called after the calibration (when the number of photons is // +// available for all pixels) // +// // +///////////////////////////////////////////////////////////////////////////// + +#include "MFCT1SelBasic.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MCerPhotEvt.h" +#include "MMcEvt.hxx" +#include "MRawRunHeader.h" +#include "MGeomCam.h" +#include "MPedestalCam.h" +#include "MPedestalPix.h" +#include "MGeomPix.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MFCT1SelBasic); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MFCT1SelBasic::MFCT1SelBasic(const char *name, const char *title) +{ + fName = name ? name : "MFCT1SelBasic"; + fTitle = title ? title : "Filter to evaluate basic cuts"; + + // default values of cuts + SetCuts(13.0, 0.0, 60.0); +} + +// -------------------------------------------------------------------------- +// +// Set the cut values +// +// +void MFCT1SelBasic::SetCuts(Float_t minphotons, + Float_t thetamin, Float_t thetamax) +{ + fMinPhotons = minphotons; + fThetaMin = thetamin; + fThetaMax = thetamax; + + *fLog << inf << "MFCT1SelBasic cut values : fMinPhotons, fThetaMin, fThetaMax = " + << fMinPhotons <<", " << fThetaMin << ", " << fThetaMax << endl; +} + +// -------------------------------------------------------------------------- +// +// Set the pointers +// +// +Bool_t MFCT1SelBasic::PreProcess(MParList *pList) +{ + fRawRun = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); + if (!fRawRun) + { + *fLog << dbginf << "MRawRunHeader not found... aborting." << endl; + return kFALSE; + } + + fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); + if (!fMcEvt) + { + *fLog << dbginf << "MMcEvt not found... aborting." << endl; + return kFALSE; + } + + fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + + fCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fCam) + { + *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; + return kFALSE; + } + + fPed = (MPedestalCam*)pList->FindObject("MPedestalCam"); + if (!fPed) + { + *fLog << dbginf << "MPedestalCam missing in Parameter List... aborting." << endl; + return kFALSE; + } + + memset(fCut, 0, sizeof(fCut)); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Evaluate basic cuts +// +// if cuts are fulfilled : fResult = kTRUE; +// if they are not fullfilled : fResult = kFALSE; +// +Bool_t MFCT1SelBasic::Process() +{ + Int_t rc = 0; + fResult = kFALSE; + + //if ( fRawRun->GetRunNumber() < 16279 ) + //{ + // rc = 1; + // fResult = kTRUE; + //} + + Double_t theta = kRad2Deg*fMcEvt->GetTelescopeTheta(); + if ( theta < fThetaMin ) + { + rc = 1; + fResult = kTRUE; + } + + else if ( theta > fThetaMax ) + { + rc = 2; + fResult = kTRUE; + } + + else if ( !SwTrigger() ) + { + //*fLog << "MFCT1SelBasic::Process; SwTrigger = " << SwTrigger << endl; + rc = 3; + fResult = kTRUE; + } + + fCut[rc]++; + + return kTRUE; +} +// -------------------------------------------------------------------------- +// +// Software trigger +// +// require 2 neighboring pixels (which are not in the outermost ring), +// each having at least 'fMinPhotons' photons +// +// +Bool_t MFCT1SelBasic::SwTrigger() +{ + const Int_t entries = fEvt->GetNumPixels(); + + for (Int_t i=0; iIsPixelUsed(id1) ) continue; + + MGeomPix &gpix1 = (*fCam)[id1]; + if ( gpix1.IsInOutermostRing() ) continue; + + MCerPhotPix &pix1 = *(fEvt->GetPixById(id1)); + + Double_t photons1 = pix1.GetNumPhotons(); + if (photons1 >= fMinPhotons) return kTRUE; + } + } + return kFALSE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the Basic selections. +// +Bool_t MFCT1SelBasic::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) + << (int)(fCut[1]*100/GetNumExecutions()) + << "%) Evts skipped due to: Zenith angle < " << fThetaMin << endl; + + *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) + << (int)(fCut[2]*100/GetNumExecutions()) + << "%) Evts skipped due to: Zenith angle > " << fThetaMax << endl; + + *fLog << " " << setw(7) << fCut[3] << " (" << setw(3) + << (int)(fCut[3]*100/GetNumExecutions()) + << "%) Evts skipped due to: Software trigger not fullfilled" + << " (with fMinPhotons = " << fMinPhotons << ")" << endl; + + *fLog << " " << fCut[0] << " (" << (int)(fCut[0]*100/GetNumExecutions()) + << "%) Evts survived Basic selections!" << endl; + *fLog << endl; + + return kTRUE; +} + + + + + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.h =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.h (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelBasic.h (revision 1905) @@ -0,0 +1,66 @@ +#ifndef MARS_MFCT1SelBasic +#define MARS_MFCT1SelBasic + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1SelBasic // +// // +// Class to evaluate basic cuts // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MFilter +#include "MFilter.h" +#endif + +class MGeomCam; +class MPedestalCam; +class MCerPhotEvt; +class MHillas; +class MMcEvt; +class MRawRunHeader; + +class MFCT1SelBasic : public MFilter +{ +private: + const MPedestalCam *fPed; // Pedestal information + const MGeomCam *fCam; // Camera Geometry + const MCerPhotEvt *fEvt; // Cerenkov Photon Event + const MMcEvt *fMcEvt; + const MRawRunHeader *fRawRun; + + Float_t fMinPhotons; + Float_t fThetaMin; + Float_t fThetaMax; + + Bool_t fResult; + + Int_t fCut[4]; + +public: + MFCT1SelBasic(const char *name=NULL, const char *title=NULL); + + Bool_t IsExpressionTrue() const { return fResult; } + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + Bool_t PostProcess(); + + Bool_t SwTrigger(); + void SetCuts(Float_t minphotons, Float_t thetamin, Float_t thetamax); + + ClassDef(MFCT1SelBasic, 0) // Class to evaluate basic cuts +}; + +#endif + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.cc =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.cc (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.cc (revision 1905) @@ -0,0 +1,186 @@ +/* ======================================================================== *\ +! +! * +! * 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): Wolfgang Wittek 04/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1SelFinal // +// // +// This is a class to evaluate the Final Cuts // +// (these cuts define the final sample of gammas; // +// relevant for the calculation of the effective collection areas) // +// // +// to be called after the calculation of the hadroness // +// // +///////////////////////////////////////////////////////////////////////////// + +#include "math.h" +#include "MFCT1SelFinal.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MHillasExt.h" +#include "MHillasSrc.h" +#include "MCerPhotEvt.h" +#include "MMcEvt.hxx" +#include "MGeomCam.h" +#include "MGeomPix.h" +#include "MHadronness.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MFCT1SelFinal); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MFCT1SelFinal::MFCT1SelFinal(const char *HilName, const char *HilSrcName, + const char *name, const char *title) +{ + fName = name ? name : "MFCT1SelFinal"; + fTitle = title ? title : "Class to evaluate the Final Cuts"; + + fHilName = HilName; + fHilSrcName = HilSrcName; + + // default values of cuts + SetCuts(1.0, 100.0); +} + +// -------------------------------------------------------------------------- +// +// Set cut values +// +// +void MFCT1SelFinal::SetCuts(Float_t hadmax, Float_t alphamax) +{ + fHadronnessMax = hadmax; + fAlphaMax = alphamax; + + *fLog << inf << "MFCT1SelFinal cut values : fHadronnessMax, fAlphaMax = " + << fHadronnessMax << ", " << fAlphaMax << endl; +} + +// -------------------------------------------------------------------------- +// +// Set the pointers +// +// +Bool_t MFCT1SelFinal::PreProcess(MParList *pList) +{ + fHil = (MHillasExt*)pList->FindObject(fHilName, "MHillasExt"); + if (!fHil) + { + *fLog << dbginf << "MHillasExt object " << fHilName << " not found... aborting." << endl; + return kFALSE; + } + + fHilSrc = (MHillasSrc*)pList->FindObject(fHilSrcName, "MHillasSrc"); + if (!fHilSrc) + { + *fLog << dbginf << "MHillasSrc object " << fHilSrcName << " not found... aborting." << endl; + return kFALSE; + } + + fHadronness = (MHadronness*)pList->FindObject("MHadronness"); + if (!fHadronness) + { + *fLog << dbginf << "MHadronness not found... aborting." << endl; + return kFALSE; + } + + fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); + if (!fMcEvt) + { + *fLog << dbginf << "MMcEvt not found... aborting." << endl; + return kFALSE; + } + + memset(fCut, 0, sizeof(fCut)); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Evaluate final cuts +// +// if cuts are fulfilled set fResult = kTRUE +// +Bool_t MFCT1SelFinal::Process() +{ + Int_t rc = 0; + fResult = kFALSE; + + Double_t modalpha = fabs( fHilSrc->GetAlpha() ); + + Double_t h = fHadronness->GetHadronness(); + + if ( h>fHadronnessMax ) + { + rc = 1; + fResult = kTRUE; + } + + else if ( modalpha > fAlphaMax ) + { + rc = 2; + fResult = kTRUE; + } + + fCut[rc]++; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the Final selections. +// +Bool_t MFCT1SelFinal::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) + << (int)(fCut[1]*100/GetNumExecutions()) + << "%) Evts skipped due to: hadroness > "<< fHadronnessMax << endl; + + *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) + << (int)(fCut[2]*100/GetNumExecutions()) + << "%) Evts skipped due to: |ALPHA| > " << fAlphaMax + << " [degrees]" << endl; + + *fLog << " " << fCut[0] << " (" + << (int)(fCut[0]*100/GetNumExecutions()) + << "%) Evts survived Final selections!" << endl; + *fLog << endl; + + return kTRUE; +} Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.h =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.h (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelFinal.h (revision 1905) @@ -0,0 +1,70 @@ +#ifndef MARS_MFCT1SelFinal +#define MARS_MFCT1SelFinal + +///////////////////////////////////////////////////////////////////////////// +// // +// MSelFinal // +// // +// Class to evaluate final cuts // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MFilter +#include "MFilter.h" +#endif + +class MGeomCam; +class MCerPhotEvt; +class MHillas; +class MHillasSrc; +class MMcEvt; +class MHadronness; + +class MFCT1SelFinal : public MFilter +{ +private: + MGeomCam *fCam; // Camera Geometry + MCerPhotEvt *fEvt; // Cerenkov Photon Event + MMcEvt *fMcEvt; + MHillas *fHil; + MHillasSrc *fHilSrc; + MHadronness *fHadronness; + + Double_t fMm2Deg; // conversion mm to degrees in camera + Int_t fCut[3]; + TString fHilName; + TString fHilSrcName; + + Float_t fHadronnessMax; + Float_t fAlphaMax; + + Bool_t fResult; + + +public: + MFCT1SelFinal(const char *HilName="MHillas", const char *HilSrcName="MHillasSrc", + const char *name=NULL, const char *title=NULL); + + Bool_t IsExpressionTrue() const { return fResult; } + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + Bool_t PostProcess(); + + void SetCuts(Float_t hadmax, Float_t alphamax); + + ClassDef(MFCT1SelFinal, 0) // Class to evaluate final cuts +}; + +#endif + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.cc =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.cc (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.cc (revision 1905) @@ -0,0 +1,234 @@ +/* ======================================================================== *\ +! +! * +! * 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): Wolfgang Wittek 04/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1SelStandard // +// // +// This is a class to evaluate the Standard Cuts // +// // +// to be called after the calculation of the image parameters // +// before the g/h separation // +// // +///////////////////////////////////////////////////////////////////////////// + +#include "MFCT1SelStandard.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MHillasExt.h" +#include "MHillasSrc.h" +#include "MCerPhotEvt.h" +#include "MMcEvt.hxx" +#include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MFCT1SelStandard); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MFCT1SelStandard::MFCT1SelStandard(const char *HilName, const char *HilSrcName, + const char *name, const char *title) +{ + fName = name ? name : "MFCT1SelStandard"; + fTitle = title ? title : "Class to evaluate the Standard Cuts"; + + fHilName = HilName; + fHilSrcName = HilSrcName; + + // default values of cuts + SetCuts(92, 4, 60, 0.4, 1.05, 0.0, 0.0); +} + +// -------------------------------------------------------------------------- +// +// Set the values for the cuts +// +// +void MFCT1SelStandard::SetCuts(Float_t usedpixelsmax, Float_t corepixelsmin, + Float_t sizemin, Float_t distmin, Float_t distmax, + Float_t lengthmin, Float_t widthmin) +{ + fUsedPixelsMax = usedpixelsmax; + fCorePixelsMin = corepixelsmin; + fSizeMin = sizemin; + fDistMin = distmin; + fDistMax = distmax; + fLengthMin = lengthmin; + fWidthMin = widthmin; + + *fLog << inf << "MFCT1SelStandard cut values : fUsedPixelsMax, fCorePixelsMin = " + << fUsedPixelsMax << ", " << fCorePixelsMin << endl; + *fLog << inf << " fSizeMin, fDistMin, fDistMax = " << fSizeMin + << ", " << fDistMin << ", " << fDistMax << endl; + *fLog << inf << " fLengthMin, fWidthMin = " << fLengthMin + << ", " << fWidthMin << endl; +} + +// -------------------------------------------------------------------------- +// +// +// +// +// +Bool_t MFCT1SelStandard::PreProcess(MParList *pList) +{ + fHil = (MHillasExt*)pList->FindObject(fHilName, "MHillasExt"); + if (!fHil) + { + *fLog << dbginf << "MHillasExt object " << fHilName << " not found... aborting." << endl; + return kFALSE; + } + + fHilSrc = (MHillasSrc*)pList->FindObject(fHilSrcName, "MHillasSrc"); + if (!fHilSrc) + { + *fLog << dbginf << "MHillasSrc object " << fHilSrcName << " not found... aborting." << endl; + return kFALSE; + } + + + fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); + if (!fMcEvt) + { + *fLog << dbginf << "MMcEvt not found... aborting." << endl; + return kFALSE; + } + + fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + + fCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fCam) + { + *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; + return kFALSE; + } + fMm2Deg = fCam->GetConvMm2Deg(); + + //*fLog << "fMm2Deg = " << fMm2Deg << endl; + + memset(fCut, 0, sizeof(fCut)); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Evaluate standard cuts +// +// if selections are fulfilled set fResult = kTRUE; +// +// +Bool_t MFCT1SelStandard::Process() +{ + Int_t rc = 0; + fResult = kFALSE; + + Double_t length = fHil->GetLength() * fMm2Deg; + Double_t width = fHil->GetWidth() * fMm2Deg; + Double_t dist = fHilSrc->GetDist()* fMm2Deg; + //Double_t delta = fHil->GetDelta() * kRad2Deg; + Double_t size = fHil->GetSize(); + Int_t numusedpixels = fHil->GetNumUsedPixels(); + Int_t numcorepixels = fHil->GetNumCorePixels(); + + if ( numusedpixels >= fUsedPixelsMax || numcorepixels <= fCorePixelsMin ) + { + rc = 1; + fResult = kTRUE; + } + + else if ( size <= fSizeMin ) + { + rc = 2; + fResult = kTRUE; + } + + else if ( dist< fDistMin || dist > fDistMax ) + { + rc = 3; + fResult = kTRUE; + } + + else if ( length <= fLengthMin || width <= fWidthMin ) + { + rc = 4; + fResult = kTRUE; + } + + fCut[rc]++; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the Standard selections. +// +Bool_t MFCT1SelStandard::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << setw(7) << fCut[1] << " (" << setw(3) + << (int)(fCut[1]*100/GetNumExecutions()) + << "%) Evts skipped due to: Used pixels >= " << fUsedPixelsMax + << " or Core pixels <= " << fCorePixelsMin << endl; + + *fLog << " " << setw(7) << fCut[2] << " (" << setw(3) + << (int)(fCut[2]*100/GetNumExecutions()) + << "%) Evts skipped due to: SIZE <= " << fSizeMin << endl; + + *fLog << " " << setw(7) << fCut[3] << " (" << setw(3) + << (int)(fCut[3]*100/GetNumExecutions()) + << "%) Evts skipped due to: DIST < " << fDistMin + << " or DIST > " << fDistMax << endl; + + *fLog << " " << setw(7) << fCut[4] << " (" << setw(3) + << (int)(fCut[4]*100/GetNumExecutions()) + << "%) Evts skipped due to: LENGTH <= " << fLengthMin + << " or WIDTH <= " << fWidthMin << endl; + + *fLog << " " << fCut[0] << " (" + << (int)(fCut[0]*100/GetNumExecutions()) + << "%) Evts survived Standard selections!" << endl; + *fLog << endl; + + return kTRUE; +} Index: trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.h =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.h (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1SelStandard.h (revision 1905) @@ -0,0 +1,74 @@ +#ifndef MARS_MFCT1SelStandard +#define MARS_MFCT1SelStandard + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1SelStandard // +// // +// Class to evaluate standard cuts // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MFilter +#include "MFilter.h" +#endif + +class MGeomCam; +class MCerPhotEvt; +class MHillas; +class MHillasSrc; +class MMcEvt; + +class MFCT1SelStandard : public MFilter +{ +private: + MGeomCam *fCam; // Camera Geometry + MCerPhotEvt *fEvt; // Cerenkov Photon Event + MMcEvt *fMcEvt; + MHillas *fHil; + MHillasSrc *fHilSrc; + + Double_t fMm2Deg; // conversion mm to degrees in camera + Int_t fCut[5]; + TString fHilName; + TString fHilSrcName; + + Float_t fUsedPixelsMax; + Float_t fCorePixelsMin; + Float_t fSizeMin; + Float_t fDistMin; + Float_t fDistMax; + Float_t fLengthMin; + Float_t fWidthMin; + + Bool_t fResult; + +public: + MFCT1SelStandard(const char *HilName="MHillas", const char *HilSrcName="MHillasSrc", + const char *name=NULL, const char *title=NULL); + + Bool_t IsExpressionTrue() const { return fResult; } + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + Bool_t PostProcess(); + + void SetCuts(Float_t usedpixelsmax, Float_t corepixelsmin, + Float_t sizemin, Float_t distmin, Float_t distmax, + Float_t lengthmin, Float_t widthmin); + + ClassDef(MFCT1SelStandard, 0) // Class to evaluate standard cuts +}; + +#endif + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.cc =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.cc (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.cc (revision 1905) @@ -0,0 +1,280 @@ +/* ======================================================================== *\ +! +! * +! * 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): Wolfgang Wittek 04/2003 + +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1Supercuts // +// // +///////////////////////////////////////////////////////////////////////////// +#include "MFCT1Supercuts.h" + +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" +#include "MHillasExt.h" +#include "MHillasSrc.h" +#include "MMcEvt.hxx" +#include "MCerPhotEvt.h" +#include "MGeomCam.h" + +ClassImp(MFCT1Supercuts); + + +// -------------------------------------------------------------------------- +// +MFCT1Supercuts::MFCT1Supercuts(const char *hilname, const char *hilsrcname, + const char *name, const char *title) +{ + fName = name ? name : "MFCT1Supercuts"; + fTitle = title ? title : "Class to evaluate the Supercuts"; + + fHilName = hilname; + fHilSrcName = hilsrcname; +} + +// -------------------------------------------------------------------------- +// +Bool_t MFCT1Supercuts::PreProcess(MParList *pList) +{ + fHil = (MHillasExt*)pList->FindObject(fHilName, "MHillasExt"); + if (!fHil) + { + *fLog << dbginf << "MHillasExt object " << fHilName << " not found... aborting." << endl; + return kFALSE; + } + + fHilSrc = (MHillasSrc*)pList->FindObject(fHilSrcName, "MHillasSrc"); + if (!fHilSrc) + { + *fLog << dbginf << "MHillasSrc object " << fHilSrcName << " not found... aborting." << endl; + return kFALSE; + } + + fMcEvt = (MMcEvt*)pList->FindObject("MMcEvt"); + if (!fMcEvt) + { + *fLog << dbginf << "MMcEvt not found... aborting." << endl; + return kFALSE; + } + + fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + fCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fCam) + { + *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; + return kFALSE; + } + fMm2Deg = fCam->GetConvMm2Deg(); + + //*fLog << "fMm2Deg = " << fMm2Deg << endl; + + + //--------------------------------- + // cut parameters + + lengthup[0] = 0.315585; + lengthup[1] = 0.001455; + lengthup[2] = 0.203198; + lengthup[3] = 0.005532; + lengthup[4] =-0.001670; + lengthup[5] =-0.020362; + lengthup[6] = 0.007388; + lengthup[7] =-0.013463; + + widthup[0] = 0.145412; + widthup[1] =-0.001771; + widthup[2] = 0.054462; + widthup[3] = 0.022280; + widthup[4] =-0.009893; + widthup[5] = 0.056353; + widthup[6] = 0.020711; + widthup[7] =-0.016703; + + distup[0] = 1.787943; + distup[1] = 0.; + distup[2] = 2.942310; + distup[3] = 0.199815; + distup[4] = 0.; + distup[5] = 0.249909; + distup[6] = 0.189697; + distup[7] = 0.; + + lengthlow[0] = 0.151530; + lengthlow[1] = 0.028323; + lengthlow[2] = 0.510707; + lengthlow[3] = 0.053089; + lengthlow[4] = 0.013708; + lengthlow[5] = 2.357993; + lengthlow[6] = 0.000080; + lengthlow[7] =-0.007157; + + widthlow[0] = 0.089187; + widthlow[1] =-0.006430; + widthlow[2] = 0.074442; + widthlow[3] = 0.003738; + widthlow[4] =-0.004256; + widthlow[5] =-0.014101; + widthlow[6] = 0.006126; + widthlow[7] =-0.002849; + + distlow[0] = 0.589406; + distlow[1] = 0.; + distlow[2] =-0.083964; + distlow[3] =-0.007975; + distlow[4] = 0.; + distlow[5] = 0.045374; + distlow[6] =-0.001750; + distlow[7] = 0.; + + asymup[0] = 0.061267; + asymup[1] = 0.014462; + asymup[2] = 0.014327; + asymup[3] = 0.014540; + asymup[4] = 0.013391; + asymup[5] = 0.012319; + asymup[6] = 0.010444; + asymup[7] = 0.008328; + + asymlow[0] =-0.012055; + asymlow[1] = 0.009157; + asymlow[2] = 0.005441; + asymlow[3] = 0.000399; + asymlow[4] = 0.001433; + asymlow[5] =-0.002050; + asymlow[6] =-0.000104; + asymlow[7] =-0.001188; + + alphaup[0] = 13.123440; + alphaup[1] = 0.; + alphaup[2] = 0.; + alphaup[3] = 0.; + alphaup[4] = 0.; + alphaup[5] = 0.; + alphaup[6] = 0.; + alphaup[7] = 0.; + //--------------------------------- + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Calculation of upper and lower limits +// +Double_t MFCT1Supercuts::CtsMCut(Double_t *a, Double_t ls, Double_t ct, + Double_t ls2, Double_t dd2) +{ + // define cut-function + // + // dNOMLOGSIZE = 4.1 (=log(60.0) + // dNOMCOSZA = 1.0 + // + // a: array of cut parameters + // ls: log(SIZE) - dNOMLOGSIZE + // ls2: ls^2 + // ct: Cos(ZA.) - dNOMCOSZA + // dd2: DIST^2 + + Double_t limit = *a + *(a+1) * dd2 + *(a+2) * ct + + ls * (*(a+3) + *(a+4) * dd2 + *(a+5) * ct) + + ls2 * (*(a+6) + *(a+7) * dd2); + + *fLog << "MFCT1Supercuts::Process; *a = " + << *a << ", " << *(a+1) << ", " << *(a+2) << ", " + << *(a+3) << ", " << *(a+4) << ", " << *(a+5) << ", " + << *(a+6) << ", " << *(a+7) << endl; + + *fLog << "MFCT1Supercuts::Process; ls, ls2, ct, dd2, limit = " << ls + << ", " << ls2 << ", " << ct << ", " << dd2 << ", " + << limit << endl; + + return limit; +} + +// --------------------------------------------------------------------------- +// +// Evaluate supercuts for CT1 Mkn421 2001 (Daniel Kranich) +// +Bool_t MFCT1Supercuts::Process() +{ + // apply dynamical Scuts + // optimized for mkn 421 2001 data + + Double_t dNOMLOGSIZE = 4.1; + Double_t dNOMCOSZA = 1.0; + + + Double_t newdist = fHilSrc->GetDist() * fMm2Deg; + + Double_t xbar = fHil->GetMeanX(); + Double_t ybar = fHil->GetMeanY(); + Double_t dist = sqrt(xbar*xbar + ybar*ybar) * fMm2Deg; + Double_t dd2 = dist * dist; + + Double_t dmls = log(fHil->GetSize()) - dNOMLOGSIZE; + Double_t dmls2 = dmls * dmls; + + Double_t dmcza = cos(fMcEvt->GetTelescopeTheta()) - dNOMCOSZA; + + Double_t length = fHil->GetLength() *fMm2Deg; + Double_t width = fHil->GetWidth() *fMm2Deg; + + *fLog << "MFCT1Supercuts::Process; dmls, dmcza, dmls2, dd2 = " + << dmls << ", " << dmcza << ", " << dmls2 << ", " + << dd2 << endl; + + *fLog << "MFCT1Supercuts::Process; newdist, dist, length, width = " + << newdist << ", " << dist << ", " << length << ", " + << width << endl; + + if ( newdist < 1.05 + && newdist < CtsMCut (distup, dmls, dmcza, dmls2, dd2) + && newdist > CtsMCut (distlow, dmls, dmcza, dmls2, dd2) + && dist < 1.05 + && length < CtsMCut (lengthup, dmls, dmcza, dmls2, dd2) + && length > CtsMCut (lengthlow, dmls, dmcza, dmls2, dd2) + && width < CtsMCut (widthup, dmls, dmcza, dmls2, dd2) + && width > CtsMCut (widthlow, dmls, dmcza, dmls2, dd2) + //&& asym < CtsMCut (asymup, dmls, dmcza, dmls2, dd2) + //&& asym > CtsMCut (asymlow, dmls, dmcza, dmls2, dd2) + && dist < CtsMCut (distup, dmls, dmcza, dmls2, dd2) + && dist > CtsMCut (distlow, dmls, dmcza, dmls2, dd2) ) + fResult = kTRUE; + else + fResult = kFALSE; + + return kTRUE; +} +// --------------------------------------------------------------------------- + Index: trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.h =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.h (revision 1905) +++ trunk/MagicSoft/Mars/mfilter/MFCT1Supercuts.h (revision 1905) @@ -0,0 +1,73 @@ +#ifndef MARS_MFCT1Supercuts +#define MARS_MFCT1Supercuts + +///////////////////////////////////////////////////////////////////////////// +// // +// MFCT1Supercuts // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MFilter +#include "MFilter.h" +#endif + +class MParList; +class MHillas; +class MHillasSrc; +class MMcEvt; +class MCerPhotEvt; +class MGeomCam; + + +class MFCT1Supercuts : public MFilter +{ +private: + MHillas *fHil; + MHillasSrc *fHilSrc; + MMcEvt *fMcEvt; + MCerPhotEvt *fEvt; + MGeomCam *fCam; + + TString fHilName; + TString fHilSrcName; + + Double_t fMm2Deg; + + Bool_t fResult; //! + + //--------------------------------- + // cut parameters + + double lengthup[8]; + double widthup[8]; + double distup[8]; + double lengthlow[8]; + double widthlow[8]; + double distlow[8]; + double asymup[8]; + double asymlow[8]; + double alphaup[8]; + //--------------------------------- + + +public: + MFCT1Supercuts(const char *hilname="MHillas", + const char *hilsrcname="MHillasSrc", + const char *name=NULL, const char *title=NULL); + + Double_t CtsMCut(Double_t *a, Double_t ls, Double_t ct, + Double_t ls2, Double_t dd2); + + Bool_t IsExpressionTrue() const { return fResult; } + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + + ClassDef(MFCT1Supercuts, 1) // A class to evaluate the Supercuts +}; + +#endif + + + + + Index: trunk/MagicSoft/Mars/mfilter/Makefile =================================================================== --- trunk/MagicSoft/Mars/mfilter/Makefile (revision 1904) +++ trunk/MagicSoft/Mars/mfilter/Makefile (revision 1905) @@ -20,5 +20,6 @@ # @endcode -INCLUDES = -I. -I../mbase -I../mmc -I../manalysis -I../mdata -I../mfileio +INCLUDES = -I. -I../mbase -I../mraw -I../mmc -I../mdata -I../manalysis \ + -I../mfileio -I../mgeom # @code @@ -39,5 +40,9 @@ MFDataMember.cc \ MFParticleId.cc \ - MFAlpha.cc + MFAlpha.cc \ + MFCT1SelBasic.cc \ + MFCT1SelStandard.cc \ + MFCT1SelFinal.cc \ + MFCT1Supercuts.cc # MFDataChain.cc \