diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationFileIOTool.h b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationFileIOTool.h new file mode 100644 index 0000000000000000000000000000000000000000..958562f1d3d0462e37860fcf27b9dddba297d8e6 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationFileIOTool.h @@ -0,0 +1,58 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef CalibrationFileIOTool_H +#define CalibrationFileIOTool_H + + +// MuonCalibStandAloneBase +#include "MuonCalibStandAloneBase/CalibrationIOTool.h" + +#include "MdtCalibUtils/RtDataFromFile.h" + +namespace MuonCalib { + +class RtRelation; + +class CalibrationFileIOTool : public AlgTool, virtual public CalibrationIOTool + { + public: + /** constructor*/ + CalibrationFileIOTool(const std::string& t, const std::string& n, const IInterface* p); + /** initialisation */ + inline StatusCode initialize() + { + return StatusCode::SUCCESS; + } + /** finalisation */ + inline StatusCode finalize() + { + std::cout<<"CalibrationFileIOTool Finalize"<<std::endl; + return StatusCode::SUCCESS; + } + /** write out t0 */ + StatusCode WriteT0(MdtTubeFitContainer* t0_output, const NtupleStationId & station_id, int /*iov_start*/, int /*iov_end*/); + /** write rt*/ + StatusCode WriteRt(const RtCalibrationOutput *rt_relation, const IRtResolution * resolution, const NtupleStationId & station_id, int /*iov_start*/, int /*iov_end*/, bool /*real_rt*/, bool /*real resolution*/); + /** load t0s*/ + StatusCode LoadT0(std::map<NtupleStationId, MdtStationT0Container *> &t0s, int /*iov_id*/); + /** load rts */ + StatusCode LoadRt(std::map<NtupleStationId, IRtRelation *> & rts, std::map<NtupleStationId, IRtResolution *> &res, int /*iov_id*/); + private: + //! path to calibration directory - job option + std::string m_calib_dir; + //! use fixed identifier if set tot true - job option + bool m_use_fixed_id; + //! create rt relation as lookup table if set tot true - job option + bool m_rt_lookup; + //! fill rt relation + inline bool fill_rt(RtDataFromFile::RtRelation *rt, const IRtRelation *new_rt, const MuonCalib::IRtResolution * resolut); + //! extract station identifier from file name + inline bool interpret_chamber_name(const std::string &nm, const char *prefix, std::string & station, int &eta, int & phi, int & ml) const; + //! create the rt relation and resolution + inline void read_rt_relation(const std::string & fname, std::map<NtupleStationId, IRtRelation *> & rts, std::map<NtupleStationId, IRtResolution *> &res,const MuonCalib::NtupleStationId & id); + }; + +} +#endif diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationOracleFileIOTool.h b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationOracleFileIOTool.h new file mode 100644 index 0000000000000000000000000000000000000000..7de65f55f91bca4831c1e04f992bc6a446dcbbf7 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/CalibrationOracleFileIOTool.h @@ -0,0 +1,46 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef CalibrationOracleFileIOTool_H +#define CalibrationOracleFileIOTool_H + + +// MuonCalibStandAloneBase +#include "MuonCalibStandAloneBase/CalibrationIOTool.h" + +#include "MdtCalibUtils/RtDataFromFile.h" + +namespace MuonCalib { + +class RtRelation; + +class CalibrationOracleFileIOTool : public AlgTool, virtual public CalibrationIOTool + { + public: + /** constructor*/ + CalibrationOracleFileIOTool(const std::string& t, const std::string& n, const IInterface* p); + /** initialisation */ + inline StatusCode initialize() + { + return StatusCode::SUCCESS; + } + /** finalisation */ + inline StatusCode finalize() + { + std::cout<<"CalibrationOracleFileIOTool Finalize"<<std::endl; + return StatusCode::SUCCESS; + } + /** write out t0 */ + StatusCode WriteT0(MdtTubeFitContainer* t0_output, const NtupleStationId & station_id, int iov_start, int iov_end); + /** write rt*/ + StatusCode WriteRt(const RtCalibrationOutput *rt_relation, const IRtResolution * resolution, const NtupleStationId & station_id, int iov_start, int iov_end, bool /*real_rt*/, bool /*real_resolution*/); + private: + //! path to calibration directory - job option + std::string m_calib_dir; + //! fill rt relation + inline bool fill_rt(RtDataFromFile::RtRelation *rt, const IRtRelation *new_rt, const MuonCalib::IRtResolution * resolut); + }; + +} +#endif diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibInputSvc.h b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibInputSvc.h new file mode 100644 index 0000000000000000000000000000000000000000..2d059575474a5c5434265498bb3617ebda3a2396 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibInputSvc.h @@ -0,0 +1,114 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef MdtCalibInputSvc_H +#define MdtCalibInputSvc_H + + +// c - c++ +#include "map" +#include "string" + + +//gaudi +#include "GaudiKernel/Service.h" +#include "GaudiKernel/IInterface.h" +#include "GaudiKernel/ToolHandle.h" + +#include "MuonCalibStandAloneBase/CalibrationIOTool.h" + +//MuonCalib +namespace MuonCalib { +class MdtStationT0Container; +class IRtRelation; +class BFieldCorFunc; +class IRtResolution; +} +#include "MuonCalibStandAloneBase/NtupleStationId.h" +class RegionSelectionSvc; + + +// interface to enable retrieving of a pointer to the singleton // +const InterfaceID IID_IMdtCalibInputSvc("MdtCalibInputSvc", 1, 0); + + +/** @calss MdtCalibInputSvc +Athena service which read calibration from text files and sorts them by station +id +@author rauscher@cern.ch +*/ + + +class MdtCalibInputSvc: public Service + { + public: +//============================================================================= + /** Service constructor */ + MdtCalibInputSvc(const std::string & name, ISvcLocator *svc_locator); + /** destructor */ + virtual ~MdtCalibInputSvc(); + static const InterfaceID& interfaceID() { return IID_IMdtCalibInputSvc; } + /** just some crazy atheta function */ + virtual StatusCode queryInterface(const InterfaceID& riid, + void** ppvUnknown); + /** service initalizer - reads files */ + virtual StatusCode initialize(void); + /** service finalize function - does nothing */ + virtual StatusCode finalize(void) ; + /** Get t0 container for Station */ + const MuonCalib::MdtStationT0Container * GetT0(const MuonCalib::NtupleStationId & id) const; + /** Get rt relation container */ + const MuonCalib::IRtRelation * GetRtRelation(const MuonCalib::NtupleStationId & id) const; + /** Get B field correction */ + const MuonCalib::BFieldCorFunc * GetBCorr(const MuonCalib::NtupleStationId & id) const; + /** Get Resolution */ + const MuonCalib::IRtResolution * GetResolution(const MuonCalib::NtupleStationId & id) const; + /** Get rt-Relation for calibration region */ + inline const MuonCalib::IRtRelation * GetRtRelation() const + { + return p_sel_region_rt; + } + /** Get b-field correction for calibratino region */ + inline const MuonCalib::BFieldCorFunc * GetBCorr() const + { + if(p_sel_region_b == NULL) + p_sel_region_b = GetBCorr(m_mean_station_id); + return p_sel_region_b; + } + /** Get resolution for calibration region */ + inline const MuonCalib::IRtResolution * GetResolution() const + { + return p_sel_region_res; + } +//============================================================================== + private: + //!calibration io tool to be used + ToolHandle<MuonCalib::CalibrationIOTool> m_calib_input_tool; + //! calibration data sorted by station id + std::map<MuonCalib::NtupleStationId, MuonCalib::MdtStationT0Container *> m_t0; + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtRelation *> m_rt_relation; + mutable std::map<MuonCalib::NtupleStationId, MuonCalib::BFieldCorFunc *> m_B_corr; + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtResolution *> m_spat_res; + /** pointer to region selection service */ + RegionSelectionSvc *p_reg_sel_svc; + /** create the b-field correction */ + inline bool create_b_field_correction(const MuonCalib::NtupleStationId &id) const; + inline const MuonCalib::BFieldCorFunc * findbfieldfun(const MuonCalib::NtupleStationId & id) const; + /** create mean rt relations, and resolutions for the selected calibration region */ + inline void create_mean_rts(); + inline StatusCode read_calib_input(); + /** rt relation - resolution - and correction function for the selected region - is averaged out of all matching rt relations*/ + const MuonCalib::IRtRelation * p_sel_region_rt; + mutable const MuonCalib::BFieldCorFunc * p_sel_region_b; + const MuonCalib::IRtResolution * p_sel_region_res; + /** statino id for mean rt */ + MuonCalib::NtupleStationId m_mean_station_id; + /** give warnings about missing calibratino only once per chamber */ + mutable std::set<MuonCalib::NtupleStationId> m_t0_warned; + mutable std::set<MuonCalib::NtupleStationId> m_rt_warned; + + }; + + +#endif diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibOutputDbSvc.h b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibOutputDbSvc.h new file mode 100644 index 0000000000000000000000000000000000000000..9bc6c403eba44ff83fe723ed1cdfc4e1c1c24b45 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/MdtCalibIOSvc/MdtCalibOutputDbSvc.h @@ -0,0 +1,191 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ +// 07.10.2006, AUTHOR: OLIVER KORTNER +// Modifications: 13.01.2007 by O. Kortner, singleton-like construction added +// to simplify use in stand-alon mode. +// 11.04.2007 by O. Kortner, new method which allows the user +// store a resolution function +// associated with an r-t function. +//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + +#ifndef MdtCalibOutputDbScvH +#define MdtCalibOutputDbScvH + +//::::::::::::::::::::::::::::::: +//:: CLASS MdtCalibOutputDvSvc :: +//::::::::::::::::::::::::::::::: + +/// \class MdtCalibOutputDbSvc +/// Service which holds the results of the calibration algorithms. +/// In its initial version the class can only write out the calibration to text +/// files as it is currently performed in the calibration algorithms. The +/// location of the output file is given in the job options file. +/// +/// The name of the service is MdtCalibOutputDbSvc. MdtCalibOutputDbSvc is a +/// singleton. +/// +/// \author Oliver.Kortner@CERN.CH +/// +/// \date 07.10.2006 + +//:::::::::::::::::: +//:: HEADER FILES :: +//:::::::::::::::::: + +// STL // +#include <vector> +#include <string> + +// Gaudi // +#include "GaudiKernel/Service.h" +#include "GaudiKernel/IInterface.h" +#include "GaudiKernel/ToolHandle.h" + +#include "MuonCalibStandAloneBase/NtupleStationId.h" +#include "MuonCalibStandAloneBase/CalibrationIOTool.h" +class RegionSelectionSvc; +class MdtCalibInputSvc; +class MdtIdHelper; +namespace MuonGM { +class MuonDetectorManager; +} + +// MuonCalib // +//#include "MuonCalibIdentifier/MdtRegion.h" +//#include "MuonCalibIdentifier/MdtHashTable.h" +namespace MuonCalib { +class IMdtCalibrationOutput; +class IRtRelation; +class IRtResolution; +class MdtTubeFitContainer; +} +//class CalibDBCoral; +#include "MdtCalibUtils/RtDataFromFile.h" + + +// interface to enable retrieving of a pointer to the singleton // +const InterfaceID IID_IMdtCalibOutputDbSvc("MdtCalibOutputDbSvc", 1, 0); + + +class MdtCalibOutputDbSvc : public Service { + +public: +// Constructor // + MdtCalibOutputDbSvc(const std::string & name, + ISvcLocator *svc_locator); + ///< Constructor + + virtual ~MdtCalibOutputDbSvc(void) {}; + ///< Virtual destructor + +// Methods // +/// Methods required as defined in the base class "Service" + /** interface */ + static const InterfaceID& interfaceID() { return IID_IMdtCalibOutputDbSvc; } + virtual StatusCode queryInterface(const InterfaceID& riid, + void** ppvUnknown); + ///< method required by the base class + ///< which is need to obtain a pointer + ///< to the service in the standard way + virtual StatusCode initialize(void); + ///< initialize method as required by + ///< the base class + + virtual StatusCode finalize(void); + ///< finalize method as required by + ///< the base class; + ///< the finalize method calls the + ///< method "save_calibration_results" + void AddRunNumber(int run_number); + ///< add a run number to the iov + ///< interval. The interval will begin + ///< with the smallest run number, and + ///< end with the largest run number + bool memorize(const MuonCalib::IMdtCalibrationOutput * result); + ///< memorize the result of a particular + ///< calibration (given in result) for + ///< the calibration region "regionKey"; + ///< a previous calibration result + ///< of the same type (e.g. r-t relation + ///< calibration) will be overwritten + ///< internally if overwrite is true; + ///< the calibration which is memorized + ///< will only be saved for ever after + ///< a call to the method + ///< "saved_calibration_results"; + ///< method return true in case of + ///< success, false otherwise + bool memorize(const MuonCalib::IMdtCalibrationOutput * result, + const MuonCalib::IRtResolution * resol); + ///< memorize the result of a particular + ///< calibration (given in result) for + ///< the calibration region "regionKey"; + ///< a previous calibration result + ///< of the same type (e.g. r-t relation + ///< calibration) will be overwritten + ///< internally if overwrite is true; + ///< the calibration which is memorized + ///< will only be saved for ever after + ///< a call to the method + ///< "saved_calibration_results"; + ///< method return true in case of + ///< success, false otherwise; + ///< the user can pass a pointer to a + ///< resolution function which should + ///< be associated with the calibration + ///< output + void reset(void); + ///< reset, clear memory of results + +private: + +// calibration outputs // + const MuonCalib::IMdtCalibrationOutput *m_results; + // vector holding the + // results of the + // calibration + // algorithms +// associated resolution functions // + const MuonCalib::IRtResolution * m_resolution; + // vector of pointers to + // associated resolution + // functions + +// postprocess calibration data - job option +//NOTE: Detector geometry has to be loaded to do this + bool m_postprocess_calibration; +//flat default resolution - job option + double m_flat_default_resolution; +//use default resolution even if a resolution was loaded by the input service + bool m_force_default_resolution; +//calibration io tool to be used + ToolHandle<MuonCalib::CalibrationIOTool> m_calib_output_tool; +// iov range in run numbers// + int m_iov_start, m_iov_end; +//access to geomodel + const MdtIdHelper* m_mdtIdHelper; + const MuonGM::MuonDetectorManager* m_detMgr; +//region selection service + RegionSelectionSvc *p_reg_sel_svc; +//calibration input service + MdtCalibInputSvc *p_input_service; + std::vector<MuonCalib::NtupleStationId> region_ids; +// private methods // + StatusCode saveCalibrationResults(void); + ///< write out all memorized results + ///< to text files (location specified + ///< in the job options) which can be + ///< uploaded to the calibration + ///< database; + ///< method returns true in case of + ///< success, false otherwise + +//postprocess t0 + MuonCalib::MdtTubeFitContainer * postprocess_t0s(MuonCalib::MdtTubeFitContainer * new_t0, const MuonCalib::NtupleStationId &id); + inline void create_default_resolution(const MuonCalib::IRtRelation *rt); +}; + +#endif diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/cmt/requirements b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/cmt/requirements new file mode 100644 index 0000000000000000000000000000000000000000..fb8be05ce21f502791f7a429e8f6d587446fe394 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/cmt/requirements @@ -0,0 +1,42 @@ +package MdtCalibIOSvc + +manager Felix Rauscher <rauscher@cern.ch> +author Felix Rauscher <rauscher@cern.ch> + +use AtlasPolicy * + +use GaudiInterface * External + + + +use MuonCalibStandAloneBase MuonCalibStandAloneBase-* MuonSpectrometer/MuonCalib/MuonCalibStandAlone +branches MdtCalibIOSvc src share + +use MuonCalibStandAloneBase MuonCalibStandAloneBase-* MuonSpectrometer/MuonCalib/MuonCalibStandAlone +branches MdtCalibIOSvc src share + +use MdtCalibUtils MdtCalibUtils-* MuonSpectrometer/MuonCalib/MdtCalib + + + +apply_pattern declare_joboptions files="*.py" +apply_pattern dual_use_library files= *.cxx + +private + +use Identifier Identifier-* DetectorDescription +use MdtCalibData MdtCalibData-* MuonSpectrometer/MuonCalib/MdtCalib +use MdtCalibInterfaces MdtCalibInterfaces-* MuonSpectrometer/MuonCalib/MdtCalib +use MdtCalibT0 MdtCalibT0-* MuonSpectrometer/MuonCalib/MdtCalib +use MdtCalibRt MdtCalibRt-* MuonSpectrometer/MuonCalib/MdtCalib +use MuonCalibIdentifier MuonCalibIdentifier-* MuonSpectrometer/MuonCalib +use MuonCalibMath MuonCalibMath-* MuonSpectrometer/MuonCalib/MuonCalibUtils +use MuonIdHelpers MuonIdHelpers-* MuonSpectrometer +use StoreGate StoreGate-* Control + +use MuonReadoutGeometry MuonReadoutGeometry-* MuonSpectrometer/MuonDetDescr + + +#private +#macro cppdebugflags '$(cppdebugflags_s)' +#macro_remove componentshr_linkopts "-Wl,-s" diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/share/MdtCalibIOSvc.py b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/share/MdtCalibIOSvc.py new file mode 100644 index 0000000000000000000000000000000000000000..b948ef588b24724a0c06639c6429f2875bf084ca --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/share/MdtCalibIOSvc.py @@ -0,0 +1,14 @@ +from AthenaCommon.AppMgr import ServiceMgr +from MdtCalibIOSvc.MdtCalibIOSvcConf import MdtCalibInputSvc +from MdtCalibIOSvc.MdtCalibIOSvcConf import MdtCalibOutputDbSvc +from MdtCalibIOSvc.MdtCalibIOSvcConf import MuonCalib__CalibrationFileIOTool +from MdtCalibIOSvc.MdtCalibIOSvcConf import MuonCalib__CalibrationOracleFileIOTool +MdtCalibInputSvc = MdtCalibInputSvc() +MdtCalibOutputDbSvc = MdtCalibOutputDbSvc() +CalibrationFileIOTool = MuonCalib__CalibrationFileIOTool() +CalibrationOracleFileIOTool = MuonCalib__CalibrationOracleFileIOTool() +from AthenaCommon.AppMgr import ToolSvc +ToolSvc += CalibrationFileIOTool +ToolSvc += CalibrationOracleFileIOTool +ServiceMgr += MdtCalibInputSvc +ServiceMgr += MdtCalibOutputDbSvc diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationFileIOTool.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationFileIOTool.cxx new file mode 100644 index 0000000000000000000000000000000000000000..3f025974dd29c027bf6a1d83a58bfa67b2797209 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationFileIOTool.cxx @@ -0,0 +1,337 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +//this +#include "MdtCalibIOSvc/CalibrationFileIOTool.h" + +//gaudi +#include "GaudiKernel/MsgStream.h" + +//MdtCalibUtils +#include "MdtCalibUtils/TubeDataFromFile.h" +#include "MdtCalibUtils/RtDataFromFile.h" + +//MdtCalibData +#include "MdtCalibData/MdtTubeFitContainer.h" +#include "MdtCalibData/IRtResolution.h" +#include "MdtCalibData/RtChebyshev.h" +#include "MdtCalibData/RtRelationLookUp.h" +#include "MdtCalibData/RtFromPoints.h" +#include "MdtCalibData/RtResolutionFromPoints.h" + +//MdtCalibRt +#include "MdtCalibRt/RtCalibrationOutput.h" + +//MuonCalibStandAloneBase +#include "MuonCalibStandAloneBase/NtupleStationId.h" +#include "MuonCalibStandAloneBase/MdtStationT0Container.h" + +//MuonCalibMath +#include "MuonCalibMath/SamplePoint.h" + + +//c - c++ +#include "string" +#include "iostream" +#include "fstream" +#include <sys/types.h> +#include <dirent.h> + +namespace MuonCalib { + +CalibrationFileIOTool :: CalibrationFileIOTool(const std::string& t, const std::string& n, const IInterface* p) : AlgTool(t, n, p), m_calib_dir("calibration"), m_use_fixed_id(false), m_rt_lookup(true) + { + declareInterface< CalibrationIOTool >(this) ; + declareProperty("outputLocation", m_calib_dir); + declareProperty("UseFixedId", m_use_fixed_id); + declareProperty("RtRelationLookup", m_rt_lookup); + } + + + +StatusCode CalibrationFileIOTool :: WriteT0(MdtTubeFitContainer* t0_output, const NtupleStationId & station_id, int /*iov_start*/, int /*iov_end*/) + { + MsgStream log(msgSvc(), name()); +//create directory + system(("mkdir -p "+m_calib_dir+"/t0s").c_str()); +//open file + std::string t0_file_name(m_calib_dir + "/t0s/T0" + station_id.regionId() + ".dat"); + std::ofstream output_file(t0_file_name.c_str()); + if(output_file.fail()) + { + log << MSG::ERROR << "Could not open t0 file " + << t0_file_name << "!" << endreq; + return StatusCode::FAILURE; + } +//write + TubeDataFromFile t0_file; + t0_file.setNRegions(1); + t0_file.addTubes(0, t0_output); + t0_file.write(output_file); + return StatusCode::SUCCESS; + } + + +StatusCode CalibrationFileIOTool :: WriteRt(const RtCalibrationOutput *rt_relation, const IRtResolution * resolution, const NtupleStationId & station_id, int /*iov_start*/, int /*iov_end*/, bool /*real_rt*/, bool /*real_resolution*/) + { + MsgStream log(msgSvc(), name()); + RtDataFromFile rt_data; + const IRtRelation *new_rt = rt_relation->rt(); + RtDataFromFile::RtRelation *rt(new RtDataFromFile::RtRelation()); + if(m_use_fixed_id) + { + rt->setRegionId(station_id.FixedId()); + rt_data.setVersion(2, 0); + } + else + { + if(station_id.RegionHash()!=0) + { + rt->setRegionId(station_id.RegionHash()); + } + } + int rt_region_id(0); + rt_data.setNRts(1); + rt_data.addRt(rt_region_id, rt, rt_relation->fullInfo()); + if(!fill_rt(rt, new_rt, resolution)) + return StatusCode::FAILURE; +// write out the r-t relationship // + system(("mkdir -p "+m_calib_dir+"/rts").c_str()); + std::string rt_file_name(m_calib_dir + "/rts/Rt_" + station_id.regionId() + ".dat"); + std::ofstream output_file(rt_file_name.c_str()); + rt_data.write(output_file, rt_region_id); + return StatusCode::SUCCESS; + } + + + +StatusCode CalibrationFileIOTool :: LoadT0(std::map<NtupleStationId, MdtStationT0Container *> &t0s, int /*iov_id*/) + { + MsgStream log(msgSvc(), name()); + log<< MSG::INFO << "Reading calibration from '" << (m_calib_dir + "/t0s") << "'" << endreq; + t0s.clear(); + DIR *directory(opendir((m_calib_dir + "/t0s").c_str())); + if(directory==NULL) return StatusCode::SUCCESS; + struct dirent *dent; +//loop on all files in directory + while((dent=readdir(directory))!=NULL) + { + std::string nm(dent->d_name); + std::string station_str; + int eta, phi, ml; + //interpret filename and check if it is an t0 file + if(!interpret_chamber_name(nm, "T0", station_str, eta, phi, ml)) continue; + //store + MuonCalib::NtupleStationId id; + if(!id.Initialize(station_str, eta, phi, ml)) continue; + log << MSG::INFO << "Reading t0 for " << station_str << "_" << eta << "_" << phi << endreq; + MuonCalib::MdtStationT0Container *t0_container = new MuonCalib::MdtStationT0Container((m_calib_dir + "/t0s/" + nm).c_str()); + if(t0_container->t0_loaded()) + { + t0s[id] = t0_container; + log<< MSG::INFO << "Read t0 for "<<id.regionId()<<endreq; + } + else + { + log<<MSG::ERROR<< "t0 file for "<<id.regionId()<<" is corrupted!"<<endreq; + delete t0_container; + } + } + closedir(directory); + return StatusCode::SUCCESS; + } + + + +StatusCode CalibrationFileIOTool :: LoadRt(std::map<NtupleStationId, IRtRelation *> & rts, std::map<NtupleStationId, IRtResolution *> & res, int /*iov_id*/) + { + rts.clear(); + res.clear(); + MsgStream log(msgSvc(), name()); + DIR *directory(opendir((m_calib_dir + "/rts").c_str())); + if(directory==NULL) return StatusCode::SUCCESS; + struct dirent *dent; +//loop on all files in directory + while((dent=readdir(directory))!=NULL) + { + std::string nm(dent->d_name); + std::string station_str; + int eta, phi, ml; + //interpret filename and check if it is an t0 file + if(!interpret_chamber_name(nm, "Rt_", station_str, eta, phi, ml)) continue; + MuonCalib::NtupleStationId id; + if(!id.Initialize(station_str, eta, phi, ml)) continue; + log << MSG::INFO << "Reading rt for " << station_str << "_" << eta << "_" << phi << "_" << ml << endreq; + //read rt ralation + read_rt_relation((m_calib_dir + "/rts/" + nm), rts, res, id); + log << MSG::INFO << "Read rt for " << id.regionId() << endreq; + } + closedir(directory); + return StatusCode::SUCCESS; + } + + +inline bool CalibrationFileIOTool::fill_rt(RtDataFromFile::RtRelation *rt, const IRtRelation *new_rt, const MuonCalib::IRtResolution * resolut) + { + MsgStream log(msgSvc(), name()); +/////////////// +// VARIABLES // +/////////////// + + const CalibFunc::ParVec & rt_param = new_rt->parameters(); + // parameters of the r-t relationship to be + // copied + const RtChebyshev *rt_Chebyshev(dynamic_cast<const RtChebyshev *>(new_rt)); + const RtRelationLookUp *rt_lookup(dynamic_cast<const RtRelationLookUp *>(new_rt)); + + +//////////////////////// +// FILL THE r-t CLASS // +//////////////////////// + +// case 1: r-t relationship is stored in the class RtChebyshev // + if (rt_Chebyshev) { + unsigned int nb_points(100); + double t_length(rt_Chebyshev->tUpper()-rt_Chebyshev->tLower()); + double bin_size(t_length/static_cast<double>(nb_points-1)); + for (unsigned int k=0; k<nb_points; k++) { + double time(rt_Chebyshev->tLower()+k*bin_size); + double radius(rt_Chebyshev->radius(time)); + double resol(0.0); + resol = resolut->resolution(time); + if(std::isnan(time) || std::isnan(radius) || std::isnan(resol)) + { + log<<MSG::FATAL<<"Filling nan into rt relation!"<<endreq; +// return false; + } + rt->addEntry(time, radius, resol); + } + } + +// case 2: r-t relationship is stored in the class RtRelationLookUp // + if (rt_lookup) + { + double t_min(rt_param[0]); + double bin_size = rt_param[1]; + unsigned int nb_points(rt_lookup->nPar()-2); + for (unsigned int k=0; k<nb_points; k++) + { + double radius(rt_param[k+2]); + double resol(0.0); + + resol = resolut->resolution(t_min+bin_size*k); + if(std::isnan(radius) || std::isnan(resol)) + { + log<<MSG::FATAL<<"Filling nan into rt relation!"<<endreq; +// return false; + } + rt->addEntry(t_min+bin_size*k, rt_param[k+2], resol); + } + } + +//////////////////////////////////////////// +// DECLARE THE NEW r-t RELATIONSHIP VALID // +//////////////////////////////////////////// + if (new_rt->HasTmaxDiff()) + rt->addEntry(new_rt->GetTmaxDiff(), -9999, 0); + rt->isValid(1); + + return true; + +} + +inline bool CalibrationFileIOTool :: interpret_chamber_name(const std::string &nm, const char *prefix, std::string & station, int &eta, int & phi, int & ml) const + { +//check if name begins with the prefix + std::string prefix_st(prefix); + if(nm.find(prefix_st) !=0) return false; +//check if filename ends in .dat + if(static_cast<int>(nm.find(".dat")) <0 || nm.find(".dat")!=nm.size()-4) return false; +//cut prefix and suffix from filename + std::string cutout(nm, prefix_st.size(), nm.size()-4-prefix_st.size()); +//extrect station name + int uscore_pos(cutout.find('_')); + if(uscore_pos<=0) return false; + station=std::string(cutout, 0, uscore_pos); +//get eta ,phi and multilayer + std::string cutout2(cutout, uscore_pos+1); + int count_items(sscanf(cutout2.c_str(), "%d_%d_%d", &phi, &eta, &ml)); + if(count_items<2) return false; + if(count_items!=3) ml=0; + return true; + } + + +inline void CalibrationFileIOTool :: read_rt_relation(const std::string & fname, std::map<NtupleStationId, IRtRelation *> & rts, std::map<NtupleStationId, IRtResolution *> &res_map,const MuonCalib::NtupleStationId & id) + { + MsgStream log(msgSvc(), name()); +//open file + std::ifstream infile(fname.c_str()); + if(infile.fail()) + { + log << MSG::ERROR << "Cannot open file '" << fname << "' for reading!" << endreq; + return; + } +//sample points + std::vector<double> r, t, res; + std::string sdummy; + double dummy; +//ignore the first line + getline(infile, sdummy); +//read table + float multilayer_diff(9e9); + while (!(infile.eof() || infile.fail())) + { + infile >> dummy; + if (infile.eof() || infile.fail() || dummy > 15.0) break; + if(dummy<-9998) + { + infile >> multilayer_diff; + infile >> dummy; + } + else + { + r.push_back(dummy); + infile >> dummy; + t.push_back(dummy); + infile >> dummy; res.push_back(dummy); + } + } + if(r.size()<3) + { + log<<MSG::WARNING<< "Not enough good rt points for "<<id.regionId()<<"!"<<endreq; + return; + } + std::vector<MuonCalib::SamplePoint> point(r.size()); + for(unsigned int k=0; k<point.size(); k++) + { + point[k].set_x1(t[k]); + point[k].set_x2(r[k]); + point[k].set_error(1.0); + } +//convert rt relation + MuonCalib::RtFromPoints rt_from_points; + if(m_rt_lookup) + { + rts[id] = new MuonCalib::RtRelationLookUp(rt_from_points.getRtRelationLookUp(point)); + } + else + { + rts[id] = new MuonCalib::RtChebyshev((rt_from_points.getRtChebyshev(point, 10))); + } + if(multilayer_diff<8e8) + rts[id]->SetTmaxDiff(multilayer_diff); +// create resolution function + MuonCalib::RtResolutionFromPoints res_from_points; + for (unsigned int k=0; k<point.size(); k++) + { + point[k].set_x1(t[k]); + point[k].set_x2(res[k]); + point[k].set_error(1.0); + } + res_map[id] = new MuonCalib::RtResolutionChebyshev(res_from_points.getRtResolutionChebyshev(point, 8)); + } + + +} diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationOracleFileIOTool.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationOracleFileIOTool.cxx new file mode 100644 index 0000000000000000000000000000000000000000..896f94d9c5fb13357c631042ad23c5645aae7556 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/CalibrationOracleFileIOTool.cxx @@ -0,0 +1,199 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +//this +#include "MdtCalibIOSvc/CalibrationOracleFileIOTool.h" + +//gaudi +#include "GaudiKernel/MsgStream.h" + +//MdtCalibUtils +#include "MdtCalibUtils/TubeDataFromFile.h" +#include "MdtCalibUtils/RtDataFromFile.h" + +//MdtCalibData +#include "MdtCalibData/MdtTubeFitContainer.h" +#include "MdtCalibData/IRtResolution.h" +#include "MdtCalibData/RtChebyshev.h" +#include "MdtCalibData/RtRelationLookUp.h" +#include "MdtCalibData/RtFromPoints.h" +#include "MdtCalibData/RtResolutionFromPoints.h" + +//MdtCalibRt +#include "MdtCalibRt/RtCalibrationOutput.h" + +//MuonCalibStandAloneBase +#include "MuonCalibStandAloneBase/NtupleStationId.h" +#include "MuonCalibStandAloneBase/MdtStationT0Container.h" + +//MuonCalibMath +#include "MuonCalibMath/SamplePoint.h" + + +//c - c++ +#include "string" +#include "iostream" +#include "fstream" +#include <sys/types.h> +#include <dirent.h> + +namespace MuonCalib { + +CalibrationOracleFileIOTool :: CalibrationOracleFileIOTool(const std::string& t, const std::string& n, const IInterface* p) : AlgTool(t, n, p), m_calib_dir("oracle_calib") + { + declareInterface< CalibrationIOTool >(this) ; + declareProperty("outputLocation", m_calib_dir); + } + + + +StatusCode CalibrationOracleFileIOTool :: WriteT0(MdtTubeFitContainer* t0_output, const NtupleStationId & station_id, int iov_start, int iov_end) + { + MsgStream log(msgSvc(), name()); +//create directory + system(("mkdir -p "+m_calib_dir+"/t0s").c_str()); +//open file + std::string t0_dbfile_name(m_calib_dir + "/t0s/DB_t0" + station_id.regionId() + ".dat"); + FILE * db_file; + db_file=fopen(t0_dbfile_name.c_str(),"w"); + if(db_file==NULL) + { + log<< MSG::FATAL <<"Cannot open output file!" << endreq; + return StatusCode::FAILURE; + } + TubeDataFromFile t0_file; + t0_file.setNRegions(1); + t0_file.addTubes(0, t0_output); + t0_file.write_forDB(db_file, -9999, iov_start, iov_end); + fclose(db_file); + return StatusCode::SUCCESS; + } + + +StatusCode CalibrationOracleFileIOTool :: WriteRt(const RtCalibrationOutput *rt_relation, const IRtResolution * resolution, const NtupleStationId & station_id, int iov_start, int iov_end, bool /*real_rt*/, bool /*real_resolution*/) + { + MsgStream log(msgSvc(), name()); + RtDataFromFile rt_data; + const IRtRelation *new_rt = rt_relation->rt(); + RtDataFromFile::RtRelation *rt(new RtDataFromFile::RtRelation()); + int rt_region_id(0); + rt_data.setNRts(1); + rt_data.addRt(rt_region_id, rt, rt_relation->fullInfo()); + if(!fill_rt(rt, new_rt, resolution)) + return StatusCode::FAILURE; +// write out the r-t relationship // + std::string location="location"; + int mdt_rt_id=99; + std::string mdt_rt_map_t_id="mdt_rt_map_t_id"; + int mdt_rt_map_r_id=99; + int mdt_rt_map_s_id=99; + system(("mkdir -p "+m_calib_dir+"/rts").c_str()); + std::string rt_dbfile_name(m_calib_dir + "/rts/DBrt_" + station_id.regionId() + ".dat"); + std::string rtt_dbfile_name(m_calib_dir + "/rts/DBrt_" + station_id.regionId() + "_t.dat"); + std::string rtr_dbfile_name(m_calib_dir + "/rts/DBrt_" + station_id.regionId() + "_r.dat"); + std::string rts_dbfile_name(m_calib_dir + "/rts/DBrt_" + station_id.regionId() + "_s.dat"); + FILE * db_rt_file; + db_rt_file=fopen(rt_dbfile_name.c_str(),"w"); + FILE * db_rtt_file; + db_rtt_file=fopen(rtt_dbfile_name.c_str(),"w"); + FILE * db_rtr_file; + db_rtr_file=fopen(rtr_dbfile_name.c_str(),"w"); + FILE * db_rts_file; + db_rts_file=fopen(rts_dbfile_name.c_str(),"w"); + if(db_rt_file==NULL || db_rtt_file==NULL || db_rtr_file==NULL || db_rts_file==NULL) + { + log << MSG::FATAL << "Cannot open output files" <<endreq; + if (db_rt_file) delete db_rt_file; + if (db_rtt_file) delete db_rtt_file; + if (db_rtr_file) delete db_rtr_file; + if (db_rts_file) delete db_rts_file; + return StatusCode::FAILURE; + } + log << MSG::INFO << "Writing out r-t relationships in the files for calibration db." << endreq; + // for the time being the variables: + // mdt_rt_id,mdt_rt_map_t_id,mdt_rt_map_r_id, + // mdt_rt_map_s_id and location + // are dummy. + fprintf(db_rt_file," %d,%d,%d,%d,%d,%s,", mdt_rt_id, station_id.FixedId(), -9999, iov_start, iov_end, location.c_str()); + fprintf(db_rtt_file," %s,%d,%d,", mdt_rt_map_t_id.c_str(), mdt_rt_id, mdt_rt_id); + fprintf(db_rtr_file," %d,%d,%s,", mdt_rt_map_r_id, mdt_rt_id, mdt_rt_map_t_id.c_str()); + fprintf(db_rts_file," %d,%d,",mdt_rt_map_s_id,mdt_rt_id); + rt_data.write_forDB(db_rt_file, db_rtt_file, db_rtr_file, db_rts_file, rt_region_id); + log << MSG::INFO << "r-t relationships wrote in the files for calibration db." << endreq; + fclose(db_rt_file); + fclose(db_rtt_file); + fclose(db_rtr_file); + fclose(db_rts_file); + return StatusCode::SUCCESS; + } + + +inline bool CalibrationOracleFileIOTool::fill_rt(RtDataFromFile::RtRelation *rt, const IRtRelation *new_rt, const MuonCalib::IRtResolution * resolut) + { + MsgStream log(msgSvc(), name()); +/////////////// +// VARIABLES // +/////////////// + + const CalibFunc::ParVec & rt_param = new_rt->parameters(); + // parameters of the r-t relationship to be + // copied + const RtChebyshev *rt_Chebyshev(dynamic_cast<const RtChebyshev *>(new_rt)); + const RtRelationLookUp *rt_lookup(dynamic_cast<const RtRelationLookUp *>(new_rt)); + + +//////////////////////// +// FILL THE r-t CLASS // +//////////////////////// + +// case 1: r-t relationship is stored in the class RtChebyshev // + if (rt_Chebyshev) { + unsigned int nb_points(100); + double t_length(rt_Chebyshev->tUpper()-rt_Chebyshev->tLower()); + double bin_size(t_length/static_cast<double>(nb_points-1)); + for (unsigned int k=0; k<nb_points; k++) { + double time(rt_Chebyshev->tLower()+k*bin_size); + double radius(rt_Chebyshev->radius(time)); + double resol(0.0); + resol = resolut->resolution(time); + if(std::isnan(time) || std::isnan(radius) || std::isnan(resol)) + { + log<<MSG::FATAL<<"Filling nan into rt relation!"<<endreq; +// return false; + } + rt->addEntry(time, radius, resol); + } + } + +// case 2: r-t relationship is stored in the class RtRelationLookUp // + if (rt_lookup) + { + double t_min(rt_param[0]); + double bin_size = rt_param[1]; + unsigned int nb_points(rt_lookup->nPar()-2); + for (unsigned int k=0; k<nb_points; k++) + { + double radius(rt_param[k+2]); + double resol(0.0); + + resol = resolut->resolution(t_min+bin_size*k); + if(std::isnan(radius) || std::isnan(resol)) + { + log<<MSG::FATAL<<"Filling nan into rt relation!"<<endreq; +// return false; + } + rt->addEntry(t_min+bin_size*k, rt_param[k+2], resol); + } + } + +//////////////////////////////////////////// +// DECLARE THE NEW r-t RELATIONSHIP VALID // +//////////////////////////////////////////// + + rt->isValid(1); + + return true; + +} +} diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibInputSvc.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibInputSvc.cxx new file mode 100644 index 0000000000000000000000000000000000000000..3dbced444381fe74176d1eb15d25973ae5471af6 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibInputSvc.cxx @@ -0,0 +1,285 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + + +//c - c++ +#include "fstream" +#include <sys/types.h> +#include <dirent.h> +#include "iostream" +#include "list" + +// Gaudi // +#include "GaudiKernel/MsgStream.h" + +//MdtCalibData +#include "MdtCalibData/IRtRelation.h" +#include "MdtCalibData/IRtResolution.h" +#include "MdtCalibData/BFieldCorFunc.h" +#include "MdtCalibData/RtFromPoints.h" + +//MuonCalibMath +#include "MuonCalibMath/SamplePoint.h" +#include "MdtCalibData/RtResolutionFromPoints.h" + +//this +#include "MdtCalibIOSvc/MdtCalibInputSvc.h" + +//MuonCalibStandAloneBase +#include "MuonCalibStandAloneBase/MdtStationT0Container.h" +#include "MuonCalibStandAloneBase/NtupleStationId.h" +#include "MuonCalibStandAloneBase/RegionSelectionSvc.h" + +#include "MuonCalibIdentifier/MuonFixedId.h" + +MdtCalibInputSvc ::MdtCalibInputSvc(const std::string & name, ISvcLocator *svc_locator) : Service(name, svc_locator), m_calib_input_tool("MuonCalib::CalibrationDummyIOTool") + { + declareProperty("CalibrationInputTool", m_calib_input_tool); + p_sel_region_res=NULL; + p_sel_region_rt=NULL; + p_sel_region_b=NULL; + p_reg_sel_svc=NULL; + } + + +MdtCalibInputSvc ::~MdtCalibInputSvc() + { + std::map<MuonCalib::NtupleStationId, MuonCalib::MdtStationT0Container *> :: iterator it1; + for(it1=m_t0.begin(); it1!=m_t0.end(); it1++) + delete it1->second; + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtRelation *> :: iterator it2; + for(it2=m_rt_relation.begin(); it2!=m_rt_relation.end(); it2++) + delete it2->second; + std::map<MuonCalib::NtupleStationId, MuonCalib::BFieldCorFunc *> :: iterator it3; + for(it3=m_B_corr.begin(); it3!=m_B_corr.end(); it3++) + delete it3->second; + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtResolution *> :: iterator it4; + for(it4=m_spat_res.begin(); it4!=m_spat_res.end(); it4++) + delete it4->second; + } + + +StatusCode MdtCalibInputSvc :: initialize() + { + MsgStream log(messageService(), name()); + if(!Service::initialize().isSuccess()) + { + log << MSG::FATAL <<"Failed to initialize service"<<endreq; + return StatusCode::FAILURE; + } + log<< MSG::INFO << "initialize()" << endreq; +//get region selection service + StatusCode sc=service("RegionSelectionSvc", p_reg_sel_svc); + if(!sc.isSuccess() || p_reg_sel_svc==NULL) + { + log << MSG::ERROR <<"Cannot retrieve RegionSelectionSvc!" <<endreq; + return sc; + } + sc=m_calib_input_tool.retrieve(); + if(!sc.isSuccess()) + { + log << MSG::ERROR <<"Cannot retrieve inpput tool!" <<endreq; + return sc; + } + sc=read_calib_input(); + return sc; + } + +StatusCode MdtCalibInputSvc :: finalize(void) + { + return StatusCode::SUCCESS; + } + +StatusCode MdtCalibInputSvc::queryInterface(const InterfaceID& riid, + void** ppvUnknown) { + std::cout<<"StatusCode MdtCalibInputSvc::queryInterface"<<std::endl; + + if (IID_IMdtCalibInputSvc.versionMatch(riid)) { + *ppvUnknown = (MdtCalibInputSvc *)this; + } else { + return Service::queryInterface(riid, ppvUnknown); + } + + return StatusCode::SUCCESS; +} + + +const MuonCalib::MdtStationT0Container * MdtCalibInputSvc :: GetT0(const MuonCalib::NtupleStationId & id) const + { + std::map<MuonCalib::NtupleStationId, MuonCalib::MdtStationT0Container *> :: const_iterator it; + if((it=m_t0.find(id)) == m_t0.end()) + { + MuonCalib::NtupleStationId chamber_id=id; + chamber_id.SetMultilayer(0); + if((it=m_t0.find(chamber_id)) == m_t0.end()) + { + if(m_t0_warned.find(chamber_id) == m_t0_warned.end()) + { + MsgStream log(messageService(), name()); + log<< MSG::WARNING << "T0 not loaded for station " << id.regionId() <<endreq; + m_t0_warned.insert(chamber_id); + } + return NULL; + } + } + return it->second; + } + + +const MuonCalib::IRtRelation * MdtCalibInputSvc :: GetRtRelation(const MuonCalib::NtupleStationId & id) const + { + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtRelation *> :: const_iterator it; + if((it=m_rt_relation.find(id)) == m_rt_relation.end()) + { + MuonCalib::NtupleStationId chamber_id=id; + chamber_id.SetMultilayer(0); + if((it=m_rt_relation.find(chamber_id)) == m_rt_relation.end()) + { + if(m_rt_warned.find(chamber_id) == m_rt_warned.end()) + { + MsgStream log(messageService(), name()); + log<< MSG::WARNING << "Rt relation not loaded for station" << chamber_id.regionId() <<endreq; + m_rt_warned.insert(chamber_id); + } + return NULL; + } + } + return it->second; + } + + +const MuonCalib::BFieldCorFunc * MdtCalibInputSvc :: GetBCorr(const MuonCalib::NtupleStationId & id) const + { + const MuonCalib::BFieldCorFunc * fun; + if((fun=findbfieldfun(id))!=NULL) + return fun; + MuonCalib::NtupleStationId chamber_id(id); + chamber_id.SetMultilayer(0); + return findbfieldfun(chamber_id); + } + +const MuonCalib::IRtResolution * MdtCalibInputSvc :: GetResolution(const MuonCalib::NtupleStationId & id) const + { + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtResolution *>:: const_iterator it; + if((it=m_spat_res.find(id)) == m_spat_res.end()) + { + MuonCalib::NtupleStationId chamber_id=id; + chamber_id.SetMultilayer(0); + if((it=m_spat_res.find(chamber_id)) == m_spat_res.end()) + { + if(m_rt_warned.find(chamber_id) == m_rt_warned.end()) + { + MsgStream log(messageService(), name()); + log<< MSG::FATAL << "Rt relation not loaded for station" <<endreq; + m_rt_warned.insert(chamber_id); + } + return NULL; + } + } + return it->second; + } + + + +inline StatusCode MdtCalibInputSvc :: read_calib_input() + { + MsgStream log(messageService(), name()); + StatusCode sc=m_calib_input_tool->LoadT0(m_t0, -1); + if(!sc.isSuccess()) + { + log << MSG::ERROR << "Cannot read t0s" <<endreq; + return sc; + } + sc=m_calib_input_tool->LoadRt(m_rt_relation, m_spat_res, -1); + if(!sc.isSuccess()) + { + log << MSG::ERROR << "Cannot read rts" <<endreq; + return sc; + } + create_mean_rts(); + return sc; + } + + +inline bool MdtCalibInputSvc::create_b_field_correction(const MuonCalib::NtupleStationId &id) const + { + std::map<MuonCalib::NtupleStationId, MuonCalib::IRtRelation *> :: const_iterator it(m_rt_relation.find(id)); + if (it==m_rt_relation.end()) return false; + MsgStream log(messageService(), name()); + log<< MSG::INFO << "Initiailizing B-Field correction for " << id.regionId() <<endreq; + const MuonCalib::IRtRelation *rt_rel(it->second); +// create magnetic field correction + std::vector<double> corr_params(2); + corr_params[0] = 3080.0; // high voltage + corr_params[1] = 0.11; // epsilon parameter + m_B_corr[id] = new MuonCalib::BFieldCorFunc(corr_params, rt_rel); + return true; + } + +inline void MdtCalibInputSvc::create_mean_rts() + { + std::cout<<"MdtCalibInputSvc::create_mean_rts()"<<std::endl; + MsgStream log(messageService(), name()); + std::list<const MuonCalib::IRtRelation *> matching_relations; + std::list<MuonCalib::NtupleStationId> matching_ids; +//loop over all rts + for(std::map<MuonCalib::NtupleStationId, MuonCalib::IRtRelation *>::const_iterator it = m_rt_relation.begin(); it!=m_rt_relation.end(); it++) + { + MuonCalib::MuonFixedId id; + //if the rt relation is stored per chamber, check both multilayers + if(it->first.GetMl()==0) + { + id.set_mdt(); + for(int i=1; i<=2; i++) + { + id.setStationName(it->first.GetStation()); + id.setStationEta(it->first.GetEta()); + id.setStationPhi(it->first.GetPhi()); + id.setMdtMultilayer(i); + if(p_reg_sel_svc->isInRegion(id)) + { + matching_relations.push_back( it->second); + matching_ids.push_back(it->first); + break; + } + } + } + else + { + id.setStationName(it->first.GetStation()); + id.setStationEta(it->first.GetEta()); + id.setStationPhi(it->first.GetPhi()); + id.setMdtMultilayer(it->first.GetMl()); + if(p_reg_sel_svc->isInRegion(id)) + { + matching_relations.push_back( it->second); + matching_ids.push_back(it->first); + } + } + } + log << MSG::INFO <<"Found "<<matching_relations.size()<< " rt-relations for calibration region"<<endreq; + if(matching_relations.size()==0) return; +//averageing over rt relations is not yet implemented - take the first found + if(matching_relations.size()>1) + { + log << MSG::WARNING <<"More than one rt relation for this region loaded! Taking first!" << endreq; + } + p_sel_region_rt = *(matching_relations.begin()); + m_mean_station_id = *(matching_ids.begin()); +// p_sel_region_b = GetBCorr(*(matching_ids.begin())); + p_sel_region_b = NULL; + p_sel_region_res = GetResolution(*(matching_ids.begin())); + } + +inline const MuonCalib::BFieldCorFunc * MdtCalibInputSvc :: findbfieldfun(const MuonCalib::NtupleStationId & id) const + { + std::map<MuonCalib::NtupleStationId, MuonCalib::BFieldCorFunc *> :: const_iterator it(m_B_corr.find(id)); + if(it==m_B_corr.end()) + { + if(!create_b_field_correction(id)) + return NULL; + return findbfieldfun(id); + } + return it->second; + } diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibOutputDbSvc.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibOutputDbSvc.cxx new file mode 100644 index 0000000000000000000000000000000000000000..cfb14aeae0526a83cc8e9f828b9a2ca0eab9512c --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/MdtCalibOutputDbSvc.cxx @@ -0,0 +1,500 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ +// 07.10.2006, AUTHOR: OLIVER KORTNER +//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + +//:::::::::::::::::: +//:: HEADER FILES :: +//:::::::::::::::::: + +// standard C++ // +#include <iostream> +#include <fstream> +#include <typeinfo> +#include "cmath" + +// Gaudi // +#include "GaudiKernel/MsgStream.h" +#include "Identifier/Identifier.h" +#include "StoreGate/StoreGateSvc.h" + +//geo model +#include "MuonIdHelpers/MdtIdHelper.h" +#include "MuonReadoutGeometry/MuonDetectorManager.h" + + +// MuonCalib // +#include "MdtCalibIOSvc/MdtCalibOutputDbSvc.h" +#include "MdtCalibIOSvc/MdtCalibInputSvc.h" +#include "MdtCalibT0/T0CalibrationOutput.h" +#include "MdtCalibData/MdtTubeFitContainer.h" +#include "MdtCalibRt/RtCalibrationOutput.h" +#include "MdtCalibData/IRtRelation.h" +//#include "MdtCalibData/RtChebyshev.h" +//#include "MdtCalibData/RtRelationLookUp.h" +#include "MdtCalibData/RtResolutionLookUp.h" +//#include "MuonCalibIdentifier/MuonFixedId.h" + +#include "MdtCalibInterfaces/IMdtCalibrationOutput.h" +//#include "MuonCalibDbOperations/CalibDBCoral.h" +#include "MdtCalibData/IRtResolution.h" +#include "MdtCalibData/RtResolutionFromPoints.h" + +#include "MuonCalibStandAloneBase/RegionSelectionSvc.h" + +//:::::::::::::::::::::::: +//:: NAMESPACE SETTINGS :: +//:::::::::::::::::::::::: + +using namespace std; +using namespace MuonCalib; + +//:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: +//:: IMPLEMENTATION OF METHODS DEFINED IN THE CLASS MdtCalibOutputDbSvc :: +//:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: + +////////////////// +// SET POINTERS // +////////////////// + +//MdtCalibOutputDbSvc *MdtCalibOutputDbSvc::m_MdtCalibOutputDbSvc_pointer=0; + +//***************************************************************************** + + +//:::::::::::::::::::::::: +//:: PUBLIC CONSTRUCTOR :: +//:::::::::::::::::::::::: + +//MdtCalibOutputDbSvc::MdtCalibOutputDbSvc(const std::string & name, +// ISvcLocator *svc_locator) : Service(name, svc_locator), m_results(NULL), m_postprocess_calibration(false), m_write_root_database(false), m_db(NULL), m_iov_start(-1), m_iov_end(-1) { +MdtCalibOutputDbSvc::MdtCalibOutputDbSvc(const std::string & name,ISvcLocator *svc_locator) : Service(name, svc_locator), m_results(NULL), m_postprocess_calibration(false), m_calib_output_tool("MuonCalib::CalibrationFileIOTool"), m_iov_start(-1), m_iov_end(-1) { + + declareProperty("PostprocessCalibration", m_postprocess_calibration); + m_flat_default_resolution=-1; + declareProperty("FlatDefaultResolution", m_flat_default_resolution); + declareProperty("OutputTool", m_calib_output_tool); + m_force_default_resolution=false; + declareProperty("ForceDefaultResolution", m_force_default_resolution); + + //for the sake of coverity + p_reg_sel_svc=NULL; + m_detMgr=NULL; + m_mdtIdHelper=NULL; + m_resolution=NULL; + p_input_service=NULL; + + return; + +} + +//***************************************************************************** + +//::::::::::::::::::::::::::: +//:: METHOD queryInterface :: +//::::::::::::::::::::::::::: + +StatusCode MdtCalibOutputDbSvc::queryInterface(const InterfaceID& riid, + void** ppvUnknown) { + + if (IID_IMdtCalibOutputDbSvc.versionMatch(riid)) { + *ppvUnknown = (MdtCalibOutputDbSvc *)this; + } else { + return Service::queryInterface(riid, ppvUnknown); + } + + return StatusCode::SUCCESS; + +} + +//***************************************************************************** + +//::::::::::::::::::::::: +//:: METHOD initialize :: +//::::::::::::::::::::::: + +StatusCode MdtCalibOutputDbSvc::initialize(void) { + + MsgStream log(messageService(), name()); + if(!Service::initialize().isSuccess()) + { + log << MSG::FATAL <<"Failed to initialize service"<<endreq; + return StatusCode::FAILURE; + } + +//to get service via static function +// m_MdtCalibOutputDbSvc_pointer = this; + log << MSG::INFO <<"initialize MdtCalibOutputDbSvc"<<endreq; + + +//get id helper and detector manager if postprocessing of calibration is selected + if(m_postprocess_calibration) + { + //retrieve detector store + StoreGateSvc* m_detStore; + StatusCode sc = serviceLocator()->service("DetectorStore", m_detStore); + if ( sc.isSuccess() ) + { + log << MSG::DEBUG << "Retrieved DetectorStore" << endreq; + } + else + { + log << MSG::ERROR << "Failed to retrieve DetectorStore" << endreq; + return sc; + } + //retrieve mdt id helper + sc = m_detStore->retrieve(m_mdtIdHelper, "MDTIDHELPER" ); + if (!sc.isSuccess()) + { + log << MSG::ERROR << "Can't retrieve MdtIdHelper" << endreq; + return sc; + } + //retrieve detector manager + sc = m_detStore->retrieve( m_detMgr ); + if (!sc.isSuccess()) + { + log << MSG::ERROR << "Can't retrieve MuonDetectorManager" << endreq; + return sc; + } + } +//get region selection service + StatusCode sc=service("RegionSelectionSvc", p_reg_sel_svc); + if(!sc.isSuccess()) + { + log << MSG::ERROR <<"Cannot retrieve RegionSelectionSvc!" <<endreq; + return sc; + } + region_ids=p_reg_sel_svc->GetStationsInRegions(); + log<< MSG::INFO << "Regions selected: "<<region_ids.size()<<endreq; +//retrieve tool + sc=m_calib_output_tool.retrieve(); + if(sc.isFailure()) + { + log << MSG::FATAL << "Cannot retrieve IO tool!" <<endreq; + return sc; + } + sc = service("MdtCalibInputSvc", p_input_service); + if(sc.isFailure()) + { + log << MSG::FATAL << "Cannot retrieve calibration input service" <<endreq; + } + return StatusCode::SUCCESS; + +} + +//***************************************************************************** + +//::::::::::::::::::::: +//:: METHOD finalize :: +//::::::::::::::::::::: + +StatusCode MdtCalibOutputDbSvc::finalize(void) { + StatusCode sc=saveCalibrationResults(); + MsgStream log(messageService(), name()); + log<< MSG::INFO<<"Results saved!"<<endreq; + return sc; +} + + + +//***************************************************************************** + +//:::::::::::::::::::::::::::::: +//:: METHOD AddRunNumber :: +//:::::::::::::::::::::::::::::: + + +void MdtCalibOutputDbSvc:: AddRunNumber(int run_number) + { + if(run_number < m_iov_start || m_iov_start < 0) + { + m_iov_start = run_number; + } + if(run_number > m_iov_end || m_iov_end < 0) + { + m_iov_end = run_number; + } + } + +//***************************************************************************** + +//:::::::::::::::::::::::::::::: +//:: METHOD memorize(., ., .) :: +//:::::::::::::::::::::::::::::: + +bool MdtCalibOutputDbSvc::memorize( const MuonCalib::IMdtCalibrationOutput * result) { + + +// no overwriting is allowed // + m_results=result; + m_resolution=0; + + return true; + +} + +//***************************************************************************** + +//::::::::::::::::::::::::::::::::: +//:: METHOD memorize(., ., ., .) :: +//::::::::::::::::::::::::::::::::: + +bool MdtCalibOutputDbSvc::memorize(const MuonCalib::IMdtCalibrationOutput * result, const MuonCalib::IRtResolution * resol) { + + m_results=result; + m_resolution=resol; + + return true; + +} + +//***************************************************************************** + +//::::::::::::::::::::::::::::::::::: +//:: METHOD saveCalibrationResults :: +//::::::::::::::::::::::::::::::::::: + +StatusCode MdtCalibOutputDbSvc::saveCalibrationResults(void) { + MsgStream log(messageService(), name()); + + if(m_results==NULL) return true; + + + StatusCode sc; + +/////////////// +// VARIABLES // +/////////////// + + + +///////////////////////////////////////////////////////////////////////////// +// LOOP OVER THE STORED RESULTS AND WRITE THEM OUT ACCORDING TO THEIR TYPE // +///////////////////////////////////////////////////////////////////////////// + +//----------------------------------------------------------------------------- + for (std::vector<MuonCalib::NtupleStationId>::const_iterator it=region_ids.begin(); it!=region_ids.end(); it++) + { +//----------------------------------------------------------------------------- + MuonCalib::NtupleStationId the_id(*it); +// get region geometry if required + if (m_postprocess_calibration) + { + if(!the_id.InitializeGeometry(m_mdtIdHelper, m_detMgr)) + { + log << MSG::ERROR << "Faild to get geometry for " << the_id.regionId()<<endreq; + } + } +// t0 calibration // + const T0CalibrationOutput* t0_output( + dynamic_cast<const T0CalibrationOutput*>(m_results)); + if (t0_output!=0) { + + + log << MSG::INFO << "Writing out t0s." << endreq; + MdtTubeFitContainer *new_t0s = t0_output->t0s(); + if(t0_output->GetMap().size()) + { + std::map<NtupleStationId, MdtTubeFitContainer*> :: const_iterator mit(t0_output->GetMap().find(*it)); + if(mit!=t0_output->GetMap().end()) + { + new_t0s=mit->second; + } + } + if(!new_t0s) continue; + if(m_postprocess_calibration) new_t0s = postprocess_t0s(new_t0s, the_id); + if(new_t0s == NULL) cerr<<"new_t0s == NULL"<<endl; + + sc=m_calib_output_tool->WriteT0(new_t0s, the_id, m_iov_start, m_iov_end); + if(sc.isFailure()) return sc; + } + + // r-t calibration // + bool real_resolution(true), real_rt(true); + const RtCalibrationOutput *rt_output( + dynamic_cast<const RtCalibrationOutput*>(m_results)); + if (rt_output!=0) { + if(rt_output->rt()) + { + if(!rt_output->rt()->HasTmaxDiff()) + { + const IRtRelation *old_rel=p_input_service->GetRtRelation(); + if(old_rel && old_rel->HasTmaxDiff()) + { + const_cast<IRtRelation *>(rt_output->rt())->SetTmaxDiff(old_rel->GetTmaxDiff()); + } + } + } + log << MSG::INFO + << "Writing out r-t relationships." << endreq; + + + if(m_resolution==NULL) + { + real_resolution=false; + create_default_resolution(rt_output->rt()); + } + else + { + real_rt=false; + } + sc=m_calib_output_tool->WriteRt(rt_output, m_resolution, the_id, m_iov_start, m_iov_end, real_rt, real_resolution); + + if(sc.isFailure()) return sc; + + } + +//----------------------------------------------------------------------------- + } +//----------------------------------------------------------------------------- + log << MSG::INFO << "Finished writing"<<endreq; + return StatusCode::SUCCESS; + +} + +//***************************************************************************** + +//:::::::::::::::::: +//:: METHOD reset :: +//:::::::::::::::::: + +void MdtCalibOutputDbSvc::reset() { + + m_results=0; + + return; + +} + + + +MuonCalib::MdtTubeFitContainer * MdtCalibOutputDbSvc::postprocess_t0s(MuonCalib::MdtTubeFitContainer * t0, const MuonCalib::NtupleStationId &id) + { + int n_layers(-1), n_tubes(-1); + if(id.NMultilayers() <0) return t0; +//get the maximum number of layers and tubes per layer for this station + for(int i=0; i<id.NMultilayers(); i++) + { + if(n_layers<id.NLayers(i)) n_layers=id.NLayers(i); + if(n_tubes<id.NTubes(i)) n_tubes=id.NTubes(i); + } + MdtTubeFitContainer *new_t0=new MdtTubeFitContainer(id.regionId(), id.NMultilayers(), n_layers, n_tubes); + new_t0->setImplementation(t0->implementation()); + new_t0->setGroupBy(t0->GroupBy()); +//loop over tubes + double n_goods(0.0); + double mean_t0(0.0); + bool has_bad_fits(false); + for(unsigned int ml=0; ml<t0->numMultilayers(); ml++) + { + if(static_cast<int>(ml)==id.NMultilayers()) break; + for(unsigned int ly=0; ly<t0->numLayers(); ly++) + { + if(static_cast<int>(ly)==n_layers) break; + for(unsigned int tb=0; tb<t0->numTubes(); tb++) + { + if(static_cast<int>(tb)==n_tubes) break; + new_t0->setFit(ml, ly, tb, *(t0->getFit(ml, ly, tb))); + new_t0->setCalib(ml, ly, tb, *(t0->getCalib(ml, ly, tb))); + //if tube was ok use for mean t0 + if(t0->getCalib(ml, ly, tb)->statusCode == 0) + { + n_goods++; + mean_t0+=t0->getCalib(ml, ly, tb)->t0; + } + else + has_bad_fits=true; + } + } + } +//check for bad fits and insert mean t0; + if(has_bad_fits) + { + mean_t0/=n_goods; + for(unsigned int ml=0; ml<new_t0->numMultilayers(); ml++) + for(unsigned int ly=0; ly<new_t0->numLayers(); ly++) + for(unsigned int tb=0; tb<new_t0->numTubes(); tb++) + { + if(t0->getCalib(ml, ly, tb)->statusCode != 0) + { + MdtTubeCalibContainer::SingleTubeCalib calib(*(new_t0->getCalib(ml, ly, tb))); + calib.t0=mean_t0; + new_t0->setCalib(ml, ly, tb, calib); + } + } + } + return new_t0; + } + + + +inline void MdtCalibOutputDbSvc :: create_default_resolution(const MuonCalib::IRtRelation *rt) + { + MsgStream log(messageService(), name()); +//check if resolution is saved in input service + const IRtResolution *old_res = p_input_service->GetResolution(); + const IRtRelation *old_rel=p_input_service->GetRtRelation(); + if(old_res!=NULL && old_rel!=NULL && !m_force_default_resolution) + { + //scale the old resolution to the new rt relation + log<< MSG::INFO << "Taken old resolution" << endreq; + std::vector<SamplePoint> res_points(100); + for (unsigned int i=0; i<100; i++) + { + double di=static_cast<double>(i); + double t_new=rt->tLower() + (di/100.0) * (rt->tUpper() - rt->tLower()); + double t_old=old_rel->tLower() + (di/100.0) * (old_rel->tUpper() - old_rel->tLower()); + res_points[i].set_x1(t_new); + res_points[i].set_x2(old_res->resolution(t_old)); + if(res_points[i].x2()>100 && i>0) res_points[i].set_x2(res_points[i-1].x2()); + res_points[i].set_error(1); + } + RtResolutionFromPoints respoints; + m_resolution=new RtResolutionLookUp(respoints.getRtResolutionLookUp(res_points)); + return; + } + log<< MSG::INFO << "Creating default resolution" << endreq; +//parameters for default resolution curve + double alpha[9] = { + 0.31476, + -0.217661, + 0.118337, + -0.0374466, + 0.00692553, + -0.000764969, + 4.97305E-05, + -1.75457E-06, + 2.59E-08 + }; + MuonCalib::CalibFunc::ParVec par_vec; + + if(m_flat_default_resolution<0.0) + { + double bin_width=(rt->tUpper() - rt->tLower())/100.0; + par_vec.resize(102); + par_vec[0]=rt->tLower(); + par_vec[1]=bin_width; + for(unsigned int i=0; i<100; i++) + { + double t=rt->tLower() + i*bin_width; + double r=rt->radius(t); + double resol(0.0); + for (int l=0; l<9; l++) { + resol = resol+alpha[l]*pow(r, l); + } + par_vec[i+2]=resol; + } + } + else + { + par_vec.resize(4); + par_vec[0]=rt->tLower(); + par_vec[1]=(rt->tUpper() - rt->tLower())/2.0; + par_vec[2]=m_flat_default_resolution; + par_vec[3]=m_flat_default_resolution; + } + m_resolution = new MuonCalib::RtResolutionLookUp(par_vec); + } diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_entries.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_entries.cxx new file mode 100644 index 0000000000000000000000000000000000000000..90da297a0671c4525ad8a70ad4ffcee7fc173ebc --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_entries.cxx @@ -0,0 +1,20 @@ +#include "GaudiKernel/DeclareFactoryEntries.h" +#include "MdtCalibIOSvc/MdtCalibInputSvc.h" +#include "MdtCalibIOSvc/MdtCalibOutputDbSvc.h" +#include "MdtCalibIOSvc/CalibrationFileIOTool.h" +#include "MdtCalibIOSvc/CalibrationOracleFileIOTool.h" + +using namespace MuonCalib; + +DECLARE_SERVICE_FACTORY ( MdtCalibOutputDbSvc ) +DECLARE_SERVICE_FACTORY ( MdtCalibInputSvc ) +DECLARE_TOOL_FACTORY ( CalibrationFileIOTool ) +DECLARE_TOOL_FACTORY ( CalibrationOracleFileIOTool ) + +DECLARE_FACTORY_ENTRIES(MdtCalibIOSvc) { + DECLARE_SERVICE(MdtCalibOutputDbSvc) + DECLARE_SERVICE(MdtCalibInputSvc) + DECLARE_TOOL(CalibrationFileIOTool) + DECLARE_TOOL(CalibrationOracleFileIOTool) +} + diff --git a/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_load.cxx b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_load.cxx new file mode 100644 index 0000000000000000000000000000000000000000..1f241329ff845c9efa14a2e43444d5c21648a159 --- /dev/null +++ b/MuonSpectrometer/MuonCalib/MuonCalibStandAlone/MdtCalibIOSvc/src/components/MdtCalibIOSvc_load.cxx @@ -0,0 +1,3 @@ +#include "GaudiKernel/LoadFactoryEntries.h" + +LOAD_FACTORY_ENTRIES(MdtCalibIOSvc)