收藏 分享(赏)

数学建模——声音识别模型的建立与评价.doc

上传人:精品资料 文档编号:10197167 上传时间:2019-10-18 格式:DOC 页数:20 大小:206.59KB
下载 相关 举报
数学建模——声音识别模型的建立与评价.doc_第1页
第1页 / 共20页
数学建模——声音识别模型的建立与评价.doc_第2页
第2页 / 共20页
数学建模——声音识别模型的建立与评价.doc_第3页
第3页 / 共20页
数学建模——声音识别模型的建立与评价.doc_第4页
第4页 / 共20页
数学建模——声音识别模型的建立与评价.doc_第5页
第5页 / 共20页
点击查看更多>>
资源描述

1、1数学建模论文报告题目:A 题 声音识别模型的建立与评价参赛队员(数学与统计学院):学号: 姓名: 联系电话:2声音识别模型的建立与评价【摘要】本文针对正常非正常开门(指盗窃开门等声音)的声音进行识别的问题,通过 matlab 的 sound 和 plot 采集到了正常和非正常开门的声音信号和声音波形图,附件中有正常开门声音(如正 1.mat),非正常开门声音(如非 1.mat),各 40 次开门,共 80 次开门声音数据。将这些数据利用 matlab 的 load 函数载入到计算机内存,内存中变量有 Fs 和 y 等变量,其中 Fs 为采用频率,y 为采用数据,通过建立数学模型将函数关系表示

2、出来,并利用合适的时域或(和)频域特征表达声音信号,建立出特征向量。针对问题一:利用 matlab 中的 sound 函数,播放出声音信号,试听并比较正常和非正常开门声音的差别,利用 plot 函数绘制出具体的声音波形图,总结出差别在哪些方面,正常开门声音(如正 1.mat)较短暂,波形图上各点分布较为分散;而非正常开门声音(如非 1.mat)连续,波形图上各点分布也相对集中。针对问题二:利用合适的时域或(和)频域特征表达声音信号,建立特征向量,写出提取特征向量的具体方法和程序代码。针对问题三:建立声音识别模型(二分类模型) ,利用模型区分正常和非正常声音,评价模型的好坏。针对问题四:利用特征

3、选择或变换,对特征向量进行优化,并利用参数优化技术优化模型的参数,使识别模型的准确率提高。针对问题五:若原始声音信号中有环境噪声(如白噪声) ,设 y 中叠加了一定幅值的白噪声(利用 y1=y+(-0.15+0.3*rand(size(y)*max(y)叠加噪声) ,对声音进行前期处理。本文充分运用了高等代数、数理统计等知识,并通过 MATLAB 软件模拟的方法对理论数据进行了分析整合。最后对模型的优缺点进行了评价,并给出了改进方向。关键词:声音识别、特征向量、二分类模型、MATLAB 软件、神经网络31 问题重述A 题 声音识别模型的建立与评价随着家居智能化逐渐普及,智能冰箱、智能清洁机器人

4、、智能电视等已步入平常老百姓家庭,但智能化的防盗门还处于研发阶段,未进入市场。随着人们对家居安全意识的不断增强,对防盗、防抢和防砸的门禁系统的智能性提出更高的要求。基于此,对正常和非正常开门(指盗窃开门等声音)的声音进行识别是智能防盗门的关键问题和技术,其具有广泛的应用前景和实用价值。为了进行声音识别模型的建立,我们采集到了正常和非正常开门的声音,附件中有正常开门声音(如正 1.mat),非正常开门声音(如非 1.mat),各 40次开门,共 80 次开门声音数据。该数据可利用 matlab 的 load 函数载入到计算机内存,内存中变量有 Fs 和 y 等变量,其中 Fs 为采用频率,y 为

5、采用数据。利用这些数据要求完成以下工作:1. 利用 matlab 中的 sound 函数,播放出声音信号,试听并比较正常和非正常开门声音的差别,利用 plot 函数绘制出具体的声音波形图,总结差别在哪些方面?2. 利用合适的时域或(和)频域特征表达个声音信号,建立特征向量,写出提取特征向量的具体方法和程序代码。3. 建立声音识别模型(二分类模型) ,利用模型区分正常和非正常声音,评价模型的好坏。4. 试利用特征选择或变换,对特征向量进行优化,并利用参数优化技术优化模型的参数,使识别模型的准确率提高。5. 若原始声音信号中有环境噪声(如白噪声) ,设 y 中叠加了一定幅值的白噪声(利用 y1=y

6、+(-0.15+0.3*rand(size(y)*max(y)叠加噪声) ,4如何对声音进行前期处理?2 模型假设1 假设原始声音信号中有环境噪声(白噪声)单一,无其他噪音干扰;2 假设采用的频率都为 11025Hz 来利用合适的时域或(和)频域特征表达声音信号;3 假设实验环境保持不变。3 符号说明Fs 采用频率y 采用数据 4 问题的分析本文针对正常非正常开门(指盗窃开门等声音)的声音进行识别的问题,通过 matlab 的 sound 和 plot 采集到了正常和非正常开门的声音信号和声音波形图,附件中有正常开门声音(如正 1.mat),非正常开门声音(如非 1.mat),各 40 次开门

7、,共 80 次开门声音数据。将这些数据利用 matlab 的 load 函数载入到计算机内存,内存中变量有 Fs 和 y 等变量,其中 Fs 为采用频率,y 为采用数据,通过建立数学模型将函数关系表示出来,并利用合适的时域或(和)频域特征表达声音信号,建立出特征向量。若原始声音信号中有环境噪声(如白噪声) ,设 y 中叠加了一定幅值的白噪声(利用 y1=y+(-0.15+0.3*rand(size(y)*max(y)叠加噪声) ,对声音进行前期处理。55 模型的建立与求解5.1 利用 matlab 中的 sound 函数,播放出声音信号,试听并比较正常和非正常开门声音的差别,利用 plot 函

8、数绘制出具体的波形图load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集正 1.mat)sound(y,Fs)plot(y)图 1 正常开门声音(正 1.mat)图load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集非 1.mat)6sound(y,Fs)plot(y)图 2 非正常开门声音(非 1.mat)图由图像可总结出差别在于:正常开门声与非正常开门声的区别在于,正常开门声是连续的短暂的,而非正常开门声是持续的开锁声音,用人的听觉很容易就能判断是否为正常的开门声音。当我们运用

9、matlab 来分析声音时刻看出在频域内,语音信号的频谱分量主要集中在 3003400Hz 的范围内。如图 1 正常开门时域波形图等。在时域内,两种声音信号都具有“短时性”的特点,即在总体上,开门声的特征是随着时间而变化的,但在一段较短的时间间隔内,语音信号保持平稳。但从时域波形图的对比可看出两种开门声的区别,正常开门声有一组的峰值比其他的峰值明显要大很多,而非正常开门声有多组峰值相似。5.2 利用合适的时域或(和)频域特征表达个声音信号,建立特征向量,写出提取特征向量的具体方法和程序代码。7由于已假设采用的频率都为 11025Hz 所以需要利用合适的时域或(和)频域特征表达个声音信号。在将开

10、门声音进行数字化前,必须先进行防混叠预滤波,预滤波的目的有两个:抑制输入信导各领域分量中频率超出 fs/2 的所有分量(fs 为采样频率),以防止混叠干扰。抑制 50Hz 的电源工频干扰。这样,预滤波器必须是一个带通滤波器,设其上、下截止颜率分别是 fH 和 fL,则对于绝人多数语音编译码器,fH=3400Hz、fL60100Hz、采样率为 fs8kHz;而对于开门声音的识别而言,当判断是否为正常开门声有重要意义。当使用要求较高或很高的场合时fH4500Hz 或 8000Hz、fL60Hz、fs10kHz 或 20kHz。为了建立和提取不同开门声音的特征向量,必须经过采样和量化两个步骤,从而得

11、到时间和幅度上均为离散的数字语音信号。采样也称抽样,是信号在时间上的离散化,即按照一定时间间隔t 在模拟信号 x(t)上逐点采取其瞬时值。采样时必须要注意满足奈奎斯特定理,即采样频率 fs 必须以高于受测信号的最高频率两倍以上的速度进行取样,才能正确地重建波它是通过采样脉冲和模拟信号相乘来实现的。下图时一段语音信号在采样频率 11025Hz 情况下的频谱图。%语音信号时域频域分析y=y(:,1);%我这里假设你的声音是双声道,我只取单声道作分析,如果你想分析另外一个声道,请改成 y=y(:,2) sigLength=length(y); Y = fft(y,sigLength); Pyy =

12、Y.* conj(Y) / sigLength; halflength=floor(sigLength/2); f=Fs*(0:halflength)/sigLength; figure;plot(f,Pyy(1:halflength+1);xlabel(Frequency(Hz); t=(0:sigLength-1)/Fs; figure;plot(t,y);xlabel(Time(s);y=y(:,1);lg=length(y);8yx=(0:(lg-1)*Fs/lg;yx=yx(1:lg/2);% sound(y,fs);figure(1);subplot(211);plot(y);ti

13、tle(原始图像);subplot(212);yf=abs(fft(y);yf=yf(1:lg/2);plot(yx,yf);ycq=zeros(1,lg); %先将抽取后的值全设为零 for i=1:80:lg; %通过循环,每隔 80 个点将抽取后的值赋值为原函数的 ycq(i)=y(i); %采样值 end sound(ycq,Fs);figure(2);subplot(211);plot(ycq);title(提取特征向量后的图像);subplot(212);ycqf=abs(fft(ycq);ycqf=ycqf(1:lg/2);plot(yx,ycqf);0 2 4 6 8 10 1

14、2x 104-0.200.20.40.6 像像像像0 100 200 300 400 500 60001020300 2 4 6 8 10 12x 104-0.1-0.0500.050.1像像像像像像像像像像0 100 200 300 400 500 60000.20.40.60.8图 3 正常开门声提取的特征向量图90 1 2 3 4 5 6x 104-0.1-0.0500.050.1 像像像像0 100 200 300 400 500 60001230 1 2 3 4 5 6x 104-0.04-0.0200.020.04 像像像像像像像像像像0 100 200 300 400 500 6

15、0000.050.10.150.2图 4 非正常开门声提取的特征向量图5.3 建立声音识别模型(二分类模型) ,利用模型区分正常和非正常声音,评价模型的好坏。对于分类资料的分析,当要考察的影响因素较少,且也为分类变量时,常用列联表(contingency Table)的形式对这种资料进行整理,并使用 2检验来进行分析,汉存在分类的混杂因素时,还可应用 Mantel-Haenszel 检验进行统计学检验,这种方法可以很好地控制混杂因素的影响。但是这种经典分析方法也存在局限性,首先,它虽然可以控制若干个因素的作用,但无法描述其作用大小及方向,更不能考察各因素间是否存在交互任用;其次,该方法对样本含

16、量的要求较大,当控制的分层因素较多时,单元格被划分的越来越细,列联表的格子中频数可能很小甚至为 0,将导致检验结果的不可靠。最后, 2检验无法对连续性自变量的影响进行分析,而这将大大限制其应用范围,无疑是其致使的缺陷。标准的线性回归模型: A1mYx如果对分类变量直接拟合,则实质上拟合的是发生概率,参照前面线性回归方程 ,很自然地会想到是否可以建立下面形式的回归模型: 1mPx显然,该模型可以描述当各自变量变化时,因变量的发生概率会怎样变化,可以满足分析的基本要求,对计算方法进行了改进,用加权最小二乘法来对该模型进行拟合。以 Logit()为因变量,建立包含 p 个自变量的 logistic

17、回归模型如下:01logit()pPx以上即为 logistic 回归模型。由上式可推得:01exp()pPx 01e()px上面三个方程式相互等价。通过大量的分析实践,发现 logistic 回归模型10可以很好地满足对分类数据的建模需求,因此目前它已经成为了分类因变量的标准建模方法。优势比:如前所述,把出现某种结果的概率与不出现的概率之比称为比值( odds) ,即 1Pods。两个比值之比称为优势比( odds Ratio,简称 OR) 。首先考察OR 的特性:若 2,则 122Podsods若 1P,则 122若 2,则 122Podsods显然, OR 是否大于 1 可以用作两种情形

18、下发生概率大小的比较。5.4 试利用特征选择或变换,对特征向量进行优化,并利用参数优化技术优化模型的参数,使识别模型的准确率提高。1 用 MATLAB 实现神经网络建模MATLAB 提供了许多工具箱,其中神经网络工具箱是其中之一,该工具箱由许多子程序组成,能完成神经网络的构造和学习及所涉及的许多运算和操作。本文采用 BP 网络来建立神经网络模型来模拟实际系统间的关系。神经网络建模程序如下:Inputp=;Inputt=;%神经网络训练输入、输出矩阵pn,meanp,stdp,tn,meant,stdt=prestd(p,t);%神经网络输入、输出数归一化处理net=newff(minmax(p

19、n),4101,logsiglogsigpurelin,trainlm);%建立一个三层 BP 网络net.trainparam.show=25;net.trainparam.e2pochs=300;net.trainparam.goal=1e-10;%设置训练参数net,tr,epochs=train(net,pn,t);%神经网络训练11a=sim(net,pn);%装载训练结果m,b,r=postreg(a,t);%线性回归分析save(result,net);%保存网络结果在上述程序中,隐含层的神经元数可通过试凑法得到。神经网络训练数据进行了归一化预处理,归一化采用了工具箱中的 Pre

20、std函数。神经网络训练成功后保存结果,以便优化时载入已训练好的网络结构,调入过程模型。为了验证所建立的神经网络模型的有效性,将测试集数据输入训练完的神经网络模型,运用 Sim 函数进行仿真,将神经网络的输出与现场实测数据作比较,运用 Postreg 函数,将输出结果与目标值作线性回归分析,以此来检验网络泛化能力。2 用 MATLAB 实现遗传算法寻优在用神经网络建立数学模型后,采用遗传算法寻优,其程序框图如图 5 所示。遗传算法优化程序包括编码子程序(Encoding),解码子程序(Decoding)、选择子程序(Selection)、交叉子程序(Crossover)、变异子程序(Mutat

21、ion)组成。(1)编码子程序。本设计采用最常用的二进制编码方案,即用二进制数构成的符号串来表示个体,用 Encoding 函数来实现编码并产生初始种群。程序设计中,首先根据各决策变量的上、下界及其搜索精度来确定表示各决策变量的二进制串的长度,然后随机产生一个种群大小为 Popsize 的初始种群。(2)解码子程序。它将编码后的个体构成的种群通过 Decoding 函数解码以转换成原问题空间的决策变量,并通过已建好的神经网络模型求得各个体对应的输出值,并进行适应值的计算。12图 5 程序框图(3)选择子程序。它将解码后求得的各个体适应值采用最优保存策略和比例选择法,将适应值差的淘汰掉,选出一些

22、优良的个体以进行下一步的交叉和变异操作。该子程序首先找出当前群体中适应值最高和最低的个体,将最佳个体保留并用它替换最差个体。为保证当前最佳个体不被交叉、变异操作破坏掉,允许最佳个体不参与交叉、变异操作而直接进入下一代,然后将剩下的个体按比例选择法进行操作。(4)交叉子程序。采用单点交叉的方法来实现交叉算子,即按交叉概率 pc在两两配对的个体编码串中随机设置一个交叉点,然后在该点相互交换两个配对个体的部分基因,从而形成两个新的个体。(5)变异子程序。变异操作按照变异概率 pm 随机选择变异点,对于二进制基因串来说在变异点处将其位取反即可。5.5 若原始声音信号中有环境噪声(如白噪声) ,设 y

23、中叠加了一定幅值的白噪声(利用 y1=y+(-0.15+0.3*rand(size(y)*max(y)叠加噪声) ,如何对声音进行前期处理?白噪声:白噪声是功率谱密度为常数的随机过程或者或者随机信号,信号在各13个频段上的功率是一样的。 【白噪声的时域信号中任意两个不同时刻是不相关的,白噪声的自相关函数为冲击函数,白噪声的功率谱密度为常数。 (自相关函数和功率谱密度是傅立叶变换对) 】注:理想的白噪声具有无限带宽,因而其能量是无限大,这在现实世界是不可能存在的。实际上,我们常 常将有限带宽的平整讯号视为白噪音,因为这让我们在数学分析上更加方便。然而,白噪声在数学处 理上比较方便,因此它是系统分

24、析的有力工具。一般,只要一个噪声过程所具有的频谱宽度远远大于 它所作用系统的带宽,并且在该带宽中其频谱密度基本上可以作为常数来考虑,就可以把它作为白噪 声来处理。例如,热噪声和散弹噪声在很宽的频率范围内具有均匀的功率谱密度,通常可以认为它们 是白噪声。在信号处理中经常需要把噪声叠加到信号上去,在叠加噪声时往往需要满足一定的信噪比,这样产生二个问题,其一噪声是否按指定的信噪比叠加,其二怎么样检验带噪信号中信噪比满足指定的信噪比。在 MATLAB 中可以用 randn 产生均值为0方差为1的正态分布白噪声,但在任意长度下 x=randn(1,N),x 不一定是均值为0方差为1(有些小小的偏差) ,

25、这样对后续的计算会产生影响。在这里提供3个函数用于按一定的信噪比把噪声叠加到信号上去,同时可检验带噪信号中信噪比。1,把白噪声叠加到信号上去:function Y,NOISE = noisegen(X,SNR)% noisegen add white Gaussian noise to a signal.% Y, NOISE = NOISEGEN(X,SNR) adds white Gaussian NOISE to X. The SNR is in dB.NOISE=randn(size(X);NOISE=NOISE-mean(NOISE);signal_power = 1/length(X

26、)*sum(X.*X);noise_variance = signal_power / ( 10(SNR/10) );NOISE=sqrt(noise_variance)/std(NOISE)*NOISE;Y=X+NOISE;其中 X 是纯信号,SNR 是要求的信噪比,Y 是带噪信号,NOISE 是叠加在信号上的噪声。2,把指定的噪声叠加到信号上去有标准噪声库 NOISEX-92,其中带有白噪声、办公室噪声、工厂噪声、汽车噪声、坦克噪声等等,在信号处理中往往需要把库中的噪声叠加到信号中去,而噪声的采样频率与纯信号的采样频率往往不一致,需要采样频率的校准。function Y,NOISE = a

27、dd_noisem(X,filepath_name,SNR,fs)% add_noisem add determinated noise to a signal.% X is signal, and its sample frequency is fs;% filepath_name is NOISEs path and name, and the SNR is signal to noise ratio in dB.wavin,fs1,nbits=wavread(filepath_name);14if fs1=fswavin1=resample(wavin,fs,fs1);endnx=siz

28、e(X,1);NOISE=wavin1(1:nx);NOISE=NOISE-mean(NOISE);signal_power = 1/nx*sum(X.*X);noise_variance = signal_power / ( 10(SNR/10) );NOISE=sqrt(noise_variance)/std(NOISE)*NOISE;Y=X+NOISE;其中 X 是纯信号,filepath_name 是指定噪声文件(.wav)的路径和文件名,SNR 是要求的信噪比,fs 是信号 X 的采样频率,Y 是带噪信号,NOISE 是叠加在信号上的噪声。3,检验带噪信号的信噪比信噪比的定义为信号能

29、量 (纯信号)2SNR=-=-噪声能量load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集正 1.mat)y1=y+(-0.15+0.3*rand(size(y)*max(y);plot(y1)15图 6 原始声音信号 y 中叠加了一定幅值的白噪声图(正 1.mat)load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集非 1.mat)y1=y+(-0.15+0.3*rand(size(y)*max(y);plot(y1)16图 7 原始声音信号 y 中叠加了一定幅值的白噪声图(非 1

30、.mat)7 模型的评价模型的优点:(1)本文充分运用了高等代数、数理统计等知识,并通过 MATLAB 软件模拟的方法对理论数据进行了分析整合。(2)利用 matlab 中的 sound 函数,播放出声音信号,试听并比较正常和非正常开门声音的差别,利用 plot 函数绘制出具体的声音波形图,总结出差别。(3)假设了采用的频率都为 11025Hz 来利用合适的时域或(和)频域特征表达声音信号;原始声音信号中有环境噪声(白噪声) ,对声音进行了前期处理。模型的不足:本文考虑因素不足,未做充分分析和所有模型建立。建立出的模型也较为简单,模型不可能达到实际效果,因为很多数据和因素无法处理等,因此总体1

31、7来看模型简单,收到效果甚微。模型仍需要修正和完善。8 模型的改进与推广本文假设实验环境保持不变,但实际应用中环境是不断改变的,因此原始声音信号 y 中叠加了一定幅值的白噪声在不同时间会随着其它未经量化的因素的影响而不同。8.1 模型的改进将环境的影响关系引入,从而减少理论数据与实际数据的误差。8.2 模型的推广本文通过建立声音识别模型较为有效的模拟出正常和非正常开门的声音信号以及原始声音信号 y 中叠加了一定幅值的白噪声之间关系的数学模型,具有很好的预测和拟合效果,因此有较强的现实意义,可以在现实生活中广泛应用。9 参考文献1刘卫国主编.MATLAB 程序设计与应用.北京:高等教育出版社,2

32、006 年第二版.2姜启源 谢金星 叶俊编.数学模型.北京:高等教育出版社,2011 年第四版.3茆诗松 程依明 濮晓龙编.概率论与数理统计.北京:高等教育出版社,2011年第二版.4http:/.18程序:附录一:%正常开门声音及波形图(正 1.mat)load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集正 1.mat)sound(y,Fs)plot(y)附录二:%非正常开门声音及波形图(非 1.mat)load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集非 1.mat)sound

33、(y,Fs)plot(y)附录三:%语音信号时域频域分析y=y(:,1);%我这里假设你的声音是双声道,我只取单声道作分析,如果你想分析另外一个声道,请改成 y=y(:,2) sigLength=length(y); Y = fft(y,sigLength); Pyy = Y.* conj(Y) / sigLength; halflength=floor(sigLength/2); f=Fs*(0:halflength)/sigLength; figure;plot(f,Pyy(1:halflength+1);xlabel(Frequency(Hz); t=(0:sigLength-1)/Fs

34、; figure;plot(t,y);xlabel(Time(s);y=y(:,1);lg=length(y);yx=(0:(lg-1)*Fs/lg;19yx=yx(1:lg/2);% sound(y,fs);figure(1);subplot(211);plot(y);title(原始图像);subplot(212);yf=abs(fft(y);yf=yf(1:lg/2);plot(yx,yf);ycq=zeros(1,lg); %先将抽取后的值全设为零 for i=1:80:lg; %通过循环,每隔 80 个点将抽取后的值赋值为原函数的 ycq(i)=y(i); %采样值 end soun

35、d(ycq,Fs);figure(2);subplot(211);plot(ycq);title(提取特征向量后的图像);subplot(212);ycqf=abs(fft(ycq);ycqf=ycqf(1:lg/2);plot(yx,ycqf);附录四:%参数优化Inputp=;Inputt=;%神经网络训练输入、输出矩阵pn,meanp,stdp,tn,meant,stdt=prestd(p,t);%神经网络输入、输出数据归化处理net=newff(minmax(pn),4101,logsiglogsigpurelin,trainlm);%建立一个三层 BP 网络net.trainpara

36、m.show=25;net.trainparam.e2pochs=300;net.trainparam.goal=1e-10;%设置训练参数net,tr,epochs=train(net,pn,t);%神经网络训练a=sim(net,pn);%装载训练结果20m,b,r=postreg(a,t);%线性回归分析save(result,net);%保存网络结果在上述程序中,隐含层的神经元数可通过试凑法得到附录五:%原始声音信号 y 中叠加了一定幅值的白噪声图(正 1.mat)load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集正 1.mat)y1=y+(-0.15+0.3*rand(size(y)*max(y);plot(y1)附录六:%原始声音信号 y 中叠加了一定幅值的白噪声图(非 1.mat)load(C:UsersAdministratorDesktop2014 校级数模竞赛赛题A 题开锁声音采集非 1.mat)y1=y+(-0.15+0.3*rand(size(y)*max(y);plot(y1)

展开阅读全文
相关资源
猜你喜欢
相关搜索
资源标签

当前位置:首页 > 企业管理 > 管理学资料

本站链接:文库   一言   我酷   合作


客服QQ:2549714901微博号:道客多多官方知乎号:道客多多

经营许可证编号: 粤ICP备2021046453号世界地图

道客多多©版权所有2020-2025营业执照举报