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

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

diadct.F90 in branches/UKMO/CO6_shelfclimate_fabm_noos/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

Context Navigation

source: branches/UKMO/CO6_shelfclimate_fabm_noos/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90 @ 10492

Last change on this file since 10492 was 10492, checked in by hadjt, 5 years ago
Updated CO6_shelfclimate_fabm for new NOOS Cross sections
Property svn:executable set to ``*
File size: 133.3 KB

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

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

Download in other formats: