*IF DEF,A15_1A DYNDIA1A.2
C ******************************COPYRIGHT****************************** GTS2F400.2359
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.2360
C GTS2F400.2361
C Use, duplication or disclosure of this code is subject to the GTS2F400.2362
C restrictions as set forth in the contract. GTS2F400.2363
C GTS2F400.2364
C Meteorological Office GTS2F400.2365
C London Road GTS2F400.2366
C BRACKNELL GTS2F400.2367
C Berkshire UK GTS2F400.2368
C RG12 2SZ GTS2F400.2369
C GTS2F400.2370
C If no contract has been raised with this copy of the code, the use, GTS2F400.2371
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.2372
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.2373
C Modelling at the above address. GTS2F400.2374
C ******************************COPYRIGHT****************************** GTS2F400.2375
C GTS2F400.2376
CLL SUBROUTINE DYN_DIAG--------------------------------------- DYNDIA1A.3
CLL DYNDIA1A.4
CLL PURPOSE: Calculate various diagnostics required for operational DYNDIA1A.5
CLL and climate oputput. U and V compnts of wind on P levels DYNDIA1A.6
CLL , clear air turbulence, DYNDIA1A.7
CLL and MAXIMUM WIND Compnts and LEVEL. DYNDIA1A.8
CLL and Potential Vorticty on isentropic surfaces DYNDIA1A.9
CLL DYNDIA1A.10
CLL JH, DR, RS <- programmer of some or all of previous code or changes DYNDIA1A.11
CLL DYNDIA1A.12
CLL Model Modification history from model version 3.0: DYNDIA1A.13
CLL version Date DYNDIA1A.14
CLL 3.0 30/12/92 Two references to the former deckname DYNDIAG1 DYNDIA1A.15
CLL in CMESSAGE changed to subroutine name DYN_DIAG. MJH DYNDIA1A.16
CLL 3.1 25/01/93 Include test diagnostic (a simple analytic RR250193.1
CLL function), items 231,232,233,234. R. Rawlins RR250193.2
CLL 3.1 14/01/93 Add routines to calculate potential vorticity on MM180193.1
CLL a pressure surface and theta on a pv surface. MM180193.2
CLL 3.4 26/05/94 Argument LLINTS added and passed to CALC_PV, GSS1F304.186
CLL CALC_PV_P, THETA_PV S.J.Swarbrick GSS1F304.187
CLL 4.1 31/05/96 The number of v points to be processed on a UIE2F401.408
CLL C grid differs from u by row_length. u,v UIE2F401.409
CLL dimensioned separately in call to WLLTOEQ. UIE2F401.410
CLL Requirement for VAR. UIE2F401.411
CLL Author I.Edmond Reviewer D. Goddard UIE2F401.412
!LL 4.2 08/01/97 Initialise PUV array to remove any NaNs in GPB1F403.260
!LL halo regions. P.Burton GPB1F403.261
!LL 4.4 09/04/97 : Add new diagnositics 235 qw, 236 heavyside ARS1F404.1
!LL function (1 if pressure surface above land ARS1F404.2
!LL zero if below) ARS1F404.3
!LL and 237 total kinetic energy per unit area. ARS1F404.4
!LL 30/07/97 : Also geopotential height on u grid, 238 Z, ARS1F404.5
!LL 239 uZ, 240 VZ. ARS1F404.6
!LL 19/08/97 : 241 mountain torque per unit area. ARS1F404.7
!LL R. A. Stratton. ARS1F404.8
!LL 4.5 15/04/98 Added start-end arguments to V_INT, V_INT_T and GSM1F405.615
!LL V_INT_Z routines, and also to a lot of loops GSM1F405.616
!LL over fields. Consequently, NS halos of diagnostics GSM1F405.617
!LL are not set in this routine - instead STASHWORK is GSM1F405.618
!LL initialised in ST_DIAG1. S.D.Mullerworth GSM1F405.619
!LL 23/09/98 Allow 50m winds to be above second model level GPB0F405.162
!LL P.Burton GPB0F405.163
CLL DYNDIA1A.17
CLL Programming standard: U M DOC Paper NO. 4, DYNDIA1A.18
CLL DYNDIA1A.19
CLL System components covered : D16D DYNDIA1A.20
CLL DYNDIA1A.21
CLL Project task: DYNDIA1A.22
CLL DYNDIA1A.23
CLL External documentation: DYNDIA1A.24
!LL Unified Model Documentation Paper no D4 ARS1F404.9
!LL describes the diagnostics u*v, v*T etc. ARS1F404.10
CLL DYNDIA1A.25
CLLEND------------------------------------------------------------- DYNDIA1A.26
DYNDIA1A.27
C DYNDIA1A.28
C*L ARGUMENTS:--------------------------------------------------- DYNDIA1A.29
SUBROUTINE DYN_DIAG( 2,40DYNDIA1A.30
*CALL ARGFLDPT
GSM1F405.620
C primary data in DYNDIA1A.31
& PSTAR,U,V,Q, DYNDIA1A.32
& THETA,OROG,P_EXNER_HALF,PSTAR_OLD, DYNDIA1A.33
C primary data constants DYNDIA1A.34
& U_ROWS,P_ROWS,ROW_LENGTH,P_LEVELS,Q_LEVELS,P_FIELD, DYNDIA1A.35
& U_FIELD,AK,BK,AKH,BKH,DELTA_AK,DELTA_BK, DYNDIA1A.36
& NMOST_LAT,WMOST_LONG,NS_SPACE,EW_SPACE,PHI_POLE, DYNDIA1A.37
& LAMBDA_POLE,SEC_U_LATITUDE,ROTATE_UV,ROTATE_MAX_UV, DYNDIA1A.38
& ELF,ETA_MATRIX_INV,MATRIX_P_O,LATITUDE_STEP_INVERSE, RR250193.3
& LONGITUDE_STEP_INVERSE,ADVECTION_TIMESTEP,SEC_P_LATITUDE, DYNDIA1A.40
& COS_U_LATITUDE,F3,FORECAST_HRS, RR250193.4
C Required Theta values DYNDIA1A.43
& DESIRED_THETA,PV_PRESS,DESIRED_PV,REQ_THETA_PV_LEVS, MM180193.3
& n_levels, TD141293.101
C Required pressures RR250193.5
& UCOMP_PRESS,VCOMP_PRESS,CAT_PROB_PRESS,T_PRESS,W_PRESS, DYNDIA1A.45
& Q_PRESS,TESTD_PRESS,HEAVY_PRESS,Z_PRESS, ARS1F404.11
C Required Model levels RR250193.7
& TESTD_MODEL, RR250193.8
C Indices for product fields DYNDIA1A.47
& UV_IND,UT_IND,VT_IND,T2_IND,U2_IND,V2_IND,WT_IND,WU_IND, DYNDIA1A.48
& WV_IND,QU_IND,QV_IND,QW_IND,UZ_IND,VZ_IND, ARS1F404.12
C DIAGNOSTICS OUT DYNDIA1A.50
& UCOMP_P,VCOMP_P,MAX_CAT_PROB,MAX_CAT_LEVEL, DYNDIA1A.51
& CAT_PROB_SINGLE,MAX_WIND_HEIGHT, DYNDIA1A.52
& MAX_WIND_ICAO_HEIGHT,MAX_WIND_PRESSURE ,UCOMP_MAX_WIND, DYNDIA1A.53
& VCOMP_MAX_WIND,CAT_PROB_MEAN,UCOMP50_WIND,VCOMP50_WIND, DYNDIA1A.54
& POTN_VORT_THETA, MM180193.4
& UV_P,T_P,UT_P,VT_P,T2_P,U2_P,V2_P,W_P,WT_P,WU_P,WV_P,Q_P, DYNDIA1A.56
& UQ_P,VQ_P, DYNDIA1A.57
& POTN_VORT_ON_P,THETA_ON_PV, MM180193.5
& TESTDIAG1,TESTDIAG2,TESTDIAG3,TESTDIAG4, RR250193.9
& WQ_P,HEAVYSIDE_P,TOTAL_KE, ARS1F404.13
& Z_P,UZ_P,VZ_P,M_TORQUE, ARS1F404.14
C diagnostic lengths DYNDIA1A.58
& UCOMP_P_LEVS,VCOMP_P_LEVS,CAT_PROB_LEVS,POTN_VORT_THETA_LEVS, MM180193.6
& POTN_VORT_P_LEVS,THETA_PV_LEVS,THETA_PV_P_LEVS, MM180193.7
& UV_P_LEVS,T_P_LEVS, DYNDIA1A.60
& UT_P_LEVS,VT_P_LEVS,T2_P_LEVS,U2_P_LEVS,V2_P_LEVS,W_P_LEVS, DYNDIA1A.61
& WT_P_LEVS,WU_P_LEVS,WV_P_LEVS,Q_P_LEVS,QU_P_LEVS,QV_P_LEVS, DYNDIA1A.62
& TESTD_P_LEVS,TESTD_M_LEVS, QW_P_LEVS,HEAVY_P_LEVS, ARS1F404.15
& Z_P_LEVS,UZ_P_LEVS,VZ_P_LEVS, ARS1F404.16
C diagnostic logical indicators DYNDIA1A.63
& QUCOMP_P, QVCOMP_P, QMAX_CAT_PROB, DYNDIA1A.64
& QMAX_CAT_LEVEL, QCAT_PROB_SINGLE, QMAX_WIND_HEIGHT, DYNDIA1A.65
& QMAX_WIND_ICAO_HEIGHT, QMAX_WIND_PRESSURE,QUCOMP_MAX_WIND, DYNDIA1A.66
& QVCOMP_MAX_WIND, QCAT_PROB_MEAN, QUCOMP50_WIND, DYNDIA1A.67
& QVCOMP50_WIND, QPOTN_VORT_THETA, MM180193.8
& QUV_P, QT_P, QUT_P, QVT_P, QT2_P, MM180193.9
& QU2_P, QV2_P, QW_P, QWT_P, QWU_P, QWV_P, QQ_P, QUQ_P, QVQ_P, DYNDIA1A.69
& QPOTN_VORT_PRESS, QTHETA_ON_PV, MM180193.10
& QDIA1,QDIA2,QDIA3,QDIA4, RR250193.11
& QWQ_P,QHEAVY_P,QTOTAL_KE,QZ_P,QUZ_P,QVZ_P,Q_MT,Z_REF, ARS1F404.17
C diagnostic rerun code and message RR250193.12
& ICODE,CMESSAGE, GSS1F304.188
C Logical switch LLINTS - passed to other routines GSS1F304.189
& LLINTS) GSS1F304.190
C DYNDIA1A.72
IMPLICIT NONE DYNDIA1A.73
LOGICAL LLINTS GSS1F304.191
*CALL C_R_CP
DYNDIA1A.74
*CALL C_G
DYNDIA1A.75
*CALL C_MDI
DYNDIA1A.76
*CALL C_PI
DYNDIA1A.77
*IF DEF,MPP GSM3F403.1
*CALL PARVARS
GSM3F403.2
*ENDIF GSM3F403.3
*CALL TYPFLDPT
GSM3F403.4
INTEGER DYNDIA1A.78
* P_FIELD !IN 1ST DIMENSION OF FIELD OF PSTAR DYNDIA1A.79
*, U_FIELD !IN 1ST DIMENSION OF FIELD OF U,V DYNDIA1A.80
*, U_ROWS !IN NUMBER OF ROWS FOR U,V FIELD DYNDIA1A.81
*, P_ROWS !IN NUMBER OF ROWS FOR P,T FIELD DYNDIA1A.82
*, ROW_LENGTH !IN NUMBER OF POINTS PER ROW DYNDIA1A.83
*, LEVELS !IN NUMBER OF MODEL LEVELS DYNDIA1A.84
*, P_LEVELS !IN NUMBER OF PRESSURE LEVELS DYNDIA1A.85
*, Q_LEVELS !IN NUMBER OF WET LEVELS DYNDIA1A.86
*, MATRIX_P_O !IN Order of polynomial used in calculation RR250193.13
* ! of ETA_HALF inverse matrix RR250193.14
*, FORECAST_HRS !IN FORECAST HOURS AFTER ANALYSIS T+nn RR250193.15
*, ICODE ! RETURN CODE : IRET=0 NORMAL EXIT DYNDIA1A.89
INTEGER DYNDIA1A.90
* UCOMP_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP U_P DYNDIA1A.91
*, VCOMP_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP V_P DYNDIA1A.92
*, CAT_PROB_LEVS !IN NO OF LEVS ON WHICH TO CALC/INTERP CAT DYNDIA1A.93
& ,n_levels ! number of levels for dtheta/dp spline TD141293.102
&, POTN_VORT_THETA_LEVS !IN no of theta levs on which to calc pv MM180193.11
&, POTN_VORT_P_LEVS !IN no of p levs on which to calc pv MM180193.12
&, THETA_PV_LEVS !IN no of pv levs on which to calc theta MM180193.13
&, THETA_PV_P_LEVS !IN no of p levs to calculate pv at,which MM180193.14
& ! are then used to calculate theta on p. MM180193.15
*, UV_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP UV_P DYNDIA1A.95
*, T_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP T_P DYNDIA1A.96
*, UT_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP UT_P DYNDIA1A.97
*, VT_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP VT_P DYNDIA1A.98
*, T2_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP T2_P DYNDIA1A.99
*, U2_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP U2_P DYNDIA1A.100
*, V2_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP V2_P DYNDIA1A.101
*, W_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP W_P DYNDIA1A.102
*, WT_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP WT_P DYNDIA1A.103
*, WU_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP WU_P DYNDIA1A.104
*, WV_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP WV_P DYNDIA1A.105
*, Q_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP Q_P DYNDIA1A.106
*, QU_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP UQ_P DYNDIA1A.107
*, QV_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP VQ_P DYNDIA1A.108
*, QW_P_LEVS !IN NO OF LEVS ON WHICH TO INTERP WQ_P ARS1F404.18
*, TESTD_P_LEVS !IN NO OF PRESS LEVELS TO CALC TESTDIAG RR250193.16
*, TESTD_M_LEVS !IN NO OF MODEL LEVELS TO CALC TESTDIAG RR250193.17
&, HEAVY_P_LEVS !IN NO OF PRESS levels to calc HEAVYSIDE ARS1F404.19
&, Z_P_LEVS !IN No of press levels for geopotential ARS1F404.20
&, UZ_P_LEVS !IN No of press levels for UZ ARS1F404.21
&, VZ_P_LEVS !IN No of press levels for VZ ARS1F404.22
&, Z_REF !IN Level of model used to calculate PMSL ARS1F404.23
ARS1F404.24
INTEGER DYNDIA1A.109
* UV_IND(UV_P_LEVS,2) !IN index for pressure levels in u and v DYNDIA1A.110
*, UT_IND(UT_P_LEVS,2) !IN index for pressure levels in U and T DYNDIA1A.111
*, VT_IND(VT_P_LEVS,2) !IN index for pressure levels in V and T DYNDIA1A.112
*, T2_IND(T2_P_LEVS) !IN index for pressure levels in T DYNDIA1A.113
*, U2_IND(U2_P_LEVS) !IN index for pressure levels in U DYNDIA1A.114
*, V2_IND(V2_P_LEVS) !IN index for pressure levels in V DYNDIA1A.115
*, WT_IND(WT_P_LEVS,2) !IN index for pressure levels in W and T DYNDIA1A.116
*, WU_IND(WU_P_LEVS,2) !IN index for pressure levels in W and U DYNDIA1A.117
*, WV_IND(WV_P_LEVS,2) !IN index for pressure levels in W and V DYNDIA1A.118
*, QU_IND(QU_P_LEVS,2) !IN index for pressure levels in q and U DYNDIA1A.119
*, QV_IND(QV_P_LEVS,2) !IN index for pressure levels in Q and V DYNDIA1A.120
*, QW_IND(QW_P_LEVS,2) !IN index for pressure levels in Q and W ARS1F404.25
*, UZ_IND(UZ_P_LEVS,2) !IN index for pressure levels in u and z ARS1F404.26
*, VZ_IND(VZ_P_LEVS,2) !IN index for pressure levels in v and z ARS1F404.27
C DYNDIA1A.121
LOGICAL DYNDIA1A.122
* QUCOMP_P !IN LOGICAL FLAG FOR PRESS INTER U COMPONENTS DYNDIA1A.123
*,QVCOMP_P !IN " " " " " V " DYNDIA1A.124
*,QMAX_CAT_PROB !IN " " " MAXIMUM CAT PROBABILITY DYNDIA1A.125
*,QMAX_CAT_LEVEL !IN " " " MAX CAT PROB LEVEL DYNDIA1A.126
*,QCAT_PROB_SINGLE !IN " " " CAT PROB ON PRESS SFCE DYNDIA1A.127
&,QMAX_WIND_HEIGHT !IN " " " MAX_WIND HEIGHT DYNDIA1A.128
&,QMAX_WIND_ICAO_HEIGHT !IN " " MAX_WIND ICAO HEIGHT DYNDIA1A.129
&,QMAX_WIND_PRESSURE !IN " " MAX_WIND PRESSURE DYNDIA1A.130
&,QUCOMP_MAX_WIND ! IN " " " UCOMP_MAX _WIND DYNDIA1A.131
&,QVCOMP_MAX_WIND ! IN " " " VCOMP_MAX _WIND DYNDIA1A.132
*,QCAT_PROB_MEAN ! IN " " " CAT PROB MEAN 300/250/200m DYNDIA1A.133
&,ELF ! IN True if ELF i.e rotated LAM grid RR250193.18
&,ROTATE_UV ! IN True if winds to be rotated DYNDIA1A.134
&,ROTATE_MAX_UV ! IN True if winds to be rotated DYNDIA1A.135
&,QUCOMP50_WIND ! IN Logical flag for QUCOMP50 DYNDIA1A.136
&,QVCOMP50_WIND ! IN " " " " DYNDIA1A.137
&,QPOTN_VORT_THETA ! IN " " " computing pv on theta MM180193.16
&,QPOTN_VORT_PRESS ! IN " " " computing pv on MM180193.17
& ! pressure. MM180193.18
&,QTHETA_ON_PV ! IN " " " computing theta on PV MM180193.19
&,QUV_P ! IN " " " UV on pressure levels DYNDIA1A.139
&,QT_P ! IN " " " T on pressure levels DYNDIA1A.140
&,QUT_P ! IN " " " UT on pressure levels DYNDIA1A.141
&,QVT_P ! IN " " " VT on pressure levels DYNDIA1A.142
&,QT2_P ! IN " " " T2 on pressure levels DYNDIA1A.143
&,QU2_P ! IN " " " U2 on pressure levels DYNDIA1A.144
&,QV2_P ! IN " " " V2 on pressure levels DYNDIA1A.145
&,Qw_P ! IN " " " w on pressure levels DYNDIA1A.146
&,QwT_P ! IN " " " wT on pressure levels DYNDIA1A.147
&,QwU_P ! IN " " " wU on pressure levels DYNDIA1A.148
&,QwV_P ! IN " " " wV on pressure levels DYNDIA1A.149
&,QQ_P ! IN " " " q on pressure levels DYNDIA1A.150
&,Quq_P ! IN " " " uq on pressure levels DYNDIA1A.151
&,Qvq_P ! IN " " " vq on pressure levels DYNDIA1A.152
&,Qwq_P ! IN " " " wq on pressure levels ARS1F404.28
&,QDIA1 ! IN " " " test diagnostic 1 RR250193.19
&,QDIA2 ! IN " " " test diagnostic 2 RR250193.20
&,QDIA3 ! IN " " " test diagnostic 3 RR250193.21
&,QDIA4 ! IN " " " test diagnostic 4 RR250193.22
&,QHEAVY_P ! IN " " " heavyside function p lev ARS1F404.29
&,QTOTAL_KE ! IN " " " Total KE ARS1F404.30
&,QZ_P ! IN " " " Z on pressure levels ARS1F404.31
&,QUZ_P ! IN " " " UZ on pressure levels ARS1F404.32
&,QVZ_P ! IN " " " VZ on pressure levels ARS1F404.33
&,Q_MT ! IN mountain torque per unit area ARS1F404.34
C DYNDIA1A.153
CHARACTER CMESSAGE*(*) DYNDIA1A.154
DYNDIA1A.155
REAL DYNDIA1A.156
* PSTAR(P_FIELD) !IN PRIMARY MODEL ARRAY FOR PSTAR FIELD DYNDIA1A.157
*,PSTAR_OLD(P_FIELD) !IN Pstar before dynamics. DYNDIA1A.158
*,OROG(P_FIELD) !IN PRIMARY MODEL OROGRAPHY DYNDIA1A.159
*,P_EXNER_HALF(P_FIELD,P_LEVELS+1) !IN EXNER PRESS ON 1/2 LVLS DYNDIA1A.160
*,THETA(P_FIELD,P_LEVELS)!IN PRIMARY MODEL ARRAY FOR THETA FIELD DYNDIA1A.161
*,U(U_FIELD,P_LEVELS) !INT PRIMARY MODEL ARRAY FOR U FIELD DYNDIA1A.162
*,V(U_FIELD,P_LEVELS) !IN PRIMARY MODEL ARRAY FOR V FIELD DYNDIA1A.163
*,Q(P_FIELD,Q_LEVELS) !IN PRIMARY MODEL ARRAY FOR HUMIDITY DYNDIA1A.164
C DIAGNOSTICS OUT RR250193.23
*,UCOMP_P(U_FIELD,UCOMP_P_LEVS) !OUT UCOMP ON ANY PRESSURE SFCE DYNDIA1A.165
*,VCOMP_P(U_FIELD,VCOMP_P_LEVS) !OUT VCOMP ON ANY PRESSURE SFCE DYNDIA1A.166
*,UCOMP50_WIND(U_FIELD) !OUT 50 M wind zonal cmpnt. DYNDIA1A.167
*,VCOMP50_WIND(U_FIELD) !OUT 50 M wind zonal cmpnt. DYNDIA1A.168
*,MAX_CAT_PROB(U_FIELD) !OUT MAX CAT PROB FROM LEVELS 300/250/200mb DYNDIA1A.169
*,MAX_CAT_LEVEL(U_FIELD)!OUT LEVEL OF MAX CAT PROB DYNDIA1A.170
*,CAT_PROB_SINGLE(U_FIELD,CAT_PROB_LEVS)!OUT CAT PROB ON PRESS SFCE DYNDIA1A.171
*,MAX_WIND_HEIGHT(U_FIELD) !OUT HEIGHT LEVEL OF MAX WIND DYNDIA1A.172
*,MAX_WIND_ICAO_HEIGHT(U_FIELD) !OUT ICAO HEIGHT LEVEL OF MAX WIND DYNDIA1A.173
*,MAX_WIND_PRESSURE(U_FIELD) !OUT PRESSURE LEVEL OF MAX WIND DYNDIA1A.174
*,UCOMP_MAX_WIND(U_FIELD) !OUT U COMPONENT OF MAX WIND DYNDIA1A.175
*,VCOMP_MAX_WIND(U_FIELD) !OUT V COMPONENT OF MAX WIND DYNDIA1A.176
*,CAT_PROB_MEAN(U_FIELD)!OUT CAT PROB MEAN OVER LEVELS 300/250/200m DYNDIA1A.177
&,POTN_VORT_THETA(P_FIELD,POTN_VORT_THETA_LEVS) !OUT pv on theta MM180193.20
&,POTN_VORT_ON_P(P_FIELD,POTN_VORT_P_LEVS) !OUT pv on pressure MM180193.21
&,THETA_ON_PV(P_FIELD,THETA_PV_LEVS) !OUT Pot. temp on a pv surface MM180193.22
&,e_levels(n_levels) ! Model half-levels for dtheta/dp. TD141293.103
&,dthe_dph(p_field,n_levels) ! dtheta/dp for potential vorticity TD141293.104
& ! on half-levels. TD141293.105
*,UV_P(U_FIELD,UV_P_LEVS) ! UV on pressure levels, wind grid DYNDIA1A.179
*,T_P(U_FIELD,T_P_LEVS) ! T on pressure levels, wind grid DYNDIA1A.180
*,UT_P(U_FIELD,UT_P_LEVS) ! UT on pressure levels, wind grid DYNDIA1A.181
*,VT_P(U_FIELD,VT_P_LEVS) ! VT on pressure levels, wind grid DYNDIA1A.182
*,T2_P(U_FIELD,T2_P_LEVS) ! T2 on pressure levels, wind grid DYNDIA1A.183
*,U2_P(U_FIELD,U2_P_LEVS) ! U2 on pressure levels, wind grid DYNDIA1A.184
*,V2_P(U_FIELD,V2_P_LEVS) ! V2 on pressure levels, wind grid DYNDIA1A.185
*,W_P(U_FIELD,W_P_LEVS) ! w on pressure levels, wind grid DYNDIA1A.186
*,WT_P(U_FIELD,WT_P_LEVS) ! wT on pressure levels, wind grid DYNDIA1A.187
*,WU_P(U_FIELD,WU_P_LEVS) ! wU on pressure levels, wind grid DYNDIA1A.188
*,WV_P(U_FIELD,WV_P_LEVS) ! wV on pressure levels, wind grid DYNDIA1A.189
*,Q_P(U_FIELD,Q_P_LEVS) ! q on pressure levels, wind grid DYNDIA1A.190
*,UQ_P(U_FIELD,QU_P_LEVS) ! qu on pressure levels, wind grid DYNDIA1A.191
*,VQ_P(U_FIELD,QV_P_LEVS) ! qv on pressure levels, wind grid DYNDIA1A.192
*,WQ_P(U_FIELD,QW_P_LEVS) ! qw on pressure levels, wind grid ARS1F404.35
*,TESTDIAG1(U_FIELD) ! OUT Diag 1 single lev, u grid RR250193.24
*,TESTDIAG2(P_FIELD) ! OUT Diag 2 single lev, p grid RR250193.25
*,TESTDIAG3(P_FIELD,TESTD_P_LEVS)! OUT Diag 3 press levs, p grid RR250193.26
*,TESTDIAG4(P_FIELD,TESTD_M_LEVS)! OUT Diag 4 model levs, p grid RR250193.27
&,HEAVYSIDE_P(U_FIELD,HEAVY_P_LEVS) ! OUT heavyside on p levs ARS1F404.36
&,TOTAL_KE(U_FIELD) ! total KE per unit area, u grid ARS1F404.37
&,Z_P(U_FIELD,Z_P_LEVS) ! z on pressure levels, u grid ARS1F404.38
&,UZ_P(U_FIELD,UZ_P_LEVS) ! Uz on pressure levels, u grid ARS1F404.39
&,VZ_P(U_FIELD,VZ_P_LEVS) ! Vz on pressure levels, u grid ARS1F404.40
&,M_TORQUE(U_FIELD) ! mountain torque per unit area, u ARS1F404.41
C DYNDIA1A.193
C AK,BK DEFINE HYBRID VERTICAL COORDINATES P=A+BP*, DYNDIA1A.194
C DELTA_AK,DELTA_BK DEFINE LAYER PRESSURE THICKNESS PD=AD+BDP*, DYNDIA1A.195
REAL DYNDIA1A.196
* AKH(P_LEVELS+1) !IN LAYER THICKNESS DYNDIA1A.197
*,BKH(P_LEVELS+1) !IN LAYER THICKNESS DYNDIA1A.198
*,AK (P_LEVELS) !IN VALUE AT LAYER CENTRE DYNDIA1A.199
*,BK (P_LEVELS) !IN VALUE AT LAYER CENTRE DYNDIA1A.200
*,DELTA_AK (P_LEVELS) !IN DYNDIA1A.201
*,DELTA_BK (P_LEVELS) !IN DYNDIA1A.202
*,NMOST_LAT !Northern most latitude of grid DYNDIA1A.203
*,WMOST_LONG !Western most longitude DYNDIA1A.204
*,EW_SPACE !Delta longitude DYNDIA1A.205
*,NS_SPACE !Delta latitude DYNDIA1A.206
*,PHI_POLE !Latitude of the pseudo pole DYNDIA1A.207
*,LAMBDA_POLE !Longitude of the pseudo pole DYNDIA1A.208
*,SEC_U_LATITUDE(U_FIELD)!IN 1/COS(LAT) AT U POINTS DYNDIA1A.209
*,ETA_MATRIX_INV(MATRIX_P_O,MATRIX_P_O,P_LEVELS)!IN Inverse matrix DYNDIA1A.210
* ! of ETA_HALF DYNDIA1A.211
*,UCOMP_PRESS(UCOMP_P_LEVS) !IN Required pressure surface DYNDIA1A.212
*,VCOMP_PRESS(VCOMP_P_LEVS) !IN Required pressure surface DYNDIA1A.213
*,CAT_PROB_PRESS(CAT_PROB_LEVS)!IN " " " DYNDIA1A.214
&,DESIRED_THETA(POTN_VORT_THETA_LEVS) !IN required theta surfaces MM180193.23
&,PV_PRESS(POTN_VORT_P_LEVS) !IN required pressure surfaces MM180193.24
&,DESIRED_PV(THETA_PV_LEVS) !IN required pv surfaces MM180193.25
&,REQ_THETA_PV_LEVS(THETA_PV_P_LEVS) !IN required p surfaces. MM180193.26
*,T_PRESS(T_P_LEVS) !IN " " " DYNDIA1A.216
*,W_PRESS(W_P_LEVS) !IN " " " DYNDIA1A.217
*,Q_PRESS(Q_P_LEVS) !IN " " " DYNDIA1A.218
*,TESTD_PRESS(TESTD_P_LEVS) !IN Required pressures for test diag RR250193.28
&,HEAVY_PRESS(HEAVY_P_LEVS) !IN Required pressures for heavyside ARS1F404.42
&,Z_PRESS(Z_P_LEVS) !IN Required pressures for geopot ARS1F404.43
*,TESTD_MODEL(TESTD_M_LEVS) !IN Required model lvs for test diag RR250193.29
*,LATITUDE_STEP_INVERSE !IN 1/latitude increment DYNDIA1A.219
*,LONGITUDE_STEP_INVERSE !IN 1/longitude increment DYNDIA1A.220
*,ADVECTION_TIMESTEP !IN advection timestep DYNDIA1A.221
*,SEC_P_LATITUDE(P_FIELD) !IN 1/cos(lat) p points DYNDIA1A.222
*,COS_U_LATITUDE(U_FIELD) !IN cos(lat) u points DYNDIA1A.223
&,F3(U_FIELD) !IN Coriolis term. DYNDIA1A.224
DYNDIA1A.225
C Local variables DYNDIA1A.226
DYNDIA1A.227
LOGICAL GPB0F405.164
& found_levels ! TRUE if level search is successful GPB0F405.165
GPB0F405.166
INTEGER GPB0F405.167
& level1,level2 ! Model levels either side of 50m GPB0F405.168
GPB0F405.169
C*--------------------------------------------------------------------- DYNDIA1A.228
DYNDIA1A.229
C*L WORKSPACE USAGE:------------------------------------------------- DYNDIA1A.230
C DEFINE LOCAL WORKSPACE ARRAYS: DYNDIA1A.231
C REAL ARRAYS REQUIRED AT FULL FIELD LENGTH DYNDIA1A.232
C 1 INTEGER INDEX ARRAY DYNDIA1A.233
C DYNDIA1A.234
C*--------------------------------------------------------------------- DYNDIA1A.235
C DYNDIA1A.236
C*L EXTERNAL SUBROUTINES CALLED--------------------------------------- DYNDIA1A.237
EXTERNAL TROP,V_INT,WINDMAX,ICAO_HT,V_INT_T,OMEGA_DIAG, DYNDIA1A.238
& P_TO_UV,CAT,CALC_PV,CALC_PV_P,THETA_PV MM180193.27
& ,V_INT_ZH,V_INT_Z ARS1F404.44
C*------------------------------------------------------------------ DYNDIA1A.240
CL MAXIMUM VECTOR LENGTH ASSUMED IS (ROWS+1) * ROWLENGTH DYNDIA1A.241
CL--------------------------------------------------------------------- DYNDIA1A.242
C---------------------------------------------------------------------- DYNDIA1A.243
C DEFINE LOCAL VARIABLES DYNDIA1A.244
LOGICAL DYNDIA1A.245
* TEST ! NUMBER OF P POINTS NEEDED DYNDIA1A.246
REAL DYNDIA1A.247
* PHI_STAR(P_FIELD) DYNDIA1A.248
*, P(P_FIELD,P_LEVELS) ! PRESSURE ARRAY DYNDIA1A.249
*, PUV(U_FIELD,P_LEVELS) ! PRESSURE ARRAY ON U,V POINTS DYNDIA1A.250
*, PZ(P_FIELD) ! PRESSURE SURFACE ON WHICH RESULTS REQD DYNDIA1A.251
*, WORK1(U_FIELD) ! Work array DYNDIA1A.252
*, WORK5(P_FIELD) ! Work array DYNDIA1A.253
*, ETA1,ETA2,ETA50,C1,C2 ! Used in the calculation of 50 M winds DYNDIA1A.254
*, OMEGA(P_FIELD,P_LEVELS) ! Omega array DYNDIA1A.255
&, THETA_ON_P(P_FIELD,POTN_VORT_P_LEVS) !holds Pot. temperature on MM180193.28
& !on a pressure surface. MM180193.29
& ,model_half_height(p_field,p_levels+1) !heights on model half lev ARS1F404.45
& ,PSTAR_UV(U_FIELD) ! pstar on uv grid. ARS1F404.46
& ,FACTOR ! factor for KE calculation ARS1F404.47
& ,PLEV ! pressure level for Heavyside calculation ARS1F404.48
*IF -DEF,GLOBAL DYNDIA1A.256
*, WORK2(U_FIELD) ! Work array DYNDIA1A.257
*, WORK3(U_FIELD) ! Work array DYNDIA1A.258
*, WORK4(U_FIELD) ! Work array DYNDIA1A.259
*, COEFF1(U_FIELD) ! Rotation coeffs DYNDIA1A.260
*, COEFF2(U_FIELD) ! Rotation coeffs DYNDIA1A.261
*ENDIF DYNDIA1A.262
C DYNDIA1A.263
C R IS GAS CONSTANT FOR DRY AIR DYNDIA1A.264
C CP IS SPECIFIC HEAT OF DRY AIR AT CONSTANT PRESSURE DYNDIA1A.265
C PREF IS REFERENCE SURFACE PRESSURE DYNDIA1A.266
INTEGER K,I,II,IK,J,LEVEL! LOOP COUNTERS IN ROUTINE DYNDIA1A.267
*, T_REF ! reference level for below surface T extrapolation. DYNDIA1A.268
&, U_FLD_VALID ! Set to No of points in U-field excluding GSM1F405.621
! unused rows and halos GSM1F405.622
C DYNDIA1A.269
ICODE=0 DYNDIA1A.271
U_FLD_VALID=LAST_U_FLD_PT-FIRST_FLD_PT+1 GSM1F405.623
*IF -DEF,GLOBAL DYNDIA1A.272
CL------------------Calculate rotation coeffs if required ------------- DYNDIA1A.273
IF(QVCOMP_P.OR.QUCOMP_P.OR.QMAX_WIND_PRESSURE) THEN DYNDIA1A.274
IF(ROTATE_UV.OR.ROTATE_MAX_UV) THEN DYNDIA1A.275
K=0 DYNDIA1A.278
DO I=1,U_ROWS DYNDIA1A.279
DO J=1,ROW_LENGTH RR250193.30
K=K+1 RR250193.31
*IF DEF,MPP GSM3F403.6
C Calculate relative to NW corner of LAM grid GSM3F403.7
WORK1(K)=NMOST_LAT - GSM3F403.8
& (I-1-NS_Halo+datastart(2)-0.5)*NS_SPACE ! Equatorial Phi GSM3F403.9
WORK2(K)=WMOST_LONG+ GSM3F403.10
& (J-1-EW_Halo+datastart(1)-0.5)*EW_SPACE GSM3F403.11
*ELSE GSM3F403.12
WORK1(K)=NMOST_LAT - (I-0.5)*NS_SPACE ! Equatorial Phi RR250193.32
WORK2(K)=WMOST_LONG+ (J-0.5)*EW_SPACE RR250193.33
*ENDIF GSM3F403.13
IF(WORK2(K).GT.360.0) THEN RR250193.34
WORK2(K)=WORK2(K)-360.0 RR250193.35
ENDIF RR250193.36
ENDDO RR250193.37
ENDDO DYNDIA1A.288
C WORK3 holds true PHI and WORK4 holds true LAMBDA DYNDIA1A.289
CALL EQTOLL
(WORK1,WORK2,WORK3,WORK4,PHI_POLE,LAMBDA_POLE, DYNDIA1A.290
* U_FIELD) DYNDIA1A.291
DO I=1,U_FIELD DYNDIA1A.292
IF(WORK4(I).GT.180.0) WORK4(I)=WORK4(I)-360.0 DYNDIA1A.293
ENDDO DYNDIA1A.294
CALL W_COEFF
(COEFF1,COEFF2,WORK4,WORK2,PHI_POLE,LAMBDA_POLE, DYNDIA1A.295
* U_FIELD) DYNDIA1A.296
ENDIF DYNDIA1A.297
ENDIF DYNDIA1A.298
*ENDIF DYNDIA1A.299
DYNDIA1A.300
C----------------------------------------------------------------------- DYNDIA1A.301
CL Calculate variables required by various subroutines DYNDIA1A.302
C----------------------------------------------------------------------- DYNDIA1A.303
IF(QUCOMP_P.OR.QVCOMP_P.OR.QCAT_PROB_SINGLE.OR.QCAT_PROB_MEAN.OR. DYNDIA1A.304
& QMAX_CAT_PROB.OR.QMAX_CAT_LEVEL.OR.QT_P.OR.QW_P)THEN DYNDIA1A.305
DO K=1,P_LEVELS DYNDIA1A.306
DO I=1,P_FIELD DYNDIA1A.307
P(I,K)=AK(K)+BK(K)*PSTAR(I) DYNDIA1A.308
ENDDO DYNDIA1A.309
ENDDO DYNDIA1A.310
*IF DEF,MPP GPB1F403.262
! QAN fix : Initialise unused rows GSM1F405.624
IF (at_base_of_LPG) THEN GSM1F405.625
DO K=1,P_LEVELS GSM1F405.626
DO I=LAST_U_VALID_PT,U_FIELD GSM1F405.627
PUV(I,K)=0.0 GSM1F405.628
ENDDO GSM1F405.629
ENDDO GPB1F403.267
ENDIF GSM1F405.630
*ENDIF GPB1F403.269
DO K=1,P_LEVELS DYNDIA1A.311
CALL P_TO_UV
(P(1,K),PUV(1,K),P_FIELD,U_FIELD,ROW_LENGTH, DYNDIA1A.312
& P_ROWS) DYNDIA1A.313
ENDDO DYNDIA1A.314
*IF DEF,MPP GSM3F403.14
CALL SWAPBOUNDS
(PUV,ROW_LENGTH,U_ROWS, GSM3F403.15
& EW_Halo,NS_Halo,P_LEVELS) GSM3F403.16
*ENDIF GSM3F403.17
ENDIF DYNDIA1A.315
C----------------------------------------------------------------------- DYNDIA1A.316
IF(QUCOMP_P.OR.QVCOMP_P.OR.QZ_P) THEN ARS1F404.49
DO I=1,P_FIELD DYNDIA1A.318
PHI_STAR(I)=OROG(I)*G DYNDIA1A.319
ENDDO DYNDIA1A.320
*IF -DEF,GLOBAL DYNDIA1A.321
CL------------------Test to see if winds are to be rotated------------- DYNDIA1A.322
IF(ROTATE_UV) THEN DYNDIA1A.323
TEST=.TRUE. DYNDIA1A.324
IF(UCOMP_P_LEVS.EQ.VCOMP_P_LEVS) THEN DYNDIA1A.325
DO K=1,UCOMP_P_LEVS DYNDIA1A.326
IF(UCOMP_PRESS(K).NE.VCOMP_PRESS(K)) TEST=.FALSE. DYNDIA1A.327
ENDDO DYNDIA1A.328
ELSE DYNDIA1A.329
TEST=.FALSE. DYNDIA1A.330
ENDIF DYNDIA1A.331
IF(TEST) THEN DYNDIA1A.332
DO K=1,UCOMP_P_LEVS DYNDIA1A.333
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.631
PZ(I)=UCOMP_PRESS(K)*100.0 DYNDIA1A.335
ENDDO DYNDIA1A.336
CALL V_INT
(PUV,PZ,U,WORK1,U_FIELD,P_LEVELS, DYNDIA1A.337
& UCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.632
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.633
CALL V_INT
(PUV,PZ,V,WORK2,U_FIELD,P_LEVELS, DYNDIA1A.339
& VCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.634
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.635
CALL W_EQTOLL
(COEFF1,COEFF2,WORK1(FIRST_FLD_PT) GSM1F405.636
& ,WORK2(FIRST_FLD_PT),UCOMP_P(FIRST_FLD_PT,K) GSM1F405.637
& ,VCOMP_P(FIRST_FLD_PT,K),U_FLD_VALID,U_FLD_VALID) GSM1F405.638
ENDDO ! Levels loop DYNDIA1A.343
ELSE ! Test is false UCOMP_LEVS and VCOMP_LEVS are different DYNDIA1A.344
IF(QUCOMP_P) THEN DYNDIA1A.345
DO K=1,UCOMP_P_LEVS DYNDIA1A.346
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.639
PZ(I)=UCOMP_PRESS(K)*100.0 DYNDIA1A.348
ENDDO DYNDIA1A.349
CALL V_INT
(PUV,PZ,U,WORK1,U_FIELD,P_LEVELS, DYNDIA1A.350
& UCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.640
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.641
CALL V_INT
(PUV,PZ,V,WORK2,U_FIELD,P_LEVELS, DYNDIA1A.352
& VCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.642
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.643
CALL W_EQTOLL
(COEFF1,COEFF2,WORK1(FIRST_FLD_PT) GSM1F405.644
& ,WORK2(FIRST_FLD_PT),UCOMP_P(FIRST_FLD_PT,K),WORK3 GSM1F405.645
& ,U_FLD_VALID,U_FLD_VALID) GSM1F405.646
ENDDO ! Levels loop DYNDIA1A.356
ENDIF DYNDIA1A.357
IF(QVCOMP_P) THEN DYNDIA1A.358
DO K=1,VCOMP_P_LEVS DYNDIA1A.359
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.647
PZ(I)=VCOMP_PRESS(K)*100.0 DYNDIA1A.361
ENDDO DYNDIA1A.362
CALL V_INT
(PUV,PZ,U,WORK1,U_FIELD,P_LEVELS, DYNDIA1A.363
& UCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.648
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.649
CALL V_INT
(PUV,PZ,V,WORK2,U_FIELD,P_LEVELS, DYNDIA1A.365
& VCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.650
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.651
CALL W_EQTOLL
(COEFF1,COEFF2,WORK1(FIRST_FLD_PT) GSM1F405.652
& ,WORK2(FIRST_FLD_PT),WORK3,VCOMP_P(FIRST_FLD_PT,K) GSM1F405.653
& ,U_FLD_VALID,U_FLD_VALID) GSM1F405.654
ENDDO ! Levels loop DYNDIA1A.369
ENDIF DYNDIA1A.370
ENDIF ! End of TEST IF block DYNDIA1A.371
C DYNDIA1A.372
ELSE ! Winds not to be rotated-------------------------- DYNDIA1A.373
*ENDIF DYNDIA1A.374
C DYNDIA1A.375
CL------------------Interpolate U cmpnt of wind onto Pressure --------- DYNDIA1A.376
DYNDIA1A.377
IF(QUCOMP_P) THEN DYNDIA1A.378
DO K=1,UCOMP_P_LEVS DYNDIA1A.379
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.655
PZ(I)=UCOMP_PRESS(K)*100.0 ! convert to Pascals DYNDIA1A.381
ENDDO DYNDIA1A.382
CALL V_INT
(PUV,PZ,U,UCOMP_P(1,K),U_FIELD,P_LEVELS, DYNDIA1A.383
& UCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.656
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.657
ENDDO ! Levels loop DYNDIA1A.385
ENDIF DYNDIA1A.386
CL------------------Interpolate V cmpnt of wind onto Pressure --------- DYNDIA1A.387
DYNDIA1A.388
IF(QVCOMP_P) THEN DYNDIA1A.389
DO K=1,VCOMP_P_LEVS DYNDIA1A.390
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.658
PZ(I)=VCOMP_PRESS(K)*100.0 ! convert to Pascals DYNDIA1A.392
ENDDO DYNDIA1A.393
CALL V_INT
(PUV,PZ,V,VCOMP_P(1,K),U_FIELD,P_LEVELS, DYNDIA1A.394
& VCOMP_MAX_WIND,MAX_WIND_PRESSURE,.FALSE. GSM1F405.659
& ,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.660
ENDDO ! Levels loop DYNDIA1A.396
ENDIF DYNDIA1A.397
*IF -DEF,GLOBAL DYNDIA1A.398
ENDIF ! End of ROTATE IF block DYNDIA1A.399
*ENDIF DYNDIA1A.400
ENDIF ! End of QUCOMP or QVCOMP DYNDIA1A.401
CL------------------Calculate the maximum wind------------------------- DYNDIA1A.402
DYNDIA1A.403
IF (QUCOMP_MAX_WIND.AND. DYNDIA1A.404
* QVCOMP_MAX_WIND.AND. DYNDIA1A.405
* QMAX_WIND_PRESSURE) THEN DYNDIA1A.406
*IF -DEF,GLOBAL DYNDIA1A.407
IF(ROTATE_MAX_UV) THEN DYNDIA1A.408
CALL WINDMAX
( GSM3F403.18
*CALL ARGFLDPT
GSM3F403.19
& PSTAR,U,V,U_ROWS,P_ROWS,ROW_LENGTH,P_LEVELS, GSM3F403.20
* P_FIELD,U_FIELD,AK,BK,AKH,BKH,ETA_MATRIX_INV,MATRIX_P_O, DYNDIA1A.410
* WORK1,WORK2,MAX_WIND_PRESSURE) DYNDIA1A.411
CALL W_EQTOLL
(COEFF1,COEFF2,WORK1(FIRST_FLD_PT) GSM1F405.661
& ,WORK2(FIRST_FLD_PT),UCOMP_MAX_WIND(FIRST_FLD_PT) GSM1F405.662
& ,VCOMP_MAX_WIND(FIRST_FLD_PT),U_FLD_VALID,U_FLD_VALID) GSM1F405.663
ELSE ! Do not rotate DYNDIA1A.414
*ENDIF DYNDIA1A.415
CALL WINDMAX
( GSM3F403.21
*CALL ARGFLDPT
GSM3F403.22
& PSTAR,U,V,U_ROWS,P_ROWS,ROW_LENGTH,P_LEVELS, GSM3F403.23
* P_FIELD,U_FIELD,AK,BK,AKH,BKH,ETA_MATRIX_INV,MATRIX_P_O, DYNDIA1A.417
* UCOMP_MAX_WIND,VCOMP_MAX_WIND,MAX_WIND_PRESSURE) DYNDIA1A.418
*IF -DEF,GLOBAL DYNDIA1A.419
ENDIF DYNDIA1A.420
*ENDIF DYNDIA1A.421
CL--- ICAO HT of the MAX WIND pressure ? DYNDIA1A.422
IF(QMAX_WIND_ICAO_HEIGHT) THEN DYNDIA1A.423
CALL ICAO_HT
(MAX_WIND_PRESSURE(FIRST_FLD_PT),U_FLD_VALID GSM1F405.664
& ,MAX_WIND_ICAO_HEIGHT(FIRST_FLD_PT)) GSM1F405.665
ENDIF DYNDIA1A.425
ELSEIF(QUCOMP_MAX_WIND.NEQV.QVCOMP_MAX_WIND.OR.QUCOMP_MAX_WIND. DYNDIA1A.426
* NEQV.QMAX_WIND_PRESSURE.OR.QVCOMP_MAX_WIND.NEQV. DYNDIA1A.427
* QMAX_WIND_PRESSURE)THEN DYNDIA1A.428
WRITE(6,*)' Subroutine WINDMAX not called - U & VCOMP_MAX_WIND'
GIE0F403.123
WRITE(6,*)' and MAX_WIND_PRESSURE all must be selected' GIE0F403.124
ENDIF ! Top IF block for maxwind DYNDIA1A.431
DYNDIA1A.432
DYNDIA1A.433
C----------------------------------------------------------------------- DYNDIA1A.434
CL Section 15 Item 205 CAT PROBABILITY at pressure levels DYNDIA1A.435
C----------------------------------------------------------------------- DYNDIA1A.436
IF(QCAT_PROB_SINGLE.AND.QMAX_WIND_PRESSURE.AND.QUCOMP_MAX_WIND. DYNDIA1A.437
* AND.QVCOMP_MAX_WIND)THEN DYNDIA1A.438
DO K=1,CAT_PROB_LEVS DYNDIA1A.439
DO I=FIRST_VALID_PT,LAST_U_VALID_PT GSM1F405.666
PZ(I)=CAT_PROB_PRESS(K)*100.0 ! Convert to pascals DYNDIA1A.441
ENDDO DYNDIA1A.442
CALL CAT
( GSM3F403.24
*CALL ARGFLDPT
GSM3F403.25
*IF DEF,MPP GSM3F403.26
& glsize(2), GSM3F403.27
*ENDIF GSM3F403.28
& U,V,PUV,PSTAR,PZ,MAX_WIND_PRESSURE, GSM3F403.29
& CAT_PROB_SINGLE(1,K),P_FIELD,U_FIELD,P_LEVELS,ROW_LENGTH, DYNDIA1A.444
& P_ROWS,SEC_U_LATITUDE,AK,BK,EW_SPACE,NS_SPACE) DYNDIA1A.445
ENDDO DYNDIA1A.446
ELSEIF(QCAT_PROB_SINGLE.AND.(.NOT.(QMAX_WIND_PRESSURE.AND. DYNDIA1A.447
* QUCOMP_MAX_WIND.AND.QVCOMP_MAX_WIND)))THEN DYNDIA1A.448
WRITE(6,*)' Subroutine CAT not called - PRESSURE, U & VCOMP of'
GIE0F403.125
WRITE(6,*)' MAX WIND must be selected as well as CAT_PROB_SINGLE' GIE0F403.126
ENDIF DYNDIA1A.451
C----------------------------------------------------------------------- DYNDIA1A.452
CL Section 15 Item 211 MEAN CAT PROBABILITY over levels 300,250,200mb DYNDIA1A.453
C----------------------------------------------------------------------- DYNDIA1A.454
IF(QCAT_PROB_MEAN.AND.QMAX_WIND_PRESSURE.AND.QUCOMP_MAX_WIND. DYNDIA1A.455
* AND.QVCOMP_MAX_WIND)THEN DYNDIA1A.456
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.667
CAT_PROB_MEAN(I)=0.0 DYNDIA1A.458
ENDDO DYNDIA1A.459
C----------------------------------------------------------------------- DYNDIA1A.460
CL Call CAT for three levels - CAT probability output to WORK1 DYNDIA1A.461
C----------------------------------------------------------------------- DYNDIA1A.462
DO K=1,3 DYNDIA1A.463
IF(K.EQ.1)PZ(FIRST_VALID_PT)=30000.0 GSM1F405.668
IF(K.EQ.2)PZ(FIRST_VALID_PT)=25000.0 GSM1F405.669
IF(K.EQ.3)PZ(FIRST_VALID_PT)=20000.0 GSM1F405.670
DO I=FIRST_VALID_PT+1,LAST_U_VALID_PT GSM1F405.671
PZ(I)=PZ(FIRST_VALID_PT) GSM1F405.672
ENDDO DYNDIA1A.469
CALL CAT
( GSM3F403.30
*CALL ARGFLDPT
GSM3F403.31
*IF DEF,MPP GSM3F403.32
& glsize(2), GSM3F403.33
*ENDIF GSM3F403.34
& U,V,PUV,PSTAR,PZ,MAX_WIND_PRESSURE,WORK1,P_FIELD, GSM3F403.35
& U_FIELD,P_LEVELS,ROW_LENGTH,P_ROWS,SEC_U_LATITUDE,AK,BK, DYNDIA1A.471
& EW_SPACE,NS_SPACE) DYNDIA1A.472
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM3F403.36
CAT_PROB_MEAN(I)=CAT_PROB_MEAN(I)+WORK1(I) DYNDIA1A.474
ENDDO DYNDIA1A.475
ENDDO DYNDIA1A.476
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM3F403.37
CAT_PROB_MEAN(I)=CAT_PROB_MEAN(I)/3.0 DYNDIA1A.478
ENDDO DYNDIA1A.479
ELSEIF(QCAT_PROB_MEAN.AND.(.NOT.(QMAX_WIND_PRESSURE.AND. DYNDIA1A.480
* QUCOMP_MAX_WIND.AND.QVCOMP_MAX_WIND)))THEN DYNDIA1A.481
WRITE(6,*)' Subroutine CAT not called - PRESSURE, U & VCOMP of'
GIE0F403.127
WRITE(6,*)' MAX WIND must be selected as well as CAT_PROB_MEAN' GIE0F403.128
ENDIF DYNDIA1A.484
C----------------------------------------------------------------------- DYNDIA1A.485
CL Section 15 Items 203/204 MAXIMUM CAT PROBABILITY AND LEVEL DYNDIA1A.486
C----------------------------------------------------------------------- DYNDIA1A.487
IF(QMAX_CAT_PROB.AND.QMAX_CAT_LEVEL.AND.QMAX_WIND_PRESSURE.AND. DYNDIA1A.488
* QUCOMP_MAX_WIND.AND.QVCOMP_MAX_WIND)THEN DYNDIA1A.489
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.673
MAX_CAT_PROB(I)=0.0 DYNDIA1A.491
MAX_CAT_LEVEL(I)=30000.0 DYNDIA1A.492
ENDDO DYNDIA1A.493
C----------------------------------------------------------------------- DYNDIA1A.494
CL Call CAT for three levels - CAT probability output to WORK1 DYNDIA1A.495
C----------------------------------------------------------------------- DYNDIA1A.496
DO K=1,3 DYNDIA1A.497
IF(K.EQ.1)PZ(FIRST_VALID_PT)=30000.0 GSM1F405.674
IF(K.EQ.2)PZ(FIRST_VALID_PT)=25000.0 GSM1F405.675
IF(K.EQ.3)PZ(FIRST_VALID_PT)=20000.0 GSM1F405.676
DO I=FIRST_VALID_PT+1,LAST_U_VALID_PT GSM1F405.677
PZ(I)=PZ(FIRST_VALID_PT) GSM1F405.678
ENDDO DYNDIA1A.503
CALL CAT
( GSM3F403.38
*CALL ARGFLDPT
GSM3F403.39
*IF DEF,MPP GSM3F403.40
& glsize(2), GSM3F403.41
*ENDIF GSM3F403.42
& U,V,PUV,PSTAR,PZ,MAX_WIND_PRESSURE,WORK1,P_FIELD, GSM3F403.43
& U_FIELD,P_LEVELS,ROW_LENGTH,P_ROWS,SEC_U_LATITUDE,AK,BK, DYNDIA1A.505
& EW_SPACE,NS_SPACE) DYNDIA1A.506
IF(K.EQ.1)PZ(1)=30000.0 DYNDIA1A.507
IF(K.EQ.2)PZ(1)=25000.0 DYNDIA1A.508
IF(K.EQ.3)PZ(1)=20000.0 DYNDIA1A.509
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM3F403.44
IF(MAX_CAT_PROB(I).LT.WORK1(I))THEN DYNDIA1A.511
MAX_CAT_PROB(I)=WORK1(I) DYNDIA1A.512
MAX_CAT_LEVEL(I)=PZ(1) DYNDIA1A.513
ENDIF DYNDIA1A.514
ENDDO DYNDIA1A.515
ENDDO DYNDIA1A.516
ELSEIF(QMAX_CAT_PROB.NEQV.QMAX_CAT_LEVEL)THEN DYNDIA1A.517
WRITE(6,*)' Subroutine CAT not called - both MAX_CAT_PROB and'
GIE0F403.129
WRITE(6,*)' MAX_CAT_LEVEL must be selected' GIE0F403.130
ELSEIF(QMAX_CAT_PROB.NEQV.QMAX_CAT_LEVEL)THEN DYNDIA1A.520
WRITE(6,*)' Subroutine CAT not called - MAX_CAT_PROB and'
GIE0F403.131
WRITE(6,*)' MAX_CAT_LEVEL must both be selected' GIE0F403.132
ELSEIF(QMAX_CAT_PROB.AND.QMAX_CAT_LEVEL.AND. DYNDIA1A.523
* (.NOT.(QMAX_WIND_PRESSURE.AND.QUCOMP_MAX_WIND.AND. DYNDIA1A.524
* QVCOMP_MAX_WIND)))THEN DYNDIA1A.525
WRITE(6,*)' Subroutine CAT not called - PRESSURE, U & VCOMP of'
GIE0F403.133
WRITE(6,*)' MAX WIND must be selected as well as MAX_CAT_PROB ' GIE0F403.134
WRITE(6,*)' & MAX_CAT_LEVEL' GIE0F403.135
ENDIF DYNDIA1A.529
C----------------------------------------------------------------------- DYNDIA1A.530
CL Section 15 Items 212/213 50 M U and V components. DYNDIA1A.531
C----------------------------------------------------------------------- DYNDIA1A.532
C=====================================================================C DYNDIA1A.533
C 50 METRE WINDS C DYNDIA1A.534
C USE U(50)=C1*U(ETA2)+C2(ETA1) C DYNDIA1A.535
C C1=LOG(ETA50/ETA1)/LOG(ETA2/ETA1) = 0.135 C DYNDIA1A.536
C C2=LOG(ETA2/ETA50)/LOG(ETA2/ETA1) = 0.865 C DYNDIA1A.537
C ETA50=0.994 IE CORRESPONDS TO Z=50M AND TBAR=283K and assumes C DYNDIA1A.538
C the first 5 levels are sigma levels. C DYNDIA1A.539
C=====================================================================C DYNDIA1A.540
C First check between which levels the 50 M level lies (ETA=0.994) DYNDIA1A.541
IF(QUCOMP50_WIND.OR.QVCOMP50_WIND) THEN DYNDIA1A.542
ETA50=0.994 DYNDIA1A.543
found_levels=.FALSE. GPB0F405.170
K=1 GPB0F405.171
DO WHILE ((.NOT. found_levels) .AND. (K .LT. P_LEVELS)) GPB0F405.172
level1=K GPB0F405.173
level2=K+1 GPB0F405.174
ETA1=AK(level1)/PREF+BK(level1) GPB0F405.175
ETA2=AK(level2)/PREF+BK(level2) GPB0F405.176
GPB0F405.177
IF ((ETA1 .GE. ETA50) .AND. (ETA2 .LE. ETA50)) GPB0F405.178
& found_levels=.TRUE. GPB0F405.179
GPB0F405.180
K=K+1 GPB0F405.181
ENDDO GPB0F405.182
GPB0F405.183
IF (.NOT. found_levels) THEN GPB0F405.184
ICODE=1 GPB0F405.185
CMESSAGE='DYN_DIAG: Error in calculating 50 M winds' GPB0F405.186
RETURN GPB0F405.187
ENDIF GPB0F405.188
IF(ETA1.LT.ETA50) THEN DYNDIA1A.546
ICODE=1 DYNDIA1A.547
CMESSAGE='DYN_DIAG: Error in calculating 50 M winds' DYNDIA1A.548
RETURN DYNDIA1A.549
ENDIF DYNDIA1A.550
IF(ETA2.GT.ETA50) THEN DYNDIA1A.551
ICODE=1 DYNDIA1A.552
CMESSAGE='DYN_DIAG: Error in calculating 50 M winds' DYNDIA1A.553
RETURN DYNDIA1A.554
ENDIF DYNDIA1A.555
C1=ALOG(ETA50/ETA1)/ALOG(ETA2/ETA1) DYNDIA1A.556
C2=ALOG(ETA2/ETA50)/ALOG(ETA2/ETA1) DYNDIA1A.557
ENDIF DYNDIA1A.558
IF(QUCOMP50_WIND) THEN DYNDIA1A.559
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.679
UCOMP50_WIND(I)=C1*U(I,level2)+C2*U(I,level1) GPB0F405.189
ENDDO DYNDIA1A.562
ENDIF DYNDIA1A.563
IF(QVCOMP50_WIND) THEN DYNDIA1A.564
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.680
VCOMP50_WIND(I)=C1*V(I,level2)+C2*V(I,level1) GPB0F405.190
ENDDO DYNDIA1A.567
ENDIF DYNDIA1A.568
C want to compute pv for some theta level MM180193.30
if(qpotn_vort_theta.or.qpotn_vort_press.or.qtheta_on_pv)then TD141293.106
n_levels=p_levels-1 TD141293.107
call dthe_dp
(pstar,theta,p_field,p_levels TD141293.108
2 ,ak,bk,akh,bkh,n_levels TD141293.109
3 ,e_levels,dthe_dph) TD141293.110
endif TD141293.111
IF (QPOTN_VORT_THETA) THEN MM180193.31
DO I = 1,POTN_VORT_THETA_LEVS MM180193.32
CALL CALC_PV
DYNDIA1A.572
1 (PSTAR,THETA,U,V,P_FIELD,U_FIELD,P_LEVELS, DYNDIA1A.573
2 ROW_LENGTH, GSM3F403.45
*CALL ARGFLDPT
GSM3F403.46
& RMDI,AK,BK,DESIRED_THETA(I),F3, GSM3F403.47
& e_levels,n_levels,dthe_dph, GSM3F403.48
3 LATITUDE_STEP_INVERSE,LONGITUDE_STEP_INVERSE, DYNDIA1A.575
4 COS_U_LATITUDE,SEC_P_LATITUDE, DYNDIA1A.576
5 POTN_VORT_THETA(1,I),LLINTS) GSS1F304.192
ENDDO DYNDIA1A.579
ENDIF MM180193.34
MM180193.35
C want to compute pv for some pressure level MM180193.36
IF (QPOTN_VORT_PRESS) THEN MM180193.37
DO I = 1,POTN_VORT_P_LEVS MM180193.38
CALL CALC_PV_P
MM180193.39
1 (PSTAR,THETA,U,V,P_FIELD,U_FIELD,P_LEVELS, MM180193.40
2 ROW_LENGTH, GSM3F403.49
*CALL ARGFLDPT
GSM3F403.50
& RMDI,AK,BK,PV_PRESS(I),F3, GSM3F403.51
& e_levels,n_levels,dthe_dph, TD141293.113
3 LATITUDE_STEP_INVERSE,LONGITUDE_STEP_INVERSE, MM180193.42
4 COS_U_LATITUDE,SEC_P_LATITUDE, MM180193.43
5 POTN_VORT_ON_P(1,I),THETA_ON_P(1,I),LLINTS) GSS1F304.193
ENDDO MM180193.45
END IF MM180193.46
C want to compute theta on some pv surfaces. MM180193.47
C the loop over the surfaces is contained inside the subroutine. MM180193.48
IF (QTHETA_ON_PV) THEN MM180193.49
CALL THETA_PV
( MM180193.50
1 PSTAR,THETA,U,V,P_FIELD,U_FIELD,P_LEVELS, MM180193.51
2 ROW_LENGTH, GSM3F403.52
*CALL ARGFLDPT
GSM3F403.53
& RMDI,AK,BK,F3, GSM3F403.54
& e_levels,n_levels,dthe_dph, TD141293.114
3 THETA_PV_LEVS,DESIRED_PV, MM180193.53
4 THETA_PV_P_LEVS,REQ_THETA_PV_LEVS, MM180193.54
5 LATITUDE_STEP_INVERSE,LONGITUDE_STEP_INVERSE, MM180193.55
6 COS_U_LATITUDE,SEC_P_LATITUDE, MM180193.56
7 THETA_ON_PV,LLINTS) GSS1F304.194
ENDIF DYNDIA1A.580
DYNDIA1A.581
C----------------------------------------------------------------------- DYNDIA1A.582
CL Section 15 item 215, UV on pressure levels RR250193.38
CL Only possible if U and V have been requested on pressure levels DYNDIA1A.584
CL required for UV DYNDIA1A.585
CL DYNDIA1A.586
IF (QUV_P) THEN DYNDIA1A.587
do K=1,UV_P_LEVS DYNDIA1A.588
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.681
UV_P(I,K) = UCOMP_P(I,UV_IND(K,1))*VCOMP_P(I,UV_IND(K,2)) DYNDIA1A.590
ENDDO DYNDIA1A.591
ENDDO DYNDIA1A.592
ENDIF DYNDIA1A.593
C----------------------------------------------------------------------- DYNDIA1A.594
CL Section 15 item 216, T on pressure levels on the wind grid DYNDIA1A.595
CL T is first calculated on p-grid for the pressure level and DYNDIA1A.596
CL then interpolated to the u-grid. DYNDIA1A.597
CL DYNDIA1A.598
IF (QT_P) THEN DYNDIA1A.599
T_REF=2 !used in vertical interpolation DYNDIA1A.600
DO K=1,T_P_LEVS DYNDIA1A.601
DO I=FIRST_FLD_PT,P_FIELD GSM1F405.682
PZ(I)=T_PRESS(K)*100.0 ! convert to pascals DYNDIA1A.603
ENDDO DYNDIA1A.604
CALL V_INT_T
(WORK5,PZ,P(1,T_REF),PSTAR,P_EXNER_HALF GSM1F405.683
& ,THETA,P_FIELD,P_LEVELS,T_REF,AKH,BKH GSM1F405.684
& ,FIRST_FLD_PT,P_FIELD) GSM1F405.685
*IF DEF,MPP GSM1F405.686
CALL P_TO_UV
(WORK5(FIRST_FLD_PT),T_P(FIRST_FLD_PT,K) GSM1F405.687
& ,P_FIELD-FIRST_FLD_PT+1,U_FIELD-FIRST_FLD_PT+1 GSM1F405.688
& ,ROW_LENGTH,P_LAST_ROW) GSM1F405.689
*ELSE GSM1F405.690
CALL P_TO_UV
(WORK5,T_P(1,K),P_FIELD,U_FIELD,ROW_LENGTH GSM1F405.691
& ,P_ROWS) GSM1F405.692
*ENDIF GSM1F405.693
ENDDO DYNDIA1A.609
C----------------------------------------------------------------------- DYNDIA1A.610
CL Section 15 item 217, UT on pressure levels DYNDIA1A.611
CL DYNDIA1A.612
IF (QUT_P) THEN DYNDIA1A.613
do K=1,UT_P_LEVS DYNDIA1A.614
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.694
UT_P(I,K) = UCOMP_P(I,UT_IND(K,1))*T_P(I,UT_IND(K,2)) DYNDIA1A.616
ENDDO DYNDIA1A.617
ENDDO DYNDIA1A.618
ENDIF DYNDIA1A.619
C----------------------------------------------------------------------- DYNDIA1A.620
CL Section 15 item 218, VT on pressure levels DYNDIA1A.621
CL DYNDIA1A.622
IF (QVT_P) THEN DYNDIA1A.623
do K=1,VT_P_LEVS DYNDIA1A.624
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.695
VT_P(I,K) = VCOMP_P(I,VT_IND(K,1))*T_P(I,VT_IND(K,2)) DYNDIA1A.626
ENDDO DYNDIA1A.627
ENDDO DYNDIA1A.628
ENDIF DYNDIA1A.629
C----------------------------------------------------------------------- DYNDIA1A.630
CL Section 15 item 219, T**2 on pressure levels DYNDIA1A.631
CL DYNDIA1A.632
IF (QT2_P) THEN DYNDIA1A.633
DO K=1,T2_P_LEVS DYNDIA1A.634
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.696
T2_P(I,K) = T_P(I,T2_IND(K))*T_P(I,T2_IND(K)) DYNDIA1A.636
ENDDO DYNDIA1A.637
ENDDO DYNDIA1A.638
ENDIF DYNDIA1A.639
ENDIF DYNDIA1A.640
C----------------------------------------------------------------------- DYNDIA1A.641
CL Section 15 item 220, U2 on pressure levels DYNDIA1A.642
CL Only possible if U has been requested on the same pressure levels DYNDIA1A.643
CL DYNDIA1A.644
IF (QU2_P) THEN DYNDIA1A.645
DO K=1,U2_P_LEVS DYNDIA1A.646
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.697
U2_P(I,K) = UCOMP_P(I,U2_IND(K))*UCOMP_P(I,U2_IND(K)) DYNDIA1A.648
ENDDO DYNDIA1A.649
ENDDO DYNDIA1A.650
ENDIF DYNDIA1A.651
C----------------------------------------------------------------------- DYNDIA1A.652
CL Section 15 item 221, v2 on pressure levels DYNDIA1A.653
CL Only possible if v has been requested on the same pressure levels DYNDIA1A.654
CL DYNDIA1A.655
IF (QV2_P) THEN DYNDIA1A.656
DO K=1,V2_P_LEVS DYNDIA1A.657
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.698
V2_P(I,K) = VCOMP_P(I,V2_IND(K))*VCOMP_P(I,V2_IND(K)) DYNDIA1A.659
ENDDO DYNDIA1A.660
ENDDO DYNDIA1A.661
ENDIF DYNDIA1A.662
C----------------------------------------------------------------------- DYNDIA1A.663
CL Section 15 item 222, w on pressure levels and wind grid DYNDIA1A.664
CL DYNDIA1A.665
CL DYNDIA1A.666
IF (QW_P) THEN DYNDIA1A.667
CALL OMEGA_DIAG
( GSM3F403.59
*CALL ARGFLDPT
GSM3F403.60
& U,V,OMEGA,SEC_P_LATITUDE,COS_U_LATITUDE, GSM3F403.61
1 PSTAR,PSTAR_OLD,DELTA_AK,DELTA_BK, DYNDIA1A.669
2 AK,BK,AKH,BKH,U_FIELD,P_FIELD,P_LEVELS, DYNDIA1A.670
3 ROW_LENGTH,LATITUDE_STEP_INVERSE, DYNDIA1A.671
4 LONGITUDE_STEP_INVERSE,ADVECTION_TIMESTEP) DYNDIA1A.672
DYNDIA1A.673
C omega returned at u points on all p_levels DYNDIA1A.674
C Interpolation of omega to pressure levels DYNDIA1A.675
DO K=1,W_P_LEVS DYNDIA1A.676
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.699
PZ(I) = W_PRESS(K)*100.0 ! convert to Pascals DYNDIA1A.678
ENDDO DYNDIA1A.679
CALL V_INT
(PUV,PZ,OMEGA,W_P(1,K),U_FIELD,P_LEVELS,WORK1,WORK5, DYNDIA1A.680
& .FALSE.,FIRST_FLD_PT,LAST_U_FLD_PT) GSM1F405.700
ENDDO DYNDIA1A.682
DYNDIA1A.683
C----------------------------------------------------------------------- DYNDIA1A.684
CL Section 15 item 223-225, wT, wu, wv on pressure levels DYNDIA1A.685
CL Only possible if w and T have been requested on the same pressure DYNDIA1A.686
CL levels DYNDIA1A.687
IF (QWT_P) THEN DYNDIA1A.688
do K=1,WT_P_LEVS DYNDIA1A.689
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.701
WT_P(I,K) = W_P(I,WT_IND(K,1))*T_P(I,WT_IND(K,2)) DYNDIA1A.691
ENDDO DYNDIA1A.692
ENDDO DYNDIA1A.693
ENDIF DYNDIA1A.694
IF (QWU_P) THEN DYNDIA1A.695
do K=1,WU_P_LEVS DYNDIA1A.696
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.702
WU_P(I,K) = W_P(I,WU_IND(K,1))*UCOMP_P(I,WU_IND(K,2)) DYNDIA1A.698
ENDDO DYNDIA1A.699
ENDDO DYNDIA1A.700
ENDIF DYNDIA1A.701
IF (QWV_P) THEN DYNDIA1A.702
do K=1,WV_P_LEVS DYNDIA1A.703
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.703
WV_P(I,K) = W_P(I,WV_IND(K,1))*VCOMP_P(I,WV_IND(K,2)) DYNDIA1A.705
ENDDO DYNDIA1A.706
ENDDO DYNDIA1A.707
ENDIF DYNDIA1A.708
ENDIF DYNDIA1A.709
C----------------------------------------------------------------------- DYNDIA1A.710
CL Section 15 item 226-228, q, qu, qv on pressure levels DYNDIA1A.711
CL item 235 qw on pressure levels ARS1F404.50
CL DYNDIA1A.712
IF (QQ_P) THEN DYNDIA1A.714
DO K=1,Q_P_LEVS DYNDIA1A.715
DO I=FIRST_FLD_PT,P_FIELD GSM1F405.704
PZ(I) = Q_PRESS(K)*100.0 ! convert to Pascals DYNDIA1A.717
ENDDO DYNDIA1A.718
CALL V_INT
(P,PZ,Q,WORK5,P_FIELD,Q_LEVELS,WORK1,WORK1,.FALSE. GSM1F405.705
& ,FIRST_FLD_PT,P_FIELD) GSM1F405.706
*IF DEF,MPP GSM1F405.707
CALL P_TO_UV
(WORK5(FIRST_FLD_PT),Q_P(FIRST_FLD_PT,K) GSM1F405.708
& ,P_FIELD-FIRST_FLD_PT+1,U_FIELD-FIRST_FLD_PT+1,ROW_LENGTH GSM1F405.709
& ,P_LAST_ROW) GSM1F405.710
*ELSE GSM1F405.711
CALL P_TO_UV
(WORK5,Q_P(1,K),P_FIELD,U_FIELD,ROW_LENGTH GSM1F405.712
& ,P_ROWS) GSM1F405.713
*ENDIF GSM1F405.714
ENDDO DYNDIA1A.721
DYNDIA1A.722
IF (QuQ_P) THEN DYNDIA1A.723
do K=1,QU_P_LEVS DYNDIA1A.724
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.715
UQ_P(I,K) = Q_P(I,QU_IND(K,1))*UCOMP_P(I,QU_IND(K,2)) DYNDIA1A.726
ENDDO DYNDIA1A.727
ENDDO DYNDIA1A.728
ENDIF DYNDIA1A.729
IF (QVQ_P) THEN DYNDIA1A.730
do K=1,QV_P_LEVS DYNDIA1A.731
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.716
VQ_P(I,K) = Q_P(I,QV_IND(K,1))*VCOMP_P(I,QV_IND(K,2)) DYNDIA1A.733
ENDDO DYNDIA1A.734
ENDDO DYNDIA1A.735
ENDIF DYNDIA1A.736
! added later therefore jump in stashcode ARS1F404.51
IF (QWQ_P) THEN ARS1F404.52
do K=1,QW_P_LEVS ARS1F404.53
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.717
WQ_P(I,K) = Q_P(I,QW_IND(K,1))*W_P(I,QW_IND(K,2)) ARS1F404.55
ENDDO ARS1F404.56
ENDDO ARS1F404.57
ENDIF ARS1F404.58
ENDIF RR250193.39
C --------------------------------------------------------------------- RR250193.40
CL Section 15 items 231,232,233,234 test diagnostics RR250193.41
CL 231 single level on u grid, 232 single level on p grid RR250193.42
CL 233 press levels on p grid, 234 model levels on p grid RR250193.43
CL RR250193.44
IF (QDIA1.OR.QDIA2.OR.QDIA3.OR.QDIA4) THEN RR250193.45
CALL TESTDIAG
( RR250193.46
1 P_FIELD,U_FIELD,P_ROWS,U_ROWS,ROW_LENGTH,EW_SPACE,NS_SPACE, RR250193.47
2 NMOST_LAT,WMOST_LONG,ELF,PHI_POLE,LAMBDA_POLE, RR250193.48
3 TESTD_PRESS,TESTD_P_LEVS, RR250193.49
4 TESTD_MODEL,TESTD_M_LEVS,FORECAST_HRS, RR250193.50
5 TESTDIAG1,TESTDIAG2,TESTDIAG3,TESTDIAG4, RR250193.51
6 QDIA1,QDIA2,QDIA3,QDIA4) RR250193.52
C RR250193.53
ENDIF ARS1F404.59
! --------------------------------------------------------------------- ARS1F404.60
!L Items 236 and 237 both require pstar on the uv grid ARS1F404.61
!L ARS1F404.62
IF (QHEAVY_P.or.QTOTAL_KE) THEN ARS1F404.63
CALL P_TO_UV
(PSTAR,PSTAR_UV,P_FIELD,U_FIELD,ROW_LENGTH,P_ROWS) ARS1F404.64
ARS1F404.69
!L Section 15 item 236 Heavyside function on pressure levels for ARS1F404.70
!L u-grid. ARS1F404.71
!L The Heavyside function is defined as 1.0 if the pressure level ARS1F404.72
!L is above the surface (i.e. pstar) and 0.0 if below. A time mean of ARS1F404.73
!L this will give information on the fraction of time a pressure ARS1F404.74
!L level is above the land or sea surface. ARS1F404.75
ARS1F404.76
IF (QHEAVY_P) THEN ARS1F404.77
DO K=1,HEAVY_P_LEVS ARS1F404.78
PLEV=HEAVY_PRESS(K)*100. ! pressure in Pascals ARS1F404.79
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.718
IF (PSTAR_UV(I).LT.PLEV) THEN ARS1F404.81
HEAVYSIDE_P(I,K)=0.0 ARS1F404.82
ELSE ARS1F404.83
HEAVYSIDE_P(I,K)=1.0 ARS1F404.84
ENDIF ARS1F404.85
ENDDO ARS1F404.86
ENDDO ARS1F404.87
ENDIF ARS1F404.88
! --------------------------------------------------------------------- ARS1F404.89
!L Section 15 item 237 Total kinetic energy in a column u-grid ARS1F404.90
!L ARS1F404.91
!L KE = SUM [0.5/g (u*u + v*v)] dp over model levels. ARS1F404.92
!L ARS1F404.93
!L Output scaled by 1.0e-6 to prevent accuracy problems ARS1F404.94
ARS1F404.95
IF (QTOTAL_KE) THEN ARS1F404.96
FACTOR=0.5*1.0e-6/g ARS1F404.97
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.719
TOTAL_KE(I)=0.0 ARS1F404.99
ENDDO ARS1F404.100
DO K=1,P_LEVELS ARS1F404.101
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.720
TOTAL_KE(I)=TOTAL_KE(I) - ARS1F404.103
& FACTOR*(U(I,K)*U(I,K)+V(I,K)*V(I,K)) ARS1F404.104
& *(DELTA_AK(K)+DELTA_BK(K)*PSTAR_UV(I)) ARS1F404.105
ENDDO ARS1F404.106
ENDDO ARS1F404.107
ENDIF ARS1F404.108
ENDIF ARS1F404.109
! --------------------------------------------------------------------- ARS1F404.110
!L Section 15 item 238 Geopotential on pressure levels ARS1F404.111
IF (QZ_P) THEN ARS1F404.112
CALL V_INT_ZH
(P_EXNER_HALF,THETA,Q,PHI_STAR,MODEL_HALF_HEIGHT, ARS1F404.113
& P_FIELD,P_LEVELS,Q_LEVELS) ARS1F404.114
DO K=1,Z_P_LEVS ARS1F404.115
DO I=1,P_FIELD ARS1F404.116
PZ(I)=Z_PRESS(k)*100.0 ! convert to pascals ARS1F404.117
ENDDO ARS1F404.118
CALL V_INT_Z
(PZ,P(1,Z_REF),PSTAR,P_EXNER_HALF,THETA,Q, ARS1F404.119
& MODEL_HALF_HEIGHT,WORK5,P_FIELD,p_LEVELS,Q_LEVELS, ARS1F404.120
& Z_REF,AKH,BKH,FIRST_FLD_PT,P_FIELD) GSM1F405.721
ARS1F404.122
! put on u grid ARS1F404.123
*IF DEF,MPP GSM1F405.722
CALL P_TO_UV
(WORK5(FIRST_FLD_PT),Z_P(FIRST_FLD_PT,k) GSM1F405.723
& ,P_FIELD-FIRST_FLD_PT+1,U_FIELD-FIRST_FLD_PT+1 GSM1F405.724
& ,ROW_LENGTH,P_LAST_ROW) GSM1F405.725
*ELSE GSM1F405.726
CALL P_TO_UV
(WORK5,Z_P(1,k),P_FIELD,U_FIELD,ROW_LENGTH GSM1F405.727
& ,P_ROWS) GSM1F405.728
*ENDIF GSM1F405.729
ENDDO ARS1F404.125
ARS1F404.129
!L Section 15 item 239 U*Z ARS1F404.130
ARS1F404.131
IF (QUZ_P) THEN ARS1F404.132
DO K=1,UZ_P_LEVS ARS1F404.133
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.730
UZ_P(I,K) = UCOMP_P(I,UZ_IND(k,1))*Z_P(I,UZ_IND(K,2)) ARS1F404.135
ENDDO ARS1F404.136
ENDDO ARS1F404.137
ENDIF ARS1F404.138
ARS1F404.139
ARS1F404.140
!L Section 15 item 240 V*Z ARS1F404.141
ARS1F404.142
IF (QVZ_P) THEN ARS1F404.143
DO K=1,VZ_P_LEVS ARS1F404.144
DO I=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.731
VZ_P(I,K) = VCOMP_P(I,VZ_IND(k,1))*Z_P(I,VZ_IND(K,2)) ARS1F404.146
ENDDO ARS1F404.147
ENDDO ARS1F404.148
ENDIF ARS1F404.149
ENDIF ARS1F404.150
ARS1F404.151
!L Section 15 item 241 mountain torque per unit area ARS1F404.152
!L ARS1F404.153
!L a*orog* E-W pressure gradient = a* orography * dp /(a dlong) ARS1F404.154
ARS1F404.155
IF (Q_MT) THEN ARS1F404.156
FACTOR=LONGITUDE_STEP_INVERSE ARS1F404.157
! wrong at row ends if non MPP code ARS1F404.161
DO i=FIRST_FLD_PT,LAST_U_FLD_PT GSM1F405.732
M_TORQUE(I)=0.25*factor*( (orog(i)+orog(i+1))* ARS1F404.163
& (pstar(i+1)-pstar(i)) ARS1F404.164
& + (orog(i+row_length)+orog(i+row_length+1))* ARS1F404.165
& (pstar(i+row_length+1)-pstar(i+row_length))) ARS1F404.166
ENDDO ARS1F404.167
*IF DEF,GLOBAL ARS1F404.168
*IF -DEF,MPP ARS1F404.169
! need to correct end of rows. Not required for MPP because of halo. ARS1F404.170
DO I=FIRST_ROW,U_LAST_ROW ARS1F404.171
! last point on row need to use first point on row ARS1F404.172
J=I*ROW_LENGTH ARS1F404.173
M_TORQUE(J)=0.25*factor*( (orog(J)+orog(J-ROW_LENGTH+1)) ARS1F404.174
& *(pstar(J-row_length+1)-pstar(J)) ARS1F404.175
& + (orog(J+row_length)+orog(J+1))* ARS1F404.176
& (pstar(J+1)-pstar(J+row_length))) ARS1F404.177
ENDDO ARS1F404.178
ARS1F404.179
*ENDIF ARS1F404.180
*ENDIF ARS1F404.181
ENDIF DYNDIA1A.737
C --------------------------------------------------------------------- DYNDIA1A.738
DYNDIA1A.739
RETURN DYNDIA1A.740
END DYNDIA1A.741
*ENDIF DYNDIA1A.742