From 8e2949e7a990f3240ce9d12cd1b841d1b0a9f8b2 Mon Sep 17 00:00:00 2001 From: Andreas Suter Date: Sun, 23 Nov 2025 17:47:46 +0100 Subject: [PATCH] improve the doxygen docu of PTheory.* --- src/classes/PTheory.cpp | 353 ++++++++++++++++++++++++++++++---------- src/include/PTheory.h | 286 ++++++++++++++++++++++++++------ 2 files changed, 500 insertions(+), 139 deletions(-) diff --git a/src/classes/PTheory.cpp b/src/classes/PTheory.cpp index 3fd45c08b..b39a7c535 100644 --- a/src/classes/PTheory.cpp +++ b/src/classes/PTheory.cpp @@ -53,16 +53,18 @@ extern std::vector gGlobalUserFcn; // Constructor //-------------------------------------------------------------------------- /** - *

Constructor. + * \brief Parses the THEORY block and builds a binary expression tree. * - *

The theory block his parsed and mapped to a tree. Every line (except the '+') - * is a theory object by itself, i.e. the whole theory is build up recursivly. - * Example: + * The constructor recursively parses the theory block from the MSR file and + * builds a binary tree structure. Each theory function becomes a node, with + * multiplication creating right children and addition creating left children. + * + * Example Theory Block: * \verbatim - Theory block: - a 1 - tf 2 3 - se 4 + THEORY + a 1 # asymmetry + tf 2 3 # TFieldCos (phase, frequency) + se 4 # simpleExp (rate) + a 5 tf 6 7 @@ -72,7 +74,7 @@ extern std::vector gGlobalUserFcn; tf 10 11 \endverbatim * - *

This is mapped into the following binary tree: + * Resulting Binary Tree: * \verbatim a 1 +/ \* @@ -85,11 +87,40 @@ extern std::vector gGlobalUserFcn; se 8 \endverbatim * - * \param msrInfo msr-file handler - * \param runNo msr-file run number - * \param hasParent flag needed in order to know if PTheory is the root object or a child - * false (default) means this is the root object - * true means this is part of an already existing object tree + * The tree is evaluated as: (a1 * tf23 * se4) + (a5 * tf67 * se8) + (a9 * tf1011) + * + * Parsing Algorithm: + * -# Initialize member variables and static counters (lineNo, depth) + * -# Get current line from theory block (skip if beyond end) + * -# Remove comments: text after '(' or '#' + * -# Tokenize line into function name and parameters + * -# Search fgTheoDataBase for function (by name or abbreviation) + * -# Validate parameter count (except for polynom and userFcn) + * -# For each parameter token: + * - If "mapN": resolve via RUN block map vector + * - If "funN": get function index + MSR_PARAM_FUN_OFFSET + * - Otherwise: direct parameter number (1-based → 0-based) + * -# If next line is not '+': increment depth, create fMul child recursively + * -# If depth=0 and next line is '+': skip '+', create fAdd child recursively + * -# For root node: call MakeCleanAndTidyTheoryBlock() to format output + * -# For userFcn: load shared library and instantiate class object + * + * User Function Handling: + * User functions require special parsing: + * - Token 1: shared library name (e.g., "libMyFunc.so") + * - Token 2: class name (e.g., "TMyFunction") + * - Token 3+: parameters + * The library is loaded via gSystem->Load() and the class instantiated via + * TClass::GetClass()->New(). + * + * \param msrInfo Pointer to MSR file handler (NOT owned, must outlive PTheory) + * \param runNo Zero-based run number for map parameter resolution + * \param hasParent false (default) for root theory node, true for child nodes. + * Controls reset of static line/depth counters. + * + * \warning Sets fValid=false on errors. Always check IsValid() after construction. + * + * \see fgTheoDataBase, SearchDataBase(), IsValid() */ PTheory::PTheory(PMsrHandler *msrInfo, UInt_t runNo, const Bool_t hasParent) : fMsrInfo(msrInfo) { @@ -328,7 +359,18 @@ PTheory::PTheory(PMsrHandler *msrInfo, UInt_t runNo, const Bool_t hasParent) : f // Destructor //-------------------------------------------------------------------------- /** - *

Destructor + * \brief Destructor that recursively cleans up the expression tree. + * + * Releases all resources allocated by this theory node and its children: + * -# Clears parameter number vector (fParamNo) + * -# Clears user function parameter cache (fUserParam) + * -# Clears LF integral caches (fLFIntegral, fDynLFFuncValue) + * -# Recursively deletes child nodes via CleanUp() + * -# Deletes user function object if present + * -# Clears global user function vector + * + * \note Child nodes (fAdd, fMul) are deleted via CleanUp() which traverses + * the entire tree recursively. */ PTheory::~PTheory() { @@ -353,11 +395,22 @@ PTheory::~PTheory() // IsValid //-------------------------------------------------------------------------- /** - *

Checks if the theory tree is valid. + * \brief Recursively checks if the entire theory expression tree is valid. * - * return: - * - true if valid - * - false otherwise + * Traverses the binary tree and checks fValid for this node and all + * descendant nodes. A theory tree is valid only if ALL nodes are valid. + * + * Validation logic: + * - If fMul and fAdd both exist: return fValid && fMul->IsValid() && fAdd->IsValid() + * - If only fMul exists: return fValid && fMul->IsValid() + * - If only fAdd exists: return fValid && fAdd->IsValid() + * - If leaf node: return fValid + * + * \return true if this node and all descendants are valid, false otherwise + * + * \note Should always be called after construction to verify successful parsing. + * A single invalid node anywhere in the tree causes the entire theory + * to be marked invalid. */ Bool_t PTheory::IsValid() { @@ -379,14 +432,38 @@ Bool_t PTheory::IsValid() //-------------------------------------------------------------------------- /** - *

Evaluates the theory tree. + * \brief Evaluates the theory expression tree at a given time point. * - * return: - * - theory value for the given sets of parameters + * Recursively evaluates the binary tree, combining function values according + * to the tree structure (multiplication for consecutive lines, addition for + * '+' separated sections). * - * \param t time for single histogram, asymmetry, and mu-minus fits, x-axis value non-muSR fits - * \param paramValues vector with the parameters - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Evaluation Algorithm: + * Based on tree structure, one of four cases applies: + * -# fMul && fAdd: this_func(t) * fMul->Func(t) + fAdd->Func(t) + * -# fMul only: this_func(t) * fMul->Func(t) + * -# fAdd only: this_func(t) + fAdd->Func(t) + * -# Leaf node: this_func(t) + * + * Parameter Resolution: + * Within each function, parameters are accessed via fParamNo indices: + * - If fParamNo[i] < MSR_PARAM_FUN_OFFSET: use paramValues[fParamNo[i]] + * - If fParamNo[i] >= MSR_PARAM_FUN_OFFSET: use funcValues[fParamNo[i] - MSR_PARAM_FUN_OFFSET] + * + * LF Caching: + * Longitudinal field functions (StaticGaussKTLF, DynamicLorentzKTLF, etc.) + * use cached integral values that are recalculated only when parameters change. + * This is handled transparently within each LF function. + * + * \param t Time in microseconds for μSR fits, or x-axis value for non-μSR data + * \param paramValues Vector of current fit parameter values from FITPARAMETER block + * \param funcValues Vector of evaluated function values from FUNCTIONS block + * + * \return Evaluated polarization/asymmetry value at time t + * + * \note This method is const but uses mutable members for LF caching. + * + * \see Constant(), TFCos(), StaticGaussKT(), etc. for individual function implementations */ Double_t PTheory::Func(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -765,13 +842,22 @@ void PTheory::CleanUp(PTheory *theo) //-------------------------------------------------------------------------- /** - *

Searches the internal look up table for the theory name, and if found - * set some internal variables to proper values. + * \brief Searches fgTheoDataBase for a theory function by name or abbreviation. * - * return: - * - index of the theory function + * Performs case-insensitive search through the theory database to find a + * matching function. If found, sets fType and fNoOfParam member variables. * - * \param name theory name + * Search matches either: + * - Full function name (e.g., "simplExpo", "TFieldCos") + * - Abbreviation (e.g., "se", "tf") + * + * \param name Function name or abbreviation to search for + * + * \return Index in fgTheoDataBase (0 to THEORY_MAX-1), or THEORY_UNDEFINED (-1) if not found + * + * \note Side effects: Sets fType and fNoOfParam when function is found. + * + * \see fgTheoDataBase */ Int_t PTheory::SearchDataBase(TString name) { @@ -793,11 +879,19 @@ Int_t PTheory::SearchDataBase(TString name) // GetUserFcnIdx (private) //-------------------------------------------------------------------------- /** - *

Counts the number of user functions in the theory block up to lineNo. + * \brief Counts user function occurrences in theory block up to a given line. * - * return: to number of user functions found up to lineNo + * Scans the theory block from line 1 to lineNo (inclusive) and counts + * how many lines contain "userFcn". This index is used to: + * - Identify which global user function object to use + * - Track multiple user functions in a single theory * - * \param lineNo current line number in the theory block + * \param lineNo Current line number being parsed (1-based in theory block) + * + * \return Index of current user function (0-based count), or -1 if no userFcn + * entries found up to lineNo + * + * \note The search is case-insensitive ("userFcn", "USERFCN", etc. all match). */ Int_t PTheory::GetUserFcnIdx(UInt_t lineNo) const { @@ -823,10 +917,28 @@ Int_t PTheory::GetUserFcnIdx(UInt_t lineNo) const // MakeCleanAndTidyTheoryBlock private //-------------------------------------------------------------------------- /** - *

Takes the theory block form the msr-file, and makes a nicely formated - * theory block. + * \brief Reformats the theory block with consistent spacing and comments. * - * \param fullTheoryBlock msr-file text lines of the theory block + * Processes each line in the theory block to produce clean, well-formatted + * output suitable for writing back to MSR files. Formatting includes: + * - Function name left-aligned in 10-character field + * - Parameters right-aligned in 6-character fields + * - Comment string appended showing parameter names + * + * Example transformation: + * \verbatim + * Input: "tf 2 3" + * Output: "TFieldCos 2 3 (phase frequency)" + * \endverbatim + * + * Special handling for: + * - '+' lines: left unchanged + * - polynom: delegated to MakeCleanAndTidyPolynom() + * - userFcn: delegated to MakeCleanAndTidyUserFcn() + * + * \param fullTheoryBlock Pointer to theory block lines (modified in place) + * + * \note Called only by root theory node (hasParent=false) after successful parsing. */ void PTheory::MakeCleanAndTidyTheoryBlock(PMsrLines *fullTheoryBlock) { @@ -914,10 +1026,16 @@ void PTheory::MakeCleanAndTidyTheoryBlock(PMsrLines *fullTheoryBlock) // MakeCleanAndTidyPolynom private //-------------------------------------------------------------------------- /** - *

Prettifies a polynom theory line. + * \brief Formats a polynomial theory line with proper spacing. * - * \param i line index of the theory polynom line to be prettified - * \param fullTheoryBlock msr-file text lines of the theory block + * Polynomial functions have variable parameter count, so they require + * special formatting. Creates output in format: + * \verbatim + * polynom tshift p0 p1 p2 ... pn (tshift p0 p1 ... pn) + * \endverbatim + * + * \param i Line index (1-based) in the theory block + * \param fullTheoryBlock Pointer to theory block lines (modified in place) */ void PTheory::MakeCleanAndTidyPolynom(UInt_t i, PMsrLines *fullTheoryBlock) { @@ -971,10 +1089,16 @@ void PTheory::MakeCleanAndTidyPolynom(UInt_t i, PMsrLines *fullTheoryBlock) // MakeCleanAndTidyUserFcn private //-------------------------------------------------------------------------- /** - *

Prettifies a user function theory line. + * \brief Formats a user function theory line with proper spacing. * - * \param i line index of the user function line to be prettified - * \param fullTheoryBlock msr-file text lines of the theory block + * User functions have variable structure (library, class, parameters), + * so they require special formatting. Creates output in format: + * \verbatim + * userFcn libName.so ClassName param1 param2 ... paramN + * \endverbatim + * + * \param i Line index (1-based) in the theory block + * \param fullTheoryBlock Pointer to theory block lines (modified in place) */ void PTheory::MakeCleanAndTidyUserFcn(UInt_t i, PMsrLines *fullTheoryBlock) { @@ -1010,15 +1134,20 @@ void PTheory::MakeCleanAndTidyUserFcn(UInt_t i, PMsrLines *fullTheoryBlock) //-------------------------------------------------------------------------- /** - *

theory function: Const - * \f[ = const \f] + * \brief Returns a constant value (baseline, background). * - * meaning of paramValues: const + * Mathematical form: + * \f[ P(t) = c \f] * - * return: function value + * Used for constant offsets, fixed backgrounds, or baseline corrections. * - * \param paramValues vector with the parameters - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Parameters: + * - fParamNo[0]: Constant value c + * + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Constant value c */ Double_t PTheory::Constant(const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -1038,15 +1167,22 @@ Double_t PTheory::Constant(const PDoubleVector& paramValues, const PDoubleVector //-------------------------------------------------------------------------- /** - *

theory function: Asymmetry - * \f[ = A \f] + * \brief Returns the initial asymmetry value. * - * meaning of paramValues: \f$A\f$ + * Mathematical form: + * \f[ P(t) = A \f] * - * return: function value + * The asymmetry represents the initial polarization of the muon beam and is + * typically used as a multiplicative prefactor for oscillating/relaxing functions. + * Typical values range from 0.15 to 0.30 depending on detector geometry. * - * \param paramValues vector with the parameters - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Parameters: + * - fParamNo[0]: Asymmetry value A + * + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Asymmetry value A */ Double_t PTheory::Asymmetry(const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -1066,17 +1202,25 @@ Double_t PTheory::Asymmetry(const PDoubleVector& paramValues, const PDoubleVecto //-------------------------------------------------------------------------- /** - *

theory function: simple exponential - * \f[ = \exp\left(-\lambda t\right) \f] or - * \f[ = \exp\left(-\lambda [t-t_{\rm shift}] \right) \f] + * \brief Simple exponential relaxation function. * - * meaning of paramValues: \f$\lambda\f$ [, \f$t_{\rm shift}\f$] + * Mathematical form: + * \f[ P(t) = \exp\left(-\lambda t\right) \f] + * or with time shift: + * \f[ P(t) = \exp\left(-\lambda [t-t_{\rm shift}] \right) \f] * - * return: function value + * Represents exponential muon spin relaxation due to dynamic fluctuations + * in the fast-fluctuation limit (BPP-type relaxation). * - * \param t time in \f$(\mu\mathrm{s})\f$, or x-axis value for non-muSR fit - * \param paramValues parameter values - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Parameters: + * - fParamNo[0]: Relaxation rate λ (μs⁻¹) + * - fParamNo[1]: (optional) Time shift t_shift (μs) + * + * \param t Time in μs + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Exponential relaxation value at time t */ Double_t PTheory::SimpleExp(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -1106,17 +1250,29 @@ Double_t PTheory::SimpleExp(Double_t t, const PDoubleVector& paramValues, const //-------------------------------------------------------------------------- /** - *

theory function: generalized exponential - * \f[ = \exp\left(-[\lambda t]^\beta\right) \f] or - * \f[ = \exp\left(-[\lambda (t-t_{\rm shift})]^\beta\right) \f] + * \brief Generalized (stretched) exponential relaxation function. * - * meaning of paramValues: \f$\lambda\f$, \f$\beta\f$ [, \f$t_{\rm shift}\f$] + * Mathematical form: + * \f[ P(t) = \exp\left(-[\lambda t]^\beta\right) \f] + * or with time shift: + * \f[ P(t) = \exp\left(-[\lambda (t-t_{\rm shift})]^\beta\right) \f] * - * return: function value + * The stretched exponential (Kohlrausch function) represents relaxation + * in systems with a distribution of correlation times: + * - β = 1: Simple exponential + * - β < 1: Stretched exponential (distributed relaxation) + * - β = 2: Gaussian relaxation * - * \param t time in \f$(\mu\mathrm{s})\f$, or x-axis value for non-muSR fit - * \param paramValues parameter values - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Parameters: + * - fParamNo[0]: Relaxation rate λ (μs⁻¹) + * - fParamNo[1]: Stretching exponent β + * - fParamNo[2]: (optional) Time shift t_shift (μs) + * + * \param t Time in μs + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Stretched exponential value at time t */ Double_t PTheory::GeneralExp(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -1154,17 +1310,26 @@ Double_t PTheory::GeneralExp(Double_t t, const PDoubleVector& paramValues, const //-------------------------------------------------------------------------- /** - *

theory function: simple Gaussian - * \f[ = \exp\left(-1/2 [\sigma t]^2\right) \f] or - * \f[ = \exp\left(-1/2 [\sigma (t-t_{\rm shift})]^2\right) \f] + * \brief Simple Gaussian relaxation function. * - * meaning of paramValues: \f$\sigma\f$ [, \f$t_{\rm shift}\f$] + * Mathematical form: + * \f[ P(t) = \exp\left(-\frac{1}{2} [\sigma t]^2\right) \f] + * or with time shift: + * \f[ P(t) = \exp\left(-\frac{1}{2} [\sigma (t-t_{\rm shift})]^2\right) \f] * - * return: function value + * Represents Gaussian muon spin relaxation due to static random fields + * with a Gaussian distribution. The parameter σ = γ_μ ΔB where ΔB is the + * RMS field width. * - * \param t time in \f$(\mu\mathrm{s})\f$, or x-axis value for non-muSR fit - * \param paramValues parameter values - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Parameters: + * - fParamNo[0]: Relaxation rate σ (μs⁻¹) + * - fParamNo[1]: (optional) Time shift t_shift (μs) + * + * \param t Time in μs + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Gaussian relaxation value at time t */ Double_t PTheory::SimpleGauss(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { @@ -1194,17 +1359,33 @@ Double_t PTheory::SimpleGauss(Double_t t, const PDoubleVector& paramValues, cons //-------------------------------------------------------------------------- /** - *

theory function: static Gaussian Kubo-Toyabe in zero applied field - * \f[ = \frac{1}{3} + \frac{2}{3} \left[1-(\sigma t)^2\right] \exp\left[-1/2 (\sigma t)^2\right] \f] or - * \f[ = \frac{1}{3} + \frac{2}{3} \left[1-(\sigma \{t-t_{\rm shift}\})^2\right] \exp\left[-1/2 (\sigma \{t-t_{\rm shift}\})^2\right] \f] + * \brief Static Gaussian Kubo-Toyabe relaxation function in zero field. * - * meaning of paramValues: \f$\sigma\f$ [, \f$t_{\rm shift}\f$] + * Mathematical form: + * \f[ P(t) = \frac{1}{3} + \frac{2}{3} \left[1-(\sigma t)^2\right] \exp\left[-\frac{1}{2} (\sigma t)^2\right] \f] + * or with time shift: + * \f[ P(t) = \frac{1}{3} + \frac{2}{3} \left[1-(\sigma \{t-t_{\rm shift}\})^2\right] \exp\left[-\frac{1}{2} (\sigma \{t-t_{\rm shift}\})^2\right] \f] * - * return: function value + * The Kubo-Toyabe function describes muon spin relaxation in a powder sample + * with random static nuclear magnetic moments creating a Gaussian field + * distribution at the muon site. * - * \param t time in \f$(\mu\mathrm{s})\f$, or x-axis value for non-muSR fit - * \param paramValues parameter values - * \param funcValues vector with the functions (i.e. functions of the parameters) + * Physical interpretation: + * - 1/3 component: Muons with spin parallel to local field (unrelaxed) + * - 2/3 component: Muons with spin perpendicular to local field (relaxing) + * - Recovery at long times: Partial repolarization to 1/3 + * + * Parameters: + * - fParamNo[0]: Relaxation rate σ = γ_μ ΔB (μs⁻¹), where ΔB is RMS field + * - fParamNo[1]: (optional) Time shift t_shift (μs) + * + * \param t Time in μs + * \param paramValues Vector of fit parameter values + * \param funcValues Vector of evaluated function values + * + * \return Gaussian Kubo-Toyabe value at time t + * + * \see StaticGaussKTLF() for longitudinal field version */ Double_t PTheory::StaticGaussKT(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const { diff --git a/src/include/PTheory.h b/src/include/PTheory.h index 3c5b98ae9..3bc1331f3 100644 --- a/src/include/PTheory.h +++ b/src/include/PTheory.h @@ -344,135 +344,315 @@ static PTheoDataBase fgTheoDataBase[THEORY_MAX] = { //-------------------------------------------------------------------------------------- /** - *

Theory function evaluator and manager. + * \brief Theory function evaluator and expression tree manager. * - *

This class parses, validates, and evaluates theory functions specified - * in the THEORY block of MSR files. It handles: - * - Parsing theory syntax (function names, parameters, operators) - * - Building theory expression trees (addition/multiplication) - * - Resolving parameter mappings and user functions - * - Evaluating theory at specific time points - * - Caching longitudinal field calculations - * - Loading and managing user-defined shared library functions + * PTheory is the core class responsible for parsing, validating, and evaluating + * theory functions specified in the THEORY block of MSR files. It implements + * a binary expression tree to handle complex combinations of theory functions. * - *

Theory syntax: - * - Single function: asymmetry simpleExp 1 3 - * - Addition: func1 par1... + func2 par1... - * - Multiplication: func1 par1... * func2 par1... - * - Nested: (func1 + func2) * func3 + * \section theory_overview Overview * - *

Parameter resolution: - * Parameters can be direct numbers, parameter references (1, 2, 3, ...), - * map references (map1, map2, ...), or function references (fun1, fun2, ...). + * The theory describes the expected muon polarization as a function of time. + * Different physical phenomena (relaxation, precession, field distributions) + * are represented by different theory functions that can be combined. * - *

Example: - * @code + * \section theory_tree Expression Tree Structure + * + * Theory expressions are parsed into a binary tree where: + * - Each node represents a theory function + * - Left child (fAdd) represents addition (+) + * - Right child (fMul) represents multiplication (*) + * + * Example MSR theory block: + * \verbatim * THEORY - * asymmetry 1 - * TFieldCos 2 3 (phase, freq) - * + - * simpleExp 4 (rate) - * @endcode - * This evaluates to: par1 * cos(par2 + 2π*par3*t) + exp(-par4*t) + * a 1 # asymmetry + * tf 2 3 # TFieldCos + * se 4 # simpleExp + * + + * a 5 + * tf 6 7 + * \endverbatim + * + * Becomes: (par1 * cos(φ₂ + 2πν₃t) * exp(-λ₄t)) + (par5 * cos(φ₆ + 2πν₇t)) + * + * \section theory_syntax Syntax Details + * + * Function specification: + * \code + * function_name param1 param2 ... paramN (optional comment) + * \endcode + * + * Parameter types: + * - Direct number: \c 1, \c 2, \c 3 → parameter index (1-based) + * - Map reference: \c map1, \c map2 → indirection via RUN block map + * - Function reference: \c fun1, \c fun2 → evaluated FUNCTIONS block entry + * + * Operators: + * - \c + on separate line: Addition of following terms + * - Consecutive lines without +: Implicit multiplication + * + * \section theory_functions Available Functions + * + * Basic: + * - \c const (c): Constant value + * - \c asymmetry (a): Initial asymmetry + * - \c simplExpo (se): exp(-λt) + * - \c generExpo (ge): exp(-(λt)^β) + * - \c simpleGss (sg): exp(-σ²t²/2) + * + * Kubo-Toyabe (Gaussian): + * - \c statGssKt (stg): Static ZF Gaussian KT + * - \c statGssKTLF (sgktlf): Static LF Gaussian KT + * - \c dynGssKTLF (dgktlf): Dynamic LF Gaussian KT + * + * Kubo-Toyabe (Lorentzian): + * - \c statExpKT (sekt): Static ZF Lorentzian KT + * - \c statExpKTLF (sektlf): Static LF Lorentzian KT + * - \c dynExpKTLF (dektlf): Dynamic LF Lorentzian KT + * + * Precession: + * - \c TFieldCos (tf): cos(φ + 2πνt) + * - \c bessel (b): J₀(2πνt + φ) + * - \c internFld (ifld): Internal field distribution + * + * Special: + * - \c spinGlass (spg): Spin glass relaxation + * - \c abragam (ab): Motional narrowing + * - \c userFcn (u): User-defined external function + * + * \section theory_lf Longitudinal Field Calculations + * + * LF Kubo-Toyabe functions require numerical integration which is cached + * for efficiency. The cache is invalidated when parameters change. + * Caching variables: + * - fPrevParam: Previous parameter values for change detection + * - fLFIntegral: Cached integral values + * - fDynLFFuncValue: Dynamic KT function cache + * + * \section theory_user User Functions + * + * External functions can be loaded from shared libraries: + * \code + * userFcn libMyFunctions.so MyFunctionClass param1 param2 ... + * \endcode + * + * User functions must inherit from PUserFcnBase and implement: + * - operator()(Double_t t, const std::vector& par) + * + * \see PUserFcnBase, PMsrHandler, fgTheoDataBase */ class PTheory { public: /** - *

Constructor for theory handler. + * \brief Constructor that parses the THEORY block and builds the expression tree. * - * @param msrInfo Pointer to MSR file handler - * @param runNo Run number for this theory (0-based) - * @param hasParent True if this is a sub-theory (for + / * operators) + * Parses the theory block from the MSR file, validates function names and + * parameter counts, resolves parameter references, and recursively builds + * the expression tree for operators (+ and *). + * + * Parsing algorithm: + * -# Get current line from theory block + * -# Remove comments (text after '(' or '#') + * -# Tokenize to extract function name and parameters + * -# Look up function in fgTheoDataBase + * -# Validate parameter count + * -# Resolve parameter references (direct, map, fun) + * -# If next line is not '+', recursively create fMul child + * -# If next line is '+', skip it and recursively create fAdd child + * -# For user functions: load shared library and instantiate class + * + * \param msrInfo Pointer to MSR file handler containing theory block + * \param runNo Run number (0-based) for parameter map resolution + * \param hasParent False for root theory, true for child nodes. + * Controls static counter reset for recursive parsing. + * + * \warning If parsing fails, fValid is set to false and error messages + * are output to stderr. Always check IsValid() after construction. */ PTheory(PMsrHandler *msrInfo, UInt_t runNo, const Bool_t hasParent = false); + /** + * \brief Destructor that recursively cleans up the expression tree. + * + * Releases all allocated resources: + * - Clears parameter and function value vectors + * - Clears LF integral caches + * - Recursively deletes child theory nodes (fAdd, fMul) + * - Deletes user function object if present + * - Clears global user function vector + */ virtual ~PTheory(); - /// Returns true if theory is valid and ready for evaluation + /** + * \brief Checks if the entire theory expression tree is valid. + * + * Recursively validates all nodes in the tree. A theory is valid only + * if this node and all its children (fAdd and fMul) are valid. + * + * \return true if all nodes are valid, false if any node has errors + * + * \note Call this after construction to verify the theory was parsed + * successfully before attempting evaluation. + */ virtual Bool_t IsValid(); /** - *

Evaluates the theory function at time t. + * \brief Evaluates the theory function at a given time point. * - * @param t Time in microseconds (μs) - * @param paramValues Vector of fit parameter values - * @param funcValues Vector of evaluated user function values - * @return Theory value at time t + * Recursively evaluates the expression tree at time t using the provided + * parameter values. The evaluation follows the tree structure: + * - If both fMul and fAdd exist: this_func * fMul + fAdd + * - If only fMul exists: this_func * fMul + * - If only fAdd exists: this_func + fAdd + * - If neither exists: this_func + * + * \param t Time in microseconds (μs) for μSR fits, or x-value for non-μSR + * \param paramValues Vector of current fit parameter values (FITPARAMETER block) + * \param funcValues Vector of evaluated function values (FUNCTIONS block) + * + * \return Evaluated theory value (polarization or asymmetry) at time t + * + * \note For thread-safety, LF calculations cache results which may need + * recalculation if parameters change between calls. */ virtual Double_t Func(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; private: + /** \brief Recursively deletes child theory nodes (fAdd and fMul). */ virtual void CleanUp(PTheory *theo); + + /** \brief Searches fgTheoDataBase for a function by name or abbreviation. + * \param name Function name (e.g., "simplExpo") or abbreviation (e.g., "se") + * \return Index in fgTheoDataBase, or THEORY_UNDEFINED if not found */ virtual Int_t SearchDataBase(TString name); + + /** \brief Returns the index of user functions up to the given line. + * \param lineNo Current line number in theory block + * \return Count of userFcn entries before this line */ virtual Int_t GetUserFcnIdx(UInt_t lineNo) const; + + /** \brief Reformats the theory block for clean MSR file output. */ virtual void MakeCleanAndTidyTheoryBlock(PMsrLines* fullTheoryBlock); + + /** \brief Formats a polynomial theory line with proper spacing. */ virtual void MakeCleanAndTidyPolynom(UInt_t i, PMsrLines* fullTheoryBlock); + + /** \brief Formats a user function theory line with proper spacing. */ virtual void MakeCleanAndTidyUserFcn(UInt_t i, PMsrLines* fullTheoryBlock); + // -------------------- Theory Function Implementations -------------------- + // Each function evaluates its specific physical model at time t. + // Parameters are resolved from fParamNo using paramValues and funcValues. + + /** \brief Returns constant value. Formula: c */ virtual Double_t Constant(const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Returns asymmetry value. Formula: A */ virtual Double_t Asymmetry(const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Simple exponential relaxation. Formula: exp(-λt) */ virtual Double_t SimpleExp(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief General (stretched) exponential. Formula: exp(-(λt)^β) */ virtual Double_t GeneralExp(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Simple Gaussian relaxation. Formula: exp(-σ²t²/2) */ virtual Double_t SimpleGauss(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Gaussian Kubo-Toyabe (ZF). Formula: 1/3 + 2/3(1-σ²t²)exp(-σ²t²/2) */ virtual Double_t StaticGaussKT(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Gaussian Kubo-Toyabe (LF). Requires numerical integration. */ virtual Double_t StaticGaussKTLF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Dynamic Gaussian Kubo-Toyabe (LF). Strong collision model. */ virtual Double_t DynamicGaussKTLF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Lorentzian Kubo-Toyabe (ZF). Formula: 1/3 + 2/3(1-at)exp(-at) */ virtual Double_t StaticLorentzKT(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Lorentzian Kubo-Toyabe (LF). Requires numerical integration. */ virtual Double_t StaticLorentzKTLF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Dynamic Lorentzian Kubo-Toyabe (LF). Strong collision model. */ virtual Double_t DynamicLorentzKTLF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Fast dynamic Gaussian-Lorentzian KT (ZF). Approximate fast calculation. */ virtual Double_t DynamicGauLorKTZFFast(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Fast dynamic Gaussian-Lorentzian KT (LF). Approximate fast calculation. */ virtual Double_t DynamicGauLorKTLFFast(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Dynamic Gaussian-Lorentzian KT (LF). Full numerical calculation. */ virtual Double_t DynamicGauLorKTLF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Combined Lorentzian-Gaussian KT. Product of both relaxation types. */ virtual Double_t CombiLGKT(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Stretched Kubo-Toyabe. Formula: exp(-(σt)^β) with KT-like recovery. */ virtual Double_t StrKT(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Spin glass relaxation function. Edwards-Anderson order parameter. */ virtual Double_t SpinGlass(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Random anisotropic hyperfine coupling. Powder average of anisotropic coupling. */ virtual Double_t RandomAnisotropicHyperfine(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Abragam relaxation. Motional narrowing formula. */ virtual Double_t Abragam(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Transverse field cosine. Formula: cos(φ + 2πνt) */ virtual Double_t TFCos(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Internal field distribution. Gaussian field distribution model. */ virtual Double_t InternalField(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Internal field (Kornilov model). Vortex lattice field distribution. */ virtual Double_t InternalFieldGK(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Internal field (Larkin-Ovchinnikov model). Vortex lattice field distribution. */ virtual Double_t InternalFieldLL(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Bessel function precession. Formula: J₀(2πνt + φ) */ virtual Double_t Bessel(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Internal Bessel field distribution. Combines Bessel with relaxation. */ virtual Double_t InternalBessel(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Skewed Gaussian. Asymmetric relaxation rates before/after zero crossing. */ virtual Double_t SkewedGauss(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Nakajima-Keren (ZF). Combined nuclear and electronic relaxation. */ virtual Double_t StaticNKZF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Static Nakajima-Keren (TF). Combined nuclear and electronic relaxation with precession. */ virtual Double_t StaticNKTF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Dynamic Nakajima-Keren (ZF). With spin fluctuations. */ virtual Double_t DynamicNKZF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Dynamic Nakajima-Keren (TF). With spin fluctuations and precession. */ virtual Double_t DynamicNKTF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief F-μ-F oscillation. Muon bound between two fluorine atoms. */ virtual Double_t FmuF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief μ⁻ exponential TF. Negative muon in transverse field. */ virtual Double_t MuMinusExpTF(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief Polynomial function. Formula: Σᵢ pᵢtⁱ */ virtual Double_t Polynom(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + /** \brief User-defined function. Calls external shared library function. */ virtual Double_t UserFcn(Double_t t, const PDoubleVector& paramValues, const PDoubleVector& funcValues) const; + // -------------------- LF Calculation Helpers -------------------- + /** \brief Calculates and caches Gaussian LF integral for static KT. */ virtual void CalculateGaussLFIntegral(const Double_t *val) const; + /** \brief Calculates and caches Lorentzian LF integral for static KT. */ virtual void CalculateLorentzLFIntegral(const Double_t *val) const; + /** \brief Retrieves cached LF integral value at time t using interpolation. */ virtual Double_t GetLFIntegralValue(const Double_t t) const; + /** \brief Calculates dynamic KT in LF using integral equation approach. */ virtual void CalculateDynKTLF(const Double_t *val, Int_t tag) const; + /** \brief Retrieves cached dynamic KT LF value at time t. */ virtual Double_t GetDynKTLFValue(const Double_t t) const; + /** \brief Retrieves cached dynamic Gauss-Lorentz KT LF value at time t. */ virtual Double_t GetDyn_GL_KTLFValue(const Double_t t) const; - // variables - Bool_t fValid; ///< flag to tell if the theory is valid - UInt_t fType; ///< function tag - std::vector fParamNo; ///< holds the parameter numbers for the theory (including maps and functions, see constructor desciption) - UInt_t fNoOfParam; ///< number of parameters for the given function - PTheory *fAdd, *fMul; ///< pointers to the add-sub-function or the multiply-sub-function + // -------------------- Member Variables -------------------- + Bool_t fValid; ///< True if this theory node and its parse state are valid + UInt_t fType; ///< Theory function type (THEORY_CONST, THEORY_SIMPLE_EXP, etc.) + std::vector fParamNo; ///< Resolved parameter indices (0-based). Values >= MSR_PARAM_FUN_OFFSET are function references. + UInt_t fNoOfParam; ///< Expected number of parameters for this function type + PTheory *fAdd; ///< Pointer to addition child node (left branch of tree) + PTheory *fMul; ///< Pointer to multiplication child node (right branch of tree) - Int_t fUserFcnIdx; ///< index of the user function within the theory tree - TString fUserFcnClassName; ///< name of the user function class for within root - TString fUserFcnSharedLibName; ///< name of the shared lib to which the user function belongs - PUserFcnBase *fUserFcn; ///< pointer to the user function object - mutable PDoubleVector fUserParam; ///< vector holding the resolved user function parameters, i.e. map and function resolved. + // User function members + Int_t fUserFcnIdx; ///< Index of this user function among all userFcn entries (for global state) + TString fUserFcnClassName; ///< ROOT class name for user function (e.g., "TMyFunction") + TString fUserFcnSharedLibName; ///< Shared library path (e.g., "libMyFunctions.so") + PUserFcnBase *fUserFcn; ///< Pointer to instantiated user function object + mutable PDoubleVector fUserParam; ///< Resolved parameter values for user function calls - PMsrHandler *fMsrInfo; ///< pointer to the msr-file handler + PMsrHandler *fMsrInfo; ///< Pointer to MSR file handler (not owned) - mutable Double_t fSamplingTime; ///< needed for LF. Keeps the sampling time of the non-analytic integral - mutable Double_t fPrevParam[THEORY_MAX_PARAM]; ///< needed for LF. Keeps the previous fitting parameters - mutable PDoubleVector fLFIntegral; ///< needed for LF. Keeps the non-analytic integral values - mutable Double_t fDynLFdt; ///< needed for LF. Keeps the time step for the dynamic LF calculation, used in the integral equation approach - mutable PDoubleVector fDynLFFuncValue; ///< needed for LF. Keeps the dynamic LF KT function values - mutable PDoubleVector fDyn_GL_LFFuncValue; ///< needed for LF. Keeps the dynamic LF local Gauss/global Lorentzian KT function values + // LF calculation caching (mutable for const Func() calls) + mutable Double_t fSamplingTime; ///< Time step for LF integral calculation (default 1 ns = 0.001 μs) + mutable Double_t fPrevParam[THEORY_MAX_PARAM]; ///< Previous parameter values for cache invalidation check + mutable PDoubleVector fLFIntegral; ///< Cached static LF KT integral values + mutable Double_t fDynLFdt; ///< Time step for dynamic LF integral equation + mutable PDoubleVector fDynLFFuncValue; ///< Cached dynamic Gaussian/Lorentzian LF KT values + mutable PDoubleVector fDyn_GL_LFFuncValue; ///< Cached dynamic Gauss-Lorentz LF KT values }; #endif // _PTHEORY_H_