Transcript HarmonicAnalysis

NRCEST
國立交通大學土木工程學系
中華民國一百零一年七月一十七日
調和分析之介紹
調和分析是以牛頓所提之平衡潮(equilibrium tide)理論
為基礎，將潮汐表示為各分潮的線性相加，利用最小二
乘法讓實際潮位與調和潮位的誤差值最小而求得分潮的
待定係數，其表示式如下公式
M
Y (t )  a 0   [ai cos(i t )  bi sin(i t )]
i 1
M
 a0   Ai cos(i t   i )
i 1
Algorithm AND programming?
調和分析與Fourier series之差異?
2
能譜分析
TC-2009-tide
AP-2009-tide
180
160
20
140
120
|Sf(y)|
|Sf(y)|
15
100
80
10
60
40
5
20
0
0
0
12
24
36
48
60
0
12
24
36
48
T(hr)
T(hr)
KH-2009-tide
安平
台中
18
16
14
|Sf(y)|
12
10
8
6
高雄
4
2
0
0
12
24
36
T(hr)
48
60
60
調和分析之介紹
Emery and Thomson(2001)提出三個調和分析選用分潮
的限制：
1. 最 低 頻 率 分 潮 的 週 期 需 小 於 或 等 於 資 料 長 度 ， 即
Tmax≦NΔt，其中Tmax為分潮的最長週期，N為資料數，
Δt 為 取 樣 間 隔 。 其 中 最 低 角 頻 率 分 潮 (Sa) 的 週 期 為
365.24天。
2. 最高頻率分潮的週期需大於或等於兩倍的取樣間隔，
即 Tmin≧2 Δt ，其中Tmin為分潮的最短週期。最高頻率
分潮(MSK6)之週期為4.04小時。
3. 兩分潮的頻率差Δf的倒數需小於或等於資料長度，即
1/Δf≦NΔt。60分潮有部分兩相臨頻率差之倒數超過
365天者，即為Sa，Ssa； π 1，P1，S1，K1，ψ 1，φ 1；T2，
S2，R2，K2之間的分潮。
4
混疊(aliasing)
取樣定理:取樣頻率至少為訊號最大頻率的兩倍或兩
倍以上。
取樣頻率(sampling frequency): 每秒共取樣幾次。
取樣週期(sampling period):每次取樣所間隔的時間。
快速傅立葉轉換(FFT)限制:
1.滿足取樣定理。
2.取樣點數為2的n次方個。
3.訊號為周期性。
5
y  cos(2 * 25t )  cos(2 * 200t )
Fs=500 Hz
T=1/500 s
L=100
cos(2*pi*25*t)+cos(2*pi*200*t)
cos(2*pi*25*t)+cos(2*pi*200*t)
2
1
0
0.5
-2
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0
0.1
0
50
150
200
250
200
250
200
250
1
Amplitude
0
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0.1
0.5
0
0
50
1
0
0.5
0
0.01
0.02
0.03
0.04
0.05 0.06
Fs=500
100
150
cos(2*pi*200*t)
cos(2*pi*200*t)
2
-2
100
cos(2*pi*25*t)
cos(2*pi*25*t)
2
-2
滿足取樣定理
0.07
0.08
時間序列圖(time series)
0.09
0.1
0
0
50
150
100
Frequency(Fs=500)
頻譜圖(spectrum)
未滿足取樣定理，故無法完整描述訊號
cos(2*pi*25*t)+cos(2*pi*200*t)
cos(2*pi*25*t)+cos(2*pi*200*t)
2
1
0
0.5
-2
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
Fs=300 Hz
T=1/300 s
L=100
0.08
0.09
0
0.1
0
50
cos(2*pi*25*t)
100
150
100
150
Amplitude
1
0
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0.1
0.5
0
0
50
cos(2*pi*200*t)
cos(2*pi*200*t)
2
1
0
0.5
-2
150
cos(2*pi*25*t)
2
-2
100
0
0.01
0.02
0.03
0.04
0.05 0.06
Fs=300
0.07
0.08
時間序列圖(time series)
0.09
0.1
0
0
50
Frequency(Fs=300)
頻譜圖(spectrum)
Fs=200 Hz
T=1/200 s
L=100
cos(2*pi*25*t)+cos(2*pi*200*t)
cos(2*pi*25*t)+cos(2*pi*200*t)
2
2
0
1
-2
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0
0.1
0
10
20
30
60
70
80
90
100
70
80
90
100
70
80
90
100
1
Amplitude
0
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0.1
0.5
0
0
10
20
30
2
2
0
1
0
0.01
0.02
0.03
0.04
0.05 0.06
Fs=200
40
50
60
cos(2*pi*200*t)
cos(2*pi*200*t)
-2
50
cos(2*pi*25*t)
cos(2*pi*25*t)
2
-2
40
0.07
0.08
時間序列圖(time series)
0.09
0.1
0
0
10
20
30
40
50
60
Frequency(Fs=200)
頻譜圖(spectrum)
Fs=40Hz
T=1/40 s
L=100
cos(2*pi*25*t)+cos(2*pi*200*t)
cos(2*pi*25*t)+cos(2*pi*200*t)
2
2
0
1
-2
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0
0.5
0
10
20
30
cos(2*pi*25*t)
60
70
80
90
100
14
16
18
20
14
16
18
20
Amplitude
2
0
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
1
0
0
2
4
6
cos(2*pi*200*t)
2
0
1
0
0.05
0.1
0.15
0.2
0.25
Fs=40
0.3
8
10
12
cos(2*pi*200*t)
2
-2
50
cos(2*pi*25*t)
2
-2
40
0.35
0.4
時間序列圖(time series)
0.45
0.5
0
0
2
4
6
12
10
8
Frequency(Fs=40)
頻譜圖(spectrum)
fs=500
y(t)
2
0
-2
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0.1
0.06
0.07
0.08
0.09
0.1
0.06
0.07
0.08
0.09
0.1
fs=300
y(t)
2
0
-2
0
0.01
0.02
0.03
0.04
0.05
fs=200
y(t)
2
0
-2
0
0.01
0.02
0.03
0.04
0.05
x  cos(2 *100t )  0.5 cos(2 * 400t )  noise
Time series
STFT
dt=0.0017
1
1
1
0
-2
y(t)
2
-1
0
-1
0
-2
0.06
0
-1
0
-2
0.06
0.8
0.8
0.7
0.7
0.7
0.6
0.6
0.6
0.5
0.5
0.5
0.4
0.3
Amplitude
0.8
Amplitude
Amplitude
dt=0.0029
2
y(t)
y(t)
dt=0.001
2
0.4
0.3
0.2
0.1
0.1
0.1
100 200 300 400
Frequency (Hz)
500
0
0
100 200 300 400
Frequency (Hz)
500
spectrum
0.3
0.2
0
0.06
0.4
0.2
0
0
0
0
100 200 300 400
Frequency (Hz)
500
取樣間隔之影響
200
y(cm)
100
0
高雄
-100
-200
dt
σ
R2
1hr 7.4257 0.9375
6min 7.4325 0.9374
t(hr)
30
20
y(cm)
10
0
-10
-20
-30
0
3000
6000
9000
t(hr)
 'M2' 'K1' 'O1' 'Sa' 'P1' 'N2' 'S2' 'Q1' 'μ2' 'K2'
21.817 21.026 18.408 15.774 6.242 5.409 5.143 3.830 2.475 1.845
21.802 21.019 18.407 15.790 6.236 5.421 5.186 3.827 2.455 1.860
13
研究生 李明憲
取樣間隔之影響
 中央氣象局所提供之潮汐資料有6分鐘及1小時二種的取
樣間隔，本文首先探討二種取樣間隔對調和分析之影響。
分潮
K1
O1
M2
S2
取樣間隔
基隆
淡水
竹圍
高雄
蘭嶼
梗枋
6min
21.03
19.88
24.35
19.06
14.52
15.87
1hr
21.04
19.87
24.34
19.05
14.52
15.85
差值
0.01
-0.01
-0.01
-0.01
0.00
-0.02
6min
18.49
17.14
21.78
18.59
15.15
14.10
1hr
18.49
17.15
21.77
18.56
15.16
14.09
差值
0.00
0.01
-0.01
-0.03
0.01
-0.01
6min
22.03
87.21
116.18
17.39
43.18
26.70
1hr
22.05
87.19
116.21
17.39
43.18
26.69
差值
0.02
-0.02
0.03
0.00
0.00
-0.01
6min
5.42
25.77
34.64
6.89
19.18
12.12
1hr
5.40
25.75
34.64
6.89
19.17
12.11
差值
-0.02
-0.02
0.00
0.00
-0.01
-0.01
14
研究生 李明憲
資料補遺
KH2009
200
y(cm)
100
0
-100
-200
t(hr)
30
20
y(cm)
10
0
-10
-20
-30
0
3000
6000
9000
t(hr)
高雄:缺3筆 補遺三次之STD為6.1, 6.1, 6.1, 6.1cm
若缺資料不比較時6.1, 6.1, 6.1, 6.1cm
資料補遺
KH2009
450
18
400
16
350
300
12
250
|Sf(y)|
pdf(normal)
14
200
10
8
150
6
100
4
50
2
0
-40
0
-30
-20
-10
0
10
20
30
0
12
 y(cm)
24
36
T(hr)
Dev:6.12, R2=0.950
48
60
資料補遺
AP2009
200
y(cm)
100
0
-100
-200
t(hr)
30
20
y(cm)
10
0
-10
-20
-30
0
3000
6000
9000
t(hr)
安平:缺1924筆 補遺三次之STD為11.0, 6.4, 5.6, 5.2cm
若缺資料不比較時 7.7, 5.8, 5.4, 5.3cm
資料補遺
AP2009
25
400
350
20
250
15
|Sf(y)|
pdf(normal)
300
200
10
150
100
5
50
0
-30
0
-20
-10
0
10
20
30
0
12
 y(cm)
24
36
T(hr)
Dev:5.98, R2=0.957
48
60
資料補遺
TC2009
400
y(cm)
200
0
-200
-400
t(hr)
30
20
y(cm)
10
0
-10
-20
-30
0
6000
3000
9000
t(hr)
台中:缺463筆 補遺三次之STD為30.1, 15.3, 14.3, 14.2cm
若缺資料不比較時18.6, 14.4, 14.2, 14.2cm
資料補遺
TC2009
600
160
140
120
400
100
|Sf(y)|
pdf(normal)
500
300
80
60
200
40
100
20
0
-80
-60
-40
-20
0
20
 y(cm)
40
60
80
100
120
0
0
12
24
36
T(hr)
Dev:15.60, R2=0.989
48
60
分潮的選擇
60分潮
振幅 相位
K1 19.21 -76.07
O1 18.51 23.05
M2 18.17 0.15
Sa 15.03 48.17
S2 7.04 83.76
P1 5.76 -56.38
N2 4.45 53.54
Ssa 3.90 -11.00
Q1 3.84 63.49
K2 2.08 46.42
10分潮
振幅 相位
K1 19.21 -76.09
O1 18.53 22.98
M2 18.17 0.06
Sa 15.00 48.28
S2 7.02 83.73
P1 5.75 -56.50
N2 4.36 52.71
Ssa 3.86 -11.14
Q1 3.74 63.12
K2 2.10 46.83
54分潮
振幅 相位
K1 19.21 -76.05
O1 18.49 23.03
M2 18.16 0.15
S2 7.03 83.70
P1 5.76 -56.35
N2 4.45 53.53
Saa 3.86 -10.96
Q1 3.85 63.46
K2 2.08 46.46
Mm 1.03 54.73
分潮的選擇
潮汐
資料個數 8760 (每小時) 分潮個數: 60
平均值 -0.03 (CM)
R2: 0.950
STD(ym-yp): 6.1071 (CM)
資料個數 8760 (每小時) 分潮個數: 10 (選擇60個分潮最大10個
平均值 -0.01 (CM)
R2: 0.942
STD(ym-yp): 6.5344 (CM)
刪除6個不滿足選擇資料長度之分潮(雖然只少幾筆)
天文朝能量比0.796 (此顯示濾出之天文潮不好)
平均值
-0.022
R2:
0.796
STD(ym-yp): 12.269 (高雄若刪除Sa分潮則STD變大)
資料補遺
採用內差方式(polynomial Interpolation)進行補遺。
提出具體作業化之調和內差補遺方法：HIS (Harmonic
Interpolation Supplements )。準確且有效率地一次完
成潮汐之調和分析、預報、補遺、及相關之誤差分析。
透過人造斷續缺漏的營造，應用HIS，提供純化的天文
潮位。
Effect of length of missing data
Is it necessary for HIS?
23
壓力與溫度的影響
.
0.5
KH-2009-tide
Surg Analysis
300
Diff. of (Meas.-HIS Pred.), Suppl.
24hrs MAvg of Suppl. Diff.
0.4
HAIMA
AERE
AERE
200
0.2
100
0.1
y(cm)
Diff.(Meas.-HIS) (m)
0.3
RANANIM
(24-2160)hrs MAvg of Meas.
RANANIM
24hrs MAvg of Original Diff.
0
0
-0.1
-100
-0.2
HuaLien Tides Y2K4
HL2K4M12D365L22.dat
HIS Supplement
2
-200
Measurement
1.5
Tide Ele. (m)
Error Analyses:
Meas. RMS=0.06584 (m); Corel. Coef.=0.98687
Supple. RMS=0.18084 (m); Corel. Coef.=0.93389
1
0.5
0
-300
-0.5
0
3000
6000
t(hr)
-1
5040
5400
6120 .
5760
Time (hrs)
Storm surge: 張永欣論文
Temperature+ multi-point predication:林立青博士論文
24
9