在通信中，我们常用科斯塔斯环来提取相干载波，请给出相应的matlab代码

科斯塔斯环是一种常用的相干载波提取方法，其matlab代码如下：

fs = 100e3; % 采样频率
fc = 10e3; % 载波频率
N = 1000; % 采样点数
t = (0:N-1)/fs; % 时间数组

fIF = 1e3; % 中频频率
theta = pi/4; % 初始相位
K = 10; % 环宽
B = K/2; % 带宽

% 生成QPSK信号
data = randi([0,3],1,N/4);
modData = pskmod(data,4,pi/4);
txSig = real(modData.*exp(1j*2*pi*fc*t));

% 生成本振信号
osc = exp(1j*2*pi*(fc-fIF)*t);

% 理想低通滤波器
b = fir1(100,2*B/fs);
lpf = dsp.FIRFilter('Numerator',b);

% 科斯塔斯环
for k = 1:N
    mixSig = txSig(k)*osc(k); % 混频
    I(k) = real(mixSig); % 采样
    Q(k) = imag(mixSig);
    I_lp = step(lpf,I(k)); % 低通滤波
    Q_lp = step(lpf,Q(k));
    phase(k) = atan2(Q_lp,I_lp); % 相位提取
    osc(k+1) = exp(1j*(2*pi*fc*t(k+1) + phase(k))); % 更新本振
end

figure;
subplot(2,1,1);
plot(t,txSig);
title('发送信号');
xlabel('时间');
ylabel('幅度');
subplot(2,1,2);
plot(t,phase);
title('相位提取结果');
xlabel('时间');
ylabel('相位');

这段代码中，首先生成一个QPSK调制信号，然后通过混频的方式将其与本振信号相乘，得到带通信号。接着利用理想低通滤波器将带通信号滤波得到基带信号，再通过相位提取得到相干载波的相位。最后更新本振信号，进行下一次混频。