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/2013/dev_r3987_UKMO6_C1D/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/2013/dev_r3987_UKMO6_C1D/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90 @ 4144

Last change on this file since 4144 was 4144, checked in by rfurner, 10 years ago
Commit for 2013 changes; see #1085
Property svn:keywords set to `Id`
File size: 72.6 KB

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

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

Download in other formats: