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 @ 11264

Last change on this file since 11264 was 11264, checked in by hadjt, 5 years ago
Updated the maximum size of the NOOS cross-section size (from 150-275) in DIA/diadct.F90
Property svn:executable set to ``*
File size: 133.4 KB

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

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

Download in other formats: