给定 3d space 中的一组点,找到彼此距离内的所有点集

Given a set of points in 3d space, find all sets of points within a distance of eachother

我有一组 3d 点 S。 我需要在曼哈顿距离 d 中找到 S 中所有点集的集合 X彼此。

即对于 X 中的每个集合 Y 在 3d space 中至少存在一个点在距离 d[= 内Y

中所有点的 38=]

S 的长度永远不会 >20 但我将不得不 运行 对正在生产的集流进行此分析 ~10 新每秒设置,所以无论我使用什么解决方案都必须相当有效。

帮助可视化问题的示例,给出以下内容:

输出将是 ((A,B), (B,C,E), (B,D,E))

我们只关心最大可能的集合,所以集合 (B,C)、(B,D)、(B,E)、(C,E) 和 (D,E),同时在给定的范围内参数,不在输出中,因为它们是 X

中其他集合的子集

我也在 java 中这样做,但是任何关于算法或伪代码的指示将不胜感激,在此先感谢。

20 个点之间的详尽距离计算,即 190 个距离对于 PC 来说毫无意义。时间将以微秒为单位。您可以从矩阵中编码的 "close to" 关系中提取所需信息。

伪代码的解决方案是:

calculate_intersections(areas):
    intersections = calculate every two intersecting areas
    combinations = combine_intersections(intersections)
    reduced = remove all sets in combinations that are included in bigger sets

combine_intersections(intersections):
    do:
        combinations = new HashSet

        for s1 in intersections:
            for s2 in intersections:
                diff_1_2 = s1 \ s2
                diff_2_1 = s2 \ s1

                if diff_1_2.len == 1 && diff_2_1.len == 1:
                    union = diff_1_2 + diff_2_1

                    if union in intersections:
                        union2 = s1 + s2
                        if !union2 in intersections:
                            combinations.add(union)
    while (combinations not empty)

Java 中的实现可能如下所示:

import java.util.Arrays;
import java.util.HashSet;
import java.util.Iterator;
import java.util.Set;

import org.apache.commons.collections4.SetUtils;

public class IntersectionSetCalculation {

    private static class ManhattanDistanceArea {

        private String id;
        private Vector3D center;
        private double distance;

        public ManhattanDistanceArea(Vector3D center, double distance, String id) {
            this.center = center;
            this.distance = distance;
            this.id = id;
        }

        @Override
        public String toString() {
            return id;
        }

        @Override
        public int hashCode() {
            final int prime = 31;
            int result = 1;
            result = prime * result + ((center == null) ? 0 : center.hashCode());
            long temp;
            temp = Double.doubleToLongBits(distance);
            result = prime * result + (int) (temp ^ (temp >>> 32));
            result = prime * result + ((id == null) ? 0 : id.hashCode());
            return result;
        }

        @Override
        public boolean equals(Object obj) {
            if (this == obj)
                return true;
            if (obj == null)
                return false;
            if (getClass() != obj.getClass())
                return false;
            ManhattanDistanceArea other = (ManhattanDistanceArea) obj;
            if (center == null) {
                if (other.center != null)
                    return false;
            }
            else if (!center.equals(other.center))
                return false;
            if (Double.doubleToLongBits(distance) != Double.doubleToLongBits(other.distance))
                return false;
            if (id == null) {
                if (other.id != null)
                    return false;
            }
            else if (!id.equals(other.id))
                return false;
            return true;
        }

        public boolean intersects(ManhattanDistanceArea other) {
            double maxDist = distance + other.distance;
            return center.distance(other.center, 1) < maxDist;
        }
    }

    /**
     * Calculate the intersection of all areas (maximum of 2 areas in an intersection)
     */
    public static Set<Set<ManhattanDistanceArea>> getIntersectingAreas(Set<ManhattanDistanceArea> areas) {
        Set<Set<ManhattanDistanceArea>> intersections = new HashSet<Set<ManhattanDistanceArea>>();

        for (ManhattanDistanceArea area : areas) {
            for (ManhattanDistanceArea area2 : areas) {
                if (!area.equals(area2) && area.intersects(area2)) {
                    HashSet<ManhattanDistanceArea> intersection = new HashSet<ManhattanDistanceArea>();
                    intersection.add(area);
                    intersection.add(area2);
                    intersections.add(intersection);
                }
            }
        }

        Set<Set<ManhattanDistanceArea>> combined = combineIntersections(intersections);
        Set<Set<ManhattanDistanceArea>> reduced = reduceIntersections(combined);

        return reduced;
    }

    /**
     * Combine the small intersections (with a maximum of 2 areas in an intersection) to bigger intersections
     */
    public static Set<Set<ManhattanDistanceArea>> combineIntersections(Set<Set<ManhattanDistanceArea>> inters) {
        Set<Set<ManhattanDistanceArea>> intersections = new HashSet<Set<ManhattanDistanceArea>>(inters);
        Set<Set<ManhattanDistanceArea>> combinations;

        do {
            combinations = new HashSet<Set<ManhattanDistanceArea>>();

            for (Set<ManhattanDistanceArea> intersecting1 : intersections) {
                for (Set<ManhattanDistanceArea> intersecting2 : intersections) {
                    Set<ManhattanDistanceArea> diff_1_2 = SetUtils.difference(intersecting1, intersecting2);
                    Set<ManhattanDistanceArea> diff_2_1 = SetUtils.difference(intersecting2, intersecting1);

                    if (diff_1_2.size() == 1 && diff_2_1.size() == 1) {
                        Set<ManhattanDistanceArea> union_1_2 = SetUtils.union(diff_1_2, diff_2_1);

                        if (intersections.contains(union_1_2)) {
                            Set<ManhattanDistanceArea> union = SetUtils.union(intersecting1, intersecting2);
                            if (!intersections.contains(union)) {
                                combinations.add(union);
                            }
                        }
                    }
                }
            }

            intersections.addAll(combinations);

        } while (!combinations.isEmpty());

        return intersections;
    }

    /**
     * Remove the small intersections that are completely covered by bigger intersections
     */
    public static Set<Set<ManhattanDistanceArea>> reduceIntersections(Set<Set<ManhattanDistanceArea>> inters) {
        Set<Set<ManhattanDistanceArea>> intersections = new HashSet<Set<ManhattanDistanceArea>>(inters);
        Iterator<Set<ManhattanDistanceArea>> iter = intersections.iterator();

        while (iter.hasNext()) {
            Set<ManhattanDistanceArea> intersection = iter.next();
            for (Set<ManhattanDistanceArea> intersection2 : inters) {
                if (intersection2.size() > intersection.size() && intersection2.containsAll(intersection)) {
                    iter.remove();
                    break;
                }
            }
        }

        return intersections;
    }

    public static void main(String[] args) {
        final double dist = 2d;//the manhattan distance d

        ManhattanDistanceArea A = new ManhattanDistanceArea(new Vector3D(0, -3, 0), dist, "A");
        ManhattanDistanceArea B = new ManhattanDistanceArea(new Vector3D(0, 0, 0), dist, "B");
        ManhattanDistanceArea C = new ManhattanDistanceArea(new Vector3D(3.5, 0, 0), dist, "C");
        ManhattanDistanceArea D = new ManhattanDistanceArea(new Vector3D(0, 3.5, 0), dist, "D");
        ManhattanDistanceArea E = new ManhattanDistanceArea(new Vector3D(1, 1, 0), dist, "E");

        ManhattanDistanceArea F = new ManhattanDistanceArea(new Vector3D(-1, 1, 0), dist, "F");

        //test the example you provided
        Set<ManhattanDistanceArea> abcde = new HashSet<ManhattanDistanceArea>();
        abcde.addAll(Arrays.asList(new ManhattanDistanceArea[] {A, B, C, D, E}));

        //test another example
        Set<ManhattanDistanceArea> abcdef = new HashSet<ManhattanDistanceArea>();
        abcdef.addAll(abcde);
        abcdef.add(F);

        Set<Set<ManhattanDistanceArea>> intersectionsABCDE = getIntersectingAreas(abcde);
        Set<Set<ManhattanDistanceArea>> intersectionsABCDEF = getIntersectingAreas(abcdef);

        System.out.println(intersectionsABCDE);
        System.out.println(intersectionsABCDEF);

        //test the runntime for 1000 calculation
        double startTime = System.currentTimeMillis();
        final int calculations = 1000;
        for (int i = 0; i < calculations; i++) {
            Set<ManhattanDistanceArea> areas = new HashSet<ManhattanDistanceArea>();
            for (int j = 0; j < 20; j++) {
                areas.add(new ManhattanDistanceArea(new Vector3D(Math.random() * 10 - 5, Math.random() * 10 - 5, Math.random() * 10 - 5), dist,
                        "A" + j));
            }

            getIntersectingAreas(areas);
        }
        System.out.println("\nTime used for " + calculations + " intersection calculations (with sets of size 20): "
                + (System.currentTimeMillis() - startTime) + "ms");
    }
}

为了实现,我使用了这个 class Vector3D:

public class Vector3D {

    public double x;
    public double y;
    public double z;

    public static final Vector3D NAN_VEC = new Vector3D(Double.NaN, Double.NaN, Double.NaN);
    public static final Vector3D NULL_VEC = new Vector3D(0, 0, 0);

    public enum Axis {
        X, Y, Z;
    }

    public Vector3D() {

    }

    /**
     * Crate a new Vector2D with x and y components.
     */
    public Vector3D(double x, double y, double z) {
        this.x = x;
        this.y = y;
        this.z = z;
    }
    public Vector3D(double... val) {
        x = val[0];
        y = val[1];
        z = val[2];
    }

    /**
     * Create a Vector3D by two angles (in degree).
     * 
     * The first angle is in XY direction. The second angle is the Z direction.
     * 
     * An angle (XY) of 0° results in (x, y) = (1, 0); 90° in (x, y) = (0, 1); ... An angle (Z) of 0° results in (x, y, z) = (x, y, 0); 90° in (x, y,
     * z) = (x, y, 1); -90° in (x, y, z) = (x, y, -1)
     * 
     * The resulting vector has a length of 1.
     * 
     * @param angleXY
     *        The angle of the new vector (in degree) for the XY direction (from 0 to 360).
     * 
     * @param angleZ
     *        The angle of the new vector (in degree) for the Z direction (from -90 to 90).
     */
    public Vector3D(double angleXY, double angleZ) {
        x = Math.cos(angleXY * Math.PI / 180) * Math.cos(angleZ * Math.PI / 180);
        y = Math.sin(angleXY * Math.PI / 180) * Math.cos(angleZ * Math.PI / 180);
        z = Math.sin(angleZ * Math.PI / 180);
        double len = length();
        x /= len;
        y /= len;
        z /= len;
    }

    private Vector3D(Vector3D clone) {
        this.x = clone.x;
        this.y = clone.y;
    }

    @Override
    public Vector3D clone() {
        return new Vector3D(this);
    }

    @Override
    public String toString() {
        return "Vector3D[x: " + x + " y: " + y + " z:" + z + "]";
    }

    @Override
    public boolean equals(Object obj) {
        if (obj instanceof Vector3D) {
            Vector3D v = (Vector3D) obj;
            return Math.abs(x - v.x) < 1e-8 && Math.abs(y - v.y) < 1e-8 && Math.abs(z - v.z) < 1e-8;
        }
        return false;
    }

    /**
     * Get this vector as 3D-Array.
     */
    public double[] asArray() {
        return new double[] {x, y, z};
    }

    /**
     * The (euclidean) length of the Vector.
     */
    public double length() {
        return Math.sqrt(Math.pow(x, 2) + Math.pow(y, 2) + Math.pow(z, 2));
    }
    /**
     * The length of this vector in a given norm.
     * 
     * @param norm
     *        The norm of the vector length.
     * 
     * @return The length of this vector in the given norm.
     */
    public double length(int norm) {
        if (norm == Integer.MAX_VALUE) {
            return Math.max(Math.max(x, y), z);
        }
        return Math.pow(Math.pow(x, norm) + Math.pow(y, norm) + Math.pow(z, norm), 1.0 / norm);
    }

    /**
     * Rotate this vector an angle (in degrees) around an axis resulting in a new Vector that is returned.
     * 
     * @param degrees
     *        The angle to return the vector.
     * 
     * @param axis
     *        The axis around which the vector is rotated.
     * 
     * @return The new created vector.
     */
    public Vector3D rotate(double degrees, Axis axis) {
        double cos = Math.cos(degrees * Math.PI / 180);
        double sin = Math.sin(degrees * Math.PI / 180);
        switch (axis) {
            case X:
                return new Vector3D(x, cos * y - sin * z, sin * y + cos * z);
            case Y:
                return new Vector3D(cos * x + sin * z, y, -sin * x + cos * z);
            case Z:
                return new Vector3D(cos * x - sin * y, sin * x + cos * y, z);
            default:
                return null;
        }
    }

    /**
     * Project the vector given as parameter on this vector.
     * 
     * @param vec
     *        The vector that is to be projected on this vector.
     * 
     * @return The projected vector.
     */
    public Vector3D project(Vector3D vec) {
        return mult(scalar(vec) / Math.pow(length(), 2));
    }

    /**
     * Add another Vector3D to this vector resulting in a new Vector that is returned.
     * 
     * @param vec
     *        The vector added to this vector.
     * 
     * @return The new created vector.
     */
    public Vector3D add(Vector3D vec) {
        return new Vector3D(x + vec.x, y + vec.y, z + vec.z);
    }
    /**
     * Subtract another Vector3D from this vector resulting in a new Vector that is returned.
     * 
     * @param vec
     *        The vector subtracted from this vector.
     * 
     * @return The new created vector.
     */
    public Vector3D sub(Vector3D vec) {
        return new Vector3D(x - vec.x, y - vec.y, z - vec.z);
    }
    /**
     * Multiply this vector with a scalar resulting in a new Vector that is returned.
     * 
     * @param scalar
     *        The scalar to multiply this vector with.
     * 
     * @return The new created vector.
     */
    public Vector3D mult(double scalar) {
        return new Vector3D(x * scalar, y * scalar, z * scalar);
    }

    /**
     * Check whether this vector is linearly dependent to the parameter vector.
     * 
     * @param vec
     *        The checked vector.
     * 
     * @return True if the vectors are linearly dependent. False otherwise.
     */
    public boolean isLinearlyDependent(Vector3D vec) {
        double t1 = (x == 0 ? 0 : vec.x / x);
        double t2 = (y == 0 ? 0 : vec.y / y);
        double t3 = (z == 0 ? 0 : vec.z / z);
        return Math.abs(t1 - t2) < 1e-5 && Math.abs(t1 - t3) < 1e-5 && t1 != 0;//all parameters t are equal and != 0
    }

    /**
     * Calculate the scalar product of this vector and the parameter vector.
     * 
     * @param vec
     *        The vector to calculate the scalar with this vector.
     * 
     * @return The scalar of the vectors.
     */
    public double scalar(Vector3D vec) {
        return this.x * vec.x + this.y * vec.y + this.z * vec.z;
    }

    /**
     * Calculate the cross product of this vector with another vector (resulting vector = this X parameter vector)
     * 
     * @param vec
     *        The second vector for the cross product calculation.
     * 
     * @return The cross product vector of the two vectors.
     */
    public Vector3D cross(Vector3D vec) {
        return new Vector3D(y * vec.z - z * vec.y, z * vec.x - x * vec.z, x * vec.y - y * vec.x);
    }

    /**
     * Create a new vector with the same direction but a different length as this vector.
     * 
     * @param length
     *        The length of the new vector.
     * 
     * @return The new vector with a new length.
     */
    public Vector3D setLength(double length) {
        double len = length();
        return new Vector3D(x * length / len, y * length / len, z * length / len);
    }

    /**
     * Get the distance of this point's position vector to another point's position vector.
     * 
     * @param p
     *        The second point's position vector.
     * 
     * @return The distance between the points.
     */
    public double distance(Vector3D p) {
        return Math.sqrt((this.x - p.x) * (this.x - p.x) + (this.y - p.y) * (this.y - p.y) + (this.z - p.z) * (this.z - p.z));
    }
    /**
     * Get the distance of this point's position vector to another point's position vector in a given norm.
     * 
     * @param p
     *        The second point's position vector.
     * 
     * @param norm
     *        The norm in which the distance is calculated (1 -> manhattan, 2 -> euclide, ...)
     * 
     * @return The distance between the points in the given norm.
     */
    public double distance(Vector3D p, int norm) {
        return Math.pow((Math.pow(Math.abs(this.x - p.x), norm) + Math.pow(Math.abs(this.y - p.y), norm) + Math.pow(Math.abs(this.z - p.z), norm)),
                1d / norm);
    }

    /**
     * Change this vector to the new coordinates.
     */
    public void move(double x, double y, double z) {
        this.x = x;
        this.y = y;
        this.z = z;
    }

    /**
     * Move a point's position vector in a direction (by a vector) and a distance.
     * 
     * @param p
     *        The direction vector.
     * 
     * @param distance
     *        The distance to move the new vector
     * 
     * @return The new created vector.
     */
    public Vector3D moveTo(Vector3D p, double distance) {
        double d = distance(p);
        double dx = p.x - x;
        double dy = p.y - y;
        double dz = p.z - z;
        double coef = distance / d;
        return new Vector3D(x + dx * coef, y + dy * coef, z + dz * coef);
    }

    /**
     * Get the angle difference of this vector to another vector.
     * 
     * @param vec
     *        The other vector.
     * 
     * @return The angle difference of the two vectors (from 0° to 180°).
     */
    public double getAngleTo(Vector3D vec) {
        double angle = Math.acos(scalar(vec) / (length() * vec.length())) * 180 / Math.PI;
        if (angle > 180) {
            angle = 360 - angle;
        }
        return angle;
    }

    /**
     * Get the vector from this point to another.
     * 
     * @param vec
     *        The point to which the vector is calculated.
     * 
     * @return The vector from this points position vector to the other point.
     */
    public Vector3D vectorTo(Vector3D vec) {
        return new Vector3D(vec.x - x, vec.y - y, vec.z - z);
    }

    /**
     * Checks whether a point (by its position vector) is in a given range of this point.
     * 
     * @param p
     *        The point that is checked.
     * 
     * @param range
     *        The range used for the check.
     * 
     * @return True if the point is in the range of this point (distance <= range).
     */
    public boolean isInRange(Vector3D p, double range) {
        return p != this && distance(p) <= range;
    }
}

和来自 apache 公共库的 class SetUtils

我还添加了一些测试:

  1. 你问题的测试
  2. 另一个具有更大交叉集的测试
  3. 运行时测试

结果是:

[[A, B], [B, E, C], [B, E, D]]

[[A, B], [B, E, C], [D, E, F, B]]

Time used for 1000 intersection calculations (with sets of size 20): 791.0ms

所以结果似乎是正确的,你可以在一秒钟内计算出 1000 多个交叉点。