Sample program of shape-from-polarization. The target object is assumed to have diffuse reflection only without specular reflection. The object surface should be smooth since the surface roughness coefficient is not considered in this program. Index-of-refraction should be known, and should be written in the variable "ior" in the sourcecode. This program assumes that y-axis of surface normal is positive. Input data is degree-of-polarization "dop.npy" and angle-of-polarization "phase.npy". Angle-of-polarization is represented in radian, and is the angle of linear polarizer when the maximum brightness is observed. Open Google Colab, create new file, copy&paste the text below. Run, push upload button, select "dop.npy" and "phase.npy". Output surface normal is downloaded. The output is the numpy file with 3D vector in each pixel. The y-axis of surface normal is the upper direction.
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')
