倒排索引優化的波形激光雷達數據存儲和訪問

盧 昊，龐 勇，李增元

LU Hao,PANG Yong,LI Zengyuan

中國林業科學研究院 資源信息研究所，北京100091

Institute of Forest Resources Information Techniques,Chinese Academy of Forestry,Beijing 100091,China

1 引言

激光雷達（LiDAR）技術是目前國際上發展迅速的一種主動遙感手段，具有直接快速獲取地表三維信息的特點[1]。全波形激光雷達通過記錄激光脈沖的完整回波，提供了更豐富的地物信息[2]。為了實現不同軟硬件平臺之間的數據交換，美國攝影測量與遙感協會設計發布了LAS 格式規范作為通用的激光數據存儲格式，其中LAS1.0 至LAS1.2 僅僅支持點云信息存儲，LAS1.3 開始增加了對波形數據的支持。Andre Samberg 在2007 年ISPRS 會議上解析并比較了多種版本的LAS 文件格式，張婧等（2008）進行了LAS 格式的解析并分析了其擴展域的應用[3]，趙自明（2010）等提出了采用內存映射技術對點云數據進行快速讀取和顯示的方法[4]，秦楠楠等（2013）對LAS1.3 文件格式進行了解析并提出了讀取波形數據的方法，總結了快速讀取波形數據的高效方法[5]。但LAS1.3 格式的主要關注點仍然是點數據的表達，波形信息只是作為一種額外的擴展字段附加在原始的點文件中。要訪問波形數據，只能通過從點到波形的映射關系。對于全波形激光雷達來說，波形數據包含了激光脈沖回波的全部信息，回波點僅僅是其中的一個信息子集，從信息獲取的角度講，這種邏輯關系與原始信息流相反，并且帶來了一些其他問題（見3.1 節）。波形數據的處理主要有波形分解和反卷積等，要求直接訪問波形采樣序列，而LAS1.3 格式并不支持此種訪問方式。

倒排索引（Inverted index），也常被稱為倒排索引、置入檔案或倒排檔案[6]，是一種常用于關系型數據庫和搜索引擎檢索算法的索引方法，可以用來存儲在全文搜索下某個單詞在一個文檔或者一組文檔中的存儲位置的映射。它是文檔檢索系統中最常用的數據結構，源于實際應用中需要根據屬性值來查找記錄。倒排序索引表中的每一項都包括一個屬性值和具有該屬性值的各記錄的地址。在這種索引方式下，不是由記錄來確定屬性值，而是由屬性值來確定記錄的位置。王智強等（2005）利用倒排索引實現了一種實時更新的索引結構[6]，用于實現搜索引擎對網頁數據庫的快速檢索和實時更新。鄧攀等（2008）提出一種高效的倒排索引結構[7]，設計了一種優化的倒排索引結構，充分考慮索引表的壓縮和對磁盤IO 的降低以支持大規模的網絡信息檢索。鄭榕增等（2010）利用Lucene 實現了中文倒排索引技術[8]，基于開源的全文檢索引擎工具包Lucene 實現了對中文文檔的高效分詞和索引，較好地改進了搜索引擎對中文字符編碼的支持。董長春（2011）基于Hadoop 進行了倒排索引技術的研究[9]，在Hadoop 集群上實現了多層次倒排索引結構，在減少了集群節點間通信代價的基礎上進行大量信息的查詢。可以看出，倒排索引是互聯網與信息系統中的一種感覺常用技術，主要被應用于網絡搜索引擎與文檔處理系統，而尚未被應用于海量遙感數據處理。事實上，在全波形激光雷達大量數據的研究與后處理過程中，對于LAS1.3 格式中存在的只能從點云通過多對一映射訪問波形而不能進行快速反向檢索（詳見3.1）的問題，利用倒排索引可以很好地解決。

2 波形數據的常規訪問

根據LAS1.3 格式規范[10]可知，作為對較低版本格式的改進，它將波形數據作為一個擴展變長記錄（EXTENDED VARIABLE LENGTH RECORD）字段掛接在點云文件末尾或獨立同名文件中，一個包含了波形數據的LAS1.3 文件結構如圖1 所示。

從圖1 以及文獻[10]看出，在舊版本的基礎上，LAS1.3 格式有了如下改變：

（1）文件頭部分結構基本不變，在文件頭末尾增加了波形數據包記錄起始位置（Start of Waveform Data Packet Record）字段，表示從文件開始到波形擴展變長記錄區的偏移量。

（2）變長記錄區中有一個變長記錄為波形包描述符用戶定義記錄，其擴展域中有1～255 個描述符，其格式如表1 所示，每一個記錄表示一種波形包量化方式。

表1 波形包描述符字段格式

也就是說，一個LAS1.3 文件中的波形，從量化比特數、采樣數和采樣間隔等信息來看，最多只能有255 種不同情況，超過255 種則無法完全表示。實際應用中，由于設備在使用時并沒有修改其記錄方式，所以一般只有1 個記錄，這樣對于數據讀寫比較方便。

（3）對于包含波形數據的LAS 文件，點格式必須是FORMAT4 或FORMAT5，這時點記錄中增加了幾個字段用于關聯和表達波形信息，以FORMAT4 為例，其中與波形相關的屬性如表2 所示。

其中，Wave Packet Descriptor Index 為0 代表該點不存在波形，1～255 則表示該點存在波形，且描述該波形的量化參數保存在前面變長記錄字段的1～255 里對應的位置。Byte offset to waveform data 表示該點所屬的激光脈沖回波波形數據從evlr 字段開始到其第一個采樣的字節偏移量。

圖1 帶波形的LAS1.3 文件結構

表2 LAS1.3 格式FORMAT4 點記錄字段

綜上所述，LAS1.3 文件中波形的訪問方式可以概況為：

（1）讀入公共文件頭及相關參數。

（2）遍歷可變長度記錄：逐一判斷每個可變長度記錄頭部，若為波形包描述符用戶定義記錄，則獲取波形包描述符用戶定義記錄個數，并全部讀入數組。每個記錄的長度固定為26 字節。

（3）遍歷讀取每個點記錄，判斷是否存在波形。如存在，則首先根據波形包描述符索引（Wave Packet Descriptor Index）獲取描述符數組中對應元素，然后根據描述符中信息訪問波形數據字段。此時，內存中存在（a）該點數據；（b）波形描述數據；（c）波形采樣數據，根據這些信息可以進行對該波形的算法處理。

步驟（3）中之所以會出現不存在波形的情況，是由于某些激光雷達系統受到硬件存儲速度的限制，工作時會以一定的采樣比例對點云進行波形記錄，導致部分點云并不關聯其對應激光脈沖的波形。如Leica ALS 60系統對脈沖回波以1∶1 的比例對點云進行波形采集和丟棄，使文件中一半的點云不存在相應脈沖回波記錄，這一現象源于硬件技術的限制，不利于全波形數據的研究與應用。而在具備真正全波形記錄能力的系統中，脈沖回波將被全部記錄。如Riegl LMS Q680 系統[11]可以對脈沖的全部有效回波信號進行記錄，其后處理軟件根據波形信息進行點云解算，故該系統獲取的點云數據中全部的點都會關聯其對應激光脈沖信號而無一遺漏。LAS1.3 格式通過這種可選的存儲方式可以兼容不同的硬件設備所采集的點云和波形數據。

3 倒排索引文件格式

3.1 LAS1.3 格式存儲波形數據的缺陷

按照LAS1.3 格式定義，波形數據的存儲和訪問存在的問題包括：一、允許存在的波形量化參數組合是有限的，所以不能存儲任意個數任意長度的波形。二、波形往往不是完整的，只是對部分有效信號進行記錄。三、存儲的方式是從點到波形的多對一映射，故只能從點到波形訪問，多對一的情況下存在重復訪問，在波形分解等算法中這是不被允許的。四、訪問一個波形時，由于原文件不存在相應的映射關系，因此并不知道哪些點記錄由該波形得到。一般而言，雖然同屬于一個波形的點記錄前后相鄰，但無法保證點數據讀取的準確和無遺漏。在常規的訪問方式下，如果需要隨機訪問指定回波波形數據關聯的點云數據，會引入大量的查找計算和磁盤IO，顯著降低索引訪問的速度。

3.2 構建倒排索引

針對這些問題，本文設計了一種簡化訪問的索引文件格式，目的是為了建立從波形到點的一對多映射關系，同時簡化波形的訪問，便于進行波形分解與系統生成點記錄的比對。

首先將原始LAS1.3 文件進行分割，將除最后一個波形擴展變長記錄（evlr）字段之外的其余部分視為點文件（*.pts），將原始文件中波形包描述符用戶定義記錄與末尾的波形擴展變長記錄（evlr）字段放在一起，構成波形數據包文件（*.wdp），然后創建索引文件，實現從波形到點的映射，具體操作是：

（1）計算出波形條數N，原始文件中并沒有給出，需要根據總波形長度與波形字節長度求出。

（2）創建索引數組，類型為8 字節無符號整型（保證足夠訪問大文件）。假定每個波形最多有5 次回波，每條波形對應一個5 元數組，數組對應位置存儲pts 文件中從文件開始到某個點記錄的偏移量。若該數組中某個元素（偏移量）為0 則代表波形沒有該次回波。

（3）在原始LAS 文件中讀取一個點記錄，如果某點無波形則跳過；如果有波形，根據點記錄中相應字段計算出波形包在evlr 中的位置(0，1，2，…，N-1)，該索引為idx。索引數組中下標為idx的元素即為該點對應的波形，若idx＜0 或idx＞N-1，屬于數據異常，實際中存在數據錯誤的點記錄會導致該情況。在這5 元數組中，將該點在pts 文件中從文件開始計算的偏移量填入第一個不為0 的位置。若5 個位置都已存在數據，說明原始文件中存在5 次以上回波的波形，屬于異常。

（4）重復（3）步驟直到遍歷完所有的點，得到最終的索引數組，寫入索引文件（*.idx）。

三種文件的格式及訪問邏輯關系如圖2 所示。

圖2 倒排索引文件結構

3.3 基于倒排索引的波形數據訪問

經過分割和創建索引處理，原始點云和波形信息全部進入了新生成的文件，根據索引文件即可逐個訪問波形記錄并獲取從該脈沖獲取的點云數據。訪問方式如下所示：

（1）讀取idx 文件頭部無符號長整型4 個字節，為波形個數。

（2）讀取idx 文件中一個波形的索引記錄，每個記錄包含5 個8 字節無符號整型數據，共40 個字節。

（3）根據該索引記錄的位置，計算對應波形在wdp文件中的位置，讀取該波形數據包，長度為128 字節，每個采樣1 字節（不同激光雷達設備量化方式可能不同）。

（4）遍歷該波形的5 個索引記錄數值，直到最后一個不為0 的位置，根據該文件偏移量讀取pts 文件中的對應點記錄，正常情況下數量為1～5。

（5）重復（2）直到結束。

4 實驗結果

本文在Visual Studio 2010 環境下以C++語言實現了LAS1.3 文件的數據索引構建與波形數據訪問，以Leica ALS 60 機載激光雷達[12]獲取的帶有波形數據的LAS1.3 文件為實驗對象[13]。從實驗結果看，基于倒排索引的格式與LAS1.3 格式在波形數據訪問特性對比如表3 所示。其中，前三項反映出倒排索引文件可以完整保留原始LAS1.3 文件的屬性而不造成信息丟失。由于倒排索引文件保留了原始點云數據的序列存儲特征，所以仍然可以進行順序訪問，而由于點云信息的定長特征，所以能夠通過序號計算偏移量定位任意點記錄以實現隨機訪問。在LAS1.3 文件中，由于波形的起始位置只能通過某一點記錄的屬性進行索引定位，所以無法對波形直接進行遍歷訪問，不能支持順序或隨機訪問。而在倒排索引文件中，全部波形的位置信息都以全局偏移量的形式保存在倒排索引文件中，且倒排索引文件采用了定長字段的設計，故同時支持波形數據的順序和隨機訪問。倒排索引文件同時還保留了LAS1.3 文件內部的正排索引信息，故仍可以隨機訪問指定點記錄的波形。由于采用了倒排索引結構，由指定激光脈沖所獲取的點云記錄，可以通過波形索引直接定位，這一特性在原始LAS1.3 中是無法實現的，若進行全局搜索則會帶來巨大的磁盤IO 和計算代價。此外，對于點到波形的多對一映射，倒排索引有效地消除了波形重復訪問的弊端，而隨機訪問特性保證了倒排索引較快的遍歷速度。

從表3 可以看出，本文實現的索引文件組不僅保留了LAS1.3 格式包含的所有屬性和結構特征，還增加了一些新的特性，有利于對波形數據和點云數據結合展開處理算法[14]。讀取索引文件組的波形數據并獲取關聯的點云實例顯示如圖3 所示。波形數據遍歷速度根據應用程序循環訪問波形記錄的磁盤IO 速率定量表征，一般用record/s 或MB/s 表示。表4 顯示了針對一組測試數據進行的兩種訪問速度統計。從表中可以看出，對于不同測量環境下生成的激光雷達數據，通過倒排索引方式均能一定程度地加快波形數據尋址訪問速度。從平均速率可以看出，常規訪問方式下的速率表現出文件相關的特征，這是由于不同探測目標對于激光脈沖回波次數分布有差別，導致多脈沖回波造成的冗余訪問各不相同所致。而在倒排索引訪問方式下由于冗余訪問得到了消除，所以總體速率幾乎不受不同文件的影響趨于接近，其尋址速度主要由磁盤IO 性能決定。

表3 兩種文件格式的特征比較

圖3 倒排索引文件訪問波形及激光點云數據

表4 波形數據遍歷速度量化實驗結果

圖3 中，黑色曲線為一束激光脈沖的回波波形，紅色點為原始點云文件中由該脈沖回波獲取的激光點云位置及振幅，通過直接獲取點云數據，可以快速進行波形與系統點云結合的算法處理[15]。可以看出，該束激光脈沖獲取了2 次反射回波，其位置由紅色點表示。這兩次回波信息可以從倒排文件中通過點記錄直接獲取。而在原始LAS1.3 文件中，無法得知這兩個點記錄與該波形記錄的關系，只能在訪問一個點記錄時訪問到該波形，故這種情況下同一波形記錄將被重復訪問2 次。

5 結論

本文通過對LAS1.3 格式的波形和點云數據進行倒排索引存儲，實現了波形數據的直接快速訪問，對波形數據處理是一種比較理想的轉換格式。索引文件組既可以單獨作為點云數據載體，也可以單獨作為波形數據載體，或者作為波形與點云數據結合處理算法的載體，保持了激光脈沖回波與相應回波點云的邏輯關系，具有很好的靈活性，一定程度上消除了LAS1.3 格式規范中對波形數據支持的缺陷。但由于LAS1.3 格式并沒有嚴格限定波形采樣量化的方式，波形包的結構與實際回波情況相關，難以用一種固定的形式進行描述。索引文件結構的設計也無法很好地解決波形長度、量化比特數不一致的問題，這將成為后續研究的重點。但總體而言，本文設計的倒排索引結構體現了波形數據與激光點云的內在邏輯關系，可以作為一種激光雷達波形數據通用交換格式的設計思想。

[1] Wagner W.Gaussian decomposition and calibration of a novel small-footprint full-waveform digitising airborne laser scanner[J].ISPRS Journal of Photogrammetry & Remote Sensing，2006，60：100-112.

[2] 李奇，馬洪超.基于激光雷達波形數據的點云生產[J].測繪學報，2008，37（3）：349-354.

[3] 張靖，高偉.LAS 格式解析及其擴展域的應用[J].測繪科學，2008，33（3）：154-155.

[4] 趙自明，史兵，田喜平，等.LAS 格式解析及其數據的讀取與顯示[J].測繪技術裝備，2010：17-20.

[5] 秦楠楠，賴旭東.機載LiDAR波形數據快速讀取方法研究[C]//第一屆全國高分辨率遙感數據處理與應用研討會，西安，2011：236-240.

[6] 王智強，劉建毅.一種實時更新索引結構的設計與實現[J].計算機系統與應用，2005，41（10）：79-82.

[7] 鄧攀，劉功申.一種高效的倒排索引存儲結構[J].計算機工程與應用，2008，44（31）：149-152.

[8] 鄭榕增，林世平.基于Lucene的中文倒排索引技術的研究[J].計算機技術與發展，2010，20（3）：80-83.

[9] 董長春.基于Hadoop 的倒排索引技術的研究[D]沈陽：遼寧大學，2011.

[10] ASPRS.Las Specification Version 1.3 Spe[S].USA：ASPRS，2010.

[11] IGI GmbH.LiteMapper 6800 User Manual[S].Germany：IGI GmbH，2010.

[12] Leica Geosystems.Leica ALS60 Airborne Laser Scanner Product Specifications[S].Switzerland：Leica Geosystems，2008

[13] Lu H.Quantitative comparison between full waveform decomposition and discrete returns point clouds data from Small-footprint Airborne LiDAR[C]//13th International Conference on LiDAR Applications for Assessing Forest Ecosystems（SilviLaser2013），Beijing，China，2013：63-70.

[14] Pang Y.LiCHY：CAF’s LiDAR，CCD and Hyperspectral Airborne Imager[C]//13th International Conference on LiDAR Applications for Assessing Forest Ecosystems（SilviLaser2013），Beijing，China，2013：45-54.

[15] Isenburg M.PulseWaves-full waveform LiDAR specification（Version 0.3 Revision 8）[Z].Germany：rapidlasso GmbH，2011.