#include <MuonCosmicCompatibilityFiller.h>

Public Member Functions
reco::MuonCosmicCompatibility	fillCompatibility (const reco::Muon &muon, edm::Event &, const edm::EventSetup &)
	fill cosmic compatibility variables
	MuonCosmicCompatibilityFiller (const edm::ParameterSet &)
	~MuonCosmicCompatibilityFiller ()
Private Member Functions
unsigned int	backToBack2LegCosmic (const edm::Event &, const reco::Muon &) const
	return cosmic-likeness based on presence of a track in opp side: 0 == no matching opp tracks
bool	checkMuonID (const reco::Muon &) const
	tag a muon as cosmic based on the muonID information
bool	checkMuonSegments (const reco::Muon &muon) const
	tag a muon as cosmic based on segment compatibility and the number of segment matches
float	combinedCosmicID (const edm::Event &, const edm::EventSetup &iSetup, const reco::Muon &, bool CheckMuonID, bool checkVertex) const
	combined cosmic-likeness: 0 == not cosmic-like
unsigned int	eventActivity (const edm::Event &, const reco::Muon &) const
	returns cosmic-likeness based on the event activity information: tracker track multiplicity and vertex quality. 0 == cosmic-like
bool	isOverlappingMuon (const edm::Event &, const edm::EventSetup &iSetup, const reco::Muon &) const
	returns cosmic-likeness based on overlap with traversing cosmic muon (only muon/STA hits are used)
float	muonTiming (const edm::Event &iEvent, const reco::Muon &muon, bool isLoose) const
	check muon time (DT and CSC) information: 0 == prompt-like
unsigned int	nMuons (const edm::Event &) const
	get number of muons in the vent
unsigned int	pvMatches (const edm::Event &, const reco::Muon &, bool) const
	return cosmic-likeness based on the 2D impact parameters (dxy, dz wrt to PV). 0 == cosmic-like
Private Attributes
double	angleThreshold_
double	corrTimeNeg_
double	corrTimePos_
double	deltaPt_
double	hIpTrdxy_
double	hIpTrvProb_
edm::InputTag	inputCosmicMuonCollection_
std::vector< edm::InputTag >	inputMuonCollections_
std::vector< edm::InputTag >	inputTrackCollections_
edm::InputTag	inputVertexCollection_
double	ipThreshold_
double	largedxy_
double	largedxyMult_
double	maxdxyLoose_
double	maxdxyLooseMult_
double	maxdxyTight_
double	maxdxyTightMult_
double	maxdzLoose_
double	maxdzLooseMult_
double	maxdzTight_
double	maxdzTightMult_
double	maxvertRho_
double	maxvertZ_
double	minvProb_
int	nChamberMatches_
unsigned int	nTrackThreshold_
double	offTimeNegLoose_
double	offTimeNegLooseMult_
double	offTimeNegTight_
double	offTimeNegTightMult_
double	offTimePosLoose_
double	offTimePosLooseMult_
double	offTimePosTight_
double	offTimePosTightMult_
double	segmentComp_
MuonServiceProxy *	service_
double	sharedFrac_
int	sharedHits_

Detailed Description

Description: class for cosmic muon identification

Date:: 2010/09/13 12:44:08

Revision:: 1.2

Author:: : A. Everett, Purdue University; : A. Svyatkovskiy, Purdue University; : H.D. Yoo, Purdue University

Description: class for cosmic muon identification

Date:: 2011/10/28 22:10:34

Revision:: 1.4

Author:: : A. Everett, Purdue University; : A. Svyatkovskiy, Purdue University; : H.D. Yoo, Purdue University

Definition at line 32 of file MuonCosmicCompatibilityFiller.h.

Constructor & Destructor Documentation

MuonCosmicCompatibilityFiller::MuonCosmicCompatibilityFiller ( const edm::ParameterSet & iConfig )

Definition at line 55 of file MuonCosmicCompatibilityFiller.cc.

References angleThreshold_, corrTimeNeg_, corrTimePos_, deltaPt_, edm::ParameterSet::getParameter(), hIpTrdxy_, hIpTrvProb_, ipThreshold_, largedxy_, largedxyMult_, maxdxyLoose_, maxdxyLooseMult_, maxdxyTight_, maxdxyTightMult_, maxdzLoose_, maxdzLooseMult_, maxdzTight_, maxdzTightMult_, maxvertRho_, maxvertZ_, minvProb_, MuonServiceProxy_cff::MuonServiceProxy, nChamberMatches_, offTimeNegLoose_, offTimeNegLooseMult_, offTimeNegTight_, offTimeNegTightMult_, offTimePosLoose_, offTimePosLooseMult_, offTimePosTight_, offTimePosTightMult_, segmentComp_, service_, sharedFrac_, and sharedHits_.

                                                                                          :
  inputMuonCollections_(iConfig.getParameter<std::vector<edm::InputTag> >("InputMuonCollections")),
  inputTrackCollections_(iConfig.getParameter<std::vector<edm::InputTag> >("InputTrackCollections")),
  inputCosmicMuonCollection_(iConfig.getParameter<edm::InputTag>("InputCosmicMuonCollection")),
  inputVertexCollection_(iConfig.getParameter<edm::InputTag>("InputVertexCollection")),
  service_(0)
{
  // service parameters
  edm::ParameterSet serviceParameters = iConfig.getParameter<edm::ParameterSet>("ServiceParameters");
  service_ = new MuonServiceProxy(serviceParameters);     
  
  //kinematic vars
  angleThreshold_ = iConfig.getParameter<double>("angleCut");
  deltaPt_ = iConfig.getParameter<double>("deltaPt");
  //time
  offTimePosTightMult_ = iConfig.getParameter<double>("offTimePosTightMult");
  offTimeNegTightMult_ = iConfig.getParameter<double>("offTimeNegTightMult");
  offTimePosTight_ = iConfig.getParameter<double>("offTimePosTight");
  offTimeNegTight_ = iConfig.getParameter<double>("offTimeNegTight");
  offTimePosLooseMult_ = iConfig.getParameter<double>("offTimePosLooseMult");
  offTimeNegLooseMult_ = iConfig.getParameter<double>("offTimeNegLooseMult");
  offTimePosLoose_ = iConfig.getParameter<double>("offTimePosLoose");
  offTimeNegLoose_ = iConfig.getParameter<double>("offTimeNegLoose");
  corrTimeNeg_ = iConfig.getParameter<double>("corrTimeNeg");
  corrTimePos_ = iConfig.getParameter<double>("corrTimePos");
  //rechits
  sharedHits_ = iConfig.getParameter<int>("sharedHits");
  sharedFrac_ = iConfig.getParameter<double>("sharedFrac");
  ipThreshold_ = iConfig.getParameter<double>("ipCut");
  //segment comp, matches
  nChamberMatches_ = iConfig.getParameter<int>("nChamberMatches");
  segmentComp_ = iConfig.getParameter<double>("segmentComp");
  //ip, vertex
  maxdzLooseMult_ = iConfig.getParameter<double>("maxdzLooseMult");
  maxdxyLooseMult_ = iConfig.getParameter<double>("maxdxyLooseMult");
  maxdzTightMult_ = iConfig.getParameter<double>("maxdzTightMult");
  maxdxyTightMult_ = iConfig.getParameter<double>("maxdxyTightMult");
  maxdzLoose_ = iConfig.getParameter<double>("maxdzLoose");
  maxdxyLoose_ = iConfig.getParameter<double>("maxdxyLoose");
  maxdzTight_ = iConfig.getParameter<double>("maxdzTight");
  maxdxyTight_ = iConfig.getParameter<double>("maxdxyTight");
  largedxyMult_ = iConfig.getParameter<double>("largedxyMult");
  largedxy_ = iConfig.getParameter<double>("largedxy");
  hIpTrdxy_ = iConfig.getParameter<double>("hIpTrdxy");
  hIpTrvProb_ = iConfig.getParameter<double>("hIpTrvProb");
  minvProb_ = iConfig.getParameter<double>("minvProb");
  maxvertZ_ = iConfig.getParameter<double>("maxvertZ");
  maxvertRho_ = iConfig.getParameter<double>("maxvertRho");
//  nTrackThreshold_ = iConfig.getParameter<unsigned int>("nTrackThreshold");
}

MuonCosmicCompatibilityFiller::~MuonCosmicCompatibilityFiller ( )

Definition at line 106 of file MuonCosmicCompatibilityFiller.cc.

References service_.

                                                              {
  if (service_) delete service_;
}

Member Function Documentation

unsigned int MuonCosmicCompatibilityFiller::backToBack2LegCosmic	(	const edm::Event &	iEvent,
		const reco::Muon &	muon
	)		const `[private]`

return cosmic-likeness based on presence of a track in opp side: 0 == no matching opp tracks

Definition at line 222 of file MuonCosmicCompatibilityFiller.cc.

References angleThreshold_, deltaPt_, muonid::findOppositeTrack(), edm::Event::getByLabel(), reco::Muon::innerTrack(), inputTrackCollections_, reco::Muon::isGlobalMuon(), reco::Muon::isStandAloneMuon(), reco::Muon::isTrackerMuon(), query::result, and reco::Muon::track().

Referenced by combinedCosmicID(), and fillCompatibility().

                                                                                                      {

  unsigned int result = 0; //no partners - collision
  reco::TrackRef track;
  if ( muon.isGlobalMuon()  )            track = muon.innerTrack();
  else if ( muon.isTrackerMuon() )       track = muon.track();
  else if ( muon.isStandAloneMuon() )    return false;

  for (unsigned int iColl = 0; iColl<inputTrackCollections_.size(); ++iColl){
    edm::Handle<reco::TrackCollection> trackHandle;
    iEvent.getByLabel(inputTrackCollections_[iColl],trackHandle);
    if (muonid::findOppositeTrack(trackHandle, *track, angleThreshold_, deltaPt_).isNonnull()) { 
      result++;
     }
   } //loop over track collections

  return result;
}

bool MuonCosmicCompatibilityFiller::checkMuonID ( const reco::Muon & imuon ) const [private]

tag a muon as cosmic based on the muonID information

Definition at line 605 of file MuonCosmicCompatibilityFiller.cc.

References muon::GlobalMuonPromptTight, muon::isGoodMuon(), query::result, and muon::TMOneStationLoose.

{
  bool result = false;
  // initial set up using Jordan's study: GlobalMuonPromptTight + TMOneStationLoose
  if( muon::isGoodMuon(imuon, muon::GlobalMuonPromptTight) && muon::isGoodMuon(imuon, muon::TMOneStationLoose) ) result = true;
  
  return result;
}

bool MuonCosmicCompatibilityFiller::checkMuonSegments ( const reco::Muon & muon ) const [private]

tag a muon as cosmic based on segment compatibility and the number of segment matches

Definition at line 614 of file MuonCosmicCompatibilityFiller.cc.

References nChamberMatches_, reco::Muon::numberOfMatches(), query::result, segmentComp_, and muon::segmentCompatibility().

                                                                                 {

   bool result = false;
   if (imuon.numberOfMatches() < nChamberMatches_ &&  muon::segmentCompatibility(imuon) < segmentComp_) result = true;

   return result;
}

float MuonCosmicCompatibilityFiller::combinedCosmicID	(	const edm::Event &	iEvent,
		const edm::EventSetup &	iSetup,
		const reco::Muon &	muon,
		bool	CheckMuonID,
		bool	checkVertex
	)		const `[private]`

combined cosmic-likeness: 0 == not cosmic-like

Definition at line 481 of file MuonCosmicCompatibilityFiller.cc.

References backToBack2LegCosmic(), eventActivity(), reco::Muon::isGlobalMuon(), isOverlappingMuon(), muonTiming(), pvMatches(), and query::result.

Referenced by fillCompatibility().

                                                                                                                                {

  float result = 0.0;

  // return >=1 = identify as cosmic muon (the more like cosmics, the higher is the number) 
  // return 0.0 = identify as collision muon
  if( muon.isGlobalMuon() ) {

    unsigned int cosmicVertex = eventActivity(iEvent, muon);
    bool isOverlapping = isOverlappingMuon(iEvent, iSetup, muon);
    unsigned int looseIp = pvMatches(iEvent, muon, true);
    unsigned int tightIp = pvMatches(iEvent, muon, false);
    float looseTime = muonTiming(iEvent, muon, true);
    float tightTime = muonTiming(iEvent, muon, false);
    unsigned int backToback = backToBack2LegCosmic(iEvent,muon);
    //bool cosmicSegment = checkMuonSegments(muon);

    //short cut to reject cosmic event
    if (checkVertex && cosmicVertex == 0) return 10.0;

    // compatibility (0 - 10)
    // weight is assigned by the performance of individual module
    // btob: ~90% eff / ~0% misid
    // ip: ~90% eff / ~0% misid
    // time: ~30% eff / ~0% misid
    double weight_btob = 2.0;
    double weight_ip = 2.0;
    double weight_time = 1.0;
    double weight_overlap = 0.5;

    // collision muon should have compatibility < 4 (0 - 4)
    // cosmic muon should have compatibility >= 4 (4 - 10)

    // b-to-b (max comp.: 4.0)
    if( backToback >= 1 ) {
       //in this case it is cosmic for sure
       result += weight_btob*2.;
       if( tightIp == 1 ) {
           // check with other observables to reduce mis-id (subtract compatibilities)
           if( looseIp == 1 ) {
             if( backToback < 2 ) result -= weight_btob*0.5;
           }
       }
    }

    // ip (max comp.: 4.0)
    if( tightIp == 0 ) {
        //in this case it is cosmic for sure
        result += weight_ip*2.0;
        if( backToback == 0 ) {
           // check with other observables to reduce mis-id (subtract compatibilities)
           if( tightTime == 0 ) {
             if( looseTime == 0 && !isOverlapping ) result -= weight_ip*1.0;
           }
        }
    }
    else if( tightIp >= 2 ) {
        // in this case it is almost collision-like (reduce compatibility)
        // if multi pvs: comp = -2.0 
        if( backToback >= 1 ) result -= weight_ip*1.0;
    }

    // timing (max comp.: 2.0)
    if( tightTime > 0 ) {
        // bonus track
        if( looseTime > 0 ) {
          if( backToback >= 1 ) {
            if( tightIp == 0 ) result += weight_time*tightTime;
            else if( looseIp == 0 ) result += weight_time*0.25;
          }
        }
        else {
          if( backToback >= 1 && tightIp == 0 ) result += weight_time*0.25;
        }
    }

    // overlapping
    if( backToback == 0 && isOverlapping ) {
        // bonus track
        if( tightIp == 0 && tightTime >= 1 ) {
            result += weight_overlap*1.0;
        }
    }
 }//is global muon


//  if (CheckMuonID && cosmicSegment) result += 4;  

  return result;
}

unsigned int MuonCosmicCompatibilityFiller::eventActivity	(	const edm::Event &	iEvent,
		const reco::Muon &	muon
	)		const `[private]`

returns cosmic-likeness based on the event activity information: tracker track multiplicity and vertex quality. 0 == cosmic-like

Definition at line 579 of file MuonCosmicCompatibilityFiller.cc.

References edm::HandleBase::failedToGet(), edm::Event::getByLabel(), inputTrackCollections_, inputVertexCollection_, maxvertRho_, maxvertZ_, minvProb_, edm::Handle< T >::product(), query::result, and testEve_cfg::tracks.

Referenced by combinedCosmicID(), and fillCompatibility().

{

  unsigned int result = 0; //no good vertices - cosmic-like

  //check track activity
  edm::Handle<reco::TrackCollection> tracks;
  iEvent.getByLabel(inputTrackCollections_[0],tracks);
  if (!tracks.failedToGet() && tracks->size() < 3) return 0; 

  //cosmic event should have zero good vertices
  edm::Handle<reco::VertexCollection> pvHandle;
  if (!iEvent.getByLabel(inputVertexCollection_,pvHandle)) {return 0;} else {
      const reco::VertexCollection & vertices = *pvHandle.product();
      //check if vertex collection is empty
      if (vertices.begin() == vertices.end()) return 0;
      for(reco::VertexCollection::const_iterator it=vertices.begin() ; it!=vertices.end() ; ++it){
          if ((TMath::Prob(it->chi2(),(int)it->ndof()) > minvProb_) && (fabs(it->z()) <= maxvertZ_) && (fabs(it->position().rho()) <= maxvertRho_)) result++;
          }
       }
  return result;
}

reco::MuonCosmicCompatibility MuonCosmicCompatibilityFiller::fillCompatibility	(	const reco::Muon &	muon,
		edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)

fill cosmic compatibility variables

Definition at line 111 of file MuonCosmicCompatibilityFiller.cc.

References backToBack2LegCosmic(), reco::MuonCosmicCompatibility::backToBackCompatibility, combinedCosmicID(), reco::MuonCosmicCompatibility::cosmicCompatibility, eventActivity(), reco::MuonCosmicCompatibility::ipCompatibility, isOverlappingMuon(), muonTiming(), reco::MuonCosmicCompatibility::overlapCompatibility, pvMatches(), service_, reco::MuonCosmicCompatibility::timeCompatibility, and reco::MuonCosmicCompatibility::vertexCompatibility.

Referenced by CosmicsMuonIdProducer::produce().

{
  const std::string theCategory = "MuonCosmicCompatibilityFiller";

  reco::MuonCosmicCompatibility returnComp;
  
  service_->update(iSetup);
  
  float timeCompatibility = muonTiming(iEvent, muon, false);
  float backToBackCompatibility = backToBack2LegCosmic(iEvent,muon);    
  float overlapCompatibility = isOverlappingMuon(iEvent,iSetup,muon);
  float ipCompatibility = pvMatches(iEvent,muon,false);
  float vertexCompatibility = eventActivity(iEvent,muon);
  float combinedCompatibility = combinedCosmicID(iEvent,iSetup,muon,false,false);
  
  returnComp.timeCompatibility = timeCompatibility;
  returnComp.backToBackCompatibility = backToBackCompatibility;
  returnComp.overlapCompatibility = overlapCompatibility;
  returnComp.cosmicCompatibility = combinedCompatibility;
  returnComp.ipCompatibility = ipCompatibility;
  returnComp.vertexCompatibility = vertexCompatibility;
  
  return returnComp;
  
}

bool MuonCosmicCompatibilityFiller::isOverlappingMuon	(	const edm::Event &	iEvent,
		const edm::EventSetup &	iSetup,
		const reco::Muon &	muon
	)		const `[private]`

returns cosmic-likeness based on overlap with traversing cosmic muon (only muon/STA hits are used)

Definition at line 267 of file MuonCosmicCompatibilityFiller.cc.

References alignCSCRings::e, edm::HandleBase::failedToGet(), edm::EventSetup::get(), edm::Event::getByLabel(), reco::Muon::innerTrack(), inputCosmicMuonCollection_, inputVertexCollection_, reco::Muon::isGlobalMuon(), edm::Ref< C, T, F >::isNonnull(), mag(), reco::Muon::outerTrack(), edm::Handle< T >::product(), sharedFrac_, sharedHits_, and MuonServiceProxy::trackingGeometry().

Referenced by combinedCosmicID(), and fillCompatibility().

                                                                                                                                {

  // 4 steps in this module
  // step1 : check whether it's 1leg cosmic muon or not
  // step2 : both muons (muons and muonsFromCosmics1Leg) should have close IP
  // step3 : both muons should share very close reference point
  // step4 : check shared hits in both muon tracks

  // check if this muon is available in muonsFromCosmics collection
  bool overlappingMuon = false; //false - not cosmic-like
  if( !muon.isGlobalMuon() ) return false;
  
  // reco muons for cosmics
  edm::Handle<reco::MuonCollection> muonHandle;
  iEvent.getByLabel(inputCosmicMuonCollection_, muonHandle);
  
  // Global Tracking Geometry
  ESHandle<GlobalTrackingGeometry> trackingGeometry;
  iSetup.get<GlobalTrackingGeometryRecord>().get(trackingGeometry);
  
  // PV
  math::XYZPoint RefVtx;
  RefVtx.SetXYZ(0, 0, 0);
  
  edm::Handle<reco::VertexCollection> pvHandle;
  iEvent.getByLabel(inputVertexCollection_,pvHandle);
  const reco::VertexCollection & vertices = *pvHandle.product();
  for(reco::VertexCollection::const_iterator it=vertices.begin() ; it!=vertices.end() ; ++it) {
      RefVtx = it->position();
  }
  
  
  if( !muonHandle.failedToGet() ) {
    for ( reco::MuonCollection::const_iterator cosmicMuon = muonHandle->begin();cosmicMuon !=  muonHandle->end(); ++cosmicMuon ) {
      if ( cosmicMuon->innerTrack() == muon.innerTrack() || cosmicMuon->outerTrack() == muon.outerTrack()) return true;
      
      reco::TrackRef outertrack = muon.outerTrack();
      reco::TrackRef costrack = cosmicMuon->outerTrack();
      
      bool isUp = false;
      if( outertrack->phi() > 0 ) isUp = true; 
      
      // shared hits 
      int RecHitsMuon = outertrack->numberOfValidHits();
      int RecHitsCosmicMuon = 0;
      int shared = 0;
      // count hits for same hemisphere
      if( costrack.isNonnull() ) {
        int nhitsUp = 0;
        int nhitsDown = 0;
        // unused
        //      bool isCosmic1Leg = false;
        //      bool isCloseIP = false;
        //      bool isCloseRef = false;

        for( trackingRecHit_iterator coshit = costrack->recHitsBegin(); coshit != costrack->recHitsEnd(); coshit++ ) {
          if( (*coshit)->isValid() ) {
            DetId id((*coshit)->geographicalId());
            double hity = trackingGeometry->idToDet(id)->position().y();
            if( hity > 0 ) nhitsUp++;
            if( hity < 0 ) nhitsDown++;

            if( isUp && hity > 0 ) RecHitsCosmicMuon++;
            if( !isUp && hity < 0 ) RecHitsCosmicMuon++;
          }
        }
        // step1
        //UNUSED:       if( nhitsUp > 0 && nhitsDown > 0 ) isCosmic1Leg = true;
        //if( !isCosmic1Leg ) continue;
        
        if( outertrack.isNonnull() ) {
          // step2
          //UNUSED:          const double ipErr = (double)outertrack->d0Error();
          //UNUSED:          double ipThreshold  = max(ipThreshold_, ipErr);
          //UNUSED:       if( fabs(outertrack->dxy(RefVtx) + costrack->dxy(RefVtx)) < ipThreshold ) isCloseIP = true;
          //if( !isCloseIP ) continue;

          // step3
          GlobalPoint muonRefVtx( outertrack->vx(), outertrack->vy(), outertrack->vz() );
          GlobalPoint cosmicRefVtx( costrack->vx(), costrack->vy(), costrack->vz() );
          //UNUSED:       float dist = (muonRefVtx - cosmicRefVtx).mag();
          //UNUSED:       if( dist < 0.1 ) isCloseRef = true;
          //if( !isCloseRef ) continue;

          for( trackingRecHit_iterator trkhit = outertrack->recHitsBegin(); trkhit != outertrack->recHitsEnd(); trkhit++ ) {
            if( (*trkhit)->isValid() ) {
              for( trackingRecHit_iterator coshit = costrack->recHitsBegin(); coshit != costrack->recHitsEnd(); coshit++ ) {
                if( (*coshit)->isValid() ) {
                  if( (*trkhit)->geographicalId() == (*coshit)->geographicalId() ) {
                    if( ((*trkhit)->localPosition() - (*coshit)->localPosition()).mag()< 10e-5 ) shared++;
                  }
                  
                }
              }
            }
          }
        }
      }
      // step4
      double fraction = -1;
      if( RecHitsMuon != 0 ) fraction = shared/(double)RecHitsMuon;
    //       std::cout << "shared = " << shared << " " << fraction << " " << RecHitsMuon << " " << RecHitsCosmicMuon << std::endl;
      if( shared > sharedHits_ && fraction > sharedFrac_ ) {
        overlappingMuon = true;
        break;
      }
    }
  }
  
  return overlappingMuon;
}

float MuonCosmicCompatibilityFiller::muonTiming	(	const edm::Event &	iEvent,
		const reco::Muon &	muon,
		bool	isLoose
	)		const `[private]`

check muon time (DT and CSC) information: 0 == prompt-like

Definition at line 141 of file MuonCosmicCompatibilityFiller.cc.

References corrTimeNeg_, corrTimePos_, edm::HandleBase::failedToGet(), edm::Event::getByLabel(), inputMuonCollections_, edm::Ref< C, T, F >::isNonnull(), reco::Muon::isTimeValid(), nMuons(), offTimeNegLoose_, offTimeNegLooseMult_, offTimeNegTight_, offTimeNegTightMult_, offTimePosLoose_, offTimePosLooseMult_, offTimePosTight_, offTimePosTightMult_, reco::Muon::outerTrack(), pvMatches(), query::result, and reco::Muon::time().

Referenced by combinedCosmicID(), and fillCompatibility().

                                                                                                          {

   float offTimeNegMult, offTimePosMult, offTimeNeg, offTimePos;

   if (isLoose) {
     //use "loose" parameter set
     offTimeNegMult = offTimeNegLooseMult_; 
     offTimePosMult = offTimePosLooseMult_;
     offTimeNeg     = offTimeNegLoose_;
     offTimePos     = offTimePosLoose_;
  } else {
     offTimeNegMult = offTimeNegTightMult_;
     offTimePosMult = offTimePosTightMult_;
     offTimeNeg     = offTimeNegTight_;
     offTimePos     = offTimePosTight_;
  }

  float result = 0.0;

  if( muon.isTimeValid() ) {
    //case of multiple muon event
    if (nMuons(iEvent) > 1) {

      float positiveTime = 0;
      if ( muon.time().timeAtIpInOut < offTimeNegMult || muon.time().timeAtIpInOut > offTimePosMult) result = 1.;
      if ( muon.time().timeAtIpInOut > 0.) positiveTime = muon.time().timeAtIpInOut;

      //special case, looking for time-correlation 
      // between muons in opposite hemispheres 
      if (!isLoose && result == 0 && positiveTime > corrTimePos_) { 

          //check hemi of this muon
          bool isUp = false;
          reco::TrackRef outertrack = muon.outerTrack();
          if( outertrack.isNonnull() ) {
              if( outertrack->phi() > 0 ) isUp = true;

              //loop over muons in that event and find if there are any in the opposite hemi 
              edm::Handle<reco::MuonCollection> muonHandle;
              iEvent.getByLabel(inputMuonCollections_[1], muonHandle);

              if( !muonHandle.failedToGet() ) {
                  for ( reco::MuonCollection::const_iterator iMuon = muonHandle->begin(); iMuon !=  muonHandle->end(); ++iMuon ) {
                      if (!iMuon->isGlobalMuon()) continue;

                      reco::TrackRef checkedTrack = iMuon->outerTrack();
                      if( muon.isTimeValid() ) {

                          // from bottom up
                          if (checkedTrack->phi() < 0 && isUp) {
                              if (iMuon->time().timeAtIpInOut < corrTimeNeg_) result = 1.0;
                              break;
                          } else if (checkedTrack->phi() > 0 && !isUp) {
                               // from top down 
                               if (iMuon->time().timeAtIpInOut < corrTimeNeg_) result = 1.0; 
                               break;
                             }
                         } //muon is time valid
                     } 
                 }
             } //track is nonnull
         } //double check timing
     } else {
        //case of a single muon event
        if ( muon.time().timeAtIpInOut < offTimeNeg || muon.time().timeAtIpInOut > offTimePos) result = 1.;
      }
    } //is time valid

  
   if (!isLoose && result > 0) {
     //check loose ip
     if (pvMatches(iEvent, muon, true) == 0) result *= 2.;
   } 

   return result;
}

unsigned int MuonCosmicCompatibilityFiller::nMuons ( const edm::Event & iEvent ) const [private]

get number of muons in the vent

Definition at line 245 of file MuonCosmicCompatibilityFiller.cc.

References edm::HandleBase::failedToGet(), edm::Event::getByLabel(), and inputMuonCollections_.

Referenced by muonTiming(), and pvMatches().

                                                                  {

    unsigned int nGlb = 0;

    edm::Handle<reco::MuonCollection> muonHandle;
    iEvent.getByLabel(inputMuonCollections_[1], muonHandle);

    if( !muonHandle.failedToGet() ) {
      for ( reco::MuonCollection::const_iterator iMuon = muonHandle->begin(); iMuon !=  muonHandle->end(); ++iMuon ) {
         if (!iMuon->isGlobalMuon()) continue;
           nGlb++;
        }
      } 
 
    return nGlb;
}

unsigned int MuonCosmicCompatibilityFiller::pvMatches	(	const edm::Event &	iEvent,
		const reco::Muon &	muon,
		bool	isLoose
	)		const `[private]`

return cosmic-likeness based on the 2D impact parameters (dxy, dz wrt to PV). 0 == cosmic-like

Definition at line 383 of file MuonCosmicCompatibilityFiller.cc.

References edm::HandleBase::failedToGet(), edm::Event::getByLabel(), hIpTrdxy_, hIpTrvProb_, reco::Muon::innerTrack(), inputMuonCollections_, inputVertexCollection_, reco::Muon::isGlobalMuon(), edm::Ref< C, T, F >::isNonnull(), reco::Muon::isStandAloneMuon(), reco::Muon::isTrackerMuon(), largedxy_, largedxyMult_, maxdxyLoose_, maxdxyLooseMult_, maxdxyTight_, maxdxyTightMult_, maxdzLoose_, maxdzLooseMult_, maxdzTight_, maxdzTightMult_, patZpeak::muons, nMuons(), reco::Muon::outerTrack(), edm::Handle< T >::product(), query::result, reco::Muon::standAloneMuon(), reco::Muon::track(), and testEve_cfg::tracks.

Referenced by combinedCosmicID(), fillCompatibility(), and muonTiming().

                                                                                                         {

   float maxdxyMult, maxdzMult, maxdxy, maxdz;

   if (isLoose) {
     //use "loose" parameter set
     maxdxyMult = maxdxyLooseMult_;
     maxdzMult  = maxdzLooseMult_;
     maxdxy     = maxdxyLoose_;
     maxdz      = maxdzLoose_;
  } else {
     maxdxyMult = maxdxyTightMult_;
     maxdzMult  = maxdzTightMult_;
     maxdxy     = maxdxyTight_;
     maxdz      = maxdzTight_;
  }

  unsigned int result = 0;

  reco::TrackRef track;
  if ( muon.isGlobalMuon() )          track = muon.innerTrack();
  else if ( muon.isTrackerMuon() )    track = muon.track();
  else if ( muon.isStandAloneMuon())  track = muon.standAloneMuon();
  
  bool multipleMu = false;
  if (nMuons(iEvent) > 1) multipleMu = true;

  math::XYZPoint RefVtx;
  RefVtx.SetXYZ(0, 0, 0);

  edm::Handle<reco::VertexCollection> pvHandle;
  iEvent.getByLabel(inputVertexCollection_,pvHandle);
  const reco::VertexCollection & vertices = *pvHandle.product();
  for(reco::VertexCollection::const_iterator it=vertices.begin() ; it!=vertices.end() ; ++it){
    RefVtx = it->position();

    if ( track.isNonnull() ) {
       if (multipleMu) {
         //multiple muon event

         if ( fabs( (*track).dxy(RefVtx) ) < maxdxyMult || fabs( (*track).dz(RefVtx) ) < maxdzMult) {
           result++;

           //case of extra large dxy
           if (!isLoose && fabs( (*track).dxy(RefVtx) ) > largedxyMult_) result -= 1;  

          }
     } else {
        //single muon event

        if (fabs( (*track).dxy(RefVtx) ) < maxdxy || fabs( (*track).dz(RefVtx) ) < maxdz) {
          result++;

          //case of extra large dxy
          if (!isLoose && fabs( (*track).dxy(RefVtx) ) > largedxy_) result -= 1; 

         }
      }
   }//track is nonnull
}//loop over vertices

   //special case for non-cosmic large ip muons
   if (result == 0 && multipleMu) {
      // consider all reco muons in an event
      edm::Handle<reco::MuonCollection> muonHandle;
      iEvent.getByLabel(inputMuonCollections_[1], muonHandle);

      //cosmic event should have zero good vertices
      edm::Handle<reco::VertexCollection> pvHandle;
      iEvent.getByLabel(inputVertexCollection_,pvHandle);
      const reco::VertexCollection & vertices = *pvHandle.product();

      //find the "other" one
      if( !muonHandle.failedToGet() ) {
         for ( reco::MuonCollection::const_iterator muons = muonHandle->begin(); muons !=  muonHandle->end(); ++muons ) {
            if (!muons->isGlobalMuon()) continue;
            //skip this track  
            if ( muons->innerTrack() == muon.innerTrack() && muons->outerTrack() == muon.outerTrack()) continue;
                //check ip and vertex of the "other" muon
                reco::TrackRef tracks;
                if (muons->isGlobalMuon()) tracks = muons->innerTrack();
                if (fabs((*tracks).dxy(RefVtx)) > hIpTrdxy_) continue; 
                //check if vertex collection is empty
                if (vertices.begin() == vertices.end()) continue;
                //for(reco::VertexCollection::const_iterator it=vertices.begin() ; it!=vertices.end() ; ++it) {
                    //find matching vertex by position
                    //if (fabs(it->z() - tracks->vz()) > 0.01) continue; //means will not be untagged from cosmics 
                    if (TMath::Prob(vertices.front().chi2(),(int)(vertices.front().ndof())) > hIpTrvProb_) result = 1; 
               //}
           }
      }
 }

    return result;

}