packages/micromegas/MicromegasCalibrationData.cc

0001 /*!
0002  * \file MicromegasCalibrationData.cc
0003  * \author Hugo Pereira Da Costa <hugo.pereira-da-costa@cea.fr>
0004  */
0005
0006 #include "MicromegasCalibrationData.h"
0007 #include "MicromegasMapping.h"
0008
0009 #include <cdbobjects/CDBTTree.h>
0010
0011 #include <TFile.h>
0012 #include <fstream>
0013 #include <sstream>
0014 #include <iostream>
0015
0016 namespace
0017 {
0018   static const std::string m_pedestal_key = "pedestal";
0019   static const std::string m_rms_key = "rms";
0020 }
0021
0022 //________________________________________________________________________-
0023 std::ostream& operator << (std::ostream& out, const MicromegasCalibrationData& calib_data )
0024 {
0025   out << "MicromegasCalibrationData" << std::endl;
0026   size_t total_entries = 0;
0027   for( const auto& [fee_id, data_array] : calib_data.m_raw_calibration_map )
0028   {
0029     total_entries += data_array.size();
0030     out << "fee_id: " << fee_id << " entries: " << data_array.size() << std::endl;
0031     for( size_t i=0; i<data_array.size(); ++i )
0032     {
0033       const auto& data = data_array[i];
0034       out << "fee_id: " << fee_id << " channel: " << i << " pedestal: " << data.m_pedestal << " rms: " << data.m_rms << std::endl;
0035     }
0036
0037   }
0038   out << "total entries: " << total_entries << std::endl;
0039
0040
0041   return out;
0042 }
0043
0044 //________________________________________________________________________-
0045 void MicromegasCalibrationData::read( const std::string& filename )
0046 {
0047   std::cout << "MicromegasCalibrationData::read - filename: " << filename << std::endl;
0048
0049   // clear existing data
0050   m_raw_calibration_map.clear();
0051   m_mapped_calibration_map.clear();
0052
0053   // make sure file exists before loading, otherwise crashes
0054   if( !std::ifstream( filename.c_str() ).good() )
0055   {
0056     std::cout << "MicromegasCalibrationData::read -"
0057       << " filename: " << filename << " does not exist."
0058       << " No calibration loaded" << std::endl;
0059     return;
0060   }
0061
0062   // use generic CDBTree to load
0063   CDBTTree cdbttree( filename );
0064   cdbttree.LoadCalibrations();
0065
0066   // loop over registered fee ids
0067   MicromegasMapping mapping;
0068   const auto fee_id_list( mapping.get_fee_id_list() );
0069
0070   // loop over all possible fee ids
0071   static constexpr int m_fee_count_max = 26;
0072   for( int fee_id = 0; fee_id < m_fee_count_max; ++ fee_id )
0073   {
0074
0075     // convert to new ids, for backward compatibility, and check if listed
0076     int fee_id_new = mapping.get_new_fee_id( fee_id );
0077     if( std::find(fee_id_list.begin(), fee_id_list.end(), fee_id_new) == fee_id_list.end() )
0078     { continue; }
0079
0080     // get corresponding hitset key
0081     const auto hitsetkey = mapping.get_hitsetkey(fee_id_new);
0082
0083     // loop over channels
0084     for( int i = 0; i < m_nchannels_fee; ++i )
0085     {
0086       // read pedestal and rms
0087       int channel = fee_id*m_nchannels_fee + i;
0088       double pedestal = cdbttree.GetDoubleValue( channel, m_pedestal_key, 0 );
0089       double rms = cdbttree.GetDoubleValue( channel, m_rms_key, 0 );
0090
0091       // insert in local structure
0092       if( !std::isnan( rms ) )
0093       {
0094         // create calibration data
0095         const calibration_data_t calibration_data( pedestal, rms );
0096
0097         // insert in fee,channel structure
0098         m_raw_calibration_map[fee_id_new].at(i) = calibration_data;
0099
0100         // also insert in hitsetkey, strip structure
0101         const auto strip = mapping.get_physical_strip( fee_id_new, i );
0102         if( hitsetkey > 0 && strip >= 0 )
0103         { m_mapped_calibration_map[hitsetkey].at(strip) = calibration_data; }
0104       }
0105     }
0106   }
0107 }
0108
0109 //________________________________________________________________________-
0110 void MicromegasCalibrationData::set_pedestal( int fee, int channel, double value )
0111 { m_raw_calibration_map[fee].at(channel).m_pedestal = value; }
0112
0113 //________________________________________________________________________-
0114 void MicromegasCalibrationData::set_rms( int fee, int channel, double value )
0115 { m_raw_calibration_map[fee].at(channel).m_rms = value; }
0116
0117 //________________________________________________________________________-
0118 void MicromegasCalibrationData::write( const std::string& filename ) const
0119 {
0120   std::cout << "MicromegasCalibrationData::write - filename: " << filename << std::endl;
0121   if( m_raw_calibration_map.empty() ) { return;
0122 }
0123
0124   // use generic CDBTree to load
0125   CDBTTree cdbttree( filename );
0126   for( const auto& [fee,array]:m_raw_calibration_map )
0127   {
0128     for( size_t i = 0; i < array.size(); ++i )
0129     {
0130       int channel = static_cast<unsigned long>(fee*m_nchannels_fee) + i;
0131       const auto& pedestal =  array[i].m_pedestal;
0132       const auto& rms =  array[i].m_rms;
0133       cdbttree.SetDoubleValue( channel, m_pedestal_key, pedestal );
0134       cdbttree.SetDoubleValue( channel, m_rms_key, rms );
0135     }
0136   }
0137
0138   // commit and write
0139   cdbttree.Commit();
0140   cdbttree.WriteCDBTTree();
0141 }
0142
0143 //________________________________________________________________________-
0144 double MicromegasCalibrationData::get_pedestal( int fee, int channel ) const
0145 {
0146   const auto iter = m_raw_calibration_map.find(fee);
0147   return (iter != m_raw_calibration_map.end()) ? iter->second.at(channel).m_pedestal: -1;
0148 }
0149
0150 //________________________________________________________________________-
0151 double MicromegasCalibrationData::get_rms( int fee, int channel ) const
0152 {
0153   const auto iter = m_raw_calibration_map.find(fee);
0154   return (iter != m_raw_calibration_map.end()) ? iter->second.at(channel).m_rms: -1;
0155 }
0156
0157 //________________________________________________________________________-
0158 double MicromegasCalibrationData::get_pedestal_mapped( TrkrDefs::hitsetkey hitsetkey, int strip ) const
0159 {
0160   const auto iter = m_mapped_calibration_map.find(hitsetkey);
0161   return (iter != m_mapped_calibration_map.end()) ? iter->second.at(strip).m_pedestal: -1;
0162 }
0163
0164 //________________________________________________________________________-
0165 double MicromegasCalibrationData::get_rms_mapped( TrkrDefs::hitsetkey hitsetkey, int strip ) const
0166 {
0167   const auto iter = m_mapped_calibration_map.find(hitsetkey);
0168   return (iter != m_mapped_calibration_map.end()) ? iter->second.at(strip).m_rms: -1;
0169 }