HepMC3 event record library
GenCrossSection.h
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#ifndef HEPMC3_CROSS_SECTION_H
7#define HEPMC3_CROSS_SECTION_H
8/**
9 * @file GenCrossSection.h
10 * @brief Definition of attribute \b class GenCrossSection
11 *
12 * @class HepMC3::GenCrossSection
13 * @brief Stores additional information about cross-section
14 *
15 * This is an example of event attribute used to store cross-section information
16 *
17 * This class is meant to be used to pass, on an event by event basis,
18 * the current best guess of the total cross section.
19 * It is expected that the final cross section will be stored elsewhere.
20 *
21 * - double cross_section; // cross section in pb.
22 * - double cross_section_error; // error associated with this cross section.
23 * - long accepted_events; ///< The number of events generated so far.
24 * - long attempted_events; ///< The number of events attempted so far.
25 *
26 * In addition, several cross sections and related info can be
27 * included in case of runs with mulltiple weights.
28 *
29 * The units of cross_section and cross_section_error are expected to be pb.
30 *
31 * @ingroup attributes
32 *
33 */
34#include <iostream>
35#include <algorithm>
36#include "HepMC3/Attribute.h"
37
38namespace HepMC3 {
39/** Deprecated */
40using namespace std;
41
42class GenCrossSection : public Attribute {
43
44//
45// Fields
46//
47private:
48
49 long accepted_events; ///< The number of events generated so far.
50 long attempted_events; ///< The number of events attempted so far.
51
52 std::vector<double> cross_sections; ///< Per-weight cross-section.
53 std::vector<double> cross_section_errors; ///< Per-weight errors.
54//
55// Functions
56//
57public:
58 /** @brief Implementation of Attribute::from_string */
59 bool from_string(const std::string &att) override;
60
61 /** @brief Implementation of Attribute::to_string */
62 bool to_string(std::string &att) const override;
63 /// @name Deprecated functionality
64 /// @{
65 /// @brief Set all fields
66 /// @deprecated Use set_cross_section(const std::vector<double>& xs, const std::vector<double>& xs_err instead
67 void set_cross_section(const double& xs, const double& xs_err,const long& n_acc = -1, const long& n_att = -1);
68
69 /// @}
70 /** @brief Set all fields */
71 void set_cross_section(const std::vector<double>& xs, const std::vector<double>& xs_err,const long& n_acc = -1, const long& n_att = -1);
72
73 /** @brief Get the cross-sections
74 */
75 const std::vector<double>& xsecs() const { return cross_sections; }
76
77 /** @brief Get the cross-section errors
78 */
79 const std::vector<double>& xsec_errs() const { return cross_section_errors; }
80
81
82 /** @brief Set the number of accepted events
83 */
84 void set_accepted_events(const long& n_acc ) {
85 accepted_events=n_acc;
86 }
87
88 /** @brief Set the number of attempted events
89 */
90 void set_attempted_events(const long& n_att ) {
91 attempted_events=n_att;
92 }
93
94 /** @brief Get the number of accepted events
95 */
96 long get_accepted_events() const {
97 return accepted_events;
98 }
99
100 /** @brief Get the number of attempted events
101 */
102 long get_attempted_events() const {
103 return attempted_events;
104 }
105
106 /** @brief Set the cross section corresponding to the weight
107 named \a wName.
108 */
109 void set_xsec(const std::string& wName,const double& xs) {
110 int pos = windx(wName);
111 if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec(const std::string&,const double&): no weight with given name in this run");
112 set_xsec(pos, xs);
113 }
114
115 /** @brief Set the cross section corresponding to the weight with
116 index \a indx.
117 */
118 void set_xsec(const unsigned long& index, const double& xs) {
119 if ( index >= cross_sections.size() ) {throw std::runtime_error("GenCrossSection::set_xsec(const unsigned long&): index outside of range");}
120 cross_sections[index] = xs;
121 }
122
123 /** @brief Set the cross section error corresponding to the weight
124 named \a wName.
125 */
126 void set_xsec_err(const std::string& wName, const double& xs_err) {
127 int pos = windx(wName);
128 if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec_err(const std::string&,const double&): no weight with given name in this run");
129 set_xsec_err(pos, xs_err);
130 }
131
132 /** @brief Set the cross section error corresponding to the weight
133 with index \a indx.
134 */
135 void set_xsec_err(const unsigned long& index, const double& xs_err) {
136 if ( index >= cross_section_errors.size() ) {throw std::runtime_error("GenCrossSection::set_xsec_err(const unsigned long&): index outside of range");}
137 cross_section_errors[index] = xs_err;
138 }
139
140 /** @brief Get the cross section corresponding to the weight named
141 \a wName.
142 */
143 double xsec(const std::string& wName) const {
144 int pos = windx(wName);
145 if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec(const std::string&): no weight with given name in this run");
146 return xsec(pos);
147 }
148
149 /** @brief Get the cross section corresponding to the weight with index
150 \a indx.
151 */
152 double xsec(const unsigned long& index = 0) const {
153 if ( index < cross_sections.size() ) { return cross_sections.at(index); }
154 else { throw std::runtime_error("GenCrossSection::xsec(const unsigned long&): index outside of range");}
155 return 0.0;
156 }
157
158 /** @brief Get the cross section error corresponding to the weight
159 named \a wName.
160 */
161 double xsec_err(const std::string& wName) const {
162 int pos = windx(wName);
163 if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec_err(const std::string&): no weight with given name in this run");
164 return xsec_err(pos);
165 }
166
167 /** @brief Get the cross section error corresponding to the weight
168 with index \a indx.
169 */
170 double xsec_err(const unsigned long& index = 0) const {
171 if ( index < cross_section_errors.size() ) {return cross_section_errors.at(index);}
172 else { throw std::runtime_error("GenCrossSection::xsec_err(const unsigned long&): index outside of range");}
173 return 0.0;
174 }
175
176 bool operator==( const GenCrossSection& ) const; ///< Operator ==
177 bool operator!=( const GenCrossSection& ) const; ///< Operator !=
178 bool is_valid() const; ///< Verify that the instance contains non-zero information
179
180private:
181
182 /** @brief get the weight index given a weight name. */
183 int windx(const std::string& wName) const;
184
185};
186
187} // namespace HepMC3
188
189#endif
Definition of class Attribute, class IntAttribute and class StringAttribute.
Base attribute class.
Definition: Attribute.h:44
Stores additional information about cross-section.
double xsec(const std::string &wName) const
Get the cross section corresponding to the weight named wName.
double xsec(const unsigned long &index=0) const
Get the cross section corresponding to the weight with index indx.
long get_accepted_events() const
Get the number of accepted events.
long attempted_events
The number of events attempted so far.
void set_xsec_err(const unsigned long &index, const double &xs_err)
Set the cross section error corresponding to the weight with index indx.
std::vector< double > cross_section_errors
Per-weight errors.
void set_attempted_events(const long &n_att)
Set the number of attempted events.
void set_accepted_events(const long &n_acc)
Set the number of accepted events.
const std::vector< double > & xsecs() const
Get the cross-sections.
std::vector< double > cross_sections
Per-weight cross-section.
bool is_valid() const
Verify that the instance contains non-zero information.
const std::vector< double > & xsec_errs() const
Get the cross-section errors.
bool from_string(const std::string &att) override
Implementation of Attribute::from_string.
double xsec_err(const std::string &wName) const
Get the cross section error corresponding to the weight named wName.
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.
double xsec_err(const unsigned long &index=0) const
Get the cross section error corresponding to the weight with index indx.
long get_attempted_events() const
Get the number of attempted events.
void set_xsec(const std::string &wName, const double &xs)
Set the cross section corresponding to the weight named wName.
void set_xsec_err(const std::string &wName, const double &xs_err)
Set the cross section error corresponding to the weight named wName.
void set_xsec(const unsigned long &index, const double &xs)
Set the cross section corresponding to the weight with index indx.
bool to_string(std::string &att) const override
Implementation of Attribute::to_string.
long accepted_events
The number of events generated so far.
HepMC3 main namespace.