A13 - Photometric Stereo

August 7, 2008
We have 4 synthetic images, namely I1, I2, I3 and I4, of a sphere illuminated from 4 different directions.
The illumination matrix is a composite of 4 vectors given by

V1 = [0.085832 0.17365 0.98106]
V2 = [0.085832 -0.17365 0.98106]
V3 = [0.17365 0 0.98481]
V4 = [0.16318 -0.34202 0.92542].

The illumination is from an infinitely far away source. The relationship between the reflected intensity I, the illumination vector L and the object surface normal N is given by

I = L*N (1).

Since we know I and L, we could get the surface normals via

g = [(inv(V'V))V'] I (2).

We normalize g by dividing by its magnitude to obtain the normalized surface normal

n^ = g/|g| (3).

This components of the normal vector can be transformed to the partial derivatives by using the relation

df/dx = -nx/nz and df/dy = -ny/nz.

The derivatives are then integrated via line integration to obtain the surface profile of the hemisphere (See Figure )2. The drawback of using line integration is that it produces small spikes in the reconstruction.

code:
loadmatfile('C:\Documents and Settings\AP186user20\Desktop\august7_2008\photos.mat',['I1','I2','I3','I4']);
subplot(221)
imshow(I1,[]);
subplot(222)
imshow(I2,[]);
subplot(223)
imshow(I3,[]);
subplot(224)
imshow(I4,[]);
[m,k] = size(I1)
I = [I1(:)'; I2(:)'; I3(:)';I4(:)'];
V = [0.085832 0.17365 0.98106;
0.085832 -0.17365 0.98106
0.17365 0 0.98481
0.16318 -0.34202 0.92542];
MoorePennrose = (inv(V'*V))*V';
g = MoorePennrose*I;
[M,K] = size(g);
for i = 1:M;
for j = 1:K;
absolutevalue = sqrt( g(1,j)*g(1,j) + g(2,j)*g(2,j) + g(3,j)*g(3,j) );
n(i,j) = g(i,j)./(absolutevalue + 1.0e-015);
end
end
//derivative along the x-axis
p = -n(1,:)./(n(3,:)+1.0e-20);
//derivative along the y-axis
q = -n(2,:)./(n(3,:)+1.0e-20);

Imagexderivative = matrix(p,m,k);
Imageyderivative = matrix(q,m,k);
x=cumsum(Imagexderivative,2);
y= cumsum(Imageyderivative,1);
Image = x+y;
mesh(Image);
Reference
[1]. Dr. Soriano's lecture notes

Self-grade: 10/10