diff --git a/src/classes/PRunBase.cpp b/src/classes/PRunBase.cpp index 1568ce638..5db1e2fe7 100644 --- a/src/classes/PRunBase.cpp +++ b/src/classes/PRunBase.cpp @@ -42,7 +42,14 @@ // Constructor //-------------------------------------------------------------------------- /** - *

Constructor. Needed otherwise vector's cannot be generated ;-) + * \brief Default constructor that initializes all member variables to default values. + * + * Creates an empty, invalid run object with all pointers set to nullptr and + * values set to undefined/invalid states. This constructor is needed to allow + * creation of vectors of PRunBase-derived objects. + * + * A run created with this constructor requires initialization via the main + * constructor before it can be used for fitting. */ PRunBase::PRunBase() { @@ -64,12 +71,25 @@ PRunBase::PRunBase() // Constructor //-------------------------------------------------------------------------- /** - *

Constructor. + * \brief Main constructor that initializes a run from MSR file and raw data. * - * \param msrInfo pointer to the msr-file handler - * \param rawData pointer to the raw-data handler - * \param runNo msr-file run number - * \param tag tag telling if fit, view, or rrf representation is whished. + * Performs comprehensive initialization: + * 1. Stores operation mode (fit vs. view) and pointers to MSR and data handlers + * 2. Extracts run-specific settings from the appropriate MSR RUN block + * 3. Validates function parameter mappings to ensure FUNCTIONS block is valid + * 4. Initializes metadata structures (field, energy, temperature) + * 5. Initializes function value cache for FUNCTIONS block evaluation + * 6. Creates PTheory object for evaluating the theory function + * 7. Validates that theory initialization succeeded + * + * If any initialization step fails (e.g., invalid theory, out-of-range parameters), + * the program exits with an error message. The run object is marked as valid upon + * successful completion. + * + * \param msrInfo Pointer to MSR file handler (must remain valid for object lifetime) + * \param rawData Pointer to raw data handler (must remain valid for object lifetime) + * \param runNo Run number (0-based index into MSR file RUN blocks) + * \param tag Operation mode: kFit (fitting), kView (display/plotting) */ PRunBase::PRunBase(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag) : fHandleTag(tag), fMsrInfo(msrInfo), fRawData(rawData) @@ -118,7 +138,16 @@ PRunBase::PRunBase(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, // Destructor //-------------------------------------------------------------------------- /** - *

Destructor. + * \brief Virtual destructor that cleans up allocated resources. + * + * Frees memory allocated for: + * - t0 value vectors (fT0s) + * - Additional run t0 vectors (fAddT0s) + * - Function value cache (fFuncValues) + * + * The PTheory object is automatically deleted via unique_ptr. + * Pointers to MSR handler and raw data handler are NOT deleted as they + * are owned by the calling code. */ PRunBase::~PRunBase() { @@ -135,9 +164,25 @@ PRunBase::~PRunBase() // SetFitRange (public) //-------------------------------------------------------------------------- /** - *

Sets the current fit range, and recalculated the number of fitted bins + * \brief Sets the fit time range and recalculates the number of fitted bins. * - * \param fitRange vector with fit ranges + * Updates the fitting window for this run. This method handles two scenarios: + * 1. Single fit range: If fitRange contains one pair, it applies to all runs + * 2. Multiple fit ranges: If fitRange contains multiple pairs, it selects the + * appropriate range based on fRunNo + * + * The method validates that: + * - The fit range vector is not empty (asserts) + * - The run number is within the fit range vector bounds + * - Start time is before end time (swaps if not) + * + * This is typically called by the FIT_RANGE command to dynamically adjust the + * fitting window during optimization or range scanning. + * + * \param fitRange Vector of (start, end) time pairs in microseconds (μs) + * + * Example: fitRange = {(0.1, 10.0), (0.2, 8.0)} applies first range to run 0, + * second range to run 1 */ void PRunBase::SetFitRange(PDoublePairVector fitRange) { @@ -177,7 +222,13 @@ void PRunBase::SetFitRange(PDoublePairVector fitRange) // CleanUp (public) //-------------------------------------------------------------------------- /** - *

Clean up all locally allocate memory + * \brief Cleans up allocated resources. + * + * Releases memory used by the PTheory object via unique_ptr reset. + * This is called when the run processing is complete or when preparing + * for a new fit with different settings. + * + * Other data structures (vectors) are managed automatically by their destructors. */ void PRunBase::CleanUp() { @@ -188,14 +239,37 @@ void PRunBase::CleanUp() // CalculateKaiserFilterCoeff (protected) //-------------------------------------------------------------------------- /** - *

Calculates the Kaiser filter coefficients for a low pass filter with - * a cut off frequency wc. - * For details see "Zeitdiskrete Signalverarbeitung", A.V. Oppenheim, R.W. Schafer, J.R. Buck. Pearson 2004. + * \brief Calculates Kaiser window FIR filter coefficients for RRF smoothing. * - * \param wc cut off frequency - * \param A defined as \f$ A = -\log_{10}(\delta) \f$, where \f$\delta\f$ is the tolerance band. - * \param dw defined as \f$ \Delta\omega = \omega_{\rm S} - \omega_{\rm P} \f$, where \f$ \omega_{\rm S} \f$ is the - * stop band frequency, and \f$ \omega_{\rm P} \f$ is the pass band frequency. + * Designs a low-pass FIR filter using the Kaiser window method, which provides + * excellent control over the frequency response characteristics. The filter is + * used to smooth theory curves in rotating reference frame (RRF) fits, ensuring + * consistent comparison between filtered data and filtered theory. + * + * Algorithm (based on Oppenheim, Schafer, Buck, "Discrete-Time Signal Processing"): + * 1. Determine β parameter from attenuation requirement A + * 2. Calculate filter order m from transition width dw + * 3. Generate Kaiser window: \f$ w[n] = \frac{I_0(\beta\sqrt{1-(n/\alpha)^2})}{I_0(\beta)} \f$ + * 4. Apply window to ideal sinc filter: \f$ h[n] = \frac{\sin(\omega_c n)}{\pi n} \cdot w[n] \f$ + * 5. Normalize coefficients to unity gain at DC + * + * The β parameter is chosen based on attenuation A (in dB): + * - A > 50: β = 0.1102(A - 8.7) + * - 21 ≤ A ≤ 50: β = 0.5842(A - 21)^0.4 + 0.07886(A - 21) + * - A < 21: β = 0 + * + * Filter order: m = (A - 8) / (2.285 × Δω × π), rounded to nearest odd integer + * + * Reference: A.V. Oppenheim, R.W. Schafer, J.R. Buck, + * "Discrete-Time Signal Processing", Pearson 2004, pp. 574ff + * + * \param wc Cutoff frequency (normalized, 0 to π rad/sample) + * \param A Attenuation in dB: A = -20 log₁₀(δ) where δ is the tolerance band + * \param dw Transition width: Δω = ω_S - ω_P (stop band - pass band frequencies) + * + * The computed coefficients are stored in fKaiserFilter and normalized for unity DC gain. + * + * \see FilterTheo() for application of these coefficients */ void PRunBase::CalculateKaiserFilterCoeff(Double_t wc, Double_t A, Double_t dw) { @@ -235,7 +309,33 @@ void PRunBase::CalculateKaiserFilterCoeff(Double_t wc, Double_t A, Double_t dw) // FilterTheo (protected) //-------------------------------------------------------------------------- /** - *

Filters the theory with a Kaiser FIR filter. + * \brief Applies Kaiser FIR filter to theory values for RRF fits. + * + * Performs time-domain convolution of the theory function with the Kaiser filter + * coefficients computed by CalculateKaiserFilterCoeff(). This smooths the theory + * curve to match the filtering applied to RRF-transformed data, ensuring fair + * comparison during χ² calculation. + * + * The filtering operation is: + * \f[ y_{\rm filtered}[i] = \sum_{j=0}^{M-1} h[j] \cdot y_{\rm theory}[i-j] \f] + * + * where: + * - h[j] are the Kaiser filter coefficients from fKaiserFilter + * - M is the filter length + * - For i < j, the missing samples are treated as zero (causal filter) + * + * Additional processing: + * - The filtered theory replaces the original theory in fData + * - Time start is shifted backward by half the filter length to compensate + * for the group delay introduced by the symmetric FIR filter + * + * This method is only called by RRF-derived classes (PRunAsymmetryRRF, PRunSingleHistoRRF) + * after theory calculation and RRF transformation. + * + * \pre CalculateKaiserFilterCoeff() must have been called to initialize fKaiserFilter + * \pre fData must contain theory values (CalcTheory() must have been called) + * + * \see CalculateKaiserFilterCoeff() for filter coefficient calculation */ void PRunBase::FilterTheo() { diff --git a/src/include/PRunBase.h b/src/include/PRunBase.h index aa3918735..8ef723658 100644 --- a/src/include/PRunBase.h +++ b/src/include/PRunBase.h @@ -42,33 +42,62 @@ //------------------------------------------------------------------------------------------ /** - *

Abstract base class defining the interface for all μSR fit types. + * \brief Abstract base class defining the interface for all μSR fit types. * - *

PRunBase establishes a common API for processing and fitting different - * types of μSR data (single histogram, asymmetry, RRF, etc.). Derived classes - * implement specific data processing and χ² calculation for each fit type: - * - PRunSingleHisto: Single detector histogram fits - * - PRunAsymmetry: Asymmetry fits (forward - backward) - * - PRunSingleHistoRRF: Single histogram in rotating reference frame + * PRunBase establishes a common API for processing and fitting different types of μSR data + * (single histogram, asymmetry, RRF, etc.). This class serves as the foundation of the + * musrfit framework, providing core functionality and requiring derived classes to implement + * fit-type-specific algorithms. + * + * \section derived_classes Derived Classes + * Each derived class handles a specific type of μSR measurement: + * - PRunSingleHisto: Single detector histogram fits (basic time-differential μSR) + * - PRunAsymmetry: Asymmetry fits: \f$ A(t) = \frac{F(t) - \alpha B(t)}{F(t) + \alpha B(t)} \f$ + * - PRunSingleHistoRRF: Single histogram in rotating reference frame (high-TF analysis) * - PRunAsymmetryRRF: Asymmetry in rotating reference frame - * - PRunMuMinus: Negative muon fits - * - PRunNonMusr: General x-y data fits + * - PRunAsymmetryBNMR: β-NMR asymmetry (helicity-dependent measurements) + * - PRunMuMinus: Negative muon fits (different decay properties) + * - PRunNonMusr: General x-y data fits (non-μSR time series data) * - *

Key responsibilities: - * - Loading and preprocessing raw histogram data - * - Background subtraction and normalization - * - Time bin packing/rebinning - * - Theory evaluation at data points - * - χ² or maximum likelihood calculation - * - RRF transformations (if applicable) + * \section key_responsibilities Key Responsibilities + * - Loading raw histogram data from various file formats (ROOT, NeXus, WKM, etc.) + * - Extracting metadata (magnetic field, temperature, energy) from data files + * - Managing time-zero (t0) determination and validation + * - Background subtraction (fixed or estimated from pre-t0 region) + * - Time bin packing/rebinning for improved statistics + * - Theory function evaluation via PTheory interface + * - χ² or maximum likelihood calculation for fitting + * - RRF transformations with Kaiser filtering (for RRF-derived classes) + * - Run addition (combining multiple identical measurements) + * - Histogram grouping (combining multiple detectors) * - *

Processing workflow: - * 1. PrepareData() - Load and preprocess raw data - * 2. CalcTheory() - Evaluate theory function - * 3. CalcChiSquare() or CalcMaxLikelihood() - Compute fit metric + * \section workflow Processing Workflow + * The typical processing sequence for a fit is: + * 1. Construction: Initialize from MSR file and raw data handler + * 2. PrepareData(): Load and preprocess raw data (implemented by derived classes) + * - Load histograms, validate t0 values, subtract background + * - Calculate asymmetry or apply RRF transformation (if applicable) + * - Pack bins, propagate errors + * 3. Fitting loop (called by MINUIT): + * - CalcTheory(): Evaluate theory function at data points + * - CalcChiSquare(): Compute χ² between data and theory + * - MINUIT adjusts parameters to minimize χ² + * 4. Visualization: Theory calculated with higher resolution for plotting * - *

Design pattern: Template Method - base class defines workflow, - * derived classes implement fit-type-specific algorithms. + * \section design_pattern Design Pattern + * PRunBase implements the Template Method design pattern: + * - Base class defines the overall workflow and common operations + * - Derived classes implement fit-type-specific algorithms (pure virtual methods) + * - Ensures consistency across different fit types while allowing specialization + * + * \section thread_safety Thread Safety + * PRunBase objects are NOT thread-safe. Each thread in parallel fitting should + * create its own PRunBase-derived object. Theory evaluation may use OpenMP + * internally (see PTheory documentation). + * + * \see PTheory for theory function evaluation + * \see PMsrHandler for MSR file parsing + * \see PRunDataHandler for raw data loading */ class PRunBase { @@ -77,130 +106,252 @@ class PRunBase PRunBase(); /** - *

Constructor initializing run from MSR file and raw data. + * \brief Constructor initializing run from MSR file and raw data. * - * @param msrInfo Pointer to MSR file handler - * @param rawData Pointer to raw data handler - * @param runNo Run number (0-based index in MSR file) - * @param tag Operation mode (kFit, kView) + * Initializes the run object by: + * - Storing pointers to MSR file handler and raw data handler + * - Extracting run-specific settings from MSR RUN block + * - Validating function parameter mappings + * - Creating PTheory object for theory evaluation + * - Initializing metadata and function value vectors + * + * \param msrInfo Pointer to MSR file handler containing all fit settings + * \param rawData Pointer to raw data handler for accessing histogram data + * \param runNo Run number (0-based index in MSR file RUN blocks) + * \param tag Operation mode: kFit (fitting), kView (display only) */ PRunBase(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag); + /** + * \brief Virtual destructor. + * + * Cleans up allocated resources including: + * - t0 value vectors + * - Function value vectors + * - PTheory object (via unique_ptr) + */ virtual ~PRunBase(); /** - *

Calculates χ² between data and theory (pure virtual). + * \brief Calculates χ² between data and theory (pure virtual). * - *

χ² = Σ[(data_i - theory_i)² / σ_i²] summed over all data points. - * This is the standard least-squares metric for fitting. + * Computes the chi-squared statistic: + * \f[ \chi^2 = \sum_{i=1}^{N} \frac{(y_i^{\rm data} - y_i^{\rm theory})^2}{\sigma_i^2} \f] * - * @param par Vector of fit parameter values - * @return Chi-squared value + * This is the standard least-squares metric for fitting. It assumes Gaussian + * statistics (valid for high-count data). For low-count data, consider using + * CalcMaxLikelihood() instead. + * + * Called by MINUIT during each fit iteration to evaluate parameter quality. + * Derived classes implement the specific calculation for their data type. + * + * \param par Vector of fit parameter values from MINUIT + * \return Chi-squared value (lower is better) + * + * \see CalcMaxLikelihood for alternative fit metric */ virtual Double_t CalcChiSquare(const std::vector& par) = 0; /** - *

Calculates maximum likelihood (pure virtual). + * \brief Calculates expected chi-square for statistical analysis (pure virtual). * - *

For Poisson statistics: -2ln(L) = 2Σ[theory_i - data_i·ln(theory_i)] - * Better than χ² for low-count data where Gaussian approximation fails. + * Computes the expected χ² value based on the theory and error estimates, + * useful for evaluating the quality of error bars and fit statistics. * - * @param par Vector of fit parameter values - * @return Negative 2 times log-likelihood + * For properly estimated errors: χ²_expected ≈ number of degrees of freedom + * + * \param par Vector of fit parameter values from MINUIT + * \return Expected chi-squared value + * + * \note Implementation is optional in derived classes; many return 0.0 + */ + virtual Double_t CalcChiSquareExpected(const std::vector& par) = 0; + + /** + * \brief Calculates maximum likelihood estimator (pure virtual). + * + * Computes the negative log-likelihood for Poisson statistics: + * \f[ -2\ln L = 2\sum_{i=1}^{N} \left[y_i^{\rm theory} - y_i^{\rm data} \ln(y_i^{\rm theory})\right] \f] + * + * Maximum likelihood estimation is superior to χ² for low-count data where + * the Gaussian approximation breaks down (typically when counts < 10-20 per bin). + * It naturally handles Poisson statistics without requiring error estimates. + * + * Called by MINUIT as an alternative to χ² minimization. MINUIT minimizes + * this function just like χ². + * + * \param par Vector of fit parameter values from MINUIT + * \return Negative 2 times log-likelihood (lower is better) + * + * \see CalcChiSquare for standard least-squares metric + * \note Not implemented in all derived classes */ virtual Double_t CalcMaxLikelihood(const std::vector& par) = 0; /** - *

Sets the fit time range for this run. + * \brief Sets the fit time range for this run. * - *

Updates the fitting window, useful for FIT_RANGE command which - * scans different time windows to find the optimal range. + * Updates the fitting window, useful for the FIT_RANGE command which allows + * scanning different time windows to find the optimal range for parameter + * extraction. Can be called multiple times during a fit sequence. * - * @param fitRange Vector of (start, end) time pairs in microseconds + * The fit range is specified in microseconds (μs) from t0. Multiple ranges + * can be specified for different runs in a global fit. + * + * \param fitRange Vector of (start, end) time pairs in microseconds + * + * Example: fitRange[0] = (0.1, 10.0) means fit from 0.1 μs to 10.0 μs after t0 */ virtual void SetFitRange(PDoublePairVector fitRange); /** - *

Evaluates theory function at all data points (pure virtual). + * \brief Evaluates theory function at all data points (pure virtual). * - *

Uses current parameter values to calculate expected signal - * at each time point. Called during each fit iteration. + * Calculates the expected signal at each time point using the current parameter + * values from the MSR THEORY block. This is called: + * - During each MINUIT fit iteration + * - After fitting for visualization + * - When evaluating functions from FUNCTIONS block + * + * The theory values are stored in fData for comparison with measured data. + * Derived classes implement fit-type-specific theory calculation (e.g., + * single histogram vs. asymmetry). + * + * \see PTheory for the underlying theory evaluation engine */ virtual void CalcTheory() = 0; - /// Returns the run number (0-based index in MSR file) - /// @return Run number + /** + * \brief Returns the run number (0-based index in MSR file). + * \return Run number corresponding to position in MSR RUN blocks + */ virtual UInt_t GetRunNo() { return fRunNo; } - /// Returns pointer to processed data (background-corrected, binned) - /// @return Pointer to PRunData object + /** + * \brief Returns pointer to processed data container. + * + * The PRunData object contains: + * - Background-corrected histogram data + * - Packed/rebinned data points + * - Error bars + * - Time grid information + * - Theory values (after CalcTheory() is called) + * + * \return Pointer to PRunData object with processed data + */ virtual PRunData* GetData() { return &fData; } - /// Cleans up internal data structures + /** + * \brief Cleans up internal data structures. + * + * Releases memory used by temporary data structures. Called when + * run processing is complete or when resetting for a new fit. + */ virtual void CleanUp(); - /// Returns validity status of this run object - /// @return true if run initialized successfully + /** + * \brief Returns validity status of this run object. + * + * A run becomes invalid if: + * - Required data files cannot be loaded + * - MSR file settings are inconsistent + * - Theory initialization fails + * - Data preprocessing encounters errors + * + * \return True if run initialized successfully, false otherwise + */ virtual Bool_t IsValid() { return fValid; } protected: - Bool_t fValid; ///< Flag showing if the run object initialized successfully + Bool_t fValid; ///< Flag indicating if run object initialized successfully; false if any error occurred - EPMusrHandleTag fHandleTag; ///< Operation mode: kFit (fitting), kView (display only) + EPMusrHandleTag fHandleTag; ///< Operation mode: kFit (fitting), kView (display only), kEmpty (uninitialized) Int_t fRunNo; ///< Run number (0-based index in MSR file RUN blocks) - PMsrHandler *fMsrInfo; ///< Pointer to MSR file handler - PMsrRunBlock *fRunInfo; ///< Pointer to this run's RUN block settings - PRunDataHandler *fRawData; ///< Pointer to raw data handler for this run + PMsrHandler *fMsrInfo; ///< Pointer to MSR file handler (owned externally, not deleted here) + PMsrRunBlock *fRunInfo; ///< Pointer to this run's RUN block settings within fMsrInfo + PRunDataHandler *fRawData; ///< Pointer to raw data handler (owned externally, not deleted here) - PRunData fData; ///< Processed data ready for fitting (background-corrected, packed, etc.) - Double_t fTimeResolution; ///< Time resolution of raw data in microseconds (μs) - PMetaData fMetaData; ///< Experimental metadata (field, temperature, energy) from data file - PDoubleVector fT0s; ///< Time-zero bins for all histograms in this run - std::vector fAddT0s; ///< Time-zero bins for all addrun histograms + PRunData fData; ///< Processed data container: background-corrected, packed, with theory values + Double_t fTimeResolution; ///< Time resolution of raw histogram data in microseconds (μs), e.g., 0.01953125 μs for PSI GPS + PMetaData fMetaData; ///< Experimental metadata extracted from data file header (magnetic field, temperature, beam energy) + PDoubleVector fT0s; ///< Time-zero bin values for all histograms in this run (forward, backward, etc.) + std::vector fAddT0s; ///< Time-zero bin values for additional runs to be added to main run - Double_t fFitStartTime; ///< Fit range start time in microseconds (μs) - Double_t fFitEndTime; ///< Fit range end time in microseconds (μs) + Double_t fFitStartTime; ///< Fit range start time in microseconds (μs) relative to t0 + Double_t fFitEndTime; ///< Fit range end time in microseconds (μs) relative to t0 - PDoubleVector fFuncValues; ///< Values of user-defined functions from FUNCTIONS block - std::unique_ptr fTheory; ///< Theory function evaluator for χ² calculation + PDoubleVector fFuncValues; ///< Cached values of user-defined functions from FUNCTIONS block, evaluated at current parameters + std::unique_ptr fTheory; ///< Theory function evaluator (smart pointer, automatically deleted) - PDoubleVector fKaiserFilter; ///< Kaiser window coefficients for RRF filtering + PDoubleVector fKaiserFilter; ///< Kaiser window FIR filter coefficients for smoothing RRF theory curves /** - *

Prepares raw data for fitting (pure virtual). + * \brief Prepares raw data for fitting (pure virtual). * - *

Performs fit-type-specific preprocessing: - * - Background subtraction - * - Asymmetry calculation (if applicable) - * - Time bin packing/rebinning - * - RRF transformation (if applicable) - * - Error propagation + * This is the main data preprocessing pipeline, implemented differently by each + * derived class according to the fit type. Common operations include: + * - Loading histogram data from raw data files + * - Validating and determining t0 values + * - Subtracting background (fixed or estimated from pre-t0 region) + * - Calculating asymmetry (for asymmetry-based fits) + * - Applying RRF transformation (for RRF fits) + * - Time bin packing/rebinning to improve statistics + * - Proper error propagation through all transformations + * - Adding multiple runs together (if specified) + * - Grouping detector histograms (if specified) * - *

Called during initialization before fitting begins. + * Called once during object construction. If this returns false, the run + * object is marked as invalid. * - * @return true on success, false on error + * \return True on success, false if any preprocessing step fails + * + * \see PRunAsymmetry::PrepareData() for asymmetry-specific implementation + * \see PRunSingleHisto::PrepareData() for single histogram implementation */ virtual Bool_t PrepareData() = 0; /** - *

Calculates Kaiser window filter coefficients for RRF. + * \brief Calculates Kaiser window FIR filter coefficients for RRF smoothing. * - *

The Kaiser window reduces spectral leakage when transforming - * to the rotating reference frame. Parameters control the trade-off - * between main lobe width and side lobe suppression. + * Computes a Kaiser window finite impulse response (FIR) filter for smoothing + * theory curves in rotating reference frame (RRF) fits. The Kaiser window + * provides excellent control over the trade-off between main lobe width + * (frequency resolution) and side lobe attenuation (spectral leakage). * - * @param wc Cutoff frequency (normalized, 0 to π) - * @param A Attenuation in dB (controls side lobe level) - * @param dw Transition width (normalized) + * The filter design uses the Kaiser-Bessel formula: + * \f[ w[n] = \frac{I_0\left(\beta\sqrt{1-\left(\frac{n-\alpha}{\alpha}\right)^2}\right)}{I_0(\beta)} \f] + * + * where \f$ I_0 \f$ is the modified Bessel function of the first kind, + * \f$ \alpha = (M-1)/2 \f$, and \f$ \beta \f$ is determined from the attenuation A. + * + * \param wc Cutoff frequency (normalized, 0 to π rad/sample) + * \param A Attenuation in dB (typical: 60 dB for good side lobe suppression) + * \param dw Transition width (normalized, typical: 0.1-0.3) + * + * Coefficients are stored in fKaiserFilter for use by FilterTheo(). + * + * \see FilterTheo() for application of the filter */ virtual void CalculateKaiserFilterCoeff(Double_t wc, Double_t A, Double_t dw); /** - *

Applies Kaiser filter to theory values for RRF. + * \brief Applies Kaiser FIR filter to theory values for RRF fits. * - *

Filters the theory function in the same way data is filtered, - * ensuring consistent comparison between data and theory in RRF fits. + * Filters the theory function stored in fData using the Kaiser window + * coefficients from fKaiserFilter. This ensures the theory curve is smoothed + * in the same way as the RRF-transformed data, enabling fair comparison + * in χ² calculation. + * + * The filtering is performed via convolution in the time domain: + * \f[ y_{\rm filtered}[n] = \sum_{k=0}^{M-1} h[k] \cdot y_{\rm theory}[n-k] \f] + * + * where h[k] are the Kaiser filter coefficients and M is the filter length. + * + * Only used in RRF-derived classes (PRunAsymmetryRRF, PRunSingleHistoRRF). + * + * \pre CalculateKaiserFilterCoeff() must be called first to initialize fKaiserFilter + * \pre fData must contain theory values (CalcTheory() must have been called) */ virtual void FilterTheo(); };