Here, we describe a compact Single-Pixel Imaging LIDAR System that is based on a compressive sensing technique. The application of the compressive codes to a DMD array enables compression of the spatial information, while the collection of timing histograms correlated to the pulsed laser source ensures image reconstruction at the ranged distances.

Single-pixel cameras have been compared with raster scanning and array based counterparts in terms of noise performance, and proved to be superior. Since a single photodetector is used, a better SNR and higher reliability is expected in contrast with systems using large format photodetector arrays. Furthermore, the event of failure of one or more micromirror elements in the DMD does not prevent full reconstruction of the images. This brings additional robustness to the proposed 3D imaging LIDAR.

The prototype that was implemented has three modes of operation. Range Finder: outputs the average distance between the system and the area of the target under illumination; Attitude Meter: provides the slope of the target surface based on distance measurements in three areas of the target; 3D Imager: produces 3D ranged images of the target surface. The implemented prototype demonstrated a frame rate of 30 mHz for 16x16 pixels images, a transversal (xy) resolution of 2 cm at 10 m for images with 64x64 pixels and the range (z) resolution proved to be better than 1 cm. The experimental results obtained for the “3D imaging” mode of operation demonstrated that it was possible to reconstruct spherical smooth surfaces.

The proposed solution demonstrates a great potential for: miniaturization; increase spatial resolution without using large format detector arrays; eliminate the need for scanning mechanisms; implementing simple and robust configurations.

## I.

## INTRODUCTION

It is well-known that the Nyquist-Shannon sampling theorem has been a fundamental rule of signal processing for many years and can be found in nearly all signal acquisition protocols, being extensively used from consumer video and audio electronics to medical imaging devices or communication systems. Basically, it states that a band-limited input signal can be recovered without distortion if it is sampled at a rate of at least twice the highest frequency component of interest within the signal. For some signals, such as images that are not naturally band limited, the sampling rate is dictated not by the Nyquist-Shannon theorem but by the desired temporal or spatial resolution. However, it is common in such systems to use an anti-aliasing low-pass filter to band limit the signal before sampling it, and so the Nyquist Shannon theorem plays an implicit role [1]. In the last few years, an alternative theory has emerged, showing that super-resolved signals and images can be reconstructed from far fewer data or measurements than what is usually considered necessary. This is the main concept of compressive sensing (CS), also known as compressed sensing, compressive sampling and sparse sampling.

CS relies on the empirical observation that many types of signals or images can be well approximated by a sparse expansion in terms of a suitable basis, that is, by only a small number of non-zero coefficients. This is the key aspect of many lossy compression techniques such as JPEG and MP3, where compression is achieved by simply storing only the largest basis coefficients.

In CS, since the number of samples taken is smaller than the number of coefficients in the full image or signal, converting the information back to the intended domain would involve solving an underdetermined matrix equation. Thus, there would be a huge number of candidate solutions and, as a result, we must find a strategy to select the sparsest solution.

Different approaches to recover information from incomplete data sets have existed for several decades. However, recently, that the field has gained increasing attention, when Emmanuel J. Candès, Justin Romberg and Terence Tao [2], discovered that it was possible to reconstruct Magnetic Resonance Imaging (MRI) data from what appeared to be highly incomplete data sets in face of the Nyquist-Shannon criterion. Following, Candès et al. work, this decoding or reconstruction problem can be seen as an optimization problem and be efficiently solved using the *l*_{1}-norm [3].

As a result, CS has become a kind of revolutionary research topic that draws from diverse fields, such as mathematics, engineering, signal processing, probability and statistics, convex optimization, random matrix theory and computer science.

This paper is organized as follows. After a brief introduction (section I), some mathematical background essential to the understanding of CS is shown (section II). Then, on section III, the innovative single-pixel imaging LIDAR system based on compressive sensing is presented along with some experimental results. In the end, the main conclusions of this work are exposed and relevant future prospects are outlined.

## II.

## COMPRESSIVE SENSING BACKGROUND

In order to become possible, CS is built upon two principles: sparsity, related with the signals of interest, and incoherence, related with the sensing modality.

## A.

### K-sparse and compressible signals

Let’s consider a real-valued, finite-length, one dimensional, discrete-time signal *x*, which can be viewed as a *N*×1 column vector in ℜ^{N} with elements *x*[*n*], with *n* = 1,2,…,N. Any signal in ℜ^{N} can be represented in terms of a basis of *N*×1 vectors For simplicity, let’s assume that the basis is orthonormal. Using the *N×N* basis matrix Ψ = {*ψ*_{1},*ψ*_{2},…,*ψ*_{N}} with the vectors {*ψ _{i}*} as columns, a signal

*x*can be expressed as:

where *s* is the *N*×1 column vector of weighting coefficients *s* and *x* are equivalent representations of the signal with *x* in time or space domain and *s* in **Ψ** domain.

The signal *x* is *K*-sparse if it is a linear combination of only *K* basis vectors, which means that only *K* of the *s _{i}* coefficients in (1) are nonzero, while the remaining (

*N–K*) coefficients are zero. In addition, the signal

*x*is compressible if the representation in (1) has just a few large coefficients and many small coefficients, setting the basis of transform coding. Therefore, we can say that a signal

*x*is sparse in the ψ domain if the coefficient sequence is supported on a small set, and compressible if the sequence is concentrated near a small set. In face of the typical data acquisition paradigm, huge amounts of data are collected only to be in large part discarded at the compression stage to facilitate storage and transmission. Imagine, for example, a digital camera that captures images with millions of sensors (pixels) but eventually encodes the image in just a few hundred kilobytes. Clearly, this is a tremendously wasteful process and suffers from three principal drawbacks. First, the initial number of samples

*N*may be large, even if the desired

*K*is small. Second, the set of all

*N*transform coefficients {

*s*} must be computed even though all but

_{i}*K*of them will be discarded. Third, there is an overhead that is introduced by the encoding of the large coefficients locations.

## B.

### Recovering K-sparse signals

Following the work presented in [4], Candès and Tao developed a refined version of the *Uniform Uncertainty Principle* (UUP) [5], which has proved to be essential to the study of the general robustness of CS. This key notion was then named *Restricted Isometry Property* (RIP) and can be defined as follows:

For each integer *K* = 1,2,…, define the isometry constant *δ _{K}* of a measurement matrix

*A*as the smallest number such that

holds for all *K*-sparse vectors *x*. Therefore, we can say that a matrix *A* obeys the RIP of order *K* if *δ _{K}* differs enough from one. When this condition is verified,

*A*approximately preserves the Euclidean length of

*K*-sparse signals, which in turn implies that

*K*-sparse vectors cannot be in the null space of

*A*. An alternative description of this property is to say that all subsets of

*K*columns taken from

*A*are in fact nearly orthogonal (they cannot be exactly orthogonal since we have more columns than rows).

Let’s imagine we want to acquire *K*-sparse signals making use of matrix *A*. Suppose that *δ _{2K}* is not close to one. This indicates that all pair-wise distances between

*K*-sparse signals must be well preserved in the measurement space, which means that

## C.

### Incoherence

Let’s now consider *M < N* linear measurements of *x* and a collection of test functions such that *y*[*m*] = ⟨*x,φ _{m}*⟩. By stacking the measurements

*y[m]*into the

*M*×1 vector

*y*and the test functions as rows into an

*M*×

*N*sensing matrix Φ we can write

A condition related with RIP is *incoherence*, which requires that the rows of Φ (the measurement or sensing matrix) cannot represent the columns of Ψ in a sparse way (and vice-versa).

Incoherence extends the duality between time and frequency and expresses the idea that an object having a sparse representation in Ψ must be spread out in the domain in which it was acquired. This incoherence property is verified for many pairs of bases, including, for instance, delta spikes and sine waves of the Fourier basis, or the Fourier basis and wavelets. The coherence between the sensing basis Φ and the representation basis Ψ can be given by the following equation:

which, in simple words, is measuring the largest correlation between any two elements of Φ and Ψ. CS is essentially interested in low coherence pairs. For instance, for the previously referred delta spikes and sine waves (time-frequency) pair, *μ*(Φ,Ψ)=1, therefore, indicating maximal incoherence [1, 7].

A particular aspect of interest is that random matrices are largely incoherent with any fixed basis Ψ. This empowers the use of known fast transforms such as a Walsh, Hadamard, or Noiselet transform [8].

Furthermore, what is most remarkable about this concept is that it allows capturing information contained in a sparse signal in a very efficient way without trying to previously understand that signal.

## D.

### How compressive sensing works

Compressive sensing addresses the inefficiencies presented by the *sample-then-compress* framework by directly acquiring a compressed signal representation, avoiding the intermediate stage of acquiring *N* samples [5]. CS bypasses the sampling process and directly acquires a condensed representation *y* consisting of *M* linear measurements. Furthermore, The measurement process is non-adaptive in that Φ does not depend in any way on the signal *x*. The transformation from *x* to *y* is a dimensionality reduction and so loses information in general. In particular, since *M < N*, for a given *y*, there is an infinite number of *x*’ such that Φ*x*’= *y*. The overwhelming capacity of CS is that Φ can be designed such that sparse/compressible *x* can be recovered exactly/approximately from measurements of *y*. To recover the signal *x* from the random measurements *y*, the traditional favorite method of least squares has been shown to fail with high probability. Instead, it has been demonstrated that using the *l*_{1} optimization [4]

it is possible to exactly reconstruct *K-sparse* vectors and closely approximate compressible vectors stably with high probability using just *M ≥O(K·log(N/K))* random measurements [2, 3].

## III.

## SINGLE-PIXEL IMAGING LIDAR SYSTEM

## A.

### Principle of operation

In this work, we chose to implement an active illumination imaging LIDAR system, which compressed the spatial information through the incorporation of the random measurement codes to the light used to illuminate the scene, as described in the setup depicted in Fig. 2.

Observing the diagram presented in Fig. 2, the following principle of operation can be enunciated. The pulsed laser source (Alphalas Picosecond Pulse Diode Laser @ 670 nm) emits a beam, with the desired repetition rate, that will be spatially expanded by a beam expander (Edmund Optics Beam Expander 10x), in order to illuminate the DMD array in which the compressive codes are applied. Then, through a set of projection optics, the images representing the compressive codes applied in the DMD (LogicPD LightCommander) are projected into the field of view being acquired.

The light reflected from the scene being illuminated is then collected by a finder scope (Vixen 7x50 Finder Scope with a field-of-view of 7.2°) and a fiber collimator (Thorlabs Fiber Collimator) ultimately launches the light into an optical fiber that is connected to the single photon counter (PicoQuant tau-SPAD-50). The single photon counter is then connected to a photon counter board (Brilliant Instruments BI220) that is connected to the computer that is controlling the entire system. The software running on the control computer then gathers the timing histograms that will be used to calculate distances and reconstruct the images.

Fig. 3 illustrates the final assembly of the designed system in which it is observable the finder scope assembled on an L-shaped bracket with two rotational stages, so that its pitch and yaw can be adjusted. On the smaller end of the finder scope it is represented the collimator that is connected to the single photon counting module through an optical fiber. On a parallel direction to that of the finder scope, one can see the optical module comprising the projection optics and the DMD array, with the beam expander and the laser support using mechanical parts that have been designed and manufactured.

It should be said that the projection optics and the DMD array have been mounted based on a DMD projection kit produced by LogicPD, called LightCommander, in order to allow one to evaluate the feasibility of the proposed approach at this initial stage.

## B.

### Experimental procedure and results

For the acquisition of experimental data and testing of the system, our efforts relied on the development of the software so that the Imaging LIDAR System could illuminate the scene with the compressive random codes and acquire a timing histogram for each of the codes. These compressive random codes were based on either simple or permuted Hadamard codes. After the acquisition phase, we were then in conditions of reconstructing the ranged images using the timing histograms.

For the generation of the timing histograms the BI 220 photon counter board has been used. The BI220 board outputs the information necessary to produce a timing histogram every time that a predefined number of photons is counted. Therefore, in order to evidence the effect of using the different compressive random codes for each of the histograms a normalization process had to be implemented. For that, the number of counts of each timing histogram has been divided by the total acquisition time of that particular histogram. If this normalization had not been implemented, the effect of the spatial modulation of the light (either imposed by the compressive codes, or by the scene itself) would not be evident in the number of counts of the timing histograms, as the histograms would be considered complete only when the predefined number of photons was counted, and the time it took to complete each histogram would not impact the computation. The BI220 board did not have the possibility to perform measurements during a predefined amount of time, so an alternative strategy had to be implemented. For further clarification, imagine that the single-pixel camera was using a photodiode in forward biased mode instead of an avalanche photodiode operating as a single-photon detector and it was measuring the current passing at the photodiode with a fixed sampling interval. In that case, we could state that the current would vary with different compressive codes as the amount of light hitting the photodiode area would also vary with the compressive codes. This variation imposed by the compressive codes would then allow us to reconstruct the images being acquired [9].

In order to explain how the images are reconstructed based on the histograms, consider the example shown in Table 1.

## Table 1.

Example to explain the reconstruction of ranged images from the timing histograms.

In Table 1 for each compressive code the corresponding histogram is registered. Each histogram contains in each bin a normalized value representing the ratio between the photons counted for that bin and the total acquisition time for that histogram. Then, to reconstruct an image representing the scene at a particular distance, the software used the compressive measurements associated with that particular distance as the input to the reconstruction algorithm. In the example provided, if one wanted to reconstruct an image of a scene at 5.0 m distance from the imaging LIDAR system, one would have to use the measurements of the column highlighted in red.

Some tests to characterize the timing information when the scene was being illuminated and the conditions of acquisition were kept unchanged were also performed. For that, we acquired 100 timing histograms (with 1000 bins) for a predefined amount of detected photons when the scene (white wall) was being totally illuminated (all-white code) and when the scene was being illuminated with a compressive random code. The situations in which the scene was illuminated using an all-white code or image, or when it was being illuminated with a compressive random code can be interpreted as a range finder mode operation or as a compressive imaging mode of operation, respectively. For these experiments, for each of the 100 timing histograms the number of counts of the histogram peak, its occurrence in time and the time required to perform all the measurements have been registered. Table 2 contains the results of these experiments and presents the average of the number of photon counts for each peak normalized by the time spent to acquire all the measurements for the corresponding histogram, the average distance and the average time consumed to acquire a histogram, plus the respective standard deviations (std).

From the analysis of Table 2, one can conclude that the operation mode (range finder [all-white code/image] or compressive imager [compressive random code]) did not affect significantly the distance values obtained from the timing information of the peak in each histogram.

From the information provided in Table 2 it is also obvious that the rate of detected photons along time was higher in the case in which the scene was being illuminated with an all-white image, as the average normalized photon counts have been always higher comparatively to the cases in which the scene was being illuminated with a compressive random code. This makes perfect sense since the amount of light hitting the scene when it is being illuminated with a compressive random code is lower than in the case when an all-white image is used for illumination.

In terms of measurement time, it is clear that it increases with the increase of the number of measurements associated with each histogram. This information is particularly important once it allows one to evaluate the amount of time required to acquire the desired number of histograms that will be used to reconstruct the images. In order to provide an estimate of the compromise that needs to be made, the ratio between the standard deviation and the corresponding average normalized photon counts was calculated. This allowed us to have an insight of the level of noise present in the corresponding measurements. Based on previous work [9], we can state it is possible to reconstruct images of similar quality in a compressive imaging manner when the amplitude of the noise is below 20% of the maximum amplitude of the signal without noise (SNR = 14.54 dB). Therefore, it is possible to infer from Table 2 that a histogram built using a number of measurements above 5000 shall address the stated desired conditions.

## Table 2.

Results for the characterization of the measurements under different illumination conditions.

During operation, in order to ensure ideal conditions, one shall guarantee that the variation imposed by a compressive random code to the normalized photon count number is higher than the standard deviation associated with the current parameterization of the system, specifically the number of photon detections required for each histogram. If this is guaranteed we are in perfect conditions to reconstruct images from the resulting histograms as the effect of noise is strongly reduced.

Next, some results obtained with the system will be presented. To obtain these results, two modes of operation have been defined. One in which the system was operating as a range finder and another in which the system was operating as a 3D imager. During the tests performed the laser pulse repetition rate was set to 10 MHz.

For the case of range finding, the system was completely illuminating a flat target placed at a particular distance. Then, 400 blocks with a predefined amount of measurements were acquired and for each block a histogram of distance values was built (based on the time of flight of the photons that have been detected). Afterwards, to determine the distance to the flat target it was assumed that it corresponded to the distance value associated with the center of the Gaussian trace that best fit each histogram. In the end, 400 distance values have been collected. During these experiments two block sizes have been used – i.e., 5000 and 25000 measurements per block – and the number of bins used to build the histograms has also been varied – i.e., 1000 and 5000 bins per histogram –.

In order to evaluate how spread the 400 distance values were, they were distributed in histograms and were fit using a normal or Gaussian distribution. This time only 20 bins were used to distribute the distance values. The center of the fitting trace then corresponded to the average of the most frequent distance value.

The experiments previously described have been repeated with the target placed at three different distances (differing 1 cm from the immediately previous position) in order to assess the accuracy of the system. Table 3 presents the results obtained when the different parameterizations were used. From the values presented in Table 3, one can state that the system was able to detect the 1 cm displacement of the target. It can also be said that the resolution of the system is much better than 1 cm as the three histograms (in all the cases) are far from overlapping, i.e., their standard deviation is relatively reduced.

## Table 3.

Average distance values obtained from the fittings for the different parameterizations.

Average distance value for target at initial position (m) | Average distance value for target moved 1 cm from initial position (m) | Average distance value for target moved 2 cm from initial position (m) | |
---|---|---|---|

5k measurements 5k bins | 2.5101 | 2.5202 | 2.5308 |

5k measurements 1k bins | 2.5104 | 2.5210 | 2.5310 |

25k measurements 5k bins | 2.5101 | 2.5206 | 2.5310 |

25k measurements 1k bins | 2.5104 | 2.5206 | 2.5310 |

Average +/- Standard Deviation | 2.510 ± 0.002 | 2.521 ± 0.003 | 2.531 ± 0.001 |

To acquire 3D LIDAR images, the system illuminated the scene with the random binary compressive codes and for each code a distance histogram has been recorded. Each of these histograms had 1000 bins and was computed using 50000 measurements (detected photons). In the end, the images corresponding to the desired range of distance values were reconstructed and a datacube was created. For the imaging mode of operation, a scene consisting of a white flat wall and two boxes has been composed (see Fig. 4-left). The box that was lying closer to the wall had a thickness of 10 cm, while the other box had a thickness of 5 cm.

Fig. 4-center and Fig.4-right display a 27 x 29 pixels image and a 3D mesh, respectively, that has been reconstructed and that represent the distance information of the composed scene with color coding. Due to some reconstruction errors near the border of the image, caused by the non-uniformity of the illumination and by the acceptance cone of the collection optics, the 32 x 32 pixels image has been cropped to provide a better visual representation of its content. The lack of uniformity of the illumination derives from the projection optics and from the Gaussian nature of the laser beam profile.

Another test scene that was used for the 3D imaging mode of operation consisted of a futsal ball with 20 cm diameter (see Fig. 5). The results obtained for the test scene of Fig. 5 are displayed in Fig. 6 and were reconstructed using the measurements contained in 40 distance bins of 1024 histograms (with 1000 bins) that were acquired.

Table 4 quantifies the main operational parameters of the implemented prototype. Considering the finder scope angular aperture of 3.6°, at 10 m the system’s field-of-view is defined by a circle with a diameter of 1.25 m.

## Table 4.

Main operational parameters of the single-pixel imaging LIDAR system.

Parameter | 16x16 | 32x32 | 64x64 |
---|---|---|---|

x resolution at 10 m (cm) | 8 | 4 | 2 |

y resolution at 10 m (cm) | 8 | 4 | 2 |

z resolution (cm) | 0.5 | 0.5 | 0.5 |

frame rate (mHz) | 30 | 8 | 2 |

## IV.

## CONCLUDING REMARKS

It was demonstrated that the compressive sensing technique can be used for the effective implementation of 3D imaging LIDAR systems. The implemented prototype demonstrated a xy resolution of 2 cm at 10 m for images with 64x64 pixels, a depth z resolution of 0.5 cm and a frame rate of 30 mHz for 16x16 pixels images. The sampling time was limited by the refresh time of the control of the DMD, being possible to reduce the full datacube acquisition and processing time to values below 1 s using dedicated DMD hardware with low level control. In addition to the LIDAR imaging capabilities, the system that was developed also has range finding capabilities that can be used in two distinct applications. The prototype works in three distinct modes of operation (Range Finder, Attitude Meter and 3D Imager) that can be particularly interesting for space missions, since a single compact system can fulfill the different requirements enumerated for the different sub-phases of those missions.

In face of the performance achieved we feel strongly encouraged to further develop and improve the system towards more demanding requirements with increased technology readiness levels.

In the future, considering the use of commercially available DMD that allow more than 22000 frames to be displayed each second, one can assume that using 50 % compression rate, 15 MHz for the laser pulse repetition rate and 1000 measurements per histogram it may be possible to acquire, every second, 5.3 datacubes for the reconstruction of 16 x 16 pixels images or 1.3 datacubes for the reconstruction of 32 x 32 pixels images. This improved DMD combined with the adequate optics, further empowers the imaging LIDAR system to become even more compact and light weight.

## REFERENCES

Candès, E. J. and M. B. Wakin, “An Introduction To Compressive Sampling”. Signal Processing Magazine, IEEE, 2008. 25(2): p. 21–30.Google Scholar

Candès, E. J., J. Romberg, and T. Tao, “Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information”. Information Theory, IEEE Transactions on, 2006. 52(2): p. 489–509.Google Scholar

Candès, E., J. Romberg, and T. Tao, “Stable signal recovery from incomplete and inaccurate measurements”. Communications on Pure and Applied Mathematics, 2006. 59(8): p. 1207–1223.Google Scholar

Candès, E. J. and T. Tao, “Near-Optimal Signal Recovery From Random Projections: Universal Encoding Strategies”. Information Theory, IEEE Transactions on, 2006. 52(12): p. 5406–5425.Google Scholar

Candès, E. J. and T. Tao, “Decoding by Linear Programming”. IEEE Transactions on Information Theory, 2005. 51(12): p. 4203–4215.Google Scholar

Candès, E. J., “Compressive sampling”, in Int. Cong. Mathematicians. 2006: Madrid, Spain. p. 1433–1452.Google Scholar

Baraniuk, R. G., “Compressive Sensing [Lecture Notes]”. Signal Processing Magazine, IEEE, 2007. 24(4): p. 118–121.Google Scholar

Candès, E. J. and J. Romberg, “Sparsity and incoherence in compressive sampling”. Inverse Problems, 2007. 23(3): p. 969–985.Google Scholar

F. Magalhães, F. M. Araújo, M. V. Correia, M. Abolbashari, and F. Farahi, “Active illumination single-pixel camera based on compressive sensing”. Appl. Opt. 50, 405–414 (2011).Google Scholar