Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:40:50 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
10e456d88f Jean*0001 #include "FLT_OPTIONS.h"
c806179eb4 Alis*0002 
3992cf11bb Jean*0003 CBOP 0
                0004 C !ROUTINE: FLT_TRAJ
c806179eb4 Alis*0005 
3992cf11bb Jean*0006 C !INTERFACE:
10e456d88f Jean*0007       SUBROUTINE FLT_TRAJ (
eacecc7041 Jean*0008      I                      myTime, myIter, myThid )
c806179eb4 Alis*0009 
3992cf11bb Jean*0010 C     !DESCRIPTION:
                0011 C     *==========================================================*
                0012 C     | SUBROUTINE FLT_TRAJ
                0013 C     | o This routine samples the model state at float position
                0014 C     |   every flt_int_traj time steps and writes output.
                0015 C     *==========================================================*
c806179eb4 Alis*0016 
51ec3c32fe Jean*0017 C     !USES:
                0018       IMPLICIT NONE
                0019 C     == global variables ==
c806179eb4 Alis*0020 #include "SIZE.h"
10e456d88f Jean*0021 #include "EEPARAMS.h"
c806179eb4 Alis*0022 #include "PARAMS.h"
10e456d88f Jean*0023 #include "DYNVARS.h"
730d8469b1 Oliv*0024 #include "FLT_SIZE.h"
c806179eb4 Alis*0025 #include "FLT.h"
3992cf11bb Jean*0026 #include "FLT_BUFF.h"
ad773b031f Oliv*0027 #ifdef ALLOW_EXCH2
                0028 #include "W2_EXCH2_SIZE.h"
                0029 #include "W2_EXCH2_TOPOLOGY.h"
                0030 #endif
c806179eb4 Alis*0031 
3992cf11bb Jean*0032 C     !INPUT PARAMETERS:
                0033 C     myTime :: current time in simulation
                0034 C     myIter :: current iteration number
                0035 C     myThid :: my Thread Id number
10e456d88f Jean*0036       _RL myTime
eacecc7041 Jean*0037       INTEGER myIter, myThid
c806179eb4 Alis*0038 
3992cf11bb Jean*0039 C     !FUNCTIONS:
7fc4e95251 Jean*0040       _RL FLT_MAP_K2R
                0041       EXTERNAL FLT_MAP_K2R
51ec3c32fe Jean*0042 
3992cf11bb Jean*0043 C     !LOCAL VARIABLES:
                0044       INTEGER bi, bj, nFlds
10e456d88f Jean*0045       INTEGER ip, kp, ii
d5477ff298 Jean*0046       _RL ix, jy, i0x, j0y, xx, yy, zz
7fc4e95251 Jean*0047       _RL uu, vv, tt, ss, pp
c806179eb4 Alis*0048 
3992cf11bb Jean*0049       INTEGER imax
                0050       PARAMETER (imax=13)
10e456d88f Jean*0051       _RL tmp(imax)
3992cf11bb Jean*0052       _RL npart_read, npart_times
db913584c6 Jean*0053       _RS dummyRS(1)
0ad17d4ed9 Jean*0054       INTEGER fp, ioUnit, irecord
c806179eb4 Alis*0055       CHARACTER*(MAX_LEN_FNAM) fn
10e456d88f Jean*0056       CHARACTER*(MAX_LEN_MBUF) msgBuf
ad773b031f Oliv*0057 #ifdef ALLOW_EXCH2
                0058       INTEGER nT
                0059 #endif
3992cf11bb Jean*0060 CEOP
                0061 
                0062 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0063 
                0064 C--   set number of fields to write
                0065       nFlds = 0
                0066       IF ( flt_selectTrajOutp.GE.1 ) nFlds = nFlds + 8
                0067       IF ( flt_selectTrajOutp.GE.2 ) nFlds = nFlds + 5
                0068 
                0069 C--   check buffer size
                0070       IF ( nFlds.GT.fltBufDim ) THEN
d7e0a84259 Jean*0071          _BEGIN_MASTER(myThid)
3992cf11bb Jean*0072          WRITE(msgBuf,'(3(A,I4))') ' FLT_TRAJ: fltBufDim=', fltBufDim,
                0073      &                             ' too small (<', nFlds, ' )'
                0074          CALL PRINT_ERROR( msgBuf, myThid )
                0075          WRITE(msgBuf,'(2A)')     ' FLT_TRAJ: => increase fltBufDim',
                0076      &                            ' in "FLT_SIZE.h" & recompile'
                0077          CALL PRINT_ERROR( msgBuf, myThid )
d7e0a84259 Jean*0078          _END_MASTER(myThid)
3992cf11bb Jean*0079          CALL ALL_PROC_DIE( myThid )
                0080          STOP 'ABNORMAL END: S/R FLT_TRAJ'
                0081       ENDIF
                0082 
                0083       IF ( myIter.EQ.nIter0 .OR. flt_selectTrajOutp.LE.0 ) RETURN
                0084 
                0085 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0086 C--   Calculate position + other fields at float position and fill up IO-buffer
c806179eb4 Alis*0087 
                0088       DO bj=myByLo(myThid),myByHi(myThid)
51ec3c32fe Jean*0089        DO bi=myBxLo(myThid),myBxHi(myThid)
c806179eb4 Alis*0090 
ad773b031f Oliv*0091 #ifdef ALLOW_EXCH2
                0092          nT = W2_myTileList(bi,bj)
                0093          i0x = DFLOAT( exch2_txGlobalo(nT) - 1 )
                0094          j0y = DFLOAT( exch2_tyGlobalo(nT) - 1 )
                0095 #else
7fc4e95251 Jean*0096          i0x = DFLOAT( myXGlobalLo-1 + (bi-1)*sNx )
                0097          j0y = DFLOAT( myYGlobalLo-1 + (bj-1)*sNy )
ad773b031f Oliv*0098 #endif
10e456d88f Jean*0099          DO ip=1,npart_tile(bi,bj)
c806179eb4 Alis*0100 
d5477ff298 Jean*0101             ix = ipart(ip,bi,bj)
                0102             jy = jpart(ip,bi,bj)
                0103             CALL FLT_MAP_IJLOCAL2XY( xx, yy,
                0104      I                               ix, jy, bi,bj, myThid )
                0105             zz = FLT_MAP_K2R( kpart(ip,bi,bj),bi,bj,myThid )
51ec3c32fe Jean*0106             kp = NINT(kpart(ip,bi,bj))
3992cf11bb Jean*0107             tmp(1) = npart(ip,bi,bj)
                0108             tmp(2) = myTime
                0109             tmp(3) = xx
                0110             tmp(4) = yy
                0111             tmp(5) = zz
                0112             tmp(6) = ix + i0x
                0113             tmp(7) = jy + j0y
                0114             tmp(8) = kpart(ip,bi,bj)
                0115 
                0116             IF ( ( flt_selectTrajOutp.GE.2 )   .AND.
                0117      &           ( myTime.GE.tstart(ip,bi,bj)) .AND.
10e456d88f Jean*0118      &           ( tend(ip,bi,bj).EQ.-1. .OR. myTime.LE.tend(ip,bi,bj))
                0119      &         ) THEN
                0120               IF ( kp.LT.1 .OR. kp.GT.Nr ) THEN
                0121                 WRITE(msgBuf,'(2A,I8)') '** WARNING ** FLT_TRAJ: ',
                0122      &            ' illegal value for kp=',kp
                0123                 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
                0124      &                              SQUEEZE_RIGHT, myThid )
                0125                 WRITE(msgBuf,'(A,1P5E20.13)')
3992cf11bb Jean*0126      &            ' FLT_TRAJ: ', (flt_io_buff(ii,ip,bi,bj),ii=1,5)
10e456d88f Jean*0127                 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
                0128      &                              SQUEEZE_RIGHT, myThid )
                0129 c               CALL PRINT_ERROR( msgBuf, myThid )
                0130 c               STOP 'ABNORMAL END: S/R FLT_TRAJ'
                0131 C-- jmc: not sure if this is right but added to avoid Pb in FLT_BILINEAR:
51ec3c32fe Jean*0132                 kp = MIN( MAX(kp,1), Nr)
                0133               ENDIF
7fc4e95251 Jean*0134               CALL FLT_BILINEAR  (ix,jy,uu,uVel,  kp,1,bi,bj,myThid)
                0135               CALL FLT_BILINEAR  (ix,jy,vv,vVel,  kp,2,bi,bj,myThid)
                0136               CALL FLT_BILINEAR2D(ix,jy,pp,etaN,     0,bi,bj,myThid)
                0137               CALL FLT_BILINEAR  (ix,jy,tt,theta, kp,0,bi,bj,myThid)
                0138               CALL FLT_BILINEAR  (ix,jy,ss,salt,  kp,0,bi,bj,myThid)
                0139               tmp( 9) = pp
                0140               tmp(10) = uu
                0141               tmp(11) = vv
                0142               tmp(12) = tt
                0143               tmp(13) = ss
3992cf11bb Jean*0144             ELSEIF ( flt_selectTrajOutp.GE.2 ) THEN
7fc4e95251 Jean*0145               tmp( 9) = flt_nan
                0146               tmp(10) = flt_nan
                0147               tmp(11) = flt_nan
                0148               tmp(12) = flt_nan
                0149               tmp(13) = flt_nan
10e456d88f Jean*0150             ENDIF
                0151 
3992cf11bb Jean*0152             DO ii=1,nFlds
                0153               flt_io_buff(ii,ip,bi,bj) = tmp(ii)
                0154             ENDDO
                0155 
                0156          ENDDO
                0157 
                0158        ENDDO
                0159       ENDDO
                0160 
                0161 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0162 
                0163 C--   Write shared buffer to file
                0164 
                0165       _BARRIER
                0166       _BEGIN_MASTER(myThid)
                0167 
                0168       fn = 'float_trajectories'
                0169       fp = writeBinaryPrec
                0170 
                0171       DO bj=1,nSy
                0172        DO bi=1,nSx
                0173 
                0174 C (1) read actual number floats from file (if exists)
                0175          ioUnit = -2
                0176          CALL MDS_READVEC_LOC(  fn, fp, ioUnit, 'RL', nFlds,
                0177      O                          tmp, dummyRS,
                0178      I                          bi, bj, 1, myThid )
                0179          IF ( ioUnit.GT. 0 ) THEN
                0180             npart_read  = tmp(1)
                0181             npart_times = tmp(5)
                0182             ii = NINT(tmp(7))
                0183 C-       for backward compatibility with old trajectory files:
                0184             IF ( ii.EQ.0 ) ii = 13
                0185             IF ( ii.NE.nFlds ) THEN
                0186               WRITE(msgBuf,'(A,I4,A)')
                0187      &            'FLT_TRAJ: nFlds=', nFlds,' different from'
                0188               CALL PRINT_ERROR( msgBuf, myThid )
                0189               WRITE(msgBuf,'(3A,I4,A)')
                0190      &            'previous file (',fn(1:18),') value =',ii
                0191               CALL PRINT_ERROR( msgBuf, myThid )
d7e0a84259 Jean*0192               CALL ALL_PROC_DIE( 0 )
3992cf11bb Jean*0193               STOP 'ABNORMAL END: S/R FLT_TRAJ'
                0194             ENDIF
                0195 C-       close the read-unit (safer to use a different unit for writing)
                0196             CLOSE( ioUnit )
                0197          ELSE
                0198             npart_read  = 0.
                0199             npart_times = 0.
                0200             tmp(2)      = myTime
                0201          ENDIF
                0202 
                0203 C (2) write new actual number floats and time axis into file
                0204 C-    the standard routine mds_writevec_loc can be used here
                0205 
                0206 C     total number of records in this file
                0207          tmp(1) = DBLE(npart_tile(bi,bj))+npart_read
                0208 C     first time of writing floats (do not change when written)
                0209 c        tmp(2) = tmp(2)
                0210 C     current time
                0211          tmp(3) = myTime
                0212 C     timestep
                0213          tmp(4) = flt_int_traj
                0214 C     total number of timesteps
                0215          tmp(5) = npart_times + 1.
                0216 C     total number of floats
                0217          tmp(6) = max_npart
                0218 C     total number of fields
                0219          tmp(7) = nFlds
                0220          DO ii=8,nFlds
                0221            tmp(ii) = 0.
                0222          ENDDO
                0223          ioUnit = -1
                0224          CALL MDS_WRITEVEC_LOC( fn, fp, ioUnit, 'RL', nFlds,
                0225      &                          tmp, dummyRS,
                0226      &                          bi, bj, -1, myIter, myThid )
                0227 
                0228          DO ip=1,npart_tile(bi,bj)
                0229 C (3) write float positions into file
55f764277b Jean*0230             irecord = npart_read+ip+1
                0231             IF ( ip.NE.npart_tile(bi,bj) ) irecord = -irecord
3992cf11bb Jean*0232             CALL MDS_WRITEVEC_LOC( fn, fp, ioUnit, 'RL', nFlds,
                0233      I                             flt_io_buff(1,ip,bi,bj), dummyRS,
                0234      I                             bi, bj, irecord, myIter, myThid )
10e456d88f Jean*0235          ENDDO
55f764277b Jean*0236          CLOSE( ioUnit )
c806179eb4 Alis*0237 
                0238        ENDDO
51ec3c32fe Jean*0239       ENDDO
c806179eb4 Alis*0240 
3992cf11bb Jean*0241       _END_MASTER(myThid)
                0242       _BARRIER
                0243 
10e456d88f Jean*0244       RETURN
                0245       END