Mercurial > dive4elements > river
annotate artifacts/src/main/java/org/dive4elements/river/artifacts/sinfo/flowdepth/FlowDepthCalculation.java @ 8911:37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
| author | mschaefer |
|---|---|
| date | Fri, 23 Feb 2018 09:30:24 +0100 |
| parents | 0a900d605d52 |
| children | e3519c3e7a0a |
| rev | line source |
|---|---|
| 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; |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
157 if (doCalcTkh) { |
| 8901 | 158 bedMeasurementsFinder = loadBedMeasurements(river, calcRange, sounding.getYear().intValue()); |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
159 if (bedMeasurementsFinder == null) { |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
160 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingSoilKind", null, label); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
161 problems.addProblem(message); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
162 doCalcTkh = false; |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
163 } |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
164 } |
| 8877 | 165 |
| 166 final String bedHeightLabel = bedHeight.getDescription(); | |
| 167 final String wstLabel = wstKms.getName(); | |
| 168 | |
| 8884 | 169 final UnivariateRealFunction wstInterpolator = DoubleUtil.getLinearInterpolator(wstKms.allKms(), wstKms.allWs()); |
| 8898 | 170 UnivariateRealFunction qInterpolator = null; |
| 171 DoubleRange qRange = null; | |
| 8901 | 172 if (doCalcTkh) { |
| 173 qInterpolator = DoubleUtil.getLinearInterpolator(((QKms) wstKms).allKms(), ((QKms) wstKms).allQs()); | |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
174 if (qInterpolator != null) |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
175 qRange = new DoubleRange( ((QKms) wstKms).allQs().min(), ((QKms) wstKms).allQs().max()); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
176 else { |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
177 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingQ", null, label); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
178 problems.addProblem(message); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
179 doCalcTkh = false; |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
180 } |
| 8898 | 181 } |
| 182 | |
| 8883 | 183 // FIXME: sort by station first, but in what direction? |
| 8898 | 184 // FIXME: using river.getKmUp()? |
| 8883 | 185 final List<BedHeightValue> values = bedHeight.getValues(); |
| 8877 | 186 |
| 8883 | 187 final List<BedHeightValue> sortedValues = new ArrayList<>(values); |
| 188 Collections.sort(sortedValues, new BedHeightStationComparator()); | |
| 8882 | 189 |
| 8898 | 190 // FIXME: wie wird ggf. interpoliert? prüfung ob werte vorhanden? |
| 191 /* SoilKind lastKind = SoilKind.mobil; */ | |
| 192 SoilKindKmValueFinder soilKindFinder = null; | |
| 8901 | 193 if (doCalcTkh) { |
| 8898 | 194 soilKindFinder = new SoilKindKmValueFinder(); |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
195 if (!soilKindFinder.loadValues(river, calcRange)) { |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
196 doCalcTkh = false; |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
197 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingSoilKind", null, label); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
198 problems.addProblem(message); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
199 } |
| 8898 | 200 } |
| 201 | |
| 202 FlowVelocityModelKmValueFinder flowVelocitiesFinder = null; | |
| 8901 | 203 if (doCalcTkh) { |
| 8898 | 204 flowVelocitiesFinder = new FlowVelocityModelKmValueFinder(); |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
205 if (!flowVelocitiesFinder.loadValues(river, calcRange, qRange)) { |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
206 doCalcTkh = false; |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
207 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingVelocity", null, label); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
208 problems.addProblem(message); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
209 } |
| 8898 | 210 } |
| 8886 | 211 |
| 8883 | 212 for (final BedHeightValue bedHeightValue : sortedValues) { |
| 8877 | 213 |
| 8883 | 214 final Double station = bedHeightValue.getStation(); |
| 215 if (station == null || station.isNaN()) | |
| 216 continue; | |
| 217 | |
| 218 final Double meanBedHeightDbl = bedHeightValue.getHeight(); | |
| 219 if (meanBedHeightDbl == null || meanBedHeightDbl.isNaN()) | |
| 220 continue; | |
| 221 | |
| 222 final double km = station; | |
| 223 final double meanBedHeight = meanBedHeightDbl; | |
| 224 | |
| 8894 | 225 if (!calcRange.containsDouble(km)) |
| 226 continue; | |
| 227 | |
| 8883 | 228 try { |
| 229 // FIXME: check out of range | |
| 230 final double wst = wstInterpolator.value(km); | |
| 231 | |
| 232 final double flowDepth = wst - meanBedHeight; | |
| 233 | |
| 234 // FIXME: piecewise constant interpolation? | |
| 235 // final double discharge = wstKms instanceof QKms ? ((QKms) wstKms).getQ(i) : Double.NaN; | |
| 8898 | 236 double discharge; |
| 237 if (qInterpolator != null) | |
| 238 discharge = qInterpolator.value(km); | |
| 239 else | |
| 240 discharge = Double.NaN; | |
| 8883 | 241 |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
242 // Calculate tkh |
| 8898 | 243 double tkh = 0; |
| 8901 | 244 if (doCalcTkh) { |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
245 double d50 = Double.NaN; |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
246 try { |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
247 d50 = bedMeasurementsFinder.findD50(km); |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
248 } catch (Exception e) { |
| 8898 | 249 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingD50", null, label); |
| 250 problems.addProblem(km, message); | |
| 8901 | 251 //FIXME: cumulate problems to one message? |
| 8898 | 252 } |
| 8901 | 253 if (!Double.isNaN(d50)) { |
| 254 if (flowVelocitiesFinder.findKmQValues(km, discharge)) { | |
| 255 tkh = calculateTkh(wst - meanBedHeight, flowVelocitiesFinder.getFindVmainFound(), d50, flowVelocitiesFinder.getFindTauFound()); | |
|
8911
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
256 if (!Double.isNaN(tkh) && (tkh < 0)) |
|
37ff7f435912
SINFO Flowdepth: more error checks, d50 interpolation, avoid negative tkh
mschaefer
parents:
8901
diff
changeset
|
257 tkh = 0; |
| 8901 | 258 /* log.debug(String.format("calculateTkh km %.3f q %.0f w %.2f mbh %.2f vm %.1f tau %.1f d50(mm) %.1f tkh(cm) %.1f", |
| 259 km, discharge, wst, meanBedHeight, flowVelocitiesFinder.getFindVmainFound(), flowVelocitiesFinder.getFindTauFound(), d50*1000, tkh)); */ | |
| 260 } | |
| 261 else { | |
| 262 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingQ", null, label); | |
| 263 problems.addProblem(km, message); | |
| 264 //FIXME: cumulate problems to one message? | |
| 265 } | |
| 8898 | 266 } |
| 267 else | |
| 268 tkh = Double.NaN; | |
| 269 } | |
| 270 | |
| 271 // Soil kind | |
| 8901 | 272 SoilKind kind = SoilKind.mobil; |
| 273 if (doCalcTkh) { | |
| 8898 | 274 try { |
| 275 kind = soilKindFinder.findSoilKind(km); | |
| 276 } catch (Exception e) { | |
| 277 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.missingSoilKind", null, label); | |
| 278 problems.addProblem(km, message); | |
| 8901 | 279 //FIXME: cumulate problems to one message? |
| 8898 | 280 } |
| 281 } | |
| 282 | |
| 283 /* // REMARK: bissl spielerei zum testen damit die sohlart nicht zu schnell wechselt | |
| 284 SoilKind kind; | |
| 8895 | 285 final boolean changeKind = Math.random() > 0.95; |
| 8886 | 286 if (changeKind) { |
| 287 switch (lastKind) { | |
| 288 case starr: | |
| 289 kind = SoilKind.mobil; | |
| 290 break; | |
| 291 case mobil: | |
| 292 default: | |
| 293 kind = SoilKind.starr; | |
| 294 break; | |
| 295 } | |
| 296 } else | |
| 297 kind = lastKind; | |
| 298 lastKind = kind; | |
| 8898 | 299 */ |
| 300 | |
| 301 //tkh = 100 + 10 * (Math.random() - 0.5); | |
| 8886 | 302 |
| 303 final double flowDepthTkh; | |
| 304 final double tkhUp; | |
| 305 final double tkhDown; | |
| 306 switch (kind) { | |
| 307 case starr: | |
| 308 flowDepthTkh = wst - (meanBedHeight + tkh / 100); | |
| 309 tkhUp = tkh; | |
| 310 tkhDown = 0; | |
| 311 break; | |
| 312 | |
| 313 case mobil: | |
| 314 default: | |
| 315 flowDepthTkh = wst - (meanBedHeight + tkh / 200); | |
| 316 tkhUp = tkh / 2; | |
| 317 tkhDown = -tkh / 2; | |
| 318 break; | |
| 319 } | |
| 320 | |
| 8883 | 321 // REMARK: access the location once only during calculation |
| 322 final String location = LocationProvider.getLocation(river.getName(), km); | |
| 323 | |
| 324 // REMARK: access the gauge once only during calculation | |
| 325 final Gauge gauge = findGauge(waterlevel, refGauge, gaugeIndex, km); | |
| 326 | |
| 327 final String gaugeLabel = gauge == null ? notinrange : gauge.getName(); | |
| 328 | |
| 8898 | 329 resultData.addRow(km, flowDepth, flowDepthTkh, kind, tkh, tkhUp, tkhDown, wst, discharge, wstLabel, |
| 330 gaugeLabel, meanBedHeight, bedHeightLabel, location); | |
| 8883 | 331 } |
| 332 catch (final FunctionEvaluationException e) { | |
| 333 /* should only happen if out of range */ | |
| 334 e.printStackTrace(); | |
| 335 /* simply ignore */ | |
| 336 } | |
| 337 | |
| 8877 | 338 } |
| 339 | |
| 340 return resultData; | |
| 341 } | |
| 342 | |
| 8883 | 343 /** |
| 8891 | 344 * Sohlbeschaffenheit (D50 Korndurchmesser aus Seddb) |
| 345 * Abhängig von Peiljahr | |
| 346 */ | |
| 8898 | 347 private BedQualityD50KmValueFinder loadBedMeasurements(final River river, final DoubleRange kmRange, final int soundingYear) { |
| 8891 | 348 |
| 349 /* construct valid measurement time range */ | |
| 350 final Calendar cal = Calendar.getInstance(); | |
| 351 cal.clear(); | |
| 352 | |
| 353 cal.set(soundingYear - VALID_BED_MEASUREMENT_YEARS, 0, 1); | |
| 354 final Date startTime = cal.getTime(); | |
| 355 | |
| 356 cal.set(soundingYear + VALID_BED_MEASUREMENT_YEARS, 11, 31); | |
| 357 final Date endTime = cal.getTime(); | |
| 358 | |
| 8898 | 359 final BedQualityD50KmValueFinder finder = new BedQualityD50KmValueFinder(); |
| 360 if (finder.loadValues(river, kmRange, new DateRange(startTime, endTime))) | |
| 361 return finder; | |
| 362 else | |
| 363 return null; | |
| 8891 | 364 } |
| 365 | |
| 366 /** | |
| 8883 | 367 * Checks the year difference between waterlevels and sounding, and issues a warning if too big. |
| 368 * | |
| 369 * Zeitraum Zeitliche Differenz [a] | |
| 370 * X ≥ 1998 ± 3 | |
| 371 * 1958 ≤ X < 1998 ± 6 | |
| 372 * 1918 ≤ X < 1958 ± 12 | |
| 373 * X < 1918 ± 25 | |
| 374 */ | |
| 8884 | 375 private void checkYearDifference(final String label, final WaterlevelData waterlevel, final BedHeight sounding, final Calculation problems) { |
| 8883 | 376 |
| 377 final Integer soundingYear = sounding.getYear(); | |
| 378 if (soundingYear == null) | |
| 379 return; | |
| 380 | |
| 381 final int wstYear = waterlevel.getYear(); | |
| 382 if (wstYear < 0) | |
| 383 return; | |
| 384 | |
| 385 final int maxDifference = getMaxDifferenceYears(soundingYear); | |
| 386 | |
| 387 final int difference = Math.abs(soundingYear - wstYear); | |
| 388 if (difference > maxDifference) { | |
| 8894 | 389 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.year_difference", null, label, wstYear, |
| 390 soundingYear); | |
| 8883 | 391 problems.addProblem(message); |
| 392 } | |
| 393 } | |
| 394 | |
| 395 private int getMaxDifferenceYears(final int year) { | |
| 396 | |
| 397 if (year < 1918) | |
| 398 return 25; | |
| 399 | |
| 400 if (1918 <= year && year < 1958) | |
| 401 return 12; | |
| 402 | |
| 403 if (1958 <= year && year < 1998) | |
| 404 return 6; | |
| 405 | |
| 406 /* >= 1998 */ | |
| 407 return 3; | |
| 408 } | |
| 409 | |
| 8884 | 410 private Gauge findGauge(final WaterlevelData waterlevel, final Gauge refGauge, final GaugeIndex gaugeIndex, final double km) { |
| 8882 | 411 |
| 412 // REMARK: using same logic as in WaterlevelExporter here | |
| 413 | |
| 414 final boolean showAllGauges = waterlevel.isShowAllGauges(); | |
| 415 | |
| 416 if (showAllGauges) | |
| 417 return gaugeIndex.findGauge(km); | |
| 418 | |
| 419 if (refGauge.getRange().contains(km)) | |
| 420 return refGauge; | |
| 421 | |
| 422 return null; | |
| 423 } | |
| 424 | |
| 425 /* Checks if the discretisation of the waterlevel exceeds 1000m */ | |
| 8894 | 426 |
| 427 private void checkWaterlevelDiscretisation(final WKms wstKms, final DoubleRange calcRange, final Calculation problems) { | |
| 428 | |
| 8882 | 429 final int size = wstKms.size(); |
| 430 for (int i = 0; i < size - 2; i++) { | |
| 431 final double kmPrev = wstKms.getKm(i); | |
| 432 final double kmNext = wstKms.getKm(i + 1); | |
| 433 | |
| 8894 | 434 /* only check if we are within the calculation range */ |
| 435 if (calcRange.overlapsRange(new DoubleRange(kmPrev, kmNext))) { | |
| 436 if (Math.abs(kmPrev - kmNext) > 1) { | |
| 437 final String label = wstKms.getName(); | |
| 8882 | 438 |
| 8894 | 439 final String message = Resources.getMsg(this.context.getMeta(), "sinfo_calc_flow_depth.warning.waterlevel_discretisation", null, label); |
| 440 problems.addProblem(kmPrev, message); | |
| 441 } | |
| 8882 | 442 } |
| 443 } | |
| 444 } | |
| 445 | |
| 8894 | 446 private BedHeight loadBedHeight(final String soundingId) { |
| 8877 | 447 |
| 8883 | 448 // REMARK: absolutely unbelievable.... |
| 8898 | 449 // The way how bed-heights (and other data too) is accessed is different for nearly every calculation-type |
| 8882 | 450 // throughout flys. |
| 8877 | 451 // The knowledge on how to parse the datacage-ids is spread through the complete code-base... |
| 452 | |
| 8882 | 453 // We use here the way on how bed-heights are accessed by the BedDifferenceAccess/BedDifferenceCalculation, but |
| 454 // this is plain random | |
| 8877 | 455 final String[] parts = soundingId.split(";"); |
| 456 | |
| 457 final BedHeightsArtifact artifact = (BedHeightsArtifact) RiverUtils.getArtifact(parts[0], this.context); | |
| 458 | |
| 459 final Integer bedheightId = artifact.getDataAsInteger("height_id"); | |
| 8883 | 460 // REMARK: this only works with type 'single'; unclear on how to distinguish from epoch data (or whatever the |
| 8882 | 461 // other type means) |
| 8877 | 462 // Luckily, the requirement is to only access 'single' data here. |
| 463 // final String bedheightType = artifact.getDataAsString("type"); | |
| 464 | |
| 8883 | 465 // REMARK: BedDifferences uses this, but we also need the metadata of the BedHeight |
| 466 // REMARK: second absolutely awful thing: BedHeight is a hibernate binding class, accessing the database via | |
| 8882 | 467 // hibernate stuff |
| 8877 | 468 // BedHeightFactory uses its own (direct) way of accessing the data, with its own implemented data classes. |
| 8882 | 469 // return BedHeightFactory.getHeight(bedheightType, bedheightId, from, to); |
| 8877 | 470 |
| 471 return BedHeight.getBedHeightById(bedheightId); | |
| 472 } | |
| 8898 | 473 |
| 474 /** | |
| 475 * Calculates a transport body height | |
| 476 * @param h flow depth in m | |
| 477 * @param vm flow velocity in m | |
| 478 * @param d50 grain diameter D50 in m (!) | |
| 479 * @param tau shear stress in N/m^2 | |
| 480 * @return transport body height in cm (!) | |
| 481 */ | |
| 482 private double calculateTkh(double h, double vm, double d50, double tau) { | |
| 483 final double PHYS_G = 9.81; | |
| 484 final double PHYS_SPECGRAV_S = 2.6; | |
| 485 final double PHYS_VELOCCOEFF_N = 6; | |
| 486 final double PHYS_FORMCOEFF_ALPHA = 0.7; | |
| 487 final double PHYS_VISCOSITY_NUE = 1.3e-6; | |
| 488 final double PHYS_GRAIN_DENSITY_RHOS = 2603; | |
| 489 final double PHYS_WATER_DENSITY_RHO = 999.97; | |
| 490 | |
| 491 final double froude = vm / Math.sqrt(PHYS_G * h); | |
| 492 final double partReynolds = Math.sqrt((PHYS_SPECGRAV_S - 1) * PHYS_G * d50) / PHYS_VISCOSITY_NUE * d50; | |
| 493 final double critShields = 0.22 * Math.pow(partReynolds, -0.6) + 0.06 * Math.pow(10, 7.7 * Math.pow(partReynolds, -0.6)); | |
| 494 final double critTau = critShields * (PHYS_GRAIN_DENSITY_RHOS - PHYS_WATER_DENSITY_RHO) * PHYS_G * d50; | |
| 495 return 100 * h * (1 - Math.pow(froude, 2)) / (2 * PHYS_VELOCCOEFF_N * PHYS_FORMCOEFF_ALPHA) * (1 - critTau / tau); | |
| 496 } | |
| 8884 | 497 } |