Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 1939) +++ trunk/MagicSoft/Mars/Changelog (revision 1940) @@ -30,4 +30,26 @@ - remove all IDs in listbox - add entries sorted + + * manalysis/MImgCleanStd.[h,cc], manalysis/MCameraSmooth.[h,cc], + manalysis/MHillas.[h,cc], manalysis/MHillasSrc.[h,cc], + manalysis/MHillasSrc.[h,cc], manalysis/MHillasCalc.[h,cc], + manalysis/MNewImagePar.[h,cc], manalysis/MNewImageParCalc.[h,cc], + mhist/MHHillas.[h,cc], mhist/MHHillasSrc.[h,cc], + mhist/MHHillasExt.[h,cc], mhist/MHNewImagePar.[h,cc]: + - moved to new directory mimage + + * mimage/Makefile, mimage/ImageLinkDef.h, mimage/ImageIncl.h: + - new + + * Makefile: + - added mimage + + * manalysis/AnalysisLinkDef.h, manalysis/Makefile + mhist/HistLinkDef.h, mhist/Makefile: + - removed moved files + + * manalysis/Makefile, manalysis/AnalysisLinkDef.h, + mhist/Makefile, mhist/AnalysisLinkDef.h: + - adapted Index: trunk/MagicSoft/Mars/Makefile =================================================================== --- trunk/MagicSoft/Mars/Makefile (revision 1939) +++ trunk/MagicSoft/Mars/Makefile (revision 1940) @@ -39,4 +39,5 @@ mfileio \ mfilter \ + mimage \ manalysis \ mmc \ Index: trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1939) +++ trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 1940) @@ -9,7 +9,4 @@ #pragma link C++ class MCerPhotAnal+; #pragma link C++ class MCerPhotCalc+; - -#pragma link C++ class MImgCleanStd+; -#pragma link C++ class MCameraSmooth+; #pragma link C++ class MBlindPixels+; @@ -33,9 +30,4 @@ #pragma link C++ class MSrcPosCam+; -#pragma link C++ class MHillas+; -#pragma link C++ class MHillasSrc+; -#pragma link C++ class MHillasSrcCalc+; -#pragma link C++ class MHillasExt+; -#pragma link C++ class MHillasCalc+; #pragma link C++ class MEnergyEst+; @@ -56,7 +48,4 @@ #pragma link C++ class MPointingCorr+; -#pragma link C++ class MNewImagePar+; -#pragma link C++ class MNewImageParCalc+; - #pragma link C++ class MParameterI+; #pragma link C++ class MParameterD+; Index: trunk/MagicSoft/Mars/manalysis/MCameraSmooth.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MCameraSmooth.cc (revision 1939) +++ (revision ) @@ -1,143 +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): Thomas Bretz 08/2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MCameraSmooth // -// // -// This task fills each pixel in the camera with the average of the // -// number of cerenkov photons from its surrounding pixels. This can // -// be done using the central pixel or ignoring the central pixel. // -// // -///////////////////////////////////////////////////////////////////////////// -#include "MCameraSmooth.h" - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" - -#include "MGeomPix.h" -#include "MGeomCam.h" -#include "MCerPhotPix.h" -#include "MCerPhotEvt.h" -#include "MBlindPixels.h" - -ClassImp(MCameraSmooth); - -// -------------------------------------------------------------------------- -// -// Default constructor. You can specify the numer of loops how often a -// smoothing should be done. The default is 1. -// -MCameraSmooth::MCameraSmooth(Byte_t n, const char *name, const char *title) - : fCounts(n), fUseCentralPixel(kTRUE) -{ - fName = name ? name : "MCameraSmooth"; - fTitle = title ? title : "Task to smooth the camera"; -} - -// -------------------------------------------------------------------------- -// -// - get the MCerPhotEvt from the parlist (abort if missing) -// - get MGeomCam from the parameter list -// -Bool_t MCameraSmooth::PreProcess (MParList *pList) -{ - fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fEvt) - { - *fLog << err << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fGeomCam) - { - *fLog << warn << dbginf << "No camera geometry available... aborting." << endl; - - return kFALSE; - } - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Do the smoothing -// -Bool_t MCameraSmooth::Process() -{ - const UShort_t entries = fEvt->GetNumPixels(); - - // - // remove the pixels in fPixelsID if they are set to be used, - // (set them to 'unused' state) - // - - Double_t *photons = new Double_t[entries]; - Double_t *errors = new Double_t[entries]; - - for (int n=0; nGetPixById(nid); - if (evtpix) - { - photons[i] += evtpix->GetNumPhotons(); - errors[i] += evtpix->GetErrorPhot(); - } - num++; - } - - photons[i] /= num; - errors[i] /= num; - } - - for (UShort_t i=0; i -#endif - -class MGeomCam; -class MCerPhotEvt; -class MBlindPixels; - -class MCameraSmooth : public MTask -{ -private: - MCerPhotEvt *fEvt; //! Pointer to data to smooth - MGeomCam *fGeomCam; //! Camera geometry - - Byte_t fCounts; // number of smoothing loops - Bool_t fUseCentralPixel; // use also the central pixel - -public: - MCameraSmooth(Byte_t cnt=1, const char *name=NULL, const char *title=NULL); - - void SetUseCetralPixel(Bool_t b=kTRUE) { fUseCentralPixel=kTRUE; } - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - - ClassDef(MCameraSmooth, 0) // task to smooth the camera contants -}; - -#endif - Index: trunk/MagicSoft/Mars/manalysis/MHillas.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillas.cc (revision 1939) +++ (revision ) @@ -1,432 +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): Thomas Bretz 12/2000 -! Author(s): Harald Kornmayer 1/2001 -! Author(s): Rudolf Bock 10/2001 -! Author(s): Wolfgang Wittek 6/2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHillas -// -// Storage Container for image parameters -// -// basic image parameters -// -// Version 1: -// ---------- -// fLength [mm] major axis of ellipse -// fWidth [mm] minor axis -// fDelta [rad] angle of major axis with x-axis -// by definition the major axis is pointing into -// the hemisphere x>0, thus -pi/2 < delta < pi/2 -// fSize [#CerPhot] total sum of pixels -// fMeanX [mm] x of center of ellipse -// fMeanY [mm] y of center of ellipse -// -// Version 2: -// ---------- -// fNumCorePixels number of pixels called core -// fNumUsedPixels number of pixels which survived the cleaning -// -///////////////////////////////////////////////////////////////////////////// -#include "MHillas.h" - -#include - -#include -#include - -#include "MGeomPix.h" -#include "MGeomCam.h" - -#include "MCerPhotPix.h" -#include "MCerPhotEvt.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHillas); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MHillas::MHillas(const char *name, const char *title) : fEllipse(NULL) -{ - fName = name ? name : "MHillas"; - fTitle = title ? title : "Storage container for image parameters of one event"; - - Reset(); - - fEllipse = new TEllipse; -} - -// -------------------------------------------------------------------------- -// -// Destructor. Deletes the TEllipse if one exists. -// -MHillas::~MHillas() -{ - Clear(); -} - -// -------------------------------------------------------------------------- -// -// Initializes the values with defaults. For the default values see the -// source code. -// -void MHillas::Reset() -{ - fLength = -1; - fWidth = -1; - fDelta = 0; - - fSize = -1; - fMeanX = -1; - fMeanY = -1; - - fNumUsedPixels = -1; - fNumCorePixels = -1; - - Clear(); -} - -// -------------------------------------------------------------------------- -// -// Print the hillas Parameters to *fLog -// -void MHillas::Print(Option_t *) const -{ - Double_t atg = atan2(fMeanY, fMeanX)*kRad2Deg; - - if (atg<0) - atg += 180; - - *fLog << all; - *fLog << "Basic Image Parameters (" << GetName() << ")" << endl; - *fLog << " - Length [mm] = " << fLength << endl; - *fLog << " - Width [mm] = " << fWidth << endl; - *fLog << " - Delta [deg] = " << fDelta*kRad2Deg << endl; - *fLog << " - Size [1] = " << fSize << " #CherPhot" << endl; - *fLog << " - Meanx [mm] = " << fMeanX << endl; - *fLog << " - Meany [mm] = " << fMeanY << endl; - *fLog << " - atg(y/x) [deg] = " << atg << endl; - *fLog << " - Used Pixels [#] = " << fNumUsedPixels << " Pixels" << endl; - *fLog << " - Core Pixels [#] = " << fNumCorePixels << " Pixels" << endl; -} - -/* -// ----------------------------------------------------------- -// -// call the Paint function of the Ellipse if a TEllipse exists -// -void MHillas::Paint(Option_t *) -{ - fEllipse->SetLineWidth(2); - fEllipse->PaintEllipse(fMeanX, fMeanY, fLength, fWidth, - 0, 360, fDelta*kRad2Deg+180); -} -*/ - -// -------------------------------------------------------------------------- -// -// Instead of adding MHillas itself to the Pad -// (s. AppendPad in TObject) we create an ellipse, -// which is added to the Pad by its Draw function -// You can remove it by deleting the Ellipse Object -// (s. Clear() ) -// -void MHillas::Draw(Option_t *opt) -{ - - Clear(); - - if (fLength<0 || fWidth<0) - return; - - fEllipse = new TEllipse(fMeanX, fMeanY, fLength, fWidth, - 0, 360, fDelta*kRad2Deg+180); - - fEllipse->SetLineWidth(2); - fEllipse->Draw(); - - /* - fEllipse->SetPhimin(); - fEllipse->SetPhimax(); - fEllipse->SetR1(fLength); - fEllipse->SetR2(fWidth); - fEllipse->SetTheta(fDelta*kRad2Deg+180); - fEllipse->SetX1(fMeanX); - fEllipse->SetY1(fMeanY); - - fEllipse->SetLineWidth(2); - fEllipse->PaintEllipse(fMeanX, fMeanY, fLength, fWidth, - 0, 360, fDelta*kRad2Deg+180); - - AppendPad(opt); - - // This is from TH1 - TString opt = option; - opt.ToLower(); - if (gPad && !opt.Contains("same")) { - //the following statement is necessary in case one attempts to draw - //a temporary histogram already in the current pad - if (TestBit(kCanDelete)) gPad->GetListOfPrimitives()->Remove(this); - gPad->Clear(); - } - */ -} - -// -------------------------------------------------------------------------- -// -// If a TEllipse object exists it is deleted -// -void MHillas::Clear(Option_t *) -{ - if (!fEllipse) - return; - - delete fEllipse; - - fEllipse = NULL; -} - - -// -------------------------------------------------------------------------- -// -// Calculate the image parameters from a Cherenkov photon event -// assuming Cher.photons/pixel and pixel coordinates are given -// In case you don't call Calc from within an eventloop make sure, that -// you call the Reset member function before. -// Returns: -// 0 no error -// 1 number of pixels < 3 -// 2 size==0 -// 3 number of used pixel < 3 -// 4 CorrXY == 0 -// -Int_t MHillas::Calc(const MGeomCam &geom, const MCerPhotEvt &evt) -{ - const UInt_t npixevt = evt.GetNumPixels(); - - // - // sanity check - // - if (npixevt < 3) - { - //*fLog << warn << "MHillas::Calc: Event has less than three pixels... skipped." << endl; - return 1; - } - - // - // calculate mean value of pixel coordinates and fSize - // ----------------------------------------------------- - // - // Because this are only simple sums of roughly 600 values - // with an accuracy less than three digits there should not - // be any need to use double precision (Double_t) for the - // calculation of fMeanX, fMeanY and fSize. - // - fMeanX = 0; - fMeanY = 0; - fSize = 0; - - fNumUsedPixels = 0; - fNumCorePixels = 0; - - for (UInt_t i=0; i0) we - // cannot calculate Length and Width. The calculation failed - // and returns kFALSE - // In reallity it is almost impossible to have a distribution - // of cerenkov photons in the used pixels which is exactly symmetric - // along one of the axis. - // - if (corrxy==0) - { - //*fLog << inf << GetDescriptor() << ": Event has CorrXY==0... skipped." << endl; - return 4; - } - - // - // calculate the basic Hillas parameters: orientation and size of axes - // ------------------------------------------------------------------- - // - // fDelta is the angle between the major axis of the ellipse and - // the x axis. It is independent of the position of the ellipse - // in the camera it has values between -pi/2 and pi/2 degrees - // - const Double_t d0 = corryy - corrxx; - const Double_t d1 = corrxy*2; - const Double_t d2 = d0 + sqrt(d0*d0 + d1*d1); - const Double_t tand = d2 / d1; - - fDelta = atan(tand); - - const Double_t s2 = tand*tand+1; - const Double_t s = sqrt(s2); - - fCosDelta = 1.0/s; // need these in derived classes - fSinDelta = tand/s; // like MHillasExt - - Double_t axis1 = (tand*tand*corryy + d2 + corrxx)/s2/fSize; - Double_t axis2 = (tand*tand*corrxx - d2 + corryy)/s2/fSize; - - // - // fLength^2 is the second moment along the major axis of the ellipse - // fWidth^2 is the second moment along the minor axis of the ellipse - // - // From the algorithm we get: fWidth <= fLength is always true - // - // very small numbers can get negative by rounding - // - fLength = axis1<0 ? 0 : sqrt(axis1); // [mm] - fWidth = axis2<0 ? 0 : sqrt(axis2); // [mm] - - SetReadyToSave(); - - return 0; -} - -// -------------------------------------------------------------------------- -// -// This function is ment for special usage, please never try to set -// values via this function -// -void MHillas::Set(const TArrayF &arr) -{ - if (arr.GetSize() != 8) - return; - - fLength = arr.At(0); // [mm] major axis of ellipse - fWidth = arr.At(1); // [mm] minor axis of ellipse - fDelta = arr.At(2); // [rad] angle of major axis with x-axis - fSize = arr.At(3); // [#CerPhot] sum of content of all pixels (number of Cherenkov photons) - fMeanX = arr.At(4); // [mm] x-coordinate of center of ellipse - fMeanY = arr.At(5); // [mm] y-coordinate of center of ellipse - - fNumUsedPixels = (Short_t)arr.At(6); // Number of pixels which survived the image cleaning - fNumCorePixels = (Short_t)arr.At(7); // number of core pixels -} - - -// -------------------------------------------------------------------------- -// -/* -void MHillas::AsciiRead(ifstream &fin) -{ - fin >> fLength; - fin >> fWidth; - fin >> fDelta; - fin >> fSize; - fin >> fMeanX; - fin >> fMeanY; -} -*/ -// -------------------------------------------------------------------------- -/* -void MHillas::AsciiWrite(ofstream &fout) const -{ - fout << fLength << " "; - fout << fWidth << " "; - fout << fDelta << " "; - fout << fSize << " "; - fout << fMeanX << " "; - fout << fMeanY; -} -*/ Index: trunk/MagicSoft/Mars/manalysis/MHillas.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillas.h (revision 1939) +++ (revision ) @@ -1,73 +1,0 @@ -#ifndef MARS_MHillas -#define MARS_MHillas - -#ifndef MARS_MParContainer -#include "MParContainer.h" -#endif - -class TArrayF; -class TEllipse; - -class MGeomCam; -class MCerPhotEvt; - -class MHillas : public MParContainer -{ -private: - // for description see MHillas.cc - Float_t fLength; // [mm] major axis of ellipse - Float_t fWidth; // [mm] minor axis of ellipse - Float_t fDelta; // [rad] angle of major axis with x-axis - Float_t fSize; // [#CerPhot] sum of content of all pixels (number of Cherenkov photons) - Float_t fMeanX; // [mm] x-coordinate of center of ellipse - Float_t fMeanY; // [mm] y-coordinate of center of ellipse - - Short_t fNumUsedPixels; // Number of pixels which survived the image cleaning - Short_t fNumCorePixels; // number of core pixels - - Float_t fSinDelta; //! [1] sin of Delta (to be used in derived classes) - Float_t fCosDelta; //! [1] cos of Delta (to be used in derived classes) - - TEllipse *fEllipse; //! Graphical Object to Display Ellipse - -protected: - // - // This is only for calculations in derived classes because - // we don't want to read/write this data members - // - Float_t GetCosDelta() const { return fCosDelta; } - Float_t GetSinDelta() const { return fSinDelta; } - -public: - MHillas(const char *name=NULL, const char *title=NULL); - ~MHillas(); - - void Reset(); - - virtual Int_t Calc(const MGeomCam &geom, const MCerPhotEvt &pix); - - virtual void Print(Option_t *opt=NULL) const; - virtual void Draw(Option_t *opt=NULL); - //virtual void Paint(Option_t *); - - virtual void Clear(Option_t *opt=NULL); - - Float_t GetLength() const { return fLength; } - Float_t GetWidth() const { return fWidth; } - Float_t GetDelta() const { return fDelta; } - Float_t GetSize() const { return fSize; } - Float_t GetMeanX() const { return fMeanX; } - Float_t GetMeanY() const { return fMeanY; } - - Int_t GetNumUsedPixels() const { return fNumUsedPixels; } - Int_t GetNumCorePixels() const { return fNumCorePixels; } - - virtual void Set(const TArrayF &arr); - - //virtual void AsciiRead(ifstream &fin); - //virtual void AsciiWrite(ofstream &fout) const; - - ClassDef(MHillas, 2) // Storage Container for Hillas Parameter -}; - -#endif Index: trunk/MagicSoft/Mars/manalysis/MHillasCalc.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasCalc.cc (revision 1939) +++ (revision ) @@ -1,136 +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): Thomas Bretz 12/2000 -! Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de) -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MHillasCalc // -// // -// This is a task to calculate the Hillas parameters from each event // -// // -// Input Containers: // -// MCerPhotEvt, MGeomCam // -// // -// Output Containers: // -// MHillas // -// // -///////////////////////////////////////////////////////////////////////////// - -#include "MHillasCalc.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MCerPhotEvt.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHillasCalc); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MHillasCalc::MHillasCalc(const char *hil, const char *name, const char *title) -{ - fName = name ? name : "MHillasCalc"; - fTitle = title ? title : "Task to calculate Hillas parameters"; - - fHilName = hil; -} - -// -------------------------------------------------------------------------- -// -// Check for a MCerPhotEvt object from which the Hillas are calculated. -// Try to find the Geometry conatiner. And try to find the output -// (Hillas) container or create one. -// -Bool_t MHillasCalc::PreProcess(MParList *pList) -{ - fCerPhotEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fCerPhotEvt) - { - *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fGeomCam) - { - *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; - return kFALSE; - } - - fHillas = (MHillas*)pList->FindCreateObj("MHillas",fHilName); - if (!fHillas) - return kFALSE; - - memset(fErrors, 0, sizeof(fErrors)); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// If you want do complex descisions inside the calculations -// we must move the calculation code inside this function -// -// If the calculation wasn't sucessfull skip this event -// -Bool_t MHillasCalc::Process() -{ - const Int_t rc = fHillas->Calc(*fGeomCam, *fCerPhotEvt); - if (rc<0 || rc>4) - { - *fLog << err << dbginf << "MHillas::Calc returned unknown error code!" << endl; - return kFALSE; - } - fErrors[rc]++; - - return rc==0 ? kTRUE : kCONTINUE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the hillas calculation. The percentage -// is calculated with respect to the number of executions of this task. -// -Bool_t MHillasCalc::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << setw(7) << fErrors[1] << " (" << setw(3) << (int)(fErrors[1]*100/GetNumExecutions()) << "%) Evts skipped due to: Event has less than 3 pixels" << endl; - *fLog << " " << setw(7) << fErrors[2] << " (" << setw(3) << (int)(fErrors[2]*100/GetNumExecutions()) << "%) Evts skipped due to: Calculated Size == 0" << endl; - *fLog << " " << setw(7) << fErrors[3] << " (" << setw(3) << (int)(fErrors[3]*100/GetNumExecutions()) << "%) Evts skipped due to: Number of used pixels < 3" << endl; - *fLog << " " << setw(7) << fErrors[4] << " (" << setw(3) << (int)(fErrors[4]*100/GetNumExecutions()) << "%) Evts skipped due to: CorrXY==0" << endl; - *fLog << " " << fErrors[0] << " (" << (int)(fErrors[0]*100/GetNumExecutions()) << "%) Evts survived Hillas calculation!" << endl; - *fLog << endl; - - return kTRUE; -} Index: trunk/MagicSoft/Mars/manalysis/MHillasCalc.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasCalc.h (revision 1939) +++ (revision ) @@ -1,39 +1,0 @@ -#ifndef MARS_MHillasCalc -#define MARS_MHillasCalc - -///////////////////////////////////////////////////////////////////////////// -// // -// MHillasCalkc // -// // -// Task to calculate Hillas Parameters // -// // -///////////////////////////////////////////////////////////////////////////// - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MGeomCam; -class MCerPhotEvt; -class MHillas; - -class MHillasCalc : public MTask -{ - const MGeomCam *fGeomCam; // Camera Geometry used to calculate Hillas - const MCerPhotEvt *fCerPhotEvt; // Cerenkov Photon Event used for calculation - MHillas *fHillas; // ouput container to store result - - TString fHilName; - Int_t fErrors[5]; - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - Bool_t PostProcess(); - -public: - MHillasCalc(const char *hil="MHillas", const char *name=NULL, const char *title=NULL); - - ClassDef(MHillasCalc, 0) // Task to calculate Hillas parameters -}; - -#endif Index: trunk/MagicSoft/Mars/manalysis/MHillasExt.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasExt.cc (revision 1939) +++ (revision ) @@ -1,291 +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): Thomas Bretz 12/2000 -! Author(s): Rudolf Bock 10/2001 -! Author(s): Wolfgang Wittek 06/2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHillasExt -// -// Storage Container for extended image parameters -// -// extended image parameters -// fConc ratio of sum of two highest pixels over fSize -// fConc1 ratio of highest pixel over fSize -// fAsym distance from highest pixel to center, projected onto major axis -// fM3Long third moment along major axis -// fM3Trans third moment along minor axis -// -// WARNING: Before you can use fAsym, fM3Long and fM3Trans you must -// multiply by the sign of MHillasSrc::fCosDeltaAlpha -// -//////////////////////////////////////////////////////////////////////////// -/* - // fLeakage1 ratio : (photons in most outer ring of pixels) over fSize - // fLeakage2 ratio : (photons in the 2 outer rings of pixels) over fSize - // - // fAsymna d/(d na) of ( sum(x*q^na)/sum(q^na), sum(y*q^na)/sum(q^na) ) - // projected onto the major axis - // fAsym0 (F-B)/(F+B) along the major axis - */ -#include "MHillasExt.h" - -#include -#include - -#include "MGeomPix.h" -#include "MGeomCam.h" - -#include "MCerPhotPix.h" -#include "MCerPhotEvt.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHillasExt); - -// ------------------------------------------------------------------------- -// -// Default constructor. -// -MHillasExt::MHillasExt(const char *name, const char *title) -{ - fName = name ? name : "MHillas"; - fTitle = title ? title : "Storage container for extended parameter set of one event"; - - Reset(); - // FIXME: (intelligent) initialization of values missing -} - -// ------------------------------------------------------------------------- -// -void MHillasExt::Reset() -{ - MHillas::Reset(); - - fConc = -1; - fConc1 = -1; - fAsym = 0; - fM3Long = 0; - fM3Trans = 0; -} - -// ------------------------------------------------------------------------- -// -void MHillasExt::Print(Option_t *) const -{ - MHillas::Print(); - - *fLog << "Extended Image Parameters (" << GetName() << ")" << endl; - *fLog << " - Conc = " << fConc << " (ratio)" << endl; - *fLog << " - Conc1 = " << fConc1 << " (ratio)" << endl; - *fLog << " - Asymmetry = " << fAsym << " mm" << endl; - *fLog << " - 3rd Moment Long = " << fM3Long << " mm" << endl; - *fLog << " - 3rd Moment Trans = " << fM3Trans << " mm" << endl; -} - -// ------------------------------------------------------------------------- -// -// calculation of additional parameters based on the camera geometry -// and the cerenkov photon event -// -Int_t MHillasExt::Calc(const MGeomCam &geom, const MCerPhotEvt &evt) -{ - const Int_t rc = MHillas::Calc(geom, evt); - if (rc>0) - return rc; - - // - // calculate the additional image parameters - // -------------------------------------------- - // - // loop to get third moments along ellipse axes and two max pixels - // - // For the moments double precision is used to make sure, that - // the complex matrix multiplication and sum is evaluated correctly. - // - Double_t m3x = 0; - Double_t m3y = 0; - - Float_t maxpix1 = 0; // [#phot] - Float_t maxpix2 = 0; // [#phot] - - Int_t maxpixid = 0; - - const UInt_t npixevt = evt.GetNumPixels(); - - const Float_t A0 = geom[0].GetA(); - - for (UInt_t i=0; imaxpix1) - { - maxpix2 = maxpix1; - maxpix1 = nphot; // [1] - maxpixid = pix.GetPixId(); - continue; // [1] - } - - if (nphot>maxpix2) - maxpix2 = nphot; // [1] - - /* - // - // power na for calculating fAsymna; - // the value 1.5 was suggested by Thomas Schweizer - // - Double_t na = 1.5; - - // - // get sums for calculating fAsymna - // the outer pixels are 4 times as big (in area) - // as the inner pixels ! - // - const Double_t dummy = pow(nphot, na)/r; - - sna += dummy; - xna += dzx*dummy; - - sna1 += sna/nphot; - xna1 += xna/nphot; - - // - // forward-backward asymmetry - // - fb += dzx<0 ? -nphot: nphot; - */ - } - - const MGeomPix &maxpix = geom[maxpixid]; - - fAsym = (GetMeanX()-maxpix.GetX())*GetCosDelta() + - (GetMeanY()-maxpix.GetY())*GetSinDelta(); // [mm] - - fConc = (maxpix1+maxpix2)/GetSize(); // [ratio] - fConc1 = maxpix1/GetSize(); // [ratio] - - /* - fLeakage1 = edgepix1 / GetSize(); - fLeakage2 = edgepix2 / GetSize(); - fAsym0 = fb / GetSize(); - - fAsymna = na * (sna*xna1 - sna1*xna) / (sna*sna); - */ - - // - // Third moments along axes get normalized - // - m3x /= GetSize(); - m3y /= GetSize(); - - fM3Long = m3x<0 ? -pow(-m3x, 1./3) : pow(m3x, 1./3); // [mm] - fM3Trans = m3y<0 ? -pow(-m3y, 1./3) : pow(m3y, 1./3); // [mm] - - SetReadyToSave(); - - return 0; -} - -// -------------------------------------------------------------------------- -// -// This function is ment for special usage, please never try to set -// values via this function -// -void MHillasExt::Set(const TArrayF &arr) -{ - if (arr.GetSize() != 13) - return; - - fConc = arr.At(8); // [ratio] concentration ratio: sum of the two highest pixels / fSize - fConc1 = arr.At(9); // [ratio] concentration ratio: sum of the highest pixel / fSize - fAsym = arr.At(10); // [mm] fDist minus dist: center of ellipse, highest pixel - fM3Long = arr.At(11); // [mm] 3rd moment (e-weighted) along major axis - fM3Trans = arr.At(12); // [mm] 3rd moment (e-weighted) along minor axis - - TArrayF n(arr); - n.Set(8); - MHillas::Set(n); -} - -/* -// ------------------------------------------------------------------------- -// -void MHillasExt::AsciiRead(ifstream &fin) -{ - MHillas::AsciiRead(fin); - - fin >> fConc; - fin >> fConc1; - fin >> fAsym; - fin >> fM3Long; - fin >> fM3Trans; -} -*/ -// ------------------------------------------------------------------------- -/* -void MHillasExt::AsciiWrite(ofstream &fout) const -{ - MHillas::AsciiWrite(fout); - - fout << " "; - fout << fConc << " "; - fout << fConc1 << " "; - fout << fAsym << " "; - fout << fM3Long << " "; - fout << fM3Trans; -} -*/ Index: trunk/MagicSoft/Mars/manalysis/MHillasExt.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasExt.h (revision 1939) +++ (revision ) @@ -1,43 +1,0 @@ -#ifndef MARS_MHillasExt -#define MARS_MHillasExt - -#ifndef MARS_MHillas -#include "MHillas.h" -#endif - -class MGeomCam; -class MCerPhotEvt; - -class MHillasExt : public MHillas -{ -private: - // for description see MExtHillas.cc - Float_t fConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize - Float_t fConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize - Float_t fAsym; // [mm] fDist minus dist: center of ellipse, highest pixel - Float_t fM3Long; // [mm] 3rd moment (e-weighted) along major axis - Float_t fM3Trans; // [mm] 3rd moment (e-weighted) along minor axis - -public: - MHillasExt(const char *name=NULL, const char *title=NULL); - - void Reset(); - - Float_t GetConc() const { return fConc; } - Float_t GetConc1() const { return fConc1; } - Float_t GetAsym() const { return fAsym; } - Float_t GetM3Long() const { return fM3Long; } - Float_t GetM3Trans() const { return fM3Trans; } - - Int_t Calc(const MGeomCam &geom, const MCerPhotEvt &pix); - - void Print(Option_t *opt=NULL) const; - - void Set(const TArrayF &arr); - - //void AsciiRead(ifstream &fin); - //void AsciiWrite(ofstream &fout) const; - - ClassDef(MHillasExt, 1) // Storage Container for extended Hillas Parameter -}; -#endif Index: trunk/MagicSoft/Mars/manalysis/MHillasSrc.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasSrc.cc (revision 1939) +++ (revision ) @@ -1,181 +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): Thomas Bretz 12/2000 -! Author(s): Harald Kornmayer 1/2001 -! Author(s): Rudolf Bock 10/2001 -! Author(s): Wolfgang Wittek 06/2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHillasSrc -// -// Storage Container for image parameters -// -// source-dependent image parameters -// -// Version 1: -// ---------- -// fAlpha angle between major axis and line source-to-center -// fDist distance from source to center of ellipse -// -// Version 2: -// ---------- -// fHeadTail -// -// Version 3: -// ---------- -// fCosDeltaAlpha cosine of angle between d and a, where -// - d is the vector from the source position to the -// center of the ellipse -// - a is a vector along the main axis of the ellipse, -// defined with positive x-component -// -///////////////////////////////////////////////////////////////////////////// -#include "MHillasSrc.h" - -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MSrcPosCam.h" - -ClassImp(MHillasSrc); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MHillasSrc::MHillasSrc(const char *name, const char *title) -{ - fName = name ? name : "MHillasSrc"; - fTitle = title ? title : "Parameters depending in source position"; -} - -void MHillasSrc::Reset() -{ - fDist = -1; - fAlpha = 0; - fHeadTail = 0; - fCosDeltaAlpha = 0; -} - -// -------------------------------------------------------------------------- -// -// Calculation of source-dependent parameters -// In case you don't call Calc from within an eventloop make sure, that -// you call the Reset member function before. -// -Bool_t MHillasSrc::Calc(const MHillas *hillas) -{ - fHillas = hillas; - - const Double_t mx = GetMeanX(); // [mm] - const Double_t my = GetMeanY(); // [mm] - - const Double_t sx = mx - fSrcPos->GetX(); // [mm] - const Double_t sy = my - fSrcPos->GetY(); // [mm] - - const Double_t sd = sin(GetDelta()); // [1] - const Double_t cd = cos(GetDelta()); // [1] - - const Double_t tand = tan(GetDelta()); // [1] - - fHeadTail = cd*sx + sd*sy; // [mm] - - // - // Distance from source position to center of ellipse. - // If the distance is 0 distance, Alpha is not specified. - // The calculation has failed and returnes kFALSE. - // - Double_t dist = sqrt(sx*sx + sy*sy); // [mm] - - if (dist==0) - { - //*fLog << warn << GetDescriptor() << ": Event has Dist==0... skipped." << endl; - return kFALSE; - } - - // - // Calculate Alpha and Cosda = cos(d,a) - // The sign of Cosda will be used for quantities containing - // a head-tail information - // - const Double_t arg = (sy-tand*sx) / (dist*sqrt(tand*tand+1)); - - fAlpha = asin(arg)*kRad2Deg; // [deg] - fCosDeltaAlpha = fHeadTail/dist; // [1] - fDist = dist; // [mm] - - SetReadyToSave(); - - return kTRUE; -} - -void MHillasSrc::Print(Option_t *) const -{ - *fLog << all; - *fLog << "Source dependant Image Parameters (" << GetName() << ")" << endl; - *fLog << " - Dist [mm] = " << fDist << endl; - *fLog << " - Alpha [deg] = " << fAlpha << endl; - *fLog << " - HeadTail [mm] = " << fHeadTail << endl; - *fLog << " - CosDeltaAlpha = " << fCosDeltaAlpha << endl; -} -// -------------------------------------------------------------------------- -// -// This function is ment for special usage, please never try to set -// values via this function -// -void MHillasSrc::Set(const TArrayF &arr) -{ - if (arr.GetSize() != 4) - return; - - fAlpha = arr.At(0); // [deg] angle of major axis with vector to src - fDist = arr.At(1); // [mm] distance between src and center of ellipse - fHeadTail = arr.At(2); // [mm] - fCosDeltaAlpha = arr.At(3); // [1] cosine of angle between d and a -} - -// ----------------------------------------------------------------------- -// -// overloaded MParContainer to read MHillasSrc from an ascii file -// -/* -void MHillasSrc::AsciiRead(ifstream &fin) -{ - fin >> fAlpha; - fin >> fDist; - fin >> fHeadTail; -} -*/ -// ----------------------------------------------------------------------- -// -// overloaded MParContainer to write MHillasSrc to an ascii file -/* -void MHillasSrc::AsciiWrite(ofstream &fout) const -{ - fout << fAlpha << " " << fDist; -} -*/ Index: trunk/MagicSoft/Mars/manalysis/MHillasSrc.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasSrc.h (revision 1939) +++ (revision ) @@ -1,52 +1,0 @@ -#ifndef MARS_MHillasSrc -#define MARS_MHillasSrc - -#ifndef MARS_MHillas -#include "MHillas.h" -#endif - -class MSrcPosCam; - -class MHillasSrc : public MParContainer -{ -private: - const MHillas *fHillas; //! Input parameters - const MSrcPosCam *fSrcPos; //! Source position in the camere - - Float_t fAlpha; // [deg] angle of major axis with vector to src - Float_t fDist; // [mm] distance between src and center of ellipse - Float_t fHeadTail; // [mm] - Float_t fCosDeltaAlpha; // [1] cosine of angle between d and a - -public: - MHillasSrc(const char *name=NULL, const char *title=NULL); - - void SetSrcPos(const MSrcPosCam *pos) { fSrcPos = pos; } - const MSrcPosCam *GetSrcPos() const { return fSrcPos; } - - void Reset(); - - Float_t GetLength() const { return fHillas->GetLength(); } - Float_t GetWidth() const { return fHillas->GetWidth(); } - Float_t GetDelta() const { return fHillas->GetDelta(); } - Float_t GetSize() const { return fHillas->GetSize(); } - Float_t GetMeanX() const { return fHillas->GetMeanX(); } - Float_t GetMeanY() const { return fHillas->GetMeanY(); } - Float_t GetAlpha() const { return fAlpha; } - Float_t GetDist() const { return fDist; } - Float_t GetHeadTail() const { return fHeadTail; } - Float_t GetCosDeltaAlpha() const { return fCosDeltaAlpha; } - - void Print(Option_t *opt=NULL) const; - - virtual Bool_t Calc(const MHillas *hillas); - - void Set(const TArrayF &arr); - - //virtual void AsciiRead(ifstream &fin); - //virtual void AsciiWrite(ofstream &fout) const; - - ClassDef(MHillasSrc, 3) // Container to hold source position dependant parameters -}; - -#endif Index: trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.cc (revision 1939) +++ (revision ) @@ -1,183 +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): Thomas Bretz 12/2000 -! Author(s): Rudolf Bock 10/2001 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -////////////////////////////////////////////////////////////////////////////// -// -// MHillasSrcCalc -// -// Task to calculate the source dependant part of the hillas parameters -// -////////////////////////////////////////////////////////////////////////////// -#include "MHillasSrcCalc.h" - -#include - -#include "MParList.h" - -#include "MSrcPosCam.h" -#include "MHillasSrc.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHillasSrcCalc); - -static const TString gsDefName = "MHillasSrcCalc"; -static const TString gsDefTitle = "Calculate position dependant image parameters"; - -// ------------------------------------------------------------------------- -// -// Default constructor. The first argument is the name of a container -// containing the source position in the camera plain (MScrPosCam). -// The default is "MSrcPosCam". hil is the name of a container -// of type MHillasSrc (or derived) in which the parameters are stored -// The default is "MHillasSrc" -// -//MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, -// const char *name, const char *title) -// : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) -//{ -// fName = name ? name : gsDefName.Data(); -// fTitle = title ? title : gsDefTitle.Data(); - -// fSrcName = src; -// fHillasName = hil; -// fHillasInput = "MHillas"; -//} -// ------------------------------------------------------------------------- -// -MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, - const char *name, const char *title) - : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) -{ - fName = name ? name : gsDefName.Data(); - fTitle = title ? title : gsDefTitle.Data(); - - fSrcName = src; - fHillasName = hil; - fHillasInput = "MHillas"; -} - -// ------------------------------------------------------------------------- -// -Bool_t MHillasSrcCalc::PreProcess(MParList *pList) -{ - fHillas = (MHillas*)pList->FindObject(fHillasInput, "MHillas"); - if (!fHillas) - { - *fLog << err << dbginf << "MHillas not found... aborting." << endl; - return kFALSE; - } - - fSrcPos = (MSrcPosCam*)pList->FindObject(fSrcName, "MSrcPosCam"); - if (!fSrcPos) - { - *fLog << err << dbginf << fSrcName << " [MSrcPosCam] not found... aborting." << endl; - return kFALSE; - } - - fHillasSrc = (MHillasSrc*)pList->FindCreateObj("MHillasSrc", fHillasName); - if (!fHillasSrc) - return kFALSE; - - fHillasSrc->SetSrcPos(fSrcPos); - - fErrors = 0; - - return kTRUE; -} - -// ------------------------------------------------------------------------- -// -Bool_t MHillasSrcCalc::Process() -{ - - if (!fHillasSrc->Calc(fHillas)) - { - fErrors++; - return kCONTINUE; - - } - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the hillas calculation. The percentage -// is calculated with respect to the number of executions of this task. -// -Bool_t MHillasSrcCalc::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << fErrors << " (" << (int)(fErrors*100/GetNumExecutions()) << "%) Evts skipped due to: Dist==0" << endl; - *fLog << endl; - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Implementation of SavePrimitive. Used to write the call to a constructor -// to a macro. In the original root implementation it is used to write -// gui elements to a macro-file. -// -void MHillasSrcCalc::StreamPrimitive(ofstream &out) const -{ - if (fHillas) - fHillas->SavePrimitive(out); - - if (fSrcPos) - fSrcPos->SavePrimitive(out); - - if (fHillasSrc) - fHillasSrc->SavePrimitive(out); - - out << " MHillasSrcCalc " << GetUniqueName() << "("; - - if (fSrcPos) - out << "&" << fSrcPos->GetUniqueName(); - else - out << "\"" << fSrcName << "\""; - - out << ", "; - - if (fHillasSrc) - out << "&" << fHillasSrc->GetUniqueName(); - else - out << "\"" << fHillasName << "\""; - - if (fName!=gsDefName || fTitle!=gsDefTitle) - { - out << ", \"" << fName << "\""; - if (fTitle!=gsDefTitle) - out << ", \"" << fTitle << "\""; - } - out << ");" << endl; -} Index: trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.h (revision 1939) +++ (revision ) @@ -1,44 +1,0 @@ -#ifndef MARS_MHillasSrcCalc -#define MARS_MHillasSrcCalc - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MHillas; -class MHillasSrc; -class MSrcPosCam; - -class MHillasSrcCalc : public MTask -{ -private: - MHillas *fHillas; //! Pointer to the source independant hillas parameters - MSrcPosCam *fSrcPos; //! Pointer to the source position - MHillasSrc *fHillasSrc; //! Pointer to the output container for the source dependant parameters - - TString fSrcName; - TString fHillasName; - TString fHillasInput; - - Int_t fErrors; - - void StreamPrimitive(ofstream &out) const; - - Bool_t PreProcess(MParList *plist); - Bool_t Process(); - Bool_t PostProcess(); - -public: - // MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", - // const char *name=NULL, const char *title=NULL); - - MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", - const char *name=NULL, const char *title=NULL); - - void SetInput(TString hilname) { fHillasInput = hilname; } - - ClassDef(MHillasSrcCalc, 1) // task to calculate the source position depandant hillas parameters -}; - -#endif - Index: trunk/MagicSoft/Mars/manalysis/MImgCleanStd.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MImgCleanStd.cc (revision 1939) +++ (revision ) @@ -1,512 +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): Thomas Bretz 12/2000 -! Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de) -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MImgCleanStd // -// // -// This is the standard image cleaning. If you want to know how it works // -// Please look at the three CleanSteps and Process // -// // -// FIXME: MImgCleanStd is not yet completely optimized for speed. // -// Maybe we don't have to loop over all pixels all the time... // -// // -// Input Containers: // -// MGeomCam, MCerPhotEvt // -// // -// Output Containers: // -// -/- // -// // -///////////////////////////////////////////////////////////////////////////// -#include "MImgCleanStd.h" - -#include // atof -#include // ofstream, SavePrimitive - -#include // TGFrame -#include // TGLabel -#include // TGTextEntry - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" -#include "MGeomPix.h" -#include "MGeomCam.h" -#include "MCerPhotPix.h" -#include "MCerPhotEvt.h" - -#include "MGGroupFrame.h" // MGGroupFrame - -ClassImp(MImgCleanStd); - -enum { - kImgCleanLvl1, - kImgCleanLvl2 -}; - -static const TString gsDefName = "MImgCleanStd"; -static const TString gsDefTitle = "Task to perform a standard image cleaning"; - -// -------------------------------------------------------------------------- -// -// Default constructor. Here you can specify the cleaning levels. If you -// don't specify them the 'common standard' values 3.0 and 2.5 (sigma -// above mean) are used -// -MImgCleanStd::MImgCleanStd(const Float_t lvl1, const Float_t lvl2, - const char *name, const char *title) - : fCleanLvl1(lvl1), fCleanLvl2(lvl2) -{ - fName = name ? name : gsDefName.Data(); - fTitle = title ? title : gsDefTitle.Data(); - - Print(); -} - -// -------------------------------------------------------------------------- -// -// This method looks for all pixels with an entry (photons) -// that is three times bigger than the noise of the pixel -// (std: 3 sigma, clean level 1) -// -// -// AM 18/11/2002: now cut levels are proportional to the square root -// of the pixel area. In this way the cut corresponds to a fixed -// phe-density (otherwise, it would bias the images). -// -// Returns the maximum Pixel Id (used for ispixused in CleanStep2) -// -Int_t MImgCleanStd::CleanStep1() -{ - const Int_t entries = fEvt->GetNumPixels(); - - Int_t max = entries; - - // - // check the number of all pixels against the noise level and - // set them to 'unused' state if necessary - // - for (Int_t i=0; iGetPixRatio(id); - - // COBB: '<=' to skip entry=noise=0 - if (entry <= fCleanLvl1 * noise / ratio) - pix.SetPixelUnused(); - - if (id>max) - max = id; - } - return max; -} - -// -------------------------------------------------------------------------- -// -// check if the survived pixel have a neighbor, that also -// survived, otherwise set pixel to unused. (removes pixels without -// neighbors) -// -// takes the maximum pixel id from CleanStep1 as an argument -// -void MImgCleanStd::CleanStep2(Int_t max) -{ - const Int_t entries = fEvt->GetNumPixels(); - - // - // In the worst case we have to loop 6 times 577 times, to - // catch the behaviour of all next neighbors. Here we can gain - // much by using an array instead of checking through all pixels - // (MCerPhotEvt::IsPixelUsed) all the time. - // - Byte_t *ispixused = new Byte_t[max+1]; - memset(ispixused, 0, max+1); - - for (Int_t i=0; imax || !ispixused[id2]) - continue; - - cnt = kTRUE; - break; - } - if (cnt) - continue; -#else - Int_t cnt = 0; - for (Int_t j=0; jmax || !ispixused[id2]) - continue; - - if (cnt++>nnmax-4) - break; - } - if (cnt==nnmax-2 && nnmax>=4) - { - pix.SetPixelUsed(); - continue; - } - if (cnt>0) - continue; -#endif - - // - // check if no next neighbor has the state 'used' - // set this pixel to 'unused', too. - // - pix.SetPixelUnused(); - } - - delete ispixused; - - // - // now we declare all pixels that survive as CorePixels - // - for (Int_t i=0; iGetNumPixels(); - - for (Int_t i=0; iGetPixRatio(id)); - - if (entry <= fCleanLvl2 * noise / ratio) - continue; - - // - // check if the pixel's next neighbor is a core pixel. - // if it is a core pixel set pixel state to: used. - // - MGeomPix &gpix = (*fCam)[id]; - const Int_t nnmax = gpix.GetNumNeighbors(); - - for (Int_t j=0; jIsPixelCore(id2)) - continue; - - pix.SetPixelUsed(); - break; - } - } -} - -// -------------------------------------------------------------------------- -// -// check if MEvtHeader exists in the Parameter list already. -// if not create one and add them to the list -// -Bool_t MImgCleanStd::PreProcess (MParList *pList) -{ - fCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fCam) - { - *fLog << dbginf << "MGeomCam not found (no geometry information available)... aborting." << endl; - return kFALSE; - } - - fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fEvt) - { - *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - return kTRUE; -} - - -// -------------------------------------------------------------------------- -// -// Cleans the image. -// -Bool_t MImgCleanStd::Process() -{ - const Int_t max = CleanStep1(); - CleanStep2(max); - CleanStep3(); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Print descriptor and cleaning levels. -// -void MImgCleanStd::Print(Option_t *o) const -{ - *fLog << GetDescriptor() << " initialized with noise level "; - *fLog << fCleanLvl1 << " and " << fCleanLvl2 << endl; -} - -// -------------------------------------------------------------------------- -// -// Craete two text entry fields, one for each cleaning level and a -// describing text line. -// -void MImgCleanStd::CreateGuiElements(MGGroupFrame *f) -{ - // - // Create a frame for line 3 and 4 to be able - // to align entry field and label in one line - // - TGHorizontalFrame *f1 = new TGHorizontalFrame(f, 0, 0); - TGHorizontalFrame *f2 = new TGHorizontalFrame(f, 0, 0); - - /* - * --> use with root >=3.02 <-- - * - - TGNumberEntry *fNumEntry1 = new TGNumberEntry(frame, 3.0, 2, M_NENT_LVL1, kNESRealOne, kNEANonNegative); - TGNumberEntry *fNumEntry2 = new TGNumberEntry(frame, 2.5, 2, M_NENT_LVL1, kNESRealOne, kNEANonNegative); - - */ - TGTextEntry *entry1 = new TGTextEntry(f1, "****", kImgCleanLvl1); - TGTextEntry *entry2 = new TGTextEntry(f2, "****", kImgCleanLvl2); - - // --- doesn't work like expected (until root 3.02?) --- fNumEntry1->SetAlignment(kTextRight); - // --- doesn't work like expected (until root 3.02?) --- fNumEntry2->SetAlignment(kTextRight); - - entry1->SetText("3.0"); - entry2->SetText("2.5"); - - entry1->Associate(f); - entry2->Associate(f); - - TGLabel *l1 = new TGLabel(f1, "Cleaning Level 1"); - TGLabel *l2 = new TGLabel(f2, "Cleaning Level 2"); - - l1->SetTextJustify(kTextLeft); - l2->SetTextJustify(kTextLeft); - - // - // Align the text of the label centered, left in the row - // with a left padding of 10 - // - TGLayoutHints *laylabel = new TGLayoutHints(kLHintsCenterY|kLHintsLeft, 10); - TGLayoutHints *layframe = new TGLayoutHints(kLHintsCenterY|kLHintsLeft, 5, 0, 10); - - // - // Add one entry field and the corresponding label to each line - // - f1->AddFrame(entry1); - f2->AddFrame(entry2); - - f1->AddFrame(l1, laylabel); - f2->AddFrame(l2, laylabel); - - f->AddFrame(f1, layframe); - f->AddFrame(f2, layframe); - - f->AddToList(entry1); - f->AddToList(entry2); - f->AddToList(l1); - f->AddToList(l2); - f->AddToList(laylabel); - f->AddToList(layframe); -} - -void MImgCleanStd::SetLvl1(Float_t lvl) -{ - fCleanLvl1 = lvl; - *fLog << inf << "Cleaning level 1 set to " << lvl << " sigma." << endl; -} - -void MImgCleanStd::SetLvl2(Float_t lvl) -{ - fCleanLvl2 = lvl; - *fLog << inf << "Cleaning level 2 set to " << lvl << " sigma." << endl; -} - -// -------------------------------------------------------------------------- -// -// Process the GUI Events comming from the two text entry fields. -// -Bool_t MImgCleanStd::ProcessMessage(Int_t msg, Int_t submsg, Long_t param1, Long_t param2) -{ - if (msg!=kC_TEXTENTRY || submsg!=kTE_ENTER) - return kTRUE; - - TGTextEntry *txt = (TGTextEntry*)FindWidget(param1); - - if (!txt) - return kTRUE; - - Float_t lvl = atof(txt->GetText()); - - switch (param1) - { - case kImgCleanLvl1: - SetLvl1(lvl); - return kTRUE; - - case kImgCleanLvl2: - SetLvl2(lvl); - return kTRUE; - } - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Implementation of SavePrimitive. Used to write the call to a constructor -// to a macro. In the original root implementation it is used to write -// gui elements to a macro-file. -// -void MImgCleanStd::StreamPrimitive(ofstream &out) const -{ - out << " MImgCleanStd " << GetUniqueName() << "("; - out << fCleanLvl1 << ", " << fCleanLvl2; - - if (fName!=gsDefName || fTitle!=gsDefTitle) - { - out << ", \"" << fName << "\""; - if (fTitle!=gsDefTitle) - out << ", \"" << fTitle << "\""; - } - out << ");" << endl; -} - -// -------------------------------------------------------------------------- -// -// Read the setup from a TEnv: -// Float_t fCleanLvl1: CleaningLevel1 -// Float_t fCleanLvl2: CleaningLevel2 -// -Bool_t MImgCleanStd::ReadEnv(const TEnv &env, TString prefix, Bool_t print) -{ - Bool_t rc = kTRUE; - if (!IsEnvDefined(env, prefix, "CleaningLevel1", print)) - rc = kFALSE; - else - { - SetLvl1(GetEnvValue(env, prefix, "CleaningLevel1", fCleanLvl1)); - if (fCleanLvl1<0) - { - *fLog << err << "ERROR - Negative values for Cleaning Level 1 forbidden." << endl; - return kERROR; - } - } - - if (!IsEnvDefined(env, prefix, "CleaningLevel2", print)) - rc = kFALSE; - else - { - SetLvl2(GetEnvValue(env, prefix, "CleaningLevel2", fCleanLvl2)); - if (fCleanLvl2<0) - { - *fLog << err << "ERROR - Negative values for Cleaning Level 2 forbidden." << endl; - return kERROR; - } - } - - return rc; -} Index: trunk/MagicSoft/Mars/manalysis/MImgCleanStd.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MImgCleanStd.h (revision 1939) +++ (revision ) @@ -1,51 +1,0 @@ -#ifndef MARS_MImgCleanStd -#define MARS_MImgCleanStd - -#ifndef MARS_MGTask -#include "MGTask.h" -#endif - -class MGeomCam; -class MCerPhotEvt; - -class MGGroupFrame; - -class MImgCleanStd : public MGTask -{ -private: - const MGeomCam *fCam; //! - MCerPhotEvt *fEvt; //! - - Float_t fCleanLvl1; - Float_t fCleanLvl2; - - void SetLvl1(Float_t lvl); - void SetLvl2(Float_t lvl); - - void CreateGuiElements(MGGroupFrame *f); - void StreamPrimitive(ofstream &out) const; - -public: - MImgCleanStd(const Float_t lvl1=3.0, const Float_t lvl2=2.5, - const char *name=NULL, const char *title=NULL); - - Int_t CleanStep1(); - void CleanStep2(Int_t max); - void CleanStep3(); - - Bool_t PreProcess(MParList *pList); - Bool_t Process(); - - void Print(Option_t *o="") const; - - Float_t GetCleanLvl1() const { return fCleanLvl1; } - Float_t GetCleanLvl2() const { return fCleanLvl2; } - - Bool_t ProcessMessage(Int_t msg, Int_t submsg, Long_t param1, Long_t param2); - Bool_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); - - ClassDef(MImgCleanStd, 0) // task doing a standard image cleaning -}; - -#endif - Index: trunk/MagicSoft/Mars/manalysis/MNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MNewImagePar.cc (revision 1939) +++ (revision ) @@ -1,114 +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 03/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MNewImagePar -// -// Storage Container for new image parameters -// -///////////////////////////////////////////////////////////////////////////// -#include "MNewImagePar.h" - -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MGeomCam.h" -#include "MGeomPix.h" -#include "MCerPhotEvt.h" -#include "MCerPhotPix.h" -#include "MSrcPosCam.h" - -ClassImp(MNewImagePar); - -// -------------------------------------------------------------------------- -// -// Default constructor. -// -MNewImagePar::MNewImagePar(const char *name, const char *title) -{ - fName = name ? name : "MNewImagePar"; - fTitle = title ? title : "New image parameters"; -} - -// -------------------------------------------------------------------------- -// -void MNewImagePar::Reset() -{ - fLeakage1 = 0; - fLeakage2 = 0; -} - -// -------------------------------------------------------------------------- -// -// Calculation of new image parameters -// -// -Bool_t MNewImagePar::Calc(const MGeomCam &geom, const MCerPhotEvt &evt, - const MHillas *hillas) -{ - //fHillas = (MHillas*)hillas; - - const UInt_t npixevt = evt.GetNumPixels(); - - Double_t edgepix1 = 0.0; - Double_t edgepix2 = 0.0; - - for (UInt_t i=0; iGetSize(); - fLeakage2 = edgepix2 / hillas->GetSize(); - - SetReadyToSave(); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -void MNewImagePar::Print(Option_t *) const -{ - *fLog << all; - *fLog << "New Image Parameters (" << GetName() << ")" << endl; - *fLog << " - Leakage1 = " << fLeakage1 << endl; - *fLog << " - Leakage2 = " << fLeakage2 << endl; -} Index: trunk/MagicSoft/Mars/manalysis/MNewImagePar.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MNewImagePar.h (revision 1939) +++ (revision ) @@ -1,65 +1,0 @@ -#ifndef MARS_MNewImagePar -#define MARS_MNewImagePar - -#ifndef MARS_MHillas -#include "MHillas.h" -#endif - -class MSrcPosCam; - -class MNewImagePar : public MParContainer -{ -private: - Float_t fLeakage1; // (photons in most outer ring of pixels) over fSize - Float_t fLeakage2; // (photons in the 2 outer rings of pixels) over fSize - -public: - MNewImagePar(const char *name=NULL, const char *title=NULL); - - // void SetSrcPos(MSrcPosCam *pos) { fSrcPos = pos; } - // const MSrcPosCam *GetSrcPos() const { return fSrcPos; } - - void Reset(); - - Float_t GetLeakage1() const { return fLeakage1; } - Float_t GetLeakage2() const { return fLeakage2; } - - void Print(Option_t *opt=NULL) const; - - virtual Bool_t Calc(const MGeomCam &geom, const MCerPhotEvt &evt, - const MHillas *hillas); - - //virtual void AsciiRead(ifstream &fin); - //virtual void AsciiWrite(ofstream &fout) const; - - ClassDef(MNewImagePar, 1) // Container to hold new image parameters -}; - -#endif - - - - - - - - - - - - - - - - - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MNewImageParCalc.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MNewImageParCalc.cc (revision 1939) +++ (revision ) @@ -1,159 +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 03/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -////////////////////////////////////////////////////////////////////////////// -// -// MNewImageParCalc -// -// Task to calculate the source dependant part of the hillas parameters -// -////////////////////////////////////////////////////////////////////////////// -#include "MNewImageParCalc.h" - -#include - -#include "MParList.h" - -#include "MGeomCam.h" -#include "MSrcPosCam.h" -#include "MCerPhotEvt.h" -#include "MNewImagePar.h" -#include "MNewImagePar.h" -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MNewImageParCalc); - -static const TString gsDefName = "MNewImageParCalc"; -static const TString gsDefTitle = "Calculate new image parameters"; - -// ------------------------------------------------------------------------- -// -// Default constructor. The first argument is the name of a container -// containing the source position in the camera plain (MScrPosCam). -// The default is "MSrcPosCam". newpar is the name of a container -// of type MNewImagePar, in which the parameters are stored. -// The default is "MNewImagePar" -// -// -MNewImageParCalc::MNewImageParCalc(const char *src, const char *newpar, - const char *name, const char *title) - : fHillas(NULL), fSrcPos(NULL), fNewImagePar(NULL) -{ - fName = name ? name : gsDefName.Data(); - fTitle = title ? title : gsDefTitle.Data(); - - fSrcName = src; - fNewParName = newpar; - fHillasInput = "MHillas"; -} - -// ------------------------------------------------------------------------- -// -Bool_t MNewImageParCalc::PreProcess(MParList *pList) -{ - fHillas = (MHillas*)pList->FindObject(fHillasInput, "MHillas"); - if (!fHillas) - { - *fLog << err << dbginf << "MHillas not found... aborting." << endl; - return kFALSE; - } - - fSrcPos = (MSrcPosCam*)pList->FindObject(fSrcName, "MSrcPosCam"); - if (!fSrcPos) - { - *fLog << err << dbginf << fSrcName << " [MSrcPosCam] not found... aborting." << endl; - return kFALSE; - } - - fCerPhotEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); - if (!fCerPhotEvt) - { - *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; - return kFALSE; - } - - fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fGeomCam) - { - *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; - return kFALSE; - } - - - fNewImagePar = (MNewImagePar*)pList->FindCreateObj("MNewImagePar", fNewParName); - if (!fNewImagePar) - return kFALSE; - - fErrors = 0; - - return kTRUE; -} - -// ------------------------------------------------------------------------- -// -Bool_t MNewImageParCalc::Process() -{ - - if (!fNewImagePar->Calc(*fGeomCam, *fCerPhotEvt, fHillas)) - { - fErrors++; - return kCONTINUE; - - } - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Prints some statistics about the hillas calculation. The percentage -// is calculated with respect to the number of executions of this task. -// -Bool_t MNewImageParCalc::PostProcess() -{ - if (GetNumExecutions()==0) - return kTRUE; - - *fLog << inf << endl; - *fLog << GetDescriptor() << " execution statistics:" << endl; - *fLog << dec << setfill(' '); - *fLog << " " << fErrors << " (" << (int)(fErrors*100/GetNumExecutions()) << "%) Evts skipped due to: calculation failed" << endl; - *fLog << endl; - - return kTRUE; -} - -// -------------------------------------------------------------------------- - - - - - - - - - - - - Index: trunk/MagicSoft/Mars/manalysis/MNewImageParCalc.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MNewImageParCalc.h (revision 1939) +++ (revision ) @@ -1,47 +1,0 @@ -#ifndef MARS_MNewImageParCalc -#define MARS_MNewImageParCalc - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -class MHillas; -class MNewImagePar; -class MSrcPosCam; -class MGeomCam; -class MCerPhotEvt; - -class MNewImageParCalc : public MTask -{ -private: - MGeomCam *fGeomCam; - MCerPhotEvt *fCerPhotEvt; - - MHillas *fHillas; //! Pointer to the source independent hillas parameters - MSrcPosCam *fSrcPos; //! Pointer to the source position - MNewImagePar *fNewImagePar; //! Pointer to the output container for the new image parameters - - TString fSrcName; - TString fNewParName; - TString fHillasInput; - - Int_t fErrors; - - Bool_t PreProcess(MParList *plist); - Bool_t Process(); - Bool_t PostProcess(); - -public: - MNewImageParCalc(const char *src="MSrcPosCam", const char *newpar="MNewImagePar", - const char *name=NULL, const char *title=NULL); - - void SetInput(TString hilname) { fHillasInput = hilname; } - - ClassDef(MNewImageParCalc, 1) // task to calculate new image parameters -}; - -#endif - - - - Index: trunk/MagicSoft/Mars/manalysis/Makefile =================================================================== --- trunk/MagicSoft/Mars/manalysis/Makefile (revision 1939) +++ trunk/MagicSoft/Mars/manalysis/Makefile (revision 1940) @@ -23,5 +23,5 @@ # INCLUDES = -I. -I../mbase -I../mmc -I../mraw -I../mgeom \ - -I../mdata -I../mhist -I../mgui + -I../mdata -I../mhist -I../mgui -I../mimage #------------------------------------------------------------------------------ @@ -33,10 +33,8 @@ MMcPedestalCopy.cc \ MMcPedestalNSBAdd.cc \ - MImgCleanStd.cc \ MEnergyEst.cc \ MEnergyEstimate.cc \ MEnergyEstParam.cc \ MSrcPosCam.cc \ - MCameraSmooth.cc \ MHadronness.cc \ MMatrixLoop.cc \ @@ -48,9 +46,4 @@ MRanForestCalc.cc \ MRanForestFill.cc \ - MHillas.cc \ - MHillasSrc.cc \ - MHillasExt.cc \ - MHillasCalc.cc \ - MHillasSrcCalc.cc \ MCerPhotPix.cc \ MCerPhotEvt.cc \ @@ -66,11 +59,8 @@ MPadSchweizer.cc \ MPointingCorr.cc \ - MNewImagePar.cc \ - MNewImageParCalc.cc \ MParameters.cc \ MMcTriggerLvl2.cc \ MMcTriggerLvl2Calc.cc \ MCT1SupercutsCalc.cc - SRCS = $(SRCFILES) Index: trunk/MagicSoft/Mars/mhist/HistLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mhist/HistLinkDef.h (revision 1939) +++ trunk/MagicSoft/Mars/mhist/HistLinkDef.h (revision 1940) @@ -20,7 +20,4 @@ #pragma link C++ class MHCerPhotEvt+; -#pragma link C++ class MHHillas+; -#pragma link C++ class MHHillasSrc+; -#pragma link C++ class MHHillasExt+; #pragma link C++ class MHStarMap+; #pragma link C++ class MHEnergyTime+; @@ -59,17 +56,3 @@ #pragma link C++ class MHMcTriggerLvl2+; -#pragma link C++ class MHNewImagePar+; - #endif - - - - - - - - - - - - Index: trunk/MagicSoft/Mars/mhist/MHHillas.cc =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHillas.cc (revision 1939) +++ (revision ) @@ -1,388 +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): Thomas Bretz 2001 -! Wolfgang Wittek 2002 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MHHillas -// -// This class contains histograms for the source independent image parameters -// -///////////////////////////////////////////////////////////////////////////// - -#include "MHHillas.h" - -#include - -#include -#include -#include -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MGeomCam.h" -#include "MBinning.h" - - -ClassImp(MHHillas); - -// -------------------------------------------------------------------------- -// -// Setup four histograms for Width, Length -// -MHHillas::MHHillas(const char *name, const char *title) - : fMm2Deg(1), fUseMmScale(kTRUE) -{ - // - // set the name and title of this object - // - fName = name ? name : "MHHillas"; - fTitle = title ? title : "Source independent image parameters"; - - fLength = new TH1F("Length", "Length of Ellipse", 100, 0, 296.7); - fWidth = new TH1F("Width", "Width of Ellipse", 100, 0, 296.7); - fDistC = new TH1F("DistC", "Distance from center of camera", 100, 0, 445); - fDelta = new TH1F("Delta", "Angle (Main axis - x-axis)", 101, -90, 90); - fUsedPix = new TH1F("UsedPix", "Number of used pixels", 150, 0, 150); - fCorePix = new TH1F("CorePix", "Number of core pixels", 150, 0, 150); - - fLength->SetLineColor(kBlue); - fCorePix->SetLineColor(kRed); - fUsedPix->SetLineColor(kGreen); - - fLength->SetDirectory(NULL); - fWidth->SetDirectory(NULL); - fDistC->SetDirectory(NULL); - fDelta->SetDirectory(NULL); - fUsedPix->SetDirectory(NULL); - fCorePix->SetDirectory(NULL); - - fLength->SetXTitle("Length [mm]"); - fWidth->SetXTitle("Width [mm]"); - fDistC->SetXTitle("Distance [mm]"); - fDelta->SetXTitle("Delta [\\circ]"); - fUsedPix->SetXTitle("Number of Pixels"); - fCorePix->SetXTitle("Number of Pixels"); - - fLength->SetYTitle("Counts"); - fWidth->SetYTitle("Counts"); - fDistC->SetYTitle("Counts"); - fDelta->SetYTitle("Counts"); - fUsedPix->SetYTitle("Counts"); - fCorePix->SetYTitle("Counts"); - - MBinning bins; - bins.SetEdgesLog(50, 1, 1e7); - - fSize = new TH1F; - fSize->SetName("Size"); - fSize->SetTitle("Number of Photons"); - fSize->SetDirectory(NULL); - fSize->SetXTitle("Size"); - fSize->SetYTitle("Counts"); - fSize->GetXaxis()->SetTitleOffset(1.2); - fSize->GetXaxis()->SetLabelOffset(-0.015); - - bins.Apply(*fSize); - - fCenter = new TH2F("Center", "Center of Ellipse", 51, -445, 445, 51, -445, 445); - fCenter->SetDirectory(NULL); - fCenter->SetXTitle("x [mm]"); - fCenter->SetYTitle("y [mm]"); - fCenter->SetZTitle("Counts"); -} - -// -------------------------------------------------------------------------- -// -// Delete the histograms -// -MHHillas::~MHHillas() -{ - delete fLength; - delete fWidth; - - delete fDistC; - delete fDelta; - - delete fSize; - delete fCenter; - - delete fUsedPix; - delete fCorePix; -} - -// -------------------------------------------------------------------------- -// -// Setup the Binning for the histograms automatically if the correct -// instances of MBinning (with the names 'BinningWidth' and 'BinningLength') -// are found in the parameter list -// Use this function if you want to set the conversion factor which -// is used to convert the mm-scale in the camera plain into the deg-scale -// used for histogram presentations. The conversion factor is part of -// the camera geometry. Please create a corresponding MGeomCam container. -// -Bool_t MHHillas::SetupFill(const MParList *plist) -{ - const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); - if (!geom) - *fLog << warn << GetDescriptor() << ": No Camera Geometry available. Using mm-scale for histograms." << endl; - else - { - fMm2Deg = geom->GetConvMm2Deg(); - SetMmScale(kFALSE); - } - - ApplyBinning(*plist, "Width", fWidth); - ApplyBinning(*plist, "Length", fLength); - ApplyBinning(*plist, "Dist", fDistC); - ApplyBinning(*plist, "Delta", fDelta); - ApplyBinning(*plist, "Size", fSize); - ApplyBinning(*plist, "Pixels", fUsedPix); - ApplyBinning(*plist, "Pixels", fCorePix); - - const MBinning *bins = (MBinning*)plist->FindObject("BinningCamera"); - if (!bins) - { - float r = geom ? geom->GetMaxRadius() : 600; - r *= 0.9; - if (!fUseMmScale) - r *= fMm2Deg; - - MBinning b; - b.SetEdges(61, -r, r); - SetBinning(fCenter, &b, &b); - } - else - SetBinning(fCenter, bins, bins); - - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Use this function to setup your own conversion factor between degrees -// and millimeters. The conversion factor should be the one calculated in -// MGeomCam. Use this function with Caution: You could create wrong values -// by setting up your own scale factor. -// -void MHHillas::SetMm2Deg(Float_t mmdeg) -{ - if (mmdeg<0) - { - *fLog << warn << dbginf << "Warning - Conversion factor < 0 - nonsense. Ignored." << endl; - return; - } - - if (fMm2Deg>=0) - *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; - - fMm2Deg = mmdeg; -} - -// -------------------------------------------------------------------------- -// -// With this function you can convert the histogram ('on the fly') between -// degrees and millimeters. -// -void MHHillas::SetMmScale(Bool_t mmscale) -{ - if (fUseMmScale == mmscale) - return; - - if (fMm2Deg<0) - { - *fLog << warn << dbginf << "Warning - Sorry, no conversion factor for conversion available." << endl; - return; - } - - const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; - MH::ScaleAxis(fLength, scale); - MH::ScaleAxis(fWidth, scale); - MH::ScaleAxis(fDistC, scale); - MH::ScaleAxis(fCenter, scale, scale); - - if (mmscale) - { - fLength->SetXTitle("Length [mm]"); - fWidth->SetXTitle("Width [mm]"); - fDistC->SetXTitle("Distance [mm]"); - fCenter->SetXTitle("x [mm]"); - fCenter->SetYTitle("y [mm]"); - } - else - { - fLength->SetXTitle("Length [\\circ]"); - fWidth->SetXTitle("Width [\\circ]"); - fDistC->SetXTitle("Distance [\\circ]"); - fCenter->SetXTitle("x [\\circ]"); - fCenter->SetYTitle("y [\\circ]"); - } - - fUseMmScale = mmscale; -} - -// -------------------------------------------------------------------------- -// -// Fill the histograms with data from a MHillas-Container. -// Be careful: Only call this with an object of type MHillas -// -Bool_t MHHillas::Fill(const MParContainer *par) -{ - const MHillas &h = *(MHillas*)par; - - const Double_t d = sqrt(h.GetMeanX()*h.GetMeanX() + h.GetMeanY()*h.GetMeanY()); - const Double_t scale = fUseMmScale ? 1 : fMm2Deg; - - fLength ->Fill(scale*h.GetLength()); - fWidth ->Fill(scale*h.GetWidth()); - fDistC ->Fill(scale*d); - fCenter ->Fill(scale*h.GetMeanX(), scale*h.GetMeanY()); - fDelta ->Fill(kRad2Deg*h.GetDelta()); - fSize ->Fill(h.GetSize()); - fUsedPix->Fill(h.GetNumUsedPixels()); - fCorePix->Fill(h.GetNumCorePixels()); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Setup a inversed deep blue sea palette for the fCenter histogram. -// -void MHHillas::SetColors() const -{ - // FIXME: This must be redone each time the canvas is repainted.... - gStyle->SetPalette(51, NULL); - Int_t c[50]; - for (int i=0; i<50; i++) - c[49-i] = gStyle->GetColorPalette(i); - gStyle->SetPalette(50, c); -} - -// -------------------------------------------------------------------------- -// -// Draw clones of four histograms. So that the object can be deleted -// and the histograms are still visible in the canvas. -// The cloned object are deleted together with the canvas if the canvas is -// destroyed. If you want to handle dostroying the canvas you can get a -// pointer to it from this function -// -TObject *MHHillas::DrawClone(Option_t *opt) const -{ - TCanvas *c = MakeDefCanvas(this, 720, 810); - c->Divide(2,3); - - gROOT->SetSelectedPad(NULL); - - c->cd(1); - DrawCopy(*fWidth, *fLength, "Width / Length"); - - c->cd(2); - gPad->SetLogx(); - fSize->DrawCopy(); - - c->cd(3); - DrawCopy(*fCorePix, *fUsedPix, "Number of core/used Pixels"); - - c->cd(4); - fDelta->DrawCopy(); - - c->cd(5); - fDistC->DrawCopy(); - - c->cd(6); - SetColors(); - fCenter->DrawCopy("colz"); - - c->Modified(); - c->Update(); - - return c; -} - -// -------------------------------------------------------------------------- -// -// Creates a new canvas and draws the four histograms into it. -// Be careful: The histograms belongs to this object and won't get deleted -// together with the canvas. -// -void MHHillas::Draw(Option_t *) -{ - if (!gPad) - MakeDefCanvas(this, 720, 810); - - gPad->Divide(2,3); - - gPad->cd(1); - MH::Draw(*fWidth, *fLength, "Width / Length"); - - gPad->cd(2); - gPad->SetLogx(); - fSize->Draw(); - - gPad->cd(3); - MH::Draw(*fCorePix, *fUsedPix, "Number of core/used Pixels"); - - gPad->cd(4); - fDelta->Draw(); - - gPad->cd(5); - fDistC->Draw(); - - gPad->cd(6); - SetColors(); - fCenter->Draw("colz"); - - gPad->Modified(); - gPad->Update(); -} - -TH1 *MHHillas::GetHistByName(const TString name) -{ - if (name.Contains("Width", TString::kIgnoreCase)) - return fWidth; - if (name.Contains("Length", TString::kIgnoreCase)) - return fLength; - if (name.Contains("Size", TString::kIgnoreCase)) - return fSize; - if (name.Contains("Core", TString::kIgnoreCase)) - return fCorePix; - if (name.Contains("Used", TString::kIgnoreCase)) - return fUsedPix; - if (name.Contains("Delta", TString::kIgnoreCase)) - return fDelta; - if (name.Contains("DistC", TString::kIgnoreCase)) - return fDistC; - if (name.Contains("Center", TString::kIgnoreCase)) - return fCenter; - - return NULL; -} Index: trunk/MagicSoft/Mars/mhist/MHHillas.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHillas.h (revision 1939) +++ (revision ) @@ -1,60 +1,0 @@ -#ifndef MARS_MHHillas -#define MARS_MHHillas - -#ifndef MARS_MH -#include "MH.h" -#endif - -class TH1F; -class TH2F; -class MHillas; - -class MHHillas : public MH -{ -private: - - TH1F *fLength; //-> - TH1F *fWidth; //-> - - TH1F *fDistC; //-> - TH1F *fDelta; //-> - - TH1F *fSize; //-> - TH2F *fCenter; //-> - - TH1F *fUsedPix; //-> - TH1F *fCorePix; //-> - - void SetColors() const; - - Float_t fMm2Deg; - Bool_t fUseMmScale; - -public: - MHHillas(const char *name=NULL, const char *title=NULL); - ~MHHillas(); - - void SetMmScale(Bool_t mmscale=kTRUE); - virtual void SetMm2Deg(Float_t mmdeg); - - Bool_t SetupFill(const MParList *pList); - Bool_t Fill(const MParContainer *par); - - TH1 *GetHistByName(const TString name); - - TH1F *GetHistLength() { return fLength; } - TH1F *GetHistWidth() { return fWidth; } - - TH1F *GetHistDistC() { return fDistC; } - TH1F *GetHistDelta() { return fDelta; } - - TH1F *GetHistSize() { return fSize; } - TH2F *GetHistCenter() { return fCenter; } - - void Draw(Option_t *opt=NULL); - TObject *DrawClone(Option_t *opt=NULL) const; - - ClassDef(MHHillas, 1) // Container which holds histograms for the source independent image parameters -}; - -#endif Index: trunk/MagicSoft/Mars/mhist/MHHillasExt.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHillasExt.h (revision 1939) +++ (revision ) @@ -1,50 +1,0 @@ -#ifndef MARS_MHHillasExt -#define MARS_MHHillasExt - -#ifndef ROOT_TH1 -#include -#endif -#ifndef MARS_MH -#include "MH.h" -#endif - -class MHillasExt; - -class MHHillasExt : public MH -{ -private: - MHillasExt *fHillasExt; //! Pointer to the MHillasExt container - - TH1F fHConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize - TH1F fHConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize - TH1F fHAsym; // [mm] fDist minus dist: center of ellipse, highest pixel - TH1F fHM3Long; // [mm] 3rd moment (e-weighted) along major axis - TH1F fHM3Trans; // [mm] 3rd moment (e-weighted) along minor axis - - Float_t fMm2Deg; - Bool_t fUseMmScale; - - TString fHilName; - -public: - //MHHillasExt(const char *name=NULL, const char *title=NULL); - MHHillasExt(const char *name=NULL, const char *title=NULL, - const char *hil="MHillas"); - - ~MHHillasExt(); - - void SetMmScale(Bool_t mmscale=kTRUE); - virtual void SetMm2Deg(Float_t mmdeg); - - Bool_t SetupFill(const MParList *pList); - Bool_t Fill(const MParContainer *par); - - TH1 *GetHistByName(const TString name); - - void Draw(Option_t *opt=NULL); - TObject *DrawClone(Option_t *opt=NULL) const; - - ClassDef(MHHillasExt, 1) // Container which holds histograms for the extended hillas parameters -}; - -#endif Index: trunk/MagicSoft/Mars/mhist/MHHillasSrc.cc =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHillasSrc.cc (revision 1939) +++ (revision ) @@ -1,286 +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): Thomas Bretz 2001 -! -! Copyright: MAGIC Software Development, 2000-2001 -! -! -\* ======================================================================== */ - -/////////////////////////////////////////////////////////////////////// -// -// MHHillasSrc -// -// This class contains histograms for every Hillas parameter -// -/////////////////////////////////////////////////////////////////////// -#include "MHHillasSrc.h" - -#include - -#include -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MGeomCam.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MHillasSrc.h" - -ClassImp(MHHillasSrc); - -// -------------------------------------------------------------------------- -// -// Setup four histograms for Alpha, and Dist -// -MHHillasSrc::MHHillasSrc(const char *name, const char *title) - : fUseMmScale(kTRUE) -{ - // - // set the name and title of this object - // - fName = name ? name : "MHHillasSrc"; - fTitle = title ? title : "Container for Hillas histograms"; - - // - // loop over all Pixels and create two histograms - // one for the Low and one for the High gain - // connect all the histogram with the container fHist - // - fAlpha = new TH1F("Alpha", "Alpha of Ellipse", 181, -90, 90); - fDist = new TH1F("Dist", "Dist of Ellipse", 100, 0, 445); - fHeadTail = new TH1F("HeadTail", "HeadTail of Ellipse", 101, -445, 445); - fCosDA = new TH1F("CosDA", "cos(Delta,Alpha) of Ellipse", 101, -1, 1); - - fAlpha->SetDirectory(NULL); - fDist->SetDirectory(NULL); - fHeadTail->SetDirectory(NULL); - fCosDA->SetDirectory(NULL); - - fAlpha->SetXTitle("\\alpha [\\circ]"); - fDist->SetXTitle("Dist [mm]"); - fHeadTail->SetXTitle("Head-Tail [mm]"); - fCosDA->SetXTitle("cos(\\delta,\\alpha)"); - - fAlpha->SetYTitle("Counts"); - fDist->SetYTitle("Counts"); - fHeadTail->SetYTitle("Counts"); - fCosDA->SetYTitle("Counts"); -} - -// -------------------------------------------------------------------------- -// -// Delete the four histograms -// -MHHillasSrc::~MHHillasSrc() -{ - delete fAlpha; - delete fDist; - delete fHeadTail; - delete fCosDA; -} - -// -------------------------------------------------------------------------- -// -// Setup the Binning for the histograms automatically if the correct -// instances of MBinning (with the names 'BinningAlpha' and 'BinningDist') -// are found in the parameter list -// Use this function if you want to set the conversion factor which -// is used to convert the mm-scale in the camera plain into the deg-scale -// used for histogram presentations. The conversion factor is part of -// the camera geometry. Please create a corresponding MGeomCam container. -// -Bool_t MHHillasSrc::SetupFill(const MParList *plist) -{ - const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); - if (!geom) - *fLog << warn << dbginf << "No Camera Geometry available. Using mm-scale for histograms." << endl; - else - { - fMm2Deg = geom->GetConvMm2Deg(); - SetMmScale(kFALSE); - } - - ApplyBinning(*plist, "Alpha", fAlpha); - ApplyBinning(*plist, "Dist", fDist); - ApplyBinning(*plist, "HeadTail", fHeadTail); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Fill the four histograms with data from a MHillas-Container. -// Be careful: Only call this with an object of type MHillas -// -Bool_t MHHillasSrc::Fill(const MParContainer *par) -{ - const MHillasSrc &h = *(MHillasSrc*)par; - - fAlpha ->Fill(h.GetAlpha()); - fDist ->Fill(fUseMmScale ? h.GetDist() : fMm2Deg*h.GetDist()); - fHeadTail->Fill(fUseMmScale ? h.GetHeadTail() : fMm2Deg*h.GetHeadTail()); - fCosDA ->Fill(h.GetCosDeltaAlpha()); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Use this function to setup your own conversion factor between degrees -// and millimeters. The conversion factor should be the one calculated in -// MGeomCam. Use this function with Caution: You could create wrong values -// by setting up your own scale factor. -// -void MHHillasSrc::SetMm2Deg(Float_t mmdeg) -{ - if (mmdeg<=0) - { - *fLog << warn << dbginf << "Warning - Conversion factor <= 0 - nonsense. Ignored." << endl; - return; - } - - if (fMm2Deg>0) - *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; - - fMm2Deg = mmdeg; -} - -// -------------------------------------------------------------------------- -// -// With this function you can convert the histogram ('on the fly') between -// degrees and millimeters. -// -void MHHillasSrc::SetMmScale(Bool_t mmscale) -{ - if (fUseMmScale == mmscale) - return; - - if (fMm2Deg<0) - { - *fLog << warn << GetDescriptor() << ": Warning - Sorry, no conversion factor for conversion available." << endl; - return; - } - - const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; - MH::ScaleAxis(fDist, scale); - MH::ScaleAxis(fHeadTail, scale); - - if (mmscale) - { - fDist->SetXTitle("Dist [mm]"); - fHeadTail->SetXTitle("Head-Tail [mm]"); - } - else - { - fDist->SetXTitle("Dist [\\circ]"); - fHeadTail->SetXTitle("Head-Tail [\\circ]"); - } - - fUseMmScale = mmscale; -} - -// -------------------------------------------------------------------------- -// -// Draw clones of all two histograms. So that the object can be deleted -// and the histograms are still visible in the canvas. -// The cloned object are deleted together with the canvas if the canvas is -// destroyed. If you want to handle dostroying the canvas you can get a -// pointer to it from this function -// -TObject *MHHillasSrc::DrawClone(Option_t *opt) const -{ - TCanvas *c = MakeDefCanvas(this, 700, 500); - c->Divide(2, 2); - - // FIXME: Display Source position - - gROOT->SetSelectedPad(NULL); - - // - // This is necessary to get the expected bahviour of DrawClone - // - c->cd(1); - fAlpha->DrawCopy(); - - c->cd(2); - fDist->DrawCopy(); - - c->cd(3); - fHeadTail->DrawCopy(); - - c->cd(4); - gPad->SetLogy(); - fCosDA->DrawCopy(); - - c->Modified(); - c->Update(); - - return c; -} - -// -------------------------------------------------------------------------- -// -// Creates a new canvas and draws the two histograms into it. -// Be careful: The histograms belongs to this object and won't get deleted -// together with the canvas. -// -void MHHillasSrc::Draw(Option_t *) -{ - if (!gPad) - MakeDefCanvas(this, 700, 500); - - // FIXME: Display Source position - - gPad->Divide(2, 2); - - gPad->cd(1); - fAlpha->Draw(); - - gPad->cd(2); - fDist->Draw(); - - gPad->cd(1); - fHeadTail->Draw(); - - gPad->cd(2); - gPad->SetLogy(); - fCosDA->Draw(); - - gPad->Modified(); - gPad->Update(); -} - -TH1 *MHHillasSrc::GetHistByName(const TString name) -{ - if (name.Contains("Alpha", TString::kIgnoreCase)) - return fAlpha; - if (name.Contains("Dist", TString::kIgnoreCase)) - return fDist; - if (name.Contains("HeadTail", TString::kIgnoreCase)) - return fHeadTail; - if (name.Contains("CosDA", TString::kIgnoreCase)) - return fCosDA; - - return NULL; -} Index: trunk/MagicSoft/Mars/mhist/MHHillasSrc.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHHillasSrc.h (revision 1939) +++ (revision ) @@ -1,45 +1,0 @@ -#ifndef MARS_MHHillasSrc -#define MARS_MHHillasSrc - -#ifndef MARS_MH -#include "MH.h" -#endif - -class TH1F; -class MHillas; - -class MHHillasSrc : public MH -{ -private: - TH1F *fAlpha; //-> - TH1F *fDist; //-> - TH1F *fHeadTail; //-> - TH1F *fCosDA; //-> - - Float_t fMm2Deg; - Bool_t fUseMmScale; - -public: - MHHillasSrc(const char *name=NULL, const char *title=NULL); - ~MHHillasSrc(); - - void SetMmScale(Bool_t mmscale=kTRUE); - void SetMm2Deg(Float_t mmdeg); - - Bool_t SetupFill(const MParList *pList); - Bool_t Fill(const MParContainer *par); - - TH1 *GetHistByName(const TString name); - - TH1F *GetHistAlpha() { return fAlpha; } - TH1F *GetHistDist() { return fDist; } - TH1F *GetHistHeadTail() { return fHeadTail; } - TH1F *GetHistCosDeltaAlpha() { return fCosDA; } - - void Draw(Option_t *opt=NULL); - TObject *DrawClone(Option_t *opt=NULL) const; - - ClassDef(MHHillasSrc, 1) // Container which holds histograms for the source dependant parameters -}; - -#endif Index: trunk/MagicSoft/Mars/mhist/MHNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/mhist/MHNewImagePar.cc (revision 1939) +++ (revision ) @@ -1,256 +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 03/2003 -! -! Copyright: MAGIC Software Development, 2000-2003 -! -! -\* ======================================================================== */ - -/////////////////////////////////////////////////////////////////////// -// -// MHNewImagePar -// -// This class contains histograms for every Hillas parameter -// -/////////////////////////////////////////////////////////////////////// -#include "MHNewImagePar.h" - -#include - -#include -#include -#include - -#include "MLog.h" -#include "MLogManip.h" - -#include "MGeomCam.h" - -#include "MParList.h" - -#include "MHillas.h" -#include "MNewImagePar.h" - -ClassImp(MHNewImagePar); - -// -------------------------------------------------------------------------- -// -// Setup histograms -// -MHNewImagePar::MHNewImagePar(const char *name, const char *title) - : fUseMmScale(kTRUE) -{ - fName = name ? name : "MHNewImagePar"; - fTitle = title ? title : "Container for histograms of new image parameters"; - - fLeakage1 = new TH1F("Leakage1", "Leakage1", 100, 0.0, 1.0); - fLeakage1->SetDirectory(NULL); - fLeakage1->SetXTitle("Leakage1"); - fLeakage1->SetYTitle("Counts"); - - fLeakage2 = new TH1F("Leakage2", "Leakage2", 100, 0.0, 1.0); - fLeakage2->SetDirectory(NULL); - fLeakage2->SetXTitle("Leakage2"); - fLeakage2->SetYTitle("Counts"); -} - -// -------------------------------------------------------------------------- -// -// Delete the four histograms -// -MHNewImagePar::~MHNewImagePar() -{ - delete fLeakage1; - delete fLeakage2; -} - -// -------------------------------------------------------------------------- -// -// Setup the Binning for the histograms automatically if the correct -// instances of MBinning (with the names 'BinningAlpha' and 'BinningDist') -// are found in the parameter list -// Use this function if you want to set the conversion factor which -// is used to convert the mm-scale in the camera plain into the deg-scale -// used for histogram presentations. The conversion factor is part of -// the camera geometry. Please create a corresponding MGeomCam container. -// -Bool_t MHNewImagePar::SetupFill(const MParList *plist) -{ - const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); - if (!geom) - *fLog << warn << dbginf << "No Camera Geometry available. Using mm-scale for histograms." << endl; - else - { - fMm2Deg = geom->GetConvMm2Deg(); - SetMmScale(kFALSE); - } - - //ApplyBinning(*plist, "Alpha", fAlpha); - //ApplyBinning(*plist, "Dist", fDist); - //ApplyBinning(*plist, "HeadTail", fHeadTail); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Fill the histograms with data from a MNewImagePar container. -// -Bool_t MHNewImagePar::Fill(const MParContainer *par) -{ - const MNewImagePar &h = *(MNewImagePar*)par; - - fLeakage1->Fill(h.GetLeakage1()); - fLeakage2->Fill(h.GetLeakage2()); - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Use this function to setup your own conversion factor between degrees -// and millimeters. The conversion factor should be the one calculated in -// MGeomCam. Use this function with Caution: You could create wrong values -// by setting up your own scale factor. -// -void MHNewImagePar::SetMm2Deg(Float_t mmdeg) -{ - if (mmdeg<=0) - { - *fLog << warn << dbginf << "Warning - Conversion factor <= 0 - nonsense. Ignored." << endl; - return; - } - - if (fMm2Deg>0) - *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; - - fMm2Deg = mmdeg; -} - -// -------------------------------------------------------------------------- -// -// With this function you can convert the histogram ('on the fly') between -// degrees and millimeters. -// -void MHNewImagePar::SetMmScale(Bool_t mmscale) -{ - if (fUseMmScale == mmscale) - return; - - if (fMm2Deg<0) - { - *fLog << warn << GetDescriptor() << ": Warning - Sorry, no conversion factor for conversion available." << endl; - return; - } - - //const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; - //MH::ScaleAxis(fDist, scale); - //MH::ScaleAxis(fHeadTail, scale); - - if (mmscale) - { - // fDist->SetXTitle("Dist [mm]"); - // fHeadTail->SetXTitle("Head-Tail [mm]"); - } - else - { - // fDist->SetXTitle("Dist [\\circ]"); - // fHeadTail->SetXTitle("Head-Tail [\\circ]"); - } - - fUseMmScale = mmscale; -} - -// -------------------------------------------------------------------------- -// -// Draw clones of all two histograms. So that the object can be deleted -// and the histograms are still visible in the canvas. -// The cloned object are deleted together with the canvas if the canvas is -// destroyed. If you want to handle dostroying the canvas you can get a -// pointer to it from this function -// -TObject *MHNewImagePar::DrawClone(Option_t *opt) const -{ - TCanvas *c = MakeDefCanvas(this, 300, 600); - c->Divide(1, 2); - - gROOT->SetSelectedPad(NULL); - - c->cd(1); - fLeakage1->DrawCopy(); - - c->cd(2); - fLeakage2->DrawCopy(); - - c->Modified(); - c->Update(); - - return c; -} - -// -------------------------------------------------------------------------- -// -// Creates a new canvas and draws the two histograms into it. -// Be careful: The histograms belongs to this object and won't get deleted -// together with the canvas. -// -void MHNewImagePar::Draw(Option_t *) -{ - if (!gPad) - MakeDefCanvas(this, 300, 600); - - gPad->Divide(2, 2); - - gPad->cd(1); - fLeakage1->Draw(); - - gPad->cd(2); - fLeakage2->Draw(); - - gPad->Modified(); - gPad->Update(); -} - -TH1 *MHNewImagePar::GetHistByName(const TString name) -{ - if (name.Contains("Leakage1", TString::kIgnoreCase)) - return fLeakage1; - - if (name.Contains("Leakage2", TString::kIgnoreCase)) - return fLeakage2; - - return NULL; -} - - - - - - - - - - - - - - - - - Index: trunk/MagicSoft/Mars/mhist/MHNewImagePar.h =================================================================== --- trunk/MagicSoft/Mars/mhist/MHNewImagePar.h (revision 1939) +++ (revision ) @@ -1,43 +1,0 @@ -#ifndef MARS_MHNewImagePar -#define MARS_MHNewImagePar - -#ifndef MARS_MH -#include "MH.h" -#endif - -class TH1F; -class MHillas; - -class MHNewImagePar : public MH -{ -private: - TH1F *fLeakage1; //-> - TH1F *fLeakage2; //-> - - Float_t fMm2Deg; - Bool_t fUseMmScale; - -public: - MHNewImagePar(const char *name=NULL, const char *title=NULL); - ~MHNewImagePar(); - - void SetMmScale(Bool_t mmscale=kTRUE); - void SetMm2Deg(Float_t mmdeg); - - Bool_t SetupFill(const MParList *pList); - Bool_t Fill(const MParContainer *par); - - TH1 *GetHistByName(const TString name); - - TH1F *GetHistLeakage1() { return fLeakage1; } - TH1F *GetHistLeakage2() { return fLeakage2; } - - void Draw(Option_t *opt=NULL); - TObject *DrawClone(Option_t *opt=NULL) const; - - ClassDef(MHNewImagePar, 1) // Container which holds histograms for the new image parameters -}; - -#endif - - Index: trunk/MagicSoft/Mars/mhist/Makefile =================================================================== --- trunk/MagicSoft/Mars/mhist/Makefile (revision 1939) +++ trunk/MagicSoft/Mars/mhist/Makefile (revision 1940) @@ -23,5 +23,5 @@ # INCLUDES = -I. -I../mbase -I../mraw -I../manalysis -I../mmc \ - -I../mgui -I../mgeom -I../mdata -I../mfilter + -I../mgui -I../mgeom -I../mdata -I../mfilter -I../mimage #------------------------------------------------------------------------------ @@ -38,7 +38,4 @@ MHFadcPix.cc \ MHFadcCam.cc \ - MHHillas.cc \ - MHHillasSrc.cc \ - MHHillasExt.cc \ MHStarMap.cc \ MHEnergyTime.cc \ @@ -69,5 +66,4 @@ MHSigmaPixel.cc \ MHSigmabarTheta.cc \ - MHNewImagePar.cc \ MHSigmaTheta.cc \ MHMcCT1CollectionArea.cc \ Index: trunk/MagicSoft/Mars/mimage/ImageIncl.h =================================================================== --- trunk/MagicSoft/Mars/mimage/ImageIncl.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/ImageIncl.h (revision 1940) @@ -0,0 +1,3 @@ +#ifndef __CINT__ + +#endif // __CINT__ Index: trunk/MagicSoft/Mars/mimage/ImageLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mimage/ImageLinkDef.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/ImageLinkDef.h (revision 1940) @@ -0,0 +1,24 @@ +#ifdef __CINT__ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; + +#pragma link C++ class MImgCleanStd+; +#pragma link C++ class MCameraSmooth+; + +#pragma link C++ class MHillas+; +#pragma link C++ class MHillasSrc+; +#pragma link C++ class MHillasSrcCalc+; +#pragma link C++ class MHillasExt+; +#pragma link C++ class MHillasCalc+; + +#pragma link C++ class MNewImagePar+; +#pragma link C++ class MNewImageParCalc+; + +#pragma link C++ class MHHillas+; +#pragma link C++ class MHHillasSrc+; +#pragma link C++ class MHHillasExt+; +#pragma link C++ class MHNewImagePar+; + +#endif Index: trunk/MagicSoft/Mars/mimage/MCameraSmooth.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MCameraSmooth.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MCameraSmooth.cc (revision 1940) @@ -0,0 +1,143 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 08/2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MCameraSmooth // +// // +// This task fills each pixel in the camera with the average of the // +// number of cerenkov photons from its surrounding pixels. This can // +// be done using the central pixel or ignoring the central pixel. // +// // +///////////////////////////////////////////////////////////////////////////// +#include "MCameraSmooth.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MGeomPix.h" +#include "MGeomCam.h" +#include "MCerPhotPix.h" +#include "MCerPhotEvt.h" +#include "MBlindPixels.h" + +ClassImp(MCameraSmooth); + +// -------------------------------------------------------------------------- +// +// Default constructor. You can specify the numer of loops how often a +// smoothing should be done. The default is 1. +// +MCameraSmooth::MCameraSmooth(Byte_t n, const char *name, const char *title) + : fCounts(n), fUseCentralPixel(kTRUE) +{ + fName = name ? name : "MCameraSmooth"; + fTitle = title ? title : "Task to smooth the camera"; +} + +// -------------------------------------------------------------------------- +// +// - get the MCerPhotEvt from the parlist (abort if missing) +// - get MGeomCam from the parameter list +// +Bool_t MCameraSmooth::PreProcess (MParList *pList) +{ + fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fEvt) + { + *fLog << err << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fGeomCam) + { + *fLog << warn << dbginf << "No camera geometry available... aborting." << endl; + + return kFALSE; + } + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Do the smoothing +// +Bool_t MCameraSmooth::Process() +{ + const UShort_t entries = fEvt->GetNumPixels(); + + // + // remove the pixels in fPixelsID if they are set to be used, + // (set them to 'unused' state) + // + + Double_t *photons = new Double_t[entries]; + Double_t *errors = new Double_t[entries]; + + for (int n=0; nGetPixById(nid); + if (evtpix) + { + photons[i] += evtpix->GetNumPhotons(); + errors[i] += evtpix->GetErrorPhot(); + } + num++; + } + + photons[i] /= num; + errors[i] /= num; + } + + for (UShort_t i=0; i +#endif + +class MGeomCam; +class MCerPhotEvt; +class MBlindPixels; + +class MCameraSmooth : public MTask +{ +private: + MCerPhotEvt *fEvt; //! Pointer to data to smooth + MGeomCam *fGeomCam; //! Camera geometry + + Byte_t fCounts; // number of smoothing loops + Bool_t fUseCentralPixel; // use also the central pixel + +public: + MCameraSmooth(Byte_t cnt=1, const char *name=NULL, const char *title=NULL); + + void SetUseCetralPixel(Bool_t b=kTRUE) { fUseCentralPixel=kTRUE; } + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + + ClassDef(MCameraSmooth, 0) // task to smooth the camera contants +}; + +#endif + Index: trunk/MagicSoft/Mars/mimage/MHHillas.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHHillas.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHHillas.cc (revision 1940) @@ -0,0 +1,388 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 2001 +! Wolfgang Wittek 2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHHillas +// +// This class contains histograms for the source independent image parameters +// +///////////////////////////////////////////////////////////////////////////// + +#include "MHHillas.h" + +#include + +#include +#include +#include +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MGeomCam.h" +#include "MBinning.h" + + +ClassImp(MHHillas); + +// -------------------------------------------------------------------------- +// +// Setup four histograms for Width, Length +// +MHHillas::MHHillas(const char *name, const char *title) + : fMm2Deg(1), fUseMmScale(kTRUE) +{ + // + // set the name and title of this object + // + fName = name ? name : "MHHillas"; + fTitle = title ? title : "Source independent image parameters"; + + fLength = new TH1F("Length", "Length of Ellipse", 100, 0, 296.7); + fWidth = new TH1F("Width", "Width of Ellipse", 100, 0, 296.7); + fDistC = new TH1F("DistC", "Distance from center of camera", 100, 0, 445); + fDelta = new TH1F("Delta", "Angle (Main axis - x-axis)", 101, -90, 90); + fUsedPix = new TH1F("UsedPix", "Number of used pixels", 150, 0, 150); + fCorePix = new TH1F("CorePix", "Number of core pixels", 150, 0, 150); + + fLength->SetLineColor(kBlue); + fCorePix->SetLineColor(kRed); + fUsedPix->SetLineColor(kGreen); + + fLength->SetDirectory(NULL); + fWidth->SetDirectory(NULL); + fDistC->SetDirectory(NULL); + fDelta->SetDirectory(NULL); + fUsedPix->SetDirectory(NULL); + fCorePix->SetDirectory(NULL); + + fLength->SetXTitle("Length [mm]"); + fWidth->SetXTitle("Width [mm]"); + fDistC->SetXTitle("Distance [mm]"); + fDelta->SetXTitle("Delta [\\circ]"); + fUsedPix->SetXTitle("Number of Pixels"); + fCorePix->SetXTitle("Number of Pixels"); + + fLength->SetYTitle("Counts"); + fWidth->SetYTitle("Counts"); + fDistC->SetYTitle("Counts"); + fDelta->SetYTitle("Counts"); + fUsedPix->SetYTitle("Counts"); + fCorePix->SetYTitle("Counts"); + + MBinning bins; + bins.SetEdgesLog(50, 1, 1e7); + + fSize = new TH1F; + fSize->SetName("Size"); + fSize->SetTitle("Number of Photons"); + fSize->SetDirectory(NULL); + fSize->SetXTitle("Size"); + fSize->SetYTitle("Counts"); + fSize->GetXaxis()->SetTitleOffset(1.2); + fSize->GetXaxis()->SetLabelOffset(-0.015); + + bins.Apply(*fSize); + + fCenter = new TH2F("Center", "Center of Ellipse", 51, -445, 445, 51, -445, 445); + fCenter->SetDirectory(NULL); + fCenter->SetXTitle("x [mm]"); + fCenter->SetYTitle("y [mm]"); + fCenter->SetZTitle("Counts"); +} + +// -------------------------------------------------------------------------- +// +// Delete the histograms +// +MHHillas::~MHHillas() +{ + delete fLength; + delete fWidth; + + delete fDistC; + delete fDelta; + + delete fSize; + delete fCenter; + + delete fUsedPix; + delete fCorePix; +} + +// -------------------------------------------------------------------------- +// +// Setup the Binning for the histograms automatically if the correct +// instances of MBinning (with the names 'BinningWidth' and 'BinningLength') +// are found in the parameter list +// Use this function if you want to set the conversion factor which +// is used to convert the mm-scale in the camera plain into the deg-scale +// used for histogram presentations. The conversion factor is part of +// the camera geometry. Please create a corresponding MGeomCam container. +// +Bool_t MHHillas::SetupFill(const MParList *plist) +{ + const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); + if (!geom) + *fLog << warn << GetDescriptor() << ": No Camera Geometry available. Using mm-scale for histograms." << endl; + else + { + fMm2Deg = geom->GetConvMm2Deg(); + SetMmScale(kFALSE); + } + + ApplyBinning(*plist, "Width", fWidth); + ApplyBinning(*plist, "Length", fLength); + ApplyBinning(*plist, "Dist", fDistC); + ApplyBinning(*plist, "Delta", fDelta); + ApplyBinning(*plist, "Size", fSize); + ApplyBinning(*plist, "Pixels", fUsedPix); + ApplyBinning(*plist, "Pixels", fCorePix); + + const MBinning *bins = (MBinning*)plist->FindObject("BinningCamera"); + if (!bins) + { + float r = geom ? geom->GetMaxRadius() : 600; + r *= 0.9; + if (!fUseMmScale) + r *= fMm2Deg; + + MBinning b; + b.SetEdges(61, -r, r); + SetBinning(fCenter, &b, &b); + } + else + SetBinning(fCenter, bins, bins); + + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Use this function to setup your own conversion factor between degrees +// and millimeters. The conversion factor should be the one calculated in +// MGeomCam. Use this function with Caution: You could create wrong values +// by setting up your own scale factor. +// +void MHHillas::SetMm2Deg(Float_t mmdeg) +{ + if (mmdeg<0) + { + *fLog << warn << dbginf << "Warning - Conversion factor < 0 - nonsense. Ignored." << endl; + return; + } + + if (fMm2Deg>=0) + *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; + + fMm2Deg = mmdeg; +} + +// -------------------------------------------------------------------------- +// +// With this function you can convert the histogram ('on the fly') between +// degrees and millimeters. +// +void MHHillas::SetMmScale(Bool_t mmscale) +{ + if (fUseMmScale == mmscale) + return; + + if (fMm2Deg<0) + { + *fLog << warn << dbginf << "Warning - Sorry, no conversion factor for conversion available." << endl; + return; + } + + const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; + MH::ScaleAxis(fLength, scale); + MH::ScaleAxis(fWidth, scale); + MH::ScaleAxis(fDistC, scale); + MH::ScaleAxis(fCenter, scale, scale); + + if (mmscale) + { + fLength->SetXTitle("Length [mm]"); + fWidth->SetXTitle("Width [mm]"); + fDistC->SetXTitle("Distance [mm]"); + fCenter->SetXTitle("x [mm]"); + fCenter->SetYTitle("y [mm]"); + } + else + { + fLength->SetXTitle("Length [\\circ]"); + fWidth->SetXTitle("Width [\\circ]"); + fDistC->SetXTitle("Distance [\\circ]"); + fCenter->SetXTitle("x [\\circ]"); + fCenter->SetYTitle("y [\\circ]"); + } + + fUseMmScale = mmscale; +} + +// -------------------------------------------------------------------------- +// +// Fill the histograms with data from a MHillas-Container. +// Be careful: Only call this with an object of type MHillas +// +Bool_t MHHillas::Fill(const MParContainer *par) +{ + const MHillas &h = *(MHillas*)par; + + const Double_t d = sqrt(h.GetMeanX()*h.GetMeanX() + h.GetMeanY()*h.GetMeanY()); + const Double_t scale = fUseMmScale ? 1 : fMm2Deg; + + fLength ->Fill(scale*h.GetLength()); + fWidth ->Fill(scale*h.GetWidth()); + fDistC ->Fill(scale*d); + fCenter ->Fill(scale*h.GetMeanX(), scale*h.GetMeanY()); + fDelta ->Fill(kRad2Deg*h.GetDelta()); + fSize ->Fill(h.GetSize()); + fUsedPix->Fill(h.GetNumUsedPixels()); + fCorePix->Fill(h.GetNumCorePixels()); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Setup a inversed deep blue sea palette for the fCenter histogram. +// +void MHHillas::SetColors() const +{ + // FIXME: This must be redone each time the canvas is repainted.... + gStyle->SetPalette(51, NULL); + Int_t c[50]; + for (int i=0; i<50; i++) + c[49-i] = gStyle->GetColorPalette(i); + gStyle->SetPalette(50, c); +} + +// -------------------------------------------------------------------------- +// +// Draw clones of four histograms. So that the object can be deleted +// and the histograms are still visible in the canvas. +// The cloned object are deleted together with the canvas if the canvas is +// destroyed. If you want to handle dostroying the canvas you can get a +// pointer to it from this function +// +TObject *MHHillas::DrawClone(Option_t *opt) const +{ + TCanvas *c = MakeDefCanvas(this, 720, 810); + c->Divide(2,3); + + gROOT->SetSelectedPad(NULL); + + c->cd(1); + DrawCopy(*fWidth, *fLength, "Width / Length"); + + c->cd(2); + gPad->SetLogx(); + fSize->DrawCopy(); + + c->cd(3); + DrawCopy(*fCorePix, *fUsedPix, "Number of core/used Pixels"); + + c->cd(4); + fDelta->DrawCopy(); + + c->cd(5); + fDistC->DrawCopy(); + + c->cd(6); + SetColors(); + fCenter->DrawCopy("colz"); + + c->Modified(); + c->Update(); + + return c; +} + +// -------------------------------------------------------------------------- +// +// Creates a new canvas and draws the four histograms into it. +// Be careful: The histograms belongs to this object and won't get deleted +// together with the canvas. +// +void MHHillas::Draw(Option_t *) +{ + if (!gPad) + MakeDefCanvas(this, 720, 810); + + gPad->Divide(2,3); + + gPad->cd(1); + MH::Draw(*fWidth, *fLength, "Width / Length"); + + gPad->cd(2); + gPad->SetLogx(); + fSize->Draw(); + + gPad->cd(3); + MH::Draw(*fCorePix, *fUsedPix, "Number of core/used Pixels"); + + gPad->cd(4); + fDelta->Draw(); + + gPad->cd(5); + fDistC->Draw(); + + gPad->cd(6); + SetColors(); + fCenter->Draw("colz"); + + gPad->Modified(); + gPad->Update(); +} + +TH1 *MHHillas::GetHistByName(const TString name) +{ + if (name.Contains("Width", TString::kIgnoreCase)) + return fWidth; + if (name.Contains("Length", TString::kIgnoreCase)) + return fLength; + if (name.Contains("Size", TString::kIgnoreCase)) + return fSize; + if (name.Contains("Core", TString::kIgnoreCase)) + return fCorePix; + if (name.Contains("Used", TString::kIgnoreCase)) + return fUsedPix; + if (name.Contains("Delta", TString::kIgnoreCase)) + return fDelta; + if (name.Contains("DistC", TString::kIgnoreCase)) + return fDistC; + if (name.Contains("Center", TString::kIgnoreCase)) + return fCenter; + + return NULL; +} Index: trunk/MagicSoft/Mars/mimage/MHHillas.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHHillas.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHHillas.h (revision 1940) @@ -0,0 +1,60 @@ +#ifndef MARS_MHHillas +#define MARS_MHHillas + +#ifndef MARS_MH +#include "MH.h" +#endif + +class TH1F; +class TH2F; +class MHillas; + +class MHHillas : public MH +{ +private: + + TH1F *fLength; //-> + TH1F *fWidth; //-> + + TH1F *fDistC; //-> + TH1F *fDelta; //-> + + TH1F *fSize; //-> + TH2F *fCenter; //-> + + TH1F *fUsedPix; //-> + TH1F *fCorePix; //-> + + void SetColors() const; + + Float_t fMm2Deg; + Bool_t fUseMmScale; + +public: + MHHillas(const char *name=NULL, const char *title=NULL); + ~MHHillas(); + + void SetMmScale(Bool_t mmscale=kTRUE); + virtual void SetMm2Deg(Float_t mmdeg); + + Bool_t SetupFill(const MParList *pList); + Bool_t Fill(const MParContainer *par); + + TH1 *GetHistByName(const TString name); + + TH1F *GetHistLength() { return fLength; } + TH1F *GetHistWidth() { return fWidth; } + + TH1F *GetHistDistC() { return fDistC; } + TH1F *GetHistDelta() { return fDelta; } + + TH1F *GetHistSize() { return fSize; } + TH2F *GetHistCenter() { return fCenter; } + + void Draw(Option_t *opt=NULL); + TObject *DrawClone(Option_t *opt=NULL) const; + + ClassDef(MHHillas, 1) // Container which holds histograms for the source independent image parameters +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHHillasExt.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHHillasExt.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHHillasExt.h (revision 1940) @@ -0,0 +1,50 @@ +#ifndef MARS_MHHillasExt +#define MARS_MHHillasExt + +#ifndef ROOT_TH1 +#include +#endif +#ifndef MARS_MH +#include "MH.h" +#endif + +class MHillasExt; + +class MHHillasExt : public MH +{ +private: + MHillasExt *fHillasExt; //! Pointer to the MHillasExt container + + TH1F fHConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize + TH1F fHConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize + TH1F fHAsym; // [mm] fDist minus dist: center of ellipse, highest pixel + TH1F fHM3Long; // [mm] 3rd moment (e-weighted) along major axis + TH1F fHM3Trans; // [mm] 3rd moment (e-weighted) along minor axis + + Float_t fMm2Deg; + Bool_t fUseMmScale; + + TString fHilName; + +public: + //MHHillasExt(const char *name=NULL, const char *title=NULL); + MHHillasExt(const char *name=NULL, const char *title=NULL, + const char *hil="MHillas"); + + ~MHHillasExt(); + + void SetMmScale(Bool_t mmscale=kTRUE); + virtual void SetMm2Deg(Float_t mmdeg); + + Bool_t SetupFill(const MParList *pList); + Bool_t Fill(const MParContainer *par); + + TH1 *GetHistByName(const TString name); + + void Draw(Option_t *opt=NULL); + TObject *DrawClone(Option_t *opt=NULL) const; + + ClassDef(MHHillasExt, 1) // Container which holds histograms for the extended hillas parameters +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHHillasSrc.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHHillasSrc.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHHillasSrc.cc (revision 1940) @@ -0,0 +1,286 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 2001 +! +! Copyright: MAGIC Software Development, 2000-2001 +! +! +\* ======================================================================== */ + +/////////////////////////////////////////////////////////////////////// +// +// MHHillasSrc +// +// This class contains histograms for every Hillas parameter +// +/////////////////////////////////////////////////////////////////////// +#include "MHHillasSrc.h" + +#include + +#include +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MGeomCam.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MHillasSrc.h" + +ClassImp(MHHillasSrc); + +// -------------------------------------------------------------------------- +// +// Setup four histograms for Alpha, and Dist +// +MHHillasSrc::MHHillasSrc(const char *name, const char *title) + : fUseMmScale(kTRUE) +{ + // + // set the name and title of this object + // + fName = name ? name : "MHHillasSrc"; + fTitle = title ? title : "Container for Hillas histograms"; + + // + // loop over all Pixels and create two histograms + // one for the Low and one for the High gain + // connect all the histogram with the container fHist + // + fAlpha = new TH1F("Alpha", "Alpha of Ellipse", 181, -90, 90); + fDist = new TH1F("Dist", "Dist of Ellipse", 100, 0, 445); + fHeadTail = new TH1F("HeadTail", "HeadTail of Ellipse", 101, -445, 445); + fCosDA = new TH1F("CosDA", "cos(Delta,Alpha) of Ellipse", 101, -1, 1); + + fAlpha->SetDirectory(NULL); + fDist->SetDirectory(NULL); + fHeadTail->SetDirectory(NULL); + fCosDA->SetDirectory(NULL); + + fAlpha->SetXTitle("\\alpha [\\circ]"); + fDist->SetXTitle("Dist [mm]"); + fHeadTail->SetXTitle("Head-Tail [mm]"); + fCosDA->SetXTitle("cos(\\delta,\\alpha)"); + + fAlpha->SetYTitle("Counts"); + fDist->SetYTitle("Counts"); + fHeadTail->SetYTitle("Counts"); + fCosDA->SetYTitle("Counts"); +} + +// -------------------------------------------------------------------------- +// +// Delete the four histograms +// +MHHillasSrc::~MHHillasSrc() +{ + delete fAlpha; + delete fDist; + delete fHeadTail; + delete fCosDA; +} + +// -------------------------------------------------------------------------- +// +// Setup the Binning for the histograms automatically if the correct +// instances of MBinning (with the names 'BinningAlpha' and 'BinningDist') +// are found in the parameter list +// Use this function if you want to set the conversion factor which +// is used to convert the mm-scale in the camera plain into the deg-scale +// used for histogram presentations. The conversion factor is part of +// the camera geometry. Please create a corresponding MGeomCam container. +// +Bool_t MHHillasSrc::SetupFill(const MParList *plist) +{ + const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); + if (!geom) + *fLog << warn << dbginf << "No Camera Geometry available. Using mm-scale for histograms." << endl; + else + { + fMm2Deg = geom->GetConvMm2Deg(); + SetMmScale(kFALSE); + } + + ApplyBinning(*plist, "Alpha", fAlpha); + ApplyBinning(*plist, "Dist", fDist); + ApplyBinning(*plist, "HeadTail", fHeadTail); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Fill the four histograms with data from a MHillas-Container. +// Be careful: Only call this with an object of type MHillas +// +Bool_t MHHillasSrc::Fill(const MParContainer *par) +{ + const MHillasSrc &h = *(MHillasSrc*)par; + + fAlpha ->Fill(h.GetAlpha()); + fDist ->Fill(fUseMmScale ? h.GetDist() : fMm2Deg*h.GetDist()); + fHeadTail->Fill(fUseMmScale ? h.GetHeadTail() : fMm2Deg*h.GetHeadTail()); + fCosDA ->Fill(h.GetCosDeltaAlpha()); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Use this function to setup your own conversion factor between degrees +// and millimeters. The conversion factor should be the one calculated in +// MGeomCam. Use this function with Caution: You could create wrong values +// by setting up your own scale factor. +// +void MHHillasSrc::SetMm2Deg(Float_t mmdeg) +{ + if (mmdeg<=0) + { + *fLog << warn << dbginf << "Warning - Conversion factor <= 0 - nonsense. Ignored." << endl; + return; + } + + if (fMm2Deg>0) + *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; + + fMm2Deg = mmdeg; +} + +// -------------------------------------------------------------------------- +// +// With this function you can convert the histogram ('on the fly') between +// degrees and millimeters. +// +void MHHillasSrc::SetMmScale(Bool_t mmscale) +{ + if (fUseMmScale == mmscale) + return; + + if (fMm2Deg<0) + { + *fLog << warn << GetDescriptor() << ": Warning - Sorry, no conversion factor for conversion available." << endl; + return; + } + + const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; + MH::ScaleAxis(fDist, scale); + MH::ScaleAxis(fHeadTail, scale); + + if (mmscale) + { + fDist->SetXTitle("Dist [mm]"); + fHeadTail->SetXTitle("Head-Tail [mm]"); + } + else + { + fDist->SetXTitle("Dist [\\circ]"); + fHeadTail->SetXTitle("Head-Tail [\\circ]"); + } + + fUseMmScale = mmscale; +} + +// -------------------------------------------------------------------------- +// +// Draw clones of all two histograms. So that the object can be deleted +// and the histograms are still visible in the canvas. +// The cloned object are deleted together with the canvas if the canvas is +// destroyed. If you want to handle dostroying the canvas you can get a +// pointer to it from this function +// +TObject *MHHillasSrc::DrawClone(Option_t *opt) const +{ + TCanvas *c = MakeDefCanvas(this, 700, 500); + c->Divide(2, 2); + + // FIXME: Display Source position + + gROOT->SetSelectedPad(NULL); + + // + // This is necessary to get the expected bahviour of DrawClone + // + c->cd(1); + fAlpha->DrawCopy(); + + c->cd(2); + fDist->DrawCopy(); + + c->cd(3); + fHeadTail->DrawCopy(); + + c->cd(4); + gPad->SetLogy(); + fCosDA->DrawCopy(); + + c->Modified(); + c->Update(); + + return c; +} + +// -------------------------------------------------------------------------- +// +// Creates a new canvas and draws the two histograms into it. +// Be careful: The histograms belongs to this object and won't get deleted +// together with the canvas. +// +void MHHillasSrc::Draw(Option_t *) +{ + if (!gPad) + MakeDefCanvas(this, 700, 500); + + // FIXME: Display Source position + + gPad->Divide(2, 2); + + gPad->cd(1); + fAlpha->Draw(); + + gPad->cd(2); + fDist->Draw(); + + gPad->cd(1); + fHeadTail->Draw(); + + gPad->cd(2); + gPad->SetLogy(); + fCosDA->Draw(); + + gPad->Modified(); + gPad->Update(); +} + +TH1 *MHHillasSrc::GetHistByName(const TString name) +{ + if (name.Contains("Alpha", TString::kIgnoreCase)) + return fAlpha; + if (name.Contains("Dist", TString::kIgnoreCase)) + return fDist; + if (name.Contains("HeadTail", TString::kIgnoreCase)) + return fHeadTail; + if (name.Contains("CosDA", TString::kIgnoreCase)) + return fCosDA; + + return NULL; +} Index: trunk/MagicSoft/Mars/mimage/MHHillasSrc.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHHillasSrc.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHHillasSrc.h (revision 1940) @@ -0,0 +1,45 @@ +#ifndef MARS_MHHillasSrc +#define MARS_MHHillasSrc + +#ifndef MARS_MH +#include "MH.h" +#endif + +class TH1F; +class MHillas; + +class MHHillasSrc : public MH +{ +private: + TH1F *fAlpha; //-> + TH1F *fDist; //-> + TH1F *fHeadTail; //-> + TH1F *fCosDA; //-> + + Float_t fMm2Deg; + Bool_t fUseMmScale; + +public: + MHHillasSrc(const char *name=NULL, const char *title=NULL); + ~MHHillasSrc(); + + void SetMmScale(Bool_t mmscale=kTRUE); + void SetMm2Deg(Float_t mmdeg); + + Bool_t SetupFill(const MParList *pList); + Bool_t Fill(const MParContainer *par); + + TH1 *GetHistByName(const TString name); + + TH1F *GetHistAlpha() { return fAlpha; } + TH1F *GetHistDist() { return fDist; } + TH1F *GetHistHeadTail() { return fHeadTail; } + TH1F *GetHistCosDeltaAlpha() { return fCosDA; } + + void Draw(Option_t *opt=NULL); + TObject *DrawClone(Option_t *opt=NULL) const; + + ClassDef(MHHillasSrc, 1) // Container which holds histograms for the source dependant parameters +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc (revision 1940) @@ -0,0 +1,256 @@ +/* ======================================================================== *\ +! +! * +! * 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 03/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +/////////////////////////////////////////////////////////////////////// +// +// MHNewImagePar +// +// This class contains histograms for every Hillas parameter +// +/////////////////////////////////////////////////////////////////////// +#include "MHNewImagePar.h" + +#include + +#include +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MGeomCam.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MNewImagePar.h" + +ClassImp(MHNewImagePar); + +// -------------------------------------------------------------------------- +// +// Setup histograms +// +MHNewImagePar::MHNewImagePar(const char *name, const char *title) + : fUseMmScale(kTRUE) +{ + fName = name ? name : "MHNewImagePar"; + fTitle = title ? title : "Container for histograms of new image parameters"; + + fLeakage1 = new TH1F("Leakage1", "Leakage1", 100, 0.0, 1.0); + fLeakage1->SetDirectory(NULL); + fLeakage1->SetXTitle("Leakage1"); + fLeakage1->SetYTitle("Counts"); + + fLeakage2 = new TH1F("Leakage2", "Leakage2", 100, 0.0, 1.0); + fLeakage2->SetDirectory(NULL); + fLeakage2->SetXTitle("Leakage2"); + fLeakage2->SetYTitle("Counts"); +} + +// -------------------------------------------------------------------------- +// +// Delete the four histograms +// +MHNewImagePar::~MHNewImagePar() +{ + delete fLeakage1; + delete fLeakage2; +} + +// -------------------------------------------------------------------------- +// +// Setup the Binning for the histograms automatically if the correct +// instances of MBinning (with the names 'BinningAlpha' and 'BinningDist') +// are found in the parameter list +// Use this function if you want to set the conversion factor which +// is used to convert the mm-scale in the camera plain into the deg-scale +// used for histogram presentations. The conversion factor is part of +// the camera geometry. Please create a corresponding MGeomCam container. +// +Bool_t MHNewImagePar::SetupFill(const MParList *plist) +{ + const MGeomCam *geom = (MGeomCam*)plist->FindObject("MGeomCam"); + if (!geom) + *fLog << warn << dbginf << "No Camera Geometry available. Using mm-scale for histograms." << endl; + else + { + fMm2Deg = geom->GetConvMm2Deg(); + SetMmScale(kFALSE); + } + + //ApplyBinning(*plist, "Alpha", fAlpha); + //ApplyBinning(*plist, "Dist", fDist); + //ApplyBinning(*plist, "HeadTail", fHeadTail); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Fill the histograms with data from a MNewImagePar container. +// +Bool_t MHNewImagePar::Fill(const MParContainer *par) +{ + const MNewImagePar &h = *(MNewImagePar*)par; + + fLeakage1->Fill(h.GetLeakage1()); + fLeakage2->Fill(h.GetLeakage2()); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Use this function to setup your own conversion factor between degrees +// and millimeters. The conversion factor should be the one calculated in +// MGeomCam. Use this function with Caution: You could create wrong values +// by setting up your own scale factor. +// +void MHNewImagePar::SetMm2Deg(Float_t mmdeg) +{ + if (mmdeg<=0) + { + *fLog << warn << dbginf << "Warning - Conversion factor <= 0 - nonsense. Ignored." << endl; + return; + } + + if (fMm2Deg>0) + *fLog << warn << dbginf << "Warning - Conversion factor already set. Overwriting" << endl; + + fMm2Deg = mmdeg; +} + +// -------------------------------------------------------------------------- +// +// With this function you can convert the histogram ('on the fly') between +// degrees and millimeters. +// +void MHNewImagePar::SetMmScale(Bool_t mmscale) +{ + if (fUseMmScale == mmscale) + return; + + if (fMm2Deg<0) + { + *fLog << warn << GetDescriptor() << ": Warning - Sorry, no conversion factor for conversion available." << endl; + return; + } + + //const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg; + //MH::ScaleAxis(fDist, scale); + //MH::ScaleAxis(fHeadTail, scale); + + if (mmscale) + { + // fDist->SetXTitle("Dist [mm]"); + // fHeadTail->SetXTitle("Head-Tail [mm]"); + } + else + { + // fDist->SetXTitle("Dist [\\circ]"); + // fHeadTail->SetXTitle("Head-Tail [\\circ]"); + } + + fUseMmScale = mmscale; +} + +// -------------------------------------------------------------------------- +// +// Draw clones of all two histograms. So that the object can be deleted +// and the histograms are still visible in the canvas. +// The cloned object are deleted together with the canvas if the canvas is +// destroyed. If you want to handle dostroying the canvas you can get a +// pointer to it from this function +// +TObject *MHNewImagePar::DrawClone(Option_t *opt) const +{ + TCanvas *c = MakeDefCanvas(this, 300, 600); + c->Divide(1, 2); + + gROOT->SetSelectedPad(NULL); + + c->cd(1); + fLeakage1->DrawCopy(); + + c->cd(2); + fLeakage2->DrawCopy(); + + c->Modified(); + c->Update(); + + return c; +} + +// -------------------------------------------------------------------------- +// +// Creates a new canvas and draws the two histograms into it. +// Be careful: The histograms belongs to this object and won't get deleted +// together with the canvas. +// +void MHNewImagePar::Draw(Option_t *) +{ + if (!gPad) + MakeDefCanvas(this, 300, 600); + + gPad->Divide(2, 2); + + gPad->cd(1); + fLeakage1->Draw(); + + gPad->cd(2); + fLeakage2->Draw(); + + gPad->Modified(); + gPad->Update(); +} + +TH1 *MHNewImagePar::GetHistByName(const TString name) +{ + if (name.Contains("Leakage1", TString::kIgnoreCase)) + return fLeakage1; + + if (name.Contains("Leakage2", TString::kIgnoreCase)) + return fLeakage2; + + return NULL; +} + + + + + + + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mimage/MHNewImagePar.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHNewImagePar.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHNewImagePar.h (revision 1940) @@ -0,0 +1,43 @@ +#ifndef MARS_MHNewImagePar +#define MARS_MHNewImagePar + +#ifndef MARS_MH +#include "MH.h" +#endif + +class TH1F; +class MHillas; + +class MHNewImagePar : public MH +{ +private: + TH1F *fLeakage1; //-> + TH1F *fLeakage2; //-> + + Float_t fMm2Deg; + Bool_t fUseMmScale; + +public: + MHNewImagePar(const char *name=NULL, const char *title=NULL); + ~MHNewImagePar(); + + void SetMmScale(Bool_t mmscale=kTRUE); + void SetMm2Deg(Float_t mmdeg); + + Bool_t SetupFill(const MParList *pList); + Bool_t Fill(const MParContainer *par); + + TH1 *GetHistByName(const TString name); + + TH1F *GetHistLeakage1() { return fLeakage1; } + TH1F *GetHistLeakage2() { return fLeakage2; } + + void Draw(Option_t *opt=NULL); + TObject *DrawClone(Option_t *opt=NULL) const; + + ClassDef(MHNewImagePar, 1) // Container which holds histograms for the new image parameters +}; + +#endif + + Index: trunk/MagicSoft/Mars/mimage/MHillas.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillas.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillas.cc (revision 1940) @@ -0,0 +1,432 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Harald Kornmayer 1/2001 +! Author(s): Rudolf Bock 10/2001 +! Author(s): Wolfgang Wittek 6/2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHillas +// +// Storage Container for image parameters +// +// basic image parameters +// +// Version 1: +// ---------- +// fLength [mm] major axis of ellipse +// fWidth [mm] minor axis +// fDelta [rad] angle of major axis with x-axis +// by definition the major axis is pointing into +// the hemisphere x>0, thus -pi/2 < delta < pi/2 +// fSize [#CerPhot] total sum of pixels +// fMeanX [mm] x of center of ellipse +// fMeanY [mm] y of center of ellipse +// +// Version 2: +// ---------- +// fNumCorePixels number of pixels called core +// fNumUsedPixels number of pixels which survived the cleaning +// +///////////////////////////////////////////////////////////////////////////// +#include "MHillas.h" + +#include + +#include +#include + +#include "MGeomPix.h" +#include "MGeomCam.h" + +#include "MCerPhotPix.h" +#include "MCerPhotEvt.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MHillas); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MHillas::MHillas(const char *name, const char *title) : fEllipse(NULL) +{ + fName = name ? name : "MHillas"; + fTitle = title ? title : "Storage container for image parameters of one event"; + + Reset(); + + fEllipse = new TEllipse; +} + +// -------------------------------------------------------------------------- +// +// Destructor. Deletes the TEllipse if one exists. +// +MHillas::~MHillas() +{ + Clear(); +} + +// -------------------------------------------------------------------------- +// +// Initializes the values with defaults. For the default values see the +// source code. +// +void MHillas::Reset() +{ + fLength = -1; + fWidth = -1; + fDelta = 0; + + fSize = -1; + fMeanX = -1; + fMeanY = -1; + + fNumUsedPixels = -1; + fNumCorePixels = -1; + + Clear(); +} + +// -------------------------------------------------------------------------- +// +// Print the hillas Parameters to *fLog +// +void MHillas::Print(Option_t *) const +{ + Double_t atg = atan2(fMeanY, fMeanX)*kRad2Deg; + + if (atg<0) + atg += 180; + + *fLog << all; + *fLog << "Basic Image Parameters (" << GetName() << ")" << endl; + *fLog << " - Length [mm] = " << fLength << endl; + *fLog << " - Width [mm] = " << fWidth << endl; + *fLog << " - Delta [deg] = " << fDelta*kRad2Deg << endl; + *fLog << " - Size [1] = " << fSize << " #CherPhot" << endl; + *fLog << " - Meanx [mm] = " << fMeanX << endl; + *fLog << " - Meany [mm] = " << fMeanY << endl; + *fLog << " - atg(y/x) [deg] = " << atg << endl; + *fLog << " - Used Pixels [#] = " << fNumUsedPixels << " Pixels" << endl; + *fLog << " - Core Pixels [#] = " << fNumCorePixels << " Pixels" << endl; +} + +/* +// ----------------------------------------------------------- +// +// call the Paint function of the Ellipse if a TEllipse exists +// +void MHillas::Paint(Option_t *) +{ + fEllipse->SetLineWidth(2); + fEllipse->PaintEllipse(fMeanX, fMeanY, fLength, fWidth, + 0, 360, fDelta*kRad2Deg+180); +} +*/ + +// -------------------------------------------------------------------------- +// +// Instead of adding MHillas itself to the Pad +// (s. AppendPad in TObject) we create an ellipse, +// which is added to the Pad by its Draw function +// You can remove it by deleting the Ellipse Object +// (s. Clear() ) +// +void MHillas::Draw(Option_t *opt) +{ + + Clear(); + + if (fLength<0 || fWidth<0) + return; + + fEllipse = new TEllipse(fMeanX, fMeanY, fLength, fWidth, + 0, 360, fDelta*kRad2Deg+180); + + fEllipse->SetLineWidth(2); + fEllipse->Draw(); + + /* + fEllipse->SetPhimin(); + fEllipse->SetPhimax(); + fEllipse->SetR1(fLength); + fEllipse->SetR2(fWidth); + fEllipse->SetTheta(fDelta*kRad2Deg+180); + fEllipse->SetX1(fMeanX); + fEllipse->SetY1(fMeanY); + + fEllipse->SetLineWidth(2); + fEllipse->PaintEllipse(fMeanX, fMeanY, fLength, fWidth, + 0, 360, fDelta*kRad2Deg+180); + + AppendPad(opt); + + // This is from TH1 + TString opt = option; + opt.ToLower(); + if (gPad && !opt.Contains("same")) { + //the following statement is necessary in case one attempts to draw + //a temporary histogram already in the current pad + if (TestBit(kCanDelete)) gPad->GetListOfPrimitives()->Remove(this); + gPad->Clear(); + } + */ +} + +// -------------------------------------------------------------------------- +// +// If a TEllipse object exists it is deleted +// +void MHillas::Clear(Option_t *) +{ + if (!fEllipse) + return; + + delete fEllipse; + + fEllipse = NULL; +} + + +// -------------------------------------------------------------------------- +// +// Calculate the image parameters from a Cherenkov photon event +// assuming Cher.photons/pixel and pixel coordinates are given +// In case you don't call Calc from within an eventloop make sure, that +// you call the Reset member function before. +// Returns: +// 0 no error +// 1 number of pixels < 3 +// 2 size==0 +// 3 number of used pixel < 3 +// 4 CorrXY == 0 +// +Int_t MHillas::Calc(const MGeomCam &geom, const MCerPhotEvt &evt) +{ + const UInt_t npixevt = evt.GetNumPixels(); + + // + // sanity check + // + if (npixevt < 3) + { + //*fLog << warn << "MHillas::Calc: Event has less than three pixels... skipped." << endl; + return 1; + } + + // + // calculate mean value of pixel coordinates and fSize + // ----------------------------------------------------- + // + // Because this are only simple sums of roughly 600 values + // with an accuracy less than three digits there should not + // be any need to use double precision (Double_t) for the + // calculation of fMeanX, fMeanY and fSize. + // + fMeanX = 0; + fMeanY = 0; + fSize = 0; + + fNumUsedPixels = 0; + fNumCorePixels = 0; + + for (UInt_t i=0; i0) we + // cannot calculate Length and Width. The calculation failed + // and returns kFALSE + // In reallity it is almost impossible to have a distribution + // of cerenkov photons in the used pixels which is exactly symmetric + // along one of the axis. + // + if (corrxy==0) + { + //*fLog << inf << GetDescriptor() << ": Event has CorrXY==0... skipped." << endl; + return 4; + } + + // + // calculate the basic Hillas parameters: orientation and size of axes + // ------------------------------------------------------------------- + // + // fDelta is the angle between the major axis of the ellipse and + // the x axis. It is independent of the position of the ellipse + // in the camera it has values between -pi/2 and pi/2 degrees + // + const Double_t d0 = corryy - corrxx; + const Double_t d1 = corrxy*2; + const Double_t d2 = d0 + sqrt(d0*d0 + d1*d1); + const Double_t tand = d2 / d1; + + fDelta = atan(tand); + + const Double_t s2 = tand*tand+1; + const Double_t s = sqrt(s2); + + fCosDelta = 1.0/s; // need these in derived classes + fSinDelta = tand/s; // like MHillasExt + + Double_t axis1 = (tand*tand*corryy + d2 + corrxx)/s2/fSize; + Double_t axis2 = (tand*tand*corrxx - d2 + corryy)/s2/fSize; + + // + // fLength^2 is the second moment along the major axis of the ellipse + // fWidth^2 is the second moment along the minor axis of the ellipse + // + // From the algorithm we get: fWidth <= fLength is always true + // + // very small numbers can get negative by rounding + // + fLength = axis1<0 ? 0 : sqrt(axis1); // [mm] + fWidth = axis2<0 ? 0 : sqrt(axis2); // [mm] + + SetReadyToSave(); + + return 0; +} + +// -------------------------------------------------------------------------- +// +// This function is ment for special usage, please never try to set +// values via this function +// +void MHillas::Set(const TArrayF &arr) +{ + if (arr.GetSize() != 8) + return; + + fLength = arr.At(0); // [mm] major axis of ellipse + fWidth = arr.At(1); // [mm] minor axis of ellipse + fDelta = arr.At(2); // [rad] angle of major axis with x-axis + fSize = arr.At(3); // [#CerPhot] sum of content of all pixels (number of Cherenkov photons) + fMeanX = arr.At(4); // [mm] x-coordinate of center of ellipse + fMeanY = arr.At(5); // [mm] y-coordinate of center of ellipse + + fNumUsedPixels = (Short_t)arr.At(6); // Number of pixels which survived the image cleaning + fNumCorePixels = (Short_t)arr.At(7); // number of core pixels +} + + +// -------------------------------------------------------------------------- +// +/* +void MHillas::AsciiRead(ifstream &fin) +{ + fin >> fLength; + fin >> fWidth; + fin >> fDelta; + fin >> fSize; + fin >> fMeanX; + fin >> fMeanY; +} +*/ +// -------------------------------------------------------------------------- +/* +void MHillas::AsciiWrite(ofstream &fout) const +{ + fout << fLength << " "; + fout << fWidth << " "; + fout << fDelta << " "; + fout << fSize << " "; + fout << fMeanX << " "; + fout << fMeanY; +} +*/ Index: trunk/MagicSoft/Mars/mimage/MHillas.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillas.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillas.h (revision 1940) @@ -0,0 +1,73 @@ +#ifndef MARS_MHillas +#define MARS_MHillas + +#ifndef MARS_MParContainer +#include "MParContainer.h" +#endif + +class TArrayF; +class TEllipse; + +class MGeomCam; +class MCerPhotEvt; + +class MHillas : public MParContainer +{ +private: + // for description see MHillas.cc + Float_t fLength; // [mm] major axis of ellipse + Float_t fWidth; // [mm] minor axis of ellipse + Float_t fDelta; // [rad] angle of major axis with x-axis + Float_t fSize; // [#CerPhot] sum of content of all pixels (number of Cherenkov photons) + Float_t fMeanX; // [mm] x-coordinate of center of ellipse + Float_t fMeanY; // [mm] y-coordinate of center of ellipse + + Short_t fNumUsedPixels; // Number of pixels which survived the image cleaning + Short_t fNumCorePixels; // number of core pixels + + Float_t fSinDelta; //! [1] sin of Delta (to be used in derived classes) + Float_t fCosDelta; //! [1] cos of Delta (to be used in derived classes) + + TEllipse *fEllipse; //! Graphical Object to Display Ellipse + +protected: + // + // This is only for calculations in derived classes because + // we don't want to read/write this data members + // + Float_t GetCosDelta() const { return fCosDelta; } + Float_t GetSinDelta() const { return fSinDelta; } + +public: + MHillas(const char *name=NULL, const char *title=NULL); + ~MHillas(); + + void Reset(); + + virtual Int_t Calc(const MGeomCam &geom, const MCerPhotEvt &pix); + + virtual void Print(Option_t *opt=NULL) const; + virtual void Draw(Option_t *opt=NULL); + //virtual void Paint(Option_t *); + + virtual void Clear(Option_t *opt=NULL); + + Float_t GetLength() const { return fLength; } + Float_t GetWidth() const { return fWidth; } + Float_t GetDelta() const { return fDelta; } + Float_t GetSize() const { return fSize; } + Float_t GetMeanX() const { return fMeanX; } + Float_t GetMeanY() const { return fMeanY; } + + Int_t GetNumUsedPixels() const { return fNumUsedPixels; } + Int_t GetNumCorePixels() const { return fNumCorePixels; } + + virtual void Set(const TArrayF &arr); + + //virtual void AsciiRead(ifstream &fin); + //virtual void AsciiWrite(ofstream &fout) const; + + ClassDef(MHillas, 2) // Storage Container for Hillas Parameter +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHillasCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasCalc.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasCalc.cc (revision 1940) @@ -0,0 +1,136 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de) +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MHillasCalc // +// // +// This is a task to calculate the Hillas parameters from each event // +// // +// Input Containers: // +// MCerPhotEvt, MGeomCam // +// // +// Output Containers: // +// MHillas // +// // +///////////////////////////////////////////////////////////////////////////// + +#include "MHillasCalc.h" + +#include "MParList.h" + +#include "MHillas.h" +#include "MCerPhotEvt.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MHillasCalc); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MHillasCalc::MHillasCalc(const char *hil, const char *name, const char *title) +{ + fName = name ? name : "MHillasCalc"; + fTitle = title ? title : "Task to calculate Hillas parameters"; + + fHilName = hil; +} + +// -------------------------------------------------------------------------- +// +// Check for a MCerPhotEvt object from which the Hillas are calculated. +// Try to find the Geometry conatiner. And try to find the output +// (Hillas) container or create one. +// +Bool_t MHillasCalc::PreProcess(MParList *pList) +{ + fCerPhotEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fCerPhotEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fGeomCam) + { + *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; + return kFALSE; + } + + fHillas = (MHillas*)pList->FindCreateObj("MHillas",fHilName); + if (!fHillas) + return kFALSE; + + memset(fErrors, 0, sizeof(fErrors)); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// If you want do complex descisions inside the calculations +// we must move the calculation code inside this function +// +// If the calculation wasn't sucessfull skip this event +// +Bool_t MHillasCalc::Process() +{ + const Int_t rc = fHillas->Calc(*fGeomCam, *fCerPhotEvt); + if (rc<0 || rc>4) + { + *fLog << err << dbginf << "MHillas::Calc returned unknown error code!" << endl; + return kFALSE; + } + fErrors[rc]++; + + return rc==0 ? kTRUE : kCONTINUE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the hillas calculation. The percentage +// is calculated with respect to the number of executions of this task. +// +Bool_t MHillasCalc::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << setw(7) << fErrors[1] << " (" << setw(3) << (int)(fErrors[1]*100/GetNumExecutions()) << "%) Evts skipped due to: Event has less than 3 pixels" << endl; + *fLog << " " << setw(7) << fErrors[2] << " (" << setw(3) << (int)(fErrors[2]*100/GetNumExecutions()) << "%) Evts skipped due to: Calculated Size == 0" << endl; + *fLog << " " << setw(7) << fErrors[3] << " (" << setw(3) << (int)(fErrors[3]*100/GetNumExecutions()) << "%) Evts skipped due to: Number of used pixels < 3" << endl; + *fLog << " " << setw(7) << fErrors[4] << " (" << setw(3) << (int)(fErrors[4]*100/GetNumExecutions()) << "%) Evts skipped due to: CorrXY==0" << endl; + *fLog << " " << fErrors[0] << " (" << (int)(fErrors[0]*100/GetNumExecutions()) << "%) Evts survived Hillas calculation!" << endl; + *fLog << endl; + + return kTRUE; +} Index: trunk/MagicSoft/Mars/mimage/MHillasCalc.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasCalc.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasCalc.h (revision 1940) @@ -0,0 +1,39 @@ +#ifndef MARS_MHillasCalc +#define MARS_MHillasCalc + +///////////////////////////////////////////////////////////////////////////// +// // +// MHillasCalkc // +// // +// Task to calculate Hillas Parameters // +// // +///////////////////////////////////////////////////////////////////////////// + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MGeomCam; +class MCerPhotEvt; +class MHillas; + +class MHillasCalc : public MTask +{ + const MGeomCam *fGeomCam; // Camera Geometry used to calculate Hillas + const MCerPhotEvt *fCerPhotEvt; // Cerenkov Photon Event used for calculation + MHillas *fHillas; // ouput container to store result + + TString fHilName; + Int_t fErrors[5]; + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + Bool_t PostProcess(); + +public: + MHillasCalc(const char *hil="MHillas", const char *name=NULL, const char *title=NULL); + + ClassDef(MHillasCalc, 0) // Task to calculate Hillas parameters +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHillasExt.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasExt.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasExt.cc (revision 1940) @@ -0,0 +1,291 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Rudolf Bock 10/2001 +! Author(s): Wolfgang Wittek 06/2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHillasExt +// +// Storage Container for extended image parameters +// +// extended image parameters +// fConc ratio of sum of two highest pixels over fSize +// fConc1 ratio of highest pixel over fSize +// fAsym distance from highest pixel to center, projected onto major axis +// fM3Long third moment along major axis +// fM3Trans third moment along minor axis +// +// WARNING: Before you can use fAsym, fM3Long and fM3Trans you must +// multiply by the sign of MHillasSrc::fCosDeltaAlpha +// +//////////////////////////////////////////////////////////////////////////// +/* + // fLeakage1 ratio : (photons in most outer ring of pixels) over fSize + // fLeakage2 ratio : (photons in the 2 outer rings of pixels) over fSize + // + // fAsymna d/(d na) of ( sum(x*q^na)/sum(q^na), sum(y*q^na)/sum(q^na) ) + // projected onto the major axis + // fAsym0 (F-B)/(F+B) along the major axis + */ +#include "MHillasExt.h" + +#include +#include + +#include "MGeomPix.h" +#include "MGeomCam.h" + +#include "MCerPhotPix.h" +#include "MCerPhotEvt.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MHillasExt); + +// ------------------------------------------------------------------------- +// +// Default constructor. +// +MHillasExt::MHillasExt(const char *name, const char *title) +{ + fName = name ? name : "MHillas"; + fTitle = title ? title : "Storage container for extended parameter set of one event"; + + Reset(); + // FIXME: (intelligent) initialization of values missing +} + +// ------------------------------------------------------------------------- +// +void MHillasExt::Reset() +{ + MHillas::Reset(); + + fConc = -1; + fConc1 = -1; + fAsym = 0; + fM3Long = 0; + fM3Trans = 0; +} + +// ------------------------------------------------------------------------- +// +void MHillasExt::Print(Option_t *) const +{ + MHillas::Print(); + + *fLog << "Extended Image Parameters (" << GetName() << ")" << endl; + *fLog << " - Conc = " << fConc << " (ratio)" << endl; + *fLog << " - Conc1 = " << fConc1 << " (ratio)" << endl; + *fLog << " - Asymmetry = " << fAsym << " mm" << endl; + *fLog << " - 3rd Moment Long = " << fM3Long << " mm" << endl; + *fLog << " - 3rd Moment Trans = " << fM3Trans << " mm" << endl; +} + +// ------------------------------------------------------------------------- +// +// calculation of additional parameters based on the camera geometry +// and the cerenkov photon event +// +Int_t MHillasExt::Calc(const MGeomCam &geom, const MCerPhotEvt &evt) +{ + const Int_t rc = MHillas::Calc(geom, evt); + if (rc>0) + return rc; + + // + // calculate the additional image parameters + // -------------------------------------------- + // + // loop to get third moments along ellipse axes and two max pixels + // + // For the moments double precision is used to make sure, that + // the complex matrix multiplication and sum is evaluated correctly. + // + Double_t m3x = 0; + Double_t m3y = 0; + + Float_t maxpix1 = 0; // [#phot] + Float_t maxpix2 = 0; // [#phot] + + Int_t maxpixid = 0; + + const UInt_t npixevt = evt.GetNumPixels(); + + const Float_t A0 = geom[0].GetA(); + + for (UInt_t i=0; imaxpix1) + { + maxpix2 = maxpix1; + maxpix1 = nphot; // [1] + maxpixid = pix.GetPixId(); + continue; // [1] + } + + if (nphot>maxpix2) + maxpix2 = nphot; // [1] + + /* + // + // power na for calculating fAsymna; + // the value 1.5 was suggested by Thomas Schweizer + // + Double_t na = 1.5; + + // + // get sums for calculating fAsymna + // the outer pixels are 4 times as big (in area) + // as the inner pixels ! + // + const Double_t dummy = pow(nphot, na)/r; + + sna += dummy; + xna += dzx*dummy; + + sna1 += sna/nphot; + xna1 += xna/nphot; + + // + // forward-backward asymmetry + // + fb += dzx<0 ? -nphot: nphot; + */ + } + + const MGeomPix &maxpix = geom[maxpixid]; + + fAsym = (GetMeanX()-maxpix.GetX())*GetCosDelta() + + (GetMeanY()-maxpix.GetY())*GetSinDelta(); // [mm] + + fConc = (maxpix1+maxpix2)/GetSize(); // [ratio] + fConc1 = maxpix1/GetSize(); // [ratio] + + /* + fLeakage1 = edgepix1 / GetSize(); + fLeakage2 = edgepix2 / GetSize(); + fAsym0 = fb / GetSize(); + + fAsymna = na * (sna*xna1 - sna1*xna) / (sna*sna); + */ + + // + // Third moments along axes get normalized + // + m3x /= GetSize(); + m3y /= GetSize(); + + fM3Long = m3x<0 ? -pow(-m3x, 1./3) : pow(m3x, 1./3); // [mm] + fM3Trans = m3y<0 ? -pow(-m3y, 1./3) : pow(m3y, 1./3); // [mm] + + SetReadyToSave(); + + return 0; +} + +// -------------------------------------------------------------------------- +// +// This function is ment for special usage, please never try to set +// values via this function +// +void MHillasExt::Set(const TArrayF &arr) +{ + if (arr.GetSize() != 13) + return; + + fConc = arr.At(8); // [ratio] concentration ratio: sum of the two highest pixels / fSize + fConc1 = arr.At(9); // [ratio] concentration ratio: sum of the highest pixel / fSize + fAsym = arr.At(10); // [mm] fDist minus dist: center of ellipse, highest pixel + fM3Long = arr.At(11); // [mm] 3rd moment (e-weighted) along major axis + fM3Trans = arr.At(12); // [mm] 3rd moment (e-weighted) along minor axis + + TArrayF n(arr); + n.Set(8); + MHillas::Set(n); +} + +/* +// ------------------------------------------------------------------------- +// +void MHillasExt::AsciiRead(ifstream &fin) +{ + MHillas::AsciiRead(fin); + + fin >> fConc; + fin >> fConc1; + fin >> fAsym; + fin >> fM3Long; + fin >> fM3Trans; +} +*/ +// ------------------------------------------------------------------------- +/* +void MHillasExt::AsciiWrite(ofstream &fout) const +{ + MHillas::AsciiWrite(fout); + + fout << " "; + fout << fConc << " "; + fout << fConc1 << " "; + fout << fAsym << " "; + fout << fM3Long << " "; + fout << fM3Trans; +} +*/ Index: trunk/MagicSoft/Mars/mimage/MHillasExt.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasExt.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasExt.h (revision 1940) @@ -0,0 +1,43 @@ +#ifndef MARS_MHillasExt +#define MARS_MHillasExt + +#ifndef MARS_MHillas +#include "MHillas.h" +#endif + +class MGeomCam; +class MCerPhotEvt; + +class MHillasExt : public MHillas +{ +private: + // for description see MExtHillas.cc + Float_t fConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize + Float_t fConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize + Float_t fAsym; // [mm] fDist minus dist: center of ellipse, highest pixel + Float_t fM3Long; // [mm] 3rd moment (e-weighted) along major axis + Float_t fM3Trans; // [mm] 3rd moment (e-weighted) along minor axis + +public: + MHillasExt(const char *name=NULL, const char *title=NULL); + + void Reset(); + + Float_t GetConc() const { return fConc; } + Float_t GetConc1() const { return fConc1; } + Float_t GetAsym() const { return fAsym; } + Float_t GetM3Long() const { return fM3Long; } + Float_t GetM3Trans() const { return fM3Trans; } + + Int_t Calc(const MGeomCam &geom, const MCerPhotEvt &pix); + + void Print(Option_t *opt=NULL) const; + + void Set(const TArrayF &arr); + + //void AsciiRead(ifstream &fin); + //void AsciiWrite(ofstream &fout) const; + + ClassDef(MHillasExt, 1) // Storage Container for extended Hillas Parameter +}; +#endif Index: trunk/MagicSoft/Mars/mimage/MHillasSrc.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasSrc.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasSrc.cc (revision 1940) @@ -0,0 +1,181 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Harald Kornmayer 1/2001 +! Author(s): Rudolf Bock 10/2001 +! Author(s): Wolfgang Wittek 06/2002 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MHillasSrc +// +// Storage Container for image parameters +// +// source-dependent image parameters +// +// Version 1: +// ---------- +// fAlpha angle between major axis and line source-to-center +// fDist distance from source to center of ellipse +// +// Version 2: +// ---------- +// fHeadTail +// +// Version 3: +// ---------- +// fCosDeltaAlpha cosine of angle between d and a, where +// - d is the vector from the source position to the +// center of the ellipse +// - a is a vector along the main axis of the ellipse, +// defined with positive x-component +// +///////////////////////////////////////////////////////////////////////////// +#include "MHillasSrc.h" + +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MSrcPosCam.h" + +ClassImp(MHillasSrc); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MHillasSrc::MHillasSrc(const char *name, const char *title) +{ + fName = name ? name : "MHillasSrc"; + fTitle = title ? title : "Parameters depending in source position"; +} + +void MHillasSrc::Reset() +{ + fDist = -1; + fAlpha = 0; + fHeadTail = 0; + fCosDeltaAlpha = 0; +} + +// -------------------------------------------------------------------------- +// +// Calculation of source-dependent parameters +// In case you don't call Calc from within an eventloop make sure, that +// you call the Reset member function before. +// +Bool_t MHillasSrc::Calc(const MHillas *hillas) +{ + fHillas = hillas; + + const Double_t mx = GetMeanX(); // [mm] + const Double_t my = GetMeanY(); // [mm] + + const Double_t sx = mx - fSrcPos->GetX(); // [mm] + const Double_t sy = my - fSrcPos->GetY(); // [mm] + + const Double_t sd = sin(GetDelta()); // [1] + const Double_t cd = cos(GetDelta()); // [1] + + const Double_t tand = tan(GetDelta()); // [1] + + fHeadTail = cd*sx + sd*sy; // [mm] + + // + // Distance from source position to center of ellipse. + // If the distance is 0 distance, Alpha is not specified. + // The calculation has failed and returnes kFALSE. + // + Double_t dist = sqrt(sx*sx + sy*sy); // [mm] + + if (dist==0) + { + //*fLog << warn << GetDescriptor() << ": Event has Dist==0... skipped." << endl; + return kFALSE; + } + + // + // Calculate Alpha and Cosda = cos(d,a) + // The sign of Cosda will be used for quantities containing + // a head-tail information + // + const Double_t arg = (sy-tand*sx) / (dist*sqrt(tand*tand+1)); + + fAlpha = asin(arg)*kRad2Deg; // [deg] + fCosDeltaAlpha = fHeadTail/dist; // [1] + fDist = dist; // [mm] + + SetReadyToSave(); + + return kTRUE; +} + +void MHillasSrc::Print(Option_t *) const +{ + *fLog << all; + *fLog << "Source dependant Image Parameters (" << GetName() << ")" << endl; + *fLog << " - Dist [mm] = " << fDist << endl; + *fLog << " - Alpha [deg] = " << fAlpha << endl; + *fLog << " - HeadTail [mm] = " << fHeadTail << endl; + *fLog << " - CosDeltaAlpha = " << fCosDeltaAlpha << endl; +} +// -------------------------------------------------------------------------- +// +// This function is ment for special usage, please never try to set +// values via this function +// +void MHillasSrc::Set(const TArrayF &arr) +{ + if (arr.GetSize() != 4) + return; + + fAlpha = arr.At(0); // [deg] angle of major axis with vector to src + fDist = arr.At(1); // [mm] distance between src and center of ellipse + fHeadTail = arr.At(2); // [mm] + fCosDeltaAlpha = arr.At(3); // [1] cosine of angle between d and a +} + +// ----------------------------------------------------------------------- +// +// overloaded MParContainer to read MHillasSrc from an ascii file +// +/* +void MHillasSrc::AsciiRead(ifstream &fin) +{ + fin >> fAlpha; + fin >> fDist; + fin >> fHeadTail; +} +*/ +// ----------------------------------------------------------------------- +// +// overloaded MParContainer to write MHillasSrc to an ascii file +/* +void MHillasSrc::AsciiWrite(ofstream &fout) const +{ + fout << fAlpha << " " << fDist; +} +*/ Index: trunk/MagicSoft/Mars/mimage/MHillasSrc.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasSrc.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasSrc.h (revision 1940) @@ -0,0 +1,52 @@ +#ifndef MARS_MHillasSrc +#define MARS_MHillasSrc + +#ifndef MARS_MHillas +#include "MHillas.h" +#endif + +class MSrcPosCam; + +class MHillasSrc : public MParContainer +{ +private: + const MHillas *fHillas; //! Input parameters + const MSrcPosCam *fSrcPos; //! Source position in the camere + + Float_t fAlpha; // [deg] angle of major axis with vector to src + Float_t fDist; // [mm] distance between src and center of ellipse + Float_t fHeadTail; // [mm] + Float_t fCosDeltaAlpha; // [1] cosine of angle between d and a + +public: + MHillasSrc(const char *name=NULL, const char *title=NULL); + + void SetSrcPos(const MSrcPosCam *pos) { fSrcPos = pos; } + const MSrcPosCam *GetSrcPos() const { return fSrcPos; } + + void Reset(); + + Float_t GetLength() const { return fHillas->GetLength(); } + Float_t GetWidth() const { return fHillas->GetWidth(); } + Float_t GetDelta() const { return fHillas->GetDelta(); } + Float_t GetSize() const { return fHillas->GetSize(); } + Float_t GetMeanX() const { return fHillas->GetMeanX(); } + Float_t GetMeanY() const { return fHillas->GetMeanY(); } + Float_t GetAlpha() const { return fAlpha; } + Float_t GetDist() const { return fDist; } + Float_t GetHeadTail() const { return fHeadTail; } + Float_t GetCosDeltaAlpha() const { return fCosDeltaAlpha; } + + void Print(Option_t *opt=NULL) const; + + virtual Bool_t Calc(const MHillas *hillas); + + void Set(const TArrayF &arr); + + //virtual void AsciiRead(ifstream &fin); + //virtual void AsciiWrite(ofstream &fout) const; + + ClassDef(MHillasSrc, 3) // Container to hold source position dependant parameters +}; + +#endif Index: trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.cc (revision 1940) @@ -0,0 +1,183 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Rudolf Bock 10/2001 +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MHillasSrcCalc +// +// Task to calculate the source dependant part of the hillas parameters +// +////////////////////////////////////////////////////////////////////////////// +#include "MHillasSrcCalc.h" + +#include + +#include "MParList.h" + +#include "MSrcPosCam.h" +#include "MHillasSrc.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MHillasSrcCalc); + +static const TString gsDefName = "MHillasSrcCalc"; +static const TString gsDefTitle = "Calculate position dependant image parameters"; + +// ------------------------------------------------------------------------- +// +// Default constructor. The first argument is the name of a container +// containing the source position in the camera plain (MScrPosCam). +// The default is "MSrcPosCam". hil is the name of a container +// of type MHillasSrc (or derived) in which the parameters are stored +// The default is "MHillasSrc" +// +//MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, +// const char *name, const char *title) +// : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) +//{ +// fName = name ? name : gsDefName.Data(); +// fTitle = title ? title : gsDefTitle.Data(); + +// fSrcName = src; +// fHillasName = hil; +// fHillasInput = "MHillas"; +//} +// ------------------------------------------------------------------------- +// +MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, + const char *name, const char *title) + : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) +{ + fName = name ? name : gsDefName.Data(); + fTitle = title ? title : gsDefTitle.Data(); + + fSrcName = src; + fHillasName = hil; + fHillasInput = "MHillas"; +} + +// ------------------------------------------------------------------------- +// +Bool_t MHillasSrcCalc::PreProcess(MParList *pList) +{ + fHillas = (MHillas*)pList->FindObject(fHillasInput, "MHillas"); + if (!fHillas) + { + *fLog << err << dbginf << "MHillas not found... aborting." << endl; + return kFALSE; + } + + fSrcPos = (MSrcPosCam*)pList->FindObject(fSrcName, "MSrcPosCam"); + if (!fSrcPos) + { + *fLog << err << dbginf << fSrcName << " [MSrcPosCam] not found... aborting." << endl; + return kFALSE; + } + + fHillasSrc = (MHillasSrc*)pList->FindCreateObj("MHillasSrc", fHillasName); + if (!fHillasSrc) + return kFALSE; + + fHillasSrc->SetSrcPos(fSrcPos); + + fErrors = 0; + + return kTRUE; +} + +// ------------------------------------------------------------------------- +// +Bool_t MHillasSrcCalc::Process() +{ + + if (!fHillasSrc->Calc(fHillas)) + { + fErrors++; + return kCONTINUE; + + } + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the hillas calculation. The percentage +// is calculated with respect to the number of executions of this task. +// +Bool_t MHillasSrcCalc::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << fErrors << " (" << (int)(fErrors*100/GetNumExecutions()) << "%) Evts skipped due to: Dist==0" << endl; + *fLog << endl; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Implementation of SavePrimitive. Used to write the call to a constructor +// to a macro. In the original root implementation it is used to write +// gui elements to a macro-file. +// +void MHillasSrcCalc::StreamPrimitive(ofstream &out) const +{ + if (fHillas) + fHillas->SavePrimitive(out); + + if (fSrcPos) + fSrcPos->SavePrimitive(out); + + if (fHillasSrc) + fHillasSrc->SavePrimitive(out); + + out << " MHillasSrcCalc " << GetUniqueName() << "("; + + if (fSrcPos) + out << "&" << fSrcPos->GetUniqueName(); + else + out << "\"" << fSrcName << "\""; + + out << ", "; + + if (fHillasSrc) + out << "&" << fHillasSrc->GetUniqueName(); + else + out << "\"" << fHillasName << "\""; + + if (fName!=gsDefName || fTitle!=gsDefTitle) + { + out << ", \"" << fName << "\""; + if (fTitle!=gsDefTitle) + out << ", \"" << fTitle << "\""; + } + out << ");" << endl; +} Index: trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MHillasSrcCalc.h (revision 1940) @@ -0,0 +1,44 @@ +#ifndef MARS_MHillasSrcCalc +#define MARS_MHillasSrcCalc + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MHillas; +class MHillasSrc; +class MSrcPosCam; + +class MHillasSrcCalc : public MTask +{ +private: + MHillas *fHillas; //! Pointer to the source independant hillas parameters + MSrcPosCam *fSrcPos; //! Pointer to the source position + MHillasSrc *fHillasSrc; //! Pointer to the output container for the source dependant parameters + + TString fSrcName; + TString fHillasName; + TString fHillasInput; + + Int_t fErrors; + + void StreamPrimitive(ofstream &out) const; + + Bool_t PreProcess(MParList *plist); + Bool_t Process(); + Bool_t PostProcess(); + +public: + // MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", + // const char *name=NULL, const char *title=NULL); + + MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", + const char *name=NULL, const char *title=NULL); + + void SetInput(TString hilname) { fHillasInput = hilname; } + + ClassDef(MHillasSrcCalc, 1) // task to calculate the source position depandant hillas parameters +}; + +#endif + Index: trunk/MagicSoft/Mars/mimage/MImgCleanStd.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MImgCleanStd.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MImgCleanStd.cc (revision 1940) @@ -0,0 +1,512 @@ +/* ======================================================================== *\ +! +! * +! * 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): Thomas Bretz 12/2000 +! Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de) +! +! Copyright: MAGIC Software Development, 2000-2002 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MImgCleanStd // +// // +// This is the standard image cleaning. If you want to know how it works // +// Please look at the three CleanSteps and Process // +// // +// FIXME: MImgCleanStd is not yet completely optimized for speed. // +// Maybe we don't have to loop over all pixels all the time... // +// // +// Input Containers: // +// MGeomCam, MCerPhotEvt // +// // +// Output Containers: // +// -/- // +// // +///////////////////////////////////////////////////////////////////////////// +#include "MImgCleanStd.h" + +#include // atof +#include // ofstream, SavePrimitive + +#include // TGFrame +#include // TGLabel +#include // TGTextEntry + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" +#include "MGeomPix.h" +#include "MGeomCam.h" +#include "MCerPhotPix.h" +#include "MCerPhotEvt.h" + +#include "MGGroupFrame.h" // MGGroupFrame + +ClassImp(MImgCleanStd); + +enum { + kImgCleanLvl1, + kImgCleanLvl2 +}; + +static const TString gsDefName = "MImgCleanStd"; +static const TString gsDefTitle = "Task to perform a standard image cleaning"; + +// -------------------------------------------------------------------------- +// +// Default constructor. Here you can specify the cleaning levels. If you +// don't specify them the 'common standard' values 3.0 and 2.5 (sigma +// above mean) are used +// +MImgCleanStd::MImgCleanStd(const Float_t lvl1, const Float_t lvl2, + const char *name, const char *title) + : fCleanLvl1(lvl1), fCleanLvl2(lvl2) +{ + fName = name ? name : gsDefName.Data(); + fTitle = title ? title : gsDefTitle.Data(); + + Print(); +} + +// -------------------------------------------------------------------------- +// +// This method looks for all pixels with an entry (photons) +// that is three times bigger than the noise of the pixel +// (std: 3 sigma, clean level 1) +// +// +// AM 18/11/2002: now cut levels are proportional to the square root +// of the pixel area. In this way the cut corresponds to a fixed +// phe-density (otherwise, it would bias the images). +// +// Returns the maximum Pixel Id (used for ispixused in CleanStep2) +// +Int_t MImgCleanStd::CleanStep1() +{ + const Int_t entries = fEvt->GetNumPixels(); + + Int_t max = entries; + + // + // check the number of all pixels against the noise level and + // set them to 'unused' state if necessary + // + for (Int_t i=0; iGetPixRatio(id); + + // COBB: '<=' to skip entry=noise=0 + if (entry <= fCleanLvl1 * noise / ratio) + pix.SetPixelUnused(); + + if (id>max) + max = id; + } + return max; +} + +// -------------------------------------------------------------------------- +// +// check if the survived pixel have a neighbor, that also +// survived, otherwise set pixel to unused. (removes pixels without +// neighbors) +// +// takes the maximum pixel id from CleanStep1 as an argument +// +void MImgCleanStd::CleanStep2(Int_t max) +{ + const Int_t entries = fEvt->GetNumPixels(); + + // + // In the worst case we have to loop 6 times 577 times, to + // catch the behaviour of all next neighbors. Here we can gain + // much by using an array instead of checking through all pixels + // (MCerPhotEvt::IsPixelUsed) all the time. + // + Byte_t *ispixused = new Byte_t[max+1]; + memset(ispixused, 0, max+1); + + for (Int_t i=0; imax || !ispixused[id2]) + continue; + + cnt = kTRUE; + break; + } + if (cnt) + continue; +#else + Int_t cnt = 0; + for (Int_t j=0; jmax || !ispixused[id2]) + continue; + + if (cnt++>nnmax-4) + break; + } + if (cnt==nnmax-2 && nnmax>=4) + { + pix.SetPixelUsed(); + continue; + } + if (cnt>0) + continue; +#endif + + // + // check if no next neighbor has the state 'used' + // set this pixel to 'unused', too. + // + pix.SetPixelUnused(); + } + + delete ispixused; + + // + // now we declare all pixels that survive as CorePixels + // + for (Int_t i=0; iGetNumPixels(); + + for (Int_t i=0; iGetPixRatio(id)); + + if (entry <= fCleanLvl2 * noise / ratio) + continue; + + // + // check if the pixel's next neighbor is a core pixel. + // if it is a core pixel set pixel state to: used. + // + MGeomPix &gpix = (*fCam)[id]; + const Int_t nnmax = gpix.GetNumNeighbors(); + + for (Int_t j=0; jIsPixelCore(id2)) + continue; + + pix.SetPixelUsed(); + break; + } + } +} + +// -------------------------------------------------------------------------- +// +// check if MEvtHeader exists in the Parameter list already. +// if not create one and add them to the list +// +Bool_t MImgCleanStd::PreProcess (MParList *pList) +{ + fCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fCam) + { + *fLog << dbginf << "MGeomCam not found (no geometry information available)... aborting." << endl; + return kFALSE; + } + + fEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + return kTRUE; +} + + +// -------------------------------------------------------------------------- +// +// Cleans the image. +// +Bool_t MImgCleanStd::Process() +{ + const Int_t max = CleanStep1(); + CleanStep2(max); + CleanStep3(); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Print descriptor and cleaning levels. +// +void MImgCleanStd::Print(Option_t *o) const +{ + *fLog << GetDescriptor() << " initialized with noise level "; + *fLog << fCleanLvl1 << " and " << fCleanLvl2 << endl; +} + +// -------------------------------------------------------------------------- +// +// Craete two text entry fields, one for each cleaning level and a +// describing text line. +// +void MImgCleanStd::CreateGuiElements(MGGroupFrame *f) +{ + // + // Create a frame for line 3 and 4 to be able + // to align entry field and label in one line + // + TGHorizontalFrame *f1 = new TGHorizontalFrame(f, 0, 0); + TGHorizontalFrame *f2 = new TGHorizontalFrame(f, 0, 0); + + /* + * --> use with root >=3.02 <-- + * + + TGNumberEntry *fNumEntry1 = new TGNumberEntry(frame, 3.0, 2, M_NENT_LVL1, kNESRealOne, kNEANonNegative); + TGNumberEntry *fNumEntry2 = new TGNumberEntry(frame, 2.5, 2, M_NENT_LVL1, kNESRealOne, kNEANonNegative); + + */ + TGTextEntry *entry1 = new TGTextEntry(f1, "****", kImgCleanLvl1); + TGTextEntry *entry2 = new TGTextEntry(f2, "****", kImgCleanLvl2); + + // --- doesn't work like expected (until root 3.02?) --- fNumEntry1->SetAlignment(kTextRight); + // --- doesn't work like expected (until root 3.02?) --- fNumEntry2->SetAlignment(kTextRight); + + entry1->SetText("3.0"); + entry2->SetText("2.5"); + + entry1->Associate(f); + entry2->Associate(f); + + TGLabel *l1 = new TGLabel(f1, "Cleaning Level 1"); + TGLabel *l2 = new TGLabel(f2, "Cleaning Level 2"); + + l1->SetTextJustify(kTextLeft); + l2->SetTextJustify(kTextLeft); + + // + // Align the text of the label centered, left in the row + // with a left padding of 10 + // + TGLayoutHints *laylabel = new TGLayoutHints(kLHintsCenterY|kLHintsLeft, 10); + TGLayoutHints *layframe = new TGLayoutHints(kLHintsCenterY|kLHintsLeft, 5, 0, 10); + + // + // Add one entry field and the corresponding label to each line + // + f1->AddFrame(entry1); + f2->AddFrame(entry2); + + f1->AddFrame(l1, laylabel); + f2->AddFrame(l2, laylabel); + + f->AddFrame(f1, layframe); + f->AddFrame(f2, layframe); + + f->AddToList(entry1); + f->AddToList(entry2); + f->AddToList(l1); + f->AddToList(l2); + f->AddToList(laylabel); + f->AddToList(layframe); +} + +void MImgCleanStd::SetLvl1(Float_t lvl) +{ + fCleanLvl1 = lvl; + *fLog << inf << "Cleaning level 1 set to " << lvl << " sigma." << endl; +} + +void MImgCleanStd::SetLvl2(Float_t lvl) +{ + fCleanLvl2 = lvl; + *fLog << inf << "Cleaning level 2 set to " << lvl << " sigma." << endl; +} + +// -------------------------------------------------------------------------- +// +// Process the GUI Events comming from the two text entry fields. +// +Bool_t MImgCleanStd::ProcessMessage(Int_t msg, Int_t submsg, Long_t param1, Long_t param2) +{ + if (msg!=kC_TEXTENTRY || submsg!=kTE_ENTER) + return kTRUE; + + TGTextEntry *txt = (TGTextEntry*)FindWidget(param1); + + if (!txt) + return kTRUE; + + Float_t lvl = atof(txt->GetText()); + + switch (param1) + { + case kImgCleanLvl1: + SetLvl1(lvl); + return kTRUE; + + case kImgCleanLvl2: + SetLvl2(lvl); + return kTRUE; + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Implementation of SavePrimitive. Used to write the call to a constructor +// to a macro. In the original root implementation it is used to write +// gui elements to a macro-file. +// +void MImgCleanStd::StreamPrimitive(ofstream &out) const +{ + out << " MImgCleanStd " << GetUniqueName() << "("; + out << fCleanLvl1 << ", " << fCleanLvl2; + + if (fName!=gsDefName || fTitle!=gsDefTitle) + { + out << ", \"" << fName << "\""; + if (fTitle!=gsDefTitle) + out << ", \"" << fTitle << "\""; + } + out << ");" << endl; +} + +// -------------------------------------------------------------------------- +// +// Read the setup from a TEnv: +// Float_t fCleanLvl1: CleaningLevel1 +// Float_t fCleanLvl2: CleaningLevel2 +// +Bool_t MImgCleanStd::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + Bool_t rc = kTRUE; + if (!IsEnvDefined(env, prefix, "CleaningLevel1", print)) + rc = kFALSE; + else + { + SetLvl1(GetEnvValue(env, prefix, "CleaningLevel1", fCleanLvl1)); + if (fCleanLvl1<0) + { + *fLog << err << "ERROR - Negative values for Cleaning Level 1 forbidden." << endl; + return kERROR; + } + } + + if (!IsEnvDefined(env, prefix, "CleaningLevel2", print)) + rc = kFALSE; + else + { + SetLvl2(GetEnvValue(env, prefix, "CleaningLevel2", fCleanLvl2)); + if (fCleanLvl2<0) + { + *fLog << err << "ERROR - Negative values for Cleaning Level 2 forbidden." << endl; + return kERROR; + } + } + + return rc; +} Index: trunk/MagicSoft/Mars/mimage/MImgCleanStd.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MImgCleanStd.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MImgCleanStd.h (revision 1940) @@ -0,0 +1,51 @@ +#ifndef MARS_MImgCleanStd +#define MARS_MImgCleanStd + +#ifndef MARS_MGTask +#include "MGTask.h" +#endif + +class MGeomCam; +class MCerPhotEvt; + +class MGGroupFrame; + +class MImgCleanStd : public MGTask +{ +private: + const MGeomCam *fCam; //! + MCerPhotEvt *fEvt; //! + + Float_t fCleanLvl1; + Float_t fCleanLvl2; + + void SetLvl1(Float_t lvl); + void SetLvl2(Float_t lvl); + + void CreateGuiElements(MGGroupFrame *f); + void StreamPrimitive(ofstream &out) const; + +public: + MImgCleanStd(const Float_t lvl1=3.0, const Float_t lvl2=2.5, + const char *name=NULL, const char *title=NULL); + + Int_t CleanStep1(); + void CleanStep2(Int_t max); + void CleanStep3(); + + Bool_t PreProcess(MParList *pList); + Bool_t Process(); + + void Print(Option_t *o="") const; + + Float_t GetCleanLvl1() const { return fCleanLvl1; } + Float_t GetCleanLvl2() const { return fCleanLvl2; } + + Bool_t ProcessMessage(Int_t msg, Int_t submsg, Long_t param1, Long_t param2); + Bool_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); + + ClassDef(MImgCleanStd, 0) // task doing a standard image cleaning +}; + +#endif + Index: trunk/MagicSoft/Mars/mimage/MNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MNewImagePar.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MNewImagePar.cc (revision 1940) @@ -0,0 +1,114 @@ +/* ======================================================================== *\ +! +! * +! * 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 03/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MNewImagePar +// +// Storage Container for new image parameters +// +///////////////////////////////////////////////////////////////////////////// +#include "MNewImagePar.h" + +#include +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MGeomCam.h" +#include "MGeomPix.h" +#include "MCerPhotEvt.h" +#include "MCerPhotPix.h" +#include "MSrcPosCam.h" + +ClassImp(MNewImagePar); + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MNewImagePar::MNewImagePar(const char *name, const char *title) +{ + fName = name ? name : "MNewImagePar"; + fTitle = title ? title : "New image parameters"; +} + +// -------------------------------------------------------------------------- +// +void MNewImagePar::Reset() +{ + fLeakage1 = 0; + fLeakage2 = 0; +} + +// -------------------------------------------------------------------------- +// +// Calculation of new image parameters +// +// +Bool_t MNewImagePar::Calc(const MGeomCam &geom, const MCerPhotEvt &evt, + const MHillas *hillas) +{ + //fHillas = (MHillas*)hillas; + + const UInt_t npixevt = evt.GetNumPixels(); + + Double_t edgepix1 = 0.0; + Double_t edgepix2 = 0.0; + + for (UInt_t i=0; iGetSize(); + fLeakage2 = edgepix2 / hillas->GetSize(); + + SetReadyToSave(); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +void MNewImagePar::Print(Option_t *) const +{ + *fLog << all; + *fLog << "New Image Parameters (" << GetName() << ")" << endl; + *fLog << " - Leakage1 = " << fLeakage1 << endl; + *fLog << " - Leakage2 = " << fLeakage2 << endl; +} Index: trunk/MagicSoft/Mars/mimage/MNewImagePar.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MNewImagePar.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MNewImagePar.h (revision 1940) @@ -0,0 +1,65 @@ +#ifndef MARS_MNewImagePar +#define MARS_MNewImagePar + +#ifndef MARS_MHillas +#include "MHillas.h" +#endif + +class MSrcPosCam; + +class MNewImagePar : public MParContainer +{ +private: + Float_t fLeakage1; // (photons in most outer ring of pixels) over fSize + Float_t fLeakage2; // (photons in the 2 outer rings of pixels) over fSize + +public: + MNewImagePar(const char *name=NULL, const char *title=NULL); + + // void SetSrcPos(MSrcPosCam *pos) { fSrcPos = pos; } + // const MSrcPosCam *GetSrcPos() const { return fSrcPos; } + + void Reset(); + + Float_t GetLeakage1() const { return fLeakage1; } + Float_t GetLeakage2() const { return fLeakage2; } + + void Print(Option_t *opt=NULL) const; + + virtual Bool_t Calc(const MGeomCam &geom, const MCerPhotEvt &evt, + const MHillas *hillas); + + //virtual void AsciiRead(ifstream &fin); + //virtual void AsciiWrite(ofstream &fout) const; + + ClassDef(MNewImagePar, 1) // Container to hold new image parameters +}; + +#endif + + + + + + + + + + + + + + + + + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mimage/MNewImageParCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MNewImageParCalc.cc (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MNewImageParCalc.cc (revision 1940) @@ -0,0 +1,159 @@ +/* ======================================================================== *\ +! +! * +! * 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 03/2003 +! +! Copyright: MAGIC Software Development, 2000-2003 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MNewImageParCalc +// +// Task to calculate the source dependant part of the hillas parameters +// +////////////////////////////////////////////////////////////////////////////// +#include "MNewImageParCalc.h" + +#include + +#include "MParList.h" + +#include "MGeomCam.h" +#include "MSrcPosCam.h" +#include "MCerPhotEvt.h" +#include "MNewImagePar.h" +#include "MNewImagePar.h" +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MNewImageParCalc); + +static const TString gsDefName = "MNewImageParCalc"; +static const TString gsDefTitle = "Calculate new image parameters"; + +// ------------------------------------------------------------------------- +// +// Default constructor. The first argument is the name of a container +// containing the source position in the camera plain (MScrPosCam). +// The default is "MSrcPosCam". newpar is the name of a container +// of type MNewImagePar, in which the parameters are stored. +// The default is "MNewImagePar" +// +// +MNewImageParCalc::MNewImageParCalc(const char *src, const char *newpar, + const char *name, const char *title) + : fHillas(NULL), fSrcPos(NULL), fNewImagePar(NULL) +{ + fName = name ? name : gsDefName.Data(); + fTitle = title ? title : gsDefTitle.Data(); + + fSrcName = src; + fNewParName = newpar; + fHillasInput = "MHillas"; +} + +// ------------------------------------------------------------------------- +// +Bool_t MNewImageParCalc::PreProcess(MParList *pList) +{ + fHillas = (MHillas*)pList->FindObject(fHillasInput, "MHillas"); + if (!fHillas) + { + *fLog << err << dbginf << "MHillas not found... aborting." << endl; + return kFALSE; + } + + fSrcPos = (MSrcPosCam*)pList->FindObject(fSrcName, "MSrcPosCam"); + if (!fSrcPos) + { + *fLog << err << dbginf << fSrcName << " [MSrcPosCam] not found... aborting." << endl; + return kFALSE; + } + + fCerPhotEvt = (MCerPhotEvt*)pList->FindObject("MCerPhotEvt"); + if (!fCerPhotEvt) + { + *fLog << dbginf << "MCerPhotEvt not found... aborting." << endl; + return kFALSE; + } + + fGeomCam = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fGeomCam) + { + *fLog << dbginf << "MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl; + return kFALSE; + } + + + fNewImagePar = (MNewImagePar*)pList->FindCreateObj("MNewImagePar", fNewParName); + if (!fNewImagePar) + return kFALSE; + + fErrors = 0; + + return kTRUE; +} + +// ------------------------------------------------------------------------- +// +Bool_t MNewImageParCalc::Process() +{ + + if (!fNewImagePar->Calc(*fGeomCam, *fCerPhotEvt, fHillas)) + { + fErrors++; + return kCONTINUE; + + } + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Prints some statistics about the hillas calculation. The percentage +// is calculated with respect to the number of executions of this task. +// +Bool_t MNewImageParCalc::PostProcess() +{ + if (GetNumExecutions()==0) + return kTRUE; + + *fLog << inf << endl; + *fLog << GetDescriptor() << " execution statistics:" << endl; + *fLog << dec << setfill(' '); + *fLog << " " << fErrors << " (" << (int)(fErrors*100/GetNumExecutions()) << "%) Evts skipped due to: calculation failed" << endl; + *fLog << endl; + + return kTRUE; +} + +// -------------------------------------------------------------------------- + + + + + + + + + + + + Index: trunk/MagicSoft/Mars/mimage/MNewImageParCalc.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MNewImageParCalc.h (revision 1940) +++ trunk/MagicSoft/Mars/mimage/MNewImageParCalc.h (revision 1940) @@ -0,0 +1,47 @@ +#ifndef MARS_MNewImageParCalc +#define MARS_MNewImageParCalc + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MHillas; +class MNewImagePar; +class MSrcPosCam; +class MGeomCam; +class MCerPhotEvt; + +class MNewImageParCalc : public MTask +{ +private: + MGeomCam *fGeomCam; + MCerPhotEvt *fCerPhotEvt; + + MHillas *fHillas; //! Pointer to the source independent hillas parameters + MSrcPosCam *fSrcPos; //! Pointer to the source position + MNewImagePar *fNewImagePar; //! Pointer to the output container for the new image parameters + + TString fSrcName; + TString fNewParName; + TString fHillasInput; + + Int_t fErrors; + + Bool_t PreProcess(MParList *plist); + Bool_t Process(); + Bool_t PostProcess(); + +public: + MNewImageParCalc(const char *src="MSrcPosCam", const char *newpar="MNewImagePar", + const char *name=NULL, const char *title=NULL); + + void SetInput(TString hilname) { fHillasInput = hilname; } + + ClassDef(MNewImageParCalc, 1) // task to calculate new image parameters +}; + +#endif + + + + Index: trunk/MagicSoft/Mars/mimage/Makefile =================================================================== --- trunk/MagicSoft/Mars/mimage/Makefile (revision 1940) +++ trunk/MagicSoft/Mars/mimage/Makefile (revision 1940) @@ -0,0 +1,58 @@ +################################################################## +# +# makefile +# +# for the MARS software +# +################################################################## +include ../Makefile.conf.$(OSTYPE) +include ../Makefile.conf.general + +# +# Handling name of the Root Dictionary Files +# +CINT = Image + +# +# Library name to creatre +# +LIB = mimage.a + +# +# connect the include files defined in the config.mk file +# +INCLUDES = -I. -I../mbase -I../mhist + +#------------------------------------------------------------------------------ + +.SUFFIXES: .c .cc .cxx .h .hxx .o + +SRCFILES = MImgCleanStd.cc \ + MCameraSmooth.cc \ + MHillas.cc \ + MHillasSrc.cc \ + MHillasExt.cc \ + MHillasCalc.cc \ + MHillasSrcCalc.cc \ + MNewImagePar.cc \ + MNewImageParCalc.cc \ + MHHillas.cc \ + MHHillasSrc.cc \ + MHHillasExt.cc \ + MHNewImagePar.cc + +SRCS = $(SRCFILES) +HEADERS = $(SRCFILES:.cc=.h) +OBJS = $(SRCFILES:.cc=.o) + +############################################################ + +all: $(LIB) + +include ../Makefile.rules + +clean: rmcint rmobjs rmcore rmlib + +mrproper: clean rmbak + +# @endcode