点云处理入门--PointNetPointNet++论文与代码详解

基础知识

点云数据：
点云是一种通过三维扫描设备或计算机图形学技术获取的三维空间数据，通常由一系列点组成，每个点包含其在三维空间中的坐标（如 x,y,z），有时还可能包含颜色、强度等附加信息。
介绍几种常见的点云存储格式：
1、XYZ
最简单的点云文件格式，通常是一个纯文本文件，每行表示一个点，包含三个数值（x,y,z坐标）
2、PLY
文件可以是纯文本或二进制格式，支持存储点的坐标、颜色、法线等附加信息，支持点云和网格（mesh）数据。
3、OBJ
主要用于存储三维模型（网格），但也可用于点云数据。文件格式为纯文本，支持顶点坐标、纹理坐标、法线等信息，不支持颜色。
Mesh和点云的转换：
Mesh（网格）是一种由顶点、边和面组成的三维模型，而点云是离散的点集合。虽然它们在形式上有所不同，但可以通过一些算法相互转换：

Mesh 转点云：可以通过泊松采样（Poisson Sampling）、均匀采样（Uniform Sampling）等方法从网格表面生成点云。
点云转 Mesh：可以通过泊松重建（Poisson Surface Reconstruction）、Marching Cubes 等算法从点云生成网格。

PointNet论文讲解

论文地址：https://arxiv.org/abs/1612.00593
code：https://github.com/charlesq34/pointnet
pytorch实现：https://github.com/fxia22/pointnet.pytorch

简单来说，点云数据就是一堆坐标点。

如何处理点云数据？

体素化方法：将点云数据转换为规则的三维网格（体素），每个体素包含一定范围内的点云数据，对每个体素内的点云进行特征提取（平均值、最大值等），然后利用3D卷积神经网络（CNN）进行处理。
多视图方法：将点云从多个视角投影到二维平面上，生成多个二维图像，然后使用传统的2D卷积网络进行处理。
这两种方法都没有直接利用点云原始数据，而是经过转换后输入网络。

PointNet直接使用原始点云数据，近年来越来越多的方法这样做。
直接对点云特征进行学习，就需要解决两个问题：点云数据的无序性、旋转的不变性。

1、数据的无序性

点云中的点的排列顺序是随机的，不具有固定的结构。点云中的点可以任意交换位置，而不会影响其对物体或场景的描述。这点和图像数据是不同的，而传统的深度学习方法（如卷积神经网络）通常依赖于输入数据的固定结构（如图像的像素排列）。
这就要求模型能够不受数据排列顺序的影响，提取到准确的特征。
如何能够不受排列顺序影响呢？我们可以使用一些函数来处理，比如max，min，sum，avg等，这些函数的处理结果不受序列顺序影响。但是还有一个问题，如果说一个点云nx3，我们对每一个维度采取max处理，得到1x3的特征（这反映了这个点云在三个维度上的最大值），这样做损失的信息太多了，输出的特征仅仅继承了三个坐标轴上的最大特征。
所以我们需要更多的信息，我们首先将特征映射到高维空间，简单来说将点云数据通过MLP扩充维度，如nx1024，然后取点云序列在每一个维度的最大值1x1024，组成的向量代表整个点云序列的特征，这样一来，避免了序列顺序的影响，同时也增加了更多的特征信息。
在这里插入图片描述

2、旋转不变性

点云在经过旋转变换后，其特征表示应保持不变。
换句话说，如果我们在空间内旋转一个物体，比如一把伞，那么他的点云数据（坐标）都已经被改变，我们需要模型能够包容这种旋转操作，依然识别出这是一把伞。
PointNet保证旋转不变性的做法是–设计了T-Net来学习点云的旋转。
T-Net是一个小型的神经网络。主要目的是通过学习一个仿射变换矩阵，对输入点云或特征空间进行对齐，从而减少姿态变化对模型性能的影响。
具体做法：输入点云为nx3，先扩充维度（三次卷积操作，卷积核为1x1，通道数分别是64，128，1024）得到nx1024向量，然后对其进行最大池化（确保模型对点云无序性具有不变性），得到特征1x1024，然后先通过两层全连接层，分别映射到 512 维和 256 维，最后把256维的特征映射到kxk的矩阵（在网络架构中有两处，一处是3x3，一处是64x64）。这个映射也很简单，就是再过一个全连接，维度变成 $k^2$ ，然后reshape成kxk的矩阵就行。最后这个矩阵就是学习了点云旋转的特征，将其乘以输出的点云即可。

PointNet代码

PointNet分为两个版本，点云分类和点云分割
首先看分类，分类的网络结构图应该如下图所示，红色箭头代表跳过中间的feature_transform
在这里插入图片描述
输入为一个 $\times 3$ 的点云，代码里实际的输入格式为Bx3xN
先过一个inout_transform，也就是论文中说的T-Net(这部分怎么处理的，可以看下面的代码注释)
然后得到特征为Bx3xN，然后过一个卷积，变为Bx64xN，也就是图中的mlp(64,64)，两个64分别代表之前的通道数和之后的通道数，
然后按照我画的红色箭头，再过一个mlp(64,128,1024)，也就是三个卷积，最后输出Bx1024xN，然后做一个max pool变为Bx1024
然后用全卷积神经网络nn.Linear，将1024–>512–>256–>k,也就是图中的mlp(512,256,k)
k就是分类的类别数量，最后再过一个softmax得到分类结果。
比如分两类，最后输出的Bx2也就是每个点云2个概率值，分别代表该点云是这两类的概率。
这就是分类架构。
分割与分类不同，分类预测整个点云的概率，如果一共有五类，那么就输出五个数，代表该点云是这五类的概率。
而分割需要预测每一个点的类别，要输出Nx5个数字，代表n个点是这个五类的概率。
分割的网络结构如下图所示。
在这里插入图片描述

与分类相比有几个地方不一样。
首先就是第一次mlp(64,64)之后，需要再过一个feature_transformer，这次仿射变换矩阵是64x64的，过了之后保持大小形状不变，然后把这个Bx64xN保存下来，后面继续正常往下走，直到走到最后max pool之后输出一个全局特征Bx1024，这个时候我们使用复制操作，将其复制为Bx1024xN,也就是复制n个，然后把这个特征与前面保存下来的那个Bx64xN，做一个拼接，得到Bx1088xN（64+1024），拿这个特征去做一维卷积，将通道数由1088一直降到128，也就是得到Bx128xN，最后再来一个一维卷积，直接维度降到k，BxkxN，k代表分类的类别。然后做一个softmax，得到每个点的k个类别概率（softmax之前要转换一下，可以参考代码操作）。

代码以及注释（pytorch版本）

from __future__ import print_function
import torch
import torch.nn as nn
import torch.nn.parallel
import torch.utils.data
from torch.autograd import Variable
import numpy as np
import torch.nn.functional as Fclass STN3d(nn.Module): #空间变换网络（Spatial Transformer Network）def __init__(self):super(STN3d, self).__init__()self.conv1 = torch.nn.Conv1d(3, 64, 1)self.conv2 = torch.nn.Conv1d(64, 128, 1)self.conv3 = torch.nn.Conv1d(128, 1024, 1)self.fc1 = nn.Linear(1024, 512)self.fc2 = nn.Linear(512, 256)self.fc3 = nn.Linear(256, 9)self.relu = nn.ReLU()self.bn1 = nn.BatchNorm1d(64)self.bn2 = nn.BatchNorm1d(128)self.bn3 = nn.BatchNorm1d(1024)self.bn4 = nn.BatchNorm1d(512)self.bn5 = nn.BatchNorm1d(256)def forward(self, x): # x:BxDxN  eg:(32, 3, 2500) 32个点云，每个点云2500个点，每个点三个坐标xyzbatchsize = x.size()[0] # 输入点云数据 BxDxN  D=3x = F.relu(self.bn1(self.conv1(x))) # 一维卷积 通道数3-->64x = F.relu(self.bn2(self.conv2(x))) # 一维卷积 通道数64-->128x = F.relu(self.bn3(self.conv3(x))) # 一维卷积 通道数128-->1024x = torch.max(x, 2, keepdim=True)[0] # Bx1024xN 在维度2上取最大值 Bx1024x1 即在每个维度上对N个点取最大值x = x.view(-1, 1024) # 去除最后一列 Bx1024x = F.relu(self.bn4(self.fc1(x))) # 1024-->512x = F.relu(self.bn5(self.fc2(x))) # 512-->256x = self.fc3(x) # 256-->9# 为每一个batch生成一个单位矩阵iden = Variable(torch.from_numpy(np.array([1,0,0,0,1,0,0,0,1]).astype(np.float32))).view(1,9).repeat(batchsize,1)if x.is_cuda:iden = iden.cuda()x = x + iden # 变换矩阵与单位矩阵相加 这个操作可以理解为：防止网络乱学 引导其接近单位矩阵x = x.view(-1, 3, 3) # 调整为Bx3x3return xclass STNkd(nn.Module):#这个与上面的STN是一样的 这个输出64x64  上面输出3x3def __init__(self, k=64):super(STNkd, self).__init__()self.conv1 = torch.nn.Conv1d(k, 64, 1)self.conv2 = torch.nn.Conv1d(64, 128, 1)self.conv3 = torch.nn.Conv1d(128, 1024, 1)self.fc1 = nn.Linear(1024, 512)self.fc2 = nn.Linear(512, 256)self.fc3 = nn.Linear(256, k*k)self.relu = nn.ReLU()self.bn1 = nn.BatchNorm1d(64)self.bn2 = nn.BatchNorm1d(128)self.bn3 = nn.BatchNorm1d(1024)self.bn4 = nn.BatchNorm1d(512)self.bn5 = nn.BatchNorm1d(256)self.k = kdef forward(self, x):batchsize = x.size()[0]x = F.relu(self.bn1(self.conv1(x)))x = F.relu(self.bn2(self.conv2(x)))x = F.relu(self.bn3(self.conv3(x)))x = torch.max(x, 2, keepdim=True)[0]x = x.view(-1, 1024)x = F.relu(self.bn4(self.fc1(x)))x = F.relu(self.bn5(self.fc2(x)))x = self.fc3(x)iden = Variable(torch.from_numpy(np.eye(self.k).flatten().astype(np.float32))).view(1,self.k*self.k).repeat(batchsize,1)if x.is_cuda:iden = iden.cuda()x = x + idenx = x.view(-1, self.k, self.k)return xclass PointNetfeat(nn.Module):def __init__(self, global_feat = True, feature_transform = False):super(PointNetfeat, self).__init__()self.stn = STN3d()self.conv1 = torch.nn.Conv1d(3, 64, 1)self.conv2 = torch.nn.Conv1d(64, 128, 1)self.conv3 = torch.nn.Conv1d(128, 1024, 1)self.bn1 = nn.BatchNorm1d(64)self.bn2 = nn.BatchNorm1d(128)self.bn3 = nn.BatchNorm1d(1024)self.global_feat = global_featself.feature_transform = feature_transformif self.feature_transform:self.fstn = STNkd(k=64)def forward(self, x): # x:BxDxN D就是通道数 输入的时候D=3 代表xyz坐标n_pts = x.size()[2]trans = self.stn(x) # 做一个仿射变换  trans：Bx3x3x = x.transpose(2, 1) # BxDxN --> BxNx3x = torch.bmm(x, trans) # 与仿射变换的结果相乘 结果大小为BxNx3x = x.transpose(2, 1) # x:Bx3xNx = F.relu(self.bn1(self.conv1(x))) # 3-->64  Bx64xNif self.feature_transform:trans_feat = self.fstn(x) # 做一个仿射变换  trans_feat: 64x64x = x.transpose(2,1)x = torch.bmm(x, trans_feat) # 与仿射变换的结果相乘 x = x.transpose(2,1) # 结果大小为BxNx64else:trans_feat = Nonepointfeat = x # feature_transform=false : pointfeat Bx64xNx = F.relu(self.bn2(self.conv2(x))) # 64-->128  Bx128xNx = self.bn3(self.conv3(x)) # 128-->1024  Bx1024xNx = torch.max(x, 2, keepdim=True)[0] # 最大池化 Bx1024x1x = x.view(-1, 1024) # Bx1024if self.global_feat:return x, trans, trans_feat # x是编码后的特征 teans是仿射变换矩阵 trans_feat是else:x = x.view(-1, 1024, 1).repeat(1, 1, n_pts) # 使用复制 将x:Bx1024x1-->Bx1024xN# x(Bx1024xN)和pointfeat(Bx64xN)在维度1上拼接 结果大小为Bx1088xNreturn torch.cat([x, pointfeat], 1), trans, trans_featclass PointNetCls(nn.Module): # 分类def __init__(self, k=2, feature_transform=False):# 分类时 feature_transform=Falsesuper(PointNetCls, self).__init__()self.feature_transform = feature_transformself.feat = PointNetfeat(global_feat=True, feature_transform=feature_transform) self.fc1 = nn.Linear(1024, 512)self.fc2 = nn.Linear(512, 256)self.fc3 = nn.Linear(256, k)self.dropout = nn.Dropout(p=0.3) # 每个神经元有 30% 的概率被置零 最后一层使用self.bn1 = nn.BatchNorm1d(512)self.bn2 = nn.BatchNorm1d(256)self.relu = nn.ReLU()def forward(self, x): # BxDxNx, trans, trans_feat = self.feat(x) # 提取特征 x:Bx1024 trans:Bx3x3 trans_feat=Nonex = F.relu(self.bn1(self.fc1(x))) # Bx1024-->Bx512x = F.relu(self.bn2(self.dropout(self.fc2(x))))  # Bx512-->Bx256 x = self.fc3(x) # Bx256 --> k  要分几类 默认是2return F.log_softmax(x, dim=1), trans, trans_feat # 过一个log_softmax输出结果class PointNetDenseCls(nn.Module): # 分割def __init__(self, k = 2, feature_transform=False):super(PointNetDenseCls, self).__init__()self.k = kself.feature_transform=feature_transformself.feat = PointNetfeat(global_feat=False, feature_transform=feature_transform)self.conv1 = torch.nn.Conv1d(1088, 512, 1)self.conv2 = torch.nn.Conv1d(512, 256, 1)self.conv3 = torch.nn.Conv1d(256, 128, 1)self.conv4 = torch.nn.Conv1d(128, self.k, 1)self.bn1 = nn.BatchNorm1d(512)self.bn2 = nn.BatchNorm1d(256)self.bn3 = nn.BatchNorm1d(128)def forward(self, x):batchsize = x.size()[0]n_pts = x.size()[2]x, trans, trans_feat = self.feat(x) # 提取特征 x:Bx1088xN trans:Bx3x3 trans_feat=Bx64x64x = F.relu(self.bn1(self.conv1(x))) # Bx1088xN --> Bx512xNx = F.relu(self.bn2(self.conv2(x))) # Bx512xN --> Bx256xNx = F.relu(self.bn3(self.conv3(x))) # Bx256xN --> Bx128xNx = self.conv4(x) # Bx128xN --> BxkxNx = x.transpose(2,1).contiguous() # BxNxKx = F.log_softmax(x.view(-1,self.k), dim=-1) # BxNxK-->B*NxK 然后过一个log_softmaxx = x.view(batchsize, n_pts, self.k) # BxNxKreturn x, trans, trans_feat # x:BxNxK, trans:Bx3x3, trans_feat:Bx64x64def feature_transform_regularizer(trans):d = trans.size()[1]batchsize = trans.size()[0]I = torch.eye(d)[None, :, :]if trans.is_cuda:I = I.cuda()loss = torch.mean(torch.norm(torch.bmm(trans, trans.transpose(2,1)) - I, dim=(1,2)))return lossif __name__ == '__main__':sim_data = Variable(torch.rand(32,3,2500))trans = STN3d()out = trans(sim_data)print('stn', out.size())print('loss', feature_transform_regularizer(out))sim_data_64d = Variable(torch.rand(32, 64, 2500))trans = STNkd(k=64)out = trans(sim_data_64d)print('stn64d', out.size())print('loss', feature_transform_regularizer(out))pointfeat = PointNetfeat(global_feat=True)out, _, _ = pointfeat(sim_data)print('global feat', out.size())pointfeat = PointNetfeat(global_feat=False)out, _, _ = pointfeat(sim_data)print('point feat', out.size())cls = PointNetCls(k = 5)out, _, _ = cls(sim_data)print('class', out.size())seg = PointNetDenseCls(k = 3)out, _, _ = seg(sim_data)print('seg', out.size())