qgis-3 – Calculating Angle of View from Points to Polygon in QGIS

Question

This question Calculate field of view / angle of a point to a shapefile was attempting to achieve a similar outcome to my need (calculating the angle of view between observer locations and a polygon). The image from the question is shown below for ease, with the angle of interest being alpha:

My solution is to use a virtual layer to calculate the bearing between each viewpoint and each of the nodes on the polygon layer and then use aggregate expressions to give me the information I need.

I have the following code at present but testing it throws a syntax error (1 – near ".": syntax error) that I cannot work out (there may be other errors and I'm not convinced it's the most optimal code but I would like to use a virtual layer).

SELECT 
  vp_union.*,
  minimum(vp_union.bearing) AS bearing_min,
  median(vp_union.bearing) AS bearing_med,
  maximum(vp_union.bearing) AS bearing_max,
  range(vp_union.bearing) AS bearing_range,
  distance(
    vp_union.geometry,
    closest_point(
      st_union(redline_outer.geometry),
      vp_union.geometry
    ),
    vp_union.geometry
  ) AS distance
FROM
  (SELECT
    vp.*,
    degrees(
      azimuth(
        vp.geometry,
        redline_nodes.geometry
      )
    ) AS vp.bearing
  FROM
    "Viewpoint Layer 2D" AS vp
  UNION
    (SELECT
      nodes_to_points(
        st_union(redline.geometry)
      ) AS geometry
    FROM
      "Redline" AS redline
    ) AS redline_nodes
  ) AS vp_union,
  "Redline" AS redline_outer
GROUP BY
  vp_union.name

Answer 1

So I was well out on my original idea. In the end I realised that there is a possible case where a viewpoint could lie within the convex hull of the polygons, so I went back to testing the individual vertices. It's long, convoluted and a bit slow (I'm sure someone could optimise).

I'm aware this doesn't quite give you alpha (but it did what I needed it too) and it's relatively trivial to calculate that angle using the start and end points of the linestring for each viewpoint.

WITH RECURSIVE
  -- generate table with all the points from the polygons layer
  polygons_points AS(
  SELECT
    1 AS fid,
    1 AS n,
    point_n(polygons.geometry,1) AS geometry
  FROM
    "Redline" AS polygons
  UNION ALL
  SELECT
    ogc_fid AS fid,
    n + 1 AS n,
    point_n(polygons.geometry,n+1)
  FROM
    "Redline" AS polygons
  INNER JOIN
    polygons_points ON polygons.ogc_fid = polygons_points.fid
  WHERE
    n < ST_NPoints(polygons.geometry)
  ),
  
  -- generate table with spokes from polygons_points to the points layer
  spokes AS(
  SELECT
    points.VP AS pid,
    polygons_points.fid,
    polygons_points.n,
    MakeLine(points.geometry,polygons_points.geometry) AS geometry
  FROM
    "Viewpoint Locations 2D" AS points
  CROSS JOIN
    polygons_points
  ),
  -- select sightlines from spokes that do not cross the original polygons
  sightlines AS(
  SELECT spokes.*
  FROM
    spokes
  INNER JOIN
    "Redline" AS polygons
  ON
    polygons.ogc_fid = spokes.fid
  AND
    ST_Crosses(spokes.geometry,polygons.geometry) <> 1
  ),
  -- select and join the polygon segments that intersect with the sightlines
  polygons_segments AS(
  SELECT
    sightlines.pid,
    sightlines.fid,
    sightlines.n,
    ST_AsText(sightlines.geometry) AS point_txt,
    ST_Union(ST_GeometryN(ST_DissolveSegments(ST_ExteriorRing(polygons.geometry)),n)) AS geometry
  FROM
    sightlines
  INNER JOIN
    "Redline" AS polygons
  ON
    polygons.ogc_fid = sightlines.fid
  AND
    ST_Intersects(sightlines.geometry, polygons.geometry)
  GROUP BY
    sightlines.pid
  )

SELECT * FROM polygons_segments