# Introduction The Measurement Equation Analysis Library (MEAL) is a set of classes, implemented in C++, that may be used to perform non-linear least-squares and error propagation through matrix equations. Arbitrary functions are built by modular construction of elementary functions, which may be parameterized by any number or type of arguments. It was written by Willem van Straten ([\@straten](https://github.com/straten){.user-mention}) to support the research and development of 1. Measurement Equation Modelling (MEM) [Radio Astronomical Polarimetry and Point-Source Calibration](http://dx.doi.org/10.1086/383187), *The Astrophysical Journal Supplement Series*, 152:129 (2004) 2. Matrix Template Matching (MTM) [Radio Astronomical Polarimetry and High-Precision Pulsar Timing](http://dx.doi.org/10.1086/501001), *The Astrophysical Journal*, 642:1004 (2006) 3. Measurement Equation Template Matching (METM) [High-fidelity Radio Astronomical Polarimetry Using a Millisecond Pulsar as a Polarized Reference Source](http://dx.doi.org/10.1088/0067-0049/204/1/13), *The Astrophysical Journal Supplement Series*, 204:13 (2013) ## Function Characteristics A mathematical function or expression is characterized by its return type and input variables; a distinction is made between function parameters, $\boldsymbol{a}=(a_0, a_1, ... , a_N)$ and independent variables, or arguments, $\boldsymbol{x}=(x_0, x_1, ... , x_M)$. Parameters are double-precision floating point values that can vary freely during modeling. Arguments have arbitrary type and cannot be treated as free variables during modeling. The parameter and argument access interface is defined by the {cpp:class}`MEAL::Function` base class, and the return type is defined by *evaluation base classes* that inherit Function. The evaluation base classes define an evaluate method that returns a result, $f(\boldsymbol{x},\boldsymbol{a})$, and the partial derivatives of this result with respect to each of the parameters, $\partial f/\partial a_i$. The behaviour of function expression classes is organized into three categories: {cpp:class}`MEAL::Parameter`, {cpp:class}`MEAL::Argument`, and {cpp:class}`MEAL::Evaluation`. These behaviours are implemented by the children of three base classes, each of which inherits the {cpp:class}`MEAL::FunctionPolicy` base class. The behaviour policies are mutually exclusive, and derived classes can incorporate pre-defined behaviours by setting the appropriate policy attribute. ### Parameter Policy A function expression may have an arbitary number of scalar parameters, $\boldsymbol{a}=(a_0, a_1, ... , a_N)$. Each parameter has an associated name, estimated variance, and a flag that indicates if the parameter is free or fixed. Parameter management and access is implemented by children of the {cpp:class}`MEAL::ParameterPolicy` base class. ### Argument Policy A function expression may be further parameterized by an arbitrary number of independent variables, or arguments, $\boldsymbol{x}=(x_0, x_1, ... , x_M)$. Unlike parameters, arguments have no specified type, no estimated variance, and can never be free. Because they have no specified type, the interface between a Function and its arguments is mediated through the {cpp:class}`MEAL::Argument` and {cpp:class}`MEAL::Argument::Value` abstract base classes. Argument management and behaviour is implemented by children of the {cpp:class}`MEAL::ArgumentPolicy` base class. ### Evaluation Policy The return type of a function expression is unspecified in the {cpp:class}`MEAL::Function` base class definition. Therefore, derived classes inherit {cpp:class}`MEAL::Function` via the {cpp:class}`MEAL::Evaluable` template, which has a template argument that determines the return type, `Return`, and defines the following method. ```cpp virtual Result evaluate (std::vector* gradient = 0) const; ``` This method returns a result $M$ and, if `gradient` is not null, the partial derivative of the result with respect to each of the model parameters, $\partial M/\partial a_i$. For convenience of notation, there are intermediate classes that inherit the {cpp:class}`MEAL::Evaluable` template. For example, - {cpp:class}`MEAL::Scalar` - inherits `MEAL::Evaluable` and represents functions that return a scalar value, $f({\bf a}; {\bf x})$ - {cpp:class}`MEAL::Complex2` - inherits `MEAL::Evaluable>` and represents functions that return a $2\times2$ complex-valued matrix, $J({\bf a}; {\bf x})$ ## Model Composition A number of template classes may be used to simplify the modular construction of more complicated functions. These templates implement the following basic rules of differentiation. - {cpp:class}`MEAL::ChainRule` - an arbitrary function in which one or more parameters is set equal to the ordinate of a {cpp:class}`MEAL::Scalar` function - {cpp:class}`MEAL::BinaryRule` - an associative binary operation, such as the {cpp:class}`MEAL::SumRule` or {cpp:class}`MEAL::ProductRule`. For example, use {cpp:class}`MEAL::ChainRule` to create a feed in which the receptors have equal ellipticities. ```cpp auto feed = new MEAL::Feed; auto feed_chain = new MEAL::ChainRule; feed_chain->set_model (feed); auto ellipticity = new MEAL::ScalarParameter; feed_chain->set_constraint (0, ellipticity); feed_chain->set_constraint (2, ellipticity); ``` Use {cpp:class}`MEAL::ProductRule` create multiply a boost along the Stokes V axis and a rotation about the Stokes U axis (in a basis defined by linearly polarized receptors). ```cpp auto product = new MEAL::ProductRule; auto boost = new MEAL::Boost1 (Vector<3, double>::basis(2)); product->add_model (boost); auto rotation = new MEAL::Rotation1 (Vector<3, double>::basis(1)); product->add_model (boost); ``` ## Modular Construction TO DO: Show how new functions can be built up from more basic elements. ## Example Usage ### Non-linear Least-Squares Estimation TO DO: Document lmfit ### Error Propagation The Estimate template class is very useful for storing a value and its estimated variance. There are also operators and functions which enable the propagation of error estimates to derived quantities; for example: ```cpp Estimate y (0.596,0.0034); Estimate x (-0.83,0.0072); Estimate z = pow (y,x); ``` automatically computes the variance of the new variable, z. However, the Estimate template class fails when a variable appears more than once in an expression; e.g. Estimate ratio = x/x; should yield $1\pm0$; however, the Estimate template class does not recognize that the numerator and denominator are the same variable, and incorrectly sums the weighted variances. The problem of correctly computing the partial derivatives of an expression with respect to its variables makes use of the exact same functionality used to generate the gradient and Hessian matrix in non-linear least squares fitting. A simplified interface to this functionality is implemented by the {cpp:class}`MEAL::ScalarMath` class. {cpp:class}`MEAL::ScalarMath` objects may be conveniently initialized as a single parameter and its estimated variance using the Estimate template class. As with float and double types, ScalarMath objects may be combined using normal arithmetic operations and basic mathematical functions, creating Scalar functions of any number of parameters. For example: MEAL::ScalarMath x (Estimate (0.9587, 0.00058)); MEAL::ScalarMath y (Estimate (-0.283, 0.00034)); cerr << "Polar angle = " << atan2 (y, x) << endl; yields the output Polar angle = (-0.287039 +/- 0.0189612) As with any native type, the ScalarMath class can be used as a template argument, e.g. complex z (Estimate (0.87, 0.0041), Estimate (2.38, 0.0095)); complex w (Estimate (1.74, 0.0081), Estimate (-.63, 0.0043)); Jones jones (z, conj(z), conj(w), w); enabling error propagation through increasingly complex expressions.