=Paper=
{{Paper
|id=Vol-1900/paper3
|storemode=property
|title=An algorithm for correcting X-ray image distortions caused by central projection
|pdfUrl=https://ceur-ws.org/Vol-1900/paper3.pdf
|volume=Vol-1900
|authors=Andrey V. Ustinov,Nataly Yu. Ilyasova,Nikita S. Demin
}}
==An algorithm for correcting X-ray image distortions caused by central projection ==
https://ceur-ws.org/Vol-1900/paper3.pdf
An algorithm for correcting X-ray image distortions caused by central
projection
A.V. Ustinov2 , N.Yu. Ilyasova1,2, N.S. Demin1
1
Samara National Research University, 34 Moskovskoye shosse, 443086, Samara, Russia
2
Image Processing Systems Unstitute - Branch of the Federal Scientific Research Centre "Crystallography and Photonics" of Russian Academy of Sciences, 151
Molodogvardeyskaya street, 443001, Samara, Russia
Abstract
We propose an algorithm for correcting X-ray image distortions caused by central projection. Relationships between coordinates of X-ray image
points obtained using parallel and central projection are derived. We describe two correction techniques using which the original image can be
made similar to the one based on parallel projection. In this way, the three-dimensional heart vessel model can be essentially simplified and
diagnostic parameters can be assessed with higher accuracy, resulting in a more accurate early disease diagnosis.
Key words: X-ray image; distortion correction; central projection
1. Introduction
Roentgenology is an extensively used and dynamic branch of medicine that uses X-ray imagery for the diagnosis of a variety
of diseases. By way of illustration, X-ray angiography is utilized for diagnosis of cardiovascular diseases [1]. As a rule, the
diagnosis is made based on visual assessment of angiograms, however, its accuracy essentially depends on the projection angle.
A three-dimensional model of heart vessels serves to visualize three-dimensional geometric and topological information, thus
enabling the diagnosis to be made with higher accuracy [2-6].
The initial data is a sequence of DICOM frames [7, 8]. The projection procedure is affected by a number of technical
limitations, resulting in the imagery characterized by a variety of distortions. Various techniques for distortion correction were
described in Ref. [9,10,13]. The X-ray imagery is obtained using specialized clinical equipment, such as an operating-room X-
ray unit C-ARM. Examples of such equipment are illustrated in Figure 1. The unit is composed of an X-ray source and receiver
connected by an arc-shaped holder freely moving on a support. With such a design, the X-ray camera has two degrees of
freedom. The camera can also move relative to the holder in the longitudinal direction, enabling the image to be scaled. The
spatial position of the camera is described by two angles: primary and secondary angles of the camera position.
a) b)
c) d)
Fig. 1. Specialized clinical equipment: a) AXIOM Multista, b) AXIOM Artis BC, c) primary angle, d) secondary angle.
The primary angle (denoted as α) is provided by rotating the camera holder and the support as a whole relative to the mount
beam. The secondary angle β (analogous to the geographic latitude) is provided by sliding the arc-shaped holder on the support
guide, with the camera and X-ray source moving on a circular arc. The imaging is done by a diverging X-ray, with the
divergence angle defined by technical characteristics of the scanning device and usually found within 10-12°. As a result, the X-
3rd International conference “Information Technology and Nanotechnology 2017” 10
Computer Optics and Nanophotonics / A.V. Ustinov , N.Yu. Ilyasova, N.S. Demin
ray image is observed in the central projection. Based on such projections, it is not possible to reconstruct a model of the original
object without additional information concerning imaging conditions, which is not available. At the same time, the original
image can be reconstructed from parallel projections without use of any additional information [7,8]. Our idea is that effects
caused by central projection can partially be compensated for by making the X-ray image look as if it were built using parallel
projections.
By correcting X-ray image distortions caused by central projection, the three-dimensional heart vessel model can be
essentially simplified [12] and diagnostic parameters can be assessed with higher accuracy [5], resulting in a more accurate early
disease diagnosis.
2. Mathematical model
Because the refractive index for X-rays for all substances is very close to one, it is impossible to make a collimator that
converts a divergent beam from a source close to a point beam into a parallel beam. This leads to an object image deformation
even at the planar object. At straight rays falling (the primary and secondary angles are zero), the distortion reduces to a change
in scale - the image remains similar to the image obtained by parallel projection. At inclined falling an additional distortion
appears: the scale becomes different on image area, which leads to a violation of similarity (a change in the proportions of the
object in the case of oblique incidence also occurs in parallel projection). We propose an algorithm for compensating of central
projection influence (divergent rays). It should be noted that full compensation is possible only for a planar object, because
otherwise it is principally impossible to indicate precise values of some parameters. Nevertheless, the achieved compensation is
enough for more precision of three-dimensional trace process described in Ref. [7, 12].
In order that to compensate the central projection influence we are need to get a relation of image point coordinates at parallel
and central projection. For this purpose a calculation of world and planar coordinates of projection point is required. Let a point
S is a light source and distance OS H is given. A direction to the light source is defined as OS nH . Point A ( x0 ; y0 ; h) is
an intersection of projection plane with straight-line having a directing vector SA, where h is distance from image plane to the
object at zero angles. Now we obtain the projection of point A ( x; y; z ). Planar coordinates of point are
u OAu xux yu y zuz , v OAv xvx yvy zvz , where u , v are basic vectors of planar coordinates system. If the plane is
defined by an equation Ax By Cz D 0 and straight-line is defined by formulas x x0 lt ; y y0 mt ; z z 0 nt , then
coordinates of point A are obtained at substitution of value t derived from equation ( Al Bm Cn )t Ax0 By0 Cz 0 D 0 .
In our case:
D 0; ( A; B ; C ) ( n x ; n y ; n z ); ( x 0 ; y 0 ; z 0 ) ( x0 ; y 0 ; h ); ( l ; m ; n ) AS nH ( x0 ; y 0 ; h ) ( Hn x x0 ; Hn y y 0 ; Hn z z 0 ).
At the parallel projection: (l ; m; n ) n ( A; B; C ) . Taking in account these values we get the projection point coordinates:
x0 ( B 2 C 2 ) y0 AB hAC
x ,
1 ( Ax0 By0 Ch) / H
x0 AB y0 ( A2 C 2 ) hBC
y , (1)
1 ( Ax0 By0 Ch) / H
x0 AC y0 BC h( A2 B 2 )
z .
1 ( Ax0 By0 Ch) / H
Fig. 2. Effect of divergent beam. Left - straight falling (α=β=0), right - inclined (α=45°, β=0).The second frame has larger crosshairs
on left part in image than in right part.
3rd International conference “Information Technology and Nanotechnology 2017” 11
Computer Optics and Nanophotonics / A.V. Ustinov , N.Yu. Ilyasova, N.S. Demin
Numerators are projection point coordinates at the parallel projection and L Ax0 By0 Ch in denominators is signed
distance from point A to the projection plane.
Taking in account formulas for planar coordinates we get:
u u par 1 L / H ,
v v par 1 L / H .
In some cases value of H we can get from DICOM-file header, else we can evaluate it with aid of mira image - special etalon
image (grid in fig.2). In item 3 an algorithm invented for evaluation of H is described. About value of h: if the object is not
planar then h can not be obtain in principle, since a reconstruction of three-dimensional object is need to this. For approximate
solution of correction task we propose three approaches:
1) If object length in OZ axis direction compare with H is small and object is near to receiver then we use h=0.
2) If object length is small but object is not near to receiver there we use some average value, for example, a heart size.
3) Assuming that distortions are not large (three-dimensional tracing process is not failed) we build three-dimensional tree
ignoring the central projection influence. That h equal to z-coordinate of corresponding tree point. After correction the three-
dimensional tree building is repeated. But this variant increases calculation volume.
3. Determination of distance from the source to the receiver with usage of mira
For the determination of distance from the source to receiver with usage of mira, we preliminarily consider an auxiliary case.
Let we have straight falling rays (at neglecting distortion of central projection) and a task is planar: source, object and
coordinates origin lie in one plane perpendicular to the plane of the receiver. The distance from point source to receiver plane
(the origin) is H, the distance from point object to the origin is h. Then a following relation holds: xc H H h x p , it is
connects coordinates of the object image on the receiver plane with the parallel projection and the central one.
Now we consider our case. The task is still planar - the source and the non-point object passing through the origin lie in one
plane perpendicular to the plane of the receiver. Let he object is segment of length 2a lying horizontally in this plane and the
origin is its center, i.e. h=0 for segment center. A left point will be the first, and a right point will be the second. Rays fall bias at
angle φ to vertical (the angle is counter on clockwise). The distance of the point source to the origin is H. Under this scheme of
observation the central point of segment-object retains its zero coordinate on segment-image, but it will not be a center of
segment-image. For the parallel projection, image stays symmetrical: if we use a length (instead of coordinate with sign) then
x p1 x p 2 a cos .
Because of the segment is not lie in in receiver plane, value h 0 for its ending points. A module of this value is equal for
the first and the second points: h1 h2 a sin .
But for the first (left) point it has minus sign. For absolute values we have final result:
H
xc1 a cos ;
H a sin
(2)
H
xc 2 a cos .
H a sin
This means that if we know angle φ and the lengths of segments measured on the image which are equal on the object, then
we can calculate distance H and the length of segment by formulas:
2 xc1 xc 2
a ;
( xc1 xc 2 ) cos
(3)
2 xc1 xc 2
H tg .
xc 2 xc1
A really used mira image contains many segments. Since mira grid on object is rectangular then it is desirable to use of
images photographed when one of camera orientation angle is zero. This give two advantages: angle φ is calculated easily (it is
equal to modulus of the second orientation angle) and lines of equal scale stay be straight and parallel to image sides. If β=0,
these lines are vertical, if α=0, ones are horizontal. For definiteness, we shall consider case of β=0. At this condition we have
the following algorithm.
Step 1. We do search of mira grid crosshair nearest to the image center.
Step 2. We do step back on some squares to left (not necessary on line containing center). A difference of horizontal
coordinates of central and new crosshair is value xc1 .
3rd International conference “Information Technology and Nanotechnology 2017” 12
Computer Optics and Nanophotonics / A.V. Ustinov , N.Yu. Ilyasova, N.S. Demin
Step 3. We do step back on some squares to right (not necessary on line containing center). A difference of horizontal
coordinates of new and central crosshair is value xc 2 .
Step 4. The sought distance H is calculated by the second formula (3). If it is negative then modulus is used - minus sign means
that right segment is shorter of left one that takes place of rays falling left from vertical. At obtaining of formulas (2) we assumed
rays falling right from vertical.
4. Correction distortion of central projection with cycle on source image
The first method correction distortion of central projection consists in following: for points on source image we search
corresponding points on corrected image.
If x0 , y0 are found then correction is implemented by formulas:
Ax By0 Ch
u par u 1 0
H
.
Ax0 By0 Ch
v par v 1
H
Solution of equation set for values x0 , y0 is following:
a (u b1h ) a2 (v b2 h)
x0 4
a1a4 a2 a3
. (4)
a (v b2 h ) a3 (u b1 h )
y0 1
a1a4 a2 a3
Knowing expressions for vectors dependence on primary and secondary angles
n ( sin ;sin cos ;cos cos )
u (cos ;sin sin ;cos sin )
v (0;cos ; sin )
we get expressions for coefficients a1 , a2 , a3 , a4 , b1 , b2 :
a1 u H sin cos ,
a2 u H sin cos sin sin ,
a3 v H sin ,
(5)
a4 v H sin cos cos ,
b1 u H cos cos cos sin ,
b2 u H cos cos sin .
But if we want to form the correct image this method is undesirable: because of discretization part of output image points (in
area of stretching) will not have a prototype. There is a grid from points without the prototype on fig. 3 that usually is not suited
for us.
Fig. 3. Correction with cycle on source image. Left – source image, right – corrected one. A white grid is points without prototype.
3rd International conference “Information Technology and Nanotechnology 2017” 13
Computer Optics and Nanophotonics / A.V. Ustinov , N.Yu. Ilyasova, N.S. Demin
5. Correction distortion of central projection with cycle on corrected image
For avoidance of empty points we employ the second method correction effects of central projection. Here the standard
approach of image spatial transformation (for example, rotation or reflecting in non-planar mirror) is used - a cycle employs on
output image and the prototype of current point is found on the source image.
If x0 , y0 are found then the prototype point has coordinates:
Ax By0 Ch
u u par / 1 0 ,
H
Ax By0 Ch
v v par / 1 0 .
H
Coordinates of point x0 , y0 is calculated by formulas:
a4 (u par b1h) a2 (v par b2 h)
x0 ,
a1 a4 a2 a3
(6)
a1 (v par b2 h) a3 (u par b1h)
y0 .
a1a4 a2 a3
Coefficients are equal to a1 cos ; a2 sin sin ; a3 0; a4 cos ; b1 cos sin ; b2 sin , where , are primary
and secondary angles of camera rotation, and h is distance from the object to the projection plane.
After a substitution of these values we find simpler formulas:
a4 (u par b1h) a2 (v par b2 h)
x0 ,
a1a4
(7)
v par b2 h
y0 .
a4
On fig. 4 we see that the grid from points without prototype is absent, i.e. the correction is more precise.
Fig. 4. Correction with cycle on corrected image. Left – source image, right – corrected one.
6. Conclusion
In this work, we propose an algorithm enabling X-ray image distortions caused by central projection to be corrected for.
Considering that it is not possible to reconstruct the original object from a three-dimensional heart vessel model without
additional information concerning the imaging conditions, which is not available in most cases, the process gets more
complicated. At the same time, the reconstruction based on parallel projections requires no additional information. We propose a
relation of image point coordinates at parallel and central projection. We describe two correction techniques using which the
original image can be made similar to the one based on parallel projection. In this way, the three-dimensional heart vessel model
can be essentially simplified and diagnostic parameters can be assessed with higher accuracy, resulting in a more accurate early
disease diagnosis.
3rd International conference “Information Technology and Nanotechnology 2017” 14
Computer Optics and Nanophotonics / A.V. Ustinov , N.Yu. Ilyasova, N.S. Demin
Acknowledgements
The work was carried out with the partial support of the Ministry of Education and Science of the Russian Federation within
the framework of the activities of the Program for Enhancing the Competitiveness of the SSAU among the world's leading
scientific and educational centers for 2013-2020; Grants of the Russian Foundation for Basic Research No. 14-07- 97040, No.
15-29- 03823, No. 15-29- 07077, No. 16-57-48006, No. 16-41-630761; Program № 6 of fundamental research Department of
Nanotechnologies and Information Technologies of the Russian Academy of Sciences "Bioinformatics, modern information
technologies and mathematical methods in medicine" 2016 -2017.
References
[1] Chandra T, Pukenas B, Mohan S, Melhem E. Contrast-Enhanced Magnetic Resonance Angiography. Magnetic Resonance Imaging Clinics of North America
2012; 20(4): 687–698.
[2] Ilyasova N. Computer Systems for Geometrical Analysis of Blood Vessels Diagnostic Images. Optical Memory and Neural Networks (Information Optics)
2014; 23(4): 278–286.
[3] Ilyasova N. Methods to evaluate the three-dimensional features of blood vessels. Optical Memory and Neural Networks (Information Optics) 2015; 24(1):
36–47.
[4] Ilyasova N. Evaluation of Geometric Characteristics of the Spatial Structure of Vessels. Pattern Recognition and Image Analysis 2015; 25(4): 621–625.
[5] Ilyasova NYu. Estimating the geometric features of a 3D vascular structure. Computer Optics 2014; 38(3): 529–538.
[6] Ilyasova NYu. Methods for digital analysis of human vascular system. Literature review. Computer Optics 2013; 37(4): 511–535.
[7] Ilyasova NYu, Kupriyanov AV, Khramov AG. Information technologies of image analysis in medical diagnostics. Radio and communication, 2012.
[8] Soifer VA. Computer Image Processing. Part II: Methods and algorithms: Appendix A2. Biomedi-cal Images Processing. VDM Verlag, 2009.
[9] Hong C, Lee D-H, Han BS. Characteristics of geometric distortion correction with increasing field-of-view in open-configuration MRI. Magnetic Resonance
Imaging 2014; 32(6): 786–790.
[10] Cardoso PL, Dymerska B, Bachratá B, Fischmeister FPhS, Mahr N, Matt E, Trattnig S, Beisteiner R, Robinson SD. The clinical relevance of distortion
correction in presurgical fMRI.NeuroImage 2016.
[11] Menga C, Zhang J, Zhou F, Wang T. New method for geometric calibration and distortion correction of conventional C-arm. Computers in Biology and
Medicine 2014; 52: 49–56.
[12] Ilyasova NYu, Kazanskiy NL, Korepanov AO, Kupriyanov AV, Ustinov AV, Khramov AG, Computer technology for reconstructing the 3D structure of
coronary arteries from angiographic projections. Computer Optics 2009; 33(3): 281–318.
[13] Moravec J, Hub M. Automatic correction of barrel distorted images using a cascaded evolutionary estimator. Information Sciences 2016; 366: 70–98.
3rd International conference “Information Technology and Nanotechnology 2017” 15