Refactoring

continuous/FrenzelPompe.go

@@ -13,7 +13,7 @@ import (
 // Phys. Rev. Lett., 99:204101, Nov 2007.
 func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, eta bool) (r float64) {
 
-	hk := harmonic(k - 1)
+	hk := Harmonic(k - 1)
 
 	var bar *pb.ProgressBar
 
@@ -26,13 +26,13 @@ func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, et
 		epsilon := fpGetEpsilon(k, xyz[t], xyz, xIndices, yIndices, zIndices)
 
 		cNxz := fpCount2(epsilon, xyz[t], xyz, xIndices, zIndices)
-		hNxz := harmonic(cNxz)
+		hNxz := Harmonic(cNxz)
 
 		cNyz := fpCount2(epsilon, xyz[t], xyz, yIndices, zIndices)
-		hNyz := harmonic(cNyz)
+		hNyz := Harmonic(cNyz)
 
 		cNz := fpCount1(epsilon, xyz[t], xyz, zIndices)
-		hNz := harmonic(cNz)
+		hNz := Harmonic(cNz)
 
 		r += hNxz + hNyz - hNz
 
@@ -55,9 +55,9 @@ func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, et
 // fpMaxNorm3 computes the max-norm of two 3-dimensional vectors
 //   maxnorm(a,b) = max( |a[0] - b[0]|, |a[1] - b[1]|, |a[2] - b[2]|)
 func fpMaxNorm3(a, b []float64, xIndices, yIndices, zIndices []int) float64 {
-	xDist := distance(a, b, xIndices)
-	yDist := distance(a, b, yIndices)
-	zDist := distance(a, b, zIndices)
+	xDist := Distance(a, b, xIndices)
+	yDist := Distance(a, b, yIndices)
+	zDist := Distance(a, b, zIndices)
 	return math.Max(xDist, math.Max(yDist, zDist))
 }
 
@@ -92,8 +92,8 @@ func fpCount2(epsilon float64, xyz []float64, data [][]float64, xIndices, yIndic
 }
 
 func fpMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
-	xDist := distance(a, b, xIndices)
-	yDist := distance(a, b, yIndices)
+	xDist := Distance(a, b, xIndices)
+	yDist := Distance(a, b, yIndices)
 	return math.Max(xDist, yDist)
 }
 
@@ -102,7 +102,7 @@ func fpMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
 func fpCount1(epsilon float64, xyz []float64, data [][]float64, zIndices []int) (c int) {
 	c = -1 // because we will also count xyz[t] vs. xyz[t]
 	for t := 0; t < len(data); t++ {
-		if distance(xyz, data[t], zIndices) < epsilon {
+		if Distance(xyz, data[t], zIndices) < epsilon {
 			c++
 		}
 	}

continuous/Functions.go

@@ -2,17 +2,20 @@ package continuous
 
 import "math"
 
-func distance(a, b []float64, indices []int) float64 {
+// Distance calculates the distance between to vectors,
+// given the set of indices
+func Distance(a, b []float64, indices []int) float64 {
 	d := 0.0
 	for _, v := range indices {
 		d += (a[v] - b[v]) * (a[v] - b[v])
 	}
 	return math.Sqrt(d)
 }
 
-func harmonic(n int) (r float64) {
-	// harmonic(1) = -C, see A. Kraskov, H. Stoeogbauer, and P. Grassberger.
-	// Estimating mutual information. Phys. Rev. E, 69:066138, Jun 2004.
+// Harmonic calculates the harmonic according to
+// A. Kraskov, H. Stoeogbauer, and P. Grassberger.
+// Estimating mutual information. Phys. Rev. E, 69:066138, Jun 2004.
+func Harmonic(n int) (r float64) {
 	if n == 0 {
 		return
 	}

continuous/Functions_test.go

@@ -11,7 +11,7 @@ func TestDistanceZero(t *testing.T) {
 	a := []float64{1.0, 2.0, 3.0}
 	b := []float64{1.0, 2.0, 3.0}
 
-	if diff := distance(a, b, []int{0, 1, 2}); diff > 0.0 {
+	if diff := Distance(a, b, []int{0, 1, 2}); diff > 0.0 {
 		t.Errorf(fmt.Sprintf("Distance should be zero but it is %f", diff))
 	}
 }
@@ -21,15 +21,15 @@ func TestDistanceNotZero(t *testing.T) {
 	b := []float64{2.0, 4.0, 6.0}
 
 	dist1 := math.Sqrt(1.0 + 4.0 + 9.0)
-	dist2 := distance(a, b, []int{0, 1, 2})
+	dist2 := Distance(a, b, []int{0, 1, 2})
 
 	if math.Abs(dist1-dist2) > 0.00001 {
 		t.Errorf(fmt.Sprintf("Distance should be %f but it is %f", dist1, dist2))
 	}
 }
 
 func TestHarmonic1(t *testing.T) {
-	h := harmonic(1)
+	h := Harmonic(1)
 	r := -0.5772156649
 
 	if math.Abs(h-r) > 0.00001 {
@@ -39,14 +39,14 @@ func TestHarmonic1(t *testing.T) {
 }
 
 func TestHarmonic5(t *testing.T) {
-	h := harmonic(5)
+	h := Harmonic(5)
 	r := -1.0/2.0 - 1.0/3.0 - 1.0/4.0 - 1.0/5.0 - 0.5772156649
 
 	if math.Abs(h-r) > 0.00001 {
 		t.Errorf(fmt.Sprintf("Harmonic(5) should be %f but it is %f", r, h))
 	}
 
-	h2 := harmonic(0)
+	h2 := Harmonic(0)
 	if h2 != 0.0 {
 		t.Errorf(fmt.Sprintf("Harmonic(0) should be 0.0 but it is %f", h2))
 	}

continuous/KraskovStoegbauerGrassberger.go

@@ -15,8 +15,8 @@ func KraskovStoegbauerGrassberger1(xy [][]float64, xIndices, yIndices []int, k i
 
 	r = 0.0
 
-	hk := harmonic(k)       // h(k)
-	hN := harmonic(len(xy)) // h(N)
+	hk := Harmonic(k)       // h(k)
+	hN := Harmonic(len(xy)) // h(N)
 
 	var bar *pb.ProgressBar
 
@@ -27,10 +27,10 @@ func KraskovStoegbauerGrassberger1(xy [][]float64, xIndices, yIndices []int, k i
 		epsilon := ksgGetEpsilon(k, xy[t], xy, xIndices, yIndices)
 
 		cNx := ksgCount(epsilon, xy[t], xy, xIndices) // N_x
-		hNx := harmonic(cNx + 1)                      // h(N_x)
+		hNx := Harmonic(cNx + 1)                      // h(N_x)
 
 		cNy := ksgCount(epsilon, xy[t], xy, yIndices) // N_y
-		hNy := harmonic(cNy + 1)                      // h(N_y)
+		hNy := Harmonic(cNy + 1)                      // h(N_y)
 
 		r -= hNx + hNy
 
@@ -58,8 +58,8 @@ func KraskovStoegbauerGrassberger2(xy [][]float64, xIndices, yIndices []int, k i
 
 	r = 0.0
 
-	hk := harmonic(k)
-	hN := harmonic(len(xy))
+	hk := Harmonic(k)
+	hN := Harmonic(len(xy))
 
 	var bar *pb.ProgressBar
 
@@ -70,10 +70,10 @@ func KraskovStoegbauerGrassberger2(xy [][]float64, xIndices, yIndices []int, k i
 		epsilon := ksgGetEpsilon(k, xy[t], xy, xIndices, yIndices)
 
 		cNx := ksgCount(epsilon, xy[t], xy, xIndices)
-		hNx := harmonic(cNx)
+		hNx := Harmonic(cNx)
 
 		cNy := ksgCount(epsilon, xy[t], xy, yIndices)
-		hNy := harmonic(cNy)
+		hNy := Harmonic(cNy)
 
 		r -= hNx + hNy
 
@@ -110,8 +110,8 @@ func ksgGetEpsilon(k int, xy []float64, data [][]float64, xIndices, yIndices []i
 }
 
 func ksgMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
-	xDistance := distance(a, b, xIndices)
-	yDistance := distance(a, b, yIndices)
+	xDistance := Distance(a, b, xIndices)
+	yDistance := Distance(a, b, yIndices)
 	return math.Max(xDistance, yDistance)
 }
 
@@ -120,7 +120,7 @@ func ksgMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
 func ksgCount(epsilon float64, xy []float64, data [][]float64, indices []int) (c int) {
 
 	for t := 0; t < len(data); t++ {
-		if distance(xy, data[t], indices) < epsilon {
+		if Distance(xy, data[t], indices) < epsilon {
 			c++
 		}
 	}

continuous/state/FrenzelPompe.go

@@ -4,6 +4,8 @@ import (
 	"math"
 	"sort"
 
+	"github.com/kzahedi/goent/continuous"
+
 	pb "gopkg.in/cheggaaa/pb.v1"
 )
 
@@ -15,7 +17,7 @@ func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, et
 
 	r := make([]float64, len(xyz), len(xyz))
 
-	hk := harmonic(k - 1)
+	hk := continuous.Harmonic(k - 1)
 
 	var bar *pb.ProgressBar
 
@@ -27,13 +29,13 @@ func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, et
 		epsilon := fpGetEpsilon(k, v, xyz, xIndices, yIndices, zIndices)
 
 		cNxz := fpCount2(epsilon, v, xyz, xIndices, zIndices)
-		hNxz := harmonic(cNxz)
+		hNxz := continuous.Harmonic(cNxz)
 
 		cNyz := fpCount2(epsilon, v, xyz, yIndices, zIndices)
-		hNyz := harmonic(cNyz)
+		hNyz := continuous.Harmonic(cNyz)
 
 		cNz := fpCount1(epsilon, v, xyz, zIndices)
-		hNz := harmonic(cNz)
+		hNz := continuous.Harmonic(cNz)
 
 		r[t] = hNxz + hNyz - hNz - hk
 
@@ -53,9 +55,9 @@ func FrenzelPompe(xyz [][]float64, xIndices, yIndices, zIndices []int, k int, et
 // fpMaxNorm3 computes the max-norm of two 3-dimensional vectors
 //   maxnorm(a,b) = max( |a[0] - b[0]|, |a[1] - b[1]|, |a[2] - b[2]|)
 func fpMaxNorm3(a, b []float64, xIndices, yIndices, zIndices []int) float64 {
-	xDist := distance(a, b, xIndices)
-	yDist := distance(a, b, yIndices)
-	zDist := distance(a, b, zIndices)
+	xDist := continuous.Distance(a, b, xIndices)
+	yDist := continuous.Distance(a, b, yIndices)
+	zDist := continuous.Distance(a, b, zIndices)
 	return math.Max(xDist, math.Max(yDist, zDist))
 }
 
@@ -92,8 +94,8 @@ func fpCount2(epsilon float64, xyz []float64, data [][]float64, xIndices, yIndic
 // fpMaxNorm2 computes the max-norm of two 3-dimensional vectors
 //   maxnorm(a,b) = max( |a[0] - b[0]|, |a[1] - b[1]|)
 func fpMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
-	xDist := distance(a, b, xIndices)
-	yDist := distance(a, b, yIndices)
+	xDist := continuous.Distance(a, b, xIndices)
+	yDist := continuous.Distance(a, b, yIndices)
 	return math.Max(xDist, yDist)
 }
 
@@ -102,7 +104,7 @@ func fpMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
 func fpCount1(epsilon float64, xyz []float64, data [][]float64, zIndices []int) (c int) {
 	c = -1 // because we will also count xyz[t] vs. xyz[t]
 	for t := 0; t < len(data); t++ {
-		if distance(xyz, data[t], zIndices) < epsilon {
+		if continuous.Distance(xyz, data[t], zIndices) < epsilon {
 			c++
 		}
 	}

continuous/state/KraskovStoegbauerGrassberger.go

@@ -4,6 +4,7 @@ import (
 	"math"
 	"sort"
 
+	"github.com/kzahedi/goent/continuous"
 	pb "gopkg.in/cheggaaa/pb.v1"
 )
 
@@ -16,8 +17,8 @@ func KraskovStoegbauerGrassberger1(xy [][]float64, xIndices, yIndices []int, k i
 	N := float64(len(xy))
 	r := make([]float64, len(xy), len(xy))
 
-	hk := harmonic(k)       // h(k)
-	hN := harmonic(len(xy)) // h(N)
+	hk := continuous.Harmonic(k)       // h(k)
+	hN := continuous.Harmonic(len(xy)) // h(N)
 
 	var bar *pb.ProgressBar
 
@@ -29,10 +30,10 @@ func KraskovStoegbauerGrassberger1(xy [][]float64, xIndices, yIndices []int, k i
 		epsilon := ksgGetEpsilon(k, xy[t], xy, xIndices, yIndices)
 
 		cNx := ksgCount(epsilon, xy[t], xy, xIndices) // N_x
-		hNx := harmonic(cNx + 1)                      // h(N_x)
+		hNx := continuous.Harmonic(cNx + 1)           // h(N_x)
 
 		cNy := ksgCount(epsilon, xy[t], xy, yIndices) // N_y
-		hNy := harmonic(cNy + 1)                      // h(N_y)
+		hNy := continuous.Harmonic(cNy + 1)           // h(N_y)
 
 		r[t] = (-hNx - hNy + hk + hN) / N
 
@@ -59,8 +60,8 @@ func KraskovStoegbauerGrassberger2(xy [][]float64, xIndices, yIndices []int, k i
 
 	N := float64(len(xy))
 
-	hk := harmonic(k)
-	hN := harmonic(n)
+	hk := continuous.Harmonic(k)
+	hN := continuous.Harmonic(n)
 	k1 := 1.0 / float64(k)
 
 	var bar *pb.ProgressBar
@@ -72,10 +73,10 @@ func KraskovStoegbauerGrassberger2(xy [][]float64, xIndices, yIndices []int, k i
 		epsilon := ksgGetEpsilon(k, xy[t], xy, xIndices, yIndices)
 
 		cNx := ksgCount(epsilon, xy[t], xy, xIndices)
-		hNx := harmonic(cNx)
+		hNx := continuous.Harmonic(cNx)
 
 		cNy := ksgCount(epsilon, xy[t], xy, yIndices)
-		hNy := harmonic(cNy)
+		hNy := continuous.Harmonic(cNy)
 
 		r[t] = (-hNx - hNy + hk + hN - k1) / N
 
@@ -108,8 +109,8 @@ func ksgGetEpsilon(k int, xy []float64, data [][]float64, xIndices, yIndices []i
 }
 
 func ksgMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
-	xDistance := distance(a, b, xIndices)
-	yDistance := distance(a, b, yIndices)
+	xDistance := continuous.Distance(a, b, xIndices)
+	yDistance := continuous.Distance(a, b, yIndices)
 	return math.Max(xDistance, yDistance)
 }
 
@@ -118,7 +119,7 @@ func ksgMaxNorm2(a, b []float64, xIndices, yIndices []int) float64 {
 func ksgCount(epsilon float64, xy []float64, data [][]float64, indices []int) (c int) {
 
 	for t := 0; t < len(data); t++ {
-		if distance(xy, data[t], indices) < epsilon {
+		if continuous.Distance(xy, data[t], indices) < epsilon {
 			c++
 		}
 	}