Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Multiple Input/Multiple Output Sine Control

14Multiple Input/Multiple Output Sine Control

The MIMO Sine environment in Rattlesnake is used to generate deterministic sinusoidal excitations at one or more drive channels so that the measured response at one or more control channels matches a prescribed specification. Unlike the Random Vibration environment, which operates on power spectral density matrices, the Sine environment operates on explicitly defined amplitude, phase, frequency, and sweep-rate information over time.

At a high level, the MIMO Sine environment does three things:

  1. defines one or more sine sweeps as a specification,

  2. uses a system identification phase to estimate the transfer functions between drive and response channels,

  3. computes drive signals that should reproduce the specified response amplitudes and phases, and then updates those drives during a run based on measured tracking data.

The environment supports:

14.1Governing Equations

At any instant in a sine test, the desired response of the structure may be represented as a complex response vector

x(ω)=ax(ω)eiϕx(ω)\mathbf{x}(\omega) = \mathbf{a}_x(\omega)e^{i\boldsymbol{\phi}_x(\omega)}
(14.1)

where ax\mathbf{a}_x is the vector of desired response amplitudes and ϕx\boldsymbol{\phi}_x is the vector of desired response phases.

The measured structural dynamics are represented by a frequency response matrix

Hxv(ω)\mathbf{H}_{xv}(\omega)
(14.2)

relating excitation voltages v(ω)\mathbf{v}(\omega) to control responses x(ω)\mathbf{x}(\omega):

x(ω)=Hxv(ω)v(ω)\mathbf{x}(\omega) = \mathbf{H}_{xv}(\omega)\mathbf{v}(\omega)
(14.3)

The basic control problem is therefore to determine a drive vector v(ω)\mathbf{v}(\omega) such that the desired response is achieved. In a least-squares sense, a first estimate can be obtained using the pseudoinverse

v(ω)=Hxv+(ω)x(ω)\mathbf{v}(\omega) = \mathbf{H}_{xv}^{+}(\omega)\mathbf{x}(\omega)
(14.4)

where (⋅)+(\cdot)^+ denotes the pseudoinverse.

Rattlesnake computes these complex drive values across the active tones and then reconstructs real time-domain signals for the output channels. During the run, measured response amplitudes and phases are extracted from the acquired time histories, and the control law computes corrections to the drive amplitudes and phases.

Because the Sine environment operates on deterministic sinusoidal content rather than stationary random spectra, a great deal of the implementation revolves around:

14.2Specification Definition

The Sine environment specification is defined as one or more sine tones, each consisting of a set of frequency breakpoints and associated amplitude/phase information for each control channel.

Each tone may be a linear sweep or a logarithmic sweep depending on how the sweep type and sweep rate are specified between breakpoints. Currently, there is no “dwell” capability to pause at a specific frequency for a certain amount of time; however this can be approximated by using a slow sweep rate between two very close frequency breakpoints.

In Rattlesnake, the specifications are entered on the Environment Definition tab using one tab per sine tone. Each tab corresponds to a SineSpecification object internally.

A specification tone is defined by:

The environment can contain multiple such tones, which are combined in time to form the overall specification.

Specification information can be entered manually into the breakpoint table (see Section 14.3.7). Alternatively, the sine tone information may be loaded from a file.

14.2.1Specification File Format

The MIMO Sine environment accepts specification files in either NumPy archive format (*.npz) or MATLAB format (*.mat). These files define one sine tone at a time. A specification file loaded into one sine-tone tab of the environment must contain the information needed to construct a Sine Specification. If multiple sine tones are desired for a given environment, multiple files should be prepared and loaded individually in separate sine tone tabs in the Environment Definition tab (see Section 14.3.7).

Both NumPy and MATLAB files are structured identically, with the same field names and shapes associated with each field. Note that for 1D arrays, MATLAB can either specify n×1n \times 1 or 1×n1 \times n arrays; Rattlesnake will squeeze out the extra dimension.

For a specification where nfn_f is the number of frequency breakpoints and ncn_c is the number of control channels, the fields are:

For any array with dimension of size ncn_c, the ordering of this dimension must be identical to the ordering of the control degrees of freedom in the environment loading the file. No bookkeeping or reordering of specification data to match the channel data occurs in the Sine environment. If a transformation matrix is used, then the ordering of this dimension must be identical to the rows of the transformation matrix.

14.3Defining the MIMO Sine Environment in Rattlesnake

In addition to the specification, there are a number of sampling and signal processing parameters that are used by the MIMO Sine environment. These, along with the specification, are defined on the Environment Definition tab in the Rattlesnake controller, on a sub-tab corresponding to a MIMO Sine environment, as shown in Figure 14.1.

UI used to define a MIMO Random Vibration environment.

Figure 14.1:UI used to define a MIMO Random Vibration environment.

14.3.1Sampling Parameters

The Sampling Parameters section of the MIMO Sine definition sub-tab consists of the following parameters:

14.3.2Signal Generation Parameters

The Signal Generation Parameters section of the MIMO Sine definition sub-tab consists of the following parameters:

14.3.3Tracking Filter Parameters

The Tracking Filter Parameters section of the MIMO Sine definition sub-tab consists of properties used to specify the tracking filter, which is used to extract sine-tone information from the measured time histories. Depending on which filter type is selected, different parameters will be made available to the user.

With the digital tracking filter selected, the following parameters are available:

With the Vold-Kalman filter selected, the following parameters are available:

If the user is unsure as to which filter they should use, they can use the filter explorer to investigate filter performance, see Section 14.7.4.

14.3.4Control Parameters

The Control Parameters section of the MIMO Sine definition sub-tab consists of properties used to specify the control settings used by the controller. Like other environment types, MIMO Sine allows users to load Python scripts containing custom control laws. Unfortunately, MIMO Sine control laws are by their nature complex, so simple implementations like functions or generators are not available; only class definitions are available. See Section 14.11 for more information on defining a custom control law.

14.3.5Control Channels

The Control Channels section of the MIMO Sine definition sub-tab allows users to select the channels in the environment that will be used for control.

14.3.6Control and Drive Transforms

The Control and Drive Transforms section of the MIMO Sine definition sub-tab consists of properties used to specify the response and drive transformations, allowing virtual control degrees of freedom or drive signals to be constructed from physical channels.

Note that if transformation matrices are defined, the number of control channels ends up being the number of rows of the Response Transformation Matrix, rather than the number of physical control channels. The number of physical control channels will be equal to the number of columns of the transformation matrix. The number of control channels in the specification should be equal to the number of rows in the transformation.

See Section 12.8 for more information on specifying transformation matrices.

14.3.7Test Specification

The Test Specification section of the MIMO Sine definition sub-tab consists of properties used to specify the test specification. The main tabs of this section denote the individual sine tones present in this environment. Clicking on the tab labelled + will create a new sine tone. This section also includes four plots to display the specification data in different formats.

Within each sine tone, there are three tables behind a sub-tab interface. The Breakpoint Table tab shows the frequency breakpoints, as well as various buttons to add or remove breakpoints.

The Warning Table contains amplitude limits that will trigger warnings.

The Abort Table contains amplitude limits that will trigger aborts.

Both warning and abort limit tables allow for specifying the amplitude on the left and right side of the breakpoint independently. This can allow for, for example, looser tolerances during segments of the test that have higher sweep rates. Upper and lower tolerances can also be specified. The warning and abort limits will be shown in the plots along with the breakpoints themselves.

14.4System Identification for the MIMO Sine Environment

When all environments are defined and the Initialize Environments button is pressed, Rattlesnake will proceed to the next phase of the test, which is defined on the System Identification tab.

Before the Sine environment can predict or control the response, it must estimate the transfer function matrix between the drive channels and control channels. This is done through the shared system identification workflow used by the system-ID-based environments.

A typical system identification UI for the sine environment is shown in Figure 14.2.

System identification UI used by the MIMO Sine environment.

Figure 14.2:System identification UI used by the MIMO Sine environment.

Rattlesnake’s system identification phase will start with a noise floor check, where the data acquisition records data on all the channels without specifying an output signal. After the noise floor is computed, the system identification phase will play out the specified signals to the excitation devices, and transfer functions will be computed using the responses of the control channels to those excitation signals. Section 3.3 describes the System Identification tab and its various parameters and capabilities.

The Sine environment then uses those transfer function to compute initial drive amplitudes and phases.

14.5Test Predictions for the MIMO Sine Environment

Once the system identification is available, the Sine environment can use that information in conjunction with the defined control law to compute initial excitation signals and resulting control responses. The control responses are computed by convolving the excitation signals with the impulse response (the inverse fast Fourier transform of the transfer functions). The specified tracking filters are then used to extract the sine tone amplitudes and phases from the signals. These predictions appear on the Test Predictions tab, shown conceptually in Figure 14.3.

Prediction UI used by the MIMO Sine environment.

Figure 14.3:Prediction UI used by the MIMO Sine environment.

The prediction UI allows the user to inspect:

The predicted response amplitudes and phases are compared directly against the specification. If time histories are compared, they are compared against time histories synthesized from the specification. Channels that are predicted to violate warning or abort limits are highlighted in the response error matrix.

14.5.1Excitation Prediction

The excitation display on the top portion of the window contains the following displays:

Users can plot drive signals, amplitude, and phase quantities over time or frequency. Certain quantities must be plotted on a tone-by-tone basis, while others can be plotted as a superposition of all tones.

14.5.2Response Prediction

Similarly, the response prediction display on the bottom of the window contains the following displays:

14.6Running the MIMO Sine Environment

The Run Test tab is used to actually run the environment after system identification and prediction are complete.

A typical run page is shown in Figure 14.4.

Run GUI used by the MIMO Sine environment.

Figure 14.4:Run GUI used by the MIMO Sine environment.

The run page allows the user to:

14.6.1Drive Updates

The Drive Updates shows the updates that the controller is applying to the open-loop drive signal based on error correction. It is plotted as a complex amplitude, so the phase is encoded as the angle to the xx-axis, and the amplitude is the distance from the origin. A history of previous drive updates is plotted with lightening colors, so users can identify if the drive updates are diverging.

14.6.2Environment Control

The Environment Control allows the user to define how the Sine environment will be run. It contains the following controls:

14.6.3Response Amplitude and Phase

The Response Amplitudes and Response Phases portions of the controller show the amplitude and phase extracted from the channel (column) and sine tone (row) selected in the Response Signal Selector.

14.6.4Response and Drive Channel Selection

The Response Signal Selector and Drive Signal Selector allow the user to visualize the response amplitude error (in dB) and the drive update amplitude (in dB) respectively. Clicking on a cell in either table will display the data for that sine tone (row) and control or drive channel (column).

14.6.5Individual Tone/Channel Displays

The Data Display portion of the Run Test tab allows the user to create new windows showing specific sine tones and channels.

An example window is shown in Figure 14.5.

A window displaying the amplitude and phase response compared to the specification.

Figure 14.5:A window displaying the amplitude and phase response compared to the specification.

14.7Tracking Amplitude and Phase

A central task in the MIMO Sine environment is the extraction of the instantaneous amplitude and instantaneous phase of the measured response for each control channel and each active sine tone. These extracted quantities are used for:

Because the Sine environment is deterministic rather than stationary, Rattlesnake does not primarily control to averaged random quantities such as CPSD matrices during a run. Instead, it estimates the response of each active tone directly in the time domain and expresses that response as an amplitude and phase evolving over time.

At a conceptual level, each tracked response component is modeled as

y(t)≈A(t)cos⁡ ⁣(θ(t)+ϕ(t))y(t) \approx A(t)\cos\!\bigl(\theta(t)+\phi(t)\bigr)
(14.5)

where:

The primary filtering problem is therefore to recover A(t)A(t) and ϕ(t)\phi(t) from the measured response y(t)y(t), given the known excitation or response reference phase history.

Rattlesnake currently supports two approaches for this task:

  1. a Digital Tracking Filter (DTF), and

  2. an Overlapped Vold-Kalman (VK) Filter.

The DTF is computationally faster and relatively intuitive. The VK filter is more mathematically sophisticated and can provide better discrimination between nearby or crossing sine tones, at the cost of additional computational expense.

14.7.1Digital Tracking Filter

The digital tracking filter operates by demodulating the signal against a known sinusoidal reference and then low-pass filtering the demodulated components.

Suppose the measured signal contains a component of interest whose known instantaneous argument is θ(t)\theta(t). Then the measured signal may be projected onto in-phase and quadrature reference signals:

y0(t)=y(t)cos⁡θ(t)y_0(t) = y(t)\cos\theta(t)
(14.6)
y90(t)=−y(t)sin⁡θ(t)y_{90}(t) = -y(t)\sin\theta(t)
(14.7)

These two signals contain a slowly varying baseband component associated with the desired sine tone, along with higher-frequency content that is rejected by low-pass filtering.

After low-pass filtering, Rattlesnake obtains filtered in-phase and quadrature estimates, which may be denoted x0(t)x_0(t) and x90(t)x_{90}(t). From these, the instantaneous amplitude and phase are reconstructed as

A(t)=2x0(t)2+x90(t)2A(t) = 2\sqrt{x_0(t)^2 + x_{90}(t)^2}
(14.8)
ϕ(t)=atan2⁡ ⁣(x90(t),x0(t))\phi(t) = \operatorname{atan2}\!\bigl(x_{90}(t), x_0(t)\bigr)
(14.9)

This is exactly the logic implemented in the digital tracking filter generator used by the Sine environment.

The DTF therefore behaves like a classical synchronous detector followed by a low-pass filter:

14.7.1.1Tracking Filter Cutoff

The Tracking Filter Cutoff parameter sets the low-pass cutoff as a fraction of the instantaneous tone frequency.

A higher cutoff:

A lower cutoff:

14.7.1.2Tracking Filter Order

The Tracking Filter Order parameter sets the order of the Butterworth low-pass filter used after demodulation. Higher-order filters roll off more sharply, which can improve rejection of unwanted content, but may introduce a slower or more oscillatory transient response.

14.7.1.3Digital Tracking Filter Summary

The DTF is generally a good choice when:

14.7.2Vold-Kalman Filter

The Vold-Kalman Filter 1,2 is a more advanced order tracking method, and includes the method to track multiple sine tones simultaneously. The basic formulation of the VK filter is as follows.

For a given discretely-sampled signal y(n)y(n) containing primarily a harmonic component, the equation for y(n)y(n) can be written as

y(n)=x(n)exp⁡(jΘ(n))+η(n)y(n) = x(n)\exp(j\Theta(n)) + \eta(n)
(14.10)

Here x(n)x(n) is the instantaneous complex amplitude of the sine tone, and Θ(n)\Theta(n) is the instantaneous phase or argument of the sine tone. η(n)\eta(n) represents an error term. To use the VK filter, we therefore need to know the argument Θ(n)\Theta(n) of the sine tone over time; because we can analytically construct the sine sweep frequencies from the breakpoint table in the specification, we can solve for this quantity as the integral of frequency over time. The solution of the VK filter is therefore concerned with the identification of the complex amplitude x(n)x(n) at each sample for the sine tone, as this will contain instantaneous amplitude and phase of that sine tone, which will be compared to the specified amplitudes and phases to judge test accuracy and perform closed-loop control.

In addition to the data equation (14.10), the VK filter also contains structural equations which describes the mathematical characteristics of the sine tone to be extracted. The goal in this case is that the complex envelope varies slowly over time, at least in comparison to the frequency of the sine tones themselves. Therefore, we can write the structural equation in terms of finite difference operations on the envelope. For first-order, second-order, and third-order filters, these equations are:

∇x(n)=x(n)−x(n+1)=ε(n)\nabla x(n) = x(n) - x(n+1) = \varepsilon(n)
(14.11)
∇2x(n)=x(n)−2x(n+1)+x(n+2)=ε(n)\nabla^2 x(n) = x(n) - 2x(n+1) + x(n+2) = \varepsilon(n)
(14.12)
∇3x(n)=x(n)−3x(n+1)+3x(n+2)−x(n+3)=ε(n)\nabla^3 x(n) = x(n) - 3x(n+1) + 3x(n+2) - x(n+3) = \varepsilon(n)
(14.13)

In these equations the value ε(n)\varepsilon(n) represents a small term that allows for the slow modulation of the envelope.

All discrete variables from these equations can be arranged into vectors of length NN where NN is the total number of samples.

y=[y(1),y(2),…,y(N)]T\mathbf{y} = \left[y(1), y(2), \dots, y(N)\right]^T
(14.14)
x=[x(1),x(2),…,x(N)]T\mathbf{x} = \left[x(1), x(2), \dots, x(N)\right]^T
(14.15)
η=[η(1),η(2),…,η(N)]T\boldsymbol{\eta} = \left[\eta(1), \eta(2), \dots, \eta(N)\right]^T
(14.16)
ε=[ε(1),ε(2),…,ε(N)]T\boldsymbol{\varepsilon} = \left[\varepsilon(1), \varepsilon(2), \dots, \varepsilon(N)\right]^T
(14.17)

We can then set up matrix forms for the structural equations (14.18) and data equations (14.19).

Ax=ε\mathbf{A}\mathbf{x} = \boldsymbol{\varepsilon}
(14.18)
y−Cx=η\mathbf{y}-\mathbf{C}\mathbf{x} = \boldsymbol{\eta}
(14.19)

where the coefficient matrix A\mathbf{A} is

A=[1−10⋯0001−1⋯00001⋯00⋮⋮⋮⋱⋮⋮000⋯1−1]\mathbf{A} = \begin{bmatrix} 1 & -1 & 0 & \cdots & 0 & 0 \\ 0 & 1 & -1 & \cdots & 0 & 0 \\ 0 & 0 & 1 & \cdots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots\\ 0 & 0 & 0 & \cdots & 1 & -1 \end{bmatrix}
(14.20)

for the first-order filter or

A=[1−21⋯00001−2⋯000001⋯000⋮⋮⋮⋱⋮⋮⋮000⋯1−21]\mathbf{A} = \begin{bmatrix} 1 & -2 & 1 & \cdots & 0 & 0 & 0 \\ 0 & 1 & -2 & \cdots & 0 & 0 & 0 \\ 0 & 0 & 1 & \cdots & 0 & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots\\ 0 & 0 & 0 & \cdots & 1 & -2 & 1 \end{bmatrix}
(14.21)

for the second-order filter, etc. And, the coefficient matrix C\mathbf{C} is a diagonal matrix consisting of the complex argument phasors at each sample in time.

C=[exp⁡(jΘ(1))00⋯00exp⁡(jΘ(2))0⋯000exp⁡(jΘ(3))⋯0⋮⋮⋮⋱⋮000⋯exp⁡(jΘ(N))]\mathbf{C} = \begin{bmatrix} \exp(j\Theta(1)) & 0 & 0 & \cdots & 0 \\ 0 & \exp(j\Theta(2)) & 0 & \cdots & 0 \\ 0 & 0 & \exp(j\Theta(3)) & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & \exp(j\Theta(N)) \\ \end{bmatrix}
(14.22)

The goal is to minimize the values ε\boldsymbol{\varepsilon} and η\boldsymbol{\eta}, so we can represent the vectors as a scalar product.

εTε=xTATAx\boldsymbol{\varepsilon}^T\boldsymbol{\varepsilon} = \mathbf{x}^T\mathbf{A}^T\mathbf{A}\mathbf{x}
(14.23)
ηHη=(yT−xHCH)(y−xC)\boldsymbol{\eta}^H\boldsymbol{\eta} = \left(\mathbf{y}^T - \mathbf{x}^H\mathbf{C}^H\right)\left(\mathbf{y} - \mathbf{x}\mathbf{C}\right)
(14.24)

The weighted sum of these parameters form the loss function

J=r2εTε+ηHηJ = r^2\boldsymbol{\varepsilon}^T\boldsymbol{\varepsilon} + \boldsymbol{\eta}^H\boldsymbol{\eta}
(14.25)

where rr is a weighting parameter that determines how heavily the structural equations are weighted compared to the data equations.

The derivative of this function with respect to xx set to zero gives the minimum of this function.

∂J∂x=2r2ATAx+2(x−CHy)=0\frac{\partial J}{\partial x} = 2r^2\mathbf{A}^T\mathbf{A}\mathbf{x} + 2\left(\mathbf{x}-\mathbf{C}^H\mathbf{y}\right) = \mathbf{0}
(14.26)

The solution is then

x=(r2ATA+E)−1CHy\mathbf{x} = \left(r^2\mathbf{A}^T\mathbf{A}+\mathbf{E}\right)^{-1}\mathbf{C}^H\mathbf{y}
(14.27)

The above equations are satisfactory for a single tone. However, when multiple tones are simultaneously present in a given signal, we can solve for all tones simultaneously. In this case, the equation for y(n)y(n) becomes

y(n)=∑k=1Pxk(n)exp⁡(jΘk(n))+η(n)y(n) = \sum_{k=1}^{P} x_k(n)\exp(j\Theta_k(n)) + \eta(n)
(14.28)

where now there is a separate complex amplitude xk(n)x_k(n) and argument Θk(n)\Theta_k(n) for each of the PP sine tones in the signal.

Following the same logic as above, the function to minimize becomes

J=∑k=1Pr2εkTεk+ηHη=∑k=1Pr2xkHATAxk+(yT−∑k=1PxkHCkH)(y−∑k=1PCkxk)J = \sum_{k=1}^{P} r^2{\boldsymbol{\varepsilon}_k}^T\boldsymbol{\varepsilon}_k + \boldsymbol{\eta}^H\boldsymbol{\eta} \\ = \sum_{k=1}^{P} r^2{\mathbf{x}_k}^H\mathbf{A}^T\mathbf{A}\mathbf{x}_k + \left( \mathbf{y}^T - \sum_{k=1}^{P} {\mathbf{x}_k}^H{\mathbf{C}_k}^H \right)\left( \mathbf{y} - \sum_{k=1}^{P} {\mathbf{C}_k}{\mathbf{x}_k} \right)
(14.29)

The solution is then the xk\mathbf{x}_k for each of the PP sine tones that minimize this equation. Setting the derivatives to zero gives:

∂J∂xiH=Bixi+CiH+∑k=1,k≠iPCkxk−CiHy=0\frac{\partial J}{\partial \mathbf{x}_i^H} = \mathbf{B}_i \mathbf{x}_i + \mathbf{C}_i^H + \sum_{k=1, k\ne i}^P \mathbf{C}_k\mathbf{x}_k - \mathbf{C}_i^H\mathbf{y}=\mathbf{0}
(14.30)

where Bi=r2ATA+I\mathbf{B}_i = r^2\mathbf{A}^T\mathbf{A} + \mathbf{I} and I\mathbf{I} is the identity matrix. We also note that CiHCi=I\mathbf{C}_i^H\mathbf{C}_i=\mathbf{I} due to the phases of the Hermetian diagonal phasor matrices cancelling out when multiplied.

We can assemble this into a large system of equations of the form Bx=b\mathbf{B} \mathbf{x} = \mathbf{b} with

B=[B1C1HC2⋯C1HCPC2HC1B2⋯C2HCP⋮⋮⋱⋮CPHC1CPHC2⋯BP]\mathbf{B} = \begin{bmatrix} \mathbf{B}_1 & \mathbf{C}_1^H\mathbf{C}_2 & \cdots & \mathbf{C}_1^H\mathbf{C}_P \\ \mathbf{C}_2^H\mathbf{C}_1 & \mathbf{B}_2 & \cdots & \mathbf{C}_2^H\mathbf{C}_P \\ \vdots & \vdots & \ddots & \vdots \\ \mathbf{C}_P^H\mathbf{C}_1 & \mathbf{C}_P^H\mathbf{C}_2 & \cdots & \mathbf{B}_P \\ \end{bmatrix}
(14.31)
x=[x1x2⋮xP]\mathbf{x} = \begin{bmatrix} \mathbf{x}_1 \\ \mathbf{x}_2 \\ \vdots \\ \mathbf{x}_P \end{bmatrix}
(14.32)
b=[C1HyC2Hy⋮CPHy]\mathbf{b} = \begin{bmatrix} \mathbf{C}_1^H\mathbf{y} \\ \mathbf{C}_2^H\mathbf{y} \\ \vdots \\ \mathbf{C}_P^H\mathbf{y} \\ \end{bmatrix}
(14.33)

As the B\mathbf{B} matrix has a sparse, banded form, the equation is solved using the sparse solvers in SciPy’s linalg package. Once the envelope x\mathbf{x} is found, it can be split out into the individual complex amplitudes for each sine tone. From these complex amplitudes, the amplitude and phase of the individual sine tones can be derived.

The complex envelope contains the information of interest:

A(n)=∣x(n)∣A(n) = |x(n)|
(14.34)
ϕ(n)=∠x(n)\phi(n) = \angle x(n)
(14.35)

There are a few parameters that can be selected when setting up the VK filter.

14.7.2.1Filter order in the Vold-Kalman method

The Vold-Kalman Filter Order parameter controls the order of the finite-difference smoothness constraint.

Conceptually:

Higher order generally increases filter selectivity and mathematical smoothness of the solution at the expense of additional computation.

14.7.2.2Bandwidth in the Vold-Kalman method

The Vold-Kalman Filter Bandwidth controls how narrowly the filter tracks the desired tone.

A smaller bandwidth gives a more selective filter and improves rejection of neighboring tones and broadband noise, but it can make the response slower and may distort rapidly changing amplitude modulation. A larger bandwidth allows faster tracking of amplitude or phase changes but reduces selectivity and can admit contamination from other content. The bandwidth parameter is closely related to the weighting parameter rr in (14.25), as it controls the tradeoff between smoothness of the solution and matching of the data.

14.7.2.3Multiple simultaneous tones

One of the major advantages of the Vold-Kalman formulation is that multiple sinusoidal components can be solved simultaneously. This is especially important when multiple sine tones exist in the signal and is absolutely necessary if the sine tones cross in frequency.

In these cases, a simple tracking filter may struggle to separate the contributions of the different tones, while the VK formulation can solve for them together as a coupled estimation problem. This is one of the main reasons the Sine environment offers the VK filter as an alternative to the DTF.

14.7.2.4Online and Overlapped Vold-Kalman Filtering

The classical Vold-Kalman filter is most naturally formulated over a full signal block. However, Rattlesnake must use it in an online control environment, where the signal is still being acquired and the control law must update in a timely way.

To make this possible, Rattlesnake implements the VK filter in an overlapped blockwise form.

A blockwise VK solution can exhibit startup and ending transients on each analysis block. If consecutive blocks were simply stitched together end-to-end, these edge effects would become visible in the extracted amplitudes and phases.

To mitigate this, Rattlesnake:

  1. divides the signal into overlapping blocks,

  2. solves the VK problem on each block,

  3. rejects the responses in the overlapped regions, keeping only the “center” portion of each block where the transients are minimal

The Vold-Kalman Filter Block Size sets the number of samples used in each VK solve.

A larger block size generally gives the filter more information, which can remove the effects of startup and ending transients. However, this increases computational cost and decreases controller responsiveness, as the control must acquire more data prior to making a control decision.

The Vold-Kalman Filter Block Overlap sets the overlap fraction between consecutive blocks.

This overlap is used to reduce discontinuities between blocks and to suppress the effect of the VK startup/end transients. In practice, too little overlap can leave visible stitching artifacts due to the startup and ending transients. Increasing overlap improves continuity between blocks; however, it increases redundant computation.

14.7.3Practical Filter Tradeoffs

From a practical standpoint, the DTF and VK filters represent a tradeoff between:

Digital Tracking Filter is often preferable when:

Vold-Kalman Filter is often preferable when:

The DTF is relatively cheap computationally. The VK filter can become expensive, especially because in the Sine environment the filter is applied separately for each control channel. Thus, the total computational burden scales with:

For that reason, a filter configuration that is acceptable for a single-channel example may become too expensive for a large multi-channel test.

This is why the filter explorer described in Section 14.7.4 is important: it gives the user the ability to investigate not only the quality of the extracted amplitude and phase, but also the practical computational behavior of the chosen filter.

14.7.4Using the Filter Explorer

The Filter Explorer is accessed by clicking the Explore... button on the Environment Definition tab. This brings up the Filter Explorer dialog box, which is shown in Figure 14.6. The Filter Explorer dialog will take the defined specification and generate a time signal that exactly matches the amplitudes and phases defined for each sine tone in the specification for the selected channel. Noise can also be added to this perfect specification realization to represent realistic test conditions that the filter should reject. The user can then select from filter parameters mirroring those found on the Environment Definition tab described in Section 14.3.3. The filter explorer will then attempt to filter the known signal from the specification, extract amplitude and phase information, and present those data in comparison to the known data from the specification. This allows the user to investigate filter choices prior to running a test.

The Filter Explorer dialog box that allows users to understand the effect of their tracking filter choices.

Figure 14.6:The Filter Explorer dialog box that allows users to understand the effect of their tracking filter choices.

Descriptions of the filter selection parameters are found in Section 14.3.3. Additional widgets and displays on the Filter Explorer dialog box include:

14.8Output NetCDF File Structure

Like the other environments in Rattlesnake, the MIMO Sine environment stores its metadata in a netCDF group whose name matches the environment name. This group contains the parameters needed to reconstruct the environment definition, including the sine specifications, the control law configuration, the tracking filter settings, the shared system identification settings, and any transformation matrices.

The root netCDF dataset also contains the global hardware metadata and channel table information described is Section 3.7. The material in this section focuses only on the additional data stored inside the Sine environment’s group.

Because the Sine environment derives from the shared system-identification environment infrastructure, its netCDF group contains both sine-specific fields and shared system-identification fields.

14.8.1NetCDF Dimensions

The Sine environment creates the following dimensions in its netCDF group.

14.8.2NetCDF Attributes

The following attributes are stored directly on the Sine environment’s netCDF group.

14.8.3NetCDF Variables

The following variables are stored directly on the Sine environment’s netCDF group.

14.8.4The specifications Group

A subgroup named specifications is created inside the environment group. Within this group, one subgroup is created for each sine tone in the environment. The subgroup name is the name of the tone as shown in the UI.

For each tone subgroup, the following dimensions are defined:

For each sine tone subgroup, the following attribute is stored:

Each tone subgroup also stores the following variables:

For the warning and abort limits, the meanings of the dimensions are:

14.8.5Interpretation of Stored Specification Data

The Sine environment stores the specification in the netCDF4 file in the same logical form as used internally during the run:

This is slightly different from the external .mat and .npz sine specification files used when loading individual tones, which are designed to be convenient for exchange with MATLAB or NumPy workflows. In particular:

14.9Saving Control Data

The Sine environment allows saving current control data from the run tab. This occurs when the user clicks the Save Control Data button.

Data is stored to a NumPy Archive (*.npz) file.

This typically includes:

This data is useful for:

14.10Using Transformation Matrices

Transformation matrices in the Sine environment behave identically to the the Random Vibration environment. See Section 12.8 for more information.

14.11Writing a Custom Sine Control Law

Unlike the Random Vibration environment, where relatively simple function-based control laws can often be sufficient, the MIMO Sine environment is inherently more stateful. The control logic must manage:

For this reason, the Sine environment is built around class-based custom control laws rather than one-shot function calls. A custom Sine control law must therefore be implemented as a Python class with the methods described below. The easiest way to understand the intent of these methods is to compare them to the current default control law, implemented in rattlesnake/environment/sine_sys_id_utilities.py as DefaultSineControlLaw.

A MIMO sine test is not simply a static inversion problem at one frequency. The environment must:

The control law therefore needs to maintain internal state such as:

This is why the Sine environment expects a control-law class rather than a simple stateless function.

14.11.1Required Control-Law Lifecycle Methods

A custom Sine control law class is expected to expose the following methods:

The environment calls these methods at different stages of the test.

14.11.1.1__init__(...)

The constructor is called when the environment metadata are initialized and the control law is created. It should store all persistent information required by the control law, including:

The constructor is also a reasonable place to precompute signal representations derived from the specification, such as:

The current default control law does exactly this. In particular, it uses SineSpecification.create_combined_signals(...) to build:

These are stored as persistent arrays so they do not need to be recomputed repeatedly during the run.

14.11.1.2system_id_update(...)

This method is called when the system identification data have been obtained or updated. Its purpose is to convert the system-ID results into a usable open-loop drive estimate.

A custom implementation should use the transfer functions to compute some initial estimate of the drive amplitudes and phases required to reproduce the specified response.

In the default implementation, this method performs the following major steps:

  1. stores the transfer function matrix and related system-ID quantities,

  2. interpolates the FRF pseudoinverse onto the instantaneous frequencies of the specification,

  3. computes the desired complex response

    x(t)=ax(t)eiϕx(t)\mathbf{x}(t) = \mathbf{a}_x(t)e^{i\boldsymbol{\phi}_x(t)}
    (14.36)

  4. maps that desired response into a complex drive estimate

    v(t)=Hxv+(t)x(t)\mathbf{v}(t) = \mathbf{H}_{xv}^{+}(t)\mathbf{x}(t)
    (14.37)

  5. extracts the resulting drive amplitudes and phases,

  6. reconstructs the corresponding time-domain preshaped drive signals.

This stage is where the default control law solves for the open-loop drive signals that would produce the best estimate of the desired response, assuming the transfer functions are measured exactly.

Any replacement control law must do something equivalent, even if the exact mathematical method differs.

14.11.1.3initialize_control(...)

This method is called at the start of a particular run. It is responsible for initializing the runtime state needed for online control.

Typical responsibilities include:

The current default control law uses this method to:

This method returns the first excitation signal block that will be sent to the output process.

14.11.1.4update_control(...)

This is the key online feedback method. It is called when a new portion of the measured response has been analyzed and reduced to:

A custom control law should compare the achieved response to the target response over the current block and update its internal drive-correction state.

The default implementation does this by:

  1. storing the achieved response data,

  2. identifying the corresponding time block within the full specification,

  3. reconstructing the target complex response over that block,

  4. reconstructing the achieved complex response over that block,

  5. computing the complex error,

  6. projecting that error back through the transfer function pseudoinverse,

  7. scaling the result by a convergence factor,

  8. accumulating the resulting correction into self.control_drive_correction.

Conceptually, the default implementation is applying a correction of the form

Δv(t)∝Hxv+(t)(xtarget(t)−xachieved(t))\Delta \mathbf{v}(t) \propto \mathbf{H}_{xv}^{+}(t)\left(\mathbf{x}_{\mathrm{target}}(t)-\mathbf{x}_{\mathrm{achieved}}(t)\right)
(14.38)

and accumulating this correction across blocks.

Any custom control law that performs closed-loop sine control must implement an analogous update rule, even if the exact form of the correction differs.

14.11.1.5generate_signal(...)

This method is called whenever the signal generation process needs another block of output samples.

It should:

In the default implementation, this method forms a complex excitation signal of the form

vblock(t)=av(t)eiϕv(t)+Δv\mathbf{v}_{\mathrm{block}}(t) = \mathbf{a}_v(t)e^{i\boldsymbol{\phi}_v(t)} + \Delta\mathbf{v}
(14.39)

and then converts it into a real time signal using cosine reconstruction.

It also applies an optional drive-voltage limit by clipping the magnitude of the complex excitation vector.

This is the method that turns the control law’s internal state into actual samples to be output by the hardware.

14.11.1.6finalize_control(...)

This method is called after the run has completed. It allows the control law to perform any final processing or to prepare updated preshaped drive information for later use.

In the current default implementation, this method primarily returns the current drive signals and associated metadata:

A more advanced custom implementation could use this stage to:

14.11.2What the Default Control Law Does

The current DefaultSineControlLaw in sine_sys_id_utilities.py can be understood as having the following stages.

14.11.2.1Stage 1: Parse and store the specification

In __init__, the default control law:

This establishes the full target trajectory that the environment will attempt to realize.

14.11.2.2Stage 2: Ingest the system identification results

In system_id_update, the default law:

This stage generates the initial open-loop estimate of the excitation.

14.11.2.3Stage 3: Build preshaped drive signals

Still within system_id_update, the default implementation reconstructs the corresponding time-domain preshaped drives by applying the amplitude and phase trajectories to the tone arguments.

These preshaped drives are what the environment uses as the baseline excitation before any online correction is applied. This is the signal used to develop the test predictions shown on the Test Predictions tab.

14.11.2.4Stage 4: Initialize runtime control state

In initialize_control, the default control law:

This transitions the control law from an open-loop prediction object into a live online controller.

14.11.2.5Stage 5: Update based on measured response

In update_control, the default implementation:

This is the core online feedback stage.

14.11.2.6Stage 6: Generate the next output block

In generate_signal, the default control law:

This continues until all requested signal data have been generated.

14.11.2.7Stage 7: Finalize

In finalize_control, the default law returns the final drive predictions and associated arrays so they can be reused or displayed after the run.

14.11.3What a Custom Sine Control Law Must Do

A replacement control law should preserve the same broad responsibilities, even if the mathematics differ.

At minimum, a custom implementation should be able to:

  1. accept and store specification information,

  2. accept updated system identification results,

  3. compute an initial open-loop drive estimate,

  4. initialize per-run state,

  5. compare achieved response to target response during the run,

  6. compute and accumulate drive corrections,

  7. generate output signal blocks,

  8. optionally finalize or update future drive predictions after the run.

In other words, even if a custom control law uses a completely different control strategy, it still needs to participate correctly in the same lifecycle as the default implementation.

14.11.4Practical Guidance for Implementing a Replacement

14.11.4.1Preserve state explicitly

The Sine control problem is inherently stateful. A custom class should explicitly track:

14.11.4.2Be careful about coordinate systems

The control law may be operating on:

A custom implementation should ensure that the transfer functions, specification, and drive outputs are all interpreted in compatible coordinates.

14.11.4.3Respect blockwise execution

The environment acquires and outputs data in finite blocks. A custom control law must therefore be designed to work incrementally, not just as a single full-signal solve.

14.11.4.4Use tracked data, not raw data

By the time update_control(...) is called, the environment has already extracted:

Most custom control laws should work with those reduced quantities rather than trying to reinterpret the raw time histories from scratch.

14.11.4.5Think carefully about output limits

The default implementation optionally clips the magnitude of the complex excitation to a maximum drive voltage. A replacement control law should decide explicitly:

14.11.5Summary

A custom Sine control law in Rattlesnake is best understood as a stateful object managing the full sine-test lifecycle, not merely as a function that computes one correction.

The default implementation provides a working reference for:

When implementing a replacement control law, the most important thing is not to duplicate the exact mathematics of the default implementation, but rather to preserve the required lifecycle:

Any custom control law that satisfies that lifecycle correctly can be integrated into the Sine environment.

References
  1. Vold, H., & Leuridan, J. (1995). High Resolution Order Tracking at Extreme Slew Rates Using Kalman Tracking Filters. Shock and Vibration, 2(6), 919402. https://doi.org/10.3233/SAV-1995-2609
  2. Blough, J. R. (2007). Understanding the Kalman/Vold-Kalman Order Tracking Filters’ Formulation and Behavior. SAE Noise and Vibration Conference and Exhibition, 07NVC-254.
Rattlesnake Vibration Controller
Multiple Input/Multiple Output Transient Control
Rattlesnake Vibration Controller
Multiple Input/Multiple Output Shock Control