|
8854
|
1 /* Copyright (C) 2017 by Bundesanstalt für Gewässerkunde |
|
8877
|
2 * Software engineering by |
|
|
3 * Björnsen Beratende Ingenieure GmbH |
|
8854
|
4 * Dr. Schumacher Ingenieurbüro für Wasser und Umwelt |
|
|
5 * |
|
|
6 * This file is Free Software under the GNU AGPL (>=v3) |
|
|
7 * and comes with ABSOLUTELY NO WARRANTY! Check out the |
|
|
8 * documentation coming with Dive4Elements River for details. |
|
|
9 */ |
|
|
10 package org.dive4elements.river.artifacts.sinfo.flowdepth; |
|
|
11 |
|
8883
|
12 import java.util.ArrayList; |
|
8891
|
13 import java.util.Calendar; |
|
8854
|
14 import java.util.Collection; |
|
8883
|
15 import java.util.Collections; |
|
8891
|
16 import java.util.Date; |
|
8854
|
17 import java.util.List; |
|
|
18 |
|
8894
|
19 import org.apache.commons.lang.math.DoubleRange; |
|
8883
|
20 import org.apache.commons.math.FunctionEvaluationException; |
|
|
21 import org.apache.commons.math.analysis.UnivariateRealFunction; |
|
8898
|
22 import org.apache.log4j.Logger; |
|
8879
|
23 import org.dive4elements.artifacts.ArtifactDatabase; |
|
8854
|
24 import org.dive4elements.artifacts.CallContext; |
|
|
25 import org.dive4elements.river.artifacts.BedHeightsArtifact; |
|
|
26 import org.dive4elements.river.artifacts.model.Calculation; |
|
|
27 import org.dive4elements.river.artifacts.model.CalculationResult; |
|
8898
|
28 import org.dive4elements.river.artifacts.model.DateRange; |
|
8854
|
29 import org.dive4elements.river.artifacts.model.LocationProvider; |
|
8901
|
30 import org.dive4elements.river.artifacts.model.QKms; |
|
8854
|
31 import org.dive4elements.river.artifacts.model.WKms; |
|
|
32 import org.dive4elements.river.artifacts.resources.Resources; |
|
|
33 import org.dive4elements.river.artifacts.sinfo.SINFOArtifact; |
|
|
34 import org.dive4elements.river.artifacts.sinfo.flowdepth.FlowDepthAccess.DifferencesPair; |
|
8894
|
35 import org.dive4elements.river.artifacts.sinfo.util.BedHeightInfo; |
|
|
36 import org.dive4elements.river.artifacts.sinfo.util.RiverInfo; |
|
|
37 import org.dive4elements.river.artifacts.sinfo.util.WstInfo; |
|
8882
|
38 import org.dive4elements.river.artifacts.states.WaterlevelData; |
|
|
39 import org.dive4elements.river.artifacts.states.WaterlevelFetcher; |
|
8854
|
40 import org.dive4elements.river.model.BedHeight; |
|
8883
|
41 import org.dive4elements.river.model.BedHeightValue; |
|
8854
|
42 import org.dive4elements.river.model.Gauge; |
|
|
43 import org.dive4elements.river.model.River; |
|
8883
|
44 import org.dive4elements.river.utils.DoubleUtil; |
|
8854
|
45 import org.dive4elements.river.utils.GaugeIndex; |
|
|
46 import org.dive4elements.river.utils.RiverUtils; |
|
|
47 |
|
|
48 class FlowDepthCalculation { |
|
|
49 |
|
8898
|
50 private static Logger log = Logger.getLogger(FlowDepthCalculation.class); |
|
|
51 |
|
8891
|
52 private static final int VALID_BED_MEASUREMENT_YEARS = 20; |
|
|
53 |
|
8854
|
54 private static final String CSV_NOT_IN_GAUGE_RANGE = "export.waterlevel.csv.not.in.gauge.range"; |
|
|
55 |
|
8877
|
56 private final CallContext context; |
|
8854
|
57 |
|
8882
|
58 public FlowDepthCalculation(final CallContext context) { |
|
8877
|
59 this.context = context; |
|
8854
|
60 } |
|
|
61 |
|
8877
|
62 public CalculationResult calculate(final SINFOArtifact sinfo) { |
|
8854
|
63 |
|
8879
|
64 /* |
|
|
65 * find the user of this artifact, sadly this is not part of the calling context, so instead we determine the |
|
|
66 * owner oft the artifact |
|
|
67 */ |
|
|
68 final ArtifactDatabase database = this.context.getDatabase(); |
|
|
69 final String user = database.findArtifactUser(sinfo.identifier()); |
|
8854
|
70 |
|
8877
|
71 /* access input data */ |
|
|
72 final FlowDepthAccess access = new FlowDepthAccess(sinfo); |
|
|
73 final River river = access.getRiver(); |
|
8894
|
74 final RiverInfo riverInfo = new RiverInfo(river); |
|
8854
|
75 |
|
8877
|
76 final Collection<DifferencesPair> diffPairs = access.getDifferencePairs(); |
|
|
77 |
|
|
78 final double from = access.getFrom(); |
|
|
79 final double to = access.getTo(); |
|
8894
|
80 final DoubleRange calcRange = new DoubleRange(from, to); |
|
8877
|
81 |
|
|
82 final boolean useTkh = access.isUseTransportBodies(); |
|
|
83 |
|
|
84 /* calculate results for each diff pair */ |
|
|
85 final Calculation problems = new Calculation(); |
|
|
86 |
|
|
87 final List<Gauge> gauges = river.determineGauges(from, to); |
|
|
88 final GaugeIndex gaugeIndex = new GaugeIndex(gauges); |
|
|
89 |
|
8882
|
90 final String calcModeLabel = Resources.getMsg(this.context.getMeta(), sinfo.getCalculationMode().name()); |
|
8877
|
91 |
|
8894
|
92 final FlowDepthCalculationResults results = new FlowDepthCalculationResults(calcModeLabel, user, riverInfo, calcRange, useTkh); |
|
8877
|
93 |
|
|
94 for (final DifferencesPair diffPair : diffPairs) { |
|
8898
|
95 final FlowDepthCalculationResult result = calculateResult(river, calcRange, diffPair, problems, gaugeIndex, useTkh); |
|
8882
|
96 if (result != null) |
|
8877
|
97 results.addResult(result); |
|
|
98 } |
|
|
99 |
|
8882
|
100 return new CalculationResult(results, problems); |
|
8877
|
101 } |
|
|
102 |
|
8898
|
103 /** |
|
|
104 * Calculates one W-MSH differences pair. |
|
|
105 */ |
|
8894
|
106 private FlowDepthCalculationResult calculateResult(final River river, final DoubleRange calcRange, final DifferencesPair diffPair, |
|
8898
|
107 final Calculation problems, final GaugeIndex gaugeIndex, final boolean useTkh) { |
|
8877
|
108 |
|
|
109 /* access real input data from database */ |
|
|
110 final String soundingId = diffPair.getSoundingId(); |
|
|
111 final String wstId = diffPair.getWstId(); |
|
|
112 |
|
8894
|
113 final BedHeight bedHeight = loadBedHeight(soundingId); |
|
8882
|
114 if (bedHeight == null) { |
|
8884
|
115 final String message = Resources.format(this.context.getMeta(), "Failed to access sounding with id '{0}'", soundingId); |
|
8877
|
116 problems.addProblem(message); |
|
|
117 return null; |
|
|
118 } |
|
|
119 |
|
8882
|
120 /* REMARK: fetch ALL wst kms, because we want to determine the original reference gauge */ |
|
8884
|
121 final WaterlevelData waterlevel = new WaterlevelFetcher().findWaterlevel(this.context, wstId, Double.NaN, Double.NaN); |
|
8882
|
122 if (waterlevel == null) { |
|
8884
|
123 final String message = Resources.format(this.context.getMeta(), "Failed to access waterlevel with id '{0}'", wstId); |
|
8877
|
124 problems.addProblem(message); |
|
|
125 return null; |
|
|
126 } |
|
8882
|
127 final WKms wstKms = waterlevel.getWkms(); |
|
|
128 |
|
8883
|
129 final String wspLabel = wstKms.getName(); |
|
|
130 final String soundingLabel = bedHeight.getDescription(); |
|
|
131 final String label = String.format("%s - %s", wspLabel, soundingLabel); |
|
8877
|
132 |
|
8883
|
133 checkYearDifference(label, waterlevel, bedHeight, problems); |
|
8894
|
134 checkWaterlevelDiscretisation(wstKms, calcRange, problems); |
|
8901
|
135 // TODO: prüfen, ob sohlhöhen die calcRange abdecken/überschneiden |
|
8882
|
136 |
|
|
137 /* re-determine the reference gauge, in the same way as the WaterlevelArtifact would do it */ |
|
8884
|
138 final String notinrange = Resources.getMsg(this.context.getMeta(), CSV_NOT_IN_GAUGE_RANGE, CSV_NOT_IN_GAUGE_RANGE); |
|
8882
|
139 |
|
|
140 final Gauge refGauge = waterlevel.findReferenceGauge(river); |
|
|
141 final String refGaugeName = refGauge == null ? notinrange : refGauge.getName(); |
|
8877
|
142 |
|
|
143 final BedHeightInfo sounding = BedHeightInfo.from(bedHeight); |
|
8883
|
144 final int wspYear = waterlevel.getYear(); |
|
8882
|
145 final WstInfo wstInfo = new WstInfo(wspLabel, wspYear, refGaugeName); |
|
8877
|
146 |
|
|
147 final FlowDepthCalculationResult resultData = new FlowDepthCalculationResult(label, wstInfo, sounding); |
|
|
148 |
|
8901
|
149 boolean doCalcTkh = useTkh; |
|
|
150 if (doCalcTkh && !(wstKms instanceof QKms)) { |
|
8884
|
151 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingQ", null, label); |
|
8877
|
152 problems.addProblem(message); |
|
8901
|
153 doCalcTkh = false; |
|
8877
|
154 } |
|
|
155 |
|
8898
|
156 BedQualityD50KmValueFinder bedMeasurementsFinder = null; |
|
8901
|
157 if (doCalcTkh) |
|
|
158 bedMeasurementsFinder = loadBedMeasurements(river, calcRange, sounding.getYear().intValue()); |
|
8877
|
159 |
|
|
160 final String bedHeightLabel = bedHeight.getDescription(); |
|
|
161 final String wstLabel = wstKms.getName(); |
|
|
162 |
|
8884
|
163 final UnivariateRealFunction wstInterpolator = DoubleUtil.getLinearInterpolator(wstKms.allKms(), wstKms.allWs()); |
|
8898
|
164 UnivariateRealFunction qInterpolator = null; |
|
|
165 DoubleRange qRange = null; |
|
8901
|
166 if (doCalcTkh) { |
|
|
167 qInterpolator = DoubleUtil.getLinearInterpolator(((QKms) wstKms).allKms(), ((QKms) wstKms).allQs()); |
|
|
168 qRange = new DoubleRange( ((QKms) wstKms).allQs().min(), ((QKms) wstKms).allQs().max()); |
|
8898
|
169 } |
|
|
170 |
|
8883
|
171 // FIXME: sort by station first, but in what direction? |
|
8898
|
172 // FIXME: using river.getKmUp()? |
|
8883
|
173 final List<BedHeightValue> values = bedHeight.getValues(); |
|
8877
|
174 |
|
8883
|
175 final List<BedHeightValue> sortedValues = new ArrayList<>(values); |
|
|
176 Collections.sort(sortedValues, new BedHeightStationComparator()); |
|
8882
|
177 |
|
8898
|
178 // FIXME: wie wird ggf. interpoliert? prüfung ob werte vorhanden? |
|
|
179 /* SoilKind lastKind = SoilKind.mobil; */ |
|
|
180 SoilKindKmValueFinder soilKindFinder = null; |
|
8901
|
181 if (doCalcTkh) { |
|
8898
|
182 soilKindFinder = new SoilKindKmValueFinder(); |
|
|
183 soilKindFinder.loadValues(river, calcRange); |
|
|
184 } |
|
|
185 |
|
|
186 FlowVelocityModelKmValueFinder flowVelocitiesFinder = null; |
|
8901
|
187 if (doCalcTkh) { |
|
8898
|
188 flowVelocitiesFinder = new FlowVelocityModelKmValueFinder(); |
|
|
189 flowVelocitiesFinder.loadValues(river, calcRange, qRange); |
|
|
190 } |
|
8886
|
191 |
|
8883
|
192 for (final BedHeightValue bedHeightValue : sortedValues) { |
|
8877
|
193 |
|
8883
|
194 final Double station = bedHeightValue.getStation(); |
|
|
195 if (station == null || station.isNaN()) |
|
|
196 continue; |
|
|
197 |
|
|
198 final Double meanBedHeightDbl = bedHeightValue.getHeight(); |
|
|
199 if (meanBedHeightDbl == null || meanBedHeightDbl.isNaN()) |
|
|
200 continue; |
|
|
201 |
|
|
202 final double km = station; |
|
|
203 final double meanBedHeight = meanBedHeightDbl; |
|
|
204 |
|
8894
|
205 if (!calcRange.containsDouble(km)) |
|
|
206 continue; |
|
|
207 |
|
8883
|
208 try { |
|
|
209 // FIXME: check out of range |
|
|
210 final double wst = wstInterpolator.value(km); |
|
|
211 |
|
|
212 final double flowDepth = wst - meanBedHeight; |
|
|
213 |
|
|
214 // FIXME: piecewise constant interpolation? |
|
|
215 // final double discharge = wstKms instanceof QKms ? ((QKms) wstKms).getQ(i) : Double.NaN; |
|
8898
|
216 double discharge; |
|
|
217 if (qInterpolator != null) |
|
|
218 discharge = qInterpolator.value(km); |
|
|
219 else |
|
|
220 discharge = Double.NaN; |
|
8883
|
221 |
|
|
222 // FIXME: calculate tkh |
|
8898
|
223 double tkh = 0; |
|
8901
|
224 if (doCalcTkh) { |
|
|
225 double d50 = bedMeasurementsFinder.findD50(km); |
|
|
226 if (Double.isNaN(d50)) { |
|
8898
|
227 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingD50", null, label); |
|
|
228 problems.addProblem(km, message); |
|
8901
|
229 //FIXME: cumulate problems to one message? |
|
8898
|
230 } |
|
8901
|
231 if (!Double.isNaN(d50)) { |
|
|
232 if (flowVelocitiesFinder.findKmQValues(km, discharge)) { |
|
|
233 tkh = calculateTkh(wst - meanBedHeight, flowVelocitiesFinder.getFindVmainFound(), d50, flowVelocitiesFinder.getFindTauFound()); |
|
|
234 /* log.debug(String.format("calculateTkh km %.3f q %.0f w %.2f mbh %.2f vm %.1f tau %.1f d50(mm) %.1f tkh(cm) %.1f", |
|
|
235 km, discharge, wst, meanBedHeight, flowVelocitiesFinder.getFindVmainFound(), flowVelocitiesFinder.getFindTauFound(), d50*1000, tkh)); */ |
|
|
236 } |
|
|
237 else { |
|
|
238 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingQ", null, label); |
|
|
239 problems.addProblem(km, message); |
|
|
240 //FIXME: cumulate problems to one message? |
|
|
241 } |
|
8898
|
242 } |
|
|
243 else |
|
|
244 tkh = Double.NaN; |
|
|
245 } |
|
|
246 |
|
|
247 // Soil kind |
|
8901
|
248 SoilKind kind = SoilKind.mobil; |
|
|
249 if (doCalcTkh) { |
|
8898
|
250 try { |
|
|
251 kind = soilKindFinder.findSoilKind(km); |
|
|
252 } catch (Exception e) { |
|
|
253 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingSoilKind", null, label); |
|
|
254 problems.addProblem(km, message); |
|
8901
|
255 //FIXME: cumulate problems to one message? |
|
8898
|
256 } |
|
|
257 } |
|
|
258 |
|
|
259 /* // REMARK: bissl spielerei zum testen damit die sohlart nicht zu schnell wechselt |
|
|
260 SoilKind kind; |
|
8895
|
261 final boolean changeKind = Math.random() > 0.95; |
|
8886
|
262 if (changeKind) { |
|
|
263 switch (lastKind) { |
|
|
264 case starr: |
|
|
265 kind = SoilKind.mobil; |
|
|
266 break; |
|
|
267 case mobil: |
|
|
268 default: |
|
|
269 kind = SoilKind.starr; |
|
|
270 break; |
|
|
271 } |
|
|
272 } else |
|
|
273 kind = lastKind; |
|
|
274 lastKind = kind; |
|
8898
|
275 */ |
|
|
276 |
|
|
277 //tkh = 100 + 10 * (Math.random() - 0.5); |
|
8886
|
278 |
|
|
279 final double flowDepthTkh; |
|
|
280 final double tkhUp; |
|
|
281 final double tkhDown; |
|
|
282 switch (kind) { |
|
|
283 case starr: |
|
|
284 flowDepthTkh = wst - (meanBedHeight + tkh / 100); |
|
|
285 tkhUp = tkh; |
|
|
286 tkhDown = 0; |
|
|
287 break; |
|
|
288 |
|
|
289 case mobil: |
|
|
290 default: |
|
|
291 flowDepthTkh = wst - (meanBedHeight + tkh / 200); |
|
|
292 tkhUp = tkh / 2; |
|
|
293 tkhDown = -tkh / 2; |
|
|
294 break; |
|
|
295 } |
|
|
296 |
|
8883
|
297 // REMARK: access the location once only during calculation |
|
|
298 final String location = LocationProvider.getLocation(river.getName(), km); |
|
|
299 |
|
|
300 // REMARK: access the gauge once only during calculation |
|
|
301 final Gauge gauge = findGauge(waterlevel, refGauge, gaugeIndex, km); |
|
|
302 |
|
|
303 final String gaugeLabel = gauge == null ? notinrange : gauge.getName(); |
|
|
304 |
|
8898
|
305 resultData.addRow(km, flowDepth, flowDepthTkh, kind, tkh, tkhUp, tkhDown, wst, discharge, wstLabel, |
|
|
306 gaugeLabel, meanBedHeight, bedHeightLabel, location); |
|
8883
|
307 } |
|
|
308 catch (final FunctionEvaluationException e) { |
|
|
309 /* should only happen if out of range */ |
|
|
310 e.printStackTrace(); |
|
|
311 /* simply ignore */ |
|
|
312 } |
|
|
313 |
|
8877
|
314 } |
|
|
315 |
|
|
316 return resultData; |
|
|
317 } |
|
|
318 |
|
8883
|
319 /** |
|
8891
|
320 * Sohlbeschaffenheit (D50 Korndurchmesser aus Seddb) |
|
|
321 * Abhängig von Peiljahr |
|
|
322 */ |
|
8898
|
323 private BedQualityD50KmValueFinder loadBedMeasurements(final River river, final DoubleRange kmRange, final int soundingYear) { |
|
8891
|
324 |
|
|
325 /* construct valid measurement time range */ |
|
|
326 final Calendar cal = Calendar.getInstance(); |
|
|
327 cal.clear(); |
|
|
328 |
|
|
329 cal.set(soundingYear - VALID_BED_MEASUREMENT_YEARS, 0, 1); |
|
|
330 final Date startTime = cal.getTime(); |
|
|
331 |
|
|
332 cal.set(soundingYear + VALID_BED_MEASUREMENT_YEARS, 11, 31); |
|
|
333 final Date endTime = cal.getTime(); |
|
|
334 |
|
8898
|
335 final BedQualityD50KmValueFinder finder = new BedQualityD50KmValueFinder(); |
|
|
336 if (finder.loadValues(river, kmRange, new DateRange(startTime, endTime))) |
|
|
337 return finder; |
|
|
338 else |
|
|
339 return null; |
|
8891
|
340 } |
|
|
341 |
|
|
342 /** |
|
8883
|
343 * Checks the year difference between waterlevels and sounding, and issues a warning if too big. |
|
|
344 * |
|
|
345 * Zeitraum Zeitliche Differenz [a] |
|
|
346 * X ≥ 1998 ± 3 |
|
|
347 * 1958 ≤ X < 1998 ± 6 |
|
|
348 * 1918 ≤ X < 1958 ± 12 |
|
|
349 * X < 1918 ± 25 |
|
|
350 */ |
|
8884
|
351 private void checkYearDifference(final String label, final WaterlevelData waterlevel, final BedHeight sounding, final Calculation problems) { |
|
8883
|
352 |
|
|
353 final Integer soundingYear = sounding.getYear(); |
|
|
354 if (soundingYear == null) |
|
|
355 return; |
|
|
356 |
|
|
357 final int wstYear = waterlevel.getYear(); |
|
|
358 if (wstYear < 0) |
|
|
359 return; |
|
|
360 |
|
|
361 final int maxDifference = getMaxDifferenceYears(soundingYear); |
|
|
362 |
|
|
363 final int difference = Math.abs(soundingYear - wstYear); |
|
|
364 if (difference > maxDifference) { |
|
8894
|
365 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.year_difference", null, label, wstYear, |
|
|
366 soundingYear); |
|
8883
|
367 problems.addProblem(message); |
|
|
368 } |
|
|
369 } |
|
|
370 |
|
|
371 private int getMaxDifferenceYears(final int year) { |
|
|
372 |
|
|
373 if (year < 1918) |
|
|
374 return 25; |
|
|
375 |
|
|
376 if (1918 <= year && year < 1958) |
|
|
377 return 12; |
|
|
378 |
|
|
379 if (1958 <= year && year < 1998) |
|
|
380 return 6; |
|
|
381 |
|
|
382 /* >= 1998 */ |
|
|
383 return 3; |
|
|
384 } |
|
|
385 |
|
8884
|
386 private Gauge findGauge(final WaterlevelData waterlevel, final Gauge refGauge, final GaugeIndex gaugeIndex, final double km) { |
|
8882
|
387 |
|
|
388 // REMARK: using same logic as in WaterlevelExporter here |
|
|
389 |
|
|
390 final boolean showAllGauges = waterlevel.isShowAllGauges(); |
|
|
391 |
|
|
392 if (showAllGauges) |
|
|
393 return gaugeIndex.findGauge(km); |
|
|
394 |
|
|
395 if (refGauge.getRange().contains(km)) |
|
|
396 return refGauge; |
|
|
397 |
|
|
398 return null; |
|
|
399 } |
|
|
400 |
|
|
401 /* Checks if the discretisation of the waterlevel exceeds 1000m */ |
|
8894
|
402 |
|
|
403 private void checkWaterlevelDiscretisation(final WKms wstKms, final DoubleRange calcRange, final Calculation problems) { |
|
|
404 |
|
8882
|
405 final int size = wstKms.size(); |
|
|
406 for (int i = 0; i < size - 2; i++) { |
|
|
407 final double kmPrev = wstKms.getKm(i); |
|
|
408 final double kmNext = wstKms.getKm(i + 1); |
|
|
409 |
|
8894
|
410 /* only check if we are within the calculation range */ |
|
|
411 if (calcRange.overlapsRange(new DoubleRange(kmPrev, kmNext))) { |
|
|
412 if (Math.abs(kmPrev - kmNext) > 1) { |
|
|
413 final String label = wstKms.getName(); |
|
8882
|
414 |
|
8894
|
415 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.waterlevel_discretisation", null, label); |
|
|
416 problems.addProblem(kmPrev, message); |
|
|
417 } |
|
8882
|
418 } |
|
|
419 } |
|
|
420 } |
|
|
421 |
|
8894
|
422 private BedHeight loadBedHeight(final String soundingId) { |
|
8877
|
423 |
|
8883
|
424 // REMARK: absolutely unbelievable.... |
|
8898
|
425 // The way how bed-heights (and other data too) is accessed is different for nearly every calculation-type |
|
8882
|
426 // throughout flys. |
|
8877
|
427 // The knowledge on how to parse the datacage-ids is spread through the complete code-base... |
|
|
428 |
|
8882
|
429 // We use here the way on how bed-heights are accessed by the BedDifferenceAccess/BedDifferenceCalculation, but |
|
|
430 // this is plain random |
|
8877
|
431 final String[] parts = soundingId.split(";"); |
|
|
432 |
|
|
433 final BedHeightsArtifact artifact = (BedHeightsArtifact) RiverUtils.getArtifact(parts[0], this.context); |
|
|
434 |
|
|
435 final Integer bedheightId = artifact.getDataAsInteger("height_id"); |
|
8883
|
436 // REMARK: this only works with type 'single'; unclear on how to distinguish from epoch data (or whatever the |
|
8882
|
437 // other type means) |
|
8877
|
438 // Luckily, the requirement is to only access 'single' data here. |
|
|
439 // final String bedheightType = artifact.getDataAsString("type"); |
|
|
440 |
|
8883
|
441 // REMARK: BedDifferences uses this, but we also need the metadata of the BedHeight |
|
|
442 // REMARK: second absolutely awful thing: BedHeight is a hibernate binding class, accessing the database via |
|
8882
|
443 // hibernate stuff |
|
8877
|
444 // BedHeightFactory uses its own (direct) way of accessing the data, with its own implemented data classes. |
|
8882
|
445 // return BedHeightFactory.getHeight(bedheightType, bedheightId, from, to); |
|
8877
|
446 |
|
|
447 return BedHeight.getBedHeightById(bedheightId); |
|
|
448 } |
|
8898
|
449 |
|
|
450 /** |
|
|
451 * Calculates a transport body height |
|
|
452 * @param h flow depth in m |
|
|
453 * @param vm flow velocity in m |
|
|
454 * @param d50 grain diameter D50 in m (!) |
|
|
455 * @param tau shear stress in N/m^2 |
|
|
456 * @return transport body height in cm (!) |
|
|
457 */ |
|
|
458 private double calculateTkh(double h, double vm, double d50, double tau) { |
|
|
459 final double PHYS_G = 9.81; |
|
|
460 final double PHYS_SPECGRAV_S = 2.6; |
|
|
461 final double PHYS_VELOCCOEFF_N = 6; |
|
|
462 final double PHYS_FORMCOEFF_ALPHA = 0.7; |
|
|
463 final double PHYS_VISCOSITY_NUE = 1.3e-6; |
|
|
464 final double PHYS_GRAIN_DENSITY_RHOS = 2603; |
|
|
465 final double PHYS_WATER_DENSITY_RHO = 999.97; |
|
|
466 |
|
|
467 final double froude = vm / Math.sqrt(PHYS_G * h); |
|
|
468 final double partReynolds = Math.sqrt((PHYS_SPECGRAV_S - 1) * PHYS_G * d50) / PHYS_VISCOSITY_NUE * d50; |
|
|
469 final double critShields = 0.22 * Math.pow(partReynolds, -0.6) + 0.06 * Math.pow(10, 7.7 * Math.pow(partReynolds, -0.6)); |
|
|
470 final double critTau = critShields * (PHYS_GRAIN_DENSITY_RHOS - PHYS_WATER_DENSITY_RHO) * PHYS_G * d50; |
|
|
471 return 100 * h * (1 - Math.pow(froude, 2)) / (2 * PHYS_VELOCCOEFF_N * PHYS_FORMCOEFF_ALPHA) * (1 - critTau / tau); |
|
|
472 } |
|
8884
|
473 } |
