Single-shot mid-infrared incoherent holography using Lucy-Richardson-Rosen algorithm

In recent years, there has been a significant transformation in the field of incoherent imaging with new possibilities of compressing three-dimensional (3D) information into a two-dimensional intensity distribution without two-beam interference (TBI). Most of the incoherent 3D imagers without TBI are based on scattering by a random phase mask exhibiting sharp autocorrelation and low cross-correlation along the depth. Consequently, during reconstruction, high lateral and axial resolutions are obtained. Imaging based on scattering requires an astronomical photon budget and is therefore precluded in many power-sensitive applications. In this study, a proof-of-concept 3D imaging method without TBI using deterministic fields has been demonstrated. A new reconstruction method called the Lucy-Richardson-Rosen algorithm has been developed for this imaging concept. We believe that the proposed approach will cause a paradigm-shift in the current state-of-the-art incoherent imaging, fluorescence microscopy, mid-infrared fingerprinting, astronomical imaging, and fast object recognition applications.


Section 1: Theoretical analysis
The phase of a Cassegrain objective lens (COL) is approximated by an annular diffractive lens or diffractive equivalent COL (DE-COL) at a single wavelength λ. In addition to that COL is mounted with a cross shaped block as shown in Fig.  S1(a), which modulates the transmittance. The phase image of the COL can be approximated as shown in Fig. S1(b). A simplified optical configuration of imaging is shown in Fig. S1(c). The complex amplitude introduced by the COL can be approximated (Fresnel approximation) 1 as: where M 1 is the transmittance function of the cross shaped block, ; R 1 and R 2 are the inner and outer radii of the annulus and f is the focal length, which is designed for finite conjugate mode given as .
The theoretical analysis is from the object plane (o-plane) located at u' from DE-COL to the sensor plane (s-plane) at v from DE-COL. The proposed system is a spatially incoherent imaging system. A self-luminous point object located at emits light, which reaches the DE-COL with an intensity where u = u'+Δz and Δz is the axial shift error. The complex amplitude reaching the DE-COL plane can be given as , where C 1 is a complex constant, , and are the linear and quadratic phase factors. The complex amplitude after the DE-COL is given as . As COL does not have notable spectral aberrations, the wavelength dependent analysis is not considered. The intensity pattern obtained at the sensor plane located at v from the DE-COL is given as: ⊗ where ' ' is the 2D convolutional operator. The above equation can be expressed as: A two-dimensional chemical object p located in the object plane o consisting of M points is given as: where is the transmitted intensity. Every δ point generates an intensity pattern given as with a shift from the optical axis depending upon the acquired linear phase. The intensity distribution obtained for the object is given as: where the transverse magnification M T =(v/u). There are two cases: Δz = 0 and Δz ≠ 0. When Δz = 0, direct imaging condition is satisfied and , which is a magnified version of the object the with minimum feature given

210006-S2
by 1.22λv/D on the camera, where D is the diameter of the DE-COL. When Δz ≠ 0, , which is a distorted image of the object formed by the convolution of distorted PSF with p. As it is seen here, unlike random field based sharp autocorrelation and low cross-correlation along depth (SALCAD), deterministic SALCADs can have dual mode, i.e., both direct imaging and indirect imaging can co-exist. In any thick sample, the planes within the depth of focus ±2(u') 2 λ/D 2 can be observed without the need for any reconstruction, which is different from SALCADs based on random fields.
Section 2: Simulation of focal characteristics of DE-COL A simulative study of the DE-COL was carried out and the images of the intensity distributions obtained in the sensor plane for different values of shift errors for a regular diffractive lens are shown in Fig. S2. It demonstrates that even though the recorded intensity distribution is not a point, the autocorrelation is sharp which is the resolving power in indirect imaging mode. The cross-correlation for all images was carried out with respect to the reference image Δz (ref), except for the two planes that the cross-correlation is lower for other planes. The approaches for reconstruction by LRRA and NLR are quite different from one another. If and p are the point spread function and object function respectively, then the object intensity distribution for a linear optical system in intensity is given as . However, in practical cases, the expression is not always true due to noise σ. The noise σ can be signal dependent Poisson noise or additive noise or both. The correct expression for the intensity distribution for an object is given as . For this reason, the correlation by the matched filter, which is exactly the opposite operation of convolution does not yield the optimal solution 5 . refers to the complex conjugate of and the loop is iterated until an optimal reconstruction is obtained. In fact, the LRA has been widely used for astronomical imaging, where the recorded image is distorted or blurred. In many cases, the blurred image is not very different from the original image unlike scattering based images, where the object lntensity at the sensor plane

210006-S3
I n R ⊗ I PSF I p I PSF information is converted into speckles. For this reason, LRA's initial guess is often the recorded image itself and the final solution is a maximum-likelihood solution. As seen in the above equation, there is a forward convolution and the ratio between this and is correlated with , which is replaced by the NLR and yields a better estimation. Consequently, the process achieves a rapid convergence.
Here, a test object "Lucy Richardson Rosen" has been selected as shown in Fig. S3(a). The image of the deterministic PSF generated by COL is shown in Fig. S3(b). The distorted image of the test object is shown in Fig. S3(c). The reconstruction results using LRA (iterations = 500), NLR (α = 0, β = 0.5) and LRRA (α = 0, β = 0.6, iterations = 10) are shown in Fig. S3(d-f), respectively. The LRRA is not only more than 50 times faster than LRA, but also a significantly better estimate than both LRA and NLR.

Section 4: Synthesis of PSFs from recorded PSF
In most of the studies of scattering-based 3D imagers, it was necessary to record the PSFs at all possible axial locations mainly due to the fact that they are not deterministic 6,7 . Some studies had utilized the linear region of propagation to apply the scaling factor to synthesize the PSFs from one or two recorded PSFs 8,9 . However, this linear region is quite short. The above disadvantage does not exist with deterministic fields, where the modulation function can be generated using one or two recorded PSFs based on phase-retrieval algorithms. The schematic of the modified phase retrieval algorithm is shown in Fig. S4 10 . Once the phase is synthesized in plane -1, the complex amplitude can be propagated by any distance and the entire focal characteristics can be obtained.

Section 5. Data structure conversion
The infrared microspectrometry unit (IRM) and the Fourier transform infrared (FTIR) spectrometer are linked by OPUS software of Bruker. The output from the OPUS software is saved as data point table format (*.dpt). The spectral images (64 × 64) for 765 channels obtained from the IRM FTIR system are structured into a matrix size of 765 × 4097.
The spectral image data (64 × 64) is obtained from every row of the matrix 1∶4096 by rearrangement. The resulting cube data is of the structure (765 × 64 × 64). A single matrix is noisy and so multiple images (50 images) are averaged to obtain an image with a high signal to noise ratio. The image obtained from visible light is shown in Fig. S5(a). The im-  age of a single recording and average of 100 images are shown in Fig. S5(b) and S5(c) respectively. A MATLAB code is provided has been designed for the reformatting of data. The code is deposited online and can be downloaded here (10.5281/zenodo.5541384).