《DSP using MATLAB》Problem 7.32
2024-10-08 00:02:29
代码:
%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
%% Output Info about this m-file
fprintf('\n***********************************************************\n');
fprintf(' <DSP using MATLAB> Problem 7.32 \n\n'); banner();
%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ w1 = 0; w2 = 0.4*pi; delta1 = 0.05; gain1 = 0.4;
w3 = 0.5*pi; w4 = 0.7*pi; delta2 = 0.05; gain2 = 0.95;
w5 = 0.8*pi; w6 = pi; delta3 = 0.025; gain3 = 0.0; f = [w2, w3, w4, w5]/pi; m = [1, 1, 0]; delta = [delta1, delta2, delta3]; [N, f0, m0, weights] = firpmord(f, m, delta);
fprintf('\n-------------------------------------------\n');
fprintf('\n Results by firpmord function:\n');
N
%f0
%m0
%weights
fprintf('\n-------------------------------------------\n'); weights = [1 delta1/delta2 delta1/delta3]
deltaH = max([delta1,delta2,delta3]); deltaL = min([delta1,delta2,delta3]); %Dw = min((w3-w2), (w5-w3));
%M = ceil((-20*log10((delta1*delta2*delta3)^(1/3)) - 13) / (2.285*Dw) + 1) f = [ 0 w2 w3 w4 w5 pi]/pi;
m = [ 0.4 0.4 0.95 0.95 0.025 0.025]; h = firpm(N+1, f, m, weights); % even-number order
[db, mag, pha, grd, w] = freqz_m(h, [1]);
delta_w = 2*pi/1000;
w1i = floor(w1/delta_w)+1; w2i = floor(w2/delta_w)+1;
w3i = floor(w3/delta_w)+1; w4i = floor(w4/delta_w)+1;
w5i = floor(w5/delta_w)+1; w6i = floor(w6/delta_w)+1; Asd = -max(db(w5i:w6i)) %[Hr, ww, a, L] = Hr_Type1(h);
[Hr,omega,P,L] = ampl_res(h); M = N+2
l = 0:M-1; %% -------------------------------------------------
%% Plot
%% ------------------------------------------------- figure('NumberTitle', 'off', 'Name', 'Problem 7.32')
set(gcf,'Color','white'); subplot(2,2,1); stem(l, h); axis([-1, M, -0.3, 0.6]); grid on;
xlabel('n'); ylabel('h(n)'); title('Actual Impulse Response, M=27');
set(gca,'XTickMode','manual','XTick',[0,M-1])
set(gca,'YTickMode','manual','YTick',[-0.3:0.1:0.6]) subplot(2,2,2); plot(w/pi, db); axis([0, 1, -80, 10]); grid on;
xlabel('frequency in \pi units'); ylabel('Decibels'); title('Magnitude Response in dB ');
set(gca,'XTickMode','manual','XTick',f)
set(gca,'YTickMode','manual','YTick',[-60,-26,0]);
set(gca,'YTickLabelMode','manual','YTickLabel',['60';'26';' 0']); subplot(2,2,3); plot(omega/pi, Hr); axis([0, 1, -0.2, 1.2]); grid on;
xlabel('frequency in \pi nuits'); ylabel('Hr(w)'); title('Amplitude Response');
set(gca,'XTickMode','manual','XTick',f)
set(gca,'YTickMode','manual','YTick',[0,0.05,0.35,0.45,0.9,1]); subplot(2,2,4); axis([0, 1, -deltaH, deltaH]);
b1w = omega(w1i:w2i)/pi; b1e = (Hr(w1i:w2i)-m(1)); %b1e = (Hr(w1i:w2i)-m(1))*weights(1);
b2w = omega(w3i:w4i)/pi; b2e = (Hr(w3i:w4i)-m(3)); %b2e = (Hr(w3i:w4i)-m(3))*weights(2);
b3w = omega(w5i:w6i)/pi; b3e = (Hr(w5i:w6i)-m(5)); %b3e = (Hr(w5i:w6i)-m(5))*weights(3);
plot(b1w,b1e,b2w,b2e,b3w,b3e); grid on;
xlabel('frequency in \pi units'); ylabel('Hr(w)'); title('Error Response'); %title('Weighted Error');
set(gca,'XTickMode','manual','XTick',f);
set(gca,'YTickMode','manual','YTick',[-deltaH,0,deltaH]); figure('NumberTitle', 'off', 'Name', 'Problem 7.32 AmpRes of h(n), Parks-McClellan Method')
set(gcf,'Color','white'); plot(omega/pi, Hr); grid on; %axis([0 1 -100 10]);
xlabel('frequency in \pi units'); ylabel('Hr'); title('Amplitude Response');
set(gca,'YTickMode','manual','YTick',[0,delta2, 0.35, 0.45, 0.95-delta1, 0.95, 0.95+delta1]);
%set(gca,'YTickLabelMode','manual','YTickLabel',['90';'40';' 0']);
set(gca,'XTickMode','manual','XTick',[0,0.4,0.5,0.7,0.8,1]);
运行结果:
最新文章
- MySQL中如何插入反斜杠,反斜杠被吃掉,反斜杠转义
- MongoDB入门三:MongoDB shell
- 修改Azure Website 移动服务 默认时区
- python学习笔记13(模块、包)
- ACdream 1114(莫比乌斯反演)
- 馋-c语言的规则
- 敏捷冲刺每日报告--day1
- 不定高元素的高度transition动画实现
- css_css 盒子水平居中 垂直居中
- Vue2全家桶之二:vue-router(路由)详细教程,看这个就够了
- ApplicaitionContext妙用request解耦合
- install kubernetes cluster k8s集群安装
- maven中修改可用的版本
- js for in 获得遍历数组索引和对象属性
- html5-div布局
- C++实现的服务不可用发送报警邮件
- 去除UITableView多余的seperator
- linux CMakeLists.txt 语法
- 油田(DFS)
- Linux mount指令
热门文章
- 【Vue】组件的基础与组件间通信
- python 对excel操作
- <;scrapy爬虫>;爬取quotes.toscrape.com
- 【LGP5349】幂
- hadoop 环境下不知道yarn端口可以通过此命令查找
- Python+Selenium基础入门及实践
- C++:多线程002
- Java监控工具介绍,VisualVm ,JProfiler,Perfino,Yourkit,Perf4J,JProbe,Java微基准测试【转】
- vue之vant组件下拉加载更多
- vue + element-ui实现简洁的导入导出功能