# Atmospheric Sulfur Photochemistry on Hot Jupiters

K. Zahnle  NASA Ames Research Center, Moffett Field, CA 94035
<EMAIL_ADDRESS>M. S. Marley  NASA Ames Research Center, Moffett
Field, CA 94035<EMAIL_ADDRESS>R. S. Freedman  NASA Ames Research
Center, Moffett Field, CA 94035<EMAIL_ADDRESS>K. Lodders
Washington University - St. Louis  J. J. Fortney  Department of Astronomy and
Astrophysics, University of California - Santa Cruz

###### Abstract

We develop a new 1D photochemical kinetics code to address stratospheric
chemistry and stratospheric heating in hot Jupiters. Here we address optically
active S-containing species and CO2 at $1200\leq T\leq 2000$ K. HS (mercapto)
and S2 are highly reactive species that are generated photochemically and
thermochemically from H2S with peak abundances between 1-10 mbar. S2 absorbs
UV between 240 and 340 nm and is optically thick for metallicities $[{\rm
S}/{\rm H}]>0$ at $T\geq 1200$ K. HS is probably more important than S2, as it
is generally more abundant than S2 under hot Jupiter conditions and it absorbs
at somewhat redder wavelengths. We use molecular theory to compute an HS
absorption spectrum from sparse available data and find that HS should absorb
strongly between 300 and 460 nm, with absorption at the longer wavelengths
being temperature sensitive. When the two absorbers are combined, radiative
heating (per kg of gas) peaks at 100 $\mu$bars, with a total stratospheric
heating of $\sim\\!8\times 10^{4}$ W/m2 for a jovian planet orbiting a solar-
twin at 0.032 AU. Total heating is insensitive to metallicity. The CO2 mixing
ratio is a well-behaved quadratic function of metallicity, ranging from
$1.6\times 10^{-8}$ to $1.6\times 10^{-4}$ for $-0.3<[{\rm M}/{\rm H}]<1.7$.
CO2 is insensitive to insolation, vertical mixing, temperature
($1200<T<2000$), and gravity. The photochemical calculations confirm that CO2
should prove a useful probe of planetary metallicity.

planetary systems — stars: individual(HD 209458, HD 149026)

††slugcomment: submitted to Ap. J. Lett.

##  1 Introduction

Stratospheric temperature inversions are ubiquitous in the Solar System, and
it is beginning to look as if they are commonplace on hot Jupiters as well.
Stratospheric temperature inversions form when substantial amounts of light
are absorbed at low pressures (high altitudes) where radiative cooling is
inefficient. Hubeny et al. (2003) pointed out that efficient absorption of
visible light by gaseous TiO and VO would greatly heat the upper atmospheres
of those planets already hot enough for these molecules to be present as
vapor. Thermal inversions on transiting hot Jupiters were first seen by
Richardson et al. (2007) for HD 209458b and Harrington et al. (2007) for HD
149026b. The observed flux ratio at 8 $\mu$m for HD 149026b agreed only with
models that included a thermal inversion (Fortney et al. 2006). Temperature
inversions have since been confirmed by _Spitzer_ observations of HD 209458b
(Knutson et al. 2008a), XO-1b (Machalek et al. 2008), and TrES-4 (Knutson et
al. 2009), all of which show distinctive flux ratios in IRAC bands that
suggest inversions (Fortney et al. 2006; Burrows et al. 2007). More
circumstantial evidence exists for HD 179949b (Barnes et al. 2008).

On the other hand TrES-1, the least irradiated planet with published _Spitzer_
observations, does not appear to have a pronounced inversion (Burrows et al.
2008). Nor, seemingly, does HD 189733b, which is also modestly irradated
(Charbonneau et al. 2008; Barman et al. 2008). One suggestion is that
temperature inversions are triggered by irradiation reaching a critical level
that is hot enough to evaporate TiO and VO from grains, as discussed by
Burrows et al. (2007), Fortney et al. (2008), and Burrows et al. (2008).
However, irradiation of XO-1b and HD 189733b is within uncertainties the same
(Torres et al. 2008), which poses a challenge to the irradiation trigger.

In the Solar System, stratospheric temperature inversions are often caused by
absorption of UV light by gases or aerosols produced by photochemistry. Here
we ask if atmospheric chemistry might play a similar role in hot Jupiters.
Speculation has tended to focus on sulfur-containing species (Tinetti 2008),
as the reservoir species H2S is expected to be abundant (Visscher et al 2006)
in these atmospheres and many of its breakdown products (S2, in particular)
absorb violet and ultraviolet light.

##  2 The Photochemical Model

Previous photochemical modeling of hot Jupiters addressed the abundance of
photochemical H (Liang et al 2003) and the absence of photochemical smogs
(Liang et al 2004). Liang et al (2003) focused on the high H/H2 ratio that
arises from H2O photolysis. In their second paper, Liang et al (2004) argued
that simple hydrocarbons would not condense to form photochemical smogs in hot
solar composition atmospheres. Neither study considered sulfur.

We have developed a new general purpose 1D photochemical kinetics code
applicable to hot extrasolar planets. The code is based on the sulfur
photochemistry model for early Earth originally described by Kasting et al
(1989) and Kasting (1990), and subsequently adapted by Zahnle et al (2006) and
Claire et al (2006) to address sulfur photochemistry of Earth’s atmosphere
during the Archean, and by Zahnle et al. (2008) to address martian atmospheric
chemistry. Steady state solutions are found by integrating the system through
time using a fully implicit backward-difference method.

Our chemical network has been upgraded from that used by Zahnle et al (1995)
to address the chemistry generated when the fragments of Comet Shoemaker Levy
9 struck Jupiter. We have assembled a reasonably complete list of the
reactions that can take place between the small molecules and free radicals
that can be made from H, C, O, N, and S. The code solves 507 chemical
reactions for 49 chemical species: H, H2O, OH, O, O2, CO, CO2, HCO, H2CO, C,
CH, CH2, CH3, CH4, CH3O, C2, C2H, C2H2, C2H3, C2H4, C 2H5, C2H6, C4H, C4H2,
CN, HCN, N, N2, NO, NH, NH2, NH3, NS, H2S, HS, S, S2, S3, S4, S8, SO, HSO,
SO2, OCS, CS, HCS, H2CS, CS2, and H2. Reaction rates, when known, are selected
from the publicly available NIST database (http://kinetics.nist.gov/kinetics).
In order of decreasing priority, we choose between reported reaction rates
according to relevant temperature range, newest review, newest experiment, and
newest theory. Reverse reaction rates $k_{r}=K_{\rm eq}k_{f}$ of two-body
reactions are determined from the forward reaction rate $k_{f}$ and the
equilibrium $K_{\rm eq}=\exp{\left\\{\left(-\Delta H+T\Delta
S\right)/RT\right\\}}$ by using $H^{\circ}(T)$ and $S^{\circ}(T)$ as available
($R$ is the gas constant). Rates are not available for all reactions,
especially for reactions involving elemental sulfur. We will present a full
listing of the chemical reactions important to sulfur in a more general
followup study.

Here we use simple descriptions of atmospheric properties. The background
atmosphere is 84% H2 and 16% He. We include Rayleigh scattering by H2
(Dalgarno and Williams 1962). For our base case we assume an isothermal
atmosphere with T=1400 K; constant vertical eddy diffusivity $K_{zz}=1\times
10^{7}$ cm2/s; a surface gravity of 20 m/s2; and insolation levels $I$ by a
solar twin that are 1000$\times$ greater than at Earth. Metallicity proved to
be the most interesting parameter and was varied $-0.3\leq[{\rm M}/{\rm
H}]\leq 1.7$. In these units, solar metallicity is $[{\rm M}/{\rm H}]=0$,
Jupiter’s is $[{\rm M}/{\rm H}]=0.5$, and Saturn’s is $[{\rm M}/{\rm H}]=0.8$.
Short chemical lifetimes of S-containing species make our results insensitive
to $K_{zz}$. Model parameters are listed in Table 1.

At the upper boundary we set a zero flux lid at 1 $\mu$bar, with neither
escape nor exogenous supply. For the lower boundary we use fixed equilibrium
mixing ratios of the most abundant species at 1 bar of H2 and temperature $T$
(Lodders and Fegley 2002, Visscher et al 2006). For other species we force the
mixing ratio at 1 bar to approach zero. We scale the lower boundary conditions
such that the total mixing ratios of C, O, N, and S all scale linearly with
metallicity.

Absorption by S2 between 240 nm and $\sim$360 nm from the ground state is
analogous to the Schumann-Runge system in O2 (Okabe 1978). Strong, distinctive
S2 emission near 300 nm was observed on Jupiter after the impact with
Shoemaker-Levy 9 with Jupiter in 1994 (Noll et al 1995). Subsequent
thermochemical modeling showed that S2 readily forms as a major product in a
shock-heated ($T>1000$ K) gas of either cometary or jovian composition (Zahnle
et al 1995, Zahnle 1996). S2 has also been seen in gases vented by volcanoes
on Io (Spencer et al, 2000; Moses et al 2002). For S2, we use absorption cross
sections at 1500 K computed by van der Heijden and van der Mullen (2001).

The HS (mercapto) radical absorbs from its ground state at 324 nm (Okabe
1978). Visscher et al (2006) predicted that HS would be very abundant in
equilibrium at hot Jupiter conditions. We find the same. We therefore
calculated absorption cross sections of HS at four temperatures at 30 mbar
pressure using literature values of the molecular properties. The ground
$X^{2}\Pi$ state has been well studied (Ram et al 1995) but the upper level
$A^{2}\Sigma^{+}$ is subject to strong pre-dissociation (Resende and Ornellas
2001, Wheeler et al 1997, Schneider et al 1990, Henneker and Popkie 1971), and
only the value of the rotational constant $B$ and the spacing of the lowest
vibrational energy levels have been well measured. Using these constants and a
value for the electronic band oscillator strength of the 0-0 transition
derived from a study of HS in the solar spectra by Berdyugina and Livingston
(2002), a line list was computed using the RLS code developed by R.N. Zare and
D. Albritton (Zare et al 1973). This RLS code uses the molecular constants and
band strengths to predict line positions and strengths by fitting to an RKR
potential (Zare et al 1973). Other needed data — Franck-Condon factors,
partition functions, etc. — were derived either from the cited literature, the
program itself, or from Sauval and Tatum (1984) or Larsson (1983). The
calculations were carried out for values of $v^{\prime\prime}$(0-4) and
$v^{\prime}(012)$. Because the excited vibrational levels of the
$A^{2}\Sigma^{+}$ state are unstable with respect to predissociation, the
corresponding optical transitions are likely to be broad and shallow, or even
continuous. These uncertainties principally affect the absorption spectrum at
wavelengths shorter than 324 nm, which is in the range that is absorbed
strongly by S2. Results are shown in Figure 1. In the photochemical model we
used only the 1500 K absorption coefficients.

Other sulfur allotropes are better absorbers than S2 but less abundant. S3
absorbs strongly between 350 and 500 nm, and S4 absorbs between 450 nm and 600
nm, but more weakly (Billmers and Smith 1991). Unfortunately, the chemistries
of S3 and S4 are very uncertain, and we have had to estimate the important
reaction rates. In an earlier version of this study, we focused on the heats
of formation, and we tentatively concluded that S3 heating would be important
for metallicities $[{\rm S}/{\rm H}]>0.7$. We have since learned that
reactions of the form H + Sn $\rightarrow$ HS + Sn-1, where $n\geq 2$, are
strongly favored by entropy. The revised model predicts less S2 and much less
S3, which reduces the importance of S3 heating considerably.

Sulfanes (H2Sn, hydropolysulfides) will be present in cooler hot Jupiters. At
low temperatures sulfanes absorb VUV between 260 nm and 330 nm (Steudel and
Eckert 2003). Absorption may extend beyond 400 nm at higher temperatures as
the ground state becomes vibrationally excited, as in HS, but to first
approximation these wavelengths are covered by the more abundant S2 and HS. We
have not included sulfanes in this study.

##  3 Results

Figure 2 shows how CO2 and the abundant S-containing species vary as a
function of altitude. This particular case shows a hot Jupiter at 1400 K with
a “planetary” metallicity of $[{\rm M}/{\rm H}]=0.7$. Figure 2 is broadly
representative of all our models with $1200\leq T\leq 2000$ K and $-0.3<[{\rm
M}/{\rm H}]<1.7$. In particular, S2 and HS show well-defined peaks at
$\sim\\!2$ mbars that coincide with the altitude where H2S photolysis becomes
important. At lower altitudes H2S is the main S-containing species, and at
higher altitudes S is. It is also notable that the atmosphere becomes more
oxidizing at higher altitudes where H2O photolysis is important.

Table 1 lists some key results pertinent to sulfur for several variations of
basic model parameters. The models assume that $K_{zz}=10^{7}$ cm2/s and
$g=2000$ cm/s2 unless otherwise noted. In this temperature range the models
are insensitive to $K_{zz}$ (results not shown). Model G shows that, as
expected, column densities vary inversely with $g$.

Column densities of S2 and HS are sensitive to metallicity. To first
approximation, species with one metal atom, such as H2O and H2S, increase
linearly with metallicity, and species with two metal atoms, such as SO and
S2, increase as the square of metallicity (VIsscher et al 2006). A slight
complication is that CO and N2 increase linearly with metallicity because
these are the major reservoir species for C and N, respectively; hence CO2
increases as the square of metallicity (as CO $\times$ O), rather than as the
cube.

The models are not sensitive to temperature and insolation over the parameter
ranges ($1200\leq T\leq 2000$ K and $1\leq I\leq 1000$) presented here.
Insensitivity of the chemistry to $T$ and $I$ surprised us, and suggests that
thermochemical equilibrium is more important for sulfur than photochemistry or
kinetics. Minor differences are that HS is favored by higher temperatures and
SO and S2 are favored by high $I$. Not shown here is that the chemistry
changes markedly for $T<1100$ K: hydrocarbons, CS, and CS2 become abundant,
and the results become sensitive to $K_{zz}$. Cooler atmospheres introduce a
variety of new topics best left for another study.

Carbon dioxide, a robust molecule and a potential observable, has been
reported in HD 189733b by Swain et al (2009). CO2 is generated from CO by
reaction with OH radicals. The chief source of OH is the reaction of H2O with
atomic hydrogen; at high altitudes UV photolysis of H2O is also important. We
find that CO2 mixing ratios range from $1.6\times 10^{-8}$ to $1.6\times
10^{-4}$ for $-0.3\leq[M/H]\leq 1.7$, scaling as the square of metallicity.
Table 1 lists computed CO2 mixing ratios in the models discussed here. These
results are insensitive to insolation, vertical mixing, temperature between
1200 K and 2000 K, and gravity. The CO2/CO ratio is nearly independent of
pressure, as seen in Figure 2. Pressure independence is expected because the
controlling reactions, CO2+H $\leftrightarrow$ CO + OH and H+H2O
$\leftrightarrow$ H2 \+ OH, and the controlling equilibrium, CO2+H2
$\leftrightarrow$ CO + H2O, all leave the total pressure unchanged. (At very
high altitudes photochemistry alters the CO2/CO ratio.) The computed CO2
abundances are in good agreement with the reported observation of CO2 at the
ppmv level in HD 189733b (Swain et al. 2009). The sensitivity of CO2 to
metallicity and insensitivity to other atmospheric parameters makes CO2 a good
probe of planetary metallicity, as pointed out by Lodders and Fegley (2002).

###  3.1 Optical depth and stratospheric heating

Figure 3 shows the pressure levels where the solar and planetary metallicity
atmospheres of Models A, M, and MM become optically thick. Opacity is
dominated by HS, with some contribution by S2 at wavelengths shorter than 300
nm. The twin peaks between 300 nm and 320 nm may be fictitious, but the peak
at 324 nm could prove diagnostic of HS. A solar and a K0V stellar spectrum are
shown for comparison.

Figure 4 shows the magnitude of stratospheric heating and the pressure level
where the heating occurs for a solar-twin primary at 0.032 AU ($I=1000$) for 3
metallicities (Models A, M, and MM). Radiative heating is dominated by HS, and
is nearly saturated through the stratosphere for all these models (see also
Table 1). By contrast, peak heating at $\sim\\!100\,\mu$bars takes place where
SO and SO2 are significant. The sensitivity of SO and SO2 to metallicity is
reflected in greater heating rates at $\sim\\!100\,\mu$bars.

Cumulative stratospheric heating rates for these models are listed in Table 1.
For a solar-twin at 0.032 AU, cumulative heating above 1 mbar is typically
$4\times 10^{4}$ W/m2 and above 0.1 bars is typically $8\times 10^{4}$ W/m2,
i.e., about half the energy is absorbed in the lower stratosphere. Burrows et
al. (2008) modeled hot stratospheres by adding an unknown gray absorber. They
found that gray cross-sections of $0.05-0.6$ cm2/g, averaged over 430 to 1000
nm for altitudes above 0.03 bars, could produce the observed heating. Heating
profiles using gray opacities in this range are plotted for comparison on Fig
4 for the same planet and star. The gray opacities produce more heating in
total (indeed, the stratospheres in both these models are optically thick),
and more heating at low altitudes, but at higher altitudes sulfur generates
heating at levels quite similar to what Burrows et al find useful.

##  4 Conclusions

We develop a new 1D photochemical model for stratospheric modeling of
hydrogen-rich atmospheres of warm or hot exoplanets. This model is applicable
to any H-rich planet subject to high insolation, including hot Neptunes,
superearths, and waterworlds. Here we apply the model to sulfur chemistry,
stratospheric heating, and CO2 abundance.

We find that hot stratospheres of hot Jupiters could be explained by
absorption of UV and violet visible light by HS and S2, two highly reactive
species that are generated chemically from H2S. For a hot Jupiter orbiting a
solar-twin at 0.032 AU, for a wide range of possible planetary compositions,
HS and S2 together absorb $4\times 10^{4}$ W/m2 at altitudes above 1 mbar and
another $4\times 10^{4}$ W/m2 at altitudes between 1 mbar and 0.1 bar. This
level of heating approaches what Fortney et al (2006) and Burrows et al (2008)
use in their most successful LTE spectral models. Non-LTE mechanisms may
improve the agreement, because LTE models systematically overestimate
radiative cooling and thus underestimate the temperature. Chemiluminescence by
H2O, formed by the exothermic reaction of OH+H2, might also be expected.

Although our computed HS and S2 column densities increase with metallicity,
optically thick columns are predicted for all plausible atmospheric
compositions, which means that millibar-level temperature inversions are
expected to be commonplace. The distinctive interaction of S2 and HS with near
ultraviolet light could make these species detectable in transit by the
refurbished HST; there is evidence for a blue absorber in legacy HST data of
HD 209458b (Sing et al 2008).

On the other hand, sulfur does not give an easy answer to why some hot
Jupiters have superheated stratospheres, and others not. In an earlier draft
of this study, we speculated that S3—which is very sensitive to
metallicity—might be part of the explanation. This no longer appears likely.
We have since developed a better understanding of HS’s opacity, which turns
out to be considerable. We no longer see a strong connection between
metallicity and radiative heating, save at very low pressures
($<\\!100\mu$bars) where SO and SO2 become important. It now seems that sulfur
chemistry by itself is unlikely to explain differences between planets,
although planetary metallicity may still be key.

Heating by sulfur compounds does not preclude heating by TiO and VO on hotter
planets. Sulfur species provide considerable heating from below 1000 K to
above 2000 K, but they do not provide the spectral coverage at visible
wavelengths that TiO and VO provide. For TiO and VO to be abundant enough to
explain stratospheric heating, the temperature needs to be very high, in
excess of 2000 K, and not just in the stratosphere but also at deeper levels
in the planet where these two refractory oxides would otherwise be cold-
trapped in silicate clouds. OGLE-TR-56b (Sing and López-Morales 2009) seems to
meet the TiO-VO threshold.

CO2 is generated by the reaction of CO with OH and destroyed by the reverse
(endothermic) reaction with H, CO + OH $\leftrightarrow$ CO2 \+ H. At low
altitudes OH is generated by the reaction of H2O with atomic H, supplemented
at high altitudes by UV photolysis of H2O. As both the major source and major
sink of CO2 are proportional to atomic hydrogen densities, the kinetic
inhibition against hydrogen recombination does not disturb CO2’s
thermochemical equilibrium. We find that CO2 mixing ratios vary quadratically
with metallicity from $1.6\times 10^{-8}$ to $1.6\times 10^{-4}$ for $0<[{\rm
M}/{\rm H}]<0.7$. This result is insensitive to insolation, vertical mixing,
temperature (for $1200\leq T\leq 2000$ K), and gravity. Because the reactions
that form and destroy CO2 leave the total number of molecules unchanged, the
CO2/CO ratio is also pressure independent. The computed CO2 abundances are in
good agreement with the observation of CO2 at the ppmv level in HD 189733b
(Swain et al 2009). Therefore we confirm Lodders and Fegley’s (2002)
suggestion that CO2 is a promising probe of planetary metallicity.

##  5 Acknowledgements

We thank R. V. Yelle for discussions regarding the potential importance of S3,
and G. Tinetti for an insightful review. We thank NASA’s Exobiology and
Planetary Atmospheres Programs for support. KL was also supported by NSF Grant
AST-0707377.

## References

  * (1)
  * (2) Allard, F., Hauschildt, P. H., Alexander, D. R., Tamanai, A., & Schweitzer, A. 2001, ApJ, 556, 357
  * (3)
  * (4) Barman, T. S. 2008, ApJ, 676, L61
  * (5)
  * (6) Barnes, J. R., Barman, T. S., Jones, H. R. A., Leigh, C. J., Cameron, A. C., Barber, R. J., & Pinfield, D. J. 2008, MNRAS, 390, 1258
  * (7)
  * (8) Berdyugina, S.V. and Livingston, W.C. (2002) Astron. Astrophys. 387, L6-L9.
  * (9)
  * (10) Billmers RI and Smith AL (1991). J. Phys. Chem. 95, 4242-4245.
  * Burrows et al. (2007) Burrows, A., Hubeny, I., Budaj, J., Knutson, H. A., & Charbonneau, D. 2007, ApJ, 668, L171
  * Burrows et al. (2008) Burrows, A., Budaj, J., & Hubeny, I. 2008, ApJ, 678, 1436
  * Burrows et al. (2008) Burrows, A., Budaj, J., & Hubeny, I. 2008, ApJ, 678, 1436
  * Charbonneau et al. (2008) Charbonneau, D., Knutson, H. A., Barman, T., Allen, L. E., Mayor, M., Megeath, S. T., Queloz, D., & Udry, S. 2008, ApJ, 686, 1341
  * (15) Claire MW, Catling DC, Zahnle KJ (2006). Geobiology 4, 239-269.
  * (16) Dalgarno A and Williams DA (1962). . Astrophy. J. 136, 690-692.
  * (17) Désert, J.-M., Vidal-Madjar, A., Lecavelier Des Etangs, A., Sing, D., Ehrenreich, D., Hébrard, G., & Ferlet, R. 2008, A&A, 492, 585
  * Fegley & Lodders (1994) Fegley, B. J. & Lodders, K. 1994, Icarus, 110, 117
  * Fortney et al. (2006) Fortney, J. J., Saumon, D., Marley, M. S., Lodders, K., & Freedman, R. S. 2006, ApJ, 642, 495
  * Fortney et al. (2008) Fortney, J. J., Lodders, K., Marley, M. S., & Freedman, R. S. 2008, ApJ, 678, 1419
  * Harrington et al. (2007) Harrington, J., Luszcz, S., Seager, S., Deming, D., & Richardson, L. J. 2007, Nature, 447, 691
  * (22) Henneker, W. H. and Popkie, H.E. (1971) J. Chem. Phys. 54, 1763-1778.
  * Hubeny et al. (2003) Hubeny, I., Burrows, A., & Sudarsky, D. 2003, ApJ, 594, 1011
  * (24) Kasting JF, Zahnle KJ, Pinto JP, Young AT (1989) Origins of Life, 19, 95-108.
  * (25) Kasting JF (1990). Origins of Life 20, 199-231.
  * Knutson et al. (2008) Knutson, H. A., Charbonneau, D., Allen, L. E., Burrows, A., & Megeath, S. T. 2008, ApJ, 673, 526
  * Knutson et al. (2009) Knutson, H. A., Charbonneau, D., Burrows, A., O’Donovan, F. T., & Mandushev, G. 2009, ApJ, 691, 866
  * (28) Larsson, M. (1983) Astron. Astrophys. 128, 291-298.
  * (29) Liang M-C, Parkinson CD, Lee AYT, Yung YL, and Seager S, (2003) “Source of atomic hydrogen in the atmosphere of HD 209458b” Astrophys. J. 596, L247 L250.
  * (30) Liang M-C, Seager S, Parkinson CD, Lee AYT, and Yung YL (2004) “On the insignificance of photochemical hydrocarbon aerosols in the atmospheres of close-in extrasolar giant planets” Astrophys. J. 605, L61 L64.
  * (31) Lodders, K. and Fegley B. (1999) The Planetary Scientist’s Companion. Oxford.
  * Lodders (1999) Lodders, K. 1999, ApJ, 519, 793
  * Lodders (2002) —. 2002, ApJ, 577, 974
  * (34) Lodders, K. and Fegley, B.J. (2002). Icarus 155, 393-424.
  * (35) Machalek P, McCullough PR, Burke CJ, Valenti JA, Burrows A, and Hora JL (2008) Astrophys. J. 684, 1427-1432.
  * Marley et al. (2007) Marley, M. S., Fortney, J., Seager, S., & Barman, T. 2007, in Protostars and Planets V, ed. B. Reipurth, D. Jewitt, & K. Keil, 733–747
  * (37) Moses JI, Zolotov MY, and Fegley B (2002). Icarus 156, 76 106.
  * (38) Nicholas, J.E., Amodio, C.A., Baker, M.J. (1979) J. Chem. Soc. Faraday Trans. 1:75, 1868-1880.
  * (39) Noll KS, McGrath MA, Trafton LM, Atreya SK, Caldwell JJ, Weaver HA, Yelle RV, Barnet C, and Edgington S (1995) Science, 267, 1307\.
  * (40) Okabe H (1978) The Photochemistry of Small Molecules. Wiley-Interscience, New York, 431 pp.
  * (41) Ram, R S, Bernath P. F., Engleman R. and Brault J. W. (1995) J. Molec. Spect. 172, 34-42.
  * (42) Resende, S.M., and Ornellas, F.R. (2001) J. Chem. Phys. 115, 2178-2187.
  * Richardson et al. (2007) Richardson, L. J., Deming, D., Horning, K., Seager, S., & Harrington, J. 2007, Nature, 445, 892
  * (44) Sauval, A.J. and Tatum, J.B. (1984) Astrophs. J. Supp. 56, 193-209.
  * (45) Schneider, L., Meier W, and Weige KH (1990) J. Chem. Phys. 92, 7027-7037.
  * Showman et al. (2008) Showman, A. P., Fortney, J. J., Lian, Y., Marley, M. S., Freedman, R. S., Knutson, H. A., & Charbonneau, D. 2008, ArXiv e-prints
  * (47) Sing, D.K., Vidal-Madjar A., Désert, J.-M., Lecavelier des Etangs, A., and Ballester, G. (2008). Astrophys. J. 686, 658-666.
  * (48) Sing, D.K., and López-Morales, M. (2009). Astron. Astrophys. 493, L31-L34.
  * (49) Spencer JR., Jessup, KL, McGrath MA, Ballester GE, amd Yelle RV (2000). Science 288, 1208-1210.
  * (50) Steudel R and Eckert B (2003). Elemental sulfur and sulfur-rich compounds. Springer, 202 pp.
  * (51) Swain, M. R., Vasisht, G., Tinetti, G., Bouwman, J., Chen, P., Yung, Y., Deming, D., & Deroo, P. 2009, Astrophys. J. Lett. 690, L114-L117.
  * (52) Tinetti, G. (2008). Bull. Am. Astron. Soc. 40, 463\.
  * Torres et al. (2008) Torres, G., Winn, J. N., & Holman, M. J. 2008, ApJ, 677, 1324
  * (54) van der Heijden and van der Mullen (2001). J. Phys. B. Atom. Mol. Opt. Phys. 34, 4183-4201.
  * (55) Visscher C, Lodders K, and Fegley B. (2006). Astrophys. J. 648, 1181-1195.
  * (56) Wheeler, M.D., Orr-Ewing, AJ, and Ashfold, MNR (1997) J. Chem. Phys. 107, 7591-7600.
  * (57) Zahnle KJ, Mac Low M-M, Lodders K, and Fegly B. (1995). Geophys. Res. Lett. 22, 1593-1596.
  * (58) Zahnle KJ, Claire MW, Catling DC (2006). Geobiology 4, 271-282.
  * (59) Zahnle KJ, Haberle RM, Catling DC, Kasting JF (2008). J. Geophys. Res. 113, E11004, doi:10.1029/2008JE003160.
  * (60) Zare, R.N., Schmeltekopf AL, Harrop WJ, and Albritton DL (1973) J. Molec. Spect. 46, 37-66.

Model | $[{\rm M}/{\rm H}]^{a}$ | $I^{b}$ | $T$ |  S2 [cm-2]c |  HS [cm-2]c |  HSd |  SOd |  CO${}_{2}^{d}$ |  Heatinge |  Heatingf
---|---|---|---|---|---|---|---|---|---|---
A | $0$ | 1000 | 1400 | $4.2\times 10^{18}$ | $1.2\times 10^{20}$ | 6 | $0.07$ | $0.065$ | $6.4\times 10^{4}$ | $2.7\times 10^{4}$
M | $0.7$ | 1000 | 1400 | $2.2\times 10^{20}$ | $1.0\times 10^{21}$ | 26 | $1.7$ | $1.6$ | $8.4\times 10^{4}$ | $4.1\times 10^{4}$
MM | $1.4$ | 1000 | 1400 | $6.0\times 10^{21}$ | $5.4\times 10^{21}$ | 100 | $27$ | $41$ | $1.1\times 10^{5}$ | $5.1\times 10^{4}$
H | $0.7$ | 1000 | 1600 | $2.2\times 10^{20}$ | $2.3\times 10^{21}$ | 43 | $1.2$ | $1.4$ | $9.0\times 10^{4}$ | $4.3\times 10^{4}$
HH | $0.7$ | 1000 | 1800 | $1.9\times 10^{20}$ | $4.0\times 10^{21}$ | 52 | $1.3$ | $1.5$ | $9.5\times 10^{4}$ | $4.4\times 10^{4}$
HHH | $0.7$ | 1000 | 2000 | $1.3\times 10^{20}$ | $5.0\times 10^{21}$ | 41 | $1.4$ | $1.3$ | $9.6\times 10^{4}$ | $4.1\times 10^{4}$
C | $0.7$ | 1000 | 1200 | $1.6\times 10^{20}$ | $2.5\times 10^{20}$ | 11 | $1.9$ | $1.9$ | $7.5\times 10^{4}$ | $3.7\times 10^{4}$
G | $0.7$ | 1000 | 1400 | $4.3\times 10^{20}$ | $2.0\times 10^{21}$ | 32 | $1.3$ | $1.6$ | $9.3\times 10^{4}$ | $4.7\times 10^{4}$
I | $0.7$ | 200 | 1400 | $2.1\times 10^{20}$ | $1.0\times 10^{21}$ | 37 | $0.9$ | $1.6$ | $1.7\times 10^{4}$ | $8.8\times 10^{3}$
SSC | $0.7$ | 1 | 1200 | $1.1\times 10^{20}$ | $2.4\times 10^{20}$ | 16 | $0.06$ | $1.9$ | $72$ | $44$
$a$ – Metallicity. This notation means that the planet is $10^{[{\rm M}/{\rm
H}]}$ richer in C, S, N, and O than the Sun.
$b$ – Insolation. $I=1000$ corresponds to a solar twin primary at 0.032 AU.
$c$ – Column densities above 1 bar.
$d$ – Mixing ratio in ppmv at 1 mbar.
$e$ – Total atmospheric heating [W/m2] above 0.1 bar for a solar twin source.
$f$ – Total atmospheric heating [W/m2] above 1 mbar for a solar twin source.
Figure 1: Theoretical absorption cross sections of HS radicals at near UV,
violet and indigo wavelengths at four temperatures at 30 mbar pressure. Cross
sections were computed from the lowest five vibrational levels of the ground
electronic state $X^{2}\Pi$ to the lowest three vibrational levels of the
upper level $A^{2}\Sigma^{+}$. The excited vibrational levels of
$A^{2}\Sigma^{+}$ are strongly predissociating, which suggests that absorption
at wavelengths shorter than 324 nm is probably continuous rather than
allocated into the well-defined bands shown here.  Figure 2: Important sulfur
species, CO, and CO2 in the atmosphere of a hot Jupiter with a “planetary”
metallicity of $[{\rm M}/{\rm H}]=0.7$. The atmosphere is assumed isothermal
at 1400 K and insolated 1000$\times$ more strongly than Earth. Other model M
parameters are listed in Table 1. The prominent transition at $\sim$2 mbar —
where the S2 mixing ratio peaks — is associated with photolysis of H2S. The
bump in CO2 at 6 $\mu$bars is attributable to photochemistry. Abundance
profiles in the 1400 K atmosphere are generally representative of atmospheres
with $1200\leq T\leq 2000$ K. Figure 3: Pressure levels of the $\tau=1$
surface as a function of wavelength for three metallicities, $[{\rm S}/{\rm
H}]=0$, 0.7, and 1.4. These metallicities correspond to models A, M, and MM of
Table Atmospheric Sulfur Photochemistry on Hot Jupiters. Absorption between
250 and 300 nm is mostly by S2 and absorption between 300 and 460 nm is by HS.
Structure blueward of 324 nm is associated with transitions to predissociating
states and is probably fictitious. The $\tau=1$ surface of a pure H2 Rayleigh
scattering atmosphere and two incident stellar spectra, one for the Sun and
another for a generic K0V dwarf, are shown for comparison.  Figure 4:
Radiative heating at different altitudes for three metallicities, $[{\rm
S}/{\rm H}]=0$, 0.7, and 1.4. These correspond to models A, M, and MM of Table
Atmospheric Sulfur Photochemistry on Hot Jupiters. Heating rates are given in
W/kg, which emphasizes the potential impact on temperature. Heating peaks at
100 $\mu$bars but extends through the stratosphere. Heating with constant gray
opacities of $0.05$ and $0.6$ cm2/g for $430<\lambda<1000$ nm is shown for
comparison.


# Polarization-artifact reduction and accuracy improvement of Jones-matrix
polarization-sensitive optical coherence tomography by multi-focus averaging

Lida Zhu  1 Shuichi Makita  1 Junya Tamaoki  2 Yiqiang Zhu  1 Yiheng Lim  1
Makoto Kobayashi  2 and Yoshiaki Yasuno1  1Computational Optics Group,
University of Tsukuba, Tsukuba, Ibaraki, Japan
2Department of Molecular and Developmental Biology, Institute of Medicine,
University of Tsukuba, Tsukuba, Ibaraki, Japan
<EMAIL_ADDRESS>
###### Abstract

Polarization-sensitive optical coherence tomography (PS-OCT) is a promising
biomedical imaging tool for differentiation of various tissue properties.
However, the presence of multiple-scattering (MS) signals can degrade the
quantitative polarization measurement accuracy. We demonstrate a method to
reduce MS signals and increase the measurement accuracy of Jones matrix PS-
OCT. This method suppresses MS signals by averaging of multiple Jones matrix
volumes measured using different focal positions. The MS signals are
decorrelated among the volumes by focus position modulation and are thus
reduced by averaging. However, the single scattering signals are kept
consistent among the focus-modulated volumes by computational refocusing. We
validated the proposed method using a scattering phantom and a postmortem
medaka fish. The results showed reduced artifacts in birefringence and degree-
of-polarization uniformity measurements, particularly in deeper regions in the
samples. This method offers a practical solution to mitigate MS-induced
artifacts in PS-OCT imaging and improves quantitative polarization measurement
accuracy.

††journal: opticajournal

https://optics.bk.tsukuba.ac.jp/COG/

##  1 Introduction

Optical coherence tomography (OCT) is a non-invasive and label-free imaging
modality that visualizes the three-dimensional (3D) morphology of a sample
using the backscattering intensity as a contrast source [1]. Polarization-
sensitive OCT (PS-OCT), a functional extension of OCT, provides additional
contrasts that are sensitive to a diverse range of tissue properties[2, 3].
One example of these contrasts is the birefringence, which is effective for
identification of fibrous tissues such as muscle and collagen[4, 5, 6, 7].
Another example is the degree-of-polarization uniformity (DOPU), which
measures the polarization randomness. The DOPU is sensitive to the
polarization scrambling property of the sample and has been used in detection
and diagnosis of the retinal pigment epithelium layer [8] and the choroid [9,
10] in the human retina. These contrasts in PS-OCT enable comprehensive
investigations of tissue properties in biological studies[11, 12]

Typically, both OCT and PS-OCT rely on use of single-back-scattering (SS)
signals to reconstruct the images of a sample[13]. However, in practice, OCT
and PS-OCT images both deteriorate by the multiple-scattering (MS) signals,
especially when the sample exhibit strong scattering. In such cases, the SS
signal at a specific depth is overlaid by MS signals originating from
different depths, and this effect causes contrast degradation of the image.
Additionally, these MS signals can interfere with polarization measurements.
For example, Lichtenegger et al.measured zebrafish using a Jones-matrix PS-OCT
(JM-OCT) and noticed erroneously increased birefringence and reduced DOPU
values[14].

To mitigate the image degradation caused by the MS signals, several wavefront
shaping-based methods have been demonstrated. In these methods, the probe
beam’s wavefront is modulated using a control device such as a spatial light
modulator (SLM) [15, 16], a deformable mirror (DM) [17], or an electrical
tunable lens (ETL) [18]. Because the modulation process decorrelates the MS
and SS signals, only the MS signals can be eliminated by the post-acquisition
signal processing methods. Although the methods noted above improved the
contrast of OCT intensity images, their effects on PS-OCT imaging have not
been demonstrated to date.

In this study, we modify and adapt our previously demonstrated wavefront
shaping-based MS signal reduction method, the multi-focus averaging (MFA)
method [18], for application to polarization imaging in JM-OCT. In this
method, we simultaneously modulate four OCT images that form a measured Jones
matrix using an ETL. Specifically, multiple JM-OCT volumes are acquired
sequentially with different defocus configurations. The defocus is corrected
by a computational refocusing process, and this process decorrelates the SS
and MS signals. Subsequent complex averaging of the volumes reduces the MS
signals while retaining the SS signals. The proposed method was validated by
measurements of a phantom and a postmortem medaka fish. We found that the MS
signals caused errors in the DOPU and birefringence measurements, particularly
in the deeper regions in which the MS signals are prominent, and also found
that these errors can be mitigated using our MFA method.

##  2 Principle and implementation

###  2.1 Principle

The JM-OCT measures the cumulative Jones matrix of the light that is
backscattered from a sample and multiple contrasts can then be computed from
the measured Jones matrix. The measured Jones matrix $\mathbf{J}_{m}$ at a
depth $z$ can be represented by

$\mathbf{J}_{m}(x,y;z)={\mathbf{J}_{\mathrm{out}}}\mathbf{J}_{s}(x,y;z){\mathbf{J}_{\mathrm{in}}},$
(1)

where $x$ and $y$ are the lateral positions, and ${\mathbf{J}_{\mathrm{in}}}$
and ${\mathbf{J}_{\mathrm{out}}}$ represent the Jones matrices of the
illumination and collection paths of the system, respectively.
$\mathbf{J}_{s}$ is the round-trip cumulative Jones matrix of the sample
within a small depth region around a depth of $z$. Notably, $\mathbf{J}_{m}$
has four entries that correspond to four complex OCT images (or volumes)
acquired using the four polarization channels of the JM-OCT system.

We then apply our previously demonstrated MFA method to each of the entries
(i.e., to each of the complex OCT volumes) in the measured Jones matrix
$\mathbf{J}_{m}$. Because the principle of MFA can be found in detail
elsewhere [18], we present only a brief summary of the method here. Assuming
that the depth $z$ in the OCT image experiences a defocus amount of $z_{d}$,
the complex OCT signal of the $j$-th entry of $\mathbf{J}_{m}$ can then be
written as

$S^{j}(x,y;z,z_{d})=S^{j}_{SS}(x,y;z,0)*\exp\left[i\phi(x,y;z,z_{d})\right]+S^{j}_{MS}(x,y;z,z_{d}),$
(2)

where $*$ denotes the convolution operation, and $S^{j}_{SS}$ and $S^{j}_{MS}$
denote the SS and MS signal components, respectively. $\phi$ is a defocus-
dependent quadratic phase function and the phase term convolved with the SS
component term represents the defocus.

By applying a computational refocusing process to the complex OCT signal of
the $j-$th entry, the signal then becomes

${S^{j}}^{\prime}(x,y;z,z_{d})=S^{j}_{SS}(x,y;z,0)+S^{j}_{MS}(x,y;z,z_{d})*\exp\left[-i\phi(x,y;z,z_{d})\right].$
(3)

As shown in Eq. (3), the refocused SS signals $S^{j}_{SS}$ remain consistent,
regardless of the different defocus amounts. In contrast, the MS lights follow
practically random paths and are thus altered by the defocus. Therefore, the
MS signals cannot remain consistent with variations in the defocus amount,
even after computational refocusing. Therefore, averaging of the multiple
focus-modulated OCT signals after application of the refocusing process can
reduce the MS signals while preserving the SS signals as follows:

$\overline{{S^{j}}^{\prime}}(x,y;z)=S^{j}_{SS}(x,y;z,0)+\frac{1}{N}\sum_{i=0}^{N-1}S^{j}_{MS}(x,y;z,z_{d,i})*\exp\left[-i\phi(x,y;z,z_{d,i})\right],$
(4)

where $\overline{{S^{j}}^{\prime}}$ is the averaged OCT signal, $N$ is the
number of defocus-modulated OCT volumes for averaging, and $z_{d,i}$ denotes
the defocus amount of the $i-$th acquisition.

Because all four entries, i.e., the four complex OCT volumes contained in a
volumetric measured Jones matrix, are assumed to experience equal defocus, the
same computational refocusing process can be applied to each of the four
entries. The OCT intensity, the birefringence, and the DOPU are then computed
from a new Jones matrix constructed using the averaged OCT signals. The
implementation of this method will be described in detail in Section 2.2.2

###  2.2 Implementation

####  2.2.1 JM-OCT setup

A custom-built JM-OCT system was used in this study. The system uses a
wavelength-swept light source (AXP50124-8, Axsun Technologies, MA) with a
wavelength centered at 1,310 nm and a scan width of approximately 106 nm when
measured at -10 dB width. The axial resolution was measured to be 14 $\mu$m in
tissue. The A-line rate of this system is 50 kHz. Use of an objective with an
effective focal length of 36 mm (LSM03, Thorlabs, NJ) and an incident beam
diameter of 3.49 mm on the objective provides a lateral resolution of 18
$\mu$m and a depth-of-focus (DOF) of 0.36 mm in air. An ETL (EL-10-30-CI-NIR-
LD-MV, Optotune, Switzerland) was integrated into the sample arm to modulate
the defocus for the MFA method. Further details about the system with the
exception of the ETL part can be found in Refs. [19, 20], and details of the
ETL implementation can be found in Section 4.4 of Ref. [18].

####  2.2.2 MFA signal processing for JM-OCT

The JM-OCT-based MFA measures multiple volumetric Jones matrices with the
focus at different depth positions, i.e., with different defocus amounts.
These volumetric Jones matrices consist of four complex OCT volumes that
correspond to the entries in the Jones matrix. The volumetric Jones matrices
are processed using a signal processing flow that consists of bulk-phase error
correction, computational refocusing, inter-volume axial shift correction, and
inter-volume phase offset correction, and the resulting matrices are then
averaged to reduce the MS signals. Multi-contrast OCT images, including the
OCT intensity, DOPU, and birefringence are computed from the averaged Jones
matrix volume. This process is similar to our previous MFA method but has been
modified to deal with the Jones matrix signals. The processing details are
described in the following.

First, for each acquired volumetric Jones matrix, we estimate and correct the
bulk phase errors using a smart-integration-path method [21]. The bulk phase
error is estimated from a single entry in the Jones matrix and the errors in
all entries are then corrected using this estimate.

Second, a computational refocusing process was applied to correct the defocus.
We estimated the defocus amount at each depth from a single entry in the
volumetric Jones matrix and then corrected the defocus for all entries using
this defocus amount. Details of the refocusing process have been presented in
Section 2.2 of Ref. [22].

The deformation of the liquid lens in the ETL causes an axial shift in the
image. In the third step, this shift is estimated and corrected using a sub-
pixel linear intensity cross-correlation method [23]. The axial shift is
estimated from a single entry in the volumetric Jones matrix and this shift is
used to correct all four entries. Details of the process for estimation for a
single entry are given in Section 2.2.2 in Ref. [18].

In the fourth step, the inter-volume phase offsets are estimated using all
four entries (i.e., the volumes), while they were estimated from a single
entry in our previous MFA method. In this estimation process, one volumetric
Jones matrix among the Jones matrices that were obtained with the different
defocus amounts is designated as a reference. The product of the reference
volume and the complex conjugate of another volume, which is regarded as the
target volume, is computed for each entry. These products are averaged in
complex form along the depth direction and are then averaged over the four
entries, thus allowing a complex en face map to be obtained. The phase of the
en face map represents the relative phase offset of the target volume with
respect to the reference volume. The estimated phase offset is used to correct
the phase offsets of all four entries in the target volume.

After all the aforementioned corrections have been performed, the corrected
Jones matrix volumes are averaged in complex form, and multi-contrast OCT
images, including the OCT intensity, the DOPU, and the birefringence, are
computed from the averaged Jones matrix. The intensity image is the intensity
average of all four entries. In our implementation, the DOPU is computed by
applying Makita’s noise correction [24] using a spatial kernel of 3 $\times$ 3
pixels (17.5 $\mu$m in the lateral direction $\times$ 21.7 $\mu$m in the axial
direction). For the birefringence computation (i.e., the local phase
retardation computation), we first compute the local Jones matrix with an
axial separation of 8 pixels (57.9 $\mu$m in tissue). The birefringence is
then computed using the maximum a posteriori (MAP) birefringence estimator
[25, 26]. The estimation kernel used here is a single pixel, i.e., the
estimator morphed the single value of the measured birefringence into its
maximum likelihood estimation using a numerically obtained likelihood function
and the measured effective signal-to-noise ratio (SNR). Details of the image
formation processes for these contrasts can be found in Section 2.3 of Ref.
[22].

###  2.3 Validation study design

####  2.3.1 Samples

Figure 1:  A scattering phantom sample and a postmortem medaka fish were used
to validate the MS-signal rejection method. (a) and (b) show a photograph and
a schematic of the scattering phantom, respectively. This phantom is composed
of a scattering layer sandwiched between two glass plates, where the
scattering layer is a mixture of polystyrene microparticles with a diameter of
10 $\mu$m and ultrasound gel. The particle density is 0.2 g/mL. (c) shows a
photograph of the postmortem medaka fish, where the yellow box and the dashed
line indicate the measurement area and the location of the cross-sectional
images shown in Figs. 4–6, respectively.

To validate the proposed method, we measured two sample types. The first was a
non-birefringent scattering phantom that was similar to the phantom used in
our previous MFA study [18]. A photograph of this phantom is shown in Fig.
1(a). As the phantom schematic shows [Fig. 1(b)], the phantom consists of a
glass slip at the top, a scattering layer in the center, and another glass
slip with black tape at the bottom. The scattering medium of this phantom was
composed of a mixture of 10%-concentration polystyrene microparticles
(diameter of 10 $\mu$m , 72986-10ML-F, Sigma-Aldrich, St Louis, MO) and
transparent ultrasound gel (Pro Jelly, Jex, Japan). The particle density of
the resulting phantom is 0.2 g/mL.

The second sample is an adult postmortem medaka fish, which contains
birefringent tissues such as muscle. A photograph of the sample fish is shown
in Fig. 1(c). Before conducting the measurements, we anesthetized the fish
using tricaine and then sacrificed it by placing it on ice for 2 min. We
placed the fish in a petri dish and immersed it in a saline solution to
perform the measurements. This protocol follows the animal experiment
guidelines of the University of Tsukuba and is approved by the Institutional
Animal Care and Use Committee of the University of Tsukuba.

####  2.3.2 Measurement protocols

During the measurements of both the phantom and the medaka fish, we acquired
seven volumes with different focal positions for the averaging process. We set
the interval between the focal positions at 0.12 mm, which represents
approximately one-third of the DOF. It results in an overall focus shifting
distance of 0.72 mm for the set of seven acquisitions. We determined these
parameters based on an experimental optimization process that is described in
Section 4.2 of Ref. [18]. The images obtained by this method are referred to
as the MFA images in the later sections.

We also acquired another set of volumes without focus shifting for comparison
purposes. We processed these non-focus-shifting volumes using the same
processing flow that was described in Section 2.2.2, and the resulting images
are referred to here as the single focus averaging (SFA) images.

In addition, images that were obtained from a non-averaged but defocus-
corrected single Jones matrix were computed and are referred to here as
single-acquisition images. The single-acquisition, SFA, and MFA methods are
compared qualitatively and quantitatively.

During the measurements, the focal positions were set at depths of more than 1
mm below the sample surface. The lateral imaging area was 1.5 mm $\times$ 1.5
mm, and this area was sampled using 256 $\times$ 256 A-lines. This scan
protocol results in an isotropic lateral pixel separation of 5.86 $\mu$m
$\times$ 5.86 $\mu$m , which is approximately one-third of the lateral optical
resolution. The total acquisition time for the seven volumes was approximately
20 s.

##  3 Results

Note that all OCT intensity images presented in this paper are displayed in dB
scale, where 0 dB refers to the noise floor in the air region in the single-
acquisition images.

###  3.1 Scattering phantom

####  3.1.1 OCT intensity and DOPU of the phantom

Figure 2:  OCT intensity and DOPU measurement results for the scattering
phantom. (a)–(c) and (d)–(f) show the cross-sectional intensity and DOPU
images of the phantom, respectively. (g)–(l) show the en face DOPU images at
two depths (Z1 and Z2) indicated by the arrows in (d) and by the dashed lines
in (m). The polystyrene particles show high DOPU (red) and the surrounding
space shows low DOPU (yellow). When compared with the single-acquisition
results, the SFA and MFA images show a reduction in the low DOPU signals. The
mean DOPU depth profiles of the three methods in (m) show that the DOPU of the
MFA method is higher than that of the other methods. (n) and (o) are
normalized histograms computed from (g)–(i) and (j)–(l), respectively. The
dashed vertical lines in (n) and (o) indicate the modes of each histogram. The
focal positions are indicated by orange arrowheads in (a)–(c). The histogram
at the superficial depth of Z1 in (n) shows the reduction in the low DOPU
pixels in the SFA and MFA images when compared with the single-acquisition
image. The histogram at the deeper depth of Z2 in (o) demonstrates the
increased mode values of the DOPU for the SFA and MFA when compared with the
single-acquisition case.

Figures 2(a)–(l) show the intensity and DOPU images of the phantom obtained
using the single-acquisition, SFA, and MFA methods. The depth profiles of the
mean DOPU obtained by these three methods are shown in Fig. 2(m), and these
profiles were computed by averaging 10,000 A-lines (50 B-scans, with 200
A-lines per B-scan). The normalized histograms of the DOPU shown in Figs. 2(n)
and (o) were computed from the en face DOPU images at two depths using a
shifted average histogram method [27] with a bin size of 0.03125. An intensity
threshold for the noise floor of +5 dB was applied to exclude the pixels with
very low scattering intensities when computing the mean DOPU depth profiles
and the corresponding histograms.

Among the three cross-sectional intensity images [Figs. 2(a)–(c)], the MFA
image shows lower background intensity in most regions of the scattering layer
when compared with the other two images.

In the cross-sectional single-acquisition DOPU image [Fig. 2(d)], the
polystyrene microparticles show high DOPU (red), while the surrounding space
shows low DOPU (yellow). Because the polystyrene spheres may not scramble the
polarization, this appearance of low DOPU may be the artifact caused by the
overlaying MS signals. The mean DOPU depth profile [Fig. 2(m), gray curve]
shows a sharp reduction along the depth until approximately 0.36 mm, and the
profile becomes stably low in the deeper regions. When compared with the
single-acquisition DOPU, the SFA cross-sectional DOPU image shows slightly
fewer low-DOPU artifacts [Fig. 2(e)]. The mean-DOPU depth profile of the SFA
image [Fig. 2(m), green curve] is slightly higher than that of the single-
acquisition method, but the difference is not evident. In contrast, the MFA
image shows an evident reduction in the low-DOPU artifacts [Fig. 2(f)], and
its mean-DOPU depth profile shows much higher mean DOPU values than those
observed for the other two methods [Fig. 2(m)].

Figures 2(g), (h), and (i) show the en face DOPU images acquired at a
superficial depth [Z1, as indicated in Fig. 2(d)]. In the single-acquisition
image [Fig. 2(g)], the polystyrene particles show high DOPU and are enclosed
by the low-DOPU artifact. The low-DOPU artifact is reduced in the SFA image
[Fig. 2(h)] because of the reduction in the standard noise (i.e., shot noise,
relative intensity noise, and detection noise) caused by the complex averaging
process. Further reduction in the low-DOPU artifact is observed in the MFA
image [Fig. 2(i)]. Normalized histograms of these en face images are shown in
Fig. 2(n). Although the histograms of the SFA and MFA images are sharper than
that of the single-acquisition image, the mode values for these three images
are all high and are nearly identical, with values of 0.967 (single-
acquisition), 0.968 (SFA), and 0.968 (MFA). These consistently high mode
values indicate that there are fewer MS signals at this superficial depth,
which would also explain the nearly identical histogram shapes for the SFA and
MFA methods. In addition, the sharpening of the histogram may be accounted for
by the reduction in the standard noise.

At the deeper depth of Z2 [as indicated in Fig. 2(d)], more low-DOPU artifacts
are presented in the single-acquisition en face image [Fig. 2(j)] than in the
image at the superficial depth Z1. The SFA image [Fig. 2(k)] shows a reduction
in the low-DOPU artifacts [Fig. 2(k)] and the MFA image [Fig. 2(l)] shows a
further reduction in these artifacts. In the corresponding histograms [Fig.
2(o)], the mode value is low for the single-acquisition image (0.746), is
higher for the SFA image (0.778), and is then at its highest for the MFA image
(0.832). The low mode value of the single-acquisition image may be caused by
the standard noise and MS signals. The higher mode value of the SFA image may
be accounted for by the reduction in the standard noise, and the even higher
mode value of the MFA image can be accounted for by the reduction in the MS
signals. These findings indicate that the MFA method is beneficial,
particularly in the deeper regions.

####  3.1.2 Birefringence of the phantom

Figure 3:  Birefringence measurement results for the scattering phantom.
(a)–(c) and (d)–(i) show the cross-sectional and en face slices of the
birefringence images, respectively. These images correspond to those of Figs.
2(d)–(l). At the superficial depth of Z1, the polystyrene particles show
moderately high birefringence (green) and are surrounded by a moderately high
birefringent halo. The halo is reduced by the MFA method when compared with
the other methods. At the deeper depth (Z2), the particles are surrounded by a
halo with high birefringence, which may be artifacts caused by the dominance
of the MS signals. (j) shows the mean birefringence depth profiles for the
three methods. Here, the single-acquisition method shows slightly higher
birefringence than the other methods. (k) and (l) show the histograms computed
from (d)–(f) (at depth z1) and from (g)–(i) (at depth Z2), respectively. The
dashed vertical lines in (k) and (l) indicate the modes of each histogram. The
axial label “biref.” represents an abbreviation of birefringence. Although the
histograms at the deeper depth in (l) do not show evident differences among
the three methods, the corresponding histograms at the superficial depth in
(k) showed a reduction in the high-birefringence pixels with the SFA and MFA
methods.

Figures 3(a)–(i) show the birefringence images of the phantom, where the
positions in the sample are identical to those in Fig. 2. At the superficial
depth (Z1), the polystyrene particles that show moderate birefringence (green)
are surrounded by a moderately birefringent (green) halo. Because the
surrounding gel was not birefringent, the moderately birefringent halo may be
an artifact. In the JM-OCT method, the birefringence can be artifactually high
if the SNR is sufficiently low. This artifactual elevation of the
birefringence will be discussed in detail later in Section 4.2. When compared
with the single-acquisition image [Fig. 3(d)], the MFA image [Fig. 3(f)] shows
a less moderately birefringent (green) halo. This difference may be caused by
the reduction in the standard noise and the MS signals.

In the deep region (Z2), both the particles and the surrounding gel showed
high birefringence (indicated by colors ranging from yellow to red). This may
be caused by the predominance of the MS signals, and the MFA method could not
show evident improvements for this high-scattering phantom.

Figure 3(j) shows the mean birefringence depth profiles, where the three
profiles almost overlay each other, with the exception that the single-
acquisition profile shows slightly higher birefringence values. Figures 3(k)
and (l) are the birefringence histograms obtained at depths Z1 (superficial)
and Z2 (deep), respectively. At the superficial depth (Z1), some parts of the
high birefringence appearance (i.e., the right shoulder of the histogram) are
reduced by the SFA method and further reduced by the MFA method. This
reduction may occur because of the reduction in the standard noise (from
single-acquisition to SFA) and the reduction in the MS signals (from SFA to
MFA), and this behavior is consistent with the image appearances. At the
deeper depth (Z2), the histograms are all nearly identical to each other. The
ineffectuality of the SFA and MFA methods is also confirmed by the appearances
of the images, and it may be caused by the scattering being too strong.

###  3.2 Postmortem medaka fish

First, we compare the OCT intensity, DOPU, and birefringence images of the
postmortem medaka fish acquired using the three methods objectively.

In addition, for quantitative comparison of the measured DOPU and
birefringence, we manually selected four regions-of-interests (ROIs) within
the cross-sectional images. The size of each ROI is 15 (axial) $\times$ 35
(lateral) pixels (0.2 mm axial $\times$ 0.1 mm lateral). The ROIs are
indicated by white boxes in both Figs. 5(a) and 6(a). We selected the skin and
muscle regions in the upper and lower halves of the fish, where the terms
upper and lower refer to the upper and lower areas of the image, i.e., they
correspond to the right and left sides of the fish. Because the fish body
under investigation is symmetrical in terms of the right and left sides, the
upper and lower ROIs should have similar DOPU and birefringence values as long
as polarization artifacts do not exist. We applied an intensity threshold of
+5 dB to exclude pixels with very low scattering intensities from the
computation of the statistics, i.e., the means and percentiles.

####  3.2.1 OCT intensity images of the medaka fish

Figure 4:  OCT intensity images of the postmortem medaka fish. (a)–(c) show
the cross-sectional images and (d)–(i) show the en face intensity images. The
white arrows in (a)–(c) indicate the shadows in the deeper region of the
sample, which become darker in the MFA image than in the other images because
of the reduction in the MS signals. The white ellipses in (d)–(f) indicate the
fibrous muscle tissue. The contrast in this region is improved using the MFA
method. The yellow arrows in (g)–(i) indicate low-intensity shadows, and the
contrasts of these shadows are found to be improved by the MFA method. The red
arrows in (a) and the orange dashed line in (d) indicate the locations in
which the en face and cross-sectional images were acquired, respectively. The
focal positions are indicated by the orange arrowheads shown in (a)–(c).

Figure 4 shows the intensity images of the medaka fish for the single-
acquisition, SFA, and MFA methods. In the cross-sectional images [Fig.
4(a)–(c)], the shadow region on the lower side of the fish (white arrows)
becomes darker in the SFA image than in the single-acquisition image, and is
even darker in the MFA image. This change may be caused by the reduction in
the standard noise (from single-acquisition to SFA) and by the reduction in
the MS signals (from SFA to MFA).

In the en face images at the superficial depth (Z1) [Figs. 4(d)–(f)], an area
that shows fibrous structures can be observed (indicated by white ellipses),
and this is believed to be the muscle region. The contrast in this region is
improved using the MFA method.

At the deeper depth (Z2), the improvement in the contrast is more evident
[Figs. 4(g)–(i)]. The SFA image shows a slightly enhanced contrast between the
tissues and some low-intensity structures (indicated by the yellow arrows),
and the MFA image shows further enhancement. These low-intensity structures
are the projection shadows that resulted from hyper scattering or from the
absorption of the tissues at the shallower depths.

####  3.2.2 DOPU of the medaka fish

Figure 5:  DOPU measurement results for the medaka fish. (a)–(c) show the
cross-sectional images and (d)–(i) show the en face DOPU images at two depths.
When compared with the single-acquisition and SFA images, the MFA images [(c)
and (i)] show fewer low-DOPU artifacts in the deeper regions. The white boxes
in (a) indicate the ROIs selected for the quantitative analysis. (j) and (k)
show the box plots of the DOPU values in the skin and muscle ROIs,
respectively. In each box, the upper and lower boundaries denote the upper and
lower quartiles, respectively, the center lines denote the medians, the stars
denote the means, and the whiskers denote the minimum and maximum values. In
the deeper region (i.e., the area with the lower ROIs), the MFA produces
increased DOPU values, which may indicate the reduction of the MS signals
caused by the MFA method. (l) shows the differences in the mean DOPU between
the lower and upper ROIs, which is called the upper-lower difference, for the
skin and muscle regions. Both the SFA and MFA methods provided smaller upper-
lower differences of DOPU than the single acquisition method.

Figure 5 shows the DOPU images of the medaka fish. In the single-acquisition
images, both the skin and muscle regions in the upper half of the fish show
homogeneously high DOPU [see the upper half of Fig. 5(a) (cross-section) and
Fig. 5(d) (en face at superficial depth Z1)]. We also see a few small low-DOPU
(yellow) dots, which may represent the artifacts caused by the standard noise.
These small dots are eliminated by both the SFA and MFA methods.

The lower half of the fish shows numerous low-DOPU signals [see the lower half
of Fig. 5(a) (cross-section) and Fig. 5(g) (en face at the deeper depth Z2)].
On the basis of the symmetry of the fish body, these low-DOPU signals may be
artifacts. When compared with the single-acquisition images, the SFA images
show reduced artifacts [Figs. 5(b) and (h)], and the MFA images show a further
reduction in the artifacts [Figs. 5(c) and (i)].

The box-whisker plots shown in Figs. 5(j) and (k) compare the results for the
three methods for the upper and lower ROIs for the skin [Fig. 5(j)] and for
the muscle [Fig. 5(k)] ROIs. Here, the boxes indicate the upper and lower
quartiles, the center bar denotes the median, and the whiskers denote the
maximum and minimum values. At the upper ROIs, for both the skin and the
muscle, all three methods produced similarly high DOPU because the MS signals
are not dominant at these superficial depths. Although the mean DOPU values
obtained with the single-acquisition method are close to 1, they are improved
(i.e., they become higher) by the SFA method, and are improved further by the
MFA method.

At the lower ROIs, these improvements are more evident. For both the skin and
muscle ROIs, the single-acquisition method shows artifactually low DOPU, and
this DOPU is clearly reduced by the SFA method. For the skin ROI, the mean
DOPU was increased further (i.e., improved) by the MFA method. For the muscle
ROI, the mean DOPU was increased further by the MFA method (from 0.93 for the
SFA to 0.94 for the MFA).

The differences between the mean DOPU values for the upper and lower ROIs
(i.e., the values obtained by subtracting the lower ROI mean from the upper
ROI mean) are depicted in Fig. 5(l). Because of the symmetry of the fish’s
anatomy, a smaller difference represents a better reduction in the artifact.
When compared with the other two methods, the MFA method provided the smallest
upper-lower difference values for both the skin and muscle ROIs.

####  3.2.3 Birefringence of the medaka fish

Figure 6:  Birefringence imaging results of the medaka fish. (a)–(c) show the
cross-sectional birefringence images and (d)–(i) show the en face
birefringence images at two depths (the superficial depth Z1 and the deeper
depth Z2). At depth Z1 in (d)–(f), the en face images of the three methods
show a muscle region (indicated by white ellipses) that exhibits high
birefringence (green). In the single-acquisition image in (d), the neighboring
skin region (indicated by the arrowheads) shows moderate birefringence (blue
to green), but this birefringence is slightly reduced by the SFA method in (e)
and the MFA method in (f). At depth Z2 in (g)–(i), only the MFA image shows a
region with low birefringence (indicated by the white arrows), which is the
skin in the lower half of the fish [see also the region indicated by the white
arrow in (c)]. (j) and (k) show the box plots of the birefringence values in
the four ROIs that are indicated in (a). (l) shows the differences in the mean
birefringence between the lower and upper ROIs for both the skin and muscle
regions. The SFA and MFA methods reduce the birefringence, especially in the
lower ROIs. (m) shows the “birefringence contrasts” of the three methods,
where the birefringence contrast is defined as the difference between the mean
birefringence values in the low- and high-birefringence regions [red boxes in
(g)–(i)]. The contrast is improved by the SFA method and is improved further
by the MFA method. The axial label “biref.” is an abbreviation of
birefringence.

Figure 6 presents the birefringence images of the medaka fish. In the cross-
sectional images [Figs. 6(a)–(c)], the skin in the upper half of the fish
shows low birefringence (blue) for all three methods. In contrast, the skin in
the lower half shows high birefringence (red to green) in the single-
acquisition and SFA images [arrows in Figs. 6(a) and (b)], but the MFA image
shows lower birefringence (blue) than the other methods [arrow in Fig. 6(c)]
This reduction in the birefringence also can be observed in the en face images
at the deeper depth (Z2) [arrows in Figs. 6(g)–(i)]. The en face images at the
superficial depth (Z1) [Figs. 6(d)–(e)] also show slight reductions in the
birefringence (indicated by the white arrowheads) produced by the SFA and MFA
methods. This reduction may be less obvious because the MS signals were less.

Two regions showing high and low birefringence were selected using a 15
$\times$ 15 pixel window in the deeper-depth en face images [red boxes in
Figs. 6(g)–(i)]. The difference between the mean birefringence values of these
two regions, which was computed by subtracting the mean of the low
birefringence region from that of the high birefringence region and denotes
the birefringence contrast, was plotted in Fig. 6(m). When compared with the
single-acquisition method, the SFA method shows a higher birefringence
contrast, and the MFA method shows an even higher birefringence contrast. This
finding highlights the ability of the MFA method to enhance the birefringence
contrast between the skin and the neighboring muscle tissues in the deeper
regions.

Figures 6(j) and (k) show box-whisker plots for the skin ROIs [Fig. 6(j)] and
the muscle ROIs [Fig. 6(k)]. These ROIs are indicated in Fig. 6(a). At the
upper ROIs of both the skin and the muscle, all methods show similarly low
birefringence. In contrast, the SFA method shows reduced mean birefringence in
the lower ROIs when compared with the single-acquisition method, and this
birefringence is reduced further by the MFA method.

Because of the anatomical symmetry of the fish, the mean birefringence should
be similar for the upper and lower ROIs. The plots of the upper-lower
differences [Fig. 6(l)] demonstrate that the MFA method achieves the closest
similarity in these mean birefringence values for both the skin and the muscle
ROIs among the three methods.

In summary, the polarization metrics (i.e., the DOPU and the birefringence)
are found to be altered by the MS signals in the deeper regions. The results
indicate that the MFA method is less susceptible to these alterations than
both the single-acquisition and SFA methods. In other words, the MFA method
reduces these artifactual alterations effectively. A discussion of the
interaction between the MS signal and the estimates of the DOPU and
birefringence will be presented in Section 4.1.

##  4 Discussions

###  4.1 Interpretation of polarization artifacts associated with multiple
scattering

In this work, the DOPU and birefringence measured by PS-OCT were found to be
prone to the MS signals in the deeper regions. Adie et al.found that the
presence of MS signals increases the randomness of the measured polarization
states[28]. This randomization can then lead to low DOPU estimation.
Specifically, as the depth in the sample increases, stronger MS signals will
be present, and these signals cause the erroneously low DOPU.

Understanding of the interaction between the MS signal and the birefringence
measurement remains an open issue. In this study, we applied a MAP estimator
to perform the maximum likelihood estimation of the birefringence. However,
this method did not account for the presence of the MS signal and it was thus
not expected to yield accurate birefringence estimation in such cases.
Nonetheless, we observed a reduction in the skin birefringence in the deeper
regions of the medaka fish after applying the MFA method, which was consistent
with the appearance of the skin birefringence in the superficial regions. This
finding indicates that removal of the MS signal can improve the birefringence
estimation accuracy.

###  4.2 Artifactual elevation of birefringence in the low signal-to-noise
region

Although the SFA method cannot reduce the MS signals, the results for both the
phantom and the fish showed that the SFA method reduces the measured
birefringence values. This can be explained by the noise susceptibility of the
phase retardation and birefringence measurements of JM-OCT [29, 30]. The phase
retardation when measured by JM-OCT converges at approximately 2/3 $pi$ rad as
the effective SNR decreases. If the true phase retardation is low, this effect
results in an erroneous up-shift in the measured phase retardation. Because
the birefringence used in this study is defined as a local phase retardation,
the birefringence also suffered from this systemic error. Although the MAP
birefringence estimator used in this study can mitigate this effect, it
remains effectual at very low SNR [25, 26].

This artifactual elevation of the birefringence can explain the reduction in
the measured birefringence realized by the SFA method, because the SFA method
increases the SNR by complex averaging.

###  4.3 Limitations of the MFA implementation

The relationship between the MS behavior and the polarization artifacts is not
well understood at present. Theoretical modeling and/or numerical simulations
of the PS-OCT signal in an MS regime may provide a more comprehensive
understanding and further insights in future work.

Another limitation of the MFA method is the difficulty of performing in vivo
imaging. Because the MFA method requires multiple volumetric imaging, it is
vulnerable to the sample’s motion. In other words, the sample should remain
stationary during all the volumetric measurements, which have taken
approximately 20 s to complete in the present study.

One possible solution is to integrate the MFA method with a high-speed PS-OCT
system [31, 32, 33, 34, 35]. Another potential solution is to use B-scan-based
focus modulation, in which the depth position of the focus changes B-scan by
B-scan, which means that the sample only need to remain stationary for a few
B-scan acquisitions and not for a few volume acquisitions. In our preliminary
study, a set of B-scans with focus modulation was acquired in 86.2 ms [36]
using the same JM-OCT system.

Either of these two solutions, or a combination of them, may help achieve MS
signals reduction for in vivo imaging and enhance visualization of the
polarization features in deep tissues.

##  5 Conclusion

In this work, we demonstrate a technique called MFA to mitigate the MS signals
in Jones matrix PS-OCT. The proposed method was validated using a scattering
phantom and a postmortem medaka fish, and it was demonstrated that the MFA
method reduces the artifacts for both DOPU and birefringence effectively. As a
result, improved DOPU and birefringence image contrasts were obtained. The
proposed method is a cost-effective solution to improve the accuracy of
quantitative polarization measurements. This work may have important
implications for polarization-sensitive imaging and polarization properties
quantification of thick biological tissues.

### Acknowledgments

Lida Zhu is supported by the China Scholarship Council through the Chinese
Government Graduate Student Overseas Study Program.

## Funding

Core Research for Evolutional Science and Technology (JPMJCR2105); Japan
Society for the Promotion of Science (21H01836, 22K04962); China Scholarship
Council (201908130130).

### Disclosures

L. Zhu, S. Makita, Y. Lim, Y. Zhu, Y. Yasuno: Yokogawa Electric Corp. (F), Sky
Technology (F), Nikon (F), Kao Corp. (F), Topcon (F). J. Tamaoki, M.
Kobayashi: None.

### Data, Materials, and Code Availability

Data underlying the results presented in this paper are not publicly available
at this time but may be obtained from the authors upon reasonable request.

## References

  * [1] D. Huang, E. A. Swanson, C. P. Lin, J. S. Schuman, W. G. Stinson, W. Chang, M. R. H. a. Flotte, K. Gregory, C. A. Puliafito, and J. G. Fujimoto, “Optical coherence tomography,” Science 254, 1178–1181 (1991).
  * [2] J. F. de Boer, C. K. Hitzenberger, and Y. Yasuno, “Polarization sensitive optical coherence tomography – a review [Invited],” Biomed. Opt. Express 8, 1838–1873 (2017).
  * [3] B. Baumann, “Polarization sensitive optical coherence tomography: A review of technology and applications,” Appl. Sci. 7, 474 (2017).
  * [4] J. F. de Boer, T. E. Milner, M. J. C. van Gemert, and J. S. Nelson, “Two-dimensional birefringence imaging in biological tissue by polarization-sensitive optical coherence tomography,” Opt. Lett. 22, 934–936 (1997).
  * [5] S. K. Nadkarni, M. C. Pierce, B. H. Park, J. F. de Boer, P. Whittaker, B. E. Bouma, J. E. Bressner, E. Halpern, S. L. Houser, and G. J. Tearney, “Measurement of collagen and smooth muscle cell content in atherosclerotic plaques using polarization-sensitive optical coherence tomography,” J. Am. Coll. Cardiol. 49, 1474–1481 (2007).
  * [6] A. Miyazawa, M. Yamanari, S. Makita, M. Miura, K. Kawana, K. Iwaya, H. Goto, and Y. Yasuno, “Tissue discrimination in anterior eye using three optical parameters obtained by polarization sensitive optical coherence tomography,” Opt. Express 17, 17426–17440 (2009).
  * [7] Y. Lim, M. Yamanari, S. Fukuda, Y. Kaji, T. Kiuchi, M. Miura, T. Oshika, and Y. Yasuno, “Birefringence measurement of cornea and anterior segment by office-based polarization-sensitive optical coherence tomography,” Biomed. Opt. Express 2, 2392–2402 (2011).
  * [8] E. Götzinger, M. Pircher, W. Geitzenauer, C. Ahlers, B. Baumann, S. Michels, U. Schmidt-Erfurth, and C. K. Hitzenberger, “Retinal pigment epithelium segmentation by polarization sensitive optical coherence tomography,” Opt. Express 16, 16410–16422 (2008).
  * [9] M. Miura, S. Makita, Y. Yasuno, S. Azuma, T. Mino, T. Yamaguchi, T. Iwasaki, R. Nemoto, H. Shimizu, and H. Goto, “Objective evaluation of choroidal melanin loss in patients with Vogt–Koyanagi–Harada disease using polarization-sensitive optical coherence tomography,” Sci. Rep. 12, 3526 (2022).
  * [10] M. Miura, S. Makita, Y. Yasuno, T. Iwasaki, S. Azuma, T. Mino, and T. Yamaguchi, “Evaluation of choroidal melanin-containing tissue in healthy Japanese subjects by polarization-sensitive optical coherence tomography,” Sci. Rep. 12, 4048 (2022).
  * [11] D. Yang, Z. Yuan, M. Hu, and Y. Liang, “Zebrafish brain and skull imaging based on polarization-sensitive optical coherence tomography,” J. Biophotonics 15, e202200112 (2022).
  * [12] A. Lichtenegger, B. Baumann, and Y. Yasuno, “Optical coherence tomography is a promising tool for zebrafish-based research—a review,” Bioengineering 10, 5 (2022).
  * [13] A. F. Fercher, W. Drexler, C. K. Hitzenberger, and T. Lasser, “Optical coherence tomography - principles and applications,” Rep. Prog. Phys. 66, 239–303 (2003).
  * [14] A. Lichtenegger, P. Mukherjee, L. Zhu, R. Morishita, K. Tomita, D. Oida, K. Leskovar, I. A. El-Sadek, S. Makita, S. Kirchberger, M. Distel, B. Baumann, and Y. Yasuno, “Non-destructive characterization of adult zebrafish models using jones matrix optical coherence tomography,” Biomed. Opt. Express 13, 2202–2223 (2022).
  * [15] A. Badon, D. Li, G. Lerosey, A. C. Boccara, M. Fink, and A. Aubry, “Smart optical coherence tomography for ultra-deep imaging through highly scattering media,” Sci. Adv. 2, e1600370 (2016).
  * [16] M. Wojtkowski, P. Stremplewski, E. Auksorius, and D. Borycki, “Spatio-temporal optical coherence imaging — a new tool for in vivo microscopy,” Photonics Lett. Pol. 11, 44–49 (2019).
  * [17] S. Liu, M. R. E. Lamont, J. A. Mulligan, and S. G. Adie, “Aberration-diverse optical coherence tomography for suppression of multiple scattering and speckle,” Biomed. Opt. Express 9, 4919–4935 (2018).
  * [18] L. Zhu, S. Makita, J. Tamaoki, A. Lichtenegger, Y. Lim, Y. Zhu, M. Kobayashi, and Y. Yasuno, “Multi-focus averaging for multiple scattering suppression in optical coherence tomography,” Biomed. Opt. Express 14, 4828–4844 (2023).
  * [19] E. Li, S. Makita, Y.-J. Hong, D. Kasaragod, and Y. Yasuno, “Three-dimensional multi-contrast imaging of in vivo human skin by Jones matrix optical coherence tomography,” Biomed. Opt. Express 8, 1290–1305 (2017).
  * [20] A. Miyazawa, S. Makita, E. Li, K. Yamazaki, M. Kobayashi, S. Sakai, and Y. Yasuno, “Polarization-sensitive optical coherence elastography,” Biomed. Opt. Express 10, 5162–5181 (2019).
  * [21] K. Oikawa, D. Oida, S. Makita, and Y. Yasuno, “Bulk-phase-error correction for phase-sensitive signal processing of optical coherence tomography,” Biomed. Opt. Express 11, 5886–5902 (2020).
  * [22] L. Zhu, S. Makita, D. Oida, A. Miyazawa, K. Oikawa, P. Mukherjee, A. Lichtenegger, M. Distel, and Y. Yasuno, “Computational refocusing of jones matrix polarization-sensitive optical coherence tomography and investigation of defocus-induced polarization artifacts,” Biomed. Opt. Express 13, 2975–2994 (2022).
  * [23] D. Hillmann, T. Bonin, C. Lührs, G. Franke, M. Hagen-Eggert, P. Koch, and G. Hüttmann, “Common approach for compensation of axial motion artifacts in swept-source OCT and dispersion in Fourier-domain OCT,” Opt. Express 20, 6761–6776 (2012).
  * [24] S. Makita, Y.-J. Hong, M. Miura, and Y. Yasuno, “Degree of polarization uniformity with high noise immunity using polarization sensitive optical coherence tomography,” Opt. Lett. 39, 6783–6786 (2014).
  * [25] D. Kasaragod, S. Makita, S. Fukuda, S. Beheregaray, T. Oshika, and Y. Yasuno, “Bayesian maximum likelihood estimator of phase retardation for quantitative polarization-sensitive optical coherence tomography,” Opt. Express 22, 16472–16492 (2014).
  * [26] D. Kasaragod, S. Makita, Y.-J. Hong, and Y. Yasuno, “Noise stochastic corrected maximum a posteriori estimator for birefringence imaging using polarization-sensitive optical coherence tomography,” Biomed. Opt. Express 8, 653–669 (2017).
  * [27] D. A. Freedman and P. Diaconis, “On the histogram as a density estimator:l2 theory,” Z. Wahrscheinlichkeitstheorie verw. Gebiete 57, 453–476 (1981).
  * [28] S. G. Adie, T. R. Hillman, and D. D. Sampson, “Detection of multiple scattering in optical coherence tomography using the spatial distribution of Stokes vectors,” Opt. Express 15, 18033–18049 (2007).
  * [29] S. Makita, M. Yamanari, and Y. Yasuno, “Generalized Jones matrix optical coherence tomography: performance and local birefringence imaging,” Opt. Express 18, 854–876 (2010).
  * [30] L. Duan, S. Makita, M. Yamanari, Y. Lim, and Y. Yasuno, “Monte-Carlo-based phase retardation estimator for polarization sensitive optical coherence tomography,” Opt. Express 19, 16330–16345 (2011).
  * [31] F. A. South, E. J. Chaney, M. Marjanovic, S. G. Adie, and S. A. Boppart, “Differentiation of ex vivo human breast tissue using polarization-sensitive optical coherence tomography,” Biomed. Opt. Express 5, 3417–3426 (2014).
  * [32] T. Klein and R. Huber, “High-speed oct light sources and systems [Invited],” Biomed. Opt. Express 8, 828–859 (2017).
  * [33] J. P. Kolb, W. Draxinger, J. Klee, T. Pfeiffer, M. Eibl, T. Klein, W. Wieser, and R. Huber, “Live video rate volumetric oct imaging of the retina with multi-mhz a-scan rates,” PLoS One 14, 1–20 (2019).
  * [34] E. Auksorius, D. Borycki, P. Stremplewski, K. Liżewski, S. Tomczewski, P. Niedźwiedziuk, B. L. Sikorski, and M. Wojtkowski, “In vivo imaging of the human cornea with high-speed and high-resolution fourier-domain full-field optical coherence tomography,” Biomed. Opt. Express 11, 2849–2865 (2020).
  * [35] S. Chen, B. Potsaid, Y. Li, J. Lin, Y. Hwang, E. M. Moult, J. Zhang, D. Huang, and J. G. Fujimoto, “High speed, long range, deep penetration swept source OCT for structural and angiographic imaging of the anterior eye,” Sci. Rep. 12, 992 (2022).
  * [36] Y. Zhu, L. Zhu, S. Makita, Y. Lim, and Y. Yasuno, “B-scan-wise multi-focus averaging method to suppress multiple scattering signal,” in _Optical Coherence Tomography and Coherence Domain Optical Methods in Biomedicine XXVII,_ vol. 12367 J. A. Izatt and J. G. Fujimoto, eds. (SPIE, 2023), pp. 175–178.