PotreePicker.cpp

#include "Services/Variables.h"
#include "Systems/Core/CameraSystem.h"
#include "Systems/Core/TransformSystem.h"
#include "Utilities/FrustumClassification.h"
#include "Context.h"
 
#include "PotreeSystem.h"
 
#include "Foundation/Logging/Logger.h"
 
namespace {
 
  Cogs::Logging::Log logger = Cogs::Logging::getLogger("PotreePicker");
 
  using namespace Cogs::Core;
 
  unsigned popcount(uint8_t x)
  {
    x = (x & 0x55u) + ((x >> 1) & 0x55u);
    x = (x & 0x33u) + ((x >> 2) & 0x33u);
    x = (x & 0x0Fu) + ((x >> 4) & 0x0Fu);
    return x;
  }
 
  glm::vec3 euclidean(const glm::vec4& h)
  {
    return (1.f / h.w) * glm::vec3(h);
  }
 
  bool isBBoxInFrustum(const glm::mat4& pickFromLocal,
                       const glm::mat4& localFromPick,
                       const glm::vec3& min,
                       const glm::vec3& max)
  {
    unsigned mask = 0;
    for (unsigned i = 0; i < 8; i++) {
      glm::vec4 p = pickFromLocal * glm::vec4((i & 1) ? min.x : max.x,
                                              (i & 2) ? min.y : max.y,
                                              (i & 4) ? min.z : max.z,
                                              1.f);
      if (-p.w <= p.x) mask |= 0b000001;
      if (-p.w <= p.y) mask |= 0b000010;
      if (-p.w <= p.z) mask |= 0b000100;
      if (+p.x <= p.w) mask |= 0b001000;
      if (+p.y <= p.w) mask |= 0b010000;
      if (+p.z <= p.w) mask |= 0b100000;
    }
    if (mask != 0b111111) return false;
 
    // Inverse test as well to get even tighter estimate
    mask = 0;
    for (unsigned i = 0; i < 8; i++) {
      glm::vec3 p = euclidean(localFromPick * glm::vec4((i & 1) ? -1.f : 1.f,
                                                        (i & 2) ? -1.f : 1.f,
                                                        (i & 4) ? -1.f : 1.f,
                                                        1.f));
      if (min.x <= p.x) mask |= 0b000001;
      if (min.y <= p.y) mask |= 0b000010;
      if (min.z <= p.z) mask |= 0b000100;
      if (p.x <= max.x) mask |= 0b001000;
      if (p.y <= max.y) mask |= 0b010000;
      if (p.z <= max.z) mask |= 0b100000;
    }
    return mask == 0b111111;
  }
 
  inline const glm::mat4 createPickMatrix(const CameraData& camData, const glm::vec2& pixelSize, const glm::vec2& query_pos_norm)
  {
    // Create a frustum around ray of size maxpointsize
    const float sx = camData.viewportSize.x / pixelSize.x;
    const float sy = camData.viewportSize.y / pixelSize.y;
    return glm::mat4(sx, 0, 0, 0,
                     0, sy, 0, 0,
                     0, 0, 1, 0,
                     -sx * query_pos_norm.x, -sy * query_pos_norm.y, 0, 1);
  }
 
  // Create a ortho matrix that is sized to match a viewspace radius plus a guard band
  inline const glm::mat4 createOrthoPickMatrix(const CameraData& camData, const float guard_norm, const float viewSpaceRadius, const glm::vec2& query_norm)
  {
    const glm::mat4 inverseRawProjectionMatrix = glm::inverse(camData.rawProjectionMatrix);
    glm::vec3 viewCenterA = euclidean(inverseRawProjectionMatrix * glm::vec4(query_norm, -1.f, 1.f));
    glm::vec3 viewCenterB = euclidean(inverseRawProjectionMatrix * glm::vec4(query_norm, 1.f, 1.f));
 
    // Calculate size of guard band at front and back planes
    glm::vec3 viewCenterA_guard = euclidean(inverseRawProjectionMatrix * glm::vec4(query_norm + guard_norm, -1.f, 1.f));
    glm::vec3 viewCenterB_guard = euclidean(inverseRawProjectionMatrix * glm::vec4(query_norm + guard_norm, 1.f, 1.f));
 
    const float g_n = viewCenterA_guard.x - viewCenterA.x;    // extra guard at near plane
    const float g_f = viewCenterB_guard.x - viewCenterB.x;    // additional guard at far plane
    const float r = viewSpaceRadius + g_n;
    const float s = g_f - g_n;
    const float d = -(viewCenterB.z - viewCenterA.z);
    //const float a = (2.f * r + s) / d; //unused param
 
 
    // Shear x and y so center line aligns with perspective direction
    float shearX = -(viewCenterB.x - viewCenterA.x) / (viewCenterB.z - viewCenterA.z);
    float shearY = -(viewCenterB.y - viewCenterA.y) / (viewCenterB.z - viewCenterA.z);
    glm::mat4 orthoPickFromViewFrame =
 
      // Orthogonal projection
      glm::mat4(1.f, 0.f, 0.f, 0.f,
                0.f, 1.f, 0.f, 0.f,
                0.f, 0.f, -(s + 2.f * r) / d, -s / d,
                0.f, 0.f, -r, r) *
 
      // Shear
      glm::mat4(1.f, 0.f, 0.f, 0.f,
                0.f, 1.f, 0.f, 0.f,
                shearX, shearY, 1.f, 0.f,
                0.f, 0.f, 0.f, 1.f) *
 
      // Translate so viewCenterA is at the origin
      glm::mat4(1.f, 0.f, 0.f, 0.f,
                0.f, 1.f, 0.f, 0.f,
                0.f, 0.f, 1.f, 0.f,
                -viewCenterA.x, -viewCenterA.y, -viewCenterA.z, 1.f);
 
    return orthoPickFromViewFrame;
  }
 
  bool rayPickFixedSize(PotreeData* poData, const CameraData& camData,
                        const glm::mat4& viewFromOcttreeFrame,
                        const glm::vec2& query_norm,
                        const float guard_px,
                        const float viewspaceRadius,
                        float& viewspaceDepthBest, glm::vec3& worldspacePosBest, bool roundPoint, bool depthPoint)
  {
    const glm::mat4& worldFromOcttreeFrame = poData->worldFromOcttreeFrame;
    bool anyHits = false;
 
    // Perspective matrix that limits point size to max size
    // this is in screen space, and thus we need a perspective matrix
    glm::mat4 perspPickFromOcttreeFrame = (createPickMatrix(camData, glm::vec2(poData->pointSize.val) + guard_px, query_norm) *
                                           camData.rawProjectionMatrix * viewFromOcttreeFrame);
    glm::mat4 octtreeFrameFromPerspPick = glm::inverse(perspPickFromOcttreeFrame);
 
    float viewspace_nearest = viewspaceDepthBest;
 
    if (isBBoxInFrustum(perspPickFromOcttreeFrame,
                        octtreeFrameFromPerspPick,
                        poData->octtreeFrame.tightBBoxMin,
                        poData->octtreeFrame.tightBBoxMax))
    {
 
      // Run through all active cells
      for (auto* cell : poData->activeCells) {
        if (isBBoxInFrustum(perspPickFromOcttreeFrame,
                            octtreeFrameFromPerspPick,
                            cell->tbmin,
                            cell->tbmax))
        {
          bool hit = false;
          size_t indexBest = 0;
          const glm::vec3* points = cell->points.data();
          const size_t points_N = cell->pointCount;
          for (size_t i = 0; i < points_N; i++) {
            const glm::vec3 p = points[i];
 
 
            // Check against max point size-frustum
            glm::vec4 maxSizePointCoord_h = (perspPickFromOcttreeFrame[0] * p.x +
                                             perspPickFromOcttreeFrame[1] * p.y +
                                             perspPickFromOcttreeFrame[2] * p.z +
                                             perspPickFromOcttreeFrame[3]);
            if ((maxSizePointCoord_h.x <= maxSizePointCoord_h.w) && (-maxSizePointCoord_h.w <= maxSizePointCoord_h.x) &&
                (maxSizePointCoord_h.y <= maxSizePointCoord_h.w) && (-maxSizePointCoord_h.w <= maxSizePointCoord_h.y))
            {
              glm::vec3 maxSizePointCoord = (1.f / maxSizePointCoord_h.w) * glm::vec3(maxSizePointCoord_h);
              float maxSizeRadiusSquared = maxSizePointCoord.x * maxSizePointCoord.x + maxSizePointCoord.y * maxSizePointCoord.y;
              if (roundPoint && (1.f < maxSizeRadiusSquared)) continue;
 
              // view space position of point
              glm::vec4 viewSpacePos_h = (viewFromOcttreeFrame[0] * p.x +
                                          viewFromOcttreeFrame[1] * p.y +
                                          viewFromOcttreeFrame[2] * p.z +
                                          viewFromOcttreeFrame[3]);
              float viewspace_depth = -viewSpacePos_h.z / viewSpacePos_h.w;      // flip sign due to looking along -Z.
 
              if (depthPoint) {
 
                float r = 0.f;
                if (roundPoint) {
                  r = std::sqrt(maxSizeRadiusSquared);
                }
                else {
                  r = std::max(std::abs(maxSizePointCoord.x),
                               std::abs(maxSizePointCoord.y));
                }
 
 
                viewspace_depth += viewspaceRadius * r;
              }
 
              if (viewspace_nearest < viewspace_depth) continue;
 
              viewspace_nearest = viewspace_depth;
              hit = true;
              indexBest = i;
            }
            if (hit) {
              worldspacePosBest = euclidean(worldFromOcttreeFrame * glm::vec4(cell->points[indexBest], 1.f));
              viewspaceDepthBest = viewspace_nearest;
              anyHits = true;
            }
          }
        }
      }
    }
    return anyHits;
  }
 
 
  float getDensityScale(const PotreeData* poData, const glm::vec3& p)
  {
    const glm::vec3 invScale = glm::vec3(1.f) / (poData->octtreeFrame.fullBBoxMax - poData->octtreeFrame.fullBBoxMin);
    const glm::vec3 invShift = 0.5f * (poData->octtreeFrame.tightBBoxMax - poData->octtreeFrame.tightBBoxMin);
 
    glm::vec3 unitPos = glm::clamp(invScale * (p + invShift), glm::vec3(0.f), glm::vec3(1.f));
 
    size_t offset = 0;
    float scale = 1.f;
    for (size_t l = 0; l < 30; l++) {
      glm::vec3 t = 2.f * unitPos;
      unitPos = fract(t);
 
      size_t childIndex = ((1.0 <= t.x ? 4 : 0) +
                           (1.0 <= t.y ? 2 : 0) +
                           (1.0 <= t.z ? 1 : 0));
      size_t childBit = 1 << childIndex;
 
      glm::u8vec4 data = poData->octtreeTextureData[offset];
      uint8_t childrenBitmask = data.r;
 
      if ((childrenBitmask & childBit) == 0) {
        // Scale value as the texture sampler returns it
        float factor = (1.f / 255.f) * data.w;
        // Unpack formula
        float decoded = poData->densitySized.bias - ((poData->densitySized.scale * factor) + poData->densitySized.levelScale * float(l));
        // Note: inverse exponent on purpose to replace division with multiplication further down
        scale = glm::exp2(-glm::min(0.f, decoded));
        break;
      }
      else {
        size_t childrenOffset = (static_cast<size_t>(data.b) << 8) + data.g;
        uint8_t lesserChildren = childrenBitmask & (childBit - 1u);
        size_t childrenBeforeSelectedChild = popcount(lesserChildren);
        offset = childrenOffset + childrenBeforeSelectedChild;
        assert(offset < poData->octtreeTextureData.size());
      }
    }
    return scale;
  }
 
  bool rayPickDistanceScaled(PotreeData* poData, const CameraData& camData,
                             const glm::mat4& viewFromOcttreeFrame,
                             const glm::vec2& query_norm,
                             const float guard_px,
                             const float guard_norm,
                             const float viewSpaceRadius,
                             float& viewspaceDepthBest, glm::vec3& worldspacePosBest, bool roundPoint, bool depthPoint, bool densityScale, bool hasViewport)
  {
    const glm::mat4& worldFromOcttreeFrame = poData->worldFromOcttreeFrame;
    const float modelSpaceRadius = poData->modelSpaceRadius;
 
    bool anyHits = false;
 
    // Create an orthographic projection aligned with the perspective direction
    // of the query position, sized by the view space radius of points
    glm::mat4 orthoPickFromOcttreeFrame;
    glm::mat4 octtreeFrameFromOrthoPick;
 
    // Find start and end point for ray in view space.
    glm::mat4 perspPickFromOcttreeFrame;  // only referenced if hasViewport is true
    glm::mat4 octtreeFrameFromPerspPick;
    if(hasViewport) {
      orthoPickFromOcttreeFrame = createOrthoPickMatrix(camData, guard_norm, viewSpaceRadius, query_norm) * viewFromOcttreeFrame;
      octtreeFrameFromOrthoPick = glm::inverse(orthoPickFromOcttreeFrame);
 
      // Perspective matrix that limits point size to max size
      // this is in screen space, and thus we need a perspective matrix
      perspPickFromOcttreeFrame =
        createPickMatrix(camData, glm::vec2(poData->pointSize.max) + guard_px, query_norm) *
        camData.rawProjectionMatrix * viewFromOcttreeFrame;
      octtreeFrameFromPerspPick = glm::inverse(perspPickFromOcttreeFrame);
 
    }
    else {
      orthoPickFromOcttreeFrame = camData.rawProjectionMatrix * viewFromOcttreeFrame;
      octtreeFrameFromOrthoPick = glm::inverse(orthoPickFromOcttreeFrame);
    }
 
    // Check if bounding box of full point cloud intersects pick frustum
    float viewspace_nearest = viewspaceDepthBest;
    if (!isBBoxInFrustum(orthoPickFromOcttreeFrame,
                         octtreeFrameFromOrthoPick,
                         poData->octtreeFrame.tightBBoxMin - modelSpaceRadius,
                         poData->octtreeFrame.tightBBoxMax + modelSpaceRadius))
    {
      return false;
    }
 
    // Run through all active cells
    for (auto* cell : poData->activeCells) {
 
      // Perspective early-out-test is important when we are close to or inside the point cloud
      // as this clamps the maximum point size, which will happen a lot in that case. This test
      // only works with a viewport.
      if ((!hasViewport || isBBoxInFrustum(perspPickFromOcttreeFrame,      // Max-size test
                                           octtreeFrameFromPerspPick,
                                           cell->tbmin,
                                           cell->tbmax)) &&
          isBBoxInFrustum(orthoPickFromOcttreeFrame,      // view-space-size test
                          octtreeFrameFromOrthoPick,
                          cell->tbmin - modelSpaceRadius,
                          cell->tbmax + modelSpaceRadius))
      {
        bool hit = false;
        size_t indexBest = 0;
        const glm::vec3* points = cell->points.data();
        const size_t points_N = cell->pointCount;
        for (size_t i = 0; i < points_N; i++) {
          const glm::vec3 p = points[i];
 
          // check against world-space point size frustum
          glm::vec4 maxSizePointCoord_h;
          if (hasViewport) {
            maxSizePointCoord_h = (perspPickFromOcttreeFrame[0] * p.x +
                                   perspPickFromOcttreeFrame[1] * p.y +
                                   perspPickFromOcttreeFrame[2] * p.z +
                                   perspPickFromOcttreeFrame[3]);
            if ((maxSizePointCoord_h.w < maxSizePointCoord_h.x) || (maxSizePointCoord_h.x < -maxSizePointCoord_h.w) ||
                (maxSizePointCoord_h.w < maxSizePointCoord_h.y) || (maxSizePointCoord_h.y < -maxSizePointCoord_h.w))
            {
              continue;
            }
          }
          glm::vec4 worldSizePointCoord_h = (orthoPickFromOcttreeFrame[0] * p.x +
                                             orthoPickFromOcttreeFrame[1] * p.y +
                                             orthoPickFromOcttreeFrame[2] * p.z +
                                             orthoPickFromOcttreeFrame[3]);
          if ((worldSizePointCoord_h.w < worldSizePointCoord_h.x) || (worldSizePointCoord_h.x < -worldSizePointCoord_h.w) ||
              (worldSizePointCoord_h.w < worldSizePointCoord_h.y) || (worldSizePointCoord_h.y < -worldSizePointCoord_h.w) ||
              (worldSizePointCoord_h.w < worldSizePointCoord_h.z) || (worldSizePointCoord_h.z < -worldSizePointCoord_h.w))
          {
            continue;
          }
 
          // Check against world-space round point test
          glm::vec3 worldSizePointCoord = (1.f / worldSizePointCoord_h.w) * glm::vec3(worldSizePointCoord_h);
          float worldSizeRadiusSquared = worldSizePointCoord.x * worldSizePointCoord.x + worldSizePointCoord.y * worldSizePointCoord.y;
          if (roundPoint && (1.f < worldSizeRadiusSquared)) {
            continue;
          }
 
          // Check against max point size-frustum, only done when we have a viewport that can define a max size
          float maxSizeRadiusSquared = worldSizeRadiusSquared;
          glm::vec3 maxSizePointCoord = worldSizePointCoord;
          if (hasViewport) {
 
            // Check against max point-size round point test
            if (roundPoint) {
              maxSizePointCoord = (1.f / maxSizePointCoord_h.w) * glm::vec3(maxSizePointCoord_h);
              maxSizeRadiusSquared = maxSizePointCoord.x * maxSizePointCoord.x + maxSizePointCoord.y * maxSizePointCoord.y;
              if (1.f < maxSizeRadiusSquared) {
                continue;
              }
            }
          }
 
          // Scale point size if density scale is enabled
          if (densityScale) {
            float scale = getDensityScale(poData, p);
            worldSizePointCoord.x *= scale;
            worldSizePointCoord.y *= scale;
            worldSizeRadiusSquared *= scale * scale;
            if (roundPoint) {
              if (1.f < worldSizeRadiusSquared) continue;
            }
            else {
              if (1.f < glm::abs(worldSizePointCoord.x) || 1.f < glm::abs(worldSizePointCoord.y)) continue;
            }
          }
 
          // view space position of point
          glm::vec4 viewSpacePos_h = (viewFromOcttreeFrame[0] * p.x +
                                      viewFromOcttreeFrame[1] * p.y +
                                      viewFromOcttreeFrame[2] * p.z +
                                      viewFromOcttreeFrame[3]);
          float viewspace_depth = -viewSpacePos_h.z / viewSpacePos_h.w;      // flip sign due to looking along -Z.
 
          if (depthPoint) {
            //glm::vec2 viewSpacePointCoord = (1.f / worldSizePointCoord_h.w) * glm::vec2(worldSizePointCoord_h);
 
            float r = 0.f;
            if (roundPoint) {
              r = std::sqrt(std::min(worldSizeRadiusSquared, maxSizeRadiusSquared));
            }
            else {
              float rA = std::max(std::abs(worldSizePointCoord.x),
                                  std::abs(worldSizePointCoord.y));
              float rB = std::max(std::abs(maxSizePointCoord.x),
                                  std::abs(maxSizePointCoord.y));
              r = std::min(rA, rB);
            }
 
 
            viewspace_depth += viewSpaceRadius * r;
          }
 
          // See if we got a better hit
          if (viewspace_nearest < viewspace_depth) continue;
 
          // Yes, record
          viewspace_nearest = viewspace_depth;
          hit = true;
          indexBest = i;
        }
 
        if (hit) {
          worldspacePosBest = euclidean(worldFromOcttreeFrame * glm::vec4(cell->points[indexBest], 1.f));
          viewspaceDepthBest = viewspace_nearest;
          anyHits = true;
        }
      }
    }
    return anyHits;
  }
 
}
 
 
bool Cogs::Core::PotreePicker::pickCamera(Context* context,
                                          const CameraComponent& camera,
                                          const glm::vec2& normPosition,
                                          float /*rayLength*/,
                                          float radius,
                                          PickingFlags pickingFlags,
                                          PicksReturned returnFlag,
                                          const RayPicking::RayPickFilter& filter,
                                          std::vector<RayPicking::RayPickHit>& hits)
{
  const bool returnChildEntity = (pickingFlags & PickingFlags::ReturnChildEntity) == PickingFlags::ReturnChildEntity;
 
  bool hitSomething = false;
 
  const CameraData& cameraData = context->cameraSystem->getData(&camera);
 
  // Increase tolerance for main camera picking so you don't need to precisely hit a point
  float tolerance_norm = 2.f * radius / cameraData.viewportSize.x;
 
  for (const PotreeComponent& poComp : poSystem->pool) {
    if (!poComp.supportPicking || !poComp.isPickable() || filter.isUnwantedType(poComp)) { continue; }
 
    const RenderComponent* renderComp = poComp.getComponent<RenderComponent>();
    if (renderComp && !renderComp->isVisibleInLayer(filter.layerMask)) { continue; }
 
    PotreeDataHolder& poDataHolder = poSystem->getData(&poComp);
    PotreeData* poData = poDataHolder.poData.get();
    assert(poData);
 
    if (poData->state != PotreeState::Running || poData->octtreeFrame.tightBBoxMax.x < poData->octtreeFrame.tightBBoxMin.x) { continue; }
 
    bool roundPoint = false;
    bool depthPoint = false;
    switch (poComp.pointShape) {
    case PotreePointShape::Square:
      roundPoint = false;
      depthPoint = false;
      break;
    case PotreePointShape::Disc:
      roundPoint = true;
      depthPoint = false;
      break;
    case PotreePointShape::Paraboloid:
      roundPoint = false;
      depthPoint = true;
      break;
    case PotreePointShape::Sphere:
      roundPoint = true;
      depthPoint = true;
      break;
    default:
      assert(false && "Invalid case");
      break;
    }
 
    // Transform from octtree frame to world
    const glm::mat4& worldFromOcttreeFrame = poData->worldFromOcttreeFrame;
    const glm::mat4 viewFromOcttreeFrame = cameraData.viewMatrix * worldFromOcttreeFrame;
    const float modelSpaceRadius = poData->modelSpaceRadius;
    const float viewspaceRadius = std::cbrt(glm::determinant(glm::mat3(viewFromOcttreeFrame))) * modelSpaceRadius;
 
    const Component* poCompHit = nullptr;
    float viewspaceDepth = std::numeric_limits<float>::max();
    glm::vec3 worldspacePos;
 
    switch (poComp.pointSizing) {
    case PotreePointSizing::FixedSize:
      if (rayPickFixedSize(poData, cameraData,
                           viewFromOcttreeFrame,
                           normPosition,
                           radius,
                           viewspaceRadius,
                           viewspaceDepth, worldspacePos, roundPoint, depthPoint))
      {
        poCompHit = &poComp;
      }
      break;
    case PotreePointSizing::DistanceScaled:
      if (rayPickDistanceScaled(poData, cameraData,
                                viewFromOcttreeFrame,
                                normPosition,
                                radius,
                                tolerance_norm,
                                viewspaceRadius,
                                viewspaceDepth, worldspacePos, roundPoint, depthPoint, false, true))
      {
        poCompHit = &poComp;
      }
      break;
    case PotreePointSizing::DistanceAndDensityScaled:
      if (rayPickDistanceScaled(poData, cameraData,
                                viewFromOcttreeFrame,
                                normPosition,
                                radius,
                                tolerance_norm,
                                viewspaceRadius,
                                viewspaceDepth, worldspacePos, roundPoint, depthPoint, true, true))
      {
        poCompHit = &poComp;
      }
      break;
    default:
      assert(false && "Illegal pointsizing enum value");
      break;
    }
 
    if (poCompHit) {
      glm::vec4 viewPos = cameraData.viewMatrix * glm::vec4(glm::vec3(worldspacePos), 1.f);
      float viewDist = -viewPos.z;
 
      if (returnFlag == PicksReturned::Closest && !hits.empty()) {
        if (viewDist < hits[0].distance) {
          hits[0] = {*poCompHit, returnChildEntity, worldspacePos, viewDist};
          hitSomething = true;
        }
        // else, the intersection we found is further, so don't do anything
      }
      else {
        hits.emplace_back(*poCompHit, returnChildEntity, worldspacePos, viewDist);
        hitSomething = true;
      }
    }
  }
 
  // Sort hits based on distance
  if (returnFlag == PicksReturned::AllSorted) {
    std::sort(hits.begin(), hits.end());
  }
  return hitSomething;
}
 
bool Cogs::Core::PotreePicker::pickRay(Context* /*context*/,
                                       const glm::vec3& startPos,
                                       const glm::quat& rot,
                                       float rayLength,
                                       float radius,
                                       PickingFlags pickingFlags,
                                       PicksReturned returnFlag,
                                       const RayPicking::RayPickFilter& filter,
                                       std::vector<RayPicking::RayPickHit>& hits)
{
  const bool returnChildEntity = (pickingFlags & PickingFlags::ReturnChildEntity) == PickingFlags::ReturnChildEntity;
 
  bool hitSomething = false;
 
  for (const PotreeComponent& poComp : poSystem->pool) {
    if (!poComp.supportPicking || !poComp.isPickable() || filter.isUnwantedType(poComp)) { continue; }
 
    const RenderComponent* renderComp = poComp.getComponent<RenderComponent>();
    if (renderComp && !renderComp->isVisibleInLayer(filter.layerMask)) { continue; }
 
    PotreeDataHolder& poDataHolder = poSystem->getData(&poComp);
    PotreeData* poData = poDataHolder.poData.get();
    assert(poData);
 
    const glm::vec3& tbmin = poData->octtreeFrame.tightBBoxMin;
    const glm::vec3& tbmax = poData->octtreeFrame.tightBBoxMax;
    if (poData->state != PotreeState::Running || (tbmax.x < tbmin.x)) { continue; }
 
    bool roundPoint = false;
    bool depthPoint = false;
    switch (poComp.pointShape) {
    case PotreePointShape::Square:
      roundPoint = false;
      depthPoint = false;
      break;
    case PotreePointShape::Disc:
      roundPoint = true;
      depthPoint = false;
      break;
    case PotreePointShape::Paraboloid:
      roundPoint = false;
      depthPoint = true;
      break;
    case PotreePointShape::Sphere:
      roundPoint = true;
      depthPoint = true;
      break;
    default:
      assert(false && "Invalid case");
      break;
    }
 
    // Matrix that transforms a world position into the view space
    glm::mat4 viewMatrix = glm::toMat4(glm::conjugate(rot)) * glm::translate(glm::mat4(1.f), -startPos);
    //glm::quat cgate = glm::conjugate(rot);
 
    // Transform from octtree frame to world
    const glm::mat4& worldFromOcttreeFrame = poData->worldFromOcttreeFrame;
    const glm::mat4 viewFromOcttreeFrame = viewMatrix * worldFromOcttreeFrame;
 
    const float modelSpaceRadius = poData->modelSpaceRadius;
    const float viewSpaceRadius = std::cbrt(glm::determinant(glm::mat3(viewFromOcttreeFrame))) * modelSpaceRadius;
 
    // already simplified some operations (2/right - left) == (1/right) given that right == -left
    // term -(f+n)/(f-n); n == 0 => -f/f = -1.
    const float size = radius + poData->modelSpaceRadius;
    const glm::mat4 projectionMatrix(
      1 / size, 0, 0, 0,
      0, 1 / size, 0, 0,
      0, 0, -2 / rayLength, 0,
      0, 0, -1, 1
    );
 
    glm::vec2 query_pos_norm(0.0f, 0.0f);
    CameraData cameraData = CameraData();
    cameraData.rawProjectionMatrix = projectionMatrix;
    cameraData.viewportSize = glm::vec2(8.f);
    float tolerance = radius;
    float tolerance_norm = 2.f * tolerance / cameraData.viewportSize.x;
 
    const Component* poCompHit = nullptr;
    float viewspaceDepth = rayLength;
    glm::vec3 worldspacePos;
 
    switch (poComp.pointSizing) {
    case PotreePointSizing::FixedSize:
      LOG_TRACE(logger, "PotreePointSizing::FixedSize is currently not supported for camera-less picking, no hits detected");
      break;
    case PotreePointSizing::DistanceScaled:
      if (rayPickDistanceScaled(poData, cameraData,
        viewFromOcttreeFrame,
        query_pos_norm,
        tolerance,
        tolerance_norm,
        viewSpaceRadius,
        viewspaceDepth, worldspacePos, roundPoint, depthPoint, false, false))
      {
        poCompHit = &poComp;
      }
      break;
    case PotreePointSizing::DistanceAndDensityScaled:
      if (rayPickDistanceScaled(poData, cameraData,
        viewFromOcttreeFrame,
        query_pos_norm,
        tolerance,
        tolerance_norm,
        viewSpaceRadius,
        viewspaceDepth, worldspacePos, roundPoint, depthPoint, true, false))
      {
        poCompHit = &poComp;
      }
      break;
    default:
      assert(false && "Illegal pointsizing enum value");
      break;
    }
 
    if (poCompHit) {
      glm::vec4 viewPos = viewMatrix * glm::vec4(glm::vec3(worldspacePos), 1.f);
      float viewDist = -viewPos.z;
      if (returnFlag == PicksReturned::Closest && !hits.empty()) {
        if (viewDist < hits[0].distance) {
          hits[0] = {*poCompHit, returnChildEntity, worldspacePos, viewDist};
          hitSomething = true;
        }
        // else, the intersection we found is further, so don't do anything
      }
      else {
        hits.emplace_back(*poCompHit, returnChildEntity, worldspacePos, viewDist);
        hitSomething = true;
      }
    }
  }
 
  // Sort hits based on distance
  if (returnFlag == PicksReturned::AllSorted) {
    std::sort(hits.begin(), hits.end());
  }
  return hitSomething;
}