CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MonitorTrackResiduals.cc
Go to the documentation of this file.
27 
29  : dqmStore_( edm::Service<DQMStore>().operator->() )
30  , conf_(iConfig), m_cacheID_(0)
31  , genTriggerEventFlag_(new GenericTriggerEventFlag(iConfig)) {
32 }
33 
36 }
37 
38 
40  ModOn = conf_.getParameter<bool>("Mod_On");
41  reset_me_after_each_run = conf_.getParameter<bool>("ResetAfterRun");
42 }
43 
45 
46 
47  unsigned long long cacheID = iSetup.get<SiStripDetCablingRcd>().cacheIdentifier();
48  if (m_cacheID_ != cacheID) {
49  m_cacheID_ = cacheID;
50  this->createMEs(iSetup);
51  }
53  if(ModOn) {
54  for(std::map<int32_t, MonitorElement*>::const_iterator it = HitResidual.begin(),
55  itEnd = HitResidual.end(); it!= itEnd;++it) {
56  this->resetModuleMEs(it->first);
58  }
59  } else {
60  for(std::map< std::pair<std::string,int32_t>, MonitorElement*>::const_iterator it = m_SubdetLayerResiduals.begin(),
61  itEnd = m_SubdetLayerResiduals.end(); it!= itEnd;++it) {
62  this->resetLayerMEs(it->first);
63  }
64  } // end if-else Module level on
65  } // end reset after run
66 
67  // Initialize the GenericTriggerEventFlag
68  if ( genTriggerEventFlag_->on() ) genTriggerEventFlag_->initRun( run, iSetup );
69 }
70 
72  Parameters = conf_.getParameter<edm::ParameterSet>("TH1ResModules");
73  int32_t i_residuals_Nbins = Parameters.getParameter<int32_t>("Nbinx");
74  double d_residual_xmin = Parameters.getParameter<double>("xmin");
75  double d_residual_xmax = Parameters.getParameter<double>("xmax");
76  Parameters = conf_.getParameter<edm::ParameterSet>("TH1NormResModules");
77  int32_t i_normres_Nbins = Parameters.getParameter<int32_t>("Nbinx");
78  double d_normres_xmin = Parameters.getParameter<double>("xmin");
79  double d_normres_xmax = Parameters.getParameter<double>("xmax");
80 
81 
82  // use SistripHistoId for producing histogram id (and title)
83  SiStripHistoId hidmanager;
84  folder_organizer.setSiStripFolder(); // top SiStrip folder
85 
86  // take from eventSetup the SiStripDetCabling object
88  iSetup.get<SiStripDetCablingRcd>().get(tkmechstruct);
89 
90  // get list of active detectors from SiStripDetCabling
91  std::vector<uint32_t> activeDets;
92  activeDets.clear(); // just in case
93  tkmechstruct->addActiveDetectorsRawIds(activeDets);
94 
95  // use SiStripSubStructure for selecting certain regions
96  SiStripSubStructure substructure;
97  std::vector<uint32_t> DetIds = activeDets;
98 
99  // book histo per each detector module
100  for (std::vector<uint32_t>::const_iterator DetItr=activeDets.begin(),
101  DetItrEnd = activeDets.end(); DetItr!=DetItrEnd; ++DetItr)
102  {
103  uint ModuleID = (*DetItr);
104 
105  // is this a StripModule?
106  if( SiStripDetId(ModuleID).subDetector() != 0 ) {
107 
108  folder_organizer.setDetectorFolder(ModuleID); // detid sets appropriate detector folder
109  // Book module histogramms?
110  if (ModOn) {
111  std::string hid = hidmanager.createHistoId("HitResiduals","det",ModuleID);
112  std::string normhid = hidmanager.createHistoId("NormalizedHitResiduals","det",ModuleID);
113  HitResidual[ModuleID] = dqmStore_->book1D(hid, hid,
114  i_residuals_Nbins,d_residual_xmin,d_residual_xmax);
115  HitResidual[ModuleID]->setAxisTitle("(x_{pred} - x_{rec})' [cm]");
116  NormedHitResiduals[ModuleID] = dqmStore_->book1D(normhid, normhid,
117  i_normres_Nbins,d_normres_xmin,d_normres_xmax);
118  NormedHitResiduals[ModuleID]->setAxisTitle("(x_{pred} - x_{rec})'/#sigma");
119  }
120  // book layer level histogramms
121  std::pair<std::string,int32_t> subdetandlayer = folder_organizer.GetSubDetAndLayer(ModuleID);
122  folder_organizer.setLayerFolder(ModuleID,subdetandlayer.second);
123  if(! m_SubdetLayerResiduals[subdetandlayer ] ) {
124  // book histogramms on layer level, check for barrel for correct labeling
125  std::string histoname = Form(subdetandlayer.first.find("B") != std::string::npos ?
126  "HitResiduals_%s__Layer__%d" : "HitResiduals_%s__wheel__%d" ,
127  subdetandlayer.first.c_str(),std::abs(subdetandlayer.second));
128  std::string normhistoname =
129  Form(subdetandlayer.first.find("B") != std::string::npos ?
130  "NormalizedHitResidual_%s__Layer__%d" : "NormalizedHitResidual_%s__wheel__%d" ,
131  subdetandlayer.first.c_str(),std::abs(subdetandlayer.second));
132  m_SubdetLayerResiduals[subdetandlayer] =
133  dqmStore_->book1D(histoname.c_str(),histoname.c_str(),
134  i_residuals_Nbins,d_residual_xmin,d_residual_xmax);
135  m_SubdetLayerResiduals[subdetandlayer]->setAxisTitle("(x_{pred} - x_{rec})' [cm]");
136  m_SubdetLayerNormedResiduals[subdetandlayer] =
137  dqmStore_->book1D(normhistoname.c_str(),normhistoname.c_str(),
138  i_normres_Nbins,d_normres_xmin,d_normres_xmax);
139  m_SubdetLayerNormedResiduals[subdetandlayer]->setAxisTitle("(x_{pred} - x_{rec})'/#sigma");
140  }
141  } // end 'is strip module'
142  } // end loop over activeDets
143 }
144 
145 
147  HitResidual[modid]->Reset();
148  NormedHitResiduals[modid]->Reset();
149 }
150 
151 void MonitorTrackResiduals::resetLayerMEs(const std::pair<std::string, int32_t> &subdetandlayer) {
152  m_SubdetLayerResiduals [subdetandlayer]->Reset();
153  m_SubdetLayerNormedResiduals[subdetandlayer]->Reset();
154 }
155 
156 
157 
159 }
160 
161 
163 
164  //dqmStore_->showDirStructure();
165  bool outputMEsInRootFile = conf_.getParameter<bool>("OutputMEsInRootFile");
166  std::string outputFileName = conf_.getParameter<std::string>("OutputFileName");
167  if(outputMEsInRootFile){
168  dqmStore_->save(outputFileName);
169  }
170 }
171 
172 
174 
175  // Filter out events if Trigger Filtering is requested
176  if (genTriggerEventFlag_->on()&& ! genTriggerEventFlag_->accept( iEvent, iSetup) ) return;
177 
178  TrackerValidationVariables avalidator_(iSetup,conf_);
179  std::vector<TrackerValidationVariables::AVHitStruct> v_hitstruct;
180  avalidator_.fillHitQuantities(iEvent,v_hitstruct);
181  for (std::vector<TrackerValidationVariables::AVHitStruct>::const_iterator it = v_hitstruct.begin(),
182  itEnd = v_hitstruct.end(); it != itEnd; ++it) {
183  uint RawId = it->rawDetId;
184 
185  // fill if hit belongs to StripDetector and its error is not zero
186  if( it->resXprimeErr != 0 && SiStripDetId(RawId).subDetector() != 0 ) {
187  if (ModOn && HitResidual[RawId]) {
188  HitResidual[RawId]->Fill(it->resXprime);
189  NormedHitResiduals[RawId]->Fill(it->resXprime/it->resXprimeErr);
190  }
191  std::pair<std::string, int32_t> subdetandlayer = folder_organizer.GetSubDetAndLayer(RawId);
192  if(m_SubdetLayerResiduals[subdetandlayer]) {
193  m_SubdetLayerResiduals[subdetandlayer]->Fill(it->resXprime);
194  m_SubdetLayerNormedResiduals[subdetandlayer]->Fill(it->resXprime/it->resXprimeErr);
195  }
196  }
197  }
198 
199 }
200 
201 
202 
204 
virtual void beginRun(edm::Run const &run, edm::EventSetup const &iSetup)
T getParameter(std::string const &) const
void createMEs(const edm::EventSetup &)
void resetModuleMEs(int32_t modid)
virtual void analyze(const edm::Event &, const edm::EventSetup &)
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:717
unsigned long long m_cacheID_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
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:2113
#define abs(x)
Definition: mlp_lapack.h:159
void fillHitQuantities(const Trajectory *trajectory, std::vector< AVHitStruct > &v_avhitout)
Provides a code based selection for trigger and DCS information in order to have no failing filters i...
SiStripFolderOrganizer folder_organizer
dictionary map
Definition: Association.py:205
int iEvent
Definition: GenABIO.cc:243
vector< ParameterSet > Parameters
bool accept(const edm::Event &event, const edm::EventSetup &setup)
To be called from analyze/filter() methods.
void setDetectorFolder(uint32_t rawdetid=0)
virtual void endRun(const edm::Run &, const edm::EventSetup &)
virtual void beginJob(void)
void resetLayerMEs(const std::pair< std::string, int32_t > &)
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:17
SubDetector subDetector() const
Definition: SiStripDetId.h:114
const T & get() const
Definition: EventSetup.h:55
std::map< std::pair< std::string, int32_t >, MonitorElement * > m_SubdetLayerNormedResiduals
GenericTriggerEventFlag * genTriggerEventFlag_
std::map< std::pair< std::string, int32_t >, MonitorElement * > m_SubdetLayerResiduals
std::pair< std::string, int32_t > GetSubDetAndLayer(const uint32_t &detid, bool ring_flag=0)
MonitorTrackResiduals(const edm::ParameterSet &)
void initRun(const edm::Run &run, const edm::EventSetup &setup)
To be called from beginRun() methods.
void setLayerFolder(uint32_t rawdetid=0, int32_t layer=0, bool ring_flag=0)
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
Definition: Run.h:33