{'fraction_non_alphanumeric': 0.034357021311076895, 'fraction_numerical': 0.018272425249169437, 'mean_word_length': 4.556731202161188, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A `Little Bit' Too Much? High Speed Imaging from Sparse Photon Counts

Recent advances in photographic sensing technologies have made it possible to achieve light detection in terms of a single photon. Photon counting sensors are being increasingly used in many diverse applications. We address the problem of jointly recovering spatial and temporal scene radiance from very few photon counts. Our ConvNet-based scheme effectively combines spatial and temporal information present in measurements to reduce noise. We demonstrate that using our method one can acquire videos at a high frame rate and still achieve good quality signal-to-noise ratio. Experiments show that the proposed scheme performs quite well in different challenging scenarios while the existing approaches are unable to handle them.


Introduction
We address the problem of imaging fast dynamic scenes through single photon counting sensors [14]. Single photon avalanche diode (SPAD) detectors are endowed with the ability of photon counting and time-stamping. They are getting increasingly popular for a variety of imaging applications [7,40,16,33,20]. We attempt to the enhance the fast motion capture capability of SPAD sensors by developing a robust video recovery algorithm. Consider the example shown in Fig. 1. An analog oscilloscope with a traversing sinusoid was imaged by a SwissSPAD camera [4]. Fig.  1 (a), shows a frame from the captured video obtained at 156k frames per second (fps). The value of each pixel in this frame is a binary number wherein positive values indicate the detection of a photon. Note that in a typical DSLR camera, the number of photons captured is of the order of thousands of photons per pixel [32,40] while the frame in Fig. 1 (a) shows the detection of only one photon. Due to the extremely low photon count, one can hardly infer any structure present in the scene. One could average consec-utive frames to reduce the effect of noise while sacrificing the temporal resolution ( Fig. 1 (c)). To preserve temporal resolution as well as recover accurate scene reflectance, we develop a convolutional neural network (CNN) that takes the low-photon count sequence as input and generates a high-photon count estimate at the same frame rate. Our scheme effectively combines the spatio-temporal information present in the input sequence for video recovery. The resultant frame from our method is shown in Fig. 1 (b). Note that in Figs 1 (a) and (b), one can observe the localization of the sinusoidal wave while in Fig. 1 (c), the temporal information is lost. In Fig. 1 (b), we also observe that the details of the static regions have been recovered quite well.
Previously, in the context of time-correlated SPAD imaging, regularization-based approaches have been developed to reconstruct scene reflectivity [40,49]. For oversampled binary observations, image reconstruction algorithms such as [6] are applicable. In this paper, our objective is to jointly recover spatial and temporal variations of radiance in dynamic scenes without spatial oversampling. Existing video denoising schemes devise methods for combining local and non-local structures present across space and time [31,44,46]. In extremely noisy scenarios, explicitly determining such information does not work well. Instead of "hand-crafted" approaches to combine structural information, our method uses convolutional neural networks consisting of 3D filtering across the spatio-temporal volume [45]. By accumulating a set of binary frames, one can obtain video sequences with lesser noise and reduced frame rate. In this paper, we address different scenarios in which different number of the binary frames are combined. Although we show the application of our method on SwissS-PAD cameras, our scheme can be applied to any other photon counting or binary imaging sensors [11,18]. High speed consumer cameras typically require significantly bright illumination [1]. In contrast, we do not use any high intensity illumination and operate in normal lighting conditions.

Related work
Since our work is related to SPAD imaging, photon counting sensors and video denoising, we briefly discuss relevant prior works in these topics.
SPAD-based imaging SPAD sensors are photodetectors in which photon radiation can be detected from the resulting large avalanche currents. SPAD sensor arrays are capable of photon counting at a high speed with high timing resolution and are useful in a variety of applications such as fluorescence lifetime imaging microscopy (FLIM), positron emission tomography (PET), time-of-flight imaging etc. [40,7,26]. Recently, Burri et al. developed a SPAD array known as SwissSPAD [4]. The SwissSPAD is fabricated in a high-voltage CMOS process and features a large 512×128 array with global gating. In this paper, we use data from different SwissSPAD sensor arrays for demonstrating our high speed video recovery scheme.
Quanta Imaging Closely related to the topic of single photon counting is that of Binary/Quanta image sensors [12,14]. These sensors were developed with the objective of shrinking the pixel pitch. In Quanta imaging, the densely defined sensors oversample the scene radiance to generate binary measurements [12,50]. Image reconstruction schemes have been proposed for these sensors [6,9,37,38]. The main difference between these reconstruction methods and our scheme is that these algorithms consider the availability of a higher number of samples to predict the image intensity at a pixel location. These methods have mainly focused on recovering a static image. In contrast, our focus is on recovering videos.
Video denoising Video denoising has been widely studied for many years and different kinds of algorithms have been proposed. Because of the vastness of the topic, we restrict our discussion to some of the popular and recent works. When compared with applying image denoising on each frame, video denoising has an advantage because the high temporal coherence can be leveraged to make a better prediction. Consequently, denoising algorithms adopt different strategies such as non-local means [3,44], motion estimation [29,47], 3D dictionary representations [36] etc. Maggioni et al. [31] propose a denoising method popularly known as VBM4D and is based on the collaborative filtering scheme (VBM3D) [10]. They apply this filtering scheme to a stack of 3-D spatio-temporal volumes that are obtained through non-local grouping. Sutour et al., [44] develop a variational denoising scheme by adaptively combining nonlocal means with total variation on spatio-temporal volumes. Their method adapts to different noise levels. The authors in [22] propose a denoising method for the scenario of mixed noise (like Gaussian noise mixed with impulsive noise). They formulate the denoising problem as low-rank matrix completion problem. This work has been extended in [21] to avoid pre-detection of outliers.
Recently, techniques based on deep learning have been proposed for video denoising. Chen et al. [8] propose a deep recurrent neural network for video denoising. Their method does not assume a specific noise model and performs close to state-of-the art VBM4D scheme [31]. Xue et al. [48] propose task oriented flow to achieve a specific video processing objective such as denoising. Instead of trying to achieve a precise flow estimation, they train a model whose objective is to predict a motion field tailored for a specific task. In our experiments, we observe that video denoising methods fail in the scenario of SPAD imaging. This could be due to the fact that grouping of similar structures across the spatio-temporal volume fails when the number of photons in the observation is sparse.
Contributions The contributions of this paper can be summarized as follows: i) We address the problem of joint spatio-temporal radiance estimation from single photon counting sensors. ii) We develop a CNN-based video recovery scheme that can handle very high noise levels iii) We devise methods to obtain appropriate real and synthetic datasets for training and evaluation. This data will made available subsequently.

Image formation
In this section, we describe the image formation process for SwissSPAD cameras which have a globally gated sensor array [4]. The imaging model can also be applied to any other gated SPAD arrays [11,18].
Each pixel in the SwissSPAD array is composed of a SPAD p-n junction suitably biased for enabling photon triggered avalanche. A one-bit counter is present at all the pixels. The SPAD array has a global shutter in which all the pixels can be kept active for a duration as low as 3.8 ns. The pixel counter content is transferred from the sensor pixels via a fast digital readout which takes about 6.4 µs for transferring the contents of the whole array of 128×512 pixels. In effect, a 1-bit frame of size 128×512 indicating the detection of one photon will be read out in 6.4 µs resulting in a frame rate of 156kfps [4].
Due to the dark counts generated in a SPAD by thermal events, the number of photon counts per unit time follows a Poisson distribution [2]. The probability of k photon counts in unit time is given by where χ is the expected value of counts and is related to the impinging count rate, dark count rate and photon detection efficiency [2]. Within a particular time frame, the SwissS-PAD sensor array can only report whether one or more photons were detected. Consequently, the probability of recording a detection in one readout time is P (count > 0) = 1 − e −χ . Since the number of photons impinging at a pixel depends on the scene radiance corresponding to that pixel, one can conclude that the scene radiance is sensed nonlinearly (according to 1 − e −χ ). This non-linear mapping has been experimentally verified in [2]. Since the SwissSPAD camera records only a 'binary pixel' in each frame, at a time instant t, the intensity I t (i, j) observed at a pixel (i, j) is a Bernoulli sample whose probability depends on the intensity of the corresponding scene point. Fig. 2 (a) shows a single 1-bit image of a resolution chart obtained by the SwissSPAD camera. When 3 such frames were captured and averaged, the resultant image is shown in Fig. 2 (b). This is a 2-bit image corresponding to 2 2 − 1 frames. Similarly Figs. 2 (c), (d) and (e) show 4-bit, 8-bit and 14-bit images. As the number of samples are increased, the noise in the observation reduces. The isolated bright pixels present in all the observations of Fig. 2 correspond to the hot pixels.

Objective
When a binary image sequence I t is averaged upto b bit levels, one would get another sequence denoted by u b τ with frame rate reduced by the factor N b = 2 b −1. The sequence with bit resolution b is given by While imaging a static scene, one can afford to collect as many frames as possible and average them for obtaining an accurate estimate of the scene reflectance map. However, as seen in Fig. 1, while imaging fast dynamic scenes, if one were to average many frames, the temporal information would be lost. Our objective is to overcome this trade-off between frame rate and intensity resolution. For a particular scene, consider that one requires to have a specific frame rate and thereby the bit resolution b gets fixed to a certain level. Let u b τ denote the corresponding sequence. Our aim is to estimate a high bit intensity sequence ub τ whereb is much greater than b and also at the same time preserve the frame rate of the original sequence.
In this paper, we consider the scenarios where N b takes the values of either 1, 3, 7, and 15. i.e., either 1-bit, 2-bit, 3-bit or 4-bit input sequences, respectively. From image sequences of such low photon counts, we attempt to recover sequences corresponding to a very high bit resolutionb > 12 bits. Note that, beyond a certain value of number of bits, any further increase would not be adding new information [13]. We observed that whenb > 12, the noise becomes quite negligible.

Proposed method
We propose to learn a mapping function f for generating a sequence with a high bit resolutionû = f u b ; θ , such that it is close to the noise-free sequence ub. Note that we have dropped the time index τ to simplify the notation. The term θ denotes the network parameters.

Network architecture
Our network architecture is composed of 3D convolutional layers and residual blocks. The network design is motivated by the fact that many image restoration methods employ architectures with a similar structure (but with 2D convolutions) [25,15,23]. Following such an approach also helps to avoid gradient exploding/vanishing problems [19]. Since 2D CNNs in ResNet style have achieved significant success, we intend to use 3D convolutional layers for processing videos. Fig. 3 shows one residual block of our network. Totally, our network consists of K = 3 such residual blocks. In each residual block, all the units are composed of 3D convolutional filters of size 3×3×3. The input convolutional unit consists of one input and 60 output channels. While the three intermediate 'conv' units consist of 60 input and 60 output channels, the output unit has 1 output channel. Except the output unit, all the others are followed by  a 'Leaky' ReLU to model non-linearities. For comparison, we also train with filter size 5×5×5. The input to the network is a spatio-temporal patch. At the end of each residual block indexed by k, the input is added to the resultant of the 'output conv' layer to arrive at u k , an estimate of the spatio-temporal patch corresponding to the clean high-bit sequence. For training the network, we minimize the loss function which is composed of the loss functions of each residual block. Since it is observed in [25] that the Charbonnier penalty function leads to good performance as against the standard 2 penalty, we also use the Charbonnier penalty to define our loss function. i.e., the final loss function is given by where N is the number of training samples in a batch and the penalty term is given by ρ(x) = x 2 + η 2 . The value of η is chosen to be 10 −3 following [25].

Training data
For each bit-level, we train our model using simulated data. For the case of 4-bit sequences, we also train with real pairs of low-bit and corresponding high-bit sequences obtained from SPAD camera. For the 1-bit and 2-bit scenarios, the frame rates that can be achieved by SwissSPAD are 156k fps and 78k fps, respectively. At such high frame rates, one can expect that the temporal variations are quite low. Hence, for these two scenarios, we generate a video dataset with high temporal coherence. The raw videos used in the training of video-deblurring scheme of [42], consists of image sequences at 240fps. The spatial extent of these sequences is 1280×720. We spatially downsample these sequences by a factor of 7 to obtain sequences with reduced variations across time. We randomly crop these downsampled sequences at different temporal locations and arrive at 2500 sequences of dimension 100×100×64. These sequences are directly considered to be the high intensity resolution sequences (ub). To generate the low-bit sequences, for every frame, we average N b Bernoulli sampled (binary) instances. In each Bernoulli instance, the probability of getting a one at a pixel is equal to the corresponding normalized true intensity value at that pixel. i.e., brighter pixels are more likely to generate a 1. To generate such a sample, we generate a uniform random number at every pixel. At a particular pixel, if the randomly generated number is less than the true normalized image intensity value, then that pixel is assigned as one, and zero otherwise. We also randomly add salt-and-pepper-noise at a few points to simulate hotpixels. For 3-bit and 4-bit scenarios, we found that training with videos from UCF 101 dataset [41] that have a normal frame-rate was sufficient.
Model training We trained our network models on a NVIDIA GeForce GTX 1080 Ti using stochastic gradient descent without batch normalization. The models named 3DCNN1B, 3DCNN2B, 3DCNN3B and 3DCNN4B denote the CNNs trained with 1-bit, 2-bit, 3-bit and 4-bit sequences, respectively. The network trained using 4-bit real data is referred to as 3DCNNR. From the 2500 sequences, we use 2400 for training and the rest for testing. The input in each batch of training is obtained by randomly cropping 3D patches of size 60×60×38 from any of the sequences of the training dataset. We apply random scaling by small amounts, flipping, and rotation by multiples of 90 • for data augmentation.

Experiments
To run our algorithm on a GPU for videos of realistic size, we divide an input sequence into overlapping 3D spatio-temporal patches and subsequently merge the out-One-bit frame Two-bit frame 3DCNN1B result 3DCNN2B result Ground truth  puts. For sequences of a particular bit-level, we use its corresponding CNN. If the bit-level of an input is not known, it can be determined easily by checking the number of unique levels in the pixel intensities.

1-bit and 2-bit synthetic experiments
We used the test set of simulated data to evaluate the performance of our video reconstruction method. Table 1 shows the peak signal to noise ratio (PSNR) and structural similar index measure (SSIM) averaged over all the 100 test sequences. Note that 3DCNN1B and 3DCNN2B had different inputs, but the ground-truth sequences were the same as seen in the repre-  sentative examples shown in Fig. 4. In Table 1, the value 5 within parentheses indicates that the filter size used in the CNN was 5×5×5 instead of 3×3×3. To check if any other denoising method works for this scenario, we applied the algorithms proposed in [31], [44], [48], and [30] on five of these image sequences. None of these methods were able to restore videos for these sequences. We varied the parameters of the algorithms and searched for the optimal values. Out of these other algorithms, the best performance was obtained from [44] (when applied with Poisson noise statistics). However even this is not quite satisfactory. The SSIM values of the outputs from [44] ranged between 0.5 and 0.64. For visual comparison, we show one example output of [44] in Fig. 6. The synthetic experiments clearly demonstrate that existing video denoising methods cannot be used for observations with sparse photon counts. However, we subsequently notice improvements in their performance when the bit-level improves. We also checked if a single network can be used to recover both 1-bit and 2-bit sequences. For this purpose, we trained another CNN wherein the inputs for training consisted of both 1-bit and 2-bit sequences (with equal probability). We observe a slight reduction in the performance of this jointly-trained CNN compared to that of the specific networks seen in Table 1. For the case of 1-bit test sequences, the resultant mean PSNR and SSIM from the jointly-trained CNN are 24.97 and 0.726, respectively. For the 2-bit test sequences, the resultant mean PSNR and SSIM are 26.9 and 0.80, respectively. We have also tested the proposed networks on video sequences with downsampling factors less than 7. In these test sequences, the spatiotemporal variations will be increased when compared to the training data. These results are reported in the supplementary material.  Real experiments We initially show an example of 1-bit real sequence. Fig. 7 shows a scene wherein the SwissS-PAD camera was placed close to a rotating tool. In this particular setup, because of limitations in the data-transfer rate, 1-bit frames were captured at the rate of 42 kfps. A binary image from the input sequence is shown in Fig. 7 (a). Its corresponding output is shown in Fig. 7 (b). For comparison, we averaged 15 binary frames to generate a 4bit sequence. This 4-bit sequence was fed as input to our 3DCNN4B network. In the videos included as supplementary material, one can clearly see the loss of temporal resolution in the output of 3DCNN4B when the averaged sequence was input. Fig. 8 shows another rotating tool present in a static background. The input in this scenario was a 2-bit sequence captured at about 25 kfps. While Fig. 8 (a) shows an input frame, Fig. 8 (b) shows the corresponding result from 3DCNN2B. By averaging five frames of this sequence, one can obtain a 4-bit sequence. Figs. 8 (c) and (d) show a 4-bit input and a resultant (3DCNN4B) frame, respectively. We observe that in the static regions, the output of 3DCNN2B does come close to that of 3DCNN4B. The supplementary material contains a comparison of our output with a high-intensity resolution reference image captured when the scene was still.
We next show results on additional oscilloscope sequences of [4] (Fig. 5). The 2-bit sequence was fed as input to 3DCNN2B, and 4-bit sequence was input to 3DCNN4B. One can observe that the quality of our 2-bit output does come close to that of 4-bit and even with the high-bit observation ( Fig. 5 (f)) at static regions. This shows that our algorithm is capable of producing high-quality outputs from only 2-bit frames. The output from [44] on the 2-bit sequence (Fig. 5 (c)) looks quite inferior to our output. The supplementary material contains complete videos and comparisons.
Subsequently, we compare the performances of 3DCNN1B, 3DCNN2B and 3DCNN3B on the same scene. The image sequence corresponds to breaking of glass with a resolution chart in the background. The scene was captured at 156kfps and at 1-bit resolution. We divided the sequence into groups of seven frames. For 1-bit observations, we keep the central frame in each group and drop six other frames. For 2-bit observations, we average the middle three frames and drop the rest. For 3-bit observations, we average all the seven frames in the group. Essentially, we arrive at 1-bit, 2-bit and 3-bit sequences with similar temporal variations and corresponding to about 22kfps. With these inputs, we obtain the outputs from our restoration scheme. Fig. 9 shows inputs and outputs at two different instants of time. In the second instant, we see that the particles have been scattered after the breaking of glass. Despite the inputs being highly noisy, we see that the structural information has been recovered quite well. This shows that even with just one-bit measurements, our network model is able to reconstruct scene information quite robustly. 4-bit sequences We quantitatively evaluate different video denoising methods using our real SPAD dataset that has both 4-bit noisy sequence and the corresponding high-bit sequence. We evaluate the performance on ten randomly selected test-image sequences. The performance of different schemes are presented in Table 2. The algorithms of [31] and [44], do not handle hot pixels. For these methods, to evaluate score, we replace the intensity of a hot-pixel by the value equal to the median of the 3×3 neighborhood of that hot pixel (excluding the damaged pixels while calculating median). The other algorithms can handle outliers and this step is not necessary. The table shows that our scheme of 3DCNNR clearly outperforms other methods. A representative example from the SwissSPAD dataset is shown in Fig. 11. On close observation of different regions, we can see that the reconstruction is more faithful in the proposed CNN model.
We next show a result on a 4-bit sequence captured at 12 kfps. In Fig. 10, we have embedded some of the frames of the video sequences. From the videos, one can clearly observe that the proposed method performs quite well. There has been an improvement in the performance of video denoising schemes when compared to the 1-bit and 2-bit scenarios. However, in the results of [31], [44] and [48] (supplementary material), we do observe artifacts clearly in regions where there is no motion.

Conclusions
We proposed a video recovery scheme for single-photon counting cameras with sparse photon counts. The performance of our model is quite good in real scenarios despite the training on simulated data. The level of performance achieved by our method in 1-bit and 2-bit scenarios is not possible with any existing approach. Even for 3-bit and 4bit scenarios, our method outperforms existing video denoising schemes.
In other applications of SPAD cameras such as range imaging and time-of-flight imaging, the number photon counts is of much higher magnitude than the numbers seen (c) (d) Figure 10. Real example of a balloon bursting sequence. Frames from (a) input sequence, result of (b) [31], (c) [44] and (d) 3DCNN4B.
in this paper. Our work could serve as a template for developing more photon-efficient techniques to perform these tasks. We would explore this direction in future.