Use arbitrary precision floating-point numbers

bbox.go

@@ -1,48 +1,49 @@
 package polygol
 
-type bbox struct {
-	ll point
-	ur point
+type Bbox struct {
+	ll Vector
+	ur Vector
 }
 
-func (b bbox) isInBbox(point point) bool {
-	return b.ll.x <= point.x &&
-		point.x <= b.ur.x &&
-		b.ll.y <= point.y &&
-		point.y <= b.ur.y
+func (bbox Bbox) isInBbox(point Vector) bool {
+	return bbox.ll.x.isLessThanOrEqualTo(point.x) &&
+		point.x.isLessThanOrEqualTo(bbox.ur.x) &&
+		bbox.ll.y.isLessThanOrEqualTo(point.y) &&
+		point.y.isLessThanOrEqualTo(bbox.ur.y)
+
 }
 
-func (b bbox) getBboxOverlap(ob bbox) *bbox {
+func (b1 Bbox) getBboxOverlap(b2 Bbox) *Bbox {
 	// check if the bboxes overlap at all
-	if ob.ur.x < b.ll.x ||
-		b.ur.x < ob.ll.x ||
-		ob.ur.y < b.ll.y ||
-		b.ur.y < ob.ll.y {
+	if b2.ur.x.isLessThan(b1.ll.x) ||
+		b1.ur.x.isLessThan(b2.ll.x) ||
+		b2.ur.y.isLessThan(b1.ll.y) ||
+		b1.ur.y.isLessThan(b2.ll.y) {
 		return nil
 	}
 
 	// find the middle two X values
-	lowerX := b.ll.x
-	if b.ll.x < ob.ll.x {
-		lowerX = ob.ll.x
+	lowerX := b1.ll.x
+	if b1.ll.x.isLessThan(b2.ll.x) {
+		lowerX = b2.ll.x
 	}
-	upperX := ob.ur.x
-	if b.ur.x < ob.ur.x {
-		upperX = b.ur.x
+	upperX := b2.ur.x
+	if b1.ur.x.isLessThan(b2.ur.x) {
+		upperX = b1.ur.x
 	}
 
 	// find the middle two Y values
-	lowerY := b.ll.y
-	if b.ll.y < ob.ll.y {
-		lowerY = ob.ll.y
+	lowerY := b1.ll.y
+	if b1.ll.y.isLessThan(b2.ll.y) {
+		lowerY = b2.ll.y
 	}
-	upperY := ob.ur.y
-	if b.ur.y < ob.ur.y {
-		upperY = b.ur.y
+	upperY := b2.ur.y
+	if b1.ur.y.isLessThan(b2.ur.y) {
+		upperY = b1.ur.y
 	}
 
-	return &bbox{
-		ll: point{x: lowerX, y: lowerY},
-		ur: point{x: upperX, y: upperY},
+	return &Bbox{
+		ll: Vector{x: lowerX, y: lowerY},
+		ur: Vector{x: upperX, y: upperY},
 	}
 }

bignumber.go

@@ -0,0 +1,121 @@
+package polygol
+
+import (
+	"github.com/cockroachdb/apd/v3"
+)
+
+var context = apd.BaseContext.WithPrecision(50)
+
+func init() {
+	context.Rounding = apd.RoundHalfUp
+}
+
+type BigNumber struct {
+	f *apd.Decimal
+}
+
+func (b BigNumber) String() string {
+	return b.f.Text('f')
+}
+
+func newBigNumber(f float64) BigNumber {
+	b := BigNumber{f: new(apd.Decimal)}
+	b.f.SetFloat64(f)
+	return b
+}
+
+func bigZero() BigNumber {
+	b := BigNumber{f: new(apd.Decimal)}
+	// b.f.SetPrec(1024)
+	return b
+}
+
+func (b BigNumber) plus(other BigNumber) BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Add(r.f, b.f, other.f)
+	return r
+}
+
+func (b BigNumber) minus(other BigNumber) BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Sub(r.f, b.f, other.f)
+	return r
+}
+
+func (b BigNumber) times(other BigNumber) BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Mul(r.f, b.f, other.f)
+	return r
+}
+
+func (b BigNumber) div(other BigNumber) BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Quo(r.f, b.f, other.f)
+	return r
+
+}
+func (b BigNumber) abs() BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Abs(r.f, b.f)
+	return r
+}
+
+func (b BigNumber) Cmp(other BigNumber) int {
+	return b.f.Cmp(other.f)
+}
+
+func (b BigNumber) equalTo(other BigNumber) bool {
+	return b.f.Cmp(other.f) == 0
+}
+
+func (b BigNumber) isLessThan(other BigNumber) bool {
+	return b.f.Cmp(other.f) == -1
+}
+
+func (b BigNumber) isLessThanOrEqualTo(other BigNumber) bool {
+	return b.f.Cmp(other.f) != 1
+}
+
+func (b BigNumber) isGreaterThanOrEqualTo(other BigNumber) bool {
+	return b.f.Cmp(other.f) != -1
+}
+
+func (b BigNumber) isGreaterThan(other BigNumber) bool {
+	return b.f.Cmp(other.f) == 1
+}
+
+func (b BigNumber) notEqualTo(other BigNumber) bool {
+	return b.f.Cmp(other.f) != 0
+}
+
+func (b BigNumber) closeTo(other BigNumber) bool {
+	return b.minus(other).abs().isLessThanOrEqualTo(newBigNumber(2e-16))
+}
+
+func (b BigNumber) number() float64 {
+	f, _ := b.f.Float64()
+	return f
+}
+
+func (b BigNumber) sqrt() BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	context.Sqrt(r.f, b.f)
+	return r
+}
+
+func (b BigNumber) negated() BigNumber {
+	r := BigNumber{f: new(apd.Decimal)}
+	r.f.Neg(b.f)
+	return r
+}
+
+func (b BigNumber) isZero() bool {
+	return b.f.Cmp(new(apd.Decimal)) == 0
+}
+
+func bigInf(setToNegative bool) BigNumber {
+	if setToNegative {
+		return newBigNumber(-1e99)
+	}
+	return newBigNumber(1e99)
+}

end-to-end_test.go

@@ -4,13 +4,10 @@ import (
 	"fmt"
 	"io/ioutil"
 	"log"
-	"os"
 	"path"
 	"path/filepath"
 	"strings"
 	"testing"
-
-	geojson "github.com/engelsjk/polygol/geojson"
 )
 
 const (
@@ -34,14 +31,17 @@ type testCase struct {
 }
 
 func TestEndToEnd(t *testing.T) {
-
 	targets, err := ioutil.ReadDir(endToEndDir)
 	if err != nil {
 		t.Fatal(err)
 	}
 
 	for _, target := range targets {
 
+		if targetOnly != "" && target.Name() != targetOnly {
+			continue
+		}
+
 		if contains(targetsSkip, target.Name()) {
 			fmt.Printf("skipping target %s...\n", target.Name())
 			continue
@@ -113,7 +113,7 @@ func TestEndToEnd(t *testing.T) {
 
 			t.Run(testCase.Name, func(t *testing.T) {
 
-				t.Parallel() // run all end-to-end tests in parallel
+				// t.Parallel() // run all end-to-end tests in parallel
 
 				if contains(opsSkip, testCase.OperationType) {
 					fmt.Printf("skipping op type %s...\n", testCase.OperationType)
@@ -127,11 +127,17 @@ func TestEndToEnd(t *testing.T) {
 				expected := geoms[0]
 
 				result, err := newOperation(testCase.OperationType).run(args[0], args[1:]...)
+				resetPrecision()
 				if err != nil {
 					t.Error(err)
 				}
+				same := equalMultiPoly(expected, result)
+				if !same {
+					// d, _ := diff.Diff(expected, result)
+					// t.Fatal(d)
+					t.Fatal("resulting geometry does not match expectations")
 
-				expect(t, equalMultiPoly(expected, result))
+				}
 			})
 		}
 	}
@@ -145,50 +151,3 @@ func contains(s []string, str string) bool {
 	}
 	return false
 }
-
-func loadGeoms(filepath string) ([]Geom, error) {
-
-	fmt.Println(filepath)
-	f, err := os.Open(filepath)
-	if err != nil {
-		return nil, err
-	}
-	defer f.Close()
-
-	b, err := ioutil.ReadAll(f)
-	if err != nil {
-		return nil, err
-	}
-
-	newFeatures := unmarshalFeatureOrFeatureCollection(b)
-
-	geoms := make([]Geom, len(newFeatures))
-	for i := range newFeatures {
-		fg := newFeatures[i].Geometry
-		switch fg.Type {
-		case "Polygon":
-			geoms[i] = Geom{fg.Polygon}
-		case "MultiPolygon":
-			geoms[i] = fg.MultiPolygon
-		default:
-			return nil, fmt.Errorf("only polygon or multipolygon geometry types supported")
-		}
-	}
-
-	return geoms, nil
-}
-
-func unmarshalFeatureOrFeatureCollection(b []byte) []*geojson.Feature {
-	feature, err := geojson.UnmarshalFeature(b)
-	if err != nil {
-		return nil
-	}
-	if feature.Type != "FeatureCollection" {
-		return []*geojson.Feature{feature}
-	}
-	fc, err := geojson.UnmarshalFeatureCollection(b)
-	if err != nil {
-		return nil
-	}
-	return fc.Features
-}