Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 10:42:25

0001 
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \author G. Cerminara - INFN Torino
0006  */
0007 
0008 #include <iostream>
0009 
0010 #include "DumpFileToDB.h"
0011 #include "DTCalibrationMap.h"
0012 
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 
0015 #include "FWCore/Framework/interface/ESHandle.h"
0016 #include "FWCore/Framework/interface/EventSetup.h"
0017 #include "FWCore/Utilities/interface/Exception.h"
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 
0020 #include "CondFormats/DTObjects/interface/DTMtime.h"
0021 #include "CondFormats/DTObjects/interface/DTTtrig.h"
0022 #include "CondFormats/DataRecord/interface/DTTtrigRcd.h"
0023 #include "CondFormats/DTObjects/interface/DTT0.h"
0024 #include "CondFormats/DTObjects/interface/DTStatusFlag.h"
0025 #include "CondFormats/DTObjects/interface/DTDeadFlag.h"
0026 #include "CondFormats/DTObjects/interface/DTReadOutMapping.h"
0027 #include "CondFormats/DTObjects/interface/DTRecoConditions.h"
0028 #include "CondFormats/DataRecord/interface/DTRecoConditionsTtrigRcd.h"
0029 #include "CondFormats/DataRecord/interface/DTRecoConditionsVdriftRcd.h"
0030 #include "CondFormats/DataRecord/interface/DTRecoConditionsUncertRcd.h"
0031 #include "CalibMuon/DTCalibration/interface/DTCalibDBUtils.h"
0032 
0033 using namespace edm;
0034 using namespace std;
0035 
0036 DumpFileToDB::DumpFileToDB(const ParameterSet& pset) {
0037   dbToDump = pset.getUntrackedParameter<string>("dbToDump", "TTrigDB");
0038   format = pset.getUntrackedParameter<string>("dbFormat", "Legacy");
0039   ttrigToken_ = esConsumes<edm::Transition::BeginRun>();
0040 
0041   cout << "Writing DB: " << dbToDump << " with format: " << format << endl;
0042 
0043   if (dbToDump != "ChannelsDB")
0044     theCalibFile = new DTCalibrationMap(pset.getUntrackedParameter<ParameterSet>("calibFileConfig"));
0045 
0046   mapFileName = pset.getUntrackedParameter<ParameterSet>("calibFileConfig")
0047                     .getUntrackedParameter<string>("calibConstFileName", "dummy.txt");
0048 
0049   if (dbToDump != "VDriftDB" && dbToDump != "TTrigDB" && dbToDump != "TZeroDB" && dbToDump != "NoiseDB" &&
0050       dbToDump != "DeadDB" && dbToDump != "ChannelsDB" && dbToDump != "RecoUncertDB")
0051     throw cms::Exception("IncorrectSetup") << "Parameter dbToDump is not valid";
0052 
0053   if (format != "Legacy" && format != "DTRecoConditions")
0054     throw cms::Exception("IncorrectSetup") << "Parameter format is not valid";
0055 
0056   if (format == "DTRecoConditions" && (dbToDump != "VDriftDB" && dbToDump != "TTrigDB" && dbToDump != "RecoUncertDB"))
0057     throw cms::Exception("IncorrectSetup")
0058         << "DTRecoConditions currently implemented only for TTrigDB, VDriftDB, RecoUncertDB";
0059 
0060   diffMode = pset.getUntrackedParameter<bool>("differentialMode", false);
0061   if (diffMode) {
0062     if (dbToDump != "TTrigDB") {
0063       throw cms::Exception("IncorrectSetup") << "Differential mode currently implemented for ttrig only";
0064     }
0065     cout << "Using differential mode: mean value of txt table will be added to the current DB value" << endl;
0066   }
0067 }
0068 
0069 DumpFileToDB::~DumpFileToDB() {}
0070 
0071 void DumpFileToDB::endJob() {
0072   //---------- VDrift
0073   if (dbToDump == "VDriftDB") {
0074     if (format == "Legacy") {
0075       DTMtime mtime;
0076       for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0077            keyAndCalibs != theCalibFile->keyAndConsts_end();
0078            ++keyAndCalibs) {
0079         //  cout << "key: " << (*keyAndCalibs).first
0080         //       << " vdrift (cm/ns): " << theCalibFile->meanVDrift((*keyAndCalibs).first)
0081         //       << " hit reso (cm): " << theCalibFile->sigma_meanVDrift((*keyAndCalibs).first) << endl;
0082         // vdrift is cm/ns , resolution is cm
0083         mtime.set((*keyAndCalibs).first.superlayerId(),
0084                   theCalibFile->meanVDrift((*keyAndCalibs).first),
0085                   theCalibFile->sigma_meanVDrift((*keyAndCalibs).first),
0086                   DTVelocityUnits::cm_per_ns);
0087       }
0088       DTCalibDBUtils::writeToDB<DTMtime>("DTMtimeRcd", mtime);
0089     } else if (format == "DTRecoConditions") {
0090       DTRecoConditions conds;
0091       conds.setFormulaExpr("[0]");
0092       //      conds.setFormulaExpr("[0]*(1-[1]*x)");
0093       int version = 1;
0094       conds.setVersion(version);
0095       for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0096            keyAndCalibs != theCalibFile->keyAndConsts_end();
0097            ++keyAndCalibs) {
0098         vector<float> values = (*keyAndCalibs).second;
0099         int fversion = int(values[10] / 1000);
0100         int type = (int(values[10]) % 1000) / 100;
0101         int nfields = int(values[10]) % 100;
0102         if (type != 1)
0103           throw cms::Exception("IncorrectSetup") << "Input file is not for VDriftDB";
0104         if (values.size() != unsigned(nfields + 11))
0105           throw cms::Exception("IncorrectSetup") << "Inconsistent number of fields";
0106         if (fversion != version)
0107           throw cms::Exception("IncorrectSetup") << "Inconsistent version of file";
0108 
0109         vector<double> params(values.begin() + 11, values.begin() + 11 + nfields);
0110         conds.set((*keyAndCalibs).first, params);
0111       }
0112       DTCalibDBUtils::writeToDB<DTRecoConditions>("DTRecoConditionsVdriftRcd", conds);
0113     }
0114 
0115     //---------- TTrig
0116   } else if (dbToDump == "TTrigDB") {
0117     if (format == "Legacy") {
0118       DTTtrig tTrig;
0119       for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0120            keyAndCalibs != theCalibFile->keyAndConsts_end();
0121            ++keyAndCalibs) {
0122         if (diffMode) {  // sum the correction in the txt file (for the mean value) to what is input DB
0123           // retrieve the previous value from the DB
0124           float tmean;
0125           float trms;
0126           float kFactor;
0127           // ttrig and rms are ns
0128           tTrigMapOrig->get((*keyAndCalibs).first.rawId(), tmean, trms, kFactor, DTTimeUnits::ns);
0129           if (theCalibFile->kFactor((*keyAndCalibs).first) != 0 || kFactor != 0) {
0130             throw cms::Exception("IncorrectSetup")
0131                 << "The differentialMode can only be used with kFactor = 0, old: " << kFactor
0132                 << " new: " << theCalibFile->kFactor((*keyAndCalibs).first);
0133           }
0134           tTrig.set((*keyAndCalibs).first.superlayerId(),
0135                     theCalibFile->tTrig((*keyAndCalibs).first) + tmean,
0136                     theCalibFile->sigma_tTrig((*keyAndCalibs).first),
0137                     theCalibFile->kFactor((*keyAndCalibs).first),
0138                     DTTimeUnits::ns);
0139           //      cout << "key: " << (*keyAndCalibs).first
0140           //           << " ttrig_mean (ns): " << theCalibFile->tTrig((*keyAndCalibs).first) + tmean
0141           //           << " ttrig_sigma(ns): " << theCalibFile->sigma_tTrig((*keyAndCalibs).first)
0142           //           << " kFactor: " << theCalibFile->kFactor((*keyAndCalibs).first) << endl;
0143 
0144         } else {  // Normal mode
0145                   //      cout << "key: " << (*keyAndCalibs).first
0146                   //           << " ttrig_mean (ns): " << theCalibFile->tTrig((*keyAndCalibs).first)
0147                   //           << " ttrig_sigma(ns): " << theCalibFile->sigma_tTrig((*keyAndCalibs).first)
0148                   //           << " kFactor: " << theCalibFile->kFactor((*keyAndCalibs).first) << endl;
0149           tTrig.set((*keyAndCalibs).first.superlayerId(),
0150                     theCalibFile->tTrig((*keyAndCalibs).first),
0151                     theCalibFile->sigma_tTrig((*keyAndCalibs).first),
0152                     theCalibFile->kFactor((*keyAndCalibs).first),
0153                     DTTimeUnits::ns);
0154         }
0155       }
0156       DTCalibDBUtils::writeToDB<DTTtrig>("DTTtrigRcd", tTrig);
0157 
0158     } else if (format == "DTRecoConditions") {
0159       DTRecoConditions conds;
0160       conds.setFormulaExpr("[0]");
0161       int version = 1;
0162       conds.setVersion(version);
0163 
0164       for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0165            keyAndCalibs != theCalibFile->keyAndConsts_end();
0166            ++keyAndCalibs) {
0167         vector<float> values = (*keyAndCalibs).second;
0168         int fversion = int(values[10] / 1000);
0169         int type = (int(values[10]) % 1000) / 100;
0170         int nfields = int(values[10]) % 100;
0171 
0172         if (type != 1)
0173           throw cms::Exception("IncorrectSetup") << "Only type==1 currently supported for TTrigDB";
0174         if (values.size() != unsigned(nfields + 11))
0175           throw cms::Exception("IncorrectSetup") << "Inconsistent number of fields";
0176         if (fversion != version)
0177           throw cms::Exception("IncorrectSetup") << "Inconsistent version of file";
0178 
0179         vector<double> params(values.begin() + 11, values.begin() + 11 + nfields);
0180         conds.set((*keyAndCalibs).first, params);
0181       }
0182       DTCalibDBUtils::writeToDB<DTRecoConditions>("DTRecoConditionsTtrigRcd", conds);
0183     }
0184 
0185   } else if (dbToDump == "TZeroDB") {  // Write the T0
0186 
0187     // Create the object to be written to DB
0188     DTT0 tZeroMap;
0189 
0190     // Loop over file entries
0191     for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0192          keyAndCalibs != theCalibFile->keyAndConsts_end();
0193          ++keyAndCalibs) {
0194       float t0mean = (*keyAndCalibs).second[5];
0195       float t0rms = (*keyAndCalibs).second[6];
0196       cout << "key: " << (*keyAndCalibs).first << " T0 mean (TDC counts): " << t0mean
0197            << " T0_rms (TDC counts): " << t0rms << endl;
0198       tZeroMap.set((*keyAndCalibs).first, t0mean, t0rms, DTTimeUnits::counts);
0199     }
0200 
0201     DTCalibDBUtils::writeToDB<DTT0>("DTT0Rcd", tZeroMap);
0202 
0203   } else if (dbToDump == "NoiseDB") {  // Write the Noise
0204     DTStatusFlag statusMap;
0205 
0206     // Loop over file entries
0207     for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0208          keyAndCalibs != theCalibFile->keyAndConsts_end();
0209          ++keyAndCalibs) {
0210       cout << "key: " << (*keyAndCalibs).first << " Noisy flag: " << (*keyAndCalibs).second[7] << endl;
0211       statusMap.setCellNoise((*keyAndCalibs).first, (*keyAndCalibs).second[7]);
0212     }
0213 
0214     DTCalibDBUtils::writeToDB<DTStatusFlag>("DTStatusFlagRcd", statusMap);
0215 
0216   } else if (dbToDump == "DeadDB") {  // Write the tp-dead
0217     DTDeadFlag deadMap;
0218 
0219     // Loop over file entries
0220     for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0221          keyAndCalibs != theCalibFile->keyAndConsts_end();
0222          ++keyAndCalibs) {
0223       cout << "key: " << (*keyAndCalibs).first << " dead flag: " << (*keyAndCalibs).second[7] << endl;
0224       deadMap.setCellDead_TP((*keyAndCalibs).first, (*keyAndCalibs).second[7]);
0225     }
0226 
0227     DTCalibDBUtils::writeToDB<DTDeadFlag>("DTDeadFlagRcd", deadMap);
0228 
0229   } else if (dbToDump == "ChannelsDB") {  //Write channels map
0230 
0231     DTReadOutMapping ro_map("cmssw_ROB", "cmssw_ROS");
0232     //Loop over file entries
0233     string line;
0234     ifstream file(mapFileName.c_str());
0235     while (getline(file, line)) {
0236       if (line == "" || line[0] == '#')
0237         continue;  // Skip comments and empty lines
0238       stringstream linestr;
0239       linestr << line;
0240       vector<int> channelMap = readChannelsMap(linestr);
0241       int status = ro_map.insertReadOutGeometryLink(channelMap[0],
0242                                                     channelMap[1],
0243                                                     channelMap[2],
0244                                                     channelMap[3],
0245                                                     channelMap[4],
0246                                                     channelMap[5],
0247                                                     channelMap[6],
0248                                                     channelMap[7],
0249                                                     channelMap[8],
0250                                                     channelMap[9],
0251                                                     channelMap[10]);
0252       cout << "ddu " << channelMap[0] << " "
0253            << "ros " << channelMap[1] << " "
0254            << "rob " << channelMap[2] << " "
0255            << "tdc " << channelMap[3] << " "
0256            << "channel " << channelMap[4] << " "
0257            << "wheel " << channelMap[5] << " "
0258            << "station " << channelMap[6] << " "
0259            << "sector " << channelMap[7] << " "
0260            << "superlayer " << channelMap[8] << " "
0261            << "layer " << channelMap[9] << " "
0262            << "wire " << channelMap[10] << " "
0263            << "  -> ";
0264       cout << "insert status: " << status << std::endl;
0265     }
0266     DTCalibDBUtils::writeToDB<DTReadOutMapping>("DTReadOutMappingRcd", ro_map);
0267 
0268     //---------- Uncertainties
0269   } else if (dbToDump == "RecoUncertDB") {  // Write the Uncertainties
0270 
0271     DTRecoConditions conds;
0272     conds.setFormulaExpr("par[step]");
0273     int version = 1;  // Uniform uncertainties per SL and step; parameters 0-3 are for steps 1-4.
0274     conds.setVersion(version);
0275 
0276     for (DTCalibrationMap::const_iterator keyAndCalibs = theCalibFile->keyAndConsts_begin();
0277          keyAndCalibs != theCalibFile->keyAndConsts_end();
0278          ++keyAndCalibs) {
0279       vector<float> values = (*keyAndCalibs).second;
0280       int fversion = int(values[10] / 1000);
0281       int type = (int(values[10]) % 1000) / 100;
0282       int nfields = int(values[10]) % 100;
0283       if (type != 2)
0284         throw cms::Exception("IncorrectSetup") << "Only type==2 supported for uncertainties DB";
0285       if (values.size() != unsigned(nfields + 11))
0286         throw cms::Exception("IncorrectSetup") << "Inconsistent number of fields";
0287       if (fversion != version)
0288         throw cms::Exception("IncorrectSetup") << "Inconsistent version of file";
0289 
0290       vector<double> params(values.begin() + 11, values.begin() + 11 + nfields);
0291       conds.set((*keyAndCalibs).first, params);
0292       DTCalibDBUtils::writeToDB<DTRecoConditions>("DTRecoConditionsUncertRcd", conds);
0293     }
0294   }
0295 }
0296 
0297 vector<int> DumpFileToDB::readChannelsMap(stringstream& linestr) {
0298   //The hardware channel
0299   int ddu_id = 0;
0300   int ros_id = 0;
0301   int rob_id = 0;
0302   int tdc_id = 0;
0303   int channel_id = 0;
0304   //The software channel
0305   int wheel_id = 0;
0306   int station_id = 0;
0307   int sector_id = 0;
0308   int superlayer_id = 0;
0309   int layer_id = 0;
0310   int wire_id = 0;
0311 
0312   linestr >> ddu_id >> ros_id >> rob_id >> tdc_id >> channel_id >> wheel_id >> station_id >> sector_id >>
0313       superlayer_id >> layer_id >> wire_id;
0314 
0315   vector<int> channelMap;
0316   channelMap.push_back(ddu_id);
0317   channelMap.push_back(ros_id);
0318   channelMap.push_back(rob_id);
0319   channelMap.push_back(tdc_id);
0320   channelMap.push_back(channel_id);
0321   channelMap.push_back(wheel_id);
0322   channelMap.push_back(station_id);
0323   channelMap.push_back(sector_id);
0324   channelMap.push_back(superlayer_id);
0325   channelMap.push_back(layer_id);
0326   channelMap.push_back(wire_id);
0327   return channelMap;
0328 }
0329 
0330 void DumpFileToDB::beginRun(const edm::Run& run, const edm::EventSetup& setup) {
0331   if (diffMode) {
0332     if (dbToDump == "TTrigDB") {  // read the original DB
0333       ESHandle<DTTtrig> tTrig = setup.getHandle(ttrigToken_);
0334       tTrigMapOrig = &*tTrig;
0335       cout << "[DumpDBToFile] TTrig version: " << tTrig->version() << endl;
0336     }
0337   }
0338 }