一種并行方式柵格影像切片方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明屬于地理空間信息處理技術(shù)領(lǐng)域，涉及一種并行方式柵格影像切片方法。 技術(shù)背景
[0002] 隨著衛(wèi)星傳感器技術(shù)以及無(wú)人機(jī)航拍技術(shù)的快速發(fā)展，遙感影像的空間和時(shí)間分 辨率都有大幅度地提高，單幅遙感影像文件的數(shù)據(jù)量也急劇增加。當(dāng)前主流GIS軟件以及互 聯(lián)網(wǎng)地圖應(yīng)用在WebGIS(網(wǎng)絡(luò)地理信息系統(tǒng))解決方案中都廣泛采用地圖切片（Tile，又稱(chēng) 瓦片）的發(fā)布策略，這種方法通過(guò)預(yù)先將原始影像按照一定的切分規(guī)則切割成一張張大小 相同的瓦片，以加載小數(shù)據(jù)量瓦片的方式達(dá)到在網(wǎng)絡(luò)帶寬受限的條件下實(shí)現(xiàn)影像的高效顯 示。但是目前的商業(yè)GIS軟件其地圖切片的生成以及發(fā)布成本較高，操作較復(fù)雜，以行業(yè)內(nèi) 知名的ArcGIS系列軟件為例，要將柵格影像進(jìn)行切片不僅要安裝ArcGIS Desktop，而且還 要安裝龐大的ArcGIS Server，并且操作十分復(fù)雜繁瑣，極大增加時(shí)間和物力成本。最重要 一點(diǎn)是隨著單幅遙感影像文件的分辨率以及數(shù)據(jù)量的急劇增加，其相應(yīng)的切片數(shù)量會(huì)呈現(xiàn) 幾何級(jí)數(shù)式的增長(zhǎng)，傳統(tǒng)的串行算法以及商業(yè)GIS軟件是通過(guò)單機(jī)預(yù)先切好瓦片，再對(duì)外提 供，這種傳統(tǒng)切片技術(shù)計(jì)算資源利用低下，并且沒(méi)有錯(cuò)誤恢復(fù)機(jī)制，一旦切片過(guò)程中出現(xiàn)故 障，整個(gè)切片工作得推倒重來(lái)，無(wú)法在原有的進(jìn)度上繼續(xù)進(jìn)行。因此在硬件性能高速發(fā)展的 情況下，如何利用高性能計(jì)算技術(shù)，方便快捷并且高效快速地進(jìn)行柵格影像切片，已經(jīng)成為 WebG IS地圖應(yīng)用中快速可視化方面必須要解決的重要問(wèn)題。
[0003] 柵格影像切片是指將指定地理范圍的柵格影像，在某一級(jí)別比例尺下，基于一定 的切割準(zhǔn)則切割成若干行和列的固定尺寸的正方形柵格圖片的過(guò)程，這些規(guī)整的影像切片 又稱(chēng)為瓦片。每個(gè)瓦片在該級(jí)別比例尺下都有一個(gè)獨(dú)立的編碼，該編碼由瓦片的行列號(hào)以 及比例尺的級(jí)別構(gòu)成。當(dāng)瀏覽器前端在影像瀏覽時(shí)，根據(jù)當(dāng)前的地理范圍動(dòng)態(tài)地計(jì)算所需 顯示的影像切片的行列號(hào)，通過(guò)行列號(hào)以及當(dāng)前比例尺級(jí)別得到瓦片的編碼，而后通過(guò)編 碼向服務(wù)器端請(qǐng)求相應(yīng)瓦片進(jìn)行顯示，從而達(dá)到快速響應(yīng)目的。
[0004] 目前柵格影像并行切片方法主要有三種思路:一種是CPU(Central Processing Unit，中央處理器)+GPU(Graphic Processing Unit，圖形處理單元）的方式，利用GPU的計(jì) 算能力進(jìn)行并行加速，這種方法并行能力受限于GPU硬件的能力，并行程度有限，而且會(huì)提 高系統(tǒng)架構(gòu)成本;另一種是利用分布式集群系統(tǒng)，將切片任務(wù)劃分為多個(gè)子任務(wù)，各子任務(wù) 在多個(gè)分布式節(jié)點(diǎn)上同時(shí)進(jìn)行，這種方法可以比較方便將原有的串行方法并行化，但是當(dāng) 數(shù)據(jù)規(guī)模比較大的時(shí)候，前期的數(shù)據(jù)分布式存儲(chǔ)以及后期的結(jié)果合并都比較耗時(shí);再一種 是基于多線程并行技術(shù)，各線程將任務(wù)劃分成多個(gè)子任務(wù)，每個(gè)子任務(wù)同時(shí)進(jìn)行，這種方法 可以充分利用本機(jī)計(jì)算資源，但是由于線程間并行屬于細(xì)粒度并行，各線程共享父進(jìn)程的 內(nèi)存空間，容易造成系統(tǒng)的不穩(wěn)定，并且這種方法可擴(kuò)展性不是很好；目前基于MPI (Message Passing Interface，消息傳遞接口）的多進(jìn)程方式進(jìn)行柵格影像的并行切片研 究較少，這種方法利用共享外存的高性能集群，可擴(kuò)展性強(qiáng)，穩(wěn)定性強(qiáng)，可以充分利用多機(jī) 計(jì)算資源，由于集群之間共享外存，所以數(shù)據(jù)無(wú)需提前分布式存儲(chǔ)，可以極大地提高影像切 片的效率。

【發(fā)明內(nèi)容】

[0005] 本發(fā)明的目的在于提供了一種基于新思路的柵格數(shù)據(jù)并行切片方法。本發(fā)明利用 多進(jìn)程技術(shù)對(duì)任務(wù)進(jìn)行劃分，每個(gè)進(jìn)程獨(dú)立進(jìn)行其相應(yīng)切片工作，當(dāng)某個(gè)進(jìn)程出現(xiàn)問(wèn)題時(shí)， 其余進(jìn)程仍可正常完成切片工作，互不影響。本發(fā)明切出的結(jié)果瓦片兼容目前絕大多數(shù)互 聯(lián)網(wǎng)地圖應(yīng)用，可以直接投入生產(chǎn)，并且本發(fā)明支持?jǐn)帱c(diǎn)，一旦發(fā)生斷電等現(xiàn)象，算法可在 原有工作的基礎(chǔ)上繼續(xù)進(jìn)行切片，避免了重復(fù)工作的問(wèn)題。
[0006] 針對(duì)上述技術(shù)問(wèn)題，本發(fā)明提供了一種并行方式柵格影像切片方法。步驟如下：
[0007] 第一步:獲取原始柵格影像，設(shè)置目標(biāo)瓦片的級(jí)別和進(jìn)程總數(shù)。
[0008] 設(shè)定目標(biāo)瓦片的級(jí)別level以及MPI參與的進(jìn)程總數(shù)n，以及瓦片輸出路徑。MPI會(huì) 為每一個(gè)進(jìn)程分配一個(gè)進(jìn)程號(hào)i (0 < i〈n，i取整數(shù)），后續(xù)可通過(guò)進(jìn)程號(hào)來(lái)達(dá)到對(duì)每個(gè)進(jìn)程 單獨(dú)控制的目的。
[0009] 第二步:讀取原始柵格影像的元數(shù)據(jù)信息。
[0010] 設(shè)置其中一個(gè)進(jìn)程為主進(jìn)程，主進(jìn)程讀取原始柵格影像的元數(shù)據(jù)信息，包括影像 的長(zhǎng)寬、波段數(shù)，投影坐標(biāo)系、數(shù)據(jù)的無(wú)效值；
[0011]第三步:將原始柵格影像變換至WebMercator投影。
[0012]當(dāng)前主流WebGIS軟件以及互聯(lián)網(wǎng)地圖應(yīng)用的切分瓦片坐標(biāo)系統(tǒng)普遍采用Google 提出的WebMercator投影坐標(biāo)系，因此在進(jìn)行切片前首先需將原始柵格影像的坐標(biāo)系統(tǒng)轉(zhuǎn) 換至WebMercator目標(biāo)投影系，然后基于WebMercator投影下的坐標(biāo)進(jìn)行切分。為了加快投 影變換的效率，本發(fā)明采用以下策略:主進(jìn)程根據(jù)之前讀取的原始影像投影坐標(biāo)系，如果原 始柵格影像不是WebMercator投影系(在地理柵格影像中，所有投影坐標(biāo)系都對(duì)應(yīng)一個(gè)由開(kāi) 放地理空間聯(lián)盟制定的WKT編碼，利用GDAL庫(kù)的相關(guān)函數(shù)從原始影像中讀取該影像的WTK形 式投影信息，通過(guò)比較原始影像的WKT編碼和WebMercator投影的WTK就可以判斷該影像的 投影系是否為WebMercator投影）。主進(jìn)程利用⑶AL類(lèi)庫(kù)的GDALAutoCreateWarpedVRT函數(shù)， 重采樣方法采用最鄰近內(nèi)插法，將原始柵格影像變換至WebMercator投影坐標(biāo)系下，得到投 影變換結(jié)果影像，以vrt格式文件進(jìn)行保存，后續(xù)的所有切片操作都是基于投影變換結(jié)果影 像進(jìn)行，并將原始柵格影像的無(wú)效值、影像寬度、波段數(shù)等元數(shù)據(jù)信息寫(xiě)入vrt文件中。vrt 是一個(gè)xml格式的文件，本身不存儲(chǔ)像素的色彩信息，它通過(guò)與原始柵格影像相關(guān)聯(lián)，通過(guò) 描述信息記錄自身像素的坐標(biāo)（自身像素坐標(biāo)就是指像素點(diǎn)在投影變換結(jié)果影像像素坐標(biāo) 系下的坐標(biāo))與相關(guān)聯(lián)柵格影像的映射關(guān)系，以及投影變換結(jié)果影像的長(zhǎng)寬、波段數(shù)等元數(shù) 據(jù)信息，當(dāng)通過(guò)vrt文件讀取柵格數(shù)據(jù)時(shí)，首先基于映射關(guān)系，找到其在原始柵格影像文件 中的位置，然后從原始柵格影像文件相應(yīng)位置中把柵格數(shù)據(jù)讀出。由于vrt文件只記錄一些 描述性信息，包括投影變換結(jié)果影像的長(zhǎng)寬、波段數(shù)、像素點(diǎn)的數(shù)據(jù)類(lèi)型、投影坐標(biāo)系，以及 與原始柵格影像的映射關(guān)系。所以vrt文件數(shù)據(jù)量小，因此vrt文件用于作為中間文件的投 影變換結(jié)果影像具有很大的優(yōu)勢(shì)，可以減少大量10操作。
[0013] 第四步:計(jì)算投影變換結(jié)果影像的瓦片最大級(jí)別和最小級(jí)別，并更新目標(biāo)瓦片的 級(jí)別。
[0014] 假設(shè)以像素為單位每個(gè)瓦片的大小為tilesizeXtilesize，投影變換結(jié)果影像長(zhǎng) 為Xsize、寬為Ysize，以米(m)為單位投影變換結(jié)果影像分辨率為resit。以一張世界地圖而 言，Resolution(h)表示其第h級(jí)別的瓦片分辨率，Resolution(h)可以通過(guò)以下公式計(jì)算得 出1^8〇1111:;[011(11) = (2\31\1?)/(1：;[168丨26\211)，1?為地球半徑等于6378137111，11為整數(shù)，31為 圓周率。投影變換結(jié)果影像的瓦片最大級(jí)別tmaxz為Resolution (h)最接近res It但不大于 resltB^^hfE； j|?T^minRes= (resit Xmax(Xsize ,Ysize))/tilesize 數(shù)，貝丨】瓦片最小級(jí)別tminz為1^801111:;[011(11)最接近111；[111^8但不大于111；[111^8時(shí)的11值。如果之 前設(shè)定的目標(biāo)瓦片級(jí)別level>tmaxz則將tmaxz的值賦給level，如果level〈tminz則將 tminz的值賦給level。而后在瓦片輸出目錄下創(chuàng)建名為level的文件夾。
[0015] 第五步：計(jì)算投影變換結(jié)果影像的瓦片行列號(hào)范圍，并按照路徑為/目標(biāo)瓦片級(jí) 另IJ/瓦片列號(hào)/瓦片行號(hào).png的路徑建立空瓦片文件；
[0016] 計(jì)算當(dāng)前l(fā)evel級(jí)別下投影變換結(jié)果影像地理范圍內(nèi)瓦片的行列號(hào)范圍，具體計(jì) 算方法如下:假設(shè)投影變換結(jié)果影像地理范圍為[ominX，ominY，omaxX，omaxY]，其中該范圍 為一個(gè)矩形四邊形，omi nX，omaxY為左上角點(diǎn)坐標(biāo)，omaxX，omi nY為右下角點(diǎn)坐標(biāo);瓦片行列 號(hào)范圍為[tminX，tminY，tmaxX，tmaxY]，其中該范圍為一個(gè)矩形四邊形，tminX，tmaxY為左 上角行列號(hào)，tmaxX，tminY為右下角行列號(hào)，則其之間滿足以下映射關(guān)系：
[0017] tminX = ceil( (ominX+originShift)/ (Resolution Xtilesize))-1；
[0018] tminY = cei1((ominY+originShift)/(Resolution Xtilesize))-1；
[0019] tmaxX = ceil( (omaxX+originShift)/ (Resolution Xtilesize))-1；
[0020] tmaxY = cei1((omaxY+