CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch1/src/CalibMuon/DTCalibration/src/DTSegmentSelector.cc

Go to the documentation of this file.
00001 #include "CalibMuon/DTCalibration/interface/DTSegmentSelector.h"
00002 
00003 #include "FWCore/Framework/interface/Event.h"
00004 #include "FWCore/Framework/interface/EventSetup.h"
00005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00006 #include "FWCore/Framework/interface/ESHandle.h"
00007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00008 
00009 #include "DataFormats/DTRecHit/interface/DTRecSegment2D.h"
00010 #include "DataFormats/DTRecHit/interface/DTRecSegment4D.h"
00011 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
00012 #include "DataFormats/DTRecHit/interface/DTRecHit1D.h"
00013 #include "CondFormats/DataRecord/interface/DTStatusFlagRcd.h"
00014 #include "CondFormats/DTObjects/interface/DTStatusFlag.h"
00015 
00016 bool DTSegmentSelector::operator() (DTRecSegment4D const& segment, edm::Event const& event, edm::EventSetup const& setup){
00017 
00018   bool result = true;
00019 
00020   /*
00021   // Get the DT Geometry
00022   ESHandle<DTGeometry> dtGeom;
00023   eventSetup.get<MuonGeometryRecord>().get(dtGeom);
00024   */
00025  
00026   edm::ESHandle<DTStatusFlag> statusMap;
00027   if(checkNoisyChannels_) setup.get<DTStatusFlagRcd>().get(statusMap);
00028 
00029   // Get the Phi 2D segment
00030   int nPhiHits = -1;
00031   bool segmentNoisyPhi = false;
00032   if( segment.hasPhi() ){
00033      const DTChamberRecSegment2D* phiSeg = segment.phiSegment();  // phiSeg lives in the chamber RF
00034      //LocalPoint phiSeg2DPosInCham = phiSeg->localPosition();  
00035      //LocalVector phiSeg2DDirInCham = phiSeg->localDirection();
00036      std::vector<DTRecHit1D> phiRecHits = phiSeg->specificRecHits();
00037      nPhiHits = phiRecHits.size(); 
00038      if(checkNoisyChannels_) segmentNoisyPhi = checkNoisySegment(statusMap,phiRecHits);
00039   }
00040   // Get the Theta 2D segment
00041   int nZHits = -1;
00042   bool segmentNoisyZ = false;
00043   if( segment.hasZed() ){
00044      const DTSLRecSegment2D* zSeg = segment.zSegment();  // zSeg lives in the SL RF
00045      //const DTSuperLayer* sl = chamber->superLayer(zSeg->superLayerId());
00046      //LocalPoint zSeg2DPosInCham = chamber->toLocal(sl->toGlobal((*zSeg).localPosition())); 
00047      //LocalVector zSeg2DDirInCham = chamber->toLocal(sl->toGlobal((*zSeg).localDirection()));
00048      std::vector<DTRecHit1D> zRecHits = zSeg->specificRecHits();
00049      nZHits = zRecHits.size();
00050      if(checkNoisyChannels_) segmentNoisyZ = checkNoisySegment(statusMap,zRecHits);
00051   } 
00052 
00053   // Segment selection 
00054   // Discard segment if it has a noisy cell
00055   if(segmentNoisyPhi || segmentNoisyZ)
00056      result = false;
00057 
00058   // 2D-segment number of hits
00059   if(segment.hasPhi() && nPhiHits < minHitsPhi_)
00060      result = false;
00061 
00062   if(segment.hasZed() && nZHits < minHitsZ_)
00063      result = false;
00064 
00065   // Segment chi2
00066   double chiSquare = segment.chi2()/segment.degreesOfFreedom();
00067   if(chiSquare > maxChi2_)
00068      result = false;
00069 
00070   // Segment angle
00071   LocalVector segment4DLocalDir = segment.localDirection();
00072   double angleZ = fabs( atan(segment4DLocalDir.y()/segment4DLocalDir.z())*180./Geom::pi() ); 
00073   if( angleZ > maxAngleZ_)
00074      result = false;
00075 
00076   double anglePhi = fabs( atan(segment4DLocalDir.x()/segment4DLocalDir.z())*180./Geom::pi() );
00077   if( anglePhi > maxAnglePhi_)
00078      result = false;
00079 
00080   return result;
00081 }
00082 
00083 bool DTSegmentSelector::checkNoisySegment(edm::ESHandle<DTStatusFlag> const& statusMap, std::vector<DTRecHit1D> const& dtHits){
00084 
00085   bool segmentNoisy = false;
00086 
00087   std::vector<DTRecHit1D>::const_iterator dtHit = dtHits.begin();
00088   std::vector<DTRecHit1D>::const_iterator dtHits_end = dtHits.end();
00089   for(; dtHit != dtHits_end; ++dtHit){
00090      //DTRecHit1D const* dtHit1D = dynamic_cast<DTRecHit1D const*>(*recHit);
00091      DTWireId wireId = dtHit->wireId();
00092      // Check for noisy channels to skip them
00093      bool isNoisy = false, isFEMasked = false, isTDCMasked = false, isTrigMask = false,
00094           isDead = false, isNohv = false;
00095      statusMap->cellStatus(wireId, isNoisy, isFEMasked, isTDCMasked, isTrigMask, isDead, isNohv);
00096      if(isNoisy) {
00097         LogTrace("Calibration") << "Wire: " << wireId << " is noisy, skipping!";
00098         segmentNoisy = true; break;
00099      }
00100   }
00101   return segmentNoisy;
00102 }