An automatic method for assembling a large synthetic aperture digital hologram

A major issue so far for digital holography is the low spatial resolution generally achieved. The numerical aperture is limited by the area of currently available detectors, such as CCD sensors, which is significantly lower than that of a holographic plate. This is an even more severe constraint when IR sensors such as microbolometers are taken into account. In order to increase the numerical aperture of such systems, we developed an automatic technique which is capable of recording several holograms and of stitching them together, obtaining a digital hologram with a synthetic but larger numerical aperture. In this way we show that more detail can be resolved and a wider parallax angle can be achieved. The method is demonstrated for visible as well IR digital holography, recording and displaying large size objects. ©2012 Optical Society of America OCIS codes: (090.1995) Digital holography; (100.6640) Superresolution; (100.0100) Image processing; (090.2870) Holographic display; (110.3080) Infrared imaging. References and links 1. E. Allaria, S. Brugioni, S. De Nicola, P. Ferraro, S. Grilli, and M. Meucci, “Digital holography at 10.6 μm,” Opt. Commun. 215(4-6), 257–262 (2003). 2. A. Faridian, D. Hopp, G. Pedrini, U. Eigenthaler, M. Hirscher, and W. Osten, “Nanoscale imaging using deep ultraviolet digital holographic microscopy,” Opt. Express 18(13), 14159–14164 (2010). 3. L. Miccio, A. Finizio, S. Grilli, V. Vespini, M. Paturzo, S. De Nicola, and P. Ferraro, “Tunable liquid microlens arrays in electrode-less configuration and their accurate characterization by interference microscopy,” Opt. Express 17(4), 2487–2499 (2009). 4. M. Paturzo, P. Memmolo, A. Finizio, R. Näsänen, T. J. Naughton, and P. Ferraro, “Synthesis and display of dynamic holographic 3D scenes with real-world objects,” Opt. Express 18(9), 8806–8815 (2010). 5. S. A. Alexandrov, T. R. Hillman, T. Gutzler, and D. D. Sampson, “Synthetic aperture fourier holographic optical microscopy,” Phys. Rev. Lett. 97(16), 168102 (2006). 6. Y. Kuznetsova, A. Neumann, and S. R. Brueck, “Imaging interferometric microscopy-approaching the linear systems limits of optical resolution,” Opt. Express 15(11), 6651–6663 (2007). 7. V. Mico, Z. Zalevsky, P. García-Martínez, and J. García, “Superresolved imaging in digital holography by superposition of tilted wavefronts,” Appl. Opt. 45(5), 822–828 (2006). 8. L. Martínez-León and B. Javidi, “Synthetic aperture single-exposure on-axis digital holography,” Opt. Express 16(1), 161–169 (2008). 9. F. Le Clerc, M. Gross, and L. Collot, “Synthetic-aperture experiment in the visible with on-axis digital heterodyne holography,” Opt. Lett. 26(20), 1550–1552 (2001). 10. R. Binet, J. Colineau, and J.-C. Lehureau, “Short-range synthetic aperture imaging at 633 nm by digital holography,” Appl. Opt. 41(23), 4775–4782 (2002). 11. B. M. Hennelly, T. J. Naughton, and J. McDonald, “Digital holographic superresolution by rotating the object wavefield,” in Adaptive Optics: Analysis and Methods/Computational Optical Sensing and Imaging/Information Photonics/Signal Recovery and Synthesis Topical Meetings on CD-ROM, OSA Technical Digest (CD) (Optical Society of America, 2007), paper DTuD5. 12. C. Liu, Z. Liu, F. Bo, Y. Wang, and J. Zhu, “Super-resolution digital holographic imaging method,” Appl. Phys. Lett. 81(17), 3143–3145 (2002). 13. M. Paturzo, F. Merola, S. Grilli, S. De Nicola, A. Finizio, and P. Ferraro, “Super-resolution in digital holography by a two-dimensional dynamic phase grating,” Opt. Express 16(21), 17107–17118 (2008). 14. M. Paturzo and P. Ferraro, “Correct self-assembling of spatial frequencies in super-resolution synthetic aperture digital holography,” Opt. Lett. 34(23), 3650–3652 (2009). 15. J. H. Massig, “Digital off-axis holography with a synthetic aperture,” Opt. Lett. 27(24), 2179–2181 (2002). #158129 $15.00 USD Received 15 Nov 2011; revised 6 Dec 2011; accepted 7 Dec 2011; published 13 Feb 2012 (C) 2012 OSA 27 February 2012 / Vol. 20, No. 5 / OPTICS EXPRESS 4830 16. F. Gyímesi, Z. Füzessy, V. Borbély, B. Ráczkevi, G. Molnár, A. Czitrovszky, A. Tibor Nagy, G. Molnárka, A. Lotfi, A. Nagy, I. Harmati, and D. Szigethy, “Half-magnitude extensions of resolution and field of view in digital holography by scanning and magnification,” Appl. Opt. 48(31), 6026–6034 (2009). 17. V. Micó, L. Granero, Z. Zalevsky, and J. García, “Superresolved phase-shifting Gabor holography by CCD shift,” J. Opt. A, Pure Appl. Opt. 11(12), 125408 (2009). 18. A. Gruen, “Adaptive least square correlation: a powerful image matching technique,” S. Afr. J. Photogrammet. 14, 175–187 (1985). 19. D. Lowe, “Distinctive image features from scale-invariant keypoints,” Int. J. Comput. Vis. 60(2), 91–110 (2004). 20. H. Bay, A. Ess, T. Tuytelaars, and L. Van Gool, “Speeded-up robust features (SURF),” Comput. Vis. Image Underst. 110(3), 346–359 (2008). 21. M. Paturzo, A. Pelagotti, A. Finizio, L. Miccio, M. Locatelli, A. Gertrude, P. Poggi, R. Meucci, and P. Ferraro, “Optical reconstruction of digital holograms recorded at 10.6 microm: route for 3D imaging at long infrared wavelengths,” Opt. Lett. 35(12), 2112–2114 (2010), http://www.opticsinfobase.org/ol/abstract.cfm?URI=ol-3512-2112. 22. J. M. Fitzpatrick and D. L. G. Hill, and C. R. Maurer, Jr., “Image registration,” in Handbook of Medical Imaging—Volume 2, Medical Image Processing And Analysis, M. Sonka and J. M. Fitzpatrick, eds. (SPIE Press, 2000), Chap. 8. 23. B. Zitovà and J. Flusser, “Image registration methods: a survey,” Image Vis. Comput. 21(11), 977–1000 (2003). 24. F. Maes, D. Vandermeulen, and P. Suetens, “Medical image registration using mutual information,” Proc IEEE 91(10), 1699–1722 (2003). 25. F. Maes, A. Collignon, D. Vandermeulen, G. Marchal, and P. Suetens, “Multimodality image registration by maximization of mutual information,” IEEE Trans. Med. Imaging 16(2), 187–198 (1997). 26. A. Collignon, F. Maes, D. Delaere, D. Vandermeulen, P. Suetens, and G. Marchal, “Automated multimodality medical image registration using information theory,” in Proceedings XIVth International Conference on Information Processing in Medical Imaging—IPMI'95, Computational Imaging and Vision (Kluwer Academic, 1995), Vol. 3, pp. 263–274.


Introduction
In digital holography (DH), the film material used in analog holography is replaced by an electronic matrix detector capable of recording in a wider spectrum, from deep UV to far IR [1,2]. Together with many advantages, this implies also some drawbacks, both in microcopy configuration [3] and for 3D imaging and displaying of macroscopic objects [4]. The size of pixels restricts the field of view, while the limited size of the detector bounds the resolution of the reconstructed holographic image. Holograms in the visible range are mainly acquired with charged coupled devices (CCD) or complementary metal-oxide semiconductor (CMOS) cameras. The area of their sensor is generally no larger of some Mpixel. In the infrared range this limit is much higher since their size in generally less than 1 Mpixel.
In the study reported in this paper, we managed to increase the synthetic aperture of our system by scanning a larger hologram (the interference light field generated by the object and the reference waves in the hologram plane) with a small matrix detector and stitching together the single holograms acquired. The stitching is performed in an automatic way based on registration of the overlapping portions of neighboring holograms.
Registration is the determination of a geometrical transformation that aligns points in one picture with corresponding points in another picture. Often, the registration is performed manually or in a semi-automatic way, by a user iteratively setting the parameters of the geometrical transformation. However, this approach is time consuming and can give subjective results. In this paper an automatic registration technique, based on computation of the mutual information, is applied to obtain the synthetic hologram.
In the last years, different approaches have been tested to increase the NA of the optical system in order to get super-resolution. Alexandrov et al. were able to increase the resolution by rotating the sample and recording a digital hologram for each position in order to capture the diffraction field along different directions [5]. A different approach was proposed by Kuznetsova et al. who rotated the sample with respect to the optical axis in order to re-direct the rays scattered at wider angles into the aperture of the optical system, thus going beyond its diffraction limit [6]. Mico et al. proposed and demonstrated a method for enhancing the resolution of the aperture limited imaging systems based on the use of tilted illumination and common-path interferometric recording [7]. Martínez-León and Javidi translated the camera by a few microns in order to increase both the spatial resolution and the sampling in the recoding process [8]. Synthetic-aperture was gained with on-axis digital heterodyne holography by LeClerc et al. [9] while a combination of aperture synthesis and phase-shifting DH was developed by Binet et al. [10]. Hennely et al. superimposes different digital holograms of the same 3-D object, accomplished by rotating the input object wavefield either by rotation of the object (it is 2-D) or by rotation of a mirror that is placed between the object and the CCD [11].
Moreover, super-resolved images can be obtained simply using a 1D diffraction grating [12] or a 2D tuneable phase grating [13] that allows one to collect parts of the light diffracted from the object that otherwise would fall outside the CCD. Different digital holograms, according to the grating geometry, are spatially multiplexed onto the same CCD array, and the super-resolved images are obtained by the numerical reconstruction. Moreover, in Ref. [14] it is demonstrated that the spatial frequencies are naturally self-assembled in the reconstructed image plane when the NA is increased synthetically at its maximum extent of three times.
Finally, other methods were presented which are based on the shift of the CCD to different positions in order to improve the system's resolution [15][16][17]. However, only Gyímesi et al.
[16] obtain the synthetic aperture hologram by directly stitching the acquired holograms, while Massig and Mico et al. first transform from spatial to frequency domain and then aptly process the obtained images, by selecting only specific portions and so on, thus requiring some skilled user interaction. Gyímesi and coworkers apply a correlation approach; however, these methods, least-squares matching algorithms [18] or feature-based approaches [19,20] are not always applicable to align and stitch speckle holograms, due to the few distinctive patterns they present in the spatial domain, resulting in noisy images. By practical tests we demonstrated they will not converge to the correct solution. Moreover Gyímesi and coworkers specify that the obtained accuracy is about 1-2 pixels, which is similar to what is obtained manually and sometimes is not enough to increase the resolution, as shown later in this paper. On the contrary, here we propose a new automatic technique which is capable of stitching together several speckle holograms with the aim to construct a synthetic digital hologram, thus increasing considerable the numerical aperture and demonstrating that a wider parallax angle can be achieved. The method is applied for visible wavelengths as well as demonstrated for the first time, at best of our knowledge, in the far IR spectrum.

Experimental setup
As test objects we use different statuettes, of which holograms are acquired by means of two different laser sources, a visible and an infrared one [21]. In Fig. 1 a sketch of the recording setups is shown. The infrared laser is a CW CO 2 laser emitting at 10.6 μm. The laser operates on the fundamental TEM 00 mode at about 110 W of emission power. As visible laser we use a DPSS laser emitting at λ=0.532 μm. The beam is divided by a beam splitter (BS) into two beams with different propagation direction: an object beam, composed of the 80% of the incoming intensity, and a reference beam, composed of the remaining 20% of intensity. The object beam is expanded by a lens (L 2 ) in order to cover entirely the object surface. The reference beam is redirected toward the thermal camera and expanded by a lens (L 1 ), and then interferes with the light scattered by the object. The IR holograms are recorded by means of a thermal imaging camera (Miricle Thermoteknix 307k) with 640 × 480 pixels and 25 μm pixel pitch, while the visible holograms are acquired by a CCD camera with 1024 × 1024 pixels, 4.4 μm in size. The cameras were fixed on two translation stages (horizontal and vertical direction) in order to move them manually of relatively precise step values (shorter than the sensor dimensions) and collect a set of shifted and super-imposable images of the interferometric pattern obtained.

Description of the "image registration" algorithm
As already stated, an image registration algorithm aims at determining the parameters of a geometrical transformation mapping points in one picture, the reference (or target) one, with corresponding points in another picture, the sensed (or template) one [22,23]. There are several proposed methods and many ways of classify them. One of the most interesting classification of these algorithms is the one based on the basis, i.e. the set of points or features that are involved in the registration task [22]. According to this, registration methods can be classified in point-based methods, surface-based methods and intensity-based methods. In point-based methods, a certain number of corresponding positions are a priori identified, by manual selection or by external superimposition, over both the to-be-registered images; the registration process seeks the displacement between these points. The surface-based methods determine which are the corresponding areas (instead of points) in the two images and try to find the transformation that best aligns these surfaces; in this case, the appropriate zones need to be correctly identified. Finally, the intensity-based methods calculate the transformation using all pixel values, with no reference to distinctive points or surfaces, through an iterative optimization of a similarity measure computed over two images. One of the advantages of intensity-based algorithms is that it need less user interaction compared to other methods. Within this class, one of the most interesting methods is the one based on the maximization of the mutual information (MMI) [24], that has shown excellent results compared to previously proposed methods, such as point-based or surface-based, which can often fail if there is an inherently difference of the image structures and tone dynamics. Moreover, this method is automatic, and does not need any preprocessing.

MMI-based "image registration"
Mutual information (MI), a basic notion of information theory, represents a measure of the amount of information that one image contains about the other one. The MMI approach states that MI between two images is maximum when the images are correctly registered. Given two images X and Y related by the geometric transformation T α with parameters α such that the pixel p of X with intensity value x correspond to the pixel T α (p) of Y with intensity value y, their mutual information is 2 , where p XY (x, y) is the joint distribution, p X (x) and p Y (y) the marginal ones. The mutual information registration criterion affirms that the images are geometrically aligned by the transformation T α* where * arg max ( , ) Estimates for the joint and marginal distributions can be obtained by simple normalization of the joint and marginal histograms of the overlapping parts of both images [25]. The joint histogram h α (x, y) is obtained by binning the pixel intensity value pairs (X(p),Y(T α (p))) for all the pixels in the overlapping region of X and Y. Since very often the registered pixel position T α (p) will not coincide with a grid position, an interpolation of the reference image will be required to obtain the pixel value Y(T α (p)). Next, the following distributions can be estimated: , , , By using these values in Eq. (1) it is possible to derive the Mutual Information I(α), whose maximization will give us the optimal registration parameter α*.

Implementation issue
First of all, the implementation of the registration algorithm is strictly dependent on the geometrical transformation model we are referring to. In our case, we consider that the transformation the image can undergo is an affine transformation [22]. In general, an affine transformation is composed of linear transformations L (rotation, scaling or shear) and a translation (or "shift"): ' where the coordinate vector of a single point x  is transformed in another point ' x  by the transformation t  , that is the translational displacement, [ , ] x y t t t =  (5) and by a linear transformation that can be composed by a rotation R, with a rotation matrix (a 2 by 2 matrix) cos sen sen cos R and by a scaling S, with a diagonal matrix S = diag(s x , s y ), whose elements represent the scaling factors along the two axes. In the simpler condition of isotropic scaling, as we assume, s x = s y = s; thus, the previous equation, Eq. (4), can be simplified in the following: where s is a scalar value.
In our implementation, thus, the geometric parameters are α = [t x , t y , θ, s].
Introducing the θ and the s parameters allows us to cope also with some challenging occurrences, as reported by Gyímesi et al., showing how holograms can suffer from distortions which are excluded in the method they propose. The optimization process looking for the maximum of MI is here an heuristic search procedure, in which the 4 scalar parameters are iteratively changed by small amounts in order to find the maximum value. Another important issue to be taken into account is the fact that, in general, the modified pixel position T α (p) will not fall into a point of the pixel grid, so that an interpolation will be required to obtain the pixel intensity value Y(T α (p)). Several interpolation methods can be used. We implemented the trilinear (TRI) interpolation algorithm [26], that computes the intensity value of the interpolated pixel T α (p) as a weighted sum of the intensity values of the four nearest neighbors n i (i = 1, ⋅,4); the joint histogram h α is then incremented according to this computed value, as described:

Experimental results
In this section we show some of the experimental results obtained by applying the implemented MMI-based registration algorithm in order to stitch IR and visible holograms, in challenges cases. We acquired several neighboring holograms shifting the camera both vertically and horizontally in a meandering path, with an unknown number of overlapping pixels in both directions. No preprocessing was applied to the images so obtained. Since the developed algorithm accepts as input 2 images at a time, the images were firstly joined in couples and then the two resulting images were joined together. The whole procedure was fully automated.
In Fig. 2 we show the numerical reconstruction of a single (640 × 480 pixels) and a stitched hologram (1590 × 1500 pixels) of the Perseo statuette. The stitched hologram is obtained joining 4 × 3 single holograms. The resolution improvement in the reconstructed image is clearly visible. In order to evaluate the outcome of our algorithm with what is achievable with some others methods described in literature, we plotted the MI values obtained for a number of values of the t x and t y parameters, and we then compared their distribution with what achieved by replacing the MI with the Cross Correlation (CC) within the algorithm. In Fig. 3 we show the 3D plot of the normalized values of MI (a) and of CC (b) achieved as a result of a registration procedure, over a number of overlapping positions, obtained by changing the t x and t y parameters. In case the MI is applied, a clear best match appears, which corresponds to the optimal registration value; on the other hand, in case CC is applied, a number of local maxima, with close values are present, among which the global maximum is determined by random local noise.
In order to check whether the improvement was mainly due to the increase of number of pixel in the hologram, or to the proper assembly of the single holograms, we compared in Fig.  4 the results of the numerical reconstructions of a stitched hologram, obtained by correctly, i.e. by means of the proposed algorithm, superimposing the overlapping areas of 4  neighboring IR holograms of a "Perseus" statuette (a), and by joining them manually, i.e., by a skilled operator superimposing their overlapping areas using an image editing software tool (b). It clearly appears an improvement, in the level of detail which is achieved by using our algorithm, demonstrating that a not perfect alignment severely compromises the final output quality.
Although less severe, the same problem as far as resolution of the acquired hologram is also true for the visible range, when a digital camera is employed. In order to test the algorithm also in this case, we acquired, using a CCD camera, also 4 holograms in a 2×2 pattern both in horizontal and in vertical direction. One of the four acquired visible holograms is illustrated in Fig. 5a, while the joint hologram is illustrated in Fig. 5b. Figure 5c shows a picture of the object used in this case, that is the statuette of Adriano emperor. In Fig. 6 it is possible to compare the result of the digital reconstruction of a single hologram (Fig. 6a) and  the one of the joint hologram (Fig. 6b). A definite improvement is present in case a joint hologram is used. In principle there is no limit to the number of holograms which can be used by the algorithm. The main limit resides in the possibility of having an opportunely large hologram.
In the case of IR holograms, thanks to the high emission power of our CO 2 laser, we obtain a large hologram, so that many holograms can be acquired and stitched together. In fact, in Fig. 7 (a) it is shown the synthetic hologram obtained by means of the proposed software, stitching together 3 × 7 single holograms. In Fig. 7(b) and Fig. 7(c) its numerical reconstructions at two different distances, 515 mm and 505 mm, are shown. Looking at these images, it is clear that, thanks to the increased numerical aperture of the synthetic hologram, the depth of focus is decreased. In fact, in Fig. 7 (b) the inscription is more visible in respect to Fig. 7 (c) where, on the contrary, the neck and the ear are in good focus. If we look at the reconstruction of a single hologram at these two different distances, no focus differences appears, as shown in Fig. 8, because of the low numerical aperture and the consequent high depth of focus.
Another consequence of the increased numerical aperture is the possibility to display a wider parallax angle. In our case, the increasing in horizontal parallax is more evident, because the dimension of the stitched holograms is 3624 × 1468, 5.66 × 3.06 times the while the single one which size is 640 × 480. Media 1 shows the sequence of the numerical reconstruction of the left and right part of the stitched hologram. Then, seven 1920 × 1080 holograms has been extract from the stitched hologram and has been displayed sequentially by means of a spatial light modulator. The obtained optical reconstructions has been acquired by another CCD and are shown in Media 2 [21] (Fig. 9).

Conclusion
We presented an automatic technique capable of stitching together several speckle holograms, in order to obtain a larger synthetic hologram with higher numerical aperture. The algorithm allows high accuracy with limited computational costs. The method works well for digital hologram recorded both in the visible and in the IR range. Alongside with a perceptible improvement in the resolution, that allows to visualize more details, the numerical reconstruction of such hologram presents a decreased depth of focus. Moreover, a wider parallax angle can be displayed, and can be advantageously used with a spatial light modulator.