Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 1398) +++ /trunk/MagicSoft/Mars/Changelog (revision 1399) @@ -1,3 +1,8 @@ -*-*- END -*-*- + 2002/07/11: Abelardo Moralejo + + * manalysis/MCerPhotCalc2.[h,cc] + - new class to compute nphe weighting FADC slices. + 2002/07/10: Oscar Blanch Index: /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.cc =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.cc (revision 1399) +++ /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.cc (revision 1399) @@ -0,0 +1,218 @@ +/* ======================================================================== *\ +! +! * +! * 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): Abelardo Moralejo 7/2002 (moralejo@pd.infn.it) +! +! Copyright: MAGIC Software Development, 2002 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// // +// MCerPhotCalc2 // +// // +// This is a task which calculates the number of photons from the FADC // +// time slices. It weights the each slice according to the fraction // +// of signal which, in average, falls within it. This average has been // +// determined over a run (about 1000 triggers, 0 deg z.a. gammas). // +// // +// Input Containers: // +// MRawEvtData, MPedesdtalCam // +// // +// Output Containers: // +// MCerPhotEvt // +// // +////////////////////////////////////////////////////////////////////////////// + +#include "MCerPhotCalc2.h" + +#include "MParList.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MMcRunHeader.hxx" + +#include "MRawRunHeader.h" +#include "MRawEvtData.h" // MRawEvtData::GetNumPixels +#include "MCerPhotEvt.h" +#include "MPedestalPix.h" +#include "MPedestalCam.h" +#include "MRawEvtPixelIter.h" + +ClassImp(MCerPhotCalc2); + +const Float_t fWeight[15] ={0, 0.0809835, 0.289593, 0.366926, 0.211665, 0.0508328, 0., 0., 0., 0., 0., 0., 0., 0., 0.}; + +Float_t SumQuadWeights, SumWeights; + + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MCerPhotCalc2::MCerPhotCalc2(const char *name, const char *title) +{ + fName = name ? name : "MCerPhotCalc2"; + fTitle = title ? title : "Task to calculate Cerenkov photons from raw data"; + + AddToBranchList("MRawEvtData.fHiGainPixId"); + AddToBranchList("MRawEvtData.fLoGainPixId"); + AddToBranchList("MRawEvtData.fHiGainFadcSamples"); + AddToBranchList("MRawEvtData.fLoGainFadcSamples"); + + SumQuadWeights = SumWeights = 0.; + for (Int_t i = 0; i < 15; i++) + { + SumQuadWeights += pow(fWeight[i],2.); + SumWeights += fWeight[i]; + } + SumQuadWeights = sqrt(SumQuadWeights); + +} + +// -------------------------------------------------------------------------- +// +// The PreProcess searches for the following input containers: +// - MRawRunHeader +// - MRawEvtData +// - MPedestalCam +// +// The following output containers are also searched and created if +// they were not found: +// - MCerPhotEvt +// +Bool_t MCerPhotCalc2::PreProcess(MParList *pList) +{ + fRunHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); + if (!fRunHeader) + { + *fLog << dbginf << "MRawRunHeader not found... aborting." << endl; + return kFALSE; + } + + fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData"); + if (!fRawEvt) + { + *fLog << dbginf << "MRawEvtData not found... aborting." << endl; + return kFALSE; + } + + fPedestals = (MPedestalCam*)pList->FindCreateObj("MPedestalCam"); + if (!fPedestals) + { + *fLog << dbginf << "MPedestalCam not found... aborting." << endl; + return kFALSE; + } + + fCerPhotEvt = (MCerPhotEvt*)pList->FindCreateObj("MCerPhotEvt"); + if (!fCerPhotEvt) + return kFALSE; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Check for the run type and camera version. +// If the file is a MC file and the used camera version is <= 40 +// we enable a fix for truncated pedestal means in this version. +// +Bool_t MCerPhotCalc2::ReInit(MParList *pList) +{ + const MRawRunHeader *runheader = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); + if (!runheader) + { + *fLog << warn << dbginf << "Warning - cannot check file type, MRawRunHeader not found." << endl; + return kTRUE; + } + + if (runheader->GetRunType() != kRTMonteCarlo) + return kTRUE; + + MMcRunHeader *mcrunheader = (MMcRunHeader*)pList->FindObject("MMcRunHeader"); + if (!mcrunheader) + { + *fLog << warn << dbginf << "Warning - cannot check for camera version, MC run header not found." << endl; + return kTRUE; + } + + if (mcrunheader->GetCamVersion() <= 40) + fEnableFix = kTRUE; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Calculate the integral of the FADC time slices and store them as a new +// pixel in the MCerPhotEvt container. +// +Bool_t MCerPhotCalc2::Process() +{ + fCerPhotEvt->InitSize(fRawEvt->GetNumPixels()); + + MRawEvtPixelIter pixel(fRawEvt); + + Float_t ADCBins[15]; + + while (pixel.Next()) + { + const UInt_t pixid = pixel.GetPixelId(); + const MPedestalPix &ped = (*fPedestals)[pixid]; + + // + // sanity check (old MC files sometimes have pixids>577) + // + if (!fPedestals->CheckBounds(pixid)) + { + *fLog << inf << "Pixel ID larger than camera... skipping event." << endl; + return kCONTINUE; + } + + // Mean pedestal: + const Double_t mean = fEnableFix ? ped.GetMean()-0.5 : ped.GetMean(); + + Byte_t *ptr = pixel.GetHiGainSamples(); + + Float_t nphot = 0.; + for(Int_t i = 0; i<15;i++) + { + ADCBins[i] = (Float_t) *(ptr+i); + nphot += ADCBins[i] * fWeight[i]; + } + + // + // We check that the pixel is not empty, if it is empty + // we won't substract the pedestal. Empty means that it has + // 0 signal in all the slices. + // + if (nphot!=0) + nphot -= mean*SumWeights; + + Float_t nphoterr = ped.GetSigma()* SumQuadWeights; + + fCerPhotEvt->AddPixel(pixid, nphot, nphoterr); + + // FIXME! Handling of Lo Gains is missing! + } + + fCerPhotEvt->SetReadyToSave(); + + return kTRUE; +} Index: /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.h =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.h (revision 1399) +++ /trunk/MagicSoft/Mars/manalysis/MCerPhotCalc2.h (revision 1399) @@ -0,0 +1,43 @@ +#ifndef MARS_MCerPhotCalc2 +#define MARS_MCerPhotCalc2 + +///////////////////////////////////////////////////////////////////////////// +// // +// MCerPhotCalc2 // +// // +// Integrates certain time slices of one pixel and substracts the pedestal // +// (offset) value // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MRawEvtData; +class MPedestalCam; +class MCerPhotEvt; +class MRawRunHeader; + +class MCerPhotCalc2 : public MTask +{ + MPedestalCam *fPedestals; // Pedestals of all pixels in the camera + MRawEvtData *fRawEvt; // raw event data (time slices) + MCerPhotEvt *fCerPhotEvt; // Cerenkov Photon Event used for calculation + MRawRunHeader *fRunHeader; // RunHeader information + + Bool_t fEnableFix; // fix for a bug in files from older camera versions (<=40) + + +public: + MCerPhotCalc2(const char *name=NULL, const char *title=NULL); + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + + Bool_t ReInit(MParList *pList); + + ClassDef(MCerPhotCalc2, 0) // Task to calculate cerenkov photons from raw data +}; + +#endif