宋金秀 譚弘武 張洪雷 王 成 楊大兵

(北京煜邦電力技術(shù)股份有限公司, 北京100089)

0 引言

目前,激光點(diǎn)云數(shù)據(jù)處理包括點(diǎn)云濾波、點(diǎn)云配準(zhǔn)、點(diǎn)云分割與分類、點(diǎn)云三維重建等相關(guān)技術(shù),而在點(diǎn)云三維瓦片處理方面的應(yīng)用較少,且瓦片數(shù)據(jù)處理多為地圖切片。

地圖切片Web服務(wù)簡稱Web地圖瓦片服務(wù)(Web Map Tile Service,WMTS)，是開放地理空間信息聯(lián)盟(Open Geospatial Consortium,OGC)提出的緩存技術(shù)標(biāo)準(zhǔn),即在服務(wù)器端緩存被切割成一定大小瓦片的地圖,對客戶端只提供這些預(yù)先定義好的單個(gè)瓦片的服務(wù),將更多的數(shù)據(jù)處理操作如圖層疊加等放在客戶端,從而緩解地理信息系統(tǒng)(Geographic Information System,GIS)服務(wù)器端數(shù)據(jù)處理的壓力,改善用戶體驗(yàn)[1-3]。WMTS使用瓦片矩陣集來表示切割后的地圖。瓦片就是包含地理數(shù)據(jù)的矩形影像,一幅地圖按一定的瓦片大小被切割成多個(gè)瓦片,形成瓦片矩陣,一個(gè)或多個(gè)瓦片矩陣即組成瓦片矩陣集[4]。不同的瓦片矩陣具有不同的分辨率,每個(gè)瓦片矩陣由瓦片矩陣標(biāo)識(shí)符進(jìn)行標(biāo)識(shí)。

1 地圖切片概況

隨著科技的不斷發(fā)展,用戶對瓦片服務(wù)要求的提高,瓦片數(shù)據(jù)處理技術(shù)逐漸從地圖切片向點(diǎn)云切片、模型切片等方向不斷發(fā)展。而目前的點(diǎn)云切片處理多數(shù)針對局部數(shù)據(jù)處理,在數(shù)據(jù)追加和局部更新方面存在一定的不足,在處理海量數(shù)據(jù)時(shí)存在操作繁復(fù)、數(shù)據(jù)處理效率較低等問題[5]。而且針對某些特殊需求,如對已分類的點(diǎn)云數(shù)據(jù)進(jìn)行特定類型切片處理，需要將特定類型的數(shù)據(jù)提取出來再進(jìn)行切片處理,缺乏制定化和多功能的集成[6-7]。

因此,如何提供一種瓦片數(shù)據(jù)處理效率高、支持?jǐn)?shù)據(jù)追加和局部更新的基于全球索引的點(diǎn)云切片處理方法，是本領(lǐng)域技術(shù)人員亟須解決的問題[8-10]。

2 全球索引點(diǎn)云切片發(fā)布設(shè)計(jì)

基于全球索引的點(diǎn)云切片發(fā)布處理包括：獲取欲處理的原始點(diǎn)云數(shù)據(jù);將所有原始點(diǎn)云數(shù)據(jù)的縮放級(jí)別進(jìn)行統(tǒng)一;按照統(tǒng)一后的縮放級(jí)別對原始點(diǎn)云數(shù)據(jù)進(jìn)行四叉樹分塊或八叉樹分塊,生成多個(gè)數(shù)據(jù)塊,并對每個(gè)數(shù)據(jù)塊進(jìn)行全球索引,獲得與每個(gè)數(shù)據(jù)塊一一對應(yīng)的索引值;分別按照每個(gè)數(shù)據(jù)塊索引值遍歷原有點(diǎn)云切片成果數(shù)據(jù)庫中所有三維瓦片,并按照預(yù)設(shè)條件對原有點(diǎn)云切片成果數(shù)據(jù)庫進(jìn)行更新?；谌蛩饕衅l(fā)布處理可對激光點(diǎn)云數(shù)據(jù)進(jìn)行全球索引瓦片劃分,且支持?jǐn)?shù)據(jù)追加和局部更新,并可針對不同類型的點(diǎn)云篩選進(jìn)行瓦片劃分。

全球索引的點(diǎn)云切片處理方法包括全球索引瓦片劃分和數(shù)據(jù)追加及更新,其中全球索引瓦片劃分方法包括：

(1)確定最大級(jí)別：統(tǒng)一級(jí)別,以方便融合處理。

(2)文件切割：文件太大,需要分割處理,避免內(nèi)存溢出的問題。

(3)分塊規(guī)則：前10級(jí)按全球經(jīng)緯度進(jìn)行分塊；大于10級(jí)時(shí),在10級(jí)瓦片的基礎(chǔ)上按八叉樹分塊。

(4)緩存和融合：大文件切割后的小文件或者是和其他文件做增量處理時(shí)都需要進(jìn)行融合。

(5)保存點(diǎn)云瓦片數(shù)據(jù)：數(shù)據(jù)追加及更新時(shí)需首先判斷當(dāng)前數(shù)據(jù)所在的瓦片,遍歷成果數(shù)據(jù)庫中所有瓦片并在緩存中讀取與當(dāng)前瓦片索引相同的瓦片；其次遍歷所有點(diǎn),根據(jù)點(diǎn)的網(wǎng)格索引從緩存的瓦片中查找有沒有相同網(wǎng)格的點(diǎn),有的話舍棄緩存中的點(diǎn),緩存更新后若無相同的點(diǎn)，則將新讀入的點(diǎn)全部追加進(jìn)入緩存；最后將新的瓦片保存到成果數(shù)據(jù)庫中。

3 基于全球索引的點(diǎn)云切片發(fā)布實(shí)現(xiàn)

3.1 全球索引點(diǎn)云切片發(fā)布

基于全球索引切片發(fā)布處理流程如圖1所示。

圖1 點(diǎn)云切片流程圖

(1)預(yù)先進(jìn)行參數(shù)設(shè)置(包括空間坐標(biāo)參數(shù)、不切片類型和輸出類型),并輸入原始點(diǎn)云文件,即獲取欲處理的原始點(diǎn)云數(shù)據(jù)。

(2)統(tǒng)一級(jí)別：將所有所述原始點(diǎn)云數(shù)據(jù)的縮放級(jí)別進(jìn)行統(tǒng)一。

(3)文件切割及全球索引：按照統(tǒng)一后的縮放級(jí)別對所述原始點(diǎn)云數(shù)據(jù)進(jìn)行四叉樹分塊或八叉樹分塊,生成多個(gè)數(shù)據(jù)塊,并對每個(gè)數(shù)據(jù)塊進(jìn)行全球索引,獲得與每個(gè)數(shù)據(jù)塊一一對應(yīng)的索引值。

(4)數(shù)據(jù)追加及更新：分別按照每個(gè)數(shù)據(jù)塊索引值遍歷原有點(diǎn)云切片成果數(shù)據(jù)庫中所有三維瓦片,并按照預(yù)設(shè)條件對原有點(diǎn)云切片成果數(shù)據(jù)庫進(jìn)行更新。

3.2 全球索引瓦片劃分

全球索引瓦片劃分主要包括確定最大級(jí)別和瓦片劃分兩步,具體流程如圖2所示。

圖2 全球索引流程圖

3.2.1統(tǒng)一級(jí)別

統(tǒng)一級(jí)別,以方便融合處理,若未指定則進(jìn)行如下處理：

(1)判斷是否預(yù)設(shè)有最大縮放級(jí)別i,若是,則將最大縮放級(jí)別i作為統(tǒng)一縮放級(jí)別;若否,則執(zhí)行步驟2;

(2)取所述原始點(diǎn)云數(shù)據(jù)中間位置的點(diǎn)進(jìn)行估算;

(3)計(jì)算所取點(diǎn)的包圍盒以及包圍盒的中心點(diǎn)坐標(biāo);

(4)計(jì)算包圍盒所有點(diǎn)到中心點(diǎn)的平均距離L;若平均距離L大于預(yù)設(shè)的最大縮放級(jí)別i下的點(diǎn)云三維瓦片分辨率,則確定最大縮放級(jí)別為i;否則執(zhí)行步驟5;

(5)i自增1;

(6)重復(fù)執(zhí)行步驟4～5,直至平均距離L不大于預(yù)設(shè)的最大縮放級(jí)別i下的點(diǎn)云三維瓦片分辨率。

3.2.2全球瓦片劃分

(1)全球瓦片劃分規(guī)則為：0級(jí)時(shí),1個(gè)瓦片的分辨率為赤道周長;1級(jí)時(shí),1個(gè)瓦片的分辨率為赤道周長的1/2;2級(jí)時(shí),1個(gè)瓦片的分辨率為赤道周長的1/4;N級(jí)時(shí),1個(gè)瓦片的分辨率為赤道周長的1/2N；

(2)前10級(jí)按全球經(jīng)緯度以四叉樹分塊進(jìn)行處理,即0級(jí)時(shí),全球兩塊瓦片分辨率為180度;1級(jí)時(shí),全球4×2塊瓦片分辨率為90°;N級(jí)時(shí)全球有2×2N個(gè)瓦片,1個(gè)瓦片的分辨率為180°的1/2N。當(dāng)切片等級(jí)大于10級(jí)時(shí),在10級(jí)瓦片的基礎(chǔ)上按八叉樹分塊,八叉樹分塊是可看作是三維下的四叉樹分塊規(guī)則,即在平面xy和高程z上同時(shí)進(jìn)行四叉樹分塊。

3.3 數(shù)據(jù)追加及更新

數(shù)據(jù)追加及更新流程如圖3所示。

圖3 數(shù)據(jù)追加及更新流程圖

(1)確定當(dāng)前數(shù)據(jù)塊所對應(yīng)的索引值;

(2)按照當(dāng)前數(shù)據(jù)塊的索引值遍歷原有點(diǎn)云切片成果數(shù)據(jù)庫中所有三維瓦片,判斷原有點(diǎn)云切片成果數(shù)據(jù)庫中是否存在與當(dāng)前數(shù)據(jù)塊索引值相同的原有三維瓦片;如果不存在,則將當(dāng)前數(shù)據(jù)塊新增至原有點(diǎn)云切片成果數(shù)據(jù)庫;如果存在,則執(zhí)行下一步驟;

(3)將當(dāng)前數(shù)據(jù)塊劃分為網(wǎng)格,并確定每個(gè)網(wǎng)格點(diǎn)的索引值;

(4)按照當(dāng)前數(shù)據(jù)塊網(wǎng)格點(diǎn)索引值遍歷與當(dāng)前數(shù)據(jù)塊索引值相同的原有三維瓦片的所有網(wǎng)格點(diǎn),查找是否存在與當(dāng)前數(shù)據(jù)塊網(wǎng)格點(diǎn)索引值相同的原有網(wǎng)格點(diǎn);如果存在,則將原有網(wǎng)格點(diǎn)替換為當(dāng)前數(shù)據(jù)塊網(wǎng)格點(diǎn);如果不存在,則將當(dāng)前數(shù)據(jù)塊網(wǎng)格點(diǎn)追加至原有點(diǎn)云切片成果數(shù)據(jù)庫。

3.4 點(diǎn)云數(shù)據(jù)切片程序?qū)崿F(xiàn)

后臺(tái)算法處理采用Visual C++開發(fā)語言,基于地理空間數(shù)據(jù)抽象庫(Geospatial Data Abstraction Library,GDAL)和點(diǎn)云庫(Point Cloud Library,PCL)開源庫,實(shí)現(xiàn)基于全球索引的點(diǎn)云切片發(fā)布處理算法的代碼工作,前端展示采用Direct3D進(jìn)行三維數(shù)據(jù)加載渲染。C++點(diǎn)云切片發(fā)布處理部分代碼如下所示：

∥獲取點(diǎn)云數(shù)據(jù)所在瓦片索引號(hào)

int Xcount=int(GetData->m_PointHeader.X/Gridsize)+1;

int Ycount=int(GetData->m_PointHeader.Y/Gridsize)+1;

for(int i=0;i

{

for (int j=0;j

{

∥計(jì)算得到當(dāng)前瓦片數(shù)據(jù)范圍

std::vector GridPoints;

minx = minxm_PointHeader.MinX?GetData->m_PointHeader.MinX:minx;

miny= minym_PointHeader.MinY?GetData->m_PointHeader.MinY:miny;

maxx = maxx>GetData->m_PointHeader.MaxX?GetData->m_PointHeader.MaxX:maxx;

maxy = maxy>GetData->m_PointHeader.MaxY?GetData->m_PointHeader.MaxY:maxy;

∥構(gòu)建新的瓦片數(shù)據(jù)

C_DataIndex *m_DataIndex=new C_DataIndex();

m_DataIndex->XY0[0]=minx-0.1;

m_DataIndex->XY0[1]=miny-0.1;

m_DataIndex->XY1[0]=maxx+0.1;

m_DataIndex->XY1[1]=maxy+0.1;

m_DataIndex->GridHigh=maxy-miny;

m_DataIndex->GridWide=maxx-minx;

∥瓦片數(shù)據(jù)重新構(gòu)建索引

m_DataIndex->SetGridIndex(GridPoints);

for (size_t k=0;kBlockInput.size();++k)

{

if (m_DataIndex->BlockInput[k].size()==0) continue;

int _utm =get_utm(Prj);

∥點(diǎn)云數(shù)據(jù)平面坐標(biāo)轉(zhuǎn)換為經(jīng)緯度數(shù)據(jù)UTMxy_to_latlon_deg(m_DataIndex->GridCoordinate[k].X-m_Gridsize/2,m_DataIndex->GridCoordinate[k].Y-m_Gridsize/2,_utm, false,m_miny,m_minx);UTMxy_to_latlon_deg(m_DataIndex->GridCoordinate[k].X+m_Gridsize/2,m_DataIndex->GridCoordinate[k].Y+m_Gridsize/2,_utm, false,m_maxy,m_maxx);

std::vector vecUtmXY;

vecUtmXY.resize(m_DataIndex->BlockInput[k].size());

∥數(shù)據(jù)寫入

if (!IsWrite)

{

OutLas->m_PointHeader.DataFormatId='2';

bool b_write = true;

int isizeadd = OutGridPoints.size()+Points200.size();

if(isizeadd<200&&isizeaddm_PointHeader.PointRecordsCount)

b_write = false;

∥新的瓦片數(shù)據(jù)

if(b_write)

{

for (int ii=0;ii

OutGridPoints.push_back(Points200[i]);

Points200.clear();

std::vector().swap(Points200);

OutLas->Points=OutGridPoints;

OutLas->WriteLas(LASOutPathName);

IsWrite=true;

}

∥瓦片數(shù)據(jù)局部更新

else

{

int isi = Points200.size()+OutGridPoints.size();

for (int im=0;im

Points200.push_back(OutGridPoints[im]);

}

}

else SplitCloud::AppendLas(LASOutPathName,OutLas->m_PointHeader,OutGridPoints);∥數(shù)據(jù)更新追加

OutGridPoints.clear();

std::vector().swap(OutGridPoints);

}

}

}

經(jīng)過上述方法處理后的切片點(diǎn)云數(shù)據(jù)加載到Direct3D中,數(shù)據(jù)渲染支持海量數(shù)據(jù)無縫三維瀏覽顯示。

4 結(jié)束語

經(jīng)由上述的技術(shù)方案可知,與現(xiàn)有技術(shù)相比,本文公開提供了一種基于全球索引的點(diǎn)云切片處理方法,利用全球索引瓦片劃分的方法,解決了在處理海量點(diǎn)云數(shù)據(jù)的過程中效率較低的問題,并且可以提供數(shù)據(jù)追加和更新等方式,為后續(xù)的數(shù)據(jù)更新和發(fā)布提供極大的便利。同時(shí),可對分類后的點(diǎn)云數(shù)據(jù)進(jìn)行特定類型數(shù)據(jù)的瓦片處理,在多功能集成上有較好的適用性。