Mercurial > dive4elements > river
comparison flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/WstValueTable.java @ 655:913b52064449
Refactored version of "Berechnung 4"
flys-artifacts/trunk@2053 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Sun, 05 Jun 2011 18:24:46 +0000 |
parents | 44175d4720f8 |
children | c501f27c1f71 3dc61e00385e |
comparison
equal
deleted
inserted
replaced
654:bbc966c81809 | 655:913b52064449 |
---|---|
1 package de.intevation.flys.artifacts.model; | 1 package de.intevation.flys.artifacts.model; |
2 | 2 |
3 import java.io.Serializable; | 3 import java.io.Serializable; |
4 | 4 |
5 | 5 import de.intevation.flys.artifacts.math.Linear; |
6 import de.intevation.flys.artifacts.math.LinearRemap; | 6 import de.intevation.flys.artifacts.math.Function; |
7 | 7 |
8 import java.util.Arrays; | 8 import java.util.Arrays; |
9 import java.util.ArrayList; | 9 import java.util.ArrayList; |
10 import java.util.List; | 10 import java.util.List; |
11 import java.util.Collections; | 11 import java.util.Collections; |
107 double km, | 107 double km, |
108 double [] iqs, | 108 double [] iqs, |
109 double [] ows, | 109 double [] ows, |
110 WstValueTable table | 110 WstValueTable table |
111 ) { | 111 ) { |
112 double kmWeight = factor(km, this.km, other.km); | 112 double kmWeight = Linear.factor(km, this.km, other.km); |
113 | 113 |
114 QPosition qPosition = new QPosition(); | 114 QPosition qPosition = new QPosition(); |
115 | 115 |
116 for (int i = 0; i < iqs.length; ++i) { | 116 for (int i = 0; i < iqs.length; ++i) { |
117 if (table.getQPosition(km, iqs[i], qPosition) != null) { | 117 if (table.getQPosition(km, iqs[i], qPosition) != null) { |
118 double wt = getW(qPosition); | 118 double wt = getW(qPosition); |
119 double wo = other.getW(qPosition); | 119 double wo = other.getW(qPosition); |
120 ows[i] = weight(kmWeight, wt, wo); | 120 ows[i] = Linear.weight(kmWeight, wt, wo); |
121 } | 121 } |
122 else { | 122 else { |
123 ows[i] = Double.NaN; | 123 ows[i] = Double.NaN; |
124 } | 124 } |
125 } | 125 } |
135 | 135 |
136 if (W < 1) { | 136 if (W < 1) { |
137 return new double[2][0]; | 137 return new double[2][0]; |
138 } | 138 } |
139 | 139 |
140 double factor = factor(km, this.km, other.km); | 140 double factor = Linear.factor(km, this.km, other.km); |
141 | 141 |
142 double [] splineQ = new double[W]; | 142 double [] splineQ = new double[W]; |
143 double [] splineW = new double[W]; | 143 double [] splineW = new double[W]; |
144 | 144 |
145 double minQ = Double.MAX_VALUE; | 145 double minQ = Double.MAX_VALUE; |
146 double maxQ = -Double.MAX_VALUE; | 146 double maxQ = -Double.MAX_VALUE; |
147 | 147 |
148 for (int i = 0; i < W; ++i) { | 148 for (int i = 0; i < W; ++i) { |
149 double wws = weight(factor, ws[i], other.ws[i]); | 149 double wws = Linear.weight(factor, ws[i], other.ws[i]); |
150 double wqs = weight( | 150 double wqs = Linear.weight( |
151 factor, | 151 factor, |
152 table.getQIndex(i, km), | 152 table.getQIndex(i, km), |
153 table.getQIndex(i, other.km)); | 153 table.getQIndex(i, other.km)); |
154 | 154 |
155 splineW[i] = wws; | 155 splineW[i] = wws; |
229 int index = qPosition.index; | 229 int index = qPosition.index; |
230 double weight = qPosition.weight; | 230 double weight = qPosition.weight; |
231 | 231 |
232 return weight == 1.0 | 232 return weight == 1.0 |
233 ? ws[index] | 233 ? ws[index] |
234 : weight(weight, ws[index-1], ws[index]); | 234 : Linear.weight(weight, ws[index-1], ws[index]); |
235 } | 235 } |
236 | 236 |
237 public double getW( | 237 public double getW( |
238 Row other, | 238 Row other, |
239 double km, | 239 double km, |
240 QPosition qPosition | 240 QPosition qPosition |
241 ) { | 241 ) { |
242 double kmWeight = factor(km, this.km, other.km); | 242 double kmWeight = Linear.factor(km, this.km, other.km); |
243 | 243 |
244 int index = qPosition.index; | 244 int index = qPosition.index; |
245 double weight = qPosition.weight; | 245 double weight = qPosition.weight; |
246 | 246 |
247 double tw, ow; | 247 double tw, ow; |
249 if (weight == 1.0) { | 249 if (weight == 1.0) { |
250 tw = ws[index]; | 250 tw = ws[index]; |
251 ow = other.ws[index]; | 251 ow = other.ws[index]; |
252 } | 252 } |
253 else { | 253 else { |
254 tw = weight(weight, ws[index-1], ws[index]); | 254 tw = Linear.weight(weight, ws[index-1], ws[index]); |
255 ow = weight(weight, other.ws[index-1], other.ws[index]); | 255 ow = Linear.weight(weight, other.ws[index-1], other.ws[index]); |
256 } | 256 } |
257 | 257 |
258 return weight(kmWeight, tw, ow); | 258 return Linear.weight(kmWeight, tw, ow); |
259 } | 259 } |
260 } // class Row | 260 } // class Row |
261 | 261 |
262 protected List<Row> rows; | 262 protected List<Row> rows; |
263 | 263 |
336 Row r2 = rows.get(rowIndex); | 336 Row r2 = rows.get(rowIndex); |
337 | 337 |
338 return r1.interpolateWQ(r2, km, steps, this); | 338 return r1.interpolateWQ(r2, km, steps, this); |
339 } | 339 } |
340 | 340 |
341 public void interpolate( | 341 public boolean interpolate( |
342 double [] kms, | 342 double km, |
343 double [] ws, | 343 double [] out, |
344 double [] qs, | 344 QPosition qPosition, |
345 QPosition qPosition, | 345 Function qFunction |
346 LinearRemap remap | |
347 ) { | |
348 interpolate(kms, ws, qs, qPosition, remap, 0, kms.length); | |
349 } | |
350 | |
351 public void interpolate( | |
352 double [] kms, | |
353 double [] ws, | |
354 double [] qs, | |
355 QPosition qPosition, | |
356 LinearRemap remap, | |
357 int startIndex, | |
358 int length | |
359 ) { | 346 ) { |
360 int R1 = rows.size()-1; | 347 int R1 = rows.size()-1; |
361 | 348 |
362 Row kmKey = new Row(); | 349 out[1] = qFunction.value(getQ(qPosition, km)); |
350 | |
351 if (Double.isNaN(out[1])) { | |
352 return false; | |
353 } | |
363 | 354 |
364 QPosition nPosition = new QPosition(); | 355 QPosition nPosition = new QPosition(); |
365 | 356 if (getQPosition(km, out[1], nPosition) == null) { |
366 for (int i = startIndex, end = startIndex+length; i < end; ++i) { | 357 return false; |
367 kmKey.km = kms[i]; | 358 } |
368 | 359 |
369 qs[i] = remap.eval(kms[i], getQ(qPosition, kms[i])); | 360 int rowIndex = Collections.binarySearch(rows, new Row(km)); |
370 | 361 |
371 if (getQPosition(kms[i], qs[i], nPosition) == null) { | 362 if (rowIndex >= 0) { |
372 ws[i] = Double.NaN; | 363 // direct row match |
373 continue; | 364 out[0] = rows.get(rowIndex).getW(nPosition); |
374 } | 365 return !Double.isNaN(out[0]); |
375 | 366 } |
376 int rowIndex = Collections.binarySearch(rows, kmKey); | 367 |
377 | 368 rowIndex = -rowIndex -1; |
378 if (rowIndex >= 0) { | 369 |
379 // direct row match | 370 if (rowIndex < 1 || rowIndex >= R1) { |
380 ws[i] = rows.get(rowIndex).getW(nPosition); | 371 // do not extrapolate |
381 continue; | 372 return false; |
382 } | 373 } |
383 | 374 |
384 rowIndex = -rowIndex -1; | 375 Row r1 = rows.get(rowIndex-1); |
385 | 376 Row r2 = rows.get(rowIndex); |
386 if (rowIndex < 1 || rowIndex >= R1) { | 377 out[0] = r1.getW(r2, km, nPosition); |
387 // do not extrapolate | 378 |
388 ws[i] = Double.NaN; | 379 return !Double.isNaN(out[0]); |
389 continue; | |
390 } | |
391 Row r1 = rows.get(rowIndex-1); | |
392 Row r2 = rows.get(rowIndex); | |
393 ws[i] = r1.getW(r2, kms[i], nPosition); | |
394 } | |
395 } | 380 } |
396 | 381 |
397 public QPosition interpolate( | 382 public QPosition interpolate( |
398 double q, | 383 double q, |
399 double referenceKm, | 384 double referenceKm, |
478 double qMin, qMax; | 463 double qMin, qMax; |
479 if (qLast < qCurrent) { qMin = qLast; qMax = qCurrent; } | 464 if (qLast < qCurrent) { qMin = qLast; qMax = qCurrent; } |
480 else { qMin = qCurrent; qMax = qLast; } | 465 else { qMin = qCurrent; qMax = qLast; } |
481 | 466 |
482 if (q > qMin && q < qMax) { | 467 if (q > qMin && q < qMax) { |
483 double weight = factor(q, qLast, qCurrent); | 468 double weight = Linear.factor(q, qLast, qCurrent); |
484 return qPosition.set(i, weight); | 469 return qPosition.set(i, weight); |
485 } | 470 } |
486 qLast = qCurrent; | 471 qLast = qCurrent; |
487 } | 472 } |
488 | 473 |
500 if (weight == 1d) { | 485 if (weight == 1d) { |
501 return columns[index].getQRangeTree().findQ(km); | 486 return columns[index].getQRangeTree().findQ(km); |
502 } | 487 } |
503 double q1 = columns[index-1].getQRangeTree().findQ(km); | 488 double q1 = columns[index-1].getQRangeTree().findQ(km); |
504 double q2 = columns[index ].getQRangeTree().findQ(km); | 489 double q2 = columns[index ].getQRangeTree().findQ(km); |
505 return weight(weight, q1, q2); | 490 return Linear.weight(weight, q1, q2); |
506 } | 491 } |
507 | 492 |
508 public static final double linear( | |
509 double x, | |
510 double x1, double x2, | |
511 double y1, double y2 | |
512 ) { | |
513 // y1 = m*x1 + b | |
514 // y2 = m*x2 + b | |
515 // y2 - y1 = m*(x2 - x1) | |
516 // m = (y2 - y1)/(x2 - x1) # x2 != x1 | |
517 // b = y1 - m*x1 | |
518 | |
519 if (x1 == x2) { | |
520 return 0.5*(y1 + y2); | |
521 } | |
522 double m = (y2 - y1)/(x2 - x1); | |
523 double b = y1 - m*x1; | |
524 return x*m + b; | |
525 } | |
526 | |
527 public static final double factor(double x, double p1, double p2) { | |
528 // 0 = m*p1 + b <=> b = -m*p1 | |
529 // 1 = m*p2 + b | |
530 // 1 = m*(p2 - p1) | |
531 // m = 1/(p2 - p1) # p1 != p2 | |
532 // f(x) = x/(p2-p1) - p1/(p2-p1) <=> (x-p1)/(p2-p1) | |
533 | |
534 return p1 == p2 ? 0.0 : (x-p1)/(p2-p1); | |
535 } | |
536 | |
537 public static final double weight(double factor, double a, double b) { | |
538 //return (1.0-factor)*a + factor*b; | |
539 return a + factor*(b-a); | |
540 } | |
541 } | 493 } |
542 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : | 494 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : |