Next Article in Journal
Binaural Acoustic Scene Classification Using Wavelet Scattering, Parallel Ensemble Classifiers and Nonlinear Fusion
Previous Article in Journal
Evaluation of Isomotive Insulator-Based Dielectrophoretic Device by Measuring the Particle Velocity
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Residual Interpolation Integrated Pixel-by-Pixel Adaptive Iterative Process for Division of Focal Plane Polarimeters

MOE Key Laboratory of Optoelectronic Imaging Technology and System, Beijing Institute of Technology, Beijing 100081, China
*
Author to whom correspondence should be addressed.
Sensors 2022, 22(4), 1529; https://doi.org/10.3390/s22041529
Submission received: 22 January 2022 / Revised: 8 February 2022 / Accepted: 9 February 2022 / Published: 16 February 2022
(This article belongs to the Topic Advances in Optical Sensors)

Abstract

:
Residual interpolations are effective methods to reduce the instantaneous field-of-view error of division of focal plane (DoFP) polarimeters. However, their guide-image selection strategies are improper, and do not consider the DoFP polarimeters’ spatial sampling modes. Thus, we propose a residual interpolation method with a new guide-image selection strategy based on the spatial layout of the pixeled polarizer array to improve the sampling rate of the guide image. The interpolation performance is also improved by the proposed pixel-by-pixel, adaptive iterative process and the weighted average fusion of the results of the minimized residual and minimized Laplacian energy guide filters. Visual and objective evaluations demonstrate the proposed method’s superiority to the existing state-of-the-art methods. The proposed method proves that considering the spatial layout of the pixeled polarizer array on the physical level is vital to improving the performance of interpolation methods for DoFP polarimeters.

1. Introduction

Polarization, amplitude, wavelength, and phase are the four most important physical characteristics of light. Polarimeters can obtain the intensity (amplitude) and polarization information of the target scene to calculate polarization parameters such as the Stokes vector, the degree of linear polarization (DoLP), and the angle of polarization (AoP). Subsequently, the target contrast enhancement, stealth target detection, and deduction of characteristic information such as surface shape, roughness, and spatial attitude can be achieved. Polarization imaging technology is, therefore, extensively used in target detection and classification [1,2,3], three-dimensional shape reconstruction [4,5,6], remote sensing observation [7,8,9], and medical biological imaging [10,11].
The increasingly mature nanomanufacturing technology and the urgent need to simultaneously detect polarization information promote the rapid development of miniaturized and compact division of focal plane (DoFP) polarimeters [12,13,14,15]. Companies such as FLIR [16], 4D Technology [17], and LUCID Vision Labs [18] have successively launched DoFP polarimeter products that can be used for precision measurement. This polarimeter integrates a CCD/CMOS sensor and an aluminum nanowire polarizer filter array with a similar pixel structure, as in the imaging focal plane array (FPA) (Figure 1a). The output of this polarimeter is an incompletely sampled mosaic image with four polarization channels of 0°, 45°, 90°, and 135°. Each channel corresponds to a different instantaneous field of view (IFOV) due to the spatial dislocation arrangement. When the polarization information of these four channels is directly used to calculate the polarization parameters, the spatial resolution of the calculated polarization parameter image is reduced to 1/4 of that of FPA. Further, errors (such as the sawtooth effect) will be present in regions with abundant edge and texture. These phenomena form so-called the inherent IFOV error of DoFP polarimeters.
Reducing the IFOV error and restoring the spatial resolution is generally achieved by interpolating the output mosaic image of cameras using demosaicing methods. DoFP polarimeters shares the same principle of division of focal plane as the Bayer color camera. Therefore, the IFOV error formation mechanism of these two cameras is the same. Research into demosaicing methods for DoFP polarimeters usually refers to the earlier color image demosaicing methods [19,20,21,22,23,24,25,26,27]. In recent years, several color image demosaicing methods have been effectively transferred to DoFP polarimeters. These methods can be classified as:
(1)
Methods of independently interpolating using single-channel, which mainly include polynomial interpolation methods (bilinear [28,29,30,31], bicubic, bicubic spline [32,33], etc.) and edge directionality interpolation methods (gradient [34,35], smoothness [36], etc.). they are easy to implement, but their performance is mediocre.
(2)
Methods of interpolating using other channels as reference images, which mainly include correlation-based interpolation methods [37,38,39,40] and residual interpolation methods [41,42,43,44,45]. They are balanced in performance and stability and are the main topic of this paper. Recently, some heuristic algorithms (e.g., heuristic validation mechanisms) have been shown to find some important regions in traditional images [46], and they are expected to be combined with the residual interpolation algorithms to further improve interpolation performance.
(3)
Learning-based methods, which mainly include optimization-based methods [47], sparse representation-based methods [48], and deep learning-based methods [49,50,51]. They are considered to have the best performance on the published datasets, but their algorithm designs do not directly correspond to the DoFP polarimeter model, and the current open-access datasets contain very limited polarization scenarios.
The residual interpolation demosaicing methods can utilize the similar edge and texture features of the four channel images. Two channels are selected as the input image and the guide image to generate the initial estimate using the guide filter, and interpolation is executed in the residual domain containing less high-frequency information (where the residual is the difference between the observed image and the tentatively estimated image). This method has been proven to have a better demosaicing performance than other traditional polarization interpolation methods [41,42,43,44,45].
Two residual interpolation methods to demosaic DoFP polarimeters, with minimized residual (PRI) [41,42] and minimized Laplacian energy (MLPRI) [43], have been reported. However, these two methods did not thoroughly consider the inherent differences in the spatial sampling modes of DoFP polarimeters and Bayer color cameras. The following problems are present in the selection and preprocessing of the guide image:
(1)
The spatial layout of four channels in the pixeled polarizer array is not thoroughly considered when selecting the polarization direction of the guide image in PRI and MLPRI. In the color filter array, the sampling rate of G channel is 50%, which is twice that of the R and B channels (Figure 1b). Therefore, the G channel is usually interpolated first, and its interpolation result is also used as a reference image when interpolating the R and B channels, which makes the performance of residual interpolation methods better than the performance of traditional single-channel interpolation methods. In contrast, the sampling rates of the four channels in the pixeled polarizer array of DoFP polarimeters are equal, so there is no specific dominant direction. The existing PRI or MLPRI intuitively selects the same channel as the input image or the 0° channel as the guide image. The selected guide image does not have an advantage in terms of sampling rate, which makes the improvements in the performance of PRI and MLPRI insignificant compared with the single-channel interpolation methods.
(2)
The guide filter requires the guide image to have the same high resolution as the interpolation result. High-resolution images cannot be directly obtained during the actual polarization imaging process. Therefore, the guide image is usually generated by preprocessing the low-resolution observed image. Referring to color image demosaicing, PRI and MLPRI use basic interpolation methods, such as bilinear and bicubic interpolation, to up-sample the observed image and generate the guide image. However, when the sampling rate of the observed image is low, the guide image generated by this preprocessing may exhibit large errors in regions with an abundant edge and texture. This error will be transmitted to the tentatively estimated image, and further affect the quality of the output interpolation result.
Looking at this problem, this study proposed a residual interpolation method integrating a pixel-by-pixel adaptive iterative process for DoFP polarimeters (PAIPRI). Compared with the previously published PRI and MLPRI, innovative designs of the proposed method have been carried out, focusing on the following aspects:
(1)
We proposed a new guide-image selection strategy. We considered the spatial layout of the pixeled polarizer array, and chose different channels as the guide image for the pixels in different spatial positions. In addition, cooperating with the different sizes and directions of the filter window, the sampling rate of the adopted guide image in the filter window increased to 50%.
(2)
We designed a pixel-by-pixel adaptive iterative process based on residual interpolation. The guide image and the interpolation result were adaptively updated pixel-by-pixel through two interrelated iterative processes to improve the demosaicing performance of the output interpolation result.
(3)
We performed an adaptively weighted average fusion on the local iterative optimal results of the two guide filters, and minimized residual and minimized Laplacian energy, to make the interpolation results better.
(4)
Unlike the current mainstream learning-based methods, our algorithm is completely physical-fact-based and can explain the down-sampling process of the DoFP polarimeter. Furthermore, the focus on the improving imaging system makes our algorithm completely independent of the polarized images being processed, making it more robust to unseen scenes.
We conducted comparison experiments using both the open-access dataset images collected by a division-of-time polarimeter and the indoor and outdoor scene images collected by a real-world DoFP polarimeter to analyze and compare the demosaicing performance of the proposed method and the six previously published methods in both visual comparison and objective evaluation.
The remainder of this paper is organized as follows: Section 2 summarizes the previously published demosaicing methods for DoFP polarimeters and the basic polarization theory. Section 3 describes the selection strategy and preprocessing of the guide image. Section 4 presents the principle and process of the proposed PAIPRI method in detail. Section 5 reports the visual comparison and objective evaluation of the proposed PAIPRI method using both the open-access dataset images collected by a division-of-time (DoT) polarimeter and the indoor and outdoor scene images collected by a real-world DoFP polarimeter. Finally, Section 6 concludes the study.

2. Related Works

2.1. Demosaicing Methods for DoFP Polarimeters

Since Ratliff et al. [28] first discussed the method of reducing the IFOV error for DoFP polarimeters in 2006, more than ten demosaicing methods have been reported for DoFP polarimeters. These methods can be classified as methods of interpolating independently using single-channel methods of interpolating, using other channels as reference images and learning-based methods.

2.1.1. Methods of Interpolating Independently Using Single-Channel

Methods of interpolating independently using single-channel performs analysis and interpolation independently on each channel image, which mainly include polynomial interpolation methods and edge directionality interpolation methods.
The polynomial interpolation method is based on a spatially invariant non-adaptive linear filter. This method estimates the polarization information of un-sampled pixels using the sampling polarization information in the neighborhood by polynomial fitting. This method was first studied due to its low computational complexity. Ref. [30] used the sampling pixels of the same polarization direction in a 3 × 3 neighborhood to estimate the unsampled polarization information of the center pixels through bilinear interpolation. Ref. [33] used the sampling pixels of the same polarization direction in a 5 × 5 neighborhood to estimate the three unsampled polarization information of the center pixels through weighted bilinear, bicubic, and bicubic spline interpolation. Moreover, Ref. [33] designed an approximated bicubic spline method to reduce computational complexity. Its low computational complexity and good reconstruction performance on low-frequency information make the polynomial interpolation method easy to implement on hardware platforms. Ref. [31] implemented the real-time bilinear interpolation of 1600 × 1200 images on FPGA at a speed of 50 frames/s. The adopted window sizes were highly correlated with the PSF function of the imaging system. The polynomial interpolation method performs well with low-frequency information. Increasing the polynomial order can inspire more accurate interpolation results. However, the polynomial interpolation method is usually integrated into other demosaicing methods as a basic method because of its poor performance with high-frequency information.
The edge directionality interpolation method aims to perform polynomial interpolation along the edge instead of across the edge. The foremost task of this method is to accurately detect the edge direction in the incompletely sampled mosaic observed image. Ref. [35] used the horizontal, vertical, and diagonal gradients calculated in the 7 × 7 neighborhood as the criteria to detect the edge direction and performed bicubic convolution interpolation along the edge. Ref. [36] used the block variance, which characterizes the local smoothness, calculated in the 7 × 7 neighborhood, as the criterion to detect the edge direction, and performed bicubic interpolation along the edge. The window adopted sizes are highly correlated with the criterion calculation and the chosen interpolation method. The edge directionality interpolation method performs well in single large-scale edge. However, these criteria are extremely susceptible to IFOV errors, and are less able to discriminate complexly small-scale edges and textures.

2.1.2. Methods of Interpolating Using Other Channels as Reference Images

The method of interpolating using other channels as reference images uses the linear relationship or the similar edge and texture features of the four polarization channels as the reference information source to perform interpolation, which mainly includes correlation-based interpolation methods and residual interpolation methods.
The correlation-based interpolation method uses the linear relationship between the four polarization channels as the reference information source to interpolate un-sampled pixels. Ref. [37] took at least one parameter of the intensity, the DoLP, or the AoP, as prior information to interpolate un-sampled pixels using the linear relationship between the four polarization channels as the reference information source. Ref. [39] designed an edge classifier based on the difference between the two channels and performed Newton polynomial interpolation along the recognized edge direction. Ref. [40] used the weighted fusion of the orthogonal and non-orthogonal polarization channels to interpolate. However, these methods need the’ targets prior information, or need the difference domain between the four polarization channels to be very smooth, and assumed the polarizers to be ideal. These assumptions lead to a simple linear correlation between the four analyzer channels. However, the pixeled polarizer array of DoFP polarimeters has an obvious spatially distributed non-ideality [52]. This non-ideality means that the demosaicing performance of these methods cannot be guaranteed for the real DoFP polarimeter images, and the application scenarios for these methods are extremely limited.
The residual interpolation method can utilize the similar edge and texture features among the four channel images. This method upsamples the input image using the guide filter by referring to the edge and texture features of the guide image, and executes interpolation in the residual domain with less high-frequency information. Ref. [41] selected the low-resolution observed image and high-resolution image of the same channel as the input image and the guide image, respectively, and generated the initial estimate through the minimized residual guide filter (RI). Ref. [43] selected the 0° channel image as the guide image and generated the initial estimate through the minimized Laplacian energy guide filter (MLRI). Ref. [44] selected the edge-aware intensity image generated by an edge detector using the intensity correlation as the guide image, and generated the initial estimate through MLRI. The biggest advantage of residual interpolation is that its parameter tuning is free from training images. This method can still obtain a generally better demosaicing performance, even in the new imaging scene. Moreover, this method has good interpretability for the spatial sampling modes of DoFP polarimeters. Nevertheless, the improper guide-image selection strategy in existing methods fails to fully utilize these advantages, and the performance of these methods can be further improved.

2.1.3. Learning-Based Methods

Learning-based methods construct the sampling models [47] or demosaicing models [48,49,50,51] for DoFP polarimeters by training datasets. This can be achieved by convex optimization [47], dictionary learning [48], and convolutional neural networks (CNN) [49,50,51]. Although learning-based methods generally achieve a higher performance than traditional interpolation methods, they are highly data-dependent [44,53]. A large number of highly representative training images that cover a wide range of scenes are needed to ensure the generalization ability. However, it is very difficult to construct such a training dataset. Moreover, due to the spontaneous emission of infrared polarization devices, the images of DoT and DoFP infrared polarimeters are significantly different. The demosaicing performance of the network trained by DoT infrared images cannot be guaranteed for the real DoFP infrared polarimeters images [39].

2.2. Basic Theory of Polarization Imaging

The Stokes vector S [54] is typically used to describe the polarization characteristics of any light field and can be defined as:
S = [ S 0 S 1 S 2 S 3 ] T ,
where S0 is the total light intensity, S1 is the horizontal or vertical linear polarization component, S2 is the linear polarization component of +45° or −45° polarization directions, and S3 is the left- or right-handed circular polarization component. As the circular polarization component in natural scene radiation is extremely small, S3 is typically considered to be 0. Moreover, DoFP polarimeters only respond to linear Stokes parameters (i.e., S0, S1, and S2). Thus, S3 was omitted from the Stokes vector mentioned in this study.
DoLP and AoP are typically used to investigate the polarization states of the target scene. DoLP represents the proportion of the linearly polarized component to the total intensity of the light source, while AoP represents the angle between the polarization direction of the maximum incident light energy and the x-axis in the reference coordinate system. DoLP and AoP can be calculated using the Stokes vector as follows:
DoLP = S 1 2 + S 2 2 S 0 ,   AoP = 1 2 arctan ( S 2 S 1 ) ,
where DoLP ∊ [0, 1], and DoLP = 1 for linearly polarized light.
The process of polarization imaging and that of reconstructing the incident Stokes vector S using the output grayscale of the four polarization channels (that is, the measurement process) can be represented as [52]:
D N = ( g η ) M S , S = ( g η ) 1 M D N ,
where DN is the output grayscale vector; g is the total gain of the sensor; η is the quantum efficiency of the sensor; M is the coefficient matrix, which characterizes the modulation effect of the pixelated polarizer on the incident Stokes vector, and M is the pseudo-inverse matrix of M, M = ( M T M ) 1 M T .
When the pixelated polarizer array of DoFP polarimeters satisfies the assumption of ideal polarizers (that is, the extinction ratios ε2 of the four pixels in each super-pixel approaches +∞, and polarization direction θ is equal to 0°, 45°, 90°, and 135°, respectively), the ideal normalized coefficient matrix Mideal of a single super-pixel can be represented as follows:
M i d e a l = τ 2 [ 1 1 0 1 0 1 1 1 0 1 0 1 ] ,
where τ is the transmittance coefficient of the pixelated polarizer.

3. Discussion of the Guide Image

3.1. Framework of the Residual Interpolation Methods for DoFP Polarimeters Demosaicing

The residual interpolation method for DoFP polarimeters’ demosaicing assumes that the residual domain contains less high-frequency information, and executes interpolation in the residual domain using the simple polynomial interpolation method to generate a good demosaicing performance [41]. The previously published residual interpolation methods usually consist of three steps (Figure 2):
i
Generate the guide image: We use an up-sampling filter to interpolate the low-resolution observed image I θ 1 L R in a certain polarization direction to generate the guide image I θ 1 _ g u i d e H R . We generally select the same channel as the input image or the 0° channel as the guide image.
ii
Generate the initial estimate: We select four low-resolution observation images I θ k L R (k = 1, 2, 3, 4) as input images to get initial estimates T θ k H R through RI or MLRI guide filters.
iii
Interpolation in residual domain: We calculate a low-resolution residual image R θ k L R by making difference between the initial estimate T θ k H R and the input image I θ k L R . Then, we add the high-resolution residual image R θ k H R generated by interpolating R θ k L R and initial estimate T θ k H R to output the final high-resolution image I θ k H R .
Step (iii) is relatively standardized and fixed in the framework. Thus, we pay major attention to Step (i) and Step (ii). The closer the initial estimate generated in Step (ii) is to the ground truth, the lower the amount of high-frequency information contained in the residual domain, and the better the demosaicing performance of the final high-resolution output image. The initial estimate generated in Step (ii) is the local linear transformation of the guide image generated in Step (i). Therefore, the quality of the guide image directly affects the accuracy of the initial estimation, and further affects the demosaicing performance of the final output image. The quality of the guide image generated in Step (i) is affected by the up-sampling filter and the sampling rate of the observed low-resolution image I θ 1 L R . We simulated the actual polarization imaging process of DoFP polarimeters using an open-access dataset published in SPIE Photonics Europe 2018 [55], which includes 10 real-scene 768 × 1024 16-bits near infrared (NIR) images in 0°, 45°, 90°, and 135° polarization directions. We used these simulating images to analyze the influence of the up-sampling filter and sampling rate on the quality of the guide image and the final output image, and further demonstrated the potential of the proposed method to improve the demosaicing performance of the final output image.

3.2. Influence of the Up-Sampling Filter

The better the performance of the up-sampling filter, the higher the peak signal to noise ratio (PSNR) of the guide image I θ 1 _ g u i d e H R generated in step (i), and the higher the PSNR of the final high-resolution output image. According to the spatial sampling modes of DoFP polarimeters, we down-sampled 10 full-resolution polarization images in 0°, 45°, 90°, and 135° polarization directions of the dataset to generate the observed low-resolution images I 0 ° L R , I 45 ° L R , I 90 ° L R , and I 135 ° L R . Then, we up-sampled I 0 ° L R by operating step (i) to generate the guide image using 10 up-sampling filters: (a) bilinear interpolation; (b) bicubic interpolation; (c) gradient-based interpolation; (d) Newton polynomial interpolation; (e) up-sampling filter GF1 based on guide filter, where we interpolated I 0 ° L R using bilinear interpolation to generate the guide image I 0 ° H R , and up-sampled I 0 ° L R using RI referring to I 0 ° H R [41]; (f) up-sampling filter GF2 based on guide filter, where we interpolated I 0 ° L R using bicubic interpolation to generate the guide image I 0 ° H R , and up-sampled I 0 ° L R using MLRI referring to I 0 ° H R [43]; (g) up-sampling filter GF3 based on guide filter, where we interpolated I 45 ° L R using bilinear interpolation to generate the guide image I 45 ° H R , and up-sampled I 0 ° L R using RI referring to I 45 ° H R ; (h) up-sampling filter GF4 based on guide filter, where we interpolated I 45 ° L R using bicubic interpolation to generate the guide image I 45 ° H R , and up-sampled I 0 ° L R using MLRI referring to I 45 ° H R ; (i) up-sampling filter GF5 based on guide filter, where we interpolated I 0 ° L R using GF3 to generate the guide image I 45 ° H R , and up-sampled I 0 ° L R using RI referring to I 45 ° H R , that is, we iterated over (g); (j) up-sampling filter GF6 based on guide filter, where we interpolated I 45 ° L R using GF4 to generate the guide image I 45 ° H R , and up-sampled I 0 ° L R using MLRI referring to I 45 ° H R , that is, we iterated over (h).
The PSNR of the guide image I 0 ° _ g u i d e H R generated by the above 10 up-sampling filters is illustrated in Figure 3 Compared with other up-sampling filters, the up-sampling filters GF1–GF6, based on guide filters, perform better in the 10 images of the dataset. We compared GF1–GF6 and found that the performance of the guide filter is extremely dependent on the selected polarization direction of the guide image. The PSNR of I 0 ° _ g u i d e H R , generated by GF1 and GF2, is very close to that generated by bilinear interpolation and bicubic interpolation. That is, the PSNR of the output image generated by GF1 and GF2 is basically the same as the PSNR of the guide image, which means that the guide filter does not show a practical effect. This proves that it is inappropriate to choose the same polarization direction for the guide image and the input image of the guide filter, as in the previously published methods [41,42,43]. We have noticed that when images with different polarization directions are selected as the guide image and the input image (for example, interpolating I using I45° as the guide image in GF3–GF6), the up-sampling filter based on the guide filter can improve the PSNR of I 0 ° _ g u i d e H R .
Using the guide image I 0 ° _ g u i d e H R generated by 10 up-sampling filters, we operated steps (ii) and (iii) on I 0 ° L R , I 45 ° L R , I 90 ° L R , and I 135 ° L R to generate high-resolution output images (we selected the 7 × 7 rectangular window as the filter window) and further calculated the DoLP images (Figure 4). It can be seen that the better the performance of the up-sampling filter, the higher the PSNR of the guide image, and the higher the PSNR of the final high-resolution output image. We compared the results of GF3 and GF5, GF4 and GF6, and found that if we continuously update the guide image through an iterative process and use the output image in the previous iteration as the guide image in the next iteration, the PSNR of the final output image can be increased. Therefore, the proposed method can significantly improve the quality of the output image by selecting an appropriate polarization direction for the guide image and multi-iteration.

3.3. Influence of the Sampling Rate of the Guide Image

The higher the sampling rate of the low-resolution observed image I 0 ° L R , the higher the PSNR of the guide image I θ 1 _ g u i d e H R generated in Step (i), and the higher the PSNR of the final high-resolution output image. We performed 50% and 25% down-sampling of the full-resolution polarization image I 0 ° F R in the dataset to generate the low-resolution polarization images I 0 ° _ d 2 L R and I 0 ° _ d 4 L R (Figure 5). Then, we up-sampled I 0 ° _ d 2 L R and I 0 ° _ d 4 L R by operating Step (i) to generate I 0 ° _ d 2 H R and I 0 ° _ d 4 H R using bilinear interpolation. Using I 0 ° F R , I 0 ° _ d 2 H R and I 0 ° _ d 4 H R as the guide image, we operated Steps (ii) and (iii) on I 0 ° L R , I 45 ° L R , I 90 ° L R , and I 135 ° L R to generate high resolution output images (we selected the 7 × 7 rectangular window as the filter window) and further calculated the DoLP images (Figure 6).
It can be seen from Figure 6 that the higher the sampling rate of the low-resolution observed image I θ k L R used to generate the guide image, the higher the PSNR of the final output high-resolution image. Interpolation for Bayer color camera also follows this rule. When interpolating the R and B channels, we generally choose the interpolation result of the G channel with a higher sampling rate as the guide image. For the DoFP polarimeters, although the 0°, 45°, 90°, and 135° channels have the same sampling rate, the sampling rate of the guide image in the filter window can be increased by choosing the appropriate size and direction for the filtering window in the guide filter. When using the high-resolution I 0 ° _ g u i d e H R as the guide image to interpolate the missing 45° polarization information at the 0° channel (Figure 7), if we only operate a horizontal guide filter on odd rows and choose a 1 × h rectangular filter window, the sampling rate of the guide image can be increased to 50%. For the same reason, when interpolating the missing 45° polarization information at the 90° channel, we only operate a vertical guide filter on even columns; when interpolating the missing 45° polarization information at the 135° channel, we operate the guide filter along the two diagonal directions, and then fuse the interpolation results in these two directions. According to the spatial layout of the pixeled polarizer array, choosing different channels as the guide image for the pixels in different spatial positions, and cooperating with the different sizes and directions of the filter window, can increase the sampling rate of the guide image in the filter window, and further increase the PSNR of the final output image.
The proposed PAIPRI in this study chooses different channels as the guide image for the pixels in different spatial positions, according to the spatial layout of the pixeled polarizer array. Then, it designs the filter windows with different sizes and directions, and updates the guide image through an iterative process based on the guide filter. The analysis and discussion in this section indicate that the proposed PAIPRI in this study has great potential to improve the demosaicing performance of the final output image compared with the previously published demosaicing method for DoFP polarimeters.

4. The Proposed PAIPRI

4.1. Overall Pipeline

This study proposed a residual interpolation method, with an integrated pixel-by-pixel adaptive iterative process, for DoFP polarimeters (PAIPRI). Compared with the previously published PRI and MLPRI, the proposed PAIPRI innovatively designed a new guide-image selection strategy, and fused the local iterative optimal results of RI and MLRI. We improved the demosaicing performance of the final output image by increasing the sampling rate and PSNR of the guide image.
The overall pipeline of the proposed PAIPRI is illustrated in Figure 8. We used the up-sampling process of the low-resolution observed image I 0 ° L R as an example. The up-sampling processes of low-resolution observed images in other polarization directions followed the same principle as that of I 0 ° L R . The proposed PAIPRI consists of two steps:
I
Pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal, vertical, and two diagonal directions: We chose different channels as the guide images for pixels in different spatial positions according to the spatial layout of the pixeled polarizer array, and designed the filter windows with different sizes and directions. When using I 0 ° L R as the input image, we operated iterative RI and MLRI in horizontal, vertical, and two diagonal directions, referring to the guide images I45°, I135°, and I90°, respectively. In each iterative process, a local criterion was calculated for each reconstructed pixel to adaptively determine whether to update the interpolation result in this iteration. Until all pixels in FPA completed their update or the iterative number reached the maximum iterative number, eight sets of interpolation images, with RI and MLRI in the horizontal, vertical and two diagonal directions, could be obtained.
II
According to the spatial layout of the reconstructed pixels, we performed an adaptively weighted average fusion on the eight sets of interpolation images with RI and MLRI in the horizontal, vertical and two diagonal directions to generate the final output up-sampling image, I 0 ° H R .

4.2. Pixel-by-Pixel Adaptive Iterative Processes Based on Residual Interpolation

As an example, we used the pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal direction to interpolate the missing 0° polarization information at the 45° channel (Figure 9). We performed horizontal RI and MLRI, referring to the guide images I and I45°. The interpolation result I and the guide image I45° were adaptively updated pixel-by-pixel using two interrelated iterative processes in the primary branch. Then, the high-resolution horizontal interpolation images R I 0 ° H R H and M L R I 0 ° H R H were generated. In this step, the guide image, size, and direction of the filter window were selected according to the spatial layout of the pixeled polarizer array of DoFP polarimeters. When interpolating the missing 0° polarization information at the 45° channel, we selected I45° as the guide image to operate horizontal RI and MLRI. The sampling rate of the guide image in the filter window was increased to 50% (however, the directionless square window selected in the previously published residual interpolation methods made the sampling rate of the guide image only 25%). The increased sampling rate, cooperating with the iterative process, contributed to the increase in the PSNR of the finial output image. Interpolation of the missing 0° polarization information at the 90° and 135° channels followed the same principle as that at the 0° channel, with the guide image replaced with I90° and I135° and the filtering direction adjusted to the vertical and two diagonal directions.
The pixel-by-pixel adaptive iterative processes based on residual interpolation consists of four steps:
(i)
Calculate the initial value I 0 ° H R ( 0 ) and I 45 ° H R ( 0 ) of the iteration
We performed a horizontal linear interpolation on the observed low-resolution image I 0 ° L R and I 45 ° L R to calculate the initial value I 0 ° H R ( 0 ) and I 45 ° H R ( 0 ) of the iteration:
I 0 ° H R ( 0 ) ( 2 i 1 , j ) = { ( I 0 ° L R ( 2 i 1 , j 1 ) + I 0 ° L R ( 2 i 1 , j + 1 ) ) / 2   if   j = 2 d   I 0 ° L R ( 2 i 1 , j )                                                                                     if   j = 2 d 1 , I 45 ° H R ( 0 ) ( 2 i 1 , j ) = { ( I 45 ° L R ( 2 i 1 , j 1 ) + I 45 ° L R ( 2 i 1 , j + 1 ) ) / 2   if   j = 2 d 1 I 45 ° L R ( 2 i 1 , j )                                                                                       if   j = 2 d   ,
where dN+; (i, j) are the pixel indices; m ∈ [1, 2M], n ∈ [1, 2N], and 2M × 2N is the size of the sensor.
In the first iteration, the up-sampling filter was selected as the simple bilinear interpolation. Although the initial iterative value obtained by this simple up-sampling filter was not the best, it greatly simplifies the calculation steps and saves time. The impact of this imperfection in the initial iterative value on the final output images was almost negligible.
(ii)
Calculate the initial estimate T 0 ° k and T 45 ° k
Except for the first iteration, the up-sampling filter in each iteration was the guide filter that was proved to be optimal in Section 3 to increase the PSNR of the final output image. The initial estimate was calculated by the horizontal RI and MLRI through two interrelated iterative processes. In the primary branch, the input and the guide image of the k-th iteration, respectively, selected the output of the previous iteration result I 0 ° H R ( k 1 ) and I 45 ° H R ( k 1 ) of the primary branch and the auxiliary branch. In the auxiliary branch, the input and the guide image of the k-th iteration, respectively, selected the output of the previous iteration result I 45 ° H R ( k 1 ) and I 0 ° H R ( k 1 ) of the auxiliary branch and the primary branch. The initial estimate T 0 ° k and T 45 ° k by the horizontal RI and MLRI can be expressed as the local linear transformation of the guide image I 45 ° H R ( k 1 ) and I 0 ° H R ( k 1 ) :
T 0 ° k ( i , j ) = a 0 ° k ( m , n ) I 45 ° H R ( k 1 ) ( i , j ) + b 0 ° k ( m , n ) , ( i , j ) ω m n k , T 45 ° k ( i , j ) = a 45 ° k ( m , n ) I 0 ° H R ( k 1 ) ( i , j ) + b 45 ° k ( m , n ) , ( i , j ) ω m n k ,
where ω m n k represents the filter window selected in the k-th iteration; (m, n) is the index of the center pixel in the filter window ω m n k ; Hk × Vk is the size of filter window. In the iterative process, we empirically chose a gradually increasing window size [24]; the window size of the k-th iteration was
H k = { 1 if   k = 1   in   RI 5 if   k = 1   in   MLRI H k 1 + 2 if   k > 1   in   RI   and   MLRI , V k = { 5 if   k = 1   in   RI   and   MLRI V k 1 + 2 if   k > 1   in   RI   and   MLRI ,
( a 0 ° k ( m , n ) , b 0 ° k ( m , n ) ) and ( a 45 ° k ( m , n ) , b 45 ° k ( m , n ) ) are linear coefficients, which were assumed to be constant in the filter window with the center pixel (m, n).
The main difference between RI and MLRI is the different cost functions for solving linear coefficients. When solving linear coefficients in RI, the total difference between the initial estimate T 0 ° k and T 45 ° k and the input image I 0 ° H R ( k 1 ) and I 45 ° H R ( k 1 ) in the filter window must be minimized. When solving linear coefficients in MLRI, the total Laplacian energy of the difference between the initial estimate T 0 ° k and T 45 ° k and the input image I 0 ° H R ( k 1 ) and I 45 ° H R ( k 1 ) in the filter window must be minimized to ensure similar image smoothness between the guide image and the initial estimate. RI and MLRI calculate linear coefficients by minimizing the following cost functions in ω m n k , respectively:
E ( a 0 ° k ( m , n ) , b 0 ° k ( m , n ) ) =   { i , j ω m n k ( T 0 ° k ( i , j ) I 0 ° H R ( k 1 ) ( i , j ) ) 2 in   RI i , j ω m n k ( Δ ( T 0 ° k ( i , j ) I 0 ° H R ( k 1 ) ( i , j ) ) ) 2 in   MLRI ,
E ( a 45 ° k ( m , n ) , b 45 ° k ( m , n ) ) = { i , j ω m n k ( T 45 ° k ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) ) 2 in   RI i , j ω m n k ( Δ ( T 45 ° k ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) ) ) 2 in   MLRI ,
where Δ is the operation calculating Laplacian energy, Δ I = [ 0 1 0 1 4 1 0 1 0 ] I , and ⊗ is a convolution operation.
We solved Equations (7) and (8) through linear regression to calculate a set of solutions to linear coefficients:
a 0 ° k ( m , n ) = { 1 C m n k i , j ω m n k ( I 0 ° H R ( k 1 ) ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) ) μ 0 ° k 1 ( m , n ) μ 45 ° k 1 ( m , n ) ( σ 45 ° k 1 ( m , n ) ) 2 in   RI i , j ω m n k ( Δ I 0 ° H R ( k 1 ) ( i , j ) Δ I 45 ° H R ( k 1 ) ( i , j ) ) i , j ω m n ( Δ I 0 ° H R ( k 1 ) ( i , j ) ) 2 in   MLRI , b 0 ° k ( m , n ) = μ 0 ° k 1 ( m , n ) a 0 ° k ( m , n ) μ 45 ° k 1 ( m , n ) in   RI   and   MLRI ,
a 45 ° k ( m , n ) = { 1 C m n k i , j ω m n k ( I 0 ° H R ( k 1 ) ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) ) μ 0 ° k 1 ( m , n ) μ 45 ° k 1 ( m , n ) ( σ 0 ° k 1 ( m , n ) ) 2 in   RI i , j ω m n k ( Δ I 0 ° H R ( k 1 ) ( i , j ) Δ I 45 ° H R ( k 1 ) ( i , j ) ) i , j ω m n ( Δ I 45 ° H R ( k 1 ) ( i , j ) ) 2 in   MLRI , b 45 ° k ( m , n ) = μ 45 ° k 1 ( m , n ) a 45 ° k ( m , n ) μ 0 ° k 1 ( m , n ) in   RI   and   MLRI ,
where C m n k is the number of whole pixels in ω m n k ; μ 0 ° k 1 ( m , n ) and σ 0 ° k 1 ( m , n ) , μ 45 ° k 1 ( m , n ) and σ 45 ° k 1 ( m , n ) are the mean and standard deviation of I 0 ° H R ( k 1 ) and I 45 ° H R ( k 1 ) in ω m n k .
We can calculate a pair of linear coefficients in ω m n k using Equations (9) and (10). When the filter window traverses all pixels on the FPA, the target pixel is contained in different windows, and corresponds to different linear coefficients. Therefore, we performed a weighted average fusion on these linear coefficients to represent the composite effect of all filter windows containing the target pixel. Then, we calculated the unique pair of linear coefficients corresponding to the target pixel located at (i, j):
a ¯ 0 ° k ( i , j ) = m , n ω i j k W 0 ° ( m , n ) a 0 ° k ( m , n ) m , n ω i j k W 0 ° ( m , n ) ,   b ¯ 0 ° k ( i , j ) = m , n ω i j k W 0 ° ( m , n ) b 0 ° k ( m , n ) m , n ω i j k W 0 ° ( m , n ) ,
a ¯ 45 ° k ( i , j ) = m , n ω i j k W 45 ° ( m , n ) a 45 ° k ( m , n ) m , n ω i j k W 45 ° ( m , n ) ,   b ¯ 45 ° k ( i , j ) = m , n ω i j k W 45 ° ( m , n ) b 45 ° k ( m , n ) m , n ω i j k W 45 ° ( m , n ) ,
where W (m, n) and W45° (m, n) were the corresponding weights of the target pixel in different filter windows.
When calculating linear coefficients ( a ¯ 0 ° k ( i , j ) , b ¯ 0 ° k ( i , j ) ) and ( a ¯ 45 ° k ( i , j ) , b ¯ 45 ° k ( i , j ) ) using Equations (9)–(12), the output initial estimate of the guide filter can be expressed as:
T 0 ° k ( i , j ) = a ¯ 0 ° k ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) + b ¯ 0 ° k ( i , j ) , T 45 ° k ( i , j ) = a ¯ 45 ° k ( i , j ) I 0 ° H R ( k 1 ) ( i , j ) + b ¯ 45 ° k ( i , j ) ,
(iii)
Calculate the residual R 0 ° H R ( k ) and R 45 ° H R ( k )
The residual represents the difference between the output initial estimates T 0 ° k and T 45 ° k of the guide filter and the low-resolution observed image I 0 ° L R and I 45 ° L R , which can characterize the accuracy of the initial estimates. The low-resolution residual images R 0 ° L R ( k ) and R 45 ° L R ( k ) can be calculated as:
R 0 ° L R ( k ) ( i , j ) = { | T 0 ° k ( i , j ) I 0 ° L R ( i , j ) | if   i = 2 d - 1 ,   j = 2 e - 1 0 o t h e r w i s e ,
R 45 ° L R ( k ) ( i , j ) = { | T 45 ° k ( i , j ) I 45 ° L R ( i , j ) | if   i = 2 d - 1 ,   j = 2 e 0 o t h e r w i s e ,
where d, eN+. We calculated high-resolution residual images R 0 ° H R ( k ) and R 45 ° H R ( k ) by operating a horizontal linear interpolation on low-resolution residual images R 0 ° L R ( k ) and R 45 ° L R ( k ) .
(iv)
Pixel-by-pixel adaptively updated iterative results
We defined the following pixel-by-pixel criterion to determine whether to update the interpolation result in the current iteration. For the interpolation result located at (i, j) of the primary branch and the auxiliary branch in k-th iteration, criteria c 0 ° H ( k ) ( i , j ) and c 45 ° H ( k ) ( i , j ) were determined:
c 0 ° H ( k ) ( i , j ) = a 0 ° H ( k ) ( i , j ) + t 0 ° H ( k ) ( i , j ) , a 0 ° H ( k ) ( i , j ) = I F G a u s s i a n ( | T 0 ° k ( i , j ) I 0 ° H R ( k 1 ) ( i , j ) | ) , t 0 ° H ( k ) ( i , j ) = 1 C m n k i , j ω m n k I F G a u s s i a n ( R 0 ° L R ( k ) ( i , j ) ) ,
c 45 ° H ( k ) ( i , j ) = a 45 ° H ( k ) ( i , j ) + t 45 ° H ( k ) ( i , j ) , a 45 ° H ( k ) ( i , j ) = I F G a u s s i a n ( | T 45 ° k ( i , j ) I 45 ° H R ( k 1 ) ( i , j ) | ) , t 45 ° H ( k ) ( i , j ) = 1 C m n k i , j ω m n k I F G a u s s i a n ( R 45 ° L R ( k ) ( i , j ) ) ,
where a 0 ° H ( k ) ( i , j ) and a 45 ° H ( k ) ( i , j ) describe the convergence of the initial estimates T 0 ° k ( i , j ) and T 45 ° k ( i , j ) , obtained from step (ii) in the k-th iteration, compared to that obtained in the previous iteration, respectively; t 0 ° H ( k ) ( i , j ) and t 45 ° H ( k ) ( i , j ) describe the closeness of the initial estimates T 0 ° k ( i , j ) and T 45 ° k ( i , j ) , obtained from step (ii) in the k-th iteration, to the observed low-resolution image I 0 ° L R and I 45 ° L R ; IFGaussian is the spatial Gaussian filter. We empirically selected a 5 × 5 Gaussian kernel with the standard variation σ = 1. The guide filter was a local linear model. Therefore, we used a spatial Gaussian filter to take the influence of neighboring pixels into consideration when calculating the criterion of the target pixel (i, j), which made the proposed criterion more reliable [24].
Using the criteria calculated from Equation (16), we adaptively updated the iterative results, pixel by pixel, according to the following decision conditions:
I 0 ° H R ( k ) ( i , j ) = { T 0 ° k ( i , j ) + R 0 ° H R ( k ) ( i , j ) i f   c 0 ° H ( k ) ( i , j ) < min 1 w k 1 ( c 0 ° H ( w ) ( i , j ) ) I 0 ° H R ( k 1 ) ( i , j ) o t h e r w i s e ,
I 45 ° H R ( k ) ( i , j ) = { T 45 ° k ( i , j ) + R 45 ° H R ( k ) ( i , j ) i f   c 45 ° H ( k ) ( i , j ) < min 1 w k 1 ( c 45 ° H ( w ) ( i , j ) ) I 45 ° H R ( k 1 ) ( i , j ) o t h e r w i s e ,
where min is the operation calculating minimum value. When the criterion of the k-th iteration is smaller than that of the previous k − 1 iterations, the interpolation result located at (i, j) is updated as the sum of the initial estimate T 0 ° k ( i , j ) and T 45 ° k ( i , j ) and the high-resolution residuals R 0 ° H R ( k ) and R 45 ° H R ( k ) . Otherwise, the interpolation result located at (i, j) is not updated in the k-th iteration. When all pixels in FPA complete update, or the iterative number reaches the maximum iterative number K, the iterative process stops, and the output I 0 ° , R I H and I 0 ° , M L R I H of step (Ⅰ) are generated.

4.3. Fusion on the Iterative Results

According to the spatial layout of the reconstructed pixels, we performed an adaptively weighted average fusion on the eight sets of interpolation images with RI and MLRI in the horizontal, vertical and two diagonal directions obtained in Step (Ⅰ) to generate the finial output up-sampling image I 0 ° H R :
I 0 ° H R ( i , j ) = { I 0 ° L R ( i , j ) i f   i = 2 d 1 ,   j = 2 e 1 W 0 ° , R I H ( i , j ) I 0 ° , R I H ( i , j ) + W 0 ° , M L R I H ( i , j ) I 0 ° , M L R I H ( i , j ) i f   i = 2 d 1 ,   j = 2 e W 0 ° , R I V ( i , j ) I 0 ° , R I V ( i , j ) + W 0 ° , M L R I V ( i , j ) I 0 ° , M L R I V ( i , j ) i f   i = 2 d ,   j = 2 e 1 w = 1 2 [ W 0 ° , M L R I D w ( i , j ) I 0 ° , R I D w ( i , j ) + W 0 ° , M L R I D w ( i , j ) I 0 ° , M L R I D w ( i , j ) ] i f   i = 2 d ,   j = 2 e ,
where W is the reciprocal of the minimum value of criteria in 1~K iterations. The smaller the criteria, the greater the weight.
W 0 ° , R I H ( i , j ) = 1 / min 1 w K ( c 0 ° , R I H ( w ) ( i , j ) ) ,   W 0 ° , M L R I H ( i , j ) = 1 / min 1 w K ( c 0 ° , M L R I H ( w ) ( i , j ) ) , W 0 ° , R I V ( i , j ) = 1 / min 1 w K ( c 0 ° , R I V ( w ) ( i , j ) ) , W 0 ° , M L R I V ( i , j ) = 1 / min 1 w K ( c 0 ° , M L R I V ( w ) ( i , j ) ) , W 0 ° , R I D 1 ( i , j ) = 1 / min 1 w K ( c 0 ° , R I D 1 ( w ) ( i , j ) ) ,   W 0 ° , M L R I D 1 ( i , j ) = 1 / min 1 w K ( c 0 ° , M L R I D 1 ( w ) ( i , j ) ) , W 0 ° , R I D 2 ( i , j ) = 1 / min 1 w K ( c 0 ° , R I D 2 ( w ) ( i , j ) ) ,   W 0 ° , M L R I D 2 ( i , j ) = 1 / min 1 w K ( c 0 ° , M L R I D 2 ( w ) ( i , j ) ) ,
We take the up-sampling process of the observed low-resolution image I 0 ° L R as an example to illustrate the overall pipeline of the proposed PAIPRI. The up-sampling processes of low-resolution observed images I 45 ° L R , I 90 ° L R , and I 135 ° L R follow the same principle as that of I 0 ° L R . When up-sampling I 45 ° L R , the guide images in horizontal, vertical and diagonal directions are I, I90°, and I135°, respectively. For the same reason, when up-sampling I 90 ° L R , the guide images in horizontal, vertical and diagonal directions are I135°, I45°, and I, respectively; when up-sampling I 135 ° L R , the guide images in horizontal, vertical and diagonal directions are I90°, I, and I45°, respectively. After completing the up-sampling of the four observed low-resolution images, four high-resolution output images I 0 ° H R , I 45 ° H R , I 90 ° H R , and I 135 ° H R were generated. Then, we substituted these high-resolution output images into Equations (1)–(4), and reconstructed the high-resolution Stokes vector, DoLP, and AoP images. The whole PAIPRI procedure is presented in (Algorithm 1).
Alogrithem 1: PAIPRI
Input: Given the low-resolution observed images of the four polarization direction I 0 ° L R , I 45 ° L R , I 90 ° L R , and I 135 ° L R , the initial value of the window size, and the maximum number of iterations kmax.
Output: Four high-resolution output images I 0 ° H R , I 45 ° H R , I 90 ° H R , and I 135 ° H R .
For k = 1:kmax:
 (i) Calculate the initial iterative value using Equation (5).
 (ii) Calculate the initial estimate using RI and MLRI in horizontal, vertical, and two diagonal directions for each polarization direction. Solve the linear coefficients using Equations (9)–(12). Then, substitute the linear coefficients and the previous iteration result into Equation (13) to generate the initial estimate in current iteration.
 (iii) Calculate the residual images in horizontal, vertical, and two diagonal directions for each polarization direction. Substitute the input low-resolution observed image and the initial estimate generated by Step (ii) into Equation (14) to generate the residual images.
 (iv) Pixel-by-pixel adaptively update iterative results.
 If criteria in k-th iteration < criteria in the previous iteration:
Update iterative results in this pixel using Equations (18) and (19).
 end
end
 (v) Generate the finial output images by adaptively weighting the eight sets of interpolation images with RI and MLRI in the horizontal, vertical and two diagonal directions using Equation (20) after k reaches kmax or all the pixels complete updating.

5. Experimental Verification and Discussion

This section aims to prove that the proposed PAIPRI exhibits a better demosaicing performance compared with the existing methods for DoFP polarimeters. We compared the demosaicing performance of the proposed PAIPRI with that of the seven existing methods for DoFP polarimeters, including the bilinear interpolation (Bilinear), the bicubic spline interpolation (BS), the gradient-based interpolation (Gradient [35]), the Newton polynomial interpolation (NP [39]), the residual interpolation with minimized residual (PRI [41]), the residual interpolation with minimized Laplacian energy (MLPRI [43]) and edge-aware residual interpolation (EARI [44]). EARI and MLPRI were proven to be the state-of-the-art, non-learning-based polarization demosaicing methods [43,44]. We conducted experiments on both the open-access dataset images [44,55] collected by a division-of-time polarimeter and the indoor and outdoor scene images collected by a real-world DoFP polarimeter demonstrate to analyze and compare the demosaicing performance of the proposed PAIPRI and the seven existing methods in both visual comparison and objective evaluation. It should be noted that the learning-based methods are highly data-dependent. The CNN-based methods in [48,49] did not disclose their training datasets or pre-trained network weights, while the low sample number of open-access datasets [44,55] makes it difficult to produce a satisfactory training result. To ensure fairness, we did not compare these learning-based methods in this section, as in [39,47,48,49,50,51].

5.1. Dataset

We used two open-access datasets, published in SPIE Photonics Europe 2018 [55] and the Morimatsu dataset [44], as the image sources in experiments. According to the spatial sampling modes of DoFP polarimeters, we generated the observed low-resolution images by down-sampling the high-resolution polarization images in 0°, 45°, 90°, and 135° polarization directions of the dataset, which simulated the actual imaging process of DoFP polarimeters. Subsequently, we performed interpolation on these simulated low-resolution observed images using the proposed PAIPRI and six previously published methods. Then, the reconstructed high-resolution polarization images I, I45°, I90°, and I135° were obtained and substituted into Equations (1)–(4) to reconstruct the high-resolution Stokes vector, DoLP, and AoP images. It should be noted that, when implementing PRI [41] in this section, we used the bilinear interpolation results of the observed low-resolution image instead of the ground-truth image as the guide image [43]. The source codes of EARI [44] were downloaded from the author’s websites, and the error in calculating the Stokes vector was corrected (the pseudo-inverse instead of transpose of M should be used, as shown in Equation (3)).
The polarization demosaicing method aims to obtain an accurate estimate of the unsampled polarization information. It is not sufficient to evaluate a method solely on whether it outputs smooth, visually good results. Therefore, we also carefully analyzed the objective evaluation results. We calculated and compared the PSNR of the reconstructed results for the dataset images. The PSNRs of the reconstructed I, S2, and DoLP images for dataset [55] are illustrated in Table 1, Table 2 and Table 3 (similarly, there were reconstructed results of I45°, I90°, I135°, S0, S1, and AoP, which are not exhibited in order to save space). The highest PSNR of each row in Table 1−3 is shown in bold. The methods using the neighborhood information could not reconstruct the correct information at the boundary of the filled image. Therefore, pixels within 10 pixels from the boundary were excluded in the calculation of PSNR to eliminate the interference of the incorrect information at the boundary in the methods’ performance evaluation. Similar results for the dataset [44] are illustrated in the Supplemental Document (Supplementary Materials). Considering Table 1, Table 2 and Table 3, it can be observed that, for the 10 scene images in the tested dataset, the proposed PAIPRI performs better than the other seven comparison methods in the objective evaluation based on the index PSNR. Compared with the optimal results in the other seven comparison methods, the average PSNR of I, S2, and DoLP images reconstructed by PAIPRI are increased by 1.33 dB, 1.31 dB, and 0.78 dB, respectively.
We selected three types of representative local regions in the dataset to complete a visual comparison of the demosaicing performance for these eight methods:
(1)
Single arc-shaped edge: Due to the difference in material and surface roughness between the target and the background, the boundary in the selected local region 1 appears as continuous and sharp arc-shaped edges in both the intensity images and the polarization images. This type of edge and its neighborhoods in the reconstructed results are easily affected by the IFOV error, and further exhibit a sawtooth effect.
(2)
Multi-directional assorted edges: The selected local region 2 contains at least two of the horizontal, vertical, multi-oblique, or arc-shaped edges. This type of edge, and their neighborhoods in the reconstructed results, are easily affected by the IFOV error, and further exhibit sawtooth effect or edge artifacts.
(3)
Abundant texture features: The selected local region 3 contains a periodic hole structure. This periodic hole structure appears as a distinct texture feature in both the intensity images and the polarization images. This texture feature is easily affected by the IFOV error in the reconstructed results, and further exhibits additional error textures.
The selected three types of local region can basically cover the image edge and texture features that are susceptible to IFOV errors, and can be used to reasonably evaluate and compare the demosaicing performance for different methods.
The reconstructed I, S2, and DoLP images of the selected three types of representative local regions in the dataset are exhibited in Figure 10, Figure 11 and Figure 12. For the three selected types of representative local regions, the proposed PAIPRI can clearly and accurately reconstruct the edge and texture features, and performed better than the other seven methods in visual comparison, which is consistent with the conclusions obtained from the objective evaluation. For polarization imaging, we apparently paid more attention to the performance of the reconstructed polarization information such as the Stokes vector, the DoLP and the AoP images. However, these images were calculated by the four interpolated polarization channel images, so they are DoLP images, and also support this viewpoint. The proposed PAIPRI still exhibits good demosaicing performance in S2 and DoLP images (Figure 10a,b, Figure 11a,b and Figure 12a,b). Although the reconstructed results generated by PAIPRI still retains a small amount of mosaic effect on the edges, it did not produce obvious sawtooth effect or edge artifacts at the edges and texture features, nor does it show blurred edges due to excessive smoothing. The reconstructed results generated by PAIPRI are also very visually close to the ground-truth images.
The reconstructed I, S2, and DoLP images of the selected three types of representative local regions in the dataset are exhibited in Figure 10, Figure 11 and Figure 12. For the selected three types of representative local regions, the proposed PAIPRI can clearly and accurately reconstruct the edge and texture features, and performed better than the other seven methods in terms of visual comparison, which is consistent with the conclusions obtained from the objective evaluation. For polarization imaging, we apparently paid more attention to the performance of the reconstructed polarization information such as the Stokes vector, the DoLP and the AoP images. However, these images were calculated by the four interpolated polarization channel images, so they are DoLP images, and also support this viewpoint. The proposed PAIPRI still exhibited a good demosaicing performance in S2 and DoLP images (Figure 10a,b, Figure 11a,b and Figure 12a,b). Although the reconstructed results generated by PAIPRI still retain a small amount of mosaic effect on the edges, it did not produce an obvious sawtooth effect or edge artifacts at the edges and texture features, nor did it show blurred edges due to excessive smoothing. The reconstructed results generated by PAIPRI were also very visually close to the ground-truth images.
The computation times of the seven comparison methods and the PAIPRI with different number of iterations are illustrated in Table 4. Figure 13 shows the relationship between the PSNR of I and the number of iterations. The results shows that, with just five iterations, it is possible to obtain a PSNR that is close to the best one obtained. With more than 15 iterations, the increase in PSNR becomes almost negligible. Compared to other methods, our algorithm has more complex but highly parallelized processing in a single iteration. Thus, the better image quality in a single iteration reduces the number of iterations required for the algorithm to converge. Considering the significant improvement in demosaicing performance, the slight increase in the time needed for the proposed PAIPRI could be a good trade-off.

5.2. Images Collected by a Real-World DoFP Polarimeter

We evaluated the demosaicing performance of the proposed PAIPRI through visual comparison using images collected by a real-world DoFP polarimeter. We used the PHX050S DoFP polarimeter of LUCID Vision Labs to collect images, with a sensor size of 2048 × 2448. The adopted image format was 8 bits. The visual comparison of the seven comparison methods using images collected by a real-world DoFP polarimeter are exhibited in Figure 14 and Figure 15. It can be observed that the proposed PAIPRI can clearly and accurately reconstruct the edge and texture features, and performed better than the other six methods in a visual comparison, which is consistent with the conclusions obtained from the objective evaluation and visual comparison of the dataset images.

6. Conclusions

Looking at the problems in the selection and preprocessing of the guide image in the previously published residual interpolation methods for DoFP polarimeters, this study proposed a residual interpolation method, with an integrated pixel-by-pixel adaptive iterative process, for DoFP polarimeters. By thoroughly considering the spatial layout of the pixeled polarizer array, we proposed a new guide-image selection strategy, that is, choosing different channels for the pixels in different spatial positions as the guide image, and cooperating with the different sizes and directions of the filter window, which increased the sampling rate of the adopted guide image in the filter window to 50%. Furthermore, the pixel-by-pixel method adaptively updated the guide image and the interpolation result through two interrelated iterative processes and performed an adaptively weighted average fusion on the iterative results of RI and MLRI, which improved the demosaicing performance of the finial output images. Comparison experiments using both the open-access dataset images collected by a DoT polarimeter and the indoor and outdoor scene images collected by a real-world DoFP polarimeter demonstrated that, in a visual comparison, the proposed PAIPRI can reconstruct clear edges and texture features; in an objective evaluation, the average PSNR of I, S2, and DoLP images reconstructed by PAIPRI are increased by 1.33 dB, 1.31 dB, and 0.78 dB compared with the optimal results in the other seven comparison methods. In brief, the proposed PAIPRI is superior to the existing state-of-the-art methods in terms of both visual comparison and objective evaluation. The results of this study prove that considering the spatial layout of the pixeled polarizer array on the physical level is vital for improving the performances of interpolation methods for DoFP polarimeters.
Redundancy exists between the four polarization channels, which are designed to reduce the noise sensitivity of DoFP polarimeters. Under the condition that the coefficient matrix is known a priori, these four polarization channels satisfy a specific linear relationship. This study has many interesting perspectives, such as using this linear relationship as a reference information source for interpolation and adding more constraints in the process of using the guide filter to solve the initial estimate and the process of residual interpolation. Through these strategies, the interpolation results of the four polarization channels will also meet the linear relationship; additionally, the demosaicing performance is expected to improve.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/s22041529/s1, Table S1: PSNR of I Reconstructed by Different Methods on Dataset.

Author Contributions

Conceptualization, J.Y. and W.J.; methodology, software, formal analysis and writing—original draft preparation validation, J.Y.; writing—review and editing, J.Y., W.J. and F.X.; data curation, M.W. and S.Q.; resources, project administration and funding acquisition, W.J. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by National Natural Science Foundation of China (NSFC), grant number 62171024.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Gurton, K.; Felton, M.; Mack, R.; LeMaster, D.; Farlow, C.; Kudenov, M.; Pezzaniti, L. MidIR and LWIR polarimetric sensor comparison study. Int. Soc. Opt. Photonics 2010, 7672, 767205. [Google Scholar]
  2. Hu, Y.; Li, Y.; Pan, Z. A Dual-Polarimetric SAR Ship Detection Dataset and a Memory-Augmented Autoencoder-Based Detection Method. Sensors 2021, 21, 8478. [Google Scholar] [CrossRef]
  3. Zhou, Y.; Li, Z.; Zhou, J.; Li, N.; Zhou, X.; Chen, P.; Zheng, Y.; Chen, X.; Lu, W. High extinction ratio super pixel for long wavelength infrared polarization imaging detection based on plasmonic microcavity quantum well infrared photodetectors. Sci. Rep. 2018, 8, 15070. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  4. Reda, M.; Zhao, Y.; Chan, J. Polarization Guided Auto-Regressive Model for Depth Recovery. IEEE Photonics J. 2017, 9, 6803016. [Google Scholar] [CrossRef]
  5. Wang, S.; Liu, B.; Chen, Z.; Li, H.; Jiang, S. The Segmentation Method of Target Point Cloud for Polarization-Modulated 3D Imaging. Sensors 2020, 20, 179. [Google Scholar] [CrossRef] [Green Version]
  6. Miyazaki, D.; Shigetomi, T.; Baba, M.; Furukawa, R. Surface normal estimation of black specular objects from multiview polarization images. Opt. Eng. 2016, 56, 041303. [Google Scholar] [CrossRef] [Green Version]
  7. Xue, F.; Jin, W.; Qiu, S.; Yang, J. Polarized observations for advanced atmosphere-ocean algorithms using airborne multi-spectral hyper-angular polarimetric imager. ISPRS J. Photogramm. Remote Sens. 2021, 178, 136–154. [Google Scholar] [CrossRef]
  8. Ge, B.; Mei, X.; Li, Z.; Hou, W.; Xie, Y.; Zhang, Y.; Xu, H.; Li, K.; Wei, Y. Airborne optical polarization imaging for observation of submarine Kelvin wakes on the sea surface: Imaging chain and simulation. Remote Sens. Environ. 2020, 247, 111894. [Google Scholar] [CrossRef]
  9. Ojha, N.; Merlin, O.; Amazirh, A.; Ouaadi, N.; Rivalland, V.; Jarlan, L.; Er-Raki, S.; Escorihuela, M. A Calibration/Disaggregation Coupling Scheme for Retrieving Soil Moisture at High Spatio-Temporal Resolution: Synergy between SMAP Passive Microwave, MODIS/Landsat Optical/Thermal and Sentinel-1 Radar Data. Sensors 2021, 21, 7406. [Google Scholar] [CrossRef]
  10. Liu, X.; Shao, Y.; Liu, L.; Li, K.; Wang, J.; Li, S.; Wang, J.; Wu, X. Effects of Plant Crown Shape on Microwave Backscattering Coefficients of Vegetation Canopy. Sensors 2021, 21, 7748. [Google Scholar] [CrossRef]
  11. He, C.; He, H.; Chang, J.; Dong, Y.; Liu, S.; Zeng, N.; He, Y.; Ma, H. Characterizing microstructures of cancerous tissues using multispectral transformed Mueller matrix polarization parameters. Biomed. Opt. Express 2015, 6, 2934–2945. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  12. Gruev, V.; Perkins, R.; York, T. CCD polarization imaging sensor with aluminum nanowire optical filters. Opt. Express 2010, 18, 19087–19094. [Google Scholar] [CrossRef] [PubMed]
  13. Andreou, A.; Kalayjian, Z. Polarization imaging: Principles and integrated polarimeters. IEEE Sens. J. 2002, 2, 566–576. [Google Scholar] [CrossRef]
  14. Zhao, X.; Bermak, A.; Boussaid, F.; Chigrinov, V. Liquid-crystal micropolarimeter array for full Stokes polarization imaging in visible spectrum. Opt. Express 2010, 18, 17776–17787. [Google Scholar] [CrossRef] [PubMed]
  15. Tyo, J.; Goldstein, D.; Chenault, D.; Shaw, J. Review of passive imaging polarimetry for remote sensing applications. Appl. Opt. 2006, 45, 5453–5469. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Lucid Vision TELEDYNE FLIR. Blackfly S USB3. Available online: https://www.flir.cn/products/blackfly-s-usb3/?model=BFS-U3-51S5P-C (accessed on 16 July 2021).
  17. 4D Technology Imaging Polarimeters. Available online: https://www.4dtechnology.com/products/imaging-polarimeters (accessed on 16 July 2021).
  18. Lucid Vision Labs Phoenix Network Camera Suitable for OEM Ultra-Small and Deformable. Available online: http://thinklucid.cn/phoenix-machine-vision (accessed on 16 July 2021).
  19. Hou, H.; Andrews, H. Cubic splines for image interpolation and digital filtering. IEEE Trans. Acoust. Speech Signal Process. 1978, 26, 508–517. [Google Scholar]
  20. Li, M.; Nguyen, T. Markov Random Field Model-Based Edge-Directed Image Interpolation. IEEE Trans. Image Process. 2008, 17, 1121–1128. [Google Scholar]
  21. Pekkucuksen, I.; Altunbasak, Y. Multiscale Gradients-Based Color Filter Array Interpolation. IEEE Trans. Image Process. 2013, 22, 157–165. [Google Scholar] [CrossRef]
  22. He, K.; Sun, J.; Tang, X. Guided Image Filtering. IEEE Trans. Pattern Anal. Mach. Intell. 2013, 35, 1397–1409. [Google Scholar] [CrossRef]
  23. Kiku, D.; Monno, Y.; Tanaka, M.; Okutomi, M. Residual interpolation for color image demosaicking. In Proceedings of the IEEE International Conference on Image Processing (ICIP), Melbourne, Australia, 15–18 September 2013. [Google Scholar]
  24. Ye, W.; Ma, K. Image demosaicing by using iterative residual interpolation. In Proceedings of the IEEE International Conference on Image Processing (ICIP), Paris, France, 27–39 October 2014. [Google Scholar]
  25. Kiku, D.; Monno, Y.; Tanaka, M.; Okutomi, M. Minimized-Laplacian Residual Interpolation for Color Image Demosaicking. Proc. SPIE-Int. Soc. Opt. Eng. 2014, 9023, 90230L. [Google Scholar]
  26. Kiku, D.; Monno, Y.; Tanaka, M.; Okutomi, M. Beyond Color Difference: Residual Interpolation for Color Image Demosaicking. IEEE Trans. Image Process. 2016, 25, 1288–1300. [Google Scholar] [CrossRef] [PubMed]
  27. Kiku, D.; Monno, Y.; Tanaka, M.; Okutomi, M. Adaptive Residual Interpolation for Color and Multispectral Image Demosaicking. Sensors 2017, 17, 2787. [Google Scholar]
  28. Ratliff, B.; Boger, J.; Fetrow, M.; Tyo, J.; Black, W. Image processing methods to compensate for IFOV errors in microgrid imaging polarimeters. Proc. SPIE 2006, 6240, 139–174. [Google Scholar]
  29. Ratliff, B.; Shaw, J.; Tyo, J.; Boger, J.; Black, W.; Bowers, D.; Kumar, R. Mitigation of image artifacts in LWIR microgrid polarimeter images. Proc. SPIE-Int. Soc. Opt. Eng. 2007, 6682, 668209. [Google Scholar]
  30. Ratliff, B.; Charles, C.; Tyo, J. Interpolation strategies for reducing IFOV artifacts in microgrid polarimeter imagery. Opt. Express 2009, 17, 9112–9125. [Google Scholar] [CrossRef]
  31. York, T.; Powell, S.; Gruev, V. A comparison of polarization image processing across different platforms. Proc. SPIE 2011, 8160, 816004. [Google Scholar]
  32. Gao, S.; Gruev, V. Image interpolation methods evaluation for division of focal plane polarimeters. Proc. SPIE-Int. Soc. Opt. Eng. 2011, 8012, 150–154. [Google Scholar]
  33. Gao, S.; Gruev, V. Bilinear and bicubic interpolation methods for division of focal plane polarimeters. Opt. Express 2011, 19, 26161–26173. [Google Scholar] [CrossRef]
  34. Gao, S.; Gruev, V. Gradient based interpolation for division of focal plane polarization imaging sensor. In Proceedings of the IEEE International Symposium on Circuits and Systems (ISCAS), Seoul, Korea, 20–23 May 2012. [Google Scholar]
  35. Gao, S.; Gruev, V. Gradient-based interpolation method for division-of-focal-plane polarimeters. Opt. Express 2013, 21, 1137–1151. [Google Scholar] [CrossRef]
  36. Zhang, J.; Ye, W.; Ahmed, A.; Qiu, Z.; Cao, Y.; Zhao, X. A Novel Smoothness-Based Interpolation Algorithm for Division of Focal Plane Polarimeters. In Proceedings of the IEEE International Symposium on Circuits & Systems (ISCAS), Baltimore, MD, USA, 28–31 May 2017. [Google Scholar]
  37. Xu, X.; Kulkarni, M.; Nehorai, A.; Gruev, V. A correlation-based interpolation algorithm for division-of-focal-plane polarization sensors. Proc. SPIE 2012, 8364, 83640L. [Google Scholar]
  38. Zhang, J.; Luo, H.; Hui, B.; Chang, Z. Image interpolation for division of focal plane polarimeters with intensity correlation. Opt. Express 2016, 24, 20799–20807. [Google Scholar] [CrossRef] [PubMed]
  39. Ning, L.; Zhao, Y.; Quan, P.; Kong, S. Demosaicking DoFP images using Newton’s polynomial interpolation and polarization difference model. Opt. Express 2019, 27, 1376–1391. [Google Scholar]
  40. Wu, R.; Zhao, Y.; Li, N.; Kong, S. Polarization image demosaicking using polarization channel difference prior. Opt. Express 2021, 29, 22066–22079. [Google Scholar] [CrossRef] [PubMed]
  41. Ahmed, A.; Zhao, X.; Gruev, V.; Chang, J.; Bermak, A. Residual interpolation for division of focal plane polarization image sensors. Opt. Express 2017, 25, 10651–10662. [Google Scholar] [CrossRef]
  42. Ahmed, A.; Zhao, X.; Chang, J.; Ma, H.; Gruev, V.; Bermak, A. Four-Directional Adaptive Residual Interpolation Technique for DoFP Polarimeters with Different Micro-polarizer Patterns. IEEE Sens. J. 2018, 18, 7990–7997. [Google Scholar] [CrossRef]
  43. Jiang, T.; Wen, D.; Song, Z.; Zhang, W.; Li, Z.; Wei, X.; Liu, G. Minimized Laplacian residual interpolation for DoFP polarization image demosaicking. Appl. Opt. 2019, 58, 7367–7374. [Google Scholar] [CrossRef] [PubMed]
  44. Morimatsu, M.; Monno, Y.; Tanaka, M.; Okutomi, M. Monochrome and Color Polarization Demosaicking Using Edge-Aware Residual Interpolation. In Proceedings of the 27th IEEE International Conference on Image Processing (ICIP), Abu Dhabi, United Arab Emirates, 25–28 October 2020. [Google Scholar]
  45. Liu, S.; Chen, J.; Xun, Y.; Zhao, X.; Chang, C. A New Polarization Image Demosaicking Algorithm by Exploiting Inter-Channel Correlations with Guided Filtering. IEEE Trans. Image Process. 2020, 29, 7076–7089. [Google Scholar] [CrossRef]
  46. Połap, D.; Srivastava, G. Neural image reconstruction using a heuristic validation mechanism. Neural Comput. Appl. 2021, 33, 10787–10797. [Google Scholar] [CrossRef]
  47. Qiu, S.; Fu, Q.; Wang, C.; Heidrich, W. Polarization demosaicking for monochrome and color polarization focal plane arrays. In Proceedings of the International Symposium on Vision, Modeling, and Visualization (VMV), Rostock, Germany, 30 September 2019. [Google Scholar]
  48. Zhang, J.; Luo, H.; Liang, R.; Ahmed, A.; Zhang, X.; Hui, B.; Chang, Z. Sparse representation-based demosaicing method for microgrid polarimeter imagery. Opt. Lett. 2018, 43, 3265–3268. [Google Scholar] [CrossRef]
  49. Zhang, J.; Shao, J.; Luo, H.; Zhang, X.; Hui, B.; Chang, Z.; Liang, R. Learning a convolutional demosaicing network for microgrid polarimeter imagery. Opt. Lett. 2018, 43, 4534–4537. [Google Scholar] [CrossRef]
  50. Wen, S.; Zheng, Y.; Lu, F.; Zhao, Q. Convolutional demosaicing network for joint chromatic and polarimetric imagery. Opt. Lett. 2019, 44, 5646–5649. [Google Scholar] [CrossRef]
  51. Zeng, X.; Luo, Y.; Zhao, X.; Ye, W. An end-to-end fully-convolutional neural network for division of focal plane sensors to reconstruct S0, DoLP, and AoP. Opt. Express 2019, 27, 8566–8577. [Google Scholar] [CrossRef] [PubMed]
  52. Yang, J.; Qiu, S.; Jin, W.; Xue, F. Temporal and spatial error model for estimating the measurement precision of the division of focal plane polarimeters. Opt. Express 2021, 29, 20808–20828. [Google Scholar] [CrossRef]
  53. Mihoubi, S.; Lapray, P.; Bigué, L. Survey of Demosaicking Methods for Polarization Filter Array Images. Sensors 2018, 18, 3688. [Google Scholar] [CrossRef] [Green Version]
  54. Stokes, G. On the Composition and Resolution of Streams of Polarized Light from different Sources. Trans. Camb. Philos. Soc. 1852, 9, 399–416. [Google Scholar]
  55. Lapray, P.; Gendre, L.; Foulonneau, A.; Bigue, L. Database of polarimetric and multispectral images in the visible and NIR regions. Proc. SPIE 2018, 10677, 1067738. [Google Scholar]
Figure 1. Schematic diagrams of DoFP polarimeter and the color filter array. (a) Illustrates the structure of DoFP polarimeters; (b) presents the spatial layout of the pixeled polarizer array and the color filter array.
Figure 1. Schematic diagrams of DoFP polarimeter and the color filter array. (a) Illustrates the structure of DoFP polarimeters; (b) presents the spatial layout of the pixeled polarizer array and the color filter array.
Sensors 22 01529 g001
Figure 2. The framework of the residual interpolation methods for DoFP polarimeters demosaicing. Guided Filter represents RI or MLRI.
Figure 2. The framework of the residual interpolation methods for DoFP polarimeters demosaicing. Guided Filter represents RI or MLRI.
Sensors 22 01529 g002
Figure 3. The PSNR of the guide image generated by the 10 up-sampling filters. The abscissa represents the image number in the dataset.
Figure 3. The PSNR of the guide image generated by the 10 up-sampling filters. The abscissa represents the image number in the dataset.
Sensors 22 01529 g003
Figure 4. The variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image with the PSNR of the guide image. The results of 0°, 90°, and 135° are similar to that of 45°, so they are not repeatedly exhibited for conciseness. (a,b) represents the variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by RI with the PSNR of the guide image, respectively. (c,d) represents the variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by MLRI with the PSNR of the guide image, respectively.
Figure 4. The variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image with the PSNR of the guide image. The results of 0°, 90°, and 135° are similar to that of 45°, so they are not repeatedly exhibited for conciseness. (a,b) represents the variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by RI with the PSNR of the guide image, respectively. (c,d) represents the variation in the PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by MLRI with the PSNR of the guide image, respectively.
Sensors 22 01529 g004
Figure 5. The guide images with different sampling rates.
Figure 5. The guide images with different sampling rates.
Sensors 22 01529 g005
Figure 6. The PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by three guide images with different sampling rates. The results of 0°, 90°, and 135° are similar to that of 45°, so they are not repeatedly exhibited for conciseness. (a,b) illustrate the PSNR of I 45 ° H R and the DoLP calculated by RI, respectively. (c,d) present the PSNR of I 45 ° H R and the DoLP calculated by MLRI, respectively.
Figure 6. The PSNR of the high-resolution finial output I 45 ° H R and the DoLP image calculated by three guide images with different sampling rates. The results of 0°, 90°, and 135° are similar to that of 45°, so they are not repeatedly exhibited for conciseness. (a,b) illustrate the PSNR of I 45 ° H R and the DoLP calculated by RI, respectively. (c,d) present the PSNR of I 45 ° H R and the DoLP calculated by MLRI, respectively.
Sensors 22 01529 g006
Figure 7. The filter windows with different directions in the guided filter.
Figure 7. The filter windows with different directions in the guided filter.
Sensors 22 01529 g007
Figure 8. The overall pipeline of the proposed PAIPRI. Horizontal, Vertical, Diagonal Adaptive Iterative Filters represent pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal, vertical, and two diagonal directions, respectively. Weighted average combiner represents fusion on the interpolation images with RI and MLRI in the horizontal, vertical and two diagonal directions.
Figure 8. The overall pipeline of the proposed PAIPRI. Horizontal, Vertical, Diagonal Adaptive Iterative Filters represent pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal, vertical, and two diagonal directions, respectively. Weighted average combiner represents fusion on the interpolation images with RI and MLRI in the horizontal, vertical and two diagonal directions.
Sensors 22 01529 g008
Figure 9. The overall pipeline of the pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal direction. Primary branch represents the branch up-sampling I, whose final output interpolation results were I 0 ° , R I H and I 0 ° , M L R I H . Auxiliary branch, the branch up-sampling the guide image I45°, whose purpose is updating the guide image of the primary branch and increasing the PSNR of the guide image to further increase the PSNR of I 0 ° , R I H and I 0 ° , M L R I H . HLI represents the horizontal linear interpolation. GF represents RI or MLRI. Adaptive pixel updater represents the process of adaptively updating the iterative results, pixel by pixel.
Figure 9. The overall pipeline of the pixel-by-pixel adaptive iterative processes based on residual interpolation in horizontal direction. Primary branch represents the branch up-sampling I, whose final output interpolation results were I 0 ° , R I H and I 0 ° , M L R I H . Auxiliary branch, the branch up-sampling the guide image I45°, whose purpose is updating the guide image of the primary branch and increasing the PSNR of the guide image to further increase the PSNR of I 0 ° , R I H and I 0 ° , M L R I H . HLI represents the horizontal linear interpolation. GF represents RI or MLRI. Adaptive pixel updater represents the process of adaptively updating the iterative results, pixel by pixel.
Sensors 22 01529 g009
Figure 10. The reconstructed results of the image numbered 4 in the tested dataset, generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a knife composed of a metal body and wooden handle. The background is a rough wall and desktop. To more clearly illustrate the visual differences in the reconstructed results of the seven comparison methods, we zoomed in on a local area marked by a red box with the size of 30 × 40. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the single arc-shaped edges reconstructed by Bilinear, BS, Gradient [35], PRI [41] and EARI [44]. The similar results of PRI [41] and bilinear methods, which are also illustrated in Table 1, Table 2 and Table 3, again confirm that it is inappropriate to choose the same polarization direction for the guide image and the input image of the guided filter (as discussed in Section 3.2). NP [39], MLPRI [43] and the proposed PAIPRI can reconstruct sharp edges. However, the reconstructed results of MLPRI [43] appear as excessive smoothing in the target. The reconstructed results of NP [39] generate additional error messages. Although the reconstructed results of PAIPRI retain a small amount of mosaic effect on the edges, it is visually the closest to the ground-truth images.
Figure 10. The reconstructed results of the image numbered 4 in the tested dataset, generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a knife composed of a metal body and wooden handle. The background is a rough wall and desktop. To more clearly illustrate the visual differences in the reconstructed results of the seven comparison methods, we zoomed in on a local area marked by a red box with the size of 30 × 40. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the single arc-shaped edges reconstructed by Bilinear, BS, Gradient [35], PRI [41] and EARI [44]. The similar results of PRI [41] and bilinear methods, which are also illustrated in Table 1, Table 2 and Table 3, again confirm that it is inappropriate to choose the same polarization direction for the guide image and the input image of the guided filter (as discussed in Section 3.2). NP [39], MLPRI [43] and the proposed PAIPRI can reconstruct sharp edges. However, the reconstructed results of MLPRI [43] appear as excessive smoothing in the target. The reconstructed results of NP [39] generate additional error messages. Although the reconstructed results of PAIPRI retain a small amount of mosaic effect on the edges, it is visually the closest to the ground-truth images.
Sensors 22 01529 g010
Figure 11. The reconstructed results of the image numbered 8 in the tested dataset generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a standard color checker marked with a white brand logo. We zoomed in on a local area marked by a red box with the size of 90 × 120. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the edges reconstructed by Bilinear and PRI [41]. The reconstructed results of BS, Gradient [35], and MLPRI [43] demonstrate blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate a high number of additional error messages in flat-field regions. EARI [44] enhances the horizontal and vertical edges, but the sawtooth effect of edges in other directions is still obvious. However, the proposed PAIPRI can reconstruct clear and sharp edges. Although the reconstructed results of PAIPRI retain a small amount of mosaic effect on horizontal and vertical edges, it is visually the closest to the ground-truth images.
Figure 11. The reconstructed results of the image numbered 8 in the tested dataset generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a standard color checker marked with a white brand logo. We zoomed in on a local area marked by a red box with the size of 90 × 120. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the edges reconstructed by Bilinear and PRI [41]. The reconstructed results of BS, Gradient [35], and MLPRI [43] demonstrate blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate a high number of additional error messages in flat-field regions. EARI [44] enhances the horizontal and vertical edges, but the sawtooth effect of edges in other directions is still obvious. However, the proposed PAIPRI can reconstruct clear and sharp edges. Although the reconstructed results of PAIPRI retain a small amount of mosaic effect on horizontal and vertical edges, it is visually the closest to the ground-truth images.
Sensors 22 01529 g011
Figure 12. The reconstructed results of the image numbered 1 in the tested dataset generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a fabric with abundant texture features. We zoomed in on a local area marked by a red box with the size of 100 × 100. (a) is the I image. (b) is the DoLP image. Bilinear, BS, Gradient [35], PRI [41], MLPRI [43], and EARI [44] cannot reconstruct correct texture features in DoLP images. The reconstructed results of NP [39] demonstrate excessive smoothing. However, the reconstructed results of the proposed PAIPRI can basically reconstruct the correct texture features, and is visually the closest to the ground-truth images.
Figure 12. The reconstructed results of the image numbered 1 in the tested dataset generated by the eight demosaicing methods for DoFP polarimeters [35,39,41,43,44]. The target is a fabric with abundant texture features. We zoomed in on a local area marked by a red box with the size of 100 × 100. (a) is the I image. (b) is the DoLP image. Bilinear, BS, Gradient [35], PRI [41], MLPRI [43], and EARI [44] cannot reconstruct correct texture features in DoLP images. The reconstructed results of NP [39] demonstrate excessive smoothing. However, the reconstructed results of the proposed PAIPRI can basically reconstruct the correct texture features, and is visually the closest to the ground-truth images.
Sensors 22 01529 g012
Figure 13. The relationship between the PSNR of I and the number of iterations.
Figure 13. The relationship between the PSNR of I and the number of iterations.
Sensors 22 01529 g013
Figure 14. The reconstructed results generated by the eight demosaicing methods on the indoor scene images collected by a real-world DoFP polarimeter [35,39,41,43,44]. The target is a metal tank model with abundant high-frequency information. We zoomed in a local area marked by a red box with the size of 125 × 150. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the arc-shaped edges reconstructed by Bilinear, BS, PRI [41], and EARI [44]. The reconstructed results of Gradient [35] demonstrate blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate additional error messages. MLPRI [43] and the proposed PAIPRI can basically reconstruct clear edges, but MLPRI retains some sawtooth effects at the left edge of the wheel. Therefore, PAIPRI is visually the best compared with the other six methods.
Figure 14. The reconstructed results generated by the eight demosaicing methods on the indoor scene images collected by a real-world DoFP polarimeter [35,39,41,43,44]. The target is a metal tank model with abundant high-frequency information. We zoomed in a local area marked by a red box with the size of 125 × 150. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the arc-shaped edges reconstructed by Bilinear, BS, PRI [41], and EARI [44]. The reconstructed results of Gradient [35] demonstrate blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate additional error messages. MLPRI [43] and the proposed PAIPRI can basically reconstruct clear edges, but MLPRI retains some sawtooth effects at the left edge of the wheel. Therefore, PAIPRI is visually the best compared with the other six methods.
Sensors 22 01529 g014
Figure 15. The reconstructed results generated by the eight demosaicing methods on the indoor scene images collected by a real-world DoFP polarimeter [35,39,41,43,44]. The target is a metal tank model with abundant high-frequency information. We zoomed in on a local area marked by a red box with the size of 125 × 150. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the arc-shaped edges reconstructed by Bilinear, BS, PRI [41], and EARI [44]. The reconstructed results of Gradient [35] present blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate additional error messages. MLPRI [43] and the proposed PAIPRI can basically reconstruct clear edges, but MLPRI [43] retains some sawtooth effects at the left edge of the wheel. Therefore, PAIPRI is visually the best compared with the other six methods.
Figure 15. The reconstructed results generated by the eight demosaicing methods on the indoor scene images collected by a real-world DoFP polarimeter [35,39,41,43,44]. The target is a metal tank model with abundant high-frequency information. We zoomed in on a local area marked by a red box with the size of 125 × 150. (a) is the I image. (b) is the DoLP image. There are serious sawtooth effects at the arc-shaped edges reconstructed by Bilinear, BS, PRI [41], and EARI [44]. The reconstructed results of Gradient [35] present blurred edges due to excessive smoothing. The reconstructed results of NP [39] generate additional error messages. MLPRI [43] and the proposed PAIPRI can basically reconstruct clear edges, but MLPRI [43] retains some sawtooth effects at the left edge of the wheel. Therefore, PAIPRI is visually the best compared with the other six methods.
Sensors 22 01529 g015
Table 1. PSNR of I Reconstructed by Different Methods on Dataset.
Table 1. PSNR of I Reconstructed by Different Methods on Dataset.
Image NumberBilinearBSGradient [35]NP [39]PRI [41]MLPRI [43]EARI [44]PAIPRI
138.5639.4038.3941.2038.5639.4038.6441.40
247.7747.7647.7742.9947.7847.7639.9748.21
341.5142.4541.9338.1041.5442.5038.4442.97
442.0042.6542.4439.2142.0442.7137.5844.42
540.9340.6840.1636.9440.9340.6841.1242.14
642.4143.2442.8038.5142.4343.3140.5244.81
742.5443.4642.2737.5242.5443.4642.5345.00
841.6542.6441.4630.9541.6542.6441.4844.44
931.4533.4930.6227.9731.4533.5231.4634.38
1039.4440.4939.3434.5539.4540.5237.8042.05
Average40.8341.6340.7236.7940.8441.6538.9542.98
Table 2. PSNR of S2 Reconstructed by Different Methods on Dataset.
Table 2. PSNR of S2 Reconstructed by Different Methods on Dataset.
Image NumberBilinearBSGradient [35]NP [39]PRI [41]MLPRI [43]EARI [44]PAIPRI
139.6642.0340.6343.0739.6642.5339.7543.87
248.9549.5649.1943.2348.9637.6640.1549.59
343.6945.6045.3233.1843.7145.4339.8645.87
442.6844.3844.0838.8842.7038.8037.9245.78
545.2045.3644.8037.0645.2046.0345.1147.30
642.8644.9543.9335.7242.8745.2341.5946.49
743.9647.2244.1638.3543.9647.6943.9348.25
841.0944.1341.2437.7341.0945.1241.0546.12
932.8235.2433.0527.0432.8236.4132.8636.78
1040.0243.2840.8633.1540.0340.1138.2744.86
Average42.0944.1842.7336.7442.1042.5040.0545.49
Table 3. PSNR of DoLP Reconstructed by Different Methods on Dataset.
Table 3. PSNR of DoLP Reconstructed by Different Methods on Dataset.
Image NumberBilinearBSGradient [35]NP [39]PRI [41]MLPRI [43]EARI [44]PAIPRI
141.0142.6741.6742.9541.0143.0640.5543.62
237.6437.7037.8028.5837.6533.8435.9437.91
339.8340.9640.9234.1239.8340.9938.8041.27
437.0237.9037.8536.0837.0337.1235.7939.38
531.6631.4931.3729.2231.6631.5430.7731.93
634.3434.6634.5823.0934.3534.8233.0435.66
731.2332.5231.2620.0931.2332.4030.8033.19
831.3733.0331.5818.0331.3733.2731.1234.35
928.3429.7728.7816.3728.3430.1728.8130.72
1030.9831.5631.0425.1730.9831.3930.3732.07
Average34.3435.2334.6827.3734.3434.8633.6036.01
Table 4. Computation Times of Different Methods on Dataset Image Numbered 1.
Table 4. Computation Times of Different Methods on Dataset Image Numbered 1.
BilinearBSGradient [35]NP [39]PRI [41]MLPRI [43]EARI [44]PAIPRI(1)PAIPRI(3)PAIPRI(5)PAIPRI(10)
Time (s)0.071.041.411.371.322.982.526.8222.4840.84101.99
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Yang, J.; Jin, W.; Qiu, S.; Xue, F.; Wang, M. Residual Interpolation Integrated Pixel-by-Pixel Adaptive Iterative Process for Division of Focal Plane Polarimeters. Sensors 2022, 22, 1529. https://doi.org/10.3390/s22041529

AMA Style

Yang J, Jin W, Qiu S, Xue F, Wang M. Residual Interpolation Integrated Pixel-by-Pixel Adaptive Iterative Process for Division of Focal Plane Polarimeters. Sensors. 2022; 22(4):1529. https://doi.org/10.3390/s22041529

Chicago/Turabian Style

Yang, Jie, Weiqi Jin, Su Qiu, Fuduo Xue, and Meishu Wang. 2022. "Residual Interpolation Integrated Pixel-by-Pixel Adaptive Iterative Process for Division of Focal Plane Polarimeters" Sensors 22, no. 4: 1529. https://doi.org/10.3390/s22041529

APA Style

Yang, J., Jin, W., Qiu, S., Xue, F., & Wang, M. (2022). Residual Interpolation Integrated Pixel-by-Pixel Adaptive Iterative Process for Division of Focal Plane Polarimeters. Sensors, 22(4), 1529. https://doi.org/10.3390/s22041529

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop