Scene-level passive polarization 3D imaging

Scene-level passive polarization 3D imaging

PDF

Xin Wang¹, Pingli Han¹, Xiyuan Luo¹, Qianqian Liu¹, Tong Zhang¹, Xue Dong¹, Meng Xiang¹, Jinpeng Liu¹, Yanyan Liu², Fei Liu¹^,³^,^*

Opto-Electronic Advances | 2026, 9(2) : 250267-1 - 250267-6

Less

Opto-Electronic Advances | 2026, 9(2): 250267-1-250267-6

• Supplementary information •

Scene-level passive polarization 3D imaging

Full

Xin Wang¹, Pingli Han¹, Xiyuan Luo¹, Qianqian Liu¹, Tong Zhang¹, Xue Dong¹, Meng Xiang¹, Jinpeng Liu¹, Yanyan Liu², Fei Liu¹^,³^,^*

Affiliations

¹School of Optoelectronic Engineering, Xidian University, Xi’an 710071, China

²National Key Laboratory of Electromagnetic Space Security, Tianjin 300308, China

³State Key Laboratory of Electromechanical Integrated Manufacturing of High-Performance Electronic Equipments, Xidian University, Xi’an 710071, China

About Author:

^†

These authors contributed equally to this work.

Published: 2026-02-12 doi: 10.29026/oea.2026.250267

Outline

Cite this Article

Xin Wang, Pingli Han, Xiyuan Luo, Qianqian Liu, Tong Zhang, Xue Dong, Meng Xiang, Jinpeng Liu, Yanyan Liu, Fei Liu. Scene-level passive polarization 3D imaging[J]. Opto-Electronic Advances, 2026 , 9 (2) : 250267-1 -250267-6 . DOI: 10.29026/oea.2026.250267

Full Text

Less

Section 1:　Polarization characteristics in natural scenes

Less

Since our research focused on 3D imaging of natural scenes based on polarization. It’s essential to first test the existence of available polarization information. The tests should cover key elements of natural scenes. Vegetation, roads, and buildings constitute fundamental elements of natural scenes, spanning typical categories from natural organic to artificial objects. Data are illustrated in Fig. S1. Overall, the polarimetric characteristics of the ground objects are relatively weak. But the inverted azimuth and zenith angles still contain structural information of targets. This information will be encoded into the surface normal during subsequent reconstruction.

Section 2:　Surface normal correction

Less

Polarization 3D establishes a mapping that correlates the polarization of the light field and the 3D morphology. The essence lies in the comprehensive exploitation of the degree of polarization (DoP) and the angle of polarization (AoP) of the light field detected by the detector, thereby facilitating the determination of the zenith θ and azimuth Φ of the surface normal field. We express surface normal by Eq. (S1),

$ {\boldsymbol{N}}_{{\mathrm{P}}}=\left\{\begin{array}{c}{N}_{{\mathrm{P}}x}\\{N}_{{\mathrm{P}}y}\\{N}_{{\mathrm{P}}z}\end{array}\right\}=\left\{\begin{array}{c}\sin \theta \cos \phi \\\sin \theta \sin \phi \\\cos \theta \end{array}\right\}\;,\tag{S1} $

(1)

where N_Px, N_Py and N_Pz denote the components of the surface normal in x, y and z directions. The trigonometric functions in the inverse solution of the AoP lead to a possible

$ \phi $

phase ambiguity of π (

$ \begin{array}{ccc}\phi =\varphi & {\mathrm{or}} & \phi =\varphi +\pi \end{array} $

) and further surface normal ambiguity. Moreover, given the unit norm length of the normal vector got from Eq. (S1), the integration to surface normal only provides a measure of relative height of the surface.

To realize polarization 3D from 2D images with M×N pixels, we need to accurately determine 2MN unknowns (θ, φ). The depth of the scene derived from stereo vision, typically aligns with the consistent variation trend observed in the polarization normal field information. Assume the gradients of this depth map along the X and Y axes are represented as

$ g={[\partial Z/\partial X,\partial Z/\partial Y]}^{{\mathrm{T}}} $

. To facilitate subsequent algorithmic processing, we convert this gradient vector field into a normal vector field perpendicular to the micro-facets and perform normalization, which can be expressed as:

$ {\boldsymbol{N}}_{{\mathrm{B}}}=\left[\begin{array}{c}{N}_{{\mathrm{B}}x}\\{N}_{{\mathrm{B}}y}\\{N}_{{\mathrm{B}}z}\end{array}\right]=\left[\begin{array}{c}-g\\1\end{array}\right]=\frac{1}{\sqrt{{(\partial Z/\partial x)}^{2}+{(\partial Z/\partial y)}^{2}+{1}^{2}}}\left[\begin{array}{c}-\partial Z/\partial x\\-\partial Z/\partial y\\1\end{array}\right]\;,\tag{S2} $

(2)

where N_Bx is the normal from binocular depth, N_Bx, N_By and N_Bz are the components of N_B along three directions.

We optimize the scene's normal distribution by aligning the stereo-vision-derived normal field with the polarization data, using a tailored optimization function in Eq. (S3),

$ \begin{cases} {N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}={N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}\,\,\,\,\,\,\,\,\,\,\,\,\,{\mathrm{if}}\,\,\,{N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}\times {N}_{{{{\mathrm{B}}}_{x}}{\mathrm{or}}{{{\mathrm{B}}}_{y}}}\geq 0\\{N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}=-{N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}\,\,\,\,\,\,\,\,\,\,\,{\mathrm{else}}\,\,\,{\mathrm{if}}\,\,\,{N}_{{{{\mathrm{P}}}_{x}}{\mathrm{or}}{{{\mathrm{P}}}_{y}}}\times {N}_{{{{\mathrm{B}}}_{x}}{\mathrm{or}}{{{\mathrm{B}}}_{y}}} \lt 0\end{cases}\;.\tag{S3} $

(3)

Section 3:　Optimization by iteration

Less

As stated in the main text, the problem we seek to solve is expressed by Eq. (S4),

$ \begin{split}L({\boldsymbol{S}},{{\boldsymbol{J}}}_{1},{{\boldsymbol{J}}}_{2},{{\boldsymbol{L}}}_{1},{{\boldsymbol{L}}}_{2})&=\frac{1}{2}\left|\left|{\boldsymbol{D}}-{\boldsymbol{S}}\right|\right|_{2}^{2}+\frac{\xi }{2}\left|\left|{{\boldsymbol{N}}}_{P}-{{\boldsymbol{J}}}_{1}\right|\right|_{2}^{2}+\eta {\left|\left|{{\boldsymbol{J}}}_{2}\right|\right|}_{\ast }-\left\langle {{\boldsymbol{L}}}_{1},{{\boldsymbol{J}}}_{1}-\nabla {\boldsymbol{S}}\right\rangle -\left\langle {{\boldsymbol{L}}}_{2},{{\boldsymbol{J}}}_{2}-{\boldsymbol{FS}}\right\rangle \\&\quad+\frac{\rho }{2}(\left|\left|{{\boldsymbol{J}}}_{1}-\nabla {\boldsymbol{S}}\right|\right|_{2}^{2}+\left|\left|{{\boldsymbol{J}}}_{2}-{\boldsymbol{FS}}\right|\right|_{2}^{2})\end{split}\;.\tag{S4} $

(4)

We reformulated the original problem into an augmented form by introducing Lagrange multipliers. This was done to facilitate variable separation for alternating iterations, simplifying the original problem. The problem is solved by ADMM where we first derive the optimality conditions for each primal variable, assuming the others are fixed:

$ \left\{\begin{array}{l} {{\boldsymbol{S}}}^{k+1}=\underset{S}{\arg \min } \dfrac{1}{2}\left|\left|{\boldsymbol{D}}-{\boldsymbol{S}}\right|\right|_{2}^{2}+\left\langle {{\boldsymbol{L}}}_{1}{}^{k},{{\boldsymbol{J}}}_{1}{}^{k}-\nabla {\boldsymbol{S}}\right\rangle +\left\langle {{\boldsymbol{L}}}_{2}{}^{k},{{\boldsymbol{J}}}_{2}{}^{k}-{\boldsymbol{FS}}\right\rangle \\\quad\;\;\;\;\;\;+\dfrac{\rho }{2}(\left|\left|{{\boldsymbol{J}}}_{1}{}^{k}-\nabla {\boldsymbol{S}}\right|\right|_{2}^{2}+\left|\left|{{\boldsymbol{J}}}_{2}{}^{k}-{\boldsymbol{FS}}\right|\right|_{2}^{2})\\{\boldsymbol{J}}_{1}^{k+1}=\underset{{{\boldsymbol{J}}}_{1}}{\arg \min } \dfrac{\xi }{2}\left|\left|{{\boldsymbol{N}}}_{P}-{J}_{1}\right|\right|_{2}^{2}+\left\langle {{\boldsymbol{L}}}_{1}{}^{k},{{\boldsymbol{J}}}_{1}-\nabla {S}^{k+1}\right\rangle +\dfrac{\rho }{2}\left|\left|{{\boldsymbol{J}}}_{1}-\nabla {{\boldsymbol{S}}}^{k+1}\right|\right|_{2}^{2}\\{\boldsymbol{J}}_{2}^{k+1}=\underset{{{\boldsymbol{J}}}_{2}}{\arg \min } \eta {\left|\left|{{\boldsymbol{J}}}_{2}\right|\right|}_{\ast }+\left\langle {{\boldsymbol{L}}}_{2}{}^{k},{{\boldsymbol{J}}}_{2}-{\boldsymbol{F}}{{\boldsymbol{S}}}^{k+1}\right\rangle +\dfrac{\rho }{2}\left|\left|{{\boldsymbol{J}}}_{2}-{\boldsymbol{F}}{{\boldsymbol{S}}}^{k+1}\right|\right|_{2}^{2}\\{\boldsymbol{L}}_{1}^{k+1}={\boldsymbol{L}}_{1}^{k}+\rho ({\boldsymbol{J}}_{1}^{k+1}-\nabla {{\boldsymbol{S}}}^{k\mathrm+1})\\{\boldsymbol{L}}_{2}^{k+1}={\boldsymbol{L}}_{2}^{k}+\rho ({\boldsymbol{J}}_{2}^{k+1}-{\boldsymbol{F}}{{\boldsymbol{S}}}^{k\mathrm+1})\end{array}\right.\;.\tag{S5} $

(5)

By updating each variable, we get the solution at iteration k expressed by Eqs. (S6–S10),

$ {{\boldsymbol{S}}}^{k+1}={(1+\rho {{\nabla }^{{\mathrm{T}}}}\nabla +\rho {{{\boldsymbol{F}}}^{{\mathrm{T}}}}{\boldsymbol{F}})}^{-1}(D+{\nabla }^{{\mathrm{T}}}{\boldsymbol{L}}_{1}^{k}+{{\boldsymbol{F}}}^{{\mathrm{T}}}{\boldsymbol{L}}_{2}^{k}+\rho {\nabla }^{{\mathrm{T}}}{\boldsymbol{J}}_{1}^{k}+\rho {{\boldsymbol{F}}}^{{\mathrm{T}}}{\boldsymbol{J}}_{2}^{k})\;,\tag{S6} $

(6)

$ {\boldsymbol{J}}_{1}^{k+1}=\frac{\xi {{\boldsymbol{N}}}_{P}+\rho \nabla {{\boldsymbol{S}}}^{k+1}-{\boldsymbol{L}}_{1}^{k}}{\xi +\rho }\;,\tag{S7} $

(7)

$ {\boldsymbol{J}}_{2}^{k+1}={{\boldsymbol{S}}}_{\eta /\rho }({\boldsymbol{F}}{{\boldsymbol{S}}}^{k+1}-{\boldsymbol{L}}_{2}^{k}/\eta )\;,\tag{S8} $

(8)

$ {\boldsymbol{L}}_{1}^{k+1}={\boldsymbol{L}}_{1}^{k}+\rho ({\boldsymbol{J}}_{1}^{k+1}-\nabla {{\boldsymbol{S}}}^{k\mathrm+1})\;,\tag{S9} $

(9)

$ {\boldsymbol{L}}_{2}^{k+1}={\boldsymbol{L}}_{2}^{k}+\rho ({\boldsymbol{J}}_{2}^{k+1}-{\boldsymbol{F}}{{\boldsymbol{S}}}^{k\mathrm+1}) \;.\tag{S10}$

(10)

S is our target variable, i.e., the to-be-reconstructed 3D structure. We need to iterate the algorithm based on the analytical solution until the error falls below the set threshold. Then the algorithm is considered converged, and S is the final reconstruction.

Figure S2 provides an illustration of the optimization. We chose a simulated finely textured saddle surface and extended by adding a background to represent the discontinuous regions in Fig. S2(h). As shown in the Fig. S2(a), the initial surface is all-zero, i.e., a flat surface with no stereo depth information. It’s expressed by a zero matrix and the normal gradients in both x and y directions are 0. Figure S2(b–d) show the depth dynamic estimation by iterations. Twenty iterations give a result visually approximate to the ground truth. Error evaluation in Fig. S2(e) also proves the algorithm has converged after 20 iterations, which means Fig. S2(d) is the final reconstruction.

Figure S2(f, g) plotted the depth value variation of Fig. S2(b–d) along the red and green directions in Fig. S2(a), where we can see the algorithm iterates rapidly at the beginning and then slowly approaches convergence. Final comparison is between the reconstruction Fig. S2(d) and ground truth Fig. S2(h), from which we each selected 400 random points and listed in Fig. S2(i). The reconstruction shows high reliability compared with the ground truth.

Section 4:　Point cloud registration of large scenes

Less

The point clouds corresponding to multiple frames are not in the same coordinate system, making point cloud registration necessary. Figure S3 give a set of matched feature points:

$ \begin{array}{cc}{\mathrm{P}}=\{{p}_{1},\cdots ,{p}_{n}\}, & {P}^{\prime}=\{p_{1}^{\prime},\cdots ,p_{n}^{\prime}\}\end{array}\;.\tag{S11} $

(11)

We need to find a Euclidean transformation R and t:

$ \forall i,{p}_{i}=Rp_{i}^{\prime}+t\;.\tag{S12} $

(12)

The transformation is completed by Iterative Closest Point (ICP) in Eq. (S5) in the main text. By continuously iterating, the algorithm gives the optimal solution at the minimum value.

In Fig. S4(a), the red and blue points represent feature points detected from two sets of point clouds. Figure S4(b) to S4(d) are the results of iterations by ICP. After several iterations, the red and blue point clouds gradually align from an initial misalignment, indicating accurate registration. Figure S4(e) depicts the ground truth of the pose transformation matrix we set, while Fig. S4(f) is the pose transformation matrix from the algorithm. The difference between the two matrices is minimal, with the Mean Squared Error of the error matrix 2.672e⁻¹⁰. It means the derived pose matrix is reliable for point cloud fusion, as illustrated by Fig. S5. Figure S5(a) is an image series when imaging natural large scene. Figure S5(b) shows the estimated camera poses and the final reconstruction is presented in the main text Fig. S5(c).

Section 5:　ADMM parameters determination

Less

The augmented Lagrange function involves three parameters: ξ, η, and ρ. We performed three groups of experiments, with only one parameter varied within each group. Figure S6(a) depicts the ground truth, while Fig. S6(b) and S6(c) show the simulated coarse depth and surface normals from polarization, respectively. By varying each parameter and calculating the SSIM between the simulated results and ground truth, we plotted the curves shown in Fig. S6(d–f). The curves indicate that the algorithm exhibits low sensitivity to weight parameter η but demonstrates considerable sensitivity to parameters ξ and ρ. Based on these tests, we set the parameter values as: ξ=0.2, η=0.1, ρ=0.2.

Appendix

Less

Year 2026 volume 9 Issue 2

PDF

Cite this Article

BibTeX

Article Info

doi: 10.29026/oea.2026.250267

Online Date：2026-10-19
Published：2026-02-12

Article Data

Affiliations

History

Affiliations

¹School of Optoelectronic Engineering, Xidian University, Xi’an 710071, China

²National Key Laboratory of Electromagnetic Space Security, Tianjin 300308, China

³State Key Laboratory of Electromechanical Integrated Manufacturing of High-Performance Electronic Equipments, Xidian University, Xi’an 710071, China

Corresponding:

F Liu, E-mail: feiliu@xidian.edu.cn

References

Share

https://castjournals.cast.org.cn/joweb/oea/EN/10.29026/oea.2026.250267

Share to

Scan QR to access full text

Cite this article

BibTeX

Citations

表12种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT

Fig.S1 Polarimetric and depth information of typical objects in natural scenes.

Fig.S2 Iteration results analyses. (a) Initial. (b) Result at the iteration k=2. (c) Result at the iteration k=5. (d) Result at the iteration k=20. (e) The error after iterations. (f) The slice of (a–d) at along the red line. (g) The slice of (a–d) along the green line. (h) Ground truth. (i) comparison.

Fig.S3 Matched feature points from two images.

Fig.S4 ICP algorithm process. (a) Two sets of unaligned feature point clouds. (b, c) Iterative process. (d) Aligned two sets of feature point clouds. (e) The ground truth of the transformation matrix. (f) Iteration result of transformation matrix.

Fig.S5 Schematic diagram of camera pose estimation. (a) Captured images. (b) The estimated camera poses corresponding to a series of captured images.

Fig.S6 Parameters ξ, η, and ρ determination. (a) Ground truth; (b) Rough depth; (c) Surface normal; (d–f) SSIM sensitivity analysis of parameter ξ, η and ρ.

Articles: Latest Articles; Most Read; Collections

Updates: Events; News; Multimedia

About: About Us

Contact

No. 86 Xueyuan South Road, Haidian District, Beijing

100081

010-62662699(Operation and Technical Support)

qkjq@cast.org.cn

Copyright © 2025 China Association for Science and Technology. All rights reserved. For all open access content, the relevant licensing terms apply.
Sponsored by the Office of the Leading Group for Cybersecurity and Informatization of CAST, and supported by Science and Technology Review Publishing House