diff --git a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/MdtCalibFitters/QuasianalyticLineReconstruction.h b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/MdtCalibFitters/QuasianalyticLineReconstruction.h index 30562bcca4dfa1a1b3d72ab076818b7cb29ab208..5b79249f1f5de7fe8ce66dafcee83dd065383dd6 100644 --- a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/MdtCalibFitters/QuasianalyticLineReconstruction.h +++ b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/MdtCalibFitters/QuasianalyticLineReconstruction.h @@ -119,6 +119,9 @@ public: void setTimeOut(const double & time_out); ///< set the time-out for the track ///< finding to time_out (in seconds) + void setMaxRadius(const double& maxR); + ///< set the max innerRadius, default for MDT + ///< sMDT 7.1mm, MDT 14.6mm // methods required by the base class "IMdtSegmentFitter" // bool fit(MuonCalibSegment & r_segment) const; ///< reconstruction of the track using diff --git a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/src/QuasianalyticLineReconstruction.cxx b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/src/QuasianalyticLineReconstruction.cxx index ab6263e052121f2bfce8568db2d053ceda2bcbea..0828dcb44f58552d2d8d8ea7ab5f05beab96b0a4 100644 --- a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/src/QuasianalyticLineReconstruction.cxx +++ b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibFitters/src/QuasianalyticLineReconstruction.cxx @@ -80,7 +80,7 @@ void QuasianalyticLineReconstruction::init(const double & r_road_width) { //:: SET THE MAXIMUM RADIUS :: //:::::::::::::::::::::::::::: - m_r_max = 15.0; + m_r_max = 14.6; //:::::::::::::::::::::::: //:: SET THE ROAD WIDTH :: @@ -547,6 +547,11 @@ void QuasianalyticLineReconstruction::setTimeOut(const double & time_out) { } +void QuasianalyticLineReconstruction::setMaxRadius(const double& maxR) { + if (maxR<7 || maxR>15) throw std::runtime_error(Form("File: %s, Line: %d\nQuasianalyticLineReconstruction::setMaxRadius() - given radius %.3f not supported, neither MDT nor sMDT", __FILE__, __LINE__, maxR)); + m_r_max = maxR; +} + //***************************************************************************** //::::::::::::::::::::::: @@ -865,10 +870,8 @@ bool QuasianalyticLineReconstruction::fit(MuonCalibSegment & r_segment, m_track_hits[k]->localPosition().z()), xhat, null, null); double d(std::abs(m_track.signDistFrom(wire))); - m_chi2 = m_chi2+ - std::pow(d-m_track_hits[k]->driftRadius(), - 2)/ - m_track_hits[k]->sigma2DriftRadius(); + double r(std::abs(m_track_hits[k]->driftRadius())); + m_chi2 = m_chi2+std::pow(d-r,2)/m_track_hits[k]->sigma2DriftRadius(); } r_segment.set( m_chi2/static_cast<double>(m_nb_track_hits-2), diff --git a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/MdtCalibRt/RtCalibrationCurved.h b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/MdtCalibRt/RtCalibrationCurved.h index adc89c73aa37874a47b1b6201b68f99df3d3576d..630df6a9c8c3a14bd268fc0da1f7acd59181f09d 100644 --- a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/MdtCalibRt/RtCalibrationCurved.h +++ b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/MdtCalibRt/RtCalibrationCurved.h @@ -1,22 +1,10 @@ /* - Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration */ -//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ -// 10.05.2008, AUTHOR: OLIVER KORTNER -// Modified: 01.03.2009 by O. Kortner, parabolic extrapolation added. -// 15.03.2009 by O. Kortner, smoothing added. -// 18.03.2009 by O. Kortner, method performParabolicExtrapolation -// added. -//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ - #ifndef MuonCalib_RtCalibrationCurvedH #define MuonCalib_RtCalibrationCurvedH -//::::::::::::::::::::::::::::::: -//:: CLASS RtCalibrationCurved :: -//::::::::::::::::::::::::::::::: - /// \class RtCalibrationCurved /// /// This class performs the autocalibration of MDT chambers from the residuals @@ -28,14 +16,10 @@ /// /// \date 10.05.2008 -//:::::::::::::::::: -//:: HEADER FILES :: -//:::::::::::::::::: - -// STL // #include <vector> #include <string> #include <list> +#include <memory> // CLHEP // #include "CLHEP/Matrix/SymMatrix.h" @@ -252,7 +236,7 @@ private: IRtRelation * m_rt_new; // r-t as determined by the autocalibration RtCalibrationOutput * m_output; // class holding the results of the // autocalibration - MultilayerRtDifference * m_multilayer_rt_difference; + std::unique_ptr<MultilayerRtDifference> m_multilayer_rt_difference; // curved-segment fitting // double m_r_max; // maximum value for accepted drift radii CurvedPatRec *m_tracker; // curved segment finder (used for track fitting) @@ -290,14 +274,17 @@ private: // describing the curved line // control histograms // - TFile *m_tfile; // ROOT file - TH1F *m_cut_evolution; // cut evolution histogram - TH1F *m_nb_segment_hits; // number of hits on the segments - TH1F *m_pull_initial; // initial pull distribution - TH1F *m_pull_final; // final pull distribution after convergence - TH2F *m_residuals_initial; // initial residual distribution - TH2F *m_residuals_final; // final residual distribution after convergence - + std::unique_ptr<TFile> m_tfile; // ROOT file + std::unique_ptr<TH1F> m_cut_evolution; // cut evolution histogram + std::unique_ptr<TH1F> m_nb_segment_hits; // number of hits on the segments + std::unique_ptr<TH1F> m_pull_initial; // initial pull distribution + std::unique_ptr<TH1F> m_pull_final; // final pull distribution after convergence + std::unique_ptr<TH2F> m_residuals_initial; // initial residual distribution + std::unique_ptr<TH2F> m_residuals_initial_all; // initial residual distribution before convergence + std::unique_ptr<TH2F> m_residuals_final; // final residual distribution after convergence + std::unique_ptr<TH2F> m_driftTime_initial; // final residual distribution after convergence + std::unique_ptr<TH2F> m_driftTime_final; // final residual distribution after convergence + std::unique_ptr<TH2F> m_adc_vs_residual_final; // final residual distribution after convergence // private methods // void init(const double & rt_accuracy, const unsigned int & func_type, diff --git a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/src/RtCalibrationCurved.cxx b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/src/RtCalibrationCurved.cxx index 5db8b5eeb1a0e74e4bb659b953a90837152f1c95..009f4ca9c93334e59d2faf3ce666d8d6d4c4589a 100644 --- a/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/src/RtCalibrationCurved.cxx +++ b/MuonSpectrometer/MuonCalib/MdtCalib/MdtCalibRt/src/RtCalibrationCurved.cxx @@ -109,16 +109,8 @@ RtCalibrationCurved::~RtCalibrationCurved(void) { if (m_base_function!=0) { delete m_base_function; } - delete m_multilayer_rt_difference; - if (m_tfile!=0) { + if (m_tfile) { m_tfile->Write(); - delete m_cut_evolution; - delete m_nb_segment_hits; - delete m_pull_initial; - delete m_pull_final; - delete m_residuals_initial; - delete m_residuals_final; - delete m_tfile; } } @@ -207,32 +199,19 @@ void RtCalibrationCurved::switch_on_control_histograms(const std::string & file_ ///////////////////////////////////////////// m_control_histograms = true; - m_tfile = new TFile(file_name.c_str(), "RECREATE"); + m_tfile = std::make_unique<TFile>(file_name.c_str(), "RECREATE"); - m_cut_evolution = new TH1F("m_cut_evolution", - "CUT EVOLUTION", - 11, -0.5, 10.5); + m_cut_evolution = std::make_unique<TH1F>("m_cut_evolution","CUT EVOLUTION",11, -0.5, 10.5); + m_nb_segment_hits = std::make_unique<TH1F>("m_nb_segment_hits","NUMBER OF HITS ON THE REFITTED SEGMENTS",11, -0.5, 10.5); + m_pull_initial = std::make_unique<TH1F>("m_pull_initial","INITIAL PULL DISTRIBUTION",200, -5.05, 5.05); + m_residuals_initial = std::make_unique<TH2F>("m_residuals_initial","INITIAL OF THE REFITTED SEGMENTS",100, -0.5, 15.0, 300, -1.5, 1.5); + m_residuals_initial_all = std::make_unique<TH2F>("m_residuals_initial_all","INITIAL OF THE REFITTED SEGMENTS BEFORE CONVERGENCE",300, -15.0, 15.0, 1000, -5, 5); + m_residuals_final = std::make_unique<TH2F>("m_residuals_final","FINAL OF THE REFITTED SEGMENTS",100, -0.5, 15.0, 300, -1.5, 1.5); + m_driftTime_final = std::make_unique<TH2F>("m_driftTime_final","FINAL DRIFTTIME OF THE REFITTED SEGMENTS",300, -15.0, 15.0, 300,-15.0, 15.0); + m_driftTime_initial = std::make_unique<TH2F>("m_driftTime_initial","FINAL DRIFTTIME OF THE REFITTED SEGMENTS",300, -15.0, 15.0, 1300, -100, 1200); + m_adc_vs_residual_final = std::make_unique<TH2F>("m_adc_resi_initial","FINAL ADC VS RESIDUAL OF THE REFITTED SEGMENTS",1350, 0, 1350, 300,-15,15); - m_nb_segment_hits = new TH1F("m_nb_segment_hits", - "NUMBER OF HITS ON THE REFITTED SEGMENTS", - 11, -0.5, 10.5); - - m_pull_initial = new TH1F("m_pull_initial", - "INITIAL PULL DISTRIBUTION", - 200, -5.05, 5.05); -// m_pull_final = new TH1F("m_pull_final", -// "FINAL PULL DISTRIBUTION", -// 200, -5.05, 5.05); - - m_residuals_initial = new TH2F("m_residuals_initial", - "INITIAL OF THE REFITTED SEGMENTS", - 100, -0.5, 15.0, 300, -1.5, 1.5); - m_residuals_final = new TH2F("m_residuals_final", - "FINAL OF THE REFITTED SEGMENTS", - 100, -0.5, 15.0, 300, -1.5, 1.5); - - delete m_multilayer_rt_difference; - m_multilayer_rt_difference=new MultilayerRtDifference(10000, m_tfile); + m_multilayer_rt_difference=std::make_unique<MultilayerRtDifference>(10000, m_tfile.get()); return; } @@ -243,10 +222,9 @@ void RtCalibrationCurved::switch_on_control_histograms(const std::string & file_ //:::::::::::::::::::::::::::::::::::::::::: void RtCalibrationCurved::switch_off_control_histograms(void) { m_control_histograms = false; - if (m_tfile!=0) { + if (m_tfile) { m_tfile->Write(); - delete m_multilayer_rt_difference; - m_multilayer_rt_difference=new MultilayerRtDifference(10000); + m_multilayer_rt_difference=std::make_unique<MultilayerRtDifference>(10000); } return; @@ -309,6 +287,18 @@ void RtCalibrationCurved::noSmoothing(void) { //:::::::::::::::::::::::::::: const IMdtCalibrationOutput * RtCalibrationCurved::analyseSegments(const std::vector<MuonCalibSegment*> & seg) { const IRtRelation *tmp_rt(0); + +///////////////////// +// Initial RESIDUALS // +///////////////////// + for (unsigned int k=0; k<seg.size(); k++) { + for (unsigned int l=0; l<seg[k]->hitsOnTrack(); l++) { + double rr = (seg[k]->mdtHOT())[l]->driftRadius(); + double dd = (seg[k]->mdtHOT())[l]->signedDistanceToTrack(); + if(m_residuals_final) m_residuals_initial_all->Fill(std::abs(dd), std::abs(rr)-std::abs(dd)); + m_driftTime_initial->Fill(rr, dd); + } + } ////////////////////////// // AUTOCALIBRATION LOOP // @@ -352,13 +342,17 @@ const IMdtCalibrationOutput * RtCalibrationCurved::analyseSegments(const std::ve ////////////////////////////////////////////// if (!m_do_smoothing) { // final residuals // - double r, d; + double r(0); + double d(0); + double adc(0); for (unsigned int k=0; k<seg.size(); k++) { for (unsigned int l=0; l<seg[k]->hitsOnTrack(); l++) { - r = std::abs((seg[k]->mdtHOT())[l]->driftRadius()); - d = std::abs((seg[k]->mdtHOT())[l]->signedDistanceToTrack()); - if(m_residuals_final != NULL) - m_residuals_final->Fill(d, r-d, 1.0); + adc = (seg[k]->mdtHOT())[l]->adcCount(); + r = (seg[k]->mdtHOT())[l]->driftRadius(); + d = (seg[k]->mdtHOT())[l]->signedDistanceToTrack(); + if(m_residuals_final) m_residuals_final->Fill(std::abs(d), std::abs(r)-std::abs(d)); + m_driftTime_final->Fill(r, d); + m_adc_vs_residual_final->Fill(adc,std::abs(r)-std::abs(d)); } } return getResults(); @@ -419,10 +413,9 @@ const IMdtCalibrationOutput * RtCalibrationCurved::analyseSegments(const std::ve // final residuals // for (unsigned int k=0; k<seg.size(); k++) { for (unsigned int l=0; l<seg[k]->hitsOnTrack(); l++) { - double r = std::abs((seg[k]->mdtHOT())[l]->driftRadius()); - double d = std::abs((seg[k]->mdtHOT())[l]->signedDistanceToTrack()); - if(m_residuals_final != NULL) - m_residuals_final->Fill(d, r-d, 1.0); + double r = (seg[k]->mdtHOT())[l]->driftRadius(); + double d = (seg[k]->mdtHOT())[l]->signedDistanceToTrack(); + if(m_residuals_final) m_residuals_final->Fill(d, std::abs(r)-std::abs(d), 1.0); } } return getResults(); @@ -461,10 +454,12 @@ const IMdtCalibrationOutput * RtCalibrationCurved::analyseSegments(const std::ve ///////////////////// for (unsigned int k=0; k<seg.size(); k++) { for (unsigned int l=0; l<seg[k]->hitsOnTrack(); l++) { - double r = std::abs((seg[k]->mdtHOT())[l]->driftRadius()); - double d = std::abs((seg[k]->mdtHOT())[l]->signedDistanceToTrack()); - if(m_residuals_final != NULL) - m_residuals_final->Fill(d, r-d, 1.0); + double adc = (seg[k]->mdtHOT())[l]->adcCount(); + double r = (seg[k]->mdtHOT())[l]->driftRadius(); + double d = (seg[k]->mdtHOT())[l]->signedDistanceToTrack(); + if(m_residuals_final) m_residuals_final->Fill(d, std::abs(r)-std::abs(d)); + m_driftTime_final->Fill(r, d); + m_adc_vs_residual_final->Fill(adc,std::abs(r)-std::abs(d)); } } @@ -890,7 +885,7 @@ bool RtCalibrationCurved::analyse(const std::vector<MuonCalibSegment*> & seg) { // do not change the r-t and the endpoints // if (((k==0 || x_r[k].x2()<0.5) && m_fix_min) || - ((k==nb_points || x_r[k].x2()>14.1) && m_fix_max)) { + ((k==nb_points || x_r[k].x2()>m_r_max) && m_fix_max)) { r_corr = 0.0; x_r[k].set_error(0.01); } @@ -935,7 +930,7 @@ bool RtCalibrationCurved::analyse(const std::vector<MuonCalibSegment*> & seg) { max_k = rt_param.size()-1; } for (unsigned int k=min_k; k<max_k; k++) { - x = (rt_param[k] - 7.6)/7.6; + x = (rt_param[k] - 0.5*m_r_max)/(0.5*m_r_max); r_corr = m_alpha[0]; for (unsigned int l=1; l<m_order; l++) { r_corr = r_corr+m_alpha[l]* @@ -1057,16 +1052,9 @@ void RtCalibrationCurved::init(const double & rt_accuracy, ///////////////////////////// // RESET PRIVATE VARIABLES // ///////////////////////////// - m_rt=NULL; + m_rt=nullptr; m_r_max = 15.0*CLHEP::mm; m_control_histograms = false; - m_tfile = 0; - m_cut_evolution = 0; - m_nb_segment_hits = 0; - m_pull_initial = 0; - m_pull_final = 0; - m_residuals_initial = 0; - m_residuals_final = 0; m_nb_segments = 0; m_nb_segments_used = 0; m_iteration = 0; @@ -1081,7 +1069,7 @@ void RtCalibrationCurved::init(const double & rt_accuracy, m_fix_max = fix_max; m_max_it = std::abs(max_it); m_do_multilayer_rt_scale=do_multilayer_rt_scale; - m_multilayer_rt_difference = new MultilayerRtDifference(10000); + m_multilayer_rt_difference = std::make_unique<MultilayerRtDifference>(10000); m_tracker = new CurvedPatRec; diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MuonCalibStandAloneBase/MuonCalibStandAloneBase/MdtStationT0Container.h b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MuonCalibStandAloneBase/MuonCalibStandAloneBase/MdtStationT0Container.h index ef5752bcdefa27a62805eae7cb77fb85be4d64a6..3bc432c61db075dda98ceae755bd877156e3b3e7 100644 --- a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MuonCalibStandAloneBase/MuonCalibStandAloneBase/MdtStationT0Container.h +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MuonCalibStandAloneBase/MuonCalibStandAloneBase/MdtStationT0Container.h @@ -1,11 +1,7 @@ /* - Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration */ -//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ -// 01.02.2007, AUTHOR: OLIVER KORTNER -//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ - #ifndef MuonCalib_MdtStationT0ContainerH #define MuonCalib_MdtStationT0ContainerH @@ -24,11 +20,6 @@ /// /// \date 01.02.2007 -//:::::::::::::::::: -//:: HEADER FILES :: -//:::::::::::::::::: - -// STL // #include <vector> #include <string> @@ -39,16 +30,16 @@ class MdtStationT0Container { public: // Constructors // MdtStationT0Container(void) : -m_t0(2, std::vector< std::vector<double> >(4, std::vector<double>(72, 9e9))), -m_adc(2, std::vector< std::vector<double> >(4, std::vector<double>(72, 9e9))), +m_t0(2, std::vector< std::vector<double> >(4, std::vector<double>(78, 9e9))), +m_adc(2, std::vector< std::vector<double> >(4, std::vector<double>(78, 9e9))), m_t0_loaded(false) { } ///< Default constructor. MdtStationT0Container(const std::string & file_name): -m_t0(2, std::vector< std::vector<double> >(4, std::vector<double>(72, 9e9))), -m_adc(2, std::vector< std::vector<double> >(4, std::vector<double>(72, 9e9))), +m_t0(2, std::vector< std::vector<double> >(4, std::vector<double>(78, 9e9))), +m_adc(2, std::vector< std::vector<double> >(4, std::vector<double>(78, 9e9))), m_t0_loaded(false) { readT0File(file_name);