offset method

  1. @override
GeoLatLng offset(
  1. GeoLatLng from,
  2. double distanceInMeter,
  3. double bearing
)
override

Vincenty inverse calculation

More on Wikipedia

Implementation

@override
GeoLatLng offset(final GeoLatLng from, final double distanceInMeter,
    final double bearing) {
  final double equatorialRadius = EQUATOR_RADIUS;
  final double polarRadius = POLAR_RADIUS;
  final double flattening = FLATTENING; // WGS-84 ellipsoid params

  final double latitude = from.latitudeInRad;
  final double longitude = from.longitudeInRad;

  final double alpha1 = degToRadian(bearing);
  final double sinAlpha1 = math.sin(alpha1);
  final double cosAlpha1 = math.cos(alpha1);

  final double tanU1 = (1 - flattening) * math.tan(latitude);
  final double cosU1 = 1 / math.sqrt((1 + tanU1 * tanU1));
  final double sinU1 = tanU1 * cosU1;

  final double sigma1 = math.atan2(tanU1, cosAlpha1);
  final double sinAlpha = cosU1 * sinAlpha1;
  final double cosSqAlpha = 1 - sinAlpha * sinAlpha;
  final double dfUSq = cosSqAlpha *
      (equatorialRadius * equatorialRadius - polarRadius * polarRadius) /
      (polarRadius * polarRadius);
  final double a = 1 +
      dfUSq / 16384 * (4096 + dfUSq * (-768 + dfUSq * (320 - 175 * dfUSq)));
  final double b =
      dfUSq / 1024 * (256 + dfUSq * (-128 + dfUSq * (74 - 47 * dfUSq)));

  double sigma = distanceInMeter / (polarRadius * a);
  double sigmaP = 2 * PI;

  double sinSigma = 0.0;
  double cosSigma = 0.0;
  double cos2SigmaM = 0.0;
  double deltaSigma;
  int maxIterations = 200;

  do {
    cos2SigmaM = math.cos(2 * sigma1 + sigma);
    sinSigma = math.sin(sigma);
    cosSigma = math.cos(sigma);
    deltaSigma = b *
        sinSigma *
        (cos2SigmaM +
            b /
                4 *
                (cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM) -
                    b /
                        6 *
                        cos2SigmaM *
                        (-3 + 4 * sinSigma * sinSigma) *
                        (-3 + 4 * cos2SigmaM * cos2SigmaM)));
    sigmaP = sigma;
    sigma = distanceInMeter / (polarRadius * a) + deltaSigma;
  } while ((sigma - sigmaP).abs() > 1e-12 && --maxIterations > 0);

  if (maxIterations == 0) {
    throw new StateError("offset calculation faild to converge!");
  }

  final double tmp = sinU1 * sinSigma - cosU1 * cosSigma * cosAlpha1;
  final double lat2 = math.atan2(
      sinU1 * cosSigma + cosU1 * sinSigma * cosAlpha1,
      (1 - flattening) * math.sqrt(sinAlpha * sinAlpha + tmp * tmp));

  final double lambda = math.atan2(
      sinSigma * sinAlpha1, cosU1 * cosSigma - sinU1 * sinSigma * cosAlpha1);
  final double c =
      flattening / 16 * cosSqAlpha * (4 + flattening * (4 - 3 * cosSqAlpha));
  final double l = lambda -
      (1 - c) *
          flattening *
          sinAlpha *
          (sigma +
              c *
                  sinSigma *
                  (cos2SigmaM +
                      c * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)));

  double lon2 = longitude + l;
  // print("LA ${radianToDeg(lat2)}, LO ${radianToDeg(lon2)}");

  if (lon2 > PI) {
    lon2 = lon2 - 2 * PI;
  }
  if (lon2 < -1 * PI) {
    lon2 = lon2 + 2 * PI;
  }

  return new GeoLatLng(radianToDeg(lat2), radianToDeg(lon2));
}