CODEWITHSUNDEEP
algorithmgpscomputational-geometry
samxiao
samxiao

Reputation: 2667

Find coordinates inside a rectangular area constructed by lat/long GPS pairs

I've never deal much with location-based data, so very much new to the whole GPS coding related questions. I have a problem that I don't seem to find a very efficient way in solving it or maybe there's an algorithm that I'm not too sure.

Let said you have given 4 lat/long coordinates which construct some kind of a rectangular area: (X0, Y0), (X1, Y0), (X0, Y1), (X1, Y1)

-----------------------------------
|                   b             |
|   a                             |
|                                 | d
|                                 |
|     c                           |
-----------------------------------
             e

Is there a way to find all the point that are inside the given rectangular area : a, b, c
And all the points outside of the area? e, d

I can easily to construct a 2D matrix to do this, but that's only if the coordinates are in integer, but with lat/long pairs, the coordinates are usually in float numbers which we cannot use it to construct a 2D table.

Any cool ideas?

Edited 1:

What about this Ray-casting algorithm? Is this a good algorithm to be used for GPS coordinates which is a float number?

Upvotes: 0

Views: 2810

Answers (3)

AlexWien
AlexWien

Reputation: 28767

There is an unlimited number of points inside a rectangle, so you have to define a step with (distane between two points).

You could just iterate with two nested loops,
lat, lon coordinates can be converted to integer using a multiplication factor of:

multiply with 1E7 (10000000) to get maximum acuracy of 1cm, or 10000000: 1cm
1000000: 10cm
100000: 1m
10000: 10m
1000: 100m
100: 1km
10: 11km
1: 111km

Now iterate

// convert  to spherical integer rectangle
double toIntFact = 1E7;
int x = (int) (x0 * toIntFact);
int y = (int) (y0 * toIntFact);
int tx1 = x1 * toIntFact;
int ty1 = y1 * toIntFact;

int yStep = 100000; // about 1.11 m latitudinal span. choose desired step above in list
int xStep = (int) (yStep / cos(Math.toRadians(y0))); // longitude adaption factor depending of cos(latitude); more acurate (symetric) is to use cos of centerLatitude: (y0 + y1) / 2;

for (int px = x; px < tx1; px+= xStep) {
   for (int py = y; py < ty1; py+= yStep) {
       drawPoint(px, py); // or whatever
   }
}

This should give an point set with same distances inbetween point for about some kilometer wide rectangles. The code does not work when overlapping the Datum limit (-180 to 180 jump) or when overlapping the poles. Delivers useable results up to latitude 80° N or S.

This code uses some kind of implicit equidistant (equirectangular) projection (see the division by cos(centerLat) to correct the fact that 1 degree of latitude is another distance measured in meters than one degree of longitude.

If the size of the rectangle exceeds some ten or hundred kilomters, then depending on your requirements have to use an advanced projection: e.g convert lat, lon with an WGS84 to UTM conversion. The result are coordinates in meters, which then you iterate analog.

But are you sure that this is what you want? Nobody wants to find all atoms inside a rectangle. May all screen pixels, or a method isInsideRectangle(lat,lon, Rectangle); So think again for what you need that.

Upvotes: 0

Tim Dumol
Tim Dumol

Reputation: 584

An alternative solution to @YvesDaoust's and @EyalSchneider's is to find the winding number or the crossing number of each point (http://geomalgorithms.com/a03-_inclusion.html). This solution scales to a polygon of any number of vertices (regardless of axis-alignment).

The Crossing Number (cn) method - which counts the number of times a ray starting from the point P crosses the polygon boundary edges. The point is outside when this "crossing number" is even; otherwise, when it is odd, the point is inside. This method is sometimes referred to as the "even-odd" test.

The Winding Number (wn) method - which counts the number of times the polygon winds around the point P. The point is outside only when this "winding number" wn = 0; otherwise, the point is inside.

Incidentally, @YvesDaoust's solution effectively calculates the crossing number of the point.

Upvotes: 0

user1196549
user1196549

Reputation:

If your rectangle is axis-aligned, @Eyal's answer is the right one (and you actually don't need 8 values but 4 are enough).

If you deal with a rotated rectangle (will work for any quadrilateral), the ray-casting method is appropriate: consider the horizontal line Y=Yt through your test point and find the edges that cross it (one endpoint above, one endpoint below). There will be 0 or 2 such edges. In case 0, you are outside. Otherwise, compute the abscissas of the intersections of these edges with the line. If 0 or 2 intersection are on the left of the test point, you are outside.

Xi= Xt + (Yt - Y0) (X1 - X0) / (Y1 - Y0)

Upvotes: 1

Related Questions