File indexing completed on 2024-05-02 05:10:04
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include "TLorentzVector.h"
0021 #include <algorithm>
0022 #include <fstream>
0023 #include <iostream>
0024 #include <iterator>
0025 #include <memory>
0026 #include <string>
0027
0028
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/one/EDProducer.h"
0031
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/Framework/interface/Run.h"
0035
0036 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0037
0038 #include "DataFormats/Math/interface/LorentzVector.h"
0039 #include "DataFormats/PatCandidates/interface/Muon.h"
0040 #include "DataFormats/VertexReco/interface/Vertex.h"
0041 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0042
0043 #include "SimDataFormats/GeneratorProducts/interface/LHECommonBlocks.h"
0044 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0045 #include "SimDataFormats/GeneratorProducts/interface/LHERunInfoProduct.h"
0046 #include "SimDataFormats/GeneratorProducts/interface/LHEXMLStringProduct.h"
0047 #include "SimDataFormats/GeneratorProducts/interface/LesHouches.h"
0048
0049 #include "GeneratorInterface/LHEInterface/interface/LHEEvent.h"
0050 #include "GeneratorInterface/LHEInterface/interface/LHEReader.h"
0051 #include "GeneratorInterface/LHEInterface/interface/LHERunInfo.h"
0052
0053 #include "CLHEP/Random/RandExponential.h"
0054 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0055 #include "FWCore/ServiceRegistry/interface/Service.h"
0056 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
0057 #include "FWCore/Utilities/interface/StreamID.h"
0058
0059 namespace CLHEP {
0060 class HepRandomEngine;
0061 }
0062
0063 class EmbeddingLHEProducer : public edm::one::EDProducer<edm::BeginRunProducer, edm::EndRunProducer> {
0064 public:
0065 explicit EmbeddingLHEProducer(const edm::ParameterSet &);
0066 ~EmbeddingLHEProducer() override;
0067
0068 static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0069
0070 private:
0071 void beginJob() override;
0072 void produce(edm::Event &, const edm::EventSetup &) override;
0073 void endJob() override;
0074
0075 void beginRunProduce(edm::Run &run, edm::EventSetup const &es) override;
0076 void endRunProduce(edm::Run &, edm::EventSetup const &) override;
0077
0078 void fill_lhe_from_mumu(TLorentzVector &positiveLepton,
0079 TLorentzVector &negativeLepton,
0080 lhef::HEPEUP &outlhe,
0081 CLHEP::HepRandomEngine *engine);
0082 void fill_lhe_with_particle(lhef::HEPEUP &outlhe, TLorentzVector &particle, int pdgid, double spin, double ctau);
0083
0084 void transform_mumu_to_tautau(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton);
0085 const reco::Candidate *find_original_muon(const reco::Candidate *muon);
0086 void assign_4vector(TLorentzVector &Lepton, const pat::Muon *muon, std::string FSRmode);
0087 void mirror(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton);
0088 void InitialRecoCorrection(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton);
0089 void rotate180(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton);
0090
0091 LHERunInfoProduct::Header give_slha();
0092
0093 edm::EDGetTokenT<edm::View<pat::Muon>> muonsCollection_;
0094 edm::EDGetTokenT<reco::VertexCollection> vertexCollection_;
0095 int particleToEmbed_;
0096 bool mirror_, rotate180_, InitialRecoCorrection_;
0097 static constexpr double tauMass_ = 1.77682;
0098 static constexpr double muonMass_ = 0.1057;
0099 static constexpr double elMass_ = 0.00051;
0100 const int embeddingParticles[3]{11, 13, 15};
0101
0102 std::ofstream file;
0103 bool write_lheout;
0104
0105
0106
0107
0108 std::string studyFSRmode_;
0109 };
0110
0111
0112
0113
0114 EmbeddingLHEProducer::EmbeddingLHEProducer(const edm::ParameterSet &iConfig) {
0115
0116 produces<LHEEventProduct>();
0117 produces<LHERunInfoProduct, edm::Transition::BeginRun>();
0118 produces<math::XYZTLorentzVectorD>("vertexPosition");
0119
0120 muonsCollection_ = consumes<edm::View<pat::Muon>>(iConfig.getParameter<edm::InputTag>("src"));
0121 vertexCollection_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertices"));
0122 particleToEmbed_ = iConfig.getParameter<int>("particleToEmbed");
0123 mirror_ = iConfig.getParameter<bool>("mirror");
0124 InitialRecoCorrection_ = iConfig.getParameter<bool>("InitialRecoCorrection");
0125 rotate180_ = iConfig.getParameter<bool>("rotate180");
0126 studyFSRmode_ = iConfig.getUntrackedParameter<std::string>("studyFSRmode", "");
0127
0128 write_lheout = false;
0129 std::string lhe_ouputfile = iConfig.getUntrackedParameter<std::string>("lhe_outputfilename", "");
0130 if (!lhe_ouputfile.empty()) {
0131 write_lheout = true;
0132 file.open(lhe_ouputfile, std::fstream::out | std::fstream::trunc);
0133 }
0134
0135
0136 if (std::find(std::begin(embeddingParticles), std::end(embeddingParticles), particleToEmbed_) ==
0137 std::end(embeddingParticles)) {
0138 throw cms::Exception("Configuration") << "The given particle to embed is not in the list of allowed particles.";
0139 }
0140
0141 edm::Service<edm::RandomNumberGenerator> rng;
0142 if (!rng.isAvailable()) {
0143 throw cms::Exception("Configuration") << "The EmbeddingLHEProducer requires the RandomNumberGeneratorService\n"
0144 "which is not present in the configuration file. \n"
0145 "You must add the service\n"
0146 "in the configuration file or remove the modules that require it.";
0147 }
0148 }
0149
0150 EmbeddingLHEProducer::~EmbeddingLHEProducer() {}
0151
0152
0153
0154
0155
0156
0157 void EmbeddingLHEProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0158 using namespace edm;
0159
0160 edm::Service<edm::RandomNumberGenerator> rng;
0161 CLHEP::HepRandomEngine *engine = &rng->getEngine(iEvent.streamID());
0162
0163 edm::Handle<edm::View<pat::Muon>> muonHandle;
0164 iEvent.getByToken(muonsCollection_, muonHandle);
0165 edm::View<pat::Muon> coll_muons = *muonHandle;
0166
0167 Handle<std::vector<reco::Vertex>> coll_vertices;
0168 iEvent.getByToken(vertexCollection_, coll_vertices);
0169
0170 TLorentzVector positiveLepton, negativeLepton;
0171 bool mu_plus_found = false;
0172 bool mu_minus_found = false;
0173 lhef::HEPEUP hepeup;
0174 hepeup.IDPRUP = 0;
0175 hepeup.XWGTUP = 1;
0176 hepeup.SCALUP = -1;
0177 hepeup.AQEDUP = -1;
0178 hepeup.AQCDUP = -1;
0179
0180 for (edm::View<pat::Muon>::const_iterator muon = coll_muons.begin(); muon != coll_muons.end(); ++muon) {
0181 if (muon->charge() == 1 && !mu_plus_found) {
0182 assign_4vector(positiveLepton, &(*muon), studyFSRmode_);
0183 mu_plus_found = true;
0184 } else if (muon->charge() == -1 && !mu_minus_found) {
0185 assign_4vector(negativeLepton, &(*muon), studyFSRmode_);
0186 mu_minus_found = true;
0187 } else if (mu_minus_found && mu_plus_found)
0188 break;
0189 }
0190 InitialRecoCorrection(
0191 positiveLepton,
0192 negativeLepton);
0193 mirror(positiveLepton, negativeLepton);
0194 rotate180(positiveLepton, negativeLepton);
0195 transform_mumu_to_tautau(positiveLepton, negativeLepton);
0196 fill_lhe_from_mumu(positiveLepton, negativeLepton, hepeup, engine);
0197
0198 double originalXWGTUP_ = 1.;
0199 std::unique_ptr<LHEEventProduct> product(new LHEEventProduct(hepeup, originalXWGTUP_));
0200
0201 if (write_lheout)
0202 std::copy(product->begin(), product->end(), std::ostream_iterator<std::string>(file));
0203
0204 iEvent.put(std::move(product));
0205
0206 std::unique_ptr<math::XYZTLorentzVectorD> vertex_position(
0207 new math::XYZTLorentzVectorD(coll_vertices->at(0).x(), coll_vertices->at(0).y(), coll_vertices->at(0).z(), 0.0));
0208 iEvent.put(std::move(vertex_position), "vertexPosition");
0209 }
0210
0211
0212 void EmbeddingLHEProducer::beginJob() {}
0213
0214
0215 void EmbeddingLHEProducer::endJob() {}
0216
0217
0218
0219 void EmbeddingLHEProducer::beginRunProduce(edm::Run &run, edm::EventSetup const &) {
0220
0221 lhef::HEPRUP heprup;
0222
0223
0224 heprup.resize(1);
0225
0226
0227 heprup.IDWTUP = 3;
0228
0229
0230 heprup.XSECUP[0] = 1.;
0231 heprup.XERRUP[0] = 0;
0232 heprup.XMAXUP[0] = 1;
0233 heprup.LPRUP[0] = 1;
0234
0235 std::unique_ptr<LHERunInfoProduct> runInfo(new LHERunInfoProduct(heprup));
0236 runInfo->addHeader(give_slha());
0237
0238 if (write_lheout)
0239 std::copy(runInfo->begin(), runInfo->end(), std::ostream_iterator<std::string>(file));
0240 run.put(std::move(runInfo));
0241 }
0242
0243 void EmbeddingLHEProducer::endRunProduce(edm::Run &run, edm::EventSetup const &es) {
0244 if (write_lheout) {
0245 file << LHERunInfoProduct::endOfFile();
0246 file.close();
0247 }
0248 }
0249
0250 void EmbeddingLHEProducer::fill_lhe_from_mumu(TLorentzVector &positiveLepton,
0251 TLorentzVector &negativeLepton,
0252 lhef::HEPEUP &outlhe,
0253 CLHEP::HepRandomEngine *engine) {
0254 TLorentzVector Z = positiveLepton + negativeLepton;
0255 int leptonPDGID = particleToEmbed_;
0256
0257 static constexpr double tau_ctau0 = 8.71100e-02;
0258 double tau_ctau_p = tau_ctau0 * CLHEP::RandExponential::shoot(engine);
0259 double tau_ctau_n = tau_ctau0 * CLHEP::RandExponential::shoot(engine);
0260
0261 fill_lhe_with_particle(outlhe, Z, 23, 9.0, 0);
0262 fill_lhe_with_particle(outlhe, positiveLepton, -leptonPDGID, 1.0, tau_ctau_p);
0263 fill_lhe_with_particle(outlhe, negativeLepton, leptonPDGID, -1.0, tau_ctau_n);
0264
0265 return;
0266 }
0267
0268 void EmbeddingLHEProducer::fill_lhe_with_particle(
0269 lhef::HEPEUP &outlhe, TLorentzVector &particle, int pdgid, double spin, double ctau) {
0270
0271
0272
0273
0274 int particleindex = outlhe.NUP;
0275 outlhe.resize(outlhe.NUP + 1);
0276
0277 outlhe.PUP[particleindex][0] = particle.Px();
0278 outlhe.PUP[particleindex][1] = particle.Py();
0279 outlhe.PUP[particleindex][2] = particle.Pz();
0280 outlhe.PUP[particleindex][3] = particle.E();
0281 outlhe.PUP[particleindex][4] = particle.M();
0282 outlhe.IDUP[particleindex] = pdgid;
0283 outlhe.SPINUP[particleindex] = spin;
0284 outlhe.VTIMUP[particleindex] = ctau;
0285
0286 outlhe.ICOLUP[particleindex].first = 0;
0287 outlhe.ICOLUP[particleindex].second = 0;
0288
0289 if (std::abs(pdgid) == 23) {
0290 outlhe.MOTHUP[particleindex].first = 0;
0291 outlhe.MOTHUP[particleindex].second = 0;
0292 outlhe.ISTUP[particleindex] = 2;
0293 }
0294
0295 if (std::find(std::begin(embeddingParticles), std::end(embeddingParticles), std::abs(pdgid)) !=
0296 std::end(embeddingParticles)) {
0297 outlhe.MOTHUP[particleindex].first = 1;
0298 outlhe.MOTHUP[particleindex].second = 1;
0299
0300 outlhe.ISTUP[particleindex] = 1;
0301 }
0302
0303 return;
0304 }
0305
0306 void EmbeddingLHEProducer::transform_mumu_to_tautau(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton) {
0307
0308 double lep_mass;
0309 if (particleToEmbed_ == 11) {
0310 lep_mass = elMass_;
0311 } else if (particleToEmbed_ == 15) {
0312 lep_mass = tauMass_;
0313 } else {
0314 return;
0315 }
0316
0317 TLorentzVector Z = positiveLepton + negativeLepton;
0318
0319 TVector3 boost_from_Z_to_LAB = Z.BoostVector();
0320 TVector3 boost_from_LAB_to_Z = -Z.BoostVector();
0321
0322 positiveLepton.Boost(boost_from_LAB_to_Z);
0323 negativeLepton.Boost(boost_from_LAB_to_Z);
0324
0325
0326 double lep_mass_squared = lep_mass * lep_mass;
0327 double lep_energy_squared = 0.25 * Z.M2();
0328 double lep_3momentum_squared = lep_energy_squared - lep_mass_squared;
0329 if (lep_3momentum_squared < 0) {
0330 edm::LogWarning("TauEmbedding") << "3-Momentum squared is negative";
0331 return;
0332 }
0333
0334
0335 double scale = std::sqrt(lep_3momentum_squared / positiveLepton.Vect().Mag2());
0336 positiveLepton.SetPxPyPzE(scale * positiveLepton.Px(),
0337 scale * positiveLepton.Py(),
0338 scale * positiveLepton.Pz(),
0339 std::sqrt(lep_energy_squared));
0340 negativeLepton.SetPxPyPzE(scale * negativeLepton.Px(),
0341 scale * negativeLepton.Py(),
0342 scale * negativeLepton.Pz(),
0343 std::sqrt(lep_energy_squared));
0344
0345
0346 positiveLepton.Boost(boost_from_Z_to_LAB);
0347 negativeLepton.Boost(boost_from_Z_to_LAB);
0348
0349 return;
0350 }
0351
0352 void EmbeddingLHEProducer::assign_4vector(TLorentzVector &Lepton, const pat::Muon *muon, std::string FSRmode) {
0353 if ("afterFSR" == FSRmode && muon->genParticle() != nullptr) {
0354 const reco::GenParticle *afterFSRMuon = muon->genParticle();
0355 Lepton.SetPxPyPzE(
0356 afterFSRMuon->p4().px(), afterFSRMuon->p4().py(), afterFSRMuon->p4().pz(), afterFSRMuon->p4().e());
0357 } else if ("beforeFSR" == FSRmode && muon->genParticle() != nullptr) {
0358 const reco::Candidate *beforeFSRMuon = find_original_muon(muon->genParticle());
0359 Lepton.SetPxPyPzE(
0360 beforeFSRMuon->p4().px(), beforeFSRMuon->p4().py(), beforeFSRMuon->p4().pz(), beforeFSRMuon->p4().e());
0361 } else {
0362 Lepton.SetPxPyPzE(muon->p4().px(), muon->p4().py(), muon->p4().pz(), muon->p4().e());
0363 }
0364 return;
0365 }
0366
0367 const reco::Candidate *EmbeddingLHEProducer::find_original_muon(const reco::Candidate *muon) {
0368 if (muon->mother(0) == nullptr)
0369 return muon;
0370 if (muon->pdgId() == muon->mother(0)->pdgId())
0371 return find_original_muon(muon->mother(0));
0372 else
0373 return muon;
0374 }
0375
0376 void EmbeddingLHEProducer::rotate180(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton) {
0377 if (!rotate180_)
0378 return;
0379 edm::LogInfo("TauEmbedding") << "Applying 180<C2><B0> rotation";
0380
0381
0382 TLorentzVector Z = positiveLepton + negativeLepton;
0383
0384 edm::LogInfo("TauEmbedding") << "MuMinus before. Pt: " << negativeLepton.Pt() << " Eta: " << negativeLepton.Eta()
0385 << " Phi: " << negativeLepton.Phi() << " Mass: " << negativeLepton.M();
0386
0387 TVector3 Z3 = Z.Vect();
0388 TVector3 positiveLepton3 = positiveLepton.Vect();
0389 TVector3 negativeLepton3 = negativeLepton.Vect();
0390
0391 TVector3 p3_perp = positiveLepton3 - positiveLepton3.Dot(Z3) / Z3.Dot(Z3) * Z3;
0392 p3_perp = p3_perp.Unit();
0393
0394 positiveLepton3 -= 2 * positiveLepton3.Dot(p3_perp) * p3_perp;
0395 negativeLepton3 -= 2 * negativeLepton3.Dot(p3_perp) * p3_perp;
0396
0397 positiveLepton.SetVect(positiveLepton3);
0398 negativeLepton.SetVect(negativeLepton3);
0399
0400 edm::LogInfo("TauEmbedding") << "MuMinus after. Pt: " << negativeLepton.Pt() << " Eta: " << negativeLepton.Eta()
0401 << " Phi: " << negativeLepton.Phi() << " Mass: " << negativeLepton.M();
0402
0403 return;
0404 }
0405
0406 void EmbeddingLHEProducer::InitialRecoCorrection(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton) {
0407 if (!InitialRecoCorrection_)
0408 return;
0409 edm::LogInfo("TauEmbedding") << "Applying initial reconstruction correction";
0410 TLorentzVector Z = positiveLepton + negativeLepton;
0411 edm::LogInfo("TauEmbedding") << "MuMinus before. Pt: " << negativeLepton.Pt() << " Mass: " << negativeLepton.M();
0412 float diLeptonMass = (positiveLepton + negativeLepton).M();
0413 if (diLeptonMass > 60. && diLeptonMass < 122.) {
0414 static constexpr float zmass = 91.1876;
0415 float correction_deviation =
0416 5.;
0417 double EmbeddingCorrection =
0418 1.138;
0419 EmbeddingCorrection /=
0420 (EmbeddingCorrection -
0421 (EmbeddingCorrection - 1.) * exp(-pow((diLeptonMass - zmass), 2.) / (2. * pow(correction_deviation, 2.))));
0422 EmbeddingCorrection = ((diLeptonMass + (EmbeddingCorrection - 1.) * zmass) / (diLeptonMass * EmbeddingCorrection));
0423 double correctedpositiveLeptonEnergy =
0424 std::sqrt((pow(muonMass_, 2) / pow(EmbeddingCorrection, 2)) + pow(positiveLepton.Px(), 2) +
0425 pow(positiveLepton.Py(), 2) + pow(positiveLepton.Pz(), 2));
0426 double correctednegativeLeptonEnergy =
0427 std::sqrt((pow(muonMass_, 2) / pow(EmbeddingCorrection, 2)) + pow(negativeLepton.Px(), 2) +
0428 pow(negativeLepton.Py(), 2) + pow(negativeLepton.Pz(), 2));
0429 positiveLepton.SetE(correctedpositiveLeptonEnergy);
0430 negativeLepton.SetE(correctednegativeLeptonEnergy);
0431 positiveLepton *= EmbeddingCorrection;
0432 negativeLepton *= EmbeddingCorrection;
0433 edm::LogInfo("TauEmbedding") << "MuMinus after. Pt: " << negativeLepton.Pt() << " Mass: " << negativeLepton.M();
0434 }
0435 return;
0436 }
0437
0438 void EmbeddingLHEProducer::mirror(TLorentzVector &positiveLepton, TLorentzVector &negativeLepton) {
0439 if (!mirror_)
0440 return;
0441 edm::LogInfo("TauEmbedding") << "Applying mirroring";
0442 TLorentzVector Z = positiveLepton + negativeLepton;
0443
0444 edm::LogInfo("TauEmbedding") << "MuMinus before. Pt: " << negativeLepton.Pt() << " Eta: " << negativeLepton.Eta()
0445 << " Phi: " << negativeLepton.Phi() << " Mass: " << negativeLepton.M();
0446
0447 TVector3 Z3 = Z.Vect();
0448 TVector3 positiveLepton3 = positiveLepton.Vect();
0449 TVector3 negativeLepton3 = negativeLepton.Vect();
0450
0451 TVector3 beam(0., 0., 1.);
0452 TVector3 perpToZandBeam = Z3.Cross(beam).Unit();
0453
0454 positiveLepton3 -= 2 * positiveLepton3.Dot(perpToZandBeam) * perpToZandBeam;
0455 negativeLepton3 -= 2 * negativeLepton3.Dot(perpToZandBeam) * perpToZandBeam;
0456
0457 positiveLepton.SetVect(positiveLepton3);
0458 negativeLepton.SetVect(negativeLepton3);
0459
0460 edm::LogInfo("TauEmbedding") << "MuMinus after. Pt: " << negativeLepton.Pt() << " Eta: " << negativeLepton.Eta()
0461 << " Phi: " << negativeLepton.Phi() << " Mass: " << negativeLepton.M();
0462
0463 return;
0464 }
0465
0466 LHERunInfoProduct::Header EmbeddingLHEProducer::give_slha() {
0467 LHERunInfoProduct::Header slhah("slha");
0468
0469 slhah.addLine("######################################################################\n");
0470 slhah.addLine("## PARAM_CARD AUTOMATICALY GENERATED BY MG5 FOLLOWING UFO MODEL ####\n");
0471 slhah.addLine("######################################################################\n");
0472 slhah.addLine("## ##\n");
0473 slhah.addLine("## Width set on Auto will be computed following the information ##\n");
0474 slhah.addLine("## present in the decay.py files of the model. ##\n");
0475 slhah.addLine("## See arXiv:1402.1178 for more details. ##\n");
0476 slhah.addLine("## ##\n");
0477 slhah.addLine("######################################################################\n");
0478 slhah.addLine("\n");
0479 slhah.addLine("###################################\n");
0480 slhah.addLine("## INFORMATION FOR MASS\n");
0481 slhah.addLine("###################################\n");
0482 slhah.addLine("Block mass \n");
0483 slhah.addLine(" 6 1.730000e+02 # MT \n");
0484 slhah.addLine(" 15 1.777000e+00 # MTA \n");
0485 slhah.addLine(" 23 9.118800e+01 # MZ \n");
0486 slhah.addLine(" 25 1.250000e+02 # MH \n");
0487 slhah.addLine("## Dependent parameters, given by model restrictions.\n");
0488 slhah.addLine("## Those values should be edited following the \n");
0489 slhah.addLine("## analytical expression. MG5 ignores those values \n");
0490 slhah.addLine("## but they are important for interfacing the output of MG5\n");
0491 slhah.addLine("## to external program such as Pythia.\n");
0492 slhah.addLine(" 1 0.000000 # d : 0.0 \n");
0493 slhah.addLine(" 2 0.000000 # u : 0.0 \n");
0494 slhah.addLine(" 3 0.000000 # s : 0.0 \n");
0495 slhah.addLine(" 4 0.000000 # c : 0.0 \n");
0496 slhah.addLine(" 5 0.000000 # b : 0.0 \n");
0497 slhah.addLine(" 11 0.000000 # e- : 0.0 \n");
0498 slhah.addLine(" 12 0.000000 # ve : 0.0 \n");
0499 slhah.addLine(" 13 0.000000 # mu- : 0.0 \n");
0500 slhah.addLine(" 14 0.000000 # vm : 0.0 \n");
0501 slhah.addLine(" 16 0.000000 # vt : 0.0 \n");
0502 slhah.addLine(" 21 0.000000 # g : 0.0 \n");
0503 slhah.addLine(" 22 0.000000 # a : 0.0 \n");
0504 slhah.addLine(
0505 " 24 80.419002 # w+ : cmath.sqrt(MZ__exp__2/2. + cmath.sqrt(MZ__exp__4/4. - "
0506 "(aEW*cmath.pi*MZ__exp__2)/(Gf*sqrt__2))) \n");
0507 slhah.addLine("\n");
0508 slhah.addLine("###################################\n");
0509 slhah.addLine("## INFORMATION FOR SMINPUTS\n");
0510 slhah.addLine("###################################\n");
0511 slhah.addLine("Block sminputs \n");
0512 slhah.addLine(" 1 1.325070e+02 # aEWM1 \n");
0513 slhah.addLine(" 2 1.166390e-05 # Gf \n");
0514 slhah.addLine(" 3 1.180000e-01 # aS \n");
0515 slhah.addLine("\n");
0516 slhah.addLine("###################################\n");
0517 slhah.addLine("## INFORMATION FOR WOLFENSTEIN\n");
0518 slhah.addLine("###################################\n");
0519 slhah.addLine("Block wolfenstein \n");
0520 slhah.addLine(" 1 2.253000e-01 # lamWS \n");
0521 slhah.addLine(" 2 8.080000e-01 # AWS \n");
0522 slhah.addLine(" 3 1.320000e-01 # rhoWS \n");
0523 slhah.addLine(" 4 3.410000e-01 # etaWS \n");
0524 slhah.addLine("\n");
0525 slhah.addLine("###################################\n");
0526 slhah.addLine("## INFORMATION FOR YUKAWA\n");
0527 slhah.addLine("###################################\n");
0528 slhah.addLine("Block yukawa \n");
0529 slhah.addLine(" 6 1.730000e+02 # ymt \n");
0530 slhah.addLine(" 15 1.777000e+00 # ymtau \n");
0531 slhah.addLine("\n");
0532 slhah.addLine("###################################\n");
0533 slhah.addLine("## INFORMATION FOR DECAY\n");
0534 slhah.addLine("###################################\n");
0535 slhah.addLine("DECAY 6 1.491500e+00 # WT \n");
0536 slhah.addLine("DECAY 15 2.270000e-12 # WTau \n");
0537 slhah.addLine("DECAY 23 2.441404e+00 # WZ \n");
0538 slhah.addLine("DECAY 24 2.047600e+00 # WW \n");
0539 slhah.addLine("DECAY 25 6.382339e-03 # WH \n");
0540 slhah.addLine("## Dependent parameters, given by model restrictions.\n");
0541 slhah.addLine("## Those values should be edited following the \n");
0542 slhah.addLine("## analytical expression. MG5 ignores those values \n");
0543 slhah.addLine("## but they are important for interfacing the output of MG5\n");
0544 slhah.addLine("## to external program such as Pythia.\n");
0545 slhah.addLine("DECAY 1 0.000000 # d : 0.0 \n");
0546 slhah.addLine("DECAY 2 0.000000 # u : 0.0 \n");
0547 slhah.addLine("DECAY 3 0.000000 # s : 0.0 \n");
0548 slhah.addLine("DECAY 4 0.000000 # c : 0.0 \n");
0549 slhah.addLine("DECAY 5 0.000000 # b : 0.0 \n");
0550 slhah.addLine("DECAY 11 0.000000 # e- : 0.0 \n");
0551 slhah.addLine("DECAY 12 0.000000 # ve : 0.0 \n");
0552 slhah.addLine("DECAY 13 0.000000 # mu- : 0.0 \n");
0553 slhah.addLine("DECAY 14 0.000000 # vm : 0.0 \n");
0554 slhah.addLine("DECAY 16 0.000000 # vt : 0.0 \n");
0555 slhah.addLine("DECAY 21 0.000000 # g : 0.0 \n");
0556 slhah.addLine("DECAY 22 0.000000 # a : 0.0\n");
0557
0558 return slhah;
0559 }
0560
0561
0562 void EmbeddingLHEProducer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0563
0564
0565 edm::ParameterSetDescription desc;
0566 desc.setUnknown();
0567 descriptions.addDefault(desc);
0568 }
0569
0570
0571 DEFINE_FWK_MODULE(EmbeddingLHEProducer);