Isotropic edge-enhancement by the Hilbert-transform in optical tomography of phase objects

Areli Montes-Perez; Cruz Meneses-Fabian; Gustavo Rodriguez-Zurita

doi:10.1364/OE.19.005350

1. Introduction

Tomography is a non-invasive technique aimed at obtaining spatial distributions of some physical property lying within a given object. A probe beam, often an input wave, travels across the object under inspection, carrying away information about its interior. At the output, the probe is modified in a distribution known as projection. In general, each projection depends on the angle the probe forms with a direction fixed to the object (projection angle). Optical tomography is a special case of tomography in which the probe is a light beam [1]. Optical tomographic techniques are of particular importance in the medical field, because these techniques can provide non-invasive diagnostic images. When the beam travels through the object interior with no deviation (refractiveless, difractionless, isotropic and without scattering), the probe beam can be considered as composed of straight rays. Being such a beam composed of only parallel rays, it is spoken of parallel ray tomography. The present proposal could be a very particular case of the optical diffraction and/or diffuse optical tomography [2]. A feature of this case is that a description with Fourier methods turns out to be very useful. In this context, a parallel projection can be conveniently described with a Radon transform (RT), as it has been done in several fields, including computer aided x-ray tomography, radar imagery, and geophysics. In optical tomography, a lot of research has been carried out to study transparent objects, such as temperature in flames considering axial symmetry [3], or considering error by this approximation [4], and also refraction effects as well [5-6]. Most of the aforementioned techniques obtain information from the detected interferograms [7]. These techniques have been oriented to obtaining the distribution of phase object slices. However, in many practical situations, the phase gradient is too large, and conventional interferometry methods are not practicable. In order to overcome or minimize these deficiencies, other approaches have been proposed to obtain some modified object slice instead of the object slice itself, such as obtaining the angular derivative by using speckle interferometry [8], or the directional derivative by using a semiderivative filter in a 4f optical imaging system [9].

These two last techniques are spatial filtering applications. They thus require of proper spatial filters in the projection detection system which demand of certain corresponding spatial amplitude-transmittance distributions to be done. These distributions are not very easily fabricated. But from the practice of spatial filtering, it is rather well-known that the Hilbert transform (HT) is a phase imaging tool with some resemblances to derivation regarding directional edge detection through its enhancement. Although the resulting images achieved with the implementation of HT do not, in general, provide a direct way to quantitatively recover the phase distribution involved, they can offer related qualitative information which can suffices for many applications. One advantage of HT over derivation schemes is that its implementation requires of a phase step as the spatial filter. This filter type is less elaborated than a derivation filter. Another feature relates to the reconstructed tomographic image resulting as a consequence of using phase projection information gathered with a HT detection approach. In this report, we show that an intermediate stage between detection and reconstruction is not necessarily required to achieve an edge enhanced tomographic image. Besides, the achievement of this tomographic image needs of no further modification in the usual algorithms for reconstruction.

In the present manuscript, the optical HT is implemented by using a 4f optical imaging system, and a phase step of π radians used as a spatial filter at the Fourier plane of the projection detection system [10–12], in order to obtain the HT of the field leaving the phase object. In this proposal, the irradiance of field at the image plane is considered as a special sinogram and is demonstrated that, with this sinogram, an enhancement of borders of tomographic images is obtained, a technique that could be particularly useful for phase-only objects. A theoretical background to support the method is presented, along with some numerical and experimental results.

2. Basic considerations

Let us consider a phase object described in a rectangular coordinate system and disposed to turn around the z-axis together with the x and y axes, so that its refraction index can be described by $f (x, y, z)$ . When a non-inclined homogenous monochromatic plane wave of amplitude A crosses this object (Fig. 1 ,) the phase of field leaving the object is given by [1]

A_{ϕ} (p, z) = A \exp [i \frac{2 π}{λ} {\overset{\lor}{f}}_{ϕ} (p, z)],

where ϕ is the azimuth rotation angle with respect to the x-axis, so that the rotated coordinates

(p, p_{⊥})

obey

p = x \cos ϕ + y \sin ϕ

and

p_{⊥} = - x \sin ϕ + y \cos ϕ

.

i = \sqrt{- 1}

is the imaginary unit, λ is the wavelength of the light, and

{\overset{\lor}{f}}_{ϕ} (p, z) = \int_{L}^{} f (x, y, z) d p_{⊥}

is the optical path. For z constant, the optical path could be related with the RT by means of

ℜ {f (x, y, z)} = {\overset{\lor}{f}}_{ϕ} (p, z) = \iint_{\infty} d ξ d η f (ξ, η, z) δ (p - ξ \cos ϕ - η \sin ϕ),

where

ℜ {\dots}

denotes the RT operator [1],

{\overset{\lor}{f}}_{ϕ} (p, z)

is the profile or projection at a ϕ angle at z height, and p is the projection coordinate. As it is well known in tomographic theory, these projection data constitute the basis for the reconstruction of an object slice [1]; then, the refraction index

f (x, y, z)

can be reconstructed from a set of projections within the range

(0, π)

by using the back-projection algorithm, for instance [13].

Fig. 1 (color online) A beam of parallel rays crosses a phase object at z constant, $z = z_{0}$ . (a) 3-D view showing a rotated reference system (x, y) about z with respect to the $(p, p_{⊥})$ system. (b) 2-D view of plane $z = z_{0}$ with p axis translated.

Download Full Size | PDF

We propose to carry out the optical HT of the field leaving the object by means of a 4f _L optical imaging system, as it is depicted in Fig. 2 , where the field indicated in Eq. (1) is considered as the entrance of the Fourier system. Lens $L_{1}$ performs the FT of the wave leaving the object, located on the back focal plane (system’s frequency plane). Lens $L_{2}$ performs the inverse FT after a sign filter is placed at the frequency plane, and then, at the image plane, the resulting field is the HT of entrance field, which is symbolically expressed by ${\overset{\cap}{A}}_{ϕ} (p, z) = H {A_{ϕ} (p, z)}$ , with H{...} indicating the operator of HT, and is described by

{\overset{\cap}{A}}_{ϕ} (p, z) = ℑ^{- 1} {- i sgn (w) {\tilde{A}}_{ϕ} (w, ζ)},

where

ℑ^{- 1} {\dots}

is the operator of the inverse FT,

{\tilde{A}}_{ϕ} (w, ζ) = ℑ {A_{ϕ} (p, z)}

,

sgn (w)

is the filter function, and

(w, ζ)

are spatial frequencies. By applying the convolution property of the FT in Eq. (3) we have

{\overset{\cap}{A}}_{ϕ} (p, z) = \frac{δ (z)}{π p} \otimes A_{ϕ} (p, z),

meaning that for every z-value, the one-dimensional HT in the p-direction is carried-out.

δ (z) / π p = ℑ^{- 1} {- i sgn (w)}

is known as the system’s impulse response. The irradiance of the field at the image plane would thus be written as

I_{ϕ} (p, z) = {\overset{\cap}{A}}_{ϕ} (p, z) \overset{\cap}{A_{ϕ}^{*}} (p, z) = {| {\overset{\cap}{A}}_{ϕ} (p, z) |}^{2} = {| \frac{δ (z)}{π p} \otimes A_{ϕ} (p, z) |}^{2},

where the symbol * indicates the complex conjugate. Now, considering Eq. (5) as the projection data, the symmetry property of RT should be

I_{ϕ + π} (p, z) = {| \frac{δ (z)}{π p} \otimes A_{ϕ + π} (p, z) |}^{2} = {| \frac{δ (z)}{π (- p)} \otimes A_{ϕ} (- p, z) |}^{2} = I_{ϕ} (- p, z),

where

A_{ϕ + π} (p, z) = A_{ϕ} (- p, z)

is verifiable from Eq. (1) because

{\overset{\lor}{f}}_{ϕ + π} (p, z) = {\overset{\lor}{f}}_{ϕ} (- p, z)

is a well-known property [1]. For the zero-moment of RT, it could be established that

\int_{- \infty}^{\infty} d p I_{ϕ} (p, z) = {\tilde{I}}_{ϕ} (w, z) |_{w = 0} = {- i sgn (w) {\tilde{A}}_{ϕ} (w, z) \otimes - i sgn (w) \tilde{A_{ϕ}^{*}} (- w, z)} |_{w = 0},

where is assumed that

{\tilde{I}}_{ϕ} (w, z) = ℑ {I_{ϕ} (p, z)} = ℑ {{\overset{\cap}{A}}_{ϕ} (p, z)} \otimes ℑ {\overset{\cap}{A_{ϕ}^{*}} (p, z)}

. Substituting the integral definition of convolution, and evaluating at

w = 0

, we obtain

\int_{- \infty}^{\infty} d p I_{ϕ} (p, z) = \int_{- \infty}^{\infty} d u {sgn}^{2} (u) {\tilde{A}}_{ϕ} (u, z) \tilde{A_{ϕ}^{*}} (u, z),

where u is an auxiliary variable of integration. Now, considering

{sgn}^{2} (u) = 1

, and applying the Parseval theorem, Eq. (8) could be written as

\int_{- \infty}^{\infty} d p I_{ϕ} (p, z) = \int_{- \infty}^{\infty} d p A_{ϕ} (p, z) A_{ϕ}^{*} (p, z),

indicating that the total intensity is the same both at the object and image planes, since the used filter is an only-phase one. Substituting Eq. (1) into Eq. (9), we finally have

\int_{- \infty}^{\infty} d p I_{ϕ} (p, z) = \int_{- Γ / 2}^{Γ / 2} d p A^{2} = A^{2} Γ,

where Γ is the detector width, and the irradiance out of this region is taken to be zero. Considering the irradiance in the image plane as the projection data, the zero-moment of the RT is the total energy at the detector, finite and independent of ϕ. That way, the zero-moment of the RT is satisfied. In addition, as it has been demonstrated, the irradiance at the image plane (Eq. (5) complies both with the symmetry property and the zero-moment of the RT; therefore, it can be considered as a projection of a given object, and therefore its tomographic reconstruction must render a consistent image.

Fig. 2 (color online) 4f _L optical imaging system with a sign filter at the Fourier plane to obtain the HT of the projection data leaving the phase object.

Download Full Size | PDF

3. Numerical Simulation

Figure 3 shows a numerical simulation about the tomographic reconstruction of a phase object slice at a height z kept constant. The images are presented in 8-bit gray levels. Figure 3-a1 shows an object slice $f (x, y)$ with $200 \times 200$ pixels, consisting of a rectangle of unit height, while in Fig. 3-a2 its sinogram is shown, obtained by using Eq. (2). In this image, p has 200 data, and the number of projections is 200 for $ϕ \in (0, 2 π)$ . Thus, the sinogram has $200 \times 200$ pixels. Column 2b shows the optical field $A_{ϕ} (p, z)$ leaving an object slice $f (x, y)$ , where the sinogram shown in Fig. 3-a2 has been used as the optical path of a beam that crosses the object slice such as indicated in Eq. (1); the real and imaginary parts are shown in Fig. 3-b1 and b2 respectively. For this simulation, $A = 1$ . Figure 3c shows the HT of the field shown in Fig. (3b) obtained by using Eq. (4). The real and imaginary parts are shown in Fig. 3-c1 and 3-c2, respectively. The upper row of the next column shows the irradiance of ${\overset{\cap}{A}}_{ϕ} (p, z)$ , which is considered as a modified sonogram: the irradiance-Hilbert-sinogram (IH-sinogram); in Fig. 3-d2 the zero moment of the RT for the IH-sinogram is shown. It is important to note that the symmetry property and zero-moment of the RT can be checked as expected. Under the filtered back-projection algorithm, the tomographic reconstruction can be carried out. Figure 3e shows the obtained reconstruction both in gray-levels and 3-D plots. In this algorithm a filter rectangular was employed and the impulse’s response of this filter was used to filtering the projection data by using numerical convolution.

Fig. 3 (color online) Numerical simulation. (a) object slice (1) and sinogram (2), (b) the field leaving the phase object slice: real (1) and imaginary (2) parts (c) HT of the field leaving the phase object slice, (d) IH-sinogram (1) and the zero-moment of the RT (2), and (e) Edge-enhancement reconstruction: gray levels (1) and 3-D plot (2).

Download Full Size | PDF

In Fig. 4 , two examples are shown: on row 1, a uniform ring is considered as the object slice, and on row 2, a non-symmetric object slice is considered. On column (a), the object slices are shown; on column (b), the IH-sinograms are shown; on column (c) the zero-moments of the RT are illustrated; finally, on columns (d) and (e), both the reconstructions in gray-levels and on 3-D plots are respectively shown.

Fig. 4 (color online) Numerical simulations. (a) Two test slices: uniform ring and non-symmetric object, (b) IH-Sinograms, (c) Zero-moment of RT, d-e) reconstructions showing isotropic edge-enhancement.

Download Full Size | PDF

4. Experimental results

Experimental implementation of HT was based on the 4f _L optical system depicted in Fig. 2. The positive lenses $L_{1}$ and $L_{2}$ have the same focal distance $f_{L} = 479 m m$ ; the filter is a phase-step of π radians; and $λ = 632.8 n m$ . The phase object is kept inside of a liquid gate filled with immersion oil of nominal refraction index 1.515, thus the refractionless condition is fulfilled; the object is attached to the axis of a step motor to turn it around the z-axis, and it is put just before the entrance plane of the imaging system. The step motor has a nominal step of 1.8°, and it is driven by a visual interface in the computer via an electronic driver (not shown). A CCD camera is placed at the image plane to capture the projection data. The visual interface controls and synchronizes projection data acquisition with the step motor. Then, to obtain an experimental IH-sinogram, a data line is selected from the actual image. This line defines an object slice and is fixed for all projections. This data row is placed on a new image as the first row to begin to build the IH-sinogram. The step motor turns a step, the next projection angle is generated, and a new image is obtained; then, the aforementioned procedure is applied to obtain the second row at the IH-sinogram. This process is repeated until the object has turned 360°, and therefore an IH-sinogram of the 200 projections is generated. During reconstruction, due to the symmetry property of parallel projections for the reconstructions only are used the projections in the range $(0 °, 180 °)$ .

Two phase objects were used in this experiment: a rectangular slide piece of $1 \times 3 \times 20 m m^{3}$ in size, and a pipette of 5.6mm in external diameter and 0.25mm in wall thickness, both made of Pyrex glass with a nominal refraction index of 1.52; a photo of each object is depicted in Fig. 5a . Figure 5b shows the IH-sinograms constructed from a slice from each one of the objects; the slices are indicated on the objects in Fig. 5a with a white-dotted line, and Fig. 5c shows the corresponding tomographic reconstructions by using the back-projection algorithm.

Fig. 5 (color online) Experimental results, (a) rectangular slide glass block, (b) pipette.

Download Full Size | PDF

5. Conclusion and remarks

The irradiance at the image plane of a 4f _L optical imaging system using a phase step filter of π radians has been considered as the projection data in optical tomography of phase objects, and it has been mathematically proved that this irradiance complies with the symmetry property and with the zero moment of the RT. It was shown that, for the all possible projection angles, a modified sinogram can be obtained, and this is called irradiance-Hilbert-sinogram (IH-sinogram). As a consequence of using directly these IH-sinograms, the obtained reconstructions consist of images showing isotropic edge-enhancement for both numerical simulations and experimental results. Thus, the HT filtering approach does not only serve to detect phase projections, but it is also capable to render phase-edge enhanced tomographic images as an extra feature after a usual reconstruction using a routine algorithm. It is important to note that, for the phase objects in the numerical simulation, the obtained projection data were smaller or equal than a wavelength, so that the phase was minor or equal than $2 π$ radians. Therefore, the irradiance did not show interference fringes (Fig. 3-d1 and Fig. 4b), and this condition was sought to be met for the phase objects shown in Fig. 5a. However, some fringes can be observed in the experimental IH-sinograms, as illustrated in Fig. 5b. Maybe due to these diffraction fringes, some discontinuities in the reconstructed images can be observed (Fig. 5c). For thick-phase objects, the interference pattern will be inevitably present, and a new method will be necessary. On the other hand, when thin-phase objects are considered, the transmittance function can be approximated to the projection data, and then a new method will also be necessary. We have taken these two considerations for future work.

Acknowledgments

One of the authors (AMP) appreciates CONACYT’s grant under number 160260/160260. This work was partially supported by PROMEP under grant PROMEP/103.5/09/4544.

References and links

1. S. R. Deans, “The Radon Transform and Some of its Applications,” (Wiley, New York. 1983).

2. A. F. Fercher, W. Drexler, C. K. Hitzenberger, and T. Lasser, “Optical coherence tomography—principles and applications,” Rep. Prog. Phys. 66(2), 239–303 (2003). [CrossRef]

3. G. P. Montgomery Jr and D. L. Reuss, “Effects of refraction on axisymmetric flame temperatures measured by holographic interferometry,” Appl. Opt. 21(8), 1373–1380 (1982). [CrossRef] [PubMed]

4. J. M. Mehta and W. Z. Black, “Errors associated with interferometric measurement of convective heat transfer coefficients,” Appl. Opt. 16(6), 1720–1726 (1977). [CrossRef] [PubMed]

5. J. M. Mehta and W. M. Worek, “Analysis of refraction errors for interferometric measurements in multicomponent systems,” Appl. Opt. 23(6), 928–933 (1984). [CrossRef] [PubMed]

6. S. Cha and C. M. Vest, “Tomographic reconstruction of strongly refracting fields and its application to interferometric measurement of boundary layers,” Appl. Opt. 20(16), 2787–2794 (1981). [CrossRef] [PubMed]

7. C. Meneses-Fabian, G. Rodriguez-Zurita, and V. Arrizón, “Optical tomography of transparent objects with phase-shifting interferometry and stepwise-shifted Ronchi ruling,” J. Opt. Soc. Am. A 23(2), 298–305 (2006). [CrossRef]

8. C. Meneses-Fabian, G. Rodriguez-Zurita, R. Rodriguez-Vera, and F. Jose, “Optical tomography with parallel projection differences and Electronic Speckle Pattern Interferometry,” Opt. Commun. 228(4-6), 201–210 (2003). [CrossRef]

9. G. Rodríguez-Zurita, C. Meneses-Fabián, J.-S. Pérez-Huerta, and J.-F. Vázquez-Castillo, ““Tomographic directional derivative of phase objects slices using 1-D derivative spatial filtering of fractional order ½,” ICO20,” Proc. SPIE 6027, 410–416 (2006).

10. J. A. Davis, D. E. McNamara, and D. M. Cottrell, “Analysis of the fractional hilbert transform,” Appl. Opt. 37(29), 6911–6913 (1998). [CrossRef]

11. J. A. Davis, D. E. McNamara, D. M. Cottrell, and J. Campos, “Image processing with the radial Hilbert transform: theory and experiments,” Opt. Lett. 25(2), 99–101 (2000). [CrossRef]

12. J. A. Davis and M. D. Nowak, “Selective edge enhancement of images with an acousto-optic light modulator,” Appl. Opt. 41(23), 4835–4839 (2002). [CrossRef] [PubMed]

13. J. Hsieh, “Computed Tomography: principles, design, artifacts, and recent advances,” (SPIE PRESS, Bellingham, Washington USA, 2003).

Isotropic edge-enhancement by the Hilbert-transform in optical tomography of phase objects

Abstract

1. Introduction

2. Basic considerations

3. Numerical Simulation

4. Experimental results

5. Conclusion and remarks

Acknowledgments

References and links

Cited By

Figures (5)

Equations (10)

Optics Express