c$Author: fredfox $
c$Date: 2002-05-08 18:47:45 $
c$Revision: 1.2 $
c$Source: /weru/cvs/weps/weps.src/util/solar/dawn.for,v $

      real function   dawn(dlat,dlong,idoy,riseangle)

c     + + + PURPOSE + + +
c     This function calculates the time of sunrise (hours) for any simulation
c     site based on the global position of the site, and day of the
c     year.  The inputs for the function are day of the year, latitude
c     of the site, and longitude of the site.

c     + + + KEYWORDS + + +
c     sunrise

c     + + + ARGUMENT DECLARATIONS + + +
      integer idoy
      real dlat
      real dlong
      real riseangle

c     + + + ARGUMENT DEFINITIONS + + +
c     idoy   - Day of year
c     dlat   - Latitude of the site, degrees (north > 0, south < 0)
c     dlong  - Longitude of the site, degrees (east > 0, west < 0)
c     riseangle - angle of earths rotation where sunrise occurs
c                 this varies depending on whether you are calculating
c                 direct beam, civil twilight, nautical twilight or
c                 astronomical twilight hourangle

c     + + + LOCAL VARIABLES + + +
      real b, dec, e, h, coshr, sn
      real dlat_rad, dec_rad

c     + + + LOCAL DEFINITIONS + + +
c     b      - sub calculation (time of year, radians)
c     dec    - declination of earth wrt the sun (degrees)
c     e      - Equation of time (minutes)
c     h      - Hour angle (degrees)
c     coshr  - Cosine of hour angle at sunrise
c     sn     - Solar noon (hour of the day, midnight = 0.0)

c     + + + COMMON BLOCKS + + +

      include 'p1const.inc'
      include 'p1unconv.inc'

c     + + + FUNCTION DECLARATIONS + + +
      real declination
      real hourangle
      real equa_time

c     + + + END SPECIFICATIONS + + +

c     declination angle (dec)
      dec = declination(idoy)

c     sunset hour angle (noon is zero degrees, sunset (+), sunrise (-))
      h = hourangle(dlat, dec, riseangle)

c     equation of time (e)
      e = equa_time(idoy)

c     Calculate solar noon (sn)
      sn = 12.0-e/60.0-4.0*mod(dlong,15.0)/60.0				!h-53

c     Calculate the time of sunrise (rise)
      dawn = sn - h/15.0									!h-52

c     to prevent errors of bleed over into previous day where
c     where daylength is 24 hours, limit time of sunrise
      dawn = max(0.0, dawn)

      return
      end