CloudCompare forum

No problem, being challenged on the code is one of the point of being "open source".

In the 'flat plane' case, it's simply that Z = Z0 (the orthorectified image is a kind of flat cloud). This is clearly a quite strong approximation as complicated 3D structures in the image won't be properly re-projected in the orthophoto. Sadly without additional information such as a proper underlying 3D model we can't generate a 'clean' orthophoto (and we don't have enough keypoints to generate such a model). Hopefully It works very well for terrains and convex shapes (luckily most buildings are, etc.).

And as I already said before I also think that you we could use something else than the collinearity equation:

Oh when I was saying to use the transformation explicitly I wasn't thinking necessarily to use the collinearity equation. What we want is only to deduce the 2D pixel associated to a 3D point. So I guess that projecting the 3D point in the camera referential then using the focal and other intrinsic parameters could do. It's just a guess though ;)

I just remember that I didn't find the proper way to do this at the time I developed this (maybe I was simply tired and I missed something simple!). What we want is only a way to project a 3D point in the 2D image. If you find the way to do this without the collinearity equation I'd be happy to update the code.

From the Bundler manual, we have that

We use a pinhole camera model; the parameters we estimate for each camera are a focal length (f), two radial distortion parameters (k1 and k2), a rotation (R), and translation (t), as described in the file specification above. The formula for projecting a 3D point X into a camera (R, t, f) is:

P = R * X + t (conversion from world to camera coordinates)
p = -P / P.z (perspective division)
p' = f * r(p) * p (conversion to pixel coordinates)

where P.z is the third (z) coordinate of P. In the last equation, r(p) is a function that computes a scaling factor to undo the radial distortion:

r(p) = 1.0 + k1 * ||p||^2 + k2 * ||p||^4.

This gives a projection in pixels, where the origin of the image is the center of the image, the positive x-axis points right, and the positive y-axis points up (in addition, in the camera coordinate system, the positive z-axis points backwards, so the camera is looking down the negative z-axis, as in OpenGL).

Finally, the equations above imply that the camera viewing direction is:

R' * [0 0 -1]' (i.e., the third row of R or third column of R')

(where ' indicates the transpose of a matrix or vector).

and the 3D position of a camera is

-R' * t .

But this is only useful for the mapping of 3D points to the image and vice-versa- which is not what I pretend.

So I also miss what you are tying to do ;)

Aren't you trying to project all the images on a plane (e.g. on the ground)?

Yes. I want to project all image pixels on the same plane, but using the R and t matrices from Bundler, with intrinsics.
But I was trying to do the math I could not deduce how Z = Z0 results in those equations.

Well for the projection, if you define the extent of the projection plane (i.e. (Xmin, Ymin) and (Xmax, Ymax)) and a step (= output pixel size), this gives you a 2D grid on the ground. For each cell of the ground you back-project the corresponding "3D point" (with Z=0) in order to get the color. You'll have to do this in all images and find a proper way to select the right color as this point may be potentially seen in various images.

The result is a nice ortho image ;) I think I'll give it a try soon.

I think that is correct. And that is simlar to what is implemented in OrthorectifyAsImage.
But you use the a0, a1, a2, b0, b1, b2, c1, c2 indices estimated with your optimization routine. However I want to project those points with R and t from Bundler. Is this what you are thinking?

Hello, I managed to get some results. However, some images, on the extremes, are getting a lot of distortion (as we are projecting the images to the same origin).
Maybe we have to shift the cameras a bit as we build the mosaic. I can show you some results and I can give you my code, just tell me how.

Interesting! You can send me an email to cloudcompare [at] danielgm.net' (replace the [at] by @).