Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 23:31:32

0001 /****************************************************************************
0002  * Authors:
0003  *   Jan Kašpar
0004  ****************************************************************************/
0005 
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0010 
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 
0014 #include "CondFormats/RunInfo/interface/LHCInfo.h"
0015 #include "CondFormats/DataRecord/interface/LHCInfoRcd.h"
0016 
0017 #include "TFile.h"
0018 #include "TH1D.h"
0019 #include "TH2D.h"
0020 
0021 //----------------------------------------------------------------------------------------------------
0022 
0023 class CTPPSLHCInfoPlotter : public edm::one::EDAnalyzer<> {
0024 public:
0025   explicit CTPPSLHCInfoPlotter(const edm::ParameterSet &);
0026 
0027   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0028 
0029 private:
0030   void analyze(const edm::Event &, const edm::EventSetup &) override;
0031   void endJob() override;
0032 
0033   edm::ESGetToken<LHCInfo, LHCInfoRcd> lhcInfoESToken_;
0034 
0035   std::string outputFile_;
0036 
0037   TH1D *h_beamEnergy_;
0038   TH1D *h_xangle_;
0039   TH1D *h_betaStar_;
0040   TH2D *h2_betaStar_vs_xangle_;
0041 
0042   TH1D *h_fill_;
0043   TH1D *h_run_;
0044 };
0045 
0046 //----------------------------------------------------------------------------------------------------
0047 
0048 using namespace std;
0049 using namespace edm;
0050 
0051 //----------------------------------------------------------------------------------------------------
0052 
0053 CTPPSLHCInfoPlotter::CTPPSLHCInfoPlotter(const edm::ParameterSet &iConfig)
0054     : lhcInfoESToken_(esConsumes(ESInputTag("", iConfig.getParameter<std::string>("lhcInfoLabel")))),
0055       outputFile_(iConfig.getParameter<string>("outputFile")),
0056 
0057       h_beamEnergy_(new TH1D("h_beamEnergy", ";beam energy   (GeV)", 81, -50., 8050.)),
0058       h_xangle_(new TH1D("h_xangle", ";(half) crossing angle   (#murad)", 201, -0.5, 200.5)),
0059       h_betaStar_(new TH1D("h_betaStar", ";#beta^{*}   (m)", 101, -0.005, 1.005)),
0060       h2_betaStar_vs_xangle_(new TH2D("h2_betaStar_vs_xangle",
0061                                       ";(half) crossing angle   (#murad);#beta^{*}   (m)",
0062                                       201,
0063                                       -0.5,
0064                                       200.5,
0065                                       101,
0066                                       -0.005,
0067                                       1.005)),
0068 
0069       h_fill_(new TH1D("h_fill", ";fill", 4001, 3999.5, 8000.5)),
0070       h_run_(new TH1D("h_run", ";run", 6000, 270E3, 330E3)) {}
0071 
0072 //----------------------------------------------------------------------------------------------------
0073 
0074 void CTPPSLHCInfoPlotter::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0075   edm::ParameterSetDescription desc;
0076 
0077   desc.add<std::string>("lhcInfoLabel", "")->setComment("label of the LHCInfo record");
0078   desc.add<std::string>("outputFile", "")->setComment("output file");
0079 
0080   descriptions.add("ctppsLHCInfoPlotter", desc);
0081 }
0082 
0083 //----------------------------------------------------------------------------------------------------
0084 
0085 void CTPPSLHCInfoPlotter::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0086   const auto &lhcInfo = iSetup.getData(lhcInfoESToken_);
0087 
0088   h_beamEnergy_->Fill(lhcInfo.energy());
0089   h_xangle_->Fill(lhcInfo.crossingAngle());
0090   h_betaStar_->Fill(lhcInfo.betaStar());
0091   h2_betaStar_vs_xangle_->Fill(lhcInfo.crossingAngle(), lhcInfo.betaStar());
0092 
0093   h_fill_->Fill(lhcInfo.fillNumber());
0094   h_run_->Fill(iEvent.id().run());
0095 }
0096 
0097 //----------------------------------------------------------------------------------------------------
0098 
0099 void CTPPSLHCInfoPlotter::endJob() {
0100   auto f_out = std::make_unique<TFile>(outputFile_.c_str(), "recreate");
0101 
0102   h_beamEnergy_->Write();
0103   h_xangle_->Write();
0104   h_betaStar_->Write();
0105   h2_betaStar_vs_xangle_->Write();
0106 
0107   h_fill_->Write();
0108   h_run_->Write();
0109 }
0110 
0111 //----------------------------------------------------------------------------------------------------
0112 
0113 DEFINE_FWK_MODULE(CTPPSLHCInfoPlotter);