tutorial_global_oce_optim/code_ad/cost_weights.F

5ed655852f Jean*0001 #include "COST_OPTIONS.h"
ec31b3f7d5 Jean*0002 #ifdef ALLOW_CTRL
                0003 # include "CTRL_OPTIONS.h"
                0004 #endif
f1ca098265 Davi*0005
11c3150c71 Mart*0006 CBOP
                0007 C     !ROUTINE: COST_WEIGHTS
                0008 C     !INTERFACE:
01a8008aaa Davi*0009       SUBROUTINE COST_WEIGHTS( myThid )
11c3150c71 Mart*0010 C     !DESCRIPTION: \bv
                0011 C     *==========================================================*
                0012 C     | SUBROUTINE COST_WEIGHTS
                0013 C     | o Set weights used in the cost function
                0014 C     *==========================================================*
391bdbafb2 Jean*0015 C
11c3150c71 Mart*0016 C     \ev
                0017 C     !USES:
01a8008aaa Davi*0018       IMPLICIT NONE
11c3150c71 Mart*0019 C     == Global variables ===
f1ca098265 Davi*0020 #include "EEPARAMS.h"
                0021 #include "SIZE.h"
                0022 #include "PARAMS.h"
                0023 #include "GRID.h"
cf705a6c8e Mart*0024 #ifdef ALLOW_OPENAD
                0025 # include "FFIELDS.h"
                0026 #else
                0027 # include "cost_local.h"
11c3150c71 Mart*0028 #endif
f1ca098265 Davi*0029
11c3150c71 Mart*0030 C     !INPUT/OUTPUT PARAMETERS:
01a8008aaa Davi*0031       INTEGER  myThid
f1ca098265 Davi*0032
11c3150c71 Mart*0033 C     !FUNCTIONS:
391bdbafb2 Jean*0034       INTEGER  MDS_RECLEN
                0035       EXTERNAL MDS_RECLEN
f1ca098265 Davi*0036
11c3150c71 Mart*0037 C     !LOCAL VARIABLES:
01a8008aaa Davi*0038       INTEGER bi,bj
                0039       INTEGER i,j,k
391bdbafb2 Jean*0040       INTEGER iUnit, length_of_rec
f1ca098265 Davi*0041
                0042       _RL dummy
01a8008aaa Davi*0043       _RL wti(Nr)
                0044       REAL*8 tmpwti(Nr)
                0045       CHARACTER*(MAX_LEN_MBUF) msgBuf
11c3150c71 Mart*0046 CEOP
f1ca098265 Davi*0047
391bdbafb2 Jean*0048 C--   Initialize variance (weight) fields.
01a8008aaa Davi*0049       DO k = 1,Nr
11c3150c71 Mart*0050        wti(k) = 0. _d 0
01a8008aaa Davi*0051       ENDDO
11c3150c71 Mart*0052       DO bj = myByLo(myThid),myByHi(myThid)
                0053        DO bi = myBxLo(myThid),myBxHi(myThid)
                0054         DO j = 1-OLy,sNy+OLy
                0055          DO i = 1-OLx,sNx+OLx
                0056           whfluxm(i,j,bi,bj)= 0. _d 0
                0057          ENDDO
                0058         ENDDO
                0059         DO k = 1,Nr
                0060          wtheta(k,bi,bj) = 0. _d 0
f1ca098265 Davi*0061         ENDDO
11c3150c71 Mart*0062        ENDDO
f1ca098265 Davi*0063       ENDDO
                0064
391bdbafb2 Jean*0065 C--   Read error information and set up weight matrices.
01a8008aaa Davi*0066
f1ca098265 Davi*0067 #ifdef ALLOW_COST_TEMP
391bdbafb2 Jean*0068 C  Temperature weights for cost function
01a8008aaa Davi*0069        _BEGIN_MASTER(myThid)
                0070        CALL MDSFINDUNIT( iUnit, myThid )
391bdbafb2 Jean*0071        length_of_rec = MDS_RECLEN( precFloat64, Nr, myThid )
                0072        OPEN( iUnit, FILE='Err_levitus_15layer.bin', STATUS='OLD',
                0073      &       FORM='UNFORMATTED',ACCESS='DIRECT',RECL=length_of_rec )
                0074        READ(iUnit,rec=1) tmpwti
01a8008aaa Davi*0075        CLOSE(iUnit)
                0076 #ifdef _BYTESWAPIO
                0077        CALL MDS_BYTESWAPR8( Nr, tmpwti )
                0078 #endif
                0079        _END_MASTER(myThid)
                0080        _BARRIER
                0081
                0082        DO k=1,Nr
                0083          wti(k) = tmpwti(k)
                0084        ENDDO
                0085        WRITE(msgBuf,'(3A)') 'S/R COST_WEIGHTS:',
                0086      &  ' Temperature weights loaded from: ','Err_levitus_15layer.bin'
                0087        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0088      &                     SQUEEZE_RIGHT , myThid )
                0089
4db199adb5 Jean*0090 c     print*,'Weights for temperature: wti', (wti(k),k=1,nr)
f1ca098265 Davi*0091
11c3150c71 Mart*0092       DO bj = myByLo(myThid),myByHi(myThid)
                0093         DO bi = myBxLo(myThid),myBxHi(myThid)
01a8008aaa Davi*0094           DO k = 1, Nr
f1ca098265 Davi*0095                wtheta(k,bi,bj) = 1. _d 0/wti(k)/wti(k)
                0096           ENDDO
                0097         ENDDO
                0098       ENDDO
ec31b3f7d5 Jean*0099 #endif /* ALLOW_COST_TEMP */
391bdbafb2 Jean*0100
                0101 C--   Then the hflux weights :
                0102
11c3150c71 Mart*0103 #ifdef ALLOW_COST_HFLUXM
391bdbafb2 Jean*0104       CALL READ_REC_3D_RL( 'Err_hflux.bin', precFloat64, 1,
                0105      &                      whfluxm, 1, 0, myThid )
f1ca098265 Davi*0106       _EXCH_XY_RL(whfluxm   , myThid )
11c3150c71 Mart*0107       DO bj = myByLo(myThid),myByHi(myThid)
                0108         DO bi = myBxLo(myThid),myBxHi(myThid)
                0109           DO j = 1-OLy,sNy+OLy
                0110             DO i = 1-OLx,sNx+OLx
4db199adb5 Jean*0111 c            print*,'Uncertainties for Heat Flux',i,j,whfluxm(i,j,bi,bj)
f1ca098265 Davi*0112              IF (whfluxm(i,j,bi,bj) .NE. 0. _d 0) THEN
                0113                  whfluxm(i,j,bi,bj) = 1. _d 0 /whfluxm(i,j,bi,bj)
                0114      &                                        /whfluxm(i,j,bi,bj)
                0115              ELSE
                0116                  whfluxm(i,j,bi,bj) = 1. _d 0
                0117              ENDIF
                0118             ENDDO
                0119           ENDDO
                0120         ENDDO
                0121       ENDDO
                0122 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
391bdbafb2 Jean*0123       CALL ACTIVE_WRITE_XY('whfluxm',whfluxm,1,0,myThid,dummy)
f1ca098265 Davi*0124 #endif
11c3150c71 Mart*0125 #endif /* ALLOW_COST_HFLUXM */
391bdbafb2 Jean*0126       RETURN
                0127       END