Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:58:28

0001 /*
0002  *  See header file for a description of this class.
0003  *
0004  *  \author S. Bolognesi
0005  */
0006 
0007 #include "CalibMuon/DTCalibration/plugins/DTTTrigWriter.h"
0008 #include "CalibMuon/DTCalibration/interface/DTTimeBoxFitter.h"
0009 #include "CalibMuon/DTCalibration/interface/DTCalibDBUtils.h"
0010 
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015 #include "DataFormats/MuonDetId/interface/DTSuperLayerId.h"
0016 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0017 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0018 
0019 #include "CondFormats/DTObjects/interface/DTTtrig.h"
0020 #include "CondFormats/DataRecord/interface/DTStatusFlagRcd.h"
0021 #include "CondFormats/DTObjects/interface/DTStatusFlag.h"
0022 
0023 /* C++ Headers */
0024 #include <vector>
0025 #include <iostream>
0026 #include <fstream>
0027 #include <string>
0028 #include <sstream>
0029 #include "TFile.h"
0030 #include "TH1.h"
0031 
0032 using namespace std;
0033 using namespace edm;
0034 
0035 // Constructor
0036 DTTTrigWriter::DTTTrigWriter(const ParameterSet& pset) : dtGeomToken_(esConsumes()) {
0037   // get selected debug option
0038   debug = pset.getUntrackedParameter<bool>("debug", false);
0039 
0040   // Open the root file which contains the histos
0041   theRootInputFile = pset.getUntrackedParameter<string>("rootFileName");
0042   theFile = new TFile(theRootInputFile.c_str(), "READ");
0043   theFile->cd();
0044   theFitter = new DTTimeBoxFitter();
0045   if (debug)
0046     theFitter->setVerbosity(1);
0047 
0048   double sigmaFit = pset.getUntrackedParameter<double>("sigmaTTrigFit", 10.);
0049   theFitter->setFitSigma(sigmaFit);
0050 
0051   // the kfactor to be uploaded in the ttrig DB
0052   kFactor = pset.getUntrackedParameter<double>("kFactor", -0.7);
0053 
0054   // Create the object to be written to DB
0055   tTrig = new DTTtrig();
0056 
0057   if (debug)
0058     cout << "[DTTTrigWriter]Constructor called!" << endl;
0059 }
0060 
0061 // Destructor
0062 DTTTrigWriter::~DTTTrigWriter() {
0063   if (debug)
0064     cout << "[DTTTrigWriter]Destructor called!" << endl;
0065   theFile->Close();
0066   delete theFitter;
0067 }
0068 
0069 // Do the job
0070 void DTTTrigWriter::analyze(const Event& event, const EventSetup& eventSetup) {
0071   if (debug)
0072     cout << "[DTTTrigWriter]Analyzer called!" << endl;
0073 
0074   // Get the DT Geometry
0075   dtGeom = eventSetup.getHandle(dtGeomToken_);
0076 
0077   // Get all the sls from the setup
0078   const vector<const DTSuperLayer*> superLayers = dtGeom->superLayers();
0079 
0080   // Loop over all SLs
0081   for (auto sl = superLayers.begin(); sl != superLayers.end(); sl++) {
0082     // Get the histo from file
0083     DTSuperLayerId slId = (*sl)->id();
0084     TH1F* histo = (TH1F*)theFile->Get((getTBoxName(slId)).c_str());
0085     if (histo) {  // Check that the histo exists
0086       // Compute mean and sigma of the rising edge
0087       pair<double, double> meanAndSigma = theFitter->fitTimeBox(histo);
0088 
0089       // Write them in DB object
0090       tTrig->set(slId, meanAndSigma.first, meanAndSigma.second, kFactor, DTTimeUnits::ns);
0091       if (debug) {
0092         cout << " SL: " << slId << " mean = " << meanAndSigma.first << " sigma = " << meanAndSigma.second << endl;
0093       }
0094     }
0095   }
0096 }
0097 
0098 // Write objects to DB
0099 void DTTTrigWriter::endJob() {
0100   if (debug)
0101     cout << "[DTTTrigWriter]Writing ttrig object to DB!" << endl;
0102 
0103   // FIXME: to be read from cfg?
0104   string tTrigRecord = "DTTtrigRcd";
0105 
0106   // Write the object to DB
0107   DTCalibDBUtils::writeToDB(tTrigRecord, *tTrig);
0108 }
0109 
0110 // Compute the name of the time box histo
0111 string DTTTrigWriter::getTBoxName(const DTSuperLayerId& slId) const {
0112   string histoName = "Ch_" + std::to_string(slId.wheel()) + "_" + std::to_string(slId.station()) + "_" +
0113                      std::to_string(slId.sector()) + "_SL" + std::to_string(slId.superlayer()) + "_hTimeBox";
0114   return histoName;
0115 }