File indexing completed on 2024-04-06 12:27:41
0001
0002
0003
0004
0005
0006
0007 #ifndef RecoPPS_ProtonReconstruction_ProtonReconstructionAlgorithm_h
0008 #define RecoPPS_ProtonReconstruction_ProtonReconstructionAlgorithm_h
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "DataFormats/CTPPSReco/interface/CTPPSLocalTrackLiteFwd.h"
0013 #include "DataFormats/ProtonReco/interface/ForwardProtonFwd.h"
0014
0015 #include "CondFormats/PPSObjects/interface/LHCInterpolatedOpticalFunctionsSet.h"
0016 #include "CondFormats/PPSObjects/interface/LHCInterpolatedOpticalFunctionsSetCollection.h"
0017
0018 #include "TSpline.h"
0019 #include "Fit/Fitter.h"
0020
0021 #include <unordered_map>
0022
0023
0024
0025 class ProtonReconstructionAlgorithm {
0026 public:
0027 ProtonReconstructionAlgorithm(bool fit_vtx_y,
0028 bool improved_estimate,
0029 const std::string &multiRPAlgorithm,
0030 unsigned int verbosity);
0031 ~ProtonReconstructionAlgorithm() = default;
0032
0033 void init(const LHCInterpolatedOpticalFunctionsSetCollection &opticalFunctions);
0034 void release();
0035
0036
0037 reco::ForwardProton reconstructFromSingleRP(const CTPPSLocalTrackLiteRef &track,
0038 const float energy,
0039 std::ostream &os) const;
0040
0041
0042 reco::ForwardProton reconstructFromMultiRP(const CTPPSLocalTrackLiteRefVector &tracks,
0043 const float energy,
0044 std::ostream &os) const;
0045
0046 private:
0047 unsigned int verbosity_;
0048 bool fitVtxY_;
0049 bool useImprovedInitialEstimate_;
0050 enum { mrpaUndefined, mrpaChi2, mrpaNewton, mrpaAnalIter } multi_rp_algorithm_;
0051 bool initialized_;
0052
0053
0054 struct RPOpticsData {
0055 const LHCInterpolatedOpticalFunctionsSet *optics;
0056 std::shared_ptr<const TSpline3> s_x_d_vs_xi, s_L_x_vs_xi, s_xi_vs_x_d, s_y_d_vs_xi, s_v_y_vs_xi, s_L_y_vs_xi;
0057 double x0;
0058 double y0;
0059 double ch0;
0060 double ch1;
0061 double la0;
0062 double la1;
0063 };
0064
0065
0066 std::map<unsigned int, RPOpticsData> m_rp_optics_;
0067
0068
0069 class ChiSquareCalculator {
0070 public:
0071 ChiSquareCalculator() = default;
0072
0073 double operator()(const double *parameters) const;
0074
0075 const CTPPSLocalTrackLiteRefVector *tracks;
0076 const std::map<unsigned int, RPOpticsData> *m_rp_optics;
0077 };
0078
0079
0080 std::unique_ptr<ROOT::Fit::Fitter> fitter_;
0081
0082
0083 std::unique_ptr<ChiSquareCalculator> chiSquareCalculator_;
0084
0085 static void doLinearFit(const std::vector<double> &vx, const std::vector<double> &vy, double &b, double &a);
0086
0087 static double newtonGoalFcn(double xi, double x_N, double x_F, const RPOpticsData &i_N, const RPOpticsData &i_F);
0088 };
0089
0090 #endif