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.
bio_medusa_diag.F90 in branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA – NEMO

source: branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_diag.F90 @ 7998

Last change on this file since 7998 was 7998, checked in by marc, 7 years ago

Pull out diagnostic calculations from trcbio_medusa.F90

File size: 84.4 KB
Line 
1MODULE bio_medusa_diag_mod
2   !!======================================================================
3   !!                         ***  MODULE bio_medusa_diag_mod  ***
4   !! Calculates diagnostics
5   !!======================================================================
6   !! History :
7   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90
8   !!----------------------------------------------------------------------
9#if defined key_medusa
10   !!----------------------------------------------------------------------
11   !!                                                   MEDUSA bio-model
12   !!----------------------------------------------------------------------
13
14   IMPLICIT NONE
15   PRIVATE
16     
17   PUBLIC   bio_medusa_diag        ! Called in trcbio_medusa.F90
18
19   !!----------------------------------------------------------------------
20   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
21   !! $Id$
22   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
23   !!----------------------------------------------------------------------
24
25CONTAINS
26
27   SUBROUTINE bio_medusa_diag( kt, jk )
28      !!-------------------------------------------------------------------
29      !!                     ***  ROUTINE bio_medusa_diag  ***
30      !! This called from TRC_BIO_MEDUSA and calculates diagnostics
31      !!-------------------------------------------------------------------
32      USE bio_medusa_mod,    ONLY: dcalc3, dms_andr, dms_aran,           &
33                                   dms_hall, dms_simo, dms_surf,         &
34                                   f_benout_c, f_benout_ca, f_benout_fe, &
35                                   f_benout_lyso_ca, f_benout_n,         &
36                                   f_benout_si,                          &
37                                   f_co2flux, f_co3,                     &
38                                   f_fbenin_c, f_fbenin_ca, f_fbenin_fe, &
39                                   f_fbenin_n, f_fbenin_si,              &
40                                   f_h2co3, f_hco3,                      &
41                                   f_kw660,                              &
42                                   f_o2flux, f_o2sat, f_omarg, f_omcal,  &
43                                   f_pco2atm, f_pco2w, f_ph, f_pp0,      &
44                                   f_riv_loc_alk, f_riv_loc_c,           &
45                                   f_riv_loc_n, f_riv_loc_si,            &
46                                   f_runoff,                             &
47                                   f_sbenin_c, f_sbenin_fe,              &
48                                   f_sbenin_n,                           &
49                                   f_TALK, f_TDIC,                       &
50                                   fccd, fcomm_resp,                     &
51                                   fd_cal3, fd_car3, fd_nit3, fd_sil3,   &
52                                   fdd, fdd2d, fddc,                     &
53                                   fdpd, fdpd2, fdpd22d, fdpd2d,         &
54                                   fdpn, fdpn2, fdpn22d, fdpn2d,         &
55                                   fdzme, fdzme2, fdzme22d, fdzme2d,     &
56                                   fdzmi, fdzmi2, fdzmi22d, fdzmi2d,     &
57                                   fediss3, fescav3,                     &
58                                   ffastc, ffastca, ffastca2d,           &
59                                   ffastn, ffastsi,                      &
60                                   ffebot, ffebot2d, ffescav, ffescav2d, &
61                                   ffetop, ffetop2d,                     &
62                                   ffld, ffld2d, ffln2, ffln2d,          &
63                                   fgmed, fgmed2d, fgmedc,               &
64                                   fgmepd, fgmepd2d, fgmepn, fgmepn2d,   &
65                                   fgmezmi, fgmezmi2d,                   &
66                                   fgmid, fgmid2d, fgmidc,               &
67                                   fgmipn, fgmipn2d,                     &
68                                   ficme, ficmi, finme, finmi,           &
69                                   fjlim_pd, fjlim_pn, fjld2d, fjln2d,   &
70                                   fmeexcr, fmegrow, fmeresp,            &
71                                   fmiexcr, fmigrow, fmiresp,            &
72                                   fnld, fnld2d, fnln, fnln2d,           &
73                                   fprd, fprd_ml, fprd2d,                &
74                                   fprds, fprds2d,                       & 
75                                   fprn, fprn_ml, fprn2d,                &
76                                   fregen2d,                     &
77                                   fregenfast, fregenfastsi,             &
78                                   fregensi2d,                 &
79                                   freminc, freminc2d,                   &
80                                   freminca, freminca2d,                 & 
81                                   freminfe, freminfe2d,                 &
82                                   freminn, freminn2d,                   &
83                                   freminsi, freminsi2d,                 &
84                                   fsdiss, fsdiss2d,                     &
85                                   fsedc, fsedca, fsedfe, fsedn, fsedsi, &
86                                   fsld, fsld2, fsld2d, fsld2d2,         &
87                                   fslown, fslowc, fslown2d,             &
88                                   fslowc2d, fslowcflux, fslownflux,     &
89                                   ftempc, ftempc2d, ftempca, ftempca2d, &
90                                   ftempfe, ftempfe2d,                   &
91                                   ftempn, ftempn2d, ftempsi, ftempsi2d, &
92                                   ftot_a, ftot_c, ftot_fe, ftot_n,      &
93                                   ftot_o2, ftot_si,                     &
94                                   gmedc2d, gmidc2d,                     &
95                                   iben_c2d, iben_ca2d, iben_fe2d,       &
96                                   iben_n2d, iben_si2d,                  &
97                                   iters, lyso_ca2d,                     &
98                                   mdetc2d, megrazd3, megrazp3,          &
99                                   megrazz3,                             & 
100                                   migrazd3, migrazp3,                   &
101                                   oben_c2d, oben_ca2d, oben_fe2d,       &
102                                   oben_n2d, oben_si2d,                  &
103                                   pbsi3, pcal3,                         &
104                                   pdlimfe3, pdlimj3, pdlimn3,           &
105                                   pnlimfe3, pnlimj3, pnlimn3, pdlimsi3, &
106                                   remin3dn, remoc3,                     &
107                                   rivalk2d, rivc2d, rivn2d, rivsi2d,    &
108                                   sfr_oarg2d, sfr_ocal2d,               &
109                                   tpp3d, tppd3, xFree,                  &
110                                   zeexcr2d, zegrow2d, zemesc2d,         &
111                                   zemesd2d, zemesdc2d, zemesn2d,        &
112                                   zeresp2d,                             &
113                                   ziexcr2d, zigrow2d, zimesc2d,         &
114                                   zimesd2d, zimesdc2d, zimesn2d,        &
115                                   ziresp2d,                             &
116                                   zpds, zphd, zphn
117      USE dom_oce,           ONLY: e3t_0, e3t_n, mbathy, tmask
118      USE in_out_manager,    ONLY: lwp, numout
119      USE iom,               ONLY: lk_iomput
120      USE par_kind,          ONLY: wp
121      USE par_oce,           ONLY: jpi, jpim1, jpj, jpjm1
122      USE phycst,            ONLY: rsmall
123      USE sbc_oce,           ONLY: qsr, wndm
124      USE sms_medusa,        ONLY: f2_ccd_arg, f2_ccd_cal,               &
125                                   f3_omarg, f3_omcal, f3_pH,            &
126                                   i0100, i0150, i0200, i0500, i1000,    &
127                                   jdms, ocal_ccd,                       &
128                                   xbetac, xbetan, xpar, xphi,           &
129                                   xthetapd, xthetapn, xthetazmi, xze
130      USE trc,               ONLY: ln_diatrc, med_diag, trc2d, trc3d 
131
132   !!* Substitution
133#  include "domzgr_substitute.h90"
134
135      !! time (integer timestep)
136      INTEGER, INTENT( in ) :: kt
137      !! level
138      INTEGER, INTENT( in ) :: jk
139
140      !! Loop avariables
141      INTEGER :: ji, jj, jn
142
143      REAL(wp), DIMENSION(jpi,jpj) :: fregen,fregensi
144
145# if defined key_trc_diabio
146      !!==========================================================
147      !! LOCAL GRID CELL DIAGNOSTICS
148      !!==========================================================
149      !!
150      !!----------------------------------------------------------
151      !! Full diagnostics key_trc_diabio:
152      !! LOBSTER and PISCES support full diagnistics option
153      !! key_trc_diabio which gives an option of FULL output of
154      !! biological sourses and sinks. I cannot see any reason
155      !! for doing this. If needed, it can be done as shown
156      !! below.
157      !!----------------------------------------------------------
158      !!
159      IF(lwp) WRITE(numout,*) ' MEDUSA does not support key_trc_diabio'
160# endif
161
162      IF( lk_iomput  .AND.  .NOT.  ln_diatrc  ) THEN
163
164         DO jj = 2,jpjm1
165            DO ji = 2,jpim1
166               !! OPEN wet point IF..THEN loop
167               IF (tmask(ji,jj,jk) == 1) THEN
168                  !!-------------------------------------------------------
169                  !! Add in XML diagnostics stuff
170                  !!-------------------------------------------------------
171                  !!
172                  !! ** 2D diagnostics
173#   if defined key_debug_medusa
174                  IF (lwp) write (numout,*)                                  &
175                     'trc_bio_medusa: diag in ij-jj-jk loop'
176                  CALL flush(numout)
177#   endif
178                  IF ( med_diag%PRN%dgsave ) THEN
179                      fprn2d(ji,jj) = fprn2d(ji,jj) +                        &
180                                      (fprn(ji,jj)  * zphn(ji,jj) *          &
181                                       fse3t(ji,jj,jk)) 
182                  ENDIF
183                  IF ( med_diag%MPN%dgsave ) THEN
184                      fdpn2d(ji,jj) = fdpn2d(ji,jj) + (fdpn(ji,jj) *         &
185                                                       fse3t(ji,jj,jk))
186                  ENDIF
187                  IF ( med_diag%PRD%dgsave ) THEN
188                      fprd2d(ji,jj) = fprd2d(ji,jj) +                        &
189                                      (fprd(ji,jj)  * zphd(ji,jj) *          &
190                                       fse3t(ji,jj,jk))
191                  ENDIF
192                  IF( med_diag%MPD%dgsave ) THEN
193                      fdpd2d(ji,jj) = fdpd2d(ji,jj) + (fdpd(ji,jj) *         &
194                                                       fse3t(ji,jj,jk)) 
195                  ENDIF
196                  !  IF( med_diag%DSED%dgsave ) THEN
197                  !      CALL iom_put( "DSED"  , ftot_n )
198                  !  ENDIF
199                  IF( med_diag%OPAL%dgsave ) THEN
200                      fprds2d(ji,jj) = fprds2d(ji,jj) +                      &
201                                       (fprds(ji,jj) * zpds(ji,jj) *         &
202                                        fse3t(ji,jj,jk)) 
203                  ENDIF
204               ENDIF
205            ENDDO
206         ENDDO
207
208         DO jj = 2,jpjm1
209            DO ji = 2,jpim1
210               IF (tmask(ji,jj,jk) == 1) THEN
211                  IF( med_diag%OPALDISS%dgsave ) THEN
212                      fsdiss2d(ji,jj) = fsdiss2d(ji,jj) + (fsdiss(ji,jj) *   &
213                                                           fse3t(ji,jj,jk)) 
214                  ENDIF
215                  IF( med_diag%GMIPn%dgsave ) THEN
216                      fgmipn2d(ji,jj) = fgmipn2d(ji,jj) +                    &
217                                        (fgmipn(ji,jj)  * fse3t(ji,jj,jk)) 
218                  ENDIF
219                  IF( med_diag%GMID%dgsave ) THEN
220                      fgmid2d(ji,jj) = fgmid2d(ji,jj) + (fgmid(ji,jj) *      &
221                                                         fse3t(ji,jj,jk)) 
222                  ENDIF
223                  IF( med_diag%MZMI%dgsave ) THEN
224                      fdzmi2d(ji,jj) = fdzmi2d(ji,jj) + (fdzmi(ji,jj) *      &
225                                                         fse3t(ji,jj,jk)) 
226                  ENDIF
227               ENDIF
228            ENDDO
229         ENDDO
230
231         DO jj = 2,jpjm1
232            DO ji = 2,jpim1
233               IF (tmask(ji,jj,jk) == 1) THEN
234                  IF( med_diag%GMEPN%dgsave ) THEN
235                      fgmepn2d(ji,jj) = fgmepn2d(ji,jj) + (fgmepn(ji,jj) *   &
236                                                           fse3t(ji,jj,jk))
237                  ENDIF
238                  IF( med_diag%GMEPD%dgsave ) THEN
239                      fgmepd2d(ji,jj) = fgmepd2d(ji,jj) + (fgmepd(ji,jj) *   &
240                                                           fse3t(ji,jj,jk)) 
241                  ENDIF
242                  IF( med_diag%GMEZMI%dgsave ) THEN
243                      fgmezmi2d(ji,jj) = fgmezmi2d(ji,jj) +                  &
244                                         (fgmezmi(ji,jj) * fse3t(ji,jj,jk)) 
245                  ENDIF
246                  IF( med_diag%GMED%dgsave ) THEN
247                      fgmed2d(ji,jj) = fgmed2d(ji,jj) +                      &
248                                       (fgmed(ji,jj) * fse3t(ji,jj,jk)) 
249                  ENDIF
250                  IF( med_diag%MZME%dgsave ) THEN
251                      fdzme2d(ji,jj) = fdzme2d(ji,jj) +                      &
252                                       (fdzme(ji,jj) * fse3t(ji,jj,jk)) 
253                  ENDIF
254                  !  IF( med_diag%DEXP%dgsave ) THEN
255                  !      CALL iom_put( "DEXP"  , ftot_n )
256                  !  ENDIF
257               ENDIF
258            ENDDO
259         ENDDO
260
261         DO jj = 2,jpjm1
262            DO ji = 2,jpim1
263               IF (tmask(ji,jj,jk) == 1) THEN
264                  IF( med_diag%DETN%dgsave ) THEN
265                      fslown2d(ji,jj) = fslown2d(ji,jj) +                    &
266                                        (fslown(ji,jj) * fse3t(ji,jj,jk)) 
267                  ENDIF
268                  IF( med_diag%MDET%dgsave ) THEN
269                      fdd2d(ji,jj) = fdd2d(ji,jj) +                          &
270                                     (fdd(ji,jj) * fse3t(ji,jj,jk)) 
271                  ENDIF
272               ENDIF
273            ENDDO
274         ENDDO
275
276         DO jj = 2,jpjm1
277            DO ji = 2,jpim1
278               IF (tmask(ji,jj,jk) == 1) THEN
279                  IF( med_diag%AEOLIAN%dgsave ) THEN
280                      ffetop2d(ji,jj) = ffetop2d(ji,jj) +                    &
281                                        (ffetop(ji,jj) * fse3t(ji,jj,jk)) 
282                  ENDIF
283                  IF( med_diag%BENTHIC%dgsave ) THEN
284                      ffebot2d(ji,jj) = ffebot2d(ji,jj) +                    &
285                                        (ffebot(ji,jj) * fse3t(ji,jj,jk)) 
286                  ENDIF
287                  IF( med_diag%SCAVENGE%dgsave ) THEN
288                      ffescav2d(ji,jj) = ffescav2d(ji,jj) +                  &
289                                         (ffescav(ji,jj) * fse3t(ji,jj,jk)) 
290                  ENDIF
291               ENDIF
292            ENDDO
293         ENDDO
294
295         DO jj = 2,jpjm1
296            DO ji = 2,jpim1
297               IF (tmask(ji,jj,jk) == 1) THEN
298                  IF( med_diag%PN_JLIM%dgsave ) THEN
299                      ! fjln2d(ji,jj) = fjln2d(ji,jj) +                       &
300                      !                 (fjln(ji,jj)  * zphn(ji,jj) *         &
301                      !                  fse3t(ji,jj,jk))
302                      fjln2d(ji,jj) = fjln2d(ji,jj) +                        &
303                                      (fjlim_pn(ji,jj) * zphn(ji,jj) *       &
304                                       fse3t(ji,jj,jk)) 
305                  ENDIF
306                  IF( med_diag%PN_NLIM%dgsave ) THEN
307                      fnln2d(ji,jj) = fnln2d(ji,jj) +                        &
308                                      (fnln(ji,jj) * zphn(ji,jj) *           &
309                                       fse3t(ji,jj,jk)) 
310                  ENDIF
311                  IF( med_diag%PN_FELIM%dgsave ) THEN
312                      ffln2d(ji,jj) = ffln2d(ji,jj) +                        &
313                                      (ffln2(ji,jj) * zphn(ji,jj) *          &
314                                       fse3t(ji,jj,jk)) 
315                  ENDIF
316               ENDIF
317            ENDDO
318         ENDDO
319
320         DO jj = 2,jpjm1
321            DO ji = 2,jpim1
322               IF (tmask(ji,jj,jk) == 1) THEN
323                  IF( med_diag%PD_JLIM%dgsave ) THEN
324                      ! fjld2d(ji,jj) = fjld2d(ji,jj) +                       &
325                      !                 (fjld(ji,jj)  * zphd(ji,jj) *         &
326                      !                  fse3t(ji,jj,jk))
327                      fjld2d(ji,jj) = fjld2d(ji,jj) +                        &
328                                      (fjlim_pd(ji,jj) * zphd(ji,jj) *       &
329                                       fse3t(ji,jj,jk)) 
330                  ENDIF
331                  IF( med_diag%PD_NLIM%dgsave ) THEN
332                      fnld2d(ji,jj) = fnld2d(ji,jj) +                        &
333                                      (fnld(ji,jj) * zphd(ji,jj) *           &
334                                       fse3t(ji,jj,jk)) 
335                  ENDIF
336                  IF( med_diag%PD_FELIM%dgsave ) THEN
337                      ffld2d(ji,jj) = ffld2d(ji,jj) +                        &
338                                      (ffld(ji,jj) * zphd(ji,jj) *           &
339                                       fse3t(ji,jj,jk)) 
340                  ENDIF
341                  IF( med_diag%PD_SILIM%dgsave ) THEN
342                      fsld2d2(ji,jj) = fsld2d2(ji,jj) +                      &
343                                       (fsld2(ji,jj) * zphd(ji,jj) *         &
344                                        fse3t(ji,jj,jk)) 
345                  ENDIF
346                  IF( med_diag%PDSILIM2%dgsave ) THEN
347                      fsld2d(ji,jj) = fsld2d(ji,jj) +                        &
348                                      (fsld(ji,jj) * zphd(ji,jj) *           &
349                                       fse3t(ji,jj,jk))
350                  ENDIF
351               ENDIF
352            ENDDO
353         ENDDO
354
355         DO jj = 2,jpjm1
356            DO ji = 2,jpim1
357               IF (tmask(ji,jj,jk) == 1) THEN
358                  !!
359                  IF( med_diag%TOTREG_N%dgsave ) THEN
360                      fregen2d(ji,jj) = fregen2d(ji,jj) + fregen(ji,jj)
361                  ENDIF
362                  IF( med_diag%TOTRG_SI%dgsave ) THEN
363                      fregensi2d(ji,jj) = fregensi2d(ji,jj) + fregensi(ji,jj)
364                  ENDIF
365               ENDIF
366            ENDDO
367         ENDDO
368
369         DO jj = 2,jpjm1
370            DO ji = 2,jpim1
371               IF (tmask(ji,jj,jk) == 1) THEN
372                  !!
373                  IF( med_diag%FASTN%dgsave ) THEN
374                      ftempn2d(ji,jj) = ftempn2d(ji,jj) +                    &
375                                        (ftempn(ji,jj)  * fse3t(ji,jj,jk))
376                  ENDIF
377                  IF( med_diag%FASTSI%dgsave ) THEN
378                      ftempsi2d(ji,jj) = ftempsi2d(ji,jj) +                  &
379                                         (ftempsi(ji,jj) * fse3t(ji,jj,jk))
380                  ENDIF
381                  IF( med_diag%FASTFE%dgsave ) THEN
382                      ftempfe2d(ji,jj) = ftempfe2d(ji,jj) +                  &
383                                         (ftempfe(ji,jj) * fse3t(ji,jj,jk)) 
384                  ENDIF
385                  IF( med_diag%FASTC%dgsave ) THEN
386                      ftempc2d(ji,jj) = ftempc2d(ji,jj) +                    &
387                                        (ftempc(ji,jj) * fse3t(ji,jj,jk))
388                  ENDIF
389                  IF( med_diag%FASTCA%dgsave ) THEN
390                      ftempca2d(ji,jj) = ftempca2d(ji,jj) +                  &
391                                         (ftempca(ji,jj) * fse3t(ji,jj,jk))
392                  ENDIF
393               ENDIF
394            ENDDO
395         ENDDO
396
397         DO jj = 2,jpjm1
398            DO ji = 2,jpim1
399               IF (tmask(ji,jj,jk) == 1) THEN
400                  !!
401                  IF( med_diag%REMINN%dgsave ) THEN
402                      freminn2d(ji,jj) = freminn2d(ji,jj) +                  &
403                                         (freminn(ji,jj)  * fse3t(ji,jj,jk))
404                  ENDIF
405                  IF( med_diag%REMINSI%dgsave ) THEN
406                      freminsi2d(ji,jj) = freminsi2d(ji,jj) +                &
407                                          (freminsi(ji,jj) * fse3t(ji,jj,jk))
408                  ENDIF
409                  IF( med_diag%REMINFE%dgsave ) THEN
410                      freminfe2d(ji,jj) = freminfe2d(ji,jj) +                &
411                                          (freminfe(ji,jj) * fse3t(ji,jj,jk)) 
412                  ENDIF
413                  IF( med_diag%REMINC%dgsave ) THEN
414                      freminc2d(ji,jj) = freminc2d(ji,jj) +                  &
415                                         (freminc(ji,jj)  * fse3t(ji,jj,jk)) 
416                  ENDIF
417                  IF( med_diag%REMINCA%dgsave ) THEN
418                      freminca2d(ji,jj) = freminca2d(ji,jj) +                &
419                                          (freminca(ji,jj) * fse3t(ji,jj,jk)) 
420                  ENDIF
421                  !!
422               ENDIF
423            ENDDO
424         ENDDO
425
426# if defined key_roam
427         DO jj = 2,jpjm1
428            DO ji = 2,jpim1
429               IF (tmask(ji,jj,jk) == 1) THEN
430                  !!
431                  !! AXY (09/11/16): CMIP6 diagnostics
432                  IF( med_diag%FD_NIT3%dgsave ) THEN
433                     fd_nit3(ji,jj,jk) = ffastn(ji,jj)
434                  ENDIF
435                  IF( med_diag%FD_SIL3%dgsave ) THEN
436                     fd_sil3(ji,jj,jk) = ffastsi(ji,jj)
437                  ENDIF
438                  IF( med_diag%FD_CAR3%dgsave ) THEN
439                     fd_car3(ji,jj,jk) = ffastc(ji,jj)
440                  ENDIF
441                  IF( med_diag%FD_CAL3%dgsave ) THEN
442                     fd_cal3(ji,jj,jk) = ffastca(ji,jj)
443                  ENDIF
444               ENDIF
445            ENDDO
446         ENDDO
447
448         IF (jk.eq.i0100) THEN
449            DO jj = 2,jpjm1
450               DO ji = 2,jpim1
451                  IF (tmask(ji,jj,jk) == 1) THEN
452                     IF( med_diag%RR_0100%dgsave ) THEN
453                        ffastca2d(ji,jj) =                                   &
454                           ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
455                     ENDIF
456                  ENDIF
457               ENDDO
458            ENDDO
459         ELSE IF (jk.eq.i0500) THEN
460            DO jj = 2,jpjm1
461               DO ji = 2,jpim1
462                  IF (tmask(ji,jj,jk) == 1) THEN
463                     IF( med_diag%RR_0500%dgsave ) THEN
464                        ffastca2d(ji,jj) =                                   &
465                           ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
466                     ENDIF
467                  ENDIF
468               ENDDO
469            ENDDO
470         ELSE IF (jk.eq.i1000) THEN
471            DO jj = 2,jpjm1
472               DO ji = 2,jpim1
473                  IF (tmask(ji,jj,jk) == 1) THEN
474                     IF( med_diag%RR_1000%dgsave ) THEN
475                        ffastca2d(ji,jj) =                                   &
476                           ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
477                     ENDIF
478                  ENDIF
479               ENDDO
480            ENDDO
481         ELSE IF (jk.eq.mbathy(ji,jj)) THEN
482            DO jj = 2,jpjm1
483               DO ji = 2,jpim1
484                  IF (tmask(ji,jj,jk) == 1) THEN
485                     IF( med_diag%IBEN_N%dgsave ) THEN
486                        iben_n2d(ji,jj) = f_sbenin_n(ji,jj) +                &
487                                          f_fbenin_n(ji,jj)
488                     ENDIF
489                     IF( med_diag%IBEN_FE%dgsave ) THEN
490                        iben_fe2d(ji,jj) = f_sbenin_fe(ji,jj) +              &
491                                           f_fbenin_fe(ji,jj)
492                     ENDIF
493                     IF( med_diag%IBEN_C%dgsave ) THEN
494                        iben_c2d(ji,jj) = f_sbenin_c(ji,jj) +                &
495                                          f_fbenin_c(ji,jj)
496                     ENDIF
497                     IF( med_diag%IBEN_SI%dgsave ) THEN
498                        iben_si2d(ji,jj) = f_fbenin_si(ji,jj)
499                     ENDIF
500                     IF( med_diag%IBEN_CA%dgsave ) THEN
501                        iben_ca2d(ji,jj) = f_fbenin_ca(ji,jj)
502                     ENDIF
503                     IF( med_diag%OBEN_N%dgsave ) THEN
504                        oben_n2d(ji,jj) = f_benout_n(ji,jj)
505                     ENDIF
506                     IF( med_diag%OBEN_FE%dgsave ) THEN
507                        oben_fe2d(ji,jj) = f_benout_fe(ji,jj)
508                     ENDIF
509                     IF( med_diag%OBEN_C%dgsave ) THEN
510                        oben_c2d(ji,jj) = f_benout_c(ji,jj)
511                     ENDIF
512                     IF( med_diag%OBEN_SI%dgsave ) THEN
513                        oben_si2d(ji,jj) = f_benout_si(ji,jj)
514                     ENDIF
515                     IF( med_diag%OBEN_CA%dgsave ) THEN
516                        oben_ca2d(ji,jj) = f_benout_ca(ji,jj)
517                     ENDIF
518                     IF( med_diag%SFR_OCAL%dgsave ) THEN
519                        sfr_ocal2d(ji,jj) = f3_omcal(ji,jj,jk)
520                     ENDIF
521                     IF( med_diag%SFR_OARG%dgsave ) THEN
522                        sfr_oarg2d(ji,jj) =  f3_omarg(ji,jj,jk)
523                     ENDIF
524                     IF( med_diag%LYSO_CA%dgsave ) THEN
525                        lyso_ca2d(ji,jj) = f_benout_lyso_ca(ji,jj)
526                     ENDIF
527                  ENDIF
528               ENDDO
529            ENDDO
530         ENDIF
531         !! end bathy-1 diags
532
533         DO jj = 2,jpjm1
534            DO ji = 2,jpim1
535               IF (tmask(ji,jj,jk) == 1) THEN
536                  !!
537                  IF( med_diag%RIV_N%dgsave ) THEN
538                     rivn2d(ji,jj) = rivn2d(ji,jj) +                         &
539                                     (f_riv_loc_n(ji,jj) * fse3t(ji,jj,jk))
540                  ENDIF
541                  IF( med_diag%RIV_SI%dgsave ) THEN
542                     rivsi2d(ji,jj) = rivsi2d(ji,jj) +                       &
543                                      (f_riv_loc_si(ji,jj) * fse3t(ji,jj,jk))
544                  ENDIF
545                  IF( med_diag%RIV_C%dgsave ) THEN
546                     rivc2d(ji,jj) = rivc2d(ji,jj) +                         &
547                                     (f_riv_loc_c(ji,jj) * fse3t(ji,jj,jk))
548                  ENDIF
549                  IF( med_diag%RIV_ALK%dgsave ) THEN
550                     rivalk2d(ji,jj) = rivalk2d(ji,jj) +                     &
551                                       (f_riv_loc_alk(ji,jj) *               &
552                                        fse3t(ji,jj,jk))
553                  ENDIF
554                  IF( med_diag%DETC%dgsave ) THEN
555                     fslowc2d(ji,jj) = fslowc2d(ji,jj) +                     &
556                                       (fslowc(ji,jj)  * fse3t(ji,jj,jk))   
557                  ENDIF
558               ENDIF
559            ENDDO
560         ENDDO
561
562         DO jj = 2,jpjm1
563            DO ji = 2,jpim1
564               IF (tmask(ji,jj,jk) == 1) THEN
565                  !!
566                  IF( med_diag%PN_LLOSS%dgsave ) THEN
567                     fdpn22d(ji,jj) = fdpn22d(ji,jj) +                       &
568                                      (fdpn2(ji,jj)  * fse3t(ji,jj,jk))
569                  ENDIF
570                  IF( med_diag%PD_LLOSS%dgsave ) THEN
571                     fdpd22d(ji,jj) = fdpd22d(ji,jj) +                       &
572                                      (fdpd2(ji,jj)  * fse3t(ji,jj,jk))
573                  ENDIF
574               ENDIF
575            ENDDO
576         ENDDO
577
578         DO jj = 2,jpjm1
579            DO ji = 2,jpim1
580               IF (tmask(ji,jj,jk) == 1) THEN
581                  IF( med_diag%ZI_LLOSS%dgsave ) THEN
582                     fdzmi22d(ji,jj) = fdzmi22d(ji,jj) +                     &
583                                       (fdzmi2(ji,jj) * fse3t(ji,jj,jk))
584                  ENDIF
585                  IF( med_diag%ZE_LLOSS%dgsave ) THEN
586                     fdzme22d(ji,jj) = fdzme22d(ji,jj) +                     &
587                                       (fdzme2(ji,jj) * fse3t(ji,jj,jk))
588                  ENDIF
589               ENDIF
590            ENDDO
591         ENDDO
592
593         DO jj = 2,jpjm1
594            DO ji = 2,jpim1
595               IF (tmask(ji,jj,jk) == 1) THEN
596                  IF( med_diag%ZI_MES_N%dgsave ) THEN
597                     zimesn2d(ji,jj) = zimesn2d(ji,jj) +                     &
598                                       (xphi * (fgmipn(ji,jj) +              &
599                                                fgmid(ji,jj)) *              &
600                                        fse3t(ji,jj,jk))
601                  ENDIF
602                  IF( med_diag%ZI_MES_D%dgsave ) THEN
603                     zimesd2d(ji,jj) = zimesd2d(ji,jj) +                     & 
604                                       ((1. - xbetan) * finmi(ji,jj) *       &
605                                        fse3t(ji,jj,jk))
606                  ENDIF
607                  IF( med_diag%ZI_MES_C%dgsave ) THEN
608                     zimesc2d(ji,jj) = zimesc2d(ji,jj) +                     &
609                                       (xphi * ((xthetapn * fgmipn(ji,jj)) + &
610                                                fgmidc(ji,jj)) *             &
611                                                fse3t(ji,jj,jk))
612                  ENDIF
613                  IF( med_diag%ZI_MESDC%dgsave ) THEN
614                     zimesdc2d(ji,jj) = zimesdc2d(ji,jj) +                   &
615                                        ((1. - xbetac) * ficmi(ji,jj) *      &
616                                         fse3t(ji,jj,jk))
617                  ENDIF
618                  IF( med_diag%ZI_EXCR%dgsave ) THEN
619                     ziexcr2d(ji,jj) = ziexcr2d(ji,jj) +                     &
620                                       (fmiexcr(ji,jj) * fse3t(ji,jj,jk))
621                  ENDIF
622                  IF( med_diag%ZI_RESP%dgsave ) THEN
623                     ziresp2d(ji,jj) = ziresp2d(ji,jj) +                     &
624                                       (fmiresp(ji,jj) * fse3t(ji,jj,jk))
625                  ENDIF
626                  IF( med_diag%ZI_GROW%dgsave ) THEN
627                     zigrow2d(ji,jj) = zigrow2d(ji,jj) +                     &
628                                       (fmigrow(ji,jj) * fse3t(ji,jj,jk))
629                  ENDIF
630               ENDIF
631            ENDDO
632         ENDDO
633
634         DO jj = 2,jpjm1
635            DO ji = 2,jpim1
636               IF (tmask(ji,jj,jk) == 1) THEN
637                  IF( med_diag%ZE_MES_N%dgsave ) THEN
638                     zemesn2d(ji,jj) = zemesn2d(ji,jj) +                     &
639                                       (xphi *                               &
640                                        (fgmepn(ji,jj) + fgmepd(ji,jj) +     &
641                                         fgmezmi(ji,jj) + fgmed(ji,jj)) *    &
642                                        fse3t(ji,jj,jk))
643                  ENDIF
644                  IF( med_diag%ZE_MES_D%dgsave ) THEN
645                     zemesd2d(ji,jj) = zemesd2d(ji,jj) +                     &
646                                       ((1. - xbetan) * finme(ji,jj) *       &
647                                        fse3t(ji,jj,jk))
648                  ENDIF
649                  IF( med_diag%ZE_MES_C%dgsave ) THEN
650                     zemesc2d(ji,jj) = zemesc2d(ji,jj) +                     & 
651                                       (xphi *                               &
652                                        ((xthetapn * fgmepn(ji,jj)) +        &
653                                         (xthetapd * fgmepd(ji,jj)) +        &
654                                         (xthetazmi * fgmezmi(ji,jj)) +      &
655                                         fgmedc(ji,jj)) * fse3t(ji,jj,jk))
656                  ENDIF
657                  IF( med_diag%ZE_MESDC%dgsave ) THEN
658                     zemesdc2d(ji,jj) = zemesdc2d(ji,jj) +                   &
659                                        ((1. - xbetac) * ficme(ji,jj) *      &
660                                         fse3t(ji,jj,jk))
661                  ENDIF
662                  IF( med_diag%ZE_EXCR%dgsave ) THEN
663                     zeexcr2d(ji,jj) = zeexcr2d(ji,jj) +                     &
664                                       (fmeexcr(ji,jj) * fse3t(ji,jj,jk))
665                  ENDIF
666                  IF( med_diag%ZE_RESP%dgsave ) THEN
667                     zeresp2d(ji,jj) = zeresp2d(ji,jj) +                     &
668                                       (fmeresp(ji,jj) * fse3t(ji,jj,jk))
669                  ENDIF
670                  IF( med_diag%ZE_GROW%dgsave ) THEN
671                     zegrow2d(ji,jj) = zegrow2d(ji,jj) +                     &
672                                       (fmegrow(ji,jj) * fse3t(ji,jj,jk))
673                  ENDIF
674               ENDIF
675            ENDDO
676         ENDDO
677
678         DO jj = 2,jpjm1
679            DO ji = 2,jpim1
680               IF (tmask(ji,jj,jk) == 1) THEN
681                  IF( med_diag%MDETC%dgsave ) THEN
682                     mdetc2d(ji,jj) = mdetc2d(ji,jj) +                       &
683                                      (fddc(ji,jj) * fse3t(ji,jj,jk))
684                  ENDIF
685                  IF( med_diag%GMIDC%dgsave ) THEN
686                     gmidc2d(ji,jj) = gmidc2d(ji,jj) +                       &
687                                      (fgmidc(ji,jj) * fse3t(ji,jj,jk))
688                  ENDIF
689                  IF( med_diag%GMEDC%dgsave ) THEN
690                     gmedc2d(ji,jj) = gmedc2d(ji,jj) +                       &
691                                      (fgmedc(ji,jj) * fse3t(ji,jj,jk))
692                  ENDIF
693               ENDIF
694            ENDDO
695         ENDDO
696# endif                   
697
698         DO jj = 2,jpjm1
699            DO ji = 2,jpim1
700               IF (tmask(ji,jj,jk) == 1) THEN
701                  !!
702                  !! ** 3D diagnostics
703                  IF( med_diag%TPP3%dgsave ) THEN
704                     tpp3d(ji,jj,jk) = (fprn(ji,jj) * zphn(ji,jj)) +         &
705                                       (fprd(ji,jj) * zphd(ji,jj))
706                     !CALL iom_put( "TPP3"  , tpp3d )
707                  ENDIF
708                  IF( med_diag%TPPD3%dgsave ) THEN
709                     tppd3(ji,jj,jk) = (fprd(ji,jj) * zphd(ji,jj))
710                  ENDIF
711               ENDIF
712            ENDDO
713         ENDDO
714
715         DO jj = 2,jpjm1
716            DO ji = 2,jpim1
717               IF (tmask(ji,jj,jk) == 1) THEN
718                 
719                  IF( med_diag%REMIN3N%dgsave ) THEN
720                     !! remineralisation
721                     remin3dn(ji,jj,jk) = fregen(ji,jj) +                    &
722                                          (freminn(ji,jj) * fse3t(ji,jj,jk))
723                     !CALL iom_put( "REMIN3N"  , remin3dn )
724                  ENDIF
725                  !! IF( med_diag%PH3%dgsave ) THEN
726                  !!     CALL iom_put( "PH3"  , f3_pH )
727                  !! ENDIF
728                  !! IF( med_diag%OM_CAL3%dgsave ) THEN
729                  !!     CALL iom_put( "OM_CAL3"  , f3_omcal )
730                  !! ENDIF
731        !!
732        !! AXY (09/11/16): CMIP6 diagnostics
733        IF ( med_diag%DCALC3%dgsave   ) THEN
734                     dcalc3(ji,jj,jk) = freminca(ji,jj)
735                  ENDIF
736               ENDIF
737            ENDDO
738         ENDDO
739
740         DO jj = 2,jpjm1
741            DO ji = 2,jpim1
742               IF (tmask(ji,jj,jk) == 1) THEN
743        IF ( med_diag%FEDISS3%dgsave  ) THEN
744                     fediss3(ji,jj,jk) = ffetop(ji,jj)
745                  ENDIF
746        IF ( med_diag%FESCAV3%dgsave  ) THEN
747                     fescav3(ji,jj,jk) = ffescav(ji,jj)
748                  ENDIF
749               ENDIF
750            ENDDO
751         ENDDO
752
753         DO jj = 2,jpjm1
754            DO ji = 2,jpim1
755               IF (tmask(ji,jj,jk) == 1) THEN
756        IF ( med_diag%MIGRAZP3%dgsave ) THEN
757                     migrazp3(ji,jj,jk) = fgmipn(ji,jj) * xthetapn
758                  ENDIF
759        IF ( med_diag%MIGRAZD3%dgsave ) THEN
760                     migrazd3(ji,jj,jk) = fgmidc(ji,jj)
761                  ENDIF
762        IF ( med_diag%MEGRAZP3%dgsave ) THEN
763                     megrazp3(ji,jj,jk) = (fgmepn(ji,jj) * xthetapn) +       &
764                                          (fgmepd(ji,jj) * xthetapd)
765                  ENDIF
766        IF ( med_diag%MEGRAZD3%dgsave ) THEN
767                     megrazd3(ji,jj,jk) = fgmedc(ji,jj)
768                  ENDIF
769        IF ( med_diag%MEGRAZZ3%dgsave ) THEN
770                     megrazz3(ji,jj,jk) = (fgmezmi(ji,jj) * xthetazmi)
771                  ENDIF
772               ENDIF
773            ENDDO
774         ENDDO
775
776         DO jj = 2,jpjm1
777            DO ji = 2,jpim1
778               IF (tmask(ji,jj,jk) == 1) THEN
779        IF ( med_diag%PBSI3%dgsave    ) THEN
780                     pbsi3(ji,jj,jk)    = (fprds(ji,jj) * zpds(ji,jj))
781                  ENDIF
782        IF ( med_diag%PCAL3%dgsave    ) THEN
783                     pcal3(ji,jj,jk)    = ftempca(ji,jj)
784                  ENDIF
785        IF ( med_diag%REMOC3%dgsave   ) THEN
786                     remoc3(ji,jj,jk)   = freminc(ji,jj)
787                  ENDIF
788               ENDIF
789            ENDDO
790         ENDDO
791
792         DO jj = 2,jpjm1
793            DO ji = 2,jpim1
794               IF (tmask(ji,jj,jk) == 1) THEN
795        IF ( med_diag%PNLIMJ3%dgsave  ) THEN
796                     ! pnlimj3(ji,jj,jk)  = fjln(ji,jj)
797                     pnlimj3(ji,jj,jk)  = fjlim_pn(ji,jj)
798                  ENDIF
799        IF ( med_diag%PNLIMN3%dgsave  ) THEN
800                     pnlimn3(ji,jj,jk)  = fnln(ji,jj)
801                  ENDIF
802        IF ( med_diag%PNLIMFE3%dgsave ) THEN
803                     pnlimfe3(ji,jj,jk) = ffln2(ji,jj)
804                  ENDIF
805        IF ( med_diag%PDLIMJ3%dgsave  ) THEN
806                     ! pdlimj3(ji,jj,jk)  = fjld(ji,jj)
807                     pdlimj3(ji,jj,jk)  = fjlim_pd(ji,jj)
808                  ENDIF
809        IF ( med_diag%PDLIMN3%dgsave  ) THEN
810                     pdlimn3(ji,jj,jk)  = fnld(ji,jj)
811                  ENDIF
812        IF ( med_diag%PDLIMFE3%dgsave ) THEN
813                     pdlimfe3(ji,jj,jk) = ffld(ji,jj)
814                  ENDIF
815        IF ( med_diag%PDLIMSI3%dgsave ) THEN
816                     pdlimsi3(ji,jj,jk) = fsld2(ji,jj)
817                  ENDIF
818               ENDIF
819            ENDDO
820         ENDDO
821
822      ELSE IF( ln_diatrc ) THEN
823
824         !!
825         !! ** Without using iom_use
826#   if defined key_debug_medusa
827         IF (lwp) write (numout,*) 'trc_bio_medusa: diag in ij-jj-jk ln_diatrc'
828         CALL flush(numout)
829#   endif
830         DO jj = 2,jpjm1
831            DO ji = 2,jpim1
832               IF (tmask(ji,jj,jk) == 1) then
833                  !!-------------------------------------------------------
834                  !! Prepare 2D diagnostics
835                  !!-------------------------------------------------------
836                  !!
837                  !! if ((kt / 240*240).eq.kt) then
838                  !!    IF (lwp) write (*,*) '*******!MEDUSA DIAADD!*******',kt
839                  !! endif     
840                  !! nitrogen inventory
841                  trc2d(ji,jj,1)  =  ftot_n(ji,jj)
842                  !! silicon  inventory
843                  trc2d(ji,jj,2)  =  ftot_si(ji,jj)
844                  !! iron     inventory
845                  trc2d(ji,jj,3)  =  ftot_fe(ji,jj)
846               ENDIF
847            ENDDO
848         ENDDO
849
850         DO jj = 2,jpjm1
851            DO ji = 2,jpim1
852               IF (tmask(ji,jj,jk) == 1) THEN
853                  !! non-diatom production
854                  trc2d(ji,jj,4)  = trc2d(ji,jj,4)  +                        &
855                                    (fprn(ji,jj)  * zphn(ji,jj) *            &
856                                     fse3t(ji,jj,jk))
857                  !! non-diatom non-grazing losses
858                  trc2d(ji,jj,5)  = trc2d(ji,jj,5)  +                        &
859                                    (fdpn(ji,jj) * fse3t(ji,jj,jk))
860                  !! diatom production
861                  trc2d(ji,jj,6)  = trc2d(ji,jj,6)  +                        &
862                                    (fprd(ji,jj) * zphd(ji,jj) *             &
863                                     fse3t(ji,jj,jk))
864                  !! diatom non-grazing losses
865                  !! diagnostic field  8 is (ostensibly) supplied by trcsed.F
866                  trc2d(ji,jj,7)  = trc2d(ji,jj,7)  +                        &
867                                    (fdpd(ji,jj) * fse3t(ji,jj,jk))
868                  !! diatom silicon production
869                  trc2d(ji,jj,9)  = trc2d(ji,jj,9)  +                        &
870                                    (fprds(ji,jj) * zpds(ji,jj) *            &
871                                     fse3t(ji,jj,jk))
872                  !! diatom silicon dissolution
873                  trc2d(ji,jj,10) = trc2d(ji,jj,10) +                        &
874                                    (fsdiss(ji,jj)  * fse3t(ji,jj,jk))
875               ENDIF
876            ENDDO
877         ENDDO
878
879         DO jj = 2,jpjm1
880            DO ji = 2,jpim1
881               IF (tmask(ji,jj,jk) == 1) THEN
882                  !! microzoo grazing on non-diatoms
883                  trc2d(ji,jj,11) = trc2d(ji,jj,11) +                        &
884                                    (fgmipn(ji,jj)  * fse3t(ji,jj,jk))
885                  !! microzoo grazing on detrital nitrogen
886                  trc2d(ji,jj,12) = trc2d(ji,jj,12) +                        &
887                                    (fgmid(ji,jj) * fse3t(ji,jj,jk))
888                  !! microzoo non-grazing losses
889                  trc2d(ji,jj,13) = trc2d(ji,jj,13) +                        &
890                                    (fdzmi(ji,jj) * fse3t(ji,jj,jk))
891               ENDIF
892            ENDDO
893         ENDDO
894
895         DO jj = 2,jpjm1
896            DO ji = 2,jpim1
897               IF (tmask(ji,jj,jk) == 1) THEN
898                  !! mesozoo  grazing on non-diatoms
899                  trc2d(ji,jj,14) = trc2d(ji,jj,14) +                        &
900                                    (fgmepn(ji,jj)  * fse3t(ji,jj,jk))
901                  !! mesozoo  grazing on diatoms
902                  trc2d(ji,jj,15) = trc2d(ji,jj,15) +                        &
903                                    (fgmepd(ji,jj)  * fse3t(ji,jj,jk))
904                  !! mesozoo  grazing on microzoo
905                  trc2d(ji,jj,16) = trc2d(ji,jj,16) +                        &
906                                    (fgmezmi(ji,jj) * fse3t(ji,jj,jk))
907                  !! mesozoo  grazing on detrital nitrogen
908                  trc2d(ji,jj,17) = trc2d(ji,jj,17) +                        &
909                                    (fgmed(ji,jj) * fse3t(ji,jj,jk))
910                  !! mesozoo  non-grazing losses
911                  trc2d(ji,jj,18) = trc2d(ji,jj,18) +                        &
912                                    (fdzme(ji,jj)   * fse3t(ji,jj,jk))
913               ENDIF
914            ENDDO
915         ENDDO
916
917         DO jj = 2,jpjm1
918            DO ji = 2,jpim1
919               IF (tmask(ji,jj,jk) == 1) THEN
920                  !! diagnostic field 19 is (ostensibly) supplied by trcexp.F
921                  !! slow sinking detritus N production
922                  trc2d(ji,jj,20) = trc2d(ji,jj,20) +                        &
923                                    (fslown(ji,jj) * fse3t(ji,jj,jk))
924                  !! detrital remineralisation
925                  trc2d(ji,jj,21) = trc2d(ji,jj,21) +                        &
926                                    (fdd(ji,jj) * fse3t(ji,jj,jk))
927                  !! aeolian  iron addition
928                  trc2d(ji,jj,22) = trc2d(ji,jj,22) +                        &
929                                    (ffetop(ji,jj) * fse3t(ji,jj,jk))
930                  !! seafloor iron addition
931                  trc2d(ji,jj,23) = trc2d(ji,jj,23) +                        &
932                                    (ffebot(ji,jj) * fse3t(ji,jj,jk))
933                  !! "free" iron scavenging
934                  trc2d(ji,jj,24) = trc2d(ji,jj,24) +                        &
935                                    (ffescav(ji,jj) * fse3t(ji,jj,jk))
936               ENDIF
937            ENDDO
938         ENDDO
939
940         DO jj = 2,jpjm1
941            DO ji = 2,jpim1
942               IF (tmask(ji,jj,jk) == 1) THEN
943                  !! non-diatom J  limitation term
944                  trc2d(ji,jj,25) = trc2d(ji,jj,25) +                        &
945                                    (fjlim_pn(ji,jj) * zphn(ji,jj) *         &
946                                     fse3t(ji,jj,jk))
947                  !! non-diatom N  limitation term
948                  trc2d(ji,jj,26) = trc2d(ji,jj,26) +                        &
949                                    (fnln(ji,jj) * zphn(ji,jj) *             &
950                                     fse3t(ji,jj,jk))
951                  !! non-diatom Fe limitation term
952                  trc2d(ji,jj,27) = trc2d(ji,jj,27) +                        &
953                                    (ffln2(ji,jj) * zphn(ji,jj) *            &
954                                     fse3t(ji,jj,jk))
955                  !! diatom     J  limitation term
956                  trc2d(ji,jj,28) = trc2d(ji,jj,28) +                        &
957                                    (fjlim_pd(ji,jj) * zphd(ji,jj) *         &
958                                     fse3t(ji,jj,jk))
959                  !! diatom     N  limitation term
960                  trc2d(ji,jj,29) = trc2d(ji,jj,29) +                        &
961                                    (fnld(ji,jj) * zphd(ji,jj) *             &
962                                     fse3t(ji,jj,jk))
963                  !! diatom     Fe limitation term
964                  trc2d(ji,jj,30) = trc2d(ji,jj,30) +                        &
965                                    (ffld(ji,jj) * zphd(ji,jj) *             &
966                                     fse3t(ji,jj,jk))
967                  !! diatom     Si limitation term
968                  trc2d(ji,jj,31) = trc2d(ji,jj,31) +                        &
969                                    (fsld2(ji,jj) * zphd(ji,jj) *            &
970                                     fse3t(ji,jj,jk))
971                  !! diatom     Si uptake limitation term
972                  trc2d(ji,jj,32) = trc2d(ji,jj,32) +                        &
973                                    (fsld(ji,jj) * zphd(ji,jj) *             &
974                                     fse3t(ji,jj,jk))
975               ENDIF
976            ENDDO
977         ENDDO
978
979         IF (jk.eq.i0100) THEN
980            DO jj = 2,jpjm1
981               DO ji = 2,jpim1
982                  IF (tmask(ji,jj,jk) == 1) THEN
983                     !! slow detritus flux at  100 m
984                     trc2d(ji,jj,33) = fslownflux(ji,jj)
985                  ENDIF
986               ENDDO
987            ENDDO
988         ENDIF
989
990         IF (jk.eq.i0200) THEN
991            DO jj = 2,jpjm1
992               DO ji = 2,jpim1
993                  IF (tmask(ji,jj,jk) == 1) THEN
994                     !! slow detritus flux at  200 m
995                     trc2d(ji,jj,34) = fslownflux(ji,jj)
996                  ENDIF
997               ENDDO
998            ENDDO
999         ENDIF
1000
1001         IF (jk.eq.i0500) THEN
1002            DO jj = 2,jpjm1
1003               DO ji = 2,jpim1
1004                  IF (tmask(ji,jj,jk) == 1) THEN
1005                     !! slow detritus flux at  500 m
1006                     trc2d(ji,jj,35) = fslownflux(ji,jj)
1007                  ENDIF
1008               ENDDO
1009            ENDDO
1010         ENDIF
1011
1012         IF (jk.eq.i1000) THEN
1013            DO jj = 2,jpjm1
1014               DO ji = 2,jpim1
1015                  IF (tmask(ji,jj,jk) == 1) THEN
1016                     !! slow detritus flux at 1000 m
1017                     trc2d(ji,jj,36) = fslownflux(ji,jj)
1018                  ENDIF
1019               ENDDO
1020            ENDDO
1021         ENDIF
1022
1023         DO jj = 2,jpjm1
1024            DO ji = 2,jpim1
1025               IF (tmask(ji,jj,jk) == 1) THEN
1026                  !! non-fast N  full column regeneration
1027                  trc2d(ji,jj,37) = trc2d(ji,jj,37) + fregen(ji,jj)
1028                  !! non-fast Si full column regeneration
1029                  trc2d(ji,jj,38) = trc2d(ji,jj,38) + fregensi(ji,jj)
1030                  !! non-fast N  regeneration to  100 m
1031               ENDIF
1032            ENDDO
1033         ENDDO
1034
1035         IF (jk.eq.i0100) THEN
1036            DO jj = 2,jpjm1
1037               DO ji = 2,jpim1
1038                  IF (tmask(ji,jj,jk) == 1) THEN
1039                     trc2d(ji,jj,39) = trc2d(ji,jj,37)
1040                  ENDIF
1041               ENDDO
1042            ENDDO
1043         ENDIF
1044
1045         IF (jk.eq.i0200) THEN
1046            DO jj = 2,jpjm1
1047               DO ji = 2,jpim1
1048                  IF (tmask(ji,jj,jk) == 1) THEN
1049                     !! non-fast N  regeneration to  200 m
1050                     trc2d(ji,jj,40) = trc2d(ji,jj,37)
1051                  ENDIF
1052               ENDDO
1053            ENDDO
1054         ENDIF
1055
1056         IF (jk.eq.i0500) THEN
1057            DO jj = 2,jpjm1
1058               DO ji = 2,jpim1
1059                  IF (tmask(ji,jj,jk) == 1) THEN
1060                     !! non-fast N  regeneration to  500 m
1061                     trc2d(ji,jj,41) = trc2d(ji,jj,37)
1062                  ENDIF
1063               ENDDO
1064            ENDDO
1065         ENDIF
1066
1067         IF (jk.eq.i1000) THEN
1068            DO jj = 2,jpjm1
1069               DO ji = 2,jpim1
1070                  IF (tmask(ji,jj,jk) == 1) THEN
1071                     !! non-fast N  regeneration to 1000 m
1072                     trc2d(ji,jj,42) = trc2d(ji,jj,37)
1073                  ENDIF
1074               ENDDO
1075            ENDDO
1076         ENDIF
1077
1078         DO jj = 2,jpjm1
1079            DO ji = 2,jpim1
1080               IF (tmask(ji,jj,jk) == 1) THEN
1081                  !! fast sinking detritus N production
1082                  trc2d(ji,jj,43) = trc2d(ji,jj,43) +                        &
1083                                    (ftempn(ji,jj) * fse3t(ji,jj,jk))
1084                  !! fast sinking detritus Si production
1085                  trc2d(ji,jj,44) = trc2d(ji,jj,44) +                        &
1086                                    (ftempsi(ji,jj) * fse3t(ji,jj,jk))
1087                  !! fast sinking detritus Fe production
1088                  trc2d(ji,jj,45) = trc2d(ji,jj,45) +                        &
1089                                    (ftempfe(ji,jj) * fse3t(ji,jj,jk))
1090                  !! fast sinking detritus C production
1091                  trc2d(ji,jj,46) = trc2d(ji,jj,46) +                        &
1092                                    (ftempc(ji,jj)  * fse3t(ji,jj,jk))
1093                  !! fast sinking detritus CaCO3 production
1094                  trc2d(ji,jj,47) = trc2d(ji,jj,47) +                        &
1095                                    (ftempca(ji,jj) * fse3t(ji,jj,jk))
1096               ENDIF
1097            ENDDO
1098         ENDDO
1099
1100         IF (jk.eq.i0100) THEN
1101            DO jj = 2,jpjm1
1102               DO ji = 2,jpim1
1103                  IF (tmask(ji,jj,jk) == 1) THEN
1104                     !! fast detritus N  flux at  100 m
1105                     trc2d(ji,jj,48) = ffastn(ji,jj)
1106                  ENDIF
1107               ENDDO
1108            ENDDO
1109         ENDIF
1110
1111         IF (jk.eq.i0200) THEN
1112            DO jj = 2,jpjm1
1113               DO ji = 2,jpim1
1114                  IF (tmask(ji,jj,jk) == 1) THEN
1115                     !! fast detritus N  flux at  200 m
1116                     trc2d(ji,jj,49) = ffastn(ji,jj)
1117                  ENDIF
1118               ENDDO
1119            ENDDO
1120         ENDIF
1121
1122         IF (jk.eq.i0500) THEN
1123            DO jj = 2,jpjm1
1124               DO ji = 2,jpim1
1125                  IF (tmask(ji,jj,jk) == 1) THEN
1126                     !! fast detritus N  flux at  500 m
1127                     trc2d(ji,jj,50) = ffastn(ji,jj)
1128                  ENDIF
1129               ENDDO
1130            ENDDO
1131         ENDIF
1132
1133         IF (jk.eq.i1000) THEN
1134            DO jj = 2,jpjm1
1135               DO ji = 2,jpim1
1136                  IF (tmask(ji,jj,jk) == 1) THEN
1137                     !! fast detritus N  flux at 1000 m
1138                     trc2d(ji,jj,51) = ffastn(ji,jj)
1139                  ENDIF
1140               ENDDO
1141            ENDDO
1142         ENDIF
1143
1144         IF (jk.eq.i0100) THEN
1145            DO jj = 2,jpjm1
1146               DO ji = 2,jpim1
1147                  IF (tmask(ji,jj,jk) == 1) THEN
1148                     !! N  regeneration to  100 m
1149                     trc2d(ji,jj,52) = fregenfast(ji,jj)
1150                  ENDIF
1151               ENDDO
1152            ENDDO
1153         ENDIF
1154
1155         IF (jk.eq.i0200) THEN
1156            DO jj = 2,jpjm1
1157               DO ji = 2,jpim1
1158                  IF (tmask(ji,jj,jk) == 1) THEN
1159                     !! N  regeneration to  200 m
1160                     trc2d(ji,jj,53) = fregenfast(ji,jj)
1161                  ENDIF
1162               ENDDO
1163            ENDDO
1164         ENDIF
1165
1166         IF (jk.eq.i0500) THEN
1167            DO jj = 2,jpjm1
1168               DO ji = 2,jpim1
1169                  IF (tmask(ji,jj,jk) == 1) THEN
1170                     !! N  regeneration to  500 m
1171                     trc2d(ji,jj,54) = fregenfast(ji,jj)
1172                  ENDIF
1173               ENDDO
1174            ENDDO
1175         ENDIF
1176
1177         IF (jk.eq.i1000) THEN
1178            DO jj = 2,jpjm1
1179               DO ji = 2,jpim1
1180                  IF (tmask(ji,jj,jk) == 1) THEN
1181                     !! N  regeneration to 1000 m
1182                     trc2d(ji,jj,55) = fregenfast(ji,jj)
1183                  ENDIF
1184               ENDDO
1185            ENDDO
1186         ENDIF
1187
1188         IF (jk.eq.i0100) THEN
1189            DO jj = 2,jpjm1
1190               DO ji = 2,jpim1
1191                  IF (tmask(ji,jj,jk) == 1) THEN
1192                     !! fast detritus Si flux at  100 m
1193                     trc2d(ji,jj,56) = ffastsi(ji,jj)
1194                  ENDIF
1195               ENDDO
1196            ENDDO
1197         ENDIF
1198
1199         IF (jk.eq.i0200) THEN
1200            DO jj = 2,jpjm1
1201               DO ji = 2,jpim1
1202                  IF (tmask(ji,jj,jk) == 1) THEN
1203                     !! fast detritus Si flux at  200 m
1204                     trc2d(ji,jj,57) = ffastsi(ji,jj)
1205                  ENDIF
1206               ENDDO
1207            ENDDO
1208         ENDIF
1209
1210         IF (jk.eq.i0500) THEN
1211            DO jj = 2,jpjm1
1212               DO ji = 2,jpim1
1213                  IF (tmask(ji,jj,jk) == 1) THEN
1214                     !! fast detritus Si flux at  500 m
1215                     trc2d(ji,jj,58) = ffastsi(ji,jj)
1216                  ENDIF
1217               ENDDO
1218            ENDDO
1219         ENDIF
1220
1221         IF (jk.eq.i1000) THEN
1222            DO jj = 2,jpjm1
1223               DO ji = 2,jpim1
1224                  IF (tmask(ji,jj,jk) == 1) THEN
1225                     !! fast detritus Si flux at 1000 m
1226                     trc2d(ji,jj,59) = ffastsi(ji,jj)
1227                  ENDIF
1228               ENDDO
1229            ENDDO
1230         ENDIF
1231
1232         IF (jk.eq.i0100) THEN
1233            DO jj = 2,jpjm1
1234               DO ji = 2,jpim1
1235                  IF (tmask(ji,jj,jk) == 1) THEN
1236                     !! Si regeneration to  100 m
1237                     trc2d(ji,jj,60) = fregenfastsi(ji,jj)
1238                  ENDIF
1239               ENDDO
1240            ENDDO
1241         ENDIF
1242
1243         IF (jk.eq.i0200) THEN
1244            DO jj = 2,jpjm1
1245               DO ji = 2,jpim1
1246                  IF (tmask(ji,jj,jk) == 1) THEN
1247                     !! Si regeneration to  200 m
1248                     trc2d(ji,jj,61) = fregenfastsi(ji,jj)
1249                  ENDIF
1250               ENDDO
1251            ENDDO
1252         ENDIF
1253
1254         IF (jk.eq.i0500) THEN
1255            DO jj = 2,jpjm1
1256               DO ji = 2,jpim1
1257                  IF (tmask(ji,jj,jk) == 1) THEN
1258                     !! Si regeneration to  500 m
1259                     trc2d(ji,jj,62) = fregenfastsi(ji,jj)
1260                  ENDIF
1261               ENDDO
1262            ENDDO
1263         ENDIF
1264
1265         IF (jk.eq.i1000) THEN
1266            DO jj = 2,jpjm1
1267               DO ji = 2,jpim1
1268                  IF (tmask(ji,jj,jk) == 1) THEN
1269                     !! Si regeneration to 1000 m
1270                     trc2d(ji,jj,63) = fregenfastsi(ji,jj)
1271                  ENDIF
1272               ENDDO
1273            ENDDO
1274         ENDIF
1275
1276         DO jj = 2,jpjm1
1277            DO ji = 2,jpim1
1278               IF (tmask(ji,jj,jk) == 1) THEN
1279                  !! sum of fast-sinking N  fluxes
1280                  trc2d(ji,jj,64) = trc2d(ji,jj,64) +                        &
1281                                    (freminn(ji,jj) * fse3t(ji,jj,jk))
1282                  !! sum of fast-sinking Si fluxes
1283                  trc2d(ji,jj,65) = trc2d(ji,jj,65) +                        &
1284                                    (freminsi(ji,jj) * fse3t(ji,jj,jk))
1285                  !! sum of fast-sinking Fe fluxes
1286                  trc2d(ji,jj,66) = trc2d(ji,jj,66) +                        &
1287                                    (freminfe(ji,jj) * fse3t(ji,jj,jk))
1288                  !! sum of fast-sinking C  fluxes
1289                  trc2d(ji,jj,67) = trc2d(ji,jj,67) +                        &
1290                                    (freminc(ji,jj) * fse3t(ji,jj,jk))
1291                  !! sum of fast-sinking Ca fluxes
1292                  trc2d(ji,jj,68) = trc2d(ji,jj,68) +                        &
1293                                    (freminca(ji,jj) * fse3t(ji,jj,jk))
1294               ENDIF
1295            ENDDO
1296         ENDDO
1297
1298
1299         if (jk.eq.mbathy(ji,jj)) then
1300            DO jj = 2,jpjm1
1301               DO ji = 2,jpim1
1302                  IF (tmask(ji,jj,jk) == 1) THEN
1303                     !! N  sedimentation flux
1304                     trc2d(ji,jj,69) = fsedn(ji,jj)
1305                     !! Si sedimentation flux
1306                     trc2d(ji,jj,70) = fsedsi(ji,jj)
1307                     !! Fe sedimentation flux
1308                     trc2d(ji,jj,71) = fsedfe(ji,jj)
1309                     !! C  sedimentation flux
1310                     trc2d(ji,jj,72) = fsedc(ji,jj)
1311                     !! Ca sedimentation flux
1312                     trc2d(ji,jj,73) = fsedca(ji,jj)
1313                  ENDIF
1314               ENDDO
1315            ENDDO
1316         endif
1317
1318         if (jk.eq.1) then
1319            DO jj = 2,jpjm1
1320               DO ji = 2,jpim1
1321                  IF (tmask(ji,jj,jk) == 1) THEN
1322                     trc2d(ji,jj,74) = qsr(ji,jj)
1323                     trc2d(ji,jj,75) = xpar(ji,jj,jk)
1324                     !! trc2d(ji,jj,75) = real(iters(ji,jj))
1325                  ENDIF
1326               ENDDO
1327            ENDDO
1328         endif
1329
1330         DO jj = 2,jpjm1
1331            DO ji = 2,jpim1
1332               IF (tmask(ji,jj,jk) == 1) THEN
1333                  !! diagnostic fields 76 to 80 calculated below
1334                  !! mixed layer non-diatom production
1335                  trc2d(ji,jj,81) = trc2d(ji,jj,81) + fprn_ml(ji,jj)
1336                  !! mixed layer     diatom production
1337                  trc2d(ji,jj,82) = trc2d(ji,jj,82) + fprd_ml(ji,jj)
1338               ENDIF
1339            ENDDO
1340         ENDDO
1341
1342# if defined key_gulf_finland
1343         if (jk.eq.1) then
1344            DO jj = 2,jpjm1
1345               DO ji = 2,jpim1
1346                  IF (tmask(ji,jj,jk) == 1) THEN
1347                     !! Gulf of Finland check
1348                     trc2d(ji,jj,83) = real(ibio_switch)
1349                  ENDIF
1350               ENDDO
1351            ENDDO
1352         endif
1353# else
1354         DO jj = 2,jpjm1
1355            DO ji = 2,jpim1
1356               IF (tmask(ji,jj,jk) == 1) THEN
1357                  !! calcite CCD depth
1358                  trc2d(ji,jj,83) = ocal_ccd(ji,jj)
1359               ENDIF
1360            ENDDO
1361         ENDDO
1362# endif
1363         DO jj = 2,jpjm1
1364            DO ji = 2,jpim1
1365               IF (tmask(ji,jj,jk) == 1) THEN
1366                  !! last model level above calcite CCD depth
1367                  trc2d(ji,jj,84) = fccd(ji,jj)
1368               ENDIF
1369            ENDDO
1370         ENDDO
1371
1372         IF (jk.eq.1) THEN
1373            DO jj = 2,jpjm1
1374               DO ji = 2,jpim1
1375                  IF (tmask(ji,jj,jk) == 1) THEN
1376                     !! surface "free" iron
1377                     trc2d(ji,jj,85) = xFree(ji,jj)
1378                  ENDIF
1379               ENDDO
1380            ENDDO
1381         ENDIF
1382
1383! I'm keeping this the same as before, but it looks like it should
1384! be i0100 and not i0200 - marc 8/5/17
1385         IF (jk.eq.i0200) THEN
1386            DO jj = 2,jpjm1
1387               DO ji = 2,jpim1
1388                  IF (tmask(ji,jj,jk) == 1) THEN
1389                     !! "free" iron at  100 m
1390                     trc2d(ji,jj,86) = xFree(ji,jj)
1391                  ENDIF
1392               ENDDO
1393            ENDDO
1394         ENDIF
1395
1396
1397         IF (jk.eq.i0200) THEN
1398            DO jj = 2,jpjm1
1399               DO ji = 2,jpim1
1400                  IF (tmask(ji,jj,jk) == 1) THEN
1401                     !! "free" iron at  200 m
1402                     trc2d(ji,jj,87) = xFree(ji,jj)
1403                  ENDIF
1404               ENDDO
1405            ENDDO
1406         ENDIF
1407
1408
1409         IF (jk.eq.i0500) THEN
1410            DO jj = 2,jpjm1
1411               DO ji = 2,jpim1
1412                  IF (tmask(ji,jj,jk) == 1) THEN
1413                     !! "free" iron at  500 m
1414                     trc2d(ji,jj,88) = xFree(ji,jj)
1415                  ENDIF
1416               ENDDO
1417            ENDDO
1418         ENDIF
1419
1420
1421         IF (jk.eq.i1000) THEN
1422            DO jj = 2,jpjm1
1423               DO ji = 2,jpim1
1424                  IF (tmask(ji,jj,jk) == 1) THEN
1425                     !! "free" iron at 1000 m
1426                     trc2d(ji,jj,89) = xFree(ji,jj)
1427                  ENDIF
1428               ENDDO
1429            ENDDO
1430         ENDIF
1431
1432
1433         IF (jk.eq.1) THEN
1434            DO jj = 2,jpjm1
1435               DO ji = 2,jpim1
1436                  IF (tmask(ji,jj,jk) == 1) THEN
1437                     !! AXY (27/06/12): extract "euphotic depth"
1438                     trc2d(ji,jj,90) = xze(ji,jj)
1439                  ENDIF
1440               ENDDO
1441            ENDDO
1442         ENDIF
1443
1444# if defined key_roam
1445         if (jk .eq. 1) then
1446            DO jj = 2,jpjm1
1447               DO ji = 2,jpim1
1448                  IF (tmask(ji,jj,jk) == 1) THEN
1449                     !! ROAM provisionally has access to a further 20 2D
1450                     !! diagnostics
1451                     !! surface wind
1452                     trc2d(ji,jj,91)  = trc2d(ji,jj,91)  + wndm(ji,jj)
1453                     !! atmospheric pCO2
1454                     trc2d(ji,jj,92)  = trc2d(ji,jj,92)  + f_pco2atm(ji,jj)
1455                     !! ocean pH
1456                     trc2d(ji,jj,93)  = trc2d(ji,jj,93)  + f_ph(ji,jj)
1457                     !! ocean pCO2
1458                     trc2d(ji,jj,94)  = trc2d(ji,jj,94)  + f_pco2w(ji,jj)
1459                     !! ocean H2CO3 conc.
1460                     trc2d(ji,jj,95)  = trc2d(ji,jj,95)  + f_h2co3(ji,jj)
1461                     !! ocean HCO3 conc.
1462                     trc2d(ji,jj,96)  = trc2d(ji,jj,96)  + f_hco3(ji,jj)
1463                     !! ocean CO3 conc.
1464                     trc2d(ji,jj,97)  = trc2d(ji,jj,97)  + f_co3(ji,jj)
1465                     !! air-sea CO2 flux
1466                     trc2d(ji,jj,98)  = trc2d(ji,jj,98)  + f_co2flux(ji,jj)
1467                  ENDIF
1468               ENDDO
1469           ENDDO
1470
1471            DO jj = 2,jpjm1
1472               DO ji = 2,jpim1
1473                  IF (tmask(ji,jj,jk) == 1) THEN
1474                     !! ocean omega calcite
1475                     trc2d(ji,jj,99)  = trc2d(ji,jj,99)  + f_omcal(ji,jj)
1476                     !! ocean omega aragonite
1477                     trc2d(ji,jj,100) = trc2d(ji,jj,100) + f_omarg(ji,jj)
1478                     !! ocean TDIC
1479                     trc2d(ji,jj,101) = trc2d(ji,jj,101) + f_TDIC(ji,jj)
1480                     !! ocean TALK
1481                     trc2d(ji,jj,102) = trc2d(ji,jj,102) + f_TALK(ji,jj)
1482                     !! surface kw660
1483                     trc2d(ji,jj,103) = trc2d(ji,jj,103) + f_kw660(ji,jj)
1484                     !! surface pressure
1485                     trc2d(ji,jj,104) = trc2d(ji,jj,104) + f_pp0(ji,jj)
1486                     !! air-sea O2 flux
1487                     trc2d(ji,jj,105) = trc2d(ji,jj,105) + f_o2flux(ji,jj)
1488                     !! ocean O2 saturation
1489                     trc2d(ji,jj,106) = trc2d(ji,jj,106) + f_o2sat(ji,jj)
1490                     !! depth calcite CCD
1491                     trc2d(ji,jj,107) = f2_ccd_cal(ji,jj)
1492                     !! depth aragonite CCD
1493                     trc2d(ji,jj,108) = f2_ccd_arg(ji,jj)
1494                  ENDIF
1495               ENDDO
1496            ENDDO
1497         endif
1498
1499         if (jk .eq. mbathy(ji,jj)) then
1500            DO jj = 2,jpjm1
1501               DO ji = 2,jpim1
1502                  IF (tmask(ji,jj,jk) == 1) THEN
1503                     !! seafloor omega calcite
1504                     trc2d(ji,jj,109) = f3_omcal(ji,jj,jk)
1505                     !! seafloor omega aragonite
1506                     trc2d(ji,jj,110) = f3_omarg(ji,jj,jk)
1507                  ENDIF
1508               ENDDO
1509            ENDDO
1510         endif
1511
1512         if (jk.eq.i0100) then
1513            DO jj = 2,jpjm1
1514               DO ji = 2,jpim1
1515                  IF (tmask(ji,jj,jk) == 1) THEN
1516                     !! diagnostic fields 111 to 117 calculated below
1517                     !! rain ratio at  100 m
1518                     trc2d(ji,jj,118) =                                      &
1519                                   ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
1520                  ENDIF
1521               ENDDO
1522            ENDDO
1523         endif
1524
1525         if (jk.eq.i0500) then
1526            DO jj = 2,jpjm1
1527               DO ji = 2,jpim1
1528                  IF (tmask(ji,jj,jk) == 1) THEN
1529                     !! rain ratio at  500 m
1530                     trc2d(ji,jj,119) =                                      &
1531                                   ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
1532                  ENDIF
1533               ENDDO
1534            ENDDO
1535         endif
1536
1537         if (jk.eq.i1000) then
1538            DO jj = 2,jpjm1
1539               DO ji = 2,jpim1
1540                  IF (tmask(ji,jj,jk) == 1) THEN
1541                     !! rain ratio at 1000 m
1542                     trc2d(ji,jj,120) =                                      &
1543                                   ffastca(ji,jj)/MAX(ffastc(ji,jj), rsmall)
1544                  ENDIF
1545               ENDDO
1546            ENDDO
1547         endif
1548
1549         if (jk.eq.mbathy(ji,jj)) then
1550            DO jj = 2,jpjm1
1551               DO ji = 2,jpim1
1552                  IF (tmask(ji,jj,jk) == 1) THEN
1553                     !! AXY (18/01/12): benthic flux diagnostics
1554                     trc2d(ji,jj,121) = f_sbenin_n(ji,jj)  + f_fbenin_n(ji,jj)
1555                     trc2d(ji,jj,122) = f_sbenin_fe(ji,jj) + f_fbenin_fe(ji,jj)
1556                     trc2d(ji,jj,123) = f_sbenin_c(ji,jj)  + f_fbenin_c(ji,jj)
1557                     trc2d(ji,jj,124) = f_fbenin_si(ji,jj)
1558                     trc2d(ji,jj,125) = f_fbenin_ca(ji,jj)
1559                     trc2d(ji,jj,126) = f_benout_n(ji,jj)
1560                     trc2d(ji,jj,127) = f_benout_fe(ji,jj)
1561                     trc2d(ji,jj,128) = f_benout_c(ji,jj)
1562                     trc2d(ji,jj,129) = f_benout_si(ji,jj)
1563                     trc2d(ji,jj,130) = f_benout_ca(ji,jj)
1564                  ENDIF
1565               ENDDO
1566            ENDDO
1567         endif
1568
1569         DO jj = 2,jpjm1
1570            DO ji = 2,jpim1
1571               IF (tmask(ji,jj,jk) == 1) THEN
1572                  !! diagnostics fields 131 to 135 calculated below
1573                  trc2d(ji,jj,136) = f_runoff(ji,jj)
1574                  !! AXY (19/07/12): amended to allow for riverine
1575                  !! nutrient addition below surface
1576                  trc2d(ji,jj,137) = trc2d(ji,jj,137) +                      &
1577                                     (f_riv_loc_n(ji,jj) * fse3t(ji,jj,jk))
1578                  trc2d(ji,jj,138) = trc2d(ji,jj,138) +                      &
1579                                     (f_riv_loc_si(ji,jj) * fse3t(ji,jj,jk))
1580                  trc2d(ji,jj,139) = trc2d(ji,jj,139) +                      &
1581                                     (f_riv_loc_c(ji,jj) * fse3t(ji,jj,jk))
1582                  trc2d(ji,jj,140) = trc2d(ji,jj,140) +                      &
1583                                     (f_riv_loc_alk(ji,jj) * fse3t(ji,jj,jk))
1584                  !! slow sinking detritus C production
1585                  trc2d(ji,jj,141) = trc2d(ji,jj,141) +                      &
1586                                     (fslowc(ji,jj)  * fse3t(ji,jj,jk))
1587               ENDIF
1588            ENDDO
1589         ENDDO
1590
1591         if (jk.eq.i0100) then
1592            DO jj = 2,jpjm1
1593               DO ji = 2,jpim1
1594                  IF (tmask(ji,jj,jk) == 1) THEN
1595                     !! slow detritus flux at  100 m
1596                     trc2d(ji,jj,142) = fslowcflux(ji,jj)
1597                  ENDIF
1598               ENDDO
1599            ENDDO
1600         endif
1601
1602         if (jk.eq.i0200) then
1603            DO jj = 2,jpjm1
1604               DO ji = 2,jpim1
1605                  IF (tmask(ji,jj,jk) == 1) THEN
1606                     !! slow detritus flux at  200 m
1607                     trc2d(ji,jj,143) = fslowcflux(ji,jj)
1608                  ENDIF
1609               ENDDO
1610            ENDDO
1611         endif
1612
1613
1614         if (jk.eq.i0500) then
1615            DO jj = 2,jpjm1
1616               DO ji = 2,jpim1
1617                  IF (tmask(ji,jj,jk) == 1) THEN
1618                     !! slow detritus flux at  500 m
1619                     trc2d(ji,jj,144) = fslowcflux(ji,jj)
1620                  ENDIF
1621               ENDDO
1622            ENDDO
1623         endif
1624
1625
1626         if (jk.eq.i1000) then
1627            DO jj = 2,jpjm1
1628               DO ji = 2,jpim1
1629                  IF (tmask(ji,jj,jk) == 1) THEN
1630                     !! slow detritus flux at 1000 m
1631                     trc2d(ji,jj,145) = fslowcflux(ji,jj)
1632                  ENDIF
1633               ENDDO
1634            ENDDO
1635         endif
1636
1637         DO jj = 2,jpjm1
1638            DO ji = 2,jpim1
1639               IF (tmask(ji,jj,jk) == 1) THEN
1640                  !! carbon     inventory
1641                  trc2d(ji,jj,146)  = trc2d(ji,jj,146)  + ftot_c(ji,jj)
1642                  !! alkalinity inventory
1643                  trc2d(ji,jj,147)  = trc2d(ji,jj,147)  + ftot_a(ji,jj)
1644                  !! oxygen     inventory
1645                  trc2d(ji,jj,148)  = trc2d(ji,jj,148)  + ftot_o2(ji,jj)
1646               ENDIF
1647            ENDDO
1648         ENDDO
1649
1650         if (jk.eq.mbathy(ji,jj)) then
1651            DO jj = 2,jpjm1
1652               DO ji = 2,jpim1
1653                  IF (tmask(ji,jj,jk) == 1) THEN
1654                     trc2d(ji,jj,149) = f_benout_lyso_ca(ji,jj)
1655                  ENDIF
1656               ENDDO
1657            ENDDO
1658         endif
1659
1660         DO jj = 2,jpjm1
1661            DO ji = 2,jpim1
1662               IF (tmask(ji,jj,jk) == 1) THEN
1663                  !! community respiration
1664                  trc2d(ji,jj,150) = fcomm_resp(ji,jj) * fse3t(ji,jj,jk)
1665               ENDIF
1666            ENDDO
1667         ENDDO
1668
1669         DO jj = 2,jpjm1
1670            DO ji = 2,jpim1
1671               IF (tmask(ji,jj,jk) == 1) THEN
1672        !!
1673        !! AXY (14/02/14): a Valentines Day gift to BASIN - a
1674                  !!                 shedload of new diagnostics that
1675                  !!                 they'll most likely never need!
1676                  !!                 (actually, as with all such gifts,
1677                  !!                 I'm giving them some things I'd like
1678                  !!                 myself!)
1679                  !!
1680                  !! ------------------------------------------------------
1681                  !! linear losses
1682                  !! non-diatom
1683                  trc2d(ji,jj,151) = trc2d(ji,jj,151) +                      &
1684                                     (fdpn2(ji,jj) * fse3t(ji,jj,jk))
1685                  !! diatom
1686                  trc2d(ji,jj,152) = trc2d(ji,jj,152) +                      &
1687                                     (fdpd2(ji,jj)  * fse3t(ji,jj,jk))
1688                  !! microzooplankton
1689                  trc2d(ji,jj,153) = trc2d(ji,jj,153) +                      &
1690                                     (fdzmi2(ji,jj) * fse3t(ji,jj,jk))
1691                  !! mesozooplankton
1692                  trc2d(ji,jj,154) = trc2d(ji,jj,154) +                      &
1693                                     (fdzme2(ji,jj) * fse3t(ji,jj,jk))
1694               ENDIF
1695            ENDDO
1696         ENDDO
1697
1698         DO jj = 2,jpjm1
1699            DO ji = 2,jpim1
1700               IF (tmask(ji,jj,jk) == 1) THEN
1701                  !! ------------------------------------------------------
1702                  !! microzooplankton grazing
1703                  !! microzooplankton messy -> N
1704                  trc2d(ji,jj,155) = trc2d(ji,jj,155) +                      &
1705                                     (xphi * (fgmipn(ji,jj) +                &
1706                                              fgmid(ji,jj)) * fse3t(ji,jj,jk))
1707                  !! microzooplankton messy -> D
1708                  trc2d(ji,jj,156) = trc2d(ji,jj,156) +                      &
1709                                     ((1. - xbetan) * finmi(ji,jj) *         &
1710                                      fse3t(ji,jj,jk))
1711                  !! microzooplankton messy -> DIC
1712                  trc2d(ji,jj,157) = trc2d(ji,jj,157) +                      &
1713                                     (xphi * ((xthetapn * fgmipn(ji,jj)) +   &
1714                                              fgmidc(ji,jj)) *               &
1715                                      fse3t(ji,jj,jk))
1716                  !! microzooplankton messy -> Dc
1717                  trc2d(ji,jj,158) = trc2d(ji,jj,158) +                      &
1718                                     ((1. - xbetac) * ficmi(ji,jj) *         &
1719                                      fse3t(ji,jj,jk))
1720                  !! microzooplankton excretion
1721                  trc2d(ji,jj,159) = trc2d(ji,jj,159) +                      &
1722                                     (fmiexcr(ji,jj) * fse3t(ji,jj,jk))
1723                  !! microzooplankton respiration
1724                  trc2d(ji,jj,160) = trc2d(ji,jj,160) +                      &
1725                                     (fmiresp(ji,jj) * fse3t(ji,jj,jk))
1726                  !! microzooplankton growth
1727                  trc2d(ji,jj,161) = trc2d(ji,jj,161) +                      &
1728                                     (fmigrow(ji,jj) * fse3t(ji,jj,jk))
1729               ENDIF
1730            ENDDO
1731         ENDDO
1732
1733         DO jj = 2,jpjm1
1734            DO ji = 2,jpim1
1735               IF (tmask(ji,jj,jk) == 1) THEN
1736                  !! ------------------------------------------------------
1737                  !! mesozooplankton grazing
1738                  !! mesozooplankton messy -> N
1739                  trc2d(ji,jj,162) = trc2d(ji,jj,162) +                      &
1740                                     (xphi *                                 &
1741                                      (fgmepn(ji,jj) + fgmepd(ji,jj) +       &
1742                                       fgmezmi(ji,jj) + fgmed(ji,jj)) *      &
1743                                      fse3t(ji,jj,jk))
1744                  !! mesozooplankton messy -> D
1745                  trc2d(ji,jj,163) = trc2d(ji,jj,163) +                      &
1746                                     ((1. - xbetan) * finme(ji,jj) *         &
1747                                      fse3t(ji,jj,jk))
1748                  !! mesozooplankton messy -> DIC
1749                  trc2d(ji,jj,164) = trc2d(ji,jj,164) +                      &
1750                                     (xphi *                                 &
1751                                      ((xthetapn * fgmepn(ji,jj)) +          &
1752                                       (xthetapd * fgmepd(ji,jj)) +          &
1753                                       (xthetazmi * fgmezmi(ji,jj)) +        &
1754                                      fgmedc(ji,jj)) * fse3t(ji,jj,jk))
1755                  !! mesozooplankton messy -> Dc
1756                  trc2d(ji,jj,165) = trc2d(ji,jj,165) +                      &
1757                                     ((1. - xbetac) * ficme(ji,jj) *         &
1758                                      fse3t(ji,jj,jk))
1759                  !! mesozooplankton excretion
1760                  trc2d(ji,jj,166) = trc2d(ji,jj,166) +                      &
1761                                     (fmeexcr(ji,jj) * fse3t(ji,jj,jk))
1762                  !! mesozooplankton respiration
1763                  trc2d(ji,jj,167) = trc2d(ji,jj,167) +                      &
1764                                     (fmeresp(ji,jj) * fse3t(ji,jj,jk))
1765                  !! mesozooplankton growth
1766                  trc2d(ji,jj,168) = trc2d(ji,jj,168) +                      &
1767                                     (fmegrow(ji,jj) * fse3t(ji,jj,jk))
1768               ENDIF
1769            ENDDO
1770         ENDDO
1771
1772         DO jj = 2,jpjm1
1773            DO ji = 2,jpim1
1774               IF (tmask(ji,jj,jk) == 1) THEN
1775                  !! ------------------------------------------------------
1776                  !! miscellaneous
1777                  !! detrital C remineralisation
1778                  trc2d(ji,jj,169) = trc2d(ji,jj,169) +                      &
1779                                     (fddc(ji,jj) * fse3t(ji,jj,jk))
1780                  !! microzoo grazing on detrital carbon
1781                  trc2d(ji,jj,170) = trc2d(ji,jj,170) +                      &
1782                                     (fgmidc(ji,jj)  * fse3t(ji,jj,jk))
1783                  !! mesozoo  grazing on detrital carbon
1784                  trc2d(ji,jj,171) = trc2d(ji,jj,171) +                      &
1785                                     (fgmedc(ji,jj)  * fse3t(ji,jj,jk))
1786                  !!
1787               ENDIF
1788            ENDDO
1789         ENDDO
1790
1791         !! ------------------------------------------------------
1792    !!
1793         !! AXY (23/10/14): extract primary production related
1794         !!                 surface fields to deal with diel
1795         !!                 cycle issues; hijacking BASIN 150m
1796         !!                 diagnostics to do so (see commented
1797         !!                 out diagnostics below this section)
1798         !!
1799         !! extract relevant BASIN fields at 150m
1800         if (jk .eq. i0150) then
1801            DO jj = 2,jpjm1
1802               DO ji = 2,jpim1
1803                  IF (tmask(ji,jj,jk) == 1) THEN
1804                     !! Pn PP
1805                     trc2d(ji,jj,172) = trc2d(ji,jj,4)
1806                     !! Pn linear loss
1807                     trc2d(ji,jj,173) = trc2d(ji,jj,151)
1808                     !! Pn non-linear loss
1809                     trc2d(ji,jj,174) = trc2d(ji,jj,5)
1810                     !! Pn grazing to Zmi
1811                     trc2d(ji,jj,175) = trc2d(ji,jj,11)
1812                     !! Pn grazing to Zme
1813                     trc2d(ji,jj,176) = trc2d(ji,jj,14)
1814                     !! Pd PP
1815                     trc2d(ji,jj,177) = trc2d(ji,jj,6)
1816                     !! Pd linear loss
1817                     trc2d(ji,jj,178) = trc2d(ji,jj,152)
1818                     !! Pd non-linear loss
1819                     trc2d(ji,jj,179) = trc2d(ji,jj,7)
1820                     !! Pd grazing to Zme
1821                     trc2d(ji,jj,180) = trc2d(ji,jj,15)
1822                     !! Zmi grazing on D
1823                     trc2d(ji,jj,181) = trc2d(ji,jj,12)
1824                     !! Zmi grazing on Dc
1825                     trc2d(ji,jj,182) = trc2d(ji,jj,170)
1826                     !! Zmi messy feeding loss to N
1827                     trc2d(ji,jj,183) = trc2d(ji,jj,155)
1828                     !! Zmi messy feeding loss to D
1829                     trc2d(ji,jj,184) = trc2d(ji,jj,156)
1830                     !! Zmi messy feeding loss to DIC
1831                     trc2d(ji,jj,185) = trc2d(ji,jj,157)
1832                     !! Zmi messy feeding loss to Dc
1833                     trc2d(ji,jj,186) = trc2d(ji,jj,158)
1834                     !! Zmi excretion
1835                     trc2d(ji,jj,187) = trc2d(ji,jj,159)
1836                     !! Zmi respiration
1837                     trc2d(ji,jj,188) = trc2d(ji,jj,160)
1838                     !! Zmi growth
1839                     trc2d(ji,jj,189) = trc2d(ji,jj,161)
1840                     !! Zmi linear loss
1841                     trc2d(ji,jj,190) = trc2d(ji,jj,153)
1842                     !! Zmi non-linear loss
1843                     trc2d(ji,jj,191) = trc2d(ji,jj,13)
1844                     !! Zmi grazing to Zme
1845                     trc2d(ji,jj,192) = trc2d(ji,jj,16)
1846                     !! Zme grazing on D
1847                     trc2d(ji,jj,193) = trc2d(ji,jj,17)
1848                     !! Zme grazing on Dc
1849                     trc2d(ji,jj,194) = trc2d(ji,jj,171)
1850                     !! Zme messy feeding loss to N
1851                     trc2d(ji,jj,195) = trc2d(ji,jj,162)
1852                     !! Zme messy feeding loss to D
1853                     trc2d(ji,jj,196) = trc2d(ji,jj,163)
1854                     !! Zme messy feeding loss to DIC
1855                     trc2d(ji,jj,197) = trc2d(ji,jj,164)
1856                     !! Zme messy feeding loss to Dc
1857                     trc2d(ji,jj,198) = trc2d(ji,jj,165)
1858                     !! Zme excretion
1859                     trc2d(ji,jj,199) = trc2d(ji,jj,166)
1860                     !! Zme respiration
1861                     trc2d(ji,jj,200) = trc2d(ji,jj,167)
1862                     !! Zme growth
1863                     trc2d(ji,jj,201) = trc2d(ji,jj,168)
1864                     !! Zme linear loss
1865                     trc2d(ji,jj,202) = trc2d(ji,jj,154)
1866                     !! Zme non-linear loss
1867                     trc2d(ji,jj,203) = trc2d(ji,jj,18)
1868                     !! Slow detritus production, N
1869                     trc2d(ji,jj,204) = trc2d(ji,jj,20)
1870                     !! Slow detritus remineralisation, N
1871                     trc2d(ji,jj,205) = trc2d(ji,jj,21)
1872                     !! Slow detritus production, C
1873                     trc2d(ji,jj,206) = trc2d(ji,jj,141)
1874                     !! Slow detritus remineralisation, C
1875                     trc2d(ji,jj,207) = trc2d(ji,jj,169)
1876                     !! Fast detritus production, N
1877                     trc2d(ji,jj,208) = trc2d(ji,jj,43)
1878                     !! Fast detritus remineralisation, N
1879                     trc2d(ji,jj,209) = trc2d(ji,jj,21)
1880                     !! Fast detritus production, C
1881                     trc2d(ji,jj,210) = trc2d(ji,jj,64)
1882                     !! Fast detritus remineralisation, C
1883                     trc2d(ji,jj,211) = trc2d(ji,jj,67)
1884                     !! Community respiration
1885                     trc2d(ji,jj,212) = trc2d(ji,jj,150)
1886                     !! Slow detritus N flux at 150 m
1887                     trc2d(ji,jj,213) = fslownflux(ji,jj)
1888                     !! Slow detritus C flux at 150 m
1889                     trc2d(ji,jj,214) = fslowcflux(ji,jj)
1890                     !! Fast detritus N flux at 150 m
1891                     trc2d(ji,jj,215) = ffastn(ji,jj)
1892                     !! Fast detritus C flux at 150 m
1893                     trc2d(ji,jj,216) = ffastc(ji,jj)
1894                  ENDIF
1895               ENDDO
1896            ENDDO
1897         endif
1898
1899         !!
1900         !! Jpalm (11-08-2014)
1901         !! Add UKESM1 diagnoatics
1902         !!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1903         if ((jk .eq. 1) .and.( jdms.eq.1)) then
1904            DO jj = 2,jpjm1
1905               DO ji = 2,jpim1
1906                  IF (tmask(ji,jj,jk) == 1) THEN
1907                     !! DMS surface concentration
1908                     trc2d(ji,jj,221) = dms_surf(ji,jj)
1909                     !! AXY (13/03/15): add in other DMS estimates
1910                     !! DMS surface concentration
1911                     trc2d(ji,jj,222) = dms_andr(ji,jj)
1912                     !! DMS surface concentration
1913                     trc2d(ji,jj,223) = dms_simo(ji,jj)
1914                     !! DMS surface concentration
1915                     trc2d(ji,jj,224) = dms_aran(ji,jj)
1916                     !! DMS surface concentration
1917                     trc2d(ji,jj,225) = dms_hall(ji,jj)
1918                  ENDIF
1919               ENDDO
1920            ENDDO
1921         endif
1922# endif
1923
1924         DO jj = 2,jpjm1
1925            DO ji = 2,jpim1
1926               IF (tmask(ji,jj,jk) == 1) THEN
1927                  !! other possible future diagnostics include:
1928                  !!   - integrated tracer values (esp. biological)
1929                  !!   - mixed layer tracer values
1930                  !!   - sub-surface chlorophyll maxima (plus depth)
1931                  !!   - different mixed layer depth criteria (T, sigma,
1932                  !!     var. sigma)
1933                  !!-------------------------------------------------------
1934                  !! Prepare 3D diagnostics
1935                  !!-------------------------------------------------------
1936                  !!
1937                  !! primary production 
1938                  trc3d(ji,jj,jk,1)  = ((fprn(ji,jj) + fprd(ji,jj)) *        &
1939                                        zphn(ji,jj))
1940                  !! detrital flux
1941                  trc3d(ji,jj,jk,2)  = fslownflux(ji,jj) + ffastn(ji,jj)
1942                  !! remineralisation
1943                  trc3d(ji,jj,jk,3)  = fregen(ji,jj) +                       &
1944                                       (freminn(ji,jj) * fse3t(ji,jj,jk))
1945               ENDIF
1946            ENDDO
1947         ENDDO
1948# if defined key_roam
1949         DO jj = 2,jpjm1
1950            DO ji = 2,jpim1
1951               IF (tmask(ji,jj,jk) == 1) THEN
1952                  !! pH
1953                  trc3d(ji,jj,jk,4)  = f3_pH(ji,jj,jk)
1954                  !! omega calcite
1955                  trc3d(ji,jj,jk,5)  = f3_omcal(ji,jj,jk)
1956               ENDIF
1957            ENDDO
1958         ENDDO
1959# else
1960         DO jj = 2,jpjm1
1961            DO ji = 2,jpim1
1962               IF (tmask(ji,jj,jk) == 1) THEN
1963                  !! fast Si flux
1964                  trc3d(ji,jj,jk,4)  = ffastsi(ji,jj)
1965               ENDIF
1966            ENDDO
1967         ENDDO
1968# endif
1969
1970      ENDIF   ! end of ln_diatrc option
1971
1972   END SUBROUTINE bio_medusa_diag
1973
1974#else
1975   !!======================================================================
1976   !!  Dummy module :                                   No MEDUSA bio-model
1977   !!======================================================================
1978CONTAINS
1979   SUBROUTINE bio_medusa_diag( )                    ! Empty routine
1980      WRITE(*,*) 'bio_medusa_diag: You should not have seen this print! error?'
1981   END SUBROUTINE bio_medusa_diag
1982#endif 
1983
1984   !!======================================================================
1985END MODULE bio_medusa_diag_mod
Note: See TracBrowser for help on using the repository browser.