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.
trcrst.F90 in branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC – NEMO

source: branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/trcrst.F90 @ 6226

Last change on this file since 6226 was 6226, checked in by jpalmier, 8 years ago

JPALM -- 08-01-2016 -- debugg
File size: 32.1 KB
Line 
1MODULE trcrst
2   !!======================================================================
3   !!                         ***  MODULE trcrst  ***
4   !! TOP :   Manage the passive tracer restart
5   !!======================================================================
6   !! History :    -   !  1991-03  ()  original code
7   !!             1.0  !  2005-03 (O. Aumont, A. El Moussaoui) F90
8   !!              -   !  2005-10 (C. Ethe) print control
9   !!             2.0  !  2005-10 (C. Ethe, G. Madec) revised architecture
10   !!----------------------------------------------------------------------
11#if defined key_top
12   !!----------------------------------------------------------------------
13   !!   'key_top'                                                TOP models
14   !!----------------------------------------------------------------------
15   !!----------------------------------------------------------------------
16   !!   trc_rst :   Restart for passive tracer
17   !!----------------------------------------------------------------------
18   !!----------------------------------------------------------------------
19   !!   'key_top'                                                TOP models
20   !!----------------------------------------------------------------------
21   !!   trc_rst_opn    : open  restart file
22   !!   trc_rst_read   : read  restart file
23   !!   trc_rst_wri    : write restart file
24   !!----------------------------------------------------------------------
25   USE oce_trc
26   USE trc
27   USE trcnam_trp
28   USE iom
29   USE ioipsl, ONLY : ju2ymds    ! for calendar
30   USE daymod
31   !! AXY (05/11/13): need these for MEDUSA to input/output benthic reservoirs
32   USE sms_medusa
33   USE trcsms_medusa
34   !!
35#if defined key_idtra
36   USE trcsms_idtra
37#endif
38   !!
39#if defined key_cfc
40   USE trcsms_cfc
41#endif
42   IMPLICIT NONE
43   PRIVATE
44
45   PUBLIC   trc_rst_opn       ! called by ???
46   PUBLIC   trc_rst_read      ! called by ???
47   PUBLIC   trc_rst_wri       ! called by ???
48   PUBLIC   trc_rst_cal
49   PUBLIC   trc_rst_stat
50
51   !! * Substitutions
52#  include "top_substitute.h90"
53   
54CONTAINS
55   
56   SUBROUTINE trc_rst_opn( kt )
57      !!----------------------------------------------------------------------
58      !!                    ***  trc_rst_opn  ***
59      !!
60      !! ** purpose  :   output of sea-trc variable in a netcdf file
61      !!----------------------------------------------------------------------
62      INTEGER, INTENT(in) ::   kt       ! number of iteration
63      INTEGER             ::   iyear, imonth, iday
64      REAL (wp)           ::   zsec
65      !
66      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step define as a character
67      CHARACTER(LEN=50)   ::   clname   ! trc output restart file name
68      CHARACTER(LEN=256)  ::   clpath   ! full path to ocean output restart file
69      !!----------------------------------------------------------------------
70      !
71      IF( lk_offline ) THEN
72         IF( kt == nittrc000 ) THEN
73            lrst_trc = .FALSE.
74            IF( ln_rst_list ) THEN
75               nrst_lst = 1
76               nitrst = nstocklist( nrst_lst )
77            ELSE
78               nitrst = nitend
79            ENDIF
80         ENDIF
81
82         IF( .NOT. ln_rst_list .AND. MOD( kt - 1, nstock ) == 0 ) THEN
83            ! we use kt - 1 and not kt - nittrc000 to keep the same periodicity from the beginning of the experiment
84            nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing
85            IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run
86         ENDIF
87      ELSE
88         IF( kt == nittrc000 ) lrst_trc = .FALSE.
89      ENDIF
90
91      ! to get better performances with NetCDF format:
92      ! we open and define the tracer restart file one tracer time step before writing the data (-> at nitrst - 2*nn_dttrc + 1)
93      ! except if we write tracer restart files every tracer time step or if a tracer restart file was writen at nitend - 2*nn_dttrc + 1
94      IF( kt == nitrst - 2*nn_dttrc .OR. nstock == nn_dttrc .OR. ( kt == nitend - nn_dttrc .AND. .NOT. lrst_trc ) ) THEN
95         IF ( ln_rstdate ) THEN
96            !! JPALM -- 22-12-2015 -- modif to get the good date on restart trc file name
97            !!                     -- the condition to open the rst file is not the same than for the dynamic rst.
98            !!                     -- here it - for an obscure reason - is open 2 time-step before the restart writing process
99            !!                     instead of 1.
100            !!                     -- i am not sure if someone forgot +1 in the if loop condition as
101            !!                     it is writen in all comments nitrst -2*nn_dttrc + 1 and the condition is
102            !!                     nitrst - 2*nn_dttrc
103            !!                     -- nevertheless we didn't wanted to broke something already working
104            !!                     and just adapted the part we added.
105            !!                     -- So instead of calling ju2ymds( fjulday + (rdttra(1))
106            !!                     we call ju2ymds( fjulday + (2*rdttra(1))
107            !!--------------------------------------------------------------------     
108            CALL ju2ymds( fjulday + (2*rdttra(1)) / rday, iyear, imonth, iday, zsec )
109            WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday
110         ELSE
111            ! beware of the format used to write kt (default is i8.8, that should be large enough)
112            IF( nitrst > 1.0e9 ) THEN   ;   WRITE(clkt,*       ) nitrst
113            ELSE                        ;   WRITE(clkt,'(i8.8)') nitrst
114            ENDIF
115         ENDIF
116         ! create the file
117         IF(lwp) WRITE(numout,*)
118         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_trcrst_out)
119         clpath = TRIM(cn_trcrst_outdir)
120         IF( clpath(LEN_TRIM(clpath):) /= '/' ) clpath = TRIM(clpath) // '/'
121         IF(lwp) WRITE(numout,*) &
122             '             open trc restart.output NetCDF file: ',TRIM(clpath)//clname
123         CALL iom_open( TRIM(clpath)//TRIM(clname), numrtw, ldwrt = .TRUE., kiolib = jprstlib )
124         lrst_trc = .TRUE.
125      ENDIF
126      !
127   END SUBROUTINE trc_rst_opn
128
129   SUBROUTINE trc_rst_read
130      !!----------------------------------------------------------------------
131      !!                    ***  trc_rst_opn  ***
132      !!
133      !! ** purpose  :   read passive tracer fields in restart files
134      !!----------------------------------------------------------------------
135      INTEGER  ::  jn, jl     
136      !! AXY (05/11/13): temporary variables
137      REAL(wp) ::    fq0,fq1,fq2
138
139      !!----------------------------------------------------------------------
140      !
141      IF(lwp) WRITE(numout,*)
142      IF(lwp) WRITE(numout,*) 'trc_rst_read : read data in the TOP restart file'
143      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
144
145      ! READ prognostic variables and computes diagnostic variable
146      DO jn = 1, jptra
147         CALL iom_get( numrtr, jpdom_autoglo, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
148      END DO
149
150      DO jn = 1, jptra
151         CALL iom_get( numrtr, jpdom_autoglo, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
152      END DO
153
154      !! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
155      !!                 call to MEDUSA-2 at this point; this suggests that the FCM
156      !!                 version of NEMO date significantly earlier than the current
157      !!                 version
158
159#if defined key_medusa
160      !! AXY (13/01/12): check if the restart contains sediment fields;
161      !!                 this is only relevant for simulations that include
162      !!                 biogeochemistry and are restarted from earlier runs
163      !!                 in which there was no sediment component
164      !!
165      IF( iom_varid( numrtr, 'B_SED_N', ldstop = .FALSE. ) > 0 ) THEN
166         !! YES; in which case read them
167         !!
168         IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields present - reading in ...'
169         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_N',  zb_sed_n(:,:)  )
170         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_N',  zn_sed_n(:,:)  )
171         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_FE', zb_sed_fe(:,:) )
172         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_FE', zn_sed_fe(:,:) )
173         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_SI', zb_sed_si(:,:) )
174         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_SI', zn_sed_si(:,:) )
175         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_C',  zb_sed_c(:,:)  )
176         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_C',  zn_sed_c(:,:)  )
177         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_CA', zb_sed_ca(:,:) )
178         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_CA', zn_sed_ca(:,:) )
179      ELSE
180         !! NO; in which case set them to zero
181         !!
182         IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields absent - setting to zero ...'
183         zb_sed_n(:,:)  = 0.0   !! organic N
184         zn_sed_n(:,:)  = 0.0
185         zb_sed_fe(:,:) = 0.0   !! organic Fe
186         zn_sed_fe(:,:) = 0.0
187         zb_sed_si(:,:) = 0.0   !! inorganic Si
188         zn_sed_si(:,:) = 0.0
189         zb_sed_c(:,:)  = 0.0   !! organic C
190         zn_sed_c(:,:)  = 0.0
191         zb_sed_ca(:,:) = 0.0   !! inorganic C
192         zn_sed_ca(:,:) = 0.0
193      ENDIF
194      !!
195      !! calculate stats on these fields
196      IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
197      fq0 = MINVAL(zn_sed_n(:,:))
198      fq1 = MAXVAL(zn_sed_n(:,:))
199      fq2 = SUM(zn_sed_n(:,:))
200      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  N ', &
201         &        fq0, fq1, fq2
202      fq0 = MINVAL(zn_sed_fe(:,:))
203      fq1 = MAXVAL(zn_sed_fe(:,:))
204      fq2 = SUM(zn_sed_fe(:,:))
205      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Fe ', &
206         &        fq0, fq1, fq2
207      fq0 = MINVAL(zn_sed_si(:,:))
208      fq1 = MAXVAL(zn_sed_si(:,:))
209      fq2 = SUM(zn_sed_si(:,:))
210      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Si ', &
211         &        fq0, fq1, fq2
212      fq0 = MINVAL(zn_sed_c(:,:))
213      fq1 = MAXVAL(zn_sed_c(:,:))
214      fq2 = SUM(zn_sed_c(:,:))
215      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  C ', &
216         &        fq0, fq1, fq2
217      fq0 = MINVAL(zn_sed_ca(:,:))
218      fq1 = MAXVAL(zn_sed_ca(:,:))
219      fq2 = SUM(zn_sed_ca(:,:))
220      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Ca ', &
221         &        fq0, fq1, fq2
222      !!
223      !! AXY (07/07/15): read in temporally averaged fields for DMS
224      !!                 calculations
225      !!
226      IF( iom_varid( numrtr, 'B_DMS_CHN', ldstop = .FALSE. ) > 0 ) THEN
227         !! YES; in which case read them
228         !!
229         IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS present - reading in ...'
230         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHN',  zb_dms_chn(:,:)  )
231         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHN',  zn_dms_chn(:,:)  )
232         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHD',  zb_dms_chd(:,:)  )
233         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHD',  zn_dms_chd(:,:)  )
234         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_MLD',  zb_dms_mld(:,:)  )
235         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_MLD',  zn_dms_mld(:,:)  )
236         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_QSR',  zb_dms_qsr(:,:)  )
237         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_QSR',  zn_dms_qsr(:,:)  )
238         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_DIN',  zb_dms_din(:,:)  )
239         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_DIN',  zn_dms_din(:,:)  )
240      ELSE
241         !! NO; in which case set them to zero
242         !!
243         IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS absent - setting to zero ...'
244         zb_dms_chn(:,:)  = 0.0   !! CHN
245         zn_dms_chn(:,:)  = 0.0
246         zb_dms_chd(:,:)  = 0.0   !! CHD
247         zn_dms_chd(:,:)  = 0.0
248         zb_dms_mld(:,:)  = 0.0   !! MLD
249         zn_dms_mld(:,:)  = 0.0
250         zb_dms_qsr(:,:)  = 0.0   !! QSR
251         zn_dms_qsr(:,:)  = 0.0
252         zb_dms_din(:,:)  = 0.0   !! DIN
253         zn_dms_din(:,:)  = 0.0
254      ENDIF
255      !!
256      !! calculate stats on these fields
257      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
258      fq0 = MINVAL(zn_dms_chn(:,:))
259      fq1 = MAXVAL(zn_dms_chn(:,:))
260      fq2 = SUM(zn_dms_chn(:,:))
261      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHN ', fq0, fq1, fq2
262      fq0 = MINVAL(zn_dms_chd(:,:))
263      fq1 = MAXVAL(zn_dms_chd(:,:))
264      fq2 = SUM(zn_dms_chd(:,:))
265      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHD ', fq0, fq1, fq2
266      fq0 = MINVAL(zn_dms_mld(:,:))
267      fq1 = MAXVAL(zn_dms_mld(:,:))
268      fq2 = SUM(zn_dms_mld(:,:))
269      if (lwp) write (numout,'(a,3f15.5)') 'DMS, MLD ', fq0, fq1, fq2
270      fq0 = MINVAL(zn_dms_qsr(:,:))
271      fq1 = MAXVAL(zn_dms_qsr(:,:))
272      fq2 = SUM(zn_dms_qsr(:,:))
273      if (lwp) write (numout,'(a,3f15.5)') 'DMS, QSR ', fq0, fq1, fq2
274      fq0 = MINVAL(zn_dms_din(:,:))
275      fq1 = MAXVAL(zn_dms_din(:,:))
276      fq2 = SUM(zn_dms_din(:,:))
277      if (lwp) write (numout,'(a,3f15.5)') 'DMS, DIN ', fq0, fq1, fq2
278#endif
279      !
280#if defined key_idtra
281      !! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
282      !!                        writting here undre their key.
283      !!                        problems in CFC restart, maybe because of this...
284      !!                        and pb in idtra diag or diad-restart writing.
285      !!----------------------------------------------------------------------
286      IF( iom_varid( numrtr, 'qint_idtra', ldstop = .FALSE. ) > 0 ) THEN
287         !! YES; in which case read them
288         !!
289         IF(lwp) WRITE(numout,*) ' IDTRA averaged properties present - reading in ...'
290         CALL iom_get( numrtr, jpdom_autoglo, 'qint_IDTRA',  qint_idtra(:,:,1)  )
291      ELSE
292         !! NO; in which case set them to zero
293         !!
294         IF(lwp) WRITE(numout,*) ' IDTRA averaged properties absent - setting to zero ...'
295         qint_idtra(:,:,1)  = 0.0   !! CHN
296      ENDIF
297      !!
298      !! calculate stats on these fields
299      call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
300      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
301      fq0 = MINVAL(qint_idtra(:,:,1))
302      fq1 = MAXVAL(qint_idtra(:,:,1))
303      fq2 = SUM(qint_idtra(:,:,1))
304      if (lwp) write (numout,'(a,3f15.5)') 'qint_IDTRA ', fq0, fq1, fq2
305#endif
306      !
307#if defined key_cfc
308      DO jl = 1, jp_cfc
309         jn = jp_cfc0 + jl - 1
310         IF( iom_varid( numrtr, 'qint_'//ctrcnm(jn), ldstop = .FALSE. ) > 0 ) THEN
311            !! YES; in which case read them
312            !!
313            IF(lwp) WRITE(numout,*) ' CFC averaged properties present - reading in ...'
314            CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
315         ELSE
316            !! NO; in which case set them to zero
317            !!
318            IF(lwp) WRITE(numout,*) ' CFC averaged properties absent - setting to zero ...'
319            qint_cfc(:,:,jn)  = 0.0   !! CHN
320         ENDIF
321         !!
322         !! calculate stats on these fields
323         call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
324         IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
325         fq0 = MINVAL(qint_cfc(:,:,jl))
326         fq1 = MAXVAL(qint_cfc(:,:,jl))
327         fq2 = SUM(qint_cfc(:,:,jl))
328         if (lwp) write (numout,'(a,3f15.5)') 'qint_'//ctrcnm(jn), fq0, fq1, fq2
329      END DO
330#endif
331      !
332   END SUBROUTINE trc_rst_read
333
334   SUBROUTINE trc_rst_wri( kt )
335      !!----------------------------------------------------------------------
336      !!                    ***  trc_rst_wri  ***
337      !!
338      !! ** purpose  :   write passive tracer fields in restart files
339      !!----------------------------------------------------------------------
340      INTEGER, INTENT( in ) ::   kt    ! ocean time-step index
341      !!
342      INTEGER  :: jn, jl
343      REAL(wp) :: zarak0
344      !! AXY (05/11/13): temporary variables
345      REAL(wp) ::    fq0,fq1,fq2
346      !!----------------------------------------------------------------------
347      !
348      CALL iom_rstput( kt, nitrst, numrtw, 'rdttrc1', rdttrc(1) )   ! surface passive tracer time step
349      ! prognostic variables
350      ! --------------------
351      DO jn = 1, jptra
352         CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
353      END DO
354
355      DO jn = 1, jptra
356         CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
357      END DO
358
359      !! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
360      !!                 call to MEDUSA-2 at this point; this suggests that the FCM
361      !!                 version of NEMO date significantly earlier than the current
362      !!                 version
363
364#if defined key_medusa
365      !! AXY (13/01/12): write out "before" and "now" state of seafloor
366      !!                 sediment pools into restart; this happens
367      !!                 whether or not the pools are to be used by
368      !!                 MEDUSA (which is controlled by a switch in the
369      !!                 namelist_top file)
370      !!
371      IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields - writing out ...'
372      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_N',  zb_sed_n(:,:)  )
373      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_N',  zn_sed_n(:,:)  )
374      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_FE', zb_sed_fe(:,:) )
375      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_FE', zn_sed_fe(:,:) )
376      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_SI', zb_sed_si(:,:) )
377      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_SI', zn_sed_si(:,:) )
378      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_C',  zb_sed_c(:,:)  )
379      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_C',  zn_sed_c(:,:)  )
380      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_CA', zb_sed_ca(:,:) )
381      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_CA', zn_sed_ca(:,:) )
382      !!
383      !! calculate stats on these fields
384      IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
385      fq0 = MINVAL(zn_sed_n(:,:))
386      fq1 = MAXVAL(zn_sed_n(:,:))
387      fq2 = SUM(zn_sed_n(:,:))
388      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  N ', &
389         &        fq0, fq1, fq2
390      fq0 = MINVAL(zn_sed_fe(:,:))
391      fq1 = MAXVAL(zn_sed_fe(:,:))
392      fq2 = SUM(zn_sed_fe(:,:))
393      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Fe ', &
394         &        fq0, fq1, fq2
395      fq0 = MINVAL(zn_sed_si(:,:))
396      fq1 = MAXVAL(zn_sed_si(:,:))
397      fq2 = SUM(zn_sed_si(:,:))
398      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Si ', &
399         &        fq0, fq1, fq2
400      fq0 = MINVAL(zn_sed_c(:,:))
401      fq1 = MAXVAL(zn_sed_c(:,:))
402      fq2 = SUM(zn_sed_c(:,:))
403      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  C ', &
404         &        fq0, fq1, fq2
405      fq0 = MINVAL(zn_sed_ca(:,:))
406      fq1 = MAXVAL(zn_sed_ca(:,:))
407      fq2 = SUM(zn_sed_ca(:,:))
408      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Ca ', &
409         &        fq0, fq1, fq2
410      !!
411      !! AXY (07/07/15): write out temporally averaged fields for DMS
412      !!                 calculations
413      !!
414      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS - writing out ...'
415      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHN',  zb_dms_chn(:,:)  )
416      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHN',  zn_dms_chn(:,:)  )
417      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHD',  zb_dms_chd(:,:)  )
418      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHD',  zn_dms_chd(:,:)  )
419      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_MLD',  zb_dms_mld(:,:)  )
420      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_MLD',  zn_dms_mld(:,:)  )
421      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_QSR',  zb_dms_qsr(:,:)  )
422      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_QSR',  zn_dms_qsr(:,:)  )
423      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_DIN',  zb_dms_din(:,:)  )
424      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_DIN',  zn_dms_din(:,:)  )
425      !!
426      !! calculate stats on these fields
427      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
428      fq0 = MINVAL(zn_dms_chn(:,:))
429      fq1 = MAXVAL(zn_dms_chn(:,:))
430      fq2 = SUM(zn_dms_chn(:,:))
431      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHN ', fq0, fq1, fq2
432      fq0 = MINVAL(zn_dms_chd(:,:))
433      fq1 = MAXVAL(zn_dms_chd(:,:))
434      fq2 = SUM(zn_dms_chd(:,:))
435      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHD ', fq0, fq1, fq2
436      fq0 = MINVAL(zn_dms_mld(:,:))
437      fq1 = MAXVAL(zn_dms_mld(:,:))
438      fq2 = SUM(zn_dms_mld(:,:))
439      if (lwp) write (numout,'(a,3f15.5)') 'DMS, MLD ', fq0, fq1, fq2
440      fq0 = MINVAL(zn_dms_qsr(:,:))
441      fq1 = MAXVAL(zn_dms_qsr(:,:))
442      fq2 = SUM(zn_dms_qsr(:,:))
443      if (lwp) write (numout,'(a,3f15.5)') 'DMS, QSR ', fq0, fq1, fq2
444      fq0 = MINVAL(zn_dms_din(:,:))
445      fq1 = MAXVAL(zn_dms_din(:,:))
446      fq2 = SUM(zn_dms_din(:,:))
447      if (lwp) write (numout,'(a,3f15.5)') 'DMS, DIN ', fq0, fq1, fq2
448      !!
449#endif
450      !
451#if defined key_idtra
452      !! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
453      !!                        writting here undre their key.
454      !!                        problems in CFC restart, maybe because of this...
455      !!                        and pb in idtra diag or diad-restart writing.
456      !!----------------------------------------------------------------------
457      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties - writing out ...'
458      CALL iom_rstput( kt, nitrst, numrtw, 'qint_IDTRA',  qint_idtra(:,:,1) )
459      !!
460      !! calculate stats on these fields
461      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
462      call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
463      fq0 = MINVAL(qint_idtra(:,:,1))
464      fq1 = MAXVAL(qint_idtra(:,:,1))
465      fq2 = SUM(qint_idtra(:,:,1))
466      if (lwp) write (numout,'(a,3f15.5)') 'qint_IDTRA ', fq0, fq1, fq2
467#endif
468      !
469#if defined key_cfc
470      DO jl = 1, jp_cfc
471         jn = jp_cfc0 + jl - 1
472         IF(lwp) WRITE(numout,*) ' CFC averaged properties - writing out ...'
473         CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
474         !!
475         !! calculate stats on these fields
476         call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
477         IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
478         fq0 = MINVAL(qint_cfc(:,:,jl))
479         fq1 = MAXVAL(qint_cfc(:,:,jl))
480         fq2 = SUM(qint_cfc(:,:,jl))
481         if (lwp) write (numout,'(a,3f15.5)') 'qint_'//ctrcnm(jn), fq0, fq1, fq2
482      END DO
483#endif
484      !
485
486      IF( kt == nitrst ) THEN
487          CALL trc_rst_stat            ! statistics
488          CALL iom_close( numrtw )     ! close the restart file (only at last time step)
489#if ! defined key_trdmxl_trc
490          lrst_trc = .FALSE.
491#endif
492          IF( lk_offline .AND. ln_rst_list ) THEN
493             nrst_lst = nrst_lst + 1
494             nitrst = nstocklist( nrst_lst )
495          ENDIF
496      ENDIF
497      !
498   END SUBROUTINE trc_rst_wri 
499
500
501   SUBROUTINE trc_rst_cal( kt, cdrw )
502      !!---------------------------------------------------------------------
503      !!                   ***  ROUTINE trc_rst_cal  ***
504      !!
505      !!  ** Purpose : Read or write calendar in restart file:
506      !!
507      !!  WRITE(READ) mode:
508      !!       kt        : number of time step since the begining of the experiment at the
509      !!                   end of the current(previous) run
510      !!       adatrj(0) : number of elapsed days since the begining of the experiment at the
511      !!                   end of the current(previous) run (REAL -> keep fractions of day)
512      !!       ndastp    : date at the end of the current(previous) run (coded as yyyymmdd integer)
513      !!
514      !!   According to namelist parameter nrstdt,
515      !!       nn_rsttr = 0  no control on the date (nittrc000 is  arbitrary).
516      !!       nn_rsttr = 1  we verify that nittrc000 is equal to the last
517      !!                   time step of previous run + 1.
518      !!       In both those options, the  exact duration of the experiment
519      !!       since the beginning (cumulated duration of all previous restart runs)
520      !!       is not stored in the restart and is assumed to be (nittrc000-1)*rdt.
521      !!       This is valid is the time step has remained constant.
522      !!
523      !!       nn_rsttr = 2  the duration of the experiment in days (adatrj)
524      !!                    has been stored in the restart file.
525      !!----------------------------------------------------------------------
526      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
527      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
528      !
529      INTEGER  ::  jlibalt = jprstlib
530      LOGICAL  ::  llok
531      REAL(wp) ::  zkt, zrdttrc1
532      REAL(wp) ::  zndastp
533
534      ! Time domain : restart
535      ! ---------------------
536
537      IF( TRIM(cdrw) == 'READ' ) THEN
538
539         IF(lwp) WRITE(numout,*)
540         IF(lwp) WRITE(numout,*) 'trc_rst_cal : read the TOP restart file for calendar'
541         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
542
543         IF ( jprstlib == jprstdimg ) THEN
544           ! eventually read netcdf file (monobloc)  for restarting on different number of processors
545           ! if {cn_trcrst_in}.nc exists, then set jlibalt to jpnf90
546           INQUIRE( FILE = TRIM(cn_trcrst_indir)//'/'//TRIM(cn_trcrst_in)//'.nc', EXIST = llok )
547           IF ( llok ) THEN ; jlibalt = jpnf90  ; ELSE ; jlibalt = jprstlib ; ENDIF
548         ENDIF
549
550         IF( ln_rsttr ) THEN
551            CALL iom_open( TRIM(cn_trcrst_indir)//'/'//cn_trcrst_in, numrtr, kiolib = jlibalt )
552            CALL iom_get ( numrtr, 'kt', zkt )   ! last time-step of previous run
553
554            IF(lwp) THEN
555               WRITE(numout,*) ' *** Info read in restart : '
556               WRITE(numout,*) '   previous time-step                               : ', NINT( zkt )
557               WRITE(numout,*) ' *** restart option'
558               SELECT CASE ( nn_rsttr )
559               CASE ( 0 )   ;   WRITE(numout,*) ' nn_rsttr = 0 : no control of nittrc000'
560               CASE ( 1 )   ;   WRITE(numout,*) ' nn_rsttr = 1 : no control the date at nittrc000 (use ndate0 read in the namelist)'
561               CASE ( 2 )   ;   WRITE(numout,*) ' nn_rsttr = 2 : calendar parameters read in restart'
562               END SELECT
563               WRITE(numout,*)
564            ENDIF
565            ! Control of date
566            IF( nittrc000  - NINT( zkt ) /= nn_dttrc .AND.  nn_rsttr /= 0 )                                  &
567               &   CALL ctl_stop( ' ===>>>> : problem with nittrc000 for the restart',                 &
568               &                  ' verify the restart file or rerun with nn_rsttr = 0 (namelist)' )
569         ENDIF
570         !
571         IF( lk_offline ) THEN   
572            !                                          ! set the date in offline mode
573            IF( ln_rsttr .AND. nn_rsttr == 2 ) THEN
574               CALL iom_get( numrtr, 'ndastp', zndastp ) 
575               ndastp = NINT( zndastp )
576               CALL iom_get( numrtr, 'adatrj', adatrj  )
577             ELSE
578               ndastp = ndate0 - 1     ! ndate0 read in the namelist in dom_nam
579               adatrj = ( REAL( nittrc000-1, wp ) * rdttra(1) ) / rday
580               ! note this is wrong if time step has changed during run
581            ENDIF
582            !
583            IF(lwp) THEN
584              WRITE(numout,*) ' *** Info used values : '
585              WRITE(numout,*) '   date ndastp                                      : ', ndastp
586              WRITE(numout,*) '   number of elapsed days since the begining of run : ', adatrj
587              WRITE(numout,*)
588            ENDIF
589            !
590            IF( ln_rsttr )  THEN   ;    neuler = 1
591            ELSE                   ;    neuler = 0
592            ENDIF
593            !
594            CALL day_init          ! compute calendar
595            !
596         ENDIF
597         !
598      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
599         !
600         IF(  kt == nitrst ) THEN
601            IF(lwp) WRITE(numout,*)
602            IF(lwp) WRITE(numout,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
603            IF(lwp) WRITE(numout,*) '~~~~~~~'
604         ENDIF
605         CALL iom_rstput( kt, nitrst, numrtw, 'kt'     , REAL( kt    , wp) )   ! time-step
606         CALL iom_rstput( kt, nitrst, numrtw, 'ndastp' , REAL( ndastp, wp) )   ! date
607         CALL iom_rstput( kt, nitrst, numrtw, 'adatrj' , adatrj            )   ! number of elapsed days since
608         !                                                                     ! the begining of the run [s]
609      ENDIF
610
611   END SUBROUTINE trc_rst_cal
612
613
614   SUBROUTINE trc_rst_stat
615      !!----------------------------------------------------------------------
616      !!                    ***  trc_rst_stat  ***
617      !!
618      !! ** purpose  :   Compute tracers statistics
619      !!----------------------------------------------------------------------
620      INTEGER  :: jk, jn
621      REAL(wp) :: ztraf, zmin, zmax, zmean, zdrift
622      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol
623      !!----------------------------------------------------------------------
624
625      IF( lwp ) THEN
626         WRITE(numout,*) 
627         WRITE(numout,*) '           ----TRACER STAT----             '
628         WRITE(numout,*) 
629      ENDIF
630      !
631      DO jk = 1, jpk
632         zvol(:,:,jk) = e1e2t(:,:) * fse3t_a(:,:,jk) * tmask(:,:,jk)
633      END DO
634      !
635      DO jn = 1, jptra
636         ztraf = glob_sum( trn(:,:,:,jn) * zvol(:,:,:) )
637         zmin  = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
638         zmax  = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
639         IF( lk_mpp ) THEN
640            CALL mpp_min( zmin )      ! min over the global domain
641            CALL mpp_max( zmax )      ! max over the global domain
642         END IF
643         zmean  = ztraf / areatot
644         zdrift = ( ( ztraf - trai(jn) ) / ( trai(jn) + 1.e-12 )  ) * 100._wp
645         IF(lwp) WRITE(numout,9000) jn, TRIM( ctrcnm(jn) ), zmean, zmin, zmax, zdrift
646      END DO
647      IF(lwp) WRITE(numout,*) 
6489000  FORMAT(' tracer nb :',i2,'    name :',a10,'    mean :',e18.10,'    min :',e18.10, &
649      &      '    max :',e18.10,'    drift :',e18.10, ' %')
650      !
651   END SUBROUTINE trc_rst_stat
652
653
654   SUBROUTINE trc_rst_dia_stat( dgtr, names)
655      !!----------------------------------------------------------------------
656      !!                    ***  trc_rst_dia_stat  ***
657      !!
658      !! ** purpose  :   Compute tracers statistics
659      !!----------------------------------------------------------------------
660      REAL(wp), DIMENSION(jpi,jpj) , INTENT(in) ::   dgtr      ! 2D diag var
661      CHARACTER(len=*)             , INTENT(in) ::   names     ! 2D diag name
662      !!---------------------------------------------------------------------
663      INTEGER  :: jk, jn
664      REAL(wp) :: ztraf, zmin, zmax, zmean, areasf
665      REAL(wp), DIMENSION(jpi,jpj) :: zvol
666      !!----------------------------------------------------------------------
667
668      IF( lwp )  WRITE(numout,*) 'STAT- ', names
669      !
670      zvol(:,:) = e1e2t(:,:) * fse3t_a(:,:,1) * tmask(:,:,1)
671      ztraf = glob_sum( dgtr(:,:) * zvol(:,:) )
672      areasf = glob_sum(e1e2t(:,:) * tmask(:,:,1) )
673      zmin  = MINVAL( dgtr(:,:), mask= ((tmask(:,:,1)*SPREAD(tmask_i,DIM=2,NCOPIES=1).NE.0.)) )
674      zmax  = MAXVAL( dgtr(:,:), mask= ((tmask(:,:,1)*SPREAD(tmask_i,DIM=2,NCOPIES=1).NE.0.)) )
675      IF( lk_mpp ) THEN
676         CALL mpp_min( zmin )      ! min over the global domain
677         CALL mpp_max( zmax )      ! max over the global domain
678      END IF
679      zmean  = ztraf / areatot
680      IF(lwp) WRITE(numout,9001) TRIM( names ), zmean, zmin, zmax
681      !
682      IF(lwp) WRITE(numout,*)
6839001  FORMAT(' tracer name :',a10,'    mean :',e18.10,'    min :',e18.10, &
684      &      '    max :',e18.10 )
685      !
686   END SUBROUTINE trc_rst_dia_stat
687
688
689#else
690   !!----------------------------------------------------------------------
691   !!  Dummy module :                                     No passive tracer
692   !!----------------------------------------------------------------------
693CONTAINS
694   SUBROUTINE trc_rst_read                      ! Empty routines
695   END SUBROUTINE trc_rst_read
696   SUBROUTINE trc_rst_wri( kt )
697      INTEGER, INTENT ( in ) :: kt
698      WRITE(*,*) 'trc_rst_wri: You should not have seen this print! error?', kt
699   END SUBROUTINE trc_rst_wri   
700#endif
701
702   !!----------------------------------------------------------------------
703   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
704   !! $Id$
705   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
706   !!======================================================================
707END MODULE trcrst
Note: See TracBrowser for help on using the repository browser.