Вычисление площади полигона на сфере и на эллипсоиде
Одна-две фразы по существу.
Общие положения
Определим полигон как простой многоугольник — участок поверхности, ограниченный замкнутой полилинией без самопересечений.
Полилиния в свою очередь — ломаная, образованная отрезками геодезических линий.
Геодезическая линия на плоскости — это прямая; геодезическая линия на сфере — дуга большой окружности.
Полный поворот контура
В общем случае определение площади многоугольника на искривлённой поверхности — нетривиальная задача. Нужно интегрировать по поверхности с пределами, заданными неявно. К счастью, математика может предложить обходные пути решения задачи.
Представим себе точку, движущуюся вдоль контура полигона. Вершины являются точками поворота. Внутренний угол при вершине θ равен разности направлений α в предыдущую и следующую вершины, а поворот есть угол τ, смежный внутреннему:
На евклидовой плоскости, обойдя любой замкнутый контур без самопересечений, точка совершает поворот ровно на одну окружность — 360°, или 2π радиан. В случае многоугольника этот поворот складывается из суммы поворотов в вершинах.
На поверхности с ненулевой гауссовой кривизной общий поворот отличается от 2π на величину избытка или недостатка, пропорционального кривизне поверхности и площади фигуры.
Площадь полигона на сфере
Сферический избыток
В общем случае кривизна поверхности меняется в каждой точке, но не на сфере! Кривизна сферы постоянна, и площадь замкнутой фигуры однозначно соотносится с полным поворотом контура.
Отличие полного поворота от 2π радиан называется сферическим избытком ε, который пропорционален площади полигона S:
где R — радиус сферы.
Алгоритм вычисления площади
Пусть n-угольник задан координатами вершин φi, λi, где φi — широта i-ой вершины, λi — долгота, i = 1, … , n.
- Для каждой стороны из решения обратной геодезической задачи для её конечных вершин находим прямые и обратные азимуты αi, i+1 и αi+1, i.
- Для каждой вершины по азимутам αi, i−1 в предыдущую и αi, i+1 в последующую вершины находим поворот τi и добавляем его к полному повороту τ.
- Вычисляем сферический избыток ε.
- Вычисляем площадь полигона S.
Радиус сферы
Радиус сферы выбирается таким образом, чтобы минимизировать погрешности, вызванные тем, что кривизна земной поверхности в разных точках различна.
В геодезической и картографической практике в качестве модели Земли принимают эллипсоид вращения, характеризуемый величинами экваториального радиуса a и сжатия f. Соответствующая сфера должна иметь такую же площадь, как и принимаемый за основу эллипсоид. Такая сфера называется эквивалентной (authalic sphaera по латыни).
Ниже мы приведём формулы вычисления радиуса эквивалентной сферы по параметрам эллипсоида. Пока же поближе посмотрим на пару-тройку популярных программных пакетов.
MapInfo — единственная известная нам программа, вычисляющая площади на сфере. Несложные эксперименты («реверс-инжиниринг») показывает, что радиус сферы при этом не зависит от того, на каком эллипсоиде построена координатная система обрабатываемых данных. Радиус всегда равен 6370997 метров ровно. Анализ используемых земных эллипсоидов показывает, что это соответствует радиусу эквивалентной сферы эллипсоида Кларка 1866 года, определяемого величинами экаториального радиуса a = 6378206.4 и полярного радиуса b = 6356583.8; точное значение Rauth = 6370997.2406, округлённое до метра — 6370997 м ровно.
В большинстве своём программы GIS не вычисляют площадей на сфере ArcGIS