Experimental Application of Time-Domain Transmissibility Identification to Fault Detection and Localization in Acoustic Systems

The study of structural vibration focuses on the displacement, velocity, and acceleration response of a structure to force and torque inputs. This response can be studied in terms of either transfer functions or differential equation models. In the former case, the frequency response of the structure can be estimated by computing the ratio of the transforms of the forcing and response signals. These procedures are based on the assumption that the input and output signals are stationary, and thus, initial conditions and transient effects are either assumed to be absent or are ignored. In the case of differential equation models, time-domain methods can be used to estimate the model parameters assuming that the input is sufficiently persistent. In this case, initial conditions and transient effects can enhance rather than degrade the accuracy of the estimated model. The effect of nonzero initial conditions on time-domain and frequency-domain identification is shown in Ref. [1].

An alternative approach to force-to-motion models is the class of models known as transmissibilities. A transmissibility is a relationship between identical variables, such as force-to-force, position-to-position, velocity-to-velocity, and acceleration-to-acceleration [2,3]. The input and output of a transmissibility, referred to as the pseudo input and pseudo output, respectively, are thus outputs of the underlying system. Although force-to-motion transfer functions can be identified using either time-domain or frequency-domain methods, transmissibility estimates are traditionally obtained using only frequency-domain methods, and this remains an active area of research [4–9].

Interestingly, time-domain methods have not been traditionally used to identify transmissibilities. The reason for this is partly due to the fact that the meaning of a transmissibility in the time domain is suspect. In particular, the input (forcing) to the underlying system plays no visible role in the transmissibility, and the states of the transmissibility have no physical interpretation. Nevertheless, for applications in which response data have a significant nonstationary component, it is advantageous to use time-domain methods to estimate transmissibility functions. A step in this direction was taken in Refs. [1,10–13], where time-domain transmissibility models were obtained. These models, called pseudo transfer functions, are formulated in terms of the differential operator $p=d/dt$ rather than the Laplace s. This framework requires some care due to the fact that $p$ is not a complex number. In addition, unlike input–output transfer functions in the Laplace domain, the differentiation operator accounts for the free response in the time domain [12], but requires special attention due to issues of causality, stability, and order. In particular, a time-domain transmissibility model may be noncausal (in either continuous time or discrete time), where its numerator and denominator dynamics are the zero dynamics of the underlying transfer functions from the forcing to the respective motion sensors.

Various fault-detection techniques have been introduced in the literature [14–17]. In some cases, health monitoring can be performed by exciting the system in a controlled manner, using a plant model and observer to predict the response, and comparing the measured response to the prediction [18–20]. This approach, known as active fault detection, is based on residual generation. In contrast, passive fault detection detects faults by analyzing the sensor signals alone and searching for anomalies [21–23].

In the present paper, we focus on a technique for fault detection based on time-domain identification of transmissibilities, which is neither active nor passive as defined earlier. Instead, this technique takes advantage of freely available and unknown external (ambient) excitation to identify a sensor-to-sensor model (i.e., a transmissibility operator), which is independent of the excitation signal. In the presence of subsequent unknown external excitation, the identified transmissibility is used to compute sensor-to-sensor residuals, which are used to detect and diagnose faults in sensors and system dynamics. The sensor-to-sensor residual is the discrepancy between the predicted sensor output (based on the transmissibility) and the actual measurements.

The ability to take advantage of unknown external excitation along with the fact that the transmissibility is independent of the excitation gives the method flexibility in practice by removing the need for a known or controlled excitation. For example, freely available ambient noise can play a useful role in transmissibility identification. Most importantly, the identified transmissibility is independent of the excitation; this means that the transmissibility identified using one data set can be used for fault detection with a different data set; for both data sets, the external excitation can be completely unknown. This feature is the key benefit of this approach relative to fault-detection methods that require known external excitation. Transmissibilities are applied to fault detection and localization in Refs. [24–31]. Note that comparing transmissibility-based fault detection to other fault-detection methods such as Refs. [18–20], [32], and [33], which require knowledge of a mathematical model of the underlying system or the excitation signal, is meaningless.

Since transmissibility operators may be noncausal, unstable, and of unknown order, we consider a class of models that can approximate transmissibility operators with these properties. This class of models consists of noncausal finite impulse response (FIR) models based on a truncated Laurent expansion [34]. These models are shown to approximate the Laurent expansion inside the annulus between the asymptotically stable pole of largest modulus and the unstable pole of smallest modulus. By delaying the measured pseudo output relative to the measured pseudo input, the identified finite impulse response model is a noncausal approximation of the transmissibility operator. The causal (backward-shift) part of the Laurent expansion is asymptotically stable since all of its poles are zero, while the noncausal (forward-shift) part of the Laurent expansion captures the unstable and noncausal components of the transmissibility operator. Subspace identification is also used with the measured pseudo input and the delayed measured pseudo output to identify a noncausal state space model.

Since the transmissibility operator has the form of an input–output time-series model, methods developed for ARMAX models are candidates for this application. Some caution is needed, however, since the pseudo input and pseudo output of the transmissibility are both noisy, leading to an errors-in-variables problem with noise on both the input and output signals [35–37]. In addition, neither the pseudo input nor the pseudo output signals (irrespective of the noise) can be expected to be white, thus complicating the consistency analysis. Due to the difficulty of the errors-in-variables problem, we compare the one-step prediction error for the estimates of the transmissibility operator obtained using least squares (LS), prediction error methods (PEM), instrumental variables (IV), and subspace identification, and we compare the accuracy of the estimates based on the one-step prediction error.

Transmissibility operators were used in Ref. [38] for rate-gyro health monitoring in aircraft. Moreover, transmissibility operators were used in Ref. [39] for aircraft structural health monitoring. In addition, transmissibility operators were used in Ref. [40] for fault detection in an acoustic system consisting of a drum, two speakers, and four microphones. In the present paper, we extend the results in Ref. [40] by considering an experimental setup consisting of an acoustic system with three speakers and six microphones. Each speaker is an actuator, and each microphone is a sensor that measures the acoustic response at its location. Measurements from the six microphones are used to construct transmissibility operators, which in turn are used to detect and localize changes in the dynamics of the acoustic system or the microphones and by computing the resulting one-step residual.

The paper is organized as follows: In Sec. 2, we derive multi-input, multi-output transmissibility operators. In Sec. 3, we show that transmissibility operators can be noncausal, unstable, and of unknown order and thus noncausal FIR models can be used to approximate transmissibility operators with these properties. In Sec. 4, we use noncausal FIR models with least squares, prediction error methods, instrumental variables, and subspace methods to identify transmissibility operators. In Sec. 5, we use the identified transmissibility operators for fault detection and localization in acoustic systems. Finally, we give conclusions in Sec. 6.

where $u∈ℝm$⁠, p is the number of measurements, m is the number of pseudo inputs, and p − m is the number of pseudo outputs. The coefficient matrices have dimensions $A∈ℝn×n,B∈ℝn×m,Ci∈ℝm×n,Co∈ℝ(p−m)×n,Di∈ℝm×m,$ and $Do∈ℝ(p−m)×m$⁠.

Note that, if x₀ is zero, then $û(s)$ can be canceled in Eq. (7), and $ŷo(s)$ and $ŷi(s)$ are related by a transmissibility that is independent of the input. However, if x₀ is not zero, then $û(s)$ cannot be canceled in Eq. (7).

The transmissibility operator (11) is in continuous time. Henceforth, we assume that measurements of the output signals are obtained in discrete time, and we consider discrete-time transmissibility operators in the forward-shift operator $q$ [41].

Expression (11) shows that $T(p)$ contains information about the zeros of the system but not the poles. Therefore, a non-minimum-phase zero in the input channel of a transmissibility operator results in an unstable transmissibility operator. Moreover, if the output channel of a transmissibility operator has more zeros than the input channel, then the transmissibility operator is improper, that is, noncausal. Since this information may not be known in advance, we consider a class of models that can be used to approximate transmissibility operators that may be both unstable and improper.

Finite impulse response models can be used to describe linear asymptotically stable systems [42]. In addition, noncausal FIR models can be used to describe systems that are both noncausal and unstable [34]. A noncausal FIR model, which approximates the Laurent series of an unstable system, involves both positive and negative powers of the Z-transform variable z. The negative powers of z approximate the asymptotically stable part of the plant outside of a disk (that is, inside a punctured plane), whereas the positive powers of z approximate the unstable part of the plant inside a disk. Inside the common region, which is an annulus, the Laurent series represents a noncausal model, as evidenced by the positive powers of z [34]. Noncausal FIR models were also used in Refs. [43–45].

To identify transmissibility operators that are possibly improper, unstable, and of unknown order, we first use noncausal FIR models with LS, PEM [46], and IV [47].

We use the Matlab functions pem and iv4 to obtain PEM and IV estimates $θ̂r,d,ℓ$ of $θr,d$⁠, respectively.

where r has the same value as in Eq. (25) in order to perform a fair comparison between Eqs. (25) and (32).

Constructing a meaningful transmissibility operator requires knowledge of the number m of independent external excitation signals acting on the system. Since m may be unknown, we estimate m using the following procedure. Let $m̂∈{1,…,p−1}$ and $p̂∈{1,…,p−m̂}.$ We identify a transmissibility operator with $m̂$ pseudo inputs and $p̂$ pseudo outputs using the methods discussed earlier. For each identified transmissibility operator, we compute the residual using Eq. (24) or (31) and the norm of the residual using Eq. (25) or (32). The estimated number of independent external excitation signals is the value of $m̂$ at which a sharp drop occurs in the norm of the residual. If a sharp drop is not obvious, then the estimated number of external excitation signals is the smallest value of $m̂$ for which no sizable improvement is obtained for larger values of $m̂$⁠. Redundant sensors can then be removed or retained for possible benefit in terms of the accuracy of the identified transmissibility operators. This method will be illustrated in Sec. 5.2.

In order to investigate the ability of transmissibility operators to detect changes in the dynamics of an acoustic system, we consider the experimental setup shown in Fig. 1. The setup consists of an acoustic space with three speakers Spk1, Spk2, and Spk3 and six microphones $Mic1,…,Mic6$⁠. Each speaker is an actuator, and each microphone is a sensor that measures the acoustic response at its location. All actuator signals are generated using matlab and sent to the speakers through a data acquisition card. The sampling rate is 1000 Hz.

Let u₁, $u2,$ and u₃ be the measurements of the signals of the speakers Spk1, Spk2, and Spk3, respectively, and let $y1,…,y6$ be the measurements obtained by the sensors $Mic1,…,Mic6$⁠, respectively. We assume that data are available for $1≤k≤20,000$⁠.

To identify transmissibility operators, we delay the pseudo output d steps with respect to the pseudo input and perform the identification using the pseudo input and the delayed pseudo output. We use LS, PEM, and IV with a noncausal FIR model with r = 24, $d=25,$ and the first $ℓ=6000$ samples to obtain the identified transmissibilities. The least squares estimate is obtained using Eq. (22), where the PEM and IV estimates are obtained using the Matlab commands pem(data,'nb',(r + d + 1)* ones(1, i),'nc' , r + d + 1) and iv4(data,'nb',(r + d + 1)*ones(1,i)), respectively, where i is the number of pseudo inputs.

For subspace identification, we delay the pseudo output d steps with respect to the pseudo input, and we use the Matlab command n4sid(data,n) and the first $ℓ=6000$ samples to obtain the identified noncausal state space model of order n. Since an $nth$-order infinite impulse response model has $2n+2$ parameters and an $nth$-order FIR model has n parameters, we compare an infinite impulse response model of order n to an FIR model of order $2n+2$⁠. Since the order of the noncausal FIR model is $r+d+1,$ for subspace identification, we use the model order $n=(r+d−1)/2=24.$

Suppose that Spk1 is driven with a realization of band-limited white noise with a bandwidth of 500 Hz and Spk2 and Spk3 are not operating. Figure 2 shows the estimated Markov parameters $θ̂r,d,ℓ$ of $T(q−1,θr,d)$ from the pseudo input y₁ to the pseudo output y₆ obtained using least squares for two different runs of the experiment. Figure 3 shows y₆ and the computed one-step prediction $ŷ6=△T(q−1,θ̂r,d,ℓ)y1$ for $15,000≤k≤15,300,$ that is, for $t∈[15,15.3]$ s.

Next, suppose that Spk1 and Spk2 are driven with realizations of bandlimited white noise with a bandwidth of 500 Hz and Spk3 is not operating. Figure 4 shows the entries of the estimated Markov parameters $θ̂r,d,ℓ$ of $T(q−1,θr,d)$ from each pseudo input y₁ and y₃ to the pseudo output y₆ obtained using least squares for two different runs of the experiment. Figure 5 shows y₆ and the computed one-step prediction $ŷ6=△T(q−1,θ̂r,d,ℓ)[y1 y3]T$ for $15,000≤k≤15,300,$ that is, for $t∈[15,15.3]$ s.

Next, suppose that Spk1, Spk2, and Spk3 are driven with realizations of bandlimited white noise with a bandwidth of 500 Hz. Figure 6 shows the entries of the estimated Markov parameters $θ̂r,d,ℓ$ of $T(q−1,θr,d)$ from each pseudo input y₁, $y3,$ and y₅ to the pseudo output y₆ obtained using least squares for two different runs of the experiment. Figure 7 shows y₆ and the computed one-step prediction $ŷ6≜T(q−1,θ̂r,d,ℓ)[y1 y3 y5]T$ for $15,000≤k≤15,300,$ that is, for $t∈[15,15.3]$ s.

Suppose that the system shown in Fig. 1 is operating under healthy conditions, and suppose that Spk1 is driven with a realization of bandlimited white noise with a bandwidth of 500 Hz and Spk2 and Spk3 are not operating. Figure 9 shows the norm of the residual for each identified transmissibility operator whose pseudo input is $Yi=△[y1 ⋯ yi]T,$ for $i=1,…,5,$ and whose pseudo output is y₆ obtained using LS, PEM, IV, and subspace methods. Note from Fig. 9 that the improved residual obtained by using 2, 3, 4, or 5 pseudo inputs is not significant, which correctly suggests that there is one independent external excitation signal acting on the system.

Next, suppose that the system shown in Fig. 1 is operating under healthy conditions, and suppose that both Spk1 and Spk2 are driven with realizations of bandlimited white noise with a bandwidth of 500 Hz, and Spk3 is not operating. Figure 10 shows the norm of the residual for each identified transmissibility operator whose pseudo input is Y_i for $i=1,…,5,$ and whose pseudo output is y₆ obtained using LS, PEM, IV, and subspace methods. Figure 10 shows that the identified transmissibility operator with two pseudo inputs yields significantly lower residual than the identified transmissibility operators with one pseudo input. However, the improved residuals obtained by using 3, 4, or 5 pseudo inputs are not significant, which correctly suggests that the number of independent external excitation signals acting on the systems is two.

Next, suppose that the system shown in Fig. 1 is operating under healthy conditions, and suppose that Spk1, Spk2, and Spk3 are driven with realizations of bandlimited white noise with a bandwidth of 500 Hz. Figure 11 shows the norm of the residual for the identified transmissibility operator whose pseudo input is Y_i for $i=1,…,5,$ and whose pseudo output is y₆ obtained using LS, PEM, IV, and subspace methods. Figure 11 shows that the identified transmissibility operators with three pseudo inputs yield significantly lower residual than the identified transmissibility operators with one and two pseudo inputs. However, the improved residual obtained by using four or five pseudo input is not significant, which correctly suggests that the number of independent external excitation signals acting on the systems is three.

Figure 12 shows the norm of the residuals for the identified transmissibility operators obtained using least squares with $1,…,5$ pseudo inputs for the case of one, two, or three speakers driven with realizations of bandlimited white noise with a bandwidth of 500 Hz. Figure 12 shows that, in order to reduce the residuals, the number of pseudo inputs has to be greater than or equal to the number of independent inputs acting on the system.

To emulate change in the dynamics of the acoustic system, we remove a blanket that covers the acoustic space during the experiment. For structural health monitoring, we use the identified transmissibility $T(q−1,θ̂r,d,ℓ)$ whose pseudo input is $yi=[y1 y3 y5]T$ and whose pseudo output is $yo=y6$ and whose Markov parameters are obtained using least squares from run 1 shown in Fig. 6.

Suppose that Spk1, Spk2, and Spk3 are operating simultaneously and are driven with a band-limited white noise with a bandwidth of 500 Hz. At approximately t = 11 s, we remove the blanket that covers the acoustic system. Figure 13 shows $E(k|θ̂r,d,ℓ,w)$ for $T(q−1,θ̂r,d,ℓ),$ where w = 500 steps. Note from Fig. 13 the change in $E(k|θ̂r,d,ℓ,w)$ at approximately t = 11 due to the change in the dynamics of the system.

Next, suppose that Spk1, Spk2, and Spk3 are operating simultaneously and are driven with the harmonic signals (33)–(35). At approximately t = 11 s, we remove the blanket that covers the acoustic system. Figure 14 shows $E(k|θ̂r,d,ℓ,w)$ for $T(q−1,θ̂r,d,ℓ),$ where w = 500 steps. Note from Fig. 14 the change in $E(k|θ̂r,d,ℓ,w)$ at approximately t = 11 s due to the change in the dynamics of the system.

For sensor health monitoring, the level of the residual for a transmissibility changes as one of the sensors that are used to construct the transmissibility becomes faulty. We introduce faults in the microphones by modifying the low-pass analog filter that the microphone signal passes through before going to the dSPACE system.

Suppose that Spk1, Spk2, and Spk3 are operating simultaneously and are driven with a band-limited white noise with a bandwidth of 500 Hz. Let $T(q−1,θ̂r,d,ℓ),$ be the identified transmissibility whose pseudo input is $yi=[y1 y2 y3]T$ and whose pseudo output is $yo=y6$ obtained using least squares with a noncausal FIR model with r = 24 and d = 25. At approximately t = 10 s, a fault occurs to Mic2. Figure 15 shows a plot of $E(k|θ̂r,d,ℓ,w)$ for $T(q−1,θ̂r,d,ℓ),$ where w = 500 steps. Note the change in $E(k|θ̂r,d,ℓ,w)$ at approximately t = 10 s due to the faulty microphone.

Next, suppose that $u1,u2,u3$ are as shown in Eqs. (33)–(35). At approximately t = 10 s, a fault occurs to Mic3. Figure 16 shows a plot of $E(k|θ̂r,d,ℓ,w),$ where w = 500 steps, for the identified transmissibility $T(q−1,θ̂r,d,ℓ)$ obtained by driving the speakers with white excitations. Note the change in $E(k|θ̂r,d,ℓ,w)$ at approximately t = 10 s due to the faulty microphone.

Consider the transmissibilities $T1,…,T6$ whose pseudo inputs and pseudo outputs are as shown in Table 1. We assume that, during the experiment, at most one microphone is faulty. Figure 17 shows an algorithm to localize the fault whether it has occurred in the acoustic system or in one of the microphones.

Suppose that only Spk1 is operating and is driven with a band-limited white noise with a bandwidth of 500 Hz. We introduce a change in the dynamics of the acoustic system by removing the blanket at approximately t = 10 s. Figure 18 shows a plot of $E(k|θ̂r,d,ℓ,w)$ for the identified transmissibilities $Ti(q−1,θ̂r,d,ℓ)$ of $Ti$ for $i=1,…,6.$ Note from Fig. 18 that a fault has been detected by all transmissibilities, and thus, using Fig. 17, we correctly conclude that the fault is due to the change in the dynamics of the acoustic system.

Next, suppose that $u1,u2,u3$ are as shown in Eqs. (33)–(35). We introduce a fault in Mic2 at approximately t = 10 s. Figure 19 shows a plot of $E(k|θ̂r,d,ℓ,w)$ for the identified transmissibilities $Ti(q−1,θ̂r,d,ℓ)$ of $Ti$ for $i=1,…,6.$ Note from Fig. 19 that a fault is detected by all transmissibilities except $T5$ and $T6$⁠, and thus, using Fig. 17, we correctly conclude that the fault is in Mic2.

Next, suppose that $u1,u2,u3$ are as shown in Eqs. (33)–(35). We introduce a fault in Mic3 at approximately t = 10 s. Figure 20 shows a plot of $E(k|θ̂r,d,ℓ,w)$ for the identified transmissibilities $Ti(q−1,θ̂r,d,ℓ)$ of $Ti$ for $i=1,…,6.$ Note from Fig. 20 that a fault is detected by all transmissibilities except $T6$⁠, and thus, using Fig. 17, we correctly conclude that the fault is in Mic3.

Next, suppose that $u1,u2,u3$ are as shown in Eqs. (33)–(35). We introduce a fault in Mic5 at approximately t = 10 s. Figure 21 shows a plot of $E(k|θ̂r,d,ℓ,w)$ for the identified transmissibilities $Ti(q−1,θ̂r,d,ℓ)$ of $Ti$ for $i=1,…,6.$ Note from Fig. 21 that a fault has been detected by all transmissibilities except $T2$⁠, and thus, using Fig. 17, we correctly conclude that the fault is in Mic5.