陳長龍 劉啟方

摘要：

利用2021年日本福島地震（MJMA 7.3）中6個KiK-net臺站的觀測記錄，分別基于移動窗解卷積法、移動窗譜比法和移動窗自相關(guān)函數(shù)法對場地非線性時變過程進(jìn)行研究，分析非線性發(fā)生的閾值和程度。研究表明：三種方法識別場地非線性時變特征的穩(wěn)定性不同。移動窗解卷積法更容易獲得較為穩(wěn)定的土體非線性時變過程，但其反映的是地表至井下平均波速的變化，識別的非線性程度偏低。對于存在強(qiáng)波阻抗比的淺表層土體，移動窗譜比法和移動窗自相關(guān)函數(shù)法識別的非線性程度更強(qiáng);對于較低波阻抗比的淺表層土體，則不容易識別出穩(wěn)定的非線性時變結(jié)果。三種方法識別的非線性最強(qiáng)時刻基本一致，皆位于整條記錄加速度峰值附近。基于三種方法識別的6個臺站的非線性閾值約在40～140 gal，剪切波速下降比在16%～53%之間，場地非線性程度較明顯。

關(guān)鍵詞：

移動窗解卷積; 移動窗譜比; 自相關(guān)函數(shù); 場地非線性; 時變

中圖分類號： P631????? 文獻(xiàn)標(biāo)志碼：A?? 文章編號： 1000-0844（2023）02-0372-10

DOI：10.20000/j.1000-0844.20211110002

Comparative study on three identification methods of the nonlinear time-

varyin? process： a case study of the 2021 Fukushima earthquake in Japan

CHEN Changlong， LIU Qifang

（Suzhou University of Science and Technology， Suzhou 215011， Jiangsu， China）

Abstract：

Based on the records of six KiK-net stations in the 2021 Fukushima earthquake （MJMA7.3）， the nonlinear time-varying process of the site was studied using the moving-time-window deconvolution method， moving-time-window spectral ratio method， and moving-time-window autocorrelation function method. Then， the threshold and degree of nonlinearity were analyzed. Results show that the ability of the three methods to identify the nonlinear time-varying characteristics of the site differs. The moving-time-window deconvolution method can easily obtain a more stable nonlinear time-varying response for soil， but the degree of nonlinearity is low. The method reflects the change of average wave velocity from the surface to the underground. For shallow soils with a strong wave impedance ratio， the nonlinear degree of identification results derived from the moving-time-window spectral ratio method and moving-time-window autocorrelation function method are strong. These two methods cannot provide stable nonlinear time-varying results of shallow soils with a low wave impedance ratio. The strongest nonlinear moments identified by the three methods are basically the same and are located near the peak acceleration of the whole record. The nonlinear threshold values of the six stations identified by the three methods are about 40-140 gal， and the shear wave velocity drops between 16% and 53%， indicating that the degree of site nonlinearity is high.

Keywords：

moving-time-window deconvolution; moving-time-window spectral ratio; autocorrelation function; site nonlinearity; time-varying

0 引言

場地非線性是地震工程領(lǐng)域的研究熱點(diǎn)，土體在強(qiáng)震動作用下會產(chǎn)生損傷，從而導(dǎo)致應(yīng)變幅值增加，剪切模量降低及阻尼比增強(qiáng)［1-7］，非線性發(fā)生的閾值、程度及恢復(fù)過程均是場地非線性研究的主要問題。目前，基于強(qiáng)震動觀測數(shù)據(jù)直接分析場地非線性的方法主要有譜比法［2-3］、廣義反演法［8］、S變換法［9］、地震干涉測量法［10］和自相關(guān)函數(shù)法［11］等。這些傳統(tǒng)分析方法一般是利用主震之前的小震記錄確定場地處于線性狀態(tài)時的水平，利用主震記錄確定非線性狀態(tài)時的水平，利用余震記錄分析非線性恢復(fù)情況?；谝陨戏椒?，研究人員判斷的不同場地非線性閾值差別較大［1-3］，在100～300 gal之間，非線性的恢復(fù)過程也有很大差異，既有震后幾分鐘甚至更短時間內(nèi)恢復(fù)［10］，也有需幾年恢復(fù)［2，12］。另外，研究人員還提出了一些非線性程度指標(biāo)［5，13-14］（DNL、ADNL、PNL等）。

場地經(jīng)歷的強(qiáng)地震動一般會從較低幅值逐漸增加到較高的震動水平，然后再降低到弱震動水平。因此可根據(jù)一條強(qiáng)震記錄作用下不同時段場地的反應(yīng)，獲得場地土的非線性時變特性。近十幾年來，研究人員將移動時間窗技術(shù)與傳統(tǒng)分析方法結(jié)合，提出了移動窗譜比法、移動窗解卷積法和移動窗自相關(guān)函數(shù)法等，針對一條記錄分析不同時間窗內(nèi)非線性程度的變化，從而獲得非線性時變特性。如Wu等［15］ 利用移動窗解卷積法，針對2011 Tohoku 地震（MW9.0）KiK-net臺網(wǎng)中6個臺站主震記錄的分析表明，場地共振頻率明顯下降，等效剪切波速最大下降到主震前的70%。Sawazaki等［16］利用移動窗譜比法和2000年日本W(wǎng)estern-Tottori地震（MW6.7）的主震記錄，對KiK-net臺網(wǎng)的SMNH01臺站場地土的分析表明，在峰值加速度最大值（PGA）附近，剪切模量下降到主震前的52%，等效剪切波速下降30%。Wu等［17］利用移動窗譜比法對KiK-net 臺網(wǎng)中6個臺站1999年1月至2008年5月2 204個地震的觀測記錄分析表明，PGA在60～100 gal時場地發(fā)生非線性反應(yīng)，明顯低于傳統(tǒng)方法的結(jié)果。Federico等 ［18］ 利用移動窗解卷積法和2011 Tohoku 地震（MW9.0）主震記錄，分析了IBRH16臺站土體剪切模量的時變特性，Bonilla等［19］則進(jìn)一步利用該臺站同一記錄，分析了自相關(guān)函數(shù)法和移動窗解卷積法在識別非線性時變特性的差異，結(jié)果表明，兩種方法均可有效識別出場地非線性的發(fā)展過程，但自相關(guān)函數(shù)法獲得的非線性程度大于移動窗解卷積法，且非線性恢復(fù)速度更快。

為了分析移動窗解卷積法、譜比法和自相關(guān)函數(shù)法在識別場地非線性時變特性的差異，本文從KiK-net臺網(wǎng)中選擇6個臺站，其中4個臺站的表層與下臥基巖波阻抗比很大，另外2個臺站土層波速則由表層向下逐漸增大，無強(qiáng)波阻抗比分界面。利用2021年日本福島地震（MJMA7.3）中接收到的主震記錄，分析三種方法對于這兩類不同場地的非線性時變特性識別的差異和可行性。

1 數(shù)據(jù)選取

2021年2月13日23時8分日本福島縣東部近海地區(qū)發(fā)生7.3級地震，震源深度約為55 km，震中位于（37.728°N，141.698°E），日本防災(zāi)技術(shù)研究所（NIED）所屬的KiK-net臺網(wǎng)有超過384個臺站記錄到了完整的主震記錄。此次地震與2011年東日本大地震（MW9.0）的震中距離約75 km，是2011年東日本大地震的余震［20］。本文選擇獲得此次地震主震記錄的6個臺站識別場地的非線性時變特性，圖1為臺站和震中位置。表1給出了這6個臺站的土層地質(zhì)資料［21］，其中1個臺站為C類場地，其余臺站均為D類場地，這些臺站記錄到的地表峰值加速度在300～1 100 gal之間，易發(fā)生非線性。其中FKSH17、FKSH18、FKSH19和FKSH11臺站存在強(qiáng)波阻抗比的分界面，而FKSH20和MYGH10臺站剪切波速變化較平緩。

2 研究方法

本文采用移動窗解卷積法、譜比法和自相關(guān)函數(shù)法分析場地土的非線性時變特征。移動窗解卷積法是對固定寬度時窗內(nèi)的地表和井下記錄解卷積，獲得一個時間窗內(nèi)的等效剪切波速，然后沿時間軸按一定的窗口距離移動時間窗，獲得下一個時間窗內(nèi)結(jié)果［如圖2（a）中紅色框表示一個時間窗，黑色框為下一個時間窗］。通過計算每個窗口的PGA（亦稱為窗內(nèi)PGA），獲得不同時段的土體波速隨時間或窗內(nèi)PGA的變化。根據(jù)等效剪切波速的變化，判斷土體是否發(fā)生非線性及其程度。該方法首先將井下和地表處的地震動記錄進(jìn)行解卷積計算，解卷積結(jié)果G（ω）為：

G（ω）=U1（ω）U2（ω）=U1（ω）U*2（ω）U′2（ω）2 （1）

式中：地表和井下加速度記錄的傅里葉譜分別為U1（ω）和U2（ω）;* 表示復(fù)共軛;分子為地表與井下加速度在頻域內(nèi)的互功率譜，分母為井下加速度的功率譜。為穩(wěn)定計算，一般采用U′2（ω）對U2（ω）的進(jìn)行校正：

U′2（ω）2=U2（ω）2;當(dāng)U2（ω）≥0.1maxU2（ω）

0.1maxU2（ω）2;當(dāng)U2（ω）<0.1maxU2（ω）（2）

然后對解卷積結(jié)果進(jìn)行傅里葉逆變換，獲得井下脈沖輸入時地表的響應(yīng)。通過脈沖響應(yīng)獲得井下臺站到地表地震波的傳播時間（一般取正或負(fù)時間軸上最大值所對應(yīng)的時刻），利用式（3）計算土體等效剪切波速。

vS=hΔt （3）

式中：vS代表等效剪切波速;h代表井下臺站到地表的距離;Δt代表井下臺站到地表脈沖波的傳播時間。

移動窗譜比法中窗口的選擇與移動窗解卷積類似，不同的是通過計算每一時間窗內(nèi)地表和井下記錄傅氏譜比，獲得時間窗內(nèi)的峰值頻率，得到峰值頻率隨時間或窗內(nèi)PGA的變化。在進(jìn)行譜比前，一般需要對傅里葉譜進(jìn)行平滑，本文采用0.5 Hz的帕曾窗進(jìn)行平滑。通過峰值頻率的變化，判斷土體是否發(fā)生非線性及其程度。設(shè)v0和f0是線性狀態(tài)土體波速和固有頻率，某一時間窗內(nèi)識別的波速和峰值頻率為v和f，由于每個時間窗內(nèi)的峰值頻率為其下覆土層的固有頻率，正比與土層的波速，則土體波速的變化Δv正比峰值頻率的變化Δf［式（4）］，因此，可以通過固有頻率的變化判斷剪切波速的變化。

Δvv0=Δff0 （4）

移動窗自相關(guān)函數(shù)法［11］是利用一定時間窗內(nèi)地表記錄進(jìn)行自相關(guān)，從地震事件波形中提取土層某一界面反射波信號，觀察時間延遲，由于反射波延時正比于界面之上土層波速，因此可以通過延時的變化計算波速的變化，如式（5）

Δtt=t0-tt0=v0-vv0=-Δvv0 （5）

式中：t0和v0是土層線性狀態(tài)下延時和波速;t和v是某一時間窗自相關(guān)函數(shù)法計算的延時和波速;Δt是延時變化;Δv波速變化。

采用這三種方法進(jìn)行分析時，時間窗寬和窗距的選取是移動窗技術(shù)的關(guān)鍵，其直接決定了結(jié)果的穩(wěn)定性和分辨率。時間窗口寬度過窄，會導(dǎo)致結(jié)果的不穩(wěn)定。時間窗口寬度過長，時間分辨率太低，無法很好地反映頻譜特性隨時間的變化特征。為此，我們參考了相關(guān)研究［16-19］，并進(jìn)行了大量的窗長和窗距測試，最終確定10 s窗長、2 s窗距可獲得穩(wěn)定結(jié)果，上述三種方法采用相同的窗長和窗距。為消除截斷誤差的影響，本文對截取窗口兩端采用10%窗口長度的邊瓣余弦函數(shù)進(jìn)行處理。由于所有觀測記錄未濾波的加速度傅氏譜對應(yīng)的能量主要集中在1～13 Hz范圍，本文選擇此頻段進(jìn)行濾波。

3 分析與討論

圖2～圖7給出了利用移動窗解卷積法、譜比法和自相關(guān)函數(shù)法識別的4個存在強(qiáng)波阻抗比分界面臺站的非線性時變過程。

FKSH17臺站的鉆孔深度100 m，vS30為544 m/s，屬于美國NEHRP規(guī)范中的C類場地，鉆孔資料表明（表1），表層6 m是剪切波速為180 m/s 的砂礫，其下是剪切波速為1 100 m/s的風(fēng)化花崗巖和花崗巖，兩層之間的波阻抗比為9.9（表1），存在強(qiáng)波阻抗比分界面，強(qiáng)波阻抗比分界面之上的土層總厚度為6 m，等效剪切波速為180 m/s。由于風(fēng)化花崗巖波速較高，非線性很可能主要發(fā)生在淺表層6 m的砂礫中，地表的峰值加速度（PGA）為333 gal。三種方法在該臺站識別非線性時變過程大致相同，即非線性程度隨著窗內(nèi)PGA的增大而增強(qiáng)，當(dāng)PGA達(dá)到整條記錄的峰值時，非線性特征最明顯，隨著PGA的逐漸減小，非線性逐步恢復(fù)。從非線性發(fā)生的閾值看，三種方法識別出結(jié)果差別不大，在60～70 gal之間。非線性下降程度也基本一致，識別的波速最大下降分別為22%、24%和22%（表2），均位于PGA最大值附近。

FKSH18臺站鉆孔深度100 m，vS30為307 m/s，屬于美國NEHRP規(guī)范中的D類場地，表層4 m是剪切波速為140 m/s的風(fēng)化花崗巖，第二層8 m是剪切波速為220 m/s的風(fēng)化花崗巖，第三層18 m是剪切波速為550 m/s的風(fēng)化花崗巖。二三層之間的波阻抗比為3.8。強(qiáng)波阻抗比分界面之上的土層總厚度為12 m，等效剪切波速為185 m/s。地表的PGA為424.8 gal，三種方法識別出的非線性閾值相差不大，在56～66 gal之間。三種方法識別的非線性下降程度不同，識別的波速最大下降為29%、40%和43%（表2），譜比法和自相關(guān)法的結(jié)果較為接近，且明顯大于解卷積法的結(jié)果。

FKSH19臺站鉆孔深度100 m，vS30為338 m/s，屬于美國NEHRP規(guī)范中的D類場地。表層2 m是剪切波速為170 m/s的表土，第二層18 m是剪切波速為270 m/s的花崗巖，第三層20 m是剪切波速為970 m/s的花崗巖，二三層之間的阻抗比為5.0，強(qiáng)波阻抗比分界面之上的土層總厚度為20 m，等效剪切波速為255 m/s。地表的PGA為513.2 gal，解卷積法和譜比法識別的非線性閾值較接近，分別為141 gal和133 gal，從自相關(guān)函數(shù)法得到的延時結(jié)果看，突變點(diǎn)不明顯，不易識別準(zhǔn)確結(jié)果，本文識別為55 gal，三種方法識別的波速最大下降為16%、28%和23%（表2），仍表現(xiàn)為譜比法和自相關(guān)函數(shù)法的結(jié)果較為接近，且大于解卷積法的結(jié)果。

FKSH11臺站鉆孔深度115 m，vS30為240 m/s，屬于美國NEHRP規(guī)范中的D類場地，表層1 m是剪切波速為110 m/s 的礫質(zhì)土，第二層33 m是剪切波速為250 m/s的粉土混砂礫，第三層22 m是剪切波速為1 200 m/s的溶結(jié)凝灰?guī)r，二三層之間的波阻抗比為5.2。強(qiáng)波阻抗比分界面之上的土層總厚度為34 m，等效剪切波速為241 m/s。地表的PGV為353.9 gal，三種方法識別非線性時變特性的能力不同，解卷積法和自相關(guān)函數(shù)法識別的非線性時變特性較明顯，識別的非線性閾值分別為59 gal和63 gal，波速最大下降為18%、32%（表2），隨著PGA的逐漸減小，波速逐步恢復(fù)。但自相關(guān)函數(shù)法識別的恢復(fù)過程振蕩較大。而譜比法則幾乎識別不出非線性的時變過程，該臺站譜比結(jié)果在1 Hz、5 Hz、10 Hz處顯示3條峰值曲線，但各峰值曲線隨時間變化均不明顯。

圖8給出了利用移動窗解卷積法、譜比法和自相關(guān)函數(shù)法識別的2個無強(qiáng)波阻抗比分界面臺站的非線性時變過程。MYGH10臺站鉆孔深度205 m，vS30為348 m/s，屬于美國NEHRP規(guī)范中的D類場地。表層1m是剪切波速為110 m/s 的黏性土，第二層2 m是剪切波速為250 m/s的砂巖，第三層31 m是剪切波速為390 m/s的礫巖和砂巖，第四層80 m是剪切波速為590 m/s的凝灰質(zhì)砂巖，第五層91 m是剪切波速為770 m/s的砂巖，土層剪切波速逐漸增大，無強(qiáng)分界面。該臺站地表PGA很高，達(dá)1041.5 gal，值很大。從圖8可以看出，三種方法識別出的非線性閾值相差不大，分別識別的閾值為103、93和98 gal。三種方法識別的非線性下降程度不同，識別的波速最大下降為18%、53%和42%，均位于PGA最大值附近。譜比法和自相關(guān)函數(shù)法的結(jié)果較為接近，且顯著大于解卷積法的結(jié)果，約為2.5～3倍。

FKSH20臺站鉆孔深度109 m，vS30為350 m/s，屬于美國NEHRP規(guī)范中的D類場地。表層32 m是剪切波速為350 m/s的腐殖土和砂礫，第二層28 m是剪切波速為500 m/s的砂礫，第三層49 m是剪切波速為610 m/s的砂巖，一二層和二三層的波阻抗比為1.4和1.3，阻抗比較小，剪切波速逐漸增大，無強(qiáng)波阻抗比分界面。地表的PGA為345 gal。從圖8可以看出，解卷積法、譜比法分別識別的閾值為39 gal和35 gal，識別的波速最大下降為32%和35%，均很接近。自相關(guān)函數(shù)法則無法識別非線性時變過程。

從上述分析可見，對于存在強(qiáng)阻抗比分層的4個臺站，三種方法識別的閾值總體上差別不大，但在識別非線性程度上，譜比法和自相關(guān)函數(shù)法的結(jié)果較為接近，除FKSH17臺站，后兩種方法識別的非線性程度均且明顯大于解卷積法的結(jié)果，這很可能因為譜比法和自相關(guān)函數(shù)法更容易反映強(qiáng)阻抗比界面之上的淺層波速變化，而解卷積法則反映地表至井下平均波速的變化，由于強(qiáng)阻抗比界面之下的波速較高，非線性程度可能很低（或者不發(fā)生），所以平均波速的變化會低于強(qiáng)阻抗比界面之上土層波速的變化。Bonilla等［19］利用獲得2011年東日本大地震主震記錄的IBRH16臺站，分析了自相關(guān)函數(shù)法和移動窗解卷積法在識別非線性時變特性的差異。結(jié)果表明，自相關(guān)函數(shù)法獲得的非線性程度大于移動窗解卷積法，本文結(jié)果與其基本一致。

對于無強(qiáng)阻抗比界面的土層，波速隨深度變化不明顯，因而通過譜比法和自相關(guān)函數(shù)法不易獲得穩(wěn)定的峰值頻率和土層界面的反射波延時，不易獲得非線性時變過程。此時解卷積法仍可得到穩(wěn)定結(jié)果，主要由于該方法給出的是井下臺站到地表的平均波速，強(qiáng)震動作用下土層平均剪切波速隨震動水平的增加仍出現(xiàn)了明顯的下降，因此仍可以較好地反映出土層非線性的平均變化，如FKSH20臺站。同時我們發(fā)現(xiàn)，盡管土層中不存在明顯強(qiáng)阻抗比界面，但當(dāng)場地PGA很高，如MYGH10臺站達(dá)到1 041.5 gal，三種方法也均能識別出較穩(wěn)定非線性時變過程。

4 結(jié)論

本文利用2021年日本福島地震（MJMA 7.3）中6個KiK-net 臺站的觀測記錄，分別用移動窗解卷積法、譜比法和自相關(guān)函數(shù)法分析了場地非線性時變過程，研究表明：

（1） 三種方法識別場地非線性時變特征的穩(wěn)定性不同。解卷積法更容易獲得較為穩(wěn)定的土體非線性時變過程，但其反映的是地表至井下平均波速的變化，識別的非線性程度偏低。但對于淺表層存在強(qiáng)阻抗比的土體，譜比法和自相關(guān)函數(shù)法可以獲得更強(qiáng)的非線性程度，對于淺表層不存在強(qiáng)阻抗比的土體則不容易識別出穩(wěn)定的非線性時變過程。三種方法識別的非線性閾值基本上相差不大，對于非線性最強(qiáng)時刻的判定也基本相同，皆位于整條記錄加速度峰值附近。

（2） 6個臺站地表峰值加速度在300～1 100 gal之間，三種方法識別的非線性閾值在40～140 gal之間，剪切波速下降比在16%～53%之間，表明場地產(chǎn)生較嚴(yán)重的非線性反應(yīng)。對于PGA最大（1 041.5 gal）的MYGH10臺站，剪切波速最大下降達(dá)到53%。

本文僅用6個臺站初步分析了三種方法的優(yōu)缺點(diǎn)，由于場地非線性的復(fù)雜性，擬在今后的研究中，我們將通過更多的鉆井臺陣記錄，進(jìn)一步分析上述方法的可靠性和有效性。

參考文獻(xiàn)（References）

［1］ BERESNEV I A，WEN U K L.Nonlinear soil response：a reality？［J］.Bulletin of the Seismological Society of America，1996，86（6）：1964-1978.

［2］ WEN K L，BERESNEV I A，YEH Y T.Nonlinear soil amplification inferred from downhole strong seismic motion data［J］.Geophysical Research Letters，1994，21（24）：2625-2628.

［3］ HARTZELL S.Variability in nonlinear sediment response during the 1994 Northridge，California，earthquake［J］.Bulletin of the Seismological Society of America，1998，88（6）：1426-1437.

［4］ 王海云.基于強(qiáng)震觀測數(shù)據(jù)的土層場地反應(yīng)的研究現(xiàn)狀［J］.地震工程與工程振動，2014，34（4）：42-47.

WANG Haiyun.A review of study on soil site response estimating from strong motion data［J］.Earthquake Engineering and Engineering Dynamics，2014，34（4）：42-47.

［5］ 任葉飛.基于強(qiáng)震動記錄的汶川地震場地效應(yīng)研究［D］.哈爾濱：中國地震局工程力學(xué)研究所，2014.

REN Yefei.Study on site effect in the Wenchuan earthquake using strong-motion recordings［D］.Harbin：Institute of Engineering Mechanics，CEA，2014.

［6］ 李小軍.非線性場地地震反應(yīng)分析方法的研究［D］.哈爾濱：中國地震局工程力學(xué)研究所，1993.

LI Xiaojun.Study on seismic response analysis method of nonlinear site［D］.Harbin：Institute of Engineering Mechanics，CEA，1993.

［7］ 陳學(xué)良.土體動力特性、復(fù)雜場地非線性地震反應(yīng)及其方法研究［D］.哈爾濱：中國地震局工程力學(xué)研究所，2006.

CHEN Xueliang.Study on soil dynamic characteristics，nonlinear seismic response of complex site and its methods［D］.Harbin：Institute of Engineering Mechanics，CEA，2006.

［8］ REN Y F，WEN R Z，YAMANAKA H，et al.Site effects by generalized inversion technique using strong motion recordings of the 2008 Wenchuan earthquake［J］.Earthquake Engineering and Engineering Vibration，2013，12（2）：165-184.

［9］ BONILLA L F，TSUDA K，PULIDO N，et al.Nonlinear site response evidence of K-NET and KiK-net records from the 2011 off the Pacific coast of Tohoku earthquake［J］.Earth，Planets and Space，2011，63（7）：785-789.

［10］ YAMADA M，MORI J，OHMI S.Temporal changes of subsurface velocities during strong shaking as seen from seismic interferometry［J］.Journal of Geophysical Research：Solid Earth，2010，115（B3）：B03302.

［11］ SCHERBAUM F.Seismic imaging of the site response using microearthquake recordings.Part Ⅱ.application to the Swabian Jura，southwest Germany，seismic network［J］.Bulletin of the Seismological Society of America，1987，77（6）：1924-1944.

［12］ 王蘇陽.基于日本KiK-net地震動數(shù)據(jù)的場地反應(yīng)研究［D］.哈爾濱：中國地震局工程力學(xué)研究所，2017.

WANG Suyang.Study on site effect using ground motion data from KiK-net in Japan［D］.Harbin：Institute of Engineering Mechanics，CEA，2017.

［13］ NOGUCHI S，SASATANI T.Quantification of degree of nonlinear site response［C］// Proceedings of the 14th World Conference on Earthquake Engineering.Beijing：Seismological Press，2008.

［14］ REGNIER J，CADET H，BONILLA L F，et al.Assessing nonlinear behavior of soils in seismic site response：statistical analysis on KiK-net strong-motion data［J］.Bulletin of the Seismological Society of America，2013，103（3）：1750-1770.

［15］ WU C Q，PENG Z G.Temporal changes of site response during the 2011 MW9.0 off the Pacific coast of Tohoku earthquake［J］.Earth，Planets and Space，2011，63（7）：791-795.

［16］ SAWAZAKI K，SATO H，NAKAHARA H，et al.Time-lapse changes of seismic velocity in the shallow ground caused by strong ground motion shock of the 2000 western-Tottori earthquake，Japan，as revealed from coda deconvolution analysis［J］.Bulletin of the Seismological Society of America，2009，99（1）：352-366.

［17］ WU C Q，PENG Z G，BEN-ZION Y.Refined thresholds for non-linear ground motion and temporal changes of site response associated with medium-size earthquakes［J］.Geophysical Journal International，2010，182（3）：1567-1576.

［18］ FEDERICO D，BONILLA L F，F(xiàn)OTI S.In-situ shear modulus reduction computation using seismic interferometry by deconvolution from borehole and surface data：theory and examples［C］//Proceedings of the 5th IASPEI/IAEE International Symposium：Effects of Surface Geology in Seismic Motion.Taiwan：［s.n.］，2016：15-17.

［19］ BONILLA L F，GUGUEN P，BEN-ZION Y.Monitoring coseismic temporal changes of shallow material during strong ground motion with interferometry and autocorrelation［J］.Bulletin of the Seismological Society of America，2019，109（1）：187-198.

［20］ 張喆，許力生.2021年2月13日日本本州東海岸MW7.2地震矩心矩張量解［J］.地震學(xué)報，2021，43（2）：255-259.

ZHANG Zhe，XU Lisheng.The centroid moment tensor solution of the 13 February 2021 MW7.2 earthquake in the east coast of Honshu，Japan［J］.Acta Seismologica Sinica，2021，43（2）：255-259.

［21］ GARDNER G H F，GARDNER L W，GREGORY A R.Formation velocity and density：the diagnostic basics for stratigraphic traps［J］.Geophysics，1974，39（6）：770-780.

［22］ NIED.NIED K-NET，KiK-net，National Research Institute for Earth Science and Disaster Resilience［Z］.2019.doi：10.17598/NIED.0004.