Create an area mask using a shapefile
显示 更早的评论
I am trying to mask my study area using a shapefile But unable to do so. Here is the code that i was trying.. I'll be really grateful if someone can help or provide me with some alternative method to achieve it.
Thank you so much.
% Load the data
load Data.mat
shapefile = 'India_shapefile.shp';
S = shaperead(shapefile);
polygon = polyshape([S.X], [S.Y]);
% Create a logical mask
logical_mask = inpolygon(lon, lat, polygon.Vertices(:, 1), polygon.Vertices(:, 2));
% Use the logical mask to extract data from matrix
extracted_data = grid_data(logical_mask);
采纳的回答
Walter Roberson
2023-9-6
0 个投票
9 个评论
Vidushi Koliyan
2023-9-6
Walter Roberson
2023-9-6
You are using a shape file for India, but your lat and lon correspond mostly to the Indian Ocean up to roughly Singapore. For example Jakarta is inside your search area.
Vidushi Koliyan
2023-9-7
Vidushi Koliyan
2023-9-7
Walter Roberson
2023-9-7
I would need something more specific -- the data file and to have you point out at least two example points for me to look at.
A common problem when examining netCDF files is that netCDF arrays are transposed compared to MATLAB conventions, and ncread does not transpose them when reading them. If the netcdf file says that the first dimension is latitude and the second dimension is longitude and you ncread() then what you get in MATLAB will have a first dimension that is longitude and a second dimension that is latitude. Another way of saying this is that arrays in netCF are indexed across rows -- netcdf_A(3,5) is third column 5th row in the file.
Vidushi Koliyan
2023-9-8
Vidushi Koliyan
2023-9-8
Walter Roberson
2023-9-9
% Load the data
load data.mat
shapefile = fullfile('India Shapefile', 'India_shapefile.shp');
S = shaperead(shapefile);
%remove small triangles and squares
coord_lengths = arrayfun(@(s) numel(s.X), S);
not_useful_mask = coord_lengths < 10;
S(not_useful_mask) = [];
polygon = polyshape({S.X}, {S.Y});
% Create a logical mask
[Lon, Lat] = meshgrid(lon, lat);
logical_mask = reshape(isinterior(polygon, Lon(:), Lat(:)), size(Lon));
% Use the logical mask to extract data from matrix
extracted_data = test_data(logical_mask);
subplot(3,1,1);
imagesc(lon, lat, test_data); title('data'); xlabel('lon'); ylabel('lat'); axis equal; set(gca, 'YDir', 'normal');
subplot(3,1,2);
plot(polygon); title('borders'); xlabel('lon'); ylabel('lat'); axis equal; set(gca, 'YDir', 'normal');
subplot(3,1,3);
imagesc(lon, lat, logical_mask); title('selected points'); xlabel('lon'); ylabel('lat'); axis equal; set(gca, 'YDir', 'normal');
You have to be careful with X and Y. Y is latitude, X is longitude, but you were passing in latitude (Y) into imagesc in the place it expects X coordinates. Also by default imagesc flips images, so we account for that.
Remember that extracted_data is going to be a vector, in the order that find() would report that logical_mask entries are true, which scans down columns. All of the true locations for one column will be extracted, then all of the true locations for the next column to the right, and so on. The code does not construct a "bounding box" or a convex hull. The extracted_data vector is generally not useful if you need to make any kind of spatial analysis. If you need to do something like contour the data within the bounds, then you would be better of doing something like
selected_data = test_data;
selected_data(~logical_mask) = nan;
Vidushi Koliyan
2023-9-9
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Coordinate Reference Systems 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
