Evaluating internal condition of hardwood logs based on AR-minimum entropy deconvolution combined with wavelet based spectral kurtosis approach

2.3.1 AR model

Assuming that x_k is the zero-mean stationary signal sequence, the pth-order AR model of x_k can be defined as follows:

where, φ_i are the AR model coefficients, and ε_k is the residual error which represents the difference between the actual signal and the AR prediction. For non-stationary signal, the residual error is a mixture of noise and transient components of signal.

The order selection of the AR model is the key. The lower an order is, the smoother the estimated spectral peaks are, which causes a significant prediction deviation, however, the higher the order is, the more intense the estimated spectrum fluctuation are, which leads to pseudo spectrum peaks and statistical instability (Wang and Wong 2000). In this paper, the Akaike Information Criterion (AIC) is used to determine the optimal order of the AR model. When AIC gets the minimum, the optimal order of the model can be obtained (Bengtsson and Cavanaugh 2006).

2.3.2 Minimum entropy deconvolution (MED) technique

Assuming that ε_k is an acoustic sequence filtered out stationary components, according to the theory of signal and system it can be described as

where, g_k is the defect signal component, n_k is the noise, h_k is the impulse response function of the transmission path, and symbol * represents convolution.

The basic purpose of MED is to seek a filter that eliminates the influence of transmission path and noise and make the output approximate the defect signal by solving an optimum set of the filter coefficients. The realization process of MED filter is as follows:

1. Assuming L, ε_k and y_k is the length, input and output of a filter f_k respectively, then

1. Calculating the kurtosis of the output y_k

1. Solving the filter coefficients by maximizing the kurtosis K_f, i.e., taking the first partial derivative with respect to f_k for Eq. (7) and set to zero

1. Taking the partial derivative of f_k for Eq. (6)

1. Substituting Eqs. (6) and (9) into Eq. (8), the Eq. (8) can be rewritten as follows:

Eq. (10) can be written in matrix form:

where, the column vector c is a weighted sum of cross-correlation of the input and output of filter f_k, T is the Toeplitz autocorrelation matrix of the input, and f is the coefficient vector of the filter f_k.

1. Using an iterative algorithm for solving the filter coefficients of the Eq. (11) as follows:

   • (a) Calculating the autocorrelation matrix T.

   • (b) Initializing the filter coefficients f⁽ⁱ⁾ (i = 0) as a unit impulse sequence.

   • (c) Computing the output y⁽ⁱ⁾ using ε_k and f⁽ⁱ⁾ according to Eq. (6).

   • (d) Calculating the column vector c⁽ⁱ⁺¹⁾ using ε_k and y⁽ⁱ⁾ according to Eq. (10).

   • (e) Solving new filter coefficients: f(i+1)=T−1c(i+1).

   • (f) Computing iteration termination conditions:

The iteration ends if the expected value E(δe) is less than the set tolerance limit. Otherwise, update the filter coefficients f⁽ⁱ⁾=f^(i + 1) and repeat the next iteration loop from step (c).

2.3.3 Feature parameter extraction based on AR-MED

When the cross-section of a log is subject to unit longitudinal impact force, the acoustic signal generated reflects the internal quality of the log. If a log is in homogeneous and continuous condition, the response signal can be represented as the superposition of an infinite number of principal vibration modes according to the vibration theory, whereas if there exist internal defects (such as void, decay, crack, etc.) in the log, mutation, distortion or attenuation components caused by defects inevitably skew or modulate the response signal. Therefore, the acoustic signal x_k collected from a log can be expressed as the convolution between the response function of transmission path h_k with the combined signal of periodic components d_k (deterministic signal), defect components g_k and background noise n_k, namely

Figure 1 illustrates AR-MED based extraction process of feature parameters. The steps are as follows:

1. Predicting and subtracting the deterministic components using the AR model from the response signal to obtain the residual signal containing the defect components.

2. Using the residual signal as the input of the MED filter, and making the filtered output signals approximate to the defect components through multiple iterative operations.

3. Computing the kurtosis of the output signal, and applying them as the global parameter for quality evaluation or internal defect prediction of logs.

4. Performing complex Morlet wavelet transform on the output signal and computing the SK, and then using the SK for approximate identification of the internal defect types of logs.

Figure 1:

AR-MED based extraction process of acoustic feature parameters.

4.1.1 AR-MED based kurtosis extraction

The kurtosis for all hardwood logs are summarized in Table 2. To demonstrate the effectiveness of the AR-MED algorithm in enhancing the mutation signal, the kurtosis of the residual signal filtered by the AR model and that of the output filtered by the AR-MED were compared using acoustic signal of log no.14 (red oak) as an example (Figure 2).

Figure 2:

AR-MED based defect components separation and kurtosis extraction (log no. 14-red oak): (a) time history of the acoustic response signal; (b) frequency spectrum of the response signal; (c) distribution between the order of the AR model and AIC; (d) residual signal filtered by the AR (K=4.64); (e) quasi -defect signal filtered by the AR-MED (K=18.00); (f) frequency spectrum of the quasi-defect signal.

The response signal of log no. 14 has almost regular variation of waveform with no clutter as shown in Figure 2a, and its corresponding spectrum (Figure 2b) does not appear to contain any anomaly frequencies except for the fundamental frequency and the corresponding higher-order multiple frequency, so it is difficult to determine whether defects exist in the log or not. But from the dissected results, there do exist some defects such as relatively severe decay, void, ring shakes, etc. in the log (Figure 3) and the defect ratio calculated is about 10% according to the ground-true map (details can be found in Xu et al. 2019). Comparing to the main vibration components, mutation components caused by defects are relatively weak and almost submerged in the dominant components. AR predictive filtering was performed on the signal to remove the stationary components and obtain the residual signal containing defect components. The distribution of the order of the AR model with respect to the AIC values and the residual signal are shown in Figures 2c, d, respectively. The order of 400 was the optimal order for the AR model that separated the stationary and mutation components. However, the AR filter could not isolate noise from mutation components due to the minimum phase condition assumed in the AR model. The kurtosis of 4.64 was small and could not reflect the extent of log deterioration. To suppress noise and enhance defect-induced mutation components, the residual signal was processed further by a MED filter, and the kurtosis of the output reached 18.00 (Figure 2e). Observing its spectrum shown in Figure 2f, the mutation components were separated effectively from the response signal and enhanced significantly.

Figure 3:

Defects located in the end face and the dissected cross-section of log no.14: (a) void in the end face; (b) rot and void in the cross-section; (c) rot and ring shake in the cross-section.

4.1.2 Relationship between defect ratio and kurtosis

The kurtosis of 15 hardwood logs ranged from 10.8 to 27.5, with a coefficient of variation (COV) of 26.6%. Although a meaningful within species statistical comparison cannot be obtained from the small sample size of the logs in each species, the trend of decreasing kurtosis with the increase of defect ratio is evident. The reason could be attributed to the fact that defects of the logs were mainly rot and void. When these defects are present in logs, the occurrence probability of relatively small signal amplitudes increases, and thus the kurtosis decrease accordingly, and the more severe the defects are, the smaller the kurtosis values are (Mucciardi et al. 2011). Figure 4 shows the data plot of defect ratio and kurtosis of the hardwood logs. An excellent negative relationship was observed between log defect ratio and kurtosis, with the coefficient of determination (R²) of 0.89 for a second-order polynomial regression, which is a significant improvement comparing with the R² of 0.17 for the regression result between defect ratio and acoustic velocity (Xu et al. 2019). If the difference of the MC was ignored, the distribution of the kurtosis with respect to defect ratio seemed to be less affected by log species compared with acoustic velocity. As a result, the kurtosis extracted from defect components isolated could represent internal condition of a log more effectively than the acoustic velocity directly derived from the original signal.

Figure 4:

Relationship between log defect ratio and log kurtosis of hardwood logs.

4.2.1 Complex Morlet wavelet-based SK extraction

The SK expands the concept of the kurtosis, as a global parameter, to that of a frequency function, revealing how the mutation characteristics of a signal are distributed in the frequency domain. Therefore, the SK is a powerful tool for detecting the presence of transients in a signal. Even if the transients were hidden in deterministic components and noise, the SK can identify them by indicating which FBs they appear in.

Since the SK is related to the central frequency and bandwidth of wavelet, the Kurtogram (Antonia and Randall 2006), a fourth-order spectral analysis tool, was first made to select the optimal central frequency (for simplicity and real-time requirement, the bandwidth is fixed of 300 Hz in this paper). A variety of defect types coexist in a log, but in this study, defect types identified were mainly decay, void, crack, or other structural defects. Also, because large SK values located in some FBs were caused by the dominant components (associated with defects), while small SK values were induced by noise, and thus the small SK can be omitted for more effective analysis (Antoni and Randall 2006; De la Rosa et al. 2015). The FBs with the larger SK values were extracted and contrasted to the actual defects from the cutting-open log one by one to determine the corresponding relationship between them.

The complex Morlet wavelet transform was performed further on the AR-MED based output signal of log no.14 to calculate the SK. The distribution map of the SK (Kurtogram) of log no.14 is shown in Figure 5. Obviously, the SK values calculated are different for wavelets with different center frequencies (namely 1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, and 5) for the signal decomposition. The maximum of the SK is used as an indicator for selecting the optimal center frequency. For this signal, when the central frequency of wavelet was of 2, the maximum of SK of 19.68, corresponding to the FB of 1850–2150 Hz, can be obtained, so the number of two is the optimum frequency of wavelet for the signal decomposition, that is, the SK values were obtained based on the optimal wavelet. Because the internal defects of the hardwood log were not single, other FBs with some SK values larger than 40% of the maximum should also be the characteristic intervals of useful signal components. The SK values were sorted and then FBs corresponding to the large values were selected for analysis. It is assumed that the larger the SK value of a certain FB is, the larger the proportion of material components corresponding to the FB in the log is. In addition, some research showed that if the internal defects of a log were dominated by decay or void, the frequency of the acoustic signal transmitted in the log is relatively low (Dackermann et al. 2014; Dunlop 1983; Mucciardi et al. 2011); and if the defects were mainly cracks or splits, the high frequency components would appear in the acoustic signal (Hsiao et al. 2008). Therefore, it is speculated that there exist some defects such as decay, void, etc. in the log if larger SK values were located in low FBs of a signal. On the contrary, defects such as cracks or splits do exist if larger SK values were found in high FBs of the signal. As seen in Figure 5, the largest three SK values were located in FBs of 1850–2150, 1550–1850, and 2450–3050 Hz, respectively, accounting for about 53.7% of the sum of the SK values. Still, the remaining values were significantly smaller, and thus they could be ignored in analysis. Observing the internal condition of log no. 14 (Figure 3), it was found that internal rot, void, and ring shake were the dominant defects, other defects, such as knots, etc., were very small. The predicted defects basically conformed to the actual internal condition of log no. 14, which indicated that the conjecture mentioned above between distribution of the FBs and the corresponding defects should be correct. Therefore, the corresponding relation between the FBs of the SK and defect types could be assumed as listed in Table 3 from which the main types and ranking of severity of defects could be roughly estimated.

Figure 5:

Kurtogram for selecting the optimal center frequency of complex Morlet wavelet (log no. 14-red oak).

According to the above method, the SK values of all 15 hardwood logs were extracted, and the Kurtogram was drawn together and shown in Figure 6. The large SK values in Figure 6 are mainly distributed in the FBs of 650–4550 Hz. To determine the defect types and overcome the observation inconvenience caused by extensively large difference among the SK values of each log, the SK values were normalized by dividing SK located in each FB by the sum of the SK values of all FBs, and in the meanwhile the frequency range of Kurtogram was focused on the interval of 650–4550 Hz. Then the Kurtogram was redrawn in Figure 7 with the smaller SK values being set to zero and only the larger ones kept.

Figure 6:

Kurtogram of all 15 logs based on the complex Morlet wavelet.

Figure 7:

Kurtogram of all 15 logs with the normalized SK values and frequency bands of 650–4550 Hz.

4.2.2 SK based defect identification

The SK of log no.1 was mainly concentrated in the frequency range of 3350–4550 Hz, and the sum of SK from four FBs accounted for about 64.7% of that of all SK values. According to the hypothesis of the corresponding relation between FBs and defects listed in Table 3, there should be no internal defects in the log except end checks, as evidenced by observing the actual sawing results (Figure 8a). For log no.3, there were no apparent differences among the SK values from FBs of 650–3650 Hz, with a COV of 18.2%. It was found that there existed a variety of defects, such as void, decay, ring shake, etc., in the log, and the defects were all significant, without obvious primary and secondary by examining internal condition of the log (Figure 8b). The observed results further validated the speculation in Table 3. The major SK distribution of log no. 4 was located in frequency ranges of 1250–1850and 2750–3650 Hz, accounting for about 62.8% of the total SK values, and the SK of each frequency range was not significantly different. It could be deduced from the Table 3 there existed some defects such as rot, voids, and cracks in the log. The basic consistence between the predicted defects with the actual sawing results (Figure 8c) further verified the correctness of the assumptions in Table 3. Similar to the analysis of the above examples, the feature FBs and predicted defects from the remaining logs were listed in Table 2, no more tautology here. By comparing Table 2 with Table 1, the types and primary and secondary of the main defects predicted were roughly consistent with the sawing results, which confirmed the validity of the relation between FBs and defects presented in Table 3.

Figure 8:

Defects located in the end face and dissected cross-section of logs: (a) end checks in the log no. 1; (b) void, cracks, rot, ring shake, etc. in the log no. 3; (c) void, decay and cracks in the log no. 4; (d) ring shake and decay in the log no. 8.

It should be noticed that log no.8 (cottonwood) was a particular case, i.e., global parameters of the longitudinal acoustic signal, whether kurtosis or acoustic velocity, time centroid and damping ratio presented in Xu et al. 2019, implied that the log could be a sound log, but actually one with abnormal internal condition (Xu et al. 2019). Further examining the distribution of the larger SK of the log, it was found that they were located in FBs of 2750–3050 and 1250–1850 Hz, respectively. According to the assumption presented in Table 3, the log should contain the defects such as cracks and decay. It was demonstrated that these defects did exist in the log by checking the sawing results (Figure 8d). The research results further indicated that the SK could determine the existence and type of major defects which not be detected even by the global acoustic parameters.