/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/DPGAnalysis/Skims/src/RPCRecHitFilter.cc

Go to the documentation of this file.

// Orso Iorio, INFN Napoli 
//
//This module is a filter based on the number of RPC Hits
//

#include "DPGAnalysis/Skims/src/RPCRecHitFilter.h"
//#include "../interface/RecHitFilter.h"

#include "FWCore/MessageLogger/interface/MessageLogger.h"
//#include "DataFormats/RPCDigi/interface/RPCDigi.h"
//#include "DataFormats/RPCDigi/interface/RPCDigiCollection.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"

#include "FWCore/Framework/interface/MakerMacros.h"
#include "DataFormats/Common/interface/Handle.h"
#include <DataFormats/RPCRecHit/interface/RPCRecHit.h>
#include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"

#include <Geometry/DTGeometry/interface/DTGeometry.h>
#include <Geometry/RPCGeometry/interface/RPCGeometry.h>
#include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
#include "MagneticField/Engine/interface/MagneticField.h"
#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
#include <Geometry/Records/interface/MuonGeometryRecord.h>
#include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
#include "Geometry/CommonDetUnit/interface/GeomDet.h"
#include "DataFormats/DetId/interface/DetId.h"
#include "DataFormats/MuonDetId/interface/MuonSubdetId.h"


#include <FWCore/Utilities/interface/Exception.h>
#include "DataFormats/Provenance/interface/Timestamp.h"


#include "DataFormats/L1Trigger/interface/L1MuonParticleFwd.h"
#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuRegionalCand.h"
#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTReadoutCollection.h"

#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTCand.h"
#include "CondFormats/L1TObjects/interface/L1GtParameters.h"
#include "CondFormats/DataRecord/interface/L1GtParametersRcd.h"


/*#include "DataFormats/L1Trigger/interface/L1MuonParticleFwd.h"
#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuRegionalCand.h"

#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTReadoutCollection.h"

#include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTCand.h"
*/

using namespace reco;

typedef std::vector<L1MuRegionalCand> L1MuRegionalCandCollection;

#include <sys/time.h>
#include <algorithm>
#include <memory>
#include <cmath>
#include "math.h"
#include "TFile.h"
#include "TMath.h"
#include "TTree.h"

#include "TDirectory.h"
#include "TFile.h"
#include "TTree.h"
#include <cstdlib>
#include <stdio.h>
#include "stdlib.h"
#include <string>
#include <memory>
#include<cmath>

using namespace edm;
using namespace reco;
using namespace std;

RPCRecHitFilter::RPCRecHitFilter(const edm::ParameterSet& iConfig){


  LogTrace("RPCEffTrackExtrapolation") <<"Dentro Costruttore"<<std::endl;

  RPCDataLabel = iConfig.getUntrackedParameter<std::string>("rpcRecHitLabel");
  
  //  RPCRecHits = iConfig.getParameter< edm::InputTag >("RPCRecHits");
  
 
  centralBX_ = iConfig.getUntrackedParameter<int>("CentralBunchCrossing",0);  
  BXWindow_ = iConfig.getUntrackedParameter<int>("BunchCrossingWindow",9999);
  
  minHits_ = iConfig.getUntrackedParameter<int>("minimumNumberOfHits",0); 
  hitsInStations_ =  iConfig.getUntrackedParameter<int>("HitsInStation",0); 


  Debug_ = iConfig.getUntrackedParameter<bool>("Debug",false);
  Verbose_ = iConfig.getUntrackedParameter<bool>("Verbose",false);


  Barrel_  = iConfig.getUntrackedParameter<bool>("UseBarrel",true);
  EndcapPositive_  = iConfig.getUntrackedParameter<bool>("UseEndcapPositive",true);
  EndcapNegative_  = iConfig.getUntrackedParameter<bool>("UseEndcapNegative",true);

  cosmicsVeto_  = iConfig.getUntrackedParameter<bool>("CosmicsVeto",false);



}

void RPCRecHitFilter::beginJob()
{

}

bool RPCRecHitFilter::filter( edm::Event& iEvent, const edm::EventSetup& iSetup){

  
  edm::ESHandle<RPCGeometry> rpcGeo;
  iSetup.get<MuonGeometryRecord>().get(rpcGeo);
  
  edm::Handle<RPCRecHitCollection> rpcHits;
  iEvent.getByLabel(RPCDataLabel,rpcHits);
  
  
  ESHandle<GlobalTrackingGeometry> theTrackingGeometry;
  iSetup.get<GlobalTrackingGeometryRecord>().get(theTrackingGeometry);
    
  
  std::map <int  , int > numberOfRecHitsBarrel;
  std::map <int  , int > numberOfDigisBarrel;
  std::map <int  , int > numberOfRecHitsEndcap;
  std::map <int  , int > numberOfDigisEndcap;
  
  std::map <pair < int  , int > , std::vector<RPCDetId> > numberOfRecHitsSameWheelSameSector;
  std::map <pair < int  , int > , std::vector<RPCDetId> > numberOfDigisSameWheelSameSector;
  std::map <pair < int  , int > , std::vector<RPCDetId> > numberOfHitsSameDiskSectorPositive;
  std::map <pair < int  , int > , std::vector<RPCDetId> > numberOfHitsSameDiskSectorNegative;
  
  
  const std::vector<  RPCRoll * > rls = rpcGeo->rolls();
  
  bool condition = true;
  
  
  int nRecHits = 0;
  int nBarrel = 0;
  int nEndcap = 0;
  
  
  for (int i = 0; i < (int) rls.size() ; ++i){
    
    RPCDetId did =  rls[i]->id();
    
    
    RPCRecHitCollection::range rpcRecHitRange = rpcHits->get(did);
    RPCRecHitCollection::const_iterator RecHitsIt;                                       
    
        
    for(RecHitsIt = rpcRecHitRange.first; RecHitsIt!=rpcRecHitRange.second; ++RecHitsIt){
      
      //std::cout<< " roll is "<< did << " bx recHit is " << RecHitsIt->BunchX()<< " event number " <<eventNumber  <<std::endl;
      
      if(did.region()==0){
        
        if (cosmicsVeto_ ){
          
          for(int u = -1;u<=1;++u){
            for(int t = -1;t<=1;++t){
              numberOfRecHitsSameWheelSameSector[pair< int, int >( did.ring()+u,did.sector()+t)].push_back(did);
            }
          }
        }
        
        else{
          if( did.station()==1) {
            for(int u = -1;u<=1;++u){
              for(int t = -1;t<=1;++t){
                numberOfRecHitsSameWheelSameSector[pair< int, int >(did.ring()+u,did.sector()+t)].push_back(did);
              }
            }
          }
        }
        
        ++numberOfRecHitsBarrel[did.ring()];
        ++nBarrel;
        
      } 
      
      else {
        if(did.region()== -1) { 
          for(int t = -1;t<=1;++t){
                numberOfHitsSameDiskSectorNegative[pair<int,int>(did.ring(),did.sector()+t)].push_back(did); 
          }
        }

        if(did.region()==  1){
          for(int t = -1;t<=1;++t){
            numberOfHitsSameDiskSectorPositive[pair<int,int>(did.ring(),did.sector()+t)].push_back(did);
          }
        }           
        ++numberOfRecHitsEndcap[did.station()]; 
        
        ++nEndcap;
      }
      
      ++nRecHits;
      
      break;
    }
        
    
  }
  
  
      bool cond1 = false;
      bool cond2 = false;
      bool cond3 = false;

      std::map < int,  bool >vectorBarrelCands;
      std::map < int,  bool >vectorEndcapCandsPositive;
      std::map < int, bool >vectorEndcapCandsNegative;
 
      bool veto = false; 
      
      // barrel
      for ( std::map <pair < int , int > , std::vector < RPCDetId> >::const_iterator iter = numberOfRecHitsSameWheelSameSector.begin();iter != numberOfRecHitsSameWheelSameSector.end();++iter){

        vectorBarrelCands[1] = false; 
        vectorBarrelCands[2] = false; 

        if(iter->second.size()>1){
          for(size_t i = 0; i < iter->second.size(); ++i){
            if(iter->second[i].layer()==1 && iter->second[i].station() == 1)vectorBarrelCands[0] = true; 
            if(iter->second[i].layer()==2 && iter->second[i].station() == 1)vectorBarrelCands[1] = true; 
            if(cosmicsVeto_)if(iter->second[i].station() > 2){ veto = true; 
            vectorBarrelCands[1] = false;vectorBarrelCands[2] = false ;
            break;
            } 
          } 
        }
        
        if(( vectorBarrelCands[0] && vectorBarrelCands[1])) {
          cond1 = true;
          break;
        } 
      } 
 
 
      // endcap positive
      for ( std::map <pair < int ,int >, std::vector < RPCDetId> >::const_iterator iter = numberOfHitsSameDiskSectorPositive.begin(); iter != numberOfHitsSameDiskSectorPositive.end(); ++iter){

        vectorEndcapCandsPositive[1] = false; 
        vectorEndcapCandsPositive[2] = false; 
        vectorEndcapCandsPositive[3] = false; 

        
        if(iter->second.size()>1){
          for(size_t i = 0; i < iter->second.size(); ++i){
            if(iter->second[i].station()==1)vectorEndcapCandsPositive[1] = true; 
            if(iter->second[i].station()==2)vectorEndcapCandsPositive[2] = true; 
            if(iter->second[i].station()==3)vectorEndcapCandsPositive[3] = true; 
          } 
        }
        
        if(( (vectorEndcapCandsPositive[1] && vectorEndcapCandsPositive[2]) ||
             (vectorEndcapCandsPositive[1] && vectorEndcapCandsPositive[3]) ||
             (vectorEndcapCandsPositive[2] && vectorEndcapCandsPositive[3])  )) {
          cond2 = true;
          break;
        }       
        
      } 

      // endcap negative
      for ( std::map <pair < int ,int >, std::vector < RPCDetId> >::const_iterator iter = numberOfHitsSameDiskSectorNegative.begin(); iter != numberOfHitsSameDiskSectorNegative.end(); ++iter){

        vectorEndcapCandsNegative[1] = false; 
        vectorEndcapCandsNegative[2] = false; 
        vectorEndcapCandsNegative[3] = false; 

        
        if(iter->second.size()>1){
          for(size_t i = 0; i < iter->second.size(); ++i){
            if(iter->second[i].station()==1)vectorEndcapCandsNegative[1] = true; 
            if(iter->second[i].station()==2)vectorEndcapCandsNegative[2] = true; 
            if(iter->second[i].station()==3)vectorEndcapCandsNegative[3] = true; 
          } 
        }
        
        if(( (vectorEndcapCandsNegative[1] && vectorEndcapCandsNegative[2]) ||
             (vectorEndcapCandsNegative[1] && vectorEndcapCandsNegative[3]) ||
             (vectorEndcapCandsNegative[2] && vectorEndcapCandsNegative[3])  )) {
          cond3 = true;
          break;
        }       
        
      } 


      condition = condition && ( nBarrel + nEndcap >= minHits_ );


      cond1 = Barrel_ && cond1;
      cond2 = EndcapPositive_ && cond2;
      cond3 = EndcapNegative_ && cond3;

      bool condition2 = (cond1 || cond2 || cond3);
      if(Barrel_ || EndcapPositive_ || EndcapNegative_)condition = condition && condition2;
      

      
      return  condition;      
      
}

void RPCRecHitFilter::endJob(){

}

RPCRecHitFilter::~RPCRecHitFilter(){}

DEFINE_FWK_MODULE(RPCRecHitFilter);