CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/RecoJets/JetProducers/src/JetMuonHitsIDHelper.cc

Go to the documentation of this file.
00001 #include "RecoJets/JetProducers/interface/JetMuonHitsIDHelper.h"
00002 
00003 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 #include "FWCore/Framework/interface/ESHandle.h"
00006 
00007 #include "DataFormats/JetReco/interface/JetTracksAssociation.h"
00008 // #include "TrackingTools/TrackAssociator/interface/MuonDetIdAssociator.h"
00009 #include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h"
00010 #include "TrackingTools/TrackAssociator/interface/DetIdAssociator.h"
00011 #include "DataFormats/DTRecHit/interface/DTRecHitCollection.h"
00012 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
00013 #include "DataFormats/CSCRecHit/interface/CSCRecHit2DCollection.h"
00014 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
00015 #include "DataFormats/JetReco/interface/CaloJet.h"
00016 #include "DataFormats/TrackReco/interface/Track.h"
00017 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
00018 #include "DataFormats/MuonDetId/interface/DTWireId.h"
00019 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
00020 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
00021 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
00022 #include "DataFormats/MuonReco/interface/Muon.h"
00023 #include "DataFormats/MuonReco/interface/MuonFwd.h"
00024 #include "DataFormats/Math/interface/deltaR.h"
00025 
00026 #include "TMath.h"
00027 #include <vector>
00028 #include <iostream>
00029 
00030 using namespace std;
00031 
00032 
00033 reco::helper::JetMuonHitsIDHelper::JetMuonHitsIDHelper( edm::ParameterSet const & pset )
00034 {
00035   isRECO_ = true; // This will be "true" initially, then if the product isn't found, set to false once
00036   numberOfHits1RPC_ = 0;
00037   numberOfHits2RPC_ = 0;
00038   numberOfHits3RPC_ = 0;
00039   numberOfHits4RPC_ = 0;
00040   numberOfHitsRPC_ = 0;
00041   rpcRecHits_ = pset.getParameter<edm::InputTag>("rpcRecHits");
00042  
00043 }
00044 
00045 
00046 
00047 
00048 void reco::helper::JetMuonHitsIDHelper::calculate( const edm::Event& event, const edm::EventSetup & iSetup, 
00049                                                    const reco::Jet &jet, const int iDbg )
00050 {
00051 
00052   // initialize
00053   numberOfHits1RPC_ = 0;
00054   numberOfHits2RPC_ = 0;
00055   numberOfHits3RPC_ = 0;
00056   numberOfHits4RPC_ = 0;
00057   numberOfHitsRPC_ = 0;
00058 
00059 
00060   if ( isRECO_ ) { // This will be "true" initially, then if the product isn't found, set to false once
00061 
00062     // Get tracking geometry
00063     edm::ESHandle<GlobalTrackingGeometry> trackingGeometry;
00064     iSetup.get<GlobalTrackingGeometryRecord> ().get(trackingGeometry);
00065     
00066     //####READ RPC RecHits Collection########
00067     //#In config: RpcRecHits     = cms.InputTag("rpcRecHits")
00068     edm::Handle<RPCRecHitCollection> rpcRecHits_handle;
00069     event.getByLabel(rpcRecHits_, rpcRecHits_handle);
00070 
00071 
00072     if ( ! rpcRecHits_handle.isValid()  ) {
00073       // don't throw exception if not running on RECO
00074       edm::LogWarning("DataNotAvailable") << "JetMuonHitsIDHelper will not be run at all, this is not a RECO file.";
00075       isRECO_ = false;
00076       return;
00077     }
00078 
00079     //####calculate rpc  variables for each jet########
00080     
00081     for (  RPCRecHitCollection::const_iterator itRPC = rpcRecHits_handle->begin(),
00082              itRPCEnd = rpcRecHits_handle->end();
00083            itRPC != itRPCEnd; ++itRPC) {
00084       RPCRecHit const & hit = *itRPC;
00085       DetId detid = hit.geographicalId();
00086       LocalPoint lp = hit.localPosition();
00087       const GeomDet* gd = trackingGeometry->idToDet(detid);
00088       GlobalPoint gp = gd->toGlobal(lp);
00089       double dR2 = reco::deltaR(jet.eta(), jet.phi(), 
00090                                 static_cast<double>( gp.eta() ), static_cast<double>(gp.phi()) );
00091       if (dR2 < 0.5) {
00092         RPCDetId rpcChamberId = (RPCDetId) detid;
00093         numberOfHitsRPC_++;
00094         if (rpcChamberId.station() == 1)
00095           numberOfHits1RPC_++;
00096         if (rpcChamberId.station() == 2)
00097           numberOfHits2RPC_++;
00098         if (rpcChamberId.station() == 3)
00099           numberOfHits3RPC_++;
00100         if (rpcChamberId.station() == 4)
00101           numberOfHits4RPC_++;
00102       }
00103     }
00104   }
00105 }