Compressive optical interferometry

TL;DR

This paper demonstrates that compressive sensing can be applied directly within optical hardware to efficiently reconstruct sparse optical fields from fewer measurements without additional components.

Contribution

It introduces a method to perform compressive sensing in native optical hardware using random sub-Nyquist sampling of interferograms, avoiding extra components.

Findings

Random sub-Nyquist sampling aids in reconstructing optical fields.

Sensing matrices are provably incoherent and isotropic.

CS can be effectively implemented in optical hardware without added complexity.

Abstract

Compressive sensing (CS) combines data acquisition with compression coding to reduce the number of measurements required to reconstruct a sparse signal. In optics, this usually takes the form of projecting the field onto sequences of random spatial patterns that are selected from an appropriate random ensemble. We show here that CS can be exploited in `native' optics hardware without introducing added components. Specifically, we show that random sub-Nyquist sampling of an interferogram helps reconstruct the field modal structure. The distribution of reduced sensing matrices corresponding to random measurements is provably incoherent and isotropic, which helps us carry out CS successfully.

Equations

$$P(\alpha )\propto \int dx|E(x;\alpha )+E(x;0){|}^2=1+\sum _{n=1}^N|c_n{|}^2\cos (n\alpha ).$$

$$\underbrace{\left[\begin{array}[]{c}P(\alpha_{1})-1\\ P(\alpha_{2})-1\\ \vdots\\ P(\alpha_{M})-1\end{array}\right]}_{\mathbf{y}}\!\!\!=\!\!\!\underbrace{\left[\begin{array}[]{cccc}\cos{\alpha_{1}}&\cos{2\alpha_{1}}&\ldots&\cos{N\alpha_{1}}\\ \cos{\alpha_{2}}&\cos{2\alpha_{2}}&\ldots&\cos{N\alpha_{2}}\\ \vdots&\vdots&\ddots&\vdots\\ \cos{\alpha_{M}}&\cos{2\alpha_{M}}&\ldots&\cos{N\alpha_{M}}\end{array}\right]}_{\mathbf{\Phi}_{\mathrm{int}}}\!\!\!\underbrace{\left[\begin{array}[]{c}|c_{1}|^{2}\\ |c_{2}|^{2}\\ \vdots\\ |c_{N}|^{2}\end{array}\right]}_{\mathbf{x}}.$$

Full text

Compressive sensing (CS) is a strategy for reducing the number of measurements required to reconstruct a signal by projecting it onto a basis of random functions, thereby combining the two steps of data acquisition and compression coding Donoho (2006); Candes (2008). This procedure is particularly effective when the signal is sparse – i.e., the signal receives contributions from a small number of basis functions in some representation Donoho (2006). Underlying this approach is an under-determined linear transformation between the sparse representation and the measurement basis called a sensing matrix that is amenable to a stable pseudoinversion Candes (2006); Donoho (2006). CS has found many uses in optics such as the ‘single-pixel’ camera Duarte et al. (2008); Magalhaes et al. (2011), holography Martinez-Leon et al. (2017); Clemente et al. (2013); Deepan et al. (2014); Rivenson et al. (2013), optical imaging and tomography Liutkus et al. (2014); Graff and Sidky (2015); Süzen et al. (2010); Yao et al. (2015), optical encryption Li et al. (2015); Rawat et al. (2015), spectral analysis Wang and Yu (2014), and efficient quantum state tomography Howland et al. (2014); Mirhosseini et al. (2014); Kalev et al. (2015); Howland et al. (2016). In all these cases, the native optical system is modified by incorporating a new component that projects the field onto random patterns. CS then reduces the number of data samples acquired by the modified system that are needed to recover the target sparse vector.

In this paper, we develop a compressive approach to interferometry for optical modal analysis in arbitrary bases and realize it experimentally. We underscore four major differences from related work. First, unlike much of the prior work which introduces random masks along the optical path to obtain randomized measurements Howland et al. (2016); Martinez-Leon et al. (2017), here we exploit CS in the native optics hardware without modification to the underlying interferometer structure nor addition of new hardware components. Second, this is the first work to realize compression in the generalized interferometry framework introduced recently by some of the authors to carry out modal analysis in arbitrary domains Abouraddy et al. (2012); Martin et al. (2017). Third, we achieve compression gains both in the number of measurements used for reconstruction and the acquisition time. This is in contrast to prior work Mohan et al. (2010) employing CS in optical interferometry, which solely focused on reducing the number of measurements. Fourth, we establish analytical guarantees for successful modal reconstruction under sensing constraints set by the limited degrees of freedom of an optical interferometer.

We introduce CS into a prototypical optical system – a two-path interferometer – and demonstrate experimentally that sampling the interferogram below the Nyquist rate still enables spectral/modal recovery. We first cast the procedure of interferometric spectral/modal reconstruction as a linear measurement problem. We find that the reduced sensing matrix associated with random sub-Nyquist sampling satisfies sub-optimal conditions of isotropy and incoherence, thereby enabling compressive reconstruction from fewer samples Mardani et al. (2015). We denote this overall optical measurement scheme compressive interferometry since it combines reduced data acquisition (reduced number of interferometric points sampled) with compressive coding (replacing harmonic analysis with CS reconstruction).

For sake of generality, we employ a two-path interferometer in which the usual temporal delay is replaced by a ‘generalized delay’ that operates in the Hilbert space defined by any modal basis of interest Abouraddy et al. (2011, 2012). Such a delay corresponds to an optical transformation that has this modal basis as eigenfunctions, which reduces to the traditional temporal delay when spectral harmonics are of interest. The ‘delays’ for Hermite-Gaussian (HG) or radial Laguerre-Gaussian (LG) modal bases are optical realizations of the fractional Fourier transform (FrFT) Namias (1980a); Ozaktas et al. (2001) or fractional Hankel transform (FrHT) Namias (1980b), respectively. In all cases, the compressive interferometry strategy is to replace harmonic analysis of the Nyquist-sampled interferogram with random sub-Nyquist sampling without modifying the native optical hardware. By constructing an intrinsically stable, basis-neutral interferometer that includes realizations of fractional transforms Martin et al. (2017), we confirm that CS helps successfully reconstruct the modal weights of sparse optical beams in a basis of HG or LG modes from the random sub-Nyquist samples. This may help in delay-intolerant applications that require real-time acquisition and processing, such as real-time imaging and high-speed communications exploiting spatial multiplexing.

A typical optical implementation of CS is depicted in Fig. 1(a) Duarte et al. (2008); Howland et al. (2014); Mirhosseini et al. (2014). Consider an optical field $E(x)=\sum_{n}c_{n}\psi_{n}(x)$ represented in an $N$-dimensional orthonormal basis $\{\psi_{n}(x)\}$, and $E(x)$ is normalized such $\sum|c_{n}|^{2}=1$. Here $x$ may refer to one or more spatial dimensions, or even time. The $N\times 1$-vector $\mathbf{x}=\{c_{n}\}$ of coefficients is $s$-sparse if it receives contributions from $\leq s$ elements. In lieu of the $M\geq 2N$ Nyquist-rate measurements normally required to reconstruct $\mathbf{x}$, CS leverages sparsity to reduce $M$ by projecting the field onto a sequence of random masks selected from a Gaussian ensemble to produce a new $M\times 1$ measurement vector $\mathbf{y}$. Effectively, the field has undergone a linear transformation represented by a $M\times N$ sensing matrix $\mathbf{\Phi}$, $\mathbf{y}=\mathbf{\Phi}\mathbf{x}$ [Fig. 1(b)]. CS exploits this under-determined linear system to compressively recover $\mathbf{x}$ from $M\approx\mathcal{O}(s\log{N})\ll N$ linear measurements of $\mathbf{y}$ – provided some conditions on $\mathbf{\Phi}$ are satisfied. For example, Basis Pursuit (BP) recovers $\mathbf{x}$ if $\mathbf{\Phi}$ satisfies the Restricted Isometry Property (RIP) Candes (2008), which requires that $(1-\delta)\|\mathbf{x}\|^{2}\leq\|\mathbf{\Phi}\mathbf{x}\|^{2}\leq(1+\delta)\|\mathbf{x}\|^{2}$, for any $s$-sparse $\mathbf{x}$, where $0<\delta<1$ is a constant. Sensing matrices from Gaussian ensembles, for example, satisfy the RIP with overwhelming probability Candes (2008). The random masks introduce new components into the optical system, thus adding to its complexity while reducing the number of measurements acquired.

The compressive interferometry scheme relies on a different strategy. Consider an interferogram $P(\alpha)$ traced by scanning a ‘delay’ $\alpha$ placed in one arm of a two-path interferometer [Fig. 1(c)]. With respect to the modal basis $\{\psi_{n}(x)\}$, such a delay is an optical transformation having a continuous real control parameter $\alpha$, $\Lambda(x,x^{\prime};\alpha)=\sum_{n}e^{in\alpha}\psi_{n}(x)\psi_{n}^{*}(x^{\prime})$, which is a fractional transform having the basis functions as eigenstates with eigenvalues $e^{in\alpha}$, and $\Lambda(x,x^{\prime};0)=\delta(x-x^{\prime})$ Abouraddy et al. (2012). The parameter $\alpha$ acts as a ‘delay’ in the Hilbert space spanned by $\{\psi_{n}(x)\}$ just as a delay $\tau$ does in the time domain. That $\alpha$ is in fact a delay can be readily seen by superposing the delayed field $E(x;\alpha)=\int dx\Lambda(x,x^{\prime};\alpha)E(x^{\prime})=\sum_{n}c_{n}e^{in\alpha}\psi_{n}(x)$ and the reference $E(x;0)$ to produce the interferogram

[TABLE]

The modal coefficients $\mathbf{x}=\{|c_{n}|^{2}\}$ are revealed by taking a Fourier transform (FT) with respect to $\alpha$. This requires sampling $\alpha$ at the Nyquist rate dictated by the highest-order mode $N$ – even if only a few modes contribute significantly [Fig.1 (d)].

To show that CS algorithms may recover the sparse vector $\mathbf{x}$ from $M\ll N$ measurements, we first cast $P(\alpha)$ from Eq. 1 with selected values of $\alpha$ into a linear model. A $M\times N$ interferometric sensing matrix $\mathbf{\Phi}_{\mathrm{int}}$ relates the $N\times 1$ vector $\mathbf{x}$ of contributions from each Hilbert-space basis-element to the $M\times 1$ measurement-vector $\mathbf{y}$ that contains the sampled interferogram [Fig. 1(e)],

[TABLE]

Here $\mathbf{\Phi}_{\mathrm{int}}$ corresponds to any interferometric measurement, temporal or otherwise. This matrix differs fundamentally from those typically employed in CS applications in optics [Fig. 1(a) and 1(b)] where controllably random transformations $\mathbf{\Phi}$ are realized by an array of random patterns that are judiciously selected such that $\mathbf{\Phi}$ has the RIP, for example. In interferometry, this freedom in designing $\mathbf{\Phi}_{\mathrm{int}}$ is not available because $\mathbf{\Phi}_{\mathrm{int}}$ has a fixed structure with few controllable degrees of freedom (the values of the delay parameter $\alpha$) imposed by the constrained sensing structure (the interferometer itself).

Given the constrained form of $\mathbf{\Phi}_{\mathrm{int}}$ in Eq. 2, it can be shown that selecting random values of $\alpha$ from a uniform distribution at a sub-Nyquist rate yields $\mathbf{\Phi}_{\mathrm{int}}$ that satisfies the RIP with high probability when $M$ and $N$ are sufficiently large. For smaller values of $M$ and $N$, however, the RIP constant $\delta$ is large. The sufficient conditions for all existing reconstruction algorithms require smaller values of $\delta$ to guarantee successful reconstruction. For example, a sufficient condition for BP requires $\delta<\sqrt{2}-1$. Having the large $\delta$ associated with $\mathbf{\Phi}_{\mathrm{int}}$ therefore does not guarantee reconstruction on the basis of these requirements. As an illustration, in Fig. 2(a) we plot $\eta(\mathbf{x})=(\|(\sqrt{2/M})\mathbf{\Phi}_{\mathrm{int}}\mathbf{x}\|^{2}/\|\mathbf{x}\|^{2})-1$ with a typical realization of the sensing matrix for a subset of an ensemble of $10^{6}$ sparse-vector realizations $\mathbf{x}$ with $s=4$, where $\max_{\mathbf{x}}|\eta(\mathbf{x})|\leq\delta$. A histogram constructed from the full ensemble is plotted in Fig. 2(b). However, the sensing matrix $\mathbf{\Phi}_{\mathrm{int}}$ with randomly selected $\alpha$ indeed belongs to a random ensemble satisfying weaker incoherence and isotropy properties despite its limited degrees of freedom imposed by the interferometric process. Provably, such matrices also yield exact reconstruction with high probability via BP Candes and Plan (2011).

To clarify, let $\bm{\phi}$ denote a row of the sensing matrix selected from a random ensemble $\cal F$. The incoherence parameter $\mu(\cal{F})$ is defined as the smallest number for which $\mathrm{max}_{n}|\phi_{n}|^{2}\leq\mu({\cal F})$, $n=1,2,\ldots,N$, where $\phi_{n}$ are the entries of $\bm{\phi}$. The distribution $\cal{F}$ is said to obey the isotropy property if $\mathrm{E}[\bm{\phi}^{{\dagger}}\bm{\phi}]=a\mathbf{I}$, where $\mathrm{E}[.]$ is the expected value, $\mathbf{I}$ the identity matrix, $a$ some constant, and $\bm{\phi}^{{\dagger}}$ the conjugate transpose of $\bm{\phi}$. $\mathbf{\Phi}_{\mathrm{int}}$ satisfies the conditions of incoherence and isotropy for $\alpha$ randomly selected from a uniform distribution over $[0,2\pi]$ with $\mu=1$, as confirmed in Fig. 2(c) that shows $\mathbf{\Phi}_{\mathrm{int}}$ consisting of $M$ rows selected from the full matrix in Fig. 1(d). The matrix appears random in spite of being sampled from a structured sensing system. The matrix $\mathbf{\Phi}_{\mathrm{int}}$ is shown to satisfy the isotropy property in Fig. 2(d).

We now move on to implementing the compressive interferometric scheme experimentally. The generalized delay transformation $\Lambda(x,x^{\prime};\alpha)$ has the modal basis $\{\psi_{n}(x)\}$ as eigenfunctions and its order parameter $\alpha$ represents the delay in the Hilbert space spanned by this basis. We carry out our experiments using two modal bases, 1D HG functions and radial LG functions. The delay transformations $\Lambda(x,x^{\prime};\alpha)$ in these two modal bases correspond to the FrFT Namias (1980a); Ozaktas et al. (2001) and the FrHT Namias (1980b), respectively Abouraddy et al. (2011, 2012). The usual temporal delay is replaced with an optical realization of the appropriate fractional transform, whose order $\alpha$ is swept in the interval $[0,2\pi]$ to produce an interferogram $P(\alpha)$. We construct these fractional transforms using spatial light modulators (SLMs) that realize generalized cylindrical or spherical lenses to produce the FrFT or FrHT, respectively, of desired order Rodrigo et al. (2009). We then exploit the polarization-selectivity of liquid-crystal-based SLMs Kagalwala et al. (2013) to create an intrinsically stable common-path interferometer. One polarization – that undergoes the transformation implemented by the SLMs – is subjected to the generalized ‘delay’ while the orthogonal polarization component is left invariant after the SLMs and thus corresponds to a reference Martin et al. (2017). A polarization projection at the output can then superpose the reference and the delayed beams, thereby yielding the target interferogram $P(\alpha)$ as the delay parameter $\alpha$ is swept evenly at the Nyquist rate.

In the compressive interferometry approach, $\alpha$ is randomly sampled at a sub-Nyquist rate. We study the performance of this strategy by comparing the modal reconstruction in the HG and radial LG bases to that obtained via the FT. $M$ measurements are collected randomly by selecting generalized phases $\alpha_{j}$, $j=1,2,\ldots,M$ from a uniform random distribution $[0,2\pi]$, and we use the BP algorithm to reconstruct the vector $\tilde{\mathbf{x}}$ of modal coefficients. The reconstructed vector is then compared to the vector $\mathbf{x}_{\mathrm{FT}}$ obtained from the FT of $128$ uniformly sampled interferometric measurements by computing the scaled error metric, $||\mathbf{x}_{\mathrm{FT}}-\tilde{\mathbf{x}}||^{2}/||\mathbf{x}_{\mathrm{FT}}||^{2}$. For each example, the average reconstruction error is evaluated by averaging over $100$ runs. The number of potential modes is $N=64$. Figures 3(a) and 3(b) show the reconstructed modal coefficients using FT and CS for input beams formed from HG or LG modes. In Fig. 3(c), we depict the reconstructed modal coefficients for two different superpositions of HG modes. The performance of the BP algorithm is comparable to that of FT using $M=30$ compressive measurements ($\approx 25\%$ of the measurements in FT), leading to substantial savings in acquisition time without any added complexity to the system. The reconstruction error is plotted in Fig. 4 as a function of $M$. The average error approaches zero when $M\approx 25$. Statistical analysis of the observed error reveals that the required number of measurements lies between $M=15$ to $M=30$ to attain nearly error-free reconstruction.

In conclusion, we have devised a compressive approach to interferometry whereby the modal content of light beams can be recovered by sampling intereferograms randomly at sub-Nyquist rate using native optics hardware, which was demonstrated experimentally. Two fundamental differences distinguish our strategy from prior work employing CS in optics. First, the interferometric sensing matrix is imposed through the structure of the interferometer. As such, compression has to be carried out under sensing constraints enforced by the limited degrees of freedom of the sensing system. Second, previous work in the context of interferometry focused on reducing the number of measurements used for reconstruction, but not on compressive data acquisition. Leveraging the notion of generalized delay, our approach enables compressive modal analysis in arbitrary modal basis. Furthermore, we established provable guarantees for successful reconstruction despite the sensing constraints set by the limited degrees of freedom of the interferometer.

Funding. Office of Naval Research (ONR) contracts N00014-14-1-0260 and N00014-17-1-2458.