We all know that Discrete Fourier Transform (DFT) corresponds to circular (not linear) convolution. That is to say, if $x(n),h(n)$ and $y(n)$ is the original signal, the filter and output signal in domain and $X(\omega),H(\omega),Y(\omega)$ are their counterpart in Fourier domain, then if we want $Y(\omega) = H(\omega)X(\omega)$ we must have $y(n) = h(n) \circledast x(n)$, where $\circledast$ is circular convolution (modulo the length of x).

If $h$ is the linear operator that stands for circular convolution with $h(n)$, then I think $h^T$, the transpose of $h$, is a circular convolution plus a circular shift. My question is, is $h^Th$ also a circular convolution of some kernel?

The motivation for thinking about this is to look for a fast method to compute $(I+h^Th)^{-1}$ through FFT. I know it seems common to introduce FFT to compute this, but how exactly when h is a convolution kernel?

Thanks in advance, Bob

  • $\begingroup$ Sorry it seems a typo, the first line should read: hi, all $\endgroup$ – Yinan Hu Apr 15 at 17:55

If $ H $ is a matrix form of Circular Convolution then it is a Circulant Matrix.
Being a Circulant Matrix means it can be diagonalized by the Fourier Matrix $ {F} $:

$$ H = {F}^{H} D F $$

Where the matrix $ D $ id a Diagonal Matrix with the Fourier Coefficients of $ \mathcal{F} \left( h \right) $ on it main diagonal.
Also pay attention that we use the Conjugate Transpose $ {\left( \cdot \right)}^{H} $ operator as $ F $ and $ D $ are built by complex numbers (The space is over the complex numbers, hence the Adjoint Operator is the Complex Conjugate).

Now, it is easy to see that:

$$ {H}^{H} = {F}^{H} {D}^{H} F $$

Which is again an operator diagonlaized with the Fourier Matrix hence it is a Circulant Matrix -> Cicrcular Operation.
Basically it is a multiplication in Fourier Domain (DFT) by the conjugate of the DFT of the vector $ h $. By the DFT properties it means it is the flipped version of $ h $.

By the way:

$$ {H}^{H} H = {F}^{H} {D}^{H} F {F}^{H} D F = {F}^{H} {D}^{H} D F $$

Now, pay attention to $ {D}^{H} D $. It is a diagonal matrix multiplied by its conjugate. Namely it is, in the Fourier Domain, the Squared Magnitude.
In time domain it means it can be done by circular convolution with the correlation function of $ h $ (Or just do circular convolution in one direction and then in the other which is basically convolution with the flipped version).

  • $\begingroup$ Thanks! That solves everything. $\endgroup$ – Yinan Hu yesterday
  • $\begingroup$ Could the one who -1 this answer explain why? $\endgroup$ – Royi yesterday

If A is a "full" linear convolution operation (|signal|+|kernel|-1), A^T is a "valid" cross-correlation (shorter output, |signal| - |kernel| + 1). So I suppose

I suppose your h^T is a circular cross-correlation.

I am guessing you represented your convolution kernel h to a matrix H with Toeplitz structure and you are solving for g = (I+H^T H)^(-1) f thru iterative methods like conjugate-gradient which involves a lot of matrix-vector calculation. In particular you need to be able to calculate this H^T H v fast for any vector v ? Then you can exploit circular Toeplitz structure for fast matrix-vector calculation. Wikipedia has some reference on that.

  • $\begingroup$ It would be helpful if you could elaborate cross-correlation. Yes I'm trying to compute the inverse of (I + H^TH) by solving linear system. It's true that normally we compute through solvers like conjugate-gradient, arnold, Krylov subspace, etc. But this time things are a bit different, I'm trying to solve it through point-wise multiplication. It's not iterative method and you need H to be shift-invariant. $\endgroup$ – Yinan Hu yesterday

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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