Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:47:32

0001 #ifndef _CSCCROSSTALKDBCONDITIONS_H
0002 #define _CSCCROSSTALKDBCONDITIONS_H
0003 
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/Framework/interface/ESProducer.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "FWCore/Framework/interface/EventSetupRecordIntervalFinder.h"
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/Framework/interface/SourceFactory.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include <cmath>
0014 #include <memory>
0015 
0016 #include "CondFormats/CSCObjects/interface/CSCDBCrosstalk.h"
0017 #include "CondFormats/DataRecord/interface/CSCDBCrosstalkRcd.h"
0018 #include <DataFormats/MuonDetId/interface/CSCDetId.h>
0019 
0020 class CSCCrosstalkDBConditions : public edm::ESProducer, public edm::EventSetupRecordIntervalFinder {
0021 public:
0022   CSCCrosstalkDBConditions(const edm::ParameterSet &);
0023   ~CSCCrosstalkDBConditions() override;
0024 
0025   inline static CSCDBCrosstalk *prefillDBCrosstalk();
0026 
0027   typedef std::unique_ptr<CSCDBCrosstalk> ReturnType;
0028 
0029   ReturnType produceDBCrosstalk(const CSCDBCrosstalkRcd &);
0030 
0031 private:
0032   // ----------member data ---------------------------
0033   void setIntervalFor(const edm::eventsetup::EventSetupRecordKey &,
0034                       const edm::IOVSyncValue &,
0035                       edm::ValidityInterval &) override;
0036 };
0037 
0038 #include <fstream>
0039 #include <iostream>
0040 #include <vector>
0041 
0042 // to workaround plugin library
0043 inline CSCDBCrosstalk *CSCCrosstalkDBConditions::prefillDBCrosstalk() {
0044   // const int MAX_SIZE = 273024; //for ME1a unganged
0045   const int MAX_SIZE = 252288;
0046   const int SLOPE_FACTOR = 10000000;
0047   const int INTERCEPT_FACTOR = 100000;
0048   const int MAX_SHORT = 32767;
0049   CSCDBCrosstalk *cndbcrosstalk = new CSCDBCrosstalk();
0050 
0051   int db_index, new_index;
0052   float db_slope_right, db_slope_left, db_intercept_right;
0053   float db_intercept_left;
0054   std::vector<int> db_index_id;
0055   std::vector<float> db_slope_r;
0056   std::vector<float> db_intercept_r;
0057   std::vector<float> db_slope_l;
0058   std::vector<float> db_intercept_l;
0059   float new_slope_right, new_slope_left, new_intercept_right;
0060   float new_intercept_left;
0061   std::vector<int> new_index_id;
0062   std::vector<float> new_slope_r;
0063   std::vector<float> new_intercept_r;
0064   std::vector<float> new_slope_l;
0065   std::vector<float> new_intercept_l;
0066 
0067   int counter;
0068   int db_nrlines = 0;
0069   int new_nrlines = 0;
0070 
0071   std::ifstream dbdata;
0072   dbdata.open("old_dbxtalk.dat", std::ios::in);
0073   if (!dbdata) {
0074     std::cerr << "Error: old_dbxtalk.dat -> no such file!" << std::endl;
0075     exit(1);
0076   }
0077 
0078   while (!dbdata.eof()) {
0079     dbdata >> db_index >> db_slope_right >> db_intercept_right >> db_slope_left >> db_intercept_left;
0080     db_index_id.push_back(db_index);
0081     db_slope_r.push_back(db_slope_right);
0082     db_slope_l.push_back(db_slope_left);
0083     db_intercept_r.push_back(db_intercept_right);
0084     db_intercept_l.push_back(db_intercept_left);
0085     db_nrlines++;
0086   }
0087   dbdata.close();
0088 
0089   std::ifstream newdata;
0090   newdata.open("xtalk.dat", std::ios::in);
0091   if (!newdata) {
0092     std::cerr << "Error: xtalk.dat -> no such file!" << std::endl;
0093     exit(1);
0094   }
0095 
0096   while (!newdata.eof()) {
0097     newdata >> new_index >> new_slope_right >> new_intercept_right >> new_slope_left >> new_intercept_left;
0098     new_index_id.push_back(new_index);
0099     new_slope_r.push_back(new_slope_right);
0100     new_slope_l.push_back(new_slope_left);
0101     new_intercept_r.push_back(new_intercept_right);
0102     new_intercept_l.push_back(new_intercept_left);
0103     new_nrlines++;
0104   }
0105   newdata.close();
0106 
0107   CSCDBCrosstalk::CrosstalkContainer &itemvector = cndbcrosstalk->crosstalk;
0108   itemvector.resize(MAX_SIZE);
0109   cndbcrosstalk->factor_slope = int(SLOPE_FACTOR);
0110   cndbcrosstalk->factor_intercept = int(INTERCEPT_FACTOR);
0111 
0112   for (int i = 0; i < MAX_SIZE; ++i) {
0113     itemvector[i].xtalk_slope_right = (short int)(db_slope_r[i] * SLOPE_FACTOR + 0.5);
0114     itemvector[i].xtalk_intercept_right = (short int)(db_intercept_r[i] * INTERCEPT_FACTOR + 0.5);
0115     itemvector[i].xtalk_slope_left = (short int)(db_slope_l[i] * SLOPE_FACTOR + 0.5);
0116     itemvector[i].xtalk_intercept_left = (short int)(db_intercept_l[i] * INTERCEPT_FACTOR + 0.5);
0117   }
0118 
0119   for (int i = 0; i < MAX_SIZE; ++i) {
0120     counter = db_index_id[i];
0121     for (unsigned int k = 0; k < new_index_id.size() - 1; k++) {
0122       if (counter == new_index_id[k]) {
0123         if ((short int)(fabs(new_slope_r[k] * SLOPE_FACTOR + 0.5)) < MAX_SHORT)
0124           itemvector[counter].xtalk_slope_right = int(new_slope_r[k] * SLOPE_FACTOR + 0.5);
0125         if ((short int)(fabs(new_intercept_r[k] * INTERCEPT_FACTOR + 0.5)) < MAX_SHORT)
0126           itemvector[counter].xtalk_intercept_right = int(new_intercept_r[k] * INTERCEPT_FACTOR + 0.5);
0127         if ((short int)(fabs(new_slope_l[k] * SLOPE_FACTOR + 0.5)) < MAX_SHORT)
0128           itemvector[counter].xtalk_slope_left = int(new_slope_l[k] * SLOPE_FACTOR + 0.5);
0129         if ((short int)(fabs(new_intercept_l[k] * INTERCEPT_FACTOR + 0.5)) < MAX_SHORT)
0130           itemvector[counter].xtalk_intercept_left = int(new_intercept_l[k] * INTERCEPT_FACTOR + 0.5);
0131         itemvector[i] = itemvector[counter];
0132         // std::cout<<" counter "<<counter <<" dbindex "<<new_index_id[k]<<"
0133         // dbslope " <<db_slope_r[k]<<" new slope "<<new_slope_r[k]<<std::endl;
0134       }
0135     }
0136 
0137     if (counter > 223968) {
0138       itemvector[counter].xtalk_slope_right = int(db_slope_r[i]);
0139       itemvector[counter].xtalk_slope_left = int(db_slope_l[i]);
0140       itemvector[counter].xtalk_intercept_right = int(db_intercept_r[i]);
0141       itemvector[counter].xtalk_intercept_left = int(db_intercept_l[i]);
0142     }
0143   }
0144 
0145   return cndbcrosstalk;
0146 }
0147 
0148 #endif