CREATE FUNCTION fdistance(pnlat0 double precision, pnlon0 double precision, pnlat1 double precision, pnlon1 double precision) RETURNS double precision
    LANGUAGE plpgsql
    AS $$
DECLARE
-- расчет расстояния между двумя координатами
Distance double precision; -- расстояние в метрах

BEGIN
/*
Кратчайшее расстояние между двумя точками на земной поверхности (если принять ее за сферу) определяется зависимостью:
cos(d) = sin(fА) * sin(fB) + cos(fА) * cos(fB) * cos(lА - lB)
где
fА и fB — широты (в радианах),
lА, lB — долготы данных пунктов (в радианах),
d — расстояние между пунктами (в радианах)
перевод в метры осуществляется по формуле:
L = d * R
где R = 6371210.0 м — средний радиус земного шара в книгах русских авторов
и 6378137.0 м - в расчетах OpenStreet и OpenLayers для WGS84

ПыСы: расстояние между двумя точками на [земной поверхности] плоскости
вычисляется также по формуле
L = SQRT(POWER(x1 - x0, 2) + POWER(y1 - y0, 2)) (в метрах)
где (x0,y0) и (x1,y1) - координаты точек В МЕТРАХ, ибо плоскость
*/

with rad as (
select radians(pnLat0) as Lat0R, radians(pnLat1) as Lat1R
)
select
round((acos(sin(rad.Lat0R) * sin(rad.Lat1R) + cos(rad.Lat0R) * cos(rad.Lat1R) * cos(radians(pnLon0) - radians(pnLon1))) * na)::numeric, 0)
into Distance
from gps.tellips, rad
where nid = 0; -- WGS84

/*
gps.tellips.na = 6378137
*/

RETURN Distance;
END;
$$;