Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 1984) +++ trunk/MagicSoft/Mars/Changelog (revision 1985) @@ -1,3 +1,12 @@ -*-*- END OF LINE -*-*- + + 2003/04/23: Antonio Stamerra + + * mfilter/MFEnergySlope.[h,cc] + - added new filter class to select events according to a given + energy slope + + * mfilter/Makefile, mfilter/FilterLinkDef.h + - modified accordingly (added MFEnergySlope) 2003/04/23: Abelardo Moralejo Index: trunk/MagicSoft/Mars/mfilter/MFEnergySlope.cc~ =================================================================== --- trunk/MagicSoft/Mars/mfilter/MFEnergySlope.cc~ (revision 1984) +++ (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): Antonio Stamerra 02/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! Filter to select events with a energy slope. 02/2003 -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MFEnergySlope // -// // -// A filter to select MC events (generated with a energy slope MCSlope) // -// with a different energy slope NewSlope set by the user. // -// // -// The new slope is set by the user with SetSlope(). // -// Only negative slopes are admitted; positive ones are internally // -// converted. // -// Events are selected following the new energy power slope, and the // -// sample is normalized to the number of events at: // -// 1. McMaxEnergy, if abs(NewSlope) < abs(McSlope); // -// 2. McMinEnergy, if abs(NewSlope) > abs(McSlope); // -// Mc{Max,Min}Energy are set with SetMcMinEnergy() and SetMcMaxEnergy(); // -// with GeV values. // -// Default values are the min/max energy of the MC sample. // -// // -// With this filter the statistics of the MC sample is reduced. // -// camera ver.0.6 and reflector ver.0.6 are required to fetch // -// the correct MC information // -// // -///////////////////////////////////////////////////////////////////////////// -#include "MFEnergySlope.h" - -//#include -#include - -#include "MParList.h" - -#include "MLog.h" -#include "MLogManip.h" - - -#include "MMcEvt.hxx" -#include "MMcCorsikaRunHeader.h" - -ClassImp(MFEnergySlope); - -// -------------------------------------------------------------------------- -// -// Constructor -// -MFEnergySlope::MFEnergySlope(const char *name, const char *title): - fNumSelectedEvts(0), fNewSlope(-1), fMcMinEnergy(-1.), fMcMaxEnergy(-1.) -{ - // fContName = cname; - fName = name ? name : "MFEnergySlope"; - fTitle = title ? title : "Filter to select energy with given slope"; -} - -// -------------------------------------------------------------------------- -// -// Destructor -// -MFEnergySlope::~MFEnergySlope() -{ -} - -// -------------------------------------------------------------------------- -// -// Preprocess -// -// MC slope and min/max energy are read -// Normalization factor is computed -// -Bool_t MFEnergySlope::PreProcess(MParList *pList) -{ - - MMcCorsikaRunHeader *runheader = (MMcCorsikaRunHeader*)pList->FindObject("MMcCorsikaRunHeader"); - - if (!runheader) - { - *fLog << err << dbginf << fName << " [MMcCorsikaRunHeader] not found... aborting." << endl; - return kFALSE; - } - // - // Read info from the MC sample (it must be generated with - // reflector ver.0.6 and camera ver. 0.6) - // - fMcSlope = runheader->GetSlopeSpec(); - if (fMcMinEnergy<0) - fMcMinEnergy = runheader->GetELowLim(); - if (fMcMinEnergy<0) - fMcMaxEnergy = runheader->GetEUppLim(); - - *fLog << "MFEnergySlope::PreProcess; fetched MC info:" << endl << - " E Slope: "<< fMcSlope << endl << - " E Min: "<< fMcMinEnergy << endl << - " E Max: "<< fMcMaxEnergy << endl << - " New E Slope: "<< fNewSlope << endl; - - // Slope is used with positive values in the code - if (fNewSlope < 0) - fNewSlope *= -1; - if (fMcSlope < 0) - fMcSlope *= -1; - - - // Set normalization on energy - fN0 = pow(fNewSlope>fMcSlope?fMcMinEnergy:fMcMaxEnergy,fNewSlope-fMcSlope); - - *fLog << "Normalization factor:" <FindObject("MMcEvt"); - if (!fEvt) - { - *fLog << err << dbginf << fName << " [MMcEvt] not found... aborting." << endl; - return kFALSE; - } - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Select events randomly according to the MC ("old") and required ("new") -// energy slope. -// Old E slope is fMcSlope -// New E slope is set by the user (fval; fNewSlope) -// If old and new energy slope are the same skip the selection. -// The MC energy slope and lower and upper limits are taken from the -// run header (requires reflector ver.>0.6 and camera ver.>0.6) -// -Bool_t MFEnergySlope::Process() -{ - fResult = kTRUE; - - // Energy slopes are the same: skip it - if (fNewSlope == fMcSlope) - return kTRUE; - - // The value of the normalized spectrum is compared with a - // random value in [0,1]; - // if the latter is higher the event is accepted - const Float_t energy = fEvt->GetEnergy(); - - // - // If energy bounds different from MC ones are selected, then - // events outside these bounds are rejected, as their slope has - // not been changed. - // - if (energy > fMcMaxEnergy || energy < fMcMinEnergy) - { - fResult = kFALSE; - return kTRUE; - } - const Float_t Nexp = fN0 * pow(energy,fMcSlope-fNewSlope); - const Float_t Nrnd = gRandom->Uniform(); - - fResult = Nexp > Nrnd; - - if (fResult) - { - fNumSelectedEvts++; - return kTRUE; - } - - return kTRUE; -} - - -// -------------------------------------------------------------------------- -// -// Postprocess all filters. -// -Bool_t MFEnergySlope::PostProcess() -{ - //--------------------------------- - - return kTRUE; -} - -void MFEnergySlope::StreamPrimitive(ofstream &out) const -{ - /* - if (fEvt) - fEvt->SavePrimitive(out); - - out << " MFEnergySlope " << GetUniqueName() << "("; - - if (fEvt) - out << "&" << fEvt->GetUniqueName(); - else - out << "\"" << fContName << "\""; - - out << ", New Energy Slope= " << fValue << " Normalization factor: " - << fN0 <<" );" << endl; - */ -} - - - - - - - - -