Various improvements related to height queries against an ellipsoid

CHANGES.md

@@ -7,10 +7,13 @@
 - Added support for `EXT_accessor_additional_types` in `AccessorView`.
 - Added `EllipsoidTilesetLoader` that will generate a tileset by tesselating the surface of an ellipsoid, producing a simple globe tileset without any terrain features.
 - The `schemaUri` property in the `EXT_structural_metadata` glTF extension is now supported, allowing structural metadata schemas to be loaded from URIs rather than being embedded in the glTF itself.
+- Added `getHeightQuery` method to `TilesetContentLoader`, allowing loaders to optionally provide a custom, more efficient means of querying heights.
 
 ##### Fixes :wrench:
 
 - Updated the CMake install process to install the vcpkg-built Debug binaries in Debug builds. Previously the Release binaries were installed instead.
+- Tightened the tolerance of `IntersectionTests::rayTriangleParametric`, allowing it to find intersections with smaller triangles.
+- Fixed a bug that could cause `GltfUtilities::intersectRayGltfModel` to crash when the model contains a primitive whose position accessor does not have min/max values.
 
 ### v0.41.0 - 2024-11-01
 

Cesium3DTilesSelection/include/Cesium3DTilesSelection/EllipsoidTilesetLoader.h

@@ -8,7 +8,8 @@ namespace Cesium3DTilesSelection {
  * @brief A loader that will generate a tileset by tesselating the surface of an
  * ellipsoid, producing a simple globe tileset without any terrain features.
  */
-class EllipsoidTilesetLoader : public TilesetContentLoader {
+class EllipsoidTilesetLoader : public TilesetContentLoader,
+                               public ITilesetHeightQuery {
 public:
   /**
    * @brief Constructs a new instance.
@@ -35,6 +36,12 @@ class EllipsoidTilesetLoader : public TilesetContentLoader {
       const CesiumGeospatial::Ellipsoid& ellipsoid
           CESIUM_DEFAULT_ELLIPSOID) override;
 
+  ITilesetHeightQuery* getHeightQuery() override;
+
+  CesiumAsync::Future<SampleHeightResult> queryHeights(
+      const CesiumAsync::AsyncSystem& asyncSystem,
+      std::vector<CesiumGeospatial::Cartographic>&& positions) override;
+
 private:
   struct Geometry {
     std::vector<uint16_t> indices;

Cesium3DTilesSelection/include/Cesium3DTilesSelection/TilesetContentLoader.h

@@ -6,6 +6,7 @@
 #include "TileLoadResult.h"
 #include "TilesetOptions.h"
 
+#include <Cesium3DTilesSelection/SampleHeightResult.h>
 #include <CesiumAsync/AsyncSystem.h>
 #include <CesiumAsync/Future.h>
 #include <CesiumAsync/IAssetAccessor.h>
@@ -108,6 +109,25 @@ struct CESIUM3DTILESSELECTION_API TileChildrenResult {
   TileLoadResultState state;
 };
 
+/**
+ * @brief An interface to query heights from a tileset that can do so
+ * efficiently without necessarily downloading individual tiles.
+ */
+class CESIUM3DTILESSELECTION_API ITilesetHeightQuery {
+public:
+  /**
+   * @brief Queries the heights at a list of locations.
+   *
+   * @param asyncSystem The async system used to do work in threads.
+   * @param positions The positions at which to query heights. The height field
+   * of each {@link Cartographic} is ignored.
+   * @return A future that will be resolved when the heights have been queried.
+   */
+  virtual CesiumAsync::Future<SampleHeightResult> queryHeights(
+      const CesiumAsync::AsyncSystem& asyncSystem,
+      std::vector<CesiumGeospatial::Cartographic>&& positions) = 0;
+};
+
 /**
  * @brief The loader interface to load the tile content
  */
@@ -145,5 +165,24 @@ class CESIUM3DTILESSELECTION_API TilesetContentLoader {
       const Tile& tile,
       const CesiumGeospatial::Ellipsoid& ellipsoid
           CESIUM_DEFAULT_ELLIPSOID) = 0;
+
+  /**
+   * @brief Gets an interface that can be used to efficiently query heights from
+   * this tileset.
+   *
+   * Some loaders may be able to query heights very efficiently by using a web
+   * service or by using an analytical model, e.g., when the "terrain" is a
+   * simple ellipsoid.
+   *
+   * For loaders that have no particular way to query heights, this method will
+   * return `nullptr`, signaling that heights should be computed by downloading
+   * and sampling individual tiles.
+   *
+   * @return The interface that can be used to efficiently query heights from
+   * this loader, or `nullptr` if this loader has no particular way to do that.
+   * The returned instance must have a lifetime that is at least as long as the
+   * loader itself.
+   */
+  virtual ITilesetHeightQuery* getHeightQuery() { return nullptr; }
 };
 } // namespace Cesium3DTilesSelection

Cesium3DTilesSelection/src/EllipsoidTilesetLoader.cpp

@@ -69,7 +69,9 @@ TileChildrenResult EllipsoidTilesetLoader::createTileChildren(
   const QuadtreeTileID* pParentID =
       std::get_if<QuadtreeTileID>(&tile.getTileID());
 
-  if (pParentID) {
+  // Due to the use of uint32_t for QuadtreeTileID X and Y, we can only support
+  // through level 30.
+  if (pParentID && pParentID->level < 30) {
     std::vector<Tile> children;
     QuadtreeChildren childIDs =
         ImplicitTilingUtilities::getChildren(*pParentID);
@@ -87,6 +89,23 @@ TileChildrenResult EllipsoidTilesetLoader::createTileChildren(
   return TileChildrenResult{{}, TileLoadResultState::Failed};
 }
 
+ITilesetHeightQuery* EllipsoidTilesetLoader::getHeightQuery() { return this; }
+
+CesiumAsync::Future<SampleHeightResult> EllipsoidTilesetLoader::queryHeights(
+    const CesiumAsync::AsyncSystem& asyncSystem,
+    std::vector<CesiumGeospatial::Cartographic>&& positions) {
+  SampleHeightResult result;
+
+  result.positions = std::move(positions);
+  result.sampleSuccess.resize(result.positions.size(), true);
+
+  for (Cartographic& position : result.positions) {
+    position.height = 0.0;
+  }
+
+  return asyncSystem.createResolvedFuture(std::move(result));
+}
+
 void EllipsoidTilesetLoader::createChildTile(
     const Tile& parent,
     std::vector<Tile>& children,
@@ -120,7 +139,9 @@ EllipsoidTilesetLoader::Geometry
 EllipsoidTilesetLoader::createGeometry(const Tile& tile) const {
   static constexpr uint16_t resolution = 24;
 
-  std::vector<uint16_t> indices(6 * (resolution - 1) * (resolution - 1));
+  std::vector<uint16_t> indices;
+  indices.reserve(6 * (resolution - 1) * (resolution - 1));
+
   std::vector<glm::vec3> vertices(resolution * resolution);
   std::vector<glm::vec3> normals(vertices.size());
 

Cesium3DTilesSelection/src/Tileset.cpp

@@ -391,6 +391,7 @@ Tileset::updateView(const std::vector<ViewState>& frustums, float deltaTime) {
   }
 
   TilesetHeightRequest::processHeightRequests(
+      this->getAsyncSystem(),
       *this->_pTilesetContentManager,
       this->_options,
       this->_loadedTiles,

Cesium3DTilesSelection/src/TilesetHeightQuery.cpp

@@ -200,6 +200,7 @@ void TilesetHeightQuery::findCandidateTiles(
 }
 
 /*static*/ void TilesetHeightRequest::processHeightRequests(
+    const AsyncSystem& asyncSystem,
     TilesetContentManager& contentManager,
     const TilesetOptions& options,
     Tile::LoadedLinkedList& loadedTiles,
@@ -214,6 +215,7 @@ void TilesetHeightQuery::findCandidateTiles(
   for (auto it = heightRequests.begin(); it != heightRequests.end();) {
     TilesetHeightRequest& request = *it;
     if (!request.tryCompleteHeightRequest(
+            asyncSystem,
             contentManager,
             options,
             loadedTiles,
@@ -249,10 +251,35 @@ void Cesium3DTilesSelection::TilesetHeightRequest::failHeightRequests(
 }
 
 bool TilesetHeightRequest::tryCompleteHeightRequest(
+    const AsyncSystem& asyncSystem,
     TilesetContentManager& contentManager,
     const TilesetOptions& options,
     Tile::LoadedLinkedList& loadedTiles,
     std::set<Tile*>& tileLoadSet) {
+  // If this TilesetContentLoader supports direct height queries, use that
+  // instead of downloading tiles.
+  if (contentManager.getRootTile() &&
+      contentManager.getRootTile()->getLoader()) {
+    ITilesetHeightQuery* pQuery =
+        contentManager.getRootTile()->getLoader()->getHeightQuery();
+    if (pQuery) {
+      std::vector<Cartographic> positions;
+      positions.reserve(this->queries.size());
+      for (TilesetHeightQuery& query : this->queries) {
+        positions.emplace_back(query.inputPosition);
+      }
+
+      pQuery->queryHeights(asyncSystem, std::move(positions))
+          .thenImmediately(
+              [promise = this->promise](SampleHeightResult&& result) {
+                promise.resolve(std::move(result));
+              });
+
+      return true;
+    }
+  }
+
+  // No direct height query possible, so download and sample tiles.
   bool tileStillNeedsLoading = false;
   std::vector<std::string> warnings;
   for (TilesetHeightQuery& query : this->queries) {

Cesium3DTilesSelection/src/TilesetHeightQuery.h

@@ -126,6 +126,7 @@ struct TilesetHeightRequest {
    * @brief Process a given list of height requests. This is called by the {@link Tileset}
    * in every call to {@link Tileset::updateView}.
    *
+   * @param asyncSystem The async system used to do work in threads.
    * @param contentManager The content manager.
    * @param options Options associated with the tileset.
    * @param loadedTiles The linked list of loaded tiles, used to ensure that
@@ -137,6 +138,7 @@ struct TilesetHeightRequest {
    * height requests can complete are added to this vector.
    */
   static void processHeightRequests(
+      const CesiumAsync::AsyncSystem& asyncSystem,
       TilesetContentManager& contentManager,
       const TilesetOptions& options,
       Tile::LoadedLinkedList& loadedTiles,
@@ -160,6 +162,7 @@ struct TilesetHeightRequest {
    * @brief Tries to complete this height request. Returns false if further data
    * still needs to be loaded and thus the request cannot yet complete.
    *
+   * @param asyncSystem The async system used to do work in threads.
    * @param contentManager The content manager.
    * @param options Options associated with the tileset.
    * @param loadedTiles The linked list of loaded tiles, used to ensure that
@@ -169,6 +172,7 @@ struct TilesetHeightRequest {
    * can complete.
    */
   bool tryCompleteHeightRequest(
+      const CesiumAsync::AsyncSystem& asyncSystem,
       TilesetContentManager& contentManager,
       const TilesetOptions& options,
       Tile::LoadedLinkedList& loadedTiles,

Cesium3DTilesSelection/test/TestTilesetHeightQueries.cpp

@@ -1,4 +1,5 @@
 #include <Cesium3DTilesContent/registerAllTileContentTypes.h>
+#include <Cesium3DTilesSelection/EllipsoidTilesetLoader.h>
 #include <Cesium3DTilesSelection/Tileset.h>
 #include <CesiumGeospatial/Cartographic.h>
 #include <CesiumNativeTests/FileAccessor.h>
@@ -230,4 +231,38 @@ TEST_CASE("Tileset height queries") {
     CHECK(!results.sampleSuccess[0]);
     CHECK(results.warnings[0].find("failed to load") != std::string::npos);
   }
+
+  SECTION("ellipsoid tileset") {
+    std::unique_ptr<Tileset> pTileset =
+        EllipsoidTilesetLoader::createTileset(externals);
+
+    Future<SampleHeightResult> future = pTileset->sampleHeightMostDetailed(
+        {Cartographic::fromDegrees(-75.612559, 40.042183, 1.0)});
+
+    while (!future.isReady()) {
+      pTileset->updateView({});
+    }
+
+    SampleHeightResult results = future.waitInMainThread();
+
+    REQUIRE(results.warnings.size() == 0);
+    REQUIRE(results.positions.size() == 1);
+    REQUIRE(results.sampleSuccess.size() == 1);
+    CHECK(results.sampleSuccess[0]);
+    CHECK(Math::equalsEpsilon(
+        results.positions[0].longitude,
+        Math::degreesToRadians(-75.612559),
+        0.0,
+        Math::Epsilon4));
+    CHECK(Math::equalsEpsilon(
+        results.positions[0].latitude,
+        Math::degreesToRadians(40.042183),
+        0.0,
+        Math::Epsilon4));
+    CHECK(Math::equalsEpsilon(
+        results.positions[0].height,
+        0.0,
+        0.0,
+        Math::Epsilon4));
+  }
 }

CesiumGeometry/src/IntersectionTests.cpp

@@ -141,7 +141,7 @@ std::optional<double> IntersectionTests::rayTriangleParametric(
   glm::dvec3 p = glm::cross(direction, edge1);
   double det = glm::dot(edge0, p);
   if (cullBackFaces) {
-    if (det < Math::Epsilon6)
+    if (det < Math::Epsilon8)
       return std::nullopt;
 
     glm::dvec3 tvec = origin - p0;
@@ -158,7 +158,7 @@ std::optional<double> IntersectionTests::rayTriangleParametric(
 
   } else {
 
-    if (glm::abs(det) < Math::Epsilon6)
+    if (glm::abs(det) < Math::Epsilon8)
       return std::nullopt;
 
     double invDet = 1.0 / det;

CesiumGltfContent/src/GltfUtilities.cpp

@@ -1264,22 +1264,25 @@ std::optional<glm::dvec3> intersectRayScenePrimitive(
   glm::dmat4x4 worldToPrimitive = glm::inverse(primitiveToWorld);
   CesiumGeometry::Ray transformedRay = ray.transform(worldToPrimitive);
 
-  // Ignore primitive if ray doesn't intersect bounding box
+  // Ignore primitive if we have an AABB from the accessor min/max and the ray
+  // doesn't intersect it.
   const std::vector<double>& min = positionAccessor.min;
   const std::vector<double>& max = positionAccessor.max;
 
-  std::optional<double> boxT =
-      CesiumGeometry::IntersectionTests::rayAABBParametric(
-          transformedRay,
-          CesiumGeometry::AxisAlignedBox(
-              min[0],
-              min[1],
-              min[2],
-              max[0],
-              max[1],
-              max[2]));
-  if (!boxT)
-    return std::optional<glm::dvec3>();
+  if (min.size() >= 3 && max.size() >= 3) {
+    std::optional<double> boxT =
+        CesiumGeometry::IntersectionTests::rayAABBParametric(
+            transformedRay,
+            CesiumGeometry::AxisAlignedBox(
+                min[0],
+                min[1],
+                min[2],
+                max[0],
+                max[1],
+                max[2]));
+    if (!boxT)
+      return std::optional<glm::dvec3>();
+  }
 
   double tClosest = -1.0;