#ifndef MPDPAIRGG_H
#define MPDPAIRGG_H

#include <deque>

#include "TChain.h"
#include "TVector3.h"
#include "TLorentzVector.h"
#include "TF1.h"
#include "TH1.h"
#include "TH2.h"
#include "TH3.h"
#include "TRandom.h"

#include "MpdEvent.h"
#include "MpdPid.h"
#include "MpdHelix.h"
#include "FairMCEventHeader.h"
#include "MpdKalmanFilter.h"
#include "TpcPoint.h"
#include "MpdTpcKalmanFilter.h"
#include "MpdAnalysisTask.h"
#include "MpdParticle.h"
#include "MpdPairGGTrack.h"
#include "MpdPairGGParams.h"
#include "MpdTofMatching.h"
#include "MpdTofMatchingData.h"
#include "MpdEmcClusterKI.h"

class MpdTpcKalmanTrack;

class MpdPairGG : public MpdAnalysisTask {

public:
   MpdPairGG() {}
   MpdPairGG(const char *name, const char *outputName = "taskName");
   ~MpdPairGG() {} // TODO: normal descructor with cleaning off histos

   void UserInit();
   void ProcessEvent(MpdAnalysisEvent &event);
   void Finish();

   void setOutFile(std::string filename = "histos.root") { mOutFile = filename; }

protected:
   bool selectEvent(MpdAnalysisEvent &event);

   void selectEMCTrack(MpdAnalysisEvent &event);

   bool selectEMCPhot(MpdEmcClusterKI *clu);

   bool selectTrack(MpdTrack *tr);

   void selectPCMTrack(MpdAnalysisEvent &event);

   void processHistograms(MpdAnalysisEvent &event);

   long int IsSameParent(long int prim1, long int prim2) const;

   float EnCorr2P(float e);
   float disthSigma(float disth, float pt);
   float ChiSigma(float Chi, float pt);
   float angleSigma(float angle, float pt);
   float line1(float decay, float mass);
   float line2(float decay, float mass);
   float line3(float decay, float mass);
   float line4(float decay, float mass);
   float line5(float decay, float mass);
   float line6(float decay, float mass);

   MpdHelix MakeHelix(const MpdTpcKalmanTrack *tr);

private:
   float cen;

   // Event properties
   bool     isInitialized = false;
   bool     isMC          = true;
   int      mCenBin       = 0;
   int      mZvtxBin      = 0;
   int      mRPBin        = 0;
   TVector3 mPrimaryVertex;

   std::string     mParamConfig;
   MpdPairGGParams mParams;

   std::string mOutFile = "histos.root";

   //   BaseTpcSectorGeo   *secGeo  = nullptr;
   //   MpdTpcKalmanFilter *recoTpc = nullptr;
   //   MpdKalmanFilter      *mKF              = nullptr;
   TClonesArray *mMCTracks        = nullptr;
   TClonesArray *mKalmanTracks    = nullptr;
   TClonesArray *mMpdGlobalTracks = nullptr;
   TObjArray    *mEMCClusters     = nullptr;
   //   MpdKalmanHit          mKHit;
   TClonesArray *eventM = nullptr; // (V)

   std::vector<MpdPairGGTrack> mP2; // (V) PCM photons
   std::vector<MpdPairGGTrack> mP1; // (V) ECAL photons

   vector<MpdParticle *> vPart;

   static constexpr short nCenBinsAna = 7; // (V) number of bins in centrality for analysis

   static constexpr short nMixEventZ    = 10; // (V) number of bins in z direction for mixing. first:10
   static constexpr short nMixEventCent = 10; // (V) number of bins in centrality for mixing first:10
   static constexpr short nMixEventRP   = 1;  // (V) number of bins in event plane for mixing
   static constexpr short nMixTot       = nMixEventZ * nMixEventCent * nMixEventRP; // (V)

   int nMixed = 2; // (V) Depth of mixing
   int mixBin;
   int anaBin;
   int iv;
   int nMaxG;

   TRandom RND;

   TList *mixedEvents[nMixTot];

   // Histograms
   TList mHistoList;

   // General QA
   TH1D *mhEvents     = nullptr;
   TH1D *mhVertex     = nullptr;
   TH1D *mhCentrality = nullptr;
   TH1D *mhEvPlane    = nullptr;

   TH2D *mMaxG;
   TH2D *mEnResEMC;
   TH2D *mInvEE_before;
   TH2D *mInvEE_after;

   TH1D *mInvGenPi0;
   TH1D *mInvGenEta;
   TH2D *mInvMassGenPi0;
   TH2D *mInvMassGenEta;

   TH2D *mInvNoPID;
   TH2D *mInvPCMNoPID;
   TH2D *mInvTimePID;
   TH2D *mInvTimeMatchPID;
   TH2D *mInvNoPID_Hy;
   TH2D *mInvTimePID_Hy;
   TH2D *mInvTimeMatchPID_Hy;

   TH2D *mInvMixNoPID;
   TH2D *mInvPCMMixNoPID;
   TH2D *mInvMixTimePID;
   TH2D *mInvMixTimeMatchPID;
   TH2D *mInvMixNoPID_Hy;
   TH2D *mInvMixTimePID_Hy;
   TH2D *mInvMixTimeMatchPID_Hy;

   TH2D *mInvTrueNoPID;
   TH2D *mInvPCMTrueNoPID;
   TH2D *mInvTrueNoPIDPi0;
   TH2D *mInvPCMTrueNoPIDPi0;
   TH2D *mInvTrueNoPIDEta;
   TH2D *mInvPCMTrueNoPIDEta;
   TH2D *mInvTrueNoPID_Hy;
   TH2D *mInvTrueNoPIDPi0_Hy;
   TH2D *mInvTrueNoPIDEta_Hy;

   TH2D *mInvTrueTimePID;
   TH2D *mInvTrueTimePIDPi0;
   TH2D *mInvTrueTimePIDEta;
   TH2D *mInvTrueTimePID_Hy;
   TH2D *mInvTrueTimePIDPi0_Hy;
   TH2D *mInvTrueTimePIDEta_Hy;

   TH2D *mInvTrueTimeMatchPID;
   TH2D *mInvTrueTimeMatchPIDPi0;
   TH2D *mInvGenTimeMatchPIDPi0;
   TH2D *mInvTrueTimeMatchPIDEta;
   TH2D *mInvGenTimeMatchPIDEta;
   TH1D *mInvGenBinPi0[nCenBinsAna];
   TH1D *mInvGenBinEta[nCenBinsAna];
   TH2D *mInvTrueNoPIDPi0Bin[nCenBinsAna];
   TH2D *mInvPCMTrueNoPIDPi0Bin[nCenBinsAna];
   TH2D *mInvTrueTimePIDPi0Bin[nCenBinsAna];
   TH2D *mInvTrueTimeMatchPIDPi0Bin[nCenBinsAna];
   TH2D *mInvGenTimeMatchPIDPi0Bin[nCenBinsAna];
   TH2D *mInvTrueNoPIDEtaBin[nCenBinsAna];
   TH2D *mInvTrueTimePIDEtaBin[nCenBinsAna];
   TH2D *mInvTrueTimeMatchPIDEtaBin[nCenBinsAna];
   TH2D *mInvGenTimeMatchPIDEtaBin[nCenBinsAna];
   TH2D *mInvTrueTimeMatchPID_Hy;
   TH2D *mInvTrueTimeMatchPIDPi0_Hy;
   TH2D *mInvGenTimeMatchPIDPi0_Hy;
   TH2D *mInvTrueTimeMatchPIDEta_Hy;
   TH2D *mInvGenTimeMatchPIDEta_Hy;
   TH1D *mInvGenBinPi0_Hy[nCenBinsAna];
   TH1D *mInvGenBinEta_Hy[nCenBinsAna];
   TH2D *mInvTrueNoPIDPi0Bin_Hy[nCenBinsAna];
   TH2D *mInvTrueTimePIDPi0Bin_Hy[nCenBinsAna];
   TH2D *mInvTrueTimeMatchPIDPi0Bin_Hy[nCenBinsAna];
   TH2D *mInvGenTimeMatchPIDPi0Bin_Hy[nCenBinsAna];
   TH2D *mInvTrueNoPIDEtaBin_Hy[nCenBinsAna];
   TH2D *mInvTrueTimePIDEtaBin_Hy[nCenBinsAna];
   TH2D *mInvTrueTimeMatchPIDEtaBin_Hy[nCenBinsAna];
   TH2D *mInvGenTimeMatchPIDEtaBin_Hy[nCenBinsAna];

   TH2D *mInvNoPIDBin[nCenBinsAna];
   TH2D *mInvPCMNoPIDBin[nCenBinsAna];
   TH2D *mInvTimePIDBin[nCenBinsAna];
   TH2D *mInvTimeMatchPIDBin[nCenBinsAna];
   TH2D *mInvNoPIDBin_Hy[nCenBinsAna];
   TH2D *mInvTimePIDBin_Hy[nCenBinsAna];
   TH2D *mInvTimeMatchPIDBin_Hy[nCenBinsAna];

   TH2D *mInvMixNoPIDBin[nCenBinsAna];
   TH2D *mInvPCMMixNoPIDBin[nCenBinsAna];
   TH2D *mInvMixTimePIDBin[nCenBinsAna];
   TH2D *mInvMixTimeMatchPIDBin[nCenBinsAna];
   TH2D *mInvMixNoPIDBin_Hy[nCenBinsAna];
   TH2D *mInvMixTimePIDBin_Hy[nCenBinsAna];
   TH2D *mInvMixTimeMatchPIDBin_Hy[nCenBinsAna];

   TH2D *mAccEffGenPi0;
   TH2D *mAccEffGenEta;
   TH2D *mAccEffRecNoPIDPi0;
   TH2D *mAccPCMEffRecNoPIDPi0;
   TH2D *mAccEffRecTimePIDPi0;
   TH2D *mAccEffRecTimeMatchPIDPi0;
   TH2D *mAccEffRecNoPIDEta;
   TH2D *mAccPCMEffRecNoPIDEta;
   TH2D *mAccEffRecTimePIDEta;
   TH2D *mAccEffRecTimeMatchPIDEta;
   TH2D *mAccEffGenPi0_Hy;
   TH2D *mAccEffGenEta_Hy;
   TH2D *mAccEffRecNoPIDPi0_Hy;
   TH2D *mAccEffRecTimePIDPi0_Hy;
   TH2D *mAccEffRecTimeMatchPIDPi0_Hy;
   TH2D *mAccEffRecNoPIDEta_Hy;
   TH2D *mAccPCMEffRecNoPIDEta_Hy;
   TH2D *mAccEffRecTimePIDEta_Hy;
   TH2D *mAccEffRecTimeMatchPIDEta_Hy;

   TH2D *EE_pt_disth_before;
   TH2D *EE_pt_chi2_before;
   TH2D *EE_pt_angle_before;
   TH2D *EE_pt_decay_before;
   TH2D *EE_pt_mass_before;
   TH2D *EE_decay_mass_before;

   TH2D *EE_pt_disth_after;
   TH2D *EE_pt_chi2_after;
   TH2D *EE_pt_angle_after;
   TH2D *EE_pt_decay_after;
   TH2D *EE_pt_mass_after;
   TH2D *EE_decay_mass_after;




   TH1D *pt_MC_22; 
   TH1D *pt_Ecal;
   TH1D *pt_Ecal_22;
   TH1D *pt_tcl;
   TH1D *pt_tcl_22;
   TH1D *pt_tcl_match;
   TH1D *pt_tcl_match_22;
   TH1D *pt_PCM_after;
   TH1D *pt_PCM_22_after;



   //   TFile        *inFileSim;
   //   TTree        *inTreeSim;
   //   TClonesArray *tpcPoints;

   ClassDef(MpdPairGG, 1);
};
#endif