CMS 3D CMS Logo

DTTTrigT0SegCorrection.cc
Go to the documentation of this file.
1 /*
2  * See header file for a description of this class.
3  *
4  * \author A. Vilela Pereira
5  */
6 
14 
15 #include "TFile.h"
16 #include "TH1F.h"
17 
18 #include <string>
19 #include <sstream>
20 
21 using namespace std;
22 using namespace edm;
23 
24 namespace dtCalibration {
25 
26 DTTTrigT0SegCorrection::DTTTrigT0SegCorrection(const ParameterSet& pset) {
27  string t0SegRootFile = pset.getParameter<string>("t0SegRootFile");
28  rootFile_ = new TFile(t0SegRootFile.c_str(),"READ");
29  dbLabel = pset.getUntrackedParameter<string>("dbLabel", "");
30 }
31 
32 DTTTrigT0SegCorrection::~DTTTrigT0SegCorrection() {
33  delete rootFile_;
34 }
35 
36 void DTTTrigT0SegCorrection::setES(const EventSetup& setup) {
37  // Get tTrig record from DB
38  ESHandle<DTTtrig> tTrig;
39  setup.get<DTTtrigRcd>().get(dbLabel,tTrig);
40  tTrigMap_ = &*tTrig;
41 }
42 
43 DTTTrigData DTTTrigT0SegCorrection::correction(const DTSuperLayerId& slId) {
44  float tTrigMean,tTrigSigma,kFactor;
45  int status = tTrigMap_->get(slId,tTrigMean,tTrigSigma,kFactor,DTTimeUnits::ns);
46  if(status != 0) throw cms::Exception("[DTTTrigT0SegCorrection]") << "Could not find tTrig entry in DB for"
47  << slId << endl;
48 
49  const TH1F* t0SegHisto = getHisto(slId);
50  double corrMean = tTrigMean;
51  double corrSigma = tTrigSigma;
52  //FIXME: can we fit the t0seg histo? How do we remove the peak at 0?;
53  double corrKFact = (kFactor*tTrigSigma + t0SegHisto->GetMean())/tTrigSigma;
54  return DTTTrigData(corrMean,corrSigma,corrKFact);
55 }
56 
57 const TH1F* DTTTrigT0SegCorrection::getHisto(const DTSuperLayerId& slId) {
58  string histoName = getHistoName(slId);
59  TH1F* histo = static_cast<TH1F*>(rootFile_->Get(histoName.c_str()));
60  if(!histo) throw cms::Exception("[DTTTrigT0SegCorrection]") << "t0-seg histogram not found:"
61  << histoName << endl;
62  return histo;
63 }
64 
66  DTChamberId chId = slId.chamberId();
67 
68  // Compose the chamber name
69  std::string wheel = std::to_string(chId.wheel());
70  std::string station = std::to_string(chId.station());
71  std::string sector = std::to_string(chId.sector());
72 
73  string chHistoName =
74  "_W" + wheel +
75  "_St" + station +
76  "_Sec" + sector;
77 
78  return (slId.superLayer() != 2)?("hRPhiSegT0"+chHistoName):("hRZSegT0"+chHistoName);
79 }
80 
81 } // namespace
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
DTChamberId chamberId() const
Return the corresponding ChamberId.
def setup(process, global_tag, zero_tesla=False)
Definition: GeneralSetup.py:2
def getHistoName(wheel, station, sector)
int superLayer() const
Return the superlayer number.
HLT enums.
int sector() const
Definition: DTChamberId.h:61
T get() const
Definition: EventSetup.h:71
int station() const
Return the station number.
Definition: DTChamberId.h:51
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:45