File indexing completed on 2025-01-09 23:33:34
0001 #include <algorithm>
0002 #include <iomanip>
0003 #include <iostream>
0004 #include <memory>
0005
0006 #include <cstdio>
0007 #include <cstring>
0008 #include <fstream>
0009 #include <sstream>
0010 #include <string>
0011 #include <vector>
0012
0013 #include <xercesc/sax2/Attributes.hpp>
0014 #include <xercesc/dom/DOM.hpp>
0015
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/Utilities/interface/Exception.h"
0019 #include "FWCore/Utilities/interface/TimeOfDay.h"
0020
0021 #include "GeneratorInterface/LHEInterface/interface/LHEReader.h"
0022 #include "GeneratorInterface/LHEInterface/interface/LHERunInfo.h"
0023 #include "GeneratorInterface/LHEInterface/interface/LHEEvent.h"
0024
0025 #include "Utilities/StorageFactory/interface/IOTypes.h"
0026 #include "Utilities/StorageFactory/interface/Storage.h"
0027 #include "Utilities/StorageFactory/interface/StorageFactory.h"
0028
0029 #include "XMLUtils.h"
0030
0031 XERCES_CPP_NAMESPACE_USE
0032
0033 namespace lhef {
0034
0035 static void logFileAction(char const *msg, std::string const &fileName) {
0036 edm::LogAbsolute("fileAction") << std::setprecision(0) << edm::TimeOfDay() << msg << fileName;
0037 edm::FlushMessageLog();
0038 }
0039
0040 class LHEReader::Source {
0041 public:
0042 Source() {}
0043 virtual ~Source() {}
0044 virtual XMLDocument *createReader(XMLDocument::Handler &handler) = 0;
0045 };
0046
0047 class LHEReader::FileSource : public LHEReader::Source {
0048 public:
0049 FileSource(const std::string &fileURL) {
0050 using namespace edm::storage;
0051 auto storage = StorageFactory::get()->open(fileURL, IOFlags::OpenRead);
0052
0053 if (!storage)
0054 throw cms::Exception("FileOpenError")
0055 << "Could not open LHE file \"" << fileURL << "\" for reading" << std::endl;
0056
0057 fileStream = std::make_unique<StorageWrap>(std::move(storage));
0058 }
0059
0060 ~FileSource() override {}
0061
0062 XMLDocument *createReader(XMLDocument::Handler &handler) override { return new XMLDocument(fileStream, handler); }
0063
0064 private:
0065 std::unique_ptr<StorageWrap> fileStream;
0066 };
0067
0068 class LHEReader::StringSource : public LHEReader::Source {
0069 public:
0070 StringSource(const std::string &inputs) {
0071 if (inputs.empty())
0072 throw cms::Exception("StreamOpenError") << "Empty LHE file string name \"" << std::endl;
0073
0074 std::stringstream *tmpis = new std::stringstream(inputs);
0075 fileStream.reset(tmpis);
0076 }
0077
0078 ~StringSource() override {}
0079
0080 XMLDocument *createReader(XMLDocument::Handler &handler) override { return new XMLDocument(fileStream, handler); }
0081
0082 private:
0083 std::unique_ptr<std::istream> fileStream;
0084 };
0085
0086 class LHEReader::XMLHandler : public XMLDocument::Handler {
0087 public:
0088 typedef std::vector<std::pair<std::string, std::string> > wgt_info;
0089 XMLHandler()
0090 : impl(nullptr),
0091 gotObject(kNone),
0092 mode(kNone),
0093 xmlHeader(nullptr),
0094 xmlEvent(nullptr),
0095 headerOk(false),
0096 npLO(-99),
0097 npNLO(-99) {}
0098 ~XMLHandler() override {
0099 if (xmlHeader)
0100 xmlHeader->release();
0101 if (xmlEvent)
0102 xmlEvent->release();
0103 }
0104
0105 enum Object { kNone = 0, kHeader, kInit, kComment, kEvent };
0106
0107 void reset() {
0108 headerOk = false;
0109 weightsinevent.clear();
0110 gotObject = kNone;
0111 mode = kNone;
0112 }
0113
0114 const wgt_info &weightInfo() const { return weightsinevent; }
0115
0116 protected:
0117 void startElement(const XMLCh *const uri,
0118 const XMLCh *const localname,
0119 const XMLCh *const qname,
0120 const Attributes &attributes) override;
0121
0122 void endElement(const XMLCh *const uri, const XMLCh *const localname, const XMLCh *const qname) override;
0123
0124 void characters(const XMLCh *const chars, const XMLSize_t length) override;
0125 void comment(const XMLCh *const chars, const XMLSize_t length) override;
0126
0127 private:
0128 friend class LHEReader;
0129
0130 bool skipEvent = false;
0131 std::unique_ptr<DOMImplementation> impl;
0132 std::string buffer;
0133 Object gotObject;
0134 Object mode;
0135 DOMDocument *xmlHeader;
0136 DOMDocument *xmlEvent;
0137 std::vector<DOMElement *> xmlNodes, xmlEventNodes;
0138 bool headerOk;
0139 std::vector<LHERunInfo::Header> headers;
0140 wgt_info weightsinevent;
0141 int npLO;
0142 int npNLO;
0143 std::vector<float> scales;
0144 int evtnum = -1;
0145 };
0146
0147 static void attributesToDom(DOMElement *dom, const Attributes &attributes) {
0148 for (unsigned int i = 0; i < attributes.getLength(); i++) {
0149 const XMLCh *name = attributes.getQName(i);
0150 const XMLCh *value = attributes.getValue(i);
0151
0152 dom->setAttribute(name, value);
0153 }
0154 }
0155
0156 static void fillHeader(LHERunInfo::Header &header, const char *data, int len = -1) {
0157 const char *end = len >= 0 ? (data + len) : nullptr;
0158 while (*data && (!end || data < end)) {
0159 std::size_t len = std::strcspn(data, "\r\n");
0160 if (end && data + len > end)
0161 len = end - data;
0162 if (data[len] == '\r' && data[len + 1] == '\n')
0163 len += 2;
0164 else if (data[len])
0165 len++;
0166 header.addLine(std::string(data, len));
0167 data += len;
0168 }
0169 }
0170
0171 void LHEReader::XMLHandler::startElement(const XMLCh *const uri,
0172 const XMLCh *const localname,
0173 const XMLCh *const qname,
0174 const Attributes &attributes) {
0175 std::string name((const char *)XMLSimpleStr(qname));
0176
0177 if (!headerOk) {
0178 if (name != "LesHouchesEvents")
0179 throw cms::Exception("InvalidFormat") << "LHE file has invalid header" << std::endl;
0180 headerOk = true;
0181 return;
0182 }
0183
0184 if (mode == kHeader) {
0185 DOMElement *elem = xmlHeader->createElement(qname);
0186 attributesToDom(elem, attributes);
0187 xmlNodes.back()->appendChild(elem);
0188 xmlNodes.push_back(elem);
0189 return;
0190 } else if (mode == kEvent) {
0191 if (skipEvent) {
0192 return;
0193 }
0194
0195 DOMElement *elem = xmlEvent->createElement(qname);
0196 attributesToDom(elem, attributes);
0197
0198
0199 if (name == "rwgt") {
0200 xmlEventNodes[0]->appendChild(elem);
0201 } else if (name == "wgt") {
0202 xmlEventNodes[1]->appendChild(elem);
0203 } else if (name == "scales") {
0204 for (XMLSize_t iscale = 0; iscale < attributes.getLength(); ++iscale) {
0205 int ipart = 0;
0206 const char *scalename = XMLSimpleStr(attributes.getQName(iscale));
0207 int nmatch = sscanf(scalename, "pt_clust_%d", &ipart);
0208
0209 if (nmatch != 1) {
0210 edm::LogError("Generator|LHEInterface") << "invalid attribute in <scales> tag" << std::endl;
0211 }
0212
0213 float scaleval;
0214 const char *scalevalstr = XMLSimpleStr(attributes.getValue(iscale));
0215 sscanf(scalevalstr, "%e", &scaleval);
0216
0217 scales.push_back(scaleval);
0218 }
0219 } else if (name == "event_num") {
0220 const char *evtnumstr = XMLSimpleStr(attributes.getValue(XMLString::transcode("num")));
0221 sscanf(evtnumstr, "%d", &evtnum);
0222 }
0223 xmlEventNodes.push_back(elem);
0224 return;
0225 } else if (mode == kInit) {
0226
0227 return;
0228 } else if (mode != kNone) {
0229 throw cms::Exception("InvalidFormat") << "LHE file has invalid format" << std::endl;
0230 }
0231
0232 if (name == "header") {
0233 if (!impl)
0234 impl.reset(DOMImplementationRegistry::getDOMImplementation(XMLUniStr("Core")));
0235
0236 xmlHeader = impl->createDocument(nullptr, qname, nullptr);
0237 xmlNodes.resize(1);
0238 xmlNodes[0] = xmlHeader->getDocumentElement();
0239 mode = kHeader;
0240 }
0241 if (name == "init") {
0242 mode = kInit;
0243 } else if (name == "event") {
0244 if (!skipEvent) {
0245 if (!impl)
0246 impl.reset(DOMImplementationRegistry::getDOMImplementation(XMLUniStr("Core")));
0247
0248 if (xmlEvent)
0249 xmlEvent->release();
0250 xmlEvent = impl->createDocument(nullptr, qname, nullptr);
0251 weightsinevent.resize(0);
0252 scales.clear();
0253
0254 npLO = -99;
0255 npNLO = -99;
0256 const XMLCh *npLOval = attributes.getValue(XMLString::transcode("npLO"));
0257 if (npLOval) {
0258 const char *npLOs = XMLSimpleStr(npLOval);
0259 sscanf(npLOs, "%d", &npLO);
0260 }
0261 const XMLCh *npNLOval = attributes.getValue(XMLString::transcode("npNLO"));
0262 if (npNLOval) {
0263 const char *npNLOs = XMLSimpleStr(npNLOval);
0264 sscanf(npNLOs, "%d", &npNLO);
0265 }
0266
0267 xmlEventNodes.resize(1);
0268 xmlEventNodes[0] = xmlEvent->getDocumentElement();
0269 }
0270 mode = kEvent;
0271 }
0272
0273 if (mode == kNone)
0274 throw cms::Exception("InvalidFormat") << "LHE file has invalid format" << std::endl;
0275
0276 buffer.clear();
0277 }
0278
0279 void LHEReader::XMLHandler::endElement(const XMLCh *const uri,
0280 const XMLCh *const localname,
0281 const XMLCh *const qname) {
0282 std::string name((const char *)XMLSimpleStr(qname));
0283
0284 if (mode) {
0285 if (mode == kHeader && xmlNodes.size() > 1) {
0286 xmlNodes.resize(xmlNodes.size() - 1);
0287 return;
0288 } else if (mode == kHeader) {
0289 std::unique_ptr<DOMLSSerializer> writer(impl->createLSSerializer());
0290 std::unique_ptr<DOMLSOutput> outputDesc(impl->createLSOutput());
0291 assert(outputDesc.get());
0292 outputDesc->setEncoding(XMLUniStr("UTF-8"));
0293
0294 for (DOMNode *node = xmlNodes[0]->getFirstChild(); node; node = node->getNextSibling()) {
0295 XMLSimpleStr buffer(writer->writeToString(node));
0296
0297 std::string type;
0298 const char *p, *q;
0299 DOMElement *elem;
0300
0301 switch (node->getNodeType()) {
0302 case DOMNode::ELEMENT_NODE:
0303 elem = static_cast<DOMElement *>(node);
0304 type = (const char *)XMLSimpleStr(elem->getTagName());
0305 p = std::strchr((const char *)buffer, '>') + 1;
0306 q = std::strrchr(p, '<');
0307 break;
0308 case DOMNode::COMMENT_NODE:
0309 type = "";
0310 p = buffer + 4;
0311 q = buffer + strlen(buffer) - 3;
0312 break;
0313 default:
0314 type = "<>";
0315 p = buffer + std::strspn(buffer, " \t\r\n");
0316 if (!*p)
0317 continue;
0318 q = p + strlen(p);
0319 }
0320 LHERunInfo::Header header(type);
0321 fillHeader(header, p, q - p);
0322 headers.push_back(header);
0323 }
0324
0325 xmlHeader->release();
0326 xmlHeader = nullptr;
0327 } else if (name == "event" && mode == kEvent &&
0328 (skipEvent || (!xmlEventNodes.empty()))) {
0329
0330 if (skipEvent) {
0331 gotObject = mode;
0332 mode = kNone;
0333 return;
0334 }
0335
0336 for (DOMNode *node = xmlEventNodes[0]->getFirstChild(); node; node = node->getNextSibling()) {
0337 switch (node->getNodeType()) {
0338 case DOMNode::ELEMENT_NODE:
0339 for (DOMNode *rwgt = xmlEventNodes[1]->getFirstChild(); rwgt; rwgt = rwgt->getNextSibling()) {
0340 DOMNode *attr = rwgt->getAttributes()->item(0);
0341 XMLSimpleStr atname(attr->getNodeValue());
0342 XMLSimpleStr weight(rwgt->getFirstChild()->getNodeValue());
0343 switch (rwgt->getNodeType()) {
0344 case DOMNode::ELEMENT_NODE:
0345 weightsinevent.push_back(std::make_pair((const char *)atname, (const char *)weight));
0346 break;
0347 default:
0348 break;
0349 }
0350 }
0351 break;
0352 case DOMNode::TEXT_NODE:
0353 {
0354 XMLSimpleStr data(node->getNodeValue());
0355 buffer.append(data);
0356 } break;
0357 default:
0358 break;
0359 }
0360 }
0361 } else if (mode == kEvent) {
0362
0363 return;
0364 }
0365
0366 if (gotObject != kNone)
0367 throw cms::Exception("InvalidState") << "Unexpected pileup in"
0368 " LHEReader::XMLHandler::endElement"
0369 << std::endl;
0370
0371 gotObject = mode;
0372 mode = kNone;
0373 }
0374 }
0375
0376 void LHEReader::XMLHandler::characters(const XMLCh *const data_, const XMLSize_t length) {
0377 if (mode == kHeader) {
0378 DOMText *text = xmlHeader->createTextNode(data_);
0379 xmlNodes.back()->appendChild(text);
0380 return;
0381 }
0382
0383 if (XMLSimpleStr::isAllSpaces(data_, length))
0384 return;
0385
0386 unsigned int offset = 0;
0387 while (offset < length && XMLSimpleStr::isSpace(data_[offset]))
0388 offset++;
0389
0390 if (mode == kEvent) {
0391 if (!skipEvent) {
0392 DOMText *text = xmlEvent->createTextNode(data_ + offset);
0393 xmlEventNodes.back()->appendChild(text);
0394 }
0395 return;
0396 }
0397
0398 if (mode == kNone)
0399 throw cms::Exception("InvalidFormat") << "LHE file has invalid format" << std::endl;
0400
0401 XMLSimpleStr data(data_ + offset);
0402 buffer.append(data);
0403 }
0404
0405 void LHEReader::XMLHandler::comment(const XMLCh *const data_, const XMLSize_t length) {
0406 if (mode == kHeader) {
0407 DOMComment *comment = xmlHeader->createComment(data_);
0408 xmlNodes.back()->appendChild(comment);
0409 return;
0410 }
0411
0412 XMLSimpleStr data(data_);
0413
0414 LHERunInfo::Header header;
0415 fillHeader(header, data);
0416 headers.push_back(header);
0417 }
0418
0419 LHEReader::LHEReader(const edm::ParameterSet ¶ms)
0420 : fileURLs(params.getUntrackedParameter<std::vector<std::string> >("fileNames")),
0421 strName(""),
0422 firstEvent(params.getUntrackedParameter<unsigned int>("skipEvents", 0)),
0423 maxEvents(params.getUntrackedParameter<int>("limitEvents", -1)),
0424 curIndex(0),
0425 handler(new XMLHandler()) {}
0426
0427 LHEReader::LHEReader(const std::vector<std::string> &fileNames, unsigned int firstEvent)
0428 : fileURLs(fileNames),
0429 strName(""),
0430 firstEvent(firstEvent),
0431 maxEvents(-1),
0432 curIndex(0),
0433 handler(new XMLHandler()) {}
0434
0435 LHEReader::LHEReader(const std::string &inputs, unsigned int firstEvent)
0436 : strName(inputs), firstEvent(firstEvent), maxEvents(-1), curIndex(0), handler(new XMLHandler()) {}
0437
0438 LHEReader::~LHEReader() {
0439
0440
0441
0442
0443 handler.release();
0444 curDoc.release();
0445 curSource.release();
0446 }
0447
0448 std::shared_ptr<LHEEvent> LHEReader::next(bool *newFileOpened) {
0449 while (curDoc.get() || curIndex < fileURLs.size() || (fileURLs.empty() && !strName.empty())) {
0450 if (!curDoc.get()) {
0451 if (!platform) {
0452
0453
0454 platform = XMLDocument::platformHandle();
0455 }
0456 if (!fileURLs.empty()) {
0457 logFileAction(" Initiating request to open LHE file ", fileURLs[curIndex]);
0458 curSource = std::make_unique<FileSource>(fileURLs[curIndex]);
0459 logFileAction(" Successfully opened LHE file ", fileURLs[curIndex]);
0460 if (newFileOpened != nullptr)
0461 *newFileOpened = true;
0462 ++curIndex;
0463 } else if (!strName.empty()) {
0464 curSource = std::make_unique<StringSource>(strName);
0465 }
0466 handler->reset();
0467 curDoc.reset(curSource->createReader(*handler));
0468 curRunInfo.reset();
0469 }
0470 handler->skipEvent = firstEvent > 0;
0471
0472 XMLHandler::Object event = handler->gotObject;
0473 handler->gotObject = XMLHandler::kNone;
0474
0475 switch (event) {
0476 case XMLHandler::kNone:
0477 if (!curDoc->parse()) {
0478 curDoc.reset();
0479 logFileAction(" Closed LHE file ", fileURLs[curIndex - 1]);
0480 return std::shared_ptr<LHEEvent>();
0481 }
0482 break;
0483
0484 case XMLHandler::kHeader:
0485 break;
0486
0487 case XMLHandler::kInit: {
0488 std::istringstream data;
0489 data.str(handler->buffer);
0490 handler->buffer.clear();
0491
0492 curRunInfo = std::make_shared<LHERunInfo>(data);
0493
0494 std::for_each(handler->headers.begin(),
0495 handler->headers.end(),
0496 std::bind(&LHERunInfo::addHeader, curRunInfo.get(), std::placeholders::_1));
0497 handler->headers.clear();
0498 } break;
0499
0500 case XMLHandler::kComment:
0501 break;
0502
0503 case XMLHandler::kEvent: {
0504 if (!curRunInfo.get())
0505 throw cms::Exception("InvalidState") << "Got LHE event without"
0506 " initialization."
0507 << std::endl;
0508
0509 if (firstEvent > 0) {
0510 firstEvent--;
0511 continue;
0512 }
0513
0514 if (maxEvents == 0)
0515 return std::shared_ptr<LHEEvent>();
0516 else if (maxEvents > 0)
0517 maxEvents--;
0518
0519 std::istringstream data;
0520 data.str(handler->buffer);
0521 handler->buffer.clear();
0522
0523 std::shared_ptr<LHEEvent> lheevent;
0524 lheevent = std::make_shared<LHEEvent>(curRunInfo, data);
0525 const XMLHandler::wgt_info &info = handler->weightsinevent;
0526 for (size_t i = 0; i < info.size(); ++i) {
0527 double num = -1.0;
0528 sscanf(info[i].second.c_str(), "%le", &num);
0529 lheevent->addWeight(gen::WeightsInfo(info[i].first, num));
0530 }
0531 lheevent->setNpLO(handler->npLO);
0532 lheevent->setNpNLO(handler->npNLO);
0533 lheevent->setEvtNum(handler->evtnum);
0534 handler->evtnum = -1;
0535
0536 if (!handler->scales.empty()) {
0537 lheevent->setScales(handler->scales);
0538 }
0539 return lheevent;
0540 }
0541 }
0542 }
0543
0544 return std::shared_ptr<LHEEvent>();
0545 }
0546
0547 }