comparison gnv-artifacts/src/main/java/de/intevation/gnv/math/XYColumn.java @ 433:828df3ddb758

Added interpolation capabilities along z axis to XYColumns. gnv-artifacts/trunk@481 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Wed, 23 Dec 2009 12:20:25 +0000
parents 6a70e8883307
children 0eed5749fd63
comparison
equal deleted inserted replaced
432:6a70e8883307 433:828df3ddb758
1 package de.intevation.gnv.math; 1 package de.intevation.gnv.math;
2 2
3 import java.util.ArrayList; 3 import java.util.ArrayList;
4 import java.util.List; 4 import java.util.List;
5 import java.util.Collections;
6
7 import org.apache.commons.math.analysis.interpolation.SplineInterpolator;
8 import org.apache.commons.math.analysis.interpolation.UnivariateRealInterpolator;
9
10 import org.apache.commons.math.analysis.polynomials.PolynomialSplineFunction;
11
12 import org.apache.commons.math.analysis.UnivariateRealFunction;
13
14 import org.apache.commons.math.MathException;
15 import org.apache.commons.math.FunctionEvaluationException;
16
17 import org.apache.log4j.Logger;
5 18
6 /** 19 /**
7 * @author Ingo Weinzierl <ingo.weinzierl@intevation.de> 20 * @author Ingo Weinzierl <ingo.weinzierl@intevation.de>
8 * @author Sascha L. Teichmann <sascha.teichmann@intevation.de> 21 * @author Sascha L. Teichmann <sascha.teichmann@intevation.de>
9 */ 22 */
10 public class XYColumn 23 public class XYColumn
11 extends Point2d 24 extends Point2d
25 implements UnivariateRealFunction
12 { 26 {
27 private static Logger log = Logger.getLogger(Interpolation2D.class);
28
13 protected List<HeightValue> values; 29 protected List<HeightValue> values;
14 30
31 protected transient UnivariateRealFunction curve;
32
15 public XYColumn() { 33 public XYColumn() {
34 values = new ArrayList<HeightValue>();
16 } 35 }
17 36
18 public XYColumn(double x, double y, int i, int j) { 37 public XYColumn(double x, double y, int i, int j) {
19 super(x, y, i, j); 38 super(x, y, i, j);
39 values = new ArrayList<HeightValue>();
20 } 40 }
21 41
22 public void add(HeightValue value) { 42 public void add(HeightValue value) {
23 if (values == null) {
24 values = new ArrayList<HeightValue>();
25 }
26 values.add(value); 43 values.add(value);
27 } 44 }
28 45
29 public List<HeightValue> getValues() { 46 public List<HeightValue> getValues() {
30 return values; 47 return values;
31 } 48 }
49
50 public double value(double depth) {
51 try {
52 if (curve != null) {
53 return curve.value(depth);
54 }
55 }
56 catch (FunctionEvaluationException fee) {
57 log.error("evaluation failed", fee);
58 }
59
60 return Double.NaN;
61 }
62
63 public void prepare(XYDepth xyDepth) {
64 int N = values.size();
65 if (curve == null && N > 0) {
66 if (N == 1) {
67 // only one value -> constant function
68 curve = new ConstantFunction(values.get(0).v);
69 }
70 else { // more than on value
71 double depth = xyDepth.depth(this);
72 Collections.sort(values, HeightValue.INV_Z_COMPARATOR);
73
74 // if there is no value at 0 repeat first value
75 HeightValue first = values.get(0);
76 if (first.z < 0d) {
77 values.add(0, new HeightValue(0d, first.z, first.k-1));
78 }
79
80 // if there is no value at depth repeat last value
81 HeightValue last = values.get(N-1);
82 if (last.z > depth) {
83 values.add(new HeightValue(depth, last.z, last.k+1));
84 }
85 N = values.size();
86 if (N < 3) { // interpolate linear
87 first = values.get(0);
88 last = values.get(N-1);
89 curve = new LinearFunction.Univariate(
90 first.z, first.v,
91 last.z, last.v);
92 }
93 else { // higher degree interpolation
94 double [] z = new double[N];
95 double [] v = new double[N];
96 for (int i = 0; i < N; ++i) {
97 HeightValue h = values.get(i);
98 z[i] = h.z;
99 v[i] = h.v;
100 }
101 try {
102 curve = getInterpolator().interpolate(z, v);
103 }
104 catch (MathException me) {
105 log.error("interpolation failed", me);
106 }
107 }
108 }
109 }
110 else {
111 log.error("no points for interpolation");
112 }
113 }
114
115 protected UnivariateRealInterpolator getInterpolator() {
116 return new SplineInterpolator();
117 }
32 } 118 }
33 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8: 119 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8:

http://dive4elements.wald.intevation.org