提问者:小点点

我可以使用Postgres函数在固定大小的旋转矩形内查找点吗?


我正在使用Postgres 9.5,我刚刚为一些扩展功能安装了PostGIS。我有一个带有(x, y)点的表格,我想找到适合最大点数的矩形。约束是矩形边长是固定的。到目前为止,我正在计算盒子里有多少点没有旋转。我的点以原点(0,0)为中心。

SELECT Sum(CASE
             WHEN x > -5
                  AND x < 5
                  AND y > -10
                  AND y < 10 THEN 1
             ELSE 0
           END) AS inside_points,
       Count(1) AS total_points
FROM   track_t;  

这个查询给了我一个原点(0,0)和长度x=10和y=20的矩形内的点数。

从这里,我将创建一个旋转矩形角点的辅助表(角度,x1,y1,x2,y2),然后交叉连接到我的数据,并计算每个角度的点,而GROUP BY角度。然后我可以选择哪个角度给了我矩形内最多的点。

但是这似乎有点过时,也许是不执行的。此外,在旋转的矩形中计数点不是一个简单的计算。

有没有更有效和优雅的方法,也许使用Postgres几何数据类型或PostGIS Box2D,旋转具有固定边长的矩形,然后计算内部的点数?几何函数看起来不错,但它们似乎提供了最小边界框,而不是相反。

除了Postgreql,我还使用了一个Python框架,以防SQL无法完成这项工作。

更新:我尝试的一件事是使用几何类型,特别是BOX

  SELECT deg, Box(Point(-5, -10), Point(5, 10)) * Point(1, Radians(deg)) 
        FROM   Generate_series(0, 360, 90) AS deg

不幸的是,点的旋转函数不适用于多边形。


共1个答案

匿名用户

我最终生成了矩形顶点,旋转这些顶点,然后将矩形的面积(常数)与包含测试点的4个三角形的面积进行比较。

这种技术是基于吝啬的答案:

制作三角形。假设,abcd是矩形,x是点,那么如果区域(abx)区域(bcx)区域(cdx)区域(dax)等于区域(abcd)那么点在它里面。

矩形定义为

>

  • 左下角(-x/2,-y/2)

    b左上角(-x/2, y/2)

    C右上角(x/2, y/2)

    D右下角(x/2,-y/2)

    然后,此代码检查point(qx, qy)是否位于宽度x=10和高度y=20的矩形内,该矩形围绕原点(0,0)旋转一个角度,范围为0到180,10度。

    这是代码。检查750k点需要9分钟,所以有一定的改进空间。此外,一旦我升级到9.6,它就可以并行化

    with t as (select 10*0.5 as x, 20*0.5 as y, 17.0 as qx, -3.0 as qy)
    
    select 
        z.angle
        -- ABC area
        --,abs(0.5*(z.ax*(z.by-z.cy)+z.bx*(z.cy-z.ay)+z.cx*(z.ay-z.by)))
    
        -- CDA area
        --,abs(0.5*(z.cx*(z.dy-z.ay)+z.dx*(z.ay-z.cy)+z.ax*(z.cy-z.dy)))
    
        -- ABCD area
        ,abs(0.5*(z.ax*(z.by-z.cy)+z.bx*(z.cy-z.ay)+z.cx*(z.ay-z.by))) + abs(0.5*(z.cx*(z.dy-z.ay)+z.dx*(z.ay-z.cy)+z.ax*(z.cy-z.dy))) as abcd_area
    
        -- ABQ area
        --,abs(0.5*(z.ax*(z.by-z.qx)+z.bx*(z.qy-z.ay)+z.qx*(z.ay-z.by)))
    
        -- BCQ area
        --,abs(0.5*(z.bx*(z.cy-z.qx)+z.cx*(z.qy-z.by)+z.qx*(z.by-z.cy)))
    
        -- CDQ area
        --,abs(0.5*(z.cx*(z.dy-z.qx)+z.dx*(z.qy-z.cy)+z.qx*(z.cy-z.dy)))
    
        -- DAQ area
        --,abs(0.5*(z.dx*(z.ay-z.qx)+z.ax*(z.qy-z.dy)+z.qx*(z.dy-z.ay)))
    
        -- total area of triangles with question point (ABQ + BCQ + CDQ + DAQ)
        ,abs(0.5*(z.ax*(z.by-z.qx)+z.bx*(z.qy-z.ay)+z.qx*(z.ay-z.by)))
            + abs(0.5*(z.bx*(z.cy-z.qx)+z.cx*(z.qy-z.by)+z.qx*(z.by-z.cy)))
            + abs(0.5*(z.cx*(z.dy-z.qx)+z.dx*(z.qy-z.cy)+z.qx*(z.cy-z.dy)))
            + abs(0.5*(z.dx*(z.ay-z.qx)+z.ax*(z.qy-z.dy)+z.qx*(z.dy-z.ay))) as point_area
    
    from
    (
    SELECT 
        a.id as angle
        -- bottom left (A)
        ,(-t.x) * cos(radians(a.id)) - (-t.y) * sin(radians(a.id)) as ax
        ,(-t.x) * sin(radians(a.id)) + (-t.y) * cos(radians(a.id)) as ay
        --top left (B)
        ,(-t.x) * cos(radians(a.id)) - (t.y) * sin(radians(a.id)) as bx
        ,(-t.x) * sin(radians(a.id)) + (t.y) * cos(radians(a.id)) as by
        --top right (C)
        ,(t.x) * cos(radians(a.id)) - (t.y) * sin(radians(a.id)) as cx
        ,(t.x) * sin(radians(a.id)) + (t.y) * cos(radians(a.id)) as cy
        --bottom right (D)
        ,(t.x) * cos(radians(a.id)) - (-t.y) * sin(radians(a.id)) as dx
        ,(t.x) * sin(radians(a.id)) + (-t.y) * cos(radians(a.id)) as dy
    
        -- point to check (Q)
        ,t.qx as qx
        ,t.qy as qy
    FROM generate_series(0,180,10) AS a(id), t
    ) z
    ;
    

    结果是

    angle;abcd_area;point_area
    0;200;340
    10;200;360.6646055963
    20;200;373.409049054212
    30;200;377.846096908265
    40;200;373.84093170467
    50;200;361.515248361426
    60;200;341.243556529821
    70;200;313.641801308188
    80;200;279.548648061772
    90;200;240
    *100;200;200*
    *110;200;200*
    *120;200;200*
    *130;200;200*
    *140;200;200*
    150;200;237.846096908265
    160;200;277.643408923024
    170;200;312.04311584956
    180;200;340
    

    其中角100、110、120、130和140度的旋转包括测试点(用*表示)