Sampling Without Time: Recovering Echoes of Light via Temporal Phase Retrieval

TL;DR

This paper introduces a novel method for sampling and reconstructing sparse signals without known sampling times, leveraging phase retrieval techniques in optical time-of-flight imaging to recover light echoes from low-pass filtered measurements.

Contribution

It proposes a new approach for temporal phase retrieval in optical sensing, enabling signal reconstruction without prior knowledge of sampling instants or rates.

Findings

Successful experimental validation with optical sensors

Effective phase retrieval from low-pass filtered data

Enhanced imaging capabilities demonstrated

Abstract

This paper considers the problem of sampling and reconstruction of a continuous-time sparse signal without assuming the knowledge of the sampling instants or the sampling rate. This topic has its roots in the problem of recovering multiple echoes of light from its low-pass filtered and auto-correlated, time-domain measurements. Our work is closely related to the topic of sparse phase retrieval and in this context, we discuss the advantage of phase-free measurements. While this problem is ill-posed, cues based on physical constraints allow for its appropriate regularization. We validate our theory with experiments based on customized, optical time-of-flight imaging sensors. What singles out our approach is that our sensing method allows for temporal phase retrieval as opposed to the usual case of spatial phase retrieval. Preliminary experiments and results demonstrate a compelling capability of our phase-retrieval based imaging device.

Figures

Equations

$$m\left(\mathbf{r},t\right)=\left({\lambda }_1{\Gamma }_1\right)\left(\mathbf{r}\right)\delta \left(t-t_1(\mathbf{r})\right)+\left({\lambda }_2{\Gamma }_2\right)\left(\mathbf{r}\right)\delta \left(t-t_2(\mathbf{r})\right)$$

$${\int }_0^{\Delta \gg t_2}m\left(\mathbf{r},t\right)dt=\left({\lambda }_1{\Gamma }_1\right)\left(\mathbf{r}\right)+\left({\lambda }_2{\Gamma }_2\right)\left(\mathbf{r}\right)\equiv I\left(\mathbf{r}\right).$$

$$m\left({\boldsymbol{r},t}\right)=\sum\limits_{k=0}^{K-1}{\left({{\lambda_{k}}{\Gamma_{k}}}\right)\left(\boldsymbol{r}\right)\delta\left({t-{t_{k}{\color[rgb]{0,0,0}(\boldsymbol{r})}}}\right)}$$

$$p\to \boxed{h}\to r\to \boxed{\varphi }\to y\stackrel{\text{Sampling}}{\to }\mathbf{y}\mapsto \underset{\text{Intensity}}{\underbrace{m=\mathcal{C}\left[\mathbf{y}\right]}}$$

$$h_{\text{SI}}\left(t,\tau \right)\equiv h_K\left(t-\tau \right)=\sum _{k=0}^{K-1}{\Gamma }_k\delta \left(t-\tau -t_k\right).$$

$$m=\mathcal{C}\left[\mathbf{y}\right]=\left(\mathcal{C}\left[\phi \right]\ast \mathcal{C}\left[h_K\right]\right),\phi =\mathcal{C}\left[p\right].$$

$$p\left(t\right)\approx\widetilde{p}\left(t\right)=\sum\limits_{\left|\ell\right|\leqslant L}{{{\widehat{p}}_{\ell}}{u_{\ell}}\left({\color[rgb]{0,0,0}\omega_{0}t}\right)},$$

$$m\left(t\right)$$

$${\color[rgb]{0,0,0}m\left(t\right)}=\sum\limits_{|\ell|\leq L}\widehat{\psi}_{\ell}u_{\ell}\left({\color[rgb]{0,0,0}\omega_{0}}t\right)\underbrace{\int{\mathcal{C}\left[h_{K}\right]\left(\tau\right)u_{\ell}^{*}{\left({\color[rgb]{0,0,0}\omega_{0}}\tau\right)}}d\tau}_{\sf{Fourier\ Integral}},$$

$$\hat{s}\left(\omega \right)$$

$$=\underset{a_0}{\underbrace{\sum _{k=0}^{K-1}|{\Gamma }_k{|}^2}}+2\sum _{k=0}^{K-1}\sum _{m=k+1}^{K-1}\underset{a_{k,m}}{\underbrace{\left|{\Gamma }_k\right|\left|{\Gamma }_m\right|}}\cos \left(\omega t_{k,m}+\angle {\Gamma }_{k,m}\right)$$

$$\hat{s}\left(\omega \right)\stackrel{(\text{ref}mag)}{=}$$

$$=$$

$$\left|\left|{\Gamma }_0\right|\pm \left|{\Gamma }_1\right|\right|=\sqrt{a_0\pm a_{01}},a_0,a_{01}>0.$$

$${\left\{|{\tilde{\Gamma }}_k|\right\}}_{k=1,2}=\frac{\sqrt{{\tilde{a}}_0+{\tilde{a}}_{01}}\pm \sqrt{{\tilde{a}}_0-{\tilde{a}}_{01}}}{2}.$$

$$\hat{y}\left(\omega \right)={\left|\hat{p}\left(\omega \right)\right|}^2{h}_K^{\ast }\left(\omega \right),h_K\left(\omega \right)=\sum _{k=0}^{K-1}{\Gamma }_k{u}_{t_k}\left(\omega \right).$$

Full text

1. I Introduction

2. II Image–Formation Model for ToF Sensors

3. III Sampling Echoes of Light

4. III-A Bandlimited Approximation of Probing Function

5. III-B Sampling Theory Context

6. IV Reconstruction via Phase Retrieval

7. V Experimental Validation

8. V-A Simulation

9. V-B Practical Validation

10. VI Conclusions

I Introduction

During a recent visit to meet a collaborator, the author (who happens to be an avid photographer) saw a stark reflection of a local monument on the window panes of the Harvard science center. This is shown in Fig. 1. This is a common place phenomenon at a macro scale as well as a micro scale (for example microscopy). At the heart of this problem is a fundamental limitation, that is, all of the conventional imaging sensors are agnostic to the time information. Alternatively stated, the image–formation process is insensitive to the potentially distinct times that the photons can spend traveling between their sources and the detector.

To elaborate, consider the Gedankenexperiment version of Fig. 1 described in Fig. 2 (a). Let us assume that a light source is co-located with the imaging sensor (such as a digital camera). The reflection from the semi-reflective sheet (such as a window pane) arrives at the sensor at time $t_{1}=2d_{1}/c$ while the mannequin is only observable on or after $t_{2}=2d_{2}/c$. Here, $c$ is the speed of light. While we typically interpret images as a two dimensional, spatial representation of the scene, let us for now, consider the photograph in Fig. 1 along the time-dimension. For the pixel $\left(x_{0},y_{0}\right)$, this time-aware image is shown in Fig. 2 (b). Mathematically, our time-aware image can be written as a $2$–sparse signal (and as a $K$–sparse signal in general),

[TABLE]

where $\left\{{{\Gamma_{k}}\left(\boldsymbol{r}\right)}\right\}$ are the constituent image intensities, $\left\{{{\lambda_{k}}\left(\boldsymbol{r}\right)}\right\}$ are the reflection coefficients, $\boldsymbol{r}=\left(x,y\right)^{\top}$ is the 2D spatial coordinate, and $\delta(\cdot)$ is the Dirac distribution.

A conventional imaging sensor produces images by marginalizing time information, resulting in the 2D photograph,

[TABLE]

Recovering $\left\{{{\Gamma_{k}}}\right\},k=1,2$ given $I$, more generally, $K$ echoes of light,

[TABLE]

is an ill-posed problem. Each year, a number of papers [3, 4, 5, 6, 7, 8] attempt to address this issue by using regularization and/or measurement-acquisition diversity based on image statistics, polarization, shift, motion, color, or scene features. Unlike previous works, here, we ask the question: Can we directly estimate $\left\{{{\Gamma_{k}}}\right\}_{k=0}^{K-1}$ in (1)? In practice, sampling (3) would require exorbitant sampling rates and this is certainly not an option. Also, we are interested in $\left\{{{\Gamma_{k}}}\right\}_{k=0}^{K-1}$ only, as opposed to $\left\{{{\Gamma_{k},t_{k}}}\right\}_{k=0}^{K-1}$ [9] where $t_{k}$ is a non-linear argument in (3). As a result, our goal is to recover the intensities of a sparse signal. For this purpose, we explore the idea of sampling without time—a method for sampling a sparse signal which does not assume the knowledge of sampling instants or the sampling rate. In this context, our contributions are twofold:

1. $1)$

For the general case of $K$-echoes of light, our work relies on estimating the constituent images $\left\{{{\Gamma_{k}}}\right\}_{k=0}^{K-1}$ from the filtered, auto-correlated, time-resolved measurements. This is the distinguishing feature of our approach and is fundamentally different from methods proposed in literature which are solely based on spatio-angular information (cf. [3, 4, 5, 6, 7, 8] and references therein). 2. $2)$

As will be apparent shortly, our work is intimately linked with the problem of (sparse) phase retrieval (cf. [10, 11, 12, 13, 14] and references therein). Our “sampling without time” architecture leads to an interesting measurement system which is based on time-of-flight imaging sensors [1] and suggests that $K^{2}-K+1$ measurements suffice for estimation of $K$ echoes of light. To the best of our knowledge, neither such a measurement device nor bounds have been studied in the context of image source separation [3, 4, 5, 6, 7].

For the sake of simplicity, we will drop the dependence of $m$ and $\Gamma$ on spatial coordinates $\boldsymbol{r}$. This is particularly well suited for our case because we do not rely on any cross-spatial information or priors. Also, scaling factors $\lambda_{k}$ are assumed to be absorbed in $\Gamma_{k}$.

Significance of Phase-free or Amplitude-only measurements: It is worth emphasizing that resolving spikes from a superposition of low-pass filtered echoes, is a problem that frequently occurs in many other disciplines. This is a prototype model used in the study of multi-layered or multi-path models. Some examples include seismic imaging [15], time-delay estimation [16], channel estimation [17], optical tomography [18], ultrasound imaging [19], computational imaging [1, 20, 21, 22, 23] and source localization [24]. Almost all of these variations rely on amplitude and phase information or amplitude and time-delay information. However, recording amplitude-only data can be advantageous due to several reasons. Consider a common-place example based on pulse-echo ranging. Let $p\left(t\right)=\sin\left({{\omega}t}\right)$ be the emitted pulse. Then, the backscattered signal reads $r\left(t\right)=\Gamma_{0}\sin\left({\omega t-{\theta_{0}}}\right)$. In this setting, on-chip estimation of phase ($\theta_{0}$) or time delay ($t_{0}$)[22],

1. $\bullet$

can either be computationally expensive or challenging since $t_{k}$ is a non-linear parameter in (3), and hence, in $r\left(t\right)$. 2. $\bullet$

requires more measurements. For instance, $2$ measurements suffice for amplitude-only estimation ($\Gamma_{0}$) while amplitude-phase pair requires $4$ measurements [2]. This aspect of phase estimation is an important bottleneck as the the frame rate of an imaging device is inversely proportional to number time-domain measurements. 3. $\bullet$

is prone to errors. In many applications, multiple-frequency measurements, $\omega=k{\omega_{0}}$, are acquired [22, 25, 2] assuming that phase and frequency, ${\theta_{0}}=\omega{t_{0}}$ are linearly proportional. However, this is not the case in practice. The usual practice is to oversample [26, 27].

In all such cases, one may adopt our methodology of working with intensity-only measurements.

II Image–Formation Model for ToF Sensors

Optical ToF sensors are active devices that capture 3D scene information. We adopt the generic image–formation model used in [1] which is common to all ToF modalities such as lidar/radar/sonar, optical coherence tomography, terahertz, ultrasound, and seismic imaging. In its full generality, and dropping dependency on the spatial coordinate for convenience, one can first formalize this ToF acquisition model as:

[TABLE]

where,

1. $\rhd$

$p\left(t\right)>0$ is a $T$-periodic probing function which is used to illuminate the scene. This may be a sinusoidal waveform, or even a periodized spline, Gaussian, or Gabor pulse, for instance. 2. $\rhd$

$h\left(t,\tau\right)$ is the scene response function (SRF). This may be a filter, a shift-invariant function $h_{\sf SI}\left(t,\tau\right)=h\left(t-\tau\right)$, or a partial differential equation modeling some physical phenomenon [25]. 3. $\rhd$

$r\left(t\right)=\int{p\left(\tau\right)h\left({t,\tau}\right)d\tau}$ is the reflected signal resulting from the interaction between the probing function and the SRF. 4. $\rhd$

$y\left(t\right)=\int{r\left(\tau\right)\phi\left({t,\tau}\right)d\tau}$ is the continuous-time signal resulting from the interaction between the reflected function and the instrument response function (IRF), or $\phi$, which characterizes the electro-optical transfer function of the sensor. 5. $\rhd$

$\boldsymbol{y}$ is a set of discrete measurements of the form ${\left.{y\left(t\right)}\right|_{t=n\Delta}}$. 6. $\rhd$

{\color[rgb]{0,0,0}\mathcal{C}\left[\boldsymbol{y}\right]=\boldsymbol{y}*\overline{\boldsymbol{y}}} is the cyclic auto-correlation of $\boldsymbol{y}$, where $*$ and $\overline{\left(\cdot\right)}$ denote the convolution and time reversal operators, respectively.

The interplay between $p$, $h$, and $\phi$ results in variations on the theme of ToF imaging [1]. In this work, we will focus on an optical ToF setting. Accordingly:

1. $1)$

The probing function corresponds to a time-localized pulse. 2. $2)$

The SRF, accounting for the echoes of light, is a $K$-sparse filter,

[TABLE] 3. $3)$

The IRF is fixed by design as ${\phi_{{\sf{SI}}}}\left(t,\tau\right)=p\left(t+\tau\right)$. This implements the so-called homodyne, lock-in sensor [1].

Due to this specific shift-invariant structure of the SRF and IRF, we have, $\boldsymbol{y}=p*\overline{p}*{\overline{h}_{K}}\equiv\mathcal{C}\left[p\right]*\overline{h}_{K}$. Finally, the measurements read,

[TABLE]

III Sampling Echoes of Light

III-A Bandlimited Approximation of Probing Function

Due to physical constraints inherent to all electro-optical systems, it is reasonable to approximate the probing signal as a bandlimited function [25, 9]. We use $L$–term Fourier series with basis functions {u_{n}}\left({\color[rgb]{0,0,0}\omega_{0}t}\right)\stackrel{{\scriptstyle\mathrm{def}}}{{=}}{e^{\jmath{\omega_{0}}nt}}, $\omega_{0}T=2\pi$, to approximate $p$ with,

[TABLE]

where the $\widehat{p}_{\ell}$ are the Fourier coefficients. Note that there is no need to compute $p\left(t\right)$ explicitly as we only require the knowledge of $\mathcal{C}\left[p\right]$ (cf. (6)). In Fig. 3, we plot the calibrated $\mathcal{C}\left[p\right]$ and its Fourier coefficients. In this case, $T=232.40$ ns. We also plot its approximation, ${\widetilde{\varphi}}=\mathcal{C}\left[\widetilde{p}\right]$ with $L=25$ together with the measured $\mathcal{C}\left[p\right]$ and $\varphi$.

III-B Sampling Theory Context

The shift-invariant characterization of the ToF image–formation model allows to re-interpret the sampled version of (6) as the filtering of a sparse signal $\mathcal{C}\left[h_{K}\right]$ with a low-pass kernel $\psi=\mathcal{C}\left[\varphi\right]$ (cf. [28, 9]).

[TABLE]

Note that the properties of the auto-correlation operation imply that the sparsity of $\mathcal{C}\left[h_{K}\right]$ is $K^{2}-K+1$, unlike $h_{K}$ that is $K$-sparse and is completely specified by $(K^{2}-K)/2+1$ due to symmetry. Based on the approximation (7) and the properties of convolution and complex exponentials, and defining $\widehat{\psi}_{\ell}=|\widehat{p}_{\ell}|^{4}$, we rewrite $m(t)$ as,

[TABLE]

Finally, the properties of the Fourier transform imply that sampling the above signal $m(t)$ at time instants $t=n\Delta$, results in discrete measurements of the form $\boldsymbol{m}=\mathbf{U}\mathbf{D}_{\hat{\psi}}\hat{\boldsymbol{s}}$, which corresponds to the available acquired data in our acquisition setting. Combining all the above definitions, it follows that:

1. —$\bullet$

${\boldsymbol{m}}\in{\mathbb{R}^{N}}$ is a vector of filtered measurements (cf. (8)). 2. —$\bullet$

${\mathbf{U}}\in{\mathbb{C}^{N\times\left({2{L}+1}\right)}}$ is a DFT matrix with elements {\left[u_{n}\left({\color[rgb]{0,0,0}\omega_{0}}\ell\Delta\right)\right]_{n,\ell}}. 3. —$\bullet$

${{\mathbf{D}}_{\hat{\psi}}}\in{\mathbb{C}^{\left(2L+1\right)\times\left(2L+1\right)}}$ is a diagonal matrix with diagonal elements $\widehat{\psi}_{\ell}$. These are the Fourier-series coefficients of $\mathcal{C}\left[\varphi\right]$. 4. —$\bullet$

$\widehat{\boldsymbol{s}}\in\mathbb{R}^{2L+1}$ is a phase-less vector containing the Fourier coefficients of $\mathcal{C}\left[h_{K}\right]$, which is obtained by sampling the Fourier transform $\widehat{s}\left({{\omega}}\right)$ of $\mathcal{C}\left[h_{K}\right]$ at instants $\omega=\ell\omega_{0}$. The signal $\widehat{s}\left({{\omega}}\right)$ directly depends on the quantities $|\Gamma_{k}|$ of interest to be retrieved and is expressed as

[TABLE]

where ${t_{k,m}}={t_{k}}-{t_{m}}$ and $\angle{\Gamma_{k,m}}=\angle{\Gamma_{k}}-\angle{\Gamma_{m}}$.

It is instructive to note that the relevant unknowns $\left\{{{\Gamma_{k}}}\right\}_{k=0}^{K-1}$ can be estimated from $\widehat{\boldsymbol{s}}\in\mathbb{R}^{2L+1}$ in (10) which in turn depend only on $L$ as opposed to sampling rate and sampling instants. This aptly justifies our philosophy of sampling without time.

IV Reconstruction via Phase Retrieval

Given data model $\boldsymbol{m}=\mathbf{U}\mathbf{D}_{\hat{\psi}}\hat{\boldsymbol{s}}$ (9), we aim to recover the image parameters $|\Gamma_{k}|$ at each sensor pixel. For this purpose, we first estimate samples of $\hat{s}(\omega)$ as $\widehat{\boldsymbol{s}}=\mathbf{D}_{\psi}^{-1}\mathbf{U}^{+}\boldsymbol{m}$, where $\mathbf{U}^{+}$ is the matrix pseudo-inverse of $\mathbf{U}$. This is akin to performing a weighted deconvolution, knowing that $\mathbf{U}$ is a DFT matrix. Next, we solve the problem of estimating $|\Gamma_{k}|$ in two steps (1) First we estimate $a_{k,m}$, and then, (2) based on the estimated values, we resolve ambiguities due to $|\cdot|$.

Parameter Identification via Spectral Estimation: Based on the coefficients $\widehat{\boldsymbol{s}}$, one can then retrieve the amplitude and frequency parameters that are associated with the oscillatory terms as well as the constant value in (10). The oscillatory-term and constant-term parameters correspond to $\{a_{k,m},t_{k,m}\}$ and $a_{0}$, respectively. All parameter values are retrievable from $\widehat{\boldsymbol{s}}$ through spectral estimation [29]; details are provided in [2] for the interested reader.

Note that, given the form of (10) and our acquisition model, the sparsity level of the sequence $\widehat{\boldsymbol{s}}$—corresponding to the total amount of oscillatory and constant terms—is $(K^{2}-K)/2+1$. The magnitude values $a_{k,m}$ and $a_{0}$ can thus be retrieved if at least $L\geq(K^{2}-K)+1$ autocorrelation measurements are performed [30, 2], which is the case in the experiments described in Section V.

Resolving Ambiguities: Based on the aforementioned retrieved parameters, one wishes to then deduce $\{|\Gamma_{k}|\}_{k=0}^{K-1}$. The estimated cross terms $\{a_{k,m},t_{k,m}\}$ allow to retrieve the values of $|\Gamma_{k}|$ through simple pointwise operations. Here we will focus on the case, $K=2$. Due to space limitations, we refer the interested readers to our companion paper [2] for details on the general case ($K>2$). The case when $K=2$,

[TABLE]

may be a result of two distinct echoes (Fig. 2) or due to approximation of higher order echoes by two dominant echoes (due to inverse-square law). In this case, we effectively estimate two terms: $a_{0}$ and $a_{01}$ (10). The set of retrieved magnitude values then amount to solving[2]:

[TABLE]

Thanks to the isoperimetric property for rectangles: $a_{0}\geq a_{01}$, the r.h.s above is always positive unless there is an estimation error in which case, an exchange of variables leads to the unique estimates,

[TABLE]

V Experimental Validation

V-A Simulation

Noting that measurements $\boldsymbol{m}=\mathbf{U}\mathbf{D}_{\hat{\psi}}\hat{\boldsymbol{s}}$ and $\hat{\boldsymbol{s}}$ are linked by an invertible, linear system of equations, knowing $\boldsymbol{m}$ amounts to knowing $\hat{\boldsymbol{s}}$. We have presented detailed comparison of simulation results in [2] for the case where one directly measures $\hat{\boldsymbol{s}}$. In this paper, we will focus on practical setting where the starting point of our algorithm is (9).

V-B Practical Validation

Our experimental setup mimics the setting of Fig. 2. A placard that reads “Time–of–Flight” is covered by a semi-transparent sheet, hence $K=2$. The sampling rate is $\Delta\approx 70\times 10^{-12}$ seconds using a our custom designed ToF imaging sensor [1]. Overall, the goal of this experiment is to recover the magnitudes ${\left\{{|{\Gamma_{k}}|}\right\}_{k=0}^{k=1}}$ given auto-correlated measurements $m\left(n\Delta\right)\stackrel{{\scriptstyle(\ref{acm})}}{{=}}\mathcal{C}\left[y\right]\left(n\Delta\right),\mathbf{m}\in\mathbb{R}^{2795}$. To be able to compare with a “ground truth”, we acquire time-domain measurements $y\left(n\Delta\right)$ before autocorrelation, whose Fourier-domain phases are intact. In Fig. 4(a), we plot the non-autocorrelated measurements $\mathbf{y}$ while phase-less measurements $\mathbf{m}$ are shown in Fig. 4(c) from which we note that the samples are symmetric in time domain due to $m=\mathcal{C}\left[y\right]$ (cf. 6). In the first case (cf. Fig. 4(a)) where $y=\mathcal{C}\left[p\right]*\overline{h}_{K}$, we have [1],

[TABLE]

Similar to $\mathbf{m}$ in (9), in this case, the measurements read $\mathbf{y}=\mathbf{U}\mathbf{D}_{\widehat{\varphi}}\widehat{\mathbf{h}}$ and we can estimate the complex-valued vector $\widehat{\mathbf{h}}$. The phase information in $\widehat{\mathbf{h}}$ allows for the precise computation of $\{\Gamma_{k},t_{k}\}_{k=0}^{K-1}$ [1, 9]. These “intermediate” measurements serve as our ground truth. The spikes corresponding to the SRF (5) are also marked in Fig. 4(a). Fourier-domain measurements $\widehat{h}_{2}^{*}$ linked with $y\left(n\Delta\right)$ are shown in Fig. 4(b). Accordingly, ${\widehat{h}^{*}_{2}}\left(\omega\right),\ \omega=\ell\omega_{0},l=\{-L,\ldots,L\}$ where $\omega_{0}=2\pi/\Delta$ and $L=20$. We estimate the frequency parameters $\{\Gamma_{k},t_{k}\}_{k=0}^{K-1}$ using Cadzow’s method [31]. The resulting fit is plotted in Fig. 4(b).

In parallel to Figs. 4(a),(b), normalized, autocorrelated data $\mathbf{m}=\mathcal{C}\left[\mathbf{y}\right]$ are marked in Fig. 4(c). The signal $\mathbf{y}$ is also shown as reference in green dashed line. We also plot the estimated $\mathcal{C}\left[h_{2}\right]$ (autocorrelated SRF). In Fig. 4(d), we plot measured and deconvolved vector $\widehat{\boldsymbol{s}}=\mathbf{D}_{\psi}^{-1}\mathbf{U}^{+}\boldsymbol{m}$ from which we estimate $a_{0},a_{k,m}$ (10). The result of fitting using Cadzow’s method with $41$ samples is shown in Fig. 4(d).

The reconstructed images, $|\Gamma_{0}|$ and $|\Gamma_{1}|$, due to amplitude-phase measurements (our ground truth, Fig. 4(a)) are shown in Fig. 4(e). Alternatively, the reconstructed images with auto-correlated/intensity-only information are shown in 4(f). One can observe the great similarity between the images obtained in both cases, where only a few outliers appear in the phase-less setting. The PSNR values between the maps reconstructed with and without the phase information are of $30.25$ dB for $|\Gamma_{1}|$ and $48.88$ dB for $|\Gamma_{0}|$. These numerical results indicate that, overall, the phase loss that occurs in our autocorrelated measurements still allows for accurate reconstruction of the field of view of the scene.

Finally, in order to determine the consistency of our reconstruction approach, the PSNR between the original available measurements and their re-synthesized versions—as obtained when reintroducing our reconstructed profiles $|\Gamma_{k}|$ into our forward model—are also provided. As shown in Figs. 4(a) and (c), the PSNR values in the oracle and phase-less settings correspond to $43.68$ dB and $42.41$ dB, respectively, which confirms that our reconstruction approach accurately takes the parameters and structure of the acquisition model into account.

VI Conclusions

In this paper, we have introduced a method to satisfactorily recover the intensities of superimposed echoes of light, using a customized ToF imaging sensor for acquisition and temporal phase retrieval for reconstruction. Up to our knowledge, this is the first method that performs time-stampless sampling of a sparse ToF signal, in the sense that we only measure the amplitudes sampled at uniform instants without caring about the particular sampling times or sampling rates. This innovation can potentially lead to alternative hardware designs and mathematical simplicity as phases need not be estimated in hardware.