Vec3d indexing refactor

CMakeLists.txt

@@ -180,7 +180,6 @@ set(LIB_SOURCE_FILES
     src/h3lib/lib/polyfill.c
     src/h3lib/lib/h3Index.c
     src/h3lib/lib/vec2d.c
-    src/h3lib/lib/vec3d.c
     src/h3lib/lib/vertex.c
     src/h3lib/lib/linkedGeo.c
     src/h3lib/lib/localij.c
@@ -255,6 +254,7 @@ set(OTHER_SOURCE_FILES
     src/apps/testapps/testPolyfillInternal.c
     src/apps/testapps/testVec2dInternal.c
     src/apps/testapps/testVec3dInternal.c
+    src/apps/testapps/testVec3.c
     src/apps/testapps/testDirectedEdge.c
     src/apps/testapps/testDirectedEdgeExhaustive.c
     src/apps/testapps/testLinkedGeoInternal.c

CMakeTests.cmake

@@ -242,6 +242,7 @@ add_h3_test(testPolygonInternal src/apps/testapps/testPolygonInternal.c)
 add_h3_test(testPolyfillInternal src/apps/testapps/testPolyfillInternal.c)
 add_h3_test(testVec2dInternal src/apps/testapps/testVec2dInternal.c)
 add_h3_test(testVec3dInternal src/apps/testapps/testVec3dInternal.c)
+add_h3_test(testVec3 src/apps/testapps/testVec3.c)
 add_h3_test(testCellToLocalIj src/apps/testapps/testCellToLocalIj.c)
 add_h3_test(testCellToLocalIjInternal
             src/apps/testapps/testCellToLocalIjInternal.c)

src/apps/filters/h3.c

@@ -2642,7 +2642,7 @@ SUBCOMMAND(edgeLengthM,
     if (err) {
         return err;
     }
-    printf("%.10lf\n", length);
+    printf("%.8lf\n", length);
     return E_SUCCESS;
 }
 

src/apps/miscapps/generateFaceCenterPoint.c

@@ -1,5 +1,5 @@
 /*
- * Copyright 2018, 2020-2021 Uber Technologies, Inc.
+ * Copyright 2018, 2020-2021, 2026 Uber Technologies, Inc.
  *
  * Licensed under the Apache License, Version 2.0 (the "License");
  * you may not use this file except in compliance with the License.
@@ -56,8 +56,7 @@ static void generate(void) {
     printf("static const Vec3d faceCenterPoint[NUM_ICOSA_FACES] = {\n");
     for (int i = 0; i < NUM_ICOSA_FACES; i++) {
         LatLng centerCoords = faceCenterGeoCopy[i];
-        Vec3d centerPoint;
-        _geoToVec3d(&centerCoords, &centerPoint);
+        Vec3d centerPoint = latLngToVec3(centerCoords);
         printf("    {%.16f, %.16f, %.16f}, // face %2d\n", centerPoint.x,
                centerPoint.y, centerPoint.z, i);
     }

src/apps/testapps/testLatLngInternal.c

@@ -1,5 +1,5 @@
 /*
- * Copyright 2017-2021 Uber Technologies, Inc.
+ * Copyright 2017-2021, 2026 Uber Technologies, Inc.
  *
  * Licensed under the Apache License, Version 2.0 (the "License");
  * you may not use this file except in compliance with the License.
@@ -64,94 +64,7 @@ SUITE(latLngInternal) {
         t_assert(constrainLng(2 * M_PI) == 0, "lng 2pi");
         t_assert(constrainLng(3 * M_PI) == M_PI, "lng 2pi");
         t_assert(constrainLng(4 * M_PI) == 0, "lng 4pi");
-    }
-
-    TEST(_geoAzDistanceRads_noop) {
-        LatLng start = {15, 10};
-        LatLng out;
-        LatLng expected = {15, 10};
-
-        _geoAzDistanceRads(&start, 0, 0, &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "0 distance produces same point");
-    }
-
-    TEST(_geoAzDistanceRads_dueNorthSouth) {
-        LatLng start;
-        LatLng out;
-        LatLng expected;
-
-        // Due north to north pole
-        setGeoDegs(&start, 45, 1);
-        setGeoDegs(&expected, 90, 0);
-        _geoAzDistanceRads(&start, 0, H3_EXPORT(degsToRads)(45), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "due north to north pole produces north pole");
-
-        // Due north to south pole, which doesn't get wrapped correctly
-        setGeoDegs(&start, 45, 1);
-        setGeoDegs(&expected, 270, 1);
-        _geoAzDistanceRads(&start, 0, H3_EXPORT(degsToRads)(45 + 180), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "due north to south pole produces south pole");
-
-        // Due south to south pole
-        setGeoDegs(&start, -45, 2);
-        setGeoDegs(&expected, -90, 0);
-        _geoAzDistanceRads(&start, H3_EXPORT(degsToRads)(180),
-                           H3_EXPORT(degsToRads)(45), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "due south to south pole produces south pole");
-
-        // Due north to non-pole
-        setGeoDegs(&start, -45, 10);
-        setGeoDegs(&expected, -10, 10);
-        _geoAzDistanceRads(&start, 0, H3_EXPORT(degsToRads)(35), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "due north produces expected result");
-    }
-
-    TEST(_geoAzDistanceRads_poleToPole) {
-        LatLng start;
-        LatLng out;
-        LatLng expected;
-
-        // Azimuth doesn't really matter in this case. Any azimuth from the
-        // north pole is south, any azimuth from the south pole is north.
-
-        setGeoDegs(&start, 90, 0);
-        setGeoDegs(&expected, -90, 0);
-        _geoAzDistanceRads(&start, H3_EXPORT(degsToRads)(12),
-                           H3_EXPORT(degsToRads)(180), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "some direction to south pole produces south pole");
-
-        setGeoDegs(&start, -90, 0);
-        setGeoDegs(&expected, 90, 0);
-        _geoAzDistanceRads(&start, H3_EXPORT(degsToRads)(34),
-                           H3_EXPORT(degsToRads)(180), &out);
-        t_assert(geoAlmostEqual(&expected, &out),
-                 "some direction to north pole produces north pole");
-    }
-
-    TEST(_geoAzDistanceRads_invertible) {
-        LatLng start;
-        setGeoDegs(&start, 15, 10);
-        LatLng out;
-
-        double azimuth = H3_EXPORT(degsToRads)(20);
-        double degrees180 = H3_EXPORT(degsToRads)(180);
-        double distance = H3_EXPORT(degsToRads)(15);
-
-        _geoAzDistanceRads(&start, azimuth, distance, &out);
-        t_assert(fabs(H3_EXPORT(greatCircleDistanceRads)(&start, &out) -
-                      distance) < EPSILON_RAD,
-                 "moved distance is as expected");
-
-        LatLng start2 = out;
-        _geoAzDistanceRads(&start2, azimuth + degrees180, distance, &out);
-        // TODO: Epsilon is relatively large
-        t_assert(H3_EXPORT(greatCircleDistanceRads)(&start, &out) < 0.01,
-                 "moved back to origin");
+        t_assert(constrainLng(-2 * M_PI) == 0, "lng -2pi");
+        t_assert(constrainLng(-3 * M_PI) == -M_PI, "lng -3pi");
     }
 }

src/apps/testapps/testVec3.c

@@ -0,0 +1,148 @@
+/*
+ * Copyright 2026 Uber Technologies, Inc.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *         http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+/** @file testVec3.c
+ * @brief Tests the Vec3d helpers used by the geodesic polyfill path.
+ */
+
+#include <float.h>
+#include <math.h>
+
+#include "h3Index.h"
+#include "test.h"
+#include "vec3d.h"
+
+SUITE(Vec3d) {
+    TEST(dotProduct) {
+        Vec3d a = {.x = 1.0, .y = 0.0, .z = 0.0};
+        Vec3d b = {.x = -1.0, .y = 0.0, .z = 0.0};
+        t_assert(vec3Dot(a, b) == -1.0, "dot product matches expected value");
+    }
+
+    TEST(crossProductOrthogonality) {
+        Vec3d i = {.x = 1.0, .y = 0.0, .z = 0.0};
+        Vec3d j = {.x = 0.0, .y = 1.0, .z = 0.0};
+        Vec3d k = vec3Cross(i, j);
+        t_assert(fabs(k.x - 0.0) < DBL_EPSILON, "x component zero");
+        t_assert(fabs(k.y - 0.0) < DBL_EPSILON, "y component zero");
+        t_assert(fabs(k.z - 1.0) < DBL_EPSILON, "z component one");
+        t_assert(fabs(vec3Dot(k, i)) < DBL_EPSILON, "cross is orthogonal to i");
+        t_assert(fabs(vec3Dot(k, j)) < DBL_EPSILON, "cross is orthogonal to j");
+    }
+
+    TEST(normalizeAndMagnitude) {
+        Vec3d v = {.x = 3.0, .y = -4.0, .z = 12.0};
+        double magSq = vec3NormSq(v);
+        t_assert(fabs(magSq - 169.0) < DBL_EPSILON,
+                 "magnitude squared matches");
+        t_assert(fabs(vec3Norm(v) - 13.0) < DBL_EPSILON, "magnitude matches");
+
+        vec3Normalize(&v);
+        t_assert(fabs(vec3Norm(v) - 1.0) < DBL_EPSILON,
+                 "normalized vector is unit");
+
+        Vec3d zero = {.x = 0.0, .y = 0.0, .z = 0.0};
+        vec3Normalize(&zero);
+        t_assert(zero.x == 0.0 && zero.y == 0.0 && zero.z == 0.0,
+                 "zero vector remains unchanged when normalizing");
+    }
+
+    TEST(distance) {
+        Vec3d a = {.x = 0.0, .y = 0.0, .z = 0.0};
+        Vec3d b = {.x = 1.0, .y = 2.0, .z = 2.0};
+        t_assert(fabs(vec3DistSq(a, b) - 9.0) < DBL_EPSILON,
+                 "distance squared matches");
+    }
+
+    TEST(latLngToVec3_unitSphere) {
+        LatLng geo = {.lat = 0.5, .lng = -1.3};
+        Vec3d v = latLngToVec3(geo);
+        t_assert(fabs(vec3Norm(v) - 1.0) < DBL_EPSILON,
+                 "converted vector lives on the unit sphere");
+    }
+
+    TEST(vec3ToCell_invalidRes) {
+        Vec3d v = {.x = 1.0, .y = 0.0, .z = 0.0};
+        H3Index out;
+        t_assert(vec3ToCell(&v, -1, &out) == E_RES_DOMAIN,
+                 "negative resolution is rejected");
+        t_assert(vec3ToCell(&v, 16, &out) == E_RES_DOMAIN,
+                 "resolution above max is rejected");
+    }
+
+    TEST(cellToVec3_unitSphere) {
+        // cellToVec3 should return a point on the unit sphere.
+        LatLng p = {.lat = 0.6, .lng = -1.2};
+        H3Index h;
+        t_assertSuccess(H3_EXPORT(latLngToCell)(&p, 5, &h));
+
+        Vec3d v;
+        t_assertSuccess(cellToVec3(h, &v));
+        t_assert(fabs(vec3Norm(v) - 1.0) < DBL_EPSILON,
+                 "cellToVec3 result is on the unit sphere");
+    }
+
+    TEST(cellToVec3_matchesCellToLatLng) {
+        // vec3ToLatLng(cellToVec3(cell)) should agree with cellToLatLng.
+        LatLng p = {.lat = 0.3, .lng = 2.1};
+        H3Index h;
+        t_assertSuccess(H3_EXPORT(latLngToCell)(&p, 7, &h));
+
+        Vec3d v;
+        t_assertSuccess(cellToVec3(h, &v));
+        LatLng fromVec3 = vec3ToLatLng(v);
+
+        LatLng fromCell;
+        t_assertSuccess(H3_EXPORT(cellToLatLng)(h, &fromCell));
+
+        t_assert(fabs(fromVec3.lat - fromCell.lat) < DBL_EPSILON,
+                 "lat matches cellToLatLng");
+        t_assert(fabs(fromVec3.lng - fromCell.lng) < DBL_EPSILON,
+                 "lng matches cellToLatLng");
+    }
+
+    TEST(cellToVec3_roundTrip) {
+        // vec3ToCell(cellToVec3(cell)) should return the same cell.
+        LatLng p = {.lat = -0.4, .lng = 0.8};
+        H3Index h;
+        t_assertSuccess(H3_EXPORT(latLngToCell)(&p, 9, &h));
+
+        Vec3d v;
+        t_assertSuccess(cellToVec3(h, &v));
+
+        H3Index h2;
+        t_assertSuccess(vec3ToCell(&v, 9, &h2));
+        t_assert(h2 == h, "round-trip through Vec3d returns same cell");
+    }
+
+    TEST(cellToVec3_invalidCell) {
+        Vec3d v;
+        t_assert(cellToVec3(0x7fffffffffffffff, &v) == E_CELL_INVALID,
+                 "invalid cell gives E_CELL_INVALID");
+    }
+
+    TEST(vec3ToCell_nonFinite) {
+        H3Index out;
+        Vec3d nanX = {.x = NAN, .y = 0.0, .z = 0.0};
+        t_assert(vec3ToCell(&nanX, 0, &out) == E_DOMAIN, "NaN x is rejected");
+        Vec3d infY = {.x = 0.0, .y = INFINITY, .z = 0.0};
+        t_assert(vec3ToCell(&infY, 0, &out) == E_DOMAIN,
+                 "infinite y is rejected");
+        Vec3d infZ = {.x = 0.0, .y = 0.0, .z = -INFINITY};
+        t_assert(vec3ToCell(&infZ, 0, &out) == E_DOMAIN,
+                 "infinite z is rejected");
+    }
+}

src/apps/testapps/testVec3dInternal.c

@@ -1,5 +1,5 @@
 /*
- * Copyright 2018, 2020-2021 Uber Technologies, Inc.
+ * Copyright 2018, 2020-2021, 2026 Uber Technologies, Inc.
  *
  * Licensed under the Apache License, Version 2.0 (the "License");
  * you may not use this file except in compliance with the License.
@@ -16,50 +16,70 @@
 
 #include <float.h>
 #include <math.h>
-#include <stdlib.h>
 
 #include "test.h"
 #include "vec3d.h"
 
 SUITE(Vec3dInternal) {
-    TEST(_pointSquareDist) {
+    TEST(vec3DistSq) {
         Vec3d v1 = {0, 0, 0};
         Vec3d v2 = {1, 0, 0};
         Vec3d v3 = {0, 1, 1};
         Vec3d v4 = {1, 1, 1};
         Vec3d v5 = {1, 1, 2};
 
-        t_assert(fabs(_pointSquareDist(&v1, &v1)) < DBL_EPSILON,
+        t_assert(fabs(vec3DistSq(v1, v1)) < DBL_EPSILON,
                  "distance to self is 0");
-        t_assert(fabs(_pointSquareDist(&v1, &v2) - 1) < DBL_EPSILON,
+        t_assert(fabs(vec3DistSq(v1, v2) - 1) < DBL_EPSILON,
                  "distance to <1,0,0> is 1");
-        t_assert(fabs(_pointSquareDist(&v1, &v3) - 2) < DBL_EPSILON,
+        t_assert(fabs(vec3DistSq(v1, v3) - 2) < DBL_EPSILON,
                  "distance to <0,1,1> is 2");
-        t_assert(fabs(_pointSquareDist(&v1, &v4) - 3) < DBL_EPSILON,
+        t_assert(fabs(vec3DistSq(v1, v4) - 3) < DBL_EPSILON,
                  "distance to <1,1,1> is 3");
-        t_assert(fabs(_pointSquareDist(&v1, &v5) - 6) < DBL_EPSILON,
+        t_assert(fabs(vec3DistSq(v1, v5) - 6) < DBL_EPSILON,
                  "distance to <1,1,2> is 6");
     }
 
-    TEST(_geoToVec3d) {
+    TEST(vec3Normalize_smallNonzero) {
+        // 1e-163 squared underflows to 0, so norm == 0.
+        // vec3Normalize should produce the zero vector.
+        Vec3d v = {1e-163, 0, 0};
+
+        t_assert(v.x != 0.0, "vector is nonzero");
+        t_assert(vec3Norm(v) == 0.0, "norm underflows to zero");
+
+        vec3Normalize(&v);
+        t_assert(v.x == 0.0 && v.y == 0.0 && v.z == 0.0,
+                 "underflowed vector normalizes to zero");
+    }
+
+    TEST(vec3Normalize_dblEpsilonHalf) {
+        // DBL_EPSILON/2 is small but normalizes fine.
+        Vec3d v = {DBL_EPSILON / 2.0, 0, 0};
+
+        t_assert(vec3Norm(v) < DBL_EPSILON, "norm is small but nonzero");
+
+        vec3Normalize(&v);
+        t_assert(fabs(v.x - 1.0) < DBL_EPSILON && v.y == 0 && v.z == 0,
+                 "still normalizable to unit vector");
+    }
+
+    TEST(latLngToVec3) {
         Vec3d origin = {0};
 
         LatLng c1 = {0, 0};
-        Vec3d p1;
-        _geoToVec3d(&c1, &p1);
-        t_assert(fabs(_pointSquareDist(&origin, &p1) - 1) < EPSILON_RAD,
+        Vec3d p1 = latLngToVec3(c1);
+        t_assert(fabs(vec3DistSq(origin, p1) - 1) < EPSILON_RAD,
                  "Geo point is on the unit sphere");
 
         LatLng c2 = {M_PI_2, 0};
-        Vec3d p2;
-        _geoToVec3d(&c2, &p2);
-        t_assert(fabs(_pointSquareDist(&p1, &p2) - 2) < EPSILON_RAD,
+        Vec3d p2 = latLngToVec3(c2);
+        t_assert(fabs(vec3DistSq(p1, p2) - 2) < EPSILON_RAD,
                  "Geo point is on another axis");
 
         LatLng c3 = {M_PI, 0};
-        Vec3d p3;
-        _geoToVec3d(&c3, &p3);
-        t_assert(fabs(_pointSquareDist(&p1, &p3) - 4) < EPSILON_RAD,
+        Vec3d p3 = latLngToVec3(c3);
+        t_assert(fabs(vec3DistSq(p1, p3) - 4) < EPSILON_RAD,
                  "Geo point is the other side of the sphere");
     }
 }