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.

fldread.F90 in branches/UKMO/dev_r5518_GO6_package_text_diagnostics/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package_text_diagnostics/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90 @ 10770

Last change on this file since 10770 was 10770, checked in by andmirek, 5 years ago
GMED 450 changes for writing to numout
File size: 85.5 KB

Rev	Line
[888]	1	MODULE fldread
	2	!!======================================================================
	3	!! * MODULE fldread *
	4	!! Ocean forcing: read input field for surface boundary condition
	5	!!=====================================================================
[2715]	6	!! History : 2.0 ! 06-2006 (S. Masson, G. Madec) Original code
	7	!! ! 05-2008 (S. Alderson) Modified for Interpolation in memory
	8	!! ! from input grid to model grid
[4230]	9	!! ! 10-2013 (D. Delrosso, P. Oddo) implement suppression of
	10	!! ! land point prior to interpolation
[888]	11	!!----------------------------------------------------------------------
	12
	13	!!----------------------------------------------------------------------
	14	!! fld_read : read input fields used for the computation of the
	15	!! surface boundary condition
	16	!!----------------------------------------------------------------------
	17	USE oce ! ocean dynamics and tracers
	18	USE dom_oce ! ocean space and time domain
	19	USE phycst ! ???
	20	USE in_out_manager ! I/O manager
	21	USE iom ! I/O manager library
[1275]	22	USE geo2ocean ! for vector rotation on to model grid
[2715]	23	USE lib_mpp ! MPP library
[3294]	24	USE wrk_nemo ! work arrays
[4245]	25	USE lbclnk ! ocean lateral boundary conditions (C1D case)
[2715]	26	USE ioipsl, ONLY : ymds2ju, ju2ymds ! for calendar
[4230]	27	USE sbc_oce
	28
[888]	29	IMPLICIT NONE
	30	PRIVATE
[3294]	31
	32	PUBLIC fld_map ! routine called by tides_init
[3851]	33	PUBLIC fld_read, fld_fill ! called by sbc... modules
[6487]	34	PUBLIC fld_clopn
[888]	35
	36	TYPE, PUBLIC :: FLD_N !: Namelist field informations
[1730]	37	CHARACTER(len = 256) :: clname ! generic name of the NetCDF flux file
[4245]	38	REAL(wp) :: nfreqh ! frequency of each flux file
[1730]	39	CHARACTER(len = 34) :: clvar ! generic name of the variable in the NetCDF flux file
	40	LOGICAL :: ln_tint ! time interpolation or not (T/F)
	41	LOGICAL :: ln_clim ! climatology or not (T/F)
[2528]	42	CHARACTER(len = 8) :: cltype ! type of data file 'daily', 'monthly' or yearly'
[4663]	43	CHARACTER(len = 256) :: wname ! generic name of a NetCDF weights file to be used, blank if not
[1730]	44	CHARACTER(len = 34) :: vcomp ! symbolic component name if a vector that needs rotation
[2715]	45	! ! a string starting with "U" or "V" for each component
	46	! ! chars 2 onwards identify which components go together
[4230]	47	CHARACTER(len = 34) :: lname ! generic name of a NetCDF land/sea mask file to be used, blank if not
	48	! ! 0=sea 1=land
[888]	49	END TYPE FLD_N
	50
	51	TYPE, PUBLIC :: FLD !: Input field related variables
	52	CHARACTER(len = 256) :: clrootname ! generic name of the NetCDF file
	53	CHARACTER(len = 256) :: clname ! current name of the NetCDF file
[4245]	54	REAL(wp) :: nfreqh ! frequency of each flux file
[888]	55	CHARACTER(len = 34) :: clvar ! generic name of the variable in the NetCDF flux file
	56	LOGICAL :: ln_tint ! time interpolation or not (T/F)
[1132]	57	LOGICAL :: ln_clim ! climatology or not (T/F)
[2528]	58	CHARACTER(len = 8) :: cltype ! type of data file 'daily', 'monthly' or yearly'
[1132]	59	INTEGER :: num ! iom id of the jpfld files to be read
[1730]	60	INTEGER , DIMENSION(2) :: nrec_b ! before record (1: index, 2: second since Jan. 1st 00h of nit000 year)
	61	INTEGER , DIMENSION(2) :: nrec_a ! after record (1: index, 2: second since Jan. 1st 00h of nit000 year)
[2528]	62	REAL(wp) , ALLOCATABLE, DIMENSION(:,:,: ) :: fnow ! input fields interpolated to now time step
	63	REAL(wp) , ALLOCATABLE, DIMENSION(:,:,:,:) :: fdta ! 2 consecutive record of input fields
[1275]	64	CHARACTER(len = 256) :: wgtname ! current name of the NetCDF weight file acting as a key
[2715]	65	! ! into the WGTLIST structure
[1275]	66	CHARACTER(len = 34) :: vcomp ! symbolic name for a vector component that needs rotation
[3851]	67	LOGICAL, DIMENSION(2) :: rotn ! flag to indicate whether before/after field has been rotated
	68	INTEGER :: nreclast ! last record to be read in the current file
[4230]	69	CHARACTER(len = 256) :: lsmname ! current name of the NetCDF mask file acting as a key
[888]	70	END TYPE FLD
	71
[5132]	72	TYPE, PUBLIC :: MAP_POINTER !: Map from input data file to local domain
	73	INTEGER, POINTER, DIMENSION(:) :: ptr ! Array of integer pointers to 1D arrays
	74	LOGICAL :: ll_unstruc ! Unstructured (T) or structured (F) boundary data file
[3294]	75	END TYPE MAP_POINTER
	76
[1275]	77	!$AGRIF_DO_NOT_TREAT
	78
	79	!! keep list of all weights variables so they're only read in once
	80	!! need to add AGRIF directives not to process this structure
	81	!! also need to force wgtname to include AGRIF nest number
	82	TYPE :: WGT !: Input weights related variables
	83	CHARACTER(len = 256) :: wgtname ! current name of the NetCDF weight file
	84	INTEGER , DIMENSION(2) :: ddims ! shape of input grid
	85	INTEGER , DIMENSION(2) :: botleft ! top left corner of box in input grid containing
[2715]	86	! ! current processor grid
[1275]	87	INTEGER , DIMENSION(2) :: topright ! top right corner of box
	88	INTEGER :: jpiwgt ! width of box on input grid
	89	INTEGER :: jpjwgt ! height of box on input grid
	90	INTEGER :: numwgt ! number of weights (4=bilinear, 16=bicubic)
	91	INTEGER :: nestid ! for agrif, keep track of nest we're in
[2528]	92	INTEGER :: overlap ! =0 when cyclic grid has no overlapping EW columns
[2715]	93	! ! =>1 when they have one or more overlapping columns
	94	! ! =-1 not cyclic
[1275]	95	LOGICAL :: cyclic ! east-west cyclic or not
[2528]	96	INTEGER, DIMENSION(:,:,:), POINTER :: data_jpi ! array of source integers
	97	INTEGER, DIMENSION(:,:,:), POINTER :: data_jpj ! array of source integers
[1275]	98	REAL(wp), DIMENSION(:,:,:), POINTER :: data_wgt ! array of weights on model grid
[2528]	99	REAL(wp), DIMENSION(:,:,:), POINTER :: fly_dta ! array of values on input grid
	100	REAL(wp), DIMENSION(:,:,:), POINTER :: col ! temporary array for reading in columns
[1275]	101	END TYPE WGT
	102
	103	INTEGER, PARAMETER :: tot_wgts = 10
	104	TYPE( WGT ), DIMENSION(tot_wgts) :: ref_wgts ! array of wgts
	105	INTEGER :: nxt_wgt = 1 ! point to next available space in ref_wgts array
[4230]	106	REAL(wp), PARAMETER :: undeff_lsm = -999.00_wp
[1275]	107
	108	!$AGRIF_END_DO_NOT_TREAT
	109
[888]	110	!!----------------------------------------------------------------------
[2528]	111	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[1156]	112	!! $Id$
[2715]	113	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[888]	114	!!----------------------------------------------------------------------
	115	CONTAINS
	116
[3851]	117	SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset )
[888]	118	!!---------------------------------------------------------------------
	119	!! * ROUTINE fld_read *
	120	!!
	121	!! ** Purpose : provide at each time step the surface ocean fluxes
	122	!! (momentum, heat, freshwater and runoff)
	123	!!
	124	!! ** Method : READ each input fields in NetCDF files using IOM
	125	!! and intepolate it to the model time-step.
	126	!! Several assumptions are made on the input file:
	127	!! blahblahblah....
	128	!!----------------------------------------------------------------------
	129	INTEGER , INTENT(in ) :: kt ! ocean time step
[1132]	130	INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step)
[888]	131	TYPE(FLD), INTENT(inout), DIMENSION(:) :: sd ! input field related variables
[3851]	132	TYPE(MAP_POINTER),INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping indices
	133	INTEGER , INTENT(in ), OPTIONAL :: kit ! subcycle timestep for timesplitting option
	134	INTEGER , INTENT(in ), OPTIONAL :: kt_offset ! provide fields at time other than "now"
	135	! kt_offset = -1 => fields at "before" time level
	136	! kt_offset = +1 => fields at "after" time level
	137	! etc.
[888]	138	!!
[3851]	139	INTEGER :: itmp ! temporary variable
[2528]	140	INTEGER :: imf ! size of the structure sd
[1132]	141	INTEGER :: jf ! dummy indices
[1730]	142	INTEGER :: isecend ! number of second since Jan. 1st 00h of nit000 year at nitend
[2323]	143	INTEGER :: isecsbc ! number of seconds between Jan. 1st 00h of nit000 year and the middle of sbc time step
[3851]	144	INTEGER :: it_offset ! local time offset variable
[1628]	145	LOGICAL :: llnxtyr ! open next year file?
	146	LOGICAL :: llnxtmth ! open next month file?
	147	LOGICAL :: llstop ! stop is the file does not exist
[3294]	148	LOGICAL :: ll_firstcall ! true if this is the first call to fld_read for this set of fields
[1132]	149	REAL(wp) :: ztinta ! ratio applied to after records when doing time interpolation
	150	REAL(wp) :: ztintb ! ratio applied to before records when doing time interpolation
[1191]	151	CHARACTER(LEN=1000) :: clfmt ! write format
[3851]	152	TYPE(MAP_POINTER) :: imap ! global-to-local mapping indices
[888]	153	!!---------------------------------------------------------------------
[3851]	154	ll_firstcall = kt == nit000
	155	IF( PRESENT(kit) ) ll_firstcall = ll_firstcall .and. kit == 1
[3294]	156
[5407]	157	IF ( nn_components == jp_iam_sas ) THEN ; it_offset = nn_fsbc
	158	ELSE ; it_offset = 0
	159	ENDIF
[3851]	160	IF( PRESENT(kt_offset) ) it_offset = kt_offset
	161
	162	imap%ptr => NULL()
	163
[2528]	164	! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar
[3851]	165	IF( present(kit) ) THEN ! ignore kn_fsbc in this case
	166	isecsbc = nsec_year + nsec1jan000 + (kit+it_offset)*NINT( rdt/REAL(nn_baro,wp) )
	167	ELSE ! middle of sbc time step
	168	isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + it_offset * NINT(rdttra(1))
[3294]	169	ENDIF
[1275]	170	imf = SIZE( sd )
[2323]	171	!
[3294]	172	IF( ll_firstcall ) THEN ! initialization
[3851]	173	DO jf = 1, imf
	174	IF( PRESENT(map) ) imap = map(jf)
	175	CALL fld_init( kn_fsbc, sd(jf), imap ) ! read each before field (put them in after as they will be swapped)
	176	END DO
[2528]	177	IF( lwp ) CALL wgt_print() ! control print
	178	ENDIF
	179	! ! ====================================== !
	180	IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN ! update field at each kn_fsbc time-step !
	181	! ! ====================================== !
[888]	182	!
[2528]	183	DO jf = 1, imf ! --- loop over field --- !
	184
[3851]	185	IF( isecsbc > sd(jf)%nrec_a(2) .OR. ll_firstcall ) THEN ! read/update the after data?
[888]	186
[3851]	187	IF( PRESENT(map) ) imap = map(jf) ! temporary definition of map
[1132]	188
[3851]	189	sd(jf)%nrec_b(:) = sd(jf)%nrec_a(:) ! swap before record informations
	190	sd(jf)%rotn(1) = sd(jf)%rotn(2) ! swap before rotate informations
	191	IF( sd(jf)%ln_tint ) sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! swap before record field
	192
	193	CALL fld_rec( kn_fsbc, sd(jf), kt_offset = it_offset, kit = kit ) ! update after record informations
	194
	195	! if kn_fsbc*rdttra is larger than nfreqh (which is kind of odd),
	196	! it is possible that the before value is no more the good one... we have to re-read it
	197	! if before is not the last record of the file currently opened and after is the first record to be read
	198	! in a new file which means after = 1 (the file to be opened corresponds to the current time)
	199	! or after = nreclast + 1 (the file to be opened corresponds to a future time step)
	200	IF( .NOT. ll_firstcall .AND. sd(jf)%ln_tint .AND. sd(jf)%nrec_b(1) /= sd(jf)%nreclast &
	201	& .AND. MOD( sd(jf)%nrec_a(1), sd(jf)%nreclast ) == 1 ) THEN
	202	itmp = sd(jf)%nrec_a(1) ! temporary storage
	203	sd(jf)%nrec_a(1) = sd(jf)%nreclast ! read the last record of the file currently opened
	204	CALL fld_get( sd(jf), imap ) ! read after data
	205	sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field
	206	sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations
[4245]	207	sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)%nfreqh * 3600 ) ! assume freq to be in hours in this case
[3851]	208	sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations
	209	sd(jf)%nrec_a(1) = itmp ! move back to after record
[3294]	210	ENDIF
[1132]	211
[3851]	212	CALL fld_clopn( sd(jf) ) ! Do we need to open a new year/month/week/day file?
	213
[2528]	214	IF( sd(jf)%ln_tint ) THEN
[3851]	215
	216	! if kn_fsbc*rdttra is larger than nfreqh (which is kind of odd),
	217	! it is possible that the before value is no more the good one... we have to re-read it
	218	! if before record is not just just before the after record...
	219	IF( .NOT. ll_firstcall .AND. MOD( sd(jf)%nrec_a(1), sd(jf)%nreclast ) /= 1 &
	220	& .AND. sd(jf)%nrec_b(1) /= sd(jf)%nrec_a(1) - 1 ) THEN
	221	sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) - 1 ! move back to before record
	222	CALL fld_get( sd(jf), imap ) ! read after data
	223	sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field
	224	sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations
[4245]	225	sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)%nfreqh * 3600 ) ! assume freq to be in hours in this case
[3851]	226	sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations
	227	sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) + 1 ! move back to after record
	228	ENDIF
	229
	230	! do we have to change the year/month/week/day of the forcing field??
[2528]	231	! if we do time interpolation we will need to open next year/month/week/day file before the end of the current
	232	! one. If so, we are still before the end of the year/month/week/day when calling fld_rec so sd(jf)%nrec_a(1)
	233	! will be larger than the record number that should be read for current year/month/week/day
	234	! do we need next file data?
[3851]	235	IF( sd(jf)%nrec_a(1) > sd(jf)%nreclast ) THEN
	236
	237	sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) - sd(jf)%nreclast !
	238
	239	IF( .NOT. ( sd(jf)%ln_clim .AND. sd(jf)%cltype == 'yearly' ) ) THEN ! close/open the current/new file
	240
[2528]	241	llnxtmth = sd(jf)%cltype == 'monthly' .OR. nday == nmonth_len(nmonth) ! open next month file?
	242	llnxtyr = sd(jf)%cltype == 'yearly' .OR. (nmonth == 12 .AND. llnxtmth) ! open next year file?
[1132]	243
[2528]	244	! if the run finishes at the end of the current year/month/week/day, we will allow next
	245	! year/month/week/day file to be not present. If the run continue further than the current
	246	! year/month/week/day, next year/month/week/day file must exist
	247	isecend = nsec_year + nsec1jan000 + (nitend - kt) * NINT(rdttra(1)) ! second at the end of the run
	248	llstop = isecend > sd(jf)%nrec_a(2) ! read more than 1 record of next year
[3851]	249	! we suppose that the date of next file is next day (should be ok even for weekly files...)
[2528]	250	CALL fld_clopn( sd(jf), nyear + COUNT((/llnxtyr /)) , &
	251	& nmonth + COUNT((/llnxtmth/)) - 12 * COUNT((/llnxtyr /)), &
	252	& nday + 1 - nmonth_len(nmonth) * COUNT((/llnxtmth/)), llstop )
	253
	254	IF( sd(jf)%num <= 0 .AND. .NOT. llstop ) THEN ! next year file does not exist
	255	CALL ctl_warn('next year/month/week/day file: '//TRIM(sd(jf)%clname)// &
	256	& ' not present -> back to current year/month/day')
[3851]	257	CALL fld_clopn( sd(jf) ) ! back to the current year/month/day
	258	sd(jf)%nrec_a(1) = sd(jf)%nreclast ! force to read the last record in the current year file
[2528]	259	ENDIF
[3851]	260
[1132]	261	ENDIF
[3851]	262	ENDIF ! open need next file?
	263
	264	ENDIF ! temporal interpolation?
[1132]	265
[2528]	266	! read after data
[3851]	267	CALL fld_get( sd(jf), imap )
[2528]	268
[3851]	269	ENDIF ! read new data?
[2528]	270	END DO ! --- end loop over field --- !
[1132]	271
[3851]	272	CALL fld_rot( kt, sd ) ! rotate vector before/now/after fields if needed
[888]	273
[2528]	274	DO jf = 1, imf ! --- loop over field --- !
[888]	275	!
[2528]	276	IF( sd(jf)%ln_tint ) THEN ! temporal interpolation
[10770]	277	IF(lwp .AND. nprint > 1 ) THEN
[4245]	278	clfmt = "('fld_read: var ', a, ' kt = ', i8, ' (', f9.4,' days), Y/M/D = ', i4.4,'/', i2.2,'/', i2.2," // &
	279	& "', records b/a: ', i6.4, '/', i6.4, ' (days ', f9.4,'/', f9.4, ')')"
[3294]	280	WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, &
[1730]	281	& sd(jf)%nrec_b(1), sd(jf)%nrec_a(1), REAL(sd(jf)%nrec_b(2),wp)/rday, REAL(sd(jf)%nrec_a(2),wp)/rday
[3851]	282	WRITE(numout, *) 'it_offset is : ',it_offset
[1191]	283	ENDIF
[2528]	284	! temporal interpolation weights
[2323]	285	ztinta = REAL( isecsbc - sd(jf)%nrec_b(2), wp ) / REAL( sd(jf)%nrec_a(2) - sd(jf)%nrec_b(2), wp )
[1132]	286	ztintb = 1. - ztinta
[888]	287	!CDIR COLLAPSE
[2528]	288	sd(jf)%fnow(:,:,:) = ztintb * sd(jf)%fdta(:,:,:,1) + ztinta * sd(jf)%fdta(:,:,:,2)
	289	ELSE ! nothing to do...
[10770]	290	IF(lwp .AND. nprint > 1 ) THEN
[4245]	291	clfmt = "('fld_read: var ', a, ' kt = ', i8,' (', f9.4,' days), Y/M/D = ', i4.4,'/', i2.2,'/', i2.2," // &
	292	& "', record: ', i6.4, ' (days ', f9.4, ' <-> ', f9.4, ')')"
[2528]	293	WRITE(numout, clfmt) TRIM(sd(jf)%clvar), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, &
	294	& sd(jf)%nrec_a(1), REAL(sd(jf)%nrec_b(2),wp)/rday, REAL(sd(jf)%nrec_a(2),wp)/rday
[1191]	295	ENDIF
[888]	296	ENDIF
	297	!
[2528]	298	IF( kt == nitend - kn_fsbc + 1 ) CALL iom_close( sd(jf)%num ) ! Close the input files
[1132]	299
[2528]	300	END DO ! --- end loop over field --- !
	301	!
	302	! ! ====================================== !
	303	ENDIF ! update field at each kn_fsbc time-step !
	304	! ! ====================================== !
	305	!
[888]	306	END SUBROUTINE fld_read
	307
	308
[3294]	309	SUBROUTINE fld_init( kn_fsbc, sdjf, map )
[888]	310	!!---------------------------------------------------------------------
[1132]	311	!! * ROUTINE fld_init *
	312	!!
[4371]	313	!! ** Purpose : - first call to fld_rec to define before values
	314	!! - if time interpolation, read before data
[1132]	315	!!----------------------------------------------------------------------
[2528]	316	INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step)
	317	TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables
[3851]	318	TYPE(MAP_POINTER),INTENT(in) :: map ! global-to-local mapping indices
[1132]	319	!!
[2528]	320	LOGICAL :: llprevyr ! are we reading previous year file?
	321	LOGICAL :: llprevmth ! are we reading previous month file?
	322	LOGICAL :: llprevweek ! are we reading previous week file?
	323	LOGICAL :: llprevday ! are we reading previous day file?
	324	LOGICAL :: llprev ! llprevyr .OR. llprevmth .OR. llprevweek .OR. llprevday
	325	INTEGER :: idvar ! variable id
	326	INTEGER :: inrec ! number of record existing for this variable
	327	INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd
	328	INTEGER :: isec_week ! number of seconds since start of the weekly file
[1191]	329	CHARACTER(LEN=1000) :: clfmt ! write format
[1132]	330	!!---------------------------------------------------------------------
[2528]	331	llprevyr = .FALSE.
	332	llprevmth = .FALSE.
	333	llprevweek = .FALSE.
	334	llprevday = .FALSE.
	335	isec_week = 0
[1132]	336
	337	! define record informations
[2528]	338	CALL fld_rec( kn_fsbc, sdjf, ldbefore = .TRUE. ) ! return before values in sdjf%nrec_a (as we will swap it later)
[1132]	339
[2528]	340	! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar
[2323]	341
[1132]	342	IF( sdjf%ln_tint ) THEN ! we need to read the previous record and we will put it in the current record structure
[2528]	343
	344	IF( sdjf%nrec_a(1) == 0 ) THEN ! we redefine record sdjf%nrec_a(1) with the last record of previous year file
	345	IF ( sdjf%nfreqh == -12 ) THEN ! yearly mean
	346	IF( sdjf%cltype == 'yearly' ) THEN ! yearly file
	347	sdjf%nrec_a(1) = 1 ! force to read the unique record
	348	llprevyr = .NOT. sdjf%ln_clim ! use previous year file?
	349	ELSE
[3851]	350	CALL ctl_stop( "fld_init: yearly mean file must be in a yearly type of file: "//TRIM(sdjf%clrootname) )
[2528]	351	ENDIF
	352	ELSEIF( sdjf%nfreqh == -1 ) THEN ! monthly mean
	353	IF( sdjf%cltype == 'monthly' ) THEN ! monthly file
	354	sdjf%nrec_a(1) = 1 ! force to read the unique record
	355	llprevmth = .TRUE. ! use previous month file?
[1628]	356	llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file?
[2528]	357	ELSE ! yearly file
	358	sdjf%nrec_a(1) = 12 ! force to read december mean
[1628]	359	llprevyr = .NOT. sdjf%ln_clim ! use previous year file?
	360	ENDIF
[2528]	361	ELSE ! higher frequency mean (in hours)
	362	IF ( sdjf%cltype == 'monthly' ) THEN ! monthly file
[4245]	363	sdjf%nrec_a(1) = NINT( 24 * nmonth_len(nmonth-1) / sdjf%nfreqh ) ! last record of previous month
[2528]	364	llprevmth = .TRUE. ! use previous month file?
[1628]	365	llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file?
[2528]	366	ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ! weekly file
	367	llprevweek = .TRUE. ! use previous week file?
[4245]	368	sdjf%nrec_a(1) = NINT( 24 * 7 / sdjf%nfreqh ) ! last record of previous week
[2528]	369	isec_week = NINT(rday) * 7 ! add a shift toward previous week
	370	ELSEIF( sdjf%cltype == 'daily' ) THEN ! daily file
[4245]	371	sdjf%nrec_a(1) = NINT( 24 / sdjf%nfreqh ) ! last record of previous day
[2528]	372	llprevday = .TRUE. ! use previous day file?
[1628]	373	llprevmth = llprevday .AND. nday == 1 ! use previous month file?
	374	llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file?
[2528]	375	ELSE ! yearly file
[4245]	376	sdjf%nrec_a(1) = NINT( 24 * nyear_len(0) / sdjf%nfreqh ) ! last record of previous year
[1628]	377	llprevyr = .NOT. sdjf%ln_clim ! use previous year file?
[1132]	378	ENDIF
	379	ENDIF
	380	ENDIF
[3851]	381	!
[2528]	382	IF ( sdjf%cltype(1:4) == 'week' ) THEN
	383	isec_week = isec_week + ksec_week( sdjf%cltype(6:8) ) ! second since the beginning of the week
	384	llprevmth = isec_week > nsec_month ! longer time since the beginning of the week than the month
	385	llprevyr = llprevmth .AND. nmonth == 1
	386	ENDIF
	387	llprev = llprevyr .OR. llprevmth .OR. llprevweek .OR. llprevday
	388	!
	389	iyear = nyear - COUNT((/llprevyr /))
	390	imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /))
	391	iday = nday - COUNT((/llprevday/)) + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday)
	392	!
	393	CALL fld_clopn( sdjf, iyear, imonth, iday, .NOT. llprev )
[1132]	394
[1628]	395	! if previous year/month/day file does not exist, we switch to the current year/month/day
[1818]	396	IF( llprev .AND. sdjf%num <= 0 ) THEN
[3851]	397	CALL ctl_warn( 'previous year/month/week/day file: '//TRIM(sdjf%clrootname)// &
[2528]	398	& ' not present -> back to current year/month/week/day' )
[1628]	399	! we force to read the first record of the current year/month/day instead of last record of previous year/month/day
[2528]	400	llprev = .FALSE.
	401	sdjf%nrec_a(1) = 1
[3851]	402	CALL fld_clopn( sdjf )
[1132]	403	ENDIF
	404
[3851]	405	IF( llprev ) THEN ! check if the record sdjf%nrec_a(1) exists in the file
[1132]	406	idvar = iom_varid( sdjf%num, sdjf%clvar ) ! id of the variable sdjf%clvar
	407	IF( idvar <= 0 ) RETURN
	408	inrec = iom_file( sdjf%num )%dimsz( iom_file( sdjf%num )%ndims(idvar), idvar ) ! size of the last dim of idvar
[2528]	409	sdjf%nrec_a(1) = MIN( sdjf%nrec_a(1), inrec ) ! make sure we select an existing record
[1132]	410	ENDIF
	411
[3851]	412	! read before data in after arrays(as we will swap it later)
	413	CALL fld_get( sdjf, map )
[1132]	414
[4245]	415	clfmt = "('fld_init : time-interpolation for ', a, ' read previous record = ', i6, ' at time = ', f7.2, ' days')"
[10770]	416	IF(lwp .AND. nprint > 0) WRITE(numout, clfmt) TRIM(sdjf%clvar), sdjf%nrec_a(1), REAL(sdjf%nrec_a(2),wp)/rday
[1132]	417
	418	ENDIF
[2715]	419	!
[1132]	420	END SUBROUTINE fld_init
	421
	422
[3851]	423	SUBROUTINE fld_rec( kn_fsbc, sdjf, ldbefore, kit, kt_offset )
[1132]	424	!!---------------------------------------------------------------------
[888]	425	!! * ROUTINE fld_rec *
	426	!!
[2528]	427	!! ** Purpose : Compute
	428	!! if sdjf%ln_tint = .TRUE.
	429	!! nrec_a: record number and its time (nrec_b is obtained from nrec_a when swapping)
	430	!! if sdjf%ln_tint = .FALSE.
	431	!! nrec_a(1): record number
	432	!! nrec_b(2) and nrec_a(2): time of the beginning and end of the record (for print only)
[888]	433	!!----------------------------------------------------------------------
[2528]	434	INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step)
	435	TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables
	436	LOGICAL , INTENT(in ), OPTIONAL :: ldbefore ! sent back before record values (default = .FALSE.)
[3851]	437	INTEGER , INTENT(in ), OPTIONAL :: kit ! index of barotropic subcycle
[2528]	438	! used only if sdjf%ln_tint = .TRUE.
[3851]	439	INTEGER , INTENT(in ), OPTIONAL :: kt_offset ! Offset of required time level compared to "now"
	440	! time level in units of time steps.
[888]	441	!!
[2528]	442	LOGICAL :: llbefore ! local definition of ldbefore
	443	INTEGER :: iendrec ! end of this record (in seconds)
	444	INTEGER :: imth ! month number
	445	INTEGER :: ifreq_sec ! frequency mean (in seconds)
	446	INTEGER :: isec_week ! number of seconds since the start of the weekly file
[3851]	447	INTEGER :: it_offset ! local time offset variable
[1132]	448	REAL(wp) :: ztmp ! temporary variable
[888]	449	!!----------------------------------------------------------------------
	450	!
[2528]	451	! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar
[2323]	452	!
[2528]	453	IF( PRESENT(ldbefore) ) THEN ; llbefore = ldbefore .AND. sdjf%ln_tint ! needed only if sdjf%ln_tint = .TRUE.
	454	ELSE ; llbefore = .FALSE.
	455	ENDIF
	456	!
[5407]	457	IF ( nn_components == jp_iam_sas ) THEN ; it_offset = nn_fsbc
	458	ELSE ; it_offset = 0
	459	ENDIF
[3851]	460	IF( PRESENT(kt_offset) ) it_offset = kt_offset
	461	IF( PRESENT(kit) ) THEN ; it_offset = ( kit + it_offset ) * NINT( rdt/REAL(nn_baro,wp) )
	462	ELSE ; it_offset = it_offset * NINT( rdttra(1) )
	463	ENDIF
[3294]	464	!
[2528]	465	! ! =========== !
	466	IF ( sdjf%nfreqh == -12 ) THEN ! yearly mean
	467	! ! =========== !
[888]	468	!
[1132]	469	IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record
	470	!
	471	! INT( ztmp )
	472	! /\|\
	473	! 1 \| *----
	474	! 0 \|----(
	475	! \|----+----\|--> time
[2528]	476	! 0 /\|\ 1 (nday/nyear_len(1))
	477	! \|
	478	! \|
	479	! forcing record : 1
	480	!
[4784]	481	ztmp = REAL( nsec_year, wp ) / ( REAL( nyear_len(1), wp ) * rday ) + 0.5 &
	482	& + REAL( it_offset, wp ) / ( REAL( nyear_len(1), wp ) * rday )
[2528]	483	sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/))
	484	! swap at the middle of the year
[4784]	485	IF( llbefore ) THEN ; sdjf%nrec_a(2) = nsec1jan000 - (1 - INT(ztmp)) * NINT(0.5 * rday) * nyear_len(0) + &
	486	& INT(ztmp) * NINT( 0.5 * rday) * nyear_len(1)
	487	ELSE ; sdjf%nrec_a(2) = nsec1jan000 + (1 - INT(ztmp)) * NINT(0.5 * rday) * nyear_len(1) + &
	488	& INT(ztmp) * INT(rday) * nyear_len(1) + INT(ztmp) * NINT( 0.5 * rday) * nyear_len(2)
[2528]	489	ENDIF
	490	ELSE ! no time interpolation
	491	sdjf%nrec_a(1) = 1
	492	sdjf%nrec_a(2) = NINT(rday) * nyear_len(1) + nsec1jan000 ! swap at the end of the year
	493	sdjf%nrec_b(2) = nsec1jan000 ! beginning of the year (only for print)
	494	ENDIF
	495	!
	496	! ! ============ !
	497	ELSEIF( sdjf%nfreqh == -1 ) THEN ! monthly mean !
	498	! ! ============ !
	499	!
	500	IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record
	501	!
	502	! INT( ztmp )
	503	! /\|\
	504	! 1 \| *----
	505	! 0 \|----(
	506	! \|----+----\|--> time
[1132]	507	! 0 /\|\ 1 (nday/nmonth_len(nmonth))
	508	! \|
	509	! \|
	510	! forcing record : nmonth
	511	!
[4784]	512	ztmp = REAL( nsec_month, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday ) + 0.5 &
	513	& + REAL( it_offset, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday )
[2528]	514	imth = nmonth + INT( ztmp ) - COUNT((/llbefore/))
	515	IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/))
	516	ELSE ; sdjf%nrec_a(1) = imth
	517	ENDIF
	518	sdjf%nrec_a(2) = nmonth_half( imth ) + nsec1jan000 ! swap at the middle of the month
	519	ELSE ! no time interpolation
	520	IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nrec_a(1) = 1
	521	ELSE ; sdjf%nrec_a(1) = nmonth
	522	ENDIF
	523	sdjf%nrec_a(2) = nmonth_end(nmonth ) + nsec1jan000 ! swap at the end of the month
	524	sdjf%nrec_b(2) = nmonth_end(nmonth-1) + nsec1jan000 ! beginning of the month (only for print)
[888]	525	ENDIF
	526	!
[2528]	527	! ! ================================ !
	528	ELSE ! higher frequency mean (in hours)
	529	! ! ================================ !
[888]	530	!
[4245]	531	ifreq_sec = NINT( sdjf%nfreqh * 3600 ) ! frequency mean (in seconds)
[2528]	532	IF( sdjf%cltype(1:4) == 'week' ) isec_week = ksec_week( sdjf%cltype(6:8) ) ! since the first day of the current week
[1132]	533	! number of second since the beginning of the file
[2528]	534	IF( sdjf%cltype == 'monthly' ) THEN ; ztmp = REAL(nsec_month,wp) ! since the first day of the current month
	535	ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; ztmp = REAL(isec_week ,wp) ! since the first day of the current week
	536	ELSEIF( sdjf%cltype == 'daily' ) THEN ; ztmp = REAL(nsec_day ,wp) ! since 00h of the current day
	537	ELSE ; ztmp = REAL(nsec_year ,wp) ! since 00h on Jan 1 of the current year
[1132]	538	ENDIF
[3851]	539	ztmp = ztmp + 0.5 * REAL(kn_fsbc - 1, wp) * rdttra(1) + REAL( it_offset, wp ) ! centrered in the middle of sbc time step
	540	ztmp = ztmp + 0.01 * rdttra(1) ! avoid truncation error
[1132]	541	IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record
	542	!
[3851]	543	! INT( ztmp/ifreq_sec + 0.5 )
[1132]	544	! /\|\
	545	! 2 \| *-----(
	546	! 1 \| *-----(
	547	! 0 \|--(
	548	! \|--+--\|--+--\|--+--\|--> time
[3851]	549	! 0 /\|\ 1 /\|\ 2 /\|\ 3 (ztmp/ifreq_sec)
[1132]	550	! \| \| \|
	551	! \| \| \|
	552	! forcing record : 1 2 3
	553	!
[2528]	554	ztmp= ztmp / REAL(ifreq_sec, wp) + 0.5
	555	ELSE ! no time interpolation
[1132]	556	!
[3851]	557	! INT( ztmp/ifreq_sec )
[1132]	558	! /\|\
	559	! 2 \| *-----(
	560	! 1 \| *-----(
	561	! 0 \|-----(
	562	! \|--+--\|--+--\|--+--\|--> time
[3851]	563	! 0 /\|\ 1 /\|\ 2 /\|\ 3 (ztmp/ifreq_sec)
[1132]	564	! \| \| \|
	565	! \| \| \|
	566	! forcing record : 1 2 3
	567	!
[2528]	568	ztmp= ztmp / REAL(ifreq_sec, wp)
[1132]	569	ENDIF
[3851]	570	sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) ! record number to be read
[1132]	571
[2528]	572	iendrec = ifreq_sec * sdjf%nrec_a(1) + nsec1jan000 ! end of this record (in second)
	573	! add the number of seconds between 00h Jan 1 and the end of previous month/week/day (ok if nmonth=1)
	574	IF( sdjf%cltype == 'monthly' ) iendrec = iendrec + NINT(rday) * SUM(nmonth_len(1:nmonth -1))
	575	IF( sdjf%cltype(1:4) == 'week' ) iendrec = iendrec + ( nsec_year - isec_week )
	576	IF( sdjf%cltype == 'daily' ) iendrec = iendrec + NINT(rday) * ( nday_year - 1 )
	577	IF( sdjf%ln_tint ) THEN
	578	sdjf%nrec_a(2) = iendrec - ifreq_sec / 2 ! swap at the middle of the record
	579	ELSE
	580	sdjf%nrec_a(2) = iendrec ! swap at the end of the record
	581	sdjf%nrec_b(2) = iendrec - ifreq_sec ! beginning of the record (only for print)
	582	ENDIF
[888]	583	!
	584	ENDIF
	585	!
[1132]	586	END SUBROUTINE fld_rec
	587
	588
[3294]	589	SUBROUTINE fld_get( sdjf, map )
[2528]	590	!!---------------------------------------------------------------------
[3294]	591	!! * ROUTINE fld_get *
[2528]	592	!!
	593	!! ** Purpose : read the data
	594	!!----------------------------------------------------------------------
[2715]	595	TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables
[3851]	596	TYPE(MAP_POINTER),INTENT(in) :: map ! global-to-local mapping indices
[2528]	597	!!
[2715]	598	INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk )
	599	INTEGER :: iw ! index into wgts array
[3680]	600	INTEGER :: ipdom ! index of the domain
[4245]	601	INTEGER :: idvar ! variable ID
	602	INTEGER :: idmspc ! number of spatial dimensions
	603	LOGICAL :: lmoor ! C1D case: point data
[2528]	604	!!---------------------------------------------------------------------
[3851]	605	!
[2528]	606	ipk = SIZE( sdjf%fnow, 3 )
[3680]	607	!
[3851]	608	IF( ASSOCIATED(map%ptr) ) THEN
[5132]	609	IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map )
	610	ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), map )
[3294]	611	ENDIF
	612	ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN
[2528]	613	CALL wgt_list( sdjf, iw )
[4230]	614	IF( sdjf%ln_tint ) THEN ; CALL fld_interp( sdjf%num, sdjf%clvar, iw , ipk , sdjf%fdta(:,:,:,2), &
	615	& sdjf%nrec_a(1), sdjf%lsmname )
	616	ELSE ; CALL fld_interp( sdjf%num, sdjf%clvar, iw , ipk , sdjf%fnow(:,:,: ), &
	617	& sdjf%nrec_a(1), sdjf%lsmname )
[2528]	618	ENDIF
	619	ELSE
[3680]	620	IF( SIZE(sdjf%fnow, 1) == jpi ) THEN ; ipdom = jpdom_data
	621	ELSE ; ipdom = jpdom_unknown
	622	ENDIF
[4245]	623	! C1D case: If product of spatial dimensions == ipk, then x,y are of
	624	! size 1 (point/mooring data): this must be read onto the central grid point
	625	idvar = iom_varid( sdjf%num, sdjf%clvar )
	626	idmspc = iom_file( sdjf%num )%ndims( idvar )
	627	IF( iom_file( sdjf%num )%luld( idvar ) ) idmspc = idmspc - 1
	628	lmoor = (idmspc == 0 .OR. PRODUCT( iom_file( sdjf%num )%dimsz( 1:MAX(idmspc,1) ,idvar ) ) == ipk)
	629	!
[2528]	630	SELECT CASE( ipk )
[3680]	631	CASE(1)
[4245]	632	IF( lk_c1d .AND. lmoor ) THEN
	633	IF( sdjf%ln_tint ) THEN
	634	CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fdta(2,2,1,2), sdjf%nrec_a(1) )
	635	CALL lbc_lnk( sdjf%fdta(:,:,1,2),'Z',1. )
	636	ELSE
	637	CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fnow(2,2,1 ), sdjf%nrec_a(1) )
	638	CALL lbc_lnk( sdjf%fnow(:,:,1 ),'Z',1. )
	639	ENDIF
	640	ELSE
	641	IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_a(1) )
	642	ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,1 ), sdjf%nrec_a(1) )
	643	ENDIF
[2528]	644	ENDIF
	645	CASE DEFAULT
[4245]	646	IF (lk_c1d .AND. lmoor ) THEN
	647	IF( sdjf%ln_tint ) THEN
	648	CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fdta(2,2,:,2), sdjf%nrec_a(1) )
	649	CALL lbc_lnk( sdjf%fdta(:,:,:,2),'Z',1. )
	650	ELSE
	651	CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fnow(2,2,: ), sdjf%nrec_a(1) )
	652	CALL lbc_lnk( sdjf%fnow(:,:,: ),'Z',1. )
	653	ENDIF
	654	ELSE
	655	IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) )
	656	ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1) )
	657	ENDIF
[2528]	658	ENDIF
	659	END SELECT
	660	ENDIF
	661	!
[3851]	662	sdjf%rotn(2) = .false. ! vector not yet rotated
[2528]	663
	664	END SUBROUTINE fld_get
	665
[3294]	666	SUBROUTINE fld_map( num, clvar, dta, nrec, map )
	667	!!---------------------------------------------------------------------
[3851]	668	!! * ROUTINE fld_map *
[3294]	669	!!
	670	!! ** Purpose : read global data from file and map onto local data
	671	!! using a general mapping (for open boundaries)
	672	!!----------------------------------------------------------------------
	673	#if defined key_bdy
[3651]	674	USE bdy_oce, ONLY: dta_global, dta_global2 ! workspace to read in global data arrays
[3294]	675	#endif
	676	INTEGER , INTENT(in ) :: num ! stream number
	677	CHARACTER(LEN=*) , INTENT(in ) :: clvar ! variable name
	678	REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional)
	679	INTEGER , INTENT(in ) :: nrec ! record number to read (ie time slice)
[5132]	680	TYPE(MAP_POINTER) , INTENT(in ) :: map ! global-to-local mapping indices
[3294]	681	!!
	682	INTEGER :: ipi ! length of boundary data on local process
	683	INTEGER :: ipj ! length of dummy dimension ( = 1 )
	684	INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk )
	685	INTEGER :: ilendta ! length of data in file
	686	INTEGER :: idvar ! variable ID
[3651]	687	INTEGER :: ib, ik, ji, jj ! loop counters
[3294]	688	INTEGER :: ierr
[3651]	689	REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read ! work space for global data
[3294]	690	!!---------------------------------------------------------------------
	691
	692	ipi = SIZE( dta, 1 )
	693	ipj = 1
	694	ipk = SIZE( dta, 3 )
	695
	696	idvar = iom_varid( num, clvar )
	697	ilendta = iom_file(num)%dimsz(1,idvar)
[3651]	698
	699	#if defined key_bdy
	700	ipj = iom_file(num)%dimsz(2,idvar)
[5132]	701	IF ( map%ll_unstruc) THEN ! unstructured open boundary data file
[3651]	702	dta_read => dta_global
[5132]	703	ELSE ! structured open boundary data file
[3651]	704	dta_read => dta_global2
	705	ENDIF
	706	#endif
	707
[3294]	708	IF(lwp) WRITE(numout,*) 'Dim size for ',TRIM(clvar),' is ', ilendta
	709	IF(lwp) WRITE(numout,*) 'Number of levels for ',TRIM(clvar),' is ', ipk
	710
	711	SELECT CASE( ipk )
[3851]	712	CASE(1) ; CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1 ), nrec )
	713	CASE DEFAULT ; CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec )
[3294]	714	END SELECT
	715	!
[5132]	716	IF ( map%ll_unstruc ) THEN ! unstructured open boundary data file
[3651]	717	DO ib = 1, ipi
	718	DO ik = 1, ipk
[5132]	719	dta(ib,1,ik) = dta_read(map%ptr(ib),1,ik)
[3651]	720	END DO
[3294]	721	END DO
[5132]	722	ELSE ! structured open boundary data file
[3651]	723	DO ib = 1, ipi
[5132]	724	jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta))
	725	ji=map%ptr(ib)-(jj-1)*ilendta
[3651]	726	DO ik = 1, ipk
	727	dta(ib,1,ik) = dta_read(ji,jj,ik)
	728	END DO
	729	END DO
	730	ENDIF
[3294]	731
	732	END SUBROUTINE fld_map
	733
	734
[2528]	735	SUBROUTINE fld_rot( kt, sd )
	736	!!---------------------------------------------------------------------
[3294]	737	!! * ROUTINE fld_rot *
[2528]	738	!!
	739	!! ** Purpose : Vector fields may need to be rotated onto the local grid direction
[2715]	740	!!----------------------------------------------------------------------
[2528]	741	INTEGER , INTENT(in ) :: kt ! ocean time step
	742	TYPE(FLD), INTENT(inout), DIMENSION(:) :: sd ! input field related variables
	743	!!
[3851]	744	INTEGER :: ju,jv,jk,jn ! loop indices
[3294]	745	INTEGER :: imf ! size of the structure sd
	746	INTEGER :: ill ! character length
	747	INTEGER :: iv ! indice of V component
	748	REAL(wp), POINTER, DIMENSION(:,:) :: utmp, vtmp ! temporary arrays for vector rotation
	749	CHARACTER (LEN=100) :: clcomp ! dummy weight name
[2528]	750	!!---------------------------------------------------------------------
[2715]	751
[3294]	752	CALL wrk_alloc( jpi,jpj, utmp, vtmp )
[2715]	753
[2528]	754	!! (sga: following code should be modified so that pairs arent searched for each time
	755	!
	756	imf = SIZE( sd )
	757	DO ju = 1, imf
	758	ill = LEN_TRIM( sd(ju)%vcomp )
[3851]	759	DO jn = 2-COUNT((/sd(ju)%ln_tint/)), 2
	760	IF( ill > 0 .AND. .NOT. sd(ju)%rotn(jn) ) THEN ! find vector rotations required
	761	IF( sd(ju)%vcomp(1:1) == 'U' ) THEN ! east-west component has symbolic name starting with 'U'
	762	! look for the north-south component which has same symbolic name but with 'U' replaced with 'V'
	763	clcomp = 'V' // sd(ju)%vcomp(2:ill) ! works even if ill == 1
	764	iv = -1
	765	DO jv = 1, imf
	766	IF( TRIM(sd(jv)%vcomp) == TRIM(clcomp) ) iv = jv
	767	END DO
	768	IF( iv > 0 ) THEN ! fields ju and iv are two components which need to be rotated together
	769	DO jk = 1, SIZE( sd(ju)%fnow, 3 )
	770	IF( sd(ju)%ln_tint )THEN
	771	CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->i', utmp(:,:) )
	772	CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->j', vtmp(:,:) )
	773	sd(ju)%fdta(:,:,jk,jn) = utmp(:,:) ; sd(iv)%fdta(:,:,jk,jn) = vtmp(:,:)
	774	ELSE
	775	CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->i', utmp(:,:) )
	776	CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->j', vtmp(:,:) )
	777	sd(ju)%fnow(:,:,jk ) = utmp(:,:) ; sd(iv)%fnow(:,:,jk ) = vtmp(:,:)
	778	ENDIF
	779	END DO
	780	sd(ju)%rotn(jn) = .TRUE. ! vector was rotated
	781	IF( lwp .AND. kt == nit000 ) WRITE(numout,*) &
	782	& 'fld_read: vector pair ('//TRIM(sd(ju)%clvar)//', '//TRIM(sd(iv)%clvar)//') rotated on to model grid'
	783	ENDIF
	784	ENDIF
	785	ENDIF
	786	END DO
[2528]	787	END DO
[2715]	788	!
[3294]	789	CALL wrk_dealloc( jpi,jpj, utmp, vtmp )
[2715]	790	!
[2528]	791	END SUBROUTINE fld_rot
	792
	793
[1628]	794	SUBROUTINE fld_clopn( sdjf, kyear, kmonth, kday, ldstop )
[1132]	795	!!---------------------------------------------------------------------
	796	!! * ROUTINE fld_clopn *
	797	!!
	798	!! ** Purpose : update the file name and open the file
	799	!!----------------------------------------------------------------------
[2715]	800	TYPE(FLD) , INTENT(inout) :: sdjf ! input field related variables
[3851]	801	INTEGER, OPTIONAL, INTENT(in ) :: kyear ! year value
	802	INTEGER, OPTIONAL, INTENT(in ) :: kmonth ! month value
	803	INTEGER, OPTIONAL, INTENT(in ) :: kday ! day value
[2715]	804	LOGICAL, OPTIONAL, INTENT(in ) :: ldstop ! stop if open to read a non-existing file (default = .TRUE.)
[3851]	805	!!
	806	LOGICAL :: llprevyr ! are we reading previous year file?
	807	LOGICAL :: llprevmth ! are we reading previous month file?
	808	INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd
	809	INTEGER :: isec_week ! number of seconds since start of the weekly file
	810	INTEGER :: indexyr ! year undex (O/1/2: previous/current/next)
	811	INTEGER :: iyear_len, imonth_len ! length (days) of iyear and imonth !
	812	CHARACTER(len = 256):: clname ! temporary file name
[2715]	813	!!----------------------------------------------------------------------
[3851]	814	IF( PRESENT(kyear) ) THEN ! use given values
	815	iyear = kyear
	816	imonth = kmonth
	817	iday = kday
[6487]	818	IF ( sdjf%cltype(1:4) == 'week' ) THEN ! find the day of the beginning of the week
	819	isec_week = ksec_week( sdjf%cltype(6:8) )- (86400 * 8 )
	820	llprevmth = isec_week > nsec_month ! longer time since beginning of the week than the month
	821	llprevyr = llprevmth .AND. nmonth == 1
	822	iyear = nyear - COUNT((/llprevyr /))
	823	imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /))
	824	iday = nday + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday)
	825	ENDIF
[3851]	826	ELSE ! use current day values
	827	IF ( sdjf%cltype(1:4) == 'week' ) THEN ! find the day of the beginning of the week
	828	isec_week = ksec_week( sdjf%cltype(6:8) ) ! second since the beginning of the week
	829	llprevmth = isec_week > nsec_month ! longer time since beginning of the week than the month
	830	llprevyr = llprevmth .AND. nmonth == 1
	831	ELSE
	832	isec_week = 0
	833	llprevmth = .FALSE.
	834	llprevyr = .FALSE.
	835	ENDIF
	836	iyear = nyear - COUNT((/llprevyr /))
	837	imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /))
	838	iday = nday + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday)
	839	ENDIF
[1132]	840
	841	! build the new filename if not climatological data
[3851]	842	clname=TRIM(sdjf%clrootname)
[2528]	843	!
[3851]	844	! note that sdjf%ln_clim is is only acting on the presence of the year in the file name
[2528]	845	IF( .NOT. sdjf%ln_clim ) THEN
[3851]	846	WRITE(clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), iyear ! add year
	847	IF( sdjf%cltype /= 'yearly' ) WRITE(clname, '(a,"m" ,i2.2)' ) TRIM( clname ), imonth ! add month
[2528]	848	ELSE
	849	! build the new filename if climatological data
[3851]	850	IF( sdjf%cltype /= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), imonth ! add month
[888]	851	ENDIF
[2528]	852	IF( sdjf%cltype == 'daily' .OR. sdjf%cltype(1:4) == 'week' ) &
[3851]	853	& WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), iday ! add day
[2528]	854	!
[3851]	855	IF( TRIM(clname) /= TRIM(sdjf%clname) .OR. sdjf%num == 0 ) THEN ! new file to be open
	856
	857	sdjf%clname = TRIM(clname)
	858	IF( sdjf%num /= 0 ) CALL iom_close( sdjf%num ) ! close file if already open
	859	CALL iom_open( sdjf%clname, sdjf%num, ldstop = ldstop, ldiof = LEN(TRIM(sdjf%wgtname)) > 0 )
	860
	861	! find the last record to be read -> update sdjf%nreclast
	862	indexyr = iyear - nyear + 1
	863	iyear_len = nyear_len( indexyr )
	864	SELECT CASE ( indexyr )
	865	CASE ( 0 ) ; imonth_len = 31 ! previous year -> imonth = 12
	866	CASE ( 1 ) ; imonth_len = nmonth_len(imonth)
	867	CASE ( 2 ) ; imonth_len = 31 ! next year -> imonth = 1
	868	END SELECT
	869
	870	! last record to be read in the current file
	871	IF ( sdjf%nfreqh == -12 ) THEN ; sdjf%nreclast = 1 ! yearly mean
	872	ELSEIF( sdjf%nfreqh == -1 ) THEN ! monthly mean
	873	IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = 1
	874	ELSE ; sdjf%nreclast = 12
	875	ENDIF
	876	ELSE ! higher frequency mean (in hours)
[4245]	877	IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = NINT( 24 * imonth_len / sdjf%nfreqh )
	878	ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; sdjf%nreclast = NINT( 24 * 7 / sdjf%nfreqh )
	879	ELSEIF( sdjf%cltype == 'daily' ) THEN ; sdjf%nreclast = NINT( 24 / sdjf%nfreqh )
	880	ELSE ; sdjf%nreclast = NINT( 24 * iyear_len / sdjf%nfreqh )
[3851]	881	ENDIF
	882	ENDIF
	883
	884	ENDIF
	885	!
[1132]	886	END SUBROUTINE fld_clopn
	887
	888
	889	SUBROUTINE fld_fill( sdf, sdf_n, cdir, cdcaller, cdtitle, cdnam )
	890	!!---------------------------------------------------------------------
	891	!! * ROUTINE fld_fill *
	892	!!
	893	!! ** Purpose : fill sdf with sdf_n and control print
	894	!!----------------------------------------------------------------------
	895	TYPE(FLD) , DIMENSION(:), INTENT(inout) :: sdf ! structure of input fields (file informations, fields read)
	896	TYPE(FLD_N), DIMENSION(:), INTENT(in ) :: sdf_n ! array of namelist information structures
	897	CHARACTER(len=*) , INTENT(in ) :: cdir ! Root directory for location of flx files
	898	CHARACTER(len=*) , INTENT(in ) :: cdcaller !
	899	CHARACTER(len=*) , INTENT(in ) :: cdtitle !
	900	CHARACTER(len=*) , INTENT(in ) :: cdnam !
[888]	901	!
[1132]	902	INTEGER :: jf ! dummy indices
	903	!!---------------------------------------------------------------------
[888]	904
[1132]	905	DO jf = 1, SIZE(sdf)
	906	sdf(jf)%clrootname = TRIM( cdir )//TRIM( sdf_n(jf)%clname )
[3851]	907	sdf(jf)%clname = "not yet defined"
[1730]	908	sdf(jf)%nfreqh = sdf_n(jf)%nfreqh
[1132]	909	sdf(jf)%clvar = sdf_n(jf)%clvar
	910	sdf(jf)%ln_tint = sdf_n(jf)%ln_tint
	911	sdf(jf)%ln_clim = sdf_n(jf)%ln_clim
[2528]	912	sdf(jf)%cltype = sdf_n(jf)%cltype
[3851]	913	sdf(jf)%num = -1
	914	sdf(jf)%wgtname = " "
[1730]	915	IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )//TRIM( sdf_n(jf)%wname )
[4230]	916	sdf(jf)%lsmname = " "
	917	IF( LEN( TRIM(sdf_n(jf)%lname) ) > 0 ) sdf(jf)%lsmname = TRIM( cdir )//TRIM( sdf_n(jf)%lname )
[3851]	918	sdf(jf)%vcomp = sdf_n(jf)%vcomp
	919	sdf(jf)%rotn(:) = .TRUE. ! pretend to be rotated -> won't try to rotate data before the first call to fld_get
	920	IF( sdf(jf)%cltype(1:4) == 'week' .AND. nn_leapy == 0 ) &
	921	& CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs nn_leapy = 1')
	922	IF( sdf(jf)%cltype(1:4) == 'week' .AND. sdf(jf)%ln_clim ) &
	923	& CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs ln_clim = .FALSE.')
[4371]	924	sdf(jf)%nreclast = -1 ! Set to non zero default value to avoid errors, is updated to meaningful value during fld_clopn
[1132]	925	END DO
	926
	927	IF(lwp) THEN ! control print
	928	WRITE(numout,*)
	929	WRITE(numout,*) TRIM( cdcaller )//' : '//TRIM( cdtitle )
	930	WRITE(numout,*) (/ ('~', jf = 1, LEN_TRIM( cdcaller ) ) /)
	931	WRITE(numout,*) ' '//TRIM( cdnam )//' Namelist'
	932	WRITE(numout,*) ' list of files and frequency (>0: in hours ; <0 in months)'
	933	DO jf = 1, SIZE(sdf)
	934	WRITE(numout,*) ' root filename: ' , TRIM( sdf(jf)%clrootname ), &
	935	& ' variable name: ' , TRIM( sdf(jf)%clvar )
[1730]	936	WRITE(numout,*) ' frequency: ' , sdf(jf)%nfreqh , &
[1132]	937	& ' time interp: ' , sdf(jf)%ln_tint , &
	938	& ' climatology: ' , sdf(jf)%ln_clim , &
[1275]	939	& ' weights : ' , TRIM( sdf(jf)%wgtname ), &
	940	& ' pairing : ' , TRIM( sdf(jf)%vcomp ), &
[4230]	941	& ' data type: ' , sdf(jf)%cltype , &
	942	& ' land/sea mask:' , TRIM( sdf(jf)%lsmname )
[10770]	943	IF(nprint > 4) call flush(numout)
[1132]	944	END DO
	945	ENDIF
	946
	947	END SUBROUTINE fld_fill
	948
	949
[1275]	950	SUBROUTINE wgt_list( sd, kwgt )
	951	!!---------------------------------------------------------------------
	952	!! * ROUTINE wgt_list *
	953	!!
	954	!! ** Purpose : search array of WGTs and find a weights file
	955	!! entry, or return a new one adding it to the end
	956	!! if it is a new entry, the weights data is read in and
	957	!! restructured (fld_weight)
	958	!!----------------------------------------------------------------------
[2715]	959	TYPE( FLD ), INTENT(in ) :: sd ! field with name of weights file
	960	INTEGER , INTENT(inout) :: kwgt ! index of weights
[1275]	961	!!
[2715]	962	INTEGER :: kw, nestid ! local integer
	963	LOGICAL :: found ! local logical
[1275]	964	!!----------------------------------------------------------------------
	965	!
	966	!! search down linked list
	967	!! weights filename is either present or we hit the end of the list
	968	found = .FALSE.
	969
	970	!! because agrif nest part of filenames are now added in iom_open
	971	!! to distinguish between weights files on the different grids, need to track
	972	!! nest number explicitly
	973	nestid = 0
	974	#if defined key_agrif
	975	nestid = Agrif_Fixed()
	976	#endif
	977	DO kw = 1, nxt_wgt-1
	978	IF( TRIM(ref_wgts(kw)%wgtname) == TRIM(sd%wgtname) .AND. &
	979	ref_wgts(kw)%nestid == nestid) THEN
	980	kwgt = kw
	981	found = .TRUE.
	982	EXIT
	983	ENDIF
	984	END DO
	985	IF( .NOT.found ) THEN
	986	kwgt = nxt_wgt
	987	CALL fld_weight( sd )
	988	ENDIF
[2715]	989	!
[1275]	990	END SUBROUTINE wgt_list
	991
[2715]	992
[1275]	993	SUBROUTINE wgt_print( )
	994	!!---------------------------------------------------------------------
	995	!! * ROUTINE wgt_print *
	996	!!
	997	!! ** Purpose : print the list of known weights
	998	!!----------------------------------------------------------------------
[2715]	999	INTEGER :: kw !
[1275]	1000	!!----------------------------------------------------------------------
	1001	!
[10770]	1002	IF(nprint > 0) THEN
	1003	DO kw = 1, nxt_wgt-1
	1004	WRITE(numout,*) 'weight file: ',TRIM(ref_wgts(kw)%wgtname)
	1005	IF(nprint > 1) THEN
	1006	WRITE(numout,*) ' ddims: ',ref_wgts(kw)%ddims(1),ref_wgts(kw)%ddims(2)
	1007	WRITE(numout,*) ' numwgt: ',ref_wgts(kw)%numwgt
	1008	WRITE(numout,*) ' jpiwgt: ',ref_wgts(kw)%jpiwgt
	1009	WRITE(numout,*) ' jpjwgt: ',ref_wgts(kw)%jpjwgt
	1010	WRITE(numout,*) ' botleft: ',ref_wgts(kw)%botleft
	1011	WRITE(numout,*) ' topright: ',ref_wgts(kw)%topright
	1012	IF( ref_wgts(kw)%cyclic ) THEN
	1013	WRITE(numout,*) ' cyclical'
	1014	IF( ref_wgts(kw)%overlap > 0 ) WRITE(numout,*) ' with overlap of ', ref_wgts(kw)%overlap
	1015	ELSE
	1016	WRITE(numout,*) ' not cyclical'
	1017	ENDIF
	1018	IF( ASSOCIATED(ref_wgts(kw)%data_wgt) ) WRITE(numout,*) ' allocated'
	1019	ENDIF
	1020	END DO
	1021	ENDIF
[2715]	1022	!
[1275]	1023	END SUBROUTINE wgt_print
	1024
[2715]	1025
[1275]	1026	SUBROUTINE fld_weight( sd )
	1027	!!---------------------------------------------------------------------
	1028	!! * ROUTINE fld_weight *
	1029	!!
	1030	!! ** Purpose : create a new WGT structure and fill in data from
	1031	!! file, restructuring as required
	1032	!!----------------------------------------------------------------------
[2715]	1033	TYPE( FLD ), INTENT(in) :: sd ! field with name of weights file
	1034	!!
[3294]	1035	INTEGER :: jn ! dummy loop indices
	1036	INTEGER :: inum ! temporary logical unit
	1037	INTEGER :: id ! temporary variable id
	1038	INTEGER :: ipk ! temporary vertical dimension
	1039	CHARACTER (len=5) :: aname
[5399]	1040	INTEGER , DIMENSION(:), ALLOCATABLE :: ddims
[3294]	1041	INTEGER , POINTER, DIMENSION(:,:) :: data_src
	1042	REAL(wp), POINTER, DIMENSION(:,:) :: data_tmp
	1043	LOGICAL :: cyclical
	1044	INTEGER :: zwrap ! local integer
[1275]	1045	!!----------------------------------------------------------------------
	1046	!
[3294]	1047	CALL wrk_alloc( jpi,jpj, data_src ) ! integer
	1048	CALL wrk_alloc( jpi,jpj, data_tmp )
[2715]	1049	!
[1275]	1050	IF( nxt_wgt > tot_wgts ) THEN
[2777]	1051	CALL ctl_stop("fld_weight: weights array size exceeded, increase tot_wgts")
[1275]	1052	ENDIF
	1053	!
	1054	!! new weights file entry, add in extra information
	1055	!! a weights file represents a 2D grid of a certain shape, so we assume that the current
	1056	!! input data file is representative of all other files to be opened and processed with the
	1057	!! current weights file
	1058
	1059	!! open input data file (non-model grid)
[1319]	1060	CALL iom_open( sd%clname, inum, ldiof = LEN(TRIM(sd%wgtname)) > 0 )
[1275]	1061
	1062	!! get dimensions
[5399]	1063	IF ( SIZE(sd%fnow, 3) > 1 ) THEN
	1064	ALLOCATE( ddims(4) )
	1065	ELSE
	1066	ALLOCATE( ddims(3) )
	1067	ENDIF
[1275]	1068	id = iom_varid( inum, sd%clvar, ddims )
	1069
[2528]	1070	!! close it
	1071	CALL iom_close( inum )
[1275]	1072
[2528]	1073	!! now open the weights file
[1275]	1074
[2528]	1075	CALL iom_open ( sd%wgtname, inum ) ! interpolation weights
	1076	IF ( inum > 0 ) THEN
[1275]	1077
[2528]	1078	!! determine whether we have an east-west cyclic grid
	1079	!! from global attribute called "ew_wrap" in the weights file
	1080	!! note that if not found, iom_getatt returns -999 and cyclic with no overlap is assumed
	1081	!! since this is the most common forcing configuration
[1275]	1082
[2528]	1083	CALL iom_getatt(inum, 'ew_wrap', zwrap)
	1084	IF( zwrap >= 0 ) THEN
[1275]	1085	cyclical = .TRUE.
[2528]	1086	ELSE IF( zwrap == -999 ) THEN
[1275]	1087	cyclical = .TRUE.
[2528]	1088	zwrap = 0
	1089	ELSE
	1090	cyclical = .FALSE.
[1275]	1091	ENDIF
	1092
	1093	ref_wgts(nxt_wgt)%ddims(1) = ddims(1)
	1094	ref_wgts(nxt_wgt)%ddims(2) = ddims(2)
	1095	ref_wgts(nxt_wgt)%wgtname = sd%wgtname
[2528]	1096	ref_wgts(nxt_wgt)%overlap = zwrap
	1097	ref_wgts(nxt_wgt)%cyclic = cyclical
[1275]	1098	ref_wgts(nxt_wgt)%nestid = 0
	1099	#if defined key_agrif
	1100	ref_wgts(nxt_wgt)%nestid = Agrif_Fixed()
	1101	#endif
	1102	!! weights file is stored as a set of weights (wgt01->wgt04 or wgt01->wgt16)
	1103	!! for each weight wgtNN there is an integer array srcNN which gives the point in
	1104	!! the input data grid which is to be multiplied by the weight
	1105	!! they are both arrays on the model grid so the result of the multiplication is
	1106	!! added into an output array on the model grid as a running sum
	1107
	1108	!! two possible cases: bilinear (4 weights) or bicubic (16 weights)
	1109	id = iom_varid(inum, 'src05', ldstop=.FALSE.)
	1110	IF( id <= 0) THEN
	1111	ref_wgts(nxt_wgt)%numwgt = 4
	1112	ELSE
	1113	ref_wgts(nxt_wgt)%numwgt = 16
	1114	ENDIF
	1115
	1116	ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(jpi,jpj,4) )
	1117	ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(jpi,jpj,4) )
	1118	ALLOCATE( ref_wgts(nxt_wgt)%data_wgt(jpi,jpj,ref_wgts(nxt_wgt)%numwgt) )
	1119
	1120	DO jn = 1,4
	1121	aname = ' '
	1122	WRITE(aname,'(a3,i2.2)') 'src',jn
	1123	data_tmp(:,:) = 0
[1955]	1124	CALL iom_get ( inum, jpdom_data, aname, data_tmp(:,:) )
[1275]	1125	data_src(:,:) = INT(data_tmp(:,:))
	1126	ref_wgts(nxt_wgt)%data_jpj(:,:,jn) = 1 + (data_src(:,:)-1) / ref_wgts(nxt_wgt)%ddims(1)
	1127	ref_wgts(nxt_wgt)%data_jpi(:,:,jn) = data_src(:,:) - ref_wgts(nxt_wgt)%ddims(1)*(ref_wgts(nxt_wgt)%data_jpj(:,:,jn)-1)
	1128	END DO
	1129
	1130	DO jn = 1, ref_wgts(nxt_wgt)%numwgt
	1131	aname = ' '
	1132	WRITE(aname,'(a3,i2.2)') 'wgt',jn
[1955]	1133	ref_wgts(nxt_wgt)%data_wgt(:,:,jn) = 0.0
	1134	CALL iom_get ( inum, jpdom_data, aname, ref_wgts(nxt_wgt)%data_wgt(:,:,jn) )
[1275]	1135	END DO
	1136	CALL iom_close (inum)
	1137
	1138	! find min and max indices in grid
[1955]	1139	ref_wgts(nxt_wgt)%botleft(1) = MINVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:))
	1140	ref_wgts(nxt_wgt)%botleft(2) = MINVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:))
	1141	ref_wgts(nxt_wgt)%topright(1) = MAXVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:))
	1142	ref_wgts(nxt_wgt)%topright(2) = MAXVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:))
[1275]	1143
	1144	! and therefore dimensions of the input box
	1145	ref_wgts(nxt_wgt)%jpiwgt = ref_wgts(nxt_wgt)%topright(1) - ref_wgts(nxt_wgt)%botleft(1) + 1
	1146	ref_wgts(nxt_wgt)%jpjwgt = ref_wgts(nxt_wgt)%topright(2) - ref_wgts(nxt_wgt)%botleft(2) + 1
	1147
	1148	! shift indexing of source grid
	1149	ref_wgts(nxt_wgt)%data_jpi(:,:,:) = ref_wgts(nxt_wgt)%data_jpi(:,:,:) - ref_wgts(nxt_wgt)%botleft(1) + 1
	1150	ref_wgts(nxt_wgt)%data_jpj(:,:,:) = ref_wgts(nxt_wgt)%data_jpj(:,:,:) - ref_wgts(nxt_wgt)%botleft(2) + 1
	1151
	1152	! create input grid, give it a halo to allow gradient calculations
[1702]	1153	! SA: +3 stencil is a patch to avoid out-of-bound computation in some configuration.
	1154	! a more robust solution will be given in next release
[2528]	1155	ipk = SIZE(sd%fnow, 3)
	1156	ALLOCATE( ref_wgts(nxt_wgt)%fly_dta(ref_wgts(nxt_wgt)%jpiwgt+3, ref_wgts(nxt_wgt)%jpjwgt+3 ,ipk) )
	1157	IF( ref_wgts(nxt_wgt)%cyclic ) ALLOCATE( ref_wgts(nxt_wgt)%col(1,ref_wgts(nxt_wgt)%jpjwgt+3,ipk) )
[1275]	1158
	1159	nxt_wgt = nxt_wgt + 1
	1160
	1161	ELSE
	1162	CALL ctl_stop( ' fld_weight : unable to read the file ' )
	1163	ENDIF
	1164
[5399]	1165	DEALLOCATE (ddims )
	1166
[3294]	1167	CALL wrk_dealloc( jpi,jpj, data_src ) ! integer
	1168	CALL wrk_dealloc( jpi,jpj, data_tmp )
[2715]	1169	!
[1275]	1170	END SUBROUTINE fld_weight
	1171
[2715]	1172
[4230]	1173	SUBROUTINE apply_seaoverland(clmaskfile,zfieldo,jpi1_lsm,jpi2_lsm,jpj1_lsm, &
	1174	& jpj2_lsm,itmpi,itmpj,itmpz,rec1_lsm,recn_lsm)
[1275]	1175	!!---------------------------------------------------------------------
[4230]	1176	!! * ROUTINE apply_seaoverland *
	1177	!!
	1178	!! ** Purpose : avoid spurious fluxes in coastal or near-coastal areas
	1179	!! due to the wrong usage of "land" values from the coarse
	1180	!! atmospheric model when spatial interpolation is required
	1181	!! D. Delrosso INGV
	1182	!!----------------------------------------------------------------------
	1183	INTEGER :: inum,jni,jnj,jnz,jc ! temporary indices
	1184	INTEGER, INTENT(in) :: itmpi,itmpj,itmpz ! lengths
	1185	INTEGER, INTENT(in) :: jpi1_lsm,jpi2_lsm,jpj1_lsm,jpj2_lsm ! temporary indices
	1186	INTEGER, DIMENSION(3), INTENT(in) :: rec1_lsm,recn_lsm ! temporary arrays for start and length
	1187	REAL(wp),DIMENSION (:,:,:),INTENT(inout) :: zfieldo ! input/output array for seaoverland application
	1188	REAL(wp),DIMENSION (:,:,:),ALLOCATABLE :: zslmec1 ! temporary array for land point detection
	1189	REAL(wp),DIMENSION (:,:), ALLOCATABLE :: zfieldn ! array of forcing field with undeff for land points
	1190	REAL(wp),DIMENSION (:,:), ALLOCATABLE :: zfield ! array of forcing field
	1191	CHARACTER (len=100), INTENT(in) :: clmaskfile ! land/sea mask file name
	1192	!!---------------------------------------------------------------------
	1193	ALLOCATE ( zslmec1(itmpi,itmpj,itmpz) )
	1194	ALLOCATE ( zfieldn(itmpi,itmpj) )
	1195	ALLOCATE ( zfield(itmpi,itmpj) )
	1196
	1197	! Retrieve the land sea mask data
	1198	CALL iom_open( clmaskfile, inum )
	1199	SELECT CASE( SIZE(zfieldo(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),3) )
	1200	CASE(1)
	1201	CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), 1, rec1_lsm, recn_lsm)
	1202	CASE DEFAULT
	1203	CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), 1, rec1_lsm, recn_lsm)
	1204	END SELECT
	1205	CALL iom_close( inum )
	1206
	1207	DO jnz=1,rec1_lsm(3) !! Loop over k dimension
	1208
	1209	DO jni=1,itmpi !! copy the original field into a tmp array
	1210	DO jnj=1,itmpj !! substituting undeff over land points
	1211	zfieldn(jni,jnj) = zfieldo(jni,jnj,jnz)
	1212	IF ( zslmec1(jni,jnj,jnz) == 1. ) THEN
	1213	zfieldn(jni,jnj) = undeff_lsm
	1214	ENDIF
	1215	END DO
	1216	END DO
	1217
	1218	CALL seaoverland(zfieldn,itmpi,itmpj,zfield)
	1219	DO jc=1,nn_lsm
	1220	CALL seaoverland(zfield,itmpi,itmpj,zfield)
	1221	END DO
	1222
	1223	! Check for Undeff and substitute original values
	1224	IF(ANY(zfield==undeff_lsm)) THEN
	1225	DO jni=1,itmpi
	1226	DO jnj=1,itmpj
	1227	IF (zfield(jni,jnj)==undeff_lsm) THEN
	1228	zfield(jni,jnj) = zfieldo(jni,jnj,jnz)
	1229	ENDIF
	1230	ENDDO
	1231	ENDDO
	1232	ENDIF
	1233
	1234	zfieldo(:,:,jnz)=zfield(:,:)
	1235
	1236	END DO !! End Loop over k dimension
	1237
	1238	DEALLOCATE ( zslmec1 )
	1239	DEALLOCATE ( zfieldn )
	1240	DEALLOCATE ( zfield )
	1241
	1242	END SUBROUTINE apply_seaoverland
	1243
	1244
	1245	SUBROUTINE seaoverland(zfieldn,ileni,ilenj,zfield)
	1246	!!---------------------------------------------------------------------
	1247	!! * ROUTINE seaoverland *
	1248	!!
	1249	!! ** Purpose : create shifted matrices for seaoverland application
	1250	!! D. Delrosso INGV
	1251	!!----------------------------------------------------------------------
	1252	INTEGER,INTENT(in) :: ileni,ilenj ! lengths
	1253	REAL,DIMENSION (ileni,ilenj),INTENT(in) :: zfieldn ! array of forcing field with undeff for land points
	1254	REAL,DIMENSION (ileni,ilenj),INTENT(out) :: zfield ! array of forcing field
	1255	REAL,DIMENSION (ileni,ilenj) :: zmat1,zmat2,zmat3,zmat4 ! temporary arrays for seaoverland application
	1256	REAL,DIMENSION (ileni,ilenj) :: zmat5,zmat6,zmat7,zmat8 ! temporary arrays for seaoverland application
	1257	REAL,DIMENSION (ileni,ilenj) :: zlsm2d ! temporary arrays for seaoverland application
	1258	REAL,DIMENSION (ileni,ilenj,8) :: zlsm3d ! temporary arrays for seaoverland application
	1259	LOGICAL,DIMENSION (ileni,ilenj,8) :: ll_msknan3d ! logical mask for undeff detection
	1260	LOGICAL,DIMENSION (ileni,ilenj) :: ll_msknan2d ! logical mask for undeff detection
	1261	!!----------------------------------------------------------------------
	1262	zmat8 = eoshift(zfieldn , SHIFT=-1, BOUNDARY = (/zfieldn(:,1)/) ,DIM=2)
	1263	zmat1 = eoshift(zmat8 , SHIFT=-1, BOUNDARY = (/zmat8(1,:)/) ,DIM=1)
	1264	zmat2 = eoshift(zfieldn , SHIFT=-1, BOUNDARY = (/zfieldn(1,:)/) ,DIM=1)
	1265	zmat4 = eoshift(zfieldn , SHIFT= 1, BOUNDARY = (/zfieldn(:,ilenj)/),DIM=2)
	1266	zmat3 = eoshift(zmat4 , SHIFT=-1, BOUNDARY = (/zmat4(1,:)/) ,DIM=1)
	1267	zmat5 = eoshift(zmat4 , SHIFT= 1, BOUNDARY = (/zmat4(ileni,:)/) ,DIM=1)
	1268	zmat6 = eoshift(zfieldn , SHIFT= 1, BOUNDARY = (/zfieldn(ileni,:)/),DIM=1)
	1269	zmat7 = eoshift(zmat8 , SHIFT= 1, BOUNDARY = (/zmat8(ileni,:)/) ,DIM=1)
	1270
	1271	zlsm3d = RESHAPE( (/ zmat1, zmat2, zmat3, zmat4, zmat5, zmat6, zmat7, zmat8 /), (/ ileni, ilenj, 8 /))
	1272	ll_msknan3d = .not.(zlsm3d==undeff_lsm)
	1273	ll_msknan2d = .not.(zfieldn==undeff_lsm) ! FALSE where is Undeff (land)
	1274	zlsm2d = (SUM ( zlsm3d, 3 , ll_msknan3d ) )/(MAX(1,(COUNT( ll_msknan3d , 3 )) ))
	1275	WHERE ((COUNT( ll_msknan3d , 3 )) == 0.0_wp) zlsm2d = undeff_lsm
	1276	zfield = MERGE (zfieldn,zlsm2d,ll_msknan2d)
	1277	END SUBROUTINE seaoverland
	1278
	1279
	1280	SUBROUTINE fld_interp( num, clvar, kw, kk, dta, &
	1281	& nrec, lsmfile)
	1282	!!---------------------------------------------------------------------
[1275]	1283	!! * ROUTINE fld_interp *
	1284	!!
	1285	!! ** Purpose : apply weights to input gridded data to create data
	1286	!! on model grid
	1287	!!----------------------------------------------------------------------
[2715]	1288	INTEGER , INTENT(in ) :: num ! stream number
	1289	CHARACTER(LEN=*) , INTENT(in ) :: clvar ! variable name
	1290	INTEGER , INTENT(in ) :: kw ! weights number
	1291	INTEGER , INTENT(in ) :: kk ! vertical dimension of kk
	1292	REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: dta ! output field on model grid
	1293	INTEGER , INTENT(in ) :: nrec ! record number to read (ie time slice)
[4230]	1294	CHARACTER(LEN=*) , INTENT(in ) :: lsmfile ! land sea mask file name
[1275]	1295	!!
[6487]	1296	REAL(wp),DIMENSION(:,:,:),ALLOCATABLE :: ztmp_fly_dta ! temporary array of values on input grid
[4230]	1297	INTEGER, DIMENSION(3) :: rec1,recn ! temporary arrays for start and length
	1298	INTEGER, DIMENSION(3) :: rec1_lsm,recn_lsm ! temporary arrays for start and length in case of seaoverland
	1299	INTEGER :: ii_lsm1,ii_lsm2,ij_lsm1,ij_lsm2 ! temporary indices
	1300	INTEGER :: jk, jn, jm, jir, jjr ! loop counters
	1301	INTEGER :: ni, nj ! lengths
	1302	INTEGER :: jpimin,jpiwid ! temporary indices
	1303	INTEGER :: jpimin_lsm,jpiwid_lsm ! temporary indices
	1304	INTEGER :: jpjmin,jpjwid ! temporary indices
	1305	INTEGER :: jpjmin_lsm,jpjwid_lsm ! temporary indices
	1306	INTEGER :: jpi1,jpi2,jpj1,jpj2 ! temporary indices
	1307	INTEGER :: jpi1_lsm,jpi2_lsm,jpj1_lsm,jpj2_lsm ! temporary indices
	1308	INTEGER :: itmpi,itmpj,itmpz ! lengths
	1309
[1275]	1310	!!----------------------------------------------------------------------
	1311	!
	1312	!! for weighted interpolation we have weights at four corners of a box surrounding
	1313	!! a model grid point, each weight is multiplied by a grid value (bilinear case)
	1314	!! or by a grid value and gradients at the corner point (bicubic case)
	1315	!! so we need to have a 4 by 4 subgrid surrounding each model point to cover both cases
	1316
[2528]	1317	!! sub grid from non-model input grid which encloses all grid points in this nemo process
[1275]	1318	jpimin = ref_wgts(kw)%botleft(1)
	1319	jpjmin = ref_wgts(kw)%botleft(2)
	1320	jpiwid = ref_wgts(kw)%jpiwgt
	1321	jpjwid = ref_wgts(kw)%jpjwgt
	1322
[2528]	1323	!! when reading in, expand this sub-grid by one halo point all the way round for calculating gradients
[1275]	1324	rec1(1) = MAX( jpimin-1, 1 )
	1325	rec1(2) = MAX( jpjmin-1, 1 )
[2528]	1326	rec1(3) = 1
[1275]	1327	recn(1) = MIN( jpiwid+2, ref_wgts(kw)%ddims(1)-rec1(1)+1 )
	1328	recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 )
[2528]	1329	recn(3) = kk
[1275]	1330
[2528]	1331	!! where we need to put it in the non-nemo grid fly_dta
	1332	!! note that jpi1 and jpj1 only differ from 1 when jpimin and jpjmin are 1
	1333	!! (ie at the extreme west or south of the whole input grid) and similarly for jpi2 and jpj2
[1275]	1334	jpi1 = 2 + rec1(1) - jpimin
	1335	jpj1 = 2 + rec1(2) - jpjmin
	1336	jpi2 = jpi1 + recn(1) - 1
	1337	jpj2 = jpj1 + recn(2) - 1
	1338
	1339
[4230]	1340	IF( LEN( TRIM(lsmfile) ) > 0 ) THEN
	1341	!! indeces for ztmp_fly_dta
	1342	! --------------------------
	1343	rec1_lsm(1)=MAX(rec1(1)-nn_lsm,1) ! starting index for enlarged external data, x direction
	1344	rec1_lsm(2)=MAX(rec1(2)-nn_lsm,1) ! starting index for enlarged external data, y direction
	1345	rec1_lsm(3) = 1 ! vertical dimension
	1346	recn_lsm(1)=MIN(rec1(1)-rec1_lsm(1)+recn(1)+nn_lsm,ref_wgts(kw)%ddims(1)-rec1_lsm(1)) ! n points in x direction
	1347	recn_lsm(2)=MIN(rec1(2)-rec1_lsm(2)+recn(2)+nn_lsm,ref_wgts(kw)%ddims(2)-rec1_lsm(2)) ! n points in y direction
	1348	recn_lsm(3) = kk ! number of vertical levels in the input file
	1349
	1350	! Avoid out of bound
	1351	jpimin_lsm = MAX( rec1_lsm(1)+1, 1 )
	1352	jpjmin_lsm = MAX( rec1_lsm(2)+1, 1 )
	1353	jpiwid_lsm = MIN( recn_lsm(1)-2,ref_wgts(kw)%ddims(1)-rec1(1)+1)
	1354	jpjwid_lsm = MIN( recn_lsm(2)-2,ref_wgts(kw)%ddims(2)-rec1(2)+1)
	1355
	1356	jpi1_lsm = 2+rec1_lsm(1)-jpimin_lsm
	1357	jpj1_lsm = 2+rec1_lsm(2)-jpjmin_lsm
	1358	jpi2_lsm = jpi1_lsm + recn_lsm(1) - 1
	1359	jpj2_lsm = jpj1_lsm + recn_lsm(2) - 1
	1360
	1361
[6487]	1362	itmpi=jpi2_lsm-jpi1_lsm+1
	1363	itmpj=jpj2_lsm-jpj1_lsm+1
[4230]	1364	itmpz=kk
	1365	ALLOCATE(ztmp_fly_dta(itmpi,itmpj,itmpz))
	1366	ztmp_fly_dta(:,:,:) = 0.0
	1367	SELECT CASE( SIZE(ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),3) )
	1368	CASE(1)
	1369	CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), &
	1370	& nrec, rec1_lsm, recn_lsm)
	1371	CASE DEFAULT
	1372	CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), &
	1373	& nrec, rec1_lsm, recn_lsm)
	1374	END SELECT
	1375	CALL apply_seaoverland(lsmfile,ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), &
	1376	& jpi1_lsm,jpi2_lsm,jpj1_lsm,jpj2_lsm, &
	1377	& itmpi,itmpj,itmpz,rec1_lsm,recn_lsm)
	1378
	1379
	1380	! Relative indeces for remapping
	1381	ii_lsm1 = (rec1(1)-rec1_lsm(1))+1
	1382	ii_lsm2 = (ii_lsm1+recn(1))-1
	1383	ij_lsm1 = (rec1(2)-rec1_lsm(2))+1
	1384	ij_lsm2 = (ij_lsm1+recn(2))-1
	1385
	1386	ref_wgts(kw)%fly_dta(:,:,:) = 0.0
	1387	ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:) = ztmp_fly_dta(ii_lsm1:ii_lsm2,ij_lsm1:ij_lsm2,:)
	1388	DEALLOCATE(ztmp_fly_dta)
	1389
	1390	ELSE
	1391
	1392	ref_wgts(kw)%fly_dta(:,:,:) = 0.0
	1393	SELECT CASE( SIZE(ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:),3) )
	1394	CASE(1)
	1395	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,1), nrec, rec1, recn)
	1396	CASE DEFAULT
	1397	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:), nrec, rec1, recn)
	1398	END SELECT
	1399	ENDIF
	1400
	1401
[1275]	1402	!! first four weights common to both bilinear and bicubic
[2528]	1403	!! data_jpi, data_jpj have already been shifted to (1,1) corresponding to botleft
[1275]	1404	!! note that we have to offset by 1 into fly_dta array because of halo
[2528]	1405	dta(:,:,:) = 0.0
[1275]	1406	DO jk = 1,4
	1407	DO jn = 1, jpj
	1408	DO jm = 1,jpi
	1409	ni = ref_wgts(kw)%data_jpi(jm,jn,jk)
	1410	nj = ref_wgts(kw)%data_jpj(jm,jn,jk)
[2528]	1411	dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk) * ref_wgts(kw)%fly_dta(ni+1,nj+1,:)
[1275]	1412	END DO
	1413	END DO
	1414	END DO
	1415
	1416	IF (ref_wgts(kw)%numwgt .EQ. 16) THEN
	1417
	1418	!! fix up halo points that we couldnt read from file
	1419	IF( jpi1 == 2 ) THEN
[2528]	1420	ref_wgts(kw)%fly_dta(jpi1-1,:,:) = ref_wgts(kw)%fly_dta(jpi1,:,:)
[1275]	1421	ENDIF
	1422	IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN
[2528]	1423	ref_wgts(kw)%fly_dta(jpi2+1,:,:) = ref_wgts(kw)%fly_dta(jpi2,:,:)
[1275]	1424	ENDIF
	1425	IF( jpj1 == 2 ) THEN
[2528]	1426	ref_wgts(kw)%fly_dta(:,jpj1-1,:) = ref_wgts(kw)%fly_dta(:,jpj1,:)
[1275]	1427	ENDIF
	1428	IF( jpj2 + jpjmin - 1 == ref_wgts(kw)%ddims(2)+1 .AND. jpj2 .lt. jpjwid+2 ) THEN
[2528]	1429	ref_wgts(kw)%fly_dta(:,jpj2+1,:) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2,:) - ref_wgts(kw)%fly_dta(:,jpj2-1,:)
[1275]	1430	ENDIF
	1431
	1432	!! if data grid is cyclic we can do better on east-west edges
	1433	!! but have to allow for whether first and last columns are coincident
	1434	IF( ref_wgts(kw)%cyclic ) THEN
	1435	rec1(2) = MAX( jpjmin-1, 1 )
[2528]	1436	recn(1) = 1
[1275]	1437	recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 )
	1438	jpj1 = 2 + rec1(2) - jpjmin
	1439	jpj2 = jpj1 + recn(2) - 1
	1440	IF( jpi1 == 2 ) THEN
[2528]	1441	rec1(1) = ref_wgts(kw)%ddims(1) - ref_wgts(kw)%overlap
	1442	SELECT CASE( SIZE( ref_wgts(kw)%col(:,jpj1:jpj2,:),3) )
	1443	CASE(1)
	1444	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,1), nrec, rec1, recn)
	1445	CASE DEFAULT
	1446	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn)
	1447	END SELECT
	1448	ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:)
[1275]	1449	ENDIF
	1450	IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN
[2528]	1451	rec1(1) = 1 + ref_wgts(kw)%overlap
	1452	SELECT CASE( SIZE( ref_wgts(kw)%col(:,jpj1:jpj2,:),3) )
	1453	CASE(1)
	1454	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,1), nrec, rec1, recn)
	1455	CASE DEFAULT
	1456	CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn)
	1457	END SELECT
	1458	ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:)
[1275]	1459	ENDIF
	1460	ENDIF
	1461
	1462	! gradient in the i direction
	1463	DO jk = 1,4
	1464	DO jn = 1, jpj
	1465	DO jm = 1,jpi
	1466	ni = ref_wgts(kw)%data_jpi(jm,jn,jk)
	1467	nj = ref_wgts(kw)%data_jpj(jm,jn,jk)
[2528]	1468	dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+4) * 0.5 * &
	1469	(ref_wgts(kw)%fly_dta(ni+2,nj+1,:) - ref_wgts(kw)%fly_dta(ni,nj+1,:))
[1275]	1470	END DO
	1471	END DO
	1472	END DO
	1473
	1474	! gradient in the j direction
	1475	DO jk = 1,4
	1476	DO jn = 1, jpj
	1477	DO jm = 1,jpi
	1478	ni = ref_wgts(kw)%data_jpi(jm,jn,jk)
	1479	nj = ref_wgts(kw)%data_jpj(jm,jn,jk)
[2528]	1480	dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+8) * 0.5 * &
	1481	(ref_wgts(kw)%fly_dta(ni+1,nj+2,:) - ref_wgts(kw)%fly_dta(ni+1,nj,:))
[1275]	1482	END DO
	1483	END DO
	1484	END DO
	1485
[2715]	1486	! gradient in the ij direction
	1487	DO jk = 1,4
	1488	DO jn = 1, jpj
	1489	DO jm = 1,jpi
	1490	ni = ref_wgts(kw)%data_jpi(jm,jn,jk)
	1491	nj = ref_wgts(kw)%data_jpj(jm,jn,jk)
	1492	dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+12) * 0.25 * ( &
[2528]	1493	(ref_wgts(kw)%fly_dta(ni+2,nj+2,:) - ref_wgts(kw)%fly_dta(ni ,nj+2,:)) - &
	1494	(ref_wgts(kw)%fly_dta(ni+2,nj ,:) - ref_wgts(kw)%fly_dta(ni ,nj ,:)))
[2715]	1495	END DO
[1275]	1496	END DO
[2715]	1497	END DO
	1498	!
[1275]	1499	END IF
[2715]	1500	!
[1275]	1501	END SUBROUTINE fld_interp
[2528]	1502
	1503
	1504	FUNCTION ksec_week( cdday )
	1505	!!---------------------------------------------------------------------
	1506	!! * FUNCTION kshift_week *
	1507	!!
	1508	!! ** Purpose :
	1509	!!---------------------------------------------------------------------
	1510	CHARACTER(len=*), INTENT(in) :: cdday !3 first letters of the first day of the weekly file
	1511	!!
	1512	INTEGER :: ksec_week ! output variable
	1513	INTEGER :: ijul !temp variable
	1514	INTEGER :: ishift !temp variable
	1515	CHARACTER(len=3),DIMENSION(7) :: cl_week
	1516	!!----------------------------------------------------------------------
	1517	cl_week = (/"sun","sat","fri","thu","wed","tue","mon"/)
	1518	DO ijul = 1, 7
	1519	IF( cl_week(ijul) == TRIM(cdday) ) EXIT
[2715]	1520	END DO
[2528]	1521	IF( ijul .GT. 7 ) CALL ctl_stop( 'ksec_week: wrong day for sdjf%cltype(6:8): '//TRIM(cdday) )
	1522	!
	1523	ishift = ijul * NINT(rday)
	1524	!
	1525	ksec_week = nsec_week + ishift
	1526	ksec_week = MOD( ksec_week, 7*NINT(rday) )
	1527	!
	1528	END FUNCTION ksec_week
	1529
[2715]	1530	!!======================================================================
[888]	1531	END MODULE fldread

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

Download in other formats: