Lora/Lora-Map
1
0
Fork 0
Code Issues 8 Pull Requests Releases 23 Wiki Activity

now with MGRS Coordinates

Browse Source
This commit is contained in:
BlubbFish 2019-04-26 00:18:05 +02:00
parent 3ad491ccd5
commit b82b071d4d
26 changed files with 46577 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1317
Lora-Map/Helper/Celestial.Assistant.cs Normal file
View File

File diff suppressed because it is too large Load Diff

412
Lora-Map/Helper/Celestial.LunarEclipseCalc.cs Normal file
View File

@ -0,0 +1,412 @@
using System;
using System.Collections.Generic;
namespace CoordinateSharp
{
//CURRENT ALTITUDE IS SET CONSTANT AT 100M. POSSIBLY NEED TO ADJUST TO ALLOW USER PASS.
//Altitude adjustments appear to have minimal effect on eclipse timing. These were mainly used
//to signify eclipses that had already started during rise and set times on the NASA calculator
//SOME TIMES AND ALTS WERE RETURNED WITH COLOR AND STYLING. DETERMINE WHY AND ADJUST VALUE AS REQUIRED. SEARCH "WAS ITALIC".
//ELLIPSOID ADJUSTMENT
//6378140.0 Ellipsoid is used in the NASA Calculator
//WGS84 Ellipsoid is 6378137.0. Adjustments to the ellipsoid appear to effect eclipse seconds in fractions.
//This can be modified if need to allow users to pass custom number with the Coordinate SetDatum() functions.
//CURRENT RANGE 1601-2600.
internal class LunarEclipseCalc
{
public static List<List<string>> CalculateLunarEclipse(DateTime d, double latRad, double longRad)
{
return Calculate(d, latRad, longRad);
}
public static List<LunarEclipseDetails> CalculateLunarEclipse(DateTime d, double latRad, double longRad, double[] events)
{
List<List<string>> evs = Calculate(d, latRad, longRad, events);
List<LunarEclipseDetails> deetsList = new List<LunarEclipseDetails>();
foreach (List<string> ls in evs)
{
LunarEclipseDetails deets = new LunarEclipseDetails(ls);
deetsList.Add(deets);
}
return deetsList;
}
public static List<List<string>> CalculateLunarEclipse(DateTime d, Coordinate coord)
{
return Calculate(d, coord.Latitude.ToRadians(), coord.Longitude.ToRadians());
}
// CALCULATE!
private static List<List<string>> Calculate(DateTime d, double latRad, double longRad, double[] ev = null)
{
//DECLARE ARRAYS
double[] obsvconst = new double[6];
double[] mid = new double[41];
double[] p1 = new double[41];
double[] u1 = new double[41];
double[] u2 = new double[41];
double[] u3 = new double[41];
double[] u4 = new double[41];
double[] p4 = new double[41];
List<List<string>> events = new List<List<string>>();
double[] el;
if (ev == null)
{
el = Eclipse.LunarData.LunarDateData(d);//Get 100 year solar data;
}
else
{
el = ev;
}
events = new List<List<string>>();
ReadData(latRad, longRad, obsvconst);
for (int i = 0; i < el.Length; i += 22)
{
if (el[5 + i] <= obsvconst[5])
{
List<string> values = new List<string>();
obsvconst[4] = i;
GetAll(el, obsvconst, mid, p1, u1, u2,u3,u4,p4);
// Is there an event...
if (mid[5] != 1)
{
values.Add(GetDate(el, p1, obsvconst));
if (el[5 + i] == 1)
{
values.Add("T");
}
else if (el[5 + i] == 2)
{
values.Add("P");
}
else
{
values.Add("N");
}
// Pen. Mag
values.Add(el[3 + i].ToString());
// Umbral Mag
values.Add(el[4 + i].ToString());
// P1
values.Add(GetTime(el, p1, obsvconst));
// P1 alt
values.Add(GetAlt(p1));
if (u1[5] == 1)
{
values.Add("-");
values.Add("-");
}
else
{
// U1
values.Add(GetTime(el, u1, obsvconst));
// U1 alt
values.Add(GetAlt(u1));
}
if (u2[5] == 1)
{
values.Add("-");
values.Add("-");
}
else
{
// U2
values.Add(GetTime(el, u2, obsvconst));
// U2 alt
values.Add(GetAlt(u2));
}
// mid
values.Add(GetTime(el, mid, obsvconst));
// mid alt
values.Add(GetAlt(mid));
if (u3[5] == 1)
{
values.Add("-");
values.Add("-");
}
else
{
// u3
values.Add(GetTime(el, u3, obsvconst));
// u3 alt
values.Add(GetAlt(u3));
}
if (u4[5] == 1)
{
values.Add("-");
values.Add("-");
}
else
{
// u4
values.Add(GetTime(el, u4, obsvconst));
// u4 alt
values.Add(GetAlt(u4));
}
// P4
values.Add(GetTime(el, p4, obsvconst));
// P4 alt
values.Add(GetAlt(p4));
events.Add(values);
}
}
}
return events;
}
// Read the data that's in the form, and populate the obsvconst array
private static void ReadData(double latRad, double longRad, double[] obsvconst)
{
// Get the latitude
obsvconst[0] = latRad;
// Get the longitude
obsvconst[1] = -1 * longRad; //PASS REVERSE RADIAN.
// Get the altitude
obsvconst[2] = 100; //CHANGE TO ALLOW USER TO PASS.
// Get the time zone
obsvconst[3] = 0; //GMT TIME
obsvconst[4] = 0; //INDEX
//SET MAX ECLIPSE TYPE
obsvconst[5] = 4;//4 is ALL Eclipses
}
// Populate the p1, u1, u2, mid, u3, u4 and p4 arrays
private static void GetAll(double[] elements, double[] obsvconst, double[] mid, double[] p1, double[] u1, double[] u2, double[] u3, double[] u4, double[] p4)
{
int index = (int)obsvconst[4];
p1[1] = elements[index + 9];
PopulateCircumstances(elements, p1, obsvconst);
mid[1] = elements[index + 12];
PopulateCircumstances(elements, mid, obsvconst);
p4[1] = elements[index + 15];
PopulateCircumstances(elements, p4, obsvconst);
if (elements[index + 5] < 3)
{
u1[1] = elements[index + 10];
PopulateCircumstances(elements, u1, obsvconst);
u4[1] = elements[index + 14];
PopulateCircumstances(elements, u4, obsvconst);
if (elements[index + 5] < 2)
{
u2[1] = elements[index + 11];
u3[1] = elements[index + 13];
PopulateCircumstances(elements, u2, obsvconst);
PopulateCircumstances(elements, u3, obsvconst);
}
else
{
u2[5] = 1;
u3[5] = 1;
}
}
else
{
u1[5] = 1;
u2[5] = 1;
u3[5] = 1;
u4[5] = 1;
}
if ((p1[5] != 0) && (u1[5] != 0) && (u2[5] != 0) && (mid[5] != 0) && (u3[5] != 0) && (u4[5] != 0) && (p4[5] != 0))
{
mid[5] = 1;
}
}
// Populate the circumstances array
// entry condition - circumstances[1] must contain the correct value
private static void PopulateCircumstances(double[] elements, double[] circumstances, double[] obsvconst)
{
double t, ra, dec, h;
int index = (int)obsvconst[4];
t = circumstances[1];
ra = elements[18 + index] * t + elements[17 + index];
ra = ra * t + elements[16 + index];
dec = elements[21 + index] * t + elements[20 + index];
dec = dec * t + elements[19 + index];
dec = dec * Math.PI / 180.0;
circumstances[3] = dec;
h = 15.0 * (elements[6 + index] + (t - elements[2 + index] / 3600.0) * 1.00273791) - ra;
h = h * Math.PI / 180.0 - obsvconst[1];
circumstances[2] = h;
circumstances[4] = Math.Asin(Math.Sin(obsvconst[0]) * Math.Sin(dec) + Math.Cos(obsvconst[0]) * Math.Cos(dec) * Math.Cos(h));
circumstances[4] -= Math.Asin(Math.Sin(elements[7 + index] * Math.PI / 180.0) * Math.Cos(circumstances[4]));
if (circumstances[4] * 180.0 / Math.PI < elements[8 + index] - 0.5667)
{
circumstances[5] = 2;
}
else if (circumstances[4] < 0.0)
{
circumstances[4] = 0.0;
circumstances[5] = 0;
}
else
{
circumstances[5] = 0;
}
}
// Get the date of an event
private static string GetDate(double[] elements, double[] circumstances, double[] obsvconst)
{
string[] month = new string[] { "Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec" };//Month string array
double t, jd, a, b, c, d, e;
string ans = "";
int index = (int)obsvconst[4];
// Calculate the JD for noon (TDT) the day before the day that contains T0
jd = Math.Floor(elements[index] - (elements[1 + index] / 24.0));
// Calculate the local time (ie the offset in hours since midnight TDT on the day containing T0).
t = circumstances[1] + elements[1 + index] - obsvconst[3] - (elements[2 + index] - 30.0) / 3600.0;
if (t < 0.0)
{
jd--;
}
if (t >= 24.0)
{
jd++;
}
if (jd >= 2299160.0)
{
a = Math.Floor((jd - 1867216.25) / 36524.25);
a = jd + 1 + a - Math.Floor(a / 4.0);
}
else
{
a = jd;
}
b = a + 1525.0;
c = Math.Floor((b - 122.1) / 365.25);
d = Math.Floor(365.25 * c);
e = Math.Floor((b - d) / 30.6001);
d = b - d - Math.Floor(30.6001 * e);
if (e < 13.5)
{
e = e - 1;
}
else
{
e = e - 13;
}
double year;
if (e > 2.5)
{
ans = c - 4716 + "-";
year = c - 4716;
}
else
{
ans = c - 4715 + "-";
year = c - 4715;
}
string m = month[(int)e - 1];
ans += m+ "-";
if (d < 10)
{
ans = ans + "0";
}
ans = ans + d;
//Leap Year Integrity Check
if (m == "Feb" && d == 29 && !DateTime.IsLeapYear((int)year))
{
ans = year.ToString() + "-Mar-01";
}
return ans;
}
// Get the time of an event
private static string GetTime(double[] elements, double[] circumstances, double[] obsvconst)
{
double t;
string ans = "";
int index = (int)obsvconst[4];
t = circumstances[1] + elements[1 + index] - obsvconst[3] - (elements[2 + index] - 30.0) / 3600.0;
if (t < 0.0)
{
t = t + 24.0;
}
if (t >= 24.0)
{
t = t - 24.0;
}
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + Math.Floor(t) + ":";
t = (t * 60.0) - 60.0 * Math.Floor(t);
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + Math.Floor(t);
if (circumstances[5] == 2)
{
return ans; //RETURNED IN ITAL DETERMINE WHY
}
else
{
return ans;
}
}
// Get the altitude
private static string GetAlt(double[] circumstances)
{
double t;
string ans = "";
t = circumstances[4] * 180.0 / Math.PI;
t = Math.Floor(t + 0.5);
if (t < 0.0)
{
ans = "-";
t = -t;
}
else
{
ans = "+";
}
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + t;
if (circumstances[5] == 2)
{
return ans; //returned in italics determine why
}
else
{
return ans;
}
}
}
}

221
Lora-Map/Helper/Celestial.MeeusTables.cs Normal file
View File

@ -0,0 +1,221 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace CoordinateSharp
{
internal partial class MeeusTables
{
//Ch 47
private static double[] Table47A_Arguments = new double[]
{
0,0,1,0,
2,0,-1,0,
2,0,0,0,
0,0,2,0,
0,1,0,0,
0,0,0,2,
2,0,-2,0,
2,-1,-1,0,
2,0,1,0,
2,-1,0,0,
0,1,-1,0,
1,0,0,0,
0,1,1,0,
2,0,0,-2,
0,0,1,2,
0,0,1,-2,
4,0,-1,0,
0,0,3,0,
4,0,-2,0,
2,1,-1,0,
2,1,0,0,
1,0,-1,0,
1,1,0,0,
2,-1,1,0,
2,0,2,0,
4,0,0,0,
2,0,-3,0,
0,1,-2,0,
2,0,-1,2,
2,-1,-2,0,
1,0,1,0,
2,-2,0,0,
0,1,2,0,
0,2,0,0,
2,-2,-1,0,
2,0,1,-2,
2,0,0,2,
4,-1,-1,0,
0,0,2,2,
3,0,-1,0,
2,1,1,0,
4,-1,-2,0,
0,2,-1,0,
2,2,-1,0,
2,1,-2,0,
2,-1,0,-2,
4,0,1,0,
0,0,4,0,
4,-1,0,0,
1,0,-2,0,
2,1,0,-2,
0,0,2,-2,
1,1,1,0,
3,0,-2,0,
4,0,-3,0,
2,-1,2,0,
0,2,1,0,
1,1,-1,0,
2,0,3,0,
2,0,-1,-2
};
private static double[] Table47B_Arguments = new double[]
{
0,0,0,1,
0,0,1,1,
0,0,1,-1,
2,0,0,-1,
2,0,-1,1,
2,0,-1,-1,
2,0,0,1,
0,0,2,1,
2,0,1,-1,
0,0,2,-1,
2,-1,0,-1,
2,0,-2,-1,
2,0,1,1,
2,1,0,-1,
2,-1,-1,1,
2,-1,0,1,
2,-1,-1,-1,
0,1,-1,-1,
4,0,-1,-1,
0,1,0,1,
0,0,0,3,
0,1,-1,1,
1,0,0,1,
0,1,1,1,
0,1,1,-1,
0,1,0,-1,
1,0,0,-1,
0,0,3,1,
4,0,0,-1,
4,0,-1,1,
0,0,1,-3,
4,0,-2,1,
2,0,0,-3,
2,0,2,-1,
2,-1,1,-1,
2,0,-2,1,
0,0,3,-1,
2,0,2,1,
2,0,-3,-1,
2,1,-1,1,
2,1,0,1,
4,0,0,1,
2,-1,1,1,
2,-2,0,-1,
0,0,1,3,
2,1,1,-1,
1,1,0,-1,
1,1,0,1,
0,1,-2,-1,
2,1,-1,-1,
1,0,1,1,
2,-1,-2,-1,
0,1,2,1,
4,0,-2,-1,
4,-1,-1,-1,
1,0,1,-1,
4,0,1,-1,
1,0,-1,-1,
4,-1,0,-1,
2,-2,0,1,
};
private static double[] Table47A_El_Er = new double[]
{
//El
6288774, 1274027,658314,213618,-185116,-114332,58793,57066,53322,45758,
-40923,-34720,-30383,15327,-12528,10980,10675,10034,8548,-7888,-6766,-5163,
4987,4036,3994,3861,3665,-2689,-2602,2390,-2348,2236,-2120,-2069,2048,-1773,
-1595,1215,-1110,-892,-810,759,-713,-700,691,596,549,537,520,-487,-399,-381,
351,-340,330,327,-323,299,294,0,
//Er
-20905355,-3699111,-2955968,-569925,48888,-3149,246158,-152138,-170733,-204586,
-129620,108743,104755,10321,0,79661,-34782,-23210,-21636,24208,30824,-8379,-16675,
-12831,-10445,-11650,14403,-7003,0,10056,6322,-9884,5751,0,-4950,4130,0,-3958,0,3258,
2616,-1897,-2117,2354,0,0,-1423,-1117,-1571,-1739,0,-4421,0,0,0,0,1165,0,0,8752
};
private static double[] Table47B_Eb = new double[]
{
5128122,280602,277693,173237,55413,46271,32573,17198,9266,8822,
8216,4324,4200,-3359,2463,2211,2065,-1870,1828,-1794,-1749,-1565,-1491,
-1475,-1410,-1344,-1335,1107,1021,833,
777,671,607,596,491,-451,439,422,421,-366,-351,331,315,302,-283,-229,
223,223,-220,-220,-185,181,-177,176,166,-164,132,-119,115,107
};
private static double Get_Table47A_Values(double[] values, int l, double t, bool sine)
{
//sine true returns El
//sine false return Er
//Er values start at 60 in the Table47A_El_Er array.
int nl = l * 4;
if (sine)
{
double e = 1;
if (Table47A_Arguments[nl + 1] != 0)
{
e = 1 - .002516 * t - .0000074 * Math.Pow(t, 2);
if (Math.Abs(Table47A_Arguments[nl + 1]) == 2)
{
e *= e;
}
}
return (Table47A_El_Er[l] * e) * Math.Sin(Table47A_Arguments[nl] * values[0] + Table47A_Arguments[nl + 1] * values[1] +
Table47A_Arguments[nl + 2] * values[2] + Table47A_Arguments[nl + 3] * values[3]);
}
else
{
double e = 1;
if (Table47A_Arguments[nl + 1] != 0)
{
e = 1 - .002516 * t - .0000074 * Math.Pow(t, 2);
if (Math.Abs(Table47A_Arguments[nl + 1]) == 2)
{
e *= e;
}
}
return (Table47A_El_Er[l + 60] * e) * Math.Cos(Table47A_Arguments[nl] * values[0] + Table47A_Arguments[nl + 1] * values[1] +
Table47A_Arguments[nl + 2] * values[2] + Table47A_Arguments[nl + 3] * values[3]);
}
}
private static double Get_Table47B_Values(double[] values, int l, double t)
{
int nl = l * 4;
double e = 1;
if (Table47B_Arguments[nl + 1] != 0)
{
e = 1 - .002516 * t - .0000074 * Math.Pow(t, 2);
if (Math.Abs(Table47B_Arguments[nl + 1]) == 2)
{
e *= e;
}
}
return (Table47B_Eb[l] * e) * Math.Sin(Table47B_Arguments[nl] * values[0] + Table47B_Arguments[nl + 1] * values[1] +
Table47B_Arguments[nl + 2] * values[2] + Table47B_Arguments[nl + 3] * values[3]);
}
}
}

872
Lora-Map/Helper/Celestial.MoonCalculations.cs Normal file
View File

@ -0,0 +1,872 @@
using System;
using System.Collections.Generic;
namespace CoordinateSharp
{
internal class MoonCalc
{
static double rad = Math.PI / 180; //For converting radians
//obliquity of the ecliptic in radians based on standard equinox 2000.
static double e = rad * 23.4392911;
/// <summary>
/// Gets Moon Times, Altitude and Azimuth
/// </summary>
/// <param name="date">Date</param>
/// <param name="lat">Latitude</param>
/// <param name="lng">Longitude</param>
/// <param name="c">Celestial</param>
public static void GetMoonTimes(DateTime date, double lat, double lng, Celestial c)
{
//Get current Moon Position to populate passed Alt / Azi for user specified date
MoonPosition mp = GetMoonPosition(date, lat, lng, c);
double altRad = mp.Altitude / Math.PI*180; //Convert alt to degrees
c.moonAltitude = (altRad - mp.ParallaxCorection); //Set altitude with adjusted parallax
c.moonAzimuth = mp.Azimuth / Math.PI*180 + 180; //Azimuth in degrees + 180 for E by N.
////New Iterations for Moon set / rise
bool moonRise = false;
bool moonSet = false;
//Start at beginning of day
DateTime t = new DateTime(date.Year, date.Month, date.Day, 0, 0, 0, DateTimeKind.Utc);
//Get start of day Moon Pos
MoonPosition moonPos = GetMoonPosition(t, lat, lng, c);
double alt1 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
DateTime? setTime = null;
DateTime? riseTime = null;
double hz = -.3 * rad;//Horizon degrees at -.3 for appearant rise / set
//Iterate for each hour of the day
for(int x = 1;x<=24;x++)
{
moonPos = GetMoonPosition(t.AddHours(x), lat, lng, c);//Get the next hours altitude for comparison
double alt2 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
//If hour 1 is below horizon and hour 2 is above
if(alt1 <hz && alt2 >=hz)
{
//Moon Rise Occurred
moonRise = true;
DateTime dt1 = t.AddHours(x - 1);
moonPos = GetMoonPosition(dt1, lat, lng, c);//Get the next hours altitude for comparison
double altM1 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
//Iterate through each minute to determine at which minute the horizon is crossed.
//Interpolation is more efficient, but yielded results with deviations up to 5 minutes.
//Investigate formula efficiency
for (int y = 1;y<=60;y++)
{
DateTime dt2 = t.AddHours(x-1).AddMinutes(y);
moonPos = GetMoonPosition(dt2, lat, lng, c);//Get the next hours altitude for comparison
double altM2 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
if (altM1<hz && altM2>=hz)
{
//interpolate seconds
double p = 60 * ((hz - altM1) / (altM2 - altM1));
riseTime = dt1.AddMinutes(y-1).AddSeconds(p);
break;
}
altM1 = altM2;
}
}
//if hour 2 is above horizon and hour 1 below
if(alt1>=hz && alt2 <hz)
{
//Moon Set Occured
moonSet = true;
DateTime dt1 = t.AddHours(x - 1);
moonPos = GetMoonPosition(dt1, lat, lng, c);//Get the next hours altitude for comparison
double altM1 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
//Iterate through each minute to determine at which minute the horizon is crossed.
//Interpolation is more efficient, but yielded results with deviations up to 5 minutes.
//Investigate formula efficiency
for (int y = 1; y <= 60; y++)
{
DateTime dt2 = t.AddHours(x - 1).AddMinutes(y);
moonPos = GetMoonPosition(dt2, lat, lng, c);//Get the next hours altitude for comparison
double altM2 = moonPos.Altitude - (moonPos.ParallaxCorection * rad);
if (altM1 >= hz && altM2 < hz)
{
//Interpolate seconds
double p = 60 * ((hz - altM2) / (altM1 - altM2));
setTime = dt1.AddMinutes(y).AddSeconds(-p);
break;
}
altM1 = altM2;
}
}
alt1 = alt2;
if(moonRise && moonSet) { break; }
}
c.moonSet = setTime;
c.moonRise = riseTime;
if (moonRise && moonSet) { c.moonCondition = CelestialStatus.RiseAndSet; }
else
{
if (!moonRise && !moonSet)
{
if (alt1 >= 0) { c.moonCondition = CelestialStatus.UpAllDay; }
else { c.moonCondition = CelestialStatus.DownAllDay; }
}
if (!moonRise && moonSet) { c.moonCondition = CelestialStatus.NoRise; }
if (moonRise && !moonSet) { c.moonCondition = CelestialStatus.NoSet; }
}
}
private static MoonPosition GetMoonPosition(DateTime date, double lat, double lng, Celestial cel)
{
//Set UTC date integrity
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
double d = JulianConversions.GetJulian_Epoch2000(date);
//Ch 47
double JDE = JulianConversions.GetJulian(date);//Get julian
double T = (JDE - 2451545) / 36525; //Get dynamic time.
double[] LDMNF = Get_Moon_LDMNF(T);
CelCoords c = GetMoonCoords(d, cel, LDMNF, T);
Distance dist = GetMoonDistance(date);
double lw = rad * -lng;
double phi = rad * lat;
double H = rad * MeeusFormulas.Get_Sidereal_Time(JDE) - lw - c.ra;
double ra = c.ra; //Adjust current RA formula to avoid needless RAD conversions
double dec = c.dec; //Adjust current RA formula to avoid needless RAD conversions
//Adjust for parallax (low accuracry increases may not be worth cost)
//Investigate
double pSinE = Get_pSinE(dec, dist.Meters) * Math.PI / 180;
double pCosE = Get_pCosE(dec, dist.Meters) * Math.PI / 180;
double cRA = Parallax_RA(dist.Meters, H, pCosE, dec, ra);
double tDEC = Parallax_Dec(dist.Meters, H, pCosE, pSinE, dec, cRA);
double tRA = ra - cRA;
dec = tDEC;
ra = tRA;
//Get true altitude
double h = altitude(H, phi, dec);
// formula 14.1 of "Astronomical Algorithms" 2nd edition by Jean Meeus (Willmann-Bell, Richmond) 1998.
double pa = Math.Atan2(Math.Sin(H), Math.Tan(phi) * Math.Cos(dec) - Math.Sin(dec) * Math.Cos(H));
//altitude correction for refraction
h = h + astroRefraction(h);
MoonPosition mp = new MoonPosition();
mp.Azimuth = azimuth(H, phi, dec);
mp.Altitude = h / Math.PI * 180;
mp.Distance = dist;
mp.ParallacticAngle = pa;
double horParal = 8.794 / (dist.Meters / 149.59787E6); // horizontal parallax (arcseconds), Meeus S. 263
double p = Math.Asin(Math.Cos(h) * Math.Sin(horParal/3600)); // parallax in altitude (degrees)
p *= 1000;
mp.ParallaxCorection = p;
mp.Altitude *= rad;
return mp;
}
private static CelCoords GetMoonCoords(double d, Celestial c, double[] LDMNF, double t)
{
// Legacy function. Updated with Meeus Calcs for increased accuracy.
// geocentric ecliptic coordinates of the moon
// Meeus Ch 47
double[] cs = Get_Moon_Coordinates(LDMNF, t);
double l = cs[0]; // longitude
double b = cs[1]; // latitude
CelCoords mc = new CelCoords();
mc.ra = rightAscension(l, b);
double ra = mc.ra / Math.PI * 180;
mc.dec = declination(l, b);
double dec = mc.dec / Math.PI * 180;
return mc;
}
public static void GetMoonIllumination(DateTime date, Celestial c, double lat, double lng)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
double d = JulianConversions.GetJulian_Epoch2000(date);
CelCoords s = GetSunCoords(d);
double JDE = JulianConversions.GetJulian(date);//Get julian
double T = (JDE - 2451545) / 36525; //Get dynamic time.
double[] LDMNF = Get_Moon_LDMNF(T);
CelCoords m = GetMoonCoords(d, c,LDMNF, T);
double sdist = 149598000,
phi = Math.Acos(Math.Sin(s.dec) * Math.Sin(m.dec) + Math.Cos(s.dec) * Math.Cos(m.dec) * Math.Cos(s.ra - m.ra)),
inc = Math.Atan2(sdist * Math.Sin(phi), m.dist - sdist * Math.Cos(phi)),
angle = Math.Atan2(Math.Cos(s.dec) * Math.Sin(s.ra - m.ra), Math.Sin(s.dec) * Math.Cos(m.dec) -
Math.Cos(s.dec) * Math.Sin(m.dec) * Math.Cos(s.ra - m.ra));
MoonIllum mi = new MoonIllum();
mi.Fraction = (1 + Math.Cos(inc)) / 2;
mi.Phase = 0.5 + 0.5 * inc * (angle < 0 ? -1 : 1) / Math.PI;
mi.Angle = angle;
c.moonIllum = mi;
string moonName = "";
int moonDate = 0;
//GET PHASE NAME
//CHECK MOON AT BEGINNING AT END OF DAY TO GET DAY PHASE
DateTime dMon = new DateTime(date.Year, date.Month, 1);
for(int x = 1;x<= date.Day;x++)
{
DateTime nDate = new DateTime(dMon.Year, dMon.Month, x, 0, 0, 0, DateTimeKind.Utc);
d = JulianConversions.GetJulian_Epoch2000(nDate);
s = GetSunCoords(d);
JDE = JulianConversions.GetJulian(nDate);//Get julian
T = (JDE - 2451545) / 36525; //Get dynamic time.
LDMNF = Get_Moon_LDMNF(T);
m = GetMoonCoords(d, c,LDMNF,T);
phi = Math.Acos(Math.Sin(s.dec) * Math.Sin(m.dec) + Math.Cos(s.dec) * Math.Cos(m.dec) * Math.Cos(s.ra - m.ra));
inc = Math.Atan2(sdist * Math.Sin(phi), m.dist - sdist * Math.Cos(phi));
angle = Math.Atan2(Math.Cos(s.dec) * Math.Sin(s.ra - m.ra), Math.Sin(s.dec) * Math.Cos(m.dec) -
Math.Cos(s.dec) * Math.Sin(m.dec) * Math.Cos(s.ra - m.ra));
double startPhase = 0.5 + 0.5 * inc * (angle < 0 ? -1 : 1) / Math.PI;
nDate = new DateTime(dMon.Year, dMon.Month, x, 23, 59, 59, DateTimeKind.Utc);
d = JulianConversions.GetJulian_Epoch2000(nDate);
s = GetSunCoords(d);
JDE = JulianConversions.GetJulian(nDate);//Get julian
T = (JDE - 2451545) / 36525; //Get dynamic time.
LDMNF = Get_Moon_LDMNF(T);
m = GetMoonCoords(d, c,LDMNF,T);
phi = Math.Acos(Math.Sin(s.dec) * Math.Sin(m.dec) + Math.Cos(s.dec) * Math.Cos(m.dec) * Math.Cos(s.ra - m.ra));
inc = Math.Atan2(sdist * Math.Sin(phi), m.dist - sdist * Math.Cos(phi));
angle = Math.Atan2(Math.Cos(s.dec) * Math.Sin(s.ra - m.ra), Math.Sin(s.dec) * Math.Cos(m.dec) -
Math.Cos(s.dec) * Math.Sin(m.dec) * Math.Cos(s.ra - m.ra));
double endPhase = 0.5 + 0.5 * inc * (angle < 0 ? -1 : 1) / Math.PI;
//Determine Moon Name.
if (startPhase <= .5 && endPhase >= .5)
{
moonDate = x;
moonName = GetMoonName(dMon.Month, moonName);
}
//Get Moon Name (month, string);
//Get Moon Phase Name
if (date.Day == x)
{
if (startPhase > endPhase)
{
mi.PhaseName = "New Moon";
break;
}
if (startPhase <= .25 && endPhase >= .25)
{
mi.PhaseName = "First Quarter";
break;
}
if (startPhase <= .5 && endPhase >= .5)
{
mi.PhaseName = "Full Moon";
break;
}
if (startPhase <= .75 && endPhase >= .75)
{
mi.PhaseName = "Last Quarter";
break;
}
if (startPhase > 0 && startPhase < .25 && endPhase > 0 && endPhase < .25)
{
mi.PhaseName = "Waxing Crescent";
break;
}
if (startPhase > .25 && startPhase < .5 && endPhase > .25 && endPhase < .5)
{
mi.PhaseName = "Waxing Gibbous";
break;
}
if (startPhase > .5 && startPhase < .75 && endPhase > .5 && endPhase < .75)
{
mi.PhaseName = "Waning Gibbous";
break;
}
if (startPhase > .75 && startPhase < 1 && endPhase > .75 && endPhase < 1)
{
mi.PhaseName = "Waning Crescent";
break;
}
}
}
if (date.Day == moonDate)
{
c.AstrologicalSigns.MoonName = moonName;
}
else { c.AstrologicalSigns.MoonName = ""; }
CalculateLunarEclipse(date, lat, lng, c);
}
public static void CalculateLunarEclipse(DateTime date, double lat, double longi, Celestial c)
{
//Convert to Radian
double latR = lat * Math.PI / 180;
double longR = longi * Math.PI / 180;
List<List<string>> se = LunarEclipseCalc.CalculateLunarEclipse(date, latR, longR);
//RETURN FIRST AND LAST
if (se.Count == 0) { return; }
//FIND LAST AND NEXT ECLIPSE
int lastE = -1;
int nextE = -1;
int currentE = 0;
DateTime lastDate = new DateTime();
DateTime nextDate = new DateTime(3300, 1, 1);
//Iterate to get last and next eclipse
foreach (List<string> values in se)
{
DateTime ld = DateTime.ParseExact(values[0], "yyyy-MMM-dd", System.Globalization.CultureInfo.InvariantCulture);
if (ld < date && ld > lastDate) { lastDate = ld; lastE = currentE; }
if (ld >= date && ld < nextDate) { nextDate = ld; nextE = currentE; }
currentE++;
}
//SET ECLIPSE DATA
if (lastE >= 0)
{
c.LunarEclipse.LastEclipse = new LunarEclipseDetails(se[lastE]);
}
if (nextE >= 0)
{
c.LunarEclipse.NextEclipse = new LunarEclipseDetails(se[nextE]);
}
}
private static string GetMoonName(int month, string name)
{
if (name != "") { return "Blue Moon"; }
switch (month)
{
case 1:
return "Wolf Moon";
case 2:
return "Snow Moon";
case 3:
return "Worm Moon";
case 4:
return "Pink Moon";
case 5:
return "Flower Moon";
case 6:
return "Strawberry Moon";
case 7:
return "Buck Moon";
case 8:
return "Sturgeon Moon";
case 9:
return "Corn Moon";
case 10:
return "Hunters Moon";
case 11:
return "Beaver Moon";
case 12:
return "Cold Moon";
default:
return "";
}
}
public static void GetMoonDistance(DateTime date, Celestial c)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
c.moonDistance = GetMoonDistance(date); //Updating distance formula
}
//Moon Time Functions
private static CelCoords GetSunCoords(double d)
{
double M = solarMeanAnomaly(d),
L = eclipticLongitude(M);
CelCoords c = new CelCoords();
c.dec = declination(L, 0);
c.ra = rightAscension(L, 0);
return c;
}
private static double solarMeanAnomaly(double d) { return rad * (357.5291 + 0.98560028 * d); }
private static double eclipticLongitude(double M)
{
double C = rad * (1.9148 * Math.Sin(M) + 0.02 * Math.Sin(2 * M) + 0.0003 * Math.Sin(3 * M)), // equation of center
P = rad * 102.9372; // perihelion of the Earth
return M + C + P + Math.PI;
}
public static void GetMoonSign(DateTime date, Celestial c)
{
//Formulas taken from https://www.astrocal.co.uk/moon-sign-calculator/
double d = date.Day;
double m = date.Month;
double y = date.Year;
double hr = date.Hour;
double mi = date.Minute;
double f = hr + (mi / 60);
double im = 12 * (y + 4800) + m - 3;
double j = (2 * (im - Math.Floor(im / 12) * 12) + 7 + 365 * im) / 12;
j = Math.Floor(j) + d + Math.Floor(im / 48) - 32083;
double jd = j + Math.Floor(im / 4800) - Math.Floor(im / 1200) + 38;
double T = ((jd - 2415020) + f / 24 - .5) / 36525;
double ob = FNr(23.452294 - .0130125 * T);
double ll = 973563 + 1732564379 * T - 4 * T * T;
double g = 1012395 + 6189 * T;
double n = 933060 - 6962911 * T + 7.5 * T * T;
double g1 = 1203586 + 14648523 * T - 37 * T * T;
d = 1262655 + 1602961611 * T - 5 * T * T;
double M = 3600;
double l = (ll - g1) / M;
double l1 = ((ll - d) - g) / M;
f = (ll - n) / M;
d = d / M;
y = 2 * d;
double ml = 22639.6 * FNs(l) - 4586.4 * FNs(l - y);
ml = ml + 2369.9 * FNs(y) + 769 * FNs(2 * l) - 669 * FNs(l1);
ml = ml - 411.6 * FNs(2 * f) - 212 * FNs(2 * l - y);
ml = ml - 206 * FNs(l + l1 - y) + 192 * FNs(l + y);
ml = ml - 165 * FNs(l1 - y) + 148 * FNs(l - l1) - 125 * FNs(d);
ml = ml - 110 * FNs(l + l1) - 55 * FNs(2 * f - y);
ml = ml - 45 * FNs(l + 2 * f) + 40 * FNs(l - 2 * f);
double tn = n + 5392 * FNs(2 * f - y) - 541 * FNs(l1) - 442 * FNs(y);
tn = tn + 423 * FNs(2 * f) - 291 * FNs(2 * l - 2 * f);
g = FNu(FNp(ll + ml));
double sign = Math.Floor(g / 30);
double degree = (g - (sign * 30));
sign = sign + 1;
switch (sign.ToString())
{
case "1": c.AstrologicalSigns.MoonSign = "Aries"; break;
case "2": c.AstrologicalSigns.MoonSign = "Taurus"; break;
case "3": c.AstrologicalSigns.MoonSign = "Gemini"; break;
case "4": c.AstrologicalSigns.MoonSign = "Cancer"; break;
case "5": c.AstrologicalSigns.MoonSign = "Leo"; break;
case "6": c.AstrologicalSigns.MoonSign = "Virgo"; break;
case "7": c.AstrologicalSigns.MoonSign = "Libra"; break;
case "8": c.AstrologicalSigns.MoonSign = "Scorpio"; break;
case "9": c.AstrologicalSigns.MoonSign = "Sagitarius"; break;
case "10": c.AstrologicalSigns.MoonSign = "Capricorn"; break;
case "11": c.AstrologicalSigns.MoonSign = "Aquarius"; break;
case "12": c.AstrologicalSigns.MoonSign = "Pisces"; break;
default: c.AstrologicalSigns.MoonSign = "Pisces"; break;
}
}
private static double FNp(double x)
{
double sgn;
if (x < 0)
{ sgn = -1; }
else
{ sgn = 1; }
return sgn * ((Math.Abs(x) / 3600) / 360 - Math.Floor((Math.Abs(x) / 3600.0) / 360.0)) * 360;
}
private static double FNu(double x)
{ return x - (Math.Floor(x / 360) * 360); }
private static double FNr(double x)
{ return Math.PI / 180 * x; }
private static double FNs(double x)
{ return Math.Sin(Math.PI / 180 * x); }
//v1.1.3 Formulas
//The following formulas are either additions
//or conversions of SunCalcs formulas into Meeus
/// <summary>
/// Grabs Perigee or Apogee of Moon based on specified time.
/// Results will return event just before, or just after specified DateTime
/// </summary>
/// <param name="d">DateTime</param>
/// <param name="md">Event Type</param>
/// <returns>PerigeeApogee</returns>
private static PerigeeApogee MoonPerigeeOrApogee(DateTime d, MoonDistanceType md)
{
//Perigee & Apogee Algorithms from Jean Meeus Astronomical Algorithms Ch. 50
//50.1
//JDE = 2451534.6698 + 27.55454989 * k
// -0.0006691 * Math.Pow(T,2)
// -0.000.01098 * Math.Pow(T,3)
// -0.0000000052 * Math.Pow(T,4)
//50.2
//K approx = (yv - 1999.97)*13.2555
//yv is the year + percentage of days that have occured in the year. 1998 Oct 1 is approx 1998.75
//k ending in .0 represent perigee and .5 apogee. Anything > .5 is an error.
//50.3
//T = k/1325.55
double yt = 365; //days in year
if (DateTime.IsLeapYear(d.Year)) { yt = 366; } //days in year if leap year
double f = d.DayOfYear / yt; //Get percentage of year that as passed
double yv = d.Year + f; //add percentage of year passed to year.
double k = (yv - 1999.97) * 13.2555; //find approximate k using formula 50.2
//Set k decimal based on apogee or perigee
if (md == MoonDistanceType.Apogee)
{
k = Math.Floor(k) + .5;
}
else
{
k = Math.Floor(k);
}
//Find T using formula 50.3
double T = k / 1325.55;
//Find JDE using formula 50.1
double JDE = 2451534.6698 + 27.55454989 * k -
0.0006691 * Math.Pow(T, 2) -
0.00001098 * Math.Pow(T, 3) -
0.0000000052 * Math.Pow(T, 4);
//Find Moon's mean elongation at time JDE.
double D = 171.9179 + 335.9106046 * k -
0.0100383 * Math.Pow(T, 2) -
0.00001156 * Math.Pow(T, 3) +
0.000000055 * Math.Pow(T, 4);
//Find Sun's mean anomaly at time JDE
double M = 347.3477 + 27.1577721 * k -
0.0008130 * Math.Pow(T, 2) -
0.0000010 * Math.Pow(T, 3);
//Find Moon's argument of latitude at Time JDE
double F = 316.6109 + 364.5287911 * k -
0.0125053 * Math.Pow(T, 2) -
0.0000148 * Math.Pow(T, 3);
//Normalize DMF to a 0-360 degree number
D %= 360;
if (D < 0) { D += 360; }
M %= 360;
if (M < 0) { M += 360; }
F %= 360;
if (F < 0) { F += 360; }
//Convert DMF to radians
D = D * Math.PI / 180;
M = M * Math.PI / 180;
F = F * Math.PI / 180;
double termsA;
//Find Terms A from Table 50.A
if (md == MoonDistanceType.Apogee)
{
termsA = MeeusTables.ApogeeTermsA(D, M, F, T);
}
else
{
termsA = MeeusTables.PerigeeTermsA(D, M, F, T);
}
JDE += termsA;
double termsB;
if (md == MoonDistanceType.Apogee)
{
termsB = MeeusTables.ApogeeTermsB(D, M, F, T);
}
else
{
termsB = MeeusTables.PerigeeTermsB(D, M, F, T);
}
//Convert julian back to date
DateTime date = JulianConversions.GetDate_FromJulian(JDE).Value;
//Obtain distance
Distance dist = GetMoonDistance(date);
PerigeeApogee ap = new PerigeeApogee(date, termsB, dist);
return ap;
}
public static Perigee GetPerigeeEvents(DateTime d)
{
//Iterate in 15 day increments due to formula variations.
//Determine closest events to date.
//per1 is last date
//per2 is next date
//integrity for new date.
if (d.Year <= 0001) { return new Perigee(new PerigeeApogee(new DateTime(), 0, new Distance(0)), new PerigeeApogee(new DateTime(), 0, new Distance(0))); }
//Start at lowest increment
PerigeeApogee per1 = MoonPerigeeOrApogee(d.AddDays(-45), MoonDistanceType.Perigee);
PerigeeApogee per2 = MoonPerigeeOrApogee(d.AddDays(-45), MoonDistanceType.Perigee);
for (int x = -30; x <= 45; x+=15)
{
//used for comparison
PerigeeApogee t = MoonPerigeeOrApogee(d.AddDays(x), MoonDistanceType.Perigee);
//Find the next pergiee after specified date
if (t.Date > per2.Date && t.Date >= d)
{
per2 = t;
break;
}
//Find last perigee before specified date
if (t.Date > per1.Date && t.Date < d)
{
per1 = t;
per2 = t;
}
}
return new Perigee(per1, per2);
}
public static Apogee GetApogeeEvents(DateTime d)
{
//Iterate in 5 month increments due to formula variations.
//Determine closest events to date.
//apo1 is last date
//apo2 is next date
//integrity for new date.
if (d.Year <= 0001) { return new Apogee(new PerigeeApogee(new DateTime(), 0, new Distance(0)), new PerigeeApogee(new DateTime(), 0, new Distance(0))); }
PerigeeApogee apo1 = MoonPerigeeOrApogee(d.AddDays(-45), MoonDistanceType.Apogee);
PerigeeApogee apo2 = MoonPerigeeOrApogee(d.AddDays(-45), MoonDistanceType.Apogee);
for (int x = -30; x <= 45; x+=15)
{
PerigeeApogee t = MoonPerigeeOrApogee(d.AddDays(x), MoonDistanceType.Apogee);
//Find next apogee after specified date
if (t.Date > apo2.Date && t.Date >= d)
{
apo2 = t;
break;
}
//Find last apogee before specified date
if (t.Date > apo1.Date && t.Date < d)
{
apo1 = t;
apo2 = t;
}
}
return new Apogee(apo1, apo2);
}
/// <summary>
/// Gets moon distance (Ch 47).
/// </summary>
/// <param name="d">DateTime</param>
/// <returns>Distance</returns>
public static Distance GetMoonDistance(DateTime d)
{
//Ch 47
double JDE = JulianConversions.GetJulian(d);//Get julian
double T = (JDE - 2451545) / 36525; //Get dynamic time.
double[] values = Get_Moon_LDMNF(T);
double D = values[1];
double M = values[2];
double N = values[3];
double F = values[4];
//Ch 47 distance formula
double dist = 385000.56 + (MeeusTables.Moon_Periodic_Er(D, M, N, F, T) / 1000);
return new Distance(dist);
}
private static Distance GetMoonDistance(DateTime d, double[] values)
{
//Ch 47
double JDE = JulianConversions.GetJulian(d);//Get julian
double T = (JDE - 2451545) / 36525; //Get dynamic time.
double D = values[1];
double M = values[2];
double N = values[3];
double F = values[4];
double dist = 385000.56 + (MeeusTables.Moon_Periodic_Er(D, M, N, F, T) / 1000);
return new Distance(dist);
}
/// <summary>
/// Gets Moon L, D, M, N, F values
/// Ch. 47
/// </summary>
/// <param name="T">Dynamic Time</param>
/// <returns>double[] containing L,D,M,N,F</returns>
static double[] Get_Moon_LDMNF(double T)
{
//T = dynamic time
//Moon's mean longitude
double L = 218.316447 + 481267.88123421 * T -
.0015786 * Math.Pow(T, 2) + Math.Pow(T, 3) / 538841 -
Math.Pow(T, 4) / 65194000;
//Moon's mean elongation
double D = 297.8501921 + 445267.1114034 * T -
0.0018819 * Math.Pow(T, 2) + Math.Pow(T, 3) / 545868 - Math.Pow(T, 4) / 113065000;
//Sun's mean anomaly
double M = 357.5291092 + 35999.0502909 * T -
.0001536 * Math.Pow(T, 2) + Math.Pow(T, 3) / 24490000;
//Moon's mean anomaly
double N = 134.9633964 + 477198.8675055 * T + .0087414 * Math.Pow(T, 2) +
Math.Pow(T, 3) / 69699 - Math.Pow(T, 4) / 14712000;
//Moon's argument of latitude
double F = 93.2720950 + 483202.0175233 * T - .0036539 * Math.Pow(T, 2) - Math.Pow(T, 3) /
3526000 + Math.Pow(T, 4) / 863310000;
//Normalize DMF to a 0-360 degree number
D %= 360;
if (D < 0) { D += 360; }
M %= 360;
if (M < 0) { M += 360; }
N %= 360;
if (N < 0) { N += 360; }
F %= 360;
if (F < 0) { F += 360; }
//Convert DMF to radians
D = D * Math.PI / 180;
M = M * Math.PI / 180;
N = N * Math.PI / 180;
F = F * Math.PI / 180;
return new double[] { L, D, M, N, F };
}
/// <summary>
/// Get moons lat/long in radians (Ch 47).
/// </summary>
/// <param name="LDMNF">L,D,M,N,F</param>
/// <param name="T">Dynamic Time</param>
/// <returns>Lat[0], Long[1]</returns>
private static double[] Get_Moon_Coordinates(double[] LDMNF,double T)
{
//Refence Ch 47.
double lat = LDMNF[0] + (MeeusTables.Moon_Periodic_El(LDMNF[0], LDMNF[1], LDMNF[2], LDMNF[3], LDMNF[4],T)/1000000);
double longi = MeeusTables.Moon_Periodic_Eb(LDMNF[0], LDMNF[1], LDMNF[2], LDMNF[3], LDMNF[4], T) / 1000000;
lat %= 360;
if (lat < 0) { lat += 360; }
//Convert to radians
double l = rad * lat; // longitude
double b = rad * longi; // latitude
return new double[] { l, b };
}
/// <summary>
/// Gets right Ascension of celestial object (Ch 13 Fig 13.3)
/// </summary>
/// <param name="l">latitude in radians</param>
/// <param name="b">longitude in radian</param>
/// <returns>Right Ascension</returns>
private static double rightAscension(double l, double b)
{
//Ch 13 Fig 13.3
//tan a = ( sin(l) * cos(e) - tan(b)-sin(e) ) / cons(l)
//Converts to the following using Atan2 for 4 quadriatic regions
return Math.Atan2(Math.Sin(l) * Math.Cos(e) - Math.Tan(b) * Math.Sin(e), Math.Cos(l));
}
/// <summary>
/// Gets declination of celestial object (Ch 13 Fig 13.4)
/// </summary>
/// <param name="l">latitude in radians</param>
/// <param name="b">longitude in radian</param>
/// <returns>Declination</returns>
private static double declination(double l, double b)
{
//Ch 13 Fig 13.4
//sin o = sin(b) * cos(e) + cos(b)*sin(e) * sin(l)
//Converts to the following using Asin
return Math.Asin(Math.Sin(b) * Math.Cos(e) + Math.Cos(b) * Math.Sin(e) * Math.Sin(l));
}
static double Parallax_Dec(double distance, double H, double pCosE, double pSinE, double dec, double cRA)
{
//Ch 40 (Correction for parallax
//H - geocentric hour angle of the body (sidereal) IAW Ch 12
double pi = Math.Asin((Math.Sin(8.794 / distance))) * Math.PI / 180; // 40.1 in radians
H = H * Math.PI / 180;
//Directly to topocencric dec
double tDEC = Math.Atan2((Math.Sin(dec) - pSinE * Math.Sin(pi)) * Math.Cos(cRA), Math.Cos(dec) - pCosE * Math.Sin(pi) * Math.Cos(H));
return tDEC;
}
static double Parallax_RA(double distance, double H, double pCosE, double dec, double ra)
{
//ENSURE RADIANS
//Ch 40 (Correction for parallax
//H - geocentric hour angle of the body (sidereal) IAW Ch 12
double pi = Math.Asin((Math.Sin(8.794 / distance))) * Math.PI / 180; // 40.1
//Convert to Radian
double t = -pCosE * Math.Sin(pi) * Math.Sin(H);
double b = Math.Cos(dec) - pCosE * Math.Sin(pi) * Math.Cos(H);
double cRA = Math.Atan2(t, b);
return cRA;
//Topocencric RA = RA - cRA
}
static double Get_pSinE(double dec, double H)
{
//ASSUME WGS 84 FOR NOW
double a = 6378.14;
double f = 1 / 298.257;
double b = a * (1 - f);
double ba = .99664719; // or 1-f
double u = (ba * dec) * Math.PI / 180;
double ps = ba * Math.Sin(u) + (H / 6378140) * Math.Sin(dec);
return ps;
}
static double Get_pCosE(double dec, double H)
{
//ASSUME WGS 84 FOR NOW
double a = 6378.14;
double f = 1 / 298.257;
double b = a * (1 - f);
double ba = .99664719; // or 1-f
double u = (ba * dec) * Math.PI / 180;
double ps = Math.Cos(u) + (H / 6378140) * Math.Cos(dec);
return ps;
}
static double azimuth(double H, double phi, double dec) { return Math.Atan2(Math.Sin(H), Math.Cos(H) * Math.Sin(phi) - Math.Tan(dec) * Math.Cos(phi)); }
static double altitude(double H, double phi, double dec)
{
return Math.Asin(Math.Sin(phi) * Math.Sin(dec) + Math.Cos(phi) * Math.Cos(dec) * Math.Cos(H));
}
static double astroRefraction(double h)
{
//CH 16
double P = 1013.25; //Average pressure of earth
double T = 16; //Average temp of earth
double alt = h / Math.PI * 180;
double Ref = P * (.1594 + .0196 * alt + .00002 * Math.Pow(alt, 2)) / ((273 + T) * (1 + .505 * alt + .0845 * Math.Pow(alt, 2)));
return Ref / 60;
}
}
}

971
Lora-Map/Helper/Celestial.SolarEclipseCalc.cs Normal file
View File

@ -0,0 +1,971 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace CoordinateSharp
{
//CURRENT ALTITUDE IS SET CONSTANT AT 100M. POSSIBLY NEED TO ADJUST TO ALLOW USER PASS.
//Altitude adjustments appear to have minimal effect on eclipse timing. These were mainly used
//to signify eclipses that had already started during rise and set times on the NASA calculator
//SOME TIMES AND ALTS WERE RETURNED WITH COLOR AND STYLING. DETERMINE WHY AND ADJUST VALUE AS REQUIRED. SEARCH "WAS ITALIC".
//ELLIPSOID ADJUSTMENT
//6378140.0 Ellipsoid is used in the NASA Calculator
//WGS84 Ellipsoid is 6378137.0. Adjustments to the ellipsoid appear to effect eclipse seconds in fractions.
//This can be modified if need to allow users to pass custom number with the Coordinate SetDatum() functions.
//CURRENT RANGE 1601-2600.
internal class SolarEclipseCalc
{
public static List<List<string>> CalculateSolarEclipse(DateTime d, double latRad, double longRad)
{
return Calculate(d, latRad, longRad, null);
}
public static List<SolarEclipseDetails> CalculateSolarEclipse(DateTime d, double latRad, double longRad, double[] events)
{
List<List<string>> evs = Calculate(d, latRad, longRad, events);
List<SolarEclipseDetails> deetsList = new List<SolarEclipseDetails>();
foreach(List<string> ls in evs)
{
SolarEclipseDetails deets = new SolarEclipseDetails(ls);
deetsList.Add(deets);
}
return deetsList;
}
public static List<List<string>> CalculateSolarEclipse(DateTime d, Coordinate coord)
{
return Calculate(d, coord.Latitude.ToRadians(), coord.Longitude.ToRadians(), null);
}
private static List<List<string>> Calculate(DateTime d, double latRad, double longRad, double[] ev)
{
//Declare storage arrays
double[] obsvconst = new double[7];
double[] mid = new double[41];//Check index to see if array needs to be this size
double[] c1 = new double[41];
double[] c2 = new double[41];
double[] c3 = new double[41];
double[] c4 = new double[41];
List<List<string>> events = new List<List<string>>();
double[] el;
if (ev == null)
{
el = Eclipse.SolarData.SolarDateData(d);//Get 100 year solar data;
}
else
{
el = ev;
}
events = new List<List<string>>();
ReadData(latRad, longRad, obsvconst);
for (int i = 0; i < el.Length; i += 28)
{
obsvconst[6] = i;
GetAll(el, obsvconst, mid,c1,c2,c3,c4);
// Is there an event...
if (mid[39] > 0)
{
List<string> values = new List<string>();
values.Add(GetDate(el, mid, obsvconst));
if (mid[39] == 1)
{
values.Add("P");
}
else if (mid[39] == 2)
{
values.Add("A");
}
else
{
values.Add("T");
}
// Partial eclipse start
if (c1[40] == 4)
{
values.Add("-");
values.Add(" ");
}
else
{
// Partial eclipse start time
values.Add(GetTime(el, c1, obsvconst));
values.Add(GetAlt(c1));
}
// Central eclipse time
if ((mid[39] > 1) && (c2[40] != 4))
{
values.Add(GetTime(el, c2, obsvconst));
}
else
{
values.Add("-");
}
//Mid Time
values.Add(GetTime(el, mid, obsvconst));
// Maximum eclipse alt
values.Add(GetAlt(mid));
// Maximum eclipse azi
values.Add(GetAzi(mid));
// Central eclipse ends
if ((mid[39] > 1) && (c3[40] != 4))
{
values.Add(GetTime(el, c3, obsvconst));
}
else
{
values.Add("-");
}
// Partial eclipse ends
if (c4[40] == 4)
{
values.Add("-");
values.Add(" ");
}
else
{
// Partial eclipse ends
values.Add(GetTime(el, c4, obsvconst));
// ... sun alt
values.Add(GetAlt(c4));
}
// Eclipse magnitude
values.Add(GetMagnitude(mid));
// Coverage
values.Add(GetCoverage(mid));
// Central duration
if (mid[39] > 1)
{
values.Add(GetDuration(mid,c2,c3));
}
else
{
values.Add("-");
}
events.Add(values);
}
}
return events;
}
//Populates the obsvcont array
private static void ReadData(double latRad, double longRad, double[] obsvconst)
{
// Get the latitude
obsvconst[0] = latRad;
//// Get the longitude
obsvconst[1] = -1 * longRad; //PASS REVERSE RADIAN.
// Get the altitude
obsvconst[2] = 100; //CHANGE TO ALLOW USER TO PASS.
// Get the time zone
obsvconst[3] = 0; //ALWAYS GMT
// Get the observer's geocentric position
double tmp = Math.Atan(0.99664719 * Math.Tan(obsvconst[0]));
obsvconst[4] = 0.99664719 * Math.Sin(tmp) + (obsvconst[2] / 6378140.0) * Math.Sin(obsvconst[0]);
obsvconst[5] = Math.Cos(tmp) + (obsvconst[2] / 6378140.0 * Math.Cos(obsvconst[0]));
}
// Populate the c1, c2, mid, c3 and c4 arrays
private static void GetAll(double[] elements, double[] obsvconst, double[] mid, double[] c1, double[] c2,double[] c3, double[] c4)
{
GetMid(elements, obsvconst, mid);
MidObservational(obsvconst, mid);
if (mid[37] > 0.0)
{
Getc1c4(elements, obsvconst, mid,c1,c2,c3,c4);
if ((mid[36] < mid[29]) || (mid[36] < -mid[29]))
{
Getc2c3(elements, obsvconst, mid,c2,c3);
if (mid[29] < 0.0)
{
mid[39] = 3; // Total eclipse
}
else
{
mid[39] = 2; // Annular eclipse
}
Observational(c1, obsvconst, mid);
Observational(c2, obsvconst, mid);
Observational(c3, obsvconst, mid);
Observational(c4, obsvconst, mid);
c2[36] = 999.9;
c3[36] = 999.9;
// Calculate how much of the eclipse is above the horizon
double pattern = 0;
if (c1[40] == 0) { pattern += 10000; }
if (c2[40] == 0) { pattern += 1000; }
if (mid[40] == 0) { pattern += 100; }
if (c3[40] == 0) { pattern += 10; }
if (c4[40] == 0) { pattern += 1; }
// Now, time to make sure that all my Observational[39] and Observational[40] are OK
if (pattern == 11110)
{
GetSunset(elements, c4, obsvconst);
Observational(c4, obsvconst, mid);
c4[40] = 3;
}
else if (pattern == 11100)
{
GetSunset(elements, c3, obsvconst);
Observational(c3, obsvconst, mid);
c3[40] = 3;
CopyCircumstances(c3, c4);
}
else if (pattern == 11000)
{
c3[40] = 4;
GetSunset(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 3;
CopyCircumstances(mid, c4);
}
else if (pattern == 10000)
{
mid[39] = 1;
GetSunset(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 3;
CopyCircumstances(mid, c4);
}
else if (pattern == 1111)
{
GetSunrise(elements, c1, obsvconst);
Observational(c1, obsvconst, mid);
c1[40] = 2;
}
else if (pattern == 111)
{
GetSunrise(elements, c2, obsvconst);
Observational(c2, obsvconst, mid);
c2[40] = 2;
CopyCircumstances(c2, c1);
}
else if (pattern == 11)
{
c2[40] = 4;
GetSunrise(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 2;
CopyCircumstances(mid, c1);
}
else if (pattern == 1)
{
mid[39] = 1;
GetSunrise(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 2;
CopyCircumstances(mid, c1);
}
else if (pattern == 0)
{
mid[39] = 0;
}
// There are other patterns, but those are the only ones we're covering!
}
else
{
mid[39] = 1; // Partial eclipse
double pattern = 0;
Observational(c1, obsvconst, mid);
Observational(c4, obsvconst, mid);
if (c1[40] == 0) { pattern += 100; }
if (mid[40] == 0) { pattern += 10; }
if (c4[40] == 0) { pattern += 1; }
if (pattern == 110)
{
GetSunset(elements, c4, obsvconst);
Observational(c4, obsvconst, mid);
c4[40] = 3;
}
else if (pattern == 100)
{
GetSunset(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 3;
CopyCircumstances(mid, c4);
}
else if (pattern == 11)
{
GetSunrise(elements, c1, obsvconst);
Observational(c1, obsvconst, mid);
c1[40] = 2;
}
else if (pattern == 1)
{
GetSunrise(elements, mid, obsvconst);
MidObservational(obsvconst, mid);
mid[40] = 2;
CopyCircumstances(mid, c1);
}
else if (pattern == 0)
{
mid[39] = 0;
}
// There are other patterns, but those are the only ones we're covering!
}
}
else
{
mid[39] = 0; // No eclipse
}
// Magnitude for total and annular eclipse is moon/sun ratio
if ((mid[39] == 2) || (mid[39] == 3))
{
mid[37] = mid[38];
}
}
// Calculate mid eclipse
private static void GetMid(double[] elements, double[] obsvconst, double[] mid)
{
double iter, tmp;
mid[0] = 0;
mid[1] = 0.0;
iter = 0;
tmp = 1.0;
TimeLocDependent(elements, mid, obsvconst);
while (((tmp > 0.000001) || (tmp < -0.000001)) && (iter < 50))
{
tmp = (mid[24] * mid[26] + mid[25] * mid[27]) / mid[30];
mid[1] = mid[1] - tmp;
iter++;
TimeLocDependent(elements, mid, obsvconst);
}
}
// Populate the circumstances array with the time and location dependent circumstances
private static double[] TimeLocDependent(double[] elements, double[] circumstances, double[] obsvconst)
{
double index, type;
TimeDependent(elements, circumstances, obsvconst);
index = obsvconst[6];
// Calculate h, sin h, cos h
circumstances[16] = circumstances[7] - obsvconst[1] - (elements[(int)index + 5] / 13713.44);
circumstances[17] = Math.Sin(circumstances[16]);
circumstances[18] = Math.Cos(circumstances[16]);
// Calculate xi
circumstances[19] = obsvconst[5] * circumstances[17];
// Calculate eta
circumstances[20] = obsvconst[4] * circumstances[6] - obsvconst[5] * circumstances[18] * circumstances[5];
// Calculate zeta
circumstances[21] = obsvconst[4] * circumstances[5] + obsvconst[5] * circumstances[18] * circumstances[6];
// Calculate dxi
circumstances[22] = circumstances[13] * obsvconst[5] * circumstances[18];
// Calculate deta
circumstances[23] = circumstances[13] * circumstances[19] * circumstances[5] - circumstances[21] * circumstances[12];
// Calculate u
circumstances[24] = circumstances[2] - circumstances[19];
// Calculate v
circumstances[25] = circumstances[3] - circumstances[20];
// Calculate a
circumstances[26] = circumstances[10] - circumstances[22];
// Calculate b
circumstances[27] = circumstances[11] - circumstances[23];
// Calculate l1'
type = circumstances[0];
if ((type == -2) || (type == 0) || (type == 2))
{
circumstances[28] = circumstances[8] - circumstances[21] * elements[26 + (int)index];
}
// Calculate l2'
if ((type == -1) || (type == 0) || (type == 1))
{
circumstances[29] = circumstances[9] - circumstances[21] * elements[27 + (int)index];
}
// Calculate n^2
circumstances[30] = circumstances[26] * circumstances[26] + circumstances[27] * circumstances[27];
return circumstances;
}
// Populate the circumstances array with the time-only dependent circumstances (x, y, d, m, ...)
private static double[] TimeDependent(double[] elements, double[] circumstances, double[] obsvconst)
{
double type, t, ans;
t = circumstances[1];
int index = (int)obsvconst[6];
// Calculate x
ans = elements[9 + index] * t + elements[8 + index];
ans = ans * t + elements[7 + index];
ans = ans * t + elements[6 + index];
circumstances[2] = ans;
// Calculate dx
ans = 3.0 * elements[9 + index] * t + 2.0 * elements[8 + index];
ans = ans * t + elements[7 + index];
circumstances[10] = ans;
// Calculate y
ans = elements[13 + index] * t + elements[12 + index];
ans = ans * t + elements[11 + index];
ans = ans * t + elements[10 + index];
circumstances[3] = ans;
// Calculate dy
ans = 3.0 * elements[13 + index] * t + 2.0 * elements[12 + index];
ans = ans * t + elements[11 + index];
circumstances[11] = ans;
// Calculate d
ans = elements[16 + index] * t + elements[15 + index];
ans = ans * t + elements[14 + index];
ans = ans * Math.PI / 180.0;
circumstances[4] = ans;
// sin d and cos d
circumstances[5] = Math.Sin(ans);
circumstances[6] = Math.Cos(ans);
// Calculate dd
ans = 2.0 * elements[16 + index] * t + elements[15 + index];
ans = ans * Math.PI / 180.0;
circumstances[12] = ans;
// Calculate m
ans = elements[19 + index] * t + elements[18 + index];
ans = ans * t + elements[17 + index];
if (ans >= 360.0)
{
ans = ans - 360.0;
}
ans = ans * Math.PI / 180.0;
circumstances[7] = ans;
// Calculate dm
ans = 2.0 * elements[19 + index] * t + elements[18 + index];
ans = ans * Math.PI / 180.0;
circumstances[13] = ans;
// Calculate l1 and dl1
type = circumstances[0];
if ((type == -2) || (type == 0) || (type == 2))
{
ans = elements[22 + index] * t + elements[21 + index];
ans = ans * t + elements[20 + index];
circumstances[8] = ans;
circumstances[14] = 2.0 * elements[22 + index] * t + elements[21 + index];
}
// Calculate l2 and dl2
if ((type == -1) || (type == 0) || (type == 1))
{
ans = elements[25 + index] * t + elements[24 + index];
ans = ans * t + elements[23 + index];
circumstances[9] = ans;
circumstances[15] = 2.0 * elements[25 + index] * t + elements[24 + index];
}
return circumstances;
}
// Get the observational circumstances for mid eclipse
private static void MidObservational(double[] obsvconst, double[] mid)
{
Observational(mid, obsvconst, mid);
// Calculate m, magnitude and moon/sun
mid[36] = Math.Sqrt(mid[24] * mid[24] + mid[25] * mid[25]);
mid[37] = (mid[28] - mid[36]) / (mid[28] + mid[29]);
mid[38] = (mid[28] - mid[29]) / (mid[28] + mid[29]);
}
// Get the observational circumstances
private static void Observational(double[] circumstances, double[] obsvconst, double[] mid)
{
double contacttype, coslat, sinlat;
// We are looking at an "external" contact UNLESS this is a total eclipse AND we are looking at
// c2 or c3, in which case it is an INTERNAL contact! Note that if we are looking at mid eclipse,
// then we may not have determined the type of eclipse (mid[39]) just yet!
if (circumstances[0] == 0)
{
contacttype = 1.0;
}
else
{
if ((mid[39] == 3) && ((circumstances[0] == -1) || (circumstances[0] == 1)))
{
contacttype = -1.0;
}
else
{
contacttype = 1.0;
}
}
// Calculate p
circumstances[31] = Math.Atan2(contacttype * circumstances[24], contacttype * circumstances[25]);
// Calculate alt
sinlat = Math.Sin(obsvconst[0]);
coslat = Math.Cos(obsvconst[0]);
circumstances[32] = Math.Asin(circumstances[5] * sinlat + circumstances[6] * coslat * circumstances[18]);
// Calculate q
circumstances[33] = Math.Asin(coslat * circumstances[17] / Math.Cos(circumstances[32]));
if (circumstances[20] < 0.0)
{
circumstances[33] = Math.PI - circumstances[33];
}
// Calculate v
circumstances[34] = circumstances[31] - circumstances[33];
// Calculate azi
circumstances[35] = Math.Atan2(-1.0 * circumstances[17] * circumstances[6], circumstances[5] * coslat - circumstances[18] * sinlat * circumstances[6]);
// Calculate visibility
if (circumstances[32] > -0.00524)
{
circumstances[40] = 0;
}
else
{
circumstances[40] = 1;
}
}
// Get C1 and C4 data
// Entry conditions -
// 1. The mid array must be populated
// 2. The magnitude at mid eclipse must be > 0.0
private static void Getc1c4(double[] elements, double[] obsvconst, double[] mid, double[] c1, double[] c2, double[] c3, double[] c4)
{
double tmp, n;
n = Math.Sqrt(mid[30]);
tmp = mid[26] * mid[25] - mid[24] * mid[27];
tmp = tmp / n / mid[28];
tmp = Math.Sqrt(1.0 - tmp * tmp) * mid[28] / n;
c1[0] = -2;
c4[0] = 2;
c1[1] = mid[1] - tmp;
c4[1] = mid[1] + tmp;
c1c4iterate(elements, c1, obsvconst);
c1c4iterate(elements, c4, obsvconst);
}
// Iterate on C1 or C4
private static double[] c1c4iterate(double[] elements, double[] circumstances, double[] obsvconst)
{
double sign, iter, tmp, n;
TimeLocDependent(elements, circumstances, obsvconst);
if (circumstances[0] < 0)
{
sign = -1.0;
}
else
{
sign = 1.0;
}
tmp = 1.0;
iter = 0;
while (((tmp > 0.000001) || (tmp < -0.000001)) && (iter < 50))
{
n = Math.Sqrt(circumstances[30]);
tmp = circumstances[26] * circumstances[25] - circumstances[24] * circumstances[27];
tmp = tmp / n / circumstances[28];
tmp = sign * Math.Sqrt(1.0 - tmp * tmp) * circumstances[28] / n;
tmp = (circumstances[24] * circumstances[26] + circumstances[25] * circumstances[27]) / circumstances[30] - tmp;
circumstances[1] = circumstances[1] - tmp;
TimeLocDependent(elements, circumstances, obsvconst);
iter++;
}
return circumstances;
}
// Get C2 and C3 data
// Entry conditions -
// 1. The mid array must be populated
// 2. There must be either a total or annular eclipse at the location!
private static void Getc2c3(double[] elements, double[] obsvconst, double[] mid, double[] c2, double[] c3)
{
double tmp, n;
n = Math.Sqrt(mid[30]);
tmp = mid[26] * mid[25] - mid[24] * mid[27];
tmp = tmp / n / mid[29];
tmp = Math.Sqrt(1.0 - tmp * tmp) * mid[29] / n;
c2[0] = -1;
c3[0] = 1;
if (mid[29] < 0.0)
{
c2[1] = mid[1] + tmp;
c3[1] = mid[1] - tmp;
}
else
{
c2[1] = mid[1] - tmp;
c3[1] = mid[1] + tmp;
}
c2c3iterate(elements, c2, obsvconst, mid);
c2c3iterate(elements, c3, obsvconst, mid);
}
// Iterate on C2 or C3
private static double[] c2c3iterate(double[] elements, double[] circumstances, double[] obsvconst, double[] mid)
{
double sign, iter, tmp, n;
TimeLocDependent(elements, circumstances, obsvconst);
if (circumstances[0] < 0)
{
sign = -1.0;
}
else
{
sign = 1.0;
}
if (mid[29] < 0.0)
{
sign = -sign;
}
tmp = 1.0;
iter = 0;
while (((tmp > 0.000001) || (tmp < -0.000001)) && (iter < 50))
{
n = Math.Sqrt(circumstances[30]);
tmp = circumstances[26] * circumstances[25] - circumstances[24] * circumstances[27];
tmp = tmp / n / circumstances[29];
tmp = sign * Math.Sqrt(1.0 - tmp * tmp) * circumstances[29] / n;
tmp = (circumstances[24] * circumstances[26] + circumstances[25] * circumstances[27]) / circumstances[30] - tmp;
circumstances[1] = circumstances[1] - tmp;
TimeLocDependent(elements, circumstances, obsvconst);
iter++;
}
return circumstances;
}
// Get the date of an event
private static string GetDate(double[] elements, double[] circumstances, double[] obsvconst)
{
string[] month = new string[] { "Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec" };
double t, jd, a, b, c, d, e, index;
string ans = "";
index = obsvconst[6];
// Calculate the JD for noon (TDT) the day before the day that contains T0
jd = Math.Floor(elements[(int)index] - (elements[1 + (int)index] / 24.0));
// Calculate the local time (ie the offset in hours since midnight TDT on the day containing T0).
t = circumstances[1] + elements[1 + (int)index] - obsvconst[3] - (elements[4 + (int)index] - 0.5) / 3600.0;
if (t < 0.0)
{
jd--;
}
if (t >= 24.0)
{
jd++;
}
if (jd >= 2299160.0)
{
a = Math.Floor((jd - 1867216.25) / 36524.25);
a = jd + 1 + a - Math.Floor(a / 4.0);
}
else
{
a = jd;
}
b = a + 1525.0;
c = Math.Floor((b - 122.1) / 365.25);
d = Math.Floor(365.25 * c);
e = Math.Floor((b - d) / 30.6001);
d = b - d - Math.Floor(30.6001 * e);
if (e < 13.5)
{
e = e - 1;
}
else
{
e = e - 13;
}
double year;
if (e > 2.5)
{
ans = c - 4716 + "-";
year = c - 4716;
}
else
{
ans = c - 4715 + "-";
year = c - 4715;
}
string m = month[(int)e - 1];
ans += m + "-";
if (d < 10)
{
ans = ans + "0";
}
ans = ans + d;
//Leap Year Integrity Check
if(m =="Feb" && d ==29 && !DateTime.IsLeapYear((int)year))
{
ans = year.ToString() + "-Mar-01";
}
return ans;
}
// Calculate the time of sunset
private static void GetSunset(double[] elements, double[] circumstances, double[] obsvconst)
{
GetSunriset(elements, circumstances, 1.0, obsvconst);
}
// Calculate the time of sunrise
private static void GetSunrise(double[] elements, double[] circumstances, double[] obsvconst)
{
GetSunriset(elements, circumstances, -1.0, obsvconst);
}
// Calculate the time of sunrise or sunset
private static void GetSunriset(double[] elements, double[] circumstances, double riset, double[] obsvconst)
{
double h0, diff, iter;
diff = 1.0;
iter = 0;
while ((diff > 0.00001) || (diff < -0.00001))
{
iter++;
if (iter == 4) { return; }
h0 = Math.Acos((Math.Sin(-0.00524) - Math.Sin(obsvconst[0]) * circumstances[5]) / Math.Cos(obsvconst[0]) / circumstances[6]);
diff = (riset * h0 - circumstances[16]) / circumstances[13];
while (diff >= 12.0) { diff -= 24.0; }
while (diff <= -12.0) { diff += 24.0; }
circumstances[1] += diff;
TimeLocDependent(elements, circumstances, obsvconst);
}
}
// Copy a set of circumstances
private static void CopyCircumstances(double[] circumstancesfrom, double[] circumstancesto)
{
for (int i = 1; i < 41; i++)
{
circumstancesto[i] = circumstancesfrom[i];
}
}
// Get the local time of an event
private static string GetTime(double[] elements, double[] circumstances, double[] obsvconst)
{
string ans = "";
int index = (int)obsvconst[6];
double t = circumstances[1] + elements[1 + index] - obsvconst[3] - (elements[4 + index] - 0.5) / 3600.0;
if (t < 0.0)
{
t = t + 24.0;
}
if (t >= 24.0)
{
t = t - 24.0;
}
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + Math.Floor(t) + ":";
t = (t * 60.0) - 60.0 * Math.Floor(t);
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + Math.Floor(t);
if (circumstances[40] <= 1)
{ // not sunrise or sunset
ans = ans + ":";
t = (t * 60.0) - 60.0 * Math.Floor(t);
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + Math.Floor(t);
}
if (circumstances[40] == 1)
{
//WAS ITALIC
return ans;
}
else if (circumstances[40] == 2)
{
//Rise (CHANGED FROM NASA CALC THE INDICATES (r) WITH STRING, INVESTIGATE REMOVAL)
return ans;
}
else if (circumstances[40] == 3)
{
//Set (CHANGED FROM NASA CALC THE INDICATES (s) WITH STRING, INVESTIGATE REMOVAL)
return ans;
}
else
{
return ans;
}
}
// Get the altitude
private static string GetAlt(double[] circumstances)
{
double t;
string ans = "";
if (circumstances[40] == 2)
{
return "0(r)";
}
if (circumstances[40] == 3)
{
return "0(s)";
}
if ((circumstances[32] < 0.0) && (circumstances[32] >= -0.00524))
{
// Crude correction for refraction (and for consistency's sake)
t = 0.0;
}
else
{
t = circumstances[32] * 180.0 / Math.PI;
}
if (t < 0.0)
{
ans = "-";
t = -t;
}
else
{
ans = "";
}
t = Math.Floor(t + 0.5);
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + t;
if (circumstances[40] == 1)
{
//WAS ITALIC
return ans;
}
else
{
return ans;
}
}
// Get the azimuth
private static string GetAzi(double[] circumstances)
{
string ans = "";
double t = circumstances[35] * 180.0 / Math.PI;
if (t < 0.0)
{
t = t + 360.0;
}
if (t >= 360.0)
{
t = t - 360.0;
}
t = Math.Floor(t + 0.5);
if (t < 100.0)
{
ans = ans + "0";
}
if (t < 10.0)
{
ans = ans + "0";
}
ans = ans + t;
if (circumstances[40] == 1)
{
//WAS ITALIC
return ans;
}
else
{
return ans;
}
}
// Get the magnitude
private static string GetMagnitude(double[] mid)
{
double a = Math.Floor(1000.0 * mid[37] + 0.5) / 1000.0;
string ans = a.ToString();
if (mid[40] == 1)
{
return ans;
}
if (mid[40] == 2)
{
ans = a.ToString() + "(r)";
}
if (mid[40] == 3)
{
ans = a.ToString() + "(s)";
}
return ans;
}
// Get the coverage
private static string GetCoverage(double[] mid)
{
double a=0, b, c;
string ans = "";
if (mid[37] <= 0.0)
{
ans = "0.0";
}
else if (mid[37] >= 1.0)
{
ans = "1.000";
}
else
{
if (mid[39] == 2)
{
c = mid[38] * mid[38];
}
else
{
c = Math.Acos((mid[28] * mid[28] + mid[29] * mid[29] - 2.0 * mid[36] * mid[36]) / (mid[28] * mid[28] - mid[29] * mid[29]));
b = Math.Acos((mid[28] * mid[29] + mid[36] * mid[36]) / mid[36] / (mid[28] + mid[29]));
a = Math.PI - b - c;
c = ((mid[38] * mid[38] * a + b) - mid[38] * Math.Sin(c)) / Math.PI;
}
a = Math.Floor(1000.0 * c + 0.5) / 1000.0;
ans = a.ToString();
}
if (mid[40] == 1)
{
//WAS ITALIC
return ans;
}
if (mid[40] == 2)
{
ans = a.ToString() + "(r)";
}
if (mid[40] == 3)
{
ans = a + "(s)";
}
return ans;
}
// Get the duration in mm:ss.s format
// Adapted from code written by Stephen McCann - 27/04/2001
private static string GetDuration(double[] mid, double[] c2, double[] c3)
{
double tmp;
string ans;
if (c3[40] == 4)
{
tmp = mid[1] - c2[1];
}
else if (c2[40] == 4)
{
tmp = c3[1] - mid[1];
}
else
{
tmp = c3[1] - c2[1];
}
if (tmp < 0.0)
{
tmp = tmp + 24.0;
}
else if (tmp >= 24.0)
{
tmp = tmp - 24.0;
}
tmp = (tmp * 60.0) - 60.0 * Math.Floor(tmp) + 0.05 / 60.0;
ans = Math.Floor(tmp) + "m";
tmp = (tmp * 60.0) - 60.0 * Math.Floor(tmp);
if (tmp < 10.0)
{
ans = ans + "0";
}
ans += Math.Floor(tmp) + "s";
return ans;
}
}
}

374
Lora-Map/Helper/Celestial.SunCalculations.cs Normal file
View File

@ -0,0 +1,374 @@
using System;
using System.Collections.Generic;
namespace CoordinateSharp
{
internal class SunCalc
{
public static void CalculateSunTime(double lat, double longi, DateTime date, Celestial c,double offset = 0)
{
if (date.Year == 0001) { return; } //Return if date vaue hasn't been established.
DateTime actualDate = new DateTime(date.Year,date.Month,date.Day,0, 0, 0, DateTimeKind.Utc);
////Sun Time Calculations
//Get Julian
double lw = rad * -longi;
double phi = rad * lat;
//Rise Set
DateTime?[] evDate = Get_Event_Time(lw, phi, -.8333, actualDate);
c.sunRise = evDate[0];
c.sunSet = evDate[1];
c.sunCondition = CelestialStatus.RiseAndSet;
//Azimuth and Altitude
CalculateSunAngle(date, longi, lat, c);
// neither sunrise nor sunset
if ((!c.SunRise.HasValue) && (!c.SunSet.HasValue))
{
if (c.SunAltitude < 0)
{
c.sunCondition = CelestialStatus.DownAllDay;
}
else
{
c.sunCondition = CelestialStatus.UpAllDay;
}
}
// sunrise or sunset
else
{
if (!c.SunRise.HasValue)
{
// No sunrise this date
c.sunCondition = CelestialStatus.NoRise;
}
else if (!c.SunSet.HasValue)
{
// No sunset this date
c.sunCondition = CelestialStatus.NoSet;
}
}
//Additional Times
c.additionalSolarTimes = new AdditionalSolarTimes();
//Dusk and Dawn
//Civil
evDate = Get_Event_Time(lw, phi, -6, actualDate);
c.AdditionalSolarTimes.CivilDawn = evDate[0];
c.AdditionalSolarTimes.CivilDusk = evDate[1];
//Nautical
evDate = Get_Event_Time(lw, phi, -12, actualDate);
c.AdditionalSolarTimes.NauticalDawn = evDate[0];
c.AdditionalSolarTimes.NauticalDusk = evDate[1];
//Astronomical
evDate = Get_Event_Time(lw, phi, -18, actualDate);
c.AdditionalSolarTimes.AstronomicalDawn = evDate[0];
c.AdditionalSolarTimes.AstronomicalDusk = evDate[1];
//BottomDisc
evDate = Get_Event_Time(lw, phi, -.2998, actualDate);
c.AdditionalSolarTimes.SunriseBottomDisc = evDate[0];
c.AdditionalSolarTimes.SunsetBottomDisc = evDate[1];
CalculateSolarEclipse(date, lat, longi, c);
}
/// <summary>
/// Gets time of event based on specified degree below horizon
/// </summary>
/// <param name="lw">Observer Longitude in radians</param>
/// <param name="phi">Observer Latitude in radians</param>
/// <param name="h">Angle in Degrees</param>
/// <param name="date">Date of Event</param>
/// <returns>DateTime?[]{rise, set}</returns>
private static DateTime?[] Get_Event_Time(double lw, double phi, double h,DateTime date)
{
//Create arrays. Index 0 = Day -1, 1 = Day, 2 = Day + 1;
//These will be used to find exact day event occurs for comparison
DateTime?[] sets = new DateTime?[] { null, null, null, null, null };
DateTime?[] rises = new DateTime?[] { null, null, null,null, null };
//Iterate starting with day -1;
for (int x = 0; x < 5; x++)
{
double d = JulianConversions.GetJulian(date.AddDays(x-2)) - j2000 + .5; //LESS PRECISE JULIAN NEEDED
double n = julianCycle(d, lw);
double ds = approxTransit(0, lw, n);
double M = solarMeanAnomaly(ds);
double L = eclipticLongitude(M);
double dec = declination(L, 0);
double Jnoon = solarTransitJ(ds, M, L);
double Jset;
double Jrise;
DateTime? solarNoon = JulianConversions.GetDate_FromJulian(Jnoon);
DateTime? nadir = JulianConversions.GetDate_FromJulian(Jnoon - 0.5);
//Rise Set
Jset = GetTime(h * rad, lw, phi, dec, n, M, L);
Jrise = Jnoon - (Jset - Jnoon);
DateTime? rise = JulianConversions.GetDate_FromJulian(Jrise);
DateTime? set = JulianConversions.GetDate_FromJulian(Jset);
rises[x] = rise;
sets[x] = set;
}
//Compare and send
DateTime? tRise = null;
for(int x=0;x<5;x++)
{
if(rises[x].HasValue)
{
if(rises[x].Value.Day == date.Day)
{
tRise = rises[x];
break;
}
}
}
DateTime? tSet = null;
for (int x = 0; x < 5; x++)
{
if (sets[x].HasValue)
{
if (sets[x].Value.Day == date.Day)
{
tSet = sets[x];
break;
}
}
}
return new DateTime?[] { tRise, tSet };
}
public static void CalculateZodiacSign(DateTime date, Celestial c)
{
//Aquarius (January 20 to February 18)
//Pisces (February 19 to March 20)
//Aries (March 21-April 19)
//Taurus (April 20-May 20)
//Gemini (May 21-June 20)
//Cancer (June 21-July 22)
//Leo (July 23-August 22)
//Virgo (August 23-September 22)
//Libra (September 23-October 22)
//Scorpio (October 23-November 21)
//Sagittarius (November 22-December 21)
//Capricorn (December 22-January 19)
if (date >= new DateTime(date.Year, 1, 1) && date <= new DateTime(date.Year, 1, 19, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Capricorn";
return;
}
if (date >= new DateTime(date.Year, 1, 20) && date <= new DateTime(date.Year, 2, 18, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Aquarius";
return;
}
if (date >= new DateTime(date.Year, 2, 19) && date <= new DateTime(date.Year, 3, 20, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Pisces";
return;
}
if (date >= new DateTime(date.Year, 3, 21) && date <= new DateTime(date.Year, 4, 19, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Aries";
return;
}
if (date >= new DateTime(date.Year, 4, 20) && date <= new DateTime(date.Year, 5, 20, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Taurus";
return;
}
if (date >= new DateTime(date.Year, 5, 21) && date <= new DateTime(date.Year, 6, 20,23,59,59))
{
c.AstrologicalSigns.ZodiacSign = "Gemini";
return;
}
if (date >= new DateTime(date.Year, 6, 21) && date <= new DateTime(date.Year, 7, 22, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Cancer";
return;
}
if (date >= new DateTime(date.Year, 7, 23) && date <= new DateTime(date.Year, 8, 22, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Leo";
return;
}
if (date >= new DateTime(date.Year, 8, 23) && date <= new DateTime(date.Year, 9, 22, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Virgo";
return;
}
if (date >= new DateTime(date.Year, 9, 23) && date <= new DateTime(date.Year, 10, 22, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Libra";
return;
}
if (date >= new DateTime(date.Year, 9, 23) && date <= new DateTime(date.Year, 11, 21, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Scorpio";
return;
}
if (date >= new DateTime(date.Year, 11, 21) && date <= new DateTime(date.Year, 12, 21, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Sagittarius";
return;
}
if (date >= new DateTime(date.Year, 12, 22) && date <= new DateTime(date.Year, 12, 31, 23, 59, 59))
{
c.AstrologicalSigns.ZodiacSign = "Capricorn";
return;
}
}
public static void CalculateSolarEclipse(DateTime date, double lat, double longi, Celestial c)
{
//Convert to Radian
double latR = lat * Math.PI / 180;
double longR = longi * Math.PI / 180;
List<List<string>> se = SolarEclipseCalc.CalculateSolarEclipse(date, latR, longR);
//RETURN FIRST AND LAST
if (se.Count == 0) { return; }
//FIND LAST AND NEXT ECLIPSE
int lastE = -1;
int nextE = -1;
int currentE = 0;
DateTime lastDate = new DateTime();
DateTime nextDate = new DateTime(3300, 1, 1);
//Iterate to get last and next eclipse
foreach(List<string> values in se)
{
DateTime ld = DateTime.ParseExact(values[0], "yyyy-MMM-dd", System.Globalization.CultureInfo.InvariantCulture);
if (ld < date && ld>lastDate) { lastDate = ld;lastE = currentE; }
if(ld>= date && ld < nextDate) { nextDate = ld;nextE = currentE; }
currentE++;
}
//SET ECLIPSE DATA
if (lastE >= 0)
{
c.SolarEclipse.LastEclipse = new SolarEclipseDetails(se[lastE]);
}
if (nextE >= 0)
{
c.SolarEclipse.NextEclipse = new SolarEclipseDetails(se[nextE]);
}
}
#region Private Suntime Members
private static readonly double dayMS = 1000 * 60 * 60 * 24, j1970 = 2440588, j2000 = 2451545;
private static readonly double rad = Math.PI / 180;
private static double LocalSiderealTimeForTimeZone(double lon, double jd, double z)
{
double s = 24110.5 + 8640184.812999999 * jd / 36525 + 86636.6 * z + 86400 * lon;
s = s / 86400;
s = s - Math.Truncate(s);
double lst = s * 360 *rad;
return lst;
}
private static double SideRealTime(double d, double lw)
{
double s = rad * (280.16 + 360.9856235 * d) - lw;
return s;
}
private static double solarTransitJ(double ds, double M, double L)
{
return j2000 + ds + 0.0053 * Math.Sin(M) - 0.0069 * Math.Sin(2 * L);
}
//CH15
//Formula 15.1
//Returns Approximate Time
private static double hourAngle(double h, double phi, double d)
{
//NUMBER RETURNING > and < 1 NaN;
double d1 = Math.Sin(h) - Math.Sin(phi) * Math.Sin(d);
double d2 = Math.Cos(phi) * Math.Cos(d);
double d3 = (d1 / d2);
return Math.Acos(d3);
}
private static double approxTransit(double Ht, double lw, double n)
{
return 0.0009 + (Ht + lw) / (2 * Math.PI) + n;
}
private static double julianCycle(double d, double lw) { return Math.Round(d - 0.0009 - lw / (2 * Math.PI)); }
//Returns Time of specified event based on suns angle
private static double GetTime(double h, double lw, double phi, double dec, double n,double M, double L)
{
double approxTime = hourAngle(h, phi, dec); //Ch15 Formula 15.1
double a = approxTransit(approxTime, lw, n);
double st = solarTransitJ(a, M, L);
return st;
}
private static double declination(double l, double b)
{
double e = (Math.PI/180) * 23.4392911; // obliquity of the Earth
return Math.Asin(Math.Sin(b) * Math.Cos(e) + Math.Cos(b) * Math.Sin(e) * Math.Sin(l));
}
private static void CalculateSunAngle(DateTime date, double longi, double lat, Celestial c)
{
TimeSpan ts = date - new DateTime(1970, 1, 1,0,0,0, DateTimeKind.Utc);
double dms = (ts.TotalMilliseconds / dayMS -.5 + j1970)-j2000;
double lw = rad * -longi;
double phi = rad * lat;
double e = rad * 23.4397;
double[] sc = sunCoords(dms);
double H = SideRealTime(dms, lw) - sc[1];
c.sunAzimuth = Math.Atan2(Math.Sin(H), Math.Cos(H) * Math.Sin(phi) - Math.Tan(sc[0]) * Math.Cos(phi)) * 180 / Math.PI + 180;
c.sunAltitude = Math.Asin(Math.Sin(phi) * Math.Sin(sc[0]) + Math.Cos(phi) * Math.Cos(sc[0]) * Math.Cos(H)) * 180 / Math.PI;
}
private static double solarMeanAnomaly(double d)
{
return rad * (357.5291 + 0.98560028 * d);
}
private static double eclipticLongitude(double m)
{
double c = rad * (1.9148 * Math.Sin(m) + 0.02 * Math.Sin(2 * m) + 0.0003 * Math.Sin(3 * m)); // equation of center
double p = rad * 102.9372; // perihelion of the Earth
return m + c + p + Math.PI;
}
private static double[] sunCoords(double d)
{
double m = solarMeanAnomaly(d);
double l = eclipticLongitude(m);
double[] sc = new double[2];
double b = 0;
double e = rad * 23.4397; // obliquity of the Earth
sc[0] = Math.Asin(Math.Sin(b) * Math.Cos(e) + Math.Cos(b) * Math.Sin(e) * Math.Sin(l)); //declination
sc[1] = Math.Atan2(Math.Sin(l) * Math.Cos(e) - Math.Tan(b) * Math.Sin(e), Math.Cos(l)); //rightAscension
return sc;
}
#endregion
}
}

620
Lora-Map/Helper/Celestial.cs Normal file
View File

@ -0,0 +1,620 @@
using System;
using System.ComponentModel;
using System.Collections.Generic;
namespace CoordinateSharp
{
/// <summary>
/// The main class for handling location based celestial information.
/// </summary>
/// <remarks>
/// This class can calculate various pieces of solar and lunar data, based on location and date
/// </remarks>
[Serializable]
public class Celestial
{
//When a rise or a set does not occur, the DateTime will return null
/// <summary>
/// Creates an empty Celestial.
/// </summary>
public Celestial()
{
astrologicalSigns = new AstrologicalSigns();
lunarEclipse = new LunarEclipse();
solarEclipse = new SolarEclipse();
CalculateCelestialTime(0, 0, new DateTime(1900, 1, 1, 0, 0, 0, DateTimeKind.Utc));
}
private Celestial(bool hasCalcs)
{
astrologicalSigns = new AstrologicalSigns();
lunarEclipse = new LunarEclipse();
solarEclipse = new SolarEclipse();
if (hasCalcs) { CalculateCelestialTime(0, 0, new DateTime(1900, 1, 1, 0, 0, 0, DateTimeKind.Utc)); }
}
/// <summary>
/// Creates a Celestial based on a location and specified date
/// </summary>
/// <param name="lat">latitude</param>
/// <param name="longi">longitude</param>
/// <param name="geoDate">DateTime (UTC)</param>
public Celestial(double lat, double longi, DateTime geoDate)
{
DateTime d = new DateTime(geoDate.Year, geoDate.Month, geoDate.Day, geoDate.Hour, geoDate.Minute, geoDate.Second, DateTimeKind.Utc);
astrologicalSigns = new AstrologicalSigns();
lunarEclipse = new LunarEclipse();
solarEclipse = new SolarEclipse();
CalculateCelestialTime(lat, longi, d);
}
/// <summary>
/// Creates a Celestial based on a location and date in the provided Coordinate.
/// </summary>
/// <param name="c">Coordinate</param>
/// <returns>Celestial</returns>
public static Celestial LoadCelestial(Coordinate c)
{
DateTime geoDate = c.GeoDate;
DateTime d = new DateTime(geoDate.Year, geoDate.Month, geoDate.Day, geoDate.Hour, geoDate.Minute, geoDate.Second, DateTimeKind.Utc);
Celestial cel = new Celestial(c.Latitude.ToDouble(), c.Longitude.ToDouble(), c.GeoDate);
return cel;
}
/// <summary>
/// Converts Celestial values to local times.
/// </summary>
/// <param name="c">Coordinate</param>
/// <param name="offset">UTC offset</param>
/// <returns></returns>
public static Celestial Celestial_LocalTime(Coordinate c, double offset)
{
if(offset < -12 || offset > 12) { throw new ArgumentOutOfRangeException("Time offsets cannot be greater 12 or less than -12."); }
//Probably need to offset initial UTC date so user can op in local
//Determine best way to do this.
DateTime d = c.GeoDate.AddHours(offset);
//Get 3 objects for comparison
Celestial cel = new Celestial(c.Latitude.ToDouble(), c.Longitude.ToDouble(), c.GeoDate);
Celestial celPre = new Celestial(c.Latitude.ToDouble(), c.Longitude.ToDouble(), c.GeoDate.AddDays(-1));
Celestial celPost = new Celestial(c.Latitude.ToDouble(), c.Longitude.ToDouble(), c.GeoDate.AddDays(1));
//Slip objects for comparison. Compare with slipped date.
celPre.Local_Convert(c, offset);
cel.Local_Convert(c, offset);
celPost.Local_Convert(c, offset);
//Get SunSet
int i = Determine_Slipped_Event_Index(cel.SunSet, celPre.SunSet, celPost.SunSet, d);
cel.sunSet = Get_Correct_Slipped_Date(cel.SunSet, celPre.SunSet, celPost.SunSet, i);
cel.AdditionalSolarTimes.CivilDusk = Get_Correct_Slipped_Date(cel.AdditionalSolarTimes.CivilDusk,
celPre.AdditionalSolarTimes.CivilDusk, celPost.AdditionalSolarTimes.CivilDusk, i);
cel.AdditionalSolarTimes.NauticalDusk = Get_Correct_Slipped_Date(cel.AdditionalSolarTimes.NauticalDusk,
celPre.AdditionalSolarTimes.NauticalDusk, celPost.AdditionalSolarTimes.NauticalDusk, i);
//Get SunRise
i = Determine_Slipped_Event_Index(cel.SunRise, celPre.SunRise, celPost.SunRise, d);
cel.sunRise = Get_Correct_Slipped_Date(cel.SunRise, celPre.SunRise, celPost.SunRise, i);
cel.AdditionalSolarTimes.CivilDawn = Get_Correct_Slipped_Date(cel.AdditionalSolarTimes.CivilDawn,
celPre.AdditionalSolarTimes.CivilDawn, celPost.AdditionalSolarTimes.CivilDawn, i);
cel.AdditionalSolarTimes.NauticalDawn = Get_Correct_Slipped_Date(cel.AdditionalSolarTimes.NauticalDawn,
celPre.AdditionalSolarTimes.NauticalDawn, celPost.AdditionalSolarTimes.NauticalDawn, i);
//MoonRise
i = Determine_Slipped_Event_Index(cel.MoonRise, celPre.MoonRise, celPost.MoonRise, d);
cel.moonRise = Get_Correct_Slipped_Date(cel.MoonRise, celPre.MoonRise, celPost.MoonRise, i);
//MoonSet
i = Determine_Slipped_Event_Index(cel.MoonSet, celPre.MoonSet, celPost.MoonSet, d);
cel.moonSet = Get_Correct_Slipped_Date(cel.MoonSet, celPre.MoonSet, celPost.MoonSet, i);
//Local Conditions
CelestialStatus[] cels = new CelestialStatus[]
{
celPre.MoonCondition,cel.MoonCondition,celPost.MoonCondition
};
cel.moonCondition = Celestial.GetStatus(cel.MoonRise, cel.MoonSet, cels);
cels = new CelestialStatus[]
{
celPre.SunCondition, cel.SunCondition, celPost.SunCondition
};
cel.sunCondition = Celestial.GetStatus(cel.SunRise, cel.SunSet, cels);
//Load IsUp values based on local time with populated Celestial
Celestial.Calculate_Celestial_IsUp_Booleans(d, cel);
return cel;
}
private static CelestialStatus GetStatus(DateTime? rise, DateTime? set, CelestialStatus[] cels)
{
if (set.HasValue && rise.HasValue) { return CelestialStatus.RiseAndSet; }
if (set.HasValue && !rise.HasValue) { return CelestialStatus.NoRise; }
if (!set.HasValue && rise.HasValue) { return CelestialStatus.NoSet; }
for (int x=0; x < 3;x++)
{
if(cels[x] == CelestialStatus.DownAllDay || cels[x] == CelestialStatus.UpAllDay)
{
return cels[x];
}
}
return cels[1];
}
/// <summary>
/// In place time slip
/// </summary>
/// <param name="c">Coordinate</param>
/// <param name="offset">hour offset</param>
private void Local_Convert(Coordinate c, double offset)
{
//Find new lunar set rise times
if (MoonSet.HasValue) { moonSet = moonSet.Value.AddHours(offset); }
if (MoonRise.HasValue) { moonRise = moonRise.Value.AddHours(offset); }
//Perigee
Perigee.ConvertTo_Local_Time(offset);
//Apogee
Apogee.ConvertTo_Local_Time(offset);
//Eclipse
LunarEclipse.ConvertTo_LocalTime(offset);
////Solar
if (sunSet.HasValue) { sunSet = sunSet.Value.AddHours(offset); }
if (SunRise.HasValue) { sunRise = SunRise.Value.AddHours(offset); }
AdditionalSolarTimes.Convert_To_Local_Time(offset);
//Eclipse
SolarEclipse.ConvertTo_LocalTime(offset);
SunCalc.CalculateZodiacSign(c.GeoDate.AddHours(offset), this);
MoonCalc.GetMoonSign(c.GeoDate.AddHours(offset), this);
}
private PerigeeApogee Get_Correct_Slipped_Date(PerigeeApogee actual, PerigeeApogee pre, PerigeeApogee post, int i)
{
switch (i)
{
case 0:
return pre;
case 1:
return actual;
case 2:
return post;
default:
return actual;
}
}
private static DateTime? Get_Correct_Slipped_Date(DateTime? actual, DateTime? pre, DateTime? post, int i)
{
switch(i)
{
case 0:
return pre;
case 1:
return actual;
case 2:
return post;
default:
return null;
}
}
private static int Determine_Slipped_Event_Index(DateTime? actual, DateTime? pre, DateTime? post, DateTime d)
{
if (actual.HasValue)
{
if (actual.Value.Day != d.Day)
{
if (pre.HasValue)
{
if (pre.Value.Day == d.Day) { return 0; }
}
if (post.HasValue)
{
if (post.Value.Day == d.Day) { return 2; }
}
return 3;
}
}
else
{
if (pre.HasValue)
{
if (pre.Value.Day == d.Day) { return 0; }
}
if (post.HasValue)
{
if (post.Value.Day == d.Day) { return 2; }
}
}
return 1;
}
internal DateTime? sunSet;
internal DateTime? sunRise;
internal DateTime? moonSet;
internal DateTime? moonRise;
internal double sunAltitude;
internal double sunAzimuth;
internal double moonAltitude;
internal double moonAzimuth;
internal Distance moonDistance;
internal CelestialStatus sunCondition;
internal CelestialStatus moonCondition;
internal bool isSunUp;
internal bool isMoonUp;
internal MoonIllum moonIllum;
internal Perigee perigee;
internal Apogee apogee;
internal AdditionalSolarTimes additionalSolarTimes;
internal AstrologicalSigns astrologicalSigns;
internal SolarEclipse solarEclipse;
internal LunarEclipse lunarEclipse;
/// <summary>
/// Sunset time.
/// </summary>
public DateTime? SunSet { get { return sunSet; } }
/// <summary>
/// Sunrise time.
/// </summary>
public DateTime? SunRise { get { return sunRise; } }
/// <summary>
/// Moonset time.
/// </summary>
public DateTime? MoonSet { get { return moonSet; } }
/// <summary>
/// Moonrise time.
/// </summary>
public DateTime? MoonRise { get { return moonRise; } }
/// <summary>
/// Sun altitude in degrees (E of N).
/// </summary>
public double SunAltitude { get { return sunAltitude; } }
/// <summary>
/// Sun azimuth in degrees (E of N).
/// </summary>
public double SunAzimuth { get { return sunAzimuth; } }
/// <summary>
/// Moon altitude in degrees (corrected for parallax and refraction).
/// </summary>
public double MoonAltitude { get { return moonAltitude; } }
/// <summary>
/// Moon azimuth in degrees (E of N).
/// </summary>
public double MoonAzimuth { get { return moonAzimuth; } }
/// <summary>
/// Estimated moon distance from the earth.
/// </summary>
public Distance MoonDistance { get { return moonDistance; } }
/// <summary>
/// Sun's Condition based on the provided date.
/// </summary>
public CelestialStatus SunCondition { get { return sunCondition; } }
/// <summary>
/// Moon's condition based on the provided date.
/// </summary>
public CelestialStatus MoonCondition { get { return moonCondition; } }
/// <summary>
/// Determine if the sun is currently up, based on sunset and sunrise time at the provided location and date.
/// </summary>
public bool IsSunUp{ get { return isSunUp; } }
/// <summary>
/// Determine if the moon is currently up, based on moonset and moonrise time at the provided location and date.
/// </summary>
public bool IsMoonUp { get { return isMoonUp; } }
/// <summary>
/// Moon ilumination details based on the provided date.
/// </summary>
/// <remarks>
/// Contains phase, phase name, fraction and angle
/// </remarks>
public MoonIllum MoonIllum { get { return moonIllum; } }
/// <summary>
/// Moons perigee details based on the provided date.
/// </summary>
public Perigee Perigee { get { return perigee; } }
/// <summary>
/// Moons apogee details based on the provided date.
/// </summary>
public Apogee Apogee { get { return apogee; } }
/// <summary>
/// Additional solar event times based on the provided date and location.
/// </summary>
/// <remarks>Contains civil and nautical dawn and dusk times.</remarks>
public AdditionalSolarTimes AdditionalSolarTimes { get { return additionalSolarTimes; } }
/// <summary>
/// Astrological signs based on the provided date.
/// </summary>
/// <remarks>
/// Contains zodiac, moon sign and moon name during full moon events
/// </remarks>
public AstrologicalSigns AstrologicalSigns { get { return astrologicalSigns; } }
/// <summary>
/// Returns a SolarEclipse.
/// </summary>
public SolarEclipse SolarEclipse { get { return solarEclipse; } }
/// <summary>
/// Returns a LunarEclipse.
/// </summary>
public LunarEclipse LunarEclipse { get { return lunarEclipse; } }
/// <summary>
/// Calculates all celestial data. Coordinates will notify as changes occur
/// </summary>
/// <param name="lat">Decimal format latitude</param>
/// <param name="longi">Decimal format longitude</param>
/// <param name="date">Geographic DateTime</param>
internal void CalculateCelestialTime(double lat, double longi, DateTime date)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
SunCalc.CalculateSunTime(lat, longi, date, this);
MoonCalc.GetMoonTimes(date, lat, longi, this);
MoonCalc.GetMoonDistance(date, this);
SunCalc.CalculateZodiacSign(date, this);
MoonCalc.GetMoonSign(date, this);
MoonCalc.GetMoonIllumination(date, this,lat,longi);
perigee = MoonCalc.GetPerigeeEvents(date);
apogee = MoonCalc.GetApogeeEvents(date);
Calculate_Celestial_IsUp_Booleans(date, this);
}
/// <summary>
/// Calculate celestial data based on lat/long and date.
/// </summary>
/// <param name="lat">Decimal format latitude</param>
/// <param name="longi">Decimal format longitude</param>
/// <param name="date">Geographic DateTime</param>
/// <returns>Fully populated Celestial object</returns>
public static Celestial CalculateCelestialTimes(double lat, double longi, DateTime date)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
Celestial c = new Celestial(false);
SunCalc.CalculateSunTime(lat, longi, date, c);
MoonCalc.GetMoonTimes(date, lat, longi, c);
MoonCalc.GetMoonDistance(date, c);
SunCalc.CalculateZodiacSign(date, c);
MoonCalc.GetMoonSign(date, c);
MoonCalc.GetMoonIllumination(date, c,lat,longi);
c.perigee = MoonCalc.GetPerigeeEvents(date);
c.apogee = MoonCalc.GetApogeeEvents(date);
Calculate_Celestial_IsUp_Booleans(date, c);
return c;
}
/// <summary>
/// Calculate sun data based on lat/long and date.
/// </summary>
/// <param name="lat">latitude</param>
/// <param name="longi">longitude</param>
/// <param name="date">DateTime</param>
/// <returns>Celestial (Partially Populated)</returns>
public static Celestial CalculateSunData(double lat, double longi, DateTime date)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
Celestial c = new Celestial(false);
SunCalc.CalculateSunTime(lat, longi, date, c);
SunCalc.CalculateZodiacSign(date, c);
return c;
}
/// <summary>
/// Calculate moon data based on lat/long and date.
/// </summary>
/// <param name="lat">latitude</param>
/// <param name="longi">longitude</param>
/// <param name="date">DateTime</param>
/// <returns>Celestial (Partially Populated)</returns>
public static Celestial CalculateMoonData(double lat, double longi, DateTime date)
{
date = new DateTime(date.Year, date.Month, date.Day, date.Hour, date.Minute, date.Second, DateTimeKind.Utc);
Celestial c = new Celestial(false);
MoonCalc.GetMoonTimes(date, lat, longi, c);
MoonCalc.GetMoonDistance(date, c);
MoonCalc.GetMoonSign(date, c);
MoonCalc.GetMoonIllumination(date, c,lat,longi);
c.perigee = MoonCalc.GetPerigeeEvents(date);
c.apogee = MoonCalc.GetApogeeEvents(date);
return c;
}
/// <summary>
/// Returns a List containing solar eclipse data for the century.
/// Century return is based on the date passed.
/// </summary>
/// <param name="lat">latitude</param>
/// <param name="longi">longitude</param>
/// <param name="date">DateTime</param>
/// <returns></returns>
public static List<SolarEclipseDetails> Get_Solar_Eclipse_Table(double lat, double longi, DateTime date)
{
//Convert to Radians
double latR = lat * Math.PI / 180;
double longR = longi * Math.PI / 180;
//Get solar data based on date
double[] events = Eclipse.SolarData.SolarDateData_100Year(date);
//Return list of solar data.
return SolarEclipseCalc.CalculateSolarEclipse(date, latR, longR, events);
}
/// <summary>
/// Returns a List containing solar eclipse data for the century.
/// Century return is based on the date passed.
/// </summary>
/// <param name="lat">latitude</param>
/// <param name="longi">longitude</param>
/// <param name="date">DateTime</param>
/// <returns></returns>
public static List<LunarEclipseDetails> Get_Lunar_Eclipse_Table(double lat, double longi, DateTime date)
{
//Convert to Radians
double latR = lat * Math.PI / 180;
double longR = longi * Math.PI / 180;
//Get solar data based on date
double[] events = Eclipse.LunarData.LunarDateData_100Year(date);
//Return list of solar data.
return LunarEclipseCalc.CalculateLunarEclipse(date, latR, longR, events);
}
/// <summary>
/// Set bool SunIsUp and MoonIsUp values
/// </summary>
/// <param name="date">Coordinate GeoDate</param>
/// <param name="cel">Celestial Object</param>
private static void Calculate_Celestial_IsUp_Booleans(DateTime date, Celestial cel)
{
//SUN
switch (cel.SunCondition)
{
case CelestialStatus.DownAllDay:
cel.isSunUp = false;
break;
case CelestialStatus.UpAllDay:
cel.isSunUp = true;
break;
case CelestialStatus.NoRise:
if(date<cel.SunSet)
{
cel.isSunUp = true;
}
else { cel.isSunUp = false; }
break;
case CelestialStatus.NoSet:
if (date > cel.SunRise)
{
cel.isSunUp = true;
}
else { cel.isSunUp = false; }
break;
case CelestialStatus.RiseAndSet:
if (cel.SunRise < cel.SunSet)
{
if (date > cel.SunRise && date < cel.SunSet)
{
cel.isSunUp = true;
}
else
{
cel.isSunUp = false;
}
}
else
{
if (date > cel.SunRise || date < cel.SunSet)
{
cel.isSunUp = true;
}
else
{
cel.isSunUp = false;
}
}
break;
default:
//Should never be reached. If reached, previous calculations failed somewhere.
break;
}
//MOON
switch (cel.MoonCondition)
{
case CelestialStatus.DownAllDay:
cel.isMoonUp = false;
break;
case CelestialStatus.UpAllDay:
cel.isMoonUp = true;
break;
case CelestialStatus.NoRise:
if (date < cel.MoonSet)
{
cel.isMoonUp = true;
}
else { cel.isMoonUp = false; }
break;
case CelestialStatus.NoSet:
if (date > cel.MoonRise)
{
cel.isMoonUp = true;
}
else { cel.isMoonUp = false; }
break;
case CelestialStatus.RiseAndSet:
if (cel.MoonRise < cel.MoonSet)
{
if (date > cel.MoonRise && date < cel.MoonSet)
{
cel.isMoonUp = true;
}
else
{
cel.isMoonUp = false;
}
}
else
{
if (date > cel.MoonRise || date < cel.MoonSet)
{
cel.isMoonUp = true;
}
else
{
cel.isMoonUp = false;
}
}
break;
default:
//Should never be reached. If reached, previous calculations failed somewhere.
break;
}
}
/// <summary>
/// Returns Apogee object containing last and next apogee based on the specified date.
/// </summary>
/// <param name="d">DateTime</param>
/// <returns>Apogee</returns>
public static Apogee GetApogees(DateTime d)
{
return MoonCalc.GetApogeeEvents(d);
}
/// <summary>
/// Returns Perigee object containing last and next perigee based on the specified date.
/// </summary>
/// <param name="d">DateTime</param>
/// <returns>Perigee</returns>
public static Perigee GetPerigees(DateTime d)
{
return MoonCalc.GetPerigeeEvents(d);
}
}
}

269
Lora-Map/Helper/Coordinate.Assistant.cs Normal file
View File

@ -0,0 +1,269 @@
using System;
using System.Collections.Generic;
using System.Linq;
namespace CoordinateSharp
{
/// <summary>
/// Used for UTM/MGRS Conversions
/// </summary>
[Serializable]
internal class LatZones
{
public static List<string> longZongLetters = new List<string>(new string[]{"C", "D", "E", "F", "G", "H", "J", "K", "L", "M", "N", "P", "Q", "R", "S", "T",
"U", "V", "W", "X"});
}
/// <summary>
/// Used for handling diagraph determination
/// </summary>
[Serializable]
internal class Digraphs
{
private List<Digraph> digraph1;
private List<Digraph> digraph2;
private String[] digraph1Array = { "A", "B", "C", "D", "E", "F", "G", "H",
"J", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "U", "V", "W", "X",
"Y", "Z" };
private String[] digraph2Array = { "V", "A", "B", "C", "D", "E", "F", "G",
"H", "J", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "U", "V" };
public Digraphs()
{
digraph1 = new List<Digraph>();
digraph2 = new List<Digraph>();
digraph1.Add(new Digraph() { Zone = 1, Letter = "A" });
digraph1.Add(new Digraph() { Zone = 2, Letter = "B" });
digraph1.Add(new Digraph() { Zone = 3, Letter = "C" });
digraph1.Add(new Digraph() { Zone = 4, Letter = "D" });
digraph1.Add(new Digraph() { Zone = 5, Letter = "E" });
digraph1.Add(new Digraph() { Zone = 6, Letter = "F" });
digraph1.Add(new Digraph() { Zone = 7, Letter = "G" });
digraph1.Add(new Digraph() { Zone = 8, Letter = "H" });
digraph1.Add(new Digraph() { Zone = 9, Letter = "J" });
digraph1.Add(new Digraph() { Zone = 10, Letter = "K" });
digraph1.Add(new Digraph() { Zone = 11, Letter = "L" });
digraph1.Add(new Digraph() { Zone = 12, Letter = "M" });
digraph1.Add(new Digraph() { Zone = 13, Letter = "N" });
digraph1.Add(new Digraph() { Zone = 14, Letter = "P" });
digraph1.Add(new Digraph() { Zone = 15, Letter = "Q" });
digraph1.Add(new Digraph() { Zone = 16, Letter = "R" });
digraph1.Add(new Digraph() { Zone = 17, Letter = "S" });
digraph1.Add(new Digraph() { Zone = 18, Letter = "T" });
digraph1.Add(new Digraph() { Zone = 19, Letter = "U" });
digraph1.Add(new Digraph() { Zone = 20, Letter = "V" });
digraph1.Add(new Digraph() { Zone = 21, Letter = "W" });
digraph1.Add(new Digraph() { Zone = 22, Letter = "X" });
digraph1.Add(new Digraph() { Zone = 23, Letter = "Y" });
digraph1.Add(new Digraph() { Zone = 24, Letter = "Z" });
digraph1.Add(new Digraph() { Zone = 1, Letter = "A" });
digraph2.Add(new Digraph() { Zone = 0, Letter = "V"});
digraph2.Add(new Digraph() { Zone = 1, Letter = "A" });
digraph2.Add(new Digraph() { Zone = 2, Letter = "B" });
digraph2.Add(new Digraph() { Zone = 3, Letter = "C" });
digraph2.Add(new Digraph() { Zone = 4, Letter = "D" });
digraph2.Add(new Digraph() { Zone = 5, Letter = "E" });
digraph2.Add(new Digraph() { Zone = 6, Letter = "F" });
digraph2.Add(new Digraph() { Zone = 7, Letter = "G" });
digraph2.Add(new Digraph() { Zone = 8, Letter = "H" });
digraph2.Add(new Digraph() { Zone = 9, Letter = "J" });
digraph2.Add(new Digraph() { Zone = 10, Letter = "K" });
digraph2.Add(new Digraph() { Zone = 11, Letter = "L" });
digraph2.Add(new Digraph() { Zone = 12, Letter = "M" });
digraph2.Add(new Digraph() { Zone = 13, Letter = "N" });
digraph2.Add(new Digraph() { Zone = 14, Letter = "P" });
digraph2.Add(new Digraph() { Zone = 15, Letter = "Q" });
digraph2.Add(new Digraph() { Zone = 16, Letter = "R" });
digraph2.Add(new Digraph() { Zone = 17, Letter = "S" });
digraph2.Add(new Digraph() { Zone = 18, Letter = "T" });
digraph2.Add(new Digraph() { Zone = 19, Letter = "U" });
digraph2.Add(new Digraph() { Zone = 20, Letter = "V" });
}
internal int getDigraph1Index(String letter)
{
for (int i = 0; i < digraph1Array.Length; i++)
{
if (digraph1Array[i].Equals(letter))
{
return i + 1;
}
}
return -1;
}
internal int getDigraph2Index(String letter)
{
for (int i = 0; i < digraph2Array.Length; i++)
{
if (digraph2Array[i].Equals(letter))
{
return i;
}
}
return -1;
}
internal String getDigraph1(int longZone, double easting)
{
int a1 = longZone;
double a2 = 8 * ((a1 - 1) % 3) + 1;
double a3 = easting;
double a4 = a2 + ((int)(a3 / 100000)) - 1;
return digraph1.Where(x=>x.Zone == Math.Floor(a4)).FirstOrDefault().Letter;
}
internal String getDigraph2(int longZone, double northing)
{
int a1 = longZone;
double a2 = 1 + 5 * ((a1 - 1) % 2);
double a3 = northing;
double a4 = (a2 + ((int)(a3 / 100000)));
a4 = (a2 + ((int)(a3 / 100000.0))) % 20;
a4 = Math.Floor(a4);
if (a4 < 0)
{
a4 = a4 + 19;
}
return digraph2.Where(x => x.Zone == Math.Floor(a4)).FirstOrDefault().Letter;
}
}
/// <summary>
/// Diagraph model
/// </summary>
[Serializable]
internal class Digraph
{
public int Zone { get; set; }
public string Letter { get; set; }
}
/// <summary>
/// Used for setting whether a coordinate part is latitudinal or longitudinal.
/// </summary>
[Serializable]
public enum CoordinateType
{
/// <summary>
/// Latitude
/// </summary>
Lat,
/// <summary>
/// Longitude
/// </summary>
Long
}
/// <summary>
/// Used to set a coordinate part position.
/// </summary>
[Serializable]
public enum CoordinatesPosition :int
{
/// <summary>
/// North
/// </summary>
N,
/// <summary>
/// East
/// </summary>
E,
/// <summary>
/// South
/// </summary>
S,
/// <summary>
/// West
/// </summary>
W
}
/// <summary>
/// Coordinate type datum specification
/// </summary>
[Serializable]
[Flags]
public enum Coordinate_Datum
{
/// <summary>
/// Lat Long GeoDetic
/// </summary>
LAT_LONG = 1,
/// <summary>
/// UTM and MGRS
/// </summary>
UTM_MGRS = 2,
/// <summary>
/// ECEF
/// </summary>
ECEF = 4,
}
/// <summary>
/// Cartesian Coordinate Type
/// </summary>
public enum CartesianType
{
/// <summary>
/// Spherical Cartesian
/// </summary>
Cartesian,
/// <summary>
/// Earth Centered Earth Fixed
/// </summary>
ECEF,
}
/// <summary>
/// Used for easy read math functions
/// </summary>
[Serializable]
internal static class ModM
{
public static double Mod(double x, double y)
{
return x - y * Math.Floor(x / y);
}
public static double ModLon(double x)
{
return Mod(x + Math.PI, 2 * Math.PI) - Math.PI;
}
public static double ModCrs(double x)
{
return Mod(x, 2 * Math.PI);
}
public static double ModLat(double x)
{
return Mod(x + Math.PI / 2, 2 * Math.PI) - Math.PI / 2;
}
}
/// <summary>
/// Earth Shape for Calculations.
/// </summary>
[Serializable]
public enum Shape
{
/// <summary>
/// Calculate as sphere (less accurate, more efficient).
/// </summary>
Sphere,
/// <summary>
/// Calculate as ellipsoid (more accurate, less efficient).
/// </summary>
Ellipsoid
}
}

155
Lora-Map/Helper/Coordinate.Cartesian.cs Normal file
View File

@ -0,0 +1,155 @@
using System;
using System.ComponentModel;
namespace CoordinateSharp
{
/// <summary>
/// Cartesian (X, Y, Z) Coordinate
/// </summary>
[Serializable]
public class Cartesian : INotifyPropertyChanged
{
/// <summary>
/// Create a Cartesian Object
/// </summary>
/// <param name="c"></param>
public Cartesian(Coordinate c)
{
//formulas:
x = Math.Cos(c.Latitude.ToRadians()) * Math.Cos(c.Longitude.ToRadians());
y = Math.Cos(c.Latitude.ToRadians()) * Math.Sin(c.Longitude.ToRadians());
z = Math.Sin(c.Latitude.ToRadians());
}
/// <summary>
/// Create a Cartesian Object
/// </summary>
/// <param name="xc">X</param>
/// <param name="yc">Y</param>
/// <param name="zc">Z</param>
public Cartesian(double xc, double yc, double zc)
{
//formulas:
x = xc;
y = yc;
z = zc;
}
/// <summary>
/// Updates Cartesian Values
/// </summary>
/// <param name="c"></param>
public void ToCartesian(Coordinate c)
{
x = Math.Cos(c.Latitude.ToRadians()) * Math.Cos(c.Longitude.ToRadians());
y = Math.Cos(c.Latitude.ToRadians()) * Math.Sin(c.Longitude.ToRadians());
z = Math.Sin(c.Latitude.ToRadians());
}
private double x;
private double y;
private double z;
/// <summary>
/// X Coordinate
/// </summary>
public double X
{
get { return x; }
set
{
if(x != value)
{
x = value;
NotifyPropertyChanged("X");
}
}
}
/// <summary>
/// y Coordinate
/// </summary>
public double Y
{
get { return y; }
set
{
if (y != value)
{
y = value;
NotifyPropertyChanged("Y");
}
}
}
/// <summary>
/// Z Coordinate
/// </summary>
public double Z
{
get { return z; }
set
{
if (z != value)
{
z = value;
NotifyPropertyChanged("Z");
}
}
}
/// <summary>
/// Returns a Lat Long Coordinate object based on the provided Cartesian Coordinate
/// </summary>
/// <param name="x">X</param>
/// <param name="y">Y</param>
/// <param name="z">Z</param>
/// <returns></returns>
public static Coordinate CartesianToLatLong(double x, double y, double z)
{
double lon = Math.Atan2(y, x);
double hyp = Math.Sqrt(x * x + y * y);
double lat = Math.Atan2(z, hyp);
double Lat = lat * (180 / Math.PI);
double Lon = lon * (180 / Math.PI);
return new Coordinate(Lat, Lon);
}
/// <summary>
/// Returns a Lat Long Coordinate object based on the provided Cartesian Coordinate
/// </summary>
/// <param name="cart">Cartesian Coordinate</param>
/// <returns></returns>
public static Coordinate CartesianToLatLong(Cartesian cart)
{
double x = cart.X;
double y = cart.Y;
double z = cart.Z;
double lon = Math.Atan2(y, x);
double hyp = Math.Sqrt(x * x + y * y);
double lat = Math.Atan2(z, hyp);
double Lat = lat * (180 / Math.PI);
double Lon = lon * (180 / Math.PI);
return new Coordinate(Lat, Lon);
}
/// <summary>
/// Cartesian Default String Format
/// </summary>
/// <returns>Cartesian Formatted Coordinate String</returns>
/// <returns>Values rounded to the 8th place</returns>
public override string ToString()
{
return Math.Round(x,8).ToString() + " " + Math.Round(y, 8).ToString() + " " + Math.Round(z, 8).ToString();
}
/// <summary>
/// Property changed event
/// </summary>
public event PropertyChangedEventHandler PropertyChanged;
/// <summary>
/// Notify property changed
/// </summary>
/// <param name="propName">Property name</param>
public void NotifyPropertyChanged(string propName)
{
if (this.PropertyChanged != null)
{
PropertyChanged(this, new PropertyChangedEventArgs(propName));
}
}
}
}

566
Lora-Map/Helper/Coordinate.ECEF.cs Normal file
View File

@ -0,0 +1,566 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.ComponentModel;
namespace CoordinateSharp
{
/// <summary>
/// Earth Centered - Earth Fixed (X,Y,Z) Coordinate
/// </summary>
[Serializable]
public class ECEF : INotifyPropertyChanged
{
/// <summary>
/// Create an ECEF Object
/// </summary>
/// <param name="c">Coordinate</param>
public ECEF(Coordinate c)
{
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
WGS84();
geodetic_height = new Distance(0);
double[] ecef = LatLong_To_ECEF(c.Latitude.DecimalDegree, c.Longitude.DecimalDegree, geodetic_height.Kilometers);
x = ecef[0];
y = ecef[1];
z = ecef[2];
}
/// <summary>
/// Create an ECEF Object
/// </summary>
/// <param name="c">Coordinate</param>
/// <param name="height">Coordinate</param>
public ECEF(Coordinate c, Distance height)
{
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
WGS84();
geodetic_height = height;
double[] ecef = LatLong_To_ECEF(c.Latitude.DecimalDegree, c.Longitude.DecimalDegree, geodetic_height.Kilometers);
x = ecef[0];
y = ecef[1];
z = ecef[2];
}
/// <summary>
/// Create an ECEF Object
/// </summary>
/// <param name="xc">X</param>
/// <param name="yc">Y</param>
/// <param name="zc">Z</param>
public ECEF(double xc, double yc, double zc)
{
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
WGS84();
geodetic_height = new Distance(0);
x = xc;
y = yc;
z = zc;
}
/// <summary>
/// Updates ECEF Values
/// </summary>
/// <param name="c">Coordinate</param>
public void ToECEF(Coordinate c)
{
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
WGS84();
double[] ecef = LatLong_To_ECEF(c.Latitude.DecimalDegree, c.Longitude.DecimalDegree, geodetic_height.Kilometers);
x = ecef[0];
y = ecef[1];
z = ecef[2];
}
//Globals for calucations
private double EARTH_A;
private double EARTH_B;
private double EARTH_F;
private double EARTH_Ecc;
private double EARTH_Esq;
//ECEF Values
private double x;
private double y;
private double z;
private Distance geodetic_height;
//Datum
internal double equatorial_radius;
internal double inverse_flattening;
/// <summary>
/// Datum Equatorial Radius / Semi Major Axis
/// </summary>
public double Equatorial_Radius
{
get { return equatorial_radius; }
}
/// <summary>
/// Datum Flattening
/// </summary>
public double Inverse_Flattening
{
get { return inverse_flattening; }
}
/// <summary>
/// X Coordinate
/// </summary>
public double X
{
get { return x; }
set
{
if (x != value)
{
x = value;
NotifyPropertyChanged("X");
}
}
}
/// <summary>
/// y Coordinate
/// </summary>
public double Y
{
get { return y; }
set
{
if (y != value)
{
y = value;
NotifyPropertyChanged("Y");
}
}
}
/// <summary>
/// Z Coordinate
/// </summary>
public double Z
{
get { return z; }
set
{
if (z != value)
{
z = value;
NotifyPropertyChanged("Z");
}
}
}
/// <summary>
/// GeoDetic Height from Mean Sea Level.
/// Used for converting Lat Long / ECEF.
/// Default value is 0. Adjust as needed.
/// </summary>
public Distance GeoDetic_Height
{
get { return geodetic_height; }
internal set
{
if (geodetic_height != value)
{
geodetic_height = value;
NotifyPropertyChanged("Height");
}
}
}
/// <summary>
/// Sets GeoDetic height for ECEF conversion.
/// Recalculate ECEF Coordinate
/// </summary>
/// <param name="c">Coordinate</param>
/// <param name="dist">Height</param>
public void Set_GeoDetic_Height(Coordinate c, Distance dist)
{
geodetic_height = dist;
double[] values = LatLong_To_ECEF(c.Latitude.DecimalDegree, c.Longitude.DecimalDegree, dist.Kilometers);
x = values[0];
y = values[1];
z = values[2];
}
/// <summary>
/// Returns a Geodetic Coordinate object based on the provided ECEF Coordinate
/// </summary>
/// <param name="x">X</param>
/// <param name="y">Y</param>
/// <param name="z">Z</param>
/// <returns>Coordinate</returns>
public static Coordinate ECEFToLatLong(double x, double y, double z)
{
ECEF ecef = new ECEF(x, y, z);
double[] values = ecef.ECEF_To_LatLong(x, y, z);
ecef.geodetic_height =new Distance(values[2]);
Coordinate c = new Coordinate(values[0], values[1]);
c.ECEF = ecef;
return c;
}
/// <summary>
/// Returns a Geodetic Coordinate object based on the provided ECEF Coordinate
/// </summary>
/// <param name="ecef">ECEF Coordinate</param>
/// <returns>Coordinate</returns>
public static Coordinate ECEFToLatLong(ECEF ecef)
{
double[] values = ecef.ECEF_To_LatLong(ecef.X, ecef.Y, ecef.Z);
Coordinate c = new Coordinate(values[0], values[1]);
Distance height = new Distance(values[2]);
ecef.geodetic_height = new Distance(values[2]);
c.ECEF = ecef;
return c;
}
/// <summary>
/// ECEF Default String Format
/// </summary>
/// <returns>ECEF Formatted Coordinate String</returns>
/// <returns>Values rounded to the 3rd place</returns>
public override string ToString()
{
return Math.Round(x, 3).ToString() + " km, " + Math.Round(y, 3).ToString() + " km, " + Math.Round(z, 3).ToString() + " km";
}
/// <summary>
/// Property changed event
/// </summary>
public event PropertyChangedEventHandler PropertyChanged;
/// <summary>
/// Notify property changed
/// </summary>
/// <param name="propName">Property name</param>
public void NotifyPropertyChanged(string propName)
{
if (PropertyChanged != null)
{
PropertyChanged(this, new PropertyChangedEventArgs(propName));
}
}
//CONVERSION LOGIC
/// <summary>
/// Initialize EARTH global variables based on the Datum
/// </summary>
private void WGS84()
{
double wgs84a = equatorial_radius / 1000;
double wgs84f = 1.0 / inverse_flattening;
double wgs84b = wgs84a * (1.0 - wgs84f);
EarthCon(wgs84a, wgs84b);
}
/// <summary>
/// Sets Earth Constants as Globals
/// </summary>
/// <param name="a">a</param>
/// <param name="b">b</param>
private void EarthCon(double a, double b)
{
double f = 1 - b / a;
double eccsq = 1 - b * b / (a * a);
double ecc = Math.Sqrt(eccsq);
EARTH_A = a;
EARTH_B = b;
EARTH_F = f;
EARTH_Ecc = ecc;
EARTH_Esq = eccsq;
}
/// <summary>
/// Compute the radii at the geodetic latitude (degrees)
/// </summary>
/// <param name="lat">Latitude in degres</param>
/// <returns>double[]</returns>
private double[] radcur(double lat)
{
double[] rrnrm = new double[3];
double dtr = Math.PI / 180.0;
double a = EARTH_A;
double b = EARTH_B;
double asq = a * a;
double bsq = b * b;
double eccsq = 1 - bsq / asq;
double ecc = Math.Sqrt(eccsq);
double clat = Math.Cos(dtr * lat);
double slat = Math.Sin(dtr * lat);
double dsq = 1.0 - eccsq * slat * slat;
double d = Math.Sqrt(dsq);
double rn = a / d;
double rm = rn * (1.0 - eccsq) / dsq;
double rho = rn * clat;
double z = (1.0 - eccsq) * rn * slat;
double rsq = rho * rho + z * z;
double r = Math.Sqrt(rsq);
rrnrm[0] = r;
rrnrm[1] = rn;
rrnrm[2] = rm;
return (rrnrm);
}
/// <summary>
/// Physical radius of the Earth
/// </summary>
/// <param name="lat">Latidude in degrees</param>
/// <returns>double</returns>
private double rearth(double lat)
{
double[] rrnrm;
rrnrm = radcur(lat);
double r = rrnrm[0];
return r;
}
/// <summary>
/// Converts geocentric latitude to geodetic latitude
/// </summary>
/// <param name="flatgc">Geocentric latitude</param>
/// <param name="altkm">Altitude in KM</param>
/// <returns>double</returns>
private double gc2gd(double flatgc, double altkm)
{
var dtr = Math.PI / 180.0;
var rtd = 1 / dtr;
double ecc = EARTH_Ecc;
double esq = ecc * ecc;
//approximation by stages
//1st use gc-lat as if is gd, then correct alt dependence
double altnow = altkm;
double[] rrnrm = radcur(flatgc);
double rn = rrnrm[1];
double ratio = 1 - esq * rn / (rn + altnow);
double tlat = Math.Tan(dtr * flatgc) / ratio;
double flatgd = rtd * Math.Atan(tlat);
//now use this approximation for gd-lat to get rn etc.
rrnrm = radcur(flatgd);
rn = rrnrm[1];
ratio = 1 - esq * rn / (rn + altnow);
tlat = Math.Tan(dtr * flatgc) / ratio;
flatgd = rtd * Math.Atan(tlat);
return flatgd;
}
/// <summary>
/// Converts geodetic latitude to geocentric latitude
/// </summary>
/// <param name="flatgd">Geodetic latitude tp geocentric latitide</param>
/// <param name="altkm">Altitude in KM</param>
/// <returns>double</returns>
private double gd2gc(double flatgd, double altkm)
{
double dtr = Math.PI / 180.0;
double rtd = 1 / dtr;
double ecc = EARTH_Ecc;
double esq = ecc * ecc;
double altnow = altkm;
double[] rrnrm = radcur(flatgd);
double rn = rrnrm[1];
double ratio = 1 - esq * rn / (rn + altnow);
double tlat = Math.Tan(dtr * flatgd) * ratio;
double flatgc = rtd * Math.Atan(tlat);
return flatgc;
}
/// <summary>
/// Converts lat / long to east, north, up vectors
/// </summary>
/// <param name="flat">Latitude</param>
/// <param name="flon">Longitude</param>
/// <returns>Array[] of double[]</returns>
private Array[] llenu(double flat, double flon)
{
double clat, slat, clon, slon;
double[] ee = new double[3];
double[] en = new double[3];
double[] eu = new double[3];
Array[] enu = new Array[3];
double dtr = Math.PI / 180.0;
clat = Math.Cos(dtr * flat);
slat = Math.Sin(dtr * flat);
clon = Math.Cos(dtr * flon);
slon = Math.Sin(dtr * flon);
ee[0] = -slon;
ee[1] = clon;
ee[2] = 0.0;
en[0] = -clon * slat;
en[1] = -slon * slat;
en[2] = clat;
eu[0] = clon * clat;
eu[1] = slon * clat;
eu[2] = slat;
enu[0] = ee;
enu[1] = en;
enu[2] = eu;
return enu;
}
/// <summary>
/// Gets ECEF vector in KM
/// </summary>
/// <param name="lat">Latitude</param>
/// <param name="longi">Longitude</param>
/// <param name="altkm">Altitude in KM</param>
/// <returns>double[]</returns>
private double[] LatLong_To_ECEF(double lat, double longi, double altkm)
{
double dtr = Math.PI / 180.0;
double clat = Math.Cos(dtr * lat);
double slat = Math.Sin(dtr * lat);
double clon = Math.Cos(dtr * longi);
double slon = Math.Sin(dtr * longi);
double[] rrnrm = radcur(lat);
double rn = rrnrm[1];
double re = rrnrm[0];
double ecc = EARTH_Ecc;
double esq = ecc * ecc;
double x = (rn + altkm) * clat * clon;
double y = (rn + altkm) * clat * slon;
double z = ((1 - esq) * rn + altkm) * slat;
double[] xvec = new double[3];
xvec[0] = x;
xvec[1] = y;
xvec[2] = z;
return xvec;
}
/// <summary>
/// Converts ECEF X, Y, Z to GeoDetic Lat / Long and Height in KM
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <param name="z"></param>
/// <returns></returns>
private double[] ECEF_To_LatLong(double x, double y, double z)
{
var dtr = Math.PI / 180.0;
double[] rrnrm = new double[3];
double[] llhvec = new double[3];
double slat, tangd, flatn, dlat, clat;
double flat;
double altkm;
double esq = EARTH_Esq;
double rp = Math.Sqrt(x * x + y * y + z * z);
double flatgc = Math.Asin(z / rp) / dtr;
double flon;
double testval = Math.Abs(x) + Math.Abs(y);
if (testval < 1.0e-10)
{ flon = 0.0; }
else
{ flon = Math.Atan2(y, x) / dtr; }
if (flon < 0.0) { flon = flon + 360.0; }
double p = Math.Sqrt(x * x + y * y);
//Pole special case
if (p < 1.0e-10)
{
flat = 90.0;
if (z < 0.0) { flat = -90.0; }
altkm = rp - rearth(flat);
llhvec[0] = flat;
llhvec[1] = flon;
llhvec[2] = altkm;
return llhvec;
}
//first iteration, use flatgc to get altitude
//and alt needed to convert gc to gd lat.
double rnow = rearth(flatgc);
altkm = rp - rnow;
flat = gc2gd(flatgc, altkm);
rrnrm = radcur(flat);
double rn = rrnrm[1];
for (int kount = 0; kount < 5; kount++)
{
slat = Math.Sin(dtr * flat);
tangd = (z + rn * esq * slat) / p;
flatn = Math.Atan(tangd) / dtr;
dlat = flatn - flat;
flat = flatn;
clat = Math.Cos(dtr * flat);
rrnrm = radcur(flat);
rn = rrnrm[1];
altkm = (p / clat) - rn;
if (Math.Abs(dlat) < 1.0e-12) { break; }
}
//CONVERTER WORKS IN E LAT ONLY, IF E LAT > 180 LAT IS WEST SO IT MUCST BE CONVERTED TO Decimal
if (flon > 180) { flon = flon - 360; }
llhvec[0] = flat;
llhvec[1] = flon;
llhvec[2] = altkm;
return llhvec;
}
}
}

119
Lora-Map/Helper/Coordinate.EagerLoad.cs Normal file
View File

@ -0,0 +1,119 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace CoordinateSharp
{
/// <summary>
/// Turn on/off eager loading of certain properties.
/// </summary>
[Serializable]
public class EagerLoad
{
/// <summary>
/// Create an EagerLoad object
/// </summary>
public EagerLoad()
{
Celestial = true;
UTM_MGRS = true;
Cartesian = true;
ECEF = true;
}
/// <summary>
/// Create an EagerLoad object with all options on or off
/// </summary>
/// <param name="isOn">Turns EagerLoad on or off</param>
public EagerLoad(bool isOn)
{
Celestial = isOn;
UTM_MGRS = isOn;
Cartesian = isOn;
ECEF = isOn;
}
/// <summary>
/// Create an EagerLoad object with only the specified flag options turned on.
/// </summary>
/// <param name="et">EagerLoadType</param>
public EagerLoad(EagerLoadType et)
{
Cartesian = false;
Celestial = false;
UTM_MGRS = false;
ECEF = false;
if (et.HasFlag(EagerLoadType.Cartesian))
{
Cartesian = true;
}
if (et.HasFlag(EagerLoadType.Celestial))
{
Celestial = true;
}
if (et.HasFlag(EagerLoadType.UTM_MGRS))
{
UTM_MGRS = true;
}
if (et.HasFlag(EagerLoadType.ECEF))
{
ECEF = true;
}
}
/// <summary>
/// Creates an EagerLoad object. Only the specified flags will be set to EagerLoad.
/// </summary>
/// <param name="et">EagerLoadType</param>
/// <returns>EagerLoad</returns>
public static EagerLoad Create(EagerLoadType et)
{
EagerLoad el = new EagerLoad(et);
return el;
}
/// <summary>
/// Eager load celestial information.
/// </summary>
public bool Celestial { get; set; }
/// <summary>
/// Eager load UTM and MGRS information
/// </summary>
public bool UTM_MGRS { get; set; }
/// <summary>
/// Eager load Cartesian information
/// </summary>
public bool Cartesian { get; set; }
/// <summary>
/// Eager load ECEF information
/// </summary>
public bool ECEF { get; set; }
}
/// <summary>
/// EagerLoad Enumerator
/// </summary>
[Serializable]
[Flags]
public enum EagerLoadType
{
/// <summary>
/// UTM and MGRS
/// </summary>
UTM_MGRS = 1,
/// <summary>
/// Celestial
/// </summary>
Celestial = 2,
/// <summary>
/// Cartesian
/// </summary>
Cartesian = 4,
/// <summary>
/// ECEF
/// </summary>
ECEF = 8
}
}

107
Lora-Map/Helper/Coordinate.Formatting.cs Normal file
View File

@ -0,0 +1,107 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace CoordinateSharp
{
/// <summary>
/// Coordinate formatting options for a Coordinate object.
/// </summary>
[Serializable]
public class CoordinateFormatOptions
{
/// <summary>
/// Set default values with the constructor.
/// </summary>
public CoordinateFormatOptions()
{
Format = CoordinateFormatType.Degree_Minutes_Seconds;
Round = 3;
Display_Leading_Zeros = false;
Display_Trailing_Zeros = false;
Display_Symbols = true;
Display_Degree_Symbol = true;
Display_Minute_Symbol = true;
Display_Seconds_Symbol = true;
Display_Hyphens = false;
Position_First = true;
}
/// <summary>
/// Coordinate format type.
/// </summary>
public CoordinateFormatType Format { get; set; }
/// <summary>
/// Rounds Coordinates to the set value.
/// </summary>
public int Round { get; set; }
/// <summary>
/// Displays leading zeros.
/// </summary>
public bool Display_Leading_Zeros { get; set; }
/// <summary>
/// Display trailing zeros.
/// </summary>
public bool Display_Trailing_Zeros { get; set; }
/// <summary>
/// Allow symbols to display.
/// </summary>
public bool Display_Symbols { get; set; }
/// <summary>
/// Display degree symbols.
/// </summary>
public bool Display_Degree_Symbol { get; set; }
/// <summary>
/// Display minute symbols.
/// </summary>
public bool Display_Minute_Symbol { get; set; }
/// <summary>
/// Display secons symbol.
/// </summary>
public bool Display_Seconds_Symbol { get; set; }
/// <summary>
/// Display hyphens between values.
/// </summary>
public bool Display_Hyphens { get; set; }
/// <summary>
/// Show coordinate position first.
/// Will show last if set 'false'.
/// </summary>
public bool Position_First { get; set; }
}
/// <summary>
/// Coordinate Format Types.
/// </summary>
[Serializable]
public enum CoordinateFormatType
{
/// <summary>
/// Decimal Degree Format
/// </summary>
/// <remarks>
/// Example: N 40.456 W 75.456
/// </remarks>
Decimal_Degree,
/// <summary>
/// Decimal Degree Minutes Format
/// </summary>
/// <remarks>
/// Example: N 40º 34.552' W 70º 45.408'
/// </remarks>
Degree_Decimal_Minutes,
/// <summary>
/// Decimal Degree Minutes Format
/// </summary>
/// <remarks>
/// Example: N 40º 34" 36.552' W 70º 45" 24.408'
/// </remarks>
Degree_Minutes_Seconds,
/// <summary>
/// Decimal Format
/// </summary>
/// <remarks>
/// Example: 40.57674 -70.46574
/// </remarks>
Decimal
}
}

236
Lora-Map/Helper/Coordinate.MGRS.cs Normal file
View File

@ -0,0 +1,236 @@
using System;
using System.Linq;
using System.Diagnostics;
using System.ComponentModel;
namespace CoordinateSharp {
/// <summary>
/// Military Grid Reference System (MGRS). Uses the WGS 84 Datum.
/// Relies upon values from the UniversalTransverseMercator class
/// </summary>
[Serializable]
public class MilitaryGridReferenceSystem : INotifyPropertyChanged {
/// <summary>
/// Create an MGRS object with WGS84 datum
/// </summary>
/// <param name="latz">Lat Zone</param>
/// <param name="longz">Long Zone</param>
/// <param name="d">Digraph</param>
/// <param name="e">Easting</param>
/// <param name="n">Northing</param>
public MilitaryGridReferenceSystem(String latz, Int32 longz, String d, Double e, Double n) {
String digraphLettersE = "ABCDEFGHJKLMNPQRSTUVWXYZ";
String digraphLettersN = "ABCDEFGHJKLMNPQRSTUV";
if(longz < 1 || longz > 60) { Debug.WriteLine("Longitudinal zone out of range", "UTM longitudinal zones must be between 1-60."); }
if(!Verify_Lat_Zone(latz)) { throw new ArgumentException("Latitudinal zone invalid", "UTM latitudinal zone was unrecognized."); }
if(n < 0 || n > 10000000) { throw new ArgumentOutOfRangeException("Northing out of range", "Northing must be between 0-10,000,000."); }
if(d.Count() < 2 || d.Count() > 2) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
if(digraphLettersE.ToCharArray().ToList().Where(x => x.ToString() == d.ToUpper()[0].ToString()).Count() == 0) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
if(digraphLettersN.ToCharArray().ToList().Where(x => x.ToString() == d.ToUpper()[1].ToString()).Count() == 0) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
this.LatZone = latz;
this.LongZone = longz;
this.Digraph = d;
this.Easting = e;
this.Northing = n;
//WGS84
this.equatorialRadius = 6378137.0;
this.inverseFlattening = 298.257223563;
}
/// <summary>
/// Create an MGRS object with custom datum
/// </summary>
/// <param name="latz">Lat Zone</param>
/// <param name="longz">Long Zone</param>
/// <param name="d">Digraph</param>
/// <param name="e">Easting</param>
/// <param name="n">Northing</param>
/// <param name="rad">Equatorial Radius</param>
/// <param name="flt">Inverse Flattening</param>
public MilitaryGridReferenceSystem(String latz, Int32 longz, String d, Double e, Double n, Double rad, Double flt) {
String digraphLettersE = "ABCDEFGHJKLMNPQRSTUVWXYZ";
String digraphLettersN = "ABCDEFGHJKLMNPQRSTUV";
if(longz < 1 || longz > 60) { Debug.WriteLine("Longitudinal zone out of range", "UTM longitudinal zones must be between 1-60."); }
if(!Verify_Lat_Zone(latz)) { throw new ArgumentException("Latitudinal zone invalid", "UTM latitudinal zone was unrecognized."); }
if(n < 0 || n > 10000000) { throw new ArgumentOutOfRangeException("Northing out of range", "Northing must be between 0-10,000,000."); }
if(d.Count() < 2 || d.Count() > 2) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
if(digraphLettersE.ToCharArray().ToList().Where(x => x.ToString() == d.ToUpper()[0].ToString()).Count() == 0) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
if(digraphLettersN.ToCharArray().ToList().Where(x => x.ToString() == d.ToUpper()[1].ToString()).Count() == 0) { throw new ArgumentException("Digraph invalid", "MGRS Digraph was unrecognized."); }
this.LatZone = latz;
this.LongZone = longz;
this.Digraph = d;
this.Easting = e;
this.Northing = n;
this.equatorialRadius = rad;
this.inverseFlattening = flt;
}
private Double equatorialRadius;
private Double inverseFlattening;
private Boolean Verify_Lat_Zone(String l) {
if(LatZones.longZongLetters.Where(x => x == l.ToUpper()).Count() != 1) {
return false;
}
return true;
}
/// <summary>
/// MGRS Zone Letter
/// </summary>
public String LatZone { get; private set; }
/// <summary>
/// MGRS Zone Number
/// </summary>
public Int32 LongZone { get; private set; }
/// <summary>
/// MGRS Easting
/// </summary>
public Double Easting { get; private set; }
/// <summary>
/// MGRS Northing
/// </summary>
public Double Northing { get; private set; }
/// <summary>
/// MGRS Digraph
/// </summary>
public String Digraph { get; private set; }
/// <summary>
/// Is MGRS conversion within the coordinate system's accurate boundaries after conversion from Lat/Long.
/// </summary>
public Boolean WithinCoordinateSystemBounds { get; private set; } = true;
internal MilitaryGridReferenceSystem(UniversalTransverseMercator utm) => ToMGRS(utm);
internal void ToMGRS(UniversalTransverseMercator utm) {
Digraphs digraphs = new Digraphs();
String digraph1 = digraphs.getDigraph1(utm.LongZone, utm.Easting);
String digraph2 = digraphs.getDigraph2(utm.LongZone, utm.Northing);
this.Digraph = digraph1 + digraph2;
this.LatZone = utm.LatZone;
this.LongZone = utm.LongZone;
//String easting = String.valueOf((int)_easting);
String e = ((Int32)utm.Easting).ToString();
if(e.Length < 5) {
e = "00000" + ((Int32)utm.Easting).ToString();
}
e = e.Substring(e.Length - 5);
this.Easting = Convert.ToInt32(e);
String n = ((Int32)utm.Northing).ToString();
if(n.Length < 5) {
n = "0000" + ((Int32)utm.Northing).ToString();
}
n = n.Substring(n.Length - 5);
this.Northing = Convert.ToInt32(n);
this.equatorialRadius = utm.equatorial_radius;
this.inverseFlattening = utm.inverse_flattening;
this.WithinCoordinateSystemBounds = utm.WithinCoordinateSystemBounds;
}
/// <summary>
/// Creates a Coordinate object from an MGRS/NATO UTM Coordinate
/// </summary>
/// <param name="mgrs">MilitaryGridReferenceSystem</param>
/// <returns>Coordinate object</returns>
public static Coordinate MGRStoLatLong(MilitaryGridReferenceSystem mgrs) {
String latz = mgrs.LatZone;
String digraph = mgrs.Digraph;
Char eltr = digraph[0];
Char nltr = digraph[1];
String digraphLettersE = "ABCDEFGHJKLMNPQRSTUVWXYZ";
String digraphLettersN = "ABCDEFGHJKLMNPQRSTUV";
String digraphLettersAll = "";
for(Int32 lt = 1; lt < 25; lt++) {
digraphLettersAll += "ABCDEFGHJKLMNPQRSTUV";
}
Int32 eidx = digraphLettersE.IndexOf(eltr);
Int32 nidx = digraphLettersN.IndexOf(nltr);
if(mgrs.LongZone / 2.0 == Math.Floor(mgrs.LongZone / 2.0)) {
nidx -= 5; // correction for even numbered zones
}
Double ebase = 100000 * (1 + eidx - 8 * Math.Floor(Convert.ToDouble(eidx) / 8));
Int32 latBand = digraphLettersE.IndexOf(latz);
Int32 latBandLow = 8 * latBand - 96;
Int32 latBandHigh = 8 * latBand - 88;
if(latBand < 2) {
latBandLow = -90;
latBandHigh = -80;
} else if(latBand == 21) {
latBandLow = 72;
latBandHigh = 84;
} else if(latBand > 21) {
latBandLow = 84;
latBandHigh = 90;
}
Double lowLetter = Math.Floor(100 + 1.11 * latBandLow);
Double highLetter = Math.Round(100 + 1.11 * latBandHigh);
String latBandLetters = null;
Int32 l = Convert.ToInt32(lowLetter);
Int32 h = Convert.ToInt32(highLetter);
if(mgrs.LongZone / 2.0 == Math.Floor(mgrs.LongZone / 2.0)) {
latBandLetters = digraphLettersAll.Substring(l + 5, h + 5).ToString();
} else {
latBandLetters = digraphLettersAll.Substring(l, h).ToString();
}
Double nbase = 100000 * (lowLetter + latBandLetters.IndexOf(nltr));
//latBandLetters.IndexOf(nltr) value causing incorrect Northing below -80
Double x = ebase + mgrs.Easting;
Double y = nbase + mgrs.Northing;
if(y > 10000000) {
y = y - 10000000;
}
if(nbase >= 10000000) {
y = nbase + mgrs.Northing - 10000000;
}
Boolean southern = nbase < 10000000;
UniversalTransverseMercator utm = new UniversalTransverseMercator(mgrs.LatZone, mgrs.LongZone, x, y) {
equatorial_radius = mgrs.equatorialRadius,
inverse_flattening = mgrs.inverseFlattening
};
Coordinate c = UniversalTransverseMercator.ConvertUTMtoLatLong(utm);
c.Set_Datum(mgrs.equatorialRadius, mgrs.inverseFlattening);
return c;
}
/// <summary>
/// MGRS Default String Format
/// </summary>
/// <returns>MGRS Formatted Coordinate String</returns>
public override String ToString() {
if(!this.WithinCoordinateSystemBounds) { return ""; }//MGRS Coordinate is outside its reliable boundaries. Return empty.
return this.LongZone.ToString() + this.LatZone + " " + this.Digraph + " " + ((Int32)this.Easting).ToString("00000") + " " + ((Int32)this.Northing).ToString("00000");
}
/// <summary>
/// Property changed event
/// </summary>
public event PropertyChangedEventHandler PropertyChanged;
/// <summary>
/// Notify property changed
/// </summary>
/// <param name="propName">Property name</param>
public void NotifyPropertyChanged(String propName) => PropertyChanged?.Invoke(this, new PropertyChangedEventArgs(propName));
}
}

925
Lora-Map/Helper/Coordinate.Parser.cs Normal file
View File

@ -0,0 +1,925 @@
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Text.RegularExpressions;
namespace CoordinateSharp
{
/// <summary>
/// Type of format a Coordinate parsed from.
/// </summary>
[Serializable]
public enum Parse_Format_Type
{
/// <summary>
/// Coordinate was not initialized from a parser method.
/// </summary>
None,
/// <summary>
/// Signed Degree
/// DD.dddd
/// </summary>
Signed_Degree,
/// <summary>
/// Decimal Degree
/// P DD.dddd
/// </summary>
Decimal_Degree,
/// <summary>
/// Degree Decimal Minute
/// P DD MM.sss
/// </summary>
Degree_Decimal_Minute,
/// <summary>
/// Degree Minute Second
/// P DD MM SS.sss
/// </summary>
Degree_Minute_Second,
/// <summary>
/// Universal Transverse Mercator
/// </summary>
UTM,
/// <summary>
/// Military Grid Reference System
/// </summary>
MGRS,
/// <summary>
/// Spherical Cartesian
/// </summary>
Cartesian_Spherical,
/// <summary>
/// Earth Centered Earth Fixed
/// </summary>
Cartesian_ECEF
}
internal class FormatFinder
{
//Add main to Coordinate and tunnel to Format class. Add private methods to format.
//WHEN PARSING NO EXCPETIONS FOR OUT OF RANGE ARGS WILL BE THROWN
public static bool TryParse(string coordString, CartesianType ct, out Coordinate c)
{
//Turn of eagerload for efficiency
EagerLoad eg = new EagerLoad();
eg.Cartesian = false;
eg.Celestial = false;
eg.UTM_MGRS = false;
c = new Coordinate(eg);
string s = coordString;
s = s.Trim(); //Trim all spaces before and after string
double[] d;
//Try Signed Degree
if (TrySignedDegree(s, out d))
{
try
{
c = new Coordinate(d[0], d[1], eg);
c.Parse_Format = Parse_Format_Type.Signed_Degree;
return true;
}
catch
{//Parser failed try next method
}
}
//Try Decimal Degree
if (TryDecimalDegree(s, out d))
{
try
{
c = new Coordinate(d[0], d[1], eg);
c.Parse_Format = Parse_Format_Type.Decimal_Degree;
return true;
}
catch
{//Parser failed try next method
}
}
//Try DDM
if (TryDegreeDecimalMinute(s, out d))
{
try
{
//0 Lat Degree
//1 Lat Minute
//2 Lat Direction (0 = N, 1 = S)
//3 Long Degree
//4 Long Minute
//5 Long Direction (0 = E, 1 = W)
CoordinatesPosition latP = CoordinatesPosition.N;
CoordinatesPosition lngP = CoordinatesPosition.E;
if (d[2] != 0) { latP = CoordinatesPosition.S; }
if (d[5] != 0) { lngP = CoordinatesPosition.W; }
CoordinatePart lat = new CoordinatePart((int)d[0], d[1], latP);
CoordinatePart lng = new CoordinatePart((int)d[3], d[4], lngP);
c = new Coordinate(eg);
c.Latitude = lat;
c.Longitude = lng;
c.Parse_Format = Parse_Format_Type.Degree_Decimal_Minute;
return true;
}
catch
{//Parser failed try next method
}
}
//Try DMS
if (TryDegreeMinuteSecond(s, out d))
{
try
{
//0 Lat Degree
//1 Lat Minute
//2 Lat Second
//3 Lat Direction (0 = N, 1 = S)
//4 Long Degree
//5 Long Minute
//6 Long Second
//7 Long Direction (0 = E, 1 = W)
CoordinatesPosition latP = CoordinatesPosition.N;
CoordinatesPosition lngP = CoordinatesPosition.E;
if (d[3] != 0) { latP = CoordinatesPosition.S; }
if (d[7] != 0) { lngP = CoordinatesPosition.W; }
CoordinatePart lat = new CoordinatePart((int)d[0], (int)d[1], d[2], latP);
CoordinatePart lng = new CoordinatePart((int)d[4], (int)d[5], d[6], lngP);
c = new Coordinate(eg);
c.Latitude = lat;
c.Longitude = lng;
c.Parse_Format = Parse_Format_Type.Degree_Minute_Second;
return true;
}
catch
{//Parser failed try next method
}
}
string[] um;
//Try MGRS
if (TryMGRS(s, out um))
{
try
{
double zone = Convert.ToDouble(um[0]);
double easting = Convert.ToDouble(um[3]);
double northing = Convert.ToDouble(um[4]);
MilitaryGridReferenceSystem mgrs = new MilitaryGridReferenceSystem(um[1], (int)zone, um[2], easting, northing);
c = MilitaryGridReferenceSystem.MGRStoLatLong(mgrs);
c.Parse_Format = Parse_Format_Type.MGRS;
return true;
}
catch
{//Parser failed try next method
}
}
//Try UTM
if (TryUTM(s, out um))
{
try
{
double zone = Convert.ToDouble(um[0]);
double easting = Convert.ToDouble(um[2]);
double northing = Convert.ToDouble(um[3]);
UniversalTransverseMercator utm = new UniversalTransverseMercator(um[1], (int)zone, easting, northing);
c = UniversalTransverseMercator.ConvertUTMtoLatLong(utm);
c.Parse_Format = Parse_Format_Type.UTM;
return true;
}
catch
{//Parser failed try next method
}
}
//Try Cartesian
if (TryCartesian(s.ToUpper().Replace("KM", "").Replace("X","").Replace("Y", "").Replace("Z", ""), out d))
{
if (ct == CartesianType.Cartesian)
{
try
{
Cartesian cart = new Cartesian(d[0], d[1], d[2]);
c = Cartesian.CartesianToLatLong(cart);
c.Parse_Format = Parse_Format_Type.Cartesian_Spherical;
return true;
}
catch
{//Parser failed try next method
}
}
if (ct == CartesianType.ECEF)
{
try
{
ECEF ecef = new ECEF(d[0], d[1], d[2]);
c = ECEF.ECEFToLatLong(ecef);
c.Parse_Format = Parse_Format_Type.Cartesian_ECEF;
return true;
}
catch
{//Parser failed try next method
}
}
}
c = null;
return false;
}
private static bool TrySignedDegree(string s, out double[] d)
{
d = null;
if (Regex.Matches(s, @"[a-zA-Z]").Count != 0) { return false; } //Should contain no letters
string[] sA = SpecialSplit(s,false);
double lat;
double lng;
double degLat;
double minLat; //Minutes & MinSeconds
double secLat;
int signLat = 1;
double degLng;
double minLng; //Minutes & MinSeconds
double secLng;
int signLng = 1;
switch (sA.Count())
{
case 2:
if (!double.TryParse(sA[0], out lat))
{ return false; }
if (!double.TryParse(sA[1], out lng))
{ return false; }
d = new double[] { lat, lng };
return true;
case 4:
if (!double.TryParse(sA[0], out degLat))
{ return false; }
if (!double.TryParse(sA[1], out minLat))
{ return false; }
if (!double.TryParse(sA[2], out degLng))
{ return false; }
if (!double.TryParse(sA[3], out minLng))
{ return false; }
if (degLat < 0) { signLat = -1; }
if (degLng < 0) { signLng = -1; }
if (minLat >= 60 || minLat < 0) { return false; } //Handle in parser as degree will be incorrect.
if (minLng >= 60 || minLng < 0) { return false; } //Handle in parser as degree will be incorrect.
lat = (Math.Abs(degLat) + (minLat / 60.0)) * signLat;
lng = (Math.Abs(degLng) + (minLng / 60.0)) * signLng;
d = new double[] { lat, lng };
return true;
case 6:
if (!double.TryParse(sA[0], out degLat))
{ return false; }
if (!double.TryParse(sA[1], out minLat))
{ return false; }
if (!double.TryParse(sA[2], out secLat))
{ return false; }
if (!double.TryParse(sA[3], out degLng))
{ return false; }
if (!double.TryParse(sA[4], out minLng))
{ return false; }
if (!double.TryParse(sA[5], out secLng))
{ return false; }
if (degLat < 0) { signLat = -1; }
if (degLng < 0) { signLng = -1; }
if (minLat >= 60 || minLat < 0) { return false; } //Handle in parser as degree will be incorrect.
if (minLng >= 60 || minLng < 0) { return false; } //Handle in parser as degree will be incorrect.
if (secLat >= 60 || secLat < 0) { return false; } //Handle in parser as degree will be incorrect.
if (secLng >= 60 || secLng < 0) { return false; } //Handle in parser as degree will be incorrect.
lat = (Math.Abs(degLat) + (minLat / 60.0) + (secLat / 3600)) * signLat;
lng = (Math.Abs(degLng) + (minLng / 60.0) + (secLng / 3600)) * signLng;
d = new double[] { lat, lng };
return true;
default:
return false;
}
}
private static bool TryDecimalDegree(string s, out double[] d)
{
d = null;
if (Regex.Matches(s, @"[a-zA-Z]").Count != 2) { return false; } //Should only contain 1 letter.
string[] sA = SpecialSplit(s,true);
if (sA.Count() == 2 || sA.Count() == 4)
{
double lat;
double lng;
double latR = 1; //Sets negative if South
double lngR = 1; //Sets negative if West
//Contact get brin directional indicator together with string
if (sA.Count() == 4)
{
sA[0] += sA[1];
sA[1] = sA[2] + sA[3];
}
//Find Directions
if (!sA[0].Contains("N") && !sA[0].Contains("n"))
{
if (!sA[0].Contains("S") && !sA[0].Contains("s"))
{
return false;//No Direction Found
}
latR = -1;
}
if (!sA[1].Contains("E") && !sA[1].Contains("e"))
{
if (!sA[1].Contains("W") && !sA[1].Contains("w"))
{
return false;//No Direction Found
}
lngR = -1;
}
sA[0] = Regex.Replace(sA[0], "[^0-9.]", "");
sA[1] = Regex.Replace(sA[1], "[^0-9.]", "");
if (!double.TryParse(sA[0], out lat))
{ return false; }
if (!double.TryParse(sA[1], out lng))
{ return false; }
lat *= latR;
lng *= lngR;
d = new double[] { lat, lng };
return true;
}
return false;
}
private static bool TryDegreeDecimalMinute(string s, out double[] d)
{
d = null;
if (Regex.Matches(s, @"[a-zA-Z]").Count != 2) { return false; } //Should only contain 1 letter.
string[] sA = SpecialSplit(s,true);
if (sA.Count() == 4 || sA.Count() == 6)
{
double latD;
double latMS;
double lngD;
double lngMS;
double latR = 0; //Sets 1 if South
double lngR = 0; //Sets 1 if West
//Contact get in order to combine directional indicator together with string
//Should reduce 6 items to 4
if (sA.Count() == 6)
{
if (char.IsLetter(sA[0][0])) { sA[0] += sA[1]; sA[1] = sA[2]; }
else if (char.IsLetter(sA[1][0])) { sA[0] += sA[1]; sA[1] = sA[2]; }
else if (char.IsLetter(sA[2][0])) { sA[0] += sA[2]; }
else { return false; }
if (char.IsLetter(sA[3][0])) { sA[3] += sA[4]; sA[4] = sA[5]; }
else if (char.IsLetter(sA[4][0])) { sA[3] += sA[4]; sA[4] = sA[5]; }
else if (char.IsLetter(sA[5][0])) { sA[3] += sA[5]; }
else { return false; }
//Shift values for below logic
sA[2] = sA[3];
sA[3] = sA[4];
}
//Find Directions
if (!sA[0].Contains("N") && !sA[0].Contains("n") && !sA[1].Contains("N") && !sA[1].Contains("n"))
{
if (!sA[0].Contains("S") && !sA[0].Contains("s") && !sA[1].Contains("S") && !sA[1].Contains("s"))
{
return false;//No Direction Found
}
latR = 1;
}
if (!sA[2].Contains("E") && !sA[2].Contains("e") && !sA[3].Contains("E") && !sA[3].Contains("e"))
{
if (!sA[2].Contains("W") && !sA[2].Contains("w") && !sA[3].Contains("W") && !sA[3].Contains("w"))
{
return false;//No Direction Found
}
lngR = 1;
}
sA[0] = Regex.Replace(sA[0], "[^0-9.]", "");
sA[1] = Regex.Replace(sA[1], "[^0-9.]", "");
sA[2] = Regex.Replace(sA[2], "[^0-9.]", "");
sA[3] = Regex.Replace(sA[3], "[^0-9.]", "");
if (!double.TryParse(sA[0], out latD))
{ return false; }
if (!double.TryParse(sA[1], out latMS))
{ return false; }
if (!double.TryParse(sA[2], out lngD))
{ return false; }
if (!double.TryParse(sA[3], out lngMS))
{ return false; }
d = new double[] { latD, latMS, latR, lngD, lngMS, lngR };
return true;
}
return false;
}
private static bool TryDegreeMinuteSecond(string s, out double[] d)
{
d = null;
if (Regex.Matches(s, @"[a-zA-Z]").Count != 2) { return false; } //Should only contain 1 letter.
string[] sA = SpecialSplit(s,true);
if (sA.Count() == 6 || sA.Count() == 8)
{
double latD;
double latM;
double latS;
double lngD;
double lngM;
double lngS;
double latR = 0; //Sets 1 if South
double lngR = 0; //Sets 1 if West
//Contact get in order to combine directional indicator together with string
//Should reduce 8 items to 6
if (sA.Count() == 8)
{
if (char.IsLetter(sA[0][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else if (char.IsLetter(sA[1][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else if (char.IsLetter(sA[3][0])) { sA[0] += sA[3]; }
else { return false; }
if (char.IsLetter(sA[4][0])) { sA[4] += sA[5]; sA[5] = sA[6]; sA[6] = sA[7]; }
else if (char.IsLetter(sA[5][0])) { sA[4] += sA[5]; sA[5] = sA[6]; sA[6] = sA[7]; }
else if (char.IsLetter(sA[7][0])) { sA[4] += sA[7]; }
else { return false; }
//Shift values for below logic
sA[3] = sA[4];
sA[4] = sA[5];
sA[5] = sA[6];
}
//Find Directions
if (!sA[0].Contains("N") && !sA[0].Contains("n") && !sA[2].Contains("N") && !sA[2].Contains("n"))
{
if (!sA[0].Contains("S") && !sA[0].Contains("s") && !sA[2].Contains("S") && !sA[2].Contains("s"))
{
return false;//No Direction Found
}
latR = 1;
}
if (!sA[3].Contains("E") && !sA[3].Contains("e") && !sA[5].Contains("E") && !sA[5].Contains("e"))
{
if (!sA[3].Contains("W") && !sA[3].Contains("w") && !sA[5].Contains("W") && !sA[5].Contains("w"))
{
return false;//No Direction Found
}
lngR = 1;
}
sA[0] = Regex.Replace(sA[0], "[^0-9.]", "");
sA[1] = Regex.Replace(sA[1], "[^0-9.]", "");
sA[2] = Regex.Replace(sA[2], "[^0-9.]", "");
sA[3] = Regex.Replace(sA[3], "[^0-9.]", "");
sA[4] = Regex.Replace(sA[4], "[^0-9.]", "");
sA[5] = Regex.Replace(sA[5], "[^0-9.]", "");
if (!double.TryParse(sA[0], out latD))
{ return false; }
if (!double.TryParse(sA[1], out latM))
{ return false; }
if (!double.TryParse(sA[2], out latS))
{ return false; }
if (!double.TryParse(sA[3], out lngD))
{ return false; }
if (!double.TryParse(sA[4], out lngM))
{ return false; }
if (!double.TryParse(sA[5], out lngS))
{ return false; }
d = new double[] { latD, latM, latS, latR, lngD, lngM, lngS, lngR };
return true;
}
return false;
}
private static bool TryUTM(string s, out string[] utm)
{
utm = null;
string[] sA = SpecialSplit(s,false);
if (sA.Count() == 3 || sA.Count() == 4)
{
double zone;
string zoneL;
double easting;
double northing;
if (sA.Count() == 4)
{
if (char.IsLetter(sA[0][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else if (char.IsLetter(sA[1][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else { return false; }
}
zoneL = new string(sA[0].Where(Char.IsLetter).ToArray());
if (zoneL == string.Empty) { return false; }
sA[0] = Regex.Replace(sA[0], "[^0-9.]", "");
if (!double.TryParse(sA[0], out zone))
{ return false; }
if (!double.TryParse(sA[1], out easting))
{ return false; }
if (!double.TryParse(sA[2], out northing))
{ return false; }
utm = new string[] { zone.ToString(), zoneL, easting.ToString(), northing.ToString() };
return true;
}
return false;
}
private static bool TryMGRS(string s, out string[] mgrs)
{
mgrs = null;
string[] sA = SpecialSplit(s,false);
if (sA.Count() == 4 || sA.Count() == 5)
{
double zone;
string zoneL;
string diagraph;
double easting;
double northing;
if (sA.Count() == 5)
{
if (char.IsLetter(sA[0][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else if (char.IsLetter(sA[1][0])) { sA[0] += sA[1]; sA[1] = sA[2]; sA[2] = sA[3]; }
else { return false; }
}
zoneL = new string(sA[0].Where(Char.IsLetter).ToArray());
if (zoneL == string.Empty) { return false; }
sA[0] = Regex.Replace(sA[0], "[^0-9.]", "");
diagraph = sA[1];
if (!double.TryParse(sA[0], out zone))
{ return false; }
if (!double.TryParse(sA[2], out easting))
{ return false; }
if (!double.TryParse(sA[3], out northing))
{ return false; }
mgrs = new string[] { zone.ToString(), zoneL, diagraph, easting.ToString(), northing.ToString() };
return true;
}
return false;
}
private static bool TryCartesian(string s, out double[] d)
{
d = null;
string[] sA = SpecialSplit(s,false);
if (sA.Count() == 3)
{
double x;
double y;
double z;
if (!double.TryParse(sA[0], out x))
{ return false; }
if (!double.TryParse(sA[1], out y))
{ return false; }
if (!double.TryParse(sA[2], out z))
{ return false; }
d = new double[] { x, y, z };
return true;
}
return false;
}
//KEEP DASHES FOR SIGNED AND CARTESIAN AS THEY ARE USED FOR NEGATVE VALUES
private static string[] SpecialSplit(string s, bool removeDashes)
{
s = s.Replace("°", " ");
s = s.Replace("º", " ");
s = s.Replace("'", " ");
s = s.Replace("\"", " ");
s = s.Replace(",", " ");
s = s.Replace("mE", " ");
s = s.Replace("mN", " ");
if (removeDashes)
{
s = s.Replace("-", " ");
}
return s.Split(new char[0], StringSplitOptions.RemoveEmptyEntries);
}
}
internal class FormatFinder_CoordPart
{
//Add main to Coordinate and tunnel to Format class. Add private methods to format.
//WHEN PARSING NO EXCPETIONS FOR OUT OF RANGE ARGS WILL BE THROWN
public static bool TryParse(string coordString, out CoordinatePart cp)
{
//Turn of eagerload for efficiency
EagerLoad eg = new EagerLoad();
int type = 0; //0 = unspecifed, 1 = lat, 2 = long;
eg.Cartesian = false;
eg.Celestial = false;
eg.UTM_MGRS = false;
cp = null;
Coordinate c = new Coordinate(eg);
string s = coordString;
s = s.Trim(); //Trim all spaces before and after string
double[] d;
if (s[0] == ',')
{
type = 2;
s = s.Replace(",", "");
s = s.Trim();
}
if (s[0] == '*')
{
type = 1;
s = s.Replace("*", "");
s = s.Trim();
}
if (TrySignedDegree(s, type, out d))
{
try
{
switch (type)
{
case 0:
//Attempt Lat first (default for signed)
try
{
cp = new CoordinatePart(d[0], CoordinateType.Lat);
c.Parse_Format = Parse_Format_Type.Signed_Degree;
return true;
}
catch
{
cp = new CoordinatePart(d[0], CoordinateType.Long);
c.Parse_Format = Parse_Format_Type.Signed_Degree;
return true;
}
case 1:
//Attempt Lat
cp = new CoordinatePart(d[0], CoordinateType.Lat);
c.Parse_Format = Parse_Format_Type.Signed_Degree;
return true;
case 2:
//Attempt long
cp = new CoordinatePart(d[0], CoordinateType.Long);
c.Parse_Format = Parse_Format_Type.Signed_Degree;
return true;
}
}
catch
{
//silent fail
}
}
//SIGNED DEGREE FAILED, REMOVE DASHES FOR OTHER FORMATS
s = s.Replace("-", " ");
//All other formats should contain 1 letter.
if (Regex.Matches(s, @"[a-zA-Z]").Count != 1) { return false; } //Should only contain 1 letter.
//Get Coord Direction
int direction = Find_Position(s);
if (direction == -1)
{
return false; //No direction found
}
//If Coordinate type int specified, look for mismatch
if (type == 1 && (direction == 1 || direction == 3))
{
return false; //mismatch
}
if (type == 2 && (direction == 0 || direction == 2))
{
return false; //mismatch
}
CoordinateType t;
if (direction == 0 || direction == 2) { t = CoordinateType.Lat; }
else { t = CoordinateType.Long; }
s = Regex.Replace(s, "[^0-9. ]", ""); //Remove directional character
s = s.Trim(); //Trim all spaces before and after string
//Try Decimal Degree with Direction
if (TryDecimalDegree(s, direction, out d))
{
try
{
cp = new CoordinatePart(d[0], t);
c.Parse_Format = Parse_Format_Type.Decimal_Degree;
return true;
}
catch
{//Parser failed try next method
}
}
//Try DDM
if (TryDegreeDecimalMinute(s, out d))
{
try
{
//0 Degree
//1 Minute
//2 Direction (0 = N, 1 = E, 2 = S, 3 = W)
cp = new CoordinatePart((int)d[0], d[1], (CoordinatesPosition)direction);
c.Parse_Format = Parse_Format_Type.Degree_Decimal_Minute;
return true;
}
catch
{
//Parser failed try next method
}
}
//Try DMS
if (TryDegreeMinuteSecond(s, out d))
{
try
{
//0 Degree
//1 Minute
//2 Second
//3 Direction (0 = N, 1 = E, 2 = S, 3 = W)
cp = new CoordinatePart((int)d[0], (int)d[1], d[2], (CoordinatesPosition)direction);
c.Parse_Format = Parse_Format_Type.Degree_Minute_Second;
return true;
}
catch
{//Parser failed try next method
}
}
return false;
}
private static bool TrySignedDegree(string s, int t, out double[] d)
{
d = null;
if (Regex.Matches(s, @"[a-zA-Z]").Count != 0) { return false; } //Should contain no letters
string[] sA = SpecialSplit(s, false);
double deg;
double min; //Minutes & MinSeconds
double sec;
int sign = 1;
switch (sA.Count())
{
case 1:
if (!double.TryParse(sA[0], out deg))
{ return false; }
d = new double[] { deg };
return true;
case 2:
if (!double.TryParse(sA[0], out deg))
{ return false; }
if (!double.TryParse(sA[1], out min))
{ return false; }
if (deg < 0) { sign = -1; }
if (min >= 60 || min < 0) { return false; } //Handle in parser as degree will be incorrect.
d = new double[] { (Math.Abs(deg) + (min / 60.0)) * sign };
return true;
case 3:
if (!double.TryParse(sA[0], out deg))
{ return false; }
if (!double.TryParse(sA[1], out min))
{ return false; }
if (!double.TryParse(sA[2], out sec))
{ return false; }
if (min >= 60 || min < 0) { return false; } //Handle in parser as degree will be incorrect.
if (sec >= 60 || sec < 0) { return false; } //Handle in parser as degree will be incorrect.
if (deg < 0) { sign = -1; }
d = new double[] { (Math.Abs(deg) + (min / 60.0) + (sec / 3600.0)) * sign };
return true;
default:
return false;
}
}
private static bool TryDecimalDegree(string s, int direction, out double[] d)
{
d = null;
int sign = 1;
//S or W
if (direction == 2 || direction == 3)
{
sign = -1;
}
double coord;
string[] sA = SpecialSplit(s, true);
if (sA.Count() == 1)
{
if (!double.TryParse(s, out coord))
{ return false; }
coord *= sign;
d = new double[] { coord };
return true;
}
return false;
}
private static bool TryDegreeDecimalMinute(string s, out double[] d)
{
d = null;
double deg;
double minSec;
string[] sA = SpecialSplit(s,true);
if (sA.Count() == 2)
{
if (!double.TryParse(sA[0], out deg))
{ return false; }
if (!double.TryParse(sA[1], out minSec))
{ return false; }
d = new double[] { deg, minSec };
return true;
}
return false;
}
private static bool TryDegreeMinuteSecond(string s, out double[] d)
{
d = null;
double deg;
double min;
double sec;
string[] sA = SpecialSplit(s,true);
if (sA.Count() == 3)
{
if (!double.TryParse(sA[0], out deg))
{ return false; }
if (!double.TryParse(sA[1], out min))
{ return false; }
if (!double.TryParse(sA[2], out sec))
{ return false; }
d = new double[] { deg, min, sec };
return true;
}
return false;
}
private static int Find_Position(string s)
{
//N=0
//E=1
//S=2
//W=3
//NOPOS = -1
//Find Directions
int part = -1;
if (s.Contains("N") || s.Contains("n"))
{
part = 0;
}
if (s.Contains("E") || s.Contains("e"))
{
part = 1;
}
if (s.Contains("S") || s.Contains("s"))
{
part = 2;
}
if (s.Contains("W") || s.Contains("w"))
{
part = 3;
}
return part;
}
//KEEP DASHES FOR SIGNED AND CARTESIAN AS THEY ARE USED FOR NEGATVE VALUES
private static string[] SpecialSplit(string s, bool removeDashes)
{
s = s.Replace("°", " ");
s = s.Replace("º", " ");
s = s.Replace("'", " ");
s = s.Replace("\"", " ");
s = s.Replace(",", " ");
s = s.Replace("mE", " ");
s = s.Replace("mN", " ");
if(removeDashes)
{
s = s.Replace("-", " ");
}
return s.Split(new char[0], StringSplitOptions.RemoveEmptyEntries);
}
}
}

592
Lora-Map/Helper/Coordinate.UTM.cs Normal file
View File

@ -0,0 +1,592 @@
using System;
using System.Linq;
using System.Diagnostics;
using System.ComponentModel;
namespace CoordinateSharp
{
/// <summary>
/// Universal Transverse Mercator (UTM) coordinate system. Uses the WGS 84 Datum.
/// </summary>
[Serializable]
public class UniversalTransverseMercator : INotifyPropertyChanged
{
/// <summary>
/// Creates a UniversalTransverMercator object with a WGS84 Datum.
/// </summary>
/// <param name="latz">Latitude zone</param>
/// <param name="longz">Longitude zone</param>
/// <param name="est">Easting</param>
/// <param name="nrt">Northing</param>
public UniversalTransverseMercator(string latz, int longz, double est, double nrt)
{
if (longz < 1 || longz > 60) { Debug.WriteLine("Longitudinal zone out of range", "UTM longitudinal zones must be between 1-60."); }
if (!Verify_Lat_Zone(latz)) { Debug.WriteLine("Latitudinal zone invalid", "UTM latitudinal zone was unrecognized."); }
if (est < 160000 || est > 834000) { Debug.WriteLine("The Easting value provided is outside the max allowable range. Use with caution."); }
if (nrt < 0 || nrt > 10000000) { Debug.WriteLine("Northing out of range", "Northing must be between 0-10,000,000."); }
latZone = latz;
longZone =longz;
easting = est;
northing = nrt;
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
}
/// <summary>
/// Creates a UniversalTransverMercator object with a custom Datum.
/// </summary>
/// <param name="latz">Latitude zone</param>
/// <param name="longz">Longitude zone</param>
/// <param name="est">Easting</param>
/// <param name="nrt">Northing</param>
/// <param name="radius">Equatorial Radius</param>
/// <param name="flaten">Inverse Flattening</param>
public UniversalTransverseMercator(string latz, int longz, double est, double nrt, double radius, double flaten)
{
if (longz < 1 || longz > 60) { Debug.WriteLine("Longitudinal zone out of range", "UTM longitudinal zones must be between 1-60."); }
if (!Verify_Lat_Zone(latz)) { Debug.WriteLine("Latitudinal zone invalid", "UTM latitudinal zone was unrecognized."); }
if (est < 160000 || est > 834000) { Debug.WriteLine("The Easting value provided is outside the max allowable range. Use with caution."); }
if (nrt < 0 || nrt > 10000000) { Debug.WriteLine("Northing out of range", "Northing must be between 0-10,000,000."); }
latZone = latz;
longZone = longz;
easting = est;
northing = nrt;
equatorial_radius = radius;
inverse_flattening = flaten;
}
private Coordinate coordinate;
internal double equatorial_radius;
internal double inverse_flattening;
private string latZone;
private int longZone;
private double easting;
private double northing;
private bool withinCoordinateSystemBounds = true;
/// <summary>
/// UTM Zone Letter
/// </summary>
public string LatZone
{
get { return latZone; }
set
{
if (latZone != value)
{
latZone = value;
}
}
}
/// <summary>
/// UTM Zone Number
/// </summary>
public int LongZone
{
get { return longZone; }
set
{
if (longZone != value)
{
longZone = value;
}
}
}
/// <summary>
/// UTM Easting
/// </summary>
public double Easting
{
get { return easting; }
set
{
if (easting != value)
{
easting = value;
}
}
}
/// <summary>
/// UTM Northing
/// </summary>
public double Northing
{
get { return northing; }
set
{
if (northing != value)
{
northing = value;
}
}
}
/// <summary>
/// Datum Equatorial Radius / Semi Major Axis
/// </summary>
public double Equatorial_Radius
{
get { return equatorial_radius; }
}
/// <summary>
/// Datum Flattening
/// </summary>
public double Inverse_Flattening
{
get { return inverse_flattening; }
}
/// <summary>
/// Is the UTM conversion within the coordinate system's accurate boundaries after conversion from Lat/Long.
/// </summary>
public bool WithinCoordinateSystemBounds
{
get { return withinCoordinateSystemBounds; }
}
/// <summary>
/// Constructs a UTM object based off DD Lat/Long
/// </summary>
/// <param name="lat">DD Latitude</param>
/// <param name="longi">DD Longitide</param>
/// <param name="c">Parent Coordinate Object</param>
internal UniversalTransverseMercator(double lat, double longi, Coordinate c)
{
//validate coords
//if (lat > 180) { throw new ArgumentOutOfRangeException("Degrees out of range", "Longitudinal coordinate decimal cannot be greater than 180."); }
//if (lat < -180) { throw new ArgumentOutOfRangeException("Degrees out of range", "Longitudinal coordinate decimal cannot be less than 180."); }
//if (longi > 90) { throw new ArgumentOutOfRangeException("Degrees out of range", "Latitudinal coordinate decimal cannot be greater than 90."); }
//if (longi < -90) { throw new ArgumentOutOfRangeException("Degrees out of range", "Latitudinal coordinate decimal cannot be less than 90."); }
equatorial_radius = 6378137.0;
inverse_flattening = 298.257223563;
ToUTM(lat, longi, this);
coordinate = c;
}
/// <summary>
/// Constructs a UTM object based off DD Lat/Long
/// </summary>
/// <param name="lat">DD Latitude</param>
/// <param name="longi">DD Longitide</param>
/// <param name="c">Parent Coordinate Object</param>
/// <param name="rad">Equatorial Radius</param>
/// <param name="flt">Flattening</param>
internal UniversalTransverseMercator(double lat, double longi, Coordinate c,double rad,double flt)
{
equatorial_radius = rad;
inverse_flattening = flt;
ToUTM(lat, longi, this);
coordinate = c;
}
/// <summary>
/// Constructs a UTM object based off a UTM coordinate
/// Not yet implemented
/// </summary>
/// <param name="latz">Zone Letter</param>
/// <param name="longz">Zone Number</param>
/// <param name="e">Easting</param>
/// <param name="n">Northing</param>
/// <param name="c">Parent Coordinate Object</param>
/// <param name="rad">Equatorial Radius</param>
/// <param name="flt">Inverse Flattening</param>
internal UniversalTransverseMercator(string latz, int longz, double e, double n, Coordinate c, double rad, double flt)
{
//validate utm
if (longz < 1 || longz > 60) { Debug.WriteLine("Longitudinal zone out of range", "UTM longitudinal zones must be between 1-60."); }
if (!Verify_Lat_Zone(latz)) { throw new ArgumentException("Latitudinal zone invalid", "UTM latitudinal zone was unrecognized."); }
if (e < 160000 || e > 834000) { Debug.WriteLine("The Easting value provided is outside the max allowable range. If this is intentional, use with caution."); }
if (n < 0 || n > 10000000) { throw new ArgumentOutOfRangeException("Northing out of range", "Northing must be between 0-10,000,000."); }
equatorial_radius = rad;
inverse_flattening = flt;
latZone = latz;
longZone = longz;
easting = e;
northing = n;
coordinate = c;
if (c.Latitude.DecimalDegree <= -80 || c.Latitude.DecimalDegree >= 84) { withinCoordinateSystemBounds = false; }
else { withinCoordinateSystemBounds = true; }
}
/// <summary>
/// Verifies Lat zone when convert from UTM to DD Lat/Long
/// </summary>
/// <param name="l">Zone Letter</param>
/// <returns>boolean</returns>
private bool Verify_Lat_Zone(string l)
{
if (LatZones.longZongLetters.Where(x => x == l.ToUpper()).Count() != 1)
{
return false;
}
return true;
}
private double degreeToRadian(double degree)
{
return degree * Math.PI / 180;
}
/// <summary>
/// Assigns UTM values based of Lat/Long
/// </summary>
/// <param name="lat">DD Latitude</param>
/// <param name="longi">DD longitude</param>
/// <param name="utm">UTM Object to modify</param>
internal void ToUTM(double lat, double longi, UniversalTransverseMercator utm)
{
string letter = "";
double easting = 0;
double northing = 0;
int zone = (int)Math.Floor(longi / 6 + 31);
if (lat < -72)
letter = "C";
else if (lat < -64)
letter = "D";
else if (lat < -56)
letter = "E";
else if (lat < -48)
letter = "F";
else if (lat < -40)
letter = "G";
else if (lat < -32)
letter = "H";
else if (lat < -24)
letter = "J";
else if (lat < -16)
letter = "K";
else if (lat < -8)
letter = "L";
else if (lat < 0)
letter = "M";
else if (lat < 8)
letter = "N";
else if (lat < 16)
letter = "P";
else if (lat < 24)
letter = "Q";
else if (lat < 32)
letter = "R";
else if (lat < 40)
letter = "S";
else if (lat < 48)
letter = "T";
else if (lat < 56)
letter = "U";
else if (lat < 64)
letter = "V";
else if (lat < 72)
letter = "W";
else
letter = "X";
double a = utm.equatorial_radius;
double f = 1.0 / utm.inverse_flattening;
double b = a * (1 - f); // polar radius
double e = Math.Sqrt(1 - Math.Pow(b, 2) / Math.Pow(a, 2));
double e0 = e / Math.Sqrt(1 - Math.Pow(e, 1));
double drad = Math.PI / 180;
double k0 = 0.9996;
double phi = lat * drad; // convert latitude to radians
double lng = longi * drad; // convert longitude to radians
double utmz = 1 + Math.Floor((longi + 180) / 6.0); // longitude to utm zone
double zcm = 3 + 6.0 * (utmz - 1) - 180; // central meridian of a zone
// this gives us zone A-B for below 80S
double esq = (1 - (b / a) * (b / a));
double e0sq = e * e / (1 - Math.Pow(e, 2));
double M = 0;
double N = a / Math.Sqrt(1 - Math.Pow(e * Math.Sin(phi), 2));
double T = Math.Pow(Math.Tan(phi), 2);
double C = e0sq * Math.Pow(Math.Cos(phi), 2);
double A = (longi - zcm) * drad * Math.Cos(phi);
// calculate M (USGS style)
M = phi * (1 - esq * (1.0 / 4.0 + esq * (3.0 / 64.0 + 5.0 * esq / 256.0)));
M = M - Math.Sin(2.0 * phi) * (esq * (3.0 / 8.0 + esq * (3.0 / 32.0 + 45.0 * esq / 1024.0)));
M = M + Math.Sin(4.0 * phi) * (esq * esq * (15.0 / 256.0 + esq * 45.0 / 1024.0));
M = M - Math.Sin(6.0 * phi) * (esq * esq * esq * (35.0 / 3072.0));
M = M * a;//Arc length along standard meridian
double M0 = 0;// if another point of origin is used than the equator
// Calculate the UTM values...
// first the easting
var x = k0 * N * A * (1 + A * A * ((1 - T + C) / 6 + A * A * (5 - 18 * T + T * T + 72.0 * C - 58 * e0sq) / 120.0)); //Easting relative to CM
x = x + 500000; // standard easting
// Northing
double y = k0 * (M - M0 + N * Math.Tan(phi) * (A * A * (1 / 2.0 + A * A * ((5 - T + 9 * C + 4 * C * C) / 24.0 + A * A * (61 - 58 * T + T * T + 600 * C - 330 * e0sq) / 720.0)))); // first from the equator
double yg = y + 10000000; //yg = y global, from S. Pole
if (y < 0)
{
y = 10000000 + y; // add in false northing if south of the equator
}
easting = Math.Round(10 * (x)) / 10.0;
northing = Math.Round(10 * y) / 10.0;
utm.latZone = letter;
utm.longZone = zone;
utm.easting = easting;
utm.northing = northing;
if(lat<=-80 || lat >= 84) { withinCoordinateSystemBounds = false; }
else { withinCoordinateSystemBounds = true; }
}
/// <summary>
/// UTM Default String Format
/// </summary>
/// <returns>UTM Formatted Coordinate String</returns>
public override string ToString()
{
if (!withinCoordinateSystemBounds) { return ""; }//MGRS Coordinate is outside its reliable boundaries. Return empty.
return longZone.ToString() + LatZone + " " + (int)easting + "mE " + (int)northing + "mN";
}
private static Coordinate UTMtoLatLong(double x, double y, double zone, double equatorialRadius, double flattening)
{
//x easting
//y northing
//http://home.hiwaay.net/~taylorc/toolbox/geography/geoutm.html
double phif, Nf, Nfpow, nuf2, ep2, tf, tf2, tf4, cf;
double x1frac, x2frac, x3frac, x4frac, x5frac, x6frac, x7frac, x8frac;
double x2poly, x3poly, x4poly, x5poly, x6poly, x7poly, x8poly;
double sm_a = equatorialRadius;
double sm_b = equatorialRadius * (1 - (1.0 / flattening)); //Polar Radius
/* Get the value of phif, the footpoint latitude. */
phif = FootpointLatitude(y,equatorialRadius,flattening);
/* Precalculate ep2 */
ep2 = (Math.Pow(sm_a, 2.0) - Math.Pow(sm_b, 2.0))
/ Math.Pow(sm_b, 2.0);
/* Precalculate cos (phif) */
cf = Math.Cos(phif);
/* Precalculate nuf2 */
nuf2 = ep2 * Math.Pow(cf, 2.0);
/* Precalculate Nf and initialize Nfpow */
Nf = Math.Pow(sm_a, 2.0) / (sm_b * Math.Sqrt(1 + nuf2));
Nfpow = Nf;
/* Precalculate tf */
tf = Math.Tan(phif);
tf2 = tf * tf;
tf4 = tf2 * tf2;
/* Precalculate fractional coefficients for x**n in the equations
below to simplify the expressions for latitude and longitude. */
x1frac = 1.0 / (Nfpow * cf);
Nfpow *= Nf; /* now equals Nf**2) */
x2frac = tf / (2.0 * Nfpow);
Nfpow *= Nf; /* now equals Nf**3) */
x3frac = 1.0 / (6.0 * Nfpow * cf);
Nfpow *= Nf; /* now equals Nf**4) */
x4frac = tf / (24.0 * Nfpow);
Nfpow *= Nf; /* now equals Nf**5) */
x5frac = 1.0 / (120.0 * Nfpow * cf);
Nfpow *= Nf; /* now equals Nf**6) */
x6frac = tf / (720.0 * Nfpow);
Nfpow *= Nf; /* now equals Nf**7) */
x7frac = 1.0 / (5040.0 * Nfpow * cf);
Nfpow *= Nf; /* now equals Nf**8) */
x8frac = tf / (40320.0 * Nfpow);
/* Precalculate polynomial coefficients for x**n.
-- x**1 does not have a polynomial coefficient. */
x2poly = -1.0 - nuf2;
x3poly = -1.0 - 2 * tf2 - nuf2;
x4poly = 5.0 + 3.0 * tf2 + 6.0 * nuf2 - 6.0 * tf2 * nuf2
- 3.0 * (nuf2 * nuf2) - 9.0 * tf2 * (nuf2 * nuf2);
x5poly = 5.0 + 28.0 * tf2 + 24.0 * tf4 + 6.0 * nuf2 + 8.0 * tf2 * nuf2;
x6poly = -61.0 - 90.0 * tf2 - 45.0 * tf4 - 107.0 * nuf2
+ 162.0 * tf2 * nuf2;
x7poly = -61.0 - 662.0 * tf2 - 1320.0 * tf4 - 720.0 * (tf4 * tf2);
x8poly = 1385.0 + 3633.0 * tf2 + 4095.0 * tf4 + 1575 * (tf4 * tf2);
/* Calculate latitude */
double nLat = phif + x2frac * x2poly * (x * x)
+ x4frac * x4poly * Math.Pow(x, 4.0)
+ x6frac * x6poly * Math.Pow(x, 6.0)
+ x8frac * x8poly * Math.Pow(x, 8.0);
/* Calculate longitude */
double nLong = zone + x1frac * x
+ x3frac * x3poly * Math.Pow(x, 3.0)
+ x5frac * x5poly * Math.Pow(x, 5.0)
+ x7frac * x7poly * Math.Pow(x, 7.0);
double dLat = RadToDeg(nLat);
double dLong = RadToDeg(nLong);
if (dLat > 90) { dLat = 90; }
if (dLat < -90) { dLat = -90; }
if (dLong > 180) { dLong = 180; }
if (dLong < -180) { dLong = -180; }
Coordinate c = new Coordinate(equatorialRadius,flattening, true);
CoordinatePart cLat = new CoordinatePart(dLat, CoordinateType.Lat);
CoordinatePart cLng = new CoordinatePart(dLong, CoordinateType.Long);
c.Latitude = cLat;
c.Longitude = cLng;
return c;
}
private static double RadToDeg(double rad)
{
double pi = 3.14159265358979;
return (rad / pi * 180.0);
}
private static double DegToRad(double deg)
{
double pi = 3.14159265358979;
return (deg / 180.0 * pi);
}
private static double FootpointLatitude(double y, double equatorialRadius, double flattening)
{
double y_, alpha_, beta_, gamma_, delta_, epsilon_, n;
double result;
/* Ellipsoid model constants (actual values here are for WGS84) */
double sm_a = equatorialRadius;
double sm_b = equatorialRadius * (1 - (1.0 / flattening));
/* Precalculate n (Eq. 10.18) */
n = (sm_a - sm_b) / (sm_a + sm_b);
/* Precalculate alpha_ (Eq. 10.22) */
/* (Same as alpha in Eq. 10.17) */
alpha_ = ((sm_a + sm_b) / 2.0) * (1 + (Math.Pow(n, 2.0) / 4) + (Math.Pow(n, 4.0) / 64));
/* Precalculate y_ (Eq. 10.23) */
y_ = y / alpha_;
/* Precalculate beta_ (Eq. 10.22) */
beta_ = (3.0 * n / 2.0) + (-27.0 * Math.Pow(n, 3.0) / 32.0)
+ (269.0 * Math.Pow(n, 5.0) / 512.0);
/* Precalculate gamma_ (Eq. 10.22) */
gamma_ = (21.0 * Math.Pow(n, 2.0) / 16.0)
+ (-55.0 * Math.Pow(n, 4.0) / 32.0);
/* Precalculate delta_ (Eq. 10.22) */
delta_ = (151.0 * Math.Pow(n, 3.0) / 96.0)
+ (-417.0 * Math.Pow(n, 5.0) / 128.0);
/* Precalculate epsilon_ (Eq. 10.22) */
epsilon_ = (1097.0 * Math.Pow(n, 4.0) / 512.0);
/* Now calculate the sum of the series (Eq. 10.21) */
result = y_ + (beta_ * Math.Sin(2.0 * y_))
+ (gamma_ * Math.Sin(4.0 * y_))
+ (delta_ * Math.Sin(6.0 * y_))
+ (epsilon_ * Math.Sin(8.0 * y_));
return result;
}
/// <summary>
/// Converts UTM coordinate to Lat/Long
/// </summary>
/// <param name="utm">utm</param>
/// <returns>Coordinate object</returns>
public static Coordinate ConvertUTMtoLatLong(UniversalTransverseMercator utm)
{
bool southhemi = false;
if (utm.latZone == "A" || utm.latZone == "B" || utm.latZone == "C" || utm.latZone == "D" || utm.latZone == "E" || utm.latZone == "F" || utm.latZone == "G" || utm.latZone == "H" || utm.latZone == "J" ||
utm.latZone == "K" || utm.latZone == "L" || utm.latZone == "M")
{
southhemi = true;
}
double cmeridian;
double x = utm.Easting - 500000.0;
double UTMScaleFactor = 0.9996;
x /= UTMScaleFactor;
/* If in southern hemisphere, adjust y accordingly. */
double y = utm.Northing;
if (southhemi)
{
y -= 10000000.0;
}
y /= UTMScaleFactor;
cmeridian = UTMCentralMeridian(utm.LongZone);
Coordinate c = UTMtoLatLong(x, y, cmeridian, utm.equatorial_radius, utm.inverse_flattening);
if (c.Latitude.ToDouble() > 85 || c.Latitude.ToDouble() < -85)
{
Debug.WriteLine("UTM conversions greater than 85 degrees or less than -85 degree latitude contain major deviations and should be used with caution.");
}
return c;
}
private static double UTMCentralMeridian(double zone)
{
double cmeridian;
cmeridian = DegToRad(-183.0 + (zone * 6.0));
return cmeridian;
}
/// <summary>
/// Property changed event
/// </summary>
public event PropertyChangedEventHandler PropertyChanged;
/// <summary>
/// Notify property changed
/// </summary>
/// <param name="propName">Property name</param>
public void NotifyPropertyChanged(string propName)
{
if (this.PropertyChanged != null)
{
PropertyChanged(this, new PropertyChangedEventArgs(propName));
}
}
}
}

1979
Lora-Map/Helper/Coordinate.cs Normal file
View File

File diff suppressed because it is too large Load Diff

436
Lora-Map/Helper/Distance.cs Normal file
View File

@ -0,0 +1,436 @@
using System;
using System.Diagnostics;
namespace CoordinateSharp
{
/// <summary>
/// Contains distance values between two coordinates.
/// </summary>
[Serializable]
public class Distance
{
private double kilometers;
private double miles;
private double feet;
private double meters;
private double bearing;
private double nauticalMiles;
/// <summary>
/// Initializes a distance object using Haversine (Spherical Earth).
/// </summary>
/// <param name="c1">Coordinate 1</param>
/// <param name="c2">Coordinate 2</param>
public Distance(Coordinate c1, Coordinate c2)
{
Haversine(c1, c2);
}
/// <summary>
/// Initializes a distance object using Haversine (Spherical Earth) or Vincenty (Elliptical Earth).
/// </summary>
/// <param name="c1">Coordinate 1</param>
/// <param name="c2">Coordinate 2</param>
/// <param name="shape">Shape of earth</param>
public Distance(Coordinate c1, Coordinate c2, Shape shape)
{
if (shape == Shape.Sphere)
{
Haversine(c1, c2);
}
else
{
Vincenty(c1, c2);
}
}
/// <summary>
/// Initializes distance object based on distance in KM
/// </summary>
/// <param name="km">Kilometers</param>
public Distance(double km)
{
kilometers = km;
meters = km * 1000;
feet = meters * 3.28084;
miles = meters * 0.000621371;
nauticalMiles = meters * 0.0005399565;
bearing = 0;//None specified
}
/// <summary>
/// Initializaes distance object based on specified distance and measurement type
/// </summary>
/// <param name="distance">Distance</param>
/// <param name="type">Measurement type</param>
public Distance(double distance, DistanceType type)
{
bearing = 0;
switch (type)
{
case DistanceType.Feet:
feet = distance;
meters = feet * 0.3048;
kilometers = meters / 1000;
miles = meters * 0.000621371;
nauticalMiles = meters * 0.0005399565;
break;
case DistanceType.Kilometers:
kilometers = distance;
meters = kilometers * 1000;
feet = meters * 3.28084;
miles = meters * 0.000621371;
nauticalMiles = meters * 0.0005399565;
break;
case DistanceType.Meters:
meters = distance;
kilometers = meters / 1000;
feet = meters * 3.28084;
miles = meters * 0.000621371;
nauticalMiles = meters * 0.0005399565;
break;
case DistanceType.Miles:
miles = distance;
meters = miles * 1609.344;
feet = meters * 3.28084;
kilometers = meters / 1000;
nauticalMiles = meters * 0.0005399565;
break;
case DistanceType.NauticalMiles:
nauticalMiles = distance;
meters = nauticalMiles * 1852.001;
feet = meters * 3.28084;
kilometers = meters / 1000;
miles = meters * 0.000621371;
break;
default:
kilometers = distance;
meters = distance * 1000;
feet = meters * 3.28084;
miles = meters * 0.000621371;
nauticalMiles = meters * 0.0005399565;
break;
}
}
private void Vincenty(Coordinate c1, Coordinate c2)
{
double lat1, lat2, lon1, lon2;
double d, crs12, crs21;
lat1 = c1.Latitude.ToRadians();
lat2 = c2.Latitude.ToRadians();
lon1 = c1.Longitude.ToRadians();
lon2 = c2.Longitude.ToRadians();
//Ensure datums match between coords
if ((c1.equatorial_radius != c2.equatorial_radius) || (c1.inverse_flattening != c2.inverse_flattening))
{
throw new InvalidOperationException("The datum set does not match between Coordinate objects.");
}
double[] ellipse = new double[] { c1.equatorial_radius, c1.inverse_flattening };
// elliptic code
double[] cde = Distance_Assistant.Dist_Ell(lat1, -lon1, lat2, -lon2, ellipse); // ellipse uses East negative
crs12 = cde[1] * (180 / Math.PI); //Bearing
crs21 = cde[2] * (180 / Math.PI); //Reverse Bearing
d = cde[0]; //Distance
bearing = crs21;
//reverseBearing = crs12;
meters = d;
kilometers = d / 1000;
feet = d * 3.28084;
miles = d * 0.000621371;
nauticalMiles = d * 0.0005399565;
}
private void Haversine(Coordinate c1, Coordinate c2)
{
////RADIANS
double nLat = c1.Latitude.ToDouble() * Math.PI / 180;
double nLong = c1.Longitude.ToDouble() * Math.PI / 180;
double cLat = c2.Latitude.ToDouble() * Math.PI / 180;
double cLong = c2.Longitude.ToDouble() * Math.PI / 180;
//Calcs
double R = 6371e3; //meters
double v1 = nLat;
double v2 = cLat;
double latRad = (c2.Latitude.ToDouble() - c1.Latitude.ToDouble()) * Math.PI / 180;
double longRad = (c2.Longitude.ToDouble() - c1.Longitude.ToDouble()) * Math.PI / 180;
double a = Math.Sin(latRad / 2.0) * Math.Sin(latRad / 2.0) +
Math.Cos(nLat) * Math.Cos(cLat) * Math.Sin(longRad / 2.0) * Math.Sin(longRad / 2.0);
double cl = 2 * Math.Atan2(Math.Sqrt(a), Math.Sqrt(1 - a));
double dist = R * cl;
//Get bearing
double y = Math.Sin(cLong - nLong) * Math.Cos(cLat);
double x = Math.Cos(nLat) * Math.Sin(cLat) -
Math.Sin(nLat) * Math.Cos(cLat) * Math.Cos(cLong - nLong);
double brng = Math.Atan2(y, x) * (180 / Math.PI); //Convert bearing back to degrees.
if (brng < 0) { brng -= 180; brng = Math.Abs(brng); }
kilometers = dist / 1000;
meters = dist;
feet = dist * 3.28084;
miles = dist * 0.000621371;
nauticalMiles = dist * 0.0005399565;
bearing = brng;
}
/// <summary>
/// Distance in Kilometers
/// </summary>
public double Kilometers
{
get { return kilometers; }
}
/// <summary>
/// Distance in Statute Miles
/// </summary>
public double Miles
{
get { return miles; }
}
/// <summary>
/// Distance in Nautical Miles
/// </summary>
public double NauticalMiles
{
get { return nauticalMiles; }
}
/// <summary>
/// Distance in Meters
/// </summary>
public double Meters
{
get { return meters; }
}
/// <summary>
/// Distance in Feet
/// </summary>
public double Feet
{
get { return feet; }
}
/// <summary>
/// Initial Bearing from Coordinate 1 to Coordinate 2
/// </summary>
public double Bearing
{
get { return bearing; }
}
}
/// <summary>
/// Distance measurement type
/// </summary>
public enum DistanceType
{
/// <summary>
/// Distance in Meters
/// </summary>
Meters,
/// <summary>
/// Distance in Kilometers
/// </summary>
Kilometers,
/// <summary>
/// Distance in Feet
/// </summary>
Feet,
/// <summary>
/// Distance in Statute Miles
/// </summary>
Miles,
/// <summary>
/// Distance in Nautical Miles
/// </summary>
NauticalMiles
}
[Serializable]
internal class Distance_Assistant
{
/// <summary>
/// Returns new geodetic coordinate in radians
/// </summary>
/// <param name="glat1">Latitude in Radians</param>
/// <param name="glon1">Longitude in Radians</param>
/// <param name="faz">Bearing</param>
/// <param name="s">Distance</param>
/// <param name="ellipse">Earth Ellipse Values</param>
/// <returns>double[]</returns>
public static double[] Direct_Ell(double glat1, double glon1, double faz, double s, double[] ellipse)
{
double EPS = 0.00000000005;//Used to determine if starting at pole.
double r, tu, sf, cf, b, cu, su, sa, c2a, x, c, d, y, sy = 0, cy = 0, cz = 0, e = 0;
double glat2, glon2, f;
//Determine if near pole
if ((Math.Abs(Math.Cos(glat1)) < EPS) && !(Math.Abs(Math.Sin(faz)) < EPS))
{
Debug.WriteLine("Warning: Location is at earth's pole. Only N-S courses are meaningful at this location.");
}
double a = ellipse[0];//Equitorial Radius
f = 1 / ellipse[1];//Flattening
r = 1 - f;
tu = r * Math.Tan(glat1);
sf = Math.Sin(faz);
cf = Math.Cos(faz);
if (cf == 0)
{
b = 0.0;
}
else
{
b = 2.0 * Math.Atan2(tu, cf);
}
cu = 1.0 / Math.Sqrt(1 + tu * tu);
su = tu * cu;
sa = cu * sf;
c2a = 1 - sa * sa;
x = 1.0 + Math.Sqrt(1.0 + c2a * (1.0 / (r * r) - 1.0));
x = (x - 2.0) / x;
c = 1.0 - x;
c = (x * x / 4.0 + 1.0) / c;
d = (0.375 * x * x - 1.0) * x;
tu = s / (r * a * c);
y = tu;
c = y + 1;
while (Math.Abs(y - c) > EPS)
{
sy = Math.Sin(y);
cy = Math.Cos(y);
cz = Math.Cos(b + y);
e = 2.0 * cz * cz - 1.0;
c = y;
x = e * cy;
y = e + e - 1.0;
y = (((sy * sy * 4.0 - 3.0) * y * cz * d / 6.0 + x) *
d / 4.0 - cz) * sy * d + tu;
}
b = cu * cy * cf - su * sy;
c = r * Math.Sqrt(sa * sa + b * b);
d = su * cy + cu * sy * cf;
glat2 = ModM.ModLat(Math.Atan2(d, c));
c = cu * cy - su * sy * cf;
x = Math.Atan2(sy * sf, c);
c = ((-3.0 * c2a + 4.0) * f + 4.0) * c2a * f / 16.0;
d = ((e * cy * c + cz) * sy * c + y) * sa;
glon2 = ModM.ModLon(glon1 + x - (1.0 - c) * d * f); //Adjust for IDL
//baz = ModM.ModCrs(Math.Atan2(sa, b) + Math.PI);
return new double[] { glat2, glon2 };
}
/// <summary>
/// Returns new geodetic coordinate in radians
/// </summary>
/// <param name="lat1">Latitude in radians</param>
/// <param name="lon1">Longitude in radians</param>
/// <param name="crs12">Bearing</param>
/// <param name="d12">Distance</param>
/// <returns>double[]</returns>
public static double[] Direct(double lat1, double lon1, double crs12, double d12)
{
var EPS = 0.00000000005;//Used to determine if near pole.
double dlon, lat, lon;
d12 = d12 * 0.0005399565; //convert meter to nm
d12 = d12 / (180 * 60 / Math.PI);//Convert to Radian
//Determine if near pole
if ((Math.Abs(Math.Cos(lat1)) < EPS) && !(Math.Abs(Math.Sin(crs12)) < EPS))
{
Debug.WriteLine("Warning: Location is at earth's pole. Only N-S courses are meaningful at this location.");
}
lat = Math.Asin(Math.Sin(lat1) * Math.Cos(d12) +
Math.Cos(lat1) * Math.Sin(d12) * Math.Cos(crs12));
if (Math.Abs(Math.Cos(lat)) < EPS)
{
lon = 0.0; //endpoint a pole
}
else
{
dlon = Math.Atan2(Math.Sin(crs12) * Math.Sin(d12) * Math.Cos(lat1),
Math.Cos(d12) - Math.Sin(lat1) * Math.Sin(lat));
lon = ModM.Mod(lon1 - dlon + Math.PI, 2 * Math.PI) - Math.PI;
}
return new double[] { lat, lon };
}
public static double[] Dist_Ell(double glat1, double glon1, double glat2, double glon2, double[] ellipse)
{
double a = ellipse[0]; //Equitorial Radius
double f = 1 / ellipse[1]; //Flattening
double r, tu1, tu2, cu1, su1, cu2, s1, b1, f1;
double x = 0, sx = 0, cx = 0, sy = 0, cy = 0, y = 0, sa = 0, c2a = 0, cz = 0, e = 0, c = 0, d = 0;
double EPS = 0.00000000005;
double faz, baz, s;
double iter = 1;
double MAXITER = 100;
if ((glat1 + glat2 == 0.0) && (Math.Abs(glon1 - glon2) == Math.PI))
{
Debug.WriteLine("Warning: Course and distance between antipodal points is undefined");
glat1 = glat1 + 0.00001; // allow algorithm to complete
}
if (glat1 == glat2 && (glon1 == glon2 || Math.Abs(Math.Abs(glon1 - glon2) - 2 * Math.PI) < EPS))
{
Debug.WriteLine("Warning: Points 1 and 2 are identical- course undefined");
//D
//crs12
//crs21
return new double[] { 0, 0, Math.PI };
}
r = 1 - f;
tu1 = r * Math.Tan(glat1);
tu2 = r * Math.Tan(glat2);
cu1 = 1.0 / Math.Sqrt(1.0 + tu1 * tu1);
su1 = cu1 * tu1;
cu2 = 1.0 / Math.Sqrt(1.0 + tu2 * tu2);
s1 = cu1 * cu2;
b1 = s1 * tu2;
f1 = b1 * tu1;
x = glon2 - glon1;
d = x + 1; // force one pass
while ((Math.Abs(d - x) > EPS) && (iter < MAXITER))
{
iter = iter + 1;
sx = Math.Sin(x);
cx = Math.Cos(x);
tu1 = cu2 * sx;
tu2 = b1 - su1 * cu2 * cx;
sy = Math.Sqrt(tu1 * tu1 + tu2 * tu2);
cy = s1 * cx + f1;
y = Math.Atan2(sy, cy);
sa = s1 * sx / sy;
c2a = 1 - sa * sa;
cz = f1 + f1;
if (c2a > 0.0)
{
cz = cy - cz / c2a;
}
e = cz * cz * 2.0 - 1.0;
c = ((-3.0 * c2a + 4.0) * f + 4.0) * c2a * f / 16.0;
d = x;
x = ((e * cy * c + cz) * sy * c + y) * sa;
x = (1.0 - c) * x * f + glon2 - glon1;
}
faz = ModM.ModCrs(Math.Atan2(tu1, tu2));
baz = ModM.ModCrs(Math.Atan2(cu1 * sx, b1 * cx - su1 * cu2) + Math.PI);
x = Math.Sqrt((1 / (r * r) - 1) * c2a + 1);
x += 1;
x = (x - 2.0) / x;
c = 1.0 - x;
c = (x * x / 4.0 + 1.0) / c;
d = (0.375 * x * x - 1.0) * x;
x = e * cy;
s = ((((sy * sy * 4.0 - 3.0) * (1.0 - e - e) * cz * d / 6.0 - x) * d / 4.0 + cz) * sy * d + y) * c * a * r;
if (Math.Abs(iter - MAXITER) < EPS)
{
Debug.WriteLine("Warning: Distance algorithm did not converge");
}
return new double[] { s, faz, baz };
}
}
}

14610
Lora-Map/Helper/Eclipse/LunarData.cs Normal file
View File

File diff suppressed because it is too large Load Diff

21578
Lora-Map/Helper/Eclipse/SolarData.cs Normal file
View File

File diff suppressed because it is too large Load Diff

161
Lora-Map/Helper/GeoFence.cs Normal file
View File

@ -0,0 +1,161 @@
using System.Collections.Generic;
using System.Linq;
namespace CoordinateSharp
{
/// <summary>
/// Geo Fence class. It helps to check if points/coordinates are inside a polygon,
/// Next to a polyline, and counting...
/// </summary>
public class GeoFence
{
#region Fields
private List<Point> _points = new List<Point>();
#endregion
/// <summary>
/// Prepare GeoFence with a list of points
/// </summary>
/// <param name="points">List of points</param>
public GeoFence(List<Point> points)
{
_points = points;
}
/// <summary>
/// Prepare Geofence with a list of coordinates
/// </summary>
/// <param name="coordinates">List of coordinates</param>
public GeoFence(List<Coordinate> coordinates)
{
foreach (var c in coordinates)
{
_points.Add(new Point { Latitude = c.Latitude.ToDouble(), Longitude = c.Longitude.ToDouble() });
}
}
#region Utils
private Coordinate ClosestPointOnSegment(Point a, Point b, Coordinate p)
{
var d = new Point
{
Longitude = b.Longitude - a.Longitude,
Latitude = b.Latitude - a.Latitude,
};
double number = (p.Longitude.ToDouble() - a.Longitude) * d.Longitude + (p.Latitude.ToDouble() - a.Latitude) * d.Latitude;
if (number <= 0.0)
return new Coordinate(a.Latitude, a.Longitude);
double denom = d.Longitude * d.Longitude + d.Latitude * d.Latitude;
if (number >= denom)
return new Coordinate(b.Latitude, b.Longitude);
return new Coordinate(a.Latitude + (number / denom) * d.Latitude, a.Longitude + (number / denom) * d.Longitude);
}
#endregion
/// <summary>
/// The function will return true if the point x,y is inside the polygon, or
/// false if it is not. If the point is exactly on the edge of the polygon,
/// then the function may return true or false.
/// </summary>
/// <param name="point">The point to test</param>
/// <returns>bool</returns>
public bool IsPointInPolygon(Coordinate point)
{
if (point == null)
return false;
double latitude = point.Latitude.ToDouble();
double longitude = point.Longitude.ToDouble();
int sides = _points.Count;
int j = sides - 1;
bool pointStatus = false;
for (int i = 0; i < sides; i++)
{
if (_points[i].Latitude < latitude && _points[j].Latitude >= latitude || _points[j].Latitude < latitude && _points[i].Latitude >= latitude)
{
if (_points[i].Longitude + (latitude - _points[i].Latitude) / (_points[j].Latitude - _points[i].Latitude) * (_points[j].Longitude - _points[i].Longitude) < longitude)
{
pointStatus = !pointStatus;
}
}
j = i;
}
return pointStatus;
}
/// <summary>
/// The function will return true if the point x,y is next the given range of
/// the polyline, or false if it is not.
/// </summary>
/// <param name="point">The point to test</param>
/// <param name="range">The range in meters</param>
/// <returns>bool</returns>
public bool IsPointInRangeOfLine(Coordinate point, double range)
{
if (point == null)
return false;
for (int i = 0; i < _points.Count - 1; i++)
{
Coordinate c = ClosestPointOnSegment(_points[i], _points[i + 1], point);
if (c.Get_Distance_From_Coordinate(point).Meters <= range)
return true;
}
return false;
}
/// <summary>
/// The function will return true if the point x,y is next the given range of
/// the polyline, or false if it is not.
/// </summary>
/// <param name="point">The point to test</param>
/// <param name="range">The range is a distance object</param>
/// <returns>bool</returns>
public bool IsPointInRangeOfLine(Coordinate point, Distance range)
{
if (point == null || range == null)
return false;
return IsPointInRangeOfLine(point, range.Meters);
}
/// <summary>
/// This class is a help class to simplify GeoFence calculus
/// </summary>
public class Point
{
/// <summary>
/// Initialize empty point
/// </summary>
public Point()
{
}
/// <summary>
/// Initialize point with defined Latitude and Longitude
/// </summary>
/// <param name="lat">Latitude (signed)</param>
/// <param name="lng">Longitude (signed)</param>
public Point(double lat, double lng)
{
Latitude = lat;
Longitude = lng;
}
/// <summary>
/// The longitude in degrees
/// </summary>
public double Longitude;
/// <summary>
/// The latitude in degrees
/// </summary>
public double Latitude;
}
}
}

21
Lora-Map/Lora-Map.csproj
View File

@ -45,10 +45,31 @@
<Reference Include="System.Xml" /> <Reference Include="System.Xml" />
</ItemGroup> </ItemGroup>
<ItemGroup> <ItemGroup>
<Compile Include="Helper\Celestial.Assistant.cs" />
<Compile Include="Helper\Celestial.cs" />
<Compile Include="Helper\Celestial.LunarEclipseCalc.cs" />
<Compile Include="Helper\Celestial.MeeusTables.cs" />
<Compile Include="Helper\Celestial.MoonCalculations.cs" />
<Compile Include="Helper\Celestial.SolarEclipseCalc.cs" />
<Compile Include="Helper\Celestial.SunCalculations.cs" />
<Compile Include="Helper\Coordinate.Assistant.cs" />
<Compile Include="Helper\Coordinate.Cartesian.cs" />
<Compile Include="Helper\Coordinate.cs" />
<Compile Include="Helper\Coordinate.EagerLoad.cs" />
<Compile Include="Helper\Coordinate.ECEF.cs" />
<Compile Include="Helper\Coordinate.Formatting.cs" />
<Compile Include="Helper\Coordinate.MGRS.cs" />
<Compile Include="Helper\Coordinate.Parser.cs" />
<Compile Include="Helper\Coordinate.UTM.cs" />
<Compile Include="Helper\Distance.cs" />
<Compile Include="Helper\Eclipse\LunarData.cs" />
<Compile Include="Helper\Eclipse\SolarData.cs" />
<Compile Include="Helper\GeoFence.cs" />
<Compile Include="Model\Admin\AdminModel.cs" /> <Compile Include="Model\Admin\AdminModel.cs" />
<Compile Include="Model\Admin\AdminSession.cs" /> <Compile Include="Model\Admin\AdminSession.cs" />
<Compile Include="Model\Marker.cs" /> <Compile Include="Model\Marker.cs" />
<Compile Include="Model\AlarmItem.cs" /> <Compile Include="Model\AlarmItem.cs" />
<Compile Include="Model\UTMData.cs" />
<Compile Include="Server.cs" /> <Compile Include="Server.cs" />
<Compile Include="Model\PositionItem.cs" /> <Compile Include="Model\PositionItem.cs" />
<Compile Include="Program.cs" /> <Compile Include="Program.cs" />

2
Lora-Map/Model/AlarmItem.cs
View File

@ -10,6 +10,7 @@ namespace Fraunhofer.Fit.IoT.LoraMap.Model {
public DateTime Recievedtime { get; private set; } public DateTime Recievedtime { get; private set; }
public Double Latitude { get; private set; } public Double Latitude { get; private set; }
public Double Longitude { get; private set; } public Double Longitude { get; private set; }
public UTMData UTM { get; private set; }
public Double Hdop { get; private set; } public Double Hdop { get; private set; }
public Boolean Fix { get; private set; } public Boolean Fix { get; private set; }
public Double Height { get; private set; } public Double Height { get; private set; }
@ -30,6 +31,7 @@ namespace Fraunhofer.Fit.IoT.LoraMap.Model {
this.Latitude = (Double)json["Gps"]["LastLatitude"]; this.Latitude = (Double)json["Gps"]["LastLatitude"];
this.Longitude = (Double)json["Gps"]["LastLongitude"]; this.Longitude = (Double)json["Gps"]["LastLongitude"];
} }
this.UTM = new UTMData(this.Latitude, this.Longitude);
this.Hdop = (Double)json["Gps"]["Hdop"]; this.Hdop = (Double)json["Gps"]["Hdop"];
this.Height = (Double)json["Gps"]["Height"]; this.Height = (Double)json["Gps"]["Height"];
} }

4
Lora-Map/Model/PositionItem.cs
View File

@ -10,6 +10,7 @@ namespace Fraunhofer.Fit.IoT.LoraMap.Model {
public DateTime Recievedtime { get; private set; } public DateTime Recievedtime { get; private set; }
public Double Latitude { get; private set; } public Double Latitude { get; private set; }
public Double Longitude { get; private set; } public Double Longitude { get; private set; }
public UTMData UTM { get; private set; }
public Double Hdop { get; private set; } public Double Hdop { get; private set; }
public DateTime Lastgpspostime { get; private set; } public DateTime Lastgpspostime { get; private set; }
public Double Battery { get; private set; } public Double Battery { get; private set; }
@ -19,6 +20,8 @@ namespace Fraunhofer.Fit.IoT.LoraMap.Model {
public String Name { get; private set; } public String Name { get; private set; }
public String Icon { get; private set; } public String Icon { get; private set; }
public PositionItem(JsonData json, JsonData marker) { public PositionItem(JsonData json, JsonData marker) {
this.Update(json); this.Update(json);
this.UpdateMarker(marker, GetId(json)); this.UpdateMarker(marker, GetId(json));
@ -86,6 +89,7 @@ namespace Fraunhofer.Fit.IoT.LoraMap.Model {
this.Latitude = (Double)json["Gps"]["LastLatitude"]; this.Latitude = (Double)json["Gps"]["LastLatitude"];
this.Longitude = (Double)json["Gps"]["LastLongitude"]; this.Longitude = (Double)json["Gps"]["LastLongitude"];
} }
this.UTM = new UTMData(this.Latitude, this.Longitude);
if(DateTime.TryParse((String)json["Gps"]["LastGPSPos"], DateTimeFormatInfo.InvariantInfo, DateTimeStyles.AssumeUniversal, out DateTime lastgpstime)) { if(DateTime.TryParse((String)json["Gps"]["LastGPSPos"], DateTimeFormatInfo.InvariantInfo, DateTimeStyles.AssumeUniversal, out DateTime lastgpstime)) {
this.Lastgpspostime = lastgpstime.ToUniversalTime(); this.Lastgpspostime = lastgpstime.ToUniversalTime();
} }

22
Lora-Map/Model/UTMData.cs Normal file
View File

@ -0,0 +1,22 @@
using System;
using System.Text.RegularExpressions;
using CoordinateSharp;
namespace Fraunhofer.Fit.IoT.LoraMap.Model {
public struct UTMData {
public String MGRS;
public Int32 FieldWidth;
public Int32 FieldHeight;
public Int32 Width;
public Int32 Height;
public UTMData(Double latitude, Double longitude) {
this.MGRS = new Coordinate(latitude, longitude).MGRS.ToString();
String[] d = Regex.Split(this.MGRS, "[0-9]+[A-Z] [A-Z]+ ([0-9]{3})([0-9]{2}) ([0-9]{3})([0-9]{2})");
this.FieldWidth = Int32.Parse(d[1]);
this.Width = Int32.Parse(d[2]);
this.FieldHeight = Int32.Parse(d[3]);
this.Height = Int32.Parse(d[4]);
}
}
}

5
Lora-Map/resources/css/global.css
View File

@ -124,9 +124,10 @@ object {
} }
#pannels #pannels_info .coord { #pannels #pannels_info .coord {
font-family: monospace; font-family: monospace;
font-size: 14px; font-size: 22px;
font-weight: bold;
} }
#pannels #pannels_info .hdop { #pannels #pannels_info .lastgps {
margin-bottom: 10px; margin-bottom: 10px;
} }
#pannels #pannels_info .update { #pannels #pannels_info .update {

7
Lora-Map/resources/js/menu.js
View File

@ -50,10 +50,13 @@ function update_pannels_info() {
} else { } else {
html += "<div class=\"gps\" style=\"color: red;\">kein GPS-Empfang</div>"; html += "<div class=\"gps\" style=\"color: red;\">kein GPS-Empfang</div>";
} }
html += "<div class=\"coord\">" + positionItem["Latitude"].toFixed(7) + ", " + positionItem["Longitude"].toFixed(7) + "</div>"; html += "<div class=\"coord\">" + positionItem["UTM"]["MGRS"] + "</div>";
html += "<div class=\"lastgps\"><span class=\"bold\">Letzter Wert:</span> Vor: " + timeCalculation(positionItem["Lastgpspostime"], "difftext") + "</div>"; html += "<div class=\"planquad\"><span class=\"bold\">Planquadrat:</span> " + positionItem["UTM"]["FieldWidth"] + ", " + positionItem["UTM"]["FieldHeight"] + "</div>";
html += "<div class=\"section\"><span class=\"bold\">Ausschnitt:</span> " + positionItem["UTM"]["Width"] + ", " + positionItem["UTM"]["Height"]+"</div>";
html += "<div class=\"height\"><span class=\"bold\">Höhe:</span> " + positionItem["Height"].toFixed(1) + " m</div>"; html += "<div class=\"height\"><span class=\"bold\">Höhe:</span> " + positionItem["Height"].toFixed(1) + " m</div>";
html += "<div class=\"hdop\"><span class=\"bold\">HDOP:</span> " + positionItem["Hdop"].toFixed(1) + "</div>"; html += "<div class=\"hdop\"><span class=\"bold\">HDOP:</span> " + positionItem["Hdop"].toFixed(1) + "</div>";
html += "<div class=\"lanlot\"><span class=\"bold\">Dezimal:</span> " + positionItem["Latitude"].toFixed(5) + ", " + positionItem["Longitude"].toFixed(5) + "</div>";
html += "<div class=\"lastgps\"><span class=\"bold\">Letzter Wert:</span> Vor: " + timeCalculation(positionItem["Lastgpspostime"], "difftext") + "</div>";
html += "<div class=\"update\"><span class=\"bold\">Update:</span> " + timeCalculation(positionItem["Recievedtime"], "str") + "<br><span class=\"bold\">Vor:</span> " + timeCalculation(positionItem["Recievedtime"], "difftext") + "</div>"; html += "<div class=\"update\"><span class=\"bold\">Update:</span> " + timeCalculation(positionItem["Recievedtime"], "str") + "<br><span class=\"bold\">Vor:</span> " + timeCalculation(positionItem["Recievedtime"], "difftext") + "</div>";
html += "<div><span class=\"bold\">RSSI:</span> " + positionItem["Rssi"] + ", <span class=\"bold\">SNR:</span> " + positionItem["Snr"] + "</div>"; html += "<div><span class=\"bold\">RSSI:</span> " + positionItem["Rssi"] + ", <span class=\"bold\">SNR:</span> " + positionItem["Snr"] + "</div>";
document.getElementById("pannels_info").innerHTML = html; document.getElementById("pannels_info").innerHTML = html;