CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
DTChamberEfficiencyTask Class Reference

#include <DTChamberEfficiencyTask.h>

Inheritance diagram for DTChamberEfficiencyTask:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &event, const edm::EventSetup &setup)
 
void beginJob ()
 BeginJob. More...
 
void beginLuminosityBlock (edm::LuminosityBlock const &lumiSeg, edm::EventSetup const &context)
 To reset the MEs. More...
 
void beginRun (const edm::Run &run, const edm::EventSetup &setup)
 BeginRun. More...
 
 DTChamberEfficiencyTask (const edm::ParameterSet &pset)
 Constructor. More...
 
void endJob ()
 Endjob. More...
 
virtual ~DTChamberEfficiencyTask ()
 Destructor. More...
 
- Public Member Functions inherited from edm::EDAnalyzer
 EDAnalyzer ()
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 
- Public Member Functions inherited from edm::EDConsumerBase
 EDConsumerBase ()
 
ProductHolderIndex indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndex > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndex > &) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Private Member Functions

void bookHistos (DTChamberId chId)
 
const DTRecSegment4DgetBestSegment (const DTRecSegment4DCollection::range &segs) const
 
const DTRecSegment4DgetBestSegment (const DTRecSegment4D *s1, const DTRecSegment4D *s2) const
 
LocalPoint interpolate (const DTRecSegment4D &seg1, const DTRecSegment4D &seg3, const DTChamberId &MB2) const
 
bool isGoodSegment (const DTRecSegment4D &seg) const
 

Private Attributes

bool debug
 
bool detailedAnalysis
 
edm::ESHandle< DTGeometrydtGeom
 
std::map< DTChamberId,
std::vector< MonitorElement * > > 
histosPerCh
 
bool onlineMonitor
 
edm::ParameterSet parameters
 
edm::Handle
< DTRecSegment4DCollection
segs
 
DQMStoretheDbe
 
double theMinChi2NormSegment
 
double theMinCloseDist
 
unsigned int theMinHitsSegment
 
std::string theRecHits4DLabel
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
typedef WorkerT< EDAnalyzerWorkerType
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 
- Protected Member Functions inherited from edm::EDAnalyzer
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
CurrentProcessingContext const * currentContext () const
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

DQM Analysis of 4D DT segments, it produces plots about:

Class based on the code written by S. Lacaprara : RecoLocalMuon / DTSegment / test / DTEffAnalyzer.h

Date:
2010/01/05 10:14:40
Revision:
1.7
Author
G. Mila - INFN Torino

Definition at line 39 of file DTChamberEfficiencyTask.h.

Constructor & Destructor Documentation

DTChamberEfficiencyTask::DTChamberEfficiencyTask ( const edm::ParameterSet pset)

Constructor.

Definition at line 36 of file DTChamberEfficiencyTask.cc.

References debug, edm::ParameterSet::getUntrackedParameter(), cppFunctionSkipper::operator, and Parameters::parameters.

36  {
37 
38  debug = pset.getUntrackedParameter<bool>("debug",false);
39 
40  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask") << "[DTChamberEfficiencyTask] Constructor called!";
41 
42  // Get the DQM needed services
44 
45  theDbe->setCurrentFolder("DT/DTChamberEfficiencyTask");
46 
47  parameters = pset;
48 
49 }
T getUntrackedParameter(std::string const &, T const &) const
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:434
DTChamberEfficiencyTask::~DTChamberEfficiencyTask ( )
virtual

Destructor.

Definition at line 52 of file DTChamberEfficiencyTask.cc.

52  {
53 
54  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask") << "[DTChamberEfficiencyTask] Destructor called!";
55 }

Member Function Documentation

void DTChamberEfficiencyTask::analyze ( const edm::Event event,
const edm::EventSetup setup 
)
virtual

Implements edm::EDAnalyzer.

Definition at line 170 of file DTChamberEfficiencyTask.cc.

References DTRecSegment4D::hasPhi(), DTRecSegment4D::hasZed(), mergeVDriftHistosByStation::histos, DTRecSegment4D::localPosition(), mag(), DTChamberId::sector(), DTChamberId::station(), relativeConstraints::station, DTChamberId::wheel(), x, PV3DBase< T, PVType, FrameType >::x(), detailsBasic3DVector::y, and PV3DBase< T, PVType, FrameType >::y().

170  {
171 
172  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask") << "[DTChamberEfficiencyTask] Analyze #Run: " << event.id().run()
173  << " #Event: " << event.id().event();
174 
175  // Get the 4D rechit collection from the event
176  event.getByLabel(theRecHits4DLabel, segs);
177 
178  int bottom=0, top=0;
179 
180 
181  // Loop over all the chambers
182  vector<DTChamber*>::const_iterator ch_it = dtGeom->chambers().begin();
183  vector<DTChamber*>::const_iterator ch_end = dtGeom->chambers().end();
184  for (; ch_it != ch_end; ++ch_it) {
185 
186  DTChamberId ch = (*ch_it)->id();
187  int wheel = ch.wheel();
188  int sector = ch.sector();
189  int station = ch.station();
190 
191 
192  DTChamberId MidId(wheel, station, sector);
193 
194  // get efficiency for MB1 using MB2 and MB3
195  if( station == 1 ) {
196  bottom = 2;
197  top = 3;
198  }
199 
200  // get efficiency for MB2 using MB1 and MB3
201  if( station == 2 ) {
202  bottom = 1;
203  top = 3;
204  }
205 
206  // get efficiency for MB2 using MB2 and MB4
207  if( station == 3 ) {
208  bottom = 2;
209  top = 4;
210  }
211 
212  // get efficiency for MB4 using MB2 and MB3
213  if( station == 4 ) {
214  bottom = 2;
215  top = 3;
216  }
217 
218  // Select events with (good) segments in Bot and Top
219  DTChamberId BotId(wheel, bottom, sector);
220  DTChamberId TopId(wheel, top, sector);
221 
222  // Get segments in the bottom chambers (if any)
223  DTRecSegment4DCollection::range segsBot= segs->get(BotId);
224  int nSegsBot=segsBot.second-segsBot.first;
225  // check if any segments is there
226  if (nSegsBot==0) continue;
227 
228  vector<MonitorElement *> histos = histosPerCh[MidId];
229 
230  // Get segments in the top chambers (if any)
231  DTRecSegment4DCollection::range segsTop= segs->get(TopId);
232  int nSegsTop=segsTop.second-segsTop.first;
233 
234  // Select one segment for the bottom chamber
235  const DTRecSegment4D& bestBotSeg= getBestSegment(segsBot);
236 
237  // Select one segment for the top chamber
238  DTRecSegment4D* pBestTopSeg=0;
239  if (nSegsTop>0)
240  pBestTopSeg = const_cast<DTRecSegment4D*>(&getBestSegment(segsTop));
241  //if top chamber is MB4 sector 10, consider also sector 14
242  if (TopId.station() == 4 && TopId.sector() == 10) {
243  DTChamberId TopId14(wheel, top, 14);
244  DTRecSegment4DCollection::range segsTop14= segs->get(TopId14);
245  int nSegsTop14=segsTop14.second-segsTop14.first;
246  nSegsTop+=nSegsTop14;
247  if (nSegsTop14) {
248  DTRecSegment4D* pBestTopSeg14 = const_cast<DTRecSegment4D*>(&getBestSegment(segsTop14));
249  // get best between sector 10 and 14
250  pBestTopSeg = const_cast<DTRecSegment4D*>(getBestSegment(pBestTopSeg, pBestTopSeg14));
251  }
252  }
253  if (!pBestTopSeg) continue;
254  const DTRecSegment4D& bestTopSeg= *pBestTopSeg;
255 
256  DTRecSegment4DCollection::range segsMid= segs->get(MidId);
257  int nSegsMid=segsMid.second-segsMid.first;
258 
259  if(detailedAnalysis){
260  // very trivial efficiency, just count segments
261  histos[3]->Fill(0);
262  if (nSegsMid>0) histos[3]->Fill(1);
263  }
264 
265  // get position at Mid by interpolating the position (not direction) of best
266  // segment in Bot and Top to Mid surface
267  LocalPoint posAtMid = interpolate(bestBotSeg, bestTopSeg, MidId);
268 
269  // is best segment good enough?
270  if (isGoodSegment(bestBotSeg) && isGoodSegment(bestTopSeg)) {
271  if(detailedAnalysis)
272  histos[4]->Fill(posAtMid.x(),posAtMid.y());
273  //check if interpolation fall inside middle chamber
274  if ((dtGeom->chamber(MidId))->surface().bounds().inside(posAtMid)) {
275  histos[0]->Fill(posAtMid.x(),posAtMid.y());
276  if (nSegsMid>0) {
277 
278  if(detailedAnalysis){
279  histos[3]->Fill(2);
280  histos[5]->Fill(posAtMid.x(),posAtMid.y());
281  }
282 
283  const DTRecSegment4D& bestMidSeg= getBestSegment(segsMid);
284  // check if middle segments is good enough
285  if (isGoodSegment(bestMidSeg)) {
286 
287  if(detailedAnalysis)
288  histos[6]->Fill(posAtMid.x(),posAtMid.y());
289  LocalPoint midSegPos=bestMidSeg.localPosition();
290 
291  // check if middle segments is also close enough
292  double dist;
293  if (bestMidSeg.hasPhi()) {
294  if (bestTopSeg.hasZed() && bestBotSeg.hasZed() && bestMidSeg.hasZed()) {
295  dist = (midSegPos-posAtMid).mag();
296  } else {
297  dist = fabs((midSegPos-posAtMid).x());
298  }
299  } else {
300  dist = fabs((midSegPos-posAtMid).y());
301  }
302  if (dist < theMinCloseDist ) {
303  histos[1]->Fill(posAtMid.x(),posAtMid.y());
304  }
305  if(detailedAnalysis)
306  histos[2]->Fill(dist);
307  }
308  }
309  }
310  }
311  }// loop over stations
312 
313 }
std::pair< const_iterator, const_iterator > range
iterator range
Definition: RangeMap.h:52
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
T y() const
Definition: PV3DBase.h:63
edm::Handle< DTRecSegment4DCollection > segs
LocalPoint interpolate(const DTRecSegment4D &seg1, const DTRecSegment4D &seg3, const DTChamberId &MB2) const
bool hasPhi() const
Does it have the Phi projection?
virtual LocalPoint localPosition() const
Local position in Chamber frame.
bool hasZed() const
Does it have the Z projection?
edm::ESHandle< DTGeometry > dtGeom
bool isGoodSegment(const DTRecSegment4D &seg) const
std::map< DTChamberId, std::vector< MonitorElement * > > histosPerCh
const DTRecSegment4D & getBestSegment(const DTRecSegment4DCollection::range &segs) const
int sector() const
Definition: DTChamberId.h:63
Definition: DDAxes.h:10
int station() const
Return the station number.
Definition: DTChamberId.h:53
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:47
T x() const
Definition: PV3DBase.h:62
void DTChamberEfficiencyTask::beginJob ( void  )
virtual

BeginJob.

Reimplemented from edm::EDAnalyzer.

Definition at line 58 of file DTChamberEfficiencyTask.cc.

References Parameters::parameters.

58  {
59 
60  // the name of the 4D rec hits collection
61  theRecHits4DLabel = parameters.getParameter<string>("recHits4DLabel");
62 
63  // parameters to use for the segment quality check
64  theMinHitsSegment = static_cast<unsigned int>(parameters.getParameter<int>("minHitsSegment"));
65  theMinChi2NormSegment = parameters.getParameter<double>("minChi2NormSegment");
66  // parameter to use for the exstrapolated segment check
67  theMinCloseDist = parameters.getParameter<double>("minCloseDist");
68 
69  // the running modality
70  onlineMonitor = parameters.getUntrackedParameter<bool>("onlineMonitor");
71 
72  // the analysis mode
73  detailedAnalysis = parameters.getUntrackedParameter<bool>("detailedAnalysis");
74 
75 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
void DTChamberEfficiencyTask::beginLuminosityBlock ( edm::LuminosityBlock const &  lumiSeg,
edm::EventSetup const &  context 
)
virtual

To reset the MEs.

Reimplemented from edm::EDAnalyzer.

Definition at line 78 of file DTChamberEfficiencyTask.cc.

References timingPdfMaker::histo, i, edm::LuminosityBlockBase::id(), edm::LuminosityBlockID::luminosityBlock(), Association::map, Parameters::parameters, and findQualityFiles::size.

78  {
79 
80  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask")<<"[DTChamberEfficiencyTask]: Begin of LS transition";
81 
82  if(lumiSeg.id().luminosityBlock()%parameters.getUntrackedParameter<int>("ResetCycle", 3) == 0 && onlineMonitor) {
83  for(map<DTChamberId, vector<MonitorElement*> > ::const_iterator histo = histosPerCh.begin();
84  histo != histosPerCh.end();
85  histo++) {
86  int size = (*histo).second.size();
87  for(int i=0; i<size; i++){
88  (*histo).second[i]->Reset();
89  }
90  }
91  }
92 
93 }
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
dictionary map
Definition: Association.py:205
std::map< DTChamberId, std::vector< MonitorElement * > > histosPerCh
tuple size
Write out results.
void DTChamberEfficiencyTask::beginRun ( const edm::Run run,
const edm::EventSetup setup 
)
virtual

BeginRun.

Reimplemented from edm::EDAnalyzer.

Definition at line 96 of file DTChamberEfficiencyTask.cc.

References bookHistos(), and edm::EventSetup::get().

96  {
97 
98  // Get the DT Geometry
99  setup.get<MuonGeometryRecord>().get(dtGeom);
100 
101  // Loop over all the chambers
102  vector<DTChamber*>::const_iterator ch_it = dtGeom->chambers().begin();
103  vector<DTChamber*>::const_iterator ch_end = dtGeom->chambers().end();
104  for (; ch_it != ch_end; ++ch_it) {
105  // histo booking
106  bookHistos((*ch_it)->id());
107  }
108 
109 }
void bookHistos(DTChamberId chId)
edm::ESHandle< DTGeometry > dtGeom
const T & get() const
Definition: EventSetup.h:55
void DTChamberEfficiencyTask::bookHistos ( DTChamberId  chId)
private

Definition at line 121 of file DTChamberEfficiencyTask.cc.

References mergeVDriftHistosByStation::histos, DTChamberId::sector(), DTChamberId::station(), relativeConstraints::station, and DTChamberId::wheel().

121  {
122 
123  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask") << " Booking histos for CH : " << chId;
124 
125  // Compose the chamber name
126  stringstream wheel; wheel << chId.wheel();
127  stringstream station; station << chId.station();
128  stringstream sector; sector << chId.sector();
129 
130  string HistoName =
131  "_W" + wheel.str() +
132  "_St" + station.str() +
133  "_Sec" + sector.str();
134 
135  theDbe->setCurrentFolder("DT/01-DTChamberEfficiency/Task/Wheel" + wheel.str() +
136  "/Sector" + sector.str() +
137  "/Station" + station.str());
138 
139  // Create the monitor elements
140  vector<MonitorElement *> histos;
141 
142  //efficiency selection cuts
143  // a- number of segments of the top chamber > 0 && number of segments of the bottom chamber > 0
144  // b- number of segments of the middle chamber > 0
145  // c- check of the top and bottom segment quality
146  // d- check if interpolation falls inside the middle chamber
147  // e- check of the middle segment quality
148  // f- check if the distance between the reconstructed and the exstrapolated segments is ok
149 
150 
151  // histo for efficiency with cuts a-/c-/d-
152  histos.push_back(theDbe->book2D("hEffGoodSegVsPosDen"+HistoName,"Eff vs local position (good) ",25,-250.,250., 25,-250.,250.));
153  // histo for efficiency with cuts a-/b-/c-/d-/e-/f-
154  histos.push_back(theDbe->book2D("hEffGoodCloseSegVsPosNum"+HistoName, "Eff vs local position (good and close segs) ", 25,-250.,250., 25,-250.,250.));
155  if(detailedAnalysis){
156  histos.push_back(theDbe->book1D("hDistSegFromExtrap"+HistoName, "Distance segments from extrap position ",200,0.,200.));
157  // histo for efficiency from segment counting
158  histos.push_back(theDbe->book1D("hNaiveEffSeg"+HistoName, "Naive eff ",10,0.,10.));
159  // histo for efficiency with cuts a-/c-
160  histos.push_back(theDbe->book2D("hEffSegVsPosDen"+HistoName,"Eff vs local position (all) ",25,-250.,250., 25,-250.,250.));
161  // histo for efficiency with cuts a-/b-/c-/d-
162  histos.push_back(theDbe->book2D("hEffSegVsPosNum"+HistoName, "Eff vs local position ",25,-250.,250., 25,-250.,250.));
163  // histo for efficiency with cuts a-/b-/c-/d-/e-
164  histos.push_back(theDbe->book2D("hEffGoodSegVsPosNum"+HistoName, "Eff vs local position (good segs) ", 25,-250.,250., 25,-250.,250.));
165  }
166  histosPerCh[chId] = histos;
167 }
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:722
std::map< DTChamberId, std::vector< MonitorElement * > > histosPerCh
std::string HistoName
int sector() const
Definition: DTChamberId.h:63
int station() const
Return the station number.
Definition: DTChamberId.h:53
MonitorElement * book2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D histogram.
Definition: DQMStore.cc:850
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:47
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:434
void DTChamberEfficiencyTask::endJob ( void  )
virtual

Endjob.

Reimplemented from edm::EDAnalyzer.

Definition at line 113 of file DTChamberEfficiencyTask.cc.

113  {
114 
115  edm::LogVerbatim ("DTDQM|DTMonitorModule|DTChamberEfficiencyTask")<<"[DTChamberEfficiencyTask] endjob called!";
116 }
const DTRecSegment4D & DTChamberEfficiencyTask::getBestSegment ( const DTRecSegment4DCollection::range segs) const
private

Definition at line 319 of file DTChamberEfficiencyTask.cc.

References DTRecSegment4D::chi2().

319  {
321  unsigned int nHitBest=0;
322  double chi2Best=99999.;
324  seg!=segs.second ; ++seg ) {
325  unsigned int nHits= ((*seg).hasPhi() ? (*seg).phiSegment()->recHits().size() : 0 ) ;
326  nHits+= ((*seg).hasZed() ? (*seg).zSegment()->recHits().size() : 0 );
327 
328  if (nHits==nHitBest) {
329  if ((*seg).chi2()/(*seg).degreesOfFreedom() < chi2Best ) {
330  chi2Best=(*seg).chi2()/(*seg).degreesOfFreedom();
331  bestIter = seg;
332  }
333  }
334  else if (nHits>nHitBest) {
335  nHitBest=nHits;
336  bestIter = seg;
337  }
338  }
339  return *bestIter;
340 }
edm::Handle< DTRecSegment4DCollection > segs
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:45
const DTRecSegment4D * DTChamberEfficiencyTask::getBestSegment ( const DTRecSegment4D s1,
const DTRecSegment4D s2 
) const
private

Definition at line 342 of file DTChamberEfficiencyTask.cc.

References DTRecSegment4D::chi2(), DTRecSegment4D::degreesOfFreedom(), DTRecSegment4D::hasPhi(), DTRecSegment4D::hasZed(), DTRecSegment4D::phiSegment(), DTRecSegment2D::recHits(), indexGen::s2, and DTRecSegment4D::zSegment().

343  {
344 
345  if (!s1) return s2;
346  if (!s2) return s1;
347  unsigned int nHits1= (s1->hasPhi() ? s1->phiSegment()->recHits().size() : 0 ) ;
348  nHits1+= (s1->hasZed() ? s1->zSegment()->recHits().size() : 0 );
349 
350  unsigned int nHits2= (s2->hasPhi() ? s2->phiSegment()->recHits().size() : 0 ) ;
351  nHits2+= (s2->hasZed() ? s2->zSegment()->recHits().size() : 0 );
352 
353  if (nHits1==nHits2) {
354  if (s1->chi2()/s1->degreesOfFreedom() < s2->chi2()/s2->degreesOfFreedom() )
355  return s1;
356  else
357  return s2;
358  }
359  else if (nHits1>nHits2) return s1;
360  return s2;
361 }
const DTChamberRecSegment2D * phiSegment() const
The superPhi segment: 0 if no phi projection available.
virtual int degreesOfFreedom() const
Degrees of freedom of the segment fit.
tuple s2
Definition: indexGen.py:106
virtual std::vector< const TrackingRecHit * > recHits() const
Access to component RecHits (if any)
bool hasPhi() const
Does it have the Phi projection?
bool hasZed() const
Does it have the Z projection?
const DTSLRecSegment2D * zSegment() const
The Z segment: 0 if not zed projection available.
virtual double chi2() const
Chi2 of the segment fit.
LocalPoint DTChamberEfficiencyTask::interpolate ( const DTRecSegment4D seg1,
const DTRecSegment4D seg3,
const DTChamberId MB2 
) const
private

Definition at line 364 of file DTChamberEfficiencyTask.cc.

References DTRecSegment4D::chamberId(), dir, DTRecSegment4D::hasPhi(), DTRecSegment4D::hasZed(), DTRecSegment4D::localPosition(), toLocal(), csvLumiCalc::unit, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

366  {
367  // Get GlobalPoition of Seg in MB1
368  GlobalPoint gpos1=(dtGeom->chamber(seg1.chamberId()))->toGlobal(seg1.localPosition());
369 
370  // Get GlobalPoition of Seg in MB3
371  GlobalPoint gpos3=(dtGeom->chamber(seg3.chamberId()))->toGlobal(seg3.localPosition());
372 
373  // interpolate
374  // get all in MB2 frame
375  LocalPoint pos1=(dtGeom->chamber(id2))->toLocal(gpos1);
376  LocalPoint pos3=(dtGeom->chamber(id2))->toLocal(gpos3);
377 
378  // case 1: 1 and 3 has both projection. No problem
379 
380  // case 2: one projection is missing for one of the segments. Keep the other's segment position
381  if (!seg1.hasZed()) pos1=LocalPoint(pos1.x(),pos3.y(),pos1.z());
382  if (!seg3.hasZed()) pos3=LocalPoint(pos3.x(),pos1.y(),pos3.z());
383 
384  if (!seg1.hasPhi()) pos1=LocalPoint(pos3.x(),pos1.y(),pos1.z());
385  if (!seg3.hasPhi()) pos3=LocalPoint(pos1.x(),pos3.y(),pos3.z());
386 
387  // direction
388  LocalVector dir = (pos3-pos1).unit(); // z points inward!
389  LocalPoint pos2 = pos1+dir*pos1.z()/(-dir.z());
390 
391  return pos2;
392 }
T y() const
Definition: PV3DBase.h:63
virtual DTChamberId chamberId() const
The (specific) DetId of the chamber on which the segment resides.
LocalVector toLocal(const reco::Track::Vector &v, const Surface &s)
string unit
Definition: csvLumiCalc.py:46
T z() const
Definition: PV3DBase.h:64
bool hasPhi() const
Does it have the Phi projection?
virtual LocalPoint localPosition() const
Local position in Chamber frame.
bool hasZed() const
Does it have the Z projection?
edm::ESHandle< DTGeometry > dtGeom
Local3DPoint LocalPoint
Definition: LocalPoint.h:11
dbl *** dir
Definition: mlp_gen.cc:35
T x() const
Definition: PV3DBase.h:62
bool DTChamberEfficiencyTask::isGoodSegment ( const DTRecSegment4D seg) const
private

Definition at line 395 of file DTChamberEfficiencyTask.cc.

References DTRecSegment4D::chamberId(), DTRecSegment4D::chi2(), DTRecSegment4D::degreesOfFreedom(), DTRecSegment4D::hasPhi(), DTRecSegment4D::hasZed(), DTRecSegment4D::phiSegment(), DTRecSegment2D::recHits(), DTChamberId::station(), and DTRecSegment4D::zSegment().

395  {
396  if (seg.chamberId().station()!=4 && !seg.hasZed()) return false;
397  unsigned int nHits= (seg.hasPhi() ? seg.phiSegment()->recHits().size() : 0 ) ;
398  nHits+= (seg.hasZed() ? seg.zSegment()->recHits().size() : 0 );
399  return ( nHits >= theMinHitsSegment &&
401 }
const DTChamberRecSegment2D * phiSegment() const
The superPhi segment: 0 if no phi projection available.
virtual DTChamberId chamberId() const
The (specific) DetId of the chamber on which the segment resides.
virtual int degreesOfFreedom() const
Degrees of freedom of the segment fit.
virtual std::vector< const TrackingRecHit * > recHits() const
Access to component RecHits (if any)
bool hasPhi() const
Does it have the Phi projection?
bool hasZed() const
Does it have the Z projection?
const DTSLRecSegment2D * zSegment() const
The Z segment: 0 if not zed projection available.
virtual double chi2() const
Chi2 of the segment fit.
int station() const
Return the station number.
Definition: DTChamberId.h:53

Member Data Documentation

bool DTChamberEfficiencyTask::debug
private
bool DTChamberEfficiencyTask::detailedAnalysis
private

Definition at line 84 of file DTChamberEfficiencyTask.h.

edm::ESHandle<DTGeometry> DTChamberEfficiencyTask::dtGeom
private

Definition at line 97 of file DTChamberEfficiencyTask.h.

std::map<DTChamberId, std::vector<MonitorElement*> > DTChamberEfficiencyTask::histosPerCh
private

Definition at line 91 of file DTChamberEfficiencyTask.h.

bool DTChamberEfficiencyTask::onlineMonitor
private

Definition at line 82 of file DTChamberEfficiencyTask.h.

edm::ParameterSet DTChamberEfficiencyTask::parameters
private
edm::Handle<DTRecSegment4DCollection> DTChamberEfficiencyTask::segs
private

Definition at line 98 of file DTChamberEfficiencyTask.h.

DQMStore* DTChamberEfficiencyTask::theDbe
private

Definition at line 77 of file DTChamberEfficiencyTask.h.

double DTChamberEfficiencyTask::theMinChi2NormSegment
private

Definition at line 94 of file DTChamberEfficiencyTask.h.

double DTChamberEfficiencyTask::theMinCloseDist
private

Definition at line 95 of file DTChamberEfficiencyTask.h.

unsigned int DTChamberEfficiencyTask::theMinHitsSegment
private

Definition at line 93 of file DTChamberEfficiencyTask.h.

std::string DTChamberEfficiencyTask::theRecHits4DLabel
private

Definition at line 87 of file DTChamberEfficiencyTask.h.