# Image De-Quantization Using Plate Bending Model
## 4 Material and Methods
### Numerical Solution  
The non-continuous regions and the additional information is taken into account through the upper-lower limits, and the higher-order smoothness is guaranteed by the biharmonic equation.  
A straightforward solution can be implemented using finite differences schemes with a gradient descent method. However, a clipping step is required after every iteration step: if a value is lower than the lower envelope in the given point, or higher than the upper envelope, then approximation is set to be equal with the corresponding envelope. Laplace and biharmonic equations are part of the polyharmonic equations, and they can be described with similar equations where the difference is the order (\(p\)) of the Laplace-operator, which is 1 for Laplace and 2 for biharmonic operator. The gradient for time step is as follows:  
\[\frac{\mathrm{d}f(x,y)}{\mathrm{d}t}=\Delta^{p}f(x,y). \tag{7}\]  
While initially the solution function (\(f\)) does not satisfy the Laplace/biharmonic equations, evolving it with small time steps converges it to the solution. In this formulation, the solution of the problem is \(f(x,y,t)\) in \(t\rightarrow\inf\) limit.  
During the discretization of the partial differential equations, we assume that the pixel sizes are equal in both dimensions, and the lowest-order central difference schemes are satisfactory for the numerical solution. With these assumptions, we use the following stencils approximating the Laplace and biharmonic operators:  
\[\Delta\rightarrow\begin{array}{cccc}&0&1&0\\ 1&-4&1\\ &0&1&0\\ \end{array} \tag{8}\]  
\[\begin{array}{cccc}0&0&1&0&0\\ 0&2&-8&2&0\\ 1&-8&20&-8&1\\ 0&2&-8&2&0\\ 0&0&1&0&0\\ \end{array} \tag{9}\]  
These stencils can be used everywhere in the image, including the edges, due to the periodic boundary conditions.  
The convergence rate of this method is slow. Reducing execution time, the implementation uses a multi-grid approach where one level means a factor of 2 change in sampling points along an axis. First, both the envelopes and the initial image is downscaled. The image pixels take the average value of the higher-resolution pixels, the lower envelope takes the minimum of the higher-resolution envelope, and the upper envelope takes their maximum. This lower-resolution problem is solved first. When the approximation reaches convergence, a bicubic resampling is used to upscale the approximate solution, which serves as input for the next level. The pseudo-code of the iteration scheme is as follows:  
1, load initial image.
2, if the smallest edge is longer than 31 pixels and divisible with two:  
2a, downsample the problem with factor of 2:  
* lower-envelope: minimum of grid points  * upper-envelope: maximum of grid points  * image: average of grid points2b, solve the problem (jump to step 1)  2c, upscale the solution with a factor of two using bicubic interpolation  2d, use the upscaled image as initial image for the next step
3, calculate gradient for the time step
4, update the image
5, clip the image if it exceeds limits (envelopes)
6, continue with Step 3 if iteration number is less than specified
7, return the image to higher level or as final solution  
Alternatively, any other de-quantization method could be used as an initial input instead of Steps 2a-2d, and the algorithm should converge to the same solution.