How to get Lat/Lon from vec2mtx

2 views (last 30 days)
Pete sherer
Pete sherer on 14 Mar 2018
Answered: Pete sherer on 14 Mar 2018
I am trying to generate grid lat/lon for a given polygon. First I tried to get the grids on the boundary first using vec2mtx. Then I generate grid inside polygon. My questions are 1. How to get the right Lat/Lon from the vec2mtx 2. What would be most efficient way to generate grid inside polygons 3. These grid data are just the points, not polygon. Is there any Matlab function I can convert them into polygon object.
oriLat=[48.3860120000000;49.3686790000000;48.0145280000000;46.6825390000000;47.4787370000000;46.5073030000000;46.5176190000000;44.9243550000000;45.2906840000000;42.3015220000000;41.9182910000000;44.6990460000000;45.7893010000000;45.2741950000000;43.7351650000000;43.8312240000000;41.7326470000000;41.3863160000000;44.9774490000000;45.2345130000000;47.3537330000000;44.8174190000000;43.3218860000000;41.5511100000000;39.4791860000000;37.8949860000000;37.1004480000000;39.4531150000000;38.0365030000000;38.6352170000000;35.2258250000000;31.1253830000000;26.7715300000000;24.8661310000000;29.9185130000000;30.6849560000000;30.0184100000000;28.9338120000000;29.4307800000000;27.9083070000000;25.8403780000000;29.7380790000000;28.9821380000000;31.7520090000000;31.3322390000000;34.5540170000000;40.2609740000000;46.2591090000000;48.3860120000000];
oriLon=[-124.725839000000;-94.9521110000000;-89.4892260000000;-91.9618890000000;-87.9292690000000;-87.3667670000000;-84.1179250000000;-87.8434330000000;-86.9777800000000;-87.8347690000000;-86.5978990000000;-86.2484740000000;-84.7727650000000;-83.3851040000000;-83.9477400000000;-82.6336410000000;-83.4538320000000;-82.4605990000000;-74.9927560000000;-70.8444300000000;-68.2697100000000;-66.9498950000000;-70.5538540000000;-69.9649820000000;-75.5930680000000;-75.3386230000000;-75.9796080000000;-76.0123120000000;-76.3220930000000;-77.2467040000000;-75.5336270000000;-81.4020960000000;-80.0321200000000;-80.6511890000000;-83.6792190000000;-88.0083960000000;-89.8450650000000;-89.4009660000000;-94.6703890000000;-97.0033250000000;-97.4226360000000;-101.400636000000;-103.281190000000;-106.417940000000;-111.074825000000;-120.622575000000;-124.363414000000;-123.547659000000;-124.725839000000];
% plot( oriLon, oriLat, '-ob')
gridDensity = 10;
[inGrid, inRefVec] = vec2mtx( oriLat, oriLon, gridDensity);
% Trying to get grid Lat/Lon
[latlim, lonlim] = limitm(inGrid, inRefVec);
gridLon = [lonlim(1):1/gridDensity:lonlim(2)];
gridLat = [latlim(1):1/gridDensity:latlim(2)];
bndryLon = gridLon(1:end-1)+1/(gridDensity*2);
bndryLat = gridLat(1:end-1)+1/(gridDensity*2);
[Xbndry,Ybndry]= meshgrid( bndryLon, bndryLat);
tloc= inGrid>0;
% boundary grid
bndryGrid= struct;
bndryGrid.Lon = Xbndry(tloc); % only on the boundary
bndryGrid.Lat = Ybndry(tloc);
% making inner grid
[Xgrid,Ygrid]= meshgrid( gridLon+1/(gridDensity*2), gridLat+1/(gridDensity*2));
tloc= inpolygon( Xgrid, Ygrid, oriLon, oriLat); % This takes a really long time for large/dense gird or complex original polygon
innerGrid= struct;
innerGrid.Lon= Xgrid( tloc);
innerGrid.Lat= Ygrid( tloc);
figure; hold on
plot( oriLon, oriLat, '-k')
plot( bndryGrid.Lon, bndryGrid.Lat,'ob')
plot( innerGrid.Lon, innerGrid.Lat, '.r')
% After combining the bndryGrid and innerGrid, how to convert them into polygon?

Accepted Answer

KSSV
KSSV on 14 Mar 2018
You have two options:
1. You can generate points inside in a scattered way...i.e generate a unstructured mesh, and you will get points inside that boundary specified. To do this, you can use Mesh2D package. https://in.mathworks.com/matlabcentral/fileexchange/25555-mesh2d-delaunay-based-unstructured-mesh-generation.
Hey, don't forget to comment and rate if you like this option at the link.
2. Generate points in a structured way. For this option check the below code.
oriLat=[48.3860120000000;49.3686790000000;48.0145280000000;46.6825390000000;47.4787370000000;46.5073030000000;46.5176190000000;44.9243550000000;45.2906840000000;42.3015220000000;41.9182910000000;44.6990460000000;45.7893010000000;45.2741950000000;43.7351650000000;43.8312240000000;41.7326470000000;41.3863160000000;44.9774490000000;45.2345130000000;47.3537330000000;44.8174190000000;43.3218860000000;41.5511100000000;39.4791860000000;37.8949860000000;37.1004480000000;39.4531150000000;38.0365030000000;38.6352170000000;35.2258250000000;31.1253830000000;26.7715300000000;24.8661310000000;29.9185130000000;30.6849560000000;30.0184100000000;28.9338120000000;29.4307800000000;27.9083070000000;25.8403780000000;29.7380790000000;28.9821380000000;31.7520090000000;31.3322390000000;34.5540170000000;40.2609740000000;46.2591090000000;48.3860120000000];
oriLon=[-124.725839000000;-94.9521110000000;-89.4892260000000;-91.9618890000000;-87.9292690000000;-87.3667670000000;-84.1179250000000;-87.8434330000000;-86.9777800000000;-87.8347690000000;-86.5978990000000;-86.2484740000000;-84.7727650000000;-83.3851040000000;-83.9477400000000;-82.6336410000000;-83.4538320000000;-82.4605990000000;-74.9927560000000;-70.8444300000000;-68.2697100000000;-66.9498950000000;-70.5538540000000;-69.9649820000000;-75.5930680000000;-75.3386230000000;-75.9796080000000;-76.0123120000000;-76.3220930000000;-77.2467040000000;-75.5336270000000;-81.4020960000000;-80.0321200000000;-80.6511890000000;-83.6792190000000;-88.0083960000000;-89.8450650000000;-89.4009660000000;-94.6703890000000;-97.0033250000000;-97.4226360000000;-101.400636000000;-103.281190000000;-106.417940000000;-111.074825000000;-120.622575000000;-124.363414000000;-123.547659000000;-124.725839000000];
N = 100 ; % can be varied
lon = linspace(min(oriLon),max(oriLon),N) ;
lat = linspace(min(oriLat),max(oriLat),N) ;
%
[Lon,Lat] = meshgrid(lon,lat) ;
% pick points lying inside the given region
idx = inpolygon(Lon(:), Lat(:),oriLon, oriLat) ;
Lon(~idx) = NaN ; Lat(~idx) = NaN ;
plot( oriLon, oriLat, 'b')
hold on
plot(Lon,Lat,'.r')
  2 Comments
KSSV
KSSV on 14 Mar 2018
Edited: KSSV on 14 Mar 2018
Pete sherer commented:
Thank you for your suggestion.
I realized for large/complex polygon, this inpolygon function is extremely slow. IS there any spatial join in Matlab similar to panda in python?
idx = inpolygon(Lon(:), Lat(:), oriLon, oriLat) ;
KSSV
KSSV on 14 Mar 2018
What do you mean by spatial join?

Sign in to comment.

More Answers (1)

Pete sherer
Pete sherer on 14 Mar 2018
I realized this is a bit off track. Similar to the Table join in SQL but applying the concept to spatial data. Simply it's equivalent to a for-loop of inpolygon. Say we have 2 data set: 1) US-county polygons and 2) set of data points. If we want to determine which counties these belong to, we need to create for loop over the number of polygons. Within each loop, we check which belong to it using inpolygon function. (even much slower if each polygon has complex shape or many points). The geoSpatial you simply can "join" the full county-polygons with set of data points, it will identify county each point belong to. Simply it finds intersection of polygons and points.

Community Treasure Hunt

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

Start Hunting!