rect_bounder_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/r1"
    22  	"github.com/golang/geo/r3"
    23  	"github.com/golang/geo/s1"
    24  )
    25  
    26  func rectBoundForPoints(a, b Point) Rect {
    27  	bounder := NewRectBounder()
    28  	bounder.AddPoint(a)
    29  	bounder.AddPoint(b)
    30  	return bounder.RectBound()
    31  }
    32  
    33  func TestRectBounderMaxLatitudeSimple(t *testing.T) {
    34  	cubeLat := math.Asin(1 / math.Sqrt(3)) // 35.26 degrees
    35  	cubeLatRect := Rect{r1.IntervalFromPoint(-cubeLat).AddPoint(cubeLat),
    36  		s1.IntervalFromEndpoints(-math.Pi/4, math.Pi/4)}
    37  
    38  	tests := []struct {
    39  		a, b Point
    40  		want Rect
    41  	}{
    42  		// Check cases where the min/max latitude is attained at a vertex.
    43  		{
    44  			a:    Point{r3.Vector{1, 1, 1}},
    45  			b:    Point{r3.Vector{1, -1, -1}},
    46  			want: cubeLatRect,
    47  		},
    48  		{
    49  			a:    Point{r3.Vector{1, -1, 1}},
    50  			b:    Point{r3.Vector{1, 1, -1}},
    51  			want: cubeLatRect,
    52  		},
    53  	}
    54  
    55  	for _, test := range tests {
    56  		if got := rectBoundForPoints(test.a, test.b); !rectsApproxEqual(got, test.want, rectErrorLat, rectErrorLng) {
    57  			t.Errorf("RectBounder for points (%v, %v) near max lat failed: got %v, want %v", test.a, test.b, got, test.want)
    58  		}
    59  	}
    60  }
    61  
    62  func TestRectBounderMaxLatitudeEdgeInterior(t *testing.T) {
    63  	// Check cases where the min/max latitude occurs in the edge interior.
    64  	// These tests expect the result to be pretty close to the middle of the
    65  	// allowable error range (i.e., by adding 0.5 * kRectError).
    66  
    67  	tests := []struct {
    68  		got  float64
    69  		want float64
    70  	}{
    71  		// Max latitude, CW edge
    72  		{
    73  			math.Pi/4 + 0.5*rectErrorLat,
    74  			rectBoundForPoints(Point{r3.Vector{1, 1, 1}}, Point{r3.Vector{1, -1, 1}}).Lat.Hi,
    75  		},
    76  		// Min latitude, CW edge
    77  		{
    78  			-math.Pi/4 - 0.5*rectErrorLat,
    79  			rectBoundForPoints(Point{r3.Vector{1, -1, -1}}, Point{r3.Vector{-1, -1, -1}}).Lat.Lo,
    80  		},
    81  		// Max latitude, CCW edge
    82  		{
    83  			math.Pi/4 + 0.5*rectErrorLat,
    84  			rectBoundForPoints(Point{r3.Vector{1, -1, 1}}, Point{r3.Vector{1, 1, 1}}).Lat.Hi,
    85  		},
    86  		// Min latitude, CCW edge
    87  		{
    88  			-math.Pi/4 - 0.5*rectErrorLat,
    89  			rectBoundForPoints(Point{r3.Vector{-1, 1, -1}}, Point{r3.Vector{-1, -1, -1}}).Lat.Lo,
    90  		},
    91  
    92  		// Check cases where the edge passes through one of the poles.
    93  		{
    94  			math.Pi / 2,
    95  			rectBoundForPoints(Point{r3.Vector{.3, .4, 1}}, Point{r3.Vector{-.3, -.4, 1}}).Lat.Hi,
    96  		},
    97  		{
    98  			-math.Pi / 2,
    99  			rectBoundForPoints(Point{r3.Vector{.3, .4, -1}}, Point{r3.Vector{-.3, -.4, -1}}).Lat.Lo,
   100  		},
   101  	}
   102  
   103  	for _, test := range tests {
   104  		if !float64Eq(test.got, test.want) {
   105  			t.Errorf("RectBound for max lat on interior of edge failed; got %v want %v", test.got, test.want)
   106  		}
   107  	}
   108  }
   109  
   110  func TestRectBounderMaxLatitudeRandom(t *testing.T) {
   111  	// Check that the maximum latitude of edges is computed accurately to within
   112  	// 3 * dblEpsilon (the expected maximum error). We concentrate on maximum
   113  	// latitudes near the equator and north pole since these are the extremes.
   114  
   115  	for i := 0; i < 100; i++ {
   116  		// Construct a right-handed coordinate frame (U,V,W) such that U points
   117  		// slightly above the equator, V points at the equator, and W is slightly
   118  		// offset from the north pole.
   119  		u := randomPoint()
   120  		u.Z = dblEpsilon * 1e-6 * math.Pow(1e12, randomFloat64())
   121  
   122  		u = Point{u.Normalize()}
   123  		v := Point{PointFromCoords(0, 0, 1).PointCross(u).Normalize()}
   124  		w := Point{u.PointCross(v).Normalize()}
   125  
   126  		// Construct a line segment AB that passes through U, and check that the
   127  		// maximum latitude of this segment matches the latitude of U.
   128  		a := Point{u.Sub(v.Mul(randomFloat64())).Normalize()}
   129  		b := Point{u.Add(v.Mul(randomFloat64())).Normalize()}
   130  		abBound := rectBoundForPoints(a, b)
   131  		if !float64Near(latitude(u).Radians(), abBound.Lat.Hi, rectErrorLat) {
   132  			t.Errorf("bound for line AB not near enough to the latitude of point %v. got %v, want %v",
   133  				u, latitude(u).Radians(), abBound.Lat.Hi)
   134  		}
   135  
   136  		// Construct a line segment CD that passes through W, and check that the
   137  		// maximum latitude of this segment matches the latitude of W.
   138  		c := Point{w.Sub(v.Mul(randomFloat64())).Normalize()}
   139  		d := Point{w.Add(v.Mul(randomFloat64())).Normalize()}
   140  		cdBound := rectBoundForPoints(c, d)
   141  		if !float64Near(latitude(w).Radians(), cdBound.Lat.Hi, rectErrorLat) {
   142  			t.Errorf("bound for line CD not near enough to the lat of point %v. got %v, want %v",
   143  				v, latitude(w).Radians(), cdBound.Lat.Hi)
   144  		}
   145  	}
   146  }
   147  
   148  func TestRectBounderExpandForSubregions(t *testing.T) {
   149  	// Test the full and empty bounds.
   150  	if !ExpandForSubregions(FullRect()).IsFull() {
   151  		t.Errorf("Subregion Bound of full rect should be full")
   152  	}
   153  	if !ExpandForSubregions(EmptyRect()).IsEmpty() {
   154  		t.Errorf("Subregion Bound of empty rect should be empty")
   155  	}
   156  
   157  	tests := []struct {
   158  		xLat, xLng, yLat, yLng float64
   159  		wantFull               bool
   160  	}{
   161  		// Cases where the bound does not straddle the equator (but almost does),
   162  		// and spans nearly 180 degrees in longitude.
   163  		{3e-16, 0, 1e-14, math.Pi, true},
   164  		{9e-16, 0, 1e-14, math.Pi, false},
   165  		{1e-16, 7e-16, 1e-14, math.Pi, true},
   166  		{3e-16, 14e-16, 1e-14, math.Pi, false},
   167  		{1e-100, 14e-16, 1e-14, math.Pi, true},
   168  		{1e-100, 22e-16, 1e-14, math.Pi, false},
   169  		// Cases where the bound spans at most 90 degrees in longitude, and almost
   170  		// 180 degrees in latitude.  Note that DBL_EPSILON is about 2.22e-16, which
   171  		// implies that the double-precision value just below Pi/2 can be written as
   172  		// (math.Pi/2 - 2e-16).
   173  		{-math.Pi / 2, -1e-15, math.Pi/2 - 7e-16, 0, true},
   174  		{-math.Pi / 2, -1e-15, math.Pi/2 - 30e-16, 0, false},
   175  		{-math.Pi/2 + 4e-16, 0, math.Pi/2 - 2e-16, 1e-7, true},
   176  		{-math.Pi/2 + 30e-16, 0, math.Pi / 2, 1e-7, false},
   177  		{-math.Pi/2 + 4e-16, 0, math.Pi/2 - 4e-16, math.Pi / 2, true},
   178  		{-math.Pi / 2, 0, math.Pi/2 - 30e-16, math.Pi / 2, false},
   179  		// Cases where the bound straddles the equator and spans more than 90
   180  		// degrees in longitude.  These are the cases where the critical distance is
   181  		// between a corner of the bound and the opposite longitudinal edge.  Unlike
   182  		// the cases above, here the bound may contain nearly-antipodal points (to
   183  		// within 3.055 * DBL_EPSILON) even though the latitude and longitude ranges
   184  		// are both significantly less than (math.Pi - 3.055 * DBL_EPSILON).
   185  		{-math.Pi / 2, 0, math.Pi/2 - 1e-8, math.Pi - 1e-7, true},
   186  		{-math.Pi / 2, 0, math.Pi/2 - 1e-7, math.Pi - 1e-7, false},
   187  		{-math.Pi/2 + 1e-12, -math.Pi + 1e-4, math.Pi / 2, 0, true},
   188  		{-math.Pi/2 + 1e-11, -math.Pi + 1e-4, math.Pi / 2, 0, true},
   189  	}
   190  
   191  	for _, tc := range tests {
   192  		in := RectFromLatLng(LatLng{s1.Angle(tc.xLat), s1.Angle(tc.xLng)})
   193  		in = in.AddPoint(LatLng{s1.Angle(tc.yLat), s1.Angle(tc.yLng)})
   194  		got := ExpandForSubregions(in)
   195  
   196  		// Test that the bound is actually expanded.
   197  		if !got.Contains(in) {
   198  			t.Errorf("Subregion bound of (%f, %f, %f, %f) should contain original rect", tc.xLat, tc.xLng, tc.yLat, tc.yLng)
   199  		}
   200  		if in.Lat == validRectLatRange && in.Lat.ContainsInterval(got.Lat) {
   201  			t.Errorf("Subregion bound of (%f, %f, %f, %f) shouldn't be contained by original rect", tc.xLat, tc.xLng, tc.yLat, tc.yLng)
   202  		}
   203  
   204  		// We check the various situations where the bound contains nearly-antipodal points. The tests are organized into pairs
   205  		// where the two bounds are similar except that the first bound meets the nearly-antipodal criteria while the second does not.
   206  		if got.IsFull() != tc.wantFull {
   207  			t.Errorf("Subregion Bound of (%f, %f, %f, %f).IsFull should be %t", tc.xLat, tc.xLng, tc.yLat, tc.yLng, tc.wantFull)
   208  		}
   209  	}
   210  
   211  	rectTests := []struct {
   212  		xLat, xLng, yLat, yLng float64
   213  		wantRect               Rect
   214  	}{
   215  		{1.5, -math.Pi / 2, 1.5, math.Pi/2 - 2e-16, Rect{r1.Interval{1.5, 1.5}, s1.FullInterval()}},
   216  		{1.5, -math.Pi / 2, 1.5, math.Pi/2 - 7e-16, Rect{r1.Interval{1.5, 1.5}, s1.Interval{-math.Pi / 2, math.Pi/2 - 7e-16}}},
   217  		// Check for cases where the bound is expanded to include one of the poles
   218  		{-math.Pi/2 + 1e-15, 0, -math.Pi/2 + 1e-15, 0, Rect{r1.Interval{-math.Pi / 2, -math.Pi/2 + 1e-15}, s1.FullInterval()}},
   219  		{math.Pi/2 - 1e-15, 0, math.Pi/2 - 1e-15, 0, Rect{r1.Interval{math.Pi/2 - 1e-15, math.Pi / 2}, s1.FullInterval()}},
   220  	}
   221  
   222  	for _, tc := range rectTests {
   223  		// Now we test cases where the bound does not contain nearly-antipodal
   224  		// points, but it does contain points that are approximately 180 degrees
   225  		// apart in latitude.
   226  		in := RectFromLatLng(LatLng{s1.Angle(tc.xLat), s1.Angle(tc.xLng)})
   227  		in = in.AddPoint(LatLng{s1.Angle(tc.yLat), s1.Angle(tc.yLng)})
   228  		got := ExpandForSubregions(in)
   229  		if !rectsApproxEqual(got, tc.wantRect, rectErrorLat, rectErrorLng) {
   230  			t.Errorf("Subregion Bound of (%f, %f, %f, %f) = (%v) should be %v", tc.xLat, tc.xLng, tc.yLat, tc.yLng, got, tc.wantRect)
   231  		}
   232  	}
   233  }
   234  
   235  // TODO(roberts): TestRectBounderNearlyIdenticalOrAntipodalPoints(t *testing.T)
   236
View as plain text