A method for eliminating off-axis digital holographic zero-order terms

Abstract

The application relates to a method for eliminating off-axis digital holographic zero-order terms, which comprises the following steps: obtaining holographic data of off-axis digital holography, and obtaining spectrum information through Fourier transform. Pixel points with intensity higher than a certain value in the spectrum are selected as initial data objects, the initial data objects are divided into clustering categories through a K-means clustering method, and clustering centers are calculated. The clustering centers are taken as initial mean values of GMM, the covariance of each cluster is calculated, the mixing coefficient of GMM is calculated according to the proportion of data points in each cluster, the posterior probability of each data point belonging to each Gaussian distribution is calculated through Bayes theorem, and the GMM model is updated, so that the final clustering result is obtained. The clustering result located at the center position of the spectrum is eliminated, and the zero-order spectrum information is eliminated. Compared with the prior art, the application has the advantages of high accuracy and strong universality.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of digital holography, and particularly to a method for eliminating zero-order term of off-axis digital holography. BACKGROUND

[0002] Digital holography has been widely used in quantitative imaging of living cells, microscopic particle tracking, MEMS detection, etc. According to the tilt angle between the object beam and the reference beam, digital holography can be divided into on-axis, off-axis and micro-off-axis digital holography. Off-axis digital holography reconstruction usually involves a spatial-spectral filtering process. In the Fourier spectrum domain of the off-axis hologram, both the positive first-order spectrum and the negative first-order spectrum carry the spatial frequency of the object, and their centers are separated from the center of the zero-order spectrum. By filtering and intercepting either of the positive first-order spectrum and the negative first-order spectrum, the reconstruction imaging of the off-axis hologram can be achieved. However, due to the existence of the zero-order spectrum, the filtering and intercepting region of the first-order spectrum is limited, which means that the effective spatial bandwidth for reconstructing the complex amplitude is limited, which will affect the resolution of the reconstructed image. However, if the filtering and intercepting region of the spectrum is too large, it may contain some zero-order spectrum, which will cause some stripe-like noise in the reconstructed image.

[0003] In optical imaging and digital image processing, the effective spatial bandwidth product of the system is an important parameter affecting the imaging quality. Some multi-frame methods can improve the spatial bandwidth product of the system, thereby improving the imaging resolution or field of view size. However, these methods are not always suitable for off-axis digital holography reconstruction. For off-axis digital holography imaging, the resolution and quality of the reconstructed image can be significantly improved by suppressing the zero-order spectrum. In online digital holography technology, the phase-shifting method can effectively eliminate the zero-order term and the conjugate term of the hologram. Since the phase-shifting method usually requires a precise displacement device or a complex motion estimation algorithm, it is not used in off-axis digital holography reconstruction. By adjusting the mirror on the object arm, a phase shift can be introduced in the object wave, thereby suppressing the zero-order term in the hologram, but this method requires at least two holograms to be acquired, which is not suitable for real-time imaging. Using a nonlinear filter can effectively suppress the zero-order term and the conjugate term of the off-axis hologram, but requires the amplitude of the reference wave to be much larger than the amplitude of the object wave. In summary, it is necessary to design a digital holography zero-order term elimination method that is efficient and has high universality. SUMMARY

[0004] The purpose of the present application is to overcome the low universality and low accuracy of the prior art and provide a method for eliminating the zero-order term of off-axis digital holography.

[0005] The purpose of the present application can be achieved by the following technical solutions:

[0006] A method for eliminating the zero-order term of off-axis digital holography, comprising the following steps:

[0007] S1: running off-axis digital holographic optical system, obtaining the off-axis hologram to be processed;

[0008] S2: Fourier transform is carried out on the off-axis hologram to be processed, and a frequency spectrum of the off-axis hologram is obtained;

[0009] S3: according to the frequency spectrum, selecting a pixel point with intensity higher than a threshold value as an initial data object;

[0010] S4: selecting a center point of the frequency spectrum as an initial clustering center of the zero-order spectrum, and randomly selecting two points as a clustering center of the positive first-order spectrum and a clustering center of the negative first-order spectrum respectively;

[0011] S5: for each pixel point in the initial data object, calculating the Euclidean distance from it to all clustering centers, and classifying it to the nearest clustering center; recalculating the mean of all data points in each cluster to update the center of each cluster;

[0012] S6: repeating S5 until the clustering center no longer changes significantly or the set iteration number is reached;

[0013] S7: using the clustering center in the result of S6 as the initial mean of the GMM, calculating the distribution of each clustering data point according to the clustering result in S6, calculating the covariance of each cluster, calculating the proportion of the number of data points contained in each cluster to the total number of data points, and initializing the mixing coefficient of the GMM;

[0014] S8: calculating the posterior probability of each data point belonging to each Gaussian distribution by Bayes' theorem, recalculating the mean, covariance matrix and mixing coefficient of each Gaussian distribution according to the posterior probability, and maximizing the likelihood function of the Gaussian mixture model;

[0015] S9: repeating S8 until the parameters of the model no longer change significantly, or the increment of the likelihood function is less than the preset threshold value;

[0016] S10: the GMM calculates the probability of each pixel point belonging to each cluster, thereby obtaining the final clustering result, setting all the clustering results located at the center position of the frequency spectrum as 0, and obtaining the final result of eliminating the zero-order spectrum information.

[0017] Further, the parameter expression of the Gaussian mixture model is:

[0018]

[0019] In the formula, k is the number of the mixture model, α k is the mixing weight, 0<α k <1, is the Gaussian model of the kth mixture component, μ k and denote the mean and covariance matrix.

[0020] Further, the method further comprises the following step: S11: extracting positive first-order spectrum information and centering processing, and reconstructing the reconstruction intensity and phase of the object to be measured based on the angular spectrum method.

[0021] Further, the calculation expression of the angular spectrum method is:

[0022]

[0023] I(x,y,z) = |U(x,y,z)| 2

[0024]

[0025] In the formula, and denote the Fourier transform and inverse Fourier transform, respectively, is a filter for extracting +1 or -1 order spectrum, G(f x ,f y ,z) is the transfer function of the angular spectrum diffraction method, f x and f y are frequency domain coordinates, λ is the wavelength of light; z is the reconstruction distance, I(x,y,z) and are the reconstruction intensity and phase of the object to be measured, respectively.

[0026] Further, the calculation expression of the Euclidean distance is:

[0027]

[0028] In the formula, v i denotes the sample, center denotes the cluster center coordinates, and n denotes the observation value.

[0029] Further, the calculation expression of the mixture coefficient is:

[0030]

[0031] In the formula, k is the number of mixture models, α k is the mixture weight, 0 < α k <1, is the Gaussian model of the kth mixture component, μ k and denote the mean and covariance matrix.

[0032] Further, the calculation expression of the posterior probability is:

[0033]

[0034] wherein, α k is a mixing weight, is the kth mixture component Gaussian model,

[0035] Further, the calculation expression of the mean value, covariance matrix and mixing coefficient of each Gaussian distribution is recalculated as:

[0036]

[0037]

[0038] wherein, α k is a mixing weight, N represents the sample number, γ jk represents the posterior probability of the ith sample belonging to the kth cluster, μ k is the mean value of the kth Gaussian distribution, represents the center position of the cluster, X n is the data point of the nth sample, represents the covariance matrix of the kth Gaussian distribution, describes the expansion direction and shape of the cluster, (X n -μ k )(X n -μ k ) T is the outer product of the difference between the nth sample and the mean value, represents the dispersion degree of the sample and the mean value.

[0039] Further, the calculation expression of the likelihood function is:

[0040]

[0041] wherein, L(θ|X) is the likelihood function value, N represents the total sample number, k represents the number of Gaussian distributions, α k is a mixing weight, represents the probability density function of the kth Gaussian distribution.

[0042] Further, the calculation expression of the updated cluster center is:

[0043] Center n = ((∑(v i )) / p)

[0044] wherein, v i represents a sample, and p represents the sample number in the class after clustering.

[0045] Compared with the prior art, the present application has the following beneficial effects:

[0046] 1) This invention eliminates zero-order terms by combining K-means and GMM. After K-means clustering, GMM calculates the posterior probability of the Gaussian distribution to further subdivide the differences between the zero-order spectrum and other spectra, thus obtaining accurate clustering results. Especially when the spectrum distribution is complex, it can accurately eliminate the interference of the zero-order spectrum.

[0047] 2) When processing large amounts of holographic data, the combination of K-means and GMM can reduce the need for manual adjustment of filter parameters through automated clustering and modeling. Attached Figure Description

[0048] Figure 1 This is a schematic diagram of the off-axis digital holographic system structure in a specific embodiment.

[0049] Figure 2 This is a flowchart of a method for eliminating off-axis digital holographic zero-level terms provided in an embodiment of the present invention.

[0050] Figure 3 This is the original hologram collected in the embodiments of the present invention.

[0051] Figure 4 This is the spectrum of the original hologram after K-means clustering in this embodiment of the invention.

[0052] Figure 5 This is a spectral distribution diagram after GMM clustering in an embodiment of the present invention.

[0053] Figure 6 This is the spectrum diagram after eliminating the zero-order spectrum in this embodiment of the invention.

[0054] Figure 7 This is the phase map reconstructed by the angular spectrum method in an embodiment of the present invention.

[0055] The markings in the diagram are as follows: 1. Light source, 2. Beam expander and collimator, 3. Beam splitter, 4. First reflecting mirror, 5. Second reflecting mirror, 6. Measurement object, 7. Beam splitter, 8. CCD camera, 9. Computer. Detailed Implementation

[0056] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. These embodiments are based on the technical solution of the present invention and provide detailed implementation methods and specific operating procedures. However, the scope of protection of the present invention is not limited to the following embodiments.

[0057] Example 1

[0058] This invention provides a method for eliminating off-axis digital holographic zero-order terms, such as... Figure 2 As shown.

[0059] Specifically comprising the following steps:

[0060] S1: running an off-axis digital holographic optical system, obtaining an off-axis hologram to be processed;

[0061] The experimental environment of the present application is shown in Figure 1 , Figure 1 The device structure diagram of off-axis digital holography. The light source 1 is a helium-neon laser with a wavelength λ equal to 632.8nm. The light beam emitted by the light source 1 is divided into two beams by the beam splitter 3 after being expanded and collimated by the beam expander 2. One beam is the object light, which is reflected by the second mirror 5 and irradiates the object 6. One beam is the reference light, which is reflected by the first mirror 4 and combined with the object light on the beam splitter 7. The interference angle is adjusted by adjusting the first mirror 4 or the second mirror 5. Finally, the object light and the reference light interfere on the CCD camera 9 to form an off-axis hologram. The pixel size of the CCD camera is 4.65um×4.65um. The digitized hologram is input into the computer 10.

[0062] S2: Fourier transform is performed on the off-axis hologram to be processed to obtain the frequency spectrum information and the frequency spectrum matrix of the off-axis hologram;

[0063] S3: According to the frequency spectrum information, the pixel points with intensity higher than the threshold value are selected as the initial data objects;

[0064] S4: The center point of the frequency spectrum matrix is selected as the initial clustering center of the zero-order spectrum, and two points are randomly selected as the clustering center of the positive first-order spectrum and the clustering center of the negative first-order spectrum;

[0065] S5: For each pixel point in the initial data object, the Euclidean distance from it to all clustering centers is calculated, and it is classified into the nearest clustering center. Then for each cluster, according to the mean value of all pixel points classified into this class, the center is reselected;

[0066] The calculation expression of the Euclidean distance is:

[0067]

[0068] In the formula, v i represents the sample, center represents the clustering center coordinate, and n represents the observation value.

[0069] The sample is divided into the corresponding cluster based on the minimum distance principle mindis(n)=1,2,…k, k represents the number of clusters. The clustering center is updated Center n = ((∑(v i )) / p), in which v i represents the sample, and p represents the number of samples in the class after clustering.

[0070] S6: Repeat S5 until the cluster centers no longer change significantly or a set number of iterations is reached;

[0071] S7: Take the cluster centers obtained in S6 as the initial means of the Gaussian mixture model, calculate the covariance of each cluster, and calculate the mixing coefficients of the Gaussian mixture model according to the proportion of data points in each cluster; the calculation expression of the mixing coefficient is:

[0072]

[0073] In the formula, k is the number of the mixture model, α k is the mixing weight, 0 < α k <1, is the Gaussian model of the kth mixture component, μ k and represent the mean and the covariance matrix.

[0074] S8: Calculate the posterior probability of each data point belonging to each Gaussian distribution by Bayes' theorem, and recalculate the mean, covariance matrix and mixing coefficient of each Gaussian distribution according to the posterior probability, and maximize the likelihood function of the Gaussian mixture model;

[0075] Based on the given value, the calculation expression of the posterior probability is:

[0076]

[0077] In the formula, α k is the mixing weight, is the Gaussian model of the kth mixture component,

[0078] The calculation expression of the recalculated mean, covariance matrix and mixing coefficient of each Gaussian distribution is:

[0079]

[0080]

[0081] In the formula, α k is the mixing weight, N represents the sample number, γ jk represents the posterior probability of the ith sample belonging to the kth cluster. μ k is the mean of the kth Gaussian distribution, representing the center position of the cluster, X n is the data point of the nth sample. represents the covariance matrix of the kth Gaussian distribution, describing the expansion direction and shape of the cluster. (X n -μ k )(X n -μ k )T is the outer product of the difference between the nth sample and the mean, representing the degree of dispersion of the sample from the mean.

[0082] S9: repeat S8 until the parameters of the model no longer change significantly, or the increment of the likelihood function is less than a preset threshold value;

[0083] The calculation expression of the likelihood function is:

[0084]

[0085] In the formula, L(θ|X) is the value of the likelihood function, N represents the total number of samples, k represents the number of Gaussian distributions, α k is the mixing weight, represents the probability density function of the kth Gaussian distribution.

[0086] S10: The Gaussian mixture model calculates the probability of each pixel belonging to each cluster, thereby obtaining the final clustering result. The clustering results located at the center position of the spectrum are all set to 0, and the final result of eliminating the zero-order spectrum information is obtained.

[0087] In the present example, the measured object is a resolution test target. The obtained original hologram is as shown in Figure 3 . Fourier transform is performed on the original hologram to obtain the frequency information of the original hologram, that is, the holographic image spectrum diagram is obtained. In the spectrum diagram, the pixel points with an intensity higher than 0.003 times the maximum spectrum intensity are selected as initial data objects. The coordinate positions (u, v) of these pixel points are selected as the input data of clustering, and each data object contains two dimensions, that is, two-dimensional coordinates. The center point of the spectrum diagram is selected as the initial clustering center of the zero-order spectrum, and the other two clustering centers are randomly selected. The Euclidean distance of each data point to all clustering centers is calculated, and it is classified into the nearest clustering center. For each cluster, the center is recalculated until the clustering center no longer changes significantly or the set number of iterations is reached. The result of K-means preliminary clustering is as shown in Figure 4 . As can be seen from Figure 4 , K-means can successfully distinguish three types of spectrum information in the spectrum. The clustering centers in the K-means result are used as the initial mean of GMM. According to the clustering result of K-means, the distribution of data points in each cluster is calculated to obtain the covariance of each cluster. According to the proportion of the number of data points contained in each cluster in the total number of data points, the mixing coefficient of GMM is initialized. The posterior probability of each data point belonging to each Gaussian distribution is calculated by Bayes' theorem. According to the posterior probability, the mean, covariance matrix and mixing coefficient of each Gaussian distribution are recalculated. The result of the final clustering of GMM is as shown in Figure 5As shown in the figure, the GMM divides the spectrum into three Gaussian distribution spectrum information, and then sets the clustering data at the center position to 0, and finally obtains the spectrum information after eliminating the zero-order term as shown in the figure Figure 6 Then, the phase of the resolution plate is obtained based on the angular spectrum method reconstruction.

[0088] In order to verify the results of the present application, the following steps are also included:

[0089] S11: Extracting positive first-order spectrum information and centering processing, and obtaining the reconstruction intensity and phase of the object to be measured based on the angular spectrum method reconstruction.

[0090] The calculation expression of the angular spectrum method is:

[0091]

[0092] I(x,y,z) = |U(x,y,z)| 2

[0093]

[0094] In the formula, and respectively represent the Fourier transform and the inverse Fourier transform. is a filter for extracting +1 or -1 order spectrum, G(f x ,f y ,z) is the transfer function of the angular spectrum diffraction method, f x and f y are frequency domain coordinates, λ is the wavelength of light; z is the reconstruction distance, I(x,y,z) and are the reconstruction intensity and phase of the object to be measured, respectively.

[0095] The reconstructed image of the present embodiment is shown in the figure Figure 7 The partial enlarged view shows the 8th and 9th groups of elements in the resolution plate, and it can be seen that the reconstructed phase image obtained by the present application has the advantages of clear contour, no artifacts and high resolution. From the experimental results, the present method better restores the information of the object light, proving the effectiveness of the present application.

[0096] The preferred embodiments of the present application are described in detail above. It should be understood that those skilled in the art can make many modifications and changes to the present application without creative labor, according to the concept of the present application. Therefore, any technical solutions obtained by logical analysis, reasoning or limited experiments based on the prior art according to the concept of the present application shall be within the protection scope determined by the claims.

Claims

1. A method of eliminating off-axis zero-order terms of a digital hologram, characterized by, The method comprises the following steps: S1: running an off-axis digital holographic optical system to obtain an off-axis hologram to be processed; S2: performing Fourier transform on the off-axis hologram to be processed to obtain a frequency spectrum of the off-axis hologram; S3: selecting, according to the frequency spectrum, pixel points with intensity higher than a threshold value as initial data objects; S4: selecting a center point of the frequency spectrum as an initial clustering center of zero-order frequency spectrum, and randomly selecting two points as a clustering center of positive first-order frequency spectrum and a clustering center of negative first-order frequency spectrum respectively; S5: for each pixel point in the initial data objects, calculating the Euclidean distance from the pixel point to all clustering centers, and classifying the pixel point to the nearest clustering center; recalculating the mean value of all data points in each cluster to update the center of each cluster; S6: repeating S5 until the clustering centers no longer change significantly or a set number of iterations is reached; S7: using the clustering centers in the result of S6 as initial means of the GMM, calculating the distribution of data points in each cluster according to the clustering result in S6, calculating the covariance of each cluster, calculating the proportion of the number of data points contained in each cluster to the total number of data points, and initializing the mixing coefficients of the GMM; S8: calculating the posterior probability of each data point belonging to each Gaussian distribution by Bayes' theorem, recalculating the mean value, covariance matrix and mixing coefficient of each Gaussian distribution according to the posterior probability, and maximizing the likelihood function of the Gaussian mixture model; S9: repeating S8 until the parameters of the model no longer change significantly, or the increment of the likelihood function is less than a preset threshold value; S10: the GMM calculates the probability of each pixel point belonging to each cluster, thereby obtaining the final clustering result, setting all the clustering results located at the center position of the frequency spectrum as 0, and obtaining the final result of eliminating zero-order frequency spectrum information.

2. The method of claim 1, wherein, The parameter expression of the Gaussian mixture model is: wherein is the number of mixture models, is the mixture weight, 0 <1, =1, is the Gaussian model of the th mixture component, and denote the mean and covariance matrix.

3. The method of claim 1, wherein, The method further comprises the following steps: S11: extracting positive first-order frequency spectrum information and performing centering processing, and reconstructing the intensity and phase of the object to be measured based on the angular spectrum method.

4. The method of claim 3, wherein, The calculation expression of the angular spectrum method is: wherein and denote Fourier transform and inverse Fourier transform, respectively, is a filter extracting +1 or -1 order spectrum, is a transfer function of angular spectrum diffraction method, and are frequency domain coordinates, , λ is the wavelength of light; z is the reconstruction distance, and are the reconstructed intensity and phase of the object to be measured, respectively.

5. The method of claim 1, wherein, The calculation expression of the Euclidean distance is: In the formula, denotes the sample, denotes the cluster center coordinates, n denotes the observation.

6. The method of claim 1, wherein, The calculation expression of the posterior probability is: wherein is a mixing weight, is the Gaussian model of the th mixture component, = 1.

7. The method of claim 1, wherein, The calculation expression of the recalculated mean value, covariance matrix and mixing coefficient of each Gaussian distribution is: wherein, is a mixing weight, N represents the number of samples, represents the posterior probability that the ith sample belongs to the kth cluster, is the mean of the kth Gaussian distribution, representing the center position of the cluster, is the data point of the nth sample, represents the covariance matrix of the kth Gaussian distribution, describing the expansion direction and shape of the cluster, is the outer product of the difference between the nth sample and the mean, representing the degree of dispersion of the sample from the mean.

8. The method of claim 1, wherein, The calculation expression of the likelihood function is: In the formula, is a likelihood function value, N represents a total number of samples, k represents a number of Gaussian distributions, is a mixing weight, represents a probability density function of the kth Gaussian distribution.

9. The method of claim 1, wherein, The calculation expression of the updated clustering center is: wherein denotes the samples, denotes the number of samples in the cluster.

Images