[VoxelGrid] Downsampling produces wrong result

If my assumption that for the same input PointCloud, monotonically increasing voxel size should produce monotonically decreased downsampled PointCloud is correct then something wrong with the downsampling in VoxelGrid.

Let's take an example, the following code creates point cloud and voxelizes it with three different voxel sizes: 3.0f, 3.33f, 4.0f.

#include <pcl/test/gtest.h>
#include <pcl/point_types.h>
#include <pcl/filters/filter.h>
#include <pcl/filters/voxel_grid.h>

using namespace pcl;
using namespace pcl::io;
using namespace Eigen;

TEST (VoxelGrid, BugInFilter)
{
  const size_t size = 1000;
  float x = 0;
  float y = 360;
  float z = 720;
  PointCloud<PointXYZ> pcl {};
  for (size_t i = 0; i < size; ++i) {
    pcl.push_back(pcl::PointXYZ{x, y, z});
    x = x + 0.1f;
    y = y + 0.1f;
    z = z + 0.1f;
  }
  PointCloud<PointXYZ>::Ptr pclPtr(new PointCloud<PointXYZ>(pcl));

  VoxelGrid<PointXYZ> grid;

  std::array<float, 3> voxelSizes =  { 3.0f, 3.33f, 4.0f};
  for(const float& voxelSize: voxelSizes) {
    PointCloud<PointXYZ> output;
    grid.setLeafSize (voxelSize, voxelSize, voxelSize);
    grid.setInputCloud (pclPtr);
    grid.filter (output);

    PointXYZ& prev = output[0];
    double dx = 0;
    double dy = 0;
    double dz = 0;
    for (std::size_t i = 1; i < output.size(); ++i) {
      const PointXYZ& curr = output[i];
      dx += std::abs(curr.x - prev.x);
      dy += std::abs(curr.y - prev.y);
      dz += std::abs(curr.z - prev.z);
      prev = curr;
    }
    double avgDx = dx / output.size();
    double avgDy = dy / output.size();
    double avgDz = dz / output.size();
    std::cout << "voxelSize = " << voxelSize << ", downsampled size: " << output.size() << ", averate dx: " << avgDx
              << ", average dy: " << avgDy << ", average dz: " << avgDz<< std::endl;
  }
}

One would expect that voxelized(3.0f).size() > voxelized(3.33f).size() > voxelized(4.0f).size(), but that is not the case. (voxelized(3.0f).size() means use VoxelGrid to voxelize with leaf size 3.0 and get the size of output). It produces the following result:

voxelSize = 3, downsampled size: 67, averate dx: 1.46342, average dy: 1.46352, average dz: 1.46308
voxelSize = 3.33, downsampled size: 90, averate dx: 1.09388, average dy: 1.09396, average dz: 1.09362
voxelSize = 4, downsampled size: 49, averate dx: 1.96019, average dy: 1.96033, average dz: 1.95973

Is it expected that the average difference between x/y/z of voxelized points is around voxelSize / 2, is it a mistake to expect it to be around voxelSize?

Environment:

• OS: Ubuntu 20.04
• Compiler: c++ (Ubuntu 9.4.0-1ubuntu1~20.04.1) 9.4.0
• PCL Version: 7d52b107223ef8adec990608e8685392b48d79c6

Reproduce the reported counts using the included gtest example and inspect the VoxelGrid implementation reached through pcl/filters/voxel_grid.h. Determine whether monotonic output size is a valid invariant and whether the reported average spacing is expected; document the conclusion or add a focused regression test and correction if behavior is confirmed as a bug.