

Epipolar ConsistencyContact: Andre Aichert
Table of Contents
1 About this page
Epipolar Consistency has been shown to be a powerful tool in calibration and motion correction for flat panel detector CT (FDCT) and other Xray based applications. The epipolar consistency conditions (ECC) apply to projection data directly and can be used to correct for 3D parameters without 3D reconstruction. The computation of an epipolar consistency metric is realtime capable for small sets of images and can be parallelized for large sets of images. This document briefly describes an algorithm to compute a metric for the epipolar consistency conditions (ECC) between to Xray images. We assume that the reader is familiar with epipolar geometry [4] and has basic knowledge of Xray physics and the Radon transform [3]. We will briefly describe an improved version of the algorithm presented by Aichert et al. [1].
2 Data
Let
and
denote two Xray projections with known projection matrices
and finite source positions
and
accordingly, where
denotes equality up to positive scalar multiples. Let
be two lines in real projective twospace of the image
and
, respectively. Then they are corresponding epipolar lines if there exists an epipolar plane
with
where
denotes the pseudoinverse. An epipolar plane is any plane which contains both source positions
A line
can be characterized as the set of solution to the homogeneous equation
or by its angle
and signed distance to the origin
Every pixel in these images represents an Xray intensity, which can be transformed into an integral over absorption coefficients
where
is the function of Xray absorption along the backprojection ray of the pixel
with a distance to the source
. Please see Figure 1↑, left, center, for two example Xray images before preprocessing with several exemplary epipolar lines. Further, let
denote the Radon transform of the thus transformed image
with
compare Figure 2↑. See Figure 3↓ for radon derivatives of the preprocessed images from Figure 1↑. Note that a line bundle forms a sinusoid curve in Radon space.
3 Epipolar Consistency
The epipolar consistency conditions state that the derivative in orthogonal direction of the integral along an epipolar line is redundant with the same of the corresponding epipolar line in the other image
where
denotes the Radon transform of the image
and
is its derivative, with
and the symmetry
Since all epipolar lines intersect in the epipole, we have a oneparameter family of lines in each image with 11 correspondences. The lines are related via their epipolar planes, which in turn form a pencil of planes. Our goal is to parametrize the pencil of epipolar planes
with a single angle
. Using the relationship between epipolar lines and planes, we can define a metric for consistency between two images as the integral over the plane angle
Since all epipolar planes must contain the line connecting the two source positions
, we compute the consistency metric numerically as follows:
A large portion of the algorithm discussed in the previous section is linear. Our goal is to write it as a matrix multiplication by the vector
.
4 Implementation4.1 PlückerMatrices
This section assumes basic knowledge of about lines in two and threespace. The baseline
is the line in threespace connecting the two source positions
. We express the line as an antisymmetric matrix.
The antisymmetric Plücker matrix is defined as
where
defined by the six distinct elements
They always fullfill the GrassmannPlücker relation
. These 6 values are also called Plücker coordinates of the line. We recommend the paper by Blinn [2] about the Plücker matrix. Using Plücker coordinates, we can easily express the direction of the line
and the normal to the plane containing the line and passing through the origin
The 2Dequivalent exists as a
Plücker matrix that is related to the crossproduct via
where
We can thus write
and we introduce a dual Plücker matrix
Plücker matrices allow us to directly express the common plane of a point
and a line
4.2 ECC algorithm using Plücker coordinates
Assuming that the object is located in the world origin, we can define the reference plane directly
We can also define an orthogonal plane, which contains the baseline
If we assert that the normals of the two planes have unit length, we can thus express any epipolar plane as a linear combination of the two. Specifically, we have
with
. Projection to the two detectors is also a linear operation, so we arrive at the following compact formulation
with
and
. With
and
we have
and
, respectively. The two lines are generally not orthogonal.
5 Octave Source Code5.1 Loading Data
We provide two example data sets. A pumpkin as an example for an organic structure and a toy gun as an example of an artificial structure. We use NRRD as an image format and provide a C++ headeronly library to load and save NRRD images. We recommend Fiji/ImageJ as an image viewer. For MATLAB/Octave, there exists a nice function to load NRRD files from Mathworks FileExchange, distributed under its own BSD license. All NRRD files contain a humanreadable ASCII header. Our projection images contain a
projection matrix (MATLABlike format). In addition, there is a comment in the last line with the offset in bytes to the raw float data, followed by an empty line.
NRRD0004 dimension: 2 encoding: raw endian: little sizes: 2480 1920 type: float Projection Matrix:=[...] # Offset to raw data: 570 bytes.
By convention, we refer to the function
, in words the “derivative in
direction of the Radon transform” as dtr, compare Figure 3↑. These dtrs have to be precomputed from the projection images. We recommend to use our minimal working data set to get started quickly. For indepth investigations, we provide you with two complete data sets of 15 images each, along with a Windows/CUDA binary to compute the Radon transform and tools to convert and edit NRRD files (see Downloads↓ Section). Please see the first lines in ecc_demo.m and adjust the paths to two input images as desired.
5.2 Computing the ECC
The algorithms described in this paper are implemented in the function
% Mapping to epipolar lines
[K0 K1]=ecc_computeK01(P0,P1, n_x, n_y);
Here, P0 and P1 are the projection matrices and n_x and n_y the number of image columns and rows. The corresponding implementation located in the file ecc_demo/ecc_computeK01.m. The sampling of the ECC can be summarized in the following code
consistency=0; % resulting metric
dkappa=0.2*pi/180
for kappa=[pi*0.5:dkappa:pi*0.5]
x_kappa=[cos(kappa); sin(kappa)];
l0_kappa=K0*x_kappa;
% Sample Radon derivatives at corresponding locations
[v0 x0, y0, w0]=ecc_sample_rda(l0_kappa, rda0, range_t);
[v1 x1, y1, w1]=ecc_sample_rda(l1_kappa, rda1, range_t);
consistency+=(v0v1)*(v0v1)*dkappa;
end %
Here, dkappa is arbitrarily chosen to be
. Adjust this value to obtain a dense sampling of lines. The code in ecc_demo.m computes redundant values for both
at the same time and is able to terminate early, when the lines no longer intersect the image. It also contains a lot of code for visualization. A brief version of ECC computation can be found in the file ecc_demo/ecc_compute_consistency.m .
5.3 Heuristic for the sampling distance
The number of samples drawn from the Radon transforms depends linearly on the angular step
. On the one hand, choosing
too small results in slow runtimes, on the other hand, choosing
too large results in undersampling and less accurate ECC values. For best results,
should be chosen as large as possible, while consequent samples drawn from Radon space should be no farther apart than the size of one Radon bin. Let
and
denote the Radon bin size in angular and distance directions respectively. Further, let
and
, respectively, denote the image of the world origin in image 0 and 1. This point serves as a reference point and is thus always located inside the object and close to the center of the image. To find a good heuristic for an appropriate value for
, we want either the angle
or distance
between the epipolar lines
and
to be less than the Radon bin size
where
is the angle between two homogeneous lines and
the euclidian distance between a line and an image point. The same should hold for
,
and
, respectively. Whether
or
is the limiting dimension depends on the location of the epipoles. If they are close to the image, a differential change in
will affect
more than
, and vice versa if the epipoles are far away. One should choose
adaptively for each pair of images.
5.4 Screenshot of the demoDownloads
References[1] . Epipolar Consistency in Transmission Imaging. Medical Imaging, IEEE Transactions on, 34(11):22052219, 2015. [2] . A homogeneous formulation for lines in 3 space. ACM SIGGRAPH Computer Graphics, 11(2):237241, 1977. [3] . Computed Tomography: From Photon Statistics to Modern ConeBeam CT. Springer, 2008. [4] . Multiple View Geometry in Computer Vision. Cambridge University Press, 2004. 