Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

fldread.F90 @ 11738

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

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

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

Context Navigation

source: branches/UKMO/dev_r5518_GO6_under_ice_relax_dr_hook/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90 @ 11738

Download in other formats: