ATLAS Offline Software
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions | Private Attributes | List of all members
MuonCalib::AdaptiveResidualSmoothing Class Reference

#include <AdaptiveResidualSmoothing.h>

Collaboration diagram for MuonCalib::AdaptiveResidualSmoothing:

Public Member Functions

 AdaptiveResidualSmoothing ()
 Default constructor.
void clear ()
 clear the memory of the class
void addResidual (const double radius, const double residual)
 add the residual at the given radius
bool addResidualsFromSegment (MuonCalibSegment &seg, bool curved, double road_width)
 reconstruct the given segment and store the residuals; if curved is true a curved segment fit is performed, otherwise a straight segment is fitted to the drift radii; the user must set the road width used in the pattern recognition; returns true in case of success
RtRelationLookUp performSmoothing (const IRtRelation &rt_rel, unsigned int nb_entries_per_bin, bool fix_t_min, bool fix_t_max)
 use the stored residuals to improve the given r-t relationship to give smoother residuals; the user has to set the number of entries per radial bin; the user can request that the radii for the minimum and maximum drift time are untouched
RtRelationLookUp performSmoothing (const IRtRelation &rt_rel, const bool &fix_t_min, const bool &fix_t_max)
 use the stored residuals to improve the given r-t relationship to give smoother residuals; the user can request that the radii for the minimum and maximum drift time are untouched

Private Member Functions

double t_from_r (const IRtRelation &rt_rel, const double r)

Private Attributes

std::vector< DataPoint > m_residual_point
StraightPatRec m_sfitter
CurvedPatRec m_cfitter

Detailed Description

Definition at line 48 of file AdaptiveResidualSmoothing.h.

Constructor & Destructor Documentation

◆ AdaptiveResidualSmoothing()

AdaptiveResidualSmoothing::AdaptiveResidualSmoothing ( )

Default constructor.

Definition at line 25 of file AdaptiveResidualSmoothing.cxx.

25{}

Member Function Documentation

◆ addResidual()

void AdaptiveResidualSmoothing::addResidual ( const double radius,
const double residual )

add the residual at the given radius

Definition at line 27 of file AdaptiveResidualSmoothing.cxx.

27 {
28 // make the data point //
29 CLHEP::HepVector point(2);
30 point[0] = radius;
31 point[1] = residual;
32 DataPoint residual_point(point, 0);
33
34 // store the data point //
35 m_residual_point.push_back(residual_point);
36
37 return;
38}
ParticleGun_SamplingFraction.radius
radius
Definition ParticleGun_SamplingFraction.py:96

◆ addResidualsFromSegment()

bool AdaptiveResidualSmoothing::addResidualsFromSegment ( MuonCalibSegment & seg,
bool curved,
double road_width )

reconstruct the given segment and store the residuals; if curved is true a curved segment fit is performed, otherwise a straight segment is fitted to the drift radii; the user must set the road width used in the pattern recognition; returns true in case of success

Definition at line 39 of file AdaptiveResidualSmoothing.cxx.

39 {
41 // VARIABLES //
43
44 CLHEP::HepVector point(2); // auxiliary point
45 IMdtPatRecFitter *fitter = nullptr;
46
48 // RESIDUAL DETERMINATION AND STORAGE //
50
51 // perform fit //
52 if (curved) {
53 // set road width /
54 m_cfitter.setRoadWidth(road_width);
55
56 // set fitter pointer //
57 fitter = &m_cfitter;
58
59 } else {
60 // set road width /
61 m_sfitter.setRoadWidth(road_width);
62
63 // set fitter pointer //
64 fitter = &m_sfitter;
65 }
66
67 fitter->SetRefineSegmentFlag(false);
68
69 CurvedLine curved_track;
70 MTStraightLine straight_track;
71 IMdtSegmentFitter::HitSelection selection(seg.mdtHitsOnTrack(), 0);
72
73 // perform the track fit; return in case of failure //
74
75 // reject bad fits //
76 if (curved) {
77 if (!m_cfitter.fit(seg, selection, curved_track)) return false;
78 if (curved_track.chi2PerDegreesOfFreedom() > 50) { return false; }
79 } else {
80 if (!m_sfitter.fit(seg, selection, straight_track)) return false;
81 if (straight_track.chi2PerDegreesOfFreedom() > 50) { return false; }
82 }
83
84 // store the residuals //
85 const MuonCalibSegment::MdtHitVec &trackHits = curved ? curved_track.trackHits() : straight_track.trackHits();
86 for (const auto & trackHit : trackHits) {
87 point[0] = trackHit->driftRadius();
88 point[1] = trackHit->driftRadius() - std::abs(trackHit->signedDistanceToTrack());
89 m_residual_point.emplace_back(point, 0);
90 }
91
92 return true;
93}
MuonCalib::CurvedLine::chi2PerDegreesOfFreedom
double chi2PerDegreesOfFreedom() const
Return chi2 / number of TrackHits - 3.
Definition CurvedLine.cxx:186
MuonCalib::CurvedLine::trackHits
const MdtHitVec & trackHits() const
Definition CurvedLine.cxx:189
MuonCalib::IMdtSegmentFitter::HitSelection
std::vector< unsigned int > HitSelection
Definition IMdtSegmentFitter.h:32
MuonCalib::MTStraightLine::trackHits
const MdtHitVec & trackHits() const
Definition MTStraightLine.cxx:142
MuonCalib::MTStraightLine::chi2PerDegreesOfFreedom
double chi2PerDegreesOfFreedom() const
Return chi2 / number of TrackHits - 2.
Definition MTStraightLine.cxx:139
MuonCalib::MuonCalibSegment::mdtHitsOnTrack
unsigned int mdtHitsOnTrack() const
retrieve the number of MdtCalibHitBase s assigned to this segment
Definition MuonCalibSegment.cxx:146
MuonCalib::MuonCalibSegment::MdtHitVec
std::vector< MdtHitPtr > MdtHitVec
Definition MuonCalibSegment.h:45
selection
const std::string selection
Definition fbtTestBasics.cxx:75
LArSamples::FitterData::fitter
const ShapeFitter * fitter
Definition ShapeFitter.cxx:23

◆ clear()

void AdaptiveResidualSmoothing::clear ( )

clear the memory of the class

Definition at line 26 of file AdaptiveResidualSmoothing.cxx.

26{ m_residual_point.clear(); }

◆ performSmoothing() [1/2]

RtRelationLookUp AdaptiveResidualSmoothing::performSmoothing ( const IRtRelation & rt_rel,
const bool & fix_t_min,
const bool & fix_t_max )

use the stored residuals to improve the given r-t relationship to give smoother residuals; the user can request that the radii for the minimum and maximum drift time are untouched

Definition at line 178 of file AdaptiveResidualSmoothing.cxx.

178 {
180 // VARIABLES //
182
183 unsigned int nb_bins; // number of bins in r
184 unsigned int min_nb_entries_per_bin(20); // minimum number of
185 // entries per bin
186 double step_size; // real step size [mm]
187 double r_max(16.0); // maximum drift radius
188 double aux_rad; // auxiliary radius
189 double aux_corr; // auxiliary correction value
190 unsigned int bin_index;
191
193 // DETERMINE THE NUMBER OF BINS //
195
196 // take the sqrt from the number of residuals to determine the binsize.
197 // for low stat, a larger statistical error on the bins is allowed
198 // for example:
199 // 1000 segments -> 7.6% and 35 bins
200 // 5000 segments -> 5.1% and 64 bins
201 // 10000 segments-> 4.1% and 91 bins
202 nb_bins = static_cast<unsigned int>(std::sqrt(m_residual_point.size() / 6.));
203
204 // calculate min_nb_per_bin
205 if (nb_bins == 0) { nb_bins = 1; }
206 min_nb_entries_per_bin = m_residual_point.size() / nb_bins;
207
209 // RETURN IF THERE ARE NO RESIDUALS //
211
212 if (m_residual_point.size() == 0) {
213 MsgStream log(Athena::getMessageSvc(), "AdaptiveResidualSmoothing");
214 log << MSG::WARNING << "performSmoothing() - No residuals are stored. no correction applied to r(t)!" << endmsg;
215 }
216
217 step_size = r_max / static_cast<double>(nb_bins);
218
220 // DETERMINE THE AVERAGE RESIDUAL VALUES IN THE RADIAL BINS //
222
223 // vector for the correction function //
224 std::vector<SamplePoint> corr(nb_bins);
225 std::vector<double> radii(nb_bins);
226
227 // sort the residuals point in the residual value for a simple outlyer //
228 // rejection performed later //
229 for (auto & k : m_residual_point) { k.setReferenceComponent(1); }
230 sort(m_residual_point.begin(), m_residual_point.end());
231
232 // auxiliary data point arrays //
233 std::vector<std::vector<const DataPoint *> > sample_points_per_r_bin(nb_bins);
234
235 // group data points in r bins //
236 for (auto & k : m_residual_point) {
237 if (std::abs(k.dataVector()[0]) >= r_max) { continue; }
238 bin_index = static_cast<unsigned int>(std::abs(k.dataVector()[0]) / step_size);
239 sample_points_per_r_bin[bin_index].push_back(&k);
240 }
241
242 // loop over the r bins and determine the r-t corrections //
243 // static ofstream outfile("tmp.txt");
244 // outfile << endl;
245 for (unsigned int bin = 0; bin < nb_bins; bin++) {
246 aux_rad = 0.0;
247 aux_corr = 0.0;
248
249 // remove the first and last 1% of the residuals as simple outlyer removal //
250 unsigned int k_min(static_cast<unsigned int>(0.01 * sample_points_per_r_bin[bin].size()));
251 unsigned int k_max(static_cast<unsigned int>(0.99 * sample_points_per_r_bin[bin].size()));
252 // outfile << bin << "\t" << (0.5+bin)*step_size << "\t"
253 // << k_min << "\t" << k_max << endl;
254
255 // loop over residuals in a bin
256 for (unsigned int k = k_min; k < k_max; k++) {
257 aux_rad = aux_rad + sample_points_per_r_bin[bin][k]->dataVector()[0];
258 aux_corr = aux_corr + sample_points_per_r_bin[bin][k]->dataVector()[1];
259 }
260 if (k_max - k_min < 0.25 * min_nb_entries_per_bin) {
261 radii[bin] = (0.5 + bin) * step_size;
262 corr[bin] = SamplePoint(radii[bin], 0.0, 1.0);
263 } else {
264 radii[bin] = aux_rad / static_cast<double>(k_max - k_min);
265 corr[bin] = SamplePoint(radii[bin], aux_corr / static_cast<double>(k_max - k_min), 1.0);
266 }
267 }
268
269 // correction polygon //
270 PolygonBase polygon(radii);
271 BaseFunctionFitter fitter(nb_bins);
272 fitter.fit_parameters(corr, 1, corr.size(), polygon);
273
274 // create output r-t relationship //
275 std::vector<double> rt_params;
276 rt_params.push_back(rt_rel.tLower());
277 rt_params.push_back(0.01 * (rt_rel.tUpper() - rt_rel.tLower()));
278 // ofstream outfile("out2.txt");
279 for (double t = rt_rel.tLower(); t <= rt_rel.tUpper(); t = t + rt_params[1]) {
280 double delta_r(0.0);
281 for (unsigned int l = 0; l < radii.size(); l++) {
282 delta_r = delta_r + fitter.coefficients()[l] * polygon.value(l, rt_rel.radius(t));
283 }
284 if (rt_rel.radius(t) > r_max) { delta_r = 0.0; }
285 if (fix_t_min && (rt_rel.radius(t) < 0.5)) { delta_r = 0.0; }
286 if (fix_t_max && (rt_rel.radius(t) > 14.1)) { delta_r = 0.0; }
287 rt_params.push_back(rt_rel.radius(t) - delta_r);
288 }
289
290 RtRelationLookUp improved_rt(rt_params);
291 // for (double t=rt_rel.tLower(); t<=rt_rel.tUpper(); t=t+rt_params[1]) {
292 // outfile << rt_rel.radius(t) << "\t" << improved_rt.radius(t)
293 // << "\t" << endl;
294 // }
295
296 return improved_rt;
297}
endmsg
#define endmsg
Definition AnalysisConfig_Ntuple.cxx:63
sort
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
Definition DVL_algorithms.h:554
MuonCalib::IRtRelation::tLower
virtual double tLower() const =0
Returns the lower time covered by the r-t.
MuonCalib::IRtRelation::radius
virtual double radius(double t) const =0
returns drift radius for a given time
MuonCalib::IRtRelation::tUpper
virtual double tUpper() const =0
Returns the upper time covered by the r-t.
Athena::getMessageSvc
IMessageSvc * getMessageSvc(bool quiet=false)
Definition getMessageSvc.cxx:20
PrintTrkAnaSummary.l
l
Printing final latex table to .tex output file.
Definition PrintTrkAnaSummary.py:370

◆ performSmoothing() [2/2]

RtRelationLookUp AdaptiveResidualSmoothing::performSmoothing ( const IRtRelation & rt_rel,
unsigned int nb_entries_per_bin,
bool fix_t_min,
bool fix_t_max )

use the stored residuals to improve the given r-t relationship to give smoother residuals; the user has to set the number of entries per radial bin; the user can request that the radii for the minimum and maximum drift time are untouched

Definition at line 95 of file AdaptiveResidualSmoothing.cxx.

96 {
98 // VARIABLES //
100
101 ConstantContentBinMaker bin_maker(m_residual_point, 0.001);
102 std::vector<unsigned int> ref_coord(1);
103
105 // CHECKS //
107
108 if (m_residual_point.size() == 0) {
109 throw std::runtime_error(
110 Form("File: %s, Line: %d\nAdaptiveResidualSmoothing::performSmoothing - No residuals stored.", __FILE__, __LINE__));
111 }
112
113 if (m_residual_point.size() < nb_entries_per_bin) {
114 throw std::runtime_error(
115 Form("File: %s, Line: %d\nAdaptiveResidualSmoothing::performSmoothing - Not enough residuals stored.", __FILE__, __LINE__));
116 }
117
119 // PERFORM BINNING //
121
122 ref_coord[0] = 0;
123 bin_maker.binDataPoints(nb_entries_per_bin, ref_coord);
124
126 // DETERMINE r-t COORECTION //
128
129 // get a polygon through the bin centres //
130 std::vector<double> rad(bin_maker.getBins().size());
131 std::vector<SamplePoint> corr(bin_maker.getBins().size());
132 for (unsigned int k = 0; k < bin_maker.getBins().size(); k++) {
133 rad[k] = (bin_maker.getBins()[k])->centreOfGravity()[0];
134 corr[k].set_x1(rad[k]);
135 corr[k].set_x2((bin_maker.getBins()[k])->centreOfGravity()[1]);
136 corr[k].set_error(1.0);
137 }
138 std::sort(rad.begin(), rad.end());
139
140 std::vector<double> radi;
141 radi.push_back(rad[0]);
142 unsigned int counter(0);
143 for (unsigned int k = 1; k < rad.size(); k++) {
144 if (rad[counter] < rad[k]) {
145 radi.push_back(rad[k]);
146 counter++;
147 }
148 }
149
150 PolygonBase polygon(radi);
151 BaseFunctionFitter fitter(bin_maker.getBins().size());
152 fitter.fit_parameters(corr, 1, corr.size(), polygon);
153
154 // create an improved r-t relationship //
155 std::vector<double> rt_params;
156 rt_params.push_back(rt_rel.tLower());
157 rt_params.push_back(0.01 * (rt_rel.tUpper() - rt_rel.tLower()));
158 // ofstream outfile("out2.txt");
159 for (double t = rt_rel.tLower(); t <= rt_rel.tUpper(); t = t + rt_params[1]) {
160 double delta_r(0.0);
161 for (unsigned int l = 0; l < radi.size(); l++) {
162 delta_r = delta_r + fitter.coefficients()[l] * polygon.value(l, rt_rel.radius(t));
163 }
164 if (fix_t_min && (rt_rel.radius(t) < 0.5)) { delta_r = 0.0; }
165 if (fix_t_max && (rt_rel.radius(t) > 14.1)) { delta_r = 0.0; }
166 rt_params.push_back(rt_rel.radius(t) - delta_r);
167 }
168
169 RtRelationLookUp improved_rt(rt_params);
170 // for (double t=rt_rel.tLower(); t<=rt_rel.tUpper(); t=t+rt_params[1]) {
171 // outfile << rt_rel.radius(t) << "\t" << improved_rt.radius(t)
172 // << "\t" << endl;
173 // }
174
175 return improved_rt;
176}
std::sort
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
Definition DVL_algorithms.h:554
test_pyathena.counter
counter
Definition test_pyathena.py:15

◆ t_from_r()

double AdaptiveResidualSmoothing::t_from_r ( const IRtRelation & rt_rel,
const double r )
private

Definition at line 298 of file AdaptiveResidualSmoothing.cxx.

298 {
300 // VARIABLES //
302
303 double precision(0.001); // spatial precision of the inversion
304 double t_max(rt_rel.tUpper()); // upper time search limit
305 double t_min(rt_rel.tLower()); // lower time search limit
306
308 // SEARCH FOR THE CORRESPONDING DRIFT TIME //
310
311 while (t_max - t_min > 0.1 && std::abs(rt_rel.radius(0.5 * (t_min + t_max)) - r) > precision) {
312 if (rt_rel.radius(0.5 * (t_min + t_max)) > r) {
313 t_max = 0.5 * (t_min + t_max);
314 } else {
315 t_min = 0.5 * (t_min + t_max);
316 }
317 }
318
319 return 0.5 * (t_min + t_max);
320}
r
int r
Definition globals.cxx:22

Member Data Documentation

◆ m_cfitter

CurvedPatRec MuonCalib::AdaptiveResidualSmoothing::m_cfitter
private

Definition at line 86 of file AdaptiveResidualSmoothing.h.

◆ m_residual_point

std::vector<DataPoint> MuonCalib::AdaptiveResidualSmoothing::m_residual_point
private

Definition at line 84 of file AdaptiveResidualSmoothing.h.

◆ m_sfitter

StraightPatRec MuonCalib::AdaptiveResidualSmoothing::m_sfitter
private

Definition at line 85 of file AdaptiveResidualSmoothing.h.

The documentation for this class was generated from the following files:
Generated on for ATLAS Offline Software by doxygen 1.14.0