inpolygon with alphashape

I am trying to write some code that will determine whether a set of data points is within a polygon that is formed as a convex hull around the points. To do this I am using the inpolygon function. The code that I have so far is based around the alphavol code available on the file exchange.

    x = [0.285946856420363;0.543663232841362;0.984776236970610;0.715678067126622;0.838969597330000;0.433260561007948;0.470624715606541;0.560713411044556;0.269091543530998;0.749018467641946;0.503887773460489;0.646809666272781;0.307745582030941;0.138724636078145;0.475572933797563;0.362459280818746;0.788113428398547;0.780295820844110;0.668512214253545;0.133503859661312;0.0215558872034970;0.559840705872510;0.300819018069489;0.939409713873458;0.980903636046859;0.286620388894259;0.800820286951535;0.896111351432604;0.597526576817830;0.884016735723820;0.943731541195791;0.549158087419903;0.728386824594357;0.576758297858010;0.0258574710831396;0.446530978284814;0.646301957350656;0.521202952672418;0.372312660779512;0.937134666341562;0.829532824526515;0.849085479954455;0.372534239899540;0.593184575218490;0.872552564647559;0.933501608507105;0.668464274361376;0.206776457935105;0.653850592062588;0.0720515512206499;];
    y = [0.406726915089364;0.666931533207499;0.933725659545930;0.810950032238265;0.484548271834990;0.756749210065512;0.417047453742767;0.971785992989294;0.987974701230832;0.864147529031203;0.388883775912429;0.454741828039112;0.246687197638079;0.784423093024177;0.882837605583053;0.913711681293018;0.558284923622275;0.598868102746205;0.148876720337549;0.899713484535739;0.450393580652116;0.205672339463963;0.899650990509081;0.762585539373227;0.882486307400298;0.284950218024821;0.673225986452258;0.664279904418051;0.122814993899526;0.407318423056813;0.275286951315690;0.716669740639142;0.283384381628936;0.896198856849495;0.826578892313560;0.390026512813998;0.497902942964187;0.694805192672494;0.834369001099938;0.609629689675732;0.574737160856408;0.326042170922791;0.456424600851747;0.713795583233135;0.884405045275995;0.720855670816932;0.0186127747263861;0.674776467128286;0.438508824236513;0.437820179166270;];

    X=[x,y];

    % Dimension
    dim = size(X,2);
    if dim ~= 2;
        error('alphavol:dimension','X must have 2 columns.')
    end

    A=pdist(X);
    A=squareform(A);

    % for i=2:size(A,1)
    %     p_min(i,:)=min(A(i,1:i-1));
    % end

    for i=1:size(A,1)
        j = A(i,:) > 0;
        p_min(i,:)=min(A(i,j));    
    end

    r=max(p_min);

    T = delaunayn(X,{'Qt','Qbb','Qc','Qz'});
    % Limit circumradius of primitives
    dt = TriRep(T,X);
    [~,rcc] = circumcenters(dt);
    T = T(rcc < r,:);

    % Volume

    % Empty triangulation
    if isempty(T)
        V = [];
        return
    end

    % Local coordinates
    A = X(T(:,1),:);
    B = X(T(:,2),:) - A;
    C = X(T(:,3),:) - A;

    % 2D Area
    V = B(:,1).*C(:,2) - B(:,2).*C(:,1);
    V = abs(V)/2;

    % Plot 2D alpha shape

    % Remove inner edges
    E = [T(:,2:3); T(:,[1 3]); T(:,1:2)];
    E = sort(E,2);
    E = sortrows(E);
    ii = find(~any(diff(E),2));
    E([ii;ii+1],:) = [];

    % Coordinates
    x = X(:,1);
    y = X(:,2);
    xv = x(E)';
    yv = y(E)';

    % Plot point set and boundary edges
    plot(x,y,'k.',xv,yv,'b',xv,yv,'r.')

The x and y dataset that I have included is one which I know the alphavol code does not create a convex hull around (which contains every point). This section of code plots the original dataset and the alphashape around the points, which is fine. I am now trying to take the coordinates of the polygon formed and use them with the inpolygon function. The code below is what I have written to try and achieve this.

    % Put polygon coordinates in correct form
    xvv = (xv)'
    yvv = (yv)'

    xvv1 = xvv(:,1)
    yvv1 = yvv(:,1)
    xvv2 = xvv(:,2)
    yvv2 = yvv(:,2)

    xv1 = [xvv1 ; xvv2]
    yv1 = [yvv1 ; yvv2]

    % inpolygon
    in = inpolygon(x,y,xv1,yv1)
    plot(xv1,yv1,x(in),y(in),'.r',x(~in),y(~in),'.b');

The code that I have written does not work as it does not keep the coordinates of the polygon vertices in the correct order. How can I take the polygon coordinates as produced in the first chunk of code and use these with the inpolygon function?

Answer 1

[John D'Errico]

I don't know why you need an alpha shape for a convex domain. convhull will give you the polygon directly.

[Peter]

Ah, I should have explained that. I am trying to write code to create a concave hull. Using the method I have, with the for loop to choose a radius to create the alphashape, it often creates a concave hull, but I have found this data set as an example where not all of the points are included within the polygon.

I am therefore trying to use inpolygon to find how many points are outside of the polygon, and if there are points outside, I will then increase the radius by a factor.

Answer 2

[Sean de Wolski]

If you look at the polygon as you have it defined:

imshow(poly2mask(500*x,500*y,526,526))

it appears you need to sort your coordinates in a better way and that inpolygon appears to be working correctly.

[Peter]

I don't understand how I can take the coordinates of the polygon from the alphashape and put them into the correct form to use with inpolygon. I need to find a way to automatically sort the coordinates correctly as the program will usually run with varying datasets and therefore, a different polygon for each dataset.

Please excuse me if I am wrong, but surely the line of code that you have shown is plotting the original dataset, rather than the polygon coordinates?