File indexing completed on 2023-03-17 11:04:47
0001 #include <algorithm>
0002 #include <iostream>
0003 #include <iterator>
0004 #include <fstream>
0005 #include <string>
0006 #include <memory>
0007
0008 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/Run.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/Utilities/interface/InputTag.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include "SimDataFormats/GeneratorProducts/interface/LHERunInfoProduct.h"
0017 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0018
0019 using namespace lhef;
0020
0021 class LHEWriter : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0022 public:
0023 explicit LHEWriter(const edm::ParameterSet ¶ms);
0024 ~LHEWriter() override;
0025
0026 protected:
0027 void beginRun(const edm::Run &run, const edm::EventSetup &es) override;
0028 void endRun(const edm::Run &run, const edm::EventSetup &es) override;
0029 void analyze(const edm::Event &event, const edm::EventSetup &es) override;
0030
0031 private:
0032 std::ofstream file;
0033 std::ofstream file1;
0034
0035 edm::EDGetTokenT<LHERunInfoProduct> tokenLHERunInfo_;
0036 edm::EDGetTokenT<LHEEventProduct> tokenLHEEvent_;
0037 };
0038
0039 LHEWriter::LHEWriter(const edm::ParameterSet ¶ms)
0040 : tokenLHERunInfo_(consumes<LHERunInfoProduct, edm::InRun>(
0041 params.getUntrackedParameter<edm::InputTag>("moduleLabel", std::string("source")))),
0042 tokenLHEEvent_(consumes<LHEEventProduct>(
0043 params.getUntrackedParameter<edm::InputTag>("moduleLabel", std::string("source")))) {}
0044
0045 LHEWriter::~LHEWriter() {}
0046
0047 void LHEWriter::beginRun(const edm::Run &run, const edm::EventSetup &es) {
0048 file.open("writer.lhe", std::fstream::out | std::fstream::trunc);
0049 file1.open("writer1.lhe", std::fstream::out | std::fstream::trunc);
0050 }
0051
0052 void LHEWriter::endRun(const edm::Run &run, const edm::EventSetup &es) {
0053 const edm::Handle<LHERunInfoProduct> &product = run.getHandle(tokenLHERunInfo_);
0054
0055
0056
0057 std::copy(product->begin(), product->end(), std::ostream_iterator<std::string>(file));
0058
0059 file1 << LHERunInfoProduct::endOfFile();
0060 file.close();
0061 file1.close();
0062
0063 system("cat writer1.lhe >> writer.lhe");
0064 system("rm -rf writer1.lhe");
0065 }
0066
0067 void LHEWriter::analyze(const edm::Event &event, const edm::EventSetup &es) {
0068 const edm::Handle<LHEEventProduct> &product = event.getHandle(tokenLHEEvent_);
0069
0070
0071
0072 std::copy(product->begin(), product->end(), std::ostream_iterator<std::string>(file1));
0073 }
0074
0075 DEFINE_FWK_MODULE(LHEWriter);