4.76471

4.8 | 17 ratings Rate this file 110 downloads (last 30 days) File Size: 1.27 MB File ID: #4605

sun_position.m

by Vincent Roy

 

10 Mar 2004 (Updated 23 Feb 2009)

No BSD License  

Sun position given observer time/location.

Download Now | Watch this File

File Information
Description

sun = sun_position(time, location)

This function computes the sun position (zenith and azimuth angle at the observer location) as a function of the observer local time and position.

It is an implementation of the algorithm presented by Reda et Andreas in:
   Reda, I., Andreas, A. (2003) Solar position algorithm for solar radiation application. National Renewable Energy Laboratory (NREL) Technical report NREL/TP-560-34302, Revised January 2008.

This document is avalaible at http://rredc.nrel.gov/solar/codesandalgorithms/spa/ and is included in the .zip file

This algorithm is based on numerical approximation of the exact equations. The authors of the original paper state that this algorithm should be precise at +/- 0.0003 degrees. I have compared it to NOAA solar table (http://www.srrb.noaa.gov/highlights/sunrise/azel.html) and to USNO solar table (http://aa.usno.navy.mil/data/docs/AltAz.html) and found very good correspondance (up to the precision of those tables), except for large zenith angle, where the refraction by the atmosphere is significant (difference of about 1 degree). Note that in this code the correction for refraction in the atmosphere has been implemented for a temperature of 10C (283 kelvins) and a pressure of 1010 mbar. See the subfunction «sun_topocentric_zenith_angle_calculation» for possible modification to explicitely model the effects of temperature and pressure as describe in Reda & Andreas (2003).

Input parameters:
   time: a structure that specify the time when the sun position is calculated.
   time.year: year. Valid for [-2000, 6000]
   time.month: month [1-12]
   time.day: calendar day [1-31]
   time.hour: local hour [0-23]
   time.min: minute [0-59]
   time.sec: second [0-59]
   time.UTC: offset hour from UTC. Local time = Greenwich time + time.UTC
   This input can also be passed using the Matlab time format ('dd-mmm-yyyy HH:MM:SS'). In that case, the time has to be specified as UTC time (time.UTC = 0)

   location: a structure that specify the location of the observer
   location.latitude: latitude (in degrees, north of equator is positive)
   location.longitude: longitude (in degrees, positive for east of Greenwich) location.altitude: altitude above mean sea level (in meters)
 
Output parameters
   sun: a structure with the calculated sun position
   sun.zenith = zenith angle in degrees (angle from the vertical)
   sun.azimuth = azimuth angle in degrees, eastward from the north.

Only the sun zenith and azimuth angles are returned as output, but a lot of other parameters are calculated that could also be extracted as output of this function. See the documentation in the code.

Exemple of use

location.longitude = -105.1786;
location.latitude = 39.742476;
location.altitude = 1830.14;
time.year = 2003;
time.month = 10;
time.day = 17;
time.hour = 12;
time.min = 30;
time.sec = 30;
time.UTC = -7;

sun = sun_position(time, location);

 sun =
 
      zenith: 50.1080438859849
      azimuth: 194.341174010338

Acknowledgements
This submission has inspired the following:
Sun Azimuth Data
MATLAB release MATLAB 6.5 (R13)
Zip File Content  
Other Files 34302.pdf,
sun_position.m
Tags for This File  
Everyone's Tags
Tags I've Applied
Add New Tags Please login to tag files.
Comments and Ratings (19)
17 Mar 2004 Andrea Bergami

Very useful script!

07 Oct 2004 lu shengwang

It is a grate work and very halpful for me.Thanks a lot.

10 Jul 2005 Xin Jin

Bingo! This is exactly what I am looking for.

27 Aug 2005 Mat Saenger

Cool!

28 Aug 2005 frederik Steensen

Very high-quality. I have a translation to Delphi

18 Nov 2005 Steve Wales  
15 Feb 2006 DUNG TRAN THANH  
29 May 2006 nadeem nazar  
22 Jan 2007 George Young

The code works as advertise and is well documented. It even comes with a DOE report on the equations involved. The learning curve is extremely short. Thanks!

08 Feb 2007 Son Nguyen

Well report

03 Aug 2007 Andora Suwisto

Thanks for the program.....The reference from usno.navy is very well accepted and well known for scientist and engineer

24 Oct 2007 P Johnstone

This m-file is great! Very useful, well written.

17 Oct 2008 Lin Xinan

Thank you for your perfect program!!!

24 Oct 2008 Jeff Gullett

The algorithm has a flaw... The atmospheric refraction correction has a number of vertical asymptotes below -5 degrees true elevation to the sun. For elevation angles that hit these asymptotes, I was getting points that resulted in positive elevations, when the neighboring points would come out with the expected result of "sun well below the horizon". You need to add a conditional prior to applying the atmospheric correction to only add this in if the true elevation is greater than -5 degrees.

As for the code itself, it is somewhat disorganized and difficult to read, but it does work. I had a specific application where I needed to improved the efficiency/readability of the code, so I ended up rewriting the thing from scratch based on the paper, except I store everything in radians instead of degrees. This is still a nice application, though.

11 Nov 2008 Lu Salinas

Has anyone done this same script but in excel? I've been trying it and I still have some formula wrong... I guess I'm confuse with the t ut or the radians and degrees...
Would anyone want to give it a check?

19 Nov 2008 EDGAR RAMOS

Congrats guy!!

06 Jan 2009 Jeff Gullett

I have found a few other discrepancies in the algorithm. In the 2005 version of the document, the algorithm to compute the topocentric sun declination (eqn #39, line #825) contains the variable 'y' in the denominator, but the code presented later in the pdf document has the variable 'x' in its place. In the 2008 version of the algorithm document, this difference is resolved and the variable should be 'x' (as in the C-code, so it needs to be changed in the m-code presented on this page). Also, when computing the mean obliquity of the ecliptic (eqn #24, line #733), the documented steps have one of the terms as "-4680.93", but again, the C-code from the appendix used the term "-4680.96". Unfortunately, this discrepancy still exists in the 2008 document, so I don't know which term is correct. Using the wrong one will result in an error greater than the maximum error claimed by the algorithm designers.

23 Jan 2009 Falofolio

Great tool! Good job! However, I wonder, if there will be an update on this script, after reading Jeff Gullett's remarks on the code?

15 Jul 2009 K E

Very helpful submission

Please login to add a comment or rating.
Updates
15 Mar 2004

Small bugs fix (in julian day calculation)

18 Mar 2004

Small bugs fix (in julian day calculation)

18 Mar 2004

Changed information in the header (help display). location.elevation should have been location.altitude.

13 Apr 2004

Following a suggestion from Jody Klymak (jklymak@ucsd.edu), allowed the 'time' input to be passed as a Matlab time string.

23 Aug 2005

Modifications to the julian_calculation function. Bug was tranforming a string (time) to a structure directly. This was accepted until R14SP2. Script should now be compliant with future releases of Matlab...

23 Feb 2009

Thanks to the sharp eyes of Jeff Gullett, corrected eq. 39 as per updated Resa and Andreas paper (see line 825) and added a condition for the refractive correction (see lines 854-858, now only applied for sun elevation greater than -5 deg).

Tag Activity for this File
Tag Applied By Date/Time
sun position Vincent Roy 22 Oct 2008 07:16:03
zenith Vincent Roy 22 Oct 2008 07:16:03
elevation Vincent Roy 22 Oct 2008 07:16:03
earth sciences Vincent Roy 22 Oct 2008 07:16:04
azimuth Vincent Roy 22 Oct 2008 07:16:04
astronomy Ned Gulley 12 Nov 2008 00:32:54
elevation EDGAR RAMOS 19 Nov 2008 16:07:46
 

MATLAB Central Terms of Use

NOTICE: Any content you submit to MATLAB Central, including personal information, is not subject to the protections which may be afforded information collected under other sections of The MathWorks, Inc. Web site. You are entirely responsible for all content that you upload, post, e-mail, transmit or otherwise make available via MATLAB Central. The MathWorks does not control the content posted by visitors to MATLAB Central and, does not guarantee the accuracy, integrity, or quality of such content. Under no circumstances will The MathWorks be liable in any way for any content not authored by The MathWorks, or any loss or damage of any kind incurred as a result of the use of any content posted, e-mailed, transmitted or otherwise made available via MATLAB Central. Read the complete Terms prior to use.

Contact us at files@mathworks.com