create direction based polygons

6 views (last 30 days)
Leyon
Leyon on 24 Jun 2014
Commented: Leyon on 1 Jul 2014
I have vectors of lat/lon and want to find points that fall within an 8-point quadrant from a source point.
Was attempting to use inpolygon but couldn't figure out how to create the compass triangles based on a starting location. The quadrants are +- 22.5 degrees from the 8 cardinal points [0-360,45,90,135,180,225,270,315].
Therefore, if a point is northeast of the origin it would lie within a boundary 22.5 and 67.5. As i found out, you can't just simply subtract those values from the origin lat and lon. Any help would be great.
in = inpolygon(lat,lon,[src_lat,x_lim(1),x_lim(1)],[src_lon, src_lon+src_lon*.225, src_lon-src_lon*.225]);
Below is a sample lat/lon pair. Assume my start point is lat = 31.56, lon = 91.39.
30.433 -91.759
29.761 -90.3
29.181 -92.199
33.065 -92.344
32.923 -90.332
32.906 -94.435
30.481 -93.619
28.779 -93.766
28.496 -92.399
29.584 -91.312
31.557 -93.136
34.089 -89.848
33.121 -89.231
32.418 -90.307
29.384 -90.569
29.517 -92.767
33.018 -93.567
31.584 -92.277
30.551 -90.219
32.844 -91.866
30.675 -90.202
29.538 -89.232
29.514 -89.553
32.335 -92.06
32.331 -90.235
29.598 -91.418
33.126 -90.745
33.766 -89.968
32.497 -90.36
31.199 -88.882

Accepted Answer

Kelly Kearney
Kelly Kearney on 24 Jun 2014
I think inpolygon is overkill in this case; you really only need to know the azimuth between your center point and the target points in order to assign a quadrant. Do you have the Mapping Toolbox? If so, the azimuth function will do all the work for you:
pt = [...
30.433 -91.759
29.761 -90.3
29.181 -92.199
33.065 -92.344
32.923 -90.332
32.906 -94.435
30.481 -93.619
28.779 -93.766
28.496 -92.399
29.584 -91.312
31.557 -93.136
34.089 -89.848
33.121 -89.231
32.418 -90.307
29.384 -90.569
29.517 -92.767
33.018 -93.567
31.584 -92.277
30.551 -90.219
32.844 -91.866
30.675 -90.202
29.538 -89.232
29.514 -89.553
32.335 -92.06
32.331 -90.235
29.598 -91.418
33.126 -90.745
33.766 -89.968
32.497 -90.36
31.199 -88.882];
lat = 31.56;
lon = -91.39;
% Calculate azimuth
az = azimuth(pt(:,1), pt(:,2), lat, lon);
% Bin by quadrant
qedge = [0 22.5:45:360 360];
[n, iq] = histc(az, qedge);
iq(iq == 9 | iq == 10) = 1;
% A quick plot of the results
[ltout, lnout] = reckon(lat, lon, 5, qedge(2:end-1));
plot([ones(1,8)*lon; lnout], [ones(1,8)*lat; ltout], 'k');
hold on;
scatter(pt(:,2), pt(:,1), 20, iq, 'filled');
  3 Comments
Kelly Kearney
Kelly Kearney on 30 Jun 2014
Just reverse the order of inputs for the azimuth calculation:
az = azimuth(lat, lon, pt(:,1), pt(:,2));
Leyon
Leyon on 1 Jul 2014
That did the trick. Didn't even consider that when trying to fix it.
Thanks.

Sign in to comment.

More Answers (0)

Products

Community Treasure Hunt

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

Start Hunting!