QuadricToEllipsoid extends
	AbstractUnaryFunctionOp<Matrix4dc, Optional<Ellipsoid>> { @Override public Optional<Ellipsoid> calculate(final Matrix4dc quadricSolution) { final Vector3dc center = findCenter(quadricSolution); final Matrix4dc translated = translateToCenter(quadricSolution, center); final EigenDecomposition decomposition = solveEigenDecomposition( translated); final boolean invalidEv = Arrays.stream(decomposition.getRealEigenvalues()) .anyMatch(e -> e <= 0.0); if (invalidEv) { return Optional.empty(); } final double[] radii = Arrays.stream(decomposition.getRealEigenvalues()) .map(ev -> Math.sqrt(1.0 / ev)).toArray(); final Ellipsoid ellipsoid = new Ellipsoid(radii[0], radii[1], radii[2]); ellipsoid.setCentroid(center); final Matrix3dc orientation = toOrientationMatrix(decomposition); ellipsoid.setOrientation(orientation); return Optional.of(ellipsoid); }  @Override Optional<Ellipsoid> calculate(final Matrix4dc quadricSolution);  }
@Test public void testCone() { final Matrix4d cone = new Matrix4d( 1, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1, 0, 0, 0, 0, -1); IMAGE_J.op().run(QuadricToEllipsoid.class, cone); final Optional<Ellipsoid> result = quadricToEllipsoid.calculate(cone); assertFalse(result.isPresent()); }
@Test public void testEllipsoidBand() { final double[] radii = { 1, 2, 3 }; final Matrix4dc symmetry = new Matrix4d( 1, 0, 0, 0, 0, -1, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1); final Random rng = new Random(SEED); final List<Vector3d> points = ellipsoidPoints.calculate(radii, 1_000L); points.forEach(p -> { final double scale = (2 * rng.nextDouble() - 1) * 0.05 + 1.0; p.mul(scale); }); final Matrix4dc quadric = (Matrix4dc) IMAGE_J.op().run(Quadric.class, points); final Optional<Ellipsoid> result = quadricToEllipsoid.calculate(quadric); assertTrue(result.isPresent()); final Ellipsoid ellipsoid = result.get(); final Vector3dc centroid = new Vector3d(ellipsoid.getCentroid().x, ellipsoid .getCentroid().y, ellipsoid.getCentroid().z); assertTrue("Ellipsoid centre point is not within tolerance", epsilonEquals( new Vector3d(0, 0, 0), centroid, 0.05)); assertEquals(radii[0], ellipsoid.getA(), 0.025); assertEquals(radii[1], ellipsoid.getB(), 0.025); assertEquals(radii[2], ellipsoid.getC(), 0.025); epsilonEquals(symmetry, ellipsoid.getOrientation(), 0.025); }
@Test public void testTransformedEllipsoid() { final Vector3dc centroid = new Vector3d(1, 1, 1); final Quaterniondc q = new Quaterniond(new AxisAngle4d(Math.PI / 4.0, 0, 0, 1)); final double[] radii = { 1, 2, 3 }; final BinaryHybridCFI1<Vector3d, Quaterniondc, Vector3d> rotate = Hybrids .binaryCFI1(IMAGE_J.op(), Rotate3d.class, Vector3d.class, new Vector3d(), q); final double alpha = Math.sin(Math.PI / 4.0); final Matrix4dc orientation = new Matrix4d( alpha, alpha, 0, 0, alpha, -alpha, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1); final List<Vector3d> points = ellipsoidPoints.calculate(radii, 1_000L); points.forEach(rotate::mutate); points.forEach(p -> p.add(centroid)); final Matrix4dc quadric = (Matrix4dc) IMAGE_J.op().run(Quadric.class, points); final Optional<Ellipsoid> result = quadricToEllipsoid.calculate(quadric); assertTrue(result.isPresent()); final Ellipsoid transformedEllipsoid = result.get(); final Vector3d v = transformedEllipsoid.getCentroid(); assertTrue(epsilonEquals(centroid, new Vector3d(v.x, v.y, v.z), 1e-12)); assertEquals(radii[0], transformedEllipsoid.getA(), 1e-12); assertEquals(radii[1], transformedEllipsoid.getB(), 1e-12); assertEquals(radii[2], transformedEllipsoid.getC(), 1e-12); epsilonEquals(orientation, transformedEllipsoid.getOrientation(), 0.025); }
@Test public void testUnitSphere() { final Matrix4d expectedOrientation = new Matrix4d().identity(); final Optional<Ellipsoid> result = quadricToEllipsoid.calculate( UNIT_SPHERE); assertTrue(result.isPresent()); final Ellipsoid unitSphere = result.get(); assertEquals(1.0, unitSphere.getA(), 1e-12); assertEquals(1.0, unitSphere.getB(), 1e-12); assertEquals(1.0, unitSphere.getC(), 1e-12); final Vector3d c = unitSphere.getCentroid(); assertTrue(epsilonEquals(new Vector3d(0, 0, 0), new Vector3d(c.x, c.y, c.z), 1e-12)); epsilonEquals(expectedOrientation, unitSphere.getOrientation(), 1e-12); }