A remote sensing image cloud removal method using group sparse prior tensor model

Abstract

The present invention discloses a remote sensing image declouding method using a group sparse prior tensor model, which belongs to the field of information processing technology and comprises the following steps: step 1, obtaining remote sensing images of continuous time nodes, stacking them in sequence along the spectral direction according to the time series, and forming a third-order original tensor with spatial and spectral smoothness; step 2, based on the low-rank tensor representation, characterizing the low-rank property of the clean component, decomposing the clean component into a spatial factor and a time-spectral factor according to the set rank through the low-rank tensor representation, and writing this representation relationship as a constraint condition; step 3, based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparse constraint term based on the spatial factor. The invention is used to realize efficient and accurate cloud removal of multi-temporal remote sensing images in complex scenes.

Description

Technical Field

[0001] The present invention relates to the field of information processing technology, and in particular to a remote sensing image cloud removal method using a group sparse prior tensor model. Background Art

[0002] Current satellite technology and sensor equipment have the ability to obtain high-spatial-resolution multi-temporal remote sensing images, which have extremely broad application prospects in many fields such as land resource inventory, environmental monitoring, and urban planning and design. However, these technologies also face some challenges. The light waves emitted by the sensor are difficult to penetrate thick clouds, and thick clouds are not only unevenly distributed, but also often cause all spectral images at the same time point to be contaminated. Therefore, in a multi-temporal remote sensing image sequence, images at certain specific time points may lose a lot of valuable information due to the interference of thick clouds and their shadows. This large amount of information loss not only weakens the availability of image data, but also has a serious hindering effect on subsequent application scenarios.

[0003] In multi-temporal remote sensing images, the distribution of clouds changes over time, and the same area is not contaminated by clouds at every shooting time point. In order to reconstruct the contaminated area, researchers usually use uncontaminated images of the same area at different time points as reference images, and reconstruct them using three types of information: time, space, and spectrum. Time-based methods take advantage of the consistency of the shooting areas of multi-temporal images at different time points and the similarity of detail information, but these methods find it difficult to accurately fit the uncontaminated images to the characteristics of the target image. Space-based methods are also a common strategy, but due to the widespread distribution of clouds at certain time points, the available spatial information is greatly reduced.

[0004] Therefore, the above method is only applicable to the case where the cloud-contaminated area is relatively small and has less texture details, otherwise the accuracy of reconstructing the contaminated area will be very limited. The method based on spectral information makes full use of the similarity and correlation between different spectra. However, in real scenes, thick clouds often contaminate all bands.

[0005] Methods based on spectral information cannot obtain sufficient reference information, so using only spectral information and spatial information is still not effective in removing clouds. In summary, how to remove clouds from multi-temporal remote sensing images in complex scenes has become an urgent and important technical challenge.

[0006] Therefore, a remote sensing image declouding method using a group sparse prior tensor model is proposed to solve or alleviate the above problems. Summary of the invention

[0007] The purpose of the present invention is to solve the shortcomings of the prior art and to propose a remote sensing image declouding method using a group sparse prior tensor model.

[0008] In order to achieve the above object, the present invention adopts the following technical solutions: A remote sensing image cloud removal method using a group sparse prior tensor model comprises the following steps: Step 1: obtain remote sensing images of consecutive time nodes, stack them in sequence along the spectral direction according to the time series, and form a third-order original tensor with spatial and spectral smoothness; Step 2: Based on the low-rank tensor representation, the low-rank property of the clean component is characterized. The clean component is decomposed into spatial factors and time-spectral factors according to the set rank through the low-rank tensor representation, and this representation relationship is written as a constraint condition; Step 3 is based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparsity constraint term based on the spatial factor; Step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3; In step 5, the target clean component X is obtained by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier optimization algorithm, and the target clean component X is decomposed according to the reverse step 1 to obtain the multi-temporal remote sensing image.

[0009] Preferably, step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, comprising the following steps: Set fidelity constraints and select appropriate effective coefficients for each regularization term; Combine the regularization terms in steps 2 and 3 and write it as a multi-objective convex optimization problem.

[0010] Preferably, the step 5 obtains the target clean component X by using a proximal alternating minimization algorithm embedded with an improved alternating direction multiplier optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain a multi-temporal remote sensing image, including the following steps: The PAM algorithm is used to update the target clean component X, spatial factor and time-spectral factor in sequence; In each iteration cycle of PAM calculation, auxiliary variables are introduced to calculate the spatial factor B using the improved alternating direction multiplier method; Calculate the relative errors of the updated target clean component X, spatial factor and time-spectral factor; When the relative error is less than a preset value, the iteration is stopped and the target clean component X is output.

[0011] Preferably, the step 1 acquires remote sensing images of consecutive time nodes, stacks them in sequence along the spectral direction in a time series, and forms a third-order original tensor with spatial and spectral smoothness, comprising the following steps: Multi-temporal remote sensing images acquired , stacked in time order along the spectral direction, forming a third-order tensor with spatial and spectral smoothness ,in and is the image size, is the number of spectra, is the number of time nodes; According to the information contained in multi-temporal remote sensing images, the original tensor y is a linear combination of the target clean component X, the cloud component tensor S and the model noise tensor N, and the formula is: .

[0012] Preferably, the step 2 describes the low rank of the clean component based on the low rank tensor representation, decomposes the clean component into spatial factors and time-spectral factors according to the set rank through the low rank tensor representation, and writes this representation relationship as a constraint condition, including the following steps: The target clean component X is calculated by the formula The tensor low-rank representation shown is decomposed into a spatial factor B and a time-spectral factor A according to the set rank r, where represents the tensor-matrix multiplication along the third direction, the spatial factor , time-spectrum factor , r is the set tensor rank. Preferably, step 3 is based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparse constraint term based on the spatial factor, including the following steps: According to the formula Design based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean components; among them, represents the tensor-matrix multiplication along the kth direction, Based on l 2,1 The group sparse regularization term of the norm, D k represents the difference matrix, whose element d in the i-th row and j-th column i,j Defined as: .

[0013] Preferably, step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, comprising the following steps: Based on step 2 and step 3, a weighted group sparse regularized tensor model for multi-temporal remote sensing image cloud removal is reconstructed. The weighted group sparse regularized tensor model is: , where Ω represents the index set of the cloud-free area. represents the projection operator that maps the pixel at Ω to itself and sets all other pixels to 0, is the Frobenius time smoothing constraint of the clean component, are the internal weight factors of the regularization term of the clean component, is the model balance parameter, is element-wise multiplication.

[0014] Preferably, the step 5 obtains the target clean component X by using a proximal alternating minimization algorithm embedded with an improved alternating direction multiplier optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain a multi-temporal remote sensing image, including the following steps: The step 5 is specifically as follows: The variables A and B are updated by the alternating proximal minimization method. The specific formula is: , where β>0 is the penalty parameter and h is the number of iterations; Step 5.1: Fix other variables in the hth iteration and update variable A in the h+1th iteration. The formula is ; By taking the derivative and expanding all tensors along the third dimension, the variables can be solved by the following formula , where B （3） and X （3） are the third dimension expansion matrices of tensors B and X respectively, is the transpose operator, I is the identity matrix; Step 5.2: Fix other variables in the hth iteration and update variable B in the h+1th iteration. The formula is ; Step 5.3: Introducing auxiliary variables and To calculate the variable B using the improved alternating direction multiplier method, the formula in step 5.2 is converted into an augmented Lagrangian function, which is ,in is the Lagrange multiplier, α>0 is the penalty coefficient, is the coefficient α associated with the penalty term in the augmented Lagrangian function; Step 5.4: By formula ,in is the number of iterations, the variable can be updated as follows ; Step 5.5: comprising step 5.5.1, step 5.5.2, step 5.5.3, step 5.5.4, step 5.5.5, and step 5.5.6; Step 5.5.1: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is: ; Step 5.5.2: Derivative both sides of the equation in step 5.5.1 to obtain the following variables The solution formula , The formula is ,in Represents the sum of all invariant items in this round of iterative calculation; Step 5.5.3: Expand all tensors in the formula of K along the spectral dimension to obtain the following formula ,in, is a set of difference operators, and is the characteristic parameter of the difference operator; Step 5.5.4: Apply two-dimensional fast Fourier transform (2-DFFT) and singular value decomposition to diagonalize the variables C and in the formula in step 5.5.3 into and , where F and U are two-dimensional DFT matrices, and is a diagonal matrix. The formula in step 5.5.3 can be converted to the Sylvester tensor equation for fast solution. The Sylvester tensor equation is ,in, Indicates point division, represents the folding operator, ones (n 1 ,n 2 ,r) is of size n 1 n 2 *r is a matrix of all 1s, ; Step 5.5.5: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is ,in, The tube of row i and column j All have the same solution form, which can be solved by the formula: Solve, where the variable , the i-th row and j-th column element of the weight coefficient matrix , is a very small constant value, the scaling function ; Step 5.5.6: When the convergence condition is reached or the maximum number of iterations is reached, stop the internal ADMM algorithm loop iteration, otherwise return to step 5.5.1. The convergence condition is as follows ; Step 5.6: Calculate the clean fraction as follows: ; Step 5.7: Optimize the clean component with mask according to the following formula ; Step 5.8: When the convergence condition is reached or the maximum number of iterations is reached, stop the external PAM algorithm iteration, otherwise return to step 5.1. The convergence condition is as follows ; Get the target clean component X, and decompose the third-order tensor into a multi-temporal remote sensing image by the reverse step 1, that is, get a clean multi-temporal remote sensing image, where is the auxiliary variable introduced, is the difference matrix, and for * and* The identity matrix, and is the image size, and are the differential operation matrices corresponding to the row and column directions of the image, C is the set of differential operators, is the weight coefficient matrix, and l is the Lagrange multiplier matrix.

[0015] The present invention has the following beneficial effects: The present invention proposes a multi-temporal remote sensing image declouding method based on a group sparse regularized tensor decomposition model. The method effectively mines the common features between different spectral bands by introducing a weighted group sparse function based on spatial factors and cloud components. Compared with the traditional regularization method based on the original image, this method can more accurately characterize the group sparsity and significantly reduce the computational complexity. In addition, by introducing a weight mechanism, the efficiency and effect of the model are further improved. In order to efficiently solve the model, the present invention develops an alternating proximal minimization algorithm embedded with an alternating direction multiplier method. This algorithm can converge quickly and obtain high-quality declouding effects through alternating optimization, and is suitable for complex scenes of multi-temporal remote sensing images. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings required for use in the embodiments are briefly introduced below. It should be understood that the following drawings only show certain embodiments of the present invention and therefore should not be regarded as limiting the scope. For ordinary technicians in this field, other related drawings can be obtained based on these drawings without creative work.

[0017] Figure 1 It is a flowchart of the present invention; Figure 2 It is a result schematic diagram of the present invention; Figure 3 This is a diagram of the cloud removal effect of the simulated data experiment in the present invention; Figure 4 This is a diagram of the cloud removal effect of the real data experiment in the present invention. DETAILED DESCRIPTION

[0018] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Generally, the components of the embodiments of the present invention described and shown in the drawings here can be arranged and designed in various different configurations.

[0019] Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention claimed for protection, but merely represents selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0020] It should be noted that similar reference numerals and letters denote similar items in the following drawings, and therefore, once an item is defined in one drawing, further definition and explanation thereof is not required in subsequent drawings.

[0021] In the description of the present invention, it should be understood that the terms "center", "up", "down", "left", "right", "vertical", "horizontal", "inside", "outside" and the like indicate orientations or positional relationships based on the orientations or positional relationships shown in the accompanying drawings, or are the orientations or positional relationships in which the product of the invention is conventionally placed when in use, or are the orientations or positional relationships conventionally understood by those skilled in the art. They are only for the convenience of describing the present invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore should not be understood as a limitation on the present invention.

[0022] Furthermore, the terms “first”, “second”, “third”, etc. are merely used for distinguishing descriptions and are not to be understood as indicating or implying relative importance.

[0023] In the description of the present invention, it is also necessary to explain that, unless otherwise clearly specified and limited, the terms "set", "install", "connect", and "connect" should be understood in a broad sense, for example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be a direct connection, or it can be indirectly connected through an intermediate medium, or it can be the internal communication of two elements. For ordinary technicians in this field, the specific meanings of the above terms in the present invention can be understood according to specific circumstances.

[0024] A remote sensing image cloud removal method using a group sparse prior tensor model comprises the following steps: Step 1 obtains remote sensing images of continuous time nodes, stacks them in sequence along the spectral direction according to the time series, and forms a third-order original tensor with spatial and spectral smoothness. Step 1 obtains remote sensing images of continuous time nodes, stacks them in sequence along the spectral direction according to the time series, and forms a third-order original tensor with spatial and spectral smoothness, including the following steps: Multi-temporal remote sensing images acquired , stacked in time order along the spectral direction, forming a third-order tensor with spatial and spectral smoothness ,in and is the image size, is the number of spectra, is the number of time nodes; According to the information contained in multi-temporal remote sensing images, the original tensor y is a linear combination of the target clean component X, the cloud component tensor S and the model noise tensor N, and the formula is: ; Step 2: Based on the low-rank tensor representation, the low-rank property of the clean component is characterized. The clean component is decomposed into spatial factors and time-spectral factors according to the set rank through the low-rank tensor representation, and this representation relationship is written as a constraint condition. The low-rank tensor representation can effectively characterize the low-rank property of the target tensor to achieve the desired cloud extraction effect. Step 2 describes the low rank of the clean component based on the low rank tensor representation, decomposes the clean component into spatial factors and time-spectral factors according to the set rank through the low rank tensor representation, and writes this representation relationship into constraints, including the following steps: The target clean component X is calculated by the formula The tensor low-rank representation shown is decomposed into a spatial factor B and a time-spectral factor A according to the set rank r, where represents the tensor-matrix multiplication along the third direction, the spatial factor , time-spectrum factor , r is the set tensor rank; Step 3: In order to restore the non-smooth characteristics of the cloud contaminated area, based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparsity constraint term based on the spatial factor; Step 3 is based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparse constraint term based on the spatial factor, including the following steps: According to the formula Design based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean components; among them, represents the tensor-matrix multiplication along the kth direction, Based on l 2,1 The group sparse regularization term of the norm, D k represents the difference matrix, whose element d in the i-th row and j-th column i,j Defined as: ; Step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3. Step 4 is to simultaneously characterize the low rank and group sparsity of the target clean component so that it can accurately restore the characteristics of the ground object. The obtained low-rank tensor filling model is a linear combination of two regular terms. Step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, including the following steps: Set fidelity constraints and select appropriate effective coefficients for each regularization term; Combine the regularization terms in steps 2 and 3 and write it as a multi-objective convex optimization problem; Step 4 builds a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, including the following steps: Based on step 2 and step 3, a weighted group sparse regularized tensor model for multi-temporal remote sensing image cloud removal is reconstructed. The weighted group sparse regularized tensor model is: , where Ω represents the index set of the cloud-free area. represents the projection operator that maps the pixel at Ω to itself and sets all other pixels to 0, is the Frobenius time smoothing constraint of the clean component, are the internal weight factors of the regularization term of the clean component, is the model balance parameter, is element-wise multiplication; Step 5: obtaining the target clean component X by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier optimization algorithm, and decomposing the target clean component X by the reverse step 1 to obtain the multi-temporal remote sensing image; Step 5 obtains the target clean component X by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain a multi-temporal remote sensing image, including the following steps: The PAM algorithm is used to update the clean component X, the spatial factor and the time-spectral factor in sequence; In each iteration cycle of PAM calculation, auxiliary variables are introduced to calculate the spatial factor B using the improved alternating direction multiplier method; Calculate the relative errors of the updated target clean component X, spatial factor and time-spectral factor; When the relative error is less than a preset value, the iteration is stopped and the target clean component X is output; Wherein, the step 5 obtains the target clean component X by embedding the proximal alternating minimization algorithm with the improved alternating direction multiplier method optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain the multi-temporal remote sensing image, which includes the following steps: The step 5 is specifically as follows: The variables A and B are updated by the alternating proximal minimization method. The specific formula is: , where β>0 is the penalty parameter and h is the number of iterations; Step 5.1: Fix other variables in the hth iteration and update variable A in the h+1th iteration. The formula is ; By taking the derivative and expanding all tensors along the third dimension, the variables can be solved by the following formula , where B （3） and X （3） are the third dimension expansion matrices of tensors B and X respectively, is the transpose operator, I is the identity matrix; Step 5.2: Fix other variables in the hth iteration and update variable B in the h+1th iteration. The formula is ; Step 5.3: Introducing auxiliary variables and To calculate the variable B using the improved alternating direction multiplier method, the formula in step 5.2 is converted into an augmented Lagrangian function, which is ,in is the Lagrange multiplier, α>0 is the penalty coefficient, is the coefficient α associated with the penalty term in the augmented Lagrangian function; Step 5.4: By formula ,in is the number of iterations, the variable can be updated as follows ; Step 5.5: comprising step 5.5.1, step 5.5.2, step 5.5.3, step 5.5.4, step 5.5.5, and step 5.5.6; Step 5.5.1: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is: ; Step 5.5.2: Derivative both sides of the equation in step 5.5.1 to obtain the following variables The solution formula , The formula is ,in Represents the sum of all invariant items in this round of iterative calculation; Step 5.5.3: Expand all tensors in the formula of K along the spectral dimension to obtain the following formula ,in, is a set of difference operators, and is the characteristic parameter of the difference operator; Step 5.5.4: Apply two-dimensional fast Fourier transform (2-DFFT) and singular value decomposition to diagonalize the variables C and in the formula in step 5.5.3 into and , where F and U are two-dimensional DFT matrices, and is a diagonal matrix. The formula in step 5.5.3 can be converted to the Sylvester tensor equation for fast solution. The Sylvester tensor equation is ,in, Indicates point division, represents the folding operator, ones (n 1 ,n 2 ,r) is of size n 1 n 2 *r is a matrix of all 1s, ; Step 5.5.5: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is ,in, The tube of row i and column j All have the same solution form, which can be solved by the formula: Solve, where the variable , the i-th row and j-th column element of the weight coefficient matrix , is a very small constant value, the scaling function ; Step 5.5.6: When the convergence condition is reached or the maximum number of iterations is reached, stop the internal ADMM algorithm loop iteration, otherwise return to step 5.5.1. The convergence condition is as follows ; Step 5.6: Calculate the clean fraction as follows: ; Step 5.7: Optimize the clean component with mask according to the following formula The purpose is to make the model refer to the spatial information of the cloud-free area in each iteration, while avoiding the degradation of the cloud-free area in the process of tensor low-rank representation; Step 5.8: When the convergence condition is reached or the maximum number of iterations is reached, stop the external PAM algorithm iteration, otherwise return to step 5.1. The convergence condition is as follows The purpose is to make the variable update step in the formula of step 5 stop at a reasonable number of iterations. It is necessary to constrain the minimum relative change size of each iteration. When the relative change size reaches the expected value, exit the loop; Get the target clean component X, and decompose the third-order tensor into a multi-temporal remote sensing image by the reverse step 1, that is, get a clean multi-temporal remote sensing image, where is the auxiliary variable introduced, is the difference matrix, and for * and * The identity matrix, and is the image size, and are the differential operation matrices corresponding to the row and column directions of the image, C is the set of differential operators, is the weight coefficient matrix, and l is the Lagrange multiplier matrix.

[0025] Figure 3 The following are the cloud removal effect diagrams of the present invention in the simulated data experiment, wherein the first row is the original image, the second row is the simulated cloud-contaminated image, and the third row is the cloud removal result. Figure 4The figures are the declouding effects of the present invention in the real data experiment, where the first row is the real cloud-contaminated image, and the second row is the declouding result. It can be seen that the method proposed in the present invention shows excellent performance in thick cloud removal and detailed texture restoration, and has excellent remote sensing image declouding capability.

[0026] In summary, the multi-temporal remote sensing image declouding method based on the group sparse regularized tensor decomposition model proposed in the present invention effectively mines the common features between different spectral bands by introducing the weighted group sparse function of the spatial factor and the cloud component, thereby more accurately reconstructing the image area contaminated by clouds. Compared with the traditional regularization method based on the original image, the model of the present invention can not only more accurately characterize the group sparsity, but also significantly reduce the computational complexity.

[0027] In order to further improve the efficiency and effect of the model, the present invention adopts a weight mechanism to better fit the characteristics of the target image by optimizing the weight distribution. In addition, in order to efficiently solve the model, the present invention develops an alternating proximal minimization algorithm embedded with the alternating direction multiplier method. The algorithm combines the efficient decomposition capability of ADMM and the optimization characteristics of alternating proximal minimization, and can quickly converge and solve complex optimization problems. Through this method, the present invention performs well in declouding multi-temporal remote sensing images, providing strong support for applications in related fields.

[0028] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. A remote sensing image cloud removal method using a group sparse prior tensor model, characterized in that: The steps include: Step 1: obtain remote sensing images of consecutive time nodes, stack them in sequence along the spectral direction according to the time series, and form a third-order original tensor with spatial and spectral smoothness; Step 2: Based on the low-rank tensor representation, the low-rank property of the clean component is characterized. The clean component is decomposed into spatial factors and time-spectral factors according to the set rank through the low-rank tensor representation, and this representation relationship is written as a constraint condition; Step 3 is based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparsity constraint term based on the spatial factor; Step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3; In step 5, the target clean component X is obtained by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier optimization algorithm, and the target clean component X is decomposed according to the reverse step 1 to obtain the multi-temporal remote sensing image.

2. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, and includes the following steps: Set fidelity constraints and select appropriate effective coefficients for each regularization term; Combine the regularization terms in steps 2 and 3 and write it as a multi-objective convex optimization problem.

3. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 5 obtains the target clean component X by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier method optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain a multi-temporal remote sensing image, including the following steps: The PAM algorithm is used to update the target clean component X, spatial factor and time-spectral factor in sequence; In each iteration cycle of PAM calculation, auxiliary variables are introduced to calculate the spatial factor B using the improved alternating direction multiplier method; Calculate the relative errors of the updated target clean component X, spatial factor and time-spectral factor; When the relative error is less than a preset value, the iteration is stopped and the target clean component X is output.

4. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 1 acquires remote sensing images of continuous time nodes, stacks them in sequence along the spectral direction in time series, and forms a third-order original tensor with spatial and spectral smoothness, including the following steps: Multi-temporal remote sensing images acquired , stacked in time order along the spectral direction, forming a third-order tensor with spatial and spectral smoothness ,in and is the image size, is the number of spectra, is the number of time nodes; According to the information contained in multi-temporal remote sensing images, the original tensor y is a linear combination of the target clean component X, the cloud component tensor S and the model noise tensor N, and the formula is: .

5. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 2 describes the low rank of the clean component based on the low rank tensor representation, decomposes the clean component into a spatial factor and a time-spectral factor according to the set rank through the low rank tensor representation, and writes this representation relationship as a constraint condition, including the following steps: The target clean component X is calculated by the formula The tensor low-rank representation shown is decomposed into a spatial factor B and a time-spectral factor A according to the set rank r, where represents the tensor-matrix multiplication along the third direction, the spatial factor , time-spectrum factor , r is the tensor rank set.

6. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 3 is based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean component. The spatial factors are subjected to differential operations along two spatial directions to obtain the differential image. 2,1 The norm constructs a group sparse constraint term based on the spatial factor, including the following steps: According to the formula Design based on l 2,1 The norm regularization characterizes the group sparsity of the spatial factors of the clean components; among them, represents the tensor-matrix multiplication along the kth direction, Based on l 2,1 The group sparse regularization term of the norm, D k represents the difference matrix, whose element d in the i-th row and j-th column i,j Defined as: 。 7. The remote sensing image declouding method using a group sparse prior tensor model according to claim 2, characterized in that: The step 4 constructs a low-rank tensor filling model for remote sensing cloud removal based on steps 2 and 3, and includes the following steps: Based on step 2 and step 3, a weighted group sparse regularized tensor model for multi-temporal remote sensing image cloud removal is reconstructed. The weighted group sparse regularized tensor model is: , where Ω represents the index set of the cloud-free area. represents the projection operator that maps the pixel at Ω to itself and sets all other pixels to 0, is the Frobenius time smoothing constraint of the clean component, are the internal weight factors of the regularization term of the clean component, is the model balance parameter, is element-wise multiplication.

8. The remote sensing image declouding method using a group sparse prior tensor model according to claim 1, characterized in that: The step 5 obtains the target clean component X by using the proximal alternating minimization algorithm embedded with the improved alternating direction multiplier method optimization algorithm, and decomposes the target clean component X according to the reverse step 1 to obtain a multi-temporal remote sensing image, including the following steps: The step 5 is specifically as follows: The variables A and B are updated by the alternating proximal minimization method. The specific formula is: , where β>0 is the penalty parameter and h is the number of iterations; Step 5.1: Fix other variables in the hth iteration and update variable A in the h+1th iteration. The formula is ; By taking the derivative and expanding all tensors along the third dimension, the variables can be solved by the following formula , where B （3） and X （3） are the third dimension expansion matrices of tensors B and X respectively, is the transpose operator, I is the identity matrix; Step 5.2: Fix other variables in the hth iteration and update variable B in the h+1th iteration. The formula is ; Step 5.3: Introducing auxiliary variables and To calculate the variable B using the improved alternating direction multiplier method, the formula in step 5.2 is converted into an augmented Lagrangian function, which is ,in is the Lagrange multiplier, α>0 is the penalty coefficient, is the coefficient α associated with the penalty term in the augmented Lagrangian function; Step 5.4: By formula , update the variable ,in, is the number of iterations; Step 5.5: comprising step 5.5.1, step 5.5.2, step 5.5.3, step 5.5.4, step 5.5.5, and step 5.5.6; Step 5.5.1: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is: ; Step 5.5.2: Derivative both sides of the equation in step 5.5.1 to obtain the following variables The solution formula , The formula is ,in Represents the sum of all invariant items in this round of iterative calculation; Step 5.5.3: Expand all tensors in the formula of K along the spectral dimension to obtain the following formula ,in, is a set of difference operators, and is the characteristic parameter of the difference operator; Step 5.5.4: Apply two-dimensional fast Fourier transform (2-DFFT) and singular value decomposition to diagonalize the variables C and in the formula in step 5.5.3 into and , where F and U are two-dimensional DFT matrices, and is a diagonal matrix. The formula in step 5.5.3 can be converted to the Sylvester tensor equation for fast solution. The Sylvester tensor equation is ,in, Indicates point division, represents the folding operator, ones(n1,n2,r) is a matrix of all 1s with a size of n1n2*r, ; Step 5.5.5: Fix the external h+1th and internal Other iteration variables, update the Variables for iterations , the formula is ,in, The tube of row i and column j All have the same solution form, which can be solved by the formula: Solve, where the variable , the i-th row and j-th column element of the weight coefficient matrix , is a very small constant value, the scaling function ; Step 5.5.6: When the convergence condition is reached or the maximum number of iterations is reached, stop the internal ADMM algorithm loop iteration, otherwise return to step 5.5.

1. The convergence condition is as follows ; Step 5.6: Calculate the clean fraction as follows: ; Step 5.7: Optimize the clean component with mask according to the following formula ; Step 5.8: When the convergence condition is reached or the maximum number of iterations is reached, stop the external PAM algorithm iteration, otherwise return to step 5.

1. The convergence condition is as follows ; Get the target clean component X, and decompose the third-order tensor into a multi-temporal remote sensing image by the reverse step 1, that is, get a clean multi-temporal remote sensing image, where is the auxiliary variable introduced, is the difference matrix, and for * and * The identity matrix, and is the image size, and are the differential operation matrices corresponding to the row and column directions of the image, C is the set of differential operators, is the weight coefficient matrix, and l is the Lagrange multiplier matrix.

Images