Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:34:05

0001 /*
0002  *  See header file for a description of this class.
0003  *
0004  *  \author S. Bolognesi - INFN Torino
0005  */
0006 #include "CalibMuon/DTCalibration/plugins/DTT0CalibrationRMS.h"
0007 #include "CalibMuon/DTCalibration/interface/DTCalibDBUtils.h"
0008 
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 
0013 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0014 #include "Geometry/Records/interface/MuonNumberingRecord.h"
0015 
0016 #include "CondFormats/DTObjects/interface/DTT0.h"
0017 
0018 #include "TH1I.h"
0019 #include "TFile.h"
0020 #include "TKey.h"
0021 
0022 using namespace std;
0023 using namespace edm;
0024 // using namespace cond;
0025 
0026 // Constructor
0027 DTT0CalibrationRMS::DTT0CalibrationRMS(const edm::ParameterSet& pset) : dtGeomToken_(esConsumes()) {
0028   // Get the debug parameter for verbose output
0029   debug = pset.getUntrackedParameter<bool>("debug");
0030   if (debug)
0031     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS]Constructor called!";
0032 
0033   // Get the token to retrieve digis from the event
0034   digiToken = consumes<DTDigiCollection>(pset.getUntrackedParameter<string>("digiLabel"));
0035 
0036   // The root file which contain the histos per layer
0037   string rootFileName = pset.getUntrackedParameter<string>("rootFileName", "DTT0PerLayer.root");
0038   theFile = new TFile(rootFileName.c_str(), "RECREATE");
0039 
0040   theCalibWheel =
0041       pset.getUntrackedParameter<string>("calibWheel", "All");  //FIXME amke a vector of integer instead of a string
0042   if (theCalibWheel != "All") {
0043     stringstream linestr;
0044     int selWheel;
0045     linestr << theCalibWheel;
0046     linestr >> selWheel;
0047     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMSPerLayer] chosen wheel " << selWheel;
0048   }
0049 
0050   // Sector/s to calibrate
0051   theCalibSector =
0052       pset.getUntrackedParameter<string>("calibSector", "All");  //FIXME amke a vector of integer instead of a string
0053   if (theCalibSector != "All") {
0054     stringstream linestr;
0055     int selSector;
0056     linestr << theCalibSector;
0057     linestr >> selSector;
0058     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMSPerLayer] chosen sector " << selSector;
0059   }
0060 
0061   vector<string> defaultCell;
0062   defaultCell.push_back("None");
0063   cellsWithHistos = pset.getUntrackedParameter<vector<string> >("cellsWithHisto", defaultCell);
0064   for (vector<string>::const_iterator cell = cellsWithHistos.begin(); cell != cellsWithHistos.end(); ++cell) {
0065     if ((*cell) != "None") {
0066       stringstream linestr;
0067       int wheel, sector, station, sl, layer, wire;
0068       linestr << (*cell);
0069       linestr >> wheel >> sector >> station >> sl >> layer >> wire;
0070       wireIdWithHistos.push_back(DTWireId(wheel, station, sector, sl, layer, wire));
0071     }
0072   }
0073 
0074   hT0SectorHisto = nullptr;
0075 
0076   nevents = 0;
0077   eventsForLayerT0 = pset.getParameter<unsigned int>("eventsForLayerT0");
0078   eventsForWireT0 = pset.getParameter<unsigned int>("eventsForWireT0");
0079   rejectDigiFromPeak = pset.getParameter<unsigned int>("rejectDigiFromPeak");
0080   tpPeakWidth = pset.getParameter<double>("tpPeakWidth");
0081   //useReferenceWireInLayer_ = true;
0082   correctByChamberMean_ = pset.getParameter<bool>("correctByChamberMean");
0083 }
0084 
0085 // Destructor
0086 DTT0CalibrationRMS::~DTT0CalibrationRMS() {
0087   if (debug)
0088     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS]Destructor called!";
0089 
0090   theFile->Close();
0091 }
0092 
0093 /// Perform the real analysis
0094 void DTT0CalibrationRMS::analyze(const edm::Event& event, const edm::EventSetup& eventSetup) {
0095   if (debug || event.id().event() % 500 == 0)
0096     edm::LogVerbatim("DTCalibration") << "--- [DTT0CalibrationRMS] Analysing Event: #Run: " << event.id().run()
0097                                       << " #Event: " << event.id().event();
0098   nevents++;
0099 
0100   // Get the digis from the event
0101   const Handle<DTDigiCollection>& digis = event.getHandle(digiToken);
0102 
0103   // Get the DT Geometry
0104   dtGeom = eventSetup.getHandle(dtGeomToken_);
0105 
0106   // Iterate through all digi collections ordered by LayerId
0107   DTDigiCollection::DigiRangeIterator dtLayerIt;
0108   for (dtLayerIt = digis->begin(); dtLayerIt != digis->end(); ++dtLayerIt) {
0109     // Get the iterators over the digis associated with this LayerId
0110     const DTDigiCollection::Range& digiRange = (*dtLayerIt).second;
0111 
0112     // Get the layerId
0113     const DTLayerId layerId = (*dtLayerIt).first;  //FIXME: check to be in the right sector
0114 
0115     if ((theCalibWheel != "All") && (layerId.superlayerId().chamberId().wheel() != selWheel))
0116       continue;
0117     if ((theCalibSector != "All") && (layerId.superlayerId().chamberId().sector() != selSector))
0118       continue;
0119 
0120     // Loop over all digis in the given layer
0121     for (DTDigiCollection::const_iterator digi = digiRange.first; digi != digiRange.second; ++digi) {
0122       double t0 = (*digi).countsTDC();
0123 
0124       //Use first bunch of events to fill t0 per layer
0125       if (nevents < eventsForLayerT0) {
0126         //Get the per-layer histo from the map
0127         TH1I* hT0LayerHisto = theHistoLayerMap[layerId];
0128         //If it doesn't exist, book it
0129         if (hT0LayerHisto == nullptr) {
0130           theFile->cd();
0131           hT0LayerHisto = new TH1I(getHistoName(layerId).c_str(),
0132                                    "T0 from pulses by layer (TDC counts, 1 TDC count = 0.781 ns)",
0133                                    200,
0134                                    t0 - 100,
0135                                    t0 + 100);
0136           if (debug)
0137             edm::LogVerbatim("DTCalibration") << "  New T0 per Layer Histo: " << hT0LayerHisto->GetName();
0138           theHistoLayerMap[layerId] = hT0LayerHisto;
0139         }
0140 
0141         //Fill the histos
0142         theFile->cd();
0143         if (hT0LayerHisto != nullptr) {
0144           hT0LayerHisto->Fill(t0);
0145         }
0146       }
0147 
0148       //Use all the remaining events to compute t0 per wire
0149       if (nevents > eventsForLayerT0) {
0150         // Get the wireId
0151         const DTWireId wireId(layerId, (*digi).wire());
0152         if (debug) {
0153           edm::LogVerbatim("DTCalibration")
0154               << "   Wire: " << wireId << "\n       time (TDC counts): " << (*digi).countsTDC();
0155         }
0156 
0157         //Fill the histos per wire for the chosen cells
0158         vector<DTWireId>::iterator it_wire = find(wireIdWithHistos.begin(), wireIdWithHistos.end(), wireId);
0159         if (it_wire != wireIdWithHistos.end()) {
0160           if (theHistoWireMap.find(wireId) == theHistoWireMap.end()) {
0161             theHistoWireMap[wireId] = new TH1I(getHistoName(wireId).c_str(),
0162                                                "T0 from pulses by wire (TDC counts, 1 TDC count = 0.781 ns)",
0163                                                7000,
0164                                                0,
0165                                                7000);
0166             if (debug)
0167               edm::LogVerbatim("DTCalibration") << "  New T0 per wire Histo: " << (theHistoWireMap[wireId])->GetName();
0168           }
0169           if (theHistoWireMap_ref.find(wireId) == theHistoWireMap_ref.end()) {
0170             theHistoWireMap_ref[wireId] = new TH1I((getHistoName(wireId) + "_ref").c_str(),
0171                                                    "T0 from pulses by wire (TDC counts, 1 TDC count = 0.781 ns)",
0172                                                    7000,
0173                                                    0,
0174                                                    7000);
0175             if (debug)
0176               edm::LogVerbatim("DTCalibration")
0177                   << "  New T0 per wire Histo: " << (theHistoWireMap_ref[wireId])->GetName();
0178           }
0179 
0180           TH1I* hT0WireHisto = theHistoWireMap[wireId];
0181           //Fill the histos
0182           theFile->cd();
0183           if (hT0WireHisto)
0184             hT0WireHisto->Fill(t0);
0185         }
0186 
0187         //Check the tzero has reasonable value
0188         if (abs(hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin()) - t0) > rejectDigiFromPeak) {
0189           if (debug)
0190             edm::LogVerbatim("DTCalibration") << "digi skipped because t0 per sector "
0191                                               << hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin());
0192           continue;
0193         }
0194 
0195         //Use second bunch of events to compute a t0 reference per wire
0196         if (nevents < (eventsForLayerT0 + eventsForWireT0)) {
0197           //Fill reference wire histos
0198           if (it_wire != wireIdWithHistos.end()) {
0199             TH1I* hT0WireHisto_ref = theHistoWireMap_ref[wireId];
0200             theFile->cd();
0201             if (hT0WireHisto_ref)
0202               hT0WireHisto_ref->Fill(t0);
0203           }
0204           if (!nDigiPerWire_ref[wireId]) {
0205             mK_ref[wireId] = 0;
0206           }
0207           nDigiPerWire_ref[wireId] = nDigiPerWire_ref[wireId] + 1;
0208           mK_ref[wireId] = mK_ref[wireId] + (t0 - mK_ref[wireId]) / nDigiPerWire_ref[wireId];
0209         }
0210         //Use last all the remaining events to compute the mean and sigma t0 per wire
0211         else if (nevents > (eventsForLayerT0 + eventsForWireT0)) {
0212           if (abs(t0 - mK_ref[wireId]) > tpPeakWidth)
0213             continue;
0214           if (!nDigiPerWire[wireId]) {
0215             theAbsoluteT0PerWire[wireId] = 0;
0216             qK[wireId] = 0;
0217             mK[wireId] = 0;
0218           }
0219           nDigiPerWire[wireId] = nDigiPerWire[wireId] + 1;
0220           theAbsoluteT0PerWire[wireId] = theAbsoluteT0PerWire[wireId] + t0;
0221           //theSigmaT0PerWire[wireId] = theSigmaT0PerWire[wireId] + (t0*t0);
0222           qK[wireId] =
0223               qK[wireId] + ((nDigiPerWire[wireId] - 1) * (t0 - mK[wireId]) * (t0 - mK[wireId]) / nDigiPerWire[wireId]);
0224           mK[wireId] = mK[wireId] + (t0 - mK[wireId]) / nDigiPerWire[wireId];
0225         }
0226       }  //end if(nevents>1000)
0227     }    //end loop on digi
0228   }      //end loop on layer
0229 
0230   //Use the t0 per layer histos to have an indication about the t0 position
0231   if (nevents == eventsForLayerT0) {
0232     for (map<DTLayerId, TH1I*>::const_iterator lHisto = theHistoLayerMap.begin(); lHisto != theHistoLayerMap.end();
0233          ++lHisto) {
0234       if (debug)
0235         edm::LogVerbatim("DTCalibration") << "Reading histogram " << (*lHisto).second->GetName() << " with mean "
0236                                           << (*lHisto).second->GetMean() << " and RMS " << (*lHisto).second->GetRMS();
0237 
0238       //Take the mean as a first t0 estimation
0239       if ((*lHisto).second->GetRMS() < 5.0) {
0240         if (hT0SectorHisto == nullptr) {
0241           hT0SectorHisto = new TH1D("hT0AllLayerOfSector",
0242                                     "T0 from pulses per layer in sector",
0243                                     //20, (*lHisto).second->GetMean()-100, (*lHisto).second->GetMean()+100);
0244                                     700,
0245                                     0,
0246                                     7000);
0247         }
0248         if (debug)
0249           edm::LogVerbatim("DTCalibration") << " accepted";
0250         hT0SectorHisto->Fill((*lHisto).second->GetMean());
0251       }
0252       //Take the mean of noise + 400ns as a first t0 estimation
0253       // if((*lHisto).second->GetRMS()>10.0 && ((*lHisto).second->GetRMS()<15.0)){
0254       //    double t0_estim = (*lHisto).second->GetMean() + 400;
0255       //    if(hT0SectorHisto == 0){
0256       //      hT0SectorHisto = new TH1D("hT0AllLayerOfSector","T0 from pulses per layer in sector",
0257       //                    //20, t0_estim-100, t0_estim+100);
0258       //                    700, 0, 7000);
0259       //    }
0260       //    if(debug)
0261       //      edm::LogVerbatim("DTCalibration")<<" accepted + 400ns";
0262       //    hT0SectorHisto->Fill((*lHisto).second->GetMean() + 400);
0263       //       }
0264       //if (debug)
0265       //  edm::LogVerbatim("DTCalibration");
0266 
0267       theT0LayerMap[(*lHisto).second->GetName()] = (*lHisto).second->GetMean();
0268       theSigmaT0LayerMap[(*lHisto).second->GetName()] = (*lHisto).second->GetRMS();
0269     }
0270     if (!hT0SectorHisto) {
0271       edm::LogVerbatim("DTCalibration")
0272           << "[DTT0CalibrationRMS]: All the t0 per layer are still uncorrect: trying with greater number of events";
0273       eventsForLayerT0 = eventsForLayerT0 * 2;
0274       return;
0275     }
0276     if (debug)
0277       edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] t0 reference for this sector "
0278                                         << hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin());
0279   }
0280 }
0281 
0282 void DTT0CalibrationRMS::endJob() {
0283   DTT0* t0sAbsolute = new DTT0();
0284   DTT0* t0sRelative = new DTT0();
0285   DTT0* t0sWRTChamber = new DTT0();
0286 
0287   //if(debug)
0288   edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMSPerLayer]Writing histos to file!";
0289 
0290   theFile->cd();
0291   theFile->WriteTObject(hT0SectorHisto);
0292   //hT0SectorHisto->Write();
0293   for (map<DTWireId, TH1I*>::const_iterator wHisto = theHistoWireMap.begin(); wHisto != theHistoWireMap.end();
0294        ++wHisto) {
0295     (*wHisto).second->Write();
0296   }
0297   for (map<DTWireId, TH1I*>::const_iterator wHisto = theHistoWireMap_ref.begin(); wHisto != theHistoWireMap_ref.end();
0298        ++wHisto) {
0299     (*wHisto).second->Write();
0300   }
0301   for (map<DTLayerId, TH1I*>::const_iterator lHisto = theHistoLayerMap.begin(); lHisto != theHistoLayerMap.end();
0302        ++lHisto) {
0303     (*lHisto).second->Write();
0304   }
0305 
0306   //if(debug)
0307   edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Compute and store t0 and sigma per wire";
0308 
0309   for (map<DTWireId, double>::const_iterator wiret0 = theAbsoluteT0PerWire.begin();
0310        wiret0 != theAbsoluteT0PerWire.end();
0311        ++wiret0) {
0312     if (nDigiPerWire[(*wiret0).first]) {
0313       double t0 = (*wiret0).second / nDigiPerWire[(*wiret0).first];
0314 
0315       theRelativeT0PerWire[(*wiret0).first] = t0 - hT0SectorHisto->GetBinCenter(hT0SectorHisto->GetMaximumBin());
0316 
0317       //theSigmaT0PerWire[(*wiret0).first] = sqrt((theSigmaT0PerWire[(*wiret0).first] / nDigiPerWire[(*wiret0).first]) - t0*t0);
0318       theSigmaT0PerWire[(*wiret0).first] = sqrt(qK[(*wiret0).first] / nDigiPerWire[(*wiret0).first]);
0319 
0320       edm::LogVerbatim("DTCalibration") << "Wire " << (*wiret0).first << " has t0 " << t0 << "(absolute) "
0321                                         << theRelativeT0PerWire[(*wiret0).first] << "(relative)"
0322                                         << "    sigma " << theSigmaT0PerWire[(*wiret0).first];
0323 
0324       t0sAbsolute->set((*wiret0).first, t0, theSigmaT0PerWire[(*wiret0).first], DTTimeUnits::counts);
0325     } else {
0326       edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] ERROR: no digis in wire " << (*wiret0).first;
0327       abort();
0328     }
0329   }
0330 
0331   if (correctByChamberMean_) {
0332     ///Loop on superlayer to correct between even-odd layers (2 different test pulse lines!)
0333     // Get all the sls from the setup
0334     const vector<const DTSuperLayer*> superLayers = dtGeom->superLayers();
0335     // Loop over all SLs
0336     for (vector<const DTSuperLayer*>::const_iterator sl = superLayers.begin(); sl != superLayers.end(); sl++) {
0337       //Compute mean for odd and even superlayers
0338       double oddLayersMean = 0;
0339       double evenLayersMean = 0;
0340       double oddLayersDen = 0;
0341       double evenLayersDen = 0;
0342       for (map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
0343            wiret0 != theRelativeT0PerWire.end();
0344            ++wiret0) {
0345         if ((*wiret0).first.layerId().superlayerId() == (*sl)->id()) {
0346           if (debug)
0347             edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Superlayer " << (*sl)->id() << "layer "
0348                                               << (*wiret0).first.layerId().layer() << " with " << (*wiret0).second;
0349           if (((*wiret0).first.layerId().layer()) % 2) {
0350             oddLayersMean = oddLayersMean + (*wiret0).second;
0351             oddLayersDen++;
0352           } else {
0353             evenLayersMean = evenLayersMean + (*wiret0).second;
0354             evenLayersDen++;
0355           }
0356         }
0357       }
0358       oddLayersMean = oddLayersMean / oddLayersDen;
0359       evenLayersMean = evenLayersMean / evenLayersDen;
0360       //if(debug && oddLayersMean)
0361       edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Relative T0 mean for  odd layers " << oddLayersMean
0362                                         << "  even layers" << evenLayersMean;
0363 
0364       //Compute sigma for odd and even superlayers
0365       double oddLayersSigma = 0;
0366       double evenLayersSigma = 0;
0367       for (map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
0368            wiret0 != theRelativeT0PerWire.end();
0369            ++wiret0) {
0370         if ((*wiret0).first.layerId().superlayerId() == (*sl)->id()) {
0371           if (((*wiret0).first.layerId().layer()) % 2) {
0372             oddLayersSigma = oddLayersSigma + ((*wiret0).second - oddLayersMean) * ((*wiret0).second - oddLayersMean);
0373           } else {
0374             evenLayersSigma =
0375                 evenLayersSigma + ((*wiret0).second - evenLayersMean) * ((*wiret0).second - evenLayersMean);
0376           }
0377         }
0378       }
0379       oddLayersSigma = oddLayersSigma / oddLayersDen;
0380       evenLayersSigma = evenLayersSigma / evenLayersDen;
0381       oddLayersSigma = sqrt(oddLayersSigma);
0382       evenLayersSigma = sqrt(evenLayersSigma);
0383 
0384       //if(debug && oddLayersMean)
0385       edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Relative T0 sigma for  odd layers " << oddLayersSigma
0386                                         << "  even layers" << evenLayersSigma;
0387 
0388       //Recompute the mean for odd and even superlayers discarding fluctations
0389       double oddLayersFinalMean = 0;
0390       double evenLayersFinalMean = 0;
0391       for (map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
0392            wiret0 != theRelativeT0PerWire.end();
0393            ++wiret0) {
0394         if ((*wiret0).first.layerId().superlayerId() == (*sl)->id()) {
0395           if (((*wiret0).first.layerId().layer()) % 2) {
0396             if (abs((*wiret0).second - oddLayersMean) < (2 * oddLayersSigma))
0397               oddLayersFinalMean = oddLayersFinalMean + (*wiret0).second;
0398           } else {
0399             if (abs((*wiret0).second - evenLayersMean) < (2 * evenLayersSigma))
0400               evenLayersFinalMean = evenLayersFinalMean + (*wiret0).second;
0401           }
0402         }
0403       }
0404       oddLayersFinalMean = oddLayersFinalMean / oddLayersDen;
0405       evenLayersFinalMean = evenLayersFinalMean / evenLayersDen;
0406       //if(debug && oddLayersMean)
0407       edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Final relative T0 mean for  odd layers "
0408                                         << oddLayersFinalMean << "  even layers" << evenLayersFinalMean;
0409 
0410       for (map<DTWireId, double>::const_iterator wiret0 = theRelativeT0PerWire.begin();
0411            wiret0 != theRelativeT0PerWire.end();
0412            ++wiret0) {
0413         if ((*wiret0).first.layerId().superlayerId() == (*sl)->id()) {
0414           double t0 = -999;
0415           if (((*wiret0).first.layerId().layer()) % 2)
0416             t0 = (*wiret0).second + (evenLayersFinalMean - oddLayersFinalMean);
0417           else
0418             t0 = (*wiret0).second;
0419 
0420           edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS] Wire " << (*wiret0).first << " has t0 "
0421                                             << (*wiret0).second << " (relative, after even-odd layer corrections)  "
0422                                             << "    sigma " << theSigmaT0PerWire[(*wiret0).first];
0423 
0424           //Store the results into DB
0425           t0sRelative->set((*wiret0).first, t0, theSigmaT0PerWire[(*wiret0).first], DTTimeUnits::counts);
0426         }
0427       }
0428     }
0429 
0430     ///Change t0 absolute reference -> from sector peak to chamber average
0431     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS]Computing relative t0 wrt to chamber average";
0432     //Compute the reference for each chamber
0433     map<DTChamberId, double> sumT0ByChamber;
0434     map<DTChamberId, int> countT0ByChamber;
0435     for (DTT0::const_iterator tzero = t0sRelative->begin(); tzero != t0sRelative->end(); ++tzero) {
0436       int channelId = tzero->channelId;
0437       if (channelId == 0)
0438         continue;
0439       DTWireId wireId(channelId);
0440       DTChamberId chamberId(wireId.chamberId());
0441       //sumT0ByChamber[chamberId] = sumT0ByChamber[chamberId] + tzero->t0mean;
0442       // @@@ better DTT0 usage
0443       float t0mean_f;
0444       float t0rms_f;
0445       t0sRelative->get(wireId, t0mean_f, t0rms_f, DTTimeUnits::counts);
0446       sumT0ByChamber[chamberId] = sumT0ByChamber[chamberId] + t0mean_f;
0447       // @@@ NEW DTT0 END
0448       countT0ByChamber[chamberId]++;
0449     }
0450 
0451     //Change reference for each wire and store the new t0s in the new map
0452     for (DTT0::const_iterator tzero = t0sRelative->begin(); tzero != t0sRelative->end(); ++tzero) {
0453       int channelId = tzero->channelId;
0454       if (channelId == 0)
0455         continue;
0456       DTWireId wireId(channelId);
0457       DTChamberId chamberId(wireId.chamberId());
0458       //double t0mean = (tzero->t0mean) - (sumT0ByChamber[chamberId]/countT0ByChamber[chamberId]);
0459       //double t0rms = tzero->t0rms;
0460       // @@@ better DTT0 usage
0461       float t0mean_f;
0462       float t0rms_f;
0463       t0sRelative->get(wireId, t0mean_f, t0rms_f, DTTimeUnits::counts);
0464       double t0mean = t0mean_f - (sumT0ByChamber[chamberId] / countT0ByChamber[chamberId]);
0465       double t0rms = t0rms_f;
0466       // @@@ NEW DTT0 END
0467       t0sWRTChamber->set(wireId, t0mean, t0rms, DTTimeUnits::counts);
0468       edm::LogVerbatim("DTCalibration") << "Changing t0 of wire " << wireId << " from " << t0mean_f << " to " << t0mean;
0469     }
0470   }
0471 
0472   ///Write the t0 map into DB
0473   if (debug)
0474     edm::LogVerbatim("DTCalibration") << "[DTT0CalibrationRMS]Writing values in DB!";
0475   // FIXME: to be read from cfg?
0476   string t0Record = "DTT0Rcd";
0477   // Write the t0 map to DB
0478   if (correctByChamberMean_)
0479     DTCalibDBUtils::writeToDB(t0Record, t0sWRTChamber);
0480   else
0481     DTCalibDBUtils::writeToDB(t0Record, t0sAbsolute);
0482 }
0483 
0484 string DTT0CalibrationRMS::getHistoName(const DTWireId& wId) const {
0485   string histoName = "Ch_" + std::to_string(wId.wheel()) + "_" + std::to_string(wId.station()) + "_" +
0486                      std::to_string(wId.sector()) + "_SL" + std::to_string(wId.superlayer()) + "_L" +
0487                      std::to_string(wId.layer()) + "_W" + std::to_string(wId.wire()) + "_hT0Histo";
0488   return histoName;
0489 }
0490 
0491 string DTT0CalibrationRMS::getHistoName(const DTLayerId& lId) const {
0492   string histoName = "Ch_" + std::to_string(lId.wheel()) + "_" + std::to_string(lId.station()) + "_" +
0493                      std::to_string(lId.sector()) + "_SL" + std::to_string(lId.superlayer()) + "_L" +
0494                      std::to_string(lId.layer()) + "_hT0Histo";
0495   return histoName;
0496 }