Projections onto Convex Sets Super-Resolution Reconstruction Based on Point Spread Function Estimation of Low-Resolution Remote Sensing Images

To solve the problem on inaccuracy when estimating the point spread function (PSF) of the ideal original image in traditional projection onto convex set (POCS) super-resolution (SR) reconstruction, this paper presents an improved POCS SR algorithm based on PSF estimation of low-resolution (LR) remote sensing images. The proposed algorithm can improve the spatial resolution of the image and benefit agricultural crop visual interpolation. The PSF of the high-resolution (HR) image is unknown in reality. Therefore, analysis of the relationship between the PSF of the HR image and the PSF of the LR image is important to estimate the PSF of the HR image by using multiple LR images. In this study, the linear relationship between the PSFs of the HR and LR images can be proven. In addition, the novel slant knife-edge method is employed, which can improve the accuracy of the PSF estimation of LR images. Finally, the proposed method is applied to reconstruct airborne digital sensor 40 (ADS40) three-line array images and the overlapped areas of two adjacent GF-2 images by embedding the estimated PSF of the HR image to the original POCS SR algorithm. Experimental results show that the proposed method yields higher quality of reconstructed images than that produced by the blind SR method and the bicubic interpolation method.


Introduction
Image resolution refers to the number of pixels contained in an image per unit area. This parameter is an important factor used to evaluate the quality of remote sensing images. The limitations of imaging systems and the external circumstances in obtaining images, including inherent sensor sampling frequency, defocusing, and atmospheric disturbances [1], result in low-quality images, which are blurred, misshapen, and exhibiting random noise. To solve these problems, two methods were proposed. First, the resolution can be enhanced by increasing the chip size. However, this approach is costly and cannot significantly improve the image resolution. Second, super-resolution (SR) reconstruction can be implemented using various algorithms. In the present study, time is used to compensate for space. Certain constraints or algorithms are employed to build a high-resolution (HR) image with the initial image estimation by using the wavelet bicubic interpolation for POCS reconstruction algorithm; the experimental results are evident. Liang [20] presented a POCS algorithm based on text features; in this method, text features are added as constraints to preserve the edge details and smoothen the noise in the text images. Meanwhile, the blind SR method is proposed by Sroubek and Flusser [21] can incorporates blur estimation into SR by performing an advanced deconvolution task. The model is built by the sharpness of edge regions and the smoothness of smooth regions in the total variance of the image, as well as the prior information of the image ambiguity function. Subsequently, the cross-iteration method is used to solve the model, and then the PSF and HR image are obtained. Thus, the HR image can be reconstructed even if the degradation model and the model's parameters of the camera sensors are unknown [22,23]. Given that the methods of estimating PSF are different, two algorithms are used for reconstruction to compare the results.
In this work, an improved POCS SR algorithm based on PSF estimation of LR remote sensing images is proposed. Exploration of the relationship between the PSF of the HR image and the PSF of multiple LR images is an essential part of this algorithm. In this study, the formula is deduced and approved by the experiments conducted. The conclusions are provided in the succeeding section. Moreover, the estimated PSF of the HR is embedded to the original POCS SR algorithm, and the reconstruction results of the three different SR methods (proposed method, blind SR method, and bicubic interpolation method), in the simulated experiment and the real experiment, are compared.

Observation Model
Generally, SR image reconstruction techniques present an inverse problem of recovering the HR image by degrading the LR images. The HR image is obtained under certain conditions, such as satisfying the theory of Nyquist sampling, and is affected by some of the inherent resolution limitations of sensors in the acquisition process [24], including warping, blurring, subsampling operators, and additive noise. Therefore, an observation model (1) can be formulated, which relates to the ideal HR image f to the corresponding i-th observed LR images g i . The model can overcome the inherent resolution limitation of the LR imaging systems. The LR images display different subpixel shifts from each other because the spatial resolution is very low to capture all the details of the original scene. Finally, a goal image with denser pixel values and rich image information, called HR image, will be achieved: where f is the ideal undegraded image required to be calculated, g i is the observed LR image, and H i represents the blur matrix, including relative camera-scene motion blur, sensor blur, atmospheric turbulence, and optical blur. Generally, the blur matrix is modeled as convolution with an unknown PSF to estimate blurs [25]. B i is the warp matrix (e.g., rotation, translation scaling, and so on).
The relative movement parameters can be estimated using the subpixel shifts of the multiple LR images. D i represents a subsampling matrix, and n is the lexicographically ordered noise vector. The observation model is illustrated in Figure 1. initial image estimation by using the wavelet bicubic interpolation for POCS reconstruction algorithm; the experimental results are evident. Liang [20] presented a POCS algorithm based on text features; in this method, text features are added as constraints to preserve the edge details and smoothen the noise in the text images. Meanwhile, the blind SR method is proposed by Sroubek and Flusser [21] can incorporates blur estimation into SR by performing an advanced deconvolution task. The model is built by the sharpness of edge regions and the smoothness of smooth regions in the total variance of the image, as well as the prior information of the image ambiguity function. Subsequently, the cross-iteration method is used to solve the model, and then the PSF and HR image are obtained. Thus, the HR image can be reconstructed even if the degradation model and the model's parameters of the camera sensors are unknown [22,23]. Given that the methods of estimating PSF are different, two algorithms are used for reconstruction to compare the results. In this work, an improved POCS SR algorithm based on PSF estimation of LR remote sensing images is proposed. Exploration of the relationship between the PSF of the HR image and the PSF of multiple LR images is an essential part of this algorithm. In this study, the formula is deduced and approved by the experiments conducted. The conclusions are provided in the succeeding section. Moreover, the estimated PSF of the HR is embedded to the original POCS SR algorithm, and the reconstruction results of the three different SR methods (proposed method, blind SR method, and bicubic interpolation method), in the simulated experiment and the real experiment, are compared.

Observation Model
Generally, SR image reconstruction techniques present an inverse problem of recovering the HR image by degrading the LR images. The HR image is obtained under certain conditions, such as satisfying the theory of Nyquist sampling, and is affected by some of the inherent resolution limitations of sensors in the acquisition process [24], including warping, blurring, subsampling operators, and additive noise. Therefore, an observation model (1) can be formulated, which relates to the ideal HR image to the corresponding -th observed LR images . The model can overcome the inherent resolution limitation of the LR imaging systems. The LR images display different subpixel shifts from each other because the spatial resolution is very low to capture all the details of the original scene. Finally, a goal image with denser pixel values and rich image information, called HR image, will be achieved: where is the ideal undegraded image required to be calculated, is the observed LR image, and represents the blur matrix, including relative camera-scene motion blur, sensor blur, atmospheric turbulence, and optical blur. Generally, the blur matrix is modeled as convolution with an unknown PSF to estimate blurs [25].
is the warp matrix (e.g., rotation, translation scaling, and so on). The relative movement parameters can be estimated using the subpixel shifts of the multiple LR images.
represents a subsampling matrix, and is the lexicographically ordered noise vector. The observation model is illustrated in Figure 1.  Image degradation occurs when the acquired image is corrupted by many factors. The image egradation process can be viewed as a linear invariant system, in which the noises can be ignored. The degradation model is described by Equation (2): The model is composed of four main attributes: the original image without degradation γ, the degraded image ϕ, a PSF H, and some noises η. * is the convolution operating symbol. To restore the quality of the image, H can be estimated using some PSF estimation methods, such as the knife-edge method. When the original image is processed by downsampling, the PSF h of the downsampled image is obtained. However, h does not apply to the model (2); hence, the relationship between H and h is derived in Section 2.3. The derived formula will be applied to the SR reconstruction model.

Principle of POCS SR Algorithm
In this work, the POCS SR reconstruction algorithm is used to obtain high-quality remote sensing data, which can meet the requirements of agricultural data sources. The algorithm, which is simple and effective, is a collection theory of the image reconstruction method. Given the flexible space-domain observation model and the powerful prior knowledge embedding capability, the owned feasible region of the reconstructed image consists of an intersection consistency projective convexity set and a convex constraint set. The POCS algorithm [26] is an iterative operation; the operator of the corresponding convex constraint set projects the points in the solution space to the nearest point on the surface of the convex set. After a finite number of iterations, a solution to the intersection set that converges to the convex constraint set is finally found. The POCS SR algorithm is detailed as follows: Step 1: Estimate the image f 0 by using the linear interpolation method for LR images.
Step 2: Compute the motion compensation of the pixel of each LR image. The correspondence between the LR image and the HR image is given by Equation (3): g(m 1 , m 2 , l) = ∑ n 1 ,n 2 f (n 1 , n 2 ) h n 1 , n 2 ; m 1 , m 2 , l + n(m 1 , m 2 , l) where (m 1 , m 2 ) is the point in the LR image, and (n 1 , n 2 ) is the corresponding point in the HR image.

1.
Obtain the position of the pixel on the LR image of each frame g(m 1 , m 2 , l) and on the HR image f (n 1 , n 2 ).

2.
Calculate the parameter, h n 1 , n 2 ; m 1 , m 2 , l , which represents the range and the value of PSF according to the position of the pixel.

4.
Calculate the residuals r f (n 1 , n 2 , k) between the real image and the simulated image. The formula can be described by (6).
where h(n 1 , n 2 ; m 1 , m 2 , l) is the impulse response coefficient, δ 0 is the confidence level on the observed result. In this paper, δ 0 = c δ v , where the point δ v is the standard deviation of the noise, and c ≥ 0 is determined by an appropriate statistical confidence range. These settings define HR images that are consistent with the observed LR image frames within a certain confidence range proportional to the observed noise variation.

5.
Correct the pixel value of the HR image according to the residuals.
Step 3: Repeat from Step 2 until convergence Given a projection operator, the estimated valuef (m 1 , m 2 , l) of the HR image f (m 1 , m 2 , l) can be obtained from all the LR images g(n 1 , n 2 , k) through many iterations, such as Equation (7): where T is the combination of all the relaxation projection operators associated with C(n 1 , n 2 , k).
The initial estimate, f 0 (m 1 , m 2 , l), is obtained by bilinear interpolation of the reference frame in the super-resolution grid.

Relationship between H and h
In this section, the relationship between H and h is deduced and validated, and the simulation experiment is designed to verify the correctness of the formula. The derived formula is proven to be suitable for SR reconstruction. In the process of SR reconstruction, the PSF of the HR image can be estimated by the PSF of the LR image; and PSF high = k·PSF low , where the downsampling ratio is k.
When knife-edge areas are extracted from a remote sensing image with gray values from 0 to 1, the original signal along the gradient direction is represented by the unit step signal E. Given the PSF H and degrading and downsampling operator D, ε is finally expressed as the signal in the image. Equation (8) can be determined in the first downsampling model: where * is the downsampling operating symbol, the downsampling operator D k can be taken as a calculation of one-dimensional downsampling multiples of k, which is the equivalent of k compression from this function. The formula is shown in (9): Therefore, when the variable is less than 0, the signal value of E is 0; when the variable is greater than 0, the signal value is 1. With general downsampling using the convolution operation method, the knife-edge areas can be mathematically formulated as follows: According to the first downsampling model, Equation (11) can be deduced: Based on the principle of knife-edge method, the PSF can be obtained from the derivative function of the edge spread function (ESF). The PSF is calculated as shown in Equation (12) by the knife-edge method: We can assume that the original signal can be restored effectively by the PSF with the values calculated by the knife-edge method. Deconvolution is used in the downsampling image g and the PSF initially. The deconvolution image requires rise sampling, so that the original image F can be obtained, as shown in Equation (13): where ⊗ is the deconvolution operation. Upsampling operator U 1 k can be taken as the calculation of one-dimensional upsampling multiples of 1 k 2 , which is the equivalent of the 1/k compression from the two-dimensional function of two coordinate axes, as shown in (14): If Equation (14) is correct, Equation (15) must exist: Moreover: Hypothesis (13) can be proven as tenable because Equations (11) and (16) yield the same results; hence, Equation (15) is correct.
The original signal can be restored effectively by the PSF of the existing degraded image proof. The signal can be applied in the process of image SR reconstruction. The PSF of the LR image calculated by the knife-edge method can be applied in the SR reconstruction process.
The PSF is approximated by Gaussian functions with appropriate parameters because the PSF follows a Gaussian distribution [27][28][29]. The PSF H i of the HR image can be written as follows: The PSF h i of the LR image can be expressed using Equation (18): Therefore, the relationship between the Gaussian function parameters and the PSF of the H i and h i images can be deduced from Equation (12), as shown in Equations (20), where Σ is the Gaussian function parameter of H i , and σ is the Gaussian function parameter of h i : The simulation experiment is conducted to verify the proposed computation formula. A man-made knife-edge figure is drawn using computer language, and the figure is degraded by convolution by using the PSF, which is estimated by Gaussian functions. The selection process of the knife-edge area is shown in Figure 2. The seven Gaussian function parameters were set at 0.5, 0.75, 1.0, 1.5, 1.75, 2.0, and 2.5. The resized image by resampling is presented with different multiples of k set at 0.5, 1.5, 2.0, 2.5, and 3.0. The knife-edge area must be selected to obtain the Gaussian function parameter of the sampled image; the size of the region is 15 × 15 pixels. According to the formula derived, the result is in correspondence with the original parameter. The results are summarized in Table 1, and Figure 3 shows that the Gaussian function parameters between the original and scaled image basically meet the linear relationship. Therefore, the relationship between the Gaussian function parameters and the PSF of the and ℎ images can be deduced from Equation (12), as shown in Equations (20), where is the Gaussian function parameter of , and is the Gaussian function parameter of ℎ : 1 The simulation experiment is conducted to verify the proposed computation formula. A manmade knife-edge figure is drawn using computer language, and the figure is degraded by convolution by using the PSF, which is estimated by Gaussian functions. The selection process of the knife-edge area is shown in Figure 2. The seven Gaussian function parameters were set at 0.5, 0.75, 1.0, 1.5, 1.75, 2.0, and 2.5. The resized image by resampling is presented with different multiples of set at 0.5, 1.5, 2.0, 2.5, and 3.0. The knife-edge area must be selected to obtain the Gaussian function parameter of the sampled image; the size of the region is 15 × 15 pixels. According to the formula derived, the result is in correspondence with the original parameter. The results are summarized in Table 1, and Figure 3 shows that the Gaussian function parameters between the original and scaled image basically meet the linear relationship. The results in Figure 3 and Table 1 show the PSF Gaussian function parameters between the images before the downsampling and after satisfying the linear relationship when the image is scaled at different scales of . The ratio of the scaling parameter variation to the original parameter variation is similar to . The ratio satisfies the formula which was just deduced in Equation (20).
The size of the region is 15 × 15 The results in Figure 3 and Table 1 show the PSF Gaussian function parameters between the images before the downsampling and after satisfying the linear relationship when the image is scaled at different scales of k. The ratio of the scaling parameter variation to the original parameter variation is similar to k. The ratio satisfies the formula which was just deduced in Equation (20).  The correctness of Equation (20) in the real image can be proven by the following experiments. The experimental images with some knife-edge areas can be selected by the ADS 40 remote sensing image with the size of 200 × 200 (Example 1) and the unmanned aerial vehicle (UAV) image with the size of 800 × 800 (Example 2). The experimental data are shown in Figure 4.    The correctness of Equation (20) in the real image can be proven by the following experiments. The experimental images with some knife-edge areas can be selected by the ADS 40 remote sensing image with the size of 200 × 200 (Example 1) and the unmanned aerial vehicle (UAV) image with the size of 800 × 800 (Example 2). The experimental data are shown in Figure 4.   The correctness of Equation (20)      First, four LR images must be acquired from the downsampling model, as shown in Figure 5, in which the original image becomes a series of LR image sequences with size of 100 × 100. Second, the knife-edge areas with four LR images and the experimental image are estimated using the PSF estimation based on slant knife-edge method. The PSFs can be obtained separately.
Finally, an oversampling rate is set at 2, indicating that the original image is zoomed out in half in the experiment. Based on the relationship between the images before the downsampling and after deducing, the Gaussian function parameter σ, after downsampling, must be half of the original. The results are shown in Table 2. The value coincides with the equation, proving that SR image reconstruction based on the PSFs of LR images is possible.  First, four LR images must be acquired from the downsampling model, as shown in Figure 5, in which the original image becomes a series of LR image sequences with size of 100 × 100.  First, four LR images must be acquired from the downsampling model, as shown in Figure 5, in which the original image becomes a series of LR image sequences with size of 100 × 100. Second, the knife-edge areas with four LR images and the experimental image are estimated using the PSF estimation based on slant knife-edge method. The PSFs can be obtained separately.
Finally, an oversampling rate is set at 2, indicating that the original image is zoomed out in half in the experiment. Based on the relationship between the images before the downsampling and after deducing, the Gaussian function parameter σ, after downsampling, must be half of the original. The results are shown in Table 2. The value coincides with the equation, proving that SR image reconstruction based on the PSFs of LR images is possible.  Second, the knife-edge areas with four LR images and the experimental image are estimated using the PSF estimation based on slant knife-edge method. The PSFs can be obtained separately.
Finally, an oversampling rate is set at 2, indicating that the original image is zoomed out in half in the experiment. Based on the relationship between the images before the downsampling and after deducing, the Gaussian function parameter σ, after downsampling, must be half of the original. The results are shown in Table 2. The value coincides with the equation, proving that SR image reconstruction based on the PSFs of LR images is possible.

PSF Estimation of Low-Resolution Remote Sensing Images
The optical information of the remote sensing image is blurred in the process of the image capture because of the relative motion between the object being photographed and the satellite and the CCD or the atmosphere turbulence. The PSF of an imaging platform can represent the response of an imaging system to a point source. Therefore, the calculation of the PSF in the acquired image is a significant step to restore the ideal remote sensing images. The blurring process [30], as a convolution of an image, is shown in Equation (21): where * is the convolution operator, h(x, y) represents the PSF, f (x, y) is the original image, g(x, y) is the degradation of the image, and α and β are the blurring filters. However, even when the PSF [31] is measurable, it is influenced by some unpredictable conditions; hence, many methods are proposed to solve the problem. The most common methods in PSF estimation are the knife-edge [32], spotlight [33], and pulse [34] methods. In addition, the knife-edge method is the most commonly employed method. The principles of the typical knife-edge method [35,36] and the ideal knife-edge area are shown in Figure 6. The knife-edge area is assumed as a square area, and the knife edge goes through the center of the area. Each row of the knife-edge satisfies formula (22), in which the value of the image greater than the edge boundary line x 0 is 1, and the value less than x 0 is zero [36,37]. However, the typical knife-edge method for point spread function estimation is limited by edge slant angle. The knife-edge in the image must be parallel to the sampling direction and the slant angle should be within 8 • [32]. Under most circumstances, slant knife-edges have certain slope to the direction of the ideal ones; thus, the ideal knife-edge cannot be determined in all situations. Qin et al. [38] estimated the PSF by a robust method to solve this problem in the typical knife-edge method; they built a mathematical model of the relationship between the line spread function (LSF) estimated by the typical knife-edge method and the real PSF. Although an accurate PSF is obtained, the computed PSF still contains an error because the algorithm uses the discrete function to derive the subpixel error:

PSF Estimation of Low-Resolution Remote Sensing Images
The optical information of the remote sensing image is blurred in the process of the image capture because of the relative motion between the object being photographed and the satellite and the CCD or the atmosphere turbulence. The PSF of an imaging platform can represent the response of an imaging system to a point source. Therefore, the calculation of the PSF in the acquired image is a significant step to restore the ideal remote sensing images. The blurring process [30], as a convolution of an image, is shown in Equation (21): where * is the convolution operator, ℎ( , ) represents the PSF, ( , ) is the original image, ( , ) is the degradation of the image, and are the blurring filters. However, even when the PSF [31] is measurable, it is influenced by some unpredictable conditions; hence, many methods are proposed to solve the problem. The most common methods in PSF estimation are the knife-edge [32], spotlight [33], and pulse [34] methods. In addition, the knifeedge method is the most commonly employed method. The principles of the typical knife-edge method [35,36] and the ideal knife-edge area are shown in Figure 6. The knife-edge area is assumed as a square area, and the knife edge goes through the center of the area. Each row of the knife-edge satisfies formula (22), in which the value of the image greater than the edge boundary line is 1, and the value less than is zero [36,37]. However, the typical knife-edge method for point spread function estimation is limited by edge slant angle. The knife-edge in the image must be parallel to the sampling direction and the slant angle should be within 8 [32]. Under most circumstances, slant knife-edges have certain slope to the direction of the ideal ones; thus, the ideal knife-edge cannot be determined in all situations. Qin et al. [38] estimated the PSF by a robust method to solve this problem in the typical knife-edge method; they built a mathematical model of the relationship between the line spread function (LSF) estimated by the typical knife-edge method and the real PSF. Although an accurate PSF is obtained, the computed PSF still contains an error because the algorithm uses the discrete function to derive the subpixel error: In this study, a novel slant knife-edge method is used because this method fits the LSF directly and ensures the evenness of the edge spread function (ESF) sample to improve the accuracy of the PSF estimation. Figure 7 illustrates the process of the novel slant knife-edge method in a simplified sequence flow diagram. However, the area can be searched using several methods, and the three following requirements must be satisfied:


The area cannot be extracted from the borders of the image to avoid the noise around the borders. In this study, a novel slant knife-edge method is used because this method fits the LSF directly and ensures the evenness of the edge spread function (ESF) sample to improve the accuracy of the PSF estimation. Figure 7 illustrates the process of the novel slant knife-edge method in a simplified sequence flow diagram. However, the area can be searched using several methods, and the three following requirements must be satisfied:

•
The area cannot be extracted from the borders of the image to avoid the noise around the borders.

•
The area must be excellent in linearity to ensure the accuracy of the PSF estimation. • Evident gray value differences between two sides of the edge to reduce the influence of the noise must be obtained. The area must be excellent in linearity to ensure the accuracy of the PSF estimation.  Evident gray value differences between two sides of the edge to reduce the influence of the noise must be obtained. The original data and results of the experiment are displayed in Figure 8. The measurement area and three examples of extracted knife-edge areas are also shown. Most of the knife-edge areas extracted successfully and the process of removing the weak-related areas are displayed. Follow ups are described in detail. ESF sampling, ESF denoising, ESF resampling, and LSF sampling results are shown in Figure 9.  The original data and results of the experiment are displayed in Figure 8. The measurement area and three examples of extracted knife-edge areas are also shown. Most of the knife-edge areas extracted successfully and the process of removing the weak-related areas are displayed. Follow ups are described in detail. ESF sampling, ESF denoising, ESF resampling, and LSF sampling results are shown in Figure 9. The area must be excellent in linearity to ensure the accuracy of the PSF estimation.  Evident gray value differences between two sides of the edge to reduce the influence of the noise must be obtained. The original data and results of the experiment are displayed in Figure 8. The measurement area and three examples of extracted knife-edge areas are also shown. Most of the knife-edge areas extracted successfully and the process of removing the weak-related areas are displayed. Follow ups are described in detail. ESF sampling, ESF denoising, ESF resampling, and LSF sampling results are shown in Figure 9.

Examples of Simulated Images
The experimental results did not demonstrate the effectiveness of the SR reconstruction algorithms qualitatively because of the lack of the ideal HR image. Therefore, the design of a simulation experiment is necessary by applying three different SR approaches and comparing the difference between these reconstructed images and the original image by the peak value signal-tonoise ratio (PSNR) and the mean square error (MSE) [39].
In this section, two series of comparative experiments are designed to evaluate the correctness of the deduced formula presenting the relation of the PSF Gaussian function parameter before and after SR reconstruction, as well as to compare the differences among the three reconstruction methods, namely, the proposed method, blind SR method, and bicubic interpolation method, on the original HR image. Moreover, given PSNR and MSE as the performance evaluation indicators, the experiment verifies the efficiency of the SR image reconstructed by the modified algorithm.
The experimental process is designed as follows: 1. The experimental image with some knife-edge areas can be selected by the ADS 40 remote sensing image with the size of 314 × 314, as shown in Figure 10a. 2. Referencing to the SR observation model, the original blurred image with given low pass and downsampled with factor 2 generated four LR images with the size of 157 × 157, as shown in Figure 10b. One of the LR simulated images is shown in Figure 10c. These four images correspond to the actual transformation parameters of the reference image, as shown in Table 3, where Dx represents the actual offset in the horizontal direction of the simulated LR image, and Dy is the actual offset in the vertical direction. D is the rotation angle because this experiment mainly considered translation; hence, the rotation angle is 0. The default units are pixels and degrees.

Examples of Simulated Images
The experimental results did not demonstrate the effectiveness of the SR reconstruction algorithms qualitatively because of the lack of the ideal HR image. Therefore, the design of a simulation experiment is necessary by applying three different SR approaches and comparing the difference between these reconstructed images and the original image by the peak value signal-to-noise ratio (PSNR) and the mean square error (MSE) [39].
In this section, two series of comparative experiments are designed to evaluate the correctness of the deduced formula presenting the relation of the PSF Gaussian function parameter before and after SR reconstruction, as well as to compare the differences among the three reconstruction methods, namely, the proposed method, blind SR method, and bicubic interpolation method, on the original HR image. Moreover, given PSNR and MSE as the performance evaluation indicators, the experiment verifies the efficiency of the SR image reconstructed by the modified algorithm.
The experimental process is designed as follows: 1.
The experimental image with some knife-edge areas can be selected by the ADS 40 remote sensing image with the size of 314 × 314, as shown in Figure 10a.

2.
to the SR observation model, the original blurred image with given low pass and downsampled with factor 2 generated four LR images with the size of 157 × 157, as shown in Figure 10b. One of the LR simulated images is shown in Figure 10c. These four images correspond to the actual transformation parameters of the reference image, as shown in Table 3, where Dx represents the actual offset in the horizontal direction of the simulated LR image, and Dy is the actual offset in the vertical direction. Dθ is the rotation angle because this experiment mainly considered translation; hence, the rotation angle is 0. The default units are pixels and degrees.

3.
The knife-edge areas with four LR images are estimated using the PSF estimation based on slant knife-edge method, which is an accurate method. PSFs can then be obtained separately.

4.
According to the derived formula (20), the PSF must be multiplied by the downsampling factor 2. Afterward, the POCS method with the estimated 2*PSF is used to reconstruct the HR image from these four downsampling LR images.

5.
The blind SR and bicubic interpolation methods are used for the comparative experiments.
The similarity between the experimental image and the resulting images of the three methods are observed. 6.
Evaluation of experimental results.   (20), the PSF must be multiplied by the downsampling factor 2. Afterward, the POCS method with the estimated 2*PSF is used to reconstruct the HR image from these four downsampling LR images. 5. The blind SR and bicubic interpolation methods are used for the comparative experiments. The similarity between the experimental image and the resulting images of the three methods are observed. 6. Evaluation of experimental results. No significant detailed difference between the original image and the modified reconstructed image in Figure 11 can be observed. The results testify that the algorithm can effectively reconstruct images. However, the details of the resulting images of the blind SR reconstruction algorithm and the bicubic interpolation algorithm are inconsistent with the original image in terms of heavy noise, aliasing, excessive sharpening, and so on. Nevertheless, both algorithms can sharpen edges and receive large amounts of information.  No significant detailed difference between the original image and the modified reconstructed image in Figure 11 can be observed. The results testify that the algorithm can effectively reconstruct images. However, the details of the resulting images of the blind SR reconstruction algorithm and the bicubic interpolation algorithm are inconsistent with the original image in terms of heavy noise, aliasing, excessive sharpening, and so on. Nevertheless, both algorithms can sharpen edges and receive large amounts of information.  3. The knife-edge areas with four LR images are estimated using the PSF estimation based on slant knife-edge method, which is an accurate method. PSFs can then be obtained separately. 4. According to the derived formula (20), the PSF must be multiplied by the downsampling factor 2. Afterward, the POCS method with the estimated 2*PSF is used to reconstruct the HR image from these four downsampling LR images. 5. The blind SR and bicubic interpolation methods are used for the comparative experiments. The similarity between the experimental image and the resulting images of the three methods are observed. 6. Evaluation of experimental results. No significant detailed difference between the original image and the modified reconstructed image in Figure 11 can be observed. The results testify that the algorithm can effectively reconstruct images. However, the details of the resulting images of the blind SR reconstruction algorithm and the bicubic interpolation algorithm are inconsistent with the original image in terms of heavy noise, aliasing, excessive sharpening, and so on. Nevertheless, both algorithms can sharpen edges and receive large amounts of information.   Table 4 demonstrates that the proposed method contains less signal distortion and substantially surpasses the blind SR reconstruction and the bicubic interpolation algorithms in terms of the highest PSNR and the minimum MSE. Moreover, the results prove that the POCS method with the estimated PSF can improve the recovery image information and achieve a good performance of SR reconstruction.

Examples of Real Images
In this section, the proposed algorithm is applied in practice. The experiments are designed to investigate the effect of SR reconstruction of the POCS method with the estimated PSF in the agricultural application. Given that data sources in agriculture use HR remote sensing images, the agricultural region of the GF-2 image is selected as the experimental data with the image size of 400 × 400. The set of LR inputs are the overlapped areas of two adjacent images. The other data are the UAV images with size of 500 × 500 to prove the stability of the proposed algorithm. The experimental results are compared with the blind SR reconstruction and the bicubic interpolation reconstruction to verify the effectiveness. Because the reference image is unavailable, the quality metrics PSNR or MSE cannot be used to compare the advantages of the three algorithms. Therefore, we choose a noreference metric , which can react in a natural way to the presence of noise and blur, to provide a quantitative measure of true image. And its value drops means the variance of noise rises, and the image becomes blurry [40]. Simultaneously, manual visual interpretation also is a criterion method.
In the first set of experiments, the knife-edge areas are extracted using the novel slant knife-edge method. The experimental data and details are shown in Figure 12. However, the edge of the agriculture land is unclear, and the house is fuzzy, thereby preventing the researchers from surveying the area accurately.  Table 4 demonstrates that the proposed method contains less signal distortion and substantially surpasses the blind SR reconstruction and the bicubic interpolation algorithms in terms of the highest PSNR and the minimum MSE. Moreover, the results prove that the POCS method with the estimated PSF can improve the recovery image information and achieve a good performance of SR reconstruction.

Examples of Real Images
In this section, the proposed algorithm is applied in practice. The experiments are designed to investigate the effect of SR reconstruction of the POCS method with the estimated PSF in the agricultural application. Given that data sources in agriculture use HR remote sensing images, the agricultural region of the GF-2 image is selected as the experimental data with the image size of 400 × 400. The set of LR inputs are the overlapped areas of two adjacent images. The other data are the UAV images with size of 500 × 500 to prove the stability of the proposed algorithm. The experimental results are compared with the blind SR reconstruction and the bicubic interpolation reconstruction to verify the effectiveness. Because the reference image is unavailable, the quality metrics PSNR or MSE cannot be used to compare the advantages of the three algorithms. Therefore, we choose a no-reference metric Q, which can react in a natural way to the presence of noise and blur, to provide a quantitative measure of true image. And its value drops means the variance of noise rises, and the image becomes blurry [40]. Simultaneously, manual visual interpretation also is a criterion method.
In the first set of experiments, the knife-edge areas are extracted using the novel slant knife-edge method. The experimental data and details are shown in Figure 12. However, the edge of the agriculture land is unclear, and the house is fuzzy, thereby preventing the researchers from surveying the area accurately. The knife-edge areas are shown in Figure 13. The knife-edge region selection is based on the relevant center coordinate of the knife-edge area. Two knife-edge areas are selected (c) from a total of six knife-edge areas (b), as shown in Table 5. Experimental data show that the center coordinates [36,320] of the region is the most clear knife-edge area. Therefore, with the point as the center, a square is drawn with a radius of 7 pixels to select the area, which is shown in (d). Given these data, the PSF is calculated for subsequent reconstruction experiments.   The knife-edge areas are shown in Figure 13. The knife-edge region selection is based on the relevant center coordinate of the knife-edge area. Two knife-edge areas are selected (c) from a total of six knife-edge areas (b), as shown in Table 5. Experimental data show that the center coordinates [36,320] of the region is the most clear knife-edge area. Therefore, with the point as the center, a square is drawn with a radius of 7 pixels to select the area, which is shown in (d). Given these data, the PSF is calculated for subsequent reconstruction experiments. The knife-edge areas are shown in Figure 13. The knife-edge region selection is based on the relevant center coordinate of the knife-edge area. Two knife-edge areas are selected (c) from a total of six knife-edge areas (b), as shown in Table 5. Experimental data show that the center coordinates [36,320] of the region is the most clear knife-edge area. Therefore, with the point as the center, a square is drawn with a radius of 7 pixels to select the area, which is shown in (d). Given these data, the PSF is calculated for subsequent reconstruction experiments.    Subsequently, the two selected PSFs are taken into the POCS algorithm to obtain the reconstructed image. The nominal values of metric Q of three reconstructed images as shown in Table 6. Metric Q of the image based on the proposed algorithm are the maximum number, it shows that the algorithm has a good visual performance and detail preservation. What's more, Figure 14 shows the results and details of the HR image reconstruction. In region A, the edges of the building are clear and distinct in (e) and (f). The bicubic interpolation method can only interpolate, but the interpolation is unclear. Moreover, the images reconstructed by the blind SR reconstruction method present jagged edges. Therefore, the proposed method and the blind SR reconstruction method can effectively reconstruct the details of the building. In summary, the results indicate that all these algorithms can effectively produce robust SR images. However, the proposed method demonstrated better effect than the other algorithms.  Table 6. Metric of the image based on the proposed algorithm are the maximum number, it shows that the algorithm has a good visual performance and detail preservation. What's more, Figure 14 shows the results and details of the HR image reconstruction. In region A, the edges of the building are clear and distinct in (e) and (f). The bicubic interpolation method can only interpolate, but the interpolation is unclear. Moreover, the images reconstructed by the blind SR reconstruction method present jagged edges. Therefore, the proposed method and the blind SR reconstruction method can effectively reconstruct the details of the building. In summary, the results indicate that all these algorithms can effectively produce robust SR images. However, the proposed method demonstrated better effect than the other algorithms.  In the second set of experiments, we choose the most clear knife-edge area [114,329] in real images according to the previous experimental procedure. The UAV images, as well as the results and details of the HR image reconstruction, are shown in Figure 15 and the values of metric are shown in Table 7. The reconstruction results of the proposed method are the most natural; the image becomes clearer with the increase in the amount of image information. Although the blind SR algorithm can achieve a certain reconstruction effect, the edge of the reconstructed image is not sharper than the image obtained using the proposed method. In the second set of experiments, we choose the most clear knife-edge area [114,329] in real images according to the previous experimental procedure. The UAV images, as well as the results and details of the HR image reconstruction, are shown in Figure 15 and the values of metric Q are shown in Table 7. The reconstruction results of the proposed method are the most natural; the image becomes clearer with the increase in the amount of image information. Although the blind SR algorithm can achieve a certain reconstruction effect, the edge of the reconstructed image is not sharper than the image obtained using the proposed method.

Conclusions
In summary, the POCS method with the estimated PSF based on multiple LR images describes a number of key initiatives including the improvement of the accuracy of the PSF of LR images by a novel slant knife-edge method. The validity and reliability of the formula, which derives the relationship between the images before and after downsampling, have been proven. The value of the downsampling multiplied by the PSF of LR images is equal to the estimated PSF of the HR image. The formula can be applied in image restoration and SR reconstruction. The formula can also enhance the clarity in agricultural remote sensing images. Finally, the estimated PSF was combined with the POCS method to improve the accuracy of the SR reconstruction process. Our experimental results show that the deduced formula of PSF is accurate and significant in the development of the restoration and reconstruction processes. However, some problems remain to be solved. The quality of the knife-edge area significantly influences the estimation accuracy of the PSF, leading to some errors.

Conclusions
In summary, the POCS method with the estimated PSF based on multiple LR images describes a number of key initiatives including the improvement of the accuracy of the PSF of LR images by a novel slant knife-edge method. The validity and reliability of the formula, which derives the relationship between the images before and after downsampling, have been proven. The value of the downsampling multiplied by the PSF of LR images is equal to the estimated PSF of the HR image. The formula can be applied in image restoration and SR reconstruction. The formula can also enhance the clarity in agricultural remote sensing images. Finally, the estimated PSF was combined with the POCS method to improve the accuracy of the SR reconstruction process. Our experimental results show that the deduced formula of PSF is accurate and significant in the development of the restoration and reconstruction processes. However, some problems remain to be solved. The quality of the knife-edge area significantly influences the estimation accuracy of the PSF, leading to some errors.

Conflicts of Interest:
The authors declare no conflict of interest.