方留楊 何紅艷 張炳先

?

分塊變換和GPU并行的遙感影像快速正射校正方法

方留楊1何紅艷2,3張炳先2

（1 云南省交通規(guī)劃設(shè)計研究院有限公司 陸地交通氣象災害防治技術(shù)國家工程實驗室，昆明 650041）（2 北京空間機電研究所，北京 100094）（3 先進光學遙感技術(shù)北京市重點實驗室，北京 100094）

正射校正是整個遙感數(shù)據(jù)處理過程中計算量最大、耗時最長的步驟之一，已經(jīng)成為制約整個遙感數(shù)據(jù)處理快速完成的瓶頸。為了提高正射校正處理效率，文章系統(tǒng)地探討了基于分塊三維直接線性變換和圖形處理單元（GPU）并行的遙感影像快速正射校正方法。首先針對正射校正坐標轉(zhuǎn)換計算量過大的問題，提出了分塊三維直接線性變換策略，有效地降低了坐標轉(zhuǎn)換的計算量；在此基礎(chǔ)上，采用“漸進式”策略開展GPU并行處理，首先通過GPU并行映射（核函數(shù)任務映射、基本設(shè)置），使方法在GPU上可執(zhí)行，然后通過“兩層次”性能優(yōu)化（核函數(shù)性能優(yōu)化、整體流程性能優(yōu)化），進一步提高了方法的執(zhí)行效率。在CPU和GPU組成的實驗環(huán)境中，使用文中方法對“高分二號”衛(wèi)星全色標準景影像進行實驗，GPU執(zhí)行時間僅為5.13s，與CPU相比，相應加速比達到142.42倍，可以滿足對大數(shù)據(jù)遙感影像的快速正射校正需求。

正射校正 分塊三維直接線性變換 圖形處理單元并行映射 核函數(shù)性能優(yōu)化 整體流程性能優(yōu)化 遙感數(shù)據(jù)處理

0 引言

正射校正是遙感影像數(shù)據(jù)處理的關(guān)鍵環(huán)節(jié)，也是整個遙感數(shù)據(jù)處理過程中計算量最大、耗時最長的步驟之一。對大數(shù)據(jù)遙感影像的處理時間過長，已成為制約整個遙感數(shù)據(jù)快速處理的瓶頸，給遙感影像在災害防治和救援、應急搶險等高時效行業(yè)領(lǐng)域的應用帶來了巨大挑戰(zhàn)[1-3]。近年來，以圖形處理單元（Graphics Processing Unit，GPU）為代表的通用計算硬件逐漸成為國際上解決大數(shù)據(jù)計算和實時處理問題的主流方案，GPU專門為計算密集型、高度并行化的處理任務設(shè)計，為許多對計算速度有很高要求的問題提供了全新的解決方案[4]，其中在正射校正方面開展的研究和應用主要包括：2009年，Pro-Lines GeoImaging server軟件（由PCI Geomatics公司研發(fā)）的正射校正模塊引入了GPU并行處理技術(shù)，有效提高了其處理效率[5]，文獻[6]基于GPU頂點圖形渲染和紋理映射技術(shù)實現(xiàn)了航空影像正射校正，但該方法并未涉及GPU通用計算技術(shù)。在國內(nèi)，文獻[7]對遙感影像GPU正射校正的效率和正確性進行了驗證，說明了該技術(shù)的可行性；文獻[8]提出了基于GPU的高光譜影像幾何校正模型，但未對實現(xiàn)方法進行詳細闡述；文獻[5]探討了GPU/CPU協(xié)同的遙感影像正射糾正方法，但受實驗環(huán)境等客觀條件限制，并未對性能優(yōu)化方法進行深入分析；文獻[9]和文獻[10]先后提出了星載高分辨率SAR影像的GPU幾何校正方法；文獻[3]使用“層次性分塊”策略設(shè)計了基于CPU和GPU協(xié)同處理的正射校正方法，但上述方法僅將算法簡單移植至GPU上進行處理，并未進行充分的算法簡化和性能優(yōu)化，尚有較大的提升和完善空間。

本文首先針對傳統(tǒng)正射校正方法中逐像素坐標轉(zhuǎn)換計算量過大的問題，提出了一種分塊三維直接線性變換策略，以有效減少坐標轉(zhuǎn)換的計算量。在此基礎(chǔ)上，從GPU并行映射（核函數(shù)任務映射、基本設(shè)置）、性能優(yōu)化（核函數(shù)性能優(yōu)化、整體流程性能優(yōu)化）兩個方面出發(fā)，探討了使用GPU并行技術(shù)進行快速正射校正的方法。

1 快速正射校正方法

1.1 分塊三維直接線性變換

遙感影像正射校正分為直接法和間接法兩種方案[11]，其中間接法方案基于輸出影像（即校正后影像）進行任務劃分，根據(jù)坐標投影轉(zhuǎn)換關(guān)系、有理多項式（RPC）模型和數(shù)字高程模型（DEM），逐像素建立輸出影像物方空間坐標與輸入影像（即待校正影像）像方坐標的轉(zhuǎn)換關(guān)系，然后對輸入影像進行重采樣。該方案的并行程度較高，可以直接獲取輸出影像上每個像素的灰度值（不涉及輸出影像灰度值再分配），適合于在GPU進行并行處理，但其坐標轉(zhuǎn)換和重采樣的計算量都很大。為了有效降低坐標轉(zhuǎn)換的計算量，本文提出一種分塊三維直接線性變換策略（如圖1所示）：即對輸出影像進行格網(wǎng)劃分，僅在格網(wǎng)頂點處嚴格按照坐標投影轉(zhuǎn)換關(guān)系、RPC模型和DEM求解輸出影像和輸入影像像點的坐標對應關(guān)系，格網(wǎng)內(nèi)各像點坐標則通過近似的簡化模型建立聯(lián)系。由于描述的是輸出影像三維空間坐標和輸入影像二維像點坐標之間的關(guān)系，因此選擇三維直接線性變換模型作為簡化模型，即

式中s和l為輸入影像二維像點坐標；x、y、h為輸出影像三維空間坐標；a0~a6、b0~b6為三維直接線性變換模型參數(shù)。

三維直接線性變換模型包含14個參數(shù)，需至少列出14個方程才可求解，但每個格網(wǎng)同名頂點僅可列2個方程，4個格網(wǎng)同名頂點合計只能列8個方程，無法求解。為增加方程數(shù)，本文引入虛擬格網(wǎng)同名點的概念，即首先根據(jù)坐標投影轉(zhuǎn)換關(guān)系、RPC模型和DEM計算得到4個真實格網(wǎng)同名頂點坐標，再將DEM中真實高程分別增加和減少一個常數(shù)，保持坐標投影轉(zhuǎn)換關(guān)系和RPC模型不變，求解得到8個虛擬格網(wǎng)同名點坐標，這樣格網(wǎng)同名點數(shù)量從4個增加至12個，能列出24個方程，可以通過最小二乘法求解出14個參數(shù)，求解思路如圖2所示。

圖2 引入虛擬格網(wǎng)同名點計算三維直接線性變換模型參數(shù)

1.2 GPU并行處理

引入分塊三維直接線性變換后，雖然有效降低了坐標轉(zhuǎn)換的計算量，但卻額外引入了變換模型參數(shù)擬合的計算量，并且重采樣的計算量并未減少，因此整個正射校正的計算量仍然很大，需要使用GPU并行處理技術(shù)，進一步提高方法的執(zhí)行效率。本文選擇目前主流的統(tǒng)一設(shè)備計算架構(gòu)（CUDA）作為GPU并行處理框架[12]，采用“漸進式”策略開展GPU并行處理：首先通過GPU并行映射，明確核函數(shù)執(zhí)行流程和基本設(shè)置，確保方法在GPU上可執(zhí)行；然后通過“兩層次”性能優(yōu)化，進一步提高方法在GPU上的執(zhí)行效率。

1.2.1 GPU并行映射

（1） 核函數(shù)任務映射

首先進行核函數(shù)任務映射，以確定GPU上正射校正執(zhí)行的邏輯流程。由于格網(wǎng)同名點坐標轉(zhuǎn)換計算、三維直接線性變換模型參數(shù)計算、影像重采樣計算都具有逐格網(wǎng)/像素獨立的特點，具有很高的并行性，因此，基于最大限度的GPU端并行執(zhí)行原則[13]，將上述3個步驟都映射至GPU執(zhí)行。考慮到3個步驟的邏輯連續(xù)，為減少反復調(diào)用核函數(shù)帶來的額外開銷，進一步將3個步驟分別映射為3個獨立的核函數(shù)GC_CorPtsCal、GC_Dlt3DCal、GC_Resample實現(xiàn)，映射關(guān)系如圖3所示。

圖3 正射校正核函數(shù)任務映射

（2） 基本設(shè)置

進行核函數(shù)任務映射后，接下來需要進行基本設(shè)置，使算法在GPU上可執(zhí)行，具體包括：

1）將每個CUDA線程設(shè)置為僅對一個元素進行計算。具體來說，對于GC_CorPtsCal核函數(shù)，每個CUDA線程負責一個格網(wǎng)同名點坐標轉(zhuǎn)換的計算；對于GC_Dlt3DCal核函數(shù)，每個CUDA線程負責一組格網(wǎng)三維直接線性變換模型參數(shù)的計算；對于GC_Resample核函數(shù)，每個CUDA線程負責對一個像素進行重采樣。

2）將核函數(shù)線程塊大小設(shè)置為256像素（16像素×16像素）。

3）為減少對全局存儲器的多次訪問，使用程序默認的L1 cache緩存模式。

4）為減少使用循環(huán)指令帶來的開銷，使用程序默認的循環(huán)展開模式。

5）坐標投影轉(zhuǎn)換參數(shù)和RPC模型參數(shù)較少且為常數(shù)值，滿足常數(shù)存儲器的存儲容量和對數(shù)據(jù)只讀的要求，因此將其載入常數(shù)存儲器中，以提高參數(shù)的訪問效率。

1.2.2 GPU性能優(yōu)化

GPU性能優(yōu)化分為兩個層次進行：核函數(shù)性能優(yōu)化和整體流程性能優(yōu)化。首先對3個核函數(shù)進行性能優(yōu)化，以提高核函數(shù)的執(zhí)行效率；然后進行整體流程性能優(yōu)化，通過傳輸計算堆疊隱藏CPU和GPU的數(shù)據(jù)傳輸開銷，以進一步提高方法的整體執(zhí)行效率。

（1） 核函數(shù)性能優(yōu)化

對核函數(shù)進行性能優(yōu)化的方法包括：線程多元素重訪、線程塊大小和形狀最優(yōu)選取、存儲層次性訪問[14]。

線程多元素重訪基于GPU網(wǎng)格和線程尋址指令可以多次使用以及存儲訪問流水化的特性，同時分配多個元素（格網(wǎng)或像素）處理任務給單個GPU線程處理，以提高任務的整體處理效率[14]，但線程多元素重訪的效率受數(shù)據(jù)劃分模式、負載大小等因素的影響，對不同核函數(shù)的性能提升程度并不相同，因此需通過實驗分別確定3個核函數(shù)GC_CorPtsCal、GC_Dlt3DCalGC_Resample中各線程處理元素的具體數(shù)量，使算法執(zhí)行效率達到最優(yōu)。此外，隨著線程處理元素數(shù)量的不斷增長，若待處理數(shù)據(jù)不變，在線程塊大小固定的情況下，線程塊的數(shù)量將會逐漸減少，在該種情況下，需要保證給每個流式多處理器中至少分配1個線程塊，以充分利用GPU的計算性能。

線程塊大小選取時需考慮以下因素：1）線程塊的大小應至少在64線程以上；2）線程塊最大不能超過1 024線程；3）線程塊的大小應該是線程束數(shù)量（32線程）的整數(shù)倍，以避免在未完全展開的線程束上進行無意義的操作。此外，由于線程塊的寬度會影響全局存儲器的合并訪問效率，并且GC_CorPtsCalGC_Dlt3DCalGC_Resample三個核函數(shù)的全局存儲器訪問模式各不相同，因此將各線程塊的寬度從1線程增加到最大，并記錄3個核函數(shù)的執(zhí)行時間，以尋找核函數(shù)執(zhí)行時間最短時對應的線程塊寬度。

存儲層次性訪問利用GPU的多級緩存提高方法的訪存效率。通用計算GPU提供了兩級緩存，在基本設(shè)置中，核函數(shù)采用了默認的一級緩存（L1 cache）進行數(shù)據(jù)緩存，但二級緩存（L2 cache）同樣可進行數(shù)據(jù)緩存。此外，GC_Resample核函數(shù)進行影像重采樣時的數(shù)據(jù)內(nèi)插具有明顯的二維局部性訪存特點，該特點適合于使用紋理存儲器進行優(yōu)化。因此在實驗部分，分別記錄核函數(shù)GC_CorPtsCal、GC_Dlt3DCal采用L1 cache、L2 cache時的執(zhí)行時間，以及核函數(shù)GC_Resample采用L1 cache、L2 cache、紋理存儲器時的執(zhí)行時間，確保將核函數(shù)性能調(diào)至最優(yōu)。

（2） 整體流程性能優(yōu)化

通用計算GPU內(nèi)部固化有一個數(shù)據(jù)拷貝引擎和一個核函數(shù)執(zhí)行引擎[15]，而基于分塊變換的GPU正射校正也包含1次內(nèi)存至顯存的數(shù)據(jù)拷貝操作、3次核函數(shù)執(zhí)行操作（GC_CorPtsCal、GC_Dlt3DCal、GC_Resample）和1次顯存至內(nèi)存的數(shù)據(jù)拷貝操作。因此可使用流處理機制，將2次數(shù)據(jù)拷貝操作（內(nèi)存至顯存，顯存至內(nèi)存）分配給數(shù)據(jù)拷貝引擎，將3次核函數(shù)執(zhí)行操作分配給核函數(shù)執(zhí)行引擎，并指定多組流操作，同時進行多組數(shù)據(jù)拷貝和核函數(shù)執(zhí)行，以提高兩個引擎的使用效率。但由于對大數(shù)據(jù)量遙感影像進行處理時，有時受GPU顯存容量的限制，需要進行多輪數(shù)據(jù)拷貝操作，如果在數(shù)據(jù)拷貝過程中出現(xiàn)數(shù)據(jù)疊掩和重復拷貝現(xiàn)象，會降低流的整體執(zhí)行效率，因此需要通過具體實驗確定流的數(shù)量，使GPU上正射校正的執(zhí)行效率達到最優(yōu)。

綜上所述，分塊變換和GPU并行的遙感影像快速正射校正方法流程如圖4所示。

圖4 分塊變換和GPU并行的快速正射校正流程

2 實驗結(jié)果與分析

2.1 實驗數(shù)據(jù)及平臺

2014年8月19日，我國在太原衛(wèi)星發(fā)射中心成功發(fā)射“高分二號”衛(wèi)星，該衛(wèi)星搭載了2臺分辨率為0.8m的全色相機和1臺分辨率為3.2m的多光譜相機，可提供多種產(chǎn)品和應用服務，是迄今為止中國研制的空間分辨率最高的民用低軌遙感衛(wèi)星[16]。其中全色相機標準景影像大小為29 200像素×27 600像素，采用10bit量化，數(shù)據(jù)量較大（約1.5 Gbyte），具有明顯代表性，因此本節(jié)選擇該影像作為實驗數(shù)據(jù)。

本文實驗平臺的CPU為Intel Xeon E5650，主頻2.66 GHz，含6個處理核心；GPU為NVIDIA Tesla C2050，含448個處理核心，顯存容量為3 Gbyte，內(nèi)部含有一個數(shù)據(jù)拷貝引擎和一個核函數(shù)執(zhí)行引擎；主機內(nèi)存為32 Gbyte。

2.2 方法性能

為定量評價算法性能，本節(jié)定義加速比評價加速后算法與加速前算法相比的速度提升程度，計算公式為

式中bef為加速前算法的執(zhí)行時間；aft為加速后算法的執(zhí)行時間；為加速比。加速比越高，說明算法的執(zhí)行效率越優(yōu)。

進行正射校正時需要同時為輸入影像和輸出影像分配GPU顯存空間，但由于實驗影像數(shù)據(jù)量太大，無法一次性全部載入，因此需要分塊處理。本文按照行優(yōu)先原則進行分塊，將每次處理影像塊的寬度固定為整幅輸出影像的寬度，每次處理影像塊的高度則按照式（3）計算：

式中為算法分配給輸出影像的內(nèi)存空間，此處為512 Mbyte；d為影像數(shù)據(jù)類型需要占用的字節(jié)數(shù)，全色影像為10bit量化，因此其值為2；為輸出影像寬度，其值為34 100像素。

將式（3）各項數(shù)值代入計算，得到每次處理的輸出影像高度為7 872像素。由于整幅輸出影像的高度為35 000像素，因此為了對整幅輸出影像進行處理，正射校正需分為5次進行。此外，對算法作如下約定：輸出影像坐標系為UTM投影坐標系，因此坐標轉(zhuǎn)換通過UTM投影和RPC模型進行；格網(wǎng)大小為64像素×64像素；重采樣方式為雙線性內(nèi)插。另外，為了專注于對方法本身效率進行討論和分析，本節(jié)不涉及磁盤上數(shù)據(jù)的讀/寫（I/O）時間，僅討論方法本身的執(zhí)行時間（包括各步驟計算時間，以及CPU和GPU間的數(shù)據(jù)傳輸時間）。

首先，在Intel Xeon E5650 CPU上，記錄對“高分二號”衛(wèi)星全色標準景影像進行傳統(tǒng)逐點坐標轉(zhuǎn)換正射校正和分塊三維直接線性變換正射校正的執(zhí)行時間，結(jié)果如表1所示。由于UTM投影轉(zhuǎn)換的計算量很大，因此逐點坐標轉(zhuǎn)換正射校正中同名點坐標轉(zhuǎn)換計算的耗時很長（471.06s），占整個校正時間（730.61s）的64.5%；進行分塊三維直接線性變換優(yōu)化后，雖然額外引入了三維直接線性變換模型參數(shù)計算的時間（53.41s），且重采樣的計算時間相同（出現(xiàn)微小區(qū)別的原因在于程序執(zhí)行時硬件性能的微小波動和計時誤差），但由于同名點坐標轉(zhuǎn)換計算時間（1.38s）大幅減少，因此整體上正射校正的執(zhí)行時間（314.24s）仍顯著縮短，相應加速比達到2.33倍。

表1 兩種正射校正方法CPU逐步驟執(zhí)行時間和加速比

Tab.1 CPU step-by-step run times and speedup ratios of two ortho-rectification methods

然后，通過GPU并行處理，進一步提高分塊三維直接線性變換正射校正的執(zhí)行效率。表2列出了GPU并行映射后正射校正各個步驟的執(zhí)行時間，包括CPU→GPU數(shù)據(jù)傳輸時間、核函數(shù)執(zhí)行時間、GPU→CPU數(shù)據(jù)傳輸時間以及其余項耗時（包括了顯存分配、核函數(shù)啟動、計時開銷等操作的時間開銷）?？梢园l(fā)現(xiàn)，進行了GPU并行映射后，正射校正處理效率得到了顯著提升，3個核函數(shù)的加速比均達到了100倍以上，但由于額外增加了CPU和GPU之間雙向數(shù)據(jù)傳輸和其余耗時項的執(zhí)行時間開銷，因此整體加速比為55.91倍，相應的執(zhí)行時間5.62s。此外，對數(shù)據(jù)傳輸時間進行分析，發(fā)現(xiàn)CPU→GPU數(shù)據(jù)傳輸時間（1.07s）長于GPU→CPU數(shù)據(jù)傳輸時間（0.67s），這是由于衛(wèi)星特有的成像特點使得正射校正算法在對某一輸出影像塊進行重采樣時，需要的輸入影像塊范圍比輸出影像塊大。由于該現(xiàn)象的存在，導致正射校正多次分塊進行時，各輸出影像塊對應的輸入影像塊出現(xiàn)重疊，從而使得CPU→GPU數(shù)據(jù)傳輸量（3.58 Gbyte）大于GPU→CPU數(shù)據(jù)傳輸量（2.22 Gbyte），進而造成CPU→GPU數(shù)據(jù)傳輸時間長于GPU→CPU數(shù)據(jù)傳輸時間。

表2 GPU并行映射后方法各步驟執(zhí)行時間和加速比

Tab.2 Run times and speedup ratios of ortho-rectification after GPU parallel mapping

采用1.2.2節(jié)中設(shè)計的性能優(yōu)化策略對上述核函數(shù)進行性能優(yōu)化。首先使用線程多元素重訪策略，令一個CUDA線程處理多個元素，相應的核函數(shù)執(zhí)行時間如圖5所示。從圖5中可知，對于3個核函數(shù)GC_CorPtsCal、GC_Dlt3DCal、GC_Resample，當每個CUDA線程分別處理8、16、128個元素時，核函數(shù)執(zhí)行時間取得最小值，分別為9.8ms、366.2ms、1 331.6ms。進一步發(fā)現(xiàn)，GC_Resample的執(zhí)行時間變化較為平穩(wěn)，而GC_CorPtsCal、GC_Dlt3DCal的執(zhí)行時間在每個CUDA線程處理16個元素后開始顯著增加，造成該現(xiàn)象的原因是：在開始實驗時，已確定分塊格網(wǎng)大小為64像素×64像素，且線程塊初始大小設(shè)定為16線程×16線程，因此當每個CUDA線程負責計算4×4個（即16個）格網(wǎng)點坐標和三維直接線性變換模型參數(shù)時，每個線程塊的覆蓋范圍即可達到4 096像素×4 096像素，要覆蓋每次處理的影像塊（34 100像素×7 872像素），僅需9×2個線程塊即可；NVIDIA Tesla C2050 GPU共有14個流式多處理器，要充分發(fā)揮其計算性能，需要至少給每個流式多處理器分配1個線程塊，此時恰好可以滿足要求，若再將每個CUDA線程處理的元素量增加，線程塊數(shù)量會更少，已無法保證給每個流式多處理器分配一個線程塊，在這種情況下，GPU的硬件計算資源沒有得到充分利用，因此核函數(shù)的執(zhí)行效率明顯下降，對應的執(zhí)行時間顯著增加。

接下來討論使用不同大小、形狀和緩存模式的線程塊時，3個核函數(shù)執(zhí)行時間的變化情況，其中核函數(shù)GC_CorPtsCal的實驗結(jié)果繪制于圖6（a），核函數(shù)GC_Dlt3DCal、GC_Resample的實驗結(jié)果繪制于圖6（b）。柱狀子圖記錄核函數(shù)執(zhí)行時間（該執(zhí)行時間為線程塊寬度從1增加到最大時的最短執(zhí)行時間），折線子圖記錄對應的線程塊寬度。觀察圖6中所示的柱狀子圖及折線子圖可以發(fā)現(xiàn)：

1）改變線程塊大小對核函數(shù)GC_CorPtsCal、GC_Resample的執(zhí)行時間無明顯影響。這是因為這兩個核函數(shù)涉及的坐標轉(zhuǎn)換和重采樣的計算邏輯較為復雜，其計算指令遠多于訪存指令（將RPC模型載入常數(shù)存儲器后，GC_CorPtsCal甚至不需要訪存），計算訪存比高。因此，在線程塊逐漸增大過程中線程占有率最低時，流式多處理器中駐留的線程束的數(shù)量也足夠隱藏訪問延遲[17]。此外，發(fā)現(xiàn)當線程塊小于384線程時，改變線程塊大小對GC_Dlt3DCal的執(zhí)行時間無明顯影響，原因同上；但當線程塊大于384線程時，GC_Dlt3DCal的執(zhí)行時間明顯上升，這是因為當線程塊大于384線程時，線程塊數(shù)量已減少至無法保證給每個流式多處理器分配一個線程塊，導致GPU的硬件計算資源沒有得到充分利用，因此核函數(shù)的執(zhí)行效率明顯下降。

2）采用不同的緩存模式時，3個核函數(shù)的執(zhí)行時間無明顯差別。雖然核函數(shù)GC_Resample中存在著明顯的二維局部訪存模式，且該訪存模式非常適合于使用紋理存儲器進行加速，但在性能上紋理存儲器并未體現(xiàn)出對于L1/L2 cache的優(yōu)勢。

3）當核函數(shù)執(zhí)行時間取得最小值時，GC_Dlt3DCal和GC_Resample的線程塊寬度總體上小于GC_CorPtsCal的線程塊寬度，這是因為GC_Dlt3DCal和GC_Resample的訪存模式較為復雜，其中GC_Dlt3DCal訪問的是存儲控制點坐標的結(jié)構(gòu)體，因此每個線程訪問的顯存地址并不相鄰，而是間隔了結(jié)構(gòu)體大小的地址范圍；GC_Resample的訪存無固定模式，由于重采樣反算的不規(guī)則性，不同行的輸出影像反算到輸入影像后有可能對應同一行，使得顯存訪問被合并，因此GC_Dlt3DCal和GC_Resample的執(zhí)行時間在線程塊寬度較小時即可取到最優(yōu)值。GC_CorPtsCal不存在顯存訪問操作，因此核函數(shù)執(zhí)行時間基本不受線程塊寬度增加的影響，此處其執(zhí)行時間在線程塊寬度相對較大時取到最優(yōu)值。

圖6 改變線程塊大小、形狀和緩存模式時核函數(shù)執(zhí)行時間

表3列出了核函數(shù)性能優(yōu)化后的配置、執(zhí)行時間和性能提升比?？梢钥闯觯?jīng)過性能優(yōu)化后，3個核函數(shù)的性能都得到不同程度的提升，其中GC_Resample的性能提升比最高，達到了18.1%。

表3 核函數(shù)性能優(yōu)化后的配置、執(zhí)行時間和性能提升比

Tab.3 Configurations, run times and performance improvement ratios after kernel performance optimization

完成核函數(shù)性能的優(yōu)化配置后須進行整體流程性能優(yōu)化，通過流機制實現(xiàn)傳輸計算堆疊，進一步提高整個正射校正的執(zhí)行效率。表4為分別創(chuàng)建2、3、4個流時，方法執(zhí)行時間和性能提升比。從表中可知，創(chuàng)建2個流和3個流時，正射校正的性能得到小幅提升，其中創(chuàng)建2個流時性能提升比最高（6.45%），創(chuàng)建4個流時，性能不但沒有提升，反而下降（?7.05%）。造成該現(xiàn)象的原因是，當創(chuàng)建多個流時，CPU→GPU數(shù)據(jù)傳輸量多于不使用流時。如圖7所示，當不創(chuàng)建流時，對輸出子影像（圖7（b）中影像藍色框標識區(qū)域）進行重采樣需要從顯存載入對應的輸入子影像（圖7（a）中影像藍色框標識區(qū)域）；當使用2個流時，每個流負責對一半的輸出子影像（圖7（b）中影像黃色和綠色標識區(qū)域）進行重采樣，但此時輸入影像會發(fā)生疊掩（圖7（a）中影像黃色和綠色標識區(qū)域發(fā)生重合）。這種由于使用流而引入的冗余數(shù)據(jù)傳輸極大影響了正射校正的執(zhí)行效率，并且流的數(shù)量越多，冗余傳輸?shù)臄?shù)據(jù)量越大，因此算法性能不僅沒有進一步提升，反而出現(xiàn)下降。

表4 整體性能優(yōu)化后的方法執(zhí)行時間和性能提升比

Tab.4 Run times and performance improvement ratios after overall performance optimization

最終經(jīng)過整體性能優(yōu)化，GPU上分塊三維直接線性變換正射校正的執(zhí)行時間壓縮至5.13s，與CPU逐點坐標轉(zhuǎn)換正射校正方法、CPU分塊三維直接線性變換正射校正方法相比，加速比分別達到142.42倍、61.26倍（如表5所示），可以滿足對大數(shù)據(jù)遙感影像的快速處理需求。進一步將文獻[9]、[10]中提出的GPU并行處理方法應用至本文，由于僅進行了GPU基本映射，未進行性能優(yōu)化，因此其執(zhí)行時間僅為5.62s；文獻[3]在文獻[9]、[10]的基礎(chǔ)上進行了線程塊大小和形狀的選取優(yōu)化，因此其執(zhí)行時間提升至5.51s。綜合比較，由于本文提出的方法進行了更加充分、全面的性能優(yōu)化，因此執(zhí)行效率更優(yōu)。

表5 三種方法的執(zhí)行時間和加速比

Tab.5 Run times and speedup ratios of three ortho-rectification methods

3 結(jié)束語

本文系統(tǒng)地探討了基于分塊三維直接線性變換和GPU并行的遙感影像快速正射校正方法。首先針對正射校正坐標轉(zhuǎn)換計算量過大的問題，提出了分塊三維直接線性變換策略，有效地降低了坐標轉(zhuǎn)換的計算量，在求解三維直接線性變換模型參數(shù)時，引入了虛擬同名點的概念，通過增加方程個數(shù)確保模型參數(shù)可解。然后，采用“漸進式”策略開展GPU并行處理，提高了方法的執(zhí)行效率。在Intel Xeon E5650 CPU和NVIDIA Tesla C2050 GPU組成的實驗環(huán)境中，使用本文方法對“高分二號”衛(wèi)星全色影像進行實驗，發(fā)現(xiàn)引入分塊三維直接線性變換后，CPU上正射校正執(zhí)行時間從730.61s壓縮至314.24s，加速比為2.33倍；進一步將算法映射至GPU并進行性能優(yōu)化后，執(zhí)行時間縮短至5.13s，最終加速比達到142.42倍，可以滿足對大數(shù)據(jù)遙感影像的快速正射校正需求。

隨著高性能計算技術(shù)的快速發(fā)展，GPU計算架構(gòu)和型號會持續(xù)演變更新，新的計算架構(gòu)加入了一些新功能（例如動態(tài)并行、線程格網(wǎng)管理等），使得并行處理方式更加靈活，但GPU底層的硬件設(shè)計、執(zhí)行框架（線程組織）、層次性存儲結(jié)構(gòu)未發(fā)生明顯改變。因此，本文方法可以為基于新計算架構(gòu)的GPU進行正射校正提供借鑒。此外，近年來，眾核CPU、加速處理器（APU）等高性能處理硬件已逐漸發(fā)展成熟，因此可以研究正射校正在上述硬件上的并行處理方法，并逐步探索其余遙感數(shù)據(jù)處理瓶頸算法的并行處理技術(shù)。

[1] 楊靖宇. 攝影測量數(shù)據(jù)GPU并行處理若干關(guān)鍵技術(shù)研究[D]. 鄭州: 信息工程大學, 2011. YANG Jingyu. Study on Parallel Processing Technologies of Photogrammetry Data Based on GPU[D]. Zhengzhou: Information Engineering University, 2011. (in Chinese)

[2] 蔣艷凰, 楊學軍, 易會戰(zhàn). 衛(wèi)星遙感圖像并行幾何校正算法研究[J]. 計算機學報, 2004, 27(7): 944-951. JIANG Yanhuang, YANG Xuejun, YI Huizhan. Parallel Algorithm of Geometrical Correction for Satellite Images[J]. Chinese Journal of Computers, 2004, 27(7): 944-951. (in Chinese)

[3] 方留楊, 王密, 李德仁. CPU和GPU協(xié)同處理的光學衛(wèi)星遙感影像正射校正方法[J]. 測繪學報, 2013, 42(5): 668-675. FANG Liuyang, WANG Mi, LI Deren. A CPU-GPU Co-processing Orthographic Rectification Approach for Optical Satellite Imagery[J]. Acta Geodaetica et Cartographica Sinica, 2013, 42(5): 668-675. (in Chinese)

[4] 于英, 張元源, 薛武. 基于CUDA的核線影像并行生成技術(shù)[J]. 測繪通報, 2013(7): 27-29. YU Ying, ZHANG Yuanyuan, XUE Wu. Epipolar Image Parallel Generation Technology Based on CUDA[J]. Bulletin of Surveying and Mapping, 2013(7): 27-29. (in Chinese)

[5] 楊靖宇, 張永生, 李正國, 等. 遙感影像正射糾正的GPU-CPU協(xié)同處理研究[J]. 武漢大學學報·信息科學版, 2011, 36(9): 1043-1046. YANG Jingyu, ZHANG Yongsheng, LI Zhengguo, et al. GPU-CPU Cooperate Processing of RS Image Ortho-rectification[J]. Geomatics and Information Science of Wuhan University, 2011, 36(9): 1043-1046. (in Chinese)

[6] JAVIER R S, MARIA C C, JULIO M H. Geocorrection for Airborne Pushbroom Imagers[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(11): 4409-4419.

[7] 侯毅, 沈彥男, 王睿索, 等. 基于GPU的數(shù)字影像的正射糾正技術(shù)的研究[J]. 現(xiàn)代測繪, 2009, 32(3): 10-12. HOU Yi, SHEN Yannan, WANG Ruisuo, et al. The Discussion of GPU-based Digital Differential Rectification[J]. Modern Surveying and Mapping, 2009, 32(3): 10-12. (in Chinese)

[8] 李娜, 白勇, 趙慧潔, 等. 基于多核CPU和GPU的高光譜數(shù)據(jù)并行幾何校正[J]. 現(xiàn)代電子技術(shù), 2013, 36(2): 110-112. LI Na, BAI Yong, ZHAO Huijie, et al. Hyperspectral Data Parallel Geometry Correction Based on Multicore CPU and GPU[J]. Modern Electronics Technique, 2013, 36(2): 110-112. (in Chinese)

[9] 侯明輝, 邱虎. 基于GPU的高分辨率星載SAR系統(tǒng)幾何校正[J]. 現(xiàn)代雷達, 2013, 35(9): 40-44. HOU Minghui, QIU Hu. High Resolution Space-borne SAR Geometric Correction Based on GPU[J]. Modern Radar, 2013, 35(9): 40-44. (in Chinese)

[10] 楊景輝, 程春泉, 張繼賢, 等. GPU支持的SAR影像幾何校正大規(guī)模并行處理[J]. 中國圖像圖形學報, 2015, 20(3): 374-385. YANG Jinghui, CHENG Chunquan, ZHANG Jixian, et al. GPU Supported Massively Parallel Processing for Geometric Correction of SAR Imagery[J]. Journal of Image and Graphics, 2015, 20(3): 374-385. (in Chinese)

[11] 李德仁, 王樹根, 周月琴. 攝影測量與遙感概論[M]. 2版. 北京: 測繪出版社, 2008. LI Deren, WANG Shugen, ZHOU Yueqin. An Introduction to Photogrammetry and Remote Sensing[M]. 2nd ed. Beijing: Surveying and Mapping Press, 2008. (in Chinese)

[12] KIRK D B, HWU W M W. Programming Massively Parallel Processors[M]. Amsterdam: Elsevier, 2010.

[13] NVIDIA. CUDA C Best Practices Guide, v5.0[Z]. Santa Clara: NVIDIA Corporation, 2012.

[14] NVIDIA. Fermi Hardware and Performance Tips[Z]. Santa Clara: NVIDIA Corporation, 2011.

[15] SANDES J, KANDROT E. An Introduction to General-purpose GPU Programming[M]. Boston: Addison-wesley Professional Press, 2010.

[16] 潘騰, 關(guān)暉, 賀瑋. “高分二號”衛(wèi)星遙感技術(shù)[J]. 航天返回與遙感, 2015, 36(4): 16-24. PAN Teng, GUAN Hui, HE Wei. GF-2 Satellite Remote Sensing Technology[J]. Space Recovering & Remote Sensing, 2015, 36(4): 16-24.

[17] VOLKOV V. Better Performance at Lower Occupancy[EB/OL]. [2018-03-27]. http://www.nvidia.com/content/GTC-2010/ pdfs/2238_GTC2010.pdf.

A Fast Remote Sensing Image Ortho-rectification Method Based on Block-wise Transformation and GPU Parallel Processing

FANG Liuyang1HE Hongyan2,3ZHANG Bingxian2

（1 National Engineering Laboratory for Surface Transportation Weather Impacts Prevention, Broadvision Engineering Consultants, Kunming 650041, China）（2 Beijing Institute of Space Mechanics & Electricity, Beijing 100094, China）（3 Beijing Key Laboratory of Advanced Optical Remote Sensing Technology, Beijing 100094, China）

Ortho-rectification is one of the most computationally and time-consuming steps in remote sensing data processing, and has become a bottleneck restricting the rapid completion of remote sensing data processing. Therefore, a fast ortho-rectification method, based on block-wise three-dimensional linear transformation and GPU parallel processing, is proposed in this paper. First, a three-dimensional linear transformation is introduced to reduce the computational amount of coordinate projection conversion. Then, based on progressive strategy, the GPU parallel mapping (including kernel mapping and basic settings) is conducted. Furthermore, two performance optimization methods, i.e., kernel performance optimization and overall performance optimization, are used to further improve the efficiency. The approach is tested with GF-2 panchromatic images with CPU and GPU. The experimental results show that the processing time is only 5.13 s, with a speedup ratio of 142.42, which can meet the requirement of fast ortho-rectification for big-data-volume remote sensing images.

ortho-rectification; block-wise three-dimensional linear transformation; GPU parallel mapping; kernel performance optimization; overall performance optimization; processing of remote sensing image data

P237

A

1009-8518(2018)06-0080-11

10.3969/j.issn.1009-8518.2018.06.010

2018-03-27

國家自然科學基金（41601476）

方留楊，男，1987年生，2015年獲武漢大學測繪遙感信息工程國家重點實驗室攝影測量與遙感博士學位，高級工程師。主要從事遙感數(shù)據(jù)處理、交通災害防治應用等研究工作。E-mail：fangliuyang@stwip.com。

（編輯：夏淑密）