《随机波浪结课作业》由会员分享，可在线阅读，更多相关《随机波浪结课作业（11页珍藏版）》请在金锄头文库上搜索。
1、随机海浪理论及应用期末结课作业 姓 名： XXX 学 号： XXXXXXXXX 年级专业： XXXXXXXXX 二一八年 六月基于MATLAB的波面数据处理高志康 21170931263 水利工程本次试验造波机工作待波浪稳定后，波高仪采集120s的数据。关于试验数据的采集，波浪模型试验规程规定：试验的波浪数据采集时间间隔应小于波浪平均波周期的1/20，在波浪平稳的条件下连续采集波浪的个数应不少于10个，因此，数据的采样间隔取 0.05 s。 1. 上跨零点法1.1 原理在采集波高时常常用到上跨零点法，所谓上跨零点法就是以平均水位为零线，把波面上升与零线相交的点作为起点。波形不规则地振动降到零线
2、以下，接着又上升再次与零线相交，这一点作为该波的终点（也是下一个波的起点）。如横坐标是时间，则两个连续上跨零点间的间距就是波的周期，若横坐标是距离，则此间距是这个波的波长。把这两点间的波峰最高点到波谷最低点的垂直距离定义为波高。对于中间可能存在的小波动，只要不与零线相交就不予考虑。1.2 程序实现%不规则波clc;clear;A=textread(waveheight.txt);A=A;M=sum(A);pingjun=M/2401;A=A-pingjun;k=0;m=0;n=0;%求临界点for i=1:2400 if A(i)=0 k=k+1; B(k)=i+1; endendwave_m
3、ax=zeros(1,k-1);wave_min=zeros(1,k-1);wave_height= zeros(1,k-1);Q_3=ceil(k-1)/3); %求前三分之一波高数Q_10=ceil(k-1)/10); %求前十分之一波高数P_10=ceil(k-1)/100); %求前1%波高数P_13=ceil(13*(k-1)/100); %求前13%波高数%求波高for j=1:length(B)-1 m=B(j); n=B(j+1); wave_max(j)=max(A(m:n); wave_min(j)=min(A(m:n); wave_height(j)= wave_max(
4、j)-wave_min(j);end%排序并取值C=sort(wave_height,descend);D=C(1:Q_3);E=C(1:Q_10);F=C(1:P_1);G=C(1:P_13);Height_Q3=mean(D)Height_Q10=mean(E)Height_P1=min(F)Height_P13=min(G)Height_Max=max(D)1.3 结果输出单位：cmH1/3H1/10H1%H13%max7.09479.123311.09007.230011.30002 最大熵法估计谱2.1 原理最大熵谱法克服了协方差函数法和FFT法在采样时域长度T、采样点数受到的限制,
5、 使谱分析的旁瓣效应大为减轻, 便于寻找卓越频率成分, 对短时域的信号分析也可得到较好的效果。本法为频谱分析法中又提供一种有效的工具。熵在信息论中是衡量随机事件不定性程度的量。最大熵法只利用资料中已有的自相关函数的信息外推后求功率谱, 其中使谱的熵函数达到最大, 除此之外并不附加任何先验信息,而能使原始信息的不定性程度最大。当前通用的谱估计法有协方差函数法和快速傅里叶变换法。前者将协方差函数作傅里叶变换以求谱。熵在信息论中是衡量随机事件不定性程度的量，最大熵法的主要思路是把自相关函数外推至无穷后，再作频域变换。对于个随机事件，若每个事件出现的概率为（）,则个事件在相互独立条件下同时发生的熵为
6、如为连续型随机变量时，其熵为 式中，为的概率密度函数。若遵从数学期望为0，方差为的正态分布时 从信息论可知，当随机事件是以等概率可能性出现时，其熵值达到最大。对于连续型随机变量，其概率密度函数为正态分布时，熵值达最大。从上式可见，与有关，而又与功率谱相关联，故熵为 如离散采样时局=1时，随机变量的一组现实值为,其熵值达到最大时谱值也是最大，称之为最大熵谱。已知功率谱与自相关函数的关系为 即 按照拉格朗日乘子法，就可得到=1的最大熵谱表达式 2.2 程序实现clear allclcclear;data=textread(waveheight.txt);a,b=size(data)pinjun=s
7、um(data)/(a*b);data=(data-pinjun);S,w=pmem(data/1000,160,512,50); %采用最大熵法，采用滤波器阶数160，估计功率谱plot(w,S); %绘制功率谱xlabel(频率/Hz);ylabel(功率谱/c.s/rad);title(最大熵法 Order=512功率谱);3.3 功率谱图3 周期图法3.1 原理周期图法是根据各态历经的随机过程功率谱的定义进行的普估计。周期图法是把随机序列x(n)的N个观测数据视为一能量有限的序列，直接计算x(n)的离散傅立叶变换，得x(k)，然后再取其幅值的平方，并除以N，作为序列x(n)真实功率谱的
8、估计。3.2 程序实现%周期图法clc;clear;wave=textread(waveheight.txt);L=length(wave);Fs=20; %采样频率T=1/20;NFFT=2(nextpow2(L)-1);y=fft(wave,NFFT)/(NFFT/2);f=Fs*(0:NFFT/2)/NFFT; %只画频率小于折叠频率的部分zhenfu=abs(y(1:NFFT/2+1); %振幅Sw=(zhenfu.2)/2; %能量sigma1=sum(Sw(2:end); ff=(2*pi*f).2;Sw2=ff.*Sw;sigma2=sum(Sw2(2:end);%谱的二阶矩pl
9、ot(f(2:end),smooth(Sw(2:end)/(f(2)-f(1),20,lowess)xlabel(w)ylabel(S(w)title(periodogram)3.3 功率谱图4 相关函数法4.1 原理间接法先由序列x(n)估计出自相关函数R(n)，然后对R(n)进行傅立叶变换，便得到x(n)的功率谱估计。4.2 程序实现%相关函数法clc;clear;wave=textread(waveheight.txt);L=length(wave);Fs=20; %采样频率Nfft=2(nextpow2(L)-1);xn=wave;cxn=xcorr(xn,unbiased); %计算序列的自相关函数%求PSDcxk=fft(cxn,Nfft);pxx=abs(cxk)/2;%绘制图形index=0:round(Nfft/2-1);k=index*Fs/Nfft;plot_pxx=(pxx(index+1)/(2*pi);plot(k,plot_pxx)xlabel(w)ylabel(S(w)title(BT)4.3 功率谱
