///|
pub fn initial_bearing_degrees(
  a : Point,
  b : Point,
) -> Result[Double, RouteError] {
  match validate_point(a, 0) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  match validate_point(b, 1) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let lat1 = deg_to_rad(a.lat)
  let lat2 = deg_to_rad(b.lat)
  let dlon = deg_to_rad(b.lon - a.lon)
  let y = @math.sin(dlon) * @math.cos(lat2)
  let x = @math.cos(lat1) * @math.sin(lat2) -
    @math.sin(lat1) * @math.cos(lat2) * @math.cos(dlon)
  Ok(normalize_bearing_degrees(rad_to_deg(@math.atan2(y, x))))
}

///|
pub fn final_bearing_degrees(
  a : Point,
  b : Point,
) -> Result[Double, RouteError] {
  match initial_bearing_degrees(b, a) {
    Ok(value) => Ok(normalize_bearing_degrees(value + 180.0))
    Err(err) => Err(err)
  }
}

///|
pub fn segment_metric(
  start : Point,
  finish : Point,
  index? : Int = 0,
) -> Result[SegmentMetric, RouteError] {
  let distance_m = match haversine_meters(start, finish) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let bearing_deg = match initial_bearing_degrees(start, finish) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  Ok(SegmentMetric(index, start, finish, distance_m, bearing_deg))
}

///|
pub fn segment_metrics(
  points : ArrayView[Point],
) -> Result[Array[SegmentMetric], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = Array::new()
  if points.length() < 2 {
    return Ok(out)
  }
  for i in 1.. out.push(metric)
      Err(err) => return Err(err)
    }
  }
  Ok(out)
}

///|
pub fn route_segment_distances(
  points : ArrayView[Point],
) -> Result[Array[Double], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = Array::new()
  if points.length() < 2 {
    return Ok(out)
  }
  for i in 1.. out.push(distance)
      Err(err) => return Err(err)
    }
  }
  Ok(out)
}

///|
pub fn route_segment_spread(
  points : ArrayView[Point],
) -> Result[SegmentSpread, RouteError] {
  let distances = match route_segment_distances(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  if distances.length() == 0 {
    return Ok(SegmentSpread(0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -1))
  }
  let mut total = 0.0
  let mut shortest = 999999999999.0
  let mut longest = 0.0
  let mut longest_index = 0
  for i, distance in distances {
    total = total + distance
    shortest = min_double(shortest, distance)
    if distance > longest {
      longest = distance
      longest_index = i
    }
  }
  let average = total / distances.length().to_double()
  let mut squared_error = 0.0
  for distance in distances {
    let delta = distance - average
    squared_error = squared_error + delta * delta
  }
  let variance = squared_error / distances.length().to_double()
  let standard_deviation = variance.sqrt()
  let coefficient = if average > 0.0 {
    standard_deviation / average
  } else {
    0.0
  }
  Ok(
    SegmentSpread(
      distances.length(),
      total,
      shortest,
      longest,
      average,
      standard_deviation,
      coefficient,
      longest_index,
    ),
  )
}

///|
pub fn route_evenness_score(
  points : ArrayView[Point],
) -> Result[Int, RouteError] {
  let spread = match route_segment_spread(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  if spread.segment_count == 0 {
    return Ok(100)
  }
  let raw = (100.0 - spread.coefficient_of_variation * 100.0).round().to_int()
  Ok(max_int(0, min_int(100, raw)))
}

///|
pub fn route_is_evenly_sampled(
  points : ArrayView[Point],
  max_coefficient_of_variation? : Double = 0.35,
) -> Result[Bool, RouteError] {
  if max_coefficient_of_variation.is_nan() ||
    max_coefficient_of_variation.is_inf() ||
    max_coefficient_of_variation < 0.0 {
    return Err(InvalidFraction(fraction=max_coefficient_of_variation))
  }
  let spread = match route_segment_spread(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  Ok(spread.coefficient_of_variation <= max_coefficient_of_variation)
}

///|
pub fn route_bearings(
  points : ArrayView[Point],
) -> Result[Array[Double], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = Array::new()
  if points.length() < 2 {
    return Ok(out)
  }
  for i in 1.. out.push(value)
      Err(err) => return Err(err)
    }
  }
  Ok(out)
}

///|
pub fn route_turn_angles(
  points : ArrayView[Point],
) -> Result[Array[Double], RouteError] {
  let bearings = match route_bearings(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let out = Array::new()
  if bearings.length() < 2 {
    return Ok(out)
  }
  for i in 1.. Result[Double, RouteError] {
  let turns = match route_turn_angles(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let mut total = 0.0
  for turn in turns {
    total = total + abs_double(turn)
  }
  Ok(total)
}

///|
pub fn cumulative_distances(
  points : ArrayView[Point],
) -> Result[Array[Double], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = Array::new(capacity=points.length())
  out.push(0.0)
  if points.length() < 2 {
    return Ok(out)
  }
  let mut total = 0.0
  for i in 1.. {
        total = total + distance
        out.push(total)
      }
      Err(err) => return Err(err)
    }
  }
  Ok(out)
}

///|
pub fn cumulative_markers(
  points : ArrayView[Point],
) -> Result[Array[RouteMarker], RouteError] {
  let distances = match cumulative_distances(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let out = Array::new(capacity=distances.length())
  for i, distance in distances {
    out.push(RouteMarker(i, points[i], distance))
  }
  Ok(out)
}

///|
pub fn route_centroid(points : ArrayView[Point]) -> Result[Point, RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let mut lat = 0.0
  let mut lon = 0.0
  for point in points {
    lat = lat + point.lat
    lon = lon + point.lon
  }
  Ok(
    Point(lat / points.length().to_double(), lon / points.length().to_double()),
  )
}

///|
pub fn route_midpoint(
  points : ArrayView[Point],
) -> Result[RouteSample, RouteError] {
  let total = match route_distance_meters(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  point_at_distance(points, total / 2.0)
}

///|
pub fn route_endpoints(
  points : ArrayView[Point],
) -> Result[(Point, Point), RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  Ok((points[0], points[points.length() - 1]))
}

///|
pub fn is_closed_route(points : ArrayView[Point]) -> Result[Bool, RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  if points.length() < 2 {
    Ok(false)
  } else {
    Ok(points_equal(points[0], points[points.length() - 1]))
  }
}

///|
pub fn close_route(
  points : ArrayView[Point],
) -> Result[Array[Point], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = points.to_owned()
  if points.length() > 1 &&
    !points_equal(points[0], points[points.length() - 1]) {
    out.push(points[0])
  }
  Ok(out)
}

///|
pub fn reverse_route(
  points : ArrayView[Point],
) -> Result[Array[Point], RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let out = Array::new(capacity=points.length())
  let mut i = points.length()
  while i > 0 {
    i = i - 1
    out.push(points[i])
  }
  Ok(out)
}

///|
pub fn route_stats(points : ArrayView[Point]) -> Result[RouteStats, RouteError] {
  match require_points(points) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  let distance_m = match route_distance_meters(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let bbox = match bounding_box(points) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let direct_distance_m = if points.length() >= 2 {
    match haversine_meters(points[0], points[points.length() - 1]) {
      Ok(value) => value
      Err(err) => return Err(err)
    }
  } else {
    0.0
  }
  let center = bbox.center()
  let segment_count = max_int(points.length() - 1, 0)
  let mut min_segment_m = 0.0
  let mut max_segment_m = 0.0
  let mut average_segment_m = 0.0
  if segment_count > 0 {
    min_segment_m = 999999999999.0
    for i in 1.. value
        Err(err) => return Err(err)
      }
      min_segment_m = min_double(min_segment_m, distance)
      max_segment_m = max_double(max_segment_m, distance)
    }
    average_segment_m = distance_m / segment_count.to_double()
  }
  let sinuosity = if direct_distance_m > 0.0 {
    distance_m / direct_distance_m
  } else if distance_m > 0.0 {
    999999999999.0
  } else {
    1.0
  }
  Ok(
    RouteStats(
      points.length(),
      segment_count,
      distance_m,
      direct_distance_m,
      sinuosity,
      min_segment_m,
      max_segment_m,
      average_segment_m,
      bbox,
      center,
    ),
  )
}

///|
pub fn RouteStats::report(self : RouteStats) -> String {
  [
    "point_count: \{self.point_count}",
    "segment_count: \{self.segment_count}",
    "distance_m: \{round_meters(self.distance_m)}",
    "direct_distance_m: \{round_meters(self.direct_distance_m)}",
    "sinuosity: \{self.sinuosity}",
    "min_segment_m: \{round_meters(self.min_segment_m)}",
    "max_segment_m: \{round_meters(self.max_segment_m)}",
    "average_segment_m: \{round_meters(self.average_segment_m)}",
    "bbox: \{self.bbox.to_string()}",
    "center: \{self.center.to_string()}",
  ].join("\n")
}

///|
pub fn BBox::contains(self : BBox, point : Point) -> Bool {
  point.lat >= self.min_lat &&
  point.lat <= self.max_lat &&
  point.lon >= self.min_lon &&
  point.lon <= self.max_lon
}

///|
pub fn BBox::width_degrees(self : BBox) -> Double {
  self.max_lon - self.min_lon
}

///|
pub fn BBox::height_degrees(self : BBox) -> Double {
  self.max_lat - self.min_lat
}

///|
pub fn BBox::diagonal_meters(self : BBox) -> Result[Double, RouteError] {
  haversine_meters(
    Point(self.min_lat, self.min_lon),
    Point(self.max_lat, self.max_lon),
  )
}

///|
pub fn BBox::area_square_meters(self : BBox) -> Result[Double, RouteError] {
  let south_west = Point(self.min_lat, self.min_lon)
  let south_east = Point(self.min_lat, self.max_lon)
  let north_west = Point(self.max_lat, self.min_lon)
  let width = match haversine_meters(south_west, south_east) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  let height = match haversine_meters(south_west, north_west) {
    Ok(value) => value
    Err(err) => return Err(err)
  }
  Ok(width * height)
}

///|
pub fn BBox::expand_to_include(
  self : BBox,
  point : Point,
) -> Result[BBox, RouteError] {
  match validate_point(point, 0) {
    Ok(_) => ()
    Err(err) => return Err(err)
  }
  Ok(
    BBox(
      min_double(self.min_lat, point.lat),
      min_double(self.min_lon, point.lon),
      max_double(self.max_lat, point.lat),
      max_double(self.max_lon, point.lon),
    ),
  )
}

///|
pub fn BBox::merge(self : BBox, other : BBox) -> BBox {
  BBox(
    min_double(self.min_lat, other.min_lat),
    min_double(self.min_lon, other.min_lon),
    max_double(self.max_lat, other.max_lat),
    max_double(self.max_lon, other.max_lon),
  )
}

///|
pub fn BBox::pad_degrees(
  self : BBox,
  padding : Double,
) -> Result[BBox, RouteError] {
  if padding.is_nan() || padding.is_inf() || padding < 0.0 {
    return Err(InvalidDistance(distance_m=padding))
  }
  Ok(
    BBox(
      (self.min_lat - padding).clamp(min=-90.0, max=90.0),
      (self.min_lon - padding).clamp(min=-180.0, max=180.0),
      (self.max_lat + padding).clamp(min=-90.0, max=90.0),
      (self.max_lon + padding).clamp(min=-180.0, max=180.0),
    ),
  )
}