Source file src/github.com/golang/geo/s2/edge_distances_test.go

Documentation: github.com/golang/geo/s2

     1  // Copyright 2017 Google Inc. All rights reserved.
     2  //
     3  // Licensed under the Apache License, Version 2.0 (the "License");
     4  // you may not use this file except in compliance with the License.
     5  // You may obtain a copy of the License at
     6  //
     7  //     http://www.apache.org/licenses/LICENSE-2.0
     8  //
     9  // Unless required by applicable law or agreed to in writing, software
    10  // distributed under the License is distributed on an "AS IS" BASIS,
    11  // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    12  // See the License for the specific language governing permissions and
    13  // limitations under the License.
    14  
    15  package s2
    16  
    17  import (
    18  	"math"
    19  	"testing"
    20  
    21  	"github.com/golang/geo/r3"
    22  	"github.com/golang/geo/s1"
    23  )
    24  
    25  func TestEdgeDistancesCheckDistance(t *testing.T) {
    26  	tests := []struct {
    27  		x, a, b r3.Vector
    28  		distRad float64
    29  		want    r3.Vector
    30  	}{
    31  		{
    32  			x:       r3.Vector{1, 0, 0},
    33  			a:       r3.Vector{1, 0, 0},
    34  			b:       r3.Vector{0, 1, 0},
    35  			distRad: 0,
    36  			want:    r3.Vector{1, 0, 0},
    37  		},
    38  		{
    39  			x:       r3.Vector{0, 1, 0},
    40  			a:       r3.Vector{1, 0, 0},
    41  			b:       r3.Vector{0, 1, 0},
    42  			distRad: 0,
    43  			want:    r3.Vector{0, 1, 0},
    44  		},
    45  		{
    46  			x:       r3.Vector{1, 3, 0},
    47  			a:       r3.Vector{1, 0, 0},
    48  			b:       r3.Vector{0, 1, 0},
    49  			distRad: 0,
    50  			want:    r3.Vector{1, 3, 0},
    51  		},
    52  		{
    53  			x:       r3.Vector{0, 0, 1},
    54  			a:       r3.Vector{1, 0, 0},
    55  			b:       r3.Vector{0, 1, 0},
    56  			distRad: math.Pi / 2,
    57  			want:    r3.Vector{1, 0, 0},
    58  		},
    59  		{
    60  			x:       r3.Vector{0, 0, -1},
    61  			a:       r3.Vector{1, 0, 0},
    62  			b:       r3.Vector{0, 1, 0},
    63  			distRad: math.Pi / 2,
    64  			want:    r3.Vector{1, 0, 0},
    65  		},
    66  		{
    67  			x:       r3.Vector{-1, -1, 0},
    68  			a:       r3.Vector{1, 0, 0},
    69  			b:       r3.Vector{0, 1, 0},
    70  			distRad: 0.75 * math.Pi,
    71  			want:    r3.Vector{1, 0, 0},
    72  		},
    73  		{
    74  			x:       r3.Vector{0, 1, 0},
    75  			a:       r3.Vector{1, 0, 0},
    76  			b:       r3.Vector{1, 1, 0},
    77  			distRad: math.Pi / 4,
    78  			want:    r3.Vector{1, 1, 0},
    79  		},
    80  		{
    81  			x:       r3.Vector{0, -1, 0},
    82  			a:       r3.Vector{1, 0, 0},
    83  			b:       r3.Vector{1, 1, 0},
    84  			distRad: math.Pi / 2,
    85  			want:    r3.Vector{1, 0, 0},
    86  		},
    87  		{
    88  			x:       r3.Vector{0, -1, 0},
    89  			a:       r3.Vector{1, 0, 0},
    90  			b:       r3.Vector{-1, 1, 0},
    91  			distRad: math.Pi / 2,
    92  			want:    r3.Vector{1, 0, 0},
    93  		},
    94  		{
    95  			x:       r3.Vector{-1, -1, 0},
    96  			a:       r3.Vector{1, 0, 0},
    97  			b:       r3.Vector{-1, 1, 0},
    98  			distRad: math.Pi / 2,
    99  			want:    r3.Vector{-1, 1, 0},
   100  		},
   101  		{
   102  			x:       r3.Vector{1, 1, 1},
   103  			a:       r3.Vector{1, 0, 0},
   104  			b:       r3.Vector{0, 1, 0},
   105  			distRad: math.Asin(math.Sqrt(1.0 / 3.0)),
   106  			want:    r3.Vector{1, 1, 0},
   107  		},
   108  		{
   109  			x:       r3.Vector{1, 1, -1},
   110  			a:       r3.Vector{1, 0, 0},
   111  			b:       r3.Vector{0, 1, 0},
   112  			distRad: math.Asin(math.Sqrt(1.0 / 3.0)),
   113  			want:    r3.Vector{1, 1, 0}},
   114  		{
   115  			x:       r3.Vector{-1, 0, 0},
   116  			a:       r3.Vector{1, 1, 0},
   117  			b:       r3.Vector{1, 1, 0},
   118  			distRad: 0.75 * math.Pi,
   119  			want:    r3.Vector{1, 1, 0},
   120  		},
   121  		{
   122  			x:       r3.Vector{0, 0, -1},
   123  			a:       r3.Vector{1, 1, 0},
   124  			b:       r3.Vector{1, 1, 0},
   125  			distRad: math.Pi / 2,
   126  			want:    r3.Vector{1, 1, 0},
   127  		},
   128  		{
   129  			x:       r3.Vector{-1, 0, 0},
   130  			a:       r3.Vector{1, 0, 0},
   131  			b:       r3.Vector{1, 0, 0},
   132  			distRad: math.Pi,
   133  			want:    r3.Vector{1, 0, 0},
   134  		},
   135  	}
   136  
   137  	for _, test := range tests {
   138  		x := Point{test.x.Normalize()}
   139  		a := Point{test.a.Normalize()}
   140  		b := Point{test.b.Normalize()}
   141  		want := Point{test.want.Normalize()}
   142  
   143  		if d := DistanceFromSegment(x, a, b).Radians(); !float64Near(d, test.distRad, 1e-15) {
   144  			t.Errorf("DistanceFromSegment(%v, %v, %v) = %v, want %v", x, a, b, d, test.distRad)
   145  		}
   146  
   147  		closest := Project(x, a, b)
   148  		if !closest.ApproxEqual(want) {
   149  			t.Errorf("ClosestPoint(%v, %v, %v) = %v, want %v", x, a, b, closest, want)
   150  		}
   151  
   152  		if minDistance, ok := UpdateMinDistance(x, a, b, 0); ok {
   153  			t.Errorf("UpdateMinDistance(%v, %v, %v, %v) = %v, want %v", x, a, b, 0, minDistance, 0)
   154  		}
   155  
   156  		minDistance, ok := UpdateMinDistance(x, a, b, s1.InfChordAngle())
   157  		if !ok {
   158  			t.Errorf("UpdateMinDistance(%v, %v, %v, %v) = %v, want %v", x, a, b, s1.InfChordAngle(), minDistance, s1.InfChordAngle())
   159  		}
   160  
   161  		if !float64Near(test.distRad, minDistance.Angle().Radians(), 1e-15) {
   162  			t.Errorf("MinDistance between %v and %v,%v = %v, want %v within %v", x, a, b, minDistance.Angle().Radians(), test.distRad, 1e-15)
   163  		}
   164  	}
   165  }
   166  
   167  func TestEdgeDistancesUpdateMinInteriorDistanceLowerBoundOptimizationIsConservative(t *testing.T) {
   168  	// Verifies that alwaysUpdateMinInteriorDistance computes the lower bound
   169  	// on the true distance conservatively.  (This test used to fail.)
   170  	x := PointFromCoords(-0.017952729194524016, -0.30232422079175203, 0.95303607751077712)
   171  	a := PointFromCoords(-0.017894725505830295, -0.30229974986194175, 0.95304493075220664)
   172  	b := PointFromCoords(-0.017986591360900289, -0.30233851195954353, 0.95303090543659963)
   173  
   174  	minDistance, ok := UpdateMinDistance(x, a, b, s1.InfChordAngle())
   175  	if !ok {
   176  		t.Errorf("UpdateMinDistance(%v, %v, %v, %v) = %v, want %v", x, a, b, s1.InfChordAngle(), minDistance, s1.InfChordAngle())
   177  	}
   178  	minDistance = minDistance.Successor()
   179  	minDistance, ok = UpdateMinDistance(x, a, b, minDistance)
   180  	if !ok {
   181  		t.Errorf("UpdateMinDistance(%v, %v, %v, %v) = %v, want %v", x, a, b, s1.InfChordAngle(), minDistance, minDistance)
   182  	}
   183  
   184  }
   185  
   186  func TestEdgeDistancesUpdateMinInteriorDistanceRejectionTestIsConservative(t *testing.T) {
   187  	// This test checks several representative cases where previously
   188  	// UpdateMinInteriorDistance was failing to update the distance because a
   189  	// rejection test was not being done conservatively.
   190  	//
   191  	// Note that all of the edges AB in this test are nearly antipodal.
   192  	minDist := s1.ChordAngleFromSquaredLength(6.3897233584120815e-26)
   193  
   194  	tests := []struct {
   195  		x, a, b Point
   196  		minDist s1.ChordAngle
   197  		want    bool
   198  	}{
   199  		{
   200  
   201  			x:       Point{r3.Vector{1, -4.6547732744037044e-11, -5.6374428459823598e-89}},
   202  			a:       Point{r3.Vector{1, -8.9031850507928352e-11, 0}},
   203  			b:       Point{r3.Vector{-0.99999999999996347, 2.7030110029169596e-07, 1.555092348806121e-99}},
   204  			minDist: minDist,
   205  			want:    false,
   206  		},
   207  		{
   208  			x:       Point{r3.Vector{1, -4.7617930898495072e-13, 0}},
   209  			a:       Point{r3.Vector{-1, -1.6065916409055676e-10, 0}},
   210  			b:       Point{r3.Vector{1, 0, 9.9964883247706732e-35}},
   211  			minDist: minDist,
   212  			want:    false,
   213  		},
   214  		{
   215  			x:       Point{r3.Vector{1, 0, 0}},
   216  			a:       Point{r3.Vector{1, -8.4965026896454536e-11, 0}},
   217  			b:       Point{r3.Vector{-0.99999999999966138, 8.2297529603339328e-07, 9.6070344113320997e-21}},
   218  			minDist: minDist,
   219  			want:    false,
   220  		},
   221  	}
   222  
   223  	for _, test := range tests {
   224  		if _, ok := UpdateMinDistance(test.x, test.a, test.b, test.minDist); !ok {
   225  			t.Errorf("UpdateMinDistance(%v, %v, %v, %v) = %v, want %v", test.x, test.a, test.b, test.minDist, ok, test.want)
   226  		}
   227  	}
   228  }
   229  
   230  func TestEdgeDistancesCheckMaxDistance(t *testing.T) {
   231  	tests := []struct {
   232  		x, a, b r3.Vector
   233  		distRad float64
   234  	}{
   235  		{
   236  			x:       r3.Vector{1, 0, 1},
   237  			a:       r3.Vector{1, 0, 0},
   238  			b:       r3.Vector{0, 1, 0},
   239  			distRad: math.Pi / 2,
   240  		},
   241  		{
   242  			x:       r3.Vector{1, 0, -1},
   243  			a:       r3.Vector{1, 0, 0},
   244  			b:       r3.Vector{0, 1, 0},
   245  			distRad: math.Pi / 2,
   246  		},
   247  		{
   248  			x:       r3.Vector{0, 1, 1},
   249  			a:       r3.Vector{1, 0, 0},
   250  			b:       r3.Vector{0, 1, 0},
   251  			distRad: math.Pi / 2,
   252  		},
   253  		{
   254  			x:       r3.Vector{0, 1, -1},
   255  			a:       r3.Vector{1, 0, 0},
   256  			b:       r3.Vector{0, 1, 0},
   257  			distRad: math.Pi / 2,
   258  		},
   259  		{
   260  			x:       r3.Vector{1, 1, 1},
   261  			a:       r3.Vector{1, 0, 0},
   262  			b:       r3.Vector{0, 1, 0},
   263  			distRad: math.Asin(math.Sqrt(2. / 3)),
   264  		},
   265  		{
   266  			x:       r3.Vector{1, 1, -1},
   267  			a:       r3.Vector{1, 0, 0},
   268  			b:       r3.Vector{0, 1, 0},
   269  			distRad: math.Asin(math.Sqrt(2. / 3)),
   270  		},
   271  		{
   272  			x:       r3.Vector{1, 0, 0},
   273  			a:       r3.Vector{1, 1, 0},
   274  			b:       r3.Vector{1, -1, 0},
   275  			distRad: math.Pi / 4,
   276  		},
   277  		{
   278  			x:       r3.Vector{0, 1, 0},
   279  			a:       r3.Vector{1, 1, 0},
   280  			b:       r3.Vector{1, 1, 0},
   281  			distRad: math.Pi / 4,
   282  		},
   283  		{
   284  			x:       r3.Vector{0, 0, 1},
   285  			a:       r3.Vector{0, 1, 1},
   286  			b:       r3.Vector{0, -1, 1},
   287  			distRad: math.Pi / 4,
   288  		},
   289  		{
   290  			x:       r3.Vector{0, 0, 1},
   291  			a:       r3.Vector{1, 0, 0},
   292  			b:       r3.Vector{1, 0, -1},
   293  			distRad: 3 * math.Pi / 4,
   294  		},
   295  		{
   296  			x:       r3.Vector{0, 0, 1},
   297  			a:       r3.Vector{1, 0, 0},
   298  			b:       r3.Vector{1, 1, -math.Sqrt2},
   299  			distRad: 3 * math.Pi / 4,
   300  		},
   301  		{
   302  			x:       r3.Vector{0, 0, 1},
   303  			a:       r3.Vector{0, 0, -1},
   304  			b:       r3.Vector{0, 0, -1},
   305  			distRad: math.Pi,
   306  		},
   307  	}
   308  
   309  	for _, test := range tests {
   310  		x := Point{test.x.Normalize()}
   311  		a := Point{test.a.Normalize()}
   312  		b := Point{test.b.Normalize()}
   313  
   314  		var ok bool
   315  		maxDistance := s1.StraightChordAngle
   316  		if maxDistance, ok = UpdateMaxDistance(x, a, b, maxDistance); ok {
   317  			t.Errorf("UpdateMaxDistance(%v, %v, %v, %v) = %v, want %v", x, a, b, s1.StraightChordAngle, maxDistance, s1.StraightChordAngle)
   318  		}
   319  
   320  		maxDistance = s1.NegativeChordAngle
   321  		if maxDistance, ok = UpdateMaxDistance(x, a, b, maxDistance); !ok {
   322  			t.Errorf("UpdateMaxDistance(%v, %v, %v, %v) = %v, want > %v", x, a, b, s1.NegativeChordAngle, maxDistance, s1.NegativeChordAngle)
   323  		}
   324  
   325  		if !float64Near(test.distRad, maxDistance.Angle().Radians(), 1e-15) {
   326  			t.Errorf("MaxDistance between %v and %v, %v = %v, want %v within %v", x, a, b, maxDistance.Angle().Radians(), test.distRad, 1e-15)
   327  		}
   328  	}
   329  }
   330  
   331  func TestEdgeDistancesInterpolate(t *testing.T) {
   332  	// Choose test points designed to expose floating-point errors.
   333  	p1 := PointFromCoords(0.1, 1e-30, 0.3)
   334  	p2 := PointFromCoords(-0.7, -0.55, -1e30)
   335  	i := PointFromCoords(1, 0, 0)
   336  	j := PointFromCoords(0, 1, 0)
   337  
   338  	// Take a small fraction along the curve, 1/1000 of the way.
   339  	p := Interpolate(0.001, i, j)
   340  
   341  	tests := []struct {
   342  		a, b Point
   343  		dist float64
   344  		want Point
   345  	}{
   346  		// A zero-length edge.
   347  		{p1, p1, 0, p1},
   348  		{p1, p1, 1, p1},
   349  		// Start, end, and middle of a medium-length edge.
   350  		{p1, p2, 0, p1},
   351  		{p1, p2, 1, p2},
   352  		{p1, p2, 0.5, Point{(p1.Add(p2.Vector)).Mul(0.5)}},
   353  
   354  		// Test that interpolation is done using distances on the sphere
   355  		// rather than linear distances.
   356  		{
   357  			Point{r3.Vector{1, 0, 0}},
   358  			Point{r3.Vector{0, 1, 0}},
   359  			1.0 / 3.0,
   360  			Point{r3.Vector{math.Sqrt(3), 1, 0}},
   361  		},
   362  		{
   363  			Point{r3.Vector{1, 0, 0}},
   364  			Point{r3.Vector{0, 1, 0}},
   365  			2.0 / 3.0,
   366  			Point{r3.Vector{1, math.Sqrt(3), 0}},
   367  		},
   368  
   369  		// InterpolateCanExtrapolate checks
   370  
   371  		// Initial vectors at 90 degrees.
   372  		{i, j, 0, Point{r3.Vector{1, 0, 0}}},
   373  		{i, j, 1, Point{r3.Vector{0, 1, 0}}},
   374  		{i, j, 1.5, Point{r3.Vector{-1, 1, 0}}},
   375  		{i, j, 2, Point{r3.Vector{-1, 0, 0}}},
   376  		{i, j, 3, Point{r3.Vector{0, -1, 0}}},
   377  		{i, j, 4, Point{r3.Vector{1, 0, 0}}},
   378  
   379  		// Negative values of t.
   380  		{i, j, -1, Point{r3.Vector{0, -1, 0}}},
   381  		{i, j, -2, Point{r3.Vector{-1, 0, 0}}},
   382  		{i, j, -3, Point{r3.Vector{0, 1, 0}}},
   383  		{i, j, -4, Point{r3.Vector{1, 0, 0}}},
   384  
   385  		// Initial vectors at 45 degrees.
   386  		{i, Point{r3.Vector{1, 1, 0}}, 2, Point{r3.Vector{0, 1, 0}}},
   387  		{i, Point{r3.Vector{1, 1, 0}}, 3, Point{r3.Vector{-1, 1, 0}}},
   388  		{i, Point{r3.Vector{1, 1, 0}}, 4, Point{r3.Vector{-1, 0, 0}}},
   389  
   390  		// Initial vectors at 135 degrees.
   391  		{i, Point{r3.Vector{-1, 1, 0}}, 2, Point{r3.Vector{0, -1, 0}}},
   392  
   393  		// Test that we should get back where we started by interpolating
   394  		// the 1/1000th by 1000.
   395  		{i, p, 1000, j},
   396  	}
   397  
   398  	for _, test := range tests {
   399  		test.a = Point{test.a.Normalize()}
   400  		test.b = Point{test.b.Normalize()}
   401  		test.want = Point{test.want.Normalize()}
   402  
   403  		// We allow a bit more than the usual 1e-15 error tolerance because
   404  		// Interpolate() uses trig functions.
   405  		if got := Interpolate(test.dist, test.a, test.b); !pointsApproxEqual(got, test.want, 3e-15) {
   406  			t.Errorf("Interpolate(%v, %v, %v) = %v, want %v", test.dist, test.a, test.b, got, test.want)
   407  		}
   408  	}
   409  }
   410  
   411  func TestEdgeDistancesInterpolateOverLongEdge(t *testing.T) {
   412  	lng := math.Pi - 1e-2
   413  	a := Point{PointFromLatLng(LatLng{0, 0}).Normalize()}
   414  	b := Point{PointFromLatLng(LatLng{0, s1.Angle(lng)}).Normalize()}
   415  
   416  	for f := 0.4; f > 1e-15; f *= 0.1 {
   417  		// Test that interpolation is accurate on a long edge (but not so long that
   418  		// the definition of the edge itself becomes too unstable).
   419  		want := Point{PointFromLatLng(LatLng{0, s1.Angle(f * lng)}).Normalize()}
   420  		if got := Interpolate(f, a, b); !pointsApproxEqual(got, want, 3e-15) {
   421  			t.Errorf("long edge Interpolate(%v, %v, %v) = %v, want %v", f, a, b, got, want)
   422  		}
   423  
   424  		// Test the remainder of the dist also matches.
   425  		wantRem := Point{PointFromLatLng(LatLng{0, s1.Angle((1 - f) * lng)}).Normalize()}
   426  		if got := Interpolate(1-f, a, b); !pointsApproxEqual(got, wantRem, 3e-15) {
   427  			t.Errorf("long edge Interpolate(%v, %v, %v) = %v, want %v", 1-f, a, b, got, wantRem)
   428  		}
   429  	}
   430  }
   431  
   432  func TestEdgeDistancesInterpolateAntipodal(t *testing.T) {
   433  	p1 := PointFromCoords(0.1, 1e-30, 0.3)
   434  
   435  	// Test that interpolation on a 180 degree edge (antipodal endpoints) yields
   436  	// a result with the correct distance from each endpoint.
   437  	for dist := 0.0; dist <= 1.0; dist += 0.125 {
   438  		actual := Interpolate(dist, p1, Point{p1.Mul(-1)})
   439  		if !float64Near(actual.Distance(p1).Radians(), dist*math.Pi, 3e-15) {
   440  			t.Errorf("antipodal points Interpolate(%v, %v, %v) = %v, want %v", dist, p1, Point{p1.Mul(-1)}, actual, dist*math.Pi)
   441  		}
   442  	}
   443  }
   444  
   445  func TestEdgeDistancesRepeatedInterpolation(t *testing.T) {
   446  	// Check that points do not drift away from unit length when repeated
   447  	// interpolations are done.
   448  	for i := 0; i < 100; i++ {
   449  		a := randomPoint()
   450  		b := randomPoint()
   451  		for j := 0; j < 1000; j++ {
   452  			a = Interpolate(0.01, a, b)
   453  		}
   454  		if !a.Vector.IsUnit() {
   455  			t.Errorf("repeated Interpolate(%v, %v, %v) calls did not stay unit length for", 0.01, a, b)
   456  		}
   457  	}
   458  }
   459  
   460  func TestEdgeDistanceMinUpdateDistanceMaxError(t *testing.T) {
   461  	tests := []struct {
   462  		actual s1.Angle
   463  		maxErr s1.Angle
   464  	}{
   465  		{0, 1.5e-15},
   466  		{1e-8, 1e-15},
   467  		{1e-5, 1e-15},
   468  		{0.05, 1e-15},
   469  		{math.Pi/2 - 1e-8, 2e-15},
   470  		{math.Pi / 2, 2e-15},
   471  		{math.Pi/2 + 1e-8, 2e-15},
   472  		{math.Pi - 1e-5, 2e-10},
   473  		{math.Pi, 0},
   474  	}
   475  
   476  	// This checks that the error returned by minUpdateDistanceMaxError for
   477  	// the distance actual (measured in radians) corresponds to a distance error
   478  	// of less than maxErr (measured in radians).
   479  	//
   480  	// The reason for the awkward phraseology above is that the value returned by
   481  	// minUpdateDistanceMaxError is not a distance; it represents an error in
   482  	// the *squared* distance.
   483  	for _, test := range tests {
   484  		ca := s1.ChordAngleFromAngle(test.actual)
   485  		bound := ca.Expanded(minUpdateDistanceMaxError(ca)).Angle()
   486  
   487  		if got := s1.Angle(bound.Radians()) - test.actual; got > test.maxErr {
   488  			t.Errorf("minUpdateDistanceMaxError(%v)-%v = %v> %v, want <=", ca, got, test.actual, test.maxErr)
   489  		}
   490  	}
   491  }
   492  
   493  // TODO(roberts): TestEdgeDistanceUpdateMinInteriorDistanceMaxError once s2predicates
   494  // CompareEdgeDistance
   495  
   496  func TestEdgeDistancesEdgePairMinDistance(t *testing.T) {
   497  	var zero Point
   498  	tests := []struct {
   499  		a0, a1   Point
   500  		b0, b1   Point
   501  		distRads float64
   502  		wantA    Point
   503  		wantB    Point
   504  	}{
   505  		{
   506  			// One edge is degenerate.
   507  			a0:       PointFromCoords(1, 0, 1),
   508  			a1:       PointFromCoords(1, 0, 1),
   509  			b0:       PointFromCoords(1, -1, 0),
   510  			b1:       PointFromCoords(1, 1, 0),
   511  			distRads: math.Pi / 4,
   512  			wantA:    PointFromCoords(1, 0, 1),
   513  			wantB:    PointFromCoords(1, 0, 0),
   514  		},
   515  		{
   516  			// One edge is degenerate.
   517  			a0:       PointFromCoords(1, -1, 0),
   518  			a1:       PointFromCoords(1, 1, 0),
   519  			b0:       PointFromCoords(1, 0, 1),
   520  			b1:       PointFromCoords(1, 0, 1),
   521  			distRads: math.Pi / 4,
   522  			wantA:    PointFromCoords(1, 0, 0),
   523  			wantB:    PointFromCoords(1, 0, 1),
   524  		},
   525  		{
   526  			// Both edges are degenerate.
   527  			a0:       PointFromCoords(1, 0, 0),
   528  			a1:       PointFromCoords(1, 0, 0),
   529  			b0:       PointFromCoords(0, 1, 0),
   530  			b1:       PointFromCoords(0, 1, 0),
   531  			distRads: math.Pi / 2,
   532  			wantA:    PointFromCoords(1, 0, 0),
   533  			wantB:    PointFromCoords(0, 1, 0),
   534  		},
   535  		{
   536  			// Both edges are degenerate and antipodal.
   537  			a0:       PointFromCoords(1, 0, 0),
   538  			a1:       PointFromCoords(1, 0, 0),
   539  			b0:       PointFromCoords(-1, 0, 0),
   540  			b1:       PointFromCoords(-1, 0, 0),
   541  			distRads: math.Pi,
   542  			wantA:    PointFromCoords(1, 0, 0),
   543  			wantB:    PointFromCoords(-1, 0, 0),
   544  		},
   545  		{
   546  			// Two identical edges.
   547  			a0:       PointFromCoords(1, 0, 0),
   548  			a1:       PointFromCoords(0, 1, 0),
   549  			b0:       PointFromCoords(1, 0, 0),
   550  			b1:       PointFromCoords(0, 1, 0),
   551  			distRads: 0,
   552  			wantA:    zero,
   553  			wantB:    zero,
   554  		},
   555  		{
   556  			// Both edges are degenerate and identical.
   557  			a0:       PointFromCoords(1, 0, 0),
   558  			a1:       PointFromCoords(1, 0, 0),
   559  			b0:       PointFromCoords(1, 0, 0),
   560  			b1:       PointFromCoords(1, 0, 0),
   561  			distRads: 0,
   562  			wantA:    PointFromCoords(1, 0, 0),
   563  			wantB:    PointFromCoords(1, 0, 0),
   564  		},
   565  		// Edges that share exactly one vertex (all 4 possibilities).
   566  		{
   567  			a0:       PointFromCoords(1, 0, 0),
   568  			a1:       PointFromCoords(0, 1, 0),
   569  			b0:       PointFromCoords(0, 1, 0),
   570  			b1:       PointFromCoords(0, 1, 1),
   571  			distRads: 0,
   572  			wantA:    PointFromCoords(0, 1, 0),
   573  			wantB:    PointFromCoords(0, 1, 0),
   574  		},
   575  		{
   576  			a0:       PointFromCoords(0, 1, 0),
   577  			a1:       PointFromCoords(1, 0, 0),
   578  			b0:       PointFromCoords(0, 1, 0),
   579  			b1:       PointFromCoords(0, 1, 1),
   580  			distRads: 0,
   581  			wantA:    PointFromCoords(0, 1, 0),
   582  			wantB:    PointFromCoords(0, 1, 0),
   583  		},
   584  		{
   585  			a0:       PointFromCoords(1, 0, 0),
   586  			a1:       PointFromCoords(0, 1, 0),
   587  			b0:       PointFromCoords(0, 1, 1),
   588  			b1:       PointFromCoords(0, 1, 0),
   589  			distRads: 0,
   590  			wantA:    PointFromCoords(0, 1, 0),
   591  			wantB:    PointFromCoords(0, 1, 0),
   592  		},
   593  		{
   594  			a0:       PointFromCoords(0, 1, 0),
   595  			a1:       PointFromCoords(1, 0, 0),
   596  			b0:       PointFromCoords(0, 1, 1),
   597  			b1:       PointFromCoords(0, 1, 0),
   598  			distRads: 0,
   599  			wantA:    PointFromCoords(0, 1, 0),
   600  			wantB:    PointFromCoords(0, 1, 0),
   601  		},
   602  		{
   603  			// Two edges whose interiors cross.
   604  			a0:       PointFromCoords(1, -1, 0),
   605  			a1:       PointFromCoords(1, 1, 0),
   606  			b0:       PointFromCoords(1, 0, -1),
   607  			b1:       PointFromCoords(1, 0, 1),
   608  			distRads: 0,
   609  			wantA:    PointFromCoords(1, 0, 0),
   610  			wantB:    PointFromCoords(1, 0, 0),
   611  		},
   612  		// The closest distance occurs between two edge endpoints, but more than one
   613  		// endpoint pair is equally distant.
   614  		{
   615  			a0:       PointFromCoords(1, -1, 0),
   616  			a1:       PointFromCoords(1, 1, 0),
   617  			b0:       PointFromCoords(-1, 0, 0),
   618  			b1:       PointFromCoords(-1, 0, 1),
   619  			distRads: math.Acos(-0.5),
   620  			wantA:    zero,
   621  			wantB:    PointFromCoords(-1, 0, 1),
   622  		},
   623  		{
   624  			a0:       PointFromCoords(-1, 0, 0),
   625  			a1:       PointFromCoords(-1, 0, 1),
   626  			b0:       PointFromCoords(1, -1, 0),
   627  			b1:       PointFromCoords(1, 1, 0),
   628  			distRads: math.Acos(-0.5),
   629  			wantA:    PointFromCoords(-1, 0, 1),
   630  			wantB:    zero,
   631  		},
   632  		{
   633  			a0:       PointFromCoords(1, -1, 0),
   634  			a1:       PointFromCoords(1, 1, 0),
   635  			b0:       PointFromCoords(-1, 0, -1),
   636  			b1:       PointFromCoords(-1, 0, 1),
   637  			distRads: math.Acos(-0.5),
   638  			wantA:    zero,
   639  			wantB:    zero,
   640  		},
   641  	}
   642  
   643  	// Given two edges a0a1 and b0b1, check that the minimum distance
   644  	// between them is distRads, and that EdgePairClosestPoints returns
   645  	// wantA and wantB as the points that achieve this distance.
   646  	// Point{0, 0, 0} may be passed for wantA or wantB to indicate
   647  	// that both endpoints of the corresponding edge are equally distant,
   648  	// and therefore either one might be returned.
   649  	for _, test := range tests {
   650  		actualA, actualB := EdgePairClosestPoints(test.a0, test.a1, test.b0, test.b1)
   651  		if test.wantA == zero {
   652  			// either end point works.
   653  			if !(actualA == test.a0 || actualA == test.a1) {
   654  				t.Errorf("EdgePairClosestPoints(%v, %v, %v, %v) = %v, want %v or %v", test.a0, test.a1, test.b0, test.b1, actualA, test.a0, test.a1)
   655  			}
   656  		} else {
   657  			if !actualA.ApproxEqual(test.wantA) {
   658  				t.Errorf("EdgePairClosestPoints(%v, %v, %v, %v) = %v, want %v", test.a0, test.a1, test.b0, test.b1, actualA, test.wantA)
   659  			}
   660  		}
   661  
   662  		if test.wantB == zero {
   663  			// either end point works.
   664  			if !(actualB == test.b0 || actualB == test.b1) {
   665  				t.Errorf("EdgePairClosestPoints(%v, %v, %v, %v) = %v, want %v or %v", test.a0, test.a1, test.b0, test.b1, actualB, test.b0, test.b1)
   666  			}
   667  		} else {
   668  			if !actualB.ApproxEqual(test.wantB) {
   669  				t.Errorf("EdgePairClosestPoints(%v, %v, %v, %v) = %v, want %v", test.a0, test.a1, test.b0, test.b1, actualB, test.wantB)
   670  			}
   671  		}
   672  
   673  		var minDist s1.ChordAngle
   674  		var ok bool
   675  		minDist, ok = updateEdgePairMinDistance(test.a0, test.a1, test.b0, test.b1, minDist)
   676  		if ok {
   677  			t.Errorf("updateEdgePairMinDistance(%v, %v, %v, %v, %v) = %v, want updated to be false", test.a0, test.a1, test.b0, test.b1, 0, minDist)
   678  		}
   679  
   680  		minDist = s1.InfChordAngle()
   681  		minDist, ok = updateEdgePairMinDistance(test.a0, test.a1, test.b0, test.b1, minDist)
   682  		if !ok {
   683  			t.Errorf("updateEdgePairMinDistance(%v, %v, %v, %v, %v) = %v, want updated to be true", test.a0, test.a1, test.b0, test.b1, s1.InfChordAngle(), minDist)
   684  		}
   685  
   686  		if !float64Near(test.distRads, minDist.Angle().Radians(), epsilon) {
   687  			t.Errorf("minDist %v - %v = %v, want < %v", test.distRads, minDist.Angle().Radians(), (test.distRads - minDist.Angle().Radians()), epsilon)
   688  		}
   689  	}
   690  }
   691  
   692  func TestEdgeDistancesEdgePairMaxDistance(t *testing.T) {
   693  	tests := []struct {
   694  		a0, a1   Point
   695  		b0, b1   Point
   696  		distRads float64
   697  	}{
   698  		{
   699  			// Standard situation. Same hemisphere, not degenerate.
   700  			a0:       PointFromCoords(1, 0, 0),
   701  			a1:       PointFromCoords(0, 1, 0),
   702  			b0:       PointFromCoords(1, 1, 0),
   703  			b1:       PointFromCoords(1, 1, 1),
   704  			distRads: math.Acos(1 / math.Sqrt(3)),
   705  		},
   706  		{
   707  
   708  			// One edge is degenerate.
   709  			a0:       PointFromCoords(1, 0, 1),
   710  			a1:       PointFromCoords(1, 0, 1),
   711  			b0:       PointFromCoords(1, -1, 0),
   712  			b1:       PointFromCoords(1, 1, 0),
   713  			distRads: math.Acos(0.5),
   714  		},
   715  		{
   716  			a0:       PointFromCoords(1, -1, 0),
   717  			a1:       PointFromCoords(1, 1, 0),
   718  			b0:       PointFromCoords(1, 0, 1),
   719  			b1:       PointFromCoords(1, 0, 1),
   720  			distRads: math.Acos(0.5),
   721  		},
   722  		{
   723  			// Both edges are degenerate.
   724  			a0:       PointFromCoords(1, 0, 0),
   725  			a1:       PointFromCoords(1, 0, 0),
   726  			b0:       PointFromCoords(0, 1, 0),
   727  			b1:       PointFromCoords(0, 1, 0),
   728  			distRads: math.Pi / 2,
   729  		},
   730  		{
   731  			// Both edges are degenerate and antipodal.
   732  			a0:       PointFromCoords(1, 0, 0),
   733  			a1:       PointFromCoords(1, 0, 0),
   734  			b0:       PointFromCoords(-1, 0, 0),
   735  			b1:       PointFromCoords(-1, 0, 0),
   736  			distRads: math.Pi,
   737  		},
   738  		{
   739  			// Two identical edges.
   740  			a0:       PointFromCoords(1, 0, 0),
   741  			a1:       PointFromCoords(0, 1, 0),
   742  			b0:       PointFromCoords(1, 0, 0),
   743  			b1:       PointFromCoords(0, 1, 0),
   744  			distRads: math.Pi / 2,
   745  		},
   746  		{
   747  			// Both edges are degenerate and identical.
   748  			a0:       PointFromCoords(1, 0, 0),
   749  			a1:       PointFromCoords(1, 0, 0),
   750  			b0:       PointFromCoords(1, 0, 0),
   751  			b1:       PointFromCoords(1, 0, 0),
   752  			distRads: 0,
   753  		},
   754  		{
   755  			// Antipodal reflection of one edge crosses the other edge.
   756  			a0:       PointFromCoords(1, 0, 1),
   757  			a1:       PointFromCoords(1, 0, -1),
   758  			b0:       PointFromCoords(-1, -1, 0),
   759  			b1:       PointFromCoords(-1, 1, 0),
   760  			distRads: math.Pi,
   761  		},
   762  		{
   763  			// One vertex of one edge touches the interior of the antipodal reflection
   764  			// of the other edge.
   765  			a0:       PointFromCoords(1, 0, 1),
   766  			a1:       PointFromCoords(1, 0, 0),
   767  			b0:       PointFromCoords(-1, -1, 0),
   768  			b1:       PointFromCoords(-1, 1, 0),
   769  			distRads: math.Pi,
   770  		},
   771  	}
   772  
   773  	for _, test := range tests {
   774  		// Given two edges a0a1 and b0b1, check that the maximum distance between them
   775  		// is distancerads.
   776  		if maxDist, ok := updateEdgePairMaxDistance(test.a0, test.a1, test.b0, test.b1, s1.StraightChordAngle); ok {
   777  			t.Errorf("updateEdgePairMaxDistance(%v, %v, %v, %v, %v) = %v, want updated to be false", test.a0, test.a1, test.b0, test.b1, s1.StraightChordAngle, maxDist)
   778  		}
   779  
   780  		maxDist, ok := updateEdgePairMaxDistance(test.a0, test.a1, test.b0, test.b1, s1.NegativeChordAngle)
   781  		if !ok {
   782  			t.Errorf("updateEdgePairMaxDistance(%v, %v, %v, %v, %v) = %v, want updated to be false", test.a0, test.a1, test.b0, test.b1, s1.NegativeChordAngle, maxDist)
   783  		}
   784  		if !float64Near(test.distRads, maxDist.Angle().Radians(), epsilon) {
   785  			t.Errorf("maxDist %v - %v = %v, want < %v", test.distRads, maxDist.Angle().Radians(), (test.distRads - maxDist.Angle().Radians()), epsilon)
   786  
   787  		}
   788  	}
   789  }
   790  
   791  // TestEdgeDistancesEdgeBNearEdgeA
   792  

View as plain text