dt_BivaShrink123.m
上传用户:sla11nk8
上传日期:2013-03-09
资源大小:21k
文件大小:3k
源码类别:

其他

开发平台:

Matlab

  1. function y = dt_BivaShrink123(x)
  2. % Local Adaptive Image Denoising Algorithm
  3. % Usage :%        y = denoising_dtdwt(x)
  4. % INPUT :%        x - a noisy image
  5. % OUTPUT :%        y - the corresponding denoised image
  6. % Set the windowsize and the corresponding filter
  7. windowsize  = 7;
  8. windowfilt = ones(1,windowsize)/windowsize;
  9. windowfilt2=[1 1 1;1 0 1;1 1 1]/8;
  10. % Number of Stages
  11. J = 6;
  12. I=sqrt(-1);
  13. % symmetric extension
  14. L = length(x); % length of the original image.
  15. N = L+2^J;     % length after extension.
  16. x = symextend(x,2^(J-1));
  17. % Forward dual-tree DWT
  18. % Either FSfarras or AntonB function can be used to compute the stage 1 filters  
  19. [Faf, Fsf] = FSfarras;
  20. %[Faf, Fsf] = AntonB;
  21. [af, sf] = dualfilt1;
  22. W = cplxdual2D(x, J, Faf, af);
  23. %W = normcoef(W,J,nor);%-------------------------
  24. % Noise variance estimation using robust median estimator..
  25. tmp = W{1}{1}{1}{1};
  26. Nsig = median(abs(tmp(:)))/0.6745;
  27. for scale = 1:J-1
  28.     for dir = 1:2
  29.         for dir1 = 1:3            
  30.             % Noisy complex coefficients
  31.             %Real part
  32.             Y_coef_real = W{scale}{1}{dir}{dir1};
  33.             % imaginary part
  34.             Y_coef_imag = W{scale}{2}{dir}{dir1};
  35.             % The corresponding noisy parent coefficients
  36.             %Real part
  37.             Y_parent_real = W{scale+1}{1}{dir}{dir1};
  38.             % imaginary part
  39.             Y_parent_imag = W{scale+1}{2}{dir}{dir1};
  40.             % Extend noisy parent matrix to make the matrix size the same as the coefficient matrix.
  41.             Y_parent_real  = expand(Y_parent_real);
  42.             Y_parent_imag   = expand(Y_parent_imag); 
  43.             
  44.             %----计算邻域系数值----
  45.             Y_adjacent_real=sqrt(conv2((Y_coef_real).^2,windowfilt2,'same'));
  46.             Y_adjacent_imag=sqrt(conv2((Y_coef_imag).^2,windowfilt2,'same'));
  47.                         
  48.             % Signal variance estimation
  49.             Wsig = conv2(windowfilt,windowfilt,(Y_coef_real).^2,'same');%---用Gauss分布的ML估计系数方差----
  50.             Ssig = sqrt(max(Wsig-Nsig.^2,eps));          
  51.             % Threshold value estimation
  52.             T = sqrt(3)*Nsig^2./Ssig;            
  53.             % 
  54.             Y_coef = Y_coef_real+I*Y_coef_imag;
  55.             Y_parent = Y_parent_real + I*Y_parent_imag;
  56.             Y_adjacent = Y_adjacent_real + I*Y_adjacent_imag;            
  57.             % trivariate Shrinkage            
  58.             Y_coef = trishrink1(Y_coef,Y_adjacent,Y_parent,T);
  59.             W{scale}{1}{dir}{dir1} = real(Y_coef);
  60.             W{scale}{2}{dir}{dir1} = imag(Y_coef);
  61.         end
  62.     end
  63. end
  64. % Inverse Transform
  65. %W = unnormcoef(W,J,nor);------------------
  66. y = icplxdual2D(W, J, Fsf, sf);
  67. % Extract the image
  68. ind = 2^(J-1)+1:2^(J-1)+L;
  69. y = y(ind,ind);