Polarimetric Imaging Through Scattering Media: A Review

Imaging in scattering media has been a challenging and important subject in optical science. In scattering media, the image quality is often severely degraded by the scattering and absorption effects owing to the small particles and the resulting nonuniform distribution of the intensity or polarization properties. This study reviews the recent development in polarimetric imaging techniques that address these challenges. Specifically, based on the polarization properties of the backscattering light, polarimetric methods can estimate the intensity level of the backscattering and the transmittance of the media. They can also separate the target signal from the undesired ones to achieve high-quality imaging. In addition, the different designs of the polarimetric imaging systems offer additional metrics, for example, the degree/angle of polarization, to recover images with high fidelity. We first introduce the physical degradation models in scattering media. Secondly, we apply the models in different polarimetric imaging systems, such as polarization difference, Stokes vector, Mueller matrix, and deep learning-based systems. Lastly, we provide a model selection guideline and future research directions in polarimetric imaging.


INTRODUCTION
Optical imaging through scattering media, including turbid water [1,2], haze [3,4], fog [5], and biological tissue [6][7][8], enjoys a wide range of applications in areas such as underwater rescue [9], automatic driving [10], underwater archaeology [11], and biomedical imaging [12]. Therefore, the realization of clear visions in scattering media is of great interest and significance. However, the visibility and identifiability of the target scene are usually compromised as the radiance observed from a scene is scattered and absorbed by aerosols and particles existing in the environment [13]. The optical performance is thus limited in many practical applications [14]. In other words, the image quality captured by a camera deteriorates significantly, resulting in low image contrast [15], distorted color [16], and poor visibility [17].
Various dehazing or de-scattering techniques that have been developed to restore the image quality can be classified into two categories, non-physical and physical model-based methods, as shown in Figure 1. The non-physical methods, based on the image enhancement method, aim to highlight the target of interest and improve the contrast. The simplest non-physical method is the histogram equalization (HE) method, which enhances the overall image contrast by increasing the dynamic range of the gray value. Depending on the difference in the computing region of an image, HE can be divided into global HE (GHE) and local HE (LHE). The advantage of GHE lies in its higher efficiency and lower computation requirement, which is particularly suitable for enhancing excessively dark or bright images. This method does not fit images with high local brightness values, which often cause the "halo" effect. Such an issue can be addressed by applying the LHE method, for example, contrast limited adaptive histogram equalization (CLAHE) [18]. However, the "blocking effect" accompanied by increased computation complexity [19] cannot be avoided. Besides, HE methods may amplify noises during dehazing. Retinex-based algorithms [20,21] are models based on the color perception of human eyes, with the main concept of obtaining an object's reflection property to model the color invariance. Thus, the Retinex-based algorithms have been widely applied in the field of image enhancement to remove haze and scattered light. In addition, this method helps increase an image's contrast and brightness and regulate the dynamic range of its gray level. Nevertheless, the Retinex-based algorithms do not preserve edges well, which may lead to the halo phenomena in sharp boundary regions [21,22]. Finally, with regard to the dehazing technique, the frequency domain filtering (FDF) method proves to be another popular solution for image enhancement. More details can be found in the previous works [23][24][25].
The physical model-based dehazing and de-scattering methods are based on the knowledge related to the scene's physical features. To successfully restore an image, one of the key factors is to acquire an accurate depth of the scene [26], that is, the physical distance between the camera and target scene. However, the depth map is always unattainable for practical applications [27]. Therefore, a fundamental challenge in optical de-scattering techniques is to accurately estimate the depth map, that is, the transmission of scattering media [28]. Various methods have been proposed to overcome the challenge. For example, Fattal et al. [27] proposed a method that inferred the medium transmission map by estimating the albedo of the scene. However, such a method assumes that the transmission and surface shading are locally uncorrected and thus may fail when handling dense haze. Hautiere et al. [29] estimated the depth by determining the relationship between the road visibility and the contrast in the foggy image. Based on the analysis of the side geographical information obtained via an onboard optical sensor system, a 3D geographical model was established to remove the fog. Upon observing the property of haze-free outdoor images, He et al. [3] proposed the dark channel prior (DCP), based on the premise that "dark pixels" had a very low intensity in at least one-color channel except in the sky region. His method included three steps: air-light/scattered light estimation, transmission map estimation and refinement, and the final image reconstruction. Thanks to its effectiveness in dehazing, DCP has been adopted by most of the recent physical model-based techniques.
Almost all methods above are implemented based on the input of a single image, where certain assumptions or prior knowledge are necessary. The other physical model-based methods are based on multiple images corresponding to the same scene, that is, the images obtained under different visibility [30], images obtained with visible and near-infrared cameras [31], and images acquired with different polarization angles [15,32,33]. While images under different visibilities render the estimation of the depth map and scene structures to significantly enhance the image contrast, it remains challenging to handle real-time scenes [34]. Thanks to the excellent "long-distance transmission capacity" of the near-infrared light, the visible and near-infrared fusion methods improve the image quality by combining the rich color information of the visible image and the high visibility of the near-infrared image. However, the major obstacle is to acquire the visible and near-infrared images simultaneously, where expensive equipment and accurate optical alignment are both required. As opposed to the above methods, polarimetric imaging [35][36][37] is more effective because the scattered light is partially polarized [35] and the polarization information of the object and the turbid medium is different. Therefore, in principle, obtaining the polarization information of the scene and then processing them can effectively suppress the scattered light and extract the light coming from the object light [15,32,38]. A series of studies have shown that polarization-based imaging is a physical, low-cost, and applicable way to enhance the image quality, especially in highly scattering environments [2,39,40].
The typical polarimetric imaging systems include the polarization difference (PD) imaging [1,41,42], Stokes-based polarimetric (SP) imaging [43], and Mueller matrix (MM) imaging [44,45]. The PD imaging is based on two orthometric polarized sub-images to estimate the transmittance by analyzing the degree of linear polarization. The SP imaging, especially the full-SP imaging, leverages the robustness of the polarization angle [38] or the "memory effect" of circular polarization to achieve the backscatter removal [39,46,47]. The MM imaging benefits from its complete polarization characterization. These three basic models are built upon different optical systems that offer flexible options subject to different application requirements. Besides, the polarizationbased methods can be further improved by integrating with computer-vision-based and learning-based methods [2,15,48]. In other words, by introducing the polarization information into the traditional vision or learning-based method, greater application possibilities can be explored due to enhanced performance in image quality [2,[49][50][51][52].
This study first introduces the basic principles of the common polarimetric imaging models in scattering media and provides a comprehensive and up-to-date review for both traditional and advanced works. We explore topics that include the progress in model optimization and parameter estimation and the analysis of different methods from the perspectives of their limitations and potential solutions, with application-based recommendations for readers in optics and engineering communities. Section 2 introduces the imaging model through the scattering media and the related optical imaging systems. Section 3 demonstrates the polarimetric methods for imaging through turbid media based on different imaging systems. Section 4 provides a conclusion and an overlook for future development.

Physical Imaging Model in Scattering Media
The particles existing in the scattering media, such as in the atmosphere during hazy or foggy weather and under turbid water or sea, generally absorb and scatter light, resulting in the decay of image contrast, saturation attenuation, and color-shifting in the detected images [39,46]. Therefore, image dehazing or descattering plays an important role in various practical applications, such as traffic surveillance systems, security systems, object recognition, medical imaging, and remote sensing. Studies related to the image degradation mechanism, imaging systems, and recovery algorithms are receiving much attention.
Koschmieder et al. [53] proposed the first atmospheric scattering model, which was further modified by Narasiman and Nayar [54,55]. Based on their model, the image signal received by the camera is composed of two components: 1) direct transmission D(x, y), which represents the effect of scattering of light and the eventual decay of light before it reaches the camera, and 2) backscattered light A(x, y), which denotes the undesired backscattered lights from the particles in the object line of sight (LOS) [1,32]: As the light from the target progresses towards the camera, its energy is lost due to scattering and absorption. The fraction that does reach the camera is the direct transmission given by where z(x, y) is the distance between the object and the camera and depends on the pixel coordinates x and y; β is the attenuation coefficient; and L object (x, y) is the object radiance not scattered and absorbed along the LOS [1,35]. The attenuation coefficient is given by β, and the term e −βz(x,y) is also called the transmittance of light t (x, y).
A(x, y) denotes the undesired lights received by the camera mainly due to scattering by particles. It does not originate from the object on the LOS but varies with the horizontal distance by where A ∞ refers to the intensity value of backscattered light from infinity in the turbid medium. In most works, it is assumed to be a global constant independent of the location (x, y). Figure 2A shows the image formation and visual illumination components through the scattering medium. According to Eq. 2 and Eq. 3, we can observe that L object (x, y) and, thus, the recovered image could be obtained as far as the transmittance and backscattering are estimated accurately and the attenuation of the object light is compensated. Combining Eq. 1 and Eq. 2, one can recover L object (x, y) as follows: The currently reported polarimetric imaging-based dehazing methods are based on the above physical model and the scattered light's polarization property. The published results in various works show the high information restoration capacity and computational efficiency [56]. For the recovered image in Eq. 4, many quantitative criteria are used to characterize the quality of results, including the visibility range or distance [56,57], Michelson contrast (MC) [43], peak-to-correlation energy (PCE) [58], mean gradient [51], measure of enhancement (EME) [2,17,48], blind-reference-less image spatial quality evaluator (BRISQUE) [2], natural image quality evaluator (NIQE) [17], entropy [2], and peak signal-to-noise ratio (PSNR) [19,27,54]. Among the above criteria, PCE and PSNR describe the similarity between the restoration and clear image, EME and entropy describe the image contrast, and BRISQUE and NIQE quantify the distortion indicator of quality.
The superiority of the polarimetric methods also embodies the robustness of polarization parameters in various complex scattering conditions. Many works study the propagating light's physical characteristics and polarization properties in scattering media [59][60][61]. For example, with the aid of the Monte Carlo simulations, Xu et al. [62] demonstrated that the intensity of the polarized light after being transmitted underwater sharply decreases as the transmission distance increases, but the degree of polarization (DoP) of the transmitted lights remains above 0.75. It means that compared with the traditional imaging encoded with intensity information (in which the intensity will be significantly lost), the polarization encoding by DoP has overpowering advantages. Besides, Shen et al. [61] demonstrated that the depolarization behavior of light is sensitive to the mixing ratio or the distribution state of particles. Tao et al. [60] also found that the polarization properties provide additional information for the imaging, and the contrast of the polarization image can be significantly enhanced compared to the simplex intensity image in the turbid media. Moreover, the circular polarization images offer better contrast and higher visibility than linear ones under the same circumstance. All these reported results make the polarimetric methods and polarization control more promising for imaging in scattering media.
To date, whether it is on the basis of two images [15,32], three or four images [38,46], and nine or 16 images [44,45], various polarization recovery methods have been developed by modifying the basic model in Eq. 1. These methods are related to different polarization information, such as the PD, Stokes vector, DoP, angle of polarization (AoP), and MM. The following section introduces the related imaging systems and configurations.

Polarimetric Imaging Systems
Based on different polarization information, various polarimetric imaging systems have been developed. Early works in polarimetric imaging applications mainly focused on cases with linearly polarized light. This cross-polarization imagery, also called the PD imagery, was commonly used to enhance image contrast and minimize the blurring light in media with relatively low concentrations. The main element of this system is shown in Figure 2B1. In the setup, a light source, such as a lightemitting diode or a laser, is expanded with a beam expander. The light is polarized with a linear polarizer and hits on the target scene. A camera is positioned normal to the target scene to avoid most of the glare created by the interface. Finally, a rotating analyzer, always a linear polarizer, is placed in front of the camera to filter the back-reflected polarized light. By this configuration, two images (I , I ⊥ ) with orthogonal polarization states are acquired, from which one can calculate the degree of linear polarization (DoLP) as Stokes parameters have played a prominent part in the optical literature on polarized light [63,64]. As early as 1947, Chandrasekhar [65] used the Stokes vector to formulate the radiative transfer equations for scattering partially polarized light. Furthermore, the Stokes parameters give a complete description of any polarization state of light: where the first three parameters are linear components of the Stokes, while the last one is the circular component [49,63]. From this formalism, other parameters can be deduced, such as the DoP (P), the DoLP (P l ), the degree of circular polarization (DoCP), that is, P c , and the AoP (α). They are defined as follows: A more detailed description of Stokes vector can be found in the books about polarized light [63].
Unlike the PD imagery for linear polarized light, the Stokes vector contains the ellipticity of the beam. Hence, the complete imaging system for the Stokes vector requires an extension of the instrumentation. Optical retarders or wave plates (WP) are usually introduced into the system to generate or measure elliptical or circular states. Four intensity measurements are needed to calculate the complete Stokes vector parameters. Figure 2B2 shows a typical Stokes vector imagery, consisting of two sections: the polarization state generation (PSG) and the polarization analysis (PSA). In PSG, a WP and a linear polarizer are used to generate polarized illumination with an arbitrary state. The reflected intensity from the target scene is measured by adjusting the WP and/or polarizer's states in PSA. Based on these captured intensities, one can estimate the Stokes vector of the reflected light. In practice, the measurement of S 0 , S 1 , and S 2 is conducted by removing the WP in PSA. Only the last term, S 3 , requires this element to measure an elliptical/circular state. Moreover, with the development of nano-structures fabrication, a snapshot imaging solution has recently gained much attention using a division of focal plane (DoFP) polarization camera [66,67]. It can simultaneously capture four polarization angles from each video frame without image mismatch [68]. Besides, by adding a rotated WP, one can measure the full Stokes with only two shots [69]. The polarization camera makes real-time polarimetric imaging a reality and its applications possible.
The MM, proposed by Hans Mueller in the early 1940s [63,70], is another common parameter in addition to the Stokes vector in polarization imaging technology. The Stokes vector is a parameter describing the characteristics of the incident and the outgoing light when interacting with the materials, while the Mueller matrix is a "bridge" between the light and the material and describes the modulation of the incident light by the material [71,72]. The Stokes vector's description of a light beam requires four parameters. The modulation relationship between the incident and outgoing light can be fully described using a 4 × 4 matrix [73]. This matrix is called the MM. When a beam of light is incident on objects, the polarization properties of the reflected or transmitted light generally change [74,75]. Assuming that the Stokes vector of the incident light is S, the Stokes vector of the outgoing light after interaction with the medium is S ' , and the MM can express their relationship as follows: That is, In practice, the MM can be obtained using some possible PSG and PSA combinations in Figure 2B2. A straightforward implementation includes setting the light source at three linear states and the right circular state. Finally, 16 intensity measurements are needed to calculate the full 4 × 4 MM. By removing the WPs in PSG and PSA, the configuration becomes the same as that of the Stokes imager ( Figure 2B2). Then, by respectively rotating the directions of the two polarizers three times and obtaining a total of nine intensities, one can calculate a linear or incomplete MM with a size of 3 × 3.
In addition, according to the light sources' properties or the selected optical elements, the configurations can be categorized into different types as follows: (1) Depending on the existence of light sources, it can be categorized into "active illumination type" (with additional source) and "passive illumination type" (with nature light). (2) Depending on PSG's composition, it can be categorized into "unpolarized illumination type" and "polarized illumination type." In particular, the "polarized illumination type" can be further divided into "linearly polarized illumination" and "circularly polarized illumination." (3) To match the type of illumination, the PSA will contain a polarizer with "linearly polarized illumination" or a polarizer together with a retarder if it is "circularly polarized illumination." Choosing the configuration types depends not only on the polarization parameters being used but also on the environment in reality. For example, in the atmospheric environment, "passive illumination type" with sunlight is recommended, while in underwater, undersea, or low-light surroundings, the "polarized active illumination type" is preferred. In addition, for different target scenes, one needs to switch between linear and circular polarized illuminations. In the following section, we introduce well-established polarimetric methods in accordance with the above-mentioned configurations for imaging in scattering media.

Basic Model and Configuration
Inspired by the polarization-sensitive vision of some animals, PD imaging systems are proposed and developed to improve the visibility of objects in scattering media. This model is served as a common-mode rejection amplifier that can reduce the effects of background scattering and amplify the signal from targets where the PD magnitude is distinct from the background [41,76]. Based on the images captured for the same scene at two orthogonal linear polarization states (I (x, y) and I ⊥ (x, y)), the traditional PD system, proposed by Tyo et al. in 1995 [76], generates the PD and polarization-sum (PS) images as I PD I x, y − I ⊥ x, y , I PS I x, y + I ⊥ x, y , (10) where the PS image is equivalent to a polarization blind image obtained by a conventional imaging system. The PD image clearly depends on the choice of polarization axes, whereas the PS image does not. Such a relationship with the choice of axes can be used to minimize the effects of a partially polarized background in a PD image [42]. Notably, the scattered light is partially polarized and has the orthogonal and the same polarization components to the incident light at the same time. The performance of the PD method depends on the ratio of different components, which may be determined by the properties of scattering media, the incident polarization, the incident and observed angles, and so on.
Unlike Tyo's PD model, which is based on the theory of common-mode rejection, Schechner et al. [32,35] proposed a novel de-scattering method based on the atmospheric scattering model in Eq. 1, as shown in Figure 3. Similarly, this method captures two orthogonal polarized images composed of two unknown components (the scene radiance in the absence of haze and air-light). Because the air-light is usually partially polarized, these two images can be described by Frontiers in Physics | www.frontiersin.org March 2022 | Volume 10 | Article 815296 5 This model considers the natural polarization effects in imaging through the haze and builds the relationship of atmospheric properties, polarization properties, and imaging formalism. It does not require modeling the scattering particles' size or the precise scattering mechanisms. In general, such methods based on Schechner's model need to estimate two critical global parameters: air-light at infinity A ∞ and DoP p. Then, the air-light and the transmittance can be calculated bŷ respectively. Here, A ∞ and p rely on the choice of background or sky region Ω and are usually calculated by Finally, all these parameters are involved in the polarimetric model to obtain the recovered image as The PD-based methods share the same configuration shown in Figure 2B1. The polarizer in front of the camera or other intensity detectors is a must to capture two orthogonal polarized images. In contrast, the polarizer behind the light illumination source is optional and depends on the actual scene, active or passive illumination.

Representative Methods
As the first attempt, Tyo et al. [41,76] designed a special sample, an aluminum target containing two abraded patches with orthogonal directions, to verify the effectiveness of PD imaging systems. This sample is placed in a tank with inside dimensions of 30 × 30 × 15 cm 3 . This tank is filled with water and milk to simulate the scattering environments. PS and PD images after being transformed for optimal display are shown in   Figure 4A. We can see that the abraded patches are clearly visible in the transformed PD image but practically invisible in the transformed PS image. The key factor contributing to the enhanced visibility of the two patches in the PD image is the common-mode rejection feature intrinsic to PD. Based on a series of validation, Tyo et al. [76] found that the PD imaging model is quite sensitive to intrinsically small signals and possesses valuable qualities of being passive, simple, and potentially fast.
This method does not employ polarized illumination, which makes it constantly suffer from a reduced signal-to-noise ratio (SNR) because light reflected from the target and diffusive light backscattered from the turbid medium are superposed. Guan et al. [42] developed the traditional model and added a linear polarizer behind the light source to generate a polarized illumination. By measuring the co-polarization (parallel to the incident light's polarization state) and cross-polarization (perpendicular to the incident light's polarization state) images, they obtain an improved PD image, expressed as I PDI I − I ⊥ . Figure 4B presents the comparison between PS and PD images. The result shows that the proposed method significantly suppresses the background noise, and the image contrast is improved approximately 1.7 times [42].
However, the biggest challenge of the common-mode rejection-based method is that it cannot process scenes with complex conditions and objects because it does not consider the degradation mechanism in the scattering media. To solve this problem, Schechner et al. [32,35,77,78] combined both polarized optics and the typical atmospheric scattering model and showed how the polarization tool boosts the vision clarity in scattering media. An example in the air environment is shown in Figure 4C1. The dehazed image has much better contrast and color, especially in the distant regions of the scene noted as the green forest and the red roofs [32]. Unlike such applications in the air environment, which rely on natural sun illumination, active illumination is necessary for underwater applications. In 2004, Schechner et al. [57] developed their original model for the de-scattering in turbid water by introducing a polarized illumination, as shown in Figure 4C2. The corresponding results are shown in Figure 4C3. From the results, we may observe that the objects (i.e., the iron box) are well restored in both contrast and color. However, we must point out that their model is based on three assumptions: (1) Only the backscattering light/air-light is polarized, while the objects are unpolarized. (2) The total attenuation for objects at infinity also equals inhomogeneous haze; that is, t ∞ 0. (3) The multiple scattering (which affects the angular scattering distribution) is dominated by single scattering.
The first assumption does not apply to all practical cases because the object radiance could also contribute to polarization. When the depolarization degree of target objects is low, the light scattered or reflected by objects could contribute considerably to polarization. As a result, the previous methods may cause significant estimation errors. Huang et al. [59] found that the estimation produces negative values of t (as shown in Figures  5A-D) at the pixels corresponding to the low-depolarizing material if one assumes that the light emanated from objects in the scene is unpolarized. To solve this problem, they modified  Figure 5E shows that the recovered results reveal the details, which are not visible in the intensity images, regardless of the area corresponding to the objects having a high depolarization degree (plastic cube) or a low depolarization degree (metal coin).
To overcome the limitation caused by the second assumption and improve the recovery performance in complex conditions, such as that in the non-uniform optical field, Hu et al. [79] proposed a method of retrieving the radiance of the object based on estimating the spatial distributions of the DoP and the intensity of backscatter light by extrapolation fitting. For example, the progress of fitting spatial distribution is shown in Figure 5F, and the corresponding results are also shown for comparison. It shows that this method reveals the details of the scene (in orange and blue rectangles) that decay significantly in the intensity image, and it has better quality than the image recovered by Schechner's method because, under the non-uniform optical field, the DoP of backscatter A ∞ , represented by p, on the right side of the background region is considerably higher than that on the left side. Suppose we recover the image with Schechner's method, which considered the backscatter and its DoP to be constants; the transmittance can be considerably overestimated because of the improper estimations of A ∞ and p. In that case, it will lead to incomplete haze removal, thus a less clear recovered image [79].
The third assumption makes the traditional methods perform poorly in the case of dense haze [2,39]. To overcome this issue, Li et al. [2] proposed to combine the polarization-based model and the computational processing method, such as histogram stretching (HS). Their main idea is to reduce the density of haze computationally. Specifically, this method stretches the histograms of the orthogonal polarization images while maintaining the polarization relation in between. Based on the processed orthogonal polarization images, the recovered image with higher quality can be obtained by the traditional  Figure 6A presents an example in the dense turbid medium. The results show that the method significantly removes the scattered light and restores more details than Schechner's method. This method opens a door and verifies the feasibility and effectiveness of combining computer vision and polarimetric methods for image recovery in scattering media. Although the HS method and Schechner's polarimetric recovery method involved in this work are both old methods, the core idea of this combination has many perspectives. For example, Brosseau et al. [80] combined the low-pass polarization filter and DCP method and demonstrated the ability to significantly improve visibility and reduce runtime by a factor of about 50 for a 4K image. Figure 6B shows the natural imaging experiments at sea.
It is worth noting that the wavelength of the active illumination also affects the imaging performance because the scattered light caused by both water and particles is wavelengthdependent. Smith and Baker [81] measured the absorption coefficients and scattering coefficients of pure seawater on the wavelength ranging from 200 to 800 nm and validated the dependence of scattering on wavelength. Liu et al. [82] plotted the measured absorption and scattering coefficients of pure seawater in the visible-light range, as shown in Figure 6C. Results show that scattering decreases with wavelength increase. Based on this wavelength-dependent fact, they proposed a wavelength-selection-based polarization imaging method to image through highly turbid water with red light illumination. This method makes a good balance between range and vision and can turn targets from "undetectable" into "detectable," as shown in Figures 6D,E. In fact, as one of the earliest polarimetric imaging techniques, the PD-based methods have received a great deal of attention for image dehazing/de-scattering in scattering media. This makes them enjoy fast development in terms of scientific research and engineering applications. However, the PD-based methods only contain two polarized images, which means that the freedom degree of information is limited to two. A complete polarization characterization of the scattered light and objects is helpful in further enhancing the recovery performance. Naturally, possible solutions include capturing more images and obtaining polarization parameters with high information freedom degrees, such as Stokes vector and MM [64,83].

Basic Model and Configuration
As the information dimension (number of polarization sub-images) in Stokes vector configurations is higher (i.e., three for linear Stokes and four for complete Stokes) than two in the PD imaging system, the Stokes vector is more suited for characterizing polarization properties of scattered light. In contrast, as the two important parameters (AoP and DoP), which are highly relevant to the scattered lights' properties [66,84], can be directly deduced from the Stokes parameters, introducing Stokes analysis into the scattering removal is promising.
According to the fundamental model in Eq. 1, we know that the estimation of transmittance t(x, y) depends on the scattering section A(x, y), where the polarization property can be characterized by the Stokes vector S A (x, y). Here, S A 0 is the captured intensity corresponding to A(x, y). In this way, Eq. 1 can be rewritten as where S D (x, y) denotes the Stokes vector related to the direct transmission D(x, y). If the target objects are assumed unpolarized, S D (x, y) equals [D(x, y), 0, 0, 0] T ; conversely, S D (x, y) [D(x, y), S D 1 (x, y), S D 1 (x, y), S D 1 (x, y)] T , and the last three polarization parameters cannot be ignored. Based on the intensity measurements in different polarization states, we can estimate S A . Accordingly, we can obtain the A ∞ and t(x, y). In other words, the Stokes-based methods remove the veiling light by building the model between Stokes vectors (or the related polarization parameters) and (A ∞ , t(x, y)); that is, where the parameters in Eq. 16 (i.e., P, P l , P c , and α) are defined in Eq. 7. The basic configuration of the Stokes-based methods can refer to that in Figure 2B2, and the specific choice depends on the used algorithm and realistic environments. According to the number of the used Stokes parameters, the methods can be classified as the linear-Stokes (LS) and the circular-Stokes (CS) based methods. In particular, for the LS-based methods, the PSA contains a fixed linear polarizer to generate polarized illumination, while the PSG contains a rotating linear polarizer. By adjusting the polarizer in PSG at least three times, we can obtain polarization sub-images. Based on these images, we can obtain a 1 × 3 linear Stokes vector. In contrast, in the CS-based methods, both PSG and PSA contain a linear polarizer and a retarder (e.g., QWP or liquid crystal variable retarders). By adjusting PSG's states at least four times and obtaining the related intensity images, we can get a 1 × 4 complete Stokes vector. The difference between the two systems is whether the circular component of the Stokes vector (i.e., S 3 ) is considered. However, the CS-based methods do not depend on the scattering angle and can survive more multi-scattering events than linearly polarized light [47,85,86]. Figures 7A,B present the DoP as a function of optical depth in water with different diameter particles, from which we can observe that the CP light maintains better polarization characteristics than linear polarization (LP) [47]. This property, the so-called circular polarization (CP) memory effect, can be explained by the Mie scattering phase functions of LP and CP. The phase function in the circular case possesses a marked forward lobule that permits the photons to propagate around the beam axis with a higher probability than in the linear case [87,88]. Figure 1B presents the effects of size on the DoP. Results show that when the size parameter is large, the DOP for CP is greater than that for LP [85]. The characteristics in Figure 7B attest to the superiority of CS-based methods in dense turbid water or underwater environments with large-sized particles. Besides, Sara et al. [89] quantitively demonstrated the superiority of circularly polarized light in foggy environments. The experiments are carried out at CEREMA's 30 m fog chamber under controlled fog density conditions. Figure 7C compares the ratio of CO channels (probe the prevalence of the input polarized light through the fog) with respect to the intensity signal as a function of the radius of the intensity profile of the source image for circular (solid lines) and linear (dashed lines) polarizations at three visibilities. The results imply that circular polarization has a larger signal-to-noise ratio in transmission at deeper layers, whereas the signal from the linearly polarized light carries some noise due to its higher depolarization ratio when propagating in scattering media.

Representative Methods
To enhance the contrast, the two orthogonal polarized images in the traditional PD imaging model must be strictly selected to make the projections of the veiling light onto the two orthogonal axis directions equal [41]. However, this selection is time-consuming and inconvenient by rotating the polarizer mechanically, which is unsuitable for rapid imaging applications. Tian et al. optimized the traditional PD model to deal with this  limitation. They proposed a modified PD imaging method (M-PDI) based on the Stokes vector analysis of the veiling light [90]. The output image after removing the veiling light is expressed by where α is the polarization orientation angle of the veiling light. The optimal value of α corresponds to the highest image contrast. A linear Stokes vector is calculated by capturing three images in polarizer's directions of 0, 45, and 90°and searching the optimal α. The significant advantage of this method is that the recovery performance can be updated automatically by the computation program when obtaining the Stokes parameters, which makes the implementation of PD imaging ideal for rapid imaging. Figures 8A-D show the recovery results of this method. The result shows the background noise is significantly suppressed, and the contrast of the target "Z" is significantly improved. In 2018, Guan et al. [91] modified the above method and proposed another M-PDI method via the Stokes vector-based interpolation method. Figure 8E compares the traditional and the modified polarization filtering methods in PD imaging. The principle of this method is shown in Figure 8F [91]. From the image results, the object's contrast is significantly enhanced, and the background noise is significantly decreased. More details of this method and the comparison with the traditional PD method can be found in [91].
Although the M-PDI method has partially addressed certain inherent drawbacks of the traditional model based on Stokes analysis, the determination of crucial parameters, based on the computational searching, makes the model performance unstable and sensitive to noise [38]. To overcome this issue, Liang et al. [38] further explored the relationship between the Stokes vector and veiling light and estimated the backscattering/air-light using the AoP analysis. Based on the three captured images on different polarizer directions (i.e., 0, 45, and 90°), the Stokes vector is calculated by the expressions in [38] S 0 I(0) + I(90) where I(i) denotes the captured image when the direction of the polarizer is set to i degree. Then, we can calculate the intensity level of air-light as where θ 1/2actan(S 2 /S 1 ) denotes the AoP of air-light. Then, the output with a clear vision can be obtained by the typical physical model in Eq. 14. In this method, the noise in the sky area can be eliminated without any imaging-processing algorithm, which makes this method much more convenient and reliable in practical applications. Based on the technique, Liang et al. developed a series of algorithms to further enhance the recovery performance. For example, in 2015, they optimized the estimation of critical parameters (e.g., AoP and the intensity level of the air-light) to accommodate dense haze and achieved a 74% enhancement in the range of visibility (ROV) [43]. In 2016, as the infrared radiance has a better capacity for traveling through the haze, they modified the model by merging visible and infrared images. The ROV was thus improved by 100% [43]. In 2021, they introduced low-pass filtering into the AoP-based polarimetric imaging model and overcame the drawback of "noise sensitivity" in estimating the AoP value. The final imaging performance of these methods is shown in Figure 9. We can observe from the results that the faint information in hazy images is well preserved, and the contrast of the recovered image is increased significantly [92]. Such works based on AoP analysis accelerate the development of polarimetric imaging in specific scattering media. However, all these methods are based on the linear Stokes vector, which only includes three parameters and merely reflects the interaction of target objects with the linearly polarized light without considering the response to circularly polarized light. In particular, in some special scattering media requiring active illumination, the circularly polarized light tends to maintain its original polarization property better than the linearly polarized light, namely, the "circular polarization memory" effect [93,94]. Therefore, using circularly polarized light can improve the recovery performance in dense turbid media more than linearly polarized light. Therefore, Hu et al. [46] proposed a new polarimetric image recovery method in dense turbid media with the illumination light featuring circular polarization. In this method, the active illumination is modulated by a polarizer and an additional retarder to generate circularly polarized light. The estimated Stokes vector is decomposed into linear polarization, circular polarization, and un-polarization parts as follows: where the subscripts of l and c indicate the linearly and circularly polarized light, respectively. According to this decomposition, the linear and circular components of the veiling light are removed by solving the DoLP and DoCP of the backscattering. Figure 10A shows the difference between these two degrees of polarization, while the circular one is always missed in most methods. This is the first polarimetric imaging system and algorithm considering circularly polarized lighting. The recovered result and its comparison with the linear one are as shown in Figure 10B and further evidenced by experimental results for the scenes with different polarization properties, for example, a rough wooden board with patterns and words on its surface and the non-flat plastic toy. With the improvement in the theoretical model, the trend of Stokes-based polarimetric imaging in scattering media continues to optimize the estimation of key parameters and render more accurate values, such as the polarization information of object and scattering signals. For example, Jin et al. [52] proposed a scattering removal method from the perspective of global estimation of polarization information to realize polarimetric calculation of global pixels for automatically estimating the DoP. This global pixel calculation is accomplished by utilizing the gradient prior information of the total intensity image. Figure 10C shows the flowchart of this method, and ( Figure 10D) shows the comparison results of three classical recovery methods and this method in different turbidity [32,41,82]. Apparently, as opposed to methods with the assumption of constant DoP of target light, this method can retrieve the DoP of the target light of each pixel in the image. As a result, it has a better performance of recovering the scene's details even though the water is turbid [52].
In short, the Stokes-based methods outperform the PD-based methods because the Stokes vector renders more useful polarization information, such as DoP and AoP, which are closely related to the backscattering/veiling light caused by the existing particles. Improving the information dimension and estimating key parameters are two effective ways to boost  imaging performance. However, we must remember that, unlike MM, the Stokes vector is not a complete polarization characterization [63,95]. Therefore, the trial based on the MM has become the preferred approach to improve the information dimension, an exciting solution to enhance imaging performance in scattering media further.

Basic Model and Configuration
Unlike the Stokes vector, which is usually used to characterize the polarization properties of the light beam, the MM contains all the polarization information of the target materials. Therefore, the MM is always applied to distinguish different materials [8,[96][97][98] and 3D physical imaging [99][100][101].
According to the basic model in Section 2, the received signal or images can be classified into the target and scattered signals. The Stokes-based methods build the relationship between the Stokes vector and the scattered light. However, they fail to distinguish objects with different polarization properties. These methods are unsuitable for objects whose reflected light has a similar Stokes vector with the backscattering light.
With the information provided by MM, the operational space for polarization image processing is greatly improved, making it possible to distinguish objects in a high degree of freedom for the polarization information [102]. In addition, the MM-based configuration opens the door to modulate the incident illumination with a significant polarization space. In other words, the MM-based scattering suppression imaging technology modulates the illumination and analysis parts simultaneously. The basic configuration of MM-based methods is as shown in Figure 2B2, and one can remove the WPs in accordance with the type of MM (incomplete or complete) that is desired.

Representative Methods
The MM-based polarimetric suppression method for the imaging in scattering media is mainly based on modulating the polarization state of the active illumination. The earliest attempts focused on imaging linear MM and required nine intensity images [44,45,103]. In the related configurations, the WPs in PSG and PSA in Figure 2B2 are removed, and the polarizers in PSG and PSA are rotated to three different positions to capture nine images, respectively.
In 2019, Guan et al. [44] found that the illumination polarization angle and the MM difference between the medium and the object could affect the SNR of recovered results obtained by the rotation orthogonal polarization imaging method. They designed a linear MM-based polarimetric method to precisely control the illumination polarization angle and achieve a rapid imaging process. In 2021, Wang et al. [45] plotted backscattered light and target reflected light in the point cloud diagram in Figure 11A by establishing a differential imaging model. According to Figure 11A, on the premise that the angle between the polarizer's direction and the backscattering light's polarization direction is 45°, the backscattering light can be removed by the PD method to achieve scattering suppression. By analyzing the principle of polarization difference, the output result under ideal conditions is given by I out x, y I PB x, y − I PB⊥ x, y + I PT x, y − I PT⊥ x, y where I PB (x, y), I PB⊥ (x, y), I PT (x, y), and I PT⊥ (x, y) represent the horizontal and vertical projection of the polarized part of the backscattering and target light, respectively. It can be seen from the above formula that the output result of the differential image will be affected by the polarization angle between the target and backscattered light. Therefore, by modulating the incident light and changing θ, the performance of the traditional PD method can be improved. Figure 11B presents the comparison results of this method with the traditional PD method and other classical methods. We see that the modulation of incident light has a significant influence on polarization differential scattering suppression. In other work, the authors disregard the traditional PD method and instead directly process the image according to the MM 45 . By modulating the polarization state of the incident active illumination light, the DoP of the backscattered light is directly maximized, thus achieving the best suppression of the backscattered light. Based on the configuration in Figure 2B2, one can obtain any specific polarized incident light by adjusting the PSG. The Stokes vector of the incident light S in can be expressed by its DoP value P, azimuth α, and ellipticity ε as follows: where S in 0 denotes the intensity of the illustration light and S [ S in 1 , S in 2 , S in 3 ] is the normalized Stokes vector. When the DoP of backscattered light P back is higher, more backscattered light can be blocked by the PSA. Therefore, to suppress the backscattered light optimally, one must maximize the DoP of the backscattered light by choosing an optimal set of azimuths and ellipticities (α opt , ε opt ). This optimization problem can be expressed by where P back (M, α, ε) 1 (24) Figure 11C presents the MM image of the target scene, and the corresponding optimal set of azimuths and ellipticities can be found by a global search shown in Figure 11D. The imaging results and the comparison with the traditional PD method are shown in Figure 11E. The results demonstrate that this method is stable and can be implemented with any digital image processing to achieve a more scattering suppression performance.
In 2022, Liu et al. [104] developed an MM-based de-scattering method and introduced the depolarization (Dep) index into the de-scattering algorithm. Dep is derived from the MM and is defined by By studying the backscattering distribution in different NTU turbidities, it is found that the background intensity correlates linearly with Dep in a remarkable way, as shown in Figures  12A-B. Therefore, Dep is used to characterize scattering media. By quantifying the light attenuation with the transmittance map, a clear vision of targets can be recovered using the information of scattering media. An example of recovery performance is shown in Figures 12C1,C2. The results demonstrate that the image contrast is significantly improved after recovery. In particular, the paper stripe and metal ruler are both clearly visible. From the zoomed-in view in Figures 12C3,C4, the ruler in the intensity image blurs, especially the tick mark and edges. In contrast, the edges of the ruler are visible after recovery, even with distinguishable tick marks.
The MM-based methods have many advantages. For example, they provide more helpful information by increasing the degree of freedom for polarization and making a clear differentiation among objects with similar intensity appearance but different polarization properties. However, it should be noted that a very scarce amount of works, less than five to the best of our knowledge, has placed the focus on MM-based de-scattering. Indeed, making full use of the MM decomposition and other MM-related parameters is promising and needs more attention.

Learning-Based Polarimetric Imaging
As discussed in Section 2, one of the dehazing methods is based on a physical model where prior knowledge is applied to extract physical parameters related to the scattering media and then recover the targeted signal. In this case, the estimation accuracy of these key parameters determines the final performance. Therefore, almost every developed method strives to optimize the accuracy of parameter estimation to make it as close as possible to the physical values in the scene. However, the optimizations come at the cost of increased computation complexity and reduced universality. In contrast, the methods not directly based on a physical model aim to improve the image quality by enhancing the image contrast and the difference between different image structures. This kind of method handles different scenarios indiscriminately but can perform poorly, especially with complex conditions. The deep learning-(DL-) based method is data-driven and thus capable of extracting the hidden relationship and physical properties between/in the raw and target data. This makes it a promising choice for de-scattering. In some ways, the traditional "end-to-end" network corresponds to methods that are not based on the physical model, while the "physical embodiment" network is based on physical models. The following section introduces learning-based polarimetric imaging in scattering media from the basic concepts to wellestablished applications.

Basic Concepts of the DL and Neural Network
In 1943, psychologist Warren McCulloch and mathematical logician Walter Pitts proposed the concept of the artificial neural network and the mathematical model of the artificial neuron, thus prompting the era of artificial neural network research [105]. They abstracted the entire working process of neurons into the model shown in Figure 13A. In the model, input data x i are given to the network processed with weights and bias parameters. Then, a nonlinear activation function φ(p) is applied to obtain the final output Y. The whole process is called forward propagation. The weights w i and bias b are the parameters to be learned, which can be seen as the memory of the neural network. The final output is the prediction result obtained by the network according to the input data, which may be different from the ground truth. Therefore, it is necessary to calculate the deviation between the predicted result and the ground truth, or the network loss, to update network parameters. The process of updating weights is usually via some variant of a gradient descent algorithm. The training process of the network repeats the forward propagation and backpropagation process until the loss is minimal so that when we put in input data, we obtain an output nearly the same as the ground truth [106].
To boost the network performance and suitability for various tasks, different advanced network structures have been developed, such as LeNet [107], AlexNet [108], GoogleNet [109], ResNet [110], and DenseNet [111]. In the field of computer vision, learning-based solutions have become the hottest topic. Particularly, various learning-based works focus on improving imaging quality in scattering media. For example, Chen et al. [112] proposed an "end-to-end" dehazing network. In their review, a generative adversarial network (GAN) is used to realize end-to-end image dehazing. The work focuses on solving the problem of grid artifacts and has greatly improved the indexes of peak signal-to-noise ratio (PSNR) and Structural Similarity Index Measure (SSIM). In 2019, Pan [113] proposed a physics-based feature dehazing network for image dehazing network. In contrast to most existing end-to-end trainable network-based dehazing methods, they explicitly considered the physical model of the hazing process in the network design and removed haze in a deep feature space. However, all these descattering methods are based on a single intensity image. In recent years, DL has been successfully applied to polarimetric imaging [49,66,114]. Such works develop mainly from two parts: polarization acquisition and polarization processing. The related applications include denoising, dehazing, image fusion, targets detection and classification, and super-resolution (a brief classification is shown in Figure 13B). In the following, we focus on the learning-based de-scattering works.

Representative Methods
The currently popular neural network ResNet was first proposed by He et al. [110] in 2015 and is widely used in many scenarios of DL. This structure can solve the problem of gradient disappearance when the network is deep. Huang et al. [111] proposed DenseNet in 2017, which can alleviate gradient disappearance, strengthen feature propagation, encourage feature reuse, and significantly reduce the number of parameters through the dense connections between layers. Based on these two networks, Hu et al. [50] first proposed a polarimetric dense network (PDN) and applied it to underwater polarization image restoration. The network structure of PDN is shown in Figure 13C1. This network includes three main components (shallow feature extraction, residual dense block, and dense feature fusion) to deeply extract shallow features of polarization information from three polarization images and then fuse them. The loss function used is defined as where I pred i (x, y) and I gt i (x, y) refer to predicted and ground truth images, respectively, with their polarization information. Because the polarization is considered, the recovered image in Figure 13C2 has more detailed features than the intensity image used alone. Besides, there are more artifacts in the recovered result by "intensity-Net" than "Polarimetric-Net." These results demonstrate that embedding polarization information and constraints into the network helps improve performance.
In 2021, Zhang et al. [115] studied how to optimize the network structure and loss function to improve the suitable model performance. They found that, by adding polarized information along with the light intensity information to the model at the very front of the model structure, a betterrecovered image can be obtained. The model structure proposed can be used for image recovery in turbid water or other scattering environments. It should be noted that both the above methods are "end-to-end" and depend on paired training data. Although the application of neural networks has significantly improved image de-scattering performance compared with traditional methods, its disadvantages are also easily recognizable. DL networks, especially end-to-end networks, have poor interpretability as the neural networks are more like a black box. It is difficult to explain and understand the inner operation process. To solve this problem, Ren et al. [116] integrated the polarimetric imaging model with the DLbased method and proposed a lightweight network structure, which can restore underwater images with different turbidity. This method makes the image restoration through neural networks more in line with the physical principle and achieves good results.
To complete the physical meaning in the network and reduce the degree of freedom of the model, Zhou et al. [117] proposed a generalized physical formation model of hazy images and a robust polarization-based dehazing pipeline without the assumptions in traditional polarimetric methods. The designed network includes sub-networks to estimate different parameters, as shown in Figure 14A. The network divides the whole image de-scattering process into two steps. At the end of each step, semantic and contextual information is used to refine the output of the corresponding sub-network. This method provides a new perspective for the fusion of physical models and neural networks. However, the generalization ability on the real dataset is still limited for utilizing computer-synthesized datasets. In order to get rid of the paired data's dependence and make the learning-based methods applicable in practice, the unsupervised model-based and untrained model-based solutions are proposed. For example, Yang et al. [118] designed an end-to-end unsupervised generative network to remove the backscattering light, as shown in Figure 14B1. This method produces an adversarial loss with the discriminative network to improve the performance. In addition to using GAN to remove the backscattering light, they also modified the underwater imaging model based on several physical priors. The DoP of backscatter is the same as that of background light. This new model can be applied in a variety of nonuniform optical fields. Figure 14B2 shows its recovered results by different methods. Besides, they also verified that this unsupervised solution could adapt to the non-uniform optical field with different incident angles. In 2021, Zhu et al. [119] proposed a non-GAN unsupervised method by combining the polarization physics model and DL technology. Figure 14C1 presents the network's architecture. Rather than using atmospheric scattering model directly, they input the polarimetric hazy images into U-Net to obtain the corresponding de-scattered images, added haze to the output of the network through the model proposed by Liang et al. [43], and finally calculated the loss between the generated hazy images and corresponding captured images. Figure 14C2 presents a visual comparison among different de-scattering methods. From the results, we may observe that the background area with the homogeneous scattering effects is removed using this method, but the object information is preserved. In short, the unsupervised image de-scattering through U-Net does not need paired datasets or even haze-free images.
DL models and research largely depend on a particular dataset, and it is hard to guarantee similar performance from other datasets. This is where this method falls short if compared with other common traditional models. Nevertheless, we firmly believe that DL techniques hold a crucial place in this field.

Polarimetric Imaging Through Scattering Tissues
In the above sections, we have reviewed the polarimetric imaging through such scattering media as fog, haze, and turbid water. Meanwhile, biological tissues as another important scattering media and the related polarimetric imaging techniques have gained great attention in the biomedical field. To be more specific, the biological tissues contain fiber-like macromolecules (e.g., the collagen fibers in the skin and tendons, muscle fibers, and myofibrils in skeletal muscles), which exhibit a certain degree of structural anisotropy and anisotropy in dielectric response. These properties manifest themselves via birefringence [125,126], which can be observed using polarization measurement or polarimetric imaging [127]. For example, many works have shown that the depolarization, retardance, and diattenuation induced by the birefringent tissues can be considered indicators to assess macromolecules' microstructure, thus being conducive to diagnosis and the study of pathological alterations [128,129]. However, the scattering (especially the multiple scattering) in thick tissues often results in the depolarization of light, which makes detection of the remaining information-carrying polarization signals challenging. Therefore, various polarimetry or polarimetric imaging techniques have been developed to maximize measurement sensitivity to assist in analyzing useful tissue information [97,127,130,131].
In fact, the polarimetric imaging through the biological tissues often shares the same basic polarization configurations as mentioned in Section 2: polarization difference, Stokes vector, and Mueller matrix polarimeters. However, the focus and methods of the related research are significantly different from those in Section 3. Here, the main focus is to study the properties of scattering media, namely, the tissues themselves. For example, the emphasis will be placed on modeling the polarized light transport and the depolarization of multiply scattered light in tissues by both Monte Carlo simulation and the real experiments [132][133][134][135] or the study of the mechanism of depolarization and its dependence on the different tissue or media parameters (e.g., the density, size, distribution, shape, and refractive index) [131,[136][137][138][139][140]. On the contrary, for the reviewed physical degradation mode-based polarimetric recovery methods in Section 3, the properties of scattering media do not directly impact the recovery performance. This is because these techniques aim to remove the scattered light (i.e., A(x, y) in Eq. 1) and recover the direct transmission (i.e., D(x, y) in Eq. 1). The critical step is to calculate media transmittance (i.e., t(x, y)) and the intensity level of scattered light (i.e., A ∞ ) by estimating polarization properties (i.e., DoP and AoP). In short, the recovery performance mainly relies on the accuracy of the estimation polarization properties. More details can be found in the related works in Refs. [130,141]. In this review, we have presented an overview of the polarimetric imaging methods through scattering media from the perspectives of the basic model, imaging system, and representative works. Table 1 provides a brief summary and comparison across these methods.
Thanks to the property of polarized light propagated in scattering media, polarimetric methods outperform traditional intensity-based methods, particularly under complex conditions such as high-density turbid media, non-ideal illumination environments, and scenes with multi-material objects. We have demonstrated that the increase in polarization information dimension can constantly improve the imaging performance and the polarimetric methods. However, the complexity in both the imaging systems and computations is inevitably increased. It is worth noting that to achieve the balance between performance and complexity, advanced optical equipment and innovation in imaging systems must come into play. In addition, the demand for practical applications will certainly drive the development of polarimetric imaging methods. As such, we propose several topics of interest for future studies:

Multispectral-Polarization Systems and Fusion Algorithms
The combination of multi-spectral and polarization is often applied in the field of remote sensing [142,143]. They found that more complex and accurate indices and models can be developed to reveal more information when the polarization FIGURE 15 | (A) Images by LWI-DoFP camera (reprint from Elsevier: ISPRS Journal of Photogrammetry and Remote Sensing [146], copyright 2021). (B) Polarimetric imaging in the field of self-driving [153]. Different systems of polarimetric underwater imaging were developed for realistic applications. (C) Schechner's imaging system (reprint from IEEE: Transactions on Pattern Analysis and Machine Intelligence [35], copyright 2008). (D) FOREEEA's waterproof system (reprint from OSA: Optics Express [80], copyright 2019). (E) Full Stokes real-time dehazing system (reprint from OSA: Applied Optics [34], copyright 2017). (F) DoFP polarization camera (reprint from OSA: Optics Express [49], copyright 2020).
information is supplemented [144,145]. Therefore, combining multi-spectral and polarization into an integrated imaging system and fusing multi-parameters are two possible directions to enhance the image quality in complex scattering media.
As mentioned in Section 3, the long-wave light is more robust for transmitting through the haze than visible light. Another powerful example is using a long-wave infrared DoFP polarization in road detection, as shown in Figure 15A [146]. The fusion of short and long wavelengths effectively increases the visual range at the cost of decreased resolution. Therefore, an optimal tradeoff is helpful for this solution [147,148]. Besides, in the realistic underwater scene, seawater appears in different colors because the scattering and absorption depend on both the wavelength and the physical properties of the particles. It seems that a tunable wavelength source will make the imaging system have a wide range of application scenarios. On the contrary, the existing polarimetric methods mainly focus on recovering information only related to the intensity, that is, "to see it." However, the most significant advantage of polarization is that it can see what human eyes cannot (in the intensity condition), such as information related to DoP [149], AoP [150], or the index of polarization pure (IPP) [151,152]. These parameters help distinguish different materials and different optical responses. For example, the DoP and AoP have been successfully applied in self-driving. They can provide additional information in complicated meteorological environments (fog and rain), as shown in Figure 15B [153].
Fusing multiple polarization parameters and intensity into a frame, further increasing the information dimension, may open a new door for the imaging in scattering media and take a transform from "could be seen" to "see far" and from "see clearly" to "see more."

Real-Time, Real-Scenario, and Robust (3R) Solution for Polarimetric Imaging
Most reported polarimetric imaging methods for scattering media have been demonstrated under laboratory conditions. Although some were implemented in real-life scenarios, the imaging process was static, and the analysis was completed afterward. To achieve real-time, real-scenario, and robust (3R) solutions for polarimetric methods, advanced algorithms for optical information processing need to be developed and an improved collaborative imaging system is also required.
Schechner et al. first developed an underwater polarimetric imaging system, as shown in Figure 15C. This system is designed based on their basic PD imaging model and performs better in real-life scenarios by combining necessary post-processes [35]. In 2019, Khadidja et al. used the waterproof imaging system designed by FORSSEA Robotics Company, as shown in Figure 15D, to carry out the experiments under more realistic conditions [80]. Compared with the preceding imaging system, their setup is improved and can be integrated into underwater detectors, such as the underwater robot and autonomous underwater vehicles. However, these two systems are based on linear polarized light, and there are only two images with orthogonal polarization states. Zhang et al. designed an aperture-division polarimetric camera, as shown in Figure 15E, to capture four polarization images in the atmosphere via methods based on the full Stokes vector. It successfully achieves real-time image haze removal with an output rate of 25 fps [34]. However, the dehazing performance can be significantly affected by the registration accuracy. With the development of the DoFP polarization camera, as shown in Figure 15F, real-time processing is made possible without registration error. Our team has integrated this DoFP camera with a watertight device and an adjustable polarized illumination to create the underwater polarimetric imaging system (UPIS). The corresponding configuration is shown in Figure 2B3 to perform de-scattering (in air and undersea) tasks in real-life scenarios. Using custom dehazing algorithms, the visual range is increased by about 8.5 times, and the processing speed reaches 15-25 fps for images with a resolution of approximately 1000 × 1000 [154]. In other words, this custom-made system has fulfilled the "3R-criterion" to a much greater extent. However, compared with traditional solutions with a single image and advanced computer-vision-based algorithms, there remains room for further improvement in the 3R-polarimetric imaging system.
A possible approach going forward is to further exploit integrated polarimetric imaging systems based on the DoFP polarization camera and automatic rotating devices in PSG and PSA [155][156][157][158][159]. As such, the polarization modulations in both illumination and detector can be controlled simultaneously to handle multifunctional applications.
In the review, we have covered some, if not all, of the established works in the fields of polarimetric imaging in scattering media, and we are keen to see more studies in the area as further understanding in polarimetric imaging will undoubtedly benefit both the academic and industrial communities in a significant way.

AUTHOR CONTRIBUTIONS
XL, YH, and HW prepared the references and data. XL wrote the manuscript with input from all authors. TL reviewed and improved the writing. S-CC and HH supervised the project.