CN110309482B - Fast convergence and high-precision phase recovery method - Google Patents
Fast convergence and high-precision phase recovery method Download PDFInfo
- Publication number
- CN110309482B CN110309482B CN201910407018.6A CN201910407018A CN110309482B CN 110309482 B CN110309482 B CN 110309482B CN 201910407018 A CN201910407018 A CN 201910407018A CN 110309482 B CN110309482 B CN 110309482B
- Authority
- CN
- China
- Prior art keywords
- complex amplitude
- amplitude distribution
- object plane
- light field
- initial
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 18
- 238000011084 recovery Methods 0.000 title claims abstract description 15
- 238000004364 calculation method Methods 0.000 claims abstract description 32
- 230000003287 optical effect Effects 0.000 claims abstract description 18
- 238000001228 spectrum Methods 0.000 claims abstract description 15
- 230000003044 adaptive effect Effects 0.000 claims abstract description 10
- 238000005457 optimization Methods 0.000 claims abstract description 8
- 238000005070 sampling Methods 0.000 claims description 6
- 238000001816 cooling Methods 0.000 claims description 4
- 230000002068 genetic effect Effects 0.000 claims description 4
- 230000035772 mutation Effects 0.000 claims description 3
- 238000000137 annealing Methods 0.000 claims description 2
- 230000000977 initiatory effect Effects 0.000 claims description 2
- 230000009191 jumping Effects 0.000 claims 1
- 238000005516 engineering process Methods 0.000 description 3
- 238000002922 simulated annealing Methods 0.000 description 3
- 230000005540 biological transmission Effects 0.000 description 2
- 239000011159 matrix material Substances 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000012984 biological imaging Methods 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- ONUFESLQCSAYKA-UHFFFAOYSA-N iprodione Chemical compound O=C1N(C(=O)NC(C)C)CC(=O)N1C1=CC(Cl)=CC(Cl)=C1 ONUFESLQCSAYKA-UHFFFAOYSA-N 0.000 description 1
- 238000010845 search algorithm Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
- G01J1/00—Photometry, e.g. photographic exposure meter
- G01J1/42—Photometry, e.g. photographic exposure meter using electric radiation detectors
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/12—Computing arrangements based on biological models using genetic models
- G06N3/126—Evolutionary algorithms, e.g. genetic algorithms or genetic programming
Abstract
The invention provides a fast convergence and high-precision phase recovery method, which meets the requirements of fast convergence and high-precision solving of the complex amplitude of a light field, and not only avoids the problems that an algorithm is easily influenced by an initial assumed value and is not easy to converge, but also avoids the problem of low solving precision caused by parameterization of the complex amplitude of the light field. The method comprises the following steps: firstly, assuming a group of random distribution coefficients, obtaining the complex amplitude distribution of the object plane initial assumed measured optical field according to a light intensity calculation model and a phase calculation model, obtaining an optimized coefficient value through fast convergence based on a global convergence optimization algorithm of a simulated annealing-genetic algorithm, and obtaining a global optimal solution of the complex amplitude distribution of the measured optical field with low initial precision according to the optimized coefficient value; and finally solving to obtain the high-precision light field complex amplitude distribution of the object plane position by iterative computation by taking the global optimal solution of the measured light field complex amplitude distribution with low initial precision as an initial solution of the adaptive angular spectrum computation.
Description
Technical Field
The invention belongs to the field of optics, relates to a phase recovery method for solving complex amplitude of a light field with fast convergence and high precision, and particularly relates to a phase recovery method based on a parameterized mode method and adaptive angular spectrum calculation.
Background
The phase recovery technology utilizes a diffraction model of a light field to perform diffraction calculation on an assumed input light field to obtain the intensity distribution of the light field of an output surface, compares the calculated intensity of the light field of the output surface with the intensity of the light field generated by a real phase, and calculates the phase distribution which best accords with the intensity distribution of the real light field by using different phase recovery algorithms with the minimum error of the two as a criterion.
In 1972, Gerchberg and Saxton proposed a GS algorithm, which initiated the application of a phase recovery technology, and subsequently j.r. fienup et al proposed an Error Reduction (Error Reduction) algorithm and a Hybrid Input-Output (Hybrid Input-Output) algorithm, which can be classified as a first-class algorithm, a GS algorithm and an improved algorithm thereof, and the core is a fourier iterative transform algorithm.
The second type of algorithm is a gradient search algorithm, which adopts a gaussian function to simulate the intensity of the light field, adopts a zernike polynomial to fit the phase, sets an objective function (generally, a root-mean-square error value of the calculated light intensity of the output surface and the detected light intensity of the output surface), and solves the minimum value of the objective function along the gradient descending direction, so that the solution of the complex amplitude of the light field is converted into the solution of the parameterized coefficient.
The third type of algorithm is a modern optimization algorithm represented by a Simulated Annealing (SA) algorithm and a Genetic (GA) algorithm, which rely on probability search and only require the output value of a target function without gradient information thereof.
At present, the three algorithms for solving the complex amplitude of the light field have the following advantages and disadvantages: 1) the first type of algorithm has high solving precision, but is easily influenced by an initial assumed value of the algorithm and is extremely easy to not converge; 2) the second type of algorithm has strong local searching capability but weak global searching capability and is easy to fall into a state of local convergence or non-convergence; 3) the third type of algorithm does not depend on gradient information during optimization calculation, does not require continuous and derivable objective functions, and has strong competitiveness in the optimization problem; 4) compared with the first algorithm, the latter two algorithms are not easily affected by the initial assumed value of the algorithm, but parameterize the light field intensity and the phase, so that the solving precision is not high. When the phase recovery technology is used for solving the complex amplitude of the light field, the algorithm is required to be fast in convergence, and the algorithm is required to be high in solving precision.
Disclosure of Invention
The invention provides a fast convergence and high-precision phase recovery method, which meets the requirements of fast convergence and high-precision solving of the complex amplitude of a light field, and not only avoids the problems that an algorithm is easily influenced by an initial assumed value and is not easy to converge, but also avoids the problem of low solving precision caused by parameterization of the complex amplitude of the light field.
The solution of the invention is as follows:
the fast convergence and high-precision phase recovery method comprises the following steps:
firstly, assuming a group of random distribution coefficients, obtaining the complex amplitude distribution of the object plane initial assumed measured optical field according to a light intensity calculation model and a phase calculation model, obtaining an optimized coefficient value through fast convergence based on a global convergence optimization algorithm of a simulated annealing-genetic algorithm, and obtaining a global optimal solution of the complex amplitude distribution of the measured optical field with low initial precision according to the optimized coefficient value;
and finally solving to obtain the high-precision light field complex amplitude distribution of the object plane position by iterative computation by taking the global optimal solution of the measured light field complex amplitude distribution with low initial precision as an initial solution of the adaptive angular spectrum computation.
The invention has the following advantages:
1. neither gradient information nor the objective function is required to be continuously derivable.
2. The genetic algorithm has strong global search capability but weak local search capability, and generally only can obtain suboptimal solution of the problem, while the simulated annealing algorithm has strong local search capability, and the global optimal solution of the problem can be obtained by combining the two algorithms.
3. By utilizing the global convergence optimizing algorithm based on the simulated annealing-genetic algorithm, a solution of the light field complex amplitude with low precision can be obtained through fast convergence, and the problems that the light field complex amplitude is sensitive to an initial assumed value and is easy to not converge when the first type of algorithm is used for recovering the light field complex amplitude are solved.
4. Although the light field complex amplitude is parameterized, the light field complex amplitude can be solved with high precision.
5. The transmission process is calculated by utilizing the self-adaptive angular spectrum theory, the long-distance transmission can be met by adding a calculation window on an object plane and limiting the bandwidth of a transfer function, and the requirement on the experimental distance is not harsh.
6. The phase recovery method with fast convergence and high precision can be widely applied to the aspects of high-precision wavefront measurement, optical device surface type measurement, x-ray imaging, biological imaging and the like.
Drawings
FIG. 1 is a schematic flow chart of the present invention.
Detailed Description
The invention relates to a phase recovery method based on a parameterized mode method and adaptive angular spectrum calculation, which is based on a global convergence optimization algorithm and an adaptive angular spectrum calculation method of a simulated annealing-genetic algorithm and mainly comprises the following steps: firstly, setting the light field complex amplitude distribution as random distribution at the object plane position, rapidly converging and solving a global optimal solution with low initial precision through a designed global convergence optimization algorithm based on a simulated annealing-genetic algorithm, then taking the global optimal solution as an initial solution of adaptive angular spectrum calculation, and finally solving to obtain the high-precision light field complex amplitude distribution at the object plane position through iterative calculation.
The algorithm flow of the invention is shown in fig. 1, and the specific steps are as follows:
step 1: the initial assumption is that the object plane light field complex amplitude distribution is
Wherein the light intensity Im(x0,y0)＝aexp[-b(x0 2+y0 2)]Phase of
Fitting with zernike polynomials. z is a radical of0，z1，z2···znIs a set of Zernike polynomials, a, b, a, which are orthogonal on a unit circle0，a1，a2···anFor the coefficients to be fitted, the initial values are a set of random numbers in the range of (0, 1)]。
Step 2: initializing control parameters: size of population, maximum evolution algebra genmaxCross probability PcProbability of mutation PmAnnealing initiation temperature T0Temperature cooling coefficient T and termination temperature Tend。
Step 3: initializing the population and calculating the fitness function f of each individualiWherein i is 1,2, ·, size. The process uses the light intensity acquisition information of the measured optical field.
Step 4: the loop count variable gen is set to 0.
Step 5: and (4) carrying out genetic operations such as selection, crossing and mutation on the population individuals to generate new population individuals. For individuals in the new population, a new individual fitness function f is calculatedi', if fi'＞fiThen toThe new individual replaces the old individual; otherwise, the probability P is equal to exp [ (f)i-f'i)Ti]The new individuals are accepted and the old individuals are discarded.
Step 6: if gen < genmaxIf so, go to Step 5; otherwise go to Step 7.
Step 7: if Ti＜TendIf the algorithm is successfully ended, returning to the global optimal coefficient; otherwise, executing a cooling operation Ti+1＝tTiGo to Step 4.
Step 8: according to the light intensity model and the phase model in Step1, obtaining the global optimal coefficients a, b, a0，a1，a2…anTo obtain the global optimal solution U of the complex amplitude distribution of the measured optical field0(x0,y0)。
Step 9: the size of the object plane calculation window S is set, the calculation formula is S2N Δ x, and the extended sampling point M is 2N. Where Δ x is the sampling interval of the object plane and N is the object plane sampling point.
Step 10: initial measured light field complex amplitude distribution U obtained by inputting Step80(x0,y0) Calculating to obtain U (x) after dimension expansion according to the size of the calculation window0,y0) The matrix dimension is M multiplied by M, and the matrix dimension is used as an initial solution of the self-adaptive angular spectrum calculation method.
Step 11: transmitting along optical path, calculating light field complex amplitude distribution U of diffraction surfacec(xm,ym). The specific calculation process is as follows: (1) firstly, calculating the object plane frequency domain light field distribution A (U, v,0) ═ F [ U (x)0,y0)]F is forward fast Fourier transform; (2) calculating a transfer function after adding an adaptive constraint
Where H (u, v, z) represents the transfer function from the object plane to the diffraction plane, ulim＝[(2Δuz)2+1]-1/2λ-1，vlim＝[(2Δvz)2+1]-1/2λ-1，Δu＝Δv＝(2S)-1Z is the distance from the object plane to the diffraction plane; (3) Calculating the light field complex amplitude distribution of the diffraction surface: u shapec(xm,ym)＝F-1[A(u,v,0)H'(u,v,z)]，F-1Is an inverse fast fourier transform.
Step 12: determining an objective function
I(xm,ym) Setting a threshold value for the recorded actual light spot intensity; if the amplitude is smaller than the preset value, reversely transmitting the signal to the object plane along the light path according to the self-adaptive angular spectrum theory, and calculating the final expanded object plane light field complex amplitude distribution Uk(x0,y0) Jump Step 15; otherwise, jump to Step 13.
Step 13: using the actual spot intensity I (x)m,ym) Correcting the complex amplitude distribution of the diffraction surface light field obtained by calculation, keeping the phase unchanged, and obtaining the updated complex amplitude of the diffraction surface light field
Step 14: reversely transmitting to object plane along optical path according to self-adaptive angular spectrum theory, calculating to obtain new object plane light field complex amplitude distribution Uk(x0,y0) Let U (x)0,y0)＝Uk(x0,y0) And repeating Step 11-12.
Step 15: to Uk(x0,y0) Reducing dimension, converting into NxN dimension, and finally obtaining light field complex amplitude distribution U of object plane0(x0,y0)。
Claims (2)
1. A fast convergence and high precision phase recovery method is characterized by comprising the following steps:
firstly, assuming a group of random distribution coefficients, obtaining the complex amplitude distribution of the object plane initial assumed measured optical field according to a light intensity calculation model and a phase calculation model, obtaining an optimized coefficient value through fast convergence based on a global convergence optimization algorithm of a simulated annealing-genetic algorithm, and obtaining a global optimal solution of the complex amplitude distribution of the measured optical field with low initial precision according to the optimized coefficient value;
taking the global optimal solution of the complex amplitude distribution of the measured optical field with low initial precision as the initial solution of the adaptive angular spectrum calculation, and finally solving to obtain the complex amplitude distribution of the optical field with high precision at the object plane position through iterative calculation;
the method comprises the following specific steps:
step 1: substituting a group of randomly distributed initial coefficients to be optimized into a light intensity calculation model and a phase calculation model, and fitting the complex amplitude distribution of the initial hypothetical object plane light field;
step 2: a global convergence optimizing algorithm based on a simulated annealing-genetic algorithm initializes control parameters: setting the size of population individual and the maximum evolution algebra genmaxCross probability PcProbability of mutation PdAnnealing initiation temperature T0Temperature cooling coefficient T and termination temperature Tend；
And step 3: initializing the population and calculating the fitness function f of each individualiWherein i ═ 1,2, ·, size;
and 4, step 4: setting a loop counting variable gen to be 0;
and 5: carrying out genetic operation on the population individuals to generate new population individuals; for individuals in the new population, a new individual fitness function f is calculatedi', if fi'＞fiReplacing the old individual with the new individual; otherwise, the probability P is equal to exp [ (f)i-fi')Ti]Accepting new individuals, discarding old individuals, wherein TiIs the current temperature;
step 6: if the current evolution algebra gen is less than genmaxIf so, go to step 5; otherwise, turning to step 7;
and 7: if Ti＜TendIf the algorithm is successfully ended, returning the optimized coefficient value, and then executing the subsequent steps; otherwise, executing a cooling operation Ti+1＝tTiTurning to the step 4;
and 8: substituting the optimized coefficient value into the light intensity calculation model and the phase calculation model in the step1 to obtain the initial precisionGlobal optimal solution U of high measured optical field complex amplitude distribution0(x0,y0)；
And step 9: setting the size S of an object plane calculation window to be 2N delta x according to the number N of sampling points of the object plane and the sampling interval delta x, and solving the number M of the sampling points after dimension expansion to be 2N;
step 10: inputting the complex amplitude distribution U of the initial measured optical field obtained in step80(x0,y0) Calculating the complex amplitude distribution U (x) of the measured optical field after dimension expansion according to the size S of the calculation window0,y0) Taking the initial solution as an initial solution of an adaptive angular spectrum calculation method;
step 11: transmitting along optical path, calculating light field complex amplitude distribution U of diffraction surface according to self-adaptive angular spectrum calculation methodc(xm,ym)；
Step 12: determining an objective function
Whether or not less than E, I (x)m,ym) Setting a threshold value for the recorded actual light spot intensity; if the amplitude is smaller than the preset value, reversely transmitting the signal to the object plane along the light path according to the self-adaptive angular spectrum theory, and calculating the final expanded object plane light field complex amplitude distribution Uk(x0,y0) Skipping to step 15; otherwise, jumping to step 13;
step 13: using the actual spot intensity I (x)m,ym) Correcting the complex amplitude distribution of the diffraction surface light field obtained by calculation, keeping the phase unchanged, and obtaining the updated complex amplitude distribution of the diffraction surface light field
Step 14: reversely transmitting the signal to an object plane along an optical path according to the self-adaptive angular spectrum theory, and calculating to obtain a new object plane light field complex amplitude distribution U'k(x0,y0) Let U (x)0,y0)＝U'k(x0,y0) Repeating the steps 11-12;
step 15: to Uk(x0,y0) Performing dimensionality reduction operation, converting into NxN dimensionality, and finally obtaining the light field complex amplitude distribution U 'of the object plane'0(x0,y0)。
2. The fast convergence and high precision phase recovery method according to claim 1, wherein the specific calculation process of step11 is as follows:
(1) firstly, calculating the object plane frequency domain light field distribution A (U, v,0) ═ F [ U (x)0,y0)]F is forward fast Fourier transform;
(2) calculating a transfer function after adding an adaptive constraint
Where H (u, v, z) represents the transfer function from the object plane to the diffraction plane, ulim＝[(2Δuz)2+1]-1/2λ-1，vlim＝[(2Δvz)2+1]-1/2λ-1，Δu＝Δv＝(2S)-1Z is the distance from the object plane to the diffraction plane;
(3) calculating the light field complex amplitude distribution of the diffraction surface: u shapec(xm,ym)＝F-1[A(u,v,0)H'(u,v,z)]，F-1Is an inverse fast fourier transform.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910407018.6A CN110309482B (en) | 2019-05-16 | 2019-05-16 | Fast convergence and high-precision phase recovery method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910407018.6A CN110309482B (en) | 2019-05-16 | 2019-05-16 | Fast convergence and high-precision phase recovery method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110309482A CN110309482A (en) | 2019-10-08 |
| CN110309482B true CN110309482B (en) | 2020-11-17 |
Family
ID=68074803
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910407018.6A Active CN110309482B (en) | 2019-05-16 | 2019-05-16 | Fast convergence and high-precision phase recovery method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110309482B (en) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111580271B (en) * | 2020-06-24 | 2022-02-18 | 哈工大机器人(中山)无人装备与人工智能研究院 | Self-adaptive aberration correction method and light sheet microscopic imaging device based on same |
| CN113376448B (en) * | 2021-04-29 | 2023-02-28 | 北京邮电大学 | Method and device for quiet zone phase recovery in compact range test |
Family Cites Families (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102230827B (en) * | 2011-06-18 | 2012-11-28 | 四川大学 | Mach-Zehnder point diffraction interferometer and method for reconstructing laser complex amplitudes |
| CN107301458B (en) * | 2017-05-31 | 2020-08-11 | 重庆理工大学 | Phase optimization method based on intensity transmission equation |
| CN109060122B (en) * | 2018-07-05 | 2021-02-12 | 安徽大学 | Two-step phase recovery method, equipment and system based on single intensity measurement |
| CN108983579B (en) * | 2018-09-05 | 2020-05-01 | 南京大学 | Method and device for phase recovery and reconstruction of lensless digital holographic microscopic imaging |
-
2019
- 2019-05-16 CN CN201910407018.6A patent/CN110309482B/en active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN110309482A (en) | 2019-10-08 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109031654B (en) | Adaptive optical correction method and system based on convolutional neural network | |
| CN110309482B (en) | Fast convergence and high-precision phase recovery method | |
| CN110308649B (en) | PID parameter optimization method based on PSO-SOA fusion algorithm and applied to industrial process control | |
| CN104331635B (en) | The method of power optical fiber Communication ray power prediction | |
| CN108520325B (en) | Integrated life prediction method based on accelerated degradation data in variable environment | |
| CN108983412A (en) | A kind of no Wave-front measurement adaptive optics system and beam phase method of adjustment | |
| CN110543929A (en) | wind speed interval prediction method and system based on Lorenz system | |
| CN105954743B (en) | A kind of variable structure multi-model maneuvering target tracking method improving weights | |
| CN112597700B (en) | Aircraft trajectory simulation method based on neural network | |
| US20110008032A1 (en) | Autofocus method | |
| CN113158487B (en) | Wavefront phase difference detection method based on long-short term memory depth network | |
| CN106503861B (en) | Wind power prediction method based on probability statistics and particle swarm optimization and integrating wind speeds of multiple meteorological sources | |
| CN110895772A (en) | Electricity sales amount prediction method based on combination of grey correlation analysis and SA-PSO-Elman algorithm | |
| CN115204058B (en) | Ground wave radar flow field calculation method and device based on BP neural network | |
| CN110675278A (en) | Photovoltaic power short-term prediction method based on RBF neural network | |
| CN109932817A (en) | Means for correcting and its bearing calibration before a kind of adaptive optics simulation dynamic wave | |
| CN107894662A (en) | A kind of improvement SPGD algorithms for OAM beam phase reparations | |
| CN114035300B (en) | Self-adaptive correction method based on zoom lens under high-capacity space optical communication link | |
| Petković et al. | RETRACTED: Adaptive neuro-fuzzy prediction of modulation transfer function of optical lens system | |
| KR101762357B1 (en) | Control device and method for voltage source inverter | |
| Zirkohi et al. | Design of Radial Basis Function Network Using Adaptive Particle Swarm Optimization and Orthogonal Least Squares. | |
| CN115775231A (en) | Cascade R-CNN-based hardware defect detection method and system | |
| CN105357435B (en) | A kind of daytime environment Stellar orientation instrument focus adjustment method | |
| CN109901398B (en) | Peak upper limit estimation method of nonlinear system impulse response | |
| CN115035336A (en) | Three-level T-type inverter fault diagnosis method based on GASF and CNN |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |