plomb
Lomb-Scargle周期图
语法
描述
例子
不规则采样信号和信号丢失的样品
Lomb-Scargle方法可以处理信号,采样不均匀或有缺失的样本。
生成一个双通道正弦信号采样1 kHz大约0.5秒。正弦信号的频率是175 Hz, 400 Hz。嵌入在白噪声信号的方差 。
Fs = 1000;f0 = 175;f1 = 400;t = 0:1 / Fs: 0.5;wgn = randn(长度(t), 2) / 2;sigOrig =罪(2 *π* (f0; f1) * t)”+ wgn;
计算信号的周期图和阴谋。使用周期图使用默认设置。
周期图(sigOrig [] [], Fs) axisLim =轴;标题(“周期图”)
使用plomb使用默认设置来估计和情节的PSD信号。使用轴限制从前面的情节。
plomb (sigOrig t)轴(axisLim)标题(“Lomb-Scargle”)
假设信号丢失10%的原始样本。的地方南在随机位置模拟丢失的数据点。使用plomb估计和情节的PSD信号丢失的样品。
sinMiss = sigOrig;misfrac = 0.1;nTime =长度(t) * 2;sinMiss (randperm (nTime圆(nTime * misfrac))) =南;plomb (sinMiss t)轴(axisLim)标题(缺失的样本)
样品原始信号,但通过添加抖动使采样非均匀(不确定性)的时间测量。第一个即时仍然为零。使用plomb估计和PSD的非均匀采样信号的阴谋。
tirr = t + (1/2-rand(大小(t))) / Fs / 2;tirr (1) = 0;sinIrreg =罪(2 *π* (f0; f1) * tirr)”+ wgn;plomb (sinIrreg tirr)轴(axisLim)标题(非均匀采样的)
周期图数据集与失踪的样品
伽利略发现木星的四大卫星的运动在1610年的冬天。当天气允许,伽利略卫星的位置记录。使用他的观察来估计的一个卫星的轨道周期,木卫四。
Callisto弧角位置测量在分钟。缺失的数据由于多云天气指定使用nan。第一个观察是1月15日。生成一个datetime数组的观察时间。
yg =[10.5南南南南-5.5 -10.0 -12.0 10.5 11.5 -11.5 -12.0 -7.5…南南南南8.5 -12.5 -6.0 - -11.5 12.5 12.5 10.5南南南…-12.5 -10.5 -6.5 10.5 13.5 2.0 8.5 10.5南南南南南南…-8.5 -10.5 -10.5 -10.0 -8.0];obsv = datetime(1610, 14 +(1:长度(yg)));情节(yg、“o”甘氨胆酸)ax =;晚上= 18 32 46 [1];斧子。XTick =晚上;斧子。XTickLabel = char (obsv(晚上));网格
估计数据的功率谱plomb。指定一个过采样10倍。用逆天表达产生的频率。
[pxx f] = plomb (yg, obsv [] 10,“权力”);f = f * 86400;
使用findpeaks确定的位置只有主峰的光谱。功率谱以及峰值的阴谋。
[pk, f0] = findpeaks (pxx f“MinPeakHeight”10);情节(f, pxx f0、pk“o”)包含(频率(天^ {1}))标题(功率谱和突出的峰网格)
决定Callisto的轨道周期(天)的逆频率的最大能量。不同的结果从美国宇航局发布的值小于1%。
= 1 / f0时期
时间= 16.6454
美国国家航空航天局(NASA) = 16.6890184;PercentDiscrep = (Period-NASA) /美国国家航空航天局(NASA) * 100
PercentDiscrep = -0.2613
周期图与不规则采样数据集
伽利略发现木星的四大卫星在1610年1月,记录其位置的每一个晴朗的夜晚,直到3月。使用伽利略的数据发现木卫四的轨道周期,最外层的四个卫星。
伽利略的观测Callisto弧角位置在分钟。由于多云天气有一些差距。生成一个持续时间数组的观察时间。
t =[0 2 3 7 8 9 10 11 12 17 18 19 20 24日25日26日27 28 29日31日32 33 35 37…41 42 43 44 45]”;td =天(t);yg = [10.5 11.5 10.5 -5.5 -10.0 -12.0 -11.5 -12.0 -7.5 8.5 12.5 12.5…10.5 -6.0 -11.5 -12.5 -12.5 -10.5 -6.5 2.0 8.5 10.5 13.5 10.5 -8.5…-10.5 -10.5 -10.0 -8.0];情节(td、yg、“o”)
使用plomb计算数据的周期图。估计功率谱的频率
。指定一个过采样10倍。选择标准Lomb-Scargle正常化。
(1)一天=秒(天);[pxx f] = plomb (yg, td, 0.5 /一天10“归一化”);f = f *一天;
周期图有一个清晰的最大值。名字峰值频率 。周期图和注释峰值的阴谋。
(pmax lmax] = max (pxx);f0 = f (lmax);情节(pmax f pxx f0,“o”)包含(频率(天^ {1}))
使用线性最小二乘法适合数据的函数形式
拟合参数是振幅一个,B,C。
英国《金融时报》= 2 *π* f0 * t;ABC =[(大小(英尺))因为(英尺)罪(英尺)]\ yg
美国广播公司(ABC) =3×10.4243 10.4444 6.6137
使用拟合参数构建拟合函数在200点区间。数据和叠加适合的阴谋。
x = linspace (t (1), t(结束),200)';fx = 2 *π* f0 * x;y =[(大小(fx)因为(fx)罪(fx)] *美国广播公司(ABC);情节(td、yg、“o”天(x), y)
不规则采样和混叠
样品0.8 Hz正弦信号在1赫兹100年代。嵌入白噪声的正弦信号方差的1/100。重置的随机数字生成器可重复的结果。
f0 = 0.8;rng默认的wgn = randn (1100) / 10;ts = 1:10 0;s =罪(2 *π* f0 * ts) + wgn;
计算和绘制功率谱密度估计采样率。指定一个过采样10倍。
plomb(年代,ts 1 10)
有混叠因为正弦信号的频率高于奈奎斯特频率。
重复计算,但现在样品正弦信号的任意时间。包括1赫兹的频率。指定一个过采样因子2。绘制PSD。
tn =(100 *兰德(1100));n =罪(2 *π* f0 * tn) + wgn;外国资产控制办公室= 2;plomb (n, tn, 1, ofac)
混叠消失。不规则采样增加有效采样率缩减一些时间间隔。
放大的频率约为0.8赫兹。使用细网格间距为0.001赫兹。你不能指定一个过采样因素或最大频率。
df = 0.001;fvec = 0.7: df: 0.9;持有在plomb (n, tn, fvec)传说(外国资产控制办公室= 2的,“df = 0.001”)
指数分布
生成 样品白噪声的方差 鉴于1 Hz的采样率。计算出白噪声的功率谱。选择Lomb-Scargle正常化并指定一个过采样因素 。重置的随机数字生成器可重复的结果。
rng默认的N = 1024;t = (1: N) ';wgn = randn (N, 1);外国资产控制办公室= 15;[pwgn f] = plomb (wgn t[],外国资产控制办公室,“归一化”);
验证Lomb-Scargle白噪声的功率谱估计有一个指数分布与单位的意思。情节的值的直方图pwgn直方图的一组指数分布的随机数生成使用分布函数
。规范化的直方图,回想一下,周期图样本总数
。指定一个宽度为0.25。覆盖一个阴谋的理论分布函数。
dx = 0.25;br = 0: dx: 5;Nf = N * ofac / 2;hpwgn = histcounts (pwgn br) ';hRand = histcounts(日志(兰德(Nf, 1)), br) ';弯曲= br (1: end-1);栏(弯曲,[hpwgn hRand] / Nf / dx,“histc”)举行在情节(br、exp (br))传说(“wgn”,“经验pdf”,“理论pdf”)举行从
嵌入在噪声的正弦信号频率0.1赫兹。使用的信噪比 。指定正弦信号的振幅, ,使用的关系 。计算信号的功率谱和情节的直方图与实证和理论分布函数。
信噪比= 0.01;x0 =√2 *信噪比);sigsmall = wgn + x0 * sin(2 * 10π/ * t);[psigsmall f] = plomb (sigsmall t[],外国资产控制办公室,“归一化”);hpsigsmall = histcounts (psigsmall br) ';栏(弯曲,[hpsigsmall hRand] / Nf / dx,“histc”)举行在情节(br、exp (br))传说(“sigsmall”,“经验pdf”,“理论pdf”)举行从
重复计算使用 。现在的分布明显不同。
信噪比= 1;x0 =√2 *信噪比);siglarge = wgn + x0 * sin(2 * 10π/ * t);[psiglarge f] = plomb (siglarge t[],外国资产控制办公室,“归一化”);hpsiglarge = histcounts (psiglarge br) ';栏(弯曲,[hpsiglarge hRand] / Nf / dx,“histc”)举行在情节(br、exp (br))传说(“siglarge”,“经验pdf”,“理论pdf”)举行从
虚警概率
产生100个样本的一个正弦信号的采样率1 Hz。指定0.75的振幅和频率 。0.902信号嵌入白噪声方差。重置的随机数字生成器可重复的结果。
rng默认的X0 = 0.75;f0 = 0.6;vr = 0.902;Nsamp = 100;t = 1: Nsamp;X = X0 * cos (f0 * (1: Nsamp)) + randn (Nsamp) * sqrt (vr);
抛弃随机样本的10%。画出信号。
X (randperm (Nsamp Nsamp / 10)) =南;情节(t X,“o”)
计算归一化功率谱和情节。注释的水平对应的虚警概率为50%,10%,1%,0.01%。如果您生成许多90 -样本方差0.902白噪声信号,然后其中一半有一个或多个峰值高于50%,10%有一个或多个峰值高于10%,等等。
Pfa = (50 10 1 0.01) / 100;Pd = 1-Pfa;[pxx f,甲状旁腺素]= plomb (X, 1: Nsamp,“归一化”,“Pd”,Pd);pxx情节(f, f,甲状旁腺素*的(大小(f)))包含(“f”)文本(0.3 *[1 1 1 1],甲状旁腺素闲置,[repmat (“P_ {fa} = ',[4])num2str (Pfa)])
在这种情况下,峰值足够高,只有0.01%的可能信号可以实现它。
使用plomb没有输出参数重复计算。现在的情节是对数,水平的检测概率。
plomb (X, 1: Nsamp,“归一化”,“Pd”,Pd)
Lomb-Scargle周期图的嘈杂的正弦曲线
当给定一个数据向量作为唯一的输入,plomb使用归一化频率估计功率谱密度。
产生128的样本归一化频率的正弦信号 rad /样本嵌入在高斯白噪声方差的1/100。
t = (0:127)”;x =罪(2 *π* t / 4);x = x + randn(大小(x)) / 10;
估计使用Lomb-Scargle PSD的过程。重复的计算周期图。
[p f] = plomb (x);(pp、fp) =周期图(x);
情节的PSD估计分贝。验证结果是等价的。
情节(f /π,pow2db (p))在情节(fp /π,pow2db (pp))轴([0 1 -40年20])包含(“ω\ / \π”)ylabel (PSD的)传说(“plomb”,“周期图”)
估计一个三通道的Lomb-Scargle PSD信号正弦曲线组成。指定的频率
rad /样本,
rad /样品,
rad /样品。添加高斯白噪声方差的1/100。使用plomb没有输出参数计算和绘制PSD估计在分贝。
x3 = (sin(2 *π* t / 3)罪(2 *π* t / 4)罪(2 *π* t / 5)];x3 = x3 + randn(大小(x3)) / 10;图plomb (x3)
再次计算PSD估计,但现在删除随机数据的25%。
x3 (randperm(元素个数(x3), 0.25 *元素个数(x3))) =南;plomb (x3)
信号的功率谱密度与失踪的样品
如果你没有时间向量,使用南的指定失踪样本信号。
产生1024的样本归一化频率的正弦信号
rad /样本嵌入在白噪声方差的1/100。估计使用Lomb-Scargle过程的功率谱密度。使用plomb没有输出参数估计。
t = (0:1023)”;x =罪(2 *π* t / 5);x = x + randn(大小(x)) / 10;[pxx f] = plomb (x);plomb (x)
删除其他样本通过分配南“年代使用。plomb计算和绘制PSD。相同的周期图的峰值频率,因为时间轴不变。
xnew = x;xnew(2:2:结束)=南;plomb (xnew)
删除其他将采样样本。函数现在估计的周期性是原来的两倍的频率。这可能不是你想要的结果。
xdec = x(1:2:结束);plomb (xdec)
输入参数
输出参数
更多关于
引用
[1]霍恩,詹姆斯·H。,Sallie L. Baliunas. "A Prescription for Period Analysis of Unevenly Sampled Time Series."天体物理学杂志》上。302卷,1986年,页757 - 763。
[2]随着,尼古拉斯·R。“最小二乘频率分析不均匀地间隔的数据。”天体物理学和空间科学。39卷,1976年,页447 - 462。
[3]出版社,威廉H。,George B. Rybicki. "Fast Algorithm for Spectral Analysis of Unevenly Sampled Data."天体物理学杂志》上。338卷,1989年,页277 - 280。
[4]Scargle, Jeffrey D。“研究天文时间序列分析。二世。光谱分析的统计方面的不均匀间隔的数据。”天体物理学杂志》上。263卷,1982年,页835 - 853。
扩展功能
版本历史
介绍了R2014b
