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
16 // Author: ericv@google.com (Eric Veach)
17
18 #include "s2/s2edge_tessellator.h"
19
20 #include <iostream>
21 #include <gtest/gtest.h>
22 #include "s2/third_party/absl/strings/str_cat.h"
23 #include "s2/s2edge_distances.h"
24 #include "s2/s2loop.h"
25 #include "s2/s2pointutil.h"
26 #include "s2/s2projections.h"
27 #include "s2/s2testing.h"
28 #include "s2/s2text_format.h"
29
30 using absl::StrCat;
31 using s2textformat::ParsePointsOrDie;
32 using std::cout;
33 using std::endl;
34 using std::fabs;
35 using std::min;
36 using std::max;
37 using std::vector;
38
39 namespace {
40
41 class Stats {
42 public:
Stats()43 Stats() : max_(-std::numeric_limits<double>::infinity()),
44 sum_(0), count_(0) {
45 }
46
Tally(double v)47 void Tally(double v) {
48 if (std::isnan(v)) S2_LOG(FATAL) << "NaN";
49 max_ = std::max(v, max_);
50 sum_ += v;
51 count_ += 1;
52 }
53
max() const54 double max() const { return max_; }
avg() const55 double avg() const { return sum_ / count_; }
56
ToString() const57 string ToString() const {
58 return StrCat("avg = ", sum_ / count_, ", max = ", max_);
59 }
60
61 private:
62 double max_, sum_;
63 int count_;
64 };
65
66 // Determines whether the distance between the two edges is measured
67 // geometrically or parameterically (see algorithm description in .cc file).
68 enum class DistType { GEOMETRIC, PARAMETRIC };
69
GetMaxDistance(const S2::Projection & proj,const R2Point & px,const S2Point & x,const R2Point & py,const S2Point & y,DistType dist_type=DistType::GEOMETRIC)70 S1Angle GetMaxDistance(const S2::Projection& proj,
71 const R2Point& px, const S2Point& x,
72 const R2Point& py, const S2Point& y,
73 DistType dist_type = DistType::GEOMETRIC) {
74 // Step along the projected edge at a fine resolution and keep track of the
75 // maximum distance of any point to the current geodesic edge.
76 const int kNumSteps = 100;
77 S1ChordAngle max_dist = S1ChordAngle::Zero();
78 for (double f = 0.5 / kNumSteps; f < 1.0; f += 1.0 / kNumSteps) {
79 S1ChordAngle dist = S1ChordAngle::Infinity();
80 S2Point p = proj.Unproject(proj.Interpolate(f, px, py));
81 if (dist_type == DistType::GEOMETRIC) {
82 S2::UpdateMinDistance(p, x, y, &dist);
83 } else {
84 S2_DCHECK(dist_type == DistType::PARAMETRIC);
85 dist = S1ChordAngle(p, S2::Interpolate(f, x, y));
86 }
87 if (dist > max_dist) max_dist = dist;
88 }
89 // Ensure that the maximum distance estimate is a lower bound, not an upper
90 // bound, since we only want to record a failure of the distance estimation
91 // algorithm if the number it returns is definitely too small.
92 return S1Angle(
93 max_dist.PlusError(-S2::GetUpdateMinDistanceMaxError(max_dist)));
94 }
95
96 // When there are longitudes greater than 180 degrees due to wrapping, the
97 // combination of projecting and unprojecting an S2Point can have slightly more
98 // error than is allowed by S2::ApproxEquals.
99 const S1Angle kMaxProjError(S1Angle::Radians(2e-15));
100
101 // Converts a projected edge to a sequence of geodesic edges and verifies that
102 // the result satisfies the given tolerance.
TestUnprojected(const S2::Projection & proj,S1Angle tolerance,const R2Point & pa,const R2Point & pb_in,bool log_stats)103 Stats TestUnprojected(const S2::Projection& proj, S1Angle tolerance,
104 const R2Point& pa, const R2Point& pb_in, bool log_stats) {
105 S2EdgeTessellator tess(&proj, tolerance);
106 vector<S2Point> vertices;
107 tess.AppendUnprojected(pa, pb_in, &vertices);
108 R2Point pb = proj.WrapDestination(pa, pb_in);
109 EXPECT_LE(S1Angle(proj.Unproject(pa), vertices.front()), kMaxProjError);
110 EXPECT_LE(S1Angle(proj.Unproject(pb), vertices.back()), kMaxProjError);
111 Stats stats;
112 if (pa == pb) {
113 EXPECT_EQ(1, vertices.size());
114 return stats;
115 }
116 // Precompute the normal to the projected edge.
117 Vector2_d norm = (pb - pa).Ortho().Normalize();
118 S2Point x = vertices[0];
119 R2Point px = proj.Project(x);
120 for (int i = 1; i < vertices.size(); ++i) {
121 S2Point y = vertices[i];
122 R2Point py = proj.WrapDestination(px, proj.Project(y));
123 // Check that every vertex is on the projected edge.
124 EXPECT_LE((py - pa).DotProd(norm), 1e-14 * py.Norm());
125 stats.Tally(GetMaxDistance(proj, px, x, py, y) / tolerance);
126 x = y;
127 px = py;
128 }
129 if (log_stats) {
130 cout << pa << " to " << pb << ": " << vertices.size() << " vertices, "
131 << stats.ToString() << endl;
132 }
133 return stats;
134 }
135
136 // Converts a geodesic edge to a sequence of projected edges and verifies that
137 // the result satisfies the given tolerance.
TestProjected(const S2::Projection & proj,S1Angle tolerance,const S2Point & a,const S2Point & b,bool log_stats)138 Stats TestProjected(const S2::Projection& proj, S1Angle tolerance,
139 const S2Point& a, const S2Point& b, bool log_stats) {
140 S2EdgeTessellator tess(&proj, tolerance);
141 vector<R2Point> vertices;
142 tess.AppendProjected(a, b, &vertices);
143 EXPECT_LE(S1Angle(a, proj.Unproject(vertices.front())), kMaxProjError);
144 EXPECT_LE(S1Angle(b, proj.Unproject(vertices.back())), kMaxProjError);
145 Stats stats;
146 if (a == b) {
147 EXPECT_EQ(1, vertices.size());
148 return stats;
149 }
150 R2Point px = vertices[0];
151 S2Point x = proj.Unproject(px);
152 for (int i = 1; i < vertices.size(); ++i) {
153 R2Point py = vertices[i];
154 S2Point y = proj.Unproject(py);
155 // Check that every vertex is on the geodesic edge.
156 static S1ChordAngle kMaxInterpolationError(S1Angle::Radians(1e-14));
157 EXPECT_TRUE(S2::IsDistanceLess(y, a, b, kMaxInterpolationError));
158 stats.Tally(GetMaxDistance(proj, px, x, py, y) / tolerance);
159 x = y;
160 px = py;
161 }
162 if (log_stats) {
163 cout << vertices[0] << " to " << px << ": " << vertices.size()
164 << " vertices, " << stats.ToString() << endl;
165 }
166 return stats;
167 }
168
TEST(S2EdgeTessellator,ProjectedNoTessellation)169 TEST(S2EdgeTessellator, ProjectedNoTessellation) {
170 S2::PlateCarreeProjection proj(180);
171 S2EdgeTessellator tess(&proj, S1Angle::Degrees(0.01));
172 vector<R2Point> vertices;
173 tess.AppendProjected(S2Point(1, 0, 0), S2Point(0, 1, 0), &vertices);
174 EXPECT_EQ(2, vertices.size());
175 }
176
TEST(S2EdgeTessellator,UnprojectedNoTessellation)177 TEST(S2EdgeTessellator, UnprojectedNoTessellation) {
178 S2::PlateCarreeProjection proj(180);
179 S2EdgeTessellator tess(&proj, S1Angle::Degrees(0.01));
180 vector<S2Point> vertices;
181 tess.AppendUnprojected(R2Point(0, 30), R2Point(0, 50), &vertices);
182 EXPECT_EQ(2, vertices.size());
183 }
184
TEST(S2EdgeTessellator,UnprojectedWrapping)185 TEST(S2EdgeTessellator, UnprojectedWrapping) {
186 // This tests that a projected edge that crosses the 180 degree meridian
187 // goes the "short way" around the sphere.
188
189 S2::PlateCarreeProjection proj(180);
190 S2EdgeTessellator tess(&proj, S1Angle::Degrees(0.01));
191 vector<S2Point> vertices;
192 tess.AppendUnprojected(R2Point(-170, 0), R2Point(170, 80), &vertices);
193 for (const auto& v : vertices) {
194 EXPECT_GE(fabs(S2LatLng::Longitude(v).degrees()), 170);
195 }
196 }
197
TEST(S2EdgeTessellator,ProjectedWrapping)198 TEST(S2EdgeTessellator, ProjectedWrapping) {
199 // This tests projecting a geodesic edge that crosses the 180 degree
200 // meridian. This results in a set of vertices that may be non-canonical
201 // (i.e., absolute longitudes greater than 180 degrees) but that don't have
202 // any sudden jumps in value, which is convenient for interpolating them.
203 S2::PlateCarreeProjection proj(180);
204 S2EdgeTessellator tess(&proj, S1Angle::Degrees(0.01));
205 vector<R2Point> vertices;
206 tess.AppendProjected(S2LatLng::FromDegrees(0, -170).ToPoint(),
207 S2LatLng::FromDegrees(0, 170).ToPoint(), &vertices);
208 for (const auto& v : vertices) {
209 EXPECT_LE(v.x(), -170);
210 }
211 }
212
TEST(S2EdgeTessellator,UnprojectedWrappingMultipleCrossings)213 TEST(S2EdgeTessellator, UnprojectedWrappingMultipleCrossings) {
214 // Tests an edge chain that crosses the 180 degree meridian multiple times.
215 // Note that due to coordinate wrapping, the last vertex of one edge may not
216 // exactly match the first edge of the next edge after unprojection.
217 S2::PlateCarreeProjection proj(180);
218 S2EdgeTessellator tess(&proj, S1Angle::Degrees(0.01));
219 vector<S2Point> vertices;
220 for (double lat = 1; lat <= 60; ++lat) {
221 tess.AppendUnprojected(R2Point(180 - 0.03 * lat, lat),
222 R2Point(-180 + 0.07 * lat, lat), &vertices);
223 tess.AppendUnprojected(R2Point(-180 + 0.07 * lat, lat),
224 R2Point(180 - 0.03 * (lat + 1), lat + 1), &vertices);
225 }
226 for (const auto& v : vertices) {
227 EXPECT_GE(fabs(S2LatLng::Longitude(v).degrees()), 175);
228 }
229 }
230
TEST(S2EdgeTessellator,ProjectedWrappingMultipleCrossings)231 TEST(S2EdgeTessellator, ProjectedWrappingMultipleCrossings) {
232 // The following loop crosses the 180 degree meridian four times (twice in
233 // each direction).
234 auto loop = ParsePointsOrDie("0:160, 0:-40, 0:120, 0:-80, 10:120, "
235 "10:-40, 0:160");
236 S2::PlateCarreeProjection proj(180);
237 S1Angle tolerance(S1Angle::E7(1));
238 S2EdgeTessellator tess(&proj, tolerance);
239 vector<R2Point> vertices;
240 for (int i = 0; i + 1 < loop.size(); ++i) {
241 tess.AppendProjected(loop[i], loop[i + 1], &vertices);
242 }
243 EXPECT_EQ(vertices.front(), vertices.back());
244
245 // Note that the R2Point coordinates are in (lng, lat) order.
246 double min_lng = vertices[0].x();
247 double max_lng = vertices[0].x();
248 for (const R2Point& v : vertices) {
249 min_lng = min(min_lng, v.x());
250 max_lng = max(max_lng, v.x());
251 }
252 EXPECT_EQ(160, min_lng);
253 EXPECT_EQ(640, max_lng);
254 }
255
TEST(S2EdgeTessellator,InfiniteRecursionBug)256 TEST(S2EdgeTessellator, InfiniteRecursionBug) {
257 S2::PlateCarreeProjection proj(180);
258 S1Angle kOneMicron = S1Angle::Radians(1e-6 / 6371.0);
259 S2EdgeTessellator tess(&proj, kOneMicron);
260 vector<R2Point> vertices;
261 tess.AppendProjected(S2LatLng::FromDegrees(3, 21).ToPoint(),
262 S2LatLng::FromDegrees(1, -159).ToPoint(), &vertices);
263 EXPECT_EQ(36, vertices.size());
264 }
265
TEST(S2EdgeTessellator,UnprojectedAccuracy)266 TEST(S2EdgeTessellator, UnprojectedAccuracy) {
267 S2::MercatorProjection proj(180);
268 S1Angle tolerance(S1Angle::Degrees(1e-5));
269 R2Point pa(0, 0), pb(89.999999, 179);
270 Stats stats = TestUnprojected(proj, tolerance, pa, pb, true);
271 EXPECT_LE(stats.max(), 1.0);
272 }
273
TEST(S2EdgeTessellator,UnprojectedAccuracyCrossEquator)274 TEST(S2EdgeTessellator, UnprojectedAccuracyCrossEquator) {
275 S2::MercatorProjection proj(180);
276 S1Angle tolerance(S1Angle::Degrees(1e-5));
277 R2Point pa(-10, -10), pb(10, 10);
278 Stats stats = TestUnprojected(proj, tolerance, pa, pb, true);
279 EXPECT_LT(stats.max(), 1.0);
280 }
281
TEST(S2EdgeTessellator,ProjectedAccuracy)282 TEST(S2EdgeTessellator, ProjectedAccuracy) {
283 S2::PlateCarreeProjection proj(180);
284 S1Angle tolerance(S1Angle::E7(1));
285 S2Point a = S2LatLng::FromDegrees(-89.999, -170).ToPoint();
286 S2Point b = S2LatLng::FromDegrees(50, 100).ToPoint();
287 Stats stats = TestProjected(proj, tolerance, a, b, true);
288 EXPECT_LE(stats.max(), 1.0);
289 }
290
TEST(S2EdgeTessellator,UnprojectedAccuracyMidpointEquator)291 TEST(S2EdgeTessellator, UnprojectedAccuracyMidpointEquator) {
292 S2::PlateCarreeProjection proj(180);
293 S1Angle tolerance = S2Testing::MetersToAngle(1);
294 R2Point a(80, 50), b(-80, -50);
295 Stats stats = TestUnprojected(proj, tolerance, a, b, true);
296 EXPECT_LE(stats.max(), 1.0);
297 }
298
TEST(S2EdgeTessellator,ProjectedAccuracyMidpointEquator)299 TEST(S2EdgeTessellator, ProjectedAccuracyMidpointEquator) {
300 S2::PlateCarreeProjection proj(180);
301 S1Angle tolerance = S2Testing::MetersToAngle(1);
302 S2Point a = S2LatLng::FromDegrees(50, 80).ToPoint();
303 S2Point b = S2LatLng::FromDegrees(-50, -80).ToPoint();
304 Stats stats = TestProjected(proj, tolerance, a, b, true);
305 EXPECT_LE(stats.max(), 1.0);
306 }
307
TEST(S2EdgeTessellator,ProjectedAccuracyCrossEquator)308 TEST(S2EdgeTessellator, ProjectedAccuracyCrossEquator) {
309 S2::PlateCarreeProjection proj(180);
310 S1Angle tolerance(S1Angle::E7(1));
311 S2Point a = S2LatLng::FromDegrees(-20, -20).ToPoint();
312 S2Point b = S2LatLng::FromDegrees(20, 20).ToPoint();
313 Stats stats = TestProjected(proj, tolerance, a, b, true);
314 EXPECT_LT(stats.max(), 1.0);
315 }
316
TEST(S2EdgeTessellator,ProjectedAccuracySeattleToNewYork)317 TEST(S2EdgeTessellator, ProjectedAccuracySeattleToNewYork) {
318 S2::PlateCarreeProjection proj(180);
319 S1Angle tolerance = S2Testing::MetersToAngle(1);
320 S2Point seattle(S2LatLng::FromDegrees(47.6062, -122.3321).ToPoint());
321 S2Point newyork(S2LatLng::FromDegrees(40.7128, -74.0059).ToPoint());
322 Stats stats = TestProjected(proj, tolerance, seattle, newyork, true);
323 EXPECT_LE(stats.max(), 1.0);
324 }
325
326 // Given a projection, this function repeatedly chooses a pair of edge
327 // endpoints and measures the true distance between the geodesic and projected
328 // edges that connect those two endpoints. It then compares this to against
329 // the distance measurement algorithm used by S2EdgeTessellator, which
330 // consists of measuring the point-to-point distance between the edges at each
331 // of two fractions "t" and "1-t", computing the maximum of those two
332 // distances, and then scaling by the constant documented in the .cc file
333 // (based on the idea that the distance between the edges as a function of the
334 // interpolation fraction can be accurately modeled as a cubic polynomial).
335 //
336 // This function is used to (1) verify that the distance estimates are always
337 // conservative, (2) verify the optimality of the interpolation fraction "t",
338 // and (3) estimate the amount of overtessellation that occurs for various
339 // types of edges (e.g., short vs. long edges, edges that follow lines of
340 // latitude or longitude, etc).
TestEdgeError(const S2::Projection & proj,double t)341 void TestEdgeError(const S2::Projection& proj, double t) {
342 // Here we compute how much we need to scale the error measured at the
343 // chosen interpolation fraction "t" in order to bound the error along the
344 // entire edge, under the assumption that the error is a convex combination
345 // of E1(x) and E2(x) (see comments in the .cc file).
346 const double x = 1 - 2 * t;
347 const double dlat = sin(0.5 * M_PI_4 * (1 - x));
348 const double dlng = sin(M_PI_4 * (1 - x));
349 const double dsin2 = dlat * dlat + dlng * dlng * sin(M_PI_4 * x) * M_SQRT1_2;
350 const double dsin2_max = 0.5 * (1 - M_SQRT1_2);
351 // Note that this is the reciprocal of the value used in the .cc file!
352 const double kScaleFactor = max((2 * sqrt(3) / 9) / (x * (1 - x * x)),
353 asin(sqrt(dsin2_max)) / asin(sqrt(dsin2)));
354
355 // Keep track of the average and maximum geometric and parametric errors.
356 Stats stats_g, stats_p;
357 const int kIters = google::DEBUG_MODE ? 10000 : 100000;
358 for (int iter = 0; iter < kIters; ++iter) {
359 S2Testing::rnd.Reset(iter);
360 S2Point a = S2Testing::RandomPoint();
361 S2Point b = S2Testing::RandomPoint();
362 // Uncomment to test edges longer than 90 degrees.
363 if (a.DotProd(b) < -1e-14) continue;
364 // Uncomment to test edges than span more than 90 degrees longitude.
365 // if (a[0] * b[0] + a[1] * b[1] < 0) continue;
366 // Uncomment to only test edges of a certain length.
367 // b = S2::InterpolateAtDistance(S1Angle::Radians(1e-5), a, b);
368 // Uncomment to only test edges that stay in one hemisphere.
369 // if (a[2] * b[2] <= 0) continue;
370 R2Point pa = proj.Project(a);
371 R2Point pb = proj.WrapDestination(pa, proj.Project(b));
372 S1Angle max_dist_g = GetMaxDistance(proj, pa, a, pb, b,
373 DistType::GEOMETRIC);
374 // Ignore edges where the error is too small.
375 if (max_dist_g <= S2EdgeTessellator::kMinTolerance()) continue;
376 S1Angle max_dist_p = GetMaxDistance(proj, pa, a, pb, b,
377 DistType::PARAMETRIC);
378 if (max_dist_p <= S2EdgeTessellator::kMinTolerance()) continue;
379
380 // Compute the estimated error bound.
381 S1Angle d1(S2::Interpolate(t, a, b), proj.Unproject((1-t) * pa + t * pb));
382 S1Angle d2(S2::Interpolate(1-t, a, b), proj.Unproject(t * pa + (1-t) * pb));
383 S1Angle dist = kScaleFactor * max(S1Angle::Radians(1e-300), max(d1, d2));
384
385 // Compute the ratio of the true geometric/parametric errors to the
386 // estimate error bound.
387 double r_g = max_dist_g / dist;
388 double r_p = max_dist_p / dist;
389
390 // Our objective is to ensure that the geometric error ratio is at most 1.
391 // (The parametric ratio is computed only for analysis purposes.)
392 if (r_g > 0.99999) {
393 // Log any edges where the ratio is exceeded (or nearly so).
394 cout << pa << " to " << pb << ": ratio = " << r_g
395 << ", dist = " << max_dist_g.degrees() << endl;
396 }
397 stats_g.Tally(r_g);
398 stats_p.Tally(r_p);
399 }
400 cout << "t = " << t << ", scale = " << kScaleFactor << ", G["
401 << stats_g.ToString() << "], P[" << stats_p.ToString() << "]" << endl;
402 EXPECT_LE(stats_g.max(), kScaleFactor);
403 }
404
405 // The interpolation parameter actually used in the .cc file.
406 static constexpr double kBestFraction = 0.31215691082248312;
407
TEST(S2EdgeTessellator,MaxEdgeErrorPlateCarree)408 TEST(S2EdgeTessellator, MaxEdgeErrorPlateCarree) {
409 S2::PlateCarreeProjection proj(180);
410 // Uncomment to test some nearby parameter values.
411 // TestEdgeError(proj, 0.311);
412 TestEdgeError(proj, kBestFraction);
413 // TestEdgeError(proj, 0.313);
414 }
415
TEST(S2EdgeTessellator,MaxEdgeErrorMercator)416 TEST(S2EdgeTessellator, MaxEdgeErrorMercator) {
417 S2::MercatorProjection proj(180);
418 // Uncomment to test some nearby parameter values.
419 // TestEdgeError(proj, 0.311);
420 TestEdgeError(proj, kBestFraction);
421 // TestEdgeError(proj, 0.313);
422 }
423
424 // Tessellates random edges using the given projection and tolerance, and
425 // verifies that the expected criteria are satisfied.
TestRandomEdges(const S2::Projection & proj,S1Angle tolerance)426 void TestRandomEdges(const S2::Projection& proj, S1Angle tolerance) {
427 const int kIters = google::DEBUG_MODE ? 50 : 500;
428 double max_r2 = 0, max_s2 = 0;
429 for (int iter = 0; iter < kIters; ++iter) {
430 S2Testing::rnd.Reset(iter);
431 S2Point a = S2Testing::RandomPoint();
432 S2Point b = S2Testing::RandomPoint();
433 max_r2 = max(max_r2, TestProjected(proj, tolerance, a, b, false).max());
434 R2Point pa = proj.Project(a);
435 R2Point pb = proj.Project(b);
436 max_s2 = max(max_s2, TestUnprojected(proj, tolerance, pa, pb, false).max());
437 }
438 cout << "max_r2 = " << max_r2 << ", max_s2 = " << max_s2 << endl;
439 EXPECT_LE(max_r2, 1.0);
440 EXPECT_LE(max_s2, 1.0);
441 }
442
TEST(S2EdgeTessellator,RandomEdgesPlateCarree)443 TEST(S2EdgeTessellator, RandomEdgesPlateCarree) {
444 S2::PlateCarreeProjection proj(180);
445 S1Angle tolerance = S2Testing::MetersToAngle(100);
446 TestRandomEdges(proj, tolerance);
447 }
448
TEST(S2EdgeTessellator,RandomEdgesMercator)449 TEST(S2EdgeTessellator, RandomEdgesMercator) {
450 S2::MercatorProjection proj(180);
451 S1Angle tolerance = S2Testing::MetersToAngle(100);
452 TestRandomEdges(proj, tolerance);
453 }
454
455 // TODO(ericv): Superceded by random edge tests above, remove?
TEST(S2EdgeTessellator,UnprojectedAccuracyRandomCheck)456 TEST(S2EdgeTessellator, UnprojectedAccuracyRandomCheck) {
457 S2::PlateCarreeProjection proj(180);
458 S1Angle tolerance(S1Angle::Degrees(1e-3));
459 S2Testing::Random rand;
460 const int kIters = google::DEBUG_MODE ? 250 : 5000;
461 for (int i = 0; i < kIters; ++i) {
462 S2Testing::rnd.Reset(i);
463 double alat = rand.UniformDouble(-89.99, 89.99);
464 double blat = rand.UniformDouble(-89.99, 89.99);
465 double blon = rand.UniformDouble(0, 179);
466
467 R2Point pa(0, alat), pb(blon, blat);
468 Stats stats = TestUnprojected(proj, tolerance, pa, pb, false);
469 EXPECT_LT(stats.max(), 1.0) << pa << ", " << pb;
470 }
471 }
472
473 // XXX(ericv): Superceded by random edge tests above, remove?
TEST(S2EdgeTessellator,ProjectedAccuracyRandomCheck)474 TEST(S2EdgeTessellator, ProjectedAccuracyRandomCheck) {
475 S2::PlateCarreeProjection proj(180);
476 S1Angle tolerance(S1Angle::Degrees(1e-3));
477 S2Testing::Random rand;
478 const int kIters = google::DEBUG_MODE ? 250 : 5000;
479 for (int i = 0; i < kIters; ++i) {
480 S2Testing::rnd.Reset(i);
481 double alat = rand.UniformDouble(-89.99, 89.99);
482 double blat = rand.UniformDouble(-89.99, 89.99);
483 double blon = rand.UniformDouble(-180, 180);
484
485 S2Point a = S2LatLng::FromDegrees(alat, 0).ToPoint();
486 S2Point b = S2LatLng::FromDegrees(blat, blon).ToPoint();
487 Stats stats = TestProjected(proj, tolerance, a, b, false);
488 EXPECT_LT(stats.max(), 1.0) << a << ", " << b;
489 }
490 }
491
492 } // namespace
493