Abstract
In ultrahigh resolution (UHR) optical coherence tomography (OCT) group velocity dispersion (GVD) must be corrected for in order to approach the theoretical resolution limit. One approach promises not only compensation, but complete annihilation of even order dispersion effects, and that at all sample depths. This approach has hitherto been demonstrated with an experimentally demanding ‘balanced detection’ configuration based on using two detectors. We demonstrate intensity correlation (IC) OCT using a conventional spectral domain (SD) UHROCT system with a single detector. ICSDOCT configurations exhibit cross term ghost images and a reduced axial range, half of that of conventional SDOCT. We demonstrate that both shortcomings can be removed by applying a generic artefact reduction algorithm and using analytic interferograms. We show the superiority of ICSDOCT compared to conventional SDOCT by showing how ICSDOCT is able to image spatial structures behind a strongly dispersive silicon wafer. Finally, we question the resolution enhancement of \(\sqrt{2}\) that ICSDOCT is often believed to have compared to SDOCT. We show that this is simply the effect of squaring the reflectivity profile as a natural result of processing the product of two intensity spectra instead of a single spectrum.
Introduction
Indepth imaging of human tissue has been one of the greatest achievements of optical technologies. Optical coherence tomography (OCT) was initiated more than 25 years ago, when a crosssectional image of the human retina using a Michelson interferometer was demonstrated^{1}. The ability to display changes of the refractive index by detecting photons balistically backscattered millimetres inside tissue at the micrometre scale has let to a revolution in the field of ophthalmology and is essential for many other medical fields^{2}. Quite recently OCT has even been demonstrated for macroscopic imaging, thereby adding a new perspective in terms of its applications^{3}.
With a Gaussian spectral profile the axial (indepth) resolution limit is intrinsically given by \(\delta z=\frac{2\,\mathrm{ln}\,2}{\pi }\frac{{\lambda }_{c}^{2}}{{\rm{\Delta }}\lambda }\), where, λ_{ c } is the central wavelength and Δl is the fullwidth at half maximum (FWHM) spectral bandwidth of the light source^{2}. To maximize the penetration depth in tissue λ_{ c } is typically chosen to be in the near infrared (NIR) regime to minimize scattering^{2}, but well below the major absorption bands of water peaking at λ = 3 μm^{4}. In order to maintain δz when increasing the wavelength from the visible to the NIR, one is left to maximize Δλ. In doing so, chromatic dispersion in both optical components and sample will degrade the depth resolution. This is due to each wavelength experiencing a different optical path through the system and sample, causing the optical path difference to differ as well. The effect of the different optical path lengths in the two paths (reference arm and sample arm) for different wavelengths is commonly known as the dispersion mismatch^{2,5}.
To counter the dispersion mismatch, dispersion compensation (DC) is done hardware wise by ensuring that the two arms are constructed identically. However, this makes the setup more costly and increases complexity. Instead simple DC with a glass plate, such as BK7, is today used to balance the dispersion^{5}. Alternatively, a large variety of numerical approaches have been introduced, first for time domain (TD) OCT^{6,7,8} and later for spectral domain OCT (SDOCT). In particular, in SDOCT several new DC methods have been demonstrated^{9,10,11,12,13,14} that can achieve singleinterface DC, i.e., sharpening only one interface in the sample. Only a few methods promise multiinterface DC, which is necessary in order to maintain axial resolution throughout the imaging depth of a multilayered dispersive sample^{15,16,17,18,19,20,21,22}. One approach that can compensate only second order dispersion at multiple interfaces, is the fractional Fourier transform combined with numerical segmentation of the sample and a radon transform, posing a heavy computational load, which scales with the number of predefined sample segmentations in depth^{15}. Another simpler approach is to perform a linear interpolation of the depthdependent DC from two depths where the dispersion mismatch is known^{16}. Two alternative approaches, inspired by quantum OCT^{17,18}, are phase conjugate OCT^{19,20} and chirpedpulse interferometry OCT^{21,22}. These approaches can do evenorder dispersion cancellation, but are costly and complex hardwarewise, due to the requirement of sum frequency generation, while providing only low sensitivity. A numerical scheme exploiting a generalized autoconvolution function for depthdependent dispersion cancellation, also developed in the foot steps of quantum OCT, was proposed by Banaszek et al. and termed ‘blind dispersion compensation’^{23}. This method promises protection from GVD using a conventional SDOCT system with a single detector. Hardware implementations of Banaszek’s approach have also been proposed and termed ‘spectral intensity’ or ‘intensityinterferometric’ OCT^{24,25,26,27}, but these again require two detectors and added complexity of the experimental setup. We here consider the numerical technique of Banaszek and term it intensity correlation spectral domain OCT (ICSDOCT). All reports so far on implementing ICSDOCT, both numerical and hardwarewise, share two major drawbacks compared to conventional SDOCT: (1) Halving of the imaging depth, and (2) the appearance of IC artefacts stemming from intensity cross terms. Extending the numerical scheme of Shirai et al.^{28}, we here for the first time demonstrate ultrahigh resolution SDOCT with alldepth multiinterface sample dispersion removal with significant artefact reduction and full imaging depth using a conventional SDOCT setup. We do this by numerically implementing the IC scheme, but on the analytical signal, which we distinguish from the standard IC scheme by denoting it ICA. An artefact reduction scheme similar to what is presented in^{28} is then applied to the ICA signal, and we show that the scheme works equally well using data from a conventional OCT setup. By imaging two different silicon phantoms, we highlight the applicability of ICSDOCT with a conventional SDOCT setup, and show that GVD is intrinsically removed at all depths of the sample with no depth segmentation or conventional DC needed, while maintaining the imaging depth.
Theory
In this section the theory behind ICSDOCT is presented. First, we introduce the basic concept of ICSDOCT in the setting of conventional OCT, and we later apply the analytic signal to explain the image depthmaintaining procedure. Subsequently, the full mathematical framework of ICASDOCT is presented, including the artefact reduction technique. Finally we discuss the axial resolution in ICSDOCT, and show numerical simulations to validate the theoretical predictions, and ICASDOCT is compared to quantum OCT.
Intensity correlation spectraldomain optical coherence tomography  ICSDOCT
In SDOCT, the channelled spectrum (interferogram) is given by
where E_{ R } and E_{ S } are the electric fields returned to the spectrometer from the reference arm and sample arm, respectively, see Fig. 1(a) for a sketch of the setup. The electric field from the reference arm is
and for two scattering centres in the sample arm, the electrical field from the sample can be written as
where I_{0} is the source spectrum, r_{1}, r_{2} are the complex reflection coefficients, l_{ s }, l_{ r } are the sample and reference paths, measured as twice the distance from the beam splitter to the sample surface and reference mirror, respectively. L_{1}, L_{2} are twice the distances from the sample’s surface to each of the scattering centres, and \(\beta (\omega )=\frac{\omega }{c}n(\omega )\) is the wavenumber in the sample, with c being the vacuum speed of light, and n being the depthaveraged refractive index of the sample. In general the depthaveraged refractive index will of course be different for two scattering centres at different depths, but we assume this difference to be negligible, such that n_{1}(ω) ≈ n_{2}(ω) ≡ n(ω).
Assuming real reflection coefficients, Eqs (1–3) and are combined and the normalised interferogram, I_{ n } is obtained through \(\,{I}_{n}=\frac{2I(\omega ){I}_{0}(\omega )}{{I}_{0}(\omega )}\), yielding
where ΔL = L_{2} − L_{1} and Δl = l_{ s } − l_{ r }. To generate the ICSDOCT interferogram, I_{ n } is multiplied by itself, however flipped around the central frequency, ω_{0}, and complex conjugated:
where ω′ = ω − ω_{0}. This intraspectral product between the two optical frequency components ω_{0} + ω′ and ω_{0} − ω′ can be understood as probing the sample and the reference object at two different frequencies and seeking cross correlations between all the four electric fields involved, hence fourthorder field correlations. This classical approach is inspired by quantum OCT directly measuring fourth order correlations, which will be discussed in a later section. It is important to note that as ICSDOCT is a classical analogy to quantum OCT, the first realizations supposedly required two spectrometers to mimic the two photo detectors of quantum OCT, but Shirai showed that it is fundamentally equivalent to using two identical spectra obtained by one spectrometer instead of two different spectrometers in a ‘balanced detection’ configuration^{29}. This means that I_{ n }(ω_{0} − ω′) in eq. (5) can be obtained either experimentally or numerically from I_{ n }(ω_{0} + ω′).
Equation (5) contains multiplication of cosines from eq. (4), which creates oscillations with half the initial period. As a result, all peaks of the Fourier transform of eq. (5) are shifted to twice the optical path difference (OPD) due to the decreased period of the oscillations, as illustrated in Fig. 1(b) with the solid curve being shifted to the dashed curve. The spacing between points of discrete sampling of the spectrum in ωspace, Δω, is fixed by the spectrometer, which fixes the depth range (both positive and negative OPD) to z_{ S } = 2πc/Δω. It is therefore possible that peaks that were well below the Nyquist limit z_{ N } = z_{ S }/2 before the multiplication are above after, reaching up to twice the Nyquist limit, as illustrated in Fig. 1(b). This effectively reduces the available depth range without aliasing (the Nyquist Sampling Theorem) by a factor of two compared to conventional SDOCT^{26}. In addition, the cross terms from multiplication between different cosines cause artefacts, which deteriorate the image quality^{23,26,27,29}. For ICSDOCT to be relevant, the imaging depth must be restored, and the artefacts eliminated.
Restoring the imaging depth using the analytical signal – ICASDOCT
The coloured hatched areas in Fig. 1(b) indicate the ICSDOCT signal and aliased signal are trespassing into one another’s imaging range set by z_{ N } (vertical dashed lines). To eliminate this aliasing problem we here, for the first time to our knowledge, propose to use the complex analytic interferograms in eq. (5), instead of the realvalued interferograms. The complex analytic signal I_{ a } of a real signal I is computed by applying the Hilbert transform (HT),
where ⊗ denotes convolution. The analytic signal is zero for negative OPDs by definition, and thus also for depths between z_{ N } and z_{ S } due to the repetition of the spectrum of discretely sampled signals, as illustrated in Fig. 1(c) (solid line). The components of the ICSDOCT interferogram that are deeper than the Nyquist depth (red part in Fig. 1(b)) are, when using the analytic signal, therefore fully distinguishable, i.e., aliasing is eliminated, as seen in Fig. 1(c) (dashed line). We term this the ICA scheme. As a result of any of the IC and ICA procedures, the density of points is doubled, but by using the ICA scheme, the imaging depth is maintained because all points are utilised and not only half.
Shirai et al. have theoretically investigated the application of ICSDOCT for multiple scattering samples in the special case where dispersion originates from only a dispersive element in the sample arm of the SDOCT system^{27}, i.e., neglecting the dispersion from the sample itself. Here we present an extended derivation that is based on a conventional SDOCT setup and takes the dispersion from the sample into account. We want to derive a theory for multiple scatterers because the IC and ICASDOCT procedures cause artefact to emerge due to the multiplication in eq. (5) creating cross terms. We therefore consider the simplest case with cross terms, which is with two scatterers, without loss of generality.
Using the analytic signal of eq. (4), Taylor expanding \(\beta (\omega )=\sum _{j=0}^{\infty }\frac{{\beta }_{j}{\omega ^{\prime} }^{j}}{j!}={\beta }_{0}+{\beta }_{1}\omega ^{\prime} +\,{\beta }_{NL}^{(even)}+{\beta }_{NL}^{(odd)}\), with \({\beta }_{NL}^{(even)}=\sum _{i=1}^{\infty }\frac{{\beta }_{2i}{\omega }^{^{\prime} 2i}}{2i!}\), and \({\beta }_{NL}^{(even)}=\sum _{i=1}^{\infty }\frac{{\beta }_{2i+1}{\omega }^{^{\prime} 2i+1}}{(2i+1)!}\) containing, respectively, the even and odd nonlinear terms of the dispersion, we find
and from eq. (5)
with * denoting complex conjugates. The four first terms in the second and third lines in eq. (8) are equivalent to the four terms from conventional OCT in eq. (7), but now positioned at twice the OPD and without any GVD from the dominant dispersion term β_{2} and all other even orders of dispersion. Contrary, the odd dispersion terms are not removed, and they are even enhanced by a factor of 2, but the dominating term, 2β_{3}, often has a much weaker effect than β_{2} has in conventional SDOCT^{30}. The remaining five terms are artefacts emerging from the cross terms of the multiplication in eq. (5), and they will be treated in the following section. In order to maintain the correct physical distance, the zaxis must be scaled by a factor \(\frac{1}{2}\) as previously explained, and the point density is thus also increased by a factor of two, maintaining the imaging depth.
Artefact reduction
As discussed above I_{ n }(ω_{0} − ω′) in eq. (5) can be obtained either experimentally or numerically by mirroring I_{ n }(ω_{0} + ω′). The ‘balanced detection’ experimental configuration has been shown to supress some of the artefacts in ICSDOCT^{26,27}, i.e., some of the last 5 terms in eq. (8). It has also been shown that in the numerical configuration artefacts can also be removed, but only one at a time using a window function^{26,29}. Very recently, Shirai showed that a numerical scheme can generically remove all artefacts in the dualspectrometer configuration^{28}. Here, we briefly go through the scheme, improve it slightly by introducing a weight function, both on basis of eq. (8), showing that the algorithm works equally well using a conventional SDOCT setup, leaving the complex and expensive dualspectrometer redundant.
The five terms last artefact terms in eq. (8) all have a ~ cos ω_{0} dependence, either explicitly or implicitly through \({\beta }_{0}=\frac{{\omega }_{0}n({\omega }_{0})}{c}\). To reduce the artefacts, we employ a procedure based on varying the centre frequency. This helps identify the artefacts, as first noted by Banaszek et al.^{23}. Varying the centre frequency of the source is challenging, and therefore a numerical procedure is implemented instead. A flowchart illustrating the process in seen in Fig. 2. The process works by numerically splitting the normalised, analytic spectrum I_{n,a} of length N into M subspectra of length N − M + 1, whose centres are shifted 1 pixel relative to their neighbours’, as illustrated in panel 1 and 2 of Fig. 2. The first of the M spectra comprises the first N − M + 1 pixels of the full spectrum. The next subspectrum starts at pixel 2 of the full spectrum and so on, until subspectrum M, which is the last N − M + 1 pixels of the full spectrum. This procedure varies numerically the centre frequency at the cost of narrowing the spectrum by M − 1 pixels. The ICASDOCT procedure of eq. (8) is then applied to all M subspectra independently, giving M ICASDOCT subspectra, shown in panel 3 of Fig. 2. These spectra correspond to a span of ω_{0}′s with a fixed ω′ axis, and because the artefacts oscillate in ω_{0}, while the real OCT terms do not, the artefacts can be removed by averaging the M ICA subspectra. The ω_{0} span must be sufficiently large to ensure that the oscillations are averaged out effectively. Intuitively, this would require the span to cover at least one period of the oscillation, but that is, in fact not enough. Say the span of ω_{0}’s covers a noninteger number of periods. The fraction of a period in the end of the span will then, when all the values are summed, leave a residual, such that the artefact will still be visible, and because the artefacts oscillate with different periods in ω_{0}, it is not possible to choose the span to cover exactly an integer number of periods for all artefacts. This implies that the amount of periods needed are larger than one. However, applying a weighting function, w(ω_{0}) for this final summation, shown in panel 4 of Fig. 2, greatly reduces the Mvalue required for sufficient artefact reduction. The weighting function weighs each subspectrum, such that the first has a lower weight than the second does, and the central subspectrum has the highest weight. This reduces the influence of the fractional periods in either end of the ω_{0}span, which in turn greatly reduces the residual, i.e., the artefact for a given M. As a result, the M value required to suppress the artefact to a given level is reduced when using the weights. As a rule of thumb when choosing M, we shall require at least 5 full oscillations of every artefact to ensure complete removal of all artefacts. The slowest oscillations are, in most cases, the ones with \(\cos \,({\beta }_{0}{\rm{\Delta }}L)\approx \,\cos ({\omega }_{0}{\rm{\Delta }}L)\) in eq. (8), requiring
where Δz is the smallest OPD between two reflectors. Here we used an Mpoint Hanning window as weights.
The ICASDOCT Ascan is reached by performing a Fourier transform on the averaged ICASDOCT spectrum. However, due to the spectral multiplication, all the reflection coefficients r_{1,2} are also squared, and the resulting depth scan is thus a profile of the squared reflectivity instead of just the reflectivity. To reobtain the OCT reflectivity profile (first order in reflectivity) the square root of the depth scans are evaluated.
The implication of the narrowing of the spectrum by M − 1 pixels depends on the hardware employed. When using a source with a Gaussianlike spectrum, the narrowing will not matter much. The Gaussian shape means that the lost pixels near the edges have little amplitude. However, in this study we use an ultrabroadband supercontinuum source, which ensures the interference pattern covers the entire range of the spectrometer (see Methods). In this case, narrowing the spectrum will consequently also deteriorate the axial resolution with a factor of N/(N − M + 1), meaning that an optimal Mvalue can be determined as a tradeoff between the quality of the artefact reduction and the deterioration of the axial resolution.
Finally, the processing time must be discussed. In addition to conventional SDOCT processing, ICASDOCT requires a Hilbert transform on the full spectrum, M splits of the full spectrum, M × (N − M + 1) floating point multiplications to compute the M ICA subspectra, and as many additions for the averaging. For M ≪ N, the number of operations are proportional to M × N, i.e. the processing time scales linearly with M. Typical Mvalues are in the range 100–200 (see Results), and so, the time to process an ICSDOCT image is typically hundreds of times slower than conventional OCT making it unsuited for realtime applications. However, the processing of each Ascan is independent, and thus a heavily parallelized GPUimplementation could make a realtime imaging available.
OCT axial resolution in ICOCT
In the literature ICOCT, both TD and SD, is generally claimed to have a \(\sqrt{2}\) better axial resolution than conventional OCT^{21,23,25,26,27,29,31,32}. However, we find this to be misleading because it originates from not defining the axial resolution from the same signal, i.e., conventional OCT defines it from the reflection profile, whereas ICOCT defines it from the squared reflection profile. In ICSDOCT for example, the intensity spectrum is after the Michelson interferometer mirrored and combined with itself, whereas in for example chirpedpulse ICTDOCT two oppositely chirped pulses are combined. In other words, ICOCT in general exploits fourthorder correlations in that it combines two intensity spectra, i.e., four complex field spectra, whereas conventional OCT exploits secondorder correlations.
However, if the original signal that is about to be squared in ICOCT, cannot resolve two closely spaced reflectors but shows them as a single peak, then the squared signal will also only show a single peak. Thus, if the resolution was not defined as the FWHM of the Ascan of a single mirror, but as the distance between reflectors the system is able to resolve, then there would be no improvement in resolution with ICOCT.
We would like to note that a “true” resolution improvement of \(\sqrt{2}\) compared to standard (classical) OCT is found in socalled quantum OCT, which by nature requires two detectors and therefore inherently is ICOCT, as demonstrated in^{17,32}. This stems from the spectral entanglement shared between two photons. One photon travels the path of the reference arm and the other the arm of the sample. The two photons are subsequently mixed on a beam splitter (as in conventional OCT), after which a coincidence event is recorded varying the relative time delay (scanning the reference arm length similar to the procedure employed in TDOCT), also known as the HongOuMandel interferometer^{33}. The Ascan so obtained is in fact assimilated to a coincidence curve. For equal paths lengths a dip in the coincidence curve enabled by the unique temporal and spectral correlations between the two photons is observed. Due to the spectral entanglement between the two photons, the FWHM of this curve (assuming a Gaussian shape) is indeed a factor of \(\sqrt{2}\) smaller than the FWHM of the Ascan of a mirror in conventional (classical) OCT, even when comparing the same orders of field reflectivity. This resolution improvement can only be explained by nonclassical correlations between the two photons^{34,35}.
Numerical simulations
We carried out a proof of principle simulation to test our ICASDOCT approach. Figure 3(a1–a4) and (b1–b4) show the ICASDOCT procedure applied to simulated data with one and two reflectors, respectively. For a single reflector conventional OCT is shown in Fig. 3(a1), ICASDOCT without artefact reduction (M = 1) is shown in Fig. 3(a2), and ICASDOCT with M = 50 artefact reduction is shown in Fig. 3(a3). An artefact emerging from the cross term between the single reflector and the DC term is seen at ~400 microns in Fig. 3(a2), which is clearly suppressed by the M = 50 artefact reduction, as seen in Fig. 3(a3) and the zoom in Fig. 3(a4). \(M=50\ge 5\times \frac{2{z}_{S}}{800\,\mu \,m}\) satisfies the criterion that the artefact (term 6 in eq. (8)) oscillate 5 periods in the ω_{0} span. Figure 3(b) shows the result of a simulation of OCT imaging of a 100 micron thick silicon plate. The refractive index of silicon used for the simulation is the experimental data provided in^{36} and then interpolated to fit our spectral pixels by a standard piecewise cubic Hermite interpolating polynomial (PCHIP) routine. Figure 3(b1–b3) shows the same as 3(a1–a3), but for two reflectors. In this case Fig. 3(b2) shows nine peaks including the DC term (1 DC term, 2 reflectors, 1 cross term from conventional OCT, and 5 ICASDOCT artefacts), which corresponds directly to the nine terms in eq. (8).
Figure 3(b1) constitutes the baseline for what is possible to achieve with ICASDOCT in terms of artefact reduction. The two reflectors at 1200 microns and 1550 microns, stemming from the silicon plate (n = 3.5), and the cross term at 350 microns are the three peaks that should be left after the ICASDOCT windowing procedure. Figure 3(b3), which shows the result for ICASDOCT with M = 150 artefact reduction, demonstrates that M = 150 is enough to suppress the ICASDOCT artefacts and recover the 3 peaks from Fig. 3(b1), as expected from the limit in eq. (9), which gives M = 115.
Figure 3(b4) shows a zoomin of the back face of the silicon plate, which for the case of conventional OCT is broadened by GVD, but for ICASDOCT, both with and without artefact reduction is restored to its GVDfree width. The slight increase in FWHM observed for M = 150 is negligible, but for larger M values the broadening becomes more severe. A compromise between the artefact reduction and the resolution deterioration thus has to be established.
Methods
For imaging we used the conventional UHR SDOCT system sketched in Fig. 4, which was recently used in clinical skin studies on healthy patients to show that UHR SDOCT provides superior resolution sufficient to accurately delineate the dermalepidermal junction^{37} and how nanoparticles could improve the contrast in the OCT images^{38}. As optical source, we used a 320 MHz superK Extreme EXR9 OCT system (NKT Photonics A/S) with a longpass filter selecting light in the range 1000–1750 nm. This high repetition rate supercontinuum source is especially suited for SDOCT^{39}. A 50/50 fibre coupler customized for 1300 nm (Goosch and Housego, Netherlands) wavelength, served as the beam splitter and standard achromatic lenses collimated the light in each output arm. In the sample arm, galvanometer scanners were deployed for scanning of the sample through a microscope objective (LSM02, Thorlabs, UK). In the reference arm, a block of glass was placed before the mirror for approximate hardware DC. Interferograms were recorded with a 1300 nm spectrometer C10701470GL2KL (Wasatch, USA) providing a ~400 nm bandwidth and operating at a line rate of 76 kHz.
The spectrometer nonlinearity between wavenumber and pixel number is eliminated by resampling using two reference interferograms collected with a mirror placed at two different axial positions, as in^{11}. This technique can also be used for standard singlereflector DC, which we will compare with ICSDOCT alldepth multireflector DC in the following. With the standard DC, an axial resolution of 3–5 µm (FWHM of Gaussian fit) over the entire 2 mm image range was measured (using a mirror as sample). Laterally we found our system to be able to distinguish features down to 6 μm (USAF target 1951 phantom). For a power of 2.4 mW on the sample, the sensitivity is 89 dB. All interferograms are filtered with a 1300 nm Tukey window in ω′ with bandwidth 300 nm to smoothen the image. All Ascans and Bscans presented are single shot images with no temporal averaging applied.
Data availability
The datasets generated and/or analysed in the current study are available from the corresponding author on reasonable request
Results
To verify the theory and the results of the simulation, we imaged two phantoms. Standard DC, as described in the methods section, was applied only where mentioned explicitly. Phantom 1 is a polished silicon wafer of thickness 255 microns. The GVD of crystalline silicon is estimated to be 1100 ± 100 fs^{2}/mm^{36}, which is sufficient to cause significant broadening of the interface corresponding to the bottom surface of the wafer. Assuming a Gaussian spectrum, the relative broadening factor, p, due to GVD is calculated as:
Here L is the physical axial position relative to the surface of the sample, in this case 255 microns. β_{2} = ∂^{2}β/∂ω^{2} is the GVD parameter, c is the speed of light in vacuum, and Δλ and λ_{ c } are the FWHM and centre wavelength, respectively. From the estimated GVD parameter, we expect a relative broadening of the bottom surface by a factor of ~4.1 ± 0.3.
Cross sectional images, Bscans, of the silicon wafer are shown in Fig. 5, with 5(a) being the image collected without any DC, 5(b) the image with conventional DC of the top interface, and 5(c) the ICASDOCT image. Figure 5(d) shows the profile along the vertical dashed lines between the short horizontal, solid lines. The image with no DC in Fig. 5(a) shows the two surfaces having approximately the same thickness of ~10 microns despite the highly dispersive sample. The top interface is broadened due to the dispersion in the setup, while the bottom interface is broadened by the combined effect of the dispersion in the setup and in the sample. As the setup dispersion and sample dispersion have different signs, the accumulated dispersion for the bottom interface is in magnitude smaller than the setup dispersion, and therefore the bottom interface is thinner in the image than the top interface (but one is not always that lucky!). The image in Fig. 5(b) displays a narrow top interface with a FWHM of 4 microns and a bottom interface with a thickness that has increased by a factor of approximately 4 to 16 microns, as expected from eq. (10). The extra broadening of the bottom surface is due to the setup dispersion having been cancelled, and it highlights the major drawback of conventional DC: Not all depths can simultaneously achieve the theoretical dispersionfree axial resolution. As shown in Fig. 5(c) the ICASDOCT method allows thinning of all interfaces to about 4 microns simultaneously, irrespective of depth, by intrinsic cancellation of all even order GVD. The ICASDOCT image is created with M = 150 subspectra, which allows to obtain significant reduction of the artefacts originating from cross terms between two scatterers seen in Figure 4(a2) and (b2), with no trace of these artefacts even on a logarithmic scale. We note however that the ICASDOCT procedure introduces a weak set of artefacts in the background of each reflector peak, seen as the blur around both surfaces in Fig. 5(c), and as side lobes in Fig. 5(d), as marked by the black arrows. These noise side lobes stem from the cross term between a scatterer and the background noise, which was not included in the theoretical derivation or the numerical simulations. The width of these side lobes, that appear around every reflector peak, decreases with an increasing M number. However, since the side lobes are a direct consequence of the noise in the system, the side lobes can also be reduced by employing a lownoise source. In a simple phantom as imaged here, the side lobes do not obscure the signal, but in a complex biological sample, this not generally the case. The optimal M values does thus also require sufficiently reduced side lobes, where the level deemed sufficient will depend on the sample.
To further evaluate the performance of the ICASDOCT procedure, we created phantom 2 by placing a silicon wafer with a surface structure below phantom 1. The resulting images are seen in Fig. 6 with (a) showing a top view OCT image of the structured surface, and (b) and (c) the Bscan along the blue line in (a) with conventional DC and ICASDOCT respectively. (d) and (e) show zoomins of the structured surface of (b) and (c), respectively. This phantom imitates a sample with several interfaces and small scale features deep inside it. From Fig. 6(b,d) we see how an SDOCT system with conventional DC is not able to clearly visualise the fine details of the structure on the bottom wafer. In contrast, in Fig. 6(c,e), the ICASDOCT image with M = 200 artefact reduction provides such a good alldepth DC that the structure is clearly visible. The blurry parts in Fig. 6(d) are clearly seen in Fig. 6(e) to be elevated slightly relative to the rest of the wafer, a detail not visible in Fig. 6(d). From Fig. 6(c) we can appreciate that ICASDOCT restore the two surfaces of the top wafer and the top of the bottom wafer to their dispersion free widths.
Summary and Conclusions
In summary, we have theoretically and experimentally demonstrated a new ICASDOCT procedure that allows alldepth DC of all even order GVD. We show that we numerically can eliminate the GVD due to the sample, irrespective of the scattering depth, allowing us to maintain the theoretical axial resolution at all depths, and this using a conventional SDOCT setup with only a single spectrometer. Furthermore, our new numerical procedure maintains the axial range, which was not possible before, as well as generically removes all artefacts emerging from the multiplication of two interferograms. Numerical simulations were performed using a single and dual layer sample to investigate the dispersion compensating abilities and artefact reduction of the ICASDOCT procedure. These simulations demonstrated tolerance to GVD from the sample as well as excellent reduction of the artefacts. Two phantoms were imaged experimentally, a single polished silicon wafer, and the same polished silicon wafer with another structured silicon wafer placed underneath it. We demonstrated how a conventional SDOCT system with conventional singlereflector DC, showed a severely broadened bottom surface due to sample dispersion and was not able to clearly image the surface structure of the bottom wafer. In contrast, our experimental results showed how ICASDOCT processing can compensate dispersion indepth and image the bottom small features 260 microns into the phantom and reestablish a 4micron resolution for both top and bottom surfaces.
In the scope of increasingly applied supercontinuum sources, multilayer DC becomes gradually more relevant as the optical bandwidth is increased to improve the axial resolution of OCT systems, and to this end, ICASDOCT is ideal because the dispersion is automatically and intrinsically removed. Only the number of windowed spectra M around the central frequency ω_{0}, needs to be chosen, but this can become a constant for a given sample, such that after initial tuning of M, imaging can proceed as with conventional OCT. We expect that this procedure will be particularly useful for nondestructive testing and metrology, where highly dispersive samples are common.
Following the recipe of ICASDOCT, all present SDOCT systems can operate as ICASDOCT systems. A future necessary step towards advancing the applicability of ICASDOCT is to investigate the noise properties further and sensitivity, which we have not included in this study.
References
 1.
Huang, D. et al. Optical Coherence Tomography. Science 254 (1991).
 2.
Drexler, W. & Fujimoto, J. G. eds Optical Coherence Tomography Technology and Applications. (Springer International Publishing 2015).
 3.
Wang, Z. et al. Cubic meter volume optical coherence tomography. Optica 3, 1496–1503 (2016).
 4.
Curcio, J. A. & Petty, C. C. The near infrared absorption spectrum of liquid water. J. Opt. Soc. Am. 41, 302–304 (1951).
 5.
Hitzenberger, C. K., Baumgartner, A., Drexler, W. & Fercher, A. F. Dispersion effects in partial coherence interferometry: Implications for intraocular ranging. J. Biomed. Opt. 4, 144–151 (1999).
 6.
Fercher, A. F. et al. Numerical dispersion compensation for Partial Coherence Interferometry and Optical Coherence Tomography. Opt. Express 9, 610–615 (2001).
 7.
Fercher, A. F. et al. Dispersion compensation for optical coherence tomography depthscan signals by a numerical technique. Opt. Commun. 204, 67–74 (2002).
 8.
Marks, D. L., Oldenburg, A. L., Reynolds, J. J. & Boppart, S. A. Digital algorithm for dispersion correction in optical coherence tomography for homogeneous and stratified media. Appl. Opt. 42, 204–217 (2003).
 9.
Wojtkowski, M. et al. Ultrahighresolution, highspeed, Fourier Domain Optical Coherence Tomography and methods for dispersion compensation. Opt. Express 12, 2404–2422 (2004).
 10.
Cense, B. et al. Ultrahighresolution highspeed retinal imaging using spectraldomain optical coherence tomography. Opt. Express 12, 2435–2447 (2004).
 11.
Makita, S., Fabritius, T. & Yasuno, Y. Fullrange, highspeed, highresolution 1μm spectraldomain optical coherence tomography using BMscan for volumetric imaging of the human posterior eye. Opt. Express 16, 8406–8420 (2008).
 12.
Hillmann, D. et al. Common approach for compensation of axial motion artifacts in sweptsource OCT and dispersion in Fourierdomain OCT. Opt. Express 20, 6761–6776 (2012).
 13.
Choi, W., Baumann, B., Swanson, E. A. & Fujimoto, J. G. Extracting and compensating dispersion mismatch in ultrahighresolution Fourier domain OCT imaging of the retina. Opt. Express 20, 25357–68 (2012).
 14.
Bradu, A., Maria, M. & Podoleanu, A. G. Demonstration of tolerance to dispersion of master/slave interferometry. Opt. Express 23, 14148 (2015).
 15.
Lippok, N., Coen, S., Nielsen, P. & Vanholsbeeck, F. Dispersion compensation in Fourier domain optical coherence tomography using the fractional Fourier transform. Opt. Express 20, 23398–23413 (2012).
 16.
Pan, L. et al. Depthdependent dispersion compensation for fulldepth OCT image. Opt. Express 25, 10345–10354 (2017).
 17.
Nasr, M. B., Saleh, B. E. A., Sergienko, A. V. & Teich, M. C. Demonstration of DispersionCanceled QuantumOptical Coherence Tomography. Phys. Rev. Lett. 91, 083601 (2003).
 18.
Nasr, M. B., Saleh, B. E. A., Sergienko, A. V. & Teich, M. C. DispersionSensitive Quantum Optical Coherence Tomography. Opt. Express 12, 1353–1362 (2004).
 19.
Erkmen, B. I. & Shapiro, J. H. Phaseconjugate optical coherence tomography. Phys. Rev. A 74, 041601 (2006).
 20.
Le Gouët, J., Venkatraman, D., Wong, F. N. C. & Shapiro, J. H. Experimental realization of phaseconjugate optical coherence tomography. Opt. Lett. 35, 1001–1003 (2010).
 21.
Kaltenbaek, R., Lavoie, J., Biggerstaff, D. & Resch, K. J. Quantuminspired interferometry with chirped laser pulses. Nat. Phys. 4, 864–868 (2008).
 22.
Lavoie, J., Kaltenbaek, R. & Resch, K. J. Quantumoptical coherence tomography with classical light. 17, 3818–3825 (2009).
 23.
Banaszek, K., Radunsky, A. S. & Walmsley, I. A. Blind dispersion compensation for optical coherence tomography. Opt. Commun. 269, 152–155 (2007).
 24.
Resch, K. J., Puvanathasan, P., Lundeen, J. S., Mitchell, M. W. & Bizheva, K. Classical dispersioncancellation interferometry. Opt. Express 15, 8797–8804 (2007).
 25.
Lajunen, H., TorresCompany, V., Lancis, J. & Friberg, A. T. Resolutionenhanced optical coherence tomography based on classical intensity interferometry. J. Opt. Soc. Am. A 26, 1049–1054 (2009).
 26.
Ryczkowski, P., Turunen, J., Friberg, A. T. & Genty, G. Experimental Demonstration of Spectral Intensity Optical Coherence Tomography. Sci. Rep. 6, 22126 (2016).
 27.
Shirai, T. & Friberg, A. T. Intensityinterferometric spectraldomain optical coherence tomography with dispersion cancellation. J. Opt. Soc. Am. A 31, 258–263 (2014).
 28.
Shirai, T. & Friberg, A. T. Practical implementation of spectralintensity dispersioncanceled optical coherence tomography with artifact suppression. J. Opt. 20, 045610 (2018).
 29.
Shirai, T. Improving image quality in intensityinterferometric spectraldomain optical coherence tomography. J. Opt. 18, 075601 (2016).
 30.
Wu, X. & Gao, W. Dispersion analysis in micron resolution spectral domain optical coherence tomography. J. Opt. Soc. Am. B 34, 169–177 (2017).
 31.
Ogawa, K. & Kitano, M. Classical realization of dispersioncanceled, artifactfree, and backgroundfree optical coherence tomography. Opt. Express 24, 8280–8289 (2016).
 32.
Abouraddy, A. F., Nasr, M. B., Saleh, B. E. A., Sergienko, A. V. & Teich, M. C. Quantumoptical coherence tomography with dispersion cancellation. Phys. Rev. A 65, 053817 (2002).
 33.
Hong, C. K., Ou, Z. Y. & Mandel, L. Measurement of subpicosecond time intervals between two photons by interference. Phys. Rev. Lett. 59, 2044–2046 (1987).
 34.
Campos, R. A., Saleh, B. E. A. & Teich, M. C. Fourthorder interference of jointphoton wave packets in lossless optical systems. Phys. Rev. A 42, 4127–4137 (1990).
 35.
Rarity, J. et al. Twophoton interference in a MachZehnder interferometer. Phys. Rev. Lett. 65, 1348–1351 (1990).
 36.
Li, H. H. Refractive index of silicon and germanium and its wavelength and temperature derivatives. J. Phys. Chem. Ref. Data 9, 561–658 (1980).
 37.
Israelsen, N. M. et al. The value of ultrahigh resolution OCT in dermatology  delineating the dermoepidermal junction, capillaries in the dermal papillae and vellus hairs. Biomed. Opt. Express 9, 958–963 (2018).
 38.
Mogensen, M. et al. Two optical coherence tomography systems detect topical gold nanoshells in hair follicles, sweat ducts and measure epidermis. J. Biophotonics, https://doi.org/10.1002/jbio.201700348 (2018).
 39.
Maria, M. et al. Qswitchpumped supercontinuum for ultrahigh resolution optical coherence tomography. Opt. Lett. 42, 4744–4747 (2017).
Acknowledgements
M. Jensen, N.M. Israelsen, O. Bang and A. Podoleanu acknowledge support from Innovation Fund Denmark through the ShapeOCT grant No. 410700011A. Marie Curie EID UBAPHODESA FP7PEOPLE2013ITN 607627 supports M. Maria, T. Feuchter and A. Podoleanu, the NIHR Biomedical Research Centre at Moorfields Eye Hospital NHS Foundation Trust and the UCL Institute of Ophthalmology and the Royal Society Wolfson Research Merit Award also supports A. Podoleanu. O. Bang and M. Maria acknowledges support from the Horizon 2020 grant GALAHAD (project no. 732613).
Author information
Affiliations
Contributions
M.J. conceived and implemented the numerical schemes. N.M.I. conducted the experiments. All authors analysed the results and reviewed the manuscript.
Corresponding author
Ethics declarations
Competing Interests
The authors declare no competing interests.
Additional information
Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Jensen, M., Israelsen, N.M., Maria, M. et al. Alldepth dispersion cancellation in spectral domain optical coherence tomography using numerical intensity correlations. Sci Rep 8, 9170 (2018). https://doi.org/10.1038/s4159801827388z
Received:
Accepted:
Published:
Further reading

Intensity correlation OCT is a classical mimic of quantum OCT providing up to twofold resolution improvement
Scientific Reports (2021)

Artefactremoval algorithms for Fourier domain quantum optical coherence tomography
Scientific Reports (2021)

Realtime highresolution midinfrared optical coherence tomography
Light: Science & Applications (2019)

Recovering distance information in spectral domain interferometry
Scientific Reports (2018)
Comments
By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.