覃俞璋, 宋 燕, 王冶天, 趙 丹
(上海理工大學(xué) 光電信息與計(jì)算機(jī)工程學(xué)院, 上海200093)
隨著人類空間活動(dòng)的日漸增加,殘留在地球軌道的廢棄物也越來(lái)越多,產(chǎn)生了大量的空間撞擊碎片在地球軌道上高速運(yùn)動(dòng)著,對(duì)在軌航天器的安全構(gòu)成極大威脅。 超高速撞擊過(guò)程中的彈靶材料會(huì)發(fā)生大形變、碎裂以及應(yīng)力波作用導(dǎo)致的層裂破壞現(xiàn)象。 當(dāng)撞擊薄板結(jié)構(gòu)時(shí)通常會(huì)形成碎片云;當(dāng)碰撞速度很高時(shí),彈靶材料的壓力和溫度極高,碎片云會(huì)發(fā)生烙化、汽化甚至變成離子體等物理現(xiàn)象。 因此,根據(jù)碎片云的結(jié)構(gòu)演變、質(zhì)量/動(dòng)量分布、碎片云熱力學(xué)狀態(tài)和相分布等碎片云特性展開研究,成為超高速碰撞研究中的重要方向。 因此,采用相應(yīng)的手段對(duì)超高速撞擊碎片云圖像特征(如數(shù)量、質(zhì)量和速度分布等)進(jìn)行分析研究,對(duì)于指導(dǎo)航天器防護(hù)設(shè)計(jì),保障航天器可靠的安全在軌運(yùn)轉(zhuǎn),提升航天器在空間碎片環(huán)境中的適應(yīng)能力具有重要戰(zhàn)略意義。
本文采用2 mm 彈丸,以3.06 km/s 的速度撞擊厚度為2 mm 的靶板(板間距為80 mm),模擬太空中航天器的損傷情況。 對(duì)采集的圖像(如圖1)進(jìn)行預(yù)處理、自適應(yīng)閾值分割及canny 邊緣檢測(cè),并對(duì)碎片云進(jìn)行分析,同時(shí)根據(jù)后板的損傷情況估計(jì)碎片云中產(chǎn)生的有效損傷碎片數(shù)量比例。
圖1 彈丸撞擊靶板產(chǎn)生的碎片云圖像Fig.1 The debris cloud image produced by the projectile hitting the target
本文運(yùn)用opencv 圖像分析技術(shù)對(duì)撞擊產(chǎn)生的碎片云圖像進(jìn)行預(yù)處理(包含圖像對(duì)比度拉伸、圖像銳化、自適應(yīng)平滑濾波、小波去噪等)、分割(先kmeans 聚類,然后基于區(qū)域進(jìn)行圖像分割)和特征提?。╟anny 邊緣檢測(cè)等),分析碎片云特征分布的動(dòng)態(tài)變化過(guò)程。
首先,需要對(duì)碎片云圖像進(jìn)行高斯濾波,一方面為了去除噪點(diǎn),另一方面為圖像分割和邊緣檢測(cè)做準(zhǔn)備。
高斯濾波就是一種線性平滑濾波,實(shí)用于消除高斯噪聲,并廣泛應(yīng)用于圖像處理的減噪過(guò)程。 通俗地講,高斯濾波是對(duì)整幅圖像展開加權(quán)平均的過(guò)程,每一個(gè)像素點(diǎn)的數(shù)值,均應(yīng)由其本身及鄰域之內(nèi)的其它像素值經(jīng)加權(quán)平均之后獲得。 高斯濾波的具體實(shí)施是:用一個(gè)模板(或稱卷積、掩模)掃描圖像中的每一個(gè)像素,用模板設(shè)定的鄰域內(nèi)像素的加權(quán)平均灰度值去替代模板中心像素點(diǎn)的值。
本質(zhì)上,高斯模糊就是將(灰度)圖像I 和一個(gè)高斯核進(jìn)行卷積操作:
其中: ?表示卷積操作,Gσ是標(biāo)準(zhǔn)差為σ 的二維高斯核,定義為:
計(jì)算平均值時(shí),只需把“中心點(diǎn)”作為原點(diǎn),其余點(diǎn)根據(jù)其在正態(tài)曲線上的位置,分配權(quán)重,便能夠得到一個(gè)加權(quán)平均值,而這便是上述與二維高斯核進(jìn)行卷積的過(guò)程。 獲得權(quán)重矩陣之后,中心點(diǎn)以及附近n 個(gè)點(diǎn),每一點(diǎn)乘以自身的權(quán)重值并將這些數(shù)值相加,就是中心點(diǎn)的高斯模糊的數(shù)值。 對(duì)所有點(diǎn)重復(fù)這個(gè)過(guò)程,就得到了高斯濾波后的圖像。 高斯濾波后的超高速撞擊碎片云圖像如圖2 所示。
圖2 高斯濾波處理后的碎片云圖像Fig.2 Gaussian filtered image of debris cloud
高斯濾波后,對(duì)碎片云圖像進(jìn)行對(duì)比度拉伸。對(duì)比度拉伸是圖像增強(qiáng)的一種方法,本質(zhì)上屬于灰度變換操作。 當(dāng)圖像所有像素的灰度值大部分集中在0-255 中某個(gè)很小的區(qū)間內(nèi)時(shí),整個(gè)圖像就會(huì)顯得很暗淡,對(duì)比度不高。 如果通過(guò)灰度變換,將灰度值拉伸到整個(gè)0-255 的區(qū)間,那么其對(duì)比度顯然會(huì)大幅增強(qiáng),圖片也會(huì)顯得更加明亮清晰。 可以用如下公式將某個(gè)像素的灰度值映射到更大的灰度空間
其中: Imin,Imax是原始圖像的最小灰度值和最大灰度值,MIN 和MAX 是要拉伸到的灰度空間的灰度最小值和最大值。
除上述方法,對(duì)比度拉伸還可以使用直方圖位移法,即
在每個(gè)像素位置的灰度值增加一個(gè)偏移量offset。 注意,這個(gè)offset 可以是正數(shù),也可以是負(fù)數(shù)。 若為正數(shù),整體亮度變亮;若為負(fù)數(shù),整體亮度變暗。 需要注意的是控制offset 值的大小,不要越界。
碎片云圖像處理結(jié)果如圖3 所示。 與處理前的圖像相比,經(jīng)過(guò)對(duì)比度拉伸后的圖像更加清晰,且圖像中的無(wú)關(guān)信息、雜點(diǎn)等都被消除。
圖3 進(jìn)行對(duì)比度拉伸處理后的碎片云圖像Fig.3 Debris cloud image after contrast stretching
為了減少數(shù)據(jù)量,同時(shí)提取感興趣區(qū)域(ROI),即碎片云的主體部分,需要對(duì)圖像進(jìn)行自適應(yīng)閾值圖像分割。
根據(jù)圖像上灰度值的分布,把圖像分為背景和前景兩部分,前景是按照閾值分割出來(lái)的部分,背景和前景的分界值即為待求的閾值。 遍歷不同的閾值,計(jì)算各個(gè)閾值下對(duì)應(yīng)的背景和前景之間的類內(nèi)方差,當(dāng)類內(nèi)方差取到極大值時(shí),這時(shí)對(duì)應(yīng)的閾值即為所求的閾值。
對(duì)于圖像I(x,y), 前景(即目標(biāo))和背景的分割閾值記作T,屬于前景的像素點(diǎn)數(shù)占整幅圖像的比例記為ω0,其平均灰度μ0;背景像素點(diǎn)數(shù)占整幅圖像的比例為ω1, 其平均灰度為μ1。 圖像的總平均灰度記為μ,類間方差記為g。 假設(shè)圖像的背景較暗,并且圖像的大小為M ×N,圖像中像素的灰度值小于閾值T 的像素個(gè)數(shù)記作N0,像素灰度大于閾值T 的像素個(gè)數(shù)記作N1,則有:
將式(9)代入式(10),得到等價(jià)公式:
上述即為類間方差的公式表述。 采用遍歷的方法得到使類間方差g 最大的閾值T,即為所求的閾值。
將分割后的圖像二值化后,便得到圖4 的效果。
圖4 進(jìn)行自適應(yīng)閾值分割處理后的碎片云圖像Fig.4 Fragmentation cloud image after adaptive threshold segmentation
為了便于構(gòu)建碎片云的數(shù)學(xué)模型,進(jìn)一步除去無(wú)關(guān)信息,減少數(shù)據(jù)量,需要對(duì)圖像進(jìn)行canny 邊緣檢測(cè)。
進(jìn)行圖像邊緣檢測(cè)必須滿足2 個(gè)條件:一是能有效地抑制噪聲;二必須盡量精確確定邊緣的位置。根據(jù)信噪比及定位乘積進(jìn)行測(cè)度,獲得最優(yōu)化逼近算子。 這便是Canny 邊緣算子。 類似于Marr(LoG)邊緣檢測(cè)方法,都屬于先平滑后求導(dǎo)數(shù)的方法。
該算法首先使用高斯濾波器平滑圖像,原理如下:
令g(x,y) 為平滑后的圖像,用h(x,y,σ) 對(duì)圖像f(x,y) 的平滑可表示為:
其中,?代表卷積,用一階偏導(dǎo)的有限差分來(lái)計(jì)算梯度的幅值和方向。 已平滑的梯度可以使用2×2一階有限差分近似式來(lái)計(jì)算x 與y 偏導(dǎo)數(shù)的2 個(gè)陣列和
幅值和方位角可用直角坐標(biāo)系到極坐標(biāo)的坐標(biāo)轉(zhuǎn)化公式來(lái)計(jì)算:
M[x,y] 反映了圖像的邊緣強(qiáng)度;θ[x,y] 反應(yīng)了邊緣的方向。 使M[x,y] 取得局部最大值的方向角θ[x,y],就反映了邊緣的方向。
隨后對(duì)梯度幅值進(jìn)行非極大值抑制。 因僅僅得到全局的梯度并不足以確定邊緣,因此為確定邊緣,必須保留局部梯度最大的點(diǎn),而抑制非極大值。
圖5 邊緣梯度示意圖Fig.5 Edge gradient
圖6 圓周圖和3×3 窗口Fig.6 Circumference and 3?3 windows
利用梯度的方向?qū)⑻荻冉请x散為圓周的4 個(gè)扇區(qū)之一,以便用3×3 的窗口作抑制運(yùn)算。 4 個(gè)扇區(qū)的標(biāo)號(hào)為0 到3,對(duì)應(yīng)3×3 鄰域的4 種可能組合。在每一點(diǎn)上,鄰域的中心像素M[x,y] 與沿著梯度線的2 個(gè)像素相比。 如果M[x,y] 的梯度值小于沿梯度線的2 個(gè)相鄰像素梯度值,則令M[x,y] =0。邊緣梯度示意圖如圖5 所示。 圓周圖和3×3 窗口如圖6 所示。
最后用雙閾值算法檢測(cè)和連接邊緣,如圖7 所示。 圖像既保留了碎片云的結(jié)構(gòu)屬性,又最大程度地減少了數(shù)據(jù)量,便于后續(xù)的分析計(jì)算。
圖7 進(jìn)行Canny 算法處理后的碎片云圖像Fig.7 Debris cloud image processed by Canny algorithm
本文實(shí)驗(yàn)采用的是一塊19.8 cm×19.7 cm 后板。 經(jīng)過(guò)超高速撞擊靶板后產(chǎn)生的碎片云撞擊后板而形成的彈坑,用相機(jī)拍照得到如圖8 所示的圖像。
為保證圖像中后板彈坑部分和未受撞擊部分顏色灰度值具有明顯的區(qū)別,需要對(duì)后板圖片進(jìn)行預(yù)處理(包括小波去噪、分割ROI 圖像、增加圖像對(duì)比度等)。 將后板彈坑部分和未撞擊部分進(jìn)行圖像二值化,對(duì)比可知彈坑部分與未受撞擊部分顏色具有明顯的區(qū)別,如圖9 所示。 為了突出后板損傷部分,同時(shí)減小數(shù)據(jù)量,對(duì)后板圖像進(jìn)行自適應(yīng)閾值分割。本文實(shí)驗(yàn)自適應(yīng)閾值為154,處理后圖像如圖10 所示。
圖8 碎片云撞擊后板形成的彈坑圖像Fig.8 The crater image formed by the debris cloud impacting the back plate
圖9 進(jìn)行自適應(yīng)閾值分割后的后板圖像Fig.9 Back-plate image after adaptive threshold segmentation
圖10 選擇后板圖像的感興趣區(qū)域Fig.10 Select the ROI of theback panel image
對(duì)于碎片云來(lái)說(shuō),長(zhǎng)徑比與長(zhǎng)寬比的概念相同,即經(jīng)過(guò)碎片云內(nèi)部的最長(zhǎng)徑,和與其相垂直的最長(zhǎng)徑之比。 此參數(shù)常用來(lái)表述碎片云形貌,在碎片云的判斷中具有實(shí)用價(jià)值。 由于碎片云的顆粒較細(xì)數(shù)量又多,為了更好地觀察和分析碎片云的輪廓。 根據(jù)圖7 以及計(jì)算該碎片云的長(zhǎng)徑比數(shù)值,可以準(zhǔn)確構(gòu)建如圖11 所示的數(shù)學(xué)模型。
圖11 碎片云運(yùn)動(dòng)模型Fig.11 Model of debris cloud movement
圖11 中,D1為碎片云的橫向長(zhǎng),D2為碎片云的徑長(zhǎng)。 長(zhǎng)徑比公式為:
經(jīng)過(guò)測(cè)量得到D1=1 310,D2=648,所以碎片云徑長(zhǎng)比D =2.02。
后板內(nèi)部有彈坑,經(jīng)圖像處理后,彈坑與未撞擊后板部分的顏色不同, 因此對(duì)于圖像中的后板部分像素, 只需將該圖像所有的像素個(gè)數(shù)計(jì)算出來(lái)即可。 計(jì)算公式如下:
其中,M 表示后板像素個(gè)數(shù), n 和m 分別表示后板圖像像素的行數(shù)和列數(shù)。
后板彈坑面積與后板面積的比應(yīng)等于各自像素個(gè)數(shù)之比, 因此后板彈坑面積為N?S/M。 其中,S表示后板的面積,h1和h2分別表示圖像中彈坑部分和未撞擊部分的像素值。
對(duì)圖10 中的后板彈坑面積以及后板面積的像素值總數(shù)進(jìn)行測(cè)量,得到一個(gè)1 854×1 858 的灰度矩陣, 且h1=0,h2=255 則后板像素個(gè)數(shù): M =1 854 × 1 858 =3 444 732。 同理計(jì)算彈坑像素個(gè)數(shù):N =22 731 而S =19.8 cm×19.7 cm=390.06 cm2, 因此彈坑的面積為:
后板損傷率:
結(jié)合實(shí)際情況,本文認(rèn)為后板集中損傷部分為實(shí)際的有效損傷。 如圖12 所示,則有效損傷面積計(jì)算方式如下:后板尺寸為S =19.8 cm×19.7 cm=390.06cm2,經(jīng)測(cè)量圖12 中有效損傷面積的縱深長(zhǎng)度(即紅色矩形長(zhǎng)度) R1=2.07 cm,縱深寬度(即紅色矩形寬度) R2=1.39 cm,則有效損傷面積S′ =2.07 cm×1.39 cm=2.877 3 cm。
圖12 紅色矩形為后板實(shí)際有效損傷面積Fig.12 The red rectangle is the actual effective damage area of the rear plate
如果忽略碎片云的飛行角度、速度、運(yùn)動(dòng)軌跡等物理特性,則碎片云中的每個(gè)碎片將水平飛行撞擊后板。 因此想要計(jì)算對(duì)后板產(chǎn)生有效損傷的碎片的區(qū)域面積,則需要將圖12 中有效損傷面積的縱深長(zhǎng)度R1 映射到圖7 中的碎片云圖像,將兩者結(jié)合得到圖13,圖13 中紅線所圍矩陣部分即為碎片云中產(chǎn)生有效損傷的部分。 經(jīng)測(cè)量圖13 中碎片云的面積S1=620 819(像素值),產(chǎn)生有效損傷的碎片云面積S2=261 811 (像素值),則有效損傷碎片比例rate ≈42.17%。
圖13 有效面積的縱向長(zhǎng)映射在碎片云圖像中Fig.13 The longitudinal length of the effective area is mapped in the fragment cloud image
本文利用智能算法分析超高速撞擊碎片云圖像分析的初步探索,同時(shí)結(jié)合碎片云的物理特性、飛行角度、飛行軌跡等特性,建立了一種全新的航天器損傷估計(jì)模型。 這種模型對(duì)于指導(dǎo)航天器防護(hù)設(shè)計(jì),保障航天器可靠的安全在軌運(yùn)行,提高航天器在空間碎片環(huán)境中的生存能力具有重要意義。 首先將圖像進(jìn)行預(yù)處理消除了大量冗余信息與噪聲,大大減少了數(shù)據(jù)處理量;其次對(duì)圖像進(jìn)行自適應(yīng)閾值分割和邊緣檢測(cè)操作,很好地提取碎片云的輪廓,有利于對(duì)碎片云的建模。 根據(jù)實(shí)驗(yàn)分析結(jié)果表明了本文所提算法的有效性。 這種損傷估計(jì)模型的建立不僅考慮了彈丸的物理特征,還考慮了不同時(shí)刻碎片云的圖像特征,可以預(yù)見,其將更符合實(shí)際工程需要。 由于實(shí)驗(yàn)結(jié)果還存在一些誤差,尚需對(duì)算法等進(jìn)行進(jìn)一步調(diào)整和優(yōu)化,還需要更進(jìn)一步的學(xué)習(xí)研究。