/[zanavi_public1]/navit/navit/android/src/bpi/sdbm/illuminance/SolarPosition.java
ZANavi

Diff of /navit/navit/android/src/bpi/sdbm/illuminance/SolarPosition.java

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

Revision 30 Revision 31
42 */ 42 */
43 43
44public class SolarPosition 44public class SolarPosition
45{ 45{
46 46
47 private static final double dEarthMeanRadius = 6371.01; // in km 47 private static final double dEarthMeanRadius = 6371.01; // in km
48 private static final double dAstronomicalUnit = 149597890; // in km 48 private static final double dAstronomicalUnit = 149597890; // in km
49 49
50 private static final double pi = Math.PI; 50 private static final double pi = Math.PI;
51 private static final double twopi = (2 * pi); 51 private static final double twopi = (2 * pi);
52 private static final double rad = (pi / 180); 52 private static final double rad = (pi / 180);
53
54 53
55 /** result wrapper class */ 54 /** result wrapper class */
56 public static class SunCoordinates 55 public static class SunCoordinates
57 { 56 {
58 /** 57 /**
59 * zenith angle, in degrees 58 * zenith angle, in degrees
60 */ 59 */
61 public double zenithAngle; 60 public double zenithAngle;
62 61
63 /** 62 /**
64 * azimuth, in degrees, 0 - 360 degrees clockwise from north 63 * azimuth, in degrees, 0 - 360 degrees clockwise from north
65 */ 64 */
66 public double azimuth; 65 public double azimuth;
67 } 66 }
68 67
69 /** 68 /**
70 * Calculate sun position for a given time and location. 69 * Calculate sun position for a given time and location.
71 * 70 *
72 * @param time 71 * @param time
73 * Note that it's unclear how well the algorithm performs before the year 1990 or after the year 2015. 72 * Note that it's unclear how well the algorithm performs before the year 1990 or after the year 2015.
74 * @param latitude 73 * @param latitude
75 * (positive east of Greenwich) 74 * (positive east of Greenwich)
76 * @param longitude 75 * @param longitude
77 * (positive north of equator) 76 * (positive north of equator)
78 * @return 77 * @return
79 */ 78 */
80 public static SunCoordinates getSunPosition(final Date time, double latitude, double longitude) 79 public static SunCoordinates getSunPosition(final Date time, double latitude, double longitude)
81 { 80 {
82 81
103 { 102 {
104 long liAux1; 103 long liAux1;
105 long liAux2; 104 long liAux2;
106 double dJulianDate; 105 double dJulianDate;
107 // Calculate time of the day in UT decimal hours 106 // Calculate time of the day in UT decimal hours
108 dDecimalHours = utcTime.get(Calendar.HOUR_OF_DAY)
109 + (utcTime.get(Calendar.MINUTE) + utcTime.get(Calendar.SECOND) / 60.0) / 60.0; 107 dDecimalHours = utcTime.get(Calendar.HOUR_OF_DAY) + (utcTime.get(Calendar.MINUTE) + utcTime.get(Calendar.SECOND) / 60.0) / 60.0;
110 // Calculate current Julian Day 108 // Calculate current Julian Day
111 liAux1 = (utcTime.get(Calendar.MONTH) + 1 - 14) / 12; 109 liAux1 = (utcTime.get(Calendar.MONTH) + 1 - 14) / 12;
112 liAux2 = (1461 * (utcTime.get(Calendar.YEAR) + 4800 + liAux1)) / 4 110 liAux2 = (1461 * (utcTime.get(Calendar.YEAR) + 4800 + liAux1)) / 4 + (367 * (utcTime.get(Calendar.MONTH) + 1 - 2 - 12 * liAux1)) / 12 - (3 * ((utcTime.get(Calendar.YEAR) + 4900 + liAux1) / 100)) / 4 + utcTime.get(Calendar.DAY_OF_MONTH) - 32075;
113 + (367 * (utcTime.get(Calendar.MONTH) + 1 - 2 - 12 * liAux1)) / 12
114 - (3 * ((utcTime.get(Calendar.YEAR) + 4900 + liAux1) / 100)) / 4
115 + utcTime.get(Calendar.DAY_OF_MONTH) - 32075;
116 dJulianDate = (double) (liAux2) - 0.5 + dDecimalHours / 24.0; 111 dJulianDate = (double) (liAux2) - 0.5 + dDecimalHours / 24.0;
117 // Calculate difference between current Julian Day and JD 2451545.0 112 // Calculate difference between current Julian Day and JD 2451545.0
118 dElapsedJulianDays = dJulianDate - 2451545.0; 113 dElapsedJulianDays = dJulianDate - 2451545.0;
119 } 114 }
120 115
131 double dMeanAnomaly; 126 double dMeanAnomaly;
132 double dOmega; 127 double dOmega;
133 dOmega = 2.1429 - 0.0010394594 * dElapsedJulianDays; 128 dOmega = 2.1429 - 0.0010394594 * dElapsedJulianDays;
134 dMeanLongitude = 4.8950630 + 0.017202791698 * dElapsedJulianDays; // Radians 129 dMeanLongitude = 4.8950630 + 0.017202791698 * dElapsedJulianDays; // Radians
135 dMeanAnomaly = 6.2400600 + 0.0172019699 * dElapsedJulianDays; 130 dMeanAnomaly = 6.2400600 + 0.0172019699 * dElapsedJulianDays;
136 dEclipticLongitude = dMeanLongitude + 0.03341607 * Math.sin(dMeanAnomaly) + 0.00034894 131 dEclipticLongitude = dMeanLongitude + 0.03341607 * Math.sin(dMeanAnomaly) + 0.00034894 * Math.sin(2 * dMeanAnomaly) - 0.0001134 - 0.0000203 * Math.sin(dOmega);
137 * Math.sin(2 * dMeanAnomaly) - 0.0001134 - 0.0000203 * Math.sin(dOmega);
138 dEclipticObliquity = 0.4090928 - 6.2140e-9 * dElapsedJulianDays + 0.0000396 132 dEclipticObliquity = 0.4090928 - 6.2140e-9 * dElapsedJulianDays + 0.0000396 * Math.cos(dOmega);
139 * Math.cos(dOmega);
140 } 133 }
141 134
142 // System.err.println("ecl. longi. " + dEclipticLongitude); 135 // System.err.println("ecl. longi. " + dEclipticLongitude);
143 // System.err.println("ecl. obliq. " + dEclipticObliquity); 136 // System.err.println("ecl. obliq. " + dEclipticObliquity);
144
145 137
146 // Calculate celestial coordinates ( right ascension and declination ) 138 // Calculate celestial coordinates ( right ascension and declination )
147 // in radians 139 // in radians
148 // but without limiting the angle to be less than 2*Pi (i.e., the result 140 // but without limiting the angle to be less than 2*Pi (i.e., the result
149 // may be 141 // may be
169 double dHourAngle; 161 double dHourAngle;
170 double dCos_Latitude; 162 double dCos_Latitude;
171 double dSin_Latitude; 163 double dSin_Latitude;
172 double dCos_HourAngle; 164 double dCos_HourAngle;
173 double dParallax; 165 double dParallax;
174 dGreenwichMeanSiderealTime = 6.6974243242 + 0.0657098283 * dElapsedJulianDays 166 dGreenwichMeanSiderealTime = 6.6974243242 + 0.0657098283 * dElapsedJulianDays + dDecimalHours;
175 + dDecimalHours;
176 dLocalMeanSiderealTime = (dGreenwichMeanSiderealTime * 15 + longitude) * rad; 167 dLocalMeanSiderealTime = (dGreenwichMeanSiderealTime * 15 + longitude) * rad;
177 dHourAngle = dLocalMeanSiderealTime - dRightAscension; 168 dHourAngle = dLocalMeanSiderealTime - dRightAscension;
178 dLatitudeInRadians = latitude * rad; 169 dLatitudeInRadians = latitude * rad;
179 dCos_Latitude = Math.cos(dLatitudeInRadians); 170 dCos_Latitude = Math.cos(dLatitudeInRadians);
180 dSin_Latitude = Math.sin(dLatitudeInRadians); 171 dSin_Latitude = Math.sin(dLatitudeInRadians);
181 dCos_HourAngle = Math.cos(dHourAngle); 172 dCos_HourAngle = Math.cos(dHourAngle);
182 retval.zenithAngle = (Math.acos(dCos_Latitude * dCos_HourAngle * Math.cos(dDeclination) 173 retval.zenithAngle = (Math.acos(dCos_Latitude * dCos_HourAngle * Math.cos(dDeclination) + Math.sin(dDeclination) * dSin_Latitude));
183 + Math.sin(dDeclination) * dSin_Latitude));
184 dY = -Math.sin(dHourAngle); 174 dY = -Math.sin(dHourAngle);
185 dX = Math.tan(dDeclination) * dCos_Latitude - dSin_Latitude * dCos_HourAngle; 175 dX = Math.tan(dDeclination) * dCos_Latitude - dSin_Latitude * dCos_HourAngle;
186 retval.azimuth = Math.atan2(dY, dX); 176 retval.azimuth = Math.atan2(dY, dX);
187 if (retval.azimuth < 0.0) retval.azimuth = retval.azimuth + twopi; 177 if (retval.azimuth < 0.0) retval.azimuth = retval.azimuth + twopi;
188 retval.azimuth = retval.azimuth / rad; 178 retval.azimuth = retval.azimuth / rad;
192 } 182 }
193 183
194 return retval; 184 return retval;
195 } 185 }
196 186
197
198 public static void main(String[] args) 187 public static void main(String[] args)
199 { 188 {
200 189
201 SolarPosition.getSunPosition(new Date(), 48.2, 16.4); 190 SolarPosition.getSunPosition(new Date(), 48.2, 16.4);
202 //System.out.println("current azimuth " + sc.azimuth); 191 //System.out.println("current azimuth " + sc.azimuth);

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

   
Visit the ZANavi Wiki