把測站波形丟進 ObsPy 前，我先做了這三件小事

上個月幫朋友處理一批 2023 年從社區強震儀抓下來的 miniSEED，總共 6 個測站、每個站 3 通道、取樣率 100 Hz，時間跨度約 11 天。檔案加起來大概 4.2 GB。目标很簡單：把每天的事件窗切出來，看 P 波到時，順便算個 STA/LTA 當作輔助。一開始直接把整包丟進 ObsPy 的 read，結果跑了一陣子才回來，而且記憶體直接吃到 8 GB 以上。後來改成先做三件事，才比較順。

第一件是先把檔案按測站、按天切開，不要讓 read 一次吞整批。我用的是 obspy.io.mseed 的 _read_mseed 逐檔處理，再自己組 Stream：

from obspy import Stream, UTCDateTime
from obspy.io.mseed import _read_mseed
import glob

for f in sorted(glob.glob('data/*.mseed')):
    st = Stream([tr for tr in _read_mseed(f)])
    st.merge(method=1, fill_value=0)
    st.trim(UTCDateTime('2023-09-01'), UTCDateTime('2023-09-02'))
    st.write(f'day/{f.split("/")[-1]}', format='MSEED')

第二件是先檢查取樣率與時間碼有沒有跳。這批資料裡有兩站因為 GPS 斷過，時間軸出現約 0.4 秒的階梯，merge 不會幫你修，只會讓兩段接縫變得看起來合理。我是用 st[0].stats.sampling_rate 和 st[0].times() 的差分去掃，出現不連續就標記，不硬修。

第三件是先用 st.filter('highpass', freq=1.0) 再算 STA/LTA，而不是在原始波形上算。這批資料的長週期背景雜訊很重，直接算會把整個窗口的比值壓平，看起來什麼都沒觸發。代价是低頻事件（例如遠震的長週期波）會被濾掉，所以這條路只適合近場小事件。

為什麼要這樣做：miniSEED 是逐筆 record 存的，read 一次吃整包等於把所有 record 都攤在記憶體裡，加上 merge 會做時間軸對齊，代價是記憶體和時間都翻倍。先切檔等於把「一次處理的單位」縮到一天，壞掉的檔也只影響那一天。

哪裡會失效：如果事件窗跨日，切檔就會把波形切斷，要嘛改成按事件切，要嘛在切之前先確認事件時間。另外 _read_mseed 是私有函式，ObsPy 改版有可能會動，我自己是把它包一層 try，失敗就退回 obspy.read。還有一個坑：有些強震儀輸出的 MSEED 沒有正確的 npts，merge 之後長度會對不上，這時候要手動 trim 到共同長度，不然後面算 STA/LTA 會 index 爆掉。

你們在處理大批 miniSEED 的時候，是先用 obspy.read 一次吃，還是也會先切檔？如果是跨日事件窗，通常怎麼處理？