#ifndef __HPARTICLETOOL_H__
#define __HPARTICLETOOL_H__
#include "TObject.h"
#include "TObjArray.h"
#include "TLorentzVector.h"
#include "TH1.h"
#include "TF1.h"
#include "TCutG.h"
#include <vector>
using namespace std;
class HRichHit;
class HTofHit;
class HTofCluster;
class HRpcCluster;
class HShowerHit;
class HEmcCluster;
class HMetaMatch2;
class HMdcTrkCand;
class HMdcSeg;
class HMdcHit;
class HMdcCal1;
class HMdcClusInf;
class HMdcClusFit;
class HMdcWireFit;
class HMdcClus;
class HGeomVector;
class HVirtualCand;
class HParticleCand;
class HParticleCandSim;
class HParticlePair;
class HMdcPlane;
class HMdcLayer;
class HTofWalkPar;
class HGeantKine;
#define TOFSIGSIM 33
#define RPCSIG1 7.62
#define RPCSIG2 5.96
#define RPCSIG1SIM 8.32
#define RPCSIG2SIM 5.96
class HParticleTool : public TObject {
static Float_t rpcCellHalfWidth[192];
static Float_t tofCellHalfWidth[64] ;
static Double_t scaleDyPars[4];
static TF1* gf1;
static TF1* gf2;
static TF1* gfsum;
static TF1* gScaledy;
static TF1* fdxoffset;
static TF1* fdxsigma ;
static Double_t parsSX[6][8][8][3];
static Double_t parsDX[6][8][8][5];
static Double_t parsSX_apr12[6][8][8][3];
static Double_t parsDX_apr12[6][8][8][5];
static TString beamtime_tof;
public:
HParticleTool();
~HParticleTool();
static Float_t phiSecToLabDeg(Int_t sec, Float_t phiRad);
static Float_t thetaToLabDeg(Float_t thetaRad);
static Float_t phiLabToPhiSecDeg(Float_t phiLabDeg);
static Int_t phiLabToSec(Float_t phiLabDeg);
static Float_t getOpeningAngle(Float_t phi1,Float_t theta1,Float_t phi2,Float_t theta2);
static Float_t getOpeningAngle(TLorentzVector& vec1,TLorentzVector& vec2);
static Float_t getOpeningAngle(HParticleCand* cand1,HParticleCand* cand2);
static Float_t getOpeningAngle(HGeantKine* kine1,HGeantKine* kine2);
static Bool_t setCloseCandidates(Float_t oACut=15., Bool_t sameSector=kTRUE, Bool_t skipSameSeg=kTRUE);
static Int_t getCloseCandidates(HParticleCand* cand,vector<HParticleCand*>& vcand,vector<Float_t >& vopeninAngle,Float_t oACut=15., Bool_t sameSector=kTRUE, Bool_t skipSameSeg=kTRUE);
static Int_t getCloseCandidatesSegInd(HParticleCand* cand, vector<Int_t>& vSeg,Float_t oACut,Bool_t sameSector, Bool_t skipSameSeg);
static void getTLorentzVector(HGeantKine* kine, TLorentzVector& vec,Int_t pid=-1);
static void fillTLorentzVector(TLorentzVector& v,HVirtualCand* cand,Float_t mass);
static void fillTLorentzVector(TLorentzVector& v,HVirtualCand* cand,Int_t pid,Bool_t correctMom=kTRUE);
static Float_t getLabPhiDeg(TLorentzVector& vec);
static Float_t calcRichQA(HMdcSeg* seg, HRichHit* hit);
static Float_t calcRichQA(HMdcSeg* seg, Float_t richTheta,Float_t richPhi);
static HGeomVector getGlobalVertex(Int_t v,Bool_t warn=kFALSE);
static Double_t getMinimumDistToVertex(HParticleCand*,HGeomVector& vertex);
static HGeomVector getPointOfClosestApproachToVertex(HParticleCand*,HGeomVector& vertex);
static Double_t scaledy (Double_t* x, Double_t* par);
static Double_t getScaledDy(HParticleCand* c,Double_t dyCut=-1);
static TF1* getScaleTF1() { return gScaledy ;}
static Bool_t normDX(HParticleCand* c,TString beamtime="apr12");
static Float_t getNormDX(HParticleCand* c,TString beamtime="apr12");
static Float_t getSigmaDX(HParticleCand* c,TString beamtime="apr12");
static Bool_t normDX(HParticleCand* c,HTofWalkPar* p);
static Float_t getNormDX(HParticleCand* c,HTofWalkPar* p);
static Float_t getSigmaDX(HParticleCand* c,HTofWalkPar* p);
static TF1* getTofXOffsetTF1() {return fdxoffset;}
static TF1* getTofXSigmaTF1() {return fdxsigma ;}
static Float_t getRpcCellHalfWidth (Int_t mod,Int_t cell);
static Float_t getTofCellHalfWidth (Int_t mod,Int_t cell);
static Bool_t isGoodMetaCell (HParticleCand* c,Double_t bound=3.5,Bool_t doScaling=kTRUE);
static Float_t getCorrectedMomentum(HParticleCand* c);
static Float_t setCorrectedMomentum(HParticleCand* c);
static Bool_t isParticledEdx(Int_t PID, HParticleCand* c, Float_t& deloss, Float_t& dsigma);
static Bool_t isParticleBeta(Int_t PID, HParticleCand* pCand, Float_t nsigma, Float_t momMin, Float_t momMax,Float_t& dtime, Float_t& dsigma,TString beamtime="apr12");
static Bool_t correctPathLength(HParticleCand* pCand, HGeomVector& vertex,const HMdcPlane* planes, const HGeomVector& targetMidPoint, Double_t beamEnergy = 1230);
static Bool_t checkCropedLayer(HGeantKine* kine,HMdcLayer* mdcLayer, Bool_t* croped=0,Bool_t checkHit=kTRUE);
static Double_t beta(Int_t id,Double_t p);
static Double_t betaToP(Int_t id,Double_t beta);
static Double_t gamma(Int_t id,Double_t p);
static Double_t gammaToBeta(Int_t id,Double_t gamma);
static Double_t gammaToP(Int_t id,Double_t gamma);
static Double_t betagamma(Int_t id,Double_t p);
static Double_t betagammaToP(Int_t id,Double_t betagamma);
static Double_t kinEToMom(Int_t id = 14, Double_t Ekin = 3500);
static Double_t momToKinE(Int_t id = 14, Double_t p = 3500);
static Double_t dedxfunc(Double_t *x, Double_t *par);
static Double_t betaandgammafunc(Double_t *x, Double_t *par);
static Double_t ptyfunc (Double_t* x, Double_t* par);
static Double_t fcross(Double_t* xin, Double_t* par);
static Bool_t getIntersectionPoint(TF1* f1,TF1* f2,Double_t &xout,Double_t& yout,Double_t xlow,Double_t xup,Int_t n=500);
static TF1* energyLossTF1(Int_t id,TString name="",TString opt="p",Double_t scaleY=1,Double_t xoffset=0,Double_t scaleX=1,Double_t theta=-1,Double_t frac=0.1,Double_t xmin=20,Double_t xmax=2000,Int_t linecolor = 2,Int_t linestyle = 1,Int_t npoints=500);
static TF1* betaAndGammaTF1(Int_t id,TString name="",TString opt="beta",Double_t scaleY=1,Double_t xoffset=0,Double_t scaleX=1,Double_t theta=-1,Double_t frac=0.1,Double_t xmin=20,Double_t xmax=2000,Int_t linecolor = 2,Int_t linestyle = 1,Int_t npoints=500);
static TCutG* makeCut(TF1* lower,TF1* upper,TString name,Double_t xlow,Double_t xup,Double_t ymin=-1000,Double_t ymax=1000,Int_t npoint=500,Int_t linecolor = 2,Int_t linestyle = 2);
static TF1* ptyTF1(Int_t id,Double_t val=45,TString name="",TString opt="theta",Double_t xmin=0,Double_t xmax=2,Double_t midRap=0,Int_t linecolor = 2,Int_t linestyle = 1);
static vector<TF1*> ptyGrid(Int_t id,vector<Double_t>& vtheta,vector<Double_t>& vmom,TString name="",TString opt = "draw",Double_t xmin=0,Double_t xmax=2,Double_t midRap=0,Int_t linecolorTheta = 1,Int_t linestyleTheta = 2,Int_t linecolorMom = 1,Int_t linestyleMom = 2);
static vector<TF1*> ptyGrid(Int_t id,TString setup="@theta:8,0,10@momentum:20,10,100",TString name="",TString opt = "draw",
Double_t xmin=0,Double_t xmax=2,Double_t midRap=0,
Int_t linecolorTheta=1,Int_t linestyleTheta=2,Int_t linecolorMom=1,Int_t linestyleMom=2,
TString labels="@theta:draw=yes,format=%5.1f#circ,textsize=0.021,angle=0,align=-1@momentum:draw=yes,format=%5.1f MeV/c,textsize=0.023,angle=10,align=-1");
static void drawPtyGrid(vector<TF1*>& grid,TString opt = "draw");
static void calcSegVector(Double_t z, Double_t rho, Double_t phi, Double_t theta, HGeomVector &base, HGeomVector &dir);
static Double_t calcRMS(const Double_t* valArr, Double_t Mean,Int_t valNum);
static Double_t calcDeterminant(HGeomVector& v1, HGeomVector& v2, HGeomVector& v3);
static Double_t calculateMinimumDistanceStraightToPoint(HGeomVector &base, HGeomVector &dir,HGeomVector &point);
static HGeomVector calculatePointOfClosestApproachStraightToPoint(HGeomVector &base, HGeomVector &dir,HGeomVector &point);
static Double_t calculateMinimumDistance(HGeomVector &base1, HGeomVector &dir1, HGeomVector &base2, HGeomVector &dir2);
static HGeomVector calculatePointOfClosestApproach(HGeomVector &base1, HGeomVector &dir1,HGeomVector &base2, HGeomVector &dir2);
static HGeomVector calculateCrossPoint(HGeomVector &base1, HGeomVector &dir1,HGeomVector &base2, HGeomVector &dir2);
static HGeomVector calcVertexAnalytical(HGeomVector &base1, HGeomVector &dir1,HGeomVector &base2, HGeomVector &dir2);
static Int_t findFirstHitInTof (Int_t trackID,Int_t modeTrack = 2);
static Int_t findFirstHitShowerInTofino(Int_t trackID,Int_t modeTrack = 2);
static Int_t findFirstHitShowerInRpc (Int_t trackID,Int_t modeTrack = 2);
static Float_t getInterpolatedValue(TH1* h, Float_t xVal,Bool_t warn = kTRUE);
static Stat_t getValue(TH1* h,Float_t xVal, Float_t yVal = 0.0f, Float_t zVal = 0.0f);
static HRichHit* getRichHit (Int_t richind);
static HTofHit* getTofHit (Int_t tofind);
static HTofCluster* getTofCluster(Int_t tofind);
static HRpcCluster* getRpcCluster(Int_t rpcind);
static HShowerHit* getShowerHit (Int_t showerind);
static HEmcCluster* getEmcCluster(Int_t emcind);
static HMetaMatch2* getMetaMatch (Int_t metaind);
static HMdcTrkCand* getMdcTrkCand(Int_t metaind);
static HMdcSeg* getMdcSeg (Int_t segind);
static HMdcHit* getMdcHit (Int_t segind,Int_t nhit = 0);
static HMdcClusInf* getMdcClusInf(Int_t segind,Int_t nhit = 0);
static HMdcClusFit* getMdcClusFit(Int_t segind);
static HMdcClus* getMdcClus (Int_t segind);
static TObjArray* getMdcWireFitSeg (Int_t segind);
static TObjArray* getMdcCal1Seg (Int_t segind);
static TObjArray* getMdcCal1Cluster(Int_t segind);
static Int_t getMdcWireFitSeg (Int_t segind, vector<HMdcWireFit*>& v, Bool_t clear=kTRUE);
static Int_t getMdcCal1Seg (Int_t segind, vector<HMdcCal1*>& v, Bool_t clear=kTRUE);
static Int_t getMdcCal1Cluster(Int_t segind, vector<HMdcCal1*>& v, Bool_t clear=kTRUE);
static Bool_t printSimTracks(HParticleCandSim* c);
static Bool_t getSimTracks (HParticleCandSim* c,
vector<Int_t>&tracksMeta,
vector<Int_t>&tracksShowerEcal,
vector<Int_t>&tracksRich,
vector<Int_t>&weightRich,
vector<Int_t>&tracksInnerMdc,
vector<Int_t>&weightInnerMdc,
vector<Int_t>&tracksOuterMdc,
vector<Int_t>&weightOuterMdc,
Bool_t print=kFALSE
);
static Bool_t setPairFlags (UInt_t& flag,HParticleCand* cand2=0,HParticleCand* cand1=0);
static Bool_t evalPairsFlags(UInt_t flag,UInt_t fl);
static Bool_t isPairsFlagsBit(UInt_t flag,UInt_t fl);
static Bool_t evalPairsFlags(UInt_t flag,HParticleCand* cand1,HParticleCand* cand2);
static Bool_t evalPairsFlags(UInt_t flag,HParticlePair& pair);
static Bool_t evalPairsFlags(vector<UInt_t>& flags,HParticlePair& pair);
static Bool_t evalPairsFlags(vector<UInt_t>& flags,vector<Bool_t>& results,HParticlePair& pair);
static Bool_t isGoodClusterVertex(Float_t minZ);
static Bool_t isGoodRecoVertex (Float_t minZ);
static Bool_t isGoodSTART (Int_t minFlag);
static Bool_t isNoVETO (Float_t minStart=-15.,Float_t maxStart=15.);
static Bool_t isNoSTART (Float_t minStart=-15.,Float_t maxStart=15.,Float_t window=1.);
static Bool_t isGoodSTARTVETO (Float_t minStart=15.,Float_t maxStart = 350., Float_t window=1.) ;
static Bool_t isGoodSTARTMETA (Float_t minStart=80.,Float_t maxStart = 350.,Int_t tresh=4,Float_t window=5., Float_t offset=7.);
static Bool_t isNoSTARTPileUp ();
static Bool_t isNoMETAPileUp (Float_t ftimeTofCut,Int_t threshold);
static Bool_t isGoodTrigger (Int_t triggerbit);
static Int_t getTofHitMult (Float_t minTof=0,Float_t maxTof=35., Int_t* sector=NULL);
static Int_t getRpcHitMult (Float_t minTof=0,Float_t maxTof=25., Int_t* sector=NULL);
ClassDef(HParticleTool,0)
};
#endif //__HPARTICLETOOL_H__