/[zanavi_public1]/navit/navit/android/src/com/luckycatlabs/sunrisesunset/calculator/SolarEventCalculator.java
ZANavi

Diff of /navit/navit/android/src/com/luckycatlabs/sunrisesunset/calculator/SolarEventCalculator.java

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

Revision 30 Revision 31
27/** 27/**
28 * Parent class of the Sunrise and Sunset calculator classes. 28 * Parent class of the Sunrise and Sunset calculator classes.
29 */ 29 */
30public class SolarEventCalculator 30public class SolarEventCalculator
31{ 31{
32 private Location2 location; 32 private Location2 location;
33 private TimeZone timeZone; 33 private TimeZone timeZone;
34 34
35 /** 35 /**
36 * Constructs a new <code>SolarEventCalculator</code> using the given parameters. 36 * Constructs a new <code>SolarEventCalculator</code> using the given parameters.
37 * 37 *
38 * @param location 38 * @param location
39 * <code>Location</code> of the place where the solar event should be calculated from. 39 * <code>Location</code> of the place where the solar event should be calculated from.
40 * @param timeZoneIdentifier 40 * @param timeZoneIdentifier
41 * time zone identifier of the timezone of the location parameter. For example, 41 * time zone identifier of the timezone of the location parameter. For example,
42 * "America/New_York". 42 * "America/New_York".
43 */ 43 */
44 public SolarEventCalculator(Location2 location, String timeZoneIdentifier) 44 public SolarEventCalculator(Location2 location, String timeZoneIdentifier)
45 { 45 {
46 this.location = location; 46 this.location = location;
47 this.timeZone = TimeZone.getTimeZone(timeZoneIdentifier); 47 this.timeZone = TimeZone.getTimeZone(timeZoneIdentifier);
86 x = Math.abs(El_geometric + 0.589); 86 x = Math.abs(El_geometric + 0.589);
87 double refraction = Math.abs(a0 + a1 * x + a2 * x * x + a3 * x * x * x + a4 * x * x * x * x); 87 double refraction = Math.abs(a0 + a1 * x + a2 * x * x + a3 * x * x * x + a4 * x * x * x * x);
88 88
89 if (El_geometric > 10.2) 89 if (El_geometric > 10.2)
90 { 90 {
91 El_observed = El_geometric 91 El_observed = El_geometric + 0.01617 * (Math.cos(Radians(Math.abs(El_geometric))) / Math.sin(Radians(Math.abs(El_geometric))));
92 + 0.01617
93 * (Math.cos(Radians(Math.abs(El_geometric))) / Math.sin(Radians(Math
94 .abs(El_geometric))));
95 } 92 }
96 else 93 else
97 { 94 {
98 El_observed = El_geometric + refraction; 95 El_observed = El_geometric + refraction;
99 96
140 return (lmst); 137 return (lmst);
141 } 138 }
142 139
143 public class cart_ret 140 public class cart_ret
144 { 141 {
145 double x; 142 double x;
146 double y; 143 double y;
147 double z; 144 double z;
148 double radius; 145 double radius;
149 double lat; 146 double lat;
150 double lon; 147 double lon;
151 } 148 }
152 149
153 public cart_ret EquPolar2Cart(double lon, double lat, double distance) 150 public cart_ret EquPolar2Cart(double lon, double lat, double distance)
154 { 151 {
155 cart_ret cart = new cart_ret(); 152 cart_ret cart = new cart_ret();
194 return (cart); 191 return (cart);
195 } 192 }
196 193
197 public class moonCoor_ret 194 public class moonCoor_ret
198 { 195 {
199 public double lon; 196 public double lon;
200 public double lat; 197 public double lat;
201 public double orbitLon; 198 public double orbitLon;
202 public double distance; 199 public double distance;
203 public double diameter; 200 public double diameter;
204 public double parallax; 201 public double parallax;
205 public double raGeocentric; 202 public double raGeocentric;
206 public double decGeocentric; 203 public double decGeocentric;
207 public double ra; 204 public double ra;
208 public double dec; 205 public double dec;
209 public double raTopocentric; 206 public double raTopocentric;
210 public double decTopocentric; 207 public double decTopocentric;
211 public double distanceTopocentric; 208 public double distanceTopocentric;
212 public double moonAge; 209 public double moonAge;
213 public double phase; 210 public double phase;
214 public double az; 211 public double az;
215 public double alt; 212 public double alt;
216 public String moonPhase; 213 public String moonPhase;
217 public String sign; 214 public String sign;
218 } 215 }
219 216
220 public double Mod2Pi(double x) 217 public double Mod2Pi(double x)
221 { 218 {
222 x = Mod(x, 2. * Math.PI); 219 x = Mod(x, 2. * Math.PI);
227 { 224 {
228 double pi = Math.PI; 225 double pi = Math.PI;
229 double DEG = pi / 180.0; 226 double DEG = pi / 180.0;
230 227
231 double T = (TDT - 2451545.0) / 36525.; // Epoch 2000 January 1.5 228 double T = (TDT - 2451545.0) / 36525.; // Epoch 2000 January 1.5
232 double eps = (23. + (26 + 21.45 / 60.) / 60. + T * (-46.815 + T * (-0.0006 + T * 0.00181)) 229 double eps = (23. + (26 + 21.45 / 60.) / 60. + T * (-46.815 + T * (-0.0006 + T * 0.00181)) / 3600.) * DEG;
233 / 3600.)
234 * DEG;
235 double coseps = Math.cos(eps); 230 double coseps = Math.cos(eps);
236 double sineps = Math.sin(eps); 231 double sineps = Math.sin(eps);
237 232
238 double sinlon = Math.sin(coor.lon); 233 double sinlon = Math.sin(coor.lon);
239 coor.ra = Mod2Pi(Math.atan2((sinlon * coseps - Math.tan(coor.lat) * sineps), Math 234 coor.ra = Mod2Pi(Math.atan2((sinlon * coseps - Math.tan(coor.lat) * sineps), Math.cos(coor.lon)));
240 .cos(coor.lon)));
241 coor.dec = Math.asin(Math.sin(coor.lat) * coseps + Math.cos(coor.lat) * sineps * sinlon); 235 coor.dec = Math.asin(Math.sin(coor.lat) * coseps + Math.cos(coor.lat) * sineps * sinlon);
242 236
243 return coor; 237 return coor;
244 } 238 }
245 239
282 return coor; 276 return coor;
283 } 277 }
284 278
285 public class sunCoor_ret 279 public class sunCoor_ret
286 { 280 {
287 double lon; 281 double lon;
288 double anomalyMean; 282 double anomalyMean;
289 } 283 }
290 284
291 public moonCoor_ret MoonPosition(sunCoor_ret sunCoor, double TDT, cart_ret observer, double lmst) 285 public moonCoor_ret MoonPosition(sunCoor_ret sunCoor, double TDT, cart_ret observer, double lmst)
292 { 286 {
293 double D = TDT - 2447891.5; 287 double D = TDT - 2447891.5;
412 double DEG = pi / 180.0; 406 double DEG = pi / 180.0;
413 double RAD = 180.0 / pi; 407 double RAD = 180.0 / pi;
414 408
415 date.setTimeZone(TimeZone.getTimeZone("UTC")); 409 date.setTimeZone(TimeZone.getTimeZone("UTC"));
416 410
417 double JD0 = CalcJD((double) date.get(Calendar.DAY_OF_MONTH), (double) date 411 double JD0 = CalcJD((double) date.get(Calendar.DAY_OF_MONTH), (double) date.get(Calendar.MONTH) + 1, (double) date.get(Calendar.YEAR));
418 .get(Calendar.MONTH) + 1, (double) date.get(Calendar.YEAR));
419 double JD = JD0
420 + ((double) date.get(Calendar.HOUR_OF_DAY) + (double) date.get(Calendar.MINUTE) / 60 + (double) date 412 double JD = JD0 + ((double) date.get(Calendar.HOUR_OF_DAY) + (double) date.get(Calendar.MINUTE) / 60 + (double) date.get(Calendar.SECOND) / 3600) / 24;
421 .get(Calendar.SECOND) / 3600) / 24;
422 // UTC 413 // UTC
423 double TDT = JD; 414 double TDT = JD;
424 //System.out.println("TDT=" + TDT); 415 //System.out.println("TDT=" + TDT);
425 416
426 double lat = this.location.getLatitude().doubleValue() * DEG; // geodetic latitude of observer on WGS84 417 double lat = this.location.getLatitude().doubleValue() * DEG; // geodetic latitude of observer on WGS84
454 double Month = date.get(Calendar.MONTH) + 1; // month starts from "0" zero !! 445 double Month = date.get(Calendar.MONTH) + 1; // month starts from "0" zero !!
455 double Day = date.get(Calendar.DAY_OF_MONTH); 446 double Day = date.get(Calendar.DAY_OF_MONTH);
456 double Hour = date.get(Calendar.HOUR_OF_DAY); 447 double Hour = date.get(Calendar.HOUR_OF_DAY);
457 double Minute = date.get(Calendar.MINUTE); 448 double Minute = date.get(Calendar.MINUTE);
458 double Second = date.get(Calendar.SECOND); 449 double Second = date.get(Calendar.SECOND);
459 double d = 367 * Year - div((7 * (Year + (div((Month + 9), 12)))), 4) + div((275 * Month), 9) 450 double d = 367 * Year - div((7 * (Year + (div((Month + 9), 12)))), 4) + div((275 * Month), 9) + Day - 730530;
460 + Day - 730530;
461 //System.out.println("dd=" + d); 451 //System.out.println("dd=" + d);
462 d = d + Hour / 24 + Minute / (60 * 24) + Second / (24 * 60 * 60); // OK 452 d = d + Hour / 24 + Minute / (60 * 24) + Second / (24 * 60 * 60); // OK
463 453
464 // System.out.println("y=" + Year); 454 // System.out.println("y=" + Year);
465 // System.out.println("m=" + Month); 455 // System.out.println("m=" + Month);
484 474
485 w = Rev(w); 475 w = Rev(w);
486 M = Rev(M); 476 M = Rev(M);
487 N = Rev(N); 477 N = Rev(N);
488 478
489 double E = M + (180 / Math.PI) * e * Math.sin(Math.toRadians(M)) 479 double E = M + (180 / Math.PI) * e * Math.sin(Math.toRadians(M)) * (1 + e * Math.cos(Math.toRadians(M)));
490 * (1 + e * Math.cos(Math.toRadians(M)));
491 E = Rev(E); // OK 480 E = Rev(E); // OK
492 481
493 double Ebeforeit = E; 482 double Ebeforeit = E;
494 // now iterate until difference between E0 and E1 is less than 0.005_deg 483 // now iterate until difference between E0 and E1 is less than 0.005_deg
495 // use E0, calculate E1 484 // use E0, calculate E1
503 492
504 while ((E_error > 0.0005) && (Iterations < 20)) // ok - itererer korrekt 493 while ((E_error > 0.0005) && (Iterations < 20)) // ok - itererer korrekt
505 { 494 {
506 Iterations = Iterations + 1; 495 Iterations = Iterations + 1;
507 E0 = E; 496 E0 = E;
508 E1 = E0 - (E0 - (180 / Math.PI) * e * Math.sin(Math.toRadians(E0)) - M) 497 E1 = E0 - (E0 - (180 / Math.PI) * e * Math.sin(Math.toRadians(E0)) - M) / (1 - e * Math.cos(Math.toRadians(E0)));
509 / (1 - e * Math.cos(Math.toRadians(E0)));
510 //alert('1 E0='+E0+'\nNew E1='+E1+'\nE='+E+'\Diff='+Rev(E0-E1)); 498 //alert('1 E0='+E0+'\nNew E1='+E1+'\nE='+E+'\Diff='+Rev(E0-E1));
511 E = Rev(E1); 499 E = Rev(E1);
512 //alert(Math.abs(E-E0)); 500 //alert(Math.abs(E-E0));
513 501
514 //Eafterit = E; 502 //Eafterit = E;
547 r = Math.sqrt(x * x + y * y); 535 r = Math.sqrt(x * x + y * y);
548 v = Math.toDegrees(Math.atan2(y, x)); 536 v = Math.toDegrees(Math.atan2(y, x));
549 537
550 //alert('E='+E); 538 //alert('E='+E);
551 539
552 // ok s langt 540 // ok s langt
553 541
554 double sunlon = Rev(v + w); // trolig ok 542 double sunlon = Rev(v + w); // trolig ok
555 543
556 x = r * Math.cos(Math.toRadians(sunlon)); 544 x = r * Math.cos(Math.toRadians(sunlon));
557 y = r * Math.sin(Math.toRadians(sunlon)); 545 y = r * Math.sin(Math.toRadians(sunlon));
558 double z = 0; 546 double z = 0;
559 547
560 double xeclip = r 548 double xeclip = r * (Math.cos(Math.toRadians(N)) * Math.cos(Math.toRadians(v + w)) - Math.sin(Math.toRadians(N)) * Math.sin(Math.toRadians(v + w)) * Math.cos(Math.toRadians(i)));
561 * (Math.cos(Math.toRadians(N)) * Math.cos(Math.toRadians(v + w)) - Math.sin(Math 549 double yeclip = r * (Math.sin(Math.toRadians(N)) * Math.cos(Math.toRadians(v + w)) + Math.cos(Math.toRadians(N)) * Math.sin(Math.toRadians(v + w)) * Math.cos(Math.toRadians(i)));
562 .toRadians(N))
563 * Math.sin(Math.toRadians(v + w)) * Math.cos(Math.toRadians(i)));
564 double yeclip = r
565 * (Math.sin(Math.toRadians(N)) * Math.cos(Math.toRadians(v + w)) + Math.cos(Math
566 .toRadians(N))
567 * Math.sin(Math.toRadians(v + w)) * Math.cos(Math.toRadians(i)));
568 double zeclip = r * Math.sin(Math.toRadians(v + w)) * Math.sin(Math.toRadians(i)); 550 double zeclip = r * Math.sin(Math.toRadians(v + w)) * Math.sin(Math.toRadians(i));
569 551
570 double moon_longitude = Rev(Math.toDegrees(Math.atan2(yeclip, xeclip))); // OK 552 double moon_longitude = Rev(Math.toDegrees(Math.atan2(yeclip, xeclip))); // OK
571 double moon_latitude = Math.toDegrees(Math.atan2(zeclip, Math.sqrt(xeclip * xeclip + yeclip 553 double moon_latitude = Math.toDegrees(Math.atan2(zeclip, Math.sqrt(xeclip * xeclip + yeclip * yeclip))); // trolig OK
572 * yeclip))); // trolig OK
573 554
574 // ----------- SUN ----------- 555 // ----------- SUN -----------
575 // ----------- SUN ----------- 556 // ----------- SUN -----------
576 // date.setTimeZone(this.timeZone); 557 // date.setTimeZone(this.timeZone);
577 BigDecimal longitudeHour = getLongitudeHour(date, true); 558 BigDecimal longitudeHour = getLongitudeHour(date, true);
578 BigDecimal meanAnomaly = getMeanAnomaly(longitudeHour); 559 BigDecimal meanAnomaly = getMeanAnomaly(longitudeHour);
579 //BigDecimal sunTrueLong = getSunTrueLongitude(meanAnomaly); 560 //BigDecimal sunTrueLong = getSunTrueLongitude(meanAnomaly);
580 //BigDecimal cosineSunLocalHour = getCosineSunLocalHour(sunTrueLong, Zenith.OFFICIAL); 561 //BigDecimal cosineSunLocalHour = getCosineSunLocalHour(sunTrueLong, Zenith.OFFICIAL);
581 562
582// System.out.println("Sun MA=" + meanAnomaly); 563 // System.out.println("Sun MA=" + meanAnomaly);
583 564
584 // geschtzt!!!!! 565 // geschtzt!!!!!
585 // see -> http://www.obliquity.com/info/meaning.html 566 // see -> http://www.obliquity.com/info/meaning.html
586 double sun_Obliquity = 23.45; 567 double sun_Obliquity = 23.45;
587 // geschtzt!!!!! 568 // geschtzt!!!!!
588 569
589 // sunangles[11] ???? 570 // sunangles[11] ????
590 double w_S = 282.9404 + 4.70935E-5 * d; //OK 571 double w_S = 282.9404 + 4.70935E-5 * d; //OK
591 //double a_S = 1; 572 //double a_S = 1;
592 //a=6.6107940559473451507806351067866; 573 //a=6.6107940559473451507806351067866;
600 double GMST0_sun = (L_S + 180); 581 double GMST0_sun = (L_S + 180);
601 582
602 L_S = Rev(L_S); 583 L_S = Rev(L_S);
603 sun_Obliquity = oblecl_S; 584 sun_Obliquity = oblecl_S;
604 585
605// System.out.println("GMST0_sun=" + GMST0_sun); 586 // System.out.println("GMST0_sun=" + GMST0_sun);
606// System.out.println("L_S=" + L_S); 587 // System.out.println("L_S=" + L_S);
607// System.out.println("oblecl_S=" + oblecl_S); 588 // System.out.println("oblecl_S=" + oblecl_S);
608 // sunangles[11] ???? 589 // sunangles[11] ????
609 // ----------- SUN ----------- 590 // ----------- SUN -----------
610 // ----------- SUN ----------- 591 // ----------- SUN -----------
611
612 592
613 double Mm = Rev(M); // Moons mean anomaly 593 double Mm = Rev(M); // Moons mean anomaly
614 double Lm = Rev(N + w + M); // moon mean longitude 594 double Lm = Rev(N + w + M); // moon mean longitude
615 double Ms = meanAnomaly.doubleValue(); // sun mean anomaly 595 double Ms = meanAnomaly.doubleValue(); // sun mean anomaly
616 //double Ls = sunTrueLong.doubleValue(); // sun mean longtitude 596 //double Ls = sunTrueLong.doubleValue(); // sun mean longtitude
617 double Ls = L_S; 597 double Ls = L_S;
618 double D = Rev(Lm - Ls); //Moons mean elongation 598 double D = Rev(Lm - Ls); //Moons mean elongation
619 double F = Rev(Lm - N); //Moons argument of latitude 599 double F = Rev(Lm - N); //Moons argument of latitude
620
621 600
622 // Perbutations Moons Longitude 601 // Perbutations Moons Longitude
623 602
624 double P_lon1 = -1.274 * Math.sin(Radians(Mm - 2 * D)); // (Evection) 603 double P_lon1 = -1.274 * Math.sin(Radians(Mm - 2 * D)); // (Evection)
625 double P_lon2 = +0.658 * Math.sin(Radians(2 * D)); // (Variation) 604 double P_lon2 = +0.658 * Math.sin(Radians(2 * D)); // (Variation)
633 double P_lon10 = -0.031 * Math.sin(Radians(Mm + Ms)); 612 double P_lon10 = -0.031 * Math.sin(Radians(Mm + Ms));
634 double P_lon11 = -0.015 * Math.sin(Radians(2 * F - 2 * D)); 613 double P_lon11 = -0.015 * Math.sin(Radians(2 * F - 2 * D));
635 double P_lon12 = +0.011 * Math.sin(Radians(Mm - 4 * D)); 614 double P_lon12 = +0.011 * Math.sin(Radians(Mm - 4 * D));
636 // Perbutations Moons Latitude 615 // Perbutations Moons Latitude
637 616
638
639 double P_lat1 = -0.173 * Math.sin(Radians(F - 2 * D)); 617 double P_lat1 = -0.173 * Math.sin(Radians(F - 2 * D));
640 double P_lat2 = -0.055 * Math.sin(Radians(Mm - F - 2 * D)); 618 double P_lat2 = -0.055 * Math.sin(Radians(Mm - F - 2 * D));
641 double P_lat3 = -0.046 * Math.sin(Radians(Mm + F - 2 * D)); 619 double P_lat3 = -0.046 * Math.sin(Radians(Mm + F - 2 * D));
642 double P_lat4 = +0.033 * Math.sin(Radians(F + 2 * D)); 620 double P_lat4 = +0.033 * Math.sin(Radians(F + 2 * D));
643 double P_lat5 = +0.017 * Math.sin(Radians(2 * Mm + F)); 621 double P_lat5 = +0.017 * Math.sin(Radians(2 * Mm + F));
644 622
645 double P_lon = P_lon1 + P_lon2 + P_lon3 + P_lon4 + P_lon5 + P_lon6 + P_lon7 + P_lon8 + P_lon9 623 double P_lon = P_lon1 + P_lon2 + P_lon3 + P_lon4 + P_lon5 + P_lon6 + P_lon7 + P_lon8 + P_lon9 + P_lon10 + P_lon11 + P_lon12;
646 + P_lon10 + P_lon11 + P_lon12;
647 double P_lat = P_lat1 + P_lat2 + P_lat3 + P_lat4 + P_lat5; 624 double P_lat = P_lat1 + P_lat2 + P_lat3 + P_lat4 + P_lat5;
648 double P_moondistance = -0.58 * Math.cos(Radians(Mm - 2 * D)) - 0.46 625 double P_moondistance = -0.58 * Math.cos(Radians(Mm - 2 * D)) - 0.46 * Math.cos(Radians(2 * D));
649 * Math.cos(Radians(2 * D));
650 626
651 //alert('P_lon='+P_lon+'\nP_lat='+P_lat+'\nP_moondistance='+P_moondistance); 627 //alert('P_lon='+P_lon+'\nP_lat='+P_lat+'\nP_moondistance='+P_moondistance);
652 628
653 moon_longitude = moon_longitude + P_lon; 629 moon_longitude = moon_longitude + P_lon;
654 moon_latitude = moon_latitude + P_lat; 630 moon_latitude = moon_latitude + P_lat;
670 double Moon_Decl = Deg(Math.atan2(zh, Math.sqrt(xh * xh + yh * yh))); // trolig OK 646 double Moon_Decl = Deg(Math.atan2(zh, Math.sqrt(xh * xh + yh * yh))); // trolig OK
671 647
672 Moon_RA = Rev(Deg(Math.atan2(yequat, xequat))); // OK 648 Moon_RA = Rev(Deg(Math.atan2(yequat, xequat))); // OK
673 Moon_Decl = Deg(Math.atan2(zequat, Math.sqrt(xequat * xequat + yequat * yequat))); // trolig OK 649 Moon_Decl = Deg(Math.atan2(zequat, Math.sqrt(xequat * xequat + yequat * yequat))); // trolig OK
674 650
675// System.out.println("Moon Ra=" + Moon_RA); 651 // System.out.println("Moon Ra=" + Moon_RA);
676 652
677// System.out.println("Ls=" + Ls); 653 // System.out.println("Ls=" + Ls);
678 // war "+180" mit "-180" funkts aber besser :-) 654 // war "+180" mit "-180" funkts aber besser :-)
679 double GMST0 = (Ls - 180); 655 double GMST0 = (Ls - 180);
680// System.out.println("GMST0=" + GMST0); 656 // System.out.println("GMST0=" + GMST0);
681
682 657
683 //*********CALCULATE TIME ********************* 658 //*********CALCULATE TIME *********************
684 659
685// System.out.println("d1=" + d); 660 // System.out.println("d1=" + d);
686 double UT = d - Math.floor(d); 661 double UT = d - Math.floor(d);
687 //UT = 0.9; 662 //UT = 0.9;
688// System.out.println("d1=" + UT); 663 // System.out.println("d1=" + UT);
689 //alert("UT="+UT); 664 //alert("UT="+UT);
690 665
691 // ??????????????????? 666 // ???????????????????
692 // ??????????????????? 667 // ???????????????????
693 double SiteLon = this.location.getLatitude().doubleValue(); 668 double SiteLon = this.location.getLatitude().doubleValue();
700 // System.out.println("GMST0 + UT * 360 + SiteLon=" + GMST0 + " " + UT + " " + SiteLon); 675 // System.out.println("GMST0 + UT * 360 + SiteLon=" + GMST0 + " " + UT + " " + SiteLon);
701 // System.out.println("SIDEREALTIME - Moon_RA=" + SIDEREALTIME + " " + Moon_RA); 676 // System.out.println("SIDEREALTIME - Moon_RA=" + SIDEREALTIME + " " + Moon_RA);
702 // System.out.println("SIDEREALTIME=" + SIDEREALTIME); 677 // System.out.println("SIDEREALTIME=" + SIDEREALTIME);
703 // System.out.println("HourAngle=" + HourAngle); 678 // System.out.println("HourAngle=" + HourAngle);
704 679
705
706 // make things easier!! 680 // make things easier!!
707 double pi = Math.PI; 681 double pi = Math.PI;
708 682
709 x = Math.cos(HourAngle * Math.PI / 180) * Math.cos(Moon_Decl * Math.PI / 180); 683 x = Math.cos(HourAngle * Math.PI / 180) * Math.cos(Moon_Decl * Math.PI / 180);
710 y = Math.sin(HourAngle * Math.PI / 180) * Math.cos(Moon_Decl * Math.PI / 180); 684 y = Math.sin(HourAngle * Math.PI / 180) * Math.cos(Moon_Decl * Math.PI / 180);
713 double xhor = x * Math.sin(SiteLat * pi / 180) - z * Math.cos(SiteLat * pi / 180); 687 double xhor = x * Math.sin(SiteLat * pi / 180) - z * Math.cos(SiteLat * pi / 180);
714 //alert('sitelat='+SiteLat+'\nsitelon='+SiteLon); 688 //alert('sitelat='+SiteLat+'\nsitelon='+SiteLon);
715 double yhor = y; 689 double yhor = y;
716 double zhor = x * Math.cos(SiteLat * pi / 180) + z * Math.sin(SiteLat * pi / 180); 690 double zhor = x * Math.cos(SiteLat * pi / 180) + z * Math.sin(SiteLat * pi / 180);
717 691
718
719 double MoonElevation = Deg(Math.asin(zhor)); // ok regner ikke mne elevation helt riktig... 692 double MoonElevation = Deg(Math.asin(zhor)); // ok regner ikke mne elevation helt riktig...
720 //System.out.println("MoonElevation=" + MoonElevation); 693 //System.out.println("MoonElevation=" + MoonElevation);
721 694
722 MoonElevation = MoonElevation - Deg(Math.asin(1 / r * Math.cos(Radians(MoonElevation)))); 695 MoonElevation = MoonElevation - Deg(Math.asin(1 / r * Math.cos(Radians(MoonElevation))));
723 //System.out.println("MoonElevation=" + MoonElevation); 696 //System.out.println("MoonElevation=" + MoonElevation);
724 697
739 712
740 /** 713 /**
741 * Computes the sunrise time for the given zenith at the given date. 714 * Computes the sunrise time for the given zenith at the given date.
742 * 715 *
743 * @param solarZenith 716 * @param solarZenith
744 * <code>Zenith</code> enum corresponding to the type of sunrise to compute. 717 * <code>Zenith</code> enum corresponding to the type of sunrise to compute.
745 * @param date 718 * @param date
746 * <code>Calendar</code> object representing the date to compute the sunrise for. 719 * <code>Calendar</code> object representing the date to compute the sunrise for.
747 * @return the sunrise time, in HH:MM format (24-hour clock), 00:00 if the sun does not rise on the given 720 * @return the sunrise time, in HH:MM format (24-hour clock), 00:00 if the sun does not rise on the given
748 * date. 721 * date.
749 */ 722 */
750 public String computeSunriseTime(Zenith solarZenith, Calendar date) 723 public String computeSunriseTime(Zenith solarZenith, Calendar date)
751 { 724 {
754 727
755 /** 728 /**
756 * Computes the sunset time for the given zenith at the given date. 729 * Computes the sunset time for the given zenith at the given date.
757 * 730 *
758 * @param solarZenith 731 * @param solarZenith
759 * <code>Zenith</code> enum corresponding to the type of sunset to compute. 732 * <code>Zenith</code> enum corresponding to the type of sunset to compute.
760 * @param date 733 * @param date
761 * <code>Calendar</code> object representing the date to compute the sunset for. 734 * <code>Calendar</code> object representing the date to compute the sunset for.
762 * @return the sunset time, in HH:MM format (24-hour clock), 00:00 if the sun does not set on the given 735 * @return the sunset time, in HH:MM format (24-hour clock), 00:00 if the sun does not set on the given
763 * date. 736 * date.
764 */ 737 */
765 public String computeSunsetTime(Zenith solarZenith, Calendar date) 738 public String computeSunsetTime(Zenith solarZenith, Calendar date)
766 { 739 {
773 BigDecimal longitudeHour = getLongitudeHour(date, isSunrise); 746 BigDecimal longitudeHour = getLongitudeHour(date, isSunrise);
774 747
775 BigDecimal meanAnomaly = getMeanAnomaly(longitudeHour); 748 BigDecimal meanAnomaly = getMeanAnomaly(longitudeHour);
776 BigDecimal sunTrueLong = getSunTrueLongitude(meanAnomaly); 749 BigDecimal sunTrueLong = getSunTrueLongitude(meanAnomaly);
777 BigDecimal cosineSunLocalHour = getCosineSunLocalHour(sunTrueLong, solarZenith); 750 BigDecimal cosineSunLocalHour = getCosineSunLocalHour(sunTrueLong, solarZenith);
778 if ((cosineSunLocalHour.doubleValue() < -1.0) || (cosineSunLocalHour.doubleValue() > 1.0)) { return "99:99"; } 751 if ((cosineSunLocalHour.doubleValue() < -1.0) || (cosineSunLocalHour.doubleValue() > 1.0))
752 {
753 return "99:99";
754 }
779 755
780 BigDecimal sunLocalHour = getSunLocalHour(cosineSunLocalHour, isSunrise); 756 BigDecimal sunLocalHour = getSunLocalHour(cosineSunLocalHour, isSunrise);
781 BigDecimal localMeanTime = getLocalMeanTime(sunTrueLong, longitudeHour, sunLocalHour); 757 BigDecimal localMeanTime = getLocalMeanTime(sunTrueLong, longitudeHour, sunLocalHour);
782 BigDecimal localTime = getLocalTime(localMeanTime, date); 758 BigDecimal localTime = getLocalTime(localMeanTime, date);
783 return getLocalTimeAsString(localTime); 759 return getLocalTimeAsString(localTime);
784 } 760 }
785 761
786 /** 762 /**
787 * Computes the base longitude hour, lngHour in the algorithm. 763 * Computes the base longitude hour, lngHour in the algorithm.
788 * 764 *
789 * @return the longitude of the location of the solar event divided by 15 (deg/hour), in 765 * @return the longitude of the location of the solar event divided by 15 (deg/hour), in <code>BigDecimal</code> form.
790 * <code>BigDecimal</code> form.
791 */ 766 */
792 private BigDecimal getBaseLongitudeHour() 767 private BigDecimal getBaseLongitudeHour()
793 { 768 {
794 return divideBy(location.getLongitude(), BigDecimal.valueOf(15)); 769 return divideBy(location.getLongitude(), BigDecimal.valueOf(15));
795 } 770 }
817 * 792 *
818 * @return the suns mean anomaly, M, in <code>BigDecimal</code> form. 793 * @return the suns mean anomaly, M, in <code>BigDecimal</code> form.
819 */ 794 */
820 private BigDecimal getMeanAnomaly(BigDecimal longitudeHour) 795 private BigDecimal getMeanAnomaly(BigDecimal longitudeHour)
821 { 796 {
822 BigDecimal meanAnomaly = multiplyBy(new BigDecimal("0.9856"), longitudeHour).subtract( 797 BigDecimal meanAnomaly = multiplyBy(new BigDecimal("0.9856"), longitudeHour).subtract(new BigDecimal("3.289"));
823 new BigDecimal("3.289"));
824 return setScale(meanAnomaly); 798 return setScale(meanAnomaly);
825 } 799 }
826 800
827 /** 801 /**
828 * Computes the true longitude of the sun, L in the algorithm, at the given location, adjusted to fit in 802 * Computes the true longitude of the sun, L in the algorithm, at the given location, adjusted to fit in
829 * the range [0-360]. 803 * the range [0-360].
830 * 804 *
831 * @param meanAnomaly 805 * @param meanAnomaly
832 * the suns mean anomaly. 806 * the suns mean anomaly.
833 * @return the suns true longitude, in <code>BigDecimal</code> form. 807 * @return the suns true longitude, in <code>BigDecimal</code> form.
834 */ 808 */
835 private BigDecimal getSunTrueLongitude(BigDecimal meanAnomaly) 809 private BigDecimal getSunTrueLongitude(BigDecimal meanAnomaly)
836 { 810 {
837 BigDecimal sinMeanAnomaly = new BigDecimal(Math.sin(convertDegreesToRadians(meanAnomaly) 811 BigDecimal sinMeanAnomaly = new BigDecimal(Math.sin(convertDegreesToRadians(meanAnomaly).doubleValue()));
838 .doubleValue())); 812 BigDecimal sinDoubleMeanAnomaly = new BigDecimal(Math.sin(multiplyBy(convertDegreesToRadians(meanAnomaly), BigDecimal.valueOf(2)).doubleValue()));
839 BigDecimal sinDoubleMeanAnomaly = new BigDecimal(Math.sin(multiplyBy(
840 convertDegreesToRadians(meanAnomaly), BigDecimal.valueOf(2)).doubleValue()));
841 813
842 BigDecimal firstPart = meanAnomaly.add(multiplyBy(sinMeanAnomaly, new BigDecimal("1.916"))); 814 BigDecimal firstPart = meanAnomaly.add(multiplyBy(sinMeanAnomaly, new BigDecimal("1.916")));
843 BigDecimal secondPart = multiplyBy(sinDoubleMeanAnomaly, new BigDecimal("0.020")).add( 815 BigDecimal secondPart = multiplyBy(sinDoubleMeanAnomaly, new BigDecimal("0.020")).add(new BigDecimal("282.634"));
844 new BigDecimal("282.634"));
845 BigDecimal trueLongitude = firstPart.add(secondPart); 816 BigDecimal trueLongitude = firstPart.add(secondPart);
846 817
847 if (trueLongitude.doubleValue() > 360) 818 if (trueLongitude.doubleValue() > 360)
848 { 819 {
849 trueLongitude = trueLongitude.subtract(BigDecimal.valueOf(360)); 820 trueLongitude = trueLongitude.subtract(BigDecimal.valueOf(360));
854 /** 825 /**
855 * Computes the suns right ascension, RA in the algorithm, adjusting for the quadrant of L and turning it 826 * Computes the suns right ascension, RA in the algorithm, adjusting for the quadrant of L and turning it
856 * into degree-hours. Will be in the range [0,360]. 827 * into degree-hours. Will be in the range [0,360].
857 * 828 *
858 * @param sunTrueLong 829 * @param sunTrueLong
859 * Suns true longitude, in <code>BigDecimal</code> 830 * Suns true longitude, in <code>BigDecimal</code>
860 * @return suns right ascension in degree-hours, in <code>BigDecimal</code> form. 831 * @return suns right ascension in degree-hours, in <code>BigDecimal</code> form.
861 */ 832 */
862 private BigDecimal getRightAscension(BigDecimal sunTrueLong) 833 private BigDecimal getRightAscension(BigDecimal sunTrueLong)
863 { 834 {
864 BigDecimal tanL = new BigDecimal(Math.tan(convertDegreesToRadians(sunTrueLong).doubleValue())); 835 BigDecimal tanL = new BigDecimal(Math.tan(convertDegreesToRadians(sunTrueLong).doubleValue()));
865 836
866 BigDecimal innerParens = multiplyBy(convertRadiansToDegrees(tanL), new BigDecimal("0.91764")); 837 BigDecimal innerParens = multiplyBy(convertRadiansToDegrees(tanL), new BigDecimal("0.91764"));
867 BigDecimal rightAscension = new BigDecimal(Math.atan(convertDegreesToRadians(innerParens) 838 BigDecimal rightAscension = new BigDecimal(Math.atan(convertDegreesToRadians(innerParens).doubleValue()));
868 .doubleValue()));
869 rightAscension = setScale(convertRadiansToDegrees(rightAscension)); 839 rightAscension = setScale(convertRadiansToDegrees(rightAscension));
870 840
871 if (rightAscension.doubleValue() < 0) 841 if (rightAscension.doubleValue() < 0)
872 { 842 {
873 rightAscension = rightAscension.add(BigDecimal.valueOf(360)); 843 rightAscension = rightAscension.add(BigDecimal.valueOf(360));
893 BigDecimal sinSunDeclination = getSinOfSunDeclination(sunTrueLong); 863 BigDecimal sinSunDeclination = getSinOfSunDeclination(sunTrueLong);
894 BigDecimal cosineSunDeclination = getCosineOfSunDeclination(sinSunDeclination); 864 BigDecimal cosineSunDeclination = getCosineOfSunDeclination(sinSunDeclination);
895 865
896 BigDecimal zenithInRads = convertDegreesToRadians(zenith.degrees()); 866 BigDecimal zenithInRads = convertDegreesToRadians(zenith.degrees());
897 BigDecimal cosineZenith = BigDecimal.valueOf(Math.cos(zenithInRads.doubleValue())); 867 BigDecimal cosineZenith = BigDecimal.valueOf(Math.cos(zenithInRads.doubleValue()));
898 BigDecimal sinLatitude = BigDecimal.valueOf(Math.sin(convertDegreesToRadians( 868 BigDecimal sinLatitude = BigDecimal.valueOf(Math.sin(convertDegreesToRadians(location.getLatitude()).doubleValue()));
899 location.getLatitude()).doubleValue()));
900 BigDecimal cosLatitude = BigDecimal.valueOf(Math.cos(convertDegreesToRadians( 869 BigDecimal cosLatitude = BigDecimal.valueOf(Math.cos(convertDegreesToRadians(location.getLatitude()).doubleValue()));
901 location.getLatitude()).doubleValue()));
902 870
903 BigDecimal sinDeclinationTimesSinLat = sinSunDeclination.multiply(sinLatitude); 871 BigDecimal sinDeclinationTimesSinLat = sinSunDeclination.multiply(sinLatitude);
904 BigDecimal dividend = cosineZenith.subtract(sinDeclinationTimesSinLat); 872 BigDecimal dividend = cosineZenith.subtract(sinDeclinationTimesSinLat);
905 BigDecimal divisor = cosineSunDeclination.multiply(cosLatitude); 873 BigDecimal divisor = cosineSunDeclination.multiply(cosLatitude);
906 874
907 return setScale(divideBy(dividend, divisor)); 875 return setScale(divideBy(dividend, divisor));
908 } 876 }
909 877
910 private BigDecimal getSinOfSunDeclination(BigDecimal sunTrueLong) 878 private BigDecimal getSinOfSunDeclination(BigDecimal sunTrueLong)
911 { 879 {
912 BigDecimal sinTrueLongitude = BigDecimal.valueOf(Math 880 BigDecimal sinTrueLongitude = BigDecimal.valueOf(Math.sin(convertDegreesToRadians(sunTrueLong).doubleValue()));
913 .sin(convertDegreesToRadians(sunTrueLong).doubleValue()));
914 BigDecimal sinOfDeclination = sinTrueLongitude.multiply(new BigDecimal("0.39782")); 881 BigDecimal sinOfDeclination = sinTrueLongitude.multiply(new BigDecimal("0.39782"));
915 return setScale(sinOfDeclination); 882 return setScale(sinOfDeclination);
916 } 883 }
917 884
918 private BigDecimal getCosineOfSunDeclination(BigDecimal sinSunDeclination) 885 private BigDecimal getCosineOfSunDeclination(BigDecimal sinSunDeclination)
919 { 886 {
920 BigDecimal arcSinOfSinDeclination = BigDecimal.valueOf(Math.asin(sinSunDeclination 887 BigDecimal arcSinOfSinDeclination = BigDecimal.valueOf(Math.asin(sinSunDeclination.doubleValue()));
921 .doubleValue()));
922 BigDecimal cosDeclination = BigDecimal
923 .valueOf(Math.cos(arcSinOfSinDeclination.doubleValue())); 888 BigDecimal cosDeclination = BigDecimal.valueOf(Math.cos(arcSinOfSinDeclination.doubleValue()));
924 return setScale(cosDeclination); 889 return setScale(cosDeclination);
925 } 890 }
926 891
927 private BigDecimal getSunLocalHour(BigDecimal cosineSunLocalHour, Boolean isSunrise) 892 private BigDecimal getSunLocalHour(BigDecimal cosineSunLocalHour, Boolean isSunrise)
928 { 893 {
933 localHour = BigDecimal.valueOf(360).subtract(localHour); 898 localHour = BigDecimal.valueOf(360).subtract(localHour);
934 } 899 }
935 return divideBy(localHour, BigDecimal.valueOf(15)); 900 return divideBy(localHour, BigDecimal.valueOf(15));
936 } 901 }
937 902
938 private BigDecimal getLocalMeanTime(BigDecimal sunTrueLong, BigDecimal longitudeHour, 903 private BigDecimal getLocalMeanTime(BigDecimal sunTrueLong, BigDecimal longitudeHour, BigDecimal sunLocalHour)
939 BigDecimal sunLocalHour)
940 { 904 {
941 BigDecimal rightAscension = this.getRightAscension(sunTrueLong); 905 BigDecimal rightAscension = this.getRightAscension(sunTrueLong);
942 BigDecimal innerParens = longitudeHour.multiply(new BigDecimal("0.06571")); 906 BigDecimal innerParens = longitudeHour.multiply(new BigDecimal("0.06571"));
943 BigDecimal localMeanTime = sunLocalHour.add(rightAscension).subtract(innerParens); 907 BigDecimal localMeanTime = sunLocalHour.add(rightAscension).subtract(innerParens);
944 localMeanTime = localMeanTime.subtract(new BigDecimal("6.622")); 908 localMeanTime = localMeanTime.subtract(new BigDecimal("6.622"));
978 942
979 /** 943 /**
980 * Returns the local rise/set time in the form HH:MM. 944 * Returns the local rise/set time in the form HH:MM.
981 * 945 *
982 * @param localTime 946 * @param localTime
983 * <code>BigDecimal</code> representation of the local rise/set time. 947 * <code>BigDecimal</code> representation of the local rise/set time.
984 * @return <code>String</code> representation of the local rise/set time in HH:MM format. 948 * @return <code>String</code> representation of the local rise/set time in HH:MM format.
985 */ 949 */
986 private String getLocalTimeAsString(BigDecimal localTime) 950 private String getLocalTimeAsString(BigDecimal localTime)
987 { 951 {
988 String[] timeComponents = localTime.toPlainString().split("\\."); 952 String[] timeComponents = localTime.toPlainString().split("\\.");
994 { 958 {
995 minutes = BigDecimal.ZERO; 959 minutes = BigDecimal.ZERO;
996 hour += 1; 960 hour += 1;
997 } 961 }
998 962
999 String minuteString = minutes.intValue() < 10 ? "0" + minutes.toPlainString() : minutes 963 String minuteString = minutes.intValue() < 10 ? "0" + minutes.toPlainString() : minutes.toPlainString();
1000 .toPlainString();
1001 String hourString = (hour < 10) ? "0" + String.valueOf(hour) : String.valueOf(hour); 964 String hourString = (hour < 10) ? "0" + String.valueOf(hour) : String.valueOf(hour);
1002 return hourString + ":" + minuteString; 965 return hourString + ":" + minuteString;
1003 } 966 }
1004 967
1005 /** ******* UTILITY METHODS (Should probably go somewhere else. ***************** */ 968 /** ******* UTILITY METHODS (Should probably go somewhere else. ***************** */

Legend:
Removed from v.30  
changed lines
  Added in v.31

   
Visit the ZANavi Wiki