Three-dimensional subsurface defect shape reconstruction and visualisation by pulsed thermography

Defects detected by most thermographic inspection are represented in the form of 2D image, which might limit the understanding of where the defects initiate and how they grow over time. This paper introduces a novel technique to rapidly estimate the defect depth and thickness simultaneously based on one single-side inspection. For the first time, defects are reconstructed and visualised in the form of a 3D image using cost-effective and rapid pulsed thermography technology. The feasibility and effectiveness of the proposed solution is demonstrated through inspecting a composite specimen and a steel specimen with semi-closed airgaps. For the composite specimen, this technique can deliver comparatively low averaged percentage error of the estimated total 3D defect volume of less than 10%.


Introduction
As a highly efficient and powerful non-destructive testing (NDT) technique, pulsed thermography (PT) offers a rapid, contact-free, highefficient inspection and shows promising applicability to in-situ monitoring applications [1][2][3][4]. Beyond defect detection, quantitative characterisation of defect information including shape, size and depth, and estimation of thermal properties have been proven to be essential and decisive for through-life engineering [5][6][7][8][9][10][11]. The features within the transient temperature of inspected object have been widely used to evaluate defect depth (the distance from the inspected surface to the top surface of the defect). For example, the Peak Slope Time method (PST) [4,12] employs the peak time of the first derivative curve of temperature contrast, which is approximately proportional to the square of the defect depth. Furthermore, Absolute Peak Slope Time method (APST) [13] multiplies the square root of its time to temperature decay curve, and then estimates the defect depth by the peak time of the first derivative curve. Differently, the Logarithm Second Derivative method (LSD) [14] estimates the defect depth by fitting the logarithmic second derivative of temperature decay curve using a polynomial model. The Nonlinear System Identification method (NSI) [7], although similar to the LSD method on principles, improves depth estimation and realizes an automatic model order selection for polynomial fitting. The Least-Squares Fitting method (LSF) [15,16] fits the temperature decay curve with a theoretical heat transfer model to determine the defect depth directly.
It has been observed that the measured defects of aforementioned methods are usually presented as 2D images. However, 2D visualisation is rather limited in representing and explaining the evolution and progression of a defect. Defects are typically not two-dimensional objects but evolve in 3D space over time. It is of great significance to predict the life of in-service components and to give feedback to design and maintenance [17]. 3D visualisation of inspected parts can provide a better understanding and more details of the defects and reduce operational time and improve quality control of production in industry. For example, in the area of manufacturing, 3D images from computed tomography (CT) has been widely used to analyse the types of porosity defects occurred in the materials from castings [18][19][20][21]. In the area of nuclear and aerospace industries, 3D visualisation is reconstructed from the digital X-ray images and has been used to view the location, shape and size of the defects like corrosion, delamination and crack, and evaluate the thickness of walls in the object of titanium aerospace investment casting [22,23]. However, the X-ray technology is comparatively time consuming [24], comes with potential health risks [25] to the user and is typically limited with respect to the maximum size of the parts to be analysed [22,[26][27][28][29].
Only few publications have studied the 3D reconstruction and visualisation of subsurface defects by pulsed thermographic inspection. Plotnikov and Winfree [30] visualised a 3D tomogram of the defects by using reversed time of the peak slope masked by the binary contrast image which is constructed by the thermal contrast method. Although this method can visualise the defects in the form of a 3D image, it can visualise only the depth of the defect and cannot visualise the thickness of the defect. Ramirez-Granados et al. [31] proposed an approach for 3D reconstruction of subsurface defects by using a finite-difference model. Firstly, a non-defect nodal network is established with the properties and characteristics of the inspected object for detecting the inner defects. Next, the shape, size, depth, thickness, and location of subsurface defects are computed by means of the minimisation of a cost function and the determination of a depth and a thickness function. After that, each non-defect node of the established nodal network is replaced with the computed defect node. This method requires the knowledge of a-priori properties of the defects such as thermal conductivity, density, and specific heat capacity at constant pressure, which are typically unknown in on-site inspections. Elhassnaoui and Sahnoun [32] proposed a method for 3D reconstruction of defect shape located on the inaccessible back of a homogeneous material without the need of thermal properties such as thermal diffusivity. Based on thermal distribution, this method analyses thermal response on the object's surface and computes defect distance (defect depth and sample thickness). It modifies the APST method [13] to evaluate defect distance of the object by dividing arbitrary two points from the sample surface to eliminate the thermal diffusivity. However, this method can only be applied to the characterisation of surface defects.
3D visualisation of hidden defects based on PT is promising, attactive but challenging. A direct approach for pulsed thermography could be to conduct two inspections from both the front and back of a part. The defect depth for each inspection is calculated based on aforementioned methods, by which means the defect thickness can be computed by considering the sample thickness. However, the deployment of this approach can be limited because 1) one side of the inspected component could be inaccessible e.g. an aircraft wing or fuselage; 2) the accuracy of measurement could be compromised if the defect thickness is very thin due to extreme closed values of defect depths from two inspections; 3) if the defect is too deep one side, the defect could be missed, and 4) it introduces extra cost of inspection time.
This paper reports a new method of 3D reconstruction and visualisation of subsurface defect shape based on a single-side thermographic inspection. The proposed method with related theory background is presented in Section 2. The experimental results are presented in Section 3, while the conclusions are given in Section 4.

Principle of pulsed thermography
In a pulsed thermographic inspection in reflection mode, a short and high energy light pulse is projected onto the sample surface through one or two flash lamps (see Fig. 1(a) for typical experimental setup) [33].
Heat conduction then takes place from the heated surface to the interior of the sample, leading to a continuous decrease of the surface temperature [16]. An infrared camera controlled by a PC captures the timedependent response of the sample surface temperature. If the sample is defect-free, the time when the temperature deviation occurs can be used to estimate the thickness (if the thermal diffusivity is known) or thermal diffusivity of local materials (if the thickness is known [5]). For example, if the thermal diffusivity is known, the thickness of the sample can be estimated based on the time of temperature deviation of t 2 , collected from the point 2 on the inspected surface (see Fig. 1(b)). This approach can be extended to measure defect depth, d1 (see Fig. 2), by considering the first time of temperature deviation of t 1 , collected from the point 1. The surface temperature due to the back-wall at depth L for a homogeneous plate is given by [34] where T t ( ) is the temperature variation of the surface at time t, Q is the pulse energy, ρ is the material density, c is the heat capacity, k is the thermal conductivity of the material, R is the thermal reflection coefficient of the interface with air, and α is the thermal diffusivity. Eq. (1) can also be rewritten as where e denotes the thermal effusivity, which equals to ρck . Although Eq.
(1) only applies to 1D heat transfer, it is proposed to approximate the 3D heat transfer in the real world to simplify the problem.

Single-side inspection
As stated in Section 1, although the double-side inspection is straightforward, there are a few limitations. This paper introduces a novel single-side inspection approach to overcome these limitations. The proposed method estimates the defect thickness h x y ( , ) by establishing a predictive model between h, the defect depth d (e.g. either d1 or d2) estimated from a single-side inspection, thermal wave reflection coefficient (R), and shortest distance from the boundary of the defect to the inspected point location p x y ( , ) (D dis ), written as The defect thickness can then be inferred based on this model. The values of d x y ( , ) and D x y ( , ) dis are achievable using the existing methods, but the challenge is to measure R x y ( , ). In this paper the New Least-Squares Fitting method (NLSF) [8] is applied to estimate R and d simultaneously.
The thermal wave reflection coefficient (R) can be directly computed from the Eq. (4). The analytical model is written as fusivity, t is sampling duration, t s is the starting time of sampling, M is a large iteration number, and T NLSF is computed temperature from optimisation technique. There are five parameters to be estimated including A, W , R, t s , and s. This method employs a nonlinear least-squares solver (in Matlab called lsqnonlin) to solve this five-parameters optimisation problem. Through initially setting the lower and upper bounds for each parameter, the NLSF estimates the optimal parameters that has where T is the raw temperature on the surface. Once the optimal parameters are estimated, if α is known, the defect depth (d) can be estimated by Alternatively, if d is known, the thermal diffusivity can be estimated by

Double-side inspection
The method of double-side inspection reconstructs the 3D structure of defects by evaluating the defect depth from both sides of the sample. As seen in Fig. 2, the defect thickness (h), at the position of x and y, can be calculated by where h x y ( , ) is the defect thickness, L x y ( , ) is the thickness of the sample, d x y 1( , ) is the defect depth from the top surface, and d x y 2( , ) is defect depth from the bottom surface. The thickness of the sample can be measured by general measurement tools such as ruler or Vernier calipers, or by the thermographic inspection directly through the method introduced in Section 2.1 if an automatic approach is required. Defect depth can be estimated by methods such as PST [4,12], LSD [14], APST [13], NSI [7], LSF [15,16], or NLSF [8] method. In this paper, the NLSF method is employed to estimate the defect depth due to its fine performance.

Inspection process
The single-side inspection and the double-side inspection can be described in the flowcharts as shown in Fig. 3. For the single-side inspection, the first step is to select the region of interest (ROI), which aims to select a region for defect estimation and visualisation to reduce overall processing time. The ROI was selected manually in this paper. The second step is defect measurement to obtain subsurface information of ROI including defect depth (d), the shortest distance from the boundary of the defect to the inspecting point (D dis ), and thermal wave reflection coefficient of the defect (R), which are achieved by the NLSF method. The third step is noise reduction to enhance image quality (R image and d image), which is achieved by the application of a Median Filter [35] in this paper. The fourth step is to establish the relationship model between h, R, and D dis , shown in Eq.
(3), based on the estimated parameters in the second step. The fifth step is to estimate the defect thickness based on the established model. The last step is to visualise the defect in the form of a 3D image, where the reconstructed defect (size, depth, thickness, and location of the defect) is fused with the dimension (width, length and thickness) of the sample to produce a volume image. This volume image is then rendered and displayed in the form of a 3D image. For the double-side inspection, most steps are similar to the single-side inspection except the second, the fourth, and the fifth step. In the second step, the sample's thickness and defect depth of the both side of the sample are required, whilst the defect depth from the single-side inspection requires one side only. The fourth step is to align the position of the defect region between the front side and the back side and the fifth step is to estimate the defect thickness from Eq. (8).

Sample 1
A flat plate of carbon fibre reinforced polymer (CFRP) material was used in this study. The plate was made of unidirectional Toray 800 carbon fibres pre-impregnated with Hexcel M21 epoxy resin and manufactured in a traditional autoclaving process. The dimension of Sample 1 is 75 mm × 230 mm × 8 mm. As illustrated by Fig. 4(a), it includes five block defects with a thickness (h) of 0.5 mm, 1.0 mm, 2.0 mm, 3.0 mm, and 4.0 mm, respectively, named Defect1 to Defect5. The width (r), length, and depth (d) for all defects are 10 mm, 75 mm, and 2 mm respectively. The distance between two adjacent defects is 30 mm, which aims to reduce the influence from the adjacent defects on the thermal behaviour. A side view of the produced sample is shown in Fig. 4(b). It should be noted that these defects are not fully closed because two sides of the defects are open. This sample is aimed at studying the relationship between R, h, and D dis when d is fixed.

Sample 2
As a use case, Sample 2 was made from steel with a dimension of 210 mm × 240 mm × 35 mm, as shown in Fig. 5. There is an 'S' shape triangular air-gap through the sample. By measurement, at the top side of the sample, the defect is estimated to be triangular in shape, which has the base about 25 mm and the maximal thickness about 20 mm, illustrated in Fig. 5(c). At the bottom side of the sample, the defect is also estimated to be triangular in shape, which has the base about 16 mm and the maximal thickness about 16 mm, illustrated in Fig. 5(d). The sample was inspected at the flat side using the proposed single-side thermographic inspection, illustrated as Fig. 5(b).

Experimental setup
The experiments were conducted using the Thermoscope® II pulsedactive thermography system that comprises of two capacitor banks powered Xenon flash lamps mounted in an internally reflective hood and a desktop PC to capture and store data. The scheme of the experimental set-up is illustrated in Fig. 1(a). A FLIR SC7000 series infrared (IR) radiometer operating between 3.0 and 5.1 µm and a spatial resolution of 640 × 512 pixels was used to perform the inspection. The samples were placed with their surface perpendicular to the IR camera's line of sight at 250 mm from the lens. The apply energy was approximately 2 kJ over an inspection area of 250 mm × 200 mm. The pixel pitch is 0.32 mm. Considering the thickness of Sample 1 and the low thermal diffusivity of CFRP, a sampling rate of 10 Hz was used and totally 900 frames (equivalent to 90 s) were captured after the flash. Considering the thickness of Sample 2 and the high thermal diffusivity of steel, a sampling rate of 50 Hz was used and totally 500 frames (equivalent to 10 s) were captured.

Single-side inspection
Due to the large width of Sample 1 (230 mm), the defects of this sample suffer non-uniform heating for a single inspection, which could lead to unreliable results [36]. To reduce this effect, in this study the sample has been inspected five times, where each defect was placed on the centre of the camera's view once. An area of the centralised defect with the size of 160 × 120 pixels for five data files was then cropped and merged into one file with a size of 160 × 600 pixels for easier analysis. Fig. 6 shows a snapshot of the raw thermal image at time 23 s. The defects are labelled as "(1)", "(2)", "(3)", "(4)", and "(5)" to represent the defect thickness of 0.5 mm, 1.0 mm, 2.0 mm, 3.0 mm, and 4.0 mm, respectively. The ROI of each defect is about 10 mm × 10 mm (around 30 × 30 pixels), as highlighted in Fig. 6, which are used for evaluating the performance of the proposed technique.
The relationship between defect thickness (h) and thermal wave reflection coefficient (R) requires to be established before prediction of h. To reduce the influence of the heat leaked to each opened side on the results, only a region of 6 × 30 pixels on the middle area of each defect were sampled to establish the model. Estimated R values were averaged and are shown in Table 1. The averaged R values, denoted by R ave , in the middle area were used as the representative of R to establish the R ave vs h relationship. It is observed that the averaged value of R increases following the increase of defect thickness. This observation is expected because the heat is more difficult to be leaked through the 3D conduction around air-gap if the defect is thicker. As can be seen from the plot between R ave and h in Fig. 7, the associate is not linear but approximately exponential. An exponential fitting was applied on them, and the relationship between R ave and h can be written as To estimate the value of h using the measured R ave , the model (9) can be reversed as   Fig. 8 that R value at both sides of defect's boundary is less than the middle region. This observation could be caused by the fact that the heat around the middle of the defect is more difficult to leak than that around the boundary of the defect. This issue must be addressed if the identified model from the regions of 6 × 30 pixels is applied to the ROI of 30 × 30 pixels. To improve the accuracy of R estimation, the shortest distance from the boundary of the defect to the inspected point (D dis ), shown in Fig. 2, is taken into account to correct the R value. The re-calculated R (R new ) is a function of R and D dis new dis (11) This function can be approximated by a polynomial model, parameters of which can be estimated based on the known D dis . The estimated model is written as  (13) The average of the estimated defect thickness of each defect for the regions of 30 × 30 pixels is shown in Table 2. It is suggested that the maximum error of the estimated depth for all defect is less than 13%. For Defect 2-4, the error is less than 5%. It is interesting to observe that Defect 1 and Defect 5, which are close to the boundary of the sample, have the larger error than others, which could be caused by the manufacturing error. The average R value increases following the increase of the measured defect thickness (h). The error of the estimated defect thickness is less than 19%, particularly for Defect 1-4, the error is less than 10%. The large error of h for Defect 5 is caused by the saturated R value, which is 0.53, the same as the R value for Defect 4.
The reconstructed defects for the regions of 30 × 30 pixels are visualised in Fig. 9, which suggests that the proposed technique can effectively reconstruct the 3D structure of simple defects. The estimation of the bottom surface has an increasing error following the increment of defect thickness.

Double-side inspection
In this experiment Defect 4 and Defect 5, with the thickness of 3 mm and 4 mm, respectively, were in the focus because they can be detected from both sides of the sample. Defect 1, Defect 2 and Defect 3 were not selected because the defect depth of the back of these defects are too deep (5.5 mm, 5.0 mm, and 4.0 mm, respectively) to measure considering the depth to width ratio of a defect. The NLSF method was applied on the selected regions to estimate the defect depth where the thermal diffusivity (α) was chosen as 0.5 × 10 −6 m 2 /s. The estimated defect depth after reducing noise from inspection of two sides are shown in Fig. 10. Observation from the front of the part shows that the measured depths for Defect 4 and Defect 5 are similar. This is expected as both defects have a 2 mm depth, as illustrated in Fig. 4. The depth map from the back shows that Defect 4 (light blue) is much deeper than Defect 5 (deep blue). In the data analysis process, the defect position at the front and the back of the part was aligned and registered. After that, the defect thickness was calculated by Eq. (8). The estimated parameters including defect depths and thickness are shown in Table 3. It is observed that, on the front, the average of the estimated defect depth of Defect 4 is approximately identical to the design (2 mm) with an error of 1%, while Defect 5 has about 10% error. For the back of the part, the average error of the estimated defect depth for Defect 4 and Defect 5 are around 18% and 16%, respectively. The computed defect thickness of Defect 4 and Defect 5 are 3.55 mm and 4.50 mm, which have the percentage error of 18.33% and 12.50%, respectively. Finally, the surface information (dimension of the sample) and subsurface information (defect depth and defect thickness) were used to build the 3D volume image and visualised in the form of 3D image, as shown in Fig. 11. It is clearly visible that Defect 5 is thicker than Defect 4, and the estimated 3D volume for both defects is close to the design shown in Fig. 4. It should be noted that the z-axis is scaled up to better visualise the detail of the defect surface.

Results of Sample 2
Sample 2 was inspected by the single-side inspection because the defect was too deep to measure from the other side. To find the relationship between R and h, R values in three areas with a size of 100 × 60 pixels were sampled as highlighted in Fig. 12. It should be pointed out that the R values of 100 vertical pixels were averaged. Fig. 13 plots the cross-section at line 200 of the estimated R and d, where the three selected areas are highlighted. It can be observed that the middle area has the highest value of R, then the left one and the right now. Proven in the experiments of Sample 1, a higher value of R suggests a thicker defect. As illustrated in Fig. 5, the maximal thickness of the left defect is about 20 mm, and the right one is about 16 mm. Based on the observation of the R values, the maximum thickness of the middle area is assumed to be 24 mm. A numerical model was then established to represent the relationship between R and h, written as As the last step, the 3D volume image was reconstructed by combing the surface dimension (size and thickness of the sample) and subsurface information (defect depth and defect thickness) and result is shown in Fig. 14. The triangular 'S' shape defect can be clearly observed which matches the appearance of the defect shown in in Fig. 5(c) and Fig. 5(d). To validate the result using other NDT methods, Fig. 15 shows the comparison of the 2D image of the defect from pulsed thermography and X-ray. It can be seen that the defect shape of both images is similar. A quantitative comparison between these two imaging modalities in terms of 3D defect reconstruction was not conducted in this paper due to the distortion of thermal images caused by the lens. A further study is required to correct this distortion before conducting a comprehensive comparison.

Conclusions
This paper has developed a novel 3D reconstruction and visualisation approach for subsurface defect based on one single-side PT inspection under the reflection mode. The defect thickness, a key parameter to reconstruct defect in a 3D form, was measured through estimating thermal wave reflection coefficient value and establishing its relationship with defect thickness using an empirical model. For the comparison purpose, this paper also introduces a double-side inspection approach that estimates the defect thickness by measuring the defect depths from both sides of a part. And then the estimated thickness, depth and size of the defect, including dimension of the part were combined and visualised in the form of a 3D image.
The proposed technique has been tested on two samples with artificial defects on carbon fibre-reinforced polymer (CFRP) and steel materials, and the results show that:  a. In comparison with the double-side inspection, the proposed singleside inspection is not only faster, but also increases the applicability of 3D defect reconstruction. There are many cases where industrial components can be accessed only on one side, such as airplane wing. Furthermore, some defects are too deep to detect or measure the depth, such as Defect 1-3 in Sample 1 and Sample 2, which limits the applications of the double-side inspection. b. For the Sample 1, the single-side inspection can measure the defect thickness with an error less than 10% if the thickness is less than 3 mm. The error increases if the thickness increases. For example, the error for the defect with a thickness of 4 mm is 19%. c. The double-side inspection can produce more accurate thickness measurement if the defect is thick. For defects with thin thickness, the single-side inspection is more appropriate. d. The thermal wave reflection coefficient is not only related to the thickness of defect but also related to the shortest distance to the defect boundary. e. The proposed 3D defect reconstruction solution can effectively assess the defect or damage in CFRP and steel by offering more details of structure and friendly visualisation.
A drawback of this method is that the defect thickness is a function of multiple correlated parameters, which may constrain its application on irregular shape defects. Based on experiments in this research, it is also discovered that the model to estimate defect thickness is subject to materials. Besides, changing of material properties of tested specimen and defect shape may also have impact on the precision of method model and then the estimation of defect thickness. To fully explore its potential and improve the applicability of the proposed method, a further study is required by considering different materials with a variety of defect shape, size, depth and thickness. We will also focus on a more comprehensive comparison to other NDTs for 3D defect reconstruction and visualisation.

Declaration of Competing Interest
The authors declared that there is no conflict of interest.

Acknowledgement
This work was partly supported by the Lloyd's Register Foundation under Grant number GA\100113, and partly supported by the UK Engineering and Physical Sciences Research Council (EPSRC) Platform