ACloudViewer  3.9.4
A Modern Library for 3D Data Processing
ContinuousConv.h
Go to the documentation of this file.
1 // ----------------------------------------------------------------------------
2 // - CloudViewer: www.cloudViewer.org -
3 // ----------------------------------------------------------------------------
4 // Copyright (c) 2018-2024 www.cloudViewer.org
5 // SPDX-License-Identifier: MIT
6 // ----------------------------------------------------------------------------
7 
8 #pragma once
9 
10 #include <tbb/parallel_for.h>
11 
12 #include "ml/impl/continuous_conv/CoordinateTransformation.h"
13 
14 namespace cloudViewer {
15 namespace ml {
16 namespace impl {
17 
20 template <class TFeat,
21  class TOut,
22  class TReal,
23  class TIndex,
24  InterpolationMode INTERPOLATION,
25  CoordinateMapping MAPPING,
26  bool ALIGN_CORNERS,
27  bool INDIVIDUAL_EXTENT,
28  bool ISOTROPIC_EXTENT,
29  bool POINT_IMPORTANCE>
30 void _CConvComputeFeaturesCPU(TOut* out_features,
31  const std::vector<int>& filter_dims,
32  const TFeat* filter,
33  size_t num_out,
34  const TReal* out_positions,
35  size_t num_inp,
36  const TReal* inp_positions,
37  const TFeat* inp_features,
38  const TFeat* inp_importance,
39  size_t neighbors_index_size,
40  const TIndex* neighbors_index,
41  const TFeat* neighbors_importance,
42  const int64_t* neighbors_row_splits,
43  const TReal* extents,
44  const TReal* offsets,
45  bool normalize) {
46  const bool NEIGHBORS_IMPORTANCE = neighbors_importance != nullptr;
47  // const int VECSIZE = 32;
48 #define VECSIZE 32
49  typedef Eigen::Array<TReal, VECSIZE, 1> Vec_t;
50  typedef InterpolationVec<TReal, VECSIZE, INTERPOLATION> InterpolationVec_t;
51  InterpolationVec_t interpolation;
52 
53  const int in_channels = filter_dims[filter_dims.size() - 2];
54  const int out_channels = filter_dims[filter_dims.size() - 1];
55 
56  int spatial_filter_size = 1;
57  for (int i = 0; i < 3; ++i) spatial_filter_size *= filter_dims[i];
58  Eigen::Array<int, 3, 1> filter_size_xyz(filter_dims[2], filter_dims[1],
59  filter_dims[0]);
60 
61  memset(out_features, 0, sizeof(TOut) * num_out * out_channels);
62 
63  typedef Eigen::Array<TFeat, VECSIZE, Eigen::Dynamic> Matrix;
64  typedef Eigen::Array<TReal, VECSIZE, 3> Matrix3C;
65 
66  tbb::parallel_for(
67  tbb::blocked_range<size_t>(0, num_out, 32),
68  [&](const tbb::blocked_range<size_t>& r) {
69  int range_length = r.end() - r.begin();
70 
71  Eigen::Matrix<TOut, Eigen::Dynamic, 1> normalizers(range_length,
72  1);
73  normalizers.setZero();
74 
75  Eigen::Matrix<TFeat, Eigen::Dynamic, Eigen::Dynamic> B(
76  in_channels * spatial_filter_size, range_length);
77  B.setZero();
78 
79  Matrix infeat(VECSIZE, in_channels);
80 
81  Eigen::Array<TReal, 3, 1> offsets_(offsets[0], offsets[1],
82  offsets[2]);
83 
84  Matrix3C inv_extents;
85  if (INDIVIDUAL_EXTENT == false) {
86  if (ISOTROPIC_EXTENT) {
87  inv_extents = 1 / extents[0];
88  } else {
89  inv_extents.col(0) = 1 / extents[0];
90  inv_extents.col(1) = 1 / extents[1];
91  inv_extents.col(2) = 1 / extents[2];
92  }
93  }
94 
95  for (size_t out_idx = r.begin(); out_idx != r.end();
96  ++out_idx) {
97  const int out_col = out_idx - r.begin();
98  const size_t neighbor_start = neighbors_row_splits[out_idx];
99  const size_t neighbor_end =
100  neighbors_row_splits[out_idx + 1];
101 
102  if (INDIVIDUAL_EXTENT) {
103  if (ISOTROPIC_EXTENT) {
104  inv_extents = 1 / extents[out_idx];
105  } else {
106  inv_extents.col(0) = 1 / extents[3 * out_idx + 0];
107  inv_extents.col(1) = 1 / extents[3 * out_idx + 1];
108  inv_extents.col(2) = 1 / extents[3 * out_idx + 2];
109  }
110  }
111 
112  typename InterpolationVec_t::Weight_t interp_weights;
113  typename InterpolationVec_t::Idx_t interp_indices;
114 
115  int vec_valid_count = 0;
116  Vec_t x, y, z;
117 
118  // set to zero to avoid problems with vectors with less than
119  // VECSIZE valid entries
120  x.setZero();
121  y.setZero();
122  z.setZero();
123  for (size_t n = neighbor_start; n < neighbor_end; ++n) {
124  const size_t inp_idx = neighbors_index[n];
125  const int i = vec_valid_count;
126  x(i) = inp_positions[inp_idx * 3 + 0] -
127  out_positions[out_idx * 3 + 0];
128  y(i) = inp_positions[inp_idx * 3 + 1] -
129  out_positions[out_idx * 3 + 1];
130  z(i) = inp_positions[inp_idx * 3 + 2] -
131  out_positions[out_idx * 3 + 2];
132 
133  const TFeat n_importance =
134  (NEIGHBORS_IMPORTANCE ? neighbors_importance[n]
135  : TFeat(1));
136  normalizers(out_col) += TOut(n_importance);
137 
138  for (int ic = 0; ic < in_channels; ++ic)
139  infeat(i, ic) =
140  inp_features[inp_idx * in_channels + ic];
141 
142  TFeat importance(1.0);
143  if (POINT_IMPORTANCE)
144  importance = inp_importance[inp_idx];
145  if (NEIGHBORS_IMPORTANCE) importance *= n_importance;
146 
147  if (POINT_IMPORTANCE || NEIGHBORS_IMPORTANCE) {
148  for (int ic = 0; ic < in_channels; ++ic)
149  infeat(i, ic) *= importance;
150  }
151 
152  ++vec_valid_count;
153  if (vec_valid_count == VECSIZE) {
154  ComputeFilterCoordinates<ALIGN_CORNERS, MAPPING,
155  TReal, VECSIZE>(
156  x, y, z, filter_size_xyz, inv_extents,
157  offsets_);
158  interpolation.Interpolate(
159  interp_weights, interp_indices, x, y, z,
160  filter_size_xyz, in_channels);
161  for (int k = 0; k < VECSIZE; ++k)
162  for (int j = 0; j < InterpolationVec_t::Size();
163  ++j) {
164  for (int ic = 0; ic < in_channels; ++ic)
165  B(interp_indices(j, k) + ic, out_col) +=
166  TFeat(interp_weights(j, k)) *
167  infeat(k, ic);
168  }
169  vec_valid_count = 0;
170  }
171  }
172  if (vec_valid_count) {
173  ComputeFilterCoordinates<ALIGN_CORNERS, MAPPING, TReal,
174  VECSIZE>(
175  x, y, z, filter_size_xyz, inv_extents,
176  offsets_);
177  interpolation.Interpolate(interp_weights,
178  interp_indices, x, y, z,
179  filter_size_xyz, in_channels);
180  for (int k = 0; k < vec_valid_count; ++k)
181  for (int j = 0; j < InterpolationVec_t::Size();
182  ++j) {
183  for (int ic = 0; ic < in_channels; ++ic)
184  B(interp_indices(j, k) + ic, out_col) +=
185  TFeat(interp_weights(j, k)) *
186  infeat(k, ic);
187  }
188  }
189 
190  } // out_idx
191 
192  Eigen::Map<const Eigen::Matrix<TFeat, Eigen::Dynamic,
193  Eigen::Dynamic>>
194  A(filter, out_channels,
195  spatial_filter_size * in_channels);
196  Eigen::Map<Eigen::Matrix<TOut, Eigen::Dynamic, Eigen::Dynamic>>
197  C(out_features + (r.begin() * out_channels),
198  out_channels, range_length);
199 
200  C = (A * B).template cast<TOut>();
201  if (normalize) {
202  for (int i = 0; i < range_length; ++i) {
203  if (normalizers(i) != TOut(0))
204  C.col(i) /= normalizers(i);
205  }
206  }
207  });
208 #undef VECSIZE
209 }
210 
285 template <class TFeat, class TOut, class TReal, class TIndex>
286 void CConvComputeFeaturesCPU(TOut* out_features,
287  const std::vector<int>& filter_dims,
288  const TFeat* filter,
289  size_t num_out,
290  const TReal* out_positions,
291  size_t num_inp,
292  const TReal* inp_positions,
293  const TFeat* inp_features,
294  const TFeat* inp_importance,
295  size_t neighbors_index_size,
296  const TIndex* neighbors_index,
297  const TFeat* neighbors_importance,
298  const int64_t* neighbors_row_splits,
299  const TReal* extents,
300  const TReal* offsets,
301  InterpolationMode interpolation,
302  CoordinateMapping coordinate_mapping,
303  bool align_corners,
304  bool individual_extent,
305  bool isotropic_extent,
306  bool normalize) {
307  // Dispatch all template parameter combinations
308  bool has_importance = inp_importance;
309 
310 #define FN_PARAMETERS \
311  out_features, filter_dims, filter, num_out, out_positions, num_inp, \
312  inp_positions, inp_features, inp_importance, neighbors_index_size, \
313  neighbors_index, neighbors_importance, neighbors_row_splits, \
314  extents, offsets, normalize
315 
316 #define CALL_TEMPLATE(INTERPOLATION, MAPPING, ALIGN_CORNERS, \
317  INDIVIDUAL_EXTENT, ISOTROPIC_EXTENT, HAS_IMPORTANCE) \
318  if (INTERPOLATION == interpolation && MAPPING == coordinate_mapping && \
319  ALIGN_CORNERS == align_corners && \
320  INDIVIDUAL_EXTENT == individual_extent && \
321  ISOTROPIC_EXTENT == isotropic_extent && \
322  HAS_IMPORTANCE == has_importance) \
323  _CConvComputeFeaturesCPU<TFeat, TOut, TReal, TIndex, INTERPOLATION, \
324  MAPPING, ALIGN_CORNERS, INDIVIDUAL_EXTENT, \
325  ISOTROPIC_EXTENT, HAS_IMPORTANCE>( \
326  FN_PARAMETERS);
327 
328 #define CALL_TEMPLATE2(INTERPOLATION, MAPPING) \
329  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, true, true, true) \
330  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, true, true, false) \
331  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, true, false, true) \
332  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, true, false, false) \
333  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, false, true, true) \
334  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, false, true, false) \
335  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, false, false, true) \
336  CALL_TEMPLATE(INTERPOLATION, MAPPING, true, false, false, false) \
337  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, true, true, true) \
338  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, true, true, false) \
339  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, true, false, true) \
340  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, true, false, false) \
341  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, false, true, true) \
342  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, false, true, false) \
343  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, false, false, true) \
344  CALL_TEMPLATE(INTERPOLATION, MAPPING, false, false, false, false)
345 
346 #define CALL_TEMPLATE3(INTERPOLATION) \
347  CALL_TEMPLATE2(INTERPOLATION, CoordinateMapping::BALL_TO_CUBE_RADIAL) \
348  CALL_TEMPLATE2(INTERPOLATION, \
349  CoordinateMapping::BALL_TO_CUBE_VOLUME_PRESERVING) \
350  CALL_TEMPLATE2(INTERPOLATION, CoordinateMapping::IDENTITY)
351 
352 #define CALL_TEMPLATE4 \
353  CALL_TEMPLATE3(InterpolationMode::LINEAR) \
354  CALL_TEMPLATE3(InterpolationMode::LINEAR_BORDER) \
355  CALL_TEMPLATE3(InterpolationMode::NEAREST_NEIGHBOR)
356 
357  CALL_TEMPLATE4
358 
359 #undef CALL_TEMPLATE
360 #undef CALL_TEMPLATE2
361 #undef CALL_TEMPLATE3
362 #undef CALL_TEMPLATE4
363 
364 #undef FN_PARAMETERS
365 }
366 
367 } // namespace impl
368 } // namespace ml
369 } // namespace cloudViewer
CALL_TEMPLATE4
#define CALL_TEMPLATE4
VECSIZE
#define VECSIZE
normalize
__host__ __device__ float2 normalize(float2 v)
Definition: cutil_math.h:1179
cloudViewer::ml::impl::CConvComputeFeaturesCPU
void CConvComputeFeaturesCPU(TOut *out_features, const std::vector< int > &filter_dims, const TFeat *filter, size_t num_out, const TReal *out_positions, size_t num_inp, const TReal *inp_positions, const TFeat *inp_features, const TFeat *inp_importance, size_t neighbors_index_size, const TIndex *neighbors_index, const TFeat *neighbors_importance, const int64_t *neighbors_row_splits, const TReal *extents, const TReal *offsets, InterpolationMode interpolation, CoordinateMapping coordinate_mapping, bool align_corners, bool individual_extent, bool isotropic_extent, bool normalize)
Definition: ContinuousConv.h:286
cloudViewer::ml::impl::_CConvComputeFeaturesCPU
void _CConvComputeFeaturesCPU(TOut *out_features, const std::vector< int > &filter_dims, const TFeat *filter, size_t num_out, const TReal *out_positions, size_t num_inp, const TReal *inp_positions, const TFeat *inp_features, const TFeat *inp_importance, size_t neighbors_index_size, const TIndex *neighbors_index, const TFeat *neighbors_importance, const int64_t *neighbors_row_splits, const TReal *extents, const TReal *offsets, bool normalize)
Definition: ContinuousConv.h:30
cloudViewer::ml::impl::ComputeFilterCoordinates
void ComputeFilterCoordinates(Eigen::Array< T, VECSIZE, 1 > &x, Eigen::Array< T, VECSIZE, 1 > &y, Eigen::Array< T, VECSIZE, 1 > &z, const Eigen::Array< int, 3, 1 > &filter_size, const Eigen::Array< T, VECSIZE, 3 > &inv_extents, const Eigen::Array< T, 3, 1 > &offset)
Definition: CoordinateTransformation.h:107
cloudViewer
Generic file read and write utility for python interface.
Definition: AutoSegmentationTools.h:16
cloudViewer::ml::impl::InterpolationVec
Class for computing interpolation weights.
Definition: CoordinateTransformation.h:185