changeset 3099:f516ac26f4db direct-match-nash-sutcliffe

"Sharp" match for predicted and measured values in nash sutcliffe.
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Wed, 24 Apr 2019 22:26:16 +0200
parents 8f2ac24b0cb3
children d79e6045452e
files pkg/common/nashsutcliffe.go pkg/controllers/gauges.go
diffstat 2 files changed, 25 insertions(+), 68 deletions(-) [+]
line wrap: on
line diff
--- a/pkg/common/nashsutcliffe.go	Wed Apr 24 13:34:48 2019 +0200
+++ b/pkg/common/nashsutcliffe.go	Wed Apr 24 22:26:16 2019 +0200
@@ -32,66 +32,17 @@
 	})
 }
 
-var utc0 = time.Unix(0, 0).UTC()
-
-func (mvs TimedValues) Interpolate(when time.Time) (float64, bool) {
-
-	for lo, hi := 0, len(mvs)-2; lo <= hi; {
-		mid := lo + (hi-lo)/2
-		m1 := &mvs[mid]
-		m2 := &mvs[mid+1]
-		if m2.When.Before(when) {
-			hi = mid - 1
-			continue
-		}
-		if m1.When.After(when) {
-			lo = mid + 1
-			continue
-		}
-		if m1.When.Equal(when) {
-			//log.Printf("matches first: %f\n", m1.Value)
-			return m1.Value, true
-		}
-		if m2.When.Equal(when) {
-			//log.Printf("matches second: %f\n", m1.Value)
-			return m2.Value, true
-		}
+func epsEquals(a, b time.Time) bool {
+	d := a.Sub(b)
+	return -10*time.Millisecond < d && d < 10*time.Millisecond
+}
 
-		// f(m1.When) = m1.Value
-		// f(m2.When) = m2.Value
-		// m1.Value = m1.When*a + b <=> b = m1.Value - m1.When*a
-		// m2.Value = m2.When*a + b
-		// m1.Value - m2.Value = a*(m1.When - m2-When)
-		// a = (m1.Value - m2.Value)/(m1.When - m2.When) for m1.When != m2.When
-
-		if m1.When.Equal(m2.When) {
-			return (m1.Value + m2.Value) * 0.5, true
-		}
-
-		a := (m1.Value - m2.Value) / m1.When.Sub(m2.When).Seconds()
-		b := m1.Value - m1.When.Sub(utc0).Seconds()*a
-		m := when.Sub(utc0).Seconds()*a + b
-
-		//log.Printf("%f %f %f\n", m1.Value, m, m2.Value)
-		return m, true
-	}
-
-	if l := len(mvs); l > 0 {
-		if when.Equal(mvs[0].When) {
-			//log.Printf("at start\n")
-			return mvs[0].Value, true
-		}
-		if !when.Before(mvs[l-1].When) {
-			//log.Printf("after end\n")
-			return mvs[l-1].Value, true
+func (mvs TimedValues) Find(when time.Time) (float64, bool) {
+	for i := range mvs {
+		if epsEquals(when, mvs[i].When) {
+			return mvs[i].Value, true
 		}
 	}
-
-	//if len(mvs) > 0 {
-	//	log.Printf("does not match %v %v %v\n",
-	//		mvs[0].When, mvs[len(mvs)-1].When.Sub(mvs[0].When), when)
-	//}
-
 	return 0, false
 }
 
--- a/pkg/controllers/gauges.go	Wed Apr 24 13:34:48 2019 +0200
+++ b/pkg/controllers/gauges.go	Wed Apr 24 22:26:16 2019 +0200
@@ -491,16 +491,22 @@
 		}
 
 		if predicted {
-			// current prediction same as last two?
-			if l := len(current.predicted); l > 1 && current.predicted[l-1].Value == value && current.predicted[l-2].Value == value {
-				// only update last time.
-				current.predicted[l-1].When = issueDate
-			} else {
-				current.predicted = append(
-					current.predicted,
-					common.TimedValue{When: issueDate, Value: value},
-				)
-			}
+			current.predicted = append(
+				current.predicted,
+				common.TimedValue{When: issueDate, Value: value},
+			)
+			/*
+				// current prediction same as last two?
+				if l := len(current.predicted); l > 1 && current.predicted[l-1].Value == value && current.predicted[l-2].Value == value {
+					// only update last time.
+					current.predicted[l-1].When = issueDate
+				} else {
+					current.predicted = append(
+						current.predicted,
+						common.TimedValue{When: issueDate, Value: value},
+					)
+				}
+			*/
 		} else {
 			current.observed = value
 		}
@@ -582,7 +588,7 @@
 
 		for j := range values {
 			when := values[j].when.Add(delta)
-			if p, ok := values[j].predicted.Interpolate(when); ok {
+			if p, ok := values[j].predicted.Find(when); ok {
 				predicted = append(predicted, p)
 				observed = append(observed, values[j].observed)
 			}