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.
dtadyn.F90 in trunk/NEMO/OFF_SRC – NEMO

source: trunk/NEMO/OFF_SRC/dtadyn.F90 @ 1699

Last change on this file since 1699 was 1699, checked in by cetlod, 14 years ago

Correction in estimation of wind module with masked wind stress, see changeset:1643

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 40.3 KB
Line 
1MODULE dtadyn
2   !!======================================================================
3   !!                       ***  MODULE  dtadyn  ***
4   !! OFFLINE : interpolation of the physical fields
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dta_dyn_init : initialization, namelist read, and parameters control
9   !!   dta_dyn      : Interpolation of the fields
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and tracers variables
13   USE dom_oce         ! ocean space and time domain variables
14   USE zdf_oce         ! ocean vertical physics
15   USE in_out_manager  ! I/O manager
16   USE phycst          ! physical constants
17   USE sbc_oce
18   USE ldfslp
19   USE ldfeiv          ! eddy induced velocity coef.      (ldf_eiv routine)
20   USE ldftra_oce      ! ocean tracer   lateral physics
21   USE zdfmxl
22   USE trabbl
23   USE eosbn2
24   USE zdfddm          ! vertical  physics: double diffusion
25   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
26   USE zpshde
27   USE lib_mpp         ! distributed memory computing library
28
29   IMPLICIT NONE
30   PRIVATE
31
32   !! *  Routine accessibility
33   PUBLIC dta_dyn_init   ! called by opa.F90
34   PUBLIC dta_dyn        ! called by step.F90
35
36   LOGICAL , PUBLIC :: &
37      lperdyn = .TRUE. , & ! boolean for periodic fields or not
38      lfirdyn = .TRUE.     ! boolean for the first call or not
39
40   INTEGER , PUBLIC :: &
41      ndtadyn = 73 ,  & ! Number of dat in one year
42      ndtatot = 73 ,  & ! Number of data in the input field
43      nsptint = 1 ,   & ! type of spatial interpolation
44      nficdyn = 2       ! number of dynamical fields
45
46   INTEGER ::     &
47      ndyn1, ndyn2 , &
48      nlecoff = 0  , & ! switch for the first read
49      numfl_t, numfl_u, &
50      numfl_v, numfl_w
51
52   CHARACTER(len=45)  ::  &
53      cfile_grid_T = 'dyna_grid_T.nc', &   !: name of the grid_T file
54      cfile_grid_U = 'dyna_grid_U.nc', &   !: name of the grid_U file
55      cfile_grid_V = 'dyna_grid_V.nc', &   !: name of the grid_V file
56      cfile_grid_W = 'dyna_grid_W.nc'      !: name of the grid_W file
57
58   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
59      tdta   ,   & ! temperature at two consecutive times
60      sdta   ,   & ! salinity at two consecutive times
61      udta   ,   & ! zonal velocity at two consecutive times
62      vdta   ,   & ! meridional velocity at two consecutive times
63      wdta   ,   & ! vertical velocity at two consecutive times
64#if defined key_trc_diatrd
65      hdivdta,   & ! horizontal divergence
66#endif
67      avtdta       ! vertical diffusivity coefficient
68
69   REAL(wp), DIMENSION(jpi,jpj,2) ::       &
70      hmlddta,   & ! mixed layer depth at two consecutive times
71      wspddta,   & ! wind speed at two consecutive times
72      frlddta,   & ! sea-ice fraction at two consecutive times
73      empdta ,   & ! E-P at two consecutive times
74      qsrdta       ! short wave heat flux at two consecutive times
75
76#if defined key_ldfslp
77   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
78      uslpdta ,  & ! zonal isopycnal slopes
79      vslpdta ,  & ! meridional isopycnal slopes
80      wslpidta , & ! zonal diapycnal slopes
81      wslpjdta     ! meridional diapycnal slopes
82#endif
83
84#if ! defined key_off_degrad &&  defined key_traldf_c2d
85   REAL(wp), DIMENSION(jpi,jpj,2) ::   &
86      ahtwdta      ! Lateral diffusivity
87# if defined key_trcldf_eiv 
88   REAL(wp), DIMENSION(jpi,jpj,2) ::   &
89      aeiwdta      ! G&M coefficient
90# endif
91#endif
92
93#if defined key_off_degrad
94   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
95      ahtudta, ahtvdta, ahtwdta  !  Lateral diffusivity
96# if defined key_trcldf_eiv
97   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
98      aeiudta, aeivdta, aeiwdta  ! G&M coefficient
99# endif
100
101#endif
102
103#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
104   REAL(wp), DIMENSION(jpi,jpj,2) ::       &
105      bblxdta ,  & ! frequency of bbl in the x direction at 2 consecutive times
106      bblydta      ! frequency of bbl in the y direction at 2 consecutive times
107#endif
108
109   !! * Substitutions
110#  include "domzgr_substitute.h90"
111#  include "vectopt_loop_substitute.h90"
112   !!----------------------------------------------------------------------
113   !!   OPA 9.0 , LOCEAN-IPSL  (2005)
114   !!   $Id$
115   !!   This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
116   !!----------------------------------------------------------------------
117
118CONTAINS
119
120   SUBROUTINE dta_dyn( kt )
121      !!----------------------------------------------------------------------
122      !!                  ***  ROUTINE dta_dyn  ***
123      !!
124      !! ** Purpose : Prepares dynamics and physics fields from an
125      !!              OPA9 simulation  for an off-line simulation
126      !!               for passive tracer
127      !!
128      !! ** Method : calculates the position of DATA to read READ DATA
129      !!             (example month changement) computes slopes IF needed
130      !!             interpolates DATA IF needed
131      !!
132      !! ** History :
133      !!   ! original  : 92-01 (M. Imbard: sub domain)
134      !!   ! addition  : 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
135      !!   ! addition  : 98-05 (L. Bopp read output of coupled run)
136      !!   ! addition  : 05-03 (O. Aumont and A. El Moussaoui) F90
137      !!   ! addition  : 05-12 (C. Ethe) Adapted for DEGINT
138      !!----------------------------------------------------------------------
139      !! * Arguments
140      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
141
142      !! * Local declarations
143      INTEGER ::   iper, iperm1, iswap, izt   
144
145      REAL(wp) :: zpdtan, zpdtpe, zdemi, zt
146      REAL(wp) :: zweigh
147
148      ! 0. Initialization
149      ! -----------------
150
151      IF( lfirdyn ) THEN
152         ! first time step MUST BE nit000
153         IF( kt /= nit000 ) THEN
154            IF (lwp) THEN
155               WRITE (numout,*) ' kt MUST BE EQUAL to nit000. kt = ',kt ,' nit000 = ',nit000 
156              STOP 'dtadyn'
157            ENDIF
158          ENDIF 
159          ! Initialize the parameters of the interpolation
160          CALL dta_dyn_init
161      ENDIF
162
163
164      zpdtan = raass / rdt
165      zpdtpe = zpdtan / FLOAT( ndtadyn )
166      zdemi  = zpdtpe * 0.5
167      zt     = ( FLOAT (kt) + zdemi ) / zpdtpe
168   
169      izt      = INT( zt )
170      zweigh   = zt - FLOAT( INT(zt) )
171
172      IF( lperdyn ) THEN
173         iperm1 = MOD( izt, ndtadyn )
174      ELSE
175         iperm1 = MOD( izt, ndtatot - 1 ) + 1
176      ENDIF
177
178      iper = iperm1 + 1
179      IF( iperm1 == 0 ) THEN
180          IF( lperdyn ) THEN
181              iperm1 = ndtadyn
182          ELSE
183              IF( lfirdyn ) THEN
184                  IF (lwp) WRITE (numout,*) & 
185                      &   ' dynamic file is not periodic with or without interpolation  &
186                      &   we take the first value for the previous period iperm1 = 0  '
187              END IF
188          END IF
189      END IF
190
191      iswap  = 0
192
193      ! 1. First call lfirdyn = true
194      ! ----------------------------
195
196      IF( lfirdyn ) THEN
197         ! store the information of the period read
198         ndyn1 = iperm1
199         ndyn2 = iper
200         
201         IF (lwp) THEN
202            WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
203               &             ' and for the period ndyn2 = ',ndyn2
204            WRITE (numout,*) ' time step is : ', kt
205            WRITE (numout,*) ' we have ndtadyn = ',ndtadyn,' records in the dynamic file for one year'
206         END IF
207         !
208         IF( iperm1 /= 0 ) THEN         ! data read for the iperm1 period
209            CALL dynrea( kt, iperm1 ) 
210         ELSE
211            CALL dynrea( kt, 1 )
212         ENDIF
213         
214#if defined key_ldfslp
215         ! Computes slopes
216         ! Caution : here tn, sn and avt are used as workspace
217         tn (:,:,:) = tdta  (:,:,:,2)
218         sn (:,:,:) = sdta  (:,:,:,2)
219         avt(:,:,:) = avtdta(:,:,:,2)
220         
221         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
222         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
223         IF( ln_zps )   &
224            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
225            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
226            &                 gtv, gsv, grv )
227         CALL zdf_mxl( kt )              ! mixed layer depth
228         CALL ldf_slp( kt, rhd, rn2 )
229         
230         uslpdta (:,:,:,2) = uslp (:,:,:)
231         vslpdta (:,:,:,2) = vslp (:,:,:)
232         wslpidta(:,:,:,2) = wslpi(:,:,:)
233         wslpjdta(:,:,:,2) = wslpj(:,:,:)
234#endif
235         
236         ! swap from record 2 to 1
237         CALL swap_dyn_data
238         
239         iswap = 1        !  indicates swap
240         
241         CALL dynrea( kt, iper )    ! data read for the iper period
242         
243#if defined key_ldfslp
244         ! Computes slopes
245         ! Caution : here tn, sn and avt are used as workspace
246         tn (:,:,:) = tdta  (:,:,:,2)
247         sn (:,:,:) = sdta  (:,:,:,2)
248         avt(:,:,:) = avtdta(:,:,:,2)
249         
250         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
251         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
252         IF( ln_zps )   &
253            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
254            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
255            &                 gtv, gsv, grv )
256         CALL zdf_mxl( kt )              ! mixed layer depth
257         CALL ldf_slp( kt, rhd, rn2 )
258         
259         uslpdta (:,:,:,2) = uslp (:,:,:)
260         vslpdta (:,:,:,2) = vslp (:,:,:)
261         wslpidta(:,:,:,2) = wslpi(:,:,:)
262         wslpjdta(:,:,:,2) = wslpj(:,:,:)
263#endif
264         !
265         lfirdyn=.FALSE.    ! trace the first call
266      ENDIF
267      !
268      ! And now what we have to do at every time step
269      ! check the validity of the period in memory
270      !
271      IF( iperm1 /= ndyn1 ) THEN
272
273         IF( iperm1 == 0. ) THEN
274            IF (lwp) THEN
275               WRITE (numout,*) ' dynamic file is not periodic with periodic interpolation'
276               WRITE (numout,*) ' we take the last value for the last period '
277               WRITE (numout,*) ' iperm1 = 12,   iper = 13  '
278            ENDIF
279            iperm1 = 12
280            iper   = 13
281         ENDIF
282         !
283         ! We have to prepare a new read of data : swap from record 2 to 1
284         !
285         CALL swap_dyn_data
286
287         iswap = 1        !  indicates swap
288         
289         CALL dynrea( kt, iper )    ! data read for the iper period
290
291#if defined key_ldfslp
292         ! Computes slopes
293         ! Caution : here tn, sn and avt are used as workspace
294         tn (:,:,:) = tdta  (:,:,:,2)
295         sn (:,:,:) = sdta  (:,:,:,2)
296         avt(:,:,:) = avtdta(:,:,:,2)
297         
298         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
299         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
300         IF( ln_zps )   &
301            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
302            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
303            &                 gtv, gsv, grv )
304         CALL zdf_mxl( kt )              ! mixed layer depth
305         CALL ldf_slp( kt, rhd, rn2 )
306         
307         uslpdta (:,:,:,2) = uslp (:,:,:)
308         vslpdta (:,:,:,2) = vslp (:,:,:)
309         wslpidta(:,:,:,2) = wslpi(:,:,:)
310         wslpjdta(:,:,:,2) = wslpj(:,:,:)
311#endif
312       
313         ! store the information of the period read
314         ndyn1 = ndyn2
315         ndyn2 = iper
316
317         IF (lwp) THEN
318            WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
319               &             ' and for the period ndyn2 = ',ndyn2
320            WRITE (numout,*) ' time step is : ', kt
321         END IF
322         !
323      END IF
324      !
325      ! Compute the data at the given time step
326      !----------------------------------------     
327
328      IF( nsptint == 0 ) THEN   
329         ! No spatial interpolation, data are probably correct
330         ! We have to initialize data if we have changed the period         
331         CALL assign_dyn_data         
332      ELSE IF( nsptint == 1 ) THEN
333         ! linear interpolation
334         CALL linear_interp_dyn_data( zweigh )
335      ELSE 
336         ! other interpolation
337         WRITE (numout,*) ' this kind of interpolation do not exist at the moment : we stop'
338         STOP 'dtadyn'         
339      END IF
340     
341      ! In any case, we need rhop
342      CALL eos( tn, sn, rhd, rhop ) 
343     
344#if ! defined key_off_degrad && defined key_traldf_c2d
345      ! In case of 2D varying coefficients, we need aeiv and aeiu
346      IF( lk_traldf_eiv )   CALL ldf_eiv( kt )      ! eddy induced velocity coefficient
347#endif
348   
349   END SUBROUTINE dta_dyn
350
351   SUBROUTINE dynrea( kt, kenr )
352      !!----------------------------------------------------------------------
353      !!                  ***  ROUTINE dynrea  ***
354      !!
355      !! ** Purpose : READ dynamics fiels from OPA9 netcdf output
356      !!
357      !! ** Method : READ the kenr records of DATA and store in
358      !!             in udta(...,2), .... 
359      !!
360      !! ** History : additions : M. Levy et M. Benjelloul jan 2001
361      !!              (netcdf FORMAT)
362      !!              05-03 (O. Aumont and A. El Moussaoui) F90
363      !!              06-07 : (C. Ethe) use of iom module
364      !!----------------------------------------------------------------------
365      !! * Modules used
366      USE iom
367
368      !! * Arguments
369      INTEGER, INTENT( in ) ::   kt, kenr       ! time index
370      !! * Local declarations
371      INTEGER ::  jkenr
372
373      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
374        zu, zv, zw, zt, zs, zavt ,   &     ! 3-D dynamical fields
375        zhdiv                              ! horizontal divergence
376
377      REAL(wp), DIMENSION(jpi,jpj) :: &
378         zemp, zqsr, zmld, zice, zwspd, &
379         ztaux, ztauy
380#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
381      REAL(wp), DIMENSION(jpi,jpj) :: zbblx, zbbly
382#endif
383
384#if ! defined key_off_degrad && defined key_traldf_c2d
385      REAL(wp), DIMENSION(jpi,jpj) :: zahtw 
386#   if defined key_trcldf_eiv
387      REAL(wp), DIMENSION(jpi,jpj) :: zaeiw 
388#  endif
389#endif
390
391#if defined key_off_degrad
392   REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
393      zahtu, zahtv, zahtw  !  Lateral diffusivity
394# if defined key_trcldf_eiv
395   REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
396      zaeiu, zaeiv, zaeiw  ! G&M coefficient
397# endif
398#endif
399
400      !---------------------------------------------------------------
401      ! 0. Initialization
402     
403      ! cas d'un fichier non periodique : on utilise deux fois le premier et
404      ! le dernier champ temporel
405
406      jkenr = kenr
407
408      IF(lwp) THEN
409         WRITE(numout,*)
410         WRITE(numout,*) 'Dynrea : reading dynamical fields, kenr = ', jkenr
411         WRITE(numout,*) ' ~~~~~~~'
412#if defined key_off_degrad
413         WRITE(numout,*) ' Degraded fields'
414#endif
415         WRITE(numout,*)
416      ENDIF
417
418
419      IF( kt == nit000 .AND. nlecoff == 0 ) THEN
420         nlecoff = 1
421         CALL  iom_open ( cfile_grid_T, numfl_t )
422         CALL  iom_open ( cfile_grid_U, numfl_u )
423         CALL  iom_open ( cfile_grid_V, numfl_v )
424         CALL  iom_open ( cfile_grid_W, numfl_w )
425      ENDIF
426
427      ! file grid-T
428      !---------------
429      CALL iom_get( numfl_t, jpdom_data, 'votemper', zt   (:,:,:), jkenr )
430      CALL iom_get( numfl_t, jpdom_data, 'vosaline', zs   (:,:,:), jkenr )
431      CALL iom_get( numfl_t, jpdom_data, 'somixhgt', zmld (:,:  ), jkenr )
432      CALL iom_get( numfl_t, jpdom_data, 'sowaflcd', zemp (:,:  ), jkenr )
433      CALL iom_get( numfl_t, jpdom_data, 'soshfldo', zqsr (:,:  ), jkenr )
434      CALL iom_get( numfl_t, jpdom_data, 'soicecov', zice (:,:  ), jkenr )
435      IF( iom_varid( numfl_t, 'sowindsp', ldstop = .FALSE. ) > 0 ) THEN
436         CALL iom_get( numfl_t, jpdom_data, 'sowindsp', zwspd(:,:), jkenr ) 
437      ELSE
438         CALL iom_get( numfl_u, jpdom_data, 'sozotaux', ztaux(:,:), jkenr )
439         CALL iom_get( numfl_v, jpdom_data, 'sometauy', ztauy(:,:), jkenr )
440         CALL tau2wnd( ztaux, ztauy, zwspd )
441      ENDIF
442      ! files grid-U / grid_V
443      CALL iom_get( numfl_u, jpdom_data, 'vozocrtx', zu   (:,:,:), jkenr )
444      CALL iom_get( numfl_v, jpdom_data, 'vomecrty', zv   (:,:,:), jkenr )
445
446#if defined key_trcbbl_dif || defined key_trcbbl_adv
447      IF( iom_varid( numfl_u, 'sobblcox', ldstop = .FALSE. ) > 0  .AND. &
448      &   iom_varid( numfl_v, 'sobblcoy', ldstop = .FALSE. ) > 0 ) THEN
449         CALL iom_get( numfl_u, jpdom_data, 'sobblcox', zbblx(:,:), jkenr )
450         CALL iom_get( numfl_v, jpdom_data, 'sobblcoy', zbbly(:,:), jkenr )
451      ELSE
452         CALL bbl_sign( zt, zs, zbblx, zbbly )   
453      ENDIF
454#endif
455
456      ! file grid-W
457!!      CALL iom_get ( numfl_w, jpdom_data, 'vovecrtz', zw   (:,:,:), jkenr )
458      ! Computation of vertical velocity using horizontal divergence
459      CALL wzv( zu, zv, zw, zhdiv )
460
461# if defined key_zdfddm
462      CALL iom_get( numfl_w, jpdom_data, 'voddmavs', zavt (:,:,:), jkenr )
463#else
464      CALL iom_get( numfl_w, jpdom_data, 'votkeavt', zavt (:,:,:), jkenr )
465#endif 
466
467#if ! defined key_off_degrad && defined key_traldf_c2d
468      CALL iom_get( numfl_w, jpdom_data, 'soleahtw', zahtw (:,: ), jkenr )
469#  if   defined key_trcldf_eiv 
470      CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw (:,: ), jkenr )
471#  endif
472#endif
473
474#if defined key_off_degrad
475      CALL iom_get( numfl_u, jpdom_data, 'vozoahtu', zahtu(:,:,:), jkenr )
476      CALL iom_get( numfl_v, jpdom_data, 'vomeahtv', zahtv(:,:,:), jkenr )
477      CALL iom_get( numfl_w, jpdom_data, 'voveahtw', zahtw(:,:,:), jkenr )
478#  if defined key_trcldf_eiv
479      CALL iom_get( numfl_u, jpdom_data, 'vozoaeiu', zaeiu(:,:,:), jkenr )
480      CALL iom_get( numfl_v, jpdom_data, 'vomeaeiv', zaeiv(:,:,:), jkenr )
481      CALL iom_get( numfl_w, jpdom_data, 'voveaeiw', zaeiw(:,:,:), jkenr )
482#  endif
483#endif
484
485      udta(:,:,:,2) = zu(:,:,:) * umask(:,:,:)
486      vdta(:,:,:,2) = zv(:,:,:) * vmask(:,:,:) 
487      wdta(:,:,:,2) = zw(:,:,:) * tmask(:,:,:)
488
489#if defined key_trc_diatrd
490      hdivdta(:,:,:,2) = zhdiv(:,:,:) * tmask(:,:,:)
491#endif
492
493      tdta(:,:,:,2)   = zt  (:,:,:) * tmask(:,:,:)
494      sdta(:,:,:,2)   = zs  (:,:,:) * tmask(:,:,:)
495      avtdta(:,:,:,2) = zavt(:,:,:) * tmask(:,:,:)
496
497#if ! defined key_off_degrad && defined key_traldf_c2d
498      ahtwdta(:,:,2)  = zahtw(:,:) * tmask(:,:,1)
499#if defined key_trcldf_eiv
500      aeiwdta(:,:,2)  = zaeiw(:,:) * tmask(:,:,1)
501#endif
502#endif
503
504#if defined key_off_degrad
505        ahtudta(:,:,:,2) = zahtu(:,:,:) * umask(:,:,:)
506        ahtvdta(:,:,:,2) = zahtv(:,:,:) * vmask(:,:,:)
507        ahtwdta(:,:,:,2) = zahtw(:,:,:) * tmask(:,:,:)
508#  if defined key_trcldf_eiv
509        aeiudta(:,:,:,2) = zaeiu(:,:,:) * umask(:,:,:)
510        aeivdta(:,:,:,2) = zaeiv(:,:,:) * vmask(:,:,:)
511        aeiwdta(:,:,:,2) = zaeiw(:,:,:) * tmask(:,:,:)
512#  endif
513#endif
514
515      ! fluxes
516      !
517      wspddta(:,:,2)  = zwspd(:,:) * tmask(:,:,1)
518      frlddta(:,:,2)  = MIN( 1., zice(:,:) ) * tmask(:,:,1)
519      empdta (:,:,2)  = zemp(:,:) * tmask(:,:,1)
520      qsrdta (:,:,2)  = zqsr(:,:) * tmask(:,:,1)
521      hmlddta(:,:,2)  = zmld(:,:) * tmask(:,:,1)
522     
523#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
524      bblxdta(:,:,2) = MAX( 0., zbblx(:,:) )
525      bblydta(:,:,2) = MAX( 0., zbbly(:,:) )
526
527      WHERE( bblxdta(:,:,2) > 2. ) bblxdta(:,:,2) = 0.
528      WHERE( bblydta(:,:,2) > 2. ) bblydta(:,:,2) = 0.
529#endif
530
531      IF( kt == nitend ) THEN
532         CALL iom_close ( numfl_t )
533         CALL iom_close ( numfl_u )
534         CALL iom_close ( numfl_v )
535         CALL iom_close ( numfl_w )
536      ENDIF
537     
538   END SUBROUTINE dynrea
539
540   SUBROUTINE dta_dyn_init
541      !!----------------------------------------------------------------------
542      !!                  ***  ROUTINE dta_dyn_init  ***
543      !!
544      !! ** Purpose :   initializations of parameters for the interpolation
545      !!
546      !! ** Method :
547      !!
548      !! History :
549      !!    ! original  : 92-01 (M. Imbard: sub domain)
550      !!    ! 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
551      !!    ! 98-05 (L. Bopp read output of coupled run)
552      !!    ! 05-03 (O. Aumont and A. El Moussaoui) F90
553      !!----------------------------------------------------------------------
554      !! * Modules used
555
556      !! * Local declarations
557
558
559      NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, nficdyn, lperdyn,  &
560      &                cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W
561      !!----------------------------------------------------------------------
562
563      !  Define the dynamical input parameters
564      ! ======================================
565
566      ! Read Namelist namdyn : Lateral physics on tracers
567      REWIND( numnam )
568      READ  ( numnam, namdyn )
569
570      IF(lwp) THEN
571         WRITE(numout,*)
572         WRITE(numout,*) 'namdyn : offline dynamical selection'
573         WRITE(numout,*) '~~~~~~~'
574         WRITE(numout,*) '  Namelist namdyn : set parameters for the lecture of the dynamical fields'
575         WRITE(numout,*) 
576         WRITE(numout,*) ' number of elements in the FILE for a year  ndtadyn = ' , ndtadyn
577         WRITE(numout,*) ' total number of elements in the FILE       ndtatot = ' , ndtatot
578         WRITE(numout,*) ' type of interpolation                      nsptint = ' , nsptint
579         WRITE(numout,*) ' number of dynamics FILE                    nficdyn = ' , nficdyn
580         WRITE(numout,*) ' loop on the same FILE                      lperdyn = ' , lperdyn
581         WRITE(numout,*) '  '
582         WRITE(numout,*) ' name of grid_T file                   cfile_grid_T = ', TRIM(cfile_grid_T)   
583         WRITE(numout,*) ' name of grid_U file                   cfile_grid_U = ', TRIM(cfile_grid_U) 
584         WRITE(numout,*) ' name of grid_V file                   cfile_grid_V = ', TRIM(cfile_grid_V) 
585         WRITE(numout,*) ' name of grid_W file                   cfile_grid_W = ', TRIM(cfile_grid_W)     
586         WRITE(numout,*) ' '
587      ENDIF
588
589   END SUBROUTINE dta_dyn_init
590
591   SUBROUTINE wzv( pu, pv, pw, phdiv )
592      !!----------------------------------------------------------------------
593      !!                    ***  ROUTINE wzv  ***
594      !!
595      !! ** Purpose :   Compute the now vertical velocity after the array swap
596      !!
597      !! ** Method  :
598      !! ** Method  : - Divergence:
599      !!      - compute the now divergence given by :
600      !!         * z-coordinate
601      !!         hdiv = 1/(e1t*e2t) [ di(e2u  u) + dj(e1v  v) ]
602      !!     - Using the incompressibility hypothesis, the vertical
603      !!      velocity is computed by integrating the horizontal divergence
604      !!      from the bottom to the surface.
605      !!        The boundary conditions are w=0 at the bottom (no flux) and,
606      !!      in regid-lid case, w=0 at the sea surface.
607      !!
608      !!
609      !! History :
610      !!   9.0  !  02-07  (G. Madec)  Vector optimization
611      !!----------------------------------------------------------------------
612      !! * Arguments
613      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in  )   :: pu, pv    !:  horizontal velocities
614      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out )   :: pw        !:  verticla velocity
615      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: phdiv     !:  horizontal divergence
616
617      !! * Local declarations
618      INTEGER  ::  ji, jj, jk
619      REAL(wp) ::  zu, zu1, zv, zv1, zet
620
621
622      ! Computation of vertical velocity using horizontal divergence
623      phdiv(:,:,:) = 0.
624      DO jk = 1, jpkm1
625         DO jj = 2, jpjm1
626            DO ji = fs_2, fs_jpim1   ! vector opt.
627#if defined key_zco
628               zu  = pu(ji  ,jj  ,jk) * umask(ji  ,jj  ,jk) * e2u(ji  ,jj  )
629               zu1 = pu(ji-1,jj  ,jk) * umask(ji-1,jj  ,jk) * e2u(ji-1,jj  )
630               zv  = pv(ji  ,jj  ,jk) * vmask(ji  ,jj  ,jk) * e1v(ji  ,jj  )
631               zv1 = pv(ji  ,jj-1,jk) * vmask(ji  ,jj-1,jk) * e1v(ji  ,jj-1)
632               zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) )
633#else
634               zu  = pu(ji  ,jj  ,jk) * umask(ji  ,jj  ,jk) * e2u(ji  ,jj  ) * fse3u(ji  ,jj  ,jk)
635               zu1 = pu(ji-1,jj  ,jk) * umask(ji-1,jj  ,jk) * e2u(ji-1,jj  ) * fse3u(ji-1,jj  ,jk)
636               zv  = pv(ji  ,jj  ,jk) * vmask(ji  ,jj  ,jk) * e1v(ji  ,jj  ) * fse3v(ji  ,jj  ,jk)
637               zv1 = pv(ji  ,jj-1,jk) * vmask(ji  ,jj-1,jk) * e1v(ji  ,jj-1) * fse3v(ji  ,jj-1,jk)
638               zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
639#endif
640               phdiv(ji,jj,jk) = ( zu - zu1 + zv - zv1 ) * zet 
641            END DO
642         END DO
643      ENDDO
644
645      ! Lateral boundary conditions on phdiv
646      CALL lbc_lnk( phdiv, 'T', 1. )
647
648
649      ! computation of vertical velocity from the bottom
650      pw(:,:,jpk) = 0.
651      DO jk = jpkm1, 1, -1
652         pw(:,:,jk) = pw(:,:,jk+1) - fse3t(:,:,jk) * phdiv(:,:,jk)
653      END DO
654
655   END SUBROUTINE wzv
656
657   SUBROUTINE tau2wnd( ptaux, ptauy, pwspd )
658      !!---------------------------------------------------------------------
659      !!                    ***  ROUTINE sbc_tau2wnd  ***
660      !!
661      !! ** Purpose : Estimation of wind speed as a function of wind stress
662      !!
663      !! ** Method  : |tau|=rhoa*Cd*|U|^2
664      !!---------------------------------------------------------------------
665      !! * Arguments
666      REAL(wp), DIMENSION(jpi,jpj), INTENT( in  ) ::  &
667         ptaux, ptauy                              !: wind stress in i-j direction resp.
668      REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) ::  &
669         pwspd                                     !: wind speed
670      REAL(wp) ::   zrhoa  = 1.22         ! Air density kg/m3
671      REAL(wp) ::   zcdrag = 1.5e-3       ! drag coefficient
672      REAL(wp) ::   ztx, zty, ztau, zcoef ! temporary variables
673      INTEGER  ::   ji, jj                ! dummy indices
674      !!---------------------------------------------------------------------
675      zcoef = 1. / ( zrhoa * zcdrag )
676!CDIR NOVERRCHK
677      DO jj = 2, jpjm1
678!CDIR NOVERRCHK
679         DO ji = fs_2, fs_jpim1   ! vector opt.
680            ztx = ptaux(ji,jj) * umask(ji,jj,1) + ptaux(ji-1,jj  ) * umask(ji-1,jj  ,1)
681            zty = ptauy(ji,jj) * vmask(ji,jj,1) + ptauy(ji  ,jj-1) * vmask(ji  ,jj-1,1)
682            ztau = 0.5 * SQRT( ztx * ztx + zty * zty )
683            pwspd(ji,jj) = SQRT ( ztau * zcoef ) * tmask(ji,jj,1)
684         END DO
685      END DO
686      CALL lbc_lnk( pwspd(:,:), 'T', 1. )
687
688   END SUBROUTINE tau2wnd
689
690#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
691
692   SUBROUTINE bbl_sign( ptn, psn, pbblx, pbbly )
693      !!----------------------------------------------------------------------
694      !!                    ***  ROUTINE bbl_sign  ***
695      !!
696      !! ** Purpose :   Compute the sign of local gradient of density multiplied by the slope
697      !!                along the bottom slope gradient : grad( rho) * grad(h)
698      !!                Need to compute the diffusive bottom boundary layer
699      !!
700      !! ** Method  :   When the product grad( rho) * grad(h) < 0 (where grad
701      !!      is an along bottom slope gradient) an additional lateral diffu-
702      !!      sive trend along the bottom slope is added to the general tracer
703      !!      trend, otherwise nothing is done. See trcbbl.F90
704      !!
705      !!
706      !! History :
707      !!   9.0  !  02-07  (G. Madec)  Vector optimization
708      !!----------------------------------------------------------------------
709      !! * Arguments
710      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in  ) ::  &
711         ptn             ,  &                           !: temperature
712         psn                                            !: salinity
713      REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) ::  &
714         pbblx , pbbly                                  !: sign of bbl in i-j direction resp.
715     
716      !! * Local declarations
717      INTEGER  ::   ji, jj                   ! dummy loop indices
718      INTEGER  ::   ik
719      REAL(wp) ::   &
720         ztx, zsx, zhx, zalbetx, zgdrhox,     &  ! temporary scalars
721         zty, zsy, zhy, zalbety, zgdrhoy 
722      REAL(wp), DIMENSION(jpi,jpj) ::    &
723        ztnb, zsnb, zdep
724      REAL(wp) ::    fsalbt, pft, pfs, pfh   ! statement function
725      !!----------------------------------------------------------------------
726      ! ratio alpha/beta
727      ! ================
728      !  fsalbt: ratio of thermal over saline expension coefficients
729      !       pft :  potential temperature in degrees celcius
730      !       pfs :  salinity anomaly (s-35) in psu
731      !       pfh :  depth in meters
732
733      fsalbt( pft, pfs, pfh ) =                                              &
734         ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft                    &
735                                   - 0.203814e-03 ) * pft                    &
736                                   + 0.170907e-01 ) * pft                    &
737                                   + 0.665157e-01                            &
738         +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs   &
739         +  ( ( - 0.302285e-13 * pfh                                         &
740                - 0.251520e-11 * pfs                                         &
741                + 0.512857e-12 * pft * pft          ) * pfh                  &
742                                     - 0.164759e-06   * pfs                  &
743             +(   0.791325e-08 * pft - 0.933746e-06 ) * pft                  &
744                                     + 0.380374e-04 ) * pfh
745
746      ! 0. 2D fields of bottom temperature and salinity, and bottom slope
747      ! -----------------------------------------------------------------
748      ! mbathy= number of w-level, minimum value=1 (cf domrea.F90)
749#  if defined key_vectopt_loop
750      jj = 1
751      DO ji = 1, jpij   ! vector opt. (forced unrolling)
752#  else
753      DO jj = 1, jpj
754         DO ji = 1, jpi
755#  endif
756            ik          =  MAX( mbathy(ji,jj) - 1, 1 )    ! vertical index of the bottom ocean T-level
757            ztnb(ji,jj) = ptn(ji,jj,ik) * tmask(ji,jj,1)  ! masked T and S at ocean bottom
758            zsnb(ji,jj) = psn(ji,jj,ik) * tmask(ji,jj,1)
759            zdep(ji,jj) = fsdept(ji,jj,ik)                ! depth of the ocean bottom T-level
760#  if ! defined key_vectopt_loop
761         END DO
762#  endif
763      END DO
764
765      !!----------------------------------------------------------------------
766      ! 1. Criteria of additional bottom diffusivity: grad(rho).grad(h)<0
767      ! --------------------------------------------
768      ! Sign of the local density gradient along the i- and j-slopes
769      ! multiplied by the slope of the ocean bottom
770
771      SELECT CASE ( neos )
772
773      CASE ( 0 )                 ! Jackett and McDougall (1994) formulation
774
775#  if defined key_vectopt_loop
776      jj = 1
777      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
778#  else
779      DO jj = 1, jpjm1
780         DO ji = 1, jpim1
781#  endif
782            ! temperature, salinity anomalie and depth
783            ztx = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) )
784            zsx = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0
785            zhx = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
786            !
787            zty = 0.5 * ( ztnb(ji,jj+1) + ztnb(ji,jj) )
788            zsy = 0.5 * ( zsnb(ji,jj+1) + zsnb(ji,jj) ) - 35.0
789            zhy = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
790            ! masked ratio alpha/beta
791            zalbetx = fsalbt( ztx, zsx, zhx ) * umask(ji,jj,1)
792            zalbety = fsalbt( zty, zsy, zhy ) * vmask(ji,jj,1)
793            ! local density gradient along i-bathymetric slope
794            zgdrhox = zalbetx * ( ztnb(ji+1,jj) - ztnb(ji,jj) )   &
795                   -            ( zsnb(ji+1,jj) - zsnb(ji,jj) )
796            ! local density gradient along j-bathymetric slope
797            zgdrhoy = zalbety * ( ztnb(ji,jj+1) - ztnb(ji,jj) )   &
798                   -            ( zsnb(ji,jj+1) - zsnb(ji,jj) )
799            ! sign of local i-gradient of density multiplied by the i-slope
800            pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
801            ! sign of local j-gradient of density multiplied by the j-slope
802            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
803#  if ! defined key_vectopt_loop
804         END DO
805#  endif
806      END DO
807
808      CASE ( 1 )               ! Linear formulation function of temperature only
809                               !
810#  if defined key_vectopt_loop
811      jj = 1
812      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
813#  else
814      DO jj = 1, jpjm1
815         DO ji = 1, jpim1
816#  endif
817            ! local 'density/temperature' gradient along i-bathymetric slope
818            zgdrhox =  ztnb(ji+1,jj) - ztnb(ji,jj)
819            ! local density gradient along j-bathymetric slope
820            zgdrhoy =  ztnb(ji,jj+1) - ztnb(ji,jj)
821            ! sign of local i-gradient of density multiplied by the i-slope
822            pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
823            ! sign of local j-gradient of density multiplied by the j-slope
824            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
825#  if ! defined key_vectopt_loop
826         END DO
827#  endif
828      END DO
829
830      CASE ( 2 )               ! Linear formulation function of temperature and salinity
831
832#  if defined key_vectopt_loop
833      jj = 1
834      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
835#  else
836      DO jj = 1, jpjm1
837         DO ji = 1, jpim1
838#  endif     
839            ! local density gradient along i-bathymetric slope
840            zgdrhox = - ( rbeta*( zsnb(ji+1,jj) - zsnb(ji,jj) )   &
841                      -  ralpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) )
842            ! local density gradient along j-bathymetric slope
843            zgdrhoy = - ( rbeta*( zsnb(ji,jj+1) - zsnb(ji,jj) )   &
844                      -  ralpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) )
845            ! sign of local i-gradient of density multiplied by the i-slope
846            pbblx(ji,jj) = 0.5 - SIGN( 0.5, - zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
847            ! sign of local j-gradient of density multiplied by the j-slope
848            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
849#  if ! defined key_vectopt_loop
850         END DO
851#  endif
852      END DO
853
854      CASE DEFAULT
855
856         WRITE(ctmp1,*) '          bad flag value for neos = ', neos
857         CALL ctl_stop(ctmp1)
858
859      END SELECT
860   
861      ! Lateral boundary conditions
862      CALL lbc_lnk( pbblx, 'U', 1. )
863      CALL lbc_lnk( pbbly, 'V', 1. )
864
865   END SUBROUTINE bbl_sign
866
867#endif
868
869   SUBROUTINE swap_dyn_data
870      !!----------------------------------------------------------------------
871      !!                    ***  ROUTINE swap_dyn_data  ***
872      !!
873      !! ** Purpose :   swap array data
874      !!
875      !! History :
876      !!   9.0  !  07-09  (C. Ethe)
877      !!----------------------------------------------------------------------
878
879
880      ! swap from record 2 to 1
881      tdta   (:,:,:,1) = tdta   (:,:,:,2)
882      sdta   (:,:,:,1) = sdta   (:,:,:,2)
883      avtdta (:,:,:,1) = avtdta (:,:,:,2)
884      udta   (:,:,:,1) = udta   (:,:,:,2)
885      vdta   (:,:,:,1) = vdta   (:,:,:,2)
886      wdta   (:,:,:,1) = wdta   (:,:,:,2)
887#if defined key_trc_diatrd
888      hdivdta(:,:,:,1) = hdivdta(:,:,:,2)
889#endif
890
891#if defined key_ldfslp
892      uslpdta (:,:,:,1) = uslpdta (:,:,:,2)
893      vslpdta (:,:,:,1) = vslpdta (:,:,:,2)
894      wslpidta(:,:,:,1) = wslpidta(:,:,:,2)
895      wslpjdta(:,:,:,1) = wslpjdta(:,:,:,2)
896#endif
897      hmlddta(:,:,1) = hmlddta(:,:,2) 
898      wspddta(:,:,1) = wspddta(:,:,2) 
899      frlddta(:,:,1) = frlddta(:,:,2) 
900      empdta (:,:,1) = empdta (:,:,2) 
901      qsrdta (:,:,1) = qsrdta (:,:,2) 
902
903#if ! defined key_off_degrad && defined key_traldf_c2d
904      ahtwdta(:,:,1) = ahtwdta(:,:,2)
905#  if defined key_trcldf_eiv
906      aeiwdta(:,:,1) = aeiwdta(:,:,2)
907#  endif
908#endif
909
910#if defined key_off_degrad
911      ahtudta(:,:,:,1) = ahtudta(:,:,:,2)
912      ahtvdta(:,:,:,1) = ahtvdta(:,:,:,2)
913      ahtwdta(:,:,:,1) = ahtwdta(:,:,:,2)
914#  if defined key_trcldf_eiv
915      aeiudta(:,:,:,1) = aeiudta(:,:,:,2)
916      aeivdta(:,:,:,1) = aeivdta(:,:,:,2)
917      aeiwdta(:,:,:,1) = aeiwdta(:,:,:,2)
918#  endif
919#endif
920
921#if defined key_trcbbl_dif || defined key_trcbbl_adv
922      bblxdta(:,:,1) = bblxdta(:,:,2)
923      bblydta(:,:,1) = bblydta(:,:,2)
924#endif
925
926   END SUBROUTINE swap_dyn_data
927
928   SUBROUTINE assign_dyn_data
929      !!----------------------------------------------------------------------
930      !!                    ***  ROUTINE assign_dyn_data  ***
931      !!
932      !! ** Purpose :   Assign dynamical data to the data that have been read
933      !!                without time interpolation
934      !!
935      !!----------------------------------------------------------------------
936     
937      tn (:,:,:) = tdta  (:,:,:,2)
938      sn (:,:,:) = sdta  (:,:,:,2)
939      avt(:,:,:) = avtdta(:,:,:,2)
940     
941      un (:,:,:) = udta  (:,:,:,2) 
942      vn (:,:,:) = vdta  (:,:,:,2)
943      wn (:,:,:) = wdta  (:,:,:,2)
944
945#if defined key_trc_diatrd
946      hdivn(:,:,:) = hdivdta(:,:,:,2)
947#endif
948
949#if defined key_zdfddm
950      avs(:,:,:)   = avtdta (:,:,:,2)
951#endif
952
953     
954#if defined key_ldfslp
955      uslp (:,:,:) = uslpdta (:,:,:,2) 
956      vslp (:,:,:) = vslpdta (:,:,:,2) 
957      wslpi(:,:,:) = wslpidta(:,:,:,2) 
958      wslpj(:,:,:) = wslpjdta(:,:,:,2) 
959#endif
960
961      hmld(:,:) = hmlddta(:,:,2) 
962      wndm(:,:) = wspddta(:,:,2) 
963      fr_i(:,:) = frlddta(:,:,2) 
964      emp (:,:) = empdta (:,:,2) 
965      emps(:,:) = emp(:,:) 
966      qsr (:,:) = qsrdta (:,:,2) 
967
968#if ! defined key_off_degrad && defined key_traldf_c2d   
969      ahtw(:,:) = ahtwdta(:,:,2)
970#  if defined key_trcldf_eiv
971      aeiw(:,:) = aeiwdta(:,:,2)
972#  endif
973#endif
974     
975#if defined key_off_degrad
976      ahtu(:,:,:) = ahtudta(:,:,:,2)
977      ahtv(:,:,:) = ahtvdta(:,:,:,2)
978      ahtw(:,:,:) = ahtwdta(:,:,:,2)
979#  if defined key_trcldf_eiv
980      aeiu(:,:,:) = aeiudta(:,:,:,2)
981      aeiv(:,:,:) = aeivdta(:,:,:,2)
982      aeiw(:,:,:) = aeiwdta(:,:,:,2)
983#  endif
984     
985#endif
986     
987#if defined key_trcbbl_dif ||  defined key_trcbbl_adv
988      bblx(:,:) = bblxdta(:,:,2)
989      bbly(:,:) = bblydta(:,:,2)
990#endif
991
992   END SUBROUTINE assign_dyn_data
993
994   SUBROUTINE linear_interp_dyn_data( pweigh )
995      !!----------------------------------------------------------------------
996      !!                    ***  ROUTINE linear_interp_dyn_data  ***
997      !!
998      !! ** Purpose :   linear interpolation of data
999      !!
1000      !!----------------------------------------------------------------------
1001      !! * Argument
1002      REAL(wp), INTENT( in ) ::   pweigh       ! weigh
1003
1004      !! * Local declarations
1005      REAL(wp) :: zweighm1
1006      !!----------------------------------------------------------------------
1007
1008      zweighm1 = 1. - pweigh
1009     
1010      tn (:,:,:) = zweighm1 * tdta  (:,:,:,1) + pweigh * tdta  (:,:,:,2)
1011      sn (:,:,:) = zweighm1 * sdta  (:,:,:,1) + pweigh * sdta  (:,:,:,2)
1012      avt(:,:,:) = zweighm1 * avtdta(:,:,:,1) + pweigh * avtdta(:,:,:,2)
1013     
1014      un (:,:,:) = zweighm1 * udta  (:,:,:,1) + pweigh * udta  (:,:,:,2) 
1015      vn (:,:,:) = zweighm1 * vdta  (:,:,:,1) + pweigh * vdta  (:,:,:,2)
1016      wn (:,:,:) = zweighm1 * wdta  (:,:,:,1) + pweigh * wdta  (:,:,:,2)
1017
1018#if defined key_trc_diatrd
1019      hdivn(:,:,:) = zweighm1 * hdivdta(:,:,:,1) + pweigh * hdivdta(:,:,:,2)
1020#endif
1021
1022#if defined key_zdfddm
1023      avs(:,:,:)   = zweighm1 * avtdta (:,:,:,1) + pweigh * avtdta (:,:,:,2)
1024#endif
1025
1026     
1027#if defined key_ldfslp
1028      uslp (:,:,:) = zweighm1 * uslpdta (:,:,:,1) + pweigh * uslpdta (:,:,:,2) 
1029      vslp (:,:,:) = zweighm1 * vslpdta (:,:,:,1) + pweigh * vslpdta (:,:,:,2) 
1030      wslpi(:,:,:) = zweighm1 * wslpidta(:,:,:,1) + pweigh * wslpidta(:,:,:,2) 
1031      wslpj(:,:,:) = zweighm1 * wslpjdta(:,:,:,1) + pweigh * wslpjdta(:,:,:,2) 
1032#endif
1033
1034      hmld(:,:) = zweighm1 * hmlddta(:,:,1) + pweigh  * hmlddta(:,:,2) 
1035      wndm(:,:) = zweighm1 * wspddta(:,:,1) + pweigh  * wspddta(:,:,2) 
1036      fr_i(:,:) = zweighm1 * frlddta(:,:,1) + pweigh  * frlddta(:,:,2) 
1037      emp (:,:) = zweighm1 * empdta (:,:,1) + pweigh  * empdta (:,:,2) 
1038      emps(:,:) = emp(:,:) 
1039      qsr (:,:) = zweighm1 * qsrdta (:,:,1) + pweigh  * qsrdta (:,:,2) 
1040
1041#if ! defined key_off_degrad && defined key_traldf_c2d   
1042      ahtw(:,:) = zweighm1 * ahtwdta(:,:,1) + pweigh * ahtwdta(:,:,2)
1043#  if defined key_trcldf_eiv
1044      aeiw(:,:) = zweighm1 * aeiwdta(:,:,1) + pweigh * aeiwdta(:,:,2)
1045#  endif
1046#endif
1047     
1048#if defined key_off_degrad
1049      ahtu(:,:,:) = zweighm1 * ahtudta(:,:,:,1) + pweigh * ahtudta(:,:,:,2)
1050      ahtv(:,:,:) = zweighm1 * ahtvdta(:,:,:,1) + pweigh * ahtvdta(:,:,:,2)
1051      ahtw(:,:,:) = zweighm1 * ahtwdta(:,:,:,1) + pweigh * ahtwdta(:,:,:,2)
1052#  if defined key_trcldf_eiv
1053      aeiu(:,:,:) = zweighm1 * aeiudta(:,:,:,1) + pweigh * aeiudta(:,:,:,2)
1054      aeiv(:,:,:) = zweighm1 * aeivdta(:,:,:,1) + pweigh * aeivdta(:,:,:,2)
1055      aeiw(:,:,:) = zweighm1 * aeiwdta(:,:,:,1) + pweigh * aeiwdta(:,:,:,2)
1056#  endif
1057#endif
1058     
1059#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
1060      bblx(:,:) = zweighm1 * bblxdta(:,:,1) + pweigh * bblxdta(:,:,2)
1061      bbly(:,:) = zweighm1 * bblydta(:,:,1) + pweigh * bblydta(:,:,2)
1062#endif
1063
1064   END SUBROUTINE linear_interp_dyn_data
1065
1066END MODULE dtadyn
Note: See TracBrowser for help on using the repository browser.