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
TkAlCaRecoMonitor Class Reference

#include <DQM/TrackerMonitorTrack/src/TkAlCaRecoMonitor.cc>

Inheritance diagram for TkAlCaRecoMonitor:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

virtual void analyze (const edm::Event &, const edm::EventSetup &)
 
virtual void beginJob ()
 
virtual void endJob (void)
 
 TkAlCaRecoMonitor (const edm::ParameterSet &)
 
 ~TkAlCaRecoMonitor ()
 
- 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 fillHitmaps (const reco::Track &track, const TrackerGeometry &geometry)
 
void fillRawIdMap (const TrackerGeometry &geometry)
 

Private Attributes

MonitorElementAlCaRecoTrackEfficiency_
 
std::map< int, int > binByRawId_
 
edm::ParameterSet conf_
 
double daughterMass_
 
DQMStoredqmStore_
 
bool fillInvariantMass_
 
bool fillRawIdMap_
 
MonitorElementHits_perDetId_
 
MonitorElementHits_XvsY_
 
MonitorElementHits_ZvsR_
 
MonitorElementinvariantMass_
 
MonitorElementjetPt_
 
double maxJetPt_
 
MonitorElementminJetDeltaR_
 
MonitorElementminTrackDeltaR_
 
edm::InputTag referenceTrackProducer_
 
bool runsOnReco_
 
MonitorElementsumCharge_
 
MonitorElementTrackCurvature_
 
edm::InputTag trackProducer_
 
MonitorElementTrackPtNegative_
 
MonitorElementTrackPtPositive_
 
MonitorElementTrackQuality_
 
bool useSignedR_
 

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

Monitoring special quantities related to Tracker Alignment AlCaReco Production.

Definition at line 28 of file TkAlCaRecoMonitor.h.

Constructor & Destructor Documentation

TkAlCaRecoMonitor::TkAlCaRecoMonitor ( const edm::ParameterSet iConfig)
explicit

Definition at line 28 of file TkAlCaRecoMonitor.cc.

References conf_, dqmStore_, and cppFunctionSkipper::operator.

28  {
30  conf_ = iConfig;
31 }
edm::ParameterSet conf_
TkAlCaRecoMonitor::~TkAlCaRecoMonitor ( )

Definition at line 33 of file TkAlCaRecoMonitor.cc.

33 { }

Member Function Documentation

void TkAlCaRecoMonitor::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
virtual

Implements edm::EDAnalyzer.

Definition at line 184 of file TkAlCaRecoMonitor.cc.

References AlCaRecoTrackEfficiency_, binByRawId_, conf_, daughterMass_, deltaR(), PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, alignCSCRings::e, MonitorElement::Fill(), fillHitmaps(), fillInvariantMass_, fillRawIdMap(), fillRawIdMap_, geometry, edm::EventSetup::get(), edm::Event::getByLabel(), edm::ParameterSet::getParameter(), i, invariantMass_, edm::ESHandleBase::isValid(), edm::HandleBase::isValid(), jetPt_, fwrapper::jets, maxJetPt_, minJetDeltaR_, minTrackDeltaR_, reco::TrackBase::qualitySize, referenceTrackProducer_, runsOnReco_, mathSSE::sqrt(), sumCharge_, TrackCurvature_, trackProducer_, TrackPtNegative_, TrackPtPositive_, and TrackQuality_.

184  {
185 
186  edm::Handle<reco::TrackCollection> trackCollection;
187  iEvent.getByLabel(trackProducer_, trackCollection);
188  if (!trackCollection.isValid()){
189  edm::LogError("Alignment")<<"invalid trackcollection encountered!";
190  return;
191  }
192 
193  edm::Handle<reco::TrackCollection> referenceTrackCollection;
194  iEvent.getByLabel(referenceTrackProducer_, referenceTrackCollection);
195  if (!trackCollection.isValid()){
196  edm::LogError("Alignment")<<"invalid reference track-collection encountered!";
197  return;
198  }
199 
201  iSetup.get<TrackerDigiGeometryRecord>().get(geometry);
202  if(! geometry.isValid()){
203  edm::LogError("Alignment")<<"invalid geometry found in event setup!";
204  }
205 
206  edm::ESHandle<MagneticField> magneticField;
207  iSetup.get<IdealMagneticFieldRecord>().get(magneticField);
208  if (!magneticField.isValid()){
209  edm::LogError("Alignment")<<"invalid magnetic field configuration encountered!";
210  return;
211  }
212 
214  if(runsOnReco_){
215  edm::InputTag jetCollection = conf_.getParameter<edm::InputTag>("CaloJetCollection");
216  iEvent.getByLabel(jetCollection, jets);
217  if(! jets.isValid()){
218  edm::LogError("Alignment")<<"no jet collection found in event!";
219  }
220  }
221  // fill only once - not yet in beginJob since no access to geometry
222  if (fillRawIdMap_ && binByRawId_.empty()) this->fillRawIdMap(*geometry);
223 
224  AlCaRecoTrackEfficiency_->Fill( static_cast<double>((*trackCollection).size()) / (*referenceTrackCollection).size() );
225 
226  double sumOfCharges = 0;
227  for( reco::TrackCollection::const_iterator track = (*trackCollection).begin(); track < (*trackCollection).end(); ++track ){
228  double dR = 0;
229  if(runsOnReco_){
230  double minJetDeltaR = 10; // some number > 2pi
231  for(reco::CaloJetCollection::const_iterator itJet = jets->begin(); itJet != jets->end() ; ++itJet) {
232  jetPt_->Fill( (*itJet).pt() );
233  dR = deltaR((*track),(*itJet));
234  if((*itJet).pt() > maxJetPt_ && dR < minJetDeltaR)
235  minJetDeltaR = dR;
236 
237  //edm::LogInfo("Alignment") <<"> isolated: "<< isolated << " jetPt "<< (*itJet).pt() <<" deltaR: "<< deltaR(*(*it),(*itJet)) ;
238  }
239  minJetDeltaR_->Fill( minJetDeltaR );
240  }
241 
242  double minTrackDeltaR = 10; // some number > 2pi
243  for( reco::TrackCollection::const_iterator track2 = (*trackCollection).begin(); track2 < (*trackCollection).end(); ++track2 ){
244  dR = deltaR((*track),(*track2));
245  if(dR < minTrackDeltaR && dR > 1e-6)
246  minTrackDeltaR = dR;
247  }
248 
249  for ( int i = 0; i<reco::TrackBase::qualitySize; ++i){
250  if( (*track).quality( reco::TrackBase::TrackQuality(i) ) ){
251  TrackQuality_->Fill( i );
252  }
253  }
254 
255  GlobalPoint gPoint((*track).vx(), (*track).vy(), (*track).vz());
256  double B = magneticField->inTesla(gPoint).z();
257  double curv = -(*track).charge()*0.002998*B/(*track).pt();
258  TrackCurvature_->Fill( curv );
259 
260  if( (*track).charge() > 0 )
261  TrackPtPositive_->Fill( (*track).pt() );
262  if( (*track).charge() < 0 )
263  TrackPtNegative_->Fill( (*track).pt() );
264 
265  minTrackDeltaR_->Fill( minTrackDeltaR );
266  fillHitmaps( *track, *geometry );
267  sumOfCharges += (*track).charge();
268  }
269 
270  sumCharge_->Fill( sumOfCharges );
271 
272  if(fillInvariantMass_){
273  if((*trackCollection).size() == 2){
274  TLorentzVector track0((*trackCollection).at(0).px(),(*trackCollection).at(0).py(),(*trackCollection).at(0).pz(),
275  sqrt(((*trackCollection).at(0).p()*(*trackCollection).at(0).p())+daughterMass_*daughterMass_));
276  TLorentzVector track1((*trackCollection).at(1).px(),(*trackCollection).at(1).py(),(*trackCollection).at(1).pz(),
277  sqrt(((*trackCollection).at(1).p()*(*trackCollection).at(1).p())+daughterMass_*daughterMass_));
278  TLorentzVector mother = track0+track1;
279 
280  invariantMass_->Fill( mother.M() );
281  }else{
282  edm::LogInfo("Alignment")<<"wrong number of tracks trackcollection encountered: "<<(*trackCollection).size();
283  }
284  }
285 }
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
edm::InputTag trackProducer_
void fillRawIdMap(const TrackerGeometry &geometry)
TrackQuality
track quality
Definition: TrackBase.h:95
MonitorElement * TrackQuality_
MonitorElement * minJetDeltaR_
std::map< int, int > binByRawId_
void Fill(long long x)
MonitorElement * AlCaRecoTrackEfficiency_
MonitorElement * sumCharge_
T sqrt(T t)
Definition: SSEVec.h:48
vector< PseudoJet > jets
MonitorElement * invariantMass_
MonitorElement * TrackPtNegative_
void fillHitmaps(const reco::Track &track, const TrackerGeometry &geometry)
bool isValid() const
Definition: HandleBase.h:76
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:361
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
const T & get() const
Definition: EventSetup.h:55
MonitorElement * TrackCurvature_
MonitorElement * jetPt_
edm::ParameterSet conf_
ESHandle< TrackerGeometry > geometry
MonitorElement * minTrackDeltaR_
bool isValid() const
Definition: ESHandle.h:37
MonitorElement * TrackPtPositive_
edm::InputTag referenceTrackProducer_
void TkAlCaRecoMonitor::beginJob ( void  )
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 35 of file TkAlCaRecoMonitor.cc.

References AlCaRecoTrackEfficiency_, DQMStore::book1D(), DQMStore::book2D(), conf_, daughterMass_, dqmStore_, fillInvariantMass_, fillRawIdMap_, edm::ParameterSet::getParameter(), MonitorElement::getTH1(), Hits_perDetId_, Hits_XvsY_, Hits_ZvsR_, i, invariantMass_, jetPt_, edm::InputTag::label(), maxJetPt_, minJetDeltaR_, minTrackDeltaR_, NULL, reco::TrackBase::qualityName(), reco::TrackBase::qualitySize, referenceTrackProducer_, runsOnReco_, MonitorElement::setAxisTitle(), DQMStore::setCurrentFolder(), AlCaHLTBitMon_QueryRunRegistry::string, sumCharge_, TrackCurvature_, trackProducer_, TrackPtNegative_, TrackPtPositive_, TrackQuality_, and useSignedR_.

35  {
36 
37  std::string histname; //for naming the histograms according to algorithm used
38 
39  trackProducer_ = conf_.getParameter<edm::InputTag>("TrackProducer");
40  referenceTrackProducer_ = conf_.getParameter<edm::InputTag>("ReferenceTrackProducer");
41 
42  std::string AlgoName = conf_.getParameter<std::string>("AlgoName");
43  std::string MEFolderName = conf_.getParameter<std::string>("FolderName");
44 
45  daughterMass_ = conf_.getParameter<double>("daughterMass");
46 
47  maxJetPt_ = conf_.getParameter<double>("maxJetPt");
48 
49  dqmStore_->setCurrentFolder(MEFolderName+"/TkAlignmentSpecific");
50  fillInvariantMass_ = conf_.getParameter<bool>("fillInvariantMass");
51  runsOnReco_ = conf_.getParameter<bool>("runsOnReco");
52  useSignedR_ = conf_.getParameter<bool>("useSignedR");
53  fillRawIdMap_ = conf_.getParameter<bool>("fillRawIdMap");
54 
55  //
56  unsigned int MassBin = conf_.getParameter<unsigned int>("MassBin");
57  double MassMin = conf_.getParameter<double>("MassMin");
58  double MassMax = conf_.getParameter<double>("MassMax");
59 
61  histname = "InvariantMass_";
62  invariantMass_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, MassBin, MassMin, MassMax);
63  invariantMass_->setAxisTitle("invariant Mass / GeV");
64  }else{
65  invariantMass_ = 0;
66  }
67 
68  unsigned int TrackPtPositiveBin = conf_.getParameter<unsigned int>("TrackPtBin");
69  double TrackPtPositiveMin = conf_.getParameter<double>("TrackPtMin");
70  double TrackPtPositiveMax = conf_.getParameter<double>("TrackPtMax");
71 
72  histname = "TrackPtPositive_";
73  TrackPtPositive_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, TrackPtPositiveBin, TrackPtPositiveMin, TrackPtPositiveMax);
74  TrackPtPositive_->setAxisTitle("p_{T} of tracks charge > 0");
75 
76  unsigned int TrackPtNegativeBin = conf_.getParameter<unsigned int>("TrackPtBin");
77  double TrackPtNegativeMin = conf_.getParameter<double>("TrackPtMin");
78  double TrackPtNegativeMax = conf_.getParameter<double>("TrackPtMax");
79 
80  histname = "TrackPtNegative_";
81  TrackPtNegative_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, TrackPtNegativeBin, TrackPtNegativeMin, TrackPtNegativeMax);
82  TrackPtNegative_->setAxisTitle("p_{T} of tracks charge < 0");
83 
84  histname = "TrackQuality_";
85  TrackQuality_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName,
87  TrackQuality_->setAxisTitle("quality");
88  for ( int i = 0; i<reco::TrackBase::qualitySize; ++i){
89  TrackQuality_->getTH1()->GetXaxis()->SetBinLabel(i+1,
91  }
92 
93  unsigned int SumChargeBin = conf_.getParameter<unsigned int>("SumChargeBin");
94  double SumChargeMin = conf_.getParameter<double>("SumChargeMin");
95  double SumChargeMax = conf_.getParameter<double>("SumChargeMax");
96 
97  histname = "SumCharge_";
98  sumCharge_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, SumChargeBin, SumChargeMin, SumChargeMax);
99  sumCharge_->setAxisTitle("#SigmaCharge");
100 
101  unsigned int TrackCurvatureBin = conf_.getParameter<unsigned int>("TrackCurvatureBin");
102  double TrackCurvatureMin = conf_.getParameter<double>("TrackCurvatureMin");
103  double TrackCurvatureMax = conf_.getParameter<double>("TrackCurvatureMax");
104 
105  histname = "TrackCurvature_";
106  TrackCurvature_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, TrackCurvatureBin, TrackCurvatureMin, TrackCurvatureMax);
107  TrackCurvature_->setAxisTitle("#kappa track");
108 
109  if( runsOnReco_ ){
110 
111  unsigned int JetPtBin = conf_.getParameter<unsigned int>("JetPtBin");
112  double JetPtMin = conf_.getParameter<double>("JetPtMin");
113  double JetPtMax = conf_.getParameter<double>("JetPtMax");
114 
115  histname = "JetPt_";
116  jetPt_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, JetPtBin, JetPtMin, JetPtMax);
117  jetPt_->setAxisTitle("jet p_{T} / GeV");
118 
119  unsigned int MinJetDeltaRBin = conf_.getParameter<unsigned int>("MinJetDeltaRBin");
120  double MinJetDeltaRMin = conf_.getParameter<double>("MinJetDeltaRMin");
121  double MinJetDeltaRMax = conf_.getParameter<double>("MinJetDeltaRMax");
122 
123  histname = "MinJetDeltaR_";
124  minJetDeltaR_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, MinJetDeltaRBin, MinJetDeltaRMin, MinJetDeltaRMax);
125  minJetDeltaR_->setAxisTitle("minimal Jet #DeltaR / rad");
126  }else{
127  jetPt_ = NULL;
129  }
130 
131  unsigned int MinTrackDeltaRBin = conf_.getParameter<unsigned int>("MinTrackDeltaRBin");
132  double MinTrackDeltaRMin = conf_.getParameter<double>("MinTrackDeltaRMin");
133  double MinTrackDeltaRMax = conf_.getParameter<double>("MinTrackDeltaRMax");
134 
135  histname = "MinTrackDeltaR_";
136  minTrackDeltaR_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, MinTrackDeltaRBin, MinTrackDeltaRMin, MinTrackDeltaRMax);
137  minTrackDeltaR_->setAxisTitle("minimal Track #DeltaR / rad");
138 
139  unsigned int TrackEfficiencyBin = conf_.getParameter<unsigned int>("TrackEfficiencyBin");
140  double TrackEfficiencyMin = conf_.getParameter<double>("TrackEfficiencyMin");
141  double TrackEfficiencyMax = conf_.getParameter<double>("TrackEfficiencyMax");
142 
143  histname = "AlCaRecoTrackEfficiency_";
144  AlCaRecoTrackEfficiency_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, TrackEfficiencyBin, TrackEfficiencyMin, TrackEfficiencyMax);
145  AlCaRecoTrackEfficiency_->setAxisTitle("n("+trackProducer_.label()+") / n("+referenceTrackProducer_.label()+")");
146 
147  int zBin = conf_.getParameter<unsigned int>("HitMapsZBin"); //300
148  double zMax = conf_.getParameter<double>("HitMapZMax"); //300.0; //cm
149 
150  int rBin = conf_.getParameter<unsigned int>("HitMapsRBin");//120;
151  double rMax = conf_.getParameter<double>("HitMapRMax"); //120.0; //cm
152 
153  histname = "Hits_ZvsR_";
154  double rMin = 0.0;
155  if( useSignedR_ )
156  rMin = -rMax;
157 
158  Hits_ZvsR_ = dqmStore_->book2D(histname+AlgoName, histname+AlgoName, zBin, -zMax, zMax, rBin, rMin, rMax);
159 
160  histname = "Hits_XvsY_";
161  Hits_XvsY_ = dqmStore_->book2D(histname+AlgoName, histname+AlgoName, rBin, -rMax, rMax, rBin, -rMax, rMax);
162 
163  if( fillRawIdMap_){
164  histname = "Hits_perDetId_";
165 
166  //leads to differences in axsis between samples??
167  //int nModules = binByRawId_.size();
168  //Hits_perDetId_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, nModules, static_cast<double>(nModules) -0.5, static_cast<double>(nModules) -0.5);
169  Hits_perDetId_ = dqmStore_->book1D(histname+AlgoName, histname+AlgoName, 16601,-0.5, 16600.5 );
170  Hits_perDetId_->setAxisTitle("rawId Bins");
171 
173  // std::stringstream binLabel;
174  // for( std::map<int,int>::iterator it = binByRawId_.begin(); it != binByRawId_.end(); ++it ){
175  // binLabel.str() = "";
176  // binLabel << (*it).first;
177  // Hits_perDetId_->getTH1()->GetXaxis()->SetBinLabel( (*it).second +1, binLabel.str().c_str());
178  // }
179  }
180 }
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
static std::string qualityName(TrackQuality)
Definition: TrackBase.h:406
edm::InputTag trackProducer_
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:722
TrackQuality
track quality
Definition: TrackBase.h:95
MonitorElement * TrackQuality_
#define NULL
Definition: scimark2.h:8
MonitorElement * minJetDeltaR_
MonitorElement * Hits_XvsY_
MonitorElement * Hits_ZvsR_
MonitorElement * AlCaRecoTrackEfficiency_
MonitorElement * sumCharge_
MonitorElement * invariantMass_
MonitorElement * TrackPtNegative_
TH1 * getTH1(void) const
MonitorElement * Hits_perDetId_
MonitorElement * TrackCurvature_
MonitorElement * jetPt_
edm::ParameterSet conf_
std::string const & label() const
Definition: InputTag.h:42
MonitorElement * minTrackDeltaR_
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
MonitorElement * TrackPtPositive_
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
edm::InputTag referenceTrackProducer_
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:434
void TkAlCaRecoMonitor::endJob ( void  )
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 329 of file TkAlCaRecoMonitor.cc.

References conf_, dqmStore_, edm::ParameterSet::getParameter(), dumpDBToFile_GT_ttrig_cfg::outputFileName, DQMStore::save(), and AlCaHLTBitMon_QueryRunRegistry::string.

329  {
330  bool outputMEsInRootFile = conf_.getParameter<bool>("OutputMEsInRootFile");
332  if(outputMEsInRootFile){
333  //dqmStore_->showDirStructure();
334  dqmStore_->save(outputFileName);
335  }
336 }
T getParameter(std::string const &) const
void save(const std::string &filename, const std::string &path="", const std::string &pattern="", const std::string &rewrite="", SaveReferenceTag ref=SaveWithReference, int minStatus=dqm::qstatus::STATUS_OK, const std::string &fileupdate="RECREATE")
Definition: DQMStore.cc:2118
edm::ParameterSet conf_
void TkAlCaRecoMonitor::fillHitmaps ( const reco::Track track,
const TrackerGeometry geometry 
)
private

Definition at line 287 of file TkAlCaRecoMonitor.cc.

References binByRawId_, MonitorElement::Fill(), fillRawIdMap_, TrackingRecHit::geographicalId(), Hits_perDetId_, Hits_XvsY_, Hits_ZvsR_, TrackerGeometry::idToDet(), alignCSCRings::r, reco::Track::recHitsBegin(), reco::Track::recHitsEnd(), mathSSE::sqrt(), and useSignedR_.

Referenced by analyze().

288 {
289  for (trackingRecHit_iterator iHit = track.recHitsBegin(); iHit != track.recHitsEnd(); ++iHit) {
290  if( (*iHit)->isValid() ){
291  const TrackingRecHit *hit = (*iHit).get();
292  const DetId geoId(hit->geographicalId());
293  const GeomDet* gd = geometry.idToDet(geoId);
294  // since 2_1_X local hit positions are transient. taking center of the hit module for now.
295  // The alternative would be the coarse estimation or a refit.
296  //const GlobalPoint globP( gd->toGlobal( hit->localPosition() ) );
297  const GlobalPoint globP( gd->toGlobal( Local3DPoint(0.,0.,0.) ) );
298  double r = sqrt( globP.x()*globP.x() + globP.y()*globP.y() );
299  if( useSignedR_ )
300  r*= globP.y() / fabs( globP.y() );
301 
302  Hits_ZvsR_->Fill( globP.z(), r );
303  Hits_XvsY_->Fill( globP.x(), globP.y() );
304 
305  if( fillRawIdMap_)
306  Hits_perDetId_->Fill( binByRawId_[ geoId.rawId() ]);
307  }
308  }
309 }
std::map< int, int > binByRawId_
MonitorElement * Hits_XvsY_
MonitorElement * Hits_ZvsR_
void Fill(long long x)
T sqrt(T t)
Definition: SSEVec.h:48
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:63
virtual const GeomDet * idToDet(DetId) const
Definition: DetId.h:20
MonitorElement * Hits_perDetId_
DetId geographicalId() const
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:65
void TkAlCaRecoMonitor::fillRawIdMap ( const TrackerGeometry geometry)
private

Definition at line 311 of file TkAlCaRecoMonitor.cc.

References binByRawId_, TrackerGeometry::detUnitIds(), i, and python.multivaluedict::sort().

Referenced by analyze().

312 {
313  std::vector<int> sortedRawIds;
314  for (std::vector<DetId>::const_iterator iDetId = geometry.detUnitIds().begin();
315  iDetId != geometry.detUnitIds().end(); ++iDetId) {
316  sortedRawIds.push_back((*iDetId).rawId());
317  }
318  std::sort(sortedRawIds.begin(), sortedRawIds.end());
319 
320  int i = 0;
321  for (std::vector<int>::iterator iRawId = sortedRawIds.begin();
322  iRawId != sortedRawIds.end(); ++iRawId){
323  binByRawId_[*iRawId] = i;
324  ++i;
325  }
326 }
int i
Definition: DBlmapReader.cc:9
std::map< int, int > binByRawId_
virtual const DetIdContainer & detUnitIds() const
Returm a vector of all GeomDetUnit DetIds.

Member Data Documentation

MonitorElement* TkAlCaRecoMonitor::AlCaRecoTrackEfficiency_
private

Definition at line 54 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

std::map<int,int> TkAlCaRecoMonitor::binByRawId_
private

Definition at line 71 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), fillHitmaps(), and fillRawIdMap().

edm::ParameterSet TkAlCaRecoMonitor::conf_
private

Definition at line 43 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), beginJob(), endJob(), and TkAlCaRecoMonitor().

double TkAlCaRecoMonitor::daughterMass_
private

Definition at line 70 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

DQMStore* TkAlCaRecoMonitor::dqmStore_
private

Definition at line 42 of file TkAlCaRecoMonitor.h.

Referenced by beginJob(), endJob(), and TkAlCaRecoMonitor().

bool TkAlCaRecoMonitor::fillInvariantMass_
private

Definition at line 63 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

bool TkAlCaRecoMonitor::fillRawIdMap_
private

Definition at line 64 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), beginJob(), and fillHitmaps().

MonitorElement* TkAlCaRecoMonitor::Hits_perDetId_
private

Definition at line 55 of file TkAlCaRecoMonitor.h.

Referenced by beginJob(), and fillHitmaps().

MonitorElement* TkAlCaRecoMonitor::Hits_XvsY_
private

Definition at line 61 of file TkAlCaRecoMonitor.h.

Referenced by beginJob(), and fillHitmaps().

MonitorElement* TkAlCaRecoMonitor::Hits_ZvsR_
private

Definition at line 60 of file TkAlCaRecoMonitor.h.

Referenced by beginJob(), and fillHitmaps().

MonitorElement* TkAlCaRecoMonitor::invariantMass_
private

Definition at line 48 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::jetPt_
private

Definition at line 51 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

double TkAlCaRecoMonitor::maxJetPt_
private

Definition at line 45 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::minJetDeltaR_
private

Definition at line 52 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::minTrackDeltaR_
private

Definition at line 53 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

edm::InputTag TkAlCaRecoMonitor::referenceTrackProducer_
private

Definition at line 69 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

bool TkAlCaRecoMonitor::runsOnReco_
private

Definition at line 65 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::sumCharge_
private

Definition at line 49 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::TrackCurvature_
private

Definition at line 58 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

edm::InputTag TkAlCaRecoMonitor::trackProducer_
private

Definition at line 68 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::TrackPtNegative_
private

Definition at line 57 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::TrackPtPositive_
private

Definition at line 56 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

MonitorElement* TkAlCaRecoMonitor::TrackQuality_
private

Definition at line 50 of file TkAlCaRecoMonitor.h.

Referenced by analyze(), and beginJob().

bool TkAlCaRecoMonitor::useSignedR_
private

Definition at line 66 of file TkAlCaRecoMonitor.h.

Referenced by beginJob(), and fillHitmaps().