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.

bdydta.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

Context Navigation

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90 @ 4474

Last change on this file since 4474 was 4474, checked in by trackstand2, 10 years ago
Bug fix in bdydta - use jpkorig rather than jpk when reading boundary fields from file
Property svn:keywords set to `Id`
File size: 52.5 KB

Rev	Line
[911]	1	MODULE bdydta
[1125]	2	!!======================================================================
	3	!! * MODULE bdydta *
[911]	4	!! Open boundary data : read the data for the unstructured open boundaries.
[1125]	5	!!======================================================================
	6	!! History : 1.0 ! 2005-01 (J. Chanut, A. Sellar) Original code
	7	!! - ! 2007-01 (D. Storkey) Update to use IOM module
[2528]	8	!! - ! 2007-07 (D. Storkey) add bdy_dta_fla
[1125]	9	!! 3.0 ! 2008-04 (NEMO team) add in the reference version
[2528]	10	!! 3.3 ! 2010-09 (E.O'Dea) modifications for Shelf configurations
	11	!! 3.3 ! 2010-09 (D.Storkey) add ice boundary conditions
[1125]	12	!!----------------------------------------------------------------------
	13	#if defined key_bdy
	14	!!----------------------------------------------------------------------
	15	!! 'key_bdy' Unstructured Open Boundary Conditions
	16	!!----------------------------------------------------------------------
[2528]	17	!! bdy_dta_frs : read u, v, t, s data along open boundaries
	18	!! bdy_dta_fla : read depth-mean velocities and elevation along open boundaries
[1125]	19	!!----------------------------------------------------------------------
[911]	20	USE oce ! ocean dynamics and tracers
	21	USE dom_oce ! ocean space and time domain
	22	USE phycst ! physical constants
	23	USE bdy_oce ! ocean open boundary conditions
	24	USE bdytides ! tidal forcing at boundaries
	25	USE iom
	26	USE ioipsl
	27	USE in_out_manager ! I/O logical units
[2528]	28	#if defined key_lim2
	29	USE ice_2
	30	#endif
[911]	31
	32	IMPLICIT NONE
	33	PRIVATE
	34
[2528]	35	PUBLIC bdy_dta_frs ! routines called by step.F90
	36	PUBLIC bdy_dta_fla
[2715]	37	PUBLIC bdy_dta_alloc ! routine called by bdy_init.F90
[911]	38
[3432]	39	INTEGER :: numbdyt, numbdyu, numbdyv ! logical units for T-, U-, & V-points data file, resp.
	40	INTEGER :: ntimes_bdy ! exact number of time dumps in data files
	41	INTEGER :: nbdy_b, nbdy_a ! record of bdy data file for before and after time step
	42	INTEGER :: numbdyt_bt, numbdyu_bt, numbdyv_bt ! logical unit for T-, U- & V-points data file, resp.
	43	INTEGER :: ntimes_bdy_bt ! exact number of time dumps in data files
	44	INTEGER :: nbdy_b_bt, nbdy_a_bt ! record of bdy data file for before and after time step
[911]	45
[3432]	46	INTEGER, DIMENSION (jpbtime) :: istep, istep_bt ! time array in seconds in each data file
[911]	47
[3432]	48	REAL(wp) :: zoffset ! time offset between time origin in file & start time of model run
[911]	49
[2715]	50	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tbdydta, sbdydta ! time interpolated values of T and S bdy data
	51	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ubdydta, vbdydta ! time interpolated values of U and V bdy data
[2528]	52	REAL(wp), DIMENSION(jpbdim,2) :: ubtbdydta, vbtbdydta ! Arrays used for time interpolation of bdy data
	53	REAL(wp), DIMENSION(jpbdim,2) :: sshbdydta ! bdy data of ssh
[911]	54
[2528]	55	#if defined key_lim2
	56	REAL(wp), DIMENSION(jpbdim,2) :: frld_bdydta ! }
	57	REAL(wp), DIMENSION(jpbdim,2) :: hicif_bdydta ! } Arrays used for time interp. of ice bdy data
	58	REAL(wp), DIMENSION(jpbdim,2) :: hsnif_bdydta ! }
	59	#endif
	60
[3211]	61	!! * Control permutation of array indices
	62	# include "oce_ftrans.h90"
	63	# include "dom_oce_ftrans.h90"
	64
[1125]	65	!!----------------------------------------------------------------------
[2528]	66	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[1146]	67	!! $Id$
[2528]	68	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1125]	69	!!----------------------------------------------------------------------
[911]	70	CONTAINS
	71
[2715]	72	FUNCTION bdy_dta_alloc()
	73	!!----------------------------------------------------------------------
	74	USE lib_mpp, ONLY: ctl_warn, mpp_sum
	75	!
	76	INTEGER :: bdy_dta_alloc
	77	!!----------------------------------------------------------------------
	78	!
	79	ALLOCATE(tbdydta(jpbdim,jpk,2), sbdydta(jpbdim,jpk,2), &
	80	ubdydta(jpbdim,jpk,2), vbdydta(jpbdim,jpk,2), Stat=bdy_dta_alloc)
	81
	82	IF( lk_mpp ) CALL mpp_sum ( bdy_dta_alloc )
	83	IF(bdy_dta_alloc /= 0) CALL ctl_warn('bdy_dta_alloc: failed to allocate arrays')
	84
	85	END FUNCTION bdy_dta_alloc
	86
	87
[2528]	88	SUBROUTINE bdy_dta_frs( kt )
[1125]	89	!!----------------------------------------------------------------------
[2528]	90	!! * SUBROUTINE bdy_dta_frs *
[911]	91	!!
[1125]	92	!! ** Purpose : Read unstructured boundary data for FRS condition.
[911]	93	!!
[1125]	94	!! ** Method : At the first timestep, read in boundary data for two
	95	!! times from the file and time-interpolate. At other
	96	!! timesteps, check to see if we need another time from
	97	!! the file. If so read it in. Time interpolate.
	98	!!----------------------------------------------------------------------
[2528]	99	INTEGER, INTENT( in ) :: kt ! ocean time-step index (for timesplitting option, otherwise zero)
[911]	100	!!
[1125]	101	CHARACTER(LEN=80), DIMENSION(3) :: clfile ! names of input files
	102	CHARACTER(LEN=70 ) :: clunits ! units attribute of time coordinate
	103	LOGICAL :: lect ! flag for reading
	104	INTEGER :: it, ib, ik, igrd ! dummy loop indices
	105	INTEGER :: igrd_start, igrd_end ! start and end of loops on igrd
	106	INTEGER :: idvar ! netcdf var ID
	107	INTEGER :: iman, i15, imois ! Time variables for monthly clim forcing
	108	INTEGER :: ntimes_bdyt, ntimes_bdyu, ntimes_bdyv
[911]	109	INTEGER :: itimer, totime
[1125]	110	INTEGER :: ii, ij ! array addresses
[2528]	111	INTEGER :: ipi, ipj, ipk, inum ! local integers (NetCDF read)
[911]	112	INTEGER :: iyear0, imonth0, iday0
	113	INTEGER :: ihours0, iminutes0, isec0
[1125]	114	INTEGER :: iyear, imonth, iday, isecs
[911]	115	INTEGER, DIMENSION(jpbtime) :: istept, istepu, istepv ! time arrays from data files
[1125]	116	REAL(wp) :: dayfrac, zxy, zoffsett
	117	REAL(wp) :: zoffsetu, zoffsetv
[911]	118	REAL(wp) :: dayjul0, zdayjulini
[1125]	119	REAL(wp), DIMENSION(jpbtime) :: zstepr ! REAL time array from data files
[4474]	120	!! DCSE_NEMO: do not ftrans becase is used in calls to iom library that all use z-last ordering.
[3432]	121	!!!!! zdta :I :I :z
[4474]	122	REAL(wp), DIMENSION(jpbdta,1,jpkorig) :: zdta ! temporary array for data fields
[911]	123	!!---------------------------------------------------------------------------
	124
[1125]	125
[2528]	126	IF( ln_dyn_frs .OR. ln_tra_frs &
	127	& .OR. ln_ice_frs ) THEN ! If these are both false then this routine does nothing
	128
[911]	129	! -------------------- !
	130	! Initialization !
	131	! -------------------- !
	132
	133	lect = .false. ! If true, read a time record
	134
	135	! Some time variables for monthly climatological forcing:
	136	! *******************************************************
[2528]	137
	138	!!gm here use directely daymod calendar variables
[911]	139
[1125]	140	iman = INT( raamo ) ! Number of months in a year
[911]	141
[1713]	142	i15 = INT( 2*REAL( nday, wp ) / ( REAL( nmonth_len(nmonth), wp ) + 0.5 ) )
[911]	143	! i15=0 if the current day is in the first half of the month, else i15=1
	144
	145	imois = nmonth + i15 - 1 ! imois is the first month record
[1125]	146	IF( imois == 0 ) imois = iman
[911]	147
	148	! Time variable for non-climatological forcing:
	149	! *********************************************
	150	itimer = (kt-nit000+1)*rdt ! current time in seconds for interpolation
	151
	152
[1125]	153	! !-------------------!
	154	IF( kt == nit000 ) THEN ! First call only !
	155	! !-------------------!
	156	istep(:) = 0
[2528]	157	nbdy_b = 0
	158	nbdy_a = 0
[911]	159
[1125]	160	! Get time information from bdy data file
	161	! ***************************************
[911]	162
[1125]	163	IF(lwp) WRITE(numout,*)
[2528]	164	IF(lwp) WRITE(numout,*) 'bdy_dta_frs : Initialize unstructured boundary data'
[1125]	165	IF(lwp) WRITE(numout,*) '~~~~~~~'
[911]	166
[2528]	167	IF ( nn_dtactl == 0 ) THEN
[1125]	168	!
	169	IF(lwp) WRITE(numout,*) ' Bdy data are taken from initial conditions'
	170	!
[2528]	171	ELSEIF (nn_dtactl == 1) THEN
[1125]	172	!
	173	IF(lwp) WRITE(numout,*) ' Bdy data are read in netcdf files'
	174	!
[1713]	175	dayfrac = adatrj - REAL( itimer, wp ) / 86400. ! day fraction at time step kt-1
	176	dayfrac = dayfrac - INT ( dayfrac ) !
	177	totime = ( nitend - nit000 + 1 ) * rdt ! Total time of the run to verify that all the
	178	! ! necessary time dumps in file are included
[1125]	179	!
[2528]	180	clfile(1) = cn_dta_frs_T
	181	clfile(2) = cn_dta_frs_U
	182	clfile(3) = cn_dta_frs_V
[1125]	183	!
	184	! how many files are we to read in?
	185	igrd_start = 1
	186	igrd_end = 3
[2528]	187	IF(.NOT. ln_tra_frs .AND. .NOT. ln_ice_frs) THEN ! No T-grid file.
[1125]	188	igrd_start = 2
[2528]	189	ELSEIF ( .NOT. ln_dyn_frs ) THEN ! No U-grid or V-grid file.
[1125]	190	igrd_end = 1
	191	ENDIF
[911]	192
[1125]	193	DO igrd = igrd_start, igrd_end ! loop over T, U & V grid !
	194	! !---------------------------!
	195	CALL iom_open( clfile(igrd), inum )
	196	CALL iom_gettime( inum, zstepr, kntime=ntimes_bdy, cdunits=clunits )
[911]	197
[1125]	198	SELECT CASE( igrd )
[2528]	199	CASE (1) ; numbdyt = inum
	200	CASE (2) ; numbdyu = inum
	201	CASE (3) ; numbdyv = inum
[1125]	202	END SELECT
[911]	203
[1125]	204	! Calculate time offset
	205	READ(clunits,7000) iyear0, imonth0, iday0, ihours0, iminutes0, isec0
	206	! Convert time origin in file to julian days
	207	isec0 = isec0 + ihours060.60. + iminutes0*60.
[1713]	208	CALL ymds2ju(iyear0, imonth0, iday0, REAL(isec0, wp), dayjul0)
[1125]	209	! Compute model initialization time
	210	iyear = ndastp / 10000
	211	imonth = ( ndastp - iyear * 10000 ) / 100
	212	iday = ndastp - iyear * 10000 - imonth * 100
	213	isecs = dayfrac * 86400
[1713]	214	CALL ymds2ju(iyear, imonth, iday, REAL(isecs, wp) , zdayjulini)
[1125]	215	! offset from initialization date:
	216	zoffset = (dayjul0-zdayjulini)*86400
	217	!
[3432]	218	! ARP - Cannot have 'string literals' in an input format with the Cray compiler
	219	! so replace with corresponding number of spaces instead.
	220	!7000 FORMAT('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2)
	221	7000 FORMAT(14x, I4.4,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2)
[911]	222
[1125]	223	!! TO BE DONE... Check consistency between calendar from file
	224	!! (available optionally from iom_gettime) and calendar in model
	225	!! when calendar in model available outside of IOIPSL.
[911]	226
[1125]	227	IF(lwp) WRITE(numout,*) 'number of times: ',ntimes_bdy
	228	IF(lwp) WRITE(numout,*) 'offset: ',zoffset
	229	IF(lwp) WRITE(numout,*) 'totime: ',totime
[2528]	230	IF(lwp) WRITE(numout,*) 'zstepr: ',zstepr(1:ntimes_bdy)
[911]	231
[1125]	232	! Check that there are not too many times in the file.
	233	IF( ntimes_bdy > jpbtime ) THEN
	234	WRITE(ctmp1,*) 'Check file: ', clfile(igrd), 'jpbtime= ', jpbtime, ' ntimes_bdy= ', ntimes_bdy
	235	CALL ctl_stop( 'Number of time dumps in files exceed jpbtime parameter', ctmp1 )
	236	ENDIF
[911]	237
[1125]	238	! Check that time array increases:
	239	it = 1
[2528]	240	DO WHILE( zstepr(it+1) > zstepr(it) .AND. it /= ntimes_bdy - 1 )
	241	it = it + 1
[1125]	242	END DO
[2528]	243	!
	244	IF( it /= ntimes_bdy-1 .AND. ntimes_bdy > 1 ) THEN
[1125]	245	WRITE(ctmp1,*) 'Check file: ', clfile(igrd)
	246	CALL ctl_stop( 'Time array in unstructured boundary data files', &
	247	& 'does not continuously increase.' , ctmp1 )
	248	ENDIF
	249	!
	250	! Check that times in file span model run time:
	251	IF( zstepr(1) + zoffset > 0 ) THEN
	252	WRITE(ctmp1,*) 'Check file: ', clfile(igrd)
	253	CALL ctl_stop( 'First time dump in bdy file is after model initial time', ctmp1 )
	254	END IF
	255	IF( zstepr(ntimes_bdy) + zoffset < totime ) THEN
	256	WRITE(ctmp1,*) 'Check file: ', clfile(igrd)
	257	CALL ctl_stop( 'Last time dump in bdy file is before model final time', ctmp1 )
	258	END IF
	259	!
[2528]	260	SELECT CASE( igrd )
	261	CASE (1)
	262	ntimes_bdyt = ntimes_bdy
	263	zoffsett = zoffset
	264	istept(:) = INT( zstepr(:) + zoffset )
	265	numbdyt = inum
	266	CASE (2)
	267	ntimes_bdyu = ntimes_bdy
	268	zoffsetu = zoffset
	269	istepu(:) = INT( zstepr(:) + zoffset )
	270	numbdyu = inum
	271	CASE (3)
	272	ntimes_bdyv = ntimes_bdy
	273	zoffsetv = zoffset
	274	istepv(:) = INT( zstepr(:) + zoffset )
	275	numbdyv = inum
	276	END SELECT
[1125]	277	!
	278	END DO ! end loop over T, U & V grid
[911]	279
[1125]	280	IF (igrd_start == 1 .and. igrd_end == 3) THEN
	281	! Only test differences if we are reading in 3 files
	282	! Verify time consistency between files
	283	IF( ntimes_bdyu /= ntimes_bdyt .OR. ntimes_bdyv /= ntimes_bdyt ) THEN
	284	CALL ctl_stop( 'Bdy data files must have the same number of time dumps', &
	285	& 'Multiple time frequencies not implemented yet' )
	286	ENDIF
	287	ntimes_bdy = ntimes_bdyt
	288	!
	289	IF( zoffsetu /= zoffsett .OR. zoffsetv /= zoffsett ) THEN
	290	CALL ctl_stop( 'Bdy data files must have the same time origin', &
	291	& 'Multiple time frequencies not implemented yet' )
	292	ENDIF
	293	zoffset = zoffsett
	294	ENDIF
[911]	295
[2528]	296	IF( igrd_start == 1 ) THEN ; istep(:) = istept(:)
	297	ELSE ; istep(:) = istepu(:)
[1125]	298	ENDIF
[911]	299
[1125]	300	! Check number of time dumps:
[2528]	301	IF( ntimes_bdy == 1 .AND. .NOT. ln_clim ) THEN
[1125]	302	CALL ctl_stop( 'There is only one time dump in data files', &
[2528]	303	& 'Choose ln_clim=.true. in namelist for constant bdy forcing.' )
[1125]	304	ENDIF
[911]	305
[2528]	306	IF( ln_clim ) THEN
[1125]	307	IF( ntimes_bdy /= 1 .AND. ntimes_bdy /= 12 ) THEN
[2528]	308	CALL ctl_stop( 'For climatological boundary forcing (ln_clim=.true.),', &
[1125]	309	& 'bdy data files must contain 1 or 12 time dumps.' )
	310	ELSEIF( ntimes_bdy == 1 ) THEN
	311	IF(lwp) WRITE(numout,*)
	312	IF(lwp) WRITE(numout,*) 'We assume constant boundary forcing from bdy data files'
	313	ELSEIF( ntimes_bdy == 12 ) THEN
	314	IF(lwp) WRITE(numout,*)
	315	IF(lwp) WRITE(numout,*) 'We assume monthly (and cyclic) boundary forcing from bdy data files'
	316	ENDIF
[911]	317	ENDIF
	318
[1125]	319	! Find index of first record to read (before first model time).
	320	it = 1
	321	DO WHILE( istep(it+1) <= 0 .AND. it <= ntimes_bdy - 1 )
[2528]	322	it = it + 1
[911]	323	END DO
[1125]	324	nbdy_b = it
	325	!
[2528]	326	IF(lwp) WRITE(numout,*) 'Time offset is ',zoffset
	327	IF(lwp) WRITE(numout,*) 'First record to read is ',nbdy_b
[911]	328
[2528]	329	ENDIF ! endif (nn_dtactl == 1)
[911]	330
	331
[2528]	332	! 1.2 Read first record in file if necessary (ie if nn_dtactl == 1)
[1125]	333	! *****************************************************************
[911]	334
[2528]	335	IF( nn_dtactl == 0 ) THEN ! boundary data arrays are filled with initial conditions
[1125]	336	!
[2528]	337	IF (ln_tra_frs) THEN
	338	igrd = 1 ! T-points data
	339	DO ib = 1, nblen(igrd)
	340	ii = nbi(ib,igrd)
	341	ij = nbj(ib,igrd)
[4463]	342	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[2528]	343	tbdy(ib,ik) = tn(ii,ij,ik)
	344	sbdy(ib,ik) = sn(ii,ij,ik)
	345	END DO
	346	END DO
[911]	347	ENDIF
	348
[2528]	349	IF(ln_dyn_frs) THEN
	350	igrd = 2 ! U-points data
	351	DO ib = 1, nblen(igrd)
	352	ii = nbi(ib,igrd)
	353	ij = nbj(ib,igrd)
[4463]	354	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[2528]	355	ubdy(ib,ik) = un(ii, ij, ik)
	356	END DO
	357	END DO
	358	!
	359	igrd = 3 ! V-points data
	360	DO ib = 1, nblen(igrd)
	361	ii = nbi(ib,igrd)
	362	ij = nbj(ib,igrd)
[4463]	363	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[2528]	364	vbdy(ib,ik) = vn(ii, ij, ik)
	365	END DO
	366	END DO
[911]	367	ENDIF
[1125]	368	!
[2528]	369	#if defined key_lim2
	370	IF( ln_ice_frs ) THEN
	371	igrd = 1 ! T-points data
	372	DO ib = 1, nblen(igrd)
	373	frld_bdy (ib) = frld(nbi(ib,igrd), nbj(ib,igrd))
	374	hicif_bdy(ib) = hicif(nbi(ib,igrd), nbj(ib,igrd))
	375	hsnif_bdy(ib) = hsnif(nbi(ib,igrd), nbj(ib,igrd))
	376	END DO
	377	ENDIF
	378	#endif
	379	ELSEIF( nn_dtactl == 1 ) THEN ! Set first record in the climatological case:
[1125]	380	!
[2528]	381	IF( ln_clim .AND. ntimes_bdy == 1 ) THEN
[1125]	382	nbdy_a = 1
[2528]	383	ELSEIF( ln_clim .AND. ntimes_bdy == iman ) THEN
[1125]	384	nbdy_b = 0
	385	nbdy_a = imois
	386	ELSE
	387	nbdy_a = nbdy_b
	388	ENDIF
	389
	390	! Read first record:
	391	ipj = 1
[4474]	392	ipk = jpkorig
[1125]	393	igrd = 1
	394	ipi = nblendta(igrd)
[911]	395
[2528]	396	IF(ln_tra_frs) THEN
	397	!
[1125]	398	igrd = 1 ! Temperature
	399	IF( nblendta(igrd) <= 0 ) THEN
	400	idvar = iom_varid( numbdyt, 'votemper' )
	401	nblendta(igrd) = iom_file(numbdyt)%dimsz(1,idvar)
	402	ENDIF
[2528]	403	IF(lwp) WRITE(numout,*) 'Dim size for votemper is ', nblendta(igrd)
[1125]	404	ipi = nblendta(igrd)
[3432]	405	! DCSE_NEMO: Call iom_get and explicitly flag that array has
	406	! not been ftrans'd
	407	CALL iom_get ( numbdyt, jpdom_unknown, 'votemper', &
	408	zdta(1:ipi,1:ipj,1:ipk), nbdy_a, lzfirst=.FALSE. )
[2528]	409	!
[1125]	410	DO ib = 1, nblen(igrd)
[4463]	411
	412	ii = nbi(ib,igrd)
	413	ij = nbj(ib,igrd)
	414
	415	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	416	tbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	417	END DO
	418	END DO
	419	!
	420	igrd = 1 ! salinity
	421	IF( nblendta(igrd) .le. 0 ) THEN
	422	idvar = iom_varid( numbdyt, 'vosaline' )
	423	nblendta(igrd) = iom_file(numbdyt)%dimsz(1,idvar)
	424	ENDIF
[2528]	425	IF(lwp) WRITE(numout,*) 'Dim size for vosaline is ', nblendta(igrd)
[1125]	426	ipi = nblendta(igrd)
[3432]	427	! DCSE_NEMO: pass optional flag to signal that zdta is NOT
	428	! ftrans'd
	429	CALL iom_get ( numbdyt, jpdom_unknown, 'vosaline', &
	430	zdta(1:ipi,1:ipj,1:ipk), nbdy_a, lzfirst=.FALSE. )
[2528]	431	!
[1125]	432	DO ib = 1, nblen(igrd)
[4463]	433
	434	ii = nbi(ib,igrd)
	435	ij = nbj(ib,igrd)
	436
	437	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	438	sbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	439	END DO
	440	END DO
[2528]	441	ENDIF ! ln_tra_frs
[1125]	442
[2528]	443	IF( ln_dyn_frs ) THEN
	444	!
[1125]	445	igrd = 2 ! u-velocity
	446	IF ( nblendta(igrd) .le. 0 ) THEN
	447	idvar = iom_varid( numbdyu,'vozocrtx' )
	448	nblendta(igrd) = iom_file(numbdyu)%dimsz(1,idvar)
	449	ENDIF
[2528]	450	IF(lwp) WRITE(numout,*) 'Dim size for vozocrtx is ', nblendta(igrd)
[1125]	451	ipi = nblendta(igrd)
	452	CALL iom_get ( numbdyu, jpdom_unknown,'vozocrtx',zdta(1:ipi,1:ipj,1:ipk),nbdy_a )
	453	DO ib = 1, nblen(igrd)
	454	DO ik = 1, jpkm1
	455	ubdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	456	END DO
	457	END DO
	458	!
	459	igrd = 3 ! v-velocity
	460	IF ( nblendta(igrd) .le. 0 ) THEN
	461	idvar = iom_varid( numbdyv,'vomecrty' )
	462	nblendta(igrd) = iom_file(numbdyv)%dimsz(1,idvar)
	463	ENDIF
[2528]	464	IF(lwp) WRITE(numout,*) 'Dim size for vomecrty is ', nblendta(igrd)
[1125]	465	ipi = nblendta(igrd)
	466	CALL iom_get ( numbdyv, jpdom_unknown,'vomecrty',zdta(1:ipi,1:ipj,1:ipk),nbdy_a )
	467	DO ib = 1, nblen(igrd)
	468	DO ik = 1, jpkm1
	469	vbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	470	END DO
	471	END DO
[2528]	472	ENDIF ! ln_dyn_frs
[911]	473
[2528]	474	#if defined key_lim2
	475	IF( ln_ice_frs ) THEN
	476	!
	477	igrd=1 ! leads fraction
	478	IF(lwp) WRITE(numout,*) 'Dim size for ildsconc is ',nblendta(igrd)
	479	ipi=nblendta(igrd)
	480	CALL iom_get ( numbdyt, jpdom_unknown,'ildsconc',zdta(1:ipi,:,1),nbdy_a )
	481	DO ib=1, nblen(igrd)
	482	frld_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1)
	483	END DO
	484	!
	485	igrd=1 ! ice thickness
	486	IF(lwp) WRITE(numout,*) 'Dim size for iicethic is ',nblendta(igrd)
	487	ipi=nblendta(igrd)
	488	CALL iom_get ( numbdyt, jpdom_unknown,'iicethic',zdta(1:ipi,:,1),nbdy_a )
	489	DO ib=1, nblen(igrd)
	490	hicif_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1)
	491	END DO
	492	!
	493	igrd=1 ! snow thickness
	494	IF(lwp) WRITE(numout,*) 'Dim size for isnowthi is ',nblendta(igrd)
	495	ipi=nblendta(igrd)
	496	CALL iom_get ( numbdyt, jpdom_unknown,'isnowthi',zdta(1:ipi,:,1),nbdy_a )
	497	DO ib=1, nblen(igrd)
	498	hsnif_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1)
	499	END DO
	500	ENDIF ! just if ln_ice_frs is set
	501	#endif
[911]	502
[2528]	503	IF( .NOT.ln_clim .AND. istep(1) > 0 ) THEN ! First data time is after start of run
	504	nbdy_b = nbdy_a ! Put first value in both time levels
	505	IF( ln_tra_frs ) THEN
[1125]	506	tbdydta(:,:,1) = tbdydta(:,:,2)
	507	sbdydta(:,:,1) = sbdydta(:,:,2)
	508	ENDIF
[2528]	509	IF( ln_dyn_frs ) THEN
[1125]	510	ubdydta(:,:,1) = ubdydta(:,:,2)
	511	vbdydta(:,:,1) = vbdydta(:,:,2)
	512	ENDIF
[2528]	513	#if defined key_lim2
	514	IF( ln_ice_frs ) THEN
	515	frld_bdydta (:,1) = frld_bdydta(:,2)
	516	hicif_bdydta(:,1) = hicif_bdydta(:,2)
	517	hsnif_bdydta(:,1) = hsnif_bdydta(:,2)
	518	ENDIF
	519	#endif
[1125]	520	END IF
[2528]	521	!
	522	END IF ! nn_dtactl == 0/1
[911]	523
[1125]	524	! In the case of constant boundary forcing fill bdy arrays once for all
[2528]	525	IF( ln_clim .AND. ntimes_bdy == 1 ) THEN
	526	IF( ln_tra_frs ) THEN
	527	tbdy (:,:) = tbdydta (:,:,2)
	528	sbdy (:,:) = sbdydta (:,:,2)
[1125]	529	ENDIF
[2528]	530	IF( ln_dyn_frs) THEN
	531	ubdy (:,:) = ubdydta (:,:,2)
	532	vbdy (:,:) = vbdydta (:,:,2)
[1125]	533	ENDIF
[2528]	534	#if defined key_lim2
	535	IF( ln_ice_frs ) THEN
	536	frld_bdy (:) = frld_bdydta (:,2)
	537	hicif_bdy(:) = hicif_bdydta(:,2)
	538	hsnif_bdy(:) = hsnif_bdydta(:,2)
	539	ENDIF
	540	#endif
[911]	541
[2528]	542	IF( ln_tra_frs .OR. ln_ice_frs) CALL iom_close( numbdyt )
	543	IF( ln_dyn_frs ) CALL iom_close( numbdyu )
	544	IF( ln_dyn_frs ) CALL iom_close( numbdyv )
[1125]	545	END IF
[2528]	546	!
[1125]	547	ENDIF ! End if nit000
[911]	548
	549
[1125]	550	! !---------------------!
[2528]	551	IF( nn_dtactl == 1 .AND. ntimes_bdy > 1 ) THEN ! at each time step !
	552	! !---------------------!
[1125]	553	! Read one more record if necessary
	554	!**********************************
[911]	555
[2528]	556	IF( ln_clim .AND. imois /= nbdy_b ) THEN ! remember that nbdy_b=0 for kt=nit000
	557	nbdy_b = imois
	558	nbdy_a = imois + 1
	559	nbdy_b = MOD( nbdy_b, iman ) ; IF( nbdy_b == 0 ) nbdy_b = iman
	560	nbdy_a = MOD( nbdy_a, iman ) ; IF( nbdy_a == 0 ) nbdy_a = iman
	561	lect=.true.
	562	ELSEIF( .NOT.ln_clim .AND. itimer >= istep(nbdy_a) ) THEN
[911]	563
[2528]	564	IF( nbdy_a < ntimes_bdy ) THEN
	565	nbdy_b = nbdy_a
	566	nbdy_a = nbdy_a + 1
	567	lect =.true.
	568	ELSE
	569	! We have reached the end of the file
	570	! put the last data time into both time levels
	571	nbdy_b = nbdy_a
	572	IF(ln_tra_frs) THEN
	573	tbdydta(:,:,1) = tbdydta(:,:,2)
	574	sbdydta(:,:,1) = sbdydta(:,:,2)
	575	ENDIF
	576	IF(ln_dyn_frs) THEN
	577	ubdydta(:,:,1) = ubdydta(:,:,2)
	578	vbdydta(:,:,1) = vbdydta(:,:,2)
	579	ENDIF
	580	#if defined key_lim2
	581	IF(ln_ice_frs) THEN
	582	frld_bdydta (:,1) = frld_bdydta (:,2)
	583	hicif_bdydta(:,1) = hicif_bdydta(:,2)
	584	hsnif_bdydta(:,1) = hsnif_bdydta(:,2)
	585	ENDIF
	586	#endif
[1125]	587	END IF ! nbdy_a < ntimes_bdy
[2528]	588	!
[911]	589	END IF
	590
[2528]	591	IF( lect ) THEN ! Swap arrays
	592	IF( ln_tra_frs ) THEN
[1125]	593	tbdydta(:,:,1) = tbdydta(:,:,2)
	594	sbdydta(:,:,1) = sbdydta(:,:,2)
	595	ENDIF
[2528]	596	IF( ln_dyn_frs ) THEN
[1125]	597	ubdydta(:,:,1) = ubdydta(:,:,2)
	598	vbdydta(:,:,1) = vbdydta(:,:,2)
	599	ENDIF
[2528]	600	#if defined key_lim2
	601	IF( ln_ice_frs ) THEN
	602	frld_bdydta (:,1) = frld_bdydta (:,2)
	603	hicif_bdydta(:,1) = hicif_bdydta(:,2)
	604	hsnif_bdydta(:,1) = hsnif_bdydta(:,2)
	605	ENDIF
	606	#endif
[1125]	607	! read another set
	608	ipj = 1
[4474]	609	ipk = jpkorig
[911]	610
[2528]	611	IF( ln_tra_frs ) THEN
[1125]	612	!
	613	igrd = 1 ! temperature
	614	ipi = nblendta(igrd)
[3432]	615	! DCSE_NEMO: pass optional flag to signal that zdta is NOT
	616	! ftrans'd
	617	CALL iom_get ( numbdyt, jpdom_unknown, 'votemper', &
	618	zdta(1:ipi,1:ipj,1:ipk), nbdy_a, lzfirst=.FALSE. )
[1125]	619	DO ib = 1, nblen(igrd)
[4463]	620
	621	ii = nbi(ib,igrd)
	622	ij = nbj(ib,igrd)
	623
	624	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	625	tbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	626	END DO
	627	END DO
	628	!
	629	igrd = 1 ! salinity
	630	ipi = nblendta(igrd)
[3432]	631	! DCSE_NEMO: pass optional flag to signal that zdta is NOT
	632	! ftrans'd
	633	CALL iom_get ( numbdyt, jpdom_unknown, 'vosaline', &
	634	zdta(1:ipi,1:ipj,1:ipk), nbdy_a, lzfirst=.FALSE. )
[1125]	635	DO ib = 1, nblen(igrd)
[4463]	636	ii = nbi(ib,igrd)
	637	ij = nbj(ib,igrd)
	638	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	639	sbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	640	END DO
	641	END DO
[2528]	642	ENDIF ! ln_tra_frs
[911]	643
[2528]	644	IF(ln_dyn_frs) THEN
[1125]	645	!
	646	igrd = 2 ! u-velocity
	647	ipi = nblendta(igrd)
	648	CALL iom_get ( numbdyu, jpdom_unknown,'vozocrtx',zdta(1:ipi,1:ipj,1:ipk),nbdy_a )
	649	DO ib = 1, nblen(igrd)
	650	DO ik = 1, jpkm1
	651	ubdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	652	END DO
	653	END DO
	654	!
	655	igrd = 3 ! v-velocity
	656	ipi = nblendta(igrd)
	657	CALL iom_get ( numbdyv, jpdom_unknown,'vomecrty',zdta(1:ipi,1:ipj,1:ipk),nbdy_a )
	658	DO ib = 1, nblen(igrd)
	659	DO ik = 1, jpkm1
	660	vbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik)
	661	END DO
	662	END DO
[2528]	663	ENDIF ! ln_dyn_frs
[1125]	664	!
[2528]	665	#if defined key_lim2
	666	IF(ln_ice_frs) THEN
	667	!
	668	igrd = 1 ! ice concentration
	669	ipi=nblendta(igrd)
	670	CALL iom_get ( numbdyt, jpdom_unknown,'ildsconc',zdta(1:ipi,:,1),nbdy_a )
	671	DO ib=1, nblen(igrd)
	672	frld_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 )
	673	END DO
	674	!
	675	igrd=1 ! ice thickness
	676	ipi=nblendta(igrd)
	677	CALL iom_get ( numbdyt, jpdom_unknown,'iicethic',zdta(1:ipi,:,1),nbdy_a )
	678	DO ib=1, nblen(igrd)
	679	hicif_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 )
	680	END DO
	681	!
	682	igrd=1 ! snow thickness
	683	ipi=nblendta(igrd)
	684	CALL iom_get ( numbdyt, jpdom_unknown,'isnowthi',zdta(1:ipi,:,1),nbdy_a )
	685	DO ib=1, nblen(igrd)
	686	hsnif_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 )
	687	END DO
	688	ENDIF ! ln_ice_frs
	689	#endif
	690	!
	691	IF(lwp) WRITE(numout,*) 'bdy_dta_frs : first record file used nbdy_b ',nbdy_b
[1125]	692	IF(lwp) WRITE(numout,*) '~~~~~~~~ last record file used nbdy_a ',nbdy_a
[2528]	693	IF (.NOT.ln_clim) THEN
[1125]	694	IF(lwp) WRITE(numout,*) 'first record time (s): ', istep(nbdy_b)
	695	IF(lwp) WRITE(numout,*) 'model time (s) : ', itimer
	696	IF(lwp) WRITE(numout,*) 'second record time (s): ', istep(nbdy_a)
	697	ENDIF
	698	!
	699	ENDIF ! end lect=.true.
[911]	700
	701
[1125]	702	! Interpolate linearly
	703	! ********************
	704	!
[2528]	705	IF( ln_clim ) THEN ; zxy = REAL( nday ) / REAL( nmonth_len(nbdy_b) ) + 0.5 - i15
	706	ELSEIF( istep(nbdy_b) == istep(nbdy_a) ) THEN
	707	zxy = 0.0_wp
	708	ELSE ; zxy = REAL( istep(nbdy_b) - itimer ) / REAL( istep(nbdy_b) - istep(nbdy_a) )
[1125]	709	END IF
[911]	710
[2528]	711	IF(ln_tra_frs) THEN
[1125]	712	igrd = 1 ! temperature & salinity
	713	DO ib = 1, nblen(igrd)
[4463]	714	ii = nbi(ib,igrd)
	715	ij = nbj(ib,igrd)
	716	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	717	tbdy(ib,ik) = zxy * tbdydta(ib,ik,2) + (1.-zxy) * tbdydta(ib,ik,1)
	718	sbdy(ib,ik) = zxy * sbdydta(ib,ik,2) + (1.-zxy) * sbdydta(ib,ik,1)
	719	END DO
	720	END DO
	721	ENDIF
[911]	722
[2528]	723	IF(ln_dyn_frs) THEN
[1125]	724	igrd = 2 ! u-velocity
	725	DO ib = 1, nblen(igrd)
[4463]	726	ii = nbi(ib,igrd)
	727	ij = nbj(ib,igrd)
	728	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	729	ubdy(ib,ik) = zxy * ubdydta(ib,ik,2) + (1.-zxy) * ubdydta(ib,ik,1)
	730	END DO
	731	END DO
	732	!
	733	igrd = 3 ! v-velocity
	734	DO ib = 1, nblen(igrd)
[4463]	735	ii = nbi(ib,igrd)
	736	ij = nbj(ib,igrd)
	737	DO ik = 1, mbkmax(ii,ij)-1 ! jpkm1
[1125]	738	vbdy(ib,ik) = zxy * vbdydta(ib,ik,2) + (1.-zxy) * vbdydta(ib,ik,1)
	739	END DO
	740	END DO
	741	ENDIF
[911]	742
[2528]	743	#if defined key_lim2
	744	IF(ln_ice_frs) THEN
	745	igrd=1
	746	DO ib=1, nblen(igrd)
	747	frld_bdy(ib) = zxy * frld_bdydta(ib,2) + (1.-zxy) * frld_bdydta(ib,1)
	748	hicif_bdy(ib) = zxy * hicif_bdydta(ib,2) + (1.-zxy) * hicif_bdydta(ib,1)
	749	hsnif_bdy(ib) = zxy * hsnif_bdydta(ib,2) + (1.-zxy) * hsnif_bdydta(ib,1)
	750	END DO
	751	ENDIF ! just if ln_ice_frs is true
	752	#endif
	753
	754	END IF !end if ((nn_dtactl==1).AND.(ntimes_bdy>1))
[911]	755
	756
[1125]	757	! !---------------------!
	758	! ! last call !
	759	! !---------------------!
[911]	760	IF( kt == nitend ) THEN
[2528]	761	IF(ln_tra_frs .or. ln_ice_frs) CALL iom_close( numbdyt ) ! Closing of the 3 files
	762	IF(ln_dyn_frs) CALL iom_close( numbdyu )
	763	IF(ln_dyn_frs) CALL iom_close( numbdyv )
[911]	764	ENDIF
[1125]	765	!
[2528]	766	ENDIF ! ln_dyn_frs .OR. ln_tra_frs
	767	!
	768	END SUBROUTINE bdy_dta_frs
[911]	769
[3432]	770	!FTRANS CLEAR
[911]	771
[2528]	772	SUBROUTINE bdy_dta_fla( kt, jit, icycl )
[911]	773	!!---------------------------------------------------------------------------
[2528]	774	!! * SUBROUTINE bdy_dta_fla *
[911]	775	!!
[1125]	776	!! ** Purpose : Read unstructured boundary data for Flather condition
[911]	777	!!
[1125]	778	!! ** Method : At the first timestep, read in boundary data for two
	779	!! times from the file and time-interpolate. At other
	780	!! timesteps, check to see if we need another time from
	781	!! the file. If so read it in. Time interpolate.
[911]	782	!!---------------------------------------------------------------------------
[1125]	783	!!gm DOCTOR names : argument integer : start with "k"
[911]	784	INTEGER, INTENT( in ) :: kt ! ocean time-step index
	785	INTEGER, INTENT( in ) :: jit ! barotropic time step index
[2528]	786	INTEGER, INTENT( in ) :: icycl ! number of cycles need for final file close
[1125]	787	! ! (for timesplitting option, otherwise zero)
	788	!!
[911]	789	LOGICAL :: lect ! flag for reading
[1125]	790	INTEGER :: it, ib, igrd ! dummy loop indices
[911]	791	INTEGER :: idvar ! netcdf var ID
	792	INTEGER :: iman, i15, imois ! Time variables for monthly clim forcing
[1125]	793	INTEGER :: ntimes_bdyt, ntimes_bdyu, ntimes_bdyv
[911]	794	INTEGER :: itimer, totime
	795	INTEGER :: ipi, ipj, ipk, inum ! temporary integers (NetCDF read)
	796	INTEGER :: iyear0, imonth0, iday0
	797	INTEGER :: ihours0, iminutes0, isec0
[1125]	798	INTEGER :: iyear, imonth, iday, isecs
[911]	799	INTEGER, DIMENSION(jpbtime) :: istept, istepu, istepv ! time arrays from data files
[1125]	800	REAL(wp) :: dayfrac, zxy, zoffsett
	801	REAL(wp) :: zoffsetu, zoffsetv
[1126]	802	REAL(wp) :: dayjul0, zdayjulini
[1125]	803	REAL(wp) :: zinterval_s, zinterval_e ! First and last interval in time axis
	804	REAL(wp), DIMENSION(jpbtime) :: zstepr ! REAL time array from data files
	805	REAL(wp), DIMENSION(jpbdta,1) :: zdta ! temporary array for data fields
[2528]	806	CHARACTER(LEN=80), DIMENSION(6) :: clfile
[911]	807	CHARACTER(LEN=70 ) :: clunits ! units attribute of time coordinate
	808	!!---------------------------------------------------------------------------
	809
[2528]	810	!!gm add here the same style as in bdy_dta_frs
	811	!!gm clearly bdy_dta_fla and bdy_dta_frs can be combined...
[1125]	812	!!gm too many things duplicated in the read of data... simplification can be done
	813
[911]	814	! -------------------- !
	815	! Initialization !
	816	! -------------------- !
	817
	818	lect = .false. ! If true, read a time record
	819
	820	! Some time variables for monthly climatological forcing:
	821	! *******************************************************
[1125]	822	!!gm here use directely daymod variables
[911]	823
	824	iman = INT( raamo ) ! Number of months in a year
	825
[1713]	826	i15 = INT( 2*REAL( nday, wp ) / ( REAL( nmonth_len(nmonth), wp ) + 0.5 ) )
[911]	827	! i15=0 if the current day is in the first half of the month, else i15=1
	828
	829	imois = nmonth + i15 - 1 ! imois is the first month record
	830	IF( imois == 0 ) imois = iman
	831
	832	! Time variable for non-climatological forcing:
	833	! *********************************************
	834
[1241]	835	itimer = ((kt-1)-nit000+1)*rdt ! current time in seconds for interpolation
	836	itimer = itimer + jit*rdt/REAL(nn_baro,wp) ! in non-climatological case
[911]	837
[2528]	838	IF ( ln_tides ) THEN
[911]	839
[1125]	840	! -------------------------------------!
	841	! Update BDY fields with tidal forcing !
	842	! -------------------------------------!
	843
	844	CALL tide_update( kt, jit )
[911]	845
[1125]	846	ENDIF
[911]	847
[2528]	848	IF ( ln_dyn_fla ) THEN
[911]	849
[1125]	850	! -------------------------------------!
	851	! Update BDY fields with model data !
	852	! -------------------------------------!
[911]	853
[1125]	854	! !-------------------!
[2528]	855	IF( kt == nit000 .and. jit ==2 ) THEN ! First call only !
[1125]	856	! !-------------------!
	857	istep_bt(:) = 0
	858	nbdy_b_bt = 0
	859	nbdy_a_bt = 0
[911]	860
[1125]	861	! Get time information from bdy data file
	862	! ***************************************
[911]	863
	864	IF(lwp) WRITE(numout,*)
[2528]	865	IF(lwp) WRITE(numout,*) 'bdy_dta_fla :Initialize unstructured boundary data for barotropic variables.'
[911]	866	IF(lwp) WRITE(numout,*) '~~~~~~~'
	867
[2528]	868	IF( nn_dtactl == 0 ) THEN
[911]	869	IF(lwp) WRITE(numout,*) 'Bdy data are taken from initial conditions'
	870
[2528]	871	ELSEIF (nn_dtactl == 1) THEN
[911]	872	IF(lwp) WRITE(numout,*) 'Bdy data are read in netcdf files'
	873
[1713]	874	dayfrac = adatrj - REAL(itimer,wp)/86400. ! day fraction at time step kt-1
	875	dayfrac = dayfrac - INT (dayfrac) !
	876	totime = (nitend-nit000+1)*rdt ! Total time of the run to verify that all the
	877	! necessary time dumps in file are included
[911]	878
[2528]	879	clfile(4) = cn_dta_fla_T
	880	clfile(5) = cn_dta_fla_U
	881	clfile(6) = cn_dta_fla_V
[911]	882
[2528]	883	DO igrd = 4,6
[911]	884
[1125]	885	CALL iom_open( clfile(igrd), inum )
[2528]	886	CALL iom_gettime( inum, zstepr, kntime=ntimes_bdy_bt, cdunits=clunits )
[911]	887
[1125]	888	SELECT CASE( igrd )
[2528]	889	CASE (4)
	890	numbdyt_bt = inum
	891	CASE (5)
	892	numbdyu_bt = inum
	893	CASE (6)
	894	numbdyv_bt = inum
[1125]	895	END SELECT
	896
[911]	897	! Calculate time offset
	898	READ(clunits,7000) iyear0, imonth0, iday0, ihours0, iminutes0, isec0
	899	! Convert time origin in file to julian days
	900	isec0 = isec0 + ihours060.60. + iminutes0*60.
[1713]	901	CALL ymds2ju(iyear0, imonth0, iday0, REAL(isec0, wp), dayjul0)
[911]	902	! Compute model initialization time
[1125]	903	iyear = ndastp / 10000
	904	imonth = ( ndastp - iyear * 10000 ) / 100
	905	iday = ndastp - iyear * 10000 - imonth * 100
	906	isecs = dayfrac * 86400
[1713]	907	CALL ymds2ju(iyear, imonth, iday, REAL(isecs, wp) , zdayjulini)
[1125]	908	! zoffset from initialization date:
	909	zoffset = (dayjul0-zdayjulini)*86400
[911]	910	!
[3432]	911	! ARP - Cannot have 'string literals' in an input format with the Cray compiler
	912	! so replace with corresponding number of spaces instead.
	913	!7000 FORMAT('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2)
	914	7000 FORMAT(14x, I4.4,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2)
[911]	915
	916	!! TO BE DONE... Check consistency between calendar from file
	917	!! (available optionally from iom_gettime) and calendar in model
	918	!! when calendar in model available outside of IOIPSL.
	919
	920	! Check that there are not too many times in the file.
[1125]	921	IF (ntimes_bdy_bt > jpbtime) CALL ctl_stop( &
	922	'Number of time dumps in bdy file exceed jpbtime parameter', &
	923	'Check file:' // TRIM(clfile(igrd)) )
[911]	924
[1125]	925	! Check that time array increases (or interp will fail):
[2528]	926	DO it = 2, ntimes_bdy_bt
[1125]	927	IF ( zstepr(it-1) >= zstepr(it) ) THEN
	928	CALL ctl_stop('Time array in unstructured boundary data file', &
	929	'does not continuously increase.', &
	930	'Check file:' // TRIM(clfile(igrd)) )
	931	EXIT
	932	END IF
[911]	933	END DO
	934
[2528]	935	IF ( .NOT. ln_clim ) THEN
[1125]	936	! Check that times in file span model run time:
[911]	937
[1125]	938	! Note: the fields may be time means, so we allow nit000 to be before
	939	! first time in the file, provided that it falls inside the meaning
	940	! period of the first field. Until we can get the meaning period
	941	! from the file, use the interval between fields as a proxy.
	942	! If nit000 is before the first time, use the value at first time
	943	! instead of extrapolating. This is done by putting time 1 into
	944	! both time levels.
	945	! The same applies to the last time level: see setting of lect below.
[911]	946
[2528]	947	IF ( ntimes_bdy_bt == 1 ) CALL ctl_stop( &
[1125]	948	'There is only one time dump in data files', &
[2528]	949	'Set ln_clim=.true. in namelist for constant bdy forcing.' )
[911]	950
[1125]	951	zinterval_s = zstepr(2) - zstepr(1)
[2528]	952	zinterval_e = zstepr(ntimes_bdy_bt) - zstepr(ntimes_bdy_bt-1)
[1125]	953
[2528]	954	IF( zstepr(1) + zoffset > 0 ) THEN
	955	WRITE(ctmp1,*) 'Check file: ', clfile(igrd)
	956	CALL ctl_stop( 'First time dump in bdy file is after model initial time', ctmp1 )
[1125]	957	END IF
[2528]	958	IF( zstepr(ntimes_bdy_bt) + zoffset < totime ) THEN
	959	WRITE(ctmp1,*) 'Check file: ', clfile(igrd)
	960	CALL ctl_stop( 'Last time dump in bdy file is before model final time', ctmp1 )
[1125]	961	END IF
[2528]	962	END IF ! .NOT. ln_clim
[1125]	963
[2528]	964	IF ( igrd .EQ. 4) THEN
[1125]	965	ntimes_bdyt = ntimes_bdy_bt
	966	zoffsett = zoffset
	967	istept(:) = INT( zstepr(:) + zoffset )
[2528]	968	ELSE IF (igrd .EQ. 5) THEN
[1125]	969	ntimes_bdyu = ntimes_bdy_bt
	970	zoffsetu = zoffset
	971	istepu(:) = INT( zstepr(:) + zoffset )
[2528]	972	ELSE IF (igrd .EQ. 6) THEN
[1125]	973	ntimes_bdyv = ntimes_bdy_bt
	974	zoffsetv = zoffset
	975	istepv(:) = INT( zstepr(:) + zoffset )
[911]	976	ENDIF
	977
	978	ENDDO
	979
	980	! Verify time consistency between files
	981
[1125]	982	IF ( ntimes_bdyu /= ntimes_bdyt .OR. ntimes_bdyv /= ntimes_bdyt ) THEN
	983	CALL ctl_stop( &
	984	'Time axis lengths differ between bdy data files', &
	985	'Multiple time frequencies not implemented yet' )
	986	ELSE
	987	ntimes_bdy_bt = ntimes_bdyt
[911]	988	ENDIF
	989
[1125]	990	IF (zoffsetu.NE.zoffsett .OR. zoffsetv.NE.zoffsett) THEN
	991	CALL ctl_stop( &
	992	'Bdy data files must have the same time origin', &
	993	'Multiple time frequencies not implemented yet' )
[911]	994	ENDIF
[1125]	995	zoffset = zoffsett
[911]	996
	997	!! Check that times are the same in the three files... HERE.
	998	istep_bt(:) = istept(:)
	999
	1000	! Check number of time dumps:
[2528]	1001	IF (ln_clim) THEN
[1125]	1002	SELECT CASE ( ntimes_bdy_bt )
	1003	CASE( 1 )
[911]	1004	IF(lwp) WRITE(numout,*)
	1005	IF(lwp) WRITE(numout,*) 'We assume constant boundary forcing from bdy data files'
	1006	IF(lwp) WRITE(numout,*)
[1125]	1007	CASE( 12 )
[911]	1008	IF(lwp) WRITE(numout,*)
	1009	IF(lwp) WRITE(numout,*) 'We assume monthly (and cyclic) boundary forcing from bdy data files'
	1010	IF(lwp) WRITE(numout,*)
[1125]	1011	CASE DEFAULT
	1012	CALL ctl_stop( &
[2528]	1013	'For climatological boundary forcing (ln_clim=.true.),',&
[1125]	1014	'bdy data files must contain 1 or 12 time dumps.' )
	1015	END SELECT
[911]	1016	ENDIF
	1017
	1018	! Find index of first record to read (before first model time).
	1019
[1125]	1020	it=1
	1021	DO WHILE ( ((istep_bt(it+1)) <= 0 ).AND.(it.LE.(ntimes_bdy_bt-1)))
	1022	it=it+1
[911]	1023	END DO
[1125]	1024	nbdy_b_bt = it
[911]	1025
[2528]	1026	IF(lwp) WRITE(numout,*) 'Time offset is ',zoffset
	1027	IF(lwp) WRITE(numout,*) 'First record to read is ',nbdy_b_bt
[911]	1028
[2528]	1029	ENDIF ! endif (nn_dtactl == 1)
[911]	1030
[2528]	1031	! 1.2 Read first record in file if necessary (ie if nn_dtactl == 1)
[911]	1032	! *****************************************************************
	1033
[2528]	1034	IF ( nn_dtactl == 0) THEN
[911]	1035	! boundary data arrays are filled with initial conditions
[2528]	1036	igrd = 5 ! U-points data
[1125]	1037	DO ib = 1, nblen(igrd)
	1038	ubtbdy(ib) = un(nbi(ib,igrd), nbj(ib,igrd), 1)
[911]	1039	END DO
	1040
[2528]	1041	igrd = 6 ! V-points data
[1125]	1042	DO ib = 1, nblen(igrd)
	1043	vbtbdy(ib) = vn(nbi(ib,igrd), nbj(ib,igrd), 1)
[911]	1044	END DO
	1045
[2528]	1046	igrd = 4 ! T-points data
[1125]	1047	DO ib = 1, nblen(igrd)
	1048	sshbdy(ib) = sshn(nbi(ib,igrd), nbj(ib,igrd))
[911]	1049	END DO
	1050
[2528]	1051	ELSEIF (nn_dtactl == 1) THEN
[911]	1052
	1053	! Set first record in the climatological case:
[2528]	1054	IF ((ln_clim).AND.(ntimes_bdy_bt==1)) THEN
[1125]	1055	nbdy_a_bt = 1
[2528]	1056	ELSEIF ((ln_clim).AND.(ntimes_bdy_bt==iman)) THEN
[1125]	1057	nbdy_b_bt = 0
	1058	nbdy_a_bt = imois
[911]	1059	ELSE
[1125]	1060	nbdy_a_bt = nbdy_b_bt
[911]	1061	END IF
	1062
	1063	! Open Netcdf files:
	1064
[2528]	1065	CALL iom_open ( cn_dta_fla_T, numbdyt_bt )
	1066	CALL iom_open ( cn_dta_fla_U, numbdyu_bt )
	1067	CALL iom_open ( cn_dta_fla_V, numbdyv_bt )
[911]	1068
	1069	! Read first record:
	1070	ipj=1
[2528]	1071	igrd=4
[1125]	1072	ipi=nblendta(igrd)
[911]	1073
	1074	! ssh
[2528]	1075	igrd=4
[1125]	1076	IF ( nblendta(igrd) .le. 0 ) THEN
	1077	idvar = iom_varid( numbdyt_bt,'sossheig' )
	1078	nblendta(igrd) = iom_file(numbdyt_bt)%dimsz(1,idvar)
[911]	1079	ENDIF
[1125]	1080	WRITE(numout,*) 'Dim size for sossheig is ',nblendta(igrd)
	1081	ipi=nblendta(igrd)
[911]	1082
[1125]	1083	CALL iom_get ( numbdyt_bt, jpdom_unknown,'sossheig',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1084
[1125]	1085	DO ib=1, nblen(igrd)
	1086	sshbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1087	END DO
	1088
	1089	! u-velocity
[2528]	1090	igrd=5
[1125]	1091	IF ( nblendta(igrd) .le. 0 ) THEN
	1092	idvar = iom_varid( numbdyu_bt,'vobtcrtx' )
	1093	nblendta(igrd) = iom_file(numbdyu_bt)%dimsz(1,idvar)
[911]	1094	ENDIF
[1125]	1095	WRITE(numout,*) 'Dim size for vobtcrtx is ',nblendta(igrd)
	1096	ipi=nblendta(igrd)
[911]	1097
[1125]	1098	CALL iom_get ( numbdyu_bt, jpdom_unknown,'vobtcrtx',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1099
[1125]	1100	DO ib=1, nblen(igrd)
	1101	ubtbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1102	END DO
	1103
	1104	! v-velocity
[2528]	1105	igrd=6
[1125]	1106	IF ( nblendta(igrd) .le. 0 ) THEN
	1107	idvar = iom_varid( numbdyv_bt,'vobtcrty' )
	1108	nblendta(igrd) = iom_file(numbdyv_bt)%dimsz(1,idvar)
[911]	1109	ENDIF
[1125]	1110	WRITE(numout,*) 'Dim size for vobtcrty is ',nblendta(igrd)
	1111	ipi=nblendta(igrd)
[911]	1112
[1125]	1113	CALL iom_get ( numbdyv_bt, jpdom_unknown,'vobtcrty',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1114
[1125]	1115	DO ib=1, nblen(igrd)
	1116	vbtbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1117	END DO
	1118
	1119	END IF
	1120
	1121	! In the case of constant boundary forcing fill bdy arrays once for all
[2528]	1122	IF ((ln_clim).AND.(ntimes_bdy_bt==1)) THEN
[911]	1123
	1124	ubtbdy (:) = ubtbdydta (:,2)
	1125	vbtbdy (:) = vbtbdydta (:,2)
	1126	sshbdy (:) = sshbdydta (:,2)
	1127
[1125]	1128	CALL iom_close( numbdyt_bt )
	1129	CALL iom_close( numbdyu_bt )
	1130	CALL iom_close( numbdyv_bt )
[911]	1131
	1132	END IF
	1133
	1134	ENDIF ! End if nit000
	1135
	1136	! -------------------- !
	1137	! 2. At each time step !
	1138	! -------------------- !
	1139
[2528]	1140	IF ((nn_dtactl==1).AND.(ntimes_bdy_bt>1)) THEN
[911]	1141
	1142	! 2.1 Read one more record if necessary
	1143	!**************************************
	1144
[2528]	1145	IF ( (ln_clim).AND.(imois/=nbdy_b_bt) ) THEN ! remember that nbdy_b_bt=0 for kt=nit000
[1125]	1146	nbdy_b_bt = imois
	1147	nbdy_a_bt = imois+1
	1148	nbdy_b_bt = MOD( nbdy_b_bt, iman )
	1149	IF( nbdy_b_bt == 0 ) nbdy_b_bt = iman
	1150	nbdy_a_bt = MOD( nbdy_a_bt, iman )
	1151	IF( nbdy_a_bt == 0 ) nbdy_a_bt = iman
[911]	1152	lect=.true.
	1153
[2528]	1154	ELSEIF ((.NOT.ln_clim).AND.(itimer >= istep_bt(nbdy_a_bt))) THEN
[1125]	1155	nbdy_b_bt=nbdy_a_bt
	1156	nbdy_a_bt=nbdy_a_bt+1
[911]	1157	lect=.true.
	1158	END IF
	1159
	1160	IF (lect) THEN
	1161
	1162	! Swap arrays
	1163	sshbdydta(:,1) = sshbdydta(:,2)
	1164	ubtbdydta(:,1) = ubtbdydta(:,2)
	1165	vbtbdydta(:,1) = vbtbdydta(:,2)
	1166
	1167	! read another set
	1168
	1169	ipj=1
[4474]	1170	ipk=jpkorig
[2528]	1171	igrd=4
[1125]	1172	ipi=nblendta(igrd)
[911]	1173
	1174
	1175	! ssh
[2528]	1176	igrd=4
[1125]	1177	ipi=nblendta(igrd)
[911]	1178
[1125]	1179	CALL iom_get ( numbdyt_bt, jpdom_unknown,'sossheig',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1180
[1125]	1181	DO ib=1, nblen(igrd)
	1182	sshbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1183	END DO
	1184
	1185	! u-velocity
[2528]	1186	igrd=5
[1125]	1187	ipi=nblendta(igrd)
[911]	1188
[1125]	1189	CALL iom_get ( numbdyu_bt, jpdom_unknown,'vobtcrtx',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1190
[1125]	1191	DO ib=1, nblen(igrd)
	1192	ubtbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1193	END DO
	1194
	1195	! v-velocity
[2528]	1196	igrd=6
[1125]	1197	ipi=nblendta(igrd)
[911]	1198
[1125]	1199	CALL iom_get ( numbdyv_bt, jpdom_unknown,'vobtcrty',zdta(1:ipi,1:ipj),nbdy_a_bt )
[911]	1200
[1125]	1201	DO ib=1, nblen(igrd)
	1202	vbtbdydta(ib,2) = zdta(nbmap(ib,igrd),1)
[911]	1203	END DO
	1204
	1205
[2528]	1206	IF(lwp) WRITE(numout,*) 'bdy_dta_fla : first record file used nbdy_b_bt ',nbdy_b_bt
[1125]	1207	IF(lwp) WRITE(numout,*) '~~~~~~~~ last record file used nbdy_a_bt ',nbdy_a_bt
[2528]	1208	IF (.NOT.ln_clim) THEN
[1125]	1209	IF(lwp) WRITE(numout,*) 'first record time (s): ', istep_bt(nbdy_b_bt)
[911]	1210	IF(lwp) WRITE(numout,*) 'model time (s) : ', itimer
[1125]	1211	IF(lwp) WRITE(numout,*) 'second record time (s): ', istep_bt(nbdy_a_bt)
[911]	1212	ENDIF
	1213	END IF ! end lect=.true.
	1214
	1215
	1216	! 2.2 Interpolate linearly:
	1217	! ***************************
	1218
[2528]	1219	IF (ln_clim) THEN
[1713]	1220	zxy = REAL( nday, wp ) / REAL( nmonth_len(nbdy_b_bt), wp ) + 0.5 - i15
[911]	1221	ELSE
[1713]	1222	zxy = REAL(istep_bt(nbdy_b_bt)-itimer, wp) / REAL(istep_bt(nbdy_b_bt)-istep_bt(nbdy_a_bt), wp)
[911]	1223	END IF
	1224
[2528]	1225	igrd=4
[1125]	1226	DO ib=1, nblen(igrd)
	1227	sshbdy(ib) = zxy * sshbdydta(ib,2) + &
	1228	(1.-zxy) * sshbdydta(ib,1)
[911]	1229	END DO
	1230
[2528]	1231	igrd=5
[1125]	1232	DO ib=1, nblen(igrd)
	1233	ubtbdy(ib) = zxy * ubtbdydta(ib,2) + &
	1234	(1.-zxy) * ubtbdydta(ib,1)
[911]	1235	END DO
	1236
[2528]	1237	igrd=6
[1125]	1238	DO ib=1, nblen(igrd)
	1239	vbtbdy(ib) = zxy * vbtbdydta(ib,2) + &
	1240	(1.-zxy) * vbtbdydta(ib,1)
[911]	1241	END DO
	1242
	1243
[2528]	1244	END IF !end if ((nn_dtactl==1).AND.(ntimes_bdy_bt>1))
[911]	1245
	1246	! ------------------- !
	1247	! Last call kt=nitend !
	1248	! ------------------- !
	1249
	1250	! Closing of the 3 files
[2528]	1251	IF( kt == nitend .and. jit == icycl ) THEN
[1125]	1252	CALL iom_close( numbdyt_bt )
	1253	CALL iom_close( numbdyu_bt )
	1254	CALL iom_close( numbdyv_bt )
[911]	1255	ENDIF
	1256
[2528]	1257	ENDIF ! ln_dyn_frs
[911]	1258
[2528]	1259	END SUBROUTINE bdy_dta_fla
[911]	1260
	1261
	1262	#else
[1125]	1263	!!----------------------------------------------------------------------
	1264	!! Dummy module NO Unstruct Open Boundary Conditions
	1265	!!----------------------------------------------------------------------
[911]	1266	CONTAINS
[2528]	1267	SUBROUTINE bdy_dta_frs( kt ) ! Empty routine
	1268	WRITE(,) 'bdy_dta_frs: You should not have seen this print! error?', kt
	1269	END SUBROUTINE bdy_dta_frs
	1270	SUBROUTINE bdy_dta_fla( kt, kit, icycle ) ! Empty routine
	1271	WRITE(,) 'bdy_dta_frs: You should not have seen this print! error?', kt, kit
	1272	END SUBROUTINE bdy_dta_fla
[911]	1273	#endif
	1274
	1275	!!==============================================================================
	1276	END MODULE bdydta

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

Download in other formats: