A synthetic aperture imaging method based on SIFT algorithm

Abstract

The application discloses a synthetic aperture imaging method based on a SIFT algorithm, relates to the field of synthetic aperture digital holographic image reconstruction technology, and provides a method for synthetic aperture imaging based on the SIFT algorithm while ensuring resolution, effectively improving the space-bandwidth product of a digital holographic imaging system, and expanding an imaging field of view. The method can be used in the case that the imaging field of view is limited, and according to actual measurement requirements or imaging field of view requirements, a plurality of off-axis holographic sub-images with a certain overlapping area are collected by adjusting a precision translation stage, then feature points of the reconstructed sub-images are extracted by using the SIFT algorithm for matching, and finally, the images are weighted and fused to obtain a large-field-of-view synthetic aperture image containing information of a measured object.

Description

Technical fields:

[0001] This invention relates to the field of synthetic aperture digital holographic image reconstruction technology, specifically to a synthetic aperture imaging method based on the SIFT algorithm. Background technology:

[0002] As humanity enters the information age, the proliferation of complex image, text, and video data necessitates more advanced display technologies to visualize data in a more natural and realistic way. While many 3D display technologies exist, digital holography can fully record the complex wavefront information of a sample, including the amplitude and phase information of the sample's diffracted waves, greatly increasing the amount of information recorded. It is a powerful label-free quantitative phase imaging technology. With the rapid development of electronic and computer technologies, large-field-of-view, high-bandwidth digital holography will become possible. Furthermore, the improvement in chip computing power will enable real-time recording and reproduction of digital holograms. Therefore, digital holography is increasingly attracting research interest, and its applications are expanding from imaging and display to deformation measurement, medical diagnosis, particle field measurement, and 3D image recognition.

[0003] Digital holographic imaging technology uses photoelectric sensors such as CCDs or CMOS to acquire the complex amplitude of transparent microscopic objects such as scratches on the surface of optical elements, cells, and tissues. Quantitative phase reconstruction of the specimen is achieved by measuring this complex amplitude. In 2015, Verma et al. combined digital holography with filtering techniques to measure the surface defect width of a 6-line / mm grating as 20.1 μm and the surface scratch width of an optical glass plate as 6.7 μm. In 2019, Trivedi et al. used single-beam lensless Fourier transform digital holography to detect the refractive index distribution of objects and its application in defect detection. Off-axis holography, in particular, achieves single-shot detection of the complex amplitude image through off-axis interferograms, separating the complex amplitude image of the measured object from its twin image in the frequency domain, thus reducing the difficulty of holographic image reconstruction.

[0004] The amount of information contained in an image acquired by an imaging system is related to the size of the spatial bandwidth product of the imaging system. High-resolution, wide-field imaging can detect more complete sample defects and provide a full-view image of cell morphology, greatly increasing the imaging volume and field of view. In microscopy, the spatial bandwidth product is defined as the product of the image field of view area and the spatial frequency band area. Because the number of pixels that modern microscope sensors can provide is far lower than the spatial bandwidth product of optical microscopes, much optical information of the object being measured cannot be detected. Microscopy still has considerable room for improvement in large-scale high-throughput imaging, which is of great significance for fields such as cell detection, drug development, and tissue biology. Although off-axis holographic interferograms have a wide bandwidth, the spatial bandwidth product of their complex amplitude images is relatively small, and the acquired holographic images contain limited information about the sample due to the limitations of the field of view of digital holographic microscopy systems and the photosensitive size of CCD cameras. Summary of the Invention:

[0005] To improve the spatial bandwidth product of digital holographic imaging, this invention proposes a synthetic aperture imaging method based on the SIFT (Scale-invariant feature transform) algorithm. This method records multiple holograms with overlapping regions after object translation in an off-axis holographic manner. The diffraction reconstruction process of the holograms is then simulated numerically using a computer to reconstruct the wavefront information of the object. The SIFT algorithm is used to extract and match feature points from the reconstructed images, achieving overall weighted fusion of the reconstructed sub-images, thereby obtaining a large field-of-view synthetic aperture image with a high spatial bandwidth product.

[0006] The technical problem to be solved by this invention is achieved by the following technical solution:

[0007] A synthetic aperture imaging method based on the SIFT algorithm includes the following steps:

[0008] Step 1: Image the object under test using a transmission-type off-axis holographic optical path system based on the Mach-Zehnder interferometer principle. A CCD camera acquires holograms. The object is fixed on a 3D translation stage. The CCD camera surface and the object surface must maintain a fixed distance and be parallel to each other, ensuring that the acquired holograms only have a two-dimensional translational transformation relationship. To ensure that the acquired hologram set covers the entire area of ​​the object under test, the acquisition path is as follows: The CCD camera first acquires a row of images at equal intervals from left to right in the horizontal direction. These images are named sequentially as I... 11 I 12 ... I 1b Then, the translation stage was adjusted to move the object under test down a certain distance, and the second row of images was acquired at equal intervals from right to left. The images were named I sequentially. 2b ... I22 I 21 Then, the object being measured is moved down a certain distance, and another row of images is acquired from left to right. This process is repeated until a × b holograms are acquired (a = 2). v , v=1, 2, 3,...; b=2 u (u = 1, 2, 3, ...), which is the holographic set of the object being measured, where a is the number of rows and b is the number of columns.

[0009] Step 2: Reconstruct the holographic atlas obtained in Step 1 frame by frame using sub-holographic images to obtain the reconstructed atlas of the object under test. This atlas contains a×b images, which are named sequentially as A. 11 A 12 A ab .

[0010] The method for reconstructing the holographic subgraph in step 2 is as follows:

[0011] S1. The interference intensity P(x,y) of the holographic image can be expressed as:

[0012] P(x,y)=|O(x,y)| 2 +|R(x,y)| 2 +O * (x,y)R(x,y)+O(x,y)R * (x,y) (0)

[0013] In formula (1), |O(x,y)| 2 and |R(x,y)| 2 These represent the intensity distributions of the corresponding object beam and reference beam, respectively, which are also the zeroth order of diffraction, O. * (x,y) and R * (x,y) are the conjugates of O(x,y) and R(x,y), respectively, and the third term O * (x,y)R(x,y) is a twin image, while O(x,y)R * (x,y) represents the reconstructed image required.

[0014] If the reconstruction distance is z, and the reconstructed light wave incident on the hologram is C(x,y), then the diffracted light field from the holographic surface to the reconstruction distance z is the object light field, and its complex amplitude distribution U z (x',y') can be represented as:

[0015]

[0016] In formula (2), F and F -1 These represent the two-dimensional Fourier transform and the inverse two-dimensional Fourier transform, respectively. x ,f y () represents the frequency domain coordinates of the reconstructed image.

[0017] S2. In order to extract the phase information ψ(x',y') of the reconstructed image, it is necessary to analyze the complex amplitude distribution U of the reconstructed image. z Numerical analysis is performed on (x',y'), and the imaginary part of its complex amplitude Im[U] is obtained. z (x',y') and its real part Re[U] z The value of the ratio [x', y'] can be obtained by taking its arctangent, and the expression is:

[0018]

[0019] Step 3: For the image A(x,y) to be stitched, first perform image denoising and image enhancement to improve the accuracy and stability of subsequent steps, and then extract feature points from the image to be stitched.

[0020] The specific steps of step 3 are as follows:

[0021] S1. Establish scale space: Perform scale transformation on the images to be stitched A(x,y) to construct the Gaussian difference scale space D(x,y,σ):

[0022] D(x,y,σ)=(G(x,y,kσ)-G(x,y,σ))*A(x,y) (4)

[0023] In formula (4), G(x,y,σ) is a Gaussian function, σ is a scale parameter, k is a scaling factor, and * indicates the convolution operation of the two functions.

[0024] S2. Feature point detection and matching information search of the image to be stitched: At each scale in the constructed Gaussian difference scale space, a difference Gaussian filter is used to detect local extrema, i.e., possible feature points. From the second-order Taylor expansion of the scale space function D(x,y,σ), we can see that the offset of the feature point position... for:

[0025]

[0026] when If the condition is met, then the feature point is removed; if the condition is met, then the feature point is removed. At the same time, it is also necessary to remove the feature point.

[0027] S3. Feature Point Description: The SIFT algorithm uses the image region surrounding the feature point to calculate the descriptor. To ensure its invariance to rotation and scaling transformations, it is necessary to calculate the gradient direction histogram in the neighborhood around the feature point. The gradient magnitude m(x,y) and gradient direction θ(x,y) of the pixel in the feature region are:

[0028]

[0029]

[0030] In formulas (6) and (7), L(x,y) is the spatial scale function.

[0031] Step 4: Use the SIFT algorithm to stitch the two images together: Extract feature points from the two images to be stitched using the method in Step 3, match them, calculate the transformation matrix between the two images, perform perspective transformation on the two images to be stitched, project them into the image space after stitching, and then average and weight the gray values ​​of the two images to be stitched together before merging them to achieve smooth stitching between the two images and obtain a new image.

[0032] Step 5: Group all the row images, with each row consisting of two adjacent sub-images (e.g., A in the first row). 11 With A 12 As a group, A 13 With A 14 (As a group), the images in each group of the row image set are stitched together using the method in step 4 until all row image sets are stitched together. At this point, each row has one stitched image. All the stitched images are then rearranged in order to form a new row image set W. i (i=1, 2, 3,..., a).

[0033] Step 6: The row atlas W obtained in step 5... i The images in the image set W are grouped, with two adjacent sub-images forming a group. The method in step 4 is then used to process the image set W. i Each set of images in the image set is stitched together until the image set W is reached. i All images are stitched together to obtain a single stitched image, which is the synthetic aperture map W(x,y) of the high spatial bandwidth product of the reconstructed object.

[0034] The beneficial effects of this invention are:

[0035] 1. This invention provides a better fusion method for the sub-images after holographic sub-image reconstruction, which can realize the fusion and stitching between reconstructed sub-images and obtain a smooth stitching result.

[0036] 2. This invention provides a method for synthetic aperture imaging based on the SIFT algorithm to effectively improve the spatial bandwidth product and expand the imaging field of view of a digital holographic imaging system while ensuring resolution. This method can acquire multiple off-axis holographic sub-images with a certain overlap area by adjusting a precision translation stage when the imaging field of view is limited, according to the actual measurement needs or imaging field of view requirements. Then, the SIFT algorithm is used to extract feature points of the reconstructed sub-images for matching, and then the images are weighted and fused to obtain a large field of view synthetic aperture map containing information of the measured object. Attached image description:

[0037] Figure 1 This is a flowchart of the synthetic aperture imaging method described in this invention;

[0038] Figure 2 This is a schematic diagram of the transmission-type digital off-axis holographic micro-interference system in this invention;

[0039] Among them, Laser is the laser, NF is the attenuator, BE is the beam expander, BS1 and BS2 are the beam splitters, M1 and M2 are the mirrors, P is the polarizer, MA is the object under test, MO1 and MO2 are the microscope objectives, and CCD is the camera that receives the digital holographic interferogram.

[0040] Figure 3 This invention utilizes a transmission-type digital off-axis holographic micro-interference system to acquire holographic atlases; wherein, Figure 3 a-3d represent holographic sub-images with certain overlapping areas on the surface of the measured object;

[0041] Figure 4 This invention utilizes a transmission-type digital off-axis holographic micro-interference system to acquire a holographic atlas, which is then reconstructed into a set of images corresponding to the holographic atlas. Figure 4 a-4d are respectively Figure 3 The reconstructed subgraph corresponding to the holographic subgraph;

[0042] Figure 5 To rebuild the child Figure 4 Feature point detection results for a and 4b;

[0043] Figure 6 for Figure 5 Synthetic aperture diagram;

[0044] Figure 7 The synthetic aperture maps obtained using the method of this invention were used to reconstruct atlases 4a-4d;

[0045] Figure 8 Synthetic aperture maps obtained using a phase-operation-based method were used to reconstruct atlases 4a-4d. Detailed implementation method:

[0046] To make the technical means, creative features, objectives and effects of this invention easier to understand, the invention will be further described below with reference to specific embodiments and illustrations.

[0047] Using a standard resolution plate as the object under test, the method provided by this invention is adopted to acquire off-axis holographic sub-images with a certain overlapping area by adjusting a three-dimensional translation stage. The reconstruction method of off-axis holographic images is analyzed from the wave theory of light. Then, the diffraction and reconstruction process of the hologram is simulated by computer numerical simulation to realize the reconstruction of the wavefront information of the object. The reconstructed sub-images are weighted and fused using the SIFT algorithm to obtain a synthetic aperture map with a high spatial bandwidth product.

[0048] This invention provides a synthetic aperture imaging method based on the SIFT algorithm, such as... Figure 1 As shown, the specific steps are as follows:

[0049] Step 1: Construct a transmission-type off-axis holographic optical path system based on the Mach-Zehnder interferometer principle, such as... Figure 2 As shown, a helium-neon polarized laser with a wavelength of 632.8 nm emits a polarized beam as the light source. After passing through the attenuator NF, the beam is collimated and amplified by the collimator BE to form an extended beam. The beam splitter BS1 splits the extended beam into two parts. One part of the beam serves as the object beam. After the polarization direction of the beam is adjusted by the polarizer P, the beam is redirected by the reflector M2 and transmitted through the object MA under test. The polarizer can also filter stray light and improve the imaging quality of the optical system. Then, the microscope objective MO1 amplifies the beam containing the surface information of the object under test and then incident it on the beam splitter BS2. The other part of the beam serves as the reference beam and continues to propagate forward. Then, the propagation direction is changed by the reflector M1 and it is incident on the microscope objective MO2. After being combined with the object beam at the beam splitter BS2, holographic interference occurs. The holographic interference image of the object under test is acquired by the CCD image sensor.

[0050] The object under test is fixed on a 3D translation stage. The surface of the CCD camera and the surface of the object under test must maintain a fixed distance and be parallel to each other to ensure that there is only a two-dimensional translational transformation relationship between the acquired holograms. To ensure that the acquired hologram set covers the entire area of ​​the object under test, the acquisition path is as follows: the CCD camera first acquires a row of images at equal intervals from left to right in the horizontal direction, and the images are named sequentially as I... 11 I 12 Then, the translation stage was adjusted to move the object under test down a certain distance, and the second row of images was acquired at equal intervals from right to left. The images were named I sequentially. 22 I 21 A total of 2×2 holograms were collected, such as Figure 3 As shown, this is a holographic atlas of the object being measured.

[0051] Step 2: Reconstruct the holographic atlas obtained in Step 1 frame by frame using sub-holographic images to obtain the reconstructed atlas of the object being measured, such as... Figure 4 As shown. This atlas contains 2×2 images, which are named sequentially as A. 11 A12 A 21 A 22 .

[0052] The method for reconstructing the holographic subgraph in step 2 is as follows:

[0053] S1. The interference intensity P(x,y) of the holographic image can be expressed as:

[0054] P(x,y)=|O(x,y)| 2 +|R(x,y)| 2 +O * (x,y)R(x,y)+O(x,y)R * (x,y) (0)

[0055] In formula (1), |O(x,y)| 2 and |R(x,y)| 2 These are the intensity distributions of the object beam and the reference beam, respectively, which are also the zeroth order of diffraction, O. * (x,y) and R * (x,y) are the conjugates of O(x,y) and R(x,y), respectively, and the third term O * (x,y)R(x,y) is a twin image, while O(x,y)R * (x,y) represents the desired reconstructed image. In off-axis digital holography, due to the angle θ between the object beam and the reference beam, the interference terms of the hologram are separated in the frequency domain after Fourier transform. If the deflection angle θ is large enough, the spectrum of the twin image can be completely separated. Then, the spectrum of the reconstructed image can be extracted using filtering techniques, and the off-axis hologram can be reconstructed using the angular spectrum method to obtain the reconstructed image.

[0056] If the reconstruction distance is z, and the reconstructed light wave illuminating the hologram is C(x,y), then the diffraction field from the holographic surface to the reconstruction distance z is the object light field, and its complex amplitude distribution U z (x',y') can be represented as:

[0057]

[0058] In formula (2), F and F -1 These represent the two-dimensional Fourier transform and the inverse two-dimensional Fourier transform, respectively. x ,f y () represents the frequency domain coordinates of the reconstructed image.

[0059] S2. In order to extract the phase information ψ(x',y') of the reconstructed image, it is necessary to analyze the complex amplitude distribution U of the reconstructed image. z Numerical analysis is performed on (x',y'), and the imaginary part of its complex amplitude Im[U] is obtained. z(x',y') and its real part Re[U] z The value of the ratio [x', y'] can be obtained by taking its arctangent, and the expression is:

[0060]

[0061] Step 3: For the image A(x,y) to be stitched, first perform image denoising and image enhancement to improve the accuracy and stability of subsequent steps. Then, extract feature points from the image to be stitched. The specific steps are as follows:

[0062] S1. Establish scale space: Perform scale transformation on the images to be stitched A(x,y) to construct the Gaussian difference scale space D(x,y,σ):

[0063] D(x,y,σ)=(G(x,y,kσ)-G(x,y,σ))*A(x,y) (4)

[0064] In formula (4), G(x,y,σ) is a Gaussian function, σ is a scale parameter, k is a scaling factor, and * indicates the convolution operation of the two functions.

[0065] S2. Feature point detection and matching information search of the image to be stitched: At each scale in the constructed Gaussian difference scale space, a difference Gaussian filter is used to detect local extrema, i.e., possible feature points. From the second-order Taylor expansion of the scale space function D(x,y,σ), we can see that the offset of the feature point position... for:

[0066]

[0067] when If the condition is met, then the feature point is removed; if the condition is met, then the feature point is removed. In such cases, the feature point also needs to be removed. The feature point detection results for images 4a and 4b to be stitched together are as follows: Figure 5 As shown.

[0068] S3. Feature Point Description: The SIFT algorithm uses the image region surrounding the feature point to calculate the descriptor. To ensure its invariance to rotation and scaling transformations, it is necessary to calculate the gradient direction histogram in the neighborhood around the feature point. The gradient magnitude m(x,y) and gradient direction θ(x,y) of the pixel in the feature region are:

[0069]

[0070]

[0071] In formulas (6) and (7), L(x,y) is the spatial scale function.

[0072] Step 4: Stitching the two images using the SIFT algorithm: Feature points of the two images to be stitched are extracted and matched according to the method in Step 3. The transformation matrix between the two images is calculated. Perspective transformation is performed on the two images, which are then projected into the stitched image space. The grayscale values ​​of the two images are averaged and weighted before merging to achieve a smooth stitching of the two images, resulting in a new image. The stitching result of images 4a and 4b is shown below. Figure 6 As shown.

[0073] Step 5: Group all the row images. In each row, two adjacent sub-images form a group. For example, in the first row, A 11 With A 12 As a group, A 13 With A 14 As a group, the images in each group of the row image set are stitched together using the method in step 4 until all row image sets are stitched together. At this point, each row has one stitched image. All the stitched images are then rearranged in order to form a new row image set W. i (i = 1, 2).

[0074] Step 6: The row atlas W obtained in step 5... i The images in the image set W are grouped, with two adjacent sub-images forming a group. The method in step 4 is then used to process the image set W. i Each set of images in the image set is stitched together until the image set W is reached. i All images are stitched together to obtain a single stitched image, which is the synthetic aperture map W(x,y) of the high spatial bandwidth product after the reconstruction of the measured object. Figure 7 As shown, the four reconstructed sub-images use a synthetic aperture map T(x,y) based on phase operations, as follows: Figure 8 As shown.

[0075] To further compare the stitching effect of our method with commonly used phase stitching methods, the gray-level variance values ​​of the synthesized aperture images W(x,y) and T(x,y) are calculated respectively: the gray-level variance of an image characterizes the average degree of gray-level variation; the greater the average degree of gray-level variation, the clearer the image; the smaller the average degree of gray-level variation, the blurrier the image. The average gray-level value of all pixels in the synthesized aperture image W(x,y) is calculated. for:

[0076]

[0077] N in formula (8) x and N y Let be the number of pixels in the synthetic aperture image W(x,y) along the x and y directions, respectively. The gray-level variance s of image W(x,y) is calculated from the average gray-level of the pixels in image W(x,y):

[0078]

[0079] The following formulas (8) and (9) can be used to calculate the result. Figure 7 The grayscale variance is 4983.06. Figure 8 The grayscale variance is 4406.17. This indicates that, compared to Figure 8 , Figure 7 The grayscale transformation is more even, resulting in clearer details in the image.

[0080] The foregoing has shown and described the basic principles, main features, and advantages of the present invention. Those skilled in the art should understand that the present invention is not limited to the above embodiments. The embodiments and descriptions in the specification are merely illustrative of the principles of the invention. Various changes and modifications can be made to the invention without departing from its spirit and scope, and all such changes and modifications fall within the scope of the present invention as claimed. The scope of protection of this invention is defined by the appended claims and their equivalents.

Claims

1. A synthetic aperture imaging method based on the SIFT algorithm, characterized in that, Includes the following steps: Step 1: Image the object under test using a transmission-type off-axis holographic optical path system based on the Mach-Zehnder interferometer principle. Use a CCD camera to acquire holograms. Fix the object under test on a three-dimensional translation stage. The surface of the CCD camera and the surface of the object under test must maintain a fixed distance and be parallel to each other. Acquire multiple off-axis holographic sub-images with a certain overlapping area to ensure that there is only a two-dimensional translation transformation relationship between the acquired holograms. Step 2: Reconstruct the holographic atlas obtained in Step 1 frame by frame using sub-holographic images to obtain the reconstructed atlas of the object under test. This atlas contains a×b images, which are named sequentially. , ... ; Step 3: For the images to be stitched together First, perform image denoising and image enhancement on the image to be stitched together, and then extract feature points from the image. Step 4: Use the SIFT algorithm to stitch the two images together: Extract the feature points of the two images to be stitched together using the method in Step 3, match them, calculate the transformation matrix between the two images, perform perspective transformation on the two images to be stitched together, project them into the image space after stitching, and then average and weight the gray values ​​of the two images to be stitched together before merging them to achieve smooth stitching between the two images to be stitched together, resulting in a new image. Step 5: Group all the row image sets, with each row containing two adjacent sub-images as a group. Use the method from Step 4 to stitch together each group of images in the row image set until all row image sets are stitched together. At this point, each row will have one stitched image. Rearrange all the stitched images in order to form a new row image set. (i=1, 2, 3, ..., a); Step 6: The row atlas obtained in step 5... The images in the atlas are grouped, with two adjacent sub-images forming a group, and the method in step 4 is used to process the atlas. Each set of images in the image set is stitched together until the image set is complete. All images are stitched together to obtain a single stitched image, which is the synthetic aperture map of the high spatial bandwidth product after the reconstruction of the object under test. .

2. The synthetic aperture imaging method according to claim 1, characterized in that, In step 1, to ensure that the acquired holographic atlas covers the entire area of ​​the object under test, the acquisition path is as follows: the CCD camera first acquires a row of images at equal intervals from left to right in the horizontal direction, and the images are named sequentially. , ... Then, the translation stage is adjusted to move the object under test down a certain distance, and the second row of images is acquired at equal intervals from right to left. The images are named sequentially. ... , Then, the object being measured is moved down a certain distance, and another row of images is acquired from left to right. This process is repeated until a × b holograms are acquired (a = 2). v v=1, 2, 3, ...; b=2 u (u=1,2,3,……), which is the holographic set of the object being measured, where a is the number of rows and b is the number of columns.

Images