# Is it possible to calculate the line of sight between two points in a 2d polygon?

22 views (last 30 days)
Right Grievous on 29 Oct 2014
Edited: Matt J on 17 Oct 2018
Hi everybody,
I have a polygon defined by a number of vertex coordinates and I have a coordinate within this polygon. I would like to know if a point A is directly visible from point B (i.e. I need to know if there are polygon walls between the two points).
I can probably do this with for loops and Matlab's polyxpoly function but this function seems a bit slow and ultimately I want to run this calculation on thousands of points. The function also often gives two points of intersection when there is only one (because I think it gives the start and end point of a line segment)?
Is there any inbuilt function that can work out line of sight in 2d? I know los2 can do this in 3d with elevation maps, but I'm unsure if it can be used in 2d only.
Any help would be greatly appreciated,
Rod.
p.s. I have all toolboxes.

#### 1 Comment

Right Grievous on 29 Oct 2014
A piece of code to visualise what I'm trying to do:
environment = [-100,-100; 100,-100; 100,100; -100,100; -100,-100];
pointA = [20,20];
pointB = [-100,50];
figure
impoly(gca,environment)
hold on
plot(pointA(1),pointA(2),'ko')
plot(pointB(1),pointB(2),'ro')
plot([pointA(1) pointB(1)],[pointA(2) pointB(2)]);
axis([-110,110,-110,110])

Kelly Kearney on 29 Oct 2014
Edited: Kelly Kearney on 29 Oct 2014
Have you tried taking advantage of the fact that polyxpoly can process NaN-separated polylines? It seems to handle 1000 line segments pretty quickly... I made it to 15000 segments before I crossed the 1 sec threshold.
environment = [-100,-100; 100,-100; 100,100; -100,100; -100,-100];
pointA = [20,20];
pointB = [-100,50];
xp = environment(:,1);
yp = environment(:,2);
ntarget = 1000;
x1 = randn(ntarget,1)*100;
y1 = randn(ntarget,1)*100;
x2 = randn(ntarget,1)*100;
y2 = randn(ntarget,1)*100;
xlos = [x1 x2 nan(ntarget,1)]';
ylos = [y1 y2 nan(ntarget,1)]';
[xi, yi, ii] = polyxpoly(xp, yp, xlos(:), ylos(:));
[ir,iseg] = ind2sub(size(xlos), ii(:,2));
isseen = true(ntarget,1);
isseen(iseg) = false;
figure;
hln = plot(xlos, ylos, 'r');
hold on;
set(hln(isseen), 'color', 'g');
hp = plot(xp, yp, 'b');

William Chamberlain on 17 Oct 2018
Kelly's answer is better , but in case you're after something which operates on grids/maps/images, here's a couple of Bresenham's line algorithm implementations:

Matt J on 17 Oct 2018
Edited: Matt J on 17 Oct 2018
I'm not sure if it would meet your needs speed-wise, but another option would be to use intersectionHull in this File Exchange package. You would use it to test whether the line AB has a non-empty intersection with the polygon.
I = intersectionHull('vert', polygonVertices, 'vert', [A(:),B(:)].')