test/deprecated/CSCTriggerPrimitivesReader.cc

0001 //-------------------------------------------------
0002 //
0003 //   Class: CSCTriggerPrimitivesReader
0004 //
0005 //   Description: Basic analyzer class which accesses ALCTs, CLCTs, and
0006 //                correlated LCTs and plot various quantities.
0007 //
0008 //   Author List: S. Valuev, UCLA.
0009 //
0010 //
0011 //   Modifications:
0012 //
0013 //--------------------------------------------------
0014
0015 //-----------------------
0016 // This Class's Header --
0017 //-----------------------
0018 #include "L1Trigger/CSCTriggerPrimitives/plugins/CSCTriggerPrimitivesReader.h"
0019
0020 //-------------------------------
0021 // Collaborating Class Headers --
0022 //-------------------------------
0023
0024 #include "FWCore/Framework/interface/ConsumesCollector.h"
0025 #include "FWCore/Framework/interface/MakerMacros.h"
0026 #include "DataFormats/Common/interface/Handle.h"
0027 #include "FWCore/Framework/interface/ESHandle.h"
0028 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0029
0030 #include "CondFormats/CSCObjects/interface/CSCBadChambers.h"
0031
0032 // MC data
0033 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0034
0035 // MC tests
0036 #include "L1Trigger/CSCTriggerPrimitives/plugins/CSCAnodeLCTAnalyzer.h"
0037 #include "L1Trigger/CSCTriggerPrimitives/plugins/CSCCathodeLCTAnalyzer.h"
0038 #include "Geometry/CSCGeometry/interface/CSCLayer.h"
0039 #include "Geometry/CSCGeometry/interface/CSCGeometry.h"
0040
0041 #include "TCanvas.h"
0042 #include "TFile.h"
0043 #include "TH1.h"
0044 #include "TText.h"
0045 #include "TPaveLabel.h"
0046 #include "TPostScript.h"
0047 #include "TStyle.h"
0048 #include "TROOT.h"
0049
0050 #include "FWCore/ServiceRegistry/interface/Service.h"
0051 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0052
0053 #define MAXSTUBS 100
0054 using namespace std;
0055
0056 //-----------------
0057 // Static variables
0058 //-----------------
0059
0060 // Various useful constants
0061 const double CSCTriggerPrimitivesReader::TWOPI = 2. * M_PI;
0062 const string CSCTriggerPrimitivesReader::csc_type[CSC_TYPES] = {
0063     "ME1/1", "ME1/2", "ME1/3", "ME1/A", "ME2/1", "ME2/2", "ME3/1", "ME3/2", "ME4/1", "ME4/2"};
0064 const string CSCTriggerPrimitivesReader::csc_type_plus[CSC_TYPES] = {
0065     "ME+1/1", "ME+1/2", "ME+1/3", "ME+1/A", "ME+2/1", "ME+2/2", "ME+3/1", "ME+3/2", "ME+4/1", "ME+4/2"};
0066 const string CSCTriggerPrimitivesReader::csc_type_minus[CSC_TYPES] = {
0067     "ME-1/1", "ME-1/2", "ME-1/3", "ME-1/A", "ME-2/1", "ME-2/2", "ME-3/1", "ME-3/2", "ME-4/1", "ME-4/2"};
0068 const int CSCTriggerPrimitivesReader::NCHAMBERS[CSC_TYPES] = {36, 36, 36, 36, 18, 36, 18, 36, 18, 36};
0069 const int CSCTriggerPrimitivesReader::MAX_WG[CSC_TYPES] = {
0070     48, 64, 32, 48, 112, 64, 96, 64, 96, 64};  //max. number of wiregroups
0071 const int CSCTriggerPrimitivesReader::MAX_HS[CSC_TYPES] = {
0072     128, 160, 128, 96, 160, 160, 160, 160, 160, 160};  // max. # of halfstrips
0073 const int CSCTriggerPrimitivesReader::ptype_TMB07[CSCConstants::NUM_CLCT_PATTERNS] = {
0074     -999, -5, 4, -4, 3, -3, 2, -2, 1, -1, 0};  // "signed" pattern (== phiBend)
0075
0076 // LCT counters
0077 int CSCTriggerPrimitivesReader::numALCT = 0;
0078 int CSCTriggerPrimitivesReader::numCLCT = 0;
0079 int CSCTriggerPrimitivesReader::numLCTTMB = 0;
0080 int CSCTriggerPrimitivesReader::numLCTMPC = 0;
0081
0082 bool CSCTriggerPrimitivesReader::bookedHotWireHistos = false;
0083 bool CSCTriggerPrimitivesReader::bookedALCTHistos = false;
0084 bool CSCTriggerPrimitivesReader::bookedCLCTHistos = false;
0085 bool CSCTriggerPrimitivesReader::bookedLCTTMBHistos = false;
0086 bool CSCTriggerPrimitivesReader::bookedLCTMPCHistos = false;
0087
0088 bool CSCTriggerPrimitivesReader::bookedCompHistos = false;
0089
0090 bool CSCTriggerPrimitivesReader::bookedResolHistos = false;
0091 bool CSCTriggerPrimitivesReader::bookedEfficHistos = false;
0092
0093 bool CSCTriggerPrimitivesReader::printps = false;
0094
0095 // LUT for which ME1/1 wire group can cross which ME1/a halfstrip
0096 // 1st index: WG number
0097 // 2nd index: inclusive HS range
0098 const int CSCTriggerPrimitivesReader::lut_wg_vs_hs_me1a[48][2] = {
0099     {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},  {0, 95},
0100     {0, 95},  {0, 95},  {0, 77},  {0, 61},  {0, 39},  {0, 22},  {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0101     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0102     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0103     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}};
0104 // a modified LUT for ganged ME1a
0105 const int CSCTriggerPrimitivesReader::lut_wg_vs_hs_me1ag[48][2] = {
0106     {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},
0107     {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 31},  {0, 22},  {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0108     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0109     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0110     {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}};
0111
0112 // LUT for which ME1/1 wire group can cross which ME1/b halfstrip
0113 // 1st index: WG number
0114 // 2nd index: inclusive HS range
0115 const int CSCTriggerPrimitivesReader::lut_wg_vs_hs_me1b[48][2] = {
0116     {-1, -1},   {-1, -1},  {-1, -1},  {-1, -1},  {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1}, {-1, -1},
0117     {100, 127}, {73, 127}, {47, 127}, {22, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127},
0118     {0, 127},   {0, 127},  {0, 127},  {0, 127},  {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127},
0119     {0, 127},   {0, 127},  {0, 127},  {0, 127},  {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127}, {0, 127},
0120     {0, 127},   {0, 127},  {0, 127},  {0, 127},  {0, 105}, {0, 93},  {0, 78},  {0, 63}};
0121
0122 //init My ttree
0123 void TreePerStub::init(int run, int event) {
0124   t_EventNumberAnalyzed = -1;
0125   t_RUN = run;
0126   t_Event = event;
0127   t_nStubs = 0;
0128   t_nStubs_readout = 0;
0129
0130   t_chamber = -1;
0131   t_ring = -1;
0132   t_station = -1;
0133   t_chambertype = -1;
0134   t_endcap = -2;
0135 }
0136
0137 TTree* TreePerStub::bookTree(TTree* t, const std::string& name) {
0138   edm::Service<TFileService> fs;
0139   t = fs->make<TTree>(name.c_str(), name.c_str());
0140
0141   t->Branch("t_EventNumberAnalyzed", &t_EventNumberAnalyzed, "t_EventNumberAnalyzed/I");
0142   t->Branch("t_RUN", &t_RUN, "t_RUN/I");
0143   t->Branch("t_Event", &t_Event, "t_Event/I");
0144   t->Branch("t_nStubs", &t_nStubs, "t_nStubs/I");
0145   t->Branch("t_nStubs_readout", &t_nStubs_readout, "t_nStubs_readout/I");
0146   t->Branch("t_nStubs_ME119", &t_nStubs_ME119, "t_nStubs_ME119/I");
0147   t->Branch("t_nStubs_ME11no911", &t_nStubs_ME11no911, "t_nStubs_ME11no911/I");
0148   t->Branch("t_nStubs_noMEpm11", &t_nStubs_noMEpm11, "t_nStubs_noMEpm11/I");
0149
0150   t->Branch("t_chamber", &t_chamber, "t_chamber/I");
0151   t->Branch("t_ring", &t_ring, "t_ring/I");
0152   t->Branch("t_endcap", &t_endcap, "t_endcap/I");
0153   t->Branch("t_station", &t_station, "t_station/I");
0154   t->Branch("t_chambertype", &t_chambertype, "t_chambertype/I");
0155
0156   return t;
0157 }
0158
0159 //init stub comparison struct
0160 void MyStubComparison::init(int run, int event) {
0161   nEvents = -1;
0162   nRUN = run;
0163   nEvent = event;
0164   firstfill = false;
0165   totStubs_data = -1;
0166   totStubs_emul = -1;
0167   totStubs_emul_readout = -1;
0168   nStub_data = 0;
0169   nStub_emul = 0;
0170   chamber = -1;
0171   ring = -1;
0172   endcap = -1;
0173   station = -1;
0174   chambertype = -1;
0175   has_data = false;
0176   has_emul = false;
0177
0178   npretrig = 0;
0179   quality_pretrig = -1;
0180   maxquality_pretrig = -1;
0181   bend_pretrig = -1;
0182   bx_pretrig = -1;
0183   key_hs_pretrig = -1;
0184   pattern_pretrig = -1;
0185   maxpattern_pretrig = -1;
0186   quality_data = -1;
0187   bend_data = -1;
0188   bx_data = -1;
0189   quality_emul = -1;
0190   bend_emul = -1;
0191   pattern_data = -1;
0192   pattern_emul = -1;
0193   bx_emul = -1;
0194   bx_corr_emul = -1;  //corrected
0195   key_WG_data = -1;
0196   key_WG_emul = -1;
0197   key_hs_data = -1;
0198   key_hs_emul = -1;
0199   WGcrossHS_data = false;
0200   WGcrossHS_emul = false;
0201   trknmb_data = -1;
0202   trknmb_emul = -1;
0203   dphi_data = -1;
0204   dphi_emul = -1;
0205   eta_data = -1;
0206   eta_emul = -1;
0207   phi_data = -1;
0208   phi_emul = -1;
0209 }
0210
0211 TTree* MyStubComparison::bookTree(TTree* t, const std::string& name) {
0212   edm::Service<TFileService> fs;
0213   t = fs->make<TTree>(name.c_str(), name.c_str());
0214
0215   t->Branch("nEvents", &nEvents);
0216   t->Branch("nRUN", &nRUN);
0217   t->Branch("nEvent", &nEvent);
0218   t->Branch("firstfill", &firstfill);
0219   t->Branch("totStubs_data", &totStubs_data);
0220   t->Branch("totStubs_emul", &totStubs_emul);
0221   t->Branch("totStubs_emul_readout", &totStubs_emul_readout);
0222   t->Branch("nStub_data", &nStub_data);
0223   t->Branch("nStub_emul", &nStub_emul);
0224
0225   t->Branch("chamber", &chamber);
0226   t->Branch("ring", &ring);
0227   t->Branch("endcap", &endcap);
0228   t->Branch("station", &station);
0229   t->Branch("chambertype", &chambertype);
0230   t->Branch("has_data", &has_data);
0231   t->Branch("has_emul", &has_emul);
0232   t->Branch("quality_data", &quality_data);
0233   t->Branch("quality_emul", &quality_emul);
0234   t->Branch("npretrig", &npretrig);
0235   t->Branch("quality_pretrig", &quality_pretrig);
0236   t->Branch("maxquality_pretrig", &maxquality_pretrig);
0237   t->Branch("pattern_data", &pattern_data);
0238   t->Branch("pattern_emul", &pattern_emul);
0239   t->Branch("pattern_pretrig", &pattern_pretrig);
0240   t->Branch("maxpattern_pretrig", &maxpattern_pretrig);
0241   t->Branch("bend_data", &bend_data);
0242   t->Branch("bx_data", &bx_data);
0243   t->Branch("fullbx_data", &fullbx_data);
0244   t->Branch("bend_emul", &bend_emul);
0245   t->Branch("bx_emul", &bx_emul);
0246   t->Branch("fullbx_emul", &fullbx_emul);
0247   t->Branch("bend_pretrig", &bend_pretrig);
0248   t->Branch("bx_pretrig", &bx_pretrig);
0249   t->Branch("bx_corr_emul", &bx_corr_emul);
0250   t->Branch("WGcrossHS_data", &WGcrossHS_data);
0251   t->Branch("WGcrossHS_emul", &WGcrossHS_emul);
0252   t->Branch("key_WG_data", &key_WG_data);
0253   t->Branch("key_WG_emul", &key_WG_emul);
0254   t->Branch("key_hs_data", &key_hs_data);
0255   t->Branch("key_hs_emul", &key_hs_emul);
0256   t->Branch("key_hs_pretrig", &key_hs_pretrig);
0257   t->Branch("trknmb_data", &trknmb_data);
0258   t->Branch("trknmb_emul", &trknmb_emul);
0259   t->Branch("dphi_data", &dphi_data);
0260   t->Branch("dphi_emul", &dphi_emul);
0261   t->Branch("eta_data", &eta_data);
0262   t->Branch("eta_emul", &eta_emul);
0263   t->Branch("phi_data", &phi_data);
0264   t->Branch("phi_emul", &phi_emul);
0265
0266   return t;
0267 }
0268
0269 //----------------
0270 // Constructor  --
0271 //----------------
0272 CSCTriggerPrimitivesReader::CSCTriggerPrimitivesReader(const edm::ParameterSet& conf) : eventsAnalyzed(0) {
0273   edm::Service<TFileService> fs;
0274   //  rootFileName = conf.getUntrackedParameter<string>("rootFileName","TPEHists.root");
0275   // Create the root file for the histograms
0276   //  theFile = new TFile(rootFileName.c_str(), "RECREATE");
0277   //  theFile->cd();
0278
0279   // Various input parameters.
0280
0281   printps = conf.getParameter<bool>("printps");
0282   dataLctsIn_ = conf.getParameter<bool>("dataLctsIn");
0283   emulLctsIn_ = conf.getParameter<bool>("emulLctsIn");
0284   edm::ParameterSet commonParams = conf.getParameter<edm::ParameterSet>("commonParam");
0285
0286   // Switch for a new (2007) version of the TMB firmware.
0287   gangedME1a = commonParams.getParameter<bool>("gangedME1a");
0288
0289   // is it (non-upgrade algorithm) run along with upgrade one?
0290   plotME1A = true;
0291   plotME42 = true;
0292   lctProducerData_ = conf.getUntrackedParameter<string>("CSCLCTProducerData", "cscunpacker");
0293   mpclctProducerData_ = conf.getUntrackedParameter<string>("CSCMPCLCTProducerData", "csctfDigis");
0294   lctProducerEmul_ = conf.getUntrackedParameter<string>("CSCLCTProducerEmul", "cscTriggerPrimitiveDigis");
0295
0296   simHitProducer_ = conf.getParameter<edm::InputTag>("CSCSimHitProducer");
0297   wireDigiProducer_ = conf.getParameter<edm::InputTag>("CSCWireDigiProducer");
0298   compDigiProducer_ = conf.getParameter<edm::InputTag>("CSCComparatorDigiProducer");
0299
0300   simHit_token_ = consumes<edm::PSimHitContainer>(simHitProducer_);
0301   wireDigi_token_ = consumes<CSCWireDigiCollection>(wireDigiProducer_);
0302   compDigi_token_ = consumes<CSCComparatorDigiCollection>(compDigiProducer_);
0303
0304   alcts_d_token_ = consumes<CSCALCTDigiCollection>(edm::InputTag(lctProducerData_, "MuonCSCALCTDigi"));
0305   clcts_d_token_ = consumes<CSCCLCTDigiCollection>(edm::InputTag(lctProducerData_, "MuonCSCCLCTDigi"));
0306   lcts_tmb_d_token_ =
0307       consumes<CSCCorrelatedLCTDigiCollection>(edm::InputTag(lctProducerData_, "MuonCSCCorrelatedLCTDigi"));
0308   lcts_mpc_d_token_ = consumes<CSCCorrelatedLCTDigiCollection>(edm::InputTag(mpclctProducerData_));
0309
0310   alcts_e_token_ = consumes<CSCALCTDigiCollection>(edm::InputTag(lctProducerEmul_));
0311   clcts_e_token_ = consumes<CSCCLCTDigiCollection>(edm::InputTag(lctProducerEmul_));
0312   pretrigs_e_token_ = consumes<CSCCLCTPreTriggerDigiCollection>(edm::InputTag(lctProducerEmul_));
0313   lcts_tmb_e_token_ = consumes<CSCCorrelatedLCTDigiCollection>(edm::InputTag(lctProducerEmul_));
0314   lcts_mpc_e_token_ = consumes<CSCCorrelatedLCTDigiCollection>(edm::InputTag(lctProducerEmul_, "MPCSORTED"));
0315   cscGeomToken_ = esConsumes<CSCGeometry, MuonGeometryRecord>();
0316   pBadToken_ = esConsumes<CSCBadChambers, CSCBadChambersRcd>();
0317
0318   consumesMany<edm::HepMCProduct>();
0319   resultsFileNamesPrefix_ = conf.getUntrackedParameter<string>("resultsFileNamesPrefix", "");
0320   checkBadChambers_ = conf.getUntrackedParameter<bool>("checkBadChambers", true);
0321   debug = conf.getUntrackedParameter<bool>("debug", false);
0322   dataIsAnotherMC_ = conf.getUntrackedParameter<bool>("dataIsAnotherMC", false);
0323
0324   //rootFileName = conf.getUntrackedParameter<string>("rootFileName");
0325   // Create the root file.
0326   // Not sure we really need it - comment out for now. -Slava.
0327   //theFile = new TFile(rootFileName.c_str(), "RECREATE");
0328   //theFile->cd();
0329   stub_tree[0] = stubs_comparison[0].bookTree(stub_tree[0], "alcttree");
0330   stub_tree[1] = stubs_comparison[1].bookTree(stub_tree[1], "clcttree");
0331   stub_tree[2] = stubs_comparison[2].bookTree(stub_tree[2], "lcttree");
0332   stub_tree[3] = stubs_comparison[3].bookTree(stub_tree[3], "mpclcttree");
0333   // Per event TTree
0334   event_tree[0] = perStub[0].bookTree(event_tree[0], "Ev_alcttree");
0335   event_tree[1] = perStub[1].bookTree(event_tree[1], "Ev_emul_alcttree");
0336   event_tree[2] = perStub[2].bookTree(event_tree[2], "Ev_clcttree");
0337   event_tree[3] = perStub[3].bookTree(event_tree[3], "Ev_emul_clcttree");
0338   event_tree[4] = perStub[4].bookTree(event_tree[4], "Ev_lcttree");
0339   event_tree[5] = perStub[5].bookTree(event_tree[5], "Ev_emul_lcttree");
0340   // My favourite ROOT settings.
0341   setRootStyle();
0342 }
0343
0344 //----------------
0345 // Destructor   --
0346 //----------------
0347 CSCTriggerPrimitivesReader::~CSCTriggerPrimitivesReader() {
0348   //  histos->writeHists(theFile);
0349   //  theFile->Close();
0350   //delete theFile;
0351 }
0352
0353 int CSCTriggerPrimitivesReader::maxRing(int station) {
0354   if (station == 1) {
0355     if (plotME1A)
0356       return 4;
0357     else
0358       return 3;
0359   }
0360   return 2;
0361 }
0362
0363 void CSCTriggerPrimitivesReader::analyze(const edm::Event& ev, const edm::EventSetup& setup) {
0364   ++eventsAnalyzed;
0365   //if (ev.id().event()%10 == 0)
0366   LogTrace("CSCTriggerPrimitivesReader") << "\n** CSCTriggerPrimitivesReader: processing run #" << ev.id().run()
0367                                          << " event #" << ev.id().event() << "; events so far: " << eventsAnalyzed
0368                                          << " **";
0369   RUN_ = ev.id().run();
0370   Event_ = ev.id().event();
0371
0372   // Find the geometry for this event & cache it.  Needed in LCTAnalyzer
0373   // modules.
0374   edm::ESHandle<CSCGeometry> cscGeom = setup.getHandle(cscGeomToken_);
0375   geom_ = &*cscGeom;
0376
0377   // Find conditions data for bad chambers & cache it.  Needed for efficiency
0378   // calculations.
0379   if (checkBadChambers_) {
0380     edm::ESHandle<CSCBadChambers> pBad = setup.getHandle(pBadToken_);
0381     badChambers_ = pBad.product();
0382   }
0383
0384   // Get the collections of ALCTs, CLCTs, and correlated LCTs from event.
0385   edm::Handle<CSCALCTDigiCollection> alcts_data;
0386   edm::Handle<CSCCLCTDigiCollection> clcts_data;
0387   edm::Handle<CSCCorrelatedLCTDigiCollection> lcts_tmb_data;
0388   edm::Handle<CSCCorrelatedLCTDigiCollection> lcts_mpc_data;
0389   edm::Handle<CSCALCTDigiCollection> alcts_emul;
0390   edm::Handle<CSCCLCTDigiCollection> clcts_emul;
0391   edm::Handle<CSCCLCTPreTriggerDigiCollection> pretrigs_emul;
0392   edm::Handle<CSCCorrelatedLCTDigiCollection> lcts_tmb_emul;
0393   edm::Handle<CSCCorrelatedLCTDigiCollection> lcts_mpc_emul;
0394
0395   // Data
0396   if (dataLctsIn_) {
0397     HotWires(ev);
0398     //    ev.getByLabel(lctProducerData_,  alcts_data);
0399     //    ev.getByLabel(lctProducerData_,  clcts_data);
0400     //    ev.getByLabel(lctProducerData_,  lcts_tmb_data);
0401     ev.getByToken(alcts_d_token_, alcts_data);
0402     ev.getByToken(clcts_d_token_, clcts_data);
0403     ev.getByToken(lcts_tmb_d_token_, lcts_tmb_data);
0404     ev.getByToken(lcts_mpc_d_token_, lcts_mpc_data);
0405
0406     if (!alcts_data.isValid()) {
0407       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0408           << "+++ Warning: Collection of ALCTs with label MuonCSCALCTDigi"
0409           << " requested, but not found in the event... Skipping the rest +++\n";
0410       return;
0411     }
0412     if (!clcts_data.isValid()) {
0413       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0414           << "+++ Warning: Collection of CLCTs with label MuonCSCCLCTDigi"
0415           << " requested, but not found in the event... Skipping the rest +++\n";
0416       return;
0417     }
0418     if (!lcts_tmb_data.isValid()) {
0419       edm::LogWarning("L1CSCTPEmulatorWrongInput") << "+++ Warning: Collection of correlated LCTs with label"
0420                                                    << " MuonCSCCorrelatedLCTDigi requested, but not found in the"
0421                                                    << " event... Skipping the rest +++\n";
0422       return;
0423     }
0424     /*
0425       if (!lcts_mpc_data.isValid()) {
0426       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0427       << "+++ Warning: Collection of MPC correlated LCTs with label"
0428       << " MuonCSCCorrelatedLCTDigi + MCPSorted requested, but not found in the"
0429       << " event... Skipping the rest +++\n";
0430       //return;
0431       }*/
0432   }
0433
0434   // Emulator
0435   if (emulLctsIn_) {
0436     //    ev.getByLabel(lctProducerEmul_,              alcts_emul);
0437     //    ev.getByLabel(lctProducerEmul_,              clcts_emul);
0438     //    ev.getByLabel(lctProducerEmul_,              lcts_tmb_emul);
0439     //    ev.getByLabel(lctProducerEmul_, "MPCSORTED", lcts_mpc_emul);
0440     ev.getByToken(alcts_e_token_, alcts_emul);
0441     ev.getByToken(clcts_e_token_, clcts_emul);
0442     ev.getByToken(pretrigs_e_token_, pretrigs_emul);
0443     ev.getByToken(lcts_tmb_e_token_, lcts_tmb_emul);
0444     ev.getByToken(lcts_mpc_e_token_, lcts_mpc_emul);
0445
0446     if (!alcts_emul.isValid()) {
0447       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0448           << "+++ Warning: Collection of emulated ALCTs"
0449           << " requested, but not found in the event... Skipping the rest +++\n";
0450       return;
0451     }
0452     if (!clcts_emul.isValid()) {
0453       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0454           << "+++ Warning: Collection of emulated CLCTs"
0455           << " requested, but not found in the event... Skipping the rest +++\n";
0456       return;
0457     }
0458     if (!lcts_tmb_emul.isValid()) {
0459       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0460           << "+++ Warning: Collection of emulated correlated LCTs"
0461           << " requested, but not found in the event... Skipping the rest +++\n";
0462       return;
0463     }
0464     /*
0465       if (!lcts_mpc_emul.isValid()) {
0466       edm::LogWarning("L1CSCTPEmulatorWrongInput")
0467       << "+++ Warning: Collection of emulated correlated LCTs (MPCs)"
0468       << " requested, but not found in the event... Skipping the rest +++\n";
0469       //return;
0470       }*/
0471   }
0472
0473   // Fill histograms with reconstructed or emulated quantities.  If both are
0474   // present, plot LCTs in data.
0475   if (dataLctsIn_) {
0476     fillALCTHistos(alcts_data.product());
0477     fillCLCTHistos(clcts_data.product());
0478     fillLCTTMBHistos(lcts_tmb_data.product());
0479   } else if (emulLctsIn_) {
0480     fillALCTHistos(alcts_emul.product());
0481     fillCLCTHistos(clcts_emul.product());
0482     fillLCTTMBHistos(lcts_tmb_emul.product());
0483     //fillLCTMPCHistos(lcts_mpc_emul.product());
0484   }
0485
0486   // Compare LCTs in the data with the ones produced by the emulator.
0487   if (dataLctsIn_ && emulLctsIn_) {
0488     compare(alcts_data.product(),
0489             alcts_emul.product(),
0490             clcts_data.product(),
0491             clcts_emul.product(),
0492             pretrigs_emul.product(),
0493             lcts_tmb_data.product(),
0494             lcts_tmb_emul.product());
0495   }
0496   // Fill MC-based resolution/efficiency histograms, if needed.
0497   if (emulLctsIn_) {
0498     MCStudies(ev, alcts_emul.product(), clcts_emul.product());
0499   }
0500 }  // analyze
0501
0502 void CSCTriggerPrimitivesReader::endJob() {
0503   // Note: all operations involving ROOT should be placed here and not in the
0504   // destructor.
0505   // Plot histos if they were booked/filled.
0506   if (printps) {
0507     if (bookedALCTHistos)
0508       drawALCTHistos();
0509     if (bookedCLCTHistos)
0510       drawCLCTHistos();
0511     if (bookedLCTTMBHistos)
0512       drawLCTTMBHistos();
0513     //if (bookedLCTMPCHistos) drawLCTMPCHistos();
0514
0515     if (bookedCompHistos)
0516       drawCompHistos();
0517
0518     if (bookedResolHistos)
0519       drawResolHistos();
0520     if (bookedEfficHistos)
0521       drawEfficHistos();
0522   }
0523   //drawHistosForTalks();
0524
0525   //theFile->cd();
0526   //theFile->Write();
0527   //theFile->Close();
0528
0529   // Job summary.
0530   edm::LogInfo("CSCTriggerPrimitivesReader")
0531       << "\n  Average number of ALCTs/event = " << static_cast<float>(numALCT) / eventsAnalyzed << endl;
0532   edm::LogInfo("CSCTriggerPrimitivesReader")
0533       << "  Average number of CLCTs/event = " << static_cast<float>(numCLCT) / eventsAnalyzed << endl;
0534   edm::LogInfo("CSCTriggerPrimitivesReader")
0535       << "  Average number of TMB LCTs/event = " << static_cast<float>(numLCTTMB) / eventsAnalyzed << endl;
0536   edm::LogInfo("CSCTriggerPrimitivesReader")
0537       << "  Average number of MPC LCTs/event = " << static_cast<float>(numLCTMPC) / eventsAnalyzed << endl;
0538
0539   if (bookedEfficHistos) {
0540     {
0541       edm::LogInfo("CSCTriggerPrimitivesReader") << "\n  ALCT efficiencies:";
0542       double tot_simh = 0.0, tot_alct = 0.0;
0543       for (int idh = 0; idh < CSC_TYPES; idh++) {
0544         double simh = hEfficHitsEtaCsc[idh]->Integral();
0545         double alct = hEfficALCTEtaCsc[idh]->Integral();
0546         double eff = 0.;
0547         if (simh > 0)
0548           eff = alct / simh;
0549         edm::LogInfo("CSCTriggerPrimitivesReader")
0550             << "    " << csc_type[idh] << ": alct = " << alct << ", simh = " << simh << " eff = " << eff;
0551         tot_simh += simh;
0552         tot_alct += alct;
0553       }
0554       edm::LogInfo("CSCTriggerPrimitivesReader")
0555           << "    overall: alct = " << tot_alct << ", simh = " << tot_simh << " eff = " << tot_alct / tot_simh;
0556     }
0557
0558     {
0559       edm::LogInfo("CSCTriggerPrimitivesReader") << "\n  CLCT efficiencies:";
0560       double tot_simh = 0.0, tot_clct = 0.0;
0561       for (int idh = 0; idh < CSC_TYPES; idh++) {
0562         double simh = hEfficHitsEtaCsc[idh]->Integral();
0563         double clct = hEfficCLCTEtaCsc[idh]->Integral();
0564         double eff = 0.;
0565         if (simh > 0.)
0566           eff = clct / simh;
0567         edm::LogInfo("CSCTriggerPrimitivesReader")
0568             << "    " << csc_type[idh] << ": clct = " << clct << ", simh = " << simh << " eff = " << eff;
0569         tot_simh += simh;
0570         tot_clct += clct;
0571       }
0572       edm::LogInfo("CSCTriggerPrimitivesReader")
0573           << "    overall: clct = " << tot_clct << ", simh = " << tot_simh << " eff = " << tot_clct / tot_simh;
0574     }
0575   }
0576
0577   if (bookedResolHistos) {
0578     double cor = 0.0, tot = 0.0;
0579     cor = hResolDeltaHS->GetBinContent(hResolDeltaHS->FindBin(0.));
0580     tot = hResolDeltaHS->GetEntries();
0581     edm::LogInfo("CSCTriggerPrimitivesReader")
0582         << "\n  Correct half-strip assigned in " << cor << "/" << tot << " = " << cor / tot << " of half-strip CLCTs";
0583     //cor = hResolDeltaDS->GetBinContent(hResolDeltaDS->FindBin(0.));
0584     //tot = hResolDeltaDS->GetEntries();
0585     //edm::LogInfo("CSCTriggerPrimitivesReader")
0586     //  << "  Correct di-strip assigned in " << cor << "/" << tot
0587     //  << " = " << cor/tot << " of di-strip CLCTs";
0588     cor = hResolDeltaWG->GetBinContent(hResolDeltaWG->FindBin(0.));
0589     tot = hResolDeltaWG->GetEntries();
0590     edm::LogInfo("CSCTriggerPrimitivesReader")
0591         << "  Correct wire group assigned in " << cor << "/" << tot << " = " << cor / tot << " of ALCTs";
0592   }
0593 }
0594
0595 //---------------
0596 // ROOT settings
0597 //---------------
0598 void CSCTriggerPrimitivesReader::setRootStyle() {
0599   TH1::AddDirectory(false);
0600
0601   gROOT->SetStyle("Plain");
0602   gStyle->SetFillColor(0);
0603   gStyle->SetOptDate();
0604   gStyle->SetOptStat(111110);
0605   gStyle->SetOptFit(1111);
0606   gStyle->SetPadTickX(1);
0607   gStyle->SetPadTickY(1);
0608   gStyle->SetMarkerSize(0.5);
0609   gStyle->SetMarkerStyle(8);
0610   gStyle->SetGridStyle(3);
0611   gStyle->SetPaperSize(TStyle::kA4);
0612   gStyle->SetStatW(0.25);           // width of statistics box; default is 0.19
0613   gStyle->SetStatH(0.10);           // height of statistics box; default is 0.1
0614   gStyle->SetStatFormat("6.4g");    // leave default format for now
0615   gStyle->SetTitleSize(0.055, "");  // size for pad title; default is 0.02
0616   // Really big; useful for talks.
0617   //gStyle->SetTitleSize(0.1, "");   // size for pad title; default is 0.02
0618   gStyle->SetLabelSize(0.05, "XYZ");   // size for axis labels; default is 0.04
0619   gStyle->SetStatFontSize(0.06);       // size for stat. box
0620   gStyle->SetTitleFont(32, "XYZ");     // times-bold-italic font (p. 153) for axes
0621   gStyle->SetTitleFont(32, "");        // same for pad title
0622   gStyle->SetLabelFont(32, "XYZ");     // same for axis labels
0623   gStyle->SetStatFont(32);             // same for stat. box
0624   gStyle->SetLabelOffset(0.006, "Y");  // default is 0.005
0625 }
0626
0627 //---------------------
0628 // Histograms for LCTs
0629 //---------------------
0630 void CSCTriggerPrimitivesReader::bookHotWireHistos() {
0631   edm::Service<TFileService> fs;
0632   hHotWire1 = fs->make<TH1F>("hHotWire1", "hHotWire1", 570 * 6 * 112, 0, 570 * 6 * 112);
0633   hHotCham1 = fs->make<TH1F>("hHotCham1", "hHotCham1", 570, 0, 570);
0634   bookedHotWireHistos = true;
0635 }
0636
0637 void CSCTriggerPrimitivesReader::bookALCTHistos() {
0638   string s;
0639
0640   edm::Service<TFileService> fs;
0641   hAlctPerEvent = fs->make<TH1F>("ALCTs_per_event", "ALCTs per event", 31, -0.5, 30.5);
0642   hAlctPerChamber = fs->make<TH1F>("ALCTs_per_chamber", "ALCTs per chamber", 4, -0.5, 3.5);
0643   hAlctPerCSC = fs->make<TH1F>("ALCTs_per_CSCtype", "ALCTs per CSC type", 10, -0.5, 9.5);
0644   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0645     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0646       float csc_max = NCHAMBERS[j] + 0.5;
0647       char asdf[256];
0648       sprintf(asdf, "ALCTs_%i", i * CSC_TYPES + j);
0649       if (i == 0)
0650         s = "ALCTs, " + csc_type_plus[j];
0651       else
0652         s = "ALCTs, " + csc_type_minus[j];
0653       hAlctCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0654     }
0655   }
0656
0657   hAlctValid = fs->make<TH1F>("ALCT_validity", "ALCT validity", 3, -0.5, 2.5);
0658   hAlctQuality = fs->make<TH1F>("ALCT_quality", "ALCT quality", 5, -0.5, 4.5);
0659   hAlctAccel = fs->make<TH1F>("ALCT_accel_flag", "ALCT accel. flag", 3, -0.5, 2.5);
0660   hAlctCollis = fs->make<TH1F>("ALCT_collision_flag", "ALCT collision. flag", 3, -0.5, 2.5);
0661   hAlctKeyGroup = fs->make<TH1F>("ALCT_key_wiregroup", "ALCT key wiregroup", 120, -0.5, 119.5);
0662   hAlctBXN = fs->make<TH1F>("ALCT_bx", "ALCT bx", 20, -0.5, 19.5);
0663
0664   hAlctKeyGroupME11 = fs->make<TH1F>("hAlctKeyGroupME11", "ALCT key wiregroup ME1/1", 50, -0.5, 49.5);
0665
0666   bookedALCTHistos = true;
0667   std::cout << " bookedALCTHistos " << std::endl;
0668 }
0669
0670 void CSCTriggerPrimitivesReader::bookCLCTHistos() {
0671   string s;
0672
0673   edm::Service<TFileService> fs;
0674   hClctPerEvent = fs->make<TH1F>("CLCTs_per_event", "CLCTs per event", 31, -0.5, 30.5);
0675   hClctPerChamber = fs->make<TH1F>("CLCTs_per_chamber", "CLCTs per chamber", 3, -0.5, 2.5);
0676   hClctPerCSC = fs->make<TH1F>("CLCTs_per_CSCtype", "CLCTs per CSC type", 10, -0.5, 9.5);
0677   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0678     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0679       char asdf[256];
0680       sprintf(asdf, "CLCTs_%i", i * CSC_TYPES + j);
0681       float csc_max = NCHAMBERS[j] + 0.5;
0682       if (i == 0)
0683         s = "CLCTs, " + csc_type_plus[j];
0684       else
0685         s = "CLCTs, " + csc_type_minus[j];
0686       hClctCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0687     }
0688   }
0689
0690   hClctValid = fs->make<TH1F>("CLCT_validity", "CLCT validity", 3, -0.5, 2.5);
0691   hClctQuality = fs->make<TH1F>("CLCT_layers_hit", "CLCT layers hit", 9, -0.5, 8.5);
0692   hClctStripType = fs->make<TH1F>("CLCT_strip_type", "CLCT strip type", 3, -0.5, 2.5);
0693   hClctSign = fs->make<TH1F>("CLCT_sing_(L/R)", "CLCT sign (L/R)", 3, -0.5, 2.5);
0694   hClctCFEB = fs->make<TH1F>("CLCT_cfeb_#", "CLCT cfeb #", 6, -0.5, 5.5);
0695   hClctBXN = fs->make<TH1F>("CLCT_bx", "CLCT bx", 20, -0.5, 19.5);
0696
0697   hClctKeyStrip[0] = fs->make<TH1F>("CLCT_keystrip_distrips", "CLCT keystrip, distrips", 40, -0.5, 39.5);
0698   //hClctKeyStrip[0] = fs->make<TH1F>("","CLCT keystrip, distrips",  160, -0.5, 159.5);
0699   hClctKeyStrip[1] = fs->make<TH1F>("CLCT_keystrip_halfstrips", "CLCT keystrip, halfstrips", 160, -0.5, 159.5);
0700   hClctPattern[0] = fs->make<TH1F>("CLCT_pattern_distrips", "CLCT pattern, distrips", 13, -0.5, 12.5);
0701   hClctPattern[1] = fs->make<TH1F>("CLCT_pattern_halfstrips", "CLCT pattern, halfstrips", 13, -0.5, 12.5);
0702
0703   for (int i = 0; i < CSC_TYPES; i++) {
0704     char asdf[256];
0705     string s1 = "CLCT bend, " + csc_type[i];
0706     sprintf(asdf, "CLCT_bend0_%i", i + 1);
0707     hClctBendCsc[i][0] = fs->make<TH1F>(asdf, s1.c_str(), 5, -0.5, 4.5);
0708     sprintf(asdf, "CLCT_bend1_%i", i + 1);
0709     hClctBendCsc[i][1] = fs->make<TH1F>(asdf, s1.c_str(), 5, -0.5, 4.5);
0710
0711     sprintf(asdf, "CLCT_keystrip_%i", i + 1);
0712     string s2 = "CLCT keystrip, " + csc_type[i];
0713     int max_ds = MAX_HS[i] / 4;
0714     hClctKeyStripCsc[i] = fs->make<TH1F>(asdf, s2.c_str(), max_ds, 0., max_ds);
0715   }
0716
0717   hClctKeyStripME11 = fs->make<TH1F>("hClctKeyStripME11", "CLCT keystrip, halfstrips ME1/1", 161, -0.5, 160.5);
0718
0719   bookedCLCTHistos = true;
0720   std::cout << " bookedCLCTHistos " << std::endl;
0721 }
0722
0723 void CSCTriggerPrimitivesReader::bookLCTTMBHistos() {
0724   string s;
0725
0726   edm::Service<TFileService> fs;
0727   hLctTMBPerEvent = fs->make<TH1F>("LCTs_per_event", "LCTs per event", 31, -0.5, 30.5);
0728   hLctTMBPerChamber = fs->make<TH1F>("LCTs_per_chamber", "LCTs per chamber", 3, -0.5, 2.5);
0729   hLctTMBPerCSC = fs->make<TH1F>("LCTs_per_CSCtype", "LCTs per CSC type", 10, -0.5, 9.5);
0730   hCorrLctTMBPerCSC = fs->make<TH1F>("CorrLCTs_per_CSCtype", "Corr. LCTs per CSC type", 10, -0.5, 9.5);
0731   hLctTMBEndcap = fs->make<TH1F>("LCTS_endcap", "Endcap", 4, -0.5, 3.5);
0732   hLctTMBStation = fs->make<TH1F>("LCTS_station", "Station", 6, -0.5, 5.5);
0733   hLctTMBSector = fs->make<TH1F>("LCTS_sector", "Sector", 8, -0.5, 7.5);
0734   hLctTMBRing = fs->make<TH1F>("LCTS_ring", "Ring", 5, -0.5, 4.5);
0735   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0736     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0737       char asdf[256];
0738       float csc_max = NCHAMBERS[j] + 0.5;
0739       if (i == 0)
0740         s = "LCTs, " + csc_type_plus[j];
0741       else
0742         s = "LCTs, " + csc_type_minus[j];
0743       sprintf(asdf, "LCTs_%i", i * CSC_TYPES + j);
0744       hLctTMBCsc[i][j] = fs->make<TH1F>(s.c_str(), s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0745     }
0746   }
0747
0748   hLctTMBValid = fs->make<TH1F>("LCT_validity", "LCT validity", 3, -0.5, 2.5);
0749   hLctTMBQuality = fs->make<TH1F>("LCT_quality", "LCT quality", 17, -0.5, 16.5);
0750   hLctTMBKeyGroup = fs->make<TH1F>("LCT_key_wiregroup", "LCT key wiregroup", 120, -0.5, 119.5);
0751   hLctTMBKeyStrip = fs->make<TH1F>("LCT_key_strip", "LCT key strip", 160, -0.5, 159.5);
0752   hLctTMBStripType = fs->make<TH1F>("LCT_strip_type", "LCT strip type", 3, -0.5, 2.5);
0753   hLctTMBPattern = fs->make<TH1F>("LCT_pattern", "LCT pattern", 13, -0.5, 12.5);
0754   hLctTMBBend = fs->make<TH1F>("LCT_bend", "LCT L/R bend", 3, -0.5, 2.5);
0755   hLctTMBBXN = fs->make<TH1F>("LCT_bx", "LCT bx", 20, -0.5, 19.5);
0756
0757   // LCT quantities per station
0758   char histname[60];
0759   for (int istat = 0; istat < MAX_STATIONS; istat++) {
0760     sprintf(histname, "LCT_CSCId, station %d", istat + 1);
0761     hLctTMBChamber[istat] = fs->make<TH1F>("", histname, 10, -0.5, 9.5);
0762   }
0763
0764   hLctTMBKeyGroupME11 = fs->make<TH1F>("hLctTMBKeyGroupME11", "LCT key wiregroup ME1/1", 50, -0.5, 49.5);
0765   hLctTMBKeyStripME11 = fs->make<TH1F>("hLctTMBKeyStripME11", "LCT key strip ME1/1", 161, -0.5, 160.5);
0766
0767   bookedLCTTMBHistos = true;
0768   std::cout << " bookedLCTTMBHistos " << std::endl;
0769 }
0770
0771 int CSCTriggerPrimitivesReader::chamberIXi(CSCDetId id) {
0772   //    1/1 1/2 1/3 2/1 2/2 3/1 3/2 4/1 4/2 -1/1 -1/2 -1/3 -2/1 -2/2 -3/1 -3/2 -4/1 -4/2
0773   //ix= 0   1   2   3   4   5   6   7   8    9    10   11   12   13   14   15   16   17
0774   int ix = 0;
0775   if (id.station() != 1) {
0776     ix = (id.station() - 2) * 2 + 3;
0777   }
0778   ix += id.ring() - 1;
0779   if (id.endcap() == 2) {
0780     ix += 9;
0781   }
0782   return ix;
0783 }
0784
0785 int CSCTriggerPrimitivesReader::chamberIX(CSCDetId id) {
0786   int ix = 1;
0787   if (id.station() != 1) {
0788     ix = (id.station() - 2) * 2 + 3 + 1;
0789   }
0790   ix += id.ring() - 1;
0791   if (id.endcap() == 2) {
0792     ix *= -1;
0793   }
0794   return ix;
0795 }
0796
0797 int CSCTriggerPrimitivesReader::chamberSerial(CSCDetId id) {
0798   int st = id.station();
0799   int ri = id.ring();
0800   int ch = id.chamber();
0801   int ec = id.endcap();
0802   int kSerial = ch;
0803   if (st == 1 && ri == 1)
0804     kSerial = ch;
0805   if (st == 1 && ri == 2)
0806     kSerial = ch + 36;
0807   if (st == 1 && ri == 3)
0808     kSerial = ch + 72;
0809   if (st == 1 && ri == 4)
0810     kSerial = ch;
0811   if (st == 2 && ri == 1)
0812     kSerial = ch + 108;
0813   if (st == 2 && ri == 2)
0814     kSerial = ch + 126;
0815   if (st == 3 && ri == 1)
0816     kSerial = ch + 162;
0817   if (st == 3 && ri == 2)
0818     kSerial = ch + 180;
0819   if (st == 4 && ri == 1)
0820     kSerial = ch + 216;
0821   if (st == 4 && ri == 2)
0822     kSerial = ch + 234;  // one day...
0823   if (ec == 2)
0824     kSerial = kSerial + 300;
0825   return kSerial;
0826 }
0827
0828 void CSCTriggerPrimitivesReader::bookLCTMPCHistos() {
0829   edm::Service<TFileService> fs;
0830   hLctMPCPerEvent = fs->make<TH1F>("MPC_per_event", "LCTs per event", 31, -0.5, 30.5);
0831   hLctMPCPerCSC = fs->make<TH1F>("MPC_per_CSCtype", "LCTs per CSC type", 10, -0.5, 9.5);
0832   hCorrLctMPCPerCSC = fs->make<TH1F>("CorrMPC_per_CSCtype", "Corr. LCTs per CSC type", 10, -0.5, 9.5);
0833   hLctMPCEndcap = fs->make<TH1F>("MPC_Endcap", "Endcap", 4, -0.5, 3.5);
0834   hLctMPCStation = fs->make<TH1F>("MPC_Station", "Station", 6, -0.5, 5.5);
0835   hLctMPCSector = fs->make<TH1F>("MPC_Sector", "Sector", 8, -0.5, 7.5);
0836   hLctMPCRing = fs->make<TH1F>("MPC_Ring", "Ring", 5, -0.5, 4.5);
0837
0838   hLctMPCValid = fs->make<TH1F>("MPC_validity", "LCT validity", 3, -0.5, 2.5);
0839   hLctMPCQuality = fs->make<TH1F>("MPC_quality", "LCT quality", 17, -0.5, 16.5);
0840   hLctMPCKeyGroup = fs->make<TH1F>("MPC_key_wiregroup", "LCT key wiregroup", 120, -0.5, 119.5);
0841   hLctMPCKeyStrip = fs->make<TH1F>("MPC_key_strip", "LCT key strip", 160, -0.5, 159.5);
0842   hLctMPCStripType = fs->make<TH1F>("MPC_strip_type", "LCT strip type", 3, -0.5, 2.5);
0843   hLctMPCPattern = fs->make<TH1F>("MPC_pattern", "LCT pattern", 13, -0.5, 12.5);
0844   hLctMPCBend = fs->make<TH1F>("MPC_bend", "LCT L/R bend", 3, -0.5, 2.5);
0845   hLctMPCBXN = fs->make<TH1F>("MPC_bx", "LCT bx", 20, -0.5, 19.5);
0846
0847   // LCT quantities per station
0848   char histname[60];
0849   for (int istat = 0; istat < MAX_STATIONS; istat++) {
0850     sprintf(histname, "MPC_CSCId, station %d", istat + 1);
0851     hLctMPCChamber[istat] = fs->make<TH1F>("", histname, 10, -0.5, 9.5);
0852   }
0853
0854   hLctMPCKeyGroupME11 = fs->make<TH1F>("hLctMPCKeyGroupME11", "MPC LCT key wiregroup ME1/1", 50, -0.5, 49.5);
0855   hLctMPCKeyStripME11 = fs->make<TH1F>("hLctMPCKeyStripME11", "MPC LCT key strip ME1/1", 161, -0.5, 160.5);
0856
0857   bookedLCTMPCHistos = true;
0858   if (debug)
0859     LogTrace("CSCTriggerPrimitivesReader") << "bookedLCTMPCHistos  ";
0860   std::cout << " bookedLCTMPCHistos " << std::endl;
0861 }
0862
0863 void CSCTriggerPrimitivesReader::bookCompHistos() {
0864   string s;
0865
0866   edm::Service<TFileService> fs;
0867   //  hAlctCompMatch;
0868
0869   hAlctCompFound = fs->make<TH1F>("h_ALCT_found", "h_ALCT_found", 600, 0.5, 600.5);
0870   hAlctCompFound2 = fs->make<TH2F>("h_ALCT_found2", "h_ALCT_found2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0871   hAlctCompFound2x = fs->make<TH2F>("h_ALCT_found2x", "h_ALCT_found2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0872
0873   hAlctCompSameN = fs->make<TH1F>("h_ALCT_SameN", "h_ALCT_SameN", 600, 0.5, 600.5);
0874   hAlctCompSameN2 = fs->make<TH2F>("h_ALCT_SameN2", "h_ALCT_SameN2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0875   hAlctCompSameN2x = fs->make<TH2F>("h_ALCT_SameN2x", "h_ALCT_SameN2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0876
0877   hAlctCompMatch = fs->make<TH1F>("h_ALCT_match", "h_ALCT_match", 600, 0.5, 600.5);
0878   hAlctCompMatch2 = fs->make<TH2F>("h_ALCT_match2", "h_ALCT_match2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0879   hAlctCompMatch2x = fs->make<TH2F>("h_ALCT_match2x", "h_ALCT_match2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0880
0881   hAlctCompTotal = fs->make<TH1F>("h_ALCT_total", "h_ALCT_total", 600, 0.5, 600.5);
0882   hAlctCompTotal2 = fs->make<TH2F>("h_ALCT_total2", "h_ALCT_total2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0883   hAlctCompTotal2x = fs->make<TH2F>("h_ALCT_total2x", "h_ALCT_total2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0884
0885   hClctCompFound = fs->make<TH1F>("h_CLCT_found", "h_CLCT_found", 600, 0.5, 600.5);
0886   hClctCompFound2 = fs->make<TH2F>("h_CLCT_found2", "h_CLCT_found2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0887   hClctCompFound2x = fs->make<TH2F>("h_CLCT_found2x", "h_CLCT_found2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0888
0889   hClctCompSameN = fs->make<TH1F>("h_CLCT_SameN", "h_CLCT_SameN", 600, 0.5, 600.5);
0890   hClctCompSameN2 = fs->make<TH2F>("h_CLCT_SameN2", "h_CLCT_SameN2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0891   hClctCompSameN2x = fs->make<TH2F>("h_CLCT_SameN2x", "h_CLCT_SameN2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0892
0893   hClctCompMatch = fs->make<TH1F>("h_CLCT_match", "h_CLCT_match", 600, 0.5, 600.5);
0894   hClctCompMatch2 = fs->make<TH2F>("h_CLCT_match2", "h_CLCT_match2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0895   hClctCompMatch2x = fs->make<TH2F>("h_CLCT_match2x", "h_CLCT_match2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0896
0897   hClctCompTotal = fs->make<TH1F>("h_CLCT_total", "h_CLCT_total", 600, 0.5, 600.5);
0898   hClctCompTotal2 = fs->make<TH2F>("h_CLCT_total2", "h_CLCT_total2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0899   hClctCompTotal2x = fs->make<TH2F>("h_CLCT_total2x", "h_CLCT_total2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0900
0901   hLCTCompFound = fs->make<TH1F>("h_LCT_found", "h_LCT_found", 600, 0.5, 600.5);
0902   hLCTCompFound2 = fs->make<TH2F>("h_LCT_found2", "h_LCT_found2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0903   hLCTCompFound2x = fs->make<TH2F>("h_LCT_found2x", "h_LCT_found2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0904
0905   hLCTCompSameN = fs->make<TH1F>("h_LCT_SameN", "h_LCT_SameN", 600, 0.5, 600.5);
0906   hLCTCompSameN2 = fs->make<TH2F>("h_LCT_SameN2", "h_LCT_SameN2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0907   hLCTCompSameN2x = fs->make<TH2F>("h_LCT_SameN2x", "h_LCT_SameN2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0908
0909   hLCTCompMatch = fs->make<TH1F>("h_LCT_match", "h_LCT_match", 600, 0.5, 600.5);
0910   hLCTCompMatch2 = fs->make<TH2F>("h_LCT_match2", "h_LCT_match2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0911   hLCTCompMatch2x = fs->make<TH2F>("h_LCT_match2x", "h_LCT_match2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0912
0913   hLCTCompTotal = fs->make<TH1F>("h_LCT_total", "h_LCT_total", 600, 0.5, 600.5);
0914   hLCTCompTotal2 = fs->make<TH2F>("h_LCT_total2", "h_LCT_total2", 19, -9.5, 9.5, 36, 0.5, 36.5);
0915   hLCTCompTotal2x = fs->make<TH2F>("h_LCT_total2x", "h_LCT_total2x", 19, -9.5, 9.5, 18, 0.5, 36.5);
0916
0917   //Chad's improved historgrams
0918   hAlctCompFound2i = fs->make<TH2F>("h_ALCT_found2i", "h_ALCT_found2i", 18, 0, 18, 36, 0.5, 36.5);
0919   hAlctCompSameN2i = fs->make<TH2F>("h_ALCT_SameN2i", "h_ALCT_SameN2i", 18, 0, 18, 36, 0.5, 36.5);
0920   hAlctCompMatch2i = fs->make<TH2F>("h_ALCT_match2i", "h_ALCT_match2i", 18, 0, 18, 36, 0.5, 36.5);
0921   hAlctCompTotal2i = fs->make<TH2F>("h_ALCT_total2i", "h_ALCT_total2i", 18, 0, 18, 36, 0.5, 36.5);
0922   hClctCompFound2i = fs->make<TH2F>("h_CLCT_found2i", "h_CLCT_found2i", 18, 0, 18, 36, 0.5, 36.5);
0923   hClctCompSameN2i = fs->make<TH2F>("h_CLCT_SameN2i", "h_CLCT_SameN2i", 18, 0, 18, 36, 0.5, 36.5);
0924   hClctCompMatch2i = fs->make<TH2F>("h_CLCT_match2i", "h_CLCT_match2i", 18, 0, 18, 36, 0.5, 36.5);
0925   hClctCompTotal2i = fs->make<TH2F>("h_CLCT_total2i", "h_CLCT_total2i", 18, 0, 18, 36, 0.5, 36.5);
0926   hLCTCompFound2i = fs->make<TH2F>("h_LCT_found2i", "h_LCT_found2i", 18, 0, 18, 36, 0.5, 36.5);
0927   hLCTCompSameN2i = fs->make<TH2F>("h_LCT_SameN2i", "h_LCT_SameN2i", 18, 0, 18, 36, 0.5, 36.5);
0928   hLCTCompMatch2i = fs->make<TH2F>("h_LCT_match2i", "h_LCT_match2i", 18, 0, 18, 36, 0.5, 36.5);
0929   hLCTCompTotal2i = fs->make<TH2F>("h_LCT_total2i", "h_LCT_total2i", 18, 0, 18, 36, 0.5, 36.5);
0930
0931   //  hAlctCompFound = fs->make<TH1F>("h_ALCT_found","h_ALCT_found",600,0.5,600.5);
0932
0933   // ALCTs.
0934   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0935     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0936       char asdf[256];
0937
0938       float csc_max = NCHAMBERS[j] + 0.5;
0939       if (i == 0)
0940         s = "Comp_ALCTs, " + csc_type_plus[j];
0941       else
0942         s = "Comp_ALCTs, " + csc_type_minus[j];
0943       sprintf(asdf, "Comp_ALCTsFound_%i", i * CSC_TYPES + j);
0944       hAlctCompFoundCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0945       sprintf(asdf, "Comp_ALCTsSame_%i", i * CSC_TYPES + j);
0946       hAlctCompSameNCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0947       sprintf(asdf, "Comp_ALCTsTotal_%i", i * CSC_TYPES + j);
0948       hAlctCompTotalCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0949       sprintf(asdf, "Comp_ALCTsMatch_%i", i * CSC_TYPES + j);
0950       hAlctCompMatchCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0951       hAlctCompFoundCsc[i][j]->Sumw2();
0952       hAlctCompSameNCsc[i][j]->Sumw2();
0953       hAlctCompTotalCsc[i][j]->Sumw2();
0954       hAlctCompMatchCsc[i][j]->Sumw2();
0955     }
0956   }
0957
0958   // CLCTs.
0959   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0960     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0961       float csc_max = NCHAMBERS[j] + 0.5;
0962       char asdf[256];
0963       if (i == 0)
0964         s = "Comp_CLCTs, " + csc_type_plus[j];
0965       else
0966         s = "Comp_CLCTs, " + csc_type_minus[j];
0967       sprintf(asdf, "Comp_CLCTsFound_%i", i * CSC_TYPES + j);
0968       hClctCompFoundCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0969       sprintf(asdf, "Comp_CLCTsSame_%i", i * CSC_TYPES + j);
0970       hClctCompSameNCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0971       sprintf(asdf, "Comp_CLCTsTotal_%i", i * CSC_TYPES + j);
0972       hClctCompTotalCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0973       sprintf(asdf, "Comp_CLCTsMatch_%i", i * CSC_TYPES + j);
0974       hClctCompMatchCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0975       hClctCompFoundCsc[i][j]->Sumw2();
0976       hClctCompSameNCsc[i][j]->Sumw2();
0977       hClctCompTotalCsc[i][j]->Sumw2();
0978       hClctCompMatchCsc[i][j]->Sumw2();
0979     }
0980   }
0981
0982   // Correlated LCTs.
0983   for (int i = 0; i < MAX_ENDCAPS; i++) {  // endcaps
0984     for (int j = 0; j < CSC_TYPES; j++) {  // station/ring
0985       float csc_max = NCHAMBERS[j] + 0.5;
0986       char asdf[256];
0987       if (i == 0)
0988         s = "LCTs, " + csc_type_plus[j];
0989       else
0990         s = "LCTs, " + csc_type_minus[j];
0991       sprintf(asdf, "LCTs_CompFound_%i", i * CSC_TYPES + j);
0992       hLctCompFoundCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0993       sprintf(asdf, "LCTs_CompSame_%i", i * CSC_TYPES + j);
0994       hLctCompSameNCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0995       sprintf(asdf, "LCTs_CompTotal_%i", i * CSC_TYPES + j);
0996       hLctCompTotalCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0997       sprintf(asdf, "LCTs_CompMatch_%i", i * CSC_TYPES + j);
0998       hLctCompMatchCsc[i][j] = fs->make<TH1F>(asdf, s.c_str(), NCHAMBERS[j], 0.5, csc_max);
0999       hLctCompFoundCsc[i][j]->Sumw2();
1000       hLctCompSameNCsc[i][j]->Sumw2();
1001       hLctCompTotalCsc[i][j]->Sumw2();
1002       hLctCompMatchCsc[i][j]->Sumw2();
1003     }
1004   }
1005
1006   bookedCompHistos = true;
1007   if (debug)
1008     LogTrace("CSCTriggerPrimitivesReader") << "bookedCompHistos  ";
1009   std::cout << "bookCompHistos " << std::endl;
1010 }
1011
1012 void CSCTriggerPrimitivesReader::bookResolHistos() {
1013   edm::Service<TFileService> fs;
1014
1015   // Limits for resolution histograms
1016   const double EDMIN = -0.05;  // eta min
1017   const double EDMAX = 0.05;   // eta max
1018   const double PDMIN = -5.0;   // phi min (mrad)
1019   const double PDMAX = 5.0;    // phi max (mrad)
1020
1021   hResolDeltaWG = fs->make<TH1F>("", "Delta key wiregroup", 10, -5., 5.);
1022
1023   hResolDeltaHS = fs->make<TH1F>("", "Delta key halfstrip", 10, -5., 5.);
1024   hResolDeltaDS = fs->make<TH1F>("", "Delta key distrip", 10, -5., 5.);
1025
1026   hResolDeltaEta = fs->make<TH1F>("", "#eta_rec-#eta_sim", 100, EDMIN, EDMAX);
1027   hResolDeltaPhi = fs->make<TH1F>("", "#phi_rec-#phi_sim (mrad)", 100, -10., 10.);
1028   hResolDeltaPhiHS = fs->make<TH1F>("", "#phi_rec-#phi_sim (mrad), halfstrips", 100, -10., 10.);
1029   hResolDeltaPhiDS = fs->make<TH1F>("", "#phi_rec-#phi_sim (mrad), distrips", 100, -10., 10.);
1030
1031   hEtaRecVsSim = fs->make<TH2F>("", "#eta_rec vs #eta_sim", 64, 0.9, 2.5, 64, 0.9, 2.5);
1032   hPhiRecVsSim = fs->make<TH2F>("", "#phi_rec vs #phi_sim", 100, 0., TWOPI, 100, 0., TWOPI);
1033
1034   // LCT quantities per station
1035   char histname[60];
1036   for (int i = 0; i < MAX_STATIONS; i++) {
1037     sprintf(histname, "ALCTs vs eta, station %d", i + 1);
1038     hAlctVsEta[i] = fs->make<TH1F>("", histname, 66, 0.875, 2.525);
1039
1040     sprintf(histname, "CLCTs vs phi, station %d", i + 1);
1041     hClctVsPhi[i] = fs->make<TH1F>("", histname, 100, 0., TWOPI);
1042
1043     sprintf(histname, "#LT#eta_rec-#eta_sim#GT, station %d", i + 1);
1044     hEtaDiffVsEta[i] = fs->make<TH1F>("", histname, 66, 0.875, 2.525);
1045
1046     sprintf(histname, "#LT#phi_rec-#phi_sim#GT, station %d", i + 1);
1047     hPhiDiffVsPhi[i] = fs->make<TH1F>("", histname, 100, 0., TWOPI);
1048   }
1049
1050   for (int i = 0; i < CSC_TYPES; i++) {
1051     string t0 = "#eta_rec-#eta_sim, " + csc_type[i];
1052     hEtaDiffCsc[i][0] = fs->make<TH1F>("", t0.c_str(), 100, EDMIN, EDMAX);
1053     string t1 = t0 + ", endcap1";
1054     hEtaDiffCsc[i][1] = fs->make<TH1F>("", t1.c_str(), 100, EDMIN, EDMAX);
1055     string t2 = t0 + ", endcap2";
1056     hEtaDiffCsc[i][2] = fs->make<TH1F>("", t2.c_str(), 100, EDMIN, EDMAX);
1057
1058     string t4 = "#eta_rec-#eta_sim vs wiregroup, " + csc_type[i];
1059     hEtaDiffVsWireCsc[i] = fs->make<TH2F>("", t4.c_str(), MAX_WG[i], 0., MAX_WG[i], 100, EDMIN, EDMAX);
1060
1061     string u0 = "#phi_rec-#phi_sim, " + csc_type[i];
1062     hPhiDiffCsc[i][0] = fs->make<TH1F>("", u0.c_str(), 100, PDMIN, PDMAX);
1063     string u1 = u0 + ", endcap1";
1064     hPhiDiffCsc[i][1] = fs->make<TH1F>("", u1.c_str(), 100, PDMIN, PDMAX);
1065     string u2 = u0 + ", endcap2";
1066     hPhiDiffCsc[i][2] = fs->make<TH1F>("", u2.c_str(), 100, PDMIN, PDMAX);
1067     hPhiDiffCsc[i][3] = fs->make<TH1F>("", u0.c_str(), 100, PDMIN, PDMAX);
1068     hPhiDiffCsc[i][4] = fs->make<TH1F>("", u0.c_str(), 100, PDMIN, PDMAX);
1069
1070     int MAX_DS = MAX_HS[i] / 4;
1071     string u5 = "#phi_rec-#phi_sim (mrad) vs distrip, " + csc_type[i];
1072     hPhiDiffVsStripCsc[i][0] = fs->make<TH2F>("", u5.c_str(), MAX_DS, 0., MAX_DS, 100, PDMIN, PDMAX);
1073     string u6 = "#phi_rec-#phi_sim (mrad) vs halfstrip, " + csc_type[i];
1074     hPhiDiffVsStripCsc[i][1] = fs->make<TH2F>("", u6.c_str(), MAX_HS[i], 0., MAX_HS[i], 100, PDMIN, PDMAX);
1075
1076     string u7 = "#phi(layer 1)-#phi(layer 6), mrad, " + csc_type[i];
1077     hTrueBendCsc[i] = fs->make<TH1F>("", u7.c_str(), 100, -10., 10.);
1078   }
1079
1080   int max_patterns, phibend;
1081   max_patterns = CSCConstants::NUM_CLCT_PATTERNS;
1082   for (int i = 0; i < max_patterns; i++) {
1083     phibend = ptype_TMB07[i];
1084     sprintf(histname, "#phi_rec-#phi_sim, bend = %d", phibend);
1085     hPhiDiffPattern[i] = fs->make<TH1F>("", histname, 100, PDMIN, PDMAX);
1086   }
1087
1088   bookedResolHistos = true;
1089   if (debug)
1090     LogTrace("CSCTriggerPrimitivesReader") << "bookedResolHistos  ";
1091 }
1092
1093 void CSCTriggerPrimitivesReader::bookEfficHistos() {
1094   edm::Service<TFileService> fs;
1095
1096   // Efficiencies per station.
1097   char histname[60];
1098   for (int i = 0; i < MAX_STATIONS; i++) {
1099     sprintf(histname, "SimHits vs eta, station %d", i + 1);
1100     hEfficHitsEta[i] = fs->make<TH1F>("", histname, 66, 0.875, 2.525);
1101
1102     sprintf(histname, "ALCTs vs eta, station %d", i + 1);
1103     hEfficALCTEta[i] = fs->make<TH1F>("", histname, 66, 0.875, 2.525);
1104
1105     sprintf(histname, "CLCTs vs eta, station %d", i + 1);
1106     hEfficCLCTEta[i] = fs->make<TH1F>("", histname, 66, 0.875, 2.525);
1107   }
1108
1109   // Efficiencies per chamber type.
1110   for (int i = 0; i < CSC_TYPES; i++) {
1111     string t0 = "SimHits vs eta, " + csc_type[i];
1112     hEfficHitsEtaCsc[i] = fs->make<TH1F>("", t0.c_str(), 66, 0.875, 2.525);
1113     string t1 = "ALCTs vs eta, " + csc_type[i];
1114     hEfficALCTEtaCsc[i] = fs->make<TH1F>("", t1.c_str(), 66, 0.875, 2.525);
1115     string t2 = "CLCTs vs eta, " + csc_type[i];
1116     hEfficCLCTEtaCsc[i] = fs->make<TH1F>("", t1.c_str(), 66, 0.875, 2.525);
1117   }
1118
1119   bookedEfficHistos = true;
1120   if (debug)
1121     LogTrace("CSCTriggerPrimitivesReader") << "bookedEfficHistos  ";
1122   std::cout << " bookedEfficHistos " << std::endl;
1123 }
1124
1125 void CSCTriggerPrimitivesReader::fillALCTHistos(const CSCALCTDigiCollection* alcts) {
1126   // Book histos when called for the first time.
1127   if (!bookedALCTHistos)
1128     bookALCTHistos();
1129
1130   int nValidALCTs = 0;
1131   for (auto detUnitIt = alcts->begin(); detUnitIt != alcts->end(); detUnitIt++) {
1132     int nValidALCTsPerCSC = 0;
1133     const CSCDetId& id = (*detUnitIt).first;
1134     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
1135       continue;
1136     const auto& range = (*detUnitIt).second;
1137     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
1138       bool alct_valid = (*digiIt).isValid();
1139       hAlctValid->Fill(alct_valid);
1140       if (alct_valid) {
1141         hAlctQuality->Fill((*digiIt).getQuality());
1142         hAlctAccel->Fill((*digiIt).getAccelerator());
1143         hAlctCollis->Fill((*digiIt).getCollisionB());
1144         hAlctKeyGroup->Fill((*digiIt).getKeyWG());
1145         hAlctBXN->Fill((*digiIt).getBX());
1146
1147         int csctype = getCSCType(id);
1148         hAlctPerCSC->Fill(csctype);
1149         hAlctCsc[id.endcap() - 1][csctype]->Fill(id.chamber());
1150
1151         if (csctype == 0)
1152           hAlctKeyGroupME11->Fill((*digiIt).getKeyWG());
1153
1154         nValidALCTs++;
1155         nValidALCTsPerCSC++;
1156
1157         if (debug)
1158           LogTrace("CSCTriggerPrimitivesReader")
1159               << (*digiIt) << " found in ME" << ((id.endcap() == 1) ? "+" : "-") << id.station() << "/" << id.ring()
1160               << "/" << id.chamber() << " (sector " << id.triggerSector() << " trig id. " << id.triggerCscId() << ")";
1161         //cout << "raw id = " << id.rawId() << endl;
1162       }
1163     }
1164     hAlctPerChamber->Fill(nValidALCTsPerCSC);
1165   }
1166   hAlctPerEvent->Fill(nValidALCTs);
1167   if (debug)
1168     LogTrace("CSCTriggerPrimitivesReader") << nValidALCTs << " valid ALCTs found in this event";
1169   numALCT += nValidALCTs;
1170 }
1171
1172 void CSCTriggerPrimitivesReader::fillCLCTHistos(const CSCCLCTDigiCollection* clcts) {
1173   // Book histos when called for the first time.
1174   if (!bookedCLCTHistos)
1175     bookCLCTHistos();
1176
1177   int nValidCLCTs = 0;
1178   for (auto detUnitIt = clcts->begin(); detUnitIt != clcts->end(); detUnitIt++) {
1179     int nValidCLCTsPerCSC = 0;
1180     const CSCDetId& id = (*detUnitIt).first;
1181     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
1182       continue;
1183     const auto& range = (*detUnitIt).second;
1184     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
1185       bool clct_valid = (*digiIt).isValid();
1186       hClctValid->Fill(clct_valid);
1187       if (clct_valid) {
1188         int striptype = (*digiIt).getStripType();
1189         int keystrip = (*digiIt).getKeyStrip();  // halfstrip #
1190         if (striptype == 0)
1191           keystrip /= 4;  // distrip # for distrip ptns
1192         hClctQuality->Fill((*digiIt).getQuality());
1193         hClctStripType->Fill(striptype);
1194         hClctSign->Fill((*digiIt).getBend());
1195         hClctCFEB->Fill((*digiIt).getCFEB());
1196         hClctBXN->Fill((*digiIt).getBX());
1197         hClctKeyStrip[striptype]->Fill(keystrip);
1198         hClctPattern[striptype]->Fill((*digiIt).getPattern());
1199
1200         int csctype = getCSCType(id);
1201         hClctPerCSC->Fill(csctype);
1202         hClctCsc[id.endcap() - 1][csctype]->Fill(id.chamber());
1203
1204         if (striptype == 1 && csctype == 0)
1205           hClctKeyStripME11->Fill(keystrip);
1206
1207         int phibend;
1208         int pattern = (*digiIt).getPattern();
1209         phibend = ptype_TMB07[pattern];
1210         hClctBendCsc[csctype][striptype]->Fill(abs(phibend));
1211
1212         if (striptype == 0)  // distrips
1213           hClctKeyStripCsc[csctype]->Fill(keystrip);
1214
1215         nValidCLCTs++;
1216         nValidCLCTsPerCSC++;
1217
1218         if (debug)
1219           LogTrace("CSCTriggerPrimitivesReader")
1220               << (*digiIt) << " found in ME" << ((id.endcap() == 1) ? "+" : "-") << id.station() << "/" << id.ring()
1221               << "/" << id.chamber() << " (sector " << id.triggerSector() << " trig id. " << id.triggerCscId() << ")";
1222       }
1223     }
1224     hClctPerChamber->Fill(nValidCLCTsPerCSC);
1225   }
1226   hClctPerEvent->Fill(nValidCLCTs);
1227   if (debug)
1228     LogTrace("CSCTriggerPrimitivesReader") << nValidCLCTs << " valid CLCTs found in this event";
1229   numCLCT += nValidCLCTs;
1230 }
1231
1232 void CSCTriggerPrimitivesReader::fillLCTTMBHistos(const CSCCorrelatedLCTDigiCollection* lcts) {
1233   // Book histos when called for the first time.
1234   if (!bookedLCTTMBHistos)
1235     bookLCTTMBHistos();
1236
1237   int nValidLCTs = 0;
1238   bool alct_valid, clct_valid;
1239
1240   for (auto detUnitIt = lcts->begin(); detUnitIt != lcts->end(); detUnitIt++) {
1241     int nValidLCTsPerCSC = 0;
1242     const CSCDetId& id = (*detUnitIt).first;
1243     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
1244       continue;
1245     const auto& range = (*detUnitIt).second;
1246     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
1247       bool lct_valid = (*digiIt).isValid();
1248       hLctTMBValid->Fill(lct_valid);
1249       if (lct_valid) {
1250         hLctTMBEndcap->Fill(id.endcap());
1251         hLctTMBStation->Fill(id.station());
1252         hLctTMBSector->Fill(id.triggerSector());
1253         hLctTMBRing->Fill(id.ring());
1254         hLctTMBChamber[id.station() - 1]->Fill(id.triggerCscId());
1255
1256         int quality = (*digiIt).getQuality();
1257         hLctTMBQuality->Fill(quality);
1258         hLctTMBBXN->Fill((*digiIt).getBX());
1259
1260         alct_valid = (quality != 0 && quality != 2);
1261         if (alct_valid) {
1262           hLctTMBKeyGroup->Fill((*digiIt).getKeyWG());
1263         }
1264
1265         clct_valid = (quality != 0 && quality != 1);
1266         if (clct_valid) {
1267           hLctTMBKeyStrip->Fill((*digiIt).getStrip());
1268           hLctTMBStripType->Fill(1.);
1269           hLctTMBPattern->Fill((*digiIt).getPattern());
1270           hLctTMBBend->Fill((*digiIt).getBend());
1271         }
1272
1273         int csctype = getCSCType(id);
1274         hLctTMBPerCSC->Fill(csctype);
1275         hLctTMBCsc[id.endcap() - 1][csctype]->Fill(id.chamber());
1276         // Truly correlated LCTs; for DAQ
1277         if (alct_valid && clct_valid)
1278           hCorrLctTMBPerCSC->Fill(csctype);
1279
1280         if (alct_valid && csctype == 0) {
1281           hLctTMBKeyGroupME11->Fill((*digiIt).getKeyWG());
1282         }
1283         if (clct_valid && csctype == 0) {
1284           hLctTMBKeyStripME11->Fill((*digiIt).getStrip());
1285         }
1286
1287         nValidLCTs++;
1288         nValidLCTsPerCSC++;
1289
1290         if (debug)
1291           LogTrace("CSCTriggerPrimitivesReader")
1292               << (*digiIt) << " found in ME" << ((id.endcap() == 1) ? "+" : "-") << id.station() << "/" << id.ring()
1293               << "/" << id.chamber() << " (sector " << id.triggerSector() << " trig id. " << id.triggerCscId() << ")";
1294       }
1295     }
1296     hLctTMBPerChamber->Fill(nValidLCTsPerCSC);
1297   }
1298   hLctTMBPerEvent->Fill(nValidLCTs);
1299   if (debug)
1300     LogTrace("CSCTriggerPrimitivesReader") << nValidLCTs << " valid LCTs found in this event";
1301   numLCTTMB += nValidLCTs;
1302 }
1303
1304 void CSCTriggerPrimitivesReader::fillLCTMPCHistos(const CSCCorrelatedLCTDigiCollection* lcts) {
1305   // Book histos when called for the first time.
1306   if (!bookedLCTMPCHistos)
1307     bookLCTMPCHistos();
1308
1309   int nValidLCTs = 0;
1310   bool alct_valid, clct_valid;
1311
1312   for (auto detUnitIt = lcts->begin(); detUnitIt != lcts->end(); detUnitIt++) {
1313     const CSCDetId& id = (*detUnitIt).first;
1314     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
1315       continue;
1316     const auto& range = (*detUnitIt).second;
1317     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
1318       bool lct_valid = (*digiIt).isValid();
1319       hLctMPCValid->Fill(lct_valid);
1320       if (lct_valid) {
1321         hLctMPCEndcap->Fill(id.endcap());
1322         hLctMPCStation->Fill(id.station());
1323         hLctMPCSector->Fill(id.triggerSector());
1324         hLctMPCRing->Fill(id.ring());
1325         hLctMPCChamber[id.station() - 1]->Fill(id.triggerCscId());
1326
1327         int quality = (*digiIt).getQuality();
1328         hLctMPCQuality->Fill(quality);
1329         hLctMPCBXN->Fill((*digiIt).getBX());
1330
1331         alct_valid = (quality != 0 && quality != 2);
1332         if (alct_valid) {
1333           hLctMPCKeyGroup->Fill((*digiIt).getKeyWG());
1334         }
1335
1336         clct_valid = (quality != 0 && quality != 1);
1337         if (clct_valid) {
1338           hLctMPCKeyStrip->Fill((*digiIt).getStrip());
1339           hLctMPCStripType->Fill(1.);
1340           hLctMPCPattern->Fill((*digiIt).getPattern());
1341           hLctMPCBend->Fill((*digiIt).getBend());
1342         }
1343
1344         int csctype = getCSCType(id);
1345         hLctMPCPerCSC->Fill(csctype);
1346         // Truly correlated LCTs; for DAQ
1347         if (alct_valid && clct_valid)
1348           hCorrLctMPCPerCSC->Fill(csctype);
1349
1350         if (alct_valid && csctype == 0) {
1351           hLctMPCKeyGroupME11->Fill((*digiIt).getKeyWG());
1352         }
1353         if (clct_valid && csctype == 0) {
1354           hLctMPCKeyStripME11->Fill((*digiIt).getStrip());
1355         }
1356
1357         nValidLCTs++;
1358
1359         if (debug)
1360           LogTrace("CSCTriggerPrimitivesReader")
1361               << "MPC " << (*digiIt) << " found in ME" << ((id.endcap() == 1) ? "+" : "-") << id.station() << "/"
1362               << id.ring() << "/" << id.chamber() << " (sector " << id.triggerSector() << " trig id. "
1363               << id.triggerCscId() << ")";
1364       }
1365     }
1366   }
1367   hLctMPCPerEvent->Fill(nValidLCTs);
1368   if (debug)
1369     LogTrace("CSCTriggerPrimitivesReader") << nValidLCTs << " MPC LCTs found in this event";
1370   numLCTMPC += nValidLCTs;
1371 }
1372
1373 void CSCTriggerPrimitivesReader::compare(const CSCALCTDigiCollection* alcts_data,
1374                                          const CSCALCTDigiCollection* alcts_emul,
1375                                          const CSCCLCTDigiCollection* clcts_data,
1376                                          const CSCCLCTDigiCollection* clcts_emul,
1377                                          const CSCCLCTPreTriggerDigiCollection* pretrigs_emul,
1378                                          const CSCCorrelatedLCTDigiCollection* lcts_data,
1379                                          const CSCCorrelatedLCTDigiCollection* lcts_emul) {
1380   // Book histos when called for the first time.
1381   if (!bookedCompHistos)
1382     bookCompHistos();
1383
1384   // Comparisons
1385   compareALCTs(alcts_data, alcts_emul);
1386   compareCLCTs(clcts_data, clcts_emul, pretrigs_emul);
1387   compareLCTs(lcts_data, lcts_emul, alcts_data, clcts_data);
1388   //compareMPCLCTs(mpclcts_data,  mpclcts_emul, alcts_data, clcts_data);
1389 }
1390
1391 void CSCTriggerPrimitivesReader::compareALCTs(const CSCALCTDigiCollection* alcts_data,
1392                                               const CSCALCTDigiCollection* alcts_emul) {
1393   int emul_corr_bx;
1394
1395   // Should be taken from config. parameters.
1396   int fifo_pretrig = 10;
1397   int fpga_latency = 6;
1398   int l1a_window_width = 7;
1399   // Time offset of raw hits w.r.t. the full 12-bit BXN.
1400   int rawhit_tbin_offset = (fifo_pretrig - fpga_latency) + (l1a_window_width - 1) / 2;
1401   // Extra difference due to additional register stages; determined
1402   // empirically.
1403   int register_delay = 2;
1404   // Loop over all chambers in search for ALCTs.
1405   std::vector<CSCALCTDigi>::const_iterator pd, pe;
1406   for (int endc = 1; endc <= 2; endc++) {
1407     for (int stat = 1; stat <= 4; stat++) {
1408       for (int ring = 1; ring <= maxRing(stat); ring++) {
1409         for (int cham = 1; cham <= 36; cham++) {
1410           // Calculate DetId.  0th layer means whole chamber.
1411           CSCDetId detid(endc, stat, ring, cham, 0);
1412
1413           // Skip chambers marked as bad.
1414           if (checkBadChambers_ && badChambers_->isInBadChamber(detid))
1415             continue;
1416           std::vector<CSCALCTDigi> alctV_data, alctV_emul;
1417           const auto& drange = alcts_data->get(detid);
1418           for (auto digiIt = drange.first; digiIt != drange.second; digiIt++) {
1419             if ((*digiIt).isValid()) {
1420               alctV_data.push_back(*digiIt);
1421             }
1422           }
1423
1424           const auto& erange = alcts_emul->get(detid);
1425           std::vector<bool> bookedalctV_emul;
1426           for (auto digiIt = erange.first; digiIt != erange.second; digiIt++) {
1427             if ((*digiIt).isValid()) {
1428               alctV_emul.push_back(*digiIt);
1429               bookedalctV_emul.push_back(false);
1430             }
1431           }
1432
1433           int ndata = alctV_data.size();
1434           int nemul = alctV_emul.size();
1435
1436           if (ndata == 0 && nemul == 0)
1437             continue;
1438           bool firstfill = true;
1439
1440           if (debug or (stat == 1 and ring == 1 and (ndata != nemul))) {
1441             ostringstream strstrm;
1442             strstrm << "\n--- ME" << ((detid.endcap() == 1) ? "+" : "-") << detid.station() << "/" << detid.ring()
1443                     << "/" << detid.chamber() << " (sector " << detid.triggerSector() << " trig id. "
1444                     << detid.triggerCscId() << "):\n";
1445             strstrm << "  **** " << ndata << " valid data ALCTs found:\n";
1446             for (pd = alctV_data.begin(); pd != alctV_data.end(); pd++) {
1447               strstrm << "     " << (*pd) << " Full BX = " << (*pd).getFullBX() << "\n";
1448             }
1449             strstrm << "  **** " << nemul << " valid emul ALCTs found:\n";
1450             for (pe = alctV_emul.begin(); pe != alctV_emul.end(); pe++) {
1451               strstrm << "     " << (*pe);
1452               for (pd = alctV_data.begin(); pd != alctV_data.end(); pd++) {
1453                 if ((*pd).getTrknmb() == (*pe).getTrknmb()) {
1454                   int emul_bx = (*pe).getBX();
1455                   emul_corr_bx = emul_bx - rawhit_tbin_offset + register_delay;
1456                   strstrm << " Corr BX = " << emul_corr_bx;
1457                   break;
1458                 }
1459               }
1460               strstrm << "\n";
1461             }
1462             LogTrace("CSCTriggerPrimitivesReader") << strstrm.str();
1463             //if (stat==1 && ring==1)
1464             //  std::cout <<"ME11  CompareALCTs "<< strstrm.str()<< std::endl;
1465           }
1466
1467           //Data, add HS quality later
1468           int perEv_nStub_data = 0;
1469           int perEv_nStub_emul = 0;
1470           perStub[0].init(RUN_, Event_);
1471           perStub[1].init(RUN_, Event_);
1472           for (pd = alctV_data.begin(); pd != alctV_data.end(); pd++) {
1473             perEv_nStub_data++;
1474           }
1475           perStub[0].t_endcap = endc;
1476           perStub[0].t_station = stat;
1477           perStub[0].t_chamber = cham;
1478           perStub[0].t_ring = ring;
1479           perStub[0].t_EventNumberAnalyzed = eventsAnalyzed;
1480           perStub[0].t_nStubs = ndata;
1481           perStub[0].t_nStubs_readout = ndata;
1482           event_tree[0]->Fill();
1483           //Emul
1484           for (pe = alctV_emul.begin(); pe != alctV_emul.end(); pe++) {
1485             perEv_nStub_emul++;
1486           }
1487           perStub[1].t_endcap = endc;
1488           perStub[1].t_station = stat;
1489           perStub[1].t_chamber = cham;
1490           perStub[1].t_ring = ring;
1491           perStub[1].t_EventNumberAnalyzed = eventsAnalyzed;
1492           perStub[1].t_nStubs = nemul;
1493           perStub[1].t_nStubs_readout = nemul;
1494           event_tree[1]->Fill();
1495
1496           int csctype = getCSCType(detid);
1497           hAlctCompFoundCsc[endc - 1][csctype]->Fill(cham);
1498           int mychamber = chamberSerial(detid);
1499           hAlctCompFound->Fill(mychamber);
1500           int ix = chamberIX(detid);
1501           int ix2 = chamberIXi(detid);
1502           //      printf("station %i, ring %i, chamber %i, ix+(detid.ring()-1) %i\n",
1503           //         detid.station(),detid.ring(),detid.chamber(),ix);
1504           if (detid.station() > 1 && detid.ring() == 1) {
1505             hAlctCompFound2x->Fill(ix, detid.chamber() * 2);
1506           } else {
1507             hAlctCompFound2->Fill(ix, detid.chamber());
1508           }
1509           hAlctCompFound2i->Fill(ix2, detid.chamber());
1510
1511           if (ndata != nemul) {
1512             LogTrace("CSCTriggerPrimitivesReader")
1513                 << "   +++ Different numbers of ALCTs found in ME" << ((endc == 1) ? "+" : "-") << stat << "/" << ring
1514                 << "/" << cham << ": data = " << ndata << " emulator = " << nemul << " +++\n";
1515           } else {
1516             hAlctCompSameNCsc[endc - 1][csctype]->Fill(cham);
1517             if (detid.station() > 1 && detid.ring() == 1) {
1518               hAlctCompSameN2x->Fill(ix, detid.chamber() * 2);
1519             } else {
1520               hAlctCompSameN2->Fill(ix, detid.chamber());
1521             }
1522             hAlctCompSameN2i->Fill(ix2, detid.chamber());
1523           }
1524
1525           for (int i = 0; i < ndata; i++) {
1526             if (alctV_data[i].isValid() == 0)
1527               continue;
1528             int data_trknmb = alctV_data[i].getTrknmb();
1529             //int data_quality   = alctV_data[i].getQuality();
1530             //int data_accel     = alctV_data[i].getAccelerator();
1531             //int data_collB     = alctV_data[i].getCollisionB();
1532             int data_wiregroup = alctV_data[i].getKeyWG();
1533             int data_bx = alctV_data[i].getBX();
1534
1535             stubs_comparison[0].init(RUN_, Event_);
1536             stubs_comparison[0].firstfill = firstfill;
1537             if (firstfill)
1538               firstfill = false;
1539             stubs_comparison[0].nEvents = eventsAnalyzed;
1540             stubs_comparison[0].endcap = endc;
1541             stubs_comparison[0].station = stat;
1542             stubs_comparison[0].ring = ring;
1543             stubs_comparison[0].chamber = cham;
1544             stubs_comparison[0].chambertype = detid.iChamberType();
1545             stubs_comparison[0].totStubs_data = ndata;
1546             stubs_comparison[0].totStubs_emul = nemul;
1547             stubs_comparison[0].totStubs_emul_readout = nemul;
1548             stubs_comparison[0].nStub_data = i + 1;
1549             stubs_comparison[0].has_data = true;
1550             stubs_comparison[0].quality_data = alctV_data[i].getQuality();
1551             stubs_comparison[0].key_WG_data = alctV_data[i].getKeyWG();
1552             stubs_comparison[0].bx_data = data_bx;
1553             stubs_comparison[0].fullbx_data = alctV_data[i].getFullBX();
1554             stubs_comparison[0].trknmb_data = data_trknmb;
1555             GlobalPoint gp_alct_data(getGlobalPosition(detid.rawId(), data_wiregroup, 60));
1556             stubs_comparison[0].eta_data = gp_alct_data.eta();
1557             for (int j = 0; j < nemul; j++) {
1558               if (alctV_emul[j].isValid() == 0)
1559                 continue;
1560               if (bookedalctV_emul[j])
1561                 continue;  //used alct
1562
1563               int emul_trknmb = alctV_emul[j].getTrknmb();
1564               //int emul_quality   = alctV_emul[j].getQuality();
1565               //int emul_accel     = alctV_emul[j].getAccelerator();
1566               //int emul_collB     = alctV_emul[j].getCollisionB();
1567               int emul_wiregroup = alctV_emul[j].getKeyWG();
1568               int emul_bx = alctV_emul[j].getBX();
1569               // Emulator BX re-calculated for comparison with BX in the data.
1570               emul_corr_bx = emul_bx - rawhit_tbin_offset + register_delay;
1571               if (dataIsAnotherMC_)
1572                 emul_corr_bx = emul_bx;
1573
1574               if (ndata == nemul) {
1575                 hAlctCompTotal->Fill(mychamber);
1576                 hAlctCompTotalCsc[endc - 1][csctype]->Fill(cham);
1577                 if (detid.station() > 1 && detid.ring() == 1) {
1578                   hAlctCompTotal2x->Fill(ix, detid.chamber() * 2);
1579                 } else {
1580                   hAlctCompTotal2->Fill(ix, detid.chamber());
1581                 }
1582                 hAlctCompTotal2i->Fill(ix2, detid.chamber());
1583               }
1584               if (data_trknmb == emul_trknmb) {
1585                 if (abs(data_wiregroup - emul_wiregroup) <= 2) {
1586                   hAlctCompMatchCsc[endc - 1][csctype]->Fill(cham);
1587                   hAlctCompMatch->Fill(mychamber);
1588                   if (detid.station() > 1 && detid.ring() == 1) {
1589                     hAlctCompMatch2x->Fill(ix, detid.chamber() * 2);
1590                   } else {
1591                     hAlctCompMatch2->Fill(ix, detid.chamber());
1592                   }
1593                   hAlctCompMatch2i->Fill(ix2, detid.chamber());
1594                 }
1595                 if (debug)
1596                   LogTrace("CSCTriggerPrimitivesReader") << "       Identical ALCTs #" << data_trknmb;
1597                 stubs_comparison[0].nStub_emul = j + 1;
1598                 stubs_comparison[0].has_emul = true;
1599                 stubs_comparison[0].quality_emul = alctV_emul[j].getQuality();
1600                 stubs_comparison[0].key_WG_emul = alctV_emul[j].getKeyWG();
1601                 stubs_comparison[0].bx_emul = alctV_emul[j].getBX();
1602                 stubs_comparison[0].trknmb_emul = emul_trknmb;
1603                 GlobalPoint gp_alct_emul(getGlobalPosition(detid.rawId(), emul_wiregroup, 60));
1604                 stubs_comparison[0].eta_emul = gp_alct_emul.eta();
1605                 stubs_comparison[0].bx_corr_emul = emul_corr_bx;
1606                 bookedalctV_emul[j] = true;
1607                 break;
1608               } else {
1609                 LogTrace("CSCTriggerPrimitivesReader")
1610                     << "       Different ALCTs #" << data_trknmb << " in ME" << ((endc == 1) ? "+" : "-") << stat << "/"
1611                     << ring << "/" << cham;
1612               }
1613
1614             }  //loop emul
1615             if (debug and stubs_comparison[0].key_WG_data != stubs_comparison[0].key_WG_emul)
1616               LogTrace("CSCTriggerPrimitivesReader")
1617                   << " not matched ALCT from data " << alctV_data[i] << " stubs_comparison 0 key_WG_data "
1618                   << stubs_comparison[0].key_WG_data << " key_WG_emul " << stubs_comparison[0].key_WG_emul;
1619             //if (stat==1) std::cout <<" stub_tree filled , ring "<< stubs_comparison[0].ring << std::endl;
1620             //cout <<"ALCT data BX "<< stubs_comparison[0].bx_data <<" WG "<< stubs_comparison[0].key_WG_data <<" emul BX "<< stubs_comparison[0].bx_emul<<" emul BX corrected "<< stubs_comparison[0].bx_corr_emul <<" WG "<< stubs_comparison[0].key_WG_emul << endl;
1621             stub_tree[0]->Fill();
1622           }  //loop data
1623           for (int i = 0; i < nemul; i++) {
1624             if (alctV_emul[i].isValid() == 0 or bookedalctV_emul[i])
1625               continue;
1626             int emul_bx = alctV_emul[i].getBX();
1627             stubs_comparison[0].init(RUN_, Event_);
1628             stubs_comparison[0].firstfill = firstfill;
1629             if (firstfill)
1630               firstfill = false;
1631             stubs_comparison[0].nEvents = eventsAnalyzed;
1632             stubs_comparison[0].endcap = endc;
1633             stubs_comparison[0].station = stat;
1634             stubs_comparison[0].ring = ring;
1635             stubs_comparison[0].chamber = cham;
1636             stubs_comparison[0].chambertype = detid.iChamberType();
1637             stubs_comparison[0].totStubs_data = ndata;
1638             stubs_comparison[0].totStubs_emul = nemul;
1639             stubs_comparison[0].totStubs_emul_readout = nemul;
1640             stubs_comparison[0].nStub_data = -1;
1641             stubs_comparison[0].nStub_emul = i + 1;
1642             stubs_comparison[0].has_data = false;
1643             stubs_comparison[0].has_emul = true;
1644             stubs_comparison[0].trknmb_emul = alctV_emul[i].getTrknmb();
1645             stubs_comparison[0].quality_emul = alctV_emul[i].getQuality();
1646             stubs_comparison[0].key_WG_emul = alctV_emul[i].getKeyWG();
1647             stubs_comparison[0].bx_emul = alctV_emul[i].getBX();
1648             stubs_comparison[0].fullbx_emul = alctV_emul[i].getFullBX();
1649             GlobalPoint gp_alct_emul(getGlobalPosition(detid.rawId(), alctV_emul[i].getKeyWG(), 60));
1650             stubs_comparison[0].eta_emul = gp_alct_emul.eta();
1651             // Emulator BX re-calculated for comparison with BX in the data.
1652             emul_corr_bx = emul_bx - rawhit_tbin_offset + register_delay;
1653             if (dataIsAnotherMC_)
1654               emul_corr_bx = emul_bx;
1655             stubs_comparison[0].bx_corr_emul = emul_corr_bx;
1656             stub_tree[0]->Fill();
1657             if (debug)
1658               LogTrace("CSCTriggerPrimitivesReader") << "not matched ALCT from emulation " << alctV_emul[i];
1659
1660           }  //loop emul
1661         }
1662       }
1663     }
1664   }
1665 }  // compare ALCT
1666
1667 void CSCTriggerPrimitivesReader::compareCLCTs(const CSCCLCTDigiCollection* clcts_data,
1668                                               const CSCCLCTDigiCollection* clcts_emul,
1669                                               const CSCCLCTPreTriggerDigiCollection* pretrigs_emul) {
1670   // Number of Tbins before pre-trigger for raw cathode hits.
1671   const int tbin_cathode_offset = 7;
1672   //const int tbin_cathode_offset = 8;//in MC, it became 8, Tao
1673   const int pretrig_trig_zone = 5;  // max distance between CLCT key hs and pretrigger hs
1674
1675   // Loop over all chambers in search for CLCTs.
1676   std::vector<CSCCLCTDigi>::const_iterator pd, pe;
1677   std::vector<CSCCLCTPreTriggerDigi>::const_iterator pretrig;
1678   perStub[2].init(RUN_, Event_);
1679   perStub[3].init(RUN_, Event_);
1680   for (int endc = 1; endc <= 2; endc++) {
1681     for (int stat = 1; stat <= 4; stat++) {
1682       for (int ring = 1; ring <= maxRing(stat); ring++) {
1683         for (int cham = 1; cham <= 36; cham++) {
1684           // Calculate DetId.  0th layer means whole chamber.
1685           CSCDetId detid(endc, stat, ring, cham, 0);
1686
1687           // Skip chambers marked as bad.
1688           if (checkBadChambers_ && badChambers_->isInBadChamber(detid))
1689             continue;
1690
1691           std::vector<CSCCLCTDigi> clctV_data, clctV_emul;
1692           const auto& drange = clcts_data->get(detid);
1693           for (auto digiIt = drange.first; digiIt != drange.second; digiIt++) {
1694             if ((*digiIt).isValid()) {
1695               clctV_data.push_back(*digiIt);
1696             }
1697           }
1698
1699           const auto& erange = clcts_emul->get(detid);
1700           std::vector<bool> bookedclctV_emul;
1701           for (auto digiIt = erange.first; digiIt != erange.second; digiIt++) {
1702             if ((*digiIt).isValid()) {
1703               for (const auto& pclct : clctV_emul) {
1704                 if (digiIt->getBX() != pclct.getBX() and abs(digiIt->getBX() - pclct.getBX()) < 5)
1705                   LogTrace("CSCTriggerPrimitivesReader")
1706                       << "Two CLCTs very close in timing!!! Special event: first clct " << pclct << " second clct "
1707                       << *digiIt << std::endl;
1708               }
1709               clctV_emul.push_back(*digiIt);
1710               bookedclctV_emul.push_back(false);
1711             }
1712           }
1713
1714           std::vector<CSCCLCTPreTriggerDigi> pretrigV_emul;
1715           const auto& pretrigrange = pretrigs_emul->get(detid);
1716           for (auto pretrigIt = pretrigrange.first; pretrigIt != pretrigrange.second; pretrigIt++) {
1717             if ((*pretrigIt).isValid()) {
1718               pretrigV_emul.push_back(*pretrigIt);
1719             }
1720           }
1721
1722           int ndata = clctV_data.size();
1723           int nemul = clctV_emul.size();
1724           if (ndata == 0 && nemul == 0)
1725             continue;
1726           bool firstfill = true;
1727           int nemul_readout = 0;
1728           for (pe = clctV_emul.begin(); pe != clctV_emul.end(); pe++) {
1729             for (pd = clctV_data.begin(); pd != clctV_data.end(); pd++) {
1730               int emul_bx = (*pe).getBX();
1731               int corr_bx = ((*pd).getFullBX() + emul_bx - tbin_cathode_offset) & 0x03;
1732               int bx_data = pd->getBX();
1733               //if (corr_bx == bx_data or abs(pe->getKeyStrip() - pd->getKeyStrip())<=1){//if emulated BX after correction is same as data bx, it will be readout
1734               if (corr_bx == bx_data) {  //if emulated BX after correction is same as data bx, it will be readout
1735                 nemul_readout++;
1736                 break;
1737               }
1738             }
1739           }
1740
1741           if (debug or nemul > ndata or (ndata != nemul) or ndata != nemul_readout) {
1742             LogTrace("CSCTriggerPrimitivesReader")
1743                 << " CLCTs from data " << ndata << " CLCTs from emul " << nemul << " readout " << nemul_readout;
1744             ostringstream strstrm;
1745             strstrm << "\n--- ME" << ((detid.endcap() == 1) ? "+" : "-") << detid.station() << "/" << detid.ring()
1746                     << "/" << detid.chamber() << " (sector " << detid.triggerSector() << " trig id. "
1747                     << detid.triggerCscId() << "):\n";
1748             strstrm << "  **** " << ndata << " valid data CLCTs found:\n";
1749             for (pd = clctV_data.begin(); pd != clctV_data.end(); pd++) {
1750               strstrm << "     " << (*pd) << " Full BX = " << (*pd).getFullBX() << "\n";
1751             }
1752             strstrm << "  **** " << nemul << " valid emul CLCTs found:\n";
1753             for (pe = clctV_emul.begin(); pe != clctV_emul.end(); pe++) {
1754               strstrm << "     " << (*pe);
1755               for (pd = clctV_data.begin(); pd != clctV_data.end(); pd++) {
1756                 if ((*pd).getTrknmb() == (*pe).getTrknmb() or abs((*pe).getKeyStrip() - (*pd).getKeyStrip()) < 5) {
1757                   int emul_bx = (*pe).getBX();
1758                   int corr_bx = ((*pd).getFullBX() + emul_bx - tbin_cathode_offset) & 0x03;
1759                   strstrm << " Corr BX = " << corr_bx;
1760                   break;
1761                 }
1762               }
1763               strstrm << "\n";
1764             }
1765             LogTrace("CSCTriggerPrimitivesReader") << strstrm.str();
1766           }
1767           if (nemul_readout > 2) {  //reduce nemul_readout to 2 by hand
1768             //cout <<"CLCT matching nemul readout is "<< nemul_readout <<", larger than 2. reduce it to 2 by hand"<<endl;
1769             nemul_readout = 2;
1770           }
1771
1772           //Per event Fill, From Luca
1773           //Data, add HS quality later
1774           int perEv_nStub_data = 0;
1775           int perEv_nStub_emul = 0;
1776           perStub[2].init(RUN_, Event_);
1777           perStub[3].init(RUN_, Event_);
1778           for (pd = clctV_data.begin(); pd != clctV_data.end(); pd++) {
1779             perEv_nStub_data++;
1780           }
1781           perStub[2].t_endcap = endc;
1782           perStub[2].t_station = stat;
1783           perStub[2].t_chamber = cham;
1784           perStub[2].t_ring = ring;
1785           perStub[2].t_EventNumberAnalyzed = eventsAnalyzed;
1786           perStub[2].t_nStubs = ndata;
1787           perStub[2].t_nStubs_readout = ndata;
1788           event_tree[2]->Fill();
1789           //Emul
1790           for (pe = clctV_emul.begin(); pe != clctV_emul.end(); pe++) {
1791             perEv_nStub_emul++;
1792           }
1793           perStub[3].t_endcap = endc;
1794           perStub[3].t_station = stat;
1795           perStub[3].t_chamber = cham;
1796           perStub[3].t_ring = ring;
1797           perStub[3].t_EventNumberAnalyzed = eventsAnalyzed;
1798           perStub[3].t_nStubs = nemul;
1799           perStub[3].t_nStubs_readout = nemul_readout;
1800           event_tree[3]->Fill();
1801
1802           int csctype = getCSCType(detid);
1803           hClctCompFoundCsc[endc - 1][csctype]->Fill(cham);
1804           int ix = chamberIX(detid);
1805           int ix2 = chamberIXi(detid);
1806           if (detid.station() > 1 && detid.ring() == 1) {
1807             hClctCompFound2x->Fill(ix, detid.chamber() * 2);
1808           } else {
1809             hClctCompFound2->Fill(ix, detid.chamber());
1810           }
1811           hClctCompFound2i->Fill(ix2, detid.chamber());
1812           if (ndata != nemul) {
1813             LogTrace("CSCTriggerPrimitivesReader")
1814                 //cerr
1815                 << "   +++ Different numbers of CLCTs found in ME" << ((endc == 1) ? "+" : "-") << stat << "/" << ring
1816                 << "/" << cham << ": data = " << ndata << " emulator = " << nemul << " +++\n";
1817           } else {
1818             hClctCompSameNCsc[endc - 1][csctype]->Fill(cham);
1819             if (detid.station() > 1 && detid.ring() == 1) {
1820               hClctCompSameN2x->Fill(ix, detid.chamber() * 2);
1821             } else {
1822               hClctCompSameN2->Fill(ix, detid.chamber());
1823             }
1824             hClctCompSameN2i->Fill(ix2, detid.chamber());
1825           }
1826           int i = -1;
1827           int testwg = 20;
1828           for (pd = clctV_data.begin(); pd != clctV_data.end(); pd++) {
1829             i++;
1830             if ((*pd).isValid() == 0)
1831               continue;
1832             int data_trknmb = (*pd).getTrknmb();
1833             int data_quality = (*pd).getQuality();
1834             int data_pattern = (*pd).getPattern();
1835             int data_striptype = (*pd).getStripType();
1836             int data_bend = (*pd).getBend();
1837             int data_keystrip = (*pd).getKeyStrip();
1838             int data_cfeb = (*pd).getCFEB();
1839             int data_bx = (*pd).getBX();
1840             int fullBX = (*pd).getFullBX();  // 12-bit full BX
1841
1842             if (data_keystrip >= 128 and stat == 1 and ring == 1)
1843               testwg = 5;
1844             else
1845               testwg = 20;
1846
1847             stubs_comparison[1].init(RUN_, Event_);
1848             stubs_comparison[1].firstfill = firstfill;
1849             if (firstfill)
1850               firstfill = false;
1851             stubs_comparison[1].nEvents = eventsAnalyzed;
1852             stubs_comparison[1].endcap = endc;
1853             stubs_comparison[1].station = stat;
1854             stubs_comparison[1].ring = ring;
1855             stubs_comparison[1].chamber = cham;
1856             stubs_comparison[1].chambertype = detid.iChamberType();
1857             stubs_comparison[1].totStubs_data = ndata;
1858             stubs_comparison[1].totStubs_emul = nemul;
1859             stubs_comparison[1].totStubs_emul_readout = nemul_readout;
1860             stubs_comparison[1].nStub_data = i + 1;
1861             stubs_comparison[1].has_data = true;
1862             stubs_comparison[1].quality_data = (*pd).getQuality();
1863             stubs_comparison[1].key_hs_data = (*pd).getKeyStrip();
1864             stubs_comparison[1].bend_data = (*pd).getBend();
1865             stubs_comparison[1].pattern_data = (*pd).getPattern();
1866             stubs_comparison[1].bx_data = (*pd).getBX();
1867             stubs_comparison[1].fullbx_data = (*pd).getFullBX();
1868             stubs_comparison[1].trknmb_data = data_trknmb;
1869             GlobalPoint gp_clct_data(getGlobalPosition(detid.rawId(), testwg, (*pd).getKeyStrip()));
1870             stubs_comparison[1].phi_data = gp_clct_data.phi();
1871
1872             int j = -1;
1873             for (pe = clctV_emul.begin(); pe != clctV_emul.end(); pe++) {
1874               j++;
1875               if ((*pe).isValid() == 0)
1876                 continue;
1877               if (bookedclctV_emul[j])
1878                 continue;  //used alct
1879               int emul_trknmb = (*pe).getTrknmb();
1880               int emul_quality = (*pe).getQuality();
1881               int emul_pattern = (*pe).getPattern();
1882               int emul_striptype = (*pe).getStripType();
1883               int emul_bend = (*pe).getBend();
1884               int emul_keystrip = (*pe).getKeyStrip();
1885               int emul_cfeb = (*pe).getCFEB();
1886               int emul_bx = (*pe).getBX();
1887
1888               if (abs(data_keystrip - emul_keystrip) <= 2) {
1889                 // Emulator BX re-calculated using 12-bit full BX number.
1890                 // Used for comparison with BX in the data.
1891                 int emul_corr_bx = (fullBX + emul_bx - tbin_cathode_offset) & 0x03;
1892                 //std::cout <<"CLCT data_bx "<< data_bx <<" emul_corr_bx "<< emul_corr_bx << std::endl;
1893                 if (dataIsAnotherMC_)
1894                   emul_corr_bx = (emul_bx & 0x03);
1895                 if (ndata == nemul) {
1896                   hClctCompTotalCsc[endc - 1][csctype]->Fill(cham);
1897                   if (detid.station() > 1 && detid.ring() == 1) {
1898                     hClctCompTotal2x->Fill(ix, detid.chamber() * 2);
1899                   } else {
1900                     hClctCompTotal2->Fill(ix, detid.chamber());
1901                   }
1902                   hClctCompTotal2i->Fill(ix2, detid.chamber());
1903                 }
1904                 if (data_quality == emul_quality && data_pattern == emul_pattern && data_striptype == emul_striptype &&
1905                     data_bend == emul_bend && data_keystrip == emul_keystrip && data_cfeb == emul_cfeb) {
1906                   if (ndata == nemul) {
1907                     hClctCompMatchCsc[endc - 1][csctype]->Fill(cham);
1908                     if (detid.station() > 1 && detid.ring() == 1) {
1909                       hClctCompMatch2x->Fill(ix, detid.chamber() * 2);
1910                     } else {
1911                       hClctCompMatch2->Fill(ix, detid.chamber());
1912                     }
1913                     hClctCompMatch2i->Fill(ix2, detid.chamber());
1914                   }
1915                   if (debug)
1916                     LogTrace("CSCTriggerPrimitivesReader") << "       Identical CLCTs #" << data_trknmb;
1917                 } else {
1918                   LogTrace("CSCTriggerPrimitivesReader")
1919                       << "       Different CLCTs #" << data_trknmb << " in ME" << ((endc == 1) ? "+" : "-") << stat
1920                       << "/" << ring << "/" << cham << " data_bx " << data_bx << " emul_corr_bx " << emul_corr_bx;
1921                 }
1922                 stubs_comparison[1].nStub_emul = j + 1;
1923                 stubs_comparison[1].has_emul = true;
1924                 stubs_comparison[1].quality_emul = (*pe).getQuality();
1925                 stubs_comparison[1].key_hs_emul = (*pe).getKeyStrip();
1926                 stubs_comparison[1].bend_emul = (*pe).getBend();
1927                 stubs_comparison[1].pattern_emul = (*pe).getPattern();
1928                 stubs_comparison[1].bx_emul = (*pe).getBX();
1929                 stubs_comparison[1].bx_corr_emul = emul_corr_bx;
1930                 stubs_comparison[1].trknmb_emul = emul_trknmb;
1931                 GlobalPoint gp_clct_emul(getGlobalPosition(detid.rawId(), testwg, (*pe).getKeyStrip()));
1932                 stubs_comparison[1].phi_emul = gp_clct_emul.phi();
1933                 bookedclctV_emul[j] = true;
1934
1935                 //int mingap_trig_pretrig = pretrig_trig_zone*2;
1936                 int num_pretrig = 0;
1937                 for (pretrig = pretrigV_emul.begin(); pretrig != pretrigV_emul.end(); pretrig++) {
1938                   if ((*pretrig).getBX() != (*pe).getBX())
1939                     continue;
1940                   int hsgap = std::abs((*pretrig).getKeyStrip() - (*pe).getKeyStrip());
1941                   bool samechamber = true;
1942                   if (stat == 1 and ring == 1)
1943                     samechamber = (((*pretrig).getKeyStrip() < 128 and (*pe).getKeyStrip() < 128) ||
1944                                    ((*pretrig).getKeyStrip() >= 128 and (*pe).getKeyStrip() >= 128));
1945                   if (not samechamber)
1946                     continue;
1947
1948                   if (hsgap <= pretrig_trig_zone)
1949                     num_pretrig++;
1950                   else
1951                     continue;
1952                   if ((*pretrig).getPattern() > stubs_comparison[1].maxpattern_pretrig)
1953                     stubs_comparison[1].maxpattern_pretrig = (*pretrig).getPattern();
1954                   if ((*pretrig).getQuality() > stubs_comparison[1].maxquality_pretrig) {
1955                     stubs_comparison[1].quality_pretrig = (*pretrig).getQuality();
1956                     stubs_comparison[1].key_hs_pretrig = (*pretrig).getKeyStrip();
1957                     stubs_comparison[1].bend_pretrig = (*pretrig).getBend();
1958                     stubs_comparison[1].pattern_pretrig = (*pretrig).getPattern();
1959                     stubs_comparison[1].bx_pretrig = (*pretrig).getBX();
1960                   }
1961                 }
1962                 stubs_comparison[1].npretrig = num_pretrig;
1963                 break;
1964               }  //if (data_trknmb == emul_trknmb)
1965             }    //loop emul
1966             if (debug and stubs_comparison[1].key_hs_data != stubs_comparison[1].key_hs_emul)
1967               LogTrace("CSCTriggerPrimitivesReader")
1968                   << "CSCDetid " << detid << " not matched CLCT from data " << (*pd)
1969                   << " stubs_comparison 1 key_hs_data " << stubs_comparison[1].key_hs_data << " key_hs_emul "
1970                   << stubs_comparison[1].key_hs_emul << endl;
1971             //cout <<"CLCT data BX "<< stubs_comparison[1].bx_data <<" emul BX "<< stubs_comparison[1].bx_emul<<" emul BX corrected "<< stubs_comparison[1].bx_corr_emul << endl;
1972             stub_tree[1]->Fill();
1973           }  //loop data
1974           for (int k = 0; k < nemul; k++) {
1975             if (clctV_emul[i].isValid() == 0)
1976               continue;
1977             if (bookedclctV_emul[k])
1978               continue;
1979             stubs_comparison[1].init(RUN_, Event_);
1980             stubs_comparison[1].firstfill = firstfill;
1981             if (firstfill)
1982               firstfill = false;
1983             stubs_comparison[1].nEvents = eventsAnalyzed;
1984             stubs_comparison[1].endcap = endc;
1985             stubs_comparison[1].station = stat;
1986             stubs_comparison[1].ring = ring;
1987             stubs_comparison[1].chamber = cham;
1988             stubs_comparison[1].chambertype = detid.iChamberType();
1989             stubs_comparison[1].totStubs_data = ndata;
1990             stubs_comparison[1].totStubs_emul = nemul;
1991             stubs_comparison[1].totStubs_emul_readout = nemul_readout;
1992             stubs_comparison[1].trknmb_emul = clctV_emul[i].getTrknmb();
1993             stubs_comparison[1].nStub_data = -1;
1994             stubs_comparison[1].has_data = false;
1995             stubs_comparison[1].nStub_emul = k + 1;
1996             stubs_comparison[1].has_emul = true;
1997             stubs_comparison[1].quality_emul = clctV_emul[k].getQuality();
1998             stubs_comparison[1].key_hs_emul = clctV_emul[k].getKeyStrip();
1999             stubs_comparison[1].bend_emul = clctV_emul[k].getBend();
2000             stubs_comparison[1].pattern_emul = clctV_emul[k].getPattern();
2001             stubs_comparison[1].bx_emul = clctV_emul[k].getBX();
2002             stubs_comparison[1].fullbx_emul = clctV_emul[k].getFullBX();
2003             if (clctV_emul[k].getKeyStrip() >= 128 and stat == 1 and ring == 1)
2004               testwg = 5;
2005             else
2006               testwg = 20;
2007             // Emulator BX NOT Known from  the data.
2008             GlobalPoint gp_clct_emul(getGlobalPosition(detid.rawId(), testwg, clctV_emul[k].getKeyStrip()));
2009             stubs_comparison[1].phi_emul = gp_clct_emul.phi();
2010             bookedclctV_emul[k] = true;
2011             stub_tree[1]->Fill();
2012           }
2013         }
2014       }
2015     }
2016   }  // for all chambers
2017 }  // compareCLCTs
2018
2019 void CSCTriggerPrimitivesReader::compareLCTs(const CSCCorrelatedLCTDigiCollection* lcts_data,
2020                                              const CSCCorrelatedLCTDigiCollection* lcts_emul,
2021                                              const CSCALCTDigiCollection* alcts_data,
2022                                              const CSCCLCTDigiCollection* clcts_data) {
2023   // Need ALCT and CLCT digi collections to convert emulator bx into
2024   // hardware bx.
2025   // Loop over all chambers in search for correlated LCTs.
2026   std::vector<CSCCorrelatedLCTDigi>::const_iterator pd, pe;
2027   for (int endc = 1; endc <= 2; endc++) {
2028     for (int stat = 1; stat <= 4; stat++) {
2029       for (int ring = 1; ring <= maxRing(stat); ring++) {
2030         for (int cham = 1; cham <= 36; cham++) {
2031           // Calculate DetId.  0th layer means whole chamber.
2032           CSCDetId detid(endc, stat, ring, cham, 0);
2033
2034           // Skip chambers marked as bad.
2035           if (checkBadChambers_ && badChambers_->isInBadChamber(detid))
2036             continue;
2037
2038           std::vector<CSCCorrelatedLCTDigi> lctV_data, lctV_emul;
2039           std::vector<bool> bookedlctV_emul;
2040           const auto& drange = lcts_data->get(detid);
2041           for (auto digiIt = drange.first; digiIt != drange.second; digiIt++) {
2042             if ((*digiIt).isValid()) {
2043               lctV_data.push_back(*digiIt);
2044             }
2045           }
2046
2047           const auto& erange = lcts_emul->get(detid);
2048           for (auto digiIt = erange.first; digiIt != erange.second; digiIt++) {
2049             if ((*digiIt).isValid()) {
2050               lctV_emul.push_back(*digiIt);
2051               bookedlctV_emul.push_back(false);
2052             }
2053           }
2054
2055           int ndata = lctV_data.size();
2056           int nemul = lctV_emul.size();
2057           if (ndata == 0 && nemul == 0)
2058             continue;
2059           bool firstfill = true;
2060           int nemul_readout = 0;
2061           for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2062             for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2063               int bx_data = pd->getBX();
2064               int bx_corr = convertBXofLCT((*pe).getBX(), detid, alcts_data, clcts_data);
2065               if ((bx_data == bx_corr) or (abs(pe->getKeyWG() - pd->getKeyWG()) <= 2)) {
2066                 nemul_readout++;
2067                 break;
2068               }
2069             }
2070           }
2071
2072           if (debug or nemul > ndata or (nemul_readout != ndata) or nemul != ndata) {
2073             LogTrace("CSCTriggerPrimitivesReader")
2074                 << " LCTs from data " << ndata << " LCTs from emul " << nemul << " readout " << nemul_readout;
2075             ostringstream strstrm;
2076             strstrm << "\n--- ME" << ((detid.endcap() == 1) ? "+" : "-") << detid.station() << "/" << detid.ring()
2077                     << "/" << detid.chamber() << " (sector " << detid.triggerSector() << " trig id. "
2078                     << detid.triggerCscId() << "):\n";
2079             strstrm << "  **** " << ndata << " valid data LCTs found:\n";
2080             for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2081               strstrm << "     " << (*pd);
2082             }
2083             strstrm << "\n  **** " << nemul << " valid emul LCTs found:\n";
2084             for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2085               strstrm << "     " << (*pe);
2086               strstrm << "    corr BX = " << convertBXofLCT((*pe).getBX(), detid, alcts_data, clcts_data);
2087               strstrm << " LCT pattern = " << (*pe).getPattern();
2088               strstrm << "\n";
2089             }
2090             LogTrace("CSCTriggerPrimitivesReader") << strstrm.str();
2091             //std::cout <<"CompareLCTs "<< strstrm.str()<< std::endl;
2092           }
2093           if (nemul_readout > 2) {  //reduce nemul_readout to 2 by hand
2094             //cout <<"LCT matching nemul readout is "<< nemul_readout <<", larger than 2. reduce it to 2 by hand"<<endl;
2095             nemul_readout = 2;
2096           }
2097
2098           //
2099           //Data, add HS quality later
2100           int perEv_nStub_data = 0;
2101           int perEv_nStub_emul = 0;
2102           perStub[4].init(RUN_, Event_);
2103           perStub[5].init(RUN_, Event_);
2104           for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2105             perEv_nStub_data++;
2106           }
2107           perStub[4].t_endcap = endc;
2108           perStub[4].t_station = stat;
2109           perStub[4].t_chamber = cham;
2110           perStub[4].t_ring = ring;
2111           perStub[4].t_EventNumberAnalyzed = eventsAnalyzed;
2112           perStub[4].t_nStubs = ndata;
2113           perStub[4].t_nStubs_readout = ndata;
2114           event_tree[4]->Fill();
2115           //Emul
2116           for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2117             perEv_nStub_emul++;
2118           }
2119           perStub[5].t_endcap = endc;
2120           perStub[5].t_station = stat;
2121           perStub[5].t_chamber = cham;
2122           perStub[5].t_ring = ring;
2123           perStub[5].t_EventNumberAnalyzed = eventsAnalyzed;
2124           perStub[5].t_nStubs = nemul;
2125           perStub[5].t_nStubs_readout = nemul_readout;
2126           event_tree[5]->Fill();
2127
2128           int csctype = getCSCType(detid);
2129           hLctCompFoundCsc[endc - 1][csctype]->Fill(cham);
2130           int ix = chamberIX(detid);
2131           int ix2 = chamberIXi(detid);
2132           if (detid.station() > 1 && detid.ring() == 1) {
2133             hLCTCompFound2x->Fill(ix, detid.chamber() * 2);
2134           } else {
2135             hLCTCompFound2->Fill(ix, detid.chamber());
2136           }
2137           hLCTCompFound2i->Fill(ix2, detid.chamber());
2138           if (ndata != nemul) {
2139             LogTrace("CSCTriggerPrimitivesReader")
2140                 //cerr
2141                 << "   +++ Different numbers of LCTs found in ME" << ((endc == 1) ? "+" : "-") << stat << "/" << ring
2142                 << "/" << cham << ": data = " << ndata << " emulator = " << nemul << " +++\n";
2143           } else {
2144             hLctCompSameNCsc[endc - 1][csctype]->Fill(cham);
2145             if (detid.station() > 1 && detid.ring() == 1) {
2146               hLCTCompSameN2x->Fill(ix, detid.chamber() * 2);
2147             } else {
2148               hLCTCompSameN2->Fill(ix, detid.chamber());
2149             }
2150             hLCTCompSameN2i->Fill(ix2, detid.chamber());
2151           }
2152           int i = -1;
2153           for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2154             i++;
2155             if ((*pd).isValid() == 0)
2156               continue;
2157             int data_trknmb = (*pd).getTrknmb();
2158             int data_quality = (*pd).getQuality();
2159             int data_wiregroup = (*pd).getKeyWG();
2160             int data_keystrip = (*pd).getStrip();
2161             int data_pattern = (*pd).getCLCTPattern();
2162             int data_striptype = (*pd).getStripType();
2163             int data_bend = (*pd).getBend();
2164             int data_bx = (*pd).getBX();
2165
2166             stubs_comparison[2].init(RUN_, Event_);
2167             stubs_comparison[2].firstfill = firstfill;
2168             if (firstfill)
2169               firstfill = false;
2170             stubs_comparison[2].nEvents = eventsAnalyzed;
2171             stubs_comparison[2].endcap = endc;
2172             stubs_comparison[2].station = stat;
2173             stubs_comparison[2].ring = ring;
2174             stubs_comparison[2].chamber = cham;
2175             stubs_comparison[2].chambertype = detid.iChamberType();
2176             stubs_comparison[2].totStubs_data = ndata;
2177             stubs_comparison[2].totStubs_emul = nemul;
2178             stubs_comparison[2].totStubs_emul_readout = nemul_readout;
2179             stubs_comparison[2].nStub_data = i + 1;
2180             stubs_comparison[2].has_data = true;
2181             stubs_comparison[2].quality_data = (*pd).getQuality();
2182             stubs_comparison[2].key_WG_data = (*pd).getKeyWG();
2183             stubs_comparison[2].key_hs_data = (*pd).getStrip();
2184             //cout <<" CSC detid "<< detid <<" LCT in data "<< (*pd) <<" do doesALCTCrossCLCT "<< endl;
2185             stubs_comparison[2].WGcrossHS_data = doesALCTCrossCLCT(detid, (*pd).getKeyWG(), (*pd).getStrip());
2186             stubs_comparison[2].bend_data = (*pd).getBend();
2187             stubs_comparison[2].pattern_data = (*pd).getCLCTPattern();
2188             stubs_comparison[2].bx_data = (*pd).getBX();
2189             stubs_comparison[2].trknmb_data = data_trknmb;
2190             GlobalPoint gp_lct_data(getGlobalPosition(detid.rawId(), (*pd).getKeyWG(), (*pd).getStrip()));
2191             stubs_comparison[2].eta_data = gp_lct_data.eta();
2192             stubs_comparison[2].phi_data = gp_lct_data.phi();
2193             int j = -1;
2194             for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2195               j++;
2196               if ((*pe).isValid() == 0)
2197                 continue;
2198               if (bookedlctV_emul[j])
2199                 continue;  //used alct
2200               int emul_trknmb = (*pe).getTrknmb();
2201               int emul_quality = (*pe).getQuality();
2202               int emul_wiregroup = (*pe).getKeyWG();
2203               int emul_keystrip = (*pe).getStrip();
2204               int emul_pattern = (*pe).getCLCTPattern();
2205               int emul_striptype = (*pe).getStripType();
2206               int emul_bend = (*pe).getBend();
2207               int emul_bx = (*pe).getBX();
2208               if (abs(data_keystrip - emul_keystrip) <= 2 || abs(data_wiregroup - emul_wiregroup) <= 2) {
2209                 // Convert emulator BX into hardware BX using full 12-bit
2210                 // BX words in ALCT and CLCT digi collections.
2211                 int emul_corr_bx = convertBXofLCT(emul_bx, detid, alcts_data, clcts_data);
2212                 //std::cout <<"LCT data_bx "<< data_bx <<" emul_corr_bx "<< emul_corr_bx << std::endl;
2213                 if (dataIsAnotherMC_)
2214                   emul_corr_bx = (emul_bx & 0x01);
2215
2216                 if (ndata == nemul) {
2217                   hLctCompTotalCsc[endc - 1][csctype]->Fill(cham);
2218                   if (detid.station() > 1 && detid.ring() == 1) {
2219                     hLCTCompTotal2x->Fill(ix, detid.chamber() * 2);
2220                   } else {
2221                     hLCTCompTotal2->Fill(ix, detid.chamber());
2222                   }
2223                   hLCTCompTotal2i->Fill(ix2, detid.chamber());
2224                 }
2225                 if (data_quality == emul_quality && data_wiregroup == emul_wiregroup &&
2226                     data_keystrip == emul_keystrip && data_pattern == emul_pattern &&
2227                     data_striptype == emul_striptype && data_bend == emul_bend) {
2228                   if (ndata == nemul) {
2229                     hLctCompMatchCsc[endc - 1][csctype]->Fill(cham);
2230                     if (detid.station() > 1 && detid.ring() == 1) {
2231                       hLCTCompMatch2x->Fill(ix, detid.chamber() * 2);
2232                     } else {
2233                       hLCTCompMatch2->Fill(ix, detid.chamber());
2234                     }
2235                     hLCTCompMatch2i->Fill(ix2, detid.chamber());
2236                   }
2237                   if (debug)
2238                     LogTrace("CSCTriggerPrimitivesReader") << "       Identical LCTs #" << data_trknmb;
2239                 } else {
2240                   LogTrace("CSCTriggerPrimitivesReader")
2241                       << "       Different LCTs #" << data_trknmb << " in ME" << ((endc == 1) ? "+" : "-") << stat
2242                       << "/" << ring << "/" << cham << " data_lct bx " << data_bx << " emul bx " << emul_corr_bx;
2243                 }
2244                 stubs_comparison[2].nStub_emul = j + 1;
2245                 stubs_comparison[2].has_emul = true;
2246                 stubs_comparison[2].quality_emul = (*pe).getQuality();
2247                 stubs_comparison[2].key_WG_emul = (*pe).getKeyWG();
2248                 stubs_comparison[2].key_hs_emul = (*pe).getStrip();
2249                 stubs_comparison[2].bend_emul = (*pe).getBend();
2250                 stubs_comparison[2].pattern_emul = (*pe).getCLCTPattern();
2251                 stubs_comparison[2].bx_emul = (*pe).getBX();
2252                 //cout <<" CSC detid "<< detid <<" LCT in emul "<< (*pe) <<" do doesALCTCrossCLCT "<< endl;
2253                 stubs_comparison[2].WGcrossHS_emul = doesALCTCrossCLCT(detid, (*pe).getKeyWG(), (*pe).getStrip());
2254                 stubs_comparison[2].bx_corr_emul = emul_corr_bx;
2255                 stubs_comparison[2].trknmb_emul = emul_trknmb;
2256                 GlobalPoint gp_lct_emul(getGlobalPosition(detid.rawId(), (*pe).getKeyWG(), (*pe).getStrip()));
2257                 stubs_comparison[2].eta_emul = gp_lct_emul.eta();
2258                 stubs_comparison[2].phi_emul = gp_lct_emul.phi();
2259                 bookedlctV_emul[j] = true;
2260                 break;
2261               }  //if (data_trknmb == emul_trknmb)
2262             }    //loop emul
2263             if (debug and stubs_comparison[2].key_hs_data != stubs_comparison[2].key_hs_emul)
2264               LogTrace("CSCTriggerPrimitivesReader")
2265                   << "CSCDetid " << detid << " not matched LCT from Data " << (*pd)
2266                   << " stubs_comparison 2 key_hs_data " << stubs_comparison[2].key_hs_data << " key_hs_emul "
2267                   << stubs_comparison[2].key_hs_emul << endl;
2268             stub_tree[2]->Fill();
2269           }  //loop data
2270           for (int k = 0; k < nemul; k++) {
2271             if (bookedlctV_emul[k])
2272               continue;
2273             if (lctV_emul[i].isValid() == 0)
2274               continue;
2275             stubs_comparison[2].init(RUN_, Event_);
2276             stubs_comparison[2].nEvents = eventsAnalyzed;
2277             stubs_comparison[2].endcap = endc;
2278             stubs_comparison[2].station = stat;
2279             stubs_comparison[2].ring = ring;
2280             stubs_comparison[2].chamber = cham;
2281             stubs_comparison[2].chambertype = detid.iChamberType();
2282             stubs_comparison[2].totStubs_data = ndata;
2283             stubs_comparison[2].totStubs_emul = nemul;
2284             stubs_comparison[2].totStubs_emul_readout = nemul_readout;
2285             stubs_comparison[2].trknmb_emul = lctV_emul[i].getTrknmb();
2286             stubs_comparison[2].nStub_data = -1;
2287             stubs_comparison[2].has_data = false;
2288             stubs_comparison[2].nStub_emul = k + 1;
2289             stubs_comparison[2].has_emul = true;
2290             stubs_comparison[2].WGcrossHS_emul =
2291                 doesALCTCrossCLCT(detid, lctV_emul[k].getKeyWG(), lctV_emul[k].getStrip());
2292             stubs_comparison[2].quality_emul = lctV_emul[k].getQuality();
2293             stubs_comparison[2].key_WG_emul = lctV_emul[k].getKeyWG();
2294             stubs_comparison[2].key_hs_emul = lctV_emul[k].getStrip();
2295             stubs_comparison[2].bend_emul = lctV_emul[k].getBend();
2296             stubs_comparison[2].pattern_emul = lctV_emul[k].getCLCTPattern();
2297             stubs_comparison[2].bx_emul = lctV_emul[k].getBX();
2298             //stubs_comparison[2].fullbx_emul = lctV_emul[k].getFullBX();
2299             // Emulator BX NOT Known from  the data.
2300             GlobalPoint gp_lct_emul(getGlobalPosition(detid.rawId(), lctV_emul[k].getKeyWG(), lctV_emul[k].getStrip()));
2301             stubs_comparison[2].eta_emul = gp_lct_emul.eta();
2302             stubs_comparison[2].phi_emul = gp_lct_emul.phi();
2303             bookedlctV_emul[k] = true;
2304             stub_tree[2]->Fill();
2305           }
2306         }
2307       }
2308     }
2309   }
2310 }
2311
2312 void CSCTriggerPrimitivesReader::compareMPCLCTs(const CSCCorrelatedLCTDigiCollection* lcts_data,
2313                                                 const CSCCorrelatedLCTDigiCollection* lcts_emul,
2314                                                 const CSCALCTDigiCollection* alcts_data,
2315                                                 const CSCCLCTDigiCollection* clcts_data) {
2316   // Need ALCT and CLCT digi collections to convert emulator bx into
2317   // hardware bx.
2318   // Loop over all chambers in search for correlated LCTs.
2319   std::vector<CSCCorrelatedLCTDigi>::const_iterator pd, pe;
2320   for (int endc = 1; endc <= 2; endc++) {
2321     for (int stat = 1; stat <= 4; stat++) {
2322       for (int ring = 1; ring <= maxRing(stat); ring++) {
2323         for (int cham = 1; cham <= 36; cham++) {
2324           // Calculate DetId.  0th layer means whole chamber.
2325           CSCDetId detid(endc, stat, ring, cham, 0);
2326
2327           // Skip chambers marked as bad.
2328           if (checkBadChambers_ && badChambers_->isInBadChamber(detid))
2329             continue;
2330
2331           std::vector<CSCCorrelatedLCTDigi> lctV_data, lctV_emul;
2332           std::vector<bool> bookedlctV_emul;
2333           const auto& drange = lcts_data->get(detid);
2334           for (auto digiIt = drange.first; digiIt != drange.second; digiIt++) {
2335             if ((*digiIt).isValid()) {
2336               lctV_data.push_back(*digiIt);
2337             }
2338           }
2339
2340           const auto& erange = lcts_emul->get(detid);
2341           for (auto digiIt = erange.first; digiIt != erange.second; digiIt++) {
2342             if ((*digiIt).isValid()) {
2343               lctV_emul.push_back(*digiIt);
2344               bookedlctV_emul.push_back(false);
2345             }
2346           }
2347
2348           int ndata = lctV_data.size();
2349           int nemul = lctV_emul.size();
2350           if (ndata == 0 && nemul == 0)
2351             continue;
2352           bool firstfill = true;
2353
2354           if (debug) {
2355             ostringstream strstrm;
2356             strstrm << "\n--- ME" << ((detid.endcap() == 1) ? "+" : "-") << detid.station() << "/" << detid.ring()
2357                     << "/" << detid.chamber() << " (sector " << detid.triggerSector() << " trig id. "
2358                     << detid.triggerCscId() << "):\n";
2359             strstrm << "  **** " << ndata << " valid data LCTs found:\n";
2360             for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2361               strstrm << "     " << (*pd);
2362             }
2363             strstrm << "\n  **** " << nemul << " valid emul LCTs found:\n";
2364             for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2365               strstrm << "     " << (*pe);
2366               strstrm << "    corr BX = " << convertBXofLCT((*pe).getBX(), detid, alcts_data, clcts_data);
2367               strstrm << " LCT pattern = " << (*pe).getPattern();
2368               strstrm << "\n";
2369             }
2370             LogTrace("CSCTriggerPrimitivesReader") << strstrm.str();
2371             //std::cout <<"CompareLCTs "<< strstrm.str()<< std::endl;
2372           }
2373
2374           //int ix = chamberIX(detid);
2375           //int ix2 = chamberIXi(detid);
2376           if (ndata != nemul) {
2377             LogTrace("CSCTriggerPrimitivesReader")
2378                 << "   +++ Different numbers of MPC LCTs found in ME" << ((endc == 1) ? "+" : "-") << stat << "/"
2379                 << ring << "/" << cham << ": data = " << ndata << " emulator = " << nemul << " +++\n";
2380           }
2381           /*else {
2382             hLctCompSameNCsc[endc-1][csctype]->Fill(cham);
2383             if(detid.station()>1 && detid.ring()==1) {
2384             hLCTCompSameN2x->Fill(ix,detid.chamber()*2);
2385             }
2386             else {
2387             hLCTCompSameN2->Fill(ix,detid.chamber());
2388             }
2389             hLCTCompSameN2i->Fill(ix2,detid.chamber());
2390             }*/
2391           int i = -1;
2392           for (pd = lctV_data.begin(); pd != lctV_data.end(); pd++) {
2393             i++;
2394             if ((*pd).isValid() == 0)
2395               continue;
2396             int data_trknmb = (*pd).getTrknmb();
2397             int data_quality = (*pd).getQuality();
2398             int data_wiregroup = (*pd).getKeyWG();
2399             int data_keystrip = (*pd).getStrip();
2400             int data_pattern = (*pd).getCLCTPattern();
2401             int data_striptype = (*pd).getStripType();
2402             int data_bend = (*pd).getBend();
2403             int data_bx = (*pd).getBX();
2404
2405             stubs_comparison[3].init(RUN_, Event_);
2406             stubs_comparison[3].firstfill = firstfill;
2407             if (firstfill)
2408               firstfill = false;
2409             stubs_comparison[3].nEvents = eventsAnalyzed;
2410             stubs_comparison[3].endcap = endc;
2411             stubs_comparison[3].station = stat;
2412             stubs_comparison[3].ring = ring;
2413             stubs_comparison[3].chamber = cham;
2414             stubs_comparison[3].chambertype = detid.iChamberType();
2415             stubs_comparison[3].totStubs_data = ndata;
2416             stubs_comparison[3].totStubs_emul = nemul;
2417             stubs_comparison[3].nStub_data = i + 1;
2418             stubs_comparison[3].has_data = true;
2419             stubs_comparison[3].WGcrossHS_data = doesALCTCrossCLCT(detid, (*pd).getKeyWG(), (*pd).getStrip());
2420             stubs_comparison[3].quality_data = (*pd).getQuality();
2421             stubs_comparison[3].key_WG_data = (*pd).getKeyWG();
2422             stubs_comparison[3].key_hs_data = (*pd).getStrip();
2423             stubs_comparison[3].bend_data = (*pd).getBend();
2424             stubs_comparison[3].pattern_data = (*pd).getCLCTPattern();
2425             stubs_comparison[3].bx_data = (*pd).getBX();
2426             stubs_comparison[3].trknmb_data = data_trknmb;
2427             GlobalPoint gp_lct_data(getGlobalPosition(detid.rawId(), (*pd).getKeyWG(), (*pd).getStrip()));
2428             stubs_comparison[3].eta_data = gp_lct_data.eta();
2429             stubs_comparison[3].phi_data = gp_lct_data.phi();
2430             int j = -1;
2431             for (pe = lctV_emul.begin(); pe != lctV_emul.end(); pe++) {
2432               j++;
2433               if ((*pe).isValid() == 0)
2434                 continue;
2435               int emul_trknmb = (*pe).getTrknmb();
2436               int emul_quality = (*pe).getQuality();
2437               int emul_wiregroup = (*pe).getKeyWG();
2438               int emul_keystrip = (*pe).getStrip();
2439               int emul_pattern = (*pe).getCLCTPattern();
2440               int emul_striptype = (*pe).getStripType();
2441               int emul_bend = (*pe).getBend();
2442               int emul_bx = (*pe).getBX();
2443               if (data_trknmb == emul_trknmb) {
2444                 // Convert emulator BX into hardware BX using full 12-bit
2445                 // BX words in ALCT and CLCT digi collections.
2446                 int emul_corr_bx = convertBXofLCT(emul_bx, detid, alcts_data, clcts_data);
2447                 //std::cout <<"LCT data_bx "<< data_bx <<" emul_corr_bx "<< emul_corr_bx << std::endl;
2448                 if (dataIsAnotherMC_)
2449                   emul_corr_bx = (emul_bx & 0x01);
2450
2451                 /*if (ndata == nemul) {
2452                   hLctCompTotalCsc[endc-1][csctype]->Fill(cham);
2453                   if(detid.station()>1 && detid.ring()==1) {
2454                   hLCTCompTotal2x->Fill(ix,detid.chamber()*2);
2455                   }
2456                   else {
2457                   hLCTCompTotal2->Fill(ix,detid.chamber());
2458                   }
2459                   hLCTCompTotal2i->Fill(ix2,detid.chamber());
2460                   }*/
2461                 if (data_quality == emul_quality && data_wiregroup == emul_wiregroup &&
2462                     data_keystrip == emul_keystrip && data_pattern == emul_pattern &&
2463                     data_striptype == emul_striptype && data_bend == emul_bend) {
2464                   //data_bx        == emul_corr_bx) {
2465                   /*if (ndata == nemul) {
2466                     hLctCompMatchCsc[endc-1][csctype]->Fill(cham);
2467                     if(detid.station()>1 && detid.ring()==1) {
2468                     hLCTCompMatch2x->Fill(ix,detid.chamber()*2);
2469                     }
2470                     else {
2471                     hLCTCompMatch2->Fill(ix,detid.chamber());
2472                     }
2473                     hLCTCompMatch2i->Fill(ix2,detid.chamber());
2474                     }*/
2475                   if (debug)
2476                     LogTrace("CSCTriggerPrimitivesReader") << "       Identical LCTs #" << data_trknmb;
2477                   stubs_comparison[3].nStub_emul = j + 1;
2478                   stubs_comparison[3].has_emul = true;
2479                   stubs_comparison[3].WGcrossHS_emul = doesALCTCrossCLCT(detid, (*pe).getKeyWG(), (*pe).getStrip());
2480                   stubs_comparison[3].quality_emul = (*pe).getQuality();
2481                   stubs_comparison[3].key_WG_emul = (*pe).getKeyWG();
2482                   stubs_comparison[3].key_hs_emul = (*pe).getStrip();
2483                   stubs_comparison[3].bend_emul = (*pe).getBend();
2484                   stubs_comparison[3].pattern_emul = (*pe).getCLCTPattern();
2485                   stubs_comparison[3].bx_emul = (*pe).getBX();
2486                   stubs_comparison[3].bx_corr_emul = emul_corr_bx;
2487                   stubs_comparison[3].trknmb_emul = emul_trknmb;
2488                   GlobalPoint gp_lct_emul(getGlobalPosition(detid.rawId(), (*pe).getKeyWG(), (*pe).getStrip()));
2489                   stubs_comparison[3].eta_emul = gp_lct_emul.eta();
2490                   stubs_comparison[3].phi_emul = gp_lct_emul.phi();
2491                   bookedlctV_emul[j] = true;
2492                   break;
2493                 } else {
2494                   LogTrace("CSCTriggerPrimitivesReader")
2495                       //cerr
2496                       << "       Different LCTs #" << data_trknmb << " in ME" << ((endc == 1) ? "+" : "-") << stat
2497                       << "/" << ring << "/" << cham << " data_lct bx " << data_bx << " emul bx " << emul_corr_bx;
2498                 }
2499               }
2500             }  //loop emul
2501             if (debug and stubs_comparison[3].key_hs_data != stubs_comparison[3].key_hs_emul)
2502               LogTrace("CSCTriggerPrimitivesReader")
2503                   << "stubs_comparison 2 key_hs_data " << stubs_comparison[3].key_hs_data << " key_hs_emul "
2504                   << stubs_comparison[3].key_hs_emul;
2505             stub_tree[3]->Fill();
2506           }  //loop data
2507
2508           for (int k = 0; k < nemul; k++) {
2509             if (bookedlctV_emul[k])
2510               continue;
2511             stubs_comparison[3].init(RUN_, Event_);
2512             stubs_comparison[3].firstfill = firstfill;
2513             if (firstfill)
2514               firstfill = false;
2515             stubs_comparison[3].nEvents = eventsAnalyzed;
2516             stubs_comparison[3].endcap = endc;
2517             stubs_comparison[3].station = stat;
2518             stubs_comparison[3].ring = ring;
2519             stubs_comparison[3].chamber = cham;
2520             stubs_comparison[3].chambertype = detid.iChamberType();
2521             stubs_comparison[3].totStubs_data = ndata;
2522             stubs_comparison[3].totStubs_emul = nemul;
2523             stubs_comparison[3].trknmb_emul = lctV_emul[i].getTrknmb();
2524             stubs_comparison[3].nStub_data = -1;
2525             stubs_comparison[3].has_data = false;
2526             stubs_comparison[3].nStub_emul = k + 1;
2527             stubs_comparison[3].has_emul = true;
2528             stubs_comparison[3].WGcrossHS_emul =
2529                 doesALCTCrossCLCT(detid, lctV_emul[k].getKeyWG(), lctV_emul[k].getStrip());
2530             stubs_comparison[3].quality_emul = lctV_emul[k].getQuality();
2531             stubs_comparison[3].key_WG_emul = lctV_emul[k].getKeyWG();
2532             stubs_comparison[3].key_hs_emul = lctV_emul[k].getStrip();
2533             stubs_comparison[3].bend_emul = lctV_emul[k].getBend();
2534             stubs_comparison[3].pattern_emul = lctV_emul[k].getCLCTPattern();
2535             stubs_comparison[3].bx_emul = lctV_emul[k].getBX();
2536             //stubs_comparison[2].fullbx_emul = lctV_emul[k].getFullBX();
2537             // Emulator BX NOT Known from  the data.
2538             GlobalPoint gp_lct_emul(getGlobalPosition(detid.rawId(), lctV_emul[k].getKeyWG(), lctV_emul[k].getStrip()));
2539             stubs_comparison[3].eta_emul = gp_lct_emul.eta();
2540             stubs_comparison[3].phi_emul = gp_lct_emul.phi();
2541             bookedlctV_emul[k] = true;
2542             stub_tree[3]->Fill();
2543           }
2544         }  //end loop of chambers
2545       }
2546     }
2547   }
2548 }
2549
2550 int CSCTriggerPrimitivesReader::convertBXofLCT(const int emul_bx,
2551                                                const CSCDetId& detid,
2552                                                const CSCALCTDigiCollection* alcts_data,
2553                                                const CSCCLCTDigiCollection* clcts_data) {
2554   int full_anode_bx = -999;
2555   //int full_cathode_bx = -999;
2556   int lct_bx = -999;
2557   int tbin_anode_offset = 5;  // 2007, run 14419.
2558
2559   // Extract full 12-bit anode BX word from ALCT collections.
2560   const auto& arange = alcts_data->get(detid);
2561   for (auto digiIt = arange.first; digiIt != arange.second; digiIt++) {
2562     if ((*digiIt).isValid()) {
2563       full_anode_bx = (*digiIt).getFullBX();
2564       break;
2565     }
2566   }
2567
2568   // Extract full 12-bit cathode BX word from CLCT collections.
2569   const auto& crange = clcts_data->get(detid);
2570   for (auto digiIt = crange.first; digiIt != crange.second; digiIt++) {
2571     if ((*digiIt).isValid()) {
2572       //full_cathode_bx = (*digiIt).getFullBX();
2573       break;
2574     }
2575   }
2576
2577   // Use these 12-bit BX's to convert emulator BX into hardware BX.
2578   if (full_anode_bx == -999) {
2579     // What to do???
2580     edm::LogWarning("L1CSCTPEmulatorWrongInput") << "+++ Warning in convertBXofLCT(): full anode BX is not available!"
2581                                                  << " +++\n";
2582   } else {
2583     // LCT BX has two bits: the least-significant bit is the LSB of ALCT BX;
2584     // the MSB is 1/0 depending on whether the 12-bit full cathode BX is 0
2585     // or not.
2586     lct_bx = (full_anode_bx + emul_bx - tbin_anode_offset) & 0x01;
2587     // SV, 12/Jun/08: it looks like this bit is never set - docu must be
2588     // wrong.
2589     //lct_bx = lct_bx | ((full_cathode_bx == 0) << 1);
2590   }
2591   //std::cout <<"convertBXofLCT old emul_bx "<< emul_bx <<" new lct bx "<< lct_bx <<" ful_anode_bx "<< full_anode_bx << std::endl;
2592   return lct_bx;
2593 }
2594
2595 void CSCTriggerPrimitivesReader::HotWires(const edm::Event& iEvent) {
2596   if (!bookedHotWireHistos)
2597     bookHotWireHistos();
2598   edm::Handle<CSCWireDigiCollection> wires;
2599   //  iEvent.getByLabel(wireDigiProducer_.label(), wireDigiProducer_.instance(), wires);
2600   iEvent.getByToken(wireDigi_token_, wires);
2601
2602   int serial_old = -1;
2603   for (auto dWDiter = wires->begin(); dWDiter != wires->end(); dWDiter++) {
2604     CSCDetId id = (CSCDetId)(*dWDiter).first;
2605     int serial = chamberSerial(id) - 1;
2606     //     printf("serial %i\n",serial);
2607     std::vector<CSCWireDigi>::const_iterator wireIter = (*dWDiter).second.first;
2608     std::vector<CSCWireDigi>::const_iterator lWire = (*dWDiter).second.second;
2609     bool has_layer = false;
2610     for (; wireIter != lWire; ++wireIter) {
2611       has_layer = true;
2612       int i_layer = id.layer() - 1;
2613       int i_wire = wireIter->getWireGroup() - 1;
2614       int nbins = wireIter->getTimeBinsOn().size();
2615       int serial2 = serial * (6 * 112) + i_layer * 112 + i_wire;
2616       /*
2617         printf("endcap %i, station %i, ring %i, chamber %i, serial %i, serial2 %i, layer %i, wg %i\n",
2618         id.endcap(), id.station(), id.ring(), id.chamber(), serial, serial2, i_layer, i_wire);
2619       */
2620       hHotWire1->Fill(serial2, nbins);
2621       /*
2622         if(id.station()==1) {
2623         printf("ring %i, chamber type %i, wg %i\n",id.ring(),id.iChamberType(),i_wire);
2624         }
2625       */
2626     }
2627     if (serial_old != serial && has_layer) {
2628       //       nHotCham[serial]++;
2629       hHotCham1->Fill(serial);
2630       //    printf("serial %i filled\n",serial);
2631       serial_old = serial;
2632     }
2633   }
2634 }
2635
2636 void CSCTriggerPrimitivesReader::MCStudies(const edm::Event& ev,
2637                                            const CSCALCTDigiCollection* alcts,
2638                                            const CSCCLCTDigiCollection* clcts) {
2639   // MC particles, if any.
2640   //edm::Handle<edm::HepMCProduct> mcp;
2641   //ev.getByLabel("source", mcp);
2642   //ev.getByType(mcp);
2643   vector<edm::Handle<edm::HepMCProduct> > allhepmcp;
2644   // Use "getManyByType" to be able to check the existence of MC info.
2645   ev.getManyByType(allhepmcp);
2646
2647   //cout << "HepMC info: " << allhepmcp.size() << endl;
2648   if (!allhepmcp.empty()) {
2649     const HepMC::GenEvent& mc = allhepmcp[0]->getHepMCData();
2650     int i = 0;
2651     for (HepMC::GenEvent::particle_const_iterator p = mc.particles_begin(); p != mc.particles_end(); ++p) {
2652       int id = (*p)->pdg_id();
2653       double phitmp = (*p)->momentum().phi();
2654       if (phitmp < 0)
2655         phitmp += 2. * M_PI;
2656       if (debug)
2657         LogDebug("CSCTriggerPrimitivesReader")
2658             << "MC part #" << ++i << ": id = " << id << ", status = " << (*p)->status()
2659             << "\n   pX = " << (*p)->momentum().x() << ", pY = " << (*p)->momentum().y()
2660             << ", pT = " << (*p)->momentum().perp() << " GeV"
2661             << ", p =  " << (*p)->momentum().rho() << " GeV"
2662             << "\n   eta = " << (*p)->momentum().pseudoRapidity() << ", phi = " << phitmp << " ("
2663             << phitmp * 180. / M_PI << " deg)";
2664     }
2665
2666     // If hepMC info is there, try to get wire and comparator digis,
2667     // and SimHits.
2668     edm::Handle<CSCWireDigiCollection> wireDigis;
2669     edm::Handle<CSCComparatorDigiCollection> compDigis;
2670     edm::Handle<edm::PSimHitContainer> simHitsH;
2671     //    ev.getByLabel(wireDigiProducer_.label(), wireDigiProducer_.instance(),
2672     //        wireDigis);
2673     //    ev.getByLabel(compDigiProducer_.label(), compDigiProducer_.instance(),
2674     //        compDigis);
2675     //    ev.getByLabel(simHitProducer_.label(), simHitProducer_.instance(),
2676     //        simHits);
2677     ev.getByToken(wireDigi_token_, wireDigis);
2678     ev.getByToken(compDigi_token_, compDigis);
2679     ev.getByToken(simHit_token_, simHitsH);
2680     const edm::PSimHitContainer& simHits = *simHitsH.product();
2681
2682     if (!wireDigis.isValid()) {
2683       edm::LogWarning("L1CSCTPEmulatorWrongInput")
2684           << "+++ Warning: Collection of wire digis with label" << wireDigiProducer_.label()
2685           << " requested, but not found in the event... Skipping the rest +++\n";
2686       return;
2687     }
2688     if (!compDigis.isValid()) {
2689       edm::LogWarning("L1CSCTPEmulatorWrongInput")
2690           << "+++ Warning: Collection of comparator digis with label" << compDigiProducer_.label()
2691           << " requested, but not found in the event... Skipping the rest +++\n";
2692       return;
2693     }
2694     if (!simHitsH.isValid()) {
2695       edm::LogWarning("L1CSCTPEmulatorWrongInput")
2696           << "+++ Warning: Collection of SimHits with label" << simHitProducer_.label()
2697           << " requested, but not found in the event... Skipping the rest +++\n";
2698       return;
2699     }
2700
2701     if (debug)
2702       LogTrace("CSCTriggerPrimitivesReader") << "   #CSC SimHits: " << simHits.size();
2703
2704     // MC-based resolution studies.
2705     calcResolution(alcts, clcts, wireDigis.product(), compDigis.product(), simHits);
2706
2707     // MC-based efficiency studies.
2708     calcEfficiency(alcts, clcts, simHits);
2709   }
2710 }
2711
2712 void CSCTriggerPrimitivesReader::calcResolution(const CSCALCTDigiCollection* alcts,
2713                                                 const CSCCLCTDigiCollection* clcts,
2714                                                 const CSCWireDigiCollection* wiredc,
2715                                                 const CSCComparatorDigiCollection* compdc,
2716                                                 const edm::PSimHitContainer& allSimHits) {
2717   // Book histos when called for the first time.
2718   if (!bookedResolHistos)
2719     bookResolHistos();
2720
2721   // ALCT resolution
2722   CSCAnodeLCTAnalyzer alct_analyzer;
2723   alct_analyzer.setGeometry(geom_);
2724
2725   for (auto adetUnitIt = alcts->begin(); adetUnitIt != alcts->end(); adetUnitIt++) {
2726     const CSCDetId& id = (*adetUnitIt).first;
2727     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
2728       continue;
2729     const auto& range = (*adetUnitIt).second;
2730     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
2731       bool alct_valid = (*digiIt).isValid();
2732       if (alct_valid) {
2733         vector<CSCAnodeLayerInfo> alctInfo = alct_analyzer.getSimInfo(*digiIt, id, wiredc, &allSimHits);
2734
2735         double hitPhi = -999.0, hitEta = -999.0;
2736         int hitWG = alct_analyzer.nearestWG(alctInfo, hitPhi, hitEta);
2737         if (hitWG >= 0.) {
2738           // Key wire group and key layer id.
2739           int wiregroup = (*digiIt).getKeyWG();
2740
2741           CSCDetId layerId(id.endcap(), id.station(), id.ring(), id.chamber(), CSCConstants::KEY_ALCT_LAYER);
2742           int endc = id.endcap();
2743           int stat = id.station();
2744           int csctype = getCSCType(id);
2745
2746           double alctEta = alct_analyzer.getWGEta(layerId, wiregroup);
2747           double deltaEta = alctEta - hitEta;
2748           hResolDeltaEta->Fill(deltaEta);
2749
2750           double deltaWG = wiregroup - hitWG;
2751           if (debug)
2752             LogTrace("CSCTriggerPrimitivesReader")
2753                 << "WG: MC = " << hitWG << " rec = " << wiregroup << " delta = " << deltaWG;
2754           hResolDeltaWG->Fill(deltaWG);
2755
2756           hEtaRecVsSim->Fill(fabs(hitEta), fabs(alctEta));
2757           hEtaDiffCsc[csctype][0]->Fill(deltaEta);
2758           hEtaDiffCsc[csctype][endc]->Fill(deltaEta);
2759           hAlctVsEta[stat - 1]->Fill(fabs(alctEta));
2760           hEtaDiffVsEta[stat - 1]->Fill(fabs(alctEta), fabs(deltaEta));
2761           hEtaDiffVsWireCsc[csctype]->Fill(wiregroup, deltaEta);
2762         }
2763         // should I comment out this "else"?
2764         //else {
2765         //  edm::LogWarning("L1CSCTPEmulatorWrongInput")
2766         //    << "+++ Warning in calcResolution(): no matched SimHit"
2767         //    << " found! +++\n";
2768         //}
2769       }
2770     }
2771   }
2772
2773   // CLCT resolution
2774   static const int key_layer = CSCConstants::KEY_CLCT_LAYER;
2775   CSCCathodeLCTAnalyzer clct_analyzer;
2776   clct_analyzer.setGeometry(geom_);
2777
2778   for (auto cdetUnitIt = clcts->begin(); cdetUnitIt != clcts->end(); cdetUnitIt++) {
2779     const CSCDetId& id = (*cdetUnitIt).first;
2780     if (checkBadChambers_ && badChambers_->isInBadChamber(id))
2781       continue;
2782     const auto& range = (*cdetUnitIt).second;
2783     for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
2784       bool clct_valid = (*digiIt).isValid();
2785       if (clct_valid) {
2786         vector<CSCCathodeLayerInfo> clctInfo = clct_analyzer.getSimInfo(*digiIt, id, compdc, &allSimHits);
2787
2788         double hitPhi = -999.0, hitEta = -999.0, deltaStrip = -999.0;
2789         int hitHS = clct_analyzer.nearestHS(clctInfo, hitPhi, hitEta);
2790         if (hitHS >= 0.) {
2791           // Key strip and key layer id.
2792           int halfstrip = (*digiIt).getKeyStrip();
2793           int strip = halfstrip / 2;
2794           int distrip = halfstrip / 4;
2795           int stripType = (*digiIt).getStripType();
2796
2797           CSCDetId layerId(id.endcap(), id.station(), id.ring(), id.chamber(), key_layer);
2798           int endc = id.endcap();
2799           int stat = id.station();
2800           int csctype = getCSCType(id);
2801
2802           // 'float strip' is in the units of 'strip', i.e., angular
2803           // widths of each strip. The counting is from 0.0 at the extreme
2804           // edge of the 'first' strip at one edge of the detector.
2805           float fstrip = -999.;
2806           if (stripType == 0) {  // di-strip CLCT
2807             fstrip = strip + 1.;
2808           } else {  // half-strip CLCT
2809             fstrip = strip + 0.5 * (halfstrip % 2) + 0.25;
2810           }
2811           double clctPhi = clct_analyzer.getStripPhi(layerId, fstrip);
2812           double deltaPhi = clctPhi - hitPhi;
2813           if (deltaPhi < -M_PI)
2814             deltaPhi += 2. * M_PI;
2815           else if (deltaPhi > M_PI)
2816             deltaPhi -= 2. * M_PI;
2817           deltaPhi *= 1000;  // in mrad
2818           if (hitPhi < 0)
2819             hitPhi += 2. * M_PI;
2820           if (clctPhi < 0)
2821             clctPhi += 2. * M_PI;
2822           if (debug)
2823             LogTrace("CSCTriggerPrimitivesReader")
2824                 << " clctPhi = " << clctPhi << " hitPhi = " << hitPhi << " deltaPhi = " << deltaPhi;
2825
2826           hResolDeltaPhi->Fill(deltaPhi);
2827           if (stripType == 0) {  // di-strip CLCT
2828             deltaStrip = distrip - hitHS / 4;
2829             hResolDeltaDS->Fill(deltaStrip);
2830             hResolDeltaPhiDS->Fill(deltaPhi);
2831             hPhiDiffVsStripCsc[csctype][0]->Fill(distrip, deltaPhi);
2832           } else {  // half-strip CLCT
2833             deltaStrip = halfstrip - hitHS;
2834             hResolDeltaHS->Fill(deltaStrip);
2835             hResolDeltaPhiHS->Fill(deltaPhi);
2836             hPhiDiffVsStripCsc[csctype][1]->Fill(halfstrip, deltaPhi);
2837           }
2838           if (debug)
2839             LogTrace("CSCTriggerPrimitivesReader") << "Half-strip: MC = " << hitHS << " rec = " << halfstrip
2840                                                    << " pattern type = " << stripType << " delta = " << deltaStrip;
2841
2842           hPhiRecVsSim->Fill(hitPhi, clctPhi);
2843           hPhiDiffCsc[csctype][0]->Fill(deltaPhi);
2844           hPhiDiffCsc[csctype][endc]->Fill(deltaPhi);
2845           hPhiDiffCsc[csctype][stripType + 3]->Fill(deltaPhi);
2846           hClctVsPhi[stat - 1]->Fill(clctPhi);
2847           hPhiDiffVsPhi[stat - 1]->Fill(clctPhi, fabs(deltaPhi));
2848
2849           // Histograms to check phi offsets for various pattern types
2850           if (stripType == 1) {                      // half-strips
2851             double hsperrad = getHsPerRad(csctype);  // halfstrips-per-radian
2852             if ((endc == 1 && (stat == 1 || stat == 2)) || (endc == 2 && (stat == 3 || stat == 4))) {
2853               hPhiDiffPattern[(*digiIt).getPattern()]->Fill(deltaPhi / 1000 * hsperrad);
2854             }
2855           }
2856         }
2857         // should I comment out this "else"?
2858         //else {
2859         //  edm::LogWarning("L1CSCTPEmulatorWrongInput")
2860         //    << "+++ Warning in calcResolution(): no matched SimHit"
2861         //    << " found! +++\n";
2862         //}
2863
2864         // "True bend", defined as difference in phi between muon hit
2865         // positions in the first and in the sixth layers.
2866         double phi1 = -999.0, phi6 = -999.0;
2867         vector<CSCCathodeLayerInfo>::const_iterator pli;
2868         for (pli = clctInfo.begin(); pli != clctInfo.end(); pli++) {
2869           CSCDetId layerId = pli->getId();
2870           int layer = layerId.layer();
2871           if (layer == 1 || layer == 6) {
2872             // Get simHits in this layer.
2873             for (const auto& psh : allSimHits) {
2874               // Find detId where simHit is located.
2875               CSCDetId hitId = (CSCDetId)psh.detUnitId();
2876               if (hitId == layerId && abs(psh.particleType()) == 13) {  // muon hits only
2877                 const CSCLayer* csclayer = geom_->layer(layerId);
2878                 const auto& thisPoint = csclayer->toGlobal(psh.localPosition());
2879                 double phi = thisPoint.phi();
2880                 if (layer == 1)
2881                   phi1 = phi;
2882                 else if (layer == 6)
2883                   phi6 = phi;
2884                 break;  // simply take the first suitable hit.
2885               }
2886             }
2887           }
2888         }
2889         if (phi1 > -99. && phi6 > -99.) {
2890           double deltaPhi = phi1 - phi6;
2891           if (deltaPhi > M_PI)
2892             deltaPhi -= 2. * M_PI;
2893           else if (deltaPhi < -M_PI)
2894             deltaPhi += 2. * M_PI;
2895           int csctype = getCSCType(id);
2896           hTrueBendCsc[csctype]->Fill(deltaPhi * 1000.);  // in mrad
2897         }
2898       }
2899     }
2900   }
2901 }
2902
2903 void CSCTriggerPrimitivesReader::calcEfficiency(const CSCALCTDigiCollection* alcts,
2904                                                 const CSCCLCTDigiCollection* clcts,
2905                                                 const edm::PSimHitContainer& allSimHits) {
2906   // Book histos when called for the first time.
2907   if (!bookedEfficHistos)
2908     bookEfficHistos();
2909
2910   // Create list of chambers having SimHits.
2911   vector<CSCDetId> chamberIds;
2912   vector<CSCDetId>::const_iterator chamberIdIt;
2913   for (const auto& simHitIt : allSimHits) {
2914     // Find detId where simHit is located.
2915     bool sameId = false;
2916     CSCDetId hitId = (CSCDetId)simHitIt.detUnitId();
2917     // Skip chambers marked as bad (includes most of ME4/2 chambers).
2918     if (checkBadChambers_ && badChambers_->isInBadChamber(hitId))
2919       continue;
2920     //if (hitId.ring() == 4) continue; // skip ME1/A for now.
2921     if (!plotME1A && hitId.ring() == 4)
2922       continue;
2923     for (const auto& chamberId : chamberIds) {
2924       if (chamberId.endcap() == hitId.endcap() && chamberId.station() == hitId.station() &&
2925           chamberId.ring() == hitId.ring() && chamberId.chamber() == hitId.chamber()) {
2926         sameId = true;
2927         break;
2928       }
2929     }
2930     if (!sameId) {
2931       CSCDetId newChamberId(hitId.endcap(), hitId.station(), hitId.ring(), hitId.chamber(), 0);
2932       chamberIds.push_back(newChamberId);
2933     }
2934   }
2935   LogTrace("CSCTriggerPrimitivesReader") << "Found SimHits in " << chamberIds.size() << " CSCs";
2936
2937   bool used[CSCConstants::NUM_LAYERS];
2938   vector<PSimHit> simHitsV[CSCConstants::NUM_LAYERS];
2939   for (const auto& chamberId : chamberIds) {
2940     // Find out how many layers of this chamber have SimHits.
2941     int nLayers = 0;
2942     for (int ilayer = 0; ilayer < CSCConstants::NUM_LAYERS; ilayer++) {
2943       used[ilayer] = false;
2944       simHitsV[ilayer].clear();
2945     }
2946
2947     int endcap = chamberId.endcap();
2948     int station = chamberId.station();
2949     int ring = chamberId.ring();
2950     int chamber = chamberId.chamber();
2951     for (const auto& simHitIt : allSimHits) {
2952       CSCDetId hitId = (CSCDetId)simHitIt.detUnitId();
2953       if (hitId.endcap() == endcap && hitId.station() == station && hitId.ring() == ring &&
2954           hitId.chamber() == chamber) {
2955         int layer = hitId.layer() - 1;
2956         if (!used[layer] && abs(simHitIt.particleType()) == 13) {
2957           nLayers++;
2958           used[layer] = true;
2959           simHitsV[layer].push_back(simHitIt);
2960         }
2961       }
2962     }
2963     LogTrace("CSCTriggerPrimitivesReader") << "CSC in ME" << ((endcap == 1) ? "+" : "-") << station << "/" << ring
2964                                            << "/" << chamber << " has muon hits in " << nLayers << " layers";
2965
2966     // If the number of layers with hits is above threshold, look for
2967     // a presence of LCTs.
2968     if (nLayers > 3) {  // Should be a parameter.
2969       // Start with the key layer and take the eta of the first hit.
2970       // Really crude; should be improved.
2971       double hitEta = -999.;
2972       for (int ilayer = 2; ilayer < CSCConstants::NUM_LAYERS; ilayer++) {
2973         vector<PSimHit> layerSimHitsV = simHitsV[ilayer];
2974         if (!layerSimHitsV.empty()) {
2975           LocalPoint hitLP = layerSimHitsV[0].localPosition();
2976           CSCDetId layerId = (CSCDetId)(layerSimHitsV[0]).detUnitId();
2977           const CSCLayer* csclayer = geom_->layer(layerId);
2978           GlobalPoint hitGP = csclayer->toGlobal(hitLP);
2979           hitEta = hitGP.eta();
2980           break;
2981         }
2982       }
2983       if (hitEta < -3.) {
2984         edm::LogWarning("L1CSCTPEmulatorWrongValues") << "+++ Warning in calcEfficiency(): no SimHit found"
2985                                                       << " where there must be at least " << nLayers << "! +++\n";
2986         continue;
2987       }
2988       int csctype = getCSCType(*chamberIdIt);
2989       hEfficHitsEta[station - 1]->Fill(fabs(hitEta));
2990       hEfficHitsEtaCsc[csctype]->Fill(fabs(hitEta));
2991
2992       bool isALCT = false;
2993       for (auto adetUnitIt = alcts->begin(); adetUnitIt != alcts->end(); adetUnitIt++) {
2994         const CSCDetId& id = (*adetUnitIt).first;
2995         if (id == (*chamberIdIt)) {
2996           const auto& range = (*adetUnitIt).second;
2997           for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
2998             if (digiIt->isValid()) {
2999               // Check the distance??
3000               LogTrace("CSCTriggerPrimitivesReader") << "ALCT was found";
3001               isALCT = true;
3002               break;
3003             }
3004           }
3005         }
3006         if (isALCT)
3007           break;
3008       }
3009       if (isALCT) {
3010         hEfficALCTEta[station - 1]->Fill(fabs(hitEta));
3011         hEfficALCTEtaCsc[csctype]->Fill(fabs(hitEta));
3012       } else {
3013         LogTrace("CSCTriggerPrimitivesReader") << "ALCT was not found";
3014       }
3015
3016       bool isCLCT = false;
3017       for (auto cdetUnitIt = clcts->begin(); cdetUnitIt != clcts->end(); cdetUnitIt++) {
3018         const CSCDetId& id = (*cdetUnitIt).first;
3019         if (id == (*chamberIdIt)) {
3020           const auto& range = (*cdetUnitIt).second;
3021           for (auto digiIt = range.first; digiIt != range.second; digiIt++) {
3022             if (digiIt->isValid()) {
3023               // Check the distance??
3024               LogTrace("CSCTriggerPrimitivesReader") << "CLCT was found";
3025               isCLCT = true;
3026               break;
3027             }
3028           }
3029         }
3030         if (isCLCT)
3031           break;
3032       }
3033       if (isCLCT) {
3034         hEfficCLCTEta[station - 1]->Fill(fabs(hitEta));
3035         hEfficCLCTEtaCsc[csctype]->Fill(fabs(hitEta));
3036       } else {
3037         LogTrace("CSCTriggerPrimitivesReader") << "CLCT was not found";
3038       }
3039     }
3040   }
3041 }
3042
3043 void CSCTriggerPrimitivesReader::drawALCTHistos() {
3044   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3045   string fname = resultsFileNamesPrefix_ + "alcts.ps";
3046   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3047
3048   TPad* pad[MAXPAGES];
3049   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3050     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3051   }
3052
3053   int page = 1;
3054   TText t;
3055   t.SetTextFont(32);
3056   t.SetTextSize(0.025);
3057   char pagenum[7], titl[50];
3058   TPaveLabel* title;
3059
3060   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
3061
3062   ps->NewPage();
3063   c1->Clear();
3064   c1->cd(0);
3065   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "Number of ALCTs");
3066   title->SetFillColor(10);
3067   title->Draw();
3068   sprintf(pagenum, "- %d -", page);
3069   t.DrawText(0.9, 0.02, pagenum);
3070   gStyle->SetOptStat(111110);
3071   pad[page]->Draw();
3072   pad[page]->Divide(1, 3);
3073   pad[page]->cd(1);
3074   hAlctPerEvent->Draw();
3075   pad[page]->cd(2);
3076   hAlctPerChamber->Draw();
3077   for (int i = 0; i < CSC_TYPES; i++) {
3078     hAlctPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3079   }
3080   // Should be multiplied by 40/nevents to convert to MHz
3081   pad[page]->cd(3);
3082   hAlctPerCSC->Draw();
3083   page++;
3084   c1->Update();
3085
3086   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3087     ps->NewPage();
3088     c1->Clear();
3089     c1->cd(0);
3090     sprintf(titl, "ALCTs per chamber, endcap %d", endc + 1);
3091     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3092     title->SetFillColor(10);
3093     title->Draw();
3094     sprintf(pagenum, "- %d -", page);
3095     t.DrawText(0.9, 0.02, pagenum);
3096     gStyle->SetOptStat(10);
3097     pad[page]->Draw();
3098     pad[page]->Divide(2, 5);
3099     for (int idh = 0; idh < max_idh; idh++) {
3100       if (!plotME1A && idh == 3)
3101         continue;
3102       hAlctCsc[endc][idh]->SetMinimum(0.0);
3103       pad[page]->cd(idh + 1);
3104       hAlctCsc[endc][idh]->Draw();
3105     }
3106     page++;
3107     c1->Update();
3108   }
3109
3110   ps->NewPage();
3111   c1->Clear();
3112   c1->cd(0);
3113   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "ALCT quantities");
3114   title->SetFillColor(10);
3115   title->Draw();
3116   sprintf(pagenum, "- %d -", page);
3117   t.DrawText(0.9, 0.02, pagenum);
3118   gStyle->SetOptStat(111110);
3119   pad[page]->Draw();
3120   pad[page]->Divide(2, 3);
3121   pad[page]->cd(1);
3122   hAlctValid->Draw();
3123   pad[page]->cd(2);
3124   hAlctQuality->Draw();
3125   pad[page]->cd(3);
3126   hAlctAccel->Draw();
3127   pad[page]->cd(4);
3128   hAlctCollis->Draw();
3129   pad[page]->cd(5);
3130   hAlctKeyGroup->Draw();
3131   pad[page]->cd(6);
3132   hAlctBXN->Draw();
3133   page++;
3134   c1->Update();
3135
3136   ps->Close();
3137   delete c1;
3138 }
3139
3140 void CSCTriggerPrimitivesReader::drawCLCTHistos() {
3141   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3142   string fname = resultsFileNamesPrefix_ + "clcts.ps";
3143   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3144
3145   TPad* pad[MAXPAGES];
3146   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3147     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3148   }
3149
3150   int page = 1;
3151   TText t;
3152   t.SetTextFont(32);
3153   t.SetTextSize(0.025);
3154   char pagenum[7], titl[50];
3155   TPaveLabel* title;
3156
3157   ps->NewPage();
3158   c1->Clear();
3159   c1->cd(0);
3160   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "Number of CLCTs");
3161   title->SetFillColor(10);
3162   title->Draw();
3163   sprintf(pagenum, "- %d -", page);
3164   t.DrawText(0.9, 0.02, pagenum);
3165   gStyle->SetOptStat(111110);
3166   pad[page]->Draw();
3167   pad[page]->Divide(1, 3);
3168   pad[page]->cd(1);
3169   hClctPerEvent->Draw();
3170
3171   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
3172
3173   edm::LogInfo("CSCTriggerPrimitivesReader") << "\n";
3174   int nbins = hClctPerChamber->GetNbinsX();
3175   for (int ibin = 1; ibin <= nbins; ibin++) {
3176     double f_bin = hClctPerChamber->GetBinContent(ibin);
3177     edm::LogInfo("CSCTriggerPrimitivesReader") << "  # CLCTs/chamber: " << ibin - 1 << "; events: " << f_bin << endl;
3178   }
3179
3180   pad[page]->cd(2);
3181   hClctPerChamber->Draw();
3182   for (int i = 0; i < CSC_TYPES; i++) {
3183     hClctPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3184   }
3185   // Should be multiplied by 40/nevents to convert to MHz
3186   pad[page]->cd(3);
3187   hClctPerCSC->Draw();
3188   page++;
3189   c1->Update();
3190
3191   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3192     ps->NewPage();
3193     c1->Clear();
3194     c1->cd(0);
3195     sprintf(titl, "CLCTs per chamber, endcap %d", endc + 1);
3196     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3197     title->SetFillColor(10);
3198     title->Draw();
3199     sprintf(pagenum, "- %d -", page);
3200     t.DrawText(0.9, 0.02, pagenum);
3201     gStyle->SetOptStat(10);
3202     pad[page]->Draw();
3203     pad[page]->Divide(2, 5);
3204     for (int idh = 0; idh < max_idh; idh++) {
3205       if (!plotME1A && idh == 3)
3206         continue;
3207       hClctCsc[endc][idh]->SetMinimum(0.0);
3208       pad[page]->cd(idh + 1);
3209       hClctCsc[endc][idh]->Draw();
3210     }
3211     page++;
3212     c1->Update();
3213   }
3214
3215   ps->NewPage();
3216   c1->Clear();
3217   c1->cd(0);
3218   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT quantities");
3219   title->SetFillColor(10);
3220   title->Draw();
3221   sprintf(pagenum, "- %d -", page);
3222   t.DrawText(0.9, 0.02, pagenum);
3223   gStyle->SetOptStat(111110);
3224   pad[page]->Draw();
3225   pad[page]->Divide(2, 4);
3226   pad[page]->cd(1);
3227   hClctValid->Draw();
3228   pad[page]->cd(2);
3229   hClctQuality->Draw();
3230   pad[page]->cd(3);
3231   hClctSign->Draw();
3232   hClctPattern[1]->SetTitle("CLCT pattern #");
3233   pad[page]->cd(4);
3234   hClctPattern[1]->Draw();
3235   pad[page]->cd(5);
3236   hClctCFEB->Draw();
3237   pad[page]->cd(7);
3238   hClctKeyStrip[1]->Draw();
3239   pad[page]->cd(8);
3240   hClctBXN->Draw();
3241   page++;
3242   c1->Update();
3243
3244   ps->NewPage();
3245   c1->Clear();
3246   c1->cd(0);
3247   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT bend for various chamber types, halfstrips");
3248   title->SetFillColor(10);
3249   title->Draw();
3250   sprintf(pagenum, "- %d -", page);
3251   t.DrawText(0.9, 0.02, pagenum);
3252   gStyle->SetOptStat(110);
3253   pad[page]->Draw();
3254   pad[page]->Divide(2, 5);
3255   for (int idh = 0; idh < max_idh; idh++) {
3256     if (!plotME1A && idh == 3)
3257       continue;
3258     pad[page]->cd(idh + 1);
3259     hClctBendCsc[idh][1]->GetXaxis()->SetTitle("Pattern bend");
3260     hClctBendCsc[idh][1]->GetYaxis()->SetTitle("Number of LCTs");
3261     hClctBendCsc[idh][1]->Draw();
3262   }
3263   page++;
3264   c1->Update();
3265
3266   ps->Close();
3267   delete c1;
3268 }
3269
3270 void CSCTriggerPrimitivesReader::drawLCTTMBHistos() {
3271   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3272   string fname = resultsFileNamesPrefix_ + "lcts_tmb.ps";
3273   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3274
3275   TPad* pad[MAXPAGES];
3276   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3277     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3278   }
3279
3280   int page = 1;
3281   TText t;
3282   t.SetTextFont(32);
3283   t.SetTextSize(0.025);
3284   char pagenum[7], titl[50];
3285   TPaveLabel* title;
3286
3287   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
3288
3289   ps->NewPage();
3290   c1->Clear();
3291   c1->cd(0);
3292   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "Number of LCTs");
3293   title->SetFillColor(10);
3294   title->Draw();
3295   sprintf(pagenum, "- %d -", page);
3296   t.DrawText(0.9, 0.02, pagenum);
3297   gStyle->SetOptStat(111110);
3298   pad[page]->Draw();
3299   pad[page]->Divide(2, 2);
3300   pad[page]->cd(1);
3301   hLctTMBPerEvent->Draw();
3302   pad[page]->cd(2);
3303   hLctTMBPerChamber->Draw();
3304   c1->Update();
3305   for (int i = 0; i < CSC_TYPES; i++) {
3306     hLctTMBPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3307     hCorrLctTMBPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3308   }
3309   // Should be multiplied by 40/nevents to convert to MHz
3310   pad[page]->cd(3);
3311   hLctTMBPerCSC->Draw();
3312   pad[page]->cd(4);
3313   hCorrLctTMBPerCSC->Draw();
3314   gStyle->SetOptStat(1110);
3315   page++;
3316   c1->Update();
3317
3318   ps->NewPage();
3319   c1->Clear();
3320   c1->cd(0);
3321   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "LCT geometry");
3322   title->SetFillColor(10);
3323   title->Draw();
3324   sprintf(pagenum, "- %d -", page);
3325   t.DrawText(0.9, 0.02, pagenum);
3326   gStyle->SetOptStat(110110);
3327   pad[page]->Draw();
3328   pad[page]->Divide(2, 4);
3329   pad[page]->cd(1);
3330   hLctTMBEndcap->Draw();
3331   pad[page]->cd(2);
3332   hLctTMBStation->Draw();
3333   pad[page]->cd(3);
3334   hLctTMBSector->Draw();
3335   pad[page]->cd(4);
3336   hLctTMBRing->Draw();
3337   for (int istat = 0; istat < MAX_STATIONS; istat++) {
3338     pad[page]->cd(istat + 5);
3339     hLctTMBChamber[istat]->Draw();
3340   }
3341   page++;
3342   c1->Update();
3343
3344   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3345     ps->NewPage();
3346     c1->Clear();
3347     c1->cd(0);
3348     sprintf(titl, "LCTs per chamber, endcap %d", endc + 1);
3349     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3350     title->SetFillColor(10);
3351     title->Draw();
3352     sprintf(pagenum, "- %d -", page);
3353     t.DrawText(0.9, 0.02, pagenum);
3354     gStyle->SetOptStat(10);
3355     pad[page]->Draw();
3356     pad[page]->Divide(2, 5);
3357     for (int idh = 0; idh < max_idh; idh++) {
3358       if (!plotME1A && idh == 3)
3359         continue;
3360       hLctTMBCsc[endc][idh]->SetMinimum(0.0);
3361       pad[page]->cd(idh + 1);
3362       hLctTMBCsc[endc][idh]->Draw();
3363     }
3364     page++;
3365     c1->Update();
3366   }
3367
3368   ps->NewPage();
3369   c1->Clear();
3370   c1->cd(0);
3371   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "LCT quantities");
3372   title->SetFillColor(10);
3373   title->Draw();
3374   sprintf(pagenum, "- %d -", page);
3375   t.DrawText(0.9, 0.02, pagenum);
3376   gStyle->SetOptStat(110110);
3377   pad[page]->Draw();
3378   pad[page]->Divide(2, 4);
3379   pad[page]->cd(1);
3380   hLctTMBValid->Draw();
3381   pad[page]->cd(2);
3382   hLctTMBQuality->Draw();
3383   pad[page]->cd(3);
3384   hLctTMBKeyGroup->Draw();
3385   pad[page]->cd(4);
3386   hLctTMBKeyStrip->Draw();
3387   pad[page]->cd(5);
3388   hLctTMBStripType->Draw();
3389   pad[page]->cd(6);
3390   hLctTMBPattern->Draw();
3391   pad[page]->cd(7);
3392   hLctTMBBend->Draw();
3393   pad[page]->cd(8);
3394   hLctTMBBXN->Draw();
3395   page++;
3396   c1->Update();
3397
3398   ps->Close();
3399   delete c1;
3400 }
3401
3402 void CSCTriggerPrimitivesReader::drawLCTMPCHistos() {
3403   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3404   string fname = resultsFileNamesPrefix_ + "lcts_mpc.ps";
3405   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3406
3407   TPad* pad[MAXPAGES];
3408   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3409     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3410   }
3411
3412   int page = 1;
3413   TText t;
3414   t.SetTextFont(32);
3415   t.SetTextSize(0.025);
3416   char pagenum[7];
3417   TPaveLabel* title;
3418
3419   ps->NewPage();
3420   c1->Clear();
3421   c1->cd(0);
3422   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "Number of LCTs");
3423   title->SetFillColor(10);
3424   title->Draw();
3425   sprintf(pagenum, "- %d -", page);
3426   t.DrawText(0.9, 0.02, pagenum);
3427   gStyle->SetOptStat(111110);
3428   pad[page]->Draw();
3429   pad[page]->Divide(1, 3);
3430   pad[page]->cd(1);
3431   hLctMPCPerEvent->Draw();
3432   for (int i = 0; i < CSC_TYPES; i++) {
3433     hLctMPCPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3434     hCorrLctMPCPerCSC->GetXaxis()->SetBinLabel(i + 1, csc_type[i].c_str());
3435   }
3436   // Should be multiplied by 40/nevents to convert to MHz
3437   pad[page]->cd(2);
3438   hLctMPCPerCSC->Draw();
3439   pad[page]->cd(3);
3440   hCorrLctMPCPerCSC->Draw();
3441   page++;
3442   c1->Update();
3443
3444   ps->NewPage();
3445   c1->Clear();
3446   c1->cd(0);
3447   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "LCT geometry");
3448   title->SetFillColor(10);
3449   title->Draw();
3450   sprintf(pagenum, "- %d -", page);
3451   t.DrawText(0.9, 0.02, pagenum);
3452   gStyle->SetOptStat(110110);
3453   pad[page]->Draw();
3454   pad[page]->Divide(2, 4);
3455   pad[page]->cd(1);
3456   hLctMPCEndcap->Draw();
3457   pad[page]->cd(2);
3458   hLctMPCStation->Draw();
3459   pad[page]->cd(3);
3460   hLctMPCSector->Draw();
3461   pad[page]->cd(4);
3462   hLctMPCRing->Draw();
3463   for (int istat = 0; istat < MAX_STATIONS; istat++) {
3464     pad[page]->cd(istat + 5);
3465     hLctMPCChamber[istat]->Draw();
3466   }
3467   page++;
3468   c1->Update();
3469
3470   ps->NewPage();
3471   c1->Clear();
3472   c1->cd(0);
3473   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "LCT quantities");
3474   title->SetFillColor(10);
3475   title->Draw();
3476   sprintf(pagenum, "- %d -", page);
3477   t.DrawText(0.9, 0.02, pagenum);
3478   pad[page]->Draw();
3479   pad[page]->Divide(2, 4);
3480   pad[page]->cd(1);
3481   hLctMPCValid->Draw();
3482   pad[page]->cd(2);
3483   hLctMPCQuality->Draw();
3484   pad[page]->cd(3);
3485   hLctMPCKeyGroup->Draw();
3486   pad[page]->cd(4);
3487   hLctMPCKeyStrip->Draw();
3488   pad[page]->cd(5);
3489   hLctMPCStripType->Draw();
3490   pad[page]->cd(6);
3491   hLctMPCPattern->Draw();
3492   pad[page]->cd(7);
3493   hLctMPCBend->Draw();
3494   pad[page]->cd(8);
3495   hLctMPCBXN->Draw();
3496   page++;
3497   c1->Update();
3498
3499   ps->Close();
3500   delete c1;
3501 }
3502
3503 void CSCTriggerPrimitivesReader::drawCompHistos() {
3504   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3505   string fname = resultsFileNamesPrefix_ + "lcts_comp.ps";
3506   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3507
3508   TPad* pad[MAXPAGES];
3509   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3510     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3511   }
3512
3513   int page = 1;
3514   TText t;
3515   t.SetTextFont(32);
3516   t.SetTextSize(0.025);
3517   char pagenum[15];
3518   TPaveLabel* title;
3519   Int_t nbins;
3520
3521   TText teff;
3522   teff.SetTextFont(32);
3523   teff.SetTextSize(0.08);
3524   char eff[25], titl[60];
3525
3526   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
3527
3528   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {  // endcaps
3529     ps->NewPage();
3530     c1->Clear();
3531     c1->cd(0);
3532     sprintf(titl, "ALCT firmware-emulator: match in number found, endcap %d", endc + 1);
3533     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3534     title->SetFillColor(10);
3535     title->Draw();
3536     sprintf(pagenum, "- %d -", page);
3537     t.DrawText(0.9, 0.02, pagenum);
3538     //gStyle->SetOptStat(110010);
3539     gStyle->SetOptStat(0);
3540     pad[page]->Draw();
3541     pad[page]->Divide(2, 5);
3542     TH1F* hAlctFoundEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3543     for (int idh = 0; idh < max_idh; idh++) {
3544       if (!plotME1A && idh == 3)
3545         continue;
3546       hAlctFoundEffVsCsc[endc][idh] = (TH1F*)hAlctCompFoundCsc[endc][idh]->Clone();
3547       hAlctFoundEffVsCsc[endc][idh]->Divide(hAlctCompSameNCsc[endc][idh], hAlctCompFoundCsc[endc][idh], 1., 1., "B");
3548       nbins = hAlctCompFoundCsc[endc][idh]->GetNbinsX();
3549       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3550         if (hAlctCompFoundCsc[endc][idh]->GetBinContent(ibin) == 0) {
3551           hAlctFoundEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3552           hAlctFoundEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3553         }
3554       }
3555       gPad->Update();
3556       gStyle->SetStatX(0.65);
3557       hAlctFoundEffVsCsc[endc][idh]->SetMinimum(-0.05);
3558       hAlctFoundEffVsCsc[endc][idh]->SetMaximum(1.05);
3559       hAlctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3560       hAlctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3561       hAlctFoundEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);  // default=0.04
3562       hAlctFoundEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);  // default=0.04
3563       hAlctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);  // default=0.05
3564       hAlctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);  // default=0.05
3565       hAlctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3566       hAlctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of same number found");
3567       pad[page]->cd(idh + 1);
3568       hAlctFoundEffVsCsc[endc][idh]->Draw("e");
3569       double numer = hAlctCompSameNCsc[endc][idh]->Integral();
3570       double denom = hAlctCompFoundCsc[endc][idh]->Integral();
3571       double ratio = 0.0, error = 0.0;
3572       if (denom > 0.) {
3573         ratio = numer / denom;
3574         error = sqrt(ratio * (1. - ratio) / denom);
3575       }
3576       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3577       teff.DrawTextNDC(0.3, 0.5, eff);
3578     }
3579     page++;
3580     c1->Update();
3581   }
3582
3583   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {  // endcaps
3584     ps->NewPage();
3585     c1->Clear();
3586     c1->cd(0);
3587     sprintf(titl, "ALCT firmware-emulator: exact match, endcap %d", endc + 1);
3588     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3589     title->SetFillColor(10);
3590     title->Draw();
3591     sprintf(pagenum, "- %d -", page);
3592     t.DrawText(0.9, 0.02, pagenum);
3593     //gStyle->SetOptStat(110010);
3594     gStyle->SetOptStat(0);
3595     pad[page]->Draw();
3596     pad[page]->Divide(2, 5);
3597     TH1F* hAlctMatchEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3598     for (int idh = 0; idh < max_idh; idh++) {
3599       if (!plotME1A && idh == 3)
3600         continue;
3601       hAlctMatchEffVsCsc[endc][idh] = (TH1F*)hAlctCompTotalCsc[endc][idh]->Clone();
3602       hAlctMatchEffVsCsc[endc][idh]->Divide(hAlctCompMatchCsc[endc][idh], hAlctCompTotalCsc[endc][idh], 1., 1., "B");
3603       nbins = hAlctCompTotalCsc[endc][idh]->GetNbinsX();
3604       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3605         if (hAlctCompTotalCsc[endc][idh]->GetBinContent(ibin) == 0) {
3606           hAlctMatchEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3607           hAlctMatchEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3608         }
3609       }
3610       gPad->Update();
3611       gStyle->SetStatX(0.65);
3612       hAlctMatchEffVsCsc[endc][idh]->SetMinimum(-0.05);
3613       hAlctMatchEffVsCsc[endc][idh]->SetMaximum(1.05);
3614       hAlctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3615       hAlctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3616       hAlctMatchEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);
3617       hAlctMatchEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);
3618       hAlctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);
3619       hAlctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);
3620       hAlctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3621       hAlctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of exact match");
3622       pad[page]->cd(idh + 1);
3623       hAlctMatchEffVsCsc[endc][idh]->Draw("e");
3624       double numer = hAlctCompMatchCsc[endc][idh]->Integral();
3625       double denom = hAlctCompTotalCsc[endc][idh]->Integral();
3626       double ratio = 0.0, error = 0.0;
3627       if (denom > 0.) {
3628         ratio = numer / denom;
3629         error = sqrt(ratio * (1. - ratio) / denom);
3630       }
3631       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3632       teff.DrawTextNDC(0.3, 0.5, eff);
3633     }
3634     page++;
3635     c1->Update();
3636   }
3637
3638   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3639     ps->NewPage();
3640     c1->Clear();
3641     c1->cd(0);
3642     sprintf(titl, "CLCT firmware-emulator: match in number found, endcap %d", endc + 1);
3643     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3644     title->SetFillColor(10);
3645     title->Draw();
3646     sprintf(pagenum, "- %d -", page);
3647     t.DrawText(0.9, 0.02, pagenum);
3648     //gStyle->SetOptStat(110010);
3649     gStyle->SetOptStat(0);
3650     pad[page]->Draw();
3651     pad[page]->Divide(2, 5);
3652     TH1F* hClctFoundEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3653     for (int idh = 0; idh < max_idh; idh++) {
3654       if (!plotME1A && idh == 3)
3655         continue;
3656       hClctFoundEffVsCsc[endc][idh] = (TH1F*)hClctCompFoundCsc[endc][idh]->Clone();
3657       hClctFoundEffVsCsc[endc][idh]->Divide(hClctCompSameNCsc[endc][idh], hClctCompFoundCsc[endc][idh], 1., 1., "B");
3658       nbins = hClctCompFoundCsc[endc][idh]->GetNbinsX();
3659       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3660         if (hClctCompFoundCsc[endc][idh]->GetBinContent(ibin) == 0) {
3661           hClctFoundEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3662           hClctFoundEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3663         }
3664       }
3665       gPad->Update();
3666       gStyle->SetStatX(0.65);
3667       hClctFoundEffVsCsc[endc][idh]->SetMinimum(-0.05);
3668       hClctFoundEffVsCsc[endc][idh]->SetMaximum(1.05);
3669       hClctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3670       hClctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3671       hClctFoundEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);
3672       hClctFoundEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);
3673       hClctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);
3674       hClctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);
3675       hClctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3676       hClctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of same number found");
3677       pad[page]->cd(idh + 1);
3678       hClctFoundEffVsCsc[endc][idh]->Draw("e");
3679       double numer = hClctCompSameNCsc[endc][idh]->Integral();
3680       double denom = hClctCompFoundCsc[endc][idh]->Integral();
3681       double ratio = 0.0, error = 0.0;
3682       if (denom > 0.) {
3683         ratio = numer / denom;
3684         error = sqrt(ratio * (1. - ratio) / denom);
3685       }
3686       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3687       teff.DrawTextNDC(0.3, 0.5, eff);
3688     }
3689     page++;
3690     c1->Update();
3691   }
3692
3693   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3694     ps->NewPage();
3695     c1->Clear();
3696     c1->cd(0);
3697     sprintf(titl, "CLCT firmware-emulator: exact match, endcap %d", endc + 1);
3698     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3699     title->SetFillColor(10);
3700     title->Draw();
3701     sprintf(pagenum, "- %d -", page);
3702     t.DrawText(0.9, 0.02, pagenum);
3703     //gStyle->SetOptStat(110010);
3704     gStyle->SetOptStat(0);
3705     pad[page]->Draw();
3706     pad[page]->Divide(2, 5);
3707     TH1F* hClctMatchEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3708     for (int idh = 0; idh < max_idh; idh++) {
3709       if (!plotME1A && idh == 3)
3710         continue;
3711       hClctMatchEffVsCsc[endc][idh] = (TH1F*)hClctCompTotalCsc[endc][idh]->Clone();
3712       hClctMatchEffVsCsc[endc][idh]->Divide(hClctCompMatchCsc[endc][idh], hClctCompTotalCsc[endc][idh], 1., 1., "B");
3713       nbins = hClctCompTotalCsc[endc][idh]->GetNbinsX();
3714       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3715         if (hClctCompTotalCsc[endc][idh]->GetBinContent(ibin) == 0) {
3716           hClctMatchEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3717           hClctMatchEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3718         }
3719       }
3720       gPad->Update();
3721       gStyle->SetStatX(0.65);
3722       hClctMatchEffVsCsc[endc][idh]->SetMinimum(-0.05);
3723       hClctMatchEffVsCsc[endc][idh]->SetMaximum(1.05);
3724       hClctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3725       hClctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3726       hClctMatchEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);
3727       hClctMatchEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);
3728       hClctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);
3729       hClctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);
3730       hClctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3731       hClctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of exact match");
3732       pad[page]->cd(idh + 1);
3733       hClctMatchEffVsCsc[endc][idh]->Draw("e");
3734       double numer = hClctCompMatchCsc[endc][idh]->Integral();
3735       double denom = hClctCompTotalCsc[endc][idh]->Integral();
3736       double ratio = 0.0, error = 0.0;
3737       if (denom > 0.) {
3738         ratio = numer / denom;
3739         error = sqrt(ratio * (1. - ratio) / denom);
3740       }
3741       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3742       teff.DrawTextNDC(0.3, 0.5, eff);
3743     }
3744     page++;
3745     c1->Update();
3746   }
3747
3748   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3749     ps->NewPage();
3750     c1->Clear();
3751     c1->cd(0);
3752     sprintf(titl, "LCT firmware-emulator: match in number found, endcap %d", endc + 1);
3753     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3754     title->SetFillColor(10);
3755     title->Draw();
3756     sprintf(pagenum, "- %d -", page);
3757     t.DrawText(0.9, 0.02, pagenum);
3758     //gStyle->SetOptStat(110010);
3759     gStyle->SetOptStat(0);
3760     pad[page]->Draw();
3761     pad[page]->Divide(2, 5);
3762     TH1F* hLctFoundEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3763     for (int idh = 0; idh < max_idh; idh++) {
3764       if (!plotME1A && idh == 3)
3765         continue;
3766       hLctFoundEffVsCsc[endc][idh] = (TH1F*)hLctCompFoundCsc[endc][idh]->Clone();
3767       hLctFoundEffVsCsc[endc][idh]->Divide(hLctCompSameNCsc[endc][idh], hLctCompFoundCsc[endc][idh], 1., 1., "B");
3768       nbins = hLctCompFoundCsc[endc][idh]->GetNbinsX();
3769       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3770         if (hLctCompFoundCsc[endc][idh]->GetBinContent(ibin) == 0) {
3771           hLctFoundEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3772           hLctFoundEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3773         }
3774       }
3775       gPad->Update();
3776       gStyle->SetStatX(0.65);
3777       hLctFoundEffVsCsc[endc][idh]->SetMinimum(-0.05);
3778       hLctFoundEffVsCsc[endc][idh]->SetMaximum(1.05);
3779       hLctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3780       hLctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3781       hLctFoundEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);
3782       hLctFoundEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);
3783       hLctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);
3784       hLctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);
3785       hLctFoundEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3786       hLctFoundEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of same number found");
3787       pad[page]->cd(idh + 1);
3788       hLctFoundEffVsCsc[endc][idh]->Draw("e");
3789       double numer = hLctCompSameNCsc[endc][idh]->Integral();
3790       double denom = hLctCompFoundCsc[endc][idh]->Integral();
3791       double ratio = 0.0, error = 0.0;
3792       if (denom > 0.) {
3793         ratio = numer / denom;
3794         error = sqrt(ratio * (1. - ratio) / denom);
3795       }
3796       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3797       teff.DrawTextNDC(0.3, 0.5, eff);
3798     }
3799     page++;
3800     c1->Update();
3801   }
3802
3803   for (int endc = 0; endc < MAX_ENDCAPS; endc++) {
3804     ps->NewPage();
3805     c1->Clear();
3806     c1->cd(0);
3807     sprintf(titl, "LCT firmware-emulator: exact match, endcap %d", endc + 1);
3808     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, titl);
3809     title->SetFillColor(10);
3810     title->Draw();
3811     sprintf(pagenum, "- %d -", page);
3812     t.DrawText(0.9, 0.02, pagenum);
3813     //gStyle->SetOptStat(110010);
3814     gStyle->SetOptStat(0);
3815     pad[page]->Draw();
3816     pad[page]->Divide(2, 5);
3817     TH1F* hLctMatchEffVsCsc[MAX_ENDCAPS][CSC_TYPES];
3818     for (int idh = 0; idh < max_idh; idh++) {
3819       if (!plotME1A && idh == 3)
3820         continue;
3821       hLctMatchEffVsCsc[endc][idh] = (TH1F*)hLctCompTotalCsc[endc][idh]->Clone();
3822       hLctMatchEffVsCsc[endc][idh]->Divide(hLctCompMatchCsc[endc][idh], hLctCompTotalCsc[endc][idh], 1., 1., "B");
3823       nbins = hLctCompTotalCsc[endc][idh]->GetNbinsX();
3824       for (Int_t ibin = 1; ibin <= nbins; ibin++) {
3825         if (hLctCompTotalCsc[endc][idh]->GetBinContent(ibin) == 0) {
3826           hLctMatchEffVsCsc[endc][idh]->SetBinContent(ibin, -1.);
3827           hLctMatchEffVsCsc[endc][idh]->SetBinError(ibin, 0.);
3828         }
3829       }
3830       gPad->Update();
3831       gStyle->SetStatX(0.65);
3832       hLctMatchEffVsCsc[endc][idh]->SetMinimum(-0.05);
3833       hLctMatchEffVsCsc[endc][idh]->SetMaximum(1.05);
3834       hLctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleOffset(0.7);
3835       hLctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleOffset(0.8);
3836       hLctMatchEffVsCsc[endc][idh]->GetXaxis()->SetLabelSize(0.06);
3837       hLctMatchEffVsCsc[endc][idh]->GetYaxis()->SetLabelSize(0.06);
3838       hLctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitleSize(0.07);
3839       hLctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitleSize(0.07);
3840       hLctMatchEffVsCsc[endc][idh]->GetXaxis()->SetTitle("CSC id");
3841       hLctMatchEffVsCsc[endc][idh]->GetYaxis()->SetTitle("% of exact match");
3842       pad[page]->cd(idh + 1);
3843       hLctMatchEffVsCsc[endc][idh]->Draw("e");
3844       double numer = hLctCompMatchCsc[endc][idh]->Integral();
3845       double denom = hLctCompTotalCsc[endc][idh]->Integral();
3846       double ratio = 0.0, error = 0.0;
3847       if (denom > 0.) {
3848         ratio = numer / denom;
3849         error = sqrt(ratio * (1. - ratio) / denom);
3850       }
3851       sprintf(eff, "eff = (%4.1f +/- %4.1f)%%", ratio * 100., error * 100.);
3852       teff.DrawTextNDC(0.3, 0.5, eff);
3853     }
3854     page++;
3855     c1->Update();
3856   }
3857
3858   ps->Close();
3859   delete c1;
3860 }
3861
3862 void CSCTriggerPrimitivesReader::drawResolHistos() {
3863   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
3864   string fname = resultsFileNamesPrefix_ + "lcts_resol.ps";
3865   TPostScript* ps = new TPostScript(fname.c_str(), 111);
3866
3867   TPad* pad[MAXPAGES];
3868   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
3869     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
3870   }
3871
3872   int page = 1;
3873   TText t;
3874   t.SetTextFont(32);
3875   t.SetTextSize(0.025);
3876   char pagenum[7];
3877   TPaveLabel* title;
3878
3879   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
3880
3881   ps->NewPage();
3882   c1->Clear();
3883   c1->cd(0);
3884   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "ALCT resolution");
3885   title->SetFillColor(10);
3886   title->Draw();
3887   sprintf(pagenum, "- %d -", page);
3888   t.DrawText(0.9, 0.02, pagenum);
3889   gStyle->SetOptStat(111110);
3890   pad[page]->Draw();
3891   pad[page]->Divide(2, 2);
3892   gStyle->SetStatX(1.00);
3893   gStyle->SetStatY(0.65);
3894   pad[page]->cd(1);
3895   hEtaRecVsSim->SetMarkerSize(0.2);
3896   hEtaRecVsSim->Draw();
3897   gPad->Update();
3898   gStyle->SetStatX(1.00);
3899   gStyle->SetStatY(0.995);
3900   hResolDeltaWG->GetXaxis()->SetTitle("WG_{rec} - WG_{sim}");
3901   hResolDeltaWG->GetXaxis()->SetTitleOffset(1.2);
3902   hResolDeltaWG->GetYaxis()->SetTitle("Entries");
3903   hResolDeltaWG->GetYaxis()->SetTitleOffset(1.9);
3904   hResolDeltaWG->GetXaxis()->SetLabelSize(0.03);
3905   hResolDeltaWG->GetYaxis()->SetLabelSize(0.03);
3906   pad[page]->cd(3);
3907   hResolDeltaWG->Draw();
3908   hResolDeltaEta->GetXaxis()->SetNdivisions(505);  // twice fewer divisions
3909   hResolDeltaEta->GetXaxis()->SetLabelSize(0.04);
3910   hResolDeltaEta->GetYaxis()->SetLabelSize(0.04);
3911   pad[page]->cd(4);
3912   hResolDeltaEta->Draw();
3913   hResolDeltaEta->Fit("gaus", "Q");
3914   page++;
3915   c1->Update();
3916
3917   ps->NewPage();
3918   c1->Clear();
3919   c1->cd(0);
3920   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#eta_rec-#eta_sim");
3921   title->SetFillColor(10);
3922   title->Draw();
3923   sprintf(pagenum, "- %d -", page);
3924   t.DrawText(0.9, 0.02, pagenum);
3925   pad[page]->Draw();
3926   pad[page]->Divide(2, 5);
3927   for (int idh = 0; idh < max_idh; idh++) {
3928     if (!plotME1A && idh == 3)
3929       continue;
3930     hEtaDiffCsc[idh][0]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
3931     hEtaDiffCsc[idh][0]->GetYaxis()->SetLabelSize(0.07);
3932     pad[page]->cd(idh + 1);
3933     hEtaDiffCsc[idh][0]->Draw();
3934     if (hEtaDiffCsc[idh][0]->GetEntries() > 1)
3935       hEtaDiffCsc[idh][0]->Fit("gaus", "Q");
3936   }
3937   page++;
3938   c1->Update();
3939
3940   ps->NewPage();
3941   c1->Clear();
3942   c1->cd(0);
3943   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#eta_rec-#eta_sim, endcap1");
3944   title->SetFillColor(10);
3945   title->Draw();
3946   sprintf(pagenum, "- %d -", page);
3947   t.DrawText(0.9, 0.02, pagenum);
3948   pad[page]->Draw();
3949   pad[page]->Divide(2, 5);
3950   for (int idh = 0; idh < max_idh; idh++) {
3951     if (!plotME1A && idh == 3)
3952       continue;
3953     hEtaDiffCsc[idh][1]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
3954     hEtaDiffCsc[idh][1]->GetYaxis()->SetLabelSize(0.07);
3955     pad[page]->cd(idh + 1);
3956     hEtaDiffCsc[idh][1]->Draw();
3957     if (hEtaDiffCsc[idh][1]->GetEntries() > 1)
3958       hEtaDiffCsc[idh][1]->Fit("gaus", "Q");
3959   }
3960   page++;
3961   c1->Update();
3962
3963   ps->NewPage();
3964   c1->Clear();
3965   c1->cd(0);
3966   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#eta_rec-#eta_sim, endcap2");
3967   title->SetFillColor(10);
3968   title->Draw();
3969   sprintf(pagenum, "- %d -", page);
3970   t.DrawText(0.9, 0.02, pagenum);
3971   pad[page]->Draw();
3972   pad[page]->Divide(2, 5);
3973   for (int idh = 0; idh < max_idh; idh++) {
3974     if (!plotME1A && idh == 3)
3975       continue;
3976     hEtaDiffCsc[idh][2]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
3977     hEtaDiffCsc[idh][2]->GetYaxis()->SetLabelSize(0.07);
3978     pad[page]->cd(idh + 1);
3979     hEtaDiffCsc[idh][2]->Draw();
3980     if (hEtaDiffCsc[idh][2]->GetEntries() > 1)
3981       hEtaDiffCsc[idh][2]->Fit("gaus", "Q");
3982   }
3983   page++;
3984   c1->Update();
3985
3986   ps->NewPage();
3987   c1->Clear();
3988   c1->cd(0);
3989   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#LT#eta_rec-#eta_sim#GT vs #eta_rec");
3990   title->SetFillColor(10);
3991   title->Draw();
3992   sprintf(pagenum, "- %d -", page);
3993   t.DrawText(0.9, 0.02, pagenum);
3994   //gStyle->SetOptStat(0);
3995   pad[page]->Draw();
3996   pad[page]->Divide(2, 2);
3997   TH1F* hMeanEtaDiffVsEta[MAX_STATIONS];
3998   for (int istation = 0; istation < MAX_STATIONS; istation++) {
3999     hMeanEtaDiffVsEta[istation] = (TH1F*)hEtaDiffVsEta[istation]->Clone();
4000     hMeanEtaDiffVsEta[istation]->Divide(hEtaDiffVsEta[istation], hAlctVsEta[istation], 1., 1.);
4001     hMeanEtaDiffVsEta[istation]->GetXaxis()->SetTitleOffset(1.2);
4002     hMeanEtaDiffVsEta[istation]->GetXaxis()->SetTitle("#eta");
4003     hMeanEtaDiffVsEta[istation]->SetMaximum(0.05);
4004     pad[page]->cd(istation + 1);
4005     hMeanEtaDiffVsEta[istation]->Draw();
4006   }
4007   page++;
4008   c1->Update();
4009
4010   c1->Clear();
4011   c1->cd(0);
4012   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#eta_rec-#eta_sim vs wiregroup");
4013   title->SetFillColor(10);
4014   title->Draw();
4015   sprintf(pagenum, "- %d -", page);
4016   t.DrawText(0.9, 0.02, pagenum);
4017   pad[page]->Draw();
4018   pad[page]->Divide(2, 5);
4019   for (int idh = 0; idh < max_idh; idh++) {
4020     if (!plotME1A && idh == 3)
4021       continue;
4022     pad[page]->cd(idh + 1);
4023     hEtaDiffVsWireCsc[idh]->SetMarkerSize(0.2);
4024     hEtaDiffVsWireCsc[idh]->GetXaxis()->SetTitle("Wiregroup");
4025     hEtaDiffVsWireCsc[idh]->GetXaxis()->SetTitleSize(0.07);
4026     hEtaDiffVsWireCsc[idh]->GetXaxis()->SetTitleOffset(1.2);
4027     hEtaDiffVsWireCsc[idh]->GetYaxis()->SetTitle("#eta_rec-#eta_sim");
4028     hEtaDiffVsWireCsc[idh]->GetYaxis()->SetTitleSize(0.07);
4029     hEtaDiffVsWireCsc[idh]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4030     hEtaDiffVsWireCsc[idh]->GetYaxis()->SetLabelSize(0.07);
4031     hEtaDiffVsWireCsc[idh]->Draw();
4032   }
4033   page++;
4034   c1->Update();
4035
4036   ps->NewPage();
4037   c1->Clear();
4038   c1->cd(0);
4039   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi resolution");
4040   title->SetFillColor(10);
4041   title->Draw();
4042   sprintf(pagenum, "- %d -", page);
4043   t.DrawText(0.9, 0.02, pagenum);
4044   gStyle->SetOptStat(111110);
4045   pad[page]->Draw();
4046   pad[page]->Divide(2, 2);
4047   gStyle->SetStatX(1.00);
4048   gStyle->SetStatY(0.65);
4049   pad[page]->cd(1);
4050   hPhiRecVsSim->SetMarkerSize(0.2);
4051   hPhiRecVsSim->Draw();
4052   gPad->Update();
4053   gStyle->SetStatX(1.00);
4054   gStyle->SetStatY(0.995);
4055   hResolDeltaPhi->GetXaxis()->SetLabelSize(0.04);
4056   hResolDeltaPhi->GetYaxis()->SetLabelSize(0.04);
4057   pad[page]->cd(2);
4058   hResolDeltaPhi->Draw();
4059   hResolDeltaPhi->Fit("gaus", "Q");
4060   hResolDeltaHS->GetXaxis()->SetTitle("HS_{rec} - HS_{sim}");
4061   hResolDeltaHS->GetXaxis()->SetTitleOffset(1.2);
4062   hResolDeltaHS->GetYaxis()->SetTitle("Entries");
4063   hResolDeltaHS->GetYaxis()->SetTitleOffset(1.7);
4064   hResolDeltaHS->GetXaxis()->SetLabelSize(0.03);
4065   hResolDeltaHS->GetYaxis()->SetLabelSize(0.03);
4066   pad[page]->cd(3);
4067   hResolDeltaHS->Draw();
4068   hResolDeltaDS->GetXaxis()->SetTitle("DS_{rec} - DS_{sim}");
4069   hResolDeltaDS->GetXaxis()->SetTitleOffset(1.2);
4070   hResolDeltaDS->GetYaxis()->SetTitle("Entries");
4071   hResolDeltaDS->GetYaxis()->SetTitleOffset(1.6);
4072   hResolDeltaDS->GetXaxis()->SetLabelSize(0.04);
4073   hResolDeltaDS->GetYaxis()->SetLabelSize(0.04);
4074   pad[page]->cd(4);
4075   hResolDeltaDS->Draw();
4076   page++;
4077   c1->Update();
4078
4079   ps->NewPage();
4080   c1->Clear();
4081   c1->cd(0);
4082   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim (mrad)");
4083   title->SetFillColor(10);
4084   title->Draw();
4085   sprintf(pagenum, "- %d -", page);
4086   t.DrawText(0.9, 0.02, pagenum);
4087   pad[page]->Draw();
4088   pad[page]->Divide(2, 5);
4089   for (int idh = 0; idh < max_idh; idh++) {
4090     if (!plotME1A && idh == 3)
4091       continue;
4092     hPhiDiffCsc[idh][0]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4093     hPhiDiffCsc[idh][0]->GetYaxis()->SetLabelSize(0.07);
4094     pad[page]->cd(idh + 1);
4095     hPhiDiffCsc[idh][0]->Draw();
4096     if (hPhiDiffCsc[idh][0]->GetEntries() > 1)
4097       hPhiDiffCsc[idh][0]->Fit("gaus", "Q");
4098   }
4099   page++;
4100   c1->Update();
4101
4102   ps->NewPage();
4103   c1->Clear();
4104   c1->cd(0);
4105   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim (mrad), endcap1");
4106   title->SetFillColor(10);
4107   title->Draw();
4108   sprintf(pagenum, "- %d -", page);
4109   t.DrawText(0.9, 0.02, pagenum);
4110   pad[page]->Draw();
4111   pad[page]->Divide(2, 5);
4112   for (int idh = 0; idh < max_idh; idh++) {
4113     if (!plotME1A && idh == 3)
4114       continue;
4115     hPhiDiffCsc[idh][1]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4116     hPhiDiffCsc[idh][1]->GetYaxis()->SetLabelSize(0.07);
4117     pad[page]->cd(idh + 1);
4118     hPhiDiffCsc[idh][1]->Draw();
4119     if (hPhiDiffCsc[idh][1]->GetEntries() > 1)
4120       hPhiDiffCsc[idh][1]->Fit("gaus", "Q");
4121   }
4122   page++;
4123   c1->Update();
4124
4125   ps->NewPage();
4126   c1->Clear();
4127   c1->cd(0);
4128   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim (mrad), endcap2");
4129   title->SetFillColor(10);
4130   title->Draw();
4131   sprintf(pagenum, "- %d -", page);
4132   t.DrawText(0.9, 0.02, pagenum);
4133   pad[page]->Draw();
4134   pad[page]->Divide(2, 5);
4135   for (int idh = 0; idh < max_idh; idh++) {
4136     if (!plotME1A && idh == 3)
4137       continue;
4138     hPhiDiffCsc[idh][2]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4139     hPhiDiffCsc[idh][2]->GetYaxis()->SetLabelSize(0.07);
4140     pad[page]->cd(idh + 1);
4141     hPhiDiffCsc[idh][2]->Draw();
4142     if (hPhiDiffCsc[idh][2]->GetEntries() > 1)
4143       hPhiDiffCsc[idh][2]->Fit("gaus", "Q");
4144   }
4145   page++;
4146   c1->Update();
4147
4148   ps->NewPage();
4149   c1->Clear();
4150   c1->cd(0);
4151   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#LT#phi_rec-#phi_sim#GT (mrad) vs #phi_rec");
4152   title->SetFillColor(10);
4153   title->Draw();
4154   sprintf(pagenum, "- %d -", page);
4155   t.DrawText(0.9, 0.02, pagenum);
4156   gStyle->SetOptStat(0);
4157   pad[page]->Draw();
4158   pad[page]->Divide(2, 2);
4159   TH1F* hMeanPhiDiffVsPhi[MAX_STATIONS];
4160   for (int istation = 0; istation < MAX_STATIONS; istation++) {
4161     hMeanPhiDiffVsPhi[istation] = (TH1F*)hPhiDiffVsPhi[istation]->Clone();
4162     hMeanPhiDiffVsPhi[istation]->Divide(hPhiDiffVsPhi[istation], hClctVsPhi[istation], 1., 1.);
4163     hMeanPhiDiffVsPhi[istation]->GetXaxis()->SetTitleOffset(1.2);
4164     hMeanPhiDiffVsPhi[istation]->GetYaxis()->SetTitleOffset(1.7);
4165     hMeanPhiDiffVsPhi[istation]->GetXaxis()->SetTitle("#phi");
4166     hMeanPhiDiffVsPhi[istation]->GetYaxis()->SetTitle("#LT#phi_rec-#phi_sim#GT (mrad)");
4167     hMeanPhiDiffVsPhi[istation]->SetMaximum(5.);
4168     pad[page]->cd(istation + 1);
4169     hMeanPhiDiffVsPhi[istation]->Draw();
4170   }
4171   page++;
4172   c1->Update();
4173
4174   c1->Clear();
4175   c1->cd(0);
4176   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim (mrad) vs halfstrip #");
4177   title->SetFillColor(10);
4178   title->Draw();
4179   sprintf(pagenum, "- %d -", page);
4180   t.DrawText(0.9, 0.02, pagenum);
4181   gStyle->SetOptStat(0);
4182   pad[page]->Draw();
4183   pad[page]->Divide(2, 5);
4184   for (int idh = 0; idh < max_idh; idh++) {
4185     if (!plotME1A && idh == 3)
4186       continue;
4187     pad[page]->cd(idh + 1);
4188     hPhiDiffVsStripCsc[idh][1]->SetMarkerSize(0.2);
4189     hPhiDiffVsStripCsc[idh][1]->GetXaxis()->SetTitle("Halfstrip");
4190     hPhiDiffVsStripCsc[idh][1]->GetXaxis()->SetTitleOffset(1.4);
4191     hPhiDiffVsStripCsc[idh][1]->GetYaxis()->SetTitle("#phi_rec-#phi_sim (mrad)");
4192     hPhiDiffVsStripCsc[idh][1]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4193     hPhiDiffVsStripCsc[idh][1]->GetYaxis()->SetLabelSize(0.07);
4194     hPhiDiffVsStripCsc[idh][1]->Draw();
4195   }
4196   page++;
4197   c1->Update();
4198
4199   c1->Clear();
4200   c1->cd(0);
4201   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim, halfstrips only, different patterns");
4202   title->SetFillColor(10);
4203   title->Draw();
4204   sprintf(pagenum, "- %d -", page);
4205   t.DrawText(0.9, 0.02, pagenum);
4206   gStyle->SetOptStat(111110);
4207   pad[page]->Draw();
4208   pad[page]->Divide(3, 3);
4209   int min_pattern, max_pattern;
4210   min_pattern = 2;
4211   max_pattern = CSCConstants::NUM_CLCT_PATTERNS;
4212   for (int idh = min_pattern; idh < max_pattern; idh++) {
4213     hPhiDiffPattern[idh]->GetXaxis()->SetTitle("Halfstrip");
4214     hPhiDiffPattern[idh]->GetXaxis()->SetTitleOffset(1.2);
4215     pad[page]->cd(idh - min_pattern + 1);
4216     hPhiDiffPattern[idh]->Draw();
4217     // if (hPhiDiffPattern[idh]->GetEntries() > 1)
4218     //   hPhiDiffPattern[idh]->Fit("gaus","Q");
4219   }
4220   page++;
4221   c1->Update();
4222
4223   ps->NewPage();
4224   c1->Clear();
4225   c1->cd(0);
4226   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_1-#phi_6 (mrad), muon SimHits");
4227   title->SetFillColor(10);
4228   title->Draw();
4229   sprintf(pagenum, "- %d -", page);
4230   t.DrawText(0.9, 0.02, pagenum);
4231   pad[page]->Draw();
4232   pad[page]->Divide(2, 5);
4233   for (int idh = 0; idh < max_idh; idh++) {
4234     if (!plotME1A && idh == 3)
4235       continue;
4236     hTrueBendCsc[idh]->GetYaxis()->SetTitle("Entries");
4237     hTrueBendCsc[idh]->GetYaxis()->SetTitleSize(0.07);
4238     hTrueBendCsc[idh]->GetXaxis()->SetLabelSize(0.07);  // default=0.04
4239     hTrueBendCsc[idh]->GetYaxis()->SetLabelSize(0.07);
4240     pad[page]->cd(idh + 1);
4241     hTrueBendCsc[idh]->Draw();
4242   }
4243   page++;
4244   c1->Update();
4245
4246   ps->Close();
4247   delete c1;
4248 }
4249
4250 void CSCTriggerPrimitivesReader::drawEfficHistos() {
4251   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 700);
4252   string fname = resultsFileNamesPrefix_ + "lcts_effic.ps";
4253   TPostScript* ps = new TPostScript(fname.c_str(), 111);
4254
4255   TPad* pad[MAXPAGES];
4256   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
4257     pad[i_page] = new TPad("", "", .05, .05, .93, .93);
4258   }
4259
4260   int page = 1;
4261   TText t;
4262   t.SetTextFont(32);
4263   t.SetTextSize(0.025);
4264   char pagenum[7];
4265   TPaveLabel* title;
4266   char histtitle[60];
4267
4268   gStyle->SetOptDate(0);
4269   gStyle->SetTitleSize(0.1, "");  // size for pad title; default is 0.02
4270
4271   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
4272
4273   ps->NewPage();
4274   c1->Clear();
4275   c1->cd(0);
4276   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "ALCT efficiency vs #eta");
4277   title->SetFillColor(10);
4278   title->Draw();
4279   sprintf(pagenum, "- %d -", page);
4280   t.DrawText(0.9, 0.02, pagenum);
4281   gStyle->SetOptStat(0);
4282   pad[page]->Draw();
4283   pad[page]->Divide(2, 2);
4284   TH1F* hALCTEffVsEta[MAX_STATIONS];
4285   for (int istation = 0; istation < MAX_STATIONS; istation++) {
4286     hALCTEffVsEta[istation] = (TH1F*)hEfficHitsEta[istation]->Clone();
4287     hALCTEffVsEta[istation]->Divide(hEfficALCTEta[istation], hEfficHitsEta[istation], 1., 1., "B");
4288     hALCTEffVsEta[istation]->GetXaxis()->SetTitleOffset(1.2);
4289     hALCTEffVsEta[istation]->GetXaxis()->SetTitle("#eta");
4290     hALCTEffVsEta[istation]->SetMaximum(1.05);
4291     sprintf(histtitle, "ALCT efficiency vs #eta, station %d", istation + 1);
4292     hALCTEffVsEta[istation]->SetTitle(histtitle);
4293     pad[page]->cd(istation + 1);
4294     hALCTEffVsEta[istation]->Draw();
4295   }
4296   page++;
4297   c1->Update();
4298
4299   c1->Clear();
4300   c1->cd(0);
4301   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "ALCT efficiency vs #eta");
4302   title->SetFillColor(10);
4303   title->Draw();
4304   sprintf(pagenum, "- %d -", page);
4305   t.DrawText(0.9, 0.02, pagenum);
4306   gStyle->SetOptStat(11111);
4307   pad[page]->Draw();
4308   pad[page]->Divide(2, 5);
4309   TH1F* hALCTEffVsEtaCsc[CSC_TYPES];
4310   for (int idh = 0; idh < max_idh; idh++) {
4311     if (!plotME1A && idh == 3)
4312       continue;
4313     hALCTEffVsEtaCsc[idh] = (TH1F*)hEfficHitsEtaCsc[idh]->Clone();
4314     hALCTEffVsEtaCsc[idh]->Divide(hEfficALCTEtaCsc[idh], hEfficHitsEtaCsc[idh], 1., 1., "B");
4315     if (idh == 3 || idh == 4 || idh == 6 || idh == 8) {
4316       gPad->Update();
4317       gStyle->SetStatX(0.43);
4318     } else {
4319       gPad->Update();
4320       gStyle->SetStatX(1.00);
4321     }
4322     hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitle("#eta");
4323     hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleOffset(0.8);
4324     hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleSize(0.07);  // default=0.05
4325     hALCTEffVsEtaCsc[idh]->GetXaxis()->SetLabelSize(0.10);  // default=0.04
4326     hALCTEffVsEtaCsc[idh]->GetYaxis()->SetLabelSize(0.10);
4327     hALCTEffVsEtaCsc[idh]->SetLabelOffset(0.012, "XY");
4328     hALCTEffVsEtaCsc[idh]->SetMinimum(0.50);
4329     hALCTEffVsEtaCsc[idh]->SetMaximum(1.05);
4330     hALCTEffVsEtaCsc[idh]->SetTitle(csc_type[idh].c_str());
4331     hALCTEffVsEtaCsc[idh]->SetTitleSize(0.1, "");
4332     hALCTEffVsEtaCsc[idh]->SetLineWidth(2);
4333     hALCTEffVsEtaCsc[idh]->SetLineColor(4);
4334     pad[page]->cd(idh + 1);
4335     gPad->SetGrid(1);
4336     hALCTEffVsEtaCsc[idh]->Draw();
4337   }
4338   page++;
4339   c1->Update();
4340
4341   ps->NewPage();
4342   c1->Clear();
4343   c1->cd(0);
4344   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT efficiency vs #eta");
4345   title->SetFillColor(10);
4346   title->Draw();
4347   sprintf(pagenum, "- %d -", page);
4348   t.DrawText(0.9, 0.02, pagenum);
4349   gStyle->SetOptStat(0);
4350   pad[page]->Draw();
4351   pad[page]->Divide(2, 2);
4352   TH1F* hCLCTEffVsEta[MAX_STATIONS];
4353   for (int istation = 0; istation < MAX_STATIONS; istation++) {
4354     hCLCTEffVsEta[istation] = (TH1F*)hEfficHitsEta[istation]->Clone();
4355     hCLCTEffVsEta[istation]->Divide(hEfficCLCTEta[istation], hEfficHitsEta[istation], 1., 1., "B");
4356     hCLCTEffVsEta[istation]->GetXaxis()->SetTitleOffset(1.2);
4357     hCLCTEffVsEta[istation]->GetXaxis()->SetTitle("#eta");
4358     hCLCTEffVsEta[istation]->SetMaximum(1.05);
4359     sprintf(histtitle, "CLCT efficiency vs #eta, station %d", istation + 1);
4360     hCLCTEffVsEta[istation]->SetTitle(histtitle);
4361     pad[page]->cd(istation + 1);
4362     hCLCTEffVsEta[istation]->Draw();
4363   }
4364   page++;
4365   c1->Update();
4366
4367   c1->Clear();
4368   c1->cd(0);
4369   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT efficiency vs #eta");
4370   title->SetFillColor(10);
4371   title->Draw();
4372   sprintf(pagenum, "- %d -", page);
4373   t.DrawText(0.9, 0.02, pagenum);
4374   gStyle->SetOptStat(111110);
4375   pad[page]->Draw();
4376   pad[page]->Divide(2, 5);
4377   TH1F* hCLCTEffVsEtaCsc[CSC_TYPES];
4378   for (int idh = 0; idh < max_idh; idh++) {
4379     if (!plotME1A && idh == 3)
4380       continue;
4381     hCLCTEffVsEtaCsc[idh] = (TH1F*)hEfficHitsEtaCsc[idh]->Clone();
4382     hCLCTEffVsEtaCsc[idh]->Divide(hEfficCLCTEtaCsc[idh], hEfficHitsEtaCsc[idh], 1., 1., "B");
4383     if (idh == 3 || idh == 4 || idh == 6 || idh == 8) {
4384       gPad->Update();
4385       gStyle->SetStatX(0.43);
4386     } else {
4387       gPad->Update();
4388       gStyle->SetStatX(1.00);
4389     }
4390     hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitle("#eta");
4391     hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleOffset(0.8);
4392     hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleSize(0.07);  // default=0.05
4393     hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetLabelSize(0.10);  // default=0.04
4394     hCLCTEffVsEtaCsc[idh]->GetYaxis()->SetLabelSize(0.10);
4395     hCLCTEffVsEtaCsc[idh]->SetLabelOffset(0.012, "XY");
4396     hCLCTEffVsEtaCsc[idh]->SetMinimum(0.50);
4397     hCLCTEffVsEtaCsc[idh]->SetMaximum(1.05);
4398     hCLCTEffVsEtaCsc[idh]->SetTitle(csc_type[idh].c_str());
4399     hCLCTEffVsEtaCsc[idh]->SetLineWidth(2);
4400     hCLCTEffVsEtaCsc[idh]->SetLineColor(4);
4401     pad[page]->cd(idh + 1);
4402     gPad->SetGrid(1);
4403     hCLCTEffVsEtaCsc[idh]->Draw();
4404   }
4405   page++;
4406   c1->Update();
4407
4408   ps->Close();
4409   delete c1;
4410 }
4411
4412 void CSCTriggerPrimitivesReader::drawHistosForTalks() {
4413   TCanvas* c1 = new TCanvas("c1", "", 0, 0, 500, 640);
4414   TCanvas* c2 = new TCanvas("c2", "", 0, 0, 540, 540);
4415
4416   TPad* pad[MAXPAGES];
4417   for (int i_page = 0; i_page < MAXPAGES; i_page++) {
4418     pad[i_page] = new TPad("", "", .07, .07, .93, .93);
4419   }
4420
4421   int page = 1;
4422   TPaveLabel* title;
4423   gStyle->SetOptDate(0);
4424
4425   int max_idh = plotME42 ? CSC_TYPES : CSC_TYPES - 1;
4426
4427   TPostScript* eps1 = new TPostScript("clcts.eps", 113);
4428   eps1->NewPage();
4429   c1->Clear();
4430   c1->cd(0);
4431   title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT quantities");
4432   title->SetFillColor(10);
4433   title->Draw();
4434   pad[page]->Draw();
4435   pad[page]->Divide(2, 3);
4436   pad[page]->cd(1);
4437   hClctQuality->Draw();
4438   pad[page]->cd(2);
4439   hClctSign->Draw();
4440   TH1F* hClctPatternTot = (TH1F*)hClctPattern[0]->Clone();
4441   hClctPatternTot->SetTitle("CLCT pattern #");
4442   hClctPatternTot->Add(hClctPattern[0], hClctPattern[1], 1., 1.);
4443   pad[page]->cd(3);
4444   hClctPatternTot->Draw();
4445   hClctPattern[0]->SetLineStyle(2);
4446   hClctPattern[0]->Draw("same");
4447   hClctPattern[1]->SetLineStyle(3);
4448   hClctPattern[1]->Draw("same");
4449   pad[page]->cd(4);
4450   hClctCFEB->Draw();
4451   pad[page]->cd(5);
4452   hClctStripType->Draw();
4453   TH1F* hClctKeyStripTot = (TH1F*)hClctKeyStrip[0]->Clone();
4454   hClctKeyStripTot->SetTitle("CLCT key strip #");
4455   hClctKeyStripTot->Add(hClctKeyStrip[0], hClctKeyStrip[1], 1., 1.);
4456   pad[page]->cd(6);
4457   hClctKeyStripTot->Draw();
4458   page++;
4459   c1->Update();
4460   c1->Print("asdf.png");
4461   eps1->Close();
4462
4463   // Resolution histograms.
4464   if (bookedResolHistos) {
4465     gStyle->SetTitleSize(0.055, "");  // size for pad title; default is 0.02
4466
4467     TPostScript* eps2 = new TPostScript("alct_deltaWG.eps", 113);
4468     eps2->NewPage();
4469     c2->Clear();
4470     c2->cd(0);
4471     gStyle->SetOptStat(0);
4472     pad[page]->Draw();
4473     pad[page]->Divide(1, 1);
4474     hResolDeltaWG->GetXaxis()->SetTitle("WG_{rec} - WG_{sim}");
4475     hResolDeltaWG->GetXaxis()->SetTitleOffset(1.2);
4476     hResolDeltaWG->GetYaxis()->SetTitle("Entries");
4477     hResolDeltaWG->GetYaxis()->SetTitleOffset(1.9);
4478     hResolDeltaWG->GetXaxis()->SetLabelSize(0.04);
4479     hResolDeltaWG->GetYaxis()->SetLabelSize(0.04);
4480     pad[page]->cd(1);
4481     hResolDeltaWG->Draw();
4482     page++;
4483     c2->Update();
4484     eps2->Close();
4485
4486     TPostScript* eps3 = new TPostScript("clct_deltaHS.eps", 113);
4487     eps3->NewPage();
4488     c2->Clear();
4489     c2->cd(0);
4490     pad[page]->Draw();
4491     pad[page]->Divide(1, 1);
4492     hResolDeltaHS->GetXaxis()->SetTitle("HS_{rec} - HS_{sim}");
4493     hResolDeltaHS->GetXaxis()->SetTitleOffset(1.2);
4494     hResolDeltaHS->GetYaxis()->SetTitle("Entries");
4495     hResolDeltaHS->GetYaxis()->SetTitleOffset(1.7);
4496     hResolDeltaHS->GetXaxis()->SetLabelSize(0.04);  // default=0.04
4497     hResolDeltaHS->GetYaxis()->SetLabelSize(0.04);
4498     pad[page]->cd(1);
4499     hResolDeltaHS->Draw();
4500     page++;
4501     c2->Update();
4502     eps3->Close();
4503
4504     TPostScript* eps4 = new TPostScript("clct_deltaDS.eps", 113);
4505     eps4->NewPage();
4506     c2->Clear();
4507     c2->cd(0);
4508     gStyle->SetOptStat(0);
4509     pad[page]->Draw();
4510     pad[page]->Divide(1, 1);
4511     hResolDeltaDS->GetXaxis()->SetTitle("DS_{rec} - DS_{sim}");
4512     hResolDeltaDS->GetXaxis()->SetTitleOffset(1.2);
4513     hResolDeltaDS->GetYaxis()->SetTitle("Entries");
4514     hResolDeltaDS->GetYaxis()->SetTitleOffset(1.6);
4515     hResolDeltaDS->GetXaxis()->SetLabelSize(0.04);  // default=0.04
4516     hResolDeltaDS->GetYaxis()->SetLabelSize(0.04);
4517     pad[page]->cd(1);
4518     hResolDeltaDS->Draw();
4519     page++;
4520     c2->Update();
4521     eps4->Close();
4522
4523     TPostScript* eps5 = new TPostScript("clct_deltaPhi_hs.eps", 113);
4524     eps5->NewPage();
4525     c1->Clear();
4526     c1->cd(0);
4527     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "#phi_rec-#phi_sim (mrad), halfstrips only");
4528     title->SetFillColor(10);
4529     title->Draw();
4530     gStyle->SetTitleSize(0.1, "");  // size for pad title; default is 0.02
4531     gStyle->SetOptStat(111110);
4532     pad[page]->Draw();
4533     pad[page]->Divide(2, 5);
4534     for (int idh = 0; idh < max_idh; idh++) {
4535       if (!plotME1A && idh == 3)
4536         continue;
4537       pad[page]->cd(idh + 1);
4538       hPhiDiffCsc[idh][4]->Draw();
4539       //if (hPhiDiffCsc[idh][4]->GetEntries() > 1)
4540       //hPhiDiffCsc[idh][4]->Fit("gaus","Q");
4541       //hPhiDiffCsc[idh][4]->GetXaxis()->SetTitle("#phi_{rec} - #phi_{sim} (mrad)");
4542       //hPhiDiffCsc[idh][4]->GetXaxis()->SetTitleSize(0.06);
4543       //hPhiDiffCsc[idh][4]->GetXaxis()->SetTitleOffset(0.9);
4544       hPhiDiffCsc[idh][4]->GetYaxis()->SetTitle("Entries        ");
4545       hPhiDiffCsc[idh][4]->GetYaxis()->SetTitleSize(0.07);
4546       hPhiDiffCsc[idh][4]->GetYaxis()->SetTitleOffset(1.0);
4547       hPhiDiffCsc[idh][4]->GetXaxis()->SetLabelSize(0.10);  // default=0.04
4548       hPhiDiffCsc[idh][4]->GetYaxis()->SetLabelSize(0.10);
4549       hPhiDiffCsc[idh][4]->SetLabelOffset(0.012, "XY");
4550     }
4551     page++;
4552     c1->Update();
4553     eps5->Close();
4554   }
4555
4556   // Efficiency histograms.
4557   if (bookedEfficHistos) {
4558     TPostScript* eps6 = new TPostScript("alct_effic.eps", 113);
4559     eps6->NewPage();
4560     c1->Clear();
4561     c1->cd(0);
4562     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "ALCT efficiency vs #eta");
4563     title->SetFillColor(10);
4564     title->Draw();
4565     gStyle->SetOptStat(0);
4566     pad[page]->Draw();
4567     pad[page]->Divide(2, 5);
4568     TH1F* hALCTEffVsEtaCsc[CSC_TYPES];
4569     for (int idh = 0; idh < max_idh; idh++) {
4570       if (!plotME1A && idh == 3)
4571         continue;
4572       hALCTEffVsEtaCsc[idh] = (TH1F*)hEfficHitsEtaCsc[idh]->Clone();
4573       hALCTEffVsEtaCsc[idh]->Divide(hEfficALCTEtaCsc[idh], hEfficHitsEtaCsc[idh], 1., 1., "B");
4574       hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitle("#eta");
4575       hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleOffset(0.8);
4576       hALCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleSize(0.07);  // default=0.05
4577       hALCTEffVsEtaCsc[idh]->GetXaxis()->SetLabelSize(0.10);  // default=0.04
4578       hALCTEffVsEtaCsc[idh]->GetYaxis()->SetLabelSize(0.10);
4579       hALCTEffVsEtaCsc[idh]->SetLabelOffset(0.012, "XY");
4580       hALCTEffVsEtaCsc[idh]->SetMinimum(0.50);
4581       hALCTEffVsEtaCsc[idh]->SetMaximum(1.05);
4582       hALCTEffVsEtaCsc[idh]->SetTitle(csc_type[idh].c_str());
4583       hALCTEffVsEtaCsc[idh]->SetTitleSize(0.1, "");
4584       hALCTEffVsEtaCsc[idh]->SetLineWidth(2);
4585       hALCTEffVsEtaCsc[idh]->SetLineColor(4);
4586       hALCTEffVsEtaCsc[idh]->SetLineColor(4);
4587       pad[page]->cd(idh + 1);
4588       gPad->SetGrid(1);
4589       hALCTEffVsEtaCsc[idh]->Draw();
4590     }
4591     page++;
4592     c1->Update();
4593     eps6->Close();
4594
4595     TPostScript* eps7 = new TPostScript("clct_effic.eps", 113);
4596     eps7->NewPage();
4597     c1->Clear();
4598     c1->cd(0);
4599     title = new TPaveLabel(0.1, 0.94, 0.9, 0.98, "CLCT efficiency vs #eta");
4600     title->SetFillColor(10);
4601     title->Draw();
4602     gStyle->SetOptStat(0);
4603     pad[page]->Draw();
4604     pad[page]->Divide(2, 5);
4605     TH1F* hCLCTEffVsEtaCsc[CSC_TYPES];
4606     for (int idh = 0; idh < max_idh; idh++) {
4607       if (!plotME1A && idh == 3)
4608         continue;
4609       hCLCTEffVsEtaCsc[idh] = (TH1F*)hEfficHitsEtaCsc[idh]->Clone();
4610       hCLCTEffVsEtaCsc[idh]->Divide(hEfficCLCTEtaCsc[idh], hEfficHitsEtaCsc[idh], 1., 1., "B");
4611       hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitle("#eta");
4612       hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleOffset(0.8);
4613       hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetTitleSize(0.07);  // default=0.05
4614       hCLCTEffVsEtaCsc[idh]->GetXaxis()->SetLabelSize(0.10);  // default=0.04
4615       hCLCTEffVsEtaCsc[idh]->GetYaxis()->SetLabelSize(0.10);
4616       hCLCTEffVsEtaCsc[idh]->SetLabelOffset(0.012, "XY");
4617       hCLCTEffVsEtaCsc[idh]->SetMinimum(0.50);
4618       hCLCTEffVsEtaCsc[idh]->SetMaximum(1.05);
4619       hCLCTEffVsEtaCsc[idh]->SetTitle(csc_type[idh].c_str());
4620       hCLCTEffVsEtaCsc[idh]->SetLineWidth(2);
4621       hCLCTEffVsEtaCsc[idh]->SetLineColor(4);
4622       pad[page]->cd(idh + 1);
4623       gPad->SetGrid(1);
4624       hCLCTEffVsEtaCsc[idh]->Draw();
4625     }
4626     page++;
4627     c1->Update();
4628     eps7->Close();
4629   }
4630   delete c1;
4631   delete c2;
4632 }
4633
4634 GlobalPoint CSCTriggerPrimitivesReader::getGlobalPosition(unsigned int rawId, int keyWg, int keyHS) const {
4635   // taken from https://github.com/cms-sw/cmssw/blob/dc9f78b6af4ad56c9342cf14041b6485a60b0691/L1Trigger/CSCTriggerPrimitives/src/CSCTriggerPrimitivesReaderGEM.cc
4636   CSCDetId cscId = CSCDetId(rawId);
4637   int ring = cscId.ring();
4638   if (cscId.station() == 1 and cscId.ring() == 1 and (lut_wg_vs_hs_me1b[keyWg][0] < 0 || keyHS >= 128)) {
4639     ring = 4;
4640     if (keyHS >= 128)
4641       keyHS = keyHS - 128;
4642   }
4643
4644   CSCDetId key_id(cscId.endcap(), cscId.station(), ring, cscId.chamber(), CSCConstants::KEY_CLCT_LAYER);
4645   auto cscChamber = geom_->chamber(cscId);
4646   float fractional_strip = 0.5 * (keyHS + 1) - 0.25;
4647   auto layer_geo = cscChamber->layer(CSCConstants::KEY_CLCT_LAYER)->geometry();
4648   // LCT::getKeyWG() also starts from 0
4649   float wire = layer_geo->middleWireOfGroup(keyWg + 1);
4650   LocalPoint csc_intersect = layer_geo->intersectionOfStripAndWire(fractional_strip, wire);
4651   GlobalPoint csc_gp = geom_->idToDet(key_id)->surface().toGlobal(csc_intersect);
4652   return csc_gp;
4653 }
4654
4655 bool CSCTriggerPrimitivesReader::doesALCTCrossCLCT(CSCDetId id, int key_wg, int key_hs) const {
4656   bool isME11 = (id.station() == 1 and (id.ring() == 1 or id.ring() == 4));
4657   if (not isME11)
4658     return true;
4659
4660   int theEndcap = (id.endcap() == 1) ? 1 : 2;
4661
4662   if (key_hs > CSCConstants::MAX_HALF_STRIP_ME1B) {
4663     key_hs = key_hs - CSCConstants::MAX_HALF_STRIP_ME1B - 1;  //convert it from 128-223 -> 0-95
4664     if (!gangedME1a) {
4665       // wrap around ME11 HS number for -z endcap
4666       if (theEndcap == 2)
4667         key_hs = CSCConstants::MAX_HALF_STRIP_ME1A_UNGANGED - key_hs;
4668       if (key_hs >= lut_wg_vs_hs_me1a[key_wg][0] && key_hs <= lut_wg_vs_hs_me1a[key_wg][1])
4669         return true;
4670       return false;
4671     } else {
4672       if (theEndcap == 2)
4673         key_hs = CSCConstants::MAX_HALF_STRIP_ME1A_GANGED - key_hs;
4674       if (key_hs >= lut_wg_vs_hs_me1ag[key_wg][0] && key_hs <= lut_wg_vs_hs_me1ag[key_wg][1])
4675         return true;
4676       return false;
4677     }
4678   }
4679   if (key_hs <= CSCConstants::MAX_HALF_STRIP_ME1B) {
4680     if (theEndcap == 2)
4681       key_hs = CSCConstants::MAX_HALF_STRIP_ME1B - key_hs;
4682     if (key_hs >= lut_wg_vs_hs_me1b[key_wg][0] && key_hs <= lut_wg_vs_hs_me1b[key_wg][1])
4683       return true;
4684   }
4685   return false;
4686 }
4687
4688 // Returns chamber type (0-9) according to the station and ring number
4689 int CSCTriggerPrimitivesReader::getCSCType(const CSCDetId& id) {
4690   int type = -999;
4691
4692   if (id.station() == 1) {
4693     type = (id.triggerCscId() - 1) / 3;
4694     if (id.ring() == 4) {
4695       type = 3;
4696     }
4697   } else {  // stations 2-4
4698     type = 3 + id.ring() + 2 * (id.station() - 2);
4699   }
4700
4701   assert(type >= 0 && type < CSC_TYPES);
4702   return type;
4703 }
4704
4705 // Returns halfstrips-per-radian for different CSC types
4706 double CSCTriggerPrimitivesReader::getHsPerRad(const int idh) { return (NCHAMBERS[idh] * MAX_HS[idh] / TWOPI); }
4707
4708 DEFINE_FWK_MODULE(CSCTriggerPrimitivesReader);
4709 //-------------------------------------------------