Ecal/plugins/StoreESCondition.cc

0001 #include "CondTools/Ecal/plugins/StoreESCondition.h"
0002
0003 #include "FWCore/ServiceRegistry/interface/Service.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008
0009 #include <fstream>
0010 #include <string>
0011 #include <cstring>
0012 #include <ctime>
0013 #include <unistd.h>
0014
0015 using std::string;
0016
0017 StoreESCondition::StoreESCondition(const edm::ParameterSet& iConfig) {
0018   prog_name_ = "StoreESCondition";
0019
0020   logfile_ = iConfig.getParameter<std::string>("logfile");
0021
0022   esgain_ = iConfig.getParameter<unsigned int>("gain");
0023
0024   typedef std::vector<edm::ParameterSet> Parameters;
0025   Parameters toPut = iConfig.getParameter<Parameters>("toPut");
0026   for (Parameters::iterator itToPut = toPut.begin(); itToPut != toPut.end(); ++itToPut) {
0027     inpFileName_.push_back(itToPut->getUntrackedParameter<std::string>("inputFile"));
0028     objectName_.push_back(itToPut->getUntrackedParameter<std::string>("conditionType"));
0029     since_.push_back(itToPut->getUntrackedParameter<unsigned int>("since"));
0030   }
0031 }
0032
0033 void StoreESCondition::endJob() {
0034   edm::Service<cond::service::PoolDBOutputService> mydbservice;
0035   if (!mydbservice.isAvailable()) {
0036     edm::LogError("StoreESCondition") << "PoolDBOutputService is unavailable"
0037                                       << "\n";
0038     return;
0039   }
0040
0041   bool toAppend = false;
0042   // copy a string to the char *
0043   std::string message = "finished OK\n";
0044   size_t messageSize = message.size() + 1;
0045   char* messChar = new char[messageSize];
0046   strncpy(messChar, message.c_str(), messageSize);
0047
0048   for (unsigned int i = 0; i < objectName_.size(); i++) {
0049     cond::Time_t newTime;
0050
0051     if (mydbservice->isNewTagRequest(objectName_[i] + std::string("Rcd"))) {
0052       // This is the first object for this tag.
0053       // Append mode should be off.
0054       // newTime is the end of this new objects IOV.
0055       newTime = mydbservice->beginOfTime();
0056     } else {
0057       // There should already be an object in the DB for this tag.
0058       // Append IOV mode should be on.
0059       // newTime is the beginning of this new objects IOV.
0060       toAppend = true;
0061       newTime = (cond::Time_t)since_[i];
0062     }
0063     edm::LogInfo("StoreESCondition") << "Reading " << objectName_[i] << " from file and writing to DB with newTime "
0064                                      << newTime << "\n";
0065     if (objectName_[i] == "ESChannelStatus") {
0066       edm::LogInfo("StoreESCondition") << " ESChannelStatus file " << inpFileName_[i] << "\n";
0067       const auto mycali = readESChannelStatusFromFile(inpFileName_[i].c_str());
0068       edm::LogInfo("StoreESCondition") << " ESChannelStatus file read "
0069                                        << "\n";
0070       if (!toAppend) {
0071         edm::LogInfo("StoreESCondition") << " before create "
0072                                          << "\n";
0073         mydbservice->createOneIOV<ESChannelStatus>(*mycali, newTime, "ESChannelStatusRcd");
0074         edm::LogInfo("StoreESCondition") << " after create "
0075                                          << "\n";
0076       } else {
0077         edm::LogInfo("StoreESCondition") << " before append "
0078                                          << "\n";
0079         mydbservice->appendOneIOV<ESChannelStatus>(*mycali, newTime, "ESChannelStatusRcd");
0080         edm::LogInfo("StoreESCondition") << " after append "
0081                                          << "\n";
0082       }
0083     } else if (objectName_[i] == "ESIntercalibConstants") {
0084       const auto myintercalib = readESIntercalibConstantsFromFile(inpFileName_[i].c_str());
0085       if (!toAppend) {
0086         mydbservice->createOneIOV<ESIntercalibConstants>(*myintercalib, newTime, "ESIntercalibConstantsRcd");
0087       } else {
0088         mydbservice->appendOneIOV<ESIntercalibConstants>(*myintercalib, newTime, "ESIntercalibConstantsRcd");
0089       }
0090     } else if (objectName_[i] == "ESTimeSampleWeights") {
0091       const auto myintercalib = readESTimeSampleWeightsFromFile(inpFileName_[i].c_str());
0092       if (!toAppend) {
0093         mydbservice->createOneIOV<ESTimeSampleWeights>(*myintercalib, newTime, "ESTimeSampleWeightsRcd");
0094       } else {
0095         mydbservice->appendOneIOV<ESTimeSampleWeights>(*myintercalib, newTime, "ESTimeSampleWeightsRcd");
0096       }
0097     } else if (objectName_[i] == "ESGain") {
0098       const auto myintercalib = readESGainFromFile(inpFileName_[i].c_str());
0099       if (!toAppend) {
0100         mydbservice->createOneIOV<ESGain>(*myintercalib, newTime, "ESGainRcd");
0101       } else {
0102         mydbservice->appendOneIOV<ESGain>(*myintercalib, newTime, "ESGainRcd");
0103       }
0104     } else if (objectName_[i] == "ESMissingEnergyCalibration") {
0105       const auto myintercalib = readESMissingEnergyFromFile(inpFileName_[i].c_str());
0106       if (!toAppend) {
0107         mydbservice->createOneIOV<ESMissingEnergyCalibration>(*myintercalib, newTime, "ESMissingEnergyCalibrationRcd");
0108       } else {
0109         mydbservice->appendOneIOV<ESMissingEnergyCalibration>(*myintercalib, newTime, "ESMissingEnergyCalibrationRcd");
0110       }
0111     } else if (objectName_[i] == "ESRecHitRatioCuts") {
0112       const auto myintercalib = readESRecHitRatioCutsFromFile(inpFileName_[i].c_str());
0113       if (!toAppend) {
0114         mydbservice->createOneIOV<ESRecHitRatioCuts>(*myintercalib, newTime, "ESRecHitRatioCutsRcd");
0115       } else {
0116         mydbservice->appendOneIOV<ESRecHitRatioCuts>(*myintercalib, newTime, "ESRecHitRatioCutsRcd");
0117       }
0118     } else if (objectName_[i] == "ESThresholds") {
0119       const auto myintercalib = readESThresholdsFromFile(inpFileName_[i].c_str());
0120       if (!toAppend) {
0121         mydbservice->createOneIOV<ESThresholds>(*myintercalib, newTime, "ESThresholdsRcd");
0122       } else {
0123         mydbservice->appendOneIOV<ESThresholds>(*myintercalib, newTime, "ESThresholdsRcd");
0124       }
0125     } else if (objectName_[i] == "ESPedestals") {
0126       const auto myintercalib = readESPedestalsFromFile(inpFileName_[i].c_str());
0127       if (!toAppend) {
0128         mydbservice->createOneIOV<ESPedestals>(*myintercalib, newTime, "ESPedestalsRcd");
0129       } else {
0130         mydbservice->appendOneIOV<ESPedestals>(*myintercalib, newTime, "ESPedestalsRcd");
0131       }
0132     } else if (objectName_[i] == "ESEEIntercalibConstants") {
0133       const auto myintercalib = readESEEIntercalibConstantsFromFile(inpFileName_[i].c_str());
0134       if (!toAppend) {
0135         mydbservice->createOneIOV<ESEEIntercalibConstants>(*myintercalib, newTime, "ESEEIntercalibConstantsRcd");
0136       } else {
0137         mydbservice->appendOneIOV<ESEEIntercalibConstants>(*myintercalib, newTime, "ESEEIntercalibConstantsRcd");
0138       }
0139     } else {
0140       edm::LogError("StoreESCondition") << "Object " << objectName_[i] << " is not supported by this program.";
0141     }
0142     // if more records write here else if ....
0143
0144     writeToLogFileResults(messChar);
0145
0146     edm::LogInfo("StoreESCondition") << "Finished endJob";
0147   }
0148
0149   delete[] messChar;
0150 }
0151
0152 StoreESCondition::~StoreESCondition() {}
0153
0154 void StoreESCondition::analyze(const edm::Event& evt, const edm::EventSetup& evtSetup) {}
0155
0156 void StoreESCondition::writeToLogFile(string a, string b, unsigned long long since) {
0157   FILE* outFile;  // output log file for appending
0158   outFile = fopen(logfile_.c_str(), "a");
0159   if (!outFile) {
0160     edm::LogError("StoreESCondition") << "*** Can not open file: " << logfile_;
0161     return;
0162   }
0163   char header[256];
0164   fillHeader(header);
0165   char appendMode[10];
0166   if (since != 0)
0167     strcpy(appendMode, "append");
0168   else
0169     strcpy(appendMode, "create");
0170
0171   //fprintf(outFile, "%s %s condition from file %s written into DB for SM %d (mapped to SM %d) in %s mode (since run %u)\n",
0172   //header, a.c_str(), b.c_str(),  appendMode, (unsigned int)since);
0173
0174   fclose(outFile);  // close out file
0175 }
0176
0177 void StoreESCondition::writeToLogFileResults(char* arg) {
0178   FILE* outFile;  // output log file for appending
0179   outFile = fopen(logfile_.c_str(), "a");
0180   if (!outFile) {
0181     edm::LogError("StoreESCondition") << "*** Can not open file: " << logfile_;
0182     return;
0183   }
0184   char header[256];
0185   fillHeader(header);
0186   fprintf(outFile, "%s %s\n", header, arg);
0187   fclose(outFile);  // close out file
0188 }
0189
0190 void StoreESCondition::fillHeader(char* header) {
0191   time_t rawtime;
0192   struct tm* timeinfo;
0193   time(&rawtime);
0194   timeinfo = localtime(&rawtime);
0195   char user[50];
0196   strcpy(user, getlogin());
0197   strcpy(header, asctime(timeinfo));
0198   strcpy(header, user);
0199 }
0200
0201 std::shared_ptr<ESThresholds> StoreESCondition::readESThresholdsFromFile(const char* inputFile) {
0202   std::ifstream ESThresholdsFile(edm::FileInPath(inputFile).fullPath().c_str());
0203   float ts2, zs;  //2nd time sample, ZS threshold
0204   ESThresholdsFile >> ts2;
0205   ESThresholdsFile >> zs;
0206   auto esThresholds = std::make_shared<ESThresholds>(ts2, zs);
0207
0208   return esThresholds;
0209 }
0210
0211 std::shared_ptr<ESEEIntercalibConstants> StoreESCondition::readESEEIntercalibConstantsFromFile(const char* inputFile) {
0212   std::ifstream ESEEIntercalibFile(edm::FileInPath(inputFile).fullPath().c_str());
0213   float gammaLow0, alphaLow0, gammaHigh0, alphaHigh0, gammaLow1, alphaLow1, gammaHigh1, alphaHigh1, gammaLow2,
0214       alphaLow2, gammaHigh2, alphaHigh2, gammaLow3, alphaLow3, gammaHigh3, alphaHigh3;
0215   ESEEIntercalibFile >> gammaLow0;
0216   ESEEIntercalibFile >> alphaLow0;
0217   ESEEIntercalibFile >> gammaHigh0;
0218   ESEEIntercalibFile >> alphaHigh0;
0219   ESEEIntercalibFile >> gammaLow1;
0220   ESEEIntercalibFile >> alphaLow1;
0221   ESEEIntercalibFile >> gammaHigh1;
0222   ESEEIntercalibFile >> alphaHigh1;
0223   ESEEIntercalibFile >> gammaLow2;
0224   ESEEIntercalibFile >> alphaLow2;
0225   ESEEIntercalibFile >> gammaHigh2;
0226   ESEEIntercalibFile >> alphaHigh2;
0227   ESEEIntercalibFile >> gammaLow3;
0228   ESEEIntercalibFile >> alphaLow3;
0229   ESEEIntercalibFile >> gammaHigh3;
0230   ESEEIntercalibFile >> alphaHigh3;
0231   auto eseeIntercalibConstants = std::make_shared<ESEEIntercalibConstants>(gammaLow0,
0232                                                                            alphaLow0,
0233                                                                            gammaHigh0,
0234                                                                            alphaHigh0,
0235                                                                            gammaLow1,
0236                                                                            alphaLow1,
0237                                                                            gammaHigh1,
0238                                                                            alphaHigh1,
0239                                                                            gammaLow2,
0240                                                                            alphaLow2,
0241                                                                            gammaHigh2,
0242                                                                            alphaHigh2,
0243                                                                            gammaLow3,
0244                                                                            alphaLow3,
0245                                                                            gammaHigh3,
0246                                                                            alphaHigh3);
0247
0248   return eseeIntercalibConstants;
0249 }
0250
0251 std::shared_ptr<ESMissingEnergyCalibration> StoreESCondition::readESMissingEnergyFromFile(const char* inputFile) {
0252   std::ifstream ESMissingEnergyFile(edm::FileInPath(inputFile).fullPath().c_str());
0253   float ConstAEta0, ConstBEta0, ConstAEta1, ConstBEta1, ConstAEta2, ConstBEta2, ConstAEta3, ConstBEta3;
0254   ESMissingEnergyFile >> ConstAEta0;
0255   ESMissingEnergyFile >> ConstBEta0;
0256   ESMissingEnergyFile >> ConstAEta1;
0257   ESMissingEnergyFile >> ConstBEta1;
0258   ESMissingEnergyFile >> ConstAEta2;
0259   ESMissingEnergyFile >> ConstBEta2;
0260   ESMissingEnergyFile >> ConstAEta3;
0261   ESMissingEnergyFile >> ConstBEta3;
0262   auto esMissingEnergy = std::make_shared<ESMissingEnergyCalibration>(
0263       ConstAEta0, ConstBEta0, ConstAEta1, ConstBEta1, ConstAEta2, ConstBEta2, ConstAEta3, ConstBEta3);
0264
0265   return esMissingEnergy;
0266 }
0267
0268 std::shared_ptr<ESPedestals> StoreESCondition::readESPedestalsFromFile(const char* inputFile) {
0269   auto esPedestals = std::make_shared<ESPedestals>();
0270
0271   // int ped[2][2][40][40][32];
0272   // for (int i=0; i<2; ++i)
0273   //  for (int j=0; j<2; ++j)
0274   //    for (int k=0; k<40; ++k)
0275   //    for (int m=0; m<40; ++m)
0276   //      for (int n=0; n<32; ++n)
0277   //        ped[i][j][k][m][n] = 0;
0278
0279   int ped[ESDetId::IZ_NUM][ESDetId::PLANE_MAX][ESDetId::IX_MAX][ESDetId::IY_MAX][ESDetId::ISTRIP_MAX] = {};
0280
0281   int iz, ip, ix, iy, is, ped_, zside;
0282   std::ifstream pedestalFile(edm::FileInPath(inputFile).fullPath().c_str());
0283
0284   for (int i = 0; i < 137216; ++i) {
0285     pedestalFile >> iz >> ip >> ix >> iy >> is >> ped_;
0286
0287     zside = (iz == -1) ? 1 : 0;
0288     ped[zside][ip - 1][ix - 1][iy - 1][is - 1] = ped_;
0289   }
0290
0291   for (int iz = -1; iz <= 1; ++iz) {
0292     if (iz == 0)
0293       continue;
0294     zside = (iz == -1) ? 1 : 0;
0295
0296     for (int iplane = ESDetId::PLANE_MIN; iplane <= ESDetId::PLANE_MAX; ++iplane)
0297       for (int ix = ESDetId::IX_MIN; ix <= ESDetId::IX_MAX; ++ix)
0298         for (int iy = ESDetId::IY_MIN; iy <= ESDetId::IY_MAX; ++iy)
0299           for (int istrip = ESDetId::ISTRIP_MIN; istrip <= ESDetId::ISTRIP_MAX; ++istrip) {
0300             ESPedestals::Item ESitem;
0301             ESitem.mean = ped[zside][iplane - 1][ix - 1][iy - 1][istrip - 1];
0302             ESitem.rms = 3;  // LG : 3, HG : 6
0303
0304             if (ESDetId::validDetId(istrip, ix, iy, iplane, iz)) {
0305               ESDetId esId(istrip, ix, iy, iplane, iz);
0306               esPedestals->insert(std::make_pair(esId.rawId(), ESitem));
0307             }
0308           }
0309   }
0310
0311   return esPedestals;
0312 }
0313
0314 std::shared_ptr<ESRecHitRatioCuts> StoreESCondition::readESRecHitRatioCutsFromFile(const char* inputFile) {
0315   std::ifstream ESRecHitRatioCutsFile(edm::FileInPath(inputFile).fullPath().c_str());
0316
0317   float r12Low, r23Low, r12High, r23High;
0318   ESRecHitRatioCutsFile >> r12Low;
0319   ESRecHitRatioCutsFile >> r23Low;
0320   ESRecHitRatioCutsFile >> r12High;
0321   ESRecHitRatioCutsFile >> r23High;
0322   auto esRecHitRatioCuts = std::make_shared<ESRecHitRatioCuts>(r12Low, r23Low, r12High, r23High);
0323
0324   return esRecHitRatioCuts;
0325 }
0326
0327 std::shared_ptr<ESGain> StoreESCondition::readESGainFromFile(const char* inputFile) {
0328   std::ifstream amplFile(edm::FileInPath(inputFile).fullPath().c_str());
0329
0330   int gain;
0331   amplFile >> gain;
0332   edm::LogInfo("StoreESCondition") << "gain : " << gain << "\n";
0333
0334   auto esGain = std::make_shared<ESGain>(gain);  // 1: LG, 2: HG
0335   return esGain;
0336 }
0337
0338 std::shared_ptr<ESTimeSampleWeights> StoreESCondition::readESTimeSampleWeightsFromFile(const char* inputFile) {
0339   std::ifstream amplFile(edm::FileInPath(inputFile).fullPath().c_str());
0340
0341   float w[3];
0342   for (int k = 0; k < 3; ++k) {
0343     float ww;
0344     amplFile >> ww;
0345     w[k] = ww;
0346     edm::LogInfo("StoreESCondition") << "weight : " << k << " " << w[k] << "\n";
0347   }
0348
0349   auto esWeights = std::make_shared<ESTimeSampleWeights>(w[0], w[1], w[2]);
0350   return esWeights;
0351 }
0352
0353 std::shared_ptr<ESIntercalibConstants> StoreESCondition::readESIntercalibConstantsFromFile(const char* inputFile) {
0354   auto ical = std::make_shared<ESIntercalibConstants>();
0355
0356   std::ifstream mipFile(edm::FileInPath(inputFile).fullPath().c_str());
0357
0358   for (int i = 0; i < 137216; ++i) {
0359     int iz, ip, ix, iy, is;
0360     double mip;
0361     mipFile >> iz >> ip >> ix >> iy >> is >> mip;
0362     //if (mip <20 || mip> 70) cout<<iz<<" "<<ip<<" "<<ix<<" "<<iy<<" "<<is<<" "<<mip<<endl; // HG
0363     // LG : HG MIP/6/1.14
0364     //mip = mip/6/1.14;
0365     // LG : HG MIP/6
0366     if (esgain_ == 1)
0367       mip = mip / 6.;  // LG
0368     if (mip < 20 || mip > 70)
0369       edm::LogInfo("StoreESCondition") << iz << " " << ip << " " << ix << " " << iy << " " << is << " " << mip
0370                                        << "\n";  // LG
0371
0372     if (ESDetId::validDetId(is, ix, iy, ip, iz)) {
0373       ESDetId esId(is, ix, iy, ip, iz);
0374       ical->setValue(esId.rawId(), mip);
0375     }
0376   }
0377
0378   return ical;
0379 }
0380
0381 std::shared_ptr<ESChannelStatus> StoreESCondition::readESChannelStatusFromFile(const char* inputFile) {
0382   int z[1000], p[1000], x[1000], y[1000], nsensors;
0383   std::ifstream statusFile(edm::FileInPath(inputFile).fullPath().c_str());
0384   statusFile >> nsensors;
0385   edm::LogInfo("StoreESCondition") << " nsensors " << nsensors << "\n";
0386   if (nsensors >= 1000) {
0387     edm::LogInfo("StoreESCondition") << " *** value too high, modify the method!***"
0388                                      << "\n";
0389     exit(-1);
0390   }
0391   for (int i = 0; i < nsensors; ++i) {
0392     statusFile >> z[i] >> p[i] >> x[i] >> y[i];
0393   }
0394   auto ecalStatus = std::make_shared<ESChannelStatus>();
0395   int Nbstatus = 0, Nbstrip = 0;
0396   for (int istrip = ESDetId::ISTRIP_MIN; istrip <= ESDetId::ISTRIP_MAX; istrip++) {
0397     for (int ix = ESDetId::IX_MIN; ix <= ESDetId::IX_MAX; ix++) {
0398       for (int iy = ESDetId::IY_MIN; iy <= ESDetId::IY_MAX; iy++) {
0399         for (int iplane = 1; iplane <= 2; iplane++) {
0400           for (int izeta = -1; izeta <= 1; izeta = izeta + 2) {
0401             //      if (izeta==0) continue;
0402             //      try {
0403
0404             //ESDetId Plane iplane Zside izeta
0405             //      if(!ESDetId::validDetId(istrip,ix,iy,iplane,izeta)) cout << " Unvalid DetId" << endl;
0406             //      else {
0407             if (ESDetId::validDetId(istrip, ix, iy, iplane, izeta)) {
0408               ESDetId anESId(istrip, ix, iy, iplane, izeta);
0409               int status = 0;
0410               //          std::ifstream statusFile(edm::FileInPath(inputFile).fullPath().c_str());
0411               Nbstrip++;
0412               for (int i = 0; i < nsensors; ++i) {
0413                 if (izeta == z[i] && iplane == p[i] && ix == x[i] && iy == y[i])
0414                   status = 1;
0415               }
0416               if (status == 1) {
0417                 Nbstatus++;
0418                 if (istrip == 1)
0419                   edm::LogInfo("StoreESCondition") << " Bad channel ix " << ix << " iy " << iy << " iplane " << iplane
0420                                                    << " iz " << izeta << "\n";  // print only once
0421               }
0422               ecalStatus->setValue(anESId, status);
0423               //          statusFile.close();
0424             }  // valid DetId
0425             //      catch ( cms::Exception &e ) { }
0426           }  // loop over z
0427         }    //  loop over plane
0428       }      //   loop over y
0429     }        //    loop over x
0430   }          //     loop over strips
0431   edm::LogInfo("StoreESCondition") << " Nb of strips " << Nbstrip << " Number of bad channels " << Nbstatus << "\n";
0432   statusFile.close();
0433
0434   // overwrite the statuses which are in the file
0435   return ecalStatus;
0436 }