2

I'm trying to reconstruct surface/depthmap from normals using the frankot/chellappa algorithm. The rows and cols are the size of the img I'm trying to reconstruct the depth for.

I obtain the normal vectors like this:

rows, cols = imglist[0].shape
def find_NormAlbedo(sources, imglist, rows, cols):
    '''
    :param sources: a list of light source coordinates as [x,y,z] coordinates per light source
                    (shape (20,3) for 20 sources)
    :param imglist: a list of all images for one object
    :param rows: shape[0] of every image
    :param cols: shape[1] of every image
    :return: returns normals and albedo's for an object
    '''
    normal = np.zeros_like(imglist[0], dtype=np.ndarray)
    albedo = np.zeros_like(imglist[0])

    # for every pixel
    for x in range(rows):
        for y in range(cols):
            I = []  # intensity matrix of pixel x,y per image
            S = []  # lightsources
            for i in range(len(imglist)):
                img = imglist[i]
                I.append([img[x][y]])
                S.append(sources[i])

            # Least squares solution if S is invertible
            # pseudoinverse
            pseudoS = np.linalg.pinv(S)

            ntilde = pseudoS @ I
            p = np.linalg.norm(ntilde, ord=2)
            if p != 0.:
                n = ntilde / p
                n = n.flatten()
                # print(n)
                # print(n.shape)
            else:
                n = np.zeros_like(ntilde.flatten())

            normal[x][y] = n
            albedo[x][y] = p

    return normal, albedo

But suspect it's wrong because my albedo looks completely different from what I've seen in examples but have no clue where my mistake is...

enter image description here

Then I try to Get the surface from that using a wavepy function surface_to_grad:

def depthfromgradient(normalmap):
    '''
    :param normalmap: Previously obtained normals per pixel
    :return: Surface/Depth map from normalmap
    '''
    surfacep = np.zeros_like(normalmap)
    surfaceq = np.zeros_like(normalmap)
    for row in range(rows):
        for x in range(cols):
            #print(x)
            a, b, c = normalmap[row][x]
            #print(a, b, c)
            if c !=0:
                p = -a / c  # p=dZ/dx
                q = -b / c  # q=dZ/dy
                surfacep[row][x] = p
                surfaceq[row][x] = q
    return surface_from_grad.frankotchellappa(surfacep, surfaceq, reflec_pad=True)

My goal is to visualise the depthmap and the normalmap with cv.imshow(), but I'm not sure where I went wrong. These are my questions/ideas of where it went wrong:

-Is the albedomap plausible? If no, I think I misunderstood part of this algorithm.

-My depthmap has complex numbers, is this normal? Where do these come from?

-I looked at the shape of the normal map, the albedo map and the depth map, they all have shape (640, 500), yet I can only visualise the albedomap, the others give me the following error, what is the problem here?:

cv2.imshow('DepthMap', surface)
TypeError: Expected cv::UMat for argument 'mat'

Any help in narrowing down this problem would be welcome.

Note:I have tried converting everything to np arrays before using imshow().

10
  • To me, it looks like you have a range issue. For example, as if your math were producing values up to 1024, but you're only displaying the low-order 8 bits. Have you printed some pixel values to see what you're doing to the data? You might consider converting the pixel values to floating point in [0,1] before the computations, but that's just a guess. Jan 6, 2022 at 19:10
  • Why did you rebuild S for every pixel ? And you have one light source per image ?
    – David
    Jan 6, 2022 at 22:01
  • @David You're right, rebuilding S isn't necessary, thanks. Yes, I have 1 object, for which I have 20 images, each with a lightsource.
    – Kate S.D.
    Jan 7, 2022 at 9:09
  • @TimRoberts Thanks for your reply, My pixel values are all ints, when I start my computations, everything gets converted to float64 (ntilde is float64). When I look at the final normalmap (normal) it stores ndarrays with ndarrays(x,y,z)and albedomap (albedo), stores ndarray with uint8, but I have no issues with that one.
    – Kate S.D.
    Jan 7, 2022 at 9:23
  • I will be glad to help you more but without access to your image and the expected result this is difficult. Where did you take your example ? Can you share your image on GitHub or similar site ?
    – David
    Jan 7, 2022 at 10:59

2 Answers 2

1

The problem is exactly as I described in my very first comment. The values you get for p (which become your albedo image) range from 0 to 1998.34. When you store that in a byte, you're just getting the low-order 8 bits, which wrap around.

If you change this:

            albedo[x][y] = p

to this:

            albedo[x,y] = p/8

you'll see that the resulting image looks great.

By the way, there are several optimizations you can do. Where you have xxx[x][y] with a numpy array, do xxx[x,y] instead. When you build your lightsources, instead of

            I = []  # intensity matrix of pixel x,y per image
            S = []  # lightsources
            for i in range(len(imglist)):
                img = imglist[i]
                I.append([img[x][y]])
                S.append(sources[i])

do

            I = [img[x,y] for img in imglist]
            S = sources[:len(imglist),:]

and your obtainData function can be made more readable by doing:

    # Read all images
    paths = (
     fr".\PSData\PSData\{things[i]}\Objects\Image_01.png",
     fr".\PSData\PSData\{things[i]}\Objects\Image_02.png",
     fr".\PSData\PSData\{things[i]}\Objects\Image_03.png",
     fr".\PSData\PSData\{things[i]}\Objects\Image_04.png",
     fr".\PSData\PSData\{things[i]}\Objects\Image_05.png"
    )

    imgs = [cv2.imread(p,0) for p in paths]
...
    # Apply masks to images: cv2.bitwise
    imglist = [cv2.bitwise_or(img, img, mask=mask(img, threshold)) for img in imgs]
2
  • Dear Tim, thank you for your help. The albedo looks much better! I'm still having the issues visualising the surface however. I tried colorizing it using stackoverflow.com/questions/17044052/…, however, my image is completely lost (I simply get an image that is blue on 1 side, pink on the other side.) Do you have any idea what this might be?
    – Kate S.D.
    Jan 8, 2022 at 11:43
  • 1
    I didn't dig any further because I didn't want to install wavepy and colorsys. My best advice is to check your ranges. Print the values in surface, including mins and maxes, to make sure they really are complex and to find the ranges. See if the values returned from colorize are in the right range. Jan 9, 2022 at 4:28
0

Base on the documentation of the wavepy.surface_from_grad.frankotchellappa() the complex part can be ignored. Base on that they can be displayed with matplotlib to obtain an image like this one :enter image description here

import matplotlib.pyplot as plt
from matplotlib.ticker import LinearLocator
import numpy as np


ax = plt.figure().add_subplot(projection='3d')
surface = [[-2.19312825e+00-1.62679906e-02j,-1.76625653e+00-1.62093005e-02j ,-1.21359325e+00-1.63381200e-02j,-5.91501045e-01-1.45226036e-02j ,-9.24685295e-02-1.73373776e-03j, 1.59549135e-01+8.17785543e-05j , 2.90806910e-01-4.70409288e-05j, 3.47604091e-01+1.16492161e-05j],[-2.40280993e+00+1.62679906e-02j,-1.66061279e+00+1.62093005e-02j ,-1.11510109e+00+1.63381200e-02j,-3.08374007e-01+1.45226036e-02j , 1.57978453e-01+1.73373776e-03j, 2.59650686e-01-8.17785543e-05j , 3.63995725e-01+4.70409288e-05j, 4.01138015e-01-1.16492161e-05j],[-2.23283386e+00-1.62679906e-02j,-1.37689420e+00-1.62093005e-02j ,-7.79238609e-01-1.63381200e-02j, 6.13042608e-02-1.45226036e-02j , 4.21894421e-01-1.73373776e-03j, 4.22562434e-01+8.17785543e-05j , 4.81126228e-01-4.70409288e-05j, 4.97191034e-01+1.16492161e-05j],[-2.03380248e+00+1.62679906e-02j,-1.15214942e+00+1.62093005e-02j ,-5.31293338e-01+1.63381200e-02j, 3.06448040e-01+1.45226036e-02j ,6.43854839e-01+1.73373776e-03j, 5.97558594e-01-8.17785543e-05j ,6.18487416e-01+4.70409288e-05j, 6.16102854e-01-1.16492161e-05j],[-1.78298688e+00-1.62679906e-02j,-8.95377575e-01-1.62093005e-02j ,2.76063262e-01-1.63381200e-02j, 5.45799539e-01-1.45226036e-02j ,8.50153479e-01-1.73373776e-03j, 7.66200038e-01+8.17785543e-05j ,7.57462265e-01-4.70409288e-05j, 7.39255327e-01+1.16492161e-05j],[-1.48865585e+00+1.62679906e-02j,-6.06732343e-01+1.62093005e-02j ,2.25981613e-03+1.63381200e-02j, 7.88696843e-01+1.45226036e-02j ,1.04351159e+00+1.73373776e-03j, 9.19975122e-01-8.17785543e-05j ,8.83235485e-01+4.70409288e-05j, 8.50665864e-01-1.16492161e-05j],[-1.20317686e+00-1.62679906e-02j,-3.29863821e-01-1.62093005e-02j ,2.50911348e-01-1.63381200e-02j, 9.99471222e-01-1.45226036e-02j ,1.21789476e+00-1.73373776e-03j, 1.04762150e+00+8.17785543e-05j ,9.81522674e-01-4.70409288e-05j, 9.36002572e-01+1.16492161e-05j],[-7.74950625e-01+1.62679906e-02j, 9.04024675e-02+1.62093005e-02j
  ,6.51399438e-01+1.63381200e-02j, 1.32951041e+00+1.45226036e-02j
  ,1.36765625e+00+1.73373776e-03j, 1.12552107e+00-8.17785543e-05j
  ,1.03744750e+00+4.70409288e-05j, 9.82554438e-01-1.16492161e-05j]]

surface=np.array(surface)
X = np.arange(0,8)
Y = np.arange(0,8)
Z = [surface[x,y].real for x in X for y in Y]
xx, yy = np.meshgrid(X, Y)
Z= np.array(Z)
Z= np.reshape(Z, xx.shape)
colortuple = ('g', 'cyan')
colors = np.empty(xx.shape, dtype=str)
for y in range(len(Y)):
    for x in range(len(Y)):
        colors[y, x] = colortuple[(x + y) % len(colortuple)]

surf = ax.plot_surface(xx, yy, Z, facecolors=colors, linewidth=0)
ax.zaxis.set_major_locator(LinearLocator(6))

plt.show()

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.