#ifndef MARS_MSkyPlot
#define MARS_MSkyPlot

#ifndef MARS_MTask
#include "MTask.h"
#endif

#include <vector>

#ifndef ROOT_TH2
#include <TH2.h>
#endif

#ifndef ROOT_TH1
#include <TH1.h>
#endif

class TH2D;
class TH1D;

class MParList;
class MTime;
class MPointingPos;
class MSrcPosCam;
class MReportDrive;
class MObservatory;
class MHillas;
class MHillasExt;
class MHillasSrc;
class MNewImagePar;
class MGeomCam;

class MSkyPlot : public MTask
{
private:
    MGeomCam      *fGeomCam;        //! container to take the event time from
    MTime         *fTime;        //! container to take the event time from
    MPointingPos  *fPointPos;    //! container to take pointing position from
    MReportDrive  *fRepDrive;      
    MSrcPosCam    *fSrcPosCam;    //! container with x and y of the source
    MSrcPosCam    *fPntPosCam;    //! container with x and y of the position MReportDrive.GetRa, MReportDrive.GetDec
    MObservatory  *fObservatory; //! container to take observatory location from
    MHillas       *fHillas;      
    MHillasExt    *fHillasExt;
    MHillasSrc    *fHillasSrc;
    MNewImagePar  *fNewImagePar;

    Float_t        fMm2Deg;             // conversion factor for display in degrees
    Double_t       fGridBinning;   // degrees
    Double_t       fGridFineBin;   // degrees

//    Float_t fAlphaCut;           // Alpha cut
//    Float_t fBgMean;             // Background mean

//    Float_t fMinDist;            // Min dist
//    Float_t fMaxDist;            // Max dist

//    Float_t fMinLD;              // Minimum distance in percent of dist
//    Float_t fMaxLD;              // Maximum distance in percent of dist

    std::vector <TH1D> fHistAlpha;           // vector of histograms for alpha
    Int_t fNumalphahist;		// number of histograms for alpha
    Int_t fNumBinsAlpha;
    Float_t fHistAlphaBinWidth;
    Float_t fAlphaLeftEdge;
    Float_t fAlphaRightEdge;
    Float_t fAlphaONMax;			//  [deg] , upper cut for alpha ON region in the alpha plot 
                      // [deg], ON region in the alpha plot, maybe 5 deg is better
                      // NOTE: up to now only values of 5, 10, 15, 20 degrees are possible
    Float_t fAlphaBgLow;	   // lower limit for bg region in the ON alpha plot 
    Float_t fAlphaBgUp;		   // upper limit for bg region in the ON alpha plot
    
    TH2D     fHistSignif;                // sky plot of significance vs. x and y
    TH2D     fHistNexcess;               // sky plot of number of excess events vs. x and y
    TH2D     fHistOn;                    // sky plot of events below fAlphaONMax vs. x and y
    TH1D     fHistSignifGaus;            // distribution of significance
    Bool_t   fSetCenter;                 // used to set the center of these histograms once
    Double_t fRa0;
    Double_t fDec0;
    Bool_t kSaveAlphaPlots;
    Bool_t kSaveSkyPlots;

    Float_t fMinXGrid;			//  [deg] , left edge of the skyplot
    Float_t fMaxXGrid;			//  [deg] , right edge of the skyplot
    Float_t fMinYGrid;			//  [deg] , upper edge of the skyplot
    Float_t fMaxYGrid;			//  [deg] , lower edge of the skyplot
    Float_t fBinStepGrid;			//  [deg] , grid size
    Int_t fNumStepsX;			// number of bins in x 
    Int_t fNumStepsY;			// number of bins in y


    // some filter cuts:
    Float_t fSizeMin;		// min size in photons
    Float_t fSizeMax;		// max size in photons
    Float_t fLeakMax;		// leakmax  in per cent
    Float_t fMaxDist;          // dist max cut (ever possible)
    Float_t fMinDist;          // dist min cut (ever possible)

    // coefficients for the cuts:
    Double_t fLengthUp[8], fLengthLo[8], fWidthUp[8], fWidthLo[8], fDistUp[8], fDistLo[8];

    Double_t fRa;
    Double_t fDec;

    TString fAlphaHName;          // name for histogram with alpha plots
    TString fSkyHName;            // name for histogram with sky plots

    Int_t DistancetoPrimitive(Int_t px, Int_t py);

    TObject *GetCatalog();

    Int_t PreProcess(MParList *pList);
    Int_t Process();
    Int_t PostProcess();


public:
    MSkyPlot(const char *name=NULL, const char *title=NULL);

//    TH1 *GetHistByName(const TString name) { return &fHist; }
    TH2D *GetHistSignif() { return &fHistSignif; }
    TH2D *GetHistNexcess() { return &fHistNexcess; }
    TH2D *GetHistOn() { return &fHistOn; }
    TH1D *GetHistSignifGaus() { return &fHistSignifGaus; }

    void SetMinDist(Float_t dist) { fMinDist = dist; } // Absolute minimum distance
    void SetMaxDist(Float_t dist) { fMaxDist = dist; } // Absolute maximum distance
    void SetSizeMin(Float_t size) { fSizeMin = size; } // Absolute minimum Size
    void SetSizeMax(Float_t size) { fSizeMax = size; } // Absolute maximum Size
//    void SetMinLD(Float_t ratio)  { fMinLD = ratio; }  // Minimum ratio between length/dist
//    void SetMaxLD(Float_t ratio)  { fMaxLD = ratio; }  // Maximum ratio between length/dist
    void ReadCuts(const TString parSCinit);
    void SetSkyPlot(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t step);
    Double_t CalcLimit(Double_t *a, Double_t ls, Double_t ls2, Double_t dd2);
    //void SaveSkyPlots(const TString skyplotfilename="skyplots.root");

    void SetOutputAlphaName(TString outname1)   { fAlphaHName = outname1; }
    void SaveSkyPlots(TString stri);


    //void SaveAlphaPlots(const TString alphaplotfilename="alphaplots.root");
    void SetOutputSkyName(TString outname2)     { fSkyHName = outname2; }
    void SaveAlphaPlots(const TString stri2);

    void SetAlphaCut(Float_t alpha); 
    void SetAlphaBGLimits(Float_t alphalow, Float_t alphalup); 

    //std::vector <TH1D>::iterator  GetAlphaPlots()      { std::vector <TH1D>::iterator iter; return  iter = fHistAlpha.begin()  ; }
    std::vector <TH1D> GetAlphaPlots()      { return  fHistAlpha  ; }
 
    void Draw(Option_t *option="");

    ClassDef(MSkyPlot, 0) //2D-histogram in alpha, x and y
};

#endif
