Uniform Blur Estimation via Modified PSO and Total Variation-Based Adaptive Regularization for Blind Image Restoration

Uniform Blur Estimation via Modified PSO and Total Variation-Based Adaptive Regularization for Blind Image Restoration

Mohammed Boumechaal* Nabil Boukhennoufa

Department of Electronics, Faculty of Technology, University of Batna 2, Batna 05078, Algeria

Department of Computer Science, Faculty of Mathematics and Computer Science, University of Batna 2, Batna 05078, Algeria

Laboratory of Electrical Engineering (LGE), University of M’sila, M’sila 28000, Algeria

Corresponding Author Email: 
m.boumechaal@univ-batna2.dz
Page: 
3073-3086
|
DOI: 
https://doi.org/10.18280/ts.420604
Received: 
19 September 2025
|
Revised: 
12 November 2025
|
Accepted: 
9 December 2025
|
Available online: 
31 December 2025
| Citation

© 2025 The authors. This article is published by IIETA and is licensed under the CC BY 4.0 license (http://creativecommons.org/licenses/by/4.0/).

OPEN ACCESS

Abstract: 

Traditional image restoration methods often assume prior knowledge of the blur kernel, which is rarely available in real-world scenarios. This makes the recovery process a blind image restoration problem. In this work, we introduce an iterative space-adaptive regularization framework for blind image restoration, combining two optimization techniques: a robustly improved Particle Swarm Optimization (PSO) approach with a gradient-based steepest descent technique. The enhanced PSO reliably and efficiently estimates the parameters of a parametric point spread function (PSF) by minimizing an appropriate objective function. Meanwhile, the steepest descent method reconstructs the original image from the corrupted observation. The proposed approach minimizes a cost function comprising a restoration error term and an adaptive regularization term, where the latter incorporates Total Variation (TV) concepts to effectively preserve image edges. In addition, a fine balance between the fidelity and smoothing terms, along with an adaptive selection of weight parameters based on image and noise characteristics. Experimental results on synthetically degraded images highlight the superior performance of the proposed PSO in PSF estimation, attaining excellent quality restorations with improved accuracy and robustness compared to conventional methods.

Keywords: 

enhanced PSO, blind image restoration, space-adaptive regularization, TV, steepest descent, alternating minimization, PSF estimation

1. Introduction

Digital image restoration is an engineering discipline that explores approaches for reconstructing an unknown true scene from corrupted observations. Most traditional methods assume that the convolution operator, which models the blur, is known in advance. However, in practical imaging situations, this is rarely the case [1]. Therefore, this problem is commonly formulated as a blind image deconvolution (BID) task, with the unknown blur modeled as a PSF. Thus, BID involves reconstructing both the original image and the blur kernel from the distorted image, based on partial information about the imaging framework [2]. BID is considered a practical and important method for restoring images. The blind deconvolution problem is encountered in a variety of technical fields [3], including remote sensing, medical imaging, astronomical imaging, microscopy, photography, optics, super-resolution applications, motion tracking technology, and more. However, much of the analysis assumes that the observed image is generated by a linear spatially invariant (LSI) system corrupted by noise [3].

In the literature, there are two primary approaches to blind image deconvolution [2, 3]:

(1) Separate PSF Identification: The PSF is estimated independently of the true image, allowing it to be used later with conventional restoration methods. Here, PSF estimation and image restoration are treated as distinct steps. Often, a parametric blur model is employed in this approach.

(2) Joint Estimation: The identification process is integrated with the restoration algorithm, enabling simultaneous estimation of both the PSF and the original image. However, many practical techniques in this area adopt an alternating optimization strategy, iteratively refining the image and PSF estimates, rather than computing a true joint solution.

Blind deconvolution represents an ill-posed problem. It is often effectively managed through regularization strategies. As a fundamental component of both blind and nonblind deconvolution, regularization plays a crucial role in obtaining stable solutions [3]. However, the phenomenon of ringing effects is often produced around sharp intensity transitions when using global regularization in ill-posed image reconstruction problems [4]. This issue motivates the selection of spatially adaptive methods to mitigate such restoration artifacts. These adaptive approaches typically employ iterative algorithms [1], which have become prevalent in practice due to their ability to efficiently implement adaptive regularization schemes [5]. For comprehensive discussions on these methods, we refer readers to the studies by Campisi and Egiazarian [3] and Satish et al. [6]. Furthermore, the research of Chaudhuri et al. [7] provides an informative investigation into regularization-based blind deconvolution methods, including their convergence properties.

Numerous studies in the literature explore iterative methods for blind image restoration (e.g., [8-11]). In a foundational contribution, You and Kaveh [8] proposed a classical method based on symmetric double regularization (SDR) for space-adaptive BID. Based on the assumption that both the original image and the blur are piecewise smooth, the approach minimizes a cost function including a reconstruction error term and two regularization components—one for the image and the other for the blur. The method also addresses a scaling issue inherent in the cost function through alternating optimization. Solutions using steepest descent and conjugate gradient strategies yield effective results for uniform blur under Gaussian noise at various levels. Further extending the SDR framework, Chen and Yap [9] introduced a parametric double regularization (PDR) method to address the blind deconvolution problem under various constraints and conditions related to the support size and nature of the PSF. In this work, a new cost function and a soft parametric learning term were proposed within a double regularization procedure, and an alternating minimization technique was employed to estimate both the blurring function and the true image iteratively using the conjugate gradient method.

A detailed discussion of regularization-based restoration methods can be found in the papers by Campisi and Egiazarian [3], Chen and Ma [12], and Bertero et al. [13]. Following this, and still within the regularization framework, numerous regularization approaches based on the TV concept have been developed to solve the BID problem, for example, [11, 14-17]. The methods in the papers by Chan and Wong [14] and Money and Kang [16] performed well on specific types of blurred images, such as out-of-focus blurs with simple backgrounds, but were less effective for images containing complex structures. Later, Li et al. [18] extended the approach from the study by Chan and Wong [14] to handle images degraded by more complicated blur patterns. In the present work, the double regularization process is replaced by a split objective formulation, which is minimized using an extended split Bregman iteration.

While regularization-based iterative methods have shown promising results in blind image restoration, they often rely on handcrafted assumptions and require careful parameter tuning. Moreover, accurately modeling image priors within conventional optimization frameworks remains a persistent challenge. To overcome these limitations, recent advances in machine learning have demonstrated superior performance in image restoration tasks. These techniques can broadly be classified into two categories: approaches that directly estimate the sharp image from the degraded input [19, 20] and those that focus on estimating the blur operator [21, 22]. Within the second category, the Deep-Blur method [21] introduces a convolutional neural network (CNN) architecture coupled with a dedicated training strategy to successively estimate both the blur operator and the sharp image from a single blurry input. The framework comprises two components: the Identification Network (IN), which extracts low-dimensional blur parameters, and the Deblurring Network (DN), which utilizes these parameters alongside the degraded image to reconstruct a clear version. Notably, the DN adopts an unrolled architecture to enhance reconstruction efficiency. This design is particularly well-suited for correcting optical systems affected by slowly varying PSF. Complementary to this, another study [22] has investigated the use of inverse blur kernels. In that work, the authors proposed a CNN-based approach that learns the deconvolution process by leveraging the distinctive properties of the inverse kernel, thereby enabling robust handling of spatially varying defocus blur. In a different context, Aljadaany et al. [19] present a deep learning–based solution to the BID problem. Rather than estimating the blur kernel explicitly, the proposed network aims to learn two proximal operators simultaneously: one corresponding to the data fidelity term and the other to the image prior. These operators are modeled using two distinct CNNs and integrated within the Douglas–Rachford optimization framework, offering a principled and flexible alternative to traditional kernel estimation.

On the other hand, traditional approaches to blind image restoration often rely on alternating minimization algorithms, which iteratively estimate the blur kernel and the true image (as in the methods mentioned above). However, classical optimization methods often risk getting stuck in local minima of the cost function, limiting their effectiveness. To address this issue, modern metaheuristic algorithms [23]—such as Particle Swarm Optimization (PSO), genetic algorithms (GA), and many others—have been increasingly employed for blur PSF estimation in blind deblurring, as evidenced by the studies conducted by Dash and Majhi [24], Moghaddam [25], Trubakov and Medvedkov [26], and Lai [27]. A related approach by Sun et al. [28] utilized the histogram entropy of both blurred and restored images as an objective function within a PSO framework. The estimated Gaussian blur PSF was subsequently applied in the Richardson–Lucy algorithm to restore the original image from its distorted (blurred-only) version. Sun et al. [29] proposed a technique that uses PSO to estimate the unknown PSF of out-of-focus blur. Their approach employs the wavelet transform to formulate the objective function, followed by Wiener filtering in the Fourier domain for image restoration. The study focuses solely on blurred images. To enhance optimization robustness, Li and Li [30] introduced an improved PSO strategy that integrated a GA selection mechanism into the standard PSO framework to avoid premature convergence. This hybrid approach was applied to reconstruct images degraded by both Gaussian blur and noise, optimizing the process through a fidelity-based fitness function.

Since metaheuristic algorithms are efficient at solving optimization problems in terms of computational complexity, we adopt in this work a modified version of the standard PSO algorithm to estimate the PSF of the blur instead of the conjugate gradient method used in the works of You and Kaveh [8] and Chen and Yap [9], and the generalized radial basis function network (GRBFN) technique used in the paper by Icho et al. [10]. In this paper, we address the BID problem using a space-adaptive regularization approach. Our objective is to jointly estimate the true image and the associated blurring function. This can be expressed as a constrained optimization problem, inspired by the methods proposed by You and Kaveh [8], Chen and Yap [9], Icho et al. [10], and Sun et al. [28]. We adopt an alternating minimization strategy: the cost function is optimized with respect to the image using the steepest descent method, and with respect to the blur using a modified PSO algorithm. The cost function itself serves as the objective for the PSO.

To enhance the restoration performance, particular attention is given to the selection of parameters within the adaptive regularization framework. Specifically, we integrate the local total variation (LTV) concept into the adaptive function, which enriches the regularization process by capturing the local intrinsic characteristics of the image. A key contribution of our work lies in the modified PSO, which is used to improve the convergence path toward the optimal solution. This is especially beneficial in multidimensional cases, such as ours (uniform blur degradation), where the blur is modeled parametrically as a set of coefficients within a finite support. The improved PSO yields smoother and more accurate solution trajectories, enhancing overall deconvolution precision.

In this study, Section 2 outlines the PSO algorithm and its modified variant. Section 3 explores the components of the space-adaptive regularization process. The iterative optimization process, alternating between identification and restoration steps, as well as the proposed method's framework, are examined in Section 4. Simulation results of the proposed method are analyzed in Section 5. The positive and negative aspects of the suggested method are presented in the conclusion, along with future perspectives.

2. Particle Swarm Optimization Technique

2.1 Standard PSO and background

PSO is a metaheuristic optimization algorithm introduced by Kennedy and Eberhart in 1995 [31], inspired by the collective (social) behavior of biological swarms such as fish schools and bird flocks. As one of the most prominent swarm intelligence paradigms, PSO employs a simple yet effective mechanism to guide particles toward globally optimal solutions. It resembles genetic algorithms in particular aspects, but it is considerably simpler due to the absence of crossover and mutation operators [23]. Due to its ease of implementation and efficiency, PSO has rapidly evolved and been effectively adapted for different real-world optimization problems [32, 33], including image processing [34], controller design [35], system identification [36], robot path planning [37], and more.

A swarm in PSO consists of a number of individuals. Each individual is referred to as a “particle”, which represents a potential solution to a particular problem. Every particle is considered as a point in a D-dimensional space, and it moves through the problem space to a new position in accordance with an updated velocity, which takes into account its previous velocity, the best position known to the particle, and the best known global position. The velocity of the particles determines their direction and travel distance. Each particle's performance is evaluated according to a predetermined fitness function (objective function), which depends on the specific problem being considered.

In PSO, An individual particle i includes three vectors: its location in the D-dimensional search space $X_i= \left(x_{i 1}, x_{i 2}, \ldots, x_{i D}\right)$, the best position that it has independently found $P_i=\left(p_{i 1}, p_{i 2}, \ldots, p_{i D}\right)$, and its velocity $V_i= \left(v_{i 1}, v_{i 2}, \ldots, v_{i D}\right)$. The initial particle positions were generated randomly with a uniform distribution throughout the search space. This last is particularly important for multimodal problems. The initial velocity of a particle can be taken as zero. The particles subsequently move throughout the search space according to a simple set of update equations. At each iteration, the algorithm updates the entire swarm by adjusting the velocity and the position of each particle in every dimension based on the following rules:

$\begin{aligned} v_{i d}(t+1)=w &\ v_{i d}(t)\quad+c_1 r_1(t)\left[p_{i d}(t)\right. \left.\quad-x_{i d}(t)\right]+c_2 r_2(t)\left[G b_d(t)\right.\left.\quad-x_{i d}(t)\right]\end{aligned}$             (1)

$x_{i d}(t+1)=x_{i d}(t)+v_{i d}(t+1)$             (2)

where, Gb is the best solution of all particles, and w denotes the inertia weight of velocity [38], and it can be defined as a positive constant or even a positive linear or nonlinear function of time. The parameters $c_1$ and $c_2$ are two positive constants which denote the acceleration coefficients, $r_1$ and $r_2$ are random numbers uniformly distributed in the range [0,1] introducing stochasticity, and t is the number of generations. The second part of Eq. (1) is the “cognition” part, which denotes the private thinking of the particle itself. The third one is the “social” part, which indicates the collaboration among the particles [23, 38, 39].

For minimization, the best position of each particle is updated at every iteration as follows:

$p_i(t+1)=\left\{\begin{aligned} p_i(t), & \text { if } f\left(x_i(t+1)\right) \geq f\left(p_i(t)\right) \\ x_i(t+1), & \text { otherwise }\end{aligned}\right.$            (3)

where, $f(.)$ denotes the objective function to be optimized.

The global best position, $G b$, is determined by selecting the best among all particles and is given by:

$G b(t+1)=\arg \min _{p_i} f\left(p_i(t+1)\right)$            (4)

Due to its straightforward implementation, the original PSO algorithm has advanced rapidly in recent years, demonstrating strong capabilities to solve real-world optimization problems. However, like most population-based evolutionary computation (EC) algorithms, PSO particles can easily become trapped in local optima —leading to premature convergence—particularly when addressing large-scale or multimodal optimization problems [32, 38]. To overcome these limitations and improve performance, developing enhanced PSO variants remains a practical necessity. For more details on these advancements, see the research by Fang et al. [32].

2.2 Improved PSO algorithm

The primary objective of various PSO variants is to enhance global search capability while preventing premature convergence to local optima. To achieve this, different strategies have been developed, each targeting specific aspects of the algorithm's performance. Among these approaches, some variants focus on modifying acceleration coefficients, as demonstrated in the papers of Jordehi [40] and Liu et al. [41]. Furthermore, the sine–cosine acceleration coefficients (SCAC) updating strategy described in the work of Chen et al. [42] provides an alternative approach for coefficient adjustment. Similarly, another common approach focuses on adjusting the inertia weight, as seen in the work of Bansal et al. [43]. In these variants, specialized functions (linear or nonlinear) are employed to gradually decrease the inertia weight over time, thereby improving convergence speed and refining solutions.

Beyond parameter tuning, recent developments in PSO research have turned toward designing novel velocity updating strategies, an area that has become increasingly prominent [29], as seen in studies by Zhan et al. [33], Wei et al. [44], and Xu et al. [45]. A notable example is the Adaptive Particle Swarm Optimization (APSO) method presented by Zhan et al. [33]. This approach introduces an evolutionary factor (EF) to evaluate the swarm's state, categorizing it into exploration, exploitation, convergence, or jumping-out phases. Building on this, APSO dynamically adjusts velocity update functions based on evolutionary state estimation (ESE), while incorporating an elitist learning strategy to enhance overall performance. More recently, Xu et al. [45] introduced the Adaptive Weighted Delay Velocity (PSO-AWDV) variant, which represents a significant advancement in PSO development. Specifically, this algorithm modifies the traditional velocity update rule by incorporating delayed velocity components and dynamically adjusting inertia weights based on the swarm’s evolutionary state. These enhancements help address typical PSO challenges such as local stagnation and slow convergence. Extensive evaluations on benchmark functions demonstrate that PSO-AWDV achieves superior convergence speed, enhanced stability, and reduced computational complexity compared to other existing PSO variants.

In this study, we aim to enhance the performance of the PSO algorithm by modifying the velocity component while maintaining a balance between global exploration and local exploitation. According to Eq. (1), the magnitude of the velocity significantly influences PSO performance, particularly with regard to trajectory smoothness and convergence behavior. This influence becomes especially prominent in the final optimization stages. To address this, we apply a transformation to the velocity using a function that fulfills the following desirable properties: it should be monotonically increasing, bounded, and odd. To satisfy these properties, we utilize the hyperbolic tangent function, which naturally produces output values within the interval [−1, 1]. To manage the variability of the velocity, we scale the function by a positive constant a, which depends on the search space. Additionally, the velocity (as the function’s input variable) is scaled by another positive constant b. The modified velocity is thus defined as:

$f(x)=a \frac{e^{b x}-e^{-b x}}{e^{b x}+e^{-b x}}$             (5)

The variable $x$ denotes the original velocity and $f(x)$ its transformed.

The symmetric nature of the hyperbolic tangent function ensures that the direction (sign) of the velocity remains intact, preserving the particle's movement orientation. This means positive velocities stay positive, and negative velocities stay negative, maintaining the original search behavior while controlling magnitude. To further analyze the behavior of this transformation, especially for small velocity values where nonlinear effects are minimal, we approximate the hyperbolic tangent using a first-order Taylor expansion around zero. This yield:

$f(x) \approx a b x$            (6)

This linear approximation (linear form in terms of $x$) shows that for low velocity magnitudes, the effect of the transformation is minimal. However, as velocity increases, the function smoothly limits its value, reducing the risk of oscillatory particle trajectories, especially during later iterations. Consequently, this velocity transformation enhances trajectory smoothness, leading to more accurate optimal solutions.

With this transformation, the updated Eq. (1) and Eq. (2) become:

$\begin{array}{r}\tilde{v}_{i d}(t+1)=a \frac{e^{b v_{i d}(t+1)}-e^{-b v_{i d}(t+1)}}{e^{b v_{i d}(t+1)}+e^{-b v_{i d}(t+1)}}=a \frac{1-e^{-2 b v_{i d}(t+1)}}{1+e^{-2 b v_{i d}(t+1)}}\end{array}$             (7)

$\tilde{x}_{i d}(t+1)=x_{i d}(t)+\tilde{v}_{i d}(t+1)$              (8)

In our implementation, the parameters a and b are empirically set to 1 and 0.5, respectively.

Importantly, since this modification does not alter the core PSO strategy, it remains compatible with various PSO variants, making it a flexible enhancement without disrupting existing algorithmic structures. This approach ensures smoother optimization trajectories while maintaining the algorithm’s exploratory and exploitative capabilities.

The modified PSO algorithm is performed according to the following procedure:

Algorithm 1: Modified PSO Procedure

  1. Initialization: Generate a swarm of particles with randomly assigned positions and velocities within the solution space.
  2. Fitness Evaluation: For each particle, compute the value of the optimization fitness function.
  3. Best Positions Adjust: Update the particle’s personal best $p_i$ and the swarm’s global best $G b$, according to Eq. (3) and Eq. (4), respectively.
  4. Velocity and Position Update:
  • Compute the raw velocity using the standard PSO update rule Eq. (1).
  • Apply the velocity transformation using Eq. (7).
  • Update the particle’s position using the transformed velocity Eq. (8).
  1. Termination Check: Repeat steps 2 through 4 until a stopping criterion is met (typically when an acceptable fitness value is reached or the maximum number of iterations is exceeded).
3. Formulation of Blind Image Restoration as a Regularization Problem

3.1 Space-adaptive regularization

In image restoration tasks, it is commonly assumed that the blur can be modeled as a spatially invariant convolution with a blur kernel $h(x)$, and that both the image and noise maintain consistent statistical properties across the image domain. Under these assumptions, if $f(x)$ represents the original noise-free and blur-free image, then the observed degraded image $y(x)$ can be expressed as [5, 8]:

$y(x)=h(x) * f(x)+w(x), \quad x=\left(x_1, x_2\right) \in \Omega$            (9)

where, $\Omega \subset \mathscr{R}^2$ denotes the image domain, and (*) stands for two-dimensional convolution, defined as:

$h(x) * f(x)=\sum_{s \subset D} h(s) f(x-s)$            (10)

In this formulation, $D \subset \mathscr{R}^2$ is the support of the PSF $h(x)$, and $w(x)$ represents additive noise, typically modeled as zero-mean Gaussian white noise.

The goal of image restoration is to estimate $\hat{f}(x)$, an approximation of the original image $f(x)$, given only the degraded image $y(x)$, the blurring function $h(x)$, and some statistical knowledge of the image and noise. However, in many practical scenarios, the PSF is unknown a priori, turning the restoration into a BID problem. In such cases, the task involves estimating both $\hat{f}(x)$ and $\hat{h}(x)$ from $y(x)$, using prior information about $f(x), h(x)$ and $w(x)$. While the degradation is modeled as LSI, the deconvolution method may involve nonlinear behavior, spatial dependence, or simultaneously both [3]. Furthermore, blind restoration is inherently ill-posed, as solutions lack guaranteed uniqueness and stability. To address this, regularization techniques are employed, converting the ill-posed problem into a well-posed one with an acceptable approximate solution [5].

Assuming now that both $f(x)$ and $h(x)$ are unknown, we adopt a space-adaptive regularization method to jointly estimate the original image and the PSF. This is achieved by minimizing the cost function [8]:

$\begin{aligned} & J(\hat{h}, \hat{f})=\sum_{x \in \Omega} w_1(x)[y(x)-\hat{h}(x) *\hat{f}(x)]^2+\lambda \sum_{x \in \Omega} w_2(x)[C(x) * \hat{f}(x)]^2\end{aligned}$              (11)

Subject to the constraint:

$0 \leq \hat{f}(x) \leq 1, \quad x \in \Omega$              (12)

This cost function contains two terms. The first enforces fidelity between the restored image $\hat{f}(x)$ and the observation $y(x)$, weighted by $w_1(x)$. The second promotes smoothness in the solution through a regularization operator $C(x)$, which typically corresponds to a high-pass filter such as the Laplacian. The regularization parameter $\lambda$ balances fidelity to the observed data and the smoothness of the reconstructed image. Selecting an appropriate $\lambda$ is a crucial issue. Various techniques for estimating the regularization parameter are discussed in several works; see, for example, the papers [46, 47]. Specifically, higher values of $\lambda$ enforce greater smoothness, but may introduce ringing artifacts, while smaller values preserve fidelity but may amplify noise [1]. Space adaptivity is achieved by dynamically adjusting the weights:

$w_1(x)$ is assigned higher values in regions with low noise or near sharp transitions to emphasize fidelity.

$w_2(x)$ is increased in smoother regions or areas with higher noise to emphasize regularization.

Details of the regularization components are elaborated in the next subsection. Meanwhile, the constraint in Eq. (12) ensures the restored image remains bounded and nonnegative, reflecting the physical limits of image intensity.

Ultimately, the blind restoration problem, originally posed by Eq. (9), is reformulated as the minimization of the cost function in Eq. (11). Despite its non-convex nature and susceptibility to local minima, alternating minimization offers a viable strategy. By fixing one of the variables and optimizing with respect to the other, the problem becomes convex and solvable, providing an effective solution path for blind deconvolution [3, 8].

3.2 Enhanced design of space-adaptive regularization components

To refine the restoration process, we now delve into the design of the space-adaptive components introduced in Eq. (11). These include the regularization parameter $\lambda$, the spatially varying weights $w_1(x)$ and $w_2(x)$, and the regularization operator $C(x)$.

Traditionally, $\lambda$ is treated as a global constant, applied uniformly across the entire image. However, in adaptive regularization, it is beneficial to allow $\lambda$ to vary spatially, adapting different image regions according to local features. This dynamic approach improves restoration by fine-tuning regularization strength according to regional characteristics. For instance, regions with fine textures or edges benefit from weaker regularization to preserve detail, while smoother areas require stronger regularization to suppress noise.

To achieve this, local image statistics have been employed to guide the spatial adaptation of $\lambda$, mimicking the behavior of the human visual system and yielding promising results [48]. Building on this idea, we propose linking the concept of local total variation (LTV) to the adaptive function $w_2(x)$, thereby enriching the regularization process with more nuanced information about local image features. Notably, LTV is particularly effective at distinguishing between smooth, textured, and edge regions, making it a valuable component in adaptive weighting. Further background regarding these advancements can be found in studies by Liu et al. [11] and Bredies et al. [49].

In our formulation, spatial adaptivity is governed jointly by the local standard deviation (LSD) and the LTV. Specifically, the regularization parameter $\lambda$ (constant) is modulated by the adaptive weight $w_2(x)$, defined as:

$w_2(x)=\frac{1}{\gamma+u(x)}, \quad x \in \Omega$            (13)

where,

$u(x)=$ local_std $(x)$. local_TV $(x)$              (14)

Here, $u(x)$ captures the local adaptivity within a 3×3 neighborhood around pixel $x$ in the estimated image $\hat{f}(x)$, computed as the product of LSD and LTV. This formulation ensures that $w_2(x)$ assumes higher values in smooth regions where noise is more perceptible and lower values in textured or edge regions where detail preservation is critical.

Furthermore, the parameter $\gamma$ plays a pivotal role in shaping the dynamic range of $w_2(x)$. Ideally, so as not to disrupt the role of $w_2(x)$ in different scenarios $\gamma$ should satisfy $\gamma \ll u(x)$ in edge-rich regions and $\gamma \gg u(x)$ in flat regions. As a result, $w_2(x)$ lies within the range:

$w_2(x) \in\left[\frac{1}{\max _{x \in \Omega}\{u(x)\}}, \quad \frac{1}{\gamma}\right]$             (15)

To further enhance adaptivity, we propose relating $\gamma$ to global image statistics, such as the overall variance. Specifically, a high global variance suggests the need for a wider dynamic range in $w_2(x)$, achieved by selecting a smaller $\gamma$. Conversely, a low global variance implies a narrower dynamic range, warranting a larger $\gamma$. Thus, the design of the function $u(x)$ is central to achieving precise spatial control over regularization and the key challenge lies in optimizing $\gamma$ to align with image features.

Regarding the fidelity weight, $w_1(x)$ is set as a constant value $\alpha$ across the entire image:

$w_1(x)=\alpha, \quad x \in \Omega$            (16)

The parameter $\alpha$ works in the opposite way to $\lambda$, but together they help balance data fidelity and image smoothness. Both $\alpha$ and $\lambda$ are selected based on the signal-to-noise ratio (SNR) and are held constant for each restoration scenario. In our experiments, these parameters -along with $\gamma$ - are tuned empirically to ensure robust adaptivity across various noise levels and image types.

Finally, the smoothness requirement on the solution is enforced by the regularization operator $C$, which is chosen as a high-pass filter. Specifically, we use the two-dimensional Laplacian operator [5]:

$c=\left[\begin{array}{ccc}0 & 1 & 0 \\ 1 & -4 & 1 \\ 0 & 1 & 0\end{array}\right]$             (17)

The term associated with this operator corresponds to a high-pass filtered version of the image $\hat{f}(x)$, and thus seeks to minimize the amount of high-frequency energy in the restored image [1].

4. Optimization Scheme

Building on the regularization model presented in Section 3, we now introduce a blind image restoration framework designed to jointly estimate the original image and the unknown blur kernel (PSF) from a degraded observation. The proposed approach integrates gradient-based iterative image estimation (Subsection 4.1) with PSO-based blur identification (Subsection 4.2) through an alternating minimization strategy (Subsection 4.3).

4.1 Restoration via steepest descent

When the PSF is known, the restoration problem reduces to minimizing a cost function with respect to the image. Several optimization techniques exist for this task, including gradient-based methods, expectation-maximization (EM), prediction error-based techniques, and least squares approaches [2, 8].

In our method, we adopt the steepest descent algorithm to minimize the cost function in Eq. (11) with respect to the image estimate $\hat{f}^k(x)$, defined over an M×N grid. The update rule is given by:

$\hat{f}^{k+1}(x)=\hat{f}^k(x)-\frac{\mu}{M N} p^k(x)$             (18)

where, the gradient $p^k(x)$ is calculated by:

$\begin{gathered}p^k(x)=\frac{\partial J\left[h, \hat{f}^k(x)\right]}{\partial \hat{f}^k(x)} \\ p^k(x)=-2\left\{w_1(x)\left[y(x)-h(x) * \hat{f}^k(x)\right]\right.* h(-x)\}+2 \lambda\left\{w_2(x)\left[C(x) * \hat{f}^k(x)\right]\right\}* C(-x)\end{gathered}$             (19)

Here $\hat{f}^k(x)$ is the restoration after $k$ iterations, initialized as $\hat{f}^0(x)=y(x)$. The step size $\mu$, a small positive number, balances convergence speed and stability. Since the cost function is convex, the algorithm converges to the global minimum. This method also allows for easy incorporation of constraints on the image and PSF. Iterative methods avoid explicit inverse operations, allow real-time monitoring, and enable noise control through constraints. Spatial adaptivity and parameter updates can also be incorporated during iteration [1, 5].

4.2 Blur identification via PSO

While the previous subsection assumes the PSF $h(x)$ is known, in real-world scenarios it must be estimated from noisy, blurred data. Various methods exist for this purpose [8, 10], including metaheuristic techniques [27, 29]. We adopt a PSO-based approach to estimate the PSF by minimizing the fidelity term of the cost function in Eq. (11). A common model of PSF is a uniform 2D blur, which approximates out-of-focus degradation and is widely used in simulations [1]. The blur is defined as:

$h(i, j)=\left\{\begin{array}{cc}\frac{1}{(L)^2}, & \text { if }-l \leq i, j \leq l \\ 0, & \text { otherwise }\end{array}\right.$            (20)

where, L=2l+1 is an odd integer representing the blur kernel size.

To estimate the PSF, we define it as a square matrix of size L × L, requiring the estimation of $L^2$ coefficients. The PSO-based approach seeks the coefficient values that minimize the cost function in Eq. (11), subject to the constraints:

$\hat{h}(x) \geq 0, \quad x \in D$           (21)

and

$\sum_{x \in D} \hat{h}(x)=1$           (22)

These constraints ensure the PSF is non-negative and preserves the image’s mean intensity. Incorporating prior knowledge about the blur helps reduce ambiguity in the solution space [2].

Each particle, in PSO, represents a candidate PSF estimate, modeled as a point in an n -dimensional space. The swarm is defined as:

$X=\left(x_1, x_2, \ldots, x_m\right)^T=\left[\begin{array}{ccc}x_{11} x_{12} & \cdots & x_{1 n} \\ x_{21} x_{22} & \cdots & x_{2 n} \\ \vdots & \ddots & \vdots \\ x_{m 1} x_{m 2} & \cdots & x_{m n}\end{array}\right]$               (23)

where, m is the swarm size and n = $L^2$ is the number of PSF coefficients (dimension of the search space).

Algorithm 2: PSF Identification

  1. Input: observed image $y(x)$, initial image estimate $f(x)$.
  2. Initialize swarm with random values in [0, 1].
  3. Evaluate the fidelity term of Eq. (11) for each particle.
  4. Update personal best p, global best Gb, and velocity using Eq. (3), Eq. (4), and Eq. (1).
  5. Modify velocity according to Eq. (7).
  6. Move particles to new positions using Eq. (8).
  1. Repeat steps 3–6 until convergence or maximum iterations are achieved.

To initialize the swarm, one particle is set to a uniform value (e.g., [0.1, 0.1, …, 0.1] in the case of $L=3$) to reflect the average blur, while others are randomly sampled from the interval [0, 1]. The PSO algorithm then searches for the optimal solution in the vicinity of this reference point. Upon convergence, the final global best Gb provides the estimated PSF.

4.3 Alternate minimization for blind restoration

Having addressed image restoration and blur identification separately, we now integrate both into a unified blind restoration framework. This is achieved through an alternating minimization strategy, where each component is iteratively refined based on the other.

The blind restoration process consists of two main steps:

BI-step (Blur Identification): Given an image estimate $\hat{f}(x)$, use PSO (Algorithm 2) to estimate the blur $\widehat{h}(x)$ by minimizing the fidelity term of Eq. (11).

R-step (Image Restoration): Using the estimated blur $\widehat{h}(x)$, apply the steepest descent method (Eq. (18)) to update the estimated image.

By alternating between these steps, blind image restoration is achieved.

Algorithm 3: Main Restoration Algorithm

  1. Initialization: Set the observed image $y(x)$ as the initial estimate $\hat{f}(x)$ of the true image.
  2. Initial PSF Estimation: Given $\hat{f}(x)$ and $y(x)$, minimize the fidelity term in Eq. (11) using the PSO procedure (Algorithm 2) to obtain the initial estimated blur PSF coefficients.
  3. Intermediate Image Restoration: Using the estimated PSF $\widehat{h}(x)$ and the observed image $y(x)$, apply the update procedure (Eq. (18)) for a fixed number of iterations to compute an intermediate restored image.
  4. PSF Re-estimation: Regard the intermediate restored image as the new estimate of the true image, and re-estimate the blur PSF using the PSO procedure (Algorithm 2).
  1. Termination Check: Repeat steps 3 and 4 until either no significant improvement is observed in the image restoration step (Step 3), or the maximum number of main loop iterations is reached.

In steps 2 and 4, the swarm is initialized with uniformly distributed values in [0, 1], including the best solution from the previous stage. To ensure smoothness, we can incorporate the constraint $\operatorname{std}(\mathrm{Gb})<\varepsilon$ in step (2) (for example $\varepsilon=0.07$ ) and std(Gb) denotes the standard deviation of the global best solution coefficients. Under this alternating framework, each step optimizes a quadratic cost function, resulting in a simple and efficient algorithm. Moreover, by employing distinct optimization techniques for image and blur estimation, the scale ambiguity problem discussed in the work of You and Kaveh [8] is effectively mitigated.

5. Results and Discussions

This section presents the blind restoration of images degraded by uniform blur under various conditions, including different blur sizes and noise levels. The objective is to demonstrate the effectiveness of the proposed method through simulation examples and to analyze the quality of both the PSF and image estimations across multiple degradation scenarios. The Improved Signal-to-Noise Ratio (ISNR) is used to evaluate restoration performance. ISNR is defined as [48]:

$I S N R=10 \log _{10} \frac{\|y-f\|^2}{\|\hat{f}-f\|^2}$            (24)

where, $y, f$, and $\hat{f}$ represent the degraded, true, and recovered images respectively.

This study involves two principal simulation scenarios. The first simulation compares the performance of our proposed PSO variant with both the standard PSO [39] and the PSO-AWDV variant [45], a recent advancement in the field. The second simulation evaluates the effectiveness of the adaptive regularization model applied to image restoration, using a consistent framework to enable fair comparisons with classical established methods such as SDR [8] and PDR [9].

For all experiments, we set the following parameters based on empirical tuning: K = 30 for PSO iterations, m = 20 for swarm size, Na = 10 for main loop alternations, and Nt = 8 for image update iterations.

5.1 Simulation 1

In the first set of experiments, we use the well-known “Cameraman” image (256×256 pixels, 8-bit grayscale) as the original image (Figure 1(a)).

The original image is degraded by applying a 3×3 uniform blur and adding Gaussian noise at SNR = 30 dB, resulting in the image shown in Figure 1(b). To estimate both the PSF and the original image, we apply the proposed algorithm to the degraded image and repeat the process ten times. The only prior knowledge is the PSF size. Each simulation performs a total of 80 restoration iterations (Na = 10, Nt = 8). Moreover, the parameter settings are given as: $c_1=2, c_2=1$ for the acceleration coefficients, w=0.9 for the inertia weight, and (a,b) =(1,0.5) for the hyperbolic tangent parameters. The PSO-AWDV variant uses the parameter settings reported in [45]. And about the R-step, we put μ = 60 and γ = 0.005, which depends on the image at hand. While the others, namely α and $\lambda$, depend on the used SNR value. All parameters were heuristically selected to optimize restoration performance. The obtained results are presented in Table 1.

Where $\sum h_i$ denotes the sum of estimated PSF coefficients per run, $\sigma_h$ denotes the standard deviation of those coefficients, and MSE denotes the mean squared error between the estimated and true blur kernels.

(a) Original image

(b) Degraded image

(c) Restored image using the proposed method

Figure 1. Blind restoration results for the “Cameraman” image degraded by a 3×3 uniform blur and additive white Gaussian noise (AWGN) at 30 dB SNR

Table 1. Results of the proposed method for the “Cameraman” image degraded by a 3×3 uniform blur with 30 dB AWGN

Experiment Number

PSO VARIANTS

Standard PSO

PSO-AWDV

Our PSO Variant

ISNR

$\sum h_i$  

$\sigma_h$

ISNR

$\sum h_i$  

$\sigma_h$

ISNR

$\sum h_i$  

$\sigma_h$

1

5.8117

0.9988

0.0081

5.9331

1.0019

0.0095

6.1261

1.0012

0.0050

2

5.2095

1.0009

0.0150

3.2780

1.0028

0.0267

5.9330

1.0011

0.0092

3

3.7923

1.0022

0.0318

5.8926

1.0020

0.0087

6.1366

1.0020

0.0066

4

2.1825

1.0020

0.0567

5.0248

1.0018

0.0228

6.2215

1.0017

0.0039

5

4.5550

1.0065

0.0292

6.1582

1.0016

0.0046

6.0487

1.0018

0.0071

6

5.8734

1.0022

0.0101

5.3202

1.0013

0.0185

5.9907

1.0014

0.0057

7

2.9327

1.0045

0.0551

4.7357

1.0014

0.0214

5.9844

1.0012

0.0065

8

4.7334

1.0040

0.0264

5.6745

1.0018

0.0105

5.9214

1.0015

0.0095

9

2.5687

1.0017

0.0721

4.3074

1.0018

0.0260

5.9745

1.0017

0.0081

10

2.2550

1.0040

0.0474

6.0600

1.0019

0.0052

5.9257

1.0019

0.0080

Average

3.9914

---

---

5.2384

---

---

6.0263

---

---

Std

1.4406

---

---

0.9212

---

---

0.1034

---

---

MSE

9.40e-04

1.19e-04

3.10e-05

As illustrated in Table 1, the proposed PSO variant consistently achieves the lowest standard deviation in PSF coefficients ($\sigma_h$ ranging from 0.0039 to 0.0095), indicating superior smoothness and stability in blur kernel estimation. This smoothness directly contributes to the highest average ISNR of 6.03 dB, along with the lowest ISNR standard deviation of 0.10, confirming the method’s robustness and reliability in image restoration. Moreover, in terms of MSE, the table shows that our PSO variant achieves the lowest error, reducing the MSE by nearly 74% compared to PSO‑AWDV and by about 30 times compared to Standard PSO. These results demonstrate that our proposed PSO variant provides the most accurate blur estimation among the three tested versions, highlighting the effectiveness of the modifications introduced. In comparison, the standard PSO exhibits significantly higher variability in PSF estimation, with $\sigma_h$ reaching up to 0.0721, which reflects poor smoothness and instability. Consequently, it yields the lowest restoration performance, with an average ISNR of 3.99 dB and a high standard deviation of 1.44, indicating inconsistent results across runs. The PSO-AWDV variant shows moderate improvement over the standard PSO, with $\sigma_h$ values up to 0.0267, and achieves an average ISNR of 5.24 dB with a standard deviation of 0.92. While it performs better than the standard PSO, it remains less stable and less accurate than the proposed method.

Overall, the proposed PSO demonstrates clear superiority in both PSF estimation and restoration quality, combining smooth kernel estimation with high accuracy and exceptional consistency.

(a) Evolution of ISNR

(b) Evolution of the maximum LSD function $\max _{x \in \Omega}$ local_std $(x)$

Figure 2. Restoration dynamics for the “Cameraman” image during the iterative process

Figure 2(a) depicts the evolution of ISNR throughout the restoration process. The horizontal axis corresponds to the number of iterations performed using the steepest descent method in the R-step. The step size μ is carefully selected to ensure rapid convergence of the recursive update (Eq. (18)) while maximizing ISNR.

As shown, the results indicate that the iterative process stabilizes within approximately 70 iterations. Figure 2(b) illustrates the progression of the function $max _{x \in \Omega}$ local_std $(x)$, which corresponds to the maximum local standard deviation computed at each iteration over a 3×3 window centered at pixel $x$ in the restored image $\hat{f}(x)$. This metric reflects local contrast and texture variations. By incorporating LTV into the adaptive weight $w_2(x)$, the proposed method demonstrates its ability to recover fine structural details and improve the preservation of edges and textured regions.

As seen in the restored image in Figure 1(c), the proposed approach achieves effective restoration by recovering fine details while minimizing ringing artifacts and noise in relevant areas of the image.

5.2 Simulation 2

This set of experiments explores the trade-off between data fidelity and solution smoothness within the proposed adaptive regularization framework. The analysis focuses on the influence of two key parameters-α, which controls fidelity to the observed data and λ, which governs the strength of regularization-on restoration performance under varying noise conditions. Additionally, spatial adaptivity is achieved through the weight function $w_2(x)$, which captures local image characteristics by combining LTV and LSD, allowing the method to adapt to diverse feature types across the image domain.

To comprehensively evaluate the method’s performance, four experimental scenarios are considered:

Noise-level sensitivity: Assessing robustness across a range of SNR.

Image-type sensitivity: Evaluating restoration quality on images with varying textures and frequency content.

Time complexity: Measuring runtime performance across blur kernels and image sizes.

Comparative analysis: Benchmarking against existing methods (SDR and PDR) under severe blur conditions.

5.2.1 Sensitivity to noise levels

To assess the algorithm’s robustness to noise, the original “Cameraman” image Figure 1(a) is degraded using a 3×3 uniform blur and AWGN at SNR levels of 20, 30, 40, 50, and 60 dB. Some of corresponding degraded images are shown in Figures 3(a), 3(b), and 3(c), while the restored results from the final run of each scenario are presented in Figures 4(a), 4(b), and 4(c).

Table 2. Restoration results using the proposed method for the “Cameraman” imagedegraded by a 3×3 uniform blur under varying levels of AWGN

SNR (dB)

$\alpha\left(w_1\right)$

$\lambda$

Average ISNR (dB)

$\sigma_{\text {ISNR }}$

20

300

0.2

3.2746

0.0519

30

3000

0.2

6.0263

0.1034

40

5000

0.08

8.6135

0.1118

50

5000

0.08

9.4424

0.1515

60

5000

0.08

9.5018

0.1221

For each noise level, the algorithm is executed ten times, and the average ISNR and its standard deviation are computed. The results are provided in Table 2.

(a) SNR=20 dB

(b) SNR=40 dB

(c) SNR=60 dB

Figure 3. “Cameraman” image degraded by a 3×3 uniform blur and AWGN at different SNR levels

The data reveal a clear trend: ISNR increases as noise decreases, confirming the algorithm’s effectiveness in cleaner environments. This behavior reflects the dynamic role of the regularization parameters:

At higher SNRs, larger values of α and smaller values of λ emphasize deblurring, allowing the algorithm to recover finer details.

At lower SNRs, smaller α and larger λ shift the focus toward smoothing, helping to suppress noise and stabilize the solution.

These observations show the importance of tuning α and λ according to the noise level to achieve optimal restoration performance.

It is also noteworthy that, for SNR ≥ 40 dB, the restored images exhibit strong deblurring. However, ringing artifacts may appear in smooth regions near sharp transitions such as the boundary between the “Cameraman” and the camera support in Figure 4(c). These artifacts are likely caused by inaccuracies in blur estimation, compression effects, or residual noise [18].

(a) SNR=20 dB

(b) SNR=40 dB

(c) SNR=60 dB

Figure 4. Blind restoration results of the degraded images in Figure 3

5.2.2 Sensitivity to image characteristics

In this study, we assess the adaptability of the proposed method to images with varying structural and textural properties. Specifically, we apply the algorithm to the well-known “Lena” and “Baboon” images, each sized 256×256 pixels in 8-bit grayscale. Both images (Figure 5(a) and Figure 6(a)) are degraded using a 3×3 uniform blur and 30 dB AWGN, as shown in Figure 5(b) and Figure 6(b).

(a) Original image

(b) Degraded image

(c) Restored image using the proposed method

Figure 5. Blind restoration results for the “Lena” image degraded by a 3×3 uniform blur and AWGN at 30 dB SNR

(a) Original image

(b) Degraded image

(c) Restored image using the proposed method

Figure 6. Blind restoration results for the “Baboon” image degraded by a 3×3 uniform blur and AWGN at 30 dB SNR

The corresponding restored results are presented in Figure 5(c) and Figure 6(c), where the quantitative evaluations are summarized in Table 3.

In this table, the proposed method exhibits strong detail preservation despite the presence of prominent high-frequency components—such as hair strands in “Baboon” and feather textures in “Lena”—the proposed method demonstrates strong detail preservation. This is primarily attributed to the incorporation of LTV into the adaptive weight function $w_2(x)$, which enables spatially responsive regularization.

The parameter γ, which modulates the behavior of $w_2(x)$, is observed to be inversely proportional to the image variance, allowing the method to adapt to varying texture complexity.

Table 3. Restoration results using the proposed method on various images degraded by a 3×3 uniform blur and 30 dB AWGN

Images

$w_s$

$\gamma$

Average ISNR (dB)

Cameraman

0.0598

0.005

6.0263

Lena

0.0352

0.008

4.2805

Baboon

0.0272

0.08

4.5261

Interestingly, the “Cameraman” image yields the highest ISNR, which can be explained by its relatively simple structure and lower frequency content. Surprisingly, “Baboon” achieves a slightly higher ISNR than “Lena”, despite Lena’s smoother appearance. This counterintuitive result may be due to the presence of sharp transitions in Lena’s feather region, which pose greater challenges for accurate restoration. These plots confirm the algorithm’s ability to dynamically adapt its regularization strategy to preserve fine details in both textured and edge-rich regions.

To further illustrate the method’s responsiveness to local image features, Figures 7(a) and 7(b) visualize the evolution of the maximum LSD during the restoration process for Lena and Baboon, respectively. These plots confirm the algorithm’s ability to dynamically adapt its regularization strategy to preserve fine details in both textured and edge-rich regions.

(a) For “Lena” image

(b) For “Baboon” image

Figure 7. Evolution of the maximum LSD function, $\max _{x \in \Omega}$ local_std $(x)$, during the iterative process

5.2.3 Runtime analysis

To evaluate the performance of the method in terms of computational complexity, we report runtime statistics for the “Lena” image at resolutions of 256×256 and 512×512 pixels. Both images were degraded using blur kernels of sizes 3×3, 5×5 and 7×7. The runtime measurements were obtained on an Intel Core i3-2310M CPU (2.1 GHz, 4 GB RAM, Windows 7) and are summarized in Table 4.

Table 4. Runtime (in seconds) statistics illustrating the scaling behavior of the modified PSO algorithm and the gradient‑based method with varying image and kernel sizes

Size of image

Size of Kernel

3×3

5×5

7×7

Our PSO

Gradient method

Our PSO

Gradient method

Our PSO

Gradient method

256×256

14

25

32

26

57

27

512×512

51

110

124

117

222

125

These results clearly illustrate the computational cost and the scaling behavior of our method. From the table, we observe that image size has a considerable impact on the runtime of both optimization methods, reflecting the higher computational load required to process more pixels. In contrast, blur size influences the runtime of the PSO algorithm much more strongly than that of the gradient method, which shows only a modest increase. As the PSF size grows, the number of coefficients to be estimated also increases, thereby raising the computational complexity, as is typical for metaheuristic methods.

5.2.4 Comparison with existing methods

To benchmark the proposed method against existing approaches, we perform blind deconvolution on the “Lena” image degraded by uniform blur kernels of size 5×5 and 7×7, under noiseless and 30 dB noise conditions (see Figure 8).

We compare our method with SDR and PDR, which use similar regularization schemes. The results are summarized in Table 5, using ISNR as the evaluation metric. SDR and PDR results are taken from the paper by Chen and Yap [9].

(a) 5×5 blur with SNR = 30 dB

(b) 5×5 blur without noise

(c) 7×7 blur with SNR = 30 dB

(d) 7×7 blur without noise

Figure 8. “Lena” image degraded by uniform blur under noiseless and SNR = 30 dB

In the noiseless scenario, the proposed method achieves the highest ISNR values, significantly outperforming both SDR and PDR. Under noisy conditions, it still surpasses SDR and performs slightly better than PDR. These results highlight the method’s robustness to blur severity and its ability to maintain high restoration quality even in challenging environments.

Table 5. Blind deconvolution results for the “Lena” image degraded by severe uniform blurs of various sizes, using the SDR and PDR methods alongside the proposed algorithm

 

ISNR (Noise at 30 dB)

ISNR (Noiseless)

Size of PSF

SDR [8]

PDR [9]

Proposed method

SDR [8]

PDR [9]

Proposed method

5×5

1.67

3.64

3.97

2.00

3.98

5.69

7×7

1.02

3.39

3.53

1.25

3.60

4.45

(a) Restoration from 5×5 blur for SNR = 30 dB

(b) Restoration from 5×5 blur without noise

(c) Restoration from 7×7 blur for SNR = 30 dB

(d) Restoration from 7×7 blur without noise

Figure 9. Restored “Lena” images using the proposed method, corresponding to the degradation cases in Figure 8

Although the ISNR values obtained with our method may appear relatively low in certain cases, the visual quality of the restored images, particularly under low-noise conditions, is in fact quite satisfactory. As shown in Figure 9, the visual quality of the restored images is better than what the numerical metrics alone may indicate, showing that our method works well in real-world applications.

Overall, the experimental findings confirm the efficiency of the proposed adaptive regularization framework. Specifically, the method:

Balances fidelity and smoothness across a wide range of noise levels.

Adapts dynamically to noise and image structure, preserving fine details in textured and edge-rich regions.

Outperforms existing techniques under severe blur, particularly in low-noise conditions.

These results confirm the effectiveness of the proposed PSO variant and its suitability for high-quality blind image restoration across diverse scenarios.

However, like most iterative methods, our approach has a significant drawback in terms of execution time, which limits its applicability to large images, complex blurs, and real-time scenarios. This explains why metaheuristic methods are often employed to identify blurs modeled by only a few parameters, such as the standard deviation of a Gaussian blur or the angle and length of a motion blur.

Another limitation of this method lies in its generalization to more complex types of blurs, such as non‑uniform or spatially varying blur. Addressing this challenge remains difficult in blur identification, as highlighted in several studies (see, for example the study by Almeida and Figueiredo [46]). Methods that perform well for a specific blur type often fail for others, since they require adapting the identification strategy. Even recent methods based on deep learning encounter this limitation [21], as their performance is often restricted to the blur types seen during training. Hybrid identification methods based on metaheuristics, however, may help mitigate this issue, as they allow the straightforward incorporation of additional constraints and blur priors.

6. Conclusions

This work proposed a robust and adaptive framework for blind image restoration, combining a space-adaptive regularization model with a PSO strategy. The enhanced PSO variant demonstrated clear advantages in estimating blur kernels with smooth, accurate and stable parameters, as evidenced by consistently low standard deviations in PSF coefficients across multiple trials. This smoothness directly enhances restoration quality, yielding the highest ISNR values and the lowest variability among competing methods (standard PSO and the recent PSO-AWDV variant). Such stability in blur estimation ensures consistent and reliable performance, even under varying noise levels.

Beyond optimization, the adaptive regularization model plays a pivotal role in tailoring the restoration process to different image and noise conditions. The fidelity weight α and regularization parameter λ were effectively tuned in relation to SNR values, achieving a balanced trade-off between deblurring and noise suppression. Moreover, the spatially varying weight $w_2(x)$, governed by LTV and LSD, allowed the method to respond dynamically to diverse image features, preserving fine details in textured and edge-rich regions across a range of degradation scenarios.

Comprehensive experimental evaluations, including comparisons with classical approaches such as SDR and PDR, confirmed the effectiveness and versatility of the proposed method. The results highlight its capability to handle fairly severe blur and diverse image types. When the weighting function parameters are selected appropriately based on to the image characteristics, the method demonstrates high efficiency, making it suitable for a wide range of practical imaging applications. A key limitation of the method is its computational complexity when processing large‑scale images and under severe blur conditions.

Future research could focus on extending the framework to address non-uniform or spatially variant blur and exploring hybrid approaches that incorporate deep learning to further enhance flexibility and performance.

  References

[1] Banham, M.R., Katsaggelos, A.K. (1997). Digital image restoration. IEEE Signal Processing Magazine, 14(2): 24-41. https://doi.org/10.1109/79.581363

[2] Kundur, D., Hatzinakos, D. (1996). Blind image deconvolution. IEEE Signal Processing Magazine, 13(3): 43-64. https://doi.org/10.1109/79.489268

[3] Campisi, P., Egiazarian, K. (2017). Blind Image Deconvolution: Theory and Applications. CRC Press. https://doi.org/10.1201/9781420007299

[4] Lagendijk, R.L., Biemond, J., Boekee, D.E. (1988). Regularized iterative image restoration with ringing reduction. IEEE Transactions on Acoustics, Speech, and Signal Processing, 36(12): 1874-1888. https://doi.org/10.1109/29.9032

[5] Bovik, A.C. (2005). Handbook of Image and Video Processing (2nd ed.). Elsevier Academic Press, London. 

[6] Satish, P., Srikantaswamy, M., Ramaswamy, N.K. (2020). A comprehensive review of blind deconvolution techniques for image deblurring. Traitement du Signal, 37(3): 527-539. https://doi.org/10.18280/ts.370321

[7] Chaudhuri, S., Velmurugan, R., Rameshan, R. (2014). Blind Image Deconvolution: Methods and Convergence. Springer. https://doi.org/10.1007/978-3-319-10485-0

[8] You, Y.L., Kaveh, M. (1996). A regularization approach to joint blur identification and image restoration. IEEE Transactions on Image Processing, 5(3): 416-428. https://doi.org/10.1109/83.491316

[9] Chen, L., Yap, K.H. (2005). A soft double regularization approach to parametric blind image deconvolution. IEEE Transactions on Image Processing, 14(5): 624-633. https://doi.org/10.1109/TIP.2005.846024

[10] Icho, K., Iiguni, Y., Maeda, H. (2004). Nonlinear image restoration using a radial basis function network. EURASIP Journal on Applied Signal Processing, 2004(20): 2441-2455. https://doi.org/10.1155/S1110865704408166 

[11] Liu, H., Gu, J., Meng, M.Q.H., Lu, W.S. (2016). Fast weighted total variation regularization algorithm for blur identification and image restoration. IEEE Access, 4: 6792-6801. https://doi.org/10.1109/ACCESS.2016.2516949

[12] Chen, F., Ma, J. (2009). An empirical identification method of Gaussian blur parameter for image deblurring. IEEE Transactions on Signal Processing, 57(7): 2467-2478. http://doi.org/10.1109/TSP.2009.2018358

[13] Bertero, M., Boccacci, P., De Mol, C. (2021). Introduction to Inverse Problems in Imaging (2nd ed.). CRC Press. https://doi.org/10.1201/9781003032755

[14] Chan, T.F., Wong, C. (1998). Total variation blind deconvolution. IEEE Transactions on Image Processing, 7(3): 370–375. https://doi.org/10.1109/83.661187

[15] He, L., Marquina, A., Osher, S.J. (2005). Blind deconvolution using TV regularization and Bregman iteration. Journal of Imaging Systems and Technology, 15(1): 74-83. https://doi.org/10.1002/ima.20040

[16] Money, J.H., Kang, S.H. (2008). Total variation minimizing blind deconvolution with shock filter reference. Image and Vision Computing, 26(2): 302-314. https://doi.org/10.1016/j.imavis.2007.06.005

[17] Wang, Y., Yang, J., Yin, W., Zhang, Y. (2008). A new alternating minimization algorithm for total variation image reconstruction. SIAM Journal on Imaging Sciences, 1(3): 248-272. https://doi.org/10.1137/080724265 

[18] Li, W.H., Li, Q.L., Gong, W.G., Tang, S. (2012). Total variation blind deconvolution employing split Bregman iteration. Journal of Visual Communication and Image Representation, 23(3): 409-417. https://doi.org/10.1016/j.jvcir.2011.12.003

[19] Aljadaany, R., Pal, D.K., Savvides, M. (2019). Douglas-Rachford networks: Learning both the image prior and data fidelity terms for blind image deconvolution. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pp. 10235-10244. 

[20] Abuolaim, A., Brown, M.S. (2020). Defocus deblurring using dual-pixel data. arXiv:2005.00305 https://doi.org/10.48550/arXiv.2005.00305

[21] Debarnot, V., Weiss, P. (2024). Deep-Blur: Blind identification and deblurring with convolutional neural networks. Biological Imaging, 4: e13.

[22] Son, H., Lee, J., Cho, S., Lee, S. (2021). Single image defocus deblurring using kernel-sharing parallel atrous convolutions. In: Proceedings of the IEEE International Conference on Computer Vision (ICCV), pp. 2642–2650. http://dx.doi.org/10.1109/ICCV48922.2021.00264

[23] Yang, X.S. (2010). Engineering Optimization: An Introduction with Metaheuristic Applications. John Wiley & Sons. https://doi.org/10.1002/9780470640425

[24] Dash, R., Majhi, B. (2009). Particle swarm optimization based regularization for image restoration. In Proceedings of the World Congress on Nature and Biologically Inspired Computing (NaBIC), Coimbatore, India, pp. 1253–1257. https://doi.org/10.1109/NABIC.2009.5393754

[25] Moghaddam, M.E. (2008). Out of focus blur estimation using genetic algorithm. Journal of Computer Science, 4(4): 298-304. https://doi.org/10.3844/jcssp.2008.298.304

[26] Trubakov, A.O., Medvedkov, N.V. (2022). Image restoration using genetic algorithms. Programming and Computer Software, 48(3): 199-207. https://doi.org/10.1134/S0361768822030112

[27] Lai, Y.C. (2011). PSO-based estimation for Gaussian blur in blind image deconvolution problem. In Proceedings of the IEEE International Conference on Fuzzy Systems (FUZZ-IEEE), Taipei, Taiwan, pp. 1143-1148. https://doi.org/10.1109/FUZZY.2011.6007681

[28] Sun, T.Y., Liu, C.C., Jheng, Y.P., Jheng, J.H., Tsai, S.J., Hsieh, S.T. (2008). Blind image deconvolution via particle swarm optimization with entropy evaluation. Proceedings of the Eighth International Conference on Intelligent Systems Design and Applications, 2: 256-270. https://doi.org/10.1109/ISDA.2008.238

[29] Sun, T.Y., Ciou, S.J., Liu, C.C., Huo, C.L. (2009). Out-of-focus blur estimation for blind image deconvolution: Using particle swarm optimization. In 2009 IEEE International Conference on Systems, Man and Cybernetics, San Antonio, USA, pp. 1627-1632. https://doi.org/10.1109/ICSMC.2009.5346769

[30] Li, N., Li, Y.X. (2011). Image restoration using improved particle swarm optimization. In 2011 International Conference on Network Computing and Information Security, Guilin, China, pp. 394-397. https://doi.org/10.1109/NCIS.2011.86

[31] Kennedy, J., Eberhart, R.C. (1995). Particle swarm optimization. In Proceedings of the IEEE International Conference on Neural Networks, Perth, Australia, pp. 1942–1948. https://doi.org/10.1109/ICNN.1995.488968

[32] Fang, J., Liu, W., Chen, L., Lauria, S., Miron, A., Liu, X. (2023). A survey of algorithms, applications and trends for particle swarm optimization. International Journal of Network Dynamics and Intelligence, 2(1): 24-50. https://doi.org/10.53941/ijndi0201002

[33] Zhan, Z.H., Zhang, J., Li, Y., Chung, H.S.H. (2009). Adaptive particle swarm optimization. IEEE Transactions on Systems, Man, and Cybernetics – Part B: Cybernetics, 39(6): 1362-1381. https://doi.org/10.1109/TSMCB.2009.2015956

[34] Liu, F. (2023). Restoration method of motion blurred image based on feature fusion and particle swarm optimization algorithm. In 2023 International Conference on Distributed Computing and Electrical Circuits and Electronics (ICDCECE), Ballar, India, pp. 1-5. https://doi.org/10.1109/ICDCECE57866.2023.10151124

[35] Gaing, Z.L. (2004). A particle swarm optimization approach for optimum design of PID controller in AVR system. IEEE Transactions on Energy Conversion, 19(2): 384-391. https://doi.org/10.1109/TEC.2003.821821

[36] Liu, L., Liu, W.X., Cartes, D.A. (2008). Particle swarm optimization-based parameter identification applied to permanent magnet synchronous motors. Engineering Applications of Artificial Intelligence, 21: 1092-1100. https://doi.org/10.1016/j.engappai.2007.10.002

[37] Song, B., Wang, Z., Zou, L. (2017). On global smooth path planning for mobile robots using a novel multimodal delayed PSO algorithm. Cognitive Computation, 9(1): 5-17. https://doi.org/10.1007/s12559-016-9442-4

[38] Shi, Y., Eberhart, R. (1998). A modified particle swarm optimizer. In 1998 IEEE International Conference on Evolutionary Computation Proceedings. IEEE World Congress on Computational Intelligence (Cat. No.98TH8360), Anchorage, USA, pp. 69-73. https://doi.org/10.1109/ICEC.1998.699146

[39] Bratton, D., Kennedy, J. (2007). Defining a standard for particle swarm optimization. In Proceedings of the IEEE Swarm Intelligence Symposium, Honolulu, pp. 120-127. https://doi.org/10.1109/SIS.2007.368035

[40] Jordehi, A.R. (2016). Time varying acceleration coefficients particle swarm optimisation (TVACPSO): A new optimisation algorithm for estimating parameters of PV cells and modules. Energy Conversion and Management, 129: 262-274. https://doi.org/10.1016/j.enconman.2016.09.085

[41] Liu, W., Wang, Z., Yuan, Y., Zeng, N., Hone, K., Liu, X. (2021). A novel sigmoid-function-based adaptive weighted particle swarm optimizer. IEEE Transactions on Cybernetics, 51(2): 1085-1093. https://doi.org/10.1109/TCYB.2019.2925015

[42] Chen, K., Zhou, F., Wang, S., Wang, Y., Wan, F., Yin, L. (2018). A hybrid particle swarm optimizer with sine cosine acceleration coefficients. Information Sciences, 422: 218-241. https://doi.org/10.1016/j.ins.2017.09.015

[43] Bansal, J.C., Singh, P.K., Saraswat, M., Verma, A., Jadon, S.S., Abraham, A. (2011). Inertia weight strategies in particle swarm optimization. In Proceedings of the 2011 Third World Congress on Nature and Biologically Inspired Computing (NaBIC), Salamanca, Spain, pp. 633-640. https://doi.org/10.1109/NaBIC.2011.6089659

[44] Wei, L., Li, X., Fan, R. (2019). A new multi-objective particle swarm optimisation algorithm based on R² indicator selection mechanism. International Journal of Systems Science, 50(10): 1920-1932. https://doi.org/10.1080/00207721.2019.1645914

[45] Xu, L., Song, B., Cao, M. (2021). An improved particle swarm optimization algorithm with adaptive weighted delay velocity. Systems Science and Control Engineering, 9(1): 188-197. https://doi.org/10.1080/21642583.2021.1891153

[46] Almeida, M.S.C., Figueiredo, M.A.T. (2013). Parameter estimation for blind and non-blind deblurring using residual whiteness measures. IEEE Transactions on Image Processing, 22(7): 2751-2763. https://doi.org/10.1109/TIP.2013.2257810

[47] Chen, Z., Basarab, A., Kouamé, D. (2015). A simulation study on the choice of regularization parameter in ℓ2-norm ultrasound image restoration. In 2015 37th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pp. 6346-6349. https://doi.org/10.1109/EMBC.2015.7319844

[48] Perry, S.W., Wong, H.S., Guan, L. (2002). Adaptive Image Processing: A Computational Intelligence Perspective. CRC Press. https://doi.org/10.1201/9781315218625

[49] Bredies, K., Kunisch, K., Pock, T. (2010). Total generalized variation. SIAM Journal on Imaging Sciences, 3(3): 492-526. https://doi.org/10.1137/090769521