频谱细化之Chirp-Z变换算法

---

一、 算法简介

Chirp-Z变换的核心思想是在 Z 平面的单位圆上沿螺旋线进行等采样，若 已知离散采样序列 x(m)，那么该信号的 Z 变换可表示为：

Chirp-Z 变换算法在 Z 平面的一段螺旋线进行等角采样的示意图如下：

其中，Z0,Z1,…,ZM-1表示螺旋线上等分角的 M 个采样点；A0表示螺旋线上起始点 的半径 A0<1 表示起始采样点在圆内；θ0 表示沿着逆时针的起始采样点相角，可 正可负；φ0表示相邻采样点的角频率差，φ0<0 表示 Zk沿着顺时针顺序采样，φ0<0 表示 Zk沿着逆时针顺序采样
等角采样点 Zk可以表示为：



因此求 X(Zk)的过程转变成求 g(m)和 h(m)的卷积，即 g(m)*h(m)，线性调频 Z 变换的算法流程可以总结如下： (1)选择 FFT 点数 L,要求 L≥M1+M2-1，并且 L=2m； (2)将序列 g(m)扩展为 L 点

(3)计算 g(m)的 L 点离散傅里叶变化 G(k)=FFT[g(m)]；
(4)序列 h(m)选取 L 点；

二、算法仿真

clear;
close all;
clc;
fs = 1024; % 采样率
N = 1024; % 信号长度
nfft = N;
n = 0:1:N-1;
n1 = fs * (0:nfft/2-1) / nfft; % 时间
x = 1.5 * cos(2 * pi * 98 * n / fs) + 2 * cos(2 * pi * 99 * n / fs) + 3 * cos(2 * pi * 100. * n / fs) + 3.5 * cos(2 * pi * 101 * n / fs);
range_win = hamming(N);   %海明窗
% range_win = hanning(N);   %汉宁窗
x=x.*range_win';
figure(1);
plot(n, x, 'LineWidth', 2, 'Color', 'b');
ax1 = gca;  % 获取当前子图的坐标轴对象
ax1.FontSize = 25;  % 设置坐标轴上文字的字体大小
xlabel('采样点数','FontSize',30);
ylabel('幅值','FontSize',30);
title('原始信号','FontSize',30);

% 信号的FFT
XK = fft(x, nfft);
figure(2);
subplot(211);
plot(n1, 2 * abs(XK(1:(N/2))) / N);
grid on;
title('信号的FFT');

% CZT频谱
f1 = 93; % 细化起始范围
f2 = 106; % 细化结束范围
M = 200; % 细化倍数
w = exp(-1j * 2 * pi * (f2 - f1) / (fs * M));
a = exp(1j * 2 * pi * f1 / fs);
xk = czt(x, M, w, a);
h = 0:1:M-1;
f0 = (f2 - f1) / M * h + f1;
subplot(212);
plot(f0, 2 * abs(xk) / N);
xlabel('f');
ylabel('value');
title('CZT频谱细化后');

% 对比FFT和CZT频谱
figure(3);
subplot(211);
plot(n1, 2 * abs(XK(1:(N/2))) / N, 'b', 'LineWidth', 3);
grid on;
title('FFT和CZT仿真对比','FontSize',30);
ax1 = gca;  % 获取当前子图的坐标轴对象
ax1.FontSize = 25;  % 设置坐标轴上文字的字体大小
xlabel('Frequency (Hz)','FontSize',30);
ylabel('Amplitude','FontSize',30);
legend('FFT','FontSize',25);

subplot(212);
plot(f0, 2 * abs(xk) / N, 'r', 'LineWidth', 3);
ax1 = gca;  % 获取当前子图的坐标轴对象
ax1.FontSize = 25;  % 设置坐标轴上文字的字体大小
grid on;
% title('CZT频谱细化后');
xlabel('Frequency (Hz)','FontSize',30);
ylabel('Amplitude','FontSize',30);
legend('CZT','FontSize',25);

三、仿真结果

---