OrderedTuple.java

  1. /*
  2.  * Licensed to the Apache Software Foundation (ASF) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * The ASF licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *      http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.apache.commons.math3.geometry.partitioning.utilities;

  19. import java.util.Arrays;

  21. import org.apache.commons.math3.util.FastMath;

  23. /** This class implements an ordering operation for T-uples.
  24.  *
  25.  * <p>Ordering is done by encoding all components of the T-uple into a
  26.  * single scalar value and using this value as the sorting
  27.  * key. Encoding is performed using the method invented by Georg
  28.  * Cantor in 1877 when he proved it was possible to establish a
  29.  * bijection between a line and a plane. The binary representations of
  30.  * the components of the T-uple are mixed together to form a single
  31.  * scalar. This means that the 2<sup>k</sup> bit of component 0 is
  32.  * followed by the 2<sup>k</sup> bit of component 1, then by the
  33.  * 2<sup>k</sup> bit of component 2 up to the 2<sup>k</sup> bit of
  34.  * component {@code t}, which is followed by the 2<sup>k-1</sup>
  35.  * bit of component 0, followed by the 2<sup>k-1</sup> bit of
  36.  * component 1 ... The binary representations are extended as needed
  37.  * to handle numbers with different scales and a suitable
  38.  * 2<sup>p</sup> offset is added to the components in order to avoid
  39.  * negative numbers (this offset is adjusted as needed during the
  40.  * comparison operations).</p>
  41.  *
  42.  * <p>The more interesting property of the encoding method for our
  43.  * purpose is that it allows to select all the points that are in a
  44.  * given range. This is depicted in dimension 2 by the following
  45.  * picture:</p>
  46.  *
  47.  * <img src="doc-files/OrderedTuple.png" />
  48.  *
  49.  * <p>This picture shows a set of 100000 random 2-D pairs having their
  50.  * first component between -50 and +150 and their second component
  51.  * between -350 and +50. We wanted to extract all pairs having their
  52.  * first component between +30 and +70 and their second component
  53.  * between -120 and -30. We built the lower left point at coordinates
  54.  * (30, -120) and the upper right point at coordinates (70, -30). All
  55.  * points smaller than the lower left point are drawn in red and all
  56.  * points larger than the upper right point are drawn in blue. The
  57.  * green points are between the two limits. This picture shows that
  58.  * all the desired points are selected, along with spurious points. In
  59.  * this case, we get 15790 points, 4420 of which really belonging to
  60.  * the desired rectangle. It is possible to extract very small
  61.  * subsets. As an example extracting from the same 100000 points set
  62.  * the points having their first component between +30 and +31 and
  63.  * their second component between -91 and -90, we get a subset of 11
  64.  * points, 2 of which really belonging to the desired rectangle.</p>
  65.  *
  66.  * <p>the previous selection technique can be applied in all
  67.  * dimensions, still using two points to define the interval. The
  68.  * first point will have all its components set to their lower bounds
  69.  * while the second point will have all its components set to their
  70.  * upper bounds.</p>
  71.  *
  72.  * <p>T-uples with negative infinite or positive infinite components
  73.  * are sorted logically.</p>
  74.  *
  75.  * <p>Since the specification of the {@code Comparator} interface
  76.  * allows only {@code ClassCastException} errors, some arbitrary
  77.  * choices have been made to handle specific cases. The rationale for
  78.  * these choices is to keep <em>regular</em> and consistent T-uples
  79.  * together.</p>
  80.  * <ul>
  81.  * <li>instances with different dimensions are sorted according to
  82.  * their dimension regardless of their components values</li>
  83.  * <li>instances with {@code Double.NaN} components are sorted
  84.  * after all other ones (even after instances with positive infinite
  85.  * components</li>
  86.  * <li>instances with both positive and negative infinite components
  87.  * are considered as if they had {@code Double.NaN}
  88.  * components</li>
  89.  * </ul>
  90.  *
  91.  * @since 3.0
  92.  * @deprecated as of 3.4, this class is not used anymore and considered
  93.  * to be out of scope of Apache Commons Math
  94.  */
  95. @Deprecated
  96. public class OrderedTuple implements Comparable<OrderedTuple> {

  98.     /** Sign bit mask. */
  99.     private static final long SIGN_MASK     = 0x8000000000000000L;

  101.     /** Exponent bits mask. */
  102.     private static final long EXPONENT_MASK = 0x7ff0000000000000L;

  104.     /** Mantissa bits mask. */
  105.     private static final long MANTISSA_MASK = 0x000fffffffffffffL;

  107.     /** Implicit MSB for normalized numbers. */
  108.     private static final long IMPLICIT_ONE  = 0x0010000000000000L;

  110.     /** Double components of the T-uple. */
  111.     private double[] components;

  113.     /** Offset scale. */
  114.     private int offset;

  116.     /** Least Significant Bit scale. */
  117.     private int lsb;

  119.     /** Ordering encoding of the double components. */
  120.     private long[] encoding;

  122.     /** Positive infinity marker. */
  123.     private boolean posInf;

  125.     /** Negative infinity marker. */
  126.     private boolean negInf;

  128.     /** Not A Number marker. */
  129.     private boolean nan;

  131.     /** Build an ordered T-uple from its components.
  132.      * @param components double components of the T-uple
  133.      */
  134.     public OrderedTuple(final double ... components) {
  135.         this.components = components.clone();
  136.         int msb = Integer.MIN_VALUE;
  137.         lsb     = Integer.MAX_VALUE;
  138.         posInf  = false;
  139.         negInf  = false;
  140.         nan     = false;
  141.         for (int i = 0; i < components.length; ++i) {
  142.             if (Double.isInfinite(components[i])) {
  143.                 if (components[i] < 0) {
  144.                     negInf = true;
  145.                 } else {
  146.                     posInf = true;
  147.                 }
  148.             } else if (Double.isNaN(components[i])) {
  149.                 nan = true;
  150.             } else {
  151.                 final long b = Double.doubleToLongBits(components[i]);
  152.                 final long m = mantissa(b);
  153.                 if (m != 0) {
  154.                     final int e = exponent(b);
  155.                     msb = FastMath.max(msb, e + computeMSB(m));
  156.                     lsb = FastMath.min(lsb, e + computeLSB(m));
  157.                 }
  158.             }
  159.         }

  161.         if (posInf && negInf) {
  162.             // instance cannot be sorted logically
  163.             posInf = false;
  164.             negInf = false;
  165.             nan    = true;
  166.         }

  168.         if (lsb <= msb) {
  169.             // encode the T-upple with the specified offset
  170.             encode(msb + 16);
  171.         } else {
  172.             encoding = new long[] {
  173.                 0x0L
  174.             };
  175.         }

  177.     }

  179.     /** Encode the T-uple with a given offset.
  180.      * @param minOffset minimal scale of the offset to add to all
  181.      * components (must be greater than the MSBs of all components)
  182.      */
  183.     private void encode(final int minOffset) {

  185.         // choose an offset with some margins
  186.         offset  = minOffset + 31;
  187.         offset -= offset % 32;

  189.         if ((encoding != null) && (encoding.length == 1) && (encoding[0] == 0x0L)) {
  190.             // the components are all zeroes
  191.             return;
  192.         }

  194.         // allocate an integer array to encode the components (we use only
  195.         // 63 bits per element because there is no unsigned long in Java)
  196.         final int neededBits  = offset + 1 - lsb;
  197.         final int neededLongs = (neededBits + 62) / 63;
  198.         encoding = new long[components.length * neededLongs];

  200.         // mix the bits from all components
  201.         int  eIndex = 0;
  202.         int  shift  = 62;
  203.         long word   = 0x0L;
  204.         for (int k = offset; eIndex < encoding.length; --k) {
  205.             for (int vIndex = 0; vIndex < components.length; ++vIndex) {
  206.                 if (getBit(vIndex, k) != 0) {
  207.                     word |= 0x1L << shift;
  208.                 }
  209.                 if (shift-- == 0) {
  210.                     encoding[eIndex++] = word;
  211.                     word  = 0x0L;
  212.                     shift = 62;
  213.                 }
  214.             }
  215.         }

  217.     }

  219.     /** Compares this ordered T-uple with the specified object.

  221.      * <p>The ordering method is detailed in the general description of
  222.      * the class. Its main property is to be consistent with distance:
  223.      * geometrically close T-uples stay close to each other when stored
  224.      * in a sorted collection using this comparison method.</p>

  226.      * <p>T-uples with negative infinite, positive infinite are sorted
  227.      * logically.</p>

  229.      * <p>Some arbitrary choices have been made to handle specific
  230.      * cases. The rationale for these choices is to keep
  231.      * <em>normal</em> and consistent T-uples together.</p>
  232.      * <ul>
  233.      * <li>instances with different dimensions are sorted according to
  234.      * their dimension regardless of their components values</li>
  235.      * <li>instances with {@code Double.NaN} components are sorted
  236.      * after all other ones (evan after instances with positive infinite
  237.      * components</li>
  238.      * <li>instances with both positive and negative infinite components
  239.      * are considered as if they had {@code Double.NaN}
  240.      * components</li>
  241.      * </ul>

  243.      * @param ot T-uple to compare instance with
  244.      * @return a negative integer if the instance is less than the
  245.      * object, zero if they are equal, or a positive integer if the
  246.      * instance is greater than the object

  248.      */
  249.     public int compareTo(final OrderedTuple ot) {
  250.         if (components.length == ot.components.length) {
  251.             if (nan) {
  252.                 return +1;
  253.             } else if (ot.nan) {
  254.                 return -1;
  255.             } else if (negInf || ot.posInf) {
  256.                 return -1;
  257.             } else if (posInf || ot.negInf) {
  258.                 return +1;
  259.             } else {

  261.                 if (offset < ot.offset) {
  262.                     encode(ot.offset);
  263.                 } else if (offset > ot.offset) {
  264.                     ot.encode(offset);
  265.                 }

  267.                 final int limit = FastMath.min(encoding.length, ot.encoding.length);
  268.                 for (int i = 0; i < limit; ++i) {
  269.                     if (encoding[i] < ot.encoding[i]) {
  270.                         return -1;
  271.                     } else if (encoding[i] > ot.encoding[i]) {
  272.                         return +1;
  273.                     }
  274.                 }

  276.                 if (encoding.length < ot.encoding.length) {
  277.                     return -1;
  278.                 } else if (encoding.length > ot.encoding.length) {
  279.                     return +1;
  280.                 } else {
  281.                     return 0;
  282.                 }

  284.             }
  285.         }

  287.         return components.length - ot.components.length;

  289.     }

  291.     /** {@inheritDoc} */
  292.     @Override
  293.     public boolean equals(final Object other) {
  294.         if (this == other) {
  295.             return true;
  296.         } else if (other instanceof OrderedTuple) {
  297.             return compareTo((OrderedTuple) other) == 0;
  298.         } else {
  299.             return false;
  300.         }
  301.     }

  303.     /** {@inheritDoc} */
  304.     @Override
  305.     public int hashCode() {
  306.         // the following constants are arbitrary small primes
  307.         final int multiplier = 37;
  308.         final int trueHash   = 97;
  309.         final int falseHash  = 71;

  311.         // hash fields and combine them
  312.         // (we rely on the multiplier to have different combined weights
  313.         //  for all int fields and all boolean fields)
  314.         int hash = Arrays.hashCode(components);
  315.         hash = hash * multiplier + offset;
  316.         hash = hash * multiplier + lsb;
  317.         hash = hash * multiplier + (posInf ? trueHash : falseHash);
  318.         hash = hash * multiplier + (negInf ? trueHash : falseHash);
  319.         hash = hash * multiplier + (nan    ? trueHash : falseHash);

  321.         return hash;

  323.     }

  325.     /** Get the components array.
  326.      * @return array containing the T-uple components
  327.      */
  328.     public double[] getComponents() {
  329.         return components.clone();
  330.     }

  332.     /** Extract the sign from the bits of a double.
  333.      * @param bits binary representation of the double
  334.      * @return sign bit (zero if positive, non zero if negative)
  335.      */
  336.     private static long sign(final long bits) {
  337.         return bits & SIGN_MASK;
  338.     }

  340.     /** Extract the exponent from the bits of a double.
  341.      * @param bits binary representation of the double
  342.      * @return exponent
  343.      */
  344.     private static int exponent(final long bits) {
  345.         return ((int) ((bits & EXPONENT_MASK) >> 52)) - 1075;
  346.     }

  348.     /** Extract the mantissa from the bits of a double.
  349.      * @param bits binary representation of the double
  350.      * @return mantissa
  351.      */
  352.     private static long mantissa(final long bits) {
  353.         return ((bits & EXPONENT_MASK) == 0) ?
  354.                ((bits & MANTISSA_MASK) << 1) :          // subnormal number
  355.                (IMPLICIT_ONE | (bits & MANTISSA_MASK)); // normal number
  356.     }

  358.     /** Compute the most significant bit of a long.
  359.      * @param l long from which the most significant bit is requested
  360.      * @return scale of the most significant bit of {@code l},
  361.      * or 0 if {@code l} is zero
  362.      * @see #computeLSB
  363.      */
  364.     private static int computeMSB(final long l) {

  366.         long ll = l;
  367.         long mask  = 0xffffffffL;
  368.         int  scale = 32;
  369.         int  msb   = 0;

  371.         while (scale != 0) {
  372.             if ((ll & mask) != ll) {
  373.                 msb |= scale;
  374.                 ll >>= scale;
  375.             }
  376.             scale >>= 1;
  377.             mask >>= scale;
  378.         }

  380.         return msb;

  382.     }

  384.     /** Compute the least significant bit of a long.
  385.      * @param l long from which the least significant bit is requested
  386.      * @return scale of the least significant bit of {@code l},
  387.      * or 63 if {@code l} is zero
  388.      * @see #computeMSB
  389.      */
  390.     private static int computeLSB(final long l) {

  392.         long ll = l;
  393.         long mask  = 0xffffffff00000000L;
  394.         int  scale = 32;
  395.         int  lsb   = 0;

  397.         while (scale != 0) {
  398.             if ((ll & mask) == ll) {
  399.                 lsb |= scale;
  400.                 ll >>= scale;
  401.             }
  402.             scale >>= 1;
  403.             mask >>= scale;
  404.         }

  406.         return lsb;

  408.     }

  410.     /** Get a bit from the mantissa of a double.
  411.      * @param i index of the component
  412.      * @param k scale of the requested bit
  413.      * @return the specified bit (either 0 or 1), after the offset has
  414.      * been added to the double
  415.      */
  416.     private int getBit(final int i, final int k) {
  417.         final long bits = Double.doubleToLongBits(components[i]);
  418.         final int e = exponent(bits);
  419.         if ((k < e) || (k > offset)) {
  420.             return 0;
  421.         } else if (k == offset) {
  422.             return (sign(bits) == 0L) ? 1 : 0;
  423.         } else if (k > (e + 52)) {
  424.             return (sign(bits) == 0L) ? 0 : 1;
  425.         } else {
  426.             final long m = (sign(bits) == 0L) ? mantissa(bits) : -mantissa(bits);
  427.             return (int) ((m >> (k - e)) & 0x1L);
  428.         }
  429.     }

  431. }
Created with JaCoCo 0.7.2.201409121644