Guided Wave Phase Velocity Dispersion Reconstruction Based on Enhanced Phased Spectrum Method

Fibre-reinforced composite laminates are frequently used in various engineering structures, due to their increased weight-to-stiffness ratio, which allows to fulfil certain regulations of CO2 emissions. Limited inter-laminar strength makes composites prone to formation of various defects, which leads to progressive degradation of residual strength and fatigue life of the structure. Using ultrasonic guided waves is a common technique for assessing the structural integrity of composite laminates. Phase velocity is one of the fundamental characteristics of guided waves and can be used for defect detection, material property estimation, and evaluation of dispersion. In this paper, a phase velocity reconstruction approach, based on the phase-shift method, was proposed, which uses frequency sweep excitation to estimate velocity at specific frequency harmonics. In contrast to the conventional phase spectrum technique, the proposed approach is applicable to the narrowband piezoelectric transducers and suitable for the reconstruction of dispersion curves for direct, converted, and multiple co-existing modes with high accuracy. The proposed technique was validated with finite element simulations and experiments, both on isotropic and anisotropic structures, analysing the direct, converted, and overlapped modes. The results demonstrated that, using the proposed technique, the phase velocity dispersion can be reconstructed at −20 dB level bandwidth of the transducer, with a relative error of ±4%, compared to the theoretical velocity predictions.


Introduction
The composites market is one of the strategic development areas of the European Union, which aims to strengthen their competitiveness and extend the use of composites in the sectors of aerospace, automotive, and renewable energy [1]. The current EU demand of carbon fibre is estimated to be 35% of the global demand, and it will have an annual growth of 10-12% [1,2], while the UK market will grow from 2.5 billion up to 10 billion pounds a year by 2030 [3][4][5]. Fibre-reinforced composite laminates can offer increased strength-and stiffness-to-weight ratios, which allow for meeting the demanding requirements of CO 2 emissions. However, composites have limited interlaminar strength and are prone to formation of fibre breakage, matrix cracking, delaminations, porosity, and other structural defects. Such defects are usually hidden and progressively degrade the residual strength and fatigue life, eventually leading to sudden structural failure. Using ultrasonic guided waves is a common method for periodic inspection and monitoring of structural integrity of plate-like composite laminates, that offers large inspection areas and sensitivity to structural damage of various kinds [6][7][8]. To date, many studies are available that employ guided waves for the detection and quantification of impact damage [9][10][11][12], delaminations [13][14][15][16][17], and other defects in composite laminates. Guided wave propagation in composites is determined by many factors, including, but not limited to, multi-layered structure and anisotropy, object boundaries, dispersion, multiple co-existing modes, and mode conversion. Phase velocity is one of the fundamental properties of guided wave modes that depends on composition, structural integrity, elastic properties, and frequency-thickness product of composite.
phase velocity values can be collected at different frequencies, allowing us to achieve a wideband reconstruction. The validity of the approach is demonstrated through simulations and experiments by reconstructing the phase velocities of S 0 and converted A 0 modes, as well as identifying guided wave modes in complex multimodal signals.
In contrast to the classic phased spectrum method, the proposed approach allows to reduce the relative error of the phase velocity reconstruction from ±11% to ±4% and increase significantly the reconstruction bandwidth from −6 dB to −20 dB of the ultrasonic probe. As a result, using only two signals, measured in close proximity, the proposed phased spectrum method can achieve the reconstruction accuracy and bandwidth, which, to date, could be achieved only with techniques that include scanning of the sensor over a sufficiently large area.

Description of Proposed Phase Velocity Estimation Method
The proposed phase velocity reconstruction approach employs a classic phase-shift method to estimate the velocity values at specific frequencies that correspond to peak values of the magnitude spectra of received signal. By repeating this procedure at different excitation frequencies, velocity values can be reconstructed at wide band, covering the entire bandwidth of the transducer. Variation of the excitation frequency allow different harmonics to be enhanced or suppressed, which is the key factor if reconstruction is performed at peak values of magnitude spectra only. The algorithm of the proposed method can be summarized with the following steps: 1.
The transducer is driven by a burst at a central frequency of f 1 , and the waveforms u r1f1 (t) and u r2f1 (t) are registered with receivers r 1 and r 2 , each positioned at a distances d 1 and d 2 from the source (see Figure 1a for reference). 2.
The waveforms u r1f1 (t) and u r2f1 (t) are windowed using the tapered cosine window w(t) to isolate the wave packets of particular mode (see Figure 1b): u r 1 f 1 w (t) = u r 1 f 1 (t)·w(t − t 1 ), u r 2 f 1 w (t) = u r 2 f 1 (t)·w(t − t 2 ) (1) where u r1f1w (t) and u r2f1w (t) represent the windowed versions of the waveforms u r1f1 (t) and u r2f1 (t), respectively; t 1 and t 2 correspond to the time instances of the maximum amplitude of the wave packet.
entire band of transducer. Two sensors, positioned in close proximity, are used to record signals propagated through the structure and estimate the phase-shift between the signals. At each excitation frequency, the reconstruction of phase velocity is performed at specific frequency components only, which correspond to the peak values of the magnitude spectra. These peak frequencies depend on the frequency response of the excitation signal; hence, phase velocity values can be collected at different frequencies, allowing us to achieve a wideband reconstruction. The validity of the approach is demonstrated through simulations and experiments by reconstructing the phase velocities of S0 and converted A0 modes, as well as identifying guided wave modes in complex multimodal signals.
In contrast to the classic phased spectrum method, the proposed approach allows to reduce the relative error of the phase velocity reconstruction from ±11% to ±4% and increase significantly the reconstruction bandwidth from −6 dB to −20 dB of the ultrasonic probe. As a result, using only two signals, measured in close proximity, the proposed phased spectrum method can achieve the reconstruction accuracy and bandwidth, which, to date, could be achieved only with techniques that include scanning of the sensor over a sufficiently large area.

Description of Proposed Phase Velocity Estimation Method
The proposed phase velocity reconstruction approach employs a classic phase-shift method to estimate the velocity values at specific frequencies that correspond to peak values of the magnitude spectra of received signal. By repeating this procedure at different excitation frequencies, velocity values can be reconstructed at wide band, covering the entire bandwidth of the transducer. Variation of the excitation frequency allow different harmonics to be enhanced or suppressed, which is the key factor if reconstruction is performed at peak values of magnitude spectra only. The algorithm of the proposed method can be summarized with the following steps: 1. The transducer is driven by a burst at a central frequency of f1, and the waveforms ur 1 f 1 (t) and ur 2 f 1 (t) are registered with receivers r1 and r2, each positioned at a distances d1 and d2 from the source (see Figure 1a for reference). 2. The waveforms ur 1 f 1 (t) and ur 2 f 1 (t) are windowed using the tapered cosine window w(t) to isolate the wave packets of particular mode (see Figure 1b): where ur 1 f 1 w(t) and ur 2 f 1 w(t) represent the windowed versions of the waveforms ur 1 f 1 (t) and ur 2 f 1 (t), respectively; t1 and t2 correspond to the time instances of the maximum amplitude of the wave packet.
(a) (b) Figure 1. (a) The example of the waveform, registered with receivers r 1 and r 2 , at distances d 1 and d 2 ; (b) the illustration of waveform windowing to isolate the wave packet of single mode.

3.
Each waveform, u r1f1w (t) and u r2f1w (t), is shifted in the time domain by −t m1 and −t m2 , to avoid the uncertainties in the phase unwrapping procedure. The waveforms can be shifted according to the centroid of signal [35] or maximum value of the Hilbert envelope [36], in case of moderate dispersion: where HT denotes the Hilbert transform; t m1 and t m2 are the time instances, which corresponds to the maximum of Hilbert envelope, in such a way that the influence of the signal delay due to phase velocity is compensated. The shift in time domain is illustrated in Figure 2a. 3. Each waveform, ur 1 f 1 w(t) and ur 2 f 1 w(t), is shifted in the time domain by −tm 1 and −tm 2 , to avoid the uncertainties in the phase unwrapping procedure. The waveforms can be shifted according to the centroid of signal [35] or maximum value of the Hilbert envelope [36], in case of moderate dispersion: = arg max HT , = arg max HT , where HT denotes the Hilbert transform; tm 1 and tm 2 are the time instances, which corresponds to the maximum of Hilbert envelope, in such a way that the influence of the signal delay due to phase velocity is compensated. The

4.
The complex frequency spectra of each time-shifted waveform, u r1f1s (t) and u r2f1s (t), is obtained employing the Fourier transform: where FT represents the Fourier transform.

5.
The phase difference ∆φ(f ) between shifted signals u r1f1s (t) and u r2f1s (t), is estimated for a given frequency band f (see Figure 2b): where Im and Re represent the imaginary and real of the complex Fourier spectra. Note that the phases α r1f1 (f ) and α r2f1 (f ) are calculated in a range of [−π . . . π] radians. If the true phase of the particular frequency is less than −π radians, it will be represented below the π radians. This means that some discontinuities will appear, in case the phase goes beyond the ±π radian limit. Therefore, the phases α r1f1 (f ) and α r2f1 (f ) have to be unwrapped.

6.
The phase velocity, as a function of frequency, is calculated at particular frequencies, f 1,k1 , using a modified version of the phase spectrum method: where f 1,k1 are the frequencies that corresponds to the peak values of the magnitude spectra |U r1f1 (jf )| at excitation frequency f 1 ; k = 1 ÷ K 1 , K 1 -is the total number of detected peaks at excitation frequency f 1 , and d is the separation distance between the receivers r 1 and r 2 (d = d 2 − d 1 ). The frequency selection for phase velocity estimation is illustrated in Figure 2c.

7.
The intermediate values of the phase velocities at other frequencies are obtained by changing the excitation frequency to f 2 and repeating the whole routine described above. The final result is obtained by combining the calculations at different excitation frequencies f 1 . . . f N : where N is the number of excitation frequencies used to drive the emitter.
The method presented above is applicable to flat structures with uniform thickness, which can be multi-layered, anisotropic, or isotropic. In contrast to the conventional phase spectrum method, it provides better accuracy of velocity estimation, which will be demonstrated in the subsequent Chapter.

Experimental Validation on Isotropic Samples
In this section, the proposed phase velocity estimation approach is validated with the appropriate experiments. For this purpose, the phase velocity values, extracted with the proposed approach, are compared with the theoretical calculations, which were considered a reference. In this study, the velocities of the S 0 mode in the aluminium sample will be analysed.
The experiments were carried out on the aluminium alloy 2024 T6 plate, which was 2 mm thick, 650 mm wide, and 1250 mm long. The well-known isotropic material was deliberately selected for this study, in order to be able to compare the experimental results with the theoretically estimated values. The S 0 mode was launched into the structure by attaching the thickness mode transducer to the edge of the Al plate, as is shown in Figure 3. For the reception, two transducers, r 1 and r 2 , possessing the same characteristics, were bonded perpendicularly to the upper surface of the specimen at distances d 1 = 450 mm and d 2 = 550 mm from the source (see Figure 3).
In this paper, transducers with a central frequency of 240 kHz and bandwidth of 340 kHz at −6 dB level were used. The frequency response of the probe can be seen in Figure 4a. To reconstruct the dispersion curve under the wide band, two different scenarios employing the square pulse excitation were used, as follows: n 1 = 3 cycles, f 1 = 150 kHz; and n 2 = 3 cycles, f 2 = 200 kHz. Such excitation frequencies were deliberately selected, according to the magnitude spectrum of excitation pulse, which can be seen in Figure 4b. The results, presented in the figure, demonstrate that a minor shift of excitation frequency from 150 to 200 kHz enables peak amplitudes of the magnitude spectra to be obtained at different frequencies. Moreover, the local maximum values, in case of 200 kHz excitation, mostly correspond to the local minimum frequencies of 150 kHz excitation. Thus, excitation under the selected frequencies enables a large variety of reconstruction frequencies to be obtained. In this case, it was presumed that the selected excitation frequencies will provide a sufficient amount of velocity values. In other cases, more excitation frequencies may be used, exploiting the whole bandwidth of the transducer.
4a. To reconstruct the dispersion curve under the wide band, two different scenarios employing the square pulse excitation were used, as follows: n1 = 3 cycles, f1 = 150 kHz; and n2 = 3 cycles, f2 = 200 kHz. Such excitation frequencies were deliberately selected, according to the magnitude spectrum of excitation pulse, which can be seen in Figure 4b. The results, presented in the figure, demonstrate that a minor shift of excitation frequency from 150 to 200 kHz enables peak amplitudes of the magnitude spectra to be obtained at different frequencies. Moreover, the local maximum values, in case of 200 kHz excitation, mostly correspond to the local minimum frequencies of 150 kHz excitation. Thus, excitation under the selected frequencies enables a large variety of reconstruction frequencies to be obtained. In this case, it was presumed that the selected excitation frequencies will provide a sufficient amount of velocity values. In other cases, more excitation frequencies may be used, exploiting the whole bandwidth of the transducer.  The experimental waveforms of the S0 mode, at distances d1 and d2, under the f1 = 150 and f2 = 200 kHz excitation, are presented in Figure 5a,b, respectively. The magnitude spectra, |Ur 1 f 1 (jf)| and |Ur 2 f 2 (jf)|, of the windowed S0 mode wave packet can be seen in Figure 5c. The frequencies at which the phase velocity values were extracted are indicated with circle markers. Finally, the reconstructed dispersion curve of the phase velocity for the S0 mode, along with theoretical estimation, is shown on Figure 5d. The theoretical dispersion curve was calculated by employing the SAFE method and material properties of 4a. To reconstruct the dispersion curve under the wide band, two different scenarios employing the square pulse excitation were used, as follows: n1 = 3 cycles, f1 = 150 kHz; and n2 = 3 cycles, f2 = 200 kHz. Such excitation frequencies were deliberately selected, according to the magnitude spectrum of excitation pulse, which can be seen in Figure 4b. The results, presented in the figure, demonstrate that a minor shift of excitation frequency from 150 to 200 kHz enables peak amplitudes of the magnitude spectra to be obtained at different frequencies. Moreover, the local maximum values, in case of 200 kHz excitation, mostly correspond to the local minimum frequencies of 150 kHz excitation. Thus, excitation under the selected frequencies enables a large variety of reconstruction frequencies to be obtained. In this case, it was presumed that the selected excitation frequencies will provide a sufficient amount of velocity values. In other cases, more excitation frequencies may be used, exploiting the whole bandwidth of the transducer.  The experimental waveforms of the S0 mode, at distances d1 and d2, under the f1 = 150 and f2 = 200 kHz excitation, are presented in Figure 5a,b, respectively. The magnitude spectra, |Ur 1 f 1 (jf)| and |Ur 2 f 2 (jf)|, of the windowed S0 mode wave packet can be seen in Figure 5c. The frequencies at which the phase velocity values were extracted are indicated with circle markers. Finally, the reconstructed dispersion curve of the phase velocity for the S0 mode, along with theoretical estimation, is shown on Figure 5d. The theoretical dispersion curve was calculated by employing the SAFE method and material properties of  Figure 5a,b, respectively. The magnitude spectra, |U r1f1 (jf )| and |U r2f2 (jf )|, of the windowed S 0 mode wave packet can be seen in Figure 5c. The frequencies at which the phase velocity values were extracted are indicated with circle markers. Finally, the reconstructed dispersion curve of the phase velocity for the S 0 mode, along with theoretical estimation, is shown on Figure 5d The results in Figure 5d show that the phase velocities are reconstructed in the frequency band up to 0.8 MHz. According to the frequency response of the transducer used in this study (see Figure 4a), the technique enables the phase velocities in the −20 dB level bandwidth of the actuator to be reconstructed. In this study, a total of K = 52 velocity values were extracted at a band up to 1 MHz. This means that using two frequencies to drive the transducer, 52 reconstruction points were observed that correspond to peak frequencies of the magnitude spectra. Such a number of reconstruction points is relative and depends on the total number of excitation frequencies, N, and obtained number of peak values of magnitude spectra, in case of each excitation frequency.
It is noteworthy that the general reliability of the phase spectrum method depends on the proper selection of the time window to crop the wave packet of the single mode for FFT. The proposed method implicitly assumes that only one mode is present at the selected time window. aluminium 2024 T6 (the density: ρ = 2780 kg/m 3 ; Young's modulus: E = 72 GPa; Poisson's ratio: υ = 0.35).
The results in Figure 5d show that the phase velocities are reconstructed in the frequency band up to 0.8 MHz. According to the frequency response of the transducer used in this study (see Figure 4a), the technique enables the phase velocities in the −20 dB level bandwidth of the actuator to be reconstructed. In this study, a total of K = 52 velocity values were extracted at a band up to 1 MHz. This means that using two frequencies to drive the transducer, 52 reconstruction points were observed that correspond to peak frequencies of the magnitude spectra. Such a number of reconstruction points is relative and depends on the total number of excitation frequencies, N, and obtained number of peak values of magnitude spectra, in case of each excitation frequency.
It is noteworthy that the general reliability of the phase spectrum method depends on the proper selection of the time window to crop the wave packet of the single mode for FFT. The proposed method implicitly assumes that only one mode is present at the selected time window. In order to estimate the agreement of the results with theoretical phase velocities, the standard deviation (STD) was used as a measure of spread: In order to estimate the agreement of the results with theoretical phase velocities, the standard deviation (STD) was used as a measure of spread: where K  where K1 is number of points in reconstructed phase velocities, cp(fi) is a vector of reconstructed phase velocity values, and ct(fi) are the corresponding reference phase velocity values, calculated using the SAFE method. The estimated standard deviation of the calculated phase velocity values is σ = 161 m/s. This leads to the conclusion that 40 out of 52 velocity values (77%) are within the standard deviation range, as shown in Figure 6a. The experimental results, presented in this section, demonstrate that proposed approach reconstructs the phase velocity values at frequencies up to 800 kHz for the selected probe. At frequencies above 800 kHz, the approach starts to fail at capturing the pattern of the dispersion curve. Hence, it can be said that the phase velocity values of the S0 mode can be reconstructed at −20 dB bandwidth or 0.1 level of the transducer, according to its normalized magnitude spectra, presented at Figure 4a. The standard deviation of the reconstructed phase velocities, calculated according to Equation (7), is 161 m/s, which provides relative error of phase velocity estimation equal to ±3% for the S0 mode, calculated according to: where µc t (f) is the mean theoretical phase velocity value in the selected frequency band under analysis. In order to emphasize the achieved improvement, the signals of the S0 mode, obtained at 200 kHz, were processed using classic phased spectrum method, described in [33,34]. The reconstructed phase velocity curve is presented at Figure 6b. The results indicate that highest velocity reconstruction accuracy can be obtained at frequency band 200-340 kHz, which corresponds to −6 dB bandwidth of the sensor. The standard deviation of the S0 mode phase velocity reconstruction is estimated to be 592 m/s for the classic phase spectrum method, which gives ±11% relative phase velocity reconstruction error. In can be concluded that proposed approach allow to increase the reconstruction bandwidth, from 140 to 800 kHz, and reduce the relative velocity estimation error, from ±11% to ±3%, for the S0 mode. The experimental results, presented in this section, demonstrate that proposed approach reconstructs the phase velocity values at frequencies up to 800 kHz for the selected probe. At frequencies above 800 kHz, the approach starts to fail at capturing the pattern of the dispersion curve. Hence, it can be said that the phase velocity values of the S 0 mode can be reconstructed at −20 dB bandwidth or 0.1 level of the transducer, according to its normalized magnitude spectra, presented at Figure 4a. The standard deviation of the reconstructed phase velocities, calculated according to Equation (7), is 161 m/s, which provides relative error of phase velocity estimation equal to ±3% for the S 0 mode, calculated according to: where µ ct(f ) is the mean theoretical phase velocity value in the selected frequency band under analysis. In order to emphasize the achieved improvement, the signals of the S 0 mode, obtained at 200 kHz, were processed using classic phased spectrum method, described in [33,34]. The reconstructed phase velocity curve is presented at Figure 6b. The results indicate that highest velocity reconstruction accuracy can be obtained at frequency band 200-340 kHz, which corresponds to −6 dB bandwidth of the sensor. The standard deviation of the S 0 mode phase velocity reconstruction is estimated to be 592 m/s for the classic phase spectrum method, which gives ±11% relative phase velocity reconstruction error. In can be concluded that proposed approach allow to increase the reconstruction bandwidth, from 140 to 800 kHz, and reduce the relative velocity estimation error, from ±11% to ±3%, for the S 0 mode.

Identification of Converted Modes
In this section, the numerical validation of the proposed phase velocity reconstruction method will be presented. The major focus will be given to the method performance, in case the analysed signal is surrounded by the wave packets of other co-existing modes. To achieve the purpose of this study, the phase velocities of the converted A 0 mode will be analysed, which convert from the S 0 mode, due to the presence of notch.
To fulfil the scope of this research, the 3D linear structural mechanics finite element model of isotropic aluminium alloy 2024 T6 plate (600 × 200 × 2 mm) is considered. The top view of the analysed structure is presented on Figure 7. The S 0 mode was initially launched into the structure by applying the in-plane force to the shortest edge of the Al plate. To generate the converted A 0 mode, the vertical 36 mm wide (along x axis) crack-type defect, with a depth of 66% of the plate thickness, was introduced by duplicating the nodes of the mesh. In such way, a complete disbond was simulated, without changing the shape of finite element model. It was shown by the various researchers that, if a crack is not symmetrical to the middle plane of the plate, according to the thickness, the mode conversion takes place upon the wave interaction with the notch, and both the S 0 and A 0 modes are expected as the reflected and transmitted waves [37]. The defect was centred, with respect to the short edge of the sample, and situated at the distance of 200 mm from source of Lamb waves (see Figure 7). be analysed, which convert from the S0 mode, due to the presence of notch.
To fulfil the scope of this research, the 3D linear structural mechanics finite element model of isotropic aluminium alloy 2024 T6 plate (600 × 200 × 2 mm) is considered. The top view of the analysed structure is presented on Figure 7. The S0 mode was initially launched into the structure by applying the in-plane force to the shortest edge of the Al plate. To generate the converted A0 mode, the vertical 36 mm wide (along x axis) cracktype defect, with a depth of 66% of the plate thickness, was introduced by duplicating the nodes of the mesh. In such way, a complete disbond was simulated, without changing the shape of finite element model. It was shown by the various researchers that, if a crack is not symmetrical to the middle plane of the plate, according to the thickness, the mode conversion takes place upon the wave interaction with the notch, and both the S0 and A0 modes are expected as the reflected and transmitted waves [37]. The defect was centred, with respect to the short edge of the sample, and situated at the distance of 200 mm from source of Lamb waves (see Figure 7). Throughout the simulations, the ANSYS 17.1 implicit solver and 3D structural solid solid64 finite elements were used, which are defined by eight nodes having three degrees of freedom at each node and 2 × 2 × 2 integration points. The finite elements were hexahedrons, meshed using structured grid. Once again, two different scenarios employing the square pulse excitation were used, as it was described in the previous section. At first, the excitation pulse consisted of n1 = 3 cycles and a central frequency of f1 = 150 kHz. Meanwhile, in the second case, the Lamb waves were excited with n2 = 3 cycles at f2 = 200 kHz. The average mesh size was equal to 0.5 mm, which corresponds to 21 nodes per wavelength for the slowest A0 mode at f1 and 17 nodes per wavelength at f2. The integration steps in the time domain were 0.33 and 0.25 µs, respectively, which produces a 1/20 of the period, both at f1 and at f2. The variable monitored in this study was a vertical component of particle velocity (y) along the centreline of the sample. The waveforms for the phase velocity estimation were selected along the centreline of the sample at distances d1 = 240 mm and d2 = 360 mm. The B-scan images of the longitudinal (z) and vertical component (y) of the particle velocity, showing the S0 and converted A0 modes, are presented in Figure  8a,b. Throughout the simulations, the ANSYS 17.1 implicit solver and 3D structural solid solid64 finite elements were used, which are defined by eight nodes having three degrees of freedom at each node and 2 × 2 × 2 integration points. The finite elements were hexahedrons, meshed using structured grid. Once again, two different scenarios employing the square pulse excitation were used, as it was described in the previous section. At first, the excitation pulse consisted of n 1 = 3 cycles and a central frequency of f 1 = 150 kHz. Meanwhile, in the second case, the Lamb waves were excited with n 2 = 3 cycles at f 2 = 200 kHz. The average mesh size was equal to 0.5 mm, which corresponds to 21 nodes per wavelength for the slowest A 0 mode at f 1 and 17 nodes per wavelength at f 2 . The integration steps in the time domain were 0.33 and 0.25 µs, respectively, which produces a 1/20 of the period, both at f 1 and at f 2 . The variable monitored in this study was a vertical component of particle velocity (y) along the centreline of the sample. The waveforms for the phase velocity estimation were selected along the centreline of the sample at distances d 1 = 240 mm and d 2 = 360 mm. The B-scan images of the longitudinal (z) and vertical component (y) of the particle velocity, showing the S 0 and converted A 0 modes, are presented in Figure 8a,b.
The simulated waveforms of the converted A 0 modes, at distances d 1 and d 2 , in case of f 1 = 150 kHz and f 2 = 200 kHz excitation, are presented in Figure 9a,b. The selected time windows to cut the wave packet of single mode are indicated with vertical dashed lines. The magnitude spectra of windowed A 0 mode, at frequencies f 1 and f 2 , along with indicated reconstruction frequencies, can be seen on Figure 9c. Finally, the comparison of estimated DC, with the theoretical calculations, is shown on Figure 9d. The results demonstrate a good match between the estimated results and theoretical phase velocities, calculated with the SAFE method. The standard deviation of the reconstructed velocities is equal to σ = 47.3 m/s. Overall, the K = 32 velocity values were extracted, while 20 (63%) of them were within the range of standard deviation. Even though the number of reconstruction points is less than from the experiments present in previous section, Figure 9d suggests that its quite sufficient for the reconstruction of the segment of dispersion curve. The proposed approach is not limited with two excitation frequencies; hence, the number of reconstruction points can be increased if the segment of dispersion curve is not represented properly. The simulated waveforms of the converted A0 modes, at distances d1 and d2, in case of f1 = 150 kHz and f2 = 200 kHz excitation, are presented in Figure 9a,b. The selected time windows to cut the wave packet of single mode are indicated with vertical dashed lines. The magnitude spectra of windowed A0 mode, at frequencies f1 and f2, along with indicated reconstruction frequencies, can be seen on Figure 9c. Finally, the comparison of estimated DC, with the theoretical calculations, is shown on Figure 9d. The results demonstrate a good match between the estimated results and theoretical phase velocities, calculated with the SAFE method. The standard deviation of the reconstructed velocities is equal to σ = 47.3 m/s. Overall, the K = 32 velocity values were extracted, while 20 (63%) of them were within the range of standard deviation. Even though the number of reconstruction points is less than from the experiments present in previous section, Figure 9d suggests that its quite sufficient for the reconstruction of the segment of dispersion curve. The proposed approach is not limited with two excitation frequencies; hence, the number of reconstruction points can be increased if the segment of dispersion curve is not represented properly.
As it was mentioned previously, the time window selection ambiguity is quite essential in the success of phase velocity reconstruction using phase-shift method, especially for overlapped modes. Hence, the selected time window must hold the single mode only. For complex structures, where signals undergo many reflections, the reconstruction can be quite uncertain. On the other hand, it will be demonstrated in the next chapter that phase velocities can be estimated using part of the signal only. In such a case, the position of the time window must be optimised, i.e., by solving a minimisation problem, to get reasonable velocity reconstruction results.

Analysis of Multimodal Signals in Anisotropic Structures
In this section, the performance of the proposed phase velocity reconstruction approach is validated qualitatively by analysing the experimental multimodal signals in an anisotropic structure. For this purpose, the experiments have been carried out in a pitchcatch configuration on the 6-ply GFRP plate (biaxial: 0° and 90°/bias: ±45°/biaxial: 0° and 90°), with dimensions xo = 2000, yo = 1000, and 4 mm thickness (see Figure 10). As it was mentioned previously, the time window selection ambiguity is quite essential in the success of phase velocity reconstruction using phase-shift method, especially for overlapped modes. Hence, the selected time window must hold the single mode only. For complex structures, where signals undergo many reflections, the reconstruction can be quite uncertain. On the other hand, it will be demonstrated in the next chapter that phase velocities can be estimated using part of the signal only. In such a case, the position of the time window must be optimised, i.e., by solving a minimisation problem, to get reasonable velocity reconstruction results.

Analysis of Multimodal Signals in Anisotropic Structures
In this section, the performance of the proposed phase velocity reconstruction approach is validated qualitatively by analysing the experimental multimodal signals in an anisotropic structure. For this purpose, the experiments have been carried out in a pitch-catch configuration on the 6-ply GFRP plate (biaxial: 0 • and 90 • /bias: ±45 • /biaxial: 0 • and 90 • ), with dimensions x o = 2000, y o = 1000, and 4 mm thickness (see Figure 10).
(c) (d) Figure 9. The simulated waveforms of the converted A0 mode, at distances d1 and d2, in case of (a) 150 and (b) 200 kHz excitation; (c) the magnitude spectra of windowed A0 mode at different excitation frequencies; (d) the combined reconstruction of phase velocity dispersion curve along with the theoretical estimation.

Analysis of Multimodal Signals in Anisotropic Structures
In this section, the performance of the proposed phase velocity reconstruction approach is validated qualitatively by analysing the experimental multimodal signals in an anisotropic structure. For this purpose, the experiments have been carried out in a pitchcatch configuration on the 6-ply GFRP plate (biaxial: 0° and 90°/bias: ±45°/biaxial: 0° and 90°), with dimensions xo = 2000, yo = 1000, and 4 mm thickness (see Figure 10). The Lamb waves were generated using the MFC transducer, centred at the coordinates x e = 500 mm, y e = 250 mm. It was bonded to the surface of the specimen using a thin layer of gasket maker. The emitter was excited by a three-cycle square pulse, with a central frequency of 100 kHz, where the fundamental A 0 and S 0 modes exist in the structure. In this case, the measurements were recorded at a single excitation frequency. Two waveforms were recorded along the wave path (0 • propagation), at the distances d 1 = 773 mm and d 2 = 895 mm from the source of Lamb waves (see Figure 10). The proposed phase velocity estimation method was used to extract velocities of the four wave packets: direct S 0 , bottom reflected S 0 , left top edge reflected S 0 , and direct A 0 mode. The experimentally obtained waveforms, at the distances d 1 and d 2 , are presented in Figure 11a Table 1. Note that the reconstructed velocity values below 40 kHz were not considered in the calculations of STD. forms were recorded along the wave path (0° propagation), at the distances d1 = 773 mm and d2 = 895 mm from the source of Lamb waves (see Figure 10). The proposed phase velocity estimation method was used to extract velocities of the four wave packets: direct S0, bottom reflected S0, left top edge reflected S0, and direct A0 mode. The experimentally obtained waveforms, at the distances d1 and d2, are presented in Figure 11a The reconstructed phase velocities of different reflections can be seen in Figure 12ad. The standard deviations for each case of reconstruction are summarized in Table 1. Note that the reconstructed velocity values below 40kHz were not considered in the calculations of STD. obtained waveforms, at the distances d1 and d2, are presented in Figure 11a,b. The start and stop points of the time windows used to crop the wave packets are indicated by dashed squares.  The results presented above ( Figure 12) were found to be in quite good agreement The results presented above ( Figure 12) were found to be in quite good agreement with the theoretical calculations. It suggests that the proposed technique can be used with a certain reliability to extract the phase velocities of GW and identify modes in complex signals. The results show that the velocities of the direct modes are closer to the theoretical values, in comparison to the reflected ones. The average deviation for the direct modes (A 0 and S 0 ) is approximately 75 m/s, while for the reflected S 0 modes, it is 213 m/s. Several factors may influence the reliability of the results, though. First of all, the selected time windows in Figure 11a,b (dashed squares) may give an idea that this procedure is not very straightforward, especially for the reflected modes. As it turns out, in some cases, part of the wave packet has to be cropped to get better velocity estimation. Another important factor is the propagation distance, which varies for modes arriving at different directions. It means that the distance (d) has to be predefined for each wave packet separately. If the propagation distance is not known in advance, an additional error will be obtained. The study revealed that the proposed velocity estimation technique gives an approximate relative experimental error of ±4%, in comparison to theoretical predictions. Meanwhile, for the incident modes, the relative error is always less than a ±2.5%. For example, the 2D FFT method gives an error of approximately of 1% [38]. However, in the study above, the authors used a set of 64 time series, spatially sampled at 1 mm, to achieve such accuracy.

Conclusions
In this paper, a phase velocity reconstruction approach, based on the phased spectrum method, was developed, which exploits several excitation frequencies of the ultrasonic probe and estimates phase velocity values at peak frequencies of the magnitude spectra. The proposed approach allows us to reconstruct phase velocities with high accuracy in wide frequency bandwidth using only two waveforms measured at close proximity. In contrast to the classic phase spectrum method, the proposed technique offers an increased reconstruction bandwidth (from −6 to −20 dB) and reduced relative error of phase velocity reconstruction (from ±11% to ±4%). The main outcomes of the research can be summarized as follows: • It was found that the accuracy of the classic phase spectrum method can be improved if several frequencies are used to drive the transducer, while the phase velocities are reconstructed at peak values of Fourier spectra only. Such an approach allows us to avoid low energy frequency components, where the velocity estimation error is likely to increase.

•
The initial experiments demonstrated that the proposed phased spectrum method can increase the reconstruction bandwidth, from −6 to −20 dB, of the sensor and improve the standard deviation of velocity reconstruction, from 592 to 161 m/s. For the experimental S 0 mode, this results in a velocity estimation relative error improvement, from ±11% to ±3%.

•
The finite element simulations demonstrated the applicability of the proposed approach in detecting converted guided wave modes. It was demonstrated that the phase velocities of converted modes can be reconstructed with a standard deviation of 47.3 m/s, even if the modes are partly overlapped with direct waves. • Finally, the proposed method was demonstrated to be appropriate for the analysis of complex guided wave signals, with multiple co-existing modes. It was estimated that average deviation for the direct modes (A 0 and S 0 ) is approximately 75 m/s, while for the reflected S 0 modes, it is 213 m/s. While analysing the overlapped complex guided wave signals, the proper selection of time gate is the most important parameter for the accuracy of reconstruction. • It was estimated that the average phase velocity reconstruction error of the proposed method, including both symmetrical and asymmetrical modes, is up to ±4%. The classic phase spectrum method provides an approximate reconstruction error of ±11%.
Other techniques, reported in the literature, can achieve velocity reconstruction with an average error of ±1%; however, at least 64 signals need to be acquired to achieve such accuracy.