Context Navigation

diadct.F90 @ 11038

Last change on this file since 11038 was 11038, checked in by rrenshaw, 5 years ago
code fixes
Property svn:executable set to ``*
File size: 132.1 KB

Line
1	MODULE diadct
2	!!=====================================================================
3	!! * MODULE diadct *
4	!! Ocean diagnostics: Compute the transport trough a sec.
5	!!===============================================================
6	!! History :
7	!!
8	!! original : 02/99 (Y Drillet)
9	!! addition : 10/01 (Y Drillet, R Bourdalle Badie)
10	!! : 10/05 (M Laborie) F90
11	!! addition : 04/07 (G Garric) Ice sections
12	!! bugfix : 04/07 (C Bricaud) test on sec%nb_point
13	!! initialisation of ztransp1,ztransp2,...
14	!! nemo_v_3_4: 09/2011 (C Bricaud)
15	!!
16	!!
17	!!----------------------------------------------------------------------
18	#if defined key_diadct
19	!!----------------------------------------------------------------------
20	!! 'key_diadct' :
21	!!----------------------------------------------------------------------
22	!!----------------------------------------------------------------------
23	!! dia_dct : Compute the transport through a sec.
24	!! dia_dct_init : Read namelist.
25	!! readsec : Read sections description and pathway
26	!! removepoints : Remove points which are common to 2 procs
27	!! transport : Compute transport for each sections
28	!! dia_dct_wri : Write tranports results in ascii files
29	!! interp : Compute temperature/salinity/density at U-point or V-point
30	!!
31	!!----------------------------------------------------------------------
32	!! * Modules used
33	USE oce ! ocean dynamics and tracers
34	USE dom_oce ! ocean space and time domain
35	USE phycst ! physical constants
36	USE in_out_manager ! I/O manager
37	USE iom ! I/0 library
38	USE daymod ! calendar
39	USE dianam ! build name of file
40	USE lib_mpp ! distributed memory computing library
41	#if defined key_lim2
42	USE ice_2
43	#endif
44	#if defined key_lim3
45	USE ice
46	#endif
47	USE domvvl
48	USE timing ! preformance summary
49	USE wrk_nemo ! working arrays
50
51	IMPLICIT NONE
52	PRIVATE
53
54	!! * Routine accessibility
55	PUBLIC dia_dct ! routine called by step.F90
56	PUBLIC dia_dct_init ! routine called by opa.F90
57	PUBLIC diadct_alloc ! routine called by nemo_init in nemogcm.F90
58	PRIVATE readsec
59	PRIVATE removepoints
60	PRIVATE transport
61	PRIVATE dia_dct_wri
62	PRIVATE dia_dct_wri_NOOS
63
64	#include "domzgr_substitute.h90"
65
66	!! * Shared module variables
67	LOGICAL, PUBLIC, PARAMETER :: lk_diadct = .TRUE. !: model-data diagnostics flag
68	LOGICAL, PUBLIC :: ln_dct_calc_noos_25h !: Calcuate noos 25 h means
69	LOGICAL, PUBLIC :: ln_dct_calc_noos_hr !: Calcuate noos hourly means
70	! JT
71	LOGICAL, PUBLIC :: ln_dct_iom_cont !: Use IOM Output?
72	LOGICAL, PUBLIC :: ln_dct_ascii !: Output ascii or binary
73	LOGICAL, PUBLIC :: ln_dct_h !: Output hourly instantaneous or mean values
74	! JT
75
76	!! * Module variables
77	INTEGER :: nn_dct = 1 ! Frequency of computation
78	INTEGER :: nn_dctwri = 1 ! Frequency of output
79	INTEGER :: nn_secdebug = 0 ! Number of the section to debug
80	INTEGER :: nn_dct_h = 1 ! Frequency of computation for NOOS hourly files
81	INTEGER :: nn_dctwri_h = 1 ! Frequency of output for NOOS hourly files
82
83	INTEGER, PARAMETER :: nb_class_max = 12 ! maximum number of classes, i.e. depth levels or density classes
84	INTEGER, PARAMETER :: nb_sec_max = 100 ! maximum number of sections
85	INTEGER, PARAMETER :: nb_point_max = 375 ! maximum number of points in a single section
86	INTEGER, PARAMETER :: nb_type_class = 14 ! types of calculations, i.e. pos transport, neg transport, heat transport, salt transport
87	INTEGER, PARAMETER :: nb_3d_vars = 5
88	INTEGER, PARAMETER :: nb_2d_vars = 2
89	INTEGER :: nb_sec
90
91	TYPE POINT_SECTION
92	INTEGER :: I,J
93	END TYPE POINT_SECTION
94
95	TYPE COORD_SECTION
96	REAL(wp) :: lon,lat
97	END TYPE COORD_SECTION
98
99	TYPE SECTION
100	CHARACTER(len=60) :: name ! name of the sec
101	LOGICAL :: llstrpond ! true if you want the computation of salt and
102	! heat transports
103	LOGICAL :: ll_ice_section ! ice surface and ice volume computation
104	LOGICAL :: ll_date_line ! = T if the section crosses the date-line
105	TYPE(COORD_SECTION), DIMENSION(2) :: coordSec ! longitude and latitude of the extremities of the sec
106	INTEGER :: nb_class ! number of boundaries for density classes
107	INTEGER, DIMENSION(nb_point_max) :: direction ! vector direction of the point in the section
108	CHARACTER(len=40),DIMENSION(nb_class_max) :: classname ! characteristics of the class
109	REAL(wp), DIMENSION(nb_class_max) :: zsigi ,&! in-situ density classes (99 if you don't want)
110	zsigp ,&! potential density classes (99 if you don't want)
111	zsal ,&! salinity classes (99 if you don't want)
112	ztem ,&! temperature classes(99 if you don't want)
113	zlay ! level classes (99 if you don't want)
114	REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport ! transport output
115	REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport_h ! transport output
116	REAL(wp) :: slopeSection ! slope of the section
117	INTEGER :: nb_point ! number of points in the section
118	TYPE(POINT_SECTION),DIMENSION(nb_point_max) :: listPoint ! list of points in the sections
119	END TYPE SECTION
120
121	TYPE(SECTION),DIMENSION(nb_sec_max) :: secs ! Array of sections
122
123	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d
124	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d
125	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d_h
126	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d_h
127	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z_hr_output
128
129	!! $Id$
130	CONTAINS
131
132
133	INTEGER FUNCTION diadct_alloc()
134	!!----------------------------------------------------------------------
135	!! * FUNCTION diadct_alloc *
136	!!----------------------------------------------------------------------
137	INTEGER :: ierr(5)
138	!!----------------------------------------------------------------------
139
140	ALLOCATE(transports_3d(nb_3d_vars,nb_sec_max,nb_point_max,jpk) , STAT=ierr(1) )
141	ALLOCATE(transports_2d(nb_2d_vars,nb_sec_max,nb_point_max) , STAT=ierr(2) )
142	ALLOCATE(transports_3d_h(nb_3d_vars,nb_sec_max,nb_point_max,jpk), STAT=ierr(3) )
143	ALLOCATE(transports_2d_h(nb_2d_vars,nb_sec_max,nb_point_max) , STAT=ierr(4) )
144	ALLOCATE(z_hr_output(nb_sec_max,3,nb_class_max) , STAT=ierr(5) )
145
146	diadct_alloc = MAXVAL( ierr )
147	IF( diadct_alloc /= 0 ) CALL ctl_warn('diadct_alloc: failed to allocate arrays')
148
149	END FUNCTION diadct_alloc
150
151	SUBROUTINE dia_dct_init
152	!!---------------------------------------------------------------------
153	!! * ROUTINE diadct *
154	!!
155	!! ** Purpose: Read the namelist parameters
156	!! Open output files
157	!!
158	!!---------------------------------------------------------------------
159	NAMELIST/namdct/nn_dct,ln_dct_h,nn_dctwri,ln_dct_ascii,nn_secdebug,ln_dct_calc_noos_25h,ln_dct_calc_noos_hr,ln_dct_iom_cont
160	INTEGER :: ios,jsec ! Local integer output status for namelist read
161	CHARACTER(len=3) :: jsec_str ! name of the jsec
162
163	IF( nn_timing == 1 ) CALL timing_start('dia_dct_init')
164
165	REWIND( numnam_ref ) ! Namelist namdct in reference namelist : Diagnostic: transport through sections
166	READ ( numnam_ref, namdct, IOSTAT = ios, ERR = 901)
167	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdct in reference namelist', lwp )
168
169	REWIND( numnam_cfg ) ! Namelist namdct in configuration namelist : Diagnostic: transport through sections
170	READ ( numnam_cfg, namdct, IOSTAT = ios, ERR = 902 )
171	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdct in configuration namelist', lwp )
172	IF(lwm) WRITE ( numond, namdct )
173
174	IF( ln_NOOS ) THEN
175
176	!Do calculation for daily, 25hourly mean every hour
177	nn_dct=3600./rdt ! hard coded for NOOS transects, to give 25 hour means from hourly instantaneous values
178
179	!write out daily, 25hourly mean every day
180	nn_dctwri=86400./rdt
181
182
183	!nn_dct_h=1 ! hard coded for NOOS transects, to give hourly data
184	! If you want hourly instantaneous values, you only do the calculation every 12 timesteps (if rdt = 300)
185	! and output it every 12 time steps. For this, you set the ln_dct_h to be True, and it calcuates it automatically
186	! if you want hourly mean values, set ln_dct_h to be False, and it will do the calculate every time step.
187	!
188	!SELECT CASE( ln_dct_h )
189	! CASE(.TRUE.)
190	! nn_dct_h=3600./rdt
191	! CASE(.FALSE.)
192	! nn_dct_h=1
193	!END SELECT
194
195	IF ( ln_dct_h ) THEN
196	nn_dct_h=3600./rdt
197	ELSE
198	nn_dct_h=1.
199	ENDIF
200
201	!JT write out hourly calculation every hour
202	nn_dctwri_h=3600./rdt
203	ENDIF
204
205	IF( lwp ) THEN
206	WRITE(numout,*) " "
207	WRITE(numout,*) "diadct_init: compute transports through sections "
208	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
209	IF( ln_NOOS ) THEN
210	WRITE(numout,*) " Calculate NOOS hourly output: ln_dct_calc_noos_hr = ",ln_dct_calc_noos_hr
211	WRITE(numout,*) " Calculate NOOS 25 hour mean output: ln_dct_calc_noos_25h = ",ln_dct_calc_noos_25h
212	WRITE(numout,*) " Use IOM Output: ln_dct_iom_cont = ",ln_dct_iom_cont
213	WRITE(numout,*) " Output in ASCII (True) or Binary (False): ln_dct_ascii = ",ln_dct_ascii
214	WRITE(numout,*) " Frequency of hourly computation - instantaneous (TRUE) or hourly mean (FALSE): ln_dct_h = ",ln_dct_h
215
216	WRITE(numout,*) " Frequency of computation hard coded to be every hour: nn_dct = ",nn_dct
217	WRITE(numout,*) " Frequency of write hard coded to average 25 instantaneous hour values: nn_dctwri = ",nn_dctwri
218	WRITE(numout,*) " Frequency of hourly computation (timestep) : nn_dct_h = ",nn_dct_h
219	WRITE(numout,*) " Frequency of hourly computation Not hard coded to be every timestep, or : nn_dct_h = ",nn_dct_h
220	WRITE(numout,*) " Frequency of hourly write hard coded to every hour: nn_dctwri_h = ",nn_dctwri_h
221	ELSE
222	WRITE(numout,*) " Frequency of computation: nn_dct = ",nn_dct
223	WRITE(numout,*) " Frequency of write: nn_dctwri = ",nn_dctwri
224	ENDIF
225
226	IF ( nn_secdebug .GE. 1 .AND. nn_secdebug .LE. nb_sec_max )THEN
227	WRITE(numout,*)" Debug section number: ", nn_secdebug
228	ELSE IF ( nn_secdebug == 0 )THEN ; WRITE(numout,*)" No section to debug"
229	ELSE IF ( nn_secdebug == -1 )THEN ; WRITE(numout,*)" Debug all sections"
230	ELSE ; WRITE(numout,*)" Wrong value for nn_secdebug : ",nn_secdebug
231	ENDIF
232
233	IF(nn_dct .GE. nn_dctwri .AND. MOD(nn_dct,nn_dctwri) .NE. 0) &
234	& CALL ctl_stop( 'diadct: nn_dct should be smaller and a multiple of nn_dctwri' )
235
236	ENDIF
237
238
239	IF ( ln_NOOS ) THEN
240	IF ( ln_dct_calc_noos_25h .or. ln_dct_calc_noos_hr ) CALL readsec
241	ENDIF
242
243	!open output file
244	IF( lwp ) THEN
245	IF( ln_NOOS ) THEN
246	WRITE(numout,*) "diadct_init: Open output files. ASCII? ",ln_dct_ascii
247	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
248	IF ( ln_dct_ascii ) THEN
249	if ( ln_dct_calc_noos_25h ) CALL ctl_opn( numdct_NOOS ,'NOOS_transport' , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .TRUE. )
250	if ( ln_dct_calc_noos_hr ) CALL ctl_opn( numdct_NOOS_h,'NOOS_transport_h', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .TRUE. )
251	ELSE
252	if ( ln_dct_calc_noos_25h ) CALL ctl_opn( numdct_NOOS ,'NOOS_transport_bin' , 'REPLACE', 'UNFORMATTED', 'SEQUENTIAL', -1, numout, .TRUE. )
253	if ( ln_dct_calc_noos_hr ) CALL ctl_opn( numdct_NOOS_h,'NOOS_transport_bin_h', 'REPLACE', 'UNFORMATTED', 'SEQUENTIAL', -1, numout, .TRUE. )
254	ENDIF
255	ELSE
256	CALL ctl_opn( numdct_vol, 'volume_transport', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. )
257	CALL ctl_opn( numdct_heat, 'heat_transport' , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. )
258	CALL ctl_opn( numdct_salt, 'salt_transport' , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. )
259	ENDIF
260	ENDIF
261
262	! Initialise arrays to zero
263	transports_3d(:,:,:,:) =0._wp
264	transports_2d(:,:,:) =0._wp
265	transports_3d_h(:,:,:,:) =0._wp
266	transports_2d_h(:,:,:) =0._wp
267	z_hr_output(:,:,:) =0._wp
268
269	IF( nn_timing == 1 ) CALL timing_stop('dia_dct_init')
270
271	IF (ln_dct_iom_cont) THEN
272	IF( lwp ) THEN
273	WRITE(numout,*) " "
274	WRITE(numout,*) "diadct_init: using xios iom_put for output: field_def.xml and iodef.xml code"
275	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
276	WRITE(numout,*) ""
277	WRITE(numout,*) " field_def.xml"
278	WRITE(numout,*) " ~~~~~~~~~~~~~"
279	WRITE(numout,*) ""
280	WRITE(numout,*) ""
281
282	WRITE(numout,*) ' <field_group id="noos_cross_section" domain_ref="1point" axis_ref="noos" operation="average">'
283
284	DO jsec=1,nb_sec
285	WRITE (jsec_str, "(I3.3)") jsec
286
287	WRITE(numout,*) ' <field id="noos_'//jsec_str//'_trans" long_name="' // TRIM(secs(jsec)%name) // ' 25h mean NOOS transport cross-section number: '//jsec_str//' (total, positive, negative)" unit="m^3/s" />'
288	WRITE(numout,*) ' <field id="noos_'//jsec_str//'_heat" long_name="' // TRIM(secs(jsec)%name) // ' 25h mean NOOS heat cross-section number: '//jsec_str//' (total, positive, negative)" unit="J/s" />'
289	WRITE(numout,*) ' <field id="noos_'//jsec_str//'_salt" long_name="' // TRIM(secs(jsec)%name) // ' 25h mean NOOS salt cross-section number: '//jsec_str//' (total, positive, negative)" unit="g/s" />'
290
291	ENDDO
292
293	WRITE(numout,*) ' </field_group>'
294
295	WRITE(numout,*) ""
296	WRITE(numout,*) ""
297	WRITE(numout,*) " iodef.xml"
298	WRITE(numout,*) " ~~~~~~~~~"
299	WRITE(numout,*) ""
300	WRITE(numout,*) ""
301
302	WRITE(numout,*) ' <file_group id="1d" output_freq="1d" output_level="10" enabled=".TRUE.">'
303	WRITE(numout,*) ""
304	WRITE(numout,*) ' <file id="noos_cross_section" name="NOOS_transport">'
305	DO jsec=1,nb_sec
306	WRITE (jsec_str, "(I3.3)") jsec
307
308	WRITE(numout,*) ' <field field_ref="noos_'//jsec_str//'_trans" />'
309	WRITE(numout,*) ' <field field_ref="noos_'//jsec_str//'_heat" />'
310	WRITE(numout,*) ' <field field_ref="noos_'//jsec_str//'_salt" />'
311
312	ENDDO
313	WRITE(numout,*) ' </file>'
314	WRITE(numout,*) ""
315	WRITE(numout,*) ' </file_group>'
316
317	WRITE(numout,*) ""
318	WRITE(numout,*) ""
319	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
320	WRITE(numout,*) ""
321
322	ENDIF
323	ENDIF
324
325
326	!
327	END SUBROUTINE dia_dct_init
328
329
330	SUBROUTINE dia_dct(kt)
331	!!---------------------------------------------------------------------
332	!! * ROUTINE diadct *
333	!!
334	!! Purpose :: Compute section transports and write it in numdct files
335	!!
336	!! Method :: All arrays initialised to zero in dct_init
337	!! Each nn_dct time step call subroutine 'transports' for
338	!! each section to sum the transports.
339	!! Each nn_dctwri time step:
340	!! Divide the arrays by the number of summations to gain
341	!! an average value
342	!! Call dia_dct_sum to sum relevant grid boxes to obtain
343	!! totals for each class (density, depth, temp or sal)
344	!! Call dia_dct_wri to write the transports into file
345	!! Reinitialise all relevant arrays to zero
346	!!---------------------------------------------------------------------
347	!! * Arguments
348	INTEGER,INTENT(IN) ::kt
349
350	!! * Local variables
351	INTEGER :: jsec, &! loop on sections
352	itotal ! nb_sec_maxnb_type_classnb_class_max
353	LOGICAL :: lldebug =.FALSE. ! debug a section
354
355	INTEGER , DIMENSION(1) :: ish ! tmp array for mpp_sum
356	INTEGER , DIMENSION(3) :: ish2 ! "
357	REAL(wp), POINTER, DIMENSION(:) :: zwork ! "
358	REAL(wp), POINTER, DIMENSION(:,:,:):: zsum ! "
359
360	!!---------------------------------------------------------------------
361	IF( nn_timing == 1 ) CALL timing_start('dia_dct')
362
363	IF( lk_mpp )THEN
364	itotal = nb_sec_maxnb_type_classnb_class_max
365	CALL wrk_alloc( itotal , zwork )
366	CALL wrk_alloc( nb_sec_max,nb_type_class,nb_class_max , zsum )
367	ENDIF
368	lldebug=.TRUE.
369	! Initialise arrays
370	zwork(:) = 0.0
371	zsum(:,:,:) = 0.0
372
373	IF( lwp .AND. kt==nit000+nn_dct-1 ) THEN
374	WRITE(numout,*) " "
375	WRITE(numout,*) "diadct: compute transport"
376	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~~~~~"
377	WRITE(numout,*) "nb_sec = ",nb_sec
378	WRITE(numout,*) "nn_dct = ",nn_dct
379	WRITE(numout,*) "ln_NOOS = ",ln_NOOS
380	WRITE(numout,*) "nb_sec = ",nb_sec
381	WRITE(numout,*) "nb_sec_max = ",nb_sec_max
382	WRITE(numout,*) "nb_type_class = ",nb_type_class
383	WRITE(numout,*) "nb_class_max = ",nb_class_max
384	ENDIF
385
386
387	IF ( ln_dct_calc_noos_25h ) THEN
388
389	! Compute transport and write only at nn_dctwri
390	IF ( MOD(kt,nn_dct)==0 .or. & ! compute transport every nn_dct time steps
391	(ln_NOOS .and. kt==nn_it000 ) ) THEN ! also include first time step when calculating NOOS 25 hour averages
392
393
394
395	DO jsec=1,nb_sec
396
397	lldebug=.FALSE.
398	IF( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. kt==nit000+nn_dct-1 .AND. lwp ) lldebug=.TRUE.
399
400	!Compute transport through section
401	CALL transport(secs(jsec),lldebug,jsec)
402
403
404	ENDDO
405
406	IF( MOD(kt,nn_dctwri)==0 )THEN
407
408
409
410	IF( lwp .AND. kt==nit000+nn_dctwri-1 ) WRITE(numout,*)" diadct: average and write at kt = ",kt
411
412
413	! Not 24 values, but 25! divide by ((nn_dctwri/nn_dct) +1)
414	!! divide arrays by nn_dctwri/nn_dct to obtain average
415	transports_3d(:,:,:,:)= transports_3d(:,:,:,:)/((nn_dctwri/nn_dct)+1.)
416	transports_2d(:,:,:) = transports_2d(:,:,:) /((nn_dctwri/nn_dct)+1.)
417
418	! Sum over each class
419	DO jsec=1,nb_sec
420	CALL dia_dct_sum(secs(jsec),jsec)
421	ENDDO
422
423	!Sum on all procs
424	IF( lk_mpp )THEN
425	zsum(:,:,:)=0.0_wp
426	ish(1) = nb_sec_maxnb_type_classnb_class_max
427	ish2 = (/nb_sec_max,nb_type_class,nb_class_max/)
428	DO jsec=1,nb_sec ; zsum(jsec,:,:) = secs(jsec)%transport(:,:) ; ENDDO
429	zwork(:)= RESHAPE(zsum(:,:,:), ish )
430	CALL mpp_sum(zwork, ish(1))
431	zsum(:,:,:)= RESHAPE(zwork,ish2)
432	DO jsec=1,nb_sec ; secs(jsec)%transport(:,:) = zsum(jsec,:,:) ; ENDDO
433	ENDIF
434
435	!Write the transport
436	DO jsec=1,nb_sec
437
438	IF( lwp .and. .not. ln_NOOS )CALL dia_dct_wri(kt,jsec,secs(jsec))
439	!IF( lwp .and. ln_NOOS )CALL dia_dct_wri_NOOS(kt,jsec,secs(jsec)) ! use NOOS specific formatting
440	IF( ln_NOOS )CALL dia_dct_wri_NOOS(kt,jsec,secs(jsec)) ! use NOOS specific formatting
441
442
443	!nullify transports values after writing
444	transports_3d(:,jsec,:,:)=0.0
445	transports_2d(:,jsec,: )=0.0
446	secs(jsec)%transport(:,:)=0.
447
448
449	IF ( ln_NOOS ) CALL transport(secs(jsec),lldebug,jsec) ! reinitialise for next 25 hour instantaneous average (overlapping values)
450
451
452
453	ENDDO
454
455	ENDIF
456
457	ENDIF
458
459	ENDIF
460	IF ( ln_dct_calc_noos_hr ) THEN
461	IF ( MOD(kt,nn_dct_h)==0 ) THEN ! compute transport every nn_dct_h time steps
462
463	DO jsec=1,nb_sec
464
465	lldebug=.FALSE.
466	IF( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. kt==nit000+nn_dct_h-1 .AND. lwp ) lldebug=.TRUE.
467
468	!Compute transport through section
469	CALL transport_h(secs(jsec),lldebug,jsec)
470
471	ENDDO
472
473	IF( MOD(kt,nn_dctwri_h)==0 )THEN
474
475	IF( lwp .AND. kt==nit000+nn_dctwri_h-1 )WRITE(numout,*)" diadct: average and write hourly files at kt = ",kt
476
477	!! divide arrays by nn_dctwri/nn_dct to obtain average
478	!
479	! JT - I think this is wrong. I think it is trying to sum over 25 hours, but only dividing by 24.
480	! I think it might work for daily cycles, but not for monthly cycles,
481	!
482	transports_3d_h(:,:,:,:)=transports_3d_h(:,:,:,:)/(nn_dctwri_h/nn_dct_h)
483	transports_2d_h(:,:,:) =transports_2d_h(:,:,:) /(nn_dctwri_h/nn_dct_h)
484
485	! Sum over each class
486	DO jsec=1,nb_sec
487	CALL dia_dct_sum_h(secs(jsec),jsec)
488	ENDDO
489
490	!Sum on all procs
491	IF( lk_mpp )THEN
492	ish(1) = nb_sec_maxnb_type_classnb_class_max
493	ish2 = (/nb_sec_max,nb_type_class,nb_class_max/)
494	DO jsec=1,nb_sec ; zsum(jsec,:,:) = secs(jsec)%transport_h(:,:) ; ENDDO
495	zwork(:)= RESHAPE(zsum(:,:,:), ish )
496	CALL mpp_sum(zwork, ish(1))
497	zsum(:,:,:)= RESHAPE(zwork,ish2)
498	DO jsec=1,nb_sec ; secs(jsec)%transport_h(:,:) = zsum(jsec,:,:) ; ENDDO
499	ENDIF
500
501	!Write the transport
502	DO jsec=1,nb_sec
503
504	IF( lwp .and. ln_NOOS ) THEN
505	CALL dia_dct_wri_NOOS_h(kt/nn_dctwri_h,jsec,secs(jsec)) ! use NOOS specific formatting
506	endif
507	!nullify transports values after writing
508	transports_3d_h(:,jsec,:,:)=0.0
509	transports_2d_h(:,jsec,:)=0.0
510	secs(jsec)%transport_h(:,:)=0.0
511
512	! for hourly mean or hourly instantaneous, you don't initialise! start with zero!
513	!IF ( ln_NOOS ) CALL transport_h(secs(jsec),lldebug,jsec) ! reinitialise for next 25 hour instantaneous average (overlapping values)
514
515	ENDDO
516
517	ENDIF
518
519	ENDIF
520
521	ENDIF
522
523	IF( lk_mpp )THEN
524	itotal = nb_sec_maxnb_type_classnb_class_max
525	CALL wrk_dealloc( itotal , zwork )
526	CALL wrk_dealloc( nb_sec_max,nb_type_class,nb_class_max , zsum )
527	ENDIF
528
529	IF( nn_timing == 1 ) CALL timing_stop('dia_dct')
530	!
531	END SUBROUTINE dia_dct
532
533	SUBROUTINE readsec
534	!!---------------------------------------------------------------------
535	!! * ROUTINE readsec *
536	!!
537	!! ** Purpose:
538	!! Read a binary file(section_ijglobal.diadct)
539	!! generated by the tools "NEMOGCM/TOOLS/SECTIONS_DIADCT"
540	!!
541	!!
542	!!---------------------------------------------------------------------
543	!! * Local variables
544	INTEGER :: iptglo , iptloc ! Global and local number of points for a section
545	INTEGER :: isec, iiglo, ijglo, iiloc, ijloc,iost,i1 ,i2 ! temporary integer
546	INTEGER :: jsec, jpt ! dummy loop indices
547
548	INTEGER, DIMENSION(2) :: icoord
549	CHARACTER(len=160) :: clname !filename
550	CHARACTER(len=200) :: cltmp
551	CHARACTER(len=200) :: clformat !automatic format
552	TYPE(POINT_SECTION),DIMENSION(nb_point_max) ::coordtemp !contains listpoints coordinates
553	!read in the file
554	INTEGER, POINTER, DIMENSION(:) :: directemp !contains listpoints directions
555	!read in the files
556	LOGICAL :: llbon ,&!local logical
557	lldebug !debug the section
558	!!-------------------------------------------------------------------------------------
559	CALL wrk_alloc( nb_point_max, directemp )
560
561	! ---------------------------------------------------------------------------------------------
562	! This section seems to be needed, else later reads find no sections in section_ijglobal.diadct
563	!open input file
564	!---------------
565	!write(numout,*) 'dct low-level pre open: little endian '
566	!OPEN(UNIT=107,FILE='section_ijglobal.diadct', FORM='UNFORMATTED', ACCESS='SEQUENTIAL', STATUS='OLD',convert='LITTLE_ENDIAN')
567
568	write(numout,*) 'dct low-level pre open: big endian :',nproc,narea
569	OPEN(UNIT=107,FILE='section_ijglobal.diadct', FORM='UNFORMATTED', ACCESS='SEQUENTIAL', STATUS='OLD',convert='BIG_ENDIAN')
570
571	!write(numout,*) 'dct low-level pre open: SWAP '
572	!OPEN(UNIT=107,FILE='section_ijglobal.diadct', FORM='UNFORMATTED', ACCESS='SEQUENTIAL', STATUS='OLD',convert='SWAP')
573
574	!write(numout,*) 'dct low-level pre open: NATIVE '
575	!OPEN(UNIT=107,FILE='section_ijglobal.diadct', FORM='UNFORMATTED', ACCESS='SEQUENTIAL', STATUS='OLD',convert='NATIVE')
576
577	READ(107) isec
578	CLOSE(107)
579	! ---------------------------------------------------------------------------------------------
580
581	CALL ctl_opn( numdct_in, 'section_ijglobal.diadct', 'OLD', 'UNFORMATTED', 'SEQUENTIAL', -1, numout, .TRUE. )
582
583	!---------------
584	!Read input file
585	!---------------
586
587	DO jsec=1,nb_sec_max !loop on the nb_sec sections
588
589	IF ( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) ) &
590	& WRITE(numout,*)'debugging for section number: ',jsec
591
592	!initialization
593	!---------------
594	secs(jsec)%name=''
595	secs(jsec)%llstrpond = .FALSE.
596	secs(jsec)%ll_ice_section = .FALSE.
597	secs(jsec)%ll_date_line = .FALSE.
598	secs(jsec)%nb_class = 0
599	secs(jsec)%zsigi = 99._wp
600	secs(jsec)%zsigp = 99._wp
601	secs(jsec)%zsal = 99._wp
602	secs(jsec)%ztem = 99._wp
603	secs(jsec)%zlay = 99._wp
604	secs(jsec)%transport = 0._wp
605	secs(jsec)%transport_h = 0._wp
606	secs(jsec)%nb_point = 0
607
608	!read section's number / name / computing choices / classes / slopeSection / points number
609	!-----------------------------------------------------------------------------------------
610
611	READ(numdct_in,iostat=iost) isec
612	IF (iost .NE. 0 ) then
613	write(numout,*) 'unable to read section_ijglobal.diadct. iost = ',iost
614	EXIT !end of file
615	ENDIF
616
617
618	WRITE(cltmp,'(a,i4.4,a,i4.4)')'diadct: read sections : Problem of section number: isec= ',isec,' and jsec= ',jsec
619
620
621	IF( jsec .NE. isec ) CALL ctl_stop( cltmp )
622
623	READ(numdct_in)secs(jsec)%name
624	READ(numdct_in)secs(jsec)%llstrpond
625	READ(numdct_in)secs(jsec)%ll_ice_section
626	READ(numdct_in)secs(jsec)%ll_date_line
627	READ(numdct_in)secs(jsec)%coordSec
628	READ(numdct_in)secs(jsec)%nb_class
629	READ(numdct_in)secs(jsec)%zsigi
630	READ(numdct_in)secs(jsec)%zsigp
631	READ(numdct_in)secs(jsec)%zsal
632	READ(numdct_in)secs(jsec)%ztem
633	READ(numdct_in)secs(jsec)%zlay
634	READ(numdct_in)secs(jsec)%slopeSection
635	READ(numdct_in)iptglo
636
637	IF ( ln_NOOS ) THEN
638	WRITE(numout,*) 'Section name = ',TRIM(secs(jsec)%name)
639	WRITE(numout,*) 'coordSec = ',secs(jsec)%coordSec
640	WRITE(numout,*) 'iptglo = ',iptglo
641	ENDIF
642
643	!debug
644	!-----
645
646	IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )THEN
647
648	WRITE(clformat,'(a,i2,a)') '(A40,', nb_class_max,'(f8.3,1X))'
649
650	WRITE(numout,*) " Section name : ",TRIM(secs(jsec)%name)
651	WRITE(numout,*) " Compute temperature and salinity transports ? ",secs(jsec)%llstrpond
652	WRITE(numout,*) " Compute ice transport ? ",secs(jsec)%ll_ice_section
653	WRITE(numout,*) " Section crosses date-line ? ",secs(jsec)%ll_date_line
654	WRITE(numout,*) " Slope section : ",secs(jsec)%slopeSection
655	WRITE(numout,*) " Number of points in the section: ",iptglo
656	WRITE(numout,*) " Number of classes ",secs(jsec)%nb_class
657	WRITE(numout,clformat)" Insitu density classes : ",secs(jsec)%zsigi
658	WRITE(numout,clformat)" Potential density classes : ",secs(jsec)%zsigp
659	WRITE(numout,clformat)" Salinity classes : ",secs(jsec)%zsal
660	WRITE(numout,clformat)" Temperature classes : ",secs(jsec)%ztem
661	WRITE(numout,clformat)" Depth classes : ",secs(jsec)%zlay
662	ENDIF
663
664	IF( iptglo .NE. 0 )THEN
665
666	!read points'coordinates and directions
667	!--------------------------------------
668	IF ( ln_NOOS ) THEN
669	WRITE(numout,*) 'Coords and direction read'
670	ENDIF
671
672	coordtemp(:) = POINT_SECTION(0,0) !list of points read
673	directemp(:) = 0 !value of directions of each points
674	DO jpt=1,iptglo
675	READ(numdct_in)i1,i2
676	coordtemp(jpt)%I = i1
677	coordtemp(jpt)%J = i2
678	ENDDO
679	READ(numdct_in)directemp(1:iptglo)
680
681	!debug
682	!-----
683	IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )THEN
684	WRITE(numout,*)" List of points in global domain:"
685	DO jpt=1,iptglo
686	WRITE(numout,*)' # I J ',jpt,coordtemp(jpt),directemp(jpt)
687	ENDDO
688	ENDIF
689
690	!Now each proc selects only points that are in its domain:
691	!--------------------------------------------------------
692	iptloc = 0 !initialize number of points selected
693	DO jpt=1,iptglo !loop on listpoint read in the file
694
695	iiglo=coordtemp(jpt)%I ! global coordinates of the point
696	ijglo=coordtemp(jpt)%J ! "
697
698	IF( iiglo==jpidta .AND. nimpp==1 ) iiglo = 2
699
700	iiloc=iiglo-jpizoom+1-nimpp+1 ! local coordinates of the point
701	ijloc=ijglo-jpjzoom+1-njmpp+1 ! "
702
703	!verify if the point is on the local domain:(1,nlei)*(1,nlej)
704	IF( iiloc .GE. 1 .AND. iiloc .LE. nlei .AND. &
705	ijloc .GE. 1 .AND. ijloc .LE. nlej )THEN
706	iptloc = iptloc + 1 ! count local points
707	secs(jsec)%listPoint(iptloc) = POINT_SECTION(mi0(iiglo),mj0(ijglo)) ! store local coordinates
708	secs(jsec)%direction(iptloc) = directemp(jpt) ! store local direction
709	ENDIF
710
711	ENDDO
712
713	secs(jsec)%nb_point=iptloc !store number of section's points
714
715	!debug
716	!-----
717	IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )THEN
718	WRITE(numout,*)" List of points selected by the proc:"
719	DO jpt = 1,iptloc
720	iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1
721	ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1
722	WRITE(numout,*)' # I J : ',iiglo,ijglo
723	ENDDO
724	ENDIF
725
726	IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN
727	DO jpt = 1,iptloc
728	iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1
729	ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1
730	ENDDO
731	ENDIF
732
733	!remove redundant points between processors
734	!------------------------------------------
735	lldebug = .FALSE. ; IF ( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. lwp ) lldebug = .TRUE.
736	IF( iptloc .NE. 0 )THEN
737	CALL removepoints(secs(jsec),'I','top_list',lldebug)
738	CALL removepoints(secs(jsec),'I','bot_list',lldebug)
739	CALL removepoints(secs(jsec),'J','top_list',lldebug)
740	CALL removepoints(secs(jsec),'J','bot_list',lldebug)
741	ENDIF
742	IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN
743	DO jpt = 1,secs(jsec)%nb_point
744	iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1
745	ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1
746	ENDDO
747	ENDIF
748
749	!debug
750	!-----
751	IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )THEN
752	WRITE(numout,*)" List of points after removepoints:"
753	iptloc = secs(jsec)%nb_point
754	DO jpt = 1,iptloc
755	iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1
756	ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1
757	WRITE(numout,*)' # I J : ',iiglo,ijglo
758	CALL FLUSH(numout)
759	ENDDO
760	ENDIF
761
762	ELSE ! iptglo = 0
763	IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )&
764	WRITE(numout,*)' No points for this section.'
765	ENDIF
766
767	ENDDO !end of the loop on jsec
768
769	nb_sec = jsec-1 !number of section read in the file
770
771	IF( lwp ) WRITE(numout,*)'diadct: read sections: Finished readsec.'
772
773	CALL wrk_dealloc( nb_point_max, directemp )
774	!
775	END SUBROUTINE readsec
776
777	SUBROUTINE removepoints(sec,cdind,cdextr,ld_debug)
778	!!---------------------------------------------------------------------------
779	!! *** function removepoints
780	!!
781	!! ** Purpose :: Remove points which are common to 2 procs
782	!!
783	!----------------------------------------------------------------------------
784	!! * arguments
785	TYPE(SECTION),INTENT(INOUT) :: sec
786	CHARACTER(len=1),INTENT(IN) :: cdind ! = 'I'/'J'
787	CHARACTER(len=8),INTENT(IN) :: cdextr ! = 'top_list'/'bot_list'
788	LOGICAL,INTENT(IN) :: ld_debug
789
790	!! * Local variables
791	INTEGER :: iextr ,& !extremity of listpoint that we verify
792	iind ,& !coord of listpoint that we verify
793	itest ,& !indice value of the side of the domain
794	!where points could be redundant
795	isgn ,& ! isgn= 1 : scan listpoint from start to end
796	! isgn=-1 : scan listpoint from end to start
797	istart,iend !first and last points selected in listpoint
798	INTEGER :: jpoint !loop on list points
799	INTEGER, POINTER, DIMENSION(:) :: idirec !contains temporary sec%direction
800	INTEGER, POINTER, DIMENSION(:,:) :: icoord !contains temporary sec%listpoint
801	!----------------------------------------------------------------------------
802	CALL wrk_alloc( nb_point_max, idirec )
803	CALL wrk_alloc( 2, nb_point_max, icoord )
804
805	IF( ld_debug )WRITE(numout,*)' -------------------------'
806	IF( ld_debug )WRITE(numout,*)' removepoints in listpoint'
807
808	!iextr=extremity of list_point that we verify
809	IF ( cdextr=='bot_list' )THEN ; iextr=1 ; isgn=1
810	ELSE IF ( cdextr=='top_list' )THEN ; iextr=sec%nb_point ; isgn=-1
811	ELSE ; CALL ctl_stop("removepoints :Wrong value for cdextr")
812	ENDIF
813
814	!which coordinate shall we verify ?
815	IF ( cdind=='I' )THEN ; itest=nlei ; iind=1
816	ELSE IF ( cdind=='J' )THEN ; itest=nlej ; iind=2
817	ELSE ; CALL ctl_stop("removepoints :Wrong value for cdind")
818	ENDIF
819
820	IF( ld_debug )THEN
821	WRITE(numout,*)' case: coord/list extr/domain side'
822	WRITE(numout,*)' ', cdind,' ',cdextr,' ',itest
823	WRITE(numout,*)' Actual number of points: ',sec%nb_point
824	ENDIF
825
826	icoord(1,1:nb_point_max) = sec%listPoint%I
827	icoord(2,1:nb_point_max) = sec%listPoint%J
828	idirec = sec%direction
829	sec%listPoint = POINT_SECTION(0,0)
830	sec%direction = 0
831
832	jpoint=iextr+isgn
833	DO WHILE( jpoint .GE. 1 .AND. jpoint .LE. sec%nb_point )
834	IF( icoord( iind,jpoint-isgn ) == itest .AND. icoord( iind,jpoint ) == itest )THEN ; jpoint=jpoint+isgn
835	ELSE ; EXIT
836	ENDIF
837	ENDDO
838
839	IF( cdextr=='bot_list')THEN ; istart=jpoint-1 ; iend=sec%nb_point
840	ELSE ; istart=1 ; iend=jpoint+1
841	ENDIF
842
843	sec%listPoint(1:1+iend-istart)%I = icoord(1,istart:iend)
844	sec%listPoint(1:1+iend-istart)%J = icoord(2,istart:iend)
845	sec%direction(1:1+iend-istart) = idirec(istart:iend)
846	sec%nb_point = iend-istart+1
847
848	IF( ld_debug )THEN
849	WRITE(numout,*)' Number of points after removepoints :',sec%nb_point
850	WRITE(numout,*)' sec%direction after removepoints :',sec%direction(1:sec%nb_point)
851	ENDIF
852
853	CALL wrk_dealloc( nb_point_max, idirec )
854	CALL wrk_dealloc( 2, nb_point_max, icoord )
855	END SUBROUTINE removepoints
856
857	SUBROUTINE transport(sec,ld_debug,jsec)
858	!!-------------------------------------------------------------------------------------------
859	!! * ROUTINE transport *
860	!!
861	!! Purpose :: Compute the transport for each point in a section
862	!!
863	!! Method :: Loop over each segment, and each vertical level and add the transport
864	!! Be aware :
865	!! One section is a sum of segments
866	!! One segment is defined by 2 consecutive points in sec%listPoint
867	!! All points of sec%listPoint are positioned on the F-point of the cell
868	!!
869	!! There are two loops:
870	!! loop on the segment between 2 nodes
871	!! loop on the level jk !!
872	!!
873	!! ** Output: Arrays containing the volume,density,salinity,temperature etc
874	!! transports for each point in a section, summed over each nn_dct.
875	!!
876	!!-------------------------------------------------------------------------------------------
877	!! * Arguments
878	TYPE(SECTION),INTENT(INOUT) :: sec
879	LOGICAL ,INTENT(IN) :: ld_debug
880	INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section
881
882	!! * Local variables
883	INTEGER :: jk,jseg,jclass, &!loop on level/segment/classes
884	isgnu , isgnv !
885	REAL(wp):: zumid , zvmid , &!U/V velocity on a cell segment
886	zumid_ice , zvmid_ice , &!U/V ice velocity
887	zTnorm !transport of velocity through one cell's sides
888	REAL(wp):: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point
889
890	TYPE(POINT_SECTION) :: k
891	!!--------------------------------------------------------
892
893	IF( ld_debug )WRITE(numout,*)' Compute transport (jsec,sec%nb_point,sec%slopeSection) : ', jsec,sec%nb_point,sec%slopeSection
894
895	!---------------------------!
896	! COMPUTE TRANSPORT !
897	!---------------------------!
898	IF(sec%nb_point .NE. 0)THEN
899
900	!----------------------------------------------------------------------------------------------------
901	!----------------------------------------------------------------------------------------------------
902	!----------------------------------------------------------------------------------------------------
903	!
904	!
905	! ! ! ! JT 1/09/2018 - changing convention. Always direction + is toward left hand of section
906	!
907	! Making sign of the velocities used to calculate the volume transport a function of direction, not slopesection
908	! (isgnu, isgnv)
909	!
910	! They vary for each segment of the section.
911	!
912	!----------------------------------------------------------------------------------------------------
913	!----------------------------------------------------------------------------------------------------
914	!----------------------------------------------------------------------------------------------------
915	!Compute sign for velocities:
916	!
917	!convention:
918	! non horizontal section: direction + is toward left hand of section
919	! horizontal section: direction + is toward north of section
920	!
921	!
922	! slopeSection < 0 slopeSection > 0 slopeSection=inf slopeSection=0
923	! ---------------- ----------------- --------------- --------------
924	!
925	! isgnv=1 direction +
926	! ______ _____ ______
927	! \| //\| \| \| direction +
928	! \| isgnu=1 // \| \|isgnu=1 \|isgnu=1 /\|\
929	! \|_______ // ______\| \\ \| ---\ \|
930	! \| \| isgnv=-1 \\ \| \| ---/ direction + ____________
931	! \| \| __\\\| \|
932	! \| \| direction + \| isgnv=1
933	!
934	!----------------------------------------------------------------------------------------------------
935
936	IF( ld_debug )write(numout,*)"sec%slopeSection isgnu isgnv ",sec%slopeSection,isgnu,isgnv
937
938	!--------------------------------------!
939	! LOOP ON THE SEGMENT BETWEEN 2 NODES !
940	!--------------------------------------!
941	DO jseg=1,MAX(sec%nb_point-1,0)
942
943
944	!Compute sign for velocities:
945
946	!isgnu = 1
947	!isgnv = 1
948	!
949	!changing sign of u and v is dependent on the direction of the section.
950	!isgnu = 1
951	!isgnv = 1
952	!SELECT CASE( sec%direction(jseg) )
953	!CASE(0) ; isgnv = -1
954	!CASE(3) ; isgnu = -1
955	!END SELECT
956
957
958	SELECT CASE( sec%direction(jseg) )
959	CASE(0)
960	isgnu = 1
961	isgnv = -1
962	CASE(1)
963	isgnu = 1
964	isgnv = 1
965	CASE(2)
966	isgnu = 1
967	isgnv = 1
968	CASE(3)
969	isgnu = -1
970	isgnv = 1
971	END SELECT
972
973	!-------------------------------------------------------------------------------------------
974	! Select the appropriate coordinate for computing the velocity of the segment
975	! Corrected by JT 01/09/2018 (#)
976	!
977	! CASE(0) Case (2)
978	! ------- --------
979	! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j)
980	! F(i,j)---------#V(i,j)-------F(i+1,j) \|
981	! --------> \|
982	! \| \|
983	! \| \|
984	! Case (3) \| U(i,j)
985	! -------- \| \|
986	! V \|
987	! listPoint(jseg+1) F(i,j+1) \|
988	! \| \|
989	! \| \|
990	! \| listPoint(jseg+1) F(i,j-1)
991	! ^ \|
992	! \| \|
993	! \| U(i,j+1)
994	! \| \| Case(1)
995	! \| \| ------
996	! \|
997	! \| listPoint(jseg+1) listPoint(jseg)
998	! \| F(i-1,j)----------#V(i-1,j) ------#f(i,j)
999	! listPoint(jseg) F(i,j) <-------
1000	!
1001	!-------------------------------------------------------------------------------------------
1002
1003	SELECT CASE( sec%direction(jseg) )
1004	CASE(0) ; k = sec%listPoint(jseg)
1005	CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J)
1006	CASE(2) ; k = sec%listPoint(jseg)
1007	CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1)
1008	END SELECT
1009
1010	!---------------------------\|
1011	! LOOP ON THE LEVEL \|
1012	!---------------------------\|
1013	!Sum of the transport on the vertical
1014	DO jk=1,mbathy(k%I,k%J)
1015
1016	! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point
1017	SELECT CASE( sec%direction(jseg) )
1018	CASE(0,1)
1019	ztn = interp(k%I,k%J,jk,'V',0 )
1020	zsn = interp(k%I,k%J,jk,'V',1 )
1021	zrhop = interp(k%I,k%J,jk,'V',2)
1022	zrhoi = interp(k%I,k%J,jk,'V',3)
1023	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) vmask(k%I,k%J,1)
1024	CASE(2,3)
1025	ztn = interp(k%I,k%J,jk,'U',0)
1026	zsn = interp(k%I,k%J,jk,'U',1)
1027	zrhop = interp(k%I,k%J,jk,'U',2)
1028	zrhoi = interp(k%I,k%J,jk,'U',3)
1029	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) umask(k%I,k%J,1)
1030	END SELECT
1031
1032	zfsdep= fsdept(k%I,k%J,jk)
1033
1034	!compute velocity with the correct direction
1035	SELECT CASE( sec%direction(jseg) )
1036	CASE(0,1)
1037	zumid=0.
1038	zvmid=isgnvvn(k%I,k%J,jk)vmask(k%I,k%J,jk)
1039	CASE(2,3)
1040	zumid=isgnuun(k%I,k%J,jk)umask(k%I,k%J,jk)
1041	zvmid=0.
1042	END SELECT
1043
1044	!zTnorm=transport through one cell;
1045	!velocity* cell's length * cell's thickness
1046	zTnorm=zumide2u(k%I,k%J) fse3u(k%I,k%J,jk)+ &
1047	zvmide1v(k%I,k%J) fse3v(k%I,k%J,jk)
1048
1049	#if ! defined key_vvl
1050	!add transport due to free surface
1051	IF( jk==1 )THEN
1052	zTnorm = zTnorm + zumid* e2u(k%I,k%J) * zsshn * umask(k%I,k%J,jk) + &
1053	zvmid* e1v(k%I,k%J) * zsshn * vmask(k%I,k%J,jk)
1054	ENDIF
1055	#endif
1056	!COMPUTE TRANSPORT
1057
1058	transports_3d(1,jsec,jseg,jk) = transports_3d(1,jsec,jseg,jk) + zTnorm
1059
1060	IF ( sec%llstrpond ) THEN
1061	transports_3d(2,jsec,jseg,jk) = transports_3d(2,jsec,jseg,jk) + zTnorm * zrhoi
1062	transports_3d(3,jsec,jseg,jk) = transports_3d(3,jsec,jseg,jk) + zTnorm * zrhop
1063	transports_3d(4,jsec,jseg,jk) = transports_3d(4,jsec,jseg,jk) + zTnorm * 4.e+3_wp * (ztn+273.15) * 1026._wp
1064	transports_3d(5,jsec,jseg,jk) = transports_3d(5,jsec,jseg,jk) + zTnorm * 0.001 * zsn * 1026._wp
1065	ENDIF
1066
1067	ENDDO !end of loop on the level
1068
1069	#if defined key_lim2 \|\| defined key_lim3
1070
1071	!ICE CASE
1072	!------------
1073	IF( sec%ll_ice_section )THEN
1074	SELECT CASE (sec%direction(jseg))
1075	CASE(0)
1076	zumid_ice = 0
1077	zvmid_ice = isgnv0.5(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1))
1078	CASE(1)
1079	zumid_ice = 0
1080	zvmid_ice = isgnv0.5(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1))
1081	CASE(2)
1082	zvmid_ice = 0
1083	zumid_ice = isgnu0.5(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1))
1084	CASE(3)
1085	zvmid_ice = 0
1086	zumid_ice = isgnu0.5(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1))
1087	END SELECT
1088
1089	zTnorm=zumid_icee2u(k%I,k%J)+zvmid_icee1v(k%I,k%J)
1090
1091	transports_2d(1,jsec,jseg) = transports_2d(1,jsec,jseg) + (zTnorm)* &
1092	(1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1093	*(hsnif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J) + &
1094	hicif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1095	+zice_vol_pos
1096	transports_2d(2,jsec,jseg) = transports_2d(2,jsec,jseg) + (zTnorm)* &
1097	(1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1098	+zice_surf_pos
1099
1100	ENDIF !end of ice case
1101	#endif
1102
1103	ENDDO !end of loop on the segment
1104
1105	ENDIF !end of sec%nb_point =0 case
1106	!
1107	END SUBROUTINE transport
1108
1109	SUBROUTINE transport_h(sec,ld_debug,jsec)
1110	!!-------------------------------------------------------------------------------------------
1111	!! * ROUTINE hourly_transport *
1112	!!
1113	!! Exactly as routine transport but for data summed at
1114	!! each time step and averaged each hour
1115	!!
1116	!! Purpose :: Compute the transport for each point in a section
1117	!!
1118	!! Method :: Loop over each segment, and each vertical level and add the transport
1119	!! Be aware :
1120	!! One section is a sum of segments
1121	!! One segment is defined by 2 consecutive points in sec%listPoint
1122	!! All points of sec%listPoint are positioned on the F-point of the cell
1123	!!
1124	!! There are two loops:
1125	!! loop on the segment between 2 nodes
1126	!! loop on the level jk
1127	!!
1128	!! ** Output: Arrays containing the volume,density,salinity,temperature etc
1129	!! transports for each point in a section, summed over each nn_dct.
1130	!!
1131	!!-------------------------------------------------------------------------------------------
1132	!! * Arguments
1133	TYPE(SECTION),INTENT(INOUT) :: sec
1134	LOGICAL ,INTENT(IN) :: ld_debug
1135	INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section
1136
1137	!! * Local variables
1138	INTEGER :: jk,jseg,jclass, &!loop on level/segment/classes
1139	isgnu , isgnv !
1140	REAL(wp):: zumid , zvmid , &!U/V velocity on a cell segment
1141	zumid_ice , zvmid_ice , &!U/V ice velocity
1142	zTnorm !transport of velocity through one cell's sides
1143	REAL(wp):: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point
1144
1145	TYPE(POINT_SECTION) :: k
1146	!!--------------------------------------------------------
1147
1148	!!NIALL IF( ld_debug )WRITE(numout,*)' Compute transport'
1149
1150	!---------------------------!
1151	! COMPUTE TRANSPORT !
1152	!---------------------------!
1153	IF(sec%nb_point .NE. 0)THEN
1154
1155	!----------------------------------------------------------------------------------------------------
1156	!----------------------------------------------------------------------------------------------------
1157	!----------------------------------------------------------------------------------------------------
1158	!
1159	!
1160	! ! ! ! JT 1/09/2018 - changing convention. Always direction + is toward left hand of section
1161	!
1162	! Making sign of the velocities used to calculate the volume transport a function of direction, not slopesection
1163	! (isgnu, isgnv)
1164	!
1165	! They vary for each segment of the section.
1166	!
1167	!----------------------------------------------------------------------------------------------------
1168	!----------------------------------------------------------------------------------------------------
1169	!----------------------------------------------------------------------------------------------------
1170	!Compute sign for velocities:
1171	!
1172	!convention:
1173	! non horizontal section: direction + is toward left hand of section
1174	! horizontal section: direction + is toward north of section
1175	!
1176	!
1177	! slopeSection < 0 slopeSection > 0 slopeSection=inf slopeSection=0
1178	! ---------------- ----------------- --------------- --------------
1179	!
1180	! isgnv=1 direction +
1181	! ______ _____ ______
1182	! \| //\| \| \| direction +
1183	! \| isgnu=1 // \| \|isgnu=1 \|isgnu=1 /\|\
1184	! \|_______ // ______\| \\ \| ---\ \|
1185	! \| \| isgnv=-1 \\ \| \| ---/ direction + ____________
1186	! \| \| __\\\| \|
1187	! \| \| direction + \| isgnv=1
1188	!
1189	!----------------------------------------------------------------------------------------------------
1190
1191	IF( ld_debug )write(numout,*)"isgnu isgnv ",isgnu,isgnv
1192
1193	!--------------------------------------!
1194	! LOOP ON THE SEGMENT BETWEEN 2 NODES !
1195	!--------------------------------------!
1196	DO jseg=1,MAX(sec%nb_point-1,0)
1197
1198
1199	!Compute sign for velocities:
1200
1201	!isgnu = 1
1202	!isgnv = 1
1203	!
1204	! changing sign of u and v is dependent on the direction of the section.
1205	!isgnu = 1
1206	!isgnv = 1
1207	!SELECT CASE( sec%direction(jseg) )
1208	!CASE(0) ; isgnv = -1
1209	!CASE(3) ; isgnu = -1
1210	!END SELECT
1211
1212
1213	SELECT CASE( sec%direction(jseg) )
1214	CASE(0)
1215	isgnu = 1
1216	isgnv = -1
1217	CASE(1)
1218	isgnu = 1
1219	isgnv = 1
1220	CASE(2)
1221	isgnu = 1
1222	isgnv = 1
1223	CASE(3)
1224	isgnu = -1
1225	isgnv = 1
1226	END SELECT
1227
1228	!-------------------------------------------------------------------------------------------
1229	! Select the appropriate coordinate for computing the velocity of the segment
1230	! Corrected by JT 01/09/2018 (#)
1231	!
1232	! CASE(0) Case (2)
1233	! ------- --------
1234	! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j)
1235	! F(i,j)---------#V(i,j)-------F(i+1,j) \|
1236	! --------> \|
1237	! \| \|
1238	! \| \|
1239	! Case (3) \| U(i,j)
1240	! -------- \| \|
1241	! V \|
1242	! listPoint(jseg+1) F(i,j+1) \|
1243	! \| \|
1244	! \| \|
1245	! \| listPoint(jseg+1) F(i,j-1)
1246	! ^ \|
1247	! \| \|
1248	! \| U(i,j+1)
1249	! \| \| Case(1)
1250	! \| \| ------
1251	! \|
1252	! \| listPoint(jseg+1) listPoint(jseg)
1253	! \| F(i-1,j)----------#V(i-1,j) ------#f(i,j)
1254	! listPoint(jseg) F(i,j) <-------
1255	!
1256	!-------------------------------------------------------------------------------------------
1257
1258	SELECT CASE( sec%direction(jseg) )
1259	CASE(0) ; k = sec%listPoint(jseg)
1260	CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J)
1261	CASE(2) ; k = sec%listPoint(jseg)
1262	CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1)
1263	END SELECT
1264
1265	!---------------------------\|
1266	! LOOP ON THE LEVEL \|
1267	!---------------------------\|
1268	!Sum of the transport on the vertical
1269	DO jk=1,mbathy(k%I,k%J)
1270
1271	! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point
1272	SELECT CASE( sec%direction(jseg) )
1273	CASE(0,1)
1274	ztn = interp(k%I,k%J,jk,'V',0)
1275	zsn = interp(k%I,k%J,jk,'V',1)
1276	zrhop = interp(k%I,k%J,jk,'V',2)
1277	zrhoi = interp(k%I,k%J,jk,'V',3)
1278	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) vmask(k%I,k%J,1)
1279	CASE(2,3)
1280	ztn = interp(k%I,k%J,jk,'U',0)
1281	zsn = interp(k%I,k%J,jk,'U',1)
1282	zrhop = interp(k%I,k%J,jk,'U',2)
1283	zrhoi = interp(k%I,k%J,jk,'U',3)
1284	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) umask(k%I,k%J,1)
1285	END SELECT
1286
1287	zfsdep= fsdept(k%I,k%J,jk)
1288
1289	!compute velocity with the correct direction
1290	SELECT CASE( sec%direction(jseg) )
1291	CASE(0,1)
1292	zumid=0.
1293	zvmid=isgnvvn(k%I,k%J,jk)vmask(k%I,k%J,jk)
1294	CASE(2,3)
1295	zumid=isgnuun(k%I,k%J,jk)umask(k%I,k%J,jk)
1296	zvmid=0.
1297	END SELECT
1298
1299	!zTnorm=transport through one cell;
1300	!velocity* cell's length * cell's thickness
1301	zTnorm=zumide2u(k%I,k%J) fse3u(k%I,k%J,jk)+ &
1302	zvmide1v(k%I,k%J) fse3v(k%I,k%J,jk)
1303
1304	#if ! defined key_vvl
1305	!add transport due to free surface
1306	IF( jk==1 )THEN
1307	zTnorm = zTnorm + zumid* e2u(k%I,k%J) * zsshn * umask(k%I,k%J,jk) + &
1308	zvmid* e1v(k%I,k%J) * zsshn * vmask(k%I,k%J,jk)
1309	ENDIF
1310	#endif
1311	!COMPUTE TRANSPORT
1312
1313	transports_3d_h(1,jsec,jseg,jk) = transports_3d_h(1,jsec,jseg,jk) + zTnorm
1314
1315	IF ( sec%llstrpond ) THEN
1316	transports_3d_h(2,jsec,jseg,jk) = transports_3d_h(2,jsec,jseg,jk) + zTnorm * zrhoi
1317	transports_3d_h(3,jsec,jseg,jk) = transports_3d_h(3,jsec,jseg,jk) + zTnorm * zrhop
1318	transports_3d_h(4,jsec,jseg,jk) = transports_3d_h(4,jsec,jseg,jk) + zTnorm * 4.e+3_wp * (ztn+273.15) * 1026._wp
1319	transports_3d_h(5,jsec,jseg,jk) = transports_3d_h(5,jsec,jseg,jk) + zTnorm * 0.001 * zsn * 1026._wp
1320	ENDIF
1321
1322	ENDDO !end of loop on the level
1323
1324	#if defined key_lim2 \|\| defined key_lim3
1325
1326	!ICE CASE
1327	!------------
1328	IF( sec%ll_ice_section )THEN
1329	SELECT CASE (sec%direction(jseg))
1330	CASE(0)
1331	zumid_ice = 0
1332	zvmid_ice = isgnv0.5(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1))
1333	CASE(1)
1334	zumid_ice = 0
1335	zvmid_ice = isgnv0.5(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1))
1336	CASE(2)
1337	zvmid_ice = 0
1338	zumid_ice = isgnu0.5(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1))
1339	CASE(3)
1340	zvmid_ice = 0
1341	zumid_ice = isgnu0.5(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1))
1342	END SELECT
1343
1344	zTnorm=zumid_icee2u(k%I,k%J)+zvmid_icee1v(k%I,k%J)
1345
1346	transports_2d_h(1,jsec,jseg) = transports_2d_h(1,jsec,jseg) + (zTnorm)* &
1347	(1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1348	*(hsnif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J) + &
1349	hicif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1350	+zice_vol_pos
1351	transports_2d_h(2,jsec,jseg) = transports_2d_h(2,jsec,jseg) + (zTnorm)* &
1352	(1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) &
1353	+zice_surf_pos
1354
1355	ENDIF !end of ice case
1356	#endif
1357
1358	ENDDO !end of loop on the segment
1359
1360	ENDIF !end of sec%nb_point =0 case
1361	!
1362	END SUBROUTINE transport_h
1363
1364	SUBROUTINE dia_dct_sum(sec,jsec)
1365	!!-------------------------------------------------------------
1366	!! Purpose: Average the transport over nn_dctwri time steps
1367	!! and sum over the density/salinity/temperature/depth classes
1368	!!
1369	!! Method: Sum over relevant grid cells to obtain values
1370	!! for each class
1371	!! There are several loops:
1372	!! loop on the segment between 2 nodes
1373	!! loop on the level jk
1374	!! loop on the density/temperature/salinity/level classes
1375	!! test on the density/temperature/salinity/level
1376	!!
1377	!! Note: Transport through a given section is equal to the sum of transports
1378	!! computed on each proc.
1379	!! On each proc,transport is equal to the sum of transport computed through
1380	!! segments linking each point of sec%listPoint with the next one.
1381	!!
1382	!!-------------------------------------------------------------
1383	!! * arguments
1384	TYPE(SECTION),INTENT(INOUT) :: sec
1385	INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section
1386
1387	TYPE(POINT_SECTION) :: k
1388	INTEGER :: jk,jseg,jclass ! dummy variables for looping on level/segment/classes
1389	REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point
1390	!!-------------------------------------------------------------
1391
1392
1393	!! Sum the relevant segments to obtain values for each class
1394	IF(sec%nb_point .NE. 0)THEN
1395
1396	!--------------------------------------!
1397	! LOOP ON THE SEGMENT BETWEEN 2 NODES !
1398	!--------------------------------------!
1399	DO jseg=1,MAX(sec%nb_point-1,0)
1400
1401	!-------------------------------------------------------------------------------------------
1402	! Select the appropriate coordinate for computing the velocity of the segment
1403	!
1404	! CASE(0) Case (2)
1405	! ------- --------
1406	! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j)
1407	! F(i,j)----------V(i+1,j)-------F(i+1,j) \|
1408	! \|
1409	! \|
1410	! \|
1411	! Case (3) U(i,j)
1412	! -------- \|
1413	! \|
1414	! listPoint(jseg+1) F(i,j+1) \|
1415	! \| \|
1416	! \| \|
1417	! \| listPoint(jseg+1) F(i,j-1)
1418	! \|
1419	! \|
1420	! U(i,j+1)
1421	! \| Case(1)
1422	! \| ------
1423	! \|
1424	! \| listPoint(jseg+1) listPoint(jseg)
1425	! \| F(i-1,j)-----------V(i,j) -------f(jseg)
1426	! listPoint(jseg) F(i,j)
1427	!
1428	!-------------------------------------------------------------------------------------------
1429
1430	SELECT CASE( sec%direction(jseg) )
1431	CASE(0) ; k = sec%listPoint(jseg)
1432	CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J)
1433	CASE(2) ; k = sec%listPoint(jseg)
1434	CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1)
1435	END SELECT
1436
1437	!---------------------------\|
1438	! LOOP ON THE LEVEL \|
1439	!---------------------------\|
1440	!Sum of the transport on the vertical
1441	DO jk=1,mbathy(k%I,k%J)
1442
1443	! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point
1444	SELECT CASE( sec%direction(jseg) )
1445	CASE(0,1)
1446	ztn = interp(k%I,k%J,jk,'V',0)
1447	zsn = interp(k%I,k%J,jk,'V',1)
1448	zrhop = interp(k%I,k%J,jk,'V',2)
1449	zrhoi = interp(k%I,k%J,jk,'V',3)
1450	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) vmask(k%I,k%J,1)
1451	CASE(2,3)
1452	ztn = interp(k%I,k%J,jk,'U',0)
1453	zsn = interp(k%I,k%J,jk,'U',1)
1454	zrhop = interp(k%I,k%J,jk,'U',2)
1455	zrhoi = interp(k%I,k%J,jk,'U',3)
1456	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) umask(k%I,k%J,1)
1457	END SELECT
1458
1459	zfsdep= fsdept(k%I,k%J,jk)
1460
1461	!-------------------------------
1462	! LOOP ON THE DENSITY CLASSES \|
1463	!-------------------------------
1464	!The computation is made for each density/temperature/salinity/depth class
1465	DO jclass=1,MAX(1,sec%nb_class-1)
1466
1467	!----------------------------------------------!
1468	!TEST ON THE DENSITY/SALINITY/TEMPERATURE/LEVEL!
1469	!----------------------------------------------!
1470
1471	IF ( ( &
1472	((( zrhop .GE. (sec%zsigp(jclass)+1000. )) .AND. &
1473	( zrhop .LE. (sec%zsigp(jclass+1)+1000. ))) .OR. &
1474	( sec%zsigp(jclass) .EQ. 99.)) .AND. &
1475
1476	((( zrhoi .GE. (sec%zsigi(jclass) + 1000. )) .AND. &
1477	( zrhoi .LE. (sec%zsigi(jclass+1)+1000. ))) .OR. &
1478	( sec%zsigi(jclass) .EQ. 99.)) .AND. &
1479
1480	((( zsn .GT. sec%zsal(jclass)) .AND. &
1481	( zsn .LE. sec%zsal(jclass+1))) .OR. &
1482	( sec%zsal(jclass) .EQ. 99.)) .AND. &
1483
1484	((( ztn .GE. sec%ztem(jclass)) .AND. &
1485	( ztn .LE. sec%ztem(jclass+1))) .OR. &
1486	( sec%ztem(jclass) .EQ.99.)) .AND. &
1487
1488	((( zfsdep .GE. sec%zlay(jclass)) .AND. &
1489	( zfsdep .LE. sec%zlay(jclass+1))) .OR. &
1490	( sec%zlay(jclass) .EQ. 99. )) &
1491	)) THEN
1492
1493	!SUM THE TRANSPORTS FOR EACH CLASSES FOR THE POSITIVE AND NEGATIVE DIRECTIONS
1494	!----------------------------------------------------------------------------
1495	IF (transports_3d(1,jsec,jseg,jk) .GE. 0.0) THEN
1496	sec%transport(1,jclass) = sec%transport(1,jclass)+transports_3d(1,jsec,jseg,jk)
1497	ELSE
1498	sec%transport(2,jclass) = sec%transport(2,jclass)+transports_3d(1,jsec,jseg,jk)
1499	ENDIF
1500	IF( sec%llstrpond )THEN
1501
1502	IF( transports_3d(1,jsec,jseg,jk) .NE. 0._wp ) THEN
1503
1504	IF (transports_3d(2,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk) .GE. 0.0 ) THEN
1505	sec%transport(3,jclass) = sec%transport(3,jclass)+transports_3d(2,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk)
1506	ELSE
1507	sec%transport(4,jclass) = sec%transport(4,jclass)+transports_3d(2,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk)
1508	ENDIF
1509
1510	IF ( transports_3d(3,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk) .GE. 0.0 ) THEN
1511	sec%transport(5,jclass) = sec%transport(5,jclass)+transports_3d(3,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk)
1512	ELSE
1513	sec%transport(6,jclass) = sec%transport(6,jclass)+transports_3d(3,jsec,jseg,jk)/transports_3d(1,jsec,jseg,jk)
1514	ENDIF
1515
1516	ENDIF
1517
1518	IF ( transports_3d(4,jsec,jseg,jk) .GE. 0.0 ) THEN
1519	sec%transport(7,jclass) = sec%transport(7,jclass)+transports_3d(4,jsec,jseg,jk)
1520	ELSE
1521	sec%transport(8,jclass) = sec%transport(8,jclass)+transports_3d(4,jsec,jseg,jk)
1522	ENDIF
1523
1524	IF ( transports_3d(5,jsec,jseg,jk) .GE. 0.0 ) THEN
1525	sec%transport( 9,jclass) = sec%transport( 9,jclass)+transports_3d(5,jsec,jseg,jk)
1526	ELSE
1527	sec%transport(10,jclass) = sec%transport(10,jclass)+transports_3d(5,jsec,jseg,jk)
1528	ENDIF
1529
1530	ELSE
1531	sec%transport( 3,jclass) = 0._wp
1532	sec%transport( 4,jclass) = 0._wp
1533	sec%transport( 5,jclass) = 0._wp
1534	sec%transport( 6,jclass) = 0._wp
1535	sec%transport( 7,jclass) = 0._wp
1536	sec%transport( 8,jclass) = 0._wp
1537	sec%transport( 9,jclass) = 0._wp
1538	sec%transport(10,jclass) = 0._wp
1539	ENDIF
1540
1541	ENDIF ! end of test if point is in class
1542
1543	ENDDO ! end of loop on the classes
1544
1545	ENDDO ! loop over jk
1546
1547	#if defined key_lim2 \|\| defined key_lim3
1548
1549	!ICE CASE
1550	IF( sec%ll_ice_section )THEN
1551
1552	IF ( transports_2d(1,jsec,jseg) .GE. 0.0 ) THEN
1553	sec%transport(11,1) = sec%transport(11,1)+transports_2d(1,jsec,jseg)
1554	ELSE
1555	sec%transport(12,1) = sec%transport(12,1)+transports_2d(1,jsec,jseg)
1556	ENDIF
1557
1558	IF ( transports_2d(3,jsec,jseg) .GE. 0.0 ) THEN
1559	sec%transport(13,1) = sec%transport(13,1)+transports_2d(2,jsec,jseg)
1560	ELSE
1561	sec%transport(14,1) = sec%transport(14,1)+transports_2d(2,jsec,jseg)
1562	ENDIF
1563
1564	ENDIF !end of ice case
1565	#endif
1566
1567	ENDDO !end of loop on the segment
1568
1569	ELSE !if sec%nb_point =0
1570	sec%transport(1:2,:)=0.
1571	IF (sec%llstrpond) sec%transport(3:10,:)=0.
1572	IF (sec%ll_ice_section) sec%transport( 11:14,:)=0.
1573	ENDIF !end of sec%nb_point =0 case
1574
1575	END SUBROUTINE dia_dct_sum
1576
1577	SUBROUTINE dia_dct_sum_h(sec,jsec)
1578	!!-------------------------------------------------------------
1579	!! Exactly as dia_dct_sum but for hourly files containing data summed at each time step
1580	!!
1581	!! Purpose: Average the transport over nn_dctwri time steps
1582	!! and sum over the density/salinity/temperature/depth classes
1583	!!
1584	!! Method:
1585	!! Sum over relevant grid cells to obtain values
1586	!! for each
1587	!! There are several loops:
1588	!! loop on the segment between 2 nodes
1589	!! loop on the level jk
1590	!! loop on the density/temperature/salinity/level classes
1591	!! test on the density/temperature/salinity/level
1592	!!
1593	!! ** Method :Transport through a given section is equal to the sum of transports
1594	!! computed on each proc.
1595	!! On each proc,transport is equal to the sum of transport computed through
1596	!! segments linking each point of sec%listPoint with the next one.
1597	!!
1598	!!-------------------------------------------------------------
1599	!! * arguments
1600	TYPE(SECTION),INTENT(INOUT) :: sec
1601	INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section
1602
1603	TYPE(POINT_SECTION) :: k
1604	INTEGER :: jk,jseg,jclass !loop on level/segment/classes
1605	REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point
1606	!!-------------------------------------------------------------
1607
1608	!! Sum the relevant segments to obtain values for each class
1609	IF(sec%nb_point .NE. 0)THEN
1610
1611	!--------------------------------------!
1612	! LOOP ON THE SEGMENT BETWEEN 2 NODES !
1613	!--------------------------------------!
1614	DO jseg=1,MAX(sec%nb_point-1,0)
1615
1616	!-------------------------------------------------------------------------------------------
1617	! Select the appropriate coordinate for computing the velocity of the segment
1618	!
1619	! CASE(0) Case (2)
1620	! ------- --------
1621	! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j)
1622	! F(i,j)----------V(i+1,j)-------F(i+1,j) \|
1623	! \|
1624	! \|
1625	! \|
1626	! Case (3) U(i,j)
1627	! -------- \|
1628	! \|
1629	! listPoint(jseg+1) F(i,j+1) \|
1630	! \| \|
1631	! \| \|
1632	! \| listPoint(jseg+1) F(i,j-1)
1633	! \|
1634	! \|
1635	! U(i,j+1)
1636	! \| Case(1)
1637	! \| ------
1638	! \|
1639	! \| listPoint(jseg+1) listPoint(jseg)
1640	! \| F(i-1,j)-----------V(i,j) -------f(jseg)
1641	! listPoint(jseg) F(i,j)
1642	!
1643	!-------------------------------------------------------------------------------------------
1644
1645	SELECT CASE( sec%direction(jseg) )
1646	CASE(0) ; k = sec%listPoint(jseg)
1647	CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J)
1648	CASE(2) ; k = sec%listPoint(jseg)
1649	CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1)
1650	END SELECT
1651
1652	!---------------------------\|
1653	! LOOP ON THE LEVEL \|
1654	!---------------------------\|
1655	!Sum of the transport on the vertical
1656	DO jk=1,mbathy(k%I,k%J)
1657
1658	! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point
1659	SELECT CASE( sec%direction(jseg) )
1660	CASE(0,1)
1661	ztn = interp(k%I,k%J,jk,'V',0)
1662	zsn = interp(k%I,k%J,jk,'V',1)
1663	zrhop = interp(k%I,k%J,jk,'V',2)
1664	zrhoi = interp(k%I,k%J,jk,'V',3)
1665	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) vmask(k%I,k%J,1)
1666	CASE(2,3)
1667	ztn = interp(k%I,k%J,jk,'U',0)
1668	zsn = interp(k%I,k%J,jk,'U',1)
1669	zrhop = interp(k%I,k%J,jk,'U',2)
1670	zrhoi = interp(k%I,k%J,jk,'U',3)
1671	zsshn = 0.5( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) umask(k%I,k%J,1)
1672	END SELECT
1673
1674	zfsdep= fsdept(k%I,k%J,jk)
1675
1676	!-------------------------------
1677	! LOOP ON THE DENSITY CLASSES \|
1678	!-------------------------------
1679	!The computation is made for each density/heat/salt/... class
1680	DO jclass=1,MAX(1,sec%nb_class-1)
1681
1682	!----------------------------------------------!
1683	!TEST ON THE DENSITY/SALINITY/TEMPERATURE/LEVEL!
1684	!----------------------------------------------!
1685
1686	IF ( ( &
1687	((( zrhop .GE. (sec%zsigp(jclass)+1000. )) .AND. &
1688	( zrhop .LE. (sec%zsigp(jclass+1)+1000. ))) .OR. &
1689	( sec%zsigp(jclass) .EQ. 99.)) .AND. &
1690
1691	((( zrhoi .GE. (sec%zsigi(jclass) + 1000. )) .AND. &
1692	( zrhoi .LE. (sec%zsigi(jclass+1)+1000. ))) .OR. &
1693	( sec%zsigi(jclass) .EQ. 99.)) .AND. &
1694
1695	((( zsn .GT. sec%zsal(jclass)) .AND. &
1696	( zsn .LE. sec%zsal(jclass+1))) .OR. &
1697	( sec%zsal(jclass) .EQ. 99.)) .AND. &
1698
1699	((( ztn .GE. sec%ztem(jclass)) .AND. &
1700	( ztn .LE. sec%ztem(jclass+1))) .OR. &
1701	( sec%ztem(jclass) .EQ.99.)) .AND. &
1702
1703	((( zfsdep .GE. sec%zlay(jclass)) .AND. &
1704	( zfsdep .LE. sec%zlay(jclass+1))) .OR. &
1705	( sec%zlay(jclass) .EQ. 99. )) &
1706	)) THEN
1707
1708	!SUM THE TRANSPORTS FOR EACH CLASSES FOR THE POSITIVE AND NEGATIVE DIRECTIONS
1709	!----------------------------------------------------------------------------
1710	IF (transports_3d_h(1,jsec,jseg,jk) .GE. 0.0) THEN
1711	sec%transport_h(1,jclass) = sec%transport_h(1,jclass)+transports_3d_h(1,jsec,jseg,jk)
1712	ELSE
1713	sec%transport_h(2,jclass) = sec%transport_h(2,jclass)+transports_3d_h(1,jsec,jseg,jk)
1714	ENDIF
1715	IF( sec%llstrpond )THEN
1716
1717	IF( transports_3d_h(1,jsec,jseg,jk) .NE. 0._wp ) THEN
1718
1719	IF (transports_3d_h(2,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk) .GE. 0.0 ) THEN
1720	sec%transport_h(3,jclass) = sec%transport_h(3,jclass)+transports_3d_h(2,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk)
1721	ELSE
1722	sec%transport_h(4,jclass) = sec%transport_h(4,jclass)+transports_3d_h(2,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk)
1723	ENDIF
1724
1725	IF ( transports_3d_h(3,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk) .GE. 0.0 ) THEN
1726	sec%transport_h(5,jclass) = sec%transport_h(5,jclass)+transports_3d_h(3,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk)
1727	ELSE
1728	sec%transport_h(6,jclass) = sec%transport_h(6,jclass)+transports_3d_h(3,jsec,jseg,jk)/transports_3d_h(1,jsec,jseg,jk)
1729	ENDIF
1730
1731	ENDIF
1732
1733	IF ( transports_3d_h(4,jsec,jseg,jk) .GE. 0.0 ) THEN
1734	sec%transport_h(7,jclass) = sec%transport_h(7,jclass)+transports_3d_h(4,jsec,jseg,jk)
1735	ELSE
1736	sec%transport_h(8,jclass) = sec%transport_h(8,jclass)+transports_3d_h(4,jsec,jseg,jk)
1737	ENDIF
1738
1739	IF ( transports_3d_h(5,jsec,jseg,jk) .GE. 0.0 ) THEN
1740	sec%transport_h( 9,jclass) = sec%transport_h( 9,jclass)+transports_3d_h(5,jsec,jseg,jk)
1741	ELSE
1742	sec%transport_h(10,jclass) = sec%transport_h(10,jclass)+transports_3d_h(5,jsec,jseg,jk)
1743	ENDIF
1744
1745	ELSE
1746	sec%transport_h( 3,jclass) = 0._wp
1747	sec%transport_h( 4,jclass) = 0._wp
1748	sec%transport_h( 5,jclass) = 0._wp
1749	sec%transport_h( 6,jclass) = 0._wp
1750	sec%transport_h( 7,jclass) = 0._wp
1751	sec%transport_h( 8,jclass) = 0._wp
1752	sec%transport_h( 9,jclass) = 0._wp
1753	sec%transport_h(10,jclass) = 0._wp
1754	ENDIF
1755
1756	ENDIF ! end of test if point is in class
1757
1758	ENDDO ! end of loop on the classes
1759
1760	ENDDO ! loop over jk
1761
1762	#if defined key_lim2 \|\| defined key_lim3
1763
1764	!ICE CASE
1765	IF( sec%ll_ice_section )THEN
1766
1767	IF ( transports_2d_h(1,jsec,jseg) .GE. 0.0 ) THEN
1768	sec%transport_h(11,1) = sec%transport_h(11,1)+transports_2d_h(1,jsec,jseg)
1769	ELSE
1770	sec%transport_h(12,1) = sec%transport_h(12,1)+transports_2d_h(1,jsec,jseg)
1771	ENDIF
1772
1773	IF ( transports_2d_h(3,jsec,jseg) .GE. 0.0 ) THEN
1774	sec%transport_h(13,1) = sec%transport_h(13,1)+transports_2d_h(2,jsec,jseg)
1775	ELSE
1776	sec%transport_h(14,1) = sec%transport_h(14,1)+transports_2d_h(2,jsec,jseg)
1777	ENDIF
1778
1779	ENDIF !end of ice case
1780	#endif
1781
1782	ENDDO !end of loop on the segment
1783
1784	ELSE !if sec%nb_point =0
1785	sec%transport_h(1:2,:)=0.
1786	IF (sec%llstrpond) sec%transport_h(3:10,:)=0.
1787	IF (sec%ll_ice_section) sec%transport_h( 11:14,:)=0.
1788	ENDIF !end of sec%nb_point =0 case
1789
1790	END SUBROUTINE dia_dct_sum_h
1791
1792	SUBROUTINE dia_dct_wri_NOOS(kt,ksec,sec)
1793	!!-------------------------------------------------------------
1794	!! Write transport output in numdct using NOOS formatting
1795	!!
1796	!! Purpose: Write transports in ascii files
1797	!!
1798	!! Method:
1799	!! 1. Write volume transports in "volume_transport"
1800	!! Unit: Sv : area * Velocity / 1.e6
1801	!!
1802	!! 2. Write heat transports in "heat_transport"
1803	!! Unit: Peta W : area * Velocity * T * rhau * Cp / 1.e15
1804	!!
1805	!! 3. Write salt transports in "salt_transport"
1806	!! Unit: 10^9 g m^3 / s : area * Velocity * S / 1.e6
1807	!!
1808	!!-------------------------------------------------------------
1809	!!arguments
1810	INTEGER, INTENT(IN) :: kt ! time-step
1811	TYPE(SECTION), INTENT(INOUT) :: sec ! section to write
1812	INTEGER ,INTENT(IN) :: ksec ! section number
1813
1814	!!local declarations
1815	INTEGER :: jclass,ji ! Dummy loop
1816	CHARACTER(len=2) :: classe ! Classname
1817	REAL(wp) :: zbnd1,zbnd2 ! Class bounds
1818	REAL(wp) :: zslope ! section's slope coeff
1819	!
1820	REAL(wp), POINTER, DIMENSION(:):: zsumclasses ! 1D workspace
1821	CHARACTER(len=3) :: noos_sect_name ! Classname
1822	CHARACTER(len=25) :: noos_var_sect_name
1823	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: noos_iom_dummy
1824	INTEGER :: IERR
1825
1826	REAL(wp), DIMENSION(3) :: tmp_iom_output
1827	REAL(wp) :: max_iom_val
1828
1829	!!-------------------------------------------------------------
1830
1831
1832
1833	IF( lwp ) THEN
1834	WRITE(numout,*) " "
1835	WRITE(numout,*) "dia_dct_wri_NOOS: write transports through sections at timestep: ", kt
1836	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
1837	ENDIF
1838
1839	CALL wrk_alloc(nb_type_class , zsumclasses )
1840
1841	zsumclasses(:)=0._wp
1842	zslope = sec%slopeSection
1843
1844	IF( lwp ) THEN
1845	IF ( ln_dct_ascii ) THEN
1846	WRITE(numdct_NOOS,'(I4,a1,I2,a1,I2,a12,i3,a17,i3,a10,a25)') nyear,'.',nmonth,'.',nday,' Transect:',ksec-1,' No. of layers:',sec%nb_class-1,' Name: ',sec%name
1847	ELSE
1848	WRITE(numdct_NOOS) nyear,nmonth,nday,ksec-1,sec%nb_class-1,sec%name
1849	ENDIF
1850	ENDIF
1851
1852	! Sum all classes together, to give one values per type (pos tran, neg vol trans etc...).
1853	DO jclass=1,MAX(1,sec%nb_class-1)
1854	zsumclasses(1:nb_type_class)=zsumclasses(1:nb_type_class)+sec%transport(1:nb_type_class,jclass)
1855	ENDDO
1856
1857	classe = 'total '
1858	zbnd1 = 0._wp
1859	zbnd2 = 0._wp
1860
1861
1862
1863	write (noos_sect_name, "(I0.3)") ksec
1864
1865	IF ( ln_dct_iom_cont ) THEN
1866	max_iom_val = 1.e10
1867	ALLOCATE( noos_iom_dummy(jpi,jpj,3), STAT= ierr )
1868	IF( ierr /= 0 ) CALL ctl_stop( 'dia_dct_wri_NOOS: failed to allocate noos_iom_dummy array' )
1869	ENDIF
1870
1871	! JT I think changing the sign on the output based on the zslope value is redunant.
1872	!
1873	! IF ( zslope .gt. 0._wp .and. zslope .ne. 10000._wp ) THEN
1874	!
1875	! IF( lwp ) THEN
1876	! WRITE(numdct_NOOS,'(9e12.4E2)') -(zsumclasses( 1)+zsumclasses( 2)), -zsumclasses( 2),-zsumclasses( 1), &
1877	! -(zsumclasses( 7)+zsumclasses( 8)), -zsumclasses( 8),-zsumclasses( 7), &
1878	! -(zsumclasses( 9)+zsumclasses(10)), -zsumclasses(10),-zsumclasses( 9)
1879	! CALL FLUSH(numdct_NOOS)
1880	! endif
1881
1882	!
1883	! IF ( ln_dct_iom_cont ) THEN
1884	!
1885	! noos_iom_dummy(:,:,:) = 0.
1886	!
1887	! tmp_iom_output(:) = 0.
1888	! tmp_iom_output(1) = -(zsumclasses( 1)+zsumclasses( 2))
1889	! tmp_iom_output(2) = -zsumclasses( 2)
1890	! tmp_iom_output(3) = -zsumclasses( 1)
1891	!
1892	! ! Convert to Sv
1893	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
1894	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
1895	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
1896	!
1897	! ! limit maximum and minimum values in iom_put
1898	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
1899	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
1900	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
1901	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
1902	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
1903	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
1904	!
1905	! ! Set NaN's to Zero
1906	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
1907	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
1908	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
1909	!
1910	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
1911	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
1912	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
1913	!
1914	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_trans'
1915	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
1916	! CALL iom_put( noos_var_sect_name, noos_iom_dummy )
1917	! noos_iom_dummy(:,:,:) = 0.
1918	! tmp_iom_output(:) = 0.
1919	! tmp_iom_output(1) = -(zsumclasses( 7)+zsumclasses( 8))
1920	! tmp_iom_output(2) = -zsumclasses( 8)
1921	! tmp_iom_output(3) = -zsumclasses( 7)
1922	!
1923	! ! Convert to TJ/s
1924	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-12
1925	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-12
1926	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-12
1927	!
1928	! ! limit maximum and minimum values in iom_put
1929	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
1930	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
1931	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
1932	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
1933	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
1934	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
1935	!
1936	! ! Set NaN's to Zero
1937	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
1938	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
1939	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
1940	!
1941	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
1942	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
1943	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
1944	!
1945	! !noos_iom_dummy(:,:,1) = -(zsumclasses( 7)+zsumclasses( 8))
1946	! !noos_iom_dummy(:,:,2) = -zsumclasses( 8)
1947	! !noos_iom_dummy(:,:,3) = -zsumclasses( 7)
1948	!
1949	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_heat'
1950	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
1951	! CALL iom_put( noos_var_sect_name, noos_iom_dummy )
1952	!
1953	! noos_iom_dummy(:,:,:) = 0.
1954	! tmp_iom_output(:) = 0.
1955	! tmp_iom_output(1) = -(zsumclasses( 9)+zsumclasses( 10))
1956	! tmp_iom_output(2) = -zsumclasses( 10)
1957	! tmp_iom_output(3) = -zsumclasses( 9)
1958	!
1959	! ! Convert to MT/s
1960	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
1961	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
1962	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
1963	!
1964	! ! limit maximum and minimum values in iom_put
1965	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
1966	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
1967	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
1968	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
1969	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
1970	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
1971	!
1972	! ! Set NaN's to Zero
1973	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
1974	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
1975	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
1976	!
1977	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
1978	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
1979	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
1980	!
1981	! !noos_iom_dummy(:,:,1) = -(zsumclasses( 9)+zsumclasses( 10))
1982	! !noos_iom_dummy(:,:,2) = -zsumclasses( 10)
1983	! !noos_iom_dummy(:,:,3) = -zsumclasses( 9)
1984	!
1985	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_salt'
1986	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
1987	! CALL iom_put( noos_var_sect_name, noos_iom_dummy )
1988	! noos_iom_dummy(:,:,:) = 0.
1989	! tmp_iom_output(:) = 0.
1990	! ENDIF
1991	! ELSE
1992	! IF( lwp ) THEN
1993	! WRITE(numdct_NOOS,'(9e12.4E2)') zsumclasses( 1)+zsumclasses( 2) , zsumclasses( 1), zsumclasses( 2), &
1994	! zsumclasses( 7)+zsumclasses( 8) , zsumclasses( 7), zsumclasses( 8), &
1995	! zsumclasses( 9)+zsumclasses(10) , zsumclasses( 9), zsumclasses(10)
1996	! CALL FLUSH(numdct_NOOS)
1997	! endif
1998	!
1999	!
2000	! IF ( ln_dct_iom_cont ) THEN
2001	!
2002	! noos_iom_dummy(:,:,:) = 0.
2003	! tmp_iom_output(:) = 0.
2004	!
2005	! tmp_iom_output(1) = (zsumclasses( 1)+zsumclasses( 2))
2006	! tmp_iom_output(2) = zsumclasses( 1)
2007	! tmp_iom_output(3) = zsumclasses( 2)
2008	!
2009	! ! Convert to Sv
2010	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
2011	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
2012	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
2013	!
2014	! ! limit maximum and minimum values in iom_put
2015	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2016	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2017	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2018	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2019	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2020	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2021	!
2022	! ! Set NaN's to Zero
2023	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2024	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2025	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2026	!
2027	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2028	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2029	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2030	!
2031	! !noos_iom_dummy(:,:,1) = (zsumclasses( 1)+zsumclasses( 2))
2032	! !noos_iom_dummy(:,:,2) = zsumclasses( 1)
2033	! !noos_iom_dummy(:,:,3) = zsumclasses( 2)
2034	!
2035	!
2036	!
2037	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_trans'
2038	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2039	! CALL iom_put( noos_var_sect_name, noos_iom_dummy )
2040	! noos_iom_dummy(:,:,:) = 0.
2041	! tmp_iom_output(:) = 0.
2042	!
2043	! tmp_iom_output(1) = (zsumclasses( 7)+zsumclasses( 8))
2044	! tmp_iom_output(2) = zsumclasses( 7)
2045	! tmp_iom_output(3) = zsumclasses( 8)
2046	!
2047	! ! Convert to TJ/s
2048	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-12
2049	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-12
2050	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-12
2051	!
2052	! ! limit maximum and minimum values in iom_put
2053	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2054	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2055	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2056	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2057	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2058	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2059	!
2060	! ! Set NaN's to Zero
2061	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2062	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2063	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2064	!
2065	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2066	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2067	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2068	!
2069	! !noos_iom_dummy(:,:,1) = (zsumclasses( 7)+zsumclasses( 8))
2070	! !noos_iom_dummy(:,:,2) = zsumclasses( 7)
2071	! !noos_iom_dummy(:,:,3) = zsumclasses( 8)
2072	!
2073	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_heat'
2074	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2075	! CALL iom_put(noos_var_sect_name, noos_iom_dummy )
2076	! noos_iom_dummy(:,:,:) = 0.
2077	! tmp_iom_output(:) = 0.
2078	!
2079	! tmp_iom_output(1) = (zsumclasses( 9)+zsumclasses( 10))
2080	! tmp_iom_output(2) = zsumclasses( 9)
2081	! tmp_iom_output(3) = zsumclasses( 10)
2082	!
2083	! ! Convert to MT/s
2084	! tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
2085	! tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
2086	! tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
2087	!
2088	!
2089	! ! limit maximum and minimum values in iom_put
2090	! if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2091	! if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2092	! if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2093	! if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2094	! if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2095	! if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2096	!
2097	! ! Set NaN's to Zero
2098	! if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2099	! if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2100	! if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2101	!
2102	! noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2103	! noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2104	! noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2105	!
2106	! !noos_iom_dummy(:,:,1) = (zsumclasses( 9)+zsumclasses( 10))
2107	! !noos_iom_dummy(:,:,2) = zsumclasses( 9)
2108	! !noos_iom_dummy(:,:,3) = zsumclasses( 10)
2109	!
2110	! noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_salt'
2111	! if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2112	! CALL iom_put(noos_var_sect_name, noos_iom_dummy )
2113	! noos_iom_dummy(:,:,:) = 0.
2114	! tmp_iom_output(:) = 0.
2115	! ENDIF
2116	!
2117	! ENDIF
2118
2119
2120
2121
2122
2123
2124
2125
2126
2127
2128	IF( lwp ) THEN
2129	IF ( ln_dct_ascii ) THEN
2130	!WRITE(numdct_NOOS,'(9e12.4E2)') zsumclasses( 1)+zsumclasses( 2) , zsumclasses( 1), zsumclasses( 2), &
2131	WRITE(numdct_NOOS,'(3F18.3,6e16.8E2)') zsumclasses( 1)+zsumclasses( 2) , zsumclasses( 1), zsumclasses( 2), &
2132	zsumclasses( 7)+zsumclasses( 8) , zsumclasses( 7), zsumclasses( 8), &
2133	zsumclasses( 9)+zsumclasses(10) , zsumclasses( 9), zsumclasses(10)
2134	CALL FLUSH(numdct_NOOS)
2135	ELSE
2136	WRITE(numdct_NOOS) zsumclasses( 1)+zsumclasses( 2) , zsumclasses( 1), zsumclasses( 2), &
2137	zsumclasses( 7)+zsumclasses( 8) , zsumclasses( 7), zsumclasses( 8), &
2138	zsumclasses( 9)+zsumclasses(10) , zsumclasses( 9), zsumclasses(10)
2139	CALL FLUSH(numdct_NOOS)
2140	ENDIF
2141	ENDIF
2142
2143	IF ( ln_dct_iom_cont ) THEN
2144
2145	noos_iom_dummy(:,:,:) = 0.
2146	tmp_iom_output(:) = 0.
2147
2148	tmp_iom_output(1) = (zsumclasses( 1)+zsumclasses( 2))
2149	tmp_iom_output(2) = zsumclasses( 1)
2150	tmp_iom_output(3) = zsumclasses( 2)
2151
2152	! Convert to Sv
2153	tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
2154	tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
2155	tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
2156
2157	! limit maximum and minimum values in iom_put
2158	if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2159	if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2160	if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2161	if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2162	if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2163	if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2164
2165	! Set NaN's to Zero
2166	if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2167	if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2168	if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2169
2170	noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2171	noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2172	noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2173
2174	!noos_iom_dummy(:,:,1) = (zsumclasses( 1)+zsumclasses( 2))
2175	!noos_iom_dummy(:,:,2) = zsumclasses( 1)
2176	!noos_iom_dummy(:,:,3) = zsumclasses( 2)
2177
2178
2179
2180	noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_trans'
2181	if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2182	CALL iom_put( noos_var_sect_name, noos_iom_dummy )
2183	noos_iom_dummy(:,:,:) = 0.
2184	tmp_iom_output(:) = 0.
2185
2186	tmp_iom_output(1) = (zsumclasses( 7)+zsumclasses( 8))
2187	tmp_iom_output(2) = zsumclasses( 7)
2188	tmp_iom_output(3) = zsumclasses( 8)
2189
2190	! Convert to TJ/s
2191	tmp_iom_output(1) = tmp_iom_output(1)*1.E-12
2192	tmp_iom_output(2) = tmp_iom_output(2)*1.E-12
2193	tmp_iom_output(3) = tmp_iom_output(3)*1.E-12
2194
2195	! limit maximum and minimum values in iom_put
2196	if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2197	if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2198	if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2199	if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2200	if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2201	if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2202
2203	! Set NaN's to Zero
2204	if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2205	if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2206	if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2207
2208	noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2209	noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2210	noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2211
2212	!noos_iom_dummy(:,:,1) = (zsumclasses( 7)+zsumclasses( 8))
2213	!noos_iom_dummy(:,:,2) = zsumclasses( 7)
2214	!noos_iom_dummy(:,:,3) = zsumclasses( 8)
2215
2216	noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_heat'
2217	if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2218	CALL iom_put(noos_var_sect_name, noos_iom_dummy )
2219	noos_iom_dummy(:,:,:) = 0.
2220	tmp_iom_output(:) = 0.
2221
2222	tmp_iom_output(1) = (zsumclasses( 9)+zsumclasses( 10))
2223	tmp_iom_output(2) = zsumclasses( 9)
2224	tmp_iom_output(3) = zsumclasses( 10)
2225
2226	! Convert to MT/s
2227	tmp_iom_output(1) = tmp_iom_output(1)*1.E-6
2228	tmp_iom_output(2) = tmp_iom_output(2)*1.E-6
2229	tmp_iom_output(3) = tmp_iom_output(3)*1.E-6
2230
2231
2232	! limit maximum and minimum values in iom_put
2233	if ( tmp_iom_output(1) .gt. max_iom_val ) tmp_iom_output(1) = max_iom_val
2234	if ( tmp_iom_output(1) .lt. -max_iom_val ) tmp_iom_output(1) = -max_iom_val
2235	if ( tmp_iom_output(2) .gt. max_iom_val ) tmp_iom_output(2) = max_iom_val
2236	if ( tmp_iom_output(2) .lt. -max_iom_val ) tmp_iom_output(2) = -max_iom_val
2237	if ( tmp_iom_output(3) .gt. max_iom_val ) tmp_iom_output(3) = max_iom_val
2238	if ( tmp_iom_output(3) .lt. -max_iom_val ) tmp_iom_output(3) = -max_iom_val
2239
2240	! Set NaN's to Zero
2241	if ( tmp_iom_output(1) .ne. tmp_iom_output(1) ) tmp_iom_output(1) = max_iom_val*2
2242	if ( tmp_iom_output(2) .ne. tmp_iom_output(2) ) tmp_iom_output(1) = max_iom_val*2
2243	if ( tmp_iom_output(3) .ne. tmp_iom_output(3) ) tmp_iom_output(1) = max_iom_val*2
2244
2245	noos_iom_dummy(:,:,1) = tmp_iom_output(1)
2246	noos_iom_dummy(:,:,2) = tmp_iom_output(2)
2247	noos_iom_dummy(:,:,3) = tmp_iom_output(3)
2248
2249	!noos_iom_dummy(:,:,1) = (zsumclasses( 9)+zsumclasses( 10))
2250	!noos_iom_dummy(:,:,2) = zsumclasses( 9)
2251	!noos_iom_dummy(:,:,3) = zsumclasses( 10)
2252
2253	noos_var_sect_name = "noos_" // trim(noos_sect_name) // '_salt'
2254	if ( lwp ) WRITE(numout,*) 'dia_dct_wri_NOOS iom_put: ', kt,ksec, noos_var_sect_name
2255	CALL iom_put(noos_var_sect_name, noos_iom_dummy )
2256	noos_iom_dummy(:,:,:) = 0.
2257	tmp_iom_output(:) = 0.
2258
2259
2260	DEALLOCATE(noos_iom_dummy)
2261	ENDIF
2262
2263
2264	DO jclass=1,MAX(1,sec%nb_class-1)
2265
2266	classe = 'N '
2267	zbnd1 = 0._wp
2268	zbnd2 = 0._wp
2269
2270	!insitu density classes transports
2271	IF( ( sec%zsigi(jclass) .NE. 99._wp ) .AND. &
2272	( sec%zsigi(jclass+1) .NE. 99._wp ) )THEN
2273	classe = 'DI '
2274	zbnd1 = sec%zsigi(jclass)
2275	zbnd2 = sec%zsigi(jclass+1)
2276	ENDIF
2277	!potential density classes transports
2278	IF( ( sec%zsigp(jclass) .NE. 99._wp ) .AND. &
2279	( sec%zsigp(jclass+1) .NE. 99._wp ) )THEN
2280	classe = 'DP '
2281	zbnd1 = sec%zsigp(jclass)
2282	zbnd2 = sec%zsigp(jclass+1)
2283	ENDIF
2284	!depth classes transports
2285	IF( ( sec%zlay(jclass) .NE. 99._wp ) .AND. &
2286	( sec%zlay(jclass+1) .NE. 99._wp ) )THEN
2287	classe = 'Z '
2288	zbnd1 = sec%zlay(jclass)
2289	zbnd2 = sec%zlay(jclass+1)
2290	ENDIF
2291	!salinity classes transports
2292	IF( ( sec%zsal(jclass) .NE. 99._wp ) .AND. &
2293	( sec%zsal(jclass+1) .NE. 99._wp ) )THEN
2294	classe = 'S '
2295	zbnd1 = sec%zsal(jclass)
2296	zbnd2 = sec%zsal(jclass+1)
2297	ENDIF
2298	!temperature classes transports
2299	IF( ( sec%ztem(jclass) .NE. 99._wp ) .AND. &
2300	( sec%ztem(jclass+1) .NE. 99._wp ) ) THEN
2301	classe = 'T '
2302	zbnd1 = sec%ztem(jclass)
2303	zbnd2 = sec%ztem(jclass+1)
2304	ENDIF
2305
2306	!write volume transport per class
2307	IF( lwp ) THEN
2308
2309	IF ( ln_dct_ascii ) THEN
2310	CALL FLUSH(numdct_NOOS)
2311
2312	!WRITE(numdct_NOOS,'(9e12.4E2)') sec%transport( 1,jclass)+sec%transport( 2,jclass) , sec%transport( 1,jclass), sec%transport( 2,jclass), &
2313	! sec%transport( 7,jclass)+sec%transport( 8,jclass) , sec%transport( 7,jclass), sec%transport( 8,jclass), &
2314	! sec%transport( 9,jclass)+sec%transport(10,jclass) , sec%transport( 9,jclass), sec%transport(10,jclass)
2315	WRITE(numdct_NOOS,'(3F18.3,6e16.8E2)') sec%transport( 1,jclass)+sec%transport( 2,jclass) , sec%transport( 1,jclass), sec%transport( 2,jclass), &
2316	sec%transport( 7,jclass)+sec%transport( 8,jclass) , sec%transport( 7,jclass), sec%transport( 8,jclass), &
2317	sec%transport( 9,jclass)+sec%transport(10,jclass) , sec%transport( 9,jclass), sec%transport(10,jclass)
2318	ELSE
2319
2320	CALL FLUSH(numdct_NOOS)
2321	WRITE(numdct_NOOS) sec%transport( 1,jclass)+sec%transport( 2,jclass) , sec%transport( 1,jclass), sec%transport( 2,jclass), &
2322	sec%transport( 7,jclass)+sec%transport( 8,jclass) , sec%transport( 7,jclass), sec%transport( 8,jclass), &
2323	sec%transport( 9,jclass)+sec%transport(10,jclass) , sec%transport( 9,jclass), sec%transport(10,jclass)
2324	ENDIF
2325	ENDIF
2326
2327	ENDDO
2328
2329	!IF ( ln_dct_ascii ) THEN
2330	if ( lwp ) CALL FLUSH(numdct_NOOS)
2331	!ENDIF
2332
2333	CALL wrk_dealloc(nb_type_class , zsumclasses )
2334
2335	END SUBROUTINE dia_dct_wri_NOOS
2336
2337	SUBROUTINE dia_dct_wri_NOOS_h(hr,ksec,sec)
2338	!!-------------------------------------------------------------
2339	!! As routine dia_dct_wri_NOOS but for hourly output files
2340	!!
2341	!! Write transport output in numdct using NOOS formatting
2342	!!
2343	!! Purpose: Write transports in ascii files
2344	!!
2345	!! Method:
2346	!! 1. Write volume transports in "volume_transport"
2347	!! Unit: Sv : area * Velocity / 1.e6
2348	!!
2349	!!-------------------------------------------------------------
2350	!!arguments
2351	INTEGER, INTENT(IN) :: hr ! hour => effectively kt/12
2352	TYPE(SECTION), INTENT(INOUT) :: sec ! section to write
2353	INTEGER ,INTENT(IN) :: ksec ! section number
2354
2355	!!local declarations
2356	INTEGER :: jclass,jhr ! Dummy loop
2357	CHARACTER(len=2) :: classe ! Classname
2358	REAL(wp) :: zbnd1,zbnd2 ! Class bounds
2359	REAL(wp) :: zslope ! section's slope coeff
2360	!
2361	REAL(wp), POINTER, DIMENSION(:):: zsumclasses ! 1D workspace
2362	CHARACTER(len=3) :: noos_sect_name ! Classname
2363	CHARACTER(len=25) :: noos_var_sect_name
2364	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: noos_iom_dummy
2365	INTEGER :: IERR
2366
2367	!!-------------------------------------------------------------
2368
2369	IF( lwp ) THEN
2370	WRITE(numout,*) " "
2371	WRITE(numout,*) "dia_dct_wri_NOOS_h: write transports through section Transect:",ksec-1," at timestep: ", hr
2372	WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~"
2373	ENDIF
2374
2375	CALL wrk_alloc(nb_type_class , zsumclasses )
2376
2377
2378	write (noos_sect_name, "(I03)") ksec
2379
2380	ALLOCATE( noos_iom_dummy(jpi,jpj,3), STAT= ierr )
2381	IF( ierr /= 0 ) CALL ctl_stop( 'dia_dct_wri_NOOS_h: failed to allocate noos_iom_dummy array' )
2382
2383
2384
2385
2386
2387	zsumclasses(:)=0._wp
2388	zslope = sec%slopeSection
2389
2390	! Sum up all classes, to give the total per type (pos vol trans, neg vol trans etc...)
2391	DO jclass=1,MAX(1,sec%nb_class-1)
2392	zsumclasses(1:nb_type_class)=zsumclasses(1:nb_type_class)+sec%transport_h(1:nb_type_class,jclass)
2393	ENDDO
2394
2395
2396	! JT I think changing the sign of output according to the zslope is redundant
2397
2398	!write volume transport per class
2399	! Sum positive and vol trans for all classes in first cell of array
2400
2401	z_hr_output(ksec,1,1)= (zsumclasses(1)+zsumclasses(2))
2402	z_hr_output(ksec,2,1)= zsumclasses(1)
2403	z_hr_output(ksec,3,1)= zsumclasses(2)
2404
2405	! Sum positive and vol trans for each classes in following cell of array
2406	DO jclass=1,MAX(1,sec%nb_class-1)
2407	z_hr_output(ksec,1,jclass+1)= (sec%transport_h(1,jclass)+sec%transport_h(2,jclass))
2408	z_hr_output(ksec,2,jclass+1)= sec%transport_h(1,jclass)
2409	z_hr_output(ksec,3,jclass+1)= sec%transport_h(2,jclass)
2410	ENDDO
2411
2412
2413	IF( lwp ) THEN
2414	! JT IF ( hr .eq. 48._wp ) THEN
2415	! JT WRITE(numdct_NOOS_h,'(I4,a1,I2,a1,I2,a12,i3,a17,i3)') nyear,'.',nmonth,'.',nday,' Transect:',ksec-1,' No. of layers:',sec%nb_class-1
2416	! JT DO jhr=25,48
2417	! JT WRITE(numdct_NOOS_h,'(11F12.1)') z_hr_output(ksec,jhr,1), (z_hr_output(ksec,jhr,jclass+1), jclass=1,MAX(1,10) )
2418	! JT ENDDO
2419	! JT ENDIF
2420
2421
2422
2423	IF ( ln_dct_ascii ) THEN
2424	WRITE(numdct_NOOS_h,'(I4,a1,I2,a1,I2,a1,I2,a1,I2,a12,i3,a17,i3)') nyear,'.',nmonth,'.',nday,'.',MOD(hr,24),'.',0,' Transect:',ksec-1,' No. of layers:',sec%nb_class-1
2425	WRITE(numdct_NOOS_h,'(11F18.3)') z_hr_output(ksec,1,1), (z_hr_output(ksec,1,jclass+1), jclass=1,MAX(1,10) )
2426	WRITE(numdct_NOOS_h,'(11F18.3)') z_hr_output(ksec,2,1), (z_hr_output(ksec,2,jclass+1), jclass=1,MAX(1,10) )
2427	WRITE(numdct_NOOS_h,'(11F18.3)') z_hr_output(ksec,3,1), (z_hr_output(ksec,3,jclass+1), jclass=1,MAX(1,10) )
2428	CALL FLUSH(numdct_NOOS_h)
2429	ELSE
2430	WRITE(numdct_NOOS_h) nyear,nmonth,nday,MOD(hr,24),ksec-1,sec%nb_class-1
2431	WRITE(numdct_NOOS_h) z_hr_output(ksec,1,1), (z_hr_output(ksec,1,jclass+1), jclass=1,MAX(1,10) )
2432	WRITE(numdct_NOOS_h) z_hr_output(ksec,2,1), (z_hr_output(ksec,2,jclass+1), jclass=1,MAX(1,10) )
2433	WRITE(numdct_NOOS_h) z_hr_output(ksec,3,1), (z_hr_output(ksec,3,jclass+1), jclass=1,MAX(1,10) )
2434	CALL FLUSH(numdct_NOOS_h)
2435	ENDIF
2436
2437
2438	ENDIF
2439
2440
2441	CALL wrk_dealloc(nb_type_class , zsumclasses )
2442
2443	DEALLOCATE(noos_iom_dummy)
2444
2445
2446
2447	END SUBROUTINE dia_dct_wri_NOOS_h
2448
2449	SUBROUTINE dia_dct_wri(kt,ksec,sec)
2450	!!-------------------------------------------------------------
2451	!! Write transport output in numdct
2452	!!
2453	!! Purpose: Write transports in ascii files
2454	!!
2455	!! Method:
2456	!! 1. Write volume transports in "volume_transport"
2457	!! Unit: Sv : area * Velocity / 1.e6
2458	!!
2459	!! 2. Write heat transports in "heat_transport"
2460	!! Unit: Peta W : area * Velocity * T * rhau * Cp / 1.e15
2461	!!
2462	!! 3. Write salt transports in "salt_transport"
2463	!! Unit: 10^9 g m^3 / s : area * Velocity * S / 1.e6
2464	!!
2465	!!-------------------------------------------------------------
2466	!!arguments
2467	INTEGER, INTENT(IN) :: kt ! time-step
2468	TYPE(SECTION), INTENT(INOUT) :: sec ! section to write
2469	INTEGER ,INTENT(IN) :: ksec ! section number
2470
2471	!!local declarations
2472	INTEGER :: jclass ! Dummy loop
2473	CHARACTER(len=2) :: classe ! Classname
2474	REAL(wp) :: zbnd1,zbnd2 ! Class bounds
2475	REAL(wp) :: zslope ! section's slope coeff
2476	!
2477	REAL(wp), POINTER, DIMENSION(:):: zsumclasses ! 1D workspace
2478	!!-------------------------------------------------------------
2479	CALL wrk_alloc(nb_type_class , zsumclasses )
2480
2481	zsumclasses(:)=0._wp
2482	zslope = sec%slopeSection
2483
2484
2485	DO jclass=1,MAX(1,sec%nb_class-1)
2486
2487	classe = 'N '
2488	zbnd1 = 0._wp
2489	zbnd2 = 0._wp
2490	zsumclasses(1:nb_type_class)=zsumclasses(1:nb_type_class)+sec%transport(1:nb_type_class,jclass)
2491
2492
2493	!insitu density classes transports
2494	IF( ( sec%zsigi(jclass) .NE. 99._wp ) .AND. &
2495	( sec%zsigi(jclass+1) .NE. 99._wp ) )THEN
2496	classe = 'DI '
2497	zbnd1 = sec%zsigi(jclass)
2498	zbnd2 = sec%zsigi(jclass+1)
2499	ENDIF
2500	!potential density classes transports
2501	IF( ( sec%zsigp(jclass) .NE. 99._wp ) .AND. &
2502	( sec%zsigp(jclass+1) .NE. 99._wp ) )THEN
2503	classe = 'DP '
2504	zbnd1 = sec%zsigp(jclass)
2505	zbnd2 = sec%zsigp(jclass+1)
2506	ENDIF
2507	!depth classes transports
2508	IF( ( sec%zlay(jclass) .NE. 99._wp ) .AND. &
2509	( sec%zlay(jclass+1) .NE. 99._wp ) )THEN
2510	classe = 'Z '
2511	zbnd1 = sec%zlay(jclass)
2512	zbnd2 = sec%zlay(jclass+1)
2513	ENDIF
2514	!salinity classes transports
2515	IF( ( sec%zsal(jclass) .NE. 99._wp ) .AND. &
2516	( sec%zsal(jclass+1) .NE. 99._wp ) )THEN
2517	classe = 'S '
2518	zbnd1 = sec%zsal(jclass)
2519	zbnd2 = sec%zsal(jclass+1)
2520	ENDIF
2521	!temperature classes transports
2522	IF( ( sec%ztem(jclass) .NE. 99._wp ) .AND. &
2523	( sec%ztem(jclass+1) .NE. 99._wp ) ) THEN
2524	classe = 'T '
2525	zbnd1 = sec%ztem(jclass)
2526	zbnd2 = sec%ztem(jclass+1)
2527	ENDIF
2528
2529	!write volume transport per class
2530	WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope, &
2531	jclass,classe,zbnd1,zbnd2,&
2532	sec%transport(1,jclass),sec%transport(2,jclass), &
2533	sec%transport(1,jclass)+sec%transport(2,jclass)
2534
2535	IF( sec%llstrpond )THEN
2536
2537	!write heat transport per class:
2538	WRITE(numdct_heat,119) ndastp,kt,ksec,sec%name,zslope, &
2539	jclass,classe,zbnd1,zbnd2,&
2540	sec%transport(7,jclass)1000._wprcp/1.e15,sec%transport(8,jclass)1000._wprcp/1.e15, &
2541	( sec%transport(7,jclass)+sec%transport(8,jclass) )1000._wprcp/1.e15
2542	!write salt transport per class
2543	WRITE(numdct_salt,119) ndastp,kt,ksec,sec%name,zslope, &
2544	jclass,classe,zbnd1,zbnd2,&
2545	sec%transport(9,jclass)1000._wp/1.e9,sec%transport(10,jclass)1000._wp/1.e9,&
2546	(sec%transport(9,jclass)+sec%transport(10,jclass))*1000._wp/1.e9
2547	ENDIF
2548
2549	ENDDO
2550
2551	zbnd1 = 0._wp
2552	zbnd2 = 0._wp
2553	jclass=0
2554
2555	!write total volume transport
2556	WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope, &
2557	jclass,"total",zbnd1,zbnd2,&
2558	zsumclasses(1),zsumclasses(2),zsumclasses(1)+zsumclasses(2)
2559
2560	IF( sec%llstrpond )THEN
2561
2562	!write total heat transport
2563	WRITE(numdct_heat,119) ndastp,kt,ksec,sec%name,zslope, &
2564	jclass,"total",zbnd1,zbnd2,&
2565	zsumclasses(7)* 1000._wprcp/1.e15,zsumclasses(8) 1000._wp*rcp/1.e15,&
2566	(zsumclasses(7)+zsumclasses(8) )* 1000._wp*rcp/1.e15
2567	!write total salt transport
2568	WRITE(numdct_salt,119) ndastp,kt,ksec,sec%name,zslope, &
2569	jclass,"total",zbnd1,zbnd2,&
2570	zsumclasses(9)1000._wp/1.e9,zsumclasses(10)1000._wp/1.e9,&
2571	(zsumclasses(9)+zsumclasses(10))*1000._wp/1.e9
2572	ENDIF
2573
2574
2575	IF ( sec%ll_ice_section) THEN
2576	!write total ice volume transport
2577	WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope,&
2578	jclass,"ice_vol",zbnd1,zbnd2,&
2579	sec%transport(11,1),sec%transport(12,1),&
2580	sec%transport(11,1)+sec%transport(12,1)
2581	!write total ice surface transport
2582	WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope,&
2583	jclass,"ice_surf",zbnd1,zbnd2,&
2584	sec%transport(13,1),sec%transport(14,1), &
2585	sec%transport(13,1)+sec%transport(14,1)
2586	ENDIF
2587
2588	118 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3F12.4)
2589	119 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3E15.6)
2590
2591	CALL wrk_dealloc(nb_type_class , zsumclasses )
2592	END SUBROUTINE dia_dct_wri
2593
2594	PURE FUNCTION interp(ki, kj, kk, cd_point,var)
2595	!!----------------------------------------------------------------------
2596	!!
2597	!! Purpose: compute temperature/salinity/density at U-point or V-point
2598	!! --------
2599	!!
2600	!! Method:
2601	!! ------
2602	!!
2603	!! ====> full step and partial step
2604	!!
2605	!!
2606	!! \| I \| I+1 \| Z=temperature/salinity/density at U-poinT
2607	!! \| \| \|
2608	!! ---------------------------------------- 1. Veritcale interpolation: compute zbis
2609	!! \| \| \| interpolation between ptab(I,J,K) and ptab(I,J,K+1)
2610	!! \| \| \| zbis =
2611	!! \| \| \| [ e3w(I+1,J,K)ptab(I,J,K) + ( e3w(I,J,K) - e3w(I+1,J,K) ) ptab(I,J,K-1) ]
2612	!! \| \| \| /[ e3w(I+1,J,K) + e3w(I,J,K) - e3w(I+1,J,K) ]
2613	!! \| \| \|
2614	!! \| \| \| 2. Horizontal interpolation: compute value at U/V point
2615	!!K-1 \| ptab(I,J,K-1) \| \| interpolation between zbis and ptab(I+1,J,K)
2616	!! \| . \| \|
2617	!! \| . \| \| interp = ( 0.5zet2zbis + 0.5zet1ptab(I+1,J,K) )/(0.5zet2+0.5zet1)
2618	!! \| . \| \|
2619	!! ------------------------------------------
2620	!! \| . \| \|
2621	!! \| . \| \|
2622	!! \| . \| \|
2623	!!K \| zbis.......U...ptab(I+1,J,K) \|
2624	!! \| . \| \|
2625	!! \| ptab(I,J,K) \| \|
2626	!! \| \|------------------\|
2627	!! \| \| partials \|
2628	!! \| \| steps \|
2629	!! -------------------------------------------
2630	!! <----zet1------><----zet2--------->
2631	!!
2632	!!
2633	!! ====> s-coordinate
2634	!!
2635	!! \| \| \| 1. Compute distance between T1 and U points: SQRT( zdep1^2 + (0.5 * zet1 )^2
2636	!! \| \| \| Compute distance between T2 and U points: SQRT( zdep2^2 + (0.5 * zet2 )^2
2637	!! \| \| ptab(I+1,J,K) \|
2638	!! \| \| T2 \| 2. Interpolation between T1 and T2 values at U point
2639	!! \| \| ^ \|
2640	!! \| \| \| zdep2 \|
2641	!! \| \| \| \|
2642	!! \| ^ U v \|
2643	!! \| \| \| \|
2644	!! \| \| zdep1 \| \|
2645	!! \| v \| \|
2646	!! \| T1 \| \|
2647	!! \| ptab(I,J,K) \| \|
2648	!! \| \| \|
2649	!! \| \| \|
2650	!!
2651	!! <----zet1--------><----zet2--------->
2652	!!
2653	!!----------------------------------------------------------------------
2654	!*arguments
2655	INTEGER, INTENT(IN) :: ki, kj, kk ! coordinate of point
2656	INTEGER, INTENT(IN) :: var ! which variable
2657	CHARACTER(len=1), INTENT(IN) :: cd_point ! type of point (U, V)
2658	REAL(wp) :: interp ! interpolated variable
2659
2660	!*local declations
2661	INTEGER :: ii1, ij1, ii2, ij2 ! local integer
2662	REAL(wp):: ze3t, zfse3, zwgt1, zwgt2, zbis, zdepu ! local real
2663	REAL(wp):: zet1, zet2 ! weight for interpolation
2664	REAL(wp):: zdep1,zdep2 ! differences of depth
2665	REAL(wp):: zmsk ! mask value
2666	!!----------------------------------------------------------------------
2667
2668
2669
2670	IF( cd_point=='U' )THEN
2671	ii1 = ki ; ij1 = kj
2672	ii2 = ki+1 ; ij2 = kj
2673
2674	zet1=e1t(ii1,ij1)
2675	zet2=e1t(ii2,ij2)
2676	zmsk=umask(ii1,ij1,kk)
2677
2678
2679	ELSE ! cd_point=='V'
2680	ii1 = ki ; ij1 = kj
2681	ii2 = ki ; ij2 = kj+1
2682
2683	zet1=e2t(ii1,ij1)
2684	zet2=e2t(ii2,ij2)
2685	zmsk=vmask(ii1,ij1,kk)
2686
2687	ENDIF
2688
2689	IF( ln_sco )THEN ! s-coordinate case
2690
2691	zdepu = ( fsdept(ii1,ij1,kk) + fsdept(ii2,ij2,kk) ) /2
2692	zdep1 = fsdept(ii1,ij1,kk) - zdepu
2693	zdep2 = fsdept(ii2,ij2,kk) - zdepu
2694
2695	! weights
2696	zwgt1 = SQRT( ( 0.5 * zet1 ) * ( 0.5 * zet1 ) + ( zdep1 * zdep1 ) )
2697	zwgt2 = SQRT( ( 0.5 * zet2 ) * ( 0.5 * zet2 ) + ( zdep2 * zdep2 ) )
2698
2699	! result
2700	SELECT CASE( var )
2701	CASE(0) ; interp = zmsk * ( zwgt2 * tsn(ii1,ij1,kk,jp_tem) + zwgt1 * tsn(ii1,ij1,kk,jp_tem) ) / ( zwgt2 + zwgt1 )
2702	CASE(1) ; interp = zmsk * ( zwgt2 * tsn(ii1,ij1,kk,jp_sal) + zwgt1 * tsn(ii1,ij1,kk,jp_sal) ) / ( zwgt2 + zwgt1 )
2703	CASE(2) ; interp = zmsk * ( zwgt2 * rhop(ii1,ij1,kk) + zwgt1 * rhop(ii1,ij1,kk) ) / ( zwgt2 + zwgt1 )
2704	CASE(3) ; interp = zmsk * ( zwgt2 * (rhd(ii1,ij1,kk)rau0+rau0) + zwgt1 (rhd(ii1,ij1,kk)*rau0+rau0) ) / ( zwgt2 + zwgt1 )
2705	END SELECT
2706
2707	ELSE ! full step or partial step case
2708
2709	#if defined key_vvl
2710
2711	ze3t = fse3t_n(ii2,ij2,kk) - fse3t_n(ii1,ij1,kk)
2712	zwgt1 = ( fse3w_n(ii2,ij2,kk) - fse3w_n(ii1,ij1,kk) ) / fse3w_n(ii2,ij2,kk)
2713	zwgt2 = ( fse3w_n(ii1,ij1,kk) - fse3w_n(ii2,ij2,kk) ) / fse3w_n(ii1,ij1,kk)
2714
2715	#else
2716
2717	ze3t = fse3t(ii2,ij2,kk) - fse3t(ii1,ij1,kk)
2718	zwgt1 = ( fse3w(ii2,ij2,kk) - fse3w(ii1,ij1,kk) ) / fse3w(ii2,ij2,kk)
2719	zwgt2 = ( fse3w(ii1,ij1,kk) - fse3w(ii2,ij2,kk) ) / fse3w(ii1,ij1,kk)
2720
2721	#endif
2722
2723	IF(kk .NE. 1)THEN
2724
2725	IF( ze3t >= 0. )THEN
2726	! zbis
2727	SELECT CASE( var )
2728	CASE(0)
2729	zbis = tsn(ii2,ij2,kk,jp_tem) + zwgt1 * (tsn(ii2,ij2,kk-1,jp_tem)-tsn(ii2,ij2,kk,jp_tem) )
2730	interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_tem) + zet1 * zbis )/( zet1 + zet2 )
2731	CASE(1)
2732	zbis = tsn(ii2,ij2,kk,jp_sal) + zwgt1 * (tsn(ii2,ij2,kk-1,jp_sal)-tsn(ii2,ij2,kk,jp_sal) )
2733	interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_sal) + zet1 * zbis )/( zet1 + zet2 )
2734	CASE(2)
2735	zbis = rhop(ii2,ij2,kk) + zwgt1 * (rhop(ii2,ij2,kk-1)-rhop(ii2,ij2,kk) )
2736	interp = zmsk * ( zet2 * rhop(ii1,ij1,kk) + zet1 * zbis )/( zet1 + zet2 )
2737	CASE(3)
2738	zbis = (rhd(ii2,ij2,kk)rau0+rau0) + zwgt1 ( (rhd(ii2,ij2,kk-1)rau0+rau0)-(rhd(ii2,ij2,kk)rau0+rau0) )
2739	interp = zmsk * ( zet2 * (rhd(ii1,ij1,kk)rau0+rau0) + zet1 zbis )/( zet1 + zet2 )
2740	END SELECT
2741	! result
2742	ELSE
2743	! zbis
2744	SELECT CASE( var )
2745	CASE(0)
2746	zbis = tsn(ii1,ij1,kk,jp_tem) + zwgt2 * ( tsn(ii1,ij1,kk-1,jp_tem) - tsn(ii1,ij2,kk,jp_tem) )
2747	interp = zmsk * ( zet2 * zbis + zet1 * tsn(ii2,ij2,kk,jp_tem) )/( zet1 + zet2 )
2748	CASE(1)
2749	zbis = tsn(ii1,ij1,kk,jp_sal) + zwgt2 * ( tsn(ii1,ij1,kk-1,jp_sal) - tsn(ii1,ij2,kk,jp_sal) )
2750	interp = zmsk * ( zet2 * zbis + zet1 * tsn(ii2,ij2,kk,jp_sal) )/( zet1 + zet2 )
2751	CASE(2)
2752	zbis = rhop(ii1,ij1,kk) + zwgt2 * ( rhop(ii1,ij1,kk-1) - rhop(ii1,ij2,kk) )
2753	interp = zmsk * ( zet2 * zbis + zet1 * rhop(ii2,ij2,kk) )/( zet1 + zet2 )
2754	CASE(3)
2755	zbis = (rhd(ii1,ij1,kk)rau0+rau0) + zwgt2 ( (rhd(ii1,ij1,kk-1)rau0+rau0) - (rhd(ii1,ij2,kk)rau0+rau0) )
2756	interp = zmsk * ( zet2 * zbis + zet1 * (rhd(ii2,ij2,kk)*rau0+rau0) )/( zet1 + zet2 )
2757	END SELECT
2758	ENDIF
2759
2760	ELSE
2761	SELECT CASE( var )
2762	CASE(0)
2763	interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_tem) + zet1 * tsn(ii2,ij2,kk,jp_tem) )/( zet1 + zet2 )
2764	CASE(1)
2765	interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_sal) + zet1 * tsn(ii2,ij2,kk,jp_sal) )/( zet1 + zet2 )
2766	CASE(2)
2767	interp = zmsk * ( zet2 * rhop(ii1,ij1,kk) + zet1 * rhop(ii2,ij2,kk) )/( zet1 + zet2 )
2768	CASE(3)
2769	interp = zmsk * ( zet2 * (rhd(ii1,ij1,kk)rau0+rau0) + zet1 (rhd(ii2,ij2,kk)*rau0+rau0) )/( zet1 + zet2 )
2770	END SELECT
2771	ENDIF
2772
2773	ENDIF
2774
2775	END FUNCTION interp
2776
2777	#else
2778	!!----------------------------------------------------------------------
2779	!! Default option : Dummy module
2780	!!----------------------------------------------------------------------
2781	LOGICAL, PUBLIC, PARAMETER :: lk_diadct = .FALSE. !: diamht flag
2782	PUBLIC
2783	!! $Id$
2784	CONTAINS
2785
2786	SUBROUTINE dia_dct_init ! Dummy routine
2787	WRITE(,) 'dia_dct_init: You should not have seen this print! error?', kt
2788	END SUBROUTINE dia_dct_init
2789
2790	SUBROUTINE dia_dct( kt ) ! Dummy routine
2791	INTEGER, INTENT( in ) :: kt ! ocean time-step index
2792	WRITE(,) 'dia_dct: You should not have seen this print! error?', kt
2793	END SUBROUTINE dia_dct
2794	#endif
2795
2796	END MODULE diadct

Note: See TracBrowser for help on using the repository browser.

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

Context Navigation

source: branches/UKMO/AMM15_v3_6_STABLE_package_reanalysis2/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90 @ 11038

Download in other formats: