2#include <pybind11/pybind11.h>
3#include <pybind11/stl.h>
4#include <pybind11/numpy.h>
9namespace py = pybind11;
19 PYBIND11_OVERRIDE_NAME(
28 PYBIND11_OVERRIDE_NAME(
40 auto m_parameters = m.def_submodule(
"parameters");
42 "This is a Python binding of MaCh3s C++ parameters library.";
46 py::enum_<SystType>(m_parameters,
"SystematicType")
55 py::init<
const std::vector<std::string>&,
const char *,
double,
int,
int>(),
56 "Construct a parameters object from a set of yaml files that define the systematic parameters \n\
57 :param yaml_files: The name of the yaml file to initialise from. \n\
58 :param name: the name of this ParameterHandler object. \n\
59 :param threshold: threshold PCA threshold from 0 to 1. Default is -1 and means no PCA. \n\
60 :param first_PCA_par: FirstPCAdpar First PCA parameter that will be decomposed. \n\
61 :param last_PCA_par: LastPCAdpar First PCA parameter that will be decomposed.",
62 py::arg(
"yaml_files"),
64 py::arg(
"threshold") = -1.0,
65 py::arg(
"firs_PCA_par") = -999,
66 py::arg(
"last_PCA_par") = -999
70 "calculate_likelihood",
72 "Calculate penalty term based on inverted covariance matrix."
78 "Throw the proposed parameter by magnitude *mag* X sigma. \n\
79 :param mag: This value multiplied by the prior value of each parameter will be the width of the distribution that the parameter values are drawn from. ",
84 "get_internal_par_name",
92 "Get the internally used name of this parameter. \n\
93 :param index: The global index of the parameter",
106 "Get the name of this parameter. \n\
107 :param index: The global index of the parameter",
114 "Get the number of parameters that this ParameterHandler object knows about."
120 "Propose a step based on the covariances. Also feel free to overwrite if you want something more funky."
124 "get_proposal_array",
127 return py::memoryview::from_buffer<double>(
129 {self.GetNParameters()},
133 "Bind a python array to the parameter proposal values for this ParameterHandler object. \n\
134 This allows you to set e.g. a numpy array to 'track' the parameter proposal values. You could either use this to directly set the proposals, or to just read the values proposed by e.g. throw_par_prop() \n\
135 :warning: This should be set *AFTER* all of the parameters have been read in from the config file as it resizes the array to fit the number of parameters. \n\
136 :param array: This is the array that will be set. Size and contents don't matter as it will be changed to fit the parameters. "
145 py::init<
const std::vector<std::string>&,
const char *,
double,
int,
int>(),
146 "Construct a systematic ParameterHandler object from a set of yaml files that define the systematic parameters \n\
147 :param yaml_files: The name of the yaml file to initialise from. \n\
148 :param name: the name of this ParameterHandler object. \n\
149 :param threshold: threshold PCA threshold from 0 to 1. Default is -1 and means no PCA. \n\
150 :param first_PCA_par: FirstPCAdpar First PCA parameter that will be decomposed. \n\
151 :param last_PCA_par: LastPCAdpar First PCA parameter that will be decomposed.",
152 py::arg(
"yaml_files"),
153 py::arg(
"name") =
"xsec_cov",
154 py::arg(
"threshold") = -1.0,
155 py::arg(
"firs_PCA_par") = -999,
156 py::arg(
"last_PCA_par") = -999
162 "Get what type of systematic this parameters is (see :py:class:`pyMaCh3.m_parameters.SystematicType` for possible types). \n\
163 :param index: The global index of the parameter",
168 "get_par_spline_type",
170 "Get what type of spline this parameter is set to use (assuming that it is a spline type parameter). \n\
171 :param index: The index of the spline parameter",
176 "get_par_spline_name",
178 "Get the name of the spline associated with a spline parameter. This is generally what it is called in input spline files and can in principle be different to the parameters name. \n\
179 :param index: The index of the spline parameter",
@ kNorm
For normalisation parameters.
@ kSpline
For splined parameters (1D)
@ kSystTypes
This only enumerates.
@ kFunc
For functional parameters.
Base class responsible for handling of systematic error parameters. Capable of using PCA or using ada...
virtual void ProposeStep()
Generate a new proposed state.
void ThrowParProp(const double mag=1.)
Throw the proposed parameter by mag sigma. Should really just have the user specify this throw by hav...
double CalcLikelihood() const _noexcept_
Calc penalty term based on inverted covariance matrix.
ParameterHandlerBase(const std::vector< std::string > &YAMLFile, std::string name, double threshold=-1, int FirstPCAdpar=-999, int LastPCAdpar=-999)
ETA - constructor for a YAML file.
const std::vector< double > & GetParPropVec()
Get a reference to the proposed parameter values Can be useful if you want to track these without hav...
Class responsible for handling of systematic error parameters with different types defined in the con...
SplineInterpolation GetParSplineInterpolation(const int i) const
Get interpolation type for a given parameter.
EW: As ParameterHandlerBase is an abstract base class we have to do some gymnastics to get it to get ...
double GetLikelihood() override
Return CalcLikelihood if some params were thrown out of boundary return LARGE_LOGL
void ProposeStep() override
Generate a new proposed state.
std::string GetParName(const int i) const
Get name of parameter.
SystType GetParamType(const int i) const
Returns enum describing our param type.
int GetNParameters() const
Get number of params which will be different depending if using Eigen decomposition or not.
std::string GetParFancyName(const int i) const
Get fancy name of the Parameter.
std::string GetParSplineName(const int i) const
Get the name of the spline associated with the spline at index i.
void initParameters(py::module &m)