Mercurial > dive4elements > river
view artifacts/src/main/java/org/dive4elements/river/utils/DoubleUtil.java @ 9006:7134a4c7d1b6
U-Info S-Info general work
author | gernotbelger |
---|---|
date | Tue, 17 Apr 2018 13:20:25 +0200 |
parents | 5e38e2924c07 |
children | 8aa1605e58f1 0a5239a1e46e |
line wrap: on
line source
/* Copyright (C) 2011, 2012, 2013 by Bundesanstalt für Gewässerkunde * Software engineering by Intevation GmbH * * This file is Free Software under the GNU AGPL (>=v3) * and comes with ABSOLUTELY NO WARRANTY! Check out the * documentation coming with Dive4Elements River for details. */ package org.dive4elements.river.utils; import org.dive4elements.river.artifacts.math.Linear; import org.apache.commons.math.analysis.interpolation.LinearInterpolator; import org.apache.commons.math.analysis.polynomials.PolynomialSplineFunction; import org.apache.commons.math.exception.DimensionMismatchException; import org.apache.commons.math.exception.NumberIsTooSmallException; import org.apache.commons.math.exception.NonMonotonousSequenceException; import gnu.trove.TDoubleArrayList; import java.util.Arrays; import java.util.Comparator; import org.apache.log4j.Logger; /** Utils to deal with Double precision values. */ public class DoubleUtil { /** Private log. */ private static Logger log = Logger.getLogger(DoubleUtil.class); public static final double DEFAULT_STEP_PRECISION = 1e6; public static final Comparator<double []> DOUBLE_PAIR_CMP = new Comparator<double []>() { @Override public int compare(double [] a, double [] b) { double diff = a[0] - b[0]; if (diff < 0d) return -1; if (diff > 0d) return +1; return 0; } }; /** EPSILON for comparison of double precision values. */ public static final double EPSILON = 1e-4; private DoubleUtil() { } public static final double [] explode( double from, double to, double step ) { return explode(from, to, step, DEFAULT_STEP_PRECISION); } public static final double round(double x, double precision) { return Math.round(x * precision)/precision; } public static final double round(double x) { return Math.round(x * DEFAULT_STEP_PRECISION)/DEFAULT_STEP_PRECISION; } /** * Returns array with values from parameter from to to * with given step width. * from and to are included. */ public static final double [] explode( double from, double to, double step, double precision ) { double lower = from; double diff = to - from; double tmp = diff / step; int num = (int)Math.abs(Math.ceil(tmp)) + 1; if (num < 1) { return new double[0]; } double [] values = new double[num]; if (from > to) { step = -step; } double max = Math.max(from, to); for (int idx = 0; idx < num; idx++) { if (lower - max > EPSILON) { return Arrays.copyOfRange(values, 0, idx); } values[idx] = round(lower, precision); lower += step; } return values; } public static final double interpolateSorted( double [] xs, double [] ys, double x ) { int lo = 0, hi = xs.length-1; int mid = -1; while (lo <= hi) { mid = (lo + hi) >> 1; double mx = xs[mid]; if (x < mx) hi = mid - 1; else if (x > mx) lo = mid + 1; else return ys[mid]; } if (mid < lo) { return lo >= xs.length ? Double.NaN : Linear.linear(x, xs[mid], xs[mid+1], ys[mid], ys[mid+1]); } return hi < 0 ? Double.NaN : Linear.linear(x, xs[mid-1], xs[mid], ys[mid-1], ys[mid]); } public static final boolean isIncreasing(double [] array) { int inc = 0; int dec = 0; int sweet = (array.length-1)/2; for (int i = 1; i < array.length; ++i) { if (array[i] > array[i-1]) { if (++inc > sweet) { return true; } } else if (++dec > sweet) { return false; } } return inc > sweet; } public static final double [] swap(double [] array) { int lo = 0; int hi = array.length-1; while (hi > lo) { double t = array[lo]; array[lo] = array[hi]; array[hi] = t; ++lo; --hi; } return array; } public static final double [] swapClone(double [] in) { double [] out = new double[in.length]; for (int j = out.length-1, i = 0; j >= 0;) { out[j--] = in[i++]; } return out; } public static final double [] sumDiffs(double [] in) { double [] out = new double[in.length]; for (int i = 1; i < out.length; ++i) { out[i] = out[i-1] + Math.abs(in[i-1] - in[i]); } return out; } public static final double sum(double [] values) { double sum = 0.0; for (double value: values) { sum += value; } return sum; } public static final double [] fill(int N, double value) { double [] result = new double[N]; Arrays.fill(result, value); return result; } /** Use with parseSegments. */ public interface SegmentCallback { void newSegment(double from, double to, double [] values); } /** Call callback for every string split by colon. * Expected format FROM:TO:VALUE1,VALUE2,VALUE3*/ public static final void parseSegments( String input, SegmentCallback callback ) { TDoubleArrayList vs = new TDoubleArrayList(); for (String segmentStr: input.split(":")) { String [] parts = segmentStr.split(";"); if (parts.length < 3) { log.warn("invalid segment: '" + segmentStr + "'"); continue; } try { double from = Double.parseDouble(parts[0].trim()); double to = Double.parseDouble(parts[1].trim()); vs.resetQuick(); for (String valueStr: parts[3].split(",")) { vs.add(round(Double.parseDouble(valueStr.trim()))); } callback.newSegment(from, to, vs.toNativeArray()); } catch (NumberFormatException nfe) { log.warn("invalid segment: '" + segmentStr + "'"); } } } public static final boolean isValid(double [][] data) { for (double [] ds: data) { if (!isValid(ds)) { return false; } } return true; } public static final boolean isValid(double [] data) { for (double d: data) { if (Double.isNaN(d)) { return false; } } return true; } public static final boolean isNaN(double [] values) { for (double value: values) { if (!Double.isNaN(value)) { return false; } } return true; } /** In an array of doubles, search and return the maximum value. */ public static final double maxInArray(double[] values) { double max = - Double.MAX_VALUE; for (double d: values) { if (d > max) max = d; } return max; } /** Sort a and b with a as key. b is ordered accordingly */ public static final void sortByFirst(double [] a, double [] b) { // XXX: Not efficient but bulletproof. double [][] pairs = new double[a.length][2]; for (int i = 0; i < a.length; ++i) { double [] p = pairs[i]; p[0] = a[i]; p[1] = b[i]; } Arrays.sort(pairs, DOUBLE_PAIR_CMP); for (int i = 0; i < a.length; ++i) { double [] p = pairs[i]; a[i] = p[0]; b[i] = p[1]; } } public static void removeNaNs(TDoubleArrayList [] arrays) { int dest = 0; int A = arrays.length; int N = arrays[0].size(); OUTER: for (int i = 0; i < N; ++i) { for (int j = 0; j < A; ++j) { TDoubleArrayList a = arrays[j]; double v = a.getQuick(i); if (Double.isNaN(v)) { continue OUTER; } a.setQuick(dest, v); } ++dest; } if (dest < N) { for (int i = 0; i < A; ++i) { arrays[i].remove(dest, N-dest); } } } /** Convieniance function for results to get an interpolator. * This is basically a static wrapper to for LinearInterpolator.interpolate * with error handling. Returns null on error.*/ public static PolynomialSplineFunction getLinearInterpolator( double[] x, double[] y ) { /* We want copies here to remove NaNs but don't * change the original data. */ TDoubleArrayList tX = new TDoubleArrayList(); TDoubleArrayList tY = new TDoubleArrayList(); for (int i = 0; i < x.length; i++) { if (!Double.isNaN(y[i])) { tX.add(x[i]); tY.add(y[i]); } } LinearInterpolator lpol = new LinearInterpolator(); try { return lpol.interpolate(tX.toNativeArray(), tY.toNativeArray()); } catch (DimensionMismatchException e) { log.error("KMs and Result values have different sizes. " + "Failed to interpolate: " + e.getMessage()); } catch (NonMonotonousSequenceException e) { log.error("KMs are not monotonous. Failed to interpolate: " + e.getMessage()); } catch (NumberIsTooSmallException e) { log.error("Result is to small. Failed to interpolate: " + e.getMessage()); } return null; } public static PolynomialSplineFunction getLinearInterpolator( TDoubleArrayList x, TDoubleArrayList y ) { return getLinearInterpolator(x.toNativeArray(), y.toNativeArray()); } } // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :