偏光画像処理ライブラリを作った話

第63回 コンピュータビジョン勉強会＠関東「コンピュータビジョンでこんなプログラム作りました大LT大会5」での発表資料です。https://kantocv.connpass.com/event/344139/

偏光画像処理ライブラリPolanalyserについての紹介をします。
https://github.com/elerac/polanalyser

Transcript

1. コンピュ タビジョン勉強会＠関東 2025/02/24 前⽥ 涼汰 偏光画像処理ライブラリ
2. 偏光画像処理ライブラリ：Polanalyser • 偏光画像処理の開発に便利なPythonライブラリ • Polarization 偏光 Analyser 解析 Polanalyser • 2019年からGitHub上で公開．2025年2⽉現在165スタ ⭐． https://github.com/elerac/polanalyser pip install polanalyser 2
3. Polanalyserの主な機能 • 偏光デモザイキング • モノクロ/カラ 偏光センサに対応 • スト クスベクトルの解析 • スト クスベクトルの推定 • 偏光パラメ タの変換 AoLP，DoLPなど • ミュラ ⾏列の解析 • ミュラ ⾏列の推定 • 基本的な光学素⼦のミュラ ⾏列 • 偏光画像の可視化 • 偏光に特化したカラ マッピング AoLP，DoLP，CoP，ToP，など • 偏光の記号計算 • スト クスベクトル・ミュラ ⾏列の記号計算 3
4. 偏光 • 光には⾊ RGB だけでなく偏光の特性がある • 偏光：光の電場の振動⽅向の偏り⽅ • 直線偏光・円偏光・⾮偏光 • ⼈間の⽬は知覚できない • シャコ・ハチなどの⼀部の⽣き物は知覚できる • 偏光のふるまいは形状や材質等による物理的な関係がある • フレネルの式 が有名 • RGB画像だけでは得られない物理的制約を⽤いたシ ン認識・理解 4
5. コンピュ タビジョンにおける偏光 形状推定 拡散・鏡⾯反射分離 散乱除去 領域分割・材質推定 反射/透過光の分離 Kadambi , ICCV2015 Kalra , CVPR2020 Schechner , CVPR2004 インバ スレンダリング Baek , SIGGRAPH Asia 2018 Debevec , SIGGRAPH2000 5
6. 偏光デモザイキング 6
7. 偏光カメラ • ワンショットで偏光画像を撮影可能 • センサ 上のピクセル毎に異なる回転⾓度のマイクロ偏光板 • 通常のカラ カメラと同様にデモザイキング処理が必要 カラ 偏光センサ モノクロ偏光センサ 0° 45° 135° 90° Yamazaki , IEDM2016 センサ上の偏光素⼦ 7
8. 偏光カメラのデモザイキング RAW画像 デモザイキング後の画像 デモザイキング処理 近傍のピクセルから値を補間 8
9. Polanalyser：偏光カメラのデモザイキング • カラ カメラのデモザイキング OpenCV • 偏光カメラのデモザイキング Polanalyser import cv2 code = cv2.COLOR_BayerBG2BGR img_bgr = cv2.demosaicing(img_raw, code) import polanalyser as pa code = pa.COLOR_PolarMono imgs = pa.demosaicing(img_raw, code) code = pa.COLOR_PolarRGB imgs = pa.demosaicing(img_raw, code) カラ 偏光カメラの場合 モノクロ偏光カメラの場合 B G R 0 45 90 135 0 45 90 135 B G R 9
10. 偏光デモザイキングの結果 0 45 90 135 RAW画像 デモザイキング後の画像 10
11. スト クスベクトルの解析 11
12. 強度画像から偏光画像への変換 • 偏光板を通して得られる画像は単なる光強度画像 • 偏光のパラメ タに変換する必要 0 45 90 135 光強度画像と偏光板の回転⾓度 偏光画像 スト クスベクトル画像 12
13. スト クスベクトル • 偏光を表現するベクトル 𝐬 = 𝑠! , 𝑠" , 𝑠# , 𝑠$ % ∈ ℝ&×" • 直線偏光のみを扱う場合は 𝑠$ を除いた 𝑠! , 𝑠" , 𝑠# 光の絶対強度 直線偏光 0° vs. 90° 直線偏光 45° vs. 135° 円偏光 右回り vs. 左回り +1 -1 +1 -1 +1 -1 13
14. • 直線偏光板を回転させた場合： • 観測した複数の光強度からスト クスベクトルを推定 スト クスベクトルの計測 直線偏光のみ 偏光板の回転⾓度 𝜃 光強度 𝐼(𝜃) ０ 90 135 スト クスベクトル 𝐼 𝜃 = 0.5 (𝑠! + 𝑠" cos 2𝜃 + 𝑠# sin 2𝜃) 偏光板の回転⾓度 45 回転 𝜃 𝐬 光強度 𝐼 14
15. Polanalyser：スト クスベクトルの計測 直線偏光のみ • ⼊⼒：光強度画像と直線偏光板の回転⾓度のペア • 出⼒：スト クスベクトル画像 import polanalyser as pa # ⼊⼒の画像と⾓度 imgs = [img_000, img_045, img_090, img_135] angles = np.deg2rad([0, 45, 90, 135])) # ストークスベクトルの推定 img_stokes = pa.calcStokes(imgs, angles) 15
16. スト クスベクトルの計測結果 直線偏光のみ s! s" s# +1 -1 +1 -1 16
17. スト クスベクトルの計測 直線 円偏光 • より⼀般的にスト クスベクトルを計測するには？ • 直線偏光と円偏光を含むスト クスベクトル • 偏光光学素⼦の組み合わせによる計測 • 例：直線偏光板，波⻑板，ビ ムスプリッタ ，ミラ 17
18. 1/4波⻑板 直線偏光板による計測 • 1/4波⻑板と直線偏光板を組み合わせ，1/4波⻑板のみを回転 • 直線偏光と円偏光を同時に取得可能 • 観測した複数の光強度からスト クスベクトルを推定 直線偏光板 𝐋(0) 1/4波⻑板 𝐐(𝜃) 回転 𝐬 𝐋(0)𝐐(𝜃)𝐬 𝟎 光強度 全体のミュラ ⾏列 𝐋(0)𝐐(𝜃) 18
19. Polanalyser：スト クスベクトルの計測 直線 円偏光 • ⼊⼒：光強度画像とミュラ ⾏列のペア • 出⼒：スト クスベクトル画像 import polanalyser as pa # ⼊⼒の画像とミュラー⾏列 imgs = [img_000, img_045, img_090, img_135] thetas = np.rad2deg([0, 45, 90, 135])] muellers = [pa.polarizer(0) @ pa.qwp(theta) for theta in thetas] # ストークスベクトルの推定 img_stokes = pa.calcStokes(imgs, muellers) 19
20. スト クスベクトルの計測結果 直線 円偏光 s! s" s# s$ +1 -1 +1 -1 +1 -1 20
21. 偏光画像の可視化 21
22. スト クスベクトルの変換と可視化 • スト クスベクトル画像のみから直感的な偏光状態の理解が難しい • スト クスベクトルの変換とカラ マッピングによる可視化 スト クスベクトル画像 AoLP AoLP DoLP AoLP DoLP DoP DoLP DoCP ToP CoP 𝑠! 𝑠" 𝑠# 𝑠$ 22
23. Polanalyser: スト クスベクトルの変換 • スト クスベクトルから別の偏光表現⽅法に変換 • AoLP Angle of Linear Polarization : 振動⾯の⾓度 0, pi • DoLP Degree of Linear Polarization : 直線偏光度 0.0, 1.0 img_aolp = pa.cvtStokesToAoLP(img_stokes) img_dolp = pa.cvtStokesToDoLP(img_stokes) 0.5 tan"# 𝑠$ 𝑠# 𝑠# $ + 𝑠$ $ 𝑠% 23
24. Polanalyser: 偏光画像の可視化 AoLP, DoLP • 変換したAoLPやDoLPはfloat型 • 適切なカラ マッピングでuint8型に変換 img_aolp_vis = pa.applyColorToAoLP(img_aolp) img_dolp_vis = pa.applyColorToDoP(img_dolp) 0 90 180 0 1 AoLP DoLP 24
25. ミュラ ⾏列の解析 25
26. ミュラ ⾏列 • 偏光状態の 変化 を表現する4 4の⾏列 • ミュラ ⾏列とスト クスベクトルの積により偏光の変化を計算 𝐬!" 𝐬#$% 𝐌 𝐬#$% = 𝐌𝐬!" 𝐌 = 𝑚!! 𝑚!" 𝑚!# 𝑚!$ 𝑚"! 𝑚"" 𝑚"# 𝑚"$ 𝑚#! 𝑚#" 𝑚## 𝑚#$ 𝑚$! 𝑚$" 𝑚$# 𝑚$$ ∈ ℝ&×& 26
27. ミュラ ⾏列の例 • 直線偏光板 • 1/4波⻑板 𝐋 𝜃 = 1 2 1 cos 2𝜃 sin 2𝜃 0 cos 2𝜃 cos$ 2𝜃 sin 2𝜃 cos 2𝜃 0 sin 2𝜃 sin 2𝜃 cos 2𝜃 sin$ 2𝜃 0 0 0 0 0 𝐐 𝜃 = 1 0 0 0 0 cos$ 2𝜃 sin 2𝜃 cos 2𝜃 −sin 2𝜃 0 sin 2𝜃 cos 2𝜃 sin$ 2𝜃 cos 2𝜃 0 sin 2𝜃 − cos 2𝜃 0 𝜃：偏光素⼦の回転⾓度 27
28. Polanalyser： ミュラ ⾏列 • 基本的な光学素⼦のミュラ ⾏列をサポ ト pa.polarizer(0) # [[0.5 0.5 0. 0. ] # [0.5 0.5 0. 0. ] # [0. 0. 0. 0. ] # [0. 0. 0. 0. ]] pa.qwp(np.pi / 4) # [[ 1. 0. 0. 0.] # [ 0. 0. 0. -1.] # [ 0. 0. 1. 0.] # [ 0. 1. 0. 0.]] 直線偏光板 0度 1/4波⻑板 45度 pa.polarizer(np.pi / 4) # [[0.5 0. 0.5 0. ] # [0. 0. 0. 0. ] # [0.5 0. 0.5 0. ] # [0. 0. 0. 0. ]] pa.hwp(np.pi / 4) # [[ 1. 0. 0. 0.] # [ 0. -1. 0. 0.] # [ 0. 0. 1. 0.] # [ 0. 0. 0. -1.]] 直線偏光板 45度 1/2波⻑板 45度 28
29. エリプソメ タ • 光源とカメラ側の偏光状態を変調してミュラ ⾏列を計測 𝐌 = 𝑚%% 𝑚%# 𝑚%$ 𝑚%& 𝑚#% 𝑚## 𝑚#$ 𝑚#& 𝑚$% 𝑚$# 𝑚$$ 𝑚$& 𝑚&% 𝑚&# 𝑚&$ 𝑚&& 𝐀𝐌𝐏𝐬 ! 𝐬 = [1,0,0,0] 𝐏 = 𝐐 𝜃 𝐋(0) 𝐀 = 𝐋 0 𝐐(5𝜃) 𝐌 カメラ 無偏光光源 偏光⽣成機 1/4波⻑板 直線偏光板 偏光解析機 1/4波⻑板 直線偏光板 未知のミュラ ⾏列 光強度 29
30. Polanalyser：ミュラ ⾏列の計測 • ⼊⼒：光強度画像とミュラ ⾏列 偏光⽣成機・解析機 のペア • 出⼒：ミュラ ⾏列画像 import polanalyser as pa # ⼊⼒の画像とミュラー⾏列(偏光⽣成機・解析機) imgs = [img_00, img_01, img_02, img_03, --省略--, img_35] thetas = np.linspace(0, np.pi, num=len(imgs), endpoint=False) mm_psg = [pa.qwp(theta) @ pa.polarizer(0) for theta in thetas] mm_psa = [pa.polarizer(0) @ pa.qwp(5 * theta) for theta in thetas] # ミュラー⾏列の推定 img_mueller = pa.calcMueller(imgs, mm_psg, mm_psa) 30
31. ミュラ ⾏列の計測結果 31
32. 偏光の記号計算 32
33. Polanalyser: 偏光の記号計算 • 偏光のモデル式を記号的に計算したいことがある • しかし複数の偏光素⼦を組み合わせると複雑に • SymPyと組み合わせた偏光の記号計算のモジュ ルを提供 from sympy import Symbol, simplify import polanalyser.sympy as pas theta = Symbol("theta", real=True) f = (pas.polarizer(theta) @ pas.polarizer(0))[0, 0] f = simplify(f) マリュスの法則の計算例：2枚の直線偏光板を通した光強度 𝑓 𝜃 = 𝐋 𝜃 𝐋 0 !,! = 0.5 cos# 𝜃 0.5*cos(theta)**2 33
34. 𝜕 𝜕𝜃 log( 𝐋 0 𝐐 5𝜃 𝐌𝐐 𝜃 𝐋 0 𝐬 & ) = Polanalyser: 偏光の記号計算の活⽤事例 Ryota Maeda, Yunseong Moon, Seung Hwan Baek, Event Ellipsometer: Event based Mueller Matrix Video Imaging, arXiv 2024 イベントカメラを⽤いた⾼速なミュラ ⾏列の計測 1 1 0 0 1 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 cos! 10𝜃 sin 10𝜃 cos 10𝜃 −sin 10𝜃 0 sin 10𝜃 cos 10𝜃 sin! 10𝜃 cos 10𝜃 0 sin 10𝜃 − cos 10𝜃 0 𝑚"" 𝑚"# 𝑚"! 𝑚"$ 𝑚#" 𝑚## 𝑚#! 𝑚#$ 𝑚!" 𝑚!# 𝑚!! 𝑚!$ 𝑚$" 𝑚$# 𝑚$! 𝑚$$ 1 0 0 0 0 cos! 2𝜃 sin 2𝜃 cos 2𝜃 −sin 2𝜃 0 sin 2𝜃 cos 2𝜃 sin! 2𝜃 cos 2𝜃 0 sin 2𝜃 − cos 2𝜃 0 1 1 0 0 1 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 イベントとミュラ ⾏列に関しての式： ＊実際の式より簡略化しています ? 𝜕 𝜕𝜃 log( 𝐋 0 𝐐 5𝜃 𝐌𝐐 𝜃 𝐋 0 𝐬 & ) = 34
35. Polanalyser: 偏光の記号計算の活⽤事例 from sympy import Symbol, Matrix, ln, diff import polanalyser.sympy as pas theta = Symbol("theta", real=True) M = pas.mueller() lqmqls = (pas.polarizer(0) @ pas.qwp(5 * theta) @ M @ pas.qwp(theta) @ pas.polarizer(0) @ Matrix([1, 0, 0, 0]))[0] diff_ln = diff(ln(lqmqls), theta) 𝜕 𝜕𝜃 log( 𝐋 0 𝐐 5𝜃 𝐌𝐐 𝜃 𝐋 0 𝐬 ! ) = # (-5.0*m10*sin(10*theta)*cos(10*theta) - 2.5*m20*sin(10*theta)**2 + 2.5*m20*cos(10*theta)**2 - 2.5*m30*cos(10*theta) - 2.0*(0.5*m01 + 0.5*m11*cos(10*theta)**2 + 0.5*m21*sin(10*theta)*cos(10*theta) - 0.5*m31*sin(10*theta))*sin(2*theta)*cos(2*theta) - 1.0*(0.5*m02 + 0.5*m12*cos(10*theta)**2 + 0.5*m22*sin(10*theta)*cos(10*theta) - 0.5*m32*sin(10*theta))*sin(2*theta)**2 + 1.0*(0.5*m02 + 0.5*m12*cos(10*theta)**2 + 0.5*m22*sin(10*theta)*cos(10*theta) - 0.5*m32*sin(10*theta))*cos(2*theta)**2 + 1.0*(0.5*m03 + 0.5*m13*cos(10*theta)**2 + 0.5*m23*sin(10*theta)*cos(10*theta) - 0.5*m33*sin(10*theta))*cos(2*theta) + 0.5*(-10.0*m11*sin(10*theta)*cos(10*theta) - 5.0*m21*sin(10*theta)**2 + 5.0*m21*cos(10*theta)**2 - 5.0*m31*cos(10*theta))*cos(2*theta)**2 + 0.5*(-10.0*m12*sin(10*theta)*cos(10*theta) - 5.0*m22*sin(10*theta)**2 + 5.0*m22*cos(10*theta)**2 - 5.0*m32*cos(10*theta))*sin(2*theta)*cos(2*theta) + 0.5*(-10.0*m13*sin(10*theta)*cos(10*theta) - 5.0*m23*sin(10*theta)**2 + 5.0*m23*cos(10*theta)**2 - 5.0*m33*cos(10*theta))*sin(2*theta))/(0.25*m00 + 0.25*m10*cos(10*theta)**2 + 0.25*m20*sin(10*theta)*cos(10*theta) - 0.25*m30*sin(10*theta) + 0.5*(0.5*m01 + 0.5*m11*cos(10*theta)**2 + 0.5*m21*sin(10*theta)*cos(10*theta) - 0.5*m31*sin(10*theta))*cos(2*theta)**2 + 0.5*(0.5*m02 + 0.5*m12*cos(10*theta)**2 + 0.5*m22*sin(10*theta)*cos(10*theta) - 0.5*m32*sin(10*theta))*sin(2*theta)*cos(2*theta) + 0.5*(0.5*m03 + 0.5*m13*cos(10*theta)**2 + 0.5*m23*sin(10*theta)*cos(10*theta) - 0.5*m33*sin(10*theta))*sin(2*theta)) 35
36. GitHubで公開したメリット 1. コ ド管理・品質向上 • ブラウザからどこでもコ ドを確認 • pipで即インスト ル • リ ダブルコ ド・ドキュメントの整備 ⾃分にとっても読みやすい 2. コミュニティへの貢献 • 他の研究でも活⽤されている 偏光CVの発展 • Polanalyserの使⽤が確認できた研究 GitHub上で検索 • Dave , PANDORA: Polarization Aided Neural Decomposition of Radiance, ECCV2022 • Tiwary , ORCa: Objects as Radiance Field Cameras, CVPR2023 • Tang , Learning to Remove Wrinkled Transparent Film with Polarized Prior, CVPR2024 • Li , GNeRP: Gaussian guided Neural Reconstruction of Reflective Objects with Noisy Polarization Priors, ICLR2024 36
37. まとめ • 偏光画像処理ライブラリ：Polanalyser • 偏光デモザイキング • スト クスベクトルの解析 • ミュラ ⾏列の解析 • 偏光画像の可視化 • 偏光の記号計算 • 今後も機能追加予定 • 偏光カメラで撮影した画像も公開中 • すぐに試せます！ https://github.com/elerac/polanalyser 37