4
$\begingroup$

I have two 2D arrays: $A$ is the original matrix that contains only 0s and 1s, and $B$ is the convolved matrix. I know the size of the convolution kernel $K$. Generally, it follows $B = A*K + S$, where $S$ is the additional noise with Gaussian distribution.

I want to estimate the kernel based on $A$ and $B$, how can I do it in Python?

Here is my example code:

import numpy as np
from scipy.signal import convolve2d,oaconvolve,fftconvolve

# the dimension of the original matrix
N = 100
# create an example original matrix with around 20% 1s
src_mat = np.random.binomial(n=1, p=0.2, size=(N, N))

# a kernel function with a dimension of 2*N-1, the values decay outwards from center of the matrix, as a function of distance (here 1/r)
kernel = np.fromfunction(lambda x,y:np.sqrt((x-(N-1))**2+(y-(N-1))**2)**(-1) ,(2*N -1,2*N-1))
kernel[N-1,N-1] = 0

# the convolved matrix in the form of B=A*K+S
tgt_mat = oaconvolve(src_mat,kernel,mode='same') + np.random.normal(0,1,src_mat.shape)

```
$\endgroup$
12
  • 1
    $\begingroup$ Does this answer your question? Estimating Convolution Kernel from Input and Output Images $\endgroup$
    – Gillespie
    Aug 25, 2022 at 12:38
  • $\begingroup$ Not really, I don't understand some of the derivations in the accepted answer, I also don't have any knowledge in Matlab. I reached Royi and we had a bit discussion, then I decided to open a new question. $\endgroup$ Aug 25, 2022 at 12:53
  • $\begingroup$ I will have a look at it. But for sure don't go as is to the frequency domain :-). $\endgroup$
    – Royi
    Aug 25, 2022 at 13:41
  • $\begingroup$ @Royi, many thanks, hope this won't take too much effort. I also prefer a similar solution as you put there. $\endgroup$ Aug 25, 2022 at 14:16
  • 1
    $\begingroup$ @Royi, thanks for your help, a friend explained a bit what your matlab code is doing, I now understand your answer in another question. I just implemented it in Python and it works quite very when the kernel is smaller than the convolution matrix. For this question, I still see the potential of solving it in a similar way as you provided. I will continue exploring it. $\endgroup$ Sep 1, 2022 at 22:47

1 Answer 1

1
$\begingroup$

The issue is with the kernel being too big and creating a case that not only we're trying to solve ill poised problem, we also have more parameters to estimates than measurements.

One way to handle this would be a lower rank approximation for the kernel. By using the approach in Estimating Convolution Kernel from Input and Output Images one can chose a big support for the kernel which still have lower number of parameters than measurements. Yet the approximation is $ {L}_{2} $ based.

One might try other regularizations on the kernel which can be solved using iterative methods such as ADMM or Accelerated Proximal Gradient Descent.

One could even dare to try solving the undetermined system but probably the results won't be satisfactory.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.