Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 3125) +++ trunk/MagicSoft/Mars/Changelog (revision 3126) @@ -4,4 +4,14 @@ -*-*- END OF LINE -*-*- + 2004/02/13: Thomas Bretz + + * mbadpixels/MBadPixelsTreat.[h.cc]: + - added + + * mbadpixels/Makefile, mbadpixels/BadPixelsLinkDef.h: + - added MBadPixelsTreat + + + 2004/02/12: Markus Gaug @@ -36,4 +46,5 @@ * mtools/MFFT.[h,cc] - PowerSpectrumDensity of TArrayI implemented + Index: trunk/MagicSoft/Mars/NEWS =================================================================== --- trunk/MagicSoft/Mars/NEWS (revision 3125) +++ trunk/MagicSoft/Mars/NEWS (revision 3126) @@ -2,6 +2,9 @@ *** Version - - added filter against cosmics: MFCosmics: - To be used together with MContinue + - added new classes for the bad-pixels treatment (MBadPixels*) + which are more powerfull than the old ones (MBlindPixel*) + and will replace them. + + - added filter against cosmics: MFCosmics - added new class MArrivalTimeCalc2: Index: trunk/MagicSoft/Mars/mbadpixels/BadPixelsLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mbadpixels/BadPixelsLinkDef.h (revision 3125) +++ trunk/MagicSoft/Mars/mbadpixels/BadPixelsLinkDef.h (revision 3126) @@ -9,4 +9,5 @@ #pragma link C++ class MBadPixelsCalc+; +#pragma link C++ class MBadPixelsTreat+; #pragma link C++ class MBadPixelsMerge+; Index: trunk/MagicSoft/Mars/mbadpixels/MBadPixelsTreat.cc =================================================================== --- trunk/MagicSoft/Mars/mbadpixels/MBadPixelsTreat.cc (revision 3126) +++ trunk/MagicSoft/Mars/mbadpixels/MBadPixelsTreat.cc (revision 3126) @@ -0,0 +1,320 @@ +/* ======================================================================== *\ +! +! * +! * 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): Oscar Blanch 12/2001 +! Author(s): Thomas Bretz 08/2002 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MBadPixelsTreat +// +// You can use MBadPixelsTreat::SetUseInterpolation to replaced the +// bad pixels by the average of its neighbors instead of unmapping +// them. If you want to include the central pixel use +// MBadPixelsTreat::SetUseCentralPixel. The bad pixels are taken from +// an existing MBadPixelsCam. +// +// Input Containers: +// MCerPhotEvt +// MBadPixelsCam +// +// Output Containers: +// MCerPhotEvt +// +///////////////////////////////////////////////////////////////////////////// +#include "MBadPixelsTreat.h" + +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MGeomPix.h" +#include "MGeomCam.h" + +#include "MCerPhotPix.h" +#include "MCerPhotEvt.h" + +#include "MPedPhotPix.h" +#include "MPedPhotCam.h" + +#include "MBadPixelsPix.h" +#include "MBadPixelsCam.h" + +ClassImp(MBadPixelsTreat); + +using namespace std; + +static const TString gsDefName = "MBadPixelsTreat"; +static const TString gsDefTitle = "Task to treat bad pixels (interpolation, unmapping)"; + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MBadPixelsTreat::MBadPixelsTreat(const char *name, const char *title) + : fFlags(0) +{ + fName = name ? name : gsDefName.Data(); + fTitle = title ? title : gsDefTitle.Data(); +} + +// -------------------------------------------------------------------------- +// +// - Try to find or create MBlindPixels in parameter list. +// - get the MCerPhotEvt from the parlist (abort if missing) +// - if no pixels are given by the user try to determin the starfield +// from the monte carlo run header. +// +Int_t MBadPixelsTreat::PreProcess (MParList *pList) +{ + fBadPixels = (MBadPixelsCam*)pList->FindObject(AddSerialNumber("MBadPixelsCam")); + if (!fBadPixels) + { + *fLog << err << "MBadPixelsCam not found... aborting." << endl; + return kFALSE; + } + + fPedPhot = (MPedPhotCam*)pList->FindObject(AddSerialNumber("MPedPhotCam")); + if (!fPedPhot) + { + *fLog << err << "MPedPhotCam not found... aborting." << endl; + return kFALSE; + } + + fEvt = (MCerPhotEvt*)pList->FindObject(AddSerialNumber("MCerPhotEvt")); + if (!fEvt) + { + *fLog << err << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + fGeomCam = (MGeomCam*)pList->FindObject(AddSerialNumber("MGeomCam")); + if (!fGeomCam && TESTBIT(fFlags, kUseInterpolation)) + { + *fLog << err << "No camera geometry available... can't use interpolation." << endl; + return kFALSE; + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Replaces each pixel (signal, signal error, pedestal, pedestal rms) +// by the average of its surrounding pixels. +// If TESTBIT(fFlags, kUseCentralPixel) is set the central pixel is also +// included. +// +// NT: Informations about the interpolated pedestals are added. +// When the option Interpolated is called, the blind pixel with the new +// values of signal and fluttuation is included in the calculation of +// the Image Parameters. +// +void MBadPixelsTreat::Interpolate() const +{ + const UShort_t entries = fGeomCam->GetNumPixels(); + + // + // Create arrays + // + Double_t *nphot = new Double_t[entries]; + Double_t *perr = new Double_t[entries]; + Double_t *ped = new Double_t[entries]; + Double_t *pedrms = new Double_t[entries]; + + // + // Loop over all pixels + // + for (UShort_t i=0; iGetPixById(i); + if (!pix) + pix = fEvt->AddPixel(i, 0, 0); + + // + // Get the corresponding geometry and pedestal + // + const MGeomPix &gpix = (*fGeomCam)[i]; + const MPedPhotPix &ppix = (*fPedPhot)[i]; + + // Do Not-Use-Central-Pixel + const Bool_t nucp = !TESTBIT(fFlags, kUseCentralPixel); + + Int_t num = nucp ? 0 : 1; + + nphot[i] = nucp ? 0 : pix->GetNumPhotons(); + ped[i] = nucp ? 0 : ppix.GetMean(); + perr[i] = nucp ? 0 : Sqr(pix->GetErrorPhot()); + pedrms[i] = nucp ? 0 : Sqr(ppix.GetRms()); + + // + // The values are rescaled to the small pixels area for the right comparison + // + const Double_t ratio = fGeomCam->GetPixRatio(i); + + nphot[i] *= ratio; + ped[i] *= ratio; + perr[i] *= Sqr(ratio); + pedrms[i] *= Sqr(pedrms[i]); + + // + // Loop over all its neighbors + // + const Int_t n = gpix.GetNumNeighbors(); + for (int j=0; jGetPixById(nidx); + if (!evtpix) + continue; + + // + // Get the geometry for the neighbor + // + const Double_t nratio = fGeomCam->GetPixRatio(nidx); + MPedPhotPix &nppix = (*fPedPhot)[nidx]; + + // + //The error is calculated as the quadratic sum of the errors + // + ped[i] += nratio*nppix.GetMean(); + nphot[i] += nratio*evtpix->GetNumPhotons(); + perr[i] += Sqr(nratio*evtpix->GetErrorPhot()); + pedrms[i] += Sqr(nratio*nppix.GetRms()); + + num++; + } + + // + // Now the mean is calculated and the values rescaled back to the pixel area + // + nphot[i] /= num*ratio; + ped[i] /= num*ratio; + perr[i] = TMath::Sqrt(perr[i]/num)*ratio; + pedrms[i] = TMath::Sqrt(pedrms[i]/num)*ratio; + + } + + // + // Now the new pixel values are calculated and can be replaced in + // the corresponding containers + // + for (UShort_t i=0; iGetPixById(i)->Set(nphot[i], perr[i]); + (*fPedPhot)[i].Set(ped[i], pedrms[i]); + } + + // + // Delete the intermediat arrays + // + delete nphot; + delete perr; + delete ped; + delete pedrms; +} + +// -------------------------------------------------------------------------- +// +// Removes all blind pixels from the analysis by setting their state +// to unused. +// +void MBadPixelsTreat::Unmap() const +{ + const UShort_t entries = fEvt->GetNumPixels(); + + // + // remove the pixels in fPixelsIdx if they are set to be used, + // (set them to 'unused' state) + // + for (UShort_t i=0; i