Open Access
Add to BibSonomyAbstract
An interferometric synthetic aperture microscopy (ISAM) system design with real-time 2D cross-sectional processing is described in detail. The system can acquire, process, and display the ISAM reconstructed images at frame rates of 2.25 frames per second for 512×1024 pixel images. This system provides quantitatively meaningful structural information from previously indistinguishable scattering intensities and provides proof of feasibility for future real-time ISAM systems.
©2008 Optical Society of America
1. Introduction
Interferometric synthetic aperture microscopy [1] (ISAM) is an imaging modality using computed imaging and synthetic aperture techniques [1–5] to extend the capabilities of instrumentation derived from optical coherence tomography [6–9] (OCT) and optical coherence microscopy [10–12] (OCM). These computed imaging techniques are obtained from physics-based scattering models for OCT and formulated via a solution of the inverse scattering problem. One of the principal demonstrated advantages of ISAM over OCT is the ability to resolve features in the sample that are outside of the confocal region. Ultimately, a more quantitatively accurate and faithful representation of the sample structure is computed.
ISAM requires more computation than OCT and reconstructions have previously been computed offline. Currently, OCT systems provide real-time visualization, an ability that is useful for immediate feedback in time critical situations or for monitoring transient dynamics [13–16]. A real-time ISAM system that could provide spatially-invariant resolution has obvious benefits for the biomedical optics community. One real-time capable system produces an extended focus for OCT scans using a Bessel beam [17]. Other numerical methods for OCT have been used to numerically move the focus to areas of interest via a diffraction kernel [18, 19].
ISAM has previously been described and demonstrated with processing rates of roughly one frame per minute on a modern personal computer [1,2]. In this work, an implementation of a real-time ISAM algorithm is described. The critical functionality of the algorithm is identified where optimizations have been made. In addition, adequate compromises to a prototype algorithm are able to further improve the speed. These improvements are described in detail and the negative effects of the compromises are discussed.
The optimized system can acquire, process, and display the cross-sectional ISAM reconstructed image at frame rates of 2.25 frames per second for 512×1024 pixel images on a personal computer with two 3.0 GHz Intel Xeon processors. As previously demonstrated, this system provides quantitatively meaningful structural information from previously indistinguishable scattering intensities and now does so in real-time. The limits to speed now rely on the parallelizability of the processing hardware.
In this paper, the theoretical expression from the signal to the reconstructed object is derived for a prototype algorithm. Then, the prototype algorithm is described in a step-by-step computational procedure. Next, the real-time algorithm is derived where each of the steps is explicitly detailed. Finally, an implementation is described, and experimental results are presented.
2. Analytical model
An object has a susceptibility represented by η(r ∥, z), where r ∥ is the transverse position and z is the longitudinal position. The object may be imaged with a variety of OCT system designs to acquire complex field measurements in time-domain OCT (TD-OCT) S(r ∥, tω) or spectral-domain OCT (SD-OCT) S̃(r ∥, ω), which take the arguments transverse position of the beam r ∥, and time tω or frequency ω. It should be noted here that only real intensities are measured from which the complex representations can be computed unambiguously [20, 21]. To represent the correction of dispersion, a change of variables is introduced from ω to k, the uniformly spaced wavenumber, and from tω to t, time as if waves were propagated in a dispersionless environment. After correcting for dispersion, the TD-OCT signal S(r ∥,t) and the SD-OCT signal S̃(r ∥, k) are related by a Fourier transform with respect to k. Similarly, the 3-D Fourier transform of S(r ∥,t) is S͌(q,k), where q is the transverse wavevector. The analytical coordinates and signals are outlined in Table 1 and Table 2, respectively.
Table 1. Analytic coordinates
The algorithm for ISAM is derived as in references 2 and 3. Explicitly, a relation is derived that associates the Fourier transform of the object η with the Fourier transform of the signal S. The 2-D Fourier transform of the spectral-domain OCT signal S̃(r ∥,k), is given by the expression
where η͌ is the 3-D Fourier transform of η, k is the wavenumber, , z 0 is the fixed position of the beam focus, A(k) is the square root of the power spectral density, α=π/NA, and NA is the numerical aperture of the beam. The corresponding Tikhonov regularized solution, η͌ +, takes the form
where
β is the longitudinal frequency coordinates of the object, and λ is the regularization constant. The term z 0 will be zero when the origin of the z coordinates are at the focal plane. Implementation of this is described in more detail in Step 4 of the prototype algorithm. Thus, rearranging the terms, the pseudo inverse can be rewritten as
where
The resampling step in Eq. (4) corrects the depth-dependent defocus and is crucial for the performance of the algorithm. The filter B͌(q, k) can be highly singular which introduces noise, hence the need for regularization. Furthermore, applying the filter provides a quantitative, but not a significant qualitative change to the data. In place of Eq. (4), it is sensible to compute the unfiltered solution,
for the real-time algorithm. The model in Fig. 1(a) and the grids in Figs. 1(b) and 1(c) describe visually in 2D the ISAM resampling seen in Eq. (4).
Fig. 1. (a). The relation between the spatial frequencies of the signal space and the spatial frequencies in the object space. (b). Sampling lattice for selected (β, |q|) values on a uniform (k, |q|) grid. (c). Sampling lattice for selected (k, |q|) values on a uniform (β, |q|) grid.
Table 2. Analytic signals
3. Prototype Algorithm
The prototype algorithm was designed in Matlab, where all the calculations were done in double precision, 64-bit. Figure 2 shows a data flow diagram for the prototype ISAM algorithm for 2D imaging where in the paraxial limit the 2D coordinates r ∥ are separable in each transverse dimension. The discrete coordinates and signals are described in Tables 3 and 4, respectively. The digitized spectra Srω is a function of the M discrete positions r in the plane perpendicular to the beam axis and the N discrete frequencies ω. To account for dispersion [22], a non-uniform interpolation is needed to map the sampled frequencies ω to the sampled wavenumbers k. Similarly, the ISAM solution requires a non-uniform interpolation mapping wavenumbers k to longitudinal frequency coordinates of the object β. The cubic B-spline [23] is chosen as the non-uniform resampling interpolator; although, a windowed sinc interpolator, such as the prolate-spheroidal interpolator [24], may be chosen to select specific convergence characteristics. Unfortunately, the frequency response for many non-uniform interpolation procedures drops in performance for frequencies higher than half the Nyquist rate. To mitigate this effect, each spectrum is interpolated by performing a fast Fourier transform (FFT), padding with N zeros, and performing an inverse FFT (IFFT) [25]. The interpolated spectra are stored in a buffer with 2N rows and M columns. Thus, the new interpolated, digitized spectra S rω′ is a function of the sampled positions r and the new sampled frequencies ω′. Similarly, Srk is interpolated to twice the Nyquist rate to get S rk′ as a function of r and the new uniformly sampled wavenumbers k′.
The non-uniformly sampled spectrum in spectral-domain OCT can be corrected in a fashion similar to material dispersion correction [22] by resampling the Fourier spectrum from ω′ to k. A desired reindexing array in is based primarily on calibrated, second- and third-order correction factors α 2 and α 3, respectively.
where n is an integer between 0 and N-1, 2N is the number of samples of the interpolated spectrum S rω′, and ω ctr is the calculated centroid of the spectrum on a scale from 0 to 1. α 2 and α 3 are set through a calibration routine. A mirror models a perfect reflector with a minimized width of the point spread function (PSF). Thus, α 2 and α 3 are adjusted such that the imaged PSF is indeed minimized.
Table 3. Discrete coordinates
Table 4. Discrete signals
Figure 2 shows a data flow chart of the prototype algorithm which has been broken down into the eight steps listed below.
Step 1. The spline coefficients and interpolations are computed as described above. The result is stored in an array S rk′ where k′ is the uniformly sampled wavenumber.
Step 2. The 1-D FFT of the columns of S rk′ is taken to get the signal Srt where t is the sampled time delay of the optical signal. The Hermitian symmetric values, where t<0, are removed. There will be no ambiguity, if the reference arm is placed such that the reference light pulse arrives at the spectrometer before the sample light pulses. The complex analytic OCT signal is represented by Srt.
Step 3. A phase correction procedure is implemented to ensure phase stability for the ISAM reconstruction [26]. A glass coverslip can be placed on top of the sample to track and reduce system instabilities. The phase corrected signal is represented as Ŝrt.
Step 4. The contribution of the constant spectral intensity is removed by setting Ŝrt to zero, for values near t=0. Then, a coordinate change from t to t′ is made by circularly shifting the rows such that the focal plane lies at t′=0, which coincides with z′=0 where the solution is derived in Eq. (4). The data S rt′ is padded with 2|t 0| rows of zeros to prevent aliasing of the shifted data, where t 0 is the number of time delay samples from the focus to the center of Ŝrt. Notice that the phase-corrected results are incorporated into S rt′, yet the hat ( ) notation is not carried through. This is done for simplicity of notation because Ŝrt represents Srt when the system produces a phase stable response.
Alternatively, the diffraction kernel , as seen in Eq. (3), may be used after Step 5 to adjust the position of the reconstructed focal plane.
Step 5. The 2-D IFFT of S rt′ is taken to get Sqk.
Step 6. The 2-D ISAM resampling grid, spline coefficients, and the interpolations are calculated. Then the result is multiplied by B͌(q, k) get ηqβ.
Step 7. The 2-D FFT ofηqβ is taken to get η rz′, where z′=0 at the focal plane.
Step 8. A coordinate change from z′ to z is made by circularly shifting the rows such that the reference delay lies at the top of the image and corresponds to the z=0 plane. Then, the magnitude of the ISAM reconstruction ηrz is calculated resulting in ISAM amplitude image |ηrz|.
4. Real-time Algorithm
There are a number of operations hindering the performance of this prototype ISAM algorithm which have been overcome. Using C++ instead of Matlab allows for a number of speed improvements. The 64-bit operations such as the FFT and interpolations can be replaced with 32-bit operations with a small, but tolerable, increase of quantization noise. A speed advantage is gained by replacing the ‘for’ loops within Matlab scripts by vectorized Intel SSE (Streaming SIMD Extentions) commands and/or compiled ‘for’ loops. Time-expensive, built-in Matlab interpolation and FFT functions are replaced with hardware optimized functions as found in the Intel Math Kernel Library (MKL). An FFT of the real spectrum is implemented using the customized real-to-complex FFT in the MKL.
The focal plane is fixed such that t0=0 by imposing a restriction that the focus be placed at the center of the image. Therefore, the complex analytic signal does not need to be padded with any zeros, and thus the computation time of the FFT is reduced because the FFT is always a power of two. While the prototype functions utilized dynamically-allocated memory, the real-time program allocates memory prior to imaging time. Similarly, a table of resampling coefficients for the dispersion compensation and the ISAM reconstruction are pre-calculated and stored in memory. The prototype algorithm is interpolated to twice the Nyquist rate to mitigate the high frequency roll-off of the non-uniform interpolator. Although the interpolation has good frequency response, it is not necessary, especially when speed is an issue. The phase stabilization procedures, which might be needed for long acquisition periods, are not necessary for 2D imaging since this SD-OCT system provides good phase stability over the short duration of the scan, about 17 ms. A good estimate of maximum drift for maintaining phase stability is no more than a half of a wavelength for the length of the synthetic aperture, the depth dependent beam width.
The key computation for ISAM is the resampling in the Fourier space. The new design is an efficient routine which requires two interpolations of the columns, one one-dimensional (1D) real-to-complex (R2C) fast Fourier transform (FFT) of the columns, and two 2D FFTs. The FFT routine from the Intel Math Kernel Library (MKL) was used for all 1D and 2D transforms. The 1D 2048-point R2C FFT is calculated for every column of the 512×2048 real-valued data, while the 2D FFT and 2D IFFT are calculated for the 512×1024 complex values.
The reindexing array in and the corresponding interpolation table is pre-computed and stored in random access memory (RAM) for repeated use. The coefficients of the cubic B-spline are pre-computed. The same type of calculations could be simply done for a table of linear spline coefficients, but the cubic B-spline is described below for completeness. The integer parts of the indices used for calculation of the coefficients to the cubic B-spline are given by
The fractional coefficients are given by
and
Similarly, the Fourier resampling indices of ISAM need to be determined and precalculated. To determine the ISAM resampling indices, the constants which specify the axial and transverse bandwidths of the object are determined. These are based on system parameters and are β min, β max, q min, and q max, respectively. By approximating the longitudinal bandwidth parameters of the object, one can avoid computationally costly zero-padding for calculation of the resampled solution. In this case, β min and β max are well approximated with the range of the optical bandwidth, q min is set to zero, and q max is set to the maximum transverse scanning frequency. There may be a small loss of spectral information across the β min and β max boundary as illustrated by Figs. 1(b) and 1(c), but this has only a marginal effect on the signal-to-noise ratio. The maximum and minimum wavenumbers are calculated for the region of interest,
A range of values for q[m] and β[n] is created in the 2-D Fourier space and the resampling grid kq[m,n] is pre-calculated. Notice here that M and N′=N/2 are assigned according to number of rows and columns in the complex analytic sampled Fourier data.
The kq[m,n] grid is used to calculate a table of values for the cubic B-spline coefficients. The values are calculated as shown above except each column of resampling parameters is different and therefore must also be saved. Figure 1(b) shows the plot of the kq[m,n] grid where the curved lines represent the constant values of β[n]. This 2D grid is used to calculate another table of cubic B-spline coefficients.
Similar to the spline coefficient calculations shown above, the reindexing array kq[m,n] and the corresponding cubic B-spline interpolation coefficient table is pre-computed and stored in random access memory (RAM) for repeated use. The integer part of the index used for calculation of the in the cubic B-spline is given by
The fractional coefficients are given by
and
Using the pre-calculated tables, a flow diagram of the real-time algorithm is shown in Fig. 3.
Here Srω[m,n] is the raw interferometric data captured from the camera and has M columns and N rows. In this implementation, M=512 columns and N=2048 rows.
Fig. 3. Computational flow chart for memory allocation for each step of the inverse scattering algorithm.
Step 1. The pre-calculated table is used to perform the interpolation as follows.
for all integers 0≤n<N and 0≤m<M.
Step 2. The real-to-complex 1-D FFT routine from the Intel Math Kernel Library (MKL) is used on all the columns.
The real-to-complex FFT will compute N/2 complex values. The new number of rows of the complex data is given by N′=N/2.
Step 3. The contribution of the noise from the average spectral intensity on the detector is removed by setting Srt[m, n] equal to zero at the t=0 plane. Also, Srt[m, n] is circularly shifted by half such that the focus will be the new t′=0 plane.
Step 4. The complex 2-D inverse FFT (IFFT) of the complex analytic signal S rt′[m, n] is calculated
Step 5. The pre-calculated table is used to perform the cubic B-spline interpolation as follows.
where the calculated cubic B-spline coefficients are from the lookup table.
Step 6. The complex 2-D FFT of the Fourier transformed object ηqβ[m, n] is calculated
Step 7. η rz′[m, n] is circularly shifted such that the focus is in the middle of the image,
then the magnitude |ηrz[m, n]| is displayed.
5. Implementation
Coherence imaging measurements are made using a spectral-domain OCT system [27,28], Fig. 4. The system can be minimally reconfigured to also collect time-domain (TD) data. A femtosecond laser (Kapteyn-Murnane Laboratories, Boulder, CO) provides the broadband illumination. The bandwidth of the source is 100 nm, with a center wavelength of 800 nm. The field fluctuates too rapidly to be detected directly, thus the optical signal must be measured with an interferometer. The illuminating source is divided by a 50:50 fiber-optic coupler (splitter) for interference measurements. The 90:10 coupler is used only for time-domain OCT reference measurements, which are not used when imaging in the spectral domain. The associated reference diodes are used to remove the source intensity in TD-OCT. The reference arm contains a delay mirror and the sample arm contains a lens (objective) to focus the Gaussian beam into the sample. In the sample arm, the objective has a focal length of 12 mm, a spot size of 5.6 µm, a confocal parameter (depth-of-focus) of 239 µm, and a NA of 0.05. Since the ISAM algorithm is relatively robust to focal plane alignment, the focal plane can be placed at the center of the image either by a rough measurement or with a calibration tissue phantom.
The spectral detector collects light from the sample and reference arms. In the spectrometer, the light is collimated with a 100 mm focal length achromatic lens and dispersed from a blazed gold diffraction grating (53004BK02-460R, Spectra-Physics), with 830.3 grooves per millimeter and a blaze angle of 19.70 degrees for a blaze wavelength of 828 nm The dispersed optical spectrum is focused using a pair of achromatic lenses each having a focal length of 300 mm. The focused light is incident upon a line-scan camera (P2-2x-02K40, Dalsa.) which contains a 2048-element charge-coupled device (CCD) linear array with 10×10 µm detection elements. The spectral resolution provides 2 mm of range depth per axial scan. The camera operates at a readout rate of over 29 kHz, thus one axial scan can be captured during an exposure interval of about 34 µs.
The tissue phantom was placed on a stage and the 2D image data were acquired by scanning the incident beam in the transverse plane using a pair of computer-controlled galvanometer-scanning mirrors. The frame rate depends on the number of axial scans acquired per image or volume. In our experiment, a series of spectrally resolved images (500×2048 pixels) was each acquired in 17 ms. A frame capture card (PCI-1428, National Instruments, Inc.) in a computer collects camera data and writes it to the RAM. The frame capture card receives an external trigger from a galvanometer controller card (PCI-6711, National Instruments, Inc.) which activates the frame acquisition. The computer CPUs are dual Intel Xeon processors operating at 3.0 GHz each.
A tissue phantom consisting of TiO2 scatterers suspended in silicone with an average diameter of 1 µm is imaged. Figure 5 shows a transverse scan of 0.4 mm where the real-time cross-sectional OCT image is on the left and the cross-sectional ISAM reconstruction is on the right. The scattering effects of the beam profile in this tissue phantom are evident and shown in Fig. 5 on the left. The previously indistinguishable scatterers are easily identified in the reconstruction, shown in Fig. 5 on the right. The reconstruction better represents the tissue phantom and thus, provides better visualization and a potential for more accurate tissue classification and disease diagnosis. Note that the image boundaries on the left and right, near the top of the phantom, still have some blurring. The solution to the inverse scattering problem is underdetermined for those particular scatterers in these regions because those scatterers were not measured by the full aperture of the beam. The real-time ISAM system provides quantitatively meaningful structural information from previously indistinguishable scattering intensities and provides proof of feasibility for future real-time systems. The acquisition speed has improved from 45 seconds per image to under 450 milliseconds per image, a factor of over one-hundred. Real-time translation of the sample is shown in an online movie (movie1.avi, 6.1 MB).
When imaging 3D volumes, the real-time 2D ISAM processed measurements can be used as a partial solution to the 3D reconstruction. By computationally rotating the data in the transverse direction and re-running the algorithm in the slow scan direction, the full 3D reconstruction can be computed. Figure 6 shows an en face area (0.4 mm×0.4 mm) of excised human breast tumor tissue imaged with OCT (left) and ISAM processed with the real-time algorithm (right), following post-processing with the same algorithm in the slow scanning direction. The en face slice shown is 450 µm above the focus, such that the cellular and nuclear features in the OCT image are not apparent. Notice that the corresponding features become apparent in the ISAM processed data. A representative histological slice of the tissue is shown in Fig. 7 for comparison. It should be noted that at increasing depths, the probability of multiple scattering increases. Multiple scattering obviates the fundamental assumption of single scattering in ISAM, as it does for OCT. The real-time 2D cross-sectional acquisition movie (movie2.avi, 6.1 MB) and a representative en face movie (movie3.avi, 15 MB) are available online.
Fig. 5. Tissue phantom imaged with real-time inverse scattering OCT system. The original OCT image (left) and the inverse scattering solution (right) are computed in real-time. Real-time translation of the sample stage is shown in the associated movie (movie1.avi, 6.1 MB). [Media 1]
Fig. 6. Human breast tumor tissue imaged with OCT (left) and with the real-time 2D ISAM system (right), and post-processed in slow scan direction. The en face slice shown in the plane of the page is 450 µm above the focus. The volumetric ISAM data shows cell and nuclei data not observed in the volumetric OCT data. The scale bar is 100 µm. A portion of the real-time 2D ISAM acquisition (movie2.avi, 6.1 MB) and a representative en face fly-through (movie3.avi, 15 MB) are available online. [Media 2] [Media 3]
Fig. 7. Corresponding histological section of human breast tumor tissue showing comparable features with Fig. 6. The scale bar represents 100 µm.
6. Conclusion
Through modifications of a reconstruction algorithm and hardware upgrades, we have demonstrated the feasibility for real-time imaging of the inverse scattering solution. There are a number of modifications made to the algorithm which compromise certain aspects of the processing for speed. The largest drawback appears to be the degradation of the frequency response for the cubic B-spline interpolator, which decreases the SNR for high frequency interpolations. As described in this paper, a variety of interpolators may be used to gauge the tradeoffs in speed versus performance.
The system described achieves spatially invariant resolution for all depths in a cross-sectional scan. Such a system has a potential to provide important cellular scattering information when used to image biological tissues, particularly for real-time, ISAM image-guided surgical interventions. Furthermore, existing OCT systems may be modified to produce similar results with only minor modifications to the processing algorithm. The limits to processing speed now depend largely on the parallelizability of the processing hardware. The parallel nature of this ISAM algorithm suggests that specialized computing hardware such as Graphics Processing Units (GPU) may offer further speed advantages.
Acknowledgments
This work was supported in part by the National Institutes of Health (Roadmap Initiative/NIBIB, 1 R21 EB005321, and NIBIB, 1 R01 EB005221, S.A.B.), the National Science Foundation (CAREER Award, BES 03-47747, and BES 05-19920, and BES 06-19257, S.A.B., CAREER Award, 0239265, P.S.C.), and the Grainger Foundation (S.A.B.). Human tissue was acquired under Institutional Review Board protocols approved by the University of Illinois at Urbana-Champaign and Carle Foundation Hospital. Additional information can be found at http://www.biophotonics.uiuc.edu.
1. T. S. Ralston, D. L. Marks, P. S. Carney, and S. A. Boppart, “Interferometric synthetic aperture microscopy,” Nat. Phys. 3, 129–134 (2007). [CrossRef]
2. T. S. Ralston, D. L. Marks, P. S. Carney, and S. A. Boppart, “Inverse scattering problem for optical coherence tomography,” J. Opt. Soc. Am. A 23, 1027–1037 (2006). [CrossRef]
3. T. S. Ralston, D. L. Marks, S. A. Boppart, and P. S. Carney, “Inverse scattering for high-resolution interferometric microscopy,” Opt. Lett. 31, 3585–3587 (2006). [CrossRef] [PubMed]
4. D. L. Marks, T. S. Ralston, P. S. Carney, and S. A. Boppart, “Inverse scattering for rotationally scanned optical coherence tomography,” J. Opt. Soc. Am. A 23, 2433–2439 (2006). [CrossRef]
5. D. L. Marks, T. S. Ralston, P. S. Carney, and S. A. Boppart, “Inverse scattering for frequency-scanned full-field optical coherence tomography,” J. Opt. Soc. Am. A 24, 129–134 (2007). [CrossRef]
6. D. Huang, E. A. Swanson, C. P. Lin, J. S. Schuman, W. G. Stinson, W. Chang, M. R. Hee, T. Flotte, K. Gregory, C. A. Puliafito, and J. G. Fujimoto, “Optical Coherence Tomography,” Science 254(5035), 1178–1181 (1991). [CrossRef] [PubMed]
7. S. A. Boppart, B. E. Bouma, C. Pitris, J. F. Southern, M. E. Brezinski, and J. G. Fujimoto, “In vivo cellular optical coherence tomography imaging,” Nat. Med. 4, 861–864 (1998). [CrossRef] [PubMed]
8. J. M. Schmitt, “Optical coherence tomography (OCT): A review,” IEEE J. Select. Topics Quantum Electron. , 5, 1205–1215 (1999). [CrossRef]
9. B. E. Bouma and G. J. Tearney, The Handbook of Optical Coherence Tomography, (Marcel Dekker, New York, 2002).
10. J. A. Izatt, M. R. Hee, G. M. Owen, E. A. Swanson, and J. G. Fujimoto, “Optical coherence microscopy in scattering media,” Opt. Lett. 19, 590–592 (1994). [CrossRef] [PubMed]
11. J. M. Schmitt, M. J. Yadlowsky, and R. F. Bonner, “Subsurface imaging of living skin with optical coherence microscopy,” Dermatology 191, 93–98 (1995). [CrossRef] [PubMed]
12. J. A. Izatt, H.-W. Kulkarni, K. Wang, M. W. Kobayashi, and M. W. Sivak, “Optical coherence tomography and microscopy in gastrointestinal tissues,” IEEE J. Sel. Tops. Quantum Electron. 2, 1017–1028 (1996). [CrossRef]
13. 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]
14. A. M. Rollins, S. Yazdanfar, J. K. Barton, and J. A. Izatt, “Real-time in vivo color Doppler optical coherence tomography,” J. Biomed. Opt. 7, 123–129 (2002). [CrossRef] [PubMed]
15. R. Leitgeb, L. Schmetterer, W. Drexler, A. Fercher, R. Zawadzki, and T. Bajraszewski, “Real-time assessment of retinal blood flow with ultrafast acquisition by color Doppler Fourier domain optical coherence tomography,” Opt. Express 11, 3116–3121 (2002). [CrossRef]
16. X. Liu, M. J. Cobb, Y. Chen, M. B. Kimmey, and X. Li, “Rapid-scanning forward-imaging miniature endoscope for real-time optical coherence tomography,” Opt. Lett. 29, 1763–1765 (2004). [CrossRef] [PubMed]
17. R. A. Leitgeb, M. Villiger, A. H. Bachmann, L. Steinmann, and T. Lasser, “Extended focus depth for Fourier domain optical coherence microscopy,” Opt. Lett. 31, 2450–2452 (2006). [CrossRef] [PubMed]
18. Y. Yasuno, J. Sugisaka, Y. Sando, Y. Nakamura, S. Makita, M. Itoh, and T. Yatagai, “Non-iterative numerical method for laterally superresolving Fourier domain optical coherence tomography,” Opt. Express 14, 1006–1020 (2006). [CrossRef] [PubMed]
19. Y. Nakamura, J. Sugisaka, Y. Sando, T. Endo, M. Itoh, T. Yatagai, and Y. Yasuno, “Complex Numerical Processing for In-Focus Line-Field Spectral-Domain Optical Coherence Tomography,” Jpn. J. Appl. Phys. 46, 1774–1778 (2007). [CrossRef]
20. B. J. Davis, S. C. Schlachter, D. L. Marks, T. S. Ralston, S. A. Boppart, and P. S. Carney, “Nonparaxial vector-field modeling of optical coherence tomography and interferometric synthetic aperture microscopy,” J. Opt. Soc. Am. A 24, 2527–2542 (2007). [CrossRef]
21. B. J. Davis, T. S. Ralston, D. L. Marks, S. A. Boppart, and P. S. Carney, “Autocorrelation artifacts in optical coherence tomography and interferometric synthetic aperture microscopy,” Opt. Lett. 32, 1441–1443 (2007). [CrossRef] [PubMed]
22. D. L. Marks, A. L. Oldenburg, J. J. Reynolds, and S. A. Boppart, “Autofocus algorithm for dispersion correction in optical coherence tomography,” Appl. Opt. 42, 3038–3046 (2003). [CrossRef] [PubMed]
23. C. Pozrikidis, Numerical Computation in Science and Engineering (Oxford University Press, New York Oxford, 1998).
24. J. J. Knab, “Interpolation of band-limited functions using the approximate prolate series,” IEEE Trans. Inf. Theory IT-25, 717–720 (1979). [CrossRef]
25. L. P. Yaroslavsky, “Efficient algorithm for discrete sinc interpolation,” Appl. Opt. 36, 460–463 (1997). [CrossRef] [PubMed]
26. T. S. Ralston, D. L. Marks, P. S. Carney, and S. A. Boppart, “Phase stability technique for inverse scattering in optical coherence tomography,” 3rd IEEE International Symposium on Biomedical Imaging: Nano to Macro, 578–581 (2006).
27. L. Lepetit, G. Chériaux, and M. Joffre, “Linear techniques of phase measurement by femtosecond spectral interferometry for applications in spectroscopy,” J. Opt. Soc. Am. B 12, 2467–2474 (1995). [CrossRef]
28. B. Cense, N. A. Nassif, T. C. Chen, M. C. Pierce, S. H. Yun, B. H. Park, B. E. Bouma, G. J. Tearney, and J. F. de Boer, “Ultrahigh-resolution high-speed retinal imaging using spectral-domain optical coherence tomography,” Opt. Express 12, 2435–2447 (2004). [CrossRef] [PubMed]
