基于小波去噪法的南水北調工程沉降預警

孫 飛，翟福雷，方 琪

(1.南水北調東線江蘇水源有限責任公司宿遷分公司，江蘇 宿遷 223800;2.中科寧圖技術江蘇有限公司，江蘇 南京 211100)

南水北調東線工程是解決我國北方地區水資源嚴重短缺問題的重大戰略舉措，從長江下游取水，沿京杭運河向北送水，至山東、河北、天津，以解決京滬沿線和山東半島城市缺水問題，并為農業和生態補水。泗洪站樞紐工程是南水北調東線工程的第四梯級泵站，設計流量為120m3/s，其主要任務是將第三梯級抽入洪澤湖的江水通過運西線徐洪河繼續北送至第五梯級濉寧站，再由房亭河入駱馬湖。樞紐位于江蘇省泗洪縣朱湖鄉東南的徐洪河上，三岔河大橋下游4km、洪澤湖顧勒河口上游約16km處。該樞紐是南水北調的重要節點，自建成運行以來一直在進行安全監測，積累了較多數據，所以，對其安全監測成果進行研究很有必要。文章采用小波去噪法對安全監測數據進行處理，剔除測量誤差，然后根據去噪后的數據進行沉降預警與分析，指導樞紐的運行與維護。

1 安全監測數據

泵站工程自2013年以來，底板位置沉降9mm左右，最大累計沉降位移量為-13mm；上左翼沉降2mm左右，最大累計沉降位移量為-7.1mm；上右翼沉降7mm左右，最大累計沉降位移量為-10.8mm；下右翼沉降6mm左右，最大累計沉降位移量為-8.7mm；下左翼沉降1mm左右，最大累計沉降位移量為-7.7mm，其中2014年9—12月，發生了上升回彈現象，2016年汛后較汛前沉降明顯，12月較9月又發生了數據回彈現象，從檢測數據來看，汛期沉降在2～7mm波動。部分觀測數據見表1—2。

表1 泵站2013—2015年部分垂直位移監測數據 單位：mm

2 小波去噪

2.1 小波去噪理論

在實際工程建設中，傳統的去噪方法是利用一個濾波器除掉帶有噪聲污染的信號中的噪聲相應的頻率波段，具有很強烈的局限性。小波分析去噪作為全新的克服信號缺陷的一種分析方法，已經成功地廣泛應用于多個領域。小波就是一族小波基函，描述了信號空間和尺度的局部特性，小波的局部分析是其最大優勢，可在任意的時間域或空間域中進行信號分析。

表2 泵站2016—2018年部分垂直位移監測數據 單位：mm

小波分析方法在圖像處理和數據去噪中有著廣泛應用，常見的小波去噪方法有小波分解與重構法、小波閾值去噪法、平移不變量法和小波變換模極大值法，其中閾值去噪方法實現簡單、計算速度快，具有廣泛的適應性，是眾多小波去噪方法中應用最為廣泛的一種[1]。閾值去噪方法是通過小波對信號進行分解得到各層系數，然后構造相應的閾值，對大于或小于閾值的系數分別處理，將處理后的小波系數進行重構達到去噪的目的[2]。其去噪原理如圖1所示。其中，S表示原始信號，A表示近似信號，D為噪聲部分，下標為分解的層數。

圖1 小波去噪原理圖

泗洪站垂直位移監測數據都是離散序列，所以采用離散小波變換(DWT)。在數據序列的時頻特性上，噪聲信號的分布是全局的，對小波系數的各個尺度都有影響，目標信號的分布是局部的，主要特征集中在比較大的小波系數上，因此只要用合適的閾值對小波系數進行處理，剔除噪聲產生的小波系數，同時保留有用信號的小波系數，就可以去除數據中的噪聲。小波去噪的流程圖如圖2所示。

圖2 小波去噪過程圖

常用的閾值去噪函數有軟閾值去噪法和硬閾值去噪法。硬閾值去噪法是當小波系數的絕對值小于給定閾值時，令其為零；大于閾值時，令其保持不變，即

(1)

軟閾值去噪法是當小波系數的絕對值小于給定閾值時，令其為零；大于閾值時，令其都減去閾值，即

(2)

文章主要采用db4小波函數將其分成3層進行分解，然后確定閾值標準進行小波重構，從而去除誤差噪聲。

2.2 小波去噪過程

由于測量時間并不固定，測量間隔有長有短，因此為了消除時間序列的不平穩趨勢，需要對測量數據進行兩次差分運算，消除其趨勢性。然后利用標準化處理公式，對數據進行標準化處理，從而得到預處理的數據。通過計算自相關和偏相關系數，來檢驗預處理后的數據是否符合AR建模要求。

根據偏相關函數截尾性，初步判斷模型階次為5。用最小二乘法估計參數，計算10階以內的模型殘差方差和AIC值，應用AIC準則為模型定階。求得的所有AIC值分別為172.6632，165.4660，153.2087，145.1442，140.7898，141.6824，142.9944，144.5601，146.3067，148.7036。

取使AIC值為最小值的階次，從而判斷模型階次為5。檢驗數據是否為白噪聲，求該數據的自相關系數，看其是否趨近于零，若自相關系數接近于零，則數據可看作是白噪聲，應予以去除。

以底板的2-3測點為例，去噪后的數據與原始數據對比如圖3所示。

圖3 去噪前后數據對比

3 沉降預測與分析

3.1 回歸分析沉降預測方法

在進行預測分析方面，回歸分析法適用于單點分析情況。該方法基于對變形觀測數據和影響因子進行回歸分析以及逐步回歸計算分析，由此得到變形與影響因子的函數關系式，回歸分析需要定性地分析變量之間的相關關系，當自變量與因變量存在某種關系時，建立的多元回歸模型才有實際預測的意義。

在建立線性回歸方程后，需要對回歸方程進行統計檢驗，即根據因變量與自變量之間是否具有相關關系進行顯著性統計假設檢驗。預測模型的檢驗一般有F顯著檢驗、T顯著檢驗、擬合優度的檢驗和回歸標準差檢驗等方法。只有當他們間存在線性關系，回歸方程才具有實際意義。

文章主要以泗洪站的監測數據為依據進行沉降分析，考慮到溫度和荷載對閘站沉降的影響，選取4期同月份數據，時間間隔為1年，建立多元線性回歸預測模型，預測泗洪站沉降曲線。

3.2 回歸分析沉降差異分析與預警

3.2.1垂直位移特征值計算

安全預警的第一步是需要對垂直位移特征值進行計算，它是預警模型中的一個參數，所以文章針對泗洪站垂直位移監測數據進行特征值計算。均方差作為數據離散程度的衡量指標，3δ準則又被稱為拉依達準則，它是先假設一組檢測數據只含有隨機誤差，對其進行計算處理得到標準偏差，按一定概率確定一個區間，認為凡超過這個區間的誤差，就不屬于隨機誤差而是粗大誤差[6]。特征值的計算就是根據已有的監測數據進行統計計算，所以將均方差和3δ準則考慮其中，以2δ作為垂直位移的預警特征值，特征值的計算公式如下：

(3)

式中，λ—垂直位移特征值；xi—每個測點的監測值；r—平均值。

其中底板、上右翼、上左翼、下右翼均值及特征值曲線如圖4—7所示。

圖4 底板均值及特征值曲線

3.2.2垂直位移預警

對于水工建筑物之一的泵站而言，某測點在t時刻的垂直位移主要是由沉降相關量、溫度、時間效應這3種因素的共同作用引起的，其通用數學模型為

δ(t)=δH(t)+δT(t)+δθ(t)

(4)

式中，δ(t)—垂直位移效應量在時間t的統計估計值；δH(t)—沉降相關量；δT(t)—溫度位移；δθ(t)—時效位移。

圖5 上右翼均值及特征值曲線

圖6 上左翼均值及特征值曲線

圖7 下右翼均值及特征值曲線

選取泵站底板等水工建筑物結構連同復合地基為研究對象，對其進行單元網格劃分，并將泵站垂直位移測點的位置作為單元節點之一，根據泵站及基礎的結構形式給定相應部位的物理力學參數。

根據若干個特征值λ1，λ2，λ3，…，λn計算得出相應的沉降相關量值，利用有限元法分別計算在不同沉降相關量作用下測點的垂直位移量δ1，δ2，…，δn。由這些位移計算值采用回歸分析建立沉降相關量H(t)與垂直位移δi(t)的關系，得到一元多項式：

(5)

式中，αi—回歸系數；n值一般取3或4。

考慮到時效位移δθ(t)是一種隨時間推移而逐漸增長的變形，一般與時間呈曲線關系，由于成因較復雜，故利用數理統計方法得出其關系式：

(6)

式中，c0，c1，c2，c3，c4—由統計分析確定的系數；t—位移觀測日相對基準觀測日的日序值乘以1/365。

根據泗洪站垂直位移的日常監測資料，建立垂直位移混合模型。由前面的垂直位移模型構造式，將水壓因子的擬合式作為一個單獨的待選因子；溫度因子以觀測日前j天到前k天的日平均氣溫Tj-k作為溫度影響的待選因子，具體采用前1～10d，前11～30d，前31～60d，前61～90d，前91～120d共5個待選因子；時效因子以ln(t+1)、1-e-t和t等3個因子作為時效影響的待選因子[9]，則泗洪站垂直位移的混合模型表達式為

(7)

取顯著性水平為0.01，引入因子、剔除因子的F統計量限值為2.50，由此建立各個測點垂直位移的混合模型，代入上式得出混合模型公式：

(8)

所建的位移混合模型的復相關系數在0.96附近，剩余標準差為0.72mm，剩余標準差S占實測變幅的4.8%，小于10%，且對實際位移情況擬合較好。

根據以上原理，建立垂直位移監控預警值的公式：

[δm]+=ΔδHmax+ΔδTmax+δθ+δ0

(9)

[δm]-=ΔδHmin+ΔδTmin+δθ+δ0

(10)

式中，ΔδHmax，ΔδHmin—最高、最低閘前沉降量與初始閘前沉降量的差值；ΔδTmax，ΔδTmin—極值溫降/升日與初始溫度分量的差值；δθ—初始日至監測序列最后一天的時效分量；δ0—初始日垂直位移實測值[10]。

由此，沉降相關分量、溫度分量和時效分量的極值可根據水平位移混合模型各部分的分離結果分別計算得出，代入式(10)可擬定出泗洪站垂直位移安全監控預警值，截取2015—2017年底板向上、向下垂直位移監控預警值，詳見表3。

表3 底板垂直位移監控預警值表 單位：mm

3.3 回歸分析預測精度與結果

回歸分析法作為預測模型最基本的預測方法，且有一個好的預測效果，所以文章將基于小波去噪的對數函數沉降預測模型應用于泗洪站，并輸出三條曲線：原始數據折線、去噪后數據折線、去噪后擬合曲線。將三條線段進行對比，可以直觀判斷出預測曲線與原有數據的擬合度。預測過程中的殘差平方和與擬合優度的求解過程如下：

(1)計算殘差平方和及相關系數R2(擬合優度)。

(2)計算平滑后的數據的總和zh。

(3)計算平滑后的數據的平均數pjs。

(4)根據模型公式求每一點的預測理論值。yc[i]=ratio3[0] *(Math.Log(x[i]))+ratio3[1]，其中ratio3為擬合數據，yc為預測理論值數組，x為時間。

(5)計算每一點的殘差值。cc[i]=rec[i]-yc[i]，rec為平滑后的數據，cc為殘差數組。

(6)根據公式計算殘差平方和ccpfh。ccpfh=Math.Pow(cc[i]，2)，并保留四位小數。

(7)計算觀測值與平均數之差的絕對誤差，并求總平方和zpfh。

(8)根據公式計算相關系數R2，公式為xgxs=1-ccpfh / zpfh。

經過上述處理最終得到如圖8所示的沉降預測圖。

圖8 泗洪站沉降預測曲線圖

圖8中兩條折線分別表示泗洪泵站原始垂直位移底板1-4監測數據和經過小波去噪后的折線，然后文章根據小波去噪后的折線重采樣數據對其使用對數函數進行擬合，生成圖中平滑曲線，其中自2018年6月7日起進行預測，預測結果生成的報表如圖9所示。

圖9 泗洪站沉降預測報表

4 結 語

文章從基于小波去噪的對數函數入手，將泗洪站沉降預測與模型相結合，通過分析，使用該方法進行預測具有一個較好的效果。同時對垂直位移的特征值進行了計算，分析出差異沉降量，將其作為影響因子運用于預警模型中，計算得出不同時間的預警值。文章也存在一些不足：

(1)實驗數據量不足且數據在時間上應盡量滿足每期的時間間隔大致相同，文章運用的時間間隔后期較大，對模型影響較大。

(2)設計的預測不夠全面，只適合該測區，換其他監測數據所建立的預測模型不是很理想。

(3)文章提出的預警模型雖求出預警值，滿足實際情況，但沒有進行可靠性驗證，今后應從多個方面綜合分析以得到更理想且更加符合實際的閘站預測預警模型。