An advanced algorithm for dispersion encoded full range frequency domain optical coherence tomography

doi:10.1364/OE.20.024925

« Show journal navigation

An advanced algorithm for dispersion encoded full range frequency domain optical coherence tomography

Felix Köttig, Peter Cimalla, Maria Gärtner, and Edmund Koch »View Author Affiliations

Optics Express, Vol. 20, Issue 22, pp. 24925-24948 (2012)
http://dx.doi.org/10.1364/OE.20.024925

View Full Text Article

Acrobat PDF (4263 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

Article Navigation

▸ Article Outline
– Abstract
1. Introduction
2. Theory and DEFR implementation
3. Materials and methods
4. Results and discussion
5. Conclusion
– References and links

Article Tools
Full Text
Article Info
References
Cited By
Figures
Metrics
Download PDF

Viewing Equations in the
HTML Article

Optimal Viewing Recommendations

Full Text
Article Info
References (35)
Cited By
Figures (10)
Metrics

Abstract

Dispersion encoded full range (DEFR) optical coherence tomography (OCT) has become highly attractive as it is a simple way to increase the measurement range of OCT systems. Full range OCT is especially favorable as it does not only increase the measurement range but also shifts the highest sensitivity into the center of the measurement range. While the early versions of DEFR were highly computational expensive, new versions reduce the number of necessary Fourier transforms. Recently it has been shown that a GPU based algorithm can perform DEFR with more than 20,000 A-lines per second. We present a new version of the DEFR algorithm that requires only one Fourier transform per A-scan and uses convolution in z-space instead of multiplication in k-space, therefore reducing the computational effort considerably. While dispersion encoding has so far only been used to suppress mirror artifacts, we show that, with dispersion encoding and only one more Fourier transform, autocorrelation terms can be removed likewise. Since very high values of dispersion reduce the effective measurement range in dispersion encoded OCT, we present an estimate for a sufficient amount of dispersion for a successful image recovery, which is depending on the thickness of the scattering layers. Furthermore, we demonstrate the usability of ZnSe as a new dispersive material with a very high dispersion and describe a simple method to extract the dispersive phase from the measurement of a single reflex of a glass surface. Using a standard consumer PC, an artifact-free recovery of 1000 – 2000 A-scans per second with 2048 depth values including autocorrelation removal was achieved. The dynamic range (sensitivity) is not reduced and the suppression ratio of mirror artifacts and autocorrelation signals is more than 50dB using ZnSe.

1. Introduction

Optical coherence tomography (OCT) is a noninvasive, contactless, in vivo optical imaging technique which provides cross-sectional or three-dimensional images of transparent to highly scattering samples, for example biological tissue, with a high axial resolution of 1 – 15μm. OCT has found its main application in ophthalmology, but many other areas in biomedical research [1

S. Marschall, B. Sander, M. Mogensen, T. Jørgensen, and P. Andersen, “Optical coherence tomography – current technology and applications in clinical and biomedical research,” Anal. Bioanal. Chem. 400, 2699–2720 (2011). [CrossRef] [PubMed]

, 2

J. Walther, M. Gaertner, P. Cimalla, A. Burkhardt, L. Kirsten, S. Meissner, and E. Koch, “Optical coherence tomography in biomedical research,” Anal. Bioanal. Chem. 400, 2721–2743 (2011). [CrossRef] [PubMed]

] and nondestructive testing are emerging. During the 20 years of development, the speed of OCT has increased dramatically from about 30 A-scans per second in the early work [3

D. Huang, E. Swanson, C. Lin, J. Schuman, W. Stinson, W. Chang, M. Hee, T. Flotte, K. Gregory, and C. Puliafito “Optical coherence tomography,” Science 254, 1178–1181 (1991). [CrossRef] [PubMed]

] to more than 20 million A-scans per second using FDML lasers [4

W. Wieser, B. R. Biedermann, T. Klein, C. M. Eigenwillig, and R. Huber, “Multi-megahertz oct: High quality 3d imaging at 20 million a-scans and 4.5 gvoxels per second,” Opt. Express 18, 14685–14704 (2010). [CrossRef] [PubMed]

]. A breakthrough in speed was achieved with the swap from time domain OCT (TD OCT) to Fourier domain OCT (FD OCT) eliminating the need for a fast mechanical change of the reference arm length, with the additional benefit of an improved signal to noise ratio (SNR) [5

R. Leitgeb, C. Hitzenberger, and A. Fercher, “Performance of fourier domain vs. time domain optical coherence tomography,” Opt. Express 11, 889–894 (2003). [CrossRef] [PubMed]

–7

J. F. de Boer, B. Cense, B. H. Park, M. C. Pierce, G. J. Tearney, and B. E. Bouma, “Improved signal-to-noise ratio in spectral-domain compared with time-domain optical coherence tomography,” Opt. Lett. 28, 2067–2069 (2003). [CrossRef] [PubMed]

]. In TD OCT, the axial resolution is limited by the spectral width of the light source (including other components of the optical system) and the measurement range is limited by the stroke of the mirror in the reference arm, which are mainly independent from each other. In FD OCT, on the contrary, axial resolution and measurement range depend on the width or resolution, respectively, of the acquired wavelength spectrum and therefore are not independent from each other. Because some applications demand both higher resolution and greater measurement range, the requirements for the FD OCT spectrometer system increase and cannot be fulfilled in all cases.

The ratio of depth range to resolution is limited in classical FD OCT by the number of spectral points [2

]. While in spectrometer based FD OCT systems, called spectral domain OCT (SD OCT), this number is limited by the pixel number of the line scan camera, in swept source OCT (SS OCT) this is limited by the ratio of the sampling frequency of the analog-to-digital converter (ADC) to the repetition frequency (neglecting the limited duty cycle) of the swept source. As both values cannot simply be increased there is a demand for other methods to increase this ratio. Furthermore, both kinds of OCT systems show a sensitivity roll off to the end of the measurement range [8

T. Bajraszewski, M. Wojtkowski, M. Szkulmowski, A. Szkulmowska, R. Huber, and A. Kowalczyk, “Improved spectral optical coherence tomography using optical frequency comb,” Opt. Express 16, 4163–4176 (2008). [CrossRef] [PubMed]

, 9

M. Hagen-Eggert, P. Koch, and G. Hüttmann, “Analysis of the signal fall-off in spectral domain optical coherence tomography systems,” in Optical Coherence Tomography and Coherence Domain Optical Methods in Biomedicine XVI , Joseph A. Izatt, James G. Fujimoto, and Valery V. Tuchin, eds., Proc. SPIE 8213, 82131K (2012).

] having therefore a reduced sensitivity in the mostly interesting center of the image. Consequently, several methods to increase the measurement range have been designed. Full range OCT eliminates the complex conjugate artifacts to effectively double the measurement range around the zero delay position (ZDP), which is in most systems the position of highest sensitivity. The removal of the complex conjugate artifacts is often performed by acquiring a second spectrum at almost the same beam position and introducing a phase shift between both A-scans. This can be done by special interferometers providing quadrature signals [10

Y. Mao, S. Sherif, C. Flueraru, and S. Chang, “3x3 mach-zehnder interferometer with unbalanced differential detection for full-range swept-source optical coherence tomography,” Appl. Opt. 47, 2004–2010 (2008). [CrossRef] [PubMed]

–12

P. Meemon, K.-S. Lee, and J. P. Rolland, “Doppler imaging with dual-detection full-range frequency domain optical coherence tomography,” Biomed. Opt. Express 1, 537–552 (2010). [CrossRef]

] or by introducing a phase shift between successive A-scans [13

T. Wu, Z. Ding, C. Wang, and M. Chen, “Full-range swept source optical coherence tomography based on carrier frequency by transmissive dispersive optical delay line,” J. Biomed. Opt. 16, 126008 (2011). [CrossRef] [PubMed]

]. The latter can be done with a phase modulator in the reference arm [14

D. Y. Kim, J. S. Werner, and R. J. Zawadzki, “Comparison of phase-shifting techniques for in vivo full-range, high-speed fourier-domain optical coherence tomography,” J. Biomed. Opt. 15, 056011 (2010). [CrossRef] [PubMed]

, 15

E. Götzinger, M. Pircher, R. Leitgeb, and C. Hitzenberger, “High speed full range complex spectral domain optical coherence tomography,” Opt. Express 13, 583–594 (2005). [CrossRef] [PubMed]

] or by changing the length of the sample arm. The length of the sample arm can easily be changed by pivoting the rotation of the deflecting scanning mirror (BM-scanning [16

F. Jaillon, S. Makita, M. Yabusaki, and Y. Yasuno, “Parabolic bm-scan technique for full range doppler spectral domain optical coherence tomography,” Opt. Express 18, 1358–1372 (2010). [CrossRef] [PubMed]

–19

L. An and R. K. Wang, “Use of a scanner to modulate spatial interferograms for in vivo full-range fourier-domain optical coherence tomography,” Opt. Lett. 32, 3423–3425 (2007). [CrossRef] [PubMed]

]). In these systems dense and, in the case of moving samples, fast sampling is needed for high efficiency of complex conjugate removal. While two adjacent scans allow the reduction of the mirror images, better results are achieved by taking up to five A-scans [20

C.-T. Wu, T.-T. Chi, Y.-W. Kiang, and C. C. Yang, “Computation time-saving mirror image suppression method in fourier-domain optical coherence tomography,” Opt. Express 20, 8270–8283 (2012). [CrossRef] [PubMed]

, 21

M. Wojtkowski, A. Kowalczyk, R. Leitgeb, and A. F. Fercher, “Full range complex spectral optical coherence tomography technique in eye imaging,” Opt. Lett. 27, 1415–1417 (2002). [CrossRef]

]. Furthermore, high phase stability is essential for the algorithm. In contrast to the previously described systems, dispersion encoded full range (DEFR) OCT uses only one spectrum to calculate the full range depth signal [22

B. Hermann, B. Hofer, C. Meier, and W. Drexler, “Spectroscopic measurements with dispersion encoded full range frequency domain optical coherence tomography in single- and multilayered non- scattering phantoms,” Opt. Express 17, 24162–24174 (2009). [CrossRef]

–26

L. Wang, B. Hofer, J. A. Guggenheim, and B. Povazay, “Graphics processing unit-based dispersion encoded full-range frequency-domain optical coherence tomography,” J. Biomed. Opt. 17, 077007 (2012). [CrossRef] [PubMed]

]. While it has been demonstrated that DEFR removes the complex conjugate data very efficient in many samples, there is a more philosophical problem with the algorithm that from only N (real) spectral data points N complex valued depth values can be retrieved. As we will show later, the solution for this dilemma is that the algorithm performs only well if many of the depth values are zero or small.

The basic idea of DEFR OCT is to exploit a dispersion mismatch between the two interferometer arms, caused by a dispersive material in one arm, which leads to a broadening of signal peaks in z-space. This can be corrected by numerical dispersion compensation in k-space before the inverse Fourier transform, thereby restoring the true signal components and broadening their complex conjugate mirror artifacts even more. A peak detector reveals the true signal components and their mirror artifacts can be removed consequently. The previously described versions of the DEFR algorithm introduce one iterative step into the common OCT processing, which comprises two steps. In a first step the spectral data are corrected for detector background then linearized in k-space and spectrally shaped. In case of dispersion mismatch between the sample and reference arm a corresponding phase shift is introduced to restore the resolution [27

M. Wojtkowski, V. Srinivasan, T. Ko, J. Fujimoto, A. Kowalczyk, and J. Duker, “Ultrahigh-resolution, high-speed, fourier domain optical coherence tomography and methods for dispersion compensation,” Opt. Express 12, 2404–2422 (2004). [CrossRef] [PubMed]

, 28

B. Cense, N. Nassif, T. Chen, M. Pierce, S.-H. Yun, B. Park, B. Bouma, G. Tearney, and J. de Boer, “Ultrahigh-resolution high-speed retinal imaging using spectral-domain optical coherence tomography,” Opt. Express 12, 2435–2447 (2004). [CrossRef] [PubMed]

]. In the second step of the usual FD OCT processing the data are then Fourier transformed and displayed using only the amplitudes and logarithmic scaling. The DEFR algorithm is introduced iteratively in between these two processing steps by identifying the highest signals in the depth data and removing them from the spectral data. This takes two or even four Fourier transforms in each iteration step. While in the early versions only one component was removed in each step, what could result in as many steps as sample data points, the recent versions remove many components in one step. A high dispersion between sample and reference arm ensures a large amplitude ratio between the true signals and the complex conjugate artifacts thereby reducing the number of necessary iteration steps. The iteration can be stopped if all the complex conjugate artifacts are below the noise level. This grouping of data reduces the processing time by about two orders of magnitude to the early versions of the DEFR algorithm. Witte et al. [29

S. Witte, M. Baclayon, E. J. Peterman, R. F. Toonen, H. D. Mansvelder, and M. L. Groot, “Single-shot two-dimensional full-range optical coherence tomography achieved by dispersion control,” Opt. Express 17, 11335–11349 (2009). [CrossRef] [PubMed]

] suggested for a system with high dispersion only one iteration, limiting the number of Fourier transforms to only three, achieving more than 40dB of mirror artifact suppression, which is sufficient for the described system with a dynamic range of 45dB. Recently the DEFR algorithm was implemented on a GPU using the parallel computing capabilities of these devices. By this, the processing time could be reduced by approximately two orders of magnitude again [26

The main processing steps of the DEFR algorithm, proposed by Hofer et al. [23

B. Hofer, B. Považay, B. Hermann, A. Unterhuber, G. Matz, and W. Drexler, “Dispersion encoded full range frequency domain optical coherence tomography,” Opt. Express 17, 7–24 (2009). [CrossRef] [PubMed]

], are depicted in Fig. 1. The iterative DEFR algorithm requires two Fourier transforms for the subtraction of a single signal component. That is in each iteration, two Fourier transforms are needed to calculate from z-space to the complex conjugate spectrum in z-space again after application of a phase shift in k-space:

ℱ^{- 1} (ℱ (c (z)) e^{\pm i 2 ϕ (k)}) (z),

(1)

where c(z) denotes the complex spectrum in z-space and e^±i2^ϕ⁽^k⁾ is twice the (inverse) dispersive phase function caused by the dispersive material in the reference arm of the interferometer. Additionally, two further Fourier transforms are required to calculate from the complex conjugate spectrum in z-space to the ’original’ spectrum in z-space after applying the opposite phase shift in k-space. The main idea behind our proposed accelerated algorithm is that the phase shift, which is just a multiplication with a complex function in k-space, can be expressed as a convolution in z-space:

ℱ^{- 1} (ℱ (c (z) e^{\pm i 2 ϕ (k)}) (z) = (c * ℱ^{- 1} (e^{\pm i 2 ϕ (k)})) (z) .

(2)

This convolution of a single spectral component, which is used for subtraction of its corresponding mirror component, and the Fourier transform of the dispersive phase function is computationally efficient because the latter needs to be calculated only once and the signal component, used for subtraction, is just a single component in z-space. Assuming the spectrum c(z) contains only one component that differs from zero (c(z) = c̃) (which is the case when the subtraction of a mirror component of a true signal is performed), the convolution reduces to a simple multiplication:

(\tilde{c} * ℱ^{- 1} (e^{\pm i 2 ϕ (k)})) (z) = \tilde{c} \times ℱ^{- 1} (e^{\pm i 2 ϕ (k)}) (z - s),

(3)

where s denotes the position of c̃ in the complex spectrum in z-space. The details of this idea and improvements to reduce the autocorrelation terms in the same iteration will be explained in the following chapter.

Fig. 1 Iteration procedure in the DEFR algorithm. A peak detector finds the strongest signal component in z-space (Fig. 1(a)) and removes it from the spectrum. In order to remove its mirror artifacts, twice the dispersive phase is applied on the signal in k-space (Fig. 1(b)), which is calculated via Fourier transform. The resulting signal after the nonlinear phase shift is depicted in Fig. 1(c). The mirror component, which is just one peak in z-space after inverse Fourier transform (Fig. 1(d)), can be finally removed via subtraction of the peak corresponding to the true signal component found in Fig. 1(a) (these two signals are symmetric around the zero delay position). Again, the nonlinear phase shift in k-space is applied after a Fourier transform and the ’original’ spectrum in z-space is retrieved via inverse Fourier transform.

2. Theory and DEFR implementation

2.1. OCT signal with dispersion and autocorrelation terms

To understand the processing steps of our proposed DEFR algorithm, it is helpful to initially derive a detailed signal model of the OCT signal with dispersion and autocorrelation terms. We will discuss the contributions of these signals to the spectrum and exploit their different properties regarding the influence of the dispersion. This model is based on the general model derived in [23

, 30

W. Drexler and J. Fujimoto, Optical Coherence Tomography: Technology and Applications (Springer, 2008). [CrossRef]

], but includes the crucial terms introduced by dispersion. To evaluate the influence of an artificial dispersion mismatch, we consider a Michelson interferometer with dispersive material in the reference arm. The electric field of the incident polychromatic electromagnetic wave with the spectral amplitude of the electric field s(k) can be written as:

E_{i} (k, t) = s (k) e^{i (k z - ω t)} {\vec{e}}_{ρ},

(4)

with the propagation constant

k = \frac{ω}{c_{0}}

(equal to the real part of the angular wavenumber), the angular frequency ω and the vacuum speed of light

c_{0} = \frac{1}{\sqrt{μ_{0} ε_{0}}}

. This is a plane wave fulfilling the Helmholtz Eq. ∇²E + k²E = 0 for a harmonic excitation. When passing the beam splitter the incident electric field E_i divides into two plane waves, one in the reference arm and one in the sample arm of the interferometer. Without any dispersion the electric field after reflection at the reference mirror would be given by:

E_{r} (k, t) = E_{i} (k, t) r_{M} e^{i 2 k z_{R}},

(5)

where r_M is the reflection coefficient of the reference mirror and z_R is the length of the reference arm. The splitting ratio of the beam splitter is incorporated in the reflection coefficients and is assumed to be independent from frequency. Using a glass plate to introduce dispersion into the reference arm, the electric field after reflection at the reference mirror with dispersion reads:

E_{r, dis} (k, t) = E_{i} (k, t) r_{M} t_{G} e^{i 2 (k z_{F} + k (ω) z_{G})},

(6)

where the transmission coefficient t_G incorporates the losses due to the reflections at the surface of the glass plate and z_F and z_G are the geometrical lengths in air and glass respectively. The interferences between the light backscattered from the glass surface and from the sample can be neglected because their optical path differences and thus their modulation frequencies strongly deviate from the wanted interferences and are therefore not displayed. Note that the propagation constant

k = \frac{ω}{c_{0}}

has a linear dependence of the frequency, whereas the propagation constant k (ω) can be written as a Taylor series expansion around the center frequency ω₀ of the spectral amplitude of the electric field s(k) [31

A. Fercher, C. Hitzenberger, M. Sticker, R. Zawadzki, B. Karamata, and T. Lasser, “Numerical dispersion compensation for partial coherence interferometry and optical coherence tomography,” Opt. Express 9, 610–615 (2001). [CrossRef] [PubMed]

k (ω) = \sum_{l = 0}^{\infty} a_{l} {(ω - ω_{0})}^{l},

(7)

with nonvanishing higher order coefficients a_l that are determined as the derivatives of k (ω) with respect to frequency:

a_{l} = \frac{1}{l!} {\frac{d^{l} k (ω)}{d ω^{l}} |}_{ω = ω_{0}} .

(8)

To expand the power of the binomial (ω – ω₀)^l, we consider the binomial series identity:

{(ω - ω_{0})}^{l} = \sum_{m = 0}^{l} (\begin{matrix} l \\ m \end{matrix}) ω^{m} {(- ω_{0})}^{l - m},

(9)

where

(\begin{matrix} l \\ m \end{matrix})

is a binomial coefficient. There are three fundamental components, corresponding to m = 0, m = 1 and m ≥ 2, that can be identified in Eq. (7). The first component is not frequency dependent and is obtained using Eq. (7) and Eq. (9) as follows:

k_{abs} = \sum_{l = 0}^{\infty} a_{l} (\begin{matrix} l \\ 0 \end{matrix}) ω^{0} {(- ω_{0})}^{l - 0} = \sum_{l = 0}^{\infty} a_{l} {(- ω_{0})}^{l} .

(10)

It causes a constant phase shift which has no impact on the width of the peaks in z-space. The second component contains the linear frequency ω¹ and is calculated the same way:

k_{lin} = \sum_{l = 1}^{\infty} a_{l} (\begin{matrix} l \\ 1 \end{matrix}) ω^{1} {(- ω_{0})}^{l - 1} = \sum_{l = 1}^{\infty} a_{l} l ω {(- ω_{0})}^{l - 1} .

(11)

By this component the optical path length of the glass plate in the reference arm is determined from the geometrical path length:

z_{G, opt} = z_{G} \frac{c_{0} k_{lin}}{ω} = z_{G} c_{0} \sum_{l = 1}^{\infty} a_{l} l {(- ω_{0})}^{l - 1} .

(12)

The third component finally includes all higher powers of ω(ω^l, l = 2. . .∞):

k_{dis} = \sum_{l = 2}^{\infty} a_{l} \sum_{m = 2}^{l} (\begin{matrix} l \\ m \end{matrix}) ω^{m} {(- ω_{0})}^{l - m},

(13)

that is what we usually refer to as dispersion and causes a broadening of the peaks in z-space. Introducing the abbreviations:

ϕ_{0} = 2 k_{abs} z_{G} = 2 \sum_{l = 0}^{\infty} a_{l} {(- ω_{0})}^{l} z_{G},

(14)

which comprises the non frequency dependent component from Eq. (10) and thus causes no dispersion, and the dispersive phase term:

ϕ (k) = 2 k_{dis} z_{G} = 2 \sum_{l = 2}^{\infty} a_{l} \sum_{m = 2}^{l} (\begin{matrix} l \\ m \end{matrix}) ω^{m} {(- ω_{0})}^{l - m} z_{G},

(15)

the reflected electric field E_r,dis (Eq. (6)) can be rewritten as:

E_{r, dis} (k, t) = E_{i} (k, t) r_{M} t_{G} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G})} e^{i ϕ (k)} .

(16)

The electric field, backscattered from the sample with reflecting layers at positions z_S_,l and reflection coefficients r_S,_l, can be determined as the sum of single waves:

E_{s} (k, t) = E_{i} (k, t) \sum_{l} r_{S, l} e^{i 2 k n_{ref} z_{S, l}} .

(17)

Note that the (depth dependent) refractive index n_ref is frequency dependent if the sample consists of dispersive material, thus dispersion compensation will be necessary to avoid blurring effects. The splitting ratio of the beam splitter is again incorporated in the reflection coefficients of the sample. The superposition of the electric fields from the reference arm E_r,dis and from the sample arm E_s of the interferometer finally results in the interference signal:

E_{int} (k, t) = E_{r, dis} (k, t) + E_{s} (k, t) = E_{i} (k, t) (r_{M} t_{G} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G})} e^{i ϕ (k)} + \sum_{l} r_{S, l} e^{i 2 k n_{ref} z_{S, l}})

(18)

The intensity of light as a function of frequency is proportional to the time averaged square of the electric field:

p (k) \propto {〈 {| E (k, t) |}^{2} 〉}_{t},

(19)

P_{int} (k) = P_{i} (k) {| r_{M} t_{G} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G})} e^{i ϕ (k)} + \sum_{l} r_{S, l} e^{i 2 k n_{ref} z_{S, l}} |}^{2} .

(20)

Expanding Eq. (20) using the relation |x|² = xx^*, the different contributions to the power spectrum can be identified as follows:

P_{int} (k) = P_{i} (k) [r_{M}^{2} t_{G}^{2} + \sum_{l} r_{S, l}^{2} +]

(21a)

\sum_{l} r_{M} t_{G} r_{S, l} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, l})} e^{i ϕ (k)} +

(21b)

\sum_{l} r_{M} t_{G} r_{S, l} e^{- i ϕ_{0}} e^{- i 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, l})} e^{- i ϕ (k)} +

(21c)

\sum_{l \neq m} r_{S, l} r_{S, m} (e^{i 2 k n_{ref} (z_{S, l} - z_{S, m})} + e^{- i 2 k n_{ref} (z_{S, l} - z_{S, m})})] .

(21d)

Equation (21), term (a) describes the DC component, determined by the power spectrum of the light source and the reflection and transmission coefficients. Equation (21), term (b) and Eq. (21), term (c) describe the cross-correlation components that means the ’true’ signal we want to retrieve from the spectrum. Equation (21), term (d) are autocorrelation signals, which are usually considered as artifacts in conventional half range FD OCT systems. As the reflection coefficients of the sample r_S,_l are usually very low compared to that one of the reference mirror r_M, the autocorrelation terms are negligible in many samples. Nevertheless, they can have a disturbing impact on images of highly scattering samples as we will demonstrate later. Note that the autocorrelation signals are not corrupted by the artificial dispersion, which will later enable us to remove them from the spectrum. Numerical dispersion compensation is achieved by multiplying Eq. (21) with the compensating phase term e⁻ⁱ^ϕ⁽^k⁾:

P_{int} (k) e^{- i ϕ (k)} = P_{i} (k) [(r_{M}^{2} t_{G}^{2} + \sum_{l} r_{S, l}^{2}) e^{- i ϕ (k)} +

(22a)

\sum_{l} r_{M} t_{G} r_{S, l} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, l})} +

(22b)

\sum_{l} r_{M} t_{G} r_{S, l} e^{- i ϕ_{0}} e^{- i 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, l})} e^{- i 2 ϕ (k)} +

(22c)

(\sum_{l \neq m} r_{S, l} r_{S, m} (e^{i 2 k n_{ref} (z_{S, l} - z_{S, m})} + e^{- i 2 k n_{ref} (z_{S, l} - z_{S, m})})) e^{- i ϕ (k)}]

(22d)

In order to retrieve the depth information, we need to calculate the inverse Fourier transform of Eq. (22):

ℱ^{- 1} (P_{int} (k) e^{- i ϕ (k)}) (z) = {ℱ^{- 1} (P_{i} (k)) * [(r_{M}^{2} t_{G}^{2} + \sum_{l} r_{S, l}^{2}) ℱ^{- 1} (e^{- i ϕ (k)}) +

(23a)

\sum_{l} r_{M} t_{G} r_{S, l} e^{i ϕ_{0}} δ (z - 2 (z_{F} + k_{lin} c_{0} z_{G} - n_{ref} z_{S, l})) +

(23b)

\sum_{l} r_{M} t_{G} r_{S, l} e^{- i ϕ_{0}} δ (z + 2 (z_{F} + k_{lin} c_{0} z_{G} - n_{ref} z_{S, l})) * ℱ^{- 1} (e^{- i 2 ϕ (k)}) +

(23c)

\sum_{l \neq m} r_{S, l} r_{S, m} δ (z \pm 2 n_{ref} (z_{S, l} - z_{S, m})) * ℱ^{- 1} (e^{- i ϕ (k)})]} (z),

(23d)

where the transform variable is

k = \frac{ω}{c_{0}}

and the inverse Fourier transform is performed as follows:

ℱ^{- 1} (x (k)) (z) = \int_{- \infty}^{\infty} x (k) e^{i 2 π k z} d k .

(24)

Thus, a plane wave in air with the propagation constant

k = \frac{ω}{c_{0}}

will be transformed into a single peak, whereas a frequency dependent speed of light c = c(ω) leads to dispersion broadening the transformed peak by a convolution with the inverse Fourier transform of the dispersive phase function ℱ⁻¹(eⁱ^ϕ⁽^k⁾). The dispersion compensated spectrum Eq. (23) consists of the desired complex full range signal Eq. (23), term (b) as well as of its complex conjugate mirror components spread by convolution with the inverse Fourier transform of twice the inverse dispersive phase function ℱ⁻¹(e⁻ⁱ²^ϕ⁽^k⁾) Eq. (23), term (c) and DC and autocorrelation components, Eq. (23), term (a) and Eq. (23), term (d), respectively, which are also spread by the numerical dispersion compensation. In order to remove autocorrelation artifacts from the spectrum, we need to consider both the power spectrum with and without numerical dispersion compensation, so we also need to calculate the inverse Fourier transform of the original, non dispersion compensated spectrum Eq. (21):

ℱ^{- 1} (P_{int} (k)) (z) = {ℱ^{- 1} (P_{i} (k)) * [(r_{M}^{2} t_{G}^{2} + \sum_{l} r_{S, l}^{2}) +

(25a)

\sum_{l} r_{M} t_{G} r_{S, l} e^{i ϕ_{0}} δ (z - 2 (z_{F} + k_{lin} c_{0} z_{G} - n_{ref} z_{S, l})) * ℱ^{- 1} (e^{i ϕ (k)}) +

(25b)

\sum_{l} r_{M} t_{G} r_{S, l} e^{- i ϕ_{0}} δ (z + 2 (z_{F} + k_{lin} c_{0} z_{G} - n_{ref} z_{S, l})) * ℱ^{- 1} (e^{- i ϕ (k)}) +

(25c)

\sum_{m \neq l} r_{S, l} r_{S, m} δ (z \pm 2 n_{ref} (z_{S, l} - z_{S, m}))]} (z) .

(25d)

With the signals derived in Eq. (23) and Eq. (25), we can now develop an implementation of the DEFR algorithm that enables us to do both a direct subtraction of the complex conjugate mirror components of a true signal in z-space and to perform an effective autocorrelation suppression. Moreover, Eq. (23) and Eq. (25) mathematically prove that a subtraction of signal components directly in z-space using the corresponding dispersive phase functions equals the subtraction indirectly via k-space as proposed in the former implementations of the DEFR algorithm. The details of our approach will be explained in the following chapter using a block diagram and pseudocode.

2.2. Description of the simple DEFR algorithm without autocorrelation removal

Compared with the proposed DEFR algorithm from Hofer et al. [23

, 26

], we implement a very efficient method to remove the complex conjugate mirror components caused by each true signal. Our improved implementation of the DEFR algorithm without autocorrelation removal can be derived directly from Eq. (23) when neglecting DC and autocorrelation terms, Eq. (23), term (a) and Eq. (23), term (d), respectively. Although these components are broadened by the dispersion compensation, which leads to a convolution with the inverse Fourier transform of the inverse dispersive phase function ℱ⁻¹(e⁻ⁱ^ϕ⁽^k⁾) (z) in z-space, resulting in a lower amplitude according to Parseval’s theorem:

\int_{- \infty}^{\infty} {| x (t) |}^{2} d t = \int_{- \infty}^{\infty} {| X (f) |}^{2} d f, where X (f) = ℱ (x (t)) (f),

(26)

DC and autocorrelation terms can have a significant impact on the reconstructed spectrum as we will show with images of an optical filter. While the first term in Eq. (23), term (a) (

r_{M}^{2} t_{G}^{2}

) is removed via background subtraction, the other terms in Eq. (23), term (a) (scaling with the reflection coefficients of the sample) and Eq. (23), term (d) (autocorrelation terms) cannot be distinguished from the wanted cross-correlation signals Eq. (23), term (b). Nevertheless, we will start with a simple implementation of our proposed DEFR algorithm neglecting DC and autocorrelation terms coming from the sample.

Figure 2 shows the block diagram of our proposed DEFR algorithm without autocorrelation removal concentrating on the iteration procedure itself, which is our basic approach of subtracting the complex conjugate mirror components of a true signal directly in z-space. In other implementations of the DEFR algorithm the subtraction of the complex conjugate mirror components of a true signal is performed by transforming from z-space to the complex conjugate spectrum in z-space again after applying a phase shift k-space [23

, 26

, 29

]. In Eq. (23), term (c) the contribution of the complex conjugate mirror components of a true signal found in Eq. (23), term (b) can be directly identified and thus removed from the spectrum by a simple subtraction. Step (1) comprises the calculation of the corresponding dispersive phase function via inverse Fourier transform. This calculation has to be done only once and the results can be saved as kind of calibration for each setup with a certain fixed dispersion. We will present a simple method of calculating the dispersive phase ϕ(k) from a single measurement in Sect. 3.2. In step (2) the inverse Fourier transform of the background corrected, in frequency resampled, dispersion compensated signal is calculated (

c_{1}^{i} (n)

, note that the index 1 is used to differentiate from a second spectrum in our advanced implementation of the DEFR algorithm; i. . .iteration index). The output for the true signal components is denoted

d_{1}^{i} (n)

. In frequency domain OCT one inverse Fourier transform is anyway necessary for each A-scan so there is no additional (inverse) Fourier transform required for our proposed algorithm. The iterative process consists of step (3) and step (4). A peak detector finds the strongest signal component (highest amplitude) and provides its value

c_{1}^{i} (n_{1}^{i})

and position in the spectrum in z-space

n_{1}^{i}

, respectively. When comparing Eq. (23), term (b) and Eq. (23), term (c), we can distinguish between the true signal component and its complex conjugate mirror component spread by twice the inverse dispersive phase. According to Parseval’s theorem, the energy of a signal in z-space remains constant if it is multiplied with the inverse dispersive phase function e⁻ⁱ^ϕ⁽^k⁾ in k-space before applying the inverse Fourier transform as the magnitude of this function equals one (|e⁻ⁱ^ϕ⁽^k⁾| = 1). As a narrow peak in Eq. (23), term (b) is broadened by twice the inverse dispersive phase in Eq. (23), term (c) and the signal energy in Eq. (23), term (b) and Eq. (23), term (c) is the same, the amplitude in Eq. (23), term (c) must be smaller by a factor

\approx \sqrt{2 D}

if the peak is one channel wide in Eq. (23), term (b) and 2D channels wide in Eq. (23), term (c), assuming a nearly rectangular shape of the inverse Fourier transform of the dispersive phase function, which is realistic as we will show later (note that the energy of a signal is proportional to the square of the amplitude). The true signal component

c_{1}^{i} (n_{1}^{i})

, found by the peak detector, is added to the output

d_{1}^{i} (n)

and both the true signal component

c_{1}^{i} (n_{1}^{i})

and its corresponding complex conjugate mirror components in

c_{1}^{i} (n)

are subtracted (directly in z-space) using

c_{1}^{i} (n_{1}^{i})

and

n_{1}^{i}

. The iteration stops if a maximum number of iterations is reached or if the amplitude of the true peak falls below a predefined threshold, which results in only noise remaining in the spectrum

c_{1}^{i} (n_{1}^{i})

(i.e. consider the complex conjugate mirror components of the true signal, if they are below the noise level their removal will not improve the image quality and the iteration can be stopped). The iteration process does not need any additional Fourier transforms, therefore the computational complexity is very low – one multiplication and one subtraction for each array element of

c_{1}^{i} (n)

in each iteration. Moreover, considering only that part of the dispersive phase function that significantly differs from zero (or whose amplitude, when multiplied with a strong peak, is significantly above the noise level) (for multiplication and subtraction), further reduces the number of computations. Optionally, the residual spectrum

c_{1}^{i} (n)

(remaining after the iteration was stopped) can be added to the output

d_{1}^{i} (n)

in order to maintain weak signals.

Fig. 2 Simple implementation of the DEFR algorithm without autocorrelation removal. Independent from the actual processing steps, the inverse Fourier transform (IFFT) of twice the inverse dispersive phase function e^−i2ϕ(k_n) is calculated (denoted as

p_{2}^{-} (n)

) for the purpose of complex conjugate artifact removal (step (1)). After performing the numerical dispersion compensation in k-space, the depth signal in z-space is calculated as usual for each A-scan via inverse Fourier transform (step (2); the detector signal, resampled in frequency and background corrected, and the initial depth signal in z-space are denoted as P_int (k_n) and

c_{1}^{0} (n)

, respectively). From now on, all processing steps will be executed directly in z-space without the need of any additional Fourier transforms. The iterations are comprised of step (3) and step (4). In step (3) a peak detector finds the strongest signal component

c_{1}^{i} (n_{1}^{i})

, whose mirror component is subtracted directly in z-space using the information provided by Eq. (23), term (b) and Eq. (23), term (c). The true signal component

c_{1}^{i} (n_{1}^{i})

is added to the output

d_{1}^{i} (n)

and its complex value is also subtracted from the spectrum. The iteration finishes if a maximum number of iterations is reached or the residual spectrum

c_{1}^{i} (n)

contains only noise. Optionally, the residual spectrum

c_{1}^{i} (n)

can be added to the output

d_{1}^{i} (n)

2.3. Description of the advanced DEFR algorithm with autocorrelation removal

After having described a simple implementation of the DEFR algorithm, we now want to present our advanced approach that effectively eliminates autocorrelation artifacts from the spectrum with only minor additional computational effort. The algorithm proceeds similar to that described in Sect. 2.2, but additionally exploits the properties of the autocorrelation signals in Eq. (23), term (d) and Eq. (25), term (d). Consequently, we can reconstruct both the complex full range signal and the autocorrelation signal. While other implementations of the DEFR algorithm do not consider autocorrelation signals or simply set the signals near the zero delay position (which contain in most cases DC and autocorrelation terms) to zero [23

, 29

], thus removing any true signal in this region likewise and causing artifacts, our algorithm effectively removes them from the spectrum without disturbing true signal components in this region. This is very beneficial as autocorrelation removal is essential in some samples with highly scattering layers as we will show later.

Figure 3 shows the block diagram of our advanced DEFR algorithm including autocorrelation detection and removal. Step (1) again comprises the calculation of the corresponding dispersive phase functions via inverse Fourier transform. Besides the already in Sect. 2.2 used phase function e^−i2ϕ(k_n), we now also need to calculate the inverse Fourier transform of the phase functions e^±iϕ(k_n) (compare Eq. (23), term (a) and Eq. (23), term (d) as well as Eq. (25), term (b) and Eq. (25), term (c)). In step (2) the inverse Fourier transform of both the dispersion compensated (upper path) and the non dispersion compensated (lower path) signal is calculated (

c_{1}^{i} (n)

and

c_{2}^{i} (n)

respectively). To reduce the computational effort only one half of the spectrum for autocorrelation suppression is needed as this spectrum from a real valued signal is complex conjugate (the magnitude spectrum is symmetric around the zero delay position) and the other half does not contain any additional information. The outputs for the true signal components and autocorrelation components are denoted

d_{1}^{i} (n)

and

d_{2}^{i} (n)

, respectively (i. . .iteration index). The iterative process now consists of step (3) to step (5). A peak detector finds the strongest signal components in both spectra and provides their values

c_{1}^{i} (n_{1}^{i})

c_{2}^{i} (n_{2}^{i})

and positions

n_{1}^{i}

n_{2}^{i}

in the spectra in z-space respectively. In addition to the removal of complex conjugate mirror components of a true signal described in Sect. 2.2, the autocorrelation removal exploits the different properties of DC and autocorrelation signals in the two spectra

c_{1}^{i} (n)

and

c_{2}^{i} (n)

. When comparing Eq. (23), term (b) and Eq. (25), term (b) it is obvious that a true signal component found in

c_{1}^{i} (n)

is spread by the dispersion in Eq. (25), term (b) and results in a lower amplitude of the corresponding signal components in

c_{2}^{i} (n)

. Thus a real signal component can be distinguished from its corresponding contribution in

c_{2}^{i} (n)

. The distinction of autocorrelation signals from real signals is performed similar. By comparing Eq. (23), term (d) and Eq. (25), term (d) we can detect autocorrelation signals as their amplitude in

c_{2}^{i} (n)

is higher than the amplitude of the corresponding signal components in

c_{1}^{i} (n)

that were spread by the numerical dispersion compensation. The important step for autocorrelation suppression is to decide whether the strongest signal component in the dispersion compensated spectrum (

c_{1}^{i} (n)

) is higher than the strongest signal component in the non dispersion compensated spectrum (

c_{2}^{i} (n)

) or vice versa. In case of True, meaning

| c_{1}^{i} (n_{1}^{i}) | > | c_{2}^{i} (n_{2}^{i}) |

, we can now subtract both the complex conjugate mirror components in

c_{1}^{i} (n)

and the corresponding signal components in

c_{2}^{i} (n)

using

c_{1}^{i} (n_{1}^{i})

and

n_{1}^{i}

. In case of False, that is the strongest signal component is found in the non dispersion compensated spectrum

c_{2}^{i} (n) (| c_{1}^{i} (n_{1}^{i}) | < | c_{2}^{i} (n_{2}^{i}) |)

, this component is said to be an autocorrelation artifact. In this case we proceed similar and subtract the corresponding signal components in

c_{1}^{i} (n)

. In both cases the outputs are updated accordingly. Note that the subtraction of true signal components in

c_{2}^{i} (n)

and the subtraction of autocorrelation signal components in

c_{1}^{i} (n)

does neither affect any of the autocorrelation signals in

c_{2}^{i} (n)

nor any of the true signals in

c_{1}^{i} (n)

, so that our approach assures the correct reconstruction of the spectrum around the zero delay position even in highly scattering samples which is not the case in all other implementations of the DEFR algorithm. Although the autocorrelation removal requires one additional Fourier transform for the processing of each A-scan, the iterative procedure comprises only two multiplications and two subtractions for each array element of

c_{1}^{i} (n)

and

c_{2}^{i} (n)

in each iteration. As the inverse Fourier transform of the dispersive phase functions e^±iϕ(k_n) contains only half the number of values significantly greater than zero compared to the inverse Fourier transform of twice the inverse dispersive phase function e^−i2ϕ(k_n), the computational effort for autocorrelation removal can be further reduced by using an appropriate part of these functions. Again, the residual spectrum

c_{1}^{i} (n)

can be added to the output

d_{1}^{i} (n)

Fig. 3 Advanced implementation of the DEFR algorithm with autocorrelation removal. While the basic idea is similar to the simple version presented in Sect. 2.2, the advanced DEFR algorithm with autocorrelation removal exploits the information provided by both the non dispersion compensated and the dispersion compensated spectrum in z-space

c_{1}^{i} (n)

and

c_{2}^{i} (n)

, respectively, which are retrieved from the detector signal after resampling in frequency and background correction (step (2)). Again, the required inverse Fourier transforms (IFFT) of the dispersive phase functions arising in Eq. (23) and Eq. (25) (denoted as

p_{1}^{+} (n)

p_{1}^{-} (n)

and

p_{2}^{-} (n)

) are calculated separately from the actual processing steps in step (1). Two peak detectors determine simultaneously the strongest signal components in both spectra and provide their values and positions for the next processing steps (step (3)). Depending on whether the amplitude of

c_{1}^{i} (n_{1}^{i})

is higher than the amplitude of

c_{2}^{i} (n_{2}^{i})

or vice versa, this component is a true signal component or an autocorrelation artifact. In both cases the subtraction of the corresponding signals is performed directly in z-space and the located signals are added to the output. The iteration finishes if a maximum number of iterations is reached or if the residual spectrum

c_{1}^{i} (n)

contains only noise. Optionally, the residual spectrum

c_{1}^{i} (n)

can be added to the output

d_{1}^{i} (n)

2.4. Implementation of the advanced DEFR algorithm with autocorrelation removal

We adopt the DEFR algorithm, which uses a greedy matching pursuit algorithm for ℓ₁ optimization [23

, 32

M. Duarte, M. Davenport, M. Wakin, and R. Baraniuk, “Sparse signal detection from incoherent projections,” in Proc. Int. Conf. on Acoustics, Speech and Signal Processing (ICASSP) 3, III305–308 (2006).

], and change the method of subtracting complex conjugate mirror components. Furthermore, we add a method to remove autocorrelation artifacts by exploiting the information provided by the different spectra Eq. (23) and Eq. (25). The discrete inverse Fourier transforms of the dispersion compensated spectrum Eq. (23) and the non dispersion compensated spectrum Equation (25) as well as the necessary phase functions are calculated as follows:

p_{1}^{+} (n) = \frac{1}{N} \sum_{j = 0}^{N - 1} e^{+ i ϕ (k_{j})} e^{i 2 π \frac{n}{N} j},

p_{1}^{-} (n) = \frac{1}{N} \sum_{j = 0}^{N - 1} e^{- i ϕ (k_{j})} e^{i 2 π \frac{n}{N} j},

p_{2}^{-} (n) = \frac{1}{N} \sum_{j = 0}^{N - 1} e^{- i 2 ϕ (k_{j})} e^{i 2 π \frac{n}{N} j},

c_{1}^{0} (n) = \frac{1}{N} \sum_{j = 0}^{N - 1} P_{int} (k_{j}) e^{- i ϕ (k_{j})} e^{i 2 π \frac{n}{N} j} (dispersion compensated sprectrum),

c_{0}^{2} (n) = \frac{1}{N} \sum_{j = 0}^{N - 1} P_{int} (k_{j}) e^{i 2 π \frac{n}{N} j} (non dispersion compensated spectrum),

where n is the pixel index and N is the number of pixels in the spectrum. If the data are equidistant sampled in k-space, one will certainly use an IFFT algorithm to perform the discrete inverse Fourier transforms. While the block diagrams presented in Sect. 2.2 and Sect. 2.3 give an overview of the necessary processing steps, we will now explain the implementation of the advanced DEFR algorithm with autocorrelation removal in more detail.

Advanced DEFR Algorithm With Autocorrelation Removal

Initialization:
i = 0, $d_{1}^{0} (n) = 0$ , $d_{2}^{0} (n) = 0$ .
Determine the signal components with the highest contribution to the spectra in $c_{1}^{i} (n)$ , and $c_{2}^{i} (n)$ :
$n_{1}^{i} = arg max_{n = - N / 2 \dots N / 2 - 1} | c_{1}^{i} (n) |$ ,
$n_{2}^{i} = arg max_{n = - N / 2 \dots 0} | c_{2}^{i} (n) |$ .
Compare $c_{1}^{i} (n_{1}^{i})$ and $c_{2}^{i} (n_{2}^{i})$ with predefined threshold value T:
if max ${| c_{1}^{i} (n_{1}^{i}) |, | c_{2}^{i} (n_{2}^{i}) |} < T$ then proceed with step 6, otherwise proceed with step 4.
Update the outputs and subtract corresponding mirror artifacts:
if $| c_{1}^{i} (n_{1}^{i}) | > | c_{2}^{i} (n_{2}^{i}) |$ then

$\begin{array}{l} d_{1}^{i + 1} (n) = d_{1}^{i} (n) + c_{1}^{i} (n_{1}^{i}), \\ d_{2}^{i + 1} (n) = d_{2}^{i} (n), \\ c_{1}^{i + 1} (n) = c_{1}^{i} (n) - c_{1}^{i} (n_{1}^{i}) - {(c_{1}^{i} (n_{1}^{i}))}^{*} p_{2}^{-} (n + n_{1}^{i}), \\ c_{2}^{i + 1} (n) = c_{2}^{i} (n) - {\begin{array}{l} c_{1}^{i} (n_{1}^{i}) p_{1}^{+} (n - n_{1}^{i}), & if n_{1}^{i} \leq N / 2 \\ {(c_{1}^{i} (n_{1}^{i}))}^{*} p_{1}^{-} (n + n_{1}^{i}), & if n_{1}^{i} > N / 2 \end{array} \end{array}$

else

$\begin{array}{l} d_{1}^{i + 1} (n) = d_{1}^{i} (n), \\ d_{2}^{i + 1} (n) = d_{2}^{i} (n) + c_{2}^{i} (n_{2}^{i}), \\ c_{1}^{i + 1} (n) = c_{1}^{i} (n) - c_{2}^{i} (n_{2}^{i}) p_{1}^{-} (n - n_{2}^{i}) - {(c_{2}^{i} (n_{2}^{i}))}^{*} p_{1}^{-} (n + n_{2}^{i}), \\ c_{2}^{i + 1} (n) = c_{2}^{i} (n) - c_{2}^{i} (n_{2}^{i}), \end{array}$

where the center of the spectrum (DC term) is located at n = 0 and the corresponding dispersive phase function is shifted at the artifacts that are intended to be removed.
Check if maximum number of iterations M reached:
i = i + 1,
if i ≤ M then proceed with step 2, otherwise proceed with step 6.
(Optional) Add remaining spectrum to the output:
$d_{1} (n) = d_{1}^{i} (n) + c_{1}^{i} (n)$ ,
$d_{2} (n) = d_{2}^{i} (n)$ .

These six steps describe the implementation of the new DEFR algorithm. We found another correction helpful in some cases. Because a true signal component always causes a broadened complex conjugate mirror artifact in the spectrum, these two signals will overlap when the signal is near the zero delay position. Therefore, it is necessary to calculate the true signal component from the complex value in the spectrum using the following assumptions: A complex value

c_{1}^{i} (n_{1}^{i})

in the spectrum is the sum of the true signal component

{\hat{c}}_{1}^{i} (n_{1}^{i})

and its complex conjugate mirror component

{({\hat{c}}_{1}^{i} (n_{1}^{i}))}^{*}

spread by the dispersive phase function

p_{2}^{-} (n)

c_{1}^{i} (n_{1}^{i}) = {\hat{c}}_{1}^{i} (n_{1}^{i}) + {({\hat{c}}_{1}^{i} (n_{1}^{i}))}^{*} p_{2}^{-} (2 n_{i}^{1}),

(27)

where

2 n_{1}^{i}

is twice the distance between the signal and the DC term located at n = 0 (so the dispersive phase function is centered at the complex conjugate mirror component). Using Eq. (27) we can retrieve the true signal

{\hat{c}}_{1}^{i} (n_{1}^{i})

{\hat{c}}_{1}^{i} (n_{1}^{i}) = \frac{c_{1}^{i} (n_{1}^{i}) - {(c_{1}^{i} (n_{1}^{i}))}^{*} p_{2}^{-} (t)}{1 - {| p_{2}^{-} (t) |}^{2}} .

(28)

This correction is necessary if the exact true signal component

{\hat{c}}_{1}^{i} (n_{1}^{i})

needs to be calculated but has a less significant impact if the dispersion is very high as the amplitude of the dispersive phase function

p_{2}^{-} (n)

is accordingly low. Nevertheless, we found that the image quality of pictures acquired with little dispersion benefits from this correction.

3. Materials and methods

3.1. Experimental OCT setup

The setup is slightly modified from a dual-band OCT system which enables simultaneous imaging in the 0.8μm and 1.3μm wavelength range [34

P. Cimalla, J. Walther, M. Mehner, M. Cuevas, and E. Koch, “Simultaneous dual-band optical coherence tomography in the spectral domain for high resolution in vivo imaging,” Opt. Express 17, 19486–19500 (2009). [CrossRef] [PubMed]

]. In brief, we use only the 800nm band of the dual-band OCT system, which is illuminated by a SuperK Versa Super Continuum Source (Koheras A/S, Denmark) light source. In the used short wavelength band the spectrometer has a central wavelength of λ₀ = 791nm with a spectral width of Δλ = 166nm focused on the line scan camera. Compared to the original setup, the detector was replaced by a Basler sprint spL2048–70km line scan camera. The sensitivity was measured to be −93dB with a spectrometer roll-off of approximately 20dB over the entire depth scan range (dispersion balanced). The fiber coupled scanner head includes the Michelson interferometer and two galvanometer scanners for beam steering. Reference and sample arm of the interferometer are nearly dispersion balanced, leading to a nonlinear phase shift with a maximum mismatch less than 7rad. To introduce additional dispersion for DEFR imaging, a glass plate was placed in the reference arm of the interferometer. Several glass plates have been used starting from a single 6.35mm SF6 plate up to two SF6 plates resulting in total thickness of 12.70mm. For even higher dispersion a 6.25mm plate of ZnSe was inserted into the reference arm. ZnSe is transparent from 600nm on, has a refractive index of ≈ 2.5 at 800nm with dispersion approximately five times as high as SF6.

3.2. Dispersion measurement

There are different methods for estimating the dispersive phase. Using the Sellmeier Eq.:

n^{2} (λ) = 1 + \frac{B_{1} λ^{2}}{λ^{2} - C_{1}} + \frac{B_{2} λ^{2}}{λ^{2} - C_{2}} + \frac{B_{3} λ^{2}}{λ^{2} - C_{3}},

(29)

with given Sellmeier coefficients and a well known spectral range on the line detector, the dispersive phase can be calculated analytically. The disadvantage of this method is the fact that the spectral range must be exactly known, otherwise there will be inaccuracys. Another method is to maximize a sharpness metric applied on the measured spectrum of a sharp reflex over a given range of parameters of the Taylor series expansion of the propagation constant k(ω) [23

], which has the disadvantage that only a limited number of coefficients can be considered and that it is not sure that the global maximum is found during the iterative optimization. Moreover, many spectra with different parameter settings must be calculated in order to have a reliable set of data used for optimization. We propose a method that is similar to that one proposed by Witte et al. and Hofer et al. [24

B. Hofer, B. Považay, A. Unterhuber, L. Wang, B. Hermann, S. Rey, G. Matz, and W. Drexler, “Fast dispersion encoded full range optical coherence tomography for retinal imaging at 800 nm and 1060 nm,” Opt. Express 18, 4898–4919 (2010). [CrossRef] [PubMed]

, 29

], but has the crucial advantages that it does not affect the spectrometer calibration and does not require measurements in different depths. If the function for resampling from detector pixel number to frequency is known, then this method is favorable as it uses only one measurement of a strong reflector and does not depend on the method for resampling. We can calculate the dispersive phase for a given amount of dispersive material in the reference arm of the interferometer by simply measuring the spectrum of a strong reflector like a mirror or a glass surface. The spectrum of a glass surface in z-space is shown in Fig. 4(a). Note that there is a strong DC component as the reflection coefficient of the glass surface is comparatively high (compare Eq. (25), term (a)). We filter the dispersion spread peak in the left half-plane (compare Fig. 4(b)) and calculate the corresponding signal in k-space via Fourier transform. The physical model for dispersion measurement can be derived from Eq. (21), term (b) and Eq. (25), term (b). The filtered signal in z-space:

\bar{c} (z) = (ℱ^{- 1} (P_{i} (k)) * r_{M} t_{G} r_{S, 0} e^{i ϕ_{0}} δ (z - 2 (z_{F} + k_{lin} c_{0} z_{G} - n_{ref} z_{S, 0})) * ℱ^{- 1} (e^{i ϕ (k)})) (z),

(30)

where r_S,0 and z_S,0 are the reflection coefficient and the geometrical position of the glass surface, respectively, yields the necessary phase information in k-space:

ℱ (\bar{c} (z)) (k) = P_{i} (k) r_{M} t_{G} r_{S, 0} e^{i ϕ_{0}} e^{i 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, 0})} e^{i ϕ (k)} .

(31)

The phase of the spectrum in Eq. (31) in k-space is retrieved and unwrapped to remove 2π phase jumps:

ϕ (k) = ϕ_{0} + 2 (k z_{F} + k_{lin} z_{G} - k n_{ref} z_{S, 0}) + ϕ (k) = φ_{lin} (k) + ϕ (k) .

(32)

A linear (first order) polynomial regression φ_lin (k) is estimated and subtracted from φ(k) in order to retain the dispersive component:

ϕ (k) = φ (k) - φ_{lin} (k) .

(33)

Another method is to subtract the linear function:

φ_{lin} (k) = ϕ (k_{a}) + \frac{ϕ (k_{b}) - ϕ (k_{a})}{k_{b} - k_{a}} (k - k_{a}),

(34)

where k_a and k_b denote the left and right positions in k-space, respectively. This method has the advantage that the dispersive phase ϕ(k) is zero at the left and right end and therefore its inverse Fourier transform is narrower and converges faster to zero. Nevertheless, the differences are very small and a remaining linear part in the dispersive phase ϕ(k) will only cause a shift of the spectrum in z-space. We used the latter method and the results are shown in Fig. 4 for a 6.25mm ZnSe glass plate in the system described in Sect. 3.1.

Fig. 4 Procedure of dispersion measurement. The dispersive phase ϕ(k) can be directly extracted from the measured interference signal of a sharp reflector, i.e. the surface of a glass plate (Fig. 4(a)). Note that after resampling to 4096 spectral points in k-space and inverse Fourier transform only the middle half of the depth signal in z-space, comprising 2048 complex values, is displayed and the range of the normalized position covers these 2048 points. The dispersion broadened peak in z-space is filtered (Fig. 4(b)) and transformed into k-space via Fourier transform. The phase of the complex signal in k-space is retrieved and unwrapped (Fig. 4(c)). The dispersive phase is finally obtained after subtracting the linear part of the nonlinear phase (Fig. 4(d)).

The dispersive phases ϕ(k) for various amounts of dispersive material in the reference arm of the interferometer have been measured using the technique described in Sect. 3.2 and are depicted in Fig. 5(a). In order to suppress noise and artifacts, the average of 50 A-scans and a fourth order polynomial regression of ϕ(k) was used. While averaging of several A-scans improves the quality of the measured dispersive phases, therefore reducing the width of a numerically dispersion compensated peak in z-space, polynomial regression is not absolutely necessary but can remove existing detector inaccuracys. The group velocity dispersions

GVD = \frac{d^{2} k (ω)}{d ω^{2}}

at the center frequency λ₀ = 791nm for SF6 and ZnSe,

GVD (SF 6) = 202 \frac{{fs}^{2}}{mm}

and

GVD (ZnSe) = 1019 \frac{{fs}^{2}}{mm}

, respectively, and the third order dispersions

{k^{'}}^{'}^{'} = \frac{d^{3} k (ω)}{d ω^{3}}

{k^{'}}^{'}^{'} (SF 6) = 126 \frac{{fs}^{3}}{mm}

and

{k^{'}}^{'}^{'} (ZnSe) = 950 \frac{{fs}^{3}}{mm}

, respectively, cause the major part of the dispersion. While the ZnSe and SF6 glass plates have a comparable thickness of 6.25mm versus 6.35mm, ZnSe causes a dispersive phase that is more than five times as high as that one introduced by SF6 (compare GVD(ZnSe) and GVD(SF6)), which makes it favorable for systems with limited space or optical path length since it enables the correct reconstruction of images showing a much larger span of echoes. Figure 5(b) shows the inverse Fourier transforms of twice the inverse dispersive phase functions (ℱ⁻¹ (e⁻ⁱ²^ϕ⁽^k⁾) corresponding to the dispersion spread complex conjugate mirror artifact of a true signal component (compare Eq. (23), term (b) and Eq. (23), term (c)). As already described in Sect. 2.2, according to Parseval’s theorem, a wider range (in z-space) of the mirror components due to a larger dispersion leads to a lower amplitude. A sufficient amount of dispersion, for instance 6.25mm ZnSe, leads to an amplitude ratio of a unity-peak compared to the amplitude of the complex conjugate mirror artifacts after numerical dispersion compensation of

R \approx \frac{1}{0.03} = 33.33 ≜ 30.46 dB

Fig. 5 Measured dispersive phases and shape of the corresponding complex conjugate mirror components. The dispersive phases were extracted from the signal phase in k-space, which was averaged over 50 A-scans, and fitted to a fourth order polynomial regression (depicted in Fig. 5(a)). The corresponding inverse Fourier transforms of twice the inverse dispersive phase functions, belonging to the complex conjugate mirror components of a true signal after numerical dispersion compensation, are depicted in Fig. 5(b).

4. Results and discussion

4.1. Necessary amount of dispersion

The presented algorithm treats each signal independently and is not grouping signals as proposed by Hofer et al. [24

, 26

]. Therefore a high dispersion, leading to a large ratio between a true signal and its complex conjugate mirror artifact, does not accelerate the algorithm. Nevertheless, we found that the algorithm needs a minimum dispersion to sort the signals to the correct side. We noted that in images with only a few sharp reflexes a small amount of dispersion between the sample and the reference arm of the interferometer already leads to the expected images, while images with a larger span of echoes show large artifacts at low dispersion. False mapping can occur if neighboring signals on the false side that are distributed over a certain range add up to a signal higher than any (residual) component on the correct side. To get more insight in the amount of dispersion needed for a successful recovery of the image we made a simulation with synthetic data. We assumed a band of echoes with Rayleigh distributed amplitude [33

B. Karamata, K. Hassler, M. Laubscher, and T. Lasser, “Speckle statistics in optical coherence tomography,” J. Opt. Soc. Am. A 22, 593–596 (2005). [CrossRef]

] and arbitrary phase over a certain number N of Fourier components. On this data in one side of the full range spectrum in z-space we applied a quadratic phase shift of arbitrary amplitude, which is the largest part in most cases and applied the presented DEFR algorithm. A quadratic phase shift of an amplitude

\frac{D π}{4}

broadens a single peak to a FWHM of D Fourier components. The recovery of the depth signal was regarded as successful if no signal was found on the false side. We found that the error rate increases strongly with the ratio of the number of Fourier components N divided by the broadening D. While the error rate is not only a function of the quotient

\frac{N}{D}

but is slightly higher at larger N, for a ratio of

\frac{N}{D}

below 0.7 the error rate was below 5% and for a ratio below 0.5 it was much below 1%. This result, which should hold true if higher order terms contribute essentially to the dispersion, is in good agreement with the observations made during the recovery of different images at variable amount of dispersion. Note that this is only a statistical argument. As the amplitude of the dispersion broadened signal drops only with the square root of the dispersion D (see Eq. (26)), in rare cases a number

N \approx \sqrt{D}

of Fourier components with similar amplitude and appropriate phase may result in a false mapping of signals. We recommend to use as less dispersion as possible for reconstructing a scattering layer with a certain thickness, because a higher amount of dispersion reduces the measurement depth range due to partial aliasing of the broadened peaks.

4.2. High resolution full range imaging

In the following, we will present images of various samples using different amounts of dispersion and explain their special characteristics. All images are displayed using a dynamic range of 45dB and a Red Hot color scale. The zero delay position is marked with a white dashed line. The upper left image (a) shows the original image disturbed by the artificial dispersion. The upper right image (b) shows the dispersion compensated image, where dispersion compensation was performed by multiplication with the inverse dispersive phase function e⁻ⁱ^ϕ⁽^k⁾ in k-space prior to inverse Fourier transform. The lower left image (c) shows the image post-processed with the simple implementation of the DEFR algorithm presented in Sect. 2.2. Note that there is no difference in the results using our simple implementation of the DEFR algorithm and the implementation of Hofer et al. [23

] as both implementations are mathematically equal. The lower right image (d) finally shows our advanced implementation of the DEFR algorithm with autocorrelation removal. The iteration for a single A-scan was finished if a maximum number of iterations was reached (here: 250 iterations) or if the remaining complex conjugate mirror components of the true signals were below the noise level and consequently not causing artifacts if not removed before the residual spectrum was added to the reconstructed A-scan. As the latter condition mostly ensues first, the processing speed benefits of a larger amount of dispersion (for very high dispersion the approach proposed by Witte et al. [29

], using three Fourier transforms, is very efficient, but does not enable autocorrelation removal). If an exact signal recovery is necessary, our advanced DEFR algorithm is to the best of our knowledge the first implementation that provides full range OCT in combination with autocorrelation removal.

Figure 6 shows the images of a nanoparticle target used for resolution estimation [35

P. D. Woolliams and P. H. Tomlins, “Estimating the resolution of a commercial optical coherence tomography system with limited spatial sampling,” Meas. Sci. Technol. 22, 065502 (2011). [CrossRef]

] that was imaged using 6.35mm SF6. The axial scanning range is 3.89mm with

1.9 \frac{μ m}{Pixel}

in air and 2048 Pixel in z-direction. Figure 6(a) shows the original image that is symmetric around the zero delay position as the amplitude spectrum of the real valued signal from the detector is an even function. The effect of the artificial dispersion is clearly visible. Figure 6(b) shows that numerical dispersion compensation leads to a sharp reconstruction of the true image components and spreads the complex conjugate mirror components even more. The sharp structures in the image allow a correct reconstruction with little dispersion. Due to the tilted surface autocorrelation artifacts are not very severe (compare Fig. 6(c) and Fig. 6(d) – autocorrelation signals are hardly visible in Fig. 6(c) and removed in Fig. 6(d)). Figure 6 shows that both our simple and advanced DEFR algorithm successfully reconstruct images with sharp structures using only a small amount of dispersion, whereas images with a larger span of echoes require a larger amount of dispersion for correct reconstruction.

Fig. 6 Nanoparticle target; 6.35mm SF6; 45dB dynamic range; 3.89mm × 4.09mm (h × w); (a) – original image without dispersion compensation; (b) – dispersion compensated image; (c) – image without autocorrelation removal (simple DEFR processing); (d) – image with autocorrelation removal (advanced DEFR processing)

Figure 7 shows the images of a capillary originally used for Doppler measurements. The inner diameter is d = 300μm and the intralipid inside the capillary has a concentration c = 1% leading to a refractive index n = 1.34. The image does intentionally not cross the zero delay position as a false mapping of signals is not visible in this position. Note that the aspect ratio of the image is stretched in order to accentuate the errors due to a false signal mapping over the width of the homogeneous back scattering intralipid layer. We used this sample to experimentally investigate the relation between the amount of dispersion and the width of a correct reconstructed homogenous back scattering layer. The horizontal white dashed lines indicate the maximal width the intralipid layer may have for a proper reconstruction, which is around 200μm. This is in very good agreement with the estimations made in Sect. 4.1. The scattering layer with a width of 200μm has an optical length of 200μm × 1.34 = 268μm. According to

1.9 \frac{μ m}{Pixel}

, this results in an

N = \frac{268 μ m}{1.9 μ m / Pixel} = 141 Pixel

wide signal in z-space (number of Fourier components). Using the width of a dispersion spread peak D = 215 Pixel (12.70 mm SF6), the ratio

\frac{N}{D} = \frac{141 Pixel}{215 Pixel} = 0.656

is below 0.7, therefore implying an error ratio below 5%. As opposed to this, the scattering layer in the center of the image, having a width of around 300μm, leads to a ratio of

\frac{N}{D} = \frac{212 Pixel}{215 Pixel} = 0.986

, which is well above 0.7 and therefore indicates that there will be lots of errors in the image. These estimations hold true when comparing the number of errors outside and inside the vertical white dashed lines in the image.

Fig. 7 Capillary; 12.70mm SF6; 45dB dynamic range; 3.89mm × 0.54mm (h × w); (a) – original image without dispersion compensation; (b) – dispersion compensated image; (c) – image without autocorrelation removal (simple DEFR processing); (d) – image with autocorrelation removal (advanced DEFR processing)

Figure 8 shows the same sample used in Fig. 7, but imaged with a much higher dispersion using 6.25mm ZnSe. As can be seen from Fig. 8(c) and Fig. 8(d), the signal recovery is correct apart from some very little errors and the autocorrelation signals are not very disturbing. The horizontal white dashed lines indicate the maximal width of the properly reconstructed intralipid layer having a thickness of around 300μm. The 300μm wide scattering layer and a dispersion spread peak having a width of D = 525Pixel in z-space (6.25mm ZnSe) lead to a ratio of

\frac{N}{D} = \frac{212 Pixel}{525 Pixel} = 0.404

, which is below 0.5, therefore predicting an error ratio well below 1%. Again, the experimental results conform with the results of the simulation that there will be only few errors but deviate from a predicted width of a (for the most part) properly reconstructed layer of around 500μm. Namely, due to the increasing number of errors in the center of the capillary, we assume that a homogenous scattering layer with a width of around 300μm (geometrical width) is the limit for a correct signal recovery using 6.25mm ZnSe and our OCT system described above, although the simulation predicts a higher limit. Note that while the dispersive phase and consequently the width of the dispersion spread complex conjugate mirror components corresponding to 6.25mm ZnSe is twice as high compared to 12.70mm SF6, the maximal width of the correct reconstructed scattering layer is not twice as high. This is probably caused by the overlap of autocorrelation and mirror artifacts seen in Fig. 8(b), which is not the case in Fig. 7(b) and was not taken into account in the simulation of Sect. 4.1.

Fig. 8 Capillary; 6.25mm ZnSe; 45dB dynamic range; 3.89mm × 0.54mm (h × w); (a) – original image without dispersion compensation; (b) – dispersion compensated image; (c) – image without autocorrelation removal (simple DEFR processing); (d) – image with autocorrelation removal (advanced DEFR processing)

Figure 9 shows the images of an optical longpass filter made from RG850 Schott glass having a thickness of 1.8mm. DC and autocorrelation removal is essential as the two highly reflecting surfaces otherwise cause very strong artifacts near the zero delay position. Comparing Fig. 9(c) and Fig. 9(d) it is obvious how well our advanced DEFR algorithm works and that full range imaging of samples having strong DC and autocorrelation signals requires their removal. After numerical dispersion compensation, the front and back surface of the filter are a little shifted in z-direction, which is probably caused by a remaining linear component in the dispersive phase ϕ(k) due to a non-ideal estimation of the linear component φ_lin (k) during the phase measurement procedure described in Sect. 3.2 (compare Eq. 33).

Fig. 9 Optical filter; 12.70mm SF6; 45dB dynamic range; 3.89mm×1.79mm (h×w); (a) – original image without dispersion compensation; (b) – dispersion compensated image; (c) – image without autocorrelation removal (simple DEFR processing); (d) – image with autocorrelation removal (advanced DEFR processing)

Figure 10 shows the images of an eardrum covering a large area in z-direction and crossing the zero delay position. While autocorrelation artifacts are still visible in Fig. 10(c), our advanced DEFR algorithm removes them from the spectrum and reconstructs a very clear image without any artifacts. Even strongly disturbed areas are reconstructed properly. Comparing Fig. 10(c) and Fig. 10(d), it becomes apparent that autocorrelation removal is important in common samples and should not be neglected in full range OCT imaging. The widest scattering layer in the image has a width of around N = 140Pixel in z-space and shows no artifacts, which is accordance with the simulation in Sect. 4.1, where the corresponding ratio of

\frac{N}{D} = \frac{140 Pixel}{525 Pixel} = 0.267

is stating that there will be hardly any errors.

Fig. 10 Excised human eardrum ex vivo; 6.25mm ZnSe; 45dB dynamic range; 3.89mm× 8.18mm (h × w); (a) – original image without dispersion compensation; (b) – dispersion compensated image; (c) – image without autocorrelation removal (simple DEFR processing); (d) – image with autocorrelation removal (advanced DEFR processing)

5. Conclusion

In this paper we presented an advanced processing algorithm for dispersion encoded full range OCT, which uses convolution in z-space instead of multiplikation in k-space, thereby signifi-cantly increasing the processing speed as only one Fourier transform per A-scan is necessary for image reconstruction. Besides the efficient complex conjugate artifact removal, autocorrelation artifacts can be removed likewise with only minor additional computational effort. Furthermore, we demonstrated the usability of ZnSE as a new dispersive material for DEFR OCT with a very high dispersion and less space demand. Apart from the stop criteria, the number of iterations M and the threshold value of noise T, no additional parameters are required for the execution of our DEFR algorithm enabling an automatic stable signal recovery. With increasing dispersion, the spatial frequency span on the detector increases thereby reducing the usable depth range as the Nyquist limit should be fulfilled on the total detector to avoid aliasing artifacts. Therefore, there is some optimal dispersion depending on the image. We presented an estimation for the necessary amount of dispersion depending on the width of the scattering layers in the image. This estimation shows that arbitrary images cannot be reconstructed correctly and thereby solves the philosophical question raised in the introduction that DEFR retrieves N complex values from only N real spectral data points. As the processing speed of our DEFR algorithm only slightly benefits from a larger dispersion (which reduces the amplitude of the complex conjugate artifacts and consequently stops the iteration earlier), our implementation is favorable if the amount of dispersion is limited. We used a standard consumer PC with an Intel Core i5-2400 CPU @ 3.10GHz (4 cores), 8GB RAM and LabVIEW 2010 (64-bit). The processing of 1024 A-scans upsampled to 4096 values in k-space, reconstructed using 6.25mm ZnSe and a maximum number of M = 250 iterations (which are very seldom reached as the mirror artifacts are well below noise before), requires around 1s including autocorrelation removal and 0.9s without, leading to around 1000 A-lines per second (images depicted in Fig. 10). A lower amount of dispersion can even accelerate the processing speed of our algorithm. For example, the inverse Fourier transforms of the dispersive phase functions for 6.35mm SF6 are narrower compared to 6.25mm ZnSe, thus reducing the number of array elements for each subtraction considerably. Using 6.35mm SF6 and M = 250, the reconstruction of 1024 A-scans requires around 0.5s leading to around 2000 A-lines per second (images depicted in Fig. 6). Therefore, using a GPU based processing, which increases the processing speed by a factor of more than 90 [26

], around 100,000 – 200,000 A-lines per second could be reconstructed with our DEFR algorithm with autocorrelation removal. A possible range of applications is high resolution imaging of biological samples that require a large depth scan range, e.g. the human eardrum. Also a future extension to full range Doppler OCT is conceivable as only one measurement is required thereby eliminating the critical phase sensitivity.

Acknowledgments

The authors would like to thank Peter H Tomlins for providing the nanoparticle target. This project was supported by SAB (Sächsische Aufbaubank, project 14302/2475).

References and links

1.	S. Marschall, B. Sander, M. Mogensen, T. Jørgensen, and P. Andersen, “Optical coherence tomography – current technology and applications in clinical and biomedical research,” Anal. Bioanal. Chem. 400, 2699–2720 (2011). [CrossRef] [PubMed]
2.	J. Walther, M. Gaertner, P. Cimalla, A. Burkhardt, L. Kirsten, S. Meissner, and E. Koch, “Optical coherence tomography in biomedical research,” Anal. Bioanal. Chem. 400, 2721–2743 (2011). [CrossRef] [PubMed]
3.	D. Huang, E. Swanson, C. Lin, J. Schuman, W. Stinson, W. Chang, M. Hee, T. Flotte, K. Gregory, and C. Puliafito “Optical coherence tomography,” Science 254, 1178–1181 (1991). [CrossRef] [PubMed]
4.	W. Wieser, B. R. Biedermann, T. Klein, C. M. Eigenwillig, and R. Huber, “Multi-megahertz oct: High quality 3d imaging at 20 million a-scans and 4.5 gvoxels per second,” Opt. Express 18, 14685–14704 (2010). [CrossRef] [PubMed]
5.	R. Leitgeb, C. Hitzenberger, and A. Fercher, “Performance of fourier domain vs. time domain optical coherence tomography,” Opt. Express 11, 889–894 (2003). [CrossRef] [PubMed]
6.	M. Choma, M. Sarunic, C. Yang, and J. Izatt, “Sensitivity advantage of swept source and fourier domain optical coherence tomography,” Opt. Express 11, 2183–2189 (2003). [CrossRef] [PubMed]
7.	J. F. de Boer, B. Cense, B. H. Park, M. C. Pierce, G. J. Tearney, and B. E. Bouma, “Improved signal-to-noise ratio in spectral-domain compared with time-domain optical coherence tomography,” Opt. Lett. 28, 2067–2069 (2003). [CrossRef] [PubMed]
8.	T. Bajraszewski, M. Wojtkowski, M. Szkulmowski, A. Szkulmowska, R. Huber, and A. Kowalczyk, “Improved spectral optical coherence tomography using optical frequency comb,” Opt. Express 16, 4163–4176 (2008). [CrossRef] [PubMed]
9.	M. Hagen-Eggert, P. Koch, and G. Hüttmann, “Analysis of the signal fall-off in spectral domain optical coherence tomography systems,” in Optical Coherence Tomography and Coherence Domain Optical Methods in Biomedicine XVI , Joseph A. Izatt, James G. Fujimoto, and Valery V. Tuchin, eds., Proc. SPIE 8213, 82131K (2012).
10.	Y. Mao, S. Sherif, C. Flueraru, and S. Chang, “3x3 mach-zehnder interferometer with unbalanced differential detection for full-range swept-source optical coherence tomography,” Appl. Opt. 47, 2004–2010 (2008). [CrossRef] [PubMed]
11.	H.-C. Cheng, J.-F. Huang, and Y.-H. Hsieh, “Numerical analysis of one-shot full-range fd-oct system based on orthogonally polarized light,” Opt. Commun. 282, 3040–3045 (2009). [CrossRef]
12.	P. Meemon, K.-S. Lee, and J. P. Rolland, “Doppler imaging with dual-detection full-range frequency domain optical coherence tomography,” Biomed. Opt. Express 1, 537–552 (2010). [CrossRef]
13.	T. Wu, Z. Ding, C. Wang, and M. Chen, “Full-range swept source optical coherence tomography based on carrier frequency by transmissive dispersive optical delay line,” J. Biomed. Opt. 16, 126008 (2011). [CrossRef] [PubMed]
14.	D. Y. Kim, J. S. Werner, and R. J. Zawadzki, “Comparison of phase-shifting techniques for in vivo full-range, high-speed fourier-domain optical coherence tomography,” J. Biomed. Opt. 15, 056011 (2010). [CrossRef] [PubMed]
15.	E. Götzinger, M. Pircher, R. Leitgeb, and C. Hitzenberger, “High speed full range complex spectral domain optical coherence tomography,” Opt. Express 13, 583–594 (2005). [CrossRef] [PubMed]
16.	F. Jaillon, S. Makita, M. Yabusaki, and Y. Yasuno, “Parabolic bm-scan technique for full range doppler spectral domain optical coherence tomography,” Opt. Express 18, 1358–1372 (2010). [CrossRef] [PubMed]
17.	S. Makita, T. Fabritius, and Y. Yasuno, “Full-range, high-speed, high-resolution 1-μm spectral-domain optical coherence tomography using bm-scan for volumetric imaging of the human posterior eye,” Opt. Express 16, 8406–8420 (2008). [CrossRef] [PubMed]
18.	Y. Yasuno, S. Makita, T. Endo, G. Aoki, M. Itoh, and T. Yatagai, “Simultaneous b-m-mode scanning method for real-time full-range fourier domain optical coherence tomography,” Appl. Opt. 45, 1861–1865 (2006). [CrossRef] [PubMed]
19.	L. An and R. K. Wang, “Use of a scanner to modulate spatial interferograms for in vivo full-range fourier-domain optical coherence tomography,” Opt. Lett. 32, 3423–3425 (2007). [CrossRef] [PubMed]
20.	C.-T. Wu, T.-T. Chi, Y.-W. Kiang, and C. C. Yang, “Computation time-saving mirror image suppression method in fourier-domain optical coherence tomography,” Opt. Express 20, 8270–8283 (2012). [CrossRef] [PubMed]
21.	M. Wojtkowski, A. Kowalczyk, R. Leitgeb, and A. F. Fercher, “Full range complex spectral optical coherence tomography technique in eye imaging,” Opt. Lett. 27, 1415–1417 (2002). [CrossRef]
22.	B. Hermann, B. Hofer, C. Meier, and W. Drexler, “Spectroscopic measurements with dispersion encoded full range frequency domain optical coherence tomography in single- and multilayered non- scattering phantoms,” Opt. Express 17, 24162–24174 (2009). [CrossRef]
23.	B. Hofer, B. Považay, B. Hermann, A. Unterhuber, G. Matz, and W. Drexler, “Dispersion encoded full range frequency domain optical coherence tomography,” Opt. Express 17, 7–24 (2009). [CrossRef] [PubMed]
24.	B. Hofer, B. Považay, A. Unterhuber, L. Wang, B. Hermann, S. Rey, G. Matz, and W. Drexler, “Fast dispersion encoded full range optical coherence tomography for retinal imaging at 800 nm and 1060 nm,” Opt. Express 18, 4898–4919 (2010). [CrossRef] [PubMed]
25.	L. Wang, B. Hofer, Y.-P. Chen, J. A. Guggenheim, W. Drexler, and B. Povazay, “Highly reproducible swept-source, dispersion-encoded full-range biometry and imaging of the mouse eye,” J. Biomed. Opt. 15, 046004 (2010). [CrossRef] [PubMed]
26.	L. Wang, B. Hofer, J. A. Guggenheim, and B. Povazay, “Graphics processing unit-based dispersion encoded full-range frequency-domain optical coherence tomography,” J. Biomed. Opt. 17, 077007 (2012). [CrossRef] [PubMed]
27.	M. Wojtkowski, V. Srinivasan, T. Ko, J. Fujimoto, A. Kowalczyk, and J. Duker, “Ultrahigh-resolution, high-speed, fourier domain optical coherence tomography and methods for dispersion compensation,” Opt. Express 12, 2404–2422 (2004). [CrossRef] [PubMed]
28.	B. Cense, N. Nassif, T. Chen, M. Pierce, S.-H. Yun, B. Park, B. Bouma, G. Tearney, and J. de Boer, “Ultrahigh-resolution high-speed retinal imaging using spectral-domain optical coherence tomography,” Opt. Express 12, 2435–2447 (2004). [CrossRef] [PubMed]
29.	S. Witte, M. Baclayon, E. J. Peterman, R. F. Toonen, H. D. Mansvelder, and M. L. Groot, “Single-shot two-dimensional full-range optical coherence tomography achieved by dispersion control,” Opt. Express 17, 11335–11349 (2009). [CrossRef] [PubMed]
30.	W. Drexler and J. Fujimoto, Optical Coherence Tomography: Technology and Applications (Springer, 2008). [CrossRef]
31.	A. Fercher, C. Hitzenberger, M. Sticker, R. Zawadzki, B. Karamata, and T. Lasser, “Numerical dispersion compensation for partial coherence interferometry and optical coherence tomography,” Opt. Express 9, 610–615 (2001). [CrossRef] [PubMed]
32.	M. Duarte, M. Davenport, M. Wakin, and R. Baraniuk, “Sparse signal detection from incoherent projections,” in Proc. Int. Conf. on Acoustics, Speech and Signal Processing (ICASSP) 3, III305–308 (2006).
33.	B. Karamata, K. Hassler, M. Laubscher, and T. Lasser, “Speckle statistics in optical coherence tomography,” J. Opt. Soc. Am. A 22, 593–596 (2005). [CrossRef]
34.	P. Cimalla, J. Walther, M. Mehner, M. Cuevas, and E. Koch, “Simultaneous dual-band optical coherence tomography in the spectral domain for high resolution in vivo imaging,” Opt. Express 17, 19486–19500 (2009). [CrossRef] [PubMed]
35.	P. D. Woolliams and P. H. Tomlins, “Estimating the resolution of a commercial optical coherence tomography system with limited spatial sampling,” Meas. Sci. Technol. 22, 065502 (2011). [CrossRef]

OCIS Codes
(100.3020) Image processing : Image reconstruction-restoration
(100.5070) Image processing : Phase retrieval
(110.0110) Imaging systems : Imaging systems
(110.4500) Imaging systems : Optical coherence tomography
(170.0170) Medical optics and biotechnology : Medical optics and biotechnology
(170.3010) Medical optics and biotechnology : Image reconstruction techniques
(170.4500) Medical optics and biotechnology : Optical coherence tomography
(300.6300) Spectroscopy : Spectroscopy, Fourier transforms
(110.3010) Imaging systems : Image reconstruction techniques

ToC Category:
Image Processing

History
Original Manuscript: August 28, 2012
Revised Manuscript: October 8, 2012
Manuscript Accepted: October 9, 2012
Published: October 16, 2012

Citation
Felix Köttig, Peter Cimalla, Maria Gärtner, and Edmund Koch, "An advanced algorithm for dispersion encoded full range frequency domain optical coherence tomography," Opt. Express 20, 24925-24948 (2012)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-20-22-24925

Sort: Author | Year | Journal | Reset

References

S. Marschall, B. Sander, M. Mogensen, T. Jørgensen, and P. Andersen, “Optical coherence tomography – current technology and applications in clinical and biomedical research,” Anal. Bioanal. Chem.400, 2699–2720 (2011). [CrossRef] [PubMed]
J. Walther, M. Gaertner, P. Cimalla, A. Burkhardt, L. Kirsten, S. Meissner, and E. Koch, “Optical coherence tomography in biomedical research,” Anal. Bioanal. Chem.400, 2721–2743 (2011). [CrossRef] [PubMed]
D. Huang, E. Swanson, C. Lin, J. Schuman, W. Stinson, W. Chang, M. Hee, T. Flotte, K. Gregory, C. Puliafito, and “Optical coherence tomography,” Science254, 1178–1181 (1991). [CrossRef] [PubMed]
W. Wieser, B. R. Biedermann, T. Klein, C. M. Eigenwillig, and R. Huber, “Multi-megahertz oct: High quality 3d imaging at 20 million a-scans and 4.5 gvoxels per second,” Opt. Express18, 14685–14704 (2010). [CrossRef] [PubMed]
R. Leitgeb, C. Hitzenberger, and A. Fercher, “Performance of fourier domain vs. time domain optical coherence tomography,” Opt. Express11, 889–894 (2003). [CrossRef] [PubMed]
M. Choma, M. Sarunic, C. Yang, and J. Izatt, “Sensitivity advantage of swept source and fourier domain optical coherence tomography,” Opt. Express11, 2183–2189 (2003). [CrossRef] [PubMed]
J. F. de Boer, B. Cense, B. H. Park, M. C. Pierce, G. J. Tearney, and B. E. Bouma, “Improved signal-to-noise ratio in spectral-domain compared with time-domain optical coherence tomography,” Opt. Lett.28, 2067–2069 (2003). [CrossRef] [PubMed]
T. Bajraszewski, M. Wojtkowski, M. Szkulmowski, A. Szkulmowska, R. Huber, and A. Kowalczyk, “Improved spectral optical coherence tomography using optical frequency comb,” Opt. Express16, 4163–4176 (2008). [CrossRef] [PubMed]
M. Hagen-Eggert, P. Koch, and G. Hüttmann, “Analysis of the signal fall-off in spectral domain optical coherence tomography systems,” in Optical Coherence Tomography and Coherence Domain Optical Methods in Biomedicine XVI, Joseph A. Izatt, James G. Fujimoto, and Valery V. Tuchin, eds., Proc. SPIE 8213, 82131K (2012).
Y. Mao, S. Sherif, C. Flueraru, and S. Chang, “3x3 mach-zehnder interferometer with unbalanced differential detection for full-range swept-source optical coherence tomography,” Appl. Opt.47, 2004–2010 (2008). [CrossRef] [PubMed]
H.-C. Cheng, J.-F. Huang, and Y.-H. Hsieh, “Numerical analysis of one-shot full-range fd-oct system based on orthogonally polarized light,” Opt. Commun.282, 3040–3045 (2009). [CrossRef]
P. Meemon, K.-S. Lee, and J. P. Rolland, “Doppler imaging with dual-detection full-range frequency domain optical coherence tomography,” Biomed. Opt. Express1, 537–552 (2010). [CrossRef]
T. Wu, Z. Ding, C. Wang, and M. Chen, “Full-range swept source optical coherence tomography based on carrier frequency by transmissive dispersive optical delay line,” J. Biomed. Opt.16, 126008 (2011). [CrossRef] [PubMed]
D. Y. Kim, J. S. Werner, and R. J. Zawadzki, “Comparison of phase-shifting techniques for in vivo full-range, high-speed fourier-domain optical coherence tomography,” J. Biomed. Opt.15, 056011 (2010). [CrossRef] [PubMed]
E. Götzinger, M. Pircher, R. Leitgeb, and C. Hitzenberger, “High speed full range complex spectral domain optical coherence tomography,” Opt. Express13, 583–594 (2005). [CrossRef] [PubMed]
F. Jaillon, S. Makita, M. Yabusaki, and Y. Yasuno, “Parabolic bm-scan technique for full range doppler spectral domain optical coherence tomography,” Opt. Express18, 1358–1372 (2010). [CrossRef] [PubMed]
S. Makita, T. Fabritius, and Y. Yasuno, “Full-range, high-speed, high-resolution 1-μm spectral-domain optical coherence tomography using bm-scan for volumetric imaging of the human posterior eye,” Opt. Express16, 8406–8420 (2008). [CrossRef] [PubMed]
Y. Yasuno, S. Makita, T. Endo, G. Aoki, M. Itoh, and T. Yatagai, “Simultaneous b-m-mode scanning method for real-time full-range fourier domain optical coherence tomography,” Appl. Opt.45, 1861–1865 (2006). [CrossRef] [PubMed]
L. An and R. K. Wang, “Use of a scanner to modulate spatial interferograms for in vivo full-range fourier-domain optical coherence tomography,” Opt. Lett.32, 3423–3425 (2007). [CrossRef] [PubMed]
C.-T. Wu, T.-T. Chi, Y.-W. Kiang, and C. C. Yang, “Computation time-saving mirror image suppression method in fourier-domain optical coherence tomography,” Opt. Express20, 8270–8283 (2012). [CrossRef] [PubMed]
M. Wojtkowski, A. Kowalczyk, R. Leitgeb, and A. F. Fercher, “Full range complex spectral optical coherence tomography technique in eye imaging,” Opt. Lett.27, 1415–1417 (2002). [CrossRef]
B. Hermann, B. Hofer, C. Meier, and W. Drexler, “Spectroscopic measurements with dispersion encoded full range frequency domain optical coherence tomography in single- and multilayered non- scattering phantoms,” Opt. Express17, 24162–24174 (2009). [CrossRef]
B. Hofer, B. Považay, B. Hermann, A. Unterhuber, G. Matz, and W. Drexler, “Dispersion encoded full range frequency domain optical coherence tomography,” Opt. Express17, 7–24 (2009). [CrossRef] [PubMed]
B. Hofer, B. Považay, A. Unterhuber, L. Wang, B. Hermann, S. Rey, G. Matz, and W. Drexler, “Fast dispersion encoded full range optical coherence tomography for retinal imaging at 800 nm and 1060 nm,” Opt. Express18, 4898–4919 (2010). [CrossRef] [PubMed]
L. Wang, B. Hofer, Y.-P. Chen, J. A. Guggenheim, W. Drexler, and B. Povazay, “Highly reproducible swept-source, dispersion-encoded full-range biometry and imaging of the mouse eye,” J. Biomed. Opt.15, 046004 (2010). [CrossRef] [PubMed]
L. Wang, B. Hofer, J. A. Guggenheim, and B. Povazay, “Graphics processing unit-based dispersion encoded full-range frequency-domain optical coherence tomography,” J. Biomed. Opt.17, 077007 (2012). [CrossRef] [PubMed]
M. Wojtkowski, V. Srinivasan, T. Ko, J. Fujimoto, A. Kowalczyk, and J. Duker, “Ultrahigh-resolution, high-speed, fourier domain optical coherence tomography and methods for dispersion compensation,” Opt. Express12, 2404–2422 (2004). [CrossRef] [PubMed]
B. Cense, N. Nassif, T. Chen, M. Pierce, S.-H. Yun, B. Park, B. Bouma, G. Tearney, and J. de Boer, “Ultrahigh-resolution high-speed retinal imaging using spectral-domain optical coherence tomography,” Opt. Express12, 2435–2447 (2004). [CrossRef] [PubMed]
S. Witte, M. Baclayon, E. J. Peterman, R. F. Toonen, H. D. Mansvelder, and M. L. Groot, “Single-shot two-dimensional full-range optical coherence tomography achieved by dispersion control,” Opt. Express17, 11335–11349 (2009). [CrossRef] [PubMed]
W. Drexler and J. Fujimoto, Optical Coherence Tomography: Technology and Applications (Springer, 2008). [CrossRef]
A. Fercher, C. Hitzenberger, M. Sticker, R. Zawadzki, B. Karamata, and T. Lasser, “Numerical dispersion compensation for partial coherence interferometry and optical coherence tomography,” Opt. Express9, 610–615 (2001). [CrossRef] [PubMed]
M. Duarte, M. Davenport, M. Wakin, and R. Baraniuk, “Sparse signal detection from incoherent projections,” in Proc. Int. Conf. on Acoustics, Speech and Signal Processing (ICASSP)3, III305–308 (2006).
B. Karamata, K. Hassler, M. Laubscher, and T. Lasser, “Speckle statistics in optical coherence tomography,” J. Opt. Soc. Am. A22, 593–596 (2005). [CrossRef]
P. Cimalla, J. Walther, M. Mehner, M. Cuevas, and E. Koch, “Simultaneous dual-band optical coherence tomography in the spectral domain for high resolution in vivo imaging,” Opt. Express17, 19486–19500 (2009). [CrossRef] [PubMed]
P. D. Woolliams and P. H. Tomlins, “Estimating the resolution of a commercial optical coherence tomography system with limited spatial sampling,” Meas. Sci. Technol.22, 065502 (2011). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Click here to see a list of articles that cite this paper