偏光から法線を計算するサンプルプログラムです．鏡面反射のない，拡散反射だけがある物体の場合のプログラムです．表面粗さ係数を考慮していないプログラムのため，表面が滑らかな物体を対象としています．屈折率は既知とし，ソースコード中の変数iorに屈折率をあらかじめ記入しておいてください．このサンプルソースコードでは，法線のy軸は全てプラス方向しか存在しないことを仮定しています．入力データは，偏光度dop.npyと偏光角phase.npyです．偏光角は，偏光板を回したときに輝度が最大となる角度とし，ラジアンで表現されているものとします．Google Colabを開いて新規作成し，以下のソースコードのテキスト情報をコピペしてください．実行したら，ファイルアップロードのボタンが出てくるので，ボタンを押して"dop.npy"と"phase.npy"を指定してください．計算結果の法線がダウンロードされます．各画素に３次元ベクトルが格納されたnumpyファイルとして出力されます．法線のy軸は上向きとします．
dop.npy
phase.npy
dop.bmp
phase.bmp
inside.bmp
from google.colab import files
import numpy as np
import cv2
from scipy.optimize import minimize_scalar
def zenith2dop(theta1):
ior=1.5
theta2=np.arcsin(np.sin(theta1)/ior)
cosval=np.cos(theta1-theta2)**2
dopval=(1-cosval)/(1+cosval)
return dopval
def costfunction(zenithval,dopval):
return (zenith2dop(zenithval)-dopval)**2
def dop2zenith(dopval):
result=minimize_scalar(
costfunction,
bounds=(0,np.pi/2),
args=(dopval),
method='bounded',
)
return result.x
files.upload()
phasedata=np.load('phase.npy')
dopdata=np.load('dop.npy')
rows,cols=phasedata.shape
normaldata=np.zeros((rows,cols,3),dtype=np.float64)
for y in range(0,rows):
if y%10==9: print('Calculating...',y+1,'/',rows)
for x in range(0,cols):
phaseval=phasedata[y,x]
azimuthval=phaseval
dopval=dopdata[y,x]
zenithval=dop2zenith(dopval)
nx=np.sin(zenithval)*np.cos(azimuthval)
ny=np.sin(zenithval)*np.sin(azimuthval)
nz=np.cos(zenithval)
normaldata[y,x,0]=nx
normaldata[y,x,1]=ny
normaldata[y,x,2]=nz
np.save('normal.npy',normaldata)
files.download('normal.npy')
