如何使用liblas获取las点的“数据”字段中存储的内容？

Question

我正在工作的多脉冲激光雷达数据，收集点沿多条线在飞行路径。我正试图确定las文件中各个行的名称和数量。我在python中使用liblas模块。

我找到了这份文件，它解释了存储在las文件中的不同字段。它提到页面底部的数据字段(get_data和set_data)。

标题中的“点数据格式”和“点数据记录长度”为这个“数据”字段留出了空间。我的标题说我为数据字段预留了28个字节，数据字段中存储了28个值。第19个值(至少在来自两个不同传感器的两个数据集中)是指行号。我在单脉冲数据中有一个值，在多脉冲数据中有4个值.

我想知道是否有一个标准存储在这个领域，或者它是否是专有的。

另外，作为获取每个扫描行名称的一种方法，我编写了以下代码：

import liblas
from liblas import file as lasfile

# Get parameters
las_file = r"E:\Testing\00101.las"

f = lasfile.File(las_file, mode='r')

line_list = []
counter = 0
for p in f:
    line_num = p.data[18]
    if line_num not in line_list:
        line_list.append(line_num)
    counter += 1
print line_list

结果会出现以下错误：

Traceback (most recent call last):
  File "D:\Tools\Python_Scripts\point_info.py", line 46, in <module>
    line_num = p.data[18]
  File "C:\Python27\ArcGIS10.1\lib\site-packages\liblas\point.py", line 560, in get_data
    length = self.header.data_record_length
  File "C:\Python27\ArcGIS10.1\lib\site-packages\liblas\point.py", line 546, in get_header
    return header.Header(handle=core.las.LASPoint_GetHeader(self.handle))
WindowsError: [Error -529697949] Windows Error 0xE06D7363

是否有人对存储在las点/头中的行号有更多的了解？有人能解释一下错误吗？在我得到错误之前，它似乎分配了将近2gb的ram。我正在使用win xp，所以我猜这是一个内存错误，但我不明白为什么访问这个“数据”字段占用内存。任何帮助都是非常感谢的。

Answer 1

我没有假装自己是这方面的专家，但我对地理信息系统的数据很感兴趣，所以这引起了我的兴趣。我在Fedora 19系统上安装了liblas及其依赖项，并使用liblas附带的示例数据文件。

使用您的代码，我遇到了同样的问题，看我所有的内存被吞噬。我不知道为什么会发生这种情况--也许是不必要的引用在周围徘徊，阻止垃圾收集器如我们所希望的那样工作。这个可能是可以修复的，但我不会尝试的。

我确实注意到了liblas模块的一些有趣的特性，并决定尝试它们。我相信你能得到你想要的数据。

打开文件后，请查看标题中的XML描述。

h = f.get_header()
print(h.get_xml())

很难看到(可以随意使用xml.dom.minidom或lxml.etree)，但在我的示例文件中，它显示了点数据的字节布局(我的也有28个字节)。在我这里，偏移量18是分配给点源ID的一个短的(2个字节)。您应该能够用p.data[18:19]、p.get_data()[18:19]、p.point_source_id或p.get_point_source_id()检索它。不幸的是，data引用占用了内存，而p.point_source_id有一个错误(bug修复、拉请求提交给开发人员)。如果我们将您的代码更改为使用最后一种访问方法，那么一切似乎都正常。因此，在您的for循环中尝试以下操作：

for p in f:
    line_num = p.get_point_source_id()
    if line_num not in line_list:
        line_list.append(line_num)
    counter += 1

请注意，

counter == h.get_count()

如果你只想要一组唯一的点源ID值.

line_set = set(p.get_point_source_id() for p in f)

希望您的数据值也可以作为p.get_point_source_id()提供。在评论中让我知道它是如何为你工作的。干杯!