如 l1_ls 等算法中需要 A 的 Transpose (应该是 conjugate transpose?)
Ax, A^T*y
若 A=RW, 则 A^T = R^T W^T
1. A:DWT, DCT 等,参见 dftmtx, dftmtx2 等
A^T 就是反变换(或者正变换,如果A是反变换)
A^T *A = I,这要求正变换反变换的系数一致,比如很多文献中(包括matlab 中),DFT 定义为
\sum_m \sum_n f(m,n) \exp ...
IDFT 定义为
1/(MN) \sum_u \sum_v F(u,v) \exp ...
这样的定义不满足A^T A = I, 而是 A^T A =MN I
我们只需采用另一种形式即可
1/sqrt(MN) \sum_m \sum_n f(m,n) \exp ...
1/sqrt(MN) \sum_u \sum_v F(u,v) \exp...
注意:在采用function handle 或者 class object.时,还是可以采用传统定义(如 matlab fft2)的 DFT,只要记住A^T 是相应的反变换,即可。
2. A*z: subset subsampling, 取 z 的一部分,令选中的下标集为 I
z=A^T*w: z=zeros(n,1); z(I)=w
3. A: space-invariant circular (periodic) convolutions
参见论文 An EM Algorithm for Wavelet-Based Image Restoration
A 可以分解为
A=U^H D U
In the above equation, U is the matrix that represents the 2-D discrete Fourier transform, U^H is its inverse (U is an orthogonal matrix, that is, U^H U = I, where H denotes conjugate transpose), and D is a diagonal matrix containing the DFT coefficients of the convolution operator represented by A. This means that multiplication by A can be performed in the discrete Fourier domain with a simple point-wise multiplication (recall that D is diagonal).
Ax = U^H D Ux = U^H D \tilde{x}
\tilde{x}is the DFT of .
thus A^H z = U^H D^H \tilde{z}
If matrix A is block-Toeplitz, but not block-circulant, it is possible to embed the nonperiodic convolution that it represents in a larger periodic convolution and still work in the DFT domain [16]. Accordingly, all the results and statements made in this paper concerning circulant observation matrices (periodic convolutions) can be extended to the Toeplitz case.
% % center and normalize the blur h = fftshift(h); h = h/sum(h(:)); % definde the function handles that compute % the blur and the conjugate blur. R = @(x) real(ifft2(fft2(h).*fft2(x))); RT = @(x) real(ifft2(conj(fft2(h)).*fft2(x)));
参见 matlab dftmtx, convmtx (convolution matrix), convmtx2 (2D convolution matrix)
以及我写的 dftmtx2 (2D dft transform matrix)
function A = dftmtx2(M,N) % 2D DFT transform matrix % A = dftmtx2(M,N) % y = Ax equal to % Y = fft2(X), x = X'(:),y = Y'(:) % row wise d = M*N; A = zeros(d,d); [X Y]=meshgrid(0:N-1,0:M-1); Y = Y/M; X = X/N; for i=0:M-1 for j=0:N-1 curr = exp(-2i*pi*(Y.*i + X.*j)); curr = curr.'; A(i*N+j+1,:)=curr(:); end endA^T * A = cI (并非等于 I)
不过A^T 确实起到了 Inverse transform 的作用
0 comments:
Post a Comment