Mercurial > dive4elements > river
comparison flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/minfo/BedQualityCalculation.java @ 3769:728ecd2afa20
Implemented bed quality calculation in minfo module.
flys-artifacts/trunk@5474 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author | Raimund Renkert <raimund.renkert@intevation.de> |
---|---|
date | Fri, 14 Sep 2012 14:20:42 +0000 |
parents | 312870fded7e |
children | 5a8f8fd5310c |
comparison
equal
deleted
inserted
replaced
3768:00aafe1fedd7 | 3769:728ecd2afa20 |
---|---|
1 package de.intevation.flys.artifacts.model.minfo; | 1 package de.intevation.flys.artifacts.model.minfo; |
2 | |
3 import gnu.trove.TDoubleArrayList; | |
2 | 4 |
3 import java.util.LinkedList; | 5 import java.util.LinkedList; |
4 import java.util.List; | 6 import java.util.List; |
7 import java.util.Map; | |
5 | 8 |
6 import org.apache.log4j.Logger; | 9 import org.apache.log4j.Logger; |
7 | 10 |
8 import de.intevation.flys.artifacts.access.BedQualityAccess; | 11 import de.intevation.flys.artifacts.access.BedQualityAccess; |
9 import de.intevation.flys.artifacts.model.Calculation; | 12 import de.intevation.flys.artifacts.model.Calculation; |
18 .getLogger(BedQualityCalculation.class); | 21 .getLogger(BedQualityCalculation.class); |
19 | 22 |
20 protected String river; | 23 protected String river; |
21 protected double from; | 24 protected double from; |
22 protected double to; | 25 protected double to; |
26 protected List<String> bedDiameter; | |
27 protected List<String> bedloadDiameter; | |
23 protected List<DateRange> ranges; | 28 protected List<DateRange> ranges; |
24 | 29 |
25 public BedQualityCalculation() { | 30 public BedQualityCalculation() { |
26 } | 31 } |
27 | 32 |
29 logger.info("BedQualityCalculation.calculate"); | 34 logger.info("BedQualityCalculation.calculate"); |
30 | 35 |
31 String river = access.getRiver(); | 36 String river = access.getRiver(); |
32 Double from = access.getFrom(); | 37 Double from = access.getFrom(); |
33 Double to = access.getTo(); | 38 Double to = access.getTo(); |
39 List<String> bedDiameter = access.getBedDiameter(); | |
40 List<String> bedloadDiameter = access.getBedloadDiameter(); | |
34 List<DateRange> ranges = access.getDateRanges(); | 41 List<DateRange> ranges = access.getDateRanges(); |
35 | 42 |
36 if (river == null) { | 43 if (river == null) { |
37 // TODO: i18n | 44 // TODO: i18n |
38 addProblem("minfo.missing.river"); | 45 addProblem("minfo.missing.river"); |
56 if (!hasProblems()) { | 63 if (!hasProblems()) { |
57 this.river = river; | 64 this.river = river; |
58 this.from = from; | 65 this.from = from; |
59 this.to = to; | 66 this.to = to; |
60 this.ranges = ranges; | 67 this.ranges = ranges; |
68 this.bedDiameter = bedDiameter; | |
69 this.bedloadDiameter = bedloadDiameter; | |
61 | 70 |
62 SedDBSessionHolder.acquire(); | 71 SedDBSessionHolder.acquire(); |
63 try { | 72 try { |
64 return internalCalculate(); | 73 return internalCalculate(); |
65 } | 74 } |
74 protected CalculationResult internalCalculate() { | 83 protected CalculationResult internalCalculate() { |
75 | 84 |
76 List<BedQualityResult> results = new LinkedList<BedQualityResult>(); | 85 List<BedQualityResult> results = new LinkedList<BedQualityResult>(); |
77 // Calculate for all time periods. | 86 // Calculate for all time periods. |
78 for (DateRange dr : ranges) { | 87 for (DateRange dr : ranges) { |
79 QualityMeasurements bedMeasurements = QualityMeasurementFactory | 88 QualityMeasurements loadMeasurements = |
80 .getBedMeasurements(river, from, to, dr.getFrom(), dr.getTo()); | 89 QualityMeasurementFactory.getBedMeasurements( |
81 QualityMeasurements loadMeasurements = QualityMeasurementFactory | 90 river, |
82 .getBedMeasurements(river, from, to, dr.getFrom(), dr.getTo()); | 91 from, |
83 | 92 to, |
84 BedQualityResult bedResult = calculateBed(bedMeasurements); | 93 dr.getFrom(), |
85 BedQualityResult loadResult = calculateBedload(loadMeasurements); | 94 dr.getTo()); |
86 results.add(bedResult); | 95 QualityMeasurements bedMeasurements = |
87 results.add(loadResult); | 96 QualityMeasurementFactory.getBedMeasurements( |
97 river, | |
98 from, | |
99 to, | |
100 dr.getFrom(), | |
101 dr.getTo()); | |
102 | |
103 for (String bd : bedDiameter) { | |
104 BedQualityResult bedResult = | |
105 calculateBed(bedMeasurements, bd, dr); | |
106 | |
107 //Avoid adding empty result sets. | |
108 if (!bedResult.isEmpty()) { | |
109 results.add(bedResult); | |
110 } | |
111 } | |
112 for (String bld : bedloadDiameter) { | |
113 BedQualityResult loadResult = | |
114 calculateBedload(loadMeasurements, bld, dr); | |
115 results.add(loadResult); | |
116 } | |
88 } | 117 } |
89 | 118 |
90 return new CalculationResult( | 119 return new CalculationResult( |
91 results.toArray(new BedQualityResult[results.size()]), this); | 120 results.toArray(new BedQualityResult[results.size()]), this); |
92 } | 121 } |
93 | 122 |
94 protected BedQualityResult calculateBed(QualityMeasurements qm) { | 123 protected BedQualityResult calculateBed( |
95 // TODO | 124 QualityMeasurements qm, |
96 return new BedQualityResult(); | 125 String diameter, |
97 } | 126 DateRange range |
98 | 127 ) { |
99 protected BedQualityResult calculateBedload(QualityMeasurements qm) { | 128 List<Double> kms = qm.getKms(); |
100 // TODO | 129 TDoubleArrayList location = new TDoubleArrayList(); |
101 return new BedQualityResult(); | 130 TDoubleArrayList avDiameterCap = new TDoubleArrayList(); |
131 TDoubleArrayList avDiameterSub = new TDoubleArrayList(); | |
132 TDoubleArrayList porosityCap = new TDoubleArrayList(); | |
133 TDoubleArrayList porositySub = new TDoubleArrayList(); | |
134 TDoubleArrayList densityCap = new TDoubleArrayList(); | |
135 TDoubleArrayList densitySub = new TDoubleArrayList(); | |
136 for (double km : kms) { | |
137 //Filter cap and sub measurements. | |
138 QualityMeasurements capFiltered = filterCapMeasurements(qm); | |
139 QualityMeasurements subFiltered = filterSubMeasurements(qm); | |
140 | |
141 List<QualityMeasurement> cm = capFiltered.getMeasurements(km); | |
142 List<QualityMeasurement> sm = subFiltered.getMeasurements(km); | |
143 | |
144 double avCap = calculateAverage(cm, diameter); | |
145 double avSub = calculateAverage(sm, diameter); | |
146 double[] pCap = calculatePorosity(capFiltered, km, diameter); | |
147 double[] pSub = calculatePorosity(subFiltered, km, diameter); | |
148 double[] dCap = calculateDensity(capFiltered, pCap); | |
149 double[] dSub = calculateDensity(subFiltered, pSub); | |
150 | |
151 double pCapRes = 0d; | |
152 double pSubRes = 0d; | |
153 double dCapRes = 0d; | |
154 double dSubRes = 0d; | |
155 for (int i = 0; i < pCap.length; i++) { | |
156 pCapRes += pCap[i]; | |
157 dCapRes += dCap[i]; | |
158 } | |
159 for (int i = 0; i < pSub.length; i++) { | |
160 pSubRes += pSub[i]; | |
161 dSubRes += dSub[i]; | |
162 } | |
163 location.add(km); | |
164 avDiameterCap.add(avCap); | |
165 avDiameterSub.add(avSub); | |
166 porosityCap.add((pCapRes / pCap.length) * 100 ); | |
167 porositySub.add((pSubRes / pSub.length) * 100); | |
168 densityCap.add((dCapRes / dCap.length) / 1000); | |
169 densitySub.add((dSubRes / dSub.length) / 1000); | |
170 } | |
171 return new BedBedQualityResult( | |
172 diameter, | |
173 avDiameterCap, | |
174 avDiameterSub, | |
175 location, | |
176 range, | |
177 porosityCap, | |
178 porositySub, | |
179 densityCap, | |
180 densitySub); | |
181 } | |
182 | |
183 private double[] calculateDensity( | |
184 QualityMeasurements capFiltered, | |
185 double[] porosity | |
186 ) { | |
187 double[] density = new double[porosity.length]; | |
188 for (int i = 0; i < porosity.length; i++) { | |
189 density[i] = (1 - porosity[i]) * 2650; | |
190 } | |
191 return density; | |
192 } | |
193 | |
194 private double[] calculatePorosity( | |
195 QualityMeasurements capFiltered, | |
196 double km, | |
197 String diameter | |
198 ) { | |
199 List<QualityMeasurement> list = capFiltered.getMeasurements(km); | |
200 double[] results = new double[list.size()]; | |
201 int i = 0; | |
202 for (QualityMeasurement qm : list) { | |
203 double deviation = calculateDeviation(qm); | |
204 double p = calculateP(qm); | |
205 double porosity = 0.353 - 0.068 * deviation + 0.146 * p; | |
206 results[i] = porosity; | |
207 i++; | |
208 } | |
209 | |
210 return results; | |
211 } | |
212 | |
213 protected BedQualityResult calculateBedload( | |
214 QualityMeasurements qm, | |
215 String diameter, | |
216 DateRange range | |
217 ) { | |
218 List<Double> kms = qm.getKms(); | |
219 TDoubleArrayList location = new TDoubleArrayList(); | |
220 TDoubleArrayList avDiameter = new TDoubleArrayList(); | |
221 for (double km : kms) { | |
222 List<QualityMeasurement> measurements = qm.getMeasurements(km); | |
223 double mid = calculateAverage(measurements, diameter); | |
224 location.add(km); | |
225 avDiameter.add(mid); | |
226 } | |
227 return new BedLoadBedQualityResult( | |
228 diameter, | |
229 avDiameter, | |
230 location, | |
231 range); | |
232 } | |
233 | |
234 protected double calculateAverage( | |
235 List<QualityMeasurement> list, | |
236 String diameter | |
237 ) { | |
238 double av = 0; | |
239 for (QualityMeasurement qm : list) { | |
240 av += qm.getDiameter(diameter); | |
241 } | |
242 return av/list.size(); | |
243 } | |
244 | |
245 protected QualityMeasurements filterCapMeasurements( | |
246 QualityMeasurements qms | |
247 ) { | |
248 List<QualityMeasurement> result = new LinkedList<QualityMeasurement>(); | |
249 for (QualityMeasurement qm : qms.getMeasurements()) { | |
250 if (qm.getDepth1() == 0d && qm.getDepth2() <= 0.3) { | |
251 result.add(qm); | |
252 } | |
253 } | |
254 return new QualityMeasurements(result); | |
255 } | |
256 | |
257 protected QualityMeasurements filterSubMeasurements( | |
258 QualityMeasurements qms | |
259 ) { | |
260 List<QualityMeasurement> result = new LinkedList<QualityMeasurement>(); | |
261 for (QualityMeasurement qm : qms.getMeasurements()) { | |
262 if (qm.getDepth1() > 0d && qm.getDepth2() <= 0.5) { | |
263 result.add(qm); | |
264 } | |
265 } | |
266 return new QualityMeasurements(result); | |
267 } | |
268 | |
269 public double calculateDeviation(QualityMeasurement qm) { | |
270 Map<String, Double> dm = qm.getAllDiameter(); | |
271 double phiM = 0; | |
272 double[] phis = new double[dm.size()]; | |
273 double[] ps = new double[dm.size()]; | |
274 int i = 0; | |
275 for (String key : dm.keySet()) { | |
276 double d = dm.get(key); | |
277 double phi = -Math.log(d)/Math.log(2); | |
278 phis[i] = phi; | |
279 double p = calculateWeight(qm, key); | |
280 ps[i] = p; | |
281 phiM += phi * p; | |
282 i++; | |
283 } | |
284 | |
285 double sig = 0d; | |
286 for (i = 0; i < dm.size(); i++) { | |
287 sig += ps[i] * Math.exp(phis[i] - phiM); | |
288 } | |
289 double deviation = Math.sqrt(sig); | |
290 return deviation; | |
291 } | |
292 | |
293 protected double calculateP(QualityMeasurement qm) { | |
294 return calculateWeight(qm, "dmin"); | |
295 } | |
296 | |
297 public double calculateWeight(QualityMeasurement qm, String diameter) { | |
298 Map<String, Double> dm = qm.getAllDiameter(); | |
299 double value = qm.getDiameter(diameter); | |
300 | |
301 double sum = 0d; | |
302 for (Double d : dm.values()) { | |
303 sum =+ d.doubleValue(); | |
304 } | |
305 double weight = sum/100*value; | |
306 return weight; | |
102 } | 307 } |
103 } | 308 } |