Ⅰ MATLAB四阶五级 Runge-Kutta-Felhberg 算法计算一阶微分方程,方程系数会随着上一步计算结果变化怎么办
题主给出的这种类型微分方程,是可以用四阶五级 Runge-Kutta-Felhberg 算法和ode45()函数。
ode45()函数实际上就是变步长四阶五级 Runge-Kutta(龙格-库塔)法。为了更好地回答,请题主给出具体的微分方程,以便说明问题。
Ⅱ 如何用matlab实现最大熵谱法
例如:
具体问题:matlab中burg法谱分析语句[P,f]=pburg(x,M,nfft,fs,range),文献对M阶、nffg,fs都有说明,但
是我不太明白。假如我的记录60秒的数据,每秒100个数据点(频率为100hz)那么对应的M、nfft、fs应该是多少
可以使用
pburg(x,50,,6000,100,‘onesided’)
Ⅲ 请问matlab中如何在功率谱密度图中添加一条置信曲线
这么多,看都懒得看!~
Ⅳ 有谁知道杜宾算法原理
matlab自带有杜宾算法的函数:levinson
有自带的burg算法函数:arburg,pburg
窗口打:doc spectrum.burg
就能查到Burg spectrum
在窗口打edit levinson
就能查看函数levinson的源代码
其它类似
Ⅳ matlab中pburg函数画出来的怎么感觉是错误的
解析:
help pburg
然后,
检查,代码中此函数调用,格式是否正确
Ⅵ 利用Burg算法计算AR模型的参数:arburg函数
(1)a=arburg(x,order);返回利用Burg算法计算得到的AR模型的参数:输入参数or-der为AR模型的阶数。
(2)[a,e]=arburg(…);当输入信号为白噪声时,同时返回最后误差e。
(3)[a,e,k]=arburg(…):同时返回反射系数k。
[例4-5]用Burg算法计算AR 模型的参数。
rand(‘seed’,0);a=[1 0.1 0.2 0.3 0.4 0.5 0.6];
x=impz(1,a,20)+randn(20,1)/20;
[a,e,k]=arburg(x,5)。
Ⅶ matlab提示 function [psdviaBurg, f, p] = myBurg(x, Fs, varargin) | Error: Function definitions
1、将下面函数复制保存为myBurg.m文件。就可以计算。预测误差E 、系数a 、误差功率psdviaBurg
function [psdviaBurg, f, p,a,E] = myBurg(x, Fs, varargin)
%MYBURG 根据burg算法实现的AR模型功率谱计算
% psdviaBurg 根据burg算法求出的功率谱值
% f 频率轴参数
% p 模型阶次
% a AR模型参数
%E AR模型误差
% x 输出信号
% Fs 采样率
% varargin 若为数值型,则为AR模型阶次
% 若为字符串,则为定阶准则,AR模型阶次由程序确定
%
% $Author: lskyp
% $Date: 2010.6.26
% 解析输入参数内容
error(nargchk(3, 3, nargin)); % 该函数的输入必须为三个个
if strcmp(class(varargin{1}), 'double')
p = varargin{1};
elseif ischar(varargin{1})
criterion = varargin{1};
else
error('参数2必须为数值型或者字符串');
end
x = x(:);
N = length(x);
% 模型参数求解
if exist('p', 'var') % p变量是否存在,存在则不需要定阶,直接使用p阶
[a, E] = computeARpara(x, p);
else % p不存在,需要定阶,定阶准则即criterion
p = ceil(N/3); % 阶次一般不超过信号长度的1/3
% 计算1到p阶的误差
[a, E] = computeARpara(x, p);
% 根据误差求解目标函数最小值
kc = 1:p + 1;
switch criterion
case 'FPE'
goalF = E.*(N + (kc + 1))./(N - (kc + 1));
case 'AIC'
goalF = N.*log(E) + 2.*kc;
end
[minF, p] = min(goalF); % p就是目标函数最小的位置,也即定阶准则给出的阶次
% 使用p阶重新求解AR模型参数
[a, E] = computeARpara(x, p);
end
% 计算功率谱密度
[h, f] = freqz(1, a, [], Fs);
psdviaBurg = E(end)*abs(h).^2./Fs;
function [a, E] = computeARpara(x, p)
% 根据信号序列x和阶次p计算AR模型参数和误差
N = length(x);
% 初始值
ef = x; % 前向预测误差
eb = x; % 后向预测误差
a = 1; % 初始模型参数
E = x'*x/N; % 初始误差
k = zeros(1, p); % 为反射系数预分配空间,提高循环速度
E = [E k]; % 为误差预分配空间,提高速度
for m = 1:p
% 根据burg算法步骤,首先计算m阶的反射系数
efm = ef(2:end); % 前一阶次的前向预测误差
ebm = eb(1:end - 1); % 前一阶次的后向预测误差
num = -2.*ebm'*efm; % 反射系数的分子项
den = efm'*efm + ebm'*ebm; % 反射系数的分母项
k(m) = num./den; % 当前阶次的反射系数
% 更新前后向预测误差
ef = efm + k(m)*ebm;
eb = ebm + conj(k(m))*efm;
% 更新模型系数a
a = [a; 0] + k(m)*[0; conj(flipud(a))];
% 当前阶次的误差功率
E(m + 1) = (1 - conj(k(m))*k(m))*E(m);
end
2、例如:
参考论坛中的例子。
randn('state', 1);
x = randn(100, 1);
y = filter(1, [1 1/2 1/3 1/4 1/5], x);
pburg(y, 4, [], 1000);
[psd, f, p,a,E] = myBurg(y, 1000, 'FPE');
figure;
a,E
plot(f, 10*log10(psd));
结果:图还是一样,就补贴出来了。
a =
1.0000
0.3116
0.3647
0.2086
0.2088
0.0425
E =
1.0224 0.9859 0.9167 0.8989 0.8644 0.8629
Ⅷ matlab中AR模型的计算
[m,refl]=ar(Y,N,Approavh,Window);
%Y观察值
%N AR模型的阶数 取 实整数
%Appproach 表示计算模型的方法 'fb','ls','yw','burg','gl'
Ⅸ 在matlab中像buttord,butter,pburg这些函数是已经定义好的吗可以直接用吗
function varargout = pburg(x,p,varargin)
%PBURG
Power Spectral Density (PSD) estimate via Burg's method.
%
Pxx = PBURG(X,ORDER) returns the PSD of a discrete-time signal vector X
%
in the vector Pxx. Pxx is the distribution of power per unit frequency.
%
The frequency is expressed in units of radians/sample. ORDER is the
%
order of the autoregressive (AR) model used to proce the PSD. PBURG
%
uses a default FFT length of 256 which determines the length of Pxx.
%
%
For real signals, PBURG returns the one-sided PSD by default; for
%
complex signals, it returns the two-sided PSD. Note that a one-sided
%
PSD contains the total power of the input signal.
%
%
Pxx = PBURG(X,ORDER,NFFT) specifies the FFT length used to calculate
%
the PSD estimates. For real X, Pxx has length (NFFT/2+1) if NFFT is
%
even, and (NFFT+1)/2 if NFFT is odd. For complex X, Pxx always has
%
length NFFT. If empty, the default NFFT is 256.
%
%
[Pxx,W] = PBURG(...) returns the vector of normalized angular
%
frequencies, W, at which the PSD is estimated. W has units of
%
radians/sample. For real signals, W spans the interval [0,Pi] when
%
NFFT is even and [0,Pi) when NFFT is odd. For complex signals, W
%
always spans the interval [0,2*Pi).
%
%
[Pxx,F] = PBURG(...,Fs) specifies a sampling frequency Fs in Hz and
%
returns the power spectral density in units of power per Hz. F is a
%
vector of frequencies, in Hz, at which the PSD is estimated. For real
%
signals, F spans the interval [0,Fs/2] when NFFT is even and [0,Fs/2)
%
when NFFT is odd. For complex signals, F always spans the interval
%
[0,Fs). If Fs is empty, [], the sampling frequency defaults to 1 Hz.
%
%
[Pxx,W] = PBURG(...,'twosided') returns the PSD over the interval
%
[0,2*Pi), and [Pxx,F] = PBURG(...,Fs,'twosided') returns the PSD over
%
the interval [0,Fs). Note that 'onesided' may be optionally specified,
%
but is only valid for real X. The string 'twosided' or 'onesided' may
%
be placed in any position in the input argument list after ORDER.
%
%
PBURG(...) with no output arguments plots the PSD estimate in dB per
%
unit frequency in the current figure window.
%
%
EXAMPLE:
%
randn('state',1); x = randn(100,1);
%
y = filter(1,[1 1/2 1/3 1/4 1/5],x); % Fourth order AR filter.
%
pburg(y,4,[],1000);
% Uses the default NFFT of 256.
%
%
See also PCOV, PMCOV, PYULEAR, PMTM, PMUSIC, PEIG, PWELCH, PERIODOGRAM,
%
ARBURG, PRONY, SPECTRUM/BURG, DSPDATA/PSD.
% Author(s): R. Losada and P. Pacheco
%
Copyright 1988-2004 The MathWorks, Inc.
%
$Revision: 1.26.4.5 $ $Date: 2004/04/13 00:18:07 $
error(nargchk(2,5,nargin))
method = 'arburg';
[Pxx,freq,msg,units,Sxx,options] = arspectra(method,x,p,varargin{:});
error(msg);
if nargout==0,
freq = {freq};
if strcmpi(units,'Hz'),
freq = {freq{:},'Fs',options.Fs};
end
hpsd = dspdata.psd(Pxx,freq{:},'SpectrumType',options.range);
% Create a spectrum object to store in the PSD object's metadata.
hspec = spectrum.burg(p);
hpsd.Metadata.setsourcespectrum(hspec);
plot(hpsd);
else
% Assign output arguments.
varargout = {Pxx,freq,Sxx}; % Pxx=PSD, Sxx=PS
end
% [EOF] pburg.m
Ⅹ matlab中现代谱估计的burg法怎么编写
matlab中,现代谱估计的很多方法都可以实现。music方法用pmusic命令实现;pburg函数利用burg法实现功率谱估计;pyulear函数利用yule-walker算法实现功率谱估计等等。