#ifndef SDRECLPNHE_H #define SDRECLPNHE_H #include #include #include #include #include #include #include #include "TMinuit.h" #include "TF1.h" #include "IoSd.h" #include "Es.h" #include "TReconstruct.h" #include "ToolsLib.h" #define SIGMA_T_2 (35.*35.) ///< old time variance #define XVGD (875.5) ///< Xmax at ground #define DISTREF 1000 #define DISTREFINFILL 450 #define SINTH2MAX 1.05 #define CHI2LDFMAX 50 #define DCHI2LDFMAX 50 #define DIST_BARY_CORE_MAX 500 #define SILSIGTHRESH 3. ///< signal threshold for silent stations #define SILSIGTHRESHNEW 1. ///same with new triggers #define SIGTHRESH 0.2 namespace SdLpnhe { enum TVarianceTs {eR4, eT50, eT50Mod, eT70, eT90}; enum TGeometry {eGeomNone, ePlane, eCurve, eVCurve, eGeomOptimal}; enum TLDF {eLDFNone, eLog, eLogXmax, eNKG, eLDFOptimal}; enum TTimeRec {eTimeNone, eRise, eFall}; enum TPrimary {eProton, eIron, eOther}; enum TFit {eGeom, eLDF, eGlobal}; } using namespace SdLpnhe; /// LPNHE reconstructor. Based on the reconstruction of Er v2* class SdRecLpnhe: public TReconstruct { private: // Singleton stuff SdRecLpnhe(); ~SdRecLpnhe(); static SdRecLpnhe * _instance; // internal parameters int _errflag; ///< reconstruction flag vector fWgt; /// Weight of station [i] in time fit; default 1.0 void SetWgt(TCalibStation *st, double w); // Internal Recontruction functions int _GeometricalRec(TGeometry geom, bool allowChanges = false); int _LDFRec (TLDF ldf); int _CompleteRec(TGeometry geom, TLDF ldf); int _Reconstruct(TGeometry geom, TLDF ldf, bool global); // Internal auxiliary functions void _EstimateLDFPar(TLDF ldf = eLog); void _EstimateSRef(); void _EstimateXYZTCore(); void _CalculateXmax(); bool _Minimize(TMinuit*); bool _FitBeta(); int _CountDoublets(); int _CompactConfig(); int _EstimateUV(); int _OptimizeStationSelection(int loop=1); // int _OptimizeStationSelection2(); int _TimeRec(TTimeRec); double _StatDt(int, int, TGeometry); double _StatT(int, int, TGeometry); double _Stat2Core(int); double _MatPlanFit(int, int&); double _MatCurvedFit(int, int&); double _TimeResiduals(TGeometry, bool&); double _EstimateTCore(); double _EstimateCurvature(double signal,double costheta,bool infillrec); double _GetResults(TMinuit*, TFit, TGeometry, TLDF); TMinuit* _NewMinuit(TFit, TGeometry, TLDF); TGeometry _Best(TGeometry); TLDF _Best(TLDF); /// This function (inherited fomr TReconstruct)) is called automatically by TErEvent.Reconstruct(). /// The reconstruction behavior is driven by the global parmaters as set with the Set-functions above. bool _Reconstruct(); public: static SdRecLpnhe* Instance(); double Wgt(TCalibStation* t); /// Weight of station [i] in time fit; default 1.0 //double T50(int stat, bool doIt = false); /// Selected stations average T50 // Set Functions to drive the reconstruction void SetRec(TGeometry rec1 = eGeomOptimal, TLDF rec2 = eLDFOptimal); void SetAsymmetry(bool asym = true); void SetPrimary(TPrimary prim = eProton); void SetEnergyConverter(string e = "LOG"); void SetVerbosity(int verb = 0); void SetCore(double xc = 0.0, double yc = 0.0, double zc = 1400.0); void SetSilent(bool silent = true, double Rcut = 5000.0); void SetForce(bool force = false); void SetVarianceModel(TVarianceTs var = eT50); void SetOfficial(bool official = false); void SetOptimize(bool optimize = true); void SetGlobal(bool = true); void SetDefault(); /// Puts back everything to original default value. void SetErEvent(TErEvent *ev); /// Called by ErEvent constructor to associate ErEvent to reconstructor. int InitReconstruction();///Called before any independent reconstruction process double CorrectTime(int); double TimeVariance(int, double, double, TGeometry); // Xmax stuff int XmaxOK(); // User driven functions bool Reconstruct() { return fErEvent->Reconstruct(this); } // theLpnheRec().Reconstruct() is equivalent to TErEvent.Reconstruct(theLpnheRec()); int GeometricalRec(TGeometry geom = eGeomOptimal);///User-callable angular reconstruction. options ePlan, eCurve, eVCurve. int LDFRec(TLDF ldf = eLDFOptimal);///User-callable LDF reconstruction. Option eLog, eXmax int CompleteRec(TGeometry geom = eGeomOptimal, TLDF ldf = eLDFOptimal);///User-callable complete reconstruction; Options as above int HorizontalReconstruction(); /// dumps for the reconstruction results void StationDump(int, TGeometry); void HeaderDump(); void RecDump(TLDF); void HorizontalReconstructionDump(TGeometry = eGeomOptimal); void AngularDump(TGeometry, TLDF); void GenResPlot(TGeometry); void LDFDump(TGeometry, TLDF); void HorizontalResultsDump(); void CompleteDump(TGeometry, TLDF); double XmaxfromRiseTime;///< Xmax estimated from the risetime at r = 1000 m double XmaxfromFallTime;///< Xmax estimated from the falltime at r = 1000 m // max time residuals // unsigned int RESMAX, RESPLANMAX, RESPLANSTATMAX, RESSTATMAX; // parameters used in fits double fn_sRef; double ftCore; double fdBeta,fdGamma; double ftau, fdtau; double falpha, fdalpha; //for timeResiduals fit vector fRPlotTimes; vector fRPlotDts; vector fRPlotDist; vector fRPlotDtsErr; vector fRPlotDistErr; vector fRPlotIds; // chi2s... double angularChi2, ldfChi2, globalChi2; double residu; unsigned int LDFFitOK, GlobalFitOK, ndof; ostream & put(ostream & s, char *option); }; SdRecLpnhe& theLpnheRec(); ostream& operator<<(ostream & s, SdRecLpnhe & r); const char *ErSetFormat(char *option); #endif