!+ Subroutine to set thermodynamic properties                                      RAD_ALL.3218  
!                                                                                  RAD_ALL.3219  
! Purpose:                                                                         RAD_ALL.3220  
!   Pressures, temperatures at the centres and edges of layers                     RAD_ALL.3221  
!   and the masses in layers are set.                                              RAD_ALL.3222  
!                                                                                  RAD_ALL.3223  
! Method:                                                                          RAD_ALL.3224  
!   Straightforward.                                                               RAD_ALL.3225  
!                                                                                  RAD_ALL.3226  
! Current Owner of Code: J. M. Edwards                                             RAD_ALL.3227  
!                                                                                  RAD_ALL.3228  
! History:                                                                         RAD_ALL.3229  
!       Version         Date                    Comment                            RAD_ALL.3230  
!       4.0             27-07-95                Original Code                      RAD_ALL.3231  
!                                               (J. M. Edwards)                    RAD_ALL.3232  
!       4.1             10-06-96                Old formulation over               RAD_ALL.3233  
!                                               sea-ice removed.                   RAD_ALL.3234  
!                                               (J. M. Edwards)                    RAD_ALL.3235  
!       4.2             08-08-96                Ground temperature                 RAD_ALL.3236  
!                                               set equal to that                  RAD_ALL.3237  
!                                               in the middle of the               RAD_ALL.3238  
!                                               bottom layer.                      RAD_ALL.3239  
!                                               (J. M. Edwards)                    RAD_ALL.3240  
!                                                                                  RAD_ALL.3241  
! Description of Code:                                                             RAD_ALL.3242  
!   FORTRAN 77  with extensions listed in documentation.                           RAD_ALL.3243  
!                                                                                  RAD_ALL.3244  
!- ---------------------------------------------------------------------           RAD_ALL.3245  
      SUBROUTINE R2_SET_THERMODYNAMIC(                                             RAD_ALL.3246  
     &     N_PROFILE, NLEVS, I_GATHER, L_BOUNDARY_TEMPERATURE                      RAD_ALL.3247  
     &   , PSTAR, TSTAR, AB, BB, AC, BC, PEXNER, TAC                               RAD_ALL.3248  
     &   , P, T, T_BDY, T_SURFACE, D_MASS                                          RAD_ALL.3249  
     &   , NPD_FIELD, NPD_PROFILE, NPD_LAYER                                       RAD_ALL.3250  
     &   )                                                                         RAD_ALL.3251  
!                                                                                  RAD_ALL.3252  
!                                                                                  RAD_ALL.3253  
!                                                                                  RAD_ALL.3254  
      IMPLICIT NONE                                                                RAD_ALL.3255  
!                                                                                  RAD_ALL.3256  
!                                                                                  RAD_ALL.3257  
!     INCLUDED COMDECKS                                                            RAD_ALL.3258  
C*L------------------COMDECK C_G----------------------------------------           C_G.2     
C G IS MEAN ACCEL DUE TO GRAVITY AT EARTH'S SURFACE                                C_G.3     
      REAL G                                                                       C_G.4     
                                                                                   C_G.5     
      PARAMETER(G=9.80665)                                                         C_G.6     
C*----------------------------------------------------------------------           C_G.7     
C*L------------------COMDECK C_R_CP-------------------------------------           C_R_CP.2     
C R IS GAS CONSTANT FOR DRY AIR                                                    C_R_CP.3     
C CP IS SPECIFIC HEAT OF DRY AIR AT CONSTANT PRESSURE                              C_R_CP.4     
C PREF IS REFERENCE SURFACE PRESSURE                                               C_R_CP.5     
      REAL R,CP,KAPPA,PREF  ! PREFTOKI                                             C_R_CP.6     
                                                                                   C_R_CP.7     
      PARAMETER(R=287.05,                                                          C_R_CP.8     
     &          CP=1005.,                                                          C_R_CP.9     
     &          KAPPA=R/CP,                                                        C_R_CP.10    
     &          PREF=100000.)                                                      C_R_CP.11    
C*----------------------------------------------------------------------           C_R_CP.12    
                                                                                   C_R_CP.13    
!     ------------------------------------------------------------------           PRMCH3A.2     
!     COMDECK FOR TWO-STREAM RADIATION CODE.                                       PRMCH3A.3     
!                                                                                  PRMCH3A.4     
!     MODULE FOR SETTING MACHINE PRECISION.                                        PRMCH3A.5     
!                                                                                  PRMCH3A.6     
      REAL                                                                         PRMCH3A.7     
     &     TOL_MACHINE                                                             PRMCH3A.8     
!             MACHINE TOLERANCE                                                    PRMCH3A.9     
     &   , SQRT_TOL_MACHINE                                                        PRMCH3A.10    
!             SQRT OF MACHINE TOLERANCE                                            PRMCH3A.11    
!                                                                                  PRMCH3A.12    
!                                                                                  PRMCH3A.13    
!     THE PRECISION SHOULD BE ABOUT 2/2^(SIZE OF SIGNIFICAND)                      PRMCH3A.14    
!                                                                                  PRMCH3A.15    
!     THE IEEE-FORMAT USES 53 BITS FOR THE SIGNIFICAND                             PRMCH3A.16    
!     IN DOUBLE PRECISION                                                          PRMCH3A.17    
!                                                                                  PRMCH3A.18    
!     THE CRAY FORMAT USES 47 BITS IN SINGLE PRECISION.                            PRMCH3A.19    
!                                                                                  PRMCH3A.20    
      PARAMETER(                                                                   PRMCH3A.21    
     &     TOL_MACHINE=1.42E-14                                                    PRMCH3A.22    
     &   , SQRT_TOL_MACHINE=1.19E-7                                                PRMCH3A.23    
     &   )                                                                         PRMCH3A.24    
!                                                                                  PRMCH3A.25    
!     ------------------------------------------------------------------           PRMCH3A.26    
!     ------------------------------------------------------------------           PRECSN3A.2     
!     COMDECK FOR TWO-STREAM RADIATION CODE.                                       PRECSN3A.3     
!                                                                                  PRECSN3A.4     
!     MODULE FOR SETTING MACHINE-DEPENDENT TOLERANCES.                             PRECSN3A.5     
!     (THE COMDECK PRMCH3A MUST ALWAYS BE INCLUDED BEFORE THIS COMDECK.)           PRECSN3A.6     
!                                                                                  PRECSN3A.7     
      REAL                                                                         PRECSN3A.8     
     &     TOL_DIV                                                                 PRECSN3A.9     
!             TOLERANCE FOR DIVISION                                               PRECSN3A.10    
     &   , TOL_TEST                                                                PRECSN3A.11    
!             TOLERANCE FOR TESTING EQUALITY                                       PRECSN3A.12    
!                                                                                  PRECSN3A.13    
      PARAMETER(                                                                   PRECSN3A.14    
     &     TOL_DIV=3.2E+01*TOL_MACHINE                                             PRECSN3A.15    
     &   , TOL_TEST=1.6E+01*TOL_MACHINE                                            PRECSN3A.16    
     &   )                                                                         PRECSN3A.17    
!                                                                                  PRECSN3A.18    
!     ------------------------------------------------------------------           PRECSN3A.19    
!                                                                                  RAD_ALL.3263  
!     DUMMY ARGUMENTS.                                                             RAD_ALL.3264  
!     SIZES OF ARRAYS:                                                             RAD_ALL.3265  
      INTEGER   !, INTENT(IN)                                                      RAD_ALL.3266  
     &     NPD_FIELD                                                               RAD_ALL.3267  
!             SIZE OF ARRAY FROM UM                                                RAD_ALL.3268  
     &   , NPD_PROFILE                                                             RAD_ALL.3269  
!             MAXIMUM NUMBER OF PROFILES                                           RAD_ALL.3270  
     &   , NPD_LAYER                                                               RAD_ALL.3271  
!             MAXIMUM NUMBER OF LAYERS                                             RAD_ALL.3272  
!                                                                                  RAD_ALL.3273  
!     SIZES USED:                                                                  RAD_ALL.3274  
      INTEGER   !, INTENT(IN)                                                      RAD_ALL.3275  
     &     N_PROFILE                                                               RAD_ALL.3276  
!             NUMBER OF PROFILES                                                   RAD_ALL.3277  
     &   , NLEVS                                                                   RAD_ALL.3278  
!             NUMBER OF LEVELS                                                     RAD_ALL.3279  
     &   , I_GATHER(NPD_FIELD)                                                     RAD_ALL.3280  
!             LIST OF POINTS TO GATHER                                             RAD_ALL.3281  
!                                                                                  RAD_ALL.3282  
      REAL      !, INTENT(IN)                                                      RAD_ALL.3283  
     &     PSTAR(NPD_FIELD)                                                        RAD_ALL.3284  
!             SURFACE PRESSURES                                                    RAD_ALL.3285  
     &   , TSTAR(NPD_FIELD)                                                        RAD_ALL.3286  
!             SURFACE TEMPERSTURES                                                 RAD_ALL.3287  
     &   , AB(NLEVS+1)                                                             RAD_ALL.3288  
!             A AT EDGES OF LAYERS                                                 RAD_ALL.3289  
     &   , BB(NLEVS+1)                                                             RAD_ALL.3290  
!             B AT EDGES OF LAYERS                                                 RAD_ALL.3291  
     &   , AC(NLEVS)                                                               RAD_ALL.3292  
!             A AT CENTRES OF LAYERS                                               RAD_ALL.3293  
     &   , BC(NLEVS)                                                               RAD_ALL.3294  
!             B AT CENTRES OF LAYERS                                               RAD_ALL.3295  
     &   , TAC(NPD_FIELD, NLEVS)                                                   RAD_ALL.3296  
!             TEMPERATURES AT CENTRES OF LAYERS                                    RAD_ALL.3297  
     &   , PEXNER(NPD_FIELD, NLEVS+1)                                              RAD_ALL.3298  
!             EXNER FUNCTION AT BOUNDARIES                                         RAD_ALL.3299  
!                                                                                  RAD_ALL.3300  
      LOGICAL   !, INTENT(IN)                                                      RAD_ALL.3301  
     &     L_BOUNDARY_TEMPERATURE                                                  RAD_ALL.3302  
!             FLAG TO CALCULATE TEMPERATURES AT BOUNADRIES OF LAYERS.              RAD_ALL.3303  
!                                                                                  RAD_ALL.3304  
!                                                                                  RAD_ALL.3305  
      REAL      !, INTENT(OUT)                                                     RAD_ALL.3306  
     &     D_MASS(NPD_PROFILE, NPD_LAYER)                                          RAD_ALL.3307  
!             MASS THICKNESSES OF LAYERS                                           RAD_ALL.3308  
     &   , P(NPD_PROFILE, 0: NPD_LAYER)                                            RAD_ALL.3309  
!             PRESSURE FIELD                                                       RAD_ALL.3310  
     &   , T(NPD_PROFILE, 0: NPD_LAYER)                                            RAD_ALL.3311  
!             TEMPERATURE FIELD                                                    RAD_ALL.3312  
     &   , T_BDY(NPD_PROFILE, 0: NPD_LAYER)                                        RAD_ALL.3313  
!             TEMPERATURES AT EDGES OF LAYERS                                      RAD_ALL.3314  
     &   , T_SURFACE(NPD_PROFILE)                                                  RAD_ALL.3315  
!             GATHERED TEMPERATURE OF SURFACE                                      RAD_ALL.3316  
!                                                                                  RAD_ALL.3317  
!                                                                                  RAD_ALL.3318  
!     LOCAL VARIABLES.                                                             RAD_ALL.3319  
!                                                                                  RAD_ALL.3320  
      INTEGER                                                                      RAD_ALL.3321  
     &     I                                                                       RAD_ALL.3322  
!             LOOP VARIABLE                                                        RAD_ALL.3323  
     &   , II                                                                      RAD_ALL.3324  
!             LOOP VARIABLE                                                        RAD_ALL.3325  
     &   , L                                                                       RAD_ALL.3326  
!             LOOP VARIABLE                                                        RAD_ALL.3327  
     &   , LG                                                                      RAD_ALL.3328  
!             INDEX TO GATHER                                                      RAD_ALL.3329  
!                                                                                  RAD_ALL.3330  
      REAL                                                                         RAD_ALL.3331  
     &     PU                                                                      RAD_ALL.3332  
!             PRESSURE FOR UPPER LAYER                                             RAD_ALL.3333  
     &   , PL                                                                      RAD_ALL.3334  
!             PRESSURE FOR LOWER LAYER                                             RAD_ALL.3335  
     &   , PML1                                                                    RAD_ALL.3336  
!             PRESSURE FOR INTERPOLATION                                           RAD_ALL.3337  
     &   , WTL                                                                     RAD_ALL.3338  
!             WEIGHT FOR LOWER LAYER                                               RAD_ALL.3339  
     &   , WTU                                                                     RAD_ALL.3340  
!             WEIGHT FOR UPPER LAYER                                               RAD_ALL.3341  
!                                                                                  RAD_ALL.3342  
C*L------------------ COMDECK P_EXNERC ---------------------------------           P_EXNERC.2     
C statement function to define exner pressure at model levels                      P_EXNERC.3     
C from exner pressures at model half-levels                                        P_EXNERC.4     
      REAL P_EXNER_C                                                               P_EXNERC.5     
      REAL R_P_EXNER_C                                                             P_EXNERC.6     
      REAL P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM                             P_EXNERC.7     
C consistent with geopotential see eqn 26 DOC PAPER 10                             P_EXNERC.15    
      P_EXNER_C(P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM) =                     P_EXNERC.16    
     & (P_EXU_DUM*PU_DUM - P_EXL_DUM*PL_DUM)/                                      P_EXNERC.17    
     & ( (PU_DUM-PL_DUM)*(KAPPA_DUM + 1) )                                         P_EXNERC.18    
      R_P_EXNER_C(P_EXU_DUM,P_EXL_DUM,PU_DUM,PL_DUM,KAPPA_DUM) =                   P_EXNERC.19    
     & ( (PU_DUM-PL_DUM)*(KAPPA_DUM + 1) )/                                        P_EXNERC.20    
     & (P_EXU_DUM*PU_DUM - P_EXL_DUM*PL_DUM)                                       P_EXNERC.21    
                                                                                   P_EXNERC.23    
C*------------------- --------------------------------------------------           P_EXNERC.24    
                                                                                   P_EXNERC.25    
!                                                                                  RAD_ALL.3344  
!                                                                                  RAD_ALL.3345  
!                                                                                  RAD_ALL.3346  
!     CALCULATE PROPERTIES AT THE CENTRES OF LAYERS.                               RAD_ALL.3347  
      DO I=1, NLEVS                                                                RAD_ALL.3348  
         DO L=1, N_PROFILE                                                         RAD_ALL.3349  
            LG=I_GATHER(L)                                                         RAD_ALL.3350  
            P(L, I)=AC(NLEVS+1-I)+BC(NLEVS+1-I)*PSTAR(LG)                          RAD_ALL.3351  
            T(L, I)=TAC(LG, NLEVS+1-I)                                             RAD_ALL.3352  
            D_MASS(L, I)=(AB(NLEVS+1-I)-AB(NLEVS+2-I)                              RAD_ALL.3353  
     &         +PSTAR(LG)*(BB(NLEVS+1-I)-BB(NLEVS+2-I)))                           RAD_ALL.3354  
     &         /G                                                                  RAD_ALL.3355  
         ENDDO                                                                     RAD_ALL.3356  
      ENDDO                                                                        RAD_ALL.3357  
!                                                                                  RAD_ALL.3358  
!                                                                                  RAD_ALL.3359  
      IF (L_BOUNDARY_TEMPERATURE) THEN                                             RAD_ALL.3360  
!                                                                                  RAD_ALL.3361  
!        GATHER THE SURFACE TEMPERATURE.                                           RAD_ALL.3362  
         DO L=1, N_PROFILE                                                         RAD_ALL.3363  
            LG=I_GATHER(L)                                                         RAD_ALL.3364  
            T_SURFACE(L)=TSTAR(LG)                                                 RAD_ALL.3365  
         ENDDO                                                                     RAD_ALL.3366  
!                                                                                  RAD_ALL.3367  
!        INTERPOLATE TEMPERATURES AT THE BOUNDARIES OF LAYERS                      RAD_ALL.3368  
!        FROM THE EXNER FUNCTION.                                                  RAD_ALL.3369  
         DO L=1, N_PROFILE                                                         RAD_ALL.3370  
            LG=I_GATHER(L)                                                         RAD_ALL.3371  
!                                                                                  RAD_ALL.3372  
!           TAKE THE TEMPERATURE OF THE AIR JUST ABOVE THE SURFACE AS              RAD_ALL.3373  
!           THE TEMPERATURE AT THE MIDDLE OF THE BOTTOM LAYER.                     RAD_ALL.3374  
            T_BDY(L, NLEVS)=TAC(LG, 1)                                             RAD_ALL.3375  
!           TAKE THE TEMPERATURE AS CONSTANT ACROSS THE TOP HALF-LAYER.            RAD_ALL.3376  
            T_BDY(L, 0)=TAC(LG, NLEVS)                                             RAD_ALL.3377  
!                                                                                  RAD_ALL.3378  
         ENDDO                                                                     RAD_ALL.3379  
!                                                                                  RAD_ALL.3380  
         DO I=1, NLEVS-1                                                           RAD_ALL.3381  
            II=NLEVS-I                                                             RAD_ALL.3382  
            DO L=1, N_PROFILE                                                      RAD_ALL.3383  
               LG=I_GATHER(L)                                                      RAD_ALL.3384  
               PU=PSTAR(LG)*BB(II+2)+AB(II+2)                                      RAD_ALL.3385  
               PL=PSTAR(LG)*BB(II+1)+AB(II+1)                                      RAD_ALL.3386  
               PML1=PSTAR(LG)*BB(II)+AB(II)                                        RAD_ALL.3387  
               WTU=TAC(LG, II+1)*(PEXNER(LG, II+1)                                 RAD_ALL.3388  
     &            /P_EXNER_C(PEXNER(LG, II+2), PEXNER(LG, II+1)                    RAD_ALL.3389  
     &            , PU, PL, KAPPA)-1.0E+00)                                        RAD_ALL.3390  
               WTL=TAC(LG, II)*(PEXNER(LG, II)                                     RAD_ALL.3391  
     &            /P_EXNER_C(PEXNER(LG, II+1), PEXNER(LG, II)                      RAD_ALL.3392  
     &            , PL, PML1, KAPPA)-1.0E+00)                                      RAD_ALL.3393  
               T_BDY(L, I)=(WTU*TAC(LG, NLEVS+1-I)                                 RAD_ALL.3394  
     &            +WTL*TAC(LG, NLEVS-I))/(WTL+WTU)                                 RAD_ALL.3395  
            ENDDO                                                                  RAD_ALL.3396  
         ENDDO                                                                     RAD_ALL.3397  
!                                                                                  RAD_ALL.3398  
      ENDIF                                                                        RAD_ALL.3399  
!                                                                                  RAD_ALL.3400  
!                                                                                  RAD_ALL.3401  
!                                                                                  RAD_ALL.3402  
      RETURN                                                                       RAD_ALL.3403  
      END   
