The lossless estimation of a target’s light-intensity information in an underwater environment requires the removal of the influence of backscattered light, which is challenging in the fields of archaeology, noninvasive pathological detection, and, particularly, three-dimensional surface contouring^[1–3]. It is well known that suspended particles existing in underwater environments produce scattered light, which results in the inaccurate acquisition of target information^[4,5].

Some approaches have been suggested to solve this fundamental and practical problem over the past few years^[6–10]. The polarization imaging technique based on a physical model is one of the most promising methods to estimate a target’s intensity information, which benefits from its ability to efficiently remove backscattered light^[11–13]. Unfortunately, traditional polarization imaging methods abide by two unrealistic assumptions. Firstly, the degree of polarization (DOP) of the light reflected by the target equals zero or a constant. Secondly, the angle of polarization (AOP) of the light reflected by the target is the same as that of the backscattered light^[14,15]. However, according to Fresnel’s law, the DOP and AOP of the light reflected from a target are nonuniform because the refractive index and reflection angle are different^[11,16]. For example, when a target comprises both specular-reflecting and diffuse-reflecting materials, its DOP and AOP will be nonuniform due to its nonuniform refractive index, and the reconstructed target information obtained via traditional methods will be distorted^[17]. As shown in Fig. 1, traditional polarization imaging methods are suitable for the reconstruction of a target’s rough surfaces because the DOP of the target is approximately zero in these regions, which is consistent with the assumption stated above^[18,19]. However, it is not appropriate for the reconstruction of metal targets, because they effectively maintain polarization.

For removing backscattered light completely without losing the target information captured by a detector, this study proposes a polarization imaging method considering the nonuniform polarization characteristics of the target. Specifically, this method allows the light reflected from the target to have its own DOP and AOP for each pixel. By analyzing the variable nature of the backscattered light intensity in each polarization subimage, independent component analysis (ICA) is adopted to calculate the intensity information of the fully polarized portion of the backscattered light. Based on the polarization imaging model, which considers both DOP and AOP, the target information can be estimated accurately. The experimental results show that the error rate of the target information recovered using our method is considerably lower than that by the existing methods.

We first analyze the influence of the polarization characteristics of the target-reflected light on the reconstruction results via simulation. In the simulation, the DOP and AOP of the backscattered light are assumed to be zero for convenience. Simultaneously, different values are assigned to the DOP and AOP of the target-reflected light. Then, the maximum and minimum light intensities of the polarization images are calculated to recover the target information using the traditional polarization imaging methods^[1]. Finally, the degree of distortion is used to measure the influence of the traditional methods on information reconstruction, as shown in Eq. (1): $$R_e=\frac{T_{\text{real}}-T_{\text{reconstructed}}}{T_{\text{real}}},$$(1)where $T_{\text{real}}$ represents the real target information, and $T_{\text{reconstructed}}$ represents the reconstructed target information obtained under the unrealistic assumptions about DOP and AOP mentioned above.

The relationship between the degree of distortion and the target polarization information (DOP and AOP) is shown in Fig. 2. The degree of distortion is zero for unpolarized light reflected from the target. However, using traditional methods causes severe distortions when the polarization characteristics of the target-reflected light do not match the assumptions of traditional methods.

To address the inability of calculating the target information with nonuniform polarization characteristics, this study proposes a novel polarization descattering imaging method that takes the polarimetric characteristic of spatial variations into account, as shown in Fig. 3(a). For an image acquired by the detector, the total light intensity $I_{\text{total}}$ (dark blue arrow) is obtained by superimposing two light beams, the target-reflected light $I_{\text{target}}$ (green arrow) and the backscattered light $I_{\text{scat}}$ (red arrow)^[20], as shown in Fig. 3(b). $I_{\text{total}}$, $I_{\text{target}}$, and $I_{\text{scat}}$ comprise partially polarized light with different DOP and AOP. The horizontal direction ($x$-axis direction) is taken as the reference axis. Then, $\mathit{\Phi }$ represents the AOP of the target-reflected light, whose value is equal to the angle between the polarization vibration direction and the reference axis. Similarly, $\beta$ and $\delta$ are the AOP of the irradiance $I_{\text{total}}$ and backscattered light $I_{\text{scat}}$, respectively. The relationship among the AOP of $I_{\text{total}}$, $I_{\text{target}}$, and $I_{\text{scat}}$ is graphically shown in Fig. 3(c).

According to the Stokes vector, the linear polarization information of light waves can be calculated by 0°, 45°, 90°, and 135° polarization subimages^[21–24] using Eq. (2): $$\begin{gathered} I_i(x,y)=T_n(x,y)/2+S_n(x,y)/2 \\ +{\cos }^2(\delta -i)S_p(x,y)+{\cos }^2[\varphi (x,y)-i]T_p(x,y), \\ i=0°,45°,90°,135°, \end{gathered}$$(2)where $T_n$, $T_p$ and $S_n$, $S_p$ are the unpolarized and polarized parts of the target-reflected light and backscattered light, respectively. $(x,y)$ represents the pixel position. For the backscattered light, because the inherent parameters of its polarization characteristics are spatial constants^[25,26], its DOP and AOP can be estimated by selecting a region without the target in the scene. However, for the target-reflected light, its intensity, DOP, and AOP are different for each pixel. Therefore, $(x,y)$ is used to indicate whether the parameter is a space constant. For example, $\delta$ is a spatial constant, whereas $I_i(x,y)$, $T_n(x,y)$, $S_n(x,y)$, $S_p(x,y)$, $\mathrm{\Phi }(x,y)$, and $T_p(x,y)$ are globally variable. Equation (2) shows that an image obtained through a polarizer is a linear sum of the target light and backscattered light. The coefficients ${\cos }^2(\delta -i)$ and ${\cos }^2[\varphi (x,y)-i]$ depend on the orientation of the polarizer and the target’s AOP. If one of these parameters is unmeasurable, the value of the coefficient remains unknown. This situation corresponds to the suppositions of ICA, which is a new and powerful nonlinear data analysis method that is gaining popularity. This method can recover the original source signals from the observed mixtures of probabilistically independent source signals without knowing the manner in which the sources are mixed^[27].

According to ICA, we assume that $T_n$, $S_n$, $T_p$, and $S_p$ are the original signals and ${\cos }^2(\delta -i)$ and ${\cos }^2[\varphi (x,y)-i]$ are the mixing coefficients. After calculating the AOP of the backscattered light from the region without the target, we can calculate the mixing coefficient ${\cos }^2(\delta -i)$. In this case, estimating the target information using traditional matrix transformation methods seems impossible because $T_n$, $S_n$, $T_p$, $S_p$, and ${\cos }^2[\varphi (x,y)-i]$ are unknown. However, ICA accomplishes this using the probabilistic independence property of the original signals. Therefore, to accurately estimate the target-reflected light, this study assumes the probabilistic independence of the target-reflected light and backscattered light. Equation (3) can be obtained from Eq. (2): $$e_1=S_p(x,y)+f({\alpha }_{0\_45})T_p(x,y),$$(3)i.e., $$e_1=\frac{I_0(x,y)-I_{45}(x,y)}{{\cos }^2\delta -{\cos }^2(\delta -45)},$$(4)where $$f({\alpha }_{0\_45})=\frac{{\cos }^2\varphi (x,y)-{\cos }^2[\varphi (x,y)-45]}{{\cos }^2\delta -{\cos }^2(\delta -45)}.$$(5)

Similarly, Eq. (6) can be obtained from Eq. (2) as well: $$\begin{gathered} e_2=S_p(x,y)+f({\alpha }_{0\_90})T_p(x,y), \\ {e}_3=S_p(x,y)+f({\alpha }_{0\_135})T_p(x,y). \end{gathered}$$(6)

Equations (3) and (6) can be expressed in the form of a matrix as follows: $$E=A\left[\begin{array}{c}{S}_p(x,y)\\ T_p(x,y)\end{array}\right],$$(7)i.e., $$E=\left[\begin{array}{c}{e}_1\\ e_2\\ e_3\end{array}\right],$$(8)where $$A=\left[\begin{array}{c}\begin{array}{cc}1& f({\alpha }_{0\_45})\end{array}\\ \begin{array}{cc}1& f({\alpha }_{0\_90})\end{array}\\ \begin{array}{cc}1& f({\alpha }_{0\_135})\end{array}\end{array}\right].$$(9)

According to Eq. (7), $S_p$ and $T_p$ can be separated based on the probabilistic independence of the target-reflected light and backscattered light. Then, singular value decomposition is used to approximate $E$ via a matrix with a smaller rank, and $E$ is decomposed into the product of two matrices^[27]: $$E=UD{V}^T.$$(10)

The desired decomposition is represented as follows: $$E_{3\times n}={\hat{U}}_{3\times 2}{W}_{2\times 2}{W}_{2\times 2}^{-1}{\hat{D}}_{2\times 2}{({\hat{V}}^T)}_{2\times n},$$(11)where $\hat{D}=\mathrm{diag}({\sigma }_1,{\sigma }_2)$ when the singular values of $E$ are ${\sigma }_1>{\sigma }_2>{\sigma }_3$. $W$ is an arbitrary $2\times 2$ nonsingular matrix. $\hat{U}$ and $\hat{V}$ are submatrices of $U$ and $V$. From this expression, estimations of $A$ and $\left[\begin{array}{c}{S}_p(x,y)\\ T_p(x,y)\end{array}\right]$ are given as follows: $$A=\hat{U}W,$$(12)$$\left[\begin{array}{c}{S}_p(x,y)\\ T_p(x,y)\end{array}\right]=W^{-1}\hat{D}({\hat{V}}^T).$$(13)

Therefore, the value that remains to be determined is $W$, which can be represented as follows: $$W=\left[\begin{array}{cc}{r}_1\cos \alpha & r_2\cos \theta \\ r_1\sin \alpha & r_2\sin \theta \end{array}\right],$$(14)where $r_1$, $r_2$, $\alpha$, and $\theta$ represent the possible values for the elements of the matrix $W$, $r_1$ and $r_2$ are positive, and $r_1\cos \alpha$ and $r_1\sin \alpha$ can be calculated from the first column of matrix $A$, as shown in Eq. (15): $$\left[\begin{array}{c}{r}_1\cos \alpha \\ r_1\sin \alpha \end{array}\right]={\hat{U}}^T\left[\begin{array}{c}1\\ 1\\ 1\end{array}\right].$$(15)

The determinant of $W$ is $\mathrm{\Delta }=r_1{r}_2\sin (\theta -\alpha )$. We set $|\mathrm{\Delta }|=1$ without losing generality^[27]. According to Eq. (14), once $\theta$ is known, the target information can be calculated by solving the matrix $W$. Based on this relationship, we gradually change $\theta$ and compute the polarization component of the target and backscattered light using Eq. (11). Mutual information (MI) is used to measure the separation effect^[28,29], where the minimum MI corresponds to the best separation effect and correct $\theta$. MI can be represented as follows: $$\begin{gathered} \mathrm{MI}[S_p(\theta ),T_p(\theta )] \\ =\sum _{s\in S_p}\sum _{t\in T_p}\text{prob}(s,t)\log \left[\frac{\text{prob}(s,t)}{\text{prob}(s)\text{prob}(t)}\right], \end{gathered}$$(16)where $\text{prob}(s,t)$ represents the joint probability distribution function of pixels in $T_p$ and $S_p$. $\text{prob}(t)$ and $\text{prob}(s)$ are the marginal distribution functions of $T_p$ and $S_p$, respectively. As MI is invariant to the sign of the data, we only need to search for $\theta$ in the range of $0\le \theta \le \pi$. For example, in the scene presented in Fig. 1, the MI for different potential $\theta$ is plotted in Figs. 4(a) and 4(b), and the accurate value of $\theta$ is 2.90. After the $\theta$ is estimated, the polarization component of the backscattered light and $W$ can be calculated using Eqs. (7) and (14), respectively. At this point, the number of unknown parameters in Eq. (2) is reduced to four, which enables us to obtain the intensity information of the target exhibiting nonuniform polarization characteristics through matrix operations.

To verify the effectiveness of the proposed method, an underwater imaging experiment was designed, as shown in Fig. 5. The light from an LED (Thorlabs M660L4, LED driver model: LEDD1B T-cube) passed through a polarizer (Thorlabs LPVISE200-A) and became linearly polarized light. A glass water tank ($50\mathrm{cm}\times 70\mathrm{cm}\times 40\mathrm{cm}$) was filled with 100 L of tap water. Then, 80 mL of skim milk was mixed into the tank to simulate a scattering environment^[4]. The distance between the target and the detector was 45 cm (the calculated optical thickness was 1.512). A detector (Canon EOS 77D) and a light source were placed on the same side of the glass tank. A polarizer was installed in front of the detector and rotated to obtain polarization subimages of 0°, 45°, 90°, and 135°.

Figure 6 shows the experimental results obtained from using direct imaging and our polarization reconstruction method in the underwater environment. In Figs. 6(a)–6(d), the intensity images obtained directly by the detector experienced severe degradation. Furthermore, the presence of the backscattered light seriously reduced the image contrast and detection distance, whereas the influence of the backscattered light on imaging was significantly suppressed in the reconstruction results. In addition, the imaging field of view in Figs. 6(a) and 6(b) contained targets with different polarization characteristics. The experimental results showed that both the information sets of the two materials could be recovered simultaneously using the proposed method.

As an attempt to visualize the improvement resulting from using the proposed method, the curves in Figs. 7(a) and 7(b) were plotted, which represent the intensity profiles along the white, dotted lines crossing the middle of Figs. 6(a), 6(b), 6(e), and 6(f). In Fig. 6(a), due to the strong reflected light in the coin area and the influence of the backscattered light, the number on the coin cannot be recognized, and the curve in this area is extremely smooth with less detailed information of the target. Noticeably, the information provided in the coin area is observed in the detected image using the proposed method, and the numbers in the coin can thus be easily identified.

To analyze the degree of distortion in the reconstruction results obtained using different methods, the following simulation is conducted. A clear image without backscattered light is first acquired. Then, the total intensity image and polarization subimages are obtained by setting the intensity value of the backscattered light as well as the DOP and AOP of both the backscattered light (both were constants) and target-reflected light (both changed with the position of the pixels in the image)^[1,3,5,8], as shown in Fig. 8. Subsequently, we use these images to calculate the intensity images required by the traditional method.

Figure 9 shows the reconstruction results obtained by different methods. In order to facilitate the observation of distortions generated by different methods, we enhanced the details of the images in the second row. The traditional methods, such as image enhancement and Schechner’s method^[1], obtain good visual results; however, several distortions still exist between the restored target image and the simulated image. The images in the second line subjectively indicate that the proposed method exhibits the best fidelity. Furthermore, the degree of distortion estimated with different methods indicates that the proposed reconstruction method loses the least amount of target information.

In conclusion, this study utilized ICA to remove the backscattered light and accurately estimate the target information. Both the experimental and simulation results showed the effectiveness of our method. Compared with the other underwater polarization imaging methods, the proposed method could estimate the target information more accurately. In particular, for targets with different material properties, our method could accurately estimate their information simultaneously. Comparison of the degree of distortion estimated for the different underwater imaging methods showed that the proposed method produced the smallest distortion and could control the degree of distortion within 9%. This method is of great value in the field of underwater archaeology and shipwreck salvage where targets with nonuniform polarization characteristics need to be detected simultaneously.