Sorting of points in 2D

Dear All,

I face issue in arranging the airfoil coordinates coming from a CAD program. The coordinates are arranged from minumum x value to maximum x value. But I require to arrange the points starting from Trailing Edge ( Max X coordinate) and then proceed on top side all the way to Leading Edge ( Minimum X coordinate) and then again proceed in the bottom side and terminate in the Trailing Edge.

The input.txt ( required input file to be sorted)

Output.txt (the coordinates sorted manually).

I have several airfoil coordinates of similar nature, so a matlab function could help me out doing the job manually everytime.

Any help is really appreciated.

Thanks in Advance,

K Vijay Anand

[Ameer Hamza]

Can you create a small input and output file? Maybe with a maximum of 20 lines. That will help in understanding the pattern.

[KSSV]

One way is to use the below rocedure:

https://in.mathworks.com/matlabcentral/answers/400247-how-to-plot-closed-polygon-from-set-of-random-points-with-edges-not-intersecting

[Vijay Anand]

Dear KSSV,

Thanks a lot for your quick responce. I have gone through your algorithm and works perfect for convex polygons !!!

In my case, it is failing near the trailing edge. There should be a smart way to retify this. Kindly go through the attached code and the figures.

%% Program to sort 2D Airfoil Coordinates

clear all; close all; clc

%% Read Input File

Table1 = readtable('Input.txt','PreserveVariableNames',true);

x = Table1.x' ; y = Table1.y';

% Plot the input airfoil

h1 = figure(1); plot(x,y,'.-');

xlabel('X-Coordinate'); ylabel('Y-Coordinate'); grid on; title('Input Coordinates');

%% Sorting using Centroid Algorithm

P = [x; y]; % coordinates / points

c = mean(P,2); % mean/ central point

d = P-c ; % vectors connecting the central point and the given points

th = atan2(d(2,:),d(1,:)); % angle above x axis

[th, idx] = sort(th); % sorting the angles

P = P(:,idx); % sorting the given points

P = [P P(:,1)]; % add the first at the end to close the polygon

x1 = P(1,:);

y1 = P(2,:);

%%

% Plot the output airfoil

h2 = figure(2); plot(x1,y1,'.-');

xlabel('X-Coordinate'); ylabel('Y-Coordinate'); grid on; title('Sorted Coordinates - Centroid Algorithm');

Answer 1

[Stephen23]

Here is some simple code that that works quite nicely for your aerofoil. The code assumes that the top and bottom are both functions of x (i.e. no implicit curves in some other variable) and that they are both sampled with approximately the same step size (like your data). A convolution is used to smooth the data to create a moving average over 3 data points. The y-data are compared against this moving average to decide if the data point belongs to the top or bottom of the aerofoil, and then those two sets of points are sorted separately.

Inp = readtable('Input.txt','Delimiter','\t');
tmp = conv(Inp.y([1,1:end,end]),ones(1,3)/3,'valid'); % smoothed data
idx = Inp.y >= tmp; % true/false == top/bottom
idx(end) = true; % trailing edge -> top
Out = [sortrows(Inp(idx,:),-1);sortrows(Inp(~idx,:),+1)];
writetable(Out,'Out.txt','Delimiter','\t')

Checking against your original "Output.txt" and plotting the data:

>> Chk = readtable('Output.txt','Delimiter','\t');
>> isequal(Chk,Out)
ans =
     1
>> plot(Inp.x,Inp.y,'-*', Inp.x,tmp,'-+', Out.x,Out.y,'-x')
>> legend('original','moving average','sorted')

[Image Analyst]

What about boundary() or envelope()? Might those work?

[Stephen23]

@Image Analyst: good suggestions. The function boundary (available for >=R2014b) is certainly easy to use and provides a simple solution but is slower (it relies on the alphashape class, which is not trivial code) and requires a bit of fiddling to find a suitable shrink factor:

>> idy = boundary(T.x,T.y,0.97);
>> V = T(idy,:);
>> plot(V.x,V.y,'-d')

Below 0.97 misses out a few data points along the bottom towards the trailing edge**, above 0.98 starts to make extra connections between the top and bottom... I guess this depends on the sample step size and curvature, so most likely there is no universal solution which does not require fine-tuning.

** now that I look at the above figure, it also seems to be missing a few. 0.98 is probably about right.

Answer 2

[Vijay Anand]

Respected Stephen Cobeldick,

Fantastic Approach to the problem !! You saved my day !

I have tried with different airfoil profiles. Some profiles require moving average with a windows size of 5 or 7 instead of 3.

So, i have made use of the following function movingmean to change the window size after visual inspection of the output.

Glen (2020). Moving Average Function (https://www.mathworks.com/matlabcentral/fileexchange/41859-moving-average-function), MATLAB Central File Exchange. Retrieved October 18, 2020.

The final code :

T = readtable('Input.txt','Delimiter','\t');

vec = movingmean(T.y,5); % Moving Average with Window Size 5

idx = T.y >= vec; % true/false == top/bottom

idx(end) = true; % trailing edge -> top

U = [sortrows(T(idx,:),-1);sortrows(T(~idx,:),+1)];

writetable(U,'Out.txt','Delimiter','\t')

THANKS A LOT FOR YOUR TIME AND EFFORT !!!

[Steven Lord]

Consider using the movmean function included in MATLAB.

[Stephen23]

I did some more experimentation and got better results replacing the moving average with a polynomial fitted to the data. First I tried polyfit, but quickly realized that the end points need to be specified. Luckily this FEX submission

https://www.mathworks.com/matlabcentral/fileexchange/38926-fit-polynomial-to-data-with-constraints

allows constraints to be specified, including points that the curve has to pass through:

pts = Inp{[1,end],:}; % polynomial must pass through the leading & trailing edges
pwr = mmpolyfit(Inp.x,Inp.y, 5, 'Point',pts) % fit order 5 polynomial
tmp = polyval(pwr, Inp.x); % replaces TMP in the original answer

And checking the fit:

plot(Inp.x,Inp.y,'-*', Inp.x,tmp,'-+')