Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-07 22:41:51

0001 // Orso Iorio, INFN Napoli
0002 //
0003 //This module is a filter based on the number of RPC Hits
0004 //
0005 
0006 #include "DPGAnalysis/Skims/src/RPCRecHitFilter.h"
0007 
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "DataFormats/Common/interface/Handle.h"
0013 #include <DataFormats/RPCRecHit/interface/RPCRecHit.h>
0014 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
0015 
0016 #include <Geometry/DTGeometry/interface/DTGeometry.h>
0017 #include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
0018 #include "MagneticField/Engine/interface/MagneticField.h"
0019 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0020 #include <Geometry/Records/interface/MuonGeometryRecord.h>
0021 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
0022 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0023 #include "DataFormats/DetId/interface/DetId.h"
0024 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0025 
0026 #include <FWCore/Utilities/interface/Exception.h>
0027 #include "DataFormats/Provenance/interface/Timestamp.h"
0028 
0029 /////Trigger
0030 
0031 #include "DataFormats/L1Trigger/interface/L1MuonParticleFwd.h"
0032 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuRegionalCand.h"
0033 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTReadoutCollection.h"
0034 
0035 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTCand.h"
0036 #include "CondFormats/L1TObjects/interface/L1GtParameters.h"
0037 #include "CondFormats/DataRecord/interface/L1GtParametersRcd.h"
0038 
0039 using namespace reco;
0040 
0041 typedef std::vector<L1MuRegionalCand> L1MuRegionalCandCollection;
0042 
0043 #include <sys/time.h>
0044 #include <algorithm>
0045 #include <memory>
0046 #include <cmath>
0047 #include <cmath>
0048 #include "TFile.h"
0049 #include "TMath.h"
0050 #include "TTree.h"
0051 
0052 #include "TDirectory.h"
0053 #include "TFile.h"
0054 #include "TTree.h"
0055 #include <cstdlib>
0056 #include <cstdio>
0057 #include <cstdlib>
0058 #include <string>
0059 #include <memory>
0060 #include <cmath>
0061 
0062 using namespace edm;
0063 using namespace reco;
0064 using namespace std;
0065 
0066 RPCRecHitFilter::RPCRecHitFilter(const edm::ParameterSet& iConfig)
0067     : rpcGeomToken_(esConsumes()), trackingGeoToken_(esConsumes()) {
0068   LogTrace("RPCEffTrackExtrapolation") << "inside constructor" << std::endl;
0069 
0070   RPCDataLabel_ = iConfig.getUntrackedParameter<std::string>("rpcRecHitLabel");
0071   rpcRecHitToken_ = consumes<RPCRecHitCollection>(edm::InputTag(RPCDataLabel_));
0072 
0073   centralBX_ = iConfig.getUntrackedParameter<int>("CentralBunchCrossing", 0);
0074   BXWindow_ = iConfig.getUntrackedParameter<int>("BunchCrossingWindow", 9999);
0075 
0076   minHits_ = iConfig.getUntrackedParameter<int>("minimumNumberOfHits", 0);
0077   hitsInStations_ = iConfig.getUntrackedParameter<int>("HitsInStation", 0);
0078 
0079   Debug_ = iConfig.getUntrackedParameter<bool>("Debug", false);
0080   Verbose_ = iConfig.getUntrackedParameter<bool>("Verbose", false);
0081 
0082   Barrel_ = iConfig.getUntrackedParameter<bool>("UseBarrel", true);
0083   EndcapPositive_ = iConfig.getUntrackedParameter<bool>("UseEndcapPositive", true);
0084   EndcapNegative_ = iConfig.getUntrackedParameter<bool>("UseEndcapNegative", true);
0085 
0086   cosmicsVeto_ = iConfig.getUntrackedParameter<bool>("CosmicsVeto", false);
0087 }
0088 
0089 bool RPCRecHitFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0090   edm::ESHandle<RPCGeometry> rpcGeo = iSetup.getHandle(rpcGeomToken_);
0091   edm::ESHandle<GlobalTrackingGeometry> theTrackingGeometry = iSetup.getHandle(trackingGeoToken_);
0092 
0093   edm::Handle<RPCRecHitCollection> rpcHits = iEvent.getHandle(rpcRecHitToken_);
0094 
0095   std::map<int, int> numberOfRecHitsBarrel;
0096   std::map<int, int> numberOfDigisBarrel;
0097   std::map<int, int> numberOfRecHitsEndcap;
0098   std::map<int, int> numberOfDigisEndcap;
0099 
0100   std::map<pair<int, int>, std::vector<RPCDetId> > numberOfRecHitsSameWheelSameSector;
0101   std::map<pair<int, int>, std::vector<RPCDetId> > numberOfDigisSameWheelSameSector;
0102   std::map<pair<int, int>, std::vector<RPCDetId> > numberOfHitsSameDiskSectorPositive;
0103   std::map<pair<int, int>, std::vector<RPCDetId> > numberOfHitsSameDiskSectorNegative;
0104 
0105   const std::vector<const RPCRoll*> rls = rpcGeo->rolls();
0106 
0107   bool condition = true;
0108 
0109   int nRecHits = 0;
0110   int nBarrel = 0;
0111   int nEndcap = 0;
0112 
0113   /////BEGIN LOOP ON THE ROLLS
0114   for (int i = 0; i < (int)rls.size(); ++i) {
0115     RPCDetId did = rls[i]->id();
0116 
0117     /// LOOP OVER THE RECHITS
0118     RPCRecHitCollection::range rpcRecHitRange = rpcHits->get(did);
0119     RPCRecHitCollection::const_iterator RecHitsIt;
0120 
0121     for (RecHitsIt = rpcRecHitRange.first; RecHitsIt != rpcRecHitRange.second; ++RecHitsIt) {
0122       //std::cout<< " roll is "<< did << " bx recHit is " << RecHitsIt->BunchX()<< " event number " <<eventNumber  <<std::endl;
0123 
0124       if (did.region() == 0) {
0125         if (cosmicsVeto_) {
0126           for (int u = -1; u <= 1; ++u) {
0127             for (int t = -1; t <= 1; ++t) {
0128               numberOfRecHitsSameWheelSameSector[pair<int, int>(did.ring() + u, did.sector() + t)].push_back(did);
0129             }
0130           }
0131         }
0132 
0133         else {
0134           if (did.station() == 1) {
0135             for (int u = -1; u <= 1; ++u) {
0136               for (int t = -1; t <= 1; ++t) {
0137                 numberOfRecHitsSameWheelSameSector[pair<int, int>(did.ring() + u, did.sector() + t)].push_back(did);
0138               }
0139             }
0140           }
0141         }
0142 
0143         ++numberOfRecHitsBarrel[did.ring()];
0144         ++nBarrel;
0145 
0146       }
0147 
0148       else {
0149         if (did.region() == -1) {
0150           for (int t = -1; t <= 1; ++t) {
0151             numberOfHitsSameDiskSectorNegative[pair<int, int>(did.ring(), did.sector() + t)].push_back(did);
0152           }
0153         }
0154 
0155         if (did.region() == 1) {
0156           for (int t = -1; t <= 1; ++t) {
0157             numberOfHitsSameDiskSectorPositive[pair<int, int>(did.ring(), did.sector() + t)].push_back(did);
0158           }
0159         }
0160         ++numberOfRecHitsEndcap[did.station()];
0161 
0162         ++nEndcap;
0163       }
0164 
0165       ++nRecHits;
0166 
0167       break;
0168     }
0169   }
0170   /////END OF THE LOOP ON THE ROLLS
0171 
0172   ////CONDITION IS THAT THERE ARE TWO HITS IN RB1in and RB1out OF TWO ADJACENT SECTORS AND WHEELS
0173   ////OR THAT IN TWO DISKS SAME OF THE ENDCAP THERE ARE TWO HITS IN ADJACENT SECTORS IN THE SAME RING
0174   bool cond1 = false;
0175   bool cond2 = false;
0176   bool cond3 = false;
0177 
0178   std::map<int, bool> vectorBarrelCands;
0179   std::map<int, bool> vectorEndcapCandsPositive;
0180   std::map<int, bool> vectorEndcapCandsNegative;
0181 
0182   // barrel
0183   for (std::map<pair<int, int>, std::vector<RPCDetId> >::const_iterator iter =
0184            numberOfRecHitsSameWheelSameSector.begin();
0185        iter != numberOfRecHitsSameWheelSameSector.end();
0186        ++iter) {
0187     vectorBarrelCands[1] = false;
0188     vectorBarrelCands[2] = false;
0189 
0190     if (iter->second.size() > 1) {
0191       for (size_t i = 0; i < iter->second.size(); ++i) {
0192         if (iter->second[i].layer() == 1 && iter->second[i].station() == 1)
0193           vectorBarrelCands[0] = true;
0194         if (iter->second[i].layer() == 2 && iter->second[i].station() == 1)
0195           vectorBarrelCands[1] = true;
0196         if (cosmicsVeto_)
0197           if (iter->second[i].station() > 2) {
0198             vectorBarrelCands[1] = false;
0199             vectorBarrelCands[2] = false;
0200             break;
0201           }
0202       }
0203     }
0204 
0205     if ((vectorBarrelCands[0] && vectorBarrelCands[1])) {
0206       cond1 = true;
0207       break;
0208     }
0209   }
0210 
0211   // endcap positive
0212   for (std::map<pair<int, int>, std::vector<RPCDetId> >::const_iterator iter =
0213            numberOfHitsSameDiskSectorPositive.begin();
0214        iter != numberOfHitsSameDiskSectorPositive.end();
0215        ++iter) {
0216     vectorEndcapCandsPositive[1] = false;
0217     vectorEndcapCandsPositive[2] = false;
0218     vectorEndcapCandsPositive[3] = false;
0219 
0220     if (iter->second.size() > 1) {
0221       for (size_t i = 0; i < iter->second.size(); ++i) {
0222         if (iter->second[i].station() == 1)
0223           vectorEndcapCandsPositive[1] = true;
0224         if (iter->second[i].station() == 2)
0225           vectorEndcapCandsPositive[2] = true;
0226         if (iter->second[i].station() == 3)
0227           vectorEndcapCandsPositive[3] = true;
0228       }
0229     }
0230 
0231     if (((vectorEndcapCandsPositive[1] && vectorEndcapCandsPositive[2]) ||
0232          (vectorEndcapCandsPositive[1] && vectorEndcapCandsPositive[3]) ||
0233          (vectorEndcapCandsPositive[2] && vectorEndcapCandsPositive[3]))) {
0234       cond2 = true;
0235       break;
0236     }
0237   }
0238 
0239   // endcap negative
0240   for (std::map<pair<int, int>, std::vector<RPCDetId> >::const_iterator iter =
0241            numberOfHitsSameDiskSectorNegative.begin();
0242        iter != numberOfHitsSameDiskSectorNegative.end();
0243        ++iter) {
0244     vectorEndcapCandsNegative[1] = false;
0245     vectorEndcapCandsNegative[2] = false;
0246     vectorEndcapCandsNegative[3] = false;
0247 
0248     if (iter->second.size() > 1) {
0249       for (size_t i = 0; i < iter->second.size(); ++i) {
0250         if (iter->second[i].station() == 1)
0251           vectorEndcapCandsNegative[1] = true;
0252         if (iter->second[i].station() == 2)
0253           vectorEndcapCandsNegative[2] = true;
0254         if (iter->second[i].station() == 3)
0255           vectorEndcapCandsNegative[3] = true;
0256       }
0257     }
0258 
0259     if (((vectorEndcapCandsNegative[1] && vectorEndcapCandsNegative[2]) ||
0260          (vectorEndcapCandsNegative[1] && vectorEndcapCandsNegative[3]) ||
0261          (vectorEndcapCandsNegative[2] && vectorEndcapCandsNegative[3]))) {
0262       cond3 = true;
0263       break;
0264     }
0265   }
0266 
0267   condition = condition && (nBarrel + nEndcap >= minHits_);
0268 
0269   cond1 = Barrel_ && cond1;
0270   cond2 = EndcapPositive_ && cond2;
0271   cond3 = EndcapNegative_ && cond3;
0272 
0273   bool condition2 = (cond1 || cond2 || cond3);
0274   if (Barrel_ || EndcapPositive_ || EndcapNegative_)
0275     condition = condition && condition2;
0276 
0277   return condition;
0278 }
0279 
0280 DEFINE_FWK_MODULE(RPCRecHitFilter);