Asked by hanif hamden
on 27 Apr 2019

I have plot circle according to the grid points and I have a problem to extract the points within each circle. This is my coding..

How should I extract the points from each circle. and save it in separate files according to the grid points. anyone can help?

clc;clear;

% Number of Points

numPoints = 5475;

data = load('Sample.txt');

x = data (:,2);

y = data (:,1);

z = data (:,3);%

plot(x, y, '.r');

% Make circle.

grid = load('Grid025.asc');

[row,column]=size(grid);

x0 = grid (:,2);

y0 = grid (:,1);

R = 0.09;

% Plot circle

for i=1:row

pos = [x0(i)-R, y0(i)-R, 2*R, 2*R];

rectangle('Position', pos, 'Curvature',[1 1]);

hold on;

% Determine how many points are in the circle .

count(i) = sum(((x-x0(i)).^2+(y-y0(i)).^2<=R^2));

listC = count.';

%Extract Data in Each Circle

distances = sqrt((x-x0(i)).^ 2 + (y-y0(i)).^ 2);

InsideCircle = distances <=R;

xS = x(InsideCircle); % Selected Lon

yS = y(InsideCircle); % Selected Lat

zS = z(InsideCircle); % Selected h

end

Answer by Cedric Wannaz
on 27 Apr 2019

Edited by Cedric Wannaz
on 27 Apr 2019

Accepted Answer

We cannot test without your data, but the approach seems fine. You could use the PDIST2 function and extract everything in one shot, but as a first approach this is fine.

I will assume that the question is about saving the points, so here is small update of your code that should (if I understand your example well) generate an Excel file with 4 columns that store Circle ID, Point ID, Point X, Point Y, Point Z:

% Prealloc array of counts, and cell array for storing relevant point parameters per circle.

counts = zeros( row, 1 ) ;

output = cell( row, 1 ) ;

% Initialize figure.

figure() ;

hold on ;

% Process circles.

for circleId = 1 : row

% Plot circle.

pos = [x0(circleId)-R, y0(circleId)-R, 2*R, 2*R] ;

rectangle( 'Position', pos, 'Curvature', [1 1] ) ;

% Compute distance from all points to current circle center.

distances = sqrt( (x - x0(circleId)).^2 + (y - y0(circleId)).^2 ) ;

% Flag and count relevant points.

isInsideCircle = distances <= R ;

counts(circleId) = sum( isInsideCircle ) ;

% Extract relevant parameters.

pointsId = find( isInsideCircle ) ;

pointsX = x(isInsideCircle) ; % Selected Lon

pointsY = y(isInsideCircle) ; % Selected Lat

pointsZ = z(isInsideCircle) ; % Selected h

% Store in relevant output cell. We build the first column as a vector of elements equal

% to circleId, whose size matches the size of the other columns (e.g. pointsId).

output{circleId} = [circleId * ones( size( pointsId )), pointsId, pointsX, pointsY, pointsZ] ;

end

% Build data to export to Excel by concatenating all cells of the output cell array vertically

% into a large numeric array, converting it into a cell array, and concatenating to a first

% "header" row.

output = num2cell( vertcat( output{:} )) ;

output = [{'CircleID', 'PointID', 'PointX', 'PointY', 'PointZ'}; output] ;

% Export to Excel.

xlswrite( 'PointsPerCircle.xlsx', output ) ;

Note that I didn't test it, but it illustrates a possible approach. Let me know if you have any question.

Cedric Wannaz
on 28 Apr 2019

It was a good attempt, but you tried to concatenate a scalar (e.g. x0(circleId) which is a single number) with column vectors, hence the error. If you look, circleId is also a scalar and I multiplied it by a vector of 1s whose size is compatible with other column vectors before concatenation. You have to do the same for x0 and y0.

I am updating the last version of the code, because it seems that you have been using the former one. Here, I create a column vector of 1s that I use for expanding circleId, X0, Y0:

clc ;

clear ;

% Load points.

numPoints = 5475 ; % *** Is this useful?

data = load( 'Sample.txt' ) ;

x = data (:, 2 ) ;

y = data (:, 1) ;

z = data (:, 3) ;

plot( x, y, '.r' ) ;

% Load circles.

grid = load( 'Grid025.asc' ) ;

[nCircles, ~] = size( grid ) ;

x0 = grid (:, 2) ;

y0 = grid (:, 1) ;

R = 0.09 ;

% Prealloc array of counts, and cell array for storing relevant point parameters per circle.

countPerCircle = zeros( nCircles, 1 ) ;

pointsPerCircle = cell( nCircles, 1 ) ;

% Initialize figure.

figure() ;

hold on ;

% Process circles.

for circleId = 1 : nCircles

% Plot circle.

pos = [x0(circleId)-R, y0(circleId)-R, 2*R, 2*R] ;

rectangle( 'Position', pos, 'Curvature', [1, 1] ) ;

% Compute distance from all points to current circle center.

distances = sqrt( (x - x0(circleId)).^2 + (y - y0(circleId)).^2 ) ;

% Flag and count relevant points.

isInsideCircle = distances <= R ;

countPerCircle(circleId) = sum( isInsideCircle ) ;

% Extract relevant parameters.

pointsId = find( isInsideCircle ) ;

pointsX = x(isInsideCircle) ; % Selected Lon

pointsY = y(isInsideCircle) ; % Selected Lat

pointsZ = z(isInsideCircle) ; % Selected h

% Build vector of 1s for expanding scalar to column vectors whose size

% is compatible with other column vectors. Expand circleId, and relevant

% elements of x0 and y0.

onesVector = ones( size( pointsId )) ;

circleIdExpanded = circleId * onesVector ;

x0Expanded = x0(circleId) * onesVector ;

y0Expanded = y0(circleId) * onesVector ;

% Store in relevant output cell. We build the first column as a vector of elements equal

% to circleId, whose size matches the size of the other columns (e.g. pointsId).

pointsPerCircle{circleId} = [circleIdExpanded, x0Expanded, y0Expanded, pointsId, pointsX, pointsY, pointsZ] ;

end

% Build data to export to Excel by concatenating all cells of the output cell array vertically

% into a large numeric array, converting it into a cell array, and concatenating to a first

% "header" row.

pointsPerCircle = num2cell( vertcat( pointsPerCircle{:} )) ;

pointsPerCircle = [{'CircleID', 'X0', 'Y0', 'PointID', 'X', 'Y', 'Z'}; pointsPerCircle] ;

% Export to Excel, spreadsheet "Points".

xlswrite( 'PerCircle.xlsx', pointsPerCircle, 'Points' ) ;

% Build output buffer for "count per circle".

countPerCircle = num2cell( [(1 : nCircles).', x0, y0, countPerCircle] ) ;

countPerCircle = [{'CircleID', 'X0', 'Y0', 'Count'}; countPerCircle] ;

% Export to Excel, spreadsheet "Counts".

xlswrite( 'PerCircle.xlsx', countPerCircle, 'Counts' ) ;

hanif hamden
on 28 Apr 2019

Now I get it.Thank you very much Sir. :)

Cedric Wannaz
on 28 Apr 2019

My pleasure! :-)

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.