DTCalibration/plugins/DTVDriftSegmentCalibration.cc

0001
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \author A. Vilela Pereira
0006  */
0007
0008 #include "DTVDriftSegmentCalibration.h"
0009
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/ESHandle.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016
0017 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0018 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0019
0020 #include "CalibMuon/DTCalibration/interface/DTCalibDBUtils.h"
0021 #include "CalibMuon/DTCalibration/interface/DTSegmentSelector.h"
0022
0023 #include <string>
0024 #include <sstream>
0025 #include "TFile.h"
0026 #include "TH1F.h"
0027 #include "TH2F.h"
0028
0029 using namespace std;
0030 using namespace edm;
0031
0032 DTVDriftSegmentCalibration::DTVDriftSegmentCalibration(const ParameterSet& pset)
0033     : theRecHits4DToken_(consumes<DTRecSegment4DCollection>(pset.getParameter<InputTag>("recHits4DLabel"))),
0034       //writeVDriftDB_(pset.getUntrackedParameter<bool>("writeVDriftDB", false)),
0035       theCalibChamber_(pset.getUntrackedParameter<string>("calibChamber", "All")),
0036       dtGeomToken_(esConsumes()) {
0037   LogVerbatim("Calibration") << "[DTVDriftSegmentCalibration] Constructor called!";
0038
0039   edm::ConsumesCollector collector(consumesCollector());
0040   select_ = new DTSegmentSelector(pset, collector);
0041   // the root file which will contain the histos
0042   string rootFileName = pset.getUntrackedParameter<string>("rootFileName", "DTVDriftHistos.root");
0043   rootFile_ = new TFile(rootFileName.c_str(), "RECREATE");
0044   rootFile_->cd();
0045 }
0046
0047 void DTVDriftSegmentCalibration::beginJob() { TH1::SetDefaultSumw2(true); }
0048
0049 void DTVDriftSegmentCalibration::beginRun(const edm::Run& run, const edm::EventSetup& setup) {}
0050
0051 DTVDriftSegmentCalibration::~DTVDriftSegmentCalibration() {
0052   rootFile_->Close();
0053   LogVerbatim("Calibration") << "[DTVDriftSegmentCalibration] Destructor called!";
0054 }
0055
0056 void DTVDriftSegmentCalibration::analyze(const Event& event, const EventSetup& eventSetup) {
0057   rootFile_->cd();
0058
0059   // Get the DT Geometry
0060   ESHandle<DTGeometry> dtGeom;
0061   dtGeom = eventSetup.getHandle(dtGeomToken_);
0062
0063   // Get the rechit collection from the event
0064   const Handle<DTRecSegment4DCollection>& all4DSegments = event.getHandle(theRecHits4DToken_);
0065
0066   DTChamberId chosenChamberId;
0067   if (theCalibChamber_ != "All") {
0068     stringstream linestr;
0069     int selWheel, selStation, selSector;
0070     linestr << theCalibChamber_;
0071     linestr >> selWheel >> selStation >> selSector;
0072     chosenChamberId = DTChamberId(selWheel, selStation, selSector);
0073     LogVerbatim("Calibration") << " Chosen chamber: " << chosenChamberId << endl;
0074   }
0075   // Loop over segments by chamber
0076   DTRecSegment4DCollection::id_iterator chamberIdIt;
0077   for (chamberIdIt = all4DSegments->id_begin(); chamberIdIt != all4DSegments->id_end(); ++chamberIdIt) {
0078     // Calibrate just the chosen chamber/s
0079     if ((theCalibChamber_ != "All") && ((*chamberIdIt) != chosenChamberId))
0080       continue;
0081
0082     // Book histos
0083     if (theVDriftHistoMapTH1F_.find(*chamberIdIt) == theVDriftHistoMapTH1F_.end()) {
0084       LogTrace("Calibration") << "   Booking histos for Chamber: " << *chamberIdIt;
0085       bookHistos(*chamberIdIt);
0086     }
0087
0088     // Get the chamber from the setup
0089     const DTChamber* chamber = dtGeom->chamber(*chamberIdIt);
0090     // Get the range for the corresponding ChamberId
0091     DTRecSegment4DCollection::range range = all4DSegments->get((*chamberIdIt));
0092     // Loop over the rechits of this DetUnit
0093     for (DTRecSegment4DCollection::const_iterator segment = range.first; segment != range.second; ++segment) {
0094       LogTrace("Calibration") << "Segment local pos (in chamber RF): " << (*segment).localPosition()
0095                               << "\nSegment global pos: " << chamber->toGlobal((*segment).localPosition());
0096
0097       if (!(*select_)(*segment, event, eventSetup))
0098         continue;
0099
0100       // Fill v-drift values
0101       if ((*segment).hasPhi()) {
0102         //if( segment->phiSegment()->ist0Valid() ){
0103         double segmentVDrift = segment->phiSegment()->vDrift();
0104         if (segmentVDrift != 0.00) {
0105           (theVDriftHistoMapTH1F_[*chamberIdIt])[0]->Fill(segmentVDrift);
0106           (theVDriftHistoMapTH2F_[*chamberIdIt])[0]->Fill(segment->localPosition().x(), segmentVDrift);
0107           (theVDriftHistoMapTH2F_[*chamberIdIt])[1]->Fill(segment->localPosition().y(), segmentVDrift);
0108         }
0109       }
0110       // Probably not meaningful
0111       if ((*segment).hasZed()) {
0112         //if( segment->zSegment()->ist0Valid() ){
0113         double segmentVDrift = segment->zSegment()->vDrift();
0114         if (segmentVDrift != 0.00) {
0115           (theVDriftHistoMapTH1F_[*chamberIdIt])[1]->Fill(segmentVDrift);
0116         }
0117       }
0118     }  // DTRecSegment4DCollection::const_iterator segment
0119   }  // DTRecSegment4DCollection::id_iterator chamberIdIt
0120 }  // DTVDriftSegmentCalibration::analyze
0121
0122 void DTVDriftSegmentCalibration::endJob() {
0123   rootFile_->cd();
0124
0125   LogVerbatim("Calibration") << "[DTVDriftSegmentCalibration] Writing histos to file!" << endl;
0126
0127   for (ChamberHistosMapTH1F::const_iterator itChHistos = theVDriftHistoMapTH1F_.begin();
0128        itChHistos != theVDriftHistoMapTH1F_.end();
0129        ++itChHistos) {
0130     vector<TH1F*>::const_iterator itHistTH1F = (*itChHistos).second.begin();
0131     vector<TH1F*>::const_iterator itHistTH1F_end = (*itChHistos).second.end();
0132     for (; itHistTH1F != itHistTH1F_end; ++itHistTH1F)
0133       (*itHistTH1F)->Write();
0134
0135     vector<TH2F*>::const_iterator itHistTH2F = theVDriftHistoMapTH2F_[(*itChHistos).first].begin();
0136     vector<TH2F*>::const_iterator itHistTH2F_end = theVDriftHistoMapTH2F_[(*itChHistos).first].end();
0137     for (; itHistTH2F != itHistTH2F_end; ++itHistTH2F)
0138       (*itHistTH2F)->Write();
0139   }
0140
0141   /*if(writeVDriftDB_){
0142      // ...
0143   }*/
0144 }
0145
0146 // Book a set of histograms for a given Chamber
0147 void DTVDriftSegmentCalibration::bookHistos(DTChamberId chId) {
0148   // Compose the chamber name
0149   std::string wheel = std::to_string(chId.wheel());
0150   std::string station = std::to_string(chId.station());
0151   std::string sector = std::to_string(chId.sector());
0152
0153   string chHistoName = "_W" + wheel + "_St" + station + "_Sec" + sector;
0154
0155   vector<TH1F*> histosTH1F;
0156   histosTH1F.push_back(
0157       new TH1F(("hRPhiVDriftCorr" + chHistoName).c_str(), "v-drift corr. from Phi segments", 200, -0.4, 0.4));
0158   if (chId.station() != 4)
0159     histosTH1F.push_back(
0160         new TH1F(("hRZVDriftCorr" + chHistoName).c_str(), "v-drift corr. from Z segments", 200, -0.4, 0.4));
0161
0162   vector<TH2F*> histosTH2F;
0163   histosTH2F.push_back(new TH2F(("hRPhiVDriftCorrVsSegmPosX" + chHistoName).c_str(),
0164                                 "v-drift corr. vs. segment x position",
0165                                 250,
0166                                 -125.,
0167                                 125.,
0168                                 200,
0169                                 -0.4,
0170                                 0.4));
0171   histosTH2F.push_back(new TH2F(("hRPhiVDriftCorrVsSegmPosY" + chHistoName).c_str(),
0172                                 "v-drift corr. vs. segment y position",
0173                                 250,
0174                                 -125.,
0175                                 125.,
0176                                 200,
0177                                 -0.4,
0178                                 0.4));
0179   //if(chId.station() != 4) ...
0180
0181   theVDriftHistoMapTH1F_[chId] = histosTH1F;
0182   theVDriftHistoMapTH2F_[chId] = histosTH2F;
0183 }