Commit 30dab101 authored by David Flynn's avatar David Flynn
Browse files

refactor: replace PCC*3T typedefs with Vec3<T>

This is part of a series attempting to remove unhelpful typedefs.

The old typedefs (PCCVector3D, PCCPoint3D, PCCColor3B) are retained for
compatibility, but marked as deprecated.
parent 00bac43d
......@@ -291,7 +291,7 @@ AttributeDecoder::computeColorPredictionWeights(
int64_t minValue[3] = {0, 0, 0};
int64_t maxValue[3] = {0, 0, 0};
for (int i = 0; i < predictor.neighborCount; ++i) {
const PCCColor3B colorNeighbor =
const Vec3<uint8_t> colorNeighbor =
pointCloud.getColor(indexes[predictor.neighbors[i].predictorIndex]);
for (size_t k = 0; k < 3; ++k) {
if (i == 0 || colorNeighbor[k] < minValue[k]) {
......@@ -354,8 +354,8 @@ AttributeDecoder::decodeColorsPred(
aps, pointCloud, indexesLOD, predictor, decoder);
decoder.decode(values);
const uint32_t pointIndex = indexesLOD[predictorIndex];
PCCColor3B& color = pointCloud.getColor(pointIndex);
const PCCColor3B predictedColor =
Vec3<uint8_t>& color = pointCloud.getColor(pointIndex);
const Vec3<uint8_t> predictedColor =
predictor.predictColor(pointCloud, indexesLOD);
const int64_t quantPredAttValue = predictedColor[0];
const int64_t delta = PCCInverseQuantization(UIntToInt(values[0]), qs);
......@@ -480,7 +480,7 @@ AttributeDecoder::decodeColorsRaht(
const int r = attributes[attribCount * n].round();
const int g = attributes[attribCount * n + 1].round();
const int b = attributes[attribCount * n + 2].round();
PCCColor3B color;
Vec3<uint8_t> color;
color[0] = uint8_t(PCCClip(r, 0, clipMax));
color[1] = uint8_t(PCCClip(g, 0, clipMax));
color[2] = uint8_t(PCCClip(b, 0, clipMax));
......@@ -528,7 +528,7 @@ AttributeDecoder::decodeColorsLift(
std::vector<double> weights;
PCCComputeQuantizationWeights(predictors, weights);
const size_t lodCount = numberOfPointsPerLOD.size();
std::vector<PCCVector3D> colors;
std::vector<Vec3<double>> colors;
colors.resize(pointCount);
// decompress
for (size_t predictorIndex = 0; predictorIndex < pointCount;
......@@ -559,7 +559,7 @@ AttributeDecoder::decodeColorsLift(
const double clipMax = (1 << desc.attr_bitdepth) - 1;
for (size_t f = 0; f < pointCount; ++f) {
PCCColor3B color;
Vec3<uint8_t> color;
for (size_t d = 0; d < 3; ++d) {
color[d] = uint8_t(PCCClip(std::round(colors[f][d]), 0.0, clipMax));
}
......
......@@ -481,8 +481,8 @@ AttributeEncoder::encodeReflectancesPred(
Vec3<int64_t>
AttributeEncoder::computeColorResiduals(
const PCCColor3B color,
const PCCColor3B predictedColor,
const Vec3<uint8_t> color,
const Vec3<uint8_t> predictedColor,
const int64_t qs,
const int64_t qs2)
{
......@@ -518,7 +518,7 @@ AttributeEncoder::computeColorPredictionWeights(
int64_t minValue[3] = {0, 0, 0};
int64_t maxValue[3] = {0, 0, 0};
for (int i = 0; i < predictor.neighborCount; ++i) {
const PCCColor3B colorNeighbor =
const Vec3<uint8_t> colorNeighbor =
pointCloud.getColor(indexesLOD[predictor.neighbors[i].predictorIndex]);
for (size_t k = 0; k < 3; ++k) {
if (i == 0 || colorNeighbor[k] < minValue[k]) {
......@@ -535,11 +535,12 @@ AttributeEncoder::computeColorPredictionWeights(
if (maxDiff > aps.adaptive_prediction_threshold) {
const int qs = aps.quant_step_size_luma;
const int qs2 = aps.quant_step_size_chroma;
PCCColor3B attrValue = pointCloud.getColor(indexesLOD[predictorIndex]);
Vec3<uint8_t> attrValue =
pointCloud.getColor(indexesLOD[predictorIndex]);
// base case: weighted average of n neighbours
predictor.predMode = 0;
PCCColor3B attrPred = predictor.predictColor(pointCloud, indexesLOD);
Vec3<uint8_t> attrPred = predictor.predictColor(pointCloud, indexesLOD);
Vec3<int64_t> attrResidualQuant =
computeColorResiduals(attrValue, attrPred, qs, qs2);
......@@ -617,8 +618,8 @@ AttributeEncoder::encodeColorsPred(
aps, pointCloud, indexesLOD, predictorIndex, predictor, encoder,
context);
const auto pointIndex = indexesLOD[predictorIndex];
const PCCColor3B color = pointCloud.getColor(pointIndex);
const PCCColor3B predictedColor =
const Vec3<uint8_t> color = pointCloud.getColor(pointIndex);
const Vec3<uint8_t> predictedColor =
predictor.predictColor(pointCloud, indexesLOD);
const int64_t quantAttValue = color[0];
const int64_t quantPredAttValue = predictedColor[0];
......@@ -628,7 +629,7 @@ AttributeEncoder::encodeColorsPred(
const int64_t reconstructedQuantAttValue =
quantPredAttValue + reconstructedDelta;
values[0] = uint32_t(IntToUInt(long(delta)));
PCCColor3B reconstructedColor;
Vec3<uint8_t> reconstructedColor;
reconstructedColor[0] =
uint8_t(PCCClip(reconstructedQuantAttValue, int64_t(0), clipMax));
for (size_t k = 1; k < 3; ++k) {
......@@ -791,7 +792,7 @@ AttributeEncoder::encodeColorsTransformRaht(
const int r = attributes[attribCount * n].round();
const int g = attributes[attribCount * n + 1].round();
const int b = attributes[attribCount * n + 2].round();
PCCColor3B color;
Vec3<uint8_t> color;
color[0] = uint8_t(PCCClip(r, 0, clipMax));
color[1] = uint8_t(PCCClip(g, 0, clipMax));
color[2] = uint8_t(PCCClip(b, 0, clipMax));
......@@ -839,7 +840,7 @@ AttributeEncoder::encodeColorsLift(
std::vector<double> weights;
PCCComputeQuantizationWeights(predictors, weights);
const size_t lodCount = numberOfPointsPerLOD.size();
std::vector<PCCVector3D> colors;
std::vector<Vec3<double>> colors;
colors.resize(pointCount);
for (size_t index = 0; index < pointCount; ++index) {
......@@ -892,7 +893,7 @@ AttributeEncoder::encodeColorsLift(
const double clipMax = (1 << desc.attr_bitdepth) - 1;
for (size_t f = 0; f < pointCount; ++f) {
PCCColor3B color;
Vec3<uint8_t> color;
for (size_t d = 0; d < 3; ++d) {
color[d] = uint8_t(PCCClip(std::round(colors[f][d]), 0.0, clipMax));
}
......
......@@ -99,8 +99,8 @@ protected:
PCCResidualsEncoder& encoder);
static Vec3<int64_t> computeColorResiduals(
const PCCColor3B color,
const PCCColor3B predictedColor,
const Vec3<uint8_t> color,
const Vec3<uint8_t> predictedColor,
const int64_t qs,
const int64_t qs2);
......
......@@ -70,13 +70,13 @@ struct PCCRangeResult {
};
struct PCCNNQuery3 {
PCCPoint3D point;
Vec3<double> point;
double radius;
size_t nearestNeighborCount;
};
struct PCCRangeQuery3 {
PCCPoint3D point;
Vec3<double> point;
double radius;
size_t maxResultCount;
};
......@@ -151,7 +151,7 @@ private:
class PCCIncrementalKdTree3 {
struct PCCIncrementalKdTree3Node {
PCCPoint3D pos;
Vec3<double> pos;
uint32_t id;
uint32_t left;
uint32_t right;
......@@ -178,7 +178,7 @@ public:
root = PCC_UNDEFINED_INDEX;
}
void reserve(const size_t pointCount) { nodes.reserve(pointCount); }
uint32_t insert(const PCCPoint3D point)
uint32_t insert(const Vec3<double> point)
{
const uint32_t id = static_cast<uint32_t>(nodes.size());
nodes.resize(id + 1);
......@@ -250,7 +250,7 @@ public:
}
uint32_t
hasNeighborWithinRange(const PCCVector3D& point, const double radius2) const
hasNeighborWithinRange(const Vec3<double>& point, const double radius2) const
{
return hasNeighborWithinRange(root, point, radius2);
}
......@@ -281,10 +281,10 @@ private:
}
PCCAxis3 computeSplitAxis(const uint32_t start, const uint32_t end) const
{
PCCPoint3D minBB = nodes[start].pos;
PCCPoint3D maxBB = nodes[start].pos;
Vec3<double> minBB = nodes[start].pos;
Vec3<double> maxBB = nodes[start].pos;
for (size_t i = start + 1; i < end; ++i) {
const PCCPoint3D& pt = nodes[i].pos;
const Vec3<double>& pt = nodes[i].pos;
for (int32_t k = 0; k < 3; ++k) {
if (minBB[k] > pt[k]) {
minBB[k] = pt[k];
......@@ -293,7 +293,7 @@ private:
}
}
}
PCCPoint3D d = maxBB - minBB;
Vec3<double> d = maxBB - minBB;
if (d.x() > d.y() && d.x() > d.z()) {
return PCC_AXIS3_X;
} else if (d.y() > d.z()) {
......@@ -434,7 +434,7 @@ private:
uint32_t hasNeighborWithinRange(
const uint32_t current,
const PCCVector3D& point,
const Vec3<double>& point,
const double radius2) const
{
if (current == PCC_UNDEFINED_INDEX) {
......@@ -475,7 +475,7 @@ private:
class PCCKdTree3 {
struct PCCKdTree3Node {
PCCBox3D BB;
PCCPoint3D centd;
Vec3<double> centd;
uint32_t id;
uint32_t start;
uint32_t end;
......@@ -484,7 +484,7 @@ class PCCKdTree3 {
uint32_t medianIdx;
};
struct PointIDNode {
PCCPoint3D pos;
Vec3<double> pos;
uint32_t id;
bool isVisisted;
};
......@@ -518,7 +518,7 @@ public:
? BB.max[nodes[parentNodeIdx].axis] = nodes[parentNodeIdx].median
: BB.min[nodes[parentNodeIdx].axis] = nodes[parentNodeIdx].median;
PCCPoint3D nodeMean = computePCCMean(start, end);
Vec3<double> nodeMean = computePCCMean(start, end);
PCCAxis3 axis = computeSplitAxisVar(start, end, nodeMean);
uint32_t medianIdx = findMedian(start, end, axis);
......@@ -548,7 +548,7 @@ public:
}
PCCBox3D BB = computeBoundingBox(0, pointCount);
PCCPoint3D nodeMean = computePCCMean(0, pointCount);
Vec3<double> nodeMean = computePCCMean(0, pointCount);
PCCAxis3 axis = computeSplitAxisVar(0, pointCount, nodeMean);
uint32_t medianIdx = findMedian(0, pointCount, axis);
......@@ -564,7 +564,7 @@ public:
rootNode.medianIdx = medianIdx;
}
uint32_t searchClosestAvailablePoint(PCCPoint3D queryPoint)
uint32_t searchClosestAvailablePoint(Vec3<double> queryPoint)
{
uint32_t idToClosestPoint = -1;
uint32_t id = 0;
......@@ -579,7 +579,7 @@ public:
uint32_t closestID = 0;
for (size_t i = start; i <= end; ++i) {
if (!pointCloudTemp[i].isVisisted) {
PCCPoint3D diff = pointCloudTemp[i].pos - queryPoint;
Vec3<double> diff = pointCloudTemp[i].pos - queryPoint;
uint32_t dist =
std::sqrt(diff[0] * diff[0] + diff[1] * diff[1] + diff[2] * diff[2]);
if (dist <= closestDistThr) {
......@@ -601,10 +601,10 @@ public:
private:
PCCBox3D computeBoundingBox(const uint32_t start, const uint32_t end) const
{
PCCPoint3D minBB = pointCloudTemp[start].pos;
PCCPoint3D maxBB = pointCloudTemp[start].pos;
Vec3<double> minBB = pointCloudTemp[start].pos;
Vec3<double> maxBB = pointCloudTemp[start].pos;
for (size_t i = start + 1; i < end; ++i) {
const PCCPoint3D& pt = pointCloudTemp[i].pos;
const Vec3<double>& pt = pointCloudTemp[i].pos;
for (int32_t k = 0; k < 3; ++k) {
if (minBB[k] > pt[k]) {
minBB[k] = pt[k];
......@@ -623,7 +623,7 @@ private:
PCCAxis3 computeSplitAxis(const uint32_t start, const uint32_t end) const
{
PCCBox3D BB = computeBoundingBox(start, end);
PCCPoint3D d = BB.max - BB.min;
Vec3<double> d = BB.max - BB.min;
if (d.x() > d.y() && d.x() > d.z()) {
return PCC_AXIS3_X;
} else if (d.y() > d.z()) {
......@@ -633,7 +633,7 @@ private:
}
}
PCCAxis3 computeSplitAxisVar(
const uint32_t start, const uint32_t end, PCCPoint3D nodeMean) const
const uint32_t start, const uint32_t end, Vec3<double> nodeMean) const
{
double nodeVar[3] = {0, 0, 0};
for (size_t axis = 0; axis < 3; ++axis) {
......@@ -652,10 +652,10 @@ private:
return PCC_AXIS3_Z;
}
}
PCCPoint3D computePCCMean(const uint32_t start, const uint32_t end)
Vec3<double> computePCCMean(const uint32_t start, const uint32_t end)
{
assert(end >= start);
PCCPoint3D nodeMean;
Vec3<double> nodeMean;
for (size_t axis = 0; axis < 3; ++axis) {
uint32_t acc = 0;
for (size_t i = start; i < end; i++) {
......
......@@ -331,10 +331,17 @@ struct PCCBox3 {
}
};
typedef Vec3<double> PCCVector3D;
typedef Vec3<double> PCCPoint3D;
//---------------------------------------------------------------------------
typedef DEPRECATED_MSVC Vec3<double> PCCVector3D DEPRECATED;
typedef DEPRECATED_MSVC Vec3<double> PCCPoint3D DEPRECATED;
typedef PCCBox3<double> PCCBox3D;
typedef Vec3<uint8_t> PCCColor3B;
typedef DEPRECATED_MSVC Vec3<uint8_t> PCCColor3B DEPRECATED;
template<typename T>
using PCCVector3 DEPRECATED = Vec3<T>;
//---------------------------------------------------------------------------
template<typename T>
T
......
......@@ -43,6 +43,14 @@
#include <utility>
#include <vector>
#if _MSC_VER
# define DEPRECATED_MSVC __declspec(deprecated)
# define DEPRECATED
#else
# define DEPRECATED_MSVC
# define DEPRECATED __attribute__((deprecated))
#endif
namespace pcc {
const uint32_t PCC_UNDEFINED_INDEX = -1;
......
......@@ -53,7 +53,7 @@
namespace pcc {
class PCCPointSet3 {
public:
typedef PCCVector3D PointType;
typedef Vec3<double> PointType;
//=========================================================================
// proxy object for use with iterator, allowing handling of PCCPointSet3's
......@@ -249,27 +249,27 @@ public:
swap(withFrameIndex, other.withFrameIndex);
}
PCCPoint3D operator[](const size_t index) const
Vec3<double> operator[](const size_t index) const
{
assert(index < positions.size());
return positions[index];
}
PCCPoint3D& operator[](const size_t index)
Vec3<double>& operator[](const size_t index)
{
assert(index < positions.size());
return positions[index];
}
PCCColor3B getColor(const size_t index) const
Vec3<uint8_t> getColor(const size_t index) const
{
assert(index < colors.size() && withColors);
return colors[index];
}
PCCColor3B& getColor(const size_t index)
Vec3<uint8_t>& getColor(const size_t index)
{
assert(index < colors.size() && withColors);
return colors[index];
}
void setColor(const size_t index, const PCCColor3B color)
void setColor(const size_t index, const Vec3<uint8_t> color)
{
assert(index < colors.size() && withColors);
colors[index] = color;
......@@ -434,7 +434,7 @@ public:
std::numeric_limits<double>::lowest()};
const size_t pointCount = getPointCount();
for (size_t i = 0; i < pointCount; ++i) {
const PCCPoint3D& pt = (*this)[i];
const Vec3<double>& pt = (*this)[i];
for (int k = 0; k < 3; ++k) {
if (pt[k] > bbox.max[k]) {
bbox.max[k] = pt[k];
......@@ -526,10 +526,10 @@ public:
// fout << std::setprecision(std::numeric_limits<double>::max_digits10);
fout << std::fixed << std::setprecision(5);
for (size_t i = 0; i < pointCount; ++i) {
const PCCPoint3D& position = (*this)[i];
const Vec3<double>& position = (*this)[i];
fout << position.x() << " " << position.y() << " " << position.z();
if (hasColors()) {
const PCCColor3B& color = getColor(i);
const Vec3<uint8_t>& color = getColor(i);
fout << " " << static_cast<int>(color[0]) << " "
<< static_cast<int>(color[1]) << " "
<< static_cast<int>(color[2]);
......@@ -547,11 +547,11 @@ public:
fout.close();
fout.open(fileName, std::ofstream::binary | std::ofstream::app);
for (size_t i = 0; i < pointCount; ++i) {
const PCCPoint3D& position = (*this)[i];
const Vec3<double>& position = (*this)[i];
fout.write(
reinterpret_cast<const char* const>(&position), sizeof(double) * 3);
if (hasColors()) {
const PCCColor3B& color = getColor(i);
const Vec3<uint8_t>& color = getColor(i);
fout.write(
reinterpret_cast<const char*>(&color), sizeof(uint8_t) * 3);
}
......@@ -895,8 +895,8 @@ public:
}
private:
std::vector<PCCPoint3D> positions;
std::vector<PCCColor3B> colors;
std::vector<Vec3<double>> positions;
std::vector<Vec3<uint8_t>> colors;
std::vector<uint16_t> reflectances;
std::vector<uint8_t> frameidx;
bool withColors;
......
......@@ -68,7 +68,7 @@ quantizePositionsUniq(
std::set<Vec3<int32_t>> uniquePoints;
int numSrcPoints = src.getPointCount();
for (int i = 0; i < numSrcPoints; ++i) {
const PCCVector3D& point = src[i];
const Vec3<double>& point = src[i];
Vec3<int32_t> quantizedPoint;
for (int k = 0; k < 3; k++) {
......@@ -127,7 +127,7 @@ quantizePositions(
};
for (int i = 0; i < numSrcPoints; ++i) {
const PCCVector3D point = src[i];
const Vec3<double> point = src[i];
auto& dstPoint = (*dst)[i];
for (int k = 0; k < 3; ++k) {
double k_pos = std::round((point[k] - offset[k]) * scaleFactor);
......@@ -184,7 +184,7 @@ inline bool
PCCTransfertColors(
const PCCPointSet3& source,
double sourceToTargetScaleFactor,
PCCVector3D targetToSourceOffset,
Vec3<double> targetToSourceOffset,
PCCPointSet3& target)
{
double targetToSourceScaleFactor = 1.0 / sourceToTargetScaleFactor;
......@@ -200,8 +200,8 @@ PCCTransfertColors(
3, source, 10);
target.addColors();
std::vector<PCCColor3B> refinedColors1;
std::vector<std::vector<PCCColor3B>> refinedColors2;
std::vector<Vec3<uint8_t>> refinedColors1;
std::vector<std::vector<Vec3<uint8_t>>> refinedColors2;
refinedColors1.resize(pointCountTarget);
refinedColors2.resize(pointCountTarget);
const size_t num_results = 1;
......@@ -212,7 +212,7 @@ PCCTransfertColors(
for (size_t index = 0; index < pointCountTarget; ++index) {
resultSet.init(&indices[0], &sqrDist[0]);
PCCVector3D posInSrc =
Vec3<double> posInSrc =
target[index] * targetToSourceScaleFactor + targetToSourceOffset;
kdtreeSource.index->findNeighbors(
......@@ -221,10 +221,10 @@ PCCTransfertColors(
}
for (size_t index = 0; index < pointCountSource; ++index) {
const PCCColor3B color = source.getColor(index);
const Vec3<uint8_t> color = source.getColor(index);
resultSet.init(&indices[0], &sqrDist[0]);
PCCVector3D posInTgt =
Vec3<double> posInTgt =
(source[index] - targetToSourceOffset) * sourceToTargetScaleFactor;
kdtreeTarget.index->findNeighbors(
......@@ -233,19 +233,19 @@ PCCTransfertColors(
}
for (size_t index = 0; index < pointCountTarget; ++index) {
const PCCColor3B color1 = refinedColors1[index];
const std::vector<PCCColor3B>& colors2 = refinedColors2[index];
const Vec3<uint8_t> color1 = refinedColors1[index];
const std::vector<Vec3<uint8_t>>& colors2 = refinedColors2[index];
if (colors2.empty()) {
target.setColor(index, color1);
} else {
PCCVector3D centroid2(0.0);
Vec3<double> centroid2(0.0);
for (const auto& color2 : colors2) {
for (size_t k = 0; k < 3; ++k) {
centroid2[k] += color2[k];
}
}
centroid2 /= colors2.size();
PCCColor3B color0;
Vec3<uint8_t> color0;
for (size_t k = 0; k < 3; ++k) {
color0[k] = uint8_t(PCCClip(round(centroid2[k]), 0.0, 255.0));
}
......@@ -274,7 +274,7 @@ inline bool
PCCTransfertReflectances(
const PCCPointSet3& source,
double sourceToTargetScaleFactor,
PCCVector3D targetToSourceOffset,
Vec3<double> targetToSourceOffset,
PCCPointSet3& target)
{
double targetToSourceScaleFactor = 1.0 / sourceToTargetScaleFactor;
......@@ -302,7 +302,7 @@ PCCTransfertReflectances(
for (size_t index = 0; index < pointCountTarget; ++index) {
resultSet.init(&indices[0], &sqrDist[0]);
PCCVector3D posInSrc =
Vec3<double> posInSrc =
target[index] * targetToSourceScaleFactor + targetToSourceOffset;
kdtreeSource.index->findNeighbors(
......@@ -314,7 +314,7 @@ PCCTransfertReflectances(
const uint16_t reflectance = source.getReflectance(index);
resultSet.init(&indices[0], &sqrDist[0]);
PCCVector3D posInTgt =
Vec3<double> posInTgt =
(source[index] - targetToSourceOffset) * sourceToTargetScaleFactor;
kdtreeTarget.index->findNeighbors(
......@@ -357,7 +357,7 @@ recolour(
Vec3<int> offset,
PCCPointSet3* target)
{
PCCVector3D combinedOffset;
Vec3<double> combinedOffset;
for (int k = 0; k < 3; k++)
combinedOffset[k] =
targetToSourceOffset[k] + double(offset[k]) / sourceToTargetScaleFactor;
......
......@@ -102,21 +102,21 @@ struct PCCPredictor {
PCCNeighborInfo neighbors[kAttributePredictionMaxNeighbourCount];
int8_t predMode;
PCCColor3B predictColor(
Vec3<uint8_t> predictColor(
const PCCPointSet3& pointCloud, const std::vector<uint32_t>& indexes) const
{
PCCVector3D predicted(0.0);
Vec3<double> predicted(0.0);
if (predMode > neighborCount) {
/* nop */
} else if (predMode > 0) {
const PCCColor3B color =
const Vec3<uint8_t> color =
pointCloud.getColor(indexes[neighbors[predMode - 1].predictorIndex]);
for (size_t k = 0; k < 3; ++k) {
predicted[k] += color[k];
}
} else {
for (size_t i = 0; i < neighborCount; ++i) {
const PCCColor3B color =
const Vec3<uint8_t> color =
pointCloud.getColor(indexes[neighbors[i].predictorIndex]);
const double w = neighbors[i].weight;
for (size_t k = 0; k < 3; ++k) {
......@@ -124,7 +124,7 @@ struct PCCPredictor {
}
}
}
return PCCColor3B(
return Vec3<uint8_t>(
uint8_t(std::round(predicted[0])), uint8_t(std::round(predicted[1])),
uint8_t(std::round(predicted[2])));
}
......
......@@ -432,15 +432,15 @@ PCCTMC3Encoder3::quantization(const PCCPointSet3& inputPointCloud)
}
// Offset the point cloud to account for (preset) _sliceOrigin.
PCCVector3D sliceOriginD{double(_sliceOrigin[0]), double(_sliceOrigin[1]),
double(_sliceOrigin[2])};
Vec3<double> sliceOriginD{double(_sliceOrigin[0]), double(_sliceOrigin[1]),
double(_sliceOrigin[2])};
// The new maximum bounds of the offset cloud
Vec3<int> maxBound{0};
const size_t pointCount = pointCloud.getPointCount();
for (size_t i = 0; i < pointCount; ++i) {
const PCCVector3D point = (pointCloud[i] -= sliceOriginD);
const Vec3<double> point = (pointCloud[i] -= sl