New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
dtatsd.F90 in NEMO/trunk/tests/ISOMIP+/MY_SRC – NEMO

source: NEMO/trunk/tests/ISOMIP+/MY_SRC/dtatsd.F90

Last change on this file was 14995, checked in by mathiot, 3 years ago

ticket #2669 : merge ticket2669_isf_fluxes into trunk

File size: 14.0 KB
RevLine 
[11889]1MODULE dtatsd
2   !!======================================================================
3   !!                     ***  MODULE  dtatsd  ***
4   !! Ocean data  :  read ocean Temperature & Salinity Data from gridded data
5   !!======================================================================
6   !! History :  OPA  ! 1991-03  ()  Original code
7   !!             -   ! 1992-07  (M. Imbard)
[14995]8   !!            8.0  ! 1999-10  (M.A. Foujols, M. Imbard)  NetCDF FORMAT
9   !!   NEMO     1.0  ! 2002-06  (G. Madec)  F90: Free form and module
[11889]10   !!            3.3  ! 2010-10  (C. Bricaud, S. Masson)  use of fldread
11   !!            3.4  ! 2010-11  (G. Madec, C. Ethe) Merge of dtatem and dtasal + remove CPP keys
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   dta_tsd      : read and time interpolated ocean Temperature & Salinity Data
16   !!----------------------------------------------------------------------
17   USE oce             ! ocean dynamics and tracers
18   USE phycst          ! physical constants
19   USE dom_oce         ! ocean space and time domain
[14090]20   USE domtile
[11889]21   USE fldread         ! read input fields
22   !
23   USE in_out_manager  ! I/O manager
24   USE lib_mpp         ! MPP library
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   dta_tsd_init   ! called by opa.F90
30   PUBLIC   dta_tsd        ! called by istate.F90 and tradmp.90
31
32   !                                  !!* namtsd  namelist : Temperature & Salinity Data *
33   LOGICAL , PUBLIC ::   ln_tsd_init   !: T & S data flag
34   LOGICAL , PUBLIC ::   ln_tsd_dmp    !: internal damping toward input data flag
35
36   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tsdini ! structure of input SST (file informations, fields read)
37   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tsddmp ! structure of input SST (file informations, fields read)
38
[12905]39   !! * Substitutions
40#  include "do_loop_substitute.h90"
[11889]41   !!----------------------------------------------------------------------
42   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
43   !! $Id: dtatsd.F90 10213 2018-10-23 14:40:09Z aumont $
44   !! Software governed by the CeCILL license (see ./LICENSE)
45   !!----------------------------------------------------------------------
46CONTAINS
47
48   SUBROUTINE dta_tsd_init( ld_tradmp )
49      !!----------------------------------------------------------------------
50      !!                   ***  ROUTINE dta_tsd_init  ***
[14995]51      !!
52      !! ** Purpose :   initialisation of T & S input data
53      !!
[11889]54      !! ** Method  : - Read namtsd namelist
[14995]55      !!              - allocates T & S data structure
[11889]56      !!----------------------------------------------------------------------
57      LOGICAL, INTENT(in), OPTIONAL ::   ld_tradmp   ! force the initialization when tradp is used
58      !
59      INTEGER ::   ios, ierr0, ierr1, ierr2, ierr3   ! local integers
60      !!
61      CHARACTER(len=100)            ::   cn_dir          ! Root directory for location of ssr files
62      TYPE(FLD_N), DIMENSION( jpts) ::   slf_i           ! array of namelist informations on the fields to read
63      TYPE(FLD_N)                   ::   sn_tem, sn_sal
64      TYPE(FLD_N)                   ::   sn_dmpt, sn_dmps
65      !!
66      NAMELIST/namtsd/   ln_tsd_init, ln_tsd_dmp, cn_dir, sn_tem, sn_sal, sn_dmpt, sn_dmps
67      !!----------------------------------------------------------------------
68      !
69      !  Initialisation
70      ierr0 = 0  ;  ierr1 = 0  ;  ierr2 = 0  ;  ierr3 = 0
71      !
72      READ  ( numnam_ref, namtsd, IOSTAT = ios, ERR = 901)
[12077]73901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtsd in reference namelist' )
[11889]74      READ  ( numnam_cfg, namtsd, IOSTAT = ios, ERR = 902 )
[12077]75902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtsd in configuration namelist' )
[11889]76      IF(lwm) WRITE ( numond, namtsd )
77
78      IF( PRESENT( ld_tradmp ) )   ln_tsd_dmp = .TRUE.     ! forces the initialization when tradmp is used
[14995]79
[11889]80      IF(lwp) THEN                  ! control print
81         WRITE(numout,*)
82         WRITE(numout,*) 'dta_tsd_init : Temperature & Salinity data '
83         WRITE(numout,*) '~~~~~~~~~~~~ '
84         WRITE(numout,*) '   Namelist namtsd'
85         WRITE(numout,*) '      Initialisation of ocean T & S with T &S input data   ln_tsd_init = ', ln_tsd_init
86         WRITE(numout,*) '      damping of ocean T & S toward T &S input data        ln_tsd_dmp  = ', ln_tsd_dmp
87         WRITE(numout,*)
88         IF( .NOT.ln_tsd_init .AND. .NOT.ln_tsd_dmp ) THEN
89            WRITE(numout,*)
90            WRITE(numout,*) '   ===>>   T & S data not used'
91         ENDIF
92      ENDIF
93      !
94      IF( ln_rstart .AND. ln_tsd_init ) THEN
95         CALL ctl_warn( 'dta_tsd_init: ocean restart and T & S data intialisation, ',   &
96            &           'we keep the restart T & S values and set ln_tsd_init to FALSE' )
97         ln_tsd_init = .FALSE.
98      ENDIF
99      !
100      !                             ! allocate the arrays (if necessary)
101      IF( ln_tsd_init ) THEN
102         !
103         ALLOCATE( sf_tsdini(jpts), STAT=ierr0 )
104         IF( ierr0 > 0 ) THEN
105            CALL ctl_stop( 'dta_tsd_init: unable to allocate sf_tsd structure' )   ;   RETURN
106         ENDIF
107         !
108                                ALLOCATE( sf_tsdini(jp_tem)%fnow(jpi,jpj,jpk)   , STAT=ierr0 )
109         IF( sn_tem%ln_tint )   ALLOCATE( sf_tsdini(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 )
110                                ALLOCATE( sf_tsdini(jp_sal)%fnow(jpi,jpj,jpk)   , STAT=ierr2 )
111         IF( sn_sal%ln_tint )   ALLOCATE( sf_tsdini(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 )
112         !
113         IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN
114            CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' )   ;   RETURN
115         ENDIF
116         !                         ! fill sf_tsd with sn_tem & sn_sal and control print
117         slf_i(jp_tem) = sn_tem   ;   slf_i(jp_sal) = sn_sal
118         CALL fld_fill( sf_tsdini, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity data', 'namtsd', no_print )
119         !
120      END IF
121
122      IF( ln_tsd_dmp ) THEN
123         !
124         ALLOCATE( sf_tsddmp(jpts), STAT=ierr0 )
125         IF( ierr0 > 0 ) THEN
126            CALL ctl_stop( 'dta_tsd_init: unable to allocate sf_tsddmp structure' )   ;   RETURN
127         ENDIF
128         !
129         ! dmp file
130                                 ALLOCATE( sf_tsddmp(jp_tem)%fnow(jpi,jpj,jpk)   , STAT=ierr0 )
131         IF( sn_dmpt%ln_tint )   ALLOCATE( sf_tsddmp(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 )
132                                 ALLOCATE( sf_tsddmp(jp_sal)%fnow(jpi,jpj,jpk)   , STAT=ierr2 )
133         IF( sn_dmps%ln_tint )   ALLOCATE( sf_tsddmp(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 )
134         !
135         IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN
136            CALL ctl_stop( 'dta_tsd : unable to allocate T & S dmp data arrays' )   ;   RETURN
137         ENDIF
138         !
139         !                         ! fill sf_tsd with sn_tem & sn_sal and control print
140         slf_i(jp_tem) = sn_dmpt   ;   slf_i(jp_sal) = sn_dmps
141         CALL fld_fill( sf_tsddmp, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity dmp data', 'namtsd', no_print )
142         !
143      ENDIF
144      !
145   END SUBROUTINE dta_tsd_init
146
147
148   SUBROUTINE dta_tsd( kt, cddta, ptsd )
149      !!----------------------------------------------------------------------
150      !!                   ***  ROUTINE dta_tsd  ***
[14995]151      !!
[11889]152      !! ** Purpose :   provides T and S data at kt
[14995]153      !!
[11889]154      !! ** Method  : - call fldread routine
[14995]155      !!              - ORCA_R2: add some hand made alteration to read data
[11889]156      !!              - 'key_orca_lev10' interpolates on 10 times more levels
157      !!              - s- or mixed z-s coordinate: vertical interpolation on model mesh
158      !!              - ln_tsd_dmp=F: deallocates the T-S data structure
159      !!                as T-S data are no are used
160      !!
161      !! ** Action  :   ptsd   T-S data on medl mesh and interpolated at time-step kt
162      !!----------------------------------------------------------------------
[14995]163      INTEGER                          , INTENT(in   ) ::   kt     ! ocean time-step
164      CHARACTER(LEN=3)                 , INTENT(in   ) ::   cddta  ! dmp or ini
[13982]165      REAL(wp), DIMENSION(A2D(nn_hls),jpk,jpts), INTENT(  out) ::   ptsd   ! T & S data
[11889]166      !
167      INTEGER ::   ji, jj, jk, jl, jkk   ! dummy loop indicies
168      INTEGER ::   ik, il0, il1, ii0, ii1, ij0, ij1   ! local integers
169      REAL(wp)::   zl, zi                             ! local scalars
170      REAL(wp), DIMENSION(jpk) ::  ztp, zsp   ! 1D workspace
171      !!----------------------------------------------------------------------
172      !
[14857]173      IF( .NOT. l_istiled .OR. ntile == 1 )  THEN                                         ! Do only for the full domain
174         IF( ln_tile ) CALL dom_tile_stop( ldhold=.TRUE. )             ! Use full domain
[13982]175
176         SELECT CASE(cddta)
177         CASE('ini')
178            CALL fld_read( kt, 1, sf_tsdini ) !==   read T & S data at kt time step   ==!
179         CASE('dmp')
180            CALL fld_read( kt, 1, sf_tsddmp ) !==   read T & S data at kt time step   ==!
181         CASE DEFAULT
182            CALL ctl_stop('STOP', 'dta_tsd: cddta case unknown')
183         END SELECT
184
[14857]185         IF( ln_tile ) CALL dom_tile_start( ldhold=.TRUE. )            ! Revert to tile domain
[13982]186      ENDIF
187      !
[11889]188      SELECT CASE(cddta)
[13982]189      CASE('ini')
190         DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
191            ptsd(ji,jj,jk,jp_tem) = sf_tsdini(jp_tem)%fnow(ji,jj,jk)    ! NO mask
192            ptsd(ji,jj,jk,jp_sal) = sf_tsdini(jp_sal)%fnow(ji,jj,jk)
193         END_3D
[11889]194      CASE('dmp')
[13982]195         DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
196            ptsd(ji,jj,jk,jp_tem) = sf_tsddmp(jp_tem)%fnow(ji,jj,jk)    ! NO mask
197            ptsd(ji,jj,jk,jp_sal) = sf_tsddmp(jp_sal)%fnow(ji,jj,jk)
198         END_3D
[11889]199      CASE DEFAULT
200         CALL ctl_stop('STOP', 'dta_tsd: cddta case unknown')
201      END SELECT
202      !
203      IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==!
204         !
[14857]205         IF( .NOT. l_istiled .OR. ntile == 1 )  THEN                       ! Do only on the first tile
[13982]206            IF( kt == nit000 .AND. lwp )THEN
207               WRITE(numout,*)
208               WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh'
209            ENDIF
[11889]210         ENDIF
211         !
[14857]212         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )                  ! vertical interpolation of T & S
[12905]213            DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points
214               zl = gdept_0(ji,jj,jk)
215               IF(     zl < gdept_1d(1  ) ) THEN          ! above the first level of data
216                  ztp(jk) =  ptsd(ji,jj,1    ,jp_tem)
217                  zsp(jk) =  ptsd(ji,jj,1    ,jp_sal)
218               ELSEIF( zl > gdept_1d(jpk) ) THEN          ! below the last level of data
219                  ztp(jk) =  ptsd(ji,jj,jpkm1,jp_tem)
220                  zsp(jk) =  ptsd(ji,jj,jpkm1,jp_sal)
221               ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1
222                  DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1)
223                     IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN
224                        zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk))
[14995]225                        ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi
[12905]226                        zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi
227                     ENDIF
228                  END DO
229               ENDIF
[11889]230            END DO
[12905]231            DO jk = 1, jpkm1
232               ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord
233               ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk)
234            END DO
235            ptsd(ji,jj,jpk,jp_tem) = 0._wp
236            ptsd(ji,jj,jpk,jp_sal) = 0._wp
237         END_2D
[14995]238         !
[11889]239      ELSE                                !==   z- or zps- coordinate   ==!
[14995]240         !
[13982]241         DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
242            ptsd(ji,jj,jk,jp_tem) = ptsd(ji,jj,jk,jp_tem) * tmask(ji,jj,jk)    ! Mask
243            ptsd(ji,jj,jk,jp_sal) = ptsd(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
244         END_3D
[11889]245         !
246         IF( ln_zps ) THEN                      ! zps-coordinate (partial steps) interpolation at the last ocean level
[14857]247            DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
[14995]248               ik = mbkt(ji,jj)
[12905]249               IF( ik > 1 ) THEN
250                  zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
251                  ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik-1,jp_tem)
252                  ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik-1,jp_sal)
253               ENDIF
254               ik = mikt(ji,jj)
255               IF( ik > 1 ) THEN
[14995]256                  zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) )
[12905]257                  ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik+1,jp_tem)
258                  ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik+1,jp_sal)
259               END IF
260            END_2D
[11889]261         ENDIF
262         !
263      ENDIF
264      !
265      SELECT CASE(cddta)
266      CASE('ini') 
267         !                        !==   deallocate T & S structure   ==!
268         !                                              (data used only for initialisation)
269         IF(lwp) WRITE(numout,*) 'dta_tsd: deallocte T & S arrays as they are only use to initialize the run'
270                                        DEALLOCATE( sf_tsdini(jp_tem)%fnow )     ! T arrays in the structure
271         IF( sf_tsdini(jp_tem)%ln_tint )   DEALLOCATE( sf_tsdini(jp_tem)%fdta )
272                                        DEALLOCATE( sf_tsdini(jp_sal)%fnow )     ! S arrays in the structure
273         IF( sf_tsdini(jp_sal)%ln_tint )   DEALLOCATE( sf_tsdini(jp_sal)%fdta )
274                                        DEALLOCATE( sf_tsdini              )     ! the structure itself
275         !
276      END SELECT
277      !
278   END SUBROUTINE dta_tsd
279
280   !!======================================================================
281END MODULE dtatsd
Note: See TracBrowser for help on using the repository browser.