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

Documentation: github.com/golang/geo/s2

     1  // Copyright 2014 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/r1"
    22  	"github.com/golang/geo/r2"
    23  	"github.com/golang/geo/r3"
    24  	"github.com/golang/geo/s1"
    25  )
    26  
    27  func TestRectEmptyAndFull(t *testing.T) {
    28  	tests := []struct {
    29  		rect  Rect
    30  		valid bool
    31  		empty bool
    32  		full  bool
    33  		point bool
    34  	}{
    35  		{EmptyRect(), true, true, false, false},
    36  		{FullRect(), true, false, true, false},
    37  	}
    38  
    39  	for _, test := range tests {
    40  		if got := test.rect.IsValid(); got != test.valid {
    41  			t.Errorf("%v.IsValid() = %v, want %v", test.rect, got, test.valid)
    42  		}
    43  		if got := test.rect.IsEmpty(); got != test.empty {
    44  			t.Errorf("%v.IsEmpty() = %v, want %v", test.rect, got, test.empty)
    45  		}
    46  		if got := test.rect.IsFull(); got != test.full {
    47  			t.Errorf("%v.IsFull() = %v, want %v", test.rect, got, test.full)
    48  		}
    49  		if got := test.rect.IsPoint(); got != test.point {
    50  			t.Errorf("%v.IsPoint() = %v, want %v", test.rect, got, test.point)
    51  		}
    52  	}
    53  }
    54  
    55  func TestRectArea(t *testing.T) {
    56  	tests := []struct {
    57  		rect Rect
    58  		want float64
    59  	}{
    60  		{Rect{}, 0},
    61  		{FullRect(), 4 * math.Pi},
    62  		{Rect{r1.Interval{0, math.Pi / 2}, s1.Interval{0, math.Pi / 2}}, math.Pi / 2},
    63  	}
    64  	for _, test := range tests {
    65  		if got := test.rect.Area(); !float64Eq(got, test.want) {
    66  			t.Errorf("%v.Area() = %v, want %v", test.rect, got, test.want)
    67  		}
    68  	}
    69  }
    70  
    71  func TestRectString(t *testing.T) {
    72  	const want = "[Lo[-90.0000000, -180.0000000], Hi[90.0000000, 180.0000000]]"
    73  	if s := FullRect().String(); s != want {
    74  		t.Errorf("FullRect().String() = %q, want %q", s, want)
    75  	}
    76  }
    77  
    78  func TestRectFromLatLng(t *testing.T) {
    79  	ll := LatLngFromDegrees(23, 47)
    80  	got := RectFromLatLng(ll)
    81  	if got.Center() != ll {
    82  		t.Errorf("RectFromLatLng(%v).Center() = %v, want %v", ll, got.Center(), ll)
    83  	}
    84  	if !got.IsPoint() {
    85  		t.Errorf("RectFromLatLng(%v) = %v, want a point", ll, got)
    86  	}
    87  }
    88  
    89  func rectFromDegrees(latLo, lngLo, latHi, lngHi float64) Rect {
    90  	// Convenience method to construct a rectangle. This method is
    91  	// intentionally *not* in the S2LatLngRect interface because the
    92  	// argument order is ambiguous, but is fine for the test.
    93  	return Rect{
    94  		Lat: r1.Interval{
    95  			Lo: (s1.Angle(latLo) * s1.Degree).Radians(),
    96  			Hi: (s1.Angle(latHi) * s1.Degree).Radians(),
    97  		},
    98  		Lng: s1.IntervalFromEndpoints(
    99  			(s1.Angle(lngLo) * s1.Degree).Radians(),
   100  			(s1.Angle(lngHi) * s1.Degree).Radians(),
   101  		),
   102  	}
   103  }
   104  
   105  func TestRectFromCenterSize(t *testing.T) {
   106  	tests := []struct {
   107  		center, size LatLng
   108  		want         Rect
   109  	}{
   110  		{
   111  			LatLngFromDegrees(80, 170),
   112  			LatLngFromDegrees(40, 60),
   113  			rectFromDegrees(60, 140, 90, -160),
   114  		},
   115  		{
   116  			LatLngFromDegrees(10, 40),
   117  			LatLngFromDegrees(210, 400),
   118  			FullRect(),
   119  		},
   120  		{
   121  			LatLngFromDegrees(-90, 180),
   122  			LatLngFromDegrees(20, 50),
   123  			rectFromDegrees(-90, 155, -80, -155),
   124  		},
   125  	}
   126  	for _, test := range tests {
   127  		if got := RectFromCenterSize(test.center, test.size); !rectsApproxEqual(got, test.want, epsilon, epsilon) {
   128  			t.Errorf("RectFromCenterSize(%v,%v) was %v, want %v", test.center, test.size, got, test.want)
   129  		}
   130  	}
   131  }
   132  
   133  func TestRectAddPoint(t *testing.T) {
   134  	tests := []struct {
   135  		input Rect
   136  		point LatLng
   137  		want  Rect
   138  	}{
   139  		{
   140  			Rect{r1.EmptyInterval(), s1.EmptyInterval()},
   141  			LatLngFromDegrees(0, 0),
   142  			rectFromDegrees(0, 0, 0, 0),
   143  		},
   144  		{
   145  			rectFromDegrees(0, 0, 0, 0),
   146  			LatLng{0 * s1.Radian, (-math.Pi / 2) * s1.Radian},
   147  			rectFromDegrees(0, -90, 0, 0),
   148  		},
   149  		{
   150  			rectFromDegrees(0, -90, 0, 0),
   151  			LatLng{(math.Pi / 4) * s1.Radian, (-math.Pi) * s1.Radian},
   152  			rectFromDegrees(0, -180, 45, 0),
   153  		},
   154  		{
   155  			rectFromDegrees(0, -180, 45, 0),
   156  			LatLng{(math.Pi / 2) * s1.Radian, 0 * s1.Radian},
   157  			rectFromDegrees(0, -180, 90, 0),
   158  		},
   159  	}
   160  	for _, test := range tests {
   161  		if got, want := test.input.AddPoint(test.point), test.want; !rectsApproxEqual(got, want, epsilon, epsilon) {
   162  			t.Errorf("%v.AddPoint(%v) was %v, want %v", test.input, test.point, got, want)
   163  		}
   164  	}
   165  }
   166  func TestRectVertex(t *testing.T) {
   167  	r1 := Rect{r1.Interval{0, math.Pi / 2}, s1.IntervalFromEndpoints(-math.Pi, 0)}
   168  	tests := []struct {
   169  		r    Rect
   170  		i    int
   171  		want LatLng
   172  	}{
   173  		{r1, 0, LatLng{0, math.Pi}},
   174  		{r1, 1, LatLng{0, 0}},
   175  		{r1, 2, LatLng{math.Pi / 2, 0}},
   176  		{r1, 3, LatLng{math.Pi / 2, math.Pi}},
   177  	}
   178  
   179  	for _, test := range tests {
   180  		if got := test.r.Vertex(test.i); got != test.want {
   181  			t.Errorf("%v.Vertex(%d) = %v, want %v", test.r, test.i, got, test.want)
   182  		}
   183  	}
   184  }
   185  func TestRectVertexCCWOrder(t *testing.T) {
   186  	for i := 0; i < 4; i++ {
   187  		lat := math.Pi / 4 * float64(i-2)
   188  		lng := math.Pi/2*float64(i-2) + 0.2
   189  		r := Rect{
   190  			r1.Interval{lat, lat + math.Pi/4},
   191  			s1.Interval{
   192  				math.Remainder(lng, 2*math.Pi),
   193  				math.Remainder(lng+math.Pi/2, 2*math.Pi),
   194  			},
   195  		}
   196  
   197  		for k := 0; k < 4; k++ {
   198  			if !Sign(PointFromLatLng(r.Vertex((k-1)&3)), PointFromLatLng(r.Vertex(k)), PointFromLatLng(r.Vertex((k+1)&3))) {
   199  				t.Errorf("%v.Vertex(%v), vertices were not in CCW order", r, k)
   200  			}
   201  		}
   202  	}
   203  }
   204  
   205  func TestRectContainsLatLng(t *testing.T) {
   206  	tests := []struct {
   207  		input Rect
   208  		ll    LatLng
   209  		want  bool
   210  	}{
   211  		{
   212  			rectFromDegrees(0, -180, 90, 0),
   213  			LatLngFromDegrees(30, -45),
   214  			true,
   215  		},
   216  		{
   217  			rectFromDegrees(0, -180, 90, 0),
   218  			LatLngFromDegrees(30, 45),
   219  			false,
   220  		},
   221  		{
   222  			rectFromDegrees(0, -180, 90, 0),
   223  			LatLngFromDegrees(0, -180),
   224  			true,
   225  		},
   226  		{
   227  			rectFromDegrees(0, -180, 90, 0),
   228  			LatLngFromDegrees(90, 0),
   229  			true,
   230  		},
   231  	}
   232  	for _, test := range tests {
   233  		if got, want := test.input.ContainsLatLng(test.ll), test.want; got != want {
   234  			t.Errorf("%v.ContainsLatLng(%v) was %v, want %v", test.input, test.ll, got, want)
   235  		}
   236  	}
   237  }
   238  
   239  func TestRectExpanded(t *testing.T) {
   240  	tests := []struct {
   241  		input  Rect
   242  		margin LatLng
   243  		want   Rect
   244  	}{
   245  		{
   246  			rectFromDegrees(70, 150, 80, 170),
   247  			LatLngFromDegrees(20, 30),
   248  			rectFromDegrees(50, 120, 90, -160),
   249  		},
   250  		{
   251  			EmptyRect(),
   252  			LatLngFromDegrees(20, 30),
   253  			EmptyRect(),
   254  		},
   255  		{
   256  			FullRect(),
   257  			LatLngFromDegrees(500, 500),
   258  			FullRect(),
   259  		},
   260  		{
   261  			rectFromDegrees(-90, 170, 10, 20),
   262  			LatLngFromDegrees(30, 80),
   263  			rectFromDegrees(-90, -180, 40, 180),
   264  		},
   265  
   266  		// Negative margins.
   267  		{
   268  			rectFromDegrees(10, -50, 60, 70),
   269  			LatLngFromDegrees(-10, -10),
   270  			rectFromDegrees(20, -40, 50, 60),
   271  		},
   272  		{
   273  			rectFromDegrees(-20, -180, 20, 180),
   274  			LatLngFromDegrees(-10, -10),
   275  			rectFromDegrees(-10, -180, 10, 180),
   276  		},
   277  		{
   278  			rectFromDegrees(-20, -180, 20, 180),
   279  			LatLngFromDegrees(-30, -30),
   280  			EmptyRect(),
   281  		},
   282  		{
   283  			rectFromDegrees(-90, 10, 90, 11),
   284  			LatLngFromDegrees(-10, -10),
   285  			EmptyRect(),
   286  		},
   287  		{
   288  			rectFromDegrees(-90, 10, 90, 100),
   289  			LatLngFromDegrees(-10, -10),
   290  			rectFromDegrees(-80, 20, 80, 90),
   291  		},
   292  		{
   293  			EmptyRect(),
   294  			LatLngFromDegrees(-50, -500),
   295  			EmptyRect(),
   296  		},
   297  		{
   298  			FullRect(),
   299  			LatLngFromDegrees(-50, -50),
   300  			rectFromDegrees(-40, -180, 40, 180),
   301  		},
   302  
   303  		// Mixed margins.
   304  		{
   305  			rectFromDegrees(10, -50, 60, 70),
   306  			LatLngFromDegrees(-10, 30),
   307  			rectFromDegrees(20, -80, 50, 100),
   308  		},
   309  		{
   310  			rectFromDegrees(-20, -180, 20, 180),
   311  			LatLngFromDegrees(10, -500),
   312  			rectFromDegrees(-30, -180, 30, 180),
   313  		},
   314  		{
   315  			rectFromDegrees(-90, -180, 80, 180),
   316  			LatLngFromDegrees(-30, 500),
   317  			rectFromDegrees(-60, -180, 50, 180),
   318  		},
   319  		{
   320  			rectFromDegrees(-80, -100, 80, 150),
   321  			LatLngFromDegrees(30, -50),
   322  			rectFromDegrees(-90, -50, 90, 100),
   323  		},
   324  		{
   325  			rectFromDegrees(0, -180, 50, 180),
   326  			LatLngFromDegrees(-30, 500),
   327  			EmptyRect(),
   328  		},
   329  		{
   330  			rectFromDegrees(-80, 10, 70, 20),
   331  			LatLngFromDegrees(30, -200),
   332  			EmptyRect(),
   333  		},
   334  		{
   335  			EmptyRect(),
   336  			LatLngFromDegrees(100, -100),
   337  			EmptyRect(),
   338  		},
   339  		{
   340  			FullRect(),
   341  			LatLngFromDegrees(100, -100),
   342  			FullRect(),
   343  		},
   344  	}
   345  	for _, test := range tests {
   346  		if got, want := test.input.expanded(test.margin), test.want; !rectsApproxEqual(got, want, epsilon, epsilon) {
   347  			t.Errorf("%v.Expanded(%v) was %v, want %v", test.input, test.margin, got, want)
   348  		}
   349  	}
   350  }
   351  
   352  func TestRectPolarClosure(t *testing.T) {
   353  	tests := []struct {
   354  		r    Rect
   355  		want Rect
   356  	}{
   357  		{
   358  			rectFromDegrees(-89, 0, 89, 1),
   359  			rectFromDegrees(-89, 0, 89, 1),
   360  		},
   361  		{
   362  			rectFromDegrees(-90, -30, -45, 100),
   363  			rectFromDegrees(-90, -180, -45, 180),
   364  		},
   365  		{
   366  			rectFromDegrees(89, 145, 90, 146),
   367  			rectFromDegrees(89, -180, 90, 180),
   368  		},
   369  		{
   370  			rectFromDegrees(-90, -145, 90, -144),
   371  			FullRect(),
   372  		},
   373  	}
   374  	for _, test := range tests {
   375  		if got := test.r.PolarClosure(); !rectsApproxEqual(got, test.want, epsilon, epsilon) {
   376  			t.Errorf("%v.PolarClosure() was %v, want %v", test.r, got, test.want)
   377  		}
   378  	}
   379  }
   380  
   381  func TestRectCapBound(t *testing.T) {
   382  	tests := []struct {
   383  		r    Rect
   384  		want Cap
   385  	}{
   386  		{ // Bounding cap at center is smaller.
   387  			rectFromDegrees(-45, -45, 45, 45),
   388  			CapFromCenterHeight(Point{r3.Vector{1, 0, 0}}, 0.5),
   389  		},
   390  		{ // Bounding cap at north pole is smaller.
   391  			rectFromDegrees(88, -80, 89, 80),
   392  			CapFromCenterAngle(Point{r3.Vector{0, 0, 1}}, s1.Angle(2)*s1.Degree),
   393  		},
   394  		{ // Longitude span > 180 degrees.
   395  			rectFromDegrees(-30, -150, -10, 50),
   396  			CapFromCenterAngle(Point{r3.Vector{0, 0, -1}}, s1.Angle(80)*s1.Degree),
   397  		},
   398  	}
   399  	for _, test := range tests {
   400  		if got := test.r.CapBound(); !test.want.ApproxEqual(got) {
   401  			t.Errorf("%v.CapBound() was %v, want %v", test.r, got, test.want)
   402  		}
   403  	}
   404  }
   405  
   406  func TestRectIntervalOps(t *testing.T) {
   407  	// Rectangle that covers one-quarter of the sphere.
   408  	rect := rectFromDegrees(0, -180, 90, 0)
   409  
   410  	// Test operations where one rectangle consists of a single point.
   411  	rectMid := rectFromDegrees(45, -90, 45, -90)
   412  	rect180 := rectFromDegrees(0, -180, 0, -180)
   413  	northPole := rectFromDegrees(90, 0, 90, 0)
   414  
   415  	tests := []struct {
   416  		rect         Rect
   417  		other        Rect
   418  		contains     bool
   419  		intersects   bool
   420  		union        Rect
   421  		intersection Rect
   422  	}{
   423  		{
   424  			rect:         rect,
   425  			other:        rectMid,
   426  			contains:     true,
   427  			intersects:   true,
   428  			union:        rect,
   429  			intersection: rectMid,
   430  		},
   431  		{
   432  			rect:         rect,
   433  			other:        rect180,
   434  			contains:     true,
   435  			intersects:   true,
   436  			union:        rect,
   437  			intersection: rect180,
   438  		},
   439  		{
   440  			rect:         rect,
   441  			other:        northPole,
   442  			contains:     true,
   443  			intersects:   true,
   444  			union:        rect,
   445  			intersection: northPole,
   446  		},
   447  		{
   448  			rect:         rect,
   449  			other:        rectFromDegrees(-10, -1, 1, 20),
   450  			contains:     false,
   451  			intersects:   true,
   452  			union:        rectFromDegrees(-10, 180, 90, 20),
   453  			intersection: rectFromDegrees(0, -1, 1, 0),
   454  		},
   455  		{
   456  			rect:         rect,
   457  			other:        rectFromDegrees(-10, -1, 0, 20),
   458  			contains:     false,
   459  			intersects:   true,
   460  			union:        rectFromDegrees(-10, 180, 90, 20),
   461  			intersection: rectFromDegrees(0, -1, 0, 0),
   462  		},
   463  		{
   464  			rect:         rect,
   465  			other:        rectFromDegrees(-10, 0, 1, 20),
   466  			contains:     false,
   467  			intersects:   true,
   468  			union:        rectFromDegrees(-10, 180, 90, 20),
   469  			intersection: rectFromDegrees(0, 0, 1, 0),
   470  		},
   471  		{
   472  			rect:         rectFromDegrees(-15, -160, -15, -150),
   473  			other:        rectFromDegrees(20, 145, 25, 155),
   474  			contains:     false,
   475  			intersects:   false,
   476  			union:        rectFromDegrees(-15, 145, 25, -150),
   477  			intersection: EmptyRect(),
   478  		},
   479  		{
   480  			rect:         rectFromDegrees(70, -10, 90, -140),
   481  			other:        rectFromDegrees(60, 175, 80, 5),
   482  			contains:     false,
   483  			intersects:   true,
   484  			union:        rectFromDegrees(60, -180, 90, 180),
   485  			intersection: rectFromDegrees(70, 175, 80, 5),
   486  		},
   487  
   488  		// Check that the intersection of two rectangles that overlap in latitude
   489  		// but not longitude is valid, and vice versa.
   490  		{
   491  			rect:         rectFromDegrees(12, 30, 60, 60),
   492  			other:        rectFromDegrees(0, 0, 30, 18),
   493  			contains:     false,
   494  			intersects:   false,
   495  			union:        rectFromDegrees(0, 0, 60, 60),
   496  			intersection: EmptyRect(),
   497  		},
   498  		{
   499  			rect:         rectFromDegrees(0, 0, 18, 42),
   500  			other:        rectFromDegrees(30, 12, 42, 60),
   501  			contains:     false,
   502  			intersects:   false,
   503  			union:        rectFromDegrees(0, 0, 42, 60),
   504  			intersection: EmptyRect(),
   505  		},
   506  	}
   507  	for _, test := range tests {
   508  		if got := test.rect.Contains(test.other); got != test.contains {
   509  			t.Errorf("%v.Contains(%v) = %t, want %t", test.rect, test.other, got, test.contains)
   510  		}
   511  
   512  		if got := test.rect.Intersects(test.other); got != test.intersects {
   513  			t.Errorf("%v.Intersects(%v) = %t, want %t", test.rect, test.other, got, test.intersects)
   514  		}
   515  
   516  		if got := test.rect.Union(test.other) == test.rect; test.rect.Contains(test.other) != got {
   517  			t.Errorf("%v.Union(%v) == %v = %t, want %t",
   518  				test.rect, test.other, test.other, got, test.rect.Contains(test.other),
   519  			)
   520  		}
   521  
   522  		if got := test.rect.Intersection(test.other).IsEmpty(); test.rect.Intersects(test.other) == got {
   523  			t.Errorf("%v.Intersection(%v).IsEmpty() = %t, want %t",
   524  				test.rect, test.other, got, test.rect.Intersects(test.other))
   525  		}
   526  
   527  		if got := test.rect.Union(test.other); got != test.union {
   528  			t.Errorf("%v.Union(%v) = %v, want %v", test.rect, test.other, got, test.union)
   529  		}
   530  
   531  		if got := test.rect.Intersection(test.other); got != test.intersection {
   532  			t.Errorf("%v.Intersection(%v) = %v, want %v", test.rect, test.other, got, test.intersection)
   533  		}
   534  	}
   535  }
   536  
   537  func TestRectCellOps(t *testing.T) {
   538  	cell0 := CellFromPoint(Point{r3.Vector{1 + 1e-12, 1, 1}})
   539  	v0 := LatLngFromPoint(cell0.Vertex(0))
   540  
   541  	cell202 := CellFromCellID(CellIDFromFacePosLevel(2, 0, 2))
   542  	bound202 := cell202.RectBound()
   543  
   544  	tests := []struct {
   545  		r          Rect
   546  		c          Cell
   547  		contains   bool
   548  		intersects bool
   549  	}{
   550  		// Special cases
   551  		{
   552  			r:          EmptyRect(),
   553  			c:          CellFromCellID(CellIDFromFacePosLevel(3, 0, 0)),
   554  			contains:   false,
   555  			intersects: false,
   556  		},
   557  		{
   558  			r:          FullRect(),
   559  			c:          CellFromCellID(CellIDFromFacePosLevel(2, 0, 0)),
   560  			contains:   true,
   561  			intersects: true,
   562  		},
   563  		{
   564  			r:          FullRect(),
   565  			c:          CellFromCellID(CellIDFromFacePosLevel(5, 0, 25)),
   566  			contains:   true,
   567  			intersects: true,
   568  		},
   569  		// This rectangle includes the first quadrant of face 0.  It's expanded
   570  		// slightly because cell bounding rectangles are slightly conservative.
   571  		{
   572  			r:          rectFromDegrees(-45.1, -45.1, 0.1, 0.1),
   573  			c:          CellFromCellID(CellIDFromFacePosLevel(0, 0, 0)),
   574  			contains:   false,
   575  			intersects: true,
   576  		},
   577  		{
   578  			r:          rectFromDegrees(-45.1, -45.1, 0.1, 0.1),
   579  			c:          CellFromCellID(CellIDFromFacePosLevel(0, 0, 1)),
   580  			contains:   true,
   581  			intersects: true,
   582  		},
   583  		{
   584  			r:          rectFromDegrees(-45.1, -45.1, 0.1, 0.1),
   585  			c:          CellFromCellID(CellIDFromFacePosLevel(1, 0, 1)),
   586  			contains:   false,
   587  			intersects: false,
   588  		},
   589  		// This rectangle intersects the first quadrant of face 0.
   590  		{
   591  			r:          rectFromDegrees(-10, -45, 10, 0),
   592  			c:          CellFromCellID(CellIDFromFacePosLevel(0, 0, 0)),
   593  			contains:   false,
   594  			intersects: true,
   595  		},
   596  		{
   597  			r:          rectFromDegrees(-10, -45, 10, 0),
   598  			c:          CellFromCellID(CellIDFromFacePosLevel(0, 0, 1)),
   599  			contains:   false,
   600  			intersects: true,
   601  		},
   602  		{
   603  			r:          rectFromDegrees(-10, -45, 10, 0),
   604  			c:          CellFromCellID(CellIDFromFacePosLevel(1, 0, 1)),
   605  			contains:   false,
   606  			intersects: false,
   607  		},
   608  		// Rectangle consisting of a single point.
   609  		{
   610  			r:          rectFromDegrees(4, 4, 4, 4),
   611  			c:          CellFromCellID(CellIDFromFace(0)),
   612  			contains:   false,
   613  			intersects: true,
   614  		},
   615  		// Rectangles that intersect the bounding rectangle of a face
   616  		// but not the face itself.
   617  		{
   618  			r:          rectFromDegrees(41, -87, 42, -79),
   619  			c:          CellFromCellID(CellIDFromFace(2)),
   620  			contains:   false,
   621  			intersects: false,
   622  		},
   623  		{
   624  			r:          rectFromDegrees(-41, 160, -40, -160),
   625  			c:          CellFromCellID(CellIDFromFace(5)),
   626  			contains:   false,
   627  			intersects: false,
   628  		},
   629  		{
   630  			// This is the leaf cell at the top right hand corner of face 0.
   631  			// It has two angles of 60 degrees and two of 120 degrees.
   632  			r: rectFromDegrees(v0.Lat.Degrees()-1e-8,
   633  				v0.Lng.Degrees()-1e-8,
   634  				v0.Lat.Degrees()-2e-10,
   635  				v0.Lng.Degrees()+1e-10),
   636  			c:          cell0,
   637  			contains:   false,
   638  			intersects: false,
   639  		},
   640  		{
   641  			// Rectangles that intersect a face but where no vertex of one region
   642  			// is contained by the other region.  The first one passes through
   643  			// a corner of one of the face cells.
   644  			r:          rectFromDegrees(-37, -70, -36, -20),
   645  			c:          CellFromCellID(CellIDFromFace(5)),
   646  			contains:   false,
   647  			intersects: true,
   648  		},
   649  		{
   650  			// These two intersect like a diamond and a square.
   651  			r: rectFromDegrees(bound202.Lo().Lat.Degrees()+3,
   652  				bound202.Lo().Lng.Degrees()+3,
   653  				bound202.Hi().Lat.Degrees()-3,
   654  				bound202.Hi().Lng.Degrees()-3),
   655  			c:          cell202,
   656  			contains:   false,
   657  			intersects: true,
   658  		},
   659  		{
   660  			// from a bug report
   661  			r:          rectFromDegrees(34.2572864, 135.2673642, 34.2707907, 135.2995742),
   662  			c:          CellFromCellID(0x6007500000000000),
   663  			contains:   false,
   664  			intersects: true,
   665  		},
   666  	}
   667  
   668  	for _, test := range tests {
   669  		if got := test.r.ContainsCell(test.c); got != test.contains {
   670  			t.Errorf("%v.ContainsCell(%v) = %t, want %t", test.r, test.c, got, test.contains)
   671  		}
   672  
   673  		if got := test.r.IntersectsCell(test.c); got != test.intersects {
   674  			t.Errorf("%v.IntersectsCell(%v) = %t, want %t", test.r, test.c, got, test.intersects)
   675  		}
   676  	}
   677  
   678  }
   679  
   680  func TestRectContainsPoint(t *testing.T) {
   681  	r1 := rectFromDegrees(0, -180, 90, 0)
   682  
   683  	tests := []struct {
   684  		r    Rect
   685  		p    Point
   686  		want bool
   687  	}{
   688  		{r1, Point{r3.Vector{0.5, -0.3, 0.1}}, true},
   689  		{r1, Point{r3.Vector{0.5, 0.2, 0.1}}, false},
   690  	}
   691  	for _, test := range tests {
   692  		if got, want := test.r.ContainsPoint(test.p), test.want; got != want {
   693  			t.Errorf("%v.ContainsPoint(%v) was %v, want %v", test.r, test.p, got, want)
   694  		}
   695  	}
   696  }
   697  
   698  func TestRectIntersectsLatEdge(t *testing.T) {
   699  	tests := []struct {
   700  		a, b  Point
   701  		lat   s1.Angle
   702  		lngLo s1.Angle
   703  		lngHi s1.Angle
   704  		want  bool
   705  	}{
   706  		{
   707  			a:     Point{r3.Vector{-1, -1, 1}},
   708  			b:     Point{r3.Vector{1, -1, 1}},
   709  			lat:   41 * s1.Degree,
   710  			lngLo: -87 * s1.Degree,
   711  			lngHi: -79 * s1.Degree,
   712  			want:  false,
   713  		},
   714  		{
   715  			a:     Point{r3.Vector{-1, -1, 1}},
   716  			b:     Point{r3.Vector{1, -1, 1}},
   717  			lat:   42 * s1.Degree,
   718  			lngLo: -87 * s1.Degree,
   719  			lngHi: -79 * s1.Degree,
   720  			want:  false,
   721  		},
   722  		{
   723  			a:     Point{r3.Vector{-1, -1, -1}},
   724  			b:     Point{r3.Vector{1, 1, 0}},
   725  			lat:   -3 * s1.Degree,
   726  			lngLo: -1 * s1.Degree,
   727  			lngHi: 23 * s1.Degree,
   728  			want:  false,
   729  		},
   730  		{
   731  			a:     Point{r3.Vector{1, 0, 1}},
   732  			b:     Point{r3.Vector{1, -1, 0}},
   733  			lat:   -28 * s1.Degree,
   734  			lngLo: 69 * s1.Degree,
   735  			lngHi: 115 * s1.Degree,
   736  			want:  false,
   737  		},
   738  		{
   739  			a:     Point{r3.Vector{0, 1, 0}},
   740  			b:     Point{r3.Vector{1, -1, -1}},
   741  			lat:   44 * s1.Degree,
   742  			lngLo: 60 * s1.Degree,
   743  			lngHi: 177 * s1.Degree,
   744  			want:  false,
   745  		},
   746  		{
   747  			a:     Point{r3.Vector{0, 1, 1}},
   748  			b:     Point{r3.Vector{0, 1, -1}},
   749  			lat:   -25 * s1.Degree,
   750  			lngLo: -74 * s1.Degree,
   751  			lngHi: -165 * s1.Degree,
   752  			want:  true,
   753  		},
   754  		{
   755  			a:     Point{r3.Vector{1, 0, 0}},
   756  			b:     Point{r3.Vector{0, 0, -1}},
   757  			lat:   -4 * s1.Degree,
   758  			lngLo: -152 * s1.Degree,
   759  			lngHi: 171 * s1.Degree,
   760  			want:  true,
   761  		},
   762  		// from a bug report
   763  		{
   764  			a:     Point{r3.Vector{-0.589375791872893683986945, 0.583248451588733285433364, 0.558978908075738245564423}},
   765  			b:     Point{r3.Vector{-0.587388131301997518107783, 0.581281455376392863776402, 0.563104832905072516524569}},
   766  			lat:   34.2572864 * s1.Degree,
   767  			lngLo: 2.3608609 * s1.Radian,
   768  			lngHi: 2.3614230 * s1.Radian,
   769  			want:  true,
   770  		},
   771  	}
   772  
   773  	for _, test := range tests {
   774  		if got := intersectsLatEdge(test.a, test.b, test.lat, s1.Interval{float64(test.lngLo), float64(test.lngHi)}); got != test.want {
   775  			t.Errorf("intersectsLatEdge(%v, %v, %v, {%v, %v}) = %t, want %t",
   776  				test.a, test.b, test.lat, test.lngLo, test.lngHi, got, test.want)
   777  		}
   778  	}
   779  }
   780  
   781  func TestRectIntersectsLngEdge(t *testing.T) {
   782  	tests := []struct {
   783  		a, b  Point
   784  		latLo s1.Angle
   785  		latHi s1.Angle
   786  		lng   s1.Angle
   787  		want  bool
   788  	}{
   789  		{
   790  			a:     Point{r3.Vector{-1, -1, 1}},
   791  			b:     Point{r3.Vector{1, -1, 1}},
   792  			latLo: 41 * s1.Degree,
   793  			latHi: 42 * s1.Degree,
   794  			lng:   -79 * s1.Degree,
   795  			want:  false,
   796  		},
   797  		{
   798  			a:     Point{r3.Vector{-1, -1, 1}},
   799  			b:     Point{r3.Vector{1, -1, 1}},
   800  			latLo: 41 * s1.Degree,
   801  			latHi: 42 * s1.Degree,
   802  			lng:   -87 * s1.Degree,
   803  			want:  false,
   804  		},
   805  		{
   806  			a:     Point{r3.Vector{-1, -1, 1}},
   807  			b:     Point{r3.Vector{1, -1, 1}},
   808  			latLo: 42 * s1.Degree,
   809  			latHi: 41 * s1.Degree,
   810  			lng:   79 * s1.Degree,
   811  			want:  false,
   812  		},
   813  		{
   814  			a:     Point{r3.Vector{-1, -1, 1}},
   815  			b:     Point{r3.Vector{1, -1, 1}},
   816  			latLo: 41 * s1.Degree,
   817  			latHi: 42 * s1.Degree,
   818  			lng:   87 * s1.Degree,
   819  			want:  false,
   820  		},
   821  		{
   822  			a:     Point{r3.Vector{0, -1, -1}},
   823  			b:     Point{r3.Vector{-1, 0, -1}},
   824  			latLo: -87 * s1.Degree,
   825  			latHi: 13 * s1.Degree,
   826  			lng:   -143 * s1.Degree,
   827  			want:  true,
   828  		},
   829  		{
   830  			a:     Point{r3.Vector{1, 1, -1}},
   831  			b:     Point{r3.Vector{1, -1, 1}},
   832  			latLo: -64 * s1.Degree,
   833  			latHi: 13 * s1.Degree,
   834  			lng:   40 * s1.Degree,
   835  			want:  true,
   836  		},
   837  		{
   838  			a:     Point{r3.Vector{1, 1, 0}},
   839  			b:     Point{r3.Vector{-1, 0, -1}},
   840  			latLo: -64 * s1.Degree,
   841  			latHi: 56 * s1.Degree,
   842  			lng:   151 * s1.Degree,
   843  			want:  true,
   844  		},
   845  		{
   846  			a:     Point{r3.Vector{-1, -1, 0}},
   847  			b:     Point{r3.Vector{1, -1, -1}},
   848  			latLo: -50 * s1.Degree,
   849  			latHi: 18 * s1.Degree,
   850  			lng:   -84 * s1.Degree,
   851  			want:  true,
   852  		},
   853  	}
   854  
   855  	for _, test := range tests {
   856  		if got := intersectsLngEdge(test.a, test.b, r1.Interval{float64(test.latLo), float64(test.latHi)}, test.lng); got != test.want {
   857  			t.Errorf("intersectsLngEdge(%v, %v, {%v, %v}, %v) = %v, want %v",
   858  				test.a, test.b, test.latLo, test.latHi, test.lng, got, test.want)
   859  		}
   860  	}
   861  }
   862  
   863  // intervalDistance returns the minimum distance (in radians) from X to the latitude
   864  // line segment defined by the given latitude and longitude interval.
   865  func intervalDistance(x LatLng, lat s1.Angle, iv s1.Interval) s1.Angle {
   866  	// Is x inside the longitude interval?
   867  	if iv.Contains(float64(x.Lng)) {
   868  		return s1.Angle(math.Abs(float64(x.Lat - lat)))
   869  	}
   870  
   871  	return minAngle(
   872  		x.Distance(LatLng{lat, s1.Angle(iv.Lo)}),
   873  		x.Distance(LatLng{lat, s1.Angle(iv.Hi)}))
   874  }
   875  
   876  // Returns the minimum distance from X to the latitude line segment defined by
   877  // the given latitude and longitude interval.
   878  func bruteForceRectLatLngDistance(r Rect, ll LatLng) s1.Angle {
   879  	pt := PointFromLatLng(ll)
   880  	if r.ContainsPoint(pt) {
   881  		return 0
   882  	}
   883  
   884  	loLat := intervalDistance(ll, s1.Angle(r.Lat.Lo), r.Lng)
   885  	hiLat := intervalDistance(ll, s1.Angle(r.Lat.Hi), r.Lng)
   886  	loLng := DistanceFromSegment(PointFromLatLng(ll),
   887  		PointFromLatLng(LatLng{s1.Angle(r.Lat.Lo), s1.Angle(r.Lng.Lo)}),
   888  		PointFromLatLng(LatLng{s1.Angle(r.Lat.Hi), s1.Angle(r.Lng.Lo)}))
   889  	hiLng := DistanceFromSegment(PointFromLatLng(ll),
   890  		PointFromLatLng(LatLng{s1.Angle(r.Lat.Lo), s1.Angle(r.Lng.Hi)}),
   891  		PointFromLatLng(LatLng{s1.Angle(r.Lat.Hi), s1.Angle(r.Lng.Hi)}))
   892  
   893  	return minAngle(loLat, hiLat, loLng, hiLng)
   894  }
   895  
   896  func TestDistanceRectFromLatLng(t *testing.T) {
   897  	// Rect that spans 180.
   898  	a := RectFromLatLng(LatLngFromDegrees(-1, -1)).AddPoint(LatLngFromDegrees(2, 1))
   899  	// Rect near north pole.
   900  	b := RectFromLatLng(LatLngFromDegrees(86, 0)).AddPoint(LatLngFromDegrees(88, 2))
   901  	// Rect that touches north pole.
   902  	c := RectFromLatLng(LatLngFromDegrees(88, 0)).AddPoint(LatLngFromDegrees(90, 2))
   903  
   904  	tests := []struct {
   905  		r        Rect
   906  		lat, lng float64 // In degrees.
   907  	}{
   908  		{a, -2, -1},
   909  		{a, 1, 2},
   910  		{b, 87, 3},
   911  		{b, 87, -1},
   912  		{b, 89, 1},
   913  		{b, 89, 181},
   914  		{b, 85, 1},
   915  		{b, 85, 181},
   916  		{b, 90, 0},
   917  		{c, 89, 3},
   918  		{c, 89, 90},
   919  		{c, 89, 181},
   920  	}
   921  
   922  	for _, test := range tests {
   923  		ll := LatLngFromDegrees(test.lat, test.lng)
   924  		got := test.r.DistanceToLatLng(ll)
   925  		want := bruteForceRectLatLngDistance(test.r, ll)
   926  		if !float64Near(float64(got), float64(want), 1e-10) {
   927  			t.Errorf("dist from %v to %v = %v, want %v", test.r, ll, got, want)
   928  		}
   929  	}
   930  }
   931  
   932  func TestDistanceRectFromLatLngRandomPairs(t *testing.T) {
   933  	latlng := func() LatLng { return LatLngFromPoint(randomPoint()) }
   934  
   935  	for i := 0; i < 10000; i++ {
   936  		r := RectFromLatLng(latlng()).AddPoint(latlng())
   937  		ll := latlng()
   938  		got := r.DistanceToLatLng(ll)
   939  		want := bruteForceRectLatLngDistance(r, ll)
   940  		if !float64Near(float64(got), float64(want), 1e-10) {
   941  			t.Errorf("dist from %v to %v = %v, want %v", r, ll, got, want)
   942  		}
   943  	}
   944  }
   945  
   946  // This function assumes that DirectedHausdorffDistance() always returns
   947  // a distance from some point in a to b. So the function mainly tests whether
   948  // the returned distance is large enough, and only does a weak test on whether
   949  // it is small enough.
   950  func verifyDirectedHausdorffDistance(t *testing.T, a, b Rect) {
   951  	t.Helper()
   952  
   953  	const resolution = 0.1
   954  
   955  	// Record the max sample distance as well as the sample point realizing the
   956  	// max for easier debugging.
   957  	var maxDistance s1.Angle
   958  
   959  	sampleSizeOnLat := int(a.Lat.Length()/resolution) + 1
   960  	sampleSizeOnLng := int(a.Lng.Length()/resolution) + 1
   961  
   962  	deltaOnLat := s1.Angle(a.Lat.Length()) / s1.Angle(sampleSizeOnLat)
   963  	deltaOnLng := s1.Angle(a.Lng.Length()) / s1.Angle(sampleSizeOnLng)
   964  
   965  	ll := LatLng{Lng: s1.Angle(a.Lng.Lo)}
   966  	for i := 0; i <= sampleSizeOnLng; i++ {
   967  		ll.Lat = s1.Angle(a.Lat.Lo)
   968  
   969  		for j := 0; j <= sampleSizeOnLat; j++ {
   970  			d := b.DistanceToLatLng(ll.Normalized())
   971  			maxDistance = maxAngle(maxDistance, d)
   972  			ll.Lat += deltaOnLat
   973  		}
   974  		ll.Lng += deltaOnLng
   975  	}
   976  
   977  	got := a.DirectedHausdorffDistance(b)
   978  
   979  	if got < maxDistance-1e-10 {
   980  		t.Errorf("hausdorff(%v, %v) = %v < %v-eps, but shouldn't", a, b, got, maxDistance)
   981  	} else if got > maxDistance+resolution {
   982  		t.Errorf("DirectedHausdorffDistance(%v, %v) = %v > %v+resolution, but shouldn't", a, b, got, maxDistance)
   983  	}
   984  }
   985  
   986  func TestRectDirectedHausdorffDistanceRandomPairs(t *testing.T) {
   987  	// Test random pairs.
   988  	rnd := func() LatLng { return LatLngFromPoint(randomPoint()) }
   989  	for i := 0; i < 1000; i++ {
   990  		a := RectFromLatLng(rnd()).AddPoint(rnd())
   991  		b := RectFromLatLng(rnd()).AddPoint(rnd())
   992  		// a and b are *minimum* bounding rectangles of two random points, in
   993  		// particular, their Voronoi diagrams are always of the same topology. We
   994  		// take the "complements" of a and b for more thorough testing.
   995  		a2 := Rect{Lat: a.Lat, Lng: a.Lng.Complement()}
   996  		b2 := Rect{Lat: b.Lat, Lng: b.Lng.Complement()}
   997  
   998  		// Note that "a" and "b" come from the same distribution, so there is no
   999  		// need to test pairs such as (b, a), (b, a2), etc.
  1000  		verifyDirectedHausdorffDistance(t, a, b)
  1001  		verifyDirectedHausdorffDistance(t, a2, b)
  1002  		verifyDirectedHausdorffDistance(t, a, b2)
  1003  		verifyDirectedHausdorffDistance(t, a2, b2)
  1004  	}
  1005  }
  1006  
  1007  func TestDirectedHausdorffDistanceContained(t *testing.T) {
  1008  	// Caller rect is contained in callee rect. Should return 0.
  1009  	a := rectFromDegrees(-10, 20, -5, 90)
  1010  	tests := []Rect{
  1011  		rectFromDegrees(-10, 20, -5, 90),
  1012  		rectFromDegrees(-10, 19, -5, 91),
  1013  		rectFromDegrees(-11, 20, -4, 90),
  1014  		rectFromDegrees(-11, 19, -4, 91),
  1015  	}
  1016  	for _, test := range tests {
  1017  		got, want := a.DirectedHausdorffDistance(test), s1.Angle(0)
  1018  		if got != want {
  1019  			t.Errorf("%v.DirectedHausdorffDistance(%v) = %v, want %v", a, test, got, want)
  1020  		}
  1021  	}
  1022  }
  1023  
  1024  func TestDirectHausdorffDistancePointToRect(t *testing.T) {
  1025  	// The Hausdorff distance from a point to a rect should be the same as its
  1026  	// distance to the rect.
  1027  	a1 := LatLngFromDegrees(5, 8)
  1028  	a2 := LatLngFromDegrees(90, 10) // North pole.
  1029  
  1030  	tests := []struct {
  1031  		ll LatLng
  1032  		b  Rect
  1033  	}{
  1034  		{a1, rectFromDegrees(-85, -50, -80, 10)},
  1035  		{a2, rectFromDegrees(-85, -50, -80, 10)},
  1036  		{a1, rectFromDegrees(4, -10, 80, 10)},
  1037  		{a2, rectFromDegrees(4, -10, 80, 10)},
  1038  		{a1, rectFromDegrees(70, 170, 80, -170)},
  1039  		{a2, rectFromDegrees(70, 170, 80, -170)},
  1040  	}
  1041  	for _, test := range tests {
  1042  		a := RectFromLatLng(test.ll)
  1043  		got, want := a.DirectedHausdorffDistance(test.b), test.b.DistanceToLatLng(test.ll)
  1044  
  1045  		if !float64Eq(float64(got), float64(want)) {
  1046  			t.Errorf("hausdorff(%v, %v) = %v, want %v, as that's the closest dist", test.b, a, got, want)
  1047  		}
  1048  	}
  1049  }
  1050  
  1051  func TestDirectedHausdorffDistanceRectToPoint(t *testing.T) {
  1052  	a := rectFromDegrees(1, -8, 10, 20)
  1053  	tests := []struct {
  1054  		lat, lng float64 // Degrees.
  1055  	}{{5, 8}, {-6, -100}, {-90, -20}, {90, 0}}
  1056  	for _, test := range tests {
  1057  		verifyDirectedHausdorffDistance(t, a, RectFromLatLng(LatLngFromDegrees(test.lat, test.lng)))
  1058  	}
  1059  }
  1060  
  1061  func TestDirectedHausdorffDistanceRectToRectNearPole(t *testing.T) {
  1062  	// Tests near south pole.
  1063  	a := rectFromDegrees(-87, 0, -85, 3)
  1064  	tests := []Rect{
  1065  		rectFromDegrees(-89, 1, -88, 2),
  1066  		rectFromDegrees(-84, 1, -83, 2),
  1067  		rectFromDegrees(-88, 90, -86, 91),
  1068  		rectFromDegrees(-84, -91, -83, -90),
  1069  		rectFromDegrees(-90, 181, -89, 182),
  1070  		rectFromDegrees(-84, 181, -83, 182),
  1071  	}
  1072  	for _, test := range tests {
  1073  		verifyDirectedHausdorffDistance(t, a, test)
  1074  	}
  1075  }
  1076  
  1077  func TestDirectedHausdorffDistanceRectToRectDegenerateCases(t *testing.T) {
  1078  	// Rectangles that contain poles.
  1079  	verifyDirectedHausdorffDistance(t,
  1080  		rectFromDegrees(0, 10, 90, 20), rectFromDegrees(-4, -10, 4, 0))
  1081  	verifyDirectedHausdorffDistance(t,
  1082  		rectFromDegrees(-4, -10, 4, 0), rectFromDegrees(0, 10, 90, 20))
  1083  
  1084  	// Two rectangles share same or complement longitudinal intervals.
  1085  	a := rectFromDegrees(-50, -10, 50, 10)
  1086  	b := rectFromDegrees(30, -10, 60, 10)
  1087  	verifyDirectedHausdorffDistance(t, a, b)
  1088  
  1089  	c := Rect{Lat: a.Lat, Lng: a.Lng.Complement()}
  1090  	verifyDirectedHausdorffDistance(t, c, b)
  1091  
  1092  	// Rectangle a touches b_opposite_lng.
  1093  	verifyDirectedHausdorffDistance(t,
  1094  		rectFromDegrees(10, 170, 30, 180), rectFromDegrees(-50, -10, 50, 10))
  1095  	verifyDirectedHausdorffDistance(t,
  1096  		rectFromDegrees(10, -180, 30, -170), rectFromDegrees(-50, -10, 50, 10))
  1097  
  1098  	// Rectangle b's Voronoi diagram is degenerate (lng interval spans 180
  1099  	// degrees), and a touches the degenerate Voronoi vertex.
  1100  	verifyDirectedHausdorffDistance(t,
  1101  		rectFromDegrees(-30, 170, 30, 180), rectFromDegrees(-10, -90, 10, 90))
  1102  	verifyDirectedHausdorffDistance(t,
  1103  		rectFromDegrees(-30, -180, 30, -170), rectFromDegrees(-10, -90, 10, 90))
  1104  
  1105  	// Rectangle a touches a voronoi vertex of rectangle b.
  1106  	verifyDirectedHausdorffDistance(t,
  1107  		rectFromDegrees(-20, 105, 20, 110), rectFromDegrees(-30, 5, 30, 15))
  1108  	verifyDirectedHausdorffDistance(t,
  1109  		rectFromDegrees(-20, 95, 20, 105), rectFromDegrees(-30, 5, 30, 15))
  1110  }
  1111  
  1112  func TestRectApproxEqual(t *testing.T) {
  1113  	// s1.Interval and r1.Interval have additional testing.
  1114  
  1115  	const ε = epsilon / 10
  1116  	tests := []struct {
  1117  		a, b Rect
  1118  		want bool
  1119  	}{
  1120  		{EmptyRect(), rectFromDegrees(1, 5, 1, 5), true},
  1121  		{rectFromDegrees(1, 5, 1, 5), EmptyRect(), true},
  1122  
  1123  		{rectFromDegrees(1, 5, 1, 5), rectFromDegrees(2, 7, 2, 7), false},
  1124  		{rectFromDegrees(1, 5, 1, 5), rectFromDegrees(1+ε, 5+ε, 1+ε, 5+ε), true},
  1125  	}
  1126  
  1127  	for _, test := range tests {
  1128  		if got := test.a.ApproxEqual(test.b); got != test.want {
  1129  			t.Errorf("%v.ApproxEquals(%v) = %t, want %t", test.a, test.b, got, test.want)
  1130  		}
  1131  	}
  1132  }
  1133  
  1134  func TestRectCentroidEmptyFull(t *testing.T) {
  1135  	// Empty and full rectangles.
  1136  	if got, want := EmptyRect().Centroid(), (Point{}); !got.ApproxEqual(want) {
  1137  		t.Errorf("%v.Centroid() = %v, want %v", EmptyRect(), got, want)
  1138  	}
  1139  	if got, want := FullRect().Centroid().Norm(), epsilon; !float64Eq(got, want) {
  1140  		t.Errorf("%v.Centroid().Norm() = %v, want %v", FullRect(), got, want)
  1141  	}
  1142  }
  1143  
  1144  func testRectCentroidSplitting(t *testing.T, r Rect, leftSplits int) {
  1145  	// Recursively verify that when a rectangle is split into two pieces, the
  1146  	// centroids of the children sum to give the centroid of the parent.
  1147  	var child0, child1 Rect
  1148  	if oneIn(2) {
  1149  		lat := randomUniformFloat64(r.Lat.Lo, r.Lat.Hi)
  1150  		child0 = Rect{r1.Interval{r.Lat.Lo, lat}, r.Lng}
  1151  		child1 = Rect{r1.Interval{lat, r.Lat.Hi}, r.Lng}
  1152  	} else {
  1153  		lng := randomUniformFloat64(r.Lng.Lo, r.Lng.Hi)
  1154  		child0 = Rect{r.Lat, s1.Interval{r.Lng.Lo, lng}}
  1155  		child1 = Rect{r.Lat, s1.Interval{lng, r.Lng.Hi}}
  1156  	}
  1157  
  1158  	if got, want := r.Centroid().Sub(child0.Centroid().Vector).Sub(child1.Centroid().Vector).Norm(), 1e-15; got > want {
  1159  		t.Errorf("%v.Centroid() - %v.Centroid() - %v.Centroid = %v, want ~0", r, child0, child1, got)
  1160  	}
  1161  	if leftSplits > 0 {
  1162  		testRectCentroidSplitting(t, child0, leftSplits-1)
  1163  		testRectCentroidSplitting(t, child1, leftSplits-1)
  1164  	}
  1165  }
  1166  
  1167  func TestRectCentroidFullRange(t *testing.T) {
  1168  	// Rectangles that cover the full longitude range.
  1169  	for i := 0; i < 100; i++ {
  1170  		lat1 := randomUniformFloat64(-math.Pi/2, math.Pi/2)
  1171  		lat2 := randomUniformFloat64(-math.Pi/2, math.Pi/2)
  1172  		r := Rect{r1.Interval{lat1, lat2}, s1.FullInterval()}
  1173  		centroid := r.Centroid()
  1174  		if want := 0.5 * (math.Sin(lat1) + math.Sin(lat2)) * r.Area(); !float64Near(want, centroid.Z, epsilon) {
  1175  			t.Errorf("%v.Centroid().Z was %v, want %v", r, centroid.Z, want)
  1176  		}
  1177  		if got := (r2.Point{centroid.X, centroid.Y}.Norm()); got > epsilon {
  1178  			t.Errorf("%v.Centroid().Norm() was %v, want > %v ", r, got, epsilon)
  1179  		}
  1180  	}
  1181  
  1182  	// Rectangles that cover the full latitude range.
  1183  	for i := 0; i < 100; i++ {
  1184  		lat1 := randomUniformFloat64(-math.Pi, math.Pi)
  1185  		lat2 := randomUniformFloat64(-math.Pi, math.Pi)
  1186  		r := Rect{r1.Interval{-math.Pi / 2, math.Pi / 2}, s1.Interval{lat1, lat2}}
  1187  		centroid := r.Centroid()
  1188  
  1189  		if got, want := math.Abs(centroid.Z), epsilon; got > want {
  1190  			t.Errorf("math.Abs(%v.Centroid().Z) = %v, want <= %v", r, got, want)
  1191  		}
  1192  
  1193  		if got, want := LatLngFromPoint(centroid).Lng.Radians(), r.Lng.Center(); !float64Near(got, want, epsilon) {
  1194  			t.Errorf("%v.Lng.Radians() = %v, want %v", centroid, got, want)
  1195  		}
  1196  
  1197  		alpha := 0.5 * r.Lng.Length()
  1198  		if got, want := (r2.Point{centroid.X, centroid.Y}.Norm()), (0.25 * math.Pi * math.Sin(alpha) / alpha * r.Area()); !float64Near(got, want, epsilon) {
  1199  			t.Errorf("%v.Centroid().Norm() = %v, want ~%v", got, want, epsilon)
  1200  		}
  1201  	}
  1202  
  1203  	// Finally, verify that when a rectangle is recursively split into pieces,
  1204  	// the centroids of the pieces add to give the centroid of their parent.
  1205  	// To make the code simpler we avoid rectangles that cross the 180 degree
  1206  	// line of longitude.
  1207  	testRectCentroidSplitting(t, Rect{r1.Interval{-math.Pi / 2, math.Pi / 2}, s1.Interval{-math.Pi, math.Pi}}, 10)
  1208  }
  1209  

View as plain text