From 7b82b2c2c8dc2826f885d4e4f8e6ce7e64b1b7a8 Mon Sep 17 00:00:00 2001 From: Andreas Suter Date: Sat, 15 Nov 2025 07:43:49 +0100 Subject: [PATCH] improve the doxygen docu of PRunAsymmetry.*. --- src/classes/PRunAsymmetry.cpp | 360 ++++++++++++++++++++++++---------- src/include/PRunAsymmetry.h | 217 ++++++++++++++++++-- 2 files changed, 456 insertions(+), 121 deletions(-) diff --git a/src/classes/PRunAsymmetry.cpp b/src/classes/PRunAsymmetry.cpp index 8ffe3fdbf..3d55eabbc 100644 --- a/src/classes/PRunAsymmetry.cpp +++ b/src/classes/PRunAsymmetry.cpp @@ -50,7 +50,11 @@ // Constructor //-------------------------------------------------------------------------- /** - *

Constructor + * \brief Default constructor that initializes all member variables. + * + * Sets all counters and indices to default/invalid values. This constructor + * creates an invalid instance that requires proper initialization via the + * main constructor. */ PRunAsymmetry::PRunAsymmetry() : PRunBase() { @@ -71,12 +75,26 @@ PRunAsymmetry::PRunAsymmetry() : PRunBase() // Constructor //-------------------------------------------------------------------------- /** - *

Constructor + * \brief Main constructor that initializes μSR asymmetry fitting. * - * \param msrInfo pointer to the msr-file handler - * \param rawData raw run data - * \param runNo number of the run within the msr-file - * \param tag tag showing what shall be done: kFit == fitting, kView == viewing + * Performs comprehensive initialization: + * 1. Validates packing parameter from RUN or GLOBAL block + * 2. Validates α parameter (required for asymmetry) + * 3. Validates β parameter (optional, defaults to 1) + * 4. Determines α/β configuration tag (1-4) + * 5. Calls PrepareData() to load and process histogram data + * + * The α/β tag determines the asymmetry calculation method: + * - Tag 1: α=1, β=1 (no corrections) + * - Tag 2: α≠1, β=1 (forward/backward correction) + * - Tag 3: α=1, β≠1 (alternative correction) + * - Tag 4: α≠1, β≠1 (both corrections) + * + * \param msrInfo Pointer to MSR file handler + * \param rawData Pointer to raw run data handler + * \param runNo Run number within the MSR file + * \param tag Operation mode (kFit or kView) + * \param theoAsData If true, calculate theory only at data points */ PRunAsymmetry::PRunAsymmetry(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag, Bool_t theoAsData) : PRunBase(msrInfo, rawData, runNo, tag), fTheoAsData(theoAsData) @@ -173,7 +191,10 @@ PRunAsymmetry::PRunAsymmetry(PMsrHandler *msrInfo, PRunDataHandler *rawData, UIn // Destructor //-------------------------------------------------------------------------- /** - *

Destructor. + * \brief Destructor that cleans up histogram data. + * + * Clears all four histogram vectors (forward/backward × data/errors) + * to free memory. */ PRunAsymmetry::~PRunAsymmetry() { @@ -187,12 +208,21 @@ PRunAsymmetry::~PRunAsymmetry() // CalcChiSquare (public) //-------------------------------------------------------------------------- /** - *

Calculate chi-square. + * \brief Calculates chi-square for μSR asymmetry fit. * - * return: - * - chisq value + * Computes χ² by comparing the measured asymmetry with theory: + * χ² = Σ[(A_data - A_theory)²/σ²] * - * \param par parameter vector iterated by minuit2 + * The asymmetry calculation depends on fAlphaBetaTag: + * - Tag 1 (α=β=1): A = (F - B)/(F + B) + * - Tag 2 (α≠1, β=1): A = (F - αB)/(F + αB) + * - Tag 3 (α=1, β≠1): Uses β correction + * - Tag 4 (α≠1, β≠1): Uses both α and β corrections + * + * Supports OpenMP parallelization for faster calculation. + * + * \param par Parameter vector from MINUIT minimizer + * \return Chi-square value */ Double_t PRunAsymmetry::CalcChiSquare(const std::vector& par) { @@ -289,12 +319,17 @@ Double_t PRunAsymmetry::CalcChiSquare(const std::vector& par) // CalcChiSquareExpected (public) //-------------------------------------------------------------------------- /** - *

Calculate expected chi-square. Currently not implemented since not clear what to be done. + * \brief Calculates expected chi-square value. * - * return: - * - chisq value == 0.0 + * This function is currently not implemented for asymmetry fits because the expected + * chi-square calculation for asymmetry data requires a more complex statistical treatment + * than for single histogram fits. The asymmetry is a ratio of count rates, and the + * proper expected value calculation is non-trivial. * - * \param par parameter vector iterated by minuit2 + * \param par Parameter vector from MINUIT minimizer + * \return Always returns 0.0 (placeholder value) + * + * \todo Implement proper expected chi-square calculation for asymmetry fits */ Double_t PRunAsymmetry::CalcChiSquareExpected(const std::vector& par) { @@ -305,9 +340,16 @@ Double_t PRunAsymmetry::CalcChiSquareExpected(const std::vector& par) // CalcMaxLikelihood (public) //-------------------------------------------------------------------------- /** - *

NOT IMPLEMENTED!! + * \brief Calculates maximum likelihood estimator for asymmetry fit. * - * \param par parameter vector iterated by minuit2 + * Maximum likelihood estimation provides an alternative to χ² minimization and can be + * more appropriate for low-count data where Gaussian statistics do not apply. However, + * the proper likelihood function for μSR asymmetry data is complex and not yet implemented. + * + * \param par Parameter vector from MINUIT minimizer + * \return Placeholder value of 1.0 + * + * \todo Implement proper Poisson-based maximum likelihood for asymmetry fits */ Double_t PRunAsymmetry::CalcMaxLikelihood(const std::vector& par) { @@ -320,9 +362,13 @@ Double_t PRunAsymmetry::CalcMaxLikelihood(const std::vector& par) // GetNoOfFitBins (public) //-------------------------------------------------------------------------- /** - *

Calculate the number of fitted bins for the current fit range. + * \brief Returns the number of bins included in the fit. * - * return: number of fitted bins. + * Calculates and returns the number of data bins that fall within the current fit range. + * This value is used for determining degrees of freedom and is recalculated when the + * fit range changes (e.g., via COMMAND block). + * + * \return Number of data bins included in the fit */ UInt_t PRunAsymmetry::GetNoOfFitBins() { @@ -335,15 +381,26 @@ UInt_t PRunAsymmetry::GetNoOfFitBins() // SetFitRangeBin (public) //-------------------------------------------------------------------------- /** - *

Allows to change the fit range on the fly. Used in the COMMAND block. - * The syntax of the string is: FIT_RANGE fgb[+n00] lgb[-n01] [fgb[+n10] lgb[-n11] ... fgb[+nN0] lgb[-nN1]]. - * If only one pair of fgb/lgb is given, it is used for all runs in the RUN block section. - * If multiple fgb/lgb's are given, the number N has to be the number of RUN blocks in - * the msr-file. + * \brief Dynamically changes the fit range in bin units. * - *

nXY are offsets which can be used to shift, limit the fit range. + * Allows modification of the fit range at runtime, typically called from the COMMAND block. + * Supports bin-based fit ranges with optional offsets for fine-tuning. * - * \param fitRange string containing the necessary information. + * Syntax formats: + * - Single range: FIT_RANGE fgb[+n0] lgb[-n1] + * - Applied to all RUN blocks + * - Multiple ranges: FIT_RANGE fgb[+n00] lgb[-n01] fgb[+n10] lgb[-n11] ... fgb[+nN0] lgb[-nN1] + * - One pair per RUN block (N must equal number of RUN blocks) + * + * Parameters: + * - fgb: First good bin (start of fit range) + * - lgb: Last good bin (end of fit range) + * - +n: Positive offset to shift start forward + * - -n: Negative offset to shift end backward + * + * Example: "FIT_RANGE 10+5 200-10" uses bins [15, 190] for fitting + * + * \param fitRange String containing fit range specification */ void PRunAsymmetry::SetFitRangeBin(const TString fitRange) { @@ -425,7 +482,16 @@ void PRunAsymmetry::SetFitRangeBin(const TString fitRange) // CalcNoOfFitBins (public) //-------------------------------------------------------------------------- /** - *

Calculate the number of fitted bins for the current fit range. + * \brief Calculates the number of bins included in the current fit range. + * + * Determines fStartTimeBin and fEndTimeBin from the fit time range (fFitStartTime, fFitEndTime) + * and data time grid. Ensures that bin indices remain within valid histogram bounds. + * The result is stored in fNoOfFitBins. + * + * This calculation accounts for: + * - Data time step and start time + * - Rounding effects (ceiling for start, floor for end) + * - Boundary conditions (clips to [0, histogram size]) */ void PRunAsymmetry::CalcNoOfFitBins() { @@ -447,7 +513,20 @@ void PRunAsymmetry::CalcNoOfFitBins() // CalcTheory (protected) //-------------------------------------------------------------------------- /** - *

Calculate theory for a given set of fit-parameters. + * \brief Calculates theoretical asymmetry values for the current parameters. + * + * Computes the expected asymmetry A(t) for all data points based on the current + * parameter values and the user-defined theory function. The calculation depends + * on the α/β correction mode: + * + * - Tag 1 (α=1, β=1): A(t) = f(t) + * - Tag 2 (α≠1, β=1): A(t) = [f(t)(α+1) - (α-1)] / [(α+1) - f(t)(α-1)] + * - Tag 3 (α=1, β≠1): A(t) = f(t)(β+1) / [2 - f(t)(β-1)] + * - Tag 4 (α≠1, β≠1): Combined α and β corrections + * + * where f(t) is the raw theory function from the THEORY block. + * + * The calculated values are stored in fData for plotting and comparison with measured data. */ void PRunAsymmetry::CalcTheory() { @@ -531,23 +610,31 @@ void PRunAsymmetry::CalcTheory() // PrepareData (protected) //-------------------------------------------------------------------------- /** - *

Prepare data for fitting or viewing. What is already processed at this stage: - * - get all needed forward/backward histograms - * - get time resolution - * - get start/stop fit time - * - get t0's and perform necessary cross checks (e.g. if t0 of msr-file (if present) are consistent with t0 of the data files, etc.) - * - add runs (if addruns are present) - * - group histograms (if grouping is present) - * - subtract background + * \brief Main data preparation routine for asymmetry fitting and viewing. * - * Error propagation for \f$ A_i = (f_i^{\rm c}-b_i^{\rm c})/(f_i^{\rm c}+b_i^{\rm c})\f$: - * \f[ \Delta A_i = \pm\frac{2}{(f_i^{\rm c}+b_i^{\rm c})^2}\left[ - * (b_i^{\rm c})^2 (\Delta f_i^{\rm c})^2 + - * (\Delta b_i^{\rm c})^2 (f_i^{\rm c})^2\right]^{1/2}\f] + * Performs comprehensive data preparation including: + * - Loading forward/backward histograms from data files + * - Extracting metadata (field, energy, temperature) + * - Determining time resolution from data file + * - Validating and retrieving t0 values for all histograms + * - Adding multiple runs together (if addruns are specified) + * - Grouping multiple histograms (if grouping is specified) + * - Subtracting background (fixed or estimated) + * - Calculating asymmetry and error bars + * - Applying bin packing (if specified) * - * return: - * - true if everthing went smooth - * - false, otherwise. + * The asymmetry is calculated as: + * \f[ A_i = \frac{f_i^{\rm c} - b_i^{\rm c}}{f_i^{\rm c} + b_i^{\rm c}} \f] + * + * Error propagation (assuming Poisson statistics): + * \f[ \Delta A_i = \pm\frac{2}{(f_i^{\rm c}+b_i^{\rm c})^2}\sqrt{ + * (b_i^{\rm c})^2 (\Delta f_i^{\rm c})^2 + + * (f_i^{\rm c})^2 (\Delta b_i^{\rm c})^2} \f] + * + * where \f$ f_i^{\rm c} \f$ and \f$ b_i^{\rm c} \f$ are background-corrected forward and + * backward histograms, respectively. + * + * \return True if data preparation succeeds, false on any error */ Bool_t PRunAsymmetry::PrepareData() { @@ -942,15 +1029,24 @@ Bool_t PRunAsymmetry::SubtractEstimatedBkg() // PrepareFitData (protected) //-------------------------------------------------------------------------- /** - *

Take the pre-processed data (i.e. grouping and addrun are preformed) and form the asymmetry for fitting. - * Before forming the asymmetry, the following checks will be performed: - * -# check if data range is given, if not try to estimate one. - * -# check that if a data range is present, that it makes any sense. - * -# check that 'first good bin'-'t0' is the same for forward and backward histogram. If not adjust it. - * -# pack data (rebin). - * -# if packed forward size != backward size, truncate the longer one such that an asymmetry can be formed. - * -# calculate the asymmetry: \f$ A_i = (f_i^c-b_i^c)/(f_i^c+b_i^c) \f$ - * -# calculate the asymmetry errors: \f$ \delta A_i = 2 \sqrt{(b_i^c)^2 (\delta f_i^c)^2 + (\delta b_i^c)^2 (f_i^c)^2}/(f_i^c+b_i^c)^2\f$ + * \brief Processes pre-grouped data and calculates asymmetry for fitting. + * + * Takes forward/backward histograms (after grouping and addrun operations) and performs + * final processing steps for fitting: + * + * 1. Validates data range, estimates if not specified + * 2. Checks data range consistency and validity + * 3. Aligns forward/backward histograms (ensures 'first good bin - t0' is identical) + * 4. Applies bin packing (rebinning) if requested + * 5. Truncates longer histogram if packed sizes differ + * 6. Calculates asymmetry: \f$ A_i = (f_i^{\rm c}-b_i^{\rm c})/(f_i^{\rm c}+b_i^{\rm c}) \f$ + * 7. Propagates errors: \f$ \delta A_i = \frac{2}{(f_i^{\rm c}+b_i^{\rm c})^2}\sqrt{(b_i^{\rm c})^2 (\delta f_i^{\rm c})^2 + (f_i^{\rm c})^2 (\delta b_i^{\rm c})^2} \f$ + * + * Bin packing averages multiple bins to improve statistics: + * - Packed value normalized to per-bin counts (value/packing) + * - Error propagation: \f$ \sigma_{\rm packed} = \sqrt{\sum \sigma_i^2}/N_{\rm pack} \f$ + * + * \return True on success, false if data preparation fails */ Bool_t PRunAsymmetry::PrepareFitData() { @@ -1056,20 +1152,25 @@ Bool_t PRunAsymmetry::PrepareFitData() // PrepareViewData (protected) //-------------------------------------------------------------------------- /** - *

Take the pre-processed data (i.e. grouping and addrun are preformed) and form the asymmetry for view representation. - * Before forming the asymmetry, the following checks will be performed: - * -# check if view packing is whished. - * -# check if data range is given, if not try to estimate one. - * -# check that data range is present, that it makes any sense. - * -# check that 'first good bin'-'t0' is the same for forward and backward histogram. If not adjust it. - * -# pack data (rebin). - * -# if packed forward size != backward size, truncate the longer one such that an asymmetry can be formed. - * -# calculate the asymmetry: \f$ A_i = (\alpha f_i^c-b_i^c)/(\alpha \beta f_i^c+b_i^c) \f$ - * -# calculate the asymmetry errors: \f$ \delta A_i = 2 \sqrt{(b_i^c)^2 (\delta f_i^c)^2 + (\delta b_i^c)^2 (f_i^c)^2}/(f_i^c+b_i^c)^2\f$ - * -# calculate the theory vector. + * \brief Prepares asymmetry data for plotting and visualization. * - * \param runData raw run data needed to perform some crosschecks - * \param histoNo histogram number (within a run). histoNo[0]: forward histogram number, histNo[1]: backward histogram number + * Processes pre-grouped data for display in plots, with special handling for view packing + * and α/β corrections. Similar to PrepareFitData but includes theory calculation and + * supports separate view packing settings. + * + * Processing steps: + * 1. Checks for view-specific packing (overrides fit packing for display) + * 2. Validates and estimates data range if needed + * 3. Aligns forward/backward histogram start bins relative to t0 + * 4. Applies bin packing for improved visualization + * 5. Ensures equal bin counts between forward/backward after packing + * 6. Calculates asymmetry with α/β corrections: \f$ A_i = (\alpha f_i^{\rm c} - b_i^{\rm c})/(\alpha \beta f_i^{\rm c} + b_i^{\rm c}) \f$ + * 7. Propagates errors: \f$ \delta A_i = \frac{2}{(f_i^{\rm c}+b_i^{\rm c})^2}\sqrt{(b_i^{\rm c})^2 (\delta f_i^{\rm c})^2 + (f_i^{\rm c})^2 (\delta b_i^{\rm c})^2} \f$ + * 8. Calculates theory curve for overlay + * + * \param runData Pointer to raw run data for validation and cross-checks + * \param histoNo Array of histogram indices: [0]=forward, [1]=backward + * \return True on success, false on error */ Bool_t PRunAsymmetry::PrepareViewData(PRawRunData* runData, UInt_t histoNo[2]) { @@ -1341,19 +1442,27 @@ Bool_t PRunAsymmetry::PrepareViewData(PRawRunData* runData, UInt_t histoNo[2]) // PrepareRRFViewData (protected) //-------------------------------------------------------------------------- /** - *

Prepares the RRF data set for visual representation. This is done the following way: - * -# make all necessary checks - * -# build the asymmetry, \f$ A(t) \f$, WITHOUT packing. - * -# \f$ A_R(t) = A(t) \cdot 2 \cos(\omega_R t + \phi_R) \f$ - * -# do the packing of \f$ A_R(t) \f$ - * -# calculate theory, \f$ T(t) \f$, as close as possible to the time resolution [compatible with the RRF frequency] - * -# \f$ T_R(t) = T(t) \cdot 2 \cos(\omega_R t + \phi_R) \f$ - * -# do the packing of \f$ T_R(t) \f$ - * -# calculate the Kaiser FIR filter coefficients - * -# filter \f$ T_R(t) \f$. + * \brief Prepares rotating reference frame (RRF) data for visualization. * - * \param runData raw run data needed to perform some crosschecks - * \param histoNo array of the histo numbers form which to build the asymmetry + * Transforms asymmetry data into a rotating reference frame for analyzing high-frequency + * oscillations. The RRF technique mixes the data with a reference frequency to shift + * oscillations down to lower frequencies, making them easier to visualize and analyze. + * + * Processing sequence: + * 1. Validates data ranges and histogram alignment + * 2. Builds asymmetry \f$ A(t) \f$ without packing + * 3. Applies RRF transformation: \f$ A_R(t) = A(t) \cdot 2\cos(\omega_R t + \phi_R) \f$ + * 4. Packs the RRF asymmetry data + * 5. Calculates theory \f$ T(t) \f$ at high time resolution + * 6. Transforms theory to RRF: \f$ T_R(t) = T(t) \cdot 2\cos(\omega_R t + \phi_R) \f$ + * 7. Packs the RRF theory curve + * 8. Applies Kaiser FIR filter to smooth the theory curve + * + * The RRF frequency (\f$ \omega_R \f$) and phase (\f$ \phi_R \f$) are specified in the PLOT block. + * + * \param runData Pointer to raw run data for validation + * \param histoNo Array of histogram indices: [0]=forward, [1]=backward + * \return True on success, false on error */ Bool_t PRunAsymmetry::PrepareRRFViewData(PRawRunData* runData, UInt_t histoNo[2]) { @@ -1675,22 +1784,32 @@ Bool_t PRunAsymmetry::PrepareRRFViewData(PRawRunData* runData, UInt_t histoNo[2] // GetProperT0 (private) //-------------------------------------------------------------------------- /** - *

Get the proper t0 for the single histogram run. - * -# the t0 vector size = number of detectors (grouping) for forward + backward. - * -# initialize t0's with -1 - * -# fill t0's from RUN block - * -# if t0's are missing (i.e. t0 == -1), try to fill from the GLOBAL block. - * -# if t0's are missing, try t0's from the data file - * -# if t0's are missing, try to estimate them + * \brief Determines and validates t0 values for all forward and backward histograms. * - * \param runData pointer to the current RUN block entry from the msr-file - * \param globalBlock pointer to the GLOBLA block entry from the msr-file - * \param forwardHistoNo histogram number vector of forward; forwardHistoNo = msr-file forward + redGreen_offset - 1 - * \param backwardHistoNo histogram number vector of backwardward; backwardHistoNo = msr-file backward + redGreen_offset - 1 + * Time zero (t0) marks the arrival time of muons in the sample and is critical for + * proper time alignment. This method attempts to find t0 values from multiple sources + * with a well-defined fallback hierarchy: * - * return: - * - true if everthing went smooth - * - false, otherwise. + * Priority order: + * 1. Individual RUN block t0 values (highest priority, user-specified) + * 2. GLOBAL block t0 values (shared defaults for all runs) + * 3. Data file t0 values (from detector electronics or auto-detection) + * 4. Estimated t0 values (last resort, may be unreliable for some facilities) + * + * The t0 vector is organized as [forward_0, backward_0, forward_1, backward_1, ...], + * accommodating different numbers of forward/backward histograms in grouped data. + * + * Also handles addT0 values for runs being added together, ensuring proper time alignment + * when combining multiple datasets. + * + * \param runData Pointer to raw run data containing file-based t0 information + * \param globalBlock Pointer to GLOBAL block with default t0 settings + * \param forwardHistoNo Vector of forward histogram indices (after red/green offset adjustment) + * \param backwardHistoNo Vector of backward histogram indices (after red/green offset adjustment) + * + * \return True on success, false if critical t0 values cannot be determined + * + * \warning Estimated t0 values may be unreliable for certain facilities (e.g., LEM) */ Bool_t PRunAsymmetry::GetProperT0(PRawRunData* runData, PMsrGlobalBlock *globalBlock, PUIntVector &forwardHistoNo, PUIntVector &backwardHistoNo) { @@ -1855,17 +1974,30 @@ Bool_t PRunAsymmetry::GetProperT0(PRawRunData* runData, PMsrGlobalBlock *globalB // GetProperDataRange (private) //-------------------------------------------------------------------------- /** - *

Get the proper data range, i.e. first/last good bin (fgb/lgb). - * -# get fgb/lgb from the RUN block - * -# if fgb/lgb still undefined, try to get it from the GLOBAL block - * -# if fgb/lgb still undefined, try to estimate them. + * \brief Determines the valid data range (first/last good bins) for analysis. * - * \param runData raw run data needed to perform some crosschecks - * \param histoNo histogram number (within a run). histoNo[0]: forward histogram number, histNo[1]: backward histogram number + * The data range defines which portion of the histograms contains usable data, + * excluding initial and final bins that may be noisy or affected by detector artifacts. * - * return: - * - true if everthing went smooth - * - false, otherwise. + * Determination hierarchy: + * 1. RUN block data range settings (highest priority, run-specific) + * 2. GLOBAL block data range (shared defaults) + * 3. Estimated range (last resort: start = t0 + 10ns, end = histogram length) + * + * Performs validation checks: + * - Ensures start < end (swaps if needed) + * - Verifies bins are within histogram bounds [0, histogram size] + * - Validates t0 is within valid range + * - Clips end bin if it exceeds histogram length + * + * The good bins are stored in fGoodBins as [forward_start, forward_end, backward_start, backward_end]. + * + * \param runData Pointer to raw run data for histogram size validation + * \param histoNo Array of histogram indices: [0]=forward, [1]=backward + * + * \return True on success, false if data range is invalid or out of bounds + * + * \warning Estimated data ranges may not be appropriate for all experiments */ Bool_t PRunAsymmetry::GetProperDataRange(PRawRunData* runData, UInt_t histoNo[2]) { @@ -2000,16 +2132,26 @@ Bool_t PRunAsymmetry::GetProperDataRange(PRawRunData* runData, UInt_t histoNo[2] // GetProperFitRange (private) //-------------------------------------------------------------------------- /** - *

Get the proper fit range. There are two possible fit range commands: - * fit given in (usec), or - * fit fgb+offset_0 lgb-offset_1 given in (bins), therefore it works the following way: - * -# get fit range assuming given in time from RUN block - * -# if fit range in RUN block is given in bins, replace start/end - * -# if fit range is NOT given yet, try fit range assuming given in time from GLOBAL block - * -# if fit range in GLOBAL block is given in bins, replace start/end - * -# if still no fit range is given, use fgb/lgb. + * \brief Determines the fit range for χ² minimization. * - * \param globalBlock pointer to the GLOBAL block information form the msr-file. + * The fit range defines the time window used for parameter fitting. It can be specified + * in two formats: + * - Time-based: "fit start end" in microseconds (e.g., "fit 0.1 10.0") + * - Bin-based: "fit fgb+offset lgb-offset" using good bin boundaries (e.g., "fit fgb+5 lgb-10") + * + * Determination sequence: + * 1. Checks RUN block for time-based fit range + * 2. If RUN block specifies bin-based range, converts to time using fgb/lgb offsets + * 3. Falls back to GLOBAL block fit range (time or bin-based) if RUN block is empty + * 4. Uses full data range (fgb to lgb) if no fit range is specified + * + * Bin-based format allows fine-tuning relative to good bin boundaries: + * - fgb+n: Start n bins after first good bin + * - lgb-n: End n bins before last good bin + * + * The resulting fit range is stored as time values in fFitStartTime and fFitEndTime. + * + * \param globalBlock Pointer to GLOBAL block containing default fit range settings */ void PRunAsymmetry::GetProperFitRange(PMsrGlobalBlock *globalBlock) { diff --git a/src/include/PRunAsymmetry.h b/src/include/PRunAsymmetry.h index 4e977548e..7149d314b 100644 --- a/src/include/PRunAsymmetry.h +++ b/src/include/PRunAsymmetry.h @@ -34,57 +34,250 @@ //--------------------------------------------------------------------------- /** - *

Class handling the asymmetry fit. + * \brief Class for handling standard μSR asymmetry fits. + * + * PRunAsymmetry implements asymmetry fitting for conventional μSR experiments, + * where the asymmetry is calculated from forward and backward detector histograms: + * + * \f[ A(t) = \frac{F(t) - \alpha B(t)}{F(t) + \alpha B(t)} \f] + * + * The class supports various configurations: + * - α and β parameters (detector efficiency corrections) + * - Background subtraction (fixed or estimated) + * - Bin packing for improved statistics + * - Rotating reference frame (RRF) data handling + * - Multiple histogram grouping and run adding + * + * The α/β parameter configurations are: + * - Tag 1: α = β = 1 (no corrections) + * - Tag 2: α ≠ 1, β = 1 (forward/backward asymmetry correction) + * - Tag 3: α = 1, β ≠ 1 (alternative correction) + * - Tag 4: α ≠ 1, β ≠ 1 (both corrections active) + * + * \see PRunBase for the base class providing common functionality + * \see PRunAsymmetryBNMR for β-NMR asymmetry (helicity-dependent) + * \see PRunSingleHisto for single histogram fits */ class PRunAsymmetry : public PRunBase { public: + /// Default constructor PRunAsymmetry(); + + /** + * \brief Main constructor for μSR asymmetry fitting. + * \param msrInfo Pointer to MSR file handler + * \param rawData Pointer to raw run data handler + * \param runNo Run number within the MSR file + * \param tag Operation mode (kFit for fitting, kView for viewing) + * \param theoAsData If true, calculate theory only at data points; if false, calculate additional points for Fourier + */ PRunAsymmetry(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag, Bool_t theoAsData); + + /// Destructor virtual ~PRunAsymmetry(); + /** + * \brief Calculates chi-square for the current parameter set. + * \param par Parameter vector from MINUIT + * \return Chi-square value + */ virtual Double_t CalcChiSquare(const std::vector& par); + + /** + * \brief Calculates expected chi-square (for statistical analysis). + * \param par Parameter vector from MINUIT + * \return Expected chi-square value + */ virtual Double_t CalcChiSquareExpected(const std::vector& par); + + /** + * \brief Calculates maximum likelihood estimator. + * \param par Parameter vector from MINUIT + * \return Maximum likelihood value + */ virtual Double_t CalcMaxLikelihood(const std::vector& par); + + /** + * \brief Calculates theoretical asymmetry function. + * + * Computes the theory values for the μSR asymmetry based on the + * current parameters and fit function. + */ virtual void CalcTheory(); + /** + * \brief Returns the number of bins used in the fit. + * \return Number of fit bins + */ virtual UInt_t GetNoOfFitBins(); + /** + * \brief Sets the fit range in bins. + * \param fitRange Fit range string (format depends on configuration) + */ virtual void SetFitRangeBin(const TString fitRange); + /** + * \brief Returns the first bin used in the fit. + * \return Start time bin index + */ virtual Int_t GetStartTimeBin() { return fStartTimeBin; } + + /** + * \brief Returns the last bin used in the fit. + * \return End time bin index + */ virtual Int_t GetEndTimeBin() { return fEndTimeBin; } + + /** + * \brief Returns the packing factor. + * \return Number of bins combined (1 = no packing) + */ virtual Int_t GetPacking() { return fPacking; } + /** + * \brief Calculates the number of bins to be fitted. + * + * Determines fNoOfFitBins based on the fit range and data availability. + */ virtual void CalcNoOfFitBins(); protected: + /** + * \brief Prepares all data for fitting or viewing. + * \return True on success, false on error + * + * Main data preparation routine that handles background subtraction, + * packing, and asymmetry calculation from forward/backward histograms. + */ virtual Bool_t PrepareData(); + + /** + * \brief Prepares data specifically for fitting. + * \return True on success, false on error + * + * Sets up data structures for the fitting process, including determining + * fit ranges and calculating the number of fit bins. + */ virtual Bool_t PrepareFitData(); + + /** + * \brief Prepares data for viewing/plotting. + * \param runData Pointer to raw run data + * \param histoNo Array of histogram numbers [0]=forward, [1]=backward + * \return True on success, false on error + */ virtual Bool_t PrepareViewData(PRawRunData* runData, UInt_t histoNo[2]); + + /** + * \brief Prepares rotating reference frame (RRF) data for viewing. + * \param runData Pointer to raw run data + * \param histoNo Array of histogram numbers [0]=forward, [1]=backward + * \return True on success, false on error + * + * Special handling for RRF data where asymmetry is calculated in a + * rotating coordinate frame. + */ virtual Bool_t PrepareRRFViewData(PRawRunData* runData, UInt_t histoNo[2]); private: - UInt_t fAlphaBetaTag; ///< \f$ 1 \to \alpha = \beta = 1\f$; \f$ 2 \to \alpha \neq 1, \beta = 1\f$; \f$ 3 \to \alpha = 1, \beta \neq 1\f$; \f$ 4 \to \alpha \neq 1, \beta \neq 1\f$. - UInt_t fNoOfFitBins; ///< number of bins to be be fitted - Int_t fPacking; ///< packing for this particular run. Either given in the RUN- or GLOBAL-block. - Bool_t fTheoAsData; ///< true=only calculate the theory points at the data points, false=calculate more points for the theory as compared to data are calculated which lead to 'nicer' Fouriers + UInt_t fAlphaBetaTag; ///< Tag indicating α/β configuration: 1=both unity, 2=α free/β unity, 3=α unity/β free, 4=both free + UInt_t fNoOfFitBins; ///< Number of bins included in the fit + Int_t fPacking; ///< Bin packing factor from RUN or GLOBAL block + Bool_t fTheoAsData; ///< If true, theory calculated only at data points; if false, extra points for nicer Fourier transforms - PDoubleVector fForward; ///< forward histo data - PDoubleVector fForwardErr; ///< forward histo errors - PDoubleVector fBackward; ///< backward histo data - PDoubleVector fBackwardErr; ///< backward histo errors + PDoubleVector fForward; ///< Forward detector histogram data + PDoubleVector fForwardErr; ///< Forward detector histogram errors + PDoubleVector fBackward; ///< Backward detector histogram data + PDoubleVector fBackwardErr; ///< Backward detector histogram errors - Int_t fGoodBins[4]; ///< keep first/last good bins. 0=fgb, 1=lgb (forward); 2=fgb, 3=lgb (backward) + Int_t fGoodBins[4]; ///< Good bin boundaries: [0]=forward first, [1]=forward last, [2]=backward first, [3]=backward last - Int_t fStartTimeBin; ///< bin at which the fit starts - Int_t fEndTimeBin; ///< bin at which the fit ends + Int_t fStartTimeBin; ///< First bin index for fitting + Int_t fEndTimeBin; ///< Last bin index for fitting + /** + * \brief Subtracts fixed background from histograms. + * + * Subtracts user-specified fixed background values from forward and backward histograms. + * Background values are read from the MSR file (e.g., "backgr.fix 2 3" for forward/backward). + * + * Error propagation: + * \f[ \Delta f_i^{\rm c} = \pm\sqrt{(\Delta f_i)^2 + (\Delta \mathrm{bkg})^2} = \pm\sqrt{f_i + \mathrm{bkg}} \f] + * + * where \f$ f_i^{\rm c} \f$ is the background-corrected histogram, \f$ f_i \f$ is the raw histogram, + * and \f$ \mathrm{bkg} \f$ is the fixed background value. + * + * \return True on success, false if background values are missing + */ Bool_t SubtractFixBkg(); + + /** + * \brief Estimates and subtracts background from histograms. + * + * Calculates background from a specified bin range (typically before t0) and subtracts it + * from forward and backward histograms. The background range is adjusted to align with + * accelerator beam cycles when applicable (PSI, RAL, TRIUMF). + * + * Background calculation: + * \f[ \mathrm{bkg} = \frac{1}{N}\sum_{i=0}^N f_i \f] + * + * Error propagation: + * \f[ \Delta f_i^{\rm c} = \pm\sqrt{(\Delta f_i)^2 + (\Delta \mathrm{bkg})^2} = \pm\sqrt{f_i + (\Delta \mathrm{bkg})^2} \f] + * + * where \f$ N \f$ is the number of background bins and + * \f[ \Delta \mathrm{bkg} = \pm\frac{1}{N}\sqrt{\sum_{i=0}^N f_i} \f] + * + * \return True on success, false if background range is out of bounds + */ Bool_t SubtractEstimatedBkg(); + /** + * \brief Retrieves proper t0 values for all histograms. + * + * Determines t0 (time zero) values for forward and backward histograms from: + * 1. Individual RUN block t0 entries + * 2. GLOBAL block t0 values (fallback) + * 3. Data file t0 values (final fallback) + * + * Also handles addT0 corrections for runs to be added together. The t0 values + * are critical for proper alignment of histograms in time. + * + * \param runData Pointer to raw run data containing histogram information + * \param globalBlock Pointer to global MSR block with default t0 values + * \param forwardHisto Vector of forward histogram indices + * \param backwardHistoNo Vector of backward histogram indices + * \return True on success, false if t0 values cannot be determined + */ virtual Bool_t GetProperT0(PRawRunData* runData, PMsrGlobalBlock *globalBlock, PUIntVector &forwardHisto, PUIntVector &backwardHistoNo); + + /** + * \brief Retrieves proper data range for histograms. + * + * Determines the "good bins" range for data analysis from: + * 1. Individual RUN block data range settings + * 2. GLOBAL block data range (fallback) + * 3. Full histogram range (final fallback) + * + * The good bins define which portion of the histograms contains valid data, + * excluding noisy or problematic regions at the start/end. + * + * \param runData Pointer to raw run data + * \param histoNo Array of histogram numbers [0]=forward, [1]=backward + * \return True on success, false on error + */ virtual Bool_t GetProperDataRange(PRawRunData* runData, UInt_t histoNo[2]); + + /** + * \brief Determines the proper fit range from global block. + * + * Extracts fit range settings from the GLOBAL block if not specified in the RUN block. + * The fit range defines the time window used for χ² minimization, specified either in + * time units (μs) or bin numbers. + * + * \param globalBlock Pointer to global MSR block containing default fit range + */ virtual void GetProperFitRange(PMsrGlobalBlock *globalBlock); };