Loading [MathJax]/extensions/tex2jax.js
ATLAS Offline Software
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
RtFromPoints.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 
8 #include "GaudiKernel/MsgStream.h"
11 
15 
18 
23 
26 
29 
31 
32 
33 using namespace MuonCalib;
34 
35 
36 
37 CalibFunc::ParVec RtFromPoints::legendreFit(const std::vector<SamplePoint>& dataPoints,
38  const unsigned order) {
39  const auto [minT, maxT] = interval(dataPoints);
40 
41  CalibFunc::ParVec pars{minT, maxT};
43  LegendrePolynomial legendre{};
44  fitter.fit_parameters(normalizeDomain(dataPoints), 1, dataPoints.size(), legendre);
45  for (unsigned k = 0; k < fitter.coefficients().size(); k++) {
46  pars.emplace_back(fitter.coefficients()[k]);
47  }
48  return pars;
49 }
50 
51 CalibFunc::ParVec RtFromPoints::chebyFit(const std::vector<SamplePoint>& dataPoints,
52  const unsigned order) {
53 
54  const auto [minT, maxT] = interval(dataPoints);
55  CalibFunc::ParVec pars{minT, maxT};
57  ChebyshevPolynomial chebyshev{}; // Chebyshev polynomial
58 
59 
60  fitter.fit_parameters(normalizeDomain(dataPoints), 1, dataPoints.size(), chebyshev);
61  for (unsigned k = 0; k < fitter.coefficients().size(); k++) {
62  pars.emplace_back(fitter.coefficients()[k]);
63  }
64  return pars;
65 }
66 CalibFunc::ParVec RtFromPoints::simplePolyFit(const std::vector<SamplePoint>& dataPoints,
67  const unsigned order) {
68  const auto [minT, maxT] = interval(dataPoints);
69  CalibFunc::ParVec pars{minT, maxT};
71  SimplePolynomial simplePoly{};
72  fitter.fit_parameters(normalizeDomain(dataPoints), 1, dataPoints.size(), simplePoly);
73  for (unsigned k = 0; k < fitter.coefficients().size(); k++) {
74  pars.emplace_back(fitter.coefficients()[k]);
75  }
76  return pars;
77 
78 }
79 std::unique_ptr<IRtRelation> RtFromPoints::getRtChebyshev(const std::vector<SamplePoint> &dataPoints, const unsigned order) {
80  return std::make_unique<RtChebyshev>(chebyFit(dataPoints, order));
81 } // end RtFromPoints::getRtChebyshev
82 std::unique_ptr<ITrRelation> RtFromPoints::getTrChebyshev(const std::vector<SamplePoint>& dataPoints, const unsigned order) {
83  return std::make_unique<TrChebyshev>(chebyFit(dataPoints, order));
84 }
85 std::unique_ptr<IRtResolution> RtFromPoints::getResoChebyshev(const std::vector<SamplePoint>& dataPoints, const unsigned order){
86  return std::make_unique<RtResolutionChebyshev>(chebyFit(dataPoints, order));
87 }
88 std::unique_ptr<IRtResolution> RtFromPoints::getResoChebyshev(const std::vector<SamplePoint>& dataPoints, const IRtRelationPtr rtRelPtr,
89  const double relUnc, const unsigned order) {
90  std::vector<double> chebyCoeff = chebyFit(resoFromRadius(dataPoints, relUnc), order);
91  chebyCoeff.erase(chebyCoeff.begin(), chebyCoeff.begin() +2);
92  return std::make_unique<RadiusResolutionChebyshev>(chebyCoeff, rtRelPtr);
93 }
94 
95 std::unique_ptr<IRtRelation> RtFromPoints::getRtLegendre(const std::vector<SamplePoint>& dataPoints, const unsigned order) {
96  return std::make_unique<RtLegendre>(legendreFit(dataPoints, order));
97 }
98 std::unique_ptr<ITrRelation> RtFromPoints::getTrLegendre(const std::vector<SamplePoint>& dataPoints, const unsigned order) {
99  return std::make_unique<TrLegendre>(legendreFit(dataPoints, order));
100 }
101 std::unique_ptr<IRtRelation> RtFromPoints::getRtSimplePoly(const std::vector<SamplePoint>& dataPoints, const unsigned order) {
102  return std::make_unique<RtSimplePolynomial>(simplePolyFit(dataPoints, order));
103 }
104 std::unique_ptr<ITrRelation> RtFromPoints::getTrSimplePoly(const std::vector<SamplePoint>& dataPoints, const unsigned order) {
105  return std::make_unique<TrSimplePolynomial>(simplePolyFit(dataPoints, order));
106 }
107 //*****************************************************************************
108 
109 //::::::::::::::::::::::::::::::::
110 //:: METHOD getRtRelationLookUp ::
111 //::::::::::::::::::::::::::::::::
112 std::unique_ptr<IRtRelation> RtFromPoints::getRtRelationLookUp(const std::vector<SamplePoint> &dataPoints) {
113  // create spline rt relation
114  CalibFunc ::ParVec pars(2 * dataPoints.size());
115  for (unsigned i = 0; i < dataPoints.size(); i++) {
116  pars[2 * i] = dataPoints[i].x1();
117  pars[2 * i + 1] = dataPoints[i].x2();
118  }
119  RtSpline rt(pars);
120 
121  // variables
122  unsigned nb_points(100); // number of (r, t) points
123  double bin_width((rt.tUpper() - rt.tLower()) / static_cast<double>(nb_points - 1)); // step size
124  std::vector<double> rt_param(nb_points + 2); // r-t parameters
125 
127  // CREATE AN RtRelationLookUp OBJECT WITH THE CORRECT PARAMETERS //
129  rt_param[0] = rt.tLower();
130  rt_param[1] = bin_width;
131  for (unsigned k = 0; k < nb_points; k++) { rt_param[k + 2] = rt.radius(rt.tLower() + k * bin_width); }
132  return std::make_unique<RtRelationLookUp>(rt_param);
133 }
TrSimplePolynomial.h
make_hlt_rep.pars
pars
Definition: make_hlt_rep.py:90
MuonCalib::RtFromPoints::getTrChebyshev
static std::unique_ptr< ITrRelation > getTrChebyshev(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of r-t data points into a t(r) relation expressed as a series of chebychev polynomial...
Definition: RtFromPoints.cxx:82
LArSamples::FitterData::fitter
const ShapeFitter * fitter
Definition: ShapeFitter.cxx:23
GeoModel::TransientConstSharedPtr< IRtRelation >
MuonCalib::RtFromPoints::getRtRelationLookUp
static std::unique_ptr< IRtRelation > getRtRelationLookUp(const std::vector< SamplePoint > &sample_points)
Definition: RtFromPoints.cxx:112
getMessageSvc.h
singleton-like access to IMessageSvc via open function and helper
MuonCalib::ChebyshevPolynomial
Definition: ChebyshevPolynomial.h:17
MuonCalib::normalizeDomain
std::vector< SamplePoint > normalizeDomain(const std::vector< SamplePoint > &dataPoints)
Normalizes the domain of the samples points to the interval -1 to 1.
Definition: SamplePointUtils.cxx:76
CheckAppliedSFs.bin_width
bin_width
Definition: CheckAppliedSFs.py:242
MuonCalib::RtFromPoints::legendreFit
static CalibFunc::ParVec legendreFit(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Executes the fit of Legendre polynomials to the data points.
Definition: RtFromPoints.cxx:37
RtRelationLookUp.h
MuonCalib::interval
std::pair< double, double > interval(const std::vector< SamplePoint > &points)
Returns the interval covered by the sample points.
Definition: SamplePointUtils.cxx:92
MuonCalib::RtFromPoints::simplePolyFit
static CalibFunc::ParVec simplePolyFit(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Exectues the fit of simple monomials to the data points.
Definition: RtFromPoints.cxx:66
RtSpline.h
MuonCalib::RtSpline::tLower
virtual double tLower() const override final
get the lower drift-time bound
Definition: RtSpline.cxx:58
MuonCalib::RtFromPoints::chebyFit
static CalibFunc::ParVec chebyFit(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Executes the fit of chebychev polynomials to the data points.
Definition: RtFromPoints.cxx:51
MuonCalib::RtFromPoints::getResoChebyshev
static std::unique_ptr< IRtResolution > getResoChebyshev(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of reso - t into a reso(t) relation expressed as a series of chebychev polynomials.
Definition: RtFromPoints.cxx:85
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
MuonCalib::BaseFunctionFitter
Definition: BaseFunctionFitter.h:39
MuonCalib::RtFromPoints::getRtChebyshev
static std::unique_ptr< IRtRelation > getRtChebyshev(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of r-t data points into a r(t) relation expressed as a series of chebychev polynomial...
Definition: RtFromPoints.cxx:79
lumiFormat.i
int i
Definition: lumiFormat.py:85
MuonCalib::RtFromPoints::getRtSimplePoly
static std::unique_ptr< IRtRelation > getRtSimplePoly(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of r(t) data points into a r(t) relation expressed as a series of elementary monomoni...
Definition: RtFromPoints.cxx:101
RtFromPoints.h
mc.order
order
Configure Herwig7.
Definition: mc.Herwig7_Dijet.py:12
MuonCalib::RtFromPoints::getTrSimplePoly
static std::unique_ptr< ITrRelation > getTrSimplePoly(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of t(r) data points into a t(r) relation expressed as a series of elementary monomoni...
Definition: RtFromPoints.cxx:104
MuonCalib::RtSpline::tUpper
virtual double tUpper() const override final
get the upper drift-time bound
Definition: RtSpline.cxx:60
MuonCalib::LegendrePolynomial
Definition: LegendrePolynomial.h:19
MuonCalib
CscCalcPed - algorithm that finds the Cathode Strip Chamber pedestals from an RDO.
Definition: CscCalcPed.cxx:22
MuonCalib::CalibFunc::ParVec
std::vector< double > ParVec
Definition: CalibFunc.h:35
RtSimplePolynomial.h
TrLegendre.h
MuonCalib::RtFromPoints::getTrLegendre
static std::unique_ptr< ITrRelation > getTrLegendre(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of t(r) data points into a t(r) relation expressed as a series of legendre polynomial...
Definition: RtFromPoints.cxx:98
BaseFunctionFitter.h
MuonCalib::SimplePolynomial
Definition: SimplePolynomial.h:13
RtResolutionChebyshev.h
MuonCalib::RtSpline
Definition: RtSpline.h:34
SimplePolynomial.h
MuonCalib::RtFromPoints::getRtLegendre
static std::unique_ptr< IRtRelation > getRtLegendre(const std::vector< SamplePoint > &dataPoints, const unsigned order)
Converts a list of r-t data points into a r(t) relation expressed as a series of legendre polynomials...
Definition: RtFromPoints.cxx:95
ChebyshevPolynomial.h
RtChebyshev.h
TrChebyshev.h
SamplePointUtils.h
PolygonBase.h
RtLegendre.h
LegendrePolynomial.h
RadiusResolutionChebyshev.h
fitman.k
k
Definition: fitman.py:528
MuonCalib::resoFromRadius
std::vector< SamplePoint > resoFromRadius(const std::vector< SamplePoint > &points, const double uncert)
Creates a new vector of sample points where the x2 coordinate becomes the x1 coordinate and the resol...
Definition: SamplePointUtils.cxx:39
MuonCalib::RtSpline::radius
virtual double radius(double t) const override final
get the radius corresponding to the drift time t; 0 or 14.6 is returned if t is outside the range
Definition: RtSpline.cxx:45