firrcos.m
上传用户:loeagle
上传日期:2013-03-02
资源大小:1236k
文件大小:7k
源码类别:

通讯编程文档

开发平台:

Matlab

  1. function [b,a]=firrcos(varargin)
  2. %FIRRCOS Raised Cosine FIR Filter design.
  3. %   B=FIRRCOS(N,Fc,DF,Fs) returns an order N low pass linear phase FIR 
  4. %   filter with a raised cosine transition band.  The filter has cutoff
  5. %   frequency Fc, sampling frequency Fs and transition bandwidth DF 
  6. %   (all in Hz).
  7. %
  8. %   Fc +/- DF/2 must be in the range [0,Fs/2].    
  9. %
  10. %   The coefficients of B are normalized so that the nominal passband 
  11. %   gain is always equal to one.
  12. %
  13. %   FIRRCOS(N,Fc,DF) uses a default sampling frequency of Fs = 2.
  14. %
  15. %   B=FIRRCOS(N,Fc,R,Fs,'rolloff') interprets the third argument as the
  16. %   rolloff factor instead of as a transition bandwidth. Alternatively,
  17. %   you can specify B=FIRRCOS(N,Fc,DF,Fs,'bandwidth') which is 
  18. %   equivalent to B=FIRRCOS(N,Fc,DF,Fs).
  19. %
  20. %   R must be in the range [0,1].
  21. %
  22. %   B=FIRRCOS(N,Fc,DF,Fs,TYPE) or B=FIRRCOS(N,Fc,R,Fs,'rolloff',TYPE) 
  23. %   will design a regular FIR raised cosine filter when TYPE is 
  24. %   'normal' or set to an empty matrix. If TYPE is 'sqrt', B is the 
  25. %   square root FIR raised cosine filter.
  26. %
  27. %   B=FIRRCOS(...,TYPE,DELAY) allows for a variable integer delay to be 
  28. %   specified. When omitted or left empty, DELAY defaults to N/2 or
  29. %   (N+1)/2 depending on whether N is even or odd.
  30. %
  31. %   DELAY must be an integer in the range [0, N+1].
  32. %
  33. %   B=FIRRCOS(...,DELAY,WINDOW) applies a length N+1 window to the 
  34. %   designed filter in order to reduce the ripple in the frequency 
  35. %   response. WINDOW must be a N+1 long column vector. If no window
  36. %   is specified a boxcar (rectangular) window is used.
  37. %
  38. %   WARNING: Care must be exercised when using a window with a delay
  39. %   other than the default.
  40. %
  41. %   [B,A]=FIRRCOS(...) will always return A = 1.
  42. %
  43. %   See also FIRLS, FIR1, FIR2.
  44. %   Author(s): R. Losada and D. Orofino
  45. %   Copyright 1988-2001 The MathWorks, Inc.
  46. %   $Revision: 1.11 $  $Date: 2001/04/02 20:21:59 $
  47. error(nargchk(3,8,nargin));
  48. [n,fc,fs,R,designType,window,msg] = parse_inputs(varargin{:});
  49. error(msg);
  50. switch designType
  51. case 'normal'   %normal raised cosine design
  52.    b = normal_design(n,fc,fs,R);
  53. case 'sqrt' % square root raised cosine design
  54. b = sqrt_design(n,fc,fs,R);
  55. end
  56. if ~isempty(window),
  57. [b,msg] = apply_win(b,window);
  58. error(msg);
  59. end
  60. if nargout > 1
  61.    a = 1.0;
  62. end
  63. %-------------------------------------------------------------------------------
  64. function b = normal_design(n,fc,fs,R)
  65. ind1 = find(abs(abs(4.*R.*fc.*n) - 1.0) > sqrt(eps));
  66. if ~isempty(ind1),
  67. nind = n(ind1);
  68. b(ind1) =  sinc(2.*fc.*nind)./fs   ...
  69. .* cos(2.*pi.*R.*fc.*nind) ...
  70. ./ (1.0 - (4.*R.*fc.*nind).^2);
  71. end
  72. ind = 1:length(n);
  73. ind(ind1) = [];
  74. b(ind) = R ./ (2.*fs) .* sin(pi ./ (2.*R));
  75. b = 2.*fc.*b;
  76. %-------------------------------------------------------------------------------
  77. function b = sqrt_design(n,fc,fs,R)
  78. ind1 = find(n == 0);
  79. if ~isempty(ind1),
  80. b(ind1) = - sqrt(2.*fc) ./ (pi.*fs) .* (pi.*(R-1) - 4.*R );
  81. end
  82. ind2 = find(abs(abs(8.*R.*fc.*n) - 1.0) < sqrt(eps));
  83. if ~isempty(ind2),
  84. b(ind2) = sqrt(2.*fc) ./ (2.*pi.*fs) ...
  85. * (    pi.*(R+1)  .* sin(pi.*(R+1)./(4.*R)) ...
  86. - 4.*R     .* sin(pi.*(R-1)./(4.*R)) ...
  87. + pi.*(R-1)  .* cos(pi.*(R-1)./(4.*R)) ...
  88. );
  89. end
  90. ind = 1:length(n);
  91. ind([ind1 ind2]) = [];
  92. nind = n(ind);
  93. b(ind) = -4.*R./fs .* ( cos((1+R).*2.*pi.*fc.*nind) + ...
  94. sin((1-R).*2.*pi.*fc.*nind) ./ (8.*R.*fc.*nind) ) ...
  95. ./ (pi .* sqrt(1./(2.*fc)) .* ((8.*R.*fc.*nind).^2 - 1));
  96. b = sqrt(2.*fc) .* b;
  97. %-------------------------------------------------------------------------------
  98. function [b,msg] = apply_win(b,window)
  99. msg = '';
  100. if length(window) ~= length(b),
  101. msg = 'WINDOW must be of the same length as the filter.';
  102. return
  103. else
  104. b = b .* window(:).';
  105. end
  106. %-------------------------------------------------------------------------------
  107. function [n,fc,fs,R,designType,window,msg] = parse_inputs(varargin)
  108. % Initialize in case of early return
  109. n = [];
  110. fc = [];
  111. fs = [];
  112. R = [];
  113. designType = '';
  114. window = [];
  115. msg = '';
  116. N = varargin{1};
  117. if isempty(N) | round(N) ~= N | N < 0,
  118.    msg = 'Order must be a positive integer.';
  119.    return
  120. end
  121. L = N+1; % Length of window
  122. fc = varargin{2};
  123. R = varargin{3};  % DF or R
  124. % If optional arguments are not passed, substitute with empty:
  125. for i = nargin+1:8,
  126.    varargin{i}=[];
  127. end
  128. arg5opts = {'rolloff','sqrt','normal','bandwidth'};
  129. % map 5th arg to one of 4 possible choices:
  130. if isempty(varargin{5}),
  131.    varargin{5} = arg5opts{3};
  132. else
  133.    idx = strmatch(lower(varargin{5}), arg5opts);
  134.    if isempty(idx),
  135.       msg = 'Argument 5 is unknown - must be one of: rolloff, bandwidth, sqrt, or normal.';
  136.   return
  137.    end
  138.    varargin{5} = arg5opts{idx};
  139. end
  140. % Apply defaults as appropriate:
  141. %
  142. % Set up default values
  143. fs = 2;
  144. designType = arg5opts{3};
  145. if rem(L,2),
  146.    delay = (L-1)/2;
  147. else
  148.    delay = L/2;
  149. end
  150. % Setup arg translation:
  151. params = {'fs','designType','delay','window'};
  152. % We define a flag to indicate whether a string for the transition region type was specified
  153. isTranRegionStr = strcmp(varargin{5},'rolloff') | strcmp(varargin{5},'bandwidth');
  154. if isTranRegionStr,
  155.    xlat = [4 6:8];
  156. else
  157.    xlat = 4:7;
  158. end
  159. % Override defaults when needed:
  160. for i=1:length(xlat),
  161.    arg = varargin{xlat(i)};
  162.    if ~isempty(arg),
  163.      eval([params{i} '=arg;']);
  164.    end
  165. end
  166. % Check for validity of fs
  167. if ischar(fs),
  168.    msg = 'Fs must be a number';
  169.    return
  170. end
  171. % Check for valid cutoff frequency
  172. if (fc <= 0) | (fc >= fs./2),
  173.    msg = 'The cutoff frequency, Fc, must satisfy 0 < Fc < Fs/2.';
  174.    return
  175. end
  176. % Check for valid rolloff or bandwidth values
  177. if strcmp(varargin{5},'rolloff'),
  178.    % check if input arguments are valid 
  179.    if R < 0 | R > 1,
  180.       msg = 'The rolloff factor, R, must satisfy 0 <= R <= 1.';
  181.   return
  182.    end
  183.    % check for range of input arguments
  184.    if (fc + R.*fc) > fs/2
  185.       msg = sprintf(['The cutoff frequency, Fc, and rolloff factor, R,n',...
  186.   'must be specified such that Fc + Fc*R <= Fs/2.']);
  187.   return
  188.    end
  189. elseif strcmp(varargin{5},'bandwidth') | ~isTranRegionStr % arg5 is bandwidth, sqrt or normal
  190.    % check for range of input arguments
  191.    if fc - R/2 < 0 | fc + R/2 > fs/2
  192.       msg = sprintf(['The cutoff frequency, Fc, and the transition bandwidth, DF,n',...
  193.   'must be specified such that Fc +/- DF/2 is between zero and Fs/2.']);
  194.   return
  195.    end
  196.    % bandwidth is valid, convert to rolloff
  197.    R = R / (2*fc);
  198. end
  199. if delay < 0 | delay > L
  200.    msg = 'DELAY must be in the range [0, L+1].';
  201.    return
  202. elseif round(delay) ~= delay
  203.    msg = 'DELAY must be an integer.';
  204.    return
  205. end
  206. % R is now always a rolloff factor - DF has been converted
  207. if R == 0,
  208.    R = realmin;
  209. end
  210. %n = -delay/fs : 1/fs : (L-delay-1)/fs;
  211. n = ((0:L-1)-delay) ./ fs;
  212. if isTranRegionStr, % 6th argument, if present, is designType
  213.    arg6opts = {'sqrt','normal'};
  214.    % map 6th arg to one of 2 possible choices:
  215.    if isempty(varargin{6}),
  216.       designType = arg6opts{2};
  217.    else
  218.       idx = strmatch(lower(varargin{6}), arg6opts);
  219.       if isempty(idx),
  220.          msg = 'Argument 6 is unknown - must be one of: sqrt, normal or [].';
  221.  return
  222.       end
  223.       designType = arg6opts{idx};
  224.    end
  225. end
  226. % EOF