A novel method for impact force estimation in composite plates under simulated environmental and operational conditions

During its lifetime, an aircraft structure is subjected to various impacts from various sources such as tool drops, hail, ground service equipment, etc. In modern composite structures, these impacts have a significant chance of generating barely visible damage (BVID) which may lead to catastrophic failure of a structure if left undetected to grow. However, BVID is difficult to detect during routine visual inspection without specialised non-destructive inspection and thus there is large interest in developing monitoring systems for estimating the location and severity of impact events. Currently, most systems and methods have been developed for controlled lab conditions and do not consider the wide range of impact parameters in real life operation (environmental conditions, vibration, impactor stiffness, angle, etc) which may severely compromise the accuracy of these methods. In this study we have explored two methods for maximum impact force estimation, deconvolution and a novel gradient method, for the purpose of reliable severity assessment in composite aircraft structures under simulated environmental and operational conditions. It is shown that both methods allow accurate and robust estimation of the maximum impact force from various cases of impacts (variation of impact energy, mass, stiffness, angle, temperature, source) using minimum initial data from a single impact case. From further testing it is demonstrated that the gradient method is robust towards the effects of impact localisation errors and noise. The gradient method also has much less computational and storage requirements and is thus more feasible to integrate with current data acquisition systems being developed for structural health monitoring. Thus, we conclude that the proposed gradient method is suitable for impact force monitoring and severity assessment in composite aircraft structures in the simulated environmental and operational conditions.


Introduction
During its lifetime, an aircraft structure is subjected to various impacts from various sources such as tool drops, hail, ground service equipment, etc [1]. In modern composite structures, these impacts have a significant chance of generating barely visible damage (BVID) which may lead to catastrophic failure of a structure if left undetected [1]. However, BVID is difficult to detect during routine visual inspection without specialised non-destructive inspection and thus there is large interest in developing monitoring systems for impact damage [2,3]. Currently, most systems and methods have been developed for controlled lab conditions and do not consider the wide range of impact parameters in real life operation (vibration, impactor stiffness, mass, temperature, angle, etc) which may severely compromise the accuracy of these methods [4][5][6][7].
It is predicted that impact damage initiates at a specific threshold force [1], thus one of the main interests of impact monitoring systems is to estimate the impact force to determine the severity of the impact (estimate possibility of damage) [5,6,[8][9][10]. Many methods have been developed using passive lamb wave measurements (lamb wave source is the impact event and not active actuation from a transducer) from sensors at discrete points to reconstruct the impact force profile [5,6,[8][9][10][11][12][13]. Analytical solutions to such inverse problems exist but may be impractical to apply for complex structures [14,15].
A popular alternative approach is to find the function which maps the recorded lamb wave response directly onto the impact force for a specific location [5,[8][9][10][11]16]. This way the actual geometric complexities of the structure can be disregarded, although there is now a requirement to store characteristic functions of each location [8,17]. Many methods have been developed to find these characteristic functions, with transfer function based methods (deconvolution [7,8,12,18,19] and state space models [13,16,20]) and data driven [5,9,11,21] based works being the most common in recent studies.
Data driven based methods can directly map the recorded lamb wave to the impact force via a meta model trained with a known training data set [5,9,11,21]. They have the flexibility to map most forms of input and outputs, but the drawback is that extrapolation for inputs that are not included in the training set yield lower accuracy [4,5,22]. This is vital for application in real life applications as the variation of possible impact cases may require an impractical amount of training data to achieve robust accuracy [4].
Deconvolution and state space models use the transfer function between the recorded signal and the impact force profile [6-8, 12, 16-19]. These methods allow reconstruction of any form of impact profile using the transfer function obtained from a single impact force and lamb wave pair and are thus able to accommodate the many variations in possible impact cases. However, they can be unstable especially for noisy signals [6,7] and many regularisation methods have been developed to give more robust results [6,7,10,19]. Additionally, these methods may require large amounts of storage for transfer functions of multiple locations.
A few studies have discovered a simpler relationship between signal features (energy, amplitude, etc) with impact severity (force, energy, etc) [12,23,24] which have the potential to be used as simplified severity indicators. It was found that there is a linear relationship between signal energy and impact energy [24] that was affected by impactor stiffness. A linear relationship was also found between signal amplitude and maximum impact force [12,23] with varying gradient depending on location.
Development of methods to identify impact force need to follow the capabilities of the data acquisition (DAQ) systems being developed for operation on aircraft structures [2,3,25]. These systems tend to be simple to minimise weight and energy required for computing especially if they are to be deployed on the entire structure. In this paper a novel impact force identification method is proposed that is robust for simulated environmental and operational conditions experienced by an aircraft structure which is feasible for integration with current DAQ systems and requires minimum initial data. The method expands on previous linear relations found between the impact signal magnitude and the impact force to create a robust and simplified estimation method. We use the more established deconvolution method as a benchmark to assess the performance of the proposed method.

Experimental setup
Experimental measurements were acquired from two different setups: i) a low energy non-damaging impact test setup for collecting large amounts of impact data under various simulated environmental and operational conditions, and, ii) a high energy impact test setup for collecting data on high energy damaging impacts. Data from both setups are then used to develop and test impact force estimation methods that can identify damaging impact events as well as being robust towards simulated environmental and operational conditions. All data processing was done in MathWorks MATLAB 2018a.

Non-damaging impact setup
The non-damaging impact setup consists of a flat carbon fibre plate (M21 T800s, [0/ + 45/-45/90/0/ + 45/-45/90]s layup) and two types of impactors: an impact hammer and a small drop impactor setup as can be seen in figure 1. Both impactors were used to generate impacts from differing sources (manually actuated vs dropped mass) as well as random (hammer) vs controlled parametric impacts (drop impactor). The drop impactor, however, has no apparatus to measure the impact force directly and as such the impact force was estimated using methods discussed later in section 4.
The impact hammer (PCB Piezotronics 086C03, figure  1(b)) can measure the impact force at the tip and was used to generate impacts with random force. It has interchangeable tips (10 mm diameter) with different materials (metal, plastic, hard rubber and soft rubber) and as such can simulate impacts from different stiffness impactors. On the other side of the hammer head there is an attachment which allows the addition of a 100 g extra mass for simulating variation in impactor mass.
The drop impactor setup consists of a rail and the impactor itself which can generate repeatable impacts of specific energy by dropping the impactor from a controlled height. The impactor has two tips (steel and soft rubber, 20 mm diameter) to simulate impacts from impactors of different stiffness. The impactor itself weighs 100 g, with attachments to add extra weights to a maximum of 200 g. The rail that guides the impactor can be tilted to generate impacts with controlled angles.
The carbon fibre plate dimensions were 290 × 200 × 4 mm, had a silicone heating mat underneath to control the temperature and steel clamps on two sides to secure the specimen. The layout of the six PZT sensors (only six out of eight used as oscilloscope only had eight channels and one channel was required to measure the impact force as mentioned later on, figure 1) and impact locations are shown in figure 1(c). The PZT sensors convert the lamb waves generated by impact into a voltage signal which were recorded using an NI PXI 5105 8 channel oscilloscope at 2 MS s −1 with a signal length of 100 000 samples. Each sensor was connected to the oscilloscope via probes with 10x attenuation. Table 1 shows the various impact cases collected using the non-damaging impact setup using the hammer and impactor which encompasses a wide variety of possible cases in real life. For all cases, impacts at each location (shown in the layout, figure 1(c)) were repeated four times. All the impacts generated using this setup are non-damaging.

Damaging impact setup
The damaging impact setup (figure 2) was used to generate impacts of much larger energy which may create barely visible damage on composite specimens. The setup consists of composite plate specimens and an impact tower. The impact tower (INSTRON CEAST 9350) has a drop impactor with a mass of 2.41 kg and a steel tip with 20 mm diameter. The impactor was dropped from a specific height to generate impacts with energies of 1, 7, 9 and 11 J resulting in four impact tests. The 1 J impact acts as the non-damaging reference impact whilst the rest are the damaging impact cases. The impact force was measured from the force sensor on the impactor tip which was recorded at 500 KS s −1 and a length of 25 000 samples.
For each impact test, an individual composite plate (8552-33%-134-IM7-12 K, [0/ + 45/-45/90] 2s layup) with dimensions 250 × 250 × 2 mm was used (4 in total for all tests). The specimens were held in place by a steel fixture which clamps all four sides of the plates (figure 2(b)). Each plate has one PZT sensor (as can be seen in the layout in figure 2(c)) to measure the impact generated lamb wave which was recorded using a GW Instek GDS-2074A oscilloscope with 10x attenuation at 500 KS s −1 and a length of 25 000 samples. All plates were impacted at the same location (figure 2(c)) and after impact an ultrasound scan using a CF08 Dolphicam was done to ascertain the extent of damage.

Signal deconvolution for impact force profile reconstruction
In signal deconvolution, it is assumed that the measured response of a system (in our case the measured lamb wave) is equal to the input signal (in our case the impact force) that is convoluted as it propagates through the system [8,9,26], as can be seen in equation (1) where s is the recorded time (t) series response, c is the convoluting function and f is the   original input signal.
In order to obtain an unknown original input (f ) from a measured response (s) we need to deconvolute the input signal by finding the inverse of the convolution [8,9,26]. As the convolution is a time series integral, it is often transformed into the Fourier domain for ease. This transforms the integral into multiplication of spectral components [8,9] as seen in equation (2) where S, C and F are the spectral components of s, c and f respectively.
In this form, the inverse of the spectral convolution (C) can be expressed in terms of the transfer function (H) as seen in equation (2). If the transfer function is known, we may find the original input signal (f ) from the measured response (s) through multiplication of the spectral components of the response (S) with the transfer function (H) followed by an inverse Fourier transform. For deconvoluting impact force profiles from measured lamb waves on a structure with a specific configuration (geometry, material and boundary conditions), the transfer function depends on the location as it is affected by the propagation path of the wave source to the sensors [8]. The transfer function may be obtained experimentally by measuring the force profile and subsequent response of an impact at each location [8,18] or through estimation [26,27]. Once the transfer function is known for a specific location, impact profiles of any form (different impactor energies, masses, materials, angles) may be reconstructed from the measured response provided that the structure configuration is constant (thus constant transfer function) and the bandwidth of the input profile is covered by the transfer function [18,19]. This property makes deconvolution a promising method for robust reconstruction impact force profiles under varying operational and environmental conditions.

Reconstruction of impact force profiles under simulated environmental and operational conditions
Here we tested the deconvolution method for reconstructing the impact force profile measured from the non-damaging impact setup (figure 1) with the impact hammer under simulated environmental and operational conditions as listed previously in table 1 (H1-H6, except H1a as locations differ). We used the first impact repetition from each location (1 × 35) of the reference case (H1) to find the transfer function for these locations (1-35, figures 1(c) and 3). The rest of the impact repetitions (3 × 35) of the reference case (H1) as well as the other impact cases (H2-H6) were then used as other test cases to see the robustness of the obtained transfer functions.
As mentioned previously, deconvolution results can be unstable [10,19] as the scaling factors of the transfer function may vary largely. It is important that the frequency components are sampled with enough resolution to prevent scaling up the wrong component, resulting in noise or reconstruction errors. To achieve this, we down sampled the signal by a factor of 4 (by keeping only 1 in every 4 samples) as it was highly oversampled and used a 10 6 point Fourier transform to ensure the spectrum contains less noise (only frequency components of the signal) and has sufficient resolution. As there were six sensors on the plate, we obtained the transfer function for each sensor and reconstructed the impact force profile individually ( figure 3). Afterwards, all six reconstructed impact force profiles were averaged to obtain a robust final result. Figure 3 shows the impact force profile reconstruction results for case H2 (added mass) and H6 (soft rubber tip) at location 1 using the transfer function obtained from case H1. It can be seen that the deconvolution results match very closely to the measured force profiles even though they have significantly different shapes and values. For severity estimation, as we are mostly interested to know whether the damage threshold force is passed, we summarise the maximum impact force reconstruction results for all the test cases (H1-H6) in table 2. It can be seen that the estimated maximum force using the single transfer function obtained from case H1 has good accuracy when compared to the measured force for all cases except H3 (increased temperature).
We know from previous studies that increased temperature decreases the amplitude of the recorded response [22,28,29]. In table 2 (case H3a) and figure 4 we can see that this temperature effect causes the maximum force to be underestimated when not compensated. Using the impactor ( figure 1(a)), we generated controlled impacts at different temperatures to see the pattern of amplitude reduction. The results can be seen in figure 5 where it shows a linear trend of reduction [28] from the reference room temperature (24 • C). Using this linear trend, we obtained the correction factor (figure 5) to compensate for the reduction in amplitude which improved the maximum impact force estimation accuracy (case H3b, table 2 and figure  6) to similar levels as other cases.
From these results we conclude that for the simulated environmental and operational cases tested the deconvolution method can give a robust and accurate reconstruction of the impact force profile using transfer functions obtained from a single impact case (H1). However, the data storage (transfer functions of all locations) and computational requirements (Fourier transforms, etc) may not necessarily be suited for the capability of current data acquisition systems [2,3,25], especially if applied to the entire structure. As we are mostly concerned with the maximum impact force for severity assessment and not necessarily the whole impact force profile, most of the processing effort in deconvolution is not required. Thus, in the next section a simple yet robust novel method is proposed for estimating the maximum impact force under simulated environmental and operational conditions.

Novel gradient method for maximum impact force estimation
Previous studies [12,23,24] have found simple linear relationships between impact severity (energy or force) and impact signal features (energy, amplitude). Here, we build on these findings and establish if there is a consistent relation between the impact force and signal features which can be used to accurately estimate the maximum impact force under the simulated environmental and operational conditions. To accomplish this, we used the impacts collected on the non-damaging impact setup as shown in table 1.
We have established in section 3.1 that the deconvolution method is able to give an accurate reconstruction of the maximum impact force regardless of the original impact force profile. Thus, we used the obtained transfer functions to estimate the impact force from the impactor cases (case I1-I6, table 1) as the impactor does not have a direct force measurement system. Along with the hammer impacts (H1-H6), this adds a comparison between two types of impact actuation sources (manually actuated and free drop). It also gives a data set of controlled parametric variations (angle, energy, etc) not possible with the hammer.    Figure 7 shows the relation between the maximum absolute (MA) signal amplitude and maximum impact force for the specific specimen configuration (geometry, material and boundary condition) at a single location for all impact cases (H1-H6 and I1-I6). It can be seen that there is indeed a consistent linear relationship between the MA signal amplitude   and the maximum impact force for most cases (except H5, H6 and I6). Using this relationship, we can estimate the maximum impact force (Fmax) from the MA signal amplitude (A) provided the gradient (G) of this relationship is known, as shown in equation (3). It is known that the gradient changes with location (x,y) [12] and thus it must be obtained for all locations. Similar to the transfer functions for deconvolution, the gradient may be obtained using a known pair of signal and impact force measurements.
From figure 7, it can be seen that for impact cases using the rubber (soft and hard) impact tip (H5, H6 and I6), the gradient of the linear relationship is significantly different to that of the other cases with stiffer impact tips (steel and plastic). Other studies have found a similar phenomenon when correlating signal energy and impact energy of differing stiffness impactors [24]. The difference in gradient deters us from being able to estimate the maximum impact force accurately for a large range of impactor stiffness (steel~200 GPa, plastic 10 GPa, silicone rubber~0.1 GPa [30],). For aircraft structures, impacts may come from various sources such as tool drop (steel), hail (stiffness similar to plastic) and ground service equipment (with rubber bumpers) [1]. Thus, it may be necessary to find an alternative signal feature other than the MA signal amplitude that is consistent for a wider range of impactor stiffness.

Gradient method for impact force estimation in different stiffness impactors
The response of a structure (y) to transient loading (impact, etc) can be expressed as the components of the response of a non-homogeneous differential damped mass spring system which consists of the particular (yp) and general (yg) response [1,31] as shown in equation (4). The particular response is the transient response of the structure (with mass m, stiffness k, and damping c, dependant on the structures geometry, material and boundary conditions) towards the particular transient load (f ) along its bandwidth (ωf ). The general response is the response of the system that is mostly independent of the transient load (only implicitly affected through the change of state, s, due to the particular response) and only occurs at the structure specific (based on configuration) natural frequencies (ωn). The components of the particular and general responses that are near to the structures natural frequencies are scaled up due to resonance [32]. In a lamb wave signal from impact, the particular response is most notable when the impact load is applied causing the structure to undergo forced vibration (figure 8). When there is no more load, the structure freely vibrates at its natural frequencies (the general response) until the kinetic energy is dissipated through damping [1].
For second order systems such as simple vibrating structures, it is known that the particular response (yp) amplitude is a linear function of the force (f ) magnitude independent from the force profile [1]. Similarly, in figure 7 we find that the recorded MA signal (y) amplitude scales linearly with the impact force (f ) magnitude for all cases but show a difference in gradient (G) between stiffer and less stiff impactors. If the particular response (yp) amplitude has a consistent linear relation to the impact force (f ) magnitude irrespective of the input force profile, then the variable that may cause a change in gradient (G) is the general response (yg) which interacts implicitly with the particular response (yp).
It is known that impacts generate a wideband impact force profile (and subsequent impact signal), with a lower bound of 0 Hz and an upper bound dependent on the impactor stiffness [18]. As can be seen in figure 3, stiffer impactors trigger higher frequency responses (small peaks and short contact time) compared to less stiff impactors (smoother response with longer contact) [18]. The particular response of the structure (yp) will follow the bandwidth of the impact force profile (ωf ) whilst the general response (yg) will only occur at the frequency components that match the natural frequencies (ωn) of the structure. For stiff impacts, as the bandwidth of the impact forces (ωf ) are comparable, the generated general response (yg) profile will also have similar bandwidths. As such, it will interact with the particular response (yp) in a similar way, causing comparable gradients (G) between the amplitude of the response (y) and the magnitude of the impact force (f ). However, for less stiff impactors, the impact force profile (f ) does not trigger the higher natural frequency modes (only the 1st mode as the main flexural mode) and thus the bandwidth of the general response (yg) differs from stiffer impactors. This causes a difference in the reaction between the particular (yp) and general (yp) response as compared to stiffer impactors which ultimately leads to differing gradients (G). This is why the problem of different impactor stiffnesses is not encountered in deconvolution (as seen in section 3.1) as the transfer function captures the complex interaction between all the frequency components.
Thus, in order to have a consistent gradient between differing stiffness impactors, we need to remove the general response (yg) from the total response (y). A crude but simple way of accomplishing this is by filtering out the natural frequencies (ωn) as we know the general response (yg) only occurs at these frequency bands. However, this will also filter out components of the particular response (yp) at the same frequency bands. As we will show later on in section 4.2, for our purpose this turns out to be an acceptable loss in terms of accuracy of results.
To identify the natural frequency bands of the flat plate used in the non-damaging impact setup, we looked at the transfer functions obtained previously in section 3.1. Figure 9 shows the average magnitude for the transfer functions for all sensors due to impacts at all 35 locations of the flat plate (figure 2). The natural frequency bands were identified as the peaks in the transfer function signifying resonance [32]. Using these frequency bands, the signal was filtered (using IIR Butterworth bandpass filters) from the steel and soft rubber impact hammer cases (H1 and H6) as a test case (shown in figure 10). The 1st natural frequency band (filtering starts from the 2nd band onwards) was not filtered out as this was the main flexural mode for all impacts regardless of stiffness. Figure 10 shows that as more frequency modes are filtered, the MA signal amplitude of the steel hammer signal (H1) decreases approaching the value of the soft rubber hammer signal (H6). This happens as the filtering removes the high frequency components only present in the stiffer impactor signal (H1). The result is that the force gradient (G) of the steel   hammer impact (H1) becomes much closer to that of the soft rubber hammer impact (H6), thus removing the differences in gradient ( figure 10). However, as mentioned previously the method of filtering also removes components of the desired particular response (yp) which is why the gradient (G) values do not match up completely. By applying the filtering to the impacts shown in figure 7, we get a new linear relationship between the filtered MA signal amplitude and the maximum Figure 11. Linear relationship between filtered maximum absolute (FMA) signal amplitude and maximum impact force for impacts tested on the non-damaging impact setup. impact force which is now consistent for all cases as shown in figure 11. Thus, by substituting the MA signal amplitude in equation (3) with the filtered maximum absolute (FMA) signal amplitude we are able to estimate the maximum impact force from impacts of all cases tested (table 1).

Gradient method for impact force estimation in simulated environmental and operational conditions
Here the accuracy of the gradient method in estimating the maximum impact force under the environmental and operational conditions simulated in the impact cases taken from the non-damaging impact setup (H1-H6, I1-I6, table 1) is investigated. The first repetition of the reference impactor case (I1) was used as a reference to obtain the maximum impact force, MA signal amplitude and FMA signal amplitude for all impact locations (figure 1). From these values the reference impact force gradient for MA and FMA signal amplitudes were obtained and used to estimate the maximum impact force of all other cases collected (H1-H6, I1-I6). As there were six sensors on the flat plate, we generated an estimate from each sensor signal which is then averaged to obtain a robust final estimation. For impact cases with increased temperature (H3a, H3b, I5a and I5b), we used the correction factor obtained in section 3.1 on the signal before extracting the MA and FMA signal amplitude values. Table 3 shows the results of the maximum impact force estimation using both the MA and FMA signal amplitude accompanied by a comparison with the deconvolution method established in section 3. It can be seen that the deconvolution method in general has a higher accuracy compared to both gradient methods (MA and FMA) even though the difference is not too large. As postulated in section 4.1, the gradient  method using the MA signal amplitude is able to give an accurate estimate of the maximum impact force for all cases except for the non-stiff impactors (H4-H6, I6; table 3, figure 12) where it is severely underestimated. With the FMA signal amplitude, the gradient method is able to accurately estimate the maximum impact force of non-stiff impactors (table 3, figure 13) as well as retaining a similar level of accuracy to the MA signal amplitude for other cases. The removal of components of the particular response due to filtering mentioned in section 4.1 does not seem to significantly reduce the accuracy of estimation. This is likely due to the final estimation being an average from all sensors thus being more robust. We conclude that the gradient method with the FMA signal amplitude is able to accurately and robustly estimate the maximum impact force for impact cases under the simulated environmental and operational conditions using gradient values obtained from a single impact case.
Compared to the deconvolution method, the data storage requirement of the gradient method is much lower as we only need to store a single gradient value (as opposed to a transfer function) for each location. The process of filtering out natural frequencies may potentially be conducted using analog circuitry [2,3,25] once the frequency bands have been identified. Moreover, the calculation involved in the estimation of the maximum impact force is much simpler (only multiplication) compared to the deconvolution method which involves large matrix operations and Fourier transforms. Thus, we can conclude that the gradient method is more feasible to integrate with current DAQ systems.

Effect of localisation and noise towards gradient method force estimation
As the characteristics used for impact force estimation (transfer functions or force gradients) are location dependent, the monitoring of impact force with these methods is thus coupled and preceded by localisation of the impact event. Many methods have been developed for such purpose including data driven methods which use signal features such as time of arrival and amplitude in a reference database (for pattern matching or neural networks) as used here for impact force estimation [22,29].
The maximum impact force estimation results from the deconvolution and gradient methods in sections 3 and 4 were based on the condition where the exact impact location is known and it coincides with the location where the reference characteristics (gradient or transfer function) was measured. In reality, most cases will not follow those conditions because of two reasons: (1) we cannot feasibly measure the reference characteristics at all possible locations for a continuous structure and (2) all impact localisation methods have some degree of error such that the reference characteristic used may not be entirely correct. Additionally, there may be noise in the form of vibrations [22,29] which may interfere with the signal amplitude. Here we test the robustness of the proposed gradient method towards these effects using impacts collected on the non-damaging setup.
To address the issue of reference characteristic collection on all locations, most studies using reference databases for either localisation or impact force estimation commonly use some form of interpolation between sampled reference locations to get a continuous map of values [8,17,29]. This is based on the assumption that the characteristics of a structure should change in a relatively gradual manner and can be captured by interpolation provided there are enough sampled reference points. Here we interpolated the gradient values collected in section 4.2 from the reference impactor case (I1) for the FMA signal amplitude using cubic interpolation by a factor of 2 between the original impact locations (figure 1).
To test the accuracy of the interpolated gradients, we used them to estimate the maximum impact force of the steel hammer impacts on intermediate locations (case H1a, table 1) between the original reference impact locations (figure 1). Figure 14 shows the results where it can be seen that the interpolated gradient values are able to give an accurate estimation of maximum impact force estimation on locations outside the original reference locations.
To test the robustness of the gradient method towards localisation inaccuracy, we simulated random localisation errors (to maintain focus on impact force estimation methods) with relatively severe magnitude (−20 to 20 mm range, based on previous our previous studies on the same specimen [22,29]) on the x and y coordinates of the reference impactor case (I1a) as shown in figure 15. Using the interpolated gradients from before, we estimated the maximum impact force for the impacts of case I1 but with these new locations. Figure 16 shows the results, where we can see that the gradient method is able to retain a good level of accuracy despite significant localisation inaccuracy. With regards to noise, here we tested the robustness of the gradient method by adding random noise with a magnitude range of −20%-20% of signal amplitude to the extracted FMA amplitude values of the reference impactor case (I1b) to simulate vibration noise [22,29]. We then used the gradient values with the new FMA amplitude to estimate the maximum impact force of case I1. Figure 17 shows the results, where it can be seen that despite a quite severe noise level, the gradient method still retains a good level of accuracy for maximum impact force estimation. Table 4 summarizes the results of testing the robustness of the gradient method towards localisation and noise using gradient values obtained from case I1. We can conclude that the gradient method is able to retain a good level of accuracy despite relatively severe levels of noise and localisation error. In terms of localisation, this suggests that the gradient values   are not so sensitive towards change in location for the specimen configuration used. For the addition of noise, we postulate that the accuracy is also retained at a good level because the maximum impact force estimations were obtained from the average of all the sensors on the plate. If the noise is random (as it is here), then by averaging the estimate of multiple sensors we are effectively averaging the random errors of estimation at each sensor caused by the given noise.

Gradient method for damage assessment
In previous sections, the accuracy of the gradient method under various impact cases and scenarios (sections 4 and 5) was investigated. However, all of these impacts are low energy non-damaging impacts with impact force levels ranging around 50 to 250 N. As we are interested to use this method for damage assessment, it is necessary to test if the linear relationship that forms the base of the gradient method still holds for much larger impact forces that can cause damage and whether effects such as large non-linear deformation have any effect [9]. Thus, here we used the impact tests from the high energy impact setup (1 J, 7 J, 9 J and 11 J) to investigate the relationship between amplitude and impact force at a higher range of forces. Figure 18 shows the recorded impact force profiles, subsequent signals and the ultrasound scan of the impact area indicating whether damage had formed. For all cases except the 1 J impact, impact damage is produced in the test specimens. From previous studies [1,33,34], it is known that the threshold force to produce damage is consistent and is characterised by the first discontinuity in the impact force profile  as stiffness is lost during damage formation. In figure 18 we see for these tests that the damage threshold force (characterised by a 'dip' in the impact force profile) value for the plate specimens is consistently around 2 kN regardless of impact energy.
Plotting the MA signal amplitude vs the maximum impact force ( figure 19), we can see that there is indeed a linear relationship as seen for the impacts with the on-damaging impact setup. Thus, it can be concluded that the gradient method is valid for impact forces at least up to the damage threshold (prior to damage formation) and thus can be used for reliable severity assessment. For forces above the damage threshold, the consistency of the linear relation found here may be affected by the relatively small damage size (<30 mm radius for 200 × 200 mm plate), although larger damage may result in more than BVID and thus can be found without the need of a monitoring system. Additionally, as the boundary conditions differ between the damaging and non-damaging impact setups (2 side vs 4 side clamp), it can be said that the boundary conditions tested do not have an effect on the linear relation used in the gradient method and only changes the value of the gradient.

Impact monitoring of structures using the gradient method and future work
With knowledge from previous sections, we can thus conduct severity assessment in two ways provided the damage threshold force is known (through testing or modelling) or there is an estimate of the level of forces that should be given attention (usually in the kN range). First, the gradient method can be used to estimate the maximum force of an impact, after which the value is compared with the threshold to see if it is a critical impact which warrants structural inspection. Alternatively, the linear relationship of the gradient method can be used in reverse to determine the threshold signal amplitude level (MA or FMA, depending on range of impactor stiffnesses expected) of critical impacts. This further reduces the processing requirements as severity assessment is then done by simply comparing the signal amplitude value to the threshold amplitude.
Similar to any data driven method, the gradient method requires initial data to be collected on the structure that is to be monitored. By having a robust method towards the various possibilities of environmental and operational conditions, we can minimise the required data into that of a single impact case. However, as the characteristics of a structure (transfer function and impact force gradients) may differ when the configuration changes (boundary condition, geometry, material), there is currently still the need to collect a new data set for each structure type monitored. Thus, to further reduce initial data collection needs and increase feasibility for complex structures (with multiple substructures), future development should be aimed at either reducing the amount of data required per structure or understanding the change in characteristics as the structural configuration changes. Additionally, the effect of pre-existing damage on the structural properties (stiffness, boundary condition, etc) may have significant effect on the transfer functions or gradient values (although it may need to be severe enough that it surpasses BVID in the end). Thus, quantification of the effect of damage on the structural characteristics is also a necessary next step.

Conclusion
In this paper we have explored two new methods for maximum impact force estimation, deconvolution and novel gradient method, for the purpose of reliable severity assessment in composite aircraft structures under simulated environmental and operational conditions. It was shown that both methods allow accurate and robust estimation of the maximum impact force from various cases of impacts using minimum initial data (transfer function or gradient) from a single impact case. From further testing we found that the gradient method is robust towards the effects of impact localisation errors and noise. The gradient method also has much less computational and storage requirements and is thus more feasible to integrate with current DAQ systems being developed for structural health monitoring. Thus, we conclude that the proposed gradient method is feasible for impact force monitoring and severity assessment in composite aircraft structures in real life conditions. Future development should be aimed towards reducing initial data collection and quantifying the effects of pre-existing damage needs to increase feasibility in application in real life structures.