A multi-plane spread spectrum iterative phase recovery method based on angular spectrum diffraction theory

Abstract

The application discloses a multi-plane spread spectrum iterative phase recovery method based on an angular spectrum diffraction theory, which comprises the following steps: acquiring light intensity distribution of a plane to be recovered and a plurality of propagation distance planes; preprocessing phase distribution of the plane to be recovered, combining the preprocessed phase distribution of the plane to be recovered with the light intensity distribution of the plane to be recovered to obtain initial complex amplitude of the plane to be recovered; iteratively processing the initial complex amplitude and the light intensity distribution on the plurality of propagation distance planes by using an angular spectrum-based spread spectrum iterative phase recovery algorithm to obtain preliminary recovery phase distribution of the plurality of planes to be recovered; fusing the preliminary recovery phase of the plurality of planes to be recovered according to a multi-plane phase spatial frequency domain fusion algorithm to obtain a fused phase spectrum; and obtaining a fused phase recovery result according to the fused phase spectrum. The angular spectrum inverse transfer function fused with a dynamic adjustment factor is used to suppress high-frequency noise and release high-frequency information; and the dynamic allocation of frequency domain weights is used to improve the noise resistance of phase recovery.

Description

Technical Field

[0001] This invention relates to the field of optical measurement and computational imaging technology, specifically to a multi-plane spread spectrum iterative phase recovery method based on angular spectrum diffraction theory. Background Technology

[0002] In image processing scenarios such as terahertz band compressed field measurement, laser wavefront detection, and microscopic imaging, phase information of a certain image plane may be lost due to various reasons, making it impossible to completely reconstruct the image. This may result in distorted reconstructed images, affecting their quality and usability. To solve the above problems, existing technologies typically employ the following methods for phase recovery, including:

[0003] Traditional iterative algorithms (such as the GS algorithm): rely on the Fourier transform transmission model, are only applicable to single-plane scenarios, and have poor noise resistance; angular spectrum iterative algorithms: although they adopt a strict diffraction model, the inverse transfer function will exponentially amplify noise at high spatial frequencies, requiring the use of a window function to truncate high frequencies, which will lead to the loss of phase recovery details; multi-plane fusion spectrum methods: existing multi-plane fusion methods do not consider the differential amplification effect of high-frequency noise at different propagation distances, and the fusion weight allocation is unreasonable. Summary of the Invention

[0004] To address the aforementioned problems, this invention provides a multi-plane spread spectrum iterative phase recovery method based on angular spectrum diffraction theory. This method collects light intensity distribution information from multiple planes to recover the phase distribution of the plane to be recovered, thereby achieving higher noise immunity and overcoming the problem of poor phase recovery performance in existing technologies.

[0005] This invention adopts the following technical solution: a multi-plane spread spectrum iterative phase recovery method based on angular spectrum diffraction theory, comprising the following steps:

[0006] Step 1: Obtain the light intensity distribution on the plane to be recovered and multiple propagation distance planes;

[0007] Step 2: Preprocess the phase distribution of the plane to be recovered. The preprocessing includes zero filling or random number filling. Then, combine the preprocessed phase distribution of the plane to be recovered with the light intensity distribution of the plane to be recovered to obtain the initial complex amplitude of the plane to be recovered.

[0008] Step 3: Iterate the initial complex amplitude against the light intensity distribution on the multiple propagation distance planes using a spread spectrum-based iterative phase recovery algorithm to obtain the preliminary recovered phase distribution of the multiple planes to be recovered;

[0009] Step 4: The preliminary recovered phases of the multiple planes to be recovered are fused using a multi-plane phase space frequency domain fusion algorithm based on angular spectrum to obtain the fused phase spectrum;

[0010] Step 5: Perform an inverse Fourier transform on the fused phase spectrum to obtain the phase recovery result after fusion.

[0011] In addition, the angular spectrum-based spread spectrum iterative phase retrieval algorithm in step 3 includes:

[0012] The angular spectrum iterative forward propagation model is constructed as follows:

[0013]

[0014] Among them, U 正 To update the complex amplitude of the propagation distance plane, A0 represents the angular spectrum corresponding to the initial complex amplitude U0. This represents the two-dimensional Fourier transform and its inverse transform, where z represents the propagation distance. Let represent the angular spectrum forward transfer function, where k, λ, and f represent the wave number, wavelength, and two-dimensional spatial frequency, respectively, j is the imaginary unit, and e is the base of the natural logarithm;

[0015] The angular spectrum iterative backpropagation model is constructed as follows:

[0016]

[0017] Among them, U 逆 To recover the updated complex amplitude of the plane, U 正 The corresponding angular spectrum, H -z To reconstruct the inverse transfer function of the angular spectrum,

[0018]

[0019] in, α(Iter.Now,f) is a dynamic adjustment factor, Iter.Now is the current iteration number, and I L It is an adjustable parameter.

[0020]

[0021] Where ln is the natural logarithm function, Iter.Total is the preset total number of iterations, and f max l1 represents the maximum frequency in two-dimensional space, and l1 is an adjustable exponent.

[0022] Additionally, the iterations in step 3 include:

[0023] The initial complex amplitude of the plane to be recovered and the light intensity distribution of multiple propagation distance planes are used as input parameters;

[0024] The initial complex amplitude of the plane to be recovered is obtained by using the angular spectrum iterative forward propagation model to obtain the updated complex amplitude of the plane with a propagation distance of z. Then, the amplitude is extracted from the updated complex amplitude of the plane with a propagation distance of z and replaced with the amplitude in the original complex amplitude of the plane with a propagation distance of z to obtain the complex amplitude of the plane with a propagation distance of z after amplitude update.

[0025] The complex amplitude of the plane with a propagation distance of z is obtained by using the angular spectrum iterative inverse propagation model to obtain the updated complex amplitude of the plane to be recovered. Then, the amplitude value in the updated complex amplitude of the plane to be recovered is extracted and replaced with the amplitude value in the initial complex amplitude of the plane to be recovered to obtain the complex amplitude of the plane to be recovered after amplitude update.

[0026] The iteration ends when the number of iterations reaches Iter.Total. Preliminary recovery phase distributions of multiple planes to be recovered are extracted from the complex amplitudes of the planes to be recovered after the amplitude update. i = 1, 2...n, where angle is the complex amplitude phase function.

[0027] Additionally, step 4 includes:

[0028] Using the preliminary phase distribution of multiple planes to be recovered obtained in step 3 According to the formula The phase spectra of each of the preliminarily recovered phases were obtained.

[0029] Constructing a multi-plane phase space frequency domain fusion algorithm based on angular spectrum:

[0030] For phase spectra |f|≤1 / λ, each phase spectrum is uniformly weighted. The phase spectrum after uniform weighting;

[0031] For phase spectra where |f|>1 / λ, each phase spectrum is assigned according to the weighting formula.

[0032]

[0033] at this time,

[0034] Among them, w i This represents the weight allocation, where m is a constant and represents a measure of distrust in the high spatial frequency domain inverse transfer function value. The phase spectrum after weight allocation;

[0035] The phase spectra after uniform weighting and weight allocation are spliced ​​and fused in the frequency domain to obtain the fused phase spectrum.

[0036] Compared with existing technologies, the invention has the following beneficial effects:

[0037] This invention introduces a rigorous angular spectrum diffraction model and combines the light intensity distribution information of multiple planes to recover the phase distribution of the plane to be recovered, thus achieving higher noise resistance. This invention adds a dynamic adjustment factor to the angular spectrum iterative transfer model, enabling iterative iteration to proceed normally even without window function isolation of high frequencies, reducing the loss of detailed information. To address the problem of multi-plane fusion spectrum methods failing to consider the differential amplification effect of high-frequency noise at different propagation distances and resulting in unreasonable fusion weight allocation, this invention assigns different weights to the phase spectrum at different spatial frequencies and performs phase spectrum fusion based on these weights, reducing the differential amplification of noise. Attached Figure Description

[0038] Figure 1 This is a block diagram of the multi-plane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory;

[0039] Figure 2 This is a schematic diagram of the iterative process in step 3 of the multiplane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory;

[0040] Figure 3 This is a comparison chart of the phase recovery experiment results. Detailed Implementation

[0041] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the various embodiments of the present invention will be described in detail below with reference to the accompanying drawings. However, those skilled in the art will understand that many technical details are presented in the various embodiments of the present invention to facilitate a better understanding of this application. However, the technical solutions claimed in this application can be implemented even without these technical details and various changes and modifications based on the following embodiments. The division of the various embodiments below is for ease of description and should not constitute any limitation on the specific implementation of the present invention. The various embodiments can be combined with and referenced by each other without contradiction. The technical solutions of the application will be further described in detail below with reference to the accompanying drawings.

[0042] The multi-plane spread-spectrum iterative phase retrieval method based on angular spectrum diffraction theory in this invention can be applied to image processing scenarios such as terahertz band compressed field measurement, laser wavefront detection, and microscopic imaging. In this embodiment, it can be used to recover the phase in a microscopic imaging scenario. Before implementing this invention, the light intensity distribution I0 of the plane to be recovered in the microscopic image, as well as the light intensity distributions on multiple other propagation planes, should already be obtained. The pre-recovery plane is the target plane for phase recovery, and the multiple propagation distance planes are the planes corresponding to the propagation distances required to recover the phase of the target plane. The propagation distances can be defined by the user. In this embodiment, a near-field scanning plane wave generator is used to obtain the light intensity distributions corresponding to the two propagation distance planes z1 and z2. And the pre-recovery planar light intensity distribution I0.

[0043] Based on this, the implementation steps of this embodiment are as follows: Figure 1 As shown, Figure 1 The multi-plane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory, as shown, includes the following steps:

[0044] Step 1: Obtain the light intensity distribution on the plane to be restored and on multiple propagation distance planes. In this embodiment, the obtained light intensity distribution on the pre-restored plane is I0, and the light intensity distribution on multiple propagation distance planes...

[0045] Step 2: Preprocess the phase distribution of the plane to be recovered. The preprocessing includes zero filling or random number filling. Then, combine the preprocessed phase distribution of the plane to be recovered with the light intensity distribution of the plane to be recovered to obtain the initial complex amplitude of the plane to be recovered.

[0046] Zero-filling refers to filling the missing parts or boundaries of the phase distribution of the plane to be recovered with zero values ​​to help reconstruct the complete phase distribution; random number filling refers to filling the missing parts or boundaries of the phase distribution of the plane to be recovered with random numbers, which can introduce more randomness and help avoid boundary effects. The specific preprocessing method can be selected according to actual needs.

[0047] In this embodiment, the phase distribution of the plane to be recovered after zero-filling or random number-filling preprocessing is as follows: By combining the light intensity distribution I0 of the plane to be recovered, the initial complex amplitude of the pre-recovered plane can be obtained. Where exp() represents an exponential function with base e, and j represents the imaginary unit.

[0048] Step 3: Iterate the initial complex amplitude against the light intensity distribution on the multiple propagation distance planes using a spread spectrum-based iterative phase recovery algorithm to obtain the preliminary recovered phase distribution of the multiple planes to be recovered.

[0049] In this embodiment, the light intensity distribution on the plane with propagation distances z1 and z2 is... As input parameters, the initial complex amplitude is also included. It is also used as an input parameter.

[0050] The angular spectrum-based spread spectrum iterative phase retrieval algorithm includes: constructing an angular spectrum iterative forward transfer model, as shown in the following formula:

[0051]

[0052] Among them, U 正 To update the complex amplitude of the propagation distance plane, A0 represents the initial complex amplitude U. 01 U 02 The corresponding angular spectrum, This represents the two-dimensional Fourier transform and its inverse transform, where z represents the propagation distance. Let represent the angular spectrum forward transfer function, where k, λ, and f represent the wave number, wavelength, and two-dimensional spatial frequency, respectively, j is the imaginary unit, and e is the natural constant.

[0053] The angular spectrum iterative backpropagation model is constructed as follows:

[0054]

[0055] Among them, U 逆 To recover the updated complex amplitude of the plane, U 正 The corresponding angular spectrum, H -z To reconstruct the inverse transfer function of the angular spectrum,

[0056]

[0057] in, α(Iter.Now,f) is a dynamic adjustment factor, Iter.Now is the current iteration number, and I L is an adjustable parameter, which can be regarded as a valve that can flexibly control the exponential growth of the inverse transfer function. l is an adjustable parameter used to suppress the growth of the inverse transfer function in the early and middle stages of iteration.

[0058]

[0059] Where ln is the natural logarithm function, Iter.Total is the preset total number of iterations, and f max l1 represents the maximum frequency in two-dimensional space, and l1 is an adjustable exponent.

[0060] The iterative process involves taking the initial complex amplitude of the plane to be recovered and the light intensity distribution of multiple propagation distance planes as input parameters;

[0061] The initial complex amplitude of the pre-recovered plane is obtained by using the angular spectrum iterative forward propagation model to obtain the updated complex amplitude of the plane with a propagation distance of z. Then, the updated complex amplitude of the plane with a propagation distance of z is replaced by amplitude substitution to obtain the complex amplitude of the plane with a propagation distance of z after amplitude update.

[0062] The complex amplitude of the plane with a propagation distance of z is obtained by using the angular spectrum iterative inverse propagation model to obtain the updated complex amplitude of the plane to be recovered. Then, the updated complex amplitude of the plane to be recovered is replaced by the amplitude to obtain the complex amplitude of the plane to be recovered after amplitude update.

[0063] The above update steps are repeated iteratively until the number of iterations reaches Iter.Total, at which point the iteration ends, ultimately yielding the preliminary recovered phase distributions of multiple planes to be recovered. i = 1, 2...n, where angle is a function for obtaining the phase of the complex amplitude.

[0064] Specifically, the angular spectrum iterative forward propagation model constructed based on the input in this embodiment is as follows: This represents two corresponding angular spectrum forward transfer functions;

[0065] Constructed angular spectrum iterative backpropagation model: in To reconstruct the inverse transfer function of the angular spectrum:

[0066]

[0067] After constructing the forward and reverse propagation models of the angular spectrum iteration, the iteration is carried out.

[0068] The initial complex amplitude of the plane to be recovered The updated complex amplitude of the plane with a propagation distance of z is obtained through the angular spectrum iterative forward propagation model. Then extract The amplitudes in the equations replace the original complex amplitude U of the plane with a propagation distance of z. z The amplitude in the middle is used to obtain the updated complex amplitude of the propagation distance plane with propagation distance z. Where angle is a function for determining the phase of the complex amplitude;

[0069] The original complex amplitude U of the plane with a propagation distance of z z The updated complex amplitude U of the plane to be recovered is obtained through the angular spectrum iterative inverse propagation model. 01new U 02new Then U 01new U 02new The amplitude in the original complex amplitude of the plane to be recovered is replaced by the amplitude in the original complex amplitude of the plane to be recovered, resulting in the updated complex amplitude of the plane to be recovered. A schematic diagram of the specific iterative process is shown below. Figure 2 As shown.

[0070] The iteration ends when the number of iterations reaches Iter.Total, ultimately obtaining the preliminary recovered phase of the plane to be recovered. At this point, the spread spectrum-based iterative phase recovery is complete.

[0071] Step 4: The preliminary recovered phases of the multiple planes to be recovered are fused using a multi-plane phase spatial frequency domain fusion algorithm based on angular spectrum to obtain the fused phase spectrum.

[0072] Specifically, the preliminary phase distribution of the multiple planes to be recovered is obtained in step 3. According to the formula The phase spectra of each of the preliminarily recovered phases were obtained.

[0073] Constructing a multi-plane phase space frequency domain fusion algorithm based on angular spectrum:

[0074] For phase spectra |f|≤1 / λ, each phase spectrum is uniformly weighted. The phase spectrum after uniform weighting;

[0075] For phase spectra where |f|>1 / λ, each phase spectrum is assigned according to the weighting formula.

[0076]

[0077] at this time,

[0078] Among them, w i This represents the weight allocation, where m is a constant and represents a measure of distrust in the high spatial frequency domain inverse transfer function value. The phase spectrum after weight allocation;

[0079] The phase spectra after uniform weighting and weight allocation are spliced ​​and fused in the frequency domain to obtain the fused phase spectrum.

[0080] The frequency domain splicing here includes summing all the uniformly weighted phase spectra with |f|≤1 / λ to obtain J1;

[0081] For each phase spectrum where |f|>1 / λ, multiply it by the corresponding weight to obtain the weighted phase spectrum. Then, add the weighted phase spectra together to obtain J2.

[0082] J1 and J2 are converted to frequency domain representation and then added to obtain complete frequency domain information. The fused phase spectrum is then obtained based on the complete frequency domain information.

[0083] Specifically, the preliminary recovered phase recovered in step 3 using multiple planes with different propagation distances z1 and z2 is... Based on the formula respectively Find their respective phase spectra

[0084] The fused phase spectrum is obtained by fusing the weighted or weighted phase spectra according to the frequency domain fusion algorithm. Different weights are assigned to the phase spectrum at different spatial frequencies, which has the effect of reducing the differential amplification of noise.

[0085] Step 5: Perform an inverse Fourier transform on the fused phase spectrum to obtain the fused phase recovery result. in, This indicates the final phase recovery result. Figure 3 This is a comparison diagram of the phase recovery results of the method of the present invention. Figure 3 (a) in the figure is the original measured image. Figure 3 (b) in the image shows the corresponding recovery result. Figure 3 The horizontal and vertical axes represent the sampling point numbers of the original measured image or the corresponding restored image, and the illustration shows the color depth.

[0086] The above-described embodiments are merely descriptions of the implementation methods of the present invention and are not intended to limit the scope of the present invention. Various modifications and improvements made by those skilled in the art to the technical solutions of the present invention without departing from the spirit of the present invention should fall within the protection scope defined by the claims of the present invention.

Claims

1. A multi-plane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory, characterized in that, Includes the following steps: Step 1: Obtain the light intensity distribution on the plane to be recovered and multiple propagation distance planes; Step 2: Preprocess the phase distribution of the plane to be recovered. The preprocessing includes zero filling or random number filling. Then, combine the preprocessed phase distribution of the plane to be recovered with the light intensity distribution of the plane to be recovered to obtain the initial complex amplitude of the plane to be recovered. Step 3: Iterate the initial complex amplitude against the light intensity distribution on the multiple propagation distance planes using a spread spectrum-based iterative phase recovery algorithm to obtain the preliminary recovered phase distribution of the multiple planes to be recovered; Step 4: The preliminary recovered phases of the multiple planes to be recovered are fused using a multi-plane phase space frequency domain fusion algorithm based on angular spectrum to obtain the fused phase spectrum; Step 5: Perform an inverse Fourier transform on the fused phase spectrum to obtain the phase recovery result after fusion; The angular spectrum-based spread spectrum iterative phase recovery algorithm in step 3 includes: The angular spectrum iterative forward propagation model is constructed as follows: ， in, To update the complex amplitude of the propagation distance plane, , Indicates the initial complex amplitude The corresponding angular spectrum, This represents the two-dimensional Fourier transform and its inverse transform, where z represents the propagation distance. Represents the forward transfer function of the angular spectrum. These represent wave number, wavelength, and two-dimensional spatial frequency, respectively. The imaginary unit, It is a natural constant; The angular spectrum iterative backpropagation model is constructed as follows: ， in, To recover the updated complex amplitude of the plane, express The corresponding angular spectrum, To reconstruct the inverse transfer function of the angular spectrum, ， in, , As a dynamic adjustment factor, This represents the current iteration number. , It is an adjustable parameter. ， in, It is the natural logarithm function. This is the preset total number of iterations. The maximum frequency in two-dimensional space. It is an adjustable index.

2. The multi-plane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory according to claim 1, characterized in that, The iterations in step 3 include: The initial complex amplitude of the plane to be recovered and the light intensity distribution of multiple propagation distance planes are used as input parameters; The initial complex amplitude of the plane to be recovered is obtained by using the angular spectrum iterative forward propagation model to obtain the updated complex amplitude of the plane with a propagation distance of z. Then, the amplitude is extracted from the updated complex amplitude of the plane with a propagation distance of z and replaced with the amplitude in the original complex amplitude of the plane with a propagation distance of z to obtain the complex amplitude of the plane with a propagation distance of z after amplitude update. The complex amplitude of the plane with a propagation distance of z is obtained by using the angular spectrum iterative inverse propagation model to obtain the updated complex amplitude of the plane to be recovered. Then, the amplitude value in the updated complex amplitude of the plane to be recovered is extracted and replaced with the amplitude value in the initial complex amplitude of the plane to be recovered to obtain the complex amplitude of the plane to be recovered after amplitude update. When the number of iterations reaches The iteration ends when the amplitude is updated, and the preliminary recovery phase distribution of multiple planes to be recovered is extracted from the complex amplitude of the plane to be recovered. i = 1, 2...n To obtain the complex amplitude phase function.

3. The multi-plane spread spectrum iterative phase retrieval method based on angular spectrum diffraction theory according to claim 2, characterized in that, Step 4 includes: Using the preliminary phase distribution of multiple planes to be recovered obtained in step 3 According to the formula The phase spectra of each of the preliminarily recovered phases were obtained. ; Constructing a multi-plane phase space frequency domain fusion algorithm based on angular spectrum: for The phase spectrum is uniformly weighted. , The phase spectrum after uniform weighting; for The phase spectrum is then distributed according to a weighting formula: ， at this time, , in, This represents the weight allocation, where m is a constant and represents a measure of distrust in the high spatial frequency domain inverse transfer function value. The phase spectrum after weight allocation; The phase spectra after uniform weighting and weight allocation are spliced ​​and fused in the frequency domain to obtain the fused phase spectrum.

Images