use of net.imagej.ops.geom.geom3d.mesh.DefaultMesh in project imagej-ops by imagej.
the class DefaultConvexHull3D method calculate.
@Override
public Mesh calculate(final Mesh input) {
DefaultMesh output = new DefaultMesh();
Set<Vertex> vertices = new LinkedHashSet<>();
for (final RealLocalizable v : input.getVertices()) {
vertices.add(new Vertex(v.getDoublePosition(0), v.getDoublePosition(1), v.getDoublePosition(2)));
}
List<TriangularFacet> facets = new ArrayList<>();
List<TriangularFacet> facetsWithPointInFront = new ArrayList<>();
final double epsilon = computeHull(vertices, facets, facetsWithPointInFront);
for (TriangularFacet f : facets) {
output.addFace(f);
}
output.setEpsilon(epsilon);
return output;
}
use of net.imagej.ops.geom.geom3d.mesh.DefaultMesh in project imagej-ops by imagej.
the class DefaultVoxelization3D method calculate.
@Override
public RandomAccessibleInterval<BitType> calculate(Mesh input) {
Img<BitType> outImg = ops.create().img(new FinalInterval(width, height, depth), new BitType());
DefaultMesh dMesh = (DefaultMesh) input;
Set<RealLocalizable> verts = dMesh.getVertices();
RealPoint minPoint = new RealPoint(verts.iterator().next());
RealPoint maxPoint = new RealPoint(verts.iterator().next());
for (RealLocalizable v : verts) {
if (v.getDoublePosition(0) < minPoint.getDoublePosition(0))
minPoint.setPosition(v.getDoublePosition(0), 0);
if (v.getDoublePosition(1) < minPoint.getDoublePosition(1))
minPoint.setPosition(v.getDoublePosition(1), 1);
if (v.getDoublePosition(2) < minPoint.getDoublePosition(2))
minPoint.setPosition(v.getDoublePosition(2), 2);
if (v.getDoublePosition(0) > maxPoint.getDoublePosition(0))
maxPoint.setPosition(v.getDoublePosition(0), 0);
if (v.getDoublePosition(1) > maxPoint.getDoublePosition(1))
maxPoint.setPosition(v.getDoublePosition(1), 1);
if (v.getDoublePosition(2) > maxPoint.getDoublePosition(2))
maxPoint.setPosition(v.getDoublePosition(2), 2);
}
RealPoint dimPoint = new RealPoint((maxPoint.getDoublePosition(0) - minPoint.getDoublePosition(0)), (maxPoint.getDoublePosition(1) - minPoint.getDoublePosition(1)), (maxPoint.getDoublePosition(2) - minPoint.getDoublePosition(2)));
double[] stepSizes = new double[3];
stepSizes[0] = dimPoint.getDoublePosition(0) / width;
stepSizes[1] = dimPoint.getDoublePosition(1) / height;
stepSizes[2] = dimPoint.getDoublePosition(2) / depth;
double[] voxelHalfsize = new double[3];
for (int k = 0; k < stepSizes.length; k++) voxelHalfsize[k] = stepSizes[k] / 2.0;
for (Facet f : dMesh.getFacets()) {
TriangularFacet tri = (TriangularFacet) f;
Vector3D v1 = tri.getP0();
Vector3D v2 = tri.getP1();
Vector3D v3 = tri.getP2();
double[] minSubBoundary = new double[] { Math.min(Math.min(v1.getX(), v2.getX()), v3.getX()) - minPoint.getDoublePosition(0), Math.min(Math.min(v1.getY(), v2.getY()), v3.getY()) - minPoint.getDoublePosition(1), Math.min(Math.min(v1.getZ(), v2.getZ()), v3.getZ()) - minPoint.getDoublePosition(2) };
double[] maxSubBoundary = new double[] { Math.max(Math.max(v1.getX(), v2.getX()), v3.getX()) - minPoint.getDoublePosition(0), Math.max(Math.max(v1.getY(), v2.getY()), v3.getY()) - minPoint.getDoublePosition(1), Math.max(Math.max(v1.getZ(), v2.getZ()), v3.getZ()) - minPoint.getDoublePosition(2) };
// Should use the
RandomAccess<BitType> ra = outImg.randomAccess();
// interval
// implementation
// for speed
long[] indices = new long[3];
for (indices[0] = (long) Math.floor(minSubBoundary[0] / stepSizes[0]); indices[0] < Math.floor(maxSubBoundary[0] / stepSizes[0]); indices[0]++) {
for (indices[1] = (long) Math.floor(minSubBoundary[1] / stepSizes[1]); indices[1] < Math.floor(maxSubBoundary[1] / stepSizes[1]); indices[1]++) {
for (indices[2] = (long) Math.floor(minSubBoundary[2] / stepSizes[2]); indices[2] < Math.floor(maxSubBoundary[2] / stepSizes[2]); indices[2]++) {
ra.setPosition(indices);
if (// Don't check if voxel is already
!ra.get().get()) // filled
{
double[] voxelCenter = new double[3];
for (int k = 0; k < 3; k++) voxelCenter[k] = indices[k] * stepSizes[k] + voxelHalfsize[k];
if (triBoxOverlap(voxelCenter, voxelHalfsize, v1, v2, v3) == 1) {
ra.get().set(true);
}
}
}
}
}
}
return outImg;
}
use of net.imagej.ops.geom.geom3d.mesh.DefaultMesh in project imagej-ops by imagej.
the class QuickHull3DTest method quickhull_12_Test.
@Test
public void quickhull_12_Test() {
LinkedHashSet<RealLocalizable> points = new LinkedHashSet<>();
points.add(new Vertex(-0.03621271768232132, 0.3728502838619522, 0.4947140370446388));
points.add(new Vertex(0.3210853052521919, 0.4807189479290684, 0.4433501688235907));
points.add(new Vertex(0.07214279572678994, -0.4960366976410492, 0.1112227161519441));
points.add(new Vertex(0.2229772524190855, -0.4213242506806965, -0.1966818060695024));
points.add(new Vertex(-0.3411871756810576, -0.3328629143842151, -0.4270033635450559));
points.add(new Vertex(-0.245701439441835, 0.495905311308713, -0.3194406286994373));
points.add(new Vertex(0.458374538420117, -0.09914027349943322, -0.2505798421339875));
points.add(new Vertex(-0.4954086979808367, -0.3339869997780649, -0.3195065691317492));
points.add(new Vertex(-0.3392973838740004, 0.4288679723896719, -0.01599531622230571));
points.add(new Vertex(0.2724846394476338, -0.3506708492996831, 0.2750346518820475));
points.add(new Vertex(0.3544683273457627, -0.450828987127942, -0.0827870439577727));
points.add(new Vertex(0.1667164640191164, 0.003605551555385444, -0.4014989499947977));
DefaultMesh df = new DefaultMesh(points);
DefaultMesh convexHull = (DefaultMesh) ops.run(DefaultConvexHull3D.class, df);
assertTrue(isConvex(convexHull.getFacets(), convexHull.getEpsilon()));
assertEquals(12, convexHull.getVertices().size());
}
use of net.imagej.ops.geom.geom3d.mesh.DefaultMesh in project imagej-ops by imagej.
the class DefaultMarchingCubes method calculate.
@SuppressWarnings({ "unchecked" })
@Override
public DefaultMesh calculate(final RandomAccessibleInterval<T> input) {
DefaultMesh output = new DefaultMesh();
ExtendedRandomAccessibleInterval<T, RandomAccessibleInterval<T>> extended = Views.extendValue(input, (T) new BoolType(false));
Cursor<T> c = Views.interval(extended, new FinalInterval(new long[] { input.min(0) - 1, input.min(1) - 1, input.min(2) - 1 }, new long[] { input.max(0) + 1, input.max(1) + 1, input.max(2) + 1 })).localizingCursor();
while (c.hasNext()) {
c.next();
int cursorX = c.getIntPosition(0);
int cursorY = c.getIntPosition(1);
int cursorZ = c.getIntPosition(2);
Cursor<T> cu = getCube(extended, cursorX, cursorY, cursorZ);
int i = 0;
double[] vertex_values = new double[8];
while (cu.hasNext()) {
vertex_values[i++] = (cu.next().get()) ? 1 : 0;
}
// 6------7
// /| /|
// 2-----3 |
// | 4---|-5
// |/ |/
// 0-----1
vertex_values = mapFlatIterableToLookUpCube(vertex_values);
// 4------5
// /| /|
// 7-----6 |
// | 0---|-1
// |/ |/
// 3-----2
int cubeindex = getCubeIndex(vertex_values);
if (EDGE_TABLE[cubeindex] != 0) {
int[] p0 = new int[] { 0 + cursorX, 0 + cursorY, 1 + cursorZ };
int[] p1 = new int[] { 1 + cursorX, 0 + cursorY, 1 + cursorZ };
int[] p2 = new int[] { 1 + cursorX, 0 + cursorY, 0 + cursorZ };
int[] p3 = new int[] { 0 + cursorX, 0 + cursorY, 0 + cursorZ };
int[] p4 = new int[] { 0 + cursorX, 1 + cursorY, 1 + cursorZ };
int[] p5 = new int[] { 1 + cursorX, 1 + cursorY, 1 + cursorZ };
int[] p6 = new int[] { 1 + cursorX, 1 + cursorY, 0 + cursorZ };
int[] p7 = new int[] { 0 + cursorX, 1 + cursorY, 0 + cursorZ };
double[][] vertlist = new double[12][];
/* Find the vertices where the surface intersects the cube */
if (0 != (EDGE_TABLE[cubeindex] & 1)) {
vertlist[0] = interpolatePoint(p0, p1, vertex_values[0], vertex_values[1]);
}
if (0 != (EDGE_TABLE[cubeindex] & 2)) {
vertlist[1] = interpolatePoint(p1, p2, vertex_values[1], vertex_values[2]);
}
if (0 != (EDGE_TABLE[cubeindex] & 4)) {
vertlist[2] = interpolatePoint(p2, p3, vertex_values[2], vertex_values[3]);
}
if (0 != (EDGE_TABLE[cubeindex] & 8)) {
vertlist[3] = interpolatePoint(p3, p0, vertex_values[3], vertex_values[0]);
}
if (0 != (EDGE_TABLE[cubeindex] & 16)) {
vertlist[4] = interpolatePoint(p4, p5, vertex_values[4], vertex_values[5]);
}
if (0 != (EDGE_TABLE[cubeindex] & 32)) {
vertlist[5] = interpolatePoint(p5, p6, vertex_values[5], vertex_values[6]);
}
if (0 != (EDGE_TABLE[cubeindex] & 64)) {
vertlist[6] = interpolatePoint(p6, p7, vertex_values[6], vertex_values[7]);
}
if (0 != (EDGE_TABLE[cubeindex] & 128)) {
vertlist[7] = interpolatePoint(p7, p4, vertex_values[7], vertex_values[4]);
}
if (0 != (EDGE_TABLE[cubeindex] & 256)) {
vertlist[8] = interpolatePoint(p0, p4, vertex_values[0], vertex_values[4]);
}
if (0 != (EDGE_TABLE[cubeindex] & 512)) {
vertlist[9] = interpolatePoint(p1, p5, vertex_values[1], vertex_values[5]);
}
if (0 != (EDGE_TABLE[cubeindex] & 1024)) {
vertlist[10] = interpolatePoint(p2, p6, vertex_values[2], vertex_values[6]);
}
if (0 != (EDGE_TABLE[cubeindex] & 2048)) {
vertlist[11] = interpolatePoint(p3, p7, vertex_values[3], vertex_values[7]);
}
/* Create the triangle */
for (i = 0; TRIANGLE_TABLE[cubeindex][i] != -1; i += 3) {
TriangularFacet face = new TriangularFacet(new Vertex(vertlist[TRIANGLE_TABLE[cubeindex][i + 2]][0], vertlist[TRIANGLE_TABLE[cubeindex][i + 2]][1], vertlist[TRIANGLE_TABLE[cubeindex][i + 2]][2]), new Vertex(vertlist[TRIANGLE_TABLE[cubeindex][i + 1]][0], vertlist[TRIANGLE_TABLE[cubeindex][i + 1]][1], vertlist[TRIANGLE_TABLE[cubeindex][i + 1]][2]), new Vertex(vertlist[TRIANGLE_TABLE[cubeindex][i]][0], vertlist[TRIANGLE_TABLE[cubeindex][i]][1], vertlist[TRIANGLE_TABLE[cubeindex][i]][2]));
face.getArea();
output.addFace(face);
}
}
}
return output;
}
use of net.imagej.ops.geom.geom3d.mesh.DefaultMesh in project imagej-ops by imagej.
the class QuickHull3DTest method quickhull_100_000_Test.
@Test
public void quickhull_100_000_Test() {
DefaultMesh df = new DefaultMesh(randomPointSet(100000, 20150818));
DefaultMesh convexHull = (DefaultMesh) ops.run(DefaultConvexHull3D.class, df);
assertTrue(isConvex(convexHull.getFacets(), convexHull.getEpsilon()));
assertEquals(175, convexHull.getVertices().size());
}
Aggregations