I have implemented Stochastic Tomography (http://www.cs.ubc.ca/labs/imager/tr/2012/StochasticTomography/) into a C++ framework. It reconstructs 3D volumetric data from a number of 2D projections at different angles. My implementation already works for synthetic data, but the next logical step for me now is to test it with real data sets.
The author provides real sample data on the above homepage and luckily there also is a description for the geometric calibration and the projection model (http://www.cs.ubc.ca/labs/imager/tr/2012/StochasticTomography/files/StochasticTomography-DatasetsREADME.pdf). Since i only need to be able to project 3D points into 2D images, the most relevant part can be found in the last paragraph of section 3.2 (page 2-3). Which i also will copy into this question:
To perform the reverse process, that is, compute an image coordinate for an known world-space coordinate as is needed by Stochastic Tomography, Equation 3 must be numerically inverted. This requires finding the pixel coordinates producing a ray that passes through the world-space point, for which we use Newton’s method. We remove the rotation component of Equation 3 by multiplying the world-space point by $R^{-1}$ and use as an initial guess at the centroid of the point correspondences detected by CALTag. To perform Newton’s method, we define two functions $f_x(u, v)$ and $f_y(u, v)$ which evaluate the distance in the x and y directions respectively between the closest point on the ray through the current, estimated solution $[u, v]^T$ and the world-space point being projected into the image. We then iteratively solve the following equation to obtain an update to the solution $[∆u, ∆v]$, evaluating the differential terms with finite differences.
Equation 3 (Full version in paper):
$$\begin{bmatrix} x_i \\ y_i \\ z_i \\ 1 \end{bmatrix} = RT \begin{bmatrix} \beta^T_x\Gamma_i \\ \beta^T_y\Gamma_i \\ 0 \\ 1 \end{bmatrix}$$
Equation 4 (iterative solver): $$\begin{bmatrix} \frac{\delta f_x(u, v)}{\delta u} && \frac{\delta f_x(u, v)}{\delta v} \\ \frac{\delta f_y(u, v)}{\delta u} && \frac{\delta f_y(u, v)}{\delta v} \\ \end{bmatrix} \begin{bmatrix} \Delta u \\ \Delta v \end{bmatrix} = \begin{bmatrix} -f_x(u, v) \\ -f_y(u, v) \end{bmatrix}$$
I have different files for all the cameras which contain the point correspondences and the angles needed. I understand most of the process described, but having some problem wrapping my head around how the functions $f_x(u, v)$ and $f_y(u, v)$ are defined and used. Especially this part: "..between the closest point on the ray through the current, estimated solution $[u,v]^T$ and the world-space point being projected into the image". Does this mean i have to calculate everything described in section 3.2 for the current estimation $[u, v]^T$ to determine the ray and thus also to be able to determine the distance in x and y direction? Because to me it sounds like $f_x$ and $f_y$ are defined differently from the process described above.
Yes that is how it was was meant. The initial guess is always the centroid of the known point-correspondences and improved iteratively by calculating:
$$\begin{bmatrix} \frac{\delta f_x(u, v)}{\delta u} && \frac{\delta f_x(u, v)}{\delta v} \\ \frac{\delta f_y(u, v)}{\delta u} && \frac{\delta f_y(u, v)}{\delta v} \\ \end{bmatrix} \begin{bmatrix} \Delta u \\ \Delta v \end{bmatrix} = \begin{bmatrix} -f_x(u, v) \\ -f_y(u, v) \end{bmatrix}$$
where $f_x(u,v)$ and $f_y(u, v)$ are the distance values in $x$ and $y $ direction for the guess $(u, v)$ to the known point in 3D. This also means that we need to calculate the pierce points and the consequent ray through the volume for each iteration.