在做實驗的時候經常會遇到一些週期性現象,當實驗資料比較多的時候手動判斷週期太過麻煩,比如有一群下面這樣的沙雕陣列
最簡單粗暴的想法便是選取某個值對陣列進行等間隔取樣。如果此時的值正好是這個陣列的週期,那麼採集的資料必然相差較小,即方差較小。也就是說方差最小時的值便是該陣列的週期。
%
%獲取陣列的週期
functiont=
higett
(seq)
lenq =
length
(seq)
;tvar = realmax;
for t =
2:lenq/
4 vart =
zeros
(t-1,1
);%當前週期不同起始點的方差
for i =
1:t-
1 samp =
seq(i:t:lenq)
;%取樣序列
vart
(i)=
var(samp)
; end
ifmean
(vart)
tvar =
mean
(vart);t
= t;
endendend
結果為219,手動量取的某個值為222,效果還是不錯的,可見這種思路可行,不過速度比較慢,可能是在陣列中間隔地選取取樣點造成的。
對於samp = seq(i:t:lenq)
這條語句,其實很類似於將陣列整型成矩陣,而且,對於某一相同的取樣週期,不同初始點採集的資料個數是相同的。
則改進後的**如下
function[t
]=gett
(seq)
%獲取乙個序列的週期
lenq =
length
(seq)
;tvar =
zeros
(floor
(lenq/4)
,1);
%不同試驗週期的方差
tvar(1
)= realmax;
for t =
2:lenq/
4 cols =
floor
(lenq/t)
;%取樣點個數
mat =
reshape
(seq(1
:t*cols)
,t,cols)
;%整型為t行cols列
tvar =
mean
(var
(mat,0,
2));
%求取每一行的方差,並取平均值
end[~,
t]=min
(tvar)
;end
tictoc一下,時間從0.6秒變為0.1秒,快了不少。
然而,通過新建乙個三角函式y = sin(x/100)
來檢測該方法時卻出現了問題,得到其週期為1257,正好是2倍週期(雖說2倍週期也是週期),其tvar
的值如下圖所示
可見對於不同的陣列,還需要通過新增不同的條件來加以限定哈。
通過Matlab 使用 FFT 分析週期性資料
前言 天文學家使用蘇黎世太陽黑子相對數將幾乎 300 年的太陽黑子的數量和大小製成 對大約 1700 至 2000 年間的蘇黎世數繪圖。load sunspot.dat year sunspot 1 relnums sunspot 2 plot year,relnums xlabel year yl...
matlab陣列的建立
原文 1 直接輸入 行向量 a 1,2,3,4,5 列向量 a 1 2 3 4 5 2 用 生成向量 a j k 生成的行向量是a j,j 1,k a j d k 生成行向量a j,j d,j m d m fix k j d 3 函式linspace 用來生成資料按等差形式排列的行向量 x lins...
matlab陣列元素的引用
1.下標法 subscripts 2.索引法 index 3.布林法 boolean 在使用這三種方法之前,大家頭腦一定要清晰的記住,matlab中陣列元素是按列儲存 與fortran一樣 比如說下面的二維陣列 a 8 1 6 3 5 7 4 9 2 matlab的儲存順序是8,3,4,1,5,9,...