Optical realization of the radon transform

Tali Ilovitsh; Asaf Ilovitsh; John Sheridan; Zeev Zalevsky

doi:10.1364/OE.22.032301

1. Introduction

The Radon transform of a given object is the sum of its line integrals along all given angles [1]. The attractiveness of this transform is the existence of an inverse formula which enables the reconstruction of the object out of its projections [2,3]. The Radon transform is widely used in tomographic imaging [4], image processing [5], pattern recognition [6,7] and in motion detection [8].

In this paper we will deal with the practical computation of the Radon transform via optical means. Most of the systems that perform the Radon transform require the mechanical rotation of the object in order to acquire its projections at all angles. In addition, they rely on digital processing, which requires a large number of calculations and as such, cannot be performed in real time. This is the motivation for an all optical Radon transform system, able to operate in real time, giving rise to new potential applications [9].

The traditional setup used to perform the optical Radon transform of a two dimension (2D) input data involves the illumination of an object with a laser beam combined with the use of a cylindrical lens. The light is captured using a photo detector, which gives a signal proportional to the line integral of the reflectivity along the line. This represents a single point in the Radon space. The line of light is then swept perpendicularly to itself to generate one line and the complete 2D Radon transform is generated by varying the orientation of the line of light for all angles. This rotation is achieved using a Dove prism and a synchronous motor [10]. The resulting system is complicated and isn't applicable to real time imaging purposes, i.e. the full Radon transform can take several seconds to minutes [11]. Recent advances involve the use of only four projections of the input object. In this case the setup includes four cylindrical lenses placed in the aperture plane and four linear array sensors parallel to the cylindrical lenses [8]. Since the system uses only four projections, the results obtained are less accurate and require a special sensors array followed by computer post processing.

In this paper we present a novel optical approach for realizing the 2D Radon transform. The proposed technique yields the full transform in a single frame, using a simple 4F setup and a single optical element. The optical output is given in Cartesian coordinates, which can then be digitally transformed into polar coordinates to match the mathematical definition.

This paper is organized as follows: First, the mathematical background is introduced. Second, computer simulations of the proposed system with a Shepp-Logan head phantom object are presented, followed by a discussion of the choice of parameters used to implement the optical element and a comparison of the obtained results to the ideal mathematical transform. Third, experimental results are presented, which validate the proposed approach.

2. Theoretical background

For a continuous 2D function f(x,y), a one dimension (1D) projection is formed by performing a line integration of the image intensity f(x,y) along a line, that is at a distance L from the origin and at angle θ to the x-axis. All points on this line satisfy the following equation:

L = x \cos (θ) + y \sin (θ)

The mathematical definition of the conventional Radon transform is:

g (L, θ) = \iint^{​} f (x, y) δ (x \cos (θ) + y \sin (θ) - L) d x d y

The collection g(L,θ) at all angles θ, is called the Radon transform of the input image f(x,y), where the term δ is the Dirac delta function. The proposed setup for optically performing the 2D Radon transform is illustrated in Fig. 1.

Fig. 1 The optical setup which consists of a 4F system, where an optical element is placed in the 2F plane.

Download Full Size | PDF

The input object is denoted by s(x,y). We start by performing its 2D optical Fourier transform (FT):

\tilde{S} (x, y) = \iint s (x', y') \exp (\frac{- 2 π i (x x' + y y')}{λ F}) d x' d y' = \tilde{S} (r, θ)

The FT, S̃(r,θ), of the input object, is then multiplied by an optical element that is placed in the 2F Fourier plane.

The optical element is denoted by T(r,θ) and defined as:

T (r, θ) = \frac{\exp (2 π i β θ r)}{r}

This multiplication further propagates through another 2F system and the obtained result at the output plane is:

O (ρ, α) = \int_{0}^{2 π} \int_{0}^{\infty} \tilde{S} (r, θ) T (r, θ) \exp (\frac{- 2 π i}{λ F} r ρ \cos (θ - α)) r d r d θ

Substituting the optical element from Eq. (4) yields:

O (ρ, α) = \int_{0}^{2 π} \int_{0}^{\infty} \tilde{S} (r, θ) \exp (\frac{- 2 π i}{λ F} r ρ \cos (θ - α)) \exp (2 π i β θ r) d r d θ

By approximating the integral over the angleθ, to a sum of θ = nδθ, where n is an integer index and δθ is the angular resolution. The result becomes:

O (ρ, α) \approx \sum_{n} O_{n} (ρ, α) = \sum_{n} \int_{0}^{\infty} \tilde{S} (r, n δ θ) \exp (\frac{- 2 π i}{λ F} r ρ \cos (n δ θ - α)) \exp (2 π i β n δ θ r) d r

Note that each term of the summation, O_n(ρ,α) appears at a different spatial location due to the term of exp(2πißnδθr). Therefore, since each term is spatially separated from all other terms, every term will appear sampled in the output plane at different angular directions α. Specifically at each angular sample α = nδθ:

O (ρ, α = n δ θ) = \int_{0}^{\infty} \tilde{S} (r, n δ θ) \exp (\frac{- 2 π i}{λ F} r ρ) \exp (2 π i β n δ θ r) d r = \tilde{\tilde{S}} (\frac{ρ}{λ F} - β n δ θ, α = n δ θ)

Where S͌ is the FT of S̃ along the radial axis. The output field equals:

O (ρ, α) = \sum_{n} \tilde{\tilde{S}} (\frac{ρ}{λ F} - β n δ θ, α = n δ θ)

S͌ is the FT transform of S̃, which is the FT of the input object s. Thus S͌ is the original inspected object (up to inversion). Therefore, the result is replicas of the input object along the radial axis. Each replica corresponds to a particular discrete angular direction nδθ. When this output is transformed to polar coordinates it becomes the Radon transform of the input, as it is defined according to the central slice theorem [12]. A schematic sketch of the obtained output of the proposed technique is shown in Fig. 2. Figure 2(a) is an input object of a slit and Fig. 2(b) is the obtained output, i.e. the Radon transform of the input object in Cartesian coordinates.

Fig. 2 A schematic sketch of the obtained output of the proposed technique. (a) an input object of a slit. (b) The obtained output that is the Radon transform of the input object in Cartesian coordinates.

Download Full Size | PDF

3. Simulation results

A set of simulations were performed using MATLAB (MathWorks, Natick, MA, USA). The simulated model was the 4F system, illustrated in Fig. 1. Light illuminates an input object, which is a Shepp-Logan head phantom (Fig. 3(a)). This object was chosen in order to demonstrate the applicability of the proposed method for a grayscale object as well.

Fig. 3 Simulation results. (a) An input object of a Shepp-Logan head phantom. (b) The mathematical Radon transform of the input object. (c) The obtained output image using the proposed setup. (d) The output image transformed to polar coordinates.

Download Full Size | PDF

The Radon transform of the input object, computed using the Radon function in Matlab, is shown in Fig. 3(b). After propagating through the first 2F system, the light comes along the optical element, described in Eq. (4). The light further propagates through the second 2F system. The output image is presented in Fig. 3(c). This image is given in Cartesian coordinates and becomes the Radon transform after it is post-processed to polar coordinates (Fig. 3(d)).

An intuitive way to understand the result obtained in Fig. 3(c) is to think of taking the input object and rotating it with respect to a certain point, that is the common origin of the axes. Since the input image has an ellipsoid shape, the result has a round shape with different thicknesses that correspond to the dimensions of the main axes of the input object. The correlation coefficient between the mathematical and the optical Radon transform results is calculated to be 0.99, which supports the validity of the proposed method.

The parameter β in the numerator of the optical element is a constant that controls the radius of the rotation of the input object around the central point. Empirically it was chosen to be β = 10⁵. A larger value will yield an image that is too large to be captured by the camera, while a smaller value will result in the overlap of the rotated replicas of the input object.

4. Experimental results

The proposed method was tested using the experimental setup shown in Fig. 4.

Fig. 4 The experimental setup (right to left). A 4F system that consists of a collimated green laser beam at a wavelength of 532 nm, illuminating an input object. The light propagates through a 2F system, and illuminates an optical element realized using an amplitude slide attached to a SLM. The reflected light from the optical element propagates through a second 2F system and is captured using a standard USB camera.

Download Full Size | PDF

The experimental setup is a 4F system that consists of a collimated green laser beam at a wavelength of 532 nm (Photop DPGL-2100F), illuminating an input object. The light propagates through a 2F system, where it passes through the optical element T(r,θ).

The implementation of the optical element was divided into two separate elements, a phase element and an amplitude element. The phase term exp(2πiβθr) with parameter β = 10⁵ was realized using a SLM (Holoeye SLM device HEO 1080P) [13]. This phase-only SLM is a reflective liquid-crystal-on-silicon micro display with 1920 × 1080 pixels resolution and pixel’s pitch of 8μm. The phase used is shown in Fig. 5(a). The amplitude term 1/r was realized using an amplitude slide that was attached to the SLM. Since the SLM is based on reflection, the laser beam passes twice through the amplitude slide and therefore the amplitude of the slide varied as1/√r and was implemented with pixel of sizes of 54μm and is shown in Fig. 5(b).

Fig. 5 The proposed optical element. (a) The phase of the optical element. (b) The amplitude of the optical element.

Download Full Size | PDF

The light reflected from the optical element further propagates through the second 2F system and is captured by a standard USB camera (Thorlabs DCC1545M). Three requirements ensure the functioning of the setup; first, the laser beam must hit both the center of the SLM and the center of the amplitude slide. Second, the SLM and the amplitude slide must be aligned. Third, the reflected light from the SLM must pass perpendicularly through the second 2F system. When the three conditions aren't fulfilled, the output image contains artifacts.

In order to experimentally validate the proposed system, two different input objects were tested. The first input object was a slit of dimensions 0.25x0.5 mm, as shown in Fig. 6(a). The focal length of all the lenses used in the 4F system was 75 mm. This focal length was chosen so that the FT of the input object fills a significant part of the active area of the SLM. First, a simulation of the intended result using the proposed system was made and is shown in Fig. 6(b). The experimental results using the proposed system are shown in Fig. 6(c). The similarity of the results to those expected (see Fig. 2) is clearly visible. In addition, there is also a good agreement between the simulated results and the experimental results (Fig. 6(b) and 6(c)). Next, in order to compare the results to the mathematical Radon transform of the input object (that is shown in Fig. 6(d)), the obtained result was transformed into polar coordinates and is shown in Fig. 6(e).

Fig. 6 Experimental results.(a) The input object is a slit of dimensions of 0.25x0.5 mm. (b) The simulated output image using the proposed setup. (c) The experimental output image using the proposed setup. (d) The mathematical Radon transform of the input object. (d) The output image transformed to polar coordinates. (e) The output image of (c) transformed to polar coordinates.

Download Full Size | PDF

There is a resemblance between the mathematical Radon transform and the experimental results (Fig. 6(d) and 6(e)) is notable. The calculated correlation coefficient between these two images is 0.99. The second input object was an amplitude slide showing the outline of an airplane with dimensions of 1.5x2.0 mm, as shown in Fig. 7(a). The experimental results are shown in Fig. 7(b), the experimental results transformed to polar coordinates are shown in Fig. 7(c). In order to reconstruct the input object from our experimental Radon transform, the inverse Radon transform was computed using the iradon function in Matlab was applied to Fig. 7(c), and the result is shown in Fig. 7(d). The correlation coefficient between the reconstruction of the input object and the input object visible (Fig. 7(a) and 7(d)) is 0.993, which reinforces the proposed approach. The quality of the reconstruction is very high with compare to that of traditional optical methods.

Fig. 7 Experimental results. (a) The input object, an amplitude slide showing the outline of an airplane with dimensions of 1.5x2.0 mm. (b) The experimental results obtained using the proposed setup. (c) The experimental results transformed to polar coordinates. (d) The reconstruction of the input object using the inverse Radon transform applied on the experimental result in (c).

Download Full Size | PDF

The method that utilizes only 4 projections [8] was able to perform the transform on simple objects, mainly slits. Since their reconstruction contained many artifacts, they mainly used it for motion detection. Another work was able to reconstruct the input object using 45 projection, however with high blurriness [9]. Therefore, the proposed method is very attractive.

The proposed setup has a couple of limitation; first, the active area of the imaging detector limits the input image size to a few mm. This is because the output image is the rotation of the input object around a central point which yields a larger shape. One way to overcome this problem is to use a 4F system with a magnification smaller than 1, which will enable to work with larger input objects.

Another consideration is the SLM resolution. When working with larger pixel sizes, an aliasing occurs which results in unwanted replications of the output that distorts the image. For the proposed setup, the largest pixel size of the SLM for which the obtained image can be reconstructed is three times bigger.

5. Conclusions

A new optical system for performing the Radon transform was presented and validated experimentally. The system is relatively easy to implement, align and maintain. The proposed technique eliminates the need to capture multiple images of the rotated object in order to obtain the Radon transform. It also significantly reduces the number of calculations required in the standard implementations of the Radon transform. Since the technique is fast, it can be used for real time imaging and therefore may permit new applications in image processing, motion detection and tomography.

References and links

1. J. Radon, “Über die Bestimmung von Funktionen durch ihre Integralwerte Längs gewisser Mannigfaltigkeiten,” Ber. Sächs. Akad. Wiss 69, 262–278 (1917).

2. M. Nishimura, D. Casasent, and F. Caimi, “Optical inverse Radon transform,” Opt. Commun. 24(3), 276–280 (1978). [CrossRef]

3. F. Colonna, G. Easley, K. Guo, and D. Labate, “Radon transform inversion using the shearlet representation,” Appl. Comput. Harmon. Anal. 29(2), 232–250 (2010). [CrossRef]

4. S. R. Deans, The Radon Transform and Some of Its Applications (New York: Wiley, 1983).

5. D. Mendlovic, R. G. Dorsch, A. W. Lohmann, Z. Zalevsky, and C. Ferreira, “Optical illustration of a varied fractional Fourier-transform order and the Radon-Wigner display,” Appl. Opt. 35(20), 3925–3929 (1996). [CrossRef] [PubMed]

6. C. G. Ho, R. C. D. Young, C. D. Bradfield, and C. R. Chatwin, “A Fast Hough Transform for the Parametrisation of Straight Lines using Fourier Methods,” Real Time Imaging 6(2), 113–127 (2000). [CrossRef]

7. W. Götz and H. Druckmüller, “A fast digital Radon transform—An efficient means for evaluating the Hough transform,” Pattern Recognit. 29(4), 711–718 (1996). [CrossRef]

8. Y. Kashter, O. Levi, and A. Stern, “Optical compressive change and motion detection,” Appl. Opt. 51(13), 2491–2496 (2012). [CrossRef] [PubMed]

9. V. Farber, Y. August, and A. Stern, “Super-resolution compressive imaging with anamorphic optics,” Opt. Express 21(22), 25851–25863 (2013). [CrossRef] [PubMed]

10. W. H. Steier and R. K. Shori, “Optical Hough transform,” Appl. Opt. 25(16), 2734 (1986). [CrossRef] [PubMed]

11. S. Woolven, V. M. Ristic, and P. Chevrette, “Hybrid implementation of a real-time Radon-space image-processing system,” Appl. Opt. 32(32), 6556–6561 (1993). [CrossRef] [PubMed]

12. H. H. Barrett and W. Swindell, Radiological Imaging: The Theory of Image Formation, Detection, and Processing‏ (Elsevier, 1996).

13. O. Fixler and Z. Zalevsky, “Geometrically superresolved lensless imaging using a spatial light modulator,” Appl. Opt. 50(29), 5662–5673 (2011). [CrossRef] [PubMed]

Optical realization of the radon transform

Abstract

1. Introduction

2. Theoretical background

3. Simulation results

4. Experimental results

5. Conclusions

References and links

Cited By

Figures (7)

Equations (9)

Optics Express