impoly only approximately correct on log-scale axes

140 views Asked by At

When determining points within polygons using MATLAB's inpolygon function, I find that the results are exactly correct for polygons drawn on linear axes but only approximately correct for polygons drawn on log-scale axes. Although my suspicions lean in favor of a MATLAB bug, it's possible I've overlooked something.

The following code reproduces the issue I have been experiencing with other data. The results are shown in the following image (the bottom set of panels are zoomed views of the top panels). One can appreciate that there are unlabeled points inside the polygon and labeled points outside the polygon, neither of which should occur, in the case of a polygon drawn on log-scale axes (right). In contrast, the polygon test is exact for polygons drawn on linear axes (left).

n=2E4;
x(:,1)=rand(n,1); y(:,1)=rand(n,1);
x(:,2)=lognrnd(0.5,0.25,n,1); y(:,2)=lognrnd(0.5,0.25,n,1);
for m=1:2
    subplot(1,2,m);
    scatter(x(:,m),y(:,m),'.'); hold on;
    if(m==2)
        set(gca,'xscale','log'); set(gca,'yscale','log');
    end
    p=impoly(gca);
    pc=getPosition(p);
    in=inpolygon(x(:,m),y(:,m),pc(:,1),pc(:,2));
    scatter(x(in,m),y(in,m),20);
end

enter image description here

1

There are 1 answers

2
Ander Biguri On BEST ANSWER

I think you missed something: A line in normal scale is not a line in log scale. Your polygons are not properly drawn in the log scale, as you draw 2 points and put them together with a straight line.

Look at the real polygon in log space:

close all
clear

n=2e4;
x(:,1)=rand(n,1); y(:,1)=rand(n,1);
x(:,2)=lognrnd(0.5,0.25,n,1); y(:,2)=lognrnd(0.5,0.25,n,1);
for m=1:2
    subplot(1,2,m);
    scatter(x(:,m),y(:,m),'.'); hold on;
    if(m==2)
        set(gca,'xscale','log'); set(gca,'yscale','log');
    end
    p=impoly(gca);
    pc=getPosition(p);
    % plot polygon
    hold on
    for ii=1:size(pc,1)-1
        plot(linspace(pc(ii,1),pc(ii+1,1),100),linspace(pc(ii,2),pc(ii+1,2),100),'g')
    end
    plot(linspace(pc(end,1),pc(1,1),100),linspace(pc(end,2),pc(1,2),100),'g')
    in=inpolygon(x(:,m),y(:,m),pc(:,1),pc(:,2));
    scatter(x(in,m),y(in,m),20);
end

Look at this zoomed in result (click to enlarge):

enter image description here

This happens because the polygon is defined in euclidean space, and it is defined as points linked by lines. If you want to work in log space, things may get complicated. One way to numerically approximate it is the inverse of what I did for plotting. Create dense enough sampled straight line on log space, convert it to linear space, and define a high vertex polygon with the resulting points. Then use inpolygon.