anuranbaka / OpenDTAM

An open source implementation of DTAM
Other
286 stars 150 forks source link

How to produce point cloud using DTAM? #34

Open XiaoshuiHuang opened 8 years ago

XiaoshuiHuang commented 8 years ago

Hi , Because I want to use DTAM to conduct 3D reconstruction, I need point cloud results. Actually, I only got depth map. How can I obtain point cloud (x,y,z) for every point? I have changed the parameters value like #8 . thetaStart = 0.6 * layers thetaEnd = 0.03 * layers thetaStep = 0.99 epsilon = 0.003 * layers lambda = 0.3 / layers And when I do the following step, it got wrong. ray = ( K * R ).inv() * T; P=1/(optimizer.depthMap())*ray + cameraCenter;

It is wrong because ray is a array 3*1, while (1/optimizer.depthMap()) is the depth of whole image. It cannot directly multiple. And I do not know what cameraCenter is in your example "testprog.cpp"? Can you give me more detail information about how to code this part? Thank you for your answer! Xiaoshui

XiaoshuiHuang commented 8 years ago

Also, I have doubt if the ray is right? because in such computation way, ray is 1_3, not 3_1 and 3*3. How can I get properly 3D point based on depth?

anuranbaka commented 8 years ago

It depends on what coordinate frame you want the results in. I can't quite tell what's going on from what you have written, but if you are trying to get the point cloud in the world coordinate frame, that isn't too hard.

If your feeling really nice you might make this into a function inside of the Optimizer class, since I think people have this question a lot. I'm currently not in a position to compile the code (need to update my cuda install for gcc 4.9), but if you want to try the following and make a pull request, I'll just incorporate it.

For convenience, the CostVolume object exposes its [world_x; world_y; layer; 1] -> [pixel_x; pixel_y; layer; 1] matrix as CostVolume.projection. Inverting this and applying it to the unnormalized depth map will give a point cloud. I'm currently not in a position to compile the code (need to update my cuda install for gcc 4.9), but try the following code snippet:

//Returns a point cloud in the world coordinate frame associated to the cost volume
// Format is a Mat with [x,y,z] as the three channels
cv::Mat Optimizer::pointCloud(){
    cvStream.waitForCompletion();
    Mat tmp;
    cv.a.download(tmp); //move the raw depth map off the GPU 
                        //(could use "d" instead of "a" for more smooth biased estimate)

    Mat_<float> depth = tmp;
    Mat_<Vec3f> xyin(depth.rows,depth.cols);
    Mat_<Vec3f> xyout(depth.rows,depth.cols);
    float* pt=(float*) (xyin.data);
    float* dp=(float*) (depth.data);
    for(int i=0;i<_depth.rows;i++){
        for(int j=0;j<_depth.cols;j++,pt+=3,d++){
            pt[0] = j; //px_x
            pt[1] = i; //px_y
            pt[2] = *dp; //layer
        }
    }
    perspectiveTransform(xyin, xyout, cv.projection.inv());
    xyout=xyout.reshape(3, depth.rows);
    return xyout;
}
anuranbaka commented 8 years ago

Oh and fyi, the formula you quoted, which I see you got from #8 , is correct only if you get everything in the correct coordinate system and if the x,y,d inputs are taken elementwise. In both the original post by avanindra and my reply we were each implicitly assuming the other one knew that. We also assumed the other one could figure out when to transpose implicitly. Sorry for the confusion, but it should be equivalent to the above. Also, make sure you are on the experimental branch if you use the above, since the main branch projection matrix didn't comply with spec.

XiaoshuiHuang commented 8 years ago

Hi , I have tried two ways these days to produce point cloud, but it got the wrong results. One way is what you tell me. Add the function in optimizer and call the function. It got the following result. https://drive.google.com/file/d/0B4tp_cfWpm9QU1BEMVU5VGk5Yk0/view?usp=sharing However, this result is meaningless for me. I want point cloud to depict the scene's surface so that I can use the point cloud to reconstruct surface and so on.

Then, I read deeply into the code, I find a function 'diagnosticInfo' , it may give the results I want.After my try, it got the following results. https://drive.google.com/file/d/0B4tp_cfWpm9QdEZSZnV1RkE1Ym8/view?usp=sharing

Can you give a hand to me? May it is a bit wrong about the coordinate system.

By the way, I want to result like this video. https://www.youtube.com/watch?v=dKksb1wQmzs Also, if I got the correct result, I will pull the request to the branch. Thank you!

Xiaoshui

anuranbaka commented 8 years ago

Can you post your depth map as an image? It should evolve like https://youtu.be/Df9WhgibCQA?t=1m16s (on the right side).That should give me an idea if it's the settings or the algorithm that isn't working. For reference, your first result is not too surprising to me, since badly reconstructed points spray off to infinity. You should look and see if the reconstructed surface is buried in the tip of that point cloud.

You're correct that the diagnosticInfo function can do what you want internally, but it also does occlusion filtering and such to simulate a virtual camera.

By the way, the thing you posted is either not DTAM or has half of the algorithm turned off, which you can tell because the map doesn't optimize at all. It looks like straight cost volume reconstruction to me.