Multi-contrast Jones-matrix optical coherence tomography -- the concept, principle, implementation, and applications

Jones-matrix optical coherence tomography (JM-OCT) is an extension of OCT that provides multiple types of optical contrasts of biological and clinical samples. JM-OCT measures the spatial distribution of the Jones matrix of the sample and also its time sequence. All contrasts (i.e., multi-contrast OCT images) are then computed from the Jones matrix. The contrasts obtained from the Jones matrix include the conventional and polarization-insensitive OCT intensity, cumulative and local phase retardation (birefringence), degree-of-polarization uniformity quantifying the polarization randomness of the sample, signal attenuation coefficient, sample scatterer density, Doppler OCT, OCT angiography, and dynamic OCT. JM-OCT is a generalized version of OCT because it measures the generalized form of the sample information; i.e., the Jones matrix sequence. This review summarizes the basic conception, mathematical principle, hardware implementation, signal and image processing, and biological and clinical applications of JM-OCT. Advanced technical topics, including JM-OCT-specific noise correction and quantity estimation and JM-OCT's self-calibration nature, are also described.


I. INTRODUCTION
A. Optical coherence tomography and Jones matrix optical coherence tomography O PTICAL coherence tomography (OCT) [1] is an optical measurement methodology that primarily measures the properties of a probe beam backscattered from a sample. The properties of the sample are then obtained from the properties of the probe beam. The probe beam can be characterized in terms of its amplitude, phase, and polarization. Conventional OCT images are obtained from the amplitude of the probe beam, allowing visualization of the morphology of the sample. Using the phase information of the beam, tiny structures and alterations of the probe path length smaller than the coherence-length-limited resolution can be measured [2], [3].
This work was supported in part by the Japan Science and Technology Agency (JPMJCR2105, JPMJMI18G8) and Japan Society for Promotion of Science (21H01836, 22K04962, 21K09684, 18K09460).
Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org.
Digital Object Identifier 10.1109/JSTQE.xxxx.xxxxxxx The polarization of the probe beam is measured by hardwareextended OCT, polarization sensitive OCT (PS-OCT) [4]- [7]. Through the probe-beam polarization, the polarization property of the sample can be measured.
Time sequential measurement can provide further information. As an example, the flow of a fluid in a tissue can be measured by analyzing the temporal change in the phase of the probe beam. This modality is known as Doppler OCT [8]- [13]. Similarly, tiny internal displacements in a sample can be measured from the phase change. Such measurements have been made to quantify the mechanical properties of tissues, such as in optical coherence elastography [14]- [16] and the assessment of laser-induced dynamic processes of tissues [17], [18]. Blood flow can also be visualized by measuring the statistical temporal fluctuation of the OCT signal amplitude, phase, or both. This modality is referred to as OCT angiography (OCTA) [11], [19], [20] and is widely used in the daily clinical routine of ophthalmology [21]. Slow dynamics of the tissue, such as those induced by tissue metabolism, have been assessed using similar methods, which is referred to as dynamic OCT (D-OCT) [22]- [27].
Almost all properties of the probe beam, including its amplitude, phase, and polarization, can be represented using a Jones vector [28]. The time sequence of Jones vectors may contain the temporal properties of the probe beam. Therefore, by measuring the sequence of Jones vectors, we can conduct not only PS-OCT but also Doppler OCT, OCTA, and D-OCT.
Although OCT directly measures the properties of the probe beam, the real interests of users are the properties of the tissue. The properties of the probe beam are represented by Jones vectors, whereas the tissue properties are represented by the Jones matrix. Jones-matrix OCT (JM-OCT) is a modality that primarily measures the Jones vectors, computes the Jones matrix, and gives multi-contrast images of the sample [29]. The scattering and polarization properties of samples can be determined from the Jones matrix. By sequentially measuring the Jones matrices, JM-OCT gives temporal properties of the sample, such as the flow and intracellular dynamics.
Notably, the Jones vector and matrix are generalized representations of the probe beam and the sample's optical properties. The time sequential Jones matrix is a more generalized sample representation because it represents the time dynamics of the sample in addition to the optical properties. JM-OCT is thus the most generalized form of OCT. This conception is discussed in more detail in Section II. Intensity OCT images are obtained through coherent composition or intensity averaging. A polarization-artifact-free attenuationcoefficient image is obtained from the intensity average image. An attenuationcoefficient image and a scatterer-density image are obtained by processing the intensity OCT image. The latter is obtained using a neural-network-based scatterer density estimator. Two types of label-free tissue dynamics image (dynamic OCT; D-OCT) are obtained by processing the time sequence of the intensity images with two algorithms, namely the logarithmic-intensity variance (LIV) and OCT correlation decay speed (OCDS) algorithms. These algorithms are sensitive to the magnitude and speed of the tissue dynamics, respectively. Through eigenvalue decomposition, cumulative-phase-retardation and cumulative-diattenuation images are obtained. A local Jones matrix is computed from the original Jones matrix, and birefringence and localdiattenuation images are obtained from the eigenvalues of the local Jones matrix. The degree-of-polarization uniformity (DOPU) is defined as a circular variance of Stokes vectors, which are obtained from the Jones matrix. Doppler OCT is obtained through temporal phase analysis and OCTA through temporal signal variation analysis.

B. What is and is not covered by this review
This review summarizes the basic conception (Section II), principle (Section III), implementation (Section IV), and applications (Section V) of JM-OCT. Sections VI and VII discuss additional technical topics and further extensions of JM-OCT.
It should be noted that, in this review, JM-OCT mainly refers to the modality developed at the University of Tsukuba, which is characterized by not only the Jones-matrix-based polarization measurement but also a multi-contrast imaging capability. PS-OCT technology in its widest sense, including single-input-polarization PS-OCT [4], [30]- [34], Stokes-based PS-OCT [5], and Jones-Muller-analysis-based PS-OCT [35]- [37], is not fully covered in this review and is only summarized briefly in Sections VI-E and VII. The reader can find good reviews of PS-OCT in the literature [6], [7]. One topic that is not thoroughly covered by these PS-OCT review papers is recently developed technologies for measuring the axis orientation [34], [37], [38]. This topic is briefly summarized in Section VII-D but extensive review of this topic is awaiting.

II. BASIC CONCEPTION OF JM-OCT
The basic conception of JM-OCT is presented in Fig. 1. JM-OCT measures a sequence of Jones matrices of the sample and regards it as the generalized representation of the sample. Here, each of four entries of the Jones matrix is a complex OCT image. In general, a Jones matrix is a complex 2-×-2 matrix that transforms an incident Jones vector into an output Jones vector. In principle, the four entries of the Jones matrix can thus be determined by measuring two output Jones vectors corresponding to two incident Jones vectors.
In JM-OCT, a polarization diversity (PD) detection, which splits the probe beam into two orthogonal polarization states and measures two interference signals, is adopted to measure the output Jones vector. Here, each entry of the Jones vector is a complex OCT image corresponding to one of the two orthogonal polarizations of the backscattered probe beam. In addition to the PD detection, two orthogonal polarization states are multiplexed in the incident probe beam by some means (see Section IV for details). JM-OCT thus measures two output Jones vectors corresponding to two incident polarization states; i.e., incident Jones vectors. Each Jones vector comprises two complex OCT images, and four complex OCT images are thus obtained in a single measurement. As detailed later in Section III, these four OCT images form a "measured Jones matrix (J m )," which is mathematically similar to the Jones matrix of the sample. By repeating this measurement, a time sequence of the Jones matrices is obtained. This Jones-matrix sequence is used to compute several contrasts.
By applying a coherent composite method to the four entries of the Jones matrix, a high-sensitive OCT image is obtained (Section 3.6 of Ref. [39]). Intensity averaging of the two or four entries of the Jones matrix gives a polarization insensitive (i.e., birefringence-artifact-free) OCT image. By processing this polarization sensitive OCT image, we also obtain an attenuation coefficient image [40], which is a more direct representation of the tissue's scattering property than the raw OCT image. The scatterer density image is a more direct measure of the tissue property than the attenuation coefficient and is obtained through a neural-network-based estimator.
The Jones-matrix images are obtained as a time sequence, and the intensity images can thus be obtained as a time sequence. By processing the time-sequential OCT intensity image, two types of label-free tissue dynamics image (D-OCT) are obtained.
As the phase difference between two eigenvalues of the Jones matrix, a classical polarization sensitive OCT image (i.e., a cumulative phase retardation image) is obtained [41]. By converting the Jones matrices into the local Jones matrix before computing the eigenvalues, a local phase retardation image, or equivalently a birefringence image, is obtained [42]. Notably, the conventional cumulative phase retardation represents the polarization properties of the backscattered probe beam, whereas the birefringence directly represents the polarization properties of the sample. The birefringence image highlights fibrous tissues, including collagenous tissues, fibrosis, and nerve fibers.
In addition to the cumulative phase retardation and the birefringence, cumulative and local diattenuations can be computed from the cumulative and local Jones matrices.
The uniformity of the polarization can also be computed from the Jones matrix. For example, the DOPU can be obtained through Stokes analysis of the Jones matrix [43]. The DOPU was originally introduced to non-Jones-matrix PS-OCT Fig. 2. Model of the sample to be measured. Here the sample is modeled as having two layers with thicknesses of ∆z 01 and ∆z 12 . P 0 to P 2 are depth positions. J sij indicates the Jones matrix for a single trip from P i to P j . J mi is the raw Jones matrix measured at P i . The figure is reprinted from [48]. [44] and is sensitive to melanin [45].
A Doppler OCT image is obtained by analyzing the phase evolution of the Jones matrix entries or the global phase of the matrix [39], [46]. In addition, an OCTA image can be obtained by quantifying the temporal fluctuation of the Jones matrix through correlation analysis [47].
The details of the image generation of each contrast are described in the next section.

III. PRINCIPLE OF JM-OCT
A. Cumulative Jones matrices and cumulative-phaseretardation measurement The first step of the JM-OCT measurement is measuring a combined Jones matrix of the sample and the JM-OCT system. Subsequently, a Jones matrix having the identical eigenvalues with the sample Jones matrix, i.e., a similar matrix of the sample Jones matrix, is computed. To mathematically describe this process, we first model the sample as depicted in Fig. 2. Here, the sample is modeled as having two layers. The first layer is a superior layer from the tissue surface P 0 to a depth P 1 and has a thickness of ∆z 01 . The second layer is a thin layer from P 1 to P 2 and has a thickness of ∆z 12 . Note that the first layer can be thick and comprise multiple tissue layers whereas the second layer is assumed to be thin as it can be regarded as homogeneous along the depth. The Jones matrices for forward and backward passes of a single trip between P i and P j are denoted J sij and J sji , respectively. The raw Jones matrices measured at P i are denoted J mi .
Using Jones matrix hardware, which will be described later in Section IV, J mi is obtained as a function of space; i.e., an image of the Jones matrix. The measured Jones matrix at the depth position P i is expressed as where J out and J in are the Jones matrices of the collection and illumination optics, respectively. These matrices are not a function of the space; i.e., they are constant matrices. Hereafter, the two subscripts i and j of J rij indicate the indexes of depth positions, and J rij is the Jones matrix for a round trip from P i to P j . That is to say, J r0i is a round-trip sample Jones matrix from the surface (P 0 ) to a depth position (P i ), and J r0i = J si0 J s0i = J T s0i J s0i . Similar to J rij , J sij is the Jones matrix for a single trip from P i to P j . Here, we use the fact that J si0 = J T s0i . The superscript T denotes a transpose. In particular, the measured Jones matrix at the surface becomes a constant matrix A cumulative Jones matrix at P i is then computed as the product of the measured Jones matrix at P i and the inverse of the surface Jones matrix where J ci is the cumulative Jones matrix. As shown on the right-hand side of this equation, the cumulative Jones matrix is a similar matrix to the round-trip Jones matrix from the sample surface to P i , J r0i . Owing to this similarity, the eigenvalues of the cumulative Jones matrix are identical to those of the roundtrip sample Jones matrix. The phase retardation, which is the phase difference between two eigenvalues of a Jones matrix, of this cumulative Jones matrix is thus identical to that of the round-trip sample Jones matrix. Finally, a cumulative phase retardation, which is one of the polarization sensitive quantities measured by conventional PS-OCTs, is defined as the phase retardation computed from this cumulative Jones matrix. Note that although the cumulative Jones matrix comprises the Jones matrices of the system optics (J out ), the cumulative phase retardation solely relates to the sample Jones matrix because of the matrix similarity. A more detailed formulation of the cumulative Jones matrix and cumulative phase retardation computation can be found in Sections 3.2 and 3.3 of Ref. [39] and Section 2.2 of Ref. [48].

B. Local Jones matrices and local-phase-retardation measurement
Although cumulative phase retardation has been widely used for PS-OCT imaging, it is not really interpretable unless the tissue between the tissue surface and the depth of interest has a uniform axis orientation. Namely, if the axis orientation is not uniform, the eigenvalues of the cumulative Jones matrix cannot have a straight-forward relationship with those of the corresponding single-trip sample Jones matrix. It makes the cumulative phase retardation uninterpretable (see Section 2.2 of Ref. [48] for details).
To overcome this limitation and obtain a more depthlocalized polarization property, a local Jones matrix analysis has been introduced [42]. Assume the second layer in Fig. 2 (from P 1 to P 2 ) is the layer of interest. The round trip local Jones matrix of this layer is defined and computed as As shown in this equation, the local Jones matrix is defined from the cumulative Jones matrices mathematically (center part of the equation), but it can be directly computed from the measured Jones matrices (right hand side of the equation). The local Jones matrix is written in a decomposed form as where A is a matrix defined from the Jones matrices of the collection optics and the sample. Λ s12 is an eigenvalue matrix of the single-trip Jones matrix J s12 , which is the Jones matrix of the layer of our interest. The eigenvalue matrix is a diagonal matrix whose diagonal entries are the eigenvalues of J s12 . J L12 is a matrix similar to Λ 2 s12 as shown in this equation, and therefore no matter the values of A, the phase retardation obtained from the local Jones matrix J L12 is identical to that obtained from Λ 2 s12 . Λ s12 is a diagonal matrix by definition, and this phase retardation is thus twice the phase retardation of Λ s12 and is a round-trip phase retardation of the layer of interest. This round-trip phase retardation is referred to as the local phase retardation of the thin layer from P 1 to P 2 .
If the birefringence and optic axis orientation are uniform between P i and P j , the tissue birefringence of this layer, b i j , is obtained as where LPR i j is the local phase retardation obtained from Λ s12 , k 0 is the center wavenumber of the probe beam, and ∆z i j is the thickness of the layer. A more detailed formulation of the local Jones matrix and local phase retardation computation can be found in Section 2.3 of Ref. [48]. Note that the correct computation of the birefringence requires the uniformity of the layer of interest. This requirement is, however, weaker than that of the measurement of the cumulative phase retardation. That is to say, the measurement of cumulative phase retardation requires uniformity in the layer between the tissue surface and the depth of interest, which is notably thicker than the "layer of interest" of the local phase retardation measurement. A more detailed derivation and requirement (assumptions) for the cumulative phase retardation and the birefringence measurement can be found in Sections 2, 6.1, 6.3 and Appendix A of Ref. [48].

C. Diattenuation
Diattenuations can be computed from the eigenvalues of a Jones matrix [49], [50]. In the case of JM-OCT, cumulative and local diattenuations can be computed from the cumulative Jones matrix [Eq. (3)] and local Jones matrix [Eq. (4)], respectively, as where d is the diattenuation, λ 1 and λ 2 are two eigenvalues of the cumulative or local Jones matrices [6], [29]. Note that the diattenuation of biological tissue is negligible in general [41], [51]- [54]. And hence, the examples of diattenuation imaging are not shown in Section V (Applications).
The DOPU is defined as the variance of the Stokes vectors on a Poincaré sphere, where the Stokes vectors are those measured within a small region in the image (i.e., the spatial kernel) and is mathematically defined as where Ī ,Q,Ū,V denotes the averaged Stokes vector within the kernel [44].
In the case of JM-OCT, two Stokes vectors that correspond to two incident polarization states are acquired at each pixel. The two Stokes vectors are combined to compute a single DOPU (DOPU ) as where p represents averaging over two polarization states [43]. In our particular implementation, we use a further modified definition of the DOPU to correct for the effect of measurement noise. The details of this noise-corrected DOPU have been described elsewhere [43].

E. Intensity OCT, attenuation coefficient, and scatterer density
By combining the four OCT signals corresponding to the four Jones matrix entries, OCT intensity images are obtained as described in Section II. The attenuation coefficient image is obtained by further processing the OCT intensity. In our particular implementation, we use the method of Vermeer et al. [40].
The attenuation coefficient [61] is sensitive to tissue scattering and absorption and is suitable for highlighting tissue abnormalities. The attenuation coefficient is, however, affected by system and measurement parameters, such as the depth of focus, focus position, and aberrations. Several methods of correcting such effects have been demonstrated [62]- [64].
In our JM-OCT, we additionally use scatterer density imaging achieved by a neural-network-based scatterer density estimator [65]. This estimator is based on a convolutional neural network trained by using numerically generated OCT images. The images in the training dataset are generated from numerical samples with various scatterer densities. The OCT resolution also varies from image to image. The finally obtained estimator is thus insensitive to the system resolution and depth position of the focus. This scatterer density imaging is a good and robust alternative to attenuation coefficient imaging.

F. Doppler OCT and OCTA
In JM-OCT, the Doppler OCT (i.e., Doppler shift between two complex OCT signals acquired at two time points) is computed as whereĒ is a complex OCT signal obtained from the measured Jones matrix using the coherent composition method (Section 3.6 of Ref. [39]). r represents the spatial position, t is time, and ∆t is the time difference between two measurements. The superscript * denotes the complex conjugate. That is to say, E(r;t) andĒ(r;t + ∆t) are OCT signals obtained at the same location but at different times. They are typically B-scans obtained at the same location or A-scans obtained at very close positions. Further details of the JM-OCT-based Doppler OCT method are given in Section 3.7 of Ref. [39]. OCTA can also be computed from the Jones matrix. In principle, almost all OCTA algorithms designed for conventional OCT are applicable to JM-OCT. As the simplest example, we can compute the coherent composition of complex OCT signals or the intensity-averaged OCT signal from the Jones matrix. One of the conventional OCTA algorithms can be applied to these complex or intensity signals to yield an OCTA signal. The examples of these algorithms include intensity based methods, such as speckle variance method [66], phase based methods such as Doppler OCTA [11] and phase variance methods [67], [68], and complex methods such as optical micro-angiography [69]. By processing each entry independently, we can also use split-spectrum amplitude-decorrelation angiography method (SSADA) [70]. We can also more effectively use the four entries of each Jones matrix. In our particular implementation, we use a complex-correlation-based OCTA algorithm that is tailored to JM-OCT [47] and provides high-quality OCTA images. Notably, although it is not used in our JM-OCT implementation, Gong et al. recently presented a Jones-matrix-specific speckle decorrelation method for OCTA [71].

G. Dynamic OCT (D-OCT)
D-OCT, similar to OCTA, visualizes the temporal dynamics of the sample. In contrast to OCTA, D-OCT uses more OCT frames acquired with a longer time period and contrasts more subtle dynamics of the tissue, such as intracellular motility or metabolism. JM-OCT is compatible with several existing D-OCT methods.
In our JM-OCT modality, two D-OCT algorithms are used. One algorithm is LIV [24], which is defined as the variance of dB-scaled OCT intensities over a long acquisition time, typically around 6 s. Here, the OCT intensity is computed as the intensity average of four Jones matrix entries.
The other algorithm is OCDS [24], which is defined as the slope of the autocorrelation curve of the dB-scaled OCT signal. In our particular implementation, two kinds of OCDS are defined with different delay time ranges of the autocorrelation slope. One is early OCDS, which is defined for a short delay time, such as [12.8, 64 ms] [24], and is sensitive to a fast fluctuation of the signal, i.e., fast dynamics of the tissue. The other is late OCDS, which is defied for a long delay-time, such as [64, 627.2 ms] [24] or [204.8, 1228.8 ms] [26], and is sensitive to slow dynamics.
The LIV quantifies the magnitude of the signal fluctuation whereas the OCDS is sensitive to the speed of the fluctuation. They are thus complementary to each other.
Similar to the case of OCTA, JM-OCT is compatible with several other dynamic OCT algorithms, such as time spectral analysis [22], [23], [25], [27], correlation analysis [72], and Brownian-bridge analysis [73].  Fig. 3. Generalized depiction of JM-OCT hardware. JM-OCT hardware is characterized by two polarization components: the input-polarization multiplexer and a PD detector. The input-polarization multiplexer multiplexes two polarization components of the incident probe beam, which are demultiplexed after OCT detection . The PD detector measures two interference signals that correspond to two output polarizations and it thus gives the Jones vector of the backscattered probe beam.
IV. IMPLEMENTATION OF JM-OCT A. Generalized hardware As described in Section II, JM-OCT uses two input polarization states and Jones vector detection. Generalized hardware for such a detection scheme is depicted in Fig. 3. Here, two polarization-hardware components, namely the inputpolarization multiplexer and PD detector module, are introduced.
The input-polarization multiplexer multiplexes two incident polarization states by some means as it can be demultiplexed after interference-signal detection. For example, an electrooptic (EO) modulator can be used to modulate one of the polarization components with a particular frequency, so that two polarization states are multiplexed at different frequencies. In swept-source (SS)-OCT, the polarization components appear as two independent OCT images at two different depths [74], [75]. In spectral-domain (SD)-OCT, the EO modulator modulates the polarization states as a galvanometric mirror is performing a transversal scan (B-scan), so that the two polarization states are multiplexed into two different spatial frequencies along the transversal direction, and the states are demultiplexed later through numerical frequency filtering [76]. As an alternative option, the EO modulator switches the incident polarization states for each A-line. The two incident polarization components are then easily demultiplexed by deinterlacing the A-lines. More modern swept-source JM-OCT adopts passive-delay-based incident-polarization multiplexing [77]- [79]. Here, two different delays are applied to the two incident polarizations, and two OCT images corresponding to the two input polarizations thus appear at different depths in a single OCT image. The polarization-dependent delay can be applied using bulk optics comprising polarization beam splitters [39], [48], [77], [78] or by a long polarization maintaining fiber that has different group delays between the two polarization states [77], [80].
The second polarization unit is a PD detector, which is used not only for JM-OCT but also for several types of polarizationsensitive [32], [34], [76], [79], [81]- [88] and non-polarizationsensitive OCT [67], [89], [90]. The PD detector splits the probe and reference beam, or equivalently the interference signal, into two polarization components and detects two interference signals corresponding to the two orthogonal polarizations of the backscattered probe beam using two photodetectors. OCT can measure not only the amplitude but also the phase, and the PD detection scheme thus provides a pair of complex OCT signals corresponding to the Jones vector of the backscattered probe beam.
The input-polarization multiplexer multiplexes two incident polarizations, the PD detector measures two output polarizations, and four OCT signals (i.e., the interactions of two inputs and two outputs) are thus acquired. If both the inputpolarization pair and output-polarization pair are orthogonal pairs, the four OCT signals evidently correspond to the four entries of the measured Jones matrix. For example, if both the input-polarization and output-polarization pairs are the pairs of horizontal polarization ([1 0] T in the Jones vector representation) and vertical polarization ([0 1] T ), the measured Jones matrix is constructed as where Es are the complex OCT signals obtained with particular input and output polarization states. The first and second subscripts indicate the input and output polarization states, respectively, where H and V denote horizontal and vertical. Note that the complex OCT signals are the functions of the space, i.e., images. And hence, the measured Jones matrix is also an image with four complex channels (i.e., the four matrix entries). It is also noteworthy that the raw measured Jones matrix, J mi , in Eq. (1) is the J m of Eq. (11) obtained at the surface P i in Fig. 2.
Even if the input-and output-polarization pairs are not orthogonal, we can use the four OCT signals as the four entries of the measured Jones matrix without modification. The cumulative and local Jones matrix obtained from such measured Jones matrices is still similar to the matrix obtained with perfect orthogonality. Since a matrix and its similar matrix have identical eigenvalues, they give the same phase retardation values. The orthogonality is thus not a real requirement of JM-OCT. The diattenuation and polarizationdependent phase offset of the system also do not affect the final results. Even with these imperfections, the resulting local Jones matrix is again similar to the ideal one. The details of this strong self-calibration nature of JM-OCT are given in Section VI-C. Hz. The light is split by a 90/10 fiber coupler (C), and the 10% portion is sent to the reference delay line through a fiber Bragg grating (FBG, reflectivity > 80%). The FBG sharply reflects a narrow spectrum light at 1,266 nm, and the reflected light is detected by an amplified photodetector (APD). This triggers a digitizer for the acquisition of a spectral interference signal for each A-line. Here, the digitizer is driven by a k-clock signal generated by the light source. The probe beam passes a semi-bulk, which is a well-encased small-optics, passive polarization delay module (Optohub Co. Ltd., Saitama, Japan) comprising a linear polarizer (LP), a polarization beam splitter (PBS), and two right-angle prisms (RAPs). Here, the incident polarization is cleaned up in the 45-degree direction by the beam splitter and then split into pand s-polarizations by the PBS. Both polarized light beams are reflected by the RAPs, recombined by the PBS, and finally introduced into the probe optics through single-mode fibers and a circulator. By displacing one of the RAPs, we arbitrarily set the mutual delay between the p-and s-polarized light beams. Through this mutual polarization delay, two OCT images corresponding to these incident polarizations appear at two different depths in the imaging. That is to say, two incident polarizations are multiplexed by the mutual delay and demultiplexed into the two different depths.
The interference signal is detected by a semi-bulk encased PD detection module (Optohub). Here, the reference beam polarization is first cleaned up in the 45-degree direction. After the reference and probe beams are mixed by the nonpolarization beam splitter, the interference signal is split into horizontal and vertical polarization components by two PBSs and detected by two balanced photodetectors. The interference signals are then cleaned up by low-pass and high-pass radiofrequency filters and digitized by a digitizer synchronized with the wavelength sweeping of the light source. This PD diversity detection scheme simultaneously provides two interference signals, which correspond to two output polarizations. Accounting for the input-polarization multiplexing by the passive polarization delay module, four OCT images are obtained at once, and they form the measured Jones matrix. Multiple contrast images are obtained by processing this matrix using the method described in Section III.
More details of this particular implementation can be found in [84], [91]. Other variations but in a similar implementation of SS-JM-OCT can be found elsewhere [39], [48], [77], [78], [92]. It is also noteworthy that the passive polarization delay and the polarization diversity detection can be implemented not by using the bulk optics but using fiber components [77], [80].
2) Spectral-domain JM-OCT: JM-OCT can be implemented also with SD-OCT technology as exemplified in Fig. 5 (reprinted from Ref. [76]). Here, an EO modulator modulates the input polarization along the transversal scan (B-scan), so that two incident polarization states can be demultiplexed by numerical frequency filtering along the transversal scan. The PD detection is performed using a polarization-sensitive spectrometer in which two polarization components are split by a PBS and the interference signals corresponding to the two polarization components are detected by two line-scan cameras. The combination of input-polarization switching and PD spectrometer provides four OCT signals with different input-output polarization combinations that form the measured Jones matrix.
In some PS-OCT, the polarization multiplexing is done with an acoustic optical modulator [93]. The PD detection can also be achieved by several other configurations, such as a single-camera polarization sensitive spectrometer [81], [82] and combinations of two identical spectrometers [32], [34], [86].
Once the Jones matrix is measured, nearly identical signal processing with the SS-OCT-based JM-OCT can be conducted to yield images with multiple contrasts.

V. APPLICATIONS
JM-OCT has been applied to several clinical and biological fields, including investigation of the retina, anterior eye, skin, small animals, and cultured tissues. Some of these applications are reviewed in this section.

A. Retinal investigation
The retina is one of the first and most successful application targets of JM-OCT. For clinical retinal investigation, several clinical prototype devices [48], [87], including a functionlimited simplified version [85], have been developed.
In addition to the sole use of polarization-sensitive imaging, its combination with flow imaging has been found to be useful in the investigation of exudative macular diseases [46]. Figure  6 shows an example case of polypoidal choroidal vasculopathy (a new visualization of the case presented in Fig. 7 of Ref. [46]). The color fundus (a) exhibits a U-shaped hemorrhage (delineated by a dashed line). In the corresponding late-phase indocyanine green angiography (c), active leakage is found (arrow), and this leakage is collocated with the RPE defect found in the fundus auto-fluorescence (FAF) image (d, arrow). These clinical features can be found in the three-dimensionally reconstructed multi-functional JM-OCT (b, cut-away volume reconstruction). Here, gray represents the OCT intensity that corresponds to a conventional OCT signal. Green represents low-DOPU regions, which indicate retinal pigments. Yellow indicates flow acquired by Doppler OCT. A large detachment of RPE is seen at the center of the image. At the top of the RPE detachment (red arrow), a hyper-scattering layer corresponding to the RPE is shown, but this part of the RPE does not exhibit low-DOPU (green) signals. This indicates that retinal pigments no longer exist in this region, which corresponds to the RPE damage shown in the FAF image (d, arrow). A large cluster of flow (yellow arrow in b) is visible just beneath the RPE detachment. This may be the abnormal vessel corresponding to the leakage found in the indocyanine green angiography (ICGA) image (arrow in c). Next to the RPE detachment, a large hyper-scattering region without a flow signal is seen (blue arrow in b). This region corresponds to the hemorrhage that appears in the color fundus (a). No blood flow exists in this region as the hemorrhage is sufficiently old.
For further assessment of retinal disease, several image processing methods tailored to multi-contrast JM-OCT have been demonstrated. Azuma et al. arithmetically combined multiple contrasts to specifically highlight melanin associated with the RPE and choroid [95], which would be useful for the assessment of age-related macular degeneration [58], [60]. Kasaragod et al. demonstrated a machine-learning-based method of quantitatively assessing the tissue properties of the lamina cribrosa [100], which is expected to be useful in glaucoma diagnosis.
Notably, several retinal JM-OCT implementations exist. JM-OCT implementations were initially based on SD-OCT [94], but currently the majority are based on swept-source OCT with a 1-µm band probe beam and passive polarization multiplexing mechanism [39], [48], [77]- [79] (see also Section IV-A). In particular, the implementations of the authors started with a bulk laboratory prototype [39], [78], [94], and a compact clinical prototype based on the identical optical setup was demonstrated [48]. By introducing semi-bulk packages of the polarization multiplexer and PD detector, a more compact and stable clinical system has been demonstrated [87]. A simplified retinal JM-OCT, which uses a PD detector but does not multiplex the input polarization, has also been demonstrated [85]. Although this system cannot measure the phase retardation or birefringence, it can measure OCT, DOPU, and OCTA. This simplified system is expected to enable low-cost clinical translation. A motion-free and wide-field retinal scan method, namely the Lissajous scan method, has also been introduced into this system [101], [102]. This system might enable comprehensive retinal investigation with low device cost. Commercialization of retinal JM-OCT and its integration to a standard diagnostic flow are awaited.

B. Anterior eye
Another well-established target of JM-OCT is the anterior eye. An early investigation was performed using a bulky laboratory prototype [74], [103], and the first compact clinical prototype was then built by Lim et al. [75]. This work has been followed by several more sophisticated clinical prototypes [83].
The main targets of JM-OCT of the anterior eye are corneal diseases and glaucoma. Among corneal diseases, JM-OCT is believed to be useful for the evaluation of keratoconus [104]. JM-OCT is expected to be used in the assessment of corneal collagen cross-linking, which is a keratoconus treatment [105]. The glaucoma applications of JM-OCT can be classified as risk assessments relating to anterior angle imaging and the follow up of glaucoma surgery. For the former, the polarization imaging of JM-OCT has been found to be useful in identifying the trabecular meshwork [106], [107]. For the latter, JM-OCT is used in investigating the fibrotic change in post-surgical structures of trabeculectomy [103], [108]- [112]. In a JM-OCT based trabeculectomy study, Fukuda et al. pointed out that the phase retardation image can visualize slight fibrotic changes at early stage after the surgery, while the birefringence image can highlight localized fibrotic tissue [112]. In addition to these applications, the anterior sclera [113] and meibomian gland [114] have been investigated by JM-OCT.  [115]. Notably, their device is specialized for skin assessment because it is equipped with a dermoscopic device .

C. Skin imaging
Polarization, i.e., phase retardation and birefringence, investigations have revealed the age-related alternation of human skin [116]. The relation between the polarization properties and structures of the skin surface, such as the wrinkle structure, has also been investigated [117], [118]. The anisotropic response of dermal birefringence to mechanical stretching and compression has been found adopting JM-OCT [119]. Collagen is not only birefringent but also a strong source of second harmonic generation, and the combination of JM-OCT and second harmonic generation imaging is thus promising for the comprehensive investigation of skin [120].
In addition, the DOPU investigation [121], and ultra-highresolution birefringence and optic axis imaging [122] of skin have recently been demonstrated.

D. Zebrafish imaging
Zebrafish are suitable for the close emulation of human diseases. Additionally, zebrafish have a short life span, providing a faster experimental cycle than other experimental animals. JM-OCT has been adopted for the noninvasive and comprehensive investigation of zebrafish.
Lichtenegger et al. recently demonstrated several applications of JM-OCT to zebrafish. The investigation of postmortem  zebrafish revealed the fine structure of gills and eye globe [123]. The birefringence and depolarization features of the eye globe were clearly visualized as shown in Fig. 7 (reprinted from [123]). Furthermore, by stitching multiple (three) JM-OCT volumes, the whole structure and birefringence of a zebrafish have been visualized as shown in Fig. 8 (reprinted from [123]).
Juvenile and young zebrafish have been investigated in vivo through multi-contrast imaging, including the OCT intensity, birefringence, DOPU, and OCTA [124]. The skeletal muscles and swim bladder were found to have high birefringence. The RPE/choroid complex has a low DOPU, which may indicate melanin. The quantitative analysis of birefringence revealed an increase in spinal birefringence with the growth of the fish (from 1 month to 2 months of age).
Human tumor cells can be xenotransplanted into zebrafish, and the zebrafish is thus expected to be an important sample for drug development. The time evolution of human-derived tumor cells xenografted on zebrafish has also been shown to possibly be quantitatively evaluated through JM-OCT [125]. In the cited study, the temporal changes in the scattering, birefringence, and DOPU of a xenografted zebrafish were quantitatively investigated. The OCTA images obtained by JM-OCT were also shown to be used to quantitatively evaluate the temporal changes in the vasculature structure of the zebrafish.

E. Ex vivo and in vitro tissue imaging
In vitro and ex vivo samples are recently emerging targets of JM-OCT because they are increasingly used for drug testing and development. The temporal change and drug effect of tumor spheroids have been investigated [24], [26]. Here, the D-OCT function of JM-OCT was used to visualize the cellular activities as exemplified in Fig. 9 (adopted from Ref. [26]). The temporal degradation of a human breast cancer (MFC-7) spheroid after extraction from a culturing environment was longitudinally and three-dimensionally visualized using the late OCDS (see Section III-G for the details of the method).
In addition to the D-OCT analysis, Seesan et al. recently applied neural-network-based scatterer-density analysis to the JM-OCT signal [65]. This revealed a temporal reduction in the scatterer density after the spheroid was extracted from the culturing environment as shown in Fig. 10. The multi-contrast capability of JM-OCT might enable the comprehensive investigation of in vitro tissues.
Ex vivo mouse tissues have also been investigated adopting JM-OCT. The structure and activities of the microvascular Fig. 11. Multi-contrast imaging of mouse liver (a-f). The LIV image revealed fast dynamics (green) beneath the vessel stutterer even in the ex vivo sample. The LIV, attenuation coefficient (AC), and birefringence signals are converted into three principal component images (PC1 to PC3) in principal component analysis (g-l). By applying a tri-variate colormap (m-p), a blood vessel (purple) and active tissue beneath the vessel are clearly visualized. The images were adopted from [126]. complex of mice have been visualized by combining simultaneously acquired attenuation coefficient, LIV, and birefringence images as shown in Fig. 11 [126]. The D-OCT (LIV) investigation of normal and non-alcoholic-fatty-liver-disease (NAFLD) model mouse livers has revealed the high activity of lipid droplets, zonation structures, and inflammation that were not visible in conventional OCT imaging [127]. Here, the NAFLD models were 1-week and 2-week methionine cholinedeficient (MCD)-diet models, which were fed with a MCD diet for 1 and 2 weeks. Figure 12(a, b) shows the three-dimensional distribution of lipid droplets of the 1-week NAFLD model.
Simultaneously obtained birefringence images show that the birefringence of the liver of the 1-week NAFLD model is higher than that of the normal liver, and it becomes low again in the liver of 2-week NAFLD model mouse as shown in Fig.  12 (c-d) [127].

F. Further applications
JM-OCT has been adopted beyond the applications described above. Tang and her colleagues have intensively investigated the polarization properties of cartilage [92], [128]. Matsuzaki et al. adopted JM-OCT in assessing human-induced pluripotent-stem-cell-derived retinal pigment epithelial sheets using the polarization-sensitive function of JM-OCT [129]. Here, the measurement was made both in vitro and in vivo. These applications suggest that JM-OCT can be used in regenerative medical procedures before and after transplantation.

VI. ADDITIONAL TECHNOLOGICAL TOPICS
A. Signal to noise ratio of JM-OCT The system sensitivity of OCT is known to be proportional to the probe beam power. Since the JM-OCT splits the probe power into four polarization channels, the sensitivity of each channel is four-times lower in average than the standard OCT with the same probe power. This reduction of sensitivity can be compensated by coherently composite the four polarization channels (see Section 3.6 of Ref. [39]).
It is also noteworthy that the signal-to-noise ratio (SNR) of the cumulative phase retardation and birefringence measurements is dominated by a quantity called "effective SNR," which is a reciprocal mean of SNRs of OCT signals involved in the computation of the cumulative phase retardation or birefringence (see Eq. (39) of Ref. [42] for the mathematical definition). In typical cases, the effective SNR is maximized when the probe power is equally distributed into the four polarization channels. This gives us a simple optimization strategy of JM-OCT, i.e., as using a non-polarized mirror as a sample, adjusting the polarization controllers and polarizers in the system to have roughly equal signal strengths among the four polarization channels.

B. Signal processing methods for better accuracy and image quality
Several signal processing techniques have been demonstrated to enhance the measurement accuracy and image quality of JM-OCT. We listed examples of these technologies in this section.
To reduce the negative effect of noise, several signal enhancement and noise correction methods have been developed. The first is adaptive Jones matrix averaging (AJA). AJA is essentially complex averaging of the Jones matrix with a global phase correction. The Jones matrix is a complex matrix, and the averaging of multiple Jones matrices thus suppresses the complex noise and fundamentally enhances the signal-to-noise ratio. However, if the global phases fluctuate among the Jones matrices, the complex averaging decreases not only the noise but also the signal. In AJA, the difference in the global phase between Jones matrices is estimated and corrected before the averaging. AJA thus optimally enhances the SNR of the Jones matrix. The details of the AJA method are described in Section 2.1.2 of Ref. [75] and Section 3.4 of Ref. [39].
If the noise is additive and has a zero mean, such as in the case of complex noise in a complex OCT signal, the mean (averaging) of measured signals asymptotically approaches the true value. However, the measured phase retardation and birefringence under noise are not symmetrically distributed around the true value [42], [130]. This is because of the nonlinear signal processing used to compute the phase retardation and birefringence from the set of complex OCT signals, i.e., the Jones matrix. To accurately estimate the true birefringence from multiple measured birefringence values, maximum aposteriori estimators have been demonstrated [131], [132]. Here, Monte-Carlo simulation of JM-OCT imaging was used to compute the probability distribution of the measured birefringence under a particular true birefringence value and a particular effective SNR value. Using a vast amount of simulation results with several configurations of the true birefringence and effective SNR, we can inversely determine the a-posteriori distribution, which is the probability distribution of the true birefringence under particular values of the measured birefringence, and an effective SNR is acquired in measurements . Using this estimator, a bias-free estimation of the birefringence can be obtained from the raw JM-OCT measurements. The birefringence images shown in Figs. 7, 8, 11, and 12 were obtained using this maximum a-posteriori estimator. A similar concept of maximum a-posteriori estimation [133] has been adopted to enhance the accuracy of the polarization-insensitive OCT intensity of JM-OCT [134].
It should be noted that, for further improvement of the accuracy of phase retardation and birefringence measurement, it might be important to validate and optimize the measurement and estimation methods by using calibrated and standardized birefringent phantoms [34], [38], [135], [136].
The DOPU is also highly affected by noise. The DOPU is the polarization uniformity among multiple measurements, and it is thus low if the SNR of the raw OCT signal is low. Makita et al. presented a modified definition of the DOPU, which avoids the measurement-noise-originating reduction of the DOPU. The details of this method are described in Ref. [43].
In our JM-OCT, the OCTA signal is defined as temporal complex decorrelation of OCT signals, and it is thus also strongly affected by measurement noise. To obtain highcontrast OCTA from the measured Jones matrices, a noiseimmune OCTA method that estimates the decorrelation in correcting the phase noise effect, has been demonstrated [47].

C. Self-calibration nature of JM-OCT
The JM-OCT theory described in Sections III-A and III-B assumes that both the multiplexed incident-polarization pair and the detected polarization pair are perfectly orthogonal pairs. However, in reality, it is difficult to achieve such ideal orthogonality. For example, a small amount of cross-talk in the PBS used in the passive polarization delay module and the PD detector degrades the orthogonality. Alignment error of the PD detector causes a mutual phase shift between the two OCT signals of the two output polarizations, and it disturbs the ideal measurement of the Jones vector. Although these imperfections are unavoidable in practical implementations, they are automatically canceled in the computation of cumulative and/or local Jones matrices. This self-calibration nature is further described as follows.
Considering the above factors, the measured Jones matrix at position P i [Eq. (1)] is modified as where E in is a matrix whose columns are two polarizations incident on the polarization delay module. X is a matrix representing the imperfection of the PBS in the polarization delay module (input-polarization multiplexer), and f ( ) is a function representing the delay between two incident polarization states introduced by the delay module. The off-diagonal entries of X account for the cross-talk and the diagonal entries represent the transmittance and reflectance of the horizontally and vertically polarized light. X appears twice in the equation, where the first (right) X represents the PBS splitting the beam into two polarization states whereas the second (right) represents the combining of the two beams. R represents the interference with the reference beam as R = H * re f 0; 0 V * re f , where H * re f and V * re f are the complex conjugates of the horizontal and vertical components of the reference light. ρ is a rotation matrix accounting for the relative rotation between the polarization-delay and PD-detection modules, and X represents the imperfection of the PBS in the PD-detection module. η is [η A 0; 0 η B ] as η A and η B are the detection efficiencies of two detectors used for PD detection. That is to say, X, ρ, R, X , and η represent the imperfections.
The cumulative Jones matrix with these imperfections can be obtained by substituting Eq. (12) into Eq. (3) as Evidently, this matrix is similar to the original cumulative Jones matrix J out J r0i J −1 out shown in Eq. (3). A matrix and its similar matrix have the same eigenvalues, and the imperfections thus do not affect the phase retardation measurement. Additionally, because similar self-calibration occurs for the local Jones matrix computation, the imperfections also do not affect the local phase retardation computation. That is to say, without any special signal processing, the imperfections are automatically canceled in the derivations of the cumulative and local Jones matrices.
Although here we have discussed an example with a passivedelay based input-polarization multiplexing, the conclusion can be generalized for all types of input-polarization multiplexing. Further details of this self-calibration are described in Section 5.4 of Ref. [39] and Sections 35.5.2 and 35.5.3 of Ref. [29].

D. Resolution of birefringence imaging
It might be fair to note that the depth resolution of birefringence imaging is lower than that of OCT. As described in Sections III-A and III-B, the birefringence is computed from a local Jones matrix, and the local Jones matrix is computed from two cumulative Jones matrices at slightly distant depths. And hence, the practical depth resolution of birefringence imaging can be roughly estimated by adding the optical axial resolution and the depth-separation of the two cumulative Jones matrices.
The separation of the cumulative Jones matrices is selected from the following three perspectives. At first, the separation should be close enough as the tissue between the two depths can be considered as homogeneous, as discussed in Section III-B of this manuscript and Sections 2, 6.1, 6.3 and Appendix A of Ref. [48]. Second, the separation must be larger than the optical resolution of the OCT, otherwise the two Jones matrices are highly correlated to each other and meaningful birefringence value cannot be obtained. Third, the separation should be large enough to make the local phase retardation larger than the measurement noise of the phase retardation.
Typically, the separation is around 20 to 50 µm , while the depth resolution is a-few to ten-few micrometers. And hence, the depth-resolution of the birefringence measurement is a-few-times lower than the depth resolution of OCT. In our particular implementations, the optical resolution and the depth separation are, respectively, 6.2 µm and 24 µm (in tissue, for a retinal JM-OCT with a 1-µm probe wavelength) [132] or 14.1 µm and 57.9 µm (in tissue, for a JM-OCT with 1.3µm probe wavelength) [91]. And hence, the resolutions of the birefringence measurement of these systems are around 30 µm and 72 µm , which are around five-times lower than the optical resolution.

E. Polarization-sensitive OCT and JM-OCT
The JM-OCT described in this paper is a multi-contrast imaging modality based on OCT. However, historically, JM-OCT was a variation of PS-OCT. It is thus worth summarizing the relationship between the several PS-OCT methods and JM-OCT.
The Jones matrix analysis of PS-OCT has been demonstrated by multiple research groups including Park et al. [41]. Some PS-OCT methods measure the Jones matrix or equivalent quantities but adopt other types of polarization analysis, such as Jones-Muller analysis [35]- [37] and Stokes analysis [137].
A well-known PS-OCT method is the circular-polarization method, which illuminates the sample with a single state of circular polarization and measures two OCT signals using a PD detector [4], [57]. For this type of PS-OCT, the phase retardation and optic-axis orientation can be computed directly from the two detected OCT images. Although it is not common, this type of PS-OCT hardware can be used for JM-OCT imaging. That is to say, by assuming the unitarity of the sample, we can reconstruct a unitary Jones matrix from the two measured OCT signals, and then, we can apply the Jones matrix analyses described in Section III. Additionally, open issues, which might be addressed by further development, are discussed.

A. Combination with mechanical-property measurement
The relationship between birefringence and a mechanical property of a tissue is an intriguing research topic because both can be related to collagen. By measuring ex vivo porcine samples adopting JM-OCT and a uniaxial tension tester, Yamanari et al. and Nagase et al. showed that the polarization and mechanical properties of sclera are interrelated [138], [139]. Intraocular pressure and scleral birefringence have also been found to be correlated in human eyes in vivo [113].
One existing but not well-established JM-OCT extension is JM-OCT elastography [91]. By introducing a tissuecompression probe into JM-OCT, the simultaneous imaging of JM-OCT contrasts and tissue deformation was demonstrated. By combining multiple contrasts, including tissue deformation, birefringence, and attenuation coefficients, five tissue layers of an ex vivo porcine esophagus were revealed. Notably, these layers were not distinguishable if only one or two of the three contrasts were used.
Recently developed tissue-compression-based optical coherence elastography enables a more quantitative assessment of tissue mechanical properties (e.g., see Section 6.5 of Ref. [15]). It might be worth introducing such a quantitative optical coherence elastography method into JM-OCT.

B. JM-OCT specific image processing
Image processing methods specialized for multi-contrast JM-OCT have been demonstrated. An early demonstration was the cumulative-phase-retardation-based segmentation of the chorio-scleral interface of the human retina [140]. This method mainly uses a single contrast of the cumulative phase retardation and is thus similar to a PS-OCT-based approach demonstrated by Torzicky [141].
Other examples have more intensively used multiple contrasts obtained through JM-OCT. Myiazawa et al. first demonstrated tissue segmentation of the anterior eye using three contrasts, namely the OCT intensity, attenuation coefficient, and depth alteration ratio of cumulative phase retardation [142]. All the pixels of a JM-OCT image of the anterior eye were mapped in the three-dimensional feature space of these three contrasts, and the sclera, conjunctiva, uvea, and trabecular meshwork were segmented. Additionally, Azuma et al. [95] and Kasaragod et al. [100] demonstrated an arithmetic computation for RPE and choroidal segmentation and machine-learning-based methods for lamina cribrosa segmentation, respectively, as mentioned in Section V-A.
Miyazawa et al. demonstrated pixel grouping for JM-OCT, superpixelization [143]. In this method, the pixels in a crosssectional JM-OCT image are mapped into a six-dimensional feature space of four JM-OCT contrasts and two (lateral and axial) spatial dimensions. The four contrasts are the OCT intensity, OCTA, birefringence, and DOPU. The pixels are clustered in this six-dimensional space such that each cluster comprises similar pixels in not only its JM-OCT contrasts (i.e., optical properties) but also its spatial proximity. Each cluster forms a group of pixels, the superpixel, which, on the one hand, follows the morphological structure of the tissues and, on the other hand, establishes a tissue region with similar optical properties. The superpixels can be used as kernels for contrast averaging and potentially as spatial kernels of estimators of optical properties such as birefringence [100], [131], [144].

C. Other birefringence calculation methods
As described in Sections III-B and VI-B, the birefringence is obtained through the local Jones matrix analysis [42] and a maximum a-posteriori estimator [132] in our particular implementation. On the other hand, there are some other methods to compute the birefringence (and, similarly, the local phase retardation) that are compatible with JM-OCT hardware.
Villiger et al. used a JM-OCT hardware and a differential Mueller matrix method to compute the local phase retardation [36]. Fan et al. combined a single-input circular-polarization PS-OCT and a Jones analysis to obtain the local phase retardation [145]. In addition, Stokes-based methods such as Refs. [34], [146] also can be implemented on JM-OCT hardware with minor modifications.

D. Optic-axis-orientation measurement
Lu et al. demonstrated the measurement of the orientation of the optic axis using a single-mode-fiber-based JM-OCT system [147]. Although this method successfully determined the absolute optic-axis orientation of a sample, it requires an on-site calibration. For this calibration, two quarter-wave plates should be imaged with the sample. This will elaborate the system or sample configuration and sacrifice the depth imaging range, and the system is thus not fully portable. Recently, Li et al. demonstrated Jones matrix symmetrization method for optics axis measurement [38], [122], [148].
In addition, several sophisticated methods of measuring the axis orientation have been demonstrated adopting other types of PS-OCT or JM-OCT in its widest sense. Fan et al. demonstrated Jones-matrix-based axis-orientation analysis [149] based on a single-polarization-input bulk PS-OCT [150], which has been adopted for PS-OCT-based tractography [151]. Additionally, the measurement of the axis orientation based on Jones-Muller analysis has been demonstrated [37]. Although this method is based on Muller matrix analysis, the raw data were obtained through Jones-matrix tomography and then converted into the Muller matrix.
A Stokes-based method for axis-orientation measurement has been demonstrated [34], [152]. This method can also be applied to JM-OCT data by converting to the Stokes representation.

E. Other polarization randomness metrics
Although we considered only the DOPU in this paper, there are several other polarization randomness metrics that are compatible with the JM-OCT hardware. Yamanari et al. used the spatial entropy of a set of Jones matrices obtained in a small local region as a depolarization measure [144], and this approach has been adopted in investigating the RPE [129]. Whereas the DOPU represents the randomness of the Jones vector, i.e., the backscattered probe beam, this entropy represents the randomness of the Jones matrices. It is thus a more direct measure of the sample's polarization randomness. Lippok et al. used degree of polarization (DOP) to represent the polarization randomness [35]. They measured the Jones matrix distribution of the sample, converted it into a Jones-Muller matrix, and computed the DOP from the Jones-Muller matrix. The entropy and DOP were both computed from the Jones matrix data sets, and they are thus compatible also with the narrow-sense JM-OCT. In addition, the first demonstration of DOPU by Götzinger et al. in 2008 [44] is preceded by two similar works. In 2007, Adie et al. used spatial variations of Stokes vectors to detect multiple scattering [153]. In 2008, Lee et al. used a quantity similar to DOPU to assess the depthdependent alteration of scattering [154]. Both of the methods are applicable to JM-OCT.
The randomness of optic-axis orientation is also useful in highlighting tissue abnormality. Willemse et al. defined a metric called optic-axis uniformity for the segmentation of a retinal nerve fiber layer [155]. This method is based on Jones-Muller analysis and thus can be adapted to JM-OCT.

F. Other open issues
Although we independently demonstrated D-OCT and PS-OCT, we did not assess the dynamics of the polarization properties of tissues. If we can establish a method of visualizing the temporal dynamics of the polarization properties, as preliminary demonstrated by Mukherjee et al. [156], we might reveal additional functional properties of the tissue.
Although we discussed noise correction methods in Section VI-B, other factors may disturb the quantitative JM-OCT measurement. One example is the defocus, which has been found to degrade the accuracy of birefringence and DOPU measurements and also generates polarization artifacts [157]. This effect can be mitigated by computational refocusing.
Computational adaptive optics [158], [159] and defocus correction [160], [161] may further improve the image quality of OCT and JM-OCT. However, the aberration can be polarization dependent [162]. The combination of JM-OCT and computational aberration correction and refocusing methods may enable polarization-dependent aberration correction.

VIII. CONCLUSION
JM-OCT is a multi-contrast OCT modality based on a conception, in which all contrasts are derived from a single quantity, i.e., a time-sequence of Jones-matrix tomographies. Historically, JM-OCT began as a variation of PS-OCT, but it later acquired a multi-contrast capability. JM-OCT imaging not only relies on optical instrumentation but also exploits several mathematical and computational techniques, such as sophisticated contrast estimators. JM-OCT has been applied to several medical and biological fields, with posterior and anterior eye imaging being two of the most successful applications. JM-OCT has recently been applied to in vitro and ex vivo samples and in vivo small animals. This variation of JM-OCT, JM-OCT microscopy, is expected to be a volumetric, label-free, and quantitative microscopy. It will thus be a useful tool in biology and drug development. Several technical issues still need to be overcome, such as the accurate optic-axisorientation measurement and polarization-dependent computational aberration correction. These on-going works will further enhance the capability of JM-OCT.

ACKNOWLEDGMENTS
The JM-OCT-related research works are based on the significant contributions of previous and current laboratory members of the Computational Optics Group at the University of Tsukuba. In particular, seminal contributions from Masahiro Yamanari and Shuichi Makita are greatly acknowledged.
This review also covers some research works of non-Jones matrix PS-OCT and other multi-contrast OCT imaging technologies. It should be acknowledged that all OCT-community members have contributed heavily to establish such multicontrast OCT technologies.
The JM-OCT projects conducted at University Tsukuba have been supported by the Japan Science and Technology Agency through contracts with CREST (JPMJCR2105), MI-RAI (JPMJMI18G8), and the Development of Advanced Measurement and Analysis Systems, and by the Japan Society for the Promotion of Science (21H01836, 22K04962, 21K09684, 18K09460, 18H01893, 17H04350, 15K13371, 15K10905, 25293353, 24592682, 22390320, 21659397).