docs: def great_circle_vec(lat1, lng1, lat2, lng2, earth_radius=6371009):
      """
      Vectorized function to calculate the great-circle distance between two
      points or between vectors of points, using haversine.
  
      Parameters
      ----------
      lat1 : float or array of float
      lng1 : float or array of float
      lat2 : float or array of float
      lng2 : float or array of float
      earth_radius : numeric
          radius of earth in units in which distance will be returned (default is
          meters)
  
      Returns
      -------
      distance : float or vector of floats
          distance or vector of distances from (lat1, lng1) to (lat2, lng2) in
          units of earth_radius
      """
  
      phi1 = np.deg2rad(lat1)
      phi2 = np.deg2rad(lat2)
      d_phi = phi2 - phi1
  
      theta1 = np.deg2rad(lng1)
      theta2 = np.deg2rad(lng2)
      d_theta = theta2 - theta1
  
      h = np.sin(d_phi / 2) ** 2 + np.cos(phi1) * np.cos(phi2) * np.sin(d_theta / 2) ** 2 
      h = np.minimum(1.0, h) # protect against floating point errors
  
      arc = 2 * np.arcsin(np.sqrt(h))
  
      # return distance in units of earth_radius
      distance = arc * earth_radius
      return distance