System Identification / Modal Analysis#

Covariance-driven SSI (BRSSICovRef) and PoGER multi-setup identification (PogerSSICovRef).

class pyOMA.core.SSICovRef.BRSSICovRef(*args, **kwargs)[source]#

Bases: ModalBase

Reference-based Covariance-driven Stochastic Subspace Identification (SSI-Cov/Ref).

Builds a block-Toeplitz matrix from output cross-correlation functions, decomposes it via SVD, and identifies modal parameters at multiple model orders. The standard workflow is:

build_toeplitz_cov() — assemble and decompose the Toeplitz matrix.
compute_modal_params() — run the multi-order modal identification.
Pass the result to StabilCalc for stabilisation-diagram analysis.

Parameters:: prep_signals (PreProcessSignals) – Pre-processed signal object providing correlation functions and channel metadata.

build_toeplitz_cov(num_block_columns=None, num_block_rows=None, shift=0, num_blocks=None, training_blocks=None)[source]#

Builds a Block-Toeplitz Matrix of Covariances with varying time lags and decomposes it by a Singular Value decomposition.

  <-num_block_columns * n_r ->  _
[     R_m      R_m-1      ...      R_1    ]^
[     R_m+1    R_m        ...      R_2    ]num_block_rows * n_l
[     ...      ...        ...      ...    ]
[     R_2m-1   ...        ...      R_m    ]v

The total number of block columns and block rows should not exceed the maximum time lag of pre-computed correlation functions:

num_block_columns + num_block_rows + shift < prep_signals.m_lags

Parameters:

num_block_columns (integer, optional) – Number of block columns. By default, half the number of time lags are used
num_block_rows (integer, optional) – Number of block rows. By default it is set equal to num_block_columns
shift (integer, optional) – Allows the assembly of a shifted Block-Toeplitz matrix, s. t. the correlation function starting at shift is assembled into the block Toeplitz matrix
num_blocks (integer, optional) – The number of blocks to split the signal into for cross-validation. If given, correlation functions are (re-)estimated block-wise via prep_signals.corr_blackman_tukey(m_lags, n_segments=num_blocks, refs_only=True) and only training_blocks are averaged into the Toeplitz matrix; the remaining blocks are then available for synthesize_correlation()/compute_modal_params() via their validation_blocks argument. If not given (default), behaviour is unchanged: whatever correlation function is already cached in prep_signals (Welch or Blackman-Tukey, full signal) is used.
training_blocks (list, optional) – The selected blocks to use for system identification (=training). Only meaningful together with num_blocks. Defaults to all blocks.

Note

Unlike the block handling in PreProcessingTools.PreProcessSignals correlation estimation, num_blocks here does not persist any state on prep_signals beyond the per-block cache (corr_matrices_bt) that corr_blackman_tukey itself maintains – the training-block average is used locally and does not overwrite prep_signals.corr_matrix_bt.

compute_modal_params(max_model_order=None, max_modes=None, algo='svd', modal_contrib=True, validation_blocks=None)[source]#

Perform a multi-order computation of modal parameters. Successively calls

estimate_state(order, max_modes, algo)

modal_analysis(A,C)

synthesize_correlation(A,C, G), if modal_contrib == True

At ascending model orders, up to max_model_order. See the explanations in the the respective methods, for a detailed explanation of parameters.

Parameters:

max_model_order (integer, optional) – Maximum model order, where to interrupt the algorithm. If not given, it is min(num_channels * (num_block_rows + 1), num_reference_channels * num_block_columns)
validation_blocks (list, optional) – Only meaningful if build_toeplitz_cov() was called with num_blocks (cross-validation mode). Forwarded to synthesize_correlation() at every order when modal_contrib is True.

estimate_state(order, max_modes=None, algo='svd')[source]#

Compute the state matrix A, output matrix C and next-state-output covariance matrix G from the singular values and vectors of the block Toeplitz matrix, truncated at the requested order. Estimation of the state matrix can be performed by QR decomposition or Singular Value decomposition of the shifted observability matrix. If max_modes is specified, the singular value decomposition is truncated additionally, also known as Crystal Clear SSI.

Parameters:

order (integer, required) – Model order, at which the state matrices should be estimated
max_modes (integer, optional) – Maximum number of modes, that are known to be present in the signal, to suppress noise modes
algo (str, optional) – Algorithm to use for estimation of A. Either ‘svd’ or ‘qr’.

Returns:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)
G (numpy.ndarray) – next-state-output covariance matrix : Array of shape (order, num_ref_channels)

classmethod init_from_config(conf_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

This is a stub that must be fully reimplemented by every derived class. Derived implementations typically read analysis parameters (e.g. model order, frequency range) from conf_file, call the relevant computation methods, and return the populated object.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

classmethod load_state(fname, prep_signals)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

modal_analysis(A, C, rescale_fun=None)[source]#

Computes the modal parameters from a given state space model as described by Peeters 1999 and Döhler 2012. Mode shapes are scaled to unit modal displacements. Complex conjugate and real modes are removed prior to further processing. Typically, order // 2 modes are in the returned arrays.

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)

Returns:

modal_frequencies ((order,) numpy.ndarray) – Array holding the modal frequencies for each mode
modal_damping ((order,) numpy.ndarray) – Array holding the modal damping ratios (0,100) for each mode
mode_shapes ((n_l, order,) numpy.ndarray) – Complex array holding the mode shapes
eigenvalues ((order,) numpy.ndarray) – Complex array holding the eigenvalues for each mode

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

synthesize_correlation(A, C, G, validation_blocks=None)[source]#

Correlation function synthetization in a modal decoupled form follows Reynders-2012-SystemIdentificationMethodsFor(Operational)ModalAnalysisReviewAndComparison Eq. 161 p. 74 (24) where Lambda are the correlation functions of the identified system

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)
G (numpy.ndarray) – next-state-output covariance matrix : Array of shape (order, num_ref_channels)
validation_blocks (list, optional) – Only meaningful if build_toeplitz_cov() was called with num_blocks (cross-validation mode). The selected blocks whose (block-wise, Blackman-Tukey) correlation function is used as ground truth for computing modal contributions, instead of prep_signals.corr_matrix. Defaults to all blocks (matching the default of training_blocks in build_toeplitz_cov() – pass disjoint sets for a held-out validation).

Returns:

corr_matrix_synth ((n_l, n_r, m_lags, n_modes) numpy.ndarray) – Array holding the modally decomposed correlation functions for each channel n_l and reference channel n_r and all modes
modal_contributions ((order,) numpy.ndarray) – Array holding the contributions of each mode to the input correlation function.

synthesize_spectrum(A, C, G)[source]#

L = N*dt (duration = number_of_samples*sampling_period) P = N*df (maximal frequency = number of samples * frequency inverval)

dt * df = 1/N L * P = N

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

class pyOMA.core.SSICovRef.PogerSSICovRef[source]#

Bases: BRSSICovRef

Post-Global Estimation and Re-scaling (PoGER) multi-setup SSI-Cov/Ref.

Merges correlation functions from multiple measurement setups by stacking them into a joint subspace matrix, then identifies global modal parameters via SSI-Cov/Ref. Partial mode shapes are rescaled in a least-squares sense to the reference DOFs of the first setup.

The standard workflow is:

For each setup, create a PreProcessSignals object with pre-computed correlation functions and call add_setup().
pair_channels() — match channels across setups.
build_merged_subspace_matrix() — assemble the joint Toeplitz matrix.
compute_modal_params() — run the multi-order modal identification.

Notes

There are two distinct roles for reference channels:

Joint identification — channels shared across all setups, used to stack correlation functions (auto-determined from channel-DOF assignments).
Mode-shape rescaling — always relative to the first setup added; ensure the first setup has the best-quality measurements at the reference DOFs.

Todo

Add modal contributions
Implement PreGER merging with variance computation in a new class

References

Döhler, M. et al. “Pre- and post-identification merging for multi-setup OMA with covariance-driven SSI.” IMAC-XXVIII, 2010, pp. 57-70.

add_setup(prep_signals)[source]#: todo: check that ref_channels are equal in each setup (by number and by DOF)

build_merged_subspace_matrix(num_block_columns, num_block_rows=None)[source]#

Build and SVD-decompose the merged block-Toeplitz subspace matrix.

Stacks the correlation functions from all added setups into a joint block-Toeplitz matrix and decomposes it via SVD.

  <- num_block_columns * num_ref_channels -> _
[     R_1      R_2      ...      R_i      ]^
[     R_2      R_3      ...      R_i+1    ]num_block_rows * (n_l * num_setups)
[     ...      ...      ...      ...      ]v
[     R_i      ...      ...      R_2i-1   ]_

R_k = [ R_k^{setup_1} ]
      [ R_k^{setup_2} ]
      [ ...           ]

Parameters:

num_block_columns (int) – Number of block columns in the joint Toeplitz matrix.
num_block_rows (int, optional) – Number of block rows. Defaults to num_block_columns.

compute_modal_params(max_model_order=None, max_modes=None, algo='svd')[source]#

Perform a multi-order computation of modal parameters. Successively calls

estimate_state(order, max_modes, algo)

modal_analysis(A,C)

synthesize_correlation(A,C, G), if modal_contrib == True

At ascending model orders, up to max_model_order. See the explanations in the the respective methods, for a detailed explanation of parameters.

Parameters:

max_model_order (integer, optional) – Maximum model order, where to interrupt the algorithm. If not given, it is min(num_channels * (num_block_rows + 1), num_reference_channels * num_block_columns)
validation_blocks (list, optional) – Only meaningful if build_toeplitz_cov() was called with num_blocks (cross-validation mode). Forwarded to synthesize_correlation() at every order when modal_contrib is True.

classmethod load_state(fname)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

modal_analysis(A, C)[source]#

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)

Returns:

modal_frequencies ((order,) numpy.ndarray) – Array holding the modal frequencies for each mode
modal_damping ((order,) numpy.ndarray) – Array holding the modal damping ratios (0,100) for each mode
mode_shapes ((n_l, order,) numpy.ndarray) – Complex array holding the mode shapes
eigenvalues ((order,) numpy.ndarray) – Complex array holding the eigenvalues for each mode

pair_channels()[source]#

pairs channels from all given setups for the poger merging methods

ssi_reference channels are common to all setups rescale reference channels are common to at least two setups

finds common dofs from all setups and their respective channels generates new channel_dof_assignments with ascending channel numbers rescale reference channels are assumed to be equal to ssi_reference channels

rescale_by_references(mode_shape)[source]#

This is PoGer Rescaling

extracts each setup’s reference and roving parts of the modeshape

compute rescaling factor from all setup’s reference channels using a least-squares approach

rescales each setup’s roving channels and assembles final modeshape vector

reference channel_pairs and final channel-dof-assignments have been determined by function pair_channels note: reference channels for SSI need not necessarily be reference channels for rescaling and vice versa

$S_\phi \times \alpha = [n \times 1, 0 .. 0]$

$\phi^{ref}_i$ : Reference-sensor part of modeshape estimated from setup $i = 0 .. n$ $j_{max} = \operatorname{argmax}(\Pi_i |\phi^{ref}_i|)$ : maximal modal component in all setups → will be approximately scaled to 1, must belong to the same sensor in each setup

\[\begin{split}S_\phi = \begin{bmatrix} \phi^{ref}_{0,j_{max}}& \phi^{ref}_{1,j_{max}}& ..& ..& \phi^{ref}_{n,j_{max}} \\ \phi^{ref}_0& -\phi^{ref}_1& 0& ..& 0 \\ \phi^{ref}_0& 0& -\phi^{ref}_2& ..& 0 \\ . &. &. & . & . \\ . &. &. & . & . \\ \phi^{ref}_0& 0& 0& ..& -\phi^{ref}_n \\ 0& \phi^{ref}_1& -\phi^{ref}_2& ..& 0 \\ . &. & . & . & . \\ . &. & . & . & . \\ 0& \phi^{ref}_1& 0& ..& -\phi^{ref}_n \\ . &. & . & . & . \\ . &. & . & . & . \\ 0& 0& \phi^{ref}_2& ..& -\phi^{ref}_n \\ . &. & . & . & . \\ . &. & . & . & . \\ 0& 0& 0& \phi^{ref}_{n-1}& -\phi^{ref}_n \end{bmatrix}\end{split}\]

if references are the same in all setups

dimensions $= 1 + (n_{setups} ! )* n_{ref_{channels}} \times n_{setups}$

not quite exact, since different setups may share different references

→ list based assembly of the $S_\phi$ matrix

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

Multi-setup SSI with pre-estimation re-scaling (PreGER).

This module implements the pre-global-estimation-re-scaling (PreGER) merging strategy of Döhler, Lam & Mevel (2013, MSSP 36(2), doi:10.1016/j.ymssp.2012.06.009) for multi-setup operational modal analysis. In contrast to the post-estimation re-scaling of PogerSSICovRef, each setup’s observability is normalised to a common modal basis before the global system identification. The approach is theoretically valid under changing ambient excitation across setups, requires no per-mode mode-shape gluing, and admits a rigorous first-order covariance computation (implemented in VarPreGERSSI).

Two classes are provided:

PreGERSSI – point estimates of the modal parameters.
VarPreGERSSI – additionally propagates measurement uncertainty to first-order covariances of the modal parameters (see the class for status).

The merging and identification chain is kept agnostic to how each setup’s subspace matrix is built (covariance-driven block Hankel or, in a later phase, projection/data-driven), mirroring the subspace_method switch of VarSSIRef.

References

Döhler, M., Lam, X.-B. & Mevel, L. “Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements.” Mechanical Systems and Signal Processing 36(2), 2013, pp. 562-581.

class pyOMA.core.MultiSetupSSI.PreGERSSI[source]#

Bases: ModalBase

Pre-Global Estimation and Re-scaling (PreGER) multi-setup SSI.

Merges correlation functions from multiple measurement setups by normalising each setup’s observability matrix to the reference-DOF basis of the first setup before the global identification, then extracts the state and output matrices from the merged observability matrix (Döhler, Lam & Mevel 2013).

The standard workflow is:

For each setup, create a PreProcessSignals object with pre-computed correlation functions and call add_setup().
pair_channels() – match channels across setups and determine the common reference DOFs.
build_subspace_matrices() – build and decompose the per-setup block-Hankel subspace matrices.
compute_modal_params() – run the multi-order modal identification.

Notes

The reference channels play two roles, both auto-determined from the channel-DOF assignments:

Merging – the DOFs common to all setups (and chosen as reference channels in every setup) define the basis into which each setup’s observability is normalised. Their correlation functions form the columns of every setup’s block-Hankel matrix.
Row ordering – each setup’s output channels are physically reordered so the common-reference channels come first, which lets the merged output matrix C be assembled by plain block-row selection.

Unlike PoGER, mode shapes are globally scaled by construction; no least-squares re-scaling of partial mode shapes is performed.

add_setup(prep_signals)[source]#

Validates consistency with previously added setups (type, channel-DOF assignments, sampling rate) and stores a deep copy of the setup’s correlation functions. The correlation properties of PreProcessSignals are mutable, _last_meth-switched caches, so a deep copy is mandatory to keep each setup’s data stable.

Parameters:: prep_signals (PreProcessSignals) – Pre-processed signal object with pre-computed correlation functions (e.g. via corr_blackman_tukey(m_lags=...)) and channel-DOF assignments.

build_subspace_matrices(num_block_columns, num_block_rows=None, subspace_method='covariance', num_blocks=None)[source]#

Build and decompose the per-setup subspace matrices H^(j).

Each setup j yields a subspace matrix H^(j) of shape ((p+1)*n_l_j, cols) whose (p+1) output block rows (reordered common-references-first) carry the observability structure. For the first setup its SVD is stored; for every subsequent setup the moving-row matrix and the SVD of the reference-row matrix are stored – these feed the observability merging in _observability_parts(). The merging and downstream identification are agnostic to how H^(j) is built.

Two constructions are supported:

'covariance' – a block-Hankel of output/reference correlations, block (i, k) = corr[:, :, i+k] (cols = q*r; see notes).
'projection' – the data-driven UPC subspace matrix H^dat from the two-pass LQ decomposition of the raw data Hankels (cols = p*r, requires q == p; Döhler paper Appendix A).

Parameters:

num_block_columns (int) – Number of block columns q. For 'projection' it must equal num_block_rows (it is coerced with a warning otherwise).
num_block_rows (int, optional) – p (there are p+1 output block rows). Defaults to num_block_columns.
subspace_method ({'covariance', 'projection'}, optional)
num_blocks (int, optional) – Projection method only: number of blocks the signal is split into for the two-pass LQ (defaults to each setup’s n_segments if correlations were computed, else 50).

Notes

pyOMA’s corr_matrix[:, :, k] holds the correlation at time lag k (0-based), and BRSSICovRef treats that very k=0 entry as the paper’s R_1 slot. The covariance method follows the same convention, so the merged Hankel spans the same column space as the single-setup Toeplitz – a prerequisite for the Ns=1 equivalence with BRSSICovRef.

compute_modal_params(max_model_order=None)[source]#

Multi-order modal identification.

Successively calls estimate_state() and modal_analysis() at ascending model orders and stores the results in the ModalBase array conventions.

Parameters:: max_model_order (int, optional) – Maximum model order. Defaults to (and is capped at) the order cap min((p+1)*r, q*r) determined by build_subspace_matrices().

estimate_state(order)[source]#

Estimate the state matrix A and output matrix C at order.

Assembles the shifted-up/shifted-down merged observability matrices from _observability_parts() and solves A = pinv(O_up) @ O_down. C is the stack of first block rows of every part, already in the merged channel order.

Returns:

A (numpy.ndarray, shape (order, order))
C (numpy.ndarray, shape (merged_num_channels, order))

classmethod init_from_config(conf_file, prep_signals_list)[source]#

Build and identify a merged multi-setup model from a config file.

Unlike the single-setup init_from_config classmethods (which take one already-built PreProcessSignals object), a PreGER merger has no single setup to construct from – setups are registered one at a time via add_setup(). prep_signals_list therefore takes their place: a sequence of already pre-processed PreProcessSignals objects, one per setup, each typically built via init_from_config() beforehand (see scripts/multi_setup_analysis_preger.py). Calling this on VarPreGERSSI requires each of them to additionally carry block-wise correlations (corr_blackman_tukey(m_lags, n_segments=...)), which VarPreGERSSI.add_setup() enforces.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals_list (sequence of PreProcessSignals) – One pre-processed setup per element, in the order they should be added.

Returns:

Populated instance (cls(), so a VarPreGERSSI.init_from_config call returns a VarPreGERSSI).

Return type:

PreGERSSI

classmethod load_state(fname)[source]#: Restore a PreGERSSI from an archive written by save_state().

modal_analysis(A, C)[source]#

Extract modal parameters from a merged state-space model.

Standard SSI modal analysis: eigendecomposition of A, mode shapes C @ eigvecs, removal of complex-conjugate/real poles, frequency and damping from the discrete eigenvalues, integration to modal displacements and phase-rotation normalisation. No PoGER-style re-scaling is applied – PreGER mode shapes are globally scaled by construction.

Returns:

modal_frequencies, modal_damping ((order,) numpy.ndarray)
mode_shapes ((merged_num_channels, order) numpy.ndarray, complex)
eigenvalues ((order,) numpy.ndarray, complex)

pair_channels()[source]#

Pair channels across setups and determine common reference DOFs.

Re-uses the static channel-pairing helpers of PogerSSICovRef, with two PreGER-specific additions:

The common reference DOFs are filtered to those chosen as reference channels in every setup (so all setups contribute the same r columns), and their per-setup lengths are checked for equality.
A physical row permutation setup_row_orders is stored for every setup (common references first, in ssi_ref_dofs order, then the remaining channels in ascending order). This lets each setup’s Hankel matrix be split into reference/moving parts by plain block-row selection and lets C land in the merged channel order.

save_state(fname)[source]#: Save the current computation state to a compressed NumPy archive.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

class pyOMA.core.MultiSetupSSI.VarPreGERSSI(cache_variance_factors=False, cache='full', cache_dtype=<class 'numpy.float32'>)[source]#

Bases: PreGERSSI

PreGER multi-setup SSI with first-order modal-parameter uncertainties.

Extends PreGERSSI by propagating the (co)variance of the estimated correlation functions to first-order standard deviations of the modal frequencies, damping ratios and mode shapes.

Method (factored perturbation propagation, Döhler, Lam & Mevel 2013). Each setup j is split into n_b_j blocks; the per-block Hankel deviations ΔH^(j)_k = (H^(j)_k - H̄^(j)) / sqrt(n_b (n_b - 1)) are treated as concrete perturbation directions (columns). Every column is propagated linearly through the whole identification chain – SVD perturbation of the subspace matrix (Lemma 6), observability merging (Lemma 2/3), the least-squares state matrix, and the per-mode eigenvalue/eigenvector sensitivities – to a modal-parameter perturbation Δ(·)_k. Because the setups (and blocks) are independent, the covariance is the sum of outer products of the columns, i.e. the variance of any modal quantity is the sum of squares of its column perturbations (Theorem 4).

Two propagation branches arise:

A setup-0 column perturbs H^(0) and hence O_ref; it ripples into every setup’s normalised observability (ΔŌ^(j,mov) = H^(j,mov) pinv_n(H^(j,ref)) ΔO_ref).
A setup-``j>=1`` column perturbs only H^(j), giving ΔŌ^(j,mov) = ΔH^(j,mov) (pinv_n O_ref) + H^(j,mov) (Δpinv_n O_ref) via the truncated-pseudo-inverse perturbation (Eq. C.5).

Grouping all three product-rule terms this way sidesteps the ambiguous j=1 term of the paper’s Eq. (25)/Algorithm 1.

Block weighting. The blocks of every setup can carry non-uniform weights, in the two flavours of VarSSIRef and VarPLSCF. Weights are always given per setup – one weight vector per setup, None for an unweighted one – because each setup’s blocks are averaged into that setup’s own subspace matrix:

Build time (build_subspace_matrices() with weights=): setup j’s Hankel matrix becomes the weighted mean of its blocks and the deviation columns are re-centred on it and scaled by sqrt(w_k) / sqrt(n_eff - 1) (Kish’s n_eff = 1 / sum(w**2)). This moves the point estimate, which is the only honest way to discount a contaminated block.
Post hoc (compute_modal_params_weighted() for the recompute path, apply_block_weights() for the cached one): the point estimates and every Jacobian stay at their original linearisation and only the std_* arrays are recomputed, by right-multiplying the already-centred factors with W(w) (see _block_weight_factor()). Requires an unweighted build; it is the delta-method covariance of the reweighted estimator around the original point estimate, first-order consistent for moderate weight changes but unable to relocate a point estimate contaminated by a bad block.

Both flavours reduce to the unweighted estimator at uniform weights, and 'substitution' (the default post-hoc convention) reproduces a fresh build-time weighted run’s covariance factor exactly.

weights#

Build-time per-setup block weights (each entry a renormalised (n_b_j,) array or None for a uniformly weighted setup), or None for an entirely unweighted build.

Type:: list or None

n_eff#

Per-setup effective block count (Kish’s n_eff), equal to setup_n_b[j] for a uniformly weighted setup.

Type:: list or None

block_weights#

The post-hoc weights the current std_* arrays reflect, per setup, or None for uniform. Distinct from weights, which records the build-time weights that moved the point estimate.

Type:: list or None

block_weight_convention#

The convention the current std_* arrays were reweighted under; see _block_weight_factor().

Type:: str or None

U_fixi_cache, U_phii_cache

Per-mode uncertainty factors of every evaluated order, keyed by order and stacked as (n_modes, 2, K) / (n_modes, 2*n_l, K) over the K = sum_j n_b_j perturbation columns; populated by compute_modal_params() with cache_variance_factors=True and consumed by apply_block_weights().

Type:: dict or None

Notes

Like the point-estimate classes, the mode shapes are integrated from modal accelerations/velocities to modal displacements (integrate_quantities()); the omega-dependence of that integration is propagated into the mode-shape standard deviations via the Delta f -> Delta omega coupling. Normalisation is then the max-amplitude phase rotation (Prop. 5), as in VarSSIRef, rather than the unit re-scaling of PreGERSSI.modal_analysis(); the two therefore differ by a per-mode complex scale, while frequencies and damping ratios are identical.

add_setup(prep_signals)[source]#

For the covariance method the setup must carry block-wise correlation estimates (corr_matrices with n_segments >= 2, e.g. from corr_blackman_tukey(m_lags, n_segments=...)). For the projection method the raw signals kept by PreGERSSI.add_setup() suffice; the per-block data Hankels are formed at build time. Whichever is present is stored; the applicable requirement is enforced in build_subspace_matrices().

apply_block_weights(weights=None, convention='substitution')[source]#

Refresh the standard deviations under new block weights, from the caches.

Reweights the cached per-mode uncertainty factors (U_fixi_cache / U_phii_cache, populated by compute_modal_params() with cache_variance_factors=True) by U @ W(w) and overwrites the std_* entries in place, so a weight sweep over one identification costs a matrix product per order instead of a modal loop. This is the algebraic equivalent of compute_modal_params_weighted() – W is real and the factors are linear in the block columns – up to the cache_dtype precision (float32 by default). The point estimates are untouched.

With cache='freqdamp' only std_frequencies and std_damping are refreshed; the mode-shape standard deviations keep whatever weighting they were computed under.

There is no automatic fallback: without a cache a RuntimeError is raised, and compute_modal_params_weighted() is the recompute path.

Parameters:

weights (sequence or None) – Per-setup post-hoc block weights; None restores the uniform standard deviations. See compute_modal_params_weighted().
convention ({'substitution', 'reliability', 'precision'}, optional) – See _block_weight_factor().

build_subspace_matrices(num_block_columns, num_block_rows=None, subspace_method='covariance', num_blocks=None, weights=None)[source]#

Build the subspace matrices and the per-block deviation factors.

The per-block Hankel deviations ΔH^(j)_k = (H^(j)_k - H̄^(j)) / sqrt(n_b (n_b - 1)) are formed from block-wise correlations (covariance method) or from the per-block data Hankels of the two-pass LQ (projection method); the perturbation chain is identical thereafter.

Parameters:

num_block_columns – As in PreGERSSI.build_subspace_matrices().
num_block_rows – As in PreGERSSI.build_subspace_matrices().
subspace_method – As in PreGERSSI.build_subspace_matrices().
num_blocks – As in PreGERSSI.build_subspace_matrices().
weights (sequence, optional) –
Build-time block weights, one entry per setup: a (n_b_j,) vector of non-negative per-block weights, or None for a uniformly weighted setup. Setup j’s subspace matrix then becomes the weighted mean of its blocks – which moves the point estimate – and its deviation columns are re-centred on that mean and scaled by sqrt(w_k) / sqrt(n_eff - 1), reducing to the unweighted 1 / sqrt(n_b (n_b - 1)) at uniform weights. Pass None (default) for the classical unweighted estimator; to change the weights after identification instead, leave this out and use apply_block_weights() or compute_modal_params_weighted().

For subspace_method='projection' the weights enter only the final combination of the per-block H^dat estimates: each block’s own two-pass LQ is carried out exactly as in the unweighted build, so the block deviations stay independent and the variance chain applies unchanged (the same reading as build_subspace_mat()’s default; its experimental pre-LQ weighting has no counterpart here).

compute_modal_params(max_model_order=None, orders=None, cache_variance_factors=None, cache=None, cache_dtype=None)[source]#

Multi-order modal identification with first-order uncertainties.

Parameters:

max_model_order (int, optional) – Maximum model order (capped at min((p+1)*r, q*r)).
orders (iterable of int, optional) – Explicit model orders to evaluate; defaults to 1 .. max_model_order-1.
cache_variance_factors (bool, optional) – Keep the per-mode uncertainty factors of every evaluated order, so that apply_block_weights() can refresh the standard deviations from them without recomputing the perturbation chain. None (default) uses the constructor’s setting.
cache ({'full', 'freqdamp'}, optional) – What to cache; None uses the constructor’s setting. See __init__().
cache_dtype (numpy dtype, optional) – Storage precision of the caches; None uses the constructor’s setting.

compute_modal_params_weighted(weights, convention='substitution', max_model_order=None, orders=None)[source]#

Recompute the standard deviations under post-hoc block weights.

Reweights the (order-independent) triplet perturbations and the block deviations by the per-setup W(w) of _block_weight_factor() and re-runs the multi-order loop against them. Because every step from a block deviation to a modal-parameter perturbation is real-linear and homogeneous in the block columns, with no coupling between blocks, reweighting the factors once up front is equivalent to reweighting every downstream quantity: L(F) @ W == L(F @ W).

The expensive Lemma-6 pseudo-inverses are reused from the cache, so this costs one modal loop; apply_block_weights() is cheaper still for sweeping weights, and this method is its oracle. The per-mode caches are deliberately left alone: they hold unweighted factors and stay valid across this call.

The point estimates are recomputed identically – the weights enter only the covariance. To move the point estimates too, pass weights to build_subspace_matrices() instead.

Parameters:

weights (sequence or None) – Post-hoc block weights, one entry per setup (a (n_b_j,) vector, or None for a uniformly weighted setup); None restores the uniform result of compute_modal_params().
convention ({'substitution', 'reliability', 'precision'}, optional) – Covariance-normalisation convention; see _block_weight_factor(). The default 'substitution' reproduces a build-time weighted run.
max_model_order – As in compute_modal_params().
orders – As in compute_modal_params().

classmethod load_state(fname)[source]#: Restore a VarPreGERSSI (view only; caches are not recomputed).

save_state(fname)[source]#

Save state including standard deviations; perturbation caches are not persisted.

Neither the Lemma-6 triplet perturbations nor the block deviations hankel_devs they are built from are stored (they are large and rebuildable), so compute_modal_params_weighted() is unavailable to a loaded object until the subspace matrices are rebuilt. The much smaller per-mode uncertainty factors are stored when they exist, so that apply_block_weights() survives a round trip.

Data-driven SSI (SSIData, SSIDataMC) for operational modal analysis.

class pyOMA.core.SSIData.SSIData(*args, **kwargs)[source]#

Bases: SSIDataMC

Data-driven SSI without correlation synthesis (non-Monte-Carlo variant).

Identical workflow to SSIDataMC but skips the synthesis step (synth_sig=False), making it faster when variance estimation is not required.

compute_modal_params(max_model_order)[source]#

Perform a multi-order computation of modal parameters. Successively calls

estimate_state(order,)

modal_analysis(A,C)

at ascending model orders, up to max_model_order. See the explanations in the the respective methods, for a detailed explanation of parameters.

Parameters:: max_model_order (integer, optional) – Maximum model order, where to interrupt the algorithm. If not given, it is determined from the previously computed subspace matrix.

estimate_state(order, max_modes=None, algo='svd')[source]#

Compute the state matrix A and output matrix C from the singular values and vectors of the projection matrix, truncated at the requested order. Estimation of the state matrix can be performed by QR decomposition or Singular Value decomposition of the shifted observability matrix. If max_modes is specified, the singular value decomposition is truncated additionally, also known as Crystal Clear SSI.

Parameters:

order (integer, required) – Model order, at which the state matrices should be estimated
max_modes (integer, optional) – Maximum number of modes, that are known to be present in the signal, to suppress noise modes
algo (str, optional) – Algorithm to use for estimation of A. Either ‘svd’ or ‘qr’.

Returns:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)

class pyOMA.core.SSIData.SSIDataCV(*args, **kwargs)[source]#

Bases: SSIDataMC

build_block_hankel(num_block_rows=None, num_blocks=1, training_blocks=None, reduced_projection=True)[source]#

Builds serveral Block-Hankel Matrices of the measured time series with varying time lags and estimates the subspace matrices from their LQ decompositions. Uniqueness of the subspace estimates is ensured by an intermediate LQ decomposition, where the diagonals of the L matrices are constrained to positive values. A subspace matrix estimate is computed by the mean over the training blocks leaving any remainig blocks for validation.

Note: Blocks are not completely i.i.d. as we borrow p+q timesteps from the previous block for the projection of a full block (assembly of Hankel matrix)

Todo

investigate correct scaling of the subspace matrices [sqrt(N_b), sqrt(N_b * num_blocks), sqrt(N_b*n_training_blocks)] ?
use sparse SVD (scipy.sparse.svds) to save memory and cpu time

Parameters:

num_block_rows (integer, required) – The number of block rows of the Subspace matrix
num_blocks (integer, optional) – The number of blocks, used for cross-validation
training_blocks (list, optional) – The selected blocks to use for system identification (=training)

synthesize_signals(A, C, Q, R, S, validation_blocks=None, N_offset=None, **kwargs)[source]#

Computes the modal response signals and the contribution of each mode. The algorithm follows Peeters 1999 and the Lyapunov equation is solved as a discrete-time algebraic Riccati equation (DARE). For long signals, the computation may become time-consuming, thus only time steps up to j may be used to synthesize the signal.

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)
Q (numpy.ndarray) – state noise covariance matrix: Symmetric array of shape (order, order)
R (numpy.ndarray) – signal noise covariance matrix: Array of shape (num_analised_channels, num_analised_channels)
S (numpy.ndarray) – system noise - signal noise covariance matrix: Array of shape (order, num_analised_channels)
validation_blocks (list, optional) – The selected blocks to be synthethized and used for system validation.
N_offset (integer, optional) – The number of samples to be used from any previous block for Kalman-Filter startup.

Returns:

sig_synth (numpy.ndarray, shape (num_analised_channels, j, order // 2)) – Array holding the modally decomposed input signals for each channel n_l and all modes
modal_contributions (numpy.ndarray, shape (order, )) – Array holding the contributions of each mode to the input signals.

class pyOMA.core.SSIData.SSIDataMC(*args, **kwargs)[source]#

Bases: ModalBase

Reference-based Data-driven Stochastic Subspace Identification (SSI-Data/MC).

Identifies modal parameters from raw time-series data by constructing a block-Hankel matrix, computing its LQ decomposition, and extracting state- space models via SVD. The standard workflow is:

build_block_hankel() — build the projection matrix from raw data.
compute_modal_params() — run the multi-order modal identification.
Pass the result to StabilCalc for stabilisation-diagram analysis.

Parameters:

prep_signals (PreProcessSignals) – Pre-processed signal object providing raw time-series and channel metadata.
TODO:: (..) –
- define unit tests to check functionality after changes
- add switch to keep synthesized time-histories

build_block_hankel(num_block_rows=None, reduced_projection=True)[source]#

Builds a Block-Hankel Matrix of the measured time series with varying time lags and estimates the subspace matrix from its LQ decomposition.

  <- num_time samples - num_block_rows->      _
[     y_0      y_1      ...      y_(j-1)     ]^
[     y_1      y_2      ...      y_j         ]num_block_rows (=i)*n_l
[     ...      ...      ...      ...         ]v
[     y_(2i-1)   y_(2i)  ...     y_(2i+j-2)  ]_

The notation mostly follows Peeters 1999.

Parameters:: num_block_rows (integer, required) – The number of block rows of the Subspace matrix

compute_modal_params(max_model_order=None, j=None, validation_blocks=None, synth_sig=True)[source]#

Perform a multi-order computation of modal parameters. Successively calls

estimate_state(order,)

modal_analysis(A,C)

synthesize_signals(A, C, Q, R, S, j)

at ascending model orders, up to max_model_order. See the explanations in the the respective methods, for a detailed explanation of parameters.

Parameters:: max_model_order (integer, optional) – Maximum model order, where to interrupt the algorithm. If not given, it is determined from the previously computed subspace matrix.

estimate_state(order)[source]#

Estimate the state matrices A, C and noise covariances Q, R and S from the subspace / projection matrix. Several methods exist, e.g.

Peeters 1999 Reference Based Stochastic Subspace Identification for OMA

DeCock 2007 Subspace Identification Methods

the algorithm used in BRSSICovRef.

Here, the first algorithm, a residual-based computation of Q, R and S, is implemented.

Parameters:

order (integer, required) – The model order, at which to truncate the singular values of the projection Matrix P_i_ref

Returns:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)
Q (numpy.ndarray) – state noise covariance matrix: Symmetric array of shape (order, order)
R (numpy.ndarray) – signal noise covariance matrix: Array of shape (num_analised_channels, num_analised_channels)
S (numpy.ndarray) – system noise - signal noise covariance matrix: Array of shape (order, num_analised_channels)

classmethod init_from_config(conf_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

classmethod load_state(fname, prep_signals)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

modal_analysis(A, C, rescale_fun=None)[source]#

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)

Returns:

modal_frequencies (numpy.ndarray, shape (order,)) – Array holding the modal frequencies for each mode
modal_damping (numpy.ndarray, shape (order,)) – Array holding the modal damping ratios (0,100) for each mode
mode_shapes (numpy.ndarray, shape (n_l, order,)) – Complex array holding the mode shapes
eigenvalues (numpy.ndarray, shape (order,)) – Complex array holding the eigenvalues for each mode

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

synthesize_signals(A, C, Q, R, S, j=None, **kwargs)[source]#

Parameters:

A (numpy.ndarray) – State matrix: Array of shape (order, order)
C (numpy.ndarray) – Output matrix: Array of shape (num_analised_channels, order)
Q (numpy.ndarray) – state noise covariance matrix: Symmetric array of shape (order, order)
R (numpy.ndarray) – signal noise covariance matrix: Array of shape (num_analised_channels, num_analised_channels)
S (numpy.ndarray) – system noise - signal noise covariance matrix: Array of shape (order, num_analised_channels)
j (integer, optional) – length of signal to synthesize (number of timesteps)

Returns:

sig_synth (numpy.ndarray, shape (num_analised_channels, j, order // 2)) – Array holding the modally decomposed input signals for each channel n_l and all modes
modal_contributions (numpy.ndarray, shape (order, )) – Array holding the contributions of each mode to the input signals.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

Covariance-driven SSI with propagated parameter variances (VarSSIRef).

class pyOMA.core.VarSSIRef.VarSSIRef(prep_signals, cache_variance_factors=False)[source]#

Bases: ModalBase

Covariance-driven SSI with first-order perturbation variance estimation.

Extends BRSSICovRef with analytical uncertainty propagation from measurement noise through the correlation functions, Toeplitz matrix, SVD, and eigendecomposition to the final modal parameters. Both covariance-based and projection-based subspace estimation are supported.

The standard workflow is:

build_subspace_mat() — build the subspace matrix and its statistical properties.
compute_state_matrices() — estimate state and output matrices.
prepare_sensitivities() — pre-compute sensitivity matrices for variance propagation.
compute_modal_params() — identify modal parameters with variances.

Parameters:

prep_signals (PreProcessSignals) – Pre-processed signal object providing correlation functions and channel metadata.
TODO:: (..) –
- define unit tests to check functionality after changes
- optimize multi-order QR-based estimation routine
- add mode-shape integration with variances
- use Monte-Carlo sampling in the last step of variance propagation

apply_block_weights(weights=None, convention='substitution')[source]#

Millisecond post-hoc block reweighting from the Tier A cache.

Reweights the cached per-mode variance factors (self.U_fixi_cache/self.U_phii_cache, populated by compute_modal_params() with cache_variance_factors=True) by U @ W(w) and overwrites the std_frequencies/std_damping (and, for the 'full' cache, std_mode_shapes) entries in place. This is the algebraic equivalent of compute_modal_params_weighted() (Tier B) – U @ W equals the recomputed reweighted factor because W is real and the factors are linear in the block columns – but costs a handful of small matrix products rather than a full modal loop, at the price of the cache memory (and the cache_dtype precision, float32 by default). The point estimates are untouched.

weights=None restores the uniform (unweighted) stds. Requires an unweighted build (self.weights is None); a warning is emitted once the effective sample size n_eff drops below ~10. See _block_weight_factor() for the convention parameter and the frozen-linearization caveat, and compute_modal_params_weighted() for the weaker block-independence caveat under subspace_method='projection'.

There is no automatic fallback: if no cache is present (Tier A disabled, or an archive saved without it) a RuntimeError is raised; use compute_modal_params_weighted() (Tier B) for the recompute path.

build_subspace_mat(num_block_columns, num_block_rows=None, num_blocks=None, subspace_method='covariance', weights=None, corr_matrices=None, hankel_matrices=None, hankel_provider=None, experimental_weighted_projection=False)[source]#

Builds a Block-Hankel Matrix of Covariances with varying time lags

R_1 R_2 … R_q |

R_2 R_3 … R_q+1 |

… … … … |

R_p+1 … … R_p+q |

Parameters weights, corr_matrices (covariance method) and hankel_matrices (projection method) turn the estimator into a weighted statistic over independent per-block estimates:

weights(num_blocks,) array, optional: Non-negative probability weights, one per block; renormalized to sum to one. None (default) keeps the classical unweighted (uniform) averaging. For subspace_method='projection', the weights enter only the final block combination by default: each block’s own per-block/joint-LQ projection is built exactly as in the unweighted case, then combined as a weighted mean instead of a plain one – the same structure as the covariance method’s weighted mean, so variance_algo='fast' works normally. See experimental_weighted_projection for the alternative that weights the per-block normalization itself.
corr_matrices(num_blocks, n_l, n_r, >=m_lags) array, optional: Externally computed correlation estimates, one per block (e.g. one per aleatory sample), replacing the correlation functions of the attached prep_signals. Unlike the internal path, no block-count inflation is applied, since each entry is a standalone estimate rather than a fragment of one shared signal. subspace_method='covariance' only.
hankel_matricessequence of (n_r*p + n_l*(p+1), block_length) arrays, optional: Externally provided raw past/future block-Hankel matrices (stacked [Y_minus; Y_plus], p = num_block_rows), one per block (e.g. one per aleatory sample), bypassing prep_signals.signals; analogous to corr_matrices for the covariance method. Blocks may have different block_length (each is scaled by its own length). subspace_method='projection' only.
hankel_providercallable, optional: Memory-frugal alternative to hankel_matrices: a callable hankel_provider(n_block) -> array returning the raw block-Hankel matrix for block n_block on demand. Each block is built, LQ-reduced and freed one at a time, so the full-size Hankels never coexist (peak memory ~one block instead of num_blocks of them); the result is identical to passing the same matrices via hankel_matrices. Requires num_blocks; mutually exclusive with hankel_matrices. subspace_method='projection' only.
experimental_weighted_projectionbool, optional: Enable the experimental weighted-least-squares reading of the projection (subspace_method='projection' only): each block Hankel matrix is scaled by its weight before the per-block LQ decompositions, so both the per-block and the joint R11 normalization see the weights. Point estimates only – variance computation is rejected for this build, because the per-block weighting is entangled through the joint R11 LQ step and breaks the independent-per-block-deviation assumption the fast/slow Hankel covariance relies on. False (default) gives the weighted-mean build described under weights instead, which does support variance. Uniform weights reproduce the unweighted build either way.

compute_modal_params(max_model_order=None, debug=False, qr=True, orders=None, cache_variance_factors=None, cache='full', cache_dtype=<class 'numpy.float32'>)[source]#

Compute modal parameters with variance estimation.

Parameters:

orders (list of int, optional) – Restrict the evaluation to these model orders (each in 1..max_model_order-1); result rows for all other orders stay zero. Saves most of the per-order eigendecomposition/Jacobian cost when only a single sampled order is of interest.
cache_variance_factors (bool, optional) – Whether to cache the per-mode variance factors U_fixi/U_phii (Tier A) for millisecond post-hoc reweighting via apply_block_weights(). None (default) uses the object default set in the constructor.
cache ({'full', 'freqdamp'}, optional) – What to cache when caching is on: 'full' stores both the frequency/damping factor U_fixi and the mode-shape factor U_phii; 'freqdamp' stores only U_fixi (so apply_block_weights reweights frequency/damping stds but leaves the mode-shape stds untouched), saving the bulk of the memory.
cache_dtype (numpy dtype, optional) – Storage dtype for the caches (default np.float32; use np.float64 for exact Tier A == Tier B agreement).

compute_modal_params_weighted(weights, convention='substitution', max_model_order=None, orders=None, debug=False)[source]#

Recompute modal-parameter variances for post-hoc block weights (Tier B).

Reruns the modal-order loop with the cached fast variance factors (Q1..Q4 / J_OHT) right-multiplied by the block-weight matrix W(w) (see _block_weight_factor()), so the std_* arrays are those of the reweighted estimator. No SVD, sensitivity preparation, or extra memory is needed – the cost is a single modal loop (seconds). The point estimates (frequencies, damping, mode shapes) are recomputed identically to compute_modal_params(); only the variances change.

Works for both subspace_method='covariance' and 'projection' (the fast Q-pipeline is identical). Requires an unweighted build: the frozen-linearization identity is defined relative to the uniform-mean cached factors, so self.weights must be None.

Note

For subspace_method='projection' the per-block samples are defined after the joint LQ normalization (the lq_decomp of the stacked R11 matrices couples the blocks). Holding that normalization fixed is consistent with the frozen-linearization semantics below, but the block independence the reweighting assumes is weaker than in the covariance method, so projection reweighting is the more approximate of the two.

Warning

This is a frozen-linearization (delta-method) reweighting: the point estimates and Jacobians stay at their original weighting. It is first-order consistent for moderate weight changes but does NOT relocate a point estimate contaminated by a bad block – for that use the build-time weighted path. A warning is emitted once the effective sample size n_eff drops below ~10.

Parameters:

weights ((num_blocks,) array or None) – Non-negative block weights (renormalized to sum to one). None restores the uniform (unweighted) result.
convention ({'substitution', 'reliability', 'precision'}) – Covariance-normalization convention for the reweighting scalar; see _block_weight_factor(). Default 'substitution' matches a fresh build-time-weighted run.
max_model_order – As in compute_modal_params().
orders – As in compute_modal_params().
debug – As in compute_modal_params().

compute_state_matrices(max_model_order=None, lsq_method='pinv')[source]#: computes the state and output matrix of the state-space-model by applying a singular value decomposition to the block-hankel-matrix of covariances the state space model matrices are obtained by appropriate truncation of the svd matrices at max_model_order the decision whether to take merged covariances is taken automatically

classmethod init_from_config(conf_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

classmethod load_state(fname, prep_signals)[source]#: Load a previously saved state from a compressed NumPy archive.

prepare_sensitivities(variance_algo='fast', debug=False)[source]#

Prepare Jacobians and covariance matrices for variance propagation.

variance_algo='none' skips all variance preparation and marks the object ready for a point-estimates-only modal run (all std_* arrays stay zero). This is the only mode available for a projection build made with experimental_weighted_projection=True; the default (False) weighted-mean projection build supports variance_algo='fast' normally.

save_state(fname)[source]#: Save the current object state to a compressed NumPy archive.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

pyOMA.core.VarSSIRef.vectorize(matrix)[source]#

\[A=\begin{bmatrix} 1 & 2 & 3 \ 4 & 5 & 6 \ 7 & 8 & 9 \ \end{bmatrix}\]

returns vertically stacked columns of matrix A

..math:

\begin{bmatrix}
\
\
\
\
\
\
\
\
\
\end{bmatrix}

Eigensystem Realisation Algorithm (ERA) for operational modal analysis.

class pyOMA.core.ERA.ERA(prep_signals)[source]#

Bases: object

Eigensystem Realisation Algorithm (ERA) for operational modal analysis.

Identifies state-space models and modal parameters from impulse-response or free-decay data by constructing a Hankel matrix and decomposing it via SVD. When forced-response data are available, call CalculateFRF() first to convert them to impulse responses before calling build_hankel_matrix().

Parameters:: prep_signals (PreProcessSignals) – Pre-processed signal object providing signals, sampling_rate, and channel metadata.

CalculateFRF()[source]#: function by anil FRF(Frequency response function) is convertion of signal from time to frequency domain. The following function performs this conversion.

build_hankel_matrix(num_block_columns)[source]#

Construct the shifted Hankel matrix from the impulse-response functions.

Parameters:: num_block_columns (int) – Number of block columns in the Hankel matrix. The number of block rows is set to num_block_columns + 1.

compute_state_matrices(max_model_order=None)[source]#

Decompose the Hankel matrix and compute the observability matrix.

Parameters:: max_model_order (int, optional) – Maximum model order to retain. When None, the full rank of the Hankel matrix is used.

classmethod load_state(fname, prep_signals)[source]#

Restore an ERA object from a previously saved archive.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored object with all previously computed results.

Return type:

ERA

static remove_conjugates_new(eigval, eigvec_r, eigvec_l=None)[source]#

removes conjugates

eigvec_l.shape = [order+1, order+1] eigval.shape = [order+1,1]

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Parameters:: fname (str) – Destination file path (without .npz extension).

Poly-reference Least-Squares Complex Frequency (pLSCF) identification method.

class pyOMA.core.PLSCF.LSFDContext(A: ndarray, X: ndarray, pinv_A: ndarray, residual: ndarray, Df1: ndarray, Df2: ndarray, eigenvalues: ndarray, participation_vectors: ndarray, mode_order: ndarray = None)[source]#

Bases: object

Assembly of the least-squares frequency-domain fit for the mode shapes.

Holds the intermediate quantities that PLSCF._fit_mode_shapes_ls() discards, so that uncertainty propagation can differentiate the fit without duplicating it.

A#

Real design matrix, shape (num_omega * 2 * n_r, 2 * n_modes + 4 * n_r).

Type:: np.ndarray

X#

Least-squares solution pinv(A) @ h, shape (2 * n_modes + 4 * n_r, n_l).

Type:: np.ndarray

pinv_A#

Pseudo-inverse of A. Note pinv_A @ pinv_A.T is (A^T A)^-1 for a full-column-rank A, which saves forming and inverting the normal equations of this stage a second time.

Type:: np.ndarray

residual#

h - A @ X, shape (num_omega * 2 * n_r, n_l). Only the residual is retained; h itself is a reordering of pos_half_spectra.

Type:: np.ndarray

Df1, Df2

Pole factors 1 / (1j * omega - lambda) and its conjugate-pole counterpart, for every frequency line, shape (num_omega, n_modes).

Type:: np.ndarray

eigenvalues#

Continuous-time poles the fit was built on, shape (n_modes,).

Type:: np.ndarray

participation_vectors#

Normalised participation vectors the fit was built on, shape (n_r, n_modes). Held together with eigenvalues so that the context fully describes the fit; both are in mode-column order, not returned order.

Type:: np.ndarray

mode_order#

Permutation relating the mode columns of A to the returned modes: returned mode k occupies mode column mode_order[k]. The fit is run in the order the poles were passed in, which is not the returned frequency-sorted order. Assigned by PLSCF.modal_analysis_residuals().

Type:: np.ndarray

class pyOMA.core.PLSCF.ModalContext(A_c: ndarray, eigvals_z: ndarray, eigvecs_l: ndarray, eigvecs_r: ndarray, mode_indices: ndarray)[source]#

Bases: object

Companion-matrix eigendecomposition and the mode selection applied to it.

A_c#

Transposed companion matrix, shape (order * n_r, order * n_r).

Type:: np.ndarray

eigvals_z#

Discrete-time eigenvalues remaining after remove_conjugates().

Type:: np.ndarray

eigvecs_l, eigvecs_r

Left and right eigenvectors of A_c, column-aligned with eigvals_z.

Type:: np.ndarray

mode_indices#

Columns of the eigen-arrays corresponding to the returned modes, in the returned (in-band, frequency-sorted) order. Lets derived quantities reproduce the identical mode selection and ordering.

Type:: np.ndarray

class pyOMA.core.PLSCF.NormalEquationsContext(order: int, X_o: ndarray, RS_solutions: ndarray, M: ndarray, M_aa: ndarray, alpha: ndarray, beta_l_i: ndarray)[source]#

Bases: object

Assembly of the reduced normal equations at a single model order.

Holds the intermediate quantities that PLSCF.estimate_model() discards, so that uncertainty propagation can differentiate the assembly without duplicating it.

order#

Model order this assembly was built at.

Type:: int

X_o#

Polynomial basis, shape (num_omega, order + 1).

Type:: np.ndarray

RS_solutions#

Per-output-channel R_o^-1 S_o, shape (order + 1, (order + 1) * n_r, n_l).

Type:: np.ndarray

M#

Reduced normal equations, shape ((order + 1) * n_r, (order + 1) * n_r).

Type:: np.ndarray

M_aa#

The block of M that is solved for the free denominator coefficients, M[:order * n_r, :order * n_r] (a view).

Type:: np.ndarray

alpha, beta_l_i

Denominator and numerator coefficients; see PLSCF.estimate_model().

Type:: np.ndarray

class pyOMA.core.PLSCF.PLSCF(*args, **kwargs)[source]#

Bases: ModalBase

Poly-reference Least-Squares Complex Frequency (pLSCF) method.

Also known as PolyMAX. Identifies modal parameters from positive half-spectra derived from correlation functions. The standard workflow is:

build_half_spectra() — construct the positive half-spectra.
compute_modal_params() — run the multi-order identification.
Pass the result to StabilCalc for stabilisation-diagram analysis.

Parameters:

prep_signals (PreProcessSignals) – Pre-processed signal object providing correlation functions and channel metadata.
TODO:: (..) –
- Test functions should be added to the test package

build_half_spectra(nperseg=None, begin_frequency=None, end_frequency=None, window_decay=0.001, num_blocks=None, training_blocks=None, weights=None, corr_matrices=None, **kwargs)[source]#

Extracts an array of positive half spectra between begin_frequency and end_frequency from a spectrum of nperseg frequency lines. If begin_frequency > 0.0 or end_frequency<nyquist freqeuncy, the resulting array has less than nperseg lines.

Positive power spectra are constructed from positive correlation functions, that are windowed by an exponential window and transformed to frequency domain by and (R)FFT. Correlation functions are computed in prep_signals by either Welch’s or Blackman-Tukey’s method, though, Welch’s method is not recommmended, because the artificial damping introduced by windowing can not be corrected.

See: Cauberghe-2004-Applied Frequency-Domain System … : Sections 3.4ff

Note: The previous implementation contained severe mistakes in the computation of positive power spectra, e.g. doubled squaring of spectral values, lazy handling of array dimensions and therefore effectively only a quarter of nperseg being used as well as numerical inefficiencies.

Todo

Move spectral estimation into prep_signals.pds_blackman_tukey and only keep bandwidth selection and argument checking here
Allow other windows than exponential

Parameters:

nperseg (integer, optional) – Number of (positive) frequency lines to consider (rfft)
begin_frequency (float, optional) – Frequency range to restrict the identified system.
end_frequency (float, optional) – Frequency range to restrict the identified system.
window_decay (float, (0,1)) – Final value of the exponential window, that is applied to the correlation functions.
num_blocks (integer, optional) – The number of blocks to split the signal into for cross-validation. If given, correlation functions are (re-)estimated block-wise via prep_signals.corr_blackman_tukey(nperseg, n_segments=num_blocks, refs_only=True) and only training_blocks are averaged into the half-spectrum; the remaining blocks are then available for synthesize_spectrum()/compute_modal_params() via their validation_blocks argument. If not given (default), behaviour is unchanged: whatever correlation function is already cached in prep_signals (Welch or Blackman-Tukey, full signal) is used.
training_blocks (list, optional) – The selected blocks to use for system identification (=training). Only meaningful together with num_blocks. Defaults to all blocks.
weights ((len(training_blocks),) array_like, optional) – Non-negative weights, one per training block, in the order of training_blocks – not per block of num_blocks. The two coincide in the default case, where all blocks train. The half-spectrum is then built from the weighted mean sum_j w_j R_j of the block correlation functions rather than their plain mean. Weights are renormalised to sum to one, so their scale carries no meaning, and uniform weights reproduce the unweighted estimate exactly. Requires num_blocks.
corr_matrices ((num_blocks, n_l, n_r, >= nperseg) array_like, optional) – Externally supplied block correlation functions, bypassing prep_signals.corr_blackman_tukey. Lets each block be an independent realization rather than a segment of one record. Requires num_blocks. Note that save_state() does not persist them.
kwargs – Additional kwargs are passed to prep_signals.correlation

compute_modal_params(max_model_order, complex_coefficients=False, algo='residuals', modal_contrib=None, validation_blocks=None)[source]#

Perform a multi-order computation of modal parameters. Successively calls

estimate_model(order, complex_coefficients)

modal_analysis_residuals(alpha, beta_l_i) or modal_analysis_state_space(alpha, beta_l_i)

synthesize_spectrum(alpha, beta_l_i), if modal_contrib == True

At ascending model orders, up to max_model_order. See the explanations in the the respective methods, for a detailed explanation of parameters.

Parameters:

max_model_order (integer) – Maximum model order, where to interrupt the algorithm.
complex_coefficients (bool, optional) – Whether to estimate a real or complex RMFD model
algo (str, optional) – Algorithm to use for modal analysis. Either ‘state-space’ or ‘residuals’ Both algorithms are approximately equally fast. The state space based algorithm seems to yield less complex mode shapes.
modal_contrib (bool, optional) – Synthesize modal spectra and estimate modal contributions. Only to be used with residual-based modal analysis algorithm.
validation_blocks (list, optional) – Only meaningful if build_half_spectra() was called with num_blocks (cross-validation mode). Forwarded to synthesize_spectrum() at every order when modal_contrib is True.

estimate_model(order, complex_coefficients=False)[source]#

Estimate a right matrix-fraction model from positive half-spectra, by constructing a set of reduced normal equations as shown in Peeters 2004. The polynomial is identified following Cauberghe 2004. Sec. 5.2.1

Verboven 2002: Sect. 5.3.3 has a discussion on the use of real or complex valued coefficients, favoring complex ones. Guillaume 2003, Peeters 2004 just assume real coefficients, while later references, e.g. Cauberghe 2004, Reynders 2012 use complex coefficients. However, with complex coefficients, stabilization diagrams seem to become corrupted.

Note: The previous implementation was wrong in the estimation of alpha coefficients and led to “bad” stabilization. Additionally there was a wrong sign in the assembly of the C_c matrix, which led to corrupted mode shapes.

Todo

implement weighting function; c.p. Peeters 2004 Sect. 2.2
improve assembly by exploiting the Toeplitz structure of S, R, T; c.p. Cauberghe 2004 Eq. 5.17ff
Investigate LS-TLS solution by using a SVD
estimate polynomial once at highest order and construct all lower order models from these coefficients; c.p. Peeters 2004 Sect. 2.4
Check, if alternative solution for lpha in Reynders 2012. Sec. 5.2.4 leads to clearer stabilization, or it it is actually equivalent to the current implementation

Parameters:

order (integer, required) – Model order, at which the RMF model should be estimated
complex_coefficients (bool, optional) – Whether to assume real or complex coefficients

Returns:

alpha (numpy.ndarray) – Denominator coefficients: Array of shape ((order + 1) * n_r, n_r)
beta_l_i (numpy.ndarray) – Numerator coefficients: Array of shape (order + 1, n_r, n_l)

classmethod init_from_config(conf_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

classmethod load_state(fname, prep_signals)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

modal_analysis_residuals(alpha, *args)[source]#

Perform a modal analysis of the identified polyomial with the least-squares residual-based method as outlined in Steffensen-2025-VarianceEstimation… Sect. 2.1 Mode shapes are scaled to unit modal displacements. Complex conjugate and real modes are removed prior to further processing. Damping values are corrected, if half-spectra were constructed with an exponential window.

Todo

numerical optimization to increase speed

Parameters:

alpha (numpy.ndarray) – Denominator coefficients: Array of shape ((order + 1) * n_r, n_r)

Returns:

modal_frequencies ((order * n_r,) numpy.ndarray) – Array holding the modal frequencies for each mode
modal_damping ((order * n_r,) numpy.ndarray) – Array holding the modal damping ratios (0,100) for each mode
mode_shapes ((n_l, order * n_r,) numpy.ndarray) – Complex array holding the mode shapes
eigenvalues ((order * n_r,) numpy.ndarray) – Complex array holding the _eigenvalues for each mode

modal_analysis_state_space(alpha, beta_l_i)[source]#

Perform a modal analysis of the identified polyomial by converting it into a state-space model, as outlined in Reynders-2012: Lemma 2.2, followed by an eigendecomposition. Mode shapes are scaled to unit modal displacements. Complex conjugate and real modes are removed prior to further processing. Damping values are corrected, if half-spectra were constructed with an exponential window.

Todo

numerical optimization to increase speed

Parameters:

alpha (numpy.ndarray) – Denominator coefficients: Array of shape ((order + 1) * n_r, n_r)
beta_l_i (numpy.ndarray) – Numerator coefficients: Array of shape (order + 1, n_r, n_l)

Returns:

modal_frequencies ((order * n_r,) numpy.ndarray) – Array holding the modal frequencies for each mode
modal_damping ((order * n_r,) numpy.ndarray) – Array holding the modal damping ratios (0,100) for each mode
mode_shapes ((n_l, order * n_r,) numpy.ndarray) – Complex array holding the mode shapes
eigenvalues ((order * n_r,) numpy.ndarray) – Complex array holding the eigenvalues for each mode

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

synthesize_spectrum(alpha, beta_l_i, modal=True, validation_blocks=None)[source]#

Spectral synthetization in a modal decoupled form follows Steffensen-2025-VarianceEstimation… Sect. 2.1.2 The spectral synthetization without modal decomposition follows Peeters-2004-ThePolyMAX…

Todo

numerical optimization to increase speed

Parameters:

alpha (numpy.ndarray) – Denominator coefficients: Array of shape ((order + 1) * n_r, n_r)
beta_l_i (numpy.ndarray) – Numerator coefficients: Array of shape (order + 1, n_r, n_l)
modal (bool, optional) – Synthesize a spectrum for each mode and its modal contribution to the full spectrum
validation_blocks (list, optional) – Only meaningful if build_half_spectra() was called with num_blocks (cross-validation mode) and modal is True. The selected blocks whose (block-wise, Blackman-Tukey) half-spectrum is used as ground truth for computing modal contributions, instead of self.pos_half_spectra. Defaults to all blocks (matching the default of training_blocks in build_half_spectra() – pass disjoint sets for a held-out validation).

Returns:

half_spec_modal ((n_l, n_r, num_omega, n_modes) numpy.ndarray) – Array holding the (modally decomposed) synthesized positive half spectra for each channel n_l and reference channel n_r and all modes
modal_contributions ((order,) numpy.ndarray) – Array holding the contributions of each mode to the input spectrum

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

Variance estimation for the pLSCF method.

Propagates the uncertainty of the positive half-spectra onto the identified modal parameters by the delta method, following Steffensen, Döhler, Tcherniak & Thomsen (MSSP 223, 2025), “Variance estimation of modal parameters from the poly-reference least-squares complex frequency-domain algorithm”.

The paper enters the propagation with the covariance of measured FRFs, obtained from the H1 estimator and the multiple coherence (their Eq. 26-27). pyOMA is output-only, so no input spectrum and hence no coherence exists, and that entry point has no counterpart here. Instead the half-spectra are re-estimated block-wise from the block correlation functions that corr_blackman_tukey() already provides, and their sample covariance over blocks is used. This is the one deliberate deviation from the paper’s derivation; all Jacobians downstream are the paper’s.

The covariance is carried as a factor rather than a matrix: the block deviations are scaled such that cov(h) = B B^T for the re/im-stacked slice B of VarPLSCF.spec_block_factor, at the scale of the covariance of the estimate (the role the paper’s 1/N_avg plays). Every downstream perturbation is then evaluated per block and variances follow as einsum('ij,ij->i', U, U), mirroring VarSSIRef. A consequence worth noting: this retains the correlations across frequency lines and output channels that the paper’s Eq. 40 discards for tractability – with an empirical covariance, imposing that independence would be extra work rather than less.

class pyOMA.core.VarPLSCF.Stage1Scores(U_theta: ndarray, d_lambda: ndarray, d_frequency: ndarray, d_damping: ndarray, d_participation: ndarray)[source]#

Bases: object

First-order perturbations of the Stage-1 quantities, one column per block.

Every array carries the block index as its last axis: column j is the estimator’s linearisation applied to the deviation of block j, scaled such that contracting the array with itself over that axis yields a variance.

U_theta#

Free denominator coefficients, shape (order * n_r, n_r, n_b), real.

Type:: np.ndarray

d_lambda#

Continuous-time eigenvalues, shape (n_modes, n_b), complex.

Type:: np.ndarray

d_frequency, d_damping

Modal frequencies and damping ratios, shape (n_modes, n_b), real.

Type:: np.ndarray

d_participation#

Normalised participation vectors, shape (n_r, n_modes, n_b), complex.

Type:: np.ndarray

Modes are ordered as :meth:`~pyOMA.core.PLSCF.PLSCF.modal_analysis_residuals`

returns them.

class pyOMA.core.VarPLSCF.VarPLSCF(*args, cache_variance_factors=False, cache='full', cache_dtype=<class 'numpy.float32'>, **kwargs)[source]#

Bases: PLSCF

pLSCF identification with uncertainty quantification.

Extends PLSCF with standard deviations of the identified modal parameters. The workflow matches the base class, except that build_half_spectra() must be given num_blocks:

build_half_spectra() — half-spectra plus their block-covariance factor.
compute_modal_params() — modal parameters and their standard deviations.

spec_block_factor#

Covariance factor of the positive half-spectra, complex, shape (n_l, n_r, num_omega, n_b) with the block index last. Column j holds the deviation of block j from the block mean, scaled by 1 / sqrt(n_b * (n_b - 1)), so that contracting the re/im-stacked slice with itself yields the covariance of the mean half-spectrum.

Not persisted by save_state(): it is reconstructible from prep_signals and would dominate the archive size.

Type:: np.ndarray or None

std_frequencies, std_damping

Standard deviations of the identified modal frequencies and damping ratios, shaped like modal_frequencies.

Type:: np.ndarray or None

std_participation_vectors#

Standard deviations of the participation vectors, shaped like participation_vectors; the standard deviation of the real part is stored in .real and that of the imaginary part in .imag, as in VarSSIRef.

Type:: np.ndarray or None

block_weights#

The post-hoc block weights the current std_* arrays reflect, renormalised to sum to one, or None for uniform. Distinct from weights, which records the build-time weights that moved the point estimate: post-hoc reweighting (apply_block_weights(), compute_modal_params_weighted()) leaves the point estimates alone and requires an unweighted build.

Type:: np.ndarray or None

block_weight_convention#

The convention the current std_* arrays were reweighted under; see _block_weight_factor().

Type:: str or None

Notes

Post-hoc reweighting keeps the point estimates and every Jacobian at their original linearisation: it is the delta-method covariance of the reweighted estimator around the original point estimate, first-order consistent for moderate weight changes. It does not relocate a point estimate contaminated by a bad block – passing weights to build_half_spectra() is the honest route for that.

apply_block_weights(weights=None, convention='substitution')[source]#

Refresh the standard deviations under new block weights, from the caches.

Reweights the cached per-mode uncertainty factors instead of re-identifying, so a weight sweep over one identification costs a matrix product per order. Requires compute_modal_params(cache_variance_factors =True) to have run; there is no automatic fallback to compute_modal_params_weighted(), which is this method’s oracle.

The point estimates and every Jacobian stay at their original linearisation – see the class notes. weights=None restores the uniform standard deviations.

With cache='freqdamp' only std_frequencies and std_damping are refreshed; the participation and mode-shape standard deviations keep whatever weighting they were computed under.

Parameters:

weights ((n_b,) array_like or None) – Non-negative per-training-block weights; see _block_weight_factor(). None restores uniform weighting.
convention (str, optional) – One of 'substitution' (default), 'reliability' or 'precision'; see _block_weight_factor().

Construct the positive half-spectra and their block-covariance factor.

Behaves as build_half_spectra(), but num_blocks is required: the variance of the half-spectra is estimated from the scatter of the individual blocks about their mean, and that mean is the half-spectrum the model is identified from.

Passing weights propagates them into the covariance factor as well as the point estimate; the covariance then scales as 1 / n_eff, which generalises the 1 / N_avg of Steffensen-2025 Eq. 26 (with N_avg the number of equally weighted Welch averages) and recovers it at uniform weights. To change the weights after identification instead, leave them out here and use apply_block_weights().

Parameters:

nperseg (integer, optional) – Number of (positive) frequency lines to consider (rfft)
begin_frequency (float, optional) – Frequency range to restrict the identified system.
end_frequency (float, optional) – Frequency range to restrict the identified system.
window_decay (float, (0,1)) – Final value of the exponential window, that is applied to the correlation functions.
num_blocks (integer, required) – The number of blocks to split the signal into. At least two blocks are needed for a variance estimate; the estimate becomes meaningful only for substantially more.
training_blocks (list, optional) – The selected blocks to use for system identification (=training). Defaults to all blocks. The covariance factor is built from these blocks only, i.e. from the blocks that entered the point estimate.
weights ((len(training_blocks),) array_like, optional) – Non-negative per-training-block weights; see build_half_spectra(). At least two blocks must carry weight.
corr_matrices ((num_blocks, n_l, n_r, >= nperseg) array_like, optional) – Externally supplied block correlation functions; see build_half_spectra(). Lets each block be an independent realization, which is the data model the block covariance assumes anyway.
kwargs – Additional kwargs are passed to prep_signals.correlation

compute_modal_params(max_model_order, complex_coefficients=False, algo='residuals', modal_contrib=None, validation_blocks=None, cache_variance_factors=None)[source]#

Perform a multi-order computation of modal parameters and their standard deviations.

Behaves as compute_modal_params(), but additionally propagates the block-covariance of the half-spectra onto the identified modal parameters at every order.

Parameters:

max_model_order (integer) – Maximum model order, where to interrupt the algorithm.
complex_coefficients (bool, optional) – Must be False: the propagation is derived for real coefficients.
algo (str, optional) – Must be ‘residuals’: the propagation differentiates the residual-based modal analysis.
modal_contrib (bool, optional) – Synthesize modal spectra and estimate modal contributions.
validation_blocks (list, optional) – Forwarded to synthesize_spectrum().
cache_variance_factors (bool, optional) – Cache the per-mode uncertainty factors for apply_block_weights(). Defaults to the constructor’s setting.

compute_modal_params_weighted(max_model_order, weights, convention='substitution', complex_coefficients=False, algo='residuals', modal_contrib=None, validation_blocks=None)[source]#

Recompute the standard deviations under block weights, from scratch.

Reweights the cached covariance factor and re-runs the whole multi-order loop against it. Because every step from the block factor to the mode shapes is real-linear and homogeneous in the block deviations, applied per block column with no coupling between blocks, reweighting the factor once up front is equivalent to reweighting every downstream quantity: L(F) @ W == L(F @ W).

This is the memory-free reference path, not a fast one. Unlike VarSSIRef, pLSCF has no separate sensitivity-preparation stage, so re-running the loop repeats the whole identification and saves only the block-spectrum FFTs. Prefer apply_block_weights() for sweeping weights; this method is its oracle.

The point estimates are recomputed identically – the weights enter only the covariance. To move the point estimate too, pass weights to build_half_spectra() instead.

Parameters:

max_model_order (integer) – Maximum model order, where to interrupt the algorithm.
weights ((n_b,) array_like or None) – Non-negative per-training-block weights; None is uniform, which reproduces compute_modal_params(). See _block_weight_factor().
convention (str, optional) – One of 'substitution' (default), 'reliability' or 'precision'; see _block_weight_factor().
complex_coefficients – As in compute_modal_params().
algo – As in compute_modal_params().
modal_contrib – As in compute_modal_params().
validation_blocks – As in compute_modal_params().

classmethod init_from_config(conf_file, prep_signals)[source]#

Build and identify a VarPLSCF from a config file.

Overrides init_from_config(): that version never supplies num_blocks, which build_half_spectra() requires here to estimate the half-spectrum covariance. Cache options (cache_variance_factors, cache, cache_dtype) are constructor-time kwargs, left at their defaults, exactly as init_from_config() leaves its own cache_variance_factors kwarg unexposed.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

Poly-reference Complex Exponential (PRCE) identification method.

class pyOMA.core.PRCE.PRCE(*args, **kwargs)[source]#

Bases: ModalBase

Poly-reference Complex Exponential (PRCE) identification method.

Identifies modal parameters from a 3-D tensor of cross-correlation functions using the Complex Exponential approach. The standard workflow is:

build_corr_tensor() — assemble the correlation tensor.
compute_modal_params() — run the multi-order identification.
Pass the result to StabilCalc for stabilisation-diagram analysis.

Parameters:: prep_signals (PreProcessSignals) – Pre-processed signal object.

build_corr_tensor(num_corr_samples)[source]#: Builds a 3D Tensor of cross correlation functions with the following directions: 1 - related to reference channels 2 - all channels 3 - time

compute_modal_params(max_model_order)[source]#: Compute modal parameters for all model orders up to max_model_order.

classmethod init_from_config(mod_ID_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

classmethod load_state(fname, prep_signals)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.

Base class shared by all pyOMA system-identification methods.

class pyOMA.core.ModalBase.ModalBase(prep_signals=None)[source]#

Bases: object

Base class from which all pyOMA system-identification classes inherit.

Provides shared functionality (conjugate removal, mode-shape integration, rescaling, persistence) so that derived classes only implement the method-specific identification steps. Post-processing tools (stabilization diagram, mode-shape plot) accept any ModalBase subclass instance.

prep_signals#

The signal object from which this analysis was created.

Type:: PreProcessSignals or None

setup_name#

Human-readable label for the measurement setup.

Type:: str

start_time#

Timestamp of the measurement.

Type:: datetime.datetime or None

num_analised_channels#

Total number of analysis channels.

Type:: int or None

num_ref_channels#

Number of reference channels.

Type:: int or None

max_model_order#

Maximum model order used in the identification.

Type:: int or None

modal_frequencies#

Identified natural frequencies (Hz), shape (max_model_order, n_modes).

Type:: np.ndarray or None

modal_damping#

Identified modal damping ratios (%), same shape as modal_frequencies.

Type:: np.ndarray or None

mode_shapes#

Identified mode shapes, shape (n_channels, n_modes, max_model_order).

Type:: np.ndarray or None

eigenvalues#

Identified (complex) eigenvalues.

Type:: np.ndarray or None

classmethod init_from_config(conf_file, prep_signals)[source]#

Initialise a modal analysis object from a text configuration file.

Parameters:

conf_file (str) – Path to a tab-separated key-value configuration file compatible with ConfigFile.
prep_signals (PreProcessSignals) – Pre-processed signal object for this setup.

Returns:

Populated subclass instance.

Return type:

ModalBase

static integrate_quantities(vector, accel_channels, velo_channels, omega)[source]#

Rescales mode shapes from modal accelerations / velocities to modal displacements, by multiplication of the relevant modal coordinates (where accelerometers, or velocimeters were used, with $-1 omega^2$ or $i omega$, respectively,

Parameters:

vector ((n_channels,) numpy.ndarray) – Complex modeshape for all n_channels
accel_channels (list) – A list containing the channel numbers of all acceleration channels
velo_channels (list) – A list containing the channel numbers of all velocity channels
omega (float) – The circular frequency of the corresponding mode ($omega = 2 pi f$)

Returns:

vector – Rescaled complex modeshape for all n_channels

Return type:

(n_channels,) numpy.ndarray

classmethod load_state(fname, prep_signals)[source]#

Restore a modal-analysis object from a previously saved archive.

Must be fully reimplemented by every derived class.

Parameters:

fname (str) – Path to the .npz archive written by save_state().
prep_signals (PreProcessSignals) – Signal object for the same setup; used to validate the archive.

Returns:

Restored subclass instance.

Return type:

ModalBase

Raises:

NotImplementedError – Always, unless overridden by a derived class.

static remove_conjugates(eigval, eigvec_r=None, eigvec_l=None, inds_only=False)[source]#

This method finds complex conjugate modes, and removes unstable and overdamped poles.

A complex conjugate is defined as: $\lambda_i = \overline{\lambda_j} \text{ for } i \neq j$

Unstable poles, i.e. negatively damped poles, are defined by: $[\ln(|\lambda|)<0]: |\lambda_i|> 1$

Overdamped poles, are purely real poles: $[\operatorname{atan}(\Im/\Re)=0]: \Im(\lambda_i)=0$

The method keeps the second occurance of a conjugate pair (usually the one with the negative imaginary part) and either returns a truncated set of eigenvalues and eigenvectors or a list of (physical) poles that can be iterated.

Parameters:

eigval ((order,) numpy.ndarray) – Complex array of all eigenvalues
eigvec_r ((order, n_channels) numpy.ndarray, optional) – Complex array(s) of all right (left) eigenvectors
eigvec_l ((order, n_channels) numpy.ndarray, optional) – Complex array(s) of all right (left) eigenvectors
inds_only (bool, optional) – Whether to return a list of pole indices, or a reduced set of eigenvalues and eigenvectors

Returns:

conj_indices (list) – list of (physical) pole indices
eigval ((order,) numpy.ndarray) – Complex array of reduced (physical) eigenvalues
eigvec_l, eigvec_r ((order, n_channels) numpy.ndarray, optional) – Complex array(s) of reduced (physical) left (right) eigenvectors

static rescale_mode_shape(modeshape, rotate_only=False)[source]#

Rescales and rotates modeshapes in the complex plane. Default behaviour is to scale the larges component to unit modal displacement. If argument rotate_only is provided, the method given in Appendix C2 of Doehler 2013 (doi:0.1016/j.ymssp.2012.11.011) is used to rotate but not rescale the mode shape. Note: The scale of identified mode shapes is arbitrary in most OMA methods.

Parameters:

modeshape ((n_channels,) numpy.ndarray) – Complex modeshape for all n_channels
rotate_only (bool, optional) – Whether to rotate, but not rescale, the mode shape.

Returns:

modeshape – Rescaled complex modeshape for all n_channels

Return type:

(n_channels,) numpy.ndarray

save_state(fname)[source]#

Save the current computation state to a compressed NumPy archive.

Must be fully reimplemented by every derived class.

Parameters:: fname (str) – Destination file path (without .npz extension).
Raises:: NotImplementedError – Always, unless overridden by a derived class.

write_config(conf_file)[source]#

Write the analysis parameters used to compute this object to a text configuration file readable by init_from_config().

This is a stub that must be fully reimplemented by every derived class (mirrors init_from_config()): write the same keys that class’s own init_from_config() reads, via write().

Parameters:: conf_file (str) – Path to write the configuration file to.