Three Dimensional Widefield Imaging with Coherent Nonlinear Scattering Optical Tomography

Abstract

A full derivation of the recently introduced technique of Harmonic Optical Tomography (HOT), which is based on a sequence of nonlinear harmonic holographic field measurements, is presented. The rigorous theory of harmonic holography is developed and the image transfer theory used for HOT is demonstrated. A novel treatment of phase matching of homogeneous and in-homogeneous samples is presented. This approach provides a simple and intuitive interpretation of coherent nonlinear scattering. This detailed derivation is aimed at an introductory level to allow anyone with an optics background to be able to understand the details of coherent imaging of linear and nonlinear scattered fields, holographic image transfer models, and harmonic optical tomography.

1. Introduction

Optical microscopy permits the noninvasive acquisition of information that is revealed through light-matter interactions. These light-matter interactions are generally referred to as contrast mechanisms and come in many forms. The information carried by an optical contrast mechanism depends on the properties of the illumination light, the properties of the light produced by the contrast mechanism, and the details of the light detection. Most optical imaging systems rely on light that can be described approximately as classical, although there is a steadily growing body of work describing microscopy methods that exploit quantum correlations to enable the extraction of new information from objects.

The coherence properties of the light used for illumination and detection are also critical drivers of the properties of an optical imaging system. The classical theory of optical coherence is concerned with the statistical properties of classical fields that are treated as random variables. Although light is, in general, partially coherent, it is often suitable to describe light in the limiting case of either fully coherent or fully incoherent. Heuristically, we can describe coherent light as a field that is statistically similar across either temporal or spatial points on the field, whereas fully incoherent light lacks any correlation either along temporal or spatial displacements. Consideration of optical coherence is critical for understanding the broader context of optical microscopy and optical tomographic imaging.

While light propagation is not directly modified by optical coherence, field coherence strongly impacts the observed signal from a detector. For our purposes, a semiclassical model of light detection is suitable, where we consider generated photocurrents in a photodiode or photo-generated electrons detected by a camera chip. In all instances treated herein, we assume detector integration times are longer than the temporal coherence times of the light fluctuations. As such, detected light intensities are inferred from a long-time average of the incident optical field’s instantaneous intensity. To provide a consistent framework for our discussion of coherent tomographic imaging, we briefly review optical imaging theory to ensure that the reader is familiar with the notation used in this treatment.

In this chapter, we focus on imaging systems that can be described in a classical optical formalism that uses coherent nonlinear scattering as a contrast mechanism. Coherent nonlinear scattering exploits the microscopic properties of materials that exhibit a nonlinear dipole in response to a sufficiently strong incident field [1, 2, 3, 4, 5]. While materials can produce higher-order responses than a nonlinear dipole, the high field strength required generally precludes the use of higher-order terms to prevent damage to the object under study. It is suitable to describe the sensitivity of the nonlinear response as a nonlinear susceptibility tensor that is obtained from a Taylor expansion of the dipole response to the applied electric field. While the tensorial nature of nonlinear response depends on the incident fields and the distribution, we will suppress the vectorial dependence of both the coherent nonlinear light-matter interactions and the light scattering [3]. Within this approximation, coherent nonlinear scattering is described by a scalar field. For the purposes of imaging and tomography, we must then build a model for the propagation of the scattered scalar field through an imaging system, the detection of that light, and the processing required to obtain a microscopic or tomographic image.

This chapter is organized as follows. Coherent imaging theory is outlined and the application to tomographic imaging using coherent scattering is described. This section will show that while the imaging system permits spatial magnification of the field to enable the observation of small features, this magnification comes at the cost of low-pass filtering of the spatial frequency span of the collected coherent scattered light. Next, the scalar model of the nonlinear scattered field is developed to produce the working equations for the contrast signal that is collected by the imaging system. Then, the physical implications of the scattering and image formation models are discussed. Finally, the implications of the image transfer model for holographic optical tomography (HOT) and the tomographic reconstruction algorithm for both second harmonic generation (SHG) and third harmonic generation (THG) are discussed. We conclude by discussing prospects for widefield HOT imaging.

2. Optical holography

Optical detectors respond to the incident optical intensity rather than the field. As a result, all phase information is lost when making any direct optical measurement of a field. However, we know that when the measured light has suitable coherence, we may convert phase differences into intensity modulations through optical interference. Generally, we may have the desired signal field, $U_s$, so that if we have a well-characterized reference field, $U_r$, we may then recover the desired signal field through suitable processing of the interference intensity. This imaging method that recovers the complex wave, i.e., amplitude and phase, is referred to as holography [6].

Holography was first described as a linear scattering model where the scattering object is much smaller than the extent of the incident wave [6]. In this in-line Gabor holography, the unscattered (ballistic) part of the incident wave is treated as the reference wave. The interference between the scattered and unscattered portions of the field constitutes the hologram. In-line holography generated limited excitement initially because the desired scattered field was contaminated by an unwanted conjugate (twin-image) field. This contamination significantly degraded the utility of early holograms. Around the same time that Gabor was working on in-line holography, Leith and co-workers were working in a secret US government program to process synthetic aperture radar films optically. Leith and Upatnieks independently discovered holography, but with a communications theory perspective that employed a spatial frequency carrier for off-axis holography [7, 8, 9]. In off-axis holography, with a suitably large incident reference beam angle, the complete complex signal field may be recovered, thus solving the twin-image problem. In 1997, Yamaguchi demonstrated another elegant solution to the twin-image problem by taking a sequence of holograms where the relative phase of the signal and reference field was shifted, allowing for unique extraction of the complex signal field from a series of in-line holograms [10]. While early holographic work made use of photographic plates, modern holography makes use of digital cameras and numerical processing algorithms [11, 12].

Off-axis holography is hailed as a 3D imaging technique. Indeed, when a hologram produced from an exposed photographic plate is illuminated by a duplicate of the reference wave, one will observe the signal wave as if the object were still present. An observer will see a 3D image of the object. Off-axis holography was first revealed in dramatic fashion with holograms of trains produced in Leith’s laboratory. However, it must be remembered that human visual perception is stereoscopic and not truly 3D. Thus, an observer perceives depth, but does not truly resolve the 3D spatial structure of an object!

Emil Wolf analyzed optical holography from the perspective of linear scattering in the first Born approximation. He developed what is now referred to as the Fourier diffraction theorem, which shows that while a holographic field can be propagated or refocused, the field does not include any axially localized (optically sectioned) information about the object under study [13]. The 3D holography observations can be explained as surface scattering from an object with varied depth so that upon viewing, the observed perceives depth through stereoscopic processing.

The focus of this chapter is on nonlinear holographic imaging. In nonlinear microscopy, signals are recorded from coherent nonlinear scattering, which arises from a distortion of the induced oscillating dipole response of an atom or molecule subjected to a suitably strong illumination (fundamental) field. This nonlinear response can scatter light to new frequencies at harmonics of the incident fundamental field frequency, ${\omega }_m=m{\omega }_1$, where ${\omega }_1$ is the fundamental frequency and $m>1$ is an integer. Both second harmonic generation (SHG, $m=2$) [3] and third harmonic generation (THG, $m=3$) [14, 15, 16] are routinely used for optical microscopy [2]. The first SHG images demonstrated the approach, but routine use had to await the arrival of reliable ultrafast laser oscillators. Today, SHG microscopy is routinely used in biological imaging [4, 5] to look primarily at muscle fibers [17] and collagen [18, 19, 20]. Beyond the study of morphology [21], one may obtain macromolecular structure [22].

Due to the weak strength of the nonlinear optical susceptibility, coherent nonlinear holography had to await the development of more powerful ultrafast sources. Roughly 60 years after the first reports of linear holographic imaging, Demetri Psaltis’ group described SHG holography of SHG-active nanoparticles using a 10 Hz laser amplifier system [23] in a special issue of applied optics that was dedicated to the memory of Emmett Leith [24]. Subsequently, Psaltis demonstrated focusing and imaging of point scatterers in biological tissues and phase conjugation to improve image quality [25, 26, 27, 28]. Shortly after this initial demonstration of SHG holography, imaging in biological systems with oscillators was demonstrated [29, 30, 31, 32, 33, 34]. By optimizing the experimental configuration, quasi-3D imaging at rates of nearly 1500 volumes per second was demonstrated [34]. This early holography work was still limited in the ability to produce detailed 3D imaging and this problem was only recently solved with the introduction of harmonic optical tomography (HOT) [35].

3. Optical diffraction tomography

A measured optical field that has scattered from an object because of a spatial inhomogeneity in either the linear or nonlinear optical susceptibility reveals information about the spatial distribution of the susceptibility. The goal of any coherent imaging system is to uncover quantitative data on the spatial distribution of susceptibility variations, $\mathit{\delta \chi }\left(\mathbf{r}\right)$. The information that is transferred from the incident light field to the scattered field—whether this involves linear or nonlinear scattering—depends on both the properties of the incident light and the nature of the optical physics exploited for the contrast mechanism. In this chapter, we focus our discussion on tomographic imaging obtained through coherent nonlinear optical scattering holographic measurements. To put this technique in context, we will recite the key properties of tomographic imaging under the range of coherent properties of the illumination light and for various contrast mechanisms.

The concept of diffraction tomography was introduced by Emil Wolf in his seminal optics communications paper analyzing optical holography for the case of linear scattering [13]. Wolf’s treatment is reproduced in this chapter and the extension to coherent nonlinear scattering [35] is developed.

In classic holography [7, 8, 9, 11, 12, 36], we consider the illumination of the object by a spatially coherent, monochromatic plane wave with a particular incident propagation wavevector ${\mathbf{k}}_i$. The key observation that follows from the Fourier diffraction theorem is that very limited information is transferred from the scattered field in the process of optical scattering of the incident light from the linear susceptibility perturbation. In a holographic imaging scenario, the transverse spatial frequency span is limited by the numerical aperture of the objective lens, which sets the transverse imaging resolution [8, 37, 38]. However, for the spatially coherent illumination case, the axial frequency support is identically zero. This means that a coherent scattered field recovered from holography does not permit optical sectioning in the imaging process. To fully resolve the object, all the spatial frequency information (or equivalently the spatial information) must be adequately sampled. The methods for fully capturing both the transverse and axial spatial frequency information, called optical diffraction tomography (ODT), make use of either object or beam rotation [37, 39, 40].

The opposite extreme of spatially coherent illumination is the case of spatially incoherent illumination [41, 42]. We may still assume a long coherence time for the case of quasi-monochromatic light. However, the temporally random variations of the illumination field sample the full range of possible incident spatial frequency that is supported by the numerical aperture of the condenser lens. The image transfer model for incoherent illumination exhibits a finite thickness at intermediate spatial frequencies, yet still exhibits no spatial frequency support near zero transverse spatial frequencies [43]. As a result, spatially incoherent imaging also lacks optical sectioning capabilities. However, when an object is placed within the depth of focus of a microscope with spatially incoherent illumination, an absorption tomographic image can be reconstructed [41] with only a slight modification to the standard computed tomography filtered back projection algorithm from a sequence of transillumination intensity images taken over a full rotation of the object [40]. This imaging modality is referred to as optical projection tomography.

Partially coherent illumination fits somewhere in between fully spatially coherent and incoherent imaging. The image transfer model was derived by Streibl for linear scattering under quasi-monochromatic, partially coherent illumination [43]. This transfer function falls between that found for fully coherent and fully incoherent illumination. While in the case of fully spatially coherent light, scattering is accumulated from all depths, the image transfer function of partially coherent light acts as a low-pass filter that rejects the defocused contributions [44]. Streibl showed that by acquiring a stack of images in axial steps of the depth of focus, the set of data can be deconvolved to obtain a three-dimensional image [45, 46, 47]. This same strategy—that of a 3D deconvolution of a stack of images—is a productive approach for 3D imaging of incoherent fluorescent emission [48, 49].

In recent years, Streibl’s approach has been expanded to other partially coherent illumination sources. A method called white light diffraction tomography makes use of spatially coherent light with very broad bandwidth, and thus short temporal coherence [50]. When this very broad bandwidth light is used to illuminate a specimen with a very high NA objective, measurement of the complex field produced by linear scattering through a variant of holography broadens the imaging transverse function axially to permit optical sectioning. Three-dimensional images of the inhomogeneity of the linear susceptibility are then obtained through a 3D deconvolution from a sequence of images taken as the object is displaced in the in the axial direction.

Another strategy that avoids the use of interferometry for extracting the complex field is the use of asymmetric illumination apertures with partially coherent light [51]. When this illumination strategy is coupled with a rigorous model of the imaging transfer function, again a 3D deconvolution can be applied to an axial image stack to obtain a 3D image of the spatial variations in linear optical susceptibility.

The forms of tomography that we have discussed so far are primarily based on optical scattering. As we discuss in detail in the later section, the reliance on scattering with spatially coherent illumination allows for the measured field to relate the input and output scattering directions, which pinpoints the spatial frequency component of the object spatial susceptibility perturbation distribution. However, the short duration of the excited state lifetime and rapid dephasing of fluorescent emitters renders them spatially incoherent. While one might expect this spatial coherence to prevent interference, an individual fluorescent emitter will interfere with itself even though the lack of spatial coherence prevents interference between emitters. As a result, diffractive optical structures can be designed to enable depth-dependent interference intensity structures that allow for holography of incoherent emitters [52, 53]. An alternate strategy can be deployed for mimicking coherent scatting and holographic imaging with incoherent emitters based on the interference of spatially coherent illumination light, either between a plane wave and a point focus [54, 55] or with a pair of plane waves [56, 57, 58, 59, 60, 61]. By using the interference between two spatially coherent illuminating plane waves, one may perform tomographic imaging with fluorescent emitters that exactly mimics ODT [62, 63, 64, 65].

Widefield coherent nonlinear scattering enables the ability to form holograms when a coherent reference beam is directed to interfere with the light produced from the nonlinear scattering process [23, 29]. In the case of illumination with a plane wave fundamental beam, the scattering picture for nonlinear scattering is nearly equivalent to that of linear scattering, with a few modifications. These similarities and differences in the scattering picture will be discussed in later sections of this chapter. The key observation is that we will not obtain any optical sectioning with strictly plane wave illumination. However, due to the weak interaction strength for nonlinear scattering, the fundamental excitation beam is generally weakly focused to provide a balance between field of view activated in the nonlinear scattering process and signal strength that is driven by suitably large field strengths. In such a weak excitation case, moderate 3D imaging resolution is observed [34].

The fact that nonlinear scattering is driven by multiple input fields allows for a completely new form of optical tomography that we called holographic optical tomography (HOT) [35]. In HOT, we employ a high NA condenser to illuminate the object with a broad range of input fundamental spatial frequencies. To ensure widefield illumination, the object is illuminated at defocused plane where the beam is spread out spatially. Because the nonlinear scattering process draws from a broad distribution of illumination spatial frequencies from the full transverse spatial frequency support of the condenser NA, the coherent transfer function for this widefield coherent nonlinear scattering imaging process gains axial spatial frequency support, and thus allows for optical sectioning. Now the strategy first demonstrated by Streibl may be deployed so that 3D tomographic imaging can be obtained from the deconvolution of an axial image stack using a model of the HOT coherent transfer function.

4. Description of optical imaging systems

Optical microscopy can be modeled as a two-dimensional or three-dimensional image collection system. As our focus here is the treatment of tomographic imaging with coherently scattered light, we will provide a discussion of the imaging of spatially coherent light. In the case of coherent nonlinear scattering, however, the weak nonlinear light-matter interaction strength necessitates the use of pulsed light fields. Because the light propagation [66, 67] is linear and shift-invariant, after the coherent nonlinear scatter has occurred (as described in the following section), we may treat the light propagation for each temporal and spatial frequency independently, so that the total field may be obtained from the superposition of the imaged fields. A schematic of the optical imaging system is shown in Figure 1.

Figure 1.

The schematic of the optical 4-F imaging system. We use $U^i$ to denote the incident field, $U^t$ to denote the total field, $U^s$ to denote the field scattered from the object that is imaged as the image field $U^{\mathrm{im}}$, and $U_r$ is the reference field that mixes with the images field to form the hologram.

The spatio-temporal variation of the scalar field is denoted by

$$U_j^{\mathrm{sc}}\left(\mathbf{r},t\right)=a_j\left(t\right)e^{-i{\omega }_{j}t}{u}_j\left(\mathbf{r}\right)$$E1

Here, $j$ is the order of the incident pulse with a complex slowly-varying temporal envelope $a_j\left(t\right)$, that is centered on optical frequency ${\omega }_j=j{\omega }_1$. The object is illuminated by the incident fundamental pulsed field denoted with $j=1$ having a field spectrum centered at ${\omega }_1$. The full spectrum of the complex analytic scalar field is simply $U_j\left(\mathbf{r},\omega \right)={\mathcal{F}}_t\left\{U_j^{\mathrm{sc}}\left(\mathbf{r},t\right)\right\}=A_j\left({\mathrm{\Omega }}_j\right)u_j\left(\mathbf{r}\right)$, where $\mathcal{F}\left\{\cdot \right\}$ is the Fourier transform operator with its subscript $t$ denoting a transform with respect to the time variable $t$ and the relative frequency is ${\mathrm{\Omega }}_j=\omega -j{\omega }_1$. The complex temporal envelope is obtained through the inverse Fourier transform of the field spectral amplitude as $a_j\left(t\right)={\mathcal{F}}_t^{-1}\left\{A_j\left({\mathrm{\Omega }}_j\right)\right\}$. The time, $t$, and optical frequency, variables are conjugate as are the three-dimensional spatial vector, $\mathbf{r}=\left({\mathbf{r}}_{\perp },z\right)$, and the spatial frequencies, $\mathbf{k}=2\pi \left({\mathbf{f}}_{\perp },f_z\right)$. The transverse spatial coordinates are ${\mathbf{r}}_{\perp }=\left(x,y\right)$, with the corresponding transverse spatial frequencies ${\mathbf{f}}_{\perp }=\left(f_x,f_y\right)$.

When imaging a coherent field from an object plane to an imaging plane, we can use a simple shift-invariant model for each optical frequency component. The imaging system will be described by an ideal telecentric 4-F imaging system as shown in Figure 1 [66, 68]. Green’s function for a coherent 4-F imaging system is referred to as the amplitude (or coherent) spread function (CSF), $h\left({\mathbf{r}}_{\perp }\right)$, which is expressed as the spatial convolution

$$u^{\mathrm{im}}\left({\mathbf{r}}_{\perp ,\mathrm{im}},\omega \right)={\int }_{\pm \infty }{u}^{\mathrm{o}}\left({\mathbf{r}}_{\perp ,\mathrm{o}},\omega \right)h\left({\mathbf{r}}_{\perp ,\mathrm{im}}-{\mathbf{r}}_{\perp ,\mathrm{o}},\omega \right)d^2{r}_{\perp ,\mathrm{o}}.$$E2

Here, $u^{\mathrm{im}}=u\left(z=z_{\mathrm{im}}\right)$ is the field at the image plane, $u^{\mathrm{o}}=u\left(z=z_{\mathrm{o}}\right)$ is the field at the object plane, ${\mathbf{r}}_{\perp ,\mathrm{im}}$ is at the image plane and ${\mathbf{r}}_{\perp ,\mathrm{o}}$ is at the object plane. The spatial frequency representation is quite compact and elucidates the low-pass transverse spatial filtering behavior of coherent optical imaging systems through the expression as follows:

$$u^{\mathrm{im}}\left({\mathbf{f}}_{\perp },\omega \right)=H\left({\mathbf{f}}_{\perp },\omega \right)u^{\mathrm{o}}\left({\mathbf{f}}_{\perp },\omega \right).$$E3

The coherent transfer function (CTF) is the Fourier transform of the CSF, $H\left({\mathbf{f}}_{\perp },\omega \right)={\mathcal{F}}_{{\mathbf{r}}_{\perp }}\left\{h\left({\mathbf{r}}_{\perp },\omega \right)\right\}$. Here, we have assumed that the 4-F imaging system has unity magnification for simplicity of notation. The expressions are easily generalized to non-unity magnification [66].

As the propagation of coherent fields through a source-free region can be readily described using the angular spectral propagator, the object field can thus, be propagated from any reference plane to the conjugate object plane of the imaging system. Similarly, the field in the imaging region can also be propagated from one plane to another. The angular spectral propagator for fields propagating in the positive $z$ direction, which we denote as our 4-F imaging system optical axis, is given simply by

$$u\left({\mathbf{f}}_{\perp },z+\mathrm{\Delta }z,\omega \right)=exp\left(i2\mathit{\pi \gamma }\left({\mathbf{f}}_{\perp },\omega \right)\mathrm{\Delta }z\right)u\left({\mathbf{f}}_{\perp },z,\omega \right),z=z^{\mathrm{im}},z^{\mathrm{o}},$$E4

where the axial spatial frequency, ${\gamma }_j\left({\mathbf{f}}_{\perp },{\omega }_j\right)=\sqrt{{\left(n\omega /2\pi c\right)}^2-{‖{\mathbf{f}}_{\perp }‖}^2}$, $c$ is the speed of light in vacuum, and $n\left(\omega \right)\to n$ is the refractive index of the background medium at optical frequency $\omega$. Only non-evanescent spatial frequencies, $\sqrt{{‖\mathbf{f}\perp ‖}^2}<n\omega /2\pi c$, propagate to the far field to be detected.

Within the theory of semiclassical light detection, we may describe the signal recorded by a camera of an incident optical field as the time average of the zero-delay field autocorrelation,

$$I^{\det }\left({\mathbf{r}}_{\perp },\tau \right)={⟨U_j^{\mathrm{s}}\left({\mathbf{r}}_{\perp },z^{\mathrm{im}},t\right)U_j^{\mathrm{r}{,}^{\ast }}\left({\mathbf{r}}_{\perp },z^{\mathrm{im}},t+\tau \right)⟩}_t.$$E5

The angle brackets, $<\cdot >$, denote a time average determined by the detector timescale (e.g., the camera integration time), which for practical nonlinear holographic imaging is orders of magnitude longer than the coherent time of the light fields. For $\tau =0$, this signal can be equivalently represented by the weighted contributions by the cross-spectral density of the light $W^{\mathrm{sr}}\left(\omega \right)=A^{\mathrm{s}}\left(\omega \right)A^{\mathrm{r}{,}^{\ast }}\left(\omega \right)$, leading to the expression

$$I^{\det }\left({\mathbf{r}}_{\perp }\right)=\int W^{\mathrm{sr}}\left(\omega \right)I\left({\mathbf{r}}_{\perp },\omega \right)\mathit{d\omega }.$$E6

The intensity for spatially coherent fields, as we assume here, is defined as:

$$I\left({\mathbf{r}}_{\perp },\omega \right)={\left|u\left({\mathbf{f}}_{\perp },z^{\mathrm{im}},\omega \right)\right|}^2.$$E7

Coherent tomographic imaging requires access to the field directly. This field can be approximately retrieved experimentally through holography that relies on interference with a reference field, which we denote as $u_{\mathrm{r}}=A_{\mathrm{r}}exp\left(i{\varphi }_{\mathrm{r}}\right)$. For simplicity, we have assumed that the reference field is unity amplitude and exhibits a relative phase shift ${\varphi }_{\mathrm{r}}$. This phase shift can vary linearly as in off-axis holography [7, 8, 36] or relative phase shifts can be imparted in a series of measurements as is in phase shifting holography [10]. In either case, the field intensity for the total field given by the sum of the reference field and the images scattered field, $u^{\mathrm{t}}=u^{\mathrm{r}}+u^{\mathrm{im}}$ leads to four terms in the intensity that read $I^{\mathrm{t}}=I^{\mathrm{r}}+I^{\mathrm{im}}+u^{\mathrm{r}{,}^{\ast }}{u}^{\mathrm{im}}+u^{\mathrm{r}}{u}^{\mathrm{im}{,}^{\ast }}$. With suitable numerical processing, we may then isolate the scattered field from the measurement, leading to

$$I^{\text{holo}}\left({\mathbf{r}}_{\perp }\right)=\int W^{\mathrm{sr}}\left(\omega \right)u^{\mathrm{r}{,}^{\ast }}\left({\mathbf{r}}_{\perp },z^{\mathrm{im}},\omega \right)u^{\mathrm{im}}\left({\mathbf{r}}_{\perp },\omega \right)\mathit{d\omega }$$E8

for the case of a unity amplitude reference field. In the transverse spatial frequency domain, we may write this expression as:

$$I^{\text{holo}}\left({\mathbf{r}}_{\perp }\right)=e^{i{\beta }_r{z}^{\mathrm{im}}}\int W^{\mathrm{sr}}\left(\omega \right)u^{\mathrm{im}}\left({\mathbf{r}}_{\perp },\omega \right)\mathit{d\omega }.$$E9

Here, we note that the image of the scattered field, $u_s\left({\mathbf{f}}_{\perp },\omega \right)$, is low-pass filtered by the imaging system CTF.

With a model of coherent imaging of the scattered field, we now need a description of the scattered field to proceed. In the following section, we derive the scattered coherent nonlinear field the the $m^{\mathrm{th}}$order harmonic driven by the fundamental field pulse with an incident fundamental center frequency of ${\omega }_1$.

5. Coherent nonlinear scattering of scalar field

Our goal is to understand the imaging properties, capabilities, and limitations of coherent nonlinear optical holography and tomography. While a full description of coherent nonlinear scattering requires a vector treatment [17, 32, 34], we will restrict our discussion to scalar fields. Such a treatment may provide an understanding of holographic and imaging properties without loss of generality, as the measurement at the camera always involves a projection of the nonlinear scattered field polarization onto the reference field polarization [34]. Thus, we post-select a particular polarization component that can then be regarded as a scalar nonlinear field.

In the scalar description that follows, we begin with the wave equation, where we have made the usual assumptions for optical propagation. Explicitly, these assumptions are that we consider a region devoid of free charges and associated free-charge current densities. Moreover, we assume nonmagnetic media, so the magnetic permeability used is simply that of free space. To simplify the wave equation, we assume that both the linear, ${\chi }^{\left(1\right)}$, and nonlinear, ${\chi }^{\left(m\right)}$, optical susceptibilities are scalar quantities to facilitate a scalar treatment. Finally, we assume that any spatial variation in optical susceptibility is weak compared to the mean optical susceptibility, i.e., ${\mathit{\delta \chi }}^{\left(m\right)}\left(\mathbf{r}\right)\ll {\overline{\chi }}^{\left(m\right)}={⟨{\chi }^{\left(m\right)}\left(\mathbf{r}\right)⟩}_{\mathbf{r}}$, and the angle brackets denote a spatial average. This assumption means that to first order we may treat the medium as spatially homogeneous, which allows for simplification of the wave equation. The inverse scattering problem for imaging the spatial variations of optical susceptibility, ${\mathit{\delta \chi }}^{\left(m\right)}\left(\mathbf{r}\right)$, are treated as a perturbation to the driven homogeneous wave equation.

5.1 Scattering model of scalar field

Making explicit use of the assumptions stated above, we may combine two of Maxwell’s equations to obtain the optical wave equation:

$$\nabla \times \nabla \times \mathcal{U}\left(\mathbf{r},t\right)+{\mu }_0\frac{{\partial }^2\mathcal{D}\left(\mathbf{r},t\right)}{\partial t^2}=0.$$E10

Here, $\mathcal{U}$ is the electric field which is assumed to be scalar. The right-hand side of the equation above is $0$ because the external source of the field is excluded from the interested region and we assume there is no internal source, i.e., the material is not self-luminous. Here, $\mathcal{D}$ is the scalar displacement vector that describes the displacement current, including the linear and nonlinear response from bound charges in the material. The nonlinear contributions to this displacement constitute the quantity that we wish to image. Here, $\mathcal{U}$ and $\mathcal{D}$ denote real fields. Below, we will drive coupled wave equations for the propagation of the fundamental and nonlinear fields, where we assume that these fields may be described as complex analytic functions.

For isotropic media that can, to the first approximation, be treated as spatially homogeneous, in the absence of free charges, Gauss’ law, $\nabla \cdot \mathcal{D}=0$, allows us to make the approximation $\nabla \cdot \mathcal{U}\approx 0$. Making use of this simplification and a standard vector identity, the first term in Eq. (10) simplifies to

$$\nabla \times \nabla \times \mathcal{U}=\nabla \left(\nabla \cdot \mathcal{U}\right)-{\nabla }^2\mathcal{U}=-{\nabla }^2\mathcal{U}.$$E11

Now, Eq. (10) becomes

$${\nabla }^2\mathcal{U}-{\mu }_0\frac{{\partial }^2\mathcal{D}}{\partial t^2}=0.$$E12

Our interest lies in the nonlinear response, which is encapsulated in the real displacement field, which is written as

$$\mathcal{D}\left(\mathbf{r},t\right)={\epsilon }_0\mathcal{U}\left(\mathbf{r},t\right)+\mathcal{P}\left(\mathbf{r},t\right).$$E13

The total real polarization density of the form of a superposition of the linear and the nonlinear response is given by

$$\mathcal{P}\left(\mathbf{r},t\right)={\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},t\right)+{\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},t\right).$$E14

The linear polarization density follows a convolution of the linear optical response of the medium

$${\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},t\right)={\epsilon }_0{\int }_{-\infty }^t{R}^{\left(1\right)}\left(\mathbf{r},\tau \right)\mathcal{U}\left(\mathbf{r},\tau \right)\mathit{d\tau }.$$E15

The linear, causal optical response function of the medium is $R^{\left(1\right)}\left(\mathbf{r},t\right)$ and is related to the frequency-dependent optical susceptibility through the temporal Fourier transform relationship $R^{\left(1\right)}\left(\mathbf{r},t\right)={\mathcal{F}}_{\omega }\left\{{\chi }^{\left(1\right)}\left(\mathbf{r},\omega \right)\right\}$. Here, ${\mathcal{F}}_{\omega }\left\{\cdot \right\}={\mathcal{F}}_t^{-1}\left\{\cdot \right\}={\left(2\pi \right)}^{-1}\int F\left(\omega \right)exp\left(-i\omega t\right)\mathit{d\omega }$ denotes an inverse Fourier transform.

For nonresonant interactions, the nonlinear polarization density may generally be expanded as a power series of the form

$${\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},t\right)={\mathcal{P}}^{\left(2\right)}\left(\mathbf{r},t\right)+{\mathcal{P}}^{\left(3\right)}\left(\mathbf{r},t\right)+\cdots +{\mathcal{P}}^{\left(m\right)}\left(\mathbf{r},t\right)+\cdots .$$E16

Details of the nonlinear polarization density will be deferred to a later section. These nonlinear polarization density terms drive a wide range of nonlinear optical processes. For our purposes, we will focus on $m=2,3$ and only consider the processes that drive second and third harmonic generation (SHG and THG).

Combining all of these expressions, we arrive at our wave equation for coherent nonlinear scattering that reads the equation:

$${\nabla }^2\mathcal{U}\left(\mathbf{r},t\right)-{\mu }_0{\epsilon }_0\frac{{\partial }^2\mathcal{U}\left(\mathbf{r},t\right)}{\partial t^2}-{\mu }_0\frac{{\partial }^2{\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},t\right)}{\partial t^2}-{\mu }_0\frac{{\partial }^2{\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},t\right)}{\partial t^2}=0.$$E17

5.2 Wave equation in the frequency domain

The time-domain equation may easily be represented in the frequency domain by noting that the fields can be represented through an inverse Fourier transform as $\mathcal{U}\left(\mathbf{r},t\right)={\mathcal{F}}_{\omega }\left\{\mathcal{U}\left(\mathbf{r},\omega \right)\right\}$, ${\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},t\right)={\mathcal{F}}_{\omega }\left\{{\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},\omega \right)\right\}$, and ${\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},t\right)={\mathcal{F}}_{\omega }\left\{{\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},\omega \right)\right\}$. The spectra of the real fields are denoted with the argument $\omega$.

Applying the second-order temporal partial derivatives to the inverse Fourier transform in Eq. (17) produces the frequency-domain wave equation:

$${\nabla }^2\mathcal{U}\left(\mathbf{r},\omega \right)+{\omega }^2{\mu }_0{\epsilon }_0\mathcal{U}\left(\mathbf{r},\omega \right)+{\omega }^2{\mu }_0{\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},\omega \right)+{\omega }^2{\mu }_0{\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},\omega \right)=0.$$E18

Making use of the time-domain linear response function in Eq. (15), we may write the equation:

$${\mathcal{P}}^{\mathrm{L}}\left(\mathbf{r},\omega \right)={\epsilon }_0{\chi }^{\left(1\right)}\left(\mathbf{r},\omega \right)\mathcal{U}\left(\mathbf{r},\omega \right).$$E19

Defining the refractive index in the usual way as $n^2\left(\mathbf{r},\omega \right)=1+{\chi }^{\left(1\right)}\left(\mathbf{r},\omega \right)$ and the wavenumber as $\beta \left(\mathbf{r},\omega \right)=\mathit{\omega n}\left(\mathbf{r},\omega \right)/c$ and where the phase velocity of light in a vacuum is $c={\left({\mu }_0{\epsilon }_0\right)}^{-1/2}$, then the wave equation may be written as:

$${\nabla }^2\mathcal{U}\left(\mathbf{r},\omega \right)+{\beta }^2\left(\mathbf{r},\omega \right)\mathcal{U}\left(\mathbf{r},\omega \right)=-{\omega }^2{\mu }_0{\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},\omega \right).$$E20

This equation is now a forced Helmholtz equation, where the LHS describes linear scattering and the RHS is the nonlinear forcing function.

5.3 Slowly varying envelope approximation

The wave equation in Eq. (20) contains the spectrum of the real fields and polarization densities, which includes the complex conjugate of the positive frequencies in the negative frequency region. In addition, these spectra include all optical frequencies, including the fundamental and the nonlinear scattered fields. To simplify these expressions, we assume that each spectral region can be written as a separate spectral envelope, so that we consider, in general, a set of optical pulses (or cw fields) with frequencies centered at ${\omega }_j$ so that we may decompose the total (real) field as the superposition

$$\mathcal{U}\left(\mathbf{r},t\right)=\sum _j{\mathcal{U}}_j\left(\mathbf{r},t\right)$$E21

with the distinct spectral bands centered about ${\omega }_j$. We now write the complex analytic scalar field $U$ for the j${}^{\mathrm{th}}$ frequency band as:

$${\mathcal{U}}_j\left(\mathbf{r},t\right)=\frac{1}{2}U\left(\mathbf{r},t\right)+\frac{1}{2}{U}^{\ast }\left(\mathbf{r},t\right),$$E22

where ${}^{\ast }$ denotes the complex conjugate.

We assume that we have pulses well described by a slowly varying envelope in time relative to the rapid oscillations of a carrier (center) frequency, ${\omega }_j$. Thus, the total field can be written as a slowly varying temporal envelope, $a_j\left(t\right)$, and spatial envelope, $u_j\left(\mathbf{r}\right)$, multiplied by a rapidly varying carrier that is nominally propagating along the direction $z$, giving

$$U_j\left(\mathbf{r},t\right)=a_j\left(t\right)u_j\left(\mathbf{r}\right)exp\left(i{\beta }_{j}z-i{\omega }_{j}t\right).$$E23

The wavenumber at frequency ${\omega }_j$ is defined by ${\beta }_j={\omega }_j{n}_j/c$, where $n_j={⟨n\left(\mathbf{r},{\omega }_j\right)⟩}_{\mathbf{r}}$. The positive frequency (complex analytic field), $U_j\left(\mathbf{r},\omega \right)={\mathcal{F}}_{\omega }\left\{U_j\left(\mathbf{r},t\right)\right\}$, then reads the following equation:

$$U_j\left(\mathbf{r},\omega \right)=A_j\left(\omega -{\omega }_j\right)u_j\left(\mathbf{r}\right)e^{i{\beta }_{j}z}.$$E24

For the convenience of notation, by defining ${\mathrm{\Omega }}_j=\omega -{\omega }_j$, we may write $a_j\left(t\right)={\mathcal{F}}_{\mathrm{\Omega }j}\left\{A_j\left({\mathrm{\Omega }}_j\right)\right\}$. With this definition, we use the standard approach of describing a pulse in terms of a center of mass (of the power spectral density of the light field), ${\omega }_j$, and a slowly varying envelope in time, $a_j\left(t\right)$. We assume that these fields have a temporal envelope that varies slowly with respect to the oscillation of the carrier, ${\omega }_j$. This is the standard slowly varying envelope approximation.

5.4 Harmonic generation

The nonlinear polarization density, e.g., for SHG scattering with $m=2$, (assuming scalar interactions for simplicity) is given by

$${\mathcal{P}}^{\mathrm{NL}}\left(\mathbf{r},t\right)={\epsilon }_0{\chi }^{\left(m\right)}\left(\mathbf{r}\right){\mathcal{U}}^m\left(\mathbf{r},t\right).$$E25

We have explicitly ignored spectral dispersion of the second-order nonlinear coefficient in this expression, and as such we do not need the second-order time response integral nor the second-order response function. Physically, these assumptions equate to assuming that the second-order polarization density responds instantly. Note also that the field and polarization density are described by real quantities in this expression. While many nonlinear interactions can occur, we focus our discussion on coherent nonlinear scattering where we scatter to new frequencies at ${\omega }_m=m{\omega }_1$ due to a nonlinear dipole that oscillates at $m$ times the input fundamental center frequency ${\omega }_1$.

The many interaction terms are considered by taking the m${}^{\mathrm{th}}$ power of the total field that we consider as a superposition of the fundamental, ${\mathcal{U}}_1$, and the m${}^{\mathrm{th}}$ harmonic field, ${\mathcal{U}}_m$. This expansion provides two complex analytic terms leading to a polarization density that generates the m${}^{\mathrm{th}}$ harmonic frequency ${\omega }_m$, given by

$$P^{\left(\mathrm{mHG}\right)}\left(\mathbf{r},t\right)=\frac{1}{2^{m-1}}{\epsilon }_0{\chi }^{\left(m\right)}\left(\mathbf{r}\right)U_1^m\left(\mathbf{r},t\right).$$E26

This term drives coherent nonlinear scattering from the fundamental optical frequency centered at ${\omega }_1$ to the harmonic frequency centered at ${\omega }_n$. We may also consider the complementary process in which the nonlinear field is back-converted to the fundamental through the polarization density term that oscillates at the center frequency ${\omega }_1$

$$P^{\left(\mathrm{bc}\right)}\left(\mathbf{r},t\right)={\epsilon }_0{\chi }^{\left(m\right)}\left(\mathbf{r}\right)U_1^{\ast }\left(\mathbf{r},t\right)U_m\left(\mathbf{r},t\right).$$E27

In an imaging scenario, we may assume that very little coherent linearly scattered power is generated. Thus, we may assume that $U_m\ll U_1$, and that $U_1$ is constant throughout the interaction region. The mathematical approximation for this condition is referred to as the undepleted pump approximation. In this approximation, we may drop the back-conversion term. In addition, we will see that this approximation naturally leads to a nonlinear scattering equation that is homologous to linear scattering in the first Born approximation.

As a specific example, consider SHG, where the time-domain polarization density for the SHG source term reads

$$P^{\left(\mathrm{SHG}\right)}\left(\mathbf{r},t\right)=\frac{1}{2}{\epsilon }_0{\chi }^{\left(2\right)}\left(\mathbf{r}\right)u_1^2\left(\mathbf{r}\right)a_1^2\left(t\right)e^{i\left(2{\beta }_{1}z-2{\omega }_{1}t\right)}.$$E28

The polarization density for SHG oscillates with a center frequency of ${\omega }_2=2{\omega }_1$ and is described by the fundamental pulse spectral autocorrelation that appears from the Fourier transform of the square of the slowly varying fundamental pulse temporal envelope as given by

$$P^{\left(\mathrm{SHG}\right)}\left(\mathbf{r},{\mathrm{\Omega }}_2\right)=\frac{1}{2}{\epsilon }_0{\chi }^{\left(2\right)}\left(\mathbf{r}\right)e^{i2{\beta }_{1}z}{u}_1^2\left(\mathbf{r}\right){\mathcal{F}}_{\omega }\left\{a_1^2\left(t\right)e^{-i2{\omega }_{1}t}\right\}.$$E29

5.5 Coupled wave equations

Making use of the nonlinear scattering assumptions noted in the previous section we are now in a position to write the coupled wave equations for the coherent nonlinear scattering process.

$$\left[{\nabla }^2+{\beta }_{10}^2{n}^2\left(\mathbf{r},\omega \right)\right]U_1\left(\mathbf{r},\omega \right)=0,$$E30

$$\left[{\nabla }^2+{\beta }_{m0}^2{n}^2\left(\mathbf{r},\omega \right)\right]U_m\left(\mathbf{r},\omega \right)=-{\omega }^2{\mu }_0{P}^{\left(\mathrm{mHG}\right)}\left(\mathbf{r},\omega \right),$$E31

and where the nonlinear polarization density reads

$$P^{\left(\mathrm{mHG}\right)}\left(\mathbf{r},\omega \right)=\frac{1}{2^{m-1}}{\epsilon }_0{\chi }^{\left(m\right)}\left(\mathbf{r}\right)u_1^m\left(\mathbf{r},\omega \right)e^{im{\beta }_{1}z}{\mathcal{F}}_{\omega }\left\{a_1^m\left(t\right)e^{-im{\omega }_{1}t}\right\}.$$E32

The first equation describes linear scattering, while the second is the nonlinear scattering at the m${}^{\mathrm{th}}$ harmonic. Critically here, we have used the undepleted pump approximation because the nonlinear scattered field is assumed to ever gain enough strength to drive the back conversion process. In addition, we assume zero input coherent nonlinear field at the input boundary. The free-space wavenumber for the j${}^{\mathrm{th}}$ frequency term is ${\beta }_{j0}={\omega }_j/c$.

By defining the m${}^{\mathrm{th}}$-order autocorrelation function as

$${\mathcal{A}}_m\left({\mathrm{\Omega }}_m\right)={\mathcal{F}}_{\omega }\left\{a_1^m\left(t\right)e^{-im{\omega }_1}\right\},$$E33

we may write the forced equation governing linear and coherent nonlinear scattering as:

$$\left[{\nabla }^2+{\beta }_{10}^2{n}^2\left(\mathbf{r},\omega \right)\right]u_1\left(\mathbf{r},\omega \right){\mathcal{A}}_1\left({\mathrm{\Omega }}_1\right)=0,$$E34

$$\left[{\nabla }^2+{\beta }_{m0}^2{n}^2\left(\mathbf{r},\omega \right)\right]u_m\left(\mathbf{r},\omega \right){\mathcal{A}}_m\left({\mathrm{\Omega }}_m\right)=-\frac{{\beta }_{m0}^2}{2^{m-1}}{\chi }^{\left(m\right)}\left(\mathbf{r}\right)u_1^m\left(\mathbf{r},\omega \right){\mathcal{A}}_m\left({\mathrm{\Omega }}_m\right).$$E35

This form of the equations admits the construction of solutions from the free space Green’s functions. Here, we have assumed that spectral width is sufficiently narrow than the multiplicative ${\omega }^2\approx {\omega }_m^2$.

The equations above allow for a general spectrally-varying treatment of coherent nonlinear holography and tomography. However, the effects of the spectral variation on propagation and on the interpretation of scattering are not strongly dependant on the pulse spectrum. In order to simplify the following interpretation of the imaging transfer function, we will assume that we have a narrow enough spectrum so we may make a continuous wave (cw) approximation, where ${\mathcal{A}}_m\left({\mathrm{\Omega }}_m\right)\approx A_m\delta \left({\mathrm{\Omega }}_m\right)$. Invoking this approximation and integrating over ${\mathrm{\Omega }}_m$ and assuming unity amplitude fields leads to the simplified form of the coupled wave equations given by

$$\left[{\nabla }^2+{\beta }_{10}^2{n}^2\left(\mathbf{r},{\omega }_1\right)\right]u_1\left(\mathbf{r},{\omega }_1\right)=0,$$E36

for the fundamental and

$$\left[{\nabla }^2+{\beta }_{m0}^2{n}^2\left(\mathbf{r},\omega \right)\right]u_m\left(\mathbf{r},\omega \right)=-\frac{{\beta }_{m0}^2}{2^{m-1}}{\chi }^{\left(m\right)}\left(\mathbf{r}\right)u_1^m\left(\mathbf{r},\omega \right)$$E37

for the nonlinear harmonic field (Figure 2).

Figure 2.

The scattering model. The object is illuminated by a fundamental wave, for example, a plane wave with wavevector ${\mathbf{k}}_i$. The field scattered from the object, with a range of detected wavevectors $\mathbf{k}$, is measured at the plane perpendicular to $z$ axis with a size limited by the objective. The origin of the coordinate is denoted by the black dot. For nonlinear scattering, we have ${\mathbf{k}}_i\to {\mathbf{k}}_1^i$, $k\to {\mathbf{k}}_m$.

5.6 Holography with a linear scattering

Inspection of the coupled wave equations in Eqs. (36) and (37) makes it clear that within the undepleted pump approximation the fundamental field solution is independent of the nonlinear scattering. Thus, it is fruitful to first obtain a solution to the linear field propagation, and we will consider the general case where the linear susceptibility varies in space. We may rewrite Eq. (36) in the form of the equation:

$$\left[{\nabla }^2+{\beta }_{10}^2{n}_1^2\right]u_1\left(\mathbf{r}\right)=-{\beta }_{10}^2{\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)u_1\left(\mathbf{r}\right).$$E38

Our goal is to solve for ${\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)={\chi }^{\left(1\right)}\left(\mathbf{r}\right)-{\overline{\chi }}^{\left(1\right)}$ to produce an image of image the susceptibility variation, which constitutes our object. The background linear optical susceptibility, ${\overline{\chi }}^{\left(1\right)}$, is chosen so that ${\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)$ lies in a compact domain, i.e., so that it is contained within some volume $V$. The background refractive index at the fundamental frequency ${\omega }_1$ is then given by $n_1^2=1+{\overline{\chi }}^{\left(1\right)}$.

Solutions to Eq. (38) in the first Born approximation may be constructed using Green’s theorem with the formula as follows:

$$u_1\left(\mathbf{r}\right)=u_1^{\left(i\right)}\left(\mathbf{r}\right)-{\beta }_{10}^2{\int }_{V}g\left(\mathbf{r},{\mathbf{r}}^{\text{'}}\right){\mathit{\delta \chi }}^{\left(1\right)}\left({\mathbf{r}}^{\text{'}}\right)u_1^{\left(i\right)}\left({\mathbf{r}}^{\text{'}}\right)d^3{\mathbf{r}}^{\text{'}}.$$E39

Such solutions in the domain outside of the compactly supported susceptibility perturbation, i.e., $\mathbf{r}\notin V$, can be constructed because we have defined the susceptibility perturbation so that it is contained inside of the volume $V$. The free-space Green’s function is defined by

$${\nabla }^{2}g\left(\mathbf{r},{\mathbf{r}}^{\text{'}}\right)+{\beta }_1^{2}g\left(\mathbf{r},{\mathbf{r}}^{\text{'}}\right)={\delta }^{\left(3\right)}\left(\mathbf{r}-{\mathbf{r}}^{\text{'}}\right),$$E40

where ${\beta }_1=n_1{\beta }_{10}$. The incident fundamental wave, $u_1^{\left(i\right)}\left(\mathbf{r}\right)$, is a solution to the homogeneous wave equation, given by Eq. (38) when the susceptibility term vanishes with ${\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)=0$.

Making use of a three-dimensional spatial frequency decomposition, where ${\mathbf{k}}_j=\left(k_{\mathit{jx}},k_{\mathit{jy}},k_{\mathit{jz}}\right)$, and where the norm of the wavevector gives the wavenumber ${\mathbf{k}}_j\cdot {\mathbf{k}}_j={\beta }_j^2$, Green’s function is written as

$$g\left(\mathbf{r},{\mathbf{r}}^{\text{'}}\right)=\frac{1}{{\left(2\pi \right)}^3}{\int }_{-\infty }^{\infty }g\left({\mathbf{k}}_1\right)e^{i{\mathbf{k}}_1\cdot \left(\mathbf{r}-\mathbf{r}\text{'}\right)}{d}^3{\mathbf{k}}_1.$$E41

Making use of this expansion in Eq. (39) produces the spatial frequency spectrum of the free space Green’s function as:

$$g\left({\mathbf{k}}_1\right)=\frac{1}{\left(k_{1z}-{\mathrm{\Gamma }}_1\right)\left(k_{1z}+{\mathrm{\Gamma }}_1\right)}.$$E42

Here, we have defined the transverse spatial frequency vector ${\mathbf{k}}_{j\perp }=\left(k_{\mathit{jx}},k_{\mathit{jy}}\right)$ as well as ${\mathrm{\Gamma }}_j=\sqrt{{\beta }_j^2-{‖{\mathbf{k}}_{j\perp }‖}^2}$.

Computation of the inverse Fourier transform along the $z$ direction by making use of contour integration leads to the Weyl expression for Green’s function, that is:

$$g\left({\mathbf{r}}_{\perp },z,\mathbf{r}{\text{'}}_{\perp },z^{\prime }\right)=\frac{1}{i4\pi }{\int }_{-\infty }^{\infty }\frac{e^{i{\mathbf{k}}_{1\perp }\cdot \left({\mathbf{r}}_{\perp }-\mathbf{r}{\text{'}}_{\perp }\right)}{e}^{i{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)\mid z-z^{\prime }\mid }}{{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)}{d}^2{\mathbf{k}}_{1\perp }.$$E43

By defining the source term in the RHS of Eq. (38) as $S_1\left(\mathbf{r}\right)={\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)u_1^{\left(i\right)}\left(\mathbf{r}\right)$ in the first Born approximation, we may write the solution for the scattered wave as:

$$u_1\left({\mathbf{k}}_1\right)=u_1^{\left(i\right)}\left({\mathbf{k}}_1\right)-{\beta }_{10}^{2}g\left({\mathbf{k}}_1\right)S_1\left({\mathbf{k}}_1\right).$$E44

Here, the spectrum of the source term reads $S_1\left({\mathbf{k}}_1\right)={\mathit{\delta \chi }}^{\left(1\right)}\left({\mathbf{k}}_1\right)⊛{\mathbf{k}}_1{u}_1^{\left(i\right)}\left({\mathbf{k}}_1\right)$ and the spectrum of the input fundamental wave that is propagating along $+z$ is $u_1^{\left(i\right)}\left({\mathbf{k}}_1\right)=H_i\left({\mathbf{k}}_{1\perp }\right)\delta \left(k_{1z}-{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)\right)$. The operator $⊛$ represents a convolution.

The forward and backward scattered fields are obtained by applying contour integration and selecting the suitable residue from the simple roots of Eq. (42). In the forward direction ($z>0$), we can obtain the field as:

$$u_1^f\left({\mathbf{k}}_{1\perp },z\right)=u_1^{\left(i\right)}\left({\mathbf{k}}_{1\perp },z\right)-i\frac{{\beta }_{10}^2{e}^{i{\mathrm{\Gamma }}_{1}z}}{2{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)}{S}_1\left({\mathbf{k}}_{1\perp },{\mathrm{\Gamma }}_1\right).$$E45

In the backward direction, $z<0$, we obtain

$$u_1^b\left({\mathbf{k}}_{1\perp },z\right)=i\frac{{\beta }_{10}^2{e}^{-i{\mathrm{\Gamma }}_{1}z}}{2{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)}{S}_1\left({\mathbf{k}}_{1\perp },-{\mathrm{\Gamma }}_1\right).$$E46

We will restrict our discussion to the forward-propagating wave collected by the imaging system with a coherent transfer function given by $H_o\left({\mathbf{k}}_{1\perp }\right)$. Thus, the imaged field reads as:

$$u_1^{\mathrm{im}}\left({\mathbf{k}}_{1\perp },z\right)=-i{H}_o\left({\mathbf{k}}_{1\perp }\right)\frac{{\beta }_{10}^2{e}^{i{\mathrm{\Gamma }}_{1}z}}{2{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)}{S}_1\left({\mathbf{k}}_{1\perp },{\mathrm{\Gamma }}_1\right).$$E47

Here, we have dropped the unscattered portion of the incident field to focus on the scattered field and simplify the discussion that follows.

In the process of recording a hologram, we multiply by a reference field, $u_1^{\mathrm{ref}}\left({\mathbf{k}}_{1\perp },z\right)={\delta }^{\left(2\right)}\left({\mathbf{k}}_{1\perp }\right)exp\left(i{\beta }_{1}z\right)$, so that our cw hologram signal, $I_1^{\text{holo}}\left({\mathbf{r}}_{\perp },z\right)=u_1^{\mathrm{ref}{,}^{\ast }}\left({\mathbf{r}}_{\perp },z\right)u_1^{\mathrm{im}}\left({\mathbf{r}}_{\perp },z\right)$, which leads to the transverse spectrum of the hologram is given by the cross-correlation operation

$$I_1^{\text{holo}}\left({\mathbf{k}}_{1\perp },z\right)=u_1^{\mathrm{im}}\left({\mathbf{k}}_{\perp },z\right){\otimes }_{{\mathbf{k}}_{1\perp }}{u}_1^{\mathrm{ref}}\left({\mathbf{k}}_{1\perp },z\right).$$E48

The operator $\otimes$ represents a correlation. We have also assumed that the holographic interference term has been shifted to baseband. After evaluation of the cross-correlation integral, the hologram transverse spectrum now reads

$$I_1^{\text{holo}}\left({\mathbf{k}}_{1\perp },z\right)=u_1^{\mathrm{im}}\left({\mathbf{k}}_{\perp },z\right)e^{-i{\beta }_{1}z}.$$E49

To simplify our analysis of the hologram, we will first consider the special case of the fundamental incident wave as a plane wave incident along the direction ${\mathbf{k}}_1^i$, with amplitude $H_i\left({\mathbf{k}}_{1\perp }^i\right)$, so that $u_1^i=H_i\left({\mathbf{k}}_{1\perp }^i\right)exp\left(i{\mathbf{k}}_1^i\cdot \mathbf{r}\right)$ and the corresponding spectrum reads $u_1^i\left({\mathbf{k}}_1\right)={\left(2\pi \right)}^3{H}_i\left({\mathbf{k}}_{1\perp }^i\right){\delta }^{\left(3\right)}\left({\mathbf{k}}_1-{\mathbf{k}}_1^i\right)$.

For the plane illumination case, we may specifically write out our source term convolution integral as follows:

$$S_1\left({\mathbf{k}}_1\right)=\int {\mathit{\delta \chi }}^{\left(1\right)}\left({\mathbf{k}}_1-{\mathbf{k}}_{1^{\prime }}\right)u_1^i\left({\mathbf{k}}_{1^{\prime }}\right)d^3{\mathbf{k}}_{1^{\prime }}$$E50

as

$$S_1\left({\mathbf{k}}_1\right)=H_i\left({\mathbf{k}}_{1\perp }^i\right)e^{i{\beta }_{1}z}{\mathit{\delta \chi }}^{\left(1\right)}\left({\mathbf{k}}_1-{\mathbf{k}}_1^i\right).$$E51

Now the imaged scattered field hologram for a single incident fundamental plane wave illumination reads as:

$$I_1^{\text{holo}}\left({\mathbf{k}}_{1\perp },z\right)=-i{H}_o\left({\mathbf{k}}_{1\perp }\right)H_i\left({\mathbf{k}}_{1\perp }^i\right)\frac{{\beta }_{10}^2{e}^{i\left({\mathrm{\Gamma }}_1-{\beta }_1\right)z}}{2{\mathrm{\Gamma }}_1\left({\mathbf{k}}_{1\perp }\right)}{\mathit{\delta \chi }}^{\left(1\right)}\left({\mathbf{k}}_{1\perp }-{\mathbf{k}}_{1\perp }^i,k_{1z}-{\mathrm{\Gamma }}_1^i\right).$$E52

Note that ${\mathrm{\Gamma }}_j^i\left({\mathbf{k}}_{j\perp }^i\right)=\sqrt{{\beta }_j^2-{‖{\mathbf{k}}_{j\perp }^i‖}^2}$ and the spatial vector is decomposed as $\mathbf{r}=\left({\mathbf{r}}_{\perp },z\right)$. We see that the 2D Fourier transform of the measured field is related to the 3D Fourier transform of the susceptibility distribution of the object.

Making a coordinate transform into the spatial frequency space of the object by defining the scattering vector $\mathbf{Q}=\left({\mathbf{Q}}_{\perp },Q_z\right)\equiv {\mathbf{k}}_1-{\mathbf{k}}_1^i$, then we may rewrite the hologram transverse spectrum for a single plane wave illumination as:

$$I_1^{\text{holo}}\left({\mathbf{Q}}_{\perp },z\right)=-i{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)H_i\left({\mathbf{k}}_{1\perp }^i\right)\frac{{\beta }_{10}^2{e}^{i\left({\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)-{\beta }_1\right)z}}{2{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)}{\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{Q}\right).$$E53

Now, we may take the Fourier transform along z, giving us

$$I_1^{\text{holo}}\left({\mathbf{Q}}_{\perp },Q_z\right)=-i{\beta }_{10}^2{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)H_i\left({\mathbf{k}}_{1\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_1-{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)\right]}{2{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)}{\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{Q}\right).$$E54

The Dirac delta function has selected the portions of the Ewald sphere that are supported by the illumination and collection optical system transfer functions given as $H_i$ and $H_o$, respectively. The hologram field may now be written as a simple linear shift-invariant model with

$$I_1^{\text{holo}}\left(\mathbf{Q}\right)=H_{\mathrm{lin}}\left(\mathbf{Q}\right){\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{Q}\right).$$E55

By inspection of Eq. (54), we may identify the transfer function for a single illumination plane wave, which is given by:

$$H_{\mathrm{lin}}^i\left(\mathbf{Q}\right)\equiv -i{\beta }_{10}^2{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)H_i\left({\mathbf{k}}_{1\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_1-{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)\right]}{2{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)}.$$E56

When using a non-negligible illumination condenser optic NA, then the super-position of all of the illumination k-vectors must be considered so that we can get:

$$H_{\mathrm{lin}}\left(\mathbf{Q}\right)\equiv -i{\beta }_{10}^2\int H_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)H_i\left({\mathbf{k}}_{1\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_1-{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)\right]}{2{\mathrm{\Gamma }}_1\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{1\perp }^i\right)}{d}^2{\mathbf{k}}_{1\perp }^i.$$E57

Here, we have suppressed the explicit optical frequency dependence. For a short pulsed illumination, we would make use of an effective transfer function weighted by the cross-spectral density of the scattered and reference waves:

$$H_{\mathrm{lin}}\left(\mathbf{Q}\right)=\int H_{\mathrm{lin}}\left(\mathbf{Q}\right)W_{\mathrm{sr}}\left(\omega \right)\mathit{d\omega }.$$E58

Notice that only the terms ${\beta }_1$ and ${\beta }_{10}$ exhibit spectral dependence, provided that we can neglect dispersion of ${\mathit{\delta \chi }}^{\left(1\right)}$.

5.7 Coherent nonlinear scattering, holography, and tomography

The linear scattering case can be easily extended to nonlinear scattering. For m${}^{\mathrm{th}}$-order scattering, we first assume that we have zero-incident harmonic field, so that only the scattered fields appear in the expressions. We assume a mean refractive index at the harmonic of ${\overline{n}}_m^2={⟨n_m^2\left(\mathbf{r}\right)⟩}_{\mathbf{r}}$, so that we may define ${\beta }_m={\overline{n}}_m{\beta }_{m0}$ and ${\beta }_{m0}=m{\omega }_1/c$. Furthermore, we define the magnitude of the phase mismatch at $\mathrm{\Delta }\beta ={\beta }_m-m{\beta }_1=m{\beta }_{10}\left({\overline{n}}_m-{\overline{n}}_1\right)$. This parameter is defined to be positive for material that exhibits normal dispersion. This presents a scattering exclusion zone near the origin.

Consider Eq. (37), which is a driven Helmholtz equation analogous to Eq. (38), but where we make the substitutions ${\beta }_{10}^2\to {\beta }_{m0}^2/2^{m-1}$, ${\mathit{\delta \chi }}^{\left(1\right)}\left(\mathbf{r}\right)\to {\chi }^{\left(m\right)}\left(\mathbf{r}\right)$, and $u_1\left(\mathbf{r}\right)\to u_1^m\left(\mathbf{r}\right)$. Now we obtain a nonlinear scattering version of Eq. (44), but where no incident field is present due to the zero boundary condition identified above and we make the further substitutions ${\mathbf{k}}_{1\perp }\to {\mathbf{k}}_{m\perp }$, ${\mathrm{\Gamma }}_1\to {\mathrm{\Gamma }}_m$, and modify the source term to read

$$S_m\left({\mathbf{k}}_m\right)={\chi }^{\left(m\right)}\left({\mathbf{k}}_m\right){⊛}_{{\mathbf{k}}_m}{u}_i^{\left(m\right)}\left(\mathbf{k}\right).$$E59

where $u_i^{\left(m\right)}\left(\mathbf{k}\right)={\mathcal{F}}_{3\mathrm{D}}\left\{u_1^m\left(\mathbf{r}\right)\right\}$.

Defining a generalized nonlinear scattering vector $\mathbf{Q}=\left({\mathbf{Q}}_{\perp },Q_z\right)\equiv {\mathbf{k}}_m-{\mathbf{k}}_m^i$, where ${\mathbf{k}}_m^i={\sum }_{j=1}^m{\mathbf{k}}_1^{\left(j\right)}$, and explicitly we have ${\mathbf{k}}_{m\perp }^i={\sum }_{j=1}^m{\mathbf{k}}_{1\perp }^{\left(j\right)}$ and ${\mathrm{\Gamma }}_m^i={\sum }_{j=1}^m{\mathrm{\Gamma }}_1^{\left(j\right)}$. Then, by following the arguments in the previous section, we may again obtain a linear shift-imaging imaging model given by

$$I_m^{\text{holo}}\left(\mathbf{Q}\right)=H_{\mathrm{mHG}}\left(\mathbf{Q}\right){\chi }^{\left(m\right)}\left(\mathbf{Q}\right).$$E60

Here, the transfer function for the harmonic holography reads:

$$H_{\mathrm{mHG}}\left(\mathbf{Q}\right)=\int H_{\mathrm{mHG}}^i\left(\mathbf{Q}\right)d^2{\mathbf{k}}_{m\perp }^i.$$E61

And where the integrand is given by:

$$H_{\mathrm{mHG}}^i\left(\mathbf{Q}\right)\equiv -\frac{i}{2^{m-1}}{\beta }_{m0}^2{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{m\perp }^i\right)u_i^{\left(m\right)}\left({\mathbf{k}}_{m\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_m-{\mathrm{\Gamma }}_m\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{m\perp }^i\right)\right]}{2{\mathrm{\Gamma }}_m\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{m\perp }^i\right)}.$$E62

A similar extension to illumination with a short optical pulse can be applied to this transfer function as was applied in the linear scattering case (Figure 3).

Figure 3.

Visualization of the construction of the SHG coherent transfer function (CTF) using Eqs. (64) and (65). Panel (a) and (b) show the fundamental pupil function of the condenser lens with an NA of .3. Panel (c) shows the resultant distribution for $u_i^{\left(2\right)}\left({\mathbf{k}}_2\right)$. Panel (d) shows the pupil function for the objective lens collecting SHG light (NA = .9). Panel (e) shows panel (c) again to graphically illustrate the convolution between (d) and (e) that appears at the SHG CTF in panel (f). The final result in panel (f) is the total contribution of all possible plane wave angles and collected SHG allowed by the optics of the system as shown in Eq. (76).

5.8 Example: second harmonic generation holography

Coherent nonlinear holography offers new possibilities for expanded spatial frequency support due to the effect of noncollinear mixing of fundamental spatial frequencies in the nonlinear mixing process [35]. The key difference lies in the source term, which for SHG reads:

$$S_2\left({\mathbf{k}}_2\right)=\int {\chi }^{\left(2\right)}\left({\mathbf{k}}_2-{\mathbf{k}}_{2^{\prime }}\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2^{\prime }}\right)d^3{\mathbf{k}}_{2^{\prime }}.$$E63

The spectrum of the square of the field is the autoconvolution of the spectrum, $u_2^i\left({\mathbf{k}}_2\right)=u_1^i\left({\mathbf{k}}_1\right){⊛}_{{\mathbf{k}}_1}{u}_1^i\left({\mathbf{k}}_1\right)$, which is given by the integral

$$u_i^{\left(2\right)}\left({\mathbf{k}}_2\right)=\int u_1^i\left({\mathbf{k}}_1\right)u_1^i\left({\mathbf{k}}_2-{\mathbf{k}}_1\right)d^3{\mathbf{k}}_1.$$E64

The autoconvolution of the fundamental field spectrum also appears in the problem of modeling reflectance confocal microscopy [69] and for an illumination objective with a half-opening angle $\alpha$ that is defined by $sin\alpha =\mathrm{NA}/n$ and reads as:

$$u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp },k_{2z}\right)=\frac{4{\beta }_2}{\pi K}{sin}^{-1}\left[\frac{1}{p}\left(1-\frac{2{\beta }_{2}cos\alpha }{\mid k_{2z}\mid }\right)\right],\text{for}\mid k_{2z}\mid \ge 2{\beta }_{2}cos\alpha .$$E65

Here, we use the parameters $p=\left(2{\beta }_2|k_{2\perp }|\right)\sqrt{1-{\left(K/2{\beta }_2\right)}^2}/\left(K|k_{2z}|\right)$ and $K=\sqrt{k_{2\perp }^2+k_{2z}^2}$.

Now, we consider the case where we have two input fundamental plane waves $u_{1^{\prime }}\left(\mathbf{r}\right)=a^{\prime }exp\left(i{\mathbf{k}}_{1^{\prime }}\cdot \mathbf{r}\right)$ and ${u_1}^{"}\left(\mathbf{r}\right)=a^{"}exp\left(i{{\mathbf{k}}_1}^{"}\cdot \mathbf{r}\right)$. Now, the integral in Eq. (64)

$$u_i^{\left(2\right)}\left({\mathbf{k}}_2\right)=a^{\prime }{a}^{"}\int {\delta }^{\left(3\right)}\left({\mathbf{k}}_1-{\mathbf{k}}_{1^{\prime }}\right){\delta }^{\left(3\right)}\left({\mathbf{k}}_2-{\mathbf{k}}_1-{{\mathbf{k}}_1}^{"}\right)d^3{\mathbf{k}}_1,$$E66

which simplifies to

$$u_i^{\left(2\right)}\left({\mathbf{k}}_2\right)=a^{\prime }{a}^{"}{\delta }^{\left(3\right)}\left({\mathbf{k}}_2-{\mathbf{k}}_{1^{\prime }}-{{\mathbf{k}}_1}^{"}\right).$$E67

So that now the source term evaluates to read

$$S_2\left({\mathbf{k}}_2\right)=a^{\prime }{a}^{"}{\chi }^{\left(2\right)}\left({\mathbf{k}}_{2\perp }-{{\mathbf{k}}_{1\perp }}^{\text{'}}-{{\mathbf{k}}_{1\perp }}^{"},k_{2z}-{\mathrm{\Gamma }}_{1^{\prime }}-{{\mathrm{\Gamma }}_1}^{"}\right).$$E68

The imaged SHG field reads

$$u_2^{\mathrm{im}}\left({\mathbf{k}}_{2\perp },z\right)=-i{H}_o\left({\mathbf{k}}_{2\perp }\right)\frac{{\beta }_{20}^2{e}^{i{\mathrm{\Gamma }}_2\left({\mathbf{k}}_{2\perp }\right)z}}{4{\mathrm{\Gamma }}_2\left({\mathbf{k}}_{2\perp }\right)}{S}_2\left({\mathbf{k}}_{2\perp },{\mathrm{\Gamma }}_2\right).$$E69

The SHG hologram reads

$$I_2^{\text{holo}}\left({\mathbf{k}}_{2\perp },z\right)=u_2^{\mathrm{im}}\left({\mathbf{k}}_{2\perp },z\right)e^{-i{\beta }_{2}z}.$$E70

For the case with a pair of fundamental plane waves, the SHG hologram term

$$I_2^{\text{holo}}\left({\mathbf{k}}_{2\perp },z\right)=-i{H}_o\left({\mathbf{k}}_{2\perp }\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp }^i\right)\frac{{\beta }_{20}^2{e}^{i\left({\mathrm{\Gamma }}_2\left({\mathbf{k}}_{2\perp }\right)-{\beta }_2\right)z}}{4{\mathrm{\Gamma }}_2\left({\mathbf{k}}_{2\perp }\right)}{\chi }^{\left(2\right)}\left({\mathbf{k}}_{2\perp }-{\mathbf{k}}_{2\perp }^i,k_{2z}-{\mathrm{\Gamma }}_1^{i\left(2\right)}\right).$$E71

We have defined in incident SHG vector as the sum of the two incident fundamental k-vectors, ${\mathbf{k}}_2^i={\mathbf{k}}_{1^{\prime }}+{{\mathbf{k}}_1}^{"}$, and explicitly we have ${\mathbf{k}}_{2\perp }^i={{\mathbf{k}}_{1\perp }}^{\text{'}}+{{\mathbf{k}}_{1\perp }}^{"}$ and ${\mathrm{\Gamma }}_1^{i\left(2\right)}={{\mathrm{\Gamma }}_1}^{\text{'}}+{{\mathrm{\Gamma }}_1}^{"}$.

Defining a new SHG scattering vector $\mathbf{Q}=\left({\mathbf{Q}}_{\perp },Q_z\right)\equiv {\mathbf{k}}_2-{\mathbf{k}}_2^i$, then we may rewrite the SHG hologram as

$$I_2^{\text{holo}}\left({\mathbf{Q}}_{\perp },z\right)=-i{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp }^i\right)\frac{{\beta }_{20}^2{e}^{i\left({\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)-{\beta }_2\right)z}}{4{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)}{\chi }^{\left(2\right)}\left(\mathbf{Q}\right).$$E72

Now we may take the Fourier transform along z, giving us (Figure 4)

Figure 4.

An illustration of the contributions to the SHG CTF as described in Eq. (76). The green arrows represent the scattered SHG direction, the red arrows show the pair of incident fundamental plane wave direction that generate SHG scattering. The angle of the red arrows is limited by the condenser lens whereas the angle of the green arrow is limited by the objective lens. The color map is the contribution of all fundamental input plane waves given in Eq. (75). The dashed lines show the possible SHG scattering contributions of the detected scattering angle.

$$I_2^{\text{holo}}\left({\mathbf{Q}}_{\perp },z\right)=-i{\beta }_{20}^2{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_2-{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)\right]}{4{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)}{\chi }^{\left(2\right)}\left(\mathbf{Q}\right).$$E73

The SHG hologram field may now be written as a simple linear shift invariant model with

$$I_2^{\text{holo}}\left(\mathbf{Q}\right)=H_{\mathrm{SHG}}\left(\mathbf{Q}\right){\chi }^{\left(2\right)}\left(\mathbf{Q}\right).$$E74

The transfer function for a pair of fundamental illumination plane waves, is given by

$$H_{\mathrm{SHG}}^i\left(\mathbf{Q}\right)=-i{\beta }_{20}^2{H}_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_2-{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)\right]}{4{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)}.$$E75

When using a non-negligible illumination condenser optic NA, then the super-position of all of the illumination k-vectors must be considered, so that

$$H_{\mathrm{SHG}}\left(\mathbf{Q}\right)=-i{\beta }_{20}^2\int H_o\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)u_i^{\left(2\right)}\left({\mathbf{k}}_{2\perp }^i\right)\frac{\delta \left[Q_z+{\beta }_2-{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)\right]}{4{\mathrm{\Gamma }}_2\left({\mathbf{Q}}_{\perp }+{\mathbf{k}}_{2\perp }^i\right)}{d}^2{\mathbf{k}}_{2\perp }^i.$$E76

6. Harmonic optical tomography (HOT) conclusions

The conventional approach to harmonic holography, that is imaging of nonlinear scattering with holographic detection using a coherent harmonic reference beam, cannot provide depth information, known as optical sectioning. While one can rotate the illumination beam or the object, neither of these strategies are very practical. In the case of object rotation, the mechanical positioning errors introduced by the translation and rotation stages make high-resolution imaging all but impossible. While illumination beam scanning is easier, one is left with the classic missing cone problem, and thus estimation of the object is a difficult inverse problem that is prone to distortion.

A few years ago, we introduced a completely new strategy that takes advantage of the fact that coherent nonlinear scattering mixes all possible pairs of incident fundamental plane waves to produce a vast array of scattering directions. The result is that with a suitably large condenser NA for focusing the fundamental light, optical sectioning is admitted into the imaging process. Clearly, point scanning nonlinear scattering takes advantage of this very feature, but in that case, the total power of the scattered harmonic field is detected. As a result, one cannot obtain direct access to the desired nonlinear susceptibility, ${\chi }^{\left(m\right)}\left(\mathbf{r}\right)$. Detection of the field in such a point scanning approach would allow for an identical information transfer from the object to the image as we demonstrate in HOT.

However, HOT is able to exploit cameras, which provides several advantages. First, we have increased speed because we capture a widefield holographic image from which the mHG field is extracted. Second, we benefit from heterodyne amplification of the field because we can bring a strong reference field to detect a weak harmonic field and push to very high imaging speeds [31]. Third, because the CTF exhibits broadening along the direction of propagation ($k_{2z}$), the SHG field image is in focus over a finite depth of field. This means that we may take an image stack by either translating the object axially (along $z$), or by imparting a defocus phase to scan the depth. To produce a 3D tomographic image of ${\chi }^{\left(2\right)}\left(\mathbf{r}\right)$, the extracted field stack is deconvolved with the CTF (see Figure 5). While the low NA example in Figure 5 shows negligible optical sectioning, increases in the condenser NA rapidly expand the axial spatial frequency support. There is ample opportunity to further optimize the resolution, speed, and sensitivity of HOT. In addition, there are opportunities to implement HOT for other coherent nonlinear scattering mechanisms, such as third harmonic generation (THG). Moreover, the polarization behavior can be exploited, and the fact that we detect the fields, rather than the intensity means that the sign, and thus the orientation of the susceptibility tensors are accessible.

Figure 5.

SHG CTF for varying illumination condenser numerical aperture values of 0.1, 0.2 and 0.3 (a), (b) and (c) respectively). The imaging objective numerical aperture is fixed at 0.9.

This chapter focused on the background to introduce the uninitiated to the topics of coherent holography and tomography. We provided a full and rigorous scalar treatment of coherent nonlinear scattering for holographic and nonlinear imaging and tomography. Because most readers will be more familiar with linear scattering, we reviewed linear scattering and then demonstrated the homology to coherent nonlinear scattering through variable substitution to convert from the linear to the nonlinear scattering formulae. The critical difference between linear and nonlinear cases is that the source term in the nonlinear case provides vastly increased spatial frequency support. We demonstrated that this spatial frequency support could be related to a linear shift-invariant imaging model for coherent nonlinear scattering when holographic detection is used. As a result, the entire imaging process can be characterized by a coherent transfer function (e.g., Figure 5). Expressions for computing the CTFs for coherent nonlinear holographic imaging and tomography are derived. We hope that the theory introduced here will inspire new researchers to investigate the use of powerful coherent nonlinear holographic imaging and tomography tools.

Acknowledgments

We gratefully acknowledge funding from the Chan Zuckerberg Initiative’s Deep Tissue Imaging Program.

Conflict of interest

The authors declare no conflict of interest.

Dedication

We dedicate this chapter to the memory of Gabriel (‘Gabi’) Poposecu with whom we collaborated on the original HOT paper. Gabi is a dear friend and colleague and he is sorely missed.