Keywords: acousto-optic imaging, ultrasound, scattering media, photorefractive holography, medical imaging, inverse problems. Résumé imagerie acousto-optique, ultrasons, milieux diffusants, holographie photoréfractive, imagerie médicale, problèmes inverses, PET: Positron Emission Tomography, PRC: Photorefractive Crystals, RTE: Radiation Transfer Equation, SNR: Signal-to-Noise Ratio, TWM: Two-wave Mixing, US: Ultrasound Acousto-optic imaging, ultrasound, scattering media, photorefractive holography, medical imaging, inverse problems

Image reconstruction for acousto-optics: towards quantitative imaging

Medical imaging is a crucial tool for medical diagnosis. Current imaging techniques such as magnetic resonance imaging or ultrasound scans are efficient and sensitive. However, to detect some pathologies such as early tumours, other methods are needed. Since those tumours strongly absorb light, optical methods would be suited to detect them but in biological tissues light is highly scattered. This restricts the imaging depth of classical optical techniques to a few millimetres in such media.

Acousto-optic imaging is a technique that couples light and ultrasound to obtain optical images deep inside scattering samples. As the ultrasound propagates in the sample, it modulates light and creates tagged photons. The number of these tagged photons is proportional to the local light fluence at the focus of the ultrasound. This method can thus map the light fluence up to a few centimetres deep inside scattering media by scanning the ultrasound. The critical step is the detection of these photons and several approaches can be used such as ultra narrow band filtering or self-adaptive interferometry.

In this thesis we present an acousto-optic setup using photorefractive crystals to perform selfadaptive wave-front holography for the detection of tagged photons. Several crystals can be used and we present here an experimental characterisation of two crystals: SPS and ZnTe. The imaging setup was also improved to operate in a reflection configuration, with source and detector on the same side of the imaged phantom and some experimental results are presented in this manuscript. We also propose to combine this experiment with a model-based reconstruction algorithm to reconstruct maps of the absorption coefficient from experimental data in the pursuit of quantitative acousto-optic imaging.

Aknowledgments

I would like to thank all the people who trusted and helped me carry out this PhD from the beginning to the end. First of all thanks to the reviewers, Stephen Morgan and Jean-Marc Dinten, who accepted to read my manuscript and made very relevant comments that helped me prepare the defense. I would also like to thank all the jury members for hearing my defense and for the very interesting discussion that followed. Thanks also to Philippe Adam, from the DGA, who followed my PhD from the beginning and accepted to come to my defense as an invited member of the jury. I also would like to aknowledge the Direction Generale de l'Armement for accepting to fund my PhD in collaboration with the DSTL. This was a great opportunity for me to work both in France and in the UK.

I would now like to thank my three supervisors, François, Terence and Samuel who trusted me and without whom nothing would have been possible. Thank you, François, for those countless hours spent in the lab helping me whenever something was weird with the setup. Thank you for the trust you put in me and the freedom you gave me to explore the things I wanted to explore but most of all, thank you for always being there when I struggled with the experiment or when I lost confidence in myself. This was the perfect dosage of supervision I needed and I am glad I chose to do this PhD with you. Thanks also for your moral integrity all the discussions about music and bass guitar. Samuel, I'd like to thank you for everything you taught me about inverse problems and numerical reconstruction and for your infinite patience when I came to you with the same questions over and over. It was a slow process for me but I learned a lot during this PhD. Thanks also for all the moral support when I was feeling down and pessimistic and for the discussions unrelated to work. It was a pleasure to work with you. Finaly, thank you Terence for accepting to supervise my work and for being a front between me and UCL administration. Thanks also for your advice.

I would also like to thank every person that has been part of the acousto-optic group in Institut Langevin. First Jean-Pierre Huignard, who has been a great help for everything related to photorefractive crystals and whose contribution is greatly appreciated. Secondly, all the PhD students who laid the stones on which my PhD was build and especially Baptiste Jayet, Emilie Benoît and Jean-Baptiste Laudereau without whom I would not have able to produce such work. Thanks Baptiste and JB for always having been there to answer my dumb questions even after you were long gone from the lab. Thanks for teaching me everything you knew about acousto-optics and for the strong moral support. Thanks also to the different post-docs who worked with me, Kevin Contreras and Maimouna Bocoum.

Many thanks also to all the collaborators and the people from Institut Langevin who helped me with the experimental work. I am thinking about Jean-Luc Genisson, who was our ultrasound expert; Vincent Servois who offered us the possibility to work in Institut Curie to try some in vivo experiment on mice; Marc Gesnik, for all his help regarding the programing of Aixplorer machines; Jean-Michel Tualle with whom I tried different experiments in the lab; Simon Arridge who gave me precious advice about inverse problems and countless others who gave a bit of much appreciated help now and then.

The success of a PhD also depends a lot on the environment and on the atmosphere in the lab. For that, I'm glad to have done my PhD in Institut Langevin and in Medical Physics for both labs are home to a fantastic crowd and this contributed a lot to my sanity. Many thanks to all the former or present members of the R31 office for all the discussions, games, drinks and more general moral support that they offered me during my PhD. Thanks also to all the PhD students and post doc in Institut Langevin for all the coffee breaks and interesting talks. Thanks also to the BORL buddies from UCL and for all those tea breaks and silly debates.

I will finish by thanking all my friends and familly who supported me during those three years and especially to my brother who has always been there for me. My thoughts also go to my late mother who I'm sure would have been proud of what I accomplished.

List of abbreviations

• AO: Acousto-Optic • MOPA: Master Oscillator Power Amplifier

• NIR: Near Infra-Red

• NIRS: Near Infra-Red Spectroscopy

• OCT: Optical Coherence Tomography

• PAI: PhotoAcoustic Imaging

Introduction

Medical diagnosis relies on the information that a physician can obtain about the body. This information can take several forms such as description of the symptoms like sensations or pain, chemical analysis of blood samples or urine for example, or medical imaging. The latter is of prime importance as it can be used in several stages of a medical process from diagnosis, to surgery and therapy monitoring. Medical imaging can provide knowledge about the structure and metabolism of a living body and can be used to better understand the different processes involved.

an optical resolution while NIRS provides spectral information, in this case about blood to assess its oxygenation. On the other hand, these techniques also highlight the difficulties of using light in scattering media. The imaging depth of microscopy is limited to a few millimetres whereas NIRS gives spatially unresolved measurements of an optical property at a few centimetre depth. One of the options to overcome these limitations is to couple light to another type of wave such as ultrasound. In addition to the fact that it is widely used, inexpensive and non-ionising, ultrasound has the advantage that it is weakly scattered by biological tissues. Two proposed techniques combining optical and acoustic waves are photoacoustics and acousto-optics, and both of them use the ultrasound to localise the optical information coming from deep areas. In photoacoustics, the sample is illuminated by a pulsed laser and light is absorbed by the tissues. The absorbers undergo thermal expansion and generate an acoustic wave that can be detected on the boundary of the sample. On the other hand, in acousto-optics, both waves are sent inside the medium. As the ultrasound propagates through an illuminated area, it locally modulates light through the acousto-optic effect and generates tagged photons. The detection of these tagged photons provides information about the local optical fluence with an ultrasound resolution.

The detection of the tagged photons is the main challenge of acousto-optic imaging as only a limited number of photons is tagged and the modulation corresponds to an extremely small frequency shift compared to the light frequency. This issue has been studied since the mid 90s and several detection techniques have been developed such as ultra narrow band spectral filtering or self-adaptive interferometry. Recent works at the Institut Langevin led to the development of a multi-modal platform combining a commercial ultrasound scanner to an acousto-optic imaging system. The detection of tagged photons is performed by an adaptive wave-front holography technique in photorefractive crystals. This platform already produced interesting ex vivo images and is the first step towards clinical imaging. Independently in the Medical Physics and Bioengineering department of University College London (UCL), a numerical method was developed to quantitatively reconstruct maps of the absorption and scattering coefficients from simulated acousto-optic images. The combination of these two approaches could be the first step towards quantitative acousto-optic imaging and is the subject of the work presented here. This manuscript is composed of five different chapters. The first chapter provides the theoretical bases needed to understand light propagation in highly scattering media and aims to show the difficulties of imaging such media with light. A quick overview of the existing optical imaging techniques will be presented as well as orders of magnitude of the optical properties of biological tissues.

The second chapter introduces the principle of acousto-optic imaging and the associated theoretical background. We will see that the measurement of the tagged photons flux is proportional to the local optical fluence and the quadratic pressure of the ultrasound. After explaining the formation of an image, an overview of the different tagged photons detection techniques will be presented.

The third chapter will focus on one such detection methods, self adaptive interferometry in photorefractive crystals. After a description of the photorefractive effect, we will see how it can be used for the detection of the modulation of the light by the ultrasound. The configuration of the Institut Langevin imaging system will be presented as well as experimental characterisation of two photorefractive crystals.

In the fourth chapter, we will see how Monte-Carlo algorithms can be used to simulate the acousto-optic signal and investigate an imaging configuration where the light source and detector are on the same side of the sample. The numerical results will be compared to experimental images obtained with the aforementioned acousto-optic imaging system. The fifth and last chapter concerns the application of the numerical reconstruction developed in UCL to experimental data obtained in Institut Langevin. After briefly describing the theory of inverse problems and their application to imaging, the algorithm for acousto-optic reconstruction will be presented and we will see that some model corrections are needed to accommodate the experimental data. The results of the reconstruction will be presented and discussed and the chapter will end with considerations on the next steps for quantitative acousto-imaging.

CHAPTER 1

Light propagation in thick scattering media

Of the interest of optical imaging

Biological tissues are made of complex formations (cells, fibre, membranes...) with a wide variety of structure and composition. To establish a diagnosis, a physician will often need to acquire information about both metabolism, through the chemical composition, and structure. Most medical imaging modalities are sensitive to only one of the two.

Light can interact with biological media through two processes: absorption and scattering. The first is influenced by the chemical properties of matter while the latter is caused by structural changes in the medium. Consequently, optical imaging can give a lot of relevant information to establish a diagnosis. Unfortunately, because biological tissues strongly scatter light, optical imaging is not always straightforward, especially at large depths.

In this chapter we will describe the properties of light-matter interaction, the different regimes of light propagation in tissues and we will discuss the different optical techniques which can be used for medical imaging.

Of the interest of optical imaging

The wide variety of medical imaging techniques developed during the 20 th century gives a lot of information on the biological tissue of interest. However, most of these techniques only provide structural information or metabolic information, rarely both.

On one hand, we have structural technique such as CT-scans (X-rays), which gives information about tissues density, or conventional ultrasound imaging which is sensitive to mechanical contrast through changes in acoustic impedance [START_REF] Analoui | Medical imaging: principles and practices[END_REF]. On the other hand, techniques such as PET (Positron Emission Tomography) have very specific contrast that give accurate functional information [2]. However, PET scans need isotopes as contrast agents, and usually give poor information about structure in tissues. The only technique which can combine structural and functional imagining is MRI (Magnetic Resonance Imaging) [START_REF] Paul | Image Formation by Induced Local Interactions: Examples Employing Nuclear Magnetic Resonance[END_REF][START_REF] Luyten | Metabolic imaging of patients with intracranial tumors: H-1 MR spectroscopic imaging and PET[END_REF]. Yet MRI has the disadvantage of being expensive and quite slow. Recent advances in ultrasound imaging lead to efficient Doppler imaging which can image blood flows [5] and gave promising first result in functional imaging [START_REF] Macé | Functional ultrasound imaging of the brain[END_REF] which enhances the potential of ultrasound imaging.

Optical imaging has the advantage to be sensitive to both structural and metabolic changes in tissues. Since it is also non invasive and non ionising and it does not always necessitate contrast agents, it can be a very powerful tool for medical imaging. However, if the interaction of light with tissue through absorption and scattering provides sensitivity, it can also reduce our capability to perform in-depth imaging.

The most straightforward way to perform optical imaging in biological tissues is to place a powerful light behind a patient's organ and to measure the transmitted light. The first optical technique using this principle, diaphanography [START_REF] Marshall | Diaphanography as a means of detecting breast cancer[END_REF], based on the idea developed by Cutler in 1931 [START_REF] Cutler | Transillumination Of The Breast[END_REF], was used to image breast lesions using near-infrared (NIR) light. Unfortunately, it was demonstrated a few years later that such technique was inefficient for breast tumour diagnosis [START_REF] Monsees | Light scan evaluation of nonpalpable breast lesions[END_REF]. The main issue with transillumination techniques such as diaphanography is that the resolution is limited to one fifth of the imaging depth because of scattering [START_REF] Moon | Resolution limits for imaging through turbid media with diffuse light[END_REF]. Consequently, if one wants to image a centimetre wide organ such as breast, it is impossible to resolve millimetre size objects such as tumours or blood vessels. [START_REF] Beard | Biomedical photoacoustic imaging[END_REF]). Hb0 2 and HHb denote respectively oxygenated and de-oxygenated blood. The window of low absorption between 700 and 1200 is commonly called the Optical Therapeutic Window

Another issue with optical imaging is light absorption by tissues. Most constituents of biological tissues absorb visible and non visible light. Water strongly absorbs non visible light and chromophores (coloured species) absorb in the visible domain as shown on figure 1.1 extracted from [START_REF] Beard | Biomedical photoacoustic imaging[END_REF]. Blood, for example, is strongly absorbing below 600 nm. Biological samples are thus strongly absorbing in general, except within a particular range of the spectrum, between 700 nm and 1.2 µm. This region is called the optical therapeutic window and it explains why light transmitted by tissues appears red. For imaging purposes, people usually work between 700 nm and 900 nm.

Light-matter interaction

Light can interact with matter through two main phenomena: absorption and scattering. Depending on their relative magnitude, they can strongly influence the behaviour of light propagation and have an important impact on our ability to obtain optical information from inside the medium.

Absorption

Absorption results from the interaction of light with the transition levels of atoms (electronic transition) or molecules (electronic, vibrational...) that compose the medium. If the energy of one photon corresponds to that of a transition, it can be absorbed. The absorbed energy is then dissipated through radiative processes (fluorescence, phosphorescence...) or non-radiative processes (heat). Absorption is characterised by the absorption coefficient µ a which quantifies the global amount of energy lost along with light propagation. This coefficient is defined for a non-scattering medium by Beer-Lambert law, Φ(x) = Φ 0 e -µax , (1.1) where x is the distance travelled by light, Φ 0 is the optical incident flux and Φ(x) is the flux after propagation. µ a is usually expressed in cm -1 or mm -1 . It is also possible to describe absorption by the means of the absorption mean free path, l a which represents the average distance travelled by a photon before being absorbed. This length is the inverse of the absorption coefficient

l a = 1/µ a . (1.2)
Absorption is strongly dependent on the wavelength of light. For biological tissues, the absorption coefficient is typically 0.5 cm -1 on average within the optical therapeutic window around 800 nm. [START_REF] Vo-Dinh | Biomedical Photonics Handbook[END_REF].

Scattering

In a homogeneous medium, light travels in a straight line and interacts with the medium through absorption only. However, most media are not homogeneous at a microscopic scale and have a complex micro-structure associated with refractive index changes. Consequently, light travelling through such medium will encounter scatterers which will change its direction through refraction. This phenomenon is called scattering and has been extensively studied and was the subject of several books by Ishimaru [START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF][START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF] or Bohren and Huffman [START_REF] Bohren | Absorption and scattering of light by small particles[END_REF]. This can be described by two theories, depending on the size of the scatterers.

Rayleigh Scattering

When interacting with a particle, an electromagnetic field induces charges oscillations and creates local oscillating dipoles which will then radiate in another direction. When the size of the particle is small compared to the wavelength, the field is seen as homogeneous by the particle and all the local dipoles will oscillate in phase. This means that the particle will behave as a unique dipole. Far away from the particle, the scattering is isotropic in first approximation, meaning that each directions is as likely as another. This regime is known as Rayleigh scattering.

When several scattering events occur, it is possible to characterize scattering by a length called the scattering mean free path, l s which represents the average distance between two scattering events. One can then define a macroscopic scattering coefficient, µ s by taking the inverse of the scattering mean free path µ s = 1/l s .

(1.3) This is the scattering equivalent of the absorption coefficient and is usually expressed in cm -1 or mm -1 . The value of µ s depends on the concentration of scatterers and their nature (shape, composition...).

Mie Scattering

When the scatterers size is large compared to the wavelength of the electromagnetic wave, the field cannot be considered homogeneous over the particle any more. In this case, scattering is no longer isotropic and mostly occurs in the forward direction. This is described by Mie theory [START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF] which is a general theory for scattering. The Rayleigh regime is just a particular case of Mie theory. The degree of anisotropy can be defined by taking the average value of the projection of the new direction on the incident one. This anisotropy factor, g, is given by g = cos(θ) , (1.4) where θ is the angle of scattering as defined on figure 1.2 and < . > denotes the statistical average. The anisotropy coefficient varies between 0 and 1. The case g = 0 corresponds to isotropic scattering as described by Rayleigh theory while in the case g = 1 all photons are scattered in the incident direction.

In the case of Mie scattering, a memory of the incident direction is preserved over several scattering events because of the anisotropy. It is then possible to define a characteristic distance over which light can propagate before losing the memory of its incident direction and polarization. This distance is called the transport mean free path, l * . Figure 1.3 shows a representation of this transport mean free path. two scattering events. The transport mean free path corresponds to the average of the blue arrows lengths. This represents the average distance after which the light has lost the memory of its initial direction and polarization.

We can once again define a characteristic coefficient, called the reduced scattering coefficient by taking the inverse of the transport mean free path µ ′ s = 1/l * . (1.5) This coefficient is also expressed in cm -1 or mm -1 .

µ s and l * quantify the ability of scattering to scramble the optical information. These coefficients are linked to the scattering mean free path and the scattering coefficient by

       µ ′ s = µ s (1 -g) l * = l s 1 -g . (1.6)

Orders of magnitude in biological tissues

Biological tissues are very inhomogeneous media, composed of a wide variety of structures ranging from a few nanometres for proteins to over ten microns for cells. Their optical properties thus vary by orders of magnitude depending on the tissues. Several reviews [START_REF] Cheong | A Review of the Optical Properties of Biological Tissues[END_REF][START_REF] Tuchin | Tissue Optics and Photonics: Light-Tissue Interaction[END_REF][START_REF] Jacques | Optical properties of biological tissues: a review[END_REF] and handbooks [START_REF] Vo-Dinh | Biomedical Photonics Handbook[END_REF] gather a lot of measurements of these properties. Table 1.1 gathers typical values found in [START_REF] Cheong | A Review of the Optical Properties of Biological Tissues[END_REF] for human tissues.

n Refractive index 1.35 -1.45 µ a Absorption coefficient 1 cm -1 -10 cm -1 at 630 nm µ s Scattering coefficient 100 cm-1 -500 cm -1 in the visible range g Anisotropy factor 0.9 -0.99 µ ′ s Reduced scattering coefficient 2 cm -1 -70 cm -1 in the visible range. ∼ 10 cm -1 at 800 nm Table 1.1 -Table of the typical optical properties of human tissues (taken from [START_REF] Cheong | A Review of the Optical Properties of Biological Tissues[END_REF])

These orders of magnitude concern a large variety of biological tissues such as blood, brain or bladder. The typical absorption coefficient value for tissues is 1 cm -1 around 800nm with the exception of epidermis (of the order of several tens of cm -1 ) and blood (a few cm -1 ) [START_REF] Vo-Dinh | Biomedical Photonics Handbook[END_REF]. Red blood cells are very anisotropic scatterers in the near infrared range with anisotropy factor around 0.99 but the rest of the tissues have an anisotropy factor of 0.9. With a typical scattering coefficient of 100 cm -1 , the usual value for the reduced scattering coefficient is 10 cm -1 .

Light propagation in scattering media

If absorption reduces the number of photons travelling through a biological sample, it is scattering that makes deep optical imaging challenging. Because of scattering, transmitted photons experience random trajectories that are impossible to reconstruct from only the output direction. Depending on the scattering strength, this can be an important problem for imaging. From ballistic photons which travel in a straight line to multiply scattered ones which have a completely random trajectories, several regimes can be considered.

Propagation regimes

When light travels through a scattering medium, it is possible to define two types of photons : the ballistic photons that did not experience any scattering event and the scattered photons that experienced one or several scattering events. Three scattering regimes can be defined depending on the thickness, L, of the medium. Since the scattering mean free path defines the average distance between two scattering events, it is useful to define these regimes. When L ≪ l s , most photons transmit through the sample before being scattered, this is called the ballistic regime. If L ∼ l s , it is likely that a few scattering events occur, it is called the single-scattering regime. When L ≫ l s , a great number of scattering events can happen before the photon exits the medium, this is called the multiple-scattering regime. (3)( 4) and ( 5) represents two particular cases of multiple scattering: (3) has a snake-like trajectory which appears in highly anisotropic scattering media and (4) and ( 5) are diffuse photons. (5) represents the particular case of back-scattered photons that are sent back in the direction of the source.

Ballistic light

Ballistic photons do not experience any scattering events. For these, it is possible to define a generalized Beer-Lambert law which takes into account scattering:

Φ(L) = Φ 0 exp[-(µ a + µ s )L), (1.7)
where Φ 0 is the incident flux and Φ(L) is the outgoing flux. The coefficient (µ a + µ s ) is called the extinction coefficient and is denoted µ ext . It is important to note that a photon that experiences a scattering event is no longer considered as a ballistic photon, which reduces the flux of ballistic photons. If L ≪ l s , most photons exit the medium before undergoing any scattering event and only ballistic photons are collected at the output. In this case, the medium is not turbid and the laws of geometrical optics can be applied. The only attenuation is due to absorption. The number of ballistic photons decreases as exp(-µ ext L), but in biological samples µ s ≫ µ a within the optical therapeutic window so µ ext ∼ µ s . This means that ballistic photons completely vanish after a several scattering mean free path i.e, a few hundred microns.

Single scattering

When the size of the medium increases, it becomes more likely for a photon to be scattered. An intermediate regime appears when L ∼ l s where most photons experience only one or two scattering events. In this single-scattering regime, ballistic and scattered photons coexist but very few photons experience multiple scattering. In this regime it is still possible to do imaging with ballistic photons.

Multiple scattering

When the scattering medium becomes thick compared to the scattering mean free path, L ≫ l s , a new regime, called multiple-scattering regime, appears where most photons experience multiple scattering events. In biological tissues, the scattering is generally not isotropic and the optical information is preserved over several scattering events. As explained before, this is no longer true after a few l * . For l s < L < l * , it is considered that simply scattered and multiply scattered photons coexist. In this case, the image quality for conventional optical imaging is driven by the ratio between ballistic or single-scattered photons and multiply scattered photons.

When L ≫ l * , the number of single scattered photons is negligible and the radiance is quasiisotropic. This regime is called the diffusive regime and the scattered light power flux obeys a diffusion equation. One can show that for a point source, the light fluence at a distance r larger than l * , Ψ s (r) is written [START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF]

Ψ s (r) ∝ Φ 0 exp[-3µ a (µ a + µ ′ s )r] r , (1.8)
where the proportionality is function of µ a and µ ′ s [START_REF] Van Staveren | Light scattering in Intralipid-10% in the wavelength range of 400-1100 nm[END_REF]. From this expression, we can define the effective extinction coefficient

µ eff = 3µ a (µ a + µ ′ s ) ∼ 3µ a µ ′ s
for strong scattering.

(1.9)

In strongly anisotropic medium, a particular case can be observed. Some photons, so-called snake-like or quasi ballistic photons, are multiply scattered but retain their original direction and propagate almost in a straight line as shown in figure 1.4. The power flux of these photons follows a generalised Beer-Lambert law

Φ qb (L) = Φ 0 exp[-(µ a + µ ′ s )].
(1.10)

Speckle Expression of the multiply scattered electric field

Optical medical imaging aims at visualizing entire organs over several centimetres. Given the reduced scattering mean free path of tissues (a few millimetres) it means that in most cases, light will be in the diffusive regime. Let us consider a static highly scattering medium illuminated by an electromagnetic plane wave propagating in the z-direction with a angular frequency ω 0 and a wave vector k 0 = 2π λ u z where u z is the unit vector of the z-axis. In the following, we will not consider the effect of scattering on the polarization of light and assume that it is random. In this scalar approximation the incident field is written .11) Given the ratio between absorption and scattering in biological sample, we will neglect absorp-

E 0 = E 0 exp[i(ω 0 t -k 0 • zu z )]. ( 1
E s (r ′ , t) = E 0 e iω 0 t p a p exp 2iπn 0 s p λ exp(ik 0 • (r ′ -r p )), (1.12)
where a p , s p and r p are respectively the amplitude, the length and the output position of path number p. The scattered field is here the coherent sum of all the different paths taken by photons inside the sample, see figure 1.5. If the coherence length of the optical source is high enough, the sum over random path creates a random interference pattern called speckle pattern.

Speckle

A speckle pattern appears when a highly scattering medium is illuminated with coherent light. This pattern is composed of bright spots -resulting from constructive interference -and dark spots -resulting from destructive interference -as shown in figure 1.6. The spatial distribution of amplitude and phase is random, consequently, the auto-correlation of a speckle pattern is non zero only on a characteristic area which size is equal to the size of a speckle grain.

A given speckle pattern corresponds to a given disordered medium. Consequently, it is extremely difficult to describe accurately. Speckle properties and statistics have been widely studied by Joseph W. Goodman [START_REF] Goodman | Some fundamental properties of speckle[END_REF][START_REF] Goodman | Speckle phenomena in optics: theory and applications[END_REF]. We will not extensively describe the physics of speckle here but we will point out two properties which are important for acousto-optic imaging.

The first one comes from the fact that the relative phase of speckle grains is random. Therefore, the integration of several speckle grains on a detector does not always increase the signal to noise ratio. This will be discussed in more details in chapter 2.

The second one comes from the fact that a given realisation of the medium disorder creates a given speckle pattern. If the medium is not static, the change in disorder will thus change the speckle pattern. Since this change is continuous, the pattern at time t 1 is correlated to the pattern at time t 2 as long as ∆t = t 2 -t 1 is short enough. This characteristic time is called the speckle decorrelation time of the medium and quantifies the time over which the speckle pattern changes.

Optical imaging in scattering media

Because of multiple scattering, optical imaging techniques can be divided in two categories. The first one uses ballistic light or single-scattered photons. This was previously referred to as conventional optical imaging. This regime can be interpreted in terms of geometrical optics but given the transport mean free path of biological tissues (l * ∼ 1mm), it does not allow for deep imaging. The second category relies on multiply scattered photons. Since the trajectories of such photons are random, the information they carry is more difficult to interpret. We will show here a brief overview of the different optical imaging techniques. More detailed information can be found in reviews like [START_REF] Ntziachristos | Going deeper than microscopy: the optical imaging frontier in biology[END_REF].

Imaging with ballistic light

Imaging with ballistic light is the most natural form of imaging because it is the way our eyes form images. The main advantage of such techniques is that they benefit from optical resolution. However, in biological tissues, their depth of imaging is seriously limited because of multiple scattering which will rapidly degrade the image quality.

Microscopy

The archetypal optical imaging technique is microscopy. Conventional transmission microscopes map the transmission of ballistic light through a sample and are thus sensitive to absorption while reflection gives a contrast of scattering by recording the single-scattered photons. Consequently, they are limited in terms of depth because multiply scattered light adds a lot of parasitic noise for depths higher than the scattering mean free path.

In order to increase the depth of imaging, one solution is to spatially filter out the multiply scattered photons. This is the case of confocal microscopy [START_REF] Inoue | Foundations of Confocal Scanned Imaging in Light Microscopy[END_REF] in which pinholes are used to limit the field of view to collect only the ballistic or single scattered photons. With such techniques, the imaging depth can be increased up to the transport mean free path, above which ballistic photons no longer exist. The main disadvantage of this technique is that it requires a mechanical scan of the focus to obtain images

Time Gating

Another approach to filter multiply scattered photons is to make use of the time of flight. Ballistic photons, single-scattered and snake-like photons travel in straight or almost straight lines whereas multiply scattered photons have complex trajectories. Consequently, it takes multiply scattered photons more time to travel through the medium. One can exploit this by sending a short pulse of light from a laser and by time gating the detection. Ballistic photons will usually go through the medium in a few hundreds of picoseconds while it takes nanoseconds for multiply scattered photons. The time gates thus need to be extremely small. This can be achieved with streak cameras, Time-Correlated single photon counting [START_REF] Desmond | Time-correlated single photon counting[END_REF] or by activating non linear effect such as Kerr effect [START_REF] Wang | Kerr-fourier imaging of hidden objects in thick turbid media[END_REF][START_REF] Martin | Picosecond laser stereometry light scattering measurements on biological material[END_REF], stimulated Raman scattering [START_REF] Duncan | Time-gated imaging through scattering media using stimulated raman amplification[END_REF] or Optical Parametric Amplification [START_REF] Gaëlle | Two-dimensional time-resolved direct imaging through thick biological tissues:?a new step toward noninvasive medical imaging[END_REF].

(a) (c) (d) (b)

Optical Coherence Tomography

Optical Coherence Tomography (OCT) is a technique that uses yet another method to filter the multiply scattered photons. Rather than using time gated detectors, OCT uses the coherence of the light source to localize information at a given depth. From the same observation we made before we can see that multiply scattered photons will have much longer illuminated optical paths inside the scattering sample than the ballistic ones. Consequently, if the sample is shone with a low coherence source, only the ballistic or single scattered photons will stay coherent with the source. This is the principle of OCT as proposed by Fujimoto's team in 1991 [START_REF] Huang | Optical Coherence Tomography[END_REF].

Two different OCT techniques exist. The first one, in-line OCT, images a longitudinal section of the sample and needs mechanical scanning of the optical beam to recover a full image. A typical in-line setup is shown on figure 1.7(a) along with a typical image 1.7(b). The second technique is called full-field OCT and performs en-face images of the sample. Figure 1.7(c) and (d) show a typical setup of full field OCT as well as a typical image. Since these techniques work with single-scattered photons, they are also limited to the transport mean free path.

Imaging with multiply scattered light

Imaging techniques working with ballistic or single-scattered photons have the great advantage of having optical resolution. Yet they are limited to a few transport mean free paths which corresponds to a few millimetres in tissues. If one wants to image deeper in biological sample, several techniques using multiply scattered light exist. However, because of the nature of multiply-scattered photons, extracting the information from the outgoing light is not as straightforward. The imaging techniques roughly fall into two categories: the techniques that use detailed model of the light propagation to reconstruct optical contrast inside a scattering medium (NIRS, DOT...) and methods that combine light with another type of wave to extract the information (Photoacoustics, Acousto-optics...).

Near Infra-Red Spectroscopy

Near Infra-Red Spectroscopy (NIRS) is a technique which can probe the optical absorption at a few centimetres depth inside multiply scattering medium. A sample is illuminated by a near infra-ref source and light is collected on the same side by a detector. The most likely trajectories form a "banana shape" between the source and the detector which probes the medium at a depth related to the distance between source and detector, as shown on 1.8(a). NIRS is not properly speaking an imaging technique since it sensitive to the average absorption on a wide area, but it is possible to probe at different depths by changing the source-detector distance. NIRS is also used to perform spectroscopy by changing the wavelength of the input source. Since it has been proposed in 1977 [START_REF] Ff Jobsis | Noninvasive, infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters[END_REF] is it commonly used to measure the oxygenation of blood, and a variety of modalities have been developed (e.g. NIRS-CW (continuous wave), NIRS-TD (time domain) or NIRS-FD (frequency domain)).

Diffuse Optical Tomography

The recent advances in modelling light propagation in highly scattering media lead to the development of Diffuse Optical Tomography (DOT) [START_REF] Boas | Imaging the body with diffuse optical tomography[END_REF]. The basic principle is to surround the medium with a high number of near infra-red sources and detectors, as shown on figure 1.8(b), extracted from [START_REF] Elizabeth | Optical brain imaging in vivo: techniques and applications from animal to man[END_REF] and to detect the multiply scattered photons on the boundaries of the medium. Then the propagation of light between sources and detectors is simulated using the Radiative Transfer Equation or one of its approximations and the collected data are used with a model based inversion to reconstruct the optical properties [35].

DOT has been implemented in several commercial devices and can be used for breast cancer imaging [START_REF] Pogue | Quantitative Hemoglobin Tomography with Diffuse Near-Infrared Spectroscopy: Pilot Results in the Breast[END_REF] or mapping of the brain haemodynamics for functional imaging [START_REF] Elizabeth | Optical brain imaging in vivo: techniques and applications from animal to man[END_REF]. DOT can give 3D images of the absorption or scattering coefficient with a resolution which is usually of the order of 5 to 10 mm which is sometimes too large to detect early stage tumours. Extracted from [START_REF] Elizabeth | Optical brain imaging in vivo: techniques and applications from animal to man[END_REF] Due to the physics of the propagation the inverse problem is ill-posed, which means that a large change in the parameters of interest tend to result in small changes in the measurements. Inverse solution will thus tend to amplify these small differences and prior knowledge or assumptions are needed to obtain a physical solution. Moreover, if one wish to reconstruct a 3D volume with a useful spatial resolution, the inverse problem will also generally be under-determined. The number of voxels in the image will be much greater than the number of measurement which is another source of non uniqueness in the solution. [START_REF] Boas | Imaging the body with diffuse optical tomography[END_REF].

Multiwave imaging

Imaging with multiply-scattered photons is challenging because each photon explores a wide area of the sample before being collected and it is impossible to know its trajectory. One solution to overcome this problem is to couple another type of wave to light [37]. In biological samples, ultrasonic waves with frequencies around ten MHz are barely scattered compared to light [START_REF] Sehgal | Scattering of Ultrasound by Tissues[END_REF] which makes ultrasound a prime candidate for multi-wave imaging. It is thus possible to use the non scattered ultrasound to localize the optical information at several centimetres depth [START_REF] Beard | Biomedical photoacoustic imaging[END_REF][START_REF] Laudereau | Multi-modal acousto-optic/ultrasound imaging of ex vivo liver tumors at 790 nm using a Sn 2 P 2 S 6 wavefront adaptive holographic setup[END_REF]. Two light-ultrasound coupling techniques exist: photoacoustic and acousto-optic imaging.

Photoacoustic imaging Photoacoustic imaging (PAI), also called opto-acoustic imaging is based on the photoacoustic effect discovered by Alexander Graham Bell at the end of the 19 th century. Bell noticed that absorbed light is able to produce an acoustic wave. Photoacoustics has been widely studied [START_REF] Beard | Biomedical photoacoustic imaging[END_REF] and the theory can be found in many works such as [START_REF] Akhmanov | Pulsed laser optoacoustics: Achievements and perspective[END_REF]. We will here give a brief description of the technique as a comparison to acousto-optic imaging.

In PAI, the sample is illuminated with nanosecond pulsed laser that is scattered inside. When light is absorbed by a chromophore inside the sample such as melanin or blood, the energy is transferred as heat and produces a local temperature rise. This increase in temperature is then dissipated through a thermal expansion which creates a broadband acoustic wave. The detection of the ultrasound wave on the boundary of the sample can be performed by an ultrasound transducer array [START_REF] Daoudi | Handheld probe integrating laser diode and ultrasound transducer array for ultrasound/photoacoustic dual modality imaging[END_REF] or by an optical mean such as a Fabry Perot cavity [START_REF] Jathoul | Deep in vivo photoacoustic imaging of mammalian tissues using a tyrosinase-based genetic reporter[END_REF] or a non contact holography setup [START_REF] Jayet | Acousto-optic and photoacoustic imaging of scattering media using wavefront adaptive holography techniques in NdYO 4[END_REF][START_REF] Jayet | Refractive index and gain grating in nd:yvo4: application to speckle vibrometry and photoacoustic detection[END_REF][START_REF] Hochreiner | Photoacoustic imaging using an adaptive interferometer with a photorefractive crystal[END_REF]. The collected data are then beamformed numerically the same way as for conventional ultrasound imaging, and a map of the absorption coefficient is produced with an ultrasonic resolution. One property of PAI is that the frequency of the ultrasound wave is inversely proportional to the size of the absorber. The resolution is limited thus by the bandwidth of the detector. Transducer arrays are easy to couple to commercial ultrasound scanners [START_REF] Kim | Deeply penetrating in vivo photoacoustic imaging using a clinical ultrasound array system[END_REF] and can image deep into tissues, but their limited bandwidth acts as a low pass filter. This prevents from imaging small absorbers and makes large absorbers appear through their contour [START_REF] Zhou | Handheld photoacoustic probe to detect both melanoma depth and volume at high speed in vivo[END_REF]. Optical detection schemes have the advantage of a much higher bandwidth and can potentially be non contact but they can be difficult to implement [START_REF] Jayet | Acousto-optic and photoacoustic imaging of scattering media using wavefront adaptive holography techniques in NdYO 4[END_REF][START_REF] Hochreiner | Photoacoustic imaging using an adaptive interferometer with a photorefractive crystal[END_REF].

Ultrasound transducer

Pulsed input light

Absorbing objects

The photoacoustic signal is classically detected on the same side the light input, which is an advantage for clinical imaging. Yet for this reason, it suffers from a problem which exists in conventional ultrasound. Elongated objects will mainly radiate acoustic waves orthogonally to their long axes as shown on figure 1.9 which makes them hard to detect if the long axis is orthogonal to the plane of detection. Several ideas have been proposed to overcome this difficulty such as ring-shaped [START_REF] Li | Real-time photoacoustic tomography of cortical hemodynamics in small animals[END_REF] or hemispheric [START_REF] Xiang | 4-D Photoacoustic Tomography[END_REF] transducer arrays coupled to a tomographic reconstruction. Another interesting modality is photoacoustic microscopy [START_REF] Yao | High-speed label-free functional photoacoustic microscopy of mouse brain in action[END_REF] which uses microscopes objectives to achieve optical resolution with PAI. This technique produces the same kind of images as conventional microscopy but with an absorption contrast.

Photoacoustic imaging has been a hot topic for the past 20 years and has shown its potential for clinical applications such as, for example, breast cancer detection [START_REF] Heijblom | Clinical Photoacoustic Breast Imaging: The Twente experience[END_REF][START_REF] Mallidi | Photoacoustic imaging in cancer detection, diagnosis, and treatment guidance[END_REF].

Acousto-optic imaging Acousto-Optic imaging (AOI), also called Ultrasound Modulated Optical Tomography (UOT) is another technique which combines light and ultrasound. In this case, the acoustic wave is generated by a transducer, the medium is illuminated by a laser and an optical signal is detected on a photo-detector. This technique will be extensively described in the next chapters so we will only give a brief overview here as a comparison to photoacoustic imaging. AOI relies on the acousto-optic effect. To image inside a scattering medium, the sample is illuminated by a continuous monochromatic light. At the same time, a focused ultrasound pulse is also sent. The ultrasound pulse modulates the refractive index of the medium as well as the scatterers positions. This results in a phase modulation of the light and the generation of side-bands shifted from the laser frequency by the frequency of the ultrasound. The number of modulated photons -so-called tagged photons -is proportional to the integrated optical fluence in the volume of the ultrasound pulse. This volume thus acts as a virtual source of tagged photons which will be collected on a photo detector. Despite their complex trajectories inside the sample, the position of the tagged photons' source is known, and they provide optical information about this region. The focus of the ultrasound can then be electronically or mechanically scanned to reconstruct a 2D or 3D map of the local fluence. Figure 1.10 shows a schematic of the AOI principle.

Similar to photoacoustics, the resolution of AOI is defined by the ultrasound. and corresponds here to the volume of the ultrasound focus. The main difference is that AOI is sensitive to both absorption and scattering where PAI is only sensitive to absorption.

The main difficulty of AOI is the detection of the tagged photons. It is important to discriminate the untagged photons in order to count the tagged photons, but the frequency shift between the two is equal to the ultrasound frequency, a few MHz. Compared to the 10 14 or 10 15 Hz optical frequency, this difference is too small to be resolved by any classical optical filter. This difficulty strongly delayed the development of acousto-optic imaging and explains why it is still at its early stages despite the fact that it was proposed at the same time as PAI by Leutz and Maret in the 90s [START_REF] Leutz | Ultrasonic modulation of multiply scattered light[END_REF].

Conclusion

Optical contrast can be very useful for medical imaging as it caused by both absorption and scattering and can thus give information about structure and metabolism. However, because of multiple scattering, the optical information is scrambled and after a few transport mean free paths, the ballistic photons vanish and only multiply scattered photons remain. This limits purely optical techniques to a few millimetres depth in biological tissues. In order to achieve deeper imaging, several techniques exist amongst which diffuse optical tomography, which relies on computationally expensive model based inversion or multi-wave imaging techniques which exploit the interaction between light and another type of wave. In biological samples, ultrasound propagation is well understood and controlled, and two main techniques use a light-ultrasound coupling to achieve centimetre depth imaging with an optical contrast and a millimetre resolution. Photoacoustic imaging uses light pulses to generate acoustic waves whereas acousto-optic imaging uses pulsed ultrasound to modulate light and extract optical information. The latter is the subject of this manuscript and will be explained in detail in the following chapters.

Chapter's Keypoints

• Light interacts with matter through two phenomena: absorption and scattering.

• Absorption is quantified by the absorption mean free path, l a , or equivalently by the absorption coefficient, µ a .

• Scattering is quantified by the scattering mean free path, l s , and the transport mean free path, l * , or equivalently by the scattering coefficient, µ s , and the reduced scattering coefficient, µ ′ s .

• Optical imaging is sensitive to absorption and scattering and thus to metabolism and structure.

• In scattering media, several propagation regimes exist: ballistic, single scattering and multiple scattering.

• When the propagation length, L, is smaller than l s , photons are ballistic.

• When L ∼ l s , photons experience a few scattering events, it is the single scattering regime

• When L ≫ l * , photons are multiply scattered and their trajectories are randomized.

• Multiple scattering with coherent light creates speckle patterns which are the result of random interferences between multiply scattered photons.

• The transport mean free path in biological tissues is typically a few millimetres.

• The imaging depth of optical imaging techniques using ballistic light is limited to a few millimetres.

• Optical imaging techniques using multiply scattered light require detailed models of the light propagation or combination to another type of wave.

• Photoacoustics and Acousto-optics are two imaging techniques that couple light and ultrasound to image at centimetre depths in scattering media.

Optical imaging in biological tissues is very challenging because of multiple scattering which scrambles the optical information. Imaging deeper than the transport mean free path -a few millimetres in biological samples -is almost impossible for classical optical methods. To perform imaging at centimetre depths, it is possible to use an intricate model of light propagation (DOT) but it comes with high computational times. Light can also be coupled to another type of wave which is not scattered in the medium of interest (photoacoustic or acousto-optic imaging).

In this chapter, we will focus on one of these multi-wave techniques, acousto-optic imaging (AOI), which uses emitted ultrasound to map the local fluence inside scattering media. The two crucial components of this method are the modulation of light by the ultrasound and detection of this modulated light.

We will start with a short description of the principle behind the technique and its potential applications, then we will theoretically describe the phenomena implicated in the modulation of scattered light. We will then discuss the image formation and the resolution of this technique. After detailing the difficulties of detecting this modulated light we will finish by giving an overview of the different possible detection schemes.

Principle of acousto-optic imaging

The acousto-optic effect

The acousto-optic effect was demonstrated in transparent crystals at the same time in 1932 by two independent teams, Debye and Sears [START_REF] Debye | On the Scattering of Light by Supersonic Waves[END_REF], and Lucas and Biquard [START_REF] Lucas | Propriétés optiques des milieux solides et liquides soumis aux vibrations élastiques ultra sonores[END_REF]. They showed that when an ultrasonic wave propagates inside a crystal, it modulates the density and thus the refractive index of the medium at the frequency of the ultrasound, f US . Since the wave is periodic, it will create a refractive index grating with a period equal to the wavelength of the ultrasound, Λ. Consequently, if the crystal is illuminated, light will be diffracted on this grating. If the grating can be considered as infinitely thin, the diffraction directions are given by:

sin θ p = p λ Λ , p = 0, ±1, ±2, ..., (2.1)
where λ is the wavelength of light, and p is an integer called the diffraction order. This regime is called the Raman-Nath diffraction regime and in addition to the light diffraction, the frequency of each diffracted beam is shifted by harmonics of the ultrasound frequency:

f p = f L + pf US , (2.2)
where f L is the frequency of the light before diffraction. Figure 2.1 summarizes both effects. In a scattering medium, this effect is more complex and we will detail it later, but the principle remains the same. If monochromatic light illuminates the medium, the ultrasound will generate side-bands shifted by the frequency of the ultrasound. In practice, the energy of the photons shifted by ±f US is much higher than the energy of higher harmonics, figure 2.2 shows a schematics Pressure wave of this. The photons that have been shifted in frequency are called tagged photons and the ones that did not interact with the ultrasound are called untagged photons. The most important properties of the tagged photons for acousto-optic imaging are the following:

f L + 2 f US f L -2 f US f L -f US f L + f US f L f L
• Tagged photons are generated by the modulation of light by the ultrasound. Consequently, they can be generated only in a region were both light and ultrasound are present.

• The number of tagged photons is proportional to the local light fluence in the region where they are created. Counting the number of tagged photons thus gives an information on the local light fluence at the focus of the ultrasound. If the acoustic wave is tightly focused, tagged photons will give information about this specific region of space, whatever the trajectories of the photons before or after tagging.

Energy density

f L f L -f US f L +f US f L

Energy density

If there is an absorbing object in the scattering medium, it will absorb light and reduce the local light fluence. Consequently, the number of tagged photons in this region will be smaller and so to will be the AO signal. One can also see that a scattering inclusion will have an influence on the AO signal because the inclusion will scatter the photons in different directions. In practice, in the diffusive regime, AOI gives information about the extinction coefficient µ eff = 3µ a (µ a + µ ′ s ) and it is difficult to separate the contributions of absorption and scattering.

From these observations we can see that if we are able to count the tagged photons, we can probe the local fluence in the volume where the ultrasound is focused and make an image of the fluence by scanning the focus of the ultrasound. However, to count the number of tagged photons, one needs to be able to filter out the untagged photons. As mentioned previously, this is a difficult task for several reasons:

• the frequency shift between tagged and untagged photons is extremely small. The tagged photons are shifted by the frequency of the ultrasound which is usually a few MHz (∼ 10 -6 nm), compared to the 10 14 or 10 15 MHz of the light (at 800 nm). No conventional optical filter has a finesse small enough to filter the untagged photons;

• the number of tagged photons is extremely small. Given the size of the ultrasound focus compared to the surrounding medium and the relatively low tagging efficiency very few photons are tagged by the ultrasound compared to the untagged -around 1 in 1000 photons;

• speckle grains have a random phase and amplitude distribution. This means that averaging over a high number of speckle grains gives the same signal-to-noise ratio (SNR) as looking at one grain;

• speckle decorrelates in a few milliseconds in dynamic biological samples. This means that to perform in vivo measurements, the detection system needs to be faster than 1 ms, which is not always possible given all the issues raised above.

All these difficulties limit the number of different techniques that can be used to detect the tagged photons but still, several approaches are possible and we will give an overview at the end of this chapter.

Applications of the technique

Acousto-optic imaging is sensitive to both absorption and scattering which makes it interesting to detect colour changes in biological tissues. Even though it is still at the proof-of-concept stage, several recent works show promising possibilities for AOI.

Since tumours are highly vascularized and contain important quantities of melanin, they strongly absorb light and appear dark. AOI was first thought as a tool for breast cancer diagnosis [START_REF] Wang | Continuous-wave ultrasonic modulation of scattered laser light to image objects in turbid media[END_REF] and a lot of work has been done towards this goal. Laudereau et al. used acousto-optic imaging to detect ex vivo liver tumours embedded in a scattering sample [START_REF] Laudereau | Multi-modal acousto-optic/ultrasound imaging of ex vivo liver tumors at 790 nm using a Sn 2 P 2 S 6 wavefront adaptive holographic setup[END_REF] as shown on figure 2.3 taken from the same publication.

Two different teams, at ESPCI Paris [START_REF] Benoit | Imagerie acousto-optique dans les milieux diffusants épais : de l'amélioration technique à l'application pré-clinique ex vivo[END_REF] and at Boston University [START_REF] Lai | Real-Time Monitoring of High-Intensity Focused Ultrasound Lesion Formation Using Acousto-Optic Sensing[END_REF] independently showed that AOI could be used for monitoring of High Intensity Focused Ultrasound (HIFU) treatments [START_REF] Ter | Ultrasound focal beam surgery[END_REF]. This process consists of sending a high intensity ultrasound beam to heat up and destroy tumorous cells at the focal point of the ultrasound. It was noted that the colour of burnt tissues changes, suggesting that optical imaging technique could be used to monitor the evolution of the lesion. Acousto-optic imaging can also be used as a spectroscopic technique to measure blood oxygenation deep inside tissues as shown by Honeysett et al. [START_REF] Honeysett | An algorithm for sensing venous oxygenation using ultrasound-modulated light enhanced by microbubbles[END_REF]. Using acoustic contrast agents, they produced quantitative measurement of blood oxygenation in deep vessels where NIR techniques fail to obtain quantitative measurements if superficial saturation varies. Kobayashi et al. also showed that acousto-optics can be used to locate fluorescent [START_REF] Kobayashi | Fluorescence tomography in turbid media based on acousto-optic modulation imaging[END_REF] or chemiluminescent [START_REF] Kobayashi | Optical tomography of fluorophores in dense scattering media based on ultrasound-enhanced chemiluminescence[END_REF] objets inside highly scattering media. In a paper by Huynh et al. , the authors even show that it is possible to perform acousto-optic detection of fluorescence without using an external light source [START_REF] Huynh | Ultrasound modulated imaging of luminescence generated within a scattering medium[END_REF].

AOI can also be used as a tool for photoactoustic imaging which is, as stated before, a more developed technique. In PAI, the intensity of the acoustic wave generated by absorbing objects depends on the light fluence inside the medium. To obtain quantitative images, the fluence has to be corrected. Daoudi et al. showed that AOI can be used for this purpose [START_REF] Daoudi | Correcting photoacoustic signals for fluence variations using acousto-optic modulation[END_REF].

Acousto-optics could also be used for non-biomedical applications such as non-contact security imaging to detect object hidden behind clothes, as suggested by the work of Lev and Sfez [START_REF] Lev | Acousto-optical detection of hidden objects via speckle based imaging[END_REF]. In applications where the ultrasound propagates in air, frequencies need to be significantly lower (50 kHz).

Even though it is not technically an imaging technique, it is worth mentioning the use of acoustooptics to refocus light inside a scattering medium. By time reversing the tagged photons -through phase conjugation [START_REF] Liu | Time-reversed ultrasonically encoded optical focusing into tissue-mimicking media with thickness up to 70 mean free paths[END_REF] or spatial light modulation (SLM) iterative optimization [START_REF] Ruan | Iterative Time-Reversed Ultrasonically Encoded Light Focusing in Backscattering Mode[END_REF] -it is possible to refocus light in the volume of the ultrasound focus. Though such techniques have many names depending on the refocusing means, the most emblematic is called Time-Reversed Ultrasonically Encoded focusing (TRUE focusing). These techniques could allow for detection of absorbers by Chapter 2. Acousto-optic imaging moving the ultrasound focus, and thus the refocussing point.

Light modulation in scattering media

In a scattering medium, the interaction between light and ultrasound is different from a case without scattering. Since the trajectories of photons are complex and random, it doesn't make sense to talk about beam deflection any more. However, phase modulation of light still occurs and is driven by two effects which both contribute to the modulation of the optical path: the modulation of the scatterers' position and the modulation of the refractive index of the medium.

The study and modelling of the acousto-optic modulation only started in the 1990s. At first only one of the phenomena was taken into account, usually the vibration of the scatterers. The first theory was published by Leutz and Maret in 1995 [START_REF] Leutz | Ultrasonic modulation of multiply scattered light[END_REF]. They modelled the light modulation by considering Brownian motion of the scatterers and a collective motion due to ultrasound. This model was completed by Kempe et al. in 1997 [68] by taking into account the inhomogeneity of the ultrasound field. In 2001, Wang analytically modelled the modulation of light by considering both the scatterers movements and the modulation of the refractive index [START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF]. For a given path p, the position of the j th scatterer is denoted r j and the length of the j th scattering event is l j . K US represents the wave vector of the acoustic wave.

In the following, we will consider a scattering medium as shown on figure 2.4 and we will consider that a photon travels in a straight line between two scatterers. This sample will be illuminated by an electromagnetic plane wave

E 0 = E 0 exp[i(ω 0 t -k 0 • zu z )], (2.3)
where E 0 represents the scalar amplitude, ω 0 is the light angular frequency and k 0 = k 0 u z = 2π λ u z is the light wave vector.

For simplification purposes, the scattering will be assumed isotropic and the scattering mean free path, l s , will be considered very long compared to the optical wavelength.

Modulation of the scatterers positions

The first phenomenon contributing to the modulation of light is the modulation of the scatterers' positions. When an acoustic wave propagates in a scattering medium, the scatterers will oscillate around their rest position, thus shifting the frequency of the light by Doppler effect. Let us consider an ultrasonic wave propagating in the sample. For simplification purposes, this wave is assumed to be a continuous plane pressure wave insonifying the whole medium.

P US (r, t) = A US sin(K US • r -ω US t), ( 2.4)
where A US represents the amplitude of the acoustic wave, K US is the acoustic wave vector, r is a vector representing the spatial coordinates, ω US is the angular frequency of the ultrasound, and t is the time variable.

The phase difference associated to the movement of the j th scatterer is then written [69]

Φ j s (t) = -n 0 k 0 (k j+1 -k j )A US sin(K US • r j -ω US t), (2.5)
where n 0 is the optical index of the medium, k j and k j+1 are respectively the unitary wave vectors of the light after scatterers j and j + 1, and r j is the position of the j th scatterer.

Modulation of the refractive index

The other phenomenon which impacts the tagging of light is the modulation of the refractive index. As a pressure wave propagates, the medium is compressed and rarefied depending on the location and time and this strain creates a refractive index modulation through photo-elasticity. The local changes of the refractive index, ∆n, creates a phase shift which can be written for the j th scattering event as

Φ j n (t) = l j 0 k 0 ∆n(r j , s j , K US • r j , t)ds j , (2.6)
where l j is the length of the j th scattering event and s j represents the coordinates along this path. The term ∆n(r j , s j , K US • r j , t) can be expressed as [START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF] ∆n(r j , s j ,

K US • r j , t) = n 0 ηK US A US sin K US • r j-1 -ω US t + K US • s j k j k 0 , ( 2.7)
where K US is the norm of K US , and η is a coefficient that depends on the adiabatic piezo-optical coefficient ∂n ∂P , the density of the scattering sample, ρ 0 , and the sound velocity in the sample, V US ,

η = ρ 0 V 2 US ∂n ∂P (2.8)

Coherent acousto-optic modulation

When the sample is illuminated by a coherent light, the light exiting the sample will form a speckle pattern because of random interferences between the different paths inside the sample. This speckle will be modulated by the ultrasound but because of the random walk inside the medium is it difficult to deduce any further information. To better understand the effect of Chapter 2. Acousto-optic imaging the acoustic modulation on light exiting the medium, several models have been developed and perfected over the years [START_REF] Leutz | Ultrasonic modulation of multiply scattered light[END_REF][START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF][START_REF] Sakadžić | Modulation of multiply scattered coherent light by ultrasonic pulses: An analytical model[END_REF].

Using the auto-correlation function of the electric field and the Wiener-Khinchin theorem, the amplitude of the light modulation can be retrieved. This theorem states that the power spectral density is equal to the Fourier transform of the auto-correlation function of the electric field, G(τ ). The power of the component modulated at mω US can thus be written:

Ψ p = ω US 2 2 ω US 0 G(τ ) cos(mω US τ )dτ m ∈ Z (2.9)
The autocorrelation function of the electric field is the summation over every possible path of the autocorrelation of the field E p (t), weighted by the probability of each given path p, Ω(p),

G(τ ) = ∞ 0 Ω(p) E p (t)E * p (t + τ ) ds p , ( 2.10)
where . * denotes the complex conjugate, . stands for the average value and s p represents the length of path p.

Taking into account both modulation of the refractive index and scatterers' displacement, and assuming their effects on the phase are independent, the autocorrelation of the field can be written [START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF]

E p (t)E * p (t + τ ) = exp -i   j ∆Φ j s (t, τ ) + j ∆Φ j n (t, τ )   (2.11)
where ∆Φ s (t, τ ) = Φ s (t + τ ) -Φ s (t) is the phase variation induced by the movement of the scatterers and ∆Φ n (t, τ ) = Φ n (t + τ ) -Φ n (t) is the phase variation induced by the modulation of the refractive index. Detailed calculation of G(τ ) can be found in the previously cited article [START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF]. Under the assumptions that the phase variations are much less than unity, the diffusion limit is verified (s p /l s ≪ 1) and the phase variation due to two phenomenons are independent, it can be written:

G(τ ) = ∞ 0 Ω(p) exp - 2s p l s (n 0 k 0 A US ) 2 (δ s + δ n )(1 -cos(ω US τ )) ds p , ( 2.12)
where δ s and δ n are two scalars that respectively represent the contribution of scatterers displacement and refractive index modulation. Their exact values can be found in [START_REF] Lihong | Mechanisms of Ultrasonic Modulation of Multiply Scattered Coherent Light: An Analytic Model[END_REF]. Using this formula it is possible to assess the relative contribution of both phenomena. In his article, Wang shows that both effects have a similar contribution while K US l s < 0.56 but that above this value the modulation of the refractive index rapidly becomes predominant. Indeed, when l s increases the distance between two scattering events becomes higher so the accumulated phase due to the refractive index modulation increases while the number of scattering events remains the same.

By developing the exponential in 2.12 under a small modulation depth hypothesis, the power spectral density of scattered light can be derived. This quantity contains a static component -which corresponds to untagged photons -and oscillating parts at multiples of the ultrasound frequency -which corresponds to tagged photons. By keeping only the first frequencies, one obtains:

        
Untagged photons power:

Ψ ut = Ψ 0 1 -(n 0 k 0 A US ) 2 (δ s + δ n ) 2 l s ∞ 0 Ω(p)s p ds p ;
Tagged photons power:

Ψ t = Ψ 0 (n 0 k 0 A US ) 2 (δ s + δ n ) 2 l s ∞ 0 Ω(p)s p ds p ;
(2.13)

where Ψ 0 represents the light power without ultrasound. The same calculation holds for ultrasonic pulses and was derived by Kempe et al. [START_REF] Kempe | Acousto-optic tomography with multiply scattered light[END_REF] and Sakadžić and Wang [START_REF] Sakadžić | Modulation of multiply scattered coherent light by ultrasonic pulses: An analytical model[END_REF]. The tagged photon power at a position r ′ outside the sample can then be expressed as:

Ψ t (r ′ ) = βΨ 0 Ω(r, r ′ )I(r)P US (r) 2 dr, (2.14)
where Ω(r, r ′ ) represents the probability for a photon tagged in r to reach the position r ′ outside the sample, I(r) represents the light fluence at position r, P US (r) represents the acoustic pressure at position r and β is a proportionality factor. This expression shows several interesting properties of the tagged photons. First, the tagged photons power is proportional to the square of the ultrasound pressure which will have an important impact on the acousto-optic signal. Second, the number of tagged photons depends on the light fluence integrated over the whole acoustic volume as qualitatively mentioned before, and finally it shows that the amount of signal depends on the probability of a photon to reach the detector. This means that the spatial shape of the acousto-optic signal will depend on the detector position. This will be of particular importance in chapter 4.

Image formation and resolution

We have seen that the number of tagged photons is proportional to the local light fluence in the volume of the ultrasound. The shape of this volume will thus define the resolution of the technique. Several approaches have been considered to construct AO images by sending different acoustic sequences. Each comes with a different image formation, resolution and signal to noise ratio. We will here give a brief overview of the different techniques.

Single-element transducer and continuous wave

The first proof-of-concept were obtained with single-element focused transducers and continuous wave ultrasound [START_REF] Wang | Continuous-wave ultrasonic modulation of scattered laser light to image objects in turbid media[END_REF][START_REF] Yao | Theoretical and experimental studies of ultrasound-modulated optical tomography in biological tissue[END_REF]. In this configuration, photons are tagged in the whole acoustic column along the direction of ultrasound propagation, meaning that the resolution in this direction is equal to the entire medium depth. The resolution in the transverse directions correspond approximately to the size of the focus (a few acoustic wavelengths). To construct an image, the transducer needs to be mechanically translated to probe different regions of the medium. A solution to obtain a longitudinal resolution is to rotate the transducer around the sample to form a tomographic image [START_REF] Li | Ultrasound-modulated optical computed tomography of biological tissues[END_REF] Pulsed ultrasound The most straightforward way to recover the longitudinal resolution is to use pulsed ultrasound instead of continuous wave. By sending a few-cycles-long pulse, the tagging volume is greatly reduced. The lateral resolution is still given by diffraction (a few acoustic wavelengths) but the longitudinal resolution is now the length of a pulse. Given the difference between the speed of light and the speed of sound, it is possible with a fast detector to follow optically the propagation of the pulse. To construct a line of the image, the number of tagged photons is recorded while the pulse travels through the medium as shown on figure 2.5. The pulse is then translated to form an new line and an image is constructed by stacking the lines together. With such an approach, the number of tagged photons is significantly reduced because of the smaller volume occupied by the ultrasound. This thus requires more sensitive detection and a trade-off has to be found between resolution and amount of signal. With an acoustic frequency of a few MHz, typical longitudinal resolutions are of the order of a few millimetres.

Pulsed ultrasound was used in combination with photorefractive-based detection -which will be the subject of the next chapter -and broadband single detectors to acquire time-resolved AO signals over the pulse propagation [START_REF] Sui | Imaging in diffuse media with pulsedultrasound-modulated light and the photorefractive effect[END_REF].

Acousto-Optical Coherent Tomography

To increase the number of tagged photons while keeping the longitudinal resolution other approaches were developed by encoding information in the ultrasound sequence. In 2009, Lesaffre et al. developed a technique called Acousto-Optical Coherent Tomography (AOCT) which was the subject of several publications [START_REF] Lesaffre | Acousto-optical coherence tomography using random phase jumps on ultrasound and light[END_REF][START_REF] Lesaffre | Theoretical study of acoustooptical coherence tomography using random phase jumps on ultrasound and light[END_REF][START_REF] Lesaffre | Experimental study of z resolution in acousto-optical coherence tomography using random phase jumps on ultrasound and light[END_REF][START_REF] Lesaffre | Imagerie acousto-optique de milieux diffusants épais par détection photoréfractive[END_REF]. This technique uses long ultrasound pulses and long exposure detectors to accumulate a high number of tagged photons. In order to obtain a longitudinal resolution, it relies on an interferometric approach and encoding of the ultrasound. The long pulses are modulated by a pattern of short random phase jumps and a reference beam, modulated with the same pattern, is used to recover the resolution. Only the photons which have exactly the same phase as the reference beam can interfere constructively while the other interfere randomly. The temporal coherence volume of the pulse thus defines an area of interest in which the phase state of the tagged photon always matches the one of the reference. The position of this area can be chosen by tuning the delay between the ultrasound and the reference beam modulation. The size of this volume and thus the resolution is linked to the characteristic time between two phase jumpsi.e. the temporal autocorrelation of the modulation pattern. This techniques allows for a high number of measured tagged photons while keeping a spatial resolution of the order of a few millimetres.

Fourier-based imaging

Another approach developed more recently by Barjean et al. [78,[START_REF] Barjean | Theoretical study of Fourier-transform acoustooptic imaging[END_REF] uses long ultrasound pulses and a Fourier-based modulation to obtain depth-resolved acoustooptic images. This technique, called Fourier-Transform Acousto-Optic Imaging (FTAOI) consists of sending an amplitude modulated (AM) ultrasound beam and measuring the tagged photons signal through demodulation at the ultrasound frequency. It can be shown [START_REF] Barjean | Theoretical study of Fourier-transform acoustooptic imaging[END_REF] that this signal contains information about the Fourier coefficient of the acoustic line at the AM frequency. This information is recovered by demodulating the tagged photons signal a second time at the AM frequency. After recording the Fourier components corresponding to a set of AM frequencies, an AO line can be reconstructed by taking the inverse Fourier transform. In this case, the resolution is linked to the highest AM frequency that can be measured.

Transducers arrays

To form an acousto-optic image requires the probing of a 2D plane of the medium. Doing so with a single-element transducers necessitates a mechanical translation of the transducer. This process can be very slow and often lacks accuracy. The use of commercial ultrasound scanners in acousto-optic experiments allows for a more robust scanning of the sample by using a transducer array and an electronic steering of the acoustic beam. All the elements of the transducer are simultaneously controlled electronically and delays can be applied between them to focus an ultrasonic wave anywhere in the 2D plane below the transducer. Similarly as with single-element transducers, an image is produced by stacking acousto-optic lines, but this time there is no movement of the transducer, only changes in the electronic sequences. This process will be explained in more details in chapter 3.

Plane wave acousto-optic imaging

The use of ultrasound scanners for AOI led to the development of a new technique at Institut Langevin by J.B Laudereau in 2016 called Ultrafast AOI with plane waves [START_REF] Laudereau | Ultrafast acousto-optic imaging with ultrasonic plane waves[END_REF][START_REF] Laudereau | Acousto-optic imaging : challenges of in vivo imaging[END_REF]. This approach, inspired by ultrafast ultrasound imaging [START_REF] Tanter | Ultrafast imaging in biomedical ultrasound[END_REF], consists of sending plane waves ultrasound instead of a focused beam. The pulses are still a few-cycles-long but the tagging area is now a whole line of the medium. The longitudinal resolution will still correspond to the length of a pulse but if only one plane wave is sent, the transverse resolution will be the whole width of the scattering medium. In order to obtain a decent lateral resolution, a tomographic approach is used and a series of plane waves are sent each with a different angle of propagation as seen on figure 2.6. This is possible thanks to controllable delays between the elements of the transducer array which allow for a good control of the ultrasound wave-front shape. A tomographic reconstruction is then used to form the image. This technique has the main advantage that it does not require any scanning of the beam and is thus considerably faster. This method also benefits from a higher number of photons since the ultrasound power sent in the medium is significantly higher. On the other hand, the lateral resolution is limited by the number of angles that can be sent. A perfect tomographic reconstruction would require to send a set of plane waves with angles ranging from 0 to 360 • . In practice, this is not possible because of the shape of the probe and a restricted number of angles can be sent which creates a distortion of the image in the lateral direction and affects the resolution.

Detection of tagged photons

We have seen in the previous sections that the ultrasound modulation of light carries information about the local light fluence in the ultrasound volume and that this volume can be manipulated to probe the scattering medium and obtain a map of the light fluence. Since AOI works with multiply scattered photons, it is not limited by the transport mean free path and images can be obtained through several centimetres of biological tissues. Yet to access this information, it is necessary to quantify the amount of modulated light. To do so, the tagged photons need to be separated from the untagged photons. This constitutes the main challenge for acousto-optic imaging for the reasons already mentioned in section 2.1. Not only is the number of tagged photons very small compared to the total number of photons but the spectral distance between modulated and unmodulated light is also extremely small (10 -5 nm). This prevents one from using classical optical filters and a less straightforward detection scheme has to be found. When using a coherent detection method, another problem arises. The random phase and amplitude of speckle grains means that averaging the signal over several speckle grain does not increase the SNR. Finally, in dynamic media, the speckle can have a temporal decorrelation because of Brownian motion and intrinsic movement of biological tissues.

Over the years, several approaches have been proposed to overcome these difficulties and produce acousto-optic images. We will here give a brief overview of these techniques which fall into two categories: incoherent techniques which aim at spectrally filtering the tagged photons and coherent techniques which aim at measuring the modulation of light induced by the ultrasound.

Tagged or untagged photons

Before describing the different techniques that can be used for the detection of the acousto-optic signal, a remark has to be made on tagged and untagged photons. In the previous sections we assumed, for the sake of clarity, that the information about the light fluence was only carried by the modulated light and that the unmodulated light needed to be filtered out. However, light modulation can be a transient phenomenon -when using pulsed ultrasound for exampleand it can add another layer of difficulty to a problem that is already complex. Fortunately, conservation of energy means that if a tagged photon is created, an untagged photon disappears. This means that the quantity of unmodulated light is closely related to the number of tagged photons. Consequently, it is also possible to obtain information about the local light fluence by measuring the component at the frequency f L and filtering out the components at f L ± mf US . In the following, the type of detected photons will not always be specified but this does not impact the reasoning.

Incoherent methods

Fabry-Perot interferometer

Since AOI aims at imaging modulation of light created by the vibration of the ultrasound, it is very similar to laser vibrometry which measures the vibration of a surface on which a laser is reflected. Consequently, it is not surprising that a lot of detection schemes are adapted from the field of laser vibrometry. Such methods include the use of the Fabry-Perot interferometer which was used in non destructive control to detect waves propagating on the surface of materials [START_REF] Monchalin | Optical detection of ultrasound at a distance using a confocal Fabry-Perot interferometer[END_REF]. A Fabry-Perot interferometer acts as an interferometric filter with a central frequency tunable by changing the distance between the two partially transmitting mirrors. The filtering power of such a cavity is determined by the finesse, a coefficient equal to the ratio between the free spectral range -inversly proportional to the cavity length -and a coefficient which depends on the reflectivity and the surface state of the mirrors. For acousto-optic applications, the transmission of the interferometer is tuned on the tagged photons (f L ± mf US ) so that only the modulated light is transmitted while the rest is reflected. A simple single-element photodetector can then be used to record the power flux of the tagged photons. This was one of the first detection techniques suggested for acousto-optics [START_REF] Leutz | Ultrasonic modulation of multiply scattered light[END_REF] but the optical etendue was too small. More recently, Sakadžić et al. [START_REF] Sakadžić | High-resolution ultrasound-modulated optical tomography in biological tissues[END_REF] as well as Rousseau et al. [START_REF] Rousseau | Ultrasound-modulated optical imaging using a high-power pulsed laser and a double-pass confocal Fabry-Perot interferometer[END_REF] used confocal Fabry-Perot interferometers which have a better etendue and a filtering resolution that can be increased by imposing a double-pass of light in the interferometer.

Spectral holeburning

Because of their limited etendue, Fabry-Perot cavities are inconvenient for speckle filtering. Recently, Li et al. suggested the use of spectral holeburning to filter modulated light [START_REF] Li | Pulsed ultrasoundmodulated optical tomography using spectral-hole burning as a narrowband spectral filter[END_REF]. Spectral holeburning is a technique which consists in burning a narrow transparency window in the frequency spectrum of an absorbing material. This method concerns material which have an inhomogeneous broadening of their absorption line. Rare-earth doped inorganic crystals are a medium of choice for this application [START_REF] Roger M Macfarlane | Direct process thermal line broadening in tm : {YAG}[END_REF]. Because of inhomogeneities, each different site of the crystal absorbs at different frequencies which correspond to transitions between atomic levels. Spectral holeburning consists of exciting one of these transitions by illuminating the medium with a very monochromatic light at a frequency ω P . This will selectively excite the atoms with a transition energy corresponding to ω P , meaning that they will not be able to absorb any more photons. This creates a gap of in the absorption spectrum of the crystal with a width equal to the homogenous linewidth ∆Γ H , as represented on figure 2.7(a). If the frequency of the pump is tuned on the frequency of the tagged photons, only these will be transmitted by the crystal while the unmodulated light will be absorbed. In normal conditions, the homogeneous linewidth is very large (some GHz) and can be broader than the inhomogeneous bandwidth ∆Γ I which corresponds to a broad and not selective filter. Yet it is possible to reduce the homogeneous linewidth to some kHz by cooling the crystal down to low temperatures (≤ 4 K). The depth of the hole is proportional to the energy of the pump beam, up to a saturation limit. To obtain an efficient filter it thus requires an energetic pump beam at the same frequency as the tagged photons. Consequently, the pump needs to be switched off in order to detect the tagged photons. The lifetime of the hole is thus of prime importance. Since the excited atoms will spontaneously relax to lower energy states through several processes, it is important to select a material in which the characteristic relaxation time is long enough to perform detection of the tagged photons. Usually, the selected crystal possesses intermediate states with long relaxation times as shown on figure 2.7 (b). This technique has a high optical etendue, only limited by the size of the cryostat window, and it has been shown that it can image through several centimetres of tissues by Xu et al. [START_REF] Xu | Spectral hole burning for ultrasound-modulated optical tomography of thick tissue[END_REF], who increased the efficiency of the spectral holeburning by double-passing of the pump in the crystal.

The main limitation of this technique, in addition to the constraint of maintaining the crystal at a cryogenic temperature, is the linewidth and the stability of the laser. Recently, our team worked in collaboration with J.L. Le Gouët and A. Chauvet from Aimé Cotton laboratory to increase the filtering efficiency of spectral holeburning by putting a T m 3+ : YAG crystal in a magnetic field. The magnetic field will split some energy states by Zeeman effect and increase the lifetime of the transitions states thus increasing the time persistence of the spectral hole. This has been reported in J.B. Laudereau PhD thesis [START_REF] Laudereau | Acousto-optic imaging : challenges of in vivo imaging[END_REF]

Coherent methods

Coherent detection methods use another approach to detect the acousto-optic signal. Rather than filtering the tagged photons, they aim at detecting the modulation of light either through synchronous detection or excitation or through interferometry. Because they are sensitive to the coherence of light, these techniques suffer from the random properties of the speckle pattern.

Fast mono-detector

The first and most straightforward approach at detecting the modulated light was proposed in the mid 1990s and consists of using a fast single element photodetector associated with an electronic filtering system to extract the tagged photon component [START_REF] Wang | Continuous-wave ultrasonic modulation of scattered laser light to image objects in turbid media[END_REF][START_REF] Kempe | Acousto-optic tomography with multiply scattered light[END_REF][START_REF] Wang | Ultrasound-modulated optical tomography of absorbing objects buried in dense tissue-simulating turbid media[END_REF]. The imaging speed of this technique is considerably smaller than the decorrelation time of the speckle which allows for in vivo imaging. However, the random phase of the speckle is a major problem for this approach which integrates the signal over a large single detector. Since both modulated and unmodulated light arrive on the detector a trade-off has to be found between the amount of light collected and the amplitude of modulation. Wang et al. have shown that by detecting N g speckle grains on a photo-multiplier, the degree of modulation is divided by N g . Since the amplitude of noise is also reduced by a factor N g , they show that the signal-to-noise ratio for N g speckle grains is identical to the one obtained by detecting a single speckle grain [START_REF] Wang | Ultrasound-modulated optical tomography of absorbing objects buried in dense tissue-simulating turbid media[END_REF].

Parallel speckle processing

In order to overcome the difficulty of the random phase of speckle grains the first techniques that were used were camera based techniques which sampled speckle patterns to process each grain individually at the same time. The first use of this approach to detect the acousto-optic signal was reported by Lévêque et al. [START_REF] Lévêque | Ultrasonic tagging of photon paths in scattering media:?parallel speckle modulation processing[END_REF]. In this paper the authors suggest the use of parallel speckle detection with a 10 5 pixels CCD camera in order to process a great number of photons through a lock-in approach. They show that adding the signals detected in parallel on N g speckle grains increases the SNR by N g . Because of the low acquisition speed of the camera compared to the period of the ultrasound, it is no longer possible to directly record the modulation of the light. Lévêque et al. propose to use synchronous excitation rather than synchronous detection. Instead of using a local oscillator at the frequency f L ± f US on the detector, a modulation is applied to the laser at the frequency f US . This creates a low frequency beating that can be detected on the camera. Four images are then recorded with a different phase difference between the ultrasound excitation and the modulation of light in order to access amplitude and phase of the acousto-optic signal on each pixel through a linear combination of these images. This was later used in other papers such as [START_REF] Yao | Theoretical and experimental studies of ultrasound-modulated optical tomography in biological tissue[END_REF]. This procedure requires the recording of four images which increases the acquisition length and can deteriorate the signal if it takes longer than the speckle decorrelation time.

With the recent progresses in CMOS technologies came the development of smart-sensors arrays with a lock-in detection associated with each pixel allowing for real-time parallel demodulation of all speckle grains [START_REF] Barjean | Fourier transform acousto-optic imaging with a custom-designed CMOS smart-pixels array[END_REF][START_REF] Morgan | CMOS detection methods in ultrasound modulated optical tomography[END_REF]. The ultrasound modulation is extracted from the interference between scattered light and a reference tuned on the frequency of the tagged photons. The beating component is then demodulated independently for each pixel and it only requires a single image compared to the four-phases technique previously mentioned. The main inconvenient of this technique was the small number of pixels (64 × 64 and 24 × 24 for [START_REF] Morgan | CMOS detection methods in ultrasound modulated optical tomography[END_REF] and [START_REF] Barjean | Fourier transform acousto-optic imaging with a custom-designed CMOS smart-pixels array[END_REF] respectively) until recent advances were made. In 2016, Liu et al. demonstrated the use of a 300 × 300 pixels camera for acousto-optic measurement in decorrelating media [START_REF] Liu | Lock-in camera based heterodyne holography for ultrasound-modulated optical tomography inside dynamic scattering media[END_REF].

Speckle contrast

Another method for single shot imaging with cameras consists in an parallel analysis of the statistical properties of the speckle grains and its contrast in particular. This approach makes use of an inconvenience mentioned above: if the acquisition time of the camera is greater than the acoustic modulation period, the speckle recorded by the camera will be blurred and contrast will be degraded. The speckle contrast is defined as the ratio between the standard deviation of a speckle grain intensity and the standard deviation of the mean intensity for the whole image [START_REF] Goodman | Some fundamental properties of speckle[END_REF]. By measuring the amount of blurring it is possible to recover the modulation depth as shown in [START_REF] Zemp | Stochastic explanation of speckle contrast detection in ultrasound-modulated optical tomography[END_REF]. The first proof of concept was demonstrated by Li et al. [START_REF] Li | Ultrasound-modulated optical tomography of biological tissue by use of contrast of laser speckles[END_REF] in chicken breast.

The main problem of this technique is that the acquisition time of the camera needs to be smaller than the decorrelation time of the speckle in order to be sure that the blurring comes from the modulation and not from movements in the scattering medium. Recently, Resink et al. showed that such a detection scheme could be used with nanosecond lasers pulses to solve this speckle decorrelation issue or to couple acousto-optic and photoacoustic imaging [START_REF] Resink | Towards acousto-optic tissue imaging with nanosecond laser pulses[END_REF][START_REF] Resink | Tandem-pulsed acousto-optics: an analytical framework of modulated high-contrast speckle patterns[END_REF][START_REF] Resink | Tandem pulsed acousto-optics: obtaining the tagged light fraction from modulated non-ideal speckle patterns[END_REF].

Digital holography

Digital holography is a technique which consists of directly recording the electro-magnetic field corresponding to the tagged photons on a camera through the use of a reference beam that does not go through the scattering medium. This method, suggested for acousto-optic imaging by Gross et al. in 2003 [99], is similar to lock-in detection in the sense that it uses a local oscillator on the detector. However, rather than performing a temporal demodulation of the collected signal, digital holography relies on a the spatial demodulation of the interference pattern between the scattered light and the reference to recover the tagged photons component. This is generally used in an off-axis and heterodyne configuration [START_REF] Gross | Digital holography with ultimate sensitivity[END_REF][START_REF] Benoit A La Guillaume | Acousto-optical coherence tomography with a digital holographic detection scheme[END_REF] meaning that the reference is tilted and slightly frequency shifted by a few Hz in order to extract the modulated light in the plane of the camera. A simple 2D Fourier transform then gives the amplitude of the tagged photons in the output plane of the scattering sample, provided that the scattering sample is far enough from the camera. It has been shown that this technique is shot noise limited [START_REF] Verpillat | Digital Holography at Shot Noise Level[END_REF]

Wave-front adaptive holography

Fast single-element photo-detectors are interesting for acousto-optic imaging because of their high bandwidth that makes the time-resolved detection of the ultrasound propagation possible. However, as mentioned before, because of the random phase of the speckle grains, the signal-to-noise ratio is not increased by integrating the signal over several speckle grains. One solution to increase the SNR is to use holographic techniques that adapt a reference wave-front to the speckle wave-front in order to add them coherently on a photo-detector such as two-wave mixing (TWM).

TWM consist in recording the interference pattern between a reference beam and scattered light in a holographic medium and reading this hologram with the same reference beam. A great diversity of holographic media exist but in all cases, one of the properties of the medium is affected by light. Consequently, the interference pattern creates a grating inside the medium on which the reference can be diffracted. Since the grating was created by the reference itself and the signal beam, the diffracted reference is proportional to the scattered electric field and propagates in the same direction as the transmitted signal beam as shown in figure 2.8. The wave-front of the reference is "adapted" to the wave-front of the signal and both beam can interfere constructively on a large detector. The recording of the hologram is not an instantaneous process but a dynamic one which usually refreshes over a characteristic time which depends on the holographic medium. This means that if the interference pattern changes, the hologram will adapt over a certain time, hence the name "self-adaptive" holography that often refers to these techniques. In the case of acousto-optic detection, the reference can be tuned to the frequency of the tagged photons so that its wave-front is adapted to the wave-front of the modulated light only. This process will be explained in more detail in the next chapter.

Photorefractive crystals

The first materials that were studied for TWM were photorefractive crystals (PRC) which were investigated for imaging vibrations on rough surfaces [START_REF] Delaye | High-amplitude vibrations detection on rough surfaces using a photorefractive velocimeter[END_REF][START_REF] Campagne | Compact and fast response ultrasonic detection device based on two-wave mixing in a gallium arsenide photorefractive crystal[END_REF]. These crystals were independently suggested for AO detection by Murray et al. and Ramaz et al. , in 2004 [105,106]. TWM in photorefractive crystals will be further studied in the next chapter, we will only briefly describe its principle. Photorefractive crystals are both photoconductive and electrooptic materials which means that an inhomogeneous illumination will create an inhomogeneous charge distribution which will in turn create an inhomogeneous spatial repartition of the refractive index. PRCs can be characterised by two quantities, the photorefractive gain γ P R and the response time

τ P R .
The first one represents the amplification of the signal beam and is defined by exp(γ P R L) which is equal to the ratio between the signal beam power with and without reference beam. The second characteristic represents the time it takes to record an hologram in the crystal. The SNR depends, amongst other factors, on the number of speckle grains that can be collected on the crystal surface. With a 1 cm 2 crystal, it is possible to process up to 10 8 speckle grains, about two orders of magnitude more than cameras, which usually have up to 10 6 pixels.

Photorefractive crystal are commonly used for AOI because of their relatively high gains and their low response times. The most used crystals are Bi 12 SiO 20 (BSO) at 532 nm, Sn 2 P 2 S 6 (SPS) at 780 nm and GaAs at 1024 nm which have response times of the order of 5 to 10 ms for SPS and 0.3 to 100 ms for BSO and GaAs.

Photorefractive polymers and liquid crystals Recently, Suzuki et al. suggested the use of photorefractive polymers [START_REF] Suzuki | High-sensitivity ultrasound-modulated optical tomography with a photorefractive polymer[END_REF] which possess an optical etendue that can be up to one order of magnitude higher than for PRC. However, such materials are impractical for imaging dynamic scattering medium because of their high response time that is of the order of several seconds.

Another investigated solution is the use of liquid crystals in Liquid Crystals Light Valves (LCLV) which relies on their high birefringence. A nematic liquid crystal is combined to a photoconductive crystal and stimulated by electrodes. When light interacts with the photoconductive crystal, it creates a local electric field which will partially reorientate the liquid crystal and create a photorefractive effect through birefringence [START_REF] Bortolozzo | Pattern formation and direct measurement of the spatial resolution in a photorefractive liquid crystal light valve[END_REF]. The main advantage is that the photoconductive and the electrooptic effects are uncoupled so both can be optimized independently. This technique was suggested by Bortolozzo et al. to detect small displacement of a vibrating surface [START_REF] Bortolozzo | Picometer detection by adaptive holographic interferometry in a liquid-crystal light valve[END_REF]. More recently, the same team showed that it is possible to detect phase modulation through long optical fibers by combining a CMOS camera to a phase-only LCOS SLM (Liquid Crystal on Silicon Spacial Light Modulator). The camera records the hologram and it is numerically printed on the SML to control the gain and phase of the hologram or to remove the DC-component.

Gain media

To reduce the response time of the holographic material, Jayet et al. investigated the use of gain media for AOI detection [START_REF] Jayet | Fast wavefront adaptive holography in Nd:YVO 4 for ultrasound optical tomography imaging[END_REF]. The method consists in using a pump beam to create a population inversion in the gain medium. When this medium is illuminated by the interference pattern, some atoms will relax to their ground state through stimulated emission. This results in a spatial modulation of the gain on which light can diffract. The great advantage of this approach is that the response time, given by the lifetime of the storing level, can potentially be much shorter than speckle decorrelation (some tens of µs for the transition at 1064 nm in Nd:YVO 4 crystals [START_REF] Jayet | Fast wavefront adaptive holography in Nd:YVO 4 for ultrasound optical tomography imaging[END_REF]).

Conclusion on acousto-optic imaging

In this chapter we have seen that ultrasound can be coupled to light in order to recover the local light fluence deep inside scattering medium through light modulation. Two main processes are responsible for the modulation of light: modulation of the refractive index and of the scatterers positions. Through these phenomena, tagged photons are created with some interesting properties. Tagged photons are modulated by the ultrasound and their frequency is shifted from the optical frequency by this amount. Their number is proportional to the local light fluence integrated over the ultrasound volume and the modulated light power is proportional to the square of the acoustic pressure.

This means that a system that can filter and count these tagged photons will be able to recover the local light fluence inside the scattering medium. However, we have also seen that this detection is not straightforward because of four reasons. First the frequency shift is extremely small compared to light frequency (∼ MHz ≡ 10 -6 nm at 800 nm). Second, the ratio of tagged and untagged photons is weak (typically 1 over 1000). Third, the amplitude and phase of the speckle wave-front are random which means that the ultrasound modulation adds incoherently on a large detector. Fourth, the speckle decorrelates with the motion of the scattering medium.

Several approaches can be used for this detection and they can be sorted into two categories: incoherent techniques that aim at filtering the tagged photons and coherent techniques that measure the amount of modulation. Incoherent techniques have the advantage that they are not sensitive to the random wave-front of the scattered light but they suffer from weak optical etendue or need to be kept at cryogenic temperature. Coherent methods, on the other hand, can increase the modulation depth through interferometric processes but the randomness of the wave-front prevents from integrating a lot of speckle grains on a large detector as it does not increase the SNR. To overcome this difficulty it is possible to process each speckle grain individually on a camera or to use self-adaptive holography techniques which adapt the wave-front of a reference beam to that of the scattered light.

Each method has its advantages and drawbacks and there is no consensus in the acousto-optic community regarding the best technique. The method of choice depends on the application, the equipment and the technical know-how of the team. In the following manuscript, the detection technique will be wave-front adaptive holography using photorefractive crystals, which will be described in more details in the next chapter.

Chapter's Keypoints

• When ultrasound propagates in an illuminated medium, light can be phase modulated at the frequency of the ultrasound.

• In a scattering medium, ultrasound modulates the refractive index as well as the scatterers positions. Both effects contribute to the modulation of light.

• The so-called tagged photons are only generated where light and ultrasound are present.

• The amount of tagged photons is proportional to the local light fluence integrated over the ultrasound volume and to the square of the ultrasound pressure.

• The detection of the tagged photons is the critical part of acousto-optic imaging:

• Few tagged photons are generated compared to the total number of photons (typically 1 over 1000;

• The frequency shift of tagged photons is extremely small compared to the optical frequency (a few MHz compared to 10 14 Hz at 800nm);

• The speckle wave-front has random amplitude and phase. The integration of several speckle grains thus does not increase the signal to noise ratio;

• The speckle decorrelates with the motion of the scattering medium. In biological tissues, the speckle decorrelation time is typically a millisecond.

• The detection of tagged photons can be performed using coherent or incoherent methods.

The former uses interferometric processes to recover the ultrasound modulation while the latter aims at filtering the tagged photons.

• Several techniques exist to reconstruct an image. The approach using focused ultrasound is the one used in this manuscript.

CHAPTER 3

Photorefractive detection of tagged photons Acousto-optic imaging allows for the recovering of optical contrast deep inside highly scattering medium. Through the interaction with ultrasound in a limited volume, some photons are frequency shifted along their propagation. The detection of these so-called tagged photons provides information about the local optical fluence at the focus of the ultrasound. However, the tagged photons filtering step is the most critical part of the technique. Several methods exist and they all have advantages and limitations. Because they allow the use of large area single detectors, adaptive wave-front interferometric techniques can time-resolve the propagation of the ultrasound. This feature makes them good candidates for the detection of the acousto-optic signal generated by pulsed ultrasound. Such methods require the formation of a hologram and a variety of holographic medium can be used.

In this chapter, I will describe the use of photorefractive crystals as such a medium to perform the detection of tagged photons. After describing the physics of the photorefractive effect, I will explain how such crystals can be used for AOI before presenting the experimental characterisation of some of these crystals.

The photorefractive effect

The photorefractive effect was discovered in 1966 at Bell Laboratories by Ashkin et al. who noticed that a focused laser beam could change the refractive index of LiNbO 3 and LiTaO 3 crystals [START_REF] Ashkin | Optically-induced refractive index inhomogeneities in LiNbO 3 and LiTaO 3[END_REF]. This non linear phenomenon was first considered an unwanted effect, but a better understanding of its mechanism revealed its potential for the recording of holograms or for optical memories. In the end of the 1990s, it was suggested that this effect could be used to detect vibrations of a rough surface. A detailed and extensive explanation of the physics of photorefractive materials can be founds in books from Yeh [START_REF] Yeh | Introduction to photorefractive nonlinear optics[END_REF] or Günter and Huignard [START_REF] Günter | Photorefractive materials and their applications[END_REF]. We will here give some theoretical bases to understand how this effect can be used for the detection of the acousto-optic signal.

Principle

The photorefractive effect can happen in materials that are both photoconductive and electrooptic. Several materials present these characteristics such as crystals or polymers [START_REF] Suzuki | High-sensitivity ultrasound-modulated optical tomography with a photorefractive polymer[END_REF]. Amongst photorefractive crystals, several families exist:

• ferroelectric crystals such as lithium niobate (LiNbO 3 ), barium titanate (BaTiO

Photoconduction

Most of the crystals mentioned above are transparent in a given spectral range. For photoconduction to occur, photons need to be absorbed by the material. Consequently, photorefractive crystals are often doped or present intrinsic defects so that some energy levels can be found in their forbidden band. These defects are often called acceptor or donor centres as they can create a charge carrier by accepting or donating an electron. When a photon has a energy which matches that of one the centres, it can induce a photoionization and a charge carrier is released in a conductive band. Through a diffusion process, this carrier can travel in the crystal and recombine in another location. If the illumination is inhomogeneous only the centres in regions of high intensity will be excited and they will recombine elsewhere. Eventually, the charge carriers will accumulate in regions of low light intensity. The time needed to arrive in a steady-state is called the photorefractive response time. The spatial distribution of light thus creates a spatial distribution of charges which will then create a local electric field distribution. This simple model is called the band transport model [START_REF] Kukhtarev | Holographic storage in electrooptic crystals. I. steady state[END_REF].

Electrooptic effect

The photoinduced electric field will in turn generate a modulation of the refractive index through the Pockels effect. The interference pattern which illuminates the crystal is thus recorded as a refractive index grating inside the crystal. Figure 3.1 summarizes the photorefractive effect.

The band-transport model

Space-charge field

The refractive index variation in photorefractive crystals is usually derived from the band transport model. This model was initially developed by Kukhtarev et al. [START_REF] Kukhtarev | Holographic storage in electrooptic crystals. I. steady state[END_REF] and is inspired by photoconduction models in semiconductors. The movement of charges in the material is considered possible because of the diffusion of free charge carriers in the conduction band or the valence band as shown on figure 3.2. With this model it is possible to write four electrostatic equations that drive the creation of a space-charge field:

• the creation-recombination equation of electron donors and acceptors centres which describes the density of these centres;

• the conduction equation which describes the movement of electrons in the conduction band;

• the equation of charge conservation;

• the Poisson equation which express the electric potential created by a given distribution of charges.

For a given illumination, this system allows for the derivation of the space-charge field as a function of space and time. There is no analytical solution in the most general case but a simple expression exists for a quasi-continuous sinusoidal illumination with a weak modulation corresponding to the interference of two plane waves with very different intensities. These approximations do not represent our problem as the scattered light has a random wave-front, but since we can decompose the speckle field as a multitude of plane waves, it is interesting to study this simple case.

Let us consider two plane waves with amplitudes E R and E S and wave vectors k R and k Swhere R subscript corresponds to the reference beam and S to the signal beam. Let us assume that the modulation is weak and that the signal beam has a very low amplitude compared to the reference beam. We can then define the modulation depth, m,

m = E * R E S |E R | 2 + |E S | 2 ∼ E S E R ≪ 1. (3.1)
The fluence resulting from the interference of the two beams is written:

I(r) = I 0 1 + ℜ (me i∆k•r ) , ( 3.2)
where ℜ denotes the real part, ∆k = k R -k S is the difference of wave vectors and I 0 is the average value of the fluence:

I 0 = |E R | 2 + |E S | 2 ∼ |E R | 2 . (3.3)
The photorefractive effect is a macroscopic effect which involves a great number of charge carriers. Consequently, the establishment of the space-charge fields is not instantaneous and a steady-state is reached after a characteristic time called the photorefractive response time, τ PR . We will here calculate the steady-state space-charge field and thus ignore the time dependency. The photorefractive response time will be presented later. By inserting the fluence formula into the aforementioned system of equations it is possible to obtain the expression of the space-charge field E sc , which is also sinusoidal

E sc = E 0 + ℜ E 1 e i∆kr , ( 3.4)
where E 0 is an external electric field potentially applied to the crystal and E 1 results from the inhomogeneous illumination. In the diffusive regime, when no external field is applied, E 1 is given by

E 1 = -im k B T e ∆k 1 + ∆k 2 k 2 D , ( 3.5)
where k B is the Boltzman constant, T the temperature, e is the electron charge and k D is the Debye wave-vector and is expressed as

k 2 D = e 2 N A ǫk B T , (3.6)
where ǫ is the permittivity of the material and N A is the effective density of electrons acceptors.

In the theory of electrodynamics in plasmas, the electrons screen out the Coulomb field of a test charge over a distance of Γ d = 2π/k D . This distance is known as the Debye screening radius. It can be noted here that the local electric field is maximum for ∆k = k D , i.e. for a step between fringes equal to the Debye screening radius. One can also notice that the space-charge field is purely imaginary which means that it is phase shifted by ±π/2 with respect to the illumination pattern.

In practice, if the interference pattern is moving, because of the crystal response time, the photoinduced electric field will be proportional to the average value of the modulation depth over this time:

E 1 = -i m τ PR k B T e ∆k 1 + ∆k 2 k 2 D , ( 3.7)
where . τ PR stands for the temporal average over a time τ PR .

Photoinduced variations of the refractive index

Since the crystal is also electrooptic, the photoinduced electric field will create a modulation of the refractive index of the crystal which will strongly depend on the characteristics of the crystal and the illumination. To describe the optical properties of a crystal, it is common to use the permittivity tensor ǫ. In the most general case, there exist three axes denoted a,b and c, so-called the principal axes of the crystal. In the basis formed by these axes, the permittivity tensor can be written as:

ǫ =     ǫ a 0 0 0 ǫ b 0 0 0 ǫ c     =     n 2 a 0 0 0 n 2 b 0 0 0 n 2 c     (3.8)
The impermeability tensor, η is defined as the inverse of the permittivity tensor such as ηǫ = ǫ 0 where ǫ 0 is the permittivity of vacuum. The diagonal elements of the impermeability tensor are related to the diagonal elements of the refractive index by η ii = 1/n ii where i stands for one of the three principal axes of the crystal.

In the following we will only consider non-centrosymmetric crystals in which the Kerr effect can be neglected. In this case, the effect of a local electric field on the impermeability tensor can be written

∆η ij (E sc ) = k r ijk E sc k , ( 3.9)
where i,j and k stand for either one of the principal axes of the crystal, r ij are the elements of the electrooptic tensor and E sc k is the component of the space-charge field along the principal axis k. Since the tensor η is symmetrical (η ij = η ji ) so is the electrooptic tensor and it is possible to use a different notation with only two indices for r. The coefficients of the tensor are usually denoted r lk where k = 1, 2 or 3 for a, b and c and l varies between 1 and 6 so that:

1 = (aa); 2 = (bb); 3 = (cc) 4 = (bc) = (cb); 5 = (ac) = (ca); 6 = (ab) = (ba) (3.10)
The electrooptic tensor may then be written

r =            
            (3.11)
The electrooptic tensors for different groups of symmetry are given in the first chapter of Yeh [START_REF] Yeh | Introduction to photorefractive nonlinear optics[END_REF]. I will present here two crystals that were used for the work described in this manuscript: SPS and ZnTe. The former has been studied in detail by Grabar et al. in the book of Günter and Huignard [START_REF] Günter | Photorefractive materials and their applications[END_REF]. SPS (Sn 2 P 2 S 6 ) is a monoclinic crystal with an electrooptic tensor containing 10 independent coefficients:

r SPS =             r 11 0 r 13 r 21 0 r 23 r 31 0 r 33 0 r 42 0 r 51 0 r 53 0 r 62 0             (3.12)
These coefficients were measured by Grabar [START_REF] Günter | Photorefractive materials and their applications[END_REF] and their values range from 25 to 174 pm.V -1 in the wavelength range of 633-1300 nm, the strongest coefficients being r 11 and r 31 .

ZnTe is a crystal with a cubic symmetry and only one electrooptic coefficient:

r ZnTe =             0 0 0 0 0 0 0 0 0 r 41 0 0 0 r 41 0 0 0 r 41             (3.13)
This coefficient has a typical value of 4 pm.V -1 in the optical therapeutic window [START_REF] Ziari | Photorefractivity in vanadium-doped znte[END_REF]. Since this value is notably smaller than for SPS, we can expect that the efficiency of the grating formation will be similarly weaker.

The modulation of the refractive index is given by:

∆n ij = - 1 2 n 3 0 ∆η ij = - 1 2 n 3 0 k r ijk E sc k . (3.14)
In general it can be shown that the two beams see a modulation of the refractive index that can be expressed as [START_REF] Montemezzani | Light diffraction at mixed phase and absorption gratings in anisotropic media for arbitrary geometries[END_REF]:

∆n = - 1 2 n 3 0 r eff ℜ E 1 e i∆k•r (3.15)
where r eff is the effective electrooptic coefficient seen by the beams. This depends on the crystal characteristics and its orientation, on the fringe wave-vector direction and on the polarization of the two beams. The crystals used in this work are cut so that the sign of the effective electrooptic coefficient can be changed by π-rotating the crystal. The refractive index variation can then be expressed as:

∆n = n 1 sin (∆k • r) , (3.16)
where

n 1 = -i r eff n 3 0 2 |E 1 | = -m τ P R r eff n 3 0 k B T 2e ∆k 1 + ∆k 2 k 2 D . (3.17)
As mentioned previously, this effect is maximum when ∆k ∼ k D which corresponds to a given step length of the interference fringes and thus to a given angle between the signal and reference beam. In SPS and ZnTe, the Debye length is of the order of 1 µm which corresponds to a 40 • angle between the two beams.

Characteristics of the photorefractive effect

With this simple description of the photorefractive effect, it is possible to outline several main properties.

Memory effect reversibility

The recorded hologram can be erased by illuminating the crystal with a spatially homogeneous light. On the other hand, if the crystal is kept away from light after the formation of the refractive index grating, the hologram can persist for a time which depends on the conductivity of the crystal in the dark and can reach up to several months in the case of some ferroelectric crystals like LiNbO 3 [START_REF] Brouard | Cristaux photoréfractifs à long effet memoire pour l'interférométrie holographique aux longueurs d'onde des diodes laser[END_REF].

Spectral bandwidth

The photorefractive efficiency depends on the wavelength and the characteristics of the material but it is a non resonant effect. A broad range of wavelengths can thus be used.

Phase shift Equation 3.7 shows that the photoinduced electric field is purely imaginary and thus phase shifted by π/2 with respect to the illumination pattern. This impacts the refractive index modulation which is phase shifted by ±π/2 depending on the sign of r eff .

Response time As mentioned before, the movement of the great number of charge carriers is not instantaneous and the establishment of the refractive index modulation takes a time, τ PR , called the photorefractive response time. In the case of a low modulation (m ≪ 1), the system of four equations mentioned in section 3.1.2 gives a solution in the transition regime of the form [START_REF] Yeh | Introduction to photorefractive nonlinear optics[END_REF]:

E 1 (t) = E 1 1 -e -t/τ PR , ( 3.18)
where E 1 is the space-charge field in the steady-state regime. The photorefractive response time is given by:

τ PR = ǫ eµN 0 1 + k B T µk 2 eγ R N A 1 + ∆k 2 k 2 D , ( 3.19)
where µ is the effective mobility of charge carriers, N 0 is the density of electrons, N A is the density of electron acceptors and γ R the electron-hole recombination rate.

This characteristic time thus strongly depends on the material but it also depends on the illumination through the wave-vector of the interference fringes and through the density of electrons, N 0 which is proportional to the intensity of the illumination. The response time can thus be reduced by increasing the energy of the illumination.

Detection of the acousto-optic signal

Two-wave mixing

When the signal and reference beams interfere in a photorefractive crystal, they create a thick index grating of wave vector ∆k = k S -k R on which they will instantly diffract. Since these beams created the grating, the Bragg conditions are automatically verified and only one diffraction order exist such as k d = k i + ∆k. This means that the reference will diffract in the direction of the transmitted signal and reciprocally as shown on figure 3 In the case of small modulation, the diffraction on an index grating will phase shift the beams by π/2 with respect to the incident beam. In addition, since the index grating is phase-shifted by ±π/2 with respect to light, it means that the total phase shift between the diffracted and transmitted beams propagating outside of the crystal in the same direction is either 0 or π depending on the sign of the phase shift. These two beams will thus interfere constructively or destructively. This results in an energy transfer from a beam to the other that can be modelled through coupled equations and solved by a perturbation approach as developed by Kogelnik [START_REF] Kogelnik | Coupled Wave Theory for Thick Hologram Gratings[END_REF]. The phenomenon is called Two-Wave Mixing (TWM) and can be used to detect the acousto-optic signal.

Let us calculate the diffraction efficiency using the model of Kogelnik. If we consider the two previous beams propagating along the z-direction with a symmetrical incidence, the two coupled equations that drive the evolution of the beams inside the crystal are given by [START_REF] Yeh | Introduction to photorefractive nonlinear optics[END_REF]:

         ∂E S ∂t = πn 3 0 r eff 2iλ cos θ E sc m τ PR E R - α 2 E S , ∂E R ∂t = πn 3 0 r eff 2iλ cos θ E sc * m τ PR E S - α 2 E R , (3.20)
where 2θ is the angle between the beams and α is the absorption coefficient of the crystal. These equations define the photorefractive amplitude gain as:

γ = πn 3 0 r eff 2iλ cos θ E sc . (3.21)
This gain is a priori complex in the general case where an external field is applied but it appears that in the particular case of the diffusion regime, it is a real variable. If we assume that the temporal evolution of the interference pattern is slow compared to the response time of the crystal we can remove the mean value. In the case of slow modulation, when the reference is much more powerful than the signal (m ≪ 1), we have m ∼ E S E R and:

         ∂E S ∂t = γ - α 2 E S , ∂E R ∂t = - α 2 E R . (3.22)
In terms of fluence we can write:

     I S (z) = I S (0)e (2γ-α)z , I R (z) = I R (0)e αz , ( 3.23)
where I S (z) and I R (z) represent the fluence of the signal and reference beam respectively. From this we can define the photorefractive gain in energy, γ PR = 2γ. For high values of the average fluence on the crystal, this gain is a constant called the saturated gain. For lower values, it evolves with the average fluence in the same manner as the photoinduced electric field. We can see here that, by controlling the sign of γ PR = 2γ through the sign of r eff , it is possible to either amplify or attenuate the signal fluence. Since the signal is assumed to be weak compared to the reference, the change in the reference intensity is neglected.

Acousto-optic signal detection

In the previous section we only considered a static interference pattern in order to highlight the amplification of the signal beam. In acousto-optic imaging, the information is transmitted through a phase modulation of the signal beam. The complete study of the transmission of time-modulated signal through a photorefractive crystal has been published by Delaye et al. [START_REF] Delaye | Transmission of time modulated optical signals through an absorbing photorefractive crystal[END_REF], we will only describe the main steps of the derivation here. Let us call, ϕ(t) the phase modulation induced by the ultrasound. We can then write the signal beam at the entrance of the crystal (z = 0) as:

E S (x, y, 0, t) = E S (x, y, 0)e iϕ(t) (3.24)
To understand the evolution of this field in the crystal we will use the dynamic equations that drive the evolution of the signal and reference beam. The dynamic photoinduced electric field obeys three equations. The temporal evolution derived from the band transport model:

∂E 1 ∂t + E 1 τ PR = m E sc τ PR (3.25)
and the coupled equations 3.20 reformulated:

∂E S ∂t = γ E sc E 1 E R - α 2 E S (3.26) ∂E R ∂t = γ E sc E * 1 E S - α 2 E R (3.27)
Delaye et al. [START_REF] Delaye | Transmission of time modulated optical signals through an absorbing photorefractive crystal[END_REF] solved these three equations under four hypotheses:

• the reference beam is not time-modulated;

• the attenuation of the reference beam by TWM is neglected, which assumes that the photorefractive gain is weak and that the ratio between the signal and reference amplitude is very small (m ≪ 1);

• the reference beam is attenuated by the absorption of the crystal;

• at t = 0, before any time-modulation, a steady-state index grating already exists and the beams can be diffracted already.

In the case of a phase modulation that is fast compared to the photorefractive time, by deriving equation 3.27 with respect to time and injecting equations 3.26 and 3.25 before taking the Laplace transform we can write: E S (x, y, L, t) = e -αL/2 E S (x, y, 0)e iϕ(t) + (e γL -1)E S (x, y, 0) .

(3.28)

This equation means that, after the crystal, the beam propagating in the direction of the transmitted signal is the sum of two fields with the same spatial variations: the transmitted signal and the diffracted reference. This phenomenon is called wave-front adaptation and it allows for a coherent sum, without the random phase of the speckle. This gives photorefractive detection a great optical etendue only limited by the size of the crystal.

A single-detector placed in the direction of this beam will record an fluence given by: I S (x, y, L, t) = e -αL I S (x, y, 0) 1 + (e γL -1) 2 + 2(e γL -1) cos(ϕ(t)) (3.29) which, in the case of a weak phase modulation, ϕ(t) ≪ 1 can be approximated to [START_REF] Delaye | Transmission of time modulated optical signals through an absorbing photorefractive crystal[END_REF]:

I S (x, y, L, t) ≃ e -αL I S (x, y, 0) e 2γL -(e γL -1)ϕ(t) 2 (3.30)
This shows that the intensity directly depends on the quadratic phase of the signal wave. It is important to highlight that the phase modulation created by the ultrasound needs to be faster than the response time of the crystal. Since the refractive index modulation is proportional to m τ PR , the modulation depth averaged over the response time of the crystal, the grating does not carry any information about the ultrasound modulation in this case. This information is only present in the signal beam and not in the diffracted reference so that when they interfere, the cross term contains the ultrasound modulation. In the opposite case, when the ultrasound modulation is slow compared to τ PR , the hologram has time to adapt and the diffracted reference thus also contains the ultrasound modulation. Consequently, when the transmitted signal and diffracted reference interfere, the cross term disappears and the phase modulation information is lost.

Experimental configuration

The theoretical analysis of the photorefractive effect in crystals shows that, provided that the ultrasound modulation is faster than the response time of the crystal, the acousto-optic signal can be detected using this effect. However, it also shows that this detection can be influenced by a lot of parameters, some intrinsic to the crystal properties and others which only depend on the illumination. It follows that several configurations can be used to detect the acoustic modulation of light. I will here describe the configuration used in the work presented in the following of this manuscript.

Detection of untagged photons

Ultrasound propagates in tissues at a speed of approximately 1500 m.s -1 , which means that it takes only a few tens of microseconds to travel through several centimetres. Most photorefractive crystals have a response time of some milliseconds or more. Consequently, if the reference is tuned on the tagged photons, the hologram buid-up is too slow to follow the propagation of the ultrasound in the medium. As stated before, by energy conservation, it is possible to infer the number of tagged photons from the number of untagged photons and a measurement of the untagged photons without the ultrasound. Since the untagged photons exist without the ultrasound, the corresponding hologram can be recorded before turning them on. When the acoustic wave travels through the medium, the number of untagged photons decreases and the amount of modulated light can be deduced.

Negative gain configuration

When propagating through the photorefractive crystal, part of the reference beam can be scattered by defects, impurities or other inhomogeneities of the crystal and create a spatial distribution of light in the crystal. Since this light also undergoes TWM, it can be strongly amplified in the case of a positive gain and lead to the creation of self-induced index grating in the medium which results in an amplification of the noise. This phenomenon, called photorefractive beam fanning [START_REF] Yeh | Introduction to photorefractive nonlinear optics[END_REF][START_REF] Günter | Photorefractive materials and their applications[END_REF] is the main source of noise in photorefractive-based AOI and for this reason it is interesting to work in a negative gain configuration in order to reduce the influence of the parasitic self-induced holograms.

Equation 3.30 shows that the signal detected on the photo-diode contains a static part, e 2γL , and a time-varying one, (e γL -1)ϕ(t) 2 , which contains the information of interest. In the case of a positive gain, the static term is strong and the modulation will decrease the overall signal through the creation of tagged photons. With a negative gain, the static term will be attenuated and the modulation will result in an increase of the overall signal detected. Figure 3.4 shows a representation of both cases. In the following manuscript, because of the high photorefractive gain of SPS, a negative gain configuration was used to attenuate beam fanning.

Choice of the crystal

Photorefractive crystals exist with a wide range of properties. The choice of the crystal will then depend on the application. In the case of acousto-optic imaging, the four main properties of a photorefractive crystal are its optical etendue, which only depends on its dimensions; its photorefractive gain; its response time; and its absorption at the wavelength of interest. The latter three depend on the nature of the crystal and its doping. A high photorefractive gain enhances the energy transfer in the two-wave mixing process and thus increases the signal to noise ratio. The response time, in the scope of medical imaging, needs to be as short as possible and smaller than the speckle decorrelation time (∼ 1 ms) in order to image dynamic tissues. Looking for the best crystal for AOI, several photorefractive crystals were studied by Farahi et al. [START_REF] Farahi | Spectroscopic and holographic filtering applied to acousto-optic imaging of thick scattering media[END_REF]. The conclusion was that the best crystal working around 800 nm was SPS (Sn 2 P 2 S 6 ) because of its high photorefractive gain at this wavelength and a response time smaller than 10 ms. These crystals exhibit a high photorefractive gain but their response time is still greater than the speckle decorrelation time in biological media. Whilst it has been shown that it is possible to image through a living mouse ear using SPS [START_REF] Liu | Optical focusing deep inside dynamic scattering media with near-infrared time-reversed ultrasonically encoded (true) light[END_REF], it is interesting to investigate faster crystals in the pursuit of in vivo imaging.

A crystal of zinc telluride doped with manganese and vanadium (ZnTe:V,Mg) was studied in order to assess its potential use for fast acousto-optic imaging. After a description of the experimental measurement method, the obtained results will be presented and compared to similar measurements made on SPS crystal by Laudereau et al. [START_REF] Laudereau | Acousto-optic imaging : challenges of in vivo imaging[END_REF]

Measurement of the crystals characteristics

The two principal characteristics of the crystal regarding TWM are the gain and the response time. Both depend on the average fluence of the crystal as shown by equations 3.19 and 3.21. In order to measure their values against the average fluence, the same protocol was repeated for different values of light power.

Photorefractive Gain

Recalling equation 3.30, we can see that in the case of TWM without any ultrasound modulation (ϕ(t) = 0) the outgoing signal fluence is equal to e γ PR L-αL times the input fluence. It follows that to determine γ PR , one can measure I s (L) with and without TWM (γ PR = 0) and calculate the ratio.

To do so, a simple interference set-up with two plane waves recombining on the photorefractive crystal was used. In order to respect the hypotheses used for the derivations in the previous sections, the power of the signal beam was very weak compared to the reference one. A square amplitude modulation was applied to the reference so that its fluence ranged from 0 to E R and the signal beam intensity was measured on a photodiode for the two extreme values.

Photorefractive response time

Measuring the response time of the crystal is less straightforward. The measurement is possible in the time domain by evaluating the build-up of the photorefractive grating through time resolved measurements of the variation of the detected signal. This method is convenient for slow crystals but when fast response times have to be determined, it requires the use of high-speed optical shutters or pulsed lasers. While this is theoretically simple, it can be difficult to implement experimentally. Consequently, I used a technique similar to [START_REF] Lesaffre | In situ monitoring of the photorefractive response time in a self-adaptive wavefront holography setup developed for acousto-optic imaging[END_REF] to measure the response time in the frequency domain.

Lesaffre et al. developed the method reported in [START_REF] Lesaffre | In situ monitoring of the photorefractive response time in a self-adaptive wavefront holography setup developed for acousto-optic imaging[END_REF] to measure the response time in situ with a scattering medium and ultrasound modulation. For the sake of simplicity, plane waves and no ultrasound were used but the reasoning is similar. This method relies on the creation of a moving hologram inside the crystal through the interference between two beams with slightly different frequencies. When the beating of the interference pattern is too fast, the crystal cannot follow and the hologram disappears but when the detuning is small enough, the hologram is blurred and loses contrast. By varying the detuning it is thus possible to measure the frequency of the grating formation inside the crystal. Figure 3.5 shows the experimental setup used for this measurement. Let us consider two plane waves E S and E R which will be referred as signal and reference beam respectively. The intensity of the signal beam will be assumed weak compared to the reference. A frequency detuning ∆ω is applied to the reference beam as well as an amplitude or phase modulation. Both modulations lead to similar results and we will only consider a phase modulation φ(t) here.

E S = |E S | exp(iω 0 t), (3.31) E R = |E R | exp[i(ω 0 + ∆ω)t + iφ(t)]. (3.32)
As we have seen before, the interference of these beams will create an index grating inside the crystal. For small diffraction efficiency and negligible absorption, De Montmorillon et al. [START_REF] De Montmorillon | Novel theoretical aspects on photorefractive ultrasonic detection and implementation of a sensor with an optimum sensitivity[END_REF] show that the transfer function H of the crystal in terms of amplitude can be written

H(L, t) = γL exp(-t/τ PR ) τ PR , for t > 0, (3.33)
where γ is the photorefractive gain, L is the crystal thickness and τ PR is the photorefractive response time. This means that the crystal acts as a low-pass filter with a cut-off frequency of 1/τ PR .

The time varying index modulation inside the PRC is then proportional to the modulation depth of the interference pattern convoluted by the transfer function of the crystal. We thus have

∆n ∝ E S E * R |E S | 2 + |E R | 2 ⊗ γL exp(-t/τ PR ) τ PR , ( 3.34)
where ⊗ denotes the convolution product with respect to time. Since the signal is assumed weak compared to the reference, this can be simplified as

∆n ∝ E S E R ⊗ γL exp(-t/τ PR ) τ PR = |E S | |E R | e -i∆ωt-iφ(t) ⊗ γL exp(-t/τ PR ) τ PR . ( 3.35)
The reference can then diffract on the grating and give the diffracted beam E D :

E D = E R ∆n ∝ |E S | e -i∆ωt-iφ(t) ⊗ γL exp(-t/τ PR )
τ PR e i(ω 0 +∆ω)t+iφ(t) .

(3.36)

In order to express the convolution product, let us suppose that the phase modulation is a sinusoid of pulsation Ω:

φ(t) = ϕ m sin(Ωt). (3.37)
The results would be similar with a periodic square modulation but the derivation is more complicated. The choice of a sine function allows for the decomposition of the modulation term as a Bessel series:

e -iϕm sin(Ωt) = n=+∞ n=-∞ J n (ϕ m )e -inΩt , (3.38)
where the J n are Bessel functions of the first kind and n is the order of the Bessel function.

If we use the fact that the Fourier transform of a convolution product is equal to the product of the Fourier transforms, equation 3.40 can be written as:

E D ∝ |E S |F T -1 F T n=+∞ n=-∞ J n (ϕ m )e -i(nΩ+∆ω)t F T γL exp(-t/τ PR )
τ PR e i(ω 0 +∆ω)t+iϕm sin(Ωt)

(3.39)

∝ |E S |F T -1 n=+∞ n=-∞ J n (ϕ m )δ(ω + (nΩ + ∆ω))) γL 1 + iωτ PR e i(ω 0 +∆ω)t+iϕm sin(Ωt) , ( 3.40)
where F T and F T -1 denotes the Fourier transform and the inverse Fourier transform respectively and δ(ω -ω ′ ) is the Dirac functional at position ω = ω ′ . If we then note that:

F T -1 [δ(ω -ω ′ )f (ω)] = δ(ω -ω ′ )f (ω)e iωt dω = f (ω ′ )e iω ′ t , ( 3.41)
equation 3.40 can be written

E D ∝ |E S |γL n=+∞ n=-∞ J n (ϕ m ) e -i(nΩ+∆ω)t 1 -i(nΩ + ∆ω)τ PR e i(ω 0 +∆ω)t+iϕm sin(Ωt) , (3.42) ∝ |E S |γL n=+∞ n=-∞ J n (ϕ m ) e -inΩt 1 -i(nΩ + ∆ω)τ PR e iω 0 t+iϕm sin(Ωt) . (3.43)
The diffracted reference beam thus contains an infinite number of terms oscillating at harmonics of the modulation frequency Ω. To measure the time response of the crystal we only need to extract one of these components, for example by using a lock-in detection tuned on the modulation frequency. In this case we can simplify the expression above by developing the term e iϕm sin(Ωt) using equation 3.38 and by keeping only the terms oscillating at the frequency Ω:

E D ∝ |E S |γL J 1 (ϕ m )J 0 (ϕ m ) e -iΩt 1 + i(Ω -∆ω)τ PR + J 0 (ϕ m )J 1 (ϕ m ) e iΩt 1 -i(∆ω)τ PR + (3.44) J 0 (ϕ m )J -1 (ϕ m ) e -iΩt 1 -i(∆ω)τ PR + J -1 (ϕ m )J 0 (ϕ m ) e iΩt 1 -i(Ω + ∆ω)τ PR   e iω 0 t .
This diffracted beam will interfere with the transmitted reference and the signal measured on the photodiode, S P D , can then be written

S PD = |E S + E D | 2 = |E S | 2 + |E D | 2 + 2ℜ(E * S E D ), (3.45)
where ℜ stands for the real part. The terms |E S | 2 and |E D | 2 are constants while the cross term oscillates at the modulation frequency (amongst other frequencies). The above expression can be written using the P and Q quadrature components of the signal which respectively correspond to the terms proportional to cos(Ωt) and sin(Ωt)

S Ω = P cos(Ωt) + Q sin(Ωt). (3.46)
Using the property of Bessel functions of the first kind: J -n (x) = (-1) n J n (x), one can write after calculations:

P (∆ω) = 2A 1 + (∆ωτ PR ) 2 + A 1 + (Ω + ∆ω) 2 τ 2 PR + A 1 + (Ω -∆ω) 2 τ 2 PR , (3.47) Q(∆ω) = A(Ω + ∆ω)τ PR 1 + (Ω + ∆ω) 2 τ 2 PR + A(Ω -∆ω)τ PR 1 + (Ω -∆ω) 2 τ 2 PR , (3.48)
where A is a constant such as A ∝ γL|E S | 2 J 0 (ϕ m )J 1 (ϕ m ). The lock-in detection can extract the magnitude, R(∆ω) = P (∆ω) 2 + Q(∆ω) 2 . If the modulation frequency, Ω is high enough compared with 1/τ PR , R(∆ω) exhibits three separated resonances for ∆ω = 0 and ±Ω, which can be used to measure τ PR . The central peak (∆ω = 0) has a Lorentzian shape with a full width at half maximum (FWHM) given by FWHM ω = 2/τ PR in terms of angular frequency, which correspond to a FWHM FWHM f = 1/πτ PR in terms of frequency. The other resonances are square root lorentizan with a FWHM given by FWHM ω √ 3/τ PR . Figure 3.6 shows an experimental example of R(∆f ) after a ZnTe crystal for modulation frequency of 20 kHz.

Frequency detuning (kHz)

SPS vs. ZnTe

Laudereau et al. used the experimental setup presented on Figure 3.5 to characterize three SPS crystals in order to choose the best one. These crystals have a square section of dimensions 1 cm × 1 cm and a thickness of 0.7 cm. The results of the measurements are shown on figure 3.7 extracted from [START_REF] Laudereau | Acousto-optic imaging : challenges of in vivo imaging[END_REF]. On can see that the photorefractive bandwidth (inverse of the response time) can vary a lot between crystals because of the variability in the doping and of the intrinsic defects.

The best crystal exhibits a high gain and a short response time which can be as short as 3 ms for a strong fluence. Despite the good performances of this SPS crystal, a response time of several milliseconds is still too high to perform in vivo imaging because of the speckle decorrelation time of the order of 1 ms in tissues. Several types of PRC with a high gain and short response time exist but most of them are efficient at wavelength that are not suitable for imaging biological tissues. According to literature, zinc telluride crystals are known to exhibit response times smaller than a millisecond [START_REF] Ziari | Photorefractivity in vanadium-doped znte[END_REF] for a range of wavelengths in the optical therapeutic window. We thus decided to analyse a ZnTe crystal doped with vandadium and manganese provided by Cradley Crystals Corporation to assess its potential use for in vivo imaging.

This crystal, cut along its 0 -10 , 001 and 110 crystallographic axes, measures 10×10×5 mm as shown on figure 3.8(a). The gain and the response time of this crystal were measured using the aforementioned protocols and the setup shown on figure 3.5. It consists in a semiconductor laser diode emitting at 790 nm, amplified by a 2 W tapered amplifier (Sacher laser technik, GmbH), forming a Master Oscillator Power Amplifier (MOPA). The amplified light is the split into two beams, the reference beam and the signal beam with a Polarizing Beam Splitter (PBS) and a half-wave plate (λ/2) to control the relative power of each beam. Both beams are sent through acousto-optic modulators (AOM) controlled by a function generator (Agilent -33250A). The signal beam is shifted by 79 MHz while the reference beam is shifted by 79 MHz plus a varying frequency ∆f and phase modulated at an angular frequency Ω far greater than the bandwidth of the crystal. After recombination on the crystal, the signal is collected on a photodiode (PD, Thorlabs PDA 36A-EC) connected to a lock-in detector (EG&G instruments, inc. 7260 DSP Lock-in Amplifier). The signal and reference beam propagated along the <110> axis of the crystal as it is the direction that exhibits the higher gain. shows the results of the gain measurement where we can clearly see the saturation regime which is already reached for an average fluence of 0.4 W.cm -2 . The photorefractive gain is noticeably weaker than that of the SPS crystals with a maximum value of 0.25 cm -1 . The results of the bandwidth and response time measurements are shown on figure 3.8(c). These results confirm that ZnTe crystals can have a response time smaller than the millisecond for reasonable fluence. For higher values of the fluence, the build-up time of the hologram in this crystal can be as low as 0.2 ms.

ZnTe crystals are thus very promising in terms of response time since it can be shorter than the speckle decorrelation time in biological tissues. Yet their photorefractive gain is much weaker than that of the SPS which makes imaging inside thick scattering medium more difficult because of the low SNR. However, M. Ziari et al. showed that some PRC such as vanadium-doped cadmium telluride, CdTe:V, exhibit an strongly enhanced photorefractive gain when an external electric field is applied [START_REF] Ziari | Enhancement of the photorefractive gain at 1.3-1.5 µm in cdte using alternating electric fields[END_REF]. With a voltage as high as 20 kV.cm -2 , they show that the gain can be multiplied by up to ten times its value without external field. Since ZnTe and CdTe crystals belong to the same crystallographic family, it is expected that such behaviour could be observed in ZnTe. These crystals can thus prove quite interesting for AOI purposes but in the absence of external field, SPS crystal are still better suited four our application.

Conclusion

The photorefractive effect is a phenomenon that appears in materials that are both photoconductive and electro-optic. Under an inhomogeneous illumination a spatial distribution of charges is created in such materials which then creates a local electric field distribution. This photo-induced electric field in turn creates a refractive index variation which can act as a grating for light. This effect can be used to create holograms in photorefractive crystals using a speckle beam and a reference so that the reference can then diffract on this grating with a wavefront adapted to that of the speckle beam. This process is called two-wave mixing and it allows for the detection of an acoustic modulation on the speckle beam which is at the heart of acousto-optic imaging using PRC. The main characteristics of this detection scheme are the photorefractive gain and the response time of these crystals. A high gain is needed to increase the SNR of the acousto-optic detection while a short response time is necessary to image dynamic media. The characterization of two types of crystals, Sn 2 P 2 S 6 and ZnTe:V,Mg showed that despite the short response time of the millisecond, which could prove useful for dynamic imaging of tissues, SPS crystals remain the most suited for thick scattering medium imaging because of their high photorefractive gain. Such crystals were thus used for all the experiments presented in the following chapters.

Chapter's Keypoints

• Photorefractive materials are both photoconductive and electrooptic. The refractive index of such materials can be modulated by light.

• A spatial variation of light irradiance induces a movement of charges which creates a local electric field and modulates the refractive index of the material.

• An interference pattern illuminating a photorefractive material can create a refractive index hologram.

• The refractive index modulation is optimal for a step between interference fringes equal to the Debye length. The efficiency of the modulation thus depends on the angle between the beams used to create the interference pattern.

• The recorded hologram can be erased by illuminating the material with a spatially homogeneous light.

• The hologram is recorded over a characteristic time called the photorefractive response time, τ PR . This time decreases with the average light intensity.

• When such an hologram is created by a signal beam and a reference beam, the reference is diffracted in the direction of the transmitted signal.

• An energy transfer can take place between the signal and reference beam. The sign of this transfer depends on the sign of the effective electro optic coefficient. This phenomenon is called two-wave mixing.

• This energy transfer is quantified by a quantity called the photorefractive gain. The gain increases with the average light intensity.

• The wave-front of the diffracted reference is adapted to the wave-front of the signal which allows for a coherent summation of the two fields.

• Two wave-mixing can be used to detect the modulation of light by the ultrasound.

• Two photorefractive crystals efficient around 800 nm have been characterized:

• SPS crystals have a high photorefractive gain and a response time of several milliseconds;

• ZnTe crystals can have response times as low as 0.2 ms but their gain is weak.

• Despite a response time greater than the speckle decorrelation time in vivo, SPS crystals are currently more suited for acousto-optic imaging because of their high gain. Such crystals will be used in the rest of this manuscript.

CHAPTER 4

Acousto-optic imaging in reflection mode Acousto-optic imaging has been mainly studied in transmission mode, where the source and the detector are on opposite sides of the medium. However, this configuration can be impractical for in vivo imaging when the biological tissues can be accessed only on a single side. Some work has been published in reflection mode, with the source and detector on the same side [START_REF] Gunadi | Optimization of the acousto-optic signal detection in cylindrical geometry[END_REF], but this configuration needs to be studied in more detail to be of use. In this chapter, I will present some results obtained numerically and experimentally in a reflection configuration in scattering media. After introducing the use of a conventional ultrasound scanner for AOI, I will present a Monte-Carlo algorithm which was used to obtain preliminary results. I will then present the experimental setup and the result obtained in reflection mode in the pursuit of in vivo imaging.

Fusion between conventional ultrasound and Acousto-Optic Imaging

Ultrasound imaging is a wide spread technique used routinely for a broad range of purposes, from medical imaging to non-destructive control. This is a fast, low-cost, non-invasive and portable technique and because of these many advantages, it is broadly available and is used routinely by physicians for diagnosis. This technique is used notably for prenatal monitoring or heart imaging but also for functional imaging [START_REF] Mace | Functional ultrasound imaging of the brain: theory and basic principles[END_REF][START_REF] Gesnik | 3d functional ultrasound imaging of the cerebral visual system in rodents[END_REF] and elastography [START_REF] Tanter | Ultrafast imaging in biomedical ultrasound[END_REF][START_REF] Bercoff | In vivo breast tumor detection using transient elastography[END_REF]. Commercial ultrasound scanners use multi-element transducers arrays of several tens or hundreds of small piezoelectric transducers [START_REF] Daoudi | Correcting photoacoustic signals for fluence variations using acousto-optic modulation[END_REF]192,256 or 512 elements for instance). Transducer arrays can have different shapes (linear, curved, 1D, 2D,...) depending on the application. Each element can be independently controlled electronically to transmit or receive at a chosen time. Using the Huygens-Fresnel law, a delay line can be derived to control the shape of the wavefront under the probe. For instance, if all elements transmit in phase, a plane wave will be generated whereas applying a linear delay will create a plane wave propagating with a tilted angle.

To form a conventional ultrasound image, parabolic delays are used to focus the energy to a given point in the domain. This can be done when emitting or receiving to create an image.

Let us consider a probe with N elements. To focus at a position [x f , y f ] a delay will be applied to the m elements centred around the position y f such that the delay for the j th element is:

∆t(j) = (y c -y m ) 2 + x 2 f -(y c -y j ) 2 + x 2 f c US (4.1)
where y c is the position on the probe of the element closest to y f , y m and y j are respectively the position on the probe of the last used element and the j th element and c US is the speed of the ultrasound in the medium. Elements on the edge emit first whereas the element closest to the center emits last. This ensures that all the wavelets arrive at the same time at the focus. See The same probes and the same kind of sequences can be used to perform acousto-optic imaging, but in this case, the transducers are only used for transmission since the detection is optical. Because the speed of ultrasound is small compared to the speed of light, given a fast detector, it is possible to follow the propagation of an ultrasound pulse travelling through the medium. This means that a single pulse can give information about a line in the medium; as the pulse travels, the number of tagged photons changes and is recorded in real time. Knowing the speed of ultrasound, it is then easy to reconstruct a 1D line. The longitudinal resolution will then be given by the most limiting factor between the pulse length and the detector bandwidth as explained in chapter 2.

Since it allows for 2D AOI without any movement of the probe, fusion between conventional ultrasound and AOI is of prime interest. Moreover, since it uses the same kind of acoustic sequences, it is possible to perform conventional ultrasound imaging and AOI at the same time thus combining the advantages of the two techniques.

All the work in this manuscript was performed with a commercial ultrasound scanner (Aixplorer, SuperSonic Imagine, Aix-en-Provence, France) coupled to a photorefractive-based acousto-optic imaging set-up using the Te:SPS crystal described in chapter 3.

This set-up was used for several purposes, including imaging ex-vivo liver tumors [START_REF] Laudereau | Multi-modal acousto-optic/ultrasound imaging of ex vivo liver tumors at 790 nm using a Sn 2 P 2 S 6 wavefront adaptive holographic setup[END_REF] and the development of ultrafast AOI with ultrasonic plane waves [START_REF] Laudereau | Ultrafast acousto-optic imaging with ultrasonic plane waves[END_REF] and proved to be quite powerful. However, all these experiments were performed in transmission mode, meaning that optical source and detector were placed on both sides of the medium of interest. Whilst this approach gives good results, one can easily see the impracticality of such a configuration when it comes to in-vivo imaging.

In order to image inside thick biological tissues, a more practical configuration would be to have all sources and detectors on the same side of the medium. With such a configuration, one can envision a multi-modal probe containing a transducer array as well as optical sources and detectors via the use of optical fibres. This handheld probe could allow physicians to add an AOI modality when performing conventional ultrasound imaging without the need to change hardware.

In order to study the feasibility of such a device, simulations and experiments were performed in a reflection mode, with both optical source and detector on one side of the sample. The AO signal was simulated by a custom Monte Carlo algorithm adapted to the detection of tagged photons. Experiments in the same configurations were then performed on phantoms to study the penetration depth and some in-vivo measurements were performed to assess the feasibility of the technique.

Monte Carlo Simulations

Before performing an experiment in a reflection configuration we used simulations to develop a more precise idea of how AOI could work in a reflection mode. To simulate the propagation of light in a scattering medium, several approaches can be considered. One can derive exact equations from the Maxwell equations or use the Radiative Transfer Equation (RTE) or the Diffusion Approximation (DA) to the RTE. However, finding an analytical solution to the RTE or to Maxwell's Equations is not always possible and it can be computationally expensive to find a numerical solution. In the case of the Diffusion Approximation to the RTE, it often requires too many approximations to give a good enough representation of the AO signal. This will be discussed in chapter 5. An alternative is to use a probabilistic method such as Monte Carlo algorithms.

Monte Carlo methods are a broad range of algorithms that rely on repeated random sampling to obtain numerical results. They are particularly useful for simulating systems with many coupled degrees of freedom. Rather than deriving an analytical solution for a system with such a high number of parameters, a great number of random events are simulated to obtain an empirical mean which will represent the solution.

In the case of light propagation inside a highly scattering medium, the repeated event is the movement of a single photon as it travels through the medium. Monte Carlo methods are particularly interesting for such a problem because of the random nature of scattering. When a photon propagates in a highly scattering medium, it encounters randomly positioned scatterers which will change its direction. It is impossible to know the position of every scatterer and to derive the trajectory of a photon but it is possible to randomly propagate a photon thousands of time in order to obtain a statistical estimate of the mean of the trajectories.

Monte Carlo methods have been extensively studied for light propagation [START_REF] Zhu | Review of monte carlo modeling of light transport in tissues[END_REF] and several open source algorithms are available. In the context of AOI, we decided to create our own based on the work of Romain Pierrat, written in Fortran 90.

The algorithm

The propagation of a photon is a classic random walk in three dimensions. After defining the size of the 3D domain, a photon packet is generated at a random position, r 0 on the surface of the source, with a given direction d 0 . A random step length is then generated with a Poisson distribution centred around the scattering mean free path and the new position is calculated.

As the medium is not isotropic, the generation of a new propagation direction needs to take into account the anisotropy factor, g. A classical approach is to use the Henyey-Greenstein phase function as proposed by Henyey and Greenstein in 1941 [START_REF] Henyey | Diffuse radiation in the Galaxy[END_REF].

P HG (θ) = 1 4π 1 -g 2 (1 + g 2 -2g cos(θ) 3/2 ) (4.2)
where P HG (θ) is a probability density function and θ ∈ [0, π] is the angle between the direction of a photon before and after a scattering event. In 1977, Witt showed that it is possible to generate a random angle θ which follows the Henyey Greenstein phase function using a uniformly distributed random number ξ ∈ [0, 1] [START_REF] Witt | Multiple scattering in reflection nebulae. i. a monte carlo approach[END_REF]. After calculation an analytical solution can be written:

cos(θ) =      1 2g 1 + g 2 -1-g 2 1-g+2gξ 2 , f org > 0, 2ξ -1, f org = 0, (4.3)
This derivation of a new direction is commonly used in Monte Carlo algorithms [START_REF] Toublanc | Henyey-greenstein and mie phase functions in monte carlo radiative transfer computations[END_REF][START_REF] Binzoni | The use of the henyey-greenstein phase function in monte carlo simulations in biomedical optics[END_REF]. After generating the new direction, d i , the photon packet is moved with a random step length following a Poisson distribution, and the process is then repeated until the photon packet exits the medium. The new position, r i+1 can be written:

r i+1 = r i -l s log(1 -α)d i α ∈ [0, 1] (4.4)
where r i is the previous position, l s is the scattering mean free path and α is a random scalar between 0 and 1.

To simulate the acousto-optic tagging of the photons, a parallel method was used. Rather than doing a full simulation for each ultrasound tagging location, a given number of tagging volumes is evenly distributed in the medium and a counter is created for each of them. The volumes are 3D ellipsoids whose three dimensions can be chosen independently. For each scattering event, the algorithm checks if the path of the photon packet crossed any of these areas and updates the counter corresponding to the crossed areas. If the same photon packet crosses an area twice, the counter is switched back to zero. To take absorption into account, each photon packet starts with a given energy. For each scattering event, the energy of the photon packet is updated by taking into account the different areas it travelled through. If a photon packet moves by a distance s i including a distance s incl inside an absorbing inclusion, its new energy, E i+1 , would be:

E i+1 = E i exp(-(s i -s incl )µ 0 a -s incl µ incl a ) (4.5)
where E i is the previous energy of the photon packet, µ 0 a is the absorption coefficient of the medium and µ incl a is the absorption coefficient of the inclusion. If the photon packet exits on the detector, its remaining energy is added to the signal corresponding to every area it crossed. If the energy of a photon packet falls below a chosen limit, ǫ, this photon packet is considered absorbed by the medium and a new photon packet is created.

When all photons packets have been simulated, the signal matrix is normalized by the total incident energy. Figure 4.3 summarizes this whole process.

Influence of source-detector distance

The algorithm was used to simulate the distribution of light in a reflection mode configuration to study the influence of the distance between source and detector. The physical properties of the medium were chosen to be close to experimental values. An 80 × 80 × 40 mm homogeneous domain was defined with a reduced scattering coefficient µ ′ s0 = 10 cm -1 , an absorption coefficient µ a0 = 0.02 cm -1 and an anisotropy factor g = 0.9. To reproduce the experimental set-up that will A photon packet is generated on the source and propagates through the medium in a random walk fashion. The average step length is the scattering mean free path. Only the data relative to photons exiting the medium on the detector will be saved. Small light ellipsoids represent the ultrasound tagging areas. Every time a photon packet crosses one of these areas, a counter is incremented. Each pixel of the final image will correspond to one such area, the amount of signal for each pixel being the total energy of photons which crossed the corresponding area. For the sake of clarity, only a 3x3 pixel grid was represented here. In this image, the counter of the blue ellipsoid will be incremented whereas the counter of the red one will not because the red photon packet did not reach the detector.

be presented in the next section, the source had a 1 cm diameter and the detector diameter was 0.8 and a numerical aperture of 0.4. The ultrasound tagging areas were chosen to be spheres with a 1 mm diameter which is comparable to experimental focus size. By using 81 × 41 spheres over an 80 × 40 cm square, the whole 2D space was covered. Simulations were performed on a cluster of 192 CPUS and for each simulation, 96 10 6 photons were sent. The source and detector were then moved away from each other to obtain the images presented in figure 4.4.

As expected, the region explored by the photons inside the medium has a characteristic "banana" shape which gets bigger as the distance between the source and detector increases. Red arrows indicate the source and detector positions and it appears that the maximum intensity of the AO signal is not exactly in front of the arrows. This can be explained by the spatial extension of source and detector. Both source and detector are large and when the distance between the two becomes closer to the transport mean free path, more photons appear to be taking the shortest path. These images can be seen as probability maps, the more intense locations represent the areas more likely to be travelled through by photons which exit on the detector. This means that a photon which enters the medium on the edge of the source closer to the detector is morel likely to exit on the detector. By reciprocal propagation of light, it also means that it is more likely to exit on the edge of the detector. This appears to be particularly true when source and detector are close. Since source and detector diameters are slightly different the images are not symmetrical. As the distance increases, the intensity of the signal decreases a lot. A photon is less and less likely to come back to the detector. More noise also appears because of the smaller number of photons reaching the detector.

In order to evaluate the penetration depth of light, profiles of the images were taken on the A moving average over ten points was applied and because the signal intensity varies drastically with the source-detector distance, the profiles were normalized by their maximum values for representation purposes. One can seen that the explored region, as well as the penetration depth, increases with the distance between source and detector but the signal gets noisier because of the decrease in the number of photons reaching the detector.

To explore the relation between penetration depth and source-detector distance, the maximum of these profiles were plotted on figure 4.6(a). A linear relationship can be observed with a coefficient 0.35. In a paper by Feng et al. [START_REF] Feng | Photon migration in the presence of a single defect: a perturbation analysis[END_REF], the authors derive a theoretical model of the banana shape in an infinite and semi-infinite medium. In the low absorption regime, they find that the maximum intensity of the banana is reached for a depth equal to √ 2 4 d ≃ 0.36d, where d is the distance between source and detector, which correspond to the value obtained with simulation. An empirical value of d/2 is mentioned in several papers such as [START_REF] Shi | Non-contact reflection photoplethysmography towards effective human physiological monitoring[END_REF] but the origin of this value is often not explained.

If it seems possible to image at a few centimetres depth by increasing the distance between source and detector, it comes with a trade-off in the intensity of the signal. Figure 4.6(b) show the integral of the AO signal over the whole space for different source-detector distance on a log scale. This represents the total energy measured by the detector normalized by the total incident energy. Since one simulated photon can contribute to the signal in different tagging areas, the total normalized energy can be bigger than one. The measured energy decreases exponentially with the source-detector distance. This means that a compromise between imaging depth and signal to noise ratio will need to be found for experimental imaging.

Influence of incidence angle

Light scattering in biological tissues is strongly anisotropic and photons are mainly scattered in the forward direction (g = 0.9). In a reflection configuration, on can thus expect that the direction of photons reaching the detector is not isotropic. To investigate whether or not a preferred angle exist, the Monte-Carlo algorithm was modified in order to record the direction of the photons reaching the detector as well as their impact position on its surface. A discrete 2D matrix representing all possible combinations of polar and azimuthal angles was created with a step equal to 1 • and the energy reaching this pixel was recorded and normalized by the total incident energy. The same was done with a 100 × 100 pixel matrix representing the 8 × 8 mm 2 around the center of the detector. I then simulated the propagation of photons for different source-detector distances.

Figure 4.7 shows the distribution of the positions of photons reaching the detector. For a distance D between the source and the detector inferior to 3 cm, it appears that the photons mainly reach the boundary of the detector closest to the source whereas for larger distances, all positions seems as probable.

Figure 4.8 shows the angle distribution for three different distances between source and detector. Each vector of the image stands a given combination of polar and azimuthal angles (discrete representation with a step length of 1 • for each angle). Its length, as well as its colour, represents the amount of energy exiting the medium at in this specific direction. In this representation, the scattering medium is positioned in the positive z area and the source is located towards the positive-x. The fact that broad angles never occurs is due to total internal reflection. The Since one photon can contribute to the signal of several tagging areas, the total normalized energy can be higher than one. The energy exponentially decreases with the source-detector distance.

numerical aperture of the detection fibre was not taken into account here. It appears that when the source and detector are closer than 2 cm, the photons reaching the detector are more likely to arrive with a direction oriented towards negative x, opposite to the source. When the distance is greater, the distribution of angles tends towards homogeneity. These results seem to indicate that two regimes exit. One in which most photons take the shortest path between the source and the detector and arrive with an angle which retains some memory of the photons origin and one in which this information is lost and all angles and exit positions tend to be as likely. The limit between the two regimes is hard to evaluate but it seems to be of the order of 3 cm, which in this case corresponds to 30 times the transport mean free path. The second regime is likely to be the so-called diffusion regime in which light has totally lost all memory of its initial direction. This regime appears for propagation distances far greater than the transport mean free path which is consistent with our results.

This means that for distances greater than 3 cm between the source and the detector, the collection angle is not a crucial parameter as the light direction is almost isotropic. In the prospect of AOI in reflection mode, this is very interesting as it removes one degree of freedom for the optimisation of light collection. These numerical results are encouraging for acousto-optic imaging in reflection and give a preliminary idea of what is achievable in this configuration. between the source and the detector. Each vector originates from 0,0,0 and its length represents the total energy exiting the medium at this specific angle. Units are arbitrary and only serve to represent the angle distribution. In this representation, the scattering medium is located in the positive-z area and the source is located towards the positive-x. For D smaller or equal to 2 cm, some angles are more likely than others whereas for greater distances, it becomes almost isotropic.

Imaging in reflection mode

To confirm the results obtained with simulations, the PRC-based imaging setup was used to image phantoms in reflection with two different configurations.

The imaging setup

As mentioned previously, our setup consists of a photorefractive based acousto-optic imaging system coupled to a commercial ultrasound scanner (Aixplorer, SuperSonic Imagine, Aix-en-Provence, France). The scanner drives a linear probe -6 MHz central frequency, 192 elements, 0.2 mm pitch (SL10-2, SuperSonic Imagine, Aix-en-Provence, France) -for multi-modal imaging purposes. The same scanner and probe can be used to produce the acousto-optic image and a conventional ultrasound image, directly provided by the commercial mode of the Aixplorer, called B-mode. A 1D acousto-optic line is obtained by following the propagation of a focused ultrasound pulse, as described in chapter 2. The focusing is done by applying a delay line to the elements of the probe as described previously. To obtain an image, several consecutive 1D lines are recorded by electronically scanning the focus in front of the different elements, just as it is the case in a conventional B-mode imaging.

In the elevation direction, perpendicularly to the array (x-axis), the size of the focal spot is determined by the width of the transducers. This size is typically 200 µm at a distance 20 mm of the probe, thus defining the x-axis resolution (200 µm). In the lateral direction (z-axis), the aperture can be chosen by using a given number of elements. This number is chosen to depend on the focusing depth so that the lateral dimension of the focal spot is also 200 µm. The resolution in the longitudinal direction (y-axis) is given by the length of the pulse. We use 4-cycle long pulses which give a 0.67µs-long pulse at 6 MHz, corresponding to a spatial length of 1 mm. The setup is shown on figure 4.9. It consists of a semiconductor laser diode emitting at 790 nm, amplified by a 2 W tapered amplifier (Sacher laser technik, GmbH), forming a Master Oscillator Power Amplifier (MOPA). The amplified light is the split into two beams, the reference beam and the signal beam with a Polarizing Beam Splitter (PBS) and a half wave plate (λ/2) to control the polarization input and thus relative power of each arm of the interferometer. The reference beam is sent directly to the photorefractive crystal (PRC) through a multimode optical fibre (MMOF) whereas the signal beam is sent to the sample, through another MMOF. Scattered light is collected on the sample via a liquid core optical fibre (LCOF, diameter 4 mm, NA 0.4) and sent to interfere with the reference beam on the PRC to build the hologram and perform self-adaptive holography. The time varying signal coming out of the crystal is measured on a Si-photodiode (PD, Thorlabs PDA 36A-EC) with a detector surface of 10 mm 2 , amplified by a high bandwidth pre-amplifier (DHPVA-100, FEMTO Messtechnik -GmbH), filtered by a band-pass filter (200 Hz -5 MHz) (Multi channel filter -Model 3945 Khrone-Hite Gmbh) and acquired by an acquisition card (Express CompuScope 4344 -50 MS/s, 14-Bits, GaGe). The acquisition card is synchronized with the ultrasonic pulse via the ultrasound scanner.

Images are obtained by scanning the ultrasound focus inside the sample and stacking the acousto-optic lines one after the other. The lines are separated one from another by the probe pitch but since the lateral size of the focal spot is equal to the probe pitch, the whole medium is imaged. Each line is then repeated 1000 times for averaging purposes. Taking into account the initialisation time, the propagation, the cool-down time of the transducers and the transfer time from the buffer to the computer, an acousto-optic image is obtained in 50 to 80 seconds depending on the image size. For example, a 38 × 60 mm image is typically obtained in 70 seconds.

In order to image in a reflection configuration, we needed to be able to easily move the source and detector. For that, three optical fibres were added to the set-up as shown on figure 4.9. On the reference arm, a multi-mode fibre with a 200 µm core and a small collimator was used illuminate the crystal. To illuminate the sample, another multi-mode fibre was used and the incident beam had a 1 cm diameter. To collect the light exiting the sample, a liquid core fibre with an inner diameter of 0.8 cm and a 0.4 numerical aperture was chosen.

The first experiments were done on an homogeneous gel matrix mimicking tissues properties. This matrix is an agar gel mixed with Intralipid-10%, whose optical properties have been studied in [START_REF] Hugo | Light scattering in lntralipid-10% in the wavelength range of 400-1100 nm[END_REF]. To obtain the chosen reduced scattering coefficient (µ ′ s = 10cm -1 in the following), the proportions of Intralipid-10% and agar were calculated using the equation derived by the authors and Mie scattering theory. The scattering of Agar was assumed to be weak compared to Intralipid-10%. Absorption of such gels is considered negligible (usually µ a = 0.05 cm -1 @790 nm).

Agar was dissolved in water at a temperature of 90 • C and the solution was cooled down to 60 • C before adding Intralipid-10%. The gel was then poured into a mould where it solidifies around 40 • C.

When absorbing inclusions are needed, half of the gel was poured into the mould and left to cool down with small plastic rods dipped in it. When the half-gel was solidified, the rods were taken out and a gel-Indian ink mix was poured into the holes left by the rods. The other half of the gel was then poured on top of it. The proportion of gel and Indian ink can be chosen to obtain a desired absorption coefficient.

In the following, the gels' dimensions are 7 × 7 × 4 cm -where the smallest direction is placed along the optical axis -and they are homogeneous. The protocol for gel fabrication will be the same for the remainder of the manuscript. As the distance increases, the AO signal decreases but not at the same rate as in MC simulations. When the distance increases the shape of the signal changes and seems to be more intense in front of the detector. Banding on the first images is due to parasitic light and is not significant.

Influence of source-detector distance

The illumination and detection fibres were placed on the same side of the aforementioned gel and images were recorded for different positions of the detection fibre while the source stayed at the same position. Resulting images are shown on figure 4.10. As expected, the signal decreases significantly with the source-detector distance but not at the same rate as in Monte-Carlo simulations. This might be due to a mismatch between the experimental and numerical parameters such as absorption or scattering coefficients. The penetration depth increases with the source-detector distance but the signal contains increasing noise. For a separation of 40 mm, some signal is still visible but the signal to noise ratio (SNR) is very low. This suggests that it would be difficult to detect something with such a separation. Another interesting point is the shape of the signal. As the source-detector distance increases, the "banana" gets less and less symmetrical and a wider signal can be observed directly in front of the source. This might be explained by a very high anisotropy factor in such gels which favours forward scattering. Since the incident beam is collimated, photons will have a tendency to go mostly forward before going back to the detector. This seems to indicate that the system might be more sensitive in the region in front of the source rather than between the source and detector. To have a better idea of the penetration depth, profiles of the images were extracted and drawn in figure 4.11(a). The profiles were taken on the median plane between source and detector and a moving average over 10 points was applied. The profiles were then normalized by their maximum for representation purposes. The increase in penetration depth can be observed and despite the moving average, the profiles corresponding to the higher source-detector distances are significantly noisier. As we did before, we plotted the maximum of these profiles versus penetration depth as shown on figure 4.11(b). Penetration depth increases linearly with the source-detector distance and a linear fit was drawn. Because of the low SNR the last too points were excluded for the regression. A factor 0.25 was found which seems coherent with the factor 0.35 in the case of the Monte-Carlo simulations. The difference could be explained by a mismatch between experiment and simulation or by the fact that the extraction of the maximum value of data containing noise can lack precision.

These results show that it is possible to image inside highly scattering media in a reflection configuration and that different depth of imaging can be achieved by changing the distance between the source and detector. This also shows that the SNR is strongly affected by this distance and that for a source-detector separation higher than 4 cm, it becomes drops significantly. The penetration depth seems to be smaller than expected and it appears difficult to image deeper than 2 cm.

Towards a handheld probe

The results of the previous section lead us to try another configuration where both source and detector were on the same side of the sample. Given the spatial distribution of light within the medium, we decided to place the ultrasound probe between source and detector as shown on figure 4.12. Given the width of the probe and the size of the optical fibres holders, we decided to tilt the optical source and detector towards the center of the probe as can be seen on figure 4.12(b). This allowed us to reduce the source-detector distance and was consistent with the observed fact that a large proportion of the signal is located directly in front of the source.

New phantoms containing absorbing inclusions were made to test the sensitivity of such configuration, with the same scattering coefficient as before (10 cm -1 ). The first gel contained a cylindrical inclusion (diameter: 3 mm, oriented along the x-axis) positioned 5 mm deep inside the gel. It was imaged in the handheld configuration and the resulting image is shown on figure 4.13(a). The inclusion is clearly visible as it is positioned in the area of highest sensitivity. The Another phantom was made with two inclusions (diameters: 3 mm and 2 mm oriented along the x-axis) positioned 1 cm deep inside the gel. The image obtained in the handheld configuration is shown on figure 4.14. The two inclusions are visible but since they are out of the high signal areas, it is more difficult to see them. The two profiles show a much smaller CNR than before but the inclusions can still be detected. To increase the imaging depth the source-detector distance was increased but it only resulted in a decrease of the SNR and did not improved the imaging depth. Other trials were done on phantoms with inclusions located 2 cm from the surface but the imaging depth was insufficient to image these inclusions.

A multiple detector approach

We have seen in the previous sections that the reflection configuration is severely limited by the shape of the illuminated region inside the sample. The maximum of signal is obtained for relatively close source and detector but in this case the illuminated region is limited. On the other hand, when the distance between the two is too big, the SNR drops rapidly. In order to have more control over the illuminated region, we decided to try a new approach with multiple sources or detectors.

Using a multiplexing liquid core optical fibre with four inputs (Lumatech Gmbh) combined in one output, we slightly modified the setup by placing the four entries of the fibre around the light source as shown on figure 4.15(a). The photons collected by these four entries are combined in a larger fibre and sent to the PRC as before. The four detectors were placed on the corners of a square centred around the optical source in order to increase the size of the illuminated area within the sample.

Images were made on a phantom slightly less scattering than before (µ ′ s = 7 cm -1 )containing an highly absorbing cylindrical inclusion (diameter 3 mm, length 20 mm) orientated along the z-axis. Images were made in the xz-plane, orthogonally to the optical axis at different depths. The resulting images are shown on figure 4.15(b) for depths 7, 12, 17, and 22 mm. For small depths, one can clearly see the position of the source in the center and the four collection fibres in the corners forming an 'X'. As we go deeper into the sample, the shape of the 'X' blurs and the illuminated area becomes more even but the level of noises increases as well. The inclusion can clearly be seen in all images, including the image performed 2 cm deep within the sample even though the contrast is poor at this depth. In all images, the illuminated area is consequently bigger than it was with only one source and one detector and this could be useful for in-vivo imaging as in increases the field of view and could make for a clearer image. On the other hand, with multiple collecting fibres, the collected signal can vary from one to the other as can be seen on the images where the detector at the top left corner seems to collect more light. This could lead to more difficult interpretations since a higher signal could be due to structural differences in the medium or a better collection from a specific fibre.

This work was preliminary and a lot is left to investigate regarding the positioning of the detectors, the collection angles and the depth of imaging but it could prove a good asset for in-vivo imaging in back scattering mode.

Towards In Vivo imaging

As stated in chapter 3, the response time of our SPS photorefractive crystals is of the order of several milliseconds. This time is very close to the decorrelation time in biological tissues (typically 1 ms), which represents the characteristic time at which a speckle decorrelates because of motion (blood flow, breathing...). Unfortunately, because of the electronics and the need to let the piezoelectric transducer cool down between emissions, the imaging time of our setup is higher than this time. This means that our system is in principle not fast enough to follow the speckle decorrelation in vivo.

Despite this observation, we decided to try to obtain acousto-optic images in vivo with the setup in reflection mode. For this, the handheld probe configuration was used, as shown in figure 4.12(a) on a volunteer's wrist. The power of the ultrasound and the laser were set within the safety limits and the wrist was tied to the optical breadboard to limit the movements during acquisition. However, since the laser power had to be reduced, it was necessary to accumulate more for averaging purposes. The resulting sequence lasted several minutes, a time over which it is very difficult to stay still. Some images were acquired, an example of which is shown on figure 4.16. The acousto-optic signal was superimposed to the B-mode acoustic image . Biological information cannot be extracted and it was not the aim of this experiment but it is interesting to see that it is possible to engrave a blurry hologram and obtain images in vivo despite the decorrelation.

A high level of signal can be observed in the first millimetres but the depth of penetration was extremely small. The reason for this limited penetration is unknown but the motion of the wrist is likely to have had an impact. A lot of work remains to be done to obtain a definite proof of principle of in vivo acousto-optic imaging but this is very encouraging.

Conclusion

Imaging in reflection mode is a requisite step to apply AOI to medical applications in vivo. Because scattering is highly anisotropic in biological tissues, this configurations is less efficient than the transmission mode as fewer photons reach the detector. Using a custom Monte-Carlo algorithm, I simulated this imaging configuration and it appears that it is possible to probe a medium at a depth equal to 0.35 times the distance between the source and the detector. I also shown that when this distance is great compared to the transport mean free path, the angular distribution of the photons exiting the medium is close to isotropic. I then used our experimental set-up with optical fibres to show similar results experimentally. Showing that, despite a penetration depth slightly shorter than the simulated one, it is possible to image up to a few centimetres in highly scattering media. I investigated the use of several detection fibres at the same time and preliminary results indicate that this could be a good solution to increase the field of view of AOI. Our photorefractive crystal detection is still too slow to perform in vivo imaging but a short preliminary trial showed that despite the speckle decorrelation in tissues being faster than our detection, it might be possible to obtain optical information in living tissues. More work remains to be done to have a good proof of concept but this is a promising perspective for acousto-optic

Chapter's Keypoints

• Acousto-optic imaging can be performed with a commercial ultrasound scanner. With a transducer array, ultrasound can be focused by adjusting the delays between the transducer's elements.

• The transverse imaging resolution depends on characteristics of the transducer array. In this manuscript the transverse resolution is 0.5 mm.

• With pulsed ultrasound, the longitudinal resolution is set by the length of the pulses. In this manuscript it ranges between 0.5 and 2 mm.

• Acousto-optic imaging can be performed in a transmission or a reflection configuration.

• The reflection configuration was studied using Monte-Carlo simulations and experiments.

• The imaging depth in this configuration depends on the distance between the source and the detector, d. The maximum intensity of the AO signal was measured at 0.36d with simulations and 0.24d experimentally.

• In this configuration, the imaging depth can go up to 2 cm. For higher depths, the signal to noise ratio decreases significantly.

• Images were obtained in a 'hand-held probe' configuration where the transducer array is located between the light source and detector. Inclusions were detected at depths as high as 1 cm.

• Some images were obtained in vivo in a volunteer's wrist. Imaging depth was limited to a few millimetres.

CHAPTER 5

Inverse Problems for Quantitative Acousto Optic Imaging

Table of contents

We have seen in the previous chapters that acousto-optic imaging can be used to probe a scattering medium at centimetre depth with an optical contrast and an ultrasound resolution. This technique provides optical information about the medium as it is sensitive to the local light fluence. However, the shape of the illuminated area can sometimes affect the interpretation of the image since a variation of the signal can be due to a variation of the medium properties or a variation of illumination. Some imaging techniques, like X-ray tomography or diffuse optical tomography (DOT) use a numerical reconstruction to obtain the parameters of interest from the measured data. These reconstructions often require the solution of an inverse problem which links the space of the parameter to the space of the measurement. Whilst the use of such methods is not needed for AOI, a reconstruction algorithm could produce a map of absorption coefficient from AO signal, thus giving a clearer view of the medium properties. Moreover, this could also allow for quantitative reconstruction of the optical parameters inside a highly scattering medium which could be of significant value in a clinical context

In this chapter, we will first give a general description of an inverse problem and give some theoretical information on how inverse problems can be solved for the purpose of medical imaging. We will then focus on a reconstruction algorithm we developed which uses a simple forward model and experimental AO data to reconstruct the absorption coefficient in a scattering medium. After describing the method, I will discuss the results and the potential improvements of the technique.

Inverse problems for medical imaging

Inverse problems have been widely studied and they represent an important part of the field of numerical optimisation. The aim of this section is not to describe inverse problems in their broad generality but to give enough theoretical information to understand their use for acousto-optic imaging. We will thus focus on some specific points in order to introduce the necessary tools to understand the following section. More information about inverse problems and their application to imaging can be found in several books such as [START_REF] Bertero | Introduction to inverse problems in imaging[END_REF][START_REF] Kaipio | Statistical and computational inverse problems[END_REF][START_REF] Curtis | Computational methods for inverse problems[END_REF].

Definition of an inverse problem

The definition of an inverse problem resides in the correspondence between an object of interest, which we will call a parameter, and measured information about this parameter, which will be called the data or the measurements. This correspondence is described by a true forward model which links the parameter space to the space of the data. If F denotes the true forward model, x represents the parameter and y the data, we can write y = F(x).

(

Solving an inverse problem comes down to finding the parameter which correspond to the measured data. Since the true forward model is unknown, one needs an estimate of it in order to be able to simulate data from a given parameter. This estimate, which will be called the forward model, will be denoted by F in the following.

If the forward model is an invertible operator, then the best approach to solving the problem is to compute the operator inverse of the forward model, F -1 , and apply it to the data to find the best estimate of the parameter, x * , x * = F -1 (y).

(

In practice, forward models are often non invertible and the data can be corrupted by noise. In these cases, the solution of an inverse problem becomes a model fitting problem in which the aim is to find the parameter which minimizes the difference between measurements and a simulation of them. This difference is called the residual and the norm of the residual is called the objective function or the error functional and will be denoted E(x). The minimization problem is written:

x * := arg min x E(x) = arg min x ||y -F (x)||, (5.3)
where ||.|| represents a chosen norm.

Ill posed problems and regularization

If the inverse problem is well-posed, meaning that a unique solution exists and that the behaviour of the solution changes continuously with the initial condition, then the inverse problem y = F (x) can be solved using equation 5.3. However, inverse problems are often somewhat ill-posed. This means that a small change in the data, y, can result in a much larger change in the estimated parameter x * . Consequently, if the data is corrupted by noise or if there are some small numerical errors, the solution might be strongly affected.

To prevent this from happening, the inverse problem needs to be reformulated by including some hypotheses or prior knowledge. This process is known as regularization and can take many forms depending on the specifics of the problem. Equation 5.3 then becomes

x λ := arg min x E λ (x) = arg min x ||y -F (x)|| + λR(x), (5.4)
where R is an operator called the regularization operator and λ is a scalar know as the regularization hyper parameter which serves to control the relative contribution of the data term and the regularisation term to the objective function.

A commonly used method is called Tikhonov regularization. This consists of using the square norm of a spatial derivative of the parameter as a regularization operator. The case R(x) = ||x|| 2 is known as the zero-order Tikhonov while the first-order Tikhonov uses the gradient of the parameter (first-order derivative). The regularization method is chosen with respect to the prior knowledge on the parameter. For example, zero-order Tikhonov imposes small solutions while first-order Tikhonov enforces smooth solutions.

Linear problems and linearisation

In general, the forward model F is a non-linear operator, but in some specific cases, this operator can be linear and the inverse problem itself is called linear if the regularization is also linear. For example, the inverse problem defined by y = Ax, (

with a matrix A, and vectors x and y, is a linear inverse problem. If the matrix A is invertible, then the solution is given by x * = A -1 y. If it is not, some approximations can be made, such as using the pseudo-inverse of A. In these cases, the regularization is often applied by the means of a positive definite operator Y such that the operator (A * A + λY ) is invertible with a bounded inverse, where A * is the adjoint operator to A. The solution can then be written

x λ = (A * A + λY ) -1 A * y. ( 5.6)
When the inverse problem is not linear, it is sometimes possible to linearise it under some hypotheses in order to solve it with a linear method. Let us consider a noiseless inverse problem with a non linear forward model: y = F (x). Under the assumption of a small perturbation around a baseline represented by x 0 for the parameter and y 0 for the data, it is possible to expand the forward model in a Taylor series around x 0 to the first order

y = y 0 + F ′ (x 0 )(x -x 0 ) + O(x 2 ).
(5.7)

Consequently, (x

-x 0 ) ≈ (F ′ (x 0 )) -1 (y -y 0 ), (5.8)
and the problem can be solved approximately by calculating the first derivative of the forward model. In a discrete case, this corresponds to the Jacobian matrix [START_REF] Simon | Optical tomography in medical imaging[END_REF]. Linear inversion is the simplest and fastest way of solving inverse problems but it is of limited use in most practical cases because of the assumptions needed to linearise. However, when considering perturbations around a baseline it can be useful to have a preliminary approximate idea of the solution before using other methods.

Reconstruction algorithms

A great number of algorithms exist to solve non-linear inverse problems and the choice of the method is dependent on the particular problem one needs to solve. Most of them are iterative processes in which the estimate of the parameter is updated over iterations in order to refine the solution. We will here describe three of the most commonly used algorithms: the Damped Gauss-Newton algorithm, the Gradient Descent and the Conjugate Gradients method [START_REF] Bertero | Introduction to inverse problems in imaging[END_REF][START_REF] Lloyd N Trefethen | Numerical linear algebra[END_REF].

Damped Gauss-Newton

Let us consider a discrete non-linear inverse problem y = F (x) where x = [x 1 , ..., x n ] and y = [y 1 , ..., y m ] with n the dimension of the parameter and m the number of measurements. To find the optimal parameter, we will try to minimize the norm of the residual, r(x) = [r 1 (x), ..., r m (x)], here given by r(x) = F (x) -y. We can then express an objective function as follows

E(x) = 1 2 ||r(x)|| 2 = 1 2 ||F (x) -y|| 2 .
(5.9)

The Gauss Newton algorithm is an iterative algorithm based on the Newton method which can be used to locate roots or extrema of a function. Newton's method is derived from the Taylor expansion of a function f (x) around a point x (0)

f (x (0) + δ) = f (x (0) ) + δ.f ′ (x (0) ) + δ 2 2! .f ′′ (x (0) ) + O(δ 3 )
If we apply this equality to the first derivative of f we have, to the first order

f ′ (x (0) + δ) ≈ f ′ (x (0) ) + δ.f ′′ (x (0) ), (5.10)
which is the equation of the tangent. If we call x (1) = x (0) + δ (0) the point at which f is minimal, we can easily derive

x (1) = x (0) + δ (0) = x (0) - f ′ (x (0) ) f ′′ (x (0) ) . ( 5.11)
This x (1) is a first-order guess of the root of f ′ , to get a better estimation of this root, Newton's method uses the iterative process:

x (i+1) = x (i) - f ′ (x (i) ) f ′′ (x (i) ) . ( 5.12)
If we now consider our function E(x) of multiple variables {r i }, the first derivative is the gradient and the second derivative is the Hessian. To calculate the gradient and Hessian of E we will need to calculate the gradient of the residual vector r(x) = [r 1 (x), r 2 (x), ...r m (x)]. Since the measured data is a fixed quantity which does not change with the parameter we are trying to reconstruct, the gradient of the residual only involves the forward model. The Jacobian, J, is a matrix which represents the gradient of a vector function. The Jacobian of the residual thus only involves the forward model:

J(x) = ∂r j
∂x i j=1,...,m;i=1,...,n = ∂F (x) j ∂x i j=1,...,m;i=1,...,n .

(

The gradient and Hessian of our objective function can be expressed in terms of the Jacobian:

∇E(x) = m j=1 r j (x)∇r j (x) = J(x) T r(x), (5.14)
∇ 2 E(x) = m j=1 ∇r j (x) T ∇r j (x) + m j=1 r j (x)∇ 2 r j (x) = J(x) T J(x) + m j=1 r j (x)∇ 2 r j (x).
(5.15)

In cases where the residual is close to the solution, the second term of the Hessian can be neglected which eliminates a lengthy calculation of ∇ 2 r j (x). Depending on the problem, it is sometimes easier to calculate the Jacobian rather than calculating the gradient and Hessian independently.

Damped Gauss Newton

The Damped Gauss Newton algorithm is based on the basic equation (5.12) except that it uses a step length that can be different from one.

x (i+1) = x (i) + s (i) d (i) = x (i) -s (i) (∇ 2 E(x i )) -1 ∇E(x i ).
(5.16)

The vector d (i) is called the search direction and s (i) is called the step length. The convergence of the Gauss Newton algorithm can be improved by looking for the optimal step length rather than always using a step of one.

In order to find this optimal step length, one can use a line search method which will minimize the objective function ||F (x (i) + s (i) d (i) ) -y|| with respect to the variable s (i) . More information about line search methods can be found in chapter three of [START_REF] Stephen | Numerical optimization[END_REF].

The Gauss-Newton algorithm has the advantage that it does not require the calculation of second order derivative that can be challenging to compute. However, depending on the initial parameter, convergence is not guaranteed [START_REF] Curtis | Computational methods for inverse problems[END_REF] and it requires the storage of the Jacobian, J T J, and the gradient which can be an issue for large matrices.

Gradient descent

To avoid this potentially huge memory usage, other methods exist amongst which is a family of algorithms based on the gradient of the objective function. Gradient descent and conjugate gradients are two methods belonging to this family.

Gradient descent, or the steepest descent algorithm, is the simplest of the gradient based methods. This is an iterative method which relies on the fact that at any given point, the gradient of a function points away from the local minimum. Consequently, by updating the current estimate of the parameter in the direction opposite to the gradient, the estimation gets closer to the minimum.

Let us consider the simple linear problem Ax = y, (5.17) where A and y are known.

One can show with some calculation, that if A is positive definite, the following function f (x), called the quadratic form of the equation, is minimal when Ax = y:

f (x) = 1 2
x T Ax -y T , (5.18) where the superscript T denotes the transpose. Thus it can be used as the objective function in this specific case. Two representations of this quadratic form are shown in figure 5.1 for given values of A and y. In this example, the optimal solution is x = [-1, 2] T For a given point x (i) , the steepest direction towards the minimum is given by the opposite of the gradient calculated at this point. In the case of the linear problem, the gradient is given by ∇f (x) = Ax (i) -y, the direction of update will then be d (i) = b -Ax (i) . The optimal step length, s (i) , can then be found by doing a line search to minimize f (x (i) + s (i) d (i) ) with respect to s (i) . When A is known, this step length can be explicitly derived [START_REF] Richard | An introduction to the conjugate gradient method without the agonizing pain[END_REF]

s (i) = d T (i) d (i) d T (i) Ad (i)
.

(5.19)

When the matrix A is unknown and for non linear problems, the line search is performed by numerically minimizing ||F (x (i) + s (i) d (i) ) -y||. The algorithm proceeds as follows

Step 1 Calculate the objective function

E n (x) = 1 2 ||F (x n ) -y||
Step 2 Calculate the gradient of the objective function ∇E n (x)

Step 3 Find the best step length, s n by doing a line search

Step 4 Take a step in the direction of the negative gradient

x n+1 = x n -s n ∇E n (x)
Step 5 Go back to step 1 Figure 5.2 shows the evolution of the estimate of the parameter using a Gradient Descent algorithm for an initial guess

x (0) = [-2, -2] T .

Conjugate Gradients

Gradient descent will always converge to the solution but it has a tendency to take many steps to do so because their direction is not exactly toward the solution. This makes gradient descent slow. To avoid this behaviour, the Conjugate Gradient (CG) method imposes the directions of each step to be orthogonal to one another. This assures that the algorithm does not go twice in the same direction but rather takes an optimal step for each of the independent dimensions.

Theoretical considerations If we write {x 1 , x 2 , ..., x n } the coordinates of vector x, then the first step will find an optimal value for x 1 , the next one the value for x 2 and so on until we have reconstructed x in n steps, where n is the dimension of x. Suppose that we want to solve the linear system:

Ax = y, ( 5.20)
where A and y are known. Then we can define the objective function that we want to minimize as:

E(x) = 1 2 ||y -Ax|| 2 2 = 1 2 ||r(x)|| 2 2 , ( 5.21)
where we have chosen to minimize a Euclidean norm in the sense of the least squares. Let us call

x (i) our estimation of x at step i, e (i) = x (i) -x the error which quantifies the distance from the exact solution and D n = {d (1) , d (2) , ..., d (n) } the family of n orthogonal directions we are going to use to find x. At each step we choose a point x (i+1) as follows:

x (i+1) = x (i) + s (i) d (i) . ( 5.22)
To find s (i) , we will use the fact that we never want to go in the direction d(i) after this step. Which means that e (i+1) should be orthogonal to d (i) . An easy way to understand this is to consider e (i+1) as the distance between our position, x (i+1) and our destination x. If e (i+1) is not orthogonal to d (i) , it means that to go from x (i+1) to x, we will need to go in the direction d (i) again at some point to reach x. We thus have:

d T (i) e (i+1) = 0 d T (i) (e (i) + s (i) d (i) ) = 0 (by equation 5.22)
s (i) = - d T (i) e (i) d T (i) d (i)
.

(5.23) However, this is not helpful because we need e (i) to find s (i) and e (i) was known, the problem would have been solved already. The solution used in the Conjugate Gradient method is to use search directions which are A-orthogonal instead of orthogonal. Two vectors u and v are A-orthogonal or conjugate if:

u T Av = 0.
(5.24)

We now want e (i+1) to be A-orthogonal to d (i)

d T (i) Ae (i+1) = 0 d T (i) A(e (i) + s (i) d (i) ) = 0 s (i) = - d T (i) Ae (i) d T (i) Ad (i) s (i) = d T (i) r (i) d T (i) Ad (i) , ( 5.25)
since

r (i) = y -Ax (i) = -Ae (i) .
Unlike equation (5.23), this s (i) can be calculated. All we need now is to find a set of Aorthogonal directions {d (i) }. This can be done by applying the conjugate Gram-Schmidt process to a set of n linearly independent vectors {u 0 , u 1 , ..., u n }. To construct the vector d (i) , we take the vector u (i) and we subtract everything which is not A-orthogonal to the previous search directions:

d (i) = u i - i-i k=0 β ik d (k) , (5.26)
where the β ik are defined for i > k. To find them, we use the same method we used to find s

d T (i) Ad (j) = u T i Ad (j) - i-i k=0 β ik d T (k) Ad (j) , 0 = u T i Ad (j) -β ij d T (j) Ad (j)
, i > j (by A-orthogonality of d vectors),

β (ij) = u T i Ad (j) d T (j) Ad (j)
.

(5.27)

The difficulty with this technique is that the old search directions have to be kept in memory to calculate the next one, which is computationally expensive. In the Conjugate Gradient method the search directions are constructed by conjugation of the residuals instead of an arbitrary set of linearly independent vectors. (That is, by setting u i = r (i) ). This choice has several implications. First, since r (i) = -Ae (i) and e (i) is A-orthogonal to the previous search direction we see that each residual is orthogonal to the previous search directions.

r (i+1) = -Ae (i+1) , d T (i) r (i+1) = -d T (i) Ae (i+1) = 0.
Thus, if we take the inner product of equation 5.26:

d T (i) r (j) = u T i r (j) + i-i k=0 β ik d T (k) r(j), 0 = u T i r (j) , i < j, 0 = r T (i) r (j) , i = j (since u i = r (i) ), (5.28)
we see that each residual is also orthogonal to all the previous residuals. We can also note that:

d T (i) r (i) = r T i r (i) , (5.29)
which will be useful later.

Then, by deriving:

r (i+1) = -Ae (i+1) , = -A(e (i) + s (i) d (i) ), = r (i) -s (i) Ad (i) , (5.30)
we see that each residual is a linear combination of the previous residuals and Ad (i) . By taking the inner product of r i+1 and equation 5.30 we can write:

r T (i+1) r (j+1) = r T (i+1) r (j) -s (j) r T (i+1) Ad (j) ,
(5.31) 0 = s (j) r T (i+1) Ad (j) , j < i (by equation 5.28). (5.32) This means that Gram-Schmidt will be easy since r (i+1) is already A-orthogonal to all the previous search direction except d i . Let us recall that the Gram-Schmidt constants are now:

β (ij) = r T i Ad (j) d T (j) Ad (j)
. If we now rewrite equation 5.31 in a more general case, we can write:

s (j) r T (i) Ad (j) = r T (i+1) r (j) -r T (i+1) r (j+1) ; r T (i) Ad (j) =            1 s (i) r T (i) r (i) if j = 1, -1 s (i-1) r T (i) r (i) if j = i -1, 0 otherwise; β ij =      - 1
s (i-1) r T (i) r (i) d T (i-1) Ad (i-1) if j = i -1, 0 if j < i -1.
(5.33)

We now see that we only need to store the previous direction to find the new one, which is a significant advantage of Conjugate Gradient over other methods. We can thus simplify β ij by using the abbreviation β (i) = β i,i-1 . We have:

β i = - r T (i) r (i) d T (i-1) r (i-1)
(by equation 5.25),

= - r T (i) r (i) r T (i-1) r (i-1)
(by equation 5.29).

(5.34)

The Conjugate Gradient method can thus be written:

d (0) = r (0) = y -Ax(0), (5.35)
r (i) = r (i-1) -s (i-1) Ad (i-1) , (5.36)

β (i) = - r T (i) r (i) r T (i-1) r (i-1)
, (5.37)

d (i) = r (i) -β (i) d (i-1) , ( 5.38
)

s (i) = r T (i) r (i) d T (i) Ad (i)
(by equations 5.25 and 5.29), (5.39)

x (i+1) = x (i) + s (i) d (i) .
(5.40)

Conjugate Gradients

x 1 -4 -3 -2 -1 0 1 2 3 4 5 6 x 2 -6 -5 -4 -3 -2 -

Non linear Algorithm

In a more general case, we are trying to solve a non linear problem which cannot be expressed as Ax = y. Consequently we do not have access to A and we will have to use a line search at each step to find the optimal value of s (i) . For the same reason, the calculation of r (i) as described by equation 5.36 will not be possible. We will thus use the fact that the residual is indeed equal to the opposite of the gradient of the objective function:

r (i) = -∇E(x (i)
). We can then use the following algorithm:

Step 1 Calculate the initial steepest direction ∆x (0) = -∇E(x (0) )

Step 2 Find s (0) by doing a line search

Step 3 Take a step in steepest direction

x (1) = x 0 + s (0) ∆x (0)
Step 4 Calculate the new steepest direction ∆x (i) = -∇E(x (i) )

Step 5 Calculate β (i) with the Polak-Ribiere formula

β (i) = ∆x T (i) (∆x (i) -∆x (i-1) ) ∆x T (i-1) ∆x (i-1)
Step 6 Update the search direction

d (i) = ∆x (i) + β (i) d (i-1)
Step 7 Find s (i) by doing a line search

Step 8 Take a step in the direction d (i)

x (i+1) = x (i) + s (i) d (i)
Step 9 Go back to step 4

Convergence In ideal conditions, the Conjugate Gradient algorithm will converge in n steps, where n is the dimension of the parameter. However, as soon as the measurement contains noise, round-off errors will cause the residual to gradually lose accuracy and the search vector to lose Aorthogonality over time. Consequently, CG can be very sensitive to noise in the measurements. For this reasons, CG was discarded in the 1960s by the mathematical society before it was considered again in the 1970s as an iterative method. Since then, CG has be extensively studied and it can been shown that for two dimensional problems, CG has a time complexity of O(n 3/2 ) versus O(n 2 ) for Gradient Descent; for three dimensional problems CG has a time complexity of O(n 4/2 ) versus O(n 5/3 ) for Gradient Descent. For this reason CG is now one of the most commonly used methods.

Application to imaging

Inverse problems occur in a wide variety of fields, from geophysics to oceanography or medical imaging. Most of the time, an inverse problem needs to be solved when measurements are performed on the boundary of a medium with no access to the inside. This is notably the case of tomography techniques such as X-ray computerized tomography (CT), positron emission tomography (PET) or diffuse optical tomography (DOT) [START_REF] Simon | Optical tomography in medical imaging[END_REF][START_REF] John | Positron-emission tomography[END_REF][START_REF] Simon R Arridge | Optical tomography: forward and inverse problems[END_REF].

CT and PET scans

In CT, one measures the attenuation of a beam of X-rays travelling through a medium. X-rays propagate in a straight line and the obtained information is thus integrated over a line. The medium is then rotated over 360 • to obtain a full scan of a 2D plane. In PET, positrons are emitted by a radioactive tracer isotope [START_REF] John | Positron-emission tomography[END_REF]. After a short travel, positrons recombine with electrons and it creates a pair of gamma rays propagating in straight lines in opposite directions that are measured by scintillators. Some PET-scanners have a short enough resolving time (under 500 picoseconds) to localize the recombination position by time-of-flight measurement but it is not the case of most PET-scanners which only obtain information about a line in the medium.

To obtain an image from these measurements thus requires the solution of an inverse problem. Mathematically, each measurement represents a line integral of the parameter along the path of the beam. This mathematical transform, called the Radon Transform, has been described in 1917 by Johann Radon [START_REF] Radon | Uber die bestimmug von funktionen durch ihre integralwerte laengs geweisser mannigfaltigkeiten[END_REF][START_REF] Radon | On the determination of functions from their integral values along certain manifolds[END_REF]. This is a linear transform which can be inverted as shown in the same paper. Despite the ill-posedness of the problem, this inverse transform, also called filtered back projection, is usually used to form the image when the amount of noise is reasonable. When the measurements contain too much noise, iterative methods can be used but the computational expense if often considered too big and corrected versions of the filtered back projection method are used.

Diffuse Optical Tomography

Diffuse Optical Tomography uses near infra-red light and several sources and detectors to reconstruct an optical image of scattering medium. Sources and detectors are positioned on the boundary of the medium and each measurement represents the amount of light going from a given source to a given detector. Because of multiple scattering, a single measurement will give information about a wide area of the medium without possibility to localize this information. Consequently, a lot of measurements are needed and an inverse problem needs to be solved to reconstruct an image. In this case, the inverse problem is non-linear and an iterative method such as Gauss-Newton or Conjugate Gradients needs to be used with an appropriate forward model.

Inverse problems for acousto-optic imaging

Acousto-optic imaging provides optical information about scattering samples through the mapping of the local light fluence deep inside the medium. This is thus a technique sensitive to both absorption and scattering and it is sometimes difficult to differentiate the influence of each phenomenon. A contrast in an AO image can then come from a change in either the absorption or scattering coefficient or both. Moreover, since the illumination is not always homogeneous in the medium, it could lead to misinterpretation of the AO image. One solution to refine the information provided by AOI is to solve an inverse problem for which the measurement is the AO data and the parameters are the absorption and scattering coefficient. More than separating the contribution of these two coefficients to the AO image, such a process could also provide quantitative information about the optical properties of the medium which would be of great value for medical diagnosis. Powell et al. showed that such an approach is viable for reconstruction of both absorption and scattering coefficient from simulated acousto-optic measurements [START_REF] Powell | Quantitative reconstruction of absorption and scattering coefficients in ultrasound-modulated optical tomography[END_REF]. We developed a reconstruction algorithm to apply this method to experimental data acquired with the photorefractive crystal-based AOI setup described in the previous chapters. After describing the specificities of our inverse problem and forward model, I will present and discuss the results we obtained on experimental data.

The inverse problem and model corrections

The approach we used is based on previous work by Powell et al. [START_REF] Powell | Quantitative reconstruction of absorption and scattering coefficients in ultrasound-modulated optical tomography[END_REF]. We adapted the technique to be able to use experimental data instead of simulated AO signals. For that we defined a new inverse problem by restricting the space of the parameter and we chose a forward model which we had to correct in order to take into account the differences between the experiment and the simulation.

The inverse problem

Reconstructing both absorption and scattering coefficients at the same time appears to be possible but in order to obtain a unique solution, several measurements are needed with different source and detector positions. When applying this technique to experimental data, we decided to start by reconstructing only one of these coefficients. We chose the absorption coefficient because the influence of absorption on the AO signal is simpler to understand and has been more studied than that of scattering.

Consequently, in the following, the scattering coefficient will be assumed constant over the whole medium and will be denoted µ ′ s while the absorption coefficient will be a function of space µ a (r). The same work could theoretically be done for scattering by inverting the role the two coefficients.

We will use the same formalism as in the previous section and call y our experimental data and F (x) = F (µ a ) our forward model. The inverse problem we want to solve can then be written:

E(µ a ) = ||y -F (µ a )|| 2 2 + λR(µ a ) (5.41) µ * a = arg min µa E(µ a ), (5.42)
where E(µ a ) is our objective function, R(µ a ) represents a regularization function, λ is the regularization parameter, µ * a represents our estimate of µ a and we have chosen to use the square of the Euclidean norm to calculate the distance between the data and the simulation.

The forward model

A wide range of forward models can be used to solve inverse problems related to light transport in a scattering medium. Monte Carlo (MC) schemes can be very accurate models but they are often computationally expensive, and even with a fast GPU calculation [START_REF] Powell | Highly parallel monte-carlo simulations of the acoustooptic effect in heterogeneous turbid media[END_REF], the computation requirement are such that they can only serve as the 'gold-standard' to validate more approximate techniques. Among these, methods relying on approximations to the Radiation Transfer Equation (RTE) are commonly used. For very simple problems, it is possible to solve the RTE analytically but for more realistic media, approximations are needed. A common approximation is the Diffusion Approximation (DA). This method can be computationally efficient and we decided to use such a method to solve the Acousto-Optic forward model.

Coupled diffusion equations

In radiative transport theory, the propagation of light through a medium is expressed in terms of conservation laws which account for losses and gains of photons due to absorption and scattering [START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF]. The quantity of interest is the specific intensity, I(r, ŝ) which represents the intensity in the domain Ω at the position r in the direction ŝ. The specific intensity obeys the RTE:

1 c ∂I ∂t + ŝ • ∇I + (µ a + µ s )I = µ s p(ŝ ′ , ŝ)I(r, ŝ′ , t)dŝ ′ + q(r, ŝ, t), r ∈ Ω (5.43)
where c is the light velocity in the medium and p(ŝ ′ , ŝ) is the phase function, which defines how light arriving at position r in the direction ŝ is scattered in the direction ŝ′ . The phase function is symmetric with respect to interchange of its arguments and obeys a normalization condition:

p(ŝ ′ , ŝ)dŝ ′ = 1.
(5.44)

One often assumes that the phase function only depends upon the angle between ŝ and ŝ′ , which holds for spherical particles. For isotropic scattering, we can write p = 1/(4π), in a more general case, we can define the anisotropy factor by:

g = (ŝ ′ • ŝ)p(ŝ ′ • ŝ)dŝ ′ .
(5.45)

By making appropriate assumptions about the propagation of photons inside a scattering medium, it is possible to simplify the RTE to reduce the number of independent variables. The diffusion approximation to the RTE is widely used to simulate light transport. This is valid in the regime where the scattering coefficient is large, the absorption coefficient is small, the point of observation is far from the boundary of the medium and the time scale is sufficiently long.

With these assumptions and a development of the specific intensity on spherical harmonics, one can obtain the diffusion equation for the energy density Φ(r) = 1 c I(r, ŝ)dŝ [START_REF] Simon R Arridge | Optical tomography: forward and inverse problems[END_REF]:

[-∇ • D(r)∇ + µ a (r)] Φ(r) = q(r), r ∈ Ω, (5.46)
where D = 1/(3(1 -g)µ s (r)) = 1/(3µ ′ s (r)) is the diffusion coefficient and q(r) describes the light source.

This approximated equation is commonly used to simulate light propagation in infinite domains. In order to apply it to finite size domains, one must add a boundary condition. Several approaches can be envisioned [START_REF] Schweiger | The finite element method for the propagation of light in scattering media: Boundary and source conditions[END_REF], we used a Robin (impedance) style boundary condition [START_REF] Schweiger | The finite element method for the propagation of light in scattering media: Boundary and source conditions[END_REF]:

Φ(r) + 2An • D(r)∇Φ(r) = 0, r ∈ ∂Ω, ( 5.47)
where A is a scalar linked to the refraction mismatch between the domain and the surrounding medium, n is the unit outward normal to the boundary and ∂Ω represents the boundary. This boundary condition represents the physical model of a non scattering medium surrounding Ω with an index mismatch and associated reflections and refractions.

To simulate the acousto-optic tagging, another diffusion equation was needed:

[-∇ • D(r)∇ + µ a (r)] Φ a (r) = q a (r), r ∈ Ω, (5.48)
where Φ a (r) represents the energy density of the tagged photons and q a (r) the acousto-optic "virtual source" given by q a (r) = η a (r)Φ 0 (r), (5.49) where η a (r) represent the acousto-optic tagging efficiency and Φ 0 (r) is the solution of equation 5.47 with a source term q 0 (r) representing an isotropic source of monochromatic coherent light. The same boundary condition was used with this second equation. To represent a finite detection aperture, the measurement was obtained by integrating the flux over the area of the detector.

These two equations, with their boundary conditions, constituted our forward model. After defining the domain and the position and type of the source, equation 5.47 is used to calculate the propagation of photons inside the medium. For each ultrasound position, a virtual source term is derived and equation 5.48 is used to propagate the tagged photons emanating from this source. The energy density is then integrated over the surface of the detector and the value is attributed to the pixel corresponding to the ultrasound location. on the left hand side. 5.4(b) shows the energy density of tagged photons for one specific position of the ultrasound focus and 5.4(c) shows the final acousto-optic image for a set of ultrasound locations and a detector on the right hand side of the domain.

Model Corrections

In the forward model, we take the tagging efficiency of the ultrasound as a known parameter. Yet in practice, it is extremely difficult to measure this coefficient experimentally. To measure it, one would need to calculate the ratio of tagged photons and untagged photons which went through the ultrasound focus. Since the localisation of light is only possible through the tagging of the ultrasound, it is almost not possible to directly measure the tagging efficiency. Moreover, several other unknown experimental parameters such as the coupling of the ultrasound or some variations of the laser mode can affect the experimental data. Consequently, a mismatch exist between the experimental image and the simulation produced by the forward model which could prevent the algorithm from converging.

To circumvent this problem, we decided to acquire a baseline measurement in addition to the image of interest. This measurement is performed in an homogeneous phantom with the same properties as the medium of interest. The objective function is then modified to take into account this baseline measurement, y 0 , and the forward model simulated for an homogeneous medium, .50) Another source of mismatch between experimental data and forward model is a scaling difference. The incident flux used for the forward model is chosen to match the experimental one. However, the forward model is chosen to be quite simple and because of the experimental light collection and other experimental uncertainties, the scaling between experimental and simulated images is often very different. To correct this difference, a scaling factor was inferred from the baseline measurements. A simple scalar rather than a pixel-wise scaling factor was chosen in order to keep the influence of these model corrections to a minimum. This scalar, α was obtained by calculating the ratio between the mean value of experimental and simulated baseline measurements α = y 0 /F (µ a0 ), (5.51) where the over-line denotes the mean value.

F (µ a0 ), E(µ a ) = ||y -[F (µ a ) -F (µ a0 ) + y 0 ] || 2 2 + λR(µ a ). ( 5
The objective function was then modified as follows to take into account this scaling factor

E(µ a ) = ||y -[α (F (µ a ) -F (µ a0 )) + y 0 ] || 2 2 + λR(µ a ).
(5.52)

The inverse problem was regularized using a zero-order or first-order Tikhonov regularization:

R 0 (µ a ) = ||µ a || 2 2 , (5.53) R 1 (µ a ) = ||∇µ a || 2 2 , ( 5.54)
to either impose small or smooth solutions.

Numerical implementation

The light propagation was numerically implemented using a Finite Element Method (FEM) included in a software suite for forward and inverse modelling called Toast++, developed by Schweiger and Arridge for optical tomography [START_REF] Schweiger | The toast++ software suite for forward and inverse modeling in optical tomography[END_REF]. This open-source suite contains a selection of tools and functions for modelling the propagation of light in scattering media.

These functions are used to create a triangular mesh representing the imaging domain. To reduce the computational time, and because the experimental data were acquired on a 2D plane of the medium, we decided to restrict the simulation of light propagation to a 2D plane. This plane typically contained 2×10 4 nodes linked by a triangular mesh. The source and detector positions are defined on the boundary of the mesh.

To represent the ultrasound tagging areas, we decided to use a square basis in which one ultrasound location correspond to one pixel of the basis. This does not represent exactly the experimental conditions in which the ultrasound intensity is not evenly distributed in a square but in the AO data, a pixel of the image corresponds to a specific ultrasound location. Using this basis representation thus combines the advantage that it represents the way the experimental data is mapped and that it makes the calculation of the objective function easier. Because the experimental data contains a great number of pixels, it is down-sampled in order to reduce the size of the ultrasound basis to make the computation faster. We use a Toast++ function to move back and forth between the mesh and the basis representations.

The calculation of the Jacobian and gradient is made according to the derivations proposed by S. Powell in his paper [START_REF] Powell | Gradient-based quantitative image reconstruction in ultrasound-modulated optical tomography: first harmonic measurement type in a linearised diffusion formulation[END_REF]. After calculating the energy density of the untagged and tagged photons, φ and φ a respectively, by solving equations 5.47 and 5.48 for each position of the ultrasound, their adjoint fields, φ * and φ * a are derived. The Jacobian and gradient are then constructed:

J(µ

a ) = -      φ T 1 φ * a,1 + φ T a,1 φ * 1 . . . φ T Na φ * a,Na + φ T a,Na φ * Na      (5.55) ∇E(µ a ) = - ρ φ T ρ φ * a,ρ + φ T a,ρ φ * ρ (5.56)
where the subscript index ρ represents a given position of the ultrasound and N a is the total number of acoustic locations.

At each iteration of the reconstruction algorithms (Gauss Newton, Gradient Descent or Conjugate Gradients) the line search used for the determination of the step lengths is performed using a function implemented in Toast++. The algorithm stops when the stopping criterion is reached: if the evolution of the objective function between two steps is less than 10 -3 times the initial value. Figure 5.5 summarize the functioning of our algorithm using the conjugate gradients method.

Reconstruction of the absorption coefficient

Experimental configuration

The phantoms were prepared using the protocol described in chapter 4 with several absorbing inclusions. The absorption coefficient of a mother India ink solution was measured with an absorption spectrometer and a dilution was performed in order to obtain the desired absorption coefficient. The diluted ink is then mixed with the gel to create the inclusion. Given the small size of the inclusion (typically a cylinder with a 3 mm diameter and a height of 3 mm), the gel-ink mixture cools down and solidify very fast. Controlling the mixing of the two is thus difficult. In order to assess our experimental error regarding the absorption of the inclusion, I measured the absorption of a set of 10 mixtures. From this, I estimated an uncertainty of 10 % on the value of µ a . In order to perform the baseline measurement required for the solution of the inverse problem, another identical phantom without inclusion was needed. Because of experimental variations, creating two phantoms with the same optical and mechanical properties would have been challenging. To ensure that the physical properties were identical for the both measurements, I decided to create wider phantoms and to position the absorbing inclusions on one side of the phantom so that the other side remained homogeneous. This allowed for both measurements to be performed on the same phantom by simply sliding it in order to image the area of interest as shown on figure 5.6. This strongly reduced the unwanted variations between the two images' acquisition.

The experimental data were acquired using our PRC-based AOI setup in Institut Langevin as described in the previous chapters in a transmission configuration. The plane of imaging was aligned on the optical axis as shown in figure 5.6 so that the source and detector were located on positions belonging to the boundary of the plane. This constraint is imposed by the choice of a 2D simulation of the forward model which imposes the source and detector to be on the boundary to the imaging plane. The length of the ultrasound pulse was chosen to have a 1 mm longitudinal resolution. The lateral resolution is fixed by the properties of the acoustic probe and is of the order of 0.5 mm. In order to reduce the amount of noise, the images were averaged 4000 times. inclusions with absorption coefficients of 0.8 and 0.5 cm -1 . The choice of the absorption coefficients was made to respect as much as possible the hypothesis µ a ≪ µ ′ s used for the diffusion approximation while keeping a high enough contrast to be able to detect the inclusions experimentally. Given the low absorption, the inclusions are barely visible on figure 5.7(a) but they can be observed on the difference data shown in figure 5.7(c) in which the baseline image was subtracted from the image. The spatial inhomogeneities in the experimental images are not fully understood but could be explained by a combination of several experimental factors such as slight variations in the ultrasound coupling, variation of the ultrasound pressure due to small defects of the acoustic transducer or inhomogeneous concentration of Intralipid-10% in the phantom due to the fabrication protocol. These spatial variations justify the use of a baseline measurement which removes these unknown spatial parameters. Figure 5.7(d),(e) and (f) show the equivalent images as simulated for the same parameters using our forward model. The contrast appears to be much greater in the simulated image than in the experimental one and that the average signal is significantly reduced between the baseline and the measurement, likely because of the light absorption by the inclusions. This results in a difference data in which the light shape is more visible than in the experimental equivalent.

Choice of the regularization parameter

The choice of the correct regularization parameter can be crucial for ill-posed inverse problem solution. This parameter quantifies the importance of the regularization term with respect to the norm of the residual. If it is too small, the noise in the data can be amplified and the algorithm may converge towards a incorrect solution or not converge at all. If it is too high, the solution can be over-smoothed and some information present in the data can be lost. To choose the optimal parameter, several techniques exist such as the L-curve [START_REF] Per | The use of the l-curve in the regularization of discrete ill-posed problems[END_REF]. This method consists in calculating the value of the regularization term and the residual norm (the two terms of the objective function) for different values of λ and to plot one as a function of the other. The best value of the regularization parameter is the one that minimizes both quantities and for typical ill-posed problem, the curve exhibits a 'L' shape, hence the name of the method. For ill-posed problems, the curve usually exhibits a 'L' shape and the optimal value for the regularization corresponds to the corner of the 'L'. Here, no 'L' can be observed.

Figure 5.8 shows an example of such a function for our problem in the case of zero order Tikhonov regularization. In our case, it can be seen that the characteristic 'L' shape is not present which prevented us from using this method to find the parameter. The reason why the 'L' is absent from our curve is not clear, it could be because our inverse problem is close enough to a well posed problem for this method not to work. As a consequence, this method could not be used to find the parameter. An estimation of the parameter using the variance of the noise [START_REF] Nikolas | Methods for choosing the regularization parameter and estimating the noise variance in image restoration and their relation[END_REF] was also investigated by gave no conclusive result. In the absence of any applicable method the regularization parameter had to be set manually. A rough estimate was found by running the algorithm on simulated data and the value was refined on the experimental data. This solution is not ideal and will be discussed in the next sections.

Reconstruction results

To reconstruct the absorption of the phantom, I used a conjugate gradient algorithm as described in figure 5.5 with a first order Tikhonov regularization. The number of ultrasound positions in the forward model was set to 38 × 70 so that one pixel corresponds to 1 mm. This is almost twice as large as the experimental ultrasound but it greatly reduced the computational time. The data were down-sampled accordingly and value of the regularization parameter, λ, was manually set to 10 -6 . Figure 5.9 shows different iterations of the algorithm. which starts by building low frequencies in the image and the finer details appear later. In this case, the algorithm converged in 32 steps but it was kept running for more iterations to show the stability of the solution. Figure 5.10 shows the evolution of the objective function on a log scale. After iteration 32, the stable solution is found and the reconstructed µ a remains the same. The computation of the 200 steps took 30 hours but in a practical case, the algorithm would have reached the stopping criterion after step 32 which corresponds to a processing time of 4.5 hours.

The reconstructed absorption values are almost four times smaller than the real values but the ratio between the two inclusions appears correct as can be seen on figure 5.11 which shows profiles along two directions of the image. With peak values of 0.16 and 0.098 cm -1 the reconstructed µ a ratio is equal to 0.61 while the experimental ratio is equal to 0.5/0.8 = 0.63. This is very interesting as it gives information about the relative absorption strength of the inclusions. On the raw experimental data, it is possible to distinguish the two inclusions but given their low absorption, the contrast is weak. The difference data, shown on figure 5.9, distinctly shows the two inclusions but it is impossible to assess their relative absorption or even to determine which inclusion absorbs more. With the reconstruction, this becomes possible even though the exact values are still unknown. However, it could be possible to use an absorber with a known µ a in order to deduce the other unknown absorption coefficients, opening a door towards quantitative acousto-optic imaging.

Moreover, the size of the inclusions (3 × 5 mm) is correctly reconstructed by the algorithm as can be seen on the profiles on figure 5.11. The longer dimension of the less absorbing one is slightly underestimated (4 mm instead of 5) but this gives a more precise idea about the size of the absorbers. Despite the rectangular shape of the inclusions, their appear as Gaussian ellipsoids in the reconstruction. This is mostly due to the shadowing effect which is inherent to acousto-optic imaging. Since the technique measures the local light intensity by measuring the number of tagged photons that reach the detector for a given acoustic location, the presence of a small absorber will affect the AO signal on a much wider area. If the ultrasound focus is located in close proximity to the absorber, the path of the photons tagged in this area is likely to cross the inclusion, thus reducing the measured light fluence in this area. This effect will thus blur the boundaries of the inclusion and make them appear like Gaussian ellipsoids. It is known in acousto-optics and it is modelled by the diffusion equations but it appears that the rectangular shape of the inclusions is not recovered with the inverse problem solution.

The reconstructed image also contains a strong artefact near the light detector. This was already visible in the difference data image on figure 5.9 and might come from a slight movement of the phantom between the measurement and the baseline acquisition. Since it is located on the boundary of the medium, it is easy to discard it as an artefact because it would have been visible with naked eye when looking at the sample.

Overall, these reconstructions are very promising and we believe that they are the first step towards quantitative acousto-optic imaging. Following these results, we tried to improve our algorithm in order to obtain the real values of the absorption coefficient. Numerous reconstructions were ran with different regularization parameters on several set of experimental data but despite our efforts, we could only obtain qualitative results. I will discuss our hypotheses on the issues preventing us from obtaining quantitative measurements in the next section.

Reconstruction outside of the model hypotheses

To respect the approximations inherent in the diffusion approximation, the absorption coefficient of the inclusions need to be very small compared to the scattering coefficient of the medium. Nevertheless, in practice, absorbing objects in biological tissues like melanomas or tumour are strongly absorbing light. To test if our method could work in these conditions, we applied the algorithm to a set of data acquired on a phantom containing three strongly absorbing inclusions of µ a equal to 4.8, 3.2 and 1.6 cm -1 . The result obtained with a conjugate gradient algorithm, first order Tikhonov regularization and a hyper parameter λ equal to 10 -6 is presented on figure 5.12. As expected, the values of the reconstructed µ a are far from the ground truth but once again, the ratio between the different inclusions is respected and the dimensions are correctly recovered. It is quite interesting since it shows that despite the incorrect forward model for these values, it is possible to obtain qualitative results and to distinguish the absorption of the three inclusions. This could prove useful for medical interpretation of biological samples with several absorbing inclusions even though some work remains to be done to obtain quantitative results.

Limitations of the current algorithm

As shown in figure 5.11, the reconstruction algorithm is able to reconstruct the inclusions but the values of the reconstructed µ a are far from the values measured experimentally. The reasons for this difference are not fully understood and the problem needs to be further studied. Moreover, some requirements of our approach, such as the acquisition of a baseline measurement, can be a bit restrictive. I will here present some hypotheses we formed regarding the issue of quantitativeness and discuss some of the limitations of our current method.

Towards quantitative reconstruction

One of the objectives when solving the acousto-optic inverse problem is to quantitatively reconstruct the optical parameters to provide more information about the scattering medium. However, in its current state, the algorithm is not able to produce such results. The reasons for this have not been extensively studied due to lack of time but we formed several hypotheses which will be discussed here. The first hypothesis that can be made is that the simulated AO data generated by the forward model are too different from the experimental data. The two coupled diffusions equations were chosen to model the data because of the computational simplicity of this approach. Few experimental parameters are needed to simulate images that are a good estimate of the experimental acquisition and the computational time is significantly smaller than other approaches such as Monte-Carlo. The previous results show that it was a sensible approach for qualitative reconstruction but this method could lack accuracy to obtain quantitative results. As shown on figure 5.7, the experimental and simulated difference data appear quite different from each other. The overall scaling difference is compensated by the model correction but the difference in the light distribution can influence the efficiency of the reconstruction.

Several factors could be responsible for this difference. The ultrasound distribution, for example, is simulated as a set of homogeneous square pixels. In practice, the ultrasound wave propagates in the sample and the pressure is not evenly distributed. Along the longitudinal direction, the ultrasound extension is given by the number of sinusoidal cycles and it can be approximated as a spatial Heavyside function even though the pressure oscillate in this zone. Along the transverse direction, however, the restricted width of the spot is created by the interference of the different waves emitted by the piezoelectric elements of the probe as explained in chapter 2. Consequently, approximating the spatial distribution as a sharp square might create differences in the image. This could be investigated by changing the numerical ultrasound distribution to make is more similar to the experiment.

Another source of mismatch could come from the tagged photons detection. In our forward model, only the photons shifted by the ultrasound frequency are considered whereas experimentally, all the different harmonics are detected. Since the setup actually detects untagged photons, the contribution of all frequency shifts are added. The number of photons shifted by higher harmonics of the light is extremely small compared to the photons tagged by the first harmonic but this could have an effect on the image formation.

The diffusion approximation to the RTE imposes several approximations, among which the fact that absorption is negligible compared to scattering. To respect this hypothesis, the absorption coefficient of the experimental inclusions were chosen to be as low as possible. However, the experimental phantom is made of Agar and Intralipid-10% which weakly absorb light at this wavelength. Consequently, if the absorption of the inclusions is too low, it is challenging to detect them experimentally. Similarly, if the scattering coefficient of the phantom is too high, the number of detected photons is reduced and the image quality deteriorates. A compromise was found and the values were set to 12 cm -1 for the reduced scattering coefficient and 0.8 and 0.5 cm -1 for the absorption of the inclusions. The absorption was thus small compared to scattering but cannot be considered negligible and this may affect the simulation accuracy.

The last point of our forward model that might strongly affect the reconstruction is the fact that the light propagation is simulated on a 2D plane. This 2D plane was chosen to match the experimental plane of imaging which corresponds to the plane of the ultrasound. If it is true that the photons are only tagged in this plane because of the ultrasound, restricting their propagation to this plane affects the simulated image. We chose to use a 2D simulation because of the computational time needed to simulate the data in 3D. Considering that the forward model is calculated between 3 and 15 times per iterations during the line search, a 3D simulation would greatly increase the reconstruction time. However, the data simulated with a 3D domain look closer to the experimental data and it might be necessary to simulate in 3D to obtain quantitative results.

These hypotheses need to be further studied to correctly identify the issues and correct the model but we believe that this method is the first step towards quantitative acousto-optic imaging. Nevertheless, a few points of the technique might be limiting its application to medical imaging.

Necessity of a baseline measurement

Our reconstruction method relies on the acquisition of a baseline measurement in a homogeneous phantom to compensate for the unknown tagging efficiency and the spatial inhomogeneities of the experimental image. This baseline was straightforward when imaging custom made phantoms but might be more difficult to implement for more complex media. In the perspective of in vivo imaging, the acquisition of a baseline with the exact same properties would be almost impossible. If the variation in the image only comes from the imaging setup, however, this baseline could be acquired on an homogeneous phantom that reproduces the properties of the biological medium.

Low absorption hypothesis

The diffusion approximation requires that absorption is small compared to the scattering. In addition to the difficulty it raises for quantitative reconstruction on phantoms, this hypothesis can be very restrictive in the prospect of in vivo imaging. Some absorbing objects in biological tissues, such as tumours or melanomas, strongly absorb light. To quantitatively determine the absorption coefficient of such objects, the diffusion approximation would not be ideal. For this purpose, other forward modelling options can be thought of, such as a non approximated radiation transfer equation or a correlation transport based approach [START_REF] Sakadžić | Correlation transfer and diffusion of ultrasoundmodulated multiply scattered light[END_REF], solved using Monte-Carlo methods. The latter can be implemented on a GPU to accelerate the computing [START_REF] Powell | Highly parallel monte-carlo simulations of the acoustooptic effect in heterogeneous turbid media[END_REF].

Regularization

The last limitation of the current technique comes from the application of regularization. As mentioned previously, no conclusive way of choosing the regularization parameter was found and it was manually set. A given value of the parameter gives similar results when the problem is solved with simulated data but it is not always the case. Moreover, finding the correct parameter using simulated data requires numerous runs of the algorithm which takes a long time. Qualitative reconstruction are moderately affected by the regularization parameter as the ratio between the different inclusions absorption is quite stable. However, if the objective is to reconstruct the absorption quantitatively, the regularization will have a strong influence on the end result. For this reason, it is important to find a better way to choose this parameter.

Next steps for quantitative acousto-optic imaging

Despite the current limitations of the algorithm we developed, the solving of the acousto-optic inverse problem represents an important step towards quantitative AOI. Once the aforementioned issues preventing from obtaining quantitative images have been overcome, several options can be considered to improve the method.

3D reconstruction

The simulation of light propagation in a three dimensional space might be the solution to reach quantitative reconstruction but it can also increase the versatility of the technique. Not only will it allow for the reconstruction of a plane orthogonal to the optical axis but it can also open the way to a full 3D reconstruction of a sample.

At the moment, the source and the detector are required to be positioned on the boundary of the domain. Since we are restricting the simulation to a plane, this plane needs to contain the source and the detector. This limits the application of the reconstruction to a plane which contains the optical axis. To apply the reconstruction to other planes of the domain, a 3D forward model needs to be implemented. Such a model could also allow for a reconstruction of a 3D domain. For this, a 3D image needs to be acquired experimentally by translating the acoustic probe and stacking the 2D images.

A forward model simulating the propagation of light in three dimensions has been used by Powell et al. [START_REF] Powell | Quantitative reconstruction of absorption and scattering coefficients in ultrasound-modulated optical tomography[END_REF] on simulated data and similar implementation has been started for experimental data using a GPU. However, it will requires more time to be properly working.

Multiple pressure measurements

To correct the unknown experimental parameters such as the tagging efficiency, we had to apply some corrections to our model. These corrections help the reconstruction but it would be preferable to measure these parameters. The number of photons tagged by the ultrasound is proportional to the local light fluence but also to the square of the acoustic pressure. Consequently, by acquiring two AO images with two different acoustic pressures, it could be possible to retrieve the experimental value of the tagging efficiency as a function of space. This could potentially reduce the need for model corrections and increase the efficiency of the reconstruction.

Reconstruction of scattering coefficients

Reconstruction of the absorption coefficient is the first step but as mentioned previously, it is theoretically possible to reconstruct both absorption and scattering coefficient. Powell et al. [START_REF] Powell | Quantitative reconstruction of absorption and scattering coefficients in ultrasound-modulated optical tomography[END_REF] showed that this is possible with simulated data but that a cross talk could exist between the influence of µ a and µ ′ s on the AO data. The authors show that it is possible to reduce the illposedness of the problem by using several sources and detectors to obtain more information about the medium. Experimentally, this could be done by moving the source and the detector around the sample or by rotating the sample to image several angles.

This could prove to be a invaluable asset for acousto-optic imaging as it would allow for the distinction between absorption and scattering contrasts in complex media.

Conclusion

The solution of inverse problem has been widely studied by the numerical optimisation community to calculate a parameter of interest from a set of indirect measurement. This process is used in many fields, including medical imaging for techniques such as diffuse optical tomography or CT-scans. The application of this technique to the acousto-optic problem could be a way to obtain quantitative information about a medium of interest. In this chapter, I described the algorithm that we implemented to reconstruct a map of the absorption coefficient in a sample using experimental measurements. This algorithm uses a custom forward model consisting in two coupled diffusion equations. To compensate for several unknown experimental parameters, such as the ultrasound tagging efficiency, corrections have been applied to our model. I showed that with the current algorithm it is possible to reconstruct the size of absorbing inclusions as well as their relative absorption coefficient. Exact measurement of the absorption was not achieved but I proposed several hypotheses that can be studied to achieve this goal. We believe that this algorithm is the first step towards quantitative acousto-optic imaging which could greatly improve the interpretation of AO images in the prospect of medical diagnosis.

Chapter's Keypoints

• Inverse problems link parameters of interest to a set of measurements through a forward model.

• To recover the parameters from the measurements, an approximated forward model and a reconstruction algorithm are needed.

• Inverse problems are often ill-posed, which means that a small change in the measurement can result in a large change in the estimated parameter. This can be prevented by using regularization.

• Reconstruction algorithms minimize the difference between the measurements and a simulation of the measurements created by the forward model.

• Such model based reconstruction algorithms are used for CT-scans or diffuse optical tomography. They can be used for acousto-optic imaging in order to recover the absorption and scattering coefficients.

• The algorithm presented in this manuscript uses a forward model consisting of two coupled diffusion equations. These equations are approximations of the radiative transfer equation and are valid for µ a ≪ µ ′ s .

• Model corrections were required to use this algorithm with experimental data such as the acquisition of a baseline measurement.

• AO images obtained on experimental phantoms containing absorbing inclusions were used with the algorithm.

• The ratio between the absorption coefficients of the different inclusions as well as the size of the inclusions was recovered.

• Similar results were obtained for inclusions with absorption coefficients outside of the diffusion approximation hypotheses.

• The current algorithm needs to be improved to perform quantitative reconstruction but it is a first step towards acousto-optic imaging.

Conclusions and perspectives

Medical imaging is a crucial tool to help physicians establish a diagnosis. Current imaging technique such as magnetic resonance imaging (MRI) or ultrasound scans have been greatly improved over the last decades in terms of efficiency and sensitivity. Yet some pathologies, such as early tumours, remain difficult to detect. Optical imaging is sensitive to both absorption and scattering and can thus detect structural and metabolic changes. Optical contrast could then be an asset for medical diagnosis. However, in-depth imaging of biological tissues with light is challenging because of multiple scattering. This phenomenon affects the propagation of light and after a few transport mean free paths, the origin of light is lost. This distance corresponds to a few millimetres in biological samples, and classical optical techniques relying on geometrical optics fail to image deeper than this length. Some approaches have been developed to image at greater depths such as diffuse optical tomography which relies on a model based inversion to reconstruct an image from measurement on the boundary. Other methods consist in combining light with a different type of wave.

Acousto-optic imaging is a technique that couples light and ultrasound to image highly scattering samples at centimetre depths. As the ultrasound propagates in the illuminated area, some photons can be phase modulated at the frequency of the ultrasound. We have seen in this manuscript that the number of these tagged photons is proportional to the local light fluence integrated over the focus of the ultrasound. Consequently, measurement of the tagged photons flux gives access to the light fluence in the medium with an ultrasound resolution of the order of the millimetre. We have seen that several approaches exist to form an image from these measurements and that, if the detector is fast enough, the propagation of the ultrasound can be followed in real time and a line of the image be formed.

The detection of the tagged photons is the critical part of acousto-optic imaging for several reasons. First, the number of tagged photons is small compared to the total number of photons because of the limited dimensions of the tagging area. Second, the frequency shift created by the ultrasound is extremely small compared to the optical frequency. Third, because of multiple scattering, a speckle pattern is formed at the output of the sample. Due to speckle properties, each grain possess a random phase and the integration of numerous speckle grains does not increase the signal to noise ratio. We have seen that several methods can be used to detect those photons such as ultra narrow band filtering techniques or self-adaptive interferometry. The latter requires an holographic medium that can be of diverse nature. In this manuscript, self-adaptive wave-front interferometry in photorefractive crystals was described in more detail.

Photorefractive crystals are both photoconductive and electro-optic which means that the refractive index of the medium can be modulated by light. This effect is used to build an hologram from the interference of the speckle pattern and a reference beam. The reference is then diffracted on the refractive index grating and its wave-front is adapted to that of the speckle. This process is called two-wave mixing and we saw that it can be used to extract the tagged photons component from the background light. We have also seen that such crystals possess two main characteristics: the photorefractive response time and the photorefractive gain. The former corresponds to the build-up time of the hologram while the latter quantifies the energy transfer between the reference and the signal beams. Those quantities vary from one crystal to another and the characterisations of two different crystals, SPS and ZnTe, were presented in this manuscript. SPS crystals have been used in Institut Langevin for several years and proved to be efficient for applications such as ex vivo imaging. Whilst their photorefractive gain is high, their response time is above the millisecond, which corresponds to the speckle decorrelation time in vivo. For this reason, a ZnTe crystal was studied by despite its short response time that can be as low as 0.2 ms, its gain was not sufficient for acousto-optic imaging in the current state. SPS crystals were thus used for the rest of the work presented in this thesis. However, an electric field could potentially be applied on a ZnTe crystal to enhance its gain which could increase the potential of such crystals.

The majority of published work on acousto-optic imaging was performed in transmission mode, with the source and the detector on opposite sides of the sample. Whilst this is useful for imaging phantoms, ex vivo tissues or compressed breast, this configuration is not optimal in the pursuit of in vivo imaging. For this reason, a reflection configuration, with both source and detector on the same side of the sample, was studied. Using a custom Monte-Carlo algorithm, this configuration was first analysed numerically. We have seen that the imaging depth depends on the distance between the source and detector and that it is possible to image a medium as deep as two centimetres. Above this depth, the signal to noise ratio is too low. Numerical results were completed by experimental images obtained with a photorefractive-based detection setup coupled to a commercial ultrasound scanner. Experimental conclusions were similar and a new configuration was presented in which the ultrasound probe is located between the light source and detector. This so-called 'hand-held probe' configuration was tested on phantoms and in vivo on a volunteer's wrist. Despite the long imaging time of several minutes which implied unwanted movements of the arm, images were obtained. The penetration depth was limited to half a centimetre but it opened the way towards in vivo acousto-optic imaging.

Acousto-optic imaging maps the local light fluence inside scattering samples. For this reason, it is sensitive to absorption and scattering but the shape of the illuminated area can lead to misinterpretation of the image. Moreover, in the presence of several absorbers, for instance, it is not always possible to assess their relative strength and impossible to measure their absolute absorption coefficient. One solution to measure these values is to apply a model-based reconstruction to acousto-optic data in order to reconstruct the absorption and scattering coefficient. To do so, a forward model is needed to simulate light propagation as well as an algorithm to minimize the difference between the experimental data and the output of the forward model. This approach has been demonstrated on simulated data and it has been shown in this thesis that it can be applied to experimental data to reconstruct the absorption coefficient. We have seen that some model corrections are needed to accommodate the differences between experiments and simulation and that, even with a simple forward model, it was possible to deduce the relative absorption strength of several inclusions. Whilst more work is needed to obtain quantitative measurement and to reconstruct both absorption and scattering, these reconstructions represent the first step towards quantitative acousto-optic imaging.

Future work in acousto-optic imaging may lead to the development of a bimodal platform combining conventional ultrasound with acousto-imaging for in vivo imaging using a hand-held probe containing a transducer array and optical fibres. This platform could be used for the detection of early tumours in the breast or liver once the penetration depth will be improved. On the other hand, the improvement of the model-based algorithm could lead to the quantitative measurement of both absorption and scattering coefficients in scattering medium. The combination of in vivo imaging and quantitative measurements would be a considerable asset for medical diagnosis.

Résumé

L'imagerie médicale est un outil crucial pour le diagnostic médical. Les techniques actuelles comme l'IRM ou l'échographie sont très performantes mais pour la détection de certaines pathologies comme les tumeurs précoces, ces méthodes ne sont pas adaptées. Etant donné que les tumeurs absorbent fortement la lumière, les techniques optiques semblent adaptées pour leur détection. Cependant, la diffusion multiple dans les milieux biologique empêche d'utiliser des méthodes optiques classiques au-delà de quelques millimètres. L'imagerie acousto-optique est une technique qui couple lumière et ultrasons pour imager en profondeur dans les milieux diffusants. Lors de leur propagation, les ultrasons modulent la lumière et génèrent des photons marqués. La quantité de photons marqués est proportionnelle à l'intensité lumineuse locale au point de focalisation des ultrasons. Cette méthode permet donc de reconstruire une carte de l'intensité lumineuse à une profondeur de quelques centimètres en balayant les ultrasons. L'étape difficile est la détection de ces photons et plusieurs techniques existent comme le filtrage hyper fin ou l'interférométrie auto-adaptative. Dans cette thèse, nous présentons un dispositif d'imagerie acousto-optique utilisant des cristaux photorefractifs pour détecter les photons marqués par adaptation de front d'onde holographique. Plusieurs types de cristaux peuvent être utilisés et nous présentons ici la caractérisation expérimentale de deux cristaux : le SPS et le ZnTe. La plateforme d'imagerie a aussi été améliorée pour fonctionner en configuration de réflexion dans laquelle la source et le détecteur de lumière se trouvent du même côté de l'échantillon à imager. Des résultats numériques et expérimentaux sont présentés dans ce manuscrit. Nous proposons aussi de combiner ce dispositif avec un algorithme de reconstruction basé sur une modélisation de la propagation de la lumière afin de reconstruire des cartes de coefficients d'absorption. Les résultats expérimentaux présentés représentent la première étape vers l'imagerie acousto-optique quantitative.

Mots Clés

Imagerie acousto-optique, ultrasons, milieux diffusants, holographie photoréfractive, imagerie médicale, problèmes inverses.

Abstract

Medical imaging is a crucial tool for medical diagnosis. Current imaging techniques such as magnetic resonance imaging or ultrasound scans are efficient and sensitive. However, to detect some pathologies such as early tumours, other methods are needed. Since those tumours strongly absorb light, optical methods would be suited to detect them but in biological tissues light is highly scattered. This restricts the imaging depth of classical optical techniques to a few millimetres in such media. Acousto-optic imaging is a technique that couples light and ultrasound to obtain optical images deep inside scattering samples. As the ultrasound propagates in the sample, it modulates light and creates tagged photons. The number of these tagged photons is proportional to the local light fluence at the focus of the ultrasound. This method can thus map the light fluence up to a few centimetres deep inside scattering media by scanning the ultrasound. The critical step is the detection of these photons and several approaches can be used such as ultranarrow band filtering or self-adaptive interferometry. In this thesis we present an acousto-optic setup using photorefractive crystals to perform selfadaptive wave-front holography for the detection of tagged photons. Several crystals can be used and we present here an experimental characterisation of two crystals: SPS and ZnTe. The imaging setup was also improved to function in a reflection configuration, with source and detector on the same side of the imaged phantom and some experimental results are presented in this manuscript. We also propose to combine this experiment with a model-based reconstruction algorithm to reconstruct maps of the absorption coefficient from experimental data in the pursuit of quantitative acousto-optic imaging.

•

  AOI: Acousto-Optic Imaging • CG: Conjugate Gradients • CT: Computerized Tomography • DOT: Diffuse Optical Tomography • GPU: Graphic Processing Unit • HIFU: High Intensity Focused Ultrasound • MC: Monte Carlo • MO: Measurement Operator

Figure 1 . 1 -

 11 Figure 1.1 -Absorption coefficient spectra of endogenous tissues chromophores (taken from [11]). Hb0 2 and HHb denote respectively oxygenated and de-oxygenated blood. The window of low absorption between 700 and 1200 is commonly called the Optical Therapeutic Window

Figure 1 . 2 -

 12 Figure 1.2 -Isotropic and anisotropic scattering

Figure 1 . 3 -

 13 Figure 1.3 -Schematic of scattering. The scattering mean free path represent the average distance betweentwo scattering events. The transport mean free path corresponds to the average of the blue arrows lengths. This represents the average distance after which the light has lost the memory of its initial direction and polarization.

Figure 1 .Figure 1 . 4 -

 114 Figure 1.4 -Schematics of the different regimes undergone by light in a scattering sample. (0) represents the incident light. (1) represents the ballistic photons travelling in a straight line. (2) corresponds to a photon that experiences only a single scattering event.(3)(4) and (5) represents two particular cases of multiple scattering: (3) has a snake-like trajectory which appears in highly anisotropic scattering media and (4) and (5) are diffuse photons. (5) represents the particular case of back-scattered photons that are sent back in the direction of the source.

Figure 1 . 5 -

 15 Figure 1.5 -The scattering sample is illuminated by an electromagnetic plane wave. When propagating along a given optical path, the electric field accumulates a phase shift which depends on the length of the path, s p

Figure 1 . 6 -

 16 Figure 1.6 -Picture of a speckle pattern. Taken from Wikipedia

Figure 1 . 7 -

 17 Figure 1.7 -(a) Schematic of a typical in-line OCT setup and (b) typical image obtained with this setup (longitudinal section of a human retina and optic disk in vitro)(both taken from [30]). (c) Full field OCT setup and (d) Typical image obtained with this setup (en face image of a fixed human oesophagus epithelium at 100 µm below the surface) (both extracted from [31])

Figure 1 . 8 -

 18 Figure 1.8 -Schematics of Near Infra-red spectroscopy and Diffuse Optical Tomography. (a) Principle of NIRS. The path most travelled by light forms a banana shape. The depth of the banana is linked to the distance between source and detector. (b) Schematics of a DOT setup for brain imaging.Extracted from[START_REF] Elizabeth | Optical brain imaging in vivo: techniques and applications from animal to man[END_REF]

Figure 1 . 9 -

 19 Figure 1.9 -Schematics of photoacoustics. The pulsed light heats the absorbers up and they emit an ultrasonic wave detected by the transducers. Elongated objects only emit orthogonally to their long axis

Figure 1 . 10 -

 110 Photodetector

Figure 2 . 1 -

 21 Figure 2.1 -Schematics of the acousto-optic effect. The ultrasonic wave creates a refractive index grating on which light is diffracted. The light frequency, f L , is shifted by multiples of the ultrasound frequency, f US , depending on the diffraction order.

Figure 2 . 2 -

 22 Figure 2.2 -Schematic of the acousto optic effect in a scattering medium. The photons that interact with the ultrasound are shifted in frequency and two side-bands appear at the laser frequency, f L , plus or minus the ultrasound frequency f US . Higher harmonics exist but have considerably less energy.

Figure 2 . 3 -

 23 Figure 2.3 -Acousto-optic imaging of an ex vivo human liver biopsy containing two tumours (circled in white). Taken from [39]. (a) is the ultrasound image where the tumours cannot be seen. (b) is the acousto-optic image. (c) represents the optical contrast which shows the two tumours. (d) is a picture of the liver sample showing the black tumours.

Figure 2 . 4 -

 24 Figure 2.[START_REF] Luyten | Metabolic imaging of patients with intracranial tumors: H-1 MR spectroscopic imaging and PET[END_REF] -Schematic of the formalism used for the light modulation derivation. For a given path p, the position of the j th scatterer is denoted r j and the length of the j th scattering event is l j . K US represents the wave vector of the acoustic wave.

Figure 2 . 5 -

 25 Figure 2.5 -Schematic of the image formation in acousto-optic imaging with pulsed ulrasound. As the ultrasound pulse travels through the medium, the number of tagged photons is recorded to construct a line of the image. The focus is then mechanically or electronically translated and the lines are stacked to form an image.

Figure 2 . 6 -

 26 Figure 2.6 -Schematic of AOI with plane-wave ultrasound. A set of pulsed plane waves propagating with different angles are sent and the amount of tagged photons is recorded along propagation. A tomographic reconstruction is then used to combine the data in an acousto-optic image. The number of angles can be limited by the characteristics of the probe; in this case the image might be distorted in the transverse direction.

Figure 2 . 7 -

 27 Figure 2.7 -(a) Absorption spectrum of a crystal with a spectral hole burned at ω P . ∆Γ I is the inhomogeneous bandwidth and ∆Γ H is the homogeneous linewidth. Taken from [88]. (b) Schematics of the atomic energy levels of Thulium in a YAG matrix. Dotted line represent different relaxation processes.

Figure 2 . 8 -

 28 Figure 2.8 -Schematic of the Two-Wave Mixing process. The interference between the plane wave reference beam and the scattered signal beam creates a grating in the holographic medium on which the reference beam diffracts in the direction of the transmitted signal. The wave-fronts of the diffracted reference and transmitted signal are identical and can interfere constructively.

Figure 3 . 1 -

 31 Figure 3.1 -Representation of the photorefractive effect. Under a spatially inhomogeneous illumination, a displacement of charges occurs which in turn creates a local electric field. Through the electrooptic effect, the electric field modulates the refractive index.

Figure 3 . 2 -

 32 PhotoexcitationThermal excitation Recombination Diffusion

Figure 3 . 3 -

 33 Figure 3.3 -Schematic of the two-waves mixing process. The interference between the plane wave reference beam and the scattered signal beam creates a grating in the holographic medium on which the reference beam diffracts in the direction of the transmitted signal. The wave-fronts of the diffracted reference and transmitted signal are identical and can interfere constructively.

Figure 3 . 4 -

 34 Figure 3.4 -Representation of the acousto-optic signal in positive or negative gain configurations. The dotted lines represent the static signal, without ultrasound modulation and the solid line represent the overall signal detected on the photo-diode in a noiseless case.

Figure 3 . 5 -

 35 Figure 3.5 -Schematic of the setup used for the measurement of the spectral bandwidth of photorefractive crystals. The laser beam exiting the Master Oscillator Power Amplifier (MOPA) is split into two by a polarizing beam splitter (PBS) and a half-wave plate (λ/2) to control the relative intensities. Both beams sent through two Acousto-optic modulators (AOM).A small detuning and a phase modulation is applied to the reference beam which recombines with the signal beam on the Photorefractive crystal (PRC). The signal is measured on a photo-diode (PD) connected to a lock-in detection.

Figure 3 . 6 -

 36 Figure 3.6 -Example of the lock-in detection signal for a phase modulation frequency of Ω = 20 kHz on a ZnTe crystal.

Figure 3 . 7 -

 37 Photorefractive bandwidth (Hz)

Figure 3 .

 3 Figure 3.8(b) shows the results of the gain measurement where we can clearly see the saturation regime which is already reached for an average fluence of 0.4 W.cm -2 . The photorefractive gain is noticeably weaker than that of the SPS crystals with a maximum value of 0.25 cm -1 . The results of the bandwidth and response time measurements are shown on figure3.8(c). These results confirm that ZnTe crystals can have a response time smaller than the millisecond for reasonable fluence. For higher values of the fluence, the build-up time of the hologram in this crystal can be as low as 0.2 ms.

Figure 3 . 8 -

 38 Photorefractive bandwidth (Hz)

  figure 4.1(a). The number of elements used, m, defines the finite aperture and thus the size of the focal spot. This process, called beamforming, can also be used in reception; after recording all the echoes coming from the whole medium, they are numerically refocused to their rightful place in order to obtain a line of the image. The focus can then be scanned electronically along the probe to reconstruct a 2D image by stacking all the lines as shown on figure 4.1(b).

Figure 4 . 1 -

 41 Figure 4.1 -Principle of conventional ultrasound imaging. (a) Schematic of a linear transducer array and a parabolic delay line of focusing at position [x f , y c ]. The same kind of delay line can be applied in reception to record signal coming from a specific area. (b) Principle of image formation. The focus is electronically scanned over the probe to reconstruct a 2D ultrasound image.

Figure 4 .

 4 2 shows an example of photons trajectories inside such a medium.

Figure 4 . 2 -

 42 Figure 4.2 -Schematic of the Monte Carlo simulation. Dark ellipses represent the source and detector surfaces.A photon packet is generated on the source and propagates through the medium in a random walk fashion. The average step length is the scattering mean free path. Only the data relative to photons exiting the medium on the detector will be saved. Small light ellipsoids represent the ultrasound tagging areas. Every time a photon packet crosses one of these areas, a counter is incremented. Each pixel of the final image will correspond to one such area, the amount of signal for each pixel being the total energy of photons which crossed the corresponding area. For the sake of clarity, only a 3x3 pixel grid was represented here. In this image, the counter of the blue ellipsoid will be incremented whereas the counter of the red one will not because the red photon packet did not reach the detector.

Figure 4 . 3 -

 43 Figure 4.3 -Workflow of the Monte Carlo algorithm. Epsilon represent the energy value below which the photon packet is considered absorbed by the medium

Figure 4 . 4 -

 44 Figure 4.4 -Monte Carlo simulation of the acousto-optic images in the YZ plane for different source-detector distances (1 cm to 6 cm). The intensity represents the energy of the photons tagged in a given pixel which reached the detector normalized by the total incident energy. Red arrows represent the source and detector positions. Since both source and detector are wide the maximum intensity does not correspond to the source or detector central position.

Figure 4 . 5 -Figure 4 . 6 -

 4546 Figure 4.5 -(a) Monte Carlo simulation of the acousto-optic images in the YZ plane for a source-detector distance of 3 cm (b) Normalized profiles of the simulated acousto-optic signal along the dotted line. A moving average over 10 points was applied to each profile and it was then normalized by its maximum value for representation purposes.

Figure 4 . 7 -

 47 Figure 4.7 -Simulation of the exit position of photons on the detector for different distances D between the source and the detector. The 8 × 8 mm 2 area was divided into 100 × 100 pixels and the energy reaching each pixel was recorded and normalized by the total incident energy. For a distance D smaller than 3 cm, most photons are collected on the side of the detector closest to the source (positive-x) whereas after 3 cm, the distributions is more homogeneous.

Figure 4 . 8 -

 48 Figure 4.8 -Representation of the exiting directions of photons reaching the detector for different distances Dbetween the source and the detector. Each vector originates from 0,0,0 and its length represents the total energy exiting the medium at this specific angle. Units are arbitrary and only serve to represent the angle distribution. In this representation, the scattering medium is located in the positive-z area and the source is located towards the positive-x. For D smaller or equal to 2 cm, some angles are more likely than others whereas for greater distances, it becomes almost isotropic.

Figure 4 . 9 -

 49 Figure 4.9 -Experimental set-up with optical fibres in reflection mode. Light goes in and out of the sample on the same side and ultrasound are sent from the top. MOPA: Master Oscillator Power Amplifier, PBS: polarizing beam splitter, MMOF: Multi-mode optical fibre, LCOF: Liquid core optical fibre, PRC: Photorefractive crystal, PD: Photodiode

Figure 4 . 10 -

 410 Figure 4.10 -Acousto-optic images in reflection mode for different source-detector distances. The source position remained the same and the optical fibre used for detection was moved with a translation stage. Red arrows represent the source and detector position.As the distance increases, the AO signal decreases but not at the same rate as in MC simulations. When the distance increases the shape of the signal changes and seems to be more intense in front of the detector. Banding on the first images is due to parasitic light and is not significant.

Figure 4 .

 4 Figure 4.11 -(a) Normalized profiles of the simulated acousto-optic signal along the mediator line between source and detector. A moving average over 10 points was applied to each profile and it was then normalized by its maximum value for representation purposes. (b) Depth of the maximum intensity of the acousto-optic signal for different source-detector distances

Figure 4 .Figure 4 .

 44 Figure 4.12 -(a) Experimental set-up in a hand-held probe configuration. Light goes in and out of the sample on both sides of the ultrasound probe. MOPA: Master Oscillator Power Amplifier, PBS: polarizing beam splitter, MMOF: Multi-mode optical fibre, LCOF: Liquid core optical fibre, PRC: Photorefractive crystal, PD: Photodiode (b) Picture of the two configurations. Source and detector are on the same side of the sample. US1 is the probe used to obtain images in figure 4.10. US2 is located between source and detector fibres. This configuration is the first trial of a handheld probe.

Figure 4 .

 4 Figure 4.14 -(a) 2D AO image obtained in the handheld probe configuration shown on figure 4.12 on a gel containing two absorbing inclusions (diameters: 3 mm and 2 mm) positioned 1 centimetre deep inside the sample. Source-detector distance: 3 cm. Both inclusions are visible but the contrast is weak because the sensitivity is low in this region due to low penetration of the light. Red arrows indicate the positions of the inclusions (b) Profile along the black dotted line. The dip caused by the inclusion is visible but hard to distinguish from noise. (c) Profile along the white dotted line. The dips caused by the inclusions are visible but the contrast is limited.

Figure 4 .

 4 Figure 4.15 -(a) Schematic of the setup with multiple detectors. A 1 cm wide beam of light illuminates the sample through fibre coupling and a collimator and the four ends of a liquid core optical fibre are placed on the four corners of a square centred on the input light. Photons collected by these fibres are then combined and sent to the photorefractive crystal. The ultrasound is sent from above the sample to image a plane orthogonal to the optical axis (b) Images obtained with this configuration at different depths in a phantom containing a cylindrical inclusion (diameter 3 mm, length 20 mm). The inclusion is clearly visible in each image but the contrast is poor for the last image.

Figure 4 . 16 -

 416 Figure 4.16 -In vivo images of a volunteer's wrist. (a) Acousto-optic image obtained after averaging 4000 times (b) Superposition of B-mode image and acousto-optic signal. Dotted lines represent the boundary of the AO image.

Chapter 4 .

 4 Acousto-optic imaging in reflection modeimaging.

Figure 5 . 1 -

 51 Figure 5.1 -Two representations of the quadratic form for A = [4, 2; 2, 8] and y = [6, -4] T . The minimum is found for x = [-1, 2] T in this case (a) Surface plot of the quadratic form (b) Contour of the quadratic form. Lines represent iso-value curves

Figure 5 . 2 -

 52 Figure 5.2 -With an initial guess of [-2,-2] the algorithm converges towards the optimal solution by going towards the local steepest direction

Figure 5 . 3 -

 53 Figure 5.3 -With an initial guess of [-2,-2] the algorithm converges towards the optimal solution in two A-orthogonal steps Figure 5.3 shows the evolution of the estimate of the parameter using a Conjugate Gradients algorithm for an initial guess x (0) = [-2, 2] T .

Figure 5 .

 5 [START_REF] Luyten | Metabolic imaging of patients with intracranial tumors: H-1 MR spectroscopic imaging and PET[END_REF] shows an example of this process.

Figure 5 .Figure 5 . 4 -

 554 Figure 5.4 -(a) Energy density of untagged photons for an isotropic source on the left hand side of the domain. (b) Energy density of the tagged photons for one given position of the ultrasound. (c) Simulated acousto-optic signal. Pixel number i represents the integration of Φ i 1 on the surface of the detector, here located on the right hand side of the domain.

Figure 5 . 5 -

 55 Figure 5.5 -Workflow of the reconstruction algorithm with the conjugate gradients method. N a is the total number of ultrasound positions.

Figure 5 . 6 -

 56 Figure 5.6 -Schematics of the experimental configuration. Light is sent inside the medium from the left and detected on the right while ultrasound are sent from the top to image a plane along the optical axis (a) Acquisition of the image of interest in the region with inclusions. (b) Acquisition of the baseline image in a region without inclusions.

Figure 5 .Figure 5 . 7 -

 557 Figure 5.7 -Experimental acousto-optic images (top) and forward model simulations (bottom) for a medium with a reduced scattering coefficient of 12 cm -1 and two inclusions with absorption coefficients of 0.8 and 0.5 cm -1 . The light source is located on the left hand side and the detector on the right hand side. (a) AO image in the plane of the two inclusions. (b) AO image in a plane without inclusions. (c) Difference between the two images. (d) Simulation of the AO image in the plane of the inclusions. (e) Simulation of the AO baseline image. (f) Difference between the two simulations.

Figure 5 . 8 -

 58 Figure 5.8 -Example of L-curve for a zero-order Tikhonov regularization.For ill-posed problems, the curve usually exhibits a 'L' shape and the optimal value for the regularization corresponds to the corner of the 'L'. Here, no 'L' can be observed.

Figure 5 . 9 -

 59 Figure 5.9 -Reconstructed µ a at different stages of the iterative conjugate gradients algorithm obtained using data shown on figure 5.7(a) and (b). The inverse problem was regularized with a zero-order Tikhonov regularization with an hyper parameter equal to 10 -6 .

Figure 5 . 10 -

 510 Figure 5.10 -Evolution of the objective function along iterations of the algorithm on a log scale. The stable solution is reached at iteration 32.

Figure 5 . 11 -

 511 Figure 5.11 -Reconstructed µ a at final stage of the algorithm. (a) 2D map of µ a in the medium. (b) Profile along the vertical dotted line. (c) Profile along the horizontal dotted line. (d) Picture of the scattering phantom cut along the plane containing the inclusions.

Figure 5 . 12 -

 512 Figure 5.12 -Reconstructed µ a at final stage of the algorithm for strongly absorbing inclusions. (a) 2D map of µ a in the medium. (b) Profile along the vertical dotted line. (c) Profile along the horizontal dotted line. (d) Picture of the scattering phantom cut along the plane containing the inclusions.

Table of contents

 of 1.1 Of the interest of optical imaging . . . . . . . . . . . . . . . . . . . . . . . . . Imaging with ballistic light . . . . . . . . . . . . . . . . . . . . . . . . . 15 Microscopy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 Time Gating . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

7 1.2 Light-matter interaction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 1.2.1 Absorption . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 1.2.2 Scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 Rayleigh Scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 Mie Scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.2.3 Orders of magnitude in biological tissues . . . . . . . . . . . . . . . . . 11 1.3 Light propagation in scattering media . . . . . . . . . . . . . . . . . . . . . . . 11 1.3.1 Propagation regimes . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 Ballistic light . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Single scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 Multiple scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 1.3.2 Speckle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 Expression of the multiply scattered electric field . . . . . . . . . . . . 13 Speckle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 1.4 Optical imaging in scattering media . . . . . . . . . . . . . . . . . . . . . . . . 15 1.4.1 Optical Coherence Tomography . . . . . . . . . . . . . . . . . . . . . . 16 1.4.2 Imaging with multiply scattered light . . . . . . . . . . . . . . . . . . . 17 Near Infra-Red Spectroscopy . . . . . . . . . . . . . . . . . . . . . . . . 17 Diffuse Optical Tomography . . . . . . . . . . . . . . . . . . . . . . . . 17 Multiwave imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18

Table of contents 3

 of3

.1 The photorefractive effect . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.1.1 Principle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 Photoconduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 Electrooptic effect . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 3.1.2 The band-transport model . . . . . . . . . . . . . . . . . . . . . . . . . 46 Space-charge field . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 Photoinduced variations of the refractive index . . . . . . . . . . . . . 48 3.1.3 Characteristics of the photorefractive effect . . . . . . . . . . . . . . . 50 3.2 Detection of the acousto-optic signal . . . . . . . . . . . . . . . . . . . . . . . . 51 3.2.1 Two-wave mixing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.2.2 Acousto-optic signal detection . . . . . . . . . . . . . . . . . . . . . . . 52 3.2.3 Experimental configuration . . . . . . . . . . . . . . . . . . . . . . . . 53 Detection of untagged photons . . . . . . . . . . . . . . . . . . . . . . . 54 Negative gain configuration . . . . . . . . . . . . . . . . . . . . . . . . 54 3.3 Choice of the crystal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 3.3.1 Measurement of the crystals characteristics . . . . . . . . . . . . . . . . 55 Photorefractive Gain . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 Photorefractive response time . . . . . . . . . . . . . . . . . . . . . . . 56 3.3.2 SPS vs. ZnTe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 3.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

  r 11 r 12 r 13 r 21 r 22 r 23 r 31 r 32 r 33 r 41 r 42 r 43 r 51 r 52 r 53 r 61 r 62 r 63

Table of contents 4

 of4

.1 Fusion between conventional ultrasound and Acousto-Optic Imaging . . . . . . 64 4.2 Monte Carlo Simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 4.2.1 The algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 4.2.2 Influence of source-detector distance . . . . . . . . . . . . . . . . . . . 67 Influence of incidence angle . . . . . . . . . . . . . . . . . . . . . . . . 69 4.3 Imaging in reflection mode . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 4.3.1 The imaging setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 4.3.2 Influence of source-detector distance . . . . . . . . . . . . . . . . . . . 76 4.3.3 Towards a handheld probe . . . . . . . . . . . . . . . . . . . . . . . . . 77 4.3.4 A multiple detector approach . . . . . . . . . . . . . . . . . . . . . . . 79 4.4 Towards In Vivo imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 4.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81

1.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21