1. Introduction

Optical coherence tomography (OCT) is a promising new imaging modality for clinical visualization of the middle ear [1]. Middle ear OCT (ME-OCT) enables cross-sectional (2D) and volumetric (3D) visualization of the tympanic membrane (TM) and middle ear space non-invasively and without exposing patients to ionizing radiation. In addition to structural imaging, phase-sensitive OCT Doppler vibrometry (OCT-DV) techniques can measure vibrations of middle ear structures in response to sound with sub-nanometer sensitivity [2,3], providing new diagnostic information highly relevant for distinguishing pathologies of the middle ear that result in conductive hearing loss [4]. ME-OCT has been used successfully to identify stapes fixation in a patient population [4], to assess the state of a post-operative ear [5], to detect distortion product otoacoustic emissions as vibrations along the ossicular chain [3], to image cholesteatoma [6], to detect cochlear implant retraction [7] and to image dimeric tympanic membrane segments [8].

Apart from its use in measuring TM thickness [3,8–10], ME-OCT imaging has so far only been used as a qualitative tool to assess the presence and relative location of middle ear structures without reference to or validation against an absolute coordinate space. For many metrology-based applications such as prosthesis sizing, monitoring anatomical changes in healing or in progressive disease, assessing congenital malformations of the middle ear, measurement of TM perforations and quantifying the degree of prosthesis migration, an ability to produce images with high geometric fidelity to the patient anatomy would enhance ME-OCT’s diagnostic value.

However, obtaining geometrically accurate OCT images of the middle ear is challenging. The middle ear lies at the end of a $3\textrm{cm}$ long ear canal that is often narrow and curved. Because of this anatomical limitation, middle ear OCT lends itself well to an entocentric imaging geometry where a wide angular field of view (FOV) is imaged from a small aperture [2,3,10]. Such an imaging geometry produces images exhibiting fan-beam distortion [11] when OCT data is displayed without correcting for the fact that the data was collected along non-parallel scan lines. The geometric accuracy of OCT images may also be affected by distortions introduced from imperfect imaging optics or a non-linear relationship between scanning mirror drive signals and scanning mirror angle. If geometrically accurate images of anatomy are to be displayed to clinicians, all these distortions must be corrected so that image data is displayed on an equally spaced Cartesian grid of points and the images represent a scaled version of the true anatomy.

In point-of-care applications, it is desirable to have imaging data collected, processed, and rendered in real time in both 2D and 3D, including any transformations or corrections applied to mitigate distortion effects. Continuous, live 3D middle ear imaging, sometimes referred to as 4D imaging (3D over time) [12], could enable intuitive clinical exploration of the middle ear space and the capturing of diagnostically relevant dynamic processes such as middle ear muscle contraction, middle ear conformation changes in response to pressurization and the insertion of surgical instruments into the middle ear space. 4D visualization necessitates a capability to process and render OCT data at real-time rates as fast as it can be produced from the OCT imaging hardware including the application of any image transformations or corrections applied to achieve a geometrically accurate image. It also entails a data collection strategy capable of rapid scanning of the OCT beam across the FOV so as to make optimal use of the available scanning mirror bandwidth and avoid dead time where the mirrors are scanning but not collecting useful imaging data [12].

Volumetric OCT imaging has previously been reported for ophthalmic OCT imaging systems employing a constant linear velocity spiral scanning (CVL-SC) pattern to drive a galvanometric scanner [12–15]. Scanning in a spiral pattern allows continuous, uninterrupted data collection and optimizes the available mirror bandwidth by scanning both mirror axes at similar speeds, resulting in efficient interrogation of a volume. In previous studies employing spiral scanning, data collected along the spiral beam trajectory was interpolated onto a uniform Cartesian grid for display either in post-processing using bi-linear interpolation [12] or in real time using nearest-neighbor interpolation [15].

Geometric distortion in OCT images is well documented, and many approaches to correcting for it have been reported [16–21]. However, to our knowledge, these algorithms have never been previously integrated into a single system capable of acquiring and rendering geometrically accurate OCT volumes, in real-time, in-vivo. Nor have they been applied to image the middle ear [21] or assessed for geometric accuracy against other gold standard volumetric imaging modalities such as computed tomography (CT).

In this study we present a novel system for real-time, geometrically accurate 4D ME-OCT imaging that images over a $10.9\,mm\, \times \,{30^\circ }\, \times \,{30^\circ }$ sector from a handheld, entocentric probe. OCT images are acquired using a two-axis, gimbal-less micro electromechanical system (MEMS) mirror driven to scan in a discretized spiral pattern, scan converted to correct fan-beam distortion, corrected for optical and mirror-related distortions, and displayed to the user in real-time using a custom-volumetric ray-cast renderer. We verify the system by quantitative measurement of geometric accuracy against benchtop phantom imaging targets and validate its use for quantitative imaging of the middle ear by comparing 3D OCT images obtained in human cadaveric middle ears to micro-CT images of the same ear. A cadaver model is also used to demonstrate 4D OCT visualization of a middle ear implant in a simulated diagnostic assessment scenario. Finally, we demonstrate the real-time imaging capabilities of the system by performing continuous, live, in-vivo imaging of a subject who performs a Valsalva maneuver that causes dynamic pressurization and conformational changes to the middle ear which are captured in real time 4D OCT.

2. Methods

2.1 ME-OCT imaging system

Image data in this study was collected using a custom OCT handpiece designed for middle ear imaging from the ear canal (Fig. 1(c–e)). OCT imaging optics consisted of a two-axis, $3.6\,mm$ diameter MEMS mirror (Mirrorcle A8L1.1, USA) [22] used to scan an OCT beam over ${\pm} {6.2^ \circ }$ and a 4-lens pupil relay that imaged the mirror axis onto the entrance pupil of a graded index (GRIN) rod lens that served as the objective. The GRIN lens imaged the mirror axis onto its exit pupil so that the exit pupil acted as a point of convergence for OCT scan lines. The entocentric design of the optics enabled scanning over ${30^\circ }$ in both lateral dimensions. Following the manufacturer’s recommendations [23], the $3\textrm{}dB$ bandwidth of the signal used to drive the MEMS mirror was limited to $200\,\textrm{Hz}$ on each axis. The control voltages for the MEMS mirrors were generated using a 10-bit digital-to-analog converter (DAC) and amplified by a high voltage amplifier to deliver a drive signal to the mirrors of $155\,{V_p}$. This imaging handpiece was coupled to a custom OCT engine (Figs. 1(a) and 1(b)) consisting of an akinetic swept source laser (Insight Photonics SLE-1001) with a center wavelength of $1550\,nm$, a sweep range of $35.4\; nm$, a nominal sweep rate of $200\,kHz$, and a 3-dB axial and lateral resolution in air of $40\,\mu m$ and $40\,\mu m$ respectively. While the coherence length of the swept laser permitted scanning over as much as $200\,\textrm{mm}$ [24], we limited the scan range to $10.9\,\textrm{mm}$ [25] which was adequate to image the full depth range of the human middle ear from the most lateral aspects of the tympanic membrane to the most medial visible parts of the cochlear promontory. A motorized delay line (MDL) was inserted into the reference arm to provide a variable delay that could move the image plane axially to accommodate for differences in patient ear canal length.

 

Fig. 1. Middle ear optical coherence tomography (ME-OCT) imaging system. (a) ME-OCT system schematic. (b) ME-OCT clinical system. (c) ME-OCT imaging at the point-of-care using a handheld imaging probe with an otoscopic form factor. (d) Computer-aided design (CAD) cross-sectional view of the handpiece. (e) Optical model of the handpiece’s entocentric imaging optics. IS, image sensor; MDL, motorized delay line; DAQ, data acquisition card; FPGA, field programmable gate array; DAC, digital to analog converter; BD, balance detector. C, three-way fiberoptic circulator.

Download Full Size | PDF

2.2 Discretized spiral scanning

Given the cylindrical symmetry of the imaging optics and the limited bandwidth of the MEMS mirror, a spiral scanning technique was developed to maximize the volumetric scanning rate and to match the extent of scanning to the exit pupil dimensions. This ensured that no scanning time was wasted imaging outside the system’s aperture or during times when no data collection occurred.

Two-axis, gimbal-less MEMS mirrors [26] are appealing for handheld OCT systems due to their small size and weight. While they lack position sensing capabilities, their position is highly repeatable when driven open-loop [27]. With an open-loop driving strategy a one-time calibration must be used to map drive signals to correct angular locations. The high repeatability of the mirrors then ensures that this calibration remains valid over time. In this section we present our approach to using a spiral scan pattern to drive a two-axis MEMS mirror to collect volumetric data.

The flowchart in Fig. 2(a) depicts the proposed DC-SC algorithm. Our spiral scan pattern is a discretized approximation to a constant linear velocity (CLV) spiral scan [12] as shown in Figs. 3(a) and 3(b) where the discretized spiral was generated using Bresenham’s algorithm [28]. Starting on the outer edge of the image FOV, we generate a discrete approximation of a circle of a fixed radius one quadrant at a time. After completion of the circular scan, the radius is decremented and scanning continues along a new circle with a radius one pixel smaller than the previous circle. This process continues until the center of the FOV is reached. The radius is then reset to acquire the next spiral volume. Each location of the discretized spiral scan (Fig. 3(b)) lies on a grid of uniformly spaced angles allowing the data to be organized into a discrete $({r,\theta ,\phi } )$ array. However, due to the spiral discretization, not all lateral locations in the volume are covered by scan lines resulting in the final constructed volume having “holes” where the density of scan lines is lower than in surrounding areas (e.g. along the diagonals in Figs. 2(b) and 3(b)). This effect can be addressed by applying a 2D median filter laterally across $({\theta ,\phi } )$ at each depth to homogenize the density of image lines. Figure 2(c) shows the result of applying such a median filter on the volume of Fig. 2(b).

 

Fig. 2. Discretized spiral scanning (DC-SC) algorithm and OCT volume smoothing. (a) Flowchart to generate a discretized spiral scan pattern based on the Bresenham’s circle drawing algorithm. Discretized spiral OCT volume before (b) and after (c) median filtering. $maxRadius$, control variable specified by the user to control the FOV of the spiral scan; $radius$, index variable tracking the spiral radius currently being generated and determine when the spiral is finished; $center$, center point of the spiral scan along the discretized grid; $quadrant$, index variable of the current circle quadrant being generated; x and y, circle indices along the x-axis and y-axis respectively; d, decision variable used to increment x and y such that circular error is minimized along the discretized grid; r, temporary variable to store the radius of the circle being generated; ${X_{DAC}}$ and ${Y_{DAC}}$, DAC output signals driving the MEMS mirror.

Download Full Size | PDF

 

Fig. 3. Discretized spiral scanning (DC-SC) pattern and waveforms. (a) MEMS mirror driving voltage waveforms for spiral scan. (b) En-face view of the discretized spiral scan over the angular FOV. For illustration, only a fraction of the spiral scan’s FOV was used to produce the waveforms (a) and (b). ${X_{DAC}}$ and ${Y_{DAC}}$, FPGA’s 10-bit DAC signals.

Download Full Size | PDF

As the Bresenham’s circle algorithm only makes use of integer math, it can be implemented using economical hardware without the need for lookup tables (LUTs). In our system we implemented the algorithm on a field programmable gate array (MachXO2, Lattice Semiconductor; USA) driving a 10-bit DAC to generate volumes containing $512 \times 512$ image lines corresponding to a ${30^\circ } \times {30^\circ }$ FOV in $4.5$ seconds ($0.22\; Hz$). A speed-up factor of 4 can be achieved for every factor of 2 decrease in acquired FOV diameter. For example, halving the FOV to ${15^\circ } \times \; {15^\circ }$ we can acquire volumes containing $256 \times 256$ image lines in $1.13$ seconds ($0.89\; Hz)$ and quartering the FOV to ${7.5^\circ } \times \; {7.5^\circ }$ we can acquire volumes containing $128 \times 128$ image lines in $0.28$ seconds ($3.56\; Hz)$. For our MEMS mirror, the manufacturer’s maximum recommended drive signal bandwidth was $200\; Hz$. With this bandwidth limitation we collect $N = 5$ A-lines per image line for all acquired volumes and average them together to form the line.

A limitation of our spiral scanning technique is that as the beam gets closer to the center of the spiral (i.e., to the center of the FOV) and the scanned radius gets smaller, the frequency of the drive signal increases as can be seen in Fig. 3(a). At some point, the frequency of the drive signal will approach the $200\; Hz$ bandwidth limit of the drive electronics, resulting in an inability of the MEMS mirror to follow the spiral path. This results in an artefact at the center of the spiral. The severity of the artefact increases as the volume scanning rate increases and the artefact takes up an increasing fraction of the image volume.

2.3 Real-time volumetric scan conversion using GPU texture memory

Because the optical system shown in Figs. 1(d) and 1(e) is designed to image the MEMS mirror rotation axis to the exit pupil, all OCT beamlines intersect at the pupil. Volumes are constructed by scanning the angles $({\theta ,\phi } )$ that the beamlines make about this common origin and OCT data is collected in spherical coordinates $({r,\theta ,\phi } )$ where r is the radial distance along the A-line and $({\theta ,\phi } )$ are the spherical angles about the exit pupil center. To display data with geometrically correct proportions, a transformation must be made from the spherical coordinates in which it was measured (Fig. 4(a)) to the Cartesian coordinates in which it is displayed (Fig. 4(b)). This type of image conversion from spherical to Cartesian coordinates is commonplace in phased array ultrasound imaging where it is called scan conversion [29], and we will use that terminology here.

 

Fig. 4. Spherical measurement space (a) in which OCT volume data is acquired and the Cartesian display space (b) in which it is displayed. $\hat{\boldsymbol{n}}$ MEMS mirror surface normal vector; $\hat{\boldsymbol{r}}$ handpiece’s optical imaging axis; $\alpha ,\; \beta $ tilt angles between $\hat{\boldsymbol{n}}$ and $\hat{\boldsymbol{r}}$ along $\theta $ and $\phi $ respectively; ${r_0},{r_1}$ radial extent of the spherical measurement space; ${\theta _{0,}},{\theta _1}$ and ${\phi _0},{\phi _1}$ are the lateral angular extents of the spherical measurement space; ${x_\tau },{y_\tau },{z_\tau }$ tilted Cartesian axis; ${x_0},{x_1}$ and ${y_{0,}},{y_1}$ the lateral extents of the Cartesian display space; ${z_0},{z_1}$ depth extent of the Cartesian display space.

Download Full Size | PDF

To perform real-time volumetric scan conversion for ME-OCT we extended the real-time 2D ultrasound scan conversion kernels in the open-source FAST toolkit [30]. FAST provides a set of common medical image processing algorithms including 2D scan conversion implemented as graphic processing unit (GPU) kernels to achieve very fast processing times by exploiting the massive parallelism of GPU architectures. GPUs provide dedicated hardware for interpolation between voxel values stored in texture memory resulting in a significant increase in processing throughput as compared to CPU-based approaches.

In our scan conversion method, the extent of the spherical measurement space $\{[{{r_0},{r_1}} ],[{{\theta_{0,}},{\theta_1}} ],$ $[{{\phi_0},{\phi_1}} ] \}$ is converted into a Cartesian cuboid of extents $\{{[{{x_0},{x_1}} ],[{{y_{0,}},{y_1}} ],[{{z_0},{z_1}} ]} \}$ that bounds the displayed data, where:

The display space bounded by these Cartesian extents is divided into an integer number of voxels ${D_x} \times {D_y} \times {D_z}$ where ${D_x},\; {D_y}$ and ${D_z}$ are the number of voxels along each axis. The integer voxel coordinates $({x,y,z} )$ of the display space are then normalized to the range $[{0,1} ]$ to form the normalized Cartesian coordinates $({{x_f},{y_f},{z_f}} )$:

Due to imperfect optical alignment, there may be a non-zero tilt between the normal of the MEMS mirror when the mirror drive signals are set to zero and the image axis. To account for this, we introduce a set of tilted floating-point coordinates $({{x_\tau },{y_\tau },{z_\tau }} )$ given by:

Finally, the measured data is interpolated onto the point in spherical space $({r,\theta ,\phi } )$ that corresponds to each set of voxel coordinates $({x,y,z} )$ of the display space [21]. All these steps are performed continuously in real-time using GPU textures at a rate of $52.8\; Hz$ ($18.9\; ms$/vol) for an OCT volume containing ${D_x} \times {D_y} \times {D_z} = 512 \times 512 \times 330$ voxels, more than enough to enable real-time 4D imaging.

A motorized delay line in the reference arm of our OCT engine was used to adjust the image depth to accommodate different ear canal lengths (see Fig. 1(a)). Changes to the reference arm delay effectively change the radial extent of the spherical measurement space from $[{{r_0},{r_1}} ]$ to $[{{r_0} + \mathrm{\Delta }/2,\; {r_1} + \mathrm{\Delta }/2} ]$ where $\mathrm{\Delta }$ is the amount of additional reference arm delay. The factor of $\frac{1}{2}$ arises because reference arm is one-way delay while sample arm delay is two-way. Our scan conversion algorithm is set up to automatically adjust for these changes by scaling the apparent width and radius of curvature of the proximal and distal boundaries of the image window while keeping it centered on the screen.

2.4 Custom-built GPU-based volumetric renderer

Generated 3D volumes are displayed using a custom real-time OpenGL based rendering engine utilizing a ray-casting technique similar to that described by Stegmaier et al. [31]. Ray-casting is a data-agnostic and computationally efficient rendering method in which each voxel within a volume is modeled as emitting and absorbing light rays [32]. The image displayed to the screen is generated by casting rays originating from each pixel on the screen, accumulating the optical spectra of each ray as it traverses a virtual volume containing the acquired OCT data, and finally compositing the accumulated spectra as a color to the originating pixel on the screen. Using our ray-casting engine, OCT volumes could be rendered at 60 frames per second (FPS) on a Microsoft Windows 10 PC with 16 GB of RAM, a Nvidia Titan XP GPU, and an Intel i7-6800 K processor.

2.5 Geometric calibration through minimization of residual surface fitting error

To apply the scan conversion method described in the previous section, we first calibrated the system to obtain estimates of the extent of the spherical measurement space, i.e., to determine $\{{[{{r_0},{r_1}} ],[{{\theta_{0,}},{\theta_1}} ],[{{\phi_0},{\phi_1}} ]} \}$ for a particular reference arm delay. We did this by imaging a flat, white acrylic plate placed near the distal end of the imaging window (Dataset 1 [33])$.$ A 3D image of the plate was then acquired, and the proximal plate surface was segmented within each image line using a peak-finding algorithm followed by filtering to smooth out the extracted surface in the lateral dimensions. The output from the segmentation were the radial indices ${S_{mn}}$ of the surface location for each discrete angle $({m,n} )$ in the spiral scan grid of Fig. 3(b). Points near the center of the image that were affected by the bandwidth limit artefact were excluded from the calibration volume. The radial extent of the imaging window $\mathrm{\Delta }R = ({{r_1} - {r_0}} )$ was measured by translating the reflector from the proximal to distal end of the window using a graduated translation stage until the image of the plate surface intersected the proximal and distal bounds of the image window along the central line. $\mathrm{\Delta }R$ was found to be $10.89\; mm$ in our system.

Defining the value of ${S_{mn}}$ along $\hat{\boldsymbol{r}}$ to be ${S_0}$, we construct a plane in Cartesian space perpendicular to the $z$-axis at a depth $z = {S_0}$. A given a set of parameters $\rho = \{[{{r_0},{r_1} = {r_0} + \mathrm{\Delta }R} ],[{{\theta_{0,}},{\theta_1}} ],$ $[{{\phi_0},{\phi_1}} ],[{\alpha ,\beta } ] \}$ describes a transformation of this plane into spherical coordinates. Applying this transformation according to the Eq. (2–4) of the previous section, we calculate the intersection of this plane in spherical coordinates with each scan line $({{\theta_m},\; {\phi_n}} )= (m\mathrm{\Delta }\theta $, $n\mathrm{\Delta }\phi )$ in the spiral scan. Here the angular grid spacing is $\mathrm{\Delta \theta } = \frac{{{\theta _1} - {\theta _0}}}{M}$ and $\mathrm{\Delta }\phi = \frac{{{\phi _1} - {\phi _0}}}{N}$ where M and N are the number of grid points in the $\theta $ and $\phi $ directions. The radial location at which the plane intersects each scan line is ${I_{mn}}(\rho )$. We then perform a least-squares minimization to determine the set of parameters $\rho $ that minimizes the error between ${I_{mn}}(\rho )$ and ${S_{mn}}$ over the scan to solve the minimization problem in Eq. (5).

This minimization process is illustrated in Fig. 5. The minimization was performed using the minimize function in the Python SciPy optimization toolbox [34] using the Sequential Least Squares Programming (SLSQP) method with default option values. The parameters $\rho $ obtained from this minimization were stored and applied to all subsequent real-time scan conversions performed in this study. The calibration was found to remain stable across the six-month period of the study and did not require updating or modification.

 

Fig. 5. Diagram describing the geometric calibration procedure for finding the scan conversion parameters $\rho = \{{[{{r_0},{r_1}} ],[{{\theta_{0,}},{\theta_1}} ],[{{\phi_0},{\phi_1}} ],[{\alpha ,\beta } ]} \}$ that optimize the fidelity of the OCT volumes with the known geometry of a flat surface. I is an ideal flat Cartesian surface reference; $I(\rho )$ is the calculated spherical surface after Cartesian-to-spherical coordinate transformation using parameters $\rho $; S is the measured spherical surface segmented and extracted from an OCT volume of a flat white acrylic plate; $\mathop \sum \limits_{m,n} \; |{S_{mn}} - {I_{mn}}{(\rho )|^2}$ is the sum of squared error between the surface model and measured surface.

Download Full Size | PDF

2.6 Lateral angular distortion correction

MEMS mirrors are known to exhibit inter-axis crosstalk in which the sensitivity of mirror tilt angle to applied voltage in the $x$ axis depends on the tilt angle in the $y$-axis and vice versa [34–36]. Using a similar MEMS scanner to the one in our system, Kim et al. determined angular error during dual axis scanning due to crosstalk to be as high as $6.13\%$ depending on MEMS scanner alignment [37]. Furthermore, the dependence of MEMS mirror angle on drive voltage is non-linear with mirror angle tending to saturate at the ends of the drive voltage range. The effect of mirror crosstalk and nonlinear distortion must be corrected to obtain geometrically accurate OCT volumes.

Our method for correcting distortion due to MEMS scanner nonlinearity and crosstalk is illustrated in Fig. 6. Taking a similar approach to Izatt et al. [12], we captured an en-face image of a uniform dot grid. Our dot grid consisted of an ultraviolet (UV)-printed black acrylic plate with a rectangular grid of circular dots (UJF-6042 Mk II, Mimaki Engineering, Japan) with a center-to-center spacing of $0.5\; \textrm{mm}$ and a dot diameter of $0.25\; \textrm{mm}$ printed at a resolution of $21\; \mu \textrm{m}\; \times \; 21\; \mu \textrm{m}$. The dot grid was imaged near the distal edge of the imaging window to maximize the number of dots visible in the image and an en-face image was generated from a scan-converted OCT volume (Dataset 2 [38]). The en-face image of the dot grid captures the observed 2D map of scan angles corresponding to the driving voltages of the spiral scan pattern after undergoing lateral distortion effects caused by the MEMS mirror inter-axis crosstalk and nonlinearity. For each dot in the en-face image, we extracted its apparent angular coordinate on the distorted grid using the Discorpy library [39]. With the number of dots in the grid known, a uniform dot grid consisting of the same number of dots was then created. A one-to-one relationship between the coordinates of the extracted and uniform grid was established to form an array of uncorrected/corrected coordinate pairs. Finally, a 2D correction map was created by interpolating this coordinate pair array across a normalized grid containing the same number of pixels as the original en-face image. With the lateral correction map stored as a GPU texture, each pair of corrected angular values is read from the correction map using its corresponding pair of drive voltage values $[{\theta ({{X_{DAC}},{Y_{DAC}}} ),\; \phi ({{X_{DAC}},{Y_{DAC}}} )} ]$ (Fig. 3) as texture coordinates and used during volumetric scan conversion to correct lateral distortion in real-time. The angular calibration only needed to be performed at a single depth to correct lateral distortions at all depths.

 

Fig. 6. Correction of lateral angular distortion caused by MEMS nonlinearity and interaxis crosstalk. An en-face image of a reference dot grid exhibiting lateral distortion was generated from a scan converted OCT volume. The dots were segmented, extracted, and matched to a uniform grid creating a correction map for correction of the observed lateral distortions in real-time during scan conversion. (Center-to-center dot spacing of $0.5\; \textrm{mm}$, dot diameter of $0.25\; \textrm{mm}$)

Download Full Size | PDF

2.7 Geometrical comparison of ME-OCT and co-registered micro-CT volumes

The geometric accuracy of corrected OCT volumetric data was compared to nominal design dimensions of a 3D printed phantom target and to co-registered micro-CT images of a cadaveric middle ear. Phantoms were designed in Fusion 360 and printed in photopolymer resin with a UV stereolithographic 3D printer (Form 3, Formlabs, USA) using a $50\; \mu m$ layer thickness. The phantom consisted of a tube with holes in the circumferential wall and a four-step spiral staircase pattern at the distal end. The STL file for the phantom target can be found in Code 1 [40].

For cadaveric middle ear imaging, a fresh-frozen human cadaveric temporal bone, the right ear from a 30-year-old female was cut en-bloc, thawed, mounted, and prepared with the otic capsule, middle ear, ear canal, tympanic membrane and soft tissues of the ear canal and pinna left intact. The tympanic membrane was debrided and wax, hair, and other debris were removed from the ear canal to provide a clear view of the tympanic membrane. The OCT view of the tympanic membrane in this preparation is typical of that obtained in live subjects. During the experiment and immediately before micro-CT imaging, the ear canal was irrigated with saline and drained to maintain hydration of the tympanic membrane. Procedures used in this study were approved by the Dalhousie University Research Ethics Board under protocol 2021-5486.

Ground truth anatomy for the temporal bone and phantoms was obtained using a Triumph LabPET4/CT (Trifoil Imaging, USA) micro-CT imaging system (volume shape: $512 \times 512 \times 512$) with a slice thickness and FOV of $54\; \mu \textrm{m}$ and $27.65\; \textrm{mm}$, respectively for the phantom target and $155\; \mu \textrm{m}$ and $79.6\; \textrm{mm}$ for the cadaveric middle ear. Both micro-CT and OCT volumetric data was exported to DICOM format files for offline co-registration in 3DSlicer [41]. Artefacts were removed from the OCT volume using a previously described convolutional dictionary learning filter [42] and the micro-CT volume was cropped to a region of interest (ROI) that overlapped with the OCT volume. The two datasets were aligned using semi-manual, rigid, fiducial registration using eight boney anatomical landmarks visible in both datasets. The geometric accuracy of the OCT volume was evaluated using root mean squared error (RMSE) by taking the micro-CT images as ground truth.

3. Results

3.1 Verification of geometric and lateral angular distortion correction

The calibration process of Section 2.5 gave spherical imaging extents for our system of $\mathrm{\Delta }R = ({{r_1} - {r_0}} )= 10.89\; \textrm{mm},\; \; ({{\theta_1} - {\theta_0}} )= {30.6^ \circ }$ and $({{\phi_1} - {\phi_0}} )= {31.9^ \circ }$ and tilt angles $({\alpha ,\beta } )= ( - {0.4^\circ },\; - {0.4^\circ })$. To assess the accuracy of the correction, it was applied to a set of OCT volumes of the flat acrylic plate translated to 19 different depths covering the axial range of the image window from ${r_0}$ to ${r_1}$. After scan conversion, the RMSE surface flatness across all depths was $69\; \mu m$. When the lateral distortion correction described in Section 2.6 was also applied, as shown in Fig. 7(e), the RMSE surface flatness improved to $42\; \mu \textrm{m}$. The improvement in RMSE error from lateral distortion correction arises because of coupling between angular and axial errors towards the edges of the FOV.

 

Fig. 7. (a) En-face image of a reference dot grid ($0.5\; \textrm{mm}$ center-to-center spacing, dot diameter of $0.25\; \textrm{mm}$) generated from a scan-converted discretized spiral OCT volume before lateral distortion correction [38]. (b) Heatmap of angular error (degrees) between the measured grid points and the true grid point locations. (c) En-face image of the dot grid after lateral distortion correction (Dataset 3 [43]). (d) Residual angular error (degrees) between the measured grid points in (c) and the true grid point locations. (e) Residual radial error $({\mu m} )$ after scan conversion and lateral distortion correction (Dataset 4 [44]). (f) Residual total lateral and radial error $\delta = \sqrt {\mathrm{\Delta }{r^2} + r\mathrm{\Delta }{\theta ^2} + r\mathrm{\Delta }{\phi ^2}} $ $({\mu m} )$ after scan conversion and lateral distortion correction.

Download Full Size | PDF

Figures 7(a) and 7(c) show scan-converted OCT volumes of the reference dot grid acquired using our discretized spiral scanning algorithm before and after lateral distortion correction. In addition to the distortion effects that we correct for, Figs. 7(a) and 7(c) also show clear evidence of the astigmatism of the GRIN rod objective in the diagonal smearing of points into streaks near the corners of the FOV.

Figures 7(b) and 7(d) show heat maps of the absolute error between the measured dot grid before and after lateral distortion correction. Prior to correction, angular error was highest around the corners of the FOV where the effects of both mirror crosstalk and non-linearity are strongest. Distortion was also high at the center of the image due to the mirror bandwidth limitation artefact. The RMSE angular error before correction was found to be ${1.15^\circ }$ ($383\; \mu m$ at the distal end of the image). After lateral distortion correction, as Fig. 7(d) shows, the RMSE residual angular error was found to be ${0.036^\circ }$ ($12\; \mu \textrm{m}$ at the distal end of the image).

Total residual angular and radial error shown in Fig. 7(f), was largest around the center and the diagonal corners of the FOV. These residual errors arise mainly from inaccuracies in correct identification of dot coordinates in the image due to astigmatic blurring and the presence of the central artefact. Even with these error sources, the worst-case total error after scan conversion and lateral distortion correction was $70\,\mu m$ as shown in Fig. 7(f). However, $92\%$ of pixels had an accuracy of $40\; \mu m$ or better indicating that the residual geometric error was less than the lateral resolution over most of the FOV. Figures 7(c) and 7(d) also show that after lateral distortion correction the FOV of the OCT volume is no longer circular since the crosstalk introduces anisotropic angular range limits with a larger angular range along the diagonals than along the mirror axes. A summary of the ME-OCT imaging system’s performance for various FOVs after geometrical and lateral distortion correction can be found in Table 1.

Table 1. Benchmarks of the ME-OCT imaging system for various field of views (FOV) including FOV dimension, volumetric acquisition rate, 3-dB lateral and axial resolution in air, and residual error after geometric and lateral distortion correction. RMSE, root mean squared error.

View Table | View all tables in this article

We verified the geometric accuracy of scan-converted and lateral distortion corrected OCT volumes by comparison against nominal design dimensions of a 3D printed phantom target [40] as shown in Fig. 8 and by the similarity of the OCT volume to a co-registered CT volume. Figure 8(a) shows overlaid OCT and micro-CT images of the phantom. Figure 8(b)–(d) show measurements made between various landmarks within the OCT volume in each view plane which were compared to the nominal dimensions in the design file used to 3D print the phantom. The comparison was made against the phantom’s design file [40] instead of the micro-CT volume because the 3D printer manufacturing tolerance of $50\; \mu \textrm{m}$ was less than the resolution and segmentation error of the micro-CT images. The design dimensions and the corresponding measurement obtained from the OCT volume are reported in Table 2. The dimensions derived from the corrected OCT images agree with the design dimensions to better than the $40\; \mu m$ lateral resolution near the central FOV of the system indicating that scan conversion and lateral distortion correction have effectively eliminated technical sources of geometric inaccuracy within the region of interest. The mean and standard deviation given in Table 2 are calculated across multiple independent attempts by a single assessor to identify dimensional boundaries in the OCT images in order to calculate the measured dimensions.

 

Fig. 8. Micro-CT (greyscale) and OCT (color overlay) of a 3D printed phantom target [40] viewed in 3DSlicer. (a) 3D render of micro-CT (green) and OCT (yellow) co-registered volume segmented contours. (b), (c), and (d) show three orthogonal cross-sectional views (labeled axial, coronal and sagittal) of the co-registered volume shown in (a). White arrows and outlines are the dimensional measurements in Table 2 made to assess the geometric accuracy of the OCT volume. A flythrough video of the co-registered volume’s image stack, seen from the coronal plane (c), and 3DSlicer scene file of the co-registered volumes can be found in the supplementary materials Visualization 1 and Code 2 [45] respectively.

Download Full Size | PDF

Table 2. Measured dimensions in the co-registered OCT images of the phantom target [45] compared against nominal dimensions taken from the phantom’s design file [40]. RMSE: root mean squared error.

View Table | View all tables in this article

3.2 Validation of geometrical accuracy in a cadaveric middle ear

The ability of our scan conversion and lateral distortion correction procedure to produce geometrically accurate OCT volumetric images of the middle ear was validated by comparison of OCT and micro-CT volumetric images in the temporal bone preparation described in Section 2.7. Figure 9 shows slices through the OCT volume overlaid on corresponding micro-CT slices. Four anatomical dimensions, measured from the incudostapedial joint (ISJ) to the umbo, the ISJ to the promontory, the umbo of the malleus to the promontory, and the ISJ to the TM along the lateral-medial axis were measured on the OCT and micro-CT volumes separately as shown in Fig. 9(c) and are compared in Table 3. These dimensions are clinically relevant because they quantify the capacity of the middle ear to accommodate middle ear micro manipulations during surgery such as prothesis placement. Each measurement was taken independently five times in each modality by the same assessor. The mean and standard deviation reported in Table 3 were calculated across these five independent measurements on the OCT dataset. The RMSE was calculated between OCT and micro-CT datasets across the five independent measurements. The comparison shows that scan-converted and distortion-corrected OCT can achieve geometric accuracy exceeding the micro-CT resolution of $155\; \mu m$. The lack of ISJ to TM measurement for micro-CT is due to poor soft tissue contrast of micro-CT which made it impossible to accurately identify the TM surface in the micro-CT data. By comparison, OCT has excellent soft tissue contrast that makes the TM highly visible in the fused images. A flythrough video of the fused OCT and micro-CT image stack is included as supplementary material (Visualization 2).

 

Fig. 9. Selected slices from a fused micro-CT (greyscale) and OCT (color overlay) volume of a cadaveric human middle ear. (a) 3D render of the OCT volume. (b) Axial plane view. (c) Coronal plane view with measurement locations 1–4 marked. (d) Sagittal plane view. CTN, chorda tympani nerve; TM, tympanic membrane; ISJ, incudostapedial joint; S, stapes; P, promontory; I, incus; U, umbo; M, malleus. Measured dimensions: (1) ISJ to umbo, (2) ISJ to promontory, (3) Umbo to promontory, and (4) ISJ to TM. A flythrough video of the co-registered volume’s image stack, seen from the coronal plane (c), and 3DSlicer scene file of the co-registered volumes can be found in the supplementary materials Visualization 2 and Code 3 [46] respectively.

Download Full Size | PDF

Table 3. Clinically relevant dimensions obtained from measurements on the micro-CT and OCT volumes of the cadaveric middle ear [46]. ISJ, incudostapedial joint; TM, tympanic membrane; RMSE, root mean squared error.

View Table | View all tables in this article

The ability of corrected OCT to produce geometrically accurate measurements was further validated using a second cadaveric right middle ear from a 68-year-old male. The facial recess in the bone was opened through a mastoidectomy, the incus was removed and a $2.5\; \textrm{mm}$ long Dresden clip partial ossicular replacement prosthesis (PORP) (Kurz, Germany) with a $2.6\; mm$ diameter head was implanted between the stapes superstructure and a cartilage graft placed against the medial side of the TM, simulating an ossiculoplasty, a common surgical treatment for conductive hearing loss. The cartilage was optically cleared using the technique described in [47] and imaged with OCT. Figure 10 shows an OCT volume of the prepared ear imaged transtympanically following scan conversion and lateral distortion correction. The PORP’s diameter as measured in the OCT volume using the distance measurement tool in 3DSlicer [41] was found to be $2.61\; \textrm{mm} \pm 34\; \mu \textrm{m}$, consistent with the specified prosthesis head diameter of $2.6\; \textrm{mm}$ to within the system’s optical resolution. OCT imaging of prostheses is of interest to verify coupling and placement of the prosthesis following middle ear surgery, especially in cases where hearing results are suboptimal [5].

 

Fig. 10. Cadaveric middle ear containing a 2.6 mm diameter Dresden clip PORP crimped onto the stapes head in a simulated ossiculoplasty. (a) OCT volume of the PORP obtained transtympanically through an optically cleared cartilage graft (Dataset 5 [48]). Annotations show the measured PORP diameter mean and standard deviation (b) Manufacturer’s diagram of the PORP specifying the nominal head diameter. PORP, partial ossicular replacement prosthesis.

Download Full Size | PDF

To validate the ability of scan conversion and lateral geometric correction to be applied to 4D OCT imaging in real-time, the right middle ear of a healthy 44-year-old male volunteer was imaged while the volunteer performed a Valsalva maneuver. The Valsalva maneuver consists of pinching the nose and attempting to exhale through the nose to pressurize the middle ear space through the Eustachian tube. Patients are often asked to perform the Valsalva maneuver during clinical ear examinations to qualitatively assess tympanic membrane compliance and middle ear ventilation. 4D datasets of the Valsalva maneuver were acquired at full, half and quarter FOVs corresponding to three different volume rates. Videos of the acquired 4D datasets (Visualization 3) can be found in the supplementary materials. In all three acquisitions, the three-dimensional displacement of the TM caused by the Valsalva maneuver is clearly visualized in real time.

4. Discussion

We have introduced a novel system for live, continuous, and geometrically accurate 4D ME-OCT imaging at the point-of-care using the novel DC-SC spiral scanning technique, real-time volumetric scan conversion and real-time lateral angular distortion correction. DC-SC offers a promising method for scanning over large volumes like the middle ear, by making efficient use of available scanning bandwidth and enabling higher volume rates than is possible with sequential B-mode scanning. The imaging technique is computationally efficient and on standard hardware can run many times faster than data is acquired in even the fastest reported 4D OCT systems [15].

For our system, the combination of scan conversion and 2D lateral distortion correction was adequate to achieve resolution-limited geometric accuracy over $92\%$ of a ${30^\circ } \times {30^\circ } \times 10.9\; \textrm{mm}$ FOV. Geometric distortion along the edges of the FOV arises due to the difficulty in correctly segmenting calibration grid points in the presence of optical astigmatism. Methods for correcting astigmatism in OCT volumes have previously been reported [49–51] and their combination with the techniques of this study could potentially improve image accuracy at the edges of the FOV. Accuracy was also low in the central few pixels of our image due to mirror bandwidth limitation and central artefacts from internal reflections in our handpiece optics. The bandwidth-related artefact could be improved by slowing down the mirror scan at the expense of a reduced volume rate and the artefact from internal reflections could likely be reduced through optical design improvements in which case the center of the image volume should have the same accuracy as other nearby locations.

Through co-registration of OCT and micro-CT volumes we have demonstrated that the geometric accuracy of our corrected OCT volumes for measuring anatomical distances is $52\; \mu \textrm{m}$ to $75\; \mu \textrm{m}$. This estimate of OCT’s accuracy is likely an underestimate since the $155\; \mu \textrm{m}$ resolution of the micro-CT system likely contributed significantly to the measurement RMSE between OCT and micro-CT images. Even a $52\; \mu m$ to $75\; \mu m$ accuracy is substantially better than the resolution of typical clinical temporal bone CT ($400\; \mu m$) or MRI ($300\; \mu m$) systems [1] and similar to the resolution of binocular surgical microscopes (${\sim} 35\; \mu m$) [52]. The limited FOV of the ME-OCT and its inability to see through bone, means that OCT cannot replace clinical CT for middle ear visualization. However, this study shows that OCT can offer images of the middle ear areas visible through the TM with excellent geometric fidelity. We also note that the techniques we used to compare OCT and micro-CT could also be applied to fuse corrected OCT images with clinical CT images. The higher resolution and soft tissue contrast of OCT could potentially enhance CT images or reduce the number of CT images needed to track longitudinal changes to the middle ear over time.

Geometrically accurate 4D OCT imaging is an appealing capability for the middle ear whose complex 3D structure cannot be fully appreciated with 2D B-mode cross-sectional imaging alone. In addition to clinical diagnostics, there may be a role for 4D ME-OCT in surgical guidance including robotic surgeries [53–55]. Surgical guidance applications have been the primary motivator of recent developments in ophthalmic 4D OCT imaging technology [15].

Using DC-SC, we have demonstrated applying corrected 4D ME-OCT in vivo for visualization of conformational changes in a healthy normal ear during a Valsalva maneuver at three different FOVs (Visualization 3), achieving volume rates of $0.22\; \textrm{vol}/s$ over the full ${30^\circ } \times {30^\circ } \times 10.9\; mm$ FOV and $3.56\; \textrm{vol}/s$ over a limited ${7.5^\circ } \times {7.5^\circ } \times 10.9\; mm$ FOV. The achievable volume rate with DC-SC is limited by the bandwidth of our scanning mirror. In our system design we use a large mirror with a large scanning range to achieve a large number of resolvable spots at the expense of scanning rate, but alternative designs could employ a smaller mirror or reduced scan range to achieve higher volume rates. More advanced driving waveforms could also potentially increase volume rates by operating the MEMS mirror closer to its resonance frequency ($450\; Hz$ for the A8L1.1). The processing volume rate of $52.8\; Hz$ that we demonstrated for our technique provides substantial headroom to support faster systems than the one presented here.

5. Conclusion

We have demonstrated techniques for real-time spiral scanning, scan conversion and correction of lateral distortion and crosstalk from a MEMS scanning mirror in a handheld, entocentric OCT imaging system designed for imaging the middle ear. Analysis of residual errors in corrected OCT images of grid targets, 3D phantom targets and cadaveric temporal bones consistently show geometric accuracy similar to the lateral resolution over most of the FOV. The techniques presented are simple and of low computational complexity, making them suitable for application in real-time 4D imaging. We demonstrate this 4D capability by performing real-time volumetric imaging on the ear of a live subject during a dynamic pressurization maneuver.

Having geometrically accurate images is important for many applications in surgical planning, multimodal image fusion, objective diagnostics, and longitudinal monitoring of disease progression. The techniques presented in this study may help to make OCT a competitive imaging modality for many of these applications.