diff --git a/lib/great_circle.dart b/lib/great_circle.dart new file mode 100644 index 00000000..4d82af0b --- /dev/null +++ b/lib/great_circle.dart @@ -0,0 +1,3 @@ +library turf_great_circle; + +export 'src/great_circle.dart'; \ No newline at end of file diff --git a/lib/src/great_circle.dart b/lib/src/great_circle.dart new file mode 100644 index 00000000..fc9ae50a --- /dev/null +++ b/lib/src/great_circle.dart @@ -0,0 +1,164 @@ +import 'dart:math' as math; +import 'package:turf/turf.dart'; + +/// Calculates the great circle route between two points on a sphere +/// +/// Useful link: https://en.wikipedia.org/wiki/Great-circle_distance + +Feature greatCircle(Position start, Position end, + {Map properties = const {}, + int npoints = + 100, // npoints = number of intermediate points less one (e.g if you want 5 intermediate points, set npoints = 6) + int offset = 10}) { + if (start.length != 2 || end.length != 2) { + /// Coordinate checking + throw ArgumentError( + "Both start and end coordinates should have two values - a latitude and longitude"); + } + + // If start and end points are the same, + if (start[0] == end[0] && start[1] == end[1]) { + return Feature(geometry: LineString(coordinates: [])); + } + // Coordinate checking for valid values + if (start[0]! < -90) { + throw ArgumentError( + "Starting latitude (vertical) coordinate is less than -90. This is not a valid coordinate."); + } + + if (start[0]! > 90) { + throw ArgumentError( + "Starting latitude (vertical) coordinate is greater than 90. This is not a valid coordinate."); + } + + if (start[1]! < -180) { + throw ArgumentError( + 'Starting longitude (horizontal) coordinate is less than -180. This is not a valid coordinate.'); + } + + if (start[1]! > 180) { + throw ArgumentError( + 'Starting longitude (horizontal) coordinate is greater than 180. This is not a valid coordinate.'); + } + + if (end[0]! < -90) { + throw ArgumentError( + "Ending latitude (vertical) coordinate is less than -90. This is not a valid coordinate."); + } + + if (end[0]! > 90) { + throw ArgumentError( + "Ending latitude (vertical) coordinate is greater than 90. This is not a valid coordinate."); + } + + if (end[1]! < -180) { + throw ArgumentError( + 'Ending longitude (horizontal) coordinate is less than -180. This is not a valid coordinate.'); + } + + if (end[1]! > 180) { + throw ArgumentError( + 'Ending longitude (horizontal) coordinate is greater than 180. This is not a valid coordinate.'); + } + + // Creates a list to store points for the great circle arc + List line = []; + + num lat1 = degreesToRadians(start[0]!); + num lng1 = degreesToRadians(start[1]!); + num lat2 = degreesToRadians(end[0]!); + num lng2 = degreesToRadians(end[1]!); + + // Harvesine formula + for (int i = 0; i <= npoints; i++) { + double f = i / npoints; + double delta = 2 * + math.asin(math.sqrt(math.pow(math.sin((lat2 - lat1) / 2), 2) + + math.cos(lat1) * + math.cos(lat2) * + math.pow(math.sin((lng2 - lng1) / 2), 2))); + double A = math.sin((1 - f) * delta) / math.sin(delta); + double B = math.sin(f * delta) / math.sin(delta); + double x = A * math.cos(lat1) * math.cos(lng1) + + B * math.cos(lat2) * math.cos(lng2); + double y = A * math.cos(lat1) * math.sin(lng1) + + B * math.cos(lat2) * math.sin(lng2); + double z = A * math.sin(lat1) + B * math.sin(lat2); + + double lat = math.atan2(z, math.sqrt(x * x + y * y)); + double lng = math.atan2(y, x); + + Position point = Position(lng, lat); + line.add(point); + } + + /// Check for multilinestring if path crosses anti-meridian + bool crossAntiMeridian = (lng1 - lng2).abs() > 180; + + /// If it crossed antimeridian, we need to split our lines + if (crossAntiMeridian) { + List> multiLine = []; + List currentLine = []; + + for (var point in line) { + if ((point[0]! - line[0][0]!).abs() > 180) { + multiLine.addAll([currentLine]); + currentLine = []; + } + currentLine.add(point); + } + multiLine.addAll([currentLine]); + return Feature( + geometry: MultiLineString(coordinates: multiLine)); + } + return Feature(geometry: LineString(coordinates: line)); +} + +Feature debugGreatCircle(Position start, Position end, + {int npoints = 2}) { + print("Input start: Position(${start[0]}, ${start[1]})"); + print("Input end: Position(${end[0]}, ${end[1]})"); + + // Current assignment (what you have) + num lng1 = degreesToRadians(start[0]!); // longitude + num lat1 = degreesToRadians(start[1]!); // latitude + num lng2 = degreesToRadians(end[0]!); // longitude + num lat2 = degreesToRadians(end[1]!); // latitude + + print("After assignment:"); + print("lng1 (from start[0]): ${radiansToDegrees(lng1)}"); + print("lat1 (from start[1]): ${radiansToDegrees(lat1)}"); + print("lng2 (from end[0]): ${radiansToDegrees(lng2)}"); + print("lat2 (from end[1]): ${radiansToDegrees(lat2)}"); + + List line = []; + + // Just add the start and end points to see what happens + for (int i = 0; i <= npoints; i++) { + double f = i / npoints; + + if (f == 0) { + // Start point + Position point = Position(lng1, lat1); + line.add(point); + print( + "Start point created: Position(${lng1}, ${lat1}) = [${radiansToDegrees(lng1)}, ${radiansToDegrees(lat1)}]"); + } else if (f == 1) { + // End point + Position point = Position(lng2, lat2); + line.add(point); + print( + "End point created: Position(${lng2}, ${lat2}) = [${radiansToDegrees(lng2)}, ${radiansToDegrees(lat2)}]"); + } else { + // For simplicity, just add a midpoint + double midLng = (lng1 + lng2) / 2; + double midLat = (lat1 + lat2) / 2; + Position point = Position(midLng, midLat); + line.add(point); + print( + "Mid point created: Position(${midLng}, ${midLat}) = [${radiansToDegrees(midLng)}, ${radiansToDegrees(midLat)}]"); + } + } + + return Feature(geometry: LineString(coordinates: line)); +} diff --git a/test/components/great_circle_test.dart b/test/components/great_circle_test.dart new file mode 100644 index 00000000..9f7a03ec --- /dev/null +++ b/test/components/great_circle_test.dart @@ -0,0 +1,131 @@ +import 'package:turf/great_circle.dart'; +import 'package:test/test.dart'; +import 'package:turf/helpers.dart'; + +/* Note: This test function could be worked on further especially when dealing with precision. +Due to the general nature of greatCircle as a visual linestring, this should not be a big issue. +However, having further precision (testing changes from 1 decimal place to 4-6 can make it more useful for when npoints is large OR for shorter distances) +A +*/ +void main() { + //First test - simple coordinates + + final start = Position(0, -90); + final end = Position(0, -80); + + List> resultsFirstTest1 = [ + [-90.0, 0.0], + [-88.0, 0.0], + [-86.0, 0.0], + [-84.0, 0.0], + [-82.0, 0.0], + [-80.0, 0.0] + ]; + List> resultsFirstTest2 = [ + [-90.0, 0.0], + [-89.0, 0.0], + [-88.0, 0.0], + [-87.0, 0.0], + [-86.0, 0.0], + [-85.0, 0.0], + [-84.0, 0.0], + [-83.0, 0.0], + [-82.0, 0.0], + [-81.0, 0.0], + [-80.0, 0.0] + ]; + + test('Great circle simple tests:', () { + var resultFirst1 = greatCircle(start, end, npoints: 5); + var convertedResultFirst1 = resultFirst1.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultFirst1, resultsFirstTest1); + + var resultFirst2 = greatCircle(start, end, npoints: 10); + var convertedResultFirst2 = resultFirst2.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultFirst2, resultsFirstTest2); + }); + + // Second test - intermediate coordiantes (non-straight lines) + final start2 = Position(48, -122); + final end2 = Position(39, -77); + + List> resultsSecondTest1 = [ + [-122.0, 48.0], + [-97.7, 45.8], + [-77.0, 39.0] + ]; + List> resultsSecondTest2 = [ + [-122.0, 48.0], + [-109.6, 47.5], + [-97.7, 45.8], + [-86.8, 42.8], + [-77.0, 39.0] + ]; + + test('Great circle intermediate tests:', () { + var resultSecond1 = greatCircle(start2, end2, npoints: 2); + var convertedResultSecond1 = resultSecond1.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultSecond1, resultsSecondTest1); + + var resultSecond2 = greatCircle(start2, end2, npoints: 4); + var convertedResultSecond2 = resultSecond2.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultSecond2, resultsSecondTest2); + }); + + // Third test - complex coordinates (crossing anti-meridian) + + final start3 = Position(-21, 143); + final end3 = Position(41, -140); + + List> resultsThirdTest1 = [ + [143.0, -21.0], + [176.7, 12.7], + [-140, 41] + ]; + List> resultsThirdTest2 = [ + [143.0, -21.0], + [160.2, -4.4], + [176.7, 12.7], + [-164.6, 28.5], + [-140, 41] + ]; + test('Great circle complex tests:', () { + var resultThird1 = greatCircle(start3, end3, npoints: 2); + var convertedResultThird1 = resultThird1.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultThird1, resultsThirdTest1); + + var resultThird2 = greatCircle(start3, end3, npoints: 4); + var convertedResultThird2 = resultThird2.geometry?.coordinates + .map((pos) => [ + double.parse(radiansToDegrees(pos[0]).toStringAsFixed(1)), + double.parse(radiansToDegrees(pos[1]).toStringAsFixed(1)) + ]) + .toList(); + expect(convertedResultThird2, resultsThirdTest2); + }); +}