Cloud Optimized Point Cloud(云优化点云)格式介绍构建解析及可视化（1）

提示：文章写完后，目录可以自动生成，如何生成可参考右边的帮助文档

---

前言

在现代的数字化世界中，点云数据成为了许多领域中重要的数据源，例如地理信息系统、遥感、地质勘探和三维建模等。为了更有效地存储和处理点云数据，出现了许多不同的点云格式。本文将重点介绍一种近两年发布的一种大数据组织的点云格式——Cloud Optimized Point Cloud（COPC 云优化点云）。该格式是2021年发布的一种符合LAS标准的以八叉树形式组织的一种大数据点云格式，COPC的数据组织与EPT点云（Potree所使用格式）类似。不同的是EPT格式的元数据、八叉树结构及各分块点云位于不同文件，而COPC的头信息，树结构描述及数据存储在单个.laz文件中，COPC仅支持LAZ点格式，不支持bin二进制格式。

一、COPC格式实现

COPC可以被看作是对传统的LAS（LiDAR数据交换标准）格式的改进和扩展。相比于LAS格式，COPC提供了更高的数据压缩率和更快的数据访问性能。这得益于COPC的层次结构以及数据压缩和索引技术的应用。八叉树组织的COPC点云与LAZ1.4标准点云的区别有以下方面：

1. copc格式必须包含LAS PDRF 6、7或8格式数据
2. copc格式必须包含信息VLR
3. copc格式必须包含层次结构VLR

信息VLR

必须是文件中的第一个VLR，必须从文件头偏移量375开始，类实现如下：

	class CopcInfo
	{
	public:
	public:
		static const int VLR_OFFSET = 375;     // COPC信息VLR必须在偏移量为375的LAS头块后的第一个VLR.
		static const int VLR_SIZE_BYTES = 160; // COPC信息VLR的数据结构大小为160字节, https://copc.io/

		CopcInfo() = default;

		CopcInfo(const lazperf::copc_info_vlr& copc_info_vlr);

		lazperf::copc_info_vlr ToLazPerf(const CopcExtent& gps_time) const;

		std::string ToString() const;
		friend std::ostream& operator<<(std::ostream& os, CopcInfo const& value)
		{
			os << value.ToString();
			return os;
		}

		double center_x{ 0 };
		double center_y{ 0 };
		double center_z{ 0 };
		double halfsize{ 0 };
		double spacing{ 0 };
		uint64_t root_hier_offset{ 0 };
		uint64_t root_hier_size{ 0 };
	};

层级结构VLR

类似与EPT格式，COPC存储的层次结构信息可以在读取时定位到特定的八叉树节点，与EPT类似，该层次结构可能为分页存储，应至少包含一个层次结构页。
VLR数据由一个或多个层次结构页组成，可有如下类表示：

/// @class VoxelKey
	/// @brief 节点键值，深度，x，y，z
	/// @note  
	class DATA_EXPORT VoxelKey
	{
	public:
		using Vector3 = lax::Vector3d;
		using Box = lax::BoundingBoxd;
		VoxelKey(int32_t d, int32_t x, int32_t y, int32_t z) : d(d), x(x), y(y), z(z) {}
		VoxelKey() : VoxelKey(-1, -1, -1, -1) {}
		VoxelKey(const std::string& key) { fromString(key); }
		VoxelKey(const std::vector<int32_t>& vec)
		{
			if (vec.size() != 4)
				throw std::runtime_error("Vector size must be 4.");
			d = vec[0];
			x = vec[1];
			y = vec[2];
			z = vec[3];
		};

		static VoxelKey InvalidKey() { return VoxelKey(); }
		static VoxelKey RootKey() { return VoxelKey(0, 0, 0, 0); }

		bool IsValid() const { return d >= 0 && x >= 0 && y >= 0 && z >= 0; }

		std::string ToString() const;
		void fromString(const std::string& str);
 
		/// @brief  根据方位[0,7]返回对应的八叉树节点键值
		VoxelKey Bisect(const uint64_t& direction) const;
		std::vector<VoxelKey> GetChildren() const;
 
		/// @brief  当前体素的在层次结构中的父节点
		VoxelKey GetParent() const;
 
		/// @brief  当前提取在指定深度的父节点
		VoxelKey GetParentAtDepth(int32_t depth) const;
 
		/// @brief  从当前节点到根节点的所有父节点,可选是否包含当前节点自身
		std::vector<VoxelKey> GetParents(bool include_self = false) const;
 
		/// @brief  检测当前节点是否为指定体素的子节点
		bool ChildOf(VoxelKey parent_key) const;

		/// @brief  根据key和las头信息创建包围盒   
		Box BoxFrom(const VoxelKey& key, const las::LasHeader& header) const;

		/*空间查找功能*/
		// Definitions taken from https://shapely.readthedocs.io/en/stable/manual.html#binary-predicates
		bool Intersects(const las::LasHeader& header, const Box& box) const;
		bool Contains(const las::LasHeader& header, const Box& vec) const;
		bool Contains(const las::LasHeader& header, const Vector3& point) const;
		bool Within(const las::LasHeader& header, const Box& box) const;
		bool Crosses(const las::LasHeader& header, const Box& box) const;

		double Resolution(const las::LasHeader& header, const CopcInfo& copc_info) const;
		static double GetResolutionAtDepth(int32_t d, const las::LasHeader& header, const CopcInfo& copc_info);

		int32_t d; ///< 深度
		int32_t x; ///< x方向位置
		int32_t y; ///< y方向位置
		int32_t z; ///< z方向位置
}；

    /// @class Entry
	/// @brief Entry类是Node、Page的基类
	/// @note  
	class Entry
	{
	public:
		static const int ENTRY_SIZE = 32;

		Entry() : offset(-1), byte_size(-1), key(VoxelKey::InvalidKey()), point_count(-1) {};
		Entry(VoxelKey key, uint64_t offset, int32_t size, int32_t point_count)
			: offset(offset), byte_size(size), key(key), point_count(point_count) {};
		Entry(const Entry& other)
			: offset(other.offset), byte_size(other.byte_size), key(other.key), point_count(other.point_count) {};
		virtual bool IsValid() const { return offset >= 0 && byte_size >= 0 && key.IsValid(); }
		virtual bool IsPage() const { return IsValid() && point_count == -1; }

		std::string ToString() const
		{
			std::stringstream ss;
			ss << "Entry " << key.ToString() << ": off=" << offset << ", size=" << byte_size << ", count=" << point_count
				<< ", is_valid=" << IsValid();
			return ss.str();
		}
		friend std::ostream& operator<<(std::ostream& os, Entry const& value)
		{
			os << value.key.ToString();
			os << value.ToString();
			return os;
		}

		void Pack(std::ostream& out_stream)
		{
			out_stream.write(reinterpret_cast<char*>(&key.d), sizeof(key.d));
			out_stream.write(reinterpret_cast<char*>(&key.x), sizeof(key.x));
			out_stream.write(reinterpret_cast<char*>(&key.y), sizeof(key.y));
			out_stream.write(reinterpret_cast<char*>(&key.z), sizeof(key.z));

			out_stream.write(reinterpret_cast<char*>(&offset), sizeof(offset));
			out_stream.write(reinterpret_cast<char*>(&byte_size), sizeof(byte_size));
			out_stream.write(reinterpret_cast<char*>(&point_count), sizeof(point_count));
		}

		static Entry Unpack(std::istream& in_stream)
		{
			VoxelKey key;
			in_stream.read(reinterpret_cast<char*>(&key.d), sizeof(key.d));
			in_stream.read(reinterpret_cast<char*>(&key.x), sizeof(key.x));
			in_stream.read(reinterpret_cast<char*>(&key.y), sizeof(key.y));
			in_stream.read(reinterpret_cast<char*>(&key.z), sizeof(key.z));

			uint64_t offset;
			in_stream.read(reinterpret_cast<char*>(&offset), sizeof(offset));
			int32_t size;
			in_stream.read(reinterpret_cast<char*>(&size), sizeof(size));
			int32_t point_count;
			in_stream.read(reinterpret_cast<char*>(&point_count), sizeof(point_count));

			return Entry(key, offset, size, point_count);
		}

		VoxelKey key;
		uint64_t offset;
		int32_t byte_size;
		int32_t point_count;

	protected:
		bool IsEqual(const Entry& rhs) const
		{
			return offset == rhs.offset && byte_size == rhs.byte_size && point_count == rhs.point_count && key == rhs.key;
		}
	};

class Node : public Entry
	{
	public:
		Node(const Entry& e, const VoxelKey& page_key = VoxelKey::InvalidKey()) : Entry(e), page_key(page_key) {};
		Node() : Entry() {};

		bool operator==(const Node& rhs) { return IsEqual(rhs); }
		VoxelKey page_key{};
	};

class Page : public Entry
	{
	public:
		Page(Entry e) : Entry(std::move(e)) {};
		Page(VoxelKey key, int64_t offset, int32_t byte_size) : Entry(key, offset, byte_size, -1) {};
 
		bool IsValid() const override { return (loaded || (offset >= 0 && byte_size >= 0)) && key.IsValid(); }
		bool IsPage() const override { return IsValid() && point_count == -1; }

		bool loaded = false;

		friend bool operator==(const Page& lhs, const Page& rhs) { return lhs.IsEqual(rhs) && lhs.loaded == rhs.loaded; }
	};

二、可视化效果预览

本文简单先介绍COPC格式的组成，有时间对COPC格式文件的构建及解析进行讲解：
先放一个COPC点云可视化效果，实验点云sofi1.88GB，3亿6千多万点，点云LAS头信息如下图：

copc格式点云渲染

使用osg三维引擎点云的调度渲染，基本可稳定在100帧左右。