daubcqf.m
上传用户:speoil
上传日期:2022-06-23
资源大小:224k
文件大小:4k
源码类别:

波变换

开发平台:

Matlab

  1. function [h_0,h_1] = daubcqf(N,TYPE)
  2. %    [h_0,h_1] = daubcqf(N,TYPE); 
  3. %
  4. %    Function computes the Daubechies' scaling and wavelet filters
  5. %    (normalized to sqrt(2)).
  6. %
  7. %    Input: 
  8. %       N    : Length of filter (must be even)
  9. %       TYPE : Optional parameter that distinguishes the minimum phase,
  10. %              maximum phase and mid-phase solutions ('min', 'max', or
  11. %              'mid'). If no argument is specified, the minimum phase
  12. %              solution is used.
  13. %
  14. %    Output: 
  15. %       h_0 : Minimal phase Daubechies' scaling filter 
  16. %       h_1 : Minimal phase Daubechies' wavelet filter 
  17. %
  18. %    Example:
  19. %       N = 4;
  20. %       TYPE = 'min';
  21. %       [h_0,h_1] = daubcqf(N,TYPE)
  22. %       h_0 = 0.4830 0.8365 0.2241 -0.1294
  23. %       h_1 = 0.1294 0.2241 -0.8365 0.4830
  24. %
  25. %    Reference: "Orthonormal Bases of Compactly Supported Wavelets",
  26. %                CPAM, Oct.89 
  27. %
  28. %File Name: daubcqf.m
  29. %Last Modification Date: 01/02/96 15:12:57
  30. %Current Version: daubcqf.m 2.4
  31. %File Creation Date: 10/10/88
  32. %Author: Ramesh Gopinath  <ramesh@dsp.rice.edu>
  33. %
  34. %Copyright (c) 2000 RICE UNIVERSITY. All rights reserved.
  35. %Created by Ramesh Gopinath, Department of ECE, Rice University. 
  36. %
  37. %This software is distributed and licensed to you on a non-exclusive 
  38. %basis, free-of-charge. Redistribution and use in source and binary forms, 
  39. %with or without modification, are permitted provided that the following 
  40. %conditions are met:
  41. %
  42. %1. Redistribution of source code must retain the above copyright notice, 
  43. %   this list of conditions and the following disclaimer.
  44. %2. Redistribution in binary form must reproduce the above copyright notice, 
  45. %   this list of conditions and the following disclaimer in the 
  46. %   documentation and/or other materials provided with the distribution.
  47. %3. All advertising materials mentioning features or use of this software 
  48. %   must display the following acknowledgment: This product includes 
  49. %   software developed by Rice University, Houston, Texas and its contributors.
  50. %4. Neither the name of the University nor the names of its contributors 
  51. %   may be used to endorse or promote products derived from this software 
  52. %   without specific prior written permission.
  53. %
  54. %THIS SOFTWARE IS PROVIDED BY WILLIAM MARSH RICE UNIVERSITY, HOUSTON, TEXAS, 
  55. %AND CONTRIBUTORS AS IS AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, 
  56. %BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS 
  57. %FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL RICE UNIVERSITY 
  58. %OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 
  59. %EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 
  60. %PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; 
  61. %OR BUSINESS INTERRUPTIONS) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, 
  62. %WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR 
  63. %OTHERWISE), PRODUCT LIABILITY, OR OTHERWISE ARISING IN ANY WAY OUT OF THE 
  64. %USE OF THIS SOFTWARE,  EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  65. %
  66. %For information on commercial licenses, contact Rice University's Office of 
  67. %Technology Transfer at techtran@rice.edu or (713) 348-6173
  68. if(nargin < 2),
  69.   TYPE = 'min';
  70. end;
  71. if(rem(N,2) ~= 0),
  72.   error('No Daubechies filter exists for ODD length');
  73. end;
  74. K = N/2;
  75. a = 1;
  76. p = 1;
  77. q = 1;
  78. h_0 = [1 1];
  79. for j  = 1:K-1,
  80.   a = -a * 0.25 * (j + K - 1)/j;
  81.   h_0 = [0 h_0] + [h_0 0];
  82.   p = [0 -p] + [p 0];
  83.   p = [0 -p] + [p 0];
  84.   q = [0 q 0] + a*p;
  85. end;
  86. q = sort(roots(q));
  87. qt = q(1:K-1);
  88. if TYPE=='mid',
  89.   if rem(K,2)==1,  
  90.     qt = q([1:4:N-2 2:4:N-2]);
  91.   else
  92.     qt = q([1 4:4:K-1 5:4:K-1 N-3:-4:K N-4:-4:K]);
  93.   end;
  94. end;
  95. h_0 = conv(h_0,real(poly(qt)));
  96. h_0 = sqrt(2)*h_0/sum(h_0);  %Normalize to sqrt(2);
  97. if(TYPE=='max'),
  98.   h_0 = fliplr(h_0);
  99. end;
  100. if(abs(sum(h_0 .^ 2))-1 > 1e-4) 
  101.   error('Numerically unstable for this value of "N".');
  102. end;
  103. h_1 = rot90(h_0,2);
  104. h_1(1:2:N)=-h_1(1:2:N);