CMS 3D CMS Logo
List of all members | Public Member Functions | Private Member Functions | Private Attributes
DTT0Calibration Class Reference

#include <DTT0Calibration.h>

Inheritance diagram for DTT0Calibration:
edm::one::EDAnalyzer<> edm::one::EDAnalyzerBase edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &event, const edm::EventSetup &eventSetup) override
 Fill the maps with t0 (by channel) More...
 
 DTT0Calibration (const edm::ParameterSet &pset)
 Constructor. More...
 
void endJob () override
 Compute the mean and the RMS of the t0 from the maps and write them to the DB with channel granularity. More...
 
 ~DTT0Calibration () override
 Destructor. More...
 
- Public Member Functions inherited from edm::one::EDAnalyzer<>
 EDAnalyzer ()=default
 
SerialTaskQueue * globalLuminosityBlocksQueue () final
 
SerialTaskQueue * globalRunsQueue () final
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
- Public Member Functions inherited from edm::one::EDAnalyzerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzerBase ()
 
ModuleDescription const & moduleDescription () const
 
bool wantsStreamLuminosityBlocks () const
 
bool wantsStreamRuns () const
 
 ~EDAnalyzerBase () override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo > consumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase &&)=default
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESProxyIndex > const & esGetTokenIndicesVector (edm::Transition iTrans) const
 
std::vector< ESRecordIndex > const & esGetTokenRecordIndicesVector (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase & operator= (EDConsumerBase &&)=default
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

std::string getHistoName (const DTLayerId &lId) const
 
std::string getHistoName (const DTWireId &wId) const
 

Private Attributes

bool debug
 
edm::EDGetTokenT< DTDigiCollection > digiToken
 
edm::ESHandle< DTGeometry > dtGeom
 
unsigned int eventsForLayerT0
 
unsigned int eventsForWireT0
 
TH1D hLayerPeaks
 
std::vector< DTLayerId > layerIdWithWireHistos
 
std::map< DTWireId, double > mK
 
std::map< DTWireId, double > mK_ref
 
std::map< DTWireId, int > nDigiPerWire
 
std::map< DTWireId, int > nDigiPerWire_ref
 
unsigned int nevents
 
std::map< DTWireId, double > qK
 
unsigned int rejectDigiFromPeak
 
int selSector
 
int selWheel
 
TSpectrum spectrum
 
std::map< DTWireId, double > theAbsoluteT0PerWire
 
std::string theCalibSector
 
std::string theCalibWheel
 
std::map< DTChamberId, int > theCountT0ByChamber
 
TFile theFile
 
std::map< DTLayerId, TH1I > theHistoLayerMap
 
std::map< DTWireId, TH1I > theHistoWireMap
 
std::map< DTChamberId, double > theMeanT0ByChamber
 
std::map< DTChamberId, double > theRefT0ByChamber
 
std::map< DTWireId, double > theRelativeT0PerWire
 
std::map< DTChamberId, double > theSigmaT0ByChamber
 
std::map< std::string, double > theSigmaT0LayerMap
 
std::map< DTWireId, double > theSigmaT0PerWire
 
std::map< DTChamberId, double > theSumT0ByChamber
 
std::map< std::string, double > theT0LayerMap
 
std::map< DTLayerId, double > theTPPeakMap
 
double tpPeakWidth
 
double tpPeakWidthPerLayer
 
std::vector< DTWireId > wireIdWithHistos
 

Additional Inherited Members

- Public Types inherited from edm::one::EDAnalyzerBase
typedef EDAnalyzerBase ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::one::EDAnalyzerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< B > consumes (edm::InputTag tag) noexcept
 
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 ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<Transition Tr = Transition::Event>
constexpr auto esConsumes () noexcept
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<Transition Tr = Transition::Event>
auto esConsumes (ESInputTag tag) noexcept
 
template<Transition Tr = Transition::Event>
ESGetTokenGeneric esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
 Used with EventSetupRecord::doGet. More...
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 

Detailed Description

Analyzer class computes the mean and RMS of t0 from pulses. Those values are written in the DB with cell granularity. The mean value for each channel is normalized to a reference time common to all the sector. The t0 of wires in odd layers are corrected for the relative difference between odd and even layers

Author
S. Bolognesi - INFN Torino

Definition at line 36 of file DTT0Calibration.h.

Constructor & Destructor Documentation

◆ DTT0Calibration()

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

Constructor.

Definition at line 34 of file DTT0Calibration.cc.

35  : debug(pset.getUntrackedParameter<bool>("debug")),
36  digiToken(consumes<DTDigiCollection>(pset.getUntrackedParameter<string>("digiLabel"))),
37  theFile(pset.getUntrackedParameter<string>("rootFileName", "DTT0PerLayer.root").c_str(), "RECREATE"),
38  nevents(0),
39  eventsForLayerT0(pset.getParameter<unsigned int>("eventsForLayerT0")),
40  eventsForWireT0(pset.getParameter<unsigned int>("eventsForWireT0")),
41  tpPeakWidth(pset.getParameter<double>("tpPeakWidth")),
42  tpPeakWidthPerLayer(pset.getParameter<double>("tpPeakWidthPerLayer")),
43  rejectDigiFromPeak(pset.getParameter<unsigned int>("rejectDigiFromPeak")),
44  hLayerPeaks("hLayerPeaks", "", 3000, 0, 3000),
45  spectrum(20)
46 
47 {
48  // Get the debug parameter for verbose output
49  if (debug)
50  cout << "[DTT0Calibration]Constructor called!" << endl;
51 
52  theCalibWheel =
53  pset.getUntrackedParameter<string>("calibWheel", "All"); //FIXME amke a vector of integer instead of a string
54  if (theCalibWheel != "All") {
55  stringstream linestr;
56  int selWheel;
57  linestr << theCalibWheel;
58  linestr >> selWheel;
59  cout << "[DTT0CalibrationPerLayer] chosen wheel " << selWheel << endl;
60  }
61 
62  // Sector/s to calibrate
63  theCalibSector =
64  pset.getUntrackedParameter<string>("calibSector", "All"); //FIXME amke a vector of integer instead of a string
65  if (theCalibSector != "All") {
66  stringstream linestr;
67  int selSector;
68  linestr << theCalibSector;
69  linestr >> selSector;
70  cout << "[DTT0CalibrationPerLayer] chosen sector " << selSector << endl;
71  }
72 
73  vector<string> defaultCell;
74  const auto& cellsWithHistos = pset.getUntrackedParameter<vector<string> >("cellsWithHisto", defaultCell);
75  for (const auto& cell : cellsWithHistos) {
76  stringstream linestr;
77  int wheel, sector, station, sl, layer, wire;
78  linestr << cell;
79  linestr >> wheel >> sector >> station >> sl >> layer >> wire;
80  wireIdWithHistos.push_back(DTWireId(wheel, station, sector, sl, layer, wire));
81  }
82 }

References gather_cfg::cout, debug, muonDTDigis_cfi::pset, selSector, selWheel, relativeConstraints::station, theCalibSector, theCalibWheel, makeMuonMisalignmentScenario::wheel, and wireIdWithHistos.

◆ ~DTT0Calibration()

DTT0Calibration::~DTT0Calibration ( )
override

Destructor.

Definition at line 85 of file DTT0Calibration.cc.

85  {
86  if (debug)
87  cout << "[DTT0Calibration]Destructor called!" << endl;
88 
89  theFile.Close();
90 }

References gather_cfg::cout, debug, and theFile.

Member Function Documentation

◆ analyze()

void DTT0Calibration::analyze ( const edm::Event &  event,
const edm::EventSetup &  eventSetup 
)
overridevirtual

Fill the maps with t0 (by channel)

Perform the real analysis.

Implements edm::one::EDAnalyzerBase.

Definition at line 93 of file DTT0Calibration.cc.

93  {
94  nevents++;
95 
96  // Get the digis from the event
97  Handle<DTDigiCollection> digis;
98  event.getByToken(digiToken, digis);
99 
100  // Get the DT Geometry
101  if (nevents == 1)
102  eventSetup.get<MuonGeometryRecord>().get(dtGeom);
103 
104  // Iterate through all digi collections ordered by LayerId
105  for (const auto& digis_per_layer : *digis) {
106  //std::cout << __LINE__ << std::endl;
107  // Get the iterators over the digis associated with this LayerId
108  const DTDigiCollection::Range& digiRange = digis_per_layer.second;
109 
110  // Get the layerId
111  const DTLayerId layerId = digis_per_layer.first;
112  //const DTChamberId chamberId = layerId.superlayerId().chamberId();
113 
114  if ((theCalibWheel != "All") && (layerId.superlayerId().chamberId().wheel() != selWheel))
115  continue;
116  if ((theCalibSector != "All") && (layerId.superlayerId().chamberId().sector() != selSector))
117  continue;
118 
119  // Loop over all digis in the given layer
120  for (DTDigiCollection::const_iterator digi = digiRange.first; digi != digiRange.second; ++digi) {
121  const double t0 = digi->countsTDC();
122  const DTWireId wireIdtmp(layerId, (*digi).wire());
123 
124  // Use first bunch of events to fill t0 per layer
125  if (nevents <= eventsForLayerT0) {
126  // If it doesn't exist, book it
127  if (not theHistoLayerMap.count(layerId)) {
128  theHistoLayerMap[layerId] = TH1I(getHistoName(layerId).c_str(),
129  "T0 from pulses by layer (TDC counts, 1 TDC count = 0.781 ns)",
130  3000,
131  0,
132  3000);
133  if (debug)
134  cout << " New T0 per Layer Histo: " << theHistoLayerMap[layerId].GetName() << endl;
135  }
136  theHistoLayerMap[layerId].Fill(t0);
137  }
138 
139  // Use all the remaining events to compute t0 per wire
140  if (nevents > eventsForLayerT0) {
141  // Get the wireId
142  const DTWireId wireId(layerId, (*digi).wire());
143  if (debug) {
144  cout << " Wire: " << wireId << endl << " time (TDC counts): " << (*digi).countsTDC() << endl;
145  }
146 
147  //Fill the histos per wire for the chosen cells
148  if (std::find(layerIdWithWireHistos.begin(), layerIdWithWireHistos.end(), layerId) !=
149  layerIdWithWireHistos.end() or
150  std::find(wireIdWithHistos.begin(), wireIdWithHistos.end(), wireId) != wireIdWithHistos.end()) {
151  //If it doesn't exist, book it
152  if (theHistoWireMap.count(wireId) == 0) {
153  theHistoWireMap[wireId] = TH1I(getHistoName(wireId).c_str(),
154  "T0 from pulses by wire (TDC counts, 1 TDC count = 0.781 ns)",
155  7000,
156  0,
157  7000);
158  if (debug)
159  cout << " New T0 per wire Histo: " << theHistoWireMap[wireId].GetName() << endl;
160  }
161  theHistoWireMap[wireId].Fill(t0);
162  }
163 
164  //Select per layer
165  if (fabs(theTPPeakMap[layerId] - t0) > rejectDigiFromPeak) {
166  if (debug)
167  cout << "digi skipped because t0 too far from peak " << theTPPeakMap[layerId] << endl;
168  continue;
169  }
170 
171  //Use second bunch of events to compute a t0 reference per wire
172  if (nevents <= (eventsForLayerT0 + eventsForWireT0)) {
173  if (!nDigiPerWire_ref[wireId]) {
174  mK_ref[wireId] = 0;
175  }
176  nDigiPerWire_ref[wireId] = nDigiPerWire_ref[wireId] + 1;
177  mK_ref[wireId] = mK_ref[wireId] + (t0 - mK_ref[wireId]) / nDigiPerWire_ref[wireId];
178  }
179  //Use last all the remaining events to compute the mean and sigma t0 per wire
180  else if (nevents > (eventsForLayerT0 + eventsForWireT0)) {
181  if (abs(t0 - mK_ref[wireId]) > tpPeakWidth)
182  continue;
183  if (!nDigiPerWire[wireId]) {
184  theAbsoluteT0PerWire[wireId] = 0;
185  qK[wireId] = 0;
186  mK[wireId] = 0;
187  }
188  nDigiPerWire[wireId] = nDigiPerWire[wireId] + 1;
189  theAbsoluteT0PerWire[wireId] = theAbsoluteT0PerWire[wireId] + t0;
190  qK[wireId] =
191  qK[wireId] + ((nDigiPerWire[wireId] - 1) * (t0 - mK[wireId]) * (t0 - mK[wireId]) / nDigiPerWire[wireId]);
192  mK[wireId] = mK[wireId] + (t0 - mK[wireId]) / nDigiPerWire[wireId];
193  }
194  } //end if(nevents>1000)
195  } //end loop on digi
196  } //end loop on layer
197 
198  //Use the t0 per layer histos to have an indication about the t0 position
199  if (nevents == eventsForLayerT0) {
200  for (const auto& lHisto : theHistoLayerMap) {
201  const auto& layerId = lHisto.first;
202  const auto& hist = lHisto.second;
203  if (debug)
204  cout << "Reading histogram " << hist.GetName() << " with mean " << hist.GetMean() << " and RMS "
205  << hist.GetRMS() << endl;
206 
207  //Find peaks
208  int npeaks = spectrum.Search(&hist, (tpPeakWidthPerLayer / 2.), "", 0.3);
209 
210  double* peaks = spectrum.GetPositionX();
211  //Put in a std::vector<float>
212  vector<double> peakMeans(peaks, peaks + npeaks);
213  //Sort the peaks in ascending order
214  sort(peakMeans.begin(), peakMeans.end());
215 
216  if (peakMeans.empty()) {
217  theTPPeakMap[layerId] = hist.GetMaximumBin();
218  std::cout << "No peaks found by peakfinder in layer " << layerId << ". Taking maximum bin at "
219  << theTPPeakMap[layerId] << ". Please check!" << std::endl;
220  layerIdWithWireHistos.push_back(layerId);
221  } else if (fabs(hist.GetXaxis()->FindBin(peakMeans.front()) - hist.GetXaxis()->FindBin(peakMeans.back())) <
222  rejectDigiFromPeak) {
223  theTPPeakMap[layerId] = peakMeans[peakMeans.size() / 2];
224  } else {
225  bool peak_set = false;
226  for (const auto& peak : peakMeans) {
227  // Skip if at low edge
228  if (peak - tpPeakWidthPerLayer <= 0)
229  continue;
230  // Get integral of peak
231  double sum = 0;
232  for (int ibin = peak - tpPeakWidthPerLayer; ibin < peak + tpPeakWidthPerLayer; ibin++) {
233  sum += hist.GetBinContent(ibin);
234  }
235  // Skip if peak too small
236  if (sum < hist.GetMaximum() / 2)
237  continue;
238 
239  // Passed all cuts
240  theTPPeakMap[layerId] = peak;
241  peak_set = true;
242  break;
243  }
244  if (peak_set) {
245  std::cout << "Peaks to far away from each other in layer " << layerId
246  << ". Maybe cross talk? Taking first good peak at " << theTPPeakMap[layerId] << ". Please check!"
247  << std::endl;
248  layerIdWithWireHistos.push_back(layerId);
249  } else {
250  theTPPeakMap[layerId] = hist.GetMaximumBin();
251  std::cout << "Peaks to far away from each other in layer " << layerId
252  << " and no good peak found. Taking maximum bin at " << theTPPeakMap[layerId] << ". Please check!"
253  << std::endl;
254  layerIdWithWireHistos.push_back(layerId);
255  }
256  }
257  if (peakMeans.size() > 5) {
258  std::cout << "Found more than 5 peaks in layer " << layerId << ". Please check!" << std::endl;
259  if (std::find(layerIdWithWireHistos.begin(), layerIdWithWireHistos.end(), layerId) ==
260  layerIdWithWireHistos.end())
261  layerIdWithWireHistos.push_back(layerId);
262  }
263  // Check for noise
264  int nspikes = 0;
265  for (int ibin = 0; ibin < hist.GetNbinsX(); ibin++) {
266  if (hist.GetBinContent(ibin + 1) > hist.GetMaximum() * 0.001)
267  nspikes++;
268  }
269  if (nspikes > 50) {
270  std::cout << "Found a lot of (>50) small spikes in layer " << layerId
271  << ". Please check if all wires are functioning as expected!" << std::endl;
272  if (std::find(layerIdWithWireHistos.begin(), layerIdWithWireHistos.end(), layerId) ==
273  layerIdWithWireHistos.end())
274  layerIdWithWireHistos.push_back(layerId);
275  }
276  hLayerPeaks.Fill(theTPPeakMap[layerId]);
277  }
278  }
279 }

References funct::abs(), DTSuperLayerId::chamberId(), gather_cfg::cout, debug, digiToken, dtGeom, eventsForLayerT0, eventsForWireT0, spr::find(), edm::EventSetup::get(), get, getHistoName(), compare::hist, hLayerPeaks, layerIdWithWireHistos, mK, mK_ref, nDigiPerWire, nDigiPerWire_ref, nevents, or, qK, rejectDigiFromPeak, DTChamberId::sector(), selSector, selWheel, spectrum, DTLayerId::superlayerId(), FrontierCondition_GT_autoExpress_cfi::t0, theAbsoluteT0PerWire, theCalibSector, theCalibWheel, theHistoLayerMap, theHistoWireMap, theTPPeakMap, tpPeakWidth, tpPeakWidthPerLayer, DTChamberId::wheel(), and wireIdWithHistos.

◆ endJob()

void DTT0Calibration::endJob ( void  )
overridevirtual

Compute the mean and the RMS of the t0 from the maps and write them to the DB with channel granularity.

Write the t0 map into DB

Reimplemented from edm::one::EDAnalyzerBase.

Definition at line 281 of file DTT0Calibration.cc.

281  {
282  std::cout << "Analyzed " << nevents << " events" << std::endl;
283 
284  DTT0* t0sWRTChamber = new DTT0();
285 
286  if (debug)
287  cout << "[DTT0CalibrationPerLayer]Writing histos to file!" << endl;
288 
289  theFile.cd();
290  //hT0SectorHisto->Write();
291  hLayerPeaks.Write();
292  for (const auto& wHisto : theHistoWireMap) {
293  wHisto.second.Write();
294  }
295  for (const auto& lHisto : theHistoLayerMap) {
296  lHisto.second.Write();
297  }
298 
299  if (debug)
300  cout << "[DTT0Calibration] Compute and store t0 and sigma per wire" << endl;
301 
302  // Calculate uncertainties per wire (counting experiment)
303  for (auto& wiret0 : theAbsoluteT0PerWire) {
304  auto& wireId = wiret0.first;
305  if (nDigiPerWire[wireId] > 1)
306  theSigmaT0PerWire[wireId] = qK[wireId] / (nDigiPerWire[wireId] - 1);
307  else
308  theSigmaT0PerWire[wireId] = 999.; // Only one measurement: uncertainty -> infinity
309  // syst uncert
310  //theSigmaT0PerWire[wireId] += pow(0.5, 2));
311  // Every time the same measurement. Use Laplace estimator as estimation how propable it is to measure another value due to limited size of sample
312  if (theSigmaT0PerWire[wireId] == 0) {
313  theSigmaT0PerWire[wireId] += pow(1. / (nDigiPerWire[wireId] + 1), 2);
314  }
315  }
316 
317  // function to calculate unweighted means
318  auto unweighted_mean_function = [](const std::list<double>& values, const std::list<double>& sigmas) {
319  double mean = 0;
320  for (auto& value : values) {
321  mean += value;
322  }
323  mean /= values.size();
324 
325  double uncertainty = 0;
326  for (auto& value : values) {
327  uncertainty += pow(value - mean, 2);
328  }
329  uncertainty /= values.size();
330  uncertainty = sqrt(uncertainty);
331  return std::make_pair(mean, uncertainty);
332  };
333 
334  // correct for odd-even effect in each super layer
335  std::map<DTSuperLayerId, std::pair<double, double> > mean_sigma_even;
336  std::map<DTSuperLayerId, std::pair<double, double> > mean_sigma_odd;
337  for (const auto& superlayer : dtGeom->superLayers()) {
338  const auto superlayer_id = superlayer->id();
339  std::list<double> values_even;
340  std::list<double> sigmas_even;
341  std::list<double> values_odd;
342  std::list<double> sigmas_odd;
343 
344  for (const auto& wiret0 : theAbsoluteT0PerWire) {
345  const auto& wireId = wiret0.first;
346  if (wireId.layerId().superlayerId() == superlayer_id) {
347  const auto& t0 = wiret0.second / nDigiPerWire[wireId];
348  if (wireId.layerId().layer() % 2) {
349  values_odd.push_back(t0);
350  sigmas_odd.push_back(sqrt(theSigmaT0PerWire[wireId]));
351  } else {
352  values_even.push_back(t0);
353  sigmas_even.push_back(sqrt(theSigmaT0PerWire[wireId]));
354  }
355  }
356  }
357  // get mean and uncertainty
358  mean_sigma_even.emplace(superlayer_id, unweighted_mean_function(values_even, sigmas_even));
359  mean_sigma_odd.emplace(superlayer_id, unweighted_mean_function(values_odd, sigmas_odd));
360  }
361 
362  // filter outliers
363  for (const auto& superlayer : dtGeom->superLayers()) {
364  const auto superlayer_id = superlayer->id();
365  std::list<double> values_even;
366  std::list<double> sigmas_even;
367  std::list<double> values_odd;
368  std::list<double> sigmas_odd;
369 
370  for (const auto& wiret0 : theAbsoluteT0PerWire) {
371  const auto& wireId = wiret0.first;
372  if (wireId.layerId().superlayerId() == superlayer_id) {
373  const auto& t0 = wiret0.second / nDigiPerWire[wireId];
374  if (wireId.layerId().layer() % 2 and
375  abs(t0 - mean_sigma_odd[superlayer_id].first) < 2 * mean_sigma_odd[superlayer_id].second) {
376  values_odd.push_back(t0);
377  sigmas_odd.push_back(sqrt(theSigmaT0PerWire[wireId]));
378  } else {
379  if (abs(t0 - mean_sigma_even[superlayer_id].first) < 2 * mean_sigma_even[superlayer_id].second) {
380  values_even.push_back(t0);
381  sigmas_even.push_back(sqrt(theSigmaT0PerWire[wireId]));
382  }
383  }
384  }
385  }
386  // get mean and uncertainty
387  mean_sigma_even[superlayer_id] = unweighted_mean_function(values_even, sigmas_even);
388  mean_sigma_odd[superlayer_id] = unweighted_mean_function(values_odd, sigmas_odd);
389  }
390 
391  // apply correction
392  for (auto& wiret0 : theAbsoluteT0PerWire) {
393  const auto& wire_id = wiret0.first;
394  const auto& superlayer_id = wiret0.first.layerId().superlayerId();
395  const auto& layer = wiret0.first.layerId().layer();
396  auto& t0 = wiret0.second;
397  t0 /= nDigiPerWire[wire_id];
398  if (not layer % 2)
399  continue;
400  // t0 is reference. Changing it changes the map
401  t0 += mean_sigma_even[superlayer_id].first - mean_sigma_odd[superlayer_id].first;
402  theSigmaT0PerWire[wire_id] +=
403  pow(mean_sigma_odd[superlayer_id].second, 2) + pow(mean_sigma_even[superlayer_id].second, 2);
404  }
405 
406  // get chamber mean
407  std::map<DTChamberId, std::list<double> > values_per_chamber;
408  std::map<DTChamberId, std::list<double> > sigmas_per_chamber;
409  for (const auto& wire_t0 : theAbsoluteT0PerWire) {
410  const auto& wire_id = wire_t0.first;
411  const auto& chamber_id = wire_id.chamberId();
412  const auto& t0 = wire_t0.second;
413  values_per_chamber[chamber_id].push_back(t0);
414  sigmas_per_chamber[chamber_id].push_back(sqrt(theSigmaT0PerWire[wire_id]));
415  }
416 
417  std::map<DTChamberId, std::pair<double, double> > mean_per_chamber;
418  for (const auto& chamber_mean : values_per_chamber) {
419  const auto& chamber_id = chamber_mean.first;
420  const auto& means = chamber_mean.second;
421  const auto& sigmas = sigmas_per_chamber[chamber_id];
422  mean_per_chamber.emplace(chamber_id, unweighted_mean_function(means, sigmas));
423  }
424 
425  // calculate relative values
426  for (const auto& wire_t0 : theAbsoluteT0PerWire) {
427  const auto& wire_id = wire_t0.first;
428  const auto& chamber_id = wire_id.chamberId();
429  const auto& t0 = wire_t0.second;
430  theRelativeT0PerWire.emplace(wire_id, t0 - mean_per_chamber[chamber_id].first);
431  cout << "[DTT0Calibration] Wire " << wire_id << " has t0 " << theRelativeT0PerWire[wire_id]
432  << " (relative, after even-odd layer corrections) "
433  << " sigma " << sqrt(theSigmaT0PerWire[wire_id]) << endl;
434  }
435 
436  for (const auto& wire_t0 : theRelativeT0PerWire) {
437  const auto& wire_id = wire_t0.first;
438  const auto& t0 = wire_t0.second;
439  t0sWRTChamber->set(wire_id, t0, sqrt(theSigmaT0PerWire[wire_id]), DTTimeUnits::counts);
440  }
441 
443  if (debug)
444  cout << "[DTT0Calibration]Writing values in DB!" << endl;
445  // FIXME: to be read from cfg?
446  string t0Record = "DTT0Rcd";
447  // Write the t0 map to DB
448  DTCalibDBUtils::writeToDB(t0Record, t0sWRTChamber);
449  delete t0sWRTChamber;
450 }

References funct::abs(), DTTimeUnits::counts, gather_cfg::cout, debug, dtGeom, dqmdumpme::first, hLayerPeaks, SiStripPI::mean, nDigiPerWire, nevents, funct::pow(), qK, edm::second(), DTT0::set(), mathSSE::sqrt(), DTGeometry::superLayers(), FrontierCondition_GT_autoExpress_cfi::t0, theAbsoluteT0PerWire, theFile, theHistoLayerMap, theHistoWireMap, theRelativeT0PerWire, theSigmaT0PerWire, relativeConstraints::value, contentValuesCheck::values, and DTCalibDBUtils::writeToDB().

Referenced by o2olib.O2ORunMgr::executeJob().

◆ getHistoName() [1/2]

string DTT0Calibration::getHistoName ( const DTLayerId &  lId) const
private

Definition at line 461 of file DTT0Calibration.cc.

461  {
462  string histoName;
463  stringstream theStream;
464  theStream << "Ch_" << lId.wheel() << "_" << lId.station() << "_" << lId.sector() << "_SL" << lId.superlayer() << "_L"
465  << lId.layer() << "_hT0Histo";
466  theStream >> histoName;
467  return histoName;
468 }

References HltBtagPostValidation_cff::histoName, DTLayerId::layer(), DTChamberId::sector(), DTChamberId::station(), DTSuperLayerId::superlayer(), and DTChamberId::wheel().

◆ getHistoName() [2/2]

string DTT0Calibration::getHistoName ( const DTWireId &  wId) const
private

Definition at line 452 of file DTT0Calibration.cc.

452  {
453  string histoName;
454  stringstream theStream;
455  theStream << "Ch_" << wId.wheel() << "_" << wId.station() << "_" << wId.sector() << "_SL" << wId.superlayer() << "_L"
456  << wId.layer() << "_W" << wId.wire() << "_hT0Histo";
457  theStream >> histoName;
458  return histoName;
459 }

References HltBtagPostValidation_cff::histoName, DTLayerId::layer(), DTChamberId::sector(), DTChamberId::station(), DTSuperLayerId::superlayer(), DTChamberId::wheel(), and DTWireId::wire().

Referenced by analyze().

Member Data Documentation

◆ debug

bool DTT0Calibration::debug
private

Definition at line 59 of file DTT0Calibration.h.

Referenced by analyze(), rrapi.RRApi::dprint(), util.rrapi.RRApi::dprint(), DTT0Calibration(), endJob(), pkg.AbstractPkg::generate(), rrapi.RRApi::get(), util.rrapi.RRApi::get(), pkg.AbstractPkg::get_kwds(), runTauIdMVA.TauIDEmbedder::loadMVA_WPs_run2_2017(), runTauIdMVA.TauIDEmbedder::runTauID(), crabFunctions.CrabController::submit(), pkg.AbstractPkg::write(), and ~DTT0Calibration().

◆ digiToken

edm::EDGetTokenT<DTDigiCollection> DTT0Calibration::digiToken
private

Definition at line 62 of file DTT0Calibration.h.

Referenced by analyze().

◆ dtGeom

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

Definition at line 124 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ eventsForLayerT0

unsigned int DTT0Calibration::eventsForLayerT0
private

Definition at line 70 of file DTT0Calibration.h.

Referenced by analyze().

◆ eventsForWireT0

unsigned int DTT0Calibration::eventsForWireT0
private

Definition at line 72 of file DTT0Calibration.h.

Referenced by analyze().

◆ hLayerPeaks

TH1D DTT0Calibration::hLayerPeaks
private

Definition at line 93 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ layerIdWithWireHistos

std::vector<DTLayerId> DTT0Calibration::layerIdWithWireHistos
private

Definition at line 99 of file DTT0Calibration.h.

Referenced by analyze().

◆ mK

std::map<DTWireId, double> DTT0Calibration::mK
private

Definition at line 107 of file DTT0Calibration.h.

Referenced by analyze().

◆ mK_ref

std::map<DTWireId, double> DTT0Calibration::mK_ref
private

Definition at line 108 of file DTT0Calibration.h.

Referenced by analyze().

◆ nDigiPerWire

std::map<DTWireId, int> DTT0Calibration::nDigiPerWire
private

Definition at line 105 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ nDigiPerWire_ref

std::map<DTWireId, int> DTT0Calibration::nDigiPerWire_ref
private

Definition at line 106 of file DTT0Calibration.h.

Referenced by analyze().

◆ nevents

unsigned int DTT0Calibration::nevents
private

Definition at line 68 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ qK

std::map<DTWireId, double> DTT0Calibration::qK
private

Definition at line 109 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ rejectDigiFromPeak

unsigned int DTT0Calibration::rejectDigiFromPeak
private

Definition at line 81 of file DTT0Calibration.h.

Referenced by analyze().

◆ selSector

int DTT0Calibration::selSector
private

Definition at line 87 of file DTT0Calibration.h.

Referenced by analyze(), and DTT0Calibration().

◆ selWheel

int DTT0Calibration::selWheel
private

Definition at line 85 of file DTT0Calibration.h.

Referenced by analyze(), and DTT0Calibration().

◆ spectrum

TSpectrum DTT0Calibration::spectrum
private

Definition at line 95 of file DTT0Calibration.h.

Referenced by analyze().

◆ theAbsoluteT0PerWire

std::map<DTWireId, double> DTT0Calibration::theAbsoluteT0PerWire
private

Definition at line 102 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ theCalibSector

std::string DTT0Calibration::theCalibSector
private

Definition at line 86 of file DTT0Calibration.h.

Referenced by analyze(), and DTT0Calibration().

◆ theCalibWheel

std::string DTT0Calibration::theCalibWheel
private

Definition at line 84 of file DTT0Calibration.h.

Referenced by analyze(), and DTT0Calibration().

◆ theCountT0ByChamber

std::map<DTChamberId, int> DTT0Calibration::theCountT0ByChamber
private

Definition at line 118 of file DTT0Calibration.h.

◆ theFile

TFile DTT0Calibration::theFile
private

Definition at line 65 of file DTT0Calibration.h.

Referenced by endJob(), and ~DTT0Calibration().

◆ theHistoLayerMap

std::map<DTLayerId, TH1I> DTT0Calibration::theHistoLayerMap
private

Definition at line 90 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ theHistoWireMap

std::map<DTWireId, TH1I> DTT0Calibration::theHistoWireMap
private

Definition at line 111 of file DTT0Calibration.h.

Referenced by analyze(), and endJob().

◆ theMeanT0ByChamber

std::map<DTChamberId, double> DTT0Calibration::theMeanT0ByChamber
private

Definition at line 120 of file DTT0Calibration.h.

◆ theRefT0ByChamber

std::map<DTChamberId, double> DTT0Calibration::theRefT0ByChamber
private

Definition at line 121 of file DTT0Calibration.h.

◆ theRelativeT0PerWire

std::map<DTWireId, double> DTT0Calibration::theRelativeT0PerWire
private

Definition at line 103 of file DTT0Calibration.h.

Referenced by endJob().

◆ theSigmaT0ByChamber

std::map<DTChamberId, double> DTT0Calibration::theSigmaT0ByChamber
private

Definition at line 119 of file DTT0Calibration.h.

◆ theSigmaT0LayerMap

std::map<std::string, double> DTT0Calibration::theSigmaT0LayerMap
private

Definition at line 114 of file DTT0Calibration.h.

◆ theSigmaT0PerWire

std::map<DTWireId, double> DTT0Calibration::theSigmaT0PerWire
private

Definition at line 104 of file DTT0Calibration.h.

Referenced by endJob().

◆ theSumT0ByChamber

std::map<DTChamberId, double> DTT0Calibration::theSumT0ByChamber
private

Definition at line 117 of file DTT0Calibration.h.

◆ theT0LayerMap

std::map<std::string, double> DTT0Calibration::theT0LayerMap
private

Definition at line 113 of file DTT0Calibration.h.

◆ theTPPeakMap

std::map<DTLayerId, double> DTT0Calibration::theTPPeakMap
private

Definition at line 115 of file DTT0Calibration.h.

Referenced by analyze().

◆ tpPeakWidth

double DTT0Calibration::tpPeakWidth
private

Definition at line 75 of file DTT0Calibration.h.

Referenced by analyze().

◆ tpPeakWidthPerLayer

double DTT0Calibration::tpPeakWidthPerLayer
private

Definition at line 78 of file DTT0Calibration.h.

Referenced by analyze().

◆ wireIdWithHistos

std::vector<DTWireId> DTT0Calibration::wireIdWithHistos
private

Definition at line 98 of file DTT0Calibration.h.

Referenced by analyze(), and DTT0Calibration().

DTT0Calibration::theFile
TFile theFile
Definition: DTT0Calibration.h:65
DTT0Calibration::layerIdWithWireHistos
std::vector< DTLayerId > layerIdWithWireHistos
Definition: DTT0Calibration.h:99
DTWireId::wire
int wire() const
Return the wire number.
Definition: DTWireId.h:42
SiStripPI::mean
Definition: SiStripPayloadInspectorHelper.h:169
relativeConstraints.station
station
Definition: relativeConstraints.py:67
gather_cfg.cout
cout
Definition: gather_cfg.py:144
DTT0Calibration::theSigmaT0PerWire
std::map< DTWireId, double > theSigmaT0PerWire
Definition: DTT0Calibration.h:104
DTT0Calibration::theCalibWheel
std::string theCalibWheel
Definition: DTT0Calibration.h:84
DTT0Calibration::mK
std::map< DTWireId, double > mK
Definition: DTT0Calibration.h:107
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:222
DTSuperLayerId::superlayer
int superlayer() const
Return the superlayer number (deprecated method name)
Definition: DTSuperLayerId.h:42
DTT0Calibration::tpPeakWidthPerLayer
double tpPeakWidthPerLayer
Definition: DTT0Calibration.h:78
DTT0Calibration::theRelativeT0PerWire
std::map< DTWireId, double > theRelativeT0PerWire
Definition: DTT0Calibration.h:103
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
edm::Handle< DTDigiCollection >
dqmdumpme.first
first
Definition: dqmdumpme.py:55
MuonDigiCollection::const_iterator
std::vector< DigiType >::const_iterator const_iterator
Definition: MuonDigiCollection.h:94
edm::EventSetup::get
T get() const
Definition: EventSetup.h:80
DTT0Calibration::theHistoLayerMap
std::map< DTLayerId, TH1I > theHistoLayerMap
Definition: DTT0Calibration.h:90
compare.hist
hist
Definition: compare.py:376
DTWireId
Definition: DTWireId.h:12
FrontierCondition_GT_autoExpress_cfi.t0
t0
Definition: FrontierCondition_GT_autoExpress_cfi.py:148
DTT0Calibration::theAbsoluteT0PerWire
std::map< DTWireId, double > theAbsoluteT0PerWire
Definition: DTT0Calibration.h:102
DTT0Calibration::theHistoWireMap
std::map< DTWireId, TH1I > theHistoWireMap
Definition: DTT0Calibration.h:111
contentValuesCheck.values
values
Definition: contentValuesCheck.py:38
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
DTT0Calibration::nDigiPerWire
std::map< DTWireId, int > nDigiPerWire
Definition: DTT0Calibration.h:105
DTT0::set
int set(int wheelId, int stationId, int sectorId, int slId, int layerId, int cellId, float t0mean, float t0rms, DTTimeUnits::type unit)
Definition: DTT0.cc:97
DTT0Calibration::eventsForWireT0
unsigned int eventsForWireT0
Definition: DTT0Calibration.h:72
DTT0Calibration::getHistoName
std::string getHistoName(const DTWireId &wId) const
Definition: DTT0Calibration.cc:452
DTLayerId
Definition: DTLayerId.h:12
DTT0Calibration::selSector
int selSector
Definition: DTT0Calibration.h:87
DTT0Calibration::qK
std::map< DTWireId, double > qK
Definition: DTT0Calibration.h:109
DTT0Calibration::dtGeom
edm::ESHandle< DTGeometry > dtGeom
Definition: DTT0Calibration.h:124
DTT0Calibration::debug
bool debug
Definition: DTT0Calibration.h:59
DTT0Calibration::eventsForLayerT0
unsigned int eventsForLayerT0
Definition: DTT0Calibration.h:70
DTT0Calibration::nevents
unsigned int nevents
Definition: DTT0Calibration.h:68
DTT0Calibration::spectrum
TSpectrum spectrum
Definition: DTT0Calibration.h:95
DTSuperLayerId::chamberId
DTChamberId chamberId() const
Return the corresponding ChamberId.
Definition: DTSuperLayerId.h:45
DTT0Calibration::digiToken
edm::EDGetTokenT< DTDigiCollection > digiToken
Definition: DTT0Calibration.h:62
makeMuonMisalignmentScenario.wheel
wheel
Definition: makeMuonMisalignmentScenario.py:319
value
Definition: value.py:1
DTGeometry::superLayers
const std::vector< const DTSuperLayer * > & superLayers() const
Return a vector of all SuperLayer.
Definition: DTGeometry.cc:86
DTT0Calibration::hLayerPeaks
TH1D hLayerPeaks
Definition: DTT0Calibration.h:93
get
#define get
DTT0Calibration::tpPeakWidth
double tpPeakWidth
Definition: DTT0Calibration.h:75
DTChamberId::sector
int sector() const
Definition: DTChamberId.h:49
DTTimeUnits::counts
Definition: DTTimeUnits.h:32
relativeConstraints.value
value
Definition: relativeConstraints.py:53
DTLayerId::superlayerId
DTSuperLayerId superlayerId() const
Return the corresponding SuperLayerId.
Definition: DTLayerId.h:45
DTT0
Definition: DTT0.h:48
DTT0Calibration::mK_ref
std::map< DTWireId, double > mK_ref
Definition: DTT0Calibration.h:108
HltBtagPostValidation_cff.histoName
histoName
Definition: HltBtagPostValidation_cff.py:17
MuonDigiCollection::Range
std::pair< const_iterator, const_iterator > Range
Definition: MuonDigiCollection.h:95
or
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
DTT0Calibration::theTPPeakMap
std::map< DTLayerId, double > theTPPeakMap
Definition: DTT0Calibration.h:115
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
DTT0Calibration::selWheel
int selWheel
Definition: DTT0Calibration.h:85
DTT0Calibration::nDigiPerWire_ref
std::map< DTWireId, int > nDigiPerWire_ref
Definition: DTT0Calibration.h:106
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
DTT0Calibration::theCalibSector
std::string theCalibSector
Definition: DTT0Calibration.h:86
DTLayerId::layer
int layer() const
Return the layer number.
Definition: DTLayerId.h:42
MuonGeometryRecord
Definition: MuonGeometryRecord.h:34
DTT0Calibration::wireIdWithHistos
std::vector< DTWireId > wireIdWithHistos
Definition: DTT0Calibration.h:98
DTChamberId::wheel
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:39
DTCalibDBUtils::writeToDB
static void writeToDB(std::string record, T *payload)
Definition: DTCalibDBUtils.h:28
DTChamberId::station
int station() const
Return the station number.
Definition: DTChamberId.h:42
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
DTT0Calibration::rejectDigiFromPeak
unsigned int rejectDigiFromPeak
Definition: DTT0Calibration.h:81
Generated for CMSSW Reference Manual by doxygen 1.8.16