版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领
文档简介
实验七_____________数字滤波器设计_____________
实脸室名称:信息学院2204实验时间:2015年11月26日
姓名:蒋逸恒学号:20131120038专业:通信工程指导教师:陶大鹏
成绩
教师签名:年月日
、实验目的
1、学会使用matlab设计符合要求的数字滤波器的参数。
2、学会使用基本的matlab命令,把获得的数字滤波器参数带到传输函数
中写出传输函数,之后利用之前的实验知识,画出滤波器级联框图,实现滤波器。
二、实验内容
Q7.1用MATLAB确定一个数字无限冲击响应低通滤波器所有四种类型的最
低阶数。指标如下:40kHz的抽样率,4kHz的通带边界频率,8kHz的阻带边
界频率,0。5dB的通带波纹,40dB的最小阻带衰减.评论你的结果。
Q7o2用MATLAB确定一个数字无限冲击响应高通滤波器所有四种类型的最
低阶数。指标如下:3500Hz的抽样率,1050Hz的通带边界频率,600Hz的阻带
边界频率,1dB的通带波纹,50dB的最小阻带衰减。评论你的结果。
Q7.5通过运行程序P7。1来设计巴特沃兹带阻滤波器。写出所产生的传输
函数的准确表达式。滤波器的指标是什么?你的设计符合指标吗?使用
MATLAB,计算并绘制滤波器的未畸变的相位响应及群延迟响应。
Q7o6修改程序P7.1来设计符合习题Q7.1所给指标的切比雪夫1型低通
滤波器。写出所产生的传输函数的准确表达式。你的设计符合指标吗?使用
MATLAB,计算并绘制滤波器的未畸变的相位响应及群延迟响应.
Q7o20使用函数fir1,设计一个线性相位有限冲击响应低通滤波器,使其
满足习题Q7.13给出的指标,并画出其增益和相位响应.使用习题Q7.13中用凯
泽公式估计出的阶数.用表格形式显示滤波器系数.你的设计满足指标吗?若不
满足,涧整滤波器阶数直到设计满足指标.满足指标的滤波器阶数是多少?
Q7o23用凯泽窗设计一个有限冲击响应低通滤波器。滤波器的指标是:
wp=0.31,ws=0.41,As=50dBo注意,函数kaiser需要参数。及阶数N的值,
他们必须先用式(7.36)和式(7。37)分别算出.你的设计满足指标吗?
Q7.25用fir2设计一个95阶有限冲击响应滤波器,它具有三个不同的常
熟幅度级:在频率范围0到0.25中为0.4,在频率范围0.3到0。45中为1.0,
在频率范围0.5到1.0中为0。8o画出所设计的滤波器的幅度响应。你的设
计满足指标吗?
Q7o27用remez设计具有如下指标的有限冲击响应带通滤波器:通带边界
为1O8kHz和3.0kHz,阻带边界为1.5kHz和4。2kHz,通带波纹bp=0.1,阻带
•波纹bs=0。02,抽样率为12kHz。用kaiserord估计滤波器的阶数。你的设计
是一个最优有限冲击响应滤波器吗?你的设计满足指标吗?若不满足,增加滤波
器阶数在满足指标方面有用吗?指标由一个较低阶数的滤波器来满足而不是由
kaiserord得到的来满足吗?在不等过渡带的情形下,用设计的滤波器可能在较
大的过度带宽中以增益响应表现不满意的行为.改进该行为的一种方法是:通过
移动阻带边界减少过渡带宽,直到使设计在过渡带中以平滑的下降来满足指标。
在通带边界保持固定的情况下,尝试这种方法并确定新的指标,它在过渡带中提
供平滑的下降。
三、实验器材及软件
1.微型计算机1台
2.MATLAB7.0软件
四、实脸原理
IIR数字滤波器采用递归型结构,即结构上带有反馈环路。IIR滤波器
运算结构通常由延时、乘以系数和相加等基本运算组成,可以组合成直接型、正
准型、级联型、并联型四种结构形式,都具有反馈回路。由于运算中的舍入处
理,使误差不断累积,有时会产生微弱的寄生振荡.
相对的,FIR滤波器是直接型结构,没有递归结构,也就没有反馈环路,
FIR滤波器运算结构通常也是由由延时、乘以系数和相加等基本运算组成,可以
组合成各种结构形式的滤波器组来实现滤波。
IIR滤波器的单位脉冲响应为无限长,网络中有反溃回路。
FIR(FiniteImpuIseResponse)滤波器的单位脉冲响应是有限长的,一
般网络中没有反馈回路。
IIR滤波器的系统函数一般是一个有理分式,分母多项式决定滤波器的反
馈网络。
FIR滤波器的系统函数一般是一个分母为1的多项式,所以没有反馈网络。
IIR数字滤波器在设计上可以借助成熟的模拟滤波器的成果,如巴特沃斯、
契比雪夫和椭圆滤波器等,有现成的设计数据或图表可查,其设计工作量比较
小,对计算工具的要求不高。在设计一个IIR数字滤波器时,我们根据指标先
写出模拟滤波器的公式,然后通过一定的变换,将模拟滤波器的公式转换成数字
滤波器的公式。
FIR的设计在matlab中其实也同样很简单,其主要思想是通过先设计出无
限长的IIR滤波器再通过截尾的方式获得一个FIR滤波落但是简单的截尾会出
现吉布斯现象,所以此时还需要对截尾之后的结果加窗;来减小吉布斯震荡现
象。
IIR数字滤波器幅频特性精度很高,不是线性相位的,可以应用于对相位
信息不敏感的音频信号上;FIR数字滤波器的幅频特性精度较之于IIR数字滤波
器低,但是线性相位,就是不同频率分量的信号经过FIR滤波器后他们的时间差
不变,这是很好的性质.FIR数字滤波器是有限的单位响应也有利于对数字信号
的处理,便于编程.
五、实验步骤
1、进行本实验,首先必须熟悉matIab的运用,所以第一步是学会使用matIab。
2、学习相关基础知识,根据《数字信号处理》课程的学习理解实险内容和目的。
3、在充分熟悉基础知识的情况下进行实脸,利用matIab完成各种简单的波形
产生和观察,理解各种波形产生的原理和方法。
4、从产生的图形中学习新的知识,掌握实验的目的,充分学习数字滤波器的设
计和运用。
5、最后需要思考各种滤波器的联系和区别,以及如何根据要求设计出最符合要
求且最容易实现的滤波器
7.5
输入:
Ws-[0.40.6];
Wp=[0.20。8];
Rp=0。4;Rs=50;
[N1,Wnl]=buttord(Wp,Ws,Rp,Rs);
[num,den]=butter(N1,Wn1,'stop');
disp('分子系数:');disp(num);
disp('分母系数:’);disp(den);
[g,w]=gain(num,den);
plot(w/pi,g);grid
axis([01-605]);
xIabeI['\omega/\pi');ylabel('增益,dB*);
title。巴特沃斯带阻滤波器的增益响应’);
输出:
分子系数:
0o04930.00000.24650.00000o49300.00000.4930
0o00000o24650.00000.0493
分母系教:
1.00000o0000—0o08500o00000.63600.0000-0.0288
0.00000.05610.0000—0.0008
相位响应和群延迟响应:
w=0:pi/255:pi;
Ws=[0o40o6];Wp=[0o2Oo8]:Rp=0o4:Rs=50;
[N1,Wnl]=buttord(Wp,Ws,Rp,Rs);
[num,den]=butter(N1,Wn1,'stop*);
h=freqz(num,den,w);
subplot(1,2,1);
pIot(w/zpi,unwrap(angIe(h)));
gridon;
k=[diff(angIe(h))。/diff(w)0];
subplot(1,2,2);
pIot(w/pi,k);
axis([0,1,3,10]);
gridon;
7.6
输入:
Ws=0.4;Wp=0o2;Rp=0o5;Rs=
[N1,Wn1]=buttord(Wp,Ws,Rp,Rs)
[num,den]=butter(N1,Wn1);
disp(?分子系数:');disp(num);
disp('分母系数:’);disp(den);
[g,w]=gain(num,den):
plot(w/pi,g);grid
axis(.01—605]);
xlabel]'\omega/\pi');
ylabelC增益,dB,);
title('切比雪夫1型低通滤波器的增益响幺
输出:
分子系数:
0o00040o00200.00400.00400。00200.0004
分母系教:
1.0000-3.82696。2742-5.44642.4915一0。4797
相位响应和群延迟响应:
Ws=0。4;Wp=0.2;Rp=0o5;Rs=40;
%估计滤波器阶数
[N1,Wn1]=cheblord(Wp,Ws,Rp,Rs):
%设计成波器
[num,den]=chebyl(N1,Rp,Wn1);
h=freqz(num,den,w);
subplot(1,2,1);
pIot(w/zpi,unwrap(angle(h)));
gridon;
k=[diff(angIe(h))./diff(w)0];
subplot(1,2,2);
plot(w/pi,k);
axis([0,1,-5,20]);
gridon;
15
w=0:pi/255:pi;
Fs=2500;Fp=2000;
Rp=0.005;Rs=0o005;
Ft=10000;
N=kaiord(Fp,Fs,Rp,Rs,Ft);
b=fir1(N,2*Fp/Ft);
[g,w]=gain(b,1):
subplot(2,1,1);plot(w/pi,g);
axis([01-1255])
xIabeI('\omega/\pi');
ylabel('增益,dB');
titIe('FIR低通滤波器的增益响应'
h=freqz(b,1,w);
subplot(2,1,2);
pIot(w/pi,unwrap(angle(h)));gridon;
xIabeI('\omega/\pi');yIabeI('相位');
title('FIR低通滤波器的相频响应'):
输出:
b=
。
-0o00070o00070o0014-00000—0o0023—0o00190.0025
0o0052—0.0000—0.0083—0.00640o00790.0157-0.0000
—0。0233-0«01760o02150o0431-0o0000-0.0706—0o0600
0o09190.30130o39980.30130.0919-0.0600—0.0706
-0o00000o04310.0215-0.0176—0.0233-0o0000Oo0157
0.0079—0.0064—0o0083—0.00000o00520.0025-0.0019
—0。0023-0o00000.00140o0007-0o0007
7o23
输入:
w2=0:pi/255:pi:
Wp=0.31;Ws=0o41;
Wn=Wp+(Ws-Wp)/2;
As=50;Rs=1(T(—As/20):
Rp=Rs;
ifAs>21
N-round((AG-7O95)*2/(14.36*(abo(WP-WG)))+1);
eIse
N=round(0.9222*2/abs(Wp-V/s)+1);
end
ifAs>50
B=0o1102*(As—8o7);
eIseifAs>=21
B=0o5842*(As-21)"0o4+0o07886*(As-21);
eIse
B=0;
end
s=kaser(N+1,B);
h=fir1(N,Wn,s);
disp('系数:');
disp(h):
[g,w]=gain(h,1);
subpIot(2,1,1):
plot(w/pi,g);解缠绕的相位响应
grid:
xlabel('\omega/\pi');
yIabeI('增益,dB');
title。增益响应’);
Hz=freqz(h,1,w2);
subplot(2,1,2);
pIot(w2/pi,unwrap(angle(Hz)));“。o-20.30.40.50.60.70.80.91
co/n
grid:
xIabeIV\omega/\pi');yIabeI「以弧度为单位的相位');
title('解缠绕的相位响应’):
输出:
系数:
0.00030o00080.0003-0.0011■-0.00170.00000o0026
0o0027—0.0010-0o0049-0.00350o00330o00800o0034
-0.0074—Oo0119—0.00180.01400o0161-0.0027—0.0241
—0。02010o01270.04060o0236-0o0354—0o0754-0.0258
0o12140o28710o35970o28710o1214-0o0258-0.0754
一0。03540o02360.04060o0127—0.0201-0.0241-0.0027
0.01610o0140—0.0018-0.0119-0o00740.00340。0080
0.0033—0.0035-0o0049-0.00100.0027Oo00260。0000
-0.0017—Oo00110.00030.00080.0003
7.25
解缠绕的相位响应
b=fir2(95,w,a);
disp('系数:');
diGP(b);
[g,w]=gain(b,1);
subplot(2,1,1);
plot(w/pi,g);
grid;
xlabel('\omega/\pi');
ylabelC增益,dB');
title('增益响应’);
h=freqz(b,1,w2);
subpIot(2,1,2);
pIot(w2/pi,unwrap(angIe(h)));grid;
xlabel('\omega/\pi'):yIabeI('以弧度为单位的相位');
title(?解缠绕的相位响应');
输出:
—0.00000.00000o00000o0000—Oo0000-0o0000Oo0000
0。0000—0.0000—Oo0000—0o0001-0.00010.00010.0003
0o0003-0.0001—0.0005-0.0005—0.0004-0.0000Oo0008
0.00160.0010-Oo0013—0.0030—0.00190o00080o0025
0o00290.00280.0004-0o0053-Oo0089-0o00340。0076
0o01210o0060—0o0024-0.0078—0o0136—0o0161-0.0003
0o03160.04440o0084—0o0537—0o0885-0.08100o7307
一0。0810-0.0885—0o05370o00840o04440。0316—0.0003
-0.0161—Oo0136-0.0078-0.00240.00600.01210o0076
—0.0034-0o0089-0.00530.00040.0028Oo00290.0025
0.0008—0.0019—0.0030-0o00130o00100o00160o0008
—0。0000—0.0004-0.0005—0.0005-0o00010o00030.0003
0.0001—Oo0001-0.0001-0o0000—Oo00000.00000o0000
-0.0000-0o00000.00000o00000.0000-0.0000
7O27kaiserord获得阶数计算出的结
果(73阶)
输入:增益.应
LAI»6I4I5t
w2=0:pi/255:pi;
Fs1=1500;
Fp1=1800;
0.10.20.30.40.50.60.70.80.9
Fp2=3000;Fs2=4200:(o/x
解缎线的相位明应
Ft=12000;
ws1=2*Fs1/Ft:
wp1=2*Fp1/Ft;
wp2=2*Fp2/Ft;
ws2=2*Fs2/Ft;0.10.20.30.40.50.60.70.80.9
3区
Dp=0.1;Ds=0.02;
F=[Fs1Fp1Fp2Fs2];比kaiserord获得的阶数更低的
A=[010];阶数(20阶)得到的输出结果:
DEV=[DsDpDsl:增益响应
3,太
解绅经的相位峋应
[N,Wn,B,Ftap]=kaiserord(F,A,DEV,Ft);
F2=[0ws1wp1wp2ws21];
A2-[001100];
h=renez(N,F2,A2);
disp('系数:');
disp(h);
[g,w]=gain(h,1);
subplot(2,1,1);
plot(w/pi,g);
grid:
axis([01—8060]);
xlabel('\omega/\pi');移动阻带边界(将4.2KHz移动到3.4KHz)来减
少
ylabel「增益,dB');过渡带宽后kaiserord获得阶数后计算的结果
(73阶)增苗明应
title('增益响应’);
*
Hz=freqz(h,1,w2);
S3
subplot(2,1,2);a
plot(w2/pi,unwrap(angle(Hz)));0.40.50.6
grid;留级税的相位峋应
xlabel('\omega/\pi');0
邛
⑤
ylabel('以弧度为单位的相位’);0
&
title('解缠绕的相位响应');依
诬
乐
不
输出:
不使用kaiserord来状得阶数(73阶),而是用比kaiserorc狄得的阶数更低的一个阶
数(20阶)得到的揄出结果(图像记录是右边从上至下的第二个):
系数:
0o0058—0.1330—0.02100.0709—0。00150.01460.1225
—0.0810—Oo29530o04670.37890。0467-0。2953—0.0810
0.12250o0146—0.00150.0709-0.0210—0。13300o0058
移动阻带边界(<4.2KHz移动到3.4KHz)来减少过渡带宽,获得的输出,此时,阶数仍
然是73阶(图像记录是右边从上至下的第三个):
系数:
-0o00470o00660o00760.0007—0.0033-0.0002-0.0015
—0.00640.00030o01200.0069—Qo0071-0。00610o0002
-0o0068—0o00730o01400o0219—0.0034一0。0194—C.0047
-0.0011-0o01450o00560o04270.0205—0。0368—0.0337
0o0046—0o0102-0o01760.07050。1129-0.0493-0.2139
-0o07750o20490.2049-0.0775-0.2139-0.04930.1129
0.0705—0.0176-0.01020.0046—0。0337—0。03680o0205
0.04270.0056—0o0145—0o0011-0.0047-0。0194-0o0034
0.02190o0140—0.0073—0o00680。0002-0。0061—0.0071
0.00690.01200o0003-0.0064—0.0015-0.0002-0o0033
0。00070.00760.0066-0.0047
七、实验思考题及解答
7O1归一化频率吗,=等=0.2%,这里角频率取0.2;归一化频率
吗=管=。4万,这里角频率取0.4;理想通带波纹Rp=0.5dB;最小阻带衰减
Rs二04dB.
巴特沃斯低通滤波器最低阶数N=8,截止频率Wn=0。2469
切比雪夫1型低通滤波器最低阶数N=5,截止频率Wn=0.2
切比雪夫2型低通滤波器最低阶数N=5,截止频率Wn=0.4
椭圆滤波器低通滤波器最低阶数N=4,裁止频率Wn=0。2
从以上结果中观察到椭圆滤波器的阶数最低,阶数为4,截止频率也最低,
为0。2,并且符合要求;此时的切比雪夫1型低通滤波器的截止频率也是最小
的0.2;而巴特沃斯滤波器的阶数最高,实现最困难;截止频率最大的是切比雪
夫2型,为0。4.
1
7.2归一化频率w=苧=0.6万,这里角频率取0。6;归一化频率
PrT
吗=等=0.34),这里角频率取0。34:理想通带波纹Rp=1dB;最小阻带衰
减Rs=50dBo
巴特沃斯低通滤波器最低阶数N=8,截止频率W『0.5615
切比雪夫1型低通滤波器最低阶数N=5,截止频率Wn=0.6
切比雪夫2型低通滤波器最低阶数N=5,截止频率Wn=0.34
椭圆滤波器低通滤波器最低阶数N=4,截止频率W产0。6
从以上结果中观察到椭圆滤波器的阶数最低,阶数为4,并且符合要求;
此时的切比雪夫2型低通滤波器的裁止频率是最小的0。34;而巴特沃斯滤波
器的阶数最高,实现最困难;栈止频率最大的是切比雪夫1型和椭圆滤波器,
为0o6o
7O5传输函数的表达式:
0.0493-h0.2465Z2+0.493z4+0.493z6+0,2465z8+O.(M93z-10
1-0.085z-2+0.636z4+0.0288z6+0.561z8+0.0008z-10
滤波器参数是:通带归一化裁至频率Wp1=0o2n,Wp2=0o8r[,阻带归一化
频率Ws1=0.4n,Ws2=0o6n,通带波纹Rp=O。4dB,最小阻带衰减Rs=50dB.
设计出来的滤波器的通带归一化截至频率Wp1=0。25n,Wp2=0o75n,很接
近指标参数,所以符合指标。
7.6根据7.1可知归一化频率%,=符=0.2〃,这里角频率取0.2;归一化频
2死
率叱==0.4万,这里角频率取0。4;理想通带波纹Rp=0o5dB;最小阻带
FT
衰减Rs=40dBo
传输函数的表达式,“㈠二0001H+S002z"+0-®z-2+SOOfe-3+O.QQ2Z7+Q(XXHz,
“1-3.8269z-1+6.2742z-2-5.4464z-3+2.4915z4-0.4797z
设计出来的滤波器的通带归一化截至频率Wp=0o2n,等于指标参数,所以符
合指标。
7.20根据7。13可知频率Fp=2000Hz;Fs=2500Hz;理想通带波纹Rp=00005dB;
阻带波纹Rp=0.005dB;抽样频率Fk10000Hz。
通过使用kaiord函数求的阶数是N=46,滤波器系数如下(这里没有使用
kaiserord来求取阶数,因为kaiserord会改变通带截止频率,从而使估算出的
阶数对应的截止频率不一定满足要求,而kaiord可以直接得出效果很好且满足
要求的阶数):
-0o00070.00070o0014—0.0000-0o0023—0o00190o0025
0«0052—0.0000-0o0083-0«00640o00790o0157—0o0000
—0.0233-0o01760o02150。0431—0.0000—0.0706—0.0600
0o09190o30130o39980.30130o0919-0.0600-0o0706
-0o00000o04310.0215—0.0176—0。0233-0o00000o0157
0.0079-0o0064-0o0083-0o00000o00520.0025—0o0019
—0.0023—0«00000.00140o0007—0o0007
通过计算可以看到增益响应和相位响应图的结果能很好的符合指标要求。
7.23从图中可以看出设计的滤波器满足指标。此时N=61。
7.25由题意可得频率点分布为w=[00o250o30.450.51.0];对应的
幅值分布a=[0.40。41.01.00.80。8].设计的滤波器见实验记录.从幅
度响应中可以看出,此滤波器满足指标.
7.27设计出来的滤波器应该是在该阶数(73阶)下最优的滤波器设计,因为使
用的remez函数就是采用优化方法设计最优的标准多频带FIR数字滤波器,但并
不是真正最优的滤波器,因为在过渡带中出现了上升再下降的情况,而且相位
在通带内也并不是完全线性的,这并不是FIR滤波器应该具有的特性。说明设
计不能完全满足指标C
增加阶数也并不能消除出现的过渡带中的意外上升再下降的状况和相位非
线性的状况,即并不能通过增m阶数来改善滤波器指标。
但是当阶数变为20阶时,我们可以看到输出图像中得到的结果却满足了指
标要求,图像见实验记录。
通过移动阻带边界来减少过渡带宽,最终获得了较为满意的输出结果,确定
的新的指标是:通带边界为1。8kHz和3。0kHz,阻带边界为1。5kHz和3.4kHz,
通带波纹bp=O。1,阻带波纹bs=0.02,抽样率为12kHz。此时的阶数为73阶。
八、实脸结果分析与总结
结果分析:
7017.2:这两道题是IIR滤波器的设计时获取最低阶数,通过结果我们
可以看到对于不同类型的滤波器,输出的结果是不一样的,即每一种滤波器对
于同一组指标都有着差别,但是从实现的结果可以看到,虽然图像上差别会比较
大,但是对于通带内的差别是很小的,所以在实际应用中,这几种滤波器都可以
被应用到。
7057O6:这两道题是对滤波器的设计,设计出来的IIR滤波器不但性
能好,而且阶数低,在数字滤波领域非常易于实现,从这里我们可以看到IIR
滤波器的优秀之处,即设计简单,实现容易。但是通过学习,我们知道IIR的
相位并不是线性的,所以数字信号处理后可能引入群延迟和相位延迟等问
温馨提示
- 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
- 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
- 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
- 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
- 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
- 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
- 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。
最新文档
- 2025版风湿性关节炎常见症状及护理技巧
- 安全教育不碰电
- 语文中的观察方法
- 综合评价方法与决策数学建模
- 日本建筑工程介绍
- AEO企业高级认证介绍
- 白板课件制作方法
- 智能康复技术科普讲解
- 心血管内科进修护士工作汇报
- 呼吸内科护理工作汇报总结
- 2025全国导游资格证考试《导游业务》真题库(含答案)
- 中国软件行业协会:2025中国软件行业基准数据报告 SSM-BK-202509
- GB/T 4026-2025人机界面标志标识的基本和安全规则设备端子、导体终端和导体的标识
- 2025年陕西延安旅游集团有限公司招聘笔试参考题库含答案解析
- 减少我们的碳排放-课件(17张)
- 体能训练概论(NSCA)
- Q∕SY 1736-2014 评标方法选择和评标标准编制规范
- 食品风味化学-6食品风味的调整和香味料
- 国家开放大学电大专科《美学与美育》简答题综合论述题题库及答案
- 闵行区加强住宅小区综合治理三年行动计划(2015-2017)
- ZZG22高频开关整流器 - 图文-
评论
0/150
提交评论