HepMC3 event record library
GenCrossSection.cc
Go to the documentation of this file.
1// -*- C++ -*-
2//
3// This file is part of HepMC
4// Copyright (C) 2014-2023 The HepMC collaboration (see AUTHORS for details)
5//
6/**
7 * @file GenCrossSection.cc
8 * @brief Implementation of \b class GenCrossSection
9 *
10 */
11#include <cstdlib> // atoi
12#include <cstring> // memcmp
13#include <iomanip>
14#include <sstream>
15
17#include "HepMC3/GenEvent.h"
18
19
20namespace HepMC3 {
21
22
23int GenCrossSection::windx(const std::string& wName) const {
24 if ( !event() || !event()->run_info() ) {return 0;}
25 return event()->run_info()->weight_index(wName);
26}
27
28void GenCrossSection::set_cross_section(const double& xs, const double& xs_err, const long& n_acc, const long& n_att) {
29 double cross_section = xs;
30 double cross_section_error = xs_err;
31 accepted_events = n_acc;
32 attempted_events = n_att;
33 size_t N = std::max( event() ? event()->weights().size() : 0, size_t{1});
34 cross_sections = std::vector<double>(N, cross_section);
35 cross_section_errors = std::vector<double>(N, cross_section_error);
36}
37
38
39void GenCrossSection::set_cross_section(const std::vector<double>& xs, const std::vector<double>& xs_err, const long& n_acc, const long& n_att) {
40 cross_sections = xs;
41 cross_section_errors = xs_err;
42 accepted_events = n_acc;
43 attempted_events = n_att;
44}
45
46
47bool GenCrossSection::from_string(const std::string &att) {
48 const char *cursor = att.data();
49 cross_sections.clear();
51
52
53 double cross_section = atof(cursor);
54 cross_sections.emplace_back(cross_section);
55
56 if ( !(cursor = strchr(cursor+1, ' ')) ) {return false;}
57 double cross_section_error = atof(cursor);
58 cross_section_errors.emplace_back(cross_section_error);
59
60 if ( !(cursor = strchr(cursor+1, ' ')) ) {
61 accepted_events = -1;
63 } else {
64 accepted_events = atoi(cursor);
65 if ( !(cursor = strchr(cursor+1, ' ')) ) { attempted_events = -1; }
66 else { attempted_events = atoi(cursor); }
67 }
68 const size_t nweights = event() ? std::max(event()->weights().size(), size_t{1}) : size_t{1};
69 for (;;) {
70 if ( !(cursor = strchr(cursor+1, ' ')) ) break;
71 cross_sections.emplace_back(atof(cursor));
72 if ( !(cursor = strchr(cursor+1, ' ')) ) break;
73 cross_section_errors.emplace_back(atof(cursor));
74 }
75 if (cross_sections.size() != cross_section_errors.size()) {
76 HEPMC3_WARNING_LEVEL(800,"GenCrossSection::from_string: number of cross-sections and errors differ "
77 << cross_sections.size() << " vs " << cross_section_errors.size() << "). Ill-formed input:" << att)
78 }
79 // Use the default values to fill the vector to the size of N.
80 size_t oldxsecsize = cross_sections.size();
81 if (oldxsecsize > 1 && oldxsecsize != nweights) {
82 HEPMC3_WARNING_LEVEL(800,"GenCrossSection::from_string: the number of cross-sections (N = " << cross_sections.size() << ") does not match the number of weights (Nw = " << event()->weights().size() << ")")
83 }
84 for (size_t i = oldxsecsize; i < nweights; i++) {
85 cross_sections.emplace_back(cross_section);
86 cross_section_errors.emplace_back(cross_section_error);
87 }
88 return true;
89}
90
91
92bool GenCrossSection::to_string(std::string &att) const {
93 std::ostringstream os;
94
95 os << std::setprecision(8) << std::scientific
96 << (cross_sections.empty()?0.0:cross_sections.at(0)) << " "
97 << (cross_section_errors.empty()?0.0:cross_section_errors.at(0)) << " "
98 << accepted_events << " "
100 if (event() && event()->weights().size() > 0 &&
101 cross_sections.size() > 1 &&
102 event()->weights().size() != cross_sections.size() ) {
103 HEPMC3_WARNING_LEVEL(800,"GenCrossSection::to_string: the number of cross-sections (N = "<< cross_sections.size() << ") does not match the number of weights (Nw = "<< event()->weights().size() << ")")
104 }
105 for (size_t i = 1; i < cross_sections.size(); ++i ) {
106 os << " " << cross_sections.at(i) << " " << (cross_section_errors.size()>i?cross_section_errors.at(i):0.0);
107 }
108 att = os.str();
109
110 return true;
111}
112
114 return ( memcmp( static_cast<const void*>(this), static_cast<const void*>(&a), sizeof(class GenCrossSection) ) == 0 );
115}
116
118 return !( a == *this );
119}
120
122 if ( cross_sections.empty() ) { return false; }
123 if ( cross_section_errors.empty() ) { return false; }
124 if ( cross_section_errors.size() != cross_sections.size() ) { return false; }
125 if ( cross_sections.at(0) != 0 ) { return true; }
126 if ( cross_section_errors.at(0) != 0 ) { return true; }
127 return false;
128}
129
130} // namespace HepMC3
#define HEPMC3_WARNING_LEVEL(LEVEL, MESSAGE)
Macro for printing HEPMC3_HEPMC3_WARNING messages.
Definition Errors.h:34
Definition of attribute class GenCrossSection.
Definition of class GenEvent.
const GenEvent * event() const
Definition Attribute.h:108
Stores additional information about cross-section.
long attempted_events
The number of events attempted so far.
std::vector< double > cross_section_errors
Per-weight errors.
void set_cross_section(const double &xs, const double &xs_err, const long &n_acc=-1, const long &n_att=-1)
Set all fields.
std::vector< double > cross_sections
Per-weight cross-section.
bool is_valid() const
Verify that the instance contains non-zero information.
bool from_string(const std::string &att) override
Implementation of Attribute::from_string.
bool operator==(const GenCrossSection &) const
Operator ==.
bool operator!=(const GenCrossSection &) const
Operator !=.
int windx(const std::string &wName) const
get the weight index given a weight name.
bool to_string(std::string &att) const override
Implementation of Attribute::to_string.
long accepted_events
The number of events generated so far.
std::shared_ptr< GenRunInfo > run_info() const
Get a pointer to the the GenRunInfo object.
Definition GenEvent.h:144
HepMC3 main namespace.