*IF DEF,C72_1A,AND,DEF,ATMOS,AND,DEF,OCEAN GLW1F404.46
C ******************************COPYRIGHT****************************** GTS2F400.7507
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.7508
C GTS2F400.7509
C Use, duplication or disclosure of this code is subject to the GTS2F400.7510
C restrictions as set forth in the contract. GTS2F400.7511
C GTS2F400.7512
C Meteorological Office GTS2F400.7513
C London Road GTS2F400.7514
C BRACKNELL GTS2F400.7515
C Berkshire UK GTS2F400.7516
C RG12 2SZ GTS2F400.7517
C GTS2F400.7518
C If no contract has been raised with this copy of the code, the use, GTS2F400.7519
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.7520
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.7521
C Modelling at the above address. GTS2F400.7522
C ******************************COPYRIGHT****************************** GTS2F400.7523
C GTS2F400.7524
SUBROUTINE PRE_AREAVER(GAPS_LAMBDA_SRCE,LAMBDA_SRCE 4PREARAV1.3
&,GAPS_PHI_SRCE,PHI_SRCE,CYCLIC_SRCE,LROW_SRCE,WANT,MASK_SRCE PREARAV1.4
&,GAPS_LAMBDA_TARG,LAMBDA_TARG,GAPS_PHI_TARG,PHI_TARG PREARAV1.5
&,CYCLIC_TARG,SPHERICAL CJG1F401.1
&,MAXL,COUNT_TARG,BASE_TARG,INDEX_SRCE,WEIGHT,ICODE,CMESSAGE) CJG1F401.2
CLL PREARAV1.8
CLL Subroutine PRE_AREAVER PREARAV1.9
CLL PREARAV1.10
CLL Calculate weights for area-averaging data on the source grid to PREARAV1.11
CLL data on the target grid. PREARAV1.12
CLL PREARAV1.13
CLL J.Gregory <- programmer of some or all of previous code or changes PREARAV1.14
CLL PREARAV1.15
CLL Model Modification history from model version 3.0: PREARAV1.16
CLL version Date PREARAV1.17
CLL 3.2 13/07/93 Changed CHARACTER*(*) to CHARACTER*(80) for TS150793.123
CLL portability. Author Tracey Smith. TS150793.124
CLL 3.3 3.9.93 Correct out-of-bounds errors for IXL and IYT JG030993.1
CLL 4.1 22.5.96 J.M.Gregory Add argument SPHERICAL to produce CJG1F401.3
CLL correct weights for a spherical surface. CJG1F401.4
CLL PREARAV1.19
CLL Standard, paper 4 version 4 (14.12.90) PREARAV1.20
CLL PREARAV1.21
CL PREARAV1.22
CL The source grid and target grid are each specified by a lambda set PREARAV1.23
CL and a phi set of coordinates delimiting the boxes. These sets are CJG1F401.5
CL supplied in 1D arrays aa_bb for coordinate aa=LAMBDA,PHI on grid CJG1F401.6
CL bb=SRCE,TARG. The number of gaps is specified by GAPS_aa_bb, CJG1F401.7
CL which is equal to the number of lines IF (CYCLIC_bb) and aa is CJG1F401.8
CL LAMBDA, otherwise one less. (By "gap" we mean the interval CJG1F401.9
CL spanning a box in one coordinate only. The total number of boxes, CJG1F401.10
CL or grid points, is the product of the numbers of gaps in the two CJG1F401.11
CL coordinates.) Whether the axes are cyclic is not known until CJG1F401.12
CL run-time, so the dimensions of the arrays LAMBDA_bb are not known CJG1F401.13
CL at compile-time, and they are dimensioned assumed-size. There are CJG1F401.14
CL no restrictions on the base meridian of lambda, and it does not CJG1F401.15
CL have to be the same for source and target grids. The lambda CJG1F401.16
CL coordinates should be in increasing order (running from left to CJG1F401.17
CL right), the phi decreasing (top to bottom). The coordinates must CJG1F401.18
CL be given in degrees for cyclic axes, because a range of 360 is CJG1F401.19
CL assumed, or IF (SPHERICAL), when trigonometric functions are CJG1F401.20
CL used. IF (SPHERICAL), the weights computed are for a spherical CJG1F401.21
CL surface, assuming that LAMBDA is longitude and PHI latitude. CJG1F401.22
CL Otherwise, LAMBDA and PHI are treated as Cartesian axes on a CJG1F401.23
CL plane. CJG1F401.24
CL PREARAV1.37
CL The array MASK_SRCE is the land/sea mask for the source grid. The PREARAV1.38
CL logical value indicating points to be used should be supplied in PREARAV1.39
CL WANT. The first dimension of MASK_SRCE should be supplied in PREARAV1.40
CL LROW_SRCE. Specifying this separately allows for the possibility PREARAV1.41
CL of extra rows and columns in MASK_SRCE which are to be ignored. PREARAV1.42
CL PREARAV1.43
CL The arrays COUNT_TARG and BASE_TARG should be dimensioned in the PREARAV1.44
CL calling program to the number of boxes in the target array. PREARAV1.45
CL PREARAV1.46
CL The arrays INDEX_SRCE and WEIGHT are returned in the form which PREARAV1.47
CL the area-averaging routine DO_AREAVER expects. They are continuous PREARAV1.48
CL lists comprising consecutive groups of entries. There is a group PREARAV1.49
CL for each target point, for which the number of entries is spec- PREARAV1.50
CL ified by COUNT_TARG, and the groups appear in the normal order of PREARAV1.51
CL grid points. The size required for INDEX_SRCE and WEIGHT depends PREARAV1.52
CL on how many source boxes go into each target box, on average, and PREARAV1.53
CL is not known at compile-time. The maximum that could be needed is PREARAV1.54
CL (GAPS_LAMBDA_SRCE+GAPS_LAMBDA_TARG)*(GAPS_PHI_SRCE+GAPS_PHI_TARG) PREARAV1.55
CL and the size to which the arrays are actually dimensioned should PREARAV1.56
CL be supplied in MAXL. The size used is returned in MAXL. It is the PREARAV1.57
CL responsibility of the calling routine to provide enough space. PREARAV1.58
CL PREARAV1.59
IMPLICIT NONE PREARAV1.60
C*L PREARAV1.61
INTEGER PREARAV1.62
& GAPS_LAMBDA_SRCE !IN number of lambda gaps in source grid PREARAV1.63
&,GAPS_PHI_SRCE !IN number of phi gaps in source grid PREARAV1.64
&,GAPS_LAMBDA_TARG !IN number of lambda gaps in target grid PREARAV1.65
&,GAPS_PHI_TARG !IN number of phi gaps in target grid PREARAV1.66
&,LROW_SRCE !IN first dimension of MASK_SRCE PREARAV1.67
&,MAXL !INOUT maximum entries in output lists PREARAV1.68
&,COUNT_TARG(GAPS_LAMBDA_TARG,GAPS_PHI_TARG) PREARAV1.69
C !OUT no. of source boxes in target box PREARAV1.70
&,BASE_TARG(GAPS_LAMBDA_TARG,GAPS_PHI_TARG) PREARAV1.71
C !OUT first index in list for target box PREARAV1.72
&,INDEX_SRCE(MAXL) !OUT list of source box indices PREARAV1.73
&,ICODE !OUT return code PREARAV1.74
LOGICAL PREARAV1.75
& CYCLIC_SRCE !IN source grid is cyclic PREARAV1.76
&,CYCLIC_TARG !IN target grid is cyclic PREARAV1.77
&,WANT !IN indicator of wanted points in mask PREARAV1.78
&,MASK_SRCE(LROW_SRCE,*) !IN land/sea mask for source grid PREARAV1.79
&,SPHERICAL !IN calculate weights for a sphere CJG1F401.25
REAL PREARAV1.80
& LAMBDA_SRCE(*) !IN source lambda line coordinates PREARAV1.81
&,PHI_SRCE(*) !IN source phi line coordinates PREARAV1.82
&,LAMBDA_TARG(*) !IN target lambda line coordinates PREARAV1.83
&,PHI_TARG(*) !IN target phi line coordinates PREARAV1.84
&,WEIGHT(MAXL) !OUT list of weights for source boxes PREARAV1.85
CHARACTER*(80) CMESSAGE !OUT error message TS150793.125
C* PREARAV1.88
*CALL C_PI 
CJG1F401.26
C CJG1F401.27
INTEGER PREARAV1.89
& LINES_LAMBDA_SRCE ! number of source lambda lines PREARAV1.90
&,LINES_PHI_SRCE ! number of source phi lines PREARAV1.91
&,LINES_LAMBDA_TARG ! number of target lambda lines PREARAV1.92
&,LINES_PHI_TARG ! number of target phi lines PREARAV1.93
&,COUNT_LAMBDA(GAPS_LAMBDA_TARG) PREARAV1.94
C ! number of source lambda gaps per target PREARAV1.95
&,COUNT_PHI(GAPS_PHI_TARG) PREARAV1.96
C ! number of source phi gaps per target PREARAV1.97
&,INDEX_LAMBDA(GAPS_LAMBDA_SRCE+GAPS_LAMBDA_TARG) PREARAV1.98
C ! source lambda gap indices PREARAV1.99
&,INDEX_PHI(GAPS_PHI_SRCE+GAPS_PHI_TARG) PREARAV1.100
C ! source phi gap indices PREARAV1.101
&,IX1,IY1,IX2,IY2 ! working SRCE/TARG LAMBDA/PHI indices PREARAV1.102
&,IX1N,IX1W ! working indices PREARAV1.103
&,IXL(GAPS_LAMBDA_TARG+1) ! source line on the left of target line JG030993.2
&,IX2N ! target gap to the right of IX2 PREARAV1.105
&,IYT(GAPS_PHI_TARG+1) ! source line above target line JG030993.3
&,IXP,IYP,IP ! pointers into lists PREARAV1.107
&,IX,IY,I ! loop indices PREARAV1.108
REAL PREARAV1.109
& BASLAM ! minimum lambda for TEMP coordinates PREARAV1.110
&,BTARG ! width of target gap PREARAV1.111
&,BLO,BHI ! limits of gap PREARAV1.112
&,TEMP_SRCE(GAPS_LAMBDA_SRCE+1) PREARAV1.113
C ! adjusted version of LAMBDA_SRCE PREARAV1.114
&,TEMP_TARG(GAPS_LAMBDA_TARG+1) PREARAV1.115
C ! adjusted version of LAMBDA_TARG PREARAV1.116
&,FRAC_LAMBDA(GAPS_LAMBDA_SRCE+GAPS_LAMBDA_TARG) PREARAV1.117
C ! fractions of target lambda gaps PREARAV1.118
&,FRAC_PHI(GAPS_PHI_SRCE+GAPS_PHI_TARG) PREARAV1.119
C ! fractions of target phi gaps PREARAV1.120
&,SUM ! sum of weights PREARAV1.121
C PREARAV1.122
CL 1 Set up the lambda coordinates to make them easier to handle. PREARAV1.123
C PREARAV1.124
CL 1.1 Produce in TEMP_SRCE a monotonically increasing set of angles PREARAV1.125
CL equivalent to LAMBDA_SRCE i.e. equal under modulo 360. PREARAV1.126
C PREARAV1.127
IF (CYCLIC_SRCE) THEN PREARAV1.128
LINES_LAMBDA_SRCE=GAPS_LAMBDA_SRCE PREARAV1.129
ELSE PREARAV1.130
LINES_LAMBDA_SRCE=GAPS_LAMBDA_SRCE+1 PREARAV1.131
ENDIF PREARAV1.132
BASLAM=LAMBDA_SRCE(1) PREARAV1.133
DO 10 IX1=1,LINES_LAMBDA_SRCE PREARAV1.134
IF (LAMBDA_SRCE(IX1).LT.BASLAM) THEN PREARAV1.135
TEMP_SRCE(IX1)=LAMBDA_SRCE(IX1)+360. PREARAV1.136
ELSE PREARAV1.137
TEMP_SRCE(IX1)=LAMBDA_SRCE(IX1) PREARAV1.138
ENDIF PREARAV1.139
10 CONTINUE ! Labelled DO for the sake of fpp PREARAV1.140
C PREARAV1.141
CL 1.2 Produce in TEMP_TARG a set of angles equivalent to PREARAV1.142
CL LAMBDA_TARG i.e. equal under modulo 360, but all in the range PREARAV1.143
CL BASLAM to BASLAM+360, where BASLAM=min(TEMP_LAMBDA). PREARAV1.144
C PREARAV1.145
IF (CYCLIC_TARG) THEN PREARAV1.146
LINES_LAMBDA_TARG=GAPS_LAMBDA_TARG PREARAV1.147
ELSE PREARAV1.148
LINES_LAMBDA_TARG=GAPS_LAMBDA_TARG+1 PREARAV1.149
ENDIF PREARAV1.150
DO 20 IX2=1,LINES_LAMBDA_TARG PREARAV1.151
TEMP_TARG(IX2)=MOD(LAMBDA_TARG(IX2)-BASLAM,360.) PREARAV1.152
IF (TEMP_TARG(IX2).LT.0.) TEMP_TARG(IX2)=TEMP_TARG(IX2)+360. PREARAV1.153
TEMP_TARG(IX2)=TEMP_TARG(IX2)+BASLAM PREARAV1.154
20 CONTINUE ! Labelled DO for the sake of fpp PREARAV1.155
C PREARAV1.156
CL 2 For each target lambda line, find the index of the next source PREARAV1.157
CL lambda line to the left. JG030993.4
C PREARAV1.161
DO IX2=1,LINES_LAMBDA_TARG PREARAV1.162
DO IX1=1,LINES_LAMBDA_SRCE PREARAV1.163
IF (TEMP_TARG(IX2).GE.TEMP_SRCE(IX1)) IXL(IX2)=IX1 PREARAV1.164
ENDDO PREARAV1.165
ENDDO PREARAV1.166
C PREARAV1.167
CL 3 Find which source lambda gaps cover each target gap and the PREARAV1.168
CL fractions they contribute. PREARAV1.169
C PREARAV1.170
C At this point IXL(target_line) gives the index of the next source PREARAV1.171
C lambda line to the left of the target lambda line, wrapping round PREARAV1.172
C if the source grid is cyclic. This is also the index of the source PREARAV1.173
C gap in which the target line falls. Similarly, the index of the PREARAV1.174
C target line is also that of the target gap of which it is the PREARAV1.175
C left-hand limit. Therefore also IXL(target_gap+1, wrapping round JG030993.5
C if reqd.), is the index of the source gap which contains the JG030993.6
C right-hand limit of the target gap. For each target gap, we loop PREARAV1.178
C over all source gaps and find the fraction covered by each. Record PREARAV1.179
C the fraction and the source index in cumulative lists. If the PREARAV1.180
C source grid is not cyclic, parts of the target gap lying outside PREARAV1.181
C the source grid are neglected. PREARAV1.182
C PREARAV1.183
IXP=0 PREARAV1.184
DO IX2=1,GAPS_LAMBDA_TARG PREARAV1.185
IX=0 PREARAV1.186
IX2N=MOD(IX2,LINES_LAMBDA_TARG)+1 PREARAV1.187
BTARG=TEMP_TARG(IX2N)-TEMP_TARG(IX2) PREARAV1.188
IF (BTARG.LT.0.) THEN PREARAV1.189
BTARG=BTARG+360. PREARAV1.190
IX1N=IXL(IX2N)+LINES_LAMBDA_SRCE PREARAV1.191
ELSE PREARAV1.192
IX1N=IXL(IX2N) PREARAV1.193
ENDIF PREARAV1.194
DO IX1W=IXL(IX2),IX1N PREARAV1.195
IX1=MOD(IX1W-1,LINES_LAMBDA_SRCE)+1 PREARAV1.196
IF (CYCLIC_SRCE.OR.IX1.NE.LINES_LAMBDA_SRCE) THEN PREARAV1.197
IF (IX1W.EQ.IXL(IX2)) THEN PREARAV1.198
BLO=TEMP_TARG(IX2) PREARAV1.199
ELSE PREARAV1.200
BLO=TEMP_SRCE(IX1) PREARAV1.201
ENDIF PREARAV1.202
IF (IX1W.EQ.IX1N) THEN PREARAV1.203
BHI=TEMP_TARG(IX2N) PREARAV1.204
ELSE PREARAV1.205
BHI=TEMP_SRCE(MOD(IX1,LINES_LAMBDA_SRCE)+1) PREARAV1.206
ENDIF PREARAV1.207
IF (BHI.LT.BLO) BHI=BHI+360. PREARAV1.208
IF (ABS(BHI-BLO).GT.1E-7*ABS(BTARG)) THEN CJG1F401.28
IX=IX+1 PREARAV1.210
INDEX_LAMBDA(IXP+IX)=IX1 PREARAV1.211
FRAC_LAMBDA(IXP+IX)=(BHI-BLO)/BTARG PREARAV1.212
ENDIF PREARAV1.213
ENDIF PREARAV1.214
ENDDO PREARAV1.215
COUNT_LAMBDA(IX2)=IX PREARAV1.216
IXP=IXP+COUNT_LAMBDA(IX2) PREARAV1.217
ENDDO PREARAV1.218
C PREARAV1.219
CL 4 For each target phi line, find the index of the next source phi PREARAV1.220
CL line above. Comments as for section 2, without wrap-round. PREARAV1.221
C PREARAV1.222
LINES_PHI_SRCE=GAPS_PHI_SRCE+1 PREARAV1.223
LINES_PHI_TARG=GAPS_PHI_TARG+1 PREARAV1.224
DO IY2=1,LINES_PHI_TARG PREARAV1.225
IYT(IY2)=0 PREARAV1.226
DO IY1=1,LINES_PHI_SRCE PREARAV1.227
IF (PHI_TARG(IY2).LE.PHI_SRCE(IY1)) IYT(IY2)=IY1 PREARAV1.228
ENDDO PREARAV1.229
ENDDO PREARAV1.230
C PREARAV1.231
CL 5 Find which source phi gaps cover each target gap and the PREARAV1.232
CL fractions they contribute. Comments as for section 3, without PREARAV1.233
CL wrap-round. PREARAV1.234
C PREARAV1.235
IYP=0 PREARAV1.236
DO IY2=1,GAPS_PHI_TARG PREARAV1.237
IY=0 PREARAV1.238
IF (SPHERICAL) THEN CJG1F401.29
C Contain angle between +-90. There is no real area outside CJG1F401.30
C these limits on a sphere. CJG1F401.31
BTARG=SIN(MAX(MIN(PHI_TARG(IY2+1),90.),-90.)*PI_OVER_180) CJG1F401.32
& -SIN(MAX(MIN(PHI_TARG(IY2),90.),-90.)*PI_OVER_180) CJG1F401.33
ELSE CJG1F401.34
BTARG=PHI_TARG(IY2+1)-PHI_TARG(IY2) CJG1F401.35
ENDIF CJG1F401.36
DO IY1=IYT(IY2),IYT(IY2+1) PREARAV1.240
IF (.NOT.(IY1.EQ.0.OR.IY1.EQ.LINES_PHI_SRCE)) THEN PREARAV1.241
IF (IY1.EQ.IYT(IY2)) THEN PREARAV1.242
BLO=PHI_TARG(IY2) PREARAV1.243
ELSE PREARAV1.244
BLO=PHI_SRCE(IY1) PREARAV1.245
ENDIF PREARAV1.246
IF (IY1.EQ.IYT(IY2+1)) THEN PREARAV1.247
BHI=PHI_TARG(IY2+1) PREARAV1.248
ELSE PREARAV1.249
BHI=PHI_SRCE(IY1+1) PREARAV1.250
ENDIF PREARAV1.251
IF (SPHERICAL) THEN CJG1F401.37
BLO=MAX(MIN(BLO,90.),-90.) CJG1F401.38
BHI=MAX(MIN(BHI,90.),-90.) CJG1F401.39
ENDIF CJG1F401.40
IF (ABS(BHI-BLO).GT.1E-7*ABS(BTARG)) THEN CJG1F401.41
IY=IY+1 PREARAV1.253
INDEX_PHI(IYP+IY)=IY1 PREARAV1.254
C Both numerator and denominator in the following are -ve. PREARAV1.255
IF (SPHERICAL) THEN CJG1F401.42
FRAC_PHI(IYP+IY) CJG1F401.43
& =(SIN(BHI*PI_OVER_180)-SIN(BLO*PI_OVER_180))/BTARG CJG1F401.44
ELSE CJG1F401.45
FRAC_PHI(IYP+IY)=(BHI-BLO)/BTARG CJG1F401.46
ENDIF CJG1F401.47
ENDIF PREARAV1.257
ENDIF PREARAV1.258
ENDDO PREARAV1.259
COUNT_PHI(IY2)=IY PREARAV1.260
IYP=IYP+COUNT_PHI(IY2) PREARAV1.261
ENDDO PREARAV1.262
C PREARAV1.263
CL 6 For each target box, loop over contributing source boxes and PREARAV1.264
CL calculate the weights for each one, ignoring land boxes. PREARAV1.265
C PREARAV1.266
C After the first pass for each target box, go back and normalise PREARAV1.267
C the weights to compensate for land source boxes and any outside PREARAV1.268
C the source grid. Record the source box index and the weight in PREARAV1.269
C cumulative lists. PREARAV1.270
C PREARAV1.271
IP=0 PREARAV1.272
IYP=0 PREARAV1.273
DO IY2=1,GAPS_PHI_TARG PREARAV1.274
IXP=0 PREARAV1.275
DO IX2=1,GAPS_LAMBDA_TARG PREARAV1.276
I=0 PREARAV1.277
SUM=0 PREARAV1.278
DO IY=IYP+1,IYP+COUNT_PHI(IY2) PREARAV1.279
DO IX=IXP+1,IXP+COUNT_LAMBDA(IX2) PREARAV1.280
IF (MASK_SRCE(INDEX_LAMBDA(IX),INDEX_PHI(IY)) PREARAV1.281
& .EQV.WANT) THEN PREARAV1.282
I=I+1 PREARAV1.283
INDEX_SRCE(IP+I)=INDEX_LAMBDA(IX) PREARAV1.284
& +(INDEX_PHI(IY)-1)*GAPS_LAMBDA_SRCE PREARAV1.285
WEIGHT(IP+I)=FRAC_LAMBDA(IX)*FRAC_PHI(IY) PREARAV1.286
SUM=SUM+WEIGHT(IP+I) PREARAV1.287
ENDIF PREARAV1.288
ENDDO PREARAV1.289
ENDDO PREARAV1.290
COUNT_TARG(IX2,IY2)=I PREARAV1.291
BASE_TARG(IX2,IY2)=IP PREARAV1.292
DO I=1,COUNT_TARG(IX2,IY2) PREARAV1.293
WEIGHT(IP+I)=WEIGHT(IP+I)/SUM PREARAV1.294
ENDDO PREARAV1.295
IP=IP+COUNT_TARG(IX2,IY2) PREARAV1.296
IXP=IXP+COUNT_LAMBDA(IX2) PREARAV1.297
ENDDO PREARAV1.298
IYP=IYP+COUNT_PHI(IY2) PREARAV1.299
ENDDO PREARAV1.300
MAXL=IP PREARAV1.301
C PREARAV1.302
ICODE=0 PREARAV1.303
CMESSAGE=' ' PREARAV1.304
RETURN PREARAV1.305
END PREARAV1.306
*ENDIF PREARAV1.307