/* * The MIT License * * Copyright (c) 2016-2021 JOML * * Permission is hereby granted, free of charge, to any person obtaining a copy * of this software and associated documentation files (the "Software"), to deal * in the Software without restriction, including without limitation the rights * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell * copies of the Software, and to permit persons to whom the Software is * furnished to do so, subject to the following conditions: * * The above copyright notice and this permission notice shall be included in * all copies or substantial portions of the Software. * * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN * THE SOFTWARE. */ package com.jozufozu.flywheel.repack.joml.sampling; import java.nio.FloatBuffer; import java.util.ArrayList; import com.jozufozu.flywheel.repack.joml.Random; import com.jozufozu.flywheel.repack.joml.Vector2f; import com.jozufozu.flywheel.repack.joml.Vector3f; /** * Creates samples using the "Best Candidate" algorithm. * * @author Kai Burjack */ public class BestCandidateSampling { private static final class IntHolder { int value; } /** * Generates Best Candidate samples on a unit sphere. *

* References: *

* * @author Kai Burjack */ public static class Sphere { /** * Implementation of a Hierarchical Triangular Mesh structure to index the sample points on the unit sphere for accelerating 1-nearest neighbor searches. * * @author Kai Burjack */ private static final class Node { private static final int MAX_OBJECTS_PER_NODE = 32; private float v0x, v0y, v0z; private float v1x, v1y, v1z; private float v2x, v2y, v2z; private float cx, cy, cz; private float arc; private ArrayList objects; private Node[] children; Node() { this.children = new Node[8]; float s = 1f; this.arc = 2.0f * (float) Math.PI; /* * See: https://arxiv.org/ftp/cs/papers/0701/0701164.pdf */ this.children[0] = new Node(-s, 0, 0, 0, 0, s, 0, s, 0); this.children[1] = new Node(0, 0, s, s, 0, 0, 0, s, 0); this.children[2] = new Node(s, 0, 0, 0, 0, -s, 0, s, 0); this.children[3] = new Node(0, 0, -s, -s, 0, 0, 0, s, 0); this.children[4] = new Node(-s, 0, 0, 0, -s, 0, 0, 0, s); this.children[5] = new Node(0, 0, s, 0, -s, 0, s, 0, 0); this.children[6] = new Node(s, 0, 0, 0, -s, 0, 0, 0, -s); this.children[7] = new Node(0, 0, -s, 0, -s, 0, -s, 0, 0); } private Node(float x0, float y0, float z0, float x1, float y1, float z1, float x2, float y2, float z2) { this.v0x = x0; this.v0y = y0; this.v0z = z0; this.v1x = x1; this.v1y = y1; this.v1z = z1; this.v2x = x2; this.v2y = y2; this.v2z = z2; cx = (v0x + v1x + v2x) / 3.0f; cy = (v0y + v1y + v2y) / 3.0f; cz = (v0z + v1z + v2z) / 3.0f; float invCLen = Math.invsqrt(cx * cx + cy * cy + cz * cz); cx *= invCLen; cy *= invCLen; cz *= invCLen; // Compute maximum radius around triangle centroid float arc1 = greatCircleDist(cx, cy, cz, v0x, v0y, v0z); float arc2 = greatCircleDist(cx, cy, cz, v1x, v1y, v1z); float arc3 = greatCircleDist(cx, cy, cz, v2x, v2y, v2z); float dist = Math.max(Math.max(arc1, arc2), arc3); /* * Convert radius into diameter, but also take into account the linear * arccos approximation we use. * This value was obtained empirically! */ dist *= 1.7f; this.arc = dist; } private void split() { float w0x = v1x + v2x; float w0y = v1y + v2y; float w0z = v1z + v2z; float len0 = Math.invsqrt(w0x * w0x + w0y * w0y + w0z * w0z); w0x *= len0; w0y *= len0; w0z *= len0; float w1x = v0x + v2x; float w1y = v0y + v2y; float w1z = v0z + v2z; float len1 = Math.invsqrt(w1x * w1x + w1y * w1y + w1z * w1z); w1x *= len1; w1y *= len1; w1z *= len1; float w2x = v0x + v1x; float w2y = v0y + v1y; float w2z = v0z + v1z; float len2 = Math.invsqrt(w2x * w2x + w2y * w2y + w2z * w2z); w2x *= len2; w2y *= len2; w2z *= len2; children = new Node[4]; /* * See: https://arxiv.org/ftp/cs/papers/0701/0701164.pdf */ children[0] = new Node(v0x, v0y, v0z, w2x, w2y, w2z, w1x, w1y, w1z); children[1] = new Node(v1x, v1y, v1z, w0x, w0y, w0z, w2x, w2y, w2z); children[2] = new Node(v2x, v2y, v2z, w1x, w1y, w1z, w0x, w0y, w0z); children[3] = new Node(w0x, w0y, w0z, w1x, w1y, w1z, w2x, w2y, w2z); } private void insertIntoChild(Vector3f o) { for (int i = 0; i < children.length; i++) { Node c = children[i]; /* * Idea: Test whether ray from origin towards point cuts triangle * * See: http://math.stackexchange.com/questions/1244512/point-in-a-spherical-triangle-test */ if (isPointOnSphericalTriangle(o.x, o.y, o.z, c.v0x, c.v0y, c.v0z, c.v1x, c.v1y, c.v1z, c.v2x, c.v2y, c.v2z, 1E-6f)) { c.insert(o); return; } } } void insert(Vector3f object) { if (children != null) { insertIntoChild(object); return; } if (objects != null && objects.size() == MAX_OBJECTS_PER_NODE) { split(); for (int i = 0; i < MAX_OBJECTS_PER_NODE; i++) insertIntoChild((Vector3f) objects.get(i)); objects = null; insertIntoChild(object); } else { if (objects == null) objects = new ArrayList(MAX_OBJECTS_PER_NODE); objects.add(object); } } /** * This is essentially a ray cast from the origin of the sphere to the point to test and then checking whether that ray goes through the triangle. *

* Reference: Fast, * Minimum Storage Ray/Triangle Intersection */ private static boolean isPointOnSphericalTriangle(float x, float y, float z, float v0X, float v0Y, float v0Z, float v1X, float v1Y, float v1Z, float v2X, float v2Y, float v2Z, float epsilon) { float edge1X = v1X - v0X; float edge1Y = v1Y - v0Y; float edge1Z = v1Z - v0Z; float edge2X = v2X - v0X; float edge2Y = v2Y - v0Y; float edge2Z = v2Z - v0Z; float pvecX = y * edge2Z - z * edge2Y; float pvecY = z * edge2X - x * edge2Z; float pvecZ = x * edge2Y - y * edge2X; float det = edge1X * pvecX + edge1Y * pvecY + edge1Z * pvecZ; if (det > -epsilon && det < epsilon) return false; float tvecX = -v0X; float tvecY = -v0Y; float tvecZ = -v0Z; float invDet = 1.0f / det; float u = (tvecX * pvecX + tvecY * pvecY + tvecZ * pvecZ) * invDet; if (u < 0.0f || u > 1.0f) return false; float qvecX = tvecY * edge1Z - tvecZ * edge1Y; float qvecY = tvecZ * edge1X - tvecX * edge1Z; float qvecZ = tvecX * edge1Y - tvecY * edge1X; float v = (x * qvecX + y * qvecY + z * qvecZ) * invDet; if (v < 0.0f || u + v > 1.0f) return false; float t = (edge2X * qvecX + edge2Y * qvecY + edge2Z * qvecZ) * invDet; return t >= epsilon; } private int child(float x, float y, float z) { for (int i = 0; i < children.length; i++) { Node c = children[i]; if (isPointOnSphericalTriangle(x, y, z, c.v0x, c.v0y, c.v0z, c.v1x, c.v1y, c.v1z, c.v2x, c.v2y, c.v2z, 1E-5f)) return i; } // No child found. This can happen in 'nearest()' when querying possible nearby nodes return 0; } /** * Reference: https://en.wikipedia.org/ */ private float greatCircleDist(float x1, float y1, float z1, float x2, float y2, float z2) { float dot = x1 * x2 + y1 * y2 + z1 * z2; /* * Just use a linear function, because we (mostly) do less-than comparisons on the result. * We just need a linear function which: * f(-1) = PI * f(0) = PI/2 * f(1) = 0 */ return (float) (-Math.PIHalf * dot + Math.PIHalf); //return (float) Math.acos(dot); } float nearest(float x, float y, float z) { return nearest(x, y, z, Float.POSITIVE_INFINITY); } float nearest(float x, float y, float z, float n) { float gcd = greatCircleDist(x, y, z, cx, cy, cz); /* * If great-circle-distance between query point and centroid is larger than the current smallest distance 'n' plus the great circle diameter 'arc', we abort here, * because then it is not possible for any point in the triangle patch to be closer to the query point than 'n'. */ /* * Yes, we are subtracting two great-circle distances from one another here, which we did not even compute correctly * using our overly linear arccos approximation. But the 1.7 factor above will take care that we still stay conservative * enough here and not rejecting triangle patches which would contain samples nearer than 'n'. */ if (gcd - arc > n) return n; float nr = n; if (children != null) { int num = children.length, mod = num-1; for (int i = child(x, y, z), c = 0; c < num; i = (i + 1) & mod, c++) { float n1 = children[i].nearest(x, y, z, nr); nr = Math.min(n1, nr); } return nr; } for (int i = 0; objects != null && i < objects.size(); i++) { Vector3f o = (Vector3f) objects.get(i); float d = greatCircleDist(o.x, o.y, o.z, x, y, z); if (d < nr) nr = d; } return nr; } } private boolean onHemisphere; private int numSamples; private int numCandidates = 60; // <- use a reasonable default private long seed; /** * Create a new instance of {@link Sphere} to configure and generate 'best candidate' sample positions on the unit sphere. */ public Sphere() {} /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xyzs float array. *

* This method performs heap allocations, so should be used sparingly. * * @param xyzs * will hold the x, y and z coordinates of all samples in the order XYZXYZXYZ.... * This array must have a length of at least numSamples * @return this */ public Sphere generate(final float[] xyzs) { final IntHolder i = new IntHolder(); return generate(new Callback3d() { public void onNewSample(float x, float y, float z) { xyzs[3 * i.value + 0] = x; xyzs[3 * i.value + 1] = y; xyzs[3 * i.value + 2] = z; i.value++; } }); } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xyzs FloatBuffer. *

* The samples will be written starting at the current position of the FloatBuffer. The position of the FloatBuffer will not be modified. *

* This method performs heap allocations, so should be used sparingly. * * @param xyzs * will hold the x, y and z coordinates of all samples in the order XYZXYZXYZ.... * This FloatBuffer must have at least numSamples remaining elements. * The position of the buffer will not be modified by this method * @return this */ public Sphere generate(final FloatBuffer xyzs) { final IntHolder i = new IntHolder(); final int pos = xyzs.position(); return generate(new Callback3d() { public void onNewSample(float x, float y, float z) { xyzs.put(pos + 3 * i.value + 0, x); xyzs.put(pos + 3 * i.value + 1, y); xyzs.put(pos + 3 * i.value + 2, z); i.value++; } }); } /** * Set the seed to initialize the pseudo-random number generator with. * * @param seed * the seed value * @return this */ public Sphere seed(long seed) { this.seed = seed; return this; } /** * Set the number of samples to generate. * * @param numSamples * the number of samples * @return this */ public Sphere numSamples(int numSamples) { this.numSamples = numSamples; return this; } /** * Set the number of candidates to try for each generated sample. * * @param numCandidates * the number of candidates to try * @return this */ public Sphere numCandidates(int numCandidates) { this.numCandidates = numCandidates; return this; } /** * Set whether to generate samples on a hemisphere around the +Z axis. *

* The default is false, which will generate samples on the whole unit sphere. * * @param onHemisphere * whether to generate samples on the hemisphere * @return this */ public Sphere onHemisphere(boolean onHemisphere) { this.onHemisphere = onHemisphere; return this; } /** * Generate 'best candidate' sample call the given callback for each generated sample. *

* This method performs heap allocations, so should be used sparingly. * * @param callback * will be called with the coordinates of each generated sample position * @return this */ public Sphere generate(Callback3d callback) { Random rnd = new Random(seed); Node otree = new Node(); for (int i = 0; i < numSamples; i++) { float bestX = Float.NaN, bestY = Float.NaN, bestZ = Float.NaN, bestDist = 0.0f; for (int c = 0; c < numCandidates; c++) { /* * Random point on sphere * * Reference: http://mathworld.wolfram.com/ */ float x1, x2; do { x1 = rnd.nextFloat() * 2.0f - 1.0f; x2 = rnd.nextFloat() * 2.0f - 1.0f; } while (x1 * x1 + x2 * x2 > 1.0f); float sqrt = (float) Math.sqrt(1.0 - x1 * x1 - x2 * x2); float x = 2 * x1 * sqrt; float y = 2 * x2 * sqrt; float z = 1.0f - 2.0f * (x1 * x1 + x2 * x2); if (onHemisphere) { z = Math.abs(z); } float minDist = otree.nearest(x, y, z); if (minDist > bestDist) { bestDist = minDist; bestX = x; bestY = y; bestZ = z; } } callback.onNewSample(bestX, bestY, bestZ); otree.insert(new Vector3f(bestX, bestY, bestZ)); } return this; } } /** * Simple quatree that can handle points and 1-nearest neighbor search. * * @author Kai Burjack */ private static class QuadTree { private static final int MAX_OBJECTS_PER_NODE = 32; // Constants for the quadrants of the quadtree private static final int PXNY = 0; private static final int NXNY = 1; private static final int NXPY = 2; private static final int PXPY = 3; private float minX, minY, hs; private ArrayList objects; private QuadTree[] children; QuadTree(float minX, float minY, float size) { this.minX = minX; this.minY = minY; this.hs = size * 0.5f; } private void split() { children = new QuadTree[4]; children[NXNY] = new QuadTree(minX, minY, hs); children[PXNY] = new QuadTree(minX + hs, minY, hs); children[NXPY] = new QuadTree(minX, minY + hs, hs); children[PXPY] = new QuadTree(minX + hs, minY + hs, hs); } private void insertIntoChild(Vector2f o) { children[quadrant(o.x, o.y)].insert(o); } void insert(Vector2f object) { if (children != null) { insertIntoChild(object); return; } if (objects != null && objects.size() == MAX_OBJECTS_PER_NODE) { split(); for (int i = 0; i < objects.size(); i++) insertIntoChild((Vector2f) objects.get(i)); objects = null; insertIntoChild(object); } else { if (objects == null) objects = new ArrayList(MAX_OBJECTS_PER_NODE); objects.add(object); } } private int quadrant(float x, float y) { if (x < minX + hs) { if (y < minY + hs) return NXNY; return NXPY; } if (y < minY + hs) return PXNY; return PXPY; } float nearest(float x, float y, float lowerBound, float upperBound) { float ub = upperBound; if (x < minX - upperBound || x > minX + hs * 2 + upperBound || y < minY - upperBound || y > minY + hs * 2 + upperBound) return ub; if (children != null) { for (int i = quadrant(x, y), c = 0; c < 4; i = (i + 1) & 3, c++) { float n1 = children[i].nearest(x, y, lowerBound, ub); ub = Math.min(n1, ub); if (ub <= lowerBound) return lowerBound; } return ub; } float ub2 = ub * ub; float lb2 = lowerBound * lowerBound; for (int i = 0; objects != null && i < objects.size(); i++) { Vector2f o = (Vector2f) objects.get(i); float d = o.distanceSquared(x, y); if (d <= lb2) return lowerBound; if (d < ub2) ub2 = d; } return (float) Math.sqrt(ub2); } } /** * Generates Best Candidate samples on a unit disk. * * @author Kai Burjack */ public static class Disk { private int numSamples; private int numCandidates = 60; // <- use a reasonable default private long seed; /** * Create a new instance of {@link Disk} to configure and generate 'best candidate' sample positions on the unit disk. */ public Disk() {} /** * Set the seed to initialize the pseudo-random number generator with. * * @param seed * the seed value * @return this */ public Disk seed(long seed) { this.seed = seed; return this; } /** * Set the number of samples to generate. * * @param numSamples * the number of samples * @return this */ public Disk numSamples(int numSamples) { this.numSamples = numSamples; return this; } /** * Set the number of candidates to try for each generated sample. * * @param numCandidates * the number of candidates to try * @return this */ public Disk numCandidates(int numCandidates) { this.numCandidates = numCandidates; return this; } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xys float array. *

* This method performs heap allocations, so should be used sparingly. * * @param xys * will hold the x and y coordinates of all samples in the order XYXYXY.... * This array must have a length of at least numSamples * @return this */ public Disk generate(final float[] xys) { final IntHolder i = new IntHolder(); return generate(new Callback2d() { public void onNewSample(float x, float y) { xys[2 * i.value + 0] = x; xys[2 * i.value + 1] = y; i.value++; } }); } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xys FloatBuffer. *

* The samples will be written starting at the current position of the FloatBuffer. The position of the FloatBuffer will not be modified. *

* This method performs heap allocations, so should be used sparingly. * * @param xys * will hold the x and y coordinates of all samples in the order XYXYXY.... This FloatBuffer must have at least numSamples remaining elements. The * position of the buffer will not be modified by this method * @return this */ public Disk generate(final FloatBuffer xys) { final IntHolder i = new IntHolder(); final int pos = xys.position(); return generate(new Callback2d() { public void onNewSample(float x, float y) { xys.put(pos + 3 * i.value + 0, x); xys.put(pos + 3 * i.value + 1, y); i.value++; } }); } /** * Generate 'best candidate' sample positions and call the given callback for each generated sample. *

* This method performs heap allocations, so should be used sparingly. * * @param callback * will be called with the coordinates of each generated sample position * @return this */ public Disk generate(Callback2d callback) { QuadTree qtree = new QuadTree(-1, -1, 2); Random rnd = new Random(seed); for (int i = 0; i < numSamples; i++) { float bestX = 0, bestY = 0, bestDist = 0.0f; for (int c = 0; c < numCandidates; c++) { float x, y; do { x = rnd.nextFloat() * 2.0f - 1.0f; y = rnd.nextFloat() * 2.0f - 1.0f; } while (x * x + y * y > 1.0f); float minDist = qtree.nearest(x, y, bestDist, Float.POSITIVE_INFINITY); if (minDist > bestDist) { bestDist = minDist; bestX = x; bestY = y; } } callback.onNewSample(bestX, bestY); qtree.insert(new Vector2f(bestX, bestY)); } return this; } } /** * Generates Best Candidate samples on a unit quad. * * @author Kai Burjack */ public static class Quad { private int numSamples; private int numCandidates = 60; // <- use a reasonable default private long seed; /** * Create a new instance of {@link Quad} to configure and generate 'best candidate' sample positions on the unit quad. */ public Quad() {} /** * Set the seed to initialize the pseudo-random number generator with. * * @param seed * the seed value * @return this */ public Quad seed(long seed) { this.seed = seed; return this; } /** * Set the number of samples to generate. * * @param numSamples * the number of samples * @return this */ public Quad numSamples(int numSamples) { this.numSamples = numSamples; return this; } /** * Set the number of candidates to try for each generated sample. * * @param numCandidates * the number of candidates to try * @return this */ public Quad numCandidates(int numCandidates) { this.numCandidates = numCandidates; return this; } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xyzs float array. *

* This method performs heap allocations, so should be used sparingly. * * @param xyzs * will hold the x, y and z coordinates of all samples in the order XYZXYZXYZ.... * This array must have a length of at least numSamples * @return this */ public Quad generate(final float[] xyzs) { final IntHolder i = new IntHolder(); return generate(new Callback2d() { public void onNewSample(float x, float y) { xyzs[2 * i.value + 0] = x; xyzs[2 * i.value + 1] = y; i.value++; } }); } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xys FloatBuffer. *

* The samples will be written starting at the current position of the FloatBuffer. The position of the FloatBuffer will not be modified. *

* This method performs heap allocations, so should be used sparingly. * * @param xys * will hold the x and y coordinates of all samples in the order XYXYXY.... This FloatBuffer must have at least numSamples remaining elements. The position of * the buffer will not be modified by this method * @return this */ public Quad generate(final FloatBuffer xys) { final IntHolder i = new IntHolder(); final int pos = xys.position(); return generate(new Callback2d() { public void onNewSample(float x, float y) { xys.put(pos + 3 * i.value + 0, x); xys.put(pos + 3 * i.value + 1, y); i.value++; } }); } /** * Generate 'best candidate' sample positions and call the given callback for each generated sample. *

* This method performs heap allocations, so should be used sparingly. * * @param callback * will be called with the coordinates of each generated sample position * @return this */ public Quad generate(Callback2d callback) { QuadTree qtree = new QuadTree(-1, -1, 2); Random rnd = new Random(seed); for (int i = 0; i < numSamples; i++) { float bestX = 0, bestY = 0, bestDist = 0.0f; for (int c = 0; c < numCandidates; c++) { float x = rnd.nextFloat() * 2.0f - 1.0f; float y = rnd.nextFloat() * 2.0f - 1.0f; float minDist = qtree.nearest(x, y, bestDist, Float.POSITIVE_INFINITY); if (minDist > bestDist) { bestDist = minDist; bestX = x; bestY = y; } } callback.onNewSample(bestX, bestY); qtree.insert(new Vector2f(bestX, bestY)); } return this; } } /** * Simple octree for points and 1-nearest neighbor distance query. * * @author Kai Burjack */ private static class Octree { private static final int MAX_OBJECTS_PER_NODE = 32; // Constants for the octants of the octree private static final int PXNYNZ = 0; private static final int NXNYNZ = 1; private static final int NXPYNZ = 2; private static final int PXPYNZ = 3; private static final int PXNYPZ = 4; private static final int NXNYPZ = 5; private static final int NXPYPZ = 6; private static final int PXPYPZ = 7; private float minX, minY, minZ, hs; private ArrayList objects; private Octree[] children; Octree(float minX, float minY, float minZ, float size) { this.minX = minX; this.minY = minY; this.minZ = minZ; this.hs = size * 0.5f; } private void split() { children = new Octree[8]; children[NXNYNZ] = new Octree(minX, minY, minZ, hs); children[PXNYNZ] = new Octree(minX + hs, minY, minZ, hs); children[NXPYNZ] = new Octree(minX, minY + hs, minZ, hs); children[PXPYNZ] = new Octree(minX + hs, minY + hs, minZ, hs); children[NXNYPZ] = new Octree(minX, minY, minZ + hs, hs); children[PXNYPZ] = new Octree(minX + hs, minY, minZ + hs, hs); children[NXPYPZ] = new Octree(minX, minY + hs, minZ + hs, hs); children[PXPYPZ] = new Octree(minX + hs, minY + hs, minZ + hs, hs); } private void insertIntoChild(Vector3f o) { children[octant(o.x, o.y, o.z)].insert(o); } void insert(Vector3f object) { if (children != null) { insertIntoChild(object); return; } if (objects != null && objects.size() == MAX_OBJECTS_PER_NODE) { split(); for (int i = 0; i < objects.size(); i++) insertIntoChild((Vector3f) objects.get(i)); objects = null; insertIntoChild(object); } else { if (objects == null) objects = new ArrayList(MAX_OBJECTS_PER_NODE); objects.add(object); } } private int octant(float x, float y, float z) { if (x < minX + hs) if (y < minY + hs) { if (z < minZ + hs) return NXNYNZ; return NXNYPZ; } else if (z < minZ + hs) return NXPYNZ; else return NXPYPZ; else if (y < minY + hs) { if (z < minZ + hs) return PXNYNZ; return PXNYPZ; } else if (z < minZ + hs) return PXPYNZ; else return PXPYPZ; } float nearest(float x, float y, float z, float lowerBound, float upperBound) { float up = upperBound; if (x < minX - upperBound || x > minX + hs * 2 + upperBound || y < minY - upperBound || y > minY + hs * 2 + upperBound || z < minZ - upperBound || z > minZ + hs * 2 + upperBound) return up; if (children != null) { for (int i = octant(x, y, z), c = 0; c < 8; i = (i + 1) & 7, c++) { float n1 = children[i].nearest(x, y, z, lowerBound, up); up = Math.min(n1, up); if (up <= lowerBound) return lowerBound; } return up; } float up2 = up * up; float lb2 = lowerBound * lowerBound; for (int i = 0; objects != null && i < objects.size(); i++) { Vector3f o = (Vector3f) objects.get(i); float d = o.distanceSquared(x, y, z); if (d <= lb2) return lowerBound; if (d < up2) up2 = d; } return (float) Math.sqrt(up2); } } /** * Generates Best Candidate samples inside a unit cube. * * @author Kai Burjack */ public static class Cube { private int numSamples; private int numCandidates = 60; // <- use a reasonable default private long seed; /** * Create a new instance of {@link Cube} to configure and generate 'best candidate' sample positions * on the unit cube with each sample tried numCandidates number of times, and call the * given callback for each sample generate. */ public Cube() {} /** * Set the seed to initialize the pseudo-random number generator with. * * @param seed * the seed value * @return this */ public Cube seed(long seed) { this.seed = seed; return this; } /** * Set the number of samples to generate. * * @param numSamples * the number of samples * @return this */ public Cube numSamples(int numSamples) { this.numSamples = numSamples; return this; } /** * Set the number of candidates to try for each generated sample. * * @param numCandidates * the number of candidates to try * @return this */ public Cube numCandidates(int numCandidates) { this.numCandidates = numCandidates; return this; } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xyzs float array. *

* This method performs heap allocations, so should be used sparingly. * * @param xyzs * will hold the x, y and z coordinates of all samples in the order XYZXYZXYZ.... * This array must have a length of at least numSamples * @return this */ public Cube generate(final float[] xyzs) { final IntHolder i = new IntHolder(); return generate(new Callback3d() { public void onNewSample(float x, float y, float z) { xyzs[3 * i.value + 0] = x; xyzs[3 * i.value + 1] = y; xyzs[3 * i.value + 2] = z; i.value++; } }); } /** * Generate 'best candidate' sample positions and store the coordinates of all generated samples into the given xyzs FloatBuffer. *

* The samples will be written starting at the current position of the FloatBuffer. The position of the FloatBuffer will not be modified. *

* This method performs heap allocations, so should be used sparingly. * * @param xyzs * will hold the x, y and z coordinates of all samples in the order XYZXYZXYZ.... * This FloatBuffer must have at least numSamples remaining elements. * The position of the buffer will not be modified by this method * @return this */ public Cube generate(final FloatBuffer xyzs) { final IntHolder i = new IntHolder(); final int pos = xyzs.position(); return generate(new Callback3d() { public void onNewSample(float x, float y, float z) { xyzs.put(pos + 3 * i.value + 0, x); xyzs.put(pos + 3 * i.value + 1, y); xyzs.put(pos + 3 * i.value + 2, z); i.value++; } }); } /** * Generate 'best candidate' sample positions and call the given callback for each generated sample. *

* This method performs heap allocations, so should be used sparingly. * * @param callback * will be called with the coordinates of each generated sample position * @return this */ public Cube generate(Callback3d callback) { Octree octree = new Octree(-1, -1, -1, 2); Random rnd = new Random(seed); for (int i = 0; i < numSamples; i++) { float bestX = 0, bestY = 0, bestZ = 0, bestDist = 0.0f; for (int c = 0; c < numCandidates; c++) { float x = rnd.nextFloat() * 2.0f - 1.0f; float y = rnd.nextFloat() * 2.0f - 1.0f; float z = rnd.nextFloat() * 2.0f - 1.0f; float minDist = octree.nearest(x, y, z, bestDist, Float.POSITIVE_INFINITY); if (minDist > bestDist) { bestDist = minDist; bestX = x; bestY = y; bestZ = z; } } callback.onNewSample(bestX, bestY, bestZ); octree.insert(new Vector3f(bestX, bestY, bestZ)); } return this; } } }