Fast wavefront adaptive holography in Nd : YVO 4 for ultrasound optical tomography imaging

Several approaches exist to perform acousto-optic imaging of multiple-scattering media such as biological samples. Up to now, most of the coherent detection methods use holographic setup based on photorefractive crystals such as BSO or SPS. One of the issue of these techniques is the moderate response time compared to the speckle decorrelation time in biological sample. We introduce a new approach for the holographic detection based on two-wave mixing in a Nd:YVO4 gain medium enabling us to perform a fast wavefront adaption (50 μs) of the speckle field from a multiple-scattering sample. © 2014 Optical Society of America OCIS codes: (090.2880) Holographic interferometry; (110.0113) Imaging through turbid media; (110.7170) Ultrasound; (170.3660) Light propagation in tissues; (170.3880) Medical and biological imaging; (190.7070) Two-wave mixing. References and links 1. D. S. Elson, R. Li, C. Dunsby, R. Eckersley, and M.-X. Tang, “Ultrasound-mediated optical tomography: a review of current methods,” Interface Focus 1, 632–648 (2011). 2. M. Kempe, M. Larionov, D. Zaslavsky, and A. Z. Genack, “Acousto-optic tomography with multiply scattered light,” J. Opt. Soc. Am. A 14, 1151–1158 (1997). 3. G. Rousseau, A. Blouin, and J.-P. Monchalin, “Ultrasound-modulated optical imaging using a high-power pulsed laser and a double-pass confocal fabry-perot interferometer,” Opt. Lett. 34, 3445–3447 (2009). 4. X. Xu, S.-R. Kothapalli, H. Liu, and L. V. Wang, “Spectral hole burning for ultrasound-modulated optical tomography of thick tissue,” J. Biomed. Opt 15, 066018 (2010). 5. H. Zhang, M. Sabooni, L. Rippe, C. Kim, S. Kroll, L. V. Wang, and P. R. Hemmer, “Slow light for deep tissue imaging with ultrasound modulation,” Applied Physics Letters 100, 131102 (2012). 6. S. Farahi, G. Montemezzani, A. A. Grabar, J.-P. Huignard, and F. Ramaz, “Photorefractive acousto-optic imaging in thick scattering media at 790 nm with a Sn2P2S6:Te crystal,” Opt. Lett. 35, 1798–1800 (2010). 7. P. Lai, X. Xu, H. Liu, and L. V. Wang, “Time-reversed ultrasonically encoded optical focusing in biological tissue,” J. Biomed. Opt 17, 030506 (2012). 8. M. Gross, P. Goy, B. C. Forget, M. Atlan, F. Ramaz, A. C. Boccara, and A. K. Dunn, “Heterodyne detection of multiply scattered monochromatic light with a multipixel detector,” Opt. Lett. 30, 1357–1359 (2005). 9. Y. Liu, P. Lai, C. Ma, X. Xu, Y. Suzuki, A. A. Grabar, and L. V. Wang, “High-speed time-reversed ultrasonically encoded (TRUE) optical focusing inside dynamic scattering media at 793 nm,” (2014), vol. 8943, pp. 894339. 10. A. E. Siegman, Lasers (University Science Books, Mill Valley, CA, 1986). 11. A. Brignon and J.-P. Huignard, “Two-wave mixing in Nd:YAG by gain saturation,” Opt. Lett. 18, 1639–1641 (1993). 12. R. Soulard, A. Brignon, J. P. Huignard, and R. Moncorge, “Non-degenerate near-resonant two-wave mixing in diode pumped Nd3+ and Yb3+ doped crystals in the presence of athermal refractive index grating,” J. Opt. Soc. Am. B 27, 2203–2210 (2010). #214650 $15.00 USD Received 24 Jun 2014; revised 4 Aug 2014; accepted 6 Aug 2014; published 19 Aug 2014 (C) 2014 OSA 25 August 2014 | Vol. 22, No. 17 | DOI:10.1364/OE.22.020622 | OPTICS EXPRESS 20622 13. P. Sillard, A. Brignon, and J.-P. Huignard, “Loop resonators with self-pumped phase-conjugate mirrors in solidstate saturable amplifiers,” J. Opt. Soc. Am. B 14, 2049–2058 (1997). 14. M. J. Damzen, A. Boyle, and A. Minassian, “Adaptive gain interferometry: a new mechanism for optical metrology with speckle beams,” Opt. Lett. 30, 2230–2232 (2005). 15. B. Jayet, J.-P. Huignard, and F. Ramaz, “Optical phase conjugation in Nd:YVO4 for acousto-optic detection in scattering media,” Opt. Lett. 38, 1256–1258 (2013). 16. P. Yeh, “Two-wave mixing in nonlinear media,” IEEE Journal of Quantum Electronics 25, 484–519 (1989). 17. M. J. Damzen, M. Trew, E. Rosas, and G. J. Crofts, “Continuous-wave Nd:YVO4 grazing-incidence laser with 22.5 w output power and 64% conversion efficiency,” Opt. Commun. 196, 237–241 (2001). 18. M. Lesaffre, F. Jean, F. Ramaz, A. C. Boccara, M. Gross, P. Delaye, and G. Roosen, “In situ monitoring of the photorefractive response time in a self-adaptive wavefront holography setup developed for acousto-optic imaging,” Opt. Express 15, 1030–1042 (2007).


Introduction
Tagging photons with ultrasound is an attractive method to image the optical properties inside multiple-scattering media such as biological tissues [1].Indeed, conventional optical imaging techniques cannot be used on biological tissues because they rely on ballistic light, which is attenuated very fast in such media due to a very high scattering coefficient.A way to harvest information from the diffused light is to couple it with ultrasound (US) through the acousto-optic effect.As the US at some MHz are ballistic inside most biological tissues, it is possible to focus them in a small volume thus creating a small virtual source of ultrasonically tagged photons.These so called tagged photons are slightly shifted in frequency, due to the US phase modulation.Then, with an appropriate detection, it is possible to probe the local optical properties of the sample by scanning the focus spot of the US within the sample.As an increase of the tagged photons results in a decrease of the untagged photons and vice versa, the information from the acousto-optic signal can be retrieved on either component.The resolution of this acousto-optic imaging technique, given by the size of the virtual tagged photons source, is connected to the size of the acoustic focal spot (typically 1 mm at echographic frequencies), whereas the contrast of the image is given by the intensity of light at the US focal point which means that it is of optical nature.
If a single detector is used at the exit of a multiple-scattering sample, the signal to noise ratio (SNR) for the detection of the modulation created by the US is weak due to first, the spatial random phase of each speckle grain [2] and second, the small fraction of photons tagged by the acoustic field.Several methods exist in order to improve this SNR, they can be either incoherent or coherent.The incoherent methods are based on creating a very sharp spectral filter to separate the tagged from the untagged photons and do only a flux detection.They are usually based on Fabry-Perot cavities [3] or spectral hole burning [4,5].The coherent methods, on the other hand, are based on the idea of phasing each speckle grain of one component (tagged or untagged) of the speckle field in order to increase the SNR.This can be performed by creating a new wavefront that is matched with, for example, the wavefront of the untagged photons.All the speckle grains of both wavefronts will interfere coherently on the single detector resulting in a higher SNR.The wavefront adaption can be realized with dynamic two-wave mixing (TWM) holography in non linear photosensitive media.
Usually photorefractive crystals are used to perform the TWM.Photorefractive dynamic holography has enabled the detection of absorbing inclusions embedded inside scattering phantoms [6] and has also been used to perform phase conjugation of the ultrasonically tagged photons in order to focus back behind a scattering sample [7].One of the issue that can come up with in vivo experiments is that the speckle has a decorrelation time constant of about 1 ms [8,9] in biological tissues, while the required time to write the hologram inside a photorefractive crystal is usually between a few ms to 100 ms [6,7].To overcome this problem, we propose in this article a new approach of TWM based on a gain media (Nd:YVO 4 ).The response time in vanadate crystal is governed by the spontaenous emission lifetime which is typically 90 µs.Nevertheless, this lifetime, and also the response time can be reduced under certain conditions such as strong pumping or saturation effects [10].
Wave mixing in gain media has already been theoretically and experimentally investigated [11,12].It has also been used to perform phase conjugate laser cavities [13], vibrometry of rough surfaces [14] and detection of acousto-optic signals from a multiple-scattering sample by four-wave mixing (FWM, phase conjugation) [15].Here we report our investigations on the use of gain media for the detection of the acousto-optic signal by TWM.
At first, we will introduce a model of the interactions of two coherent plane waves interfering inside a diode pumped Nd:YVO 4 crystal.This model will then be compared to experimental measurements.We will also investigate the experimental response time of the holography in Nd:YVO 4 in our setup.Finally, the first experimental results of imaging a scattering phantom with 3 black inclusions will be presented before discussing the outcomes and the perspectives of the use of gain media for this imaging process.

Model for TWM in Nd:YVO 4
Gain media are traditionally used either as amplifiers or in laser cavities.Light amplification in these media is the result of a diode pumping and stimulated emission at a higher wavelength.In Nd:YVO 4 , Nd 3+ ions are responsible for the amplification.A simplified diagram of their energy levels (Fig. 1(a)) shows that they absorb light around 808 nm and re-emit photons at 1064 nm by stimulated emission.If the pumping is strong enough, it is possible to create a population inversion between the two levels of the laser transition (level 1 and 2 in Fig. 1(a)).The intensity gain per unit of length of the crystal, α(r,t), defined by the following relation I OUT = I IN e αL (where L is the length of the gain region), is proportional to the population inversion ∆N = N 2 − N 1 created by the pumping.As indicated in Fig. 1(a), the dynamic of these processes is mainly governed by the lifetime τ of the transition 2 → 1, which is slow compared to the other transitions.In the presence of an incident optical field, the dynamic of α(r,t) is given by [14] : where α 0 is the small signal intensity gain, τ the spontaneous emission lifetime of the laser transition and S(r,t) = 1 + I(r,t) I S is the saturation parameter of the amplifier where I S is the saturation intensity of the gain medium.As the saturation parameter depends on the local light intensity (in W/cm 2 ), it is possible to spatially modulate the gain using an interference pattern and thus to record the spatial phase information of a complex signal beam as shown in Fig. 1(b).Let us consider the interference of two coherent beams (A R and A S ) inside the crystal.The fringes intensity pattern can be written: I(r,t) = I R + I S + A R A * S e jK.r + c.c where I R and I S are the intensities of the optical fields A R and A S and K = k r − k s is the interference wave vector.In order to solve Eq. ( 1), we expanded α(r,t) as a Fourier series in spatial harmonics and kept only the first two terms: α(r,t) = α (0) + 1 2 α (1) e jK.r + c.c., where α (0) is the uniform gain (the gain without the spatial intensity modulation) and α (1) is the grating contribution to the gain.Under these conditions, two coupled equations for α (0) and α (1) can be derived from Eq. ( 1) by only identifying the terms for K = 0 and K.r (with terms in 2K.r neglected): ∂ α (1)  ∂t If we assume that the input signal field A S is small compared to the reference field A R and that the reference field is either small or of the order of magnitude of the saturation intensity (I S ≪ I R and A S A R < I sat ), Eqs. ( 2) and ( 3) can be simplified: In order to complete the model, we have to consider the evolution of the amplitude of the two optical fields within the crystal.The system of coupled equations for A S and A R in the crystal can be derived from Maxwell's equations [16].If we consider a light propagation along the z-axis we have (with the slow varying envelope approximation): Equations ( 4)-( 7) can be simplified by taking into account our experimental conditions which consider that: • The intensity of the signal beam (I S ) is weak compared to the intensity of the reference (I R ) and to the saturation intensity (I sat ).
• The reference beam is "ON" and is also established dI R dt = 0 .
Under these assumptions, Eqs. ( 4)-( 7) can be simplified into the following system of coupled equations that can be numerically solved using, for example, a finite difference methods (Euler): ∂ α (1) (z,t) First, Eqs. ( 8) and ( 9) can be considered independently of Eqs. ( 10) and ( 11) in order to calculate I R (z) and α (0) (z).Then, Eqs. ( 10) and ( 11) can be solved step by step using a spatio-temporal Euler's method.Another method of solving Eqs. ( 10) and ( 11) is by using the Laplace transform of each equations.This leads to an analytical formula for the amplitude of the signal field at the output of the crystal given by Eq. ( 12).The details of the calculation are shown in the appendix.
where f * g is an "incomplete" convolution operator of function f and g defined as dτ and C(z,t) is the transfer function of the crystal on the input signal beam.It is given by: where Γ is the usual Gamma function.Equation ( 12) is interesting to understand the effects of the crystal amplifier on the signal.The first term, exp 1 2 z 0 α (0) (z ′ )dz ′ , exhibits the input signal amplification in the crystal and the second term is the time dependent component of the output signal, it is composed of: • The original signal A S (0,t).
• The diffracted signal A S (0,t) * C(z,t) by the gain grating.This last term represents the input signal filtered by the crystal with a transfer function C(z,t), which expression is given in Eq. ( 13).This transfer function is a first order transfer function with corrective terms given in the summation.The response time of this filter is, at the first as expected from Eqs. ( 4) and ( 5).This means that uniform gain (α (0) ) and modulated gain (α (1) ) almost exhibit the same temporal behavior, shorter than the fluorescence lifetime if the reference beam is intense enough to saturate the gain.This makes the wavefront adaption in vanadate crystals really fast.If only the first term of the summation in C(z,t) is considered, the formula is consistent with the formalism given by Damzen et al. [14] in the case of an established modulation.

Experimental conditions
This section describes our experimental conditions.The crystal is a 5 x 5 x 10 mm 3 YVO 4 crystal doped at 1%-at with Nd 3+ ions responsible for the amplification, and corresponding to a saturation intensity around 1 kW/cm 2 (manufactured by Crylight Photonics).One of the large side is anti-reflection (AR) coated at 808 nm for the pump and the two small sides are AR coated at 1064 nm for the signal and reference beams to enter the crystal in the geometrical configuration presented in Fig. 2(a).It is pumped by a 808 nm laser diode (LIMO100-F200-DL808), the output of this diode is a multimode fiber collimated, then focused inside the crystal by a 18 mm focal length cylindrical lens (CL on Fig. 1(a)).The pump works in a quasi continuous wave (QCW) regime, it consists of 2.5 ms long pulses at a repetition rate of 40 Hz.The peak power of each pulse is 80 W, the QCW regime is used to achieve a high gain and prevent thermal damages from the pumping.The crystal is cooled in an aluminium structure with a peltier module.The signal and reference beams are generated by the same laser source, a 1W-CW single longitudinal mode Nd:YAG laser (CrystaLaser inc.), coupled to a 5W Yb-doped fiber amplifier (Keopsys inc.).Each arm of the setup includes an acousto-optic modulator (AA Optoelectronic) that is used either as a fast shutter (for the experiment in part 4) or to generate a frequency shift between the signal and the reference (for the experiment in part 5).The beams enter the crystal with a small incidence angle (10 • for the reference and 20 • for the signal) to obtain a total internal reflection [17] on the pumped surface and ensure a good overlap of the interference pattern and the gain volume, which has a typical size of 0.3 x 0.3 x 5 mm 3 .At the output of the crystal, the signal beam is collected on a fast photodiode (Thorlabs DET10A) amplified by a large bandwidth amplifier (FEMTO DHPVA) and sampled at 40MHz using an ADLINK-PCI9646D DAQ board in the computer driven by MatLab r .

Model versus experiment
In order to verify the model, we performed measurements of the temporal evolution of the gain on a signal beam when it is switched ON in the presence of a continuous reference beam.A 1.5 ms long signal pulse is switched ON with a delay of δt = 500µs after the pump pulse to have an established population inversion (Fig. 2(b)).In order to observe the temporal behavior of the gain I Sout /I Sin and compare it with the model, we measured the signal both with (I Sout ) and without (I Sin ) a pump pulse.The signal with pump I Sout is corrected of the amplified spontaneous emission (ASE), which is the main source of noise, by recording the signal on the photodiode between t = 0 and t = δt.Solving Eqs. ( 8)-( 11) requires the values of the small signal gain α 0 , the ratio between the incoming reference intensity and saturation intensity I R0 /I sat and the response time τ R .The first two are estimated using the experimental measure- ments whereas τ R is adjusted to best fit the experimental data.Figure 3 shows the experimental data (dots) and the calculated gain (solid line) for several value of I R0 /I sat .In order to have a good agreement between the experimental and the calculated data, we had to introduce two corrections.First the small signal gain α 0 , which was estimated by measuring the gain without the reference beam, had to be replaced by an effective gain equal to 0.75 α 0 .After this gain correction, there still was an offset difference between the calculated data and the experimental data which has been corrected to fit the temporal behavior.As shown in Fig. 3, after these corrections, the model predicts justly the temporal evolution of the gain hologram.Several observations can be made on this figure.Firstly, the effects of the saturation can be observed, the initial gain on the signal beam decreases as the intensity of the reference beam increases.Secondly, the effect of the wave mixing and the writing of the hologram can be seen in the temporal behavior.Over time, the gain on the signal beam decreases to reach a stationary value when the static hologram is recorded inside the medium.The calculation has been made with a value of τ R = 50 µs which gives a good estimation of the temporal dynamic of the holographic recording.The magnitude of τ R will be confirmed in section 5 where we experimentally measure the response time of our setup by using a method proposed by Lesaffre et al. [18].
The corrections on the calculated data to fit the experimental signals can be justified by the different conditions between the model and the experiment.The model supposes a perfect overlap between the signal, the reference and the pump beams and also a uniform pumping inside the crystal; whereas in reality, the signal and reference beams are not parallel, and also the pump is not uniform as it is absorbed by the crystal as described in section 3. Nevertheless, after introduction of correcting factors, this model shows a good agreement with the experimental data.

Dynamic of the TWM in Nd:YVO 4
In this part we experimentally investigate the response time by using a method proposed in [18].It allows to directly measure in situ the response time of a self-adaptive wavefront holography setup using the scattering sample and US.We used the method explained in part 4, we sent on our sample an amplitude modulated ultrasonic excitation and recorded the temporal signal at the exit of the holographic system for several frequency shifts between the signal and the reference beams.Then the amplitude of the signal at the modulation frequency is retrieved by detecting the peak on the Fourier spectrum.The method requires the modulation frequency to be large compared to the inverse of the response time τ R , as we know that τ R is about 50 µs (τ −1 R = 20 kHz) we chose to do the experiment with a modulation frequency f mod = 200kHz.The US are generated by a 25.4 mm diameter Olympus focused single element transducer (model A307S) with a focal length of 51.5 mm and a center frequency f US = 4.8MHz.Figure 4 shows the amplitude of the peak at f mod = 200kHz versus the frequency shift (∆ f ) between the signal and the reference beams.In order to estimate τ R , the experimental points have been fitted by a Lorentzian model: The fit gives τ R = 44 µs ± 6 µs with a coefficient of determination r 2 = 0.91.This result confirms the previous observations with the model that lead to a response time of 50 µs.

Acousto-optic imaging using holography in Nd:YVO 4
This holographic detection is used to perform acousto-optic imaging of a multiple-scattering sample.In order to test this detection method, we put a sample into a water tank on the signal beam.An optical system (OS1, see Fig. 5(a)) is used to collect the scattered light and focus it inside the crystal within the pumped volume.After the crystal, light is collected by OS2 and focused on the photodiode (Thorlabs DET10A).A simplified drawing of the setup is shown in Fig. 5(a).In this setup, we chose to detect the untagged photons [6] thus the reference beam is not shifted from the US frequency.The pumping conditions are the same as the one described in section 3. The sample used is a 1.5 cm thick scattering phantom made of agar and zinc oxide particle for scattering (µ * s = 3 cm −1 ).In this gel we put three circular black inclusions made with Indian ink and agar with a diameter of 1 mm, 2 mm and 3 mm.A picture of the inclusions before being covered with a layer of scattering gel is shown in Fig. 5(b).The sample is placed inside a water tank (not shown in Fig. 5(a)) and placed under the same ultrasonic transducer as the one used in part 5.It sends short bursts (4 cycles) at f US = 4.8 MHz giving a resolution of 650 µm laterally and 1.2 mm axially.The transducer is aligned on the x axis in order to be al- ways above the laser and on the y axis so that its focal spot is at the same height as the incoming laser on the sample.To get a large image, the sample is moved along the x and y axes by two PI-Micos Pollux Drive translation stages (not shown in Fig. 5 Normalized AO signal (a.u.) Normalized AO signal (a.u.) 1.9 mm 0.9 mm 2.7 mm The three inclusions the gel are clearly visible on the acousto-optic image.The full width at half maximum (FWHM) of the inclusions measured on the profile are respectively 0.9 mm, 1.9 mm and 2.7 mm.These values are in agreement with their real size.

Conclusion
As outlined in the article, adaptive wavefront holography is a fast process in gain media such as Nd:YVO 4 .The response time is a main advantage for ultrasound optical tomography imaging of biological tissues due to the problems of speckle decorrelation that arises when performing in vivo experiments.We have shown through a model and measurements that our setup, with a neodymium doped vanadate crystal, has a response time around 50 µs.This value was confirmed by a method for in situ monitoring of the response time in a self-adaptive wavefront holography experiment.We then used our setup to perform the detection of an acousto-optic signal from a multiple-scattering phantom in order to map the local optical properties inside a sample.We were able to detect inclusions as small as one mm using US at 4.8 MHz.The acquisition time is still too long because of the need for a lot of averaging due to a low SNR.
For future experiments, improving the SNR is a necessity to make data acquisition faster (less averaging) on samples with scattering properties closer to biological tissue (µ * s = 10 cm −1 ).These experiments show that gain media such as Nd:YVO 4 are attractive candidates for interferometry with speckle beams with application such as in vivo UOT, where the speckle decorrelates rapidly, or for metrology of moving or vibrating rough surfaces.Even if we focused on the use of two-wave mixing for acousto-optic imaging, other application based on four-wave mixing (phase conjugation) can also be considered as we demonstrated last year (see [15]).For example, a technique such as TRUE focusing could benefit from the speed of gain media to overcome the speckle decorrelation problems inherent to in vivo experiments.The working wavelength, 1064 nm, is at the edge of the optical therapeutic windows so it does not penetrate very deep inside biological tissues in comparison to light at 800 nm.Therefore, extending these experiments to gain media operating in the 700 nm -900 nm part of the spectrum, such as Ti:Sapphire, Alexandrite or Cr:LiSAF, would be attractive for further in vivo biomedical applications.

Appendix : Analytical solution for the signal exiting the gain medium
This appendix details the calculation to go from the system of Eqs.(10) and (11) to the analytical solution given in Eq. ( 12).As mentioned in section 2, if we take the Laplace Transform of Eqs.(10) and (11), we obtain: Considering that at t = 0 there is no interference pattern, we can write α (1) (z, 0) = 0 so Eq. ( 14) can be simplified: Then inserting Eq. ( 16) into Eq.( 15) yields: Which can be rewritten: Equation (18) can then be integrated between 0 and z: A 1 (z, s) = K 2 (s) exp 1 2 The boundary conditions for z = 0 gives: A 1 (0, s) = K 2 (s) s + 1 τ 1 + I R (0) The temporal evolution of the signal field is given by the inverse Laplace transform of Eq. ( 22).In this formula, the only term that has a dependency in the variable s is Q.This means that we have to find the inverse Laplace transform of 1/Q n .For n ≥ 1: And for n = 0: where δ is the Dirac function We can now take the inverse Laplace transform of Eq. ( 22) and have the temporal evolution of the signal beam: Where f * g is the "incomplete" convolution operator define in section 2. This equation can be rearranged:

Fig. 1 .
Fig. 1.(a) 4 levels energy model of Nd 3+ ions with time constants.(b)Gain modulation with a transverse diode pumping when two beams A R and A S are interfering inside the crystal.

Fig. 2 .
Fig. 2. (a)Geometric configuration of the signal and reference beams inside the crystal.(b) Temporal diagram of the pump, signal and reference beams.

Fig. 3 .
Fig.3.Experimental (dots) and calculated (solid lines) gain of a small signal beam for several value of the ratio between the intensity of the reference and saturation intensity

Fig. 4 .
Fig. 4. Amplitude of the component of the acousto-optic signal at f mod = 200kHz vs. the frequency shift between the signal and the reference beam (dots) and lorentzian fit (solid line).
inclusions inside the gel before being covered with another layer of gel a)Simplified setup for the acousto-optic imaging of a sample.OSi are optical system to collect light, HWi are half-wave plate to change the ratio between the intensities of the reference and signal beams and correct the polarization of each beams, AOM is an acoustooptic modulator and US is an ultrasonic transducer.(b)Drawing of the sample the inclusions (left) and picture (right) of the sample before being covered with another layer of scattering gel.
(a)) while the incoming laser and the US stay aligned.The displacement step along the x and y directions are respectively 0.1 mm and 0.75 mm.The photodiode signal is first amplified by a large bandwidth amplifier (FEMTOx