Hi all, I'm attempting to wrap up the refactoring work I've been doing with GEOMETRY-32 but I've hit a snag with the code for computing barycenters in 2D spherical space. The results I'm getting seem mostly reasonable but they aren't as consistent as I'd like, which makes me wonder if I'm doing something wrong or if I've made an invalid assumption somewhere. The previous barycenter code was not really unit tested so I don't have much to compare to.
Here is what I have currently for computing the barycenter of convex regions with 3 or more sides: - Split the region into a series of triangles using a triangle fan approach (eg, if the boundary points are p0, p1, p2, p3, split the region into the triangles p0, p1, p2 and p0, p2, p3) - Compute the area of each triangle by calculating the interior angles A, B, C and using Girard's Theorem (area = A + B + C - pi) - Compute the barycenter of each triangle using the formula for the center of mass found on the first page of the article here [1]. This is done by multiplying the pole of each side's great circle by its arc length and summing the result vectors. The sum vector is normalized to place it on the surface of the unit sphere. (This is the approach used in the previous code for the barycenter [2]) - Compute the overall barycenter for the region by scaling each triangle barycenter vector by the triangle area, adding the results, and normalizing. As mentioned above, this seems to produce mostly reasonable results. The issue I'm having is with verifying their consistency. In my unit tests, I'm doing the following to check the barycenters: - Compute a number of great circles that pass through the region barycenter - Use each circle to split the area into two smaller convex regions - a plus and a minus region - Compute the barycenters and areas for the plus and minus regions using the same logic as the larger region - Combine the two barycenters by scaling their vectors by the region's area, summing, and normalizing - Verify that the barycenter computed from the split parts matches the one computed for the region as a whole This basically checks that the region is the sum of its parts in terms of barycenter. However, I'm finding that the barycenter computed from the split regions is always slightly off from the overall region barycenter, generally by a margin of between 0.001 to 0.01 radians in both azimuth and polar. So, I'm wondering if this is a math issue, a floating point issue, or something else. Any insight would be greatly appreciated. My code can be found here [3] and here [4], although I haven't yet pushed the changes for the full algorithm I outlined above. Regards, Matt J [1] https://asmedigitalcollection.asme.org/appliedmechanics/article/42/1/239/387501/The-Inertia-Tensor-for-a-Spherical-Triangle [2] https://github.com/darkma773r/commons-geometry/blob/master/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/PropertiesComputer.java#L135 [3] https://github.com/darkma773r/commons-geometry/blob/geometry-32-working/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/ConvexArea2S.java [4] https://github.com/darkma773r/commons-geometry/blob/geometry-32-working/commons-geometry-spherical/src/test/java/org/apache/commons/geometry/spherical/twod/ConvexArea2STest.java#L188