Select points inside a polygon

17 次查看(过去 30 天)
I have 96 points (longitude and latitude); each point is a center of 0.5 x 0.5 pixel (box), on the other hand, I have a polygon, I want to select each pixel that placed in or on my polygon.
I can do it easily like this scrip below but the problem of this code is it just considers a point (center of the pixel) not an entire pixel to index.
polygon1_x = polygon1_x.'; % x of polygons
polygon1_y = polygon1_y.'; % y of polygons
lat = Points.lat; % x of points = lat
lon = Points.lon; % y of points = lon
[in,on] = inpolygon(lat,lon,polygon1_x,polygon1_y); % Logical Matrix
inon = in | on; % Combine ‘in’ And ‘on’
idx = find(inon(:)); % Linear Indices Of ‘inon’ Points
latcoord = lat(idx); % X-Coordinates Of ‘inon’ Points
loncoord = lon(idx); % Y-Coordinates Of ‘inon’ Points
clf
figure(1)
plot(lon, lat, '.') % Plot All Points
hold on
plot(polygon1_y, polygon1_x, '.') % Plot Polygon
plot(loncoord, latcoord, 'gp')
OUTPUT = idx; % the output is idx
So I would be grateful if anyone can told me how I can index the row number of points from Point.mat which placed in/on my shapefile if the points are center of 0.5 x0.5 pixels. In this way, I can effectively select pixels that are in or on my polygon. I need the output like idx.
Thank you so much.
  2 个评论
BN
BN 2020-8-4
Dear KSSV,
I think I'm not clear in the previous question while the answer is correct but I need a different thing which I explain it more clear here. So I tried to explain more clearly in this new question.
I have some points which are center of 0.5x0.5 pixels; I want to select each pixel that are inside or on the polygon but select the coordinate of center of it (my points) not coordinates of pixels. I just have the coordinate of the center of each pixel. since there are just lat and lon columns and I have another row like precipitation and temperature it is important to me to have an index in order to use it and cut my data from the original big file. I'm sorry if I didn't explain well in the previous question.
Thank you so much, I would be grateful if you can help me with this issue.

请先登录,再进行评论。

采纳的回答

Bruno Luong
Bruno Luong 2020-8-4
编辑:Bruno Luong 2020-8-4
Here is the code using POLYSHAPE
load('Points.mat')
load('polygon1_x.mat')
load('polygon1_y.mat')
lat = Points.lat; % x of points = lat
lon = Points.lon; % y of points = lon
wpixel = 0.5;
hpixel = 0.5;
pixel0 = [-1 1 1 -1;
-1 -1 1 1]' .* [wpixel, hpixel]/2;
Island = polyshape([polygon1_y(:) polygon1_x(:)]);
for k=1:length(lat)
poly = polyshape([lon(k) lat(k)]+pixel0);
Ik = intersect(poly,Island);
isin = ~isempty(Ik.Vertices);
pixel(k).poly = poly;
pixel(k).isin = isin;
end
% Index of pixel that intersect with Island
idxin = find([pixel.isin]) % <=== HERE IS THE INDEX YOU NEED
% Graphical output
close all
plot(Island)
hold on
for k=1:length(pixel)
if pixel(k).isin
color = 'r';
else
color = 'w';
end
plot(pixel(k).poly, 'FaceColor', color);
end
axis equal
  1 个评论
BN
BN 2020-8-4
I am really thank you, that is exactly what I needed. Thank you again

请先登录,再进行评论。

更多回答(2 个)

Bruno Luong
Bruno Luong 2020-8-4
But you already have the index in your code. here
...
inon = in | on; % Combine .in. And .on’
idx = find(inon(:));
...
What you think the idx are? Just use it
T(idx)
Precipitation(idx)
  7 个评论

请先登录,再进行评论。


KSSV
KSSV 2020-8-4
You can use the previous questions answer..with a slight change while saving the points.
clc; clear all ;
load("polygon1_x.mat") ;
load("polygon1_y.mat") ;
load("lon.mat") ;
load("lat.mat") ;
% remove Nan's from the data
xv = polygon1_y ; yv = polygon1_x ;
xv(isnan(xv)) = [] ;
yv(isnan(yv)) = [] ;
%
loncoord = cell([],1) ;
latcoord = cell([],1) ;
count = 0 ;
for i = 1:length(lon)
i
% make grid around (lon,lat)
x = lon(i)-0.5:0.5:lon(i)+0.5 ;
y = lat(i)-0.5:0.5:lat(i)+0.5 ;
[X,Y] = meshgrid(x,y) ;
[in,on] = inpolygon(X(:),Y(:),xv,yv); % Logical Matrix
inon = in | on; % Combine ‘in’ And ‘on’
idx = find(inon(:)); % Linear Indices Of ‘inon’ Points
if any(idx)
count = count+1 ;
loncoord{count} = X(idx); % Y-Coordinates Of ‘inon’ Point
latcoord{count} = Y(idx) ;
end
end
plot(polygon1_y,polygon1_x,'b')
hold on
for i = 1:length(loncoord)
plot(loncoord{i},latcoord{i},'*r')
end
  1 个评论
BN
BN 2020-8-4
I'm so sorry, But I have no idea what is the coordinates in "latcoord" and "loncoord" are, they are different from my originals points. I need to choose among my Points, not new points.
Also, I don't need coordinates of the points (like latcoord and loncoord) I just need the index which includes the row number of points that are in/on polygon from points.mat.
I don't know what to do. But I am appreciate your help and your patiently
Thank you

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Polygons 的更多信息

产品


版本

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by