Mercurial > dive4elements > river

comparison artifacts/src/main/java/org/dive4elements/river/artifacts/sinfo/tkhcalculation/BedQualityD50KmValueFinder.java @ 8964:45f1ad66560e

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Code cleanup concerning calculations: improved error handling; improved interpolation; bed heights are now always used for spatial discretisation
author gernotbelger
date Thu, 29 Mar 2018 15:48:17 +0200
parents 48d5812e8fd5
children 9b9f5f4ddb80
comparison
equal deleted inserted replaced
8963:b98fbd91f64a 8964:45f1ad66560e
18 import org.apache.commons.math.ArgumentOutsideDomainException; 18 import org.apache.commons.math.ArgumentOutsideDomainException;
19 import org.apache.commons.math.analysis.interpolation.LinearInterpolator; 19 import org.apache.commons.math.analysis.interpolation.LinearInterpolator;
20 import org.apache.commons.math.analysis.polynomials.PolynomialSplineFunction; 20 import org.apache.commons.math.analysis.polynomials.PolynomialSplineFunction;
21 import org.apache.log4j.Logger; 21 import org.apache.log4j.Logger;
22 import org.dive4elements.river.artifacts.math.Utils; 22 import org.dive4elements.river.artifacts.math.Utils;
23 import org.dive4elements.river.artifacts.model.Calculation;
23 import org.dive4elements.river.backend.SedDBSessionHolder; 24 import org.dive4elements.river.backend.SedDBSessionHolder;
24 import org.dive4elements.river.model.River; 25 import org.dive4elements.river.model.River;
25 import org.hibernate.SQLQuery; 26 import org.hibernate.SQLQuery;
26 import org.hibernate.Session; 27 import org.hibernate.Session;
27 import org.hibernate.type.StandardBasicTypes; 28 import org.hibernate.type.StandardBasicTypes;
49 /** 50 /**
50 * Query selecting all sub layer bed measurements with their d50 for a km range and a time period<br /> 51 * Query selecting all sub layer bed measurements with their d50 for a km range and a time period<br />
51 * <br /> 52 * <br />
52 * A km may have bed measurements for multiple dates, multiple distances from the river bank, and multiple depth layers. 53 * A km may have bed measurements for multiple dates, multiple distances from the river bank, and multiple depth layers.
53 * The query filters by km range, time period and layer (sub layer: below bed to max. 50 cm depth).<br /> 54 * The query filters by km range, time period and layer (sub layer: below bed to max. 50 cm depth).<br />
54 * 55 *
55 * If PostgreSQL would support a median aggregate function like Oracle does, the aggregation could be placed into this query. 56 * If PostgreSQL would support a median aggregate function like Oracle does, the aggregation could be placed into this
57 * query.
56 */ 58 */
57 private static final String SQL_BED_D50_SUBLAYER_MEASUREMENT = 59 private static final String SQL_BED_D50_SUBLAYER_MEASUREMENT = "SELECT t.km, t.datum, p.tiefevon, p.tiefebis, a.d50"
58 "SELECT t.km, t.datum, p.tiefevon, p.tiefebis, a.d50" 60 + " FROM sohltest t INNER JOIN station s ON t.stationid = s.stationid" + " INNER JOIN gewaesser g ON s.gewaesserid = g.gewaesserid"
59 + " FROM sohltest t INNER JOIN station s ON t.stationid = s.stationid" 61 + " INNER JOIN sohlprobe p ON t.sohltestid = p.sohltestid" + " INNER JOIN siebanalyse a ON p.sohlprobeid = a.sohlprobeid"
60 + " INNER JOIN gewaesser g ON s.gewaesserid = g.gewaesserid" 62 + " WHERE (g.name = :name) AND (s.km BETWEEN :fromkm - 0.0001 AND :tokm + 0.0001)" + " AND (p.tiefevon > 0.0) AND (p.tiefebis <= 0.5)"
61 + " INNER JOIN sohlprobe p ON t.sohltestid = p.sohltestid" 63 + " AND (t.datum BETWEEN :fromdate AND :todate)" + " ORDER BY t.km ASC, a.d50 ASC";
62 + " INNER JOIN siebanalyse a ON p.sohlprobeid = a.sohlprobeid" 64
63 + " WHERE (g.name = :name) AND (s.km BETWEEN :fromkm - 0.0001 AND :tokm + 0.0001)" 65 private Calculation problems;
64 + " AND (p.tiefevon > 0.0) AND (p.tiefebis <= 0.5)"
65 + " AND (t.datum BETWEEN :fromdate AND :todate)"
66 + " ORDER BY t.km ASC, a.d50 ASC";
67 66
68 /** 67 /**
69 * Real linear interpolator for kms and d50 values (m) 68 * Real linear interpolator for kms and d50 values (m)
70 */ 69 */
71 private final PolynomialSplineFunction interpolator; 70 private final PolynomialSplineFunction interpolator;
72 71
72 /***** CONSTRUCTORS *****/
73 73
74 /***** CONSTRUCTORS *****/ 74 private BedQualityD50KmValueFinder(final Calculation problems, final double[] kms, final double[] values) {
75 75 this.problems = problems;
76 private BedQualityD50KmValueFinder(final double[] kms, final double[] values) { 76
77 // FIXME: check: max distance prüfen? dann D4E-LinearInterpolator verwenden
77 this.interpolator = new LinearInterpolator().interpolate(kms, values); 78 this.interpolator = new LinearInterpolator().interpolate(kms, values);
78 } 79 }
79 80
80 /***** METHODS *****/ 81 /***** METHODS *****/
81 82
82 /** 83 /**
83 * Sohlbeschaffenheit (D50 Korndurchmesser aus Seddb) 84 * Sohlbeschaffenheit (D50 Korndurchmesser aus Seddb)
84 * Abhängig von Peiljahr 85 * Abhängig von Peiljahr
86 *
87 * @param problems
85 */ 88 */
86 public static BedQualityD50KmValueFinder loadBedMeasurements(final River river, final DoubleRange kmRange, final int soundingYear, final int validYears) { 89 public static BedQualityD50KmValueFinder loadBedMeasurements(final Calculation problems, final River river, final DoubleRange kmRange,
90 final int soundingYear, final int validYears) {
87 91
88 /* construct valid measurement time range */ 92 /* construct valid measurement time range */
89 final Calendar cal = Calendar.getInstance(); 93 final Calendar cal = Calendar.getInstance();
90 cal.clear(); 94 cal.clear();
91 95
96 final Date endTime = cal.getTime(); 100 final Date endTime = cal.getTime();
97 101
98 log.debug(String.format("loadValues km %.3f - %.3f %tF - %tF", kmRange.getMinimumDouble(), kmRange.getMaximumDouble(), startTime, endTime)); 102 log.debug(String.format("loadValues km %.3f - %.3f %tF - %tF", kmRange.getMinimumDouble(), kmRange.getMaximumDouble(), startTime, endTime));
99 final Session session = SedDBSessionHolder.HOLDER.get(); 103 final Session session = SedDBSessionHolder.HOLDER.get();
100 final SQLQuery sqlQuery = session.createSQLQuery(SQL_BED_D50_SUBLAYER_MEASUREMENT).addScalar("km", StandardBasicTypes.DOUBLE) 104 final SQLQuery sqlQuery = session.createSQLQuery(SQL_BED_D50_SUBLAYER_MEASUREMENT).addScalar("km", StandardBasicTypes.DOUBLE)
101 .addScalar("datum", StandardBasicTypes.DATE).addScalar("tiefevon", StandardBasicTypes.DOUBLE) 105 .addScalar("datum", StandardBasicTypes.DATE).addScalar("tiefevon", StandardBasicTypes.DOUBLE).addScalar("tiefebis", StandardBasicTypes.DOUBLE)
102 .addScalar("tiefebis", StandardBasicTypes.DOUBLE).addScalar("d50", StandardBasicTypes.DOUBLE); 106 .addScalar("d50", StandardBasicTypes.DOUBLE);
103 final String seddbRiver = river.nameForSeddb(); 107 final String seddbRiver = river.nameForSeddb();
104 sqlQuery.setString("name", seddbRiver); 108 sqlQuery.setString("name", seddbRiver);
105 sqlQuery.setDouble("fromkm", kmRange.getMinimumDouble()); 109 sqlQuery.setDouble("fromkm", kmRange.getMinimumDouble());
106 sqlQuery.setDouble("tokm", kmRange.getMaximumDouble()); 110 sqlQuery.setDouble("tokm", kmRange.getMaximumDouble());
107 sqlQuery.setDate("fromdate", startTime); 111 sqlQuery.setDate("fromdate", startTime);
109 113
110 final List<Object[]> rows = sqlQuery.list(); 114 final List<Object[]> rows = sqlQuery.list();
111 final TDoubleArrayList kms = new TDoubleArrayList(); 115 final TDoubleArrayList kms = new TDoubleArrayList();
112 final TDoubleArrayList values = new TDoubleArrayList(); 116 final TDoubleArrayList values = new TDoubleArrayList();
113 final TDoubleArrayList kmd50s = new TDoubleArrayList(); 117 final TDoubleArrayList kmd50s = new TDoubleArrayList();
118
114 for (int i = 0; i <= rows.size() - 1; i++) { 119 for (int i = 0; i <= rows.size() - 1; i++) {
115 kmd50s.add((double) rows.get(i)[4]); 120 kmd50s.add((double) rows.get(i)[4]);
116 if (((i == rows.size() - 1) || !Utils.epsilonEquals((double) rows.get(i)[0], (double) rows.get(i+1)[0], 0.0001))) { 121 if (((i == rows.size() - 1) || !Utils.epsilonEquals((double) rows.get(i)[0], (double) rows.get(i + 1)[0], 0.0001))) {
117 int k = kmd50s.size() / 2; 122 final int k = kmd50s.size() / 2;
118 values.add(((k + k < kmd50s.size()) ? kmd50s.get(k) : (kmd50s.get(k-1) + kmd50s.get(k)) / 2) / 1000); 123 values.add(((k + k < kmd50s.size()) ? kmd50s.get(k) : (kmd50s.get(k - 1) + kmd50s.get(k)) / 2) / 1000);
119 kms.add((double) rows.get(i)[0]); 124 kms.add((double) rows.get(i)[0]);
120 log.debug(String.format("loadValues km %.3f d50(mm) %.1f count %d", kms.get(kms.size()-1), values.get(values.size()-1), kmd50s.size())); 125 log.debug(String.format("loadValues km %.3f d50(mm) %.1f count %d", kms.get(kms.size() - 1), values.get(values.size() - 1), kmd50s.size()));
121 kmd50s.clear(); 126 kmd50s.clear();
127 }
122 } 128 }
129
130 if (kms.size() < 2 || values.size() < 2) {
131 problems.addProblem("bedqualityd50kmvaluefinder.empty", soundingYear);
132 return null;
123 } 133 }
134
124 try { 135 try {
125 return new BedQualityD50KmValueFinder(kms.toNativeArray(), values.toNativeArray()); 136 return new BedQualityD50KmValueFinder(problems, kms.toNativeArray(), values.toNativeArray());
126 } 137 }
127 catch (final Exception e) { 138 catch (final Exception e) {
128 e.printStackTrace(); 139 e.printStackTrace();
140 problems.addProblem("bedqualityd50kmvaluefinder.error", e.getLocalizedMessage());
129 return null; 141 return null;
130 } 142 }
131 } 143 }
132 144
133 /** 145 /**
134 * Returns the d50 value interpolated according to a km 146 * Returns the d50 value interpolated according to a km
135 * 147 *
136 * @return d50 (m) of the km, or NaN 148 * @return d50 (m) of the km, or NaN
137 */ 149 */
138 public double findD50(final double km) throws ArgumentOutsideDomainException { 150 public double findD50(final double km) {
139 return this.interpolator.value(km); 151 try {
140 /* 152 return this.interpolator.value(km);
141 * ohne interpolation: 153 }
142 * if ((kms == null) || (kms.size() == 0)) 154 catch (final ArgumentOutsideDomainException e) {
143 * return Double.NaN; 155 e.printStackTrace();
144 * int i = kms.binarySearch(km); 156
145 * if (i >= 0) 157 if (this.problems != null) {
146 * return values.get(i); 158 this.problems.addProblem(km, "bedqualityd50kmvaluefinder.missing");
147 * i = -i - 1; 159 // Report only once
148 * if ((i - 1 >= 0) && Utils.epsilonEquals(km, kms.get(i - 1), 0.0001)) 160 this.problems = null;
149 * return values.get(i - 1); 161 }
150 * else if ((i >= 0) && (i <= kms.size() - 1) && Utils.epsilonEquals(km, kms.get(i), 0.0001)) 162
151 * return values.get(i); 163 return Double.NaN;
152 * else 164 }
153 * return Double.NaN;
154 */
155 } 165 }
156 } 166 }
http://dive4elements.wald.intevation.org