Elliptical ring distribution probability-based damage imaging method for complex aircraft structures

In engineering applications, the robustness and effectiveness of damage diagnostic imaging for guided wave-based structural health monitoring could be affected by the complexity of structures. In this study, an elliptical ring distribution probability-based diagnostic imaging algorithm is proposed to mitigate this effect using the estimated wave velocity and damage index. This algorithm improves the ability of damage localization by modifying the defect distribution probability of probability-based diagnostic imaging. The elliptical ring distribution probability of the presence of defect is used for each sensing path in the algorithm. The width of the elliptical ring distribution probability is determined by the range of estimated wave velocity. The amplitude of the elliptical ring distribution probability is determined by the damage index. The damage location is determined by the cross region of different elliptical rings for different sensing paths. The capability of the algorithm is validated by identifying damages at different locations on a complex composite fuselage panel. The results show that the proposed algorithm can identify a single damage accurately and it can identify multiple damages effectively as well.


Introduction
Due to the superior properties of strength, stiffness, weight and corrosion resistance, the composite materials are increasingly used in aircraft structures.Especially in recent years, some complex composite structures have become widely applied to aircrafts.However, one of the major challenges with composite materials is that these materials suffer invisible damage with different failure modes such as matrix crack and delamination.These damages can decrease the composite structure's strength significantly, so they must be detected in time.Traditional non-destructive testing (NDT) techniques, although efficient in the damage detection, are time-consuming and expensive to be applied regularly since they require the structure to be off-service, even require component disassembling.The structural health monitoring (SHM) technology offers an alternative means to detect and monitor damages in structures.
Ultrasonic guided waves have particular advantages for the purpose of SHM because they can travel over long distances and are sensitive to many types of damage such as crack, corrosion and delamination et al. [1][2][3][4][5][6][7][8].Many damage detection algorithms have been developed for guided wave-based SHM.Among these algorithms, damage diagnostic imaging algorithms have attracted much attention because they result in an interpretable and intuitive image to reflect the location of damage.Some representative algorithms are phased array method [9,10], tomography technique [11][12][13], delay-and-sum algorithm [14][15][16], correlation-based imaging technique [17], model-based damage imaging using sparse reconstruction [18,19] and probability-based diagnostic imaging (PDI) method [20].For the phased array method, the localization accuracy highly depends on the density of transducers, because it requires a number of transducers to perform as the actuator.For the tomography technique, it requires dense wave paths for imaging reconstruction, which fairly narrows its application for SHM.For the delay-and-sum algorithm, correlation-based imaging technique, and model-based damage imaging using sparse reconstruction, these methods require accurate wave velocity.However, it is quite difficult to interpret the guided wave signals for acquiring the accurate wave velocity because the propagation mechanism of guided wave is quite complicated in engineering structures due to its geometric and structural complexity, such as overlapping structure, variable thickness structure and reinforced structure.Thus, the robustness and effectiveness of these methods could be affected by the complexity of structures.
In contrast, the probability-based diagnostic imaging (PDI) enables damage imaging with a good quality in a sparse sensor network, and does not require direct interpretation of the guided wave signal.Thus, it has been studied intensively by many researchers [21].Zhao et al. [22] introduced a correlation analysis-based PDI algorithm called as the reconstruction algorithm for probabilistic inspection of defects (RAPID) which used the signal difference coefficients of sensor pairs to reconstruct the defect distribution probability map.The algorithm showed good performance in estimating the location of cracks or corrosion and monitoring its growth in a metallic skin panel of an aircraft wing.Koduru and Rose [23] improved the PDI algorithm by using designed annular array sensors for mode controlling.It showed a remarkable improvement in the ability to distinguish a real corrosion defect from any other water traces on the structure.Wang et al. [24,25] studied the PDI algorithm for the identification of multiple notches using digital damage fingerprints (DDFs) which were used to highlight the changes in signals corresponding to the presence of damage.They also developed a correlation analysis based-PDI algorithm with a concept of virtual sensing paths (VSPs) which was proposed to enhance the performance of the algorithm by increasing the number of sensing paths in data fusion.Wu et al. [26] studied the influences of multiple factors of PDI algorithm on the damage identification, and developed a methodology to determine the parameters for the PDI algorithm.Hua et al. [27] developed a local signal difference coefficient (LSDC) based-PDI algorithm to improve the imaging performance.Liu et al. [28] developed a weight-compensated probability-based diagnostic imaging (WCPDI) to weaken the influence of the dominance effect of the weight distribution for improving the ability of damage localization.
In the aforementioned researches on the PDI algorithm, the damage imaging is based on the relatively inaccurate defect distribution probability, which is linearly decreasing elliptical distribution.In these studies, for the defect distribution probability of an actuator-sensor pair, the probability of the presence of defect is larger as the imaging grid is closer to the direct path of the actuator-sensor pair.Furthermore, the probability is maximum in the direct path of the actuator-sensor pair, no matter where the actual damage is located.According to the common sense, the probability of the presence of defect is larger as the grid closer to the actual damage.Thus, this relatively inaccurate defect distribution probability does not conform to the fact.It could reduce the damage localization accuracy.This can affect the application of the PDI algorithm in real-world practices.In order to improve the damage localization accuracy of PDI algorithm, more spatial damage information provided by the actuator-sensor paths can be used.For complex structures, the spatial damage information is often difficult to be obtained accurately because the propagation mechanism of guided wave is quite complicated in these structures.In this study, a strategy is established to solve the above-mentioned problems.
In order to improve the damage localization accuracy of PDI algorithm for complex structures, modifying the defect distribution probability is an important way.In this study, an elliptical ring distribution probability-based diagnostic imaging (ERDPDI) algorithm is developed to improve the ability of damage localization for complex structures by modifying the defect distribution probability of PDI algorithm.In the ERDPDI algorithm, the probability of the presence of defect is an elliptical ring distribution for each sensing path.For aircraft structures, the composite fuselage panel is an important and complex structure.Due to the material anisotropy and existence of fuselage frame and stiffener, propagation mechanism of guided waves is quite complicated in the composite fuselage panel.It is difficult to identify damage on the composite fuselage panel.Therefore, in this study, the capability of the proposed algorithm is validated by identifying damages at different locations on a complex composite fuselage panel.
The remainder of this paper is organized as follows.Section 2 describes the proposed ERDPDI algorithm based on the brief review of PDI algorithm.Numerical simulations study on ERDPDI algorithm is presented in Section 3. Then the experimental validation for ERDPDI algorithm is presented in Section 4. In order to show the advantage of ERDPDI algorithm, the performance of ERDPDI algorithm is compared with the performance of PDI algorithm in this experimental validation.Finally, conclusions are presented in Section 5.

Brief review of PDI approach
When ultrasonic guided waves travel through a structure, the damage of structure can lead to the scattering of the guided waves.By comparing the guided wave signals before and after damage occurred, the damage information of the structure can be obtained.Guided wave-based SHM utilizes this principle to monitor the damage of structure.It uses the distributed actuator/sensors, which are permanently attached to the structure, to generate the guided wave and measure the arriving waves at sensors.The damage diagnosis is performed through the examination of the arriving waves in comparison with a 'baseline' condition.Many damage diagnosis algorithms have been developed for guided wave-based SHM.The PDI algorithm is one of the representative algorithms.
In the PDI algorithm, the monitoring area is meshed into uniformly distributed grids.A grid, which denotes a region at given position in the monitoring area, is an image point.The probability of the presence of damage at each grid is estimated by summing the signal difference features in elliptical patterns for various pitch-catch transducer pairs.Assuming that there are sensing paths in total involved for damage imaging, the estimation of the probability of the presence of damage at a certain grid ( , ) can be calculated as: where is the correlation analysis-based damage signature for the th sensing path, which can be represented by: where , is the correlation coefficient for the th sensing path, is the covariance of signal and signal of the th sensing path, and are the standard deviations of and .[ ( , )] is the weight distribution function of the th sensing path.It is dependent on the parameter, ( , ), which is defined as the relative distance from grid ( , ) to the th sensing path: where is the distance between the actuator and the sensor for the th sensing path, , ( , ) and , ( , ) are the distances between grid ( , ) and the actuator and sensor for the th sensing path, respectively.The [ ( , )] is assumed to be linearly decreasing elliptical distribution as shown in Fig. 1.The size of the effective elliptical distribution area is controlled by a scaling parameter .The weight distribution function can be expressed as: It has been shown that the damage localization result of the PDI algorithm is highly dependent on the weight distribution.According to the above, the weight increases with a decrease in the relative distance, indicating that the grid closer to the sensing path has larger probability of the presence of defect.The probability of the presence of defect is maximum in the direct path of the actuator-sensor pair, no matter where the actual damage is located.This is a relatively inaccurate defect distribution probability.According to the common sense, the probability of the presence of defect is larger as the grid closer to the actual damage.Thus, the weight distribution in the PDI algorithm can reduce the damage localization accuracy.

ERDPDI approach
Due to the relatively inaccurate defect distribution probability of PDI algorithm, the ERDPDI algorithm improves the accuracy of defect distribution probability of PDI algorithm using more spatial damage information.For an actuator-sensor pair, the actual damage location cannot be determined, but the possible damage location can be acquired.Using the guided wave group velocity and time of flight (TOF) of scattered wave caused by damage, the possible damage location is an ellipse with actuator and sensor as two foci, as show in Fig. 2. In engineering structures, due to geometric and structural complexity such as material anisotropy, variable thickness structure and reinforced structure, the guided wave group velocity usually cannot be calculated accurately.Therefore, for an actuator-sensor pair, the position of ellipse for the possible damage location cannot be determined accurately.Thus, the region of the position of ellipse for the possible damage location is considered.In order to determine the region, the range of estimated guided wave velocity is used.It supposes that the range of estimated guided wave velocity is [ , ].The region of the possible damage location for the th sensing path is determined as: where is the region of the possible damage location for the th sensing path, 1 is the nearest estimated possible damage location for the th sensing path, which is calculated by: where is the lower limit of estimated guided wave velocity, is the TOF of scattered guided wave caused by damage for th sensing path. 2 is the farthest estimated possible damage location for the th sensing path, which is calculated by: where is the upper limit of estimated guided wave velocity.For anisotropic material structures, the and can be determined by the minimum guided wave velocity and maximum guided wave velocity of structure along different directions.For isotropic material structures, a method for determination of and is proposed as: where is the theoretical value of guided wave velocity, is the allowable damage location error, is the maximum distance between arbitrary two transducers in transducers network.In the above calculation of the region of the possible damage location, the TOF estimation is an important issue.To tackle this issue, a variety of methods have been developed.Representative methods include continuous wavelet transform (CWT), Wigner-Ville distribution (WVD), Hilbert-Huang transform (HHT), short-time Fourier transform (STFT), matching pursuit decomposition (MPD), and Hilbert transform [29][30][31][32].In this study, the Hilbert transform was used to estimate the TOF by extracting the time of the maximum of the envelope of the concerned guided wave signal.The envelope of the guided wave signal is calculated as: where ( ) is the guided wave signal, ̂ ( ) is the Hilbert transform of ( ) , ( ) is the envelope of ( ).When the region of the possible damage location for the th sensing path is acquired, the weight distribution function [ ( , )] of the th sensing path for grids of monitoring area can be determined.In the ERDPDI algorithm, the weight distribution function [ ( , )] is an elliptical ring for the th sensing path, as shown in Fig. 3.The probability of the presence of defect is maximum in the region of the possible damage location for the th sensing path, which is the region between 1 and 2 as shown in Fig. 3.When the grid departs from the region between 1 and 2 , the weight decreases quickly.The [ ( , )] is expressed as: where ( , ) is the propagation distance of scattered guided wave caused by the possible damage at grid ( , ) for the th sensing path, which is written as: ( , ) = , ( , ) + , ( , ).
(11) For above mentioned weight distribution function, the probability of the presence of defect in region of the possible damage location for all sensing paths is the same.This obviously does not accord with the fact.Thus, the probability of the presence of defect in the region of the possible damage location for each sensing path must be modified.For a damage, it usually has different effects on a different sensing path.Thus, this fact can be used to improve the probability of the presence of defect for each sensing path.Usually, the damage index is used to evaluate the effect of damage on the guided wave signal for a sensing path.Thus, the damage index is used as weight to modify the probability of the presence of defect for a sensing path.For the th sensing path, the probability of the presence of defect at grid ( , ) is expressed as: where is the damage index for the th sensing path.Finally, the probability of the presence of defect at grid ( , ) can be acquired by fusion of all sensing paths.It is expressed as: The flowchart of the ERDPDI algorithm is shown in Fig. 4.

Numerical simulations study on ERDPDI algorithm
In this section, the principle of ERDPDI algorithm was studied by the numerical simulation.Due to the principle study, the verification of ERDPDI algorithm was carried on an aluminum plate which was a simple structure.Numerical simulation has been carried out using commercial finite element method (FEM) software ABAQUS.The size of the aluminum plate model is 200 × 200 × 2 mm.The material properties of the plate are Young's modulus = 70 GPa, Poisson's ratio = 0.33 and density = 2700 kg/m -3 .Four positions were selected for actuating and receiving the guided waves.When one position was used for actuating, the others were used for receiving.Each position was used as actuator for once.In order to describe the selected four positions, a coordinate system was employed with the plate.The origin of coordinate was the center of the plate.The dimensions of plate, the four positions and coordinate system are shown in Fig. 5.The coordinates of four positions are listed in Table 1.The excitation pulse, which was along the normal direction of plate, was loaded on the node at the one of selected four positions.The displacements along the normal direction of plate for nodes at the other selected positions, were used as received guided wave.The excitation pulse was a five-cycle tone burst modulated by a Hanning window.The central frequency of excitation was 90 kHz.The time span of excitation pulse was 2×10 -4 s.An Explicit dynamic step was created for the simulation.The element used for modeling is S4 which is a 4-node shell element.Through the convergence analysis, the time step was set to be 0.02 μs and element size was set to be 0.4 mm.Attenuation was not considered in numerical modeling.After the simulation of baseline signals, a notch was modeled in the plate model, which size was 1×10×2 mm, as shown in Fig. 6.The coordinate of the notch center was (25,25).Before the damage imaging, the values of and in the ERDPDI algorithm must be determined.Due to the isotropic material of plate, actuator-sensor path 1-4 was arbitrarily selected to calculate the wave group velocity.In the group velocity calculation, the TOF calculation is shown in Fig. 7. Due to small displacements for receiving positions, in Fig. 7, the signal amplitude is the ratio of the actual signal amplitude to the maximum of baseline signal amplitude.As attenuation was not considered in numerical modeling, using the Hilbert transform, the TOF was calculated by extracting the time of the maximum of the envelope of the first wave packet for the response signal.The TOF calculation of scattered signal was the same.As an example, the TOF calculation of scattered signal for actuator-sensor path 1-4 is also shown in Fig. 7. Through calculation, the wave group velocity is 1995 m/s.The maximum distance between arbitrary two transducers in the model is 100√2 mm.If the allowable damage location error is assumed to be 5 mm, using Eq. ( 8), and in the ERDPDI algorithm are calculated as 1925 m/s and 2066 m/s, respectively.Due to reciprocity, there were totally 6 sensing paths for damage imaging in this simulation.For the th sensing path, the used in ERDPDI algorithm was the correlation analysis-based damage signature, .Fig. 8 shows the probability images of ERDPDI algorithm for notch identification in the numerical simulation using different actuator-sensor paths.The white circle dots denote the four positions selected for actuating and receiving the guided waves.The center location of the notch is marked by the white cross.The grid with the highest probability value in the image denotes the location of identified damage.The images have the dimensions of 60×60 mm 2 , that is, 60×60 uniformly distributed grids with an interspatial distance of 1 mm.Fig. 8(a)-(f) show the diagnostic images of ERDPDI algorithm using one actuator-sensor path.As shown in Fig. 8(a)-(f), if damage location is not on an actuator-sensor path, the probability of the presence of defect is an elliptical ring.Otherwise, as the TOF of scattered signal is so small, the elliptical ring distribution probability of defect changes to elliptical distribution, as shown in Fig. 8(b).According to the principle of ERDPDI algorithm, three actuator-sensor paths can identify a damage.For example, three actuator-sensor paths were arbitrarily selected to identify notch, as shown in Fig. 8(g) and (h).The notch can be accurately identified by three actuator-sensor paths.However, in engineering applications, as the number of damages is usually not known before the damage monitoring, all actuator-sensor paths are used to identify damage for the ERDPDI algorithm.Fig. 8(i) shows the diagnostic image of ERDPDI algorithm using all actuator-sensor paths.As shown in Fig. 8(i), the ERDPDI algorithm can identify damage accurately.

Specimen and experimental setup
The composite fuselage panel is a representative and important aircraft structure.In recent years, this kind of composite structure has been widely applied to aircrafts.Thus, the need for damage monitoring of this kind of composite structure is also urgent.However, due to geometric and structural complexity, it is difficult to identify damage in the composite fuselage panel.Thus, experiments were carried on a complex carbon fiber composite fuselage panel to validate the performance of ERDPDI algorithm in this study.The panel has three I-shaped stiffeners and two fuselage frames.The material of panel is T700/BA9916.Sixteen lead zirconate titanate (PZT) sensors (P51) were bonded on the panel with GLEIHOW302 adhesive.The dimensions of PZT sensors are 8 mm in diameter and 0.45 mm in thickness.The dimensions of panel and PZT sensors placement are shown in Fig. 9.An integrated SHM system was used to generate and acquire the guided waves in the experiment.The experimental setup is shown in Fig. 10.
The excitation signal was a five-cycle tone burst modulated by a Hanning window.The sampling rate for sensor data acquiring was set at 10 MHz.According to the characteristics of guide wave signals tested on the panel for excitation signal center frequencies range from 50 to 150 kHz with an interval of 20 kHz, the excitation signal center frequency for this experimental validation was selected as 90 kHz.GUOQIANG LIU, YINGCHUN XIAO, HUA ZHANG, GEXUE REN PDI algorithm, the performance of ERDPDI algorithm is compared with the performance of PDI algorithm in this experimental validation.Wu et al. [26] studied the influence of two different networks of sensing paths on PDI algorithm.The result shows that group CIR (enclose a circle monitoring area) has a better performance than group SQU (enclose a square monitoring area) for the PDI algorithm.In order to make the PDI algorithm have a good performance, this experimental validation used CIR network for damage identification.Eight PZT sensors were selected for CIR network in this experimental validation.In order to describe the location of the PZT sensors, a coordinate system was employed with the plane of the monitoring area spanned by the horizontal, , and vertical, , axes.The selected PZTs for CIR network and the coordinate system are shown in Fig. 11.The red circle dots denote the selected PZTs.The coordinates of selected PZTs are listed in Table 2. Due to reciprocity, there were totally 28 sensing paths in CIR for damage imaging.Before the start of damage identification validation, the values of and in the ERDPDI algorithm must be determined.Under the 90 kHz excitation center frequency, the mode amplitude is dominant.Thus, in this experimental validation, the mode was used to damage identification.Using positive axis as 0° direction and positive y axis as 90° direction, the calculated mode group velocity along different directions is shown in Fig. 12.The actuator-sensor paths for calculating the mode group velocity along different directions are listed in Table 3.In the group velocity calculation, the TOF determination of mode is important.As an example, the TOF calculation for actuator-sensor path 3-4 is shown in Fig. 13.The Hilbert transform was used to estimate the TOF by extracting the time of the maximum of the envelope of the response signal for actuator-sensor path 3-4.In Fig. 12, the mode group velocity anisotropy is observed clearly.According to the minimum and maximum group velocity along different directions, and in this experimental validation for the ERDPDI algorithm are 1251 m/s and 1516 m/s, respectively.To reduce the experiment cost, artificial damages were introduced with M16 bolts (16 mm in diameter) bonded on the surface of the panel in this study.The diameter of the bolt was approximately equal to the maximum wavelength of mode along different directions.Just as the effect of damages to the structure, the geometry continuum is interrupted by the bolt.The procedure of each damage identification validation case is as follows:

Results and discussion
In order to make the damage identification validation representative, according to the structure symmetry, the monitoring area of the composite fuselage panel was divided into three types of small areas, as shown in Fig. 14.Area A denotes the area which has no stiffener and is far away from the fuselage frame.Area B denotes the area with stiffener.Area C denotes the area which has no stiffener and is close to the fuselage frame.The damage identification must be verified in the three types of small areas.To validate the performance of single damage identification using ERDPDI algorithm, several single artificial damage cases were introduced on the composite fuselage panel at different locations.After the measurement of the guided wave signals, the damage images were reconstructed.In the damage imaging, the scattered signal was the difference between the current signal and baseline signal.The TOF of the scattered signal was estimated by the Hilbert transform, which was used to estimate TOF by extracting the time of the maximum of the envelope of the scattered signal caused by damage.As an example, Fig. 15 shows the signals and TOF estimation of the scattered signal from actuator-sensor path 5-6, when a single artificial damage is near actuator-sensor path 5-6.The used in ERDPDI algorithm was the correlation analysis-based damage signature, , for the th sensing path, which was the same with PDI algorithm.The scaling parameter used in PDI algorithm was determined by testing, which changed from 0.1 to 1.2 by the step of 0.1.The selected scaling parameter was 0.7.4. The point with the pixel peak value of the diagnostic images is regarded as the damage location.The error in table is defined as the distance between the peak point of the individual image and actual location of damage.As shown in Table 4, the PDI localization error is generally greater than the ERDPDI localization error.
By performing the aforementioned experiments, the accurate identification of single damage using ERDPDI algorithm is achieved.In order to validate the performance of multiple damage identification using ERDPDI algorithm, two artificial damages were simultaneously introduced on the composite fuselage panel.Fig. 17 shows the probability images of PDI algorithm and ERDPDI algorithm for two-damages identification.Fig. 17   The ERDPDI algorithm can effectively identify two damages.The two damages can be clearly distinguished in the probability images of the ERDPDI algorithm.Fig. 17(b) and (f) show that some artifacts appear when damage is near PZTs.This is because the elliptical rings for distribution probability of the presence of the damage become smaller when damage is near PZTs.When the elliptical rings become smaller, the number of cross regions of different elliptical rings for sensing paths with one same PZT may increase.Therefore, the increased cross regions may produce the artifacts.
The damage localization results of two artificial damages are listed in Table 5.As two damages cannot be identified effectively by the PDI algorithm, Table 5 only lists the damage localization results of the ERDPDI algorithm.By comparing the results of Table 5 with Table 4, it can be seen that the damage localization error of two damages is greater than the damage localization error of a single damage for the ERDPDI algorithm.This is because the superposition effect of two damages on the scattered signal.

Conclusions
This study developed an ERDPDI algorithm to mitigate the effect of structures complexity on damage imaging using estimated wave velocity and damage index.The algorithm improves the ability of damage localization by modifying the defect distribution probability of PDI algorithm.In the ERDPDI algorithm, the probability of the presence of defect is an elliptical ring distribution for each sensing path.The width of the elliptical ring distribution probability is determined by the range of estimated wave velocity.The amplitude of the elliptical ring distribution probability is determined by the damage index.The effectiveness of the algorithm was assessed by identifying artificial damages on a complex composite fuselage panel.The results show that the damage localization accuracy of ERDPDI algorithm is better than that of PDI algorithm.The proposed algorithm can identify single damage accurately and it can identify multiple damages effectively as well.The proposed algorithm showed great potential for damage identification for SHM of complex aircraft structures.

Fig. 1 .
Fig. 1.Elliptical distribution function of probability of presence of defect for PDI algorithm

Fig. 2 .
Fig. 2. Illustration of possible damage location for an actuator-sensor pair

Fig. 3 .
Fig. 3. Elliptical ring distribution function of probability of presence of defect for ERDPDI algorithm

Fig. 11 .
Fig. 11.Schematic of selected PZTs for CIR network and coordinate system

Fig. 12 .
Fig. 12. mode group velocity along different directions (a) Baseline signals were collected before the artificial damage was introduced.(b) The artificial damage was introduced with bolts bonded on the surface of the panel.(c) Current signals were collected after the artificial damage had been introduced.(d) The damage image was reconstructed.

Fig. 15 .Fig. 16 .Fig. 17 .
Fig. 15.Signals and TOF estimation of scattered signal from actuator-sensor path 5-6 (a)  and (b) show the damage identification for damages located at (170, 399) and (387, 300).Fig.17 (c) and (d) show the damage identification for damages located at (176, 240) and (210, 102).Fig.17 (e) and (f) show the damage identification for damages located at (225, 88) and (280, 352).As shown in Fig.17, PDI algorithm cannot effectively identify two damages in this experiment.The indicated damage area of PDI algorithm is a simply connected region which cannot indicate the number of damage.

Table 2 .
Locations of PZTs for CIR network

Table 3 .
Actuator-sensor paths for calculating

Table 4 .
Damage localization results of single artificial damage (units: mm)

Table 5 .
Damage localization results of two artificial damages (units: mm)