Source '''Levin, Anat, et al. "Deconvolution using natural image priors." Massachusetts Institute of Technology, Computer Science and Artificial Intelligence Laboratory 3 (2007). APA '''
Generally 2-D convolution can be written as $$y = f \circledast x $$, where $x$ is an image and $f$ represents a kernel. Since convolution is linear operation, it can be written as $$y = C_f x $$ Using a set of priors , the equation which we have to solve is $$A x = b$$ where $$A = C_f^T C_f + w \sum_k C_{gk}^T C_{gk}$$ and $$b = C_f^T y$$ While solving the equation using Conjugate gradient Algorithm, we choose an intial guess for $x (= x_0)$ and its implementation in the MATLAB is written as ( for $Ax = C_f^T C_f x_0$)
$\textbf{Ax=conv2(conv2(x,fliplr(flipud(filt1)),'same').*mask, filt1,'same'); }$
$\textbf{Question:}$ It is okay that $C_f x$ can be represented as 2d Convolution but how can
$ C_f^T$ (times a vector)
can be represented as a 2-d convolution? (where $C_f$ is a convolution matrix)
Conjugate Gradient Method is just a method to solve large scale System of Linear Equations:
$$ A x = b $$
In most cases, for really large system, the matrix $ A $ is sparse so it can use Sparse Linear Algebra.
The nice thing is in many times the Linear Operator can be represented by a much simpler (To calculate) operation. In Signal / Image Processing it happens to be, in most cases, convolution.
The Conjugate Gradient Method requires also applying the Adjoint / Transpose of the operator. In the case of convolution it is the correlation operation.
For 1D kernels it means to apply convolution with the flipped signal (Reversed order).
For 2D kernels it means to flip them bot horizontally and vertically. This is why in the code you see both
fliplr()andflipud().Simpler way to do it is
fliplr(flipud(mA)) = mA(end:-1:1, end:-1:1).Chapter 5 of Jon F. Claerbout - Earth Soundings Analysis: Processing versus Inversion has a great discussion of the topic of Adjoint operators (Also known as transpose, back projection or matched filtering).