File indexing completed on 2021-11-06 05:17:38 UTC
view on githubraw file Latest commit ebcc29af on 2021-10-24 15:53:13 UTC
9c576a35a3 An T*0001 #include "ECCO_OPTIONS.h"
0002
ebcc29af97 Jean*0003 SUBROUTINE COST_GENCOST_TRANSP( myThid )
9c576a35a3 An T*0004
ebcc29af97 Jean*0005
0006
0007
0008
9c576a35a3 An T*0009
ebcc29af97 Jean*0010 IMPLICIT NONE
9c576a35a3 An T*0011
ebcc29af97 Jean*0012
9c576a35a3 An T*0013 #include "EEPARAMS.h"
0014 #include "SIZE.h"
0015 #include "PARAMS.h"
0016 #include "GRID.h"
0017 #ifdef ALLOW_CAL
0018 # include "cal.h"
0019 #endif
0020 #ifdef ALLOW_ECCO
13d362b8c1 Ou W*0021 # include "ECCO_SIZE.h"
0022 # include "ECCO.h"
9c576a35a3 An T*0023 #endif
0024
ebcc29af97 Jean*0025
0026 INTEGER myThid
9c576a35a3 An T*0027
0028 #ifdef ALLOW_GENCOST_CONTRIBUTION
ebcc29af97 Jean*0029
0030 INTEGER ILNBLNK
0031 EXTERNAL ILNBLNK
9c576a35a3 An T*0032
ebcc29af97 Jean*0033
0034 INTEGER nrecloc, localrec
0035 INTEGER localstartdate(4)
9c576a35a3 An T*0036
ebcc29af97 Jean*0037 _RL myobs (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
0038 _RL mybar (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
0039 _RL localdif(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
0040 _RL difmask (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
0041 _RL localtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
9c576a35a3 An T*0042
ebcc29af97 Jean*0043 _RL dummy, facCost, facNum
9c576a35a3 An T*0044 _RL localperiod
0045 _RL spminloc, spmaxloc, spzeroloc
ebcc29af97 Jean*0046 _RL tmpMeanTile(nSx,nSy), tmpNumTile(nSx,nSy)
0047 _RL tmpMeanGlo, tmpNumGlo
0048
0049 INTEGER kgen(NGENCOST3D)
0050 INTEGER bi, bj
0051 INTEGER k
0052 INTEGER obsrec, irec
0053 INTEGER il,k2
0054 INTEGER icount, icount_transp
0055 LOGICAL dosumsq, dovarwei, doreadobs
0056 LOGICAL exst
0057
0058 INTEGER preproc_i(NGENPPROC)
9c576a35a3 An T*0059 _RL preproc_r(NGENPPROC)
ebcc29af97 Jean*0060 CHARACTER*(MAX_LEN_FNAM) mybarfile
0061 CHARACTER*(MAX_LEN_FNAM) preproc(NGENPPROC)
0062 CHARACTER*(MAX_LEN_FNAM) preproc_c(NGENPPROC)
0063 CHARACTER*(MAX_LEN_MBUF) msgBuf
0064 CHARACTER*(128) fname1, fname0
0065
9c576a35a3 An T*0066
ebcc29af97 Jean*0067
9c576a35a3 An T*0068
ebcc29af97 Jean*0069
0070 facNum = nSx*nPx
0071 facNum = 1. _d 0 / facNum
9c576a35a3 An T*0072
ebcc29af97 Jean*0073
8af4a14f0f An T*0074 do k=1,NGENCOST3D
0075 kgen(k)=0
0076 enddo
0077
ebcc29af97 Jean*0078
8af4a14f0f An T*0079 write(msgbuf,'(A)') 'Inside cost_gencost_transp:'
0080 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
0081 & SQUEEZE_RIGHT, myThid )
0082
0083 icount_transp=0
9c576a35a3 An T*0084 do k=1,NGENCOST
0085 if ( (gencost_name(k)(1:6).EQ.'transp').AND.
8af4a14f0f An T*0086 & (using_gencost(k)) ) then
0087 icount_transp=icount_transp+1
0088 kgen(icount_transp)=k
ebcc29af97 Jean*0089 il=ILNBLNK(gencost_barfile(kgen(icount_transp)))
8af4a14f0f An T*0090 write(msgbuf,'(A,i4,A,A)') 'Cost ',kgen(icount_transp),
0091 & ': ',gencost_barfile(kgen(icount_transp))(1:il)
0092 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
0093 & SQUEEZE_RIGHT, myThid )
0094
0095 endif
9c576a35a3 An T*0096 enddo
0097
ebcc29af97 Jean*0098
8af4a14f0f An T*0099
0100 do icount=1,icount_transp
ebcc29af97 Jean*0101 if (kgen(icount).NE.0) then
9c576a35a3 An T*0102
ebcc29af97 Jean*0103
9c576a35a3 An T*0104
ebcc29af97 Jean*0105
0106 DO bj = myByLo(myThid), myByHi(myThid)
0107 DO bi = myBxLo(myThid), myBxHi(myThid)
8af4a14f0f An T*0108 objf_gencost(bi,bj,kgen(icount))=0. _d 0
0109 num_gencost(bi,bj,kgen(icount))=0. _d 0
9c576a35a3 An T*0110 ENDDO
0111 ENDDO
0112
ebcc29af97 Jean*0113
9c576a35a3 An T*0114 nrecloc=0
8af4a14f0f An T*0115 nrecloc=gencost_nrec(kgen(icount))
9c576a35a3 An T*0116
ebcc29af97 Jean*0117
9c576a35a3 An T*0118 if(nrecloc.gt.0) then
0119
ebcc29af97 Jean*0120 facCost = nrecloc
0121 facCost = 1. _d 0 / facCost
9c576a35a3 An T*0122
0123 localperiod=0.
8af4a14f0f An T*0124 localperiod=gencost_period(kgen(icount))
0125 dummy=gencost_dummy(kgen(icount))
0126 spminloc=gencost_spmin(kgen(icount))
0127 spmaxloc=gencost_spmax(kgen(icount))
0128 spzeroloc=gencost_spzero(kgen(icount))
0129
ebcc29af97 Jean*0130
0131
9c576a35a3 An T*0132 dosumsq=.FALSE.
0133 dovarwei=.FALSE.
0134 do k2 = 1, NGENPPROC
8af4a14f0f An T*0135 preproc(k2)=gencost_preproc(k2,kgen(icount))
0136 preproc_i(k2)=gencost_preproc_i(k2,kgen(icount))
0137 preproc_c(k2)=gencost_preproc_c(k2,kgen(icount))
0138 preproc_r(k2)=gencost_preproc_r(k2,kgen(icount))
9c576a35a3 An T*0139 if (preproc(k2).EQ.'variaweight') dovarwei=.TRUE.
8af4a14f0f An T*0140 if (preproc(k2)(1:7).EQ.'dosumsq') dosumsq=.TRUE.
9c576a35a3 An T*0141 enddo
0142
ebcc29af97 Jean*0143
0144 il=ILNBLNK(gencost_name(kgen(icount)))
0145 write(msgBuf,'(3A,L5,2A)')
8af4a14f0f An T*0146 & 'Cost ',gencost_name(kgen(icount))(1:il),
0147 & ' dosumsq: ',dosumsq,' preproc(1): ',preproc(1)(1:7)
fce41e6d01 An T*0148 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
0149 & SQUEEZE_RIGHT, myThid )
9c576a35a3 An T*0150
ebcc29af97 Jean*0151
8af4a14f0f An T*0152 mybarfile=gencost_barfile(kgen(icount))
0153
ebcc29af97 Jean*0154
9c576a35a3 An T*0155 doreadobs=.FALSE.
8af4a14f0f An T*0156 if( .not. gencost_datafile(kgen(icount)).eq.' ') then
0157
ebcc29af97 Jean*0158 write(msgBuf,'(A)')
8af4a14f0f An T*0159 & '**WARNING** S/R COST_GENCOST_TRANSP: gencost_datafile '
0160 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
0161 & SQUEEZE_RIGHT, myThid )
ebcc29af97 Jean*0162 write(msgBuf,'(A)')
8af4a14f0f An T*0163 & 'are currently ignored. Adjust the S/R to add code.'
0164 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
0165 & SQUEEZE_RIGHT, myThid )
0166 endif
9c576a35a3 An T*0167
ebcc29af97 Jean*0168
9c576a35a3 An T*0169
ebcc29af97 Jean*0170
8af4a14f0f An T*0171
ebcc29af97 Jean*0172
0173
9c576a35a3 An T*0174
ebcc29af97 Jean*0175
9c576a35a3 An T*0176 do irec = 1,nrecloc
8af4a14f0f An T*0177
ebcc29af97 Jean*0178
8af4a14f0f An T*0179
ebcc29af97 Jean*0180
fce41e6d01 An T*0181 call ecco_zero(mybar,Nr,zeroRL,myThid)
9c576a35a3 An T*0182
ebcc29af97 Jean*0183
9c576a35a3 An T*0184 exst=.FALSE.
0185 call ecco_zero(localtmp,Nr,zeroRL,myThid)
736e27304c Timo*0186 call cost_gencal(mybarfile,gencost_mask(kgen(icount)),
fce41e6d01 An T*0187 & irec,localstartdate,localperiod,fname1,
ebcc29af97 Jean*0188 & fname0,localrec,obsrec,exst,myThid)
0189 call cost_genread(fname1,mybar,localtmp,irec,Nr,Nr,
9c576a35a3 An T*0190 & nrecloc,preproc,preproc_c,preproc_i,preproc_r,
ebcc29af97 Jean*0191 & dummy,myThid)
9c576a35a3 An T*0192
ebcc29af97 Jean*0193
8af4a14f0f An T*0194
ebcc29af97 Jean*0195
0196
fce41e6d01 An T*0197
0198 call ecco_zero(myobs,Nr,zeroRL,myThid)
0199
9c576a35a3 An T*0200
ebcc29af97 Jean*0201
8af4a14f0f An T*0202
ebcc29af97 Jean*0203
0204 DO bj = myByLo(myThid), myByHi(myThid)
0205 DO bi = myBxLo(myThid), myBxHi(myThid)
0206 tmpMeanTile(bi,bj) = 0. _d 0
0207 tmpNumTile(bi,bj) = 0. _d 0
0208 ENDDO
9c576a35a3 An T*0209 ENDDO
0210
ebcc29af97 Jean*0211
0212
0213 call ecco_zero(localtmp,Nr,oneRL,myThid)
0214 call ecco_zero(localdif,Nr,zeroRL,myThid)
0215 call ecco_zero(difmask,Nr,zeroRL,myThid)
9c576a35a3 An T*0216
ebcc29af97 Jean*0217
0218
9c576a35a3 An T*0219 call ecco_diffmsk(
11c3150c71 Mart*0220 I mybar, myobs, localtmp,
0221 I Nr, Nr,spminloc, spmaxloc, spzeroloc,
fce41e6d01 An T*0222 O localdif, difmask,
9c576a35a3 An T*0223 I myThid )
0224
ebcc29af97 Jean*0225 call ecco_zero(localtmp,Nr,oneRL,myThid)
9c576a35a3 An T*0226 call ecco_addcost(
ebcc29af97 Jean*0227 I localdif,localtmp,difmask,Nr,Nr,dosumsq,
9c576a35a3 An T*0228 O tmpMeanTile,tmpNumTile,
ebcc29af97 Jean*0229 I myThid)
9c576a35a3 An T*0230
ebcc29af97 Jean*0231
0232
0233
8af4a14f0f An T*0234
ebcc29af97 Jean*0235
8af4a14f0f An T*0236
ebcc29af97 Jean*0237
8af4a14f0f An T*0238
ebcc29af97 Jean*0239
0240
9c576a35a3 An T*0241 tmpMeanGlo = 0. _d 0
0242 tmpNumGlo = 0. _d 0
ebcc29af97 Jean*0243 il=ILNBLNK(gencost_barfile(kgen(icount)))
9c576a35a3 An T*0244 CALL GLOBAL_SUM_TILE_RL( tmpMeanTile, tmpMeanGlo, myThid )
0245 CALL GLOBAL_SUM_TILE_RL( tmpNumTile, tmpNumGlo, myThid )
8af4a14f0f An T*0246 WRITE(msgBuf,'(2A,I3,A,1PE21.14,1PE21.14)')
0247 & 'globalsum transp ',gencost_barfile(kgen(icount))(1:il),
0248 & irec,' ',tmpMeanGlo,tmpNumGlo
9c576a35a3 An T*0249 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
0250 & SQUEEZE_RIGHT, myThid )
0251
ebcc29af97 Jean*0252
0253 DO bj = myByLo(myThid), myByHi(myThid)
0254 DO bi = myBxLo(myThid), myBxHi(myThid)
0255 objf_gencost(bi,bj,kgen(icount))=
0256 & objf_gencost(bi,bj,kgen(icount))+tmpMeanTile(bi,bj)
0257 ENDDO
9c576a35a3 An T*0258 ENDDO
0259
0260 enddo
0261
ebcc29af97 Jean*0262
0263
0264
0265 DO bj = myByLo(myThid), myByHi(myThid)
0266 DO bi = myBxLo(myThid), myBxHi(myThid)
0267 objf_gencost(bi,bj,kgen(icount))=
0268 & objf_gencost(bi,bj,kgen(icount))*facCost
0269 num_gencost(bi,bj,kgen(icount))=nrecloc*facNum
0270 ENDDO
0271 ENDDO
9c576a35a3 An T*0272
0273 endif
ebcc29af97 Jean*0274 endif
8af4a14f0f An T*0275 enddo
9c576a35a3 An T*0276
0277 #endif /* ALLOW_GENCOST_CONTRIBUTION */
0278
0279 RETURN
ebcc29af97 Jean*0280 END