出一道MATLAB信号与系统进阶题并给出代码

题目：

设计一个FIR滤波器，将输入信号x(n)中频率为1000Hz的正弦波滤除，输出y(n)。

输入信号：x(n) = 2*sin(2*pi*1000*n/Fs) + 0.5*cos(2*pi*2000*n/Fs)，采样频率为Fs=8000Hz。

滤波器长度为L=31，通带截止频率为800Hz，阻带截止频率为1200Hz，通带最大波纹为0.5dB，阻带最小衰减为40dB。

要求：画出滤波器的幅频响应和相频响应曲线，以及输入信号和输出信号的时域波形图。

代码：

% 输入信号
Fs = 8000;
n = 0:1/Fs:0.02;
x = 2*sin(2*pi*1000*n) + 0.5*cos(2*pi*2000*n);

% FIR滤波器
L = 31;
fc1 = 800;
fc2 = 1200;
wp = fc1/(Fs/2);
ws = fc2/(Fs/2);
Rp = 0.5;
Rs = 40;
b = fir1(L-1,wp,'low',kaiser(L,Rp));
[H,F] = freqz(b,1,512,Fs);

% 滤波后输出信号
y = filter(b,1,x);

% 画出滤波器的幅频响应和相频响应曲线
subplot(2,1,1);
plot(F,20*log10(abs(H)),'b');
hold on;
plot([0,fc1,fc1,fc2,fc2,Fs/2]/(Fs/2),[-Rs,-Rs,-Rp,-Rp,-Rs,-Rs],'r--');
axis([0,1,-50,5]);
xlabel('Normalized Frequency');
ylabel('Magnitude (dB)');
title('Magnitude Response of FIR Filter');

subplot(2,1,2);
plot(F,angle(H),'b');
hold on;
plot([0,fc1,fc1,fc2,fc2,Fs/2]/(Fs/2),[-pi,-pi/2,0,0,-pi,-pi],'r--');
axis([0,1,-pi,pi]);
xlabel('Normalized Frequency');
ylabel('Phase (Radians)');
title('Phase Response of FIR Filter');

% 画出输入信号和输出信号的时域波形图
figure();
subplot(2,1,1);
plot(n,x);
xlabel('Time (s)');
ylabel('Amplitude');
title('Input Signal');

subplot(2,1,2);
plot(n,y);
xlabel('Time (s)');
ylabel('Amplitude');
title('Output Signal');