switch wtype
case 1 , [nul,psi,xval] = wavefun(wname,iter);
case 2 , [nul,psi,nul,nul,xval] = wavefun(wname,iter);
case 3 , [nul,psi,xval] = wavefun(wname,iter);
case 4 , [psi,xval] = wavefun(wname,iter);
case 5 , [psi,xval] = wavefun(wname,iter);
end
T = max(xval)-min(xval);
% T is the size of the domain of psi.
n = length(psi);
psi = psi-mean(psi);
% psi is numerically centered.
psiFT = fft(psi);
% computation of the modulus
sp = (abs(psiFT));
% of the FT.
% Compute arg max of the modulus of the FT (center frequency).
[vmax,indmax] = max(sp);
if indmax > n/2
indmax = n-indmax+2;
% indmax is always >= 2.
end
per = T/(indmax-1);
% period corresponding to the maximum.
freq = 1/per;
% associated frequency.
| 计算各种小波函数的“波形”,给出“采样”的横纵坐标(当然了,计算机只能处理有限的离散点),分别为xval、psi。
采样点数为2^iter个,iter默认取8,即默认采256点
计算时域支撑宽度,即横坐标宽度
小波函数的采样点数
函数值去“直流”
小波函数做离散傅里叶变换
计算变换后有限复序列各元素的模
计算模最大的点对应的下标
考虑到实信号FFT结果是双边谱,必须考虑下标在负频率部分的情形,此时减去indmax才是对应的正频率。这是计算机计算精度有限造成的,理论上实小波函数的FFT应该是严格对称的。所以模最大点在正负频率各一个;但计算机实现必然只有一个,可能在正也可能在负。这里也可以看出MATLAB的算法是十分可靠细致的,考虑到了方方面面。
用模最大点对应的“频率”作为信号频率;注意该函数里的freq还是归一化的频率。
这也是为什么在scalfrq中要除以一个采样间隔。
|