[586fe9]: fourier / chirpzt.m Maximize Restore History

Download this file

chirpzt.m    138 lines (124 with data), 3.6 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
function c = chirpzt(f,K,fdiff,foff,fs,dim)
%CHIRPZT Chirped Z-transform
% Usage: c = chirpzt(f,K,fdiff)
% c = chirpzt(f,K,fdiff,foff)
% c = chirpzt(f,K,fdiff,foff,fs)
%
% Input parameters:
% f : Input data.
% K : Number of values.
% fdiff : Frequency increment.
% foff : Starting frequency.
% fs : Sampling frequency.
%
% Output parameters:
% c : Coefficient vector.
%
% `c = chirpzt(f,K,fdiff,foff)` computes *K* samples of the discrete-time
% fourier transform DTFT *c* of *f* at values $c(k+1)=F(2\pi(f_{off}+kf_{diff}))$
% for $k=0,\dots,K-1$ where $F=DTFT(f)$. Values `foff` and `fdiff` should
% be in range of $0-1$. If `foff` is ommited or empty, it is considered to
% be 0. If `fdiff` is ommited or empty, *K* equidistant values
% $c(k+1)=F(2\pi k/K)$ are computed. If even *K* is ommited or empty,
% input length is used instead resulting in the same values as `fft` does.
%
% `c = chirpzt(f,K,fdiff,foff,fs)` computes coefficients using frequency
% values relative to *fs* $c(k+1)=F(2\pi(f_{off}+kf_{diff})/fs)$ for $k=0,\dots,K-1$.
%
% The input *f* is processed along the first non-singleton dimension or
% along dimension *dim* if specified.
%
% Examples:
% ---------
%
% Calculating DTFT samples of interest (aka zoom FFT):::
%
% % Generate input signal
% fs = 8000;
% L = 2^10;
% k = (0:L-1).';
% f1 = 400;
% f2 = 825;
% f = 5*sin(2*pi*k*f1/fs + pi/4) + 2*sin(2*pi*k*f2/fs - pi/3);
%
% % This is equal to fft(f)
% ck = chirpzt(f,L);
%
% %chirpzt to FFT error:
% norm(ck-fft(f))
%
% % Frequency "resolution" in Hz
% fdiff = 0.4;
% % Frequency offset in Hz
% foff = 803.9;
% % Number of frequency values
% K = 125;
% % DTFT samples. The frequency range of interest is 803.9-853.5 Hz
% ckchzt = chirpzt(f,K,fdiff,foff,fs);
%
% % Plot modulus of coefficients
% figure(1);
% fax=foff+fdiff.*(0:K-1);
% hold on;
% stem(k/L*fs,abs(ck),'k');
% stem(fax,abs(ckchzt),'r:');
% set(gca,'XLim',[foff,foff+K*fdiff]);
% set(gca,'YLim',[0 1065]);
% xlabel('f[Hz]');
% ylabel('|ck|');
%
% % Plot phase of coefficients
% figure(2);
% hold on;
% stem(k/L*fs,angle(ck),'k');
% stem(fax,angle(ckchzt),'r:');
% set(gca,'XLim',[foff,foff+K*fdiff]);
% set(gca,'YLim',[-pi pi]);
% xlabel('f[Hz]');
% ylabel('angle(ck)');
%
% See also: gga
%
% References: raschra69
%% Check the input arguments
if nargin < 1
error('%s: Not enough input arguments.',upper(mfilename))
end
if isempty(f)
error('%s: X must be a nonempty vector or a matrix.',upper(mfilename))
end
if nargin<6
dim=[];
end;
[f,~,Ls,~,dim,permutedsize,order]=assert_sigreshape_pre(f,[],dim,'CHIRPZT');
if nargin > 1 && ~isempty(K)
if ~isreal(K) || ~isscalar(K) || rem(K,1)~=0
error('%s: K must be a real integer.',upper(mfilename))
end
else
K = size(f,1);
end
if nargin > 2 && ~isempty(fdiff)
if ~isreal(K) || ~isscalar(K)
error('%s: fdiff must be a real scalar.',upper(mfilename))
end
else
fdiff = 1/K;
end
if nargin > 3 && ~isempty(foff)
if ~isreal(K) || ~isscalar(K)
error('%s: foff must be a real scalar.',upper(mfilename))
end
else
foff = 0;
end
if nargin > 4 && ~isempty(fs)
if ~isreal(fs) || ~isscalar(fs)
error('%s: fs must be a real scalar.',upper(mfilename))
end
else
fs = 1;
end
c = comp_chirpzt(f,K,2*pi*fdiff/fs,2*pi*foff/fs);
permutedsize(1)=K;
c=assert_sigreshape_post(c,dim,permutedsize,order);