CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes
DTT0CalibrationRMS Class Reference

#include <DTT0CalibrationRMS.h>

Inheritance diagram for DTT0CalibrationRMS:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

void analyze (const edm::Event &event, const edm::EventSetup &eventSetup) override
 Fill the maps with t0 (by channel) More...
 
 DTT0CalibrationRMS (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...
 
 ~DTT0CalibrationRMS () override
 Destructor. More...
 
- Public Member Functions inherited from edm::EDAnalyzer
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzer ()
 
SerialTaskQueueglobalLuminosityBlocksQueue ()
 
SerialTaskQueueglobalRunsQueue ()
 
ModuleDescription const & moduleDescription () const
 
std::string workerType () const
 
 ~EDAnalyzer () override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
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 const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

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

Private Attributes

std::vector< std::string > cellsWithHistos
 
bool correctByChamberMean_
 
bool debug
 
std::string digiLabel
 
edm::ESHandle< DTGeometrydtGeom
 
unsigned int eventsForLayerT0
 
unsigned int eventsForWireT0
 
TH1D * hT0SectorHisto
 
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
 
std::map< DTWireId, double > theAbsoluteT0PerWire
 
std::string theCalibSector
 
std::string theCalibWheel
 
TFile * theFile
 
std::map< DTLayerId, TH1I * > theHistoLayerMap
 
std::map< DTWireId, TH1I * > theHistoWireMap
 
std::map< DTWireId, TH1I * > theHistoWireMap_ref
 
TFile * theOutputFile
 
std::map< DTWireId, double > theRelativeT0PerWire
 
std::map< std::string, double > theSigmaT0LayerMap
 
std::map< DTWireId, double > theSigmaT0PerWire
 
std::map< std::string, double > theT0LayerMap
 
double tpPeakWidth
 
std::vector< DTWireIdwireIdWithHistos
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 
static bool wantsGlobalLuminosityBlocks ()
 
static bool wantsGlobalRuns ()
 
static bool wantsStreamLuminosityBlocks ()
 
static bool wantsStreamRuns ()
 
- 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

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

Definition at line 26 of file DTT0CalibrationRMS.h.

Constructor & Destructor Documentation

DTT0CalibrationRMS::DTT0CalibrationRMS ( const edm::ParameterSet pset)

Constructor.

Definition at line 27 of file DTT0CalibrationRMS.cc.

References gather_cfg::cout, debug, CastorSimpleReconstructor_cfi::digiLabel, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), DTAnalyzerDetailed_cfi::rootFileName, relativeConstraints::station, interactiveExample::theFile, makeMuonMisalignmentScenario::wheel, and mixOne_premix_on_sim_cfi::wire.

27  {
28  // Get the debug parameter for verbose output
29  debug = pset.getUntrackedParameter<bool>("debug");
30  if(debug)
31  cout << "[DTT0CalibrationRMS]Constructor called!" << endl;
32 
33  // Get the label to retrieve digis from the event
34  digiLabel = pset.getUntrackedParameter<string>("digiLabel");
35 
36  // The root file which contain the histos per layer
37  string rootFileName = pset.getUntrackedParameter<string>("rootFileName","DTT0PerLayer.root");
38  theFile = new TFile(rootFileName.c_str(), "RECREATE");
39 
40  theCalibWheel = pset.getUntrackedParameter<string>("calibWheel", "All"); //FIXME amke a vector of integer instead of a string
41  if(theCalibWheel != "All") {
42  stringstream linestr;
43  int selWheel;
44  linestr << theCalibWheel;
45  linestr >> selWheel;
46  cout << "[DTT0CalibrationRMSPerLayer] chosen wheel " << selWheel << endl;
47  }
48 
49  // Sector/s to calibrate
50  theCalibSector = pset.getUntrackedParameter<string>("calibSector", "All"); //FIXME amke a vector of integer instead of a string
51  if(theCalibSector != "All") {
52  stringstream linestr;
53  int selSector;
54  linestr << theCalibSector;
55  linestr >> selSector;
56  cout << "[DTT0CalibrationRMSPerLayer] chosen sector " << selSector << endl;
57  }
58 
59  vector<string> defaultCell;
60  defaultCell.push_back("None");
61  cellsWithHistos = pset.getUntrackedParameter<vector<string> >("cellsWithHisto", defaultCell);
62  for(vector<string>::const_iterator cell = cellsWithHistos.begin(); cell != cellsWithHistos.end(); ++cell){
63  if((*cell) != "None"){
64  stringstream linestr;
65  int wheel,sector,station,sl,layer,wire;
66  linestr << (*cell);
67  linestr >> wheel >> sector >> station >> sl >> layer >> wire;
68  wireIdWithHistos.push_back(DTWireId(wheel,station,sector,sl,layer,wire));
69  }
70  }
71 
72  hT0SectorHisto=nullptr;
73 
74  nevents=0;
75  eventsForLayerT0 = pset.getParameter<unsigned int>("eventsForLayerT0");
76  eventsForWireT0 = pset.getParameter<unsigned int>("eventsForWireT0");
77  rejectDigiFromPeak = pset.getParameter<unsigned int>("rejectDigiFromPeak");
78  tpPeakWidth = pset.getParameter<double>("tpPeakWidth");
79  //useReferenceWireInLayer_ = true;
80  correctByChamberMean_ = pset.getParameter<bool>("correctByChamberMean");
81 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
unsigned int rejectDigiFromPeak
std::vector< DTWireId > wireIdWithHistos
std::string theCalibSector
unsigned int eventsForWireT0
unsigned int eventsForLayerT0
std::vector< std::string > cellsWithHistos
DTT0CalibrationRMS::~DTT0CalibrationRMS ( )
override

Destructor.

Definition at line 84 of file DTT0CalibrationRMS.cc.

References gather_cfg::cout, debug, and interactiveExample::theFile.

84  {
85  if(debug)
86  cout << "[DTT0CalibrationRMS]Destructor called!" << endl;
87 
88  theFile->Close();
89 }

Member Function Documentation

void DTT0CalibrationRMS::analyze ( const edm::Event event,
const edm::EventSetup eventSetup 
)
override

Fill the maps with t0 (by channel)

Perform the real analysis.

Definition at line 92 of file DTT0CalibrationRMS.cc.

References funct::abs(), DTSuperLayerId::chamberId(), gather_cfg::cout, debug, CastorSimpleReconstructor_cfi::digiLabel, edm::EventID::event(), spr::find(), edm::EventSetup::get(), mergeVDriftHistosByStation::getHistoName(), edm::EventBase::id(), edm::EventID::run(), DTChamberId::sector(), DTLayerId::superlayerId(), cscNeutronWriter_cfi::t0, interactiveExample::theFile, and DTChamberId::wheel().

92  {
93  if(debug || event.id().event() % 500==0)
94  cout << "--- [DTT0CalibrationRMS] Analysing Event: #Run: " << event.id().run()
95  << " #Event: " << event.id().event() << endl;
96  nevents++;
97 
98  // Get the digis from the event
100  event.getByLabel(digiLabel, digis);
101 
102  // Get the DT Geometry
103  eventSetup.get<MuonGeometryRecord>().get(dtGeom);
104 
105  // Iterate through all digi collections ordered by LayerId
107  for (dtLayerIt = digis->begin();
108  dtLayerIt != digis->end();
109  ++dtLayerIt){
110  // Get the iterators over the digis associated with this LayerId
111  const DTDigiCollection::Range& digiRange = (*dtLayerIt).second;
112 
113  // Get the layerId
114  const DTLayerId layerId = (*dtLayerIt).first; //FIXME: check to be in the right sector
115 
116  if((theCalibWheel != "All") && (layerId.superlayerId().chamberId().wheel() != selWheel))
117  continue;
118  if((theCalibSector != "All") && (layerId.superlayerId().chamberId().sector() != selSector))
119  continue;
120 
121  //if(debug) {
122  // cout << "Layer " << layerId<<" with "<<distance(digiRange.first, digiRange.second)<<" digi"<<endl;
123  //}
124 
125  // Loop over all digis in the given layer
126  for (DTDigiCollection::const_iterator digi = digiRange.first;
127  digi != digiRange.second;
128  ++digi) {
129  double t0 = (*digi).countsTDC();
130 
131  //Use first bunch of events to fill t0 per layer
133  //Get the per-layer histo from the map
134  TH1I *hT0LayerHisto = theHistoLayerMap[layerId];
135  //If it doesn't exist, book it
136  if(hT0LayerHisto == nullptr){
137  theFile->cd();
138  hT0LayerHisto = new TH1I(getHistoName(layerId).c_str(),
139  "T0 from pulses by layer (TDC counts, 1 TDC count = 0.781 ns)",
140  200, t0-100, t0+100);
141  if(debug)
142  cout << " New T0 per Layer Histo: " << hT0LayerHisto->GetName() << endl;
143  theHistoLayerMap[layerId] = hT0LayerHisto;
144  }
145 
146  //Fill the histos
147  theFile->cd();
148  if(hT0LayerHisto != nullptr) {
149  // if(debug)
150  // cout<<"Filling histo "<<hT0LayerHisto->GetName()<<" with digi "<<t0<<" TDC counts"<<endl;
151  hT0LayerHisto->Fill(t0);
152  }
153  }
154 
155  //Use all the remaining events to compute t0 per wire
157  // Get the wireId
158  const DTWireId wireId(layerId, (*digi).wire());
159  if(debug) {
160  cout << " Wire: " << wireId << endl
161  << " time (TDC counts): " << (*digi).countsTDC()<< endl;
162  }
163 
164  //Fill the histos per wire for the chosen cells
165  vector<DTWireId>::iterator it_wire = find(wireIdWithHistos.begin(),wireIdWithHistos.end(),wireId);
166  if(it_wire != wireIdWithHistos.end()){
167  if(theHistoWireMap.find(wireId) == theHistoWireMap.end()){
168  theHistoWireMap[wireId] = new TH1I(getHistoName(wireId).c_str(),"T0 from pulses by wire (TDC counts, 1 TDC count = 0.781 ns)",7000,0,7000);
169  if(debug) cout << " New T0 per wire Histo: " << (theHistoWireMap[wireId])->GetName() << endl;
170  }
171  if(theHistoWireMap_ref.find(wireId) == theHistoWireMap_ref.end()){
172  theHistoWireMap_ref[wireId] = new TH1I((getHistoName(wireId) + "_ref").c_str(),"T0 from pulses by wire (TDC counts, 1 TDC count = 0.781 ns)",7000,0,7000);
173  if(debug) cout << " New T0 per wire Histo: " << (theHistoWireMap_ref[wireId])->GetName() << endl;
174  }
175 
176  TH1I* hT0WireHisto = theHistoWireMap[wireId];
177  //Fill the histos
178  theFile->cd();
179  if(hT0WireHisto) hT0WireHisto->Fill(t0);
180  }
181 
182  //Check the tzero has reasonable value
183  if(abs(hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin()) - t0) > rejectDigiFromPeak){
184  if(debug)
185  cout<<"digi skipped because t0 per sector "<<hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin())<<endl;
186  continue;
187  }
188 
189  //Use second bunch of events to compute a t0 reference per wire
191  //Fill reference wire histos
192  if(it_wire != wireIdWithHistos.end()){
193  TH1I* hT0WireHisto_ref = theHistoWireMap_ref[wireId];
194  theFile->cd();
195  if(hT0WireHisto_ref) hT0WireHisto_ref->Fill(t0);
196  }
197  if(!nDigiPerWire_ref[wireId]){
198  mK_ref[wireId] = 0;
199  }
200  nDigiPerWire_ref[wireId] = nDigiPerWire_ref[wireId] + 1;
201  mK_ref[wireId] = mK_ref[wireId] + (t0-mK_ref[wireId])/nDigiPerWire_ref[wireId];
202  }
203  //Use last all the remaining events to compute the mean and sigma t0 per wire
205  if(abs(t0-mK_ref[wireId]) > tpPeakWidth)
206  continue;
207  if(!nDigiPerWire[wireId]){
208  theAbsoluteT0PerWire[wireId] = 0;
209  qK[wireId] = 0;
210  mK[wireId] = 0;
211  }
212  nDigiPerWire[wireId] = nDigiPerWire[wireId] + 1;
213  theAbsoluteT0PerWire[wireId] = theAbsoluteT0PerWire[wireId] + t0;
214  //theSigmaT0PerWire[wireId] = theSigmaT0PerWire[wireId] + (t0*t0);
215  qK[wireId] = qK[wireId] + ((nDigiPerWire[wireId]-1)*(t0-mK[wireId])*(t0-mK[wireId])/nDigiPerWire[wireId]);
216  mK[wireId] = mK[wireId] + (t0-mK[wireId])/nDigiPerWire[wireId];
217  }
218  }//end if(nevents>1000)
219  }//end loop on digi
220  }//end loop on layer
221 
222  //Use the t0 per layer histos to have an indication about the t0 position
223  if(nevents == eventsForLayerT0){
224  for(map<DTLayerId, TH1I*>::const_iterator lHisto = theHistoLayerMap.begin();
225  lHisto != theHistoLayerMap.end();
226  ++lHisto){
227  if(debug)
228  cout<<"Reading histogram "<<(*lHisto).second->GetName()<<" with mean "<<(*lHisto).second->GetMean()<<" and RMS "<<(*lHisto).second->GetRMS();
229 
230  //Take the mean as a first t0 estimation
231  if((*lHisto).second->GetRMS()<5.0){
232  if(hT0SectorHisto == nullptr){
233  hT0SectorHisto = new TH1D("hT0AllLayerOfSector","T0 from pulses per layer in sector",
234  //20, (*lHisto).second->GetMean()-100, (*lHisto).second->GetMean()+100);
235  700, 0, 7000);
236  }
237  if(debug)
238  cout<<" accepted"<<endl;
239  hT0SectorHisto->Fill((*lHisto).second->GetMean());
240  }
241  //Take the mean of noise + 400ns as a first t0 estimation
242  // if((*lHisto).second->GetRMS()>10.0 && ((*lHisto).second->GetRMS()<15.0)){
243 // double t0_estim = (*lHisto).second->GetMean() + 400;
244 // if(hT0SectorHisto == 0){
245 // hT0SectorHisto = new TH1D("hT0AllLayerOfSector","T0 from pulses per layer in sector",
246 // //20, t0_estim-100, t0_estim+100);
247 // 700, 0, 7000);
248 // }
249 // if(debug)
250 // cout<<" accepted + 400ns"<<endl;
251 // hT0SectorHisto->Fill((*lHisto).second->GetMean() + 400);
252 // }
253  if(debug)
254  cout<<endl;
255 
256  theT0LayerMap[(*lHisto).second->GetName()] = (*lHisto).second->GetMean();
257  theSigmaT0LayerMap[(*lHisto).second->GetName()] = (*lHisto).second->GetRMS();
258  }
259  if(!hT0SectorHisto){
260  cout<<"[DTT0CalibrationRMS]: All the t0 per layer are still uncorrect: trying with greater number of events"<<endl;
262  return;
263  }
264  if(debug)
265  cout<<"[DTT0CalibrationRMS] t0 reference for this sector "<<
266  hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin())<<endl;
267  }
268 }
RunNumber_t run() const
Definition: EventID.h:39
EventNumber_t event() const
Definition: EventID.h:41
std::map< DTWireId, TH1I * > theHistoWireMap
std::map< std::string, double > theSigmaT0LayerMap
DTChamberId chamberId() const
Return the corresponding ChamberId.
std::map< std::string, double > theT0LayerMap
std::map< DTWireId, TH1I * > theHistoWireMap_ref
edm::ESHandle< DTGeometry > dtGeom
unsigned int rejectDigiFromPeak
DTSuperLayerId superlayerId() const
Return the corresponding SuperLayerId.
Definition: DTLayerId.h:59
std::vector< DTWireId > wireIdWithHistos
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
std::string theCalibSector
unsigned int eventsForWireT0
std::string getHistoName(const DTWireId &wId) const
unsigned int eventsForLayerT0
std::map< DTLayerId, TH1I * > theHistoLayerMap
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::map< DTWireId, double > theAbsoluteT0PerWire
std::map< DTWireId, double > mK
std::map< DTWireId, int > nDigiPerWire_ref
std::vector< DTDigi >::const_iterator const_iterator
std::map< DTWireId, int > nDigiPerWire
edm::EventID id() const
Definition: EventBase.h:60
int sector() const
Definition: DTChamberId.h:61
T get() const
Definition: EventSetup.h:62
std::pair< const_iterator, const_iterator > Range
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:45
std::map< DTWireId, double > mK_ref
std::map< DTWireId, double > qK
void DTT0CalibrationRMS::endJob ( void  )
overridevirtual

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

Loop on superlayer to correct between even-odd layers (2 different test pulse lines!)

Change t0 absolute reference -> from sector peak to chamber average

Write the t0 map into DB

Reimplemented from edm::EDAnalyzer.

Definition at line 271 of file DTT0CalibrationRMS.cc.

References funct::abs(), DTT0::begin(), DTSuperLayerId::chamberId(), DTTimeUnits::counts, gather_cfg::cout, debug, DTT0::end(), DTT0::get(), DTT0::set(), mathSSE::sqrt(), cscNeutronWriter_cfi::t0, interactiveExample::theFile, tzero, and DTCalibDBUtils::writeToDB().

Referenced by o2olib.O2ORunMgr::executeJob().

271  {
272 
273  DTT0* t0sAbsolute = new DTT0();
274  DTT0* t0sRelative = new DTT0();
275  DTT0* t0sWRTChamber = new DTT0();
276 
277  //if(debug)
278  cout << "[DTT0CalibrationRMSPerLayer]Writing histos to file!" << endl;
279 
280  theFile->cd();
281  theFile->WriteTObject(hT0SectorHisto);
282  //hT0SectorHisto->Write();
283  for(map<DTWireId, TH1I*>::const_iterator wHisto = theHistoWireMap.begin();
284  wHisto != theHistoWireMap.end();
285  ++wHisto) {
286  (*wHisto).second->Write();
287  }
288  for(map<DTWireId, TH1I*>::const_iterator wHisto = theHistoWireMap_ref.begin();
289  wHisto != theHistoWireMap_ref.end();
290  ++wHisto) {
291  (*wHisto).second->Write();
292  }
293  for(map<DTLayerId, TH1I*>::const_iterator lHisto = theHistoLayerMap.begin();
294  lHisto != theHistoLayerMap.end();
295  ++lHisto) {
296  (*lHisto).second->Write();
297  }
298 
299  //if(debug)
300  cout << "[DTT0CalibrationRMS] Compute and store t0 and sigma per wire" << endl;
301 
302  for(map<DTWireId, double>::const_iterator wiret0 = theAbsoluteT0PerWire.begin();
303  wiret0 != theAbsoluteT0PerWire.end();
304  ++wiret0){
305  if(nDigiPerWire[(*wiret0).first]){
306  double t0 = (*wiret0).second/nDigiPerWire[(*wiret0).first];
307 
308  theRelativeT0PerWire[(*wiret0).first] = t0 - hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin());
309 
310  //theSigmaT0PerWire[(*wiret0).first] = sqrt((theSigmaT0PerWire[(*wiret0).first] / nDigiPerWire[(*wiret0).first]) - t0*t0);
311  theSigmaT0PerWire[(*wiret0).first] = sqrt(qK[(*wiret0).first]/nDigiPerWire[(*wiret0).first]);
312 
313  cout << "Wire " << (*wiret0).first << " has t0 " << t0 << "(absolute) "
314  << theRelativeT0PerWire[(*wiret0).first] << "(relative)"
315  << " sigma " << theSigmaT0PerWire[(*wiret0).first] << endl;
316 
317  t0sAbsolute->set((*wiret0).first, t0, theSigmaT0PerWire[(*wiret0).first],DTTimeUnits::counts);
318  }
319  else{
320  cout<<"[DTT0CalibrationRMS] ERROR: no digis in wire "<<(*wiret0).first<<endl;
321  abort();
322  }
323  }
324 
327  // Get all the sls from the setup
328  const vector<const DTSuperLayer*> superLayers = dtGeom->superLayers();
329  // Loop over all SLs
330  for(vector<const DTSuperLayer*>::const_iterator sl = superLayers.begin();
331  sl != superLayers.end(); sl++) {
332 
333 
334  //Compute mean for odd and even superlayers
335  double oddLayersMean=0;
336  double evenLayersMean=0;
337  double oddLayersDen=0;
338  double evenLayersDen=0;
339  for(map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
340  wiret0 != theRelativeT0PerWire.end();
341  ++wiret0){
342  if((*wiret0).first.layerId().superlayerId() == (*sl)->id()){
343  if(debug)
344  cout<<"[DTT0CalibrationRMS] Superlayer "<<(*sl)->id()
345  <<"layer " <<(*wiret0).first.layerId().layer()<<" with "<<(*wiret0).second<<endl;
346  if(((*wiret0).first.layerId().layer()) % 2){
347  oddLayersMean = oddLayersMean + (*wiret0).second;
348  oddLayersDen++;
349  }
350  else{
351  evenLayersMean = evenLayersMean + (*wiret0).second;
352  evenLayersDen++;
353  }
354  }
355  }
356  oddLayersMean = oddLayersMean/oddLayersDen;
357  evenLayersMean = evenLayersMean/evenLayersDen;
358  //if(debug && oddLayersMean)
359  cout<<"[DTT0CalibrationRMS] Relative T0 mean for odd layers "<<oddLayersMean<<" even layers"<<evenLayersMean<<endl;
360 
361  //Compute sigma for odd and even superlayers
362  double oddLayersSigma=0;
363  double evenLayersSigma=0;
364  for(map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
365  wiret0 != theRelativeT0PerWire.end();
366  ++wiret0){
367  if((*wiret0).first.layerId().superlayerId() == (*sl)->id()){
368  if(((*wiret0).first.layerId().layer()) % 2){
369  oddLayersSigma = oddLayersSigma + ((*wiret0).second - oddLayersMean) * ((*wiret0).second - oddLayersMean);
370  }
371  else{
372  evenLayersSigma = evenLayersSigma + ((*wiret0).second - evenLayersMean) * ((*wiret0).second - evenLayersMean);
373  }
374  }
375  }
376  oddLayersSigma = oddLayersSigma/oddLayersDen;
377  evenLayersSigma = evenLayersSigma/evenLayersDen;
378  oddLayersSigma = sqrt(oddLayersSigma);
379  evenLayersSigma = sqrt(evenLayersSigma);
380 
381  //if(debug && oddLayersMean)
382  cout<<"[DTT0CalibrationRMS] Relative T0 sigma for odd layers "<<oddLayersSigma<<" even layers"<<evenLayersSigma<<endl;
383 
384  //Recompute the mean for odd and even superlayers discarding fluctations
385  double oddLayersFinalMean=0;
386  double evenLayersFinalMean=0;
387  for(map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
388  wiret0 != theRelativeT0PerWire.end();
389  ++wiret0){
390  if((*wiret0).first.layerId().superlayerId() == (*sl)->id()){
391  if(((*wiret0).first.layerId().layer()) % 2){
392  if(abs((*wiret0).second - oddLayersMean) < (2*oddLayersSigma))
393  oddLayersFinalMean = oddLayersFinalMean + (*wiret0).second;
394  }
395  else{
396  if(abs((*wiret0).second - evenLayersMean) < (2*evenLayersSigma))
397  evenLayersFinalMean = evenLayersFinalMean + (*wiret0).second;
398  }
399  }
400  }
401  oddLayersFinalMean = oddLayersFinalMean/oddLayersDen;
402  evenLayersFinalMean = evenLayersFinalMean/evenLayersDen;
403  //if(debug && oddLayersMean)
404  cout<<"[DTT0CalibrationRMS] Final relative T0 mean for odd layers "<<oddLayersFinalMean<<" even layers"<<evenLayersFinalMean<<endl;
405 
406  for(map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
407  wiret0 != theRelativeT0PerWire.end();
408  ++wiret0){
409  if((*wiret0).first.layerId().superlayerId() == (*sl)->id()){
410  double t0=-999;
411  if(((*wiret0).first.layerId().layer()) % 2)
412  t0 = (*wiret0).second + (evenLayersFinalMean - oddLayersFinalMean);
413  else
414  t0 = (*wiret0).second;
415 
416  cout << "[DTT0CalibrationRMS] Wire " << (*wiret0).first << " has t0 " << (*wiret0).second
417  << " (relative, after even-odd layer corrections) "
418  << " sigma " << theSigmaT0PerWire[(*wiret0).first] << endl;
419 
420  //Store the results into DB
421  t0sRelative->set((*wiret0).first, t0, theSigmaT0PerWire[(*wiret0).first],DTTimeUnits::counts);
422  }
423  }
424  }
425 
427  //if(debug)
428  cout << "[DTT0CalibrationRMS]Computing relative t0 wrt to chamber average" << endl;
429  //Compute the reference for each chamber
430  map<DTChamberId,double> sumT0ByChamber;
431  map<DTChamberId,int> countT0ByChamber;
432  for(DTT0::const_iterator tzero = t0sRelative->begin();
433  tzero != t0sRelative->end(); ++tzero) {
434  int channelId = tzero->channelId;
435  if ( channelId == 0 ) continue;
436  DTWireId wireId(channelId);
437  DTChamberId chamberId(wireId.chamberId());
438  //sumT0ByChamber[chamberId] = sumT0ByChamber[chamberId] + tzero->t0mean;
439  // @@@ better DTT0 usage
440  float t0mean_f;
441  float t0rms_f;
442  t0sRelative->get(wireId,t0mean_f,t0rms_f,DTTimeUnits::counts);
443  sumT0ByChamber[chamberId] = sumT0ByChamber[chamberId] + t0mean_f;
444  // @@@ NEW DTT0 END
445  countT0ByChamber[chamberId]++;
446  }
447 
448  //Change reference for each wire and store the new t0s in the new map
449  for(DTT0::const_iterator tzero = t0sRelative->begin();
450  tzero != t0sRelative->end(); ++tzero) {
451  int channelId = tzero->channelId;
452  if ( channelId == 0 ) continue;
453  DTWireId wireId(channelId);
454  DTChamberId chamberId(wireId.chamberId());
455  //double t0mean = (tzero->t0mean) - (sumT0ByChamber[chamberId]/countT0ByChamber[chamberId]);
456  //double t0rms = tzero->t0rms;
457  // @@@ better DTT0 usage
458  float t0mean_f;
459  float t0rms_f;
460  t0sRelative->get(wireId,t0mean_f,t0rms_f,DTTimeUnits::counts);
461  double t0mean = t0mean_f - (sumT0ByChamber[chamberId]/countT0ByChamber[chamberId]);
462  double t0rms = t0rms_f;
463  // @@@ NEW DTT0 END
464  t0sWRTChamber->set(wireId,
465  t0mean,
466  t0rms,
468  //if(debug)
469  //cout<<"Chamber "<<chamberId<<" has reference "<<(sumT0ByChamber[chamberId]/countT0ByChamber[chamberId]);
470  cout << "Changing t0 of wire " << wireId << " from " << t0mean_f
471  << " to " << t0mean << endl;
472  }
473  }
474 
476  if(debug)
477  cout << "[DTT0CalibrationRMS]Writing values in DB!" << endl;
478  // FIXME: to be read from cfg?
479  string t0Record = "DTT0Rcd";
480  // Write the t0 map to DB
481  if( correctByChamberMean_ ) DTCalibDBUtils::writeToDB(t0Record, t0sWRTChamber);
482  else DTCalibDBUtils::writeToDB(t0Record, t0sAbsolute);
483 }
const_iterator begin() const
Definition: DTT0.cc:204
int set(int wheelId, int stationId, int sectorId, int slId, int layerId, int cellId, float t0mean, float t0rms, DTTimeUnits::type unit)
Definition: DTT0.cc:140
std::map< DTWireId, TH1I * > theHistoWireMap
std::map< DTWireId, double > theSigmaT0PerWire
std::map< DTWireId, TH1I * > theHistoWireMap_ref
edm::ESHandle< DTGeometry > dtGeom
int get(int wheelId, int stationId, int sectorId, int slId, int layerId, int cellId, float &t0mean, float &t0rms, DTTimeUnits::type unit) const
Definition: DTT0.cc:67
std::vector< DTT0Data >::const_iterator const_iterator
Access methods to data.
Definition: DTT0.h:140
Definition: DTT0.h:53
std::map< DTWireId, double > theRelativeT0PerWire
T sqrt(T t)
Definition: SSEVec.h:18
std::map< DTLayerId, TH1I * > theHistoLayerMap
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::map< DTWireId, double > theAbsoluteT0PerWire
const_iterator end() const
Definition: DTT0.cc:209
std::map< DTWireId, int > nDigiPerWire
static const double tzero[3]
static void writeToDB(std::string record, T *payload)
const std::vector< const DTSuperLayer * > & superLayers() const
Return a vector of all SuperLayer.
Definition: DTGeometry.cc:107
std::map< DTWireId, double > qK
string DTT0CalibrationRMS::getHistoName ( const DTWireId wId) const
private

Definition at line 485 of file DTT0CalibrationRMS.cc.

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

485  {
486  string histoName = "Ch_" + std::to_string(wId.wheel()) + "_" + std::to_string(wId.station())
487  + "_" + std::to_string(wId.sector()) + "_SL" + std::to_string(wId.superlayer())
488  + "_L" + std::to_string(wId.layer()) + "_W" + std::to_string(wId.wire()) + "_hT0Histo";
489  return histoName;
490 }
int layer() const
Return the layer number.
Definition: DTLayerId.h:53
int wire() const
Return the wire number.
Definition: DTWireId.h:56
int superlayer() const
Return the superlayer number (deprecated method name)
int sector() const
Definition: DTChamberId.h:61
int station() const
Return the station number.
Definition: DTChamberId.h:51
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:45
string DTT0CalibrationRMS::getHistoName ( const DTLayerId lId) const
private

Definition at line 492 of file DTT0CalibrationRMS.cc.

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

492  {
493  string histoName = "Ch_" + std::to_string(lId.wheel()) + "_" + std::to_string(lId.station())
494  + "_" + std::to_string(lId.sector()) + "_SL" + std::to_string(lId.superlayer())
495  + "_L" + std::to_string(lId.layer()) + "_hT0Histo";
496  return histoName;
497 }
int layer() const
Return the layer number.
Definition: DTLayerId.h:53
int superlayer() const
Return the superlayer number (deprecated method name)
int sector() const
Definition: DTChamberId.h:61
int station() const
Return the station number.
Definition: DTChamberId.h:51
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:45

Member Data Documentation

std::vector<std::string> DTT0CalibrationRMS::cellsWithHistos
private

Definition at line 90 of file DTT0CalibrationRMS.h.

bool DTT0CalibrationRMS::correctByChamberMean_
private

Definition at line 75 of file DTT0CalibrationRMS.h.

bool DTT0CalibrationRMS::debug
private
std::string DTT0CalibrationRMS::digiLabel
private

Definition at line 54 of file DTT0CalibrationRMS.h.

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

Definition at line 109 of file DTT0CalibrationRMS.h.

unsigned int DTT0CalibrationRMS::eventsForLayerT0
private

Definition at line 64 of file DTT0CalibrationRMS.h.

unsigned int DTT0CalibrationRMS::eventsForWireT0
private

Definition at line 66 of file DTT0CalibrationRMS.h.

TH1D* DTT0CalibrationRMS::hT0SectorHisto
private

Definition at line 86 of file DTT0CalibrationRMS.h.

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

Definition at line 98 of file DTT0CalibrationRMS.h.

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

Definition at line 99 of file DTT0CalibrationRMS.h.

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

Definition at line 96 of file DTT0CalibrationRMS.h.

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

Definition at line 97 of file DTT0CalibrationRMS.h.

unsigned int DTT0CalibrationRMS::nevents
private

Definition at line 62 of file DTT0CalibrationRMS.h.

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

Definition at line 100 of file DTT0CalibrationRMS.h.

unsigned int DTT0CalibrationRMS::rejectDigiFromPeak
private

Definition at line 69 of file DTT0CalibrationRMS.h.

int DTT0CalibrationRMS::selSector
private

Definition at line 81 of file DTT0CalibrationRMS.h.

int DTT0CalibrationRMS::selWheel
private

Definition at line 79 of file DTT0CalibrationRMS.h.

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

Definition at line 93 of file DTT0CalibrationRMS.h.

std::string DTT0CalibrationRMS::theCalibSector
private

Definition at line 80 of file DTT0CalibrationRMS.h.

std::string DTT0CalibrationRMS::theCalibWheel
private

Definition at line 78 of file DTT0CalibrationRMS.h.

TFile* DTT0CalibrationRMS::theFile
private

Definition at line 57 of file DTT0CalibrationRMS.h.

std::map<DTLayerId, TH1I*> DTT0CalibrationRMS::theHistoLayerMap
private

Definition at line 84 of file DTT0CalibrationRMS.h.

std::map<DTWireId,TH1I*> DTT0CalibrationRMS::theHistoWireMap
private

Definition at line 102 of file DTT0CalibrationRMS.h.

std::map<DTWireId,TH1I*> DTT0CalibrationRMS::theHistoWireMap_ref
private

Definition at line 103 of file DTT0CalibrationRMS.h.

TFile* DTT0CalibrationRMS::theOutputFile
private

Definition at line 59 of file DTT0CalibrationRMS.h.

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

Definition at line 94 of file DTT0CalibrationRMS.h.

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

Definition at line 106 of file DTT0CalibrationRMS.h.

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

Definition at line 95 of file DTT0CalibrationRMS.h.

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

Definition at line 105 of file DTT0CalibrationRMS.h.

double DTT0CalibrationRMS::tpPeakWidth
private

Definition at line 72 of file DTT0CalibrationRMS.h.

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

Definition at line 89 of file DTT0CalibrationRMS.h.