ST_HexagonGrid几何矢量查找所有点 [英] ST_HexagonGrid geom vector to find all points

查看:130
本文介绍了ST_HexagonGrid几何矢量查找所有点的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在审查PostGis中的此功能

I was reviewing this function in PostGis

https://postgis.net/docs/manual-dev/ST_HexagonGrid.html

1)我不理解的是潜在的几何数据.如图所示,获得美国地图的来源是什么?什么是数据库架构?我想这可能是一条记录,如果我只需要美国边界,而不是每个州呢?

1) What I don't understand is what the underlying geom data would be. What's the source to get the USA map as shown? What is the DB schema? I think it could be one record, if I only need the USA boundary, not each state?

2)结果是否为点列表?或几何向量?

2) Is the result a list of points? or geom vectors?

3)如果是geom向量,如何将它们转换为lat和lng点?

3) If geom vectors, how do you convert them into points of lat and lng?

4)如何近似六边形,从一个点开始大约50英里半径?

4) How to do approximate the hexigons to approximate a 50 mile radius from a point?

更新:

根据下面的吉姆·琼斯示例,我尝试使用宽度尝试获得正确的六边形数量.不幸的是,出事了..

I played with the width to try to get the correct number of hexagons based on Jim Jones example below. Unfortunately, something went wrong..

1)长度似乎与米无关

1) the length seems to have no relation to meters

2)有多个大小的六边形,看起来很奇怪.

2) There are multiple sized hexagons, which seems weird.

postgis_test=# WITH j AS (
postgis_test(# SELECT ST_Transform((hex).geom,4326) AS hex FROM ( 
postgis_test(#   SELECT 
postgis_test(#   generate_hexgrid(
postgis_test(#     5909968.8,
postgis_test(#     ST_XMin(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_YMin(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_XMax(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_YMax(ST_Extent(ST_Transform(geom,3857))) ) AS hex
postgis_test(# FROM usa_states)i) 
postgis_test-# SELECT count(j.hex) FROM j,usa_states
postgis_test-# WHERE ST_Intersects(usa_states.geom,j.hex);
 count 
-------
   119
(1 row)

postgis_test=# WITH j AS (
postgis_test(# SELECT ST_Transform((hex).geom,4326) AS hex FROM ( 
postgis_test(#   SELECT 
postgis_test(#   generate_hexgrid(
postgis_test(#     5909968.8,
postgis_test(#     ST_XMin(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_YMin(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_XMax(ST_Extent(ST_Transform(geom,3857))) ,
postgis_test(#     ST_YMax(ST_Extent(ST_Transform(geom,3857))) ) AS hex
postgis_test(# FROM usa_states)i) 
postgis_test-# SELECT DISTINCT st_area(j.hex) FROM j,usa_states
postgis_test-# WHERE ST_Intersects(usa_states.geom,j.hex);
     st_area      
------------------
 1219.78281686003
 2089.11341619338
 2089.11341619338
 3379.93344444246
  7051.4650344734
 12076.9943663072
(6 rows)

推荐答案

根据作者 >,以下函数应根据给定的BBOX和以 meters 为单位的像元大小来创建具有范围的网格.

According to the author, the following function should create a grid with the extent based on the given BBOX and the cell size in meters.

SRID 3857单位是[非常大约]米,并使用此单位 投影将创建在网络地图上看起来正确"的十六进制单元格(大多数 其中使用网络墨卡托投影).

SRID 3857 units are [very approximately] meters, and using this projection will create hex cells that "look right" on a web map (most of which use a web mercator projection).

CREATE OR REPLACE FUNCTION generate_hexgrid(width float, xmin float, ymin float, xmax float, ymax float, srid int default 3857)
RETURNS TABLE(gid text,geom geometry(Polygon)) AS $$
DECLARE
  b float := width / 2;
  a float := tan(radians(30)) * b;
  c float := 2 * a;
  height float := 2 * (a + c);
  index_xmin int := floor(xmin / width);
  index_ymin int := floor(ymin / height);
  index_xmax int := ceil(xmax / width);
  index_ymax int := ceil(ymax / height);
  snap_xmin float := index_xmin * width;
  snap_ymin float := index_ymin * height;
  snap_xmax float := index_xmax * width;
  snap_ymax float := index_ymax * height;
  ncol int := abs(index_xmax - index_xmin);
  nrow int := abs(index_ymax - index_ymin);
  polygon_string varchar := 
    'POLYGON((' || 0 || ' ' || 0 || ' , ' || b || ' ' || a || ' , ' ||
    b || ' ' || a + c || ' , ' || 0 || ' ' || a + c + a || ' , ' ||
    -1 * b || ' ' || a + c || ' , ' || -1 * b || ' ' || a || ' , ' ||
    0 || ' ' || 0 ||'))';
BEGIN
  RETURN QUERY
  SELECT 
    format('%s %s %s', width,
    x_offset + (1 * x_series + index_xmin),
    y_offset + (2 * y_series + index_ymin)),
    ST_SetSRID(ST_Translate(two_hex.geom,
    x_series * width + snap_xmin,
    y_series * height + snap_ymin), srid)
  FROM  generate_series(0, ncol, 1) AS x_series,
        generate_series(0, nrow, 1) AS y_series,
    (SELECT 0 AS x_offset, 0 AS y_offset, polygon_string::geometry AS geom
     UNION
     SELECT 0 AS x_offset, 1 AS y_offset, ST_Translate(polygon_string::geometry, b , a + c)  AS geom
    ) AS two_hex;
END; $$ LANGUAGE plpgsql;

考虑到您有一个名为usa的表,它包含此 shapefile 您应该能够使用以下查询创建与美国地图重叠的网格:

Considering you have a table called usa and it contains the geometries of this shapefile you should be able to create a grid that overlaps the USA map with the following query:

CREATE TABLE usa_hex AS
WITH j AS (
SELECT ST_Transform((hex).geom,4326) AS hex FROM ( 
  SELECT 
  generate_hexgrid(
    80467,
    ST_XMin(ST_Extent(ST_Transform(geom,3857))) ,
    ST_YMin(ST_Extent(ST_Transform(geom,3857))) ,
    ST_XMax(ST_Extent(ST_Transform(geom,3857))) ,
    ST_YMax(ST_Extent(ST_Transform(geom,3857))) ) AS hex
FROM usa)i) 
SELECT j.hex FROM j,usa
WHERE ST_Intersects(usa.geom,j.hex);

编辑:它仍然不能解决问题,因为它没有使用仪表创建六边形,但可能会对其他用户有所帮助.以下功能(源自此 answer )以创建尺寸完全相同的 geometry type 六边形.

EDIT: It is still not an answer, as it does not create the hexagons using meters, but it might help other users. The following function (derived from this answer) creates geometry type hexagons of exact the same size in degrees.

CREATE OR REPLACE FUNCTION generate_hexagons(width FLOAT, bbox BOX2D, srid INTEGER DEFAULT 4326)
RETURNS TABLE (gid INTEGER, hexagon GEOMETRY) AS $$
DECLARE
  b FLOAT := width/2;
  a FLOAT := b/2;
  c FLOAT := 2*a;
  height FLOAT := 2*a+c;
  ncol FLOAT := ceil(abs(ST_Xmax(bbox)-ST_Xmin(bbox))/width);
  nrow FLOAT := ceil(abs(ST_Ymax(bbox)-ST_Ymin(bbox))/height);
  polygon_string VARCHAR := 'POLYGON((' || 
    0 || ' ' || 0 || ' , ' || b || ' ' || a || ' , ' || b || ' ' || a+c || ' , ' || 0 || ' ' || a+c+a || ' , ' ||
   -1*b || ' ' || a+c || ' , ' || -1*b || ' ' || a || ' , ' || 0 || ' ' || 0 || '))';
BEGIN    
  RETURN QUERY 
  SELECT 
    row_number() OVER ()::INTEGER,
    ST_SetSRID(
      ST_Translate(geom, x_series*(2*a+c)+ST_Xmin(bbox), y_series*(2*(c+a))+ST_Ymin(bbox)),srid)
  FROM generate_series(0, ncol::INTEGER, 1) AS x_series,
       generate_series(0, nrow::INTEGER,1 ) AS y_series,
       (SELECT polygon_string::GEOMETRY AS geom
        UNION
        SELECT ST_Translate(polygon_string::GEOMETRY, b, a + c) AS geom) AS two_hex;    
END;
$$ LANGUAGE plpgsql;

与上面使用的数据集重叠:

Overlapping with the data set used above:

WITH j (hex_rec) AS (
  SELECT generate_hexagons(3.0,ST_Extent(geom)) 
  FROM usa
)
SELECT (hex_rec).gid,(hex_rec).hexagon FROM j, usa 
WHERE ST_Intersects(usa.geom,(hex_rec).hexagon);

进一步阅读:

这篇关于ST_HexagonGrid几何矢量查找所有点的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆