使用Python读取las点云,写入las点云,无损坐标精度

目录

  • 1 为什么要写这个博文
  • 2 提出一些关键问题
  • 3 给出全部代码
    • 安装依赖
    • 源码(laspy v2.x)

1 为什么要写这个博文

搜索使用python读写las点云数据,可以找到很多结果。但是! 有些只是简单的demo,且没有发现/说明可能遇到的问题;有些晦涩难懂,不够实用主义;有些须付费观看,没有开源精神。

本人在使用laspy v2.3读写点云文件时,着实被坐标精度问题难到了,顺便就仔细学习了下las格式到底是什么来头。

您猜怎么着,如果能打开这个网址,都在这里面详细说明了:laspy document

我还是种个树吧,直接拔了就能用。文末贴源码。

2 提出一些关键问题

  1. 注意:
    * 实测版本 laspy v2.3,2.x版本应该都可用,但不适合1.x版本。
    * las 存储数据时,需要设置 scales 和 offsets,否则会出现精度问题。
    * las 存储颜色时,数值类型为 16 位无符号整型。rgb = (normal_rgb * 65535).astype(np.uint16)

  2. las 格式原生支持的属性字段,与las版本密切相关。官方说明:https://laspy.readthedocs.io/en/latest/intro.html#point-records

  3. 对 scales 和 offsets 的理解:
    比例scales 表明数据的准确性。 0.001 是毫米精度。这意味着如果您的坐标是 0.123456,它将被限制为 0.123。
    偏移offset 的目的是避免整数溢出。假设您要存储 123456789.123。在 LAS 文件中,您实际上将存储一个整数:123456789123,该整数将在读取时使用比例因子转换为 123456789.123。但 123456789123 比 32 位整数所能存储的要大得多。因此存储时,将该值偏移 123450000,实际存的是6789123。
    (6789123 * 0.001 + 123450000 = 123456789.123)
    段落出处

3 给出全部代码

安装依赖

使用 pip (# 选择一个安装就好):

# Install _without_ LAZ support
pip install laspy

# Install with LAZ support via lazrs
pip install laspy[lazrs]

# Install with LAZ support via laszip
pip install laspy[laszip]

# Install with LAZ support via both lazrs & laszip
pip install laspy[lazrs,laszip]

使用 conda (# 选择一个安装就好):

conda install -c conda-forge laspy
conda install -c conda-forge lazrs-python

源码(laspy v2.x)

三个工具函数:
read_las_fit() : 读取 las 文件
write_las_fit(): 保存 las 文件
get_las_header_attrs() : 获取不同 las 版本支持的固有属性
备注:_fit 的意思是,可以支持各种属性信息的读取和写入

import laspy
import numpy as np


def read_las_fit(filename, attrs=None):
    """
    读取 las 文件,获取三维坐标 xyz, 颜色 rgb, 属性 attr_dict。当文件没有 RGB 信息时,返回全0的 RGB 信息
    Args:
        filename: <str> las 文件路径
        attrs: <list> 需要额外获取的属性信息 如 ['label']

    Returns:
        xyz, rgb, attr_dict
    """
    if attrs is None:
        attrs = []

    # 默认返回 scales, offsets ,合并 ["scales", "offsets"]
    attrs = list(set(attrs + ["scales", "offsets"]))

    # 读取点云
    inFile = laspy.read(filename)
    # inFile.point_format.dimensions可以获取所有的维度信息
    N_points = len(inFile)
    x = np.reshape(inFile.x, (N_points, 1))
    y = np.reshape(inFile.y, (N_points, 1))
    z = np.reshape(inFile.z, (N_points, 1))
    xyz = np.hstack((x, y, z))
    # TODO 注意。如果是大写的 X Y Z,需要转换后才是真实坐标: real_x = scale[0] * inFile.X + offset[0]

    # 初始化 rgb 全是 0
    rgb = np.zeros((N_points, 3), dtype=np.uint16)
    if hasattr(inFile, "red") and hasattr(inFile, "green") and hasattr(inFile, "blue"):
        r = np.reshape(inFile.red, (N_points, 1))
        g = np.reshape(inFile.green, (N_points, 1))
        b = np.reshape(inFile.blue, (N_points, 1))
        # i = np.reshape(inFile.Reflectance, (N_points, 1))
        rgb = np.hstack((r, g, b))
    else:
        print(f"注意:{filename.split('/')[-1]} 没有RGB信息,返回全0的RGB信息!")

    # 组织其他属性信息
    attr_dict = {}
    for attr in attrs:
        # 先判断 header 中是否有该属性
        if hasattr(inFile.header, attr):
            value = getattr(inFile.header, attr)
            if hasattr(value, "array"):
                attr_dict[attr] = np.array(value)
            else:
                attr_dict[attr] = value
        # 再判断 是否为额外属性
        elif hasattr(inFile, attr):
            value = getattr(inFile, attr)
            if hasattr(value, "array"):
                attr_dict[attr] = np.array(value)
            else:
                attr_dict[attr] = value
        else:
            attr_dict[attr] = None
            print(f"注意:{filename.split('/')[-1]} 没有属性 {attr} 信息!")

    return xyz, rgb, attr_dict


def write_las_fit(out_file, xyz, rgb=None, attrs=None):
    """
    将点云数据写入 las 文件,支持写入 坐标xyz, 颜色rgb, 属性attrs
    Args:
        out_file: 输出文件路径
        xyz: 点云坐标 ndarray (N, 3)
        rgb: 点云颜色 ndarray (N, 3)
        attrs:
            固有属性:file_source_id, gps_time, Intensity, Number of Returns, ....
            额外属性:label, pred, ...
            注意:如果不传入 scales 和 offsets,则会自动计算
    Returns:

    """
    if attrs is None:
        attrs = {}

    # 1. 创建 las 文件头。point_format和version决定了las支持哪些固有属性
    # 详情见 https://pylas.readthedocs.io/en/latest/intro.html?highlight=red#point-records
    header = laspy.LasHeader(point_format=7, version="1.4")  # 7 支持rgb

    # 自动计算 scales 和 offsets,确保坐标精度无损
    # https://stackoverflow.com/questions/77308057/conversion-accuracy-issues-of-e57-to-las-in-python-using-pye57-and-laspy
    if "offset" not in attrs:
        min_offset = np.floor(np.min(xyz, axis=0))
        attrs["offset"] = min_offset
    if "scales" not in attrs:
        attrs["scales"] = [0.001, 0.001, 0.001]  # 0.001 是毫米精度

    # 初始化一些需要保存的属性值。如果是固有属性,直接赋值; 如果是额外属性,添加到 header 中, 后续赋值
    extra_attr = []
    for attr, value in attrs.items():
        if hasattr(header, attr):  # 设置固有的属性的值, 如 scales, offsets
            header.__setattr__(attr, value)
        else:  # 添加额外属性,在 las 初始化后赋值
            header.add_extra_dim(laspy.ExtraBytesParams(name=attr, type=np.float32))
            extra_attr.append(attr)

    # 2. 创建 las 文件
    las = laspy.LasData(header)

    # 添加xyz坐标
    las.x = xyz[:, 0]
    las.y = xyz[:, 1]
    las.z = xyz[:, 2]

    # 添加RGB颜色,如果是归一化的颜色,则需要乘以 65535,转为 uint16
    if rgb is not None:
        if np.max(rgb) <= 1:
            rgb = (rgb * 65535).astype(np.uint16)  # 65535 = 2^16 - 1, las存储颜色是16位无符号整型
        las.red = rgb[:, 0]
        las.green = rgb[:, 1]
        las.blue = rgb[:, 2]

    # 添加额外属性
    for attr in extra_attr:
        # 当 value 是 n * 1 的 ndarray 时,转换为 1 维数组
        value = attrs[attr]
        if value.ndim == 2 and value.shape[1] == 1:
            value = value.flatten()
        las[attr] = value

    # 保存LAS文件
    las.write(out_file)


def get_las_header_attrs(point_format=7, version="1.4"):
    """
    根据 point_format 和 version 获取 las 文件的 header 属性
    说明文档:https://laspy.readthedocs.io/en/latest/intro.html#point-records
    Args:
        point_format: 点格式
        version: 版本

    Returns:

    """
    dimensions = []
    header = laspy.LasHeader(point_format=point_format, version=version)  # 7 支持rgb
    for dim in header.point_format.dimensions:
        dimensions.append(dim.name)
    return dimensions


if __name__ == '__main__':
    # 测试1: 获取 las 文件头属性
    fields = get_las_header_attrs(7, "1.4")
    print(f"point_format=7, version=1.4, 头文件包含字段: {fields}")

    # 测试2: 读取LAS文件
    read_las_path = "/path_2_data/one_point_cloud.las"
    xyz_, rgb_, attrs_ = read_las_fit(read_las_path, ["scales", "offsets"])
    print(attrs_)

    # 测试3: 写入LAS文件
    save_las_path = "/path_2_data/one_point_cloud_fit.las"
    write_las_fit(save_las_path, xyz_, rgb_, {
        # "scales": attrs_["scales"],
        # "offsets": attrs_["offsets"]
    })

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mfbz.cn/a/608117.html

如若内容造成侵权/违法违规/事实不符,请联系我们进行投诉反馈qq邮箱809451989@qq.com,一经查实,立即删除!

相关文章

力扣138. 随机链表的复制

Problem: 138. 随机链表的复制 文章目录 题目描述思路及解法复杂度Code 题目描述 思路及解法 1.创建Map集合Map<Node, Node> map;创建指针cur指向head&#xff1b; 2.遍历链表将cur作为键&#xff0c;new Node(cur.val)作为值&#xff0c;存入map集合&#xff1b; 3.再次…

[机器学习系列]深入探索回归决策树:从参数选择到模型可视化

目录 一、回归决策树的参数 二、准备数据 三、构建回归决策树 (一)拟合模型 (二)预测数据 (三)查看特征重要性 (四)查看模型拟合效果 (五) 可视化回归决策树真实值和预测值 (六)可视化决策树并保存 部分结果如下&#xff1a; 一、回归决策树的参数 DecisionTreeRegress…

英特尔StoryTTS:新数据集让文本到语音(TTS)表达更具丰富性和灵感

DeepVisionary 每日深度学习前沿科技推送&顶会论文分享&#xff0c;与你一起了解前沿深度学习信息&#xff01; 英特尔StoryTTS&#xff1a;新数据集让文本到语音&#xff08;TTS&#xff09;表达更具丰富性和灵感 引言&#xff1a;探索文本表达性在语音合成中的重要性 …

【深耕 Python】Quantum Computing 量子计算机(3)重要数学公式一览

写在前面 往期量子计算机博客&#xff1a; 【深耕 Python】Quantum Computing 量子计算机&#xff08;1&#xff09;图像绘制基础 【深耕 Python】Quantum Computing 量子计算机&#xff08;2&#xff09;绘制电子运动平面波 正文 偏微分&#xff1a; 交换关系&#xff…

NOIP,CSP-J,CSP-S——图

一、图的基本概念 图是顶点和边的集合 1、无向图: 每一条边都是无方向的 2、有向图: 每一条边都是有方向的 3、完全图: 任意两个顶点都有一条边相连接; 4、结论 若n个顶点的无向图有n(n-1)/2条边,称为无向完成图; 若n个顶点的有向图有n(n-1)条边,称为有向完成图…

华为eNSP Pro模拟器下载(普通账号可用)

好消息&#xff01;华为终于开放了普通账号使用权限&#xff01; 安装教程下载后见《指导手册-eNSP Pro V100R001C00.pdf》 华为eNSP Pro模拟器下载&#xff08;普通账号可用&#xff09; 下载地址 华为eNSP Pro模拟器下载&#xff08;普通账号可用&#xff09; - 下一朵云 …

cannot import name ‘ForkProcess‘ from ‘multiprocessing.context‘问题解决

问题描述 cannot import name ForkProcess from multiprocessing.context 问题原因 ForkContext用于Unix系统。SpawnContext可以在 Windows 环境中使用 解决方案 改成SpawnProcess就可以运行了 将原来的ForkProcess修改为SpawnProcess wrappers.py脚本&#xff0c;下面的代…

Android MediaCodec 简明教程(七):使用 MediaCodec 解码到 OES 纹理上

系列文章目录 Android MediaCodec 简明教程&#xff08;一&#xff09;&#xff1a;使用 MediaCodecList 查询 Codec 信息&#xff0c;并创建 MediaCodec 编解码器Android MediaCodec 简明教程&#xff08;二&#xff09;&#xff1a;使用 MediaCodecInfo.CodecCapabilities 查…

【Linux】-Linux用户和权限[3]

一、认知root用户 1、root用户&#xff08;超级管理员&#xff09; 无论是Windows、MacOS、Linux均采用多用户的管理模式进行权限管理。 在Linux系统中&#xff0c;拥有最大权限的账户为&#xff1a;root&#xff08;超级管理员&#xff09; root用户拥有最大的系统操作权限…

python 和 MATLAB 都能绘制的母亲节花束!!

hey 母亲节快到了&#xff0c;教大家用python和MATLAB两种语言绘制花束~这段代码是我七夕节发的&#xff0c;我对代码进行了简化&#xff0c;同时自己整了个python版本 MATLAB 版本代码 function roseBouquet_M() % author : slandarer% 生成花朵数据 [xr,tr]meshgrid((0:24).…

我们的小程序每天早上都白屏,真相是。。。

大家好&#xff0c;我是程序员鱼皮。最近我们在内测一款面试刷题小程序&#xff0c;没错&#xff0c;就是之前倒下的 “面试鸭”&#xff01; 在我们的内测交流群中&#xff0c;每天早上都会有同学反馈&#xff1a;打开小程序空白&#xff0c;没任何内容且登录不上。 然后过了…

感知机简介

感知机简介 导语感知机简单逻辑电路实现权重和配置与/或/与非与门实现与非门实现或门实现 线/非线性单/多层感知机异或 总结参考文献 导语 学习感知机有助于更好的理解深度学习的神经元、权重等概念&#xff0c;感知机的结构和概念很简单&#xff0c;只要学过基本线性代数、数…

Web 安全基础理论

Web 安全基础理论 培训、环境、资料、考证 公众号&#xff1a;Geek极安云科 网络安全群&#xff1a;624032112 网络系统管理群&#xff1a;223627079 网络建设与运维群&#xff1a;870959784 移动应用开发群&#xff1a;548238632 短视频制作群&#xff1a; 744125867极安云…

Star-CCM+通过将所有部件创建一个区域的方式分配至区域后子区域的分离,子区域材料属性的赋值,以及物理连续体的创建方法介绍

前言 上次介绍了将零部件分配至区域的方法与各个方法之间的区别&#xff0c;本文将继续上次的讲解&#xff0c;将其中的“将所有部件分配至一个区域”的应用进行补充。 如下图所示&#xff0c;按照将所有部件创建一个区域的方式分配至区域后&#xff0c;在区域下就会有一个区域…

Marin说PCB之如何快速打印输出整板的丝印位号图?

当小编我辛辛苦苦加班加点的把手上的板子做到投板评审状态的时候&#xff0c;坐在我旁边的日本同事龟田小郎君说让我把板子上的丝印也要调一下&#xff0c;我当时就急了&#xff0c;这么大的板子&#xff0c;将近1W多PIN 了都&#xff0c;光调丝印都要老半天啊&#xff0c;而且…

【ytb数据采集器】按关键词批量爬取视频数据,界面软件更适合文科生!

一、背景介绍 1.1 爬取目标 用Python独立开发的爬虫工具&#xff0c;作用是&#xff1a;通过搜索关键词采集油管的搜索结果&#xff0c;包含14个关键字段&#xff1a;关键词,页码,视频标题,视频id,视频链接,发布时间,视频时长,频道名称,频道id,频道链接,播放数,点赞数,评论数…

MATLAB 点云随机赋色 (68)

MATLAB 点云随机赋色 (68) 一、算法介绍二、算法介绍1.代码2.结果三、数据链接一、算法介绍 读取的点云本身带有颜色信息,有时我们需要为每个点随机赋予一种颜色,下面是具体效果和实现代码,以及使用的数据: 二、算法介绍 1.代码 代码如下(示例): % 读取点云文件 f…

linux中进程相关概念(一)

什么是程序&#xff0c;什么是进程&#xff0c;有什么区别&#xff1f; 程序是静态的概念&#xff0c;当我们使用gcc xxx.c -o pro进行编译时&#xff0c;产生的pro文件&#xff0c;就是一个程序。 进程是程序的一次运行活动&#xff0c;通俗点就是说程序跑起来了就是进程。 …

TypeScript学习日志-第二十一天(声明文件d.ts)

声明文件d.ts 在使用 Typescript 并使用第三方库 的时候 我们会发现会有很多的提示或补全&#xff0c;这都是声明文件起的作用&#xff0c;但是有写冷门的第三方库是没有声明文件的&#xff0c;这时候引用就会报错&#xff0c;我们就使用 express 库作为例子来展示一下&#x…

马蹄集oj赛(双周赛第二十六次)

目录 斐波那契数列的组合 三国杀 数列分段 小码哥的跳棋游戏新编 能量供应 小码哥爱数字 最小串 小船过河 摘果子 泼墨淋漓 很重的枪 小码哥的布阵指挥 斐波那契数列的组合 #include<bits/stdc.h> using namespace std;// 斐波那契数列 1 1 2 3 5 8 13 21 34…
最新文章