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 branches/UKMO/test_moci_test_suite_namelist_read/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

source: branches/UKMO/test_moci_test_suite_namelist_read/NEMOGCM/NEMO/OPA_SRC/DOM/dtatsd.F90

Last change on this file was 9383, checked in by andmirek, 6 years ago

#2050 fixes and changes

File size: 17.1 KB
Line 
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)
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
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 + suppression of 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 dom_oce         ! ocean space and time domain
19   USE fldread         ! read input fields
20   USE in_out_manager  ! I/O manager
21   USE phycst          ! physical constants
22   USE lib_mpp         ! MPP library
23   USE wrk_nemo        ! Memory allocation
24   USE timing          ! Timing
25   USE iom_def, ONLY:lspr
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC   dta_tsd_init   ! called by opa.F90
31   PUBLIC   dta_tsd        ! called by istate.F90 and tradmp.90
32   PRIVATE  tsd_namelist
33
34   LOGICAL , PUBLIC ::   ln_tsd_init      !: T & S data flag
35   LOGICAL , PUBLIC ::   ln_tsd_tradmp    !: internal damping toward input data flag
36   LOGICAL , PUBLIC ::   ln_tsd_sio       !: read file using 1 processor
37
38   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tsd   ! structure of input SST (file informations, fields read)
39
40   !! * Substitutions
41#  include "domzgr_substitute.h90"
42   !!----------------------------------------------------------------------
43   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
44   !! $Id$
45   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
46   !!----------------------------------------------------------------------
47CONTAINS
48
49   SUBROUTINE dta_tsd_init( ld_tradmp )
50      !!----------------------------------------------------------------------
51      !!                   ***  ROUTINE dta_tsd_init  ***
52      !!                   
53      !! ** Purpose :   initialisation of T & S input data
54      !!
55      !! ** Method  : - Read namtsd namelist
56      !!              - allocates T & S data structure
57      !!----------------------------------------------------------------------
58      LOGICAL, INTENT(in), OPTIONAL ::   ld_tradmp   ! force the initialization when tradp is used
59      !
60      INTEGER ::   ierr0, ierr1, ierr2, ierr3   ! temporary integers
61      !
62      CHARACTER(len=100)            ::   cn_dir          ! Root directory for location of ssr files
63      TYPE(FLD_N), DIMENSION( jpts) ::   slf_i           ! array of namelist informations on the fields to read
64      TYPE(FLD_N)                   ::   sn_tem, sn_sal
65      !!
66      NAMELIST/namtsd/   ln_tsd_init, ln_tsd_tradmp, cn_dir, sn_tem, sn_sal, ln_tsd_sio
67      INTEGER  ::   ios
68      !!----------------------------------------------------------------------
69      !
70      IF( nn_timing == 1 )  CALL timing_start('dta_tsd_init')
71      !
72      !  Initialisation
73      ierr0 = 0  ;  ierr1 = 0  ;  ierr2 = 0  ;  ierr3 = 0
74      !
75      ln_tsd_sio = .FALSE.
76      IF(lwm) THEN
77         REWIND( numnam_ref )              ! Namelist namtsd in reference namelist :
78         READ  ( numnam_ref, namtsd, IOSTAT = ios, ERR = 901)
79901      CONTINUE
80      ENDIF
81      call mpp_bcast(ios)
82      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtsd in reference namelist', lwp )
83      IF(lwm) THEN
84         REWIND( numnam_cfg )              ! Namelist namtsd in configuration namelist : Parameters of the run
85         READ  ( numnam_cfg, namtsd, IOSTAT = ios, ERR = 902 )
86902      CONTINUE
87      ENDIF
88      call mpp_bcast(ios)
89      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtsd in configuration namelist', lwp )
90
91      IF(lwm) WRITE ( numond, namtsd )
92
93      call tsd_namelist(cn_dir, sn_tem, sn_sal)
94
95      IF( PRESENT( ld_tradmp ) )   ln_tsd_tradmp = .TRUE.     ! forces the initialization when tradmp is used
96     
97      IF(lwp) THEN                  ! control print
98         WRITE(numout,*)
99         WRITE(numout,*) 'dta_tsd_init : Temperature & Salinity data '
100         WRITE(numout,*) '~~~~~~~~~~~~ '
101         WRITE(numout,*) '   Namelist namtsd'
102         WRITE(numout,*) '      Initialisation of ocean T & S with T &S input data   ln_tsd_init   = ', ln_tsd_init
103         WRITE(numout,*) '      damping of ocean T & S toward T &S input data        ln_tsd_tradmp = ', ln_tsd_tradmp
104         WRITE(numout,*)
105         IF( .NOT.ln_tsd_init .AND. .NOT.ln_tsd_tradmp ) THEN
106            WRITE(numout,*)
107            WRITE(numout,*) '   T & S data not used'
108         ENDIF
109      ENDIF
110      !
111      IF( ln_rstart .AND. ln_tsd_init ) THEN
112         CALL ctl_warn( 'dta_tsd_init: ocean restart and T & S data intialisation, ',   &
113            &           'we keep the restart T & S values and set ln_tsd_init to FALSE' )
114         ln_tsd_init = .FALSE.
115      ENDIF
116      !
117      !                             ! allocate the arrays (if necessary)
118      IF(  ln_tsd_init .OR. ln_tsd_tradmp  ) THEN
119         !
120         ALLOCATE( sf_tsd(jpts), STAT=ierr0 )
121         IF( ierr0 > 0 ) THEN
122            CALL ctl_stop( 'dta_tsd_init: unable to allocate sf_tsd structure' )   ;   RETURN
123         ENDIF
124         !
125                                ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk)   , STAT=ierr0 )
126         IF( sn_tem%ln_tint )   ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 )
127                                ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk)   , STAT=ierr2 )
128         IF( sn_sal%ln_tint )   ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 )
129         !
130         IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN
131            CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' )   ;   RETURN
132         ENDIF
133         !                         ! fill sf_tsd with sn_tem & sn_sal and control print
134         slf_i(jp_tem) = sn_tem   ;   slf_i(jp_sal) = sn_sal
135         CALL fld_fill( sf_tsd, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity data', 'namtsd' )
136         !
137      ENDIF
138      !
139      IF( nn_timing == 1 )  CALL timing_stop('dta_tsd_init')
140      !
141   END SUBROUTINE dta_tsd_init
142
143
144   SUBROUTINE dta_tsd( kt, ptsd )
145      !!----------------------------------------------------------------------
146      !!                   ***  ROUTINE dta_tsd  ***
147      !!                   
148      !! ** Purpose :   provides T and S data at kt
149      !!
150      !! ** Method  : - call fldread routine
151      !!              - ORCA_R2: add some hand made alteration to read data 
152      !!              - 'key_orca_lev10' interpolates on 10 times more levels
153      !!              - s- or mixed z-s coordinate: vertical interpolation on model mesh
154      !!              - ln_tsd_tradmp=F: deallocates the T-S data structure
155      !!                as T-S data are no are used
156      !!
157      !! ** Action  :   ptsd   T-S data on medl mesh and interpolated at time-step kt
158      !!----------------------------------------------------------------------
159      INTEGER                              , INTENT(in   ) ::   kt     ! ocean time-step
160      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(  out) ::   ptsd   ! T & S data
161      !
162      INTEGER ::   ji, jj, jk, jl, jkk   ! dummy loop indicies
163      INTEGER ::   ik, il0, il1, ii0, ii1, ij0, ij1   ! local integers
164      REAL(wp)::   zl, zi
165      REAL(wp), POINTER, DIMENSION(:) ::  ztp, zsp   ! 1D workspace
166      !!----------------------------------------------------------------------
167      !
168      IF( nn_timing == 1 )  CALL timing_start('dta_tsd')
169      !
170      lspr = ln_tsd_sio
171      CALL fld_read( kt, 1, sf_tsd )      !==   read T & S data at kt time step   ==!
172      lspr = .false.
173      !
174      !
175      !                                   !==   ORCA_R2 configuration and T & S damping   ==!
176      IF( cp_cfg == "orca" .AND. jp_cfg == 2 .AND. ln_tsd_tradmp ) THEN    ! some hand made alterations
177         !
178         ij0 = 101   ;   ij1 = 109                       ! Reduced T & S in the Alboran Sea
179         ii0 = 141   ;   ii1 = 155
180         DO jj = mj0(ij0), mj1(ij1)
181            DO ji = mi0(ii0), mi1(ii1)
182               sf_tsd(jp_tem)%fnow(ji,jj,13:13) = sf_tsd(jp_tem)%fnow(ji,jj,13:13) - 0.20_wp
183               sf_tsd(jp_tem)%fnow(ji,jj,14:15) = sf_tsd(jp_tem)%fnow(ji,jj,14:15) - 0.35_wp
184               sf_tsd(jp_tem)%fnow(ji,jj,16:25) = sf_tsd(jp_tem)%fnow(ji,jj,16:25) - 0.40_wp
185               !
186               sf_tsd(jp_sal)%fnow(ji,jj,13:13) = sf_tsd(jp_sal)%fnow(ji,jj,13:13) - 0.15_wp
187               sf_tsd(jp_sal)%fnow(ji,jj,14:15) = sf_tsd(jp_sal)%fnow(ji,jj,14:15) - 0.25_wp
188               sf_tsd(jp_sal)%fnow(ji,jj,16:17) = sf_tsd(jp_sal)%fnow(ji,jj,16:17) - 0.30_wp
189               sf_tsd(jp_sal)%fnow(ji,jj,18:25) = sf_tsd(jp_sal)%fnow(ji,jj,18:25) - 0.35_wp
190            END DO
191         END DO
192         IF( nn_cla == 1 ) THEN                          ! Cross Land advection
193            il0 = 138   ;   il1 = 138                          ! set T & S profile at Gibraltar Strait
194            ij0 = 101   ;   ij1 = 102
195            ii0 = 139   ;   ii1 = 139
196            DO jl = mi0(il0), mi1(il1)
197               DO jj = mj0(ij0), mj1(ij1)
198                  DO ji = mi0(ii0), mi1(ii1)
199                     sf_tsd(jp_tem)%fnow(ji,jj,:) = sf_tsd(jp_tem)%fnow(jl,jj,:)
200                     sf_tsd(jp_sal)%fnow(ji,jj,:) = sf_tsd(jp_sal)%fnow(jl,jj,:)
201                  END DO
202               END DO
203            END DO
204            il0 = 164   ;   il1 = 164                          ! set T & S profile at Bab el Mandeb Strait
205            ij0 =  87   ;   ij1 =  88
206            ii0 = 161   ;   ii1 = 163
207            DO jl = mi0(il0), mi1(il1)
208               DO jj = mj0(ij0), mj1(ij1)
209                  DO ji = mi0(ii0), mi1(ii1)
210                     sf_tsd(jp_tem)%fnow(ji,jj,:) = sf_tsd(jp_tem)%fnow(jl,jj,:)
211                     sf_tsd(jp_sal)%fnow(ji,jj,:) = sf_tsd(jp_sal)%fnow(jl,jj,:)
212                  END DO
213               END DO
214            END DO
215         ELSE                                            ! No Cross Land advection
216            ij0 =  87   ;   ij1 =  96                          ! Reduced temperature in Red Sea
217            ii0 = 148   ;   ii1 = 160
218            sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ,  4:10 ) = 7.0_wp
219            sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 11:13 ) = 6.5_wp
220            sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 14:20 ) = 6.0_wp
221         ENDIF
222      ENDIF
223      !
224      ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:)    ! NO mask
225      ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) 
226      !
227      IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==!
228         !
229         CALL wrk_alloc( jpk, ztp, zsp )
230         !
231         IF( kt == nit000 .AND. lwp )THEN
232            WRITE(numout,*)
233            WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh'
234         ENDIF
235         !
236         DO jj = 1, jpj                         ! vertical interpolation of T & S
237            DO ji = 1, jpi
238               DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points
239                  zl = gdept_0(ji,jj,jk)
240                  IF(     zl < gdept_1d(1  ) ) THEN          ! above the first level of data
241                     ztp(jk) =  ptsd(ji,jj,1    ,jp_tem)
242                     zsp(jk) =  ptsd(ji,jj,1    ,jp_sal)
243                  ELSEIF( zl > gdept_1d(jpk) ) THEN          ! below the last level of data
244                     ztp(jk) =  ptsd(ji,jj,jpkm1,jp_tem)
245                     zsp(jk) =  ptsd(ji,jj,jpkm1,jp_sal)
246                  ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1
247                     DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1)
248                        IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN
249                           zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk))
250                           ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi 
251                           zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi
252                        ENDIF
253                     END DO
254                  ENDIF
255               END DO
256               DO jk = 1, jpkm1
257                  ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord
258                  ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk)
259               END DO
260               ptsd(ji,jj,jpk,jp_tem) = 0._wp
261               ptsd(ji,jj,jpk,jp_sal) = 0._wp
262            END DO
263         END DO
264         !
265         CALL wrk_dealloc( jpk, ztp, zsp )
266         !
267      ELSE                                !==   z- or zps- coordinate   ==!
268         !                             
269         ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:)    ! Mask
270         ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) * tmask(:,:,:)
271         !
272         IF( ln_zps ) THEN                      ! zps-coordinate (partial steps) interpolation at the last ocean level
273            DO jj = 1, jpj
274               DO ji = 1, jpi
275                  ik = mbkt(ji,jj) 
276                  IF( ik > 1 ) THEN
277                     zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
278                     ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik-1,jp_tem)
279                     ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik-1,jp_sal)
280                  ENDIF
281                  ik = mikt(ji,jj)
282                  IF( ik > 1 ) THEN
283                     zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 
284                     ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik+1,jp_tem)
285                     ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik+1,jp_sal)
286                  END IF
287               END DO
288            END DO
289         ENDIF
290         !
291      ENDIF
292      !
293      IF( lwp .AND. kt == nit000 ) THEN
294         WRITE(numout,*) ' temperature Levitus '
295         WRITE(numout,*)
296         WRITE(numout,*)'  level = 1'
297         CALL prihre( ptsd(:,:,1    ,jp_tem), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
298         WRITE(numout,*)'  level = ', jpk/2
299         CALL prihre( ptsd(:,:,jpk/2,jp_tem), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
300         WRITE(numout,*)'  level = ', jpkm1
301         CALL prihre( ptsd(:,:,jpkm1,jp_tem), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
302         WRITE(numout,*)
303         WRITE(numout,*) ' salinity Levitus '
304         WRITE(numout,*)
305         WRITE(numout,*)'  level = 1'
306         CALL prihre( ptsd(:,:,1    ,jp_sal), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
307         WRITE(numout,*)'  level = ', jpk/2
308         CALL prihre( ptsd(:,:,jpk/2,jp_sal), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
309         WRITE(numout,*)'  level = ', jpkm1
310         CALL prihre( ptsd(:,:,jpkm1,jp_sal), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
311         WRITE(numout,*)
312      ENDIF
313      !
314      IF( .NOT.ln_tsd_tradmp ) THEN                   !==   deallocate T & S structure   ==!
315         !                                              (data used only for initialisation)
316         IF(lwp) WRITE(numout,*) 'dta_tsd: deallocte T & S arrays as they are only use to initialize the run'
317                                        DEALLOCATE( sf_tsd(jp_tem)%fnow )     ! T arrays in the structure
318         IF( sf_tsd(jp_tem)%ln_tint )   DEALLOCATE( sf_tsd(jp_tem)%fdta )
319                                        DEALLOCATE( sf_tsd(jp_sal)%fnow )     ! S arrays in the structure
320         IF( sf_tsd(jp_sal)%ln_tint )   DEALLOCATE( sf_tsd(jp_sal)%fdta )
321                                        DEALLOCATE( sf_tsd              )     ! the structure itself
322      ENDIF
323      !
324      IF( nn_timing == 1 )  CALL timing_stop('dta_tsd')
325      !
326   END SUBROUTINE dta_tsd
327
328   SUBROUTINE tsd_namelist(cd_dir, sd_tem, sd_sal)
329     !!---------------------------------------------------------------------
330     !!                   ***  ROUTINE tsd_namelist  ***
331     !!                     
332     !! ** Purpose :   Broadcast namelist variables read by procesor lwm
333     !!
334     !! ** Method  :   use lib_mpp
335     !!----------------------------------------------------------------------
336     CHARACTER(len=100)           ::   cd_dir          ! Root directory for location of ssr files
337     TYPE(FLD_N)                  ::   sd_tem, sd_sal
338#if defined key_mpp_mpi
339      CALL mpp_bcast(ln_tsd_init)
340      CALL mpp_bcast(ln_tsd_tradmp)
341      CALL mpp_bcast(cd_dir, 100)
342      CALL fld_n_bcast(sd_tem)
343      CALL fld_n_bcast(sd_sal)
344      CALL mpp_bcast(ln_tsd_sio)
345#endif
346   END SUBROUTINE tsd_namelist
347
348   !!======================================================================
349END MODULE dtatsd
Note: See TracBrowser for help on using the repository browser.