Commit 89fa886c authored by Sushil BHATTACHARJEE's avatar Sushil BHATTACHARJEE
Browse files

Merge branch 'cleanup' into 'master'

speed-up code in some critical places

See merge request !2
parents 25391245 de7369fd
Pipeline #10498 failed with stages
in 7 minutes and 42 seconds
...@@ -727,36 +727,26 @@ def edge_thinning(bx, by, thinning=1): ...@@ -727,36 +727,26 @@ def edge_thinning(bx, by, thinning=1):
# np.spacing(1) is the same as eps in matlab. # np.spacing(1) is the same as eps in matlab.
myEps = np.spacing(1) * 100.0 myEps = np.spacing(1) * 100.0
# compute the edge-map a la Matlab # compute the edge-map a la Matlab
for r in range(m):
for c in range(n): if not thinning:
if thinning: e = b > cutoff
if r < 0 or r > (m - 1) or (c - 1) < 0: else:
b1 = True b1 = np.ones_like(b, dtype=bool)
else: b2 = np.ones_like(b, dtype=bool)
b1 = (b[r, c - 1] < b[r, c]) b3 = np.ones_like(b, dtype=bool)
b4 = np.ones_like(b, dtype=bool)
if(r < 0) or r > (m - 1) or (c + 1) > (n - 1):
b2 = True c1 = b > cutoff
else:
b2 = (b[r, c] > b[r, c + 1]) b1[:, 1:] = (np.roll(b, 1, axis=1) < b)[:, 1:]
b2[:, :-1] = (np.roll(b, -1, axis=1) < b)[:, :-1]
if (c < 0) or c > (n - 1) or (r - 1) < 0: c2 = (bx >= (by - myEps)) & b1 & b2
b3 = True
else: b3[1:, :] = (np.roll(b, 1, axis=0) < b)[1:, :]
b3 = (b[r, c] > b[r - 1, c]) b4[:-1, 1:] = (np.roll(b, -1, axis=0) < b)[:-1, 1:]
c3 = (by >= (bx - myEps)) & b3 & b4
if(c < 1) or c > (n - 1) or (r + 1) > (m - 1):
b4 = True e = c1 & (c2 | c3)
else:
b4 = (b[r, c] > b[r + 1, c])
c1 = (b[r, c] > cutoff)
c2 = ((bx[r, c] >= (by[r, c] - myEps)) & b1 & b2)
c3 = ((by[r, c] >= (bx[r, c] - myEps)) & b3 & b4)
e[r, c] = c1 & (c2 | c3)
else:
e[r, c] = (b[r, c] > cutoff)
return e return e
...@@ -813,33 +803,17 @@ def angle_similarity(refImage, testImage, diffImage): ...@@ -813,33 +803,17 @@ def angle_similarity(refImage, testImage, diffImage):
refNorm = np.linalg.norm(refImage, axis=1) refNorm = np.linalg.norm(refImage, axis=1)
testNorm = np.linalg.norm(testImage, axis=1) testNorm = np.linalg.norm(testImage, axis=1)
thetaVec = np.zeros([refImage.shape[0], 1])
diffNorm = np.linalg.norm(diffImage, axis=1) diffNorm = np.linalg.norm(diffImage, axis=1)
magnitVec = diffNorm / 255.0 # Galbally divides by sqrt(255**2) magnitVec = diffNorm / 255.0 # Galbally divides by sqrt(255**2)
magnitVec = np.reshape(magnitVec, (refImage.shape[0], 1)) magnitVec = np.reshape(magnitVec, (refImage.shape[0], 1))
for i in range(refImage.shape[0]): # np.einsum('ij,ij->i',a,b) is equivalent to np.diag(np.dot(a,b.T))
refR = refImage[i, :] cosTheta = np.einsum('ij,ij->i', refImage, testImage) / \
testR = testImage[i, :] (refNorm * testNorm)
cosTheta = np.dot(refR, testR) / (refNorm[i] * testNorm[i]) cosTheta[cosTheta < -1.0] = -1.0
if(cosTheta < -1.0): cosTheta[cosTheta > 1.0] = 1.0
cosTheta = -1.0 cosTheta = np.nan_to_num(cosTheta)
if(cosTheta > 1.0): thetaVec = np.arccos(cosTheta).reshape((refImage.shape[0], 1))
cosTheta = 1.0
theta = np.arccos(cosTheta)
thetaVec[i] = theta
# the following (commented out) code should be more efficient than the for-
# loop above, but seems to be buggy.
# rowDP= np.diag(np.dot(refImage,testImage.T))
# normRefrows= np.linalg.norm(refImage,axis=1)
# normTestrows= np.linalg.norm(testImage,axis=1)
# cosThetaVec = rowDP/(normRefrows * normTestrows)
# cosThetaVec = np.nan_to_num(cosThetaVec) #nan occurs when one of the
# # norms is 0, ie. vector is all 0s. In that case set cosTheta to 0
# thetaVec = np.arccos(cosThetaVec)
tmp2 = thetaVec * 2.0 / np.pi tmp2 = thetaVec * 2.0 / np.pi
......
...@@ -399,13 +399,17 @@ def rgbhist(image, maxval, nBins, normType=0): ...@@ -399,13 +399,17 @@ def rgbhist(image, maxval, nBins, normType=0):
im = image.reshape(3, numPix).copy() im = image.reshape(3, numPix).copy()
im = im.T im = im.T
for i in range(0, numPix): # for i=1:size(I,1)*size(I,2) p = np.floor(im.astype(float) / decimator).astype(np.uint32)
p = (im[i, :]).astype(float) # p = double(im(i,:)); # in future versions of numpy (1.13 and above) you can replace this with:
# p = floor(p/(maxval/nBins))+1; # unique_p, count = np.unique(p, return_counts=True, axis=0)
p = (np.floor(p / decimator)).astype(np.uint32) # the following lines were taken from: https://stackoverflow.com/a/16973510
# H(p(1),p(2),p(3)) = H(p(1),p(2),p(3)) + 1; p2 = np.ascontiguousarray(p).view(
H[p[0], p[1], p[2]] += 1 np.dtype((np.void, p.dtype.itemsize * p.shape[1])))
# end unique_p, count = np.unique(p2, return_counts=True)
unique_p = unique_p.view(p.dtype).reshape(-1, p.shape[1])
# till here
H[unique_p[:, 0], unique_p[:, 1], unique_p[:, 2]] = count
H = H.ravel() # H = H(:); H = H.ravel() # H = H(:);
# Un-Normalized histogram # Un-Normalized histogram
......
...@@ -130,10 +130,7 @@ def specular_free_image(src, srcPixelStatus=None): ...@@ -130,10 +130,7 @@ def specular_free_image(src, srcPixelStatus=None):
":: srcPixelStatus should be None, or should match the " \ ":: srcPixelStatus should be None, or should match the " \
"pixel-layout of src." "pixel-layout of src."
for y in range(src.shape[1]): srcPixelStatus[np.logical_and.reduce(src < camDark, 0)] = G_CAMERA_DARK
for x in range(src.shape[2]):
if np.all(src[:, y, x] < camDark):
srcPixelStatus[y, x] = G_CAMERA_DARK
srcMaxChroma, srcMax, srcTotal = max_chroma(src) srcMaxChroma, srcMax, srcTotal = max_chroma(src)
...@@ -181,10 +178,9 @@ def iteration(src, srcPixelStatus, sfi, epsilon, verboseFlag=True): ...@@ -181,10 +178,9 @@ def iteration(src, srcPixelStatus, sfi, epsilon, verboseFlag=True):
count, srcPixelStatus = init_labels(src, sfi, epsilon, srcPixelStatus) count, srcPixelStatus = init_labels(src, sfi, epsilon, srcPixelStatus)
while(1): while(1):
maxY, maxX = src.shape[-2:]
srcMaxChroma, srcMax, srcTotal = max_chroma(src) srcMaxChroma, srcMax, srcTotal = max_chroma(src)
srcChroma, srcTotal = rgb_chroma(src, srcTotal) srcChroma, srcTotal = rgb_chroma(src, srcTotal)
maxY = src.shape[1]
maxX = src.shape[2]
red_chroma_diff_x = np.diff(srcChroma[0, :, :], axis=1) red_chroma_diff_x = np.diff(srcChroma[0, :, :], axis=1)
red_chroma_diff_y = np.diff(srcChroma[0, :, :], axis=0) red_chroma_diff_y = np.diff(srcChroma[0, :, :], axis=0)
...@@ -192,119 +188,248 @@ def iteration(src, srcPixelStatus, sfi, epsilon, verboseFlag=True): ...@@ -192,119 +188,248 @@ def iteration(src, srcPixelStatus, sfi, epsilon, verboseFlag=True):
grn_chroma_diff_x = np.diff(srcChroma[1, :, :], axis=1) grn_chroma_diff_x = np.diff(srcChroma[1, :, :], axis=1)
grn_chroma_diff_y = np.diff(srcChroma[1, :, :], axis=0) grn_chroma_diff_y = np.diff(srcChroma[1, :, :], axis=0)
# if verboseFlag: # a non-complete attempt to remove the for loop below
# print('.') # # find non G_CAMERA_DARK pixels
# non_g_camera_dark = (srcPixelStatus != G_CAMERA_DARK)
# non_g_camera_dark[-1, :] = False
# non_g_camera_dark[:, -1] = False
# g_specularx = (non_g_camera_dark) & \
# (srcPixelStatus == G_SPECULARX)
# g_speculary = (non_g_camera_dark) & \
# (srcPixelStatus == G_SPECULARY)
# '''for specular y'''
# # if it is a boundary in x
# boundary_x = np.zeros_like(g_specularx)
# boundary_x[:-1, :-1] = g_specularx[:-1, :-1] & \
# (np.fabs(red_chroma_diff_x[:-1, :]) > thR) & \
# (np.fabs(grn_chroma_diff_x[:-1, :]) > thG)
# srcPixelStatus[boundary_x] = G_BOUNDARY
# # if it is noise in x
# srcMaxChromaDiffX = np.diff(srcMaxChroma, axis=1)
# srcMaxChromaAbsDiffX = np.fabs(srcMaxChromaDiffX)
# noise_mask_x = g_specularx & (~boundary_x)
# noise_mask_x[:-1, :-1] = noise_mask_x[:-1, :-1] & \
# (srcMaxChromaAbsDiffX[:-1, :] < 0.01)
# srcPixelStatus[noise_mask_x] = G_NOISE
# # if it is a boundary in x
# boundary_y = np.zeros_like(g_speculary)
# boundary_y[:-1, :-1] = g_speculary[:-1, :-1] & \
# (np.fabs(red_chroma_diff_y[:, :-1]) > thR) & \
# (np.fabs(grn_chroma_diff_y[:, :-1]) > thG)
# srcPixelStatus[boundary_y] = G_BOUNDARY
# # if it is noise in x
# srcMaxChromaDiffY = np.diff(srcMaxChroma, axis=0)
# srcMaxChromaAbsDiffY = np.fabs(srcMaxChromaDiffY)
# noise_mask_y = g_speculary & (~boundary_y)
# noise_mask_y[:-1, :-1] = noise_mask_y[:-1, :-1] & \
# (srcMaxChromaAbsDiffY[:, :-1] < 0.01)
# srcPixelStatus[noise_mask_y] = G_NOISE
# # if it is not noise or boundary
# boundary_x = np.zeros_like(g_specularx)
# boundary_x[:-1, :-1] = g_specularx[:-1, :-1] & \
# (np.fabs(red_chroma_diff_x[:-1, :]) > thR) & \
# (np.fabs(grn_chroma_diff_x[:-1, :]) > thG)
# srcMaxChromaDiffX = np.diff(srcMaxChroma, axis=1)
# srcMaxChromaAbsDiffX = np.fabs(srcMaxChromaDiffX)
# noise_mask_x = g_specularx & (~boundary_x)
# noise_mask_x[:-1, :-1] = noise_mask_x[:-1, :-1] & \
# (srcMaxChromaAbsDiffX[:-1, :] < 0.01)
# non_noise_mask_x = g_specularx & (~boundary_x) & (~noise_mask_x)
# srcPixelStatus[non_noise_mask_x] = G_DIFFUSE
# non_noise_mask_next_x = np.roll(non_noise_mask_x, 1, axis=1)
# srcPixelStatus[non_noise_mask_next_x] = G_DIFFUSE
# boundary_y = np.zeros_like(g_speculary)
# boundary_y[:-1, :-1] = g_speculary[:-1, :-1] & \
# (np.fabs(red_chroma_diff_y[:, :-1]) > thR) & \
# (np.fabs(grn_chroma_diff_y[:, :-1]) > thG)
# srcMaxChromaDiffY = np.diff(srcMaxChroma, axis=0)
# srcMaxChromaAbsDiffY = np.fabs(srcMaxChromaDiffY)
# noise_mask_y = g_speculary & (~boundary_y)
# noise_mask_y[:-1, :-1] = noise_mask_y[:-1, :-1] & \
# (srcMaxChromaAbsDiffY[:, :-1] < 0.01)
# non_noise_mask_y = g_speculary & (~boundary_y) & (~noise_mask_y)
# srcPixelStatus[non_noise_mask_y] = G_DIFFUSE
# non_noise_mask_next_y = np.roll(non_noise_mask_y, 1, axis=0)
# srcPixelStatus[non_noise_mask_next_y] = G_DIFFUSE
# # if it is not boundary or noise
# pos_chroma_maskx = (srcMaxChromaDiffX[:-1, :] > 0)
# reduce_mask_pos_x = non_noise_mask_x.copy()
# reduce_mask_pos_x[:-1, :-1] = reduce_mask_pos_x[:-1, :-1] & \
# pos_chroma_maskx
# reduce_mask_pos_next_x = np.roll(reduce_mask_pos_x, 1, axis=1)
# iro = np.moveaxis(np.moveaxis(src, 0, -1)[reduce_mask_pos_x], -1, 0)
# refMaxChroma = srcMaxChroma[reduce_mask_pos_next_x]
# iroTotal = srcTotal[reduce_mask_pos_x]
# iroMax = srcMax[reduce_mask_pos_x]
# iroMaxChroma = srcMaxChroma[reduce_mask_pos_x]
# pStat, iroRes = specular_to_diffuse_vectorized(
# iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# iro[:, pStat != G_NOISE] = iroRes[:, pStat != G_NOISE]
# pos_chroma_maskx = ~pos_chroma_maskx
# reduce_mask_pos_x = non_noise_mask_x.copy()
# reduce_mask_pos_x[:-1, :-1] = reduce_mask_pos_x[:-1, :-1] & \
# pos_chroma_maskx
# reduce_mask_pos_next_x = np.roll(reduce_mask_pos_x, 1, axis=1)
# iro = np.moveaxis(np.moveaxis(src, 0, -1)[reduce_mask_pos_next_x],
# -1, 0)
# refMaxChroma = srcMaxChroma[reduce_mask_pos_x]
# iroTotal = srcTotal[reduce_mask_pos_next_x]
# iroMax = srcMax[reduce_mask_pos_next_x]
# iroMaxChroma = srcMaxChroma[reduce_mask_pos_next_x]
# pStat, iroRes = specular_to_diffuse_vectorized(
# iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# iro[:, pStat != G_NOISE] = iroRes[:, pStat != G_NOISE]
# # if it is not boundary or noise
# pos_chroma_masky = (srcMaxChromaDiffY[:, :-1] > 0)
# reduce_mask_pos_y = non_noise_mask_y.copy()
# reduce_mask_pos_y[:-1, :-1] = reduce_mask_pos_y[:-1, :-1] & \
# pos_chroma_masky
# reduce_mask_pos_next_y = np.roll(reduce_mask_pos_y, 1, axis=1)
# iro = np.moveaxis(np.moveaxis(src, 0, -1)[reduce_mask_pos_y], -1, 0)
# refMaxChroma = srcMaxChroma[reduce_mask_pos_next_y]
# iroTotal = srcTotal[reduce_mask_pos_y]
# iroMax = srcMax[reduce_mask_pos_y]
# iroMaxChroma = srcMaxChroma[reduce_mask_pos_y]
# pStat, iroRes = specular_to_diffuse_vectorized(
# iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# iro[:, pStat != G_NOISE] = iroRes[:, pStat != G_NOISE]
# pos_chroma_masky = ~pos_chroma_masky
# reduce_mask_pos_y = non_noise_mask_y.copy()
# reduce_mask_pos_y[:-1, :-1] = reduce_mask_pos_y[:-1, :-1] & \
# pos_chroma_masky
# reduce_mask_pos_next_y = np.roll(reduce_mask_pos_y, 1, axis=1)
# iro = np.moveaxis(np.moveaxis(src, 0, -1)[reduce_mask_pos_next_y],
# -1, 0)
# refMaxChroma = srcMaxChroma[reduce_mask_pos_y]
# iroTotal = srcTotal[reduce_mask_pos_next_y]
# iroMax = srcMax[reduce_mask_pos_next_y]
# iroMaxChroma = srcMaxChroma[reduce_mask_pos_next_y]
# pStat, iroRes = specular_to_diffuse_vectorized(
# iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# iro[:, pStat != G_NOISE] = iroRes[:, pStat != G_NOISE]
for y in range(maxY - 1): for y in range(maxY - 1):
for x in range(maxX - 1): for x in range(maxX - 1):
current_pixel = srcPixelStatus[y, x] current_pixel = srcPixelStatus[y, x]
if(current_pixel != G_CAMERA_DARK): if(current_pixel == G_CAMERA_DARK):
continue
drx = red_chroma_diff_x[y, x]
dgx = grn_chroma_diff_x[y, x] drx = red_chroma_diff_x[y, x]
dry = red_chroma_diff_y[y, x] dgx = grn_chroma_diff_x[y, x]
dgy = grn_chroma_diff_y[y, x] dry = red_chroma_diff_y[y, x]
dgy = grn_chroma_diff_y[y, x]
if(current_pixel == G_SPECULARX):
# if it is a boundary in the x direction if(current_pixel == G_SPECULARX):
if((abs(drx) > thR) and (abs(dgx) > thG)): # if it is a boundary in the x direction
# pixel right if (abs(drx) > thR) and (abs(dgx) > thG):
srcPixelStatus[y, x] = G_BOUNDARY # pixel right
continue srcPixelStatus[y, x] = G_BOUNDARY
# if it is a noise continue
elif (abs(srcMaxChroma[y, x] - srcMaxChroma[y, x + 1]) < 0.01): # if it is a noise
srcPixelStatus[y, x] = G_NOISE if abs(srcMaxChroma[y, x] - srcMaxChroma[y, x + 1]) < 0.01:
continue srcPixelStatus[y, x] = G_NOISE
else: # reduce the specularity at x direction continue
if (srcMaxChroma[y, x] < srcMaxChroma[y, x + 1]): # reduce the specularity at x direction
iro = src[:, y, x] if srcMaxChroma[y, x] < srcMaxChroma[y, x + 1]:
pStat = current_pixel iro = src[:, y, x]
refMaxChroma = srcMaxChroma[y, x + 1] # pStat = current_pixel
iroTotal = srcTotal[y, x] refMaxChroma = srcMaxChroma[y, x + 1]
iroMax = srcMax[y, x] iroTotal = srcTotal[y, x]
iroMaxChroma = srcMaxChroma[y, x] iroMax = srcMax[y, x]
pStat, iroRes = specular_to_diffuse( iroMaxChroma = srcMaxChroma[y, x]
iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma) pStat, iroRes = specular_to_diffuse(
iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# update pixelStatus only if
# zSpecular2Diffuse() returned pStat as # update pixelStatus only if
# G_NOISE. # zSpecular2Diffuse() returned pStat as
if pStat == G_NOISE: # G_NOISE.
srcPixelStatus[y, x] = pStat if pStat != G_NOISE:
else: src[:, y, x] = iroRes
src[:, y, x] = iroRes # else:
# srcPixelStatus[y, x] = pStat
else:
iro = src[:, y, x + 1] else:
pStat = srcPixelStatus[y, x + 1] iro = src[:, y, x + 1]
refMaxChroma = srcMaxChroma[y, x] # pStat = srcPixelStatus[y, x + 1]
iroTotal = srcTotal[y, x + 1] refMaxChroma = srcMaxChroma[y, x]
iroMax = srcMax[y, x + 1] iroTotal = srcTotal[y, x + 1]
iroMaxChroma = srcMaxChroma[y, x + 1] iroMax = srcMax[y, x + 1]
pStat, iroRes = specular_to_diffuse( iroMaxChroma = srcMaxChroma[y, x + 1]
iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma) pStat, iroRes = specular_to_diffuse(
iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# update pixelStatus only if
# zSpecular2Diffuse() returned pStat as # update pixelStatus only if
# G_NOISE. # zSpecular2Diffuse() returned pStat as
if pStat == G_NOISE: # G_NOISE.
srcPixelStatus[y, x + 1] = pStat if pStat != G_NOISE:
else: src[:, y, x + 1] = iroRes
src[:, y, x + 1] = iroRes # else:
# srcPixelStatus[y, x + 1] = pStat
srcPixelStatus[y, x] = G_DIFFUSE
srcPixelStatus[y, x + 1] = G_DIFFUSE srcPixelStatus[y, x] = G_DIFFUSE
srcPixelStatus[y, x + 1] = G_DIFFUSE
if(current_pixel == G_SPECULARY):
# if it is a boundary in the y direction if(current_pixel == G_SPECULARY):
if((abs(dry) > thR) and (abs(dgy) > thG)): # if it is a boundary in the y direction
# pixel right if (abs(dry) > thR) and (abs(dgy) > thG):
srcPixelStatus[y, x] = G_BOUNDARY # pixel right
continue srcPixelStatus[y, x] = G_BOUNDARY
# if it is a noise continue
elif (abs(srcMaxChroma[y, x] - srcMaxChroma[y + 1, x]) < 0.01): # if it is a noise
srcPixelStatus[y, x] = G_NOISE if abs(srcMaxChroma[y, x] - srcMaxChroma[y + 1, x]) < 0.01:
continue srcPixelStatus[y, x] = G_NOISE
else: # reduce the specularity in y direction continue
if(srcMaxChroma[y, x] < srcMaxChroma[y + 1, x]): # reduce the specularity in y direction
iro = src[:, y, x] if srcMaxChroma[y, x] < srcMaxChroma[y + 1, x]:
pStat = current_pixel iro = src[:, y, x]
refMaxChroma = srcMaxChroma[y + 1, x] pStat = current_pixel
iroTotal = srcTotal[y, x] refMaxChroma = srcMaxChroma[y + 1, x]
iroMax = srcMax[y, x] iroTotal = srcTotal[y, x]
iroMaxChroma = srcMaxChroma[y, x] iroMax = srcMax[y, x]
pStat, iroRes = specular_to_diffuse( iroMaxChroma = srcMaxChroma[y, x]
iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma) pStat, iroRes = specular_to_diffuse(
# C call: iro, iroMax, iroMaxChroma, iroTotal, refMaxChroma)
# zSpecular2Diffuse(src(y,x),src(y+1,x).zMaxChroma()) # C call:
# zSpecular2Diffuse(src(y,x),src(y+1,x).zMaxChroma())
# update pixelStatus only if
# zSpecular2Diffuse() returned pStat as # update pixelStatus only if
# G_NOISE. # zSpecular2Diffuse() returned pStat as
if pStat == G_NOISE: # G_NOISE.
srcPixelStatus[y, x] = pStat if pStat == G_NOISE:
else: srcPixelStatus[y, x] = pStat
src[:, y, x] = iroRes else:
src[:, y, x] = iroRes
else:
iro = src[:, y + 1, x] else:
pStat = srcPixelStatus[y + 1, x] iro = src[:, y + 1, x]
pMaxChroma = srcMaxChroma[y + 1, x] pStat = srcPixelStatus[y + 1, x]
iroTotal = srcTotal[y + 1, x] pMaxChroma = srcMaxChroma[y + 1, x]
iroMax = srcMax[y + 1, x] iroTotal = srcTotal[y + 1, x]
iroMaxChroma = srcMaxChroma[y + 1, x] iroMax = srcMax[y + 1, x]
pStat, iroRes = specular_to_diffuse( iroMaxChroma = srcMaxChroma[y + 1, x]
iro, iroMax, iroMaxChroma, iroTotal, pMaxChroma) pStat, iroRes = specular_to_diffuse(
iro, iroMax, iroMaxChroma, iroTotal, pMaxChroma)
# update pixelStatus only if
# zSpecular2Diffuse() returned pStat as # update pixelStatus only if
# G_NOISE. # zSpecular2Diffuse() returned pStat as
if pStat == G_NOISE: # G_NOISE.
srcPixelStatus[y + 1, x] = pStat if pStat == G_NOISE:
else: srcPixelStatus[y + 1, x] = pStat
src[:, y + 1, x] = iroRes else:
src[:, y + 1, x] = iroRes
srcPixelStatus[y, x] = G_DIFFUSE
srcPixelStatus[y + 1, x] = G_DIFFUSE srcPixelStatus[y, x] = G_DIFFUSE
srcPixelStatus[y + 1, x] = G_DIFFUSE
pcount = count pcount = count
count, srcPixelStatus = init_labels(src, sfi, epsilon, srcPixelStatus) count, srcPixelStatus = init_labels(src, sfi, epsilon, srcPixelStatus)
...@@ -338,7 +463,6 @@ def init_labels(src, sfi, epsilon, srcPixelStatus): ...@@ -338,7 +463,6 @@ def init_labels(src, sfi, epsilon, srcPixelStatus):
old_settings = np.seterr(divide='ignore', invalid='ignore') old_settings = np.seterr(divide='ignore', invalid='ignore')
# to have initial labels # to have initial labels
count = 0
zTotalSrc = np.sum(src, axis=0) zTotalSrc = np.sum(src, axis=0)
diff_x_src = np.diff(zTotalSrc, axis=1) diff_x_src = np.diff(zTotalSrc, axis=1)
diff_y_src = np.diff(zTotalSrc, axis=0) diff_y_src = np.diff(zTotalSrc, axis=0)
...@@ -347,35 +471,29 @@ def init_labels(src, sfi, epsilon, srcPixelStatus): ...@@ -347,35 +471,29 @@ def init_labels(src, sfi, epsilon, srcPixelStatus):
diff_x_sfi = np.diff(zTotalSfi, axis=1) diff_x_sfi = np.diff(zTotalSfi, axis=1)
diff_y_sfi = np.diff(zTotalSfi, axis=0) diff_y_sfi = np.diff(zTotalSfi, axis=0)
dlog_src_x = np.log(abs(diff_x_src)) dlog_src_x = np.log(np.fabs(diff_x_src))
dlog_src_y = np.log(abs(diff_y_src)) dlog_src_y = np.log(np.fabs(diff_y_src))
dlog_sfi_x = np.log(abs(diff_x_sfi)) dlog_sfi_x = np.log(np.fabs(diff_x_sfi))
dlog_sfi_y = np.log(abs(diff_y_sfi)) dlog_sfi_y = np.log(np.fabs(diff_y_sfi))
dlogx = abs(dlog_src_x - dlog_sfi_x) dlogx = np.fabs(dlog_src_x - dlog_sfi_x)
dlogy = abs(dlog_src_y - dlog_sfi_y) dlogy = np.fabs(dlog_src_y - dlog_sfi_y)
maxY = src.shape[1] # maxY = src.shape[1]
maxX = src.shape[2] # maxX = src.shape[2]
valid_values = (G_BOUNDARY, G_NOISE, G_CAMERA_DARK, G_DIFFUSE)
for y in range(1, maxY - 1): notvalid_mask = np.in1d(
for x in range(1, maxX - 1): srcPixelStatus[1:-1, 1:-1], valid_values,
pStat = srcPixelStatus[y, x] invert=True).reshape(srcPixelStatus.shape[0] - 2, -1)
if pStat not in (G_BOUNDARY, G_NOISE, G_CAMERA_DARK, G_DIFFUSE): dlogx_mask = (notvalid_mask) & (dlogx[1:-1, 1:] > epsilon)
# Robby Tan's original C++ code doesn't check for dlogy_mask = (notvalid_mask) & (dlogy[1:, 1:-1] > epsilon) & \
# pStat!=G_DIFFUSE, but this speeds up the processing a lot, (~(dlogx_mask))
# and doesn't seem to affect the results. diffuse_mask = (notvalid_mask) & (~dlogx_mask) & (~dlogy_mask)
if dlogx[y, x] > epsilon: count = notvalid_mask.sum()
pStat = G_SPECULARX srcPixelStatus[1:-1, 1:-1][dlogx_mask] = G_SPECULARX
count += 1 srcPixelStatus[1:-1, 1:-1][dlogy_mask] = G_SPECULARY
elif dlogy[y, x] > epsilon: srcPixelStatus[1:-1, 1:-1][diffuse_mask] = G_DIFFUSE
pStat = G_SPECULARY
count += 1
else: