This code is imported direct from http://rredc.nrel.gov/solar/codesandalgorithms/solpos/ Copyrights are retained in the code for each algorithm used. The library calculates the apparent solar position and intensity (theoretical maximum solar energy) based on the date, time, and location on Earth. The calculated data can be used to predict solar radiation, to be used in meteorological, solar energy and irrigation applications.

Files at this revision

API Documentation at this revision

Comitter:
jcobb
Date:
Tue Jul 20 22:45:48 2010 +0000
Commit message:

Changed in this revision

solpos.cpp Show annotated file Show diff for this revision Revisions of this file
solpos00.h Show annotated file Show diff for this revision Revisions of this file
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/solpos.cpp	Tue Jul 20 22:45:48 2010 +0000
@@ -0,0 +1,999 @@
+/*============================================================================
+*    Contains:
+*        S_solpos     (computes solar position and intensity
+*                      from time and place)
+*
+*            INPUTS:     (via posdata struct) year, daynum, hour,
+*                        minute, second, latitude, longitude, timezone,
+*                        intervl
+*            OPTIONAL:   (via posdata struct) month, day, press, temp, tilt,
+*                        aspect, function
+*            OUTPUTS:    EVERY variable in the struct posdata
+*                            (defined in solpos.h)
+*
+*                       NOTE: Certain conditions exist during which some of
+*                       the output variables are undefined or cannot be
+*                       calculated.  In these cases, the variables are
+*                       returned with flag values indicating such.  In other
+*                       cases, the variables may return a realistic, though
+*                       invalid, value. These variables and the flag values
+*                       or invalid conditions are listed below:
+*
+*                       amass     -1.0 at zenetr angles greater than 93.0
+*                                 degrees
+*                       ampress   -1.0 at zenetr angles greater than 93.0
+*                                 degrees
+*                       azim      invalid at zenetr angle 0.0 or latitude
+*                                 +/-90.0 or at night
+*                       elevetr   limited to -9 degrees at night
+*                       etr       0.0 at night
+*                       etrn      0.0 at night
+*                       etrtilt   0.0 when cosinc is less than 0
+*                       prime     invalid at zenetr angles greater than 93.0
+*                                 degrees
+*                       sretr     +/- 2999.0 during periods of 24 hour sunup or
+*                                 sundown
+*                       ssetr     +/- 2999.0 during periods of 24 hour sunup or
+*                                 sundown
+*                       ssha      invalid at the North and South Poles
+*                       unprime   invalid at zenetr angles greater than 93.0
+*                                 degrees
+*                       zenetr    limited to 99.0 degrees at night
+*
+*        S_init       (optional initialization for all input parameters in
+*                      the posdata struct)
+*           INPUTS:     struct posdata*
+*           OUTPUTS:    struct posdata*
+*
+*                     (Note: initializes the required S_solpos INPUTS above
+*                      to out-of-bounds conditions, forcing the user to
+*                      supply the parameters; initializes the OPTIONAL
+*                      S_solpos inputs above to nominal values.)
+*
+*       S_decode      (optional utility for decoding the S_solpos return code)
+*           INPUTS:     long integer S_solpos return value, struct posdata*
+*           OUTPUTS:    text to stderr
+*
+*    Usage:
+*         In calling program, just after other 'includes', insert:
+*
+*              #include "solpos00.h"
+*
+*         Function calls:
+*              S_init(struct posdata*)  [optional]
+*              .
+*              .
+*              [set time and location parameters before S_solpos call]
+*              .
+*              .
+*              int retval = S_solpos(struct posdata*)
+*              S_decode(int retval, struct posdata*) [optional]
+*                  (Note: you should always look at the S_solpos return
+*                   value, which contains error codes. S_decode is one option
+*                   for examining these codes.  It can also serve as a
+*                   template for building your own application-specific
+*                   decoder.)
+*
+*    Martin Rymes
+*    National Renewable Energy Laboratory
+*    25 March 1998
+*
+*    27 April 1999 REVISION:  Corrected leap year in S_date.
+*    13 January 2000 REVISION:  SMW converted to structure posdata parameter
+*                               and subdivided into functions.
+*    01 February 2001 REVISION: SMW corrected ecobli calculation 
+*                               (changed sign). Error is small (max 0.015 deg
+*                               in calculation of declination angle)
+*----------------------------------------------------------------------------*/
+#include <math.h>
+#include <string.h>
+#include <stdio.h>
+#include "solpos00.h"
+
+/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+*
+* Structures defined for this module
+*
+*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
+struct trigdata /* used to pass calculated values locally */
+{
+    float cd;       /* cosine of the declination */
+    float ch;       /* cosine of the hour angle */
+    float cl;       /* cosine of the latitude */
+    float sd;       /* sine of the declination */
+    float sl;       /* sine of the latitude */
+};
+
+
+/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+*
+* Temporary global variables used only in this file:
+*
+*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
+  static int  month_days[2][13] = { { 0,   0,  31,  59,  90, 120, 151,
+                                       181, 212, 243, 273, 304, 334 },
+                                    { 0,   0,  31,  60,  91, 121, 152,
+                                       182, 213, 244, 274, 305, 335 } };
+                   /* cumulative number of days prior to beginning of month */
+
+  static float degrad = 57.295779513; /* converts from radians to degrees */
+  static float raddeg = 0.0174532925; /* converts from degrees to radians */
+
+/*============================================================================
+*    Local function prototypes
+============================================================================*/
+static long int validate ( struct posdata *pdat);
+static void dom2doy( struct posdata *pdat );
+static void doy2dom( struct posdata *pdat );
+static void geometry ( struct posdata *pdat );
+static void zen_no_ref ( struct posdata *pdat, struct trigdata *tdat );
+static void ssha( struct posdata *pdat, struct trigdata *tdat );
+static void sbcf( struct posdata *pdat, struct trigdata *tdat );
+static void tst( struct posdata *pdat );
+static void srss( struct posdata *pdat );
+static void sazm( struct posdata *pdat, struct trigdata *tdat );
+static void refrac( struct posdata *pdat );
+static void amass( struct posdata *pdat );
+static void prime( struct posdata *pdat );
+static void etr( struct posdata *pdat );
+static void tilt( struct posdata *pdat );
+static void localtrig( struct posdata *pdat, struct trigdata *tdat );
+
+/*============================================================================
+*    Long integer function S_solpos, adapted from the VAX solar libraries
+*
+*    This function calculates the apparent solar position and the
+*    intensity of the sun (theoretical maximum solar energy) from
+*    time and place on Earth.
+*
+*    Requires (from the struct posdata parameter):
+*        Date and time:
+*            year
+*            daynum   (requirement depends on the S_DOY switch)
+*            month    (requirement depends on the S_DOY switch)
+*            day      (requirement depends on the S_DOY switch)
+*            hour
+*            minute
+*            second
+*            interval  DEFAULT 0
+*        Location:
+*            latitude
+*            longitude
+*        Location/time adjuster:
+*            timezone
+*        Atmospheric pressure and temperature:
+*            press     DEFAULT 1013.0 mb
+*            temp      DEFAULT 10.0 degrees C
+*        Tilt of flat surface that receives solar energy:
+*            aspect    DEFAULT 180 (South)
+*            tilt      DEFAULT 0 (Horizontal)
+*        Function Switch (codes defined in solpos.h)
+*            function  DEFAULT S_ALL
+*
+*    Returns (via the struct posdata parameter):
+*        everything defined in the struct posdata in solpos.h.
+*----------------------------------------------------------------------------*/
+long S_solpos (struct posdata *pdat)
+{
+  long int retval;
+
+  struct trigdata trigdat, *tdat;
+
+  tdat = &trigdat;   /* point to the structure */
+
+  /* initialize the trig structure */
+  tdat->sd = -999.0; /* flag to force calculation of trig data */
+  tdat->cd =    1.0;
+  tdat->ch =    1.0; /* set the rest of these to something safe */
+  tdat->cl =    1.0;
+  tdat->sl =    1.0;
+
+  if ((retval = validate ( pdat )) != 0) /* validate the inputs */
+    return retval;
+
+
+  if ( pdat->function & L_DOY )
+    doy2dom( pdat );                /* convert input doy to month-day */
+  else
+    dom2doy( pdat );                /* convert input month-day to doy */
+
+  if ( pdat->function & L_GEOM )
+    geometry( pdat );               /* do basic geometry calculations */
+
+  if ( pdat->function & L_ZENETR )  /* etr at non-refracted zenith angle */
+    zen_no_ref( pdat, tdat );
+
+  if ( pdat->function & L_SSHA )    /* Sunset hour calculation */
+    ssha( pdat, tdat );
+
+  if ( pdat->function & L_SBCF )    /* Shadowband correction factor */
+    sbcf( pdat, tdat );
+
+  if ( pdat->function & L_TST )     /* true solar time */
+    tst( pdat );
+
+  if ( pdat->function & L_SRSS )    /* sunrise/sunset calculations */
+    srss( pdat );
+
+  if ( pdat->function & L_SOLAZM )  /* solar azimuth calculations */
+    sazm( pdat, tdat );
+
+  if ( pdat->function & L_REFRAC )  /* atmospheric refraction calculations */
+    refrac( pdat );
+
+  if ( pdat->function & L_AMASS )   /* airmass calculations */
+    amass( pdat );
+
+  if ( pdat->function & L_PRIME )   /* kt-prime/unprime calculations */
+    prime( pdat );
+
+  if ( pdat->function & L_ETR )     /* ETR and ETRN (refracted) */
+    etr( pdat );
+
+  if ( pdat->function & L_TILT )    /* tilt calculations */
+    tilt( pdat );
+
+    return 0;
+}
+
+
+/*============================================================================
+*    Void function S_init
+*
+*    This function initiates all of the input parameters in the struct
+*    posdata passed to S_solpos().  Initialization is either to nominal
+*    values or to out of range values, which forces the calling program to
+*    specify parameters.
+*
+*    NOTE: This function is optional if you initialize ALL input parameters
+*          in your calling code.  Note that the required parameters of date
+*          and location are deliberately initialized out of bounds to force
+*          the user to enter real-world values.
+*
+*    Requires: Pointer to a posdata structure, members of which are
+*           initialized.
+*
+*    Returns: Void
+*----------------------------------------------------------------------------*/
+void S_init(struct posdata *pdat)
+{
+  pdat->day       =    -99;   /* Day of month (May 27 = 27, etc.) */
+  pdat->daynum    =   -999;   /* Day number (day of year; Feb 1 = 32 ) */
+  pdat->hour      =    -99;   /* Hour of day, 0 - 23 */
+  pdat->minute    =    -99;   /* Minute of hour, 0 - 59 */
+  pdat->month     =    -99;   /* Month number (Jan = 1, Feb = 2, etc.) */
+  pdat->second    =    -99;   /* Second of minute, 0 - 59 */
+  pdat->year      =    -99;   /* 4-digit year */
+  pdat->interval  =      0;   /* instantaneous measurement interval */
+  pdat->aspect    =  180.0;   /* Azimuth of panel surface (direction it
+                                    faces) N=0, E=90, S=180, W=270 */
+  pdat->latitude  =  -99.0;   /* Latitude, degrees north (south negative) */
+  pdat->longitude = -999.0;   /* Longitude, degrees east (west negative) */
+  pdat->press     = 1013.0;   /* Surface pressure, millibars */
+  pdat->solcon    = 1367.0;   /* Solar constant, 1367 W/sq m */
+  pdat->temp      =   15.0;   /* Ambient dry-bulb temperature, degrees C */
+  pdat->tilt      =    0.0;   /* Degrees tilt from horizontal of panel */
+  pdat->timezone  =  -99.0;   /* Time zone, east (west negative). */
+  pdat->sbwid     =    7.6;   /* Eppley shadow band width */
+  pdat->sbrad     =   31.7;   /* Eppley shadow band radius */
+  pdat->sbsky     =   0.04;   /* Drummond factor for partly cloudy skies */
+  pdat->function  =  S_ALL;   /* compute all parameters */
+}
+
+
+/*============================================================================
+*    Local long int function validate
+*
+*    Validates the input parameters
+*----------------------------------------------------------------------------*/
+static long int validate ( struct posdata *pdat)
+{
+
+  long int retval = 0;  /* start with no errors */
+
+  /* No absurd dates, please. */
+  if ( pdat->function & L_GEOM )
+  {
+    if ( (pdat->year < 1950) || (pdat->year > 2050) ) /* limits of algoritm */
+      retval |= (1L << S_YEAR_ERROR);
+    if ( !(pdat->function & S_DOY) && ((pdat->month < 1) || (pdat->month > 12)))
+      retval |= (1L << S_MONTH_ERROR);
+    if ( !(pdat->function & S_DOY) && ((pdat->day < 1) || (pdat->day > 31)) )
+      retval |= (1L << S_DAY_ERROR);
+    if ( (pdat->function & S_DOY) && ((pdat->daynum < 1) || (pdat->daynum > 366)) )
+      retval |= (1L << S_DOY_ERROR);
+
+    /* No absurd times, please. */
+    if ( (pdat->hour < 0) || (pdat->hour > 24) )
+      retval |= (1L << S_HOUR_ERROR);
+    if ( (pdat->minute < 0) || (pdat->minute > 59) )
+      retval |= (1L << S_MINUTE_ERROR);
+    if ( (pdat->second < 0) || (pdat->second > 59) )
+      retval |= (1L << S_SECOND_ERROR);
+    if ( (pdat->hour == 24) && (pdat->minute > 0) ) /* no more than 24 hrs */
+      retval |= ( (1L << S_HOUR_ERROR) | (1L << S_MINUTE_ERROR) );
+    if ( (pdat->hour == 24) && (pdat->second > 0) ) /* no more than 24 hrs */
+      retval |= ( (1L << S_HOUR_ERROR) | (1L << S_SECOND_ERROR) );
+    if ( fabs (pdat->timezone) > 12.0 )
+      retval |= (1L << S_TZONE_ERROR);
+    if ( (pdat->interval < 0) || (pdat->interval > 28800) )
+      retval |= (1L << S_INTRVL_ERROR);
+
+    /* No absurd locations, please. */
+    if ( fabs (pdat->longitude) > 180.0 )
+      retval |= (1L << S_LON_ERROR);
+    if ( fabs (pdat->latitude) > 90.0 )
+      retval |= (1L << S_LAT_ERROR);
+  }
+
+  /* No silly temperatures or pressures, please. */
+  if ( (pdat->function & L_REFRAC) && (fabs (pdat->temp) > 100.0) )
+    retval |= (1L << S_TEMP_ERROR);
+  if ( (pdat->function & L_REFRAC) &&
+    (pdat->press < 0.0) || (pdat->press > 2000.0) )
+    retval |= (1L << S_PRESS_ERROR);
+
+  /* No out of bounds tilts, please */
+  if ( (pdat->function & L_TILT) && (fabs (pdat->tilt) > 180.0) )
+    retval |= (1L << S_TILT_ERROR);
+  if ( (pdat->function & L_TILT) && (fabs (pdat->aspect) > 360.0) )
+    retval |= (1L << S_ASPECT_ERROR);
+
+  /* No oddball shadowbands, please */
+  if ( (pdat->function & L_SBCF) &&
+       (pdat->sbwid < 1.0) || (pdat->sbwid > 100.0) )
+    retval |= (1L << S_SBWID_ERROR);
+  if ( (pdat->function & L_SBCF) &&
+       (pdat->sbrad < 1.0) || (pdat->sbrad > 100.0) )
+    retval |= (1L << S_SBRAD_ERROR);
+  if ( (pdat->function & L_SBCF) && ( fabs (pdat->sbsky) > 1.0) )
+    retval |= (1L << S_SBSKY_ERROR);
+
+  return retval;
+}
+
+
+/*============================================================================
+*    Local Void function dom2doy
+*
+*    Converts day-of-month to day-of-year
+*
+*    Requires (from struct posdata parameter):
+*            year
+*            month
+*            day
+*
+*    Returns (via the struct posdata parameter):
+*            year
+*            daynum
+*----------------------------------------------------------------------------*/
+static void dom2doy( struct posdata *pdat )
+{
+  pdat->daynum = pdat->day + month_days[0][pdat->month];
+
+  /* (adjust for leap year) */
+  if ( ((pdat->year % 4) == 0) &&
+         ( ((pdat->year % 100) != 0) || ((pdat->year % 400) == 0) ) &&
+         (pdat->month > 2) )
+      pdat->daynum += 1;
+}
+
+
+/*============================================================================
+*    Local void function doy2dom
+*
+*    This function computes the month/day from the day number.
+*
+*    Requires (from struct posdata parameter):
+*        Year and day number:
+*            year
+*            daynum
+*
+*    Returns (via the struct posdata parameter):
+*            year
+*            month
+*            day
+*----------------------------------------------------------------------------*/
+static void doy2dom(struct posdata *pdat)
+{
+  int  imon;  /* Month (month_days) array counter */
+  int  leap;  /* leap year switch */
+
+    /* Set the leap year switch */
+    if ( ((pdat->year % 4) == 0) &&
+         ( ((pdat->year % 100) != 0) || ((pdat->year % 400) == 0) ) )
+        leap = 1;
+    else
+        leap = 0;
+
+    /* Find the month */
+    imon = 12;
+    while ( pdat->daynum <= month_days [leap][imon] )
+        --imon;
+
+    /* Set the month and day of month */
+    pdat->month = imon;
+    pdat->day   = pdat->daynum - month_days[leap][imon];
+}
+
+
+/*============================================================================
+*    Local Void function geometry
+*
+*    Does the underlying geometry for a given time and location
+*----------------------------------------------------------------------------*/
+static void geometry ( struct posdata *pdat )
+{
+  float bottom;      /* denominator (bottom) of the fraction */
+  float c2;          /* cosine of d2 */
+  float cd;          /* cosine of the day angle or delination */
+  float d2;          /* pdat->dayang times two */
+  float delta;       /* difference between current year and 1949 */
+  float s2;          /* sine of d2 */
+  float sd;          /* sine of the day angle */
+  float top;         /* numerator (top) of the fraction */
+  int   leap;        /* leap year counter */
+
+  /* Day angle */
+      /*  Iqbal, M.  1983.  An Introduction to Solar Radiation.
+            Academic Press, NY., page 3 */
+     pdat->dayang = 360.0 * ( pdat->daynum - 1 ) / 365.0;
+
+    /* Earth radius vector * solar constant = solar energy */
+        /*  Spencer, J. W.  1971.  Fourier series representation of the
+            position of the sun.  Search 2 (5), page 172 */
+    sd     = sin (raddeg * pdat->dayang);
+    cd     = cos (raddeg * pdat->dayang);
+    d2     = 2.0 * pdat->dayang;
+    c2     = cos (raddeg * d2);
+    s2     = sin (raddeg * d2);
+
+    pdat->erv  = 1.000110 + 0.034221 * cd + 0.001280 * sd;
+    pdat->erv  += 0.000719 * c2 + 0.000077 * s2;
+
+    /* Universal Coordinated (Greenwich standard) time */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->utime =
+        pdat->hour * 3600.0 +
+        pdat->minute * 60.0 +
+        pdat->second -
+        (float)pdat->interval / 2.0;
+    pdat->utime = pdat->utime / 3600.0 - pdat->timezone;
+
+    /* Julian Day minus 2,400,000 days (to eliminate roundoff errors) */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+
+    /* No adjustment for century non-leap years since this function is
+       bounded by 1950 - 2050 */
+    delta    = pdat->year - 1949;
+    leap     = (int) ( delta / 4.0 );
+    pdat->julday =
+        32916.5 + delta * 365.0 + leap + pdat->daynum + pdat->utime / 24.0;
+
+    /* Time used in the calculation of ecliptic coordinates */
+    /* Noon 1 JAN 2000 = 2,400,000 + 51,545 days Julian Date */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->ectime = pdat->julday - 51545.0;
+
+    /* Mean longitude */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->mnlong  = 280.460 + 0.9856474 * pdat->ectime;
+
+    /* (dump the multiples of 360, so the answer is between 0 and 360) */
+    pdat->mnlong -= 360.0 * (int) ( pdat->mnlong / 360.0 );
+    if ( pdat->mnlong < 0.0 )
+        pdat->mnlong += 360.0;
+
+    /* Mean anomaly */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->mnanom  = 357.528 + 0.9856003 * pdat->ectime;
+
+    /* (dump the multiples of 360, so the answer is between 0 and 360) */
+    pdat->mnanom -= 360.0 * (int) ( pdat->mnanom / 360.0 );
+    if ( pdat->mnanom < 0.0 )
+        pdat->mnanom += 360.0;
+
+    /* Ecliptic longitude */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->eclong  = pdat->mnlong + 1.915 * sin ( pdat->mnanom * raddeg ) +
+                    0.020 * sin ( 2.0 * pdat->mnanom * raddeg );
+
+    /* (dump the multiples of 360, so the answer is between 0 and 360) */
+    pdat->eclong -= 360.0 * (int) ( pdat->eclong / 360.0 );
+    if ( pdat->eclong < 0.0 )
+        pdat->eclong += 360.0;
+
+    /* Obliquity of the ecliptic */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+
+    /* 02 Feb 2001 SMW corrected sign in the following line */
+/*  pdat->ecobli = 23.439 + 4.0e-07 * pdat->ectime;     */
+    pdat->ecobli = 23.439 - 4.0e-07 * pdat->ectime;
+
+    /* Declination */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->declin = degrad * asin ( sin (pdat->ecobli * raddeg) *
+                               sin (pdat->eclong * raddeg) );
+
+    /* Right ascension */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    top      =  cos ( raddeg * pdat->ecobli ) * sin ( raddeg * pdat->eclong );
+    bottom   =  cos ( raddeg * pdat->eclong );
+
+    pdat->rascen =  degrad * atan2 ( top, bottom );
+
+    /* (make it a positive angle) */
+    if ( pdat->rascen < 0.0 )
+        pdat->rascen += 360.0;
+
+    /* Greenwich mean sidereal time */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->gmst  = 6.697375 + 0.0657098242 * pdat->ectime + pdat->utime;
+
+    /* (dump the multiples of 24, so the answer is between 0 and 24) */
+    pdat->gmst -= 24.0 * (int) ( pdat->gmst / 24.0 );
+    if ( pdat->gmst < 0.0 )
+        pdat->gmst += 24.0;
+
+    /* Local mean sidereal time */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->lmst  = pdat->gmst * 15.0 + pdat->longitude;
+
+    /* (dump the multiples of 360, so the answer is between 0 and 360) */
+    pdat->lmst -= 360.0 * (int) ( pdat->lmst / 360.0 );
+    if ( pdat->lmst < 0.)
+        pdat->lmst += 360.0;
+
+    /* Hour angle */
+        /*  Michalsky, J.  1988.  The Astronomical Almanac's algorithm for
+            approximate solar position (1950-2050).  Solar Energy 40 (3),
+            pp. 227-235. */
+    pdat->hrang = pdat->lmst - pdat->rascen;
+
+    /* (force it between -180 and 180 degrees) */
+    if ( pdat->hrang < -180.0 )
+        pdat->hrang += 360.0;
+    else if ( pdat->hrang > 180.0 )
+        pdat->hrang -= 360.0;
+}
+
+
+/*============================================================================
+*    Local Void function zen_no_ref
+*
+*    ETR solar zenith angle
+*       Iqbal, M.  1983.  An Introduction to Solar Radiation.
+*            Academic Press, NY., page 15
+*----------------------------------------------------------------------------*/
+static void zen_no_ref ( struct posdata *pdat, struct trigdata *tdat )
+{
+  float cz;          /* cosine of the solar zenith angle */
+
+    localtrig( pdat, tdat );
+    cz = tdat->sd * tdat->sl + tdat->cd * tdat->cl * tdat->ch;
+
+    /* (watch out for the roundoff errors) */
+    if ( fabs (cz) > 1.0 ) {
+        if ( cz >= 0.0 )
+            cz =  1.0;
+        else
+            cz = -1.0;
+    }
+
+    pdat->zenetr   = acos ( cz ) * degrad;
+
+    /* (limit the degrees below the horizon to 9 [+90 -> 99]) */
+    if ( pdat->zenetr > 99.0 )
+        pdat->zenetr = 99.0;
+
+    pdat->elevetr = 90.0 - pdat->zenetr;
+}
+
+
+/*============================================================================
+*    Local Void function ssha
+*
+*    Sunset hour angle, degrees
+*       Iqbal, M.  1983.  An Introduction to Solar Radiation.
+*            Academic Press, NY., page 16
+*----------------------------------------------------------------------------*/
+static void ssha( struct posdata *pdat, struct trigdata *tdat )
+{
+  float cssha;       /* cosine of the sunset hour angle */
+  float cdcl;        /* ( cd * cl ) */
+
+    localtrig( pdat, tdat );
+    cdcl    = tdat->cd * tdat->cl;
+
+    if ( fabs ( cdcl ) >= 0.001 ) {
+        cssha = -tdat->sl * tdat->sd / cdcl;
+
+        /* This keeps the cosine from blowing on roundoff */
+        if ( cssha < -1.0  )
+            pdat->ssha = 180.0;
+        else if ( cssha > 1.0 )
+            pdat->ssha = 0.0;
+        else
+            pdat->ssha = degrad * acos ( cssha );
+    }
+    else if ( ((pdat->declin >= 0.0) && (pdat->latitude > 0.0 )) ||
+              ((pdat->declin <  0.0) && (pdat->latitude < 0.0 )) )
+        pdat->ssha = 180.0;
+    else
+        pdat->ssha = 0.0;
+}
+
+
+/*============================================================================
+*    Local Void function sbcf
+*
+*    Shadowband correction factor
+*       Drummond, A. J.  1956.  A contribution to absolute pyrheliometry.
+*            Q. J. R. Meteorol. Soc. 82, pp. 481-493
+*----------------------------------------------------------------------------*/
+static void sbcf( struct posdata *pdat, struct trigdata *tdat )
+{
+  float p, t1, t2;   /* used to compute sbcf */
+
+    localtrig( pdat, tdat );
+    p       = 0.6366198 * pdat->sbwid / pdat->sbrad * pow (tdat->cd,3);
+    t1      = tdat->sl * tdat->sd * pdat->ssha * raddeg;
+    t2      = tdat->cl * tdat->cd * sin ( pdat->ssha * raddeg );
+    pdat->sbcf = pdat->sbsky + 1.0 / ( 1.0 - p * ( t1 + t2 ) );
+
+}
+
+
+/*============================================================================
+*    Local Void function tst
+*
+*    TST -> True Solar Time = local standard time + TSTfix, time
+*      in minutes from midnight.
+*        Iqbal, M.  1983.  An Introduction to Solar Radiation.
+*            Academic Press, NY., page 13
+*----------------------------------------------------------------------------*/
+static void tst( struct posdata *pdat )
+{
+    pdat->tst    = ( 180.0 + pdat->hrang ) * 4.0;
+    pdat->tstfix =
+        pdat->tst -
+        (float)pdat->hour * 60.0 -
+        pdat->minute -
+        (float)pdat->second / 60.0 +
+        (float)pdat->interval / 120.0; /* add back half of the interval */
+
+    /* bound tstfix to this day */
+    while ( pdat->tstfix >  720.0 )
+        pdat->tstfix -= 1440.0;
+    while ( pdat->tstfix < -720.0 )
+        pdat->tstfix += 1440.0;
+
+    pdat->eqntim =
+        pdat->tstfix + 60.0 * pdat->timezone - 4.0 * pdat->longitude;
+
+}
+
+
+/*============================================================================
+*    Local Void function srss
+*
+*    Sunrise and sunset times (minutes from midnight)
+*----------------------------------------------------------------------------*/
+static void srss( struct posdata *pdat )
+{
+    if ( pdat->ssha <= 1.0 ) {
+        pdat->sretr   =  2999.0;
+        pdat->ssetr   = -2999.0;
+    }
+    else if ( pdat->ssha >= 179.0 ) {
+        pdat->sretr   = -2999.0;
+        pdat->ssetr   =  2999.0;
+    }
+    else {
+        pdat->sretr   = 720.0 - 4.0 * pdat->ssha - pdat->tstfix;
+        pdat->ssetr   = 720.0 + 4.0 * pdat->ssha - pdat->tstfix;
+    }
+}
+
+
+/*============================================================================
+*    Local Void function sazm
+*
+*    Solar azimuth angle
+*       Iqbal, M.  1983.  An Introduction to Solar Radiation.
+*            Academic Press, NY., page 15
+*----------------------------------------------------------------------------*/
+static void sazm( struct posdata *pdat, struct trigdata *tdat )
+{
+  float ca;          /* cosine of the solar azimuth angle */
+  float ce;          /* cosine of the solar elevation */
+  float cecl;        /* ( ce * cl ) */
+  float se;          /* sine of the solar elevation */
+
+    localtrig( pdat, tdat );
+    ce         = cos ( raddeg * pdat->elevetr );
+    se         = sin ( raddeg * pdat->elevetr );
+
+    pdat->azim     = 180.0;
+    cecl       = ce * tdat->cl;
+    if ( fabs ( cecl ) >= 0.001 ) {
+        ca     = ( se * tdat->sl - tdat->sd ) / cecl;
+        if ( ca > 1.0 )
+            ca = 1.0;
+        else if ( ca < -1.0 )
+            ca = -1.0;
+
+        pdat->azim = 180.0 - acos ( ca ) * degrad;
+        if ( pdat->hrang > 0 )
+            pdat->azim  = 360.0 - pdat->azim;
+    }
+}
+
+
+/*============================================================================
+*    Local Int function refrac
+*
+*    Refraction correction, degrees
+*        Zimmerman, John C.  1981.  Sun-pointing programs and their
+*            accuracy.
+*            SAND81-0761, Experimental Systems Operation Division 4721,
+*            Sandia National Laboratories, Albuquerque, NM.
+*----------------------------------------------------------------------------*/
+static void refrac( struct posdata *pdat )
+{
+  float prestemp;    /* temporary pressure/temperature correction */
+  float refcor;      /* temporary refraction correction */
+  float tanelev;     /* tangent of the solar elevation angle */
+
+    /* If the sun is near zenith, the algorithm bombs; refraction near 0 */
+    if ( pdat->elevetr > 85.0 )
+        refcor = 0.0;
+
+    /* Otherwise, we have refraction */
+    else {
+        tanelev = tan ( raddeg * pdat->elevetr );
+        if ( pdat->elevetr >= 5.0 )
+            refcor  = 58.1 / tanelev -
+                      0.07 / ( pow (tanelev,3) ) +
+                      0.000086 / ( pow (tanelev,5) );
+        else if ( pdat->elevetr >= -0.575 )
+            refcor  = 1735.0 +
+                      pdat->elevetr * ( -518.2 + pdat->elevetr * ( 103.4 +
+                      pdat->elevetr * ( -12.79 + pdat->elevetr * 0.711 ) ) );
+        else
+            refcor  = -20.774 / tanelev;
+
+        prestemp    =
+            ( pdat->press * 283.0 ) / ( 1013.0 * ( 273.0 + pdat->temp ) );
+        refcor     *= prestemp / 3600.0;
+    }
+
+    /* Refracted solar elevation angle */
+    pdat->elevref = pdat->elevetr + refcor;
+
+    /* (limit the degrees below the horizon to 9) */
+    if ( pdat->elevref < -9.0 )
+        pdat->elevref = -9.0;
+
+    /* Refracted solar zenith angle */
+    pdat->zenref  = 90.0 - pdat->elevref;
+    pdat->coszen  = cos( raddeg * pdat->zenref );
+}
+
+
+/*============================================================================
+*    Local Void function  amass
+*
+*    Airmass
+*       Kasten, F. and Young, A.  1989.  Revised optical air mass
+*            tables and approximation formula.  Applied Optics 28 (22),
+*            pp. 4735-4738
+*----------------------------------------------------------------------------*/
+static void amass( struct posdata *pdat )
+{
+    if ( pdat->zenref > 93.0 )
+    {
+        pdat->amass   = -1.0;
+        pdat->ampress = -1.0;
+    }
+    else
+    {
+        pdat->amass =
+            1.0 / ( cos (raddeg * pdat->zenref) + 0.50572 *
+            pow ((96.07995 - pdat->zenref),-1.6364) );
+
+        pdat->ampress   = pdat->amass * pdat->press / 1013.0;
+    }
+}
+
+
+/*============================================================================
+*    Local Void function prime
+*
+*    Prime and Unprime
+*    Prime  converts Kt to normalized Kt', etc.
+*       Unprime deconverts Kt' to Kt, etc.
+*            Perez, R., P. Ineichen, Seals, R., & Zelenka, A.  1990.  Making
+*            full use of the clearness index for parameterizing hourly
+*            insolation conditions. Solar Energy 45 (2), pp. 111-114
+*----------------------------------------------------------------------------*/
+static void prime( struct posdata *pdat )
+{
+    pdat->unprime = 1.031 * exp ( -1.4 / ( 0.9 + 9.4 / pdat->amass ) ) + 0.1;
+    pdat->prime   = 1.0 / pdat->unprime;
+}
+
+
+/*============================================================================
+*    Local Void function etr
+*
+*    Extraterrestrial (top-of-atmosphere) solar irradiance
+*----------------------------------------------------------------------------*/
+static void etr( struct posdata *pdat )
+{
+    if ( pdat->coszen > 0.0 ) {
+        pdat->etrn = pdat->solcon * pdat->erv;
+        pdat->etr  = pdat->etrn * pdat->coszen;
+    }
+    else {
+        pdat->etrn = 0.0;
+        pdat->etr  = 0.0;
+    }
+}
+
+
+/*============================================================================
+*    Local Void function localtrig
+*
+*    Does trig on internal variable used by several functions
+*----------------------------------------------------------------------------*/
+static void localtrig( struct posdata *pdat, struct trigdata *tdat )
+{
+/* define masks to prevent calculation of uninitialized variables */
+#define SD_MASK ( L_ZENETR | L_SSHA | S_SBCF | S_SOLAZM )
+#define SL_MASK ( L_ZENETR | L_SSHA | S_SBCF | S_SOLAZM )
+#define CL_MASK ( L_ZENETR | L_SSHA | S_SBCF | S_SOLAZM )
+#define CD_MASK ( L_ZENETR | L_SSHA | S_SBCF )
+#define CH_MASK ( L_ZENETR )
+
+    if ( tdat->sd < -900.0 )  /* sd was initialized -999 as flag */
+    {
+      tdat->sd = 1.0;  /* reflag as having completed calculations */
+      if ( pdat->function | CD_MASK )
+        tdat->cd = cos ( raddeg * pdat->declin );
+      if ( pdat->function | CH_MASK )
+        tdat->ch = cos ( raddeg * pdat->hrang );
+      if ( pdat->function | CL_MASK )
+        tdat->cl = cos ( raddeg * pdat->latitude );
+      if ( pdat->function | SD_MASK )
+        tdat->sd = sin ( raddeg * pdat->declin );
+      if ( pdat->function | SL_MASK )
+        tdat->sl = sin ( raddeg * pdat->latitude );
+    }
+}
+
+
+/*============================================================================
+*    Local Void function tilt
+*
+*    ETR on a tilted surface
+*----------------------------------------------------------------------------*/
+static void tilt( struct posdata *pdat )
+{
+  float ca;          /* cosine of the solar azimuth angle */
+  float cp;          /* cosine of the panel aspect */
+  float ct;          /* cosine of the panel tilt */
+  float sa;          /* sine of the solar azimuth angle */
+  float sp;          /* sine of the panel aspect */
+  float st;          /* sine of the panel tilt */
+  float sz;          /* sine of the refraction corrected solar zenith angle */
+
+
+    /* Cosine of the angle between the sun and a tipped flat surface,
+       useful for calculating solar energy on tilted surfaces */
+    ca      = cos ( raddeg * pdat->azim );
+    cp      = cos ( raddeg * pdat->aspect );
+    ct      = cos ( raddeg * pdat->tilt );
+    sa      = sin ( raddeg * pdat->azim );
+    sp      = sin ( raddeg * pdat->aspect );
+    st      = sin ( raddeg * pdat->tilt );
+    sz      = sin ( raddeg * pdat->zenref );
+    pdat->cosinc  = pdat->coszen * ct + sz * st * ( ca * cp + sa * sp );
+
+    if ( pdat->cosinc > 0.0 )
+        pdat->etrtilt = pdat->etrn * pdat->cosinc;
+    else
+        pdat->etrtilt = 0.0;
+
+}
+
+
+/*============================================================================
+*    Void function S_decode
+*
+*    This function decodes the error codes from S_solpos return value
+*
+*    Requires the long integer return value from S_solpos
+*
+*    Returns descriptive text to stderr
+*----------------------------------------------------------------------------*/
+void S_decode(long code, struct posdata *pdat)
+{
+  if ( code & (1L << S_YEAR_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the year: %d [1950-2050]\n",
+      pdat->year);
+  if ( code & (1L << S_MONTH_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the month: %d\n",
+      pdat->month);
+  if ( code & (1L << S_DAY_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the day-of-month: %d\n",
+      pdat->day);
+  if ( code & (1L << S_DOY_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the day-of-year: %d\n",
+      pdat->daynum);
+  if ( code & (1L << S_HOUR_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the hour: %d\n",
+      pdat->hour);
+  if ( code & (1L << S_MINUTE_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the minute: %d\n",
+      pdat->minute);
+  if ( code & (1L << S_SECOND_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the second: %d\n",
+      pdat->second);
+  if ( code & (1L << S_TZONE_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the time zone: %f\n",
+      pdat->timezone);
+  if ( code & (1L << S_INTRVL_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the interval: %d\n",
+      pdat->interval);
+  if ( code & (1L << S_LAT_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the latitude: %f\n",
+      pdat->latitude);
+  if ( code & (1L << S_LON_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the longitude: %f\n",
+      pdat->longitude);
+  if ( code & (1L << S_TEMP_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the temperature: %f\n",
+      pdat->temp);
+  if ( code & (1L << S_PRESS_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the pressure: %f\n",
+      pdat->press);
+  if ( code & (1L << S_TILT_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the tilt: %f\n",
+      pdat->tilt);
+  if ( code & (1L << S_ASPECT_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the aspect: %f\n",
+      pdat->aspect);
+  if ( code & (1L << S_SBWID_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the shadowband width: %f\n",
+      pdat->sbwid);
+  if ( code & (1L << S_SBRAD_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the shadowband radius: %f\n",
+      pdat->sbrad);
+  if ( code & (1L << S_SBSKY_ERROR) )
+    fprintf(stderr, "S_decode ==> Please fix the shadowband sky factor: %f\n",
+      pdat->sbsky);
+}
+
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/solpos00.h	Tue Jul 20 22:45:48 2010 +0000
@@ -0,0 +1,364 @@
+ /*============================================================================
+*
+*    NAME:  solpos00.h
+*
+*    Contains:
+*        S_solpos     (computes the solar position and intensity
+*                      from time and place)
+*            INPUTS:     (from posdata)
+*                          year, month, day, hour, minute, second,
+*                          latitude, longitude, timezone, interval
+*            OPTIONAL:   (from posdata; defaults from S_init function)
+*                            press   DEFAULT 1013.0 (standard pressure)
+*                            temp    DEFAULT   10.0 (standard temperature)
+*                            tilt    DEFAULT    0.0 (horizontal panel)
+*                            aspect  DEFAULT  180.0 (South-facing panel)
+*                            sbwid   DEFAULT    7.6 (shadowband width)
+*                            sbrad   DEFAULT   31.7 (shadowband radius)
+*                            sbsky   DEFAULT   0.04 (shadowband sky factor)
+*
+*            OUTPUTS:    (posdata) daynum, amass, ampress, azim, cosinc,
+*                        elevref, etr, etrn, etrtilt, prime,
+*                        sbcf, sretr, ssetr, unprime, zenref
+*
+*            RETURNS:   Long int status code (defined in solpos.h)
+*
+*    Usage:
+*         In calling program, along with other 'includes', insert:
+*
+*              #include "solpos.h"
+*
+*    Martin Rymes
+*    National Renewable Energy Laboratory
+*    25 March 1998
+*----------------------------------------------------------------------------*/
+
+/*============================================================================
+*
+*     Define the function codes
+*
+*----------------------------------------------------------------------------*/
+#define L_DOY    0x0001
+#define L_GEOM   0x0002
+#define L_ZENETR 0x0004
+#define L_SSHA   0x0008
+#define L_SBCF   0x0010
+#define L_TST    0x0020
+#define L_SRSS   0x0040
+#define L_SOLAZM 0x0080
+#define L_REFRAC 0x0100
+#define L_AMASS  0x0200
+#define L_PRIME  0x0400
+#define L_TILT   0x0800
+#define L_ETR    0x1000
+#define L_ALL    0xFFFF
+
+/*============================================================================
+*
+*     Define the bit-wise masks for each function
+*
+*----------------------------------------------------------------------------*/
+#define S_DOY    ( L_DOY                          )
+#define S_GEOM   ( L_GEOM   | S_DOY               )
+#define S_ZENETR ( L_ZENETR | S_GEOM              )
+#define S_SSHA   ( L_SSHA   | S_GEOM              )
+#define S_SBCF   ( L_SBCF   | S_SSHA              )
+#define S_TST    ( L_TST    | S_GEOM              )
+#define S_SRSS   ( L_SRSS   | S_SSHA   | S_TST    )
+#define S_SOLAZM ( L_SOLAZM | S_ZENETR            )
+#define S_REFRAC ( L_REFRAC | S_ZENETR            )
+#define S_AMASS  ( L_AMASS  | S_REFRAC            )
+#define S_PRIME  ( L_PRIME  | S_AMASS             )
+#define S_TILT   ( L_TILT   | S_SOLAZM | S_REFRAC )
+#define S_ETR    ( L_ETR    | S_REFRAC            )
+#define S_ALL    ( L_ALL                          )
+
+
+/*============================================================================
+*
+*     Enumerate the error codes
+*     (Bit positions are from least significant to most significant)
+*
+*----------------------------------------------------------------------------*/
+/*          Code          Bit       Parameter            Range
+      ===============     ===  ===================  =============   */
+enum {S_YEAR_ERROR,    /*  0   year                  1950 -  2050   */
+      S_MONTH_ERROR,   /*  1   month                    1 -    12   */
+      S_DAY_ERROR,     /*  2   day-of-month             1 -    31   */
+      S_DOY_ERROR,     /*  3   day-of-year              1 -   366   */
+      S_HOUR_ERROR,    /*  4   hour                     0 -    24   */
+      S_MINUTE_ERROR,  /*  5   minute                   0 -    59   */
+      S_SECOND_ERROR,  /*  6   second                   0 -    59   */
+      S_TZONE_ERROR,   /*  7   time zone              -12 -    12   */
+      S_INTRVL_ERROR,  /*  8   interval (seconds)       0 - 28800   */
+      S_LAT_ERROR,     /*  9   latitude               -90 -    90   */
+      S_LON_ERROR,     /* 10   longitude             -180 -   180   */
+      S_TEMP_ERROR,    /* 11   temperature (deg. C)  -100 -   100   */
+      S_PRESS_ERROR,   /* 12   pressure (millibars)     0 -  2000   */
+      S_TILT_ERROR,    /* 13   tilt                   -90 -    90   */
+      S_ASPECT_ERROR,  /* 14   aspect                -360 -   360   */
+      S_SBWID_ERROR,   /* 15   shadow band width (cm)   1 -   100   */
+      S_SBRAD_ERROR,   /* 16   shadow band radius (cm)  1 -   100   */
+      S_SBSKY_ERROR};  /* 17   shadow band sky factor  -1 -     1   */
+
+struct posdata
+{
+  /***** ALPHABETICAL LIST OF COMMON VARIABLES *****/
+                           /* Each comment begins with a 1-column letter code:
+                              I:  INPUT variable
+                              O:  OUTPUT variable
+                              T:  TRANSITIONAL variable used in the algorithm,
+                                  of interest only to the solar radiation
+                                  modelers, and available to you because you
+                                  may be one of them.
+
+                              The FUNCTION column indicates which sub-function
+                              within solpos must be switched on using the
+                              "function" parameter to calculate the desired
+                              output variable.  All function codes are
+                              defined in the solpos.h file.  The default
+                              S_ALL switch calculates all output variables.
+                              Multiple functions may be or'd to create a
+                              composite function switch.  For example,
+                              (S_TST | S_SBCF). Specifying only the functions
+                              for required output variables may allow solpos
+                              to execute more quickly.
+
+                              The S_DOY mask works as a toggle between the
+                              input date represented as a day number (daynum)
+                              or as month and day.  To set the switch (to
+                              use daynum input), the function is or'd; to
+                              clear the switch (to use month and day input),
+                              the function is inverted and and'd.
+
+                              For example:
+                                  pdat->function |= S_DOY (sets daynum input)
+                                  pdat->function &= ~S_DOY (sets month and day input)
+
+                              Whichever date form is used, S_solpos will
+                              calculate and return the variables(s) of the
+                              other form.  See the soltest.c program for
+                              other examples. */
+
+  /* VARIABLE        I/O  Function    Description */
+  /* -------------  ----  ----------  ---------------------------------------*/
+
+  int   day;       /* I/O: S_DOY      Day of month (May 27 = 27, etc.)
+                                        solpos will CALCULATE this by default,
+                                        or will optionally require it as input
+                                        depending on the setting of the S_DOY
+                                        function switch. */
+  int   daynum;    /* I/O: S_DOY      Day number (day of year; Feb 1 = 32 )
+                                        solpos REQUIRES this by default, but
+                                        will optionally calculate it from
+                                        month and day depending on the setting
+                                        of the S_DOY function switch. */
+  int   function;  /* I:              Switch to choose functions for desired
+                                        output. */
+  int   hour;      /* I:              Hour of day, 0 - 23, DEFAULT = 12 */
+  int   interval;  /* I:              Interval of a measurement period in
+                                        seconds.  Forces solpos to use the
+                                        time and date from the interval
+                                        midpoint. The INPUT time (hour,
+                                        minute, and second) is assumed to
+                                        be the END of the measurement
+                                        interval. */
+  int   minute;    /* I:              Minute of hour, 0 - 59, DEFAULT = 0 */
+  int   month;     /* I/O: S_DOY      Month number (Jan = 1, Feb = 2, etc.)
+                                        solpos will CALCULATE this by default,
+                                        or will optionally require it as input
+                                        depending on the setting of the S_DOY
+                                        function switch. */
+  int   second;    /* I:              Second of minute, 0 - 59, DEFAULT = 0 */
+  int   year;      /* I:              4-digit year (2-digit year is NOT
+                                       allowed */
+
+  /***** FLOATS *****/
+
+  float amass;      /* O:  S_AMASS    Relative optical airmass */
+  float ampress;    /* O:  S_AMASS    Pressure-corrected airmass */
+  float aspect;     /* I:             Azimuth of panel surface (direction it
+                                        faces) N=0, E=90, S=180, W=270,
+                                        DEFAULT = 180 */
+  float azim;       /* O:  S_SOLAZM   Solar azimuth angle:  N=0, E=90, S=180,
+                                        W=270 */
+  float cosinc;     /* O:  S_TILT     Cosine of solar incidence angle on
+                                        panel */
+  float coszen;     /* O:  S_REFRAC   Cosine of refraction corrected solar
+                                        zenith angle */
+  float dayang;     /* T:  S_GEOM     Day angle (daynum*360/year-length)
+                                        degrees */
+  float declin;     /* T:  S_GEOM     Declination--zenith angle of solar noon
+                                        at equator, degrees NORTH */
+  float eclong;     /* T:  S_GEOM     Ecliptic longitude, degrees */
+  float ecobli;     /* T:  S_GEOM     Obliquity of ecliptic */
+  float ectime;     /* T:  S_GEOM     Time of ecliptic calculations */
+  float elevetr;    /* O:  S_ZENETR   Solar elevation, no atmospheric
+                                        correction (= ETR) */
+  float elevref;    /* O:  S_REFRAC   Solar elevation angle,
+                                        deg. from horizon, refracted */
+  float eqntim;     /* T:  S_TST      Equation of time (TST - LMT), minutes */
+  float erv;        /* T:  S_GEOM     Earth radius vector
+                                        (multiplied to solar constant) */
+  float etr;        /* O:  S_ETR      Extraterrestrial (top-of-atmosphere)
+                                        W/sq m global horizontal solar
+                                        irradiance */
+  float etrn;       /* O:  S_ETR      Extraterrestrial (top-of-atmosphere)
+                                        W/sq m direct normal solar
+                                        irradiance */
+  float etrtilt;    /* O:  S_TILT     Extraterrestrial (top-of-atmosphere)
+                                        W/sq m global irradiance on a tilted
+                                        surface */
+  float gmst;       /* T:  S_GEOM     Greenwich mean sidereal time, hours */
+  float hrang;      /* T:  S_GEOM     Hour angle--hour of sun from solar noon,
+                                        degrees WEST */
+  float julday;     /* T:  S_GEOM     Julian Day of 1 JAN 2000 minus
+                                        2,400,000 days (in order to regain
+                                        single precision) */
+  float latitude;   /* I:             Latitude, degrees north (south negative) */
+  float longitude;  /* I:             Longitude, degrees east (west negative) */
+  float lmst;       /* T:  S_GEOM     Local mean sidereal time, degrees */
+  float mnanom;     /* T:  S_GEOM     Mean anomaly, degrees */
+  float mnlong;     /* T:  S_GEOM     Mean longitude, degrees */
+  float rascen;     /* T:  S_GEOM     Right ascension, degrees */
+  float press;      /* I:             Surface pressure, millibars, used for
+                                        refraction correction and ampress */
+  float prime;      /* O:  S_PRIME    Factor that normalizes Kt, Kn, etc. */
+  float sbcf;       /* O:  S_SBCF     Shadow-band correction factor */
+  float sbwid;      /* I:             Shadow-band width (cm) */
+  float sbrad;      /* I:             Shadow-band radius (cm) */
+  float sbsky;      /* I:             Shadow-band sky factor */
+  float solcon;     /* I:             Solar constant (NREL uses 1367 W/sq m) */
+  float ssha;       /* T:  S_SRHA     Sunset(/rise) hour angle, degrees */
+  float sretr;      /* O:  S_SRSS     Sunrise time, minutes from midnight,
+                                        local, WITHOUT refraction */
+  float ssetr;      /* O:  S_SRSS     Sunset time, minutes from midnight,
+                                        local, WITHOUT refraction */
+  float temp;       /* I:             Ambient dry-bulb temperature, degrees C,
+                                        used for refraction correction */
+  float tilt;       /* I:             Degrees tilt from horizontal of panel */
+  float timezone;   /* I:             Time zone, east (west negative).
+                                      USA:  Mountain = -7, Central = -6, etc. */
+  float tst;        /* T:  S_TST      True solar time, minutes from midnight */
+  float tstfix;     /* T:  S_TST      True solar time - local standard time */
+  float unprime;    /* O:  S_PRIME    Factor that denormalizes Kt', Kn', etc. */
+  float utime;      /* T:  S_GEOM     Universal (Greenwich) standard time */
+  float zenetr;     /* T:  S_ZENETR   Solar zenith angle, no atmospheric
+                                        correction (= ETR) */
+  float zenref;     /* O:  S_REFRAC   Solar zenith angle, deg. from zenith,
+                                        refracted */
+};
+
+/* For users that wish to access individual functions, the following table
+lists all output and transition variables, the L_ mask for the function
+that calculates them, and all the input variables required by that function.
+The function variable is set to the L_ mask, which will force S_solpos to
+only call the required function.  L_ masks may be ORed as desired.
+
+VARIABLE      Mask       Required Variables
+---------  ----------  ---------------------------------------
+ amass      L_AMASS    zenref, press
+ ampress    L_AMASS    zenref, press
+ azim       L_SOLAZM   elevetr, declin, latitude, hrang
+ cosinc     L_TILT     azim, aspect, tilt, zenref, coszen,etrn
+ coszen     L_REFRAC   elevetr, press, temp
+ dayang     L_GEOM     All date, time, and location inputs
+ declin     L_GEOM     All date, time, and location inputs
+ eclong     L_GEOM     All date, time, and location inputs
+ ecobli     L_GEOM     All date, time, and location inputs
+ ectime     L_GEOM     All date, time, and location inputs
+ elevetr    L_ZENETR   declin, latitude, hrang
+ elevref    L_REFRAC   elevetr, press, temp
+ eqntim     L_TST      hrang, hour, minute, second, interval
+ erv        L_GEOM     All date, time, and location inputs
+ etr        L_ETR      coszen, solcon, erv
+ etrn       L_ETR      coszen, solcon, erv
+ etrtilt    L_TILT     azim, aspect, tilt, zenref, coszen, etrn
+ gmst       L_GEOM     All date, time, and location inputs
+ hrang      L_GEOM     All date, time, and location inputs
+ julday     L_GEOM     All date, time, and location inputs
+ lmst       L_GEOM     All date, time, and location inputs
+ mnanom     L_GEOM     All date, time, and location inputs
+ mnlong     L_GEOM     All date, time, and location inputs
+ rascen     L_GEOM     All date, time, and location inputs
+ prime      L_PRIME    amass
+ sbcf       L_SBCF     latitude, declin, ssha, sbwid, sbrad, sbsky
+ ssha       L_SRHA     latitude, declin
+ sretr      L_SRSS     ssha, tstfix
+ ssetr      L_SRSS     ssha, tstfix
+ tst        L_TST      hrang, hour, minute, second, interval
+ tstfix     L_TST      hrang, hour, minute, second, interval
+ unprime    L_PRIME    amass
+ utime      L_GEOM     All date, time, and location inputs
+ zenetr     L_ZENETR   declination, latitude, hrang
+ zenref     L_REFRAC   elevetr, press, temp
+
+
+
+*============================================================================
+*    Long int function S_solpos, adapted from the NREL VAX solar libraries
+*
+*    This function calculates the apparent solar position and intensity
+*    (theoretical maximum solar energy) based on the date, time, and
+*    location on Earth. (DEFAULT values are from the optional S_posinit
+*    function.)
+*
+*    Requires:
+*        Date and time:
+*            year
+*            month  (optional without daynum)
+*            day    (optional without daynum)
+*            daynum
+*            hour
+*            minute
+*            second
+*        Location:
+*            latitude
+*            longitude
+*        Location/time adjuster:
+*            timezone
+*        Atmospheric pressure and temperature:
+*            press     DEFAULT 1013.0 mb
+*            temp      DEFAULT 10.0 degrees C
+*        Tilt of flat surface that receives solar energy:
+*            aspect    DEFAULT 180 (South)
+*            tilt      DEFAULT 0 (Horizontal)
+*        Shadow band parameters:
+*            sbwid     DEFAULT 7.6 cm
+*            sbrad     DEFAULT 31.7 cm
+*            sbsky     DEFAULT 0.04
+*        Functionality
+*            function  DEFAULT S_ALL (all output parameters computed)
+*
+*    Returns:
+*        everything defined at the top of this listing.
+*----------------------------------------------------------------------------*/
+long S_solpos (struct posdata *pdat);
+
+/*============================================================================
+*    Void function S_init
+*
+*    This function initiates all of the input functions to S_Solpos().
+*    NOTE: This function is optional if you initialize all input parameters
+*          in your calling code.
+*
+*    Requires: Pointer to a posdata structure, members of which are
+*           initialized.
+*
+*    Returns: Void
+*
+*----------------------------------------------------------------------------*/
+void S_init(struct posdata *pdat);
+
+
+/*============================================================================
+*    Void function S_decode
+*
+*    This function decodes the error codes from S_solpos return value
+*
+*    INPUTS: Long integer S_solpos return value, struct posdata*
+*
+*    OUTPUTS: Descriptive text of errors to stderr
+*----------------------------------------------------------------------------*/
+void S_decode(long code, struct posdata *pdat);
+