Simple and robust calibration procedure for k-linearization and dispersion compensation in optical coherence tomography

Abstract. In Fourier-domain optical coherence tomography (FD-OCT), proper signal sampling and dispersion compensation are essential steps to achieve optimal axial resolution. These calibration steps can be performed through numerical signal processing, but require calibration information about the system that may require lengthy and complex measurement protocols. We report a highly robust calibration procedure that can simultaneously determine correction vectors for nonlinear wavenumber sampling and dispersion compensation. The proposed method requires only two simple mirror measurements and no prior knowledge about the system’s illumination source or detection scheme. This method applies to both spectral domain and swept-source OCT systems. Furthermore, it may be implemented as a low-cost fail-safe to validate the proper function of calibration hardware such as k-clocks. We demonstrate the method’s simple implementation, effectiveness, and robustness on both types of OCT systems.


Introduction
Optical coherence tomography (OCT) is an interferometric imaging technique, which provides depth-resolved images up to several millimeters deep with a typical resolution in the order of micrometers. 1 In order to attain the theoretical axial resolution over the entire imaging depth, however, several conditions must be met. First, in Fourier-domain OCT (FD-OCT) systems, the interference signal must be sampled linearly in wavenumber (k) prior to performing each fast Fourier transform (FFT) in order to reconstruct the depth reflectivity profile with bandwidth-limited axial resolution. 2,3 Furthermore, dispersion must be matched in both interferometer arms, as dispersion imbalance leads to an unequal, wavelength-dependent phase delay of the light. A dispersion mismatch will induce broadening of the reconstructed OCT peaks, therefore reducing SNR and resolution. 4 Various solutions have been developed to address both these requirements. However, many of these imply expensive and/ or complicated additional hardware, computationally intensive algorithms or time-consuming calibration procedures.
K-linear sampling can be achieved directly during signal acquisition with specialized hardware or recreated in postprocessing with interpolation of the acquired signal. For spectraldomain OCT (SD-OCT), hardware-driven k-linear acquisition can be achieved with k-linear spectrometers. 5 In swept-source OCT (SS-OCT) systems, linear k-sampling relies mostly on external sampling clocks, called k-clocks, coupled with highend acquisition electronics capable of nonuniform sampling frequencies. 6 If the interference signal is not directly sampled linearly in wavenumber, the most common approach is to perform numerical resampling using interpolation. 7 Interpolation requires a remapping function which relates sample index to actual wavenumber or fractional indices for interpolation. This remapping function is typically determined through an independent calibration procedure and may require further hardware. In SD-OCT, for example, the spectrometer can be wavelength-calibrated using narrow spectral lines from reference sources, such as a krypton lamp. 8 Methods to identify the remapping function directly from OCT measurements have also been proposed; 9 however, such methods can be experimentally and computationally tedious. 10,11 A more straightforward and more common method involves linearizing the unwrapped phase of the interference spectrum measured when imaging a mirror. However, this method does not decouple the effects of nonlinear sampling and dispersion mismatch and is only valid in the complete absence of dispersion, which is difficult to achieve in real systems. If applied to systems presenting a dispersion mismatch, simple phase linearization leads to remapping functions valid only for one specific depth position. 12 In a recent publication, Uribe-Patarroyo et al. 13 proposed an optimization algorithm based on mirror measurements at multiple depth locations, which successfully decouples and corrects nonlinearity in wavenumber and dispersion mismatch. While this method provides a "true" remapping function, independent of dispersion and therefore valid across the entire imaging depth, the associated measurement protocol and algorithm deployment present additional experimental complexity. By comparison, Wang and Ding 14 and Makita et al. 15 demonstrated a method which also extracts the true remapping function which relies only on two measurements and very simple algorithmic steps. In this method, a mirror is imaged at two different axial positions, and the unwrapped phase signals are subtracted from one another to remove the component due to dispersion. As the resulting phase signal must be linear with wavenumber, any nonlinearity in the subtracted signal can be attributed to nonlinear k-sampling. The remapping function can then directly be evaluated by linearizing this differential phase signal.
Dispersion compensation is a widely researched subject in the field of OCT. Many different methods to correct dispersion mismatch between the sample and reference arms have been proposed, both hardware-and software-based. Furthermore, it is important to distinguish between methods that only correct system dispersion and those that also correct for sample dispersion. The latter are particularly important when imaging structures through thick dispersive media as is often the case in ophthalmic OCT. 16 Hardware-based dispersion compensation can correct for both system and sample depending on the implementation. Solutions include physically matching dispersion in both arms by adding compensating material, 17 dispersive prisms, 18 gratings, 19 or fiber-stretchers 20 in one of the interferometer arms. However, all these methods typically only compensate up to second-order dispersion. Higher orders can be minimized by adding comparable components to both interferometer arms or through more complex setups that include several of the components mentioned above. 21 Such systems have been shown to correct dispersion up to the third order. 22,23 On top of the limitations regarding correction order, physical dispersion compensation shares the same pitfalls as hardwarebased k-linearization: increased complexity and cost. Various software methods have also been developed to address the problem of dispersion compensation. Such techniques often allow to correct dispersion mismatches to higher orders than what is possible through physical compensation and can, in some cases, simultaneously correct system and sample dispersion. Among these, some rely on an iterative adjustment of a phase correction signal to optimize image sharpness. 4,24,25 This approach is particularly useful for correcting sampleinduced dispersion which cannot be corrected prior to measurements. However, such optimization techniques are computationally intensive and, for the most part, impractical for real-time processing. Another method, presented by Uribe-Patarroyo et al., 13 relies on iterative optimization to extract the phase component due to system dispersion, which can then be corrected in future measurements at low computational cost. Recently, Singh et al. 26 reported a dispersion compensation scheme which uses differential phase measurements, much like those in Wang's k-linearization method, to exactly extract the phase component generated by system dispersion. Similarly to Uribe-Patarroyo's process, this phase component can then be used to correct dispersion in all future measurements. Singh's method requires only two mirror measurements, which must be positioned at precisely symmetrical locations about the zero-delay plane. When the phase signals from these measurements are subtracted from one another, the (linear) depth-dependent phase is removed, leaving only the dispersion-related terms.
As such, numerous methods exist for k-linearization and dispersion compensation for FD-OCT systems. However, these methods imply varying levels of technical, algorithmic and experimental complexity as well as varying accuracy and robustness. Ideally, implementation complexity should be minimized without sacrificing robustness and effectiveness. To this effect, we herein present and validate a simple and complete calibration procedure, which simultaneously extracts calibration information for k-linearization and dispersion compensation. Thorough descriptions of all experimental and algorithmic steps are provided to facilitate the method's implementation. The proposed method requires no additional hardware or a priori information, is highly robust, and relies only on two mirror measurements, obtained on either side of the zero-delay plane. Furthermore, we introduce a numerical phase shifting step allowing for the mirror measurements to be performed at arbitrary positions, as long as they are on either side of the zero-delay plane. We demonstrate that the method is equivalently applicable to SD-and SS-OCT systems, both with and without k-clock. Finally, we demonstrate that this procedure may still be applied in the case of improper configuration of the k-clock.

Theory
Our method, expanded in Sec. 7 and summarized in Fig. 1, is based on the manipulation of the phase of the complex spectral interference signal obtained from mirror measurements, using a standard FD-OCT system. We assume that there is a certain degree of dispersion mismatch due to the presence of a segment of dispersive material in one of the interferometer arms. All steps of the proposed method rely on the expression of the spectral phase signal (Δϕ), obtained experimentally, as the sum of a k-linear, depth-dependent component (ϕ lin ) and a k-nonlinear, depth-independent component (ϕ disp ): E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 3 2 6 ; 4 6 5 Δϕ In this equation, δz is the depth, k 0 represents the central wavenumber, n 0 g represents the first derivative of the group index with respect to k of the segment causing the dispersion mismatch between the two interferometer arms, and D is its thickness. A detailed mathematical derivation of Eq. (1) is given in Sec. 7. Experimentally, the phase signal is determined by computing the phase of the analytical complex interference signal, obtained from the real measured signal using the Hilbert transform (H), as described as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 3 2 6 ; 2 9 6Ĩ ðkÞ ¼ IðkÞ þ iH½IðkÞ; (2) E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 3 ; 3 2 6 ; 2 5 3 Δϕ ¼ tan −1 Im½ĨðkÞ Re½ĨðkÞ : From Eq. (1) and the properties of the Fourier transform, it is apparent that the linear term will determine the position of the peak in the processed A-line, while nonlinear, dispersioninduced terms will be responsible for peak broadening and distortion. 4 It is important to note that all the following steps assume that the phase signal has been unwrapped. From Eq. (1), we define the two calibration measurements, where Δϕ þ refers to the measurement with the mirror physically located on one (positive) side of the zero-delay plane, at δz þ , and Δϕ − refers to the measurement with the mirror located on the opposite (negative) side, at δz − : E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 4 ; 3 2 6 ; 1 0 7 The sign reversal of the dispersion component (ϕ disp ) observed in Eq. (5) is due to the experimental retrieval of the (real) interferometric signal, discussed in Sec. 7.3. We will refer to a set of two such measurements as a calibration pair throughout this paper.

k-Linearization
Using Eqs. (4) and (5), and following Wang et al., we can cancel out dispersion-related terms and isolate a linear term 10 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 6 ; 6 3 ; 3 0 4 Δϕ lin ¼ Equation (6) describes the phase term of a dispersion-free signal corresponding to a mirror located at the average absolute position of the two individual peaks. Any nonlinearity in the unwrapped phase can, therefore, be attributed to improper sampling. Fractional indices for resampling can be obtained by linearizing Δϕ lin . Linearization of the phase signal implies generating a linear vector of the same length as the vector Δϕ lin spanning from Δϕ lin;0 to Δϕ lin;N−1 where the subscripts refer to first and last samples, respectively. For each element of this newly generated, linear vector, we can identify a fractional index through interpolation corresponding to its position in the original phase vector. Thereafter, these fractional indices can be used to resample all subsequent measurements to obtain k-linear sampling. This process should be carried out on both measurements of the calibration pair before dispersion compensation can be performed.

Dispersion Compensation with Numerical Phase Shifting
Once k-linearized, the interference spectra are processed to obtain A-lines. The values of δz AE are then evaluated by finding the axial position of the peak corresponding to the mirror. The peak axial positions can be evaluated in several ways such as applying a fit (e.g., Gaussian) or directly identifying the index corresponding to the position of the peak maximum. In both cases, the interference spectra should be extensively zeropadded prior to FFT to increase the number of data points available for the measurement. Measurements based on identifying the position of the maximum are particularly sensitive to this as insufficient padding may lead to errors due to the discrete nature of the signal. If proper k-linearization has been performed, the shape of the axial PSF should be completely independent of depth. Therefore, if measured in the same way, the peak positions relative to the axial profile will be identical for both mirrors, regardless of the goodness of fit. The numerical shift computed from the peak positions [Eqs. (7)-(9) or (11)] will then always be valid, independently of how the peak positions were identified. As such, Gaussian fitting will be adequate even if the PSF deviates from a Gaussian profile due to a non-Gaussian illumination spectrum or due to the presence of significant high-order dispersion. Once the peak positions have been identified, we can numerically apply a linear phase correction, equivalent to translating both mirrors in air, which brings the two measurements Δϕ Ã þ and Δϕ Ã − ) back to perfectly symmetrical positions: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 7 ; 3 2 6 ; 1 0 7 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 9 ; 6 3 ; It is interesting to note that the applied phase shift and equivalent translation are arbitrary so long as the final positions of the peaks are symmetrical. The dispersion-induced phase term can then be isolated by finding the difference as proposed by Singh et al.: 26 ; t e m p : i n t r a l i n k -; e 0 1 0 ; 6 3 ; 5 9 2 ϕ disp ¼ This phase signal can then be fitted relative to sample index using N'th-order polynomial fitting. Thereafter, dispersion can be corrected by multiplying the complex signal by expð∓iϕ fit disp Þ for peaks on the positive side and the negative side of zero-delay, respectively. It is important to notice that it is not necessary to have any knowledge of the real values of k or δz AE . Indeed, if the interference signal is defined in index space (j ¼ ½0;1; 2; N − 1), then the processed A-lines will be located in a space defined from −π to π corresponding to the optical path length difference. The peak positions (θ peakAE ) determined previously can then directly be used to compute the phase shift vector and shifted phase signals, expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 1 ; 6 3 ; 4 1 6 Δϕ shiftAE ðjÞ ¼ j −jθ peakAE j þ jθ peakþ j þ jθ peak− j 2 ; E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 2 ; 6 3 ; 3 5 9 Δϕ Ã AE ðjÞ ¼ Δϕ AE ðjÞ þ Δϕ shiftAE ðjÞ: The shifted phase signals determined in Eq. (12) can then be used to compute the dispersion phase component described in Eq. (10). A MATLAB implementation of the full calibration process is available for download on Code Ocean: https://codeocean. com/capsule/0521437/tree/v1.

System Description
The calibration procedure was tested on two separate FD-OCT systems: one SD-OCT and one SS-OCT, shown in Fig. 2. The SD-OCT system (OCP840SR-NP, Thorlabs) utilizes a superluminescent diode source (SLD-371, Superlum, Ireland) and a 1024-pixel spectrometer with a spectral resolution of 0.06 nm. A-lines were acquired at a 5-kHz rate. The SLD spectrum was measured using a USB4000 spectrometer (Ocean Optics) and found to be centered at 850 nm with a full-width half maximum (FWHM) spectral bandwidth of 26 nm. The system also includes a 50/50 broadband fiber coupler as well as collimation and focusing optics. Through physical length measurements of the interferometer arms, we found a 1 AE 0.1 cm mismatch in fiber length between the sample and reference arms. The SS-OCT is a custom system using a 50-kHz commercial swept-source (Swept laser 1310, Axsun) with an integrated k-clock. The source spectrum was measured using a NIRQuest512 spectrometer (Ocean Optics) and found to be centered at 1313 nm with an FWHM bandwidth of 95 nm. The interference signal was detected using a balanced photodetector (PDB450C, Thorlabs) and filtered with an 80-MHz low-pass filter (VLF-80+, Mini-Circuits) prior to digitization. Signal sampling was carried out with a high-speed acquisition card (ATS9350, AlazarTech, Canada) capable of both uniform and nonuniform sampling speeds. For measurements with the k-clock, 1088 samples were recorded per acquisition at a mean acquisition rate of 130 MS∕s. Measurements without k-clock contained 3855 samples acquired at 500 MS∕s. The number of samples in measurements without the k-clock was selected to match the spectral range obtained with the k-clock based on recognizable features in the spectra. Data acquisition was controlled using custom LabView software (National Instruments) and synchronized with the swept-source laser using a custom fiber Bragg grating (FBG) sweep trigger. For the purpose of this experiment, a dispersion mismatch was induced by adding ∼1 cm of distilled water in a petri dish to the sample arm of the SS-OCT system.

Measurement Protocol
Single A-line measurements were carried out on a mirror for axial positions (δz) ranging over AE1 mm around the zerodelay plane, at 100-μm intervals. A manual micrometric stage was used to translate both the focusing lens and the sample mirror in order to avoid signal variation due to confocal gating. Raw interference signals (i.e., without any processing) were recorded from both systems using custom LabView software. All signals used for calibration were computed from the average of 64 consecutive A-line measurements. Background removal was achieved by subtracting the background spectrum from the mirror measurements as well as numerical high-pass filtering. The background was measured by successively blocking the sample arm, the reference arm, and both arms. The background signal subtracted from the other measurements was computed as The calibration described in Sec. 2 and in Fig. 1 was carried out for all possible combinations of negative and positive mirror locations. In the calibration, signals were zero-padded to 10 times their original length prior to FFT, and a fifth-order polynomial fit was used to compute the dispersion correction vector (step 5 in Fig. 1). For each combination, the extracted calibration vectors were applied to all individual measurements, and their peak position and FWHM were measured. Finally, for each mirror location, the peak position and FWHM were averaged over all calibration pairs. The negative and positive peaks located at AE100 μm were not included during the calibration as these peaks had reduced SNR due to the high-pass filtering applied previously. This combined with their position close to the edge of the A-line led to instabilities when applying fits to the data and inaccurate calibration. For the SS-OCT system, data were acquired both with and without the k-clock. Measurements were also carried out on the SS-OCT system without the added dispersion (i.e., without the 1 cm water in the sample arm) and with varying k-clock delays. K-clock delay refers to the delay between the optical trigger in the swept laser source and the electronic trigger signal at the acquisition card. Improper adjustment of this parameter may lead to incorrect sampling as discussed further below. For both systems, the spatial pixel increment (μm∕pixel) was determined from the two outermost mirror position (known separation of 2 mm).

Practical Considerations
For optimal implementation of the proposed method, it is essential to consider certain practical aspects. First, proper background removal is crucial to the method as the presence of background will strongly influence the phase signal quality. The presence of strong DC or low-frequency noise components will prevent proper phase retrieval with the Hilbert transform, hence the necessity for proper background subtraction and high-pass filtering. Furthermore, the signal of the calibration measurements should be processed to reduce random noise or contributions from system interference. Random noise can easily be attenuated by averaging multiple interference spectra. However, in systems with low phase stability, care should be taken not to wash out fringes by averaging over too many A-lines. Nonrandom system noise due to interference between various components in the optical path will appear as additional peaks in the processed A-line. The mirror peaks used for the proposed calibration procedure must be distinct from these background peaks. This background signal can be eliminated by applying spatial bandpass filtering (i.e., filtering out additional peaks in the signal after FFT and then transforming back). Finally, "null" data points, i.e., data points without any signal, should be omitted from this analysis as they will produce unstable phase values. Such points are easily identifiable by their nonsmooth appearance when plotting the unwrapped phase of the acquired signal. They are usually located at the beginning or end of data acquisitions and correspond to the extremities of the broadband illumination spectrum (dark pixels for SD-OCT systems and time-points outside of the duty cycle of the wavelength-swept laser for SS-OCTs).

Results
The results of the calibration procedure on one-dimensional measurements for both the SD and SS-OCT systems are shown in Fig. 3. The figure plots the axial resolution as a function of depth (physical mirror position) for the SD-OCT system [ Fig. 3(a)], for the SS-OCT system with k-clock [ Fig. 3

(b)]
and without k-clock [ Fig. 3(c)] acquisition. Figure 3(d) reproduces the curves of Fig. 3(c) using the same scale as 3B for simplified comparison. The error bars in Fig. 3 are obtained by calculating the standard deviation (STD) of the peak FWHM computed with all possible calibration pairs (n ¼ 81). For each case, the resolution is plotted as a function of depth for raw (blue dashed curve), k-linearized (red dashed curve), and fully corrected data (k-linearized and dispersion compensatedyellow dashed curve). The solid green line represents the optimal attainable value, often overlapping with the yellow dashed curve. It is important to note that this method corrects peak broadening due to improper sampling and dispersion mismatch. Other effects, including the spectral shape of the illumination, spectral transmission of the optical system, detector spectral efficiency, and alignment of polarization states, may lead to a poorer axial resolution compared to the expected value based only on the source spectral bandwidth. In order to properly assess the method's performance, the recovered axial resolution was compared to the bandwidth-limited, optimal value obtained by directly linearizing the phase of each peak (mean value over all acquired peaks). All numerical values are summarized in Table 1, by calculating the mean and standard deviation of the resolution over all measured axial positions. The relative variation row refers to the STD divided by the mean expressed in percent. The results of the calibration procedure applied to the SS-OCT system without any artificially added dispersion and with varying k-clock delays are shown in Fig. 4 and in Table 2. Figure 4 plots the axial resolution as a function of depth for three different k-clock delays relative to the optimal value (calculated from system fiber length): −32.7 ns [minimum, Fig. 4(a)], −16.8 ns [ Fig. 4(b)], and 3.6 ns [maximum, Fig. 4(c)]. The color scheme and the calculation of error bars and theoretical resolution limit are all identical to that for Fig. 3 and Table 1. The values in Table 2 Figure 3 shows the recovery of k-linear sampling as well as the successful compensation of any dispersion mismatch present in the systems. In all instances, raw data curves are strongly dependent on depth, while k-linearized curves are flat, demonstrating that indeed, k-linearization of the OCT signal removes the depth-dependency of the measured axial resolution, as predicted by the theoretical framework. The subsequent dispersion compensation step improves the axial resolution down to the theoretical value, allowing for the recovery of the optimal axial resolution over the entire imaging range. Table 1 shows Table 1 Statistics of the calibration procedure averaged over all axial positions. Dispersion in the SD-OCT system is due to fiber length mismatch; dispersion is artificially added in the SS-OCT measurements.

SD-OCT SS-OCT (w/ k -clock) SS-OCT (w/o k-clock)
Raw  Journal of Biomedical Optics 056001-6 May 2019 • Vol. 24 (5) that, after application of the calibration, the theoretical axial resolution is recovered to within less than 1% over the entire imaging range for all systems and configurations. The method's insensitivity to the position of the calibration pair is demonstrated by the very narrow error bars in Fig. 3 as well by the low standard deviations reported in Table 1. As such, these results indicate that the optimal axial resolution is recovered, no matter which combination of mirror positions was used. Table 1 show small differences in the recovered axial resolution for the SS-OCT system depending on whether the swept-source k-clock is used or not (13.8 and 14.2 μm, respectively). As both these measurement sets were acquired concurrently and on the same system, it would be expected that the recovered axial resolutions be identical. We propose three potential sources that may contribute to this discrepancy. First, there may have been minute differences in the spectral bandwidth between the two sets of  acquisition parameters. While in both cases the acquisitions were initiated by the spectral sweep trigger, the end of the acquisitions may have varied slightly in time and therefore in their spectral positions. Precautions were taken to minimize this issue; however, it may still account for a portion of the 400-nm difference in axial resolution. A second contributing factor might be the slightly lower SNR observed in the SS-OCT data acquired without the use of the k-clock. Indeed, the average sampling speed without the k-clock was roughly four times higher than when the k-clock was active, which would induce a higher noise equivalent power. Finally, it is possible that the discrepancy is associated with fitting instabilities due to the low number of data points contained within the reconstructed peaks in both cases. As such, the difference might be generated by a slight underestimation of the k-clocked peak width or an overestimation of the non-kclocked one. However, it is important to note that whatever the cause of this effect, its magnitude is sufficiently small to be considered negligible for all practical intents and purposes. Finally, it is particularly interesting to note that this method requires no a priori knowledge regarding the optical system, the acquisition hardware and parameters, or the illumination source. Indeed, it relies strictly on the information acquired during the two measurements. This allows for a very simple implementation, valid across all OCT platforms. Two aforementioned methods to obtain the k-linearization vector 14 and to determine the dispersion correction 26 also only require two mirror measurements. Our method maintains that advantage while generalizing to allow for mirror measurements on opposite sides of zero delay and removing the requirement of exactly symmetrical mirror positions with respect to zero-delay by numerical phase shifting.

k-Clock Sampling Issues
Another crucial point highlighted by this work is the potentially imperfect k-linearization of k-clocks. As observed in Figs. 3(b), 4(a), and 4(b), the axial resolution of the raw signal acquired with a k-clock still varies with depth, which indicates nonlinear k-sampling. This was caused by an improper clock-delay in the k-clock output. As is the case in many wavelength-swept lasers used for OCT, the k-clock signal in our setup is obtained from the zero-crossings of the interference signal of a fiber Mach-Zehnder interferometer (MZI) built into the light source. These zero-crossings are equally spaced in k-space thus producing a k-linear clock signal. However, if the total travel time of the light in the OCT system is not equivalent to the travel time in the MZI plus the electronic delay of the clock pulse, there is an offset between the clock pulse and the intended sampling point, which reintroduces nonlinearities in sampling. Such effects may, therefore, appear when the total fiber length of the OCT system is altered. While some wavelength-swept lasers offer adjustable clock-delay, the range of the delay is often limited which implies inevitable axial resolution degradation when operating outside of this range. Such situations may arise in clinical systems operating with long endoscopic probes or, inversely, in systems using photonic integrated circuits (PICs) with very short total optical paths. OCT systems based on PICs are doubly advantaged by the proposed method as they also present significant peak broadening due to dispersion mismatch. [27][28][29] The effect of an incorrect k-clock delay can be quite significant as shown in Fig. 3(b), where the axial resolution varies by 5 to 6 μm over a 1-mm imaging range. Figure 4 shows that the calibration procedure successfully recovers optimal resolution over the entire imaging range, regardless of the k-clock delay. Furthermore, our method remains stable and robust even when the signal is already linearly sampled in k-space as is the case in Fig. 4(c). In this scenario, the fractional indices for resampling are almost identical to the initial indices and do not alter the measured signal. The root-mean-square difference between the original indices and the resampling indices was found to be 0.01%, and the mean axial resolution over the imaging range was identical. As such, this method can be implemented as a low-cost fail-safe to ensure correct sampling even when a k-clock is used. Alternatively, it could be used as a validation method to ensure that the proper k-clock delay is being used. Indeed, improper delay settings would be noticeable through differences between the original indices and the fractional indices for resampling.

Calibration Using Two Measurements on the Same Side of Zero-Delay
It is possible to deploy an approximation of the procedure outlined in this paper using two mirror measurements located on the same side of zero-delay. In such a case, the phase from both measurements would be given by Eq. (4) or Eq. (5). The linear component can then be extracted by subtracting one from the other. The fractional indices can then be computed from this differential signal. 14 Once k-linearized, the dispersion component of a single measurement can be determined performing a polynomial fit on the full unwrapped phase (Δϕ AE ) and keeping only the nonlinear terms. 30 However, it was demonstrated that this method for dispersion compensation was less effective than the differential approach, particularly when higher order dispersion terms become significant. 26 This is particularly relevant for ultrabroadband or visible-light OCT systems. The use of a single measurement for dispersion compensation does not account for the presence of residual secondary peaks due to system noise, whereas in the differential approach, such additional peaks would be canceled out. As such, while this approximation remains viable, it did lead to a small degradation in axial resolution and a wider distribution for our data. It can be expected that this effect will become more significant for systems with more or higher-order dispersion.

Conclusion
In conclusion, we present a method which corrects nonlinear k-sampling and system-induced dispersion mismatch, using only two mirror measurements. Two calibration vectors are extracted which allow numerical resampling for k-linearization as well as a phase correction to account for dispersion compensation. We demonstrate the procedure's simplicity, its robustness to the position of the mirror peaks and its applicability to both SD-and SS-OCT systems. We anticipate that this method may be extended to calibrate systems which also contain sampleinduced dispersion mismatches such as in ophthalmologic imaging through the dispersive medium of the eye. In such a scenario, a strong sample reflection could be used instead of a mirror to perform the calibration.
7 Appendix A

Propagation Constant
In this section, we provide the theoretical background of the proposed calibration method. Our method relies on the manipulation of the spectral interference signal obtained Journal of Biomedical Optics 056001-8 May 2019 • Vol. 24 (5) during OCT measurements of a single reflector, which can be expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 3 ; 6 3 ; 7 3 0 I int ðkÞ ¼ I DC ðkÞ þ 2 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi I s ðkÞI r ðkÞ p Refexp½iΔϕðkÞg; (13) where I int ðkÞ is the intensity measured at the detector, I DC ðkÞ is the sum of I r ðkÞ and I s ðkÞ, the intensities in the reference and sample arms respectively, and ΔϕðkÞ is the phase difference between the fields from the reference and sample arms. The phase in each arm is given as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 4 ; 6 3 ; 6 4 2 ϕ r ¼ 2 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 5 ; 6 3 ; 5 8 8 ϕ s ¼ 2 In Eqs. (14) and (15), the variables β p;q ðkÞ denote the materialspecific propagation constants and l p;q is the physical length of the corresponding medium. The factor 2 accounts for the light traversing this distance twice. In a dispersive material, the propagation constant depends on wavenumber such that E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 6 ; 6 3 ; 5 1 2 βðkÞ ¼ knðkÞ; (16) where nðkÞ is the wavenumber-dependent complex index of refraction specific to each material. The imaginary part of βðkÞ determines the material absorption coefficient. Throughout this report, we consider only the real part of nðkÞ and therefore βðkÞ. The propagation constant can be expanded into a Taylor series around an arbitrary wavenumber, k 0 (typically the central wavenumber), as follows: where the ′ prime denotes differentiation with respect to wavenumber. Using Eq. (16) and applying the rules for chain derivation, we obtain E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 8 ; 6 3 ; 3 2 1 βðkÞ E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 9 ; 6 3 ; 2 4 2 βðkÞ To simplify the notation, we introduce the group refractive index E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 0 ; 6 3 ; 1 8 2 n g ðkÞ ¼ nðkÞ þ kn 0 ðkÞ; which we substitute into Eq. (19) to obtain E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 1 ; 6 3 ; 1 4 1 βðkÞ ¼ k 0 nðk 0 Þ þ n g ðk 0 Þðk − k 0 Þ þ 1 2 n 0 g ðk 0 Þðk − k 0 Þ 2 þ : : : ; The propagation constant, Eq. (21), can then be summarized as the infinite series E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 2 ; 3 2 6 ; 7 5 2 The first term is a constant phase term that may be neglected in subsequent analysis. First-order dispersion (m ¼ 1) accounts for group delay of the spatial interferometric signal, whereas the higher order terms lead to broadening and chirping of the OCT signal.

Interferometer with Dispersion Mismatch
First, consider a perfectly dispersion balanced OCT interferometer. We introduce a segment of glass of thickness D at z ¼ 0 and determine the phase of the interference spectrum Eq. (13), corresponding to the back surface of the glass segment. From Eqs. (14) and (15), ϕ r ¼ 2kl 0 and ϕ s ¼ 2kl 0 þ 2β glass ðkÞD.
Following the steps leading up to Eq. (22), we find that Fourier transform of Eq. (13) yields the spatial OCT signal, consisting of a peaked signal determined by the source spectrum, shifted to depth position z ¼ n g ðk 0 ÞD and possibly widened, chirped or otherwise distorted by the higher orders m ≥ 2 of the propagation constant β glass ðkÞ.
To resemble a realistic OCT system with a mismatch in fiber lengths between both arms of the interferometer, we shift the reference mirror over a distance z ¼ n g ðk 0 ÞD so that the position of zero group delay difference coincides with the back surface of the glass. We then add a further offset, δz, to the sample arm such that E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 3 ; 3 2 6 ; 4 4 3 ΔϕðkÞ ¼ 2kδz þ 2D½β glass ðkÞ − kn g ðk 0 Þ: Rearranging and expanding the propagation constant following Eq. (21) gives E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 6 ; 3 2 6 ; 3 4 6 where the last, constant term does not have any k-dependence.
Changing the k-dependence of the first two terms to ðk − k 0 Þ dependence, so that the second term in the equation and the first term within the brackets cancel, will add constant terms E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 7 ; 3 2 6 ; 2 2 8 The constant terms will have no impact on the processed A-line as they will be factored out of the Fourier transform and removed when we compute the absolute value. As such, we can simply ignore them in subsequent steps. We therefore obtain the simplified expression for the spectral phase of a single reflector: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 8 ; 3 2 6 ; 1 1 1 ΔϕðkÞ ¼ 2ðk − k 0 Þδz þ ðk − k 0 Þ 2 n 0 g ðk 0 ÞD þ : : : Journal of Biomedical Optics 056001-9 May 2019 • Vol. 24(5)

Experimental Retrieval of the Interferometric Phase
Despite the complex notation of Eq. (13), the detected signal is a real-valued modulated spectrum. In a nondispersed interferometer, signals from exact opposite locations with respect to zero-delay cannot be distinguished (the so-called complex ambiguity) unless special measures are taken: the Fourier transform on the spectral interferogram will yield an A-line with mirrorpeaks at −δz and þδz. If the mirror is located at þδz, the spectral interferogram following Eqs. (13) and (28) is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 9 ; 6 3 ; 6 3 5 I int ðkÞ þ ∝ cos½ϕ lin ðδzÞ þ ϕ disp : If the mirror is located at −δz, the spectral interferogram following Eqs. (13) and (28) is We therefore adopt the following notation of the spectral phase corresponding to peaks at positive and negative delays, respectively: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 3 1 ; 6 3 ; 4 9 9 Δϕ þ ¼ jϕ lin ðδz þ Þ þ ϕ disp j ¼ ϕ lin ðjδz þ jÞ þ ϕ disp ; E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 3 2 ; 6 3 ; 4 5 6 Δϕ − ¼ jϕ lin ðδz − Þ þ ϕ disp j ¼ ϕ lin ðjδz − jÞ − ϕ disp :

Disclosures
The authors have no relevant financial interests in this article and no potential conflicts of interest to disclose.