Towards quantum 3D imaging devices: the Qu3D project

A relevant part of the speedup that we are seeking is determined by replacing commercial high-resolution sensors, like scientific CMOS and EMCCD cameras, with sensors based on cutting-edge technology such as single-photon avalanche diode (SPAD) arrays. A SPAD is basically a photodiode which is reversely biased above its breakdown voltage, so that a single photon which impinges onto its photosensitive area can create an electron-hole pair, triggering in turn an avalanche of secondary carriers and developing a large current on a very short timescale (picoseconds) [7, 8]. This operation regime is known as Geiger mode. The SPAD output voltage is sensed by an electronic circuit and directly converted into a digital signal, further processed to store the binary information that a photon arrived, and/or the photon time of arrival. In essence, a SPAD can be seen as a photon-to-digital conversion device with exquisite temporal precision. SPADs can also be gated, in order to be sensitive only within temporal windows as short as a few nanoseconds, as shown in Fig. 3. Individual SPADs can nowadays be used as the building blocks of large arrays, with each pixel circuit containing both the SPAD and the immediate photon processing logic and interconnect. Several CMOS processes are readily available and allow to tailor both the key SPAD performance metrics and the overall sensor or imager architecture [33, 34]. Sensitivity and fill-factor have for some time lagged behind those of their scientific CMOS or EMCCD counterparts, but have been substantially catching up in recent years.

Based on the requirements of QPI, we have chosen to employ the SwisSPAD2 array developed by the AQUA laboratory group of EPFL, characterized by a 512x512 pixel resolution (see Fig. 4), which is one of the widest and most advanced SPAD arrays to date [10, 12]. The sensor is internally organized as two halves of 256x512 pixels to reduce load and skew on signal lines and enable faster operation. It is a purely binary gated imager, i.e. each pixel records either a 0 (no photon) or a 1 (one or more photons) for each frame, with basically zero readout noise. The sensor is controlled by an FPGA generating the control signals for the gating circuitry and readout sequence and collecting the pixel detection results. In the FPGA, the resulting one-bit images can be further processed, e.g. accumulated into multi-bit images, before being sent to a computer/GPU for analysis and storage. The maximum frame rate is 97.7 kfps, and the native fill factor of 10.5% can be improved by 4-5 times, for collimated light, by means of a microlens array [11] (higher values are expected from simulation after optimization); the photon detection probability is 50% (25%) at 520 nm (700 nm) and 6.5 V excess bias. The device is also characterized by low noise (typically less than 100 cps average Dark Count Rate per pixel at room temperature, with a median value about 10 times lower) and advanced circuitry for nanosecond gating. A detailed comparison of SwissSPAD2 with other large-format CMOS SPAD imaging cameras is presented in Ref. [12].

Currently, we are using the available version of SwissSPAD2, with a 512x256 pixel resolution, to generate sequences of frames and store them in an on-board 2 GB memory, before transferring them to a computer by a standard USB3 connection, which can be done using existing hardware. We are integrating this sensor in the prototype of chaotic-light base quantum plenoptic camera in a way that two disjoint halves of the sensor (of 256x256 pixels each) are used for retrieving the images of the two reference planes. The high speed of this sensor is expected to reduce the acquisition time of quantum plenoptic images by 2 order of magnitudes with respect to the first CPI experiment [6], in which the region of interest on the sensor was made of $700\times 700$ pixels for the spatial measurement side, and $600\times 600$ for the directional measurement side; $2\times 2$ binning on both sides during acquisition and a subsequent $10\times 10$ binning on the directional side led to the effective spatial resolution of $350$ pixels, and angular resolution of $25\times 25$ pixels.

We are also working towards several further optimizations of the sensor system, e.g. by developing gating for noise reduction in QPI devices. In order to employ the full sensor (512x512 pixels), we are implementing a synchronization mechanism for a pair of imagers by means of two FPGAs, so as to operate on a common time-base at the nanosecond level (to this end, two control circuits shall operate from a single clock and have a direct communication link). Finally, the SwissSPAD2 arrays are being integrated with a fast communication interface in order to speed up data transfer and make it possible to deliver, in a sustained way, full binary frame sequences to a GPU. The latter will run advanced algorithms for data pre-processing, image reconstruction and optimization, as we will discuss in more detail in the next subsection.

In a QPI device integrated with SwissSPAD2, the acquired data rate for a single frame acquisition can be estimated to 26 Gb/s, which is beyond the reach of standard data buses. This poses great challenges in both hardware and software design. Our approach is to employ the expertise of PKH to seek careful design of electronic interconnections (buses) between sensor control electronics and processing device, theoretical refinement and optimization of algorithms (e.g., compressive sensing [35, 36]), porting to an efficient computational environment, and design of a specific acquisition electronics for optimizing data flow from the light sensors to a dedicated processing system, able to guarantee the required computational performance (e.g., exploiting GPUs or FPGA) [13, 14].

The introduction of an embedded data acquisition- and processing board, integrating a GPU, aims at data pre-processing, thus significantly reducing the amount of data to be transferred to (and saved on) an external workstation. GPUs exploit a highly parallel elaboration paradigm, enabling to design algorithms that run in parallel on hundreds or thousands of cores and to make them available on embedded devices. A great advantage of GPUs is programmability: many standard tools exist (e.g. OpenCL and CUDA) that allow fast and efficient design of complex algorithms that can be injected on the fly in the GPU memory for accomplishing tasks ranging from simple filtering to advanced machine learning. Efforts are made to design a processing matrix, so that each line and column of the sensor will be managed by a dedicated portion of the heterogeneous processing platform (CPU/GPU/FPGA): the pixel series processor. Those dedicated units will be interconnected to one another to cooperate for implementing algorithms that require distributed processing on a very small scale.

The embedded acquisition-processing board is designed to best fit to SPAD array and SW application needs. The optimal system design will be evaluated by considering theoretical algorithms, engineered SW implementations, HW set-up and HW/SW trade-off. We will identify a preliminary set of possible configurations and perform a trade-off by comparing overall performances, considering the requirements in data quality, processing speed, costs, complexity, etc.

Based on the challenging objectives to be achieved, a preliminary analysis based on COTS (Commercial-Off-The-Shelf) solutions was performed, in order to identify a set of accelerating devices addressing QPI requirements in terms of computing capability and portability. The option offered by the NVIDIA Jetson Xavier AGX board shown in Fig. 5 is considered a promising candidate to achieve our goal. This device indeed offers an encouraging performance/integration ratio with low power usage and very interesting computing capabilities. Its main characteristics are reported in Table  1.

Considering the listed HW capabilities, despite the ARM processor has a limited computing power when compared to high-end desktop processors, it allows to leverage multi-core capabilities for implementing that part of the code that will feed the quite powerful GPU device on board. In addition, the foreseen optimization strategy will require a dedicated implementation able to consider the maximum amount of memory of 16GiB available, shared with the GPU. Also, the solution to be developed should take into account the bandwidth of about 136GB/s of the on-board memory, which may represent a limiting factor when GPU and CPU exchange buffers. An implementation based on the CUDA framework -over OpenCL or other technologies- will be preferred to best use the NVIDIA device. Finally, given the capability of the NVDLA Engines to perform multiplications and accumulations in a very fast way, we consider it interesting to perform an assessment of how to exploit these devices for implementing the QPI-specific correlation functions and/or other multiplication/sum intensive computations.

Along with the assessment of a dedicated HW solution, further optimizations applicable at the algorithmic level have been considered in order to enrich the engineered device with a highly customized algorithmic workflow able to exploit the peculiarities of the CPI technique and its related input data. More specifically, we are analysing those steps of QPI processing that appear as more computationally demanding, thus representing a bottleneck for performances. To this end, tailored reshaping operations applied over the original 3-dimensional multi-frame structure of the input data were explored, to facilitate the development of a parallelized elaboration paradigm for evaluating the CPI-related correlation function. Besides, the peculiar feature of the input dataset to be acquired by SPAD sensors as one-bit images will be valued through dedicated implementations able to gain from intensive multiplication/sum math among binary-valued variables.

In order to reduce the amount of required data (at present $10^{3}$–$10^{4}$ frames) by at least 1 order of magnitude, we are exploring different approaches. First, we are investigating the opportunity of implementing compression techniques for improving bus bandwidth utilization, thus acting on data transfer optimization. Data compression may also rely on manipulation of raw input data to determine only the relevant information in a sort of information bottleneck paradigm, in which software nodes in a lattice provide their contribution to a probable reconstruction of the actual decompressed data. This is very similar to artificial neural network structures, where the network stores a representation of a phenomenon and returns a response based on some similarity rating versus the observed data. Moreover, the availability of advanced processing technologies allows the investigation and implementation of alternative techniques able to exploit the sparsity of the retrieved signal to reconstruct information from a heavily sub-sampled signal, by compressive sensing techniques [35, 36].

As in other correlation imaging techniques, such as ghost imaging (GI), in the CPI protocol the object image is reconstructed by performing correlation measurements between intensities at two disjoint detectors. Katz et al. [37] demonstrated that conventional GI offers the possibility to perform compressive sensing (CS) boosting the recovered image quality.

CS theory asserts the possibility of recovering the object image from far fewer samples than required by the Nyquist–Shannon sampling theorem and it relies on two main principles: sparsity of the signal (once expressed in a proper basis) and incoherence between the sensing matrix and the sparsity basis.

In conventional GI, the transmission measured for each speckle pattern represents a projection of the object image and CS finds the sparsest image among all the possible images consistent with the projections. In practice, CS reconstruction algorithms solve a convex optimization problem, seeking for the image which minimizes the $L_{1}-$norm in the sparse basis among the ones compatible with the bucket measurements, see Refs.[38, 39, 40, 41] for a review.

We are developing a novel protocol reducing the number of measurements required for image recovery by an order of magnitude. Once properly refocused, a single acquisition can be fed into the compressive sensing algorithm several times thus exploiting the plenoptic properties of the acquired data and increasing the signal-to-noise ratio of the final refocused image. We tested the CS-CPI algorithm with numerical simulations, as summarized in Fig. 6. In panel (a) a double-slit mask is reconstructed by correlations measurements considering $N=6000$ frames, in panel (b) the standard reconstruction is repeated considering only the $10\%$ of available frames, while in panel (c) the CS reconstruction using the same reduced set of measured data. In addition to a data-fidelity term corresponding to a linear regression, we penalized the $L_{1}-$norm of the reconstructed image to account for its sparsity in the $2D-DCT$ domain. The resulting optimization problem is known as the LASSO (Least Absolute Shrinkage and Selection Operator) [42]. We employed the coordinate descent algorithm to efficiently solve it and we set the regularization parameter, controlling the degree of sparsity of the estimated coefficients, by cross-validation. In this proof-of-concept experiment, we simply use Pearson’s correlation coefficient to measure the similarity between the reconstructions obtained using the restricted dataset and the image obtained considering all the $N=6000$ frames.

Novel reconstruction algorithm with the advantage of real 3D image lack of artifacts is based on the idea of recasting plenoptic imaging as an absorption tomography. The classical refocusing algorithm is based on superimposing images of the 3D scene from different viewpoints, which necessarily results into the contribution of out-of-focus parts of the scene, the effect responsible for the existence of blurred part of 3D image reconstruction. The tomography approach is based on a different principle and provides sufficient axial and transversal resolution without artifacts. For this purpose, the object space is divided into voxels and the goal is to reconstruct absorption coefficients of each voxel. The measured correlation coefficient of a pair of points from the correlated planes is transformed to be linearly proportional to the attenuation along the ray path connecting those two points. Classical inverse Radon transform can be used in this scenario to obtain a 3D image but using the Maximum Likelihood absorption tomography algorithm [47] further enhances the quality of the tomography reconstruction and performs well even for a small range of projections angels. In fact, advanced tomographic reconstruction algorithms based on the Maximum Likelihood principle are more resistant to noise and require fewer acquisitions, for a given precision, in comparison to standard tomographic protocols. We are investigating special tools to deal with informationally incomplete detection schemes for very high resolution, and optimal methods for data analysis based on convex programming tools.

A further quantum approach to image analysis and detection schemes will be employed to achieve super-resolution (or, eventually, for maintaining the desired resolution and speeding up acquisition and elaboration of the quantum plenoptic images) and to compare correlation plenoptic detection scheme to the ultimate quantum limits. The basic concept underpinning the Fisher Information super-resolution imaging is the formal mathematical analogy between the classical wave optics and quantum theory, which makes it possible to apply the advanced tools of quantum detection and estimation theory to classical imaging and metrology [44, 19]. The dvantage of this approach lies in the ability to quantify the performance of the imaging setup based on rigorous statistical quantities. Inspired by the quantum theory of detection and estimation, quantum Fisher information, quantity connected with the ultimate limits allowed by nature, is computed for simple 3D imaging scenarios like localization and resolution of two points in the object 3D space [45, 43, 46, 20]. For example, one might be interested in measuring the separation of two point-like sources and seek the optimal detection scheme extracting the maximum amount of information about this parameter. We shall thus employ quantum Fisher information to design optimal measurement protocols within the quantum plenoptic devices, able to extract specific relevant information, for enhanced resolution, with a minimal number of acquisitions. The question of the minimal number of detections over the correlation planes which achieves acceptable reconstruction quality is of relevance because of direct connection to the amount of processed data and can lead to reduced time for data processing.