Spherical Voxelization

• 標(biāo)簽: voxelization
• AI 摘要: 文檔介紹了球面體素化的過程，包括重要的類和方法，如ConvertToSphericalVoxel和spherical_voxel_optimized，詳細(xì)解釋了參數(shù)及其作用。球面體素化通過將點云轉(zhuǎn)換為球面坐標(biāo)系，利用自適應(yīng)采樣權(quán)重來平衡不同緯度區(qū)域的點密度，從而有效捕捉幾何特征。文中還提到C++綁定的sv.compute函數(shù)，負(fù)責(zé)體素特征的計算與填充，確保在特征計算中考慮鄰近體素的信息。
• 最相關(guān)鏈接: https://github.com/CVLAB-Unibo/compass

Spherical Voxelization

參考鏈接：

• CVLAB-Unibo/compass: Repository containing the code of "Learning to Orient Surfaces by Self-supervised Spherical CNNs". (github.com)
• jonkhler/s2cnn: Spherical CNNs (github.com)
• qq456cvb/PRIN: Pointwise Rotation-Invariant Network (AAAI 2020) (github.com)

代碼組成：

ConvertToSphericalVoxel類：最高接口，實例化一個converter類，調(diào)用convert轉(zhuǎn)換局部點云
↓
spherical_voxel_optimized方法：在convert中調(diào)用，實現(xiàn)轉(zhuǎn)換，先轉(zhuǎn)換到球面坐標(biāo)系，然后進(jìn)行體素化
↓
spherical_voxel.compute方法：最終實現(xiàn)體素化，用pybind綁定C++代碼最終實現(xiàn)

ConvertToSphericalVoxel

from utils import geometry as ug

class ConvertToSphericalVoxel():

    """
    Convert point cloud to spherical voxel [beta = 2 * bandwidth,  alfa = 2 * bandwidth, num_radial_division].
    Alfa in [0, 2pi], Beta in [0, pi]
    """

    def __init__(self, bandwidth, radius_support, num_radial_division, num_points, random_sampling):

        self.bandwidth = bandwidth
        self.radius_support = radius_support
        self.num_radial_division = num_radial_division
        self.num_points = num_points
        self.random_sampling = random_sampling

    def __call__(self, point_cloud):

        features, pts_normed = ug.spherical_voxel_optimized(points=point_cloud,
                                                size_bandwidth=self.bandwidth,
                                                size_radial_divisions=self.num_radial_division,
                                                radius_support=self.radius_support,
                                                do_random_sampling=self.random_sampling,
                                                num_random_points=self.num_points)

        return features, pts_normed
……

參數(shù)解釋：

• bandwidth: 球面體素化的帶寬，通常用于定義球面信號的分辨率。它決定了角度方向上的采樣密度（球面坐標(biāo)系的 \(\beta\) ， \(\alpha\) ），影響了球面信號的分辨率，PRIN, LMVD, Compess等設(shè)置為24。
• radius_support: 支持半徑，定義了local patch的支持半徑，也就是說它確定了從關(guān)鍵點向外延伸的范圍內(nèi)哪些點將被納入local patch。
• num_radial_division: 表示徑向（從關(guān)鍵點向外輻射的方向）上的分割數(shù)目。它影響了在徑向方向上球面信號的分辨率。
• num_points: 采樣點的數(shù)量，這個值與local patch的固定點數(shù)一致（即1024點），確保輸入到轉(zhuǎn)換過程中的點數(shù)是一致的，這對于后續(xù)處理和模型輸入非常重要。
• random_sampling: 控制是否在從點云中選擇點時進(jìn)行隨機(jī)采樣，設(shè)置為 True 使得在局部區(qū)域內(nèi)的點采樣更加多樣化，避免由于局部密度過高或過低而導(dǎo)致的信息丟失。隨機(jī)采樣可以讓網(wǎng)絡(luò)更具魯棒性，適應(yīng)不同點云的分布。

spherical_voxel_optimized

def spherical_voxel_optimized(points: np.ndarray, size_bandwidth: int, size_radial_divisions: int,
                              radius_support: float, do_random_sampling: bool, num_random_points: int) \
        -> Tuple[np.ndarray, np.ndarray]:
    """Compute spherical voxel using the C++ code.

    Compute Spherical Voxel signal as defined in:
    Pointwise Rotation-Invariant Network withAdaptive Sampling and 3D Spherical Voxel Convolution.
    Yang You, Yujing Lou, Qi Liu, Yu-Wing Tai, Weiming Wang, Lizhuang Ma and Cewu Lu.
    AAAI 2020.

    :param points: the points to convert.
    :param size_bandwidth: alpha and beta bandwidth.
    :param size_radial_divisions: the number of bins along radial dimension.
    :param radius_support: the radius used to compute the points in the support.
    :param do_random_sampling: if true a subset of random points will be used to compute the spherical voxel.
    :param num_random_points: the number of points to keep if do_random_sampling is true.

    :return: A tuple containing:
        The spherical voxel, shape(size_radial_divisions, 2 * size_bandwidth, 2 * size_bandwidth).
        The points used to compute the signal normalized according the the farthest point.
    """
    if do_random_sampling:
        min_limit = 1 if points.shape[0] > 1 else 0
        indices_random = np.random.randint(min_limit, points.shape[0], num_random_points)
        points = points[indices_random]

    pts_norm = np.linalg.norm(points, axis=1)
    # Scale points to fit unit sphere
    pts_normed = points / pts_norm[:, None]
    pts_normed = np.clip(pts_normed, -1, 1)

    pts_s2_coord = S2.change_coordinates(pts_normed, p_from='C', p_to='S')
    # Convert to spherical voxel indices
    pts_s2_coord[:, 0] *= 2 * size_bandwidth / np.pi  # [0, pi]
    pts_s2_coord[:, 1] *= size_bandwidth / np.pi # raw 2*size_bandwidth/2*np.pi
    pts_s2_coord[:, 1][pts_s2_coord[:, 1] < 0] += 2 * size_bandwidth

    # Adaptive sampling factor sin{pi*[(1/2,..., 2*size_bandwidth+1/2)/(2*size_bandwidth)]}
    # 能更好的聚合點云信息，但是也會導(dǎo)致更多的形變，有得必有失
    daas_weights = np.sin(np.pi * (2 * np.arange(2 * size_bandwidth) + 1) / 4 / size_bandwidth).astype(np.float32)
    voxel = np.asarray(sv.compute(pts_on_s2=pts_s2_coord,
                                  pts_norm=pts_norm,
                                  size_bandwidth=size_bandwidth,
                                  size_radial_divisions=size_radial_divisions,
                                  radius_support=radius_support,
                                  daas_weights=daas_weights))
    pts_normed = points / np.max(pts_norm)
    return voxel.astype(np.float32), pts_normed.astype(np.float32)

• pts_norm 是local patch的點云徑向距離，所以 local patch輸入的時候最好經(jīng)過對于關(guān)鍵點的中心化操作 ，不然徑向距離會是關(guān)于坐標(biāo)系原點的。
• S2.change_coordinates 用于將點云從笛卡爾坐標(biāo)系轉(zhuǎn)換成球面坐標(biāo)系，球面坐標(biāo)系解釋見WIKI，簡單來說就是兩個坐標(biāo)，維度角度坐標(biāo)\beta，和經(jīng)度角度坐標(biāo)\alpha
• daas_weights是自適應(yīng)權(quán)重：
  • 采樣密度平衡 ：在球面坐標(biāo)系中，由于緯度（通常用 β 表示）不同區(qū)域的面積差異，不同區(qū)域的點密度會有所不同。例如，在球面的極地區(qū)域（緯度接近 0 或 π 的區(qū)域），同樣的角度變化可能覆蓋的球面面積較小，而在赤道區(qū)域，面積較大。為了避免在這些區(qū)域中出現(xiàn)過度或不足的采樣，自適應(yīng)采樣權(quán)重用于平衡不同緯度區(qū)域的影響。
  • 信息保持 ：通過在不同的緯度上使用不同的采樣權(quán)重，可以更精確地保留球面上重要的幾何特征，特別是在特定的關(guān)鍵區(qū)域。這樣可以確保球面信號在高緯度和低緯度區(qū)域都能有效地捕捉到有意義的幾何信息。
• sv.compute 用于體素轉(zhuǎn)換。

sv.compute

該函數(shù)是用pybind綁定的C++方法，文件為 spherical_voxel.cc ，代碼解釋如下：

初始化

    const float interval = radius_support / (size_radial_divisions);
    std::vector > > > > grids;
    std::vector > > features;
   
    grids.resize(size_radial_divisions);
    features.resize(size_radial_divisions);
  
    for (auto &beta: grids) {
	      beta.resize(2 * size_bandwidth);
	      for (auto &alpha: beta) {
	          alpha.resize(2 * size_bandwidth);
	      }
    }

    for (auto &beta: features) {
        beta.resize(2 * size_bandwidth);
        for (auto &alpha: beta) {
            alpha.resize(2 * size_bandwidth, 0);
        }
    }

• interval表示徑向分割下每個體素的徑向長度
• grids用來存儲每個體素覆蓋的所有點，可以通過下面的初始化看到，會初始化徑向，維度，經(jīng)度，每個體素是一種voxel
• feature用來存儲最終每個體素的特征（特征是密度特征）

grids填充

    // mapping the points to the voxel grid
    for (size_t i = 0; i < pts_on_s2.size(); i++) {
        int r_idx = int(pts_norm[i] / interval);
        // except for the points radius larger than radius_support
        if (r_idx > size_radial_divisions - 1) r_idx = size_radial_divisions - 1; 

        int beta_idx = int(pts_on_s2[i][0] + 0.5f);
        if (beta_idx > 2 * size_bandwidth - 1) beta_idx = 2 * size_bandwidth - 1;

        int alpha_idx = int(pts_on_s2[i][1] + 0.5f);
        if (alpha_idx > 2 * size_bandwidth - 1) alpha_idx = 2 * size_bandwidth - 1;

        grids[r_idx][beta_idx][alpha_idx].emplace_back(std::vector{pts_norm[i], pts_on_s2[i][0], pts_on_s2[i][1]});
    }

這里會遍歷每個點，計算每個點的徑向體素所用 r_idx ，緯度體素索引 beta_idx ，經(jīng)度體素索引 alpha_idx ，然后push到對應(yīng)的體素里面。

feature計算

首先計算每個體素的經(jīng)度左右特征計算邊界 left 、 right （也就是說每個體素的特征計算并不僅僅只考慮本體素內(nèi)部，還有一些可能出現(xiàn)的相鄰體素），這里計算左右邊界就用到自適應(yīng)權(quán)重，維度高的，左右邊界會寬一些。

之后根據(jù)左右邊界訪問對應(yīng)體素，并取出體素中所有點，基于徑向距離確定點是否靠近本體素中心，越靠近該點的特征權(quán)重越大（[0, 1]）。

然后考慮徑向相鄰體素內(nèi)部的點，用于本體素的特征計算，因為從徑向考慮，點分布相對連續(xù)，需要補(bǔ)充這樣的信息。

最后計算本體素的特征（密度特征（加過權(quán)的點個數(shù)））

    // compute the feature of each voxel
    for (size_t i = 0; i < size_radial_divisions; i++) {
        for (size_t j = 0; j < 2 * size_bandwidth; j++) {
            for (size_t k = 0; k < 2 * size_bandwidth; k++) {
                const float left = std::max(0.f, k - 0.5f / daas_weights[j]);
                const float right = std::min(2.f * size_bandwidth, k + 0.5f / daas_weights[j]);
                float sum = 0.f;
                int cnt = 0;

                for (int m = int(left + 0.5f); m < int(right + 0.5f); m++) {
                    for (int n = 0; n < grids[i][j][m].size(); n++) {
                        if (grids[i][j][m][n][2] > left && grids[i][j][m][n][2] < right) {
                            sum += 1.f - std::abs(grids[i][j][m][n][0] / interval - (i + 1)); // radial feature weight
                            cnt++;
                        }
                    }
                    
                    // 在實際情況中，點云數(shù)據(jù)可能分布在兩個相鄰的徑向分割之間，
                    // 尤其是當(dāng)點的徑向距離位于兩個徑向分割的邊界附近時。
                    // 為了防止因單純考慮當(dāng)前徑向分割而導(dǎo)致信息的丟失，
                    // 代碼會查找相鄰徑向分割中滿足條件的點，并將它們的貢獻(xiàn)也加到當(dāng)前體素單元的特征值中。
                    if (i < size_radial_divisions - 1) {
                        for (int n = 0; n < grids[i + 1][j][m].size(); n++) {
                            if (grids[i + 1][j][m][n][2] > left && grids[i + 1][j][m][n][2] < right) {
                                sum += 1.f - std::abs(grids[i + 1][j][m][n][0] / interval - (i + 1));
                                cnt++;
                            }
                        }
                    }
                }

                // 與徑向分割不同，緯度分割（即 beta 方向）代表的是球面坐標(biāo)中的角度，
                // 分割的區(qū)域代表不同的“環(huán)”或“帶”。
                // 在這種情況下，每個緯度分割對應(yīng)的球面區(qū)域是明確的，
                // 且這些分割區(qū)域之間沒有交叉，因此點不會“跨越”到另一個緯度分割。

                if (cnt > 0) {
                    features[i][j][k] = sum / cnt;
                }
            }
        }
    }