*IF DEF,A70_1A,OR,DEF,A70_1B APB4F405.27
C ******************************COPYRIGHT****************************** GASCAL1A.3
C (c) CROWN COPYRIGHT 1996, METEOROLOGICAL OFFICE, All Rights Reserved. GASCAL1A.4
C GASCAL1A.5
C Use, duplication or disclosure of this code is subject to the GASCAL1A.6
C restrictions as set forth in the contract. GASCAL1A.7
C GASCAL1A.8
C Meteorological Office GASCAL1A.9
C London Road GASCAL1A.10
C BRACKNELL GASCAL1A.11
C Berkshire UK GASCAL1A.12
C RG12 2SZ GASCAL1A.13
C GASCAL1A.14
C If no contract has been raised with this copy of the code, the use, GASCAL1A.15
C duplication or disclosure of it is strictly prohibited. Permission GASCAL1A.16
C to do so must first be obtained in writing from the Head of Numerical GASCAL1A.17
C Modelling at the above address. GASCAL1A.18
C ******************************COPYRIGHT****************************** GASCAL1A.19
CLL Subroutine GAS_CALC ---------------------------------------------- GASCAL1A.20
CLL GASCAL1A.21
CLL Purpose : GASCAL1A.22
CLL Calculates the trace gas mixing ratio (or weighting factor for GASCAL1A.23
CLL aerosol forcing fields. Rates of increase (yearly compound factors) AWI1F403.358
CLL can be supplied, or spot values (which will be linearly GASCAL1A.25
CLL interpolated) or a mixture of these. It is designed so it can be GASCAL1A.26
CLL called each time step, but when rates of increase are being used, GASCAL1A.27
CLL values are in fact only updated at New Year. GASCAL1A.28
CLL The rules are: GASCAL1A.29
CLL If rates exist (i.e. are positive) for the first & current years AWI1F403.359
CLL then all concentrations are ignored, except for the initial value. GASCAL1A.31
CLL If there is a positive rate for the current year but not for the AWI1F403.360
CLL start, the current rate & most recent concentration are used. GASCAL1A.33
CLL If rates do not exist for the current year then the concentration GASCAL1A.34
CLL is calculated by linear interpolation between the concentrations at GASCAL1A.35
CLL the given years. GASCAL1A.36
CLL The mixing ratios calculated after the last given year use the GASCAL1A.37
CLL rate for the final given year. GASCAL1A.38
CLL The 360-day year is assumed. GASCAL1A.39
CLL CARE should be taken if solitary rates are specified, as this can GASCAL1A.40
CLL result in discontinuities in the concentration time profile at GASCAL1A.41
CLL the next given year without a corresponding given rate. GASCAL1A.42
CLL GASCAL1A.43
CLL Authors : Andrew Brady, Tim Johns, William Ingram GASCAL1A.44
CLL GASCAL1A.45
CLL Version for : Cray YMP GASCAL1A.46
CLL GASCAL1A.47
CLL Model GASCAL1A.48
CLL version Date GASCAL1A.49
CLL 4.2 19/11/96 Author: William Ingram, reviewer Cath Senior. GASCAL1A.50
CLL 4.4 22/9/97 Correct code for the case when one changes AWI1F404.1
CLL from linear interpolation to a rate & then changes the rate. WJI AWI1F404.2
CLL GASCAL1A.51
CLL GASCAL1A.52
CLLEND ----------------------------------------------------------------- GASCAL1A.53
C*L Arguments GASCAL1A.54
GASCAL1A.55
SUBROUTINE GAS_CALC(GAS_NOW 10,4GASCAL1A.56
& ,GAS_INDEX_MAX GASCAL1A.57
& ,GAS_YEAR GASCAL1A.58
& ,GAS_CONC GASCAL1A.59
& ,GAS_RATE GASCAL1A.60
& ,MAX_SCENARIO_PTS AWI1F403.354
& ,ICODE GASCAL1A.61
& ,CMESSAGE) GASCAL1A.62
GASCAL1A.63
IMPLICIT NONE GASCAL1A.64
GASCAL1A.65
*CALL CMAXSIZE 
! Maximum size of common block arrays GASCAL1A.66
*CALL CSUBMODL ! Submodel parameters for array sizes GASCAL1A.67
GASCAL1A.68
REAL GAS_NOW !OUT Gas concentration at time step GASCAL1A.69
INTEGER GAS_INDEX_MAX!IN GASCAL1A.70
& , MAX_SCENARIO_PTS ! IN AWI1F403.355
INTEGER GAS_YEAR(MAX_SCENARIO_PTS) !IN GASCAL1A.71
REAL GAS_CONC(MAX_SCENARIO_PTS) !IN GASCAL1A.72
REAL GAS_RATE(MAX_SCENARIO_PTS) !IN GASCAL1A.73
INTEGER ICODE !OUT Return code: successful=0 GASCAL1A.74
CHARACTER*(*) CMESSAGE !OUT Error message if ICODE >0 GASCAL1A.75
C* GASCAL1A.76
GASCAL1A.77
! Common blocks GASCAL1A.78
GASCAL1A.79
*CALL CTIME ! Model time GASCAL1A.80
GASCAL1A.81
! Local variables GASCAL1A.82
GASCAL1A.83
INTEGER INDEX ! to subscript gas concs for NOW_TIME GASCAL1A.84
INTEGER I ! Loop over indices GASCAL1A.85
INTEGER YEAR_IN_SECS! Year length in seconds GASCAL1A.86
INTEGER NOW_TIME_DAY, NOW_TIME_SEC GASCAL1A.87
! ! Time now in days/secs from time zero GASCAL1A.88
INTEGER GAS_YR_DAY1, GAS_YR_SEC1 GASCAL1A.89
! ! Time in days/secs of current GAS_YEAR GASCAL1A.90
INTEGER TIME1 GASCAL1A.91
! ! The same converted to seconds GASCAL1A.92
INTEGER GAS_YR_DAY2, GAS_YR_SEC2 GASCAL1A.93
! ! Time in days/secs of next GAS_YEAR GASCAL1A.94
GASCAL1A.95
! Check that GASCNST namelist is defined for this year GASCAL1A.96
IF ( I_YEAR .LT. GAS_YEAR(1) ) THEN GASCAL1A.97
ICODE = 8325 GASCAL1A.98
CMESSAGE = 'GAS_CALC: no gas data for this year' GASCAL1A.99
RETURN GASCAL1A.100
ENDIF GASCAL1A.101
GASCAL1A.102
! Loop over I to find correct index for current NOW_TIME GASCAL1A.103
INDEX = 0 GASCAL1A.104
DO I=1, GAS_INDEX_MAX GASCAL1A.105
IF ( I_YEAR .GE. GAS_YEAR(I) ) INDEX = INDEX+1 GASCAL1A.106
ENDDO GASCAL1A.107
GASCAL1A.108
! Calculate time now in seconds GASCAL1A.109
CALL TIME2SEC (I_YEAR, I_MONTH, I_DAY, I_HOUR, I_MINUTE, I_SECOND, GASCAL1A.110
& 0, 0, NOW_TIME_DAY, NOW_TIME_SEC, .TRUE.) GASCAL1A.111
GASCAL1A.112
! If gas rate at current year is non zero calculate new GAS_NOW GASCAL1A.113
! by considering compound increases of GAS_RATE(1:INDEX) GASCAL1A.114
IF ( GAS_RATE(INDEX) .GT. 0. ) THEN AWI1F403.356
YEAR_IN_SECS = 360 * 86400 GASCAL1A.116
CALL TIME2SEC (GAS_YEAR(INDEX), 1, 1, 0, 0, 0, GASCAL1A.117
& 0, 0, GAS_YR_DAY1, GAS_YR_SEC1, .TRUE.) AWI1F404.3
GAS_NOW = GAS_CONC(1) AWI1F404.4
DO I=1, INDEX-1 AWI1F404.5
IF ( GAS_RATE(I) .LT. 0. ) THEN AWI1F404.6
GAS_NOW = GAS_CONC(I+1) AWI1F404.7
ELSE AWI1F404.8
GAS_NOW = GAS_NOW * AWI1F404.9
& ( GAS_RATE(I) ** REAL(GAS_YEAR(I+1)-GAS_YEAR(I)) ) AWI1F404.10
ENDIF AWI1F404.11
ENDDO AWI1F404.12
! GAS_NOW now holds the concentration in year INDEX - need only GASCAL1A.130
! update it to the current year. GASCAL1A.131
GAS_NOW=GAS_NOW*(GAS_RATE(INDEX)** GASCAL1A.132
& REAL(((NOW_TIME_DAY-GAS_YR_DAY1)*86400+ GASCAL1A.134
& NOW_TIME_SEC-GAS_YR_SEC1)/YEAR_IN_SECS)) GASCAL1A.135
GASCAL1A.136
! Otherwise calculate by linear interpolation between respective GASCAL1A.137
! GAS concentrations of given years. GASCAL1A.138
ELSE GASCAL1A.139
CALL TIME2SEC (GAS_YEAR(INDEX), 1, 1, 0, 0, 0, GASCAL1A.140
& 0, 0, GAS_YR_DAY1, GAS_YR_SEC1, .TRUE.) GASCAL1A.141
CALL TIME2SEC (GAS_YEAR(INDEX+1), 1, 1, 0, 0, 0, GASCAL1A.142
& 0, 0, GAS_YR_DAY2, GAS_YR_SEC2, .TRUE.) GASCAL1A.143
TIME1 = GAS_YR_DAY1*86400 - GAS_YR_SEC1 GASCAL1A.144
GAS_NOW = GAS_CONC(INDEX) + GASCAL1A.145
& ( GAS_CONC(INDEX+1) - GAS_CONC(INDEX) ) GASCAL1A.146
& * REAL ( NOW_TIME_DAY*86400 + NOW_TIME_SEC - TIME1 ) GASCAL1A.147
& / REAL ( GAS_YR_DAY2*86400 + GAS_YR_SEC2 - TIME1 ) GASCAL1A.148
ENDIF GASCAL1A.149
GASCAL1A.150
RETURN GASCAL1A.151
END GASCAL1A.152
*ENDIF DEF,A70_1A,OR,DEF,A70_1B APB4F405.28