File indexing completed on 2023-03-17 11:04:50
0001 #include <iostream>
0002 #include <cstring>
0003 #include <vector>
0004 #include <memory>
0005 #include <string>
0006
0007 #include <HepMC/GenEvent.h>
0008
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "GeneratorInterface/LHEInterface/interface/LHERunInfo.h"
0013 #include "GeneratorInterface/LHEInterface/interface/LHEEvent.h"
0014 #include "GeneratorInterface/PartonShowerVeto/interface/JetMatchingAlpgen.h"
0015
0016 namespace gen {
0017
0018 extern "C" {
0019
0020
0021
0022 extern void alsetp_();
0023
0024 extern void alveto_(int* ipveto);
0025
0026 extern void dbpart_();
0027 extern void pyupre_();
0028
0029 void alshcd_(char csho[3]);
0030 void alshen_();
0031 }
0032
0033
0034 JetMatchingAlpgen::JetMatchingAlpgen(const edm::ParameterSet& params)
0035 : JetMatching(params), applyMatching(params.getParameter<bool>("applyMatching")), runInitialized(false) {
0036 ahopts_.etclus = params.getParameter<double>("etMin");
0037 ahopts_.rclus = params.getParameter<double>("drMin");
0038 ahopts_.iexc = params.getParameter<bool>("exclusive");
0039 }
0040
0041
0042 JetMatchingAlpgen::~JetMatchingAlpgen() {}
0043
0044 std::set<std::string> JetMatchingAlpgen::capabilities() const {
0045 std::set<std::string> result;
0046 result.insert("psFinalState");
0047 result.insert("hepevt");
0048 result.insert("pythia6");
0049 return result;
0050 }
0051
0052
0053
0054 void JetMatchingAlpgen::init(const lhef::LHERunInfo* runInfo) {
0055
0056 std::vector<std::string> headerLines = runInfo->findHeader("AlpgenUnwParFile");
0057 if (headerLines.empty())
0058 throw cms::Exception("Generator|PartonShowerVeto") << "In order to use Alpgen jet matching, "
0059 "the input file has to contain the corresponding "
0060 "Alpgen headers."
0061 << std::endl;
0062
0063
0064 header.parse(headerLines.begin(), headerLines.end());
0065
0066
0067
0068
0069
0070
0071
0072
0073
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094 std::copy(header.masses, header.masses + AlpgenHeader::MASS_MAX, ahppara_.masses);
0095 ahppara_.ihrd = header.ihrd;
0096 ahppara_.ebeam = header.ebeam;
0097
0098 for (std::map<AlpgenHeader::Parameter, double>::const_iterator iter = header.params.begin();
0099 iter != header.params.end();
0100 ++iter) {
0101 if (iter->first <= 0 || iter->first >= (int)AHPARS::nparam - 1)
0102 continue;
0103 ahpars_.parval[(int)iter->first - 1] = iter->second;
0104 }
0105
0106
0107 alsetp_();
0108
0109
0110 runInitialized = true;
0111 }
0112
0113 void JetMatchingAlpgen::beforeHadronisation(const lhef::LHEEvent* event) {
0114
0115 if (!runInitialized)
0116 throw cms::Exception("Generator|PartonShowerVeto") << "Run not initialized in JetMatchingAlpgen" << std::endl;
0117
0118
0119
0120
0121
0122
0123
0124
0125 eventInitialized = true;
0126 }
0127
0128
0129
0130
0131
0132
0133 int JetMatchingAlpgen::match(const lhef::LHEEvent* partonLevel, const std::vector<fastjet::PseudoJet>* jetInput) {
0134
0135
0136
0137
0138
0139
0140
0141 if (!runInitialized)
0142 throw cms::Exception("Generator|PartonShowerVeto") << "Run not initialized in JetMatchingAlpgen" << std::endl;
0143
0144 if (!eventInitialized)
0145 throw cms::Exception("Generator|PartonShowerVeto") << "Event not initialized in JetMatchingAlpgen" << std::endl;
0146
0147
0148
0149 if (!applyMatching)
0150 return 0;
0151
0152
0153 int veto = 0;
0154 alveto_(&veto);
0155
0156 eventInitialized = false;
0157
0158
0159
0160
0161
0162 return veto ? 1 : 0;
0163 }
0164
0165 void alshcd_(char csho[3]) {
0166 std::memcpy(csho, "PYT", 3);
0167 }
0168
0169 void alshen_() {}
0170
0171 }