CSCSkim Class Reference

#include <RecoLocalMuon/CSCSkim/src/CSCSkim.cc>

Inheritance diagram for CSCSkim:

Public Member Functions
virtual void	beginJob ()
	CSCSkim (const edm::ParameterSet &pset)
virtual void	endJob ()
virtual bool	filter (edm::Event &event, const edm::EventSetup &eventSetup)
	~CSCSkim ()
Public Member Functions inherited from edm::EDFilter
	EDFilter ()
virtual	~EDFilter ()
Public Member Functions inherited from edm::ProducerBase
	ProducerBase ()
void	registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
boost::function< void(const BranchDescription &)>	registrationCallback () const
	used by the fwk to register list of products More...
virtual	~ProducerBase ()

Private Member Functions
int	chamberSerial (int kE, int kS, int kR, int kCh)
bool	doBFieldStudySelection (edm::Handle< reco::TrackCollection > saTracks, edm::Handle< reco::TrackCollection > Tracks, edm::Handle< reco::MuonCollection > gMuons)
bool	doCertainChamberSelection (edm::Handle< CSCWireDigiCollection > wires, edm::Handle< CSCStripDigiCollection > strips)
bool	doCSCSkimming (edm::Handle< CSCRecHit2DCollection > cscRecHits, edm::Handle< CSCSegmentCollection > cscSegments)
bool	doDTOverlap (edm::Handle< CSCSegmentCollection > cscSegments)
bool	doHaloLike (edm::Handle< CSCSegmentCollection > cscSegments)
bool	doLongSATrack (edm::Handle< reco::TrackCollection > saTracks)
bool	doMessyEventSkimming (edm::Handle< CSCRecHit2DCollection > cscRecHits, edm::Handle< CSCSegmentCollection > cscSegments)
bool	doOverlapSkimming (edm::Handle< CSCSegmentCollection > cscSegments)

Private Attributes
edm::InputTag	cscRecHitTag
edm::InputTag	cscSegmentTag
bool	demandChambersBothSides
edm::InputTag	GLBMuonTag
std::string	histogramFileName
TH1F *	hxnHitChambers
TH1F *	hxnRecHits
TH1F *	hxnRecHitsSel
TH1F *	hxnSegments
int	iEvent
int	iRun
bool	isSimulation
bool	makeHistograms
bool	makeHistogramsForMessyEvents
TH1F *	mevnChambers0
TH1F *	mevnChambers1
TH1F *	mevnRecHits0
TH1F *	mevnRecHits1
TH1F *	mevnSegments0
TH1F *	mevnSegments1
int	minimumHitChambers
int	minimumSegments
int	nCSCHitsMin
int	nEventsAnalyzed
int	nEventsCertainChamber
int	nEventsChambersBothSides
int	nEventsDTOverlap
int	nEventsForBFieldStudies
int	nEventsHaloLike
int	nEventsLongSATrack
int	nEventsMessy
int	nEventsOverlappingChambers
int	nEventsSelected
int	nLayersWithHitsMinimum
int	nTrHitsMin
int	nValidHitsMin
std::string	outputFileName
float	pMin
float	redChiSqMax
float	rExtMax
edm::InputTag	SAMuonTag
TFile *	theHistogramFile
edm::InputTag	trackTag
int	typeOfSkim
int	whichChamber
int	whichEndcap
int	whichRing
int	whichStation
TH1F *	xxnCSCHits
TH1F *	xxnTrackerHits
TH1F *	xxnValidHits
TH1F *	xxP
TH1F *	xxredChiSq
float	zInnerMax
float	zLengthMin
float	zLengthTrMin

Additional Inherited Members
Public Types inherited from edm::EDFilter
typedef EDFilter	ModuleType
typedef WorkerT< EDFilter >	WorkerType
Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList	TypeLabelList
Static Public Member Functions inherited from edm::EDFilter
static const std::string &	baseType ()
static void	fillDescriptions (ConfigurationDescriptions &descriptions)
Protected Member Functions inherited from edm::EDFilter
CurrentProcessingContext const *	currentContext () const
Protected Member Functions inherited from edm::ProducerBase
template<class TProducer , class TMethod >
void	callWhenNewProductsRegistered (TProducer *iProd, TMethod iMethod)

Detailed Description

This simple program selects minimal CSC events for output.

Michael Schmitt, Northwestern University, July 2008

Description: Offline skim module for CSC cosmic ray data

Implementation: <Notes on="" implementation>="">

Definition at line 74 of file CSCSkim.h.

Constructor & Destructor Documentation

CSCSkim::CSCSkim ( const edm::ParameterSet &  pset )

explicit

Definition at line 53 of file CSCSkim.cc.

References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), RPCNoise_example::histogramFileName, and dumpDBToFile_GT_ttrig_cfg::outputFileName.

{

  // input tags
  cscRecHitTag  = pset.getParameter<edm::InputTag>("cscRecHitTag");
  cscSegmentTag = pset.getParameter<edm::InputTag>("cscSegmentTag");
  SAMuonTag     = pset.getParameter<edm::InputTag>("SAMuonTag");
  GLBMuonTag    = pset.getParameter<edm::InputTag>("GLBMuonTag");
  trackTag      = pset.getParameter<edm::InputTag>("trackTag");

  // Get the various input parameters
  bool isSimulation = false;
  isSimulation       = pset.getUntrackedParameter<bool>("isSimulation",false);
  outputFileName     = pset.getUntrackedParameter<std::string>("outputFileName","outputSkim.root");
  histogramFileName  = pset.getUntrackedParameter<std::string>("histogramFileName","histos.root");
  typeOfSkim              = pset.getUntrackedParameter<int>("typeOfSkim",1);
  nLayersWithHitsMinimum  = pset.getUntrackedParameter<int>("nLayersWithHitsMinimum",3);
  minimumHitChambers      = pset.getUntrackedParameter<int>("minimumHitChambers",1);
  minimumSegments         = pset.getUntrackedParameter<int>("minimumSegments",3);
  demandChambersBothSides = pset.getUntrackedParameter<bool>("demandChambersBothSides",false);
  makeHistograms          = pset.getUntrackedParameter<bool>("makeHistograms",false);
  makeHistogramsForMessyEvents = pset.getUntrackedParameter<bool>("makeHistogramsForMessyEvebts",false);
  whichEndcap             = pset.getUntrackedParameter<int>("whichEndcap",2);
  whichStation            = pset.getUntrackedParameter<int>("whichStation",3);
  whichRing               = pset.getUntrackedParameter<int>("whichRing",2);
  whichChamber            = pset.getUntrackedParameter<int>("whichChamber",24);

  // for BStudy selection (skim type 9)
  pMin               = pset.getUntrackedParameter<double>("pMin",3.);
  zLengthMin         = pset.getUntrackedParameter<double>("zLengthMin",200.);
  nCSCHitsMin        = pset.getUntrackedParameter<int>("nCSCHitsMin",9);
  zInnerMax          = pset.getUntrackedParameter<double>("zInnerMax",9000.);
  nTrHitsMin         = pset.getUntrackedParameter<int>("nTrHitsMin",8);
  zLengthTrMin       = pset.getUntrackedParameter<double>("zLengthTrMin",180.);
  rExtMax            = pset.getUntrackedParameter<double>("rExtMax",3000.);
  redChiSqMax        = pset.getUntrackedParameter<double>("redChiSqMax",20.);
  nValidHitsMin      = pset.getUntrackedParameter<int>("nValidHitsMin",8);

  LogInfo("[CSCSkim] Setup")
    << "\n\t===== CSCSkim =====\n"
    << "\t\ttype of skim ...............................\t" << typeOfSkim
    << "\t\tminimum number of layers with hits .........\t" << nLayersWithHitsMinimum
    << "\n\t\tminimum number of chambers w/ hit layers..\t" << minimumHitChambers
    << "\n\t\tminimum number of segments ...............\t" << minimumSegments
    << "\n\t\tdemand chambers on both sides.............\t" << demandChambersBothSides
    << "\n\t\tmake histograms...........................\t" << makeHistograms
    << "\n\t\t..for messy events........................\t" << makeHistogramsForMessyEvents
    << "\n\t===================\n\n";

}

CSCSkim::~CSCSkim

(

)

Definition at line 107 of file CSCSkim.cc.

{
}

Member Function Documentation

void CSCSkim::beginJob ( void  )

virtual

Reimplemented from edm::EDFilter.

Definition at line 116 of file CSCSkim.cc.

References RPCNoise_example::histogramFileName, and iEvent.

{
  // set counters to zero
  nEventsAnalyzed = 0;
  nEventsSelected = 0;
  nEventsChambersBothSides = 0;
  nEventsOverlappingChambers = 0;
  nEventsMessy = 0;
  nEventsCertainChamber = 0;
  nEventsDTOverlap = 0;
  nEventsHaloLike = 0;
  nEventsLongSATrack = 0;
  nEventsForBFieldStudies = 0;
  iRun = 0;
  iEvent = 0;

  if (makeHistograms || makeHistogramsForMessyEvents) {

    // Create the root file for the histograms
    theHistogramFile = new TFile(histogramFileName.c_str(), "RECREATE");
    theHistogramFile->cd();

    if (makeHistograms) {
      // book histograms for the skimming module
      hxnRecHits     = new TH1F("hxnRecHits","n RecHits",61,-0.5,60.5);
      hxnSegments    = new TH1F("hxnSegments","n Segments",11,-0.5,10.5);
      hxnHitChambers = new TH1F("hxnHitsChambers","n chambers with hits",11,-0.5,10.5);
      hxnRecHitsSel  = new TH1F("hxnRecHitsSel","n RecHits selected",61,-0.5,60.5);

      xxP = new TH1F("xxP","P global",100,0.,200.);
      xxnValidHits = new TH1F("xxnValidHits","n valid hits global",61,-0.5,60.5);
      xxnTrackerHits = new TH1F("xxnTrackerHits","n tracker hits global",61,-0.5,60.5);
      xxnCSCHits = new TH1F("xxnCSCHits","n CSC hits global",41,-0.5,40.5);
      xxredChiSq = new TH1F("xxredChiSq","red chisq global",100,0.,100.);



    }
    if (makeHistogramsForMessyEvents) {
      // book histograms for the messy event skimming module
      mevnRecHits0     = new TH1F("mevnRecHits0","n RecHits",121,-0.5,120.5);
      mevnChambers0    = new TH1F("mevnChambers0","n chambers with hits",21,-0.5,20.5);
      mevnSegments0    = new TH1F("mevnSegments0","n Segments",21,-0.5,20.5);
      mevnRecHits1     = new TH1F("mevnRecHits1","n RecHits",100,0.,300.);
      mevnChambers1    = new TH1F("mevnChambers1","n chambers with hits",50,0.,50.);
      mevnSegments1    = new TH1F("mevnSegments1","n Segments",30,0.,30.);
    }

  }
}

int CSCSkim::chamberSerial ( int  kE, int  kS, int  kR, int  kCh  )

private

Definition at line 1276 of file CSCSkim.cc.

                                                                               {
    int kSerial = kChamber;
    if (kStation == 1 && kRing == 1) {kSerial = kChamber;}
    if (kStation == 1 && kRing == 2) {kSerial = kChamber + 36;}
    if (kStation == 1 && kRing == 3) {kSerial = kChamber + 72;}
    if (kStation == 1 && kRing == 4) {kSerial = kChamber;}
    if (kStation == 2 && kRing == 1) {kSerial = kChamber + 108;}
    if (kStation == 2 && kRing == 2) {kSerial = kChamber + 126;}
    if (kStation == 3 && kRing == 1) {kSerial = kChamber + 162;}
    if (kStation == 3 && kRing == 2) {kSerial = kChamber + 180;}
    if (kStation == 4 && kRing == 1) {kSerial = kChamber + 216;}
    if (kStation == 4 && kRing == 2) {kSerial = kChamber + 234;}  // one day...
    if (kEndcap == 2) {kSerial = kSerial + 300;}
    return kSerial;
}

bool CSCSkim::doBFieldStudySelection ( edm::Handle< reco::TrackCollection >  saTracks, edm::Handle< reco::TrackCollection >  Tracks, edm::Handle< reco::MuonCollection >  gMuons  )

private

Definition at line 1070 of file CSCSkim.cc.

References abs, MuonSubdetId::CSC, MuonSubdetId::DT, spr::goodTrack(), DetId::Muon, metsig::muon, n, reco::HitPattern::numberOfValidTrackerHits(), mathSSE::sqrt(), and PV3DBase< T, PVType, FrameType >::z().

                                                                                                                                                         {

  bool acceptThisEvent = false;

  //-----------------------------------
  // examine the stand-alone tracks
  //-----------------------------------
  int nGoodSAMuons = 0;
  for (reco::TrackCollection::const_iterator muon = saMuons->begin(); muon != saMuons->end(); ++ muon ) {
    float preco  = muon->p();

    math::XYZPoint innerPo = muon->innerPosition();
    GlobalPoint iPnt(innerPo.x(), innerPo.y(),innerPo.z());
    math::XYZPoint outerPo = muon->outerPosition();
    GlobalPoint oPnt(outerPo.x(), outerPo.y(),outerPo.z());
    float zLength = abs( iPnt.z() - oPnt.z() );

    math::XYZVector innerMom = muon->innerMomentum();
    GlobalVector iP(innerMom.x(), innerMom.y(), innerMom.z() );
    math::XYZVector outerMom = muon->outerMomentum();
    GlobalVector oP(outerMom.x(), outerMom.y(), outerMom.z() );

    const float zRef = 300.;
    float xExt = 10000.;
    float yExt = 10000.;
    if (abs(oPnt.z()) < abs(iPnt.z())) {
      float deltaZ = 0.;
      if (oPnt.z() > 0) {
    deltaZ = zRef - oPnt.z();
      } else {
    deltaZ = -zRef - oPnt.z();
      }
      xExt = oPnt.x() + deltaZ * oP.x() / oP.z();
      yExt = oPnt.y() + deltaZ * oP.y() / oP.z();
    } else {
      float deltaZ = 0.;
      if (iPnt.z() > 0) {
    deltaZ = zRef - iPnt.z();
      } else {
    deltaZ = -zRef - iPnt.z();
      }
      xExt = iPnt.x() + deltaZ * iP.x() / iP.z();
      yExt = iPnt.y() + deltaZ * iP.y() / iP.z();
    }
    float rExt = sqrt( xExt*xExt + yExt*yExt );
    
    int kHit = 0;
    int nDTHits = 0;
    int nCSCHits = 0;
    for (trackingRecHit_iterator hit = muon->recHitsBegin(); hit != muon->recHitsEnd(); ++hit ) {
      ++kHit;
      const DetId detId( (*hit)->geographicalId() );
      if (detId.det() == DetId::Muon) {
    if (detId.subdetId() == MuonSubdetId::DT) {
      nDTHits++;
    }
    else if (detId.subdetId() == MuonSubdetId::CSC) {
      nCSCHits++;
    }
      }
    } // end loop over hits

    float zInner = -1.;
    float zOuter = 20000.;
    if (nCSCHits >= nCSCHitsMin) {
      if (abs(iPnt.z()) < abs(iPnt.z())) {
    zInner = iPnt.z();
    zOuter = oPnt.z();
      } else {
    zInner = oPnt.z();
    zOuter = iPnt.z();
      }
    }

    bool goodSAMuon = (preco > pMin)
      && ( zLength > zLengthMin )
      && ( nCSCHits >= nCSCHitsMin ) 
      && ( zInner < zInnerMax )
      && ( rExt < rExtMax ) ;

    if (goodSAMuon) {nGoodSAMuons++;}
    
  } // end loop over stand-alone muon collection




  //-----------------------------------
  // examine the tracker tracks
  //-----------------------------------
  int nGoodTracks = 0;
  for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++ track ) {
    float preco  = track->p();
    int   n = track->recHitsSize();

    math::XYZPoint innerPo = track->innerPosition();
    GlobalPoint iPnt(innerPo.x(), innerPo.y(),innerPo.z());
    math::XYZPoint outerPo = track->outerPosition();
    GlobalPoint oPnt(outerPo.x(), outerPo.y(),outerPo.z());
    float zLength = abs( iPnt.z() - oPnt.z() );

    math::XYZVector innerMom = track->innerMomentum();
    GlobalVector iP(innerMom.x(), innerMom.y(), innerMom.z() );
    math::XYZVector outerMom = track->outerMomentum();
    GlobalVector oP(outerMom.x(), outerMom.y(), outerMom.z() );

    const float zRef = 300.;
    float xExt = 10000.;
    float yExt = 10000.;
    if (abs(oPnt.z()) > abs(iPnt.z())) {
      float deltaZ = 0.;
      if (oPnt.z() > 0) {
    deltaZ = zRef - oPnt.z();
      } else {
    deltaZ = -zRef - oPnt.z();
      }
      xExt = oPnt.x() + deltaZ * oP.x() / oP.z();
      yExt = oPnt.y() + deltaZ * oP.y() / oP.z();
    } else {
      float deltaZ = 0.;
      if (iPnt.z() > 0) {
    deltaZ = zRef - iPnt.z();
      } else {
    deltaZ = -zRef - iPnt.z();
      }
      xExt = iPnt.x() + deltaZ * iP.x() / iP.z();
      yExt = iPnt.y() + deltaZ * iP.y() / iP.z();
    }
    float rExt = sqrt( xExt*xExt + yExt*yExt );

    float zInner = -1.;
    float zOuter = 20000.;
    if (n >= nTrHitsMin) {
      if (abs(iPnt.z()) < abs(iPnt.z())) {
    zInner = iPnt.z();
    zOuter = oPnt.z();
      } else {
    zInner = oPnt.z();
    zOuter = iPnt.z();
      }
    }

    bool goodTrack = (preco > pMin)
      && (n >= nTrHitsMin)
      && (zLength > zLengthTrMin)
      && ( rExt < rExtMax ) ;

    if (goodTrack) {nGoodTracks++;}

  } // end loop over tracker tracks


  //-----------------------------------
  // examine the global muons
  //-----------------------------------
  int nGoodGlobalMuons = 0;
  for (reco::MuonCollection::const_iterator global = gMuons->begin(); global != gMuons->end(); ++global ) {

    if (global->isGlobalMuon()) {

      float pDef  = global->p();
      float redChiSq = global->globalTrack()->normalizedChi2();
      const reco::HitPattern& hp = (global->globalTrack())->hitPattern();
      // int nTotalHits = hp.numberOfHits();
      //    int nValidHits = hp.numberOfValidHits();
      int nTrackerHits = hp.numberOfValidTrackerHits();
      // int nPixelHits   = hp.numberOfValidPixelHits();
      // int nStripHits   = hp.numberOfValidStripHits();
      
      int nDTHits = 0;
      int nCSCHits = 0;
      for (trackingRecHit_iterator hit = (global->globalTrack())->recHitsBegin(); hit != (global->globalTrack())->recHitsEnd(); ++hit ) {
    const DetId detId( (*hit)->geographicalId() );
    if (detId.det() == DetId::Muon) {
      if (detId.subdetId() == MuonSubdetId::DT) {
        nDTHits++;
      }
      else if (detId.subdetId() == MuonSubdetId::CSC) {
        nCSCHits++;
      }
    }
      } // end loop over hits
      
      bool goodGlobalMuon = (pDef > pMin)
    && ( nTrackerHits >= nValidHitsMin )
    && ( nCSCHits >= nCSCHitsMin ) 
    && ( redChiSq < redChiSqMax );

      if (goodGlobalMuon) {nGoodGlobalMuons++;}
      
    } // this is a global muon
  } // end loop over stand-alone muon collection

  //-----------------------------------
  // do we accept this event?
  //-----------------------------------

  acceptThisEvent = ( (nGoodSAMuons > 0) && (nGoodTracks > 0) ) || (nGoodGlobalMuons > 0) ;

  return acceptThisEvent;
}

bool CSCSkim::doCertainChamberSelection ( edm::Handle< CSCWireDigiCollection >  wires, edm::Handle< CSCStripDigiCollection >  strips  )

private

Definition at line 695 of file CSCSkim.cc.

References CSCDetId, and CSCDetId::endcap().

                                                              {

  // Loop through the wire DIGIs, looking for a match
  bool certainChamberIsPresentInWires = false;
  for (CSCWireDigiCollection::DigiRangeIterator jw=wires->begin(); jw!=wires->end(); jw++) {
    CSCDetId id = (CSCDetId)(*jw).first;
    int kEndcap  = id.endcap();
    int kRing    = id.ring();
    int kStation = id.station();
    int kChamber = id.chamber();
    if ( (kEndcap     == whichEndcap) &&
         (kStation    == whichStation) &&
         (kRing       == whichRing) &&
         (kChamber    == whichChamber) )
      {certainChamberIsPresentInWires = true;}
  } // end wire loop


  // Loop through the strip DIGIs, looking for a match
  bool certainChamberIsPresentInStrips = false;
  for (CSCStripDigiCollection::DigiRangeIterator js=strips->begin(); js!=strips->end(); js++) {
    CSCDetId id = (CSCDetId)(*js).first;
    int kEndcap  = id.endcap();
    int kRing    = id.ring();
    int kStation = id.station();
    int kChamber = id.chamber();
    if ( (kEndcap     == whichEndcap) &&
         (kStation    == whichStation) &&
         (kRing       == whichRing) &&
         (kChamber    == whichChamber) )
      {certainChamberIsPresentInStrips = true;}
  }

  bool certainChamberIsPresent = certainChamberIsPresentInWires || certainChamberIsPresentInStrips;

  return certainChamberIsPresent;
}

bool CSCSkim::doCSCSkimming ( edm::Handle< CSCRecHit2DCollection >  cscRecHits, edm::Handle< CSCSegmentCollection >  cscSegments  )

private

Definition at line 354 of file CSCSkim.cc.

References CSCDetId::chamber(), CSCDetId, CSCDetId::endcap(), i, j, kLayer(), CSCDetId::layer(), LogDebug, funct::pow(), CSCDetId::ring(), and CSCDetId::station().

                                                                                                                   {

  // how many RecHits in the collection?
  int nRecHits = cscRecHits->size();

  // zero the recHit counter
  int cntRecHit[600];
  for (int i = 0; i < 600; i++) {
    cntRecHit[i] = 0;
  }

  // ---------------------
  // Loop over rechits 
  // ---------------------

  CSCRecHit2DCollection::const_iterator recIt;
  for (recIt = cscRecHits->begin(); recIt != cscRecHits->end(); recIt++) {

    // which chamber is it?
    CSCDetId idrec = (CSCDetId)(*recIt).cscDetId();
    int kEndcap  = idrec.endcap();
    int kRing    = idrec.ring();
    int kStation = idrec.station();
    int kChamber = idrec.chamber();
    int kLayer   = idrec.layer();

    // compute chamber serial number
    int kSerial = chamberSerial( kEndcap, kStation, kRing, kChamber ) ;

    // increment recHit counter
    //     (each layer is represented by a different power of 10)
    int kDigit = (int) pow((float)10.,(float)(kLayer-1));
    cntRecHit[kSerial] += kDigit;

  } //end rechit loop


  // ------------------------------------------------------
  // Are there chambers with the minimum number of hits?
  // ------------------------------------------------------

  int nChambersWithMinimalHits = 0;
  int nChambersWithMinimalHitsPOS = 0;
  int nChambersWithMinimalHitsNEG = 0;
  if (nRecHits > 0) {
    for (int i = 0; i < 600; i++) {
      if (cntRecHit[i] > 0) {
    int nLayersWithHits = 0;
    float dummy = (float) cntRecHit[i];
    for (int j = 5; j > -1; j--) {
      float digit = dummy / pow( (float)10., (float)j );
      int kCount = (int) digit;
      if (kCount > 0) nLayersWithHits++;
      dummy = dummy - ( (float) kCount) * pow( (float)10., (float)j );
    }
    if (nLayersWithHits > nLayersWithHitsMinimum) {
      if (i < 300) {nChambersWithMinimalHitsPOS++;}
      else         {nChambersWithMinimalHitsNEG++;}
    }
      }
    }
    nChambersWithMinimalHits = nChambersWithMinimalHitsPOS + nChambersWithMinimalHitsNEG;
  }

  // how many Segments?
  int nSegments = cscSegments->size();

  // ----------------------
  // fill histograms
  // ----------------------

  if (makeHistograms) {
    hxnRecHits->Fill(nRecHits);
    if (nRecHits > 0) {
      hxnSegments->Fill(nSegments);
      hxnHitChambers->Fill(nChambersWithMinimalHits);
    }
    if (nChambersWithMinimalHits > 0) {
      hxnRecHitsSel->Fill(nRecHits);
    }
  }

  // ----------------------
  // set the filter flag
  // ----------------------
  bool basicEvent = ( nChambersWithMinimalHits >= minimumHitChambers ) && ( nSegments >= minimumSegments );

  bool chambersOnBothSides = ((nChambersWithMinimalHitsPOS >= minimumHitChambers) && (nChambersWithMinimalHitsNEG >= minimumHitChambers));

  if (chambersOnBothSides) {nEventsChambersBothSides++;}

  bool selectEvent = false;
  if (typeOfSkim == 1) {selectEvent = basicEvent;}
  if (typeOfSkim == 2) {selectEvent = chambersOnBothSides;}


  // debug
  LogDebug("[CSCSkim]") << "----- nRecHits = " << nRecHits
               << "\tnChambersWithMinimalHits = " << nChambersWithMinimalHits
               << "\tnSegments = " << nSegments 
               << "\tselect? " << selectEvent << std::endl;

  /*
  if ((nChambersWithMinimalHitsPOS >= minimumHitChambers) && (nChambersWithMinimalHitsNEG >= minimumHitChambers)) {
    std::cout << "\n==========================================================================\n"
     << "\tinteresting event - chambers hit on both sides\n"
     << "\t  " <<  nEventsAnalyzed
     << "\trun " << iRun << "\tevent " << iEvent << std::endl;
    std::cout << "----- nRecHits = " << nRecHits
     << "\tnChambersWithMinimalHits = " << nChambersWithMinimalHits
     << "\tnSegments = " << nSegments 
     << "\tselect? " << selectEvent << std::endl;
    for (int i = 0; i < 600; i++) {
      if (cntRecHit[i] > 0) {
    cout << "\t\t" << i << "\tcntRecHit= " << cntRecHit[i] << std::endl;
      }
    }
    std::cout << "==========================================================================\n\n" ;
  }
  */

  return selectEvent;
}

bool CSCSkim::doDTOverlap ( edm::Handle< CSCSegmentCollection >  cscSegments )

private

Definition at line 741 of file CSCSkim.cc.

References CSCDetId, CSCDetId::endcap(), and i.

                                                                  {
  const float chisqMax = 100.;
  const int nhitsMin = 5;
  const int maxNSegments = 3;

  // initialize
  bool DTOverlapCandidate = false;
  int cntMEP13[36];
  int cntMEN13[36];
  int cntMEP22[36];
  int cntMEN22[36];
  int cntMEP32[36];
  int cntMEN32[36];
  for (int i=0; i<36; ++i) {
    cntMEP13[i] = 0;
    cntMEN13[i] = 0;
    cntMEP22[i] = 0;
    cntMEN22[i] = 0;
    cntMEP32[i] = 0;
    cntMEN32[i] = 0;
  }

  // -----------------------
  // loop over segments
  // -----------------------

  int nSegments = cscSegments->size();
  if (nSegments < 2) return DTOverlapCandidate;

  for(CSCSegmentCollection::const_iterator it=cscSegments->begin(); it != cscSegments->end(); it++) {
    // which chamber?
    CSCDetId id  = (CSCDetId)(*it).cscDetId();
    int kEndcap  = id.endcap();
    int kStation = id.station();
    int kRing    = id.ring();
    int kChamber = id.chamber();
    // segment information
    float chisq    = (*it).chi2();
    int nhits      = (*it).nRecHits();
    bool goodSegment = (chisq < chisqMax) && (nhits >= nhitsMin) ;
    if (goodSegment) {
      if ( (kStation == 1) && (kRing == 3) ) {
    if (kEndcap == 1) {cntMEP13[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN13[kChamber-1]++;}
      }
      if ( (kStation == 2) && (kRing == 2) ) {
    if (kEndcap == 1) {cntMEP22[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN22[kChamber-1]++;}
      }
      if ( (kStation == 3) && (kRing == 2) ) {
    if (kEndcap == 1) {cntMEP32[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN32[kChamber-1]++;}
      }
    } // this is a good segment
  } // end loop over segments

  // ---------------------------------------------
  // veto messy events
  // ---------------------------------------------
  bool tooManySegments = false;
  for (int i=0; i<36; ++i) {
    if ( (cntMEP13[i] > maxNSegments) ||
         (cntMEN13[i] > maxNSegments) ||
         (cntMEP22[i] > maxNSegments) ||
         (cntMEN22[i] > maxNSegments) ||
         (cntMEP32[i] > maxNSegments) ||
         (cntMEN32[i] > maxNSegments) ) tooManySegments = true;
  }
  if (tooManySegments) {
    return DTOverlapCandidate;
  }

  // ---------------------------------------------
  // check for relevant matchup of segments
  // ---------------------------------------------
  bool matchup = false;
  for (int i=0; i<36; ++i) {
    if ( (cntMEP13[i] > 0) && (cntMEP22[i]+cntMEP32[i] > 0) ) {matchup = true;}
    if ( (cntMEN13[i] > 0) && (cntMEN22[i]+cntMEN32[i] > 0) ) {matchup = true;}
  }
  /*
  if (matchup) {
    std::cout << "\tYYY looks like a good event.  Select!\n";
    std::cout << "-- pos endcap --\n"
     << "ME1/3: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP13[k];}
    std::cout << "\nME2/2: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP22[k];}
    std::cout << "\nME3/2: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP32[k];}
    std::cout << std::endl;
  }
  */

  // set the selection flag
  DTOverlapCandidate = matchup;
  return DTOverlapCandidate;
}

bool CSCSkim::doHaloLike ( edm::Handle< CSCSegmentCollection >  cscSegments )

private

Definition at line 849 of file CSCSkim.cc.

References CSCDetId, CSCDetId::endcap(), and i.

                                                                 {
  const float chisqMax = 100.;
  const int nhitsMin = 5; // on a segment
  const int maxNSegments = 3; // in a chamber

  // initialize
  bool HaloLike = false;
  int cntMEP11[36];
  int cntMEN11[36];
  int cntMEP12[36];
  int cntMEN12[36];
  int cntMEP21[36];
  int cntMEN21[36];
  int cntMEP31[36];
  int cntMEN31[36];
  int cntMEP41[36];
  int cntMEN41[36];
  for (int i=0; i<36; ++i) {
    cntMEP11[i] = 0;
    cntMEN11[i] = 0;
    cntMEP12[i] = 0;
    cntMEN12[i] = 0;
    cntMEP21[i] = 0;
    cntMEN21[i] = 0;
    cntMEP31[i] = 0;
    cntMEN31[i] = 0;
    cntMEP41[i] = 0;
    cntMEN41[i] = 0;
  }

  // -----------------------
  // loop over segments
  // -----------------------
  int nSegments = cscSegments->size();
  if (nSegments < 4) return HaloLike;

  for(CSCSegmentCollection::const_iterator it=cscSegments->begin(); it != cscSegments->end(); it++) {
    // which chamber?
    CSCDetId id  = (CSCDetId)(*it).cscDetId();
    int kEndcap  = id.endcap();
    int kStation = id.station();
    int kRing    = id.ring();
    int kChamber = id.chamber();
    // segment information
    float chisq    = (*it).chi2();
    int nhits      = (*it).nRecHits();
    bool goodSegment = (chisq < chisqMax) && (nhits >= nhitsMin) ;
    if (goodSegment) {
      if ( (kStation == 1) && (kRing == 1) ) {
    if (kEndcap == 1) {cntMEP11[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN11[kChamber-1]++;}
      }
      if ( (kStation == 1) && (kRing == 2) ) {
    if (kEndcap == 1) {cntMEP12[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN12[kChamber-1]++;}
      }
      if ( (kStation == 2) && (kRing == 1) ) {
    if (kEndcap == 1) {cntMEP21[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN21[kChamber-1]++;}
      }
      if ( (kStation == 3) && (kRing == 1) ) {
    if (kEndcap == 1) {cntMEP31[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN31[kChamber-1]++;}
      }
      if ( (kStation == 4) && (kRing == 1) ) {
    if (kEndcap == 1) {cntMEP41[kChamber-1]++;}
    if (kEndcap == 2) {cntMEN41[kChamber-1]++;}
      }
    } // this is a good segment
  } // end loop over segments

  // ---------------------------------------------
  // veto messy events
  // ---------------------------------------------
  bool tooManySegments = false;
  for (int i=0; i<36; ++i) {
    if ( (cntMEP11[i] > 3*maxNSegments) ||
         (cntMEN11[i] > 3*maxNSegments) ||
     (cntMEP12[i] > maxNSegments) ||
         (cntMEN12[i] > maxNSegments) ||
         (cntMEP21[i] > maxNSegments) ||
         (cntMEN21[i] > maxNSegments) ||
         (cntMEP31[i] > maxNSegments) ||
         (cntMEN31[i] > maxNSegments) ||
         (cntMEP41[i] > maxNSegments) ||
         (cntMEN41[i] > maxNSegments) ) tooManySegments = true;
  }
  if (tooManySegments) {
    return HaloLike;
  }

  // ---------------------------------------------
  // check for relevant matchup of segments
  // ---------------------------------------------
  bool matchup = false;
  for (int i=0; i<36; ++i) {
    if ( (cntMEP11[i]+cntMEP12[i] > 0) && 
         (cntMEP21[i] > 0) &&
         (cntMEP31[i] > 0) &&
         (cntMEP41[i] > 0) ) {matchup = true;}
    if ( (cntMEN11[i]+cntMEN12[i] > 0) && 
         (cntMEN21[i] > 0) &&
         (cntMEN31[i] > 0) &&
         (cntMEN41[i] > 0) ) {matchup = true;}
  }
  /*
  if (matchup) {
    std::cout << "\tYYY looks like a good event.  Select!\n";
    std::cout << "-- pos endcap --\n"
     << "ME1/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP11[k];}
    std::cout << "\nME1/2: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP12[k];}
    std::cout << "\nME2/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP21[k];}
    std::cout << "\nME3/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP31[k];}
    std::cout << "\nME4/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEP41[k];}
    std::cout << std::endl;
    std::cout << "-- neg endcap --\n"
     << "ME1/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEN11[k];}
    std::cout << "\nME1/2: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEN12[k];}
    std::cout << "\nME2/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEN21[k];}
    std::cout << "\nME3/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEN31[k];}
    std::cout << "\nME4/1: ";
    for (int k=0; k<36; ++k) {std::cout << " " << setw(3) << cntMEN41[k];}
    std::cout << std::endl;
    std::cout << "\tn Analyzed = " << nEventsAnalyzed << "\tn Halo-like = " << nEventsHaloLike << std::endl;
  }
  */

  // set the selection flag
  HaloLike = matchup;
  return HaloLike;
}

bool CSCSkim::doLongSATrack ( edm::Handle< reco::TrackCollection >  saTracks )

private

Definition at line 994 of file CSCSkim.cc.

References MuonSubdetId::CSC, MuonSubdetId::DT, min, DetId::Muon, metsig::muon, benchmark_cfg::select, and PV3DBase< T, PVType, FrameType >::z().

                                                                  {

  const float zDistanceMax = 2500.;
  const float zDistanceMin =  700.;
  const int nCSCHitsMin = 25;
  const int nCSCHitsMax = 50;
  const float zInnerMax = 80000.;

  const int nNiceMuonsMin = 1;

  //
  // Loop through the track collection and test each one
  //

  int nNiceMuons = 0;

  for(reco::TrackCollection::const_iterator muon = saMuons->begin(); muon != saMuons->end(); ++ muon ) {

    // basic information
    math::XYZVector innerMo = muon->innerMomentum();
    GlobalVector im(innerMo.x(),innerMo.y(),innerMo.z());
    math::XYZPoint innerPo = muon->innerPosition();
    GlobalPoint ip(innerPo.x(), innerPo.y(),innerPo.z());
    math::XYZPoint outerPo = muon->outerPosition();
    GlobalPoint op(outerPo.x(), outerPo.y(),outerPo.z());
    float zInner = ip.z();
    float zOuter = op.z();
    float zDistance = fabs(zOuter-zInner);



    // loop over hits
    int nDTHits = 0;
    int nCSCHits = 0;
    for (trackingRecHit_iterator hit = muon->recHitsBegin(); hit != muon->recHitsEnd(); ++hit ) {
      const DetId detId( (*hit)->geographicalId() );
      if (detId.det() == DetId::Muon) {
    if (detId.subdetId() == MuonSubdetId::DT) {
      //DTChamberId dtId(detId.rawId());
      //int chamberId = dtId.sector();
      nDTHits++;
    }
    else if (detId.subdetId() == MuonSubdetId::CSC) {
      //CSCDetId cscId(detId.rawId());
      //int chamberId = cscId.chamber();
      nCSCHits++;
    }
      }
    }

    // is this a nice muon?
    if ( (zDistance < zDistanceMax) && (zDistance > zDistanceMin) 
     && (nCSCHits > nCSCHitsMin) && (nCSCHits < nCSCHitsMax)
     && (min ( fabs(zInner), fabs(zOuter) ) < zInnerMax) 
     && (fabs(innerMo.z()) > 0.000000001) ) {
      nNiceMuons++;
    }
  }

  bool select = (nNiceMuons >= nNiceMuonsMin);

  return select;
}

bool CSCSkim::doMessyEventSkimming ( edm::Handle< CSCRecHit2DCollection >  cscRecHits, edm::Handle< CSCSegmentCollection >  cscSegments  )

private

Definition at line 564 of file CSCSkim.cc.

References CSCDetId::chamber(), CSCDetId, CSCDetId::endcap(), i, j, kLayer(), CSCDetId::layer(), LogDebug, funct::pow(), CSCDetId::ring(), and CSCDetId::station().

                                                                                                                          {

  // how many RecHits in the collection?
  int nRecHits = cscRecHits->size();

  // zero the recHit counter
  int cntRecHit[600];
  for (int i = 0; i < 600; i++) {
    cntRecHit[i] = 0;
  }

  // ---------------------
  // Loop over rechits 
  // ---------------------

  CSCRecHit2DCollection::const_iterator recIt;
  for (recIt = cscRecHits->begin(); recIt != cscRecHits->end(); recIt++) {

    // which chamber is it?
    CSCDetId idrec = (CSCDetId)(*recIt).cscDetId();
    int kEndcap  = idrec.endcap();
    int kRing    = idrec.ring();
    int kStation = idrec.station();
    int kChamber = idrec.chamber();
    int kLayer   = idrec.layer();

    // compute chamber serial number
    int kSerial = chamberSerial( kEndcap, kStation, kRing, kChamber ) ;

    // increment recHit counter
    //     (each layer is represented by a different power of 10)
    int kDigit = (int) pow((float)10.,(float)(kLayer-1));
    cntRecHit[kSerial] += kDigit;

  } //end rechit loop


  // ------------------------------------------------------
  // Are there chambers with the minimum number of hits?
  // ------------------------------------------------------

  int nChambersWithMinimalHits = 0;
  int nChambersWithMinimalHitsPOS = 0;
  int nChambersWithMinimalHitsNEG = 0;
  if (nRecHits > 0) {
    for (int i = 0; i < 600; i++) {
      if (cntRecHit[i] > 0) {
    int nLayersWithHits = 0;
    float dummy = (float) cntRecHit[i];
    for (int j = 5; j > -1; j--) {
      float digit = dummy / pow( (float)10., (float)j );
      int kCount = (int) digit;
      if (kCount > 0) nLayersWithHits++;
      dummy = dummy - ( (float) kCount) * pow( (float)10., (float)j );
    }
    if (nLayersWithHits > nLayersWithHitsMinimum) {
      if (i < 300) {nChambersWithMinimalHitsPOS++;}
      else         {nChambersWithMinimalHitsNEG++;}
    }
      }
    }
    nChambersWithMinimalHits = nChambersWithMinimalHitsPOS + nChambersWithMinimalHitsNEG;
  }

  // how many Segments?
  int nSegments = cscSegments->size();

  // ----------------------
  // fill histograms
  // ----------------------

  if (makeHistogramsForMessyEvents) {
    if (nRecHits > 8) {
      mevnRecHits0->Fill(nRecHits);
      mevnChambers0->Fill(nChambersWithMinimalHits);
      mevnSegments0->Fill(nSegments);
    }
    if (nRecHits > 54) {
      double dummy = (double) nRecHits;
      if (dummy > 299.9) dummy = 299.9;
      mevnRecHits1->Fill(dummy);
      dummy = (double) nChambersWithMinimalHits;
      if (dummy > 49.9) dummy = 49.9;
      mevnChambers1->Fill(dummy);
      dummy = (double) nSegments;
      if (dummy > 29.9) dummy = 29.9;
      mevnSegments1->Fill(dummy);
    }
  }

  // ----------------------
  // set the filter flag
  // ----------------------

  bool selectEvent = false;
  if ( (nRecHits > 54) && (nChambersWithMinimalHits > 5) ) {selectEvent = true;}

  // debug
  LogDebug("[CSCSkim]") << "----- nRecHits = " << nRecHits
               << "\tnChambersWithMinimalHits = " << nChambersWithMinimalHits
               << "\tnSegments = " << nSegments 
               << "\tselect? " << selectEvent << std::endl;

  /*
  if (selectEvent) {
    std::cout << "\n==========================================================================\n"
     << "\tmessy event!\n"
     << "\t  " <<  nEventsAnalyzed
     << "\trun " << iRun << "\tevent " << iEvent << std::endl;
    std::cout << "----- nRecHits = " << nRecHits
     << "\tnChambersWithMinimalHits = " << nChambersWithMinimalHits
     << "\tnSegments = " << nSegments 
     << "\tselect? " << selectEvent << std::endl;
    for (int i = 0; i < 600; i++) {
      if (cntRecHit[i] > 0) {
    cout << "\t\t" << i << "\tcntRecHit= " << cntRecHit[i] << std::endl;
      }
    }
    std::cout << "==========================================================================\n\n" ;
  }
  */

  return selectEvent;
}

bool CSCSkim::doOverlapSkimming ( edm::Handle< CSCSegmentCollection >  cscSegments )

private

Definition at line 482 of file CSCSkim.cc.

References CSCDetId, CSCDetId::endcap(), and i.

                                                                          {

  const int nhitsMinimum = 4;
  const float chisqMaximum = 100.;
  const int nAllMaximum  = 3;

  // how many Segments?
  //  int nSegments = cscSegments->size();

  // zero arrays
  int nAll[600];
  int nGood[600];
  for (int i=0; i<600; i++) {
    nAll[i] = 0;
    nGood[i] = 0;
  }

  // -----------------------
  // loop over segments
  // -----------------------
  for(CSCSegmentCollection::const_iterator it=cscSegments->begin(); it != cscSegments->end(); it++) {

    // which chamber?
    CSCDetId id  = (CSCDetId)(*it).cscDetId();
    int kEndcap  = id.endcap();
    int kStation = id.station();
    int kRing    = id.ring();
    int kChamber = id.chamber();
    int kSerial = chamberSerial( kEndcap, kStation, kRing, kChamber);

    // segment information
    float chisq    = (*it).chi2();
    int nhits      = (*it).nRecHits();

    // is this a good segment?
    bool goodSegment = (nhits >= nhitsMinimum) && (chisq < chisqMaximum) ;

    /*
    LocalPoint localPos = (*it).localPosition();
    float segX     = localPos.x();
    float segY     = localPos.y();
    std::cout << "E/S/R/Ch: " << kEndcap << "/" << kStation << "/" << kRing << "/" << kChamber
     << "\tnhits/chisq: " << nhits << "/" << chisq
     << "\tX/Y: " << segX << "/" << segY
     << "\tgood? " << goodSegment << std::endl;
    */

    // count
    nAll[kSerial-1]++;
    if (goodSegment) nGood[kSerial]++;

  } // end loop over segments

  //----------------------
  // select the event
  //----------------------

  // does any chamber have too many segments?
  bool messyChamber = false;
  for (int i = 0; i < 600; i++) {
    if (nAll[i] > nAllMaximum) messyChamber = true;
  }

  // are there consecutive chambers with good segments
  // (This is a little sloppy but is probably fine for skimming...)
  bool consecutiveChambers = false;
  for (int i = 0; i < 599; i++) {
    if ( (nGood[i]>0) && (nGood[i+1]>0) ) consecutiveChambers = true;
  }

  bool selectThisEvent = !messyChamber && consecutiveChambers;

  return selectThisEvent;

}

void CSCSkim::endJob ( void  )

virtual

Reimplemented from edm::EDFilter.

Definition at line 171 of file CSCSkim.cc.

References LogDebug.

                {

  // Write out results

  float fraction = 0.;
  if (nEventsAnalyzed > 0) {fraction = (float)nEventsSelected / (float)nEventsAnalyzed;}

  LogInfo("[CSCSkim] Summary")
    << "\n\n\t====== CSCSkim ==========================================================\n"
    << "\t\ttype of skim ...............................\t" << typeOfSkim << "\n"
    << "\t\tevents analyzed ..............\t" << nEventsAnalyzed << "\n"
    << "\t\tevents selected ..............\t" << nEventsSelected << "\tfraction= " << fraction << std::endl
    << "\t\tevents chambers both sides ...\t" << nEventsChambersBothSides << "\n"
    << "\t\tevents w/ overlaps .......... \t" << nEventsOverlappingChambers << "\n"
    << "\t\tevents lots of hit chambers . \t" << nEventsMessy << "\n"
    << "\t\tevents from certain chamber . \t" << nEventsCertainChamber << "\n"
    << "\t\tevents in DT-CSC overlap .... \t" << nEventsDTOverlap << "\n"
    << "\t\tevents halo-like ............ \t" << nEventsHaloLike << "\n"
    << "\t\tevents w/ long SA track ..... \t" << nEventsLongSATrack << "\n"
    << "\t\tevents good for BField  ..... \t" << nEventsForBFieldStudies << "\n"
    <<     "\t=========================================================================\n\n";

  if (makeHistograms || makeHistogramsForMessyEvents) {
    // Write the histos to file
    LogDebug("[CSCSkim]") << "======= write out my histograms ====\n" ;
    theHistogramFile->cd();
    if (makeHistograms) {
      hxnRecHits->Write();
      hxnSegments->Write();
      hxnHitChambers->Write();
      hxnRecHitsSel->Write();
    }
    if (makeHistogramsForMessyEvents) {
      mevnRecHits0->Write();
      mevnChambers0->Write();
      mevnSegments0->Write();
      mevnRecHits1->Write();
      mevnChambers1->Write();
      mevnSegments1->Write();
    }
    theHistogramFile->Close();
  }
}

bool CSCSkim::filter ( edm::Event &  event, const edm::EventSetup &  eventSetup  )

virtual

Implements edm::EDFilter.

Definition at line 219 of file CSCSkim.cc.

References cscSegments_cfi::cscSegments, edm::EventSetup::get(), iEvent, LogDebug, RecoTauPiZeroBuilderPlugins_cfi::strips, and testEve_cfg::tracks.

Referenced by python.Vispa.Plugins.Browser.BrowserTabController.BrowserTabController::filter(), python.Vispa.Plugins.Browser.BrowserTabController.BrowserTabController::find(), python.Vispa.Plugins.Browser.BrowserTabController.BrowserTabController::setDataAccessor(), and python.Vispa.Plugins.Browser.BrowserTabController.BrowserTabController::switchCenterView().

{
  // increment counter
  nEventsAnalyzed++;

  iRun   = event.id().run();
  iEvent = event.id().event();

  LogDebug("[CSCSkim] EventInfo") << "Run: " << iRun << "\tEvent: " << iEvent << "\tn Analyzed: " << nEventsAnalyzed;

  // Get the CSC Geometry :
  ESHandle<CSCGeometry> cscGeom;
  eventSetup.get<MuonGeometryRecord>().get(cscGeom);

  // Get the DIGI collections
  edm::Handle<CSCWireDigiCollection> wires;
  edm::Handle<CSCStripDigiCollection> strips;

  if (isSimulation){
    event.getByLabel("simMuonCSCDigis","MuonCSCWireDigi",wires);
    event.getByLabel("simMuonCSCDigis","MuonCSCStripDigi",strips);
  }
  else {
    event.getByLabel("muonCSCDigis","MuonCSCWireDigi",wires);
    event.getByLabel("muonCSCDigis","MuonCSCStripDigi",strips);
  }

  // Get the RecHits collection :
  Handle<CSCRecHit2DCollection> cscRecHits;
  event.getByLabel(cscRecHitTag,cscRecHits);

  // get CSC segment collection
  Handle<CSCSegmentCollection> cscSegments;
  event.getByLabel(cscSegmentTag, cscSegments);

  // get the cosmic muons collection
  Handle<reco::TrackCollection> saMuons;
  if (typeOfSkim == 8) {
    event.getByLabel(SAMuonTag,saMuons);
  }

  // get the stand-alone muons collection
  Handle<reco::TrackCollection> tracks;
  Handle<reco::MuonCollection> gMuons;
  if (typeOfSkim == 9) {
    event.getByLabel(SAMuonTag,saMuons);
    event.getByLabel(trackTag,tracks);
    event.getByLabel(GLBMuonTag,gMuons);
  }


  //======================================
  // evaluate the skimming routines
  //======================================


  // basic skimming
  bool basicEvent = false;
  if (typeOfSkim == 1 || typeOfSkim == 2) {
    basicEvent = doCSCSkimming(cscRecHits,cscSegments);
  }

  // overlapping chamber skim
  bool goodOverlapEvent = false;
  if (typeOfSkim == 3) {
    goodOverlapEvent = doOverlapSkimming(cscSegments);
    if (goodOverlapEvent) {nEventsOverlappingChambers++;}
  }

  // messy events skim
  bool messyEvent = false;
  if (typeOfSkim == 4) {
    messyEvent = doMessyEventSkimming(cscRecHits,cscSegments);
    if (messyEvent) {nEventsMessy++;}
  }

  // select events with DIGIs in a certain chamber
  bool hasChamber = false;
  if (typeOfSkim == 5) {
    hasChamber = doCertainChamberSelection(wires,strips);
    if (hasChamber) {nEventsCertainChamber++;}
  }

  // select events in the DT-CSC overlap region
  bool DTOverlapCandidate = false;
  if (typeOfSkim == 6) {
    DTOverlapCandidate = doDTOverlap(cscSegments);
    if (DTOverlapCandidate) {nEventsDTOverlap++;}
  }

  // select halo-like events
  bool HaloLike = false;
  if (typeOfSkim == 7) {
    HaloLike = doHaloLike(cscSegments);
    if (HaloLike) {nEventsHaloLike++;}
  }

  // select long cosmic tracks
  bool LongSATrack = false;
  if (typeOfSkim == 8) {
    LongSATrack = doLongSATrack(saMuons);
    if (LongSATrack) {nEventsLongSATrack++;}
  }

  // select events suitable for a B-field study.  They have tracks in the tracker.
  bool GoodForBFieldStudy = false;
  if (typeOfSkim == 9) {
    GoodForBFieldStudy = doBFieldStudySelection(saMuons,tracks,gMuons);
    if (GoodForBFieldStudy) {nEventsForBFieldStudies++;}
  }


  // set filter flag
  bool selectThisEvent = false;
  if (typeOfSkim == 1 || typeOfSkim == 2) {selectThisEvent = basicEvent;}
  if (typeOfSkim == 3) {selectThisEvent = goodOverlapEvent;}
  if (typeOfSkim == 4) {selectThisEvent = messyEvent;}
  if (typeOfSkim == 5) {selectThisEvent = hasChamber;}
  if (typeOfSkim == 6) {selectThisEvent = DTOverlapCandidate;}
  if (typeOfSkim == 7) {selectThisEvent = HaloLike;}
  if (typeOfSkim == 8) {selectThisEvent = LongSATrack;}
  if (typeOfSkim == 9) {selectThisEvent = GoodForBFieldStudy;}

  if (selectThisEvent) {nEventsSelected++;}

  return selectThisEvent;
}

Member Data Documentation

edm::InputTag CSCSkim::cscRecHitTag

private

Definition at line 145 of file CSCSkim.h.

edm::InputTag CSCSkim::cscSegmentTag

private

Definition at line 146 of file CSCSkim.h.

bool CSCSkim::demandChambersBothSides

private

Definition at line 157 of file CSCSkim.h.

edm::InputTag CSCSkim::GLBMuonTag

private

Definition at line 148 of file CSCSkim.h.

std::string CSCSkim::histogramFileName

private

Definition at line 142 of file CSCSkim.h.

TH1F* CSCSkim::hxnHitChambers

private

Definition at line 179 of file CSCSkim.h.

TH1F* CSCSkim::hxnRecHits

private

Definition at line 177 of file CSCSkim.h.

TH1F* CSCSkim::hxnRecHitsSel

private

Definition at line 180 of file CSCSkim.h.

TH1F* CSCSkim::hxnSegments

private

Definition at line 178 of file CSCSkim.h.

int CSCSkim::iEvent

private

Definition at line 135 of file CSCSkim.h.

int CSCSkim::iRun

private

Definition at line 134 of file CSCSkim.h.

bool CSCSkim::isSimulation

private

Definition at line 152 of file CSCSkim.h.

bool CSCSkim::makeHistograms

private

Definition at line 158 of file CSCSkim.h.

bool CSCSkim::makeHistogramsForMessyEvents

private

Definition at line 159 of file CSCSkim.h.

TH1F* CSCSkim::mevnChambers0

private

Definition at line 182 of file CSCSkim.h.

TH1F* CSCSkim::mevnChambers1

private

Definition at line 185 of file CSCSkim.h.

TH1F* CSCSkim::mevnRecHits0

private

Definition at line 181 of file CSCSkim.h.

TH1F* CSCSkim::mevnRecHits1

private

Definition at line 184 of file CSCSkim.h.

TH1F* CSCSkim::mevnSegments0

private

Definition at line 183 of file CSCSkim.h.

TH1F* CSCSkim::mevnSegments1

private

Definition at line 186 of file CSCSkim.h.

int CSCSkim::minimumHitChambers

private

Definition at line 155 of file CSCSkim.h.

int CSCSkim::minimumSegments

private

Definition at line 156 of file CSCSkim.h.

int CSCSkim::nCSCHitsMin

private

Definition at line 167 of file CSCSkim.h.

int CSCSkim::nEventsAnalyzed

private

Definition at line 122 of file CSCSkim.h.

int CSCSkim::nEventsCertainChamber

private

Definition at line 127 of file CSCSkim.h.

int CSCSkim::nEventsChambersBothSides

private

Definition at line 124 of file CSCSkim.h.

int CSCSkim::nEventsDTOverlap

private

Definition at line 128 of file CSCSkim.h.

int CSCSkim::nEventsForBFieldStudies

private

Definition at line 131 of file CSCSkim.h.

int CSCSkim::nEventsHaloLike

private

Definition at line 129 of file CSCSkim.h.

int CSCSkim::nEventsLongSATrack

private

Definition at line 130 of file CSCSkim.h.

int CSCSkim::nEventsMessy

private

Definition at line 126 of file CSCSkim.h.

int CSCSkim::nEventsOverlappingChambers

private

Definition at line 125 of file CSCSkim.h.

int CSCSkim::nEventsSelected

private

Definition at line 123 of file CSCSkim.h.

int CSCSkim::nLayersWithHitsMinimum

private

Definition at line 154 of file CSCSkim.h.

int CSCSkim::nTrHitsMin

private

Definition at line 169 of file CSCSkim.h.

int CSCSkim::nValidHitsMin

private

Definition at line 173 of file CSCSkim.h.

std::string CSCSkim::outputFileName

private

Definition at line 141 of file CSCSkim.h.

float CSCSkim::pMin

private

Definition at line 165 of file CSCSkim.h.

float CSCSkim::redChiSqMax

private

Definition at line 172 of file CSCSkim.h.

float CSCSkim::rExtMax

private

Definition at line 171 of file CSCSkim.h.

edm::InputTag CSCSkim::SAMuonTag

private

Definition at line 147 of file CSCSkim.h.

TFile* CSCSkim::theHistogramFile

private

Definition at line 138 of file CSCSkim.h.

edm::InputTag CSCSkim::trackTag

private

Definition at line 149 of file CSCSkim.h.

int CSCSkim::typeOfSkim

private

Definition at line 153 of file CSCSkim.h.

int CSCSkim::whichChamber

private

Definition at line 163 of file CSCSkim.h.

int CSCSkim::whichEndcap

private

Definition at line 160 of file CSCSkim.h.

int CSCSkim::whichRing

private

Definition at line 162 of file CSCSkim.h.

int CSCSkim::whichStation

private

Definition at line 161 of file CSCSkim.h.

TH1F * CSCSkim::xxnCSCHits

private

Definition at line 188 of file CSCSkim.h.

TH1F * CSCSkim::xxnTrackerHits

private

Definition at line 188 of file CSCSkim.h.

TH1F * CSCSkim::xxnValidHits

private

Definition at line 188 of file CSCSkim.h.

TH1F* CSCSkim::xxP

private

Definition at line 188 of file CSCSkim.h.

TH1F * CSCSkim::xxredChiSq

private

Definition at line 188 of file CSCSkim.h.

float CSCSkim::zInnerMax

private

Definition at line 168 of file CSCSkim.h.

float CSCSkim::zLengthMin

private

Definition at line 166 of file CSCSkim.h.

float CSCSkim::zLengthTrMin

private

Definition at line 170 of file CSCSkim.h.