Серверная кластеризация маркеров на карте. От теории к практике

    Привет Хабр. История начинается с того что мы решили сделать гео сервис с возможностью размещения меток на карте самими пользователями.
    И когда решили залить в базу 1 миллион маркеров то поняли, что даже если запрашивать маркеры только в определенном радиусе то все работает очень медленно и кластеризация на клиенте тоже не вариант :)

    А где-то под этим лесом находится манхетен




    Кластеризация это очень интересная задача, есть клиентская, есть серверная и все крупные картографические сервисы ее поддерживают. Известно множество видов кластеризации, но я в своих поисках выбрал Grid-кластеризацию. Про нее расказывали сотрудники Яндекса в нескольких своих докладах на конференциях. Она используется в яндекс картах для кластеризации маркеров.
    Yandex кластеризация


    Немного теории


    Grid-кластеризация основана на проекции Меркатора, вкратце это проецирование сферической поверхности земного шара на прямоугольную плоскость. Это как раз все атласы, карты итд. Как понятно из названия в этой кластеризации мы покрываем карту сеткой и вычисляем количество маркеров в каждой ячейке это сетки.

    Термины:

    Тайл(tile, плитка) — то по сути просто разделение карты на квадраты для определенного zoom-а. На первом уровне у нас карта будет состоять из 1 квадрата(тайла) на 0 уровне, на 1 уровне из 4-х квадратов(тайлов), на втором уровне на 16 на третьем на 64 и т.д. По сути это 4 в степени уровня.


    Quadkey — это уникальный идентификатор тайла на любом уровне. На первом уровне нумеруем тайлы 0,1,2,3 справа налево, сверху вниз. Затем на следующем уровне для каждого тайла из верхнего уровня делаем такую же нумера но уже добавляя quadkey родителя т.е. 00,01,02,03,10,11,12,13 и т.д. Более понятно это становится по рисунку:

    Посути мы получили число в 4-чной системе исчисления. Теперь у нас каждый тайл имеет свой уникальный идентификатор на любом уровне zoom-а. Из zoom-а, широты и долготы можно получить уже сам quadkey тайла т.е. любую точку можно преобразовать в quadkey соответсвующего уровня zoom-а. А вот так будет выглядеть quadkey для конкретных координат:
    longitude -79.377807617187500
    latitude 43.653785705566406
    quadkey 13940830302567 (base4: 03022313122033033011213)

    Описание этого подхода на англиской языке от Microsoft:
    Bing Maps Tile System
    В этой статье реализации всех этих преобразований написана на C#, я побыстрому переписал ее на PHP:
    TileSystem.php
    class TileSystem
    {
      const EarthRadius = 6378137;
      const MinLatitude = -90;
      const MaxLatitude = 90;
      const MinLongitude = -180;
      const MaxLongitude = 180;
      const MaxZoom = 23;
    
      private static function clip( $n, $minValue, $maxValue)
      {
        return min( max( $n, $minValue), $maxValue);
      }
    
      public static function mapSize( $levelOfDetail)
      {
        return 256 << $levelOfDetail;
      }
    
      public static function GroundResolution( $latitude, $levelOfDetail)
      {
        $latitude = self::clip( $latitude, self::MinLatitude, self::MaxLatitude);
        return cos( $latitude * pi() / 180) * 2 * pi() * self::EarthRadius / self::mapSize( $levelOfDetail);
      }
    
      public static function mapScale( $latitude, $levelOfDetail, $screenDpi)
      {
        return self::groundResolution( $latitude, $levelOfDetail) * $screenDpi / 0.0254;
      }
    
      public static function latToPixelY( $latitude, $levelOfDetail)
      {
        $latitude = self::clip( $latitude, self::MinLatitude, self::MaxLatitude);
        $sinLatitude = sin( $latitude * pi() / 180);
        $y = 0.5 - log((1 + $sinLatitude) / (1 - $sinLatitude)) / (4 * pi());
        $mapSize = self::mapSize( $levelOfDetail);
        $pixelY = (int) self::clip( $y * $mapSize + 0.5, 0, $mapSize - 1);
        return $pixelY;
      }
    
      public static function longToPixelX( $longitude, $levelOfDetail)
      {
        $longitude = self::clip( $longitude, self::MinLongitude, self::MaxLongitude);
        $x = ($longitude + 180) / 360;
        $mapSize = self::mapSize( $levelOfDetail);
        $pixelX = (int) self::clip( $x * $mapSize + 0.5, 0, $mapSize - 1);
        return $pixelX;
      }
    
      public static function pixelXToLong( $pixelX, $levelOfDetail)
      {
        $mapSize = self::mapSize( $levelOfDetail);
        $x = ( self::clip( $pixelX, 0, $mapSize - 1) / $mapSize) - 0.5;
        $longitude = 360 * $x;
        return $longitude;
      }
    
      public static function pixelYToLat( $pixelY, $levelOfDetail)
      {
        $mapSize = self::mapSize( $levelOfDetail);
        $y = 0.5 - (self::clip( $pixelY, 0, $mapSize - 1) / $mapSize);
        $latitude = 90 - 360 * atan( exp(-$y * 2 * pi())) / pi();
        return $latitude;
      }
    
      public static function pixelYToTileY( $pixelY)
      {
        return $pixelY / 256;
      }
    
      public static function pixelXToTileX( $pixelX)
      {
        return $pixelX / 256;
      }
    
      public static function tileXYToPixelXY( $tileX, $tileY, &$pixelX, &$pixelY)
      {
        $pixelX = $tileX * 256;
        $pixelY = $tileY * 256;
      }
    
      public static function tileXYToQuadKey( $tileX, $tileY, $levelOfDetail)
      {
        if ($levelOfDetail == 0)
          return "0";
    
        $quadKey = "";
        for ($i = $levelOfDetail; $i > 0; $i--) {
          $digit = 0;
          $mask = 1 << ($i - 1);
          if (($tileX & $mask) != 0) {
            $digit++;
          }
          if (($tileY & $mask) != 0) {
            $digit++;
            $digit++;
          }
          $quadKey .= $digit;
        }
        return $quadKey;
      }
    
      public static function latLongToQuadKeyDec( $lat, $long, $zoom)
      {
        return self::quadKey4ToDec( self::latLongToQuadKey( $lat, $long, $zoom));
      }
    
      public static function latLongToQuadKey( $lat, $long, $zoom)
      {
        $pixelX = \foci\utils\TileSystem::longToPixelX( $long, $zoom);
        $pixelY = \foci\utils\TileSystem::latToPixelY( $lat, $zoom);
        $tileX = \foci\utils\TileSystem::pixelXToTileX( $pixelX);
        $tileY = \foci\utils\TileSystem::pixelYToTileY( $pixelY);
        return \foci\utils\TileSystem::tileXYToQuadKey( $tileX, $tileY, $zoom);
      }
    
      public static function pixelXYToQuadKey( $pixelX, $pixelY, $zoom)
      {
        $tileX = \foci\utils\TileSystem::pixelXToTileX( $pixelX);
        $tileY = \foci\utils\TileSystem::pixelYToTileY( $pixelY);
        return \foci\utils\TileSystem::tileXYToQuadKey( $tileX, $tileY, $zoom);
      }
    
      public static function pixelXYToQuadKeyDec( $pixelX, $pixelY, $zoom)
      {
        return base_convert( self::pixelXYToQuadKey( $pixelX, $pixelY, $zoom), 4, 10);
      }
    
      public static function quadKeyToTileXY( $quadKey, &$tileX, &$tileY, &$levelOfDetail)
      {
        $tileX = $tileY = 0;
        $levelOfDetail = strlen( $quadKey);
        for( $i = $levelOfDetail; $i > 0; $i--)
        {
          $mask = 1 << ($i - 1);
          switch( $quadKey[$levelOfDetail - $i])
          {
            case '0':
              break;
    
            case '1':
              $tileX |= $mask;
              break;
    
            case '2':
              $tileY |= $mask;
              break;
    
            case '3':
              $tileX |= $mask;
              $tileY |= $mask;
              break;
    
            default:;
          }
        }
      }
    
      public static function tileXYToQuadKeyDec( $tileX, $tileY, $levelOfDetail)
      {
        return base_convert( self::tileXYToQuadKey( $tileX, $tileY, $levelOfDetail), 4, 10);
      }
    
      public static function quadKey4ToDec( $quadkey)
      {
        return base_convert( $quadkey, 4, 10);
      }
    }
    


    Переходим к практике


    Область видимости


    Тут надо учесть важный момент, не нужно запрашивать все кластеры на карте на конкретном zoom, а нужно получать только те, которые входят в область видимости клиента(в нашем случае это android и web).

    Варианты получения кластеров в области видимости:
    • Умный клиент. Клиент сам определяет какой zoom ему нужен и какие кластеры войдут в область его видимости и запрашивает их все сразу например при инициализации view либо запрашиват только те которые вошли в область видимости при перемещении по карте
    • Ленивый клиент. Клиент просто отсылает свою область видимости на сервер и отображет то что будет прислано сервером.

    Мы выбрали Ленивый клиент. Почему?
    • Логика клиента очень проста. Посути клиент работает как view, что получил то и отображает.
    • Это простота удобна если как у нас есть несколько клиентов на разных платформах. Не надо дублировать хитрую логику определения фокусов в области видимости. Вся логика сосредоточина а одном месте — на сервере.
    • Как показала практика требуется малое время на реализацию.

    Какова будет логика сервера?
    • На сервер приходят координаты углов прямоугольной области видимости( например верхний левый угол, нижний правый угол).
    • Сервер зная размеры этой прямоугольной области определяет тайл какого размера в нее входит полностью(и по ширине и по высоте)
    • Дальше сервер зная размер тайла может определить zoom, поэтому его передавать вместе с координатами нет необходимости.
    • Затем надо учесть что этот zoom надо увеличить, для более удобного отображения.
    • И последним сервер опеределит quadkey-и попадающие в область видимости, сделает выборку кластеров и отдаст клиенту.

    Хранение


    Для хранения quadkey-ев нам надо выбрать уровень zoom, возьмем 23, его будет вполне достаточно. БД будем использовать MySQL. Таблица хранения маркеров будет иметь следующую структуру:
    CREATE TABLE `marker` (
      `id` int(11) NOT NULL AUTO_INCREMENT,
      `latitude` FLOAT(19,15) NOT NULL,
      `longitude` FLOAT(19,15) NOT NULL,
      `quadkey` bigint(20) NOT NULL DEFAULT 0,
      PRIMARY KEY (`id`),
      INDEX `quadkey ` (`quadkey `)
    ) ENGINE=InnoDB DEFAULT CHARSET=utf8;
    

    Индекс по quadkey позволяет нам легко делать выборки количества маркеров в тайле(кластере).

    Для хранения quadkey в MySQL можно использовать bigint(20), просто переводим из 4-чной системы в десятичную и в БД и приложении работаем с 10-чным quadkey.

    Код получения кластеров


    Код получения кластеров будет выглядеть примерно так:
      public function getClusters( $quad_key, $zoom)
      {
        $q1 = $quad_key . str_repeat("0", TileSystem::MaxZoom - $zoom);
        $q2 = $quad_key . str_repeat("3", TileSystem::MaxZoom - $zoom);
        $mask = $q1;
        $mask_plus = 1;
        $mask = $quad_key . str_repeat("3", $mask_plus);
        $mask = $mask . str_repeat("0", TileSystem::MaxZoom - $zoom - $mask_plus);
    
        $q1 = TileSystem::quadKey4ToDec( $q1);
        $q2 = TileSystem::quadKey4ToDec( $q2);
        $mask = TileSystem::quadKey4ToDec( $mask);
    
         return $this-> findBetweenQuadKeys( $q1, $q2, $mask);
      }
    
     public function findBetweenQuadKeys( $quad_key_l, $quad_key_r, $mask)
      {
        $entities = [];
    
        $sql = "SELECT `quadkey`, MIN(id) AS id, COUNT(id) AS count, AVG(longitude) AS longitude, AVG(latitude) AS latitude ";
        $sql .= "FROM `marker ` ";
        $sql .= "WHERE `quadkey` BETWEEN $quad_key_l AND $quad_key_r ";
        $sql .= "AND `is_deleted`= '0' ";
        $sql .= "GROUP BY (`quadkey` & $mask) ";
    
        if( $stmt =  $this->_connection->query( $sql))
        {
          foreach( $stmt as $row)
            $entities[] = $this->_createEntityFromRow( $row);
        }
    
        return $entities;
      }
    

    Переменная $mask позволяет нам группировать quadkey-и c определенным zoom. Переменные $quad_key_l и $quad_key_r позволяют учитывать quadkey-и всех маркеров для данного родителя т.е если $quad_key_l=03022313122033033000000 и $quad_key_r=03022313122033033033333 то мы найдем в базе все маркеры для 030223131220330330.
    Так будет выглядеть результирующий sql запрос:
    SELECT `quadkey`, MIN(id) AS id, COUNT(id) AS count, AVG(longitude) AS longitude, AVG(latitude) AS latitude
    FROM marker 
    WHERE quadkey between 15499682906112 and 15499683168255 
    GROUP BY (quadkey & 15499683151872);
    
    -- более понятная 4-чня форма
    -- WHERE quadkey between "03201203031012000000000" and "03201203031012333333333"
    -- GROUP BY (quadkey & "03201203031012330000000")
    


    Что получилось


    К сожалению скриншотов с кластерами на базе в 1 миллион маркеров не сохранилось(хотя это было нечто), но есть текущая небольшая база, и вот как это выглядит:
    web

    Android


    Координаты кластеров находятся в месте наибольшего скопления маркеров в тайле а не в его центре и все благодаря простому усреднению координат в sql запросе.

    Результаты


    Сервиc использующий данную кластеризацию запущен, и работает. Сам подход очень удобен для дальнейшего масштабирования.
    Для избежания постоянного пересчитывания количества маркеров в кластере просто добавлятся кеширование клестера.
    Если рассматривать данную тематику дальше то тайлы активно используются картографическими сервисами для отрисовки карт. Так же что-то подобное используется в разработке игр и т.д.

    Источники


    www.mapbox.com/guides/how-web-maps-work
    msdn.microsoft.com/en-us/library/bb259689.aspx
    Поделиться публикацией

    Комментарии 37

      +3
      1. Конкретно этот quadkey — Morton code.
      2. 16ый уровень разбиения, который можно сохранить в 32 битный uint описывает квадрат со стороной примерно 300 метров. Для предварительной кластеризации этого достаточно.
      3. Ленивый клиент — это не самое умное решение, но за перевод в пиксельные координаты отдельное спасибо.
      4. Маску лучше представлять как маску. Т.е. 03201203031012000000000/12.

      В общем случае эта задача называется «понижением размерности» 2Д данных через использование какой либо spatial filling curve, которая сводит задачу к 1Д варианту. Люблю эту наркоманию.
        +3
        В статье W for Wikipedia, в разделе Z for Z-code есть чуть более техническое описание работы Z кривых, плюс ссылки на sql+php и mongo+js решения для серверной кластеризации.
        Вообще об этой технологии я пытаюсь рассказывать уже лет 8, жаль что до нужных людей мои презентации не доходят :(
          0
          А где можно увидеть Ваши презентации?
            0
            Последняя — events.yandex.ru/lib/talks/2354 — 2014 год, первая нормальная была 2012.codefest.ru/lecture/77 — 2012 год, даже примерчик на гитхаб выложил — github.com/theKashey/quadLoader.
            В принципе «читабельно» в первый раз я описал Z и кластеризацию в статье про gdeetotdom(2010), которая сейчас уехала на geektimes. С тех пор в каждой второй статье на хабре я это дело упоминаю, так и или иначе — уж больно широкие применения у этих чтук.
            Это не просто «число» — одновременно это nested set, B-tree индекс и материализация Linear quad-tree. И это самая простая кривая — их на любой взыскательный вкус напридумывали.
              0
              Да действительно, ваша презентация была как раз одним из источников. Еще на HighLoad 2014 девушка из яндекс карт рассказывала по этой теме.
          0
          Отлично, мого полезной информации. Спасибо. Когда искал что-то подобное найти было сложновато и наткнулся только на обрывки информации.
          А на счет ленивого клиента — это был выбран просто самый быстрый способ реализации. Естественно в дальнейшем от него надо уходить.
            0
            На самом деле, лично использую quad коды в основном для выборки обьектов, а не для кластеризации. Это так — приятный бонус, да и не все кластеризуется. Да и на клиенте можно поработать немного :)
            В частности Z коды используются не только для загрузки данных из базы — именно они (лично у меня) используются для передачи команды на сервер.
            В первый раз такой подход был замечен на wikimapia, я его активно перенял. С тех пор особо популярен он не стал — если и грузят по тайловой сетке, то в обычных x,y,z координатах.
              0
              Нам удалось еще использовать их для посика кластеров с определенными тегами маркеров.
            0
            ммм а чем квадкоды отличаются от геохэшей?
            0
            А почему не используют K-Means кластеризацию? Или это дорого на сервере?

            Это ведь суперлениво для клиента — посылаешь на сервер область видимости и максимальное желаемое количество результирующих кластеров и получаешь ответ.
              0
              Честно говоря в тот момент положились на опыт Яндекса, они в докладах утверждали что рассматривали другие варианты кластеризации но наиболее приемлимой оказалась данная.
                0
                K-means достаточно долгий, в CPU, вариант, но для многих бэкендеров нет ничего проще, чем нагрузить кластер спарков его расчетом…
                Но применимость k-means на картах просто низка. У него есть проблема — он не «стабильный».
                Например вы смотрите на Москву. Запрашиваете данные только для Москвы. Но… вы должны произвести кластеризацию всех обьектов по всему миру на 10ом зуме. Чтобы при сдвиге у вас не изменился набор исходных данных и все «взрывообразно» не перекластеризовалось.
                Потом приходит необходимость делать иеархический k-means, и ты понимаешь что это того не стоит.
                  0
                  Потом приходит необходимость делать иеархический k-means, и ты понимаешь что это того не стоит.

                  сделали, жалеем)
                    0
                    расскажите подробности
                      0
                      Что сказать, медленно работает, агргегация данных для 6000 объектов занимает пару минут. Через k-means определяли кластера элементов для 21 уровня зума, и хоть и удалось уменьшить чуть количество итераций сложность всеравно остается O(N^2). Так как объектов планируется много — не подходит.

                      В итоге перевели агрегацию на ElasticSearch + плагин geocluster-facet, который использует алгоритм с меньшим временем выполнения + оптимизирован за счет отсутствия сложной математики, все сводится к вычислению геохэшей один раз и побитовым операциям, что-то очень напоминающие приведенные в статье квадкоды.
                        0
                        А до какого примерно количества объектов этот метод дает приемлемые результаты?
                          0
                          В смысле? Что считать приемлимым результатом. На больших объемах резальтаты не плохие, но по времени кластеризации выходит накладно, и алгоритм не позволяет перестраивать данные только частично. По точности проблем особо небыло, проблемы у него те же что и у остальных (пересчет пиксальных расстояний в физические координаты), так что тут все весьма стандартно. Но выборка по фасетам выдает примерно схожий по точности результат за гораздо меньшее время. Ну то есть по фасетам точность выборки самостовима с k-means, выше чем у гридов.
                            0
                            Я имел в виду именно время. На каких примерно объемах оно ещё работает достаточно быстро.

                            Кстати, насчет пересчета координат пиксели/градусы — это да. Не могу до конца понять, как это сделать правильнее всего — мешают нелинейности проекции Мерктора. Если смотреть в масштабе одного города, то ими можно пренебречь, но если целый континент — уже нет.
                              0
                              Скажем так, при том количестве айтемов которое дает адекватное время обработки это все можно и на клиенте кластеризовать. Потому лучше взять что-то более интересное.
                0
                Не понятно какая задача решается:
                если нужно правильно хранить, то лучше взять хранилище с поддержкой геометрических индексов. Они есть уже и в SQL и в noSQL решениях. C MySQL не работал, поэтому не могу сказать, что есть в нём, но в среде открытых есть PostGIS.
                если нужна именно кластеризация, то грид быстрый, и самый не красивый из всех алгоритмов. Когда данных много, пины кластеров образуют практически регулярную сетку. Это хорошо заметно на скрине Яндекса в статье. Раз кластеры считаются на сервере, то можно выбрать более ресурсоёмкий алгоритм с красивым результатом.
                Если посмотреть на скрины того «что получилось»: алгоритм выбора репрезентативной точки выбран не удачно. Например, маркер который, видимо, должен принадлежать Кубе просто висит над морем.
                  0
                  Задача решалась реализовать выбранный способ кластеризации.
                  Да, в других БД есть поддержка хранения гео данных, но вот на счет кластеризации я не уверен. И здесь просто вырезки из большого проекта где используется MySQL.
                  Естественно что алгоритм не идеален, понятно что мы находим средную координату в квадрате и она может посать и на вершину горы и в другю страну и в озеро и в море. Все эти недостатки можно устранить в будущем.
                    0
                    Самый простой вариант — под конечной координатой выбирайте некую реальную координату в бд, близкую к средней.
                      0
                      Да, это очень простой вариант, так мы всегда попадем в реальную точку.
                  0
                  Если не секрет… Где вы использовали данную технологию? Какие задачи решали?
                    0
                    Использовали в одном гео-сервисе, название сказать не могу(посчитают пиаром).
                    Задача была размещать на карте большое количество маркеров( миллионы), и отображать это пользователям без проблем с производительностью.
                      0
                      Не посчитают, не боись.
                      Так где же?
                        +1
                        Demo версия карты fociapp.com/map. Там же есть ссылка на приложение. А вообще планируетс статья по этому проекту в разделе «Я пиаруюсь».
                          0
                          Немного хромает реализация — тут и метки «мигают» при сдвиге карты и кластеризация немного подтупливает.
                          Я бы с сервера отдавал бы больше точек — еще бы на один уровень маску группировки сдвинул бы и «докластеризовывал» бы на клиенте.
                          Ну и лимиты при выборке из базы увеличил бы (или совсем убрал) — количество данных регулируется кластеризацией.
                          Сейчас реально бывают «дырки» — в европе данные есть, в Австралии — нет. Убираем европу, и внезапно точки появляются. От этого еще хорошо разбиение запроса на тайлы лечит.
                            0
                            Согласен, уже обдумывал эти правки, но пока немного не хватает времени. Сама клиентская часть еще достаточно срая.
                    0
                    Очень хочу увидеть сравнение этого приложения с Postgresql/PostGIS и идекса типа R-tree. Что быстрее?

                    Просто в MySQL вся поддержка геоданных сделана по остаточному принципу и to be implemented. В PostGIS она полноценная и работает очень хорошо.
                      0
                      Сделать сравнение было бы доастаочно интересно, попробую найти на это время :)
                      В MySQL рассматривал Spatial Data, но там насколько помню они еще не готовы к полноценному использованию.
                        0
                        Вот я в 2009 их пробовал, тоже не были готовы. :)
                        0
                        Из опробованных решений еще стоит упоминуть ElasticSearch и плагин geocluster-facet. Ну и эластика умеет агрегацию по сетке делать из коробки (использую геохэши).
                        0
                        Кстати на практике оказалось что использование 50% медианы куда лучше выглядит чем среднее арифметическое.
                          0
                          Возможно, надо будет попробовать.
                          0
                          Сделали кластеризацию этим же методом bookingtroll.com
                          В процессе реализации возникло ряд проблем: на стыках сетки (особенно экватор, 180 меридиан, 41-й меридиан) работает не очень хорошо.
                          Если область карты далека от квадрата, то лучше разбивать ее на два или больше квадратных кластеров и получать данные для них отдельно.
                            0
                            Да, с такими проблемама тоже сталкивался, и решение прирено такое же)

                          Только полноправные пользователи могут оставлять комментарии. Войдите, пожалуйста.

                          Самое читаемое