Mozhdeh Seifi

Loic Denis
email: loic.denis@univ-st-etienne.fr

Corinne Fournier

Fast and accurate 3D object recognition directly from digital holograms

Pattern recognition methods can be used in the context of digital holography to perform the task of object detection, classification and position extraction directly from the hologram rather than from the reconstructed optical field. These approaches may exploit the differences between the holographic signatures of objects coming from distinct object classes and/or different depth positions. Direct matching of diffraction patterns however becomes computationally intractable with increasing variability of objects due to the very high dimensionality of the dictionary of all reference diffraction patterns. We show that most of the diffraction pattern variability can be captured in a lower dimensional space. Good performance for object recognition and localization is demonstrated at a reduced computational cost using a low dimensional dictionary. The principle of the method is illustrated on a digit recognition problem and on a video of experimental holograms of particles.

Introduction

In-line digital holography is a 3D imaging method which is being increasingly used in applications requiring to track micro-objects from fluid mechanics [START_REF] Katz | Applications of holography in fluid mechanics and particle dynamics[END_REF][START_REF] Coupland | Optical tomography and digital holography[END_REF][START_REF] Hinsch | Holographic particle image velocimetry[END_REF] to biology [START_REF] Garcia-Sucerquia | Digital in-line holographic microscopy[END_REF][START_REF] Allier | Bacteria detection with thin wetting film lensless imaging[END_REF][START_REF] Zhang | Lensless imaging for simultaneous microfluidic sperm monitoring and sorting[END_REF][START_REF] Mudanyali | Wide-field optical detection of nanoparticles using on-chip microscopy and self-assembled nanolenses[END_REF]. The advantage of this lensless technique compared to its alternatives is the simple imaging setup with high video frame rate resulting in accurate reconstructions of the volume of interest. During the recording step in this method, the objects are illuminated by a laser beam and the diffracted wave intensity is recorded in a video by a camera. The recorded sequence of holograms is then processed to reconstruct the volume's evolution over time.

Among the existing hologram processing techniques, model fitting methods have shown to be the most efficient and accurate for metrological applications [START_REF] Lee | Characterizing and tracking single colloidal particles with video holographic microscopy[END_REF][START_REF] Cheong | Strategies for three-dimensional particle tracking with holographic video microscopy[END_REF][START_REF] Soulez | Inverse-problem approach for particle digital holography: accurate location based on local optimization[END_REF][START_REF] Fung | Measuring translational, rotational, and vibrational dynamics in colloids with digital holographic microscopy[END_REF]. The contribution of this paper is to show that it is possible to go beyond model fitting by performing diffractionpattern matching to identify the class and 3D location of objects directly from holograms. Our approach is inspired by the method for 3D object recognition based on hologram correlation suggested in [START_REF] Javidi | Three-dimensional object recognition by use of digital holography[END_REF][START_REF] Kim | Extraction of 3-d location of matched 3-d object using power fringe-adjusted filtering and wigner analysis[END_REF]. However, rather than matching two experimental holograms corresponding to the two objects to be compared [START_REF] Kim | Three-dimensional matching by use of phase-only holographic information and the wigner distribution[END_REF], we consider matching a hologram with a whole dictionary of diffraction patterns generated for various objects and all their possible 3D locations. As observed for the task of accurate localization of particles in holographic metrology, very good performance can be reached by searching for elements of a dictionary of diffraction patterns that best match the data. This approach contrasts with conventional hologram processing methods that analyze a coarse reconstruction obtained by hologram back-propagation using Fresnel transform or other optical propagation approximations.

In-line holography of opaque objects is restricted to mostly empty volumes (typically, the projection of all objects should cover less than 1% of the sensor, see [START_REF]An application of high-speed microholography: the metrology of fogs[END_REF]). The hologram can then be well approximated as the incoherent superimposition of the diffraction pattern of each object. Hologram reconstruction is an inverse problem that can also be seen as a linear un-mixing problem, or as a multiple-object pattern-recognition problem with overlapping patterns. In the context of particle hologram analysis, a greedy algorithm has been proposed [START_REF] Soulez | Inverse-problem approach for particle digital holography: accurate location based on local optimization[END_REF][START_REF] Soulez | Inverse problem approach in particle digital holography: out-of-field particle detection made possible[END_REF]. This algorithm builds incrementally a model of the hologram by fitting at each step the parametric diffraction-pattern model of a single (spherical) particle. When considering a large diversity of particle diameters and distances to the hologram, the collection of all diffraction patterns to be considered becomes huge and the computational cost rises for reconstruction of hologram videos. The extension to more general classes of objects is restricted by the explosion of the size of the dictionary of all diffraction-patterns that must be matched with the hologram.

To overcome this computational issue, we propose to approximate the dictionary in a low-dimensional space. We show that such approximation preserves enough variability between diffraction patterns to perform correct discrimination and accurate 3D localization. The dimensionality reduction is based on a singular value decomposition (SVD) of the dictionary to identify most significant modes. It leads to an efficient diffraction-pattern matching algorithm.

The structure of this paper is as follows: the pattern recognition problem is introduced in the next section. Low-dimensional approximation of the dictionary using SVD is described in section 3. We report numerical and experimental results in section 4 and we conclude the paper in Sec. 5.

Object recognition from a digital hologram

The majority of hologram reconstruction methods reconstruct the optical wave-field using light back-propagation techniques. After reconstruction of the wave-field at different depths, objects can be segmented and recognized in their maximum-of-focus plane. This approach is illustrated in Fig. 1, on a simulated hologram of digits. Since either the in-focus or out-offocus signature of a given object is visible on all reconstructed planes, selection of the in-focus plane and separation of overlapping in-focus and out-of-focus images from different objects can be hard. Better reconstructions can be obtained by stating hologram reconstruction as an inverse problem. Regularized reconstruction using an oscillation-preventing term (total variation penalty) [START_REF] Bourquard | A practical inverse-problem approach to digital holographic reconstruction[END_REF][START_REF] Marim | Compressed sensing with off-axis frequency-shifting holography[END_REF][START_REF] Lim | Sampling and processing for compressive holography [Invited[END_REF][START_REF] Rivenson | Compressive fresnel holography[END_REF][START_REF] Choi | Compressive holography of diffuse objects[END_REF][START_REF] Zhang | Edge-preserving sectional image reconstruction in optical scanning holography[END_REF] or a sparsity constraint (ℓ 1 norm minimization) [START_REF] Denis | Inline hologram reconstruction with sparsity constraints[END_REF] have been shown to greatly enhance the quality of the reconstructed volumes. Compared to performing object recognition on the resulting volumes, it is more robust to directly recognize diffraction patterns in the hologram. In a Bayesian perspective, this can be understood as making use of a much stronger prior when solving the inverse problem of hologram reconstruction, since we specify which objects are expected to be present in the volume. Under the approximation of a hologram as a linear super-imposition of the diffraction pattern of each object, the hologram formation model (direct problem) can be formulated as

d = Mα + ϵ , (1)
where d is the vector of data, i.e., the N pixel values of the hologram re-ordered as a vector in R N . In this formulation, M represents a dictionary with columns containing all possible diffraction patterns of the family of objects under consideration. Two diffraction patterns may differ depending on the class / shape and 3D location of the objects (e.g., for opaque and spherical particles, diffraction patterns depend on the radius r and 3D location of the particle (x, y, z)). A given hologram d may be described as the weighted sum of few diffraction patterns from the dictionary M. Let m i k be the diffraction pattern created by the k th object (m i k is column i k of dictionary M). The corresponding component i k of vector α is then non-zero and defines the amplitude of that diffraction pattern on the hologram. Vector α has mostly zero components except for the amplitude of the diffraction pattern created by each object (i.e., α is a very sparse vector). Vector ϵ accounts for modeling errors and noise in the data. In the absence of specific information on the noise statistics, it will be considered white and Gaussian. Recognition of objects directly from the hologram d amounts to determining the support of vector α (i.e., location of non-zero entries). Numerous methods have been developed in the literature of compressed sensing and sparse reconstruction to solve such problem. Given the very limited number of objects present in the volume compared to the size of the dictionary M (many different shapes and 3D locations can be considered for a given object), vector α is extremely sparse. Greedy algorithms such as matching pursuit and its variants are well suited to this problem [START_REF] Tropp | Greed is good: Algorithmic results for sparse approximation[END_REF]. Object identification and detection can then be performed directly from the hologram by iterating the following steps, as proposed in the case of particle holograms in [START_REF] Soulez | Inverse-problem approach for particle digital holography: accurate location based on local optimization[END_REF][START_REF] Soulez | Inverse problem approach in particle digital holography: out-of-field particle detection made possible[END_REF]:

1. An exhaustive search for the diffraction pattern in the discrete dictionary M that best improves current modeling of the data. This step corresponds to the identification of the column j of dictionary M that is the most correlated with residuals r representing the part of data that is yet to be modeled:

arg max j m t j Wr , (2)
where W is a weighting matrix, often referred to as the precision matrix, that can account for noise correlation (W is then the inverse of noise covariance matrix), nonstationary noise variance (with a diagonal W matrix and non-constant diagonal values) or unmeasured areas (with zero values on the corresponding diagonal entries). In Eq. 2, the dictionary is supposed to be properly normalized (∀i, m t i Wm i = 1). Residuals r are initially equal to the data d and iteratively updated in step 3 described bellow.

2. A local optimization step can possibly be performed to refine shape and/or location with sub-pixel/sub-sampling accuracy. This can be done by non-linear least-squares fit when a continuous model of the diffraction patterns is available (e.g., for spherical particles).

3. A cleaning step that updates the residuals by removing the identified diffraction pattern m j . This step plays an important role as it prevents from detecting the same object several times (its contribution being removed from the residuals once and for all) and it improves the contrast of remaining diffraction patterns.

The first step requires to perform an exhaustive search on a very large dictionary which limits the applicability of the method to strongly constrained families of objects such as spherical particles. The strong hypothesis is that the object to be identified has a representative in the dictionary. The extension to more general classes of objects requires a way to speed up this search. We describe in the next section how to approximate the dictionary in a lower dimensional space and derive a fast method to perform this diffraction-pattern matching step.

Low-dimensional approximation of the dictionary of diffraction patterns

To reduce the inter-pattern redundancy of the dictionary, we propose to approximate it in a low-dimensional space. We first distinguish variations due to the transversal location of the object (i.e., its (x, y) location) from other variations. Let C be the restriction of dictionary M to diffraction patterns that are centered on the optical axis (i.e., only columns corresponding to patterns that are geometrically centered are kept). C captures all variations other than transversal shifts. Since diffraction patterns display some correlations, matrix C can be well approximated by a lower rank matrix obtained by singular value decomposition (SVD):

C ≈ k ∑ i=1 u i s i v t i , (3)
where u i and v t i represent the i th left and right singular vectors associated with singular value s i from the SVD of matrix C, and k is the rank chosen for the approximation. In this approximation, the j th diffraction pattern is given by:

c j ≈ k ∑ i=1 u i β i,j , with β i,j = s i • v i,j . (4)
The orthonormal basis formed by {u i } 1..k that spans the k-dimensional approximation subspace is called the dictionary of diffraction-pattern modes. Any diffraction pattern c j is approximated as a linear combination of modes u i , with weights β i,j . Since the array of weights β depends only on the dictionary C, it can be computed off-line once for all.

Correlations between the dictionary and residuals must be evaluated in order to identify the best matching pattern in a hologram. Let κ i denote the correlation between mode u i and the residuals. Using our low-rank approximation, correlations in Eq. 2 can be expressed as a linear combination of few mode/residuals correlations κ i :

c t j Wr ≈ k ∑ i=1 β i,j κ i , with κ i = u t i Wr . (5)
As mentioned earlier, for a given dictionary C it is necessary to compute coefficients β i,j only once. The computation of k mode/residual correlations κ i then suffices to approximate any correlation between a dictionary element c j and the residuals, resulting in a large speedup compared to direct computation of the correlation c t j Wr as long as the number of diffraction patterns in the dictionary C is much larger than the number of retained modes k.

Identification of the best-matching diffraction pattern in Eq. 2 not only involves geometrically-centered patterns c j but also all possible (x, y) shifts. Correlation with all shifted versions of a given geometrically centered pattern can be efficiently computed using fast Fourier transforms (FFT) with zero-padding to prevent from periodization artifacts. Two shifted versions m 1 and m 2 of a pattern c generally differ in norm (m t 1 Wm 1 ̸ = m t 2 Wm 2 ), so that a proper normalization must be computed. The pattern matching step given in Eq. 2 is then replaced by: arg max where normalization of the dictionary M is no longer required.

j m t j Wr m t j Wm j , (6)
Fast computation of the numerator in Eq. 6 for all integer pixel shifts of a pattern c can be performed with FFT by computing correlations between c and Wr in Fourier domain. The normalization term corresponding to the denominator depends only on the dictionary and can either be computed off-line or be approximated using another SVD. When weighting matrix W is diagonal (W = diag(w)), the denominator can be computed for all integer shifts by computing with FFT the correlations between w and diag(m)m, as done in [START_REF] Soulez | Inverse problem approach in particle digital holography: out-of-field particle detection made possible[END_REF].

Using the approximation derived in Eq. 5, the correlation between any diffraction-pattern from the dictionary M and the residuals r can be approximated as a linear combination of mode/residuals correlations. It is necessary to compute those mode/residuals correlation maps for all integer shifts only once for each mode. Each of these k maps can benefit from the fast FFT-based implementation described.

Fig. 2 compares the computationally expansive matching of diffraction patterns taken from the whole dictionary M in (a) with the use of the low-dimensional approximation based on diffraction-pattern modes in (b). For a dictionary made of n different geometrically-centered patterns (i.e., C has n columns), computation of the correlations for all integer shifts has a complexity of the order of n × N log(N ) whereas using the proposed approximation reduces the complexity to k × N log(N ) where the number of modes k can be generally chosen much smaller than n, as illustrated next section. Note that the size of matrix C can be very large. Rather than computing its full SVD, it is preferable to use (approximate) iterative methods to compute a truncated SVD. We used a method [START_REF] Simon | Low rank matrix approximation using the lanczos bidiagonalization process with applications[END_REF] based on Lanczos bi-diagonalization [START_REF] Larsen | Lanczos Bidiagonalization with Partial Reorthogonalization[END_REF] of C to simplify the calculation of the SVD of a large sparse matrix.

Experimental results

In this section we report the results of our method on both simulated and real holograms. The implemented method uses FFTW [START_REF]FFTW[END_REF] library and OpenMP [START_REF]OpenMP[END_REF] to exploit multi-threading on a six-core Xeon CPU for the calculation of the forward and backward Fourier transforms and the models. PROPACK [START_REF]PROPACK[END_REF] is used to calculate the truncated SVD of the dictionary employing an iterative Lanczos method.

4.A. Object discrimination and location: illustration with holograms of digits

We illustrate, using simulations, the capability of the proposed method to correctly identify and locate several digits by matching diffraction patterns directly on digit holograms. We aim to find the minimum rank of dictionary resulting in correct recognition and location of all digits.

First, we illustrate our reconstruction results of the hologram of Fig. 1-a which contains overlapping and non-overlapping digit patterns. Inspired by the decay rate of the singular values of the dictionary (see the green curve in Fig. 3), we verified that a rank 5 approximation of the dictionary results in correct object recognition and depth estimation for non-overlapping objects (i.e., the digits 5 and 0 in Fig. 1-b were recognized and located within pixel accuracy). However, the overlapping patterns (i.e., patterns of objects 1,8 and 6 in Fig. 1-b) are not completely recognizable using only the first 5 modes. Our experiment showed that at least 16 modes are necessary to obtain an approximation of the dictionary that provides enough information for all the digits to be recognized and located in the 3D volume. In this case, given the relatively small size of the dictionary, a time gain of a factor 3 was observed compared to matching the full dictionary. The setup parameters are chosen to simulate a hologram captured by a 400 × 400-pixel camera with the pixel size of 20 µm and the fill-factor of 0.7. the laser wavelength is set to 0.532 µm. The depth position of digits is chosen to be in the range of [0.55 0.85] m. The dictionary contained 300 patterns, 30 Singular values of digit dictionary Singular values of spherical particles dictionary Singular values of Fresnel dictionary Singular values of digit dictionary with rotation Fig. 3. First 100 singular values of the dictionary in descending order. The dictionary of digits (in green) contains 300 centered patterns for digits placed at different depth positions with the application parameters in Sec. 4.A. The dictionary of spherical particles (in red) consists of the diffraction models for spherical objects changing depth and radius with the application parameters in Sec. 4.B. The Fresnel dictionary (in blue) corresponds to Fresnel functions with changing depth values. The dictionary of rotated digits (in mustard) contains 900 centered patterns (10 (digits)×90 (angle samples) ). In the case of spherical objects, 5 first modes can approximate the dictionary accurately. The vanishing rate of the singular values of the Fresnel dictionary is very low which implies that it is not possible to deduce an accurate low rank approximation of Fresnel dictionary. patterns for different depth values per digit class. The dimensions of this relatively small dictionary is therefore 160000 × 300.

Fig. 4-a shows the hologram of digits and Fig. 4-b shows the residuals after cleaning the signature of the objects. The model that is used to clean the patterns contains the second order term of the object signature. However, the inter-pattern second order terms (i.e., the interferences between the object diffractions) are not known, since all the objects are not yet detected when the signature of one object is being removed (these interferences belong to the signatures of 1,8 and 6). These interferences (with very low intensities) are not cleaned from the residuals and are visible in Fig. 4-b.

Next, we study the discrimination power of the method as a function of the noise level and the approximation rank. The goal is to (i) show the effect of noise magnitude on the required rank of dictionary for correct object recognition and location, and (ii) to estimate the minimum required rank for successful discrimination of every digit from the others. The time costs calculated in these experiments are finally used to evaluate an order of magnitude for the time gain of the proposed algorithm in this application. Our simulations contain a set of 10000 holograms of digits. 10 noise level are selected and for each noise level and for each digit 100 holograms are simulated. Changing the rank of the approximation, object recognition is performed to calculate the rate of error of classification. One hologram and the results for one digit (a "2") are shown in Fig. 5. As shown in Fig. 5-(a), a rank 4 approximation of dictionary results in correct recognition of "2" for every level of noise. The object is placed at 0.17 m and one of the holograms with noise of σ = 0.01 is shown inFig. 5-(b). The setup parameters are the same as the previous study except for the depth range which is changed to [0.15 0.2] m.

Using the same configuration, the second goal is to illustrate the discrimination error of every digit as a function of approximation rank. Fig. 5-(c) shows the error rate of object recognition for all 10 digits for the noise with the standard deviation of 0.31. This plot shows that a rank 6 approximation of the dictionary can correctly classify all digits without error. This plot can be used in conjunction with the plot of singular values (see the green curve in Fig. 3) to choose the approximation rank of dictionary for an assured object recognition task. One of the holograms containing a '2' is shown in Fig. 5-(d).

Finally, we performed a primitive experiment to test the robustness of object identifica- tion to 2D orientation changes of the objects. Our dictionary contained 90 orientations (i.e., sample step size of 4 degrees) and 10 digits. The resulted dictionary contained 900 diffraction patterns. Truncated SVD method was then employed to obtain the closest rank 500 approximation of this dictionary. The first 100 singular values of this dictionary are plotted in Fig. 3. For the Monte Carlo study, 100 holograms of a digit "2" with an arbitrary orientation angle were simulated and degraded by white and Gaussian noise with standard deviation 0.3. Our results showed that a rank 12 approximation of the dictionary was accurate enough to estimate the orientation angle and recognize the object accurately. Fig. 6-a shows the object 2 that is rotated counter-clock wise by 210.5 degrees. This rotation value falls between two dictionary samples of 208 and 212 degrees. Fig. 6-b shows one hologram of such an object placed at depth of 0.75 m. This hologram is degraded by white and Gaussian noise with the standard deviation of 0.31. Fig. 6-c shows the error rates of the object identification and orientation estimation for this object. A full study of robustness to rotation and scale variations is left as future work.

4.B. Simulations of spherical particle holograms

In this section, we show the results of our proposed method on the reconstruction of spherical object holograms. Monte Carlo simulations are performed to illustrate the accuracy of object sizing and location, and the time gain. This experiment is performed to check the biases and standard deviations of the estimated radii and 3D position using the low rank approximation of the dictionary. The results show standard deviations of the same order as the Cramer-Rao lower bounds (CRLBs) on the estimation of the radii and 3D position (see [START_REF] Fournier | On the single point resolution of on-axis digital holography[END_REF] for details on the calculation of CRLBs) and negligible bias using the first 50 modes, which validates the use of an approximated dictionary in the exhaustive search step of the greedy algorithm (see Sec. 2).

The details of this experiment is as followed: 50 holograms are simulated for a particle placed at the center of the hologram and another set of 50 holograms are simulated for a particle at the bottom right corner of the sensor, with the depth distance of 0.2 m and radii of 80 µm. The SNR of the signal is set to approximately 19. The camera is considered to have 512 × 512 pixels with the pixel size of 7 µm and fill-factor of 0.7. The laser wavelength is 0.532 µm.

To estimate the bias and standard deviation of the proposed method, we compared our results with the reference location and size obtained using a huge dictionary with very fine sampling. To build that reference dictionary, we sampled the parameters with step sizes equal to the CRLBs (which are an order of 15 µm for z and 50 nm for radii), obtaining the dictionary C1 corresponding to 616 (z,r) pairs. The search space with this dictionary is limited in our experiments to [0.1999 0.2001] m for z and [79 81] µm for the radii to avoid having a huge number of samples.

In calculating the mode dictionary for one setup, wide ranges of parameters can be considered. Such mode dictionary can be used in different experiments . In our experiment, the mode dictionary is calculated for the range of parameters z as [0.1 0.3] m and the range of radii as [60 100] µm. To sample these ranges, a trade-off must be found between (i) the of the coarseness of the sampling and (ii) the dimensions of the dictionary. In our experiments, we chose the sampling steps as 500 times the CRLBs for z and 50 times the CRLBs for r, that are experimentally shown to satisfy both limits. The general dictionary C2 has 600 (z, r) pairs. The first 100 modes of the dictionary are calculated using the truncated SVD algorithm implemented in PROPACK. The obtained mode dictionary is then used to approximate the application dictionary C1 . It has been verified that the first 50 modes are enough to have estimations with standard deviations of the same order as the CRLBs and negligible bias. Fig. 7 illustrates the evolution of the mean square difference between diffraction patterns, (i.e., the opposite of their correlation) around its global minimum increasing the number of modes. According to this illustration, the curvature of the main basin with about 10 modes is accurate enough for a coarse detection of the object's coordinates.

To show the efficiency of using the proposed method, we compare the time costs of the exhaustive search step of the original greedy method to the proposed method. The exhaustive search takes more than one minute and half to obtain the results using the full dictionary. Applying the proposed method, the first 5 modes are showed to be enough for accurate estimation of parameters. The time costs of the exhaustive search using such an approximated dictionary takes less than 10 seconds. The second step of the greedy algorithm (i.e., the local optimization step) is not changed in the proposed method and it takes in both cases less than 10 seconds. The time gain of using the proposed method is then a factor of 8.

It shall be noted that in the case of light back-propagation methods using the Fresnel dictionary [START_REF] Kreis | Handbook of Holographic Interferometry: Optical and Digital Methods[END_REF] (i.e., a dictionary which contains only the Fresnel transforms and not the diffraction patterns), SVD does not necessarily speedup the reconstruction. The blue curve in Fig. 3 shows the singular values of the Fresnel dictionary in a descending order. As shown in this illustration, the vanishing rate of the singular values is very low, resulting in very inaccurate low rank approximations of the dictionary. Eckart-Young theorem states that the best rank-k approximation of the matrix is the one found by SVD. One can therefore conclude that it is not possible to have an accurate low rank approximate of the Fresnel dictionary. This is not surprising as the 3D capability of holography comes from the difference between diffraction patterns produced by point-like objects located at different depths (i.e., Fresnel transform is close to invertible). Extended objects affect the Fresnel transform as relative low-pass filters that decrease the significance of high frequencies. These filtering operations (a) (b) (c) Fig. 7. With an increased number of modes, the approximation of the dictionary improves and matching a reference diffraction pattern against the approximated dictionary produces a sharper correlation peak leading to more accurate 3D location. Plotted curves represent the mean square difference between diffraction patterns, i.e., the opposite of their correlation. Reddest curves correspond to the most accurate approximations of the dictionary (using up to 50 modes). The pattern of the dictionary minimizing the mean square difference is in best match with the reference diffraction-pattern and gives the 3D location and diameter of the detected particle. introduce correlations between the columns of the dictionary which result in fast vanishing singular values (see the red curve in Fig. 3). Accurate low-rank approximations can be obtained in these cases.

4.C. Experimental Holograms of spherical particles

The last set of experiments involves the reconstruction of captured holograms using an approximated dictionary. Our experiments show that a rank 5 dictionary is accurate enough to detect and successfully clean all the particles placed in the field of view of camera (see Fig. 8 for a captured hologram and the cleaned one). The time gain of the proposed method is then a factor of 8. These holograms (see Fig. 8-a) were captured of injected water droplets at the Department of Fluid Mechanics and Acoustics of Lyon (LMFA). The droplets were generated by a piezoelectric jetting device manufactured by MicroFab Technologies. This injector produces mono-dispersed droplets with radii of 31 µm ± 0.5 µm. The droplets were produced to be at distances ranging from 30 cm to 48 cm of a 1024 × 1280 pixel camera with pixel size of 21.7 µm and fill-factor of 0.84. This imaging setup had a magnification of 1.42 and the holograms were captured with an SNR close to 16. The search space contained 410 pairs of (z, r).

It should be noted that the detection of out of field particles [START_REF] Soulez | Inverse problem approach in particle digital holography: out-of-field particle detection made possible[END_REF] is very limited using the first modes of dictionary. The signatures of the out-of-the-field objects contain only the highest frequencies and these high frequencies are not present in the first modes of the dictionary. A future work can be directed towards finding a rigorous method of mode selection for out-of-field detection.

Conclusion

Signal processing methods can perform the object recognition task directly from the digital holograms. To reduce the computational costs of such approaches, we have proposed to obtain a lower dimensional space that captures the inter-pattern variability of the objects' holographic signatures. The calculation of the algorithm complexity shows that the time gain of this approach depends on the application parameters and is large when considering huge dictionaries of diffraction-patterns.

Our Monte Carlo studies on the application to spherical objects have shown that the parameter estimation using a low-rank approximation of the dictionary is unbiased with a standard deviation comparable to the Cramer-Rao lower bounds. We have also shown that very low rank approximations give results which are accurate enough for the greedy object reconstruction method to successfully find the accurate object parameters. To this end, we could achieve a speedup rate of 8 on our video of experimental holograms of water droplets.

We have shown in a primitive experiment the robustness of our proposed method to orientation changes. A more comprehensive study on the potential of the proposed approach to detect the orientation, scale and 3D location of objects from their diffraction patterns is a perspective of our work with potential interesting applications.

Fig. 1 .

 1 Fig. 1. Illustration of a digit recognition task from a digital hologram: (a) a simulation of a 400×400-pixels hologram of digits located at different distances (pixel size: 20 µm, fill-factor: 0.7, laser wavelength: 0.532 µm, depth range: [0.55 0.85] m); (b) Fresnel reconstruction of the volume displaying both infocus and out-of-focus images of the digits which makes digit recognition from the reconstructed volume a difficult task. Direct recognition of the diffraction patterns using the approach described in the paper avoids dealing with these numerous image artifacts.

Fig. 2 .

 2 Fig. 2. Identification of the best-matching pattern in the dictionary: (a) by correlation with each element of the dictionary, straightforward but costly; or (b) by linear combination of the correlation with few diffraction-pattern modes, much faster. Computation of correlation maps is denoted using the symbol ' * '.

Fig. 4 .

 4 Fig. 4. Digit recognition from a digital hologram: (a) the hologram of digits (b) the residuals after cleaning all the signatures from the hologram.

Fig. 5 .

 5 Fig. 5. Study of the error rate on object recognition for the toy problem of holograms of digits (see Sec. 4.A). (a) shows the error rates of object recognition for a '2' placed at 0.17m from the sensor. The curves show the error rate for different approximation ranks and different noise levels, (b) illustrates one hologram of 2 degraded by white Gaussian noise with σ = 0.01, (c) plots the error rates representing the discriminating power of approximated dictionaries for a fixed level of noise resulting in a SNR ≈ 3, (d) illustrates one of the holograms of Sec. 4.A for object '2' placed at 0.17m from the sensor.

Fig. 6 .

 6 Fig.6. Study of the robustness of our method to the orientation factor (see Sec. 4.A). (a) shows an object '2' rotated counter-clock wise on the object plane by 210.5 degrees, (b) illustrates one of the simulated holograms for object '2' placed at depth of 0.75m. This hologram is degraded by white and Gaussian noise with σ = 0.31, (c) shows the error rates of object recognition and orientation estimation for this object.

Fig. 8 .

 8 Fig. 8. (a): Experimental hologram of droplets from LMFA Lyon used in Sec. 4.C, (b): same experimental hologram cleaned from the in-field particles using the 5 first modes. The magnitude of the residuals is high due to the signature of the out-of-field particles placed close to the borders.


Yannick Bonthonneau

A NOTE ON THE RESONANCE COUNTING FUNCTION FOR SURFACES WITH CUSPS

We prove sharp upper bounds for the number of resonances in boxes of size 1 at high frequency for the Laplacian on finite volume surfaces with hyperbolic cusps. As a corollary, we obtain a Weyl asymptotic for the number of resonances in balls of size T → ∞ with remainder O(T 3/2 ).

In this short note, we intend to prove sharp bounds on resonance-counting functions for the Laplacian on finite volume surfaces with hyperbolic cusps. Let M be a complete non-compact surface, equipped with a Riemannian metric g. We assume that (M, g) can be decomposed as the union of a compact manifold with boundary and a finite number of hyperbolic cusps, each one being isometric to (a, +∞) y × S 1 θ with metric

dy 2 + dθ 2 y 2
for some a > 0. The spectral properties of the Laplacian ∆ g were first studied by Selberg [Sel89] and Lax-Phillips [START_REF] Lax | Scattering theory for automorphic functions[END_REF] in constant negative curvature, and by Colinde-Verdière [START_REF] Verdière | Pseudo-laplaciens[END_REF], Müller [START_REF] Müller | Spectral geometry and scattering theory for certain complete surfaces of finite volume[END_REF], Parnovski [START_REF] Parnovski | Spectral asymptotics of Laplace operators on surfaces with cusps[END_REF] in the non-constant curvature setting.

On such surfaces, the resolvent R(s) = (∆ g -s(1 -s)) -1 of the Laplacian admits a meromorphic extension from {ℜs > 1/2} to C as an operator mapping L 2 comp to L 2 loc and the natural discrete spectral set for ∆ g is the set of poles denoted by

R ⊂ {s ∈ C| ℜs ≤ 1/2} ∪ (1/2, 1].
The poles are called resonances and are counted with multiplicity m(s) (the multiplicity m(s) is defined below and corresponds, for all but finitely many resonances, to the rank of the residue of the resolvent at s). We shall recall in the next section how the set of resonances is built. To study their distribution in the complex plane, we define two counting functions :

N R (T ) := s∈R,|s-1/2|≤T m(s) (1) N R (T, δ) := s∈R,|s-1/2-ıT |≤δ m(s). (2)
The first result on the resonance counting function was proved by Selberg [Sel89,p. 25] for the special case of hyperbolic surfaces with finite volume : the following Weyl type asymptotic expansion holds as T → ∞

(3) N R (T ) = Vol(M) 2π T 2 + C 0 T log(T ) + C 1 T + O T log(T )
for some explicit constants C 0 , C 1 . In variable curvature, Müller gives a Weyl asymptotic [Mül92, Th. 1.3.a] of the form

N R (T ) = Vol(M) 2π T 2 + o(T 2 )
and this was improved by Parnovski [START_REF] Parnovski | Spectral asymptotics of Laplace operators on surfaces with cusps[END_REF] who showed that for all ǫ > 0

(4) N R (T ) = Vol(M) 2π T 2 + O(T 3/2+ǫ ).
Parnovski's proof relies on a Weyl type asymptotic expansion involving the scattering phase S(T ) (see next section for a precise definition) :

(5)

2πN d (T ) + S(T ) = Vol(M) 2 T 2 -2kT ln T + O(T ),
where k is the number of cusps, and N d is the counting function for the L 2 eigenvalues of ∆ g embedded in the continuous spectrum.

Using a Poisson formula proved by Müller [START_REF] Müller | Spectral geometry and scattering theory for certain complete surfaces of finite volume[END_REF] and estimate (5), we are able to improve the results of Parnovski :

Theorem 1. For T > 1, and 0 ≤ δ ≤ T /2, the following estimates hold

N R (T, δ) = O(T δ + T ), (6) N R (T ) = Vol(M) 2π T 2 + O(T 3/2 ). (7)
In the first estimate with δ = 1, the exponent in T is sharp in general, as can be seen from Selberg's result (3) which implies that there is

C > 0 such that as T → ∞ (8) N(T, 1) = CT + O T log T
In n-dimensional Euclidan scattering, upper bounds O(T n-1 ) on the number of resonances in boxes of fixed size at frequency T were obtained by Petkov-Zworski [PZ99] using Breit-Wigner approximation and the scattering phase ; our scheme of proof is inspired from their approach. Their result was extended to the case of noncompact perturbations of the Laplacian by Bony [START_REF] Bony | Résonances dans des domaines de taille h[END_REF]. In general, it is expected that the number of resonances in such boxes is controlled by the (fractal) dimension of the trapped set (see for example Zworski [Zwo99], Guillopé-Lin-Zworski [START_REF] Guillopé | The Selberg zeta function for convex co-compact Schottky groups[END_REF],

Sjöstrand-Zworski [SZ07], Datchev-Dyatlov [DD13]).
Acknowledgement. We thank M. Zworski for his suggestion which shortened significantly the argument of proof. We also thank J-F. Bony for sending us his work, and Colin Guillarmou and Nalini Anantharaman for their fruitful advice.

Preliminaries

We start by recalling well-known facts on scattering theory on surfaces with cusps, and we refer to the article of Müller [START_REF] Müller | Spectral geometry and scattering theory for certain complete surfaces of finite volume[END_REF] for details. Let (M, g) be a complete Riemannian surface that can be decomposed as follows:

M = M 0 ∪ Z 1 ∪ • • • ∪ Z k ,
where M 0 is a compact surface with smooth boundary, and Z j are hyperbolic cusps

Z j ≃ (a j , +∞) × S 1 , j = 1 . . . k,
with a j > 0 and the metric on Z j in coordinates (y, θ) ∈ (a j , +∞) × S 1 is

ds 2 = dy 2 + dθ 2 y 2 .
Notice that the surface has finite volume when equipped with this metric.

The non-negative Laplacian ∆ acting on C ∞ 0 (M) functions has a unique self-adjoint extension to L 2 (M) and its spectrum consists of (1) Absolutely continuous spectrum σ ac = [1/4, +∞) with multiplicity k (the number of cusps). (2) Discrete spectrum

σ d = {λ 0 = 0 < λ 1 ≤ • • • ≤ λ i ≤ . . . }, possibly finite,
and which may contain embedded eigenvalues in the continuous spectrum. To λ ∈ σ d , we associate a family of orthogonal eigenfunctions that generate its eigenspace (

u i λ ) i=1...d λ ∈ L 2 (M) ∩ C ∞ (M)
. The generalized eigenfunctions associated to the absolutely continuous spectrum are the Eisenstein functions, (E j (x, s)) i=1...k . Each E j is a meromorphic family (in s) of smooth functions on M. Its poles are contained in the open half-plane {ℜs < 1/2} or in (1/2, 1]. The Eisenstein functions are characterized by two properties :

(1) ∆ g E j (., s) = s(1 -s)E j (., s) (2) In the cusp Z i , i = 1 . . . k, the zeroth Fourier coefficient of E j in the θ variable equals δ ij y s i + φ ij (s)y 1-s i where y i denotes the y coordinate in the cusp Z i and φ ij (s) is a meromorphic function of s.

We can collect the scattering coefficients φ ij in a meromorphic family of matrices, φ(s) = (φ ij ) ij called scattering matrix. We denote its determinant by ϕ(s) = det φ(s). Then the following identities hold

φ(s)φ(1 -s) = Id, φ(s) = φ(s), φ(s) * = φ(s).
The line ℜs = 1/2 corresponds to the continuous spectrum. On that line, φ(s) is unitary, ϕ(s) has modulus 1. We also define the scattering phase ( 9)

S(T ) = - T 0 ϕ ′ ϕ ( 1 2 + ıt)dt
The set of poles of ϕ, φ and (E j ) j=1...k is the same, we call them them scattering poles and we shall denote Λ this set. It is contained in {ℜs < 1/2} ∪ (1/2, 1]. The union of this set with the set of s ∈ C such that s(1 -s) is an L 2 eigenvalue, is called the resonance set, and denoted R. Following [Mül92, pp.287], the multiplicities m(s) are defined as :

(1) If ℜs ≥ 1/2, s = 1/2, m(s) is the dimension of ker L 2 (∆ g -s(1 -s)). (2) If ℜs < 1/2, m(s) is the dimension of ker L 2 (∆ g -s(1 -s)) minus the order of ϕ at s. (3) m(1/2) equals (Tr(φ(1/2)) + k)/2 plus twice the dimension of ker L 2 (∆ g -1/4).
For convenience, we define two counting functions for the discrete spectrum and the poles of ϕ:

N d (T ) := |s i -1/2|≤T m(s i ), (10)
N Λ (T ) := s∈Λ,|s-1/2|≤T m(s), (11) so that N R (T ) := s∈R,|s-1/2|≤T m(s) = 2N d (T ) + N Λ (T ).

Main observation

In this Section, we explain how to obtain estimate for N R (T ) in boxes at high frequency.

From the asymptotic expansion (5), we deduce that for 0 ≤ δ ≤ T /2,

S ′ (T ) = log 1 q + ρ∈Λ 1 -2ℜρ (ℜρ -1/2) 2 + (ℑρ -T ) 2 .
where q is some constant (not necessarily < 1). Let C > 1, 0 < ǫ < 1 and Ω T,δ := {s ∈ C; |s -1/2 -ıT | ≤ δ/C and 0 ≤ 1/2 -ℜs ≤ ǫδ}.

Then, for s ∈ Ω T,δ , [T -δ,T +δ] 1 -2ℜs (ℜs -1/2) 2 + (t -ℑs) 2 dt = 2 arctan t -ℑs 1/2 -ℜs T +δ T -δ
The addition formula for arctan, with x, y > 0 and xy > 1 is given by arctan x + arctan y = π + arctan

x + y 1 -xy

thus [T -δ,T +δ] 1 -2ℜs (ℜs -1/2) 2 + (t -ℑs) 2 dt = 2π -2 arctan 2δ(1/2 -ℜs) δ 2 -|s -1/2 -ıT | 2 ≥ 2π -2 arctan Cǫ,
where C is set to be 2/(1 -1/C 2 ). For ǫ small enough, this is bigger than, say π.

Since all but a finite number of terms in (13) are positive, we have :

S(T + δ) -S(T -δ) ≥ O(δ) + ρ∈Λ∩Ω T,δ π.
Combining with (12), we deduce that

N d (T + δ) -N d (T -δ) + #Λ ∩ Ω T,δ = O(T δ) + O(T ) + O(δ).
This is the content of (6) in our main theorem.

Consequence

Now, we proceed to prove the second part of our theorem. We will follow the method of Müller [Mül92, pp. 282], which is a global and quantitative version of the argument used in the previous section. Integrating the Poisson formula over [-T, T ], we relate the scattering phase asymptotics to the poles of φ. Using the arctan addition formula, we are left with the sum of N Λ (T ) and an expression with arctan's (equation (4.9) in [START_REF] Müller | Spectral geometry and scattering theory for certain complete surfaces of finite volume[END_REF]) :

(

  12) 2π(N d (T +δ)-N d (T -δ))+S(T +δ)-S(T -δ) = 2Vol(M)T δ -4kδ ln T +O(T ). Next, we recall the Poisson formula for resonances proved by Müller [Mül92, Th. 3.32] (13)

  Johannes Sjöstrand and Maciej Zworski. Fractal upper bounds on the density of semiclassical resonances. Duke Math. J., 137(3):381-459, 2007. [Zwo99] Maciej Zworski. Dimension of the limit set and the density of resonances for convex cocompact hyperbolic surfaces. Invent. Math., 136(2):353-409, 1999. E-mail address: yannick.bonthonneau@ens.fr DMA, U.M.R. 8553 CNRS, École Normale Superieure, 45 rue d'Ulm, 75230 Paris cedex 05, France

The sum is then split between {1} the poles in {|T -|ρ -1/2|| > T 1/2 }, and {2}, the others. Müller proved that the sum {1} is O(T 3/2 ). The sum {2} can be bounded by 1 4

From [START_REF] Müller | Spectral geometry and scattering theory for certain complete surfaces of finite volume[END_REF]Cor. 3.29], we also recall that

If ñ(T ) is the counting function for Λ, we deduce that

This concludes the proof.