*IF DEF,A01_1A,OR,DEF,A01_1B,OR,DEF,A01_2A,OR,DEF,A01_2B,OR,DEF,A01_3A AWI3F402.6
C ******************************COPYRIGHT****************************** GTS2F400.9325
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.9326
C GTS2F400.9327
C Use, duplication or disclosure of this code is subject to the GTS2F400.9328
C restrictions as set forth in the contract. GTS2F400.9329
C GTS2F400.9330
C Meteorological Office GTS2F400.9331
C London Road GTS2F400.9332
C BRACKNELL GTS2F400.9333
C Berkshire UK GTS2F400.9334
C RG12 2SZ GTS2F400.9335
C GTS2F400.9336
C If no contract has been raised with this copy of the code, the use, GTS2F400.9337
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.9338
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.9339
C Modelling at the above address. GTS2F400.9340
C ******************************COPYRIGHT****************************** GTS2F400.9341
C GTS2F400.9342
CLL Subroutine SOLPOS ---------------------------------------------- SOLPOS1A.3
CLL SOLPOS1A.4
CLL Purpose : SOLPOS1A.5
CLL Calculations of the earth's orbit described in the first page of SOLPOS1A.6
CLL the "Calculation of incoming insolation" section of UMDP 23, i.e. SOLPOS1A.7
CLL from the day of the year (and, in forecast mode, whether it is a SOLPOS1A.8
CLL leap year) and the orbital "constants" (which vary over SOLPOS1A.9
CLL "Milankovitch" timescales) it calculates the sin of the solar SOLPOS1A.10
CLL declination and the inverse-square scaling factor for the solar SOLPOS1A.11
CLL "constant". It is thus intrinsically scalar. The FORTRAN code SOLPOS1A.12
CLL present depends on whether *DEF CAL360 is set during UPDATE: this SOLPOS1A.13
CLL replaces the Julian calendar with the climate-mode 360-day calendar SOLPOS1A.14
CLL SOLPOS1A.15
CLL Author: William Ingram SOLPOS1A.16
CLL SOLPOS1A.17
CLL Model Modification history from model version 3.0: SOLPOS1A.18
CLL version Date SOLPOS1A.19
CLL SOLPOS1A.20
CLL 3.4 20/06/94 DEF CAL360 replaced by LOGICAL LCAL360; GSS1F304.673
CLL PARAMETER statements duplicated for 360 and GSS1F304.674
CLL 365 day calendar. GSS1F304.675
CLL S.J.Swarbrick GSS1F304.676
!LL 4.4 27/02/97 Testing for leap years modified to deal with GMG1F404.337
!LL no leap every 100y except for every 400y GMG1F404.338
!LL Author: M.Gallani GMG1F404.339
CLL GSS1F304.677
CLL Programming standard : SOLPOS1A.21
CLL Written in FORTRAN 77, with the addition of "!" comments and SOLPOS1A.22
CLL underscores in variable names. SOLPOS1A.23
CLL Written to comply with 12/9/89 version of UMDP 4 (meteorological SOLPOS1A.24
CLL standard). SOLPOS1A.25
CLL SOLPOS1A.26
CLL Logical components covered : P233 SOLPOS1A.27
CLL SOLPOS1A.28
CLL Project task : SOLPOS1A.29
CLL SOLPOS1A.30
CLL External documentation: P23 SOLPOS1A.31
CLL SOLPOS1A.32
CLLEND ----------------------------------------------------------------- SOLPOS1A.33
C*L SOLPOS1A.34
SUBROUTINE SOLPOS (DAY, YEAR, SINDEC, SCS, LCAL360) 5GSS1F304.678
C GSS1F304.679
GSS1F304.680
LOGICAL LCAL360 !In, true if 360 day calendar in use. GSS1F304.681
C GSS1F304.682
INTEGER!, INTENT(IN) :: SOLPOS1A.36
& DAY, ! Day-number in the year SOLPOS1A.37
& YEAR ! Calendar year SOLPOS1A.38
REAL!, INTENT(OUT) :: SOLPOS1A.39
& SINDEC, ! sin(solar declination) SOLPOS1A.40
& SCS ! solar constant scaling SOLPOS1A.41
C* factor SOLPOS1A.42
CL This routine has no dynamically allocated work areas and no SOLPOS1A.43
CL significant structure. It calls the intrinsic functions FLOAT, SIN SOLPOS1A.44
CL & COS, but no user functions or subroutines. SOLPOS1A.45
CL SOLPOS1A.46
REAL GAMMA, E, TAU0, SINOBL, ! Basic orbital constants SOLPOS1A.47
& TAU1_360, TAU1_365, E1,E2,E3,E4,!Derived orbital contants GSS1F304.683
& TWOPI ! 2pi SOLPOS1A.49
REAL DINY_360, DINY_365 ! Number of days in year GSS1F304.684
REAL M, V ! Mean & true anomaly SOLPOS1A.51
*CALL C_PI 
SOLPOS1A.52
PARAMETER ( TWOPI = 2. * PI ) SOLPOS1A.53
PARAMETER (GAMMA=1.352631, E=.0167, ! Gamma, e SOLPOS1A.54
& TAU0 = 2.5, ! True date of perihelion SOLPOS1A.55
& SINOBL = .397789 ) ! Sin (obliquity) SOLPOS1A.56
PARAMETER ( E1 = E * (2.-.25*E*E), SOLPOS1A.57
& E2 = 1.25 * E*E, ! Coefficients for 3.1.2 SOLPOS1A.58
& E3 = E*E*E * 13./12., SOLPOS1A.59
& E4=( (1.+E*E*.5)/(1.-E*E) )**2 )! Constant for 3.1.4 SOLPOS1A.60
PARAMETER (DINY_360=360., TAU1_360=TAU0*DINY_360/365.25+0.71+.5) GSS1F304.685
C GSS1F304.686
PARAMETER (TAU1_365=TAU0+.5) GSS1F304.687
C GSS1F304.688
IF (.NOT. LCAL360) THEN GSS1F304.689
IF (mod(year,4) .eq. 0 .AND. ! is this a leap year? GMG1F404.340
& (mod(year,400) .eq. 0 .OR. mod(year,100) .ne. 0)) then GMG1F404.341
DINY_365 = 366. GSS1F304.691
ELSE GSS1F304.692
DINY_365 = 365. GSS1F304.693
END IF GSS1F304.694
END IF GSS1F304.695
C GSS1F304.696
! In forecast mode and in climate mode with real-year means, DINY GMG1F404.342
! depends on whether it is a leap year, otherwise DINY_36x = 360. GMG1F404.343
! GMG1F404.344
C TAU1 is modified so as to include the conversion of day-ordinal into SOLPOS1A.73
C fractional-number-of-days-into-the-year-at-12-Z-on-this-day. SOLPOS1A.74
C SOLPOS1A.75
IF (LCAL360) THEN GSS1F304.699
M = TWOPI * (FLOAT(DAY)-TAU1_360) / DINY_360 ! Eq 3.1.1 GSS1F304.700
ELSE GSS1F304.701
M = TWOPI * (FLOAT(DAY)-TAU1_365) / DINY_365 ! Eq 3.1.1 GSS1F304.702
END IF GSS1F304.703
V = M + E1*SIN(M) + E2*SIN(2.*M) + E3*SIN(3.*M) ! Eq 3.1.2 SOLPOS1A.77
SCS = E4 * ( 1. + E * COS(V) ) **2 ! Eq 3.1.4 SOLPOS1A.78
SINDEC = SINOBL * SIN (V - GAMMA) ! Eq 3.1.6 SOLPOS1A.79
RETURN SOLPOS1A.80
END SOLPOS1A.81
*ENDIF DEF,A01_1A,OR,DEF,A01_1B,OR,DEF,A01_2A,OR,DEF,A01_3A ADB1F400.385