AFCM聚類和張量不變量用于磁源多目標(biāo)定位

李青竹，李志寧，石志勇，范紅波

（陸軍工程大學(xué) 車輛與電氣工程系，河北 石家莊 050003）

1 引言

多個(gè)磁目標(biāo)同步探測是極具挑戰(zhàn)的工作。作為磁法勘探領(lǐng)域的研究熱點(diǎn)，目前多目標(biāo)磁定位研究已在山脈礦產(chǎn)勘探［1］、腦磁圖定位［2］、航天器 磁 特 征 檢 測［3］和 未 爆 彈（Unexploded Ordnance，UXO）探測［4］等領(lǐng)域中具有重要應(yīng)用。一般地，目標(biāo)定位和識別（Positioning and Recognition，PR）即是利用測量磁場數(shù)據(jù)對區(qū)域內(nèi)磁性目標(biāo)的位置和磁矩進(jìn)行估算的問題。當(dāng)磁化均勻且源中心與觀測點(diǎn)間距大于目標(biāo)最大尺寸的2.5倍時(shí)，磁源通常被認(rèn)為是磁偶極子［5］。磁偶極子物性參數(shù)可用六個(gè)自由度進(jìn)行描述，即空間位置（x，y，z）和磁矩（mx，my，mz），故多個(gè)緊湊源PR檢測即對多個(gè)磁偶極子的位置和磁矩共六個(gè)參數(shù)進(jìn)行估計(jì)。磁梯度張量（Magnetic Gradient Tensor，MGT）［6］是磁場在空間正交方向上的變化率，也是磁標(biāo)勢的二次偏微分，其信息度更高，受磁化方向影響小，已成為磁探測領(lǐng)域研究前沿。此外，MGT特征值和不變量等信息能減少系統(tǒng)旋轉(zhuǎn)誤差和背景磁場干擾，對于磁源PR更具優(yōu)勢。目前基于MGT測量的磁偶極子單點(diǎn)定位算法已較成熟，如張量衍生不變關(guān)系定位法［7］和二階磁張量歐拉反褶積定位法［8］，通過磁梯度張量系統(tǒng)（Magnetic Gradient Tensor System，MGTS）［5］進(jìn)行航線測量，可對測區(qū)內(nèi)孤立磁源具有較高的單點(diǎn)定位精度。然而，當(dāng)多個(gè)目標(biāo)出現(xiàn)在不同埋深，且其分布和磁矩均不相同時(shí)，弱磁性目標(biāo)常因多個(gè)磁異常場相互疊加被淹沒，這使得僅利用單點(diǎn)方法進(jìn)行磁源PR探測陷入瓶頸。

針對多磁源PR問題，1993年Hansen等［9］提出利用Werner反褶積方法對多源磁目標(biāo)進(jìn)行二維反演；此 后Hansen［10］又 將 多 源Werner反 褶 積方法拓展到了三維并進(jìn)行了盆地復(fù)雜數(shù)據(jù)測量。這些算法的一個(gè)主要缺點(diǎn)是需要在實(shí)際計(jì)算之前定義磁性目標(biāo)的數(shù)量，可實(shí)際應(yīng)用中目標(biāo)數(shù)目是很難提前知曉的。王林飛等人［11］通過歸一化磁場垂向?qū)?shù)與z方向解析信號的比值得到MGT解析信號來識別目標(biāo)最佳范圍和深度，李金朋等人［12］也提出利用歸一化磁源強(qiáng)度和同質(zhì)性方程對測區(qū)平面進(jìn)行邊界識別以估計(jì)平面下方目標(biāo)的位置和埋深，但這些方法需要人為對等高線圖進(jìn)行鑒別，以確定目標(biāo)的水平分布和數(shù)目，探測效率低、虛警率高且結(jié)果受噪聲干擾較大。

由此可見，對測區(qū)內(nèi)目標(biāo)數(shù)目進(jìn)行自適應(yīng)識別以適用于更成熟的定位方法并提高位置估計(jì)精度是必要的。Ugalde等［13］曾使用核密度分布算法過濾不相關(guān)目標(biāo)位置解，通過模糊c-means聚類（Fuzzy c-means Clustering，F(xiàn)CM）［14］算法的聚類質(zhì)心來估計(jì)異常源的地質(zhì)分布和走向，這種解決方法可能為磁偶極子的多目標(biāo)檢測和識別提供思路。此外，Davis［15］開發(fā)了一種基于希爾伯特變換的擴(kuò)展歐拉反褶積方法用于多個(gè)UXO目標(biāo)的自動(dòng)檢測，Davis使用不同滑動(dòng)窗口內(nèi)的反演結(jié)果來估計(jì)目標(biāo)位置。參考劃分方窗區(qū)域的思路，我們嘗試對目標(biāo)分布區(qū)域進(jìn)行預(yù)先識別，以提高多目標(biāo)PR檢測精度。然而，經(jīng)典FCM算法需要預(yù)設(shè)初始聚類質(zhì)心數(shù)目，但僅靠經(jīng)驗(yàn)法預(yù)設(shè)數(shù)目無法保證FCM收斂到最優(yōu)解，也不能為多目標(biāo)的自適應(yīng)PR檢測提供幫助。目前已有一些自適應(yīng)模糊c-means（Adaptive Fuzzy c-means Clustering，AFCM）聚類算法［16-17］，通過評價(jià)聚類效果的有效性函數(shù)（Clustering Validity Function，CVF）和試錯(cuò)迭代來確定最優(yōu)聚類數(shù)。本文中，我們設(shè)計(jì)了適用于目標(biāo)初始位置點(diǎn)云的AFCM聚類有效性函數(shù)，為實(shí)現(xiàn)多目標(biāo)自適應(yīng)PR檢測提供可能。

本文方法為基于AFCM聚類和張量不變量的目標(biāo)多源定位方法。首先，利用MGT不變量中對磁化方向相對不敏感的歸一化磁源強(qiáng)度和張量縮并［18］來對MGTS網(wǎng)格內(nèi)測量節(jié)點(diǎn)進(jìn)行目標(biāo)預(yù)識別和區(qū)域劃分。接著利用張量衍生不變關(guān)系定位方法對識別區(qū)域內(nèi)的網(wǎng)格節(jié)點(diǎn)進(jìn)行單點(diǎn)定位計(jì)算，獲得一定數(shù)量的目標(biāo)初始位置解集。最后，利用AFCM算法對這些解集形成的稠密點(diǎn)云進(jìn)行三維聚類并自動(dòng)估計(jì)簇質(zhì)心數(shù)目和質(zhì)心位置。目標(biāo)正上方處的MGT矩陣和位置坐標(biāo)可用以計(jì)算磁矩矢量。該方法具備如下優(yōu)勢：（1）分布區(qū)域的預(yù)識別能極大提高AFCM的聚類自適應(yīng)識別精度；（2）由于MGT隨目標(biāo)距離四次方衰減［18］，處于源正上方區(qū)域的節(jié)點(diǎn)處定位結(jié)果更加可靠，且該過程能排除大部分虛假解；（3）自適應(yīng)算法的目標(biāo)識別過程不需要人為對測區(qū)等高線特征進(jìn)行鑒別和判斷。

2 張量不變量用于目標(biāo)區(qū)域識別

2.1 MGT的不變量

無源靜磁場中，MGT矩陣可對角化且僅五個(gè)張量分量相對獨(dú)立［5］。以G表 示MGT矩 陣，λ1、λ2和λ3表示G的特征值，v1、v2、v3分別為特征值λ1、λ2、λ3對 應(yīng) 的 特 征 向 量，且 滿 足λ3-I0λ2+I1λ-I2=0 及 |λ1|≥|λ3|、|λ2|≥|λ3|、λ2≥λ3≥λ1，λ3稱中間特征值，λ1和λ2分別稱最小和最大特征值。若gij（i，j=x，y，z）為MGT分量，則G及其部分MGT不變量可由以下方程表示：

其中：不變量I0和I2分別為G的跡和行列式；u即為歸一化磁源強(qiáng)度（Normalized Source Strength，NSS）［18］；c為G的Frobenius范數(shù)，又稱張量縮并（Contraction of Tensor，CT）。

2.2 不變量改進(jìn)傾斜角的目標(biāo)區(qū)域識別

平面網(wǎng)格測量的優(yōu)勢是測區(qū)內(nèi)觀測點(diǎn)多，相比單點(diǎn)或航線測量信息量更大。網(wǎng)格測量可利用磁場特征量計(jì)算，獲取測區(qū)下方目標(biāo)的水平分布信息，并極大提高多磁性目標(biāo)定位的效率和精度。由于待測物體間磁異常場相互疊加，弱磁物體常被強(qiáng)磁物體磁場所掩蓋，一些常規(guī)特征量，如磁異?？倛鰪?qiáng)度（Total Magnetic Intensity，TMI）［6］、MGT分量、NSS、CT等，很難通過設(shè)置區(qū)域閾值等高線來圈定目標(biāo)區(qū)域。當(dāng)閾值設(shè)置過大對弱磁性和過埋深目標(biāo)發(fā)現(xiàn)率低，而閾值設(shè)置過小會(huì)使得目標(biāo)虛警率高，這為磁異常的自動(dòng)檢測帶來困難。

傾斜角（θTilt）［19］是位場標(biāo)勢的總水平導(dǎo)數(shù)和垂向?qū)?shù)之比的反正切函數(shù)的角度值，表達(dá)式為：

其中：U為位場標(biāo)勢，?U?x、?U?y和?U?z分別表示U在三軸正交方向上的偏導(dǎo)數(shù)。在磁場中，U為磁標(biāo)勢，其偏導(dǎo)數(shù)即為磁場矢量。磁場傾斜角對源埋深和磁化強(qiáng)度不敏感，且為固定閾值的無量綱特征量，可同時(shí)描述區(qū)域內(nèi)不同深度和強(qiáng)度的磁源，特別適合用于多目標(biāo)數(shù)目和分布區(qū)域的自動(dòng)檢測。

傾斜角的缺點(diǎn)是容易受到磁化方向的影響，因此傾斜角能一定程度識別測區(qū)中埋深和磁矩不同的多個(gè)目標(biāo)的數(shù)目和大致分布，但是分布區(qū)域受磁化方向影響將出現(xiàn)偏差。一般而言，MGT不變量中NSS和CT對磁化方向最不敏感［18］。在此基礎(chǔ)上，本文將NSS、CT引入傾斜角公式，得到兩個(gè)由不變量改進(jìn)的傾斜角，表達(dá)式為：

其中：u為NSS，c為CT，其水平導(dǎo)數(shù)可利用空間域差分法估計(jì)得到；?u?z和?c?z不服從Laplace方程，故其不是NSS和CT的垂向?qū)?shù)［18］，但可利用頻率域Fourier變換計(jì)算得到；θNSSTilt即為NSS不變量改進(jìn)傾斜角，θCTTilt即為CT不變量改進(jìn)傾斜角。θNSSTilt和θCTTilt無量綱，其取值范圍為-π/4到π/4，致密源正上方處取值為π/4，且當(dāng)觀測點(diǎn)遠(yuǎn)離源上方區(qū)域后將迅速衰減至-π/4。θNSSTilt和θCTTilt對磁化方向、源距離均不敏感，當(dāng)選定一固定等高線閾值時(shí)，該值能在區(qū)域MGTS檢測極限范圍內(nèi)［5］快速圈定目標(biāo)數(shù)目和分布位置，實(shí)現(xiàn)多目標(biāo)預(yù)識別。

3 AFCM聚類算法用于多目標(biāo)位置估計(jì)

航線測量中，目標(biāo)單點(diǎn)定位算法可對孤立磁源進(jìn)行連續(xù)追蹤定位，連續(xù)計(jì)算的位置解將在目標(biāo)空間區(qū)域附近形成稠密點(diǎn)云。將航線測量離散化，并將數(shù)條航線并列排開，即形成網(wǎng)格測量。相較于航線測量，網(wǎng)格測量觀測點(diǎn)雖不連續(xù)，但在目標(biāo)區(qū)域預(yù)識別基礎(chǔ)上，仍能利用單點(diǎn)定位算法計(jì)算目標(biāo)位置正上方識別區(qū)域內(nèi)的節(jié)點(diǎn)的位置點(diǎn)云，這為多目標(biāo)解集的AFCM聚類估計(jì)提供依據(jù)，如圖1所示。

圖1 識別區(qū)域內(nèi)磁偶極子定位結(jié)果形成的稠密點(diǎn)云Fig.1 Dense point cloud formed by magnetic dipole positioning results in the recognition area

3.1 張量衍生不變關(guān)系的目標(biāo)單點(diǎn)定位

尹等［20］提出了MGT幾何不變量的概念，李青竹等［7］針對MGTS測量過程將這種量歸納為MGT衍生不變關(guān)系。若磁偶極子用六個(gè)未知量描述：位置矢量r=（x，y，z）T和磁矩矢量m=（mx，my，mz）T，其中r為偶極子到MGTS觀測點(diǎn)間的相對坐標(biāo)，則有如下關(guān)系。

3.1.1 磁矩夾角不變關(guān)系

其中，‖·‖是2范數(shù)，即m與r間夾角θ是不變的，且等于中間特征值λ3與歸一化磁源強(qiáng)度u之比的反余弦函數(shù)。

3.1.2 特征向量不變關(guān)系

即中間特征值對應(yīng)特征向量v3與m、r所在平面是垂直的，而最大、最小特征值對應(yīng)特征向量v1、v2與m、r所在平面是共面的。

利用衍生不變關(guān)系可對區(qū)域內(nèi)孤立磁偶極子進(jìn)行定位，該方法對每個(gè)觀測點(diǎn)提供了四種可能的坐標(biāo)解，這些位置解關(guān)于觀測點(diǎn)空間象限對稱，且只有一個(gè)是真實(shí)解［7］：

3.2 磁源定位初始位置點(diǎn)云的AFCM聚類分析

當(dāng)檢測區(qū)域存在多個(gè)磁源時(shí)，衍生不變關(guān)系定位方法估計(jì)位置將近似聚集在磁源附近，圖1中每個(gè)點(diǎn)云簇的質(zhì)心可用于估計(jì)源的位置。聚類分析在解決許多問題中一直發(fā)揮著重要作用，而聚類數(shù)［16］是最重要的參數(shù)。計(jì)算的磁偶極子初始位置向量在目標(biāo)真實(shí)位置附近聚集程度很高，而零散的虛假解和失敗解可能極度偏離目標(biāo)區(qū)域，這使得根據(jù)經(jīng)驗(yàn)法則隨機(jī)選擇聚類質(zhì)心時(shí)會(huì)出現(xiàn)虛假質(zhì)心和位置偏移的可能。如果預(yù)定義的聚類數(shù)不正確，聚類結(jié)果可能會(huì)偏離給定數(shù)據(jù)集的真實(shí)結(jié)構(gòu)。因此，實(shí)際應(yīng)用中多目標(biāo)自適應(yīng)PR檢測的關(guān)鍵即是利用算法自動(dòng)找到的最優(yōu)聚類數(shù)。

首先給出FCM聚類的計(jì)算過程。假定R=｛r1，r2，…，rn｝為總數(shù)為n的有限目標(biāo)初始位置向量數(shù)據(jù)集，其中ri為第i個(gè)三維坐標(biāo)向量；V=｛v1，v2，…，vc｝為總數(shù)為c的聚類簇質(zhì)心，其中vi為第i個(gè)三維聚類中心；U=（uik）（n×c）為模糊劃分矩陣，uik表示第i個(gè)樣本屬于第k簇的隸屬（緊湊）程度。uik和vk的初始化及更新過程［17］遵循：

其中：ri-vk表示第i個(gè)向量到第k個(gè)簇質(zhì)心之間的相對距離向量，m（m∈［1，∞））是模糊加權(quán)指數(shù)，控制成員的模糊性。m值越大由此產(chǎn)生的成員隸屬便更加模糊，Pal等人［21］建議m應(yīng)該在1.5和2.5之間，若無特殊要求通常m=2。FCM算法通過更新U和V，實(shí)現(xiàn)如式（9）所示目標(biāo)函數(shù)的最小化迭代過程：

FCM聚類的具體步驟如下：

（1）初始化聚類的簇質(zhì)心數(shù)c、模糊指數(shù)m、最大迭代次數(shù)itermax和閾值ξ；

（2）根據(jù)隸屬度的約束隨機(jī)初始化模糊分區(qū)U（0）；

（3）在第t步迭代，根據(jù)式（8）計(jì)算全部c簇質(zhì)心V（t）；

（5）根據(jù)式（8）計(jì)算U（t+1），并返回步驟（3）；

（6）將當(dāng)前對象按照最大隸屬度原則排列進(jìn)簇，并完成迭代。

FCM聚類需要預(yù)設(shè)聚類數(shù)c，AFCM算法可通過設(shè)置聚類數(shù)的搜索范圍運(yùn)行FCM生成不同聚類數(shù)的聚類結(jié)果，選擇合適的聚類有效性函數(shù)評價(jià)聚類結(jié)果，最終確定FCM最優(yōu)聚類數(shù)。1991年，Xie和Beni設(shè)計(jì)了一種基于模糊緊湊度和分 離度的CVF［22］。隨后，Bensaid對其進(jìn)行 改進(jìn)，設(shè)計(jì)了對集群范圍不敏感的CVF［23］，這種方法十分適用于磁源定位中初始位置點(diǎn)云集群分布不規(guī)律且集群容積大小不一的情形。然而，Bensaid方法在聚類數(shù)與真實(shí)對象數(shù)量接近時(shí)，聚類有效性函數(shù)值單調(diào)遞減且接近于零，使簇質(zhì)心數(shù)尋優(yōu)過程會(huì)過早失去魯棒性和決策功能。因此，我們在Bensaid方法基礎(chǔ)上，通過引入懲罰函數(shù)消除簇?cái)?shù)尋優(yōu)過程的單調(diào)遞減趨勢，設(shè)計(jì)了一種新的基于模糊緊湊度和分離度的模糊聚類有效性函數(shù)實(shí)現(xiàn)AFCM聚類：

式（10）中的分子代表了第k簇的緊湊性，而分母代表了第k簇質(zhì)心到其他簇質(zhì)心的平均距離，用于分離性的度量。則此時(shí)分子分母之比代表了第k簇的聚類效果，L（U，V）即為全部聚類簇的聚類效果（比率）之和。顯然L值隨c的增長先變小后增大，而L取最小時(shí)對應(yīng)的c值即為最優(yōu)聚類數(shù)，這種決策模式能極大程度地排除磁偶極 子 定 位 的 虛 假 解 和 偏 差 解。Bensaid［23］和Sun［16］描述的迭代試錯(cuò)方法可以用于估計(jì)AFCM的最優(yōu)聚類數(shù)量。

3.3 多目標(biāo)磁矩矢量估計(jì)

AFCM聚類算法可以同步估計(jì)區(qū)域內(nèi)多個(gè)磁偶極子的三維位置坐標(biāo)。由于MGT分量隨源的相對距離四次方衰減［18］，對于一切致密源上方的張量場受區(qū)域其他源和地磁場的影響可以忽略，即各磁偶極子正上方處測得MGT可認(rèn)為僅由該源產(chǎn)生。設(shè)估計(jì)的源坐標(biāo)rc=（xc，yc，zc），則源正上方處MGTS與源間相對矢量r0=（0，0，-zc），MGTS測得該點(diǎn)處張量矩陣為Gc，則根據(jù)磁偶極子正演方程和歐蘭反褶積公式［24］可得

綜上，基于AFCM聚類的多目標(biāo)自適應(yīng)PR檢測方法技術(shù)路線如圖2所示。

圖2 磁偶極子多目標(biāo)自適應(yīng)PR檢測技術(shù)路線Fig.2 Technical route of multi-target adaptive PR detection for magnetic dipoles

4 磁偶極子多目標(biāo)的網(wǎng)格測量仿真

仿真環(huán)境為Matlab。在10 m×10 m的測區(qū)中，總共包含了20個(gè)磁偶極子，編號分別為1至20，其水平分布隨機(jī)，埋深范圍為-0.40 m至-0.95 m、磁矩強(qiáng)度為6 Am2至20 Am2，磁化傾角和偏角分別在［-90°，90°］和［-180°，180°］隨機(jī)取值。圖3顯示了20個(gè)磁偶極子的分布情況，預(yù)設(shè)的物性參數(shù)列于表1中。

表1 預(yù)設(shè)和估計(jì)的20個(gè)磁偶極子目標(biāo)的物性參數(shù)Tab.1 Preset and estimated physical parameters of the 20 magnetic dipoles

圖3 多源磁偶極子目標(biāo)示意圖Fig.3 Schematic diagram of multi-target magnetic dipoles

其中r0=（x0，y0，z0）為位置坐標(biāo)，M、Im和Dm分別為磁矩模、磁化傾角和偏角，磁矩m=（Mcos（Im）cos（Dm），Mcos（Im）sin（Dm），Msin（Im））。設(shè)置地磁場TMI為57 000 nT，設(shè)置測量平面z=0，坐標(biāo)軸x正向向東，y正向向北，z正向向上，此時(shí)地磁傾角和偏角分別為-55°和-5°。利用平面十字形MGTS［7］在網(wǎng)格各節(jié)點(diǎn)處差分測量MGT分量，系統(tǒng)基線距離為0.3 m，網(wǎng)格間隔為0.1 m。該平 面MGT分量、MMA、NSS和CT如圖4所示，其中（a）至（f）為測得MGT分量：（g）至（i）分別為計(jì)算的MMA、NSS、CT值?？梢?，由MMA、NSS和CT僅能發(fā)現(xiàn)部分強(qiáng)磁異常區(qū)域；磁梯度張量具有更高的分辨率，其中g(shù)zz分量識別效果最好但受磁化方向影響較大。此外，上述各參量無法找到有效等高線閾值判斷磁性體是否存在。

圖4 10 m×10 m測區(qū)平面內(nèi)MGT、MMA、NSS和CT的識別結(jié)果Fig.4 Recognition results of MGT，MMA，NSS and CT in the 10 m×10 m survey area

提出的兩個(gè)不變量改進(jìn)傾斜角θNSSTilt和θCTTilt識別結(jié)果分別如圖5（a）、（b）所示，圖5（c）、（d）為設(shè)置等高線閾值為0時(shí)圈定的目標(biāo)識別區(qū)域，圖中用白圈標(biāo)注了磁偶極子的真實(shí)位置。由圖5可知，θNSSTilt和θCTTilt有效增強(qiáng)了測區(qū)平面內(nèi)的較弱磁異常和深層磁異常目標(biāo)。識別區(qū)域剛好處于磁源位置正上方處，且區(qū)域內(nèi)節(jié)點(diǎn)數(shù)目充分，這為提高目標(biāo)定位精度和實(shí)現(xiàn)位置點(diǎn)云AFCM聚類提供了可能。

圖5 不變量改進(jìn)傾斜角θNSSTilt和θCTTilt圈定的磁偶極子識別區(qū)域（等高線閾值0）Fig.5 Magnetic dipole identification area delineated by improved tilt angles θNSSTilt and θCTTilt（contour threshold 0）

選擇θCTTilt識別結(jié)果作為目標(biāo)區(qū)域，式（7）可計(jì)算出識別區(qū)域內(nèi)各網(wǎng)格點(diǎn)對應(yīng)的目標(biāo)四組位置解。利用θCTTilt區(qū)域水平坐標(biāo)和埋深坐標(biāo)正負(fù)號篩選這些解，獲得的初始位置解集的點(diǎn)云如圖6所示，圖6（a）為估計(jì)的磁性目標(biāo)初始位置示意圖，6（b）為其俯視圖?？梢?，篩選后的解集二、三懸浮在源上方構(gòu)成了大部分被遺漏虛假解，而篩選后的解集一、四更傾向于真實(shí)解，并位于在目標(biāo)真實(shí)位置附近。但實(shí)際測量中目標(biāo)真實(shí)位置是未知的，我們無法再繼續(xù)篩選解集了，這導(dǎo)致了目標(biāo)埋深的估計(jì)值可能總體性地偏小。

圖6 識別區(qū)域內(nèi)目標(biāo)初始位置解集形成的稠密點(diǎn)云Fig.6 Dense point cloud formed by the solution set of the initial position of the target in the recognition area

提出的AFCM聚類算法對初始位置點(diǎn)云的聚類結(jié)果如圖7所示，圖中給出了估計(jì)的簇質(zhì)心位置和目標(biāo)真實(shí)位置。算法精確地估算了磁偶極子的數(shù)目，繼而獲得了各簇質(zhì)心處估計(jì)的目標(biāo)位置。表1中對比了估計(jì)的各目標(biāo)三維坐標(biāo)、磁矩強(qiáng)度和磁化方向，并統(tǒng)計(jì)了坐標(biāo)值和磁矩強(qiáng)度的總體估計(jì)精度（Total Estimation Accuracy，PTEA）作為性能指標(biāo)，表示參數(shù)真實(shí)值減去估計(jì)偏差量后與自身百分比的最小值：

其中：Pr（Pr≠0）為參數(shù)真實(shí)值；Pei為該參數(shù)第i（i=1，2，…，N）個(gè)估計(jì)值，N為目標(biāo)數(shù)量。為更真實(shí)模擬實(shí)際測量場景，在相同條件下給每個(gè)MGT分量加入了相互獨(dú)立的均值為0 nT/m、方差分別為2 nT/m和5 nT/m的高斯白噪聲。兩次測試的TEA分別記為TEAd1、TEAd2，見表1，同時(shí)給出了計(jì)算的目標(biāo)磁矩和坐標(biāo)估計(jì)偏差。結(jié)果顯示在題設(shè)噪聲環(huán)境下的目標(biāo)數(shù)目識別仍然準(zhǔn)確，該聚類結(jié)果仍具有較好的穩(wěn)定性和魯棒性。

5 多目標(biāo)小型磁鐵的網(wǎng)格測量實(shí)驗(yàn)

5.1 2.1 m×2.1 m測區(qū)五磁體定位

在中國石家莊某地平坦的地面進(jìn)行實(shí)驗(yàn)，測量區(qū)域?yàn)?.1 m×2.1 m，x軸正向向北，y軸正向向東，z軸正向向上。實(shí)驗(yàn)時(shí)間為夏季某日18：00，環(huán)境溫度28℃，系統(tǒng)為基線距離0.5 m的平面十字形MGTS［5-7］，傳感器型號Mag-03（英國Bartington），采集卡NI-6281，采樣率500 Hz。實(shí)驗(yàn)如圖8所示，在測區(qū)內(nèi)共放置序號為1～5的五個(gè)小型磁源，其磁矩強(qiáng)度和方向未知，采樣間隔為0.1 m，各網(wǎng)格點(diǎn)采樣時(shí)間不少于2 s。由于設(shè)備電磁噪聲頻段主要集中在50 Hz、200 Hz等頻段（國內(nèi)交流電場頻率），為盡量避免電磁場影響和保留真實(shí)磁異常信號，利用零相位低通濾波器進(jìn)行降噪，并將截止頻率保守地設(shè)置為5 Hz［5］。預(yù)先測得該區(qū)域內(nèi)地磁場TMI均值為53 394.36 nT，地磁傾角和偏角分別為56.33°（向下）和-20.17°（自北偏西）。

圖8 實(shí)驗(yàn)一：五個(gè)磁鐵的網(wǎng)格定位實(shí)驗(yàn)Fig.8 Experiment 1：grid positioning of five magnets

測得該平面處的MGT分量和常規(guī)特征量如圖9所示。圖9（a）至（f）為測得MGT分量；圖9（g）至（i）分別為計(jì)算的MMA、NSS、CT值。圖10為θNSSTilt和θCTTilt的區(qū)域識別效果，白圈標(biāo)注了磁體的真實(shí)位置?？梢姡逄柎盆F相對于其他目標(biāo)有更大的磁矩，導(dǎo)致1～4號的異常場幾乎被悉數(shù)淹沒。MGT和常規(guī)特征量已無法對全部磁鐵進(jìn)行有效識別，而設(shè)定閾值為0的兩個(gè)改進(jìn)傾斜角能準(zhǔn)確圈定出五個(gè)磁體的分布位置，其中θCTTilt識別效果更佳。

圖9 五個(gè)小型磁鐵的MGT、MMA、NSS和CT識別結(jié)果Fig.9 MGT，MMA，NSS and CT recognition results of five small magnets

圖10 θNSSTilt和θCTTilt圈定的5個(gè)磁鐵目標(biāo)的識別區(qū)域（等高線閾值0）Fig.10 Recognition area of 5 magnet objects delineated by θNSSTilt and θCTTilt（contour threshold 0）

而后，用張量衍生不變關(guān)系方法估計(jì)識別區(qū)域內(nèi)各網(wǎng)格節(jié)點(diǎn)處的目標(biāo)單點(diǎn)定位結(jié)果，對得到的四組解集進(jìn)行區(qū)域和象限篩選。為了獲得更高的聚類精度，在利用θCTTilt區(qū)域篩選定位解集時(shí)，將θCTTilt閾值上調(diào)至0.5，此時(shí)各區(qū)域內(nèi)單點(diǎn)定位解集形成的稠密點(diǎn)云如圖11（a）和（b）所示，AFCM聚類的結(jié)果如圖11（c）和（d）。由結(jié)果可知，經(jīng)區(qū)域識別、定位篩選的目標(biāo)解集聚集在真實(shí)位置附近，AFCM準(zhǔn)確估計(jì)了磁體數(shù)目為5，簇質(zhì)心的位置估計(jì)列于表2中，并給出了計(jì)算的目標(biāo)坐標(biāo)估計(jì)偏差。該實(shí)驗(yàn)對區(qū)域內(nèi)的五個(gè)小型磁體的位置坐標(biāo)估計(jì)偏差小于0.061 m，獲得了令人滿意的結(jié)果。

表2 定位實(shí)驗(yàn)中對九個(gè)磁性體的真實(shí)位置和估計(jì)的物性參數(shù)Tab.2 Real positions and estimated physical parameters of the nine magnets in the positioning experiment

圖11 識別區(qū)域內(nèi)的單點(diǎn)定位密集點(diǎn)云與AFCM聚類結(jié)果Fig.11 Single-point positioning dense point cloud and AFCM clustering results in the recognition area

5.2 1.2 m×1.2 m測區(qū)四磁體定位

在同地域進(jìn)行第二次實(shí)驗(yàn)，測量區(qū)域大小為1.2 m×1.2 m，時(shí)間為春季某日16：00，環(huán)境溫度17℃，測量儀器更換為五個(gè)HSF-112型（中國Huashun）磁通門傳感器搭建的平面十字形MGTS，采集卡為NI-9239，采樣率為1 kHz。實(shí)驗(yàn)如圖12所示，測器內(nèi)共放置序號6～9的四個(gè)小型磁源（6～9號磁鐵磁矩強(qiáng)度、方向和形狀均不相同［5］），磁矩強(qiáng)度和方向未知，預(yù)設(shè)坐標(biāo)列于表3中。采樣間隔為0.05 m，各節(jié)點(diǎn)采樣時(shí)間和降噪方法不變。測得該區(qū)域內(nèi)地磁總場強(qiáng)度 均 值53 407.02 nT，磁 傾 角-56.17°，磁 偏角-19.98°。

圖12 實(shí)驗(yàn)二：四個(gè)磁體的網(wǎng)格定位實(shí)驗(yàn)Fig.12 Experiment 2：grid positioning of four magnets

為了獲得更加清晰分明的區(qū)域識別結(jié)果，對張量數(shù)據(jù)網(wǎng)格進(jìn)行0.01 m間隔的三次樣條插值處理，因磁源埋深增加且各目標(biāo)間磁矩有明顯差異，將傾斜角閾值設(shè)定為0.3，MGT分量及其特征量如圖13所示，θNSSTilt和θCTTilt的區(qū)域識別如圖14所示。

圖13表明6和9號磁鐵的磁矩強(qiáng)度明顯弱于7和8號，且6和9號的磁張量場已完全被強(qiáng)磁場淹沒；7和8號的磁矩方向具有明顯差異。此時(shí)傳統(tǒng)方法已經(jīng)失去識別目標(biāo)邊界的能力，而圖14中的θNSSTilt和θCTTilt仍能清晰地識別出四個(gè)磁體的大致分布區(qū)域。根據(jù)識別結(jié)果，將等高線上調(diào)至0.4并選擇θNSSTilt識別區(qū)域篩選定位解集。各區(qū)域內(nèi)單點(diǎn)定位解集的AFCM聚類結(jié)果如圖15所示。

圖15 識別區(qū)域內(nèi)4個(gè)磁鐵的AFCM聚類結(jié)果Fig.15 AFCM clustering results of four magents in the recognition area

結(jié)果顯示在較弱磁性體磁場被淹沒情況下，兩種不變量改進(jìn)傾斜角仍能對磁體進(jìn)行區(qū)域預(yù)識別；盡管4個(gè)磁體磁矩、形狀差距較為明顯，經(jīng)區(qū)域識別、定位篩選的目標(biāo)解集仍聚集在磁體真實(shí)位置附近，AFCM準(zhǔn)確估計(jì)了磁體數(shù)目為4，估計(jì)的簇質(zhì)心位置和磁矩列于表2中；磁體8的估計(jì)磁矩強(qiáng)度數(shù)倍強(qiáng)于磁體6和9，且磁體間磁化方向差距較大，此時(shí)實(shí)驗(yàn)對區(qū)域內(nèi)的四個(gè)小型磁體的位置坐標(biāo)估計(jì)偏差仍小于0.091 m，同樣令人滿意。

6 結(jié)論

本文提出了一種利用AFCM聚類算法和張量不變量改進(jìn)傾斜角定位多目標(biāo)磁偶極子的方法。通過引入懲罰函數(shù)消除簇?cái)?shù)尋優(yōu)過程的單調(diào)遞減趨勢，設(shè)計(jì)了一種新的基于模糊緊湊度和分離度的模糊聚類有效性函數(shù)實(shí)現(xiàn)AFCM聚類。由歸一化磁源強(qiáng)度和張量縮并兩個(gè)張量不變量改進(jìn)的傾斜角可用于目標(biāo)分布區(qū)域識別和初始位置坐標(biāo)篩選，而張量衍生不變關(guān)系定位方法用于計(jì)算識別區(qū)域內(nèi)的各網(wǎng)格節(jié)點(diǎn)處對應(yīng)的四組位置解集，這些解集將在目標(biāo)真實(shí)位置附近形成稠密點(diǎn)云。設(shè)計(jì)的AFCM算法對這些點(diǎn)云進(jìn)行三維聚類并自動(dòng)估計(jì)簇質(zhì)心數(shù)目和質(zhì)心位置。

仿真結(jié)果證明了該方法的有效性和魯棒性；在2.1 m×2.1 m和1.2 m×1.2 m的網(wǎng)格測量實(shí)驗(yàn)中，該方法對埋深不超過1 m、磁矩強(qiáng)度和方向差異較大的九個(gè)小型磁鐵的坐標(biāo)估計(jì)偏差小于0.091 m。然而，由于衍生不變關(guān)系單點(diǎn)定位算法的虛假解不能完全篩除，目標(biāo)估計(jì)精度受到了一定的限制，發(fā)展更具優(yōu)勢和更高精度的磁偶極子單點(diǎn)定位方法將進(jìn)一步提高本文方法的定位性能和適用性。