Deep Learning-Based Wrapped Phase Denoising Method for Application in Digital Holographic Speckle Pattern Interferometry

: This paper presents a new processing method for denoising interferograms obtained by digital holographic speckle pattern interferometry (DHSPI) to serve in the structural diagnosis of artworks. DHSPI is a non-destructive and non-contact imaging method that has been successfully applied to the structural diagnosis of artworks by detecting hidden subsurface defects and quantifying the deformation directly from the surface illuminated by coherent light. The spatial information of structural defects is mostly delivered as local distortions interrupting the smooth distribution of intensity during the phase-shifted formation of fringe patterns. Distortions in fringe patterns are recorded and observed from the estimated wrapped phase map, but the inevitable electronic speckle noise directly a ﬀ ects the quality of the image and consequently the assessment of defects. An e ﬀ ective method for denoising DHSPI wrapped phase based on deep learning is presented in this paper. Although a related method applied to interferometry for reducing Gaussian noise has been introduced, it is not suitable for application in DHSPI to reduce speckle noise. Thus, the paper proposes a new method to remove speckle noise in the wrapped phase. Simulated data and experimental captured data from samples prove that the proposed method can e ﬀ ectively reduce the speckle noise of the DHSPI wrapped phase to extract the desired information. The proposed method is helpful for accurately detecting defects in complex defect topography maps and may help to accelerate defect detection and characterization procedures. This paper presents a deep learning wrapped phase denoising method for application in digital holographic speckle pattern interferometry for the structural diagnosis of artworks. Based on simulation results, it is further veriﬁed that deep learning can be used to reduce speckle noise in the numerator (sine fringes) and denominator (cosine fringes) of the arctangent function, and an enhanced DHSPI wrapped phase using the arctangent function is provided from the denoising numerator and denominator. The captured data show that this method can e ﬀ ectively process the high-frequency regions of the wrapped phases compared with the DHSPI post-processing, thereby improving the accuracy of defect estimation in artworks. In future work, we will explore the defect prediction problem of artworks using the deep learning method. Meanwhile, we will consider applying deep learning to deformation detection of biological specimens in biomedicine Inhomogeneous surface properties signiﬁcantly hinder reliable phase and may produce a wrong quantitative for phase data. Deep learning may useful these problems.


Introduction
Protecting cultural relics and advocating for the protection of artworks with high artistic value is a mandatory aspect of cultural heritage management to secure the very existence of valuable cultural wealth for future generations. The timely inspection and repair of potential defects in cultural relics are of great significance for long-term preservation. In the last decades, many studies based on holographic interference have been proposed as a means for the implementation of new structural diagnosis tools, skills, and practices in heritage conservation [1][2][3][4][5][6][7][8][9][10][11][12]. Digital holographic speckle pattern interferometry (DHSPI) geometry has been proved to be an important technique applied in contactless art diagnostics non-destructive testing (NDT), since it provides high-resolution imaging of surface optical displacement, which can be used out of the laboratory for direct on-field non-destructive visualization and diagnosis of invisible structural alterations. Correspondingly, a fully portable DHSPI

Computer Simulation for Speckle Fringe and Phase Map
The schematic of the DHSPI system is shown in Figure 1, and the laser is divided into an object beam and a reference beam by a beam splitter. The object beam is driven by mirrors, and eventually illuminates the surface of the rough object. After the reference beam passes through a series of optical devices, the reference beam is combined with the object beam by a collimated beam combiner. The combined beam is recorded by a charge-coupled device (CCD) sensor. In the reference beam path, a mirror on a piezoelectric transducer is used for generating the optic path difference. devices, the reference beam is combined with the object beam by a collimated beam combiner. The combined beam is recorded by a charge-coupled device (CCD) sensor. In the reference beam path, a mirror on a piezoelectric transducer is used for generating the optic path difference.
and we can further derive that  , The phase difference between   v u,  Assume that the object beam U O (u, v) and the reference beam U R (u, v) are described by where (u, v) are spatial coordinates, a O (u, v) and a R (u, v) are the corresponding amplitudes of beam, and ϕ O (u, v) and ϕ R (u, v) are the phase terms. Before the surface of the object is deformed, the resulting intensity I 1 (u, v) of the speckle pattern on the detector can be written as and we can further derive that where . When the surface of the object is deformed, the resulting optical intensity I 2 (u, v) on the detector is defined as The phase difference between I 1 (u, v) and I 2 (u, v) can be expressed as and the subtraction operation between I 1 (u, v) and I 2 (u, v) can be expressed as a speckle fringe pattern, The establishment of the database is necessary to train the model. In this work, the method of manipulating the spectrum of the light field by imitating the 4f system is used to generate a Appl. Sci. 2020, 10, 4044 4 of 13 speckle pattern, and a low-pass filter d(k 1 , k 2 ) is simulated to be applied in object beams to imitate an aperture [32,33]. An object beam can be written as where F represents the Fourier transform and F −1 represents the inverse Fourier transform. For simulation, the phase ϕ O (u, v) as a random variable is uniformly distributed in (−π, π], and the amplitude a O (u, v) is set to a constant of 1. The undeformed intensity I 1 (u, v) can be generated by Equation (3), where U R (u, v) = 1 for simulation. After the object deformation, the phase term ). Similarly, the deformed intensity I 2 (u, v) can also be simulated by Equation (3). According to Equation (7), the speckle fringe pattern can be computed by the subtraction operation between I 1 (u, v) and I 2 (u, v).
The speckle fringe pattern with the different speckle sizes can be generated by adjusting the aperture size d(k 1 , k 2 ) according to Equation (8). An example with simulated data is shown in Figure 2, and the peaks function as the phase difference ∆ϕ(u, v) is simulated. The speckle fringes (500 × 500 pixels) with the different aperture sizes generate different speckle sizes, where the smaller selected aperture produces a larger speckle size.
The establishment of the database is necessary to train the model. In this work, the method of manipulating the spectrum of the light field by imitating the 4f system is used to generate a speckle pattern, and a low-pass filter is simulated to be applied in object beams to imitate an aperture [32,33]. An object beam can be written as represents the Fourier transform and F −1 represents the inverse Fourier transform. For can also be simulated by Equation (3). According to Equation (7), the speckle fringe pattern can be computed by the subtraction operation between   v u, The speckle fringe pattern with the different speckle sizes can be generated by adjusting the aperture size according to Equation (8). An example with simulated data is shown in Figure 2, and the peaks function as the phase difference v) (u, Δ is simulated. The speckle fringes (500 × 500 pixels) with the different aperture sizes generate different speckle sizes, where the smaller selected aperture produces a larger speckle size. The wrapped phase can be estimated from the undeformed and deformed state of the intensity maps by a phase-shifting algorithm [34], and the wrapped phase is obtained by Figure 3 indicates the estimated wrapped phases under different aperture sizes. Here, the wrapped phase maps reveal the evident noise, which may produce errors in assessing defects of artworks. Therefore, effective filtering of the noisy wrapped phase map is important to help the structural diagnosis of artworks. The wrapped phase can be estimated from the undeformed and deformed state of the intensity maps by a phase-shifting algorithm [34], and the wrapped phase is obtained by where S N (u, v) are sine fringe patterns, and C N (u, v) are cosine fringe patterns. ϕ N (u, v) is the wrapped phase containing speckle noise, which is distributed in (−π, π]. Figure 3 indicates the estimated wrapped phases under different aperture sizes. Here, the wrapped phase maps reveal the evident noise, which may produce errors in assessing defects of artworks. Therefore, effective filtering of the noisy wrapped phase map is important to help the structural diagnosis of artworks.

Denoising Method
Convolutional neural networks (CNN) have been successful in processing various visual tasks and have demonstrated effective performance in handling Gaussian denoising [35]. In this work, the

Denoising Method
Convolutional neural networks (CNN) have been successful in processing various visual tasks and have demonstrated effective performance in handling Gaussian denoising [35]. In this work, the proposed network for DHSPI wrapped phase denoising is based on improved denoising convolutional neural networks (DnCNNs) [35]. The deeper neural network architecture for DHSPI wrapped phase denoising is illustrated in Figure 4.

Denoising Method
Convolutional neural networks (CNN) have been successful in processing various visual tasks and have demonstrated effective performance in handling Gaussian denoising [35]. In this work, the proposed network for DHSPI wrapped phase denoising is based on improved denoising convolutional neural networks (DnCNNs) [35]. The deeper neural network architecture for DHSPI wrapped phase denoising is illustrated in Figure 4. In the proposed method, the noisy sine and cosine fringe patterns are denoising by CNN for solving the  2 discontinuities in the noisy wrapped phase. As can be seen from Figure where i W represents the weights, and the * represents convolution operation. In the proposed method, the noisy sine and cosine fringe patterns are denoising by CNN for solving the 2π discontinuities in the noisy wrapped phase. As can be seen from Figure 4, the noisy sine S N (u, v) and cosine C N (u, v) fringe patterns are input data. The input data are first computed by 35 convolutional layers, and the rectified linear unit (ReLU) [36] for non-linear mapping and batch normalization (BN) layers [37] for easier training are used after each convolutional layer. Each convolutional layer (from the 1st to 35th layer) has 32 filters of size 3 × 3 × 32. Supposing that the input data of each layer (from 1st to 35th) is β i , the learned features Z i (β i ) on each layer are expressed mathematically as [38] Zi where W i represents the weights, and the * represents convolution operation. The 36th layer is the convolutional layer, which outputs the denoised sine S(u, v) and cosine C(u, v) fringe patterns. Then, an enhanced wrapped phase map ϕ(u, v) can be provided by arctangent function Network training is performed by finding the minimum loss function between the output data O(k, r) and the ground truth G(k, r). The adaptive moment estimation (Adam) optimizer [39] is chosen as the optimization method. A mini-batch size is expressed in b, and the used loss function L is as follows In this work, the learning rate is set to 1 × 10 −4 . The weights are automatically updated after each iteration, and the training is completed after 100 epochs. Pytorch framework is used to perform network training, and the training process is processed on Ubuntu 14.04 with an Intel ® Core TM i7-7820X CPU@3.66GHz × 13 using the NVIDIA GTX 1080.

Data Set Preparation
To generate a data set for training the model, a simulation method of phase change ∆φ(u, v) is further presented here. The phase change ∆φ(u, v) is related to the out-of-plane displacement H(u, v) as shown in Equation (13) [34], where λ is wavelength (532 nm). First, the displacement H(u, v) with size 80 × 80 pixels is simulated based on Gaussian function as shown in Equation (14). The arithmetic operations of addition and subtraction on 8 sets of Gaussian surfaces with the different mean value η (from 1 to 80) and the different standard deviation ε (9, 13, and 18) are executed to obtain the displacement H(u, v). Subsequently, phase changes with size 80 × 80 pixels are produced according to Equation (13), and 30,000 wrapped phases are generated using the simulation method mentioned in Section 2. Finally, the noisy sine and cosine fringe patterns and the corresponding noise-free data are obtained to train the model.

Evaluation in Simulated Data
The denoising results of simulated data are demonstrated to evaluate the performance. In Figure 5a, the noisy data with the size of 500 × 500 pixels are simulated under the speckle radius of 1 pixel, and the enhanced wrapped phase using the proposed method manifests a significant improvement on the noisy data as shown in Figure 5b. The phase unwrapping of Figure 5b using the least-squares algorithm by the discrete cosine transform (DCT) [40] is implemented in Figure 5c, which has a low error distribution as shown in Figure 5d. The Mean Squared Error (MSE) comparison between the noisy data and denoised data is provided by Table 1, and it can be seen that the noise is greatly reduced from the noisy wrapped phase.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 7 of 15 By adjusting the size of the aperture, the simulated data in the speckle radius of 2 pixels are analyzed as shown in Figure 6a. The performance is demonstrated from the denoising result, and the denoising effect is illustrated by visualization in Figure 6b. In Figure 6c, the phase unwrapping is implemented by DCT for illustrating the denoising effect, and Figure 6d shows the low error distribution for illustrating the feasibility. In Table 1, the MSE comparison between the noisy data and denoising data indicates the good noise suppression ability of the proposed method. Besides, the calculation time for a wrapped phase with the size of 500 × 500 pixels is 3.45 s.

Figure 5a
3.9881 0.5823 Figure 6a 4.3608 0.8654 By adjusting the size of the aperture, the simulated data in the speckle radius of 2 pixels are analyzed as shown in Figure 6a. The performance is demonstrated from the denoising result, and the denoising effect is illustrated by visualization in Figure 6b. In Figure 6c, the phase unwrapping is implemented by DCT for illustrating the denoising effect, and Figure 6d shows the low error distribution for illustrating the feasibility. In Table 1, the MSE comparison between the noisy data and denoising data indicates the good noise suppression ability of the proposed method. Besides, the calculation time for a wrapped phase with the size of 500 × 500 pixels is 3.45 s. denoising result of (a) by the proposed method; (c) unwrapped phase of (b); (d) error distribution of (c).

Data Acquisition
The precise monitoring and visualization of micro-deformation are important for understanding the structure of the defects [41]. To show the denoising effect to help the structural diagnosis, the captured data from the portable DHSPI system [13] are used to verify the performance. Figure 7 shows the schematic diagram of the measuring device [42], where IR lamps are used to provide the thermal excitation for the object. The object response after the thermal excitation is measured and recorded by the DHSPI system, and the personal computer (PC) can visualize results from the recorded data.

Data Acquisition
The precise monitoring and visualization of micro-deformation are important for understanding the structure of the defects [41]. To show the denoising effect to help the structural diagnosis, the captured data from the portable DHSPI system [13] are used to verify the performance. Figure 7 shows the schematic diagram of the measuring device [42], where IR lamps are used to provide the thermal excitation for the object. The object response after the thermal excitation is measured and recorded by the DHSPI system, and the personal computer (PC) can visualize results from the recorded data. The object to be measured is marquetry as shown in Figure 8a. Figure 8b-f shows the estimated wrapped phase from the portable DHSPI system, and every wrapped phase represents a different state of the sample while the sample is cooling down after thermal excitation with IR lamps. The image in Figure 8b is the wrapped phase with a time interval of 12 s after the excitation. Subsequently, there are the wrapped phase maps with defect reaction at 53 s, 103 s, 151 s, and 218 s in a cooling down process after thermal excitation, respectively. However, we can observe obvious noise from these wrapped phase maps. In the high-frequency regions of the phase denoted by the red arrow, noise causes phase distortion. In such cases, it is difficult to directly estimate the depth information of the defect.

Denoising for Captured Data
Post-processing software of the portable DHSPI system can process the recorded results to improve the visualization for extracting useful information. The basic post-processing tools include denoising interferogram, phase unwrapping, displaying phase unwrapped image, and 3D surface The object to be measured is marquetry as shown in Figure 8a. Figure 8b-f shows the estimated wrapped phase from the portable DHSPI system, and every wrapped phase represents a different state of the sample while the sample is cooling down after thermal excitation with IR lamps. The image in Figure 8b is the wrapped phase with a time interval of 12 s after the excitation. Subsequently, there are the wrapped phase maps with defect reaction at 53 s, 103 s, 151 s, and 218 s in a cooling down process after thermal excitation, respectively. However, we can observe obvious noise from these wrapped phase maps. In the high-frequency regions of the phase denoted by the red arrow, noise causes phase distortion. In such cases, it is difficult to directly estimate the depth information of the defect. The object to be measured is marquetry as shown in Figure 8a. Figure 8b-f shows the estimated wrapped phase from the portable DHSPI system, and every wrapped phase represents a different state of the sample while the sample is cooling down after thermal excitation with IR lamps. The image in Figure 8b is the wrapped phase with a time interval of 12 s after the excitation. Subsequently, there are the wrapped phase maps with defect reaction at 53 s, 103 s, 151 s, and 218 s in a cooling down process after thermal excitation, respectively. However, we can observe obvious noise from these wrapped phase maps. In the high-frequency regions of the phase denoted by the red arrow, noise causes phase distortion. In such cases, it is difficult to directly estimate the depth information of the defect.

Denoising for Captured Data
Post-processing software of the portable DHSPI system can process the recorded results to improve the visualization for extracting useful information. The basic post-processing tools include denoising interferogram, phase unwrapping, displaying phase unwrapped image, and 3D surface

Denoising for Captured Data
Post-processing software of the portable DHSPI system can process the recorded results to improve the visualization for extracting useful information. The basic post-processing tools include denoising interferogram, phase unwrapping, displaying phase unwrapped image, and 3D surface deformation map, etc. [13]. The noisy data of Figure 8b-f are processed by the post-processing of the portable DHSPI system as shown in column (a) of Figure 9. Then, the proposed method is further used for denoising as shown in column (b) of Figure 9. From the visual inspection, the clean phase maps are obtained by the proposed method.  To clearly compare performance, in Figure 9 the enlarged phases in the red regions of column (a) and blue regions of column (b) are shown in columns (c) and (d), and the defect state can be observed more intuitively from the enhanced wrapped phase map using the proposed method. It can be seen that the proposed method can handle the high-frequency regions of the phase well, which shows it to be a useful post-processing method to estimate the depth of the defects in artworks. In the proposed method, the calculation time is 3.60 s for a captured wrapped phase map with size 1017 × 1030 pixels. Figure 10a,c shows the phase unwrapping results of the noisy data ( Figure 8) and filtered data (Figure 9b) using DCT. The phase unwrapping effect is further illustrated by visualization of the contours in Figure 10b,d, and it further proves that the smooth phase can be obtained after denoising of the proposed method.

Conclusions
This paper presents a deep learning wrapped phase denoising method for application in digital holographic speckle pattern interferometry for the structural diagnosis of artworks. Based on simulation results, it is further verified that deep learning can be used to reduce speckle noise in the numerator (sine fringes) and denominator (cosine fringes) of the arctangent function, and an

Conclusions
This paper presents a deep learning wrapped phase denoising method for application in digital holographic speckle pattern interferometry for the structural diagnosis of artworks. Based on simulation results, it is further verified that deep learning can be used to reduce speckle noise in the numerator (sine fringes) and denominator (cosine fringes) of the arctangent function, and an enhanced DHSPI wrapped phase using the arctangent function is provided from the denoising numerator and denominator. The captured data show that this method can effectively process the high-frequency regions of the wrapped phases compared with the DHSPI post-processing, thereby improving the accuracy of defect estimation in artworks. In future work, we will explore the defect prediction problem of artworks using the deep learning method. Meanwhile, we will consider applying deep learning to deformation detection of biological specimens in biomedicine [43][44][45][46]. Inhomogeneous surface properties significantly hinder reliable phase unwrapping and may produce a wrong quantitative evaluation for phase data. Deep learning may be a useful technique to solve these problems.