npstat is hosted by Hepforge, IPPP Durham
NPStat  5.10.0
emsunfold::SparseUnfoldingBandwidthScannerND< Matrix > Class Template Reference

Public Types

typedef AbsSparseUnfoldND< Matrix > unfolding_type
 
typedef unfolding_type::response_matrix_type response_matrix_type
 
typedef unfolding_type::input_covariance_type input_covariance_type
 
typedef unfolding_type::output_covariance_type output_covariance_type
 
typedef unfolding_type::filter_type filter_type
 

Public Member Functions

 SparseUnfoldingBandwidthScannerND (unfolding_type &unfold, const std::vector< npstat::SymbetaParams1D > &filterParameters, const EigenParameters &trlanParameters, double minAbsoluteCorrelation, const npstat::ArrayND< double > &observed, double nDoFCorrectionFactor, const input_covariance_type *observationCovariance=0, const npstat::ArrayND< double > *oracle=0)
 
unsigned unfoldedDim () const
 
unsigned observedDim () const
 
virtual std::vector< std::string > variableNames () const
 
virtual unsigned variableCount () const
 
virtual unsigned ntuplize (double *buf, unsigned len) const
 
bool process (const std::vector< double > &bandwidthValues)
 
const std::vector< double > & lastBandwidthValues () const
 
virtual void setObservedData (const npstat::ArrayND< double > &observed, const input_covariance_type *observationCovarianceMatrix)
 
const npstat::ArrayND< double > & getObservedData () const
 
virtual void setBias (const npstat::ArrayND< double > &bias)
 
const npstat::ArrayND< double > & getBias () const
 
virtual void clearBias ()
 
virtual bool lastSparseUnfoldingStatus () const
 
const npstat::ArrayND< double > & unfoldedResult () const
 
const output_covariance_type & unfoldedCovariance () const
 
const response_matrix_type & responseMatrix () const
 
const npstat::ArrayND< double > & smoothedOracle () const
 
const std::vector< double > & eigenDeltas () const
 
const std::vector< double > & covEigenValues () const
 
bool usingConvolutions () const
 
void setInitialApproximation (const npstat::ArrayND< double > &a)
 
void clearInitialApproximation ()
 
const npstat::ArrayND< double > & getInitialApproximation () const
 
double nDoFCorrectionFactor () const
 
void setNDoFCorrectionFactor (const double f)
 

Protected Member Functions

virtual bool performUnfolding ()
 

Protected Attributes

unfolding_typeunfold_
 
std::vector< npstat::SymbetaParams1DfilterParameters_
 
EigenParameters trlanParameters_
 
double minAbsoluteCorrelation_
 
npstat::ArrayND< double > observed_
 
npstat::ArrayND< double > oracle_
 
input_covariance_type observationCovariance_
 
bool haveObservationCovariance_
 
double binVolume_
 
double nDoFCorr_
 
double nObserved_
 
double obsNonZeroFraction_
 
std::vector< double > bandwidthValues_
 
npstat::ArrayND< double > unfolded_
 
output_covariance_type unfoldedCovmat_
 
npstat::ArrayND< double > smoothedOracle_
 
std::vector< double > eigenDeltas_
 
std::vector< double > covEigenValues_
 
TrlanDiagnostics diagnose_
 

Constructor & Destructor Documentation

◆ SparseUnfoldingBandwidthScannerND()

template<class Matrix >
emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::SparseUnfoldingBandwidthScannerND ( unfolding_type unfold,
const std::vector< npstat::SymbetaParams1D > &  filterParameters,
const EigenParameters trlanParameters,
double  minAbsoluteCorrelation,
const npstat::ArrayND< double > &  observed,
double  nDoFCorrectionFactor,
const input_covariance_type *  observationCovariance = 0,
const npstat::ArrayND< double > *  oracle = 0 
)

The constructor arguments are as follows:

unfold – An instance of AbsSparseUnfoldND class.

filterParameters – Specifications how to build the smoothing filters for each dimension of the unfolded data. The corresponding values will be eventually passed to the "symbetaLOrPEFilter1D" call. The number of elements in this vector must be equal to the rank of the response matrix used by the "unfold" argument.

trlanParameters – Parameters which steer determination of covariance matrix eigenspectrum by TRLAN.

minAbsoluteCorrelation – Correlation coefficients of the unfolded covariance matrix with absolute values below this cutoff will be set to 0 in an attempt to reduce the matrix size. If this argument is 0 or negative, the unfolded covariance matrix will not be pruned.

observed – Observed data. Can be modified later using "setObservedData" method. Array shape, however, can not change.

nDoFCorrectionFactor – The correction factor to use for determining the number of effective parameters in the fit. If this argument is positive, it will be multiplied by the number of parameters determined by the standard procedure (so set it to 1.0 in order not to apply any correction). If this argument is 0 or negative, the correction factor will be calculated as the fraction of bins filled in the observed data.

observationCovariance – Covariance matrix for the observed values. If not NULL, this matrix should be symmetric and positive-definite. An internal copy will be made.

oracle – The "correct" unfolded distribution for use in various comparisons with unfolded results. If oracle is NULL, this distribution will be considered unknown, and the corresponding comparisons will not be made.

Member Function Documentation

◆ clearBias()

template<class Matrix >
virtual void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::clearBias ( )
inlinevirtual

Clear bias data

◆ clearInitialApproximation()

template<class Matrix >
void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::clearInitialApproximation ( )
inline

Clear the initial approximation to the unfolded solution

◆ covEigenValues()

template<class Matrix >
const std::vector<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::covEigenValues ( ) const
inline

Return covariance matrix eigenvalues (in the decreasing order)

◆ eigenDeltas()

template<class Matrix >
const std::vector<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::eigenDeltas ( ) const
inline

Return eigenvector differences divided by sigma

◆ getBias()

template<class Matrix >
const npstat::ArrayND<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::getBias ( ) const
inline

Return the bias data provided by the last "setBias" call

◆ getInitialApproximation()

template<class Matrix >
const npstat::ArrayND<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::getInitialApproximation ( ) const
inline

Return the initial approximation to the unfolded solution

◆ getObservedData()

template<class Matrix >
const npstat::ArrayND<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::getObservedData ( ) const
inline

Return the observed data

◆ lastBandwidthValues()

template<class Matrix >
const std::vector<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::lastBandwidthValues ( ) const
inline

Last bandwidth values processed

◆ lastSparseUnfoldingStatus()

template<class Matrix >
virtual bool emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::lastSparseUnfoldingStatus ( ) const
inlinevirtual

Status of the last unfolding call ("true" means success)

◆ nDoFCorrectionFactor()

template<class Matrix >
double emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::nDoFCorrectionFactor ( ) const
inline

Correction factor for the number of degrees of freedom

◆ ntuplize()

template<class Matrix >
virtual unsigned emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::ntuplize ( double *  buf,
unsigned  len 
) const
virtual

Write out produced variables into a common buffer. The order of the values will be consistent with the names returned by the "variableNames" method. The function returns the number of variables filled.

◆ observedDim()

template<class Matrix >
unsigned emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::observedDim ( ) const
inline

Dimensionality of the observed space

◆ process()

template<class Matrix >
bool emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::process ( const std::vector< double > &  bandwidthValues)

Perform unfolding with the given bandwidth values. This method returns "true" on success and "false" on failure.

◆ responseMatrix()

template<class Matrix >
const response_matrix_type& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::responseMatrix ( ) const
inline

Response matrix

◆ setBias()

template<class Matrix >
virtual void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::setBias ( const npstat::ArrayND< double > &  bias)
virtual

Set the bias (this is useful for various studies of uncertainties). If set, this bias will be added to the smoothed oracle data. The dimensions must be compatible with those of the unfolded result.

◆ setInitialApproximation()

template<class Matrix >
void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::setInitialApproximation ( const npstat::ArrayND< double > &  a)
inline

Set the initial approximation to the unfolded solution

◆ setNDoFCorrectionFactor()

template<class Matrix >
void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::setNDoFCorrectionFactor ( const double  f)
inline

Set the correction factor for the number of degrees of freedom

◆ setObservedData()

template<class Matrix >
virtual void emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::setObservedData ( const npstat::ArrayND< double > &  observed,
const input_covariance_type *  observationCovarianceMatrix 
)
virtual

Change the array of observations. Array dimensions must be compatible with those given in the constructor. The pointer to the covariance matrix of observations can be NULL. If provided, this matrix should be symmetric and positive-definite.

◆ smoothedOracle()

template<class Matrix >
const npstat::ArrayND<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::smoothedOracle ( ) const
inline

Return oracle data smoothed with the last processed bandwidth. If the bias was set, it was added to this data.

◆ unfoldedCovariance()

template<class Matrix >
const output_covariance_type& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::unfoldedCovariance ( ) const
inline

Last unfolded covariance matrix

◆ unfoldedDim()

template<class Matrix >
unsigned emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::unfoldedDim ( ) const
inline

Dimensionality of the unfolded space

◆ unfoldedResult()

template<class Matrix >
const npstat::ArrayND<double>& emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::unfoldedResult ( ) const
inline

Last unfolded distribution

◆ usingConvolutions()

template<class Matrix >
bool emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::usingConvolutions ( ) const
inline

Are we using convolutions with our filters?

◆ variableCount()

template<class Matrix >
virtual unsigned emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::variableCount ( ) const
virtual

Number of variables calculated by the "process" method. See the description of "variableNames" method for details.

◆ variableNames()

template<class Matrix >
virtual std::vector<std::string> emsunfold::SparseUnfoldingBandwidthScannerND< Matrix >::variableNames ( ) const
virtual

Names of the variables calculated by calling the "process" method. The meaning of the variables is described below. Variables which have * in front of their names will have meaningful values only in case the known "oracle" distribution was provided in the constructor.

bandwidth_N – With N = 0, 1, ..., up to filterParameters.size() - 1. Last bandwidth set processed. All these values will be set to -2.0 if "ntuplize" was called before the first call to "process".

unprunedNonZeroFrac – Fraction of non-zero entries in the unfolding covariance matrix before pruning of small correlations.

nonZeroFrac – Fraction of non-zero entries in the unfolding covariance matrix after pruning of small correlations.

foldedSum – Sum of the Poisson means fitted to the observations.

unfoldedSum – Sum of the unfolded values.

*smoothedOracleSum – Sum of the smoothed oracle data.

foldedLogli – Poisson log-likelihood for the observations assuming that the Poisson means are given by the fit. If the covariance matrix of the observations was provided by the user, this variable will be set to -chi^2/2 instead.

*unfoldedLogli – Poisson log-likelihood for the unfolded distribution assuming that the Poisson means are given by the oracle data.

*unfoldedISE – Integrated squared error for the unfolded distribution w.r.t. the oracle (both distributions are normalized to 1). This variable is the easiest one to check in order to find out whether an oracle was provided (its value will be negative if it wasn't).

*unfoldedDiagChisq – Chi-squared of the unfolded distribution w.r.t. the oracle using only the diagonal elements of the unfolding covariance matrix.

*smoothedUnfoldedLogli – Poisson log-likelihood for the unfolded distribution assuming that the Poisson means are given by the smoothed racle data.

*smoothedUnfoldedISE – Integrated squared error for the unfolded distribution w.r.t. the smoothed oracle (both distributions are normalized to 1).

*smoothedUnfoldedDiagChisq – Chi-squared of the unfolded distribution w.r.t. the smoothed oracle using only the diagonal elements of the covariance matrix.

modelNDoFEntropic – Effective NDoF for the matrix modelNDoFTrace H = K*E*(K*E)^T which plays a role similar to the hat matrix in regression. Here, E is the error propagation matrix for unfolding uniform observed distribution. "Entropic" means that the exponent of the eigenspectrum entropy is used to define effective NDoF, and "Trace" means that the ratio of the matrix trace to the largest eigenvalue is used.

AICcEntropic – AIC (Akaike information criterion) AICcTrace with a correction for the finite sample size. Calculated using "foldedLogli" and corresponding "modelNDoF".

smoothingNormfactor – Normalization factor applied during the last smoothing procedure.

integratedVariance – Product of the covariance matrix trace and the bin width in the unfolded space.

nIterations – Number of iterations used to process this bandwidth (e.g., by the expectation-maximization method).

unfoldingStatus – Status returned by the "unfold" call of the unfolding object.

trlan_* – Multiple variables which provide TRLAN diagnostic information. See the description of TrlanDiagnostics class.


The documentation for this class was generated from the following file: