Frame localisation optical projection tomography

We present a tomographic reconstruction algorithm (flOPT), which is applied to Optical Projection Tomography (OPT) images, that is robust to mechanical jitter and systematic angular and spatial drift. OPT relies on precise mechanical rotation and is less mechanically stable than large-scale computer tomography (CT) scanning systems, leading to reconstruction artefacts. The algorithm uses multiple (5+) tracked fiducial beads to recover the sample pose and the image rays are then back-projected at each orientation. The quality of the image reconstruction using the proposed algorithm shows an improvement when compared to the Radon transform. Moreover, when adding a systematic spatial and angular mechanical drift, the reconstruction shows a significant improvement over the Radon transform.


Stereoscopic imaging
When the features or fiducial markers in one view are uniquely identifiable, the stereoscopic imaging of scenes allows for the triangulation of individual features in three dimensional space (known as world points), see Figs. 1 and 2 for the coordinate system which describes this geometry. Triangulation requires that each feature is detected in both images of a stereo imaging system and for these detections to be correctly associated with one another. This is known as the correspondence problem. Various methods exist to ensure that features are detected from image data and accurately associated between two cameras or views 6 and the properties of scaleindependent features and their surrounding pixel environment in one image can thus be matched to a similar feature in a second image.
Coordinates in two adjacent views with a common epi-pole (the vector connecting the O and O ′ , see Fig. 2) are related by the essential matrix (E) for uncalibrated cameras and the fundamental matrix (F) for calibrated cameras. Their properties are described by: where K is a matrix that converts image plane coordinates to camera pixel coordinates and where p refers to a point in the image plane.

The proposed algorithm (flOPT)
The motion of a rotating sample, as in an OPT acquisition, with a transformation matrix ( [R | T] ) in view of a fixed camera is analogous to the motion of a camera around the scene with the inverse transformation matrix.
During an ideal OPT acquisition, a marker will appear to follow an elliptical path in the xy image plane. For the volume reconstruction procedure, there is a fitting step to recover the path of the fiducial marker, which is used to correct the sinogram before applying the inverse Radon transform. This type of reconstruction not only ignores any mechanical jitter of the sample, but also any affine, systematic, mechanical drift (in X, Y , Z, θ, φ, ψ ). This can be rectified by recovering the complete non-scaling transformation for every projection. Now, using two adjacent images of a scene (separated by some rotation and translation) world points in 3D space may be triangulated within the scene given the rotational and translational matrices of the respective camera views. Once a sufficient amount of fiducial markers are reliably tracked from the first to the second image, either one of the fundamental or essential matrices can be computed. Using the factorisation of one of these matrices, between each adjacent view of a rotating scene, the translation and rotational matrices can be recovered. www.nature.com/scientificreports/ To reconstruct the image, we compute F for the current image and the first image using 5 or more fiducial markers; having additional beads helps to remove ambiguity and increase confidence in F. Once F is calculated, it is decomposed into R n and T n between each view n and n + 1 . The image at view n + 1 is then back projected along the virtual optical axis within a virtual volume where the sample will be reconstructed. The size of this back projection and virtual volume is chosen to be suitably large, preventing the loss of important data. The recovered transformation matrices are then matrix inverted and applied to the back projection of the image to realign the rays in the volume to their respective source positions as shown in Fig. 3.
In both cases, a decomposed F matrix will produce four possible transformation pairs (R, T; R, − T; − R, T; − R, − T). Once the transformation matrix between the current view (n) and the first view is calculated, the Figure 2. Epi-polar geometry described for two adjacent views (or cameras of a scene). Coordinates as expressed in Fig. 1a with prime notation ( ′ ) denoting the additional right camera view. Transforming from right to left camera-centered coordinates ( X ′ c to X c ) requires a rotation (R) and a translation (T). To find the correct matrix between the n = 0 and n = 1 orientations, each of the four matrices are compared to an ideal matrix which is composed using a priori knowledge of the likely angle of rotation of the system's imaging properties.

Verification of the proposed algorithm
To verify the validity and quality of the proposed reconstruction algorithm, the image of Zelda, superposed with an orthogonal image of Cameraman, is used as a testcard volume. Virtual fiducial beads are dispersed in the volume to track the rotation and translation of the image. The reference image is then rotated through 128  www.nature.com/scientificreports/ angles over 2π radians and projected along the Y axis, then an image slice in (X, Y) is taken to create a single line projection, shown three dimensionally in Fig. 4. This is repeated for each angle, with each line projection stacked to create a sinogram.
In the standard approach for OPT reconstruction, the sinogram undergoes the inverse Radon transform, as shown in Fig. 4j, followed by post-filtering. This step is substituted for the proposed algorithm; in Fig. 5a the two techniques are compared for ideal conditions of smooth, predictable rotation. The proposed algorithm produces a faithful reconstruction on the original image, as shown in Fig. 6d. Fig. 5b illustrates the strong overlap of the images produced by the new algorithm and the Radon transform when considering the histogram of the absolute pixel-wise difference between the original source image and the respective reconstructions. The proposed algorithm generates lower deviance from the source image than the Radon transform. The mean square errors (MSE, see Eq. (4)) of the new algorithm and the Radon transform are 15.01% and 14.84%, respectively, see Fig. 5b for a histogram of a pixel-wise comparison.
where Y is the vector of observed values and Ŷ i is mean of the ith value of the predicted values The more challenging case of a sample drifting systematically along the X axis, with a constant velocity, was then considered. This drift produced a helical path of a single fiducial within the sample, see Fig. 6b. In Fig. 6c, the Radon transform fails to produce a recognisable reproduction of the test image with the addition of a slight helicity to the rotation. The proposed algorithm produces an equivalent result to that of a sample rotating without any systematic drift, see Fig. 6c. In Fig. 5c the respective reconstructions from each algorithm were compared, as before, while the helical shift was incremented. See Fig. 6b for a sinogram of a sample (shown in Fig. 6a) wherein a helical shift has been induced. When using correlation as a metric of reproduction quality, the new algorithm fares slightly worse at zero helicity, with 94% correlation compared to the Radon transform at 96%. As expected,

Recovery of R and T using matrix decomposition.
To quantitatively verify that the matrix decomposition technique was valid and robust, the accuracy of the reproduction of R and T was tested directly. The original R and T matrices were computed and compared to R and T generated from matrix decomposition. This absolute difference was computed element-wise in each matrix and then an average for each matrix was taken.
Overall, the worst-case scenario produced a percentage error of 2% (see Fig. 7 for full statistics). The accuracy of the calculated R and T deteriorated when adding in additional degrees of combined movement, but with no correlation between the degree of helicity and the error produced. The translation matrix (T) was consistently more accurately reproduced, which is likely due to it having fewer available degrees of freedom.

Discussion
A new algorithm for reconstructing OPT data has been demonstrated. The new algorithm uses multiple fiducial markers to recover the matrix which describes the rotation and translation of the sample. The quality of the reconstructions shows a slight improvement when compared to the standard Radon transform, with a great effect when a systematic drift is introduced. The accuracy of the decomposition of F into R and T was compared to the ground truth matrices. The element-wise absolute difference x−y 2(x+y) of each matrix was averaged across the matrix for R and T. In the worst-case scenario, a maximum of 2% average absolute difference was found between ground truth and recovered matrices, suggesting that the technique is robust to various forms of drift in all dimensions and general instability. Such an algorithm could be used to minimise ghosting effects seen in real samples, particularly in samples where slipping is likely to occur, such as in gels or in cheaper OPT systems which tend to be more mechanically unstable and imprecise. In particular the imaging of large mobile gels is set to become more prevalent given the surge of new techniques in Expansion Microscopy 8 , whereby fragile expanded samples embedded in thin lubricious gels.

Future work
The proposed algorithm relies on triangulation between two view points. However, it is possible to use three separate views to reconstruct a scene, one such approach being quaternion tensors 9 . Working with tensors is more complex, but a future iteration of the algorithm presented here may benefit from using three views to provide a more accurate transformation matrix. Beyond three views, there is currently no mathematical framework for four or more views. If such tools were to be developed, it may be possible to have the algorithm described above be a non-iterative, single-shot reconstruction from pixels to voxels. Fiducial markers could also be extracted from the image texture alone, circumventing the need for the additional beads embedded in the sample. To find such correspondences, points with similar local texture are found and matched in between each image using standard algorithms such as SIFT 10 and RANSAC 11 . This was attempted in this work, however, the errors introduced into the transformation matrices make this approach currently unviable; and so by requiring bright punctuate fiducial markers the burden of collecting the fiducial coordinates is shifted to well established curve fitting algorithms that are robust to noise.

Figure 7.
Box plots demonstrating that the rotational and translations matrices can be recovered accurately from fiducial marker positions. Panels (a,b) introduce an angular drift during rotation, to an observer at the detector this would appear as a tip of the sample towards them, causing precession. Panels (c,d) introduce a lateral drift in X causing a helical path to be drawn out. In all cases, the percentage error introduced by the the addition of undesirable additional movements was on the order of < 2%.(note that errors in recovering translation are much larger given a smaller helical shift as the percentage error of the recovery of the translation matrix is broadly constant).