Multi-frame blind restoration method for images corrected by adaptive optical systems

Abstract

The present invention discloses a multi-frame blind restoration method for an adaptive optical system corrected image, which specifically comprises the following steps: step 1: using a gradient filter to extract high-frequency information of an f-th frame degraded image g to obtain the extracted high-frequency information; establishing an objective function based on the extracted high-frequency information to estimate a point spread function corresponding to the f-frame degraded image; dividing the objective function into a first sub-problem and a second sub-problem based on an alternating iterative algorithm; using an iterative algorithm to perform a k-th iterative calculation on the first sub-problem to obtain a high-frequency restored image x k , based on x k The iterative algorithm is used to perform the k-th iteration calculation on the second subproblem to obtain the point spread function h k ; According to the point spread function obtained by calculation, the f-th frame original image is restored by using a multi-frame image restoration algorithm containing Poisson noise. The restored image of the present invention has better visual effects and evaluation indicators.

Description

Technical Field

[0001] The invention belongs to the technical field of image processing. Background Art

[0002] Since Galileo first used a telescope for astronomical observation in 1609, human observation of space targets has made great progress, among which ground-based large-aperture optical telescopes are one of the most important observation tools. The main performance of optical telescopes includes light-gathering power and resolution, both of which increase with the increase of the effective aperture of the telescope, so large-aperture telescopes can observe darker targets and can resolve more details of the target. However, for ground-based telescopes, there are various static and dynamic aberrations inside and outside the system, which seriously affect the imaging quality of the telescope. Among them, static aberrations mainly include the aberrations of the optical system itself, and dynamic aberrations mainly refer to structural deformation caused by temperature, gravity, etc. and wavefront distortion caused by atmospheric turbulence. These aberrations reduce the resolution of large-aperture telescopes to the level of small telescopes with an aperture of 20 cm, which greatly limits the observation capabilities and applications of large-aperture telescopes.

[0003] In order to solve this problem, adaptive optics came into being in the 1950s. The adaptive optical system uses a wavefront detector to detect the wavefront phase in real time, and controls the wavefront corrector to apply an anti-phase wavefront signal to compensate for the wavefront distortion and improve the imaging quality of the system. The adaptive optical system greatly reduces the impact of atmospheric turbulence on the imaging system of large-aperture telescopes, making it possible to use large-aperture telescopes for high-resolution observations. However, the current adaptive optical system is not perfect. Due to wavefront detection errors, corrector reconstruction wavefront errors, control cycle delays, and insufficient effective bandwidth, it cannot perfectly correct aberrations. The image obtained after correction still has correction residuals, especially when the atmospheric turbulence is strong, it is difficult to achieve high-resolution imaging, and the target image obtained at this time still has a lot of room for improvement.

[0004] The RL algorithm only processes a single image and requires the known point spread function of the image. However, when the traditional blind restoration algorithm that estimates the point spread function is applied to the image corrected by the adaptive optical system, it usually converges to the Delta function and cannot accurately estimate the point spread function. Therefore, it is impossible to restore and improve the quality of the image corrected by the adaptive optical system. Summary of the invention

[0005] Purpose of the invention: In order to solve the problems existing in the above-mentioned prior art, the present invention provides a multi-frame blind restoration method of an adaptive optical system corrected image.

[0006] Technical solution: The present invention provides a multi-frame blind restoration method for an adaptive optical system corrected image, which specifically includes the following steps:

[0007] Step 1: Use a gradient filter and Extract high-frequency information from the f-th frame of degraded image g to obtain the extracted high-frequency information T is the matrix transpose, and the degraded image g is the degraded image corrected by the adaptive optical system;

[0008] Step 2: Establish the following objective function to estimate the point spread function h corresponding to the f-th frame degraded image:

[0009]

[0010] Among them, x is the high-frequency information part of the original image, h i is the value of the i-th element of h, and h i ≥0,∑ i h i =1; λ and δ are regularization weight variables, * represents convolution operation;

[0011] Step 3: Divide the objective function of step 2 into the first and second sub-problems based on the alternating iterative algorithm;

[0012] The first sub-problem is:

[0013]

[0014] The second sub-problem is:

[0015]

[0016] Step 4: Use the iterative algorithm to perform the k-th iteration calculation on the first sub-problem to obtain the high-frequency restored image x at the k-th iteration calculation k , based on the solved x k The iterative algorithm is used to perform the k-th iteration calculation on the second subproblem, and the point spread function h at the k-th iteration calculation is obtained. k ,k=1,2,…K, where K is the total number of iterative calculations;

[0017] Step 5: Based on the point spread function obtained by calculation, use the multi-frame image restoration algorithm containing Poisson noise to restore the original image of the f-th frame, and determine whether f is less than NF. If so, then f+1, go to step 1 to restore the original image of the f+1th frame, otherwise stop the calculation.

[0018] Furthermore, in step 4, the first sub-problem is converted into the following form based on the Bregman iterative algorithm:

[0019]

[0020] Where b1 and b2 are Bregman variables, d1 = Hx, d2 = x, λ' = λ||x l ||2, where H is the matrix form of h after Fourier transformation, γ is the regularization weight variable; x l is the high-frequency restored image obtained in the last iterative calculation when the above formula 4 is iteratively calculated;

[0021] When using the Bregman iterative algorithm to solve Formula 4, Formula 4 is divided into the following forms for solution:

[0022]

[0023] Where I is the identity matrix.

[0024] Furthermore, when using an iterative algorithm to solve Formula 3, Formula 3 is divided into the following formulas for solution:

[0025]

[0026] Where N is the total number of elements in h; div(.) is the divergence calculation, is the gradient operator.

[0027] Furthermore, the step 5 is specifically as follows:

[0028] Step 5.1: Establish the following objective function:

[0029]

[0030] in, represents the estimated value of the original image of the fth frame, u f represents the actual f-th frame of the original image, f is the image frame number, and NF is the total number of image frames;

[0031] Step 5.2: Use the following RL iterative algorithm to solve Equation 7:

[0032]

[0033] Among them, j represents the jth iteration calculation, Indicates that the original image is obtained by the jth iteration calculation; when the first iteration calculation is performed Where J represents the maximum number of iterations of the RL iterative algorithm;

[0034] Step 5.3: Use the BM3D algorithm as a low-pass filter to perform low-pass filtering on the solved original image to obtain the final original image.

[0035] Beneficial effects: The image restored by the present invention has better visual effects and evaluation indicators, and the point spread function calculated by the present invention is closer to the real point spread function. The image restored by the present invention has better noise removal ability and clearer edges. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1 It is a flow chart of the present invention;

[0037] Figure 2 are the original experimental images and simulated degraded images with different noise levels; (a) is the original image; (b) is the degraded image when the noise level k1=1000; (c) is the degraded image when k1=100; (d) is the degraded image when k1=40;

[0038] Figure 3 The restored images of each algorithm when k1=1000; (a) is the restored image based on the Framelet algorithm; (b) is the restored image based on the OBD algorithm; (c) is the restored image based on the method proposed in the present invention;

[0039] Figure 4 The restored images of each algorithm when k1=100; (a) is the restored image based on the Framelet algorithm; (b) is the restored image based on the OBD algorithm; (c) is the restored image based on the method proposed in the present invention;

[0040] Figure 5 The restored images of each algorithm when k=40; (a) is the restored image based on the Framelet algorithm; (b) is the restored image based on the OBD algorithm; (c) is the restored image based on the method proposed in the present invention;

[0041] Figure 6 Point spread functions estimated by different methods: (a) is the real point spread function; (b) is the point spread function obtained based on the Framelet algorithm, RMSE = 0.742; (c) is the point spread function obtained based on the OBD algorithm, RMSE = 0.720; (d) is the point spread function obtained based on the present invention, RMSE = 0.376;

[0042] Figure 7 This is a graph showing the change of RMSE of the restored image with the number of image frames;

[0043] Figure 8 The RMSE of the restored images of the OBD algorithm and the algorithm of the present invention varies with the number of image frames when k1=1000;

[0044] Fig. 9 This is the optical path diagram of the LC-AO system;

[0045] Fig.10 The original and restored star images; Figure (a) is the original image collected; (b) is the image restored based on the OBD algorithm; (c) is the image restored using the method of the present invention;

[0046] Fig.11 This is a comparison of the height and width of the whole class for different images;

[0047] Fig.12 The ALOS satellite image and the restored image are collected; Figure (a) is the original image collected; Figure (b) is the image restored based on the OBD algorithm; Figure (c) is the image restored using the method of the present invention;

[0048] Fig.13 The captured COSMOS1461 satellite image and the restored image; Figure (a) is the captured original image; Figure (b) is the image restored based on the OBD algorithm; and Figure (c) is the image restored using the method of the present invention. DETAILED DESCRIPTION

[0049] Astronomical image data is usually collected by CCD cameras, and its statistical model can usually be expressed as:

[0050] g=g p +n g (9)

[0051] Among them, g p is the number of photoelectrons generated by the object and background emission. It is a Poisson random variable with the expected value of h*u+b, where h is the point spread function of the imaging system, u is the original image, * is the convolution operation, and b is the background term. g represents the readout noise, which can be approximated as Gaussian noise, and g is the degraded image acquired by imaging.

[0052] A classic method in the field of Poisson noise image restoration is the Richardson–Lucy algorithm [DAFish, AMBrinicombe, ERPike, and JG Walker, "Blind deconvolution by means of the Richardson–Lucy algorithm," J.Opt.Soc.Am.A 12(1), 58-65(1995)] (hereinafter referred to as the RL algorithm). The RL algorithm is developed from the Bayesian theory, and its theory can be expressed as follows:

[0053]

[0054] Among them, P(u) represents the probability of event u, and P(u|g) represents the conditional probability. Using the maximum likelihood estimation (MLE) to solve equation (2) converts the problem into a problem of minimizing the energy equation:

[0055]

[0056] Where E represents the energy function, represents u when the energy equation takes the minimum value, Indicates the value of u when the energy equation E takes the minimum value;

[0057] E(u)=h*u-glog(h*u) (12)

[0058] Converting Equation 4 into iterative form, the RL algorithm can be expressed as follows:

[0059]

[0060] Where j is the number of iterations. The point spread function h of the image needs to be known, but in the restoration of the image after the adaptive optical system correction, it is very difficult to know the true point spread function, that is, the wavefront distortion residual. Therefore, in practical applications, optimization algorithms such as minimum mean square error are usually used to estimate the point spread function. The image restoration problem in which both the original image and the point spread function are unknown is also called the blind restoration problem.

[0061] Existing image blind restoration algorithms usually use the L1 norm or L2 norm of the point spread function as a constraint when estimating the point spread function. Therefore, the blind restoration algorithm can be expressed as follows:

[0062]

[0063] However, if only the L1 norm or L2 norm is used as a constraint, for astronomical images with fewer details and more noise, the final estimated point spread function will approach the Delta function, and it is impossible to obtain an accurate point spread function to estimate the restored image.

[0064] like Figure 1 As shown, in view of the above problems, this embodiment proposes a multi-frame image blind restoration method that is more effective in estimating the point spread function, thereby improving the estimation accuracy of the point spread function, and the ultimate goal is to improve the image quality after the adaptive optical system is corrected. The specific method is:

[0065] In order to solve the problem of inaccurate estimation of the point spread function of the image corrected by the adaptive optical system, if the original image is used as a constraint to restore the point spread function, the change of the constraint is not obvious during the iterative solution process, which is not conducive to the estimation of the point spread function. Therefore, this embodiment first uses only the high-frequency part of the image to be restored to estimate the point spread function. That is, gradient filtering is first used and Extract high-frequency information from images T is the matrix transpose. g is the f-th frame degraded image acquired by imaging. If the noise is ignored, the point spread function estimated using only the high-frequency information of the image is unchanged. Therefore, y is directly used to estimate the point spread function h of the f-th frame degraded image. The objective function of the algorithm proposed in the present invention can be expressed as follows:

[0066]

[0067] x is the high-frequency information of the original image. When this objective function takes the minimum value, the value of h is the estimated point spread function. In the objective function, And i h i =1,h i is the value of the i-th element of the point spread function. The first term of the objective function is the similarity constraint of the Poisson noise model; the second term is the regularization term of the restored image, where the ratio of the L1 norm to the L2 norm of the restored image is used to ensure the scale-invariant sparsity of the restored image; the third term is the total variation constraint of the point spread function, the purpose of which is to reduce noise and maintain discontinuity. The weight variables λ and δ are regularization parameters that control the relative strength of the coordination point spread function and image constraint terms.

[0068] Equation 15 is non-convex and difficult to solve directly. We can use an alternating iterative algorithm to update the values ​​of x and h alternately. Therefore, Equation 15 can be divided into the following two sub-problems:

[0069]

[0070] Because the regularization term The first subproblem is difficult to solve directly, so the value of ||x||2 is set to the value after the previous iteration in each iteration, and the first subproblem is converted to an L1 norm solution problem. Using the Bregman iterative algorithm, the unconstrained minimization problem is converted to a constrained minimization problem by adding additional variables d1=Hx and d2=x. Among them, H is the matrix form of h after Fourier transformation. Therefore, the first subproblem is converted to the following form:

[0071]

[0072] where λ'=λ||x l ||2,xl is the high-frequency restored image obtained in the last iterative calculation when the above formula is iterated; γ is the regularization parameter, and b1 and b2 are Bregman variables.

[0073] The solution of Equation 17 can be divided into four minimization problems, minimizing x, d1 and d2 respectively. The iterative update method is as follows:

[0074]

[0075] Where I is the identity matrix, the superscript T represents the transpose operation of the matrix, and k is the number of iterations.

[0076] The second subproblem in Eq. 16 can be solved using the following iterative algorithm:

[0077]

[0078] Where N is the total number of elements in the point spread function; div(.) is the divergence calculation, is the gradient operator.

[0079] According to the process of point spread function estimation in the technical solution, it is summarized as the following algorithm 1:

[0080] Algorithm 1

[0081] Input: degraded image g, point spread function size sh

[0082] Output: point spread function h

[0083] Apply a gradient filter to the degraded image g to obtain a high-frequency degraded image y

[0084] For iterations k = 1, 2, 3…

[0085] Use formula (18) to restore the high-frequency image x k Make an update

[0086] Using equation (19) for the point spread function h k To update,

[0087] End

[0088] Return the point spread function h after the last cycle update

[0089] In terms of image restoration algorithm, the RL algorithm is extended to a multi-frame image version. Since adjacent multi-frame images are imaged in a very short period of time, it can be assumed that the original image u does not change during this period of time. Therefore, the multi-frame image restoration algorithm containing Poisson noise can be summarized as follows 20:

[0090]

[0091] in, The original image of frame f represents the estimated value, u f represents the actual f-th frame of the original image, f is the image frame number, and NF is the total number of image frames;

[0092] The RL iterative algorithm for solving Equation 20 is as follows (13):

[0093]

[0094] j is the number of iterations. In this embodiment, all frames in the image sequence are used to update the restored image at the same time, and one frame of image is used to update the restored image in each iteration. The advantage of this is that it can greatly save computing resources and computer memory, and the restored image can be obtained while the image is collected. The optimized RL iterative algorithm is as follows (14):

[0095]

[0096] Among them, j represents the jth iteration calculation, Indicates that the original image is obtained by the jth iteration calculation; when the first iteration calculation is performed Where J represents the maximum number of iterations of the RL iterative algorithm;

[0097] Combining Algorithm 1 with the above formulas (12) to (14), the multi-frame image restoration algorithm with Poisson noise proposed in this embodiment is summarized as the following Algorithm 2.

[0098] Algorithm 2

[0099] Input: f-th frame degraded image g (k ≥ 1), point spread function size sh

[0100] Output: restored image u f

[0101] For image frame sequence f=1,2,3…

[0102] Estimate the point spread function h using Algorithm 1

[0103] make

[0104] For iteration number j = 1, 2, 3...

[0105]

[0106] End

[0107] Will Through the denoising low-pass filter, the restored image u is obtained f

[0108] End

[0109] Return the u obtained in the last iteration f

[0110] When using Algorithm 2 to calculate the first frame of the original image, set the initial original image u 0 ,u 0 is the first degraded image.

[0111] In the outer loop of this embodiment, one frame of degraded image is used for processing each time, and each loop is divided into the following three steps: 1. Use algorithm 1 to estimate the point spread function of the current degraded image; 2. Use the RL iterative algorithm to obtain the deblurred image; 3. The deblurred image is subjected to low-pass filtering to obtain the final restored image. The BM3D algorithm is used as a low-pass filter to remove noise in the algorithm of this patent.

[0112] Experimental Proof

[0113] First, the verification is performed on the image corrected by the simulated adaptive optical system with artificially controlled parameters. The target image used is the 300×300 OCNR5 satellite model. This image has a high dynamic range and more details, so it is often used as a test image in astronomical image processing literature. Since there is no way to obtain the real corrected wavefront residual, this embodiment uses Kolmogorov phase theory to generate a set of random wavefronts and uses them as the corrected wavefront residual of the adaptive optical system. This embodiment generates 200 point spread functions with a D / r0 parameter value of 12 to simulate the usual environment of a 1.2m telescope in an atmospheric coherence length of r0=10cm. After the grayscale value of the image blurred by the simulated point spread function is stretched to 0-1, Poisson noise is added to simulate the detector photon shot noise. In order to simulate different imaging environments, this embodiment adds different levels of Poisson noise in the experiment for experiment. Subsequently, this embodiment re-stretches the grayscale value of the degraded image to 0-1, and adds Gaussian white noise with a mean of 0 and a variance of 0.01 to simulate the CCD readout noise. The original and degraded images are as follows Figure 2 shown.

[0114] In this experiment, this embodiment uses a point spread function of size 41*41, RL algorithm 1 iteration, and point spread function blind estimation algorithm 1 iteration. The σ parameter in the BM3D filter needs to be determined according to the noise level of the degraded image and the number of RL iterations. In the experiment of this embodiment, when the Poisson noise parameter k1=1000, σ is set to 2, and when the Poisson noise parameter k1=100 and k1=40, σ is set to 5. In the point spread function blind estimation algorithm part, λ=80, λ=20 and σ=0.03 are set. The restored image results of different algorithms and different noise levels are shown in the figure. Figure 3-Figure 5 As shown, the estimated point spread function is Figure 6 As shown, from Figure 3-6 It can be seen that in the comparison of different restoration methods, although the results of the OBD algorithm seem to have more details, most of the pixel values ​​on the sailboard have reached saturation, which is quite different from the pixel values ​​of the original image. Therefore, the RMSE and SSIM indicators are lower than those of the method in this paper. In addition, the lower RMSE also results in a lower peak signal-to-noise ratio (PSNR). The results of the OBD algorithm contain obvious noise, which affects the visual perception, especially when the sequence image is presented. The Framelet-based method has better denoising ability, but the results are seriously lacking in detail information. The results of the algorithm in this paper have better evaluation indicators and are better than the compared algorithms. In addition, the visual effect has achieved a balance between denoising and retaining details, with no obvious noise, and the target pixel value is close to the original image. At the same time, the detailed information on the satellite and surrounding antennas is retained.

[0115] In terms of point spread function restoration, the results of the OBD algorithm tend to be anisotropic Gaussian functions, which are quite different from the original point spread function. The results based on the Framelet algorithm are also quite different from the original point spread function, and the main part in the middle is not accurately estimated. The results of the algorithm in this paper are closer to the real point spread function, and the RMSE index of the point spread function also proves the above conclusion.

[0116] This embodiment uses the root mean square error (RMSE) and structural similarity (SSIM) to evaluate the effect of the restored simulated image. The RMSE and SSIM results are shown in Table 1.

[0117] Table 1

[0118]

[0119] The restored images and evaluation indicators show that the quality of the restored images decreases as the noise level of the degraded images increases. In the comparison of different restoration methods, the OBD algorithm appears to have more details, but it also causes ringing effects and pixel value saturation, so the RMSE and SSIM are lower. The Framelet-based method lacks details. The restored image of the algorithm proposed in this embodiment has better visual effects and evaluation indicators, which are better than the compared algorithms. In addition, the estimated point spread function of the algorithm proposed in this patent is closer to the real point spread function.

[0120] Due to the stochastic gradient descent strategy in the algorithm, when this embodiment is run frame by frame, the RMSE may not be strictly monotonically decreasing. The RMSE of the restored image in this embodiment changes with the number of processed image frames in the image sequence as shown in the figure. Figure 7 As shown in FIG. 1 , the RMSE of the restored image of this embodiment and the OBD algorithm is compared with the change of the number of processed image frames in the image sequence. Figure 8 As shown in the figure, it can be seen that the RMSE of the restored image of this embodiment is not strictly monotonically decreasing, but it is continuously decreasing overall, and will eventually "converge" to a range. However, the RMSE of the restored image of the OBD algorithm does not continue to decrease, and even gradually increases, and does not "converge".

[0121] In addition, this embodiment also applies the algorithm proposed in this patent to astronomical images corrected by a real adaptive optical system. The image was taken in Changchun, China using a 1.2m telescope. The adaptive optical system consists of a Shack-Hartmann wavefront detector and two liquid crystal wavefront correctors. The main parameter indicators of the imaging system are atmospheric coherence length of about 6cm, Greenwood frequency of 58Hz, telescope aperture of 1.23 meters, focal length of 84 meters, detection wavelength of 400-700nm, imaging wavelength of 700-900nm, CCD pixel size of 48μm, and diffraction limit of 0.3. The main parameter indicators of the adaptive optical system are that the Hartmann detector has 225 (15×15) microlenses, 90×90 pixels and 870Hz frame rate; the liquid crystal wavefront corrector has 256×256 pixels, a response time of 1.15 milliseconds, and wavefront reconstruction of 35 Zernike modes. The optical path diagram and the actual diagram of the adaptive optical system are as follows: Fig. 9 The experimental images were taken in Changchun, China in October.

[0122] Real star imaging experiment, the results of this embodiment are compared with the results of the OBD algorithm. Fig.10 The restored image of the OBD algorithm still has the problems of pixel saturation and severe ringing effect, while the restored image of this embodiment has better noise removal capability and visual effect. Fig.11 The degraded image and the image restored by the OBD algorithm are shown. The full width at half maximum of the stars in the restored image of this embodiment are 11.84, 8.02 and 6.16 pixels, respectively. The full width at half maximum of the restored image of this embodiment is close to the diffraction limit of the adaptive optical system.

[0123] In the real satellite imaging experiment, the results of this embodiment are compared with the results of the OBD algorithm. Fig.12 and Fig.13 shown. Fig.12 For the ALOS satellite, Fig.13 The satellites are COSMOS1461, and their orbits are 684×686km, 97.8° and 502×668km, 65.0° respectively. The restored image shows that the restored image of this embodiment has clearer edges and better noise removal effect.

[0124] It should also be noted that the various specific technical features described in the above specific embodiments can be combined in any suitable manner without contradiction. In order to avoid unnecessary repetition, the present invention will not further describe various possible combinations.

Claims

1. A multi-frame blind restoration method for an adaptive optical system corrected image, characterized in that: The specific steps include: Step 1: Use a gradient filter and Extract high-frequency information from the f-th frame degraded image g to obtain the extracted high-frequency information T is the matrix transpose, and the degraded image g is the degraded image corrected by the adaptive optical system; Step 2: Establish the following objective function to estimate the point spread function h corresponding to the f-th frame degraded image: Among them, x is the high-frequency information part of the original image, h i is the value of the i-th element of h, and h i ≥0,∑ i h i =1; λ and δ are regularization weight variables, * represents convolution operation; Step 3: Divide the objective function of step 2 into the first and second sub-problems based on the alternating iterative algorithm; The first sub-problem is: The second sub-problem is: Step 4: Use the iterative algorithm to perform the k-th iteration calculation on the first sub-problem to obtain the high-frequency restored image x at the k-th iteration calculation k , based on the solved x k The iterative algorithm is used to perform the k-th iteration calculation on the second subproblem, and the point spread function h at the k-th iteration calculation is obtained. k ,k=1,2,…K, where K is the total number of iterative calculations; Step 5: Based on the point spread function obtained by calculation, the multi-frame image restoration algorithm containing Poisson noise is used to restore the original image of the f-th frame, and it is determined whether f is less than NF. If so, f+1 is obtained, and step 1 is performed to restore the original image of the f+1th frame. Otherwise, the calculation is stopped. NF is the total number of image frames. In step 4, the first sub-problem is converted into the following form based on the Bregman iterative algorithm: Where b1 and b2 are Bregman variables, d1 = Hx, d2 = x, λ' = λ||x l ||2, where H is the matrix form of h after Fourier transformation, γ is the regularization weight variable; x l is the high-frequency restored image obtained in the last iterative calculation when the above formula 4 is iteratively calculated; When using the Bregman iterative algorithm to solve Formula 4, Formula 4 is divided into the following forms for solution: Where I is the identity matrix; When using the iterative algorithm to solve Formula 3, Formula 3 is divided into the following formulas for solution: Where N is the total number of elements in h; div(.) is the divergence calculation, is the gradient operator.

2. The method for blindly restoring a multi-frame image corrected by an adaptive optical system according to claim 1, characterized in that: The step 5 is specifically as follows: Step 5.1: Establish the following objective function: in, represents the estimated value of the original image of the fth frame, u f Represents the actual f-th frame of the original image, where f is the image frame number; Step 5.2: Use the following RL iterative algorithm to solve Equation 7: Among them, j represents the jth iteration calculation, Indicates that the original image is obtained by the jth iteration calculation; when the first iteration calculation is performed Where J represents the maximum number of iterations of the RL iterative algorithm; Step 5.3: Use the BM3D algorithm as a low-pass filter to perform low-pass filtering on the solved original image to obtain the final original image.

Images