Intersect or not Disjoint

Продолжаю работать над Mapfeatureserver.

В PostGIS есть интересный оператор «&&», вот его описание:

A && B

The "&&" operator is the "overlaps" operator. If A's bounding box overlaps B's bounding box the operator returns true.

http://postgis.net/docs/manual-1.3/ch06.html

Обратите внимание — рассматриваются охватывающие боксы фигур. Это очень удобный оператор в плане быстродействия, для поиска задействуется пространственный индекс и более ничего. Но в плане корректности результата это уже не так хорошо.

Вот пример вычисления ответа на REST запрос к featureserver ArcGIS:

запрос http://vags101.algis.com/arcgis/rest/services/seisprofiles/FeatureServer/1/query?returnGeometry=true&geometryType=esriGeometryEnvelope&geometry={"xmin":6593024.93074,"ymin":10538006.4152,"xmax":6641944.62884,"ymax":10569039.8487,"spatialReference":{"wkid":102100}}&inSR=102100&spatialRel=esriSpatialRelIntersects&outSR=102100&outFields=*&f=pjson

говорит нам о том, что есть прямоугольник для которого надо найти все пересекающиеся сейсмопрофили. Эталонный ответ содержит 16 записей.

Запрос к PostGIS с оператором «&&» вернет 20 записей

select *, st_asgeojson(1, st_transform(geom::geometry, 3857)) shape
from seisprof
where geom::geometry && ST_transform(
        ST_GeomFromText(
            'POLYGON((6593024.93074 10538006.4152, 6641944.62884 10538006.4152, 6641944.62884 10569039.8487, 6593024.93074 10569039.8487, 6593024.93074 10538006.4152))'
            , 3857)
        , st_srid(geom::geometry) );

то есть найдет дополнительно 4 записи, боксы которых перекрывают заданный бокс, но сами сейсмопрофили находятся снаружи заданного прямоугольника. Зато очень быстро.

Правильный (зато более тормозной) запрос к PostGIS, возвращающий нужные 16 записей выглядит так

select *, st_asgeojson(1, st_transform(geom::geometry, 3857)) shape
from seisprof
where not st_disjoint(
    geom::geometry
    ,ST_transform(
        ST_GeomFromText(
            'POLYGON((6593024.93074 10538006.4152, 6641944.62884 10538006.4152, 6641944.62884 10569039.8487, 6593024.93074 10569039.8487, 6593024.93074 10538006.4152))'
            , 3857)
        , st_srid(geom::geometry) ) );

причем «not st_disjoint» можно заменить на «st_intersects», что одно и то же

http://ru.convdocs.org/docs/index-9715.html?page=5#106030
http://postgis.refractions.net/docs/reference.html#Spatial_Relationships_Measurements

Конечно, нельзя не упомянуть тот факт, что разница в быстродействии будет заметна только в том случае, если правильному запросу придется отсекать гораздо больше лишних записей, нежели 4.

Чуть не забыл, оператор PostGIS «&&» более всего похож на «esriSpatialRelEnvelopeIntersects» в ArcGIS.

original post http://vasnake.blogspot.com/2013/04/intersect-or-not-disjoint.html