Home > atmlab > geographical > sun_angles.m

sun_angles

PURPOSE ^

SUN_ANGLES Finds solar conditions.

SYNOPSIS ^

function [sza,saa] = sun_angles(mjd,lat,lon)

DESCRIPTION ^

 SUN_ANGLES Finds solar conditions.
   [sza,saa] = sun_angles(mjd,lat,lon) calculates the Sun's zenith and
   azimuth angles for a given location and time, from Astronomical Almanac
   (see Michalsky, 1988).

  Refraction is not considered.

 See also: mjd2localtime, date2mjd

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

sun_angles.m

SOURCE CODE ^

0001 function [sza,saa] = sun_angles(mjd,lat,lon)
0002 % SUN_ANGLES Finds solar conditions.
0003 %   [sza,saa] = sun_angles(mjd,lat,lon) calculates the Sun's zenith and
0004 %   azimuth angles for a given location and time, from Astronomical Almanac
0005 %   (see Michalsky, 1988).
0006 %
0007 %  Refraction is not considered.
0008 %
0009 % See also: mjd2localtime, date2mjd
0010 
0011 %----------------------------------------------------------------
0012 % Craig Haley 20/09/01
0013 %   02-10-05 CSH vectorized
0014 %----------------------------------------------------------------
0015 
0016 DEGREE = pi/180.0;
0017 RADIAN = 180.0/pi;
0018 
0019 %calculate position of sun in celestial
0020 [ra,dec] = sun_celest(mjd);
0021 
0022 %calculate the Greenwich mean siderial time in hours
0023 gmst = mjd2gmst(mjd);
0024 
0025 %calculate local mean siderial time in radians
0026 lmst = gmst+lon/15.0;
0027 lmst = mod(lmst,24.0);
0028 lmst(lmst < 0) = lmst(lmst < 0)+24;
0029 lmst = lmst*15.0*DEGREE;
0030 
0031 %calculate hour angle in radians between -PI and PI
0032 %here we ignore the equation of the equinoxes factor
0033 ha = lmst-ra;
0034 ha(ha < -pi) = ha(ha < -pi)+2*pi;
0035 ha(ha > pi) = ha(ha > pi)-2*pi;
0036 
0037 %calculate altitude angle in radians
0038 a = asin(sin(dec).*sin(lat*DEGREE)+cos(dec).*cos(ha).*cos(lat*DEGREE));
0039 
0040 %calculate azimuth in radians between 0 and 2*PI
0041 az = atan2(sin(ha),cos(ha).*sin(lat*DEGREE)-tan(dec).*cos(lat*DEGREE))+pi;
0042 
0043 %     %calculate the refraction correction factor
0044 %     P = 1016.0;   %pressure in millibars
0045 %     T = 23.5;     %temperature in C
0046 %     if (a < 15.0)
0047 %         R = 0.00452*P/((273.0+T)*tan(a));
0048 %     else
0049 %         R = P*(0.1594+0.0196*a+0.00002*a*a)/((273.0+T)*(1.0+0.505*a+0.0845*a*a));
0050 %     end
0051 
0052 %calculate zenith and azimuth angles in degrees
0053 sza = 90.0-(a*RADIAN);        %add R to correct for refraction
0054 saa = az*RADIAN;

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005