Точное определение региона по GPS координатам

в 6:51, , рубрики: gps, метки:

При разработке одного приложения возникла проблема разграничения доступа для регионов.

Встала проблема определения принадлежности объекта к какому-либо региону России по его GPS координатам

Первое, что мы начали использовать — это API Google,
после того как прописали алиасы к возвращаемым строкам и оплаты доступа (чтобы убрать лимит на запросы) все заработало.
И все было нормально пока гугл не сменил выдачу, например было раньше: Moskovskaya oblast', стало Moscow oblast'
Тут то и было решено не надеяться на гугл, а определять регион своими силами.

Точное определение региона по GPS координатам

Первым, что пришло на ум было: распарсить возвращаемые данные Яндекса (в яндекс-картах версии 1.x был модуль для работы с регионами)
Так и поступили, получилось 2 таблицы.
В одной хранились названия регионов, во второй координаты многоугольника, описывающие регион. (83 — регионы и 8500 — точки )

Был найден алгоритм, для проверки лежит ли точка в многоугольнике. Код был переписан с Си на PHP, разбирать как он работает нет смысла, он просто работает.

$sx,$sy — координаты точки, которую проверяем
$coords — координаты точек выпуклого многоугольника (отсортированы в порядке обхода):
$coords = array(
array('x'=> 66.6634, 'y' => '66.4433'),
и тп.
)
$x, $y — названия ключей в массиве

class Algo{
  public static function into_poly($sx, $sy, &$coords, $x='x', $y='y')
    {
      $pj=0;
      $pk=0;
      $wrkx=0;
      $yu = 0;
      $yl = 0;
      $n = count($coords);
      for ($pj=0; $pj<$n; $pj++)
      {
        $yu = $coords[$pj][$y]>$coords[($pj+1)%$n][$y]?$coords[$pj][$y]:$coords[($pj+1)%$n][$y];
        $yl = $coords[$pj][$y]<$coords[($pj+1)%$n][$y]?$coords[$pj][$y]:$coords[($pj+1)%$n][$y];
        if ($coords[($pj+1)%$n][$y] - $coords[$pj][$y])
          $wrkx = $coords[$pj][$x] + ($coords[($pj+1)%$n][$x] - $coords[$pj][$x])*($sy - $coords[$pj][$y])/($coords[($pj+1)%$n][$y] - $coords[$pj][$y]);
        else
          $wrkx = $coords[$pj][$x];
        if ($yu >= $sy)
          if ($yl < $sy)
          {
            if ($sx > $wrkx)
              $pk++;
            if (abs($sx - $wrkx) < 0.00001) return 1;
          }
        if ((abs($sy - $yl) < 0.00001) && (abs($yu - $yl) < 0.00001) && (abs(abs($wrkx - $coords[$pj][$x]) + abs($wrkx - $coords[($pj+1)%$n][$x]) - abs($coords[$pj][$x] - $coords[($pj+1)%$n][$x])) < 0.0001))
          return 1;
      }
      if ($pk%2)
        return 1;
      else
        return 0;
    }
}

* This source code was highlighted with Source Code Highlighter.

Пример вызова:
$coords = array(
array('lng'=> 66.6634, 'lat' => '66.4433'),
array('lng'=> 66.6534, 'lat' => '66.4433'),
array('lng'=> 66.6434, 'lat' => '66.4433'),
array('lng'=> 66.6334, 'lat' => '66.4433'),
);
$in = Algo::into_poly(66.4455, 66.2255, &$coords, $x='lng', $y='lat');


Вернет true, если точка лежит внутри многоугольника.

Теперь просто перебором областей для точки мы получаем ответ.

Всё бы было хорошо, но 8500 точек границ регионов для России — это очень мало. Разброс получится +- 50 км

Чтобы получить точность больше – нам нужно больше точек границ регионов.
Поискав на просторах был найден ресурс, который нам свободно отдает эти точки: gis-lab.info/qa/rusbounds-rosreestr.html (спасибо энтузиастам за вашу работу)
Я почему-то выбрал формат KML (наверное потому, знал чем его открывать — это формат, который поддерживает Google Earth)
Распарсил данные и теперь таблицы получились пожирнее. Таблица с координатами точек границ регионов стала 620 000 точек (34 мб)

Теперь если проверять по старому алгоритму принадлежность точки к многоугольнику региона — на стареньком Core 2 Duo — около 70 секунд.
Не вариант. Алгоритм нужно оптимизировать:

Например, Москва:
Имеем 2,064 точек многоугольника (на самом деле несколько многоугольников):
Точное определение региона по GPS координатам

Но сейчас нам достаточно определить всего-лишь возможность, что точка будет в этой области. Поэтому мы воспользуемся алгоритмом обхода точек и построения выпуклого многоугольника.
algolist.manual.ru

вот его реализация на PHP:

class Algo{
    public static function sort_points($mass, $x = 'x', $y = 'y'){
      $current=0;
      $next=0;
      $p1 = $mass[0];
      $mass2 = array();
      //определяем первую точку
      for ($i=1; $i<count($mass); $i++){
        if ( ($p1[$y] > $mass[$i][$y]) || ($p1[$y] == $mass[$i][$y] && $p1[$x] < $mass[$i][$x]) ) {
          $p1 = $mass[$i];
          $current = $i;
        }
      }
      $n = count($mass);

      $p0 = $p1;
      $p0[$x]--;

      $mass2[] = $mass[$current];
      $first = $current;

      //начинаем перебор всех элементов

      do{
        $cmax_not_set=1;
        for ( $i=0; $i<$n; $i++ ) {
          $key = $i;
          //продолжаем если такой элемент уже выбран
          if ( $mass[$current][$x] == $mass[$key][$x] && $mass[$current][$y] == $mass[$key][$y] ) continue;
          //берем 1ю точку
          if ($cmax_not_set) {
            $next = $key;
            $cmax_not_set=0;
            continue;
          }

          $v1_x = $mass[$key][$x] - $mass[$current][$x];
          $v1_y = $mass[$key][$y] - $mass[$current][$y];

          $v2_x = $mass[$next][$x] - $mass[$current][$x];
          $v2_y = $mass[$next][$y] - $mass[$current][$y];

          //магия
          $c = $v1_x * $v2_y - $v1_y * $v2_x;

          if ( $c > 0 ) $next = $key;
          // if ( ($c==0) && ( self::_dist($mass[$current], $mass[$key], $x, $y) > self::_dist($mass[$current], $mass[$next], $x, $y) ) ) $next = $key;

        }
        //записываем найденный элемент в массив
        $mass2[] = $mass[$next];
        $current = $next;
      }while($first != $next);
      return $mass2;
    }
}

* This source code was highlighted with Source Code Highlighter.

На выходе мы получим массив точек выпуклого многоугольника.
Для Москвы; 19 точек:
Точное определение региона по GPS координатам

На рисунке первым слоем отмечены все 2500 точек Москвы + слой, который мы получили после расчета выпуклого многоугольника.
Сгенерируем для каждого региона такую область.

Теперь доработаем наш алгоритм и будем сначала проверять возможность, что точка может лежать в области – таких областей может быть несколько.
Теперь после прогона всех областей, мы прогоним области, в которых возможно, что лежит точка, алгоритмом описанным выше.
Итого: 0.177сек на слабеньком Core 2 Duo

На самом деле есть еще одна сложность: это – анклавы (Москва лежит внутри Московской области, как и Питер, было решено использовать приоритеты, сначала проверяется лежит ли точка в Москве, а потом уже в Московской области)

Автор: kiryam

* - обязательные к заполнению поля


https://ajax.googleapis.com/ajax/libs/jquery/3.4.1/jquery.min.js