Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 4445) +++ /trunk/MagicSoft/Mars/Changelog (revision 4446) @@ -19,4 +19,15 @@ -*-*- END OF LINE -*-*- + 2004/08/04: Thomas Bretz + + * manalysis/MBlindPixels.[h,cc], manalysis/MBlindPixelCalc.[h,cc], + manalysis/MBlindPixelCalc2.[h,cc]: + - as discussed some weeks ago the classes have been removed + + * manalysis/Makefile, manalysis/AnalysisLinkDef.h: + - removed the MBlind* entries + + + 2004/08/03: Hendrik Bartko Index: /trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h =================================================================== --- /trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 4445) +++ /trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 4446) @@ -13,8 +13,4 @@ #pragma link C++ class MCameraData+; - -#pragma link C++ class MBlindPixels+; -#pragma link C++ class MBlindPixelCalc+; -#pragma link C++ class MBlindPixelsCalc2+; #pragma link C++ class MHPedestalPix+; Index: unk/MagicSoft/Mars/manalysis/MBlindPixelCalc.cc =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MBlindPixelCalc.cc (revision 4445) +++ (revision ) @@ -1,441 +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): Oscar Blanch 12/2001 -! Author(s): Thomas Bretz 08/2002 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MBlindPixelCalc -// -// This is the specific image cleaning for a list of pixels. This task -// sets the pixels listed in fPixelsIdx to unused so they should not be -// used for analysis (eg calculation of hillas parameters). -// -// You can use MBlindPixelCalc::SetUseInterpolation to replaced the -// blind pixels by the average of its neighbors instead of unmapping -// them. If you want to include the central pixel use -// MBlindPixelCalc::SetUseCentralPixel. -// -// You have three options: -// 1) Call SetUseBlindPixels(): -// This will take an existing MBlindPixels container filled from -// elsewhere (eg. MCT1ReadPreProc) and use this pixels as blind -// pixels. -// 2) Call SetPixels(): -// This will setup an array with pixel numbers. These pixels are used -// as blind pixels. -// 3) Neither 1) nor 2) -// This options tries to identify the starfield from the -// MMcRunHeader container and tries to identifies it. If it is known -// (eg. Crab) the fixed build in pixel numbers are used as blind -// pixels. -// -// If neither an array of pixels is given (or its size is 0) and -// MMcRunHeader couldn't be found the task removes itself from the -// tasklist. -// -// Implemented star fields: -// - Crab: 400, 401, 402, 437, 438, 439 -// -// Input Containers: -// MCerPhotEvt -// [MBlindPixels] -// -// Output Containers: -// MCerPhotEvt -// MBlindPixels -// -///////////////////////////////////////////////////////////////////////////// -#include "MBlindPixelCalc.h" - -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" - -#include "MGeomPix.h" -#include "MGeomCam.h" -#include "MCerPhotPix.h" -#include "MCerPhotEvt.h" -#include "MPedestalPix.h" -#include "MPedestalCam.h" -#include "MBlindPixels.h" - -#include "MMcRunHeader.hxx" - -ClassImp(MBlindPixelCalc); - -using namespace std; - -static const TString gsDefName = "MBlindPixelCalc"; -static const TString gsDefTitle = "Find hot spots (star, broken pixels, etc)"; - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MBlindPixelCalc::MBlindPixelCalc(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 MBlindPixelCalc::PreProcess (MParList *pList) -{ - if (TESTBIT(fFlags, kUseBlindPixels)) - fPixels = (MBlindPixels*)pList->FindObject(AddSerialNumber("MBlindPixels")); - else - fPixels = (MBlindPixels*)pList->FindCreateObj(AddSerialNumber("MBlindPixels")); - if (!fPixels) - return kFALSE; - - fEvt = (MCerPhotEvt*)pList->FindObject(AddSerialNumber("MCerPhotEvt")); - if (!fEvt) - { - *fLog << err << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - fPed = (MPedestalCam*)pList->FindObject(AddSerialNumber("MPedestalCam")); - if (!fEvt) - { - *fLog << err << dbginf << "MPedestalCam not found... aborting." << endl; - return kFALSE; - } - - fGeomCam = (MGeomCam*)pList->FindObject(AddSerialNumber("MGeomCam")); - if (!fGeomCam) - *fLog << warn << dbginf << "No camera geometry available... can't use interpolation." << endl; - - if (TESTBIT(fFlags, kUseBlindPixels)) - return kTRUE; - - const UShort_t size = fPixelsIdx.GetSize(); - - if (size == 0) - { - if (!pList->FindObject("MMcRunHeader")) - { - *fLog << warn << "Warning - Neither blind pixels are given nor a MMcRunHeader was found... removing MBlindPixelCalc from list." << endl; - return kSKIP; - } - return kTRUE; - } - - // Set as blind pixels the global blind pixels, which are given - // through the macros - - UShort_t numids = fPixelsIdx.GetSize(); - - for(Int_t i = 0; iSetPixelBlind(fPixelsIdx[i]); - - 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 MBlindPixelCalc::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; iIsBlind(i)) - continue; - - // - // Get a pointer to this pixel. If it is not yet existing - // create a new entry for this pixel in MCerPhotEvt - // - const MCerPhotPix *pix = fEvt->GetPixById(i); - if (!pix) - pix = fEvt->AddPixel(i, 0, 0); - - // - // Get the corresponding geometry and pedestal - // - const MGeomPix &gpix = (*fGeomCam)[i]; - const MPedestalPix &ppix = (*fPed)[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.GetPedestal(); - perr[i] = nucp ? 0 : Sqr(pix->GetErrorPhot()); - pedrms[i] = nucp ? 0 : Sqr(ppix.GetPedestalRms()); - - // - // 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; jIsBlind(nidx)) - continue; - - // - // Check whether the neighbor has a signal stored - // - const MCerPhotPix *evtpix = fEvt->GetPixById(nidx); - if (!evtpix) - continue; - - // - // Get the geometry for the neighbor - // - const Double_t nratio = fGeomCam->GetPixRatio(nidx); - MPedestalPix &nppix = (*fPed)[nidx]; - - // - //The error is calculated as the quadratic sum of the errors - // - ped[i] += nratio*nppix.GetPedestal(); - nphot[i] += nratio*evtpix->GetNumPhotons(); - perr[i] += Sqr(nratio*evtpix->GetErrorPhot()); - pedrms[i] += Sqr(nratio*nppix.GetPedestalRms()); - - 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; iIsBlind(i)) - continue; - - // - // It must exist, we have created it in the loop before. - // - fEvt->GetPixById(i)->Set(nphot[i], perr[i]); - (*fPed)[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 MBlindPixelCalc::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; iIsBlind(pix.GetPixId())) - pix.SetPixelUnused(); - } -} - -// -------------------------------------------------------------------------- -// -// Treat the blind pixels -// -Int_t MBlindPixelCalc::Process() -{ - if (TESTBIT(fFlags, kUseInterpolation) && fGeomCam) - Interpolate(); - else - Unmap(); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// - Check whether pixels to disable are available. If pixels are -// given by the user nothing more is done. -// - Otherwise try to determin the blind pixels from the starfield -// given in MMcRunHeader. -// -Bool_t MBlindPixelCalc::ReInit(MParList *pList) -{ - if (TESTBIT(fFlags, kUseBlindPixels)) - return kTRUE; - - // - // If pixels are given by the user, we are already done - // - if (fPixelsIdx.GetSize() > 0) - return kTRUE; - - // - // Delete the old array holding the blind pixels for the last file - // - fPixels->Clear(); - - if (!fGeomCam->InheritsFrom("MGeomCamMagic")) - { - *fLog << warn << "MBlindPixelCalc::ReInit: Warning - Starfield only implemented for Magic standard Camera... no action." << endl; - return kTRUE; - } - - // - // Set as blind some particular pixels because of a particular - // Star Field of View. - // - MMcRunHeader *mcrun = (MMcRunHeader*)pList->FindObject("MMcRunHeader"); - if (!mcrun) - { - *fLog << warn << "MBlindPixelCalc::ReInit: Warning - No run header available... no action." << endl; - return kTRUE; - } - - Int_t rah, ram, ras; - Int_t ded, dem, des; - mcrun->GetStarFieldRa(&rah, &ram, &ras); - mcrun->GetStarFieldDec(&ded, &dem, &des); - - if (rah!=5 || ram!=34 || ras!=32 || ded!=22 || dem!=0 || des!=55) - { - *fLog << warn << "Warning - Starfield unknown..." << endl; - return kTRUE; - } - - // - // Case for Crab Nebula FOV - // - fPixels->SetPixelBlind(400); - fPixels->SetPixelBlind(401); - fPixels->SetPixelBlind(402); - fPixels->SetPixelBlind(437); - fPixels->SetPixelBlind(438); - fPixels->SetPixelBlind(439); - - *fLog << inf; - *fLog << "FOV is centered at CRAB NEBULA: Setting 6 blind pixels" << endl; - *fLog << "to avoid bias values of analysis due to CRAB NEBULA:" << endl; - *fLog << " Pixels: 400, 401, 402, 437, 438, 439" << endl; - - return kTRUE; -} - -void MBlindPixelCalc::StreamPrimitive(ofstream &out) const -{ - out << " MBlindPixelCalc " << GetUniqueName(); - if (fName!=gsDefName || fTitle!=gsDefTitle) - { - out << "(\"" << fName << "\""; - if (fTitle!=gsDefTitle) - out << ", \"" << fTitle << "\""; - out <<")"; - } - out << ";" << endl; - - if (TESTBIT(fFlags, kUseInterpolation)) - out << " " << GetUniqueName() << ".SetUseInterpolation();" << endl; - if (TESTBIT(fFlags, kUseCentralPixel)) - out << " " << GetUniqueName() << ".SetUseCentralPixel();" << endl; - - if (fPixelsIdx.GetSize()==0) - return; - - out << " {" << endl; - out << " TArrayS dummy;" << endl; - for (int i=0; i -#endif - -class MGeomCam; -class MPedestal; -class MCerPhotEvt; -class MBlindPixels; -class MPedestalCam; - -class MBlindPixelCalc : public MTask -{ -private: - MCerPhotEvt *fEvt; //! - MBlindPixels *fPixels; //! - MGeomCam *fGeomCam; //! - MPedestalCam *fPed; //! - - TArrayS fPixelsIdx; // Pixel Indices for blind pixels, which are entered by the user. - - Byte_t fFlags; // flag for the method which is used - - enum - { - kUseInterpolation = 1, - kUseCentralPixel = 2, - kUseBlindPixels = 3 - }; - - static Double_t Sqr(Double_t x) { return x*x; } - - void Interpolate() const; - void Unmap() const; - void StreamPrimitive(ofstream &out) const; - - Int_t PreProcess(MParList *pList); - Int_t Process(); - -public: - MBlindPixelCalc(const char *name=NULL, const char *title=NULL); - - void SetUseInterpolation(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseInterpolation) : CLRBIT(fFlags, kUseInterpolation); - } - void SetUseCentralPixel(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseCentralPixel) : CLRBIT(fFlags, kUseCentralPixel); - } - void SetUseBlindPixels(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseBlindPixels) : CLRBIT(fFlags, kUseBlindPixels); - } - - void SetPixelIndices(const TArrayS &pix) { fPixelsIdx = pix; } - virtual Bool_t ReInit(MParList *pList); - - ClassDef(MBlindPixelCalc, 1) // task to deal with hot spots (star, broken pixels, etc) -}; - -#endif - Index: unk/MagicSoft/Mars/manalysis/MBlindPixels.cc =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MBlindPixels.cc (revision 4445) +++ (revision ) @@ -1,64 +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): Oscar Blanch 12/2001 -! -! Copyright: MAGIC Software Development, 2000-2001 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MBlindPixel // -// // -// If you want to exclude pixels from the analysis (eg. hillas parameter // -// calculation) you can use this container to switch off the pixels by Id // -// // -///////////////////////////////////////////////////////////////////////////// -#include "MBlindPixels.h" - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" - -ClassImp(MBlindPixels); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MBlindPixels::MBlindPixels(const char *name, const char *title) : fPixels(577) - -{ - fName = name ? name : "MBlindPixels"; - fTitle = title ? title : "Container which have a list of blind pixels"; -} - -// ------------------------------------------------------------------------- -// -// If you don't want to use pixel with Id id call -// SetPixelBlind(id, kTRUE) -// otherwise -// SetPixelBlind(id, kFALSE) -// -void MBlindPixels::SetPixelBlind(UShort_t id, Bool_t val) -{ - fPixels.AddAt(val, id); -} - Index: unk/MagicSoft/Mars/manalysis/MBlindPixels.h =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MBlindPixels.h (revision 4445) +++ (revision ) @@ -1,35 +1,0 @@ -#ifndef MARS_MBlindPixels -#define MARS_MBlindPixels - -#ifndef MARS_MParContainer -#include "MParContainer.h" -#endif - -#ifndef ROOT_TArrayC -#include -#endif - -class MBlindPixels : public MParContainer -{ -private: - TArrayC fPixels; // Boolean value for every pixel 0=on, 1=blind - -public: - MBlindPixels(const char *name=NULL, const char *title=NULL); - - void SetPixelBlind(UShort_t id, Bool_t val=kTRUE); - void AdoptPixels(Char_t pixels[577]) { fPixels.Adopt(577, pixels); } - - void AdoptPixels(const TArrayC &pixels) { fPixels = pixels; } - - void Clear(Option_t *o="") { fPixels.Reset(); } - - Bool_t IsBlind(UShort_t id) const { return fPixels.GetSize() && ((TArrayC)fPixels)[id]; } - - void InitSize(UShort_t i) { fPixels.Set(i); } - - ClassDef(MBlindPixels, 1) // container to store blind pixels -}; - -#endif - Index: unk/MagicSoft/Mars/manalysis/MBlindPixelsCalc2.cc =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MBlindPixelsCalc2.cc (revision 4445) +++ (revision ) @@ -1,671 +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): Oscar Blanch 12/2001 -! Author(s): Thomas Bretz 08/2002 -! Author(s): Nadia Tonello 02/2004 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MBlindPixelsCalc2 -// -// NT:: WARNING: This class will change soon to fit into mbadpixels!! -// It does all MBlindPixelCalc did, plus check the pedestal Rms -// to define bad pixels. -// -// This is the specific image cleaning for a list of pixels. This task -// sets the pixels listed in fPixelsIdx to unused so they should not be -// used for analysis (eg calculation of hillas parameters). -// -// NT::You can use MBlindPixelsCalc2::SetCheckPedestalRms to see which -// pixels to blind. If the PedestalRms in one pixel is more than -// 3 times the mean or less than 1/3 the mean, they are set as blind. -// -// You can use MBlindPixelsCalc2::SetUseInterpolation to replace the -// blind pixels by the average of its neighbors instead of unmapping -// them. If you want to include the central pixel use -// MBlindPixelsCalc2::SetUseCentralPixel. -// NT:: The interpolation will be done only if there are at least -// 3 non-blind pixels. -// -// Example: -// MBadPixelCalc2 blind; -// blind.SetCheckPedestalRms(); -// blind.SetUseInterpolation(); -// -// If you do not use SetUseInterpolation(), the pixels selected with -// SetCheckPedestalRms will be Unmapped and not used for the image -// analysis. -// -// Input Containers: -// MCerPhotEvt -// [MBlindPixels] -// -// Output Containers: -// MCerPhotEvt -// MBlindPixels -// -///////////////////////////////////////////////////////////////////////////// -#include "MBlindPixelsCalc2.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 "MPedestalPix.h" -#include "MPedestalCam.h" -#include "MBlindPixels.h" - -#include "MMcRunHeader.hxx" - -ClassImp(MBlindPixelsCalc2); - -using namespace std; - -static const TString gsDefName = "MBlindPixelsCalc2"; -static const TString gsDefTitle = "Find hot spots (star, broken pixels, etc)"; - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MBlindPixelsCalc2::MBlindPixelsCalc2(const char *name, const char *title) - : fFlags(0), fErrors(3) -{ - 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) -// -Int_t MBlindPixelsCalc2::PreProcess (MParList *pList) -{ - if (TESTBIT(fFlags, kUseBlindPixels)) - fPixels = (MBlindPixels*)pList->FindObject(AddSerialNumber("MBlindPixels")); - else - fPixels = (MBlindPixels*)pList->FindCreateObj(AddSerialNumber("MBlindPixels")); - - if (!fPixels) - {*fLog << err << dbginf << "Try to remove -SetUseBlindPixels- from the task list to go on" << endl; - return kFALSE;} - - fEvt = (MCerPhotEvt*)pList->FindObject(AddSerialNumber("MCerPhotEvt")); - if (!fEvt) - { - *fLog << err << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - fPed = (MPedPhotCam*)pList->FindObject(AddSerialNumber("MPedPhotCam")); - if (!fPed) - { - *fLog << err << dbginf << "MPedPhotCam not found... aborting." << endl; - return kFALSE; - } - - fGeomCam = (MGeomCam*)pList->FindObject(AddSerialNumber("MGeomCam")); - if (!fGeomCam) - *fLog << warn << dbginf << "No camera geometry available... can't use interpolation." << endl; - - const UShort_t size = fPixelsIdx.GetSize(); - if (size == 0) - { - if (!pList->FindObject("MMcRunHeader") && !fPixels) - { - *fLog << warn << "Warning - Neither blind pixels are given nor a MMcRunHeader was found... removing MBlindPixelsCalc2 from list." << endl; - return kSKIP; - } - return kTRUE; - } - - // Set as blind pixels the global blind pixels, which are given - // through the macros - - UShort_t numids = fPixelsIdx.GetSize(); - - for(Int_t i = 0; iSetPixelBlind(fPixelsIdx[i]); - - - memset(fErrors.GetArray(), 0, sizeof(Char_t)*fErrors.GetSize()); - - - fErrors[0]= 0; - fErrors[1]= 0; - fErrors[2]= 0; - - 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. -// -void MBlindPixelsCalc2::Interpolate() const -{ - const UShort_t entries = fGeomCam->GetNumPixels(); - - // - // Create arrays - // - Double_t *nphot = new Double_t[entries]; - Double_t *perr = new Double_t[entries]; - - // - // Loop over all pixels - // - for (UShort_t i=0; iIsBlind(i)) - continue; - - // - // Get a pointer to this pixel. If it is not yet existing - // create a new entry for this pixel in MCerPhotEvt - // - - MCerPhotPix *pix = fEvt->GetPixById(i); - if (!pix) - pix = fEvt->AddPixel(i, 0, 0); - - // - // Get the corresponding geometry and pedestal - // - const MGeomPix &gpix = (*fGeomCam)[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(); - perr[i] = nucp ? 0 : Pow2(pix->GetErrorPhot()); - - // - // The values are rescaled to the small pixels area for the right comparison - // - const Double_t ratio = fGeomCam->GetPixRatio(i); - - nphot[i] *= ratio; - perr[i] *= ratio; - - // - // Loop over all its neighbors - // - const Int_t n = gpix.GetNumNeighbors(); - for (int j=0; jIsBlind(nidx)) - continue; - - // - // Check whether the neighbor has a signal stored - // - const MCerPhotPix *evtpix = fEvt->GetPixById(nidx); - if (!evtpix) - continue; - - // - // Get the geometry for the neighbor - // - const Double_t nratio = fGeomCam->GetPixRatio(nidx); - - // - //The error is calculated as the quadratic sum of the errors - // - nphot[i] += nratio*evtpix->GetNumPhotons(); - perr[i] += nratio* Pow2(evtpix->GetErrorPhot()); - - num++; - } - - // - // Now the mean is calculated and the values rescaled back to the pixel area - // - nphot[i] /= (num*ratio); - perr[i] = TMath::Sqrt(perr[i]/(num*ratio)); - - - // Check if there are enough neighbors to calculate the mean - // If not, unmap the pixel. The maximum number of blind neighbors should be 2 - - if (num < 3) - { pix->SetPixelUnmapped(); - // pix->Set(nphot[i], perr[i]); - continue; - } - - pix->Set(nphot[i], perr[i]); - } - // - // Delete the intermediat arrays - // - delete nphot; - delete perr; -} - -// -------------------------------------------------------------------------- -// - -void MBlindPixelsCalc2::InterpolatePedestals() const -{ - const Int_t entries = fPed->GetSize(); - Double_t *ped = new Double_t[entries]; - Double_t *pedrms = new Double_t[entries]; - - // - // Loop over all pixels - // - for (UShort_t i=0; iIsBlind(i)) - continue; - - // - // Get the corresponding geometry and pedestal - // - const MGeomPix &gpix = (*fGeomCam)[i]; - const MPedPhotPix &ppix = (*fPed)[i]; - - // Do Not-Use-Central-Pixel - const Bool_t nucp = !TESTBIT(fFlags, kUseCentralPixel); - - Int_t num = nucp ? 0 : 1; - - ped[i] = nucp ? 0 : ppix.GetMean(); - pedrms[i] = nucp ? 0 : Pow2(ppix.GetRms()); - - // - // The values are rescaled to the small pixels area for the right comparison - // - const Double_t ratio = fGeomCam->GetPixRatio(i); - - ped[i] *= ratio; - pedrms[i] *= ratio; - // - // Loop over all its neighbors - // - const Int_t n = gpix.GetNumNeighbors(); - for (int j=0; jIsBlind(nidx)) - continue; - - // - // Get the geometry for the neighbor - // - const Double_t nratio = fGeomCam->GetPixRatio(nidx); - MPedPhotPix &nppix = (*fPed)[nidx]; - - // - //The error is calculated as the quadratic sum of the errors - // - ped[i] += (nratio*nppix.GetMean()); - pedrms[i] += nratio*Pow2(nppix.GetRms()); - - num++; - } - - // Check if there are enough neighbors to calculate the mean - // If not, unmap the pixel. The minimum number of good neighbors - // should be 3 - - if (num < 3) - { - MCerPhotPix *pix = fEvt->GetPixById(i); - if (!pix) - pix = fEvt->AddPixel(i, 0, 0); - pix->SetPixelUnmapped(); - continue; - } - - // - // Now the mean is calculated and the values rescaled back to the pixel area - // - ped[i] /= (num*ratio); - pedrms[i] = TMath::Sqrt(pedrms[i]/(num*ratio)); - - (*fPed)[i].Set(ped[i], pedrms[i]); - } - - // - // Delete the arrays - // - delete ped; - delete pedrms; - - return; -} - -// -------------------------------------------------------------------------- -// -// Removes all blind pixels from the analysis by setting their state -// to unmapped. -// -void MBlindPixelsCalc2::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; iIsBlind(pix.GetPixId())) - pix.SetPixelUnmapped(); - } -} - -// -------------------------------------------------------------------------- -// Check the pedestal Rms of the pixels: compute with 2 iterations the mean -// for inner and outer pixels. Set as blind the pixels with too small or -// too high pedestal Rms with respect to the mean. -// -Bool_t MBlindPixelsCalc2::CheckPedestalRms() - -{ - const Int_t entries = fPed->GetSize(); - Int_t pixtoblind = 0; - Int_t pixin = 0; - Int_t pixout = 0; - - Float_t meanPedRms =0.; - Float_t meanPedRmsout = 0.; - Float_t meanPedRms2 = 0.; - Float_t meanPedRmsout2 = 0.; - - for (Int_t i=0; iGetPixRatioSqrt(i); - const Double_t pixPedRms = ppix.GetRms(); - - //Calculate the corrected means: - - if (nratio == 1 ) - { - if( pixPedRms >0. && pixPedRms <200. ) - {meanPedRms += pixPedRms; - pixin++;} - } - else - { - if( pixPedRms >0. && pixPedRms <500. ) - {meanPedRmsout += pixPedRms; - pixout++;} - } - } - - //if no pixel has a minimum signal, return - if( pixin==0 || pixout==0 ) - {fErrors[1]++; //no valid Pedestals Rms - return kFALSE; - } - - meanPedRms = meanPedRms / pixin ; - meanPedRmsout = meanPedRmsout / pixout ; - - if ( meanPedRms == 0. || meanPedRmsout ==0.) - {fErrors[1]++; //no valid Pedestals Rms - return kFALSE; - } - - pixin = 0; - pixout = 0; - - for (Int_t i=0; iGetPixRatioSqrt(i); - const Double_t pixPedRms = ppix.GetRms(); - - //Calculate the corrected means: - - if (nratio == 1) - {if (pixPedRms > 0.5 * meanPedRms && pixPedRms < 1.5 * meanPedRms) - {meanPedRms2 += pixPedRms; - pixin++;} - } - else - {if (pixPedRms > 0.5 * meanPedRmsout && pixPedRms < 1.5 * meanPedRmsout) - {meanPedRmsout2 += pixPedRms; - pixout++;} - } - } - - //if no pixel has a minimum signal, return - if( pixin==0 || pixout==0 ) - {fErrors[1]++; //no valid Pedestals Rms - return kFALSE;} - - meanPedRms = meanPedRms2 / pixin ; - meanPedRmsout = meanPedRmsout2 / pixout ; - - //Blind the Bad Pixels - for (Int_t i=0; iGetPixRatioSqrt(i); - const Double_t pixPedRms = ppix.GetRms(); - - if (nratio == 1) - { if (pixPedRms > 3.* meanPedRms || pixPedRms <= meanPedRms/3.) - { fPixels->SetPixelBlind(i); - pixtoblind++;} - } - else - { if (pixPedRms > 3.* meanPedRmsout || pixPedRms <= meanPedRmsout/3.) - { fPixels->SetPixelBlind(i); - pixtoblind++; } - } - } - - // Check if the number of pixels to blind is > 60% of total number of pixels - // - if (pixtoblind>0.6*entries) - { fErrors[2]++; - return kFALSE; } - - - return kTRUE; - -} - -// -------------------------------------------------------------------------- -// -// Treat the blind pixels -// -Int_t MBlindPixelsCalc2::Process() -{ - if (TESTBIT(fFlags, kCheckPedestalRms)) - { - // if the number of blind pixels is too high, do not interpolate - if (CheckPedestalRms()== kFALSE) - return kTRUE; - else - {if (TESTBIT(fFlags, kUseInterpolation)) - InterpolatePedestals(); - } - } - - if (TESTBIT(fFlags, kUseInterpolation) && fGeomCam) - { Interpolate(); - } - else - { Unmap(); - } - fErrors[0]++; - return kTRUE; -} - - -// -------------------------------------------------------------------------- -// -// - Check whether pixels to disable are available. If pixels are -// given by the user nothing more is done. -// -Bool_t MBlindPixelsCalc2::ReInit(MParList *pList) -{ - if (TESTBIT(fFlags, kUseBlindPixels)) - return kTRUE; - - // - // If pixels are given by the user, we are already done - // - if (fPixelsIdx.GetSize() > 0) - return kTRUE; - - // - // Delete the old array holding the blind pixels for the last file - // - fPixels->Clear(); - - if (!fGeomCam->InheritsFrom("MGeomCamMagic")) - { - *fLog << warn << "MBlindPixelsCalc2::ReInit: Warning - Starfield only implemented for Magic standard Camera... no action." << endl; - return kTRUE; - } - return kTRUE; -} - -void MBlindPixelsCalc2::StreamPrimitive(ofstream &out) const -{ - out << " MBlindPixelsCalc2 " << GetUniqueName(); - if (fName!=gsDefName || fTitle!=gsDefTitle) - { - out << "(\"" << fName << "\""; - if (fTitle!=gsDefTitle) - out << ", \"" << fTitle << "\""; - out <<")"; - } - out << ";" << endl; - - if (TESTBIT(fFlags, kUseInterpolation)) - out << " " << GetUniqueName() << ".SetUseInterpolation();" << endl; - if (TESTBIT(fFlags, kUseCentralPixel)) - out << " " << GetUniqueName() << ".SetUseCentralPixel();" << endl; - if (TESTBIT(fFlags, kCheckPedestalRms)) - out << " " << GetUniqueName() << ".SetCheckPedestalRms();" << endl; - - if (fPixelsIdx.GetSize()==0) - return; - - out << " {" << endl; - out << " TArrayS dummy;" << endl; - for (int i=0; i60% of pixels are bad, no INTERPOLATION"); - - //Number of events that have been succesfully checked and corrected - - *fLog << " " << (int)fErrors[0] << " (" << (int)(100.*fErrors[0]/GetNumExecutions()) << "%) Evts well treated!" << endl; - *fLog << endl; - - return kTRUE; -} - - - - - - Index: unk/MagicSoft/Mars/manalysis/MBlindPixelsCalc2.h =================================================================== --- /trunk/MagicSoft/Mars/manalysis/MBlindPixelsCalc2.h (revision 4445) +++ (revision ) @@ -1,95 +1,0 @@ -#ifndef MARS_MBlindPixelsCalc2 -#define MARS_MBlindPixelsCalc2 - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -#ifndef ROOT_TArrayS -#include -#endif - -#ifndef ROOT_TArrayL -#include -#endif - -class MGeomCam; -class MCerPhotEvt; -class MBlindPixels; - -class MPedPhotCam; - -class MBlindPixelsCalc2 : public MTask -{ -private: - MCerPhotEvt *fEvt; //! - MBlindPixels *fPixels; //! - MGeomCam *fGeomCam; //! - MPedPhotCam *fPed; //! - - TArrayS fPixelsIdx; // Pixel Indices for blind pixels, which are entered by the user. - Byte_t fFlags; // flag for the method which is used - - enum - { - kUseInterpolation = 1, - kUseCentralPixel = 2, - kCheckPedestalRms = 3, - kUseBlindPixels = 4 - }; - - - static Double_t Pow2(Double_t x) { return x*x; } - - void Interpolate() const; - void InterpolatePedestals() const; - void Unmap() const; - Bool_t CheckPedestalRms() ; - - - void PrintSkipped(int i, const char *str) const; - - void StreamPrimitive(ofstream &out) const; - - Int_t PreProcess(MParList *pList); - Int_t Process(); - Int_t PostProcess(); - - TArrayL fErrors; - -public: - - MBlindPixelsCalc2(const char *name=NULL, const char *title=NULL); - - void SetUseInterpolation(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseInterpolation) : CLRBIT(fFlags, kUseInterpolation); - } - void SetUseCentralPixel(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseCentralPixel) : CLRBIT(fFlags, kUseCentralPixel); - } - void SetCheckPedestalRms(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags,kCheckPedestalRms) : CLRBIT(fFlags, kCheckPedestalRms); - } - void SetUseBlindPixels(Bool_t b=kTRUE) - { - b ? SETBIT(fFlags, kUseBlindPixels) : CLRBIT(fFlags, kUseBlindPixels); - } - - void SetPixelIndices(const TArrayS &pix) { fPixelsIdx = pix; } - - virtual Bool_t ReInit(MParList *pList); - - ClassDef(MBlindPixelsCalc2, 1) // (temp!)task to deal with hot spots (star, broken pixels, etc) -}; - -#endif - - - - - - - Index: /trunk/MagicSoft/Mars/manalysis/Makefile =================================================================== --- /trunk/MagicSoft/Mars/manalysis/Makefile (revision 4445) +++ /trunk/MagicSoft/Mars/manalysis/Makefile (revision 4446) @@ -43,7 +43,4 @@ MCerPhotAnal2.cc \ MCerPhotCalc.cc \ - MBlindPixels.cc \ - MBlindPixelCalc.cc \ - MBlindPixelsCalc2.cc \ MEventRate.cc \ MEventRateCalc.cc \