I'm trying making a function that returns if a Geopoint (latitude, longitude) is inside a List<BasicGeoposition> (list of Earth's surface points). I'm showing all the code that I'm using, maybe it's a lot but if you want to make some fast testing, you have here all the code.
First, I have used the algorithm that I took from this interesting and very useful topic.
I found this topic also, with no answer, that points to the same topic.
public static bool isInPolygon(Geopoint gp, List<BasicGeoposition> polygon)
{
int polyCorners = polygon.Count; // poligonoak zenbat erpin dituen
int i;
int j = polyCorners - 1;
bool oddNodes = false;
for (i = 0; i < polyCorners; i++)
{
if ((polygon[i].Latitude < gp.Position.Latitude && polygon[j].Latitude >= gp.Position.Latitude
|| polygon[j].Latitude < gp.Position.Latitude && polygon[i].Latitude >= gp.Position.Latitude)
&& (polygon[i].Longitude <= gp.Position.Longitude || polygon[j].Longitude <= gp.Position.Longitude))
{
oddNodes ^= (polygon[i].Longitude + (gp.Position.Latitude - polygon[i].Latitude) / (polygon[j].Latitude - polygon[i].Latitude) * (polygon[j].Longitude - polygon[i].Longitude) < gp.Position.Longitude);
}
j = i;
}
return oddNodes;
}
Checking with different shapes of polygons I've found that it works ok most of the times, but there're some cases that the algorithm fails:
1-. Near pole polygons. I don't know how to handle these cases.
2-. It doesn't work with -180, 180 longitude values. So, in this cases I've made some "normalization" with the longitudes. Basically, it duplicates the points making one polygon in the negative side and another one in the positive side. Here the code:
public static List<BasicGeoposition> NormalizeByDuplicatePolygon(List<BasicGeoposition> unnormalizedPolygon)
{
var negativePolygon = new List<BasicGeoposition>();
var positivePolygon = new List<BasicGeoposition>();
if (HasPositiveAndNegativeLongitude(unnormalizedPolygon))
{
foreach (var item in unnormalizedPolygon)
{
// positiboa
if (item.Longitude > 0)
{
// Negatibo kopia sortu
BasicGeoposition bgp = new BasicGeoposition
{
Longitude = item.Longitude - 360,
Latitude = item.Latitude,
Altitude = item.Altitude
};
positivePolygon.Add(item);
negativePolygon.Add(bgp);
}
else
{
// Positibo kopia sortu
BasicGeoposition bgp = new BasicGeoposition
{
Longitude = item.Longitude + 360,
Latitude = item.Latitude,
Altitude = item.Altitude
};
positivePolygon.Add(bgp);
negativePolygon.Add(item);
}
}
var normalizedPolygon = new List<BasicGeoposition>();
normalizedPolygon.AddRange(negativePolygon);
normalizedPolygon.AddRange(positivePolygon);
return normalizedPolygon;
}
else
{
return unnormalizedPolygon;
}
}
Here again, it works ok in some cases (positive longitudes) and doesn't work when in the negative values of longitudes.
Here's the testing that I'm using to check the algorithm:
private static void TestInsideMapPolygon()
{
Debug.WriteLine("InsideMapPolygon\n-------------------------");
Debug.Indent();
Debug.Indent();
Geopoint gpIn = GlobalFunc.getGeopoint(0.5, 0.5);
Geopoint gpIn2 = GlobalFunc.getGeopoint(0.15, 0.15);
Geopoint gpOut = GlobalFunc.getGeopoint(5.5, 5.5);
// Basic Rectangle
List<BasicGeoposition> l1 = new List<BasicGeoposition>();
var bg00 = GlobalFunc.BasicGeoposition(0, 0);
var bg01 = GlobalFunc.BasicGeoposition(0, 1);
var bg11 = GlobalFunc.BasicGeoposition(1, 1);
var bg10 = GlobalFunc.BasicGeoposition(1, 0);
l1.Add(bg00);
l1.Add(bg01);
l1.Add(bg11);
l1.Add(bg10);
l1.Add(bg00);
bool inPoly = CoreEclipse.isInPolygon(gpIn, l1);
Debug.WriteLine(" Basic rectangle ");
Debug.WriteLine(" 0.5, 0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l1);
Debug.WriteLine(" 5.5, 5.5 is {0}", inPoly);
// Hollow rectangle
List<BasicGeoposition> l2 = new List<BasicGeoposition>()
{
GlobalFunc.BasicGeoposition(0, 0),
GlobalFunc.BasicGeoposition(0, 1),
GlobalFunc.BasicGeoposition(1, 1),
GlobalFunc.BasicGeoposition(1, 0),
GlobalFunc.BasicGeoposition(0.75, 0),
GlobalFunc.BasicGeoposition(0.75, 0.25),
GlobalFunc.BasicGeoposition(0.25, 0.25),
GlobalFunc.BasicGeoposition(0.25, 0.75),
GlobalFunc.BasicGeoposition(0.75, 0.75),
GlobalFunc.BasicGeoposition(0.75, 0),
//GlobalFunc.BasicGeoposition(0, 0),
};
inPoly = CoreEclipse.isInPolygon(gpIn, l2);
Debug.WriteLine(" Hollow rectangle");
Debug.WriteLine(" 0.5, 0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpIn2, l2);
Debug.WriteLine(" 0.15, 0.15 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l2);
Debug.WriteLine(" 5.5, 5.5 is {0}", inPoly);
// H rectangle
List<BasicGeoposition> l3 = new List<BasicGeoposition>()
{
GlobalFunc.BasicGeoposition(0, 0),
GlobalFunc.BasicGeoposition(1, 0),
GlobalFunc.BasicGeoposition(1, 1),
GlobalFunc.BasicGeoposition(0.75, 1),
GlobalFunc.BasicGeoposition(0.75, 0.25),
GlobalFunc.BasicGeoposition(0.25, 0.25),
GlobalFunc.BasicGeoposition(0.25, 1),
GlobalFunc.BasicGeoposition(0, 1),
//GlobalFunc.BasicGeoposition(0, 0),
};
inPoly = CoreEclipse.isInPolygon(gpIn, l3);
Debug.WriteLine(" H rectangle");
Debug.WriteLine(" 0.5, 0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpIn2, l3);
Debug.WriteLine(" 0.15, 0.15 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l3);
Debug.WriteLine(" 5.5, 5.5 is {0}", inPoly);
// -/+ rectangle
List<BasicGeoposition> l4 = new List<BasicGeoposition>()
{
GlobalFunc.BasicGeoposition(-5, 170),
GlobalFunc.BasicGeoposition(-5, -170),
GlobalFunc.BasicGeoposition(5, -170),
GlobalFunc.BasicGeoposition(5, 170),
//GlobalFunc.BasicGeoposition(-5, 170),
};
l4 = EclipseLoaderFunctions.NormalizeByDuplicatePolygon(l4);
string coord = "";
foreach (var item in l4)
{
coord += item.Latitude.ToString("0.00") + " , " + item.Longitude.ToString("0.00") + "\n";
}
Geopoint gpInExtrem = GlobalFunc.getGeopoint(0, -175);
Geopoint gpInExtrem2 = GlobalFunc.getGeopoint(0, 175);
Geopoint gpExtremOut = GlobalFunc.getGeopoint(5.5, 5.5);
inPoly = CoreEclipse.isInPolygon(gpInExtrem, l4);
Debug.WriteLine(" +/- longitude rectangle 1:");
Debug.WriteLine(" 0, -175 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpInExtrem2, l4);
Debug.WriteLine(" 0, 175 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpExtremOut, l4);
Debug.WriteLine(" 5.5, 5.5 is {0}", inPoly);
List<BasicGeoposition> l5= new List<BasicGeoposition>()
{
GlobalFunc.BasicGeoposition(-5.00 , -190.00),
GlobalFunc.BasicGeoposition(- 5.00 , -170.00),
GlobalFunc.BasicGeoposition(5.00 , -170.00),
GlobalFunc.BasicGeoposition(5.00 , -190.00),
//GlobalFunc.BasicGeoposition(-5, 170),
};
//l5.Reverse();
inPoly = CoreEclipse.isInPolygon(gpInExtrem, l4);
Debug.WriteLine(" +/- longitude rectangle 2:");
Debug.WriteLine(" 0, -175 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpInExtrem2, l4);
Debug.WriteLine(" 0, 175 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpExtremOut, l4);
Debug.WriteLine(" 5.5, 5.5 is {0}", inPoly);
List<BasicGeoposition> lNorthPole = new List<BasicGeoposition>();
for (int longitude = -180; longitude < 180; longitude+=10)
{
lNorthPole.Add(GlobalFunc.BasicGeoposition(85, longitude));
}
Geopoint gpNorthPole = GlobalFunc.getGeopoint(90, 0);
Geopoint gpNorthPoleOut = GlobalFunc.getGeopoint(80, 0);
Geopoint gpSouthPole = GlobalFunc.getGeopoint(-90, 0);
Geopoint gpSouthPoleOut = GlobalFunc.getGeopoint(-80, 0);
Debug.WriteLine(" North Pole polygon: ");
inPoly = CoreEclipse.isInPolygon(gpNorthPole, lNorthPole);
Debug.WriteLine(" 90, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpNorthPoleOut, lNorthPole);
Debug.WriteLine(" 80, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpSouthPole, lNorthPole);
Debug.WriteLine(" -90, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpSouthPoleOut, lNorthPole);
Debug.WriteLine(" -80, 0 is {0}", inPoly);
List<BasicGeoposition> lSouthPole = new List<BasicGeoposition>();
for (int longitude = -180; longitude < 180; longitude += 10)
{
lSouthPole.Add(GlobalFunc.BasicGeoposition(-85, longitude));
}
Debug.WriteLine(" South Pole polygon: ");
inPoly = CoreEclipse.isInPolygon(gpNorthPole, lSouthPole);
Debug.WriteLine(" 90, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpNorthPoleOut, lSouthPole);
Debug.WriteLine(" 80, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpSouthPole, lSouthPole);
Debug.WriteLine(" -90, 0 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpSouthPoleOut, lSouthPole);
Debug.WriteLine(" -80, 0 is {0}", inPoly);
// Negatibo bihurtuta
Debug.Unindent();
Debug.WriteLine("Negative polygons\n-------------------------");
Debug.Indent();
gpIn = GlobalFunc.getGeopoint(0.5, -0.5);
gpIn2 = GlobalFunc.getGeopoint(0.15, -0.15);
gpOut = GlobalFunc.getGeopoint(5.5, -5.5);
// Basic Rectangle
List<BasicGeoposition> l0negative = new List<BasicGeoposition>();
foreach (var item in l1)
{
l0negative.Add(GlobalFunc.BasicGeoposition(item.Latitude * 1, item.Longitude * -1));
}
inPoly = CoreEclipse.isInPolygon(gpIn, l0negative);
Debug.WriteLine(" Basic rectangle");
Debug.WriteLine(" 0.5, -0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l0negative);
Debug.WriteLine(" 5.5, -5.5 is {0}", inPoly);
// Hollow rectangle
l0negative.Clear();
foreach (var item in l2)
{
l0negative.Add(GlobalFunc.BasicGeoposition(item.Latitude * 1, item.Longitude * -1));
}
inPoly = CoreEclipse.isInPolygon(gpIn, l0negative);
Debug.WriteLine(" Hollow rectangle");
Debug.WriteLine(" 0.5, -0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpIn2, l0negative);
Debug.WriteLine(" 0.15, -0.15 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l0negative);
Debug.WriteLine(" 5.5, -5.5 is {0}", inPoly);
// H rectangle
l0negative.Clear();
foreach (var item in l3)
{
l0negative.Add(GlobalFunc.BasicGeoposition(item.Latitude * 1, item.Longitude * -1));
}
inPoly = CoreEclipse.isInPolygon(gpIn, l0negative);
Debug.WriteLine(" H rectangle ");
Debug.WriteLine(" 0.5, -0.5 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpIn2, l0negative);
Debug.WriteLine(" 0.15, -0.15 is {0}", inPoly);
inPoly = CoreEclipse.isInPolygon(gpOut, l0negative);
Debug.WriteLine(" 5.5, -5.5 is {0}", inPoly);
Debug.Unindent();
Debug.Unindent();
}
Debug output window shows:
InsideMapPolygon
-------------------------
Basic rectangle
0.5, 0.5 is True
5.5, 5.5 is False
Hollow rectangle
0.5, 0.5 is False
0.15, 0.15 is True
5.5, 5.5 is False
H rectangle
0.5, 0.5 is False
0.15, 0.15 is True
5.5, 5.5 is False
+/- longitude rectangle 1:
0, -175 is False ( wrong! should be True)
0, 175 is True
5.5, 5.5 is False
+/- longitude rectangle 2:
0, -175 is False ( wrong! should be True)
0, 175 is True
5.5, 5.5 is False
North Pole polygon:
90, 0 is False ( wrong! should be True)
80, 0 is False
-90, 0 is False
-80, 0 is False
South Pole polygon:
90, 0 is False
80, 0 is False
-90, 0 is False ( wrong! should be True)
-80, 0 is False
Negative polygons
-------------------------
Basic rectangle
0.5, -0.5 is True
5.5, -5.5 is False
Hollow rectangle
0.5, -0.5 is False
0.15, -0.15 is True
5.5, -5.5 is False
H rectangle
0.5, -0.5 is False
0.15, -0.15 is True
5.5, -5.5 is False
