盧 昊，龐 勇，李增元

LU Hao,PANG Yong,LI Zengyuan

中國林業(yè)科學(xué)研究院 資源信息研究所，北京100091

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

1 引言

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

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

2 波形數(shù)據(jù)的常規(guī)訪問

根據(jù)LAS1.3 格式規(guī)范[10]可知，作為對較低版本格式的改進(jìn)，它將波形數(shù)據(jù)作為一個(gè)擴(kuò)展變長記錄（EXTENDED VARIABLE LENGTH RECORD）字段掛接在點(diǎn)云文件末尾或獨(dú)立同名文件中，一個(gè)包含了波形數(shù)據(jù)的LAS1.3 文件結(jié)構(gòu)如圖1 所示。

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

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

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

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

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

（3）對于包含波形數(shù)據(jù)的LAS 文件，點(diǎn)格式必須是FORMAT4 或FORMAT5，這時(shí)點(diǎn)記錄中增加了幾個(gè)字段用于關(guān)聯(lián)和表達(dá)波形信息，以FORMAT4 為例，其中與波形相關(guān)的屬性如表2 所示。

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

圖1 帶波形的LAS1.3 文件結(jié)構(gòu)

表2 LAS1.3 格式FORMAT4 點(diǎn)記錄字段

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

（1）讀入公共文件頭及相關(guān)參數(shù)。

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

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

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

3 倒排索引文件格式

3.1 LAS1.3 格式存儲(chǔ)波形數(shù)據(jù)的缺陷

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

3.2 構(gòu)建倒排索引

針對這些問題，本文設(shè)計(jì)了一種簡化訪問的索引文件格式，目的是為了建立從波形到點(diǎn)的一對多映射關(guān)系，同時(shí)簡化波形的訪問，便于進(jìn)行波形分解與系統(tǒng)生成點(diǎn)記錄的比對。

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

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

（2）創(chuàng)建索引數(shù)組，類型為8 字節(jié)無符號(hào)整型（保證足夠訪問大文件）。假定每個(gè)波形最多有5 次回波，每條波形對應(yīng)一個(gè)5 元數(shù)組，數(shù)組對應(yīng)位置存儲(chǔ)pts 文件中從文件開始到某個(gè)點(diǎn)記錄的偏移量。若該數(shù)組中某個(gè)元素（偏移量）為0 則代表波形沒有該次回波。

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

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

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

圖2 倒排索引文件結(jié)構(gòu)

3.3 基于倒排索引的波形數(shù)據(jù)訪問

經(jīng)過分割和創(chuàng)建索引處理，原始點(diǎn)云和波形信息全部進(jìn)入了新生成的文件，根據(jù)索引文件即可逐個(gè)訪問波形記錄并獲取從該脈沖獲取的點(diǎn)云數(shù)據(jù)。訪問方式如下所示：

（1）讀取idx 文件頭部無符號(hào)長整型4 個(gè)字節(jié)，為波形個(gè)數(shù)。

（2）讀取idx 文件中一個(gè)波形的索引記錄，每個(gè)記錄包含5 個(gè)8 字節(jié)無符號(hào)整型數(shù)據(jù)，共40 個(gè)字節(jié)。

（3）根據(jù)該索引記錄的位置，計(jì)算對應(yīng)波形在wdp文件中的位置，讀取該波形數(shù)據(jù)包，長度為128 字節(jié)，每個(gè)采樣1 字節(jié)（不同激光雷達(dá)設(shè)備量化方式可能不同）。

（4）遍歷該波形的5 個(gè)索引記錄數(shù)值，直到最后一個(gè)不為0 的位置，根據(jù)該文件偏移量讀取pts 文件中的對應(yīng)點(diǎn)記錄，正常情況下數(shù)量為1～5。

（5）重復(fù)（2）直到結(jié)束。

4 實(shí)驗(yàn)結(jié)果

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

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

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

圖3 倒排索引文件訪問波形及激光點(diǎn)云數(shù)據(jù)

表4 波形數(shù)據(jù)遍歷速度量化實(shí)驗(yàn)結(jié)果

圖3 中，黑色曲線為一束激光脈沖的回波波形，紅色點(diǎn)為原始點(diǎn)云文件中由該脈沖回波獲取的激光點(diǎn)云位置及振幅，通過直接獲取點(diǎn)云數(shù)據(jù)，可以快速進(jìn)行波形與系統(tǒng)點(diǎn)云結(jié)合的算法處理[15]?？梢钥闯?，該束激光脈沖獲取了2 次反射回波，其位置由紅色點(diǎn)表示。這兩次回波信息可以從倒排文件中通過點(diǎn)記錄直接獲取。而在原始LAS1.3 文件中，無法得知這兩個(gè)點(diǎn)記錄與該波形記錄的關(guān)系，只能在訪問一個(gè)點(diǎn)記錄時(shí)訪問到該波形，故這種情況下同一波形記錄將被重復(fù)訪問2 次。

5 結(jié)論

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

[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] 李奇，馬洪超.基于激光雷達(dá)波形數(shù)據(jù)的點(diǎn)云生產(chǎn)[J].測繪學(xué)報(bào)，2008，37（3）：349-354.

[3] 張靖，高偉.LAS 格式解析及其擴(kuò)展域的應(yīng)用[J].測繪科學(xué)，2008，33（3）：154-155.

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

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

[6] 王智強(qiáng)，劉建毅.一種實(shí)時(shí)更新索引結(jié)構(gòu)的設(shè)計(jì)與實(shí)現(xiàn)[J].計(jì)算機(jī)系統(tǒng)與應(yīng)用，2005，41（10）：79-82.

[7] 鄧攀，劉功申.一種高效的倒排索引存儲(chǔ)結(jié)構(gòu)[J].計(jì)算機(jī)工程與應(yīng)用，2008，44（31）：149-152.

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

[9] 董長春.基于Hadoop 的倒排索引技術(shù)的研究[D]沈陽：遼寧大學(xué)，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.