*IF DEF,A02_1B LWPTSC1B.2
C ******************************COPYRIGHT****************************** GTS2F400.5653
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.5654
C GTS2F400.5655
C Use, duplication or disclosure of this code is subject to the GTS2F400.5656
C restrictions as set forth in the contract. GTS2F400.5657
C GTS2F400.5658
C Meteorological Office GTS2F400.5659
C London Road GTS2F400.5660
C BRACKNELL GTS2F400.5661
C Berkshire UK GTS2F400.5662
C RG12 2SZ GTS2F400.5663
C GTS2F400.5664
C If no contract has been raised with this copy of the code, the use, GTS2F400.5665
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.5666
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.5667
C Modelling at the above address. GTS2F400.5668
C ******************************COPYRIGHT****************************** GTS2F400.5669
C GTS2F400.5670
CLL SUBROUTINE LWPTSC LWPTSC1B.3
CLL LWPTSC1B.4
CLL PURPOSE LWPTSC1B.5
CLL It calculates scaled pathlengths of each gaseous absorber for each LWPTSC1B.6
CLL layer and returns them in DPATH for use by LWMAST, which sums them LWPTSC1B.7
CLL to get the total scaled pathlengths between each pair of layers so LWPTSC1B.8
CLL that the gaseous transmissivities can be calculated. LWPTSC1B.9
CLL Used in version 1B (gaseous effects treated as Morcrette et al, LWPTSC1B.10
CLL 1986) of the UM LW code. LWPTSC1B.11
CLL If UPDATE *DEF CRAY is off, a version is produced which except LWPTSC1B.12
CLL for the addition of ! comments is standard FORTRAN 77 (and which LWPTSC1B.13
CLL sets the "vector length" to 1) but the standard version includes LWPTSC1B.14
CLL CRAY automatic arrays also. LWPTSC1B.15
CLL Version 3, part of the alternative code giving ECMWF-like treatment LWPTSC1B.16
CLL of LW gaseous transmissivities. Almost all the calculations are LWPTSC1B.17
CLL changed: gases are generally scaled differently in each band they LWPTSC1B.18
CLL have an effect in, pressure scaling is no longer by fractional LWPTSC1B.19
CLL powers, quite elaborate pathlength-dependent temperature scaling LWPTSC1B.20
CLL is used, a diffusivity factor is included, zero pathlengths are LWPTSC1B.21
CLL permissible, and the indentation is changed. LWPTSC1B.22
CLL The increased complexity of the scaling means that the loop finding LWPTSC1B.24
CLL the water vapour scaled pathlengths does not (with current CRAY LWPTSC1B.25
CLL compilers) vectorize without compiler option "-o aggress". LWPTSC1B.26
CLL Version 3 of LWPTSC was set up from version 2.2 to be part of LWPTSC1B.28
CLL version 1B (ECMWF-like gaseous transmissivities) of the LW from LWPTSC1B.29
CLL release 2.7 of the UM. William Ingram 22 June 1992 LWPTSC1B.30
CLL LWPTSC1B.31
CLL Model Modification history from model version 3.0: LWPTSC1B.32
CLL version date LWPTSC1B.33
CLL 4.2 Sept.96 T3E migration: *DEF CRAY removed; GSS3F402.109
CLL dynamic allocation no longer *DEF controlled; GSS3F402.110
CLL cray HF functions removed. GSS3F402.111
CLL S.J.Swarbrick GSS3F402.112
CLL 4.3 Feb. 97 T3E optimisation: code restructured, cray vector GSS1F403.485
CLL library functions introduced. GSS1F403.486
CLL D.Salmond & S.J.Swarbrick GSS1F403.487
CLL 4.4 20/06/97 Add missing array indices to pbypr, dco2 ARB0F404.1
CLL & do3 in dry levels loop. RTHBarnes. ARB0F404.2
CLL LWPTSC1B.34
CLL LWPTSC1B.35
CLL It conforms to standard A of UMDP 4 (version 2, 18/1/90), and LWPTSC1B.36
CLL includes no 8X-deprecated features. LWPTSC1B.37
CLL LWPTSC1B.38
CLL It is part of component P232 (longwave radiation) which is in task LWPTSC1B.39
CLL P23 (radiation). LWPTSC1B.40
CLL LWPTSC1B.41
CLL External documentation is in UMDP 23. LWPTSC1B.42
C*L LWPTSC1B.43
SUBROUTINE LWPTSC (H2O, CO2, O3, PSTAR, AC, BC, AB, BB, TAC, 2LWPTSC1B.44
& L2, LWPTSC1B.46
& NWET, NOZONE, NLEVS, L1,DPATH) GSS3F402.113
C* LWPTSC1B.53
*CALL LWNGASES 
LWPTSC1B.54
C*L LWPTSC1B.58
INTEGER!, INTENT (IN) :: LWPTSC1B.59
& L2, ! Number of points to be treated GSS3F402.114
& NWET, ! Number of levels with moisture - LWPTSC1B.63
C ! above these zero is used. LWPTSC1B.64
& NOZONE, ! Number of levels with ozone data LWPTSC1B.65
C ! provided - below these the value in the lowest of them is used LWPTSC1B.66
& NLEVS, ! Number of levels LWPTSC1B.67
& L1 ! First dimension of input arrays LWPTSC1B.68
C ! (The different physical assumptions about water vapour and LWPTSC1B.69
C ! ozone in levels where no data is provided means that separate LWPTSC1B.70
C ! loops are used for levels with and without water vapour but LWPTSC1B.71
C ! only the indexing needs changing for levels with and without LWPTSC1B.72
C ! their own ozone data.) LWPTSC1B.73
REAL!, INTENT(IN) :: LWPTSC1B.74
& H2O(L1,NWET), CO2, ! Mass mixing ratio (mK in UMDP 23) LWPTSC1B.75
& O3(L1,NOZONE), ! of each absorbing gas LWPTSC1B.76
& TAC(L1,NLEVS), ! Mid-layer temperatures LWPTSC1B.77
& PSTAR(L1), ! Surface pressure LWPTSC1B.78
& AC(NLEVS), BC(NLEVS), ! A & B for layer centres and LWPTSC1B.79
& AB(NLEVS+1), BB(NLEVS+1) ! boundaries LWPTSC1B.80
REAL!, INTENT(OUT) :: LWPTSC1B.81
& DPATH(L2,NGASES,NLEVS) LWPTSC1B.82
C ! The scaled pathlengths are returned in DPATH, indexed by NGASES LWPTSC1B.83
C ! 1-6 are H2O line absorption in bands 1-6 respectively, 7 is CO2 LWPTSC1B.84
C ! scaled for band 2, 8 is CO2 scaled for bands 3 & 4, 9 & 10 are LWPTSC1B.85
C ! foreign & self-broadened water vapour continuum, and 11 and 12 LWPTSC1B.86
C ! are ozone without and with pressure scaling. LWPTSC1B.87
CL ! LWPTSC has no EXTERNAL calls and no significant structure LWPTSC1B.88
*CALL LWABTSAA LWPTSC1B.89
REAL RLNR10, ! lg(sqrt(e)) LWPTSC1B.90
& EPSP1, ! 1 + epsilon LWPTSC1B.91
& DIFFAC ! Diffusivity factor LWPTSC1B.92
C LWPTSC1B.93
REAL DABBYG, ! Difference of As & Bs across LWPTSC1B.94
& DBBBYG, ! model layer, divided by 10 g. LWPTSC1B.95
& DABMBP, ! Mean As & Bs across model layer, LWPTSC1B.96
& DBBMBP, ! divided by a reference pressure LWPTSC1B.97
C ! These four are used to calculate the next two quantities : LWPTSC1B.98
& DPBYGA, ! Pressure difference across model LWPTSC1B.99
C ! layer, divided by 10 g. LWPTSC1B.100
& pbypr(l2), ! Mean of layer-boundary pressure, GSS1F403.488
C ! divided by a reference pressure LWPTSC1B.102
C ! These two together give the pressure scaled pathlength - the LWPTSC1B.103
C ! integral across the layer with respect to pressure of the LWPTSC1B.104
C ! local pressure over the reference one. LWPTSC1B.105
& TN, ! Temperature diffce from TRTSAA LWPTSC1B.106
& dh2o(l2), ! Water vapour pathlength for a GSS1F403.489
C ! single layer, pressure-scaled LWPTSC1B.108
& UPH2O, ! Logarithmic function of DH2O LWPTSC1B.109
C ! used to calculate the temperature scaling LWPTSC1B.110
& dco2(l2), ! Equivalents of DH2O & UPH2O for GSS1F403.490
& UPCO2, ! CO2. LWPTSC1B.112
& do3(l2), ! Unscaled 1-layer ozone pathlength GSS1F403.491
& TSCCON, ! Temperature-scaling term for the LWPTSC1B.114
C ! self-broadened water vapour continuum pathlength LWPTSC1B.115
& SBWV, ! Water-vapour fraction of the LWPTSC1B.116
C ! atmosphere for calculating water vapour continuum pathlengths LWPTSC1B.117
& X, ! Dummy argument for statement GSS1F403.492
& exp_d(l2*10) GSS1F403.493
C ! functions LWPTSC1B.119
INTEGER LEVEL, J, ! Loopers over levels & points LWPTSC1B.120
& ONETWO, ! Flipper LWPTSC1B.121
& OLEVEL ! Index for the ozone data to be LWPTSC1B.122
C ! used in the current level LWPTSC1B.123
C* LWPTSC1B.124
*CALL C_G LWPTSC1B.125
*CALL C_EPSLON LWPTSC1B.126
PARAMETER ( EPSP1 = 1. + EPSILON ) LWPTSC1B.127
PARAMETER ( DIFFAC = 1.66 ) LWPTSC1B.128
C ! Simplify the code for temperature scaling of water vapour line LWPTSC1B.129
C ! and CO2 pathlengths by defining statement functions: LWPTSC1B.130
REAL POLYTS, LWPTSC1B.131
& TSCAL LWPTSC1B.132
POLYTS(X,J,ONETWO) = ABTSAA(1,J,ONETWO) + LWPTSC1B.133
& X * ( ABTSAA(2,J,ONETWO) + X * ABTSAA(3,J,ONETWO) ) LWPTSC1B.134
TSCAL(TN,X,J) = TN * ( POLYTS(X,J,1) + TN * POLYTS(X,J,2) ) LWPTSC1B.135
C ! FORTRAN 77 will not allow the following constants to be LWPTSC1B.136
C ! defined in a PARAMETER statement, but the CRAY compiler will LWPTSC1B.137
C ! give the same effect as if they were. LWPTSC1B.138
RLNR10 = .5 / LOG (10.) LWPTSC1B.139
C LWPTSC1B.140
DO 2 LEVEL=1, NWET LWPTSC1B.141
DABBYG = ( AB(LEVEL) - AB(LEVEL+1) ) / ( G * 10. ) LWPTSC1B.142
DBBBYG = ( BB(LEVEL) - BB(LEVEL+1) ) / ( G * 10. ) LWPTSC1B.143
DABMBP = ( AB(LEVEL) + AB(LEVEL+1) ) * 0.5 / 101325. LWPTSC1B.144
DBBMBP = ( BB(LEVEL) + BB(LEVEL+1) ) * 0.5 / 101325. LWPTSC1B.145
OLEVEL = MAX (1, LEVEL+NOZONE-NLEVS) LWPTSC1B.146
GSS1F403.494
GSS1F403.495
DO 20 J=1, L2 LWPTSC1B.149
DPBYGA = DABBYG + PSTAR(J) * DBBBYG LWPTSC1B.150
pbypr(j) = DABMBP + PSTAR(J) * DBBMBP GSS1F403.496
TN = TAC(J,LEVEL) - TRTSAA LWPTSC1B.152
C LWPTSC1B.153
dh2o(j) = DPBYGA * H2O(J,LEVEL) * pbypr(j) GSS1F403.497
z=1.0 GSS1F403.498
if(dh2o(j).eq.0) dh2o(j)=sign(dh2o(j),z) GSS1F403.499
GSS1F403.500
IF ( dh2o(j) .NE. 0. ) THEN GSS1F403.501
UPH2O = LWPTSC1B.156
& AMIN1 ( AMAX1( RLNR10 * alog(dh2o(j))+5.,0.), 6.0) GSS1F403.502
ELSE LWPTSC1B.162
UPH2O = 0. LWPTSC1B.163
ENDIF LWPTSC1B.164
dh2o(j) = dh2o(j) * DIFFAC GSS1F403.503
exp_d(j+0*l2) = TSCAL (TN, UPH2O, 1) GSS1F403.504
exp_d(j+1*l2) = TSCAL (TN, UPH2O, 2) GSS1F403.505
exp_d(j+2*l2) = TSCAL (TN, UPH2O, 3) GSS1F403.506
exp_d(j+3*l2) = TSCAL (TN, UPH2O, 4) GSS1F403.507
exp_d(j+4*l2) = TSCAL (TN, UPH2O, 5) GSS1F403.508
exp_d(j+5*l2) = TSCAL (TN, UPH2O, 6) GSS1F403.509
C LWPTSC1B.181
dco2(j) = DPBYGA * CO2 * pbypr(j) GSS1F403.510
UPCO2 = AMAX1 ( RLNR10 * alog ( dco2(j) ) + 5., 0.) GSS1F403.511
dco2(j) = dco2(j) * DIFFAC GSS1F403.512
exp_d(j+6*l2) = TSCAL (TN, UPCO2, 7) GSS1F403.513
exp_d(j+7*l2) = TSCAL (TN, UPCO2, 8) GSS1F403.514
C LWPTSC1B.194
TSCCON = EXP ( 6.08 * ( 296. / TAC(J,LEVEL) - 1. ) ) LWPTSC1B.198
SBWV = EPSP1 * H2O(J,LEVEL) / ( 1. + EPSILON*H2O(J,LEVEL) ) LWPTSC1B.200
DPATH(J,9,LEVEL) = (1.-SBWV) * dh2o(j) GSS1F403.515
DPATH(J,10,LEVEL) = SBWV * dh2o(j) * TSCCON GSS1F403.516
C LWPTSC1B.203
do3(j) = DIFFAC * DPBYGA * O3(J,OLEVEL) GSS1F403.517
exp_d(j+8*l2) = TN * ( O3T1 + TN * O3T2 ) GSS1F403.518
exp_d(j+9*l2) = TN * ( O3T3 + TN * O3T4 ) GSS1F403.519
20 CONTINUE LWPTSC1B.216
*IF DEF,VECTLIB PXVECTLB.97
call exp_v(l2*10,exp_d,exp_d) GSS1F403.521
*ELSE GSS1F403.522
do j=1,l2*10 GSS1F403.523
exp_d(j)=exp(exp_d(j)) GSS1F403.524
end do GSS1F403.525
*ENDIF GSS1F403.526
do j=1,l2 GSS1F403.527
DPATH(J,1,LEVEL) = DH2O(j) * exp_d(j+0*l2) GSS1F403.528
DPATH(J,2,LEVEL) = DH2O(j) * exp_d(j+1*l2) GSS1F403.529
enddo GSS1F403.530
do j=1,l2 GSS1F403.531
DPATH(J,3,LEVEL) = DH2O(j) * exp_d(j+2*l2) GSS1F403.532
DPATH(J,4,LEVEL) = DH2O(j) * exp_d(j+3*l2) GSS1F403.533
enddo GSS1F403.534
do j=1,l2 GSS1F403.535
DPATH(J,5,LEVEL) = DH2O(j) * exp_d(j+4*l2) GSS1F403.536
DPATH(J,6,LEVEL) = DH2O(j) * exp_d(j+5*l2) GSS1F403.537
enddo GSS1F403.538
do j=1,l2 GSS1F403.539
DPATH(J,7,LEVEL) = DCO2(j) * exp_d(j+6*l2) GSS1F403.540
DPATH(J,8,LEVEL) = DCO2(j) * exp_d(j+7*l2) GSS1F403.541
enddo GSS1F403.542
do j=1,l2 GSS1F403.543
DPATH(J,11,LEVEL) = DO3(j) * exp_d(j+8*l2) GSS1F403.544
DPATH(J,12,LEVEL) = PBYPR(j) * DO3(j)* exp_d(j+9*l2) GSS1F403.545
enddo GSS1F403.546
2 CONTINUE LWPTSC1B.217
C LWPTSC1B.218
C ! for the H2O pathlength but treat CO2 and O3 the same: LWPTSC1B.220
C LWPTSC1B.221
DO 3 LEVEL=NWET+1, NLEVS LWPTSC1B.222
DABBYG = ( AB(LEVEL) - AB(LEVEL+1) ) / ( G * 10. ) LWPTSC1B.223
DBBBYG = ( BB(LEVEL) - BB(LEVEL+1) ) / ( G * 10. ) LWPTSC1B.224
DABMBP = ( AB(LEVEL) + AB(LEVEL+1) ) * 0.5 / 101325. LWPTSC1B.225
DBBMBP = ( BB(LEVEL) + BB(LEVEL+1) ) * 0.5 / 101325. LWPTSC1B.226
OLEVEL = MAX (1, LEVEL+NOZONE-NLEVS) LWPTSC1B.227
DO 30 J=1, L2 LWPTSC1B.230
DPBYGA = DABBYG + PSTAR(J) * DBBBYG LWPTSC1B.231
PBYPR(J) = DABMBP + PSTAR(J) * DBBMBP ARB0F404.3
TN = TAC(J,LEVEL) - TRTSAA LWPTSC1B.233
C LWPTSC1B.234
DCO2(J) = DPBYGA * CO2 * PBYPR(J) ARB0F404.4
UPCO2 = AMAX1 ( RLNR10 * alog ( DCO2(j) ) + 5., 0.) GSS1F403.548
DCO2(J) = DCO2(J) * DIFFAC ARB0F404.5
exp_d(j+0*l2) = TSCAL (TN, UPCO2, 7) GSS1F403.550
exp_d(j+1*l2) = TSCAL (TN, UPCO2, 8) GSS1F403.551
C GSS1F403.552
DO3(J) = DIFFAC * DPBYGA * O3(J,OLEVEL) ARB0F404.6
exp_d(j+2*l2) = TN * ( O3T1 + TN * O3T2 ) GSS1F403.554
exp_d(j+3*l2) = TN * ( O3T3 + TN * O3T4 ) GSS1F403.555
30 CONTINUE GSS1F403.556
*IF DEF,VECTLIB PXVECTLB.98
call exp_v(l2*4,exp_d,exp_d) GSS1F403.558
*ELSE GSS1F403.559
do j=1,l2*4 GSS1F403.560
exp_d(j)=exp(exp_d(j)) GSS1F403.561
end do GSS1F403.562
*ENDIF GSS1F403.563
do j=1,l2 GSS1F403.564
DPATH(J,1,LEVEL) = 0. LWPTSC1B.235
DPATH(J,2,LEVEL) = 0. LWPTSC1B.236
enddo GSS1F403.565
do j=1,l2 GSS1F403.566
DPATH(J,3,LEVEL) = 0. LWPTSC1B.237
DPATH(J,4,LEVEL) = 0. LWPTSC1B.238
enddo GSS1F403.567
do j=1,l2 GSS1F403.568
DPATH(J,5,LEVEL) = 0. LWPTSC1B.239
DPATH(J,6,LEVEL) = 0. LWPTSC1B.240
enddo GSS1F403.569
do j=1,l2 GSS1F403.570
DPATH(J,7,LEVEL) = DCO2(j) * exp_d(j+0*l2) GSS1F403.571
DPATH(J,8,LEVEL) = DCO2(j) * exp_d(j+1*l2) GSS1F403.572
enddo GSS1F403.573
do j=1,l2 GSS1F403.574
DPATH(J,9,LEVEL) = 0. LWPTSC1B.255
DPATH(J,10,LEVEL) = 0. LWPTSC1B.256
enddo GSS1F403.575
do j=1,l2 GSS1F403.576
DPATH(J,11,LEVEL) = DO3(j) * exp_d(j+2*l2) GSS1F403.577
DPATH(J,12,LEVEL) = PBYPR(j) * DO3(j) * exp_d(j+3*l2) GSS1F403.578
enddo GSS1F403.579
3 CONTINUE LWPTSC1B.269
C LWPTSC1B.270
RETURN LWPTSC1B.271
END LWPTSC1B.272
*ENDIF A02_1B LWPTSC1B.273