pkg/darwin/darwin_insol.F

8fbfd1f382 Oliv*0001 #include "CPP_OPTIONS.h"
                0002
                0003 C     Procedure name: DARWIN_INSOL
                0004 C           Function: find time scale for plankton growth as
                0005 C                     function of light
                0006 c                     based on paltridge and parson
                0007 C           Comments: swd, April 1998
                0008 C                     following code by mick
                0009 C============================================================================
                0010 CStartofinterface
                0011       SUBROUTINE darwin_insol(Time,PARWm2,bj,returnsolz)
                0012       IMPLICIT NONE
                0013 C     === Global variables ===
                0014 #include "SIZE.h"
                0015 #include "EEPARAMS.h"
                0016 #include "PARAMS.h"
                0017 #include "FFIELDS.h"
                0018 #include "GRID.h"
                0019 #include "DYNVARS.h"
                0020 C     =============== Input ================================================
                0021 C     time :: seconds since January 1, 0:00 in a 360-day year
                0022       _RL time
                0023       integer bj
                0024       logical returnsolz
                0025 C     =============== Output ===============================================
                0026 C     PARWm2 :: daily-average PAR in W/m2
                0027       _RL PARWm2(1-OLy:sNy+OLy)
                0028
                0029 C============== Local variables ============================================
                0030       _RL  solar, albedo, par
                0031       _RL  dayfrac, yday, delta
                0032       _RL  lat, sun1, dayhrs
                0033       _RL  cosz, frac, fluxi
                0034       integer j
                0035 c
989120c6c8 Oliv*0036 C     solar constant (W/m2)
                0037       solar = 1360. _d 0
                0038 C     planetary albedo (W/W PAR range)
                0039       albedo = 0.6 _d 0
                0040 C     photosynthetically reactive fraction of PAR
                0041       par = 0.4 _d 0
8fbfd1f382 Oliv*0042 c
                0043 c find day (****NOTE for year starting in 1 Jan *****)
989120c6c8 Oliv*0044 C       fraction of year
8fbfd1f382 Oliv*0045         dayfrac=mod(Time,360. _d 0*86400. _d 0)
989120c6c8 Oliv*0046      &                    /(360. _d 0*86400. _d 0)
                0047 C       convert to radius
                0048         yday = 2.0 _d 0*3.1416 _d 0*dayfrac
                0049 C       declination (min.zenith angle) (Paltridge+Platt)
                0050         delta = (0.006918 _d 0- (0.399912 _d 0*cos(yday))
                0051      &          +(0.070257 _d 0*sin(yday))
8fbfd1f382 Oliv*0052      &          -(0.006758 _d 0*cos(2.0 _d 0*yday))
                0053      &          +(0.000907 _d 0*sin(2.0 _d 0*yday))
                0054      &          -(0.002697 _d 0*cos(3.0 _d 0*yday))
                0055      &          +(0.001480 _d 0*sin(3.0 _d 0*yday)) )
                0056        do j=1-OLy,sNy+OLy
                0057 c latitude in radians
                0058           lat=YC(1,j,1,bj)/180. _d 0*3.1416 _d 0
                0059           sun1 = -sin(delta)/cos(delta) * sin(lat)/cos(lat)
                0060           if (sun1.le.-0.999 _d 0) sun1=-0.999 _d 0
                0061           if (sun1.ge. 0.999 _d 0) sun1= 0.999 _d 0
                0062           dayhrs = abs(acos(sun1))
989120c6c8 Oliv*0063 C         average zenith angle
                0064           cosz = ( sin(delta)*sin(lat)+
8fbfd1f382 Oliv*0065      &            (cos(delta)*cos(lat)*sin(dayhrs)/dayhrs) )
                0066           if (cosz.le.0.005 _d 0) cosz=0.005 _d 0
                0067           if (returnsolz) THEN
                0068             if (cosz.ge. 0.999 _d 0) cosz= 0.999 _d 0
                0069             PARWm2(j) = acos(cosz)*180. _d 0/3.1416 _d 0
                0070           else
989120c6c8 Oliv*0071 C           fraction of daylight in day
                0072             frac = dayhrs/3.1416 _d 0
8fbfd1f382 Oliv*0073 c daily average photosynthetically active solar radiation just below surface
                0074             fluxi = solar*(1.0 _d 0-albedo)*cosz*frac*par
                0075 c
                0076 c convert to PARWm2
                0077             if (fluxi.gt.0.0 _d 0) PARWm2(j)=fluxi
                0078 c very large for polar night
                0079             if (fluxi.lt.0.00001 _d 0) PARWm2(j)=0.00001 _d 0
                0080           endif
989120c6c8 Oliv*0081 C      enddo j
                0082        enddo
8fbfd1f382 Oliv*0083 c
                0084       return
                0085       end
                0086 c
                0087 C==========================================================================