电子科技大学数字信号处理DSP课程设计-钢琴音符识别

.精选文本2014级数字信号处理课程设计报告题目：钢琴音符识别姓名：邱晨曦学号：2014010909008.精选文本答辩时间：2016/12/9一.题目要求：（1）播放和记录一段钢琴音乐中的音符；（2）记录到音符以后，找到音符所对应的现代标准钢琴的钢琴键，并分析结果。二.课程设计思路：（1）涉及到的知识点：快速傅里叶变换、钢琴音频信号的时域和频域的特性、能熵比的概念、频率校正、频率与音符的转换关系。（2）方案分析：A.预处理部分：1.直接用audioread函数读出来的原始数据。优点：准确率较高；缺点：数据量较大，采样频率为44kHz，远大于奈奎斯特采样率。2.以11kHz的采样率重新采样，并转换为单声道。优点：数据量小了很多，易于处理；缺点：牺牲了部分的准确率，但对于音符的判断影响可以忽略。B.端点检测算法：1.双门限法：1.计算短时能量（高门限）和过零率（低门限）；2.选取一个较高的门限2T，语音信号的能量包络大部分都在此门限之上，进行一次初判，语音起止点位于该门限与短时能量包络交点所对应的时间间隔之外；3.根据噪声能量，确定一个较低的门限1T，并从初判起点往左，从初判终点往右搜索，分别找到能零比曲线第一次与门限1T相交的两个点，两点之间段就是用双门限方法所判定的语音段；4.以短时平均过零率为准，从低门限点往左右搜索，找到短时平均过零率低于某阈值的两点，为语音的起止点；.精选文本图1：双门限法示意图说明：算法中的阀值是根据实验过程调节的。该算法在实际应用的过程中发现：在语音信号频率分布较为集中的时候，端点检测出来的结果比较准确，但当语音信号频率分布比较分散的时候，很难通过控制固定的阀值来检测到每个音符；2.自相关法：由于两种信号的自相关函数存在极大的差异，可以利用这种差别来提取语音端点。根据噪声的情况，设置两个阈值1T和2T，当相关函数最大值大于2T时，便判定是语音；当相关函数最大值大于或小于1T时，则判定为语音信号的端点。该算法同样存在当语音信号频率分布较广的时候，阀值比较难控制的问题。3.基于谱熵的端点检测：基于谱熵语音端点检测方法是通过检测谱的平坦程度，来进行语音端点检测的，为了更好进行语音端点检测，采用语音信号的短时功率谱构造语音信息谱熵，从而对语音段和噪声段进行区分。检测思路：1.对语音信号进行分帧加窗；2.计算每一帧的谱能量；3.计算出每一帧中每个样本点的概率密度函数20()()()nnnNnnlYkYkpEYl；4.计算出每一帧的谱熵值20()ln()NnnnlHplpl（由信息论知识知道，熵值在自变量服从均匀分布的时候，熵值达到最大值，所以噪声的熵值是比较大的，而钢琴音符的熵值是比较小的，由此区别了噪声和音符）；5.设置判决门限；6.根据各帧的谱熵值进行端点检测。在实验过程中发现：依然存在当语音信号频率分布较广时，阀值不太好控制的问题。因此对该方法进行改进，引入，能熵比的概念：.精选文本谱熵值类似于过零率，能熵比的表示为1nnnLEEEFH。由于噪声和信号的能熵比差别很大。因此在能熵比的图像中，每一个“尖刺”就代表了一个特定频率的语言信号。图2：能熵比图中的“尖刺”在检测过程中，依然不能通过简单的设置阀值的办法来进行端点检测，原因是语音频率分布较广时，每个音符的能熵比变化范围差别较大，如下图所示，有的“尖刺”完全在门限之上，而有的则完全在门限之下。图3：88阶全音的能熵比图因此，采用检测能熵比中的“低谷点”（该点比左右两边的一定数目的点的能熵比都小）的方法。语音信号一定位于两个低谷点之间的部分，再对低谷点进行适当的左右移动作为语音信号的起止点。如下图所示：.精选文本图4：标记起止点的能熵比图（绿色为起始点，红色为截止点）（3）设计框架和流程：1.用audioread函数读入钢琴音乐，并用sound函数播放；2.为了方便处理，对信号以11.025kHz的频率进行重新采样，并统一转换成单声道的信号；3.因为语言信号可以在短时间内认为是平稳的，因此对语音信号进行分帧的处理，设置帧长320，为了减小误差，两帧之间设置重叠部分，因此帧移取80；4.计算每一帧的能熵比；5.找到能熵比中的“低谷点”（该点比左右两边的一定数目的点的能熵比都小）；6.如果两个低谷点之间的距离大于miniL(认为持续长度超过一定长度的为音符，最小长度miniL可自行设置)。则低谷点右移sr(即shiftright,数值可自行调节)帧作为一段信号的起始点，将低谷点左移sl（即shiftleft,数值可自行调节）帧作为截止点[注：采用该方法的优点是通过调节相关参数能适应多种情况，缺点是检测环境发生较大变化时，需要重新设置参数]；7.将找到的语音段转换成未分帧时对应坐标的语音段，并对每段做快速傅里叶变换；8.找到每段快速傅里叶变换中的最大值以及最大值所对应的横坐标（fft点），将横坐标转换成相应的频率，得到的频率即为该段音符的频率；9.利用比值法进行频率的校正，窗函数选择矩形窗；10.根据检测到的频率确定音符，计算公式为：212log()49440fn，n为第几个按键，再通过查表得到对应音符；11.分析结果。三.具体设计过程：（1）部分代码（测试部分缺省）：主函数部分：[x,fs]=audioread('钢琴音频.WAV');formatshort;wlen=320;inc=80;%分帧的帧长和帧移.精选文本overlap=wlen-inc;%帧之间的重叠部分sound(x,fs);%播放音乐x=calsample(x,fs);%为了方便处理，重新以11025Hz的频率采样，并转换成单声道x=x-mean(x);%消去直流分量x=x/max(abs(x));%幅值归一化y=Enframe(wlen,inc,x)';%分帧fn=size(y,2);%取得帧数time=(0:length(x)-1)/11025;%计算时间坐标frameTime=frame2time(fn,wlen,inc,11025);%计算各帧对应的时间坐标sr=2;sl=13;miniL=33;%配置左右移动的帧数和要求的最短帧数[voicesegment,vos,Ef]=get_segment(y,fn,sr,sl,miniL);%获得语音段[real_f,ft,ax]=get_f(x,voicesegment,vos,wlen,inc);%检测频率的结果fori=1:length(real_f)real_f(i)=roundn(real_f(i),-4);endfori=1:length(real_f)real_node(i)=get_node(real_f(i))end%**********************************绘图部分********************************subplot211;stem(real_f);.精选文本title('频率检测结果');xlabel('音符/个');ylabel('频率/Hz');subplot212;stem(real_node,'r');title('音符检测结果');xlabel('音符/个');ylabel('对应按键');figure(2);subplot211;plot(Ef);title('能熵比图及语音起止点');xlabel('帧数/个');ylabel('能熵比');fori=1:length(voicesegment)text(voicesegment(i).begin,Ef(voicesegment(i).begin),'o','color','g')text(voicesegment(i).end,Ef(voicesegment(i).end),'o','color','r')endsubplot212,plot(time,x,'k');title('语音信号端点检测结果')axis([0max(time)-11]);ylabel('幅值');fork=1:vos%标出有话段nx1=voicesegment(k).begin;nx2=voicesegment(k).end;nxl=voicesegment(k).duration;fprintf('%4d%4d%4d%4d\n',k,nx1,nx2,nxl);subplot212line([frameTime(nx1)frameTime(nx1)],[-11],'color','r','linestyle','-');line([frameTime(nx2)frameTime(nx2)],[-11],'color','r','linestyle','--');.精选文本end其中的用到的子函数：1.calsample.m（调整采样率和声道）functionsample=calsample(sampledata,FS)temp_sample=resample(sampledata,1,FS/11025);%调整采样频率[~,n]=size(temp_sample);if(n==2)%转换成单声道sample=temp_sample(:,1);elsesample=temp_sample;endend2.Enframe.m（分帧函数）functionf=Enframe(len,inc,x)%对读入的语音进行分帧，len为帧长，%inc为帧重叠样点数，x为输入语音数据fh=fix(((size(x,1)-len)/inc)+1);%计算帧数f=zeros(fh,len);%设置一个零矩阵，行为帧数，列为帧长i=1;n=1;whilei=fh%帧间循环j=1;whilej=len%帧内循环f(i,j)=x(n);.精选文本j=j+1;n=n+1;endn=n-len+inc;%下一帧开始位置i=i+1;end3.frame2time.m（坐标刻度转换）functionframeTime=frame2time(frameNum,framelen,inc,fs)frameTime=(((1:frameNum)-1)*inc+framelen/2)/fs;%求对应的时间坐标4.get_segment.m（端点检测，确定音符段）function[voicesegment,vos,Ef]=get_segment(y,fn,sr,sl,miniL)ifsize(y,2)~=fn,y=y';end%把y转换为每列数据表示一帧语音信号wlen=size(y,1);%取得帧长fori=1:fnSp=abs(fft(y(:,i)));%FFT取幅值Sp=Sp(1:wlen/2+1);%只取正频率部分Esum(i)=sum(Sp.*Sp);%计算能量值prob=Sp/(sum(Sp));%计算概率H(i)=-sum(prob.*log(prob+eps));%求谱熵值endhindex=find(H0.1);H(hindex)=max(H);.精选文本Ef=sqrt(1+abs(Esum./H));%计算能熵比Ef=Ef/max(Ef);%归一化[x1,y1]=get_max(Ef);%找到能熵比中的“高峰点”[x2,y2]=get_min(Ef);%找到能熵比中的“低估点”voicesegment(1).begin=x1(1)-8;%由于仅仅靠低谷点无法检测出第一个音符起始位置，因此将高峰点的第一个值左移8帧作为第一个音符的起始点voicesegment(1).end=x2(1)-sl;%将第一个低谷点作为第一个音符的截止点voicesegment(1).duration=voicesegment(1).end-voicesegment(1).begin+1;%一个音符的持续帧数j=1;fork=2:length(x2)%将找到的低谷点作为音符的起止点ifx2(k)-x2(k-1)=miniL%剔除持续帧长度小于miniL的音符段j=j+1;temp1=x2(k-1)+sr;%将低谷点右