定義適應函數

我們一樣依照 (3.10) 式分別來計算點到每個圖形的距離，因此可得到輸入影像

Pattern 1 Pattern 2 Pattern k

Input point (xi , yi)

particle j: …

D1= d11*d12*…*d1k

d

11

d

12

d

1k

圖 3.7. 一個點到全部圖形的距離計算。

圖 3.8. 全部點到全部圖形的距離計算。

我們在參數空間上隨機分佈 m 個粒子，每個粒子依照欲偵測影像上參數圖 形的數量設定維度。接著依上述的方法去計算每個粒子的適應值，適應值越低，

代表粒子表示的參數圖形與欲偵測影像的圖形越相似。比較適應值，找出每一 次演化的 pBset 與 gBest，並以 PSO 演算法去調整粒子位置，直到演化次數到達 停止條件，也就是最大的演化次數 N 後，得到一組最好的粒子，並以粒子的參 數去還原圖形。整個系統如圖 3.9 所示。

圖 3.9. PSO 偵測圖形參數演算法系統。

我們將 PSO 演算法偵測圖形參數系統整理如下：

Algorithm：PSO 偵測圖形參數演算法。

Input：欲偵測影像中點的座標 (

x

_i,

y

_i)，i=1, L, n，設定要偵測的圖形數目

k

，粒子數目m，學習參數

α

、 ，終止條件：最大演化次數

β N

。 Output：適應值最佳的一個粒子，代表最佳的一組圖形參數，使偵測影像中的

點到此圖形參數所代表的圖形距離最小。

Step 1. 初始設定：依照要偵測的圖形數目 設定一個粒子 ，

k

s_i

s_i=

[ c

i_,1x

c

i_,1y

a

i_,1

b

i_,1

r

i_,1

θ

i_,1, _L,

c

i_,kx

c

i_,ky

a

i_,k

b

i_,k

r

i_,k

θ

i_,k

]

^T

同理，設定s₁, s₂, K, s_m。

Step 2. 圖形參數經過正規化處理後，輸入點座標，計算一個點 到一個

圖形的距離，如 (3.10) 式：

Chapter 4 實驗與結果

PSO 參數圖形偵測系統的實驗分為三個部分，第一部分為模擬的參數圖 形，包括混合的直線、圓、橢圓以及雙曲線。第二部分為模擬的震測訊號，我 們用已知的地層深度、震波速度、傾斜角度等，去模擬出單炸點震測圖，其中 有水平地層及傾斜地層。第三部分偵測為真實的單炸點震測圖。

實驗的程式以 Matlab R2006a 執行，作業系統為 Windows XP，Intel Pentium(R) 4，3.00GHz CPU，1.5GB RAM。

4.1 模擬參數圖形偵測

訓練的影像樣本範圍皆為 50x50，每一個圖形有 50 個點，每個點的座標加 上[-1~1]的 uniform random noise。

我們設定粒子個數為 300 個，根據圖形數目 k 去設定粒子維度，每個維度 的初始分佈皆是隨機分佈於[-20~20]，實驗中沒有限制粒子移動的大小或範圍。

演化次數設定為 300 次，學習參數

α

=

β

=2。 實驗 1. 偵測圓的實驗

輸入點為兩個圓形，範圍為 50x50，共有 100 個輸入點，每個點的座標加 上(-1~1)的 uniform random noise。設定偵測的圖形數目 k = 2。表 4.1.1 為原始 輸入參數，表 4.1.2 為參數偵測結果，圖 4.1(a) 為兩個圓的偵測結果，4.1(b) 為 gBest 的適應值 vs. iterations。

表 4.1.1. 圓的原始輸入參數。

c

x

c

_y

r

Pattern 1 17 17 10 Pattern 2 33 33 15

表 4.1.2. 圓的參數偵測結果。

c

x

c

_{y a b}

r

θ ^執行時間

(sec) Pattern 1 17.04 17.28 -0.98 -1.01 -96.53 9.31°

Pattern 2 32.67 33.15 -1.00 -0.99 -225.86 4.12°

77.14

(a)

(b)

圖 4.1. (a) 圓的偵測結果。(b) gBest 的適應值 vs. iterations。

實驗 2. 偵測橢圓的實驗

輸入點為兩個橢圓形，範圍為 50x50，共有 100 個輸入點，每個點的座標 加上(-1~1)的 uniform random noise。設定偵測的圖形數目 k = 2。表 4.2.1 為原 始輸入參數，表 4.2.2 為參數偵測結果，圖 4.2(a) 為兩個橢圓的偵測結果，4.2(b) 為 gBest 的適應值 vs. iterations。

表 4.2.1. 橢圓的原始輸入參數。

c

x

c

_{y a b}

r

θ

Pattern 1 20 35 5 25 10 -30°

Pattern 2 25 35 20 10 15 10

表 4.2.2. 橢圓的參數偵測結果。

c

x

c

_{y a b}

r

θ ^執行時間

(sec) Pattern 1 34.70 20.08 5.30 0.18 113.28 29.61

Pattern 2 18.31 38.37 0.49 2.00 201.67 -11.51

69.69

(a)

(b)

圖 4.2. (a) 橢圓的偵測結果。(b) gBest 的適應值 vs. iterations。

實驗 3. 偵測一條直線的實驗

輸入點為一條直線，範圍為 50x50，共有 50 個輸入點，每個點的座標加上 (-1~1)的 uniform random noise。設定偵測的圖形數目

k = 1。表 4.3.1 為原始輸

入參數，表 4.3.2 為參數偵測結果，圖 4.3(a) 為一條直線的偵測結果，4.3(b) 為 gBest 的適應值 vs. iterations。

表 4.3.1. 一條直線的原始輸入參數。

Pattern

x

+ y=50

表 4.3.2. 一條直線的參數偵測結果。

c

x

c

_{y a b}

r

θ ^{執行時間(sec)}

Pattern 28.57 21.49 0.59 -1.68 0.04 14.13 49.69

(a)

(b)

圖 4.3. (a) 一條直線的偵測結果。(b) gBest 的適應值 vs. iterations。

實驗 4. 偵測兩條直線的實驗

輸入點為兩條直線，範圍為 50x50，共有 100 個輸入點，每個點的座標加 上(-1~1)的 uniform random noise。設定偵測的圖形數目 k =1。表 4.4.1 為原始輸 入參數，表 4.4.2 為參數偵測結果，圖 4.4(a) 為兩條直線的偵測結果，4.4(b) 為 gBest 的適應值 vs. iterations。

表 4.4.1. 兩條直線的原始輸入參數。

Pattern 1

x

+ y=50

Pattern 2

x

− y=0

表 4.4.2. 兩條直線的參數偵測結果。

c

x

c

_{y a b}

r

θ 執行時間(sec)

Pattern 25.12 24.66 0.99 -1.00 4.16 0.26 47.39

(a)

(b)

圖 4.4. (a) 兩條直線的偵測結果。(b) gBest 的適應值 vs. iterations。

實驗 5. 偵測兩條直線的實驗

輸入點為兩條直線，範圍為 50x50，共有 100 個輸入點，每個點的座標加 上(-1~1)的 uniform random noise。設定偵測的圖形數目 k =2。表 4.5.1 為原始輸 入參數，表 4.5.2 為參數偵測結果，圖 4.5(a) 為兩條直線的偵測結果，4.5(b) 為 gBest 的適應值 vs. iterations。

表 4.5.1. 兩條直線的原始輸入參數。

Pattern 1 x=25 Pattern 2

y

=25

表 4.5.2. 兩條直線的參數偵測結果。

c

x

c

_{y a b}

r

θ 執行時間(sec)

Pattern 1 24.98 22.74 3.06 -0.32 0.06 -18.55°

Pattern 2 23.73 25.19 1.55 -0.64 -3.42 -57.43°

76.35

(a)

(b)

圖 4.5. (a) 兩條直線的偵測結果。(b) gBest 的適應值 vs. iterations。

實驗 6. 偵測兩個雙曲線的實驗

輸入點為兩個雙曲線圖形，範圍皆為 50x50，共有 100 個輸入點，每個點 的座標加上(-1~1)的 uniform random noise。設定偵測的圖形數目 k = 2。表 4.6.1 為原始輸入參數，表 4.6.2 為參數偵測結果，圖 4.6(a) 為兩條直線的偵測結果，

4.6(b)為

gBest 的適應值 vs. iterations。

表 4.6.1. 兩個雙曲線的原始輸入參數。

c

x

c

_{y a b r} θ

Pattern 1 25 25 1 -4 10 0°

Pattern 2 25 25 -4 1 10 0°

表 4.6.2. 兩個雙曲線的參數偵測結果。

c

x

c

_{y a b}

r

θ ^{執行時間(sec)}

Pattern 1 25.07 25.01 1.99 -0.50 -2.09 0.06°

Pattern 2 25.05 25.04 0.49 -2.01 1.93 -0.12°

71.38

(a)

(b)

圖 4.6. (a) 兩個雙曲線的偵測結果。(b) gBest 的適應值 vs. iterations。

4.2 模擬單炸點震測訊號圖形偵測

我們將 PSO 偵測系統應用在震測訊號圖形偵測上，首先模擬出單炸點震測 訊號圖，接著將訊號作 Envelope [29]處理，接著作 thresholding，取出震幅大於 threshold 的點，取出來的點便作為 PSO 偵測系統的輸入點，偵測出參數後，依 據參數重建圖形，整個流程如圖 4.7。

圖 4.7. PSO 單炸點震測訊號圖形偵測系統。

4.2.1 偵測水平地層單炸點震測訊號

水平地層單炸點震測訊號如圖 4.8 所示 [30]，共有兩層地層，地面上 Source 為爆炸點，在爆炸點左右各有 64 個接收站，加上爆炸點處 1 個，共有 129 個接收站，每個接收站之間距離Δx為 50 meters，接收站與爆炸點之間距離

， 為接收站編號。接收站之 sampling interval 為 0.004 sec，總共有 512 個 samples。地表與反射層的垂直深度 為 2,500 meters，為爆炸點到地下地 層的垂直深度。直接波與反射波速度 皆為 2,500 meters/sec，反射係數為 0.2。

由圖 4.8 可知直接波及反射波的行進方向，由計算可知直接波訊號為直線，(4.1) Ricker Wavelet 作 convolution，得到 impulse 訊號， 再加上 random white band 10.2539 Hz 至 59.5703 Hz 的 noise，uniform distribution over interval (-0.2,0.2)。

最後得到的震測訊號如圖 4.9(a)。圖中間上方 0 點位置為爆炸點，橫軸為接收 器的編號，左右各 64 個，縱軸由上往下，為每一個接收器接收到直接波與反射

波訊號的時間。

我們將震測訊號作 Envelope 處理，得到圖 4.9(b)。接著將振幅超過 0.2 的 點取出，共取出 432 點，作為偵測圖形系統的輸入，每個點的值是[

x

,

i

]，x是 接收器的編號， i 是接收到訊號的時間除以 sampling interval，代表該接收站到 接收到的第幾個訊號。我們得到的輸入點如圖 4.9(c)。我們設定偵測圖形數目 k

= 2，偵測出兩條雙曲線，如圖 4.9(d)，將結果顯示在原始訊號上，如圖 4.9(e)，

4.9(f)為

gBest 的適應值 vs. iterations。表 4.7 為參數偵測結果。

圖 4.9. (a) 模擬水平地層單炸點震測訊號圖。

圖 4.9. (b) 模擬單炸點震測訊號圖經過 Envelope 處理。

圖 4.9. (c) 經過 thresholding 後的模擬水平地層震測訊號輸入點。

圖 4.9. (d) 模擬水平地層單炸點震測訊號偵測結果。

圖 4.9. (e) 模擬水平地層單炸點震測訊號偵測結果與原始訊號。

圖 4.9. (f) gBest 的適應值 vs. iterations。

表 4.7. 模擬水平地層單炸點震測訊號圖參數偵測結果。

c

x

c

_{y a b}

r

θ ^{執行時間(sec)}

直接波 -0.24 -33.19 -4.83 0.20 6,636.69 -0.13°

反射波 -0.42 -12.93 5.02 -0.19 -11.04 -0.19°

112.21

由震測訊號的公式，我們可以知道直接波的圖形，實際上為反射波圖形的 漸近線，因此我們可以設定直接波圖形的係數 r = 0，同時它的方向是南北開 口，可以設定θ = 0°。因此對於直接波圖形，只需要 4 個參數

[ ^c

x

^c

y

^{a b} ]

。我 們以這樣設定的參數去偵測圖形，結果如圖 4.9(g)，直接波圖形的偵測結果為 漸近線，4.9(h)為 gBest 的適應值 vs. iterations。表 4.8 為直接波圖形使用 4 個 參數的偵測結果。

圖 4.9. (g) 直接波圖形使用 4 個參數的偵測結果。

圖 4.9. (h) 直接波圖形使用 4 個參數的 gBest 的適應值 vs. iterations。

表 4.8. 直接波圖形使用 4 個參數的偵測結果。 Ricker Wavelet 作 convolution，得到 impulse 訊號， 再加上 random white band 10.2539 Hz 至 59.5703 Hz 的 noise，uniform distribution over interval (-0.2,0.2)。

最後得到的震測訊號如圖 4.11(a)。圖中間上方 0 點位置為爆炸點，橫軸為接收 器的編號，左右各 64 個，縱軸由上往下，為每一個接收器接收到直接波與反射 波訊號的時間。

我們將震測訊號作 Envelope 處理，得到圖 4.11(b)。接著將振幅超過 0.2 的 點取出，共取出 352 點，作為偵測圖形系統的輸入，每個點的值是[

x

,

i

]，x是 接收器的編號， i 是接收到訊號的時間除以 sampling interval，代表該接收站到 接收到的第幾個訊號。我們得到的輸入點如圖 4.11(c)。我們設定偵測圖形數目 為 2，偵測出兩條雙曲線，如圖 4.11(d)，將結果顯示在原始訊號上，如圖 4.11(e)，圖 4.11(f) 為

gBest 的適應值 vs. iterations。表 4.11 為參數偵測結果。

圖 4.11. (a) 模擬傾斜地層震測訊號。

圖 4.11. (b) 模擬傾斜地層震測訊號圖經過 Envelope 處理。

圖 4.11. (c) 經過 thresholding 後的模擬傾斜地層震測訊號輸入點。

圖 4.11. (d) 模擬傾斜地層震測訊號偵測結果。

圖 4.11. (e) 模擬傾斜地層單炸點震測訊號偵測結果與原始訊號。

圖 4.11. (f) gBest 的適應值 vs. iterations。

表 4.9. 模擬傾斜地層單炸點震測訊號圖參數偵測結果。

c

x

c

_{y a b}

r

θ ^{執行時間(sec)}

直接波 -0.39 -13.02 5.08 -0.20 -12.15 0.13°

反射波 -9.34 -26.74 -4.89 -0.19 6,358.34 -0.09°

94.36

同樣的，我們設定直接波的參數

r = 0，θ = 0°。因此對於直接波圖形，只

需要 4 個參數

[ ^c

x

^c

y

^{a b} ]

。我們以這樣設定的參數去偵測圖形，結果如圖 4.11(g)，直接波圖形的偵測結果為漸近線，4.11(h) 為

gBest 的適應值 vs.

iterations。表 4.10 為直接波圖形使用 4 個參數的偵測結果。

圖 4.11. (g) 直接波圖形使用 4 個參數的偵測結果。

4.11. (h) 直接波圖形使用 4 個參數的

gBest 的適應值 vs. iterations。

表 4.10. 直接波圖形使用 4 個參數的偵測結果。

c

x

c

_{y a b r} θ ^執行時間

(sec) 直接波 -0.51 -14.32 -5.34 0.22 0 (preset) 0° (preset)

反射波 -0.31 -12.14 5.12 -0.17 -12.43 -0.51°

97.23

4.3 真實單炸點震測訊號圖形偵測

有兩張真實的水平地層單炸點震測訊號圖形，由 Colorado School of Mine 所建立的 Seismic Unix System 所取得，我們以同樣方法偵測其直接波與反射波 訊號圖形。

4.3.1 偵測 Canadian Artic 的單炸點震測訊號

第一張震測訊號圖形地點來自 Canadian Artic，共有 48 個接收站，每個接 收站接收 3,100 個訊號，接收站之 sampling interval 為 0.002 sec。原始的震測訊 號圖形如圖 4.12(a)。橫軸由左到右為接收站的號碼，縱軸由上到下為接收到訊 號的時間。

我們同樣的將訊號作 envelope 處理，如圖 4.12(b)。留下震幅超過 0.2 的 samples，此外我們去除底下多餘的圖形 (1.7 秒以下)，以免對偵測結果造成太 大干擾。最後共得到 106 點，作為偵測系統的輸入點，如圖 4.12(c)。每個輸入 點的值是[

x

,

i

]， x是接收器的編號， i 是接收到訊號的時間除以 sampling interval，代表該接收站到接收到的第幾個訊號。我們設定偵測圖形數目

k = 3，

偵測出 3 條雙曲線，如圖 4.12(d)，將結果顯示在原始訊號上，如圖 4.12(e)。

表 4.11 為震測訊號參數偵測結果。圖 4.12(f) 為 gBest 的適應值與 iteration 的 關係圖。

圖 4.12. (a) Canadian Artic 的單炸點震測訊號。

圖 4.12. (b) 震測訊號經過 Envelope 處理。

圖 4.12. (c) 經過 thresholding 後的偵測系統輸入點。

圖 4.12. (d) 震測訊號偵測結果。

圖 4.12. (e) 震測訊號偵測結果與原始訊號。

表 4.11. 震測訊號參數偵測結果。

c

x

c

_{y a b}

r

θ ^執行時間

(sec) Pattern 1 23.35 -135.62 -19.11 0.05 714.82 -0.19°

Pattern 2 25.90 17.04 -26.85 0.03 0.01 0.43°

Pattern 3 24.68 -556.45 0.15 -6.33 -2.19 90.80°

219.69

圖 4.12. (f) gBest 的適應值 vs. iterations。

4.3.2 偵測 Gulf of Cadiz 的單炸點震測訊號

第二張震測訊號圖形地點來自 Gulf of Cadiz，共有 48 個接收站，每個接收 站接收 2,050 個訊號，接收站之 sampling interval 為 0.004 sec。原始的震測訊號 圖形如圖 4.13(a)，橫軸由左到右為接收站的號碼，縱軸由上到下為接收到訊號

第二張震測訊號圖形地點來自 Gulf of Cadiz，共有 48 個接收站，每個接收 站接收 2,050 個訊號，接收站之 sampling interval 為 0.004 sec。原始的震測訊號 圖形如圖 4.13(a)，橫軸由左到右為接收站的號碼，縱軸由上到下為接收到訊號

在文檔中 粒子群演算法於圖形偵測與震測圖形識別之應用 (頁 24-0)