Skip to content
Snippets Groups Projects

Extend API of GenCrossSection

Merged Andrii Verbytskyi requested to merge prevent_out_of_the_bounds_xsection into master
Files
5
@@ -64,6 +64,18 @@ public:
/** @brief Set all fields */
void set_cross_section(const double& xs, const double& xs_err,const long& n_acc = -1, const long& n_att = -1);
/** @brief Set all fields */
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);
/** @brief Get the cross-sections
*/
const std::vector<double>& xsecs() const { return cross_sections; }
/** @brief Get the cross-section errors
*/
const std::vector<double>& xsec_errs() const { return cross_section_errors; }
/** @brief Set the number of accepted events
*/
void set_accepted_events(const long& n_acc ) {
@@ -92,56 +104,70 @@ public:
named \a wName.
*/
void set_xsec(const std::string& wName,const double& xs) {
set_xsec(windx(wName), xs);
int pos = windx(wName);
if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec(const std::string&,const double&): no weight with given name in this run");
set_xsec(pos, xs);
}
/** @brief Set the cross section corresponding to the weight with
index \a indx.
*/
void set_xsec(const int& indx, const double& xs) {
cross_sections[indx] = xs;
void set_xsec(const unsigned long& index, const double& xs) {
if ( index >= cross_sections.size() ) {throw std::runtime_error("GenCrossSection::set_xsec(const unsigned long&): index outside of range");}
cross_sections[index] = xs;
}
/** @brief Set the cross section error corresponding to the weight
named \a wName.
*/
void set_xsec_err(const std::string& wName, const double& xs_err) {
set_xsec_err(windx(wName), xs_err);
int pos = windx(wName);
if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec_err(const std::string&,const double&): no weight with given name in this run");
set_xsec_err(pos, xs_err);
}
/** @brief Set the cross section error corresponding to the weight
with index \a indx.
*/
void set_xsec_err(const int& indx, const double& xs_err) {
cross_section_errors[indx] = xs_err;
void set_xsec_err(const unsigned long& index, const double& xs_err) {
if ( index >= cross_section_errors.size() ) {throw std::runtime_error("GenCrossSection::set_xsec_err(const unsigned long&): index outside of range");}
cross_section_errors[index] = xs_err;
}
/** @brief Get the cross section corresponding to the weight named
\a wName.
*/
double xsec(const std::string& wName) const {
return xsec(windx(wName));
int pos = windx(wName);
if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec(const std::string&): no weight with given name in this run");
return xsec(pos);
}
/** @brief Get the cross section corresponding to the weight with index
\a indx.
*/
double xsec(const int& indx = 0) const {
return cross_sections[indx];
double xsec(const unsigned long& index = 0) const {
if ( index < cross_sections.size() ) { return cross_sections.at(index); }
else { throw std::runtime_error("GenCrossSection::xsec(const unsigned long&): index outside of range");}
return 0.0;
}
/** @brief Get the cross section error corresponding to the weight
named \a wName.
*/
double xsec_err(const std::string& wName) const {
return xsec_err(windx(wName));
int pos = windx(wName);
if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec_err(const std::string&): no weight with given name in this run");
return xsec_err(pos);
}
/** @brief Get the cross section error corresponding to the weight
with index \a indx.
*/
double xsec_err(const int& indx = 0) const {
return cross_section_errors[indx];
double xsec_err(const unsigned long& index = 0) const {
if ( index < cross_section_errors.size() ) {return cross_section_errors.at(index);}
else { throw std::runtime_error("GenCrossSection::xsec_err(const unsigned long&): index outside of range");}
return 0.0;
}
bool operator==( const GenCrossSection& ) const; ///< Operator ==
Loading