CN111103686B

CN111103686B - Method for realizing self-adaptive optical SPGD control algorithm based on FPGA

Info

Publication number: CN111103686B
Application number: CN201911424233.3A
Authority: CN
Inventors: 杨慧珍; 吴阳; 徐奇; 龚成龙
Original assignee: Marine Resources Development Institute Of Jiangsu (lianyungang); Jiangsu Ocean University
Current assignee: Marine Resources Development Institute Of Jiangsu (lianyungang); Jiangsu Ocean University
Priority date: 2019-12-31
Filing date: 2019-12-31
Publication date: 2022-05-03
Anticipated expiration: 2039-12-31
Also published as: CN111103686A

Abstract

The invention relates to a method for realizing a wavefront-free detection self-adaptive optical SPGD control algorithm based on an FPGA, which comprises the following steps: the random disturbance module generates a plurality of groups of pseudo-random sequences by means of Matlab software (not limited to Matlab), and each group of sequences is multiplied by disturbance voltage and added with corresponding signs in the FPGA to generate positive disturbance and negative disturbance; the performance index function calculation module selects the average Radius MR-Mean Radius of the far-field light spot as a performance index function; and the computation deformable mirror voltage control module collects image data and computes a performance index function after disturbance application is finished, performs gradient estimation on the computed performance index function, and then computes a control voltage according to a formula. The method disclosed by the invention is widely applicable to the self-adaptive optical system with any unit number, improves the universality of the mode of realizing the SPGD control algorithm by the FPGA, and ensures the correction capability of the wavefront-free detection self-adaptive optical system.

Description

Method for realizing self-adaptive optical SPGD control algorithm based on FPGA

Technical Field

The invention belongs to the technical field of adaptive optics, relates to a wavefront-free detection adaptive optics system, and particularly relates to a method for realizing an adaptive optics SPGD control algorithm based on an FPGA.

Background

Vorotstov, m.a. of the american army institute, proposes a Stochastic Parallel Gradient Descent (SPGD) algorithm and applies the algorithm to the field of adaptive optics. The SPGD-based wavefront-free detection adaptive optical system does not need wavefront measurement and phase reconstruction, and complexity of the system is reduced. The application field of the adaptive optics technology is further expanded, and particularly the application field of the conventional adaptive optics system cannot be applied.

The correction speed is one of important indexes of the correction capability of the adaptive optics system, and the SPGD control algorithm searches for an optimal solution in an iterative mode. Therefore, the iteration speed of the algorithm plays a key role in the correction capability of the adaptive optics system, and the implementation mode of the algorithm has a decisive influence on the iteration speed. In addition, the diversity of the application environment requires that the implementation of the algorithm has certain versatility. Therefore, the implementation of the SPGD algorithm needs to satisfy both the requirements of real-time performance and flexibility.

The SPGD control algorithm in the adaptive optics system is implemented in two ways: firstly, the method is realized on a general PC by using a software programming method, and has the advantages of simple realization, strong flexibility and the like, but Windows operated by the PC is a micro-kernel, preemptive multitask and soft real-time operating system, and tasks are scheduled according to the priority order of threads, so that the real-time performance of the whole system is weak. And the second is realized in a hardware mode. The method has the advantages of high iteration speed and strong real-time performance, meets the requirements of a self-adaptive optical system, but has complex circuit design, high chip cost and lack of flexibility.

Zhang jin Bao, Zhang super et al use FPGA to realize SPGD control algorithm. The FPGA is a programmable integrated circuit, the functions of the FPGA are expressed by a hardware description language, the FPGA can be configured into a hardware circuit capable of finishing specific functions (such as SPGD control algorithm), the FPGA has the characteristics of repeated programming, flexible algorithm realization and the like, and the essence of a special hardware circuit can also ensure the requirement of real-time performance. However, the existing mode of realizing the SPGD control algorithm by the FPGA is only for a specific deformable mirror, and the application range is limited. In addition, in the prior art, matrix recovery is not performed on image data serially transmitted into the FPGA, which causes that part of performance index functions cannot be realized by the FPGA or the realization mode is not accurate, and reduces the universality of the SPGD control algorithm realized by the FPGA. Meanwhile, the performance index function in the prior art is not accurate, which affects the correction capability of the adaptive optical system.

The invention content is as follows:

the technical problem to be solved by the invention is to overcome the defects of the prior art and provide the implementation method of the wavefront-free detection adaptive optical SPGD control algorithm based on the FPGA, the application range of the method is wide and the method is applicable to adaptive optical systems with any unit number, the universality of the SPGD control algorithm implementation mode of the FPGA is improved, and the correction capability of the wavefront-free detection adaptive optical system is ensured.

The technical problem to be solved by the present invention is achieved by the following technical means. The invention relates to a method for realizing a wavefront-free detection self-adaptive optical SPGD control algorithm based on an FPGA, which is characterized by comprising the following steps of:

(1) the random disturbance module generates a plurality of groups of pseudo-random sequences by means of software (including but not limited to Matlab software), multiplies each group of sequences by disturbance voltage in the FPGA and adds corresponding signs to generate positive disturbance and negative disturbance;

(2) selecting the average Radius MR-Mean Radius of the far-field light spot as a performance index function, wherein the performance index function is not limited to MR;

wherein (x)₀,y₀) Is the centroid of the spot;

i (x, y) represents the pixel gray value at the corresponding coordinate (x, y); the ROI area of the spot is NXN, so x, y ∈ [1, N ∈]；

(3) A voltage control module for calculating the deformable mirror according to a formula u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-γΔu⁽ⁿ⁾ΔJ⁽ⁿ⁾And calculating control voltage, optimizing the performance index function MR in the direction of the minimum value, and taking a positive value for gamma.

The invention relates to a method for realizing a wavefront-free detection self-adaptive optical SPGD control algorithm based on an FPGA, which further adopts the preferable technical scheme that:

the characteristics of random disturbance voltage in the SPGD algorithm provide corresponding requirements for generating a plurality of groups of pseudo-random sequences by Matlab software (not limited to Matlab):

a. the number of paths is equal to the number of deformable mirror units; n paths of pseudo-random sequences are needed for an N-unit deformable mirror;

b. the pseudo-random sequences are mutually independent in pairs; the driving units are independent to each other, so that mutual coupling among the driving units of the deformable mirror can be avoided;

c. the pseudorandom sequence follows a bernoulli distribution with a probability of two sample values each of 0.5.

The SPGD control algorithm adds random parallel disturbance to the multipath control parameters, and searches the control parameters by using the variation of the performance index measured value and the variation of the control parameters.

The invention relates to a wavefront-detection-free self-adaptive optical SPGD control method based on an FPGA, which adopts the following steps to realize an SPGD control algorithm:

(1) generating a random number sequence by software and storing the random number sequence into a development board of the FPGA;

(2) FPGA reads the random sequence and generates a disturbance voltage (delta u)ⁿ},(j＝1,...,N,|Δu⁽ⁿ⁾|＝0.2,sgn(Δu⁽ⁿ⁾)＝±1)；

(3) FPGA (field programmable Gate array) applies positive disturbance u to deformable mirror⁽ⁿ⁾+Δu⁽ⁿ⁾，u⁽ⁿ⁾The control voltage is calculated in the (n-1) th iteration, image data is collected, and a performance index function is calculated

(4) FPGA (field programmable Gate array) applies negative disturbance u to deformable mirror⁽ⁿ⁾-Δu⁽ⁿ⁾Acquiring image data, calculating performance index function

(5) Calculating the variation of the performance index function in the steps (3) and (4)

Calculating the gradient gamma delta u according to the variation⁽ⁿ⁾ΔJ⁽ⁿ⁾According to the formula u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-γΔu⁽ⁿ⁾ΔJ⁽ⁿ⁾Completing the calculation of the control voltage of the deformable mirror;

and (5) repeating the steps (3) to (5) for a plurality of iterations until the optimal control voltage of the deformable mirror is found, so that the sum of the initial aberration (psi (r)) and the correction phase (m (r)) of the deformable mirror approaches to 0.

Compared with the prior art, the invention has the following beneficial effects:

1. the invention realizes the SPGD control algorithm by utilizing the FPGA, is convenient to realize and has certain flexibility, and can ensure the requirement of the algorithm on real-time property.

2. According to the invention, the image data is stored in the address bit of the ROM core in the FPGA, the judgment condition is set, the matrix recovery is carried out on the image data, the commonly used performance index function can be accurately realized by the FPGA, the universality of the SPGD control algorithm realizing mode of the FPGA is improved, and the correction capability of the wavefront-free detection adaptive optical system based on the SPGD control algorithm of the FPGA is not influenced by the calculation of the performance index function.

3. After the matrix recovery is finished, the image data is accurately positioned by using the two counters and some auxiliary signals, the relation between the image data and the matrix coordinates is specified in the FPGA, the compilation of the centroid calculation process is finished according to a centroid formula, and the centroid of the image matrix is accurately calculated. The method for realizing the MR calculation module by using the FPGA has higher accuracy, and the wavefront-free detection adaptive optical system based on the SPGD control algorithm of the FPGA has stronger correction capability.

4. Through theoretical analysis, the gain coefficient and the disturbance amplitude are both normal numbers. The invention multiplies the two in the preprocessing part to obtain a result of a normal number, and directly defines the constant in the FPGA, thereby reducing the use of multipliers, saving the internal resources of the FPGA, reducing the complexity of the algorithm and improving the calculation speed of the algorithm to a certain extent.

Drawings

FIG. 1 is a timing diagram for generating random perturbations;

FIG. 2 is a flow chart of the calculation of the performance indicator function MR;

FIG. 3 is a timing diagram of fractional computation of the performance indicator function MR molecules;

FIG. 4 is a flow chart of calculating a deformable mirror control voltage;

FIG. 5 is a timing diagram of a portion of calculating the deformable mirror control voltage;

fig. 6 is a structural block diagram of an experimental system for a wavefront-less detection adaptive optics SPGD control algorithm based on an FPGA.

Detailed Description

The embodiments of the present invention will be further described below with reference to the accompanying drawings so that those skilled in the art can further understand the present invention.

The following embodiment takes a 97-unit deformable mirror as an example:

a module for generating random disturbance voltage:

fig. 1 is a timing diagram of a random perturbation module according to the present invention, where Clk is a system clock, and the implementation is performed as follows.

1. And the preprocessing part generates a plurality of sets of 97-way pseudorandom sequences (formed by combining 0 and 1) by utilizing a rand function by means of Matlab software but not limited to Matlab. And storing the generated sets of 97-way pseudo-random sequences for later use.

2. Taking one group of pseudo-random sequences as an example to calculate positive and negative disturbances, reading and storing the group of pseudo-random sequences into a ROM core of the FPGA. A counter cnt is set and read the pseudo-random sequence stored in the ROM as address bits (when cnt is equal to 0, it corresponds to the 1 st value of the 97 way pseudo-random sequence, and so on).

3. The corresponding relationship between the value of the counter cnt and the 97-way pseudo-random sequence is shown in fig. 1, where Q represents a random sequence read from the ROM, and when cnt is equal to 0, the corresponding Q is output when cnt is equal to 2 (there is a certain delay in reading data from the ROM). The Q between the two dotted lines in the figure is the 97-way pseudo-random sequence read from the ROM.

Q equals 0, the corresponding perturbation is 1, Q equals 1, the corresponding perturbation is-1, and the specific result is shown as signal d in FIG. 1, which is a total of 97 random combinations of 1 and-1. Since the value of Q is determined first, the corresponding value of d is delayed, and d between two dot-dash lines (labeled as start of generation of disturbance and completion of 97 disturbances) is determined to be the required 97 random disturbances with a disturbance amplitude of 1 by using a certain auxiliary signal.

5. On the basis of step 4, d is multiplied by 200 to obtain dv1 (positive perturbation) in fig. 1. 200 is the perturbation amplitude, which is defined as the binary fixed point decimal with the last 10 digits being decimal places, so 200 corresponds to an actual perturbation amplitude of 0.1953 (approximately 0.2).

6. The addition of a negative sign to each perturbation on the basis of dv1 yields the data dv2 in fig. 1 (negative perturbation).

And calculating two performance indexes:

fig. 2 is a flow chart of the calculation of the performance indicator function MR. In specific implementation, the method comprises the following steps:

1. image data (gray values of a spot image, and an ROI (region of interest) is N multiplied by N), and since the image data is transmitted to the interior of the FPGA in series, matrix recovery is needed before calculation. The image data is restored to an N × N matrix form with the address bits as the determination conditions.

2. Calculating the centroid, the formula of the centroid is

The calculation of the centroid according to the formula is calculated as follows:

(1) calculating the product of the coordinate value and the corresponding gray value according to the formulas (2) and (3):

vector X in the formula₁Expressed is the product of the abscissa 1 and all the elements of the row 1 of the matrix (image grey values), and so on, the vector X₂～X_NOf the representationIs the product of the abscissa and all its corresponding elements. Vector Y₁Representing the product of ordinate 1 and all the elements of column 1 of the matrix, and so on, vector Y₂～Y_NAll that is represented is the product of the ordinate and all the elements to which it corresponds.

(2) Will vector X₁～X_NAll elements in the formula are added to obtain the formula sigma_x,yThe result of calculation of xI (X, y) is sun _ X_iWill vector Y₁～Y_NAll elements in the formula are added to obtain the formula sigma_x,yThe calculation result of yI (x, Y) is sun _ Y_i。

(3) Adding all the elements (image gray values) of the matrix to obtain a formula sigma_x,yThe calculation result of I (x, y) is sun _ q.

(4) Performing division operation twice by using division IP core, wherein dividend is sun _ X_iThe divisor sun _ q obtains the centroid coordinate x₀Dividend is sun _ Y_iThe divisor sun _ q obtains the barycenter coordinate y₀。

3. Calculation of the molecular weight of MR

In the formula, I (x, y) represents a gray scale value of (x, y). The flow chart of the performance indicator function MR molecule stepwise calculation is shown in FIG. 2, a part of the time sequence chart is shown in FIG. 3, and the specific calculation process is as follows:

(1) another clock SCLK is generated from the clock Clk, which is a division of the clock Clk by N, and the relationship between the two is shown in fig. 3, where 1 cycle of the SCLK clock contains N Clk clock cycles.

(2) And designing a y _ cnt counter by taking the Clk as a time sequence, wherein the counter starts counting from 0, the rising edge of the Clk triggers, the counter is increased by 1, and when the count reaches N, the y _ cnt starts to count again from 0. With SCLK as a timing sequence, an x _ cnt counter is designed, the counter starts counting from 0, and stops counting when counting to N. The partial counting process of two counts is shown in x _ cnt and y _ cnt of fig. 3. The two counters are mainly used to locate the matrix elements (image gray values), x _ cnt locates the abscissa, y _ cnt locates the ordinate, x of the matrix equals to the value of x _ cnt plus 1, y of the matrix equals to the value of y _ cnt plus 1, e.g. x _ cnt equals to 0, y _ cnt equals to 1 indicates the element of row 1 and column 2 of the matrix.

(3) Calculating (x-x)₀)²-(y-y₀)²First, calculate x-x₀，y-y₀Then according to the formula (x-x)₀)×(x-x₀)-(y-y₀)×(y-y₀) The calculation result is delayed by 1 clock cycle, for example, x equals 1, y equals 3 (corresponding to the positions x _ cnt equals 0 and y _ cnt equals 2), and the calculation result is output at the position where x _ cnt equals 0 and y _ cnt equals 3 bits. The calculation results of all positions in the 1 st row of the matrix are output between the two dotted lines A and B.

(4) Computing

And (4) inputting the result calculated in the step (3) into an open IP core, and delaying for one clock cycle to output an operation result sqr _ q. And (4) integrating the calculation processes of (3) and (4), and delaying the clock by 2 clock cycles. If x equals 1, y equals 1 (corresponding to the positions x _ cnt equals 0 and y _ cnt equals 0), the result of the calculation is output at the position where x _ cnt equals 0 and y _ cnt equals 2. And the calculation results of all coordinates in the 1 st row of the matrix are output between the C and D dot-and-dash lines.

(5) According to the timing sequence of the clock Clk, a counter address is set, an enable signal for starting counting is consistent with the y _ cnt counter, and a partial timing diagram of the address counter is shown in FIG. 3. The function of the address counter is to read the image grey values from the ROM core. From the issuance of the address bit address, there is a delay of two clock cycles to the output data, and the image gradation value of 0, such as address, is read at a position where the address is equal to 2. Due to the relationship of the address counter and the y _ cnt counter, address, etc. 0 is positioned where x _ cnt equals 0, y _ cnt equals 0, and the output gray values are positioned where x _ cnt equals 0 and y _ cnt equals 2. The specific correspondence is shown in fig. 3.

(6) From FIG. 3, the formula can be seen

The result of the calculation of (c) and the gray value with coordinates (x, y) are two completeSynchronous signals, the two results can be directly multiplied to complete the formula

The calculation process of (2). The result is noted as vector sqr _ x₁。

(7) All the results (sqr _ x) from line 2 to line N can be calculated according to the above calculation process₂～sqr_x_N) The result of the addition of the results of the elements within all vectors is then the numerator of the MR.

4. And (4) dividing the result of the calculation in the step (7) in the step (3) by the sum of all the elements (image gray values) of the matrix to obtain the MR.

Calculating control voltage of the three deformable mirrors:

the calculation formula of the control voltage of the deformable mirror is u⁽ⁿ⁺¹⁾＝uⁿ-γΔuⁿΔJ⁽ⁿ⁾In the formula, gamma delta uⁿΔJ⁽ⁿ⁾For gradient estimation, γ is the gain factor, Δ uⁿFor random perturbations, Δ J⁽ⁿ⁾Is the variation of the performance index

The performance index function MR is to search the optimal solution in the direction of minimum value, so when the value of delta u isⁿAnd Δ J⁽ⁿ⁾When the signs of the two signals are consistent, the gradient estimation result is subtracted, and when the sign of the two signals is consistent, the value of the gradient is increasedⁿAnd Δ J⁽ⁿ⁾When the signs of (1) are opposite, the gradient estimation result is added. The invention realizes the calculation of the control voltage of the deformable mirror according to the formula (4).

u⁽ⁿ⁺¹⁾＝uⁿ-γΔuⁿΔJ⁽ⁿ⁾ (4)

In the formula: gamma is the gain factor taking a positive value. The calculation module of the control voltage is shown in fig. 4, and the specific implementation steps are as follows:

1. after positive disturbance is applied to the deformable mirror, image data are collected, a performance index function MR is calculated by using a performance index function calculation module and is recorded as MR₁. After negative disturbance is applied to the deformable mirror, image data is collected, and a performance index function MR is calculated by using a performance index function calculation module₂。

2. ComputingVariation Δ MR of performance indicator function is MR₁-MR₂。

3. The gradient is calculated.

4. Calculating the control voltage according to formula (4), and converting the formula into u for realizing convenience⁽ⁿ⁺¹⁾＝uⁿ-Asign(Δuⁿ)ΔJ⁽ⁿ⁾Where A denotes the product of the gain factor and the random perturbation amplitude (positive value), sign (Δ u)ⁿ) Equal to plus or minus 1, and d in the timing diagram of the random perturbation in FIG. 1 represents the different timing sign (Δ u)ⁿ) The value of (c).

5. And after the control voltage is calculated, adding the (n + 1) th random disturbance voltage, and then outputting the control voltage to the deformable mirror to acquire image data. And subtracting the random disturbance voltage of the (n + 1) th time, then outputting the control voltage to the deformable mirror, and collecting the image data again.

FIG. 5 shows a timing diagram of the control voltage calculation process, where the MR _ Done signal is high, indicating that the performance indicator function module has completed its calculation, and after the MR _ Done signal is triggered, the next 1 rising edge of the clock Clk outputs MR₁ MR₂Synchronously computing MR₁And MR₂The difference Δ MR of (a). SJU _ en enable signal triggers the random perturbation module to output d signal inside the random perturbation module. In order to determine the position of the d signal output, a counter cnt of 0-98 is designed. It can be seen from fig. 5 that the signal d is output when cnt is equal to 3, and after the cnt counting is completed, the output is delayed by 1 clock cycle, and there are 97 d outputs in total. The preprocessing part completes the product A of the random disturbance amplitude and the gain coefficient. The T _ D signal (gradient) is the product of Δ MR, D, A. Synchronously finishing the calculation process of subtracting the gradient from the current voltage to obtain u⁽ⁿ⁺¹⁾And storing. SJU _ en is also the trigger signal of the next random disturbance module, the random disturbance dv starts to be output at the position where cnt is equal to 3, after cnt counting is completed, the output is finished by delaying 1 clock cycle. Will control the voltage u⁽ⁿ⁺¹⁾And adding random disturbance to the deformable mirror and subtracting the random disturbance to obtain the control voltage to be output by the iteration.

The FPGA-based wave-front-free detection adaptive optics SPGD control algorithm experiment system comprises:

the structural block diagram of the system is shown in fig. 6, the imaging device (such as cameras with various interfaces and image acquisition systems with COMS and CCD photosensitive chips) acquires light spot image data, and the acquired image data is transmitted and stored in the FPGA development board through data lines (such as USB and optical fibers). The SPGD control algorithm realized by FPGA processes the image data, calculates the control voltage of the wave-front corrector, and transmits the calculated control voltage to the wave-front corrector (such as deformable mirrors and spatial light modulators with various unit numbers) through data lines (depending on interfaces of the wave-front corrector for receiving data) such as USB or optical fibers, so as to achieve the purpose of controlling the wave-front corrector to correct the speckle image.

Claims

1. An implementation method of an adaptive optics SPGD control algorithm based on an FPGA is characterized by comprising the following steps:

(1) a random disturbance voltage generation module: generating a plurality of groups of pseudo-random sequences, multiplying each group of sequences by disturbance voltage in the FPGA and adding corresponding signs to generate positive disturbance and negative disturbance;

(2) and (3) calculating a performance index: selecting the average Radius MR-Mean Radius of the far-field light spot as a performance index function,

wherein (x)₀,y₀) Is the centroid of the spot;

(3) And a deformable mirror voltage calculation control module: according to the formula u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-γΔu⁽ⁿ⁾ΔJ⁽ⁿ⁾Calculating control voltage, optimizing the performance index function MR in the direction of a minimum value, and taking a positive value of gamma;

the SPGD control algorithm is realized by adopting the following steps:

the module for generating the random disturbance voltage comprises the following steps:

1.1. the preprocessing part generates a plurality of groups of pseudo-random sequences by utilizing a rand function by means of Matlab software, and the sequences are formed by combining 0 and 1; storing the generated plurality of groups of 97 paths of pseudo-random sequences for later use;

1.2. calculating positive and negative disturbances by taking one group of pseudo-random sequences as an example, and reading and storing the group of pseudo-random sequences into a ROM core of the FPGA; setting a counter cnt, reading the pseudorandom sequence stored in the ROM by taking the cnt as an address bit, and when the cnt is equal to 0, corresponding to the 1 st value of the plurality of paths of pseudorandom sequences, and repeating the steps;

1.3. setting Q in the corresponding relation between the numerical value of the cnt of the counter and the plurality of paths of pseudo-random sequences to represent a random sequence read from the ROM, wherein when the cnt is equal to 0, the corresponding Q is output when the cnt is equal to 2;

1.4, when Q is equal to 0, the corresponding disturbance is 1, when Q is equal to 1, the corresponding disturbance is-1, and the number of pseudo-random sequences is 1 or-1; because the value of Q is judged firstly, the corresponding value of d is delayed, and d between two dot-dash lines is determined to be all the required paths of random disturbance with the disturbance amplitude of 1 by means of a certain auxiliary signal;

1.5. multiplying d by 200 on the basis of the step 1.4 to obtain dv1, namely positive disturbance; 200 is the disturbance amplitude, which is defined as the binary fixed point decimal with the last 10 digits being decimal digits, so the actual disturbance amplitude corresponding to 200 is 0.1953;

1.6. adding a negative sign to each disturbance on the basis of dv1 to obtain data dv2, namely negative disturbance;

and (II) calculating a performance index:

2.1. acquiring image data, wherein the gray value of a light spot image is N multiplied by N, and the ROI area is N multiplied by N, so that the matrix recovery is required before calculation because the image data is transmitted to the interior of the FPGA in series; restoring the image data into an N multiplied by N matrix form by taking the address bit as a judgment condition;

2.2. calculating the centroid, the formula of the centroid is

(2.2.1) calculating the product of the coordinate value and the corresponding gray value according to the formula (2) and (3):

vector X in the formula₁Expressed is the product of the abscissa 1 and all the elements of the 1 st row of the matrix, i.e. the image grey values, and so on, the vector X₂～X_NAll the expressions are products of the abscissa and all the elements corresponding to the abscissa; vector Y₁Representing the product of ordinate 1 and all the elements of column 1 of the matrix, and so on, vector Y₂～Y_NAll the expressions are products of the ordinate and all the elements corresponding to the ordinate;

(2.2.2) vector X₁～X_NAll elements in the formula are added to obtain the formula sigma_x,yThe result of calculation of xI (X, y) is sun _ X_iWill vector Y₁～Y_NAll elements in the formula are added to obtain the formula sigma_x,yThe calculation result of yI (x, Y) is sun _ Y_i；

(2.2.3) adding all the elements of the matrix to obtain the formula sigma_x,yThe calculation result of I (x, y) is sun _ q;

(2.2.4) carrying out two division operations by using the division IP core, wherein the dividend is sun _ X_iThe divisor sun _ q obtains the centroid coordinate x₀Dividend is sun _ Y_iThe divisor sun _ q obtains the barycenter coordinate y₀；

2.3. Calculation of the molecular weight of MR

In the formula, the coordinates represented by I (x, y) are gray values of (x, y), and the specific calculation process of the performance index function MR molecule step-by-step calculation is as follows:

(2.3.1) generating another clock SCLK from the clock Clk, SCLK being a N-division of the clock Clk, 1 cycle of the SCLK clock comprising N Clk clock cycles;

(2.3.2) designing a y _ cnt counter by taking the Clk as a time sequence, starting counting by the counter from 0, triggering the rising edge of the Clk, adding 1 by the counter, and when counting to N, starting recounting by the y _ cnt from 0; designing an x _ cnt counter by taking SCLK as a time sequence, wherein the counter starts counting from 0, and stops counting when counting to N; the two counters are mainly used for positioning matrix elements, namely image gray values, x _ cnt positions an abscissa, y _ cnt positions an ordinate, the abscissa x of the matrix is equal to the value of x _ cnt plus 1, and the ordinate y of the matrix is equal to the value of y _ cnt plus 1;

(2.3.3) calculating (x-x)₀)²-(y-y₀)²First, x-x is calculated₀，y-y₀Then according to the formula (x-x)₀)×(x-x₀)-(y-y₀)×(y-y₀) Calculating, wherein the operation result can be output by delaying 1 clock cycle;

(2.3.4) calculation

Inputting the result calculated in the step (3) into an open party IP core, and delaying one clock cycle to output an operation result sqr _ q; the (2.3.3) and (2.3.4) calculation processes are integrated, and the clock is delayed by 2 clock cycles;

(2.3.5) setting a counter address according to the time sequence of the clock Clk, wherein an enabling signal for starting counting is consistent with a y _ cnt counter, and the address counter is used for reading an image gray value from a ROM core; from the sending of the address bit address, the delay of two clock periods is provided to the output data, and the image gray value of 0 such as the address can be read at the position where the address is equal to 2; due to the relationship between the address counter and the y _ cnt counter, the location of 0 for address, etc. is x _ cnt equal to 0, y _ cnt equal to 0, the gray value output is located where x _ cnt equal to 0, y _ cnt equal to 2;

(2.3.6) formula

The result of the calculation and the gray value with coordinates (x, y) are two completely synchronous signals, and the two results can be directly multiplied to complete the formula

The calculation process of (2); the result is noted as vector sqr _ x₁(ii) a (2.3.7) all the results sqr _ x from line 2 to line N can be calculated according to the above calculation procedure₂～sqr_x_NThen, the element results in all vectors are added to obtain the result of MR molecules;

2.4. dividing the result calculated in the step (2.3.7) by the sum of all elements of the matrix, namely the image gray value, and calculating to obtain MR;

(III) calculating the control voltage of the deformable mirror:

the calculation formula of the control voltage of the deformable mirror is u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-γΔu⁽ⁿ⁾ΔJ⁽ⁿ⁾In the formula, gamma delta u⁽ⁿ⁾ΔJ⁽ⁿ⁾For gradient estimation, γ is the gain factor, Δ u⁽ⁿ⁾For random perturbations, Δ J⁽ⁿ⁾Is the variation of the performance index

The performance index function MR is to search the optimal solution in the direction of minimum value, so when the value of delta u is⁽ⁿ⁾And Δ J⁽ⁿ⁾When the signs of the two signals are consistent, the gradient estimation result is subtracted, and when the sign of the two signals is consistent, the value of the gradient is increased⁽ⁿ⁾And Δ J⁽ⁿ⁾When the signs of the two signals are opposite, adding a gradient estimation result;

and (3) realizing the calculation of the control voltage of the deformable mirror according to the formula (4):

u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-γΔu⁽ⁿ⁾ΔJ⁽ⁿ⁾ (4)

in the formula: gamma is a positive value of the gain coefficient;

the method comprises the following specific steps:

3.1. applying positive perturbations to deformable mirrorsThen, collecting image data, calculating a performance index function MR by using a performance index function calculation module, and recording the performance index function MR as MR₁(ii) a After negative disturbance is applied to the deformable mirror, image data is collected, and a performance index function MR is calculated by using a performance index function calculation module₂；

3.2. Calculating the variation Δ MR of the performance indicator function as MR₁-MR₂；

3.3. Calculating a gradient;

3.4. calculating the control voltage according to formula (4), and converting the formula into u⁽ⁿ⁺¹⁾＝u⁽ⁿ⁾-Asign(Δu⁽ⁿ⁾)ΔJ⁽ⁿ⁾Where A denotes the product of the gain factor and the amplitude of the random disturbance, sign (Δ u)⁽ⁿ⁾) Equal to plus or minus 1;

3.5. after the control voltage is calculated, adding the random disturbance voltage for the (n + 1) th time, then outputting the control voltage to the deformable mirror, and collecting image data; and subtracting the random disturbance voltage of the (n + 1) th time, then outputting the control voltage to the deformable mirror, and collecting the image data again.