Numerical method for axial motion artifact correction in retinal spectral-domain optical coherence tomography

Abstract

A numerical method that compensates image distortions caused by random fluctuations of the distance to an object in spectral-domain optical coherence tomography (SD OCT) has been proposed and verified experimentally. The proposed method is based on the analysis of the phase shifts between adjacent scans that are caused by micrometer-scale displacements and the subsequent compensation for the displacements through phase-frequency correction in the spectral space. The efficiency of the method is demonstrated in model experiments with harmonic and random movements of a scattering object as well as during in vivo imaging of the retina of the human eye.

Keywords

optical coherence tomography (OCT) / motion artifact correction / retinal imaging / numerical method

Introduction

In this paper, we describe a numerical method that compensates image distortions caused by random fluctuations of the distance to an object in spectral-domain optical coherence tomography (SD OCT) imaging. SD OCT [¹] is primarily used for in-vivo, non-invasive imaging of the internal structure of tissues with a high spatial resolution of a few microns. It is based on the spectral interferometric detection of low-coherent infrared light backscattered from the internal inhomogeneities of the object being studied and successive mathematical processing of the captured optical spectrum of the interferometric signal. In ophthalmologic applications of OCT imaging, a typical problem is the undesirable axial microscale movement of the eye relative to the optical detection system. Even with a relatively high imaging speed in SD OCT (more than 20000 A-scans per second [²]), it takes several seconds to generate a 3D image of the retinal region. Under these conditions, avoiding random involuntary movement of the eye is almost impossible. This eye movement is mainly associated with muscle tremor [³], respiration [⁴], and tissue vibration caused by blood flow in large vessels [³,⁵] and have a component directed along the probing beam propagation. The presence of relative axial displacements introduces distortions in the images, which can sometimes adversely affect the interpretation of the obtained diagnostic information.

The solution for this problem requires the development of methods for suppressing the influence of the relative motion of the object during study. These methods can be divided into hardware methods, in which special tools are used to reduce or measure the displacements during scanning, and software methods, in which a previously acquired image is corrected. The most effective method is to increase the imaging speed [⁶]; however, in a real environment, this has the limitations often associated with the threshold limit of radiation exposure and an increased signal-to-noise ratio (SNR). The SNR for excess noise is proportional to the square root of the exposure time [⁷]; therefore, increasing the imaging speed results in increased levels of excess noise. This results in the decrease of the SNR starting from the threshold value of the exposure time caused by the level of excess noise exceeding the level of shot noise [⁷]. In turn, the SNR for shot noise is proportional to the square root of the number of photons registered by the photo receiver. Therefore, to achieve the required dynamic range in the OCT system, the total number of registered photons should be maintained while decreasing the exposure value. The latter requires increasing the total power of the probing light, which is limited owing to the stipulations of the applicable safety standards.

A number of approaches have been developed to compensate for the influence of the displacement of an object during ophthalmic OCT imaging. These approaches are aimed at compensating for the effects of the natural movement of the eye. The specificity of the object under study (retinal structure) implies the predominant influence of transverse displacements, which violate the integrity of en-face images and can be compensated for by either analyzing the full 3D image [⁸–¹¹] or performing additional reference measurements, such as serial scanning in orthogonal directions using a single B-scan [¹²] or an additional 3D image [¹³,¹⁴]. In addition, a number of hardware methods for the correction of eye movement are known and have been discussed in detail in a review paper [¹⁵].

The correction of axial displacements is a separate task that requires the development of fundamentally different approaches. In Refs. [¹⁶–¹⁸], for instance, when constructing angiographic images, eye movement was corrected via correlation analysis of the images themselves. However, the results of the recovery profile of the retina were not provided in these studies. In Ref. [¹⁹], the retinal displacement is determined using a similar technology that involves identifying the three-dimensional correlation of the speckles observed in OCT images. Using this method, a good agreement between the external displacement of an object and amount of displacement recovered from the OCT image was demonstrated. Similarly, the correction of physiologic movements (breathing, heartbeats) was conducted in Ref. [²⁰]. In Ref. [⁵], the presence of axial displacements in the OCT images of the retina and other human tissues was used to determine the frequency and profile of heart contractions. However, the displacements were not compensated in these OCT images. In Refs. [²¹,²²], external hardware tracking was used to correct the arbitrary movements of an object, making it possible the identify the movements and generate a correction signal. In Ref. [¹⁰], a complex mechanism was proposed to manage the manifestations of axial displacements. This mechanism involves the creation of scanning beams with a certain spectral-geometric profile, thus making it possible to accurately distinguish between geometry-induced changes in the image profile and changes caused by the axial movements of the object. In Ref. [²³], an approach in which the correction of axial movements is based on studying the local surface curvature by statistical methods within a window covering 10–80 image lines in directions orthogonal to the scanning direction was described.

In this paper, we propose a numerical method for eliminating the influence of fluctuations of the distance to an object in an SD OCT system. This method uses of SD OCT data to obtain the information of micrometer-scale displacements that occur in the time interval between individual in-depth scans. Then this information is used to remove artifacts from the image reconstructed from the same data. The method does not involve any additional measurements.

Signal generation in SD OCT imaging

To describe the proposed correction method, we first discuss the principles of signal generation in OCT. A schematic of the SD OCT system is shown in Fig. 1. Light from a broadband light source is directed to an interferometer with a reference arm and sample arm. The sample arm contains an optical system that conducts the transverse scanning of the object using a probe beam and detection of the backscattered light. The reference and scattered waves, the fields of which are hereinafter referred to as {{\displaystyle E}}_{\mathrm{R}} and {{\displaystyle E}}_{\mathrm{S}}, respectively, are delivered to a spectrometer, where the sum of the fields is decomposed into spectral components using a diffraction grating and detected using a line scan sensor. In our experimental setup, we use a prism corrector in a spectrometer to ensure the equidistant distribution of the spectral components of the interference signal over the optical frequency w on the line scan sensor [²⁴–²⁶]. The spectral intensity signal of the detected radiation is transmitted from the spectrometer output via a specialized interface [²⁷] to a personal computer for further mathematical processing. In the case of a single scatterer in the measuring arm of the interferometer, the signal \widehat{I}\left(\omega \right) can be represented as follows:

\widehat{I} \left(\omega \right)\propto {\left|{{\displaystyle E}}_{\mathrm{S}}\left(\omega \right)\right|}^{\mathrm{2}}+{\left|{{\displaystyle E}}_{\mathrm{R}}\left(\omega \right)\right|}^{\mathrm{2}}

(1)+\mathrm{2}\left|{{\displaystyle E}}_{\mathrm{S}}\left(\omega \right)\right|\left|{{\displaystyle E}}_{\mathrm{R}}\left(\omega \right)\right|\cos \left(\mathrm{2}{\displaystyle \frac{\omega }{c}}z\right),

where z is the value of the optical path difference for the reference and scattered fields, ω is the optical frequency, c is the speed of light, and the hat notation denotes the initial data. This expression can be rewritten in a more compact form as follows:

(2)\widehat{I}\left(\omega \right)={{\displaystyle S}}_{\mathrm{A}\mathrm{C}}\left(\omega \right)+{{\displaystyle S}}_{\mathrm{C}\mathrm{C}}\left(\omega \right)\cos \left(\mathrm{2}{\displaystyle \frac{\omega }{c}}z\right),

in which the autocorrelation {{\displaystyle S}}_{\mathrm{A}\mathrm{C}}\left(\omega \right) and cross-correlation {{\displaystyle S}}_{\mathrm{C}\mathrm{C}}\left(\omega \right)\cos \left(\mathrm{2}{\displaystyle \frac{\omega }{c}}z\right) components can be distinguished [²⁸]. The {{\displaystyle S}}_{\mathrm{A}\mathrm{C}}\left(\omega \right) term may be eliminated in several ways [²⁹–³²]; however, this is outside the scope of this paper. The cross-correlation component {{\displaystyle S}}_{\mathrm{C}\mathrm{C}}\left(\omega \right)\cos \left(\mathrm{2}{\displaystyle \frac{\omega }{c}}z\right) causes the registering spectral comb with a frequency (\mathrm{2}{\displaystyle \frac{z}{c}}), which is proportional to the optical path difference between the reference and scattered waves. The cosine argument of the {{\displaystyle S}}_{\mathrm{C}\mathrm{C}}\left(\omega \right)\cos \left(\mathrm{2}{\displaystyle \frac{\omega }{c}}z\right) component depends on the distance to the scatterer and thereby reflects the fluctuations in the distance to the object. The depth (z-coordinate) of the scatterer inside the object and its backscattering coefficient can be judged by the result of the Fourier transform of the spectrometer’s response \widehat{I}\left(\omega \right). An in-depth scattering profile (A-scan) is displayed as brightness. A series of A-scans comprises 2D or 3D OCT images. A 2D dataset, or a B-scan, is acquired by moving a probe beam along the surface of the object (x-coordinate) and synchronously recording A-scans. 3D OCT images are obtained by generating a set of B-scans recorded with a successive offset along the y-coordinate, which is orthogonal to the scanning direction in the B-scan (Fig. 2).

In the present study, raster scanning is used for experimental testing, as illustrated in Fig. 2(a). To save time and ensure a continuous scanning process, the beam is moved in opposite directions along the x-coordinate in even and odd B-scans. To ensure correlation between adjacent samples, a scanning step that is smaller than the diameter of the probe beam should be chosen, as illustrated in Fig. 2(e).

Numerical displacement compensation method

Figure 3(a) presents an example of a distorted image section across the slow axis (y-axis in Fig. 2(a)), which is typical for retinal imaging. The dataset used for Fig. 3 was obtained using of an experimental SD OCT system at a 1060-nm wavelength on a healthy volunteer eye. The profile of the retina appears heavily rugged owing to a physiologic tremor. The irregular displacement range of this tremor may reach hundreds of wavelengths (see experimental example in Fig. 3(b) calculated from the dataset using the proposed method).

The initial data of the A-scan \widehat{I}\left(\omega \right) after exclusion of the autocorrelation component can be represented as a discrete sequence I_w, where w is the discrete number of the spectral samples corresponding to the optical frequency ω. The number of elements in the discrete sequence I_w is equal to the number of photocells in the line scan sensor N. A 2D dataset obtained via transverse scanning is denoted I_w,x, where x is the discrete value of the transverse coordinate. After correcting the spectrum and dispersion [³³,³⁴], a two-dimensional set of complex data is calculated as an analytical signal using the numerical Hilbert transform: {{\displaystyle I^{\text{'}}}}_{w,x}={{\displaystyle Scc}}_{w,x}^{\text{'}}{\mathrm{e}}^{{{\displaystyle \mathrm{\phi }}}_{w,x}^{\text{'}}}+\mathrm{i}H\left({{\displaystyle Scc}}_{w,x}^{\text{'}}{\mathrm{e}}^{{{\displaystyle \mathrm{\phi }}}_{w,x}^{\text{'}}}\right)\left(w\right), where {{\displaystyle Scc}}_{w,x}^{\text{'}} and {{\displaystyle \mathrm{\phi }}}_{w,x}^{\text{'}} are the corrected spectrum and phase distribution in discrete space, respectively. In this paper, the numerical Hilbert transform is used as the simplest method to obtain complex-valued spectral data. Additional methods to achieve this are described in Ref. [²⁸]. The complex dataset {{\displaystyle F}}_{z,x} obtained by the Fourier transform {{\displaystyle F}}_{z,x}=\mathfrak{F}\left({{\displaystyle I}}_{w,x}^{\text{'}}\right)\left(z\right) makes it possible to analyze the relative displacement of adjacent columns.

Determining whether adjacent A-scans have undergone vertical displacement relative to each other is possible using the total phase difference for each element {{\displaystyle F}}_{z,x} and {{\displaystyle F}}_{z,x+\mathrm{1}}. This operation can be represented as follows:

(3){{\displaystyle \mathrm{\Delta }\varphi }}_x=\arg \left({{\displaystyle \sum }}_{z=\mathrm{0}}^{N/\mathrm{2}}{{\displaystyle F}}_{z,x}{{\displaystyle F}}_{z,x+\mathrm{1}}^{*}\right),

where operator arg denotes the argument function, and the asterisk denotes the complex conjugate of {{\displaystyle F}}_{z,x+\mathrm{1}}^{\text{'}}.

The argument of the sum of complex numbers {{\displaystyle \mathrm{\Delta }\varphi }}_x (Eq. (3)) represents the change in the phase of the interference signals due to displacement and is resistant to noise. The key to the method’s robustness is the requirement of a correlation between successively registered A-scans. The object displacement dz occurring in the time gap between the A-scans results in the phase shift dϕ = 2dz⋅ω/c for all scatterers. After the Fourier transformation, this spectral phase shift becomes the additional term for the phase factors of A-scans F_z,x–ϕ_x. Because successively registered A-scans overlap, as illustrated in Fig. 2(e), the intensity overlap integral for Gaussian beams is approximately 80%. Thus, all the terms in Eq. (3) provide equal-phase (i.e., caused by the movement of the entire object) deposits to the sum. The influence of non-overlapping areas is small enough to introduce distortions into the phase factor of the total {{\displaystyle \mathrm{\Delta }\varphi }}_x in Eq. (3). Similarly, the noise terms do not affect the phase factor owing to the extremely low value of noise in relation to the signal (at least 20–40 dB). In addition, the impact of blood vessels on the change in the overall phase factor is negligible owing to the combination of the small scattering signal from blood cells and the small geometric size of the vessels compared with the total depth of imaging.

A similar approach was used for, in particular, the phase alignment of a dataset to obtain angiographic information [³⁵,³⁶]. The results of this analysis can be used to correct the axial displacement of A-scans sequentially throughout the entire dataset via the phase correction of the complex signal distribution in the Fourier domain:

(4){{\displaystyle I}}_{w,x}^{\text{'}\hspace{1em}\text{'}}={{\displaystyle I}}_{w,x}^{\text{'}}{\mathrm{e}}^{-\mathrm{i}{\varphi }_x\kappa \left(w\right)},

where the cumulative phase correction factor ϕ_x is determined by the expression {{\displaystyle \varphi }}_x={{\displaystyle \varphi }}_{x-\mathrm{1}}+\mathrm{\Delta }{{\displaystyle \varphi }}_x, ϕ₀ = 0, and k(w) is the dimensionless wavenumber normalized by the central wavenumber of registered spectrum. Equation (4) subtracts the motion-generated phase shift from A-scans and corrects the motion influence on the OCT data. For example, the displacement {{\displaystyle \varphi }}_x for the initial data shown in Fig. 3(a) is illustrated in Fig. 3(b) as a function of the A-scan number in phase units, which can be considered the x-coordinate for Eq. (4).

The use of the cumulative sum for calculating object movement may lead to the possibility of accumulation of the cumulative error caused by, for instance, a violation of the applicability of the method for the shift speed. Such errors are not critical in the case of single short-term effects because they lead to easily detectable image distortions, such as visual surface tearing (step), which can be taken into account during further image processing and analysis. The location of discontinuities in the image can be corrected after detection by other methods, such as correlation processing. In the case of a systematic violation of the conditions required for the applicability of the method, the cumulative error leads to performance loss. However, it can be identified by analyzing the dynamics of the cumulative phase correction factor. The latter may be used to generate a real-time warning to the user or automation to conduct the correction of research.

The newly calculated signal {{\displaystyle F}}_{z,x}^{\text{'}}=\mathfrak{F}\left({{\displaystyle I}}_{w,x}^{\text{'}\hspace{1em}\text{'}}\right)\left(z\right) will be corrected in the vertical direction, as presented in Fig. 3(c), where the same cross-section is presented after the implementation of the drift correction procedure. It should be noted that axial movements also cause image distortions in en-face sections; these may appear as XY-motions but are not if the surface is tilted to the probe beam direction, as illustrated in Fig. 3(d). These distortions are also successfully eliminated by implementing the proposed method. Figure 3(e) presents a restored en-face image of the same section as in Fig. 3(d).

It should be noted that the displacement may be correctly calculated only if it is sufficiently small; that is, in the time interval between adjacent A-scans, if the object under study shifts by no more than a quarter of the central wavelength of the probe radiation (2dz⋅ω/c<π/2), and the scans themselves are partially correlated (they should be at least partially overlapped, as illustrated in Fig. 2(e)). In the laboratory OCT systems used, the imaging speed is 20000 A-scans per second, and the central wavelength of the probe radiation is 1300 and 1060 nm for the model experiments and eye imaging, respectively. These parameters allow for the unambiguous determination of the displacements occurring with a longitudinal velocity of up to 6.5 and 4.5 mm/s for each of the aforementioned setups.

Further, it should be borne in mind that the procedure does not add new data to the image and the missing data (top and bottom regions of the image) are generally filled in incorrectly, making it necessary to reduce the imaging depth by the value of the maximum displacement.

A peculiarity of the proposed method is the absence of calculated displacement in the case of no motion. This signifies that surface shape of an object has no influence on the result of the procedure (unlike in intensity correlation methods); therefore, the method does not tend to flatten the object surface.

When constructing 3D images, this procedure is performed sequentially for all B-scans. In this case, the use of scanning pattern which contains both forward and backward directions (Fig. 2(a)) allows for the continuous analysis of consecutive A-scans in the image. The proposed algorithm can be applied to other continuous scanning patterns. In particular, it can be used for Lissajous scanning using resonant scanners [¹¹,³⁷–⁴⁰] (Fig. 2(c)), cyclic scanning in optic-nerve examination [⁴¹,⁴²] (Fig. 2(b)), and the endoscopic assessment of coronary vessel walls [⁴³] (Fig. 2(d)).

Experimental testing of proposed method

Figure 4 presents the experimental setup used to test the proposed compensation method. An experimental SD OCT system with an imaging speed of 20000 A-scans per second was used. The central wavelength of the probe radiation in this system was 1300 nm.

The sample used in the test was an OCT phantom, which is a flat plastic plate (BioMimic, Canada) with scattering and absorption characteristics similar to those of biotissue. The sample was fixed on the mounting surface of the vibration exciter (VEB Robotron–Messelektronik ESE 201) to simulate the axial displacements of the sample relative to the OCT system. To demonstrate the ability of the method to correct distortions in the surface image caused by the object’s motion while keeping the shape of the surface itself unchanged, the phantom was placed at a small angle relative to the incidence direction of the probing radiation (z-axis in Fig. 4). Figure 5(a) illustrates the projection of a 3D image of the sample surface that is at rest relative to the probe. Hereinafter, all projections of the 3D image are presented using a modified maximum intensity projection algorithm [⁴⁴].

To simulate smooth mutual displacements, the vibration exciter produced slow harmonic oscillations with a frequency of 0.06 Hz and an amplitude of 0.4 mm. The results of the usual 3D OCT image synthesis in this case are presented in Fig. 5(b). Here a sinusoidal-like change in the image profile along the slow scan coordinate (y, as displayed in Fig. 5(a)) is clearly visible. Figure 5(c) presents the results of 3D OCT image synthesis using the proposed method for compensating for axial displacements. It can be seen that the surface profile along the slow scan axis smoothens and becomes identical to the profile presented in Fig. 5(a). It should be noted that the algorithm used to correct the axial displacement does not introduce distortions in the restored profile of the object; it remains inclined.

In the next experiment, the vibration exciter produced chaotic movements in the sample with the center frequency of the control signal spectrum of the order of several Hz. The phantom surface appears broken in Fig. 5(d), as if covered with wrinkles of arbitrary depths. The results of image recovery after compensation for the displacements by the proposed method are presented in Fig. 5(e). The reconstructed surface appears flat again, as in Figs. 5(a) and 5(c), although at some points the instantaneous velocity of the platform reached 6 mm/s, which is almost the threshold value for the experimental setup (6.5 mm/s). A comparison between the induced displacements and those calculated using the proposed method in both model experiments confirms the reliability of the proposed method.

Figure 6 presents an example of a practical application of the proposed method for compensating for random axial displacements in human retinal imaging. An experimental SD OCT system with an imaging speed of 20000 A-scans per second was used. The central wavelength of the probe radiation in this system was 1060 nm.

This study was conducted in vivo using a healthy 37-year-old volunteer as a subject. The volunteer held his breath for the duration of the study. Pulse movement and some irregular tremors were likely the cause of motion distortions (Fig. 6(a)). No intentional movements were introduced in the experiment. After applying the proposed method, the reconstructed surface became smooth in the image (Fig. 6(b)) and visible folds disappeared from all angles. The general view of the surface of the retina corresponds to the canonical images known from the OCT literature. These images were obtained from a previously acquired OCT data. The proposed method was induced in the home-made real-time OCT acquisition software (coded in C++). The application of calculation procedure did not reveal a noticeable increase in image calculation time (approximately 2%–4%).

Conclusion

In retinal OCT imaging, the presence of random eye movement caused by tremor, pulse, and breathing leads to distortions in the OCT images; these distortions can prevent the acquisition of diagnostic information.

In this study, a numerical method for compensating for the influence of mutual axial displacements between the OCT system and object, leading to fluctuations in the measured distance to the object, has been proposed and experimentally verified. The numerical compensation method requires no additional measurements.

The proposed method is based on the analysis of micrometer-scale phase shifts occurring in the time interval between adjacent A-scans and subsequent phase-frequency correction of the shifts in the spectral space. In addition, the method does not tend to flatten the object surface.

The efficiency of the proposed method was demonstrated in model experiments with harmonic and random movements of the scattering object and during the in vivo imaging of the retina.

At the imaging speed of 20000 A-scans per second used in the experiment and a 1300-nm central wavelength of probe radiation, the proposed method leads to the effective compensation for displacements occurring at a longitudinal speed of up to 6.5 mm/s. When a light source with a wavelength of 1060 nm is used, this method is capable of compensating for displacements occurring at a longitudinal velocity of up to 4.5 mm/s.

The simplicity of the method’s calculations makes it feasible to incorporate it into the asynchronous parallel processing code described in Ref. [⁴⁵] for real-time motion artifact correction and interactive visualization.

References

[1] Fercher A F, Hitzenberger C K, Kamp G, El-Zaiat S Y. Measurement of intraocular distances by backscattering spectral interferometry. Optics Communications, 1995, 117(1–2): 43–48

[2] Gelikonov V M, Gelikonov G V, Terpelov D A, Shilyagin P A. Electronic interface systems for goals of spectral domain optical coherence tomography. Instruments and Experimental Techniques, 2012, 55(3): 392–398

[3] Rajabi H, Zirak A. Speckle noise reduction and motion artifact correction based on modified statistical parameters estimation in OCT images. Biomedical Physics & Engineering Express, 2016, 2(3): e035012

[4] Kang W, Wang H, Wang Z, Jenkins M W, Isenberg G A, Chak A, Rollins A M. Motion artifacts associated with in vivo endoscopic OCT images of the esophagus. Optics Express, 2011, 19(21): 20722–20735

[5] de Kinkelder R, Kalkman J, Faber D J, Schraa O, Kok P H B, Verbraak F D, van Leeuwen T G. Heartbeat-induced axial motion artifacts in optical coherence tomography measurements of the retina. Investigative Ophthalmology & Visual Science, 2011, 52(6): 3908–3913

[6] Zawadzki R J, Miller D T. Retinal AO OCT. In: Drexler W, Fujimoto J G, eds. Optical Coherence Tomography: Technology and Applications. 2nd ed. Switzerland: Springer International Publishing, 2015, 1849–1920

[7] Gelikonov V M, Gelikonov G V, Shilyagin P A. Optimization of Fizeau-based optical coherence tomography with a reference michelson interferometer. Bulletin of the Russian Academy of Sciences. Physics, 2008, 72(1): 93–97

[8] Kraus M F, Potsaid B, Mayer M A, Bock R, Baumann B, Liu J J, Hornegger J, Fujimoto J G. Motion correction in optical coherence tomography volumes on a per A-scan basis using orthogonal scan patterns. Biomedical Optics Express, 2012, 3(6): 1182–1199

[9] Kraus M F, Liu J J, Schottenhamml J, Chen C L, Budai A, Branchini L, Ko T, Ishikawa H, Wollstein G, Schuman J, Duker J S, Fujimoto J G, Hornegger J. Quantitative 3D-OCT motion correction with tilt and illumination correction, robust similarity measure and regularization. Biomedical Optics Express, 2014, 5(8): 2591–2613

[10] Chen Z, Shen Y, Bao W, Li P, Wang X, Ding Z. Motion correction using overlapped data correlation based on a spatial-spectral encoded parallel optical coherence tomography. Optics Express, 2017, 25(6): 7069–7083

[11] Chen Y, Hong Y J, Makita S, Yasuno Y. Eye-motion-corrected optical coherence tomography angiography using Lissajous scanning. Biomedical Optics Express, 2018, 9(3): 1111–1129

[12] Potsaid B, Gorczynska I, Srinivasan V J, Chen Y, Jiang J, Cable A, Fujimoto J G. Ultrahigh speed spectral/Fourier domain OCT ophthalmic imaging at 70000 to 312500 axial scans per second. Optics Express, 2008, 16(19): 15149–15169

[13] Lezama J, Mukherjee D, McNabb R P, Sapiro G, Kuo A N, Farsiu S. Segmentation guided registration of wide field-of-view retinal optical coherence tomography volumes. Biomedical Optics Express, 2016, 7(12): 4827–4846

[14] Camino A, Zhang M, Dongye C, Pechauer A D, Hwang T S, Bailey S T, Lujan B, Wilson D J, Huang D, Jia Y. Automated registration and enhanced processing of clinical optical coherence tomography angiography. Quantitative Imaging in Medicine and Surgery, 2016, 6(4): 391–401

[15] Baghaie A, Yu Z, D’Souza R M. Involuntary eye motion correction in retinal optical coherence tomography: hardware or software solution? Medical Image Analysis, 2017, 37: 129–145

[16] Camino A, Zhang M, Gao S S, Hwang T S, Sharma U, Wilson D J, Huang D, Jia Y. Evaluation of artifact reduction in optical coherence tomography angiography with real-time tracking and motion correction technology. Biomedical Optics Express, 2016, 7(10): 3905–3915

[17] Lang A, Carass A, Al-Louzi O, Bhargava P, Solomon S D, Calabresi P A, Prince J L. Combined registration and motion correction of longitudinal retinal OCT data. In: Proceedings of SPIE, Volume 9784, Medical Imaging 2016: Image Processing. San Diego: SPIE, 2016, 97840X

[18] Watanabe Y, Takahashi Y, Numazawa H. Graphics processing unit accelerated intensity-based optical coherence tomography angiography using differential frames with real-time motion correction. Journal of Biomedical Optics, 2013, 19(2): 021105

[19] Shemonski N D, Ahn S S, Liu Y Z, South F A, Carney P S, Boppart S A. Three-dimensional motion correction using speckle and phase for in vivo computed optical interferometric tomography. Biomedical Optics Express, 2014, 5(12): 4131–4143

[20] Lee J, Srinivasan V, Radhakrishnan H, Boas D A. Motion correction for phase-resolved dynamic optical coherence tomography imaging of rodent cerebral cortex. Optics Express, 2011, 19(22): 21258–21270

[21] Carrasco-Zevallos O M, Nankivil D, Viehland C, Keller B, Izatt J A. Pupil tracking for real-time motion corrected anterior segment optical coherence tomography. PLoS One, 2016, 11(8): e0162015

[22] Braaf B, Vienola K V, Sheehy C K, Yang Q, Vermeer K A, Tiruveedhula P, Arathorn D W, Roorda A, de Boer J F. Real-time eye motion correction in phase-resolved OCT angiography with tracking SLO. Biomedical Optics Express, 2013, 4(1): 51–65

[23] Montuoro A, Wu J, Waldstein S, Gerendas B, Langs G, Simader C, Schmidt-Erfurth U. Motion artefact correction in retinal optical coherence tomography using local symmetry. In: Proceedings of MICCAI International Conference on Medical Image Computing and Computer-Assisted Intervention. Boston: Springer, 2014, 17, 130–137

[24] Hu Z, Rollins A M. Fourier domain optical coherence tomography with a linear-in-wavenumber spectrometer. Optics Letters, 2007, 32(24): 3525–3527

[25] Gelikonov V M, Gelikonov G V, Shilyagin P A. Linear-wavenumber spectrometer for high-speed spectral-domain optical coherence tomography. Optics and Spectroscopy, 2009, 106(3): 459–465

[26] Shilyagin P A, Ksenofontov S Y, Moiseev A A, Terpelov D A, Matkivsky V A, Kasatkina I V, Mamaev Y A, Gelikonov G V, Gelikonov V M. Equidistant recording of the spectral components in ultra-wideband spectral-domain optical coherence tomography. Radiophysics and Quantum Electronics, 2018, 60(10): 769–778

[27] Terpelov D A, Ksenofontov S Y, Gelikonov G V, Gelikonov V M, Shilyagin P A. A data-acquisition and control system for spectral-domain optical coherence tomography with a speed of 91 912 A-scans/s based on a USB 3.0 interface. Instruments and Experimental Techniques, 2017, 60(6): 868–874

[28] Leitgeb R A, Wojtkowski M. Complex and coherence-noise free Fourier domain optical coherence tomography. In: Drexler W, Fujimoto J G, eds. Optical Coherence Tomography: Technology and applications. 2nd ed. Switzerland: Springer International Publishing, 2015, 195–224

[29] Gelikonov V M, Gelikonov G V, Kasatkina I V, Terpelov D A, Shilyagin P A. Coherent noise compensation in spectral-domain optical coherence tomography. Optics and Spectroscopy, 2009, 106(6): 895–900

[30] Ai J, Wang L V. Synchronous self-elimination of autocorrelation interference in Fourier-domain optical coherence tomography. Optics Letters, 2005, 30(21): 2939–2941

[31] Leitgeb R A, Hitzenberger C K, Fercher A F, Bajraszewski T. Phase-shifting algorithm to achieve high-speed long-depth-range probing by frequency-domain optical coherence tomography. Optics Letters, 2003, 28(22): 2201–2203

[32] Zhang J, Nelson J S, Chen Z. Removal of a mirror image and enhancement of the signal-to-noise ratio in Fourier-domain optical coherence tomography using an electro-optic phase modulator. Optics Letters, 2005, 30(2): 147–149

[33] Matkivsky V A, Moiseev A A, Ksenofontov S Y, Kasatkina I V, Gelikonov G V, Shabanov D V, Shilyagin P A, Gelikonov V M. Medium chromatic dispersion calculation and correction in spectral-domain optical coherence tomography. Frontiers of Optoelectronics, 2017, 10(3): 323–328

[34] Gelikonov G V, Gelikonov V M. Measurement and compensation for the amplitude and phase spectral distortions of an interference signal in optical coherence tomography for the relative optical-spectrum width exceeding 10%. Radiophysics and Quantum Electronics, 2018, 61(2): 135–145

[35] Matveev L A, Zaitsev V Y, Gelikonov G V, Matveyev A L, Moiseev A A, Ksenofontov S Y, Gelikonov V M, Sirotkina M A, Gladkova N D, Demidov V, Vitkin A. Hybrid M-mode-like OCT imaging of three-dimensional microvasculature in vivo using reference-free processing of complex valued B-scans. Optics Letters, 2015, 40(7): 1472–1475

[36] Moiseev A, Ksenofontov S, Sirotkina M, Kiseleva E, Gorozhantseva M, Shakhova N, Matveev L, Zaitsev V, Matveyev A, Zagaynova E, Gelikonov V, Gladkova N, Vitkin A, Gelikonov G. Optical coherence tomography-based angiography device with real-time angiography B-scans visualization and hand-held probe for everyday clinical use. Journal of Biophotonics, 2018, 11(10): e201700292

[37] Huo L, Xi J, Wu Y, Li X. Forward-viewing resonant fiber-optic scanning endoscope of appropriate scanning speed for 3D OCT imaging. Optics Express, 2010, 18(14): 14375–14384

[38] Moon S, Lee S W, Rubinstein M, Wong B J F, Chen Z. Semi-resonant operation of a fiber-cantilever piezotube scanner for stable optical coherence tomography endoscope imaging. Optics Express, 2010, 18(20): 21183–21197

[39] Park H C, Seo Y H, Jeong K H. Lissajous fiber scanning for forward viewing optical endomicroscopy using asymmetric stiffness modulation. Optics Express, 2014, 22(5): 5818–5825

[40] Chen Y, Hong Y J, Makita S, Yasuno Y. Three-dimensional eye motion correction by Lissajous scan optical coherence tomography. Biomedical Optics Express, 2017, 8(3): 1783–1802

[41] Chauhan B C, Stevens K T, Levesque J M, Nuschke A C, Sharpe G P, O’Leary N, Archibald M L, Wang X. Longitudinal in vivo imaging of retinal ganglion cells and retinal thickness changes following optic nerve injury in mice. PLoS One, 2012, 7(6): e40352

[42] Taibbi G, Peterson G C, Syed M F, Vizzeri G. Effect of motion artifacts and scan circle displacements on Cirrus HD-OCT retinal nerve fiber layer thickness measurements. Investigative Ophthalmology & Visual Science, 2014, 55(4): 2251–2258

[43] Bezerra H G, Costa M A, Guagliumi G, Rollins A M, Simon D I. Intracoronary optical coherence tomography: a comprehensive review clinical and research applications. JACC: Cardiovascular Interventions, 2009, 2(11): 1035–1046

[44] Ksenofontov S, Vasilenkova T. Method of optimizing maximum intensity projection technique for rendering scalar three-dimensional data in static mode, in interactive mode and in real time. Patent of Russian Federation RU 2533055, 2014

[45] Ksenofontov S Y. Application of the method of multiple mutual synchronization of parallel computational threads in spectral-domain optical coherent tomography systems. Instruments and Experimental Techniques, 2019, 62(3): 317–323

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag GmbH Germany, part of Springer Nature