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.

bias.F90 in branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/ASM – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/ASM/bias.F90 @ 7089

Last change on this file since 7089 was 7089, checked in by isabella, 7 years ago
simplifiyng the code
File size: 43.9 KB

Line
1	MODULE bias
2	!! Is used by OPA and STEP
3	!!======================================================================
4	!! * Module bias *
5	!! Code to estimate and apply bias correction.
6	!! The bias is in T/S and Pressure. It follows the generalized
7	!! bias algorithm presented in Balmaseda et al 2007.
8	!!
9	!! It can be read from a file offline, estimated online from relaxation
10	!! terms or from assimilation increments (this latter estimtd in inner loop)
11	!!
12	!! The parameters for the time evolution of the bias can be different
13	!! from the relaxation terms and for the assim increments. Only the
14	!! parameter for the relaxtion terms are considered here.
15	!!
16	!! The offline bias term can contain the seasonal cycle.
17	!!
18	!! The time evolution of the bias for relaxtion is estimated as followed
19	!! bias_rlx(t+1)=t_rlx_membias_rlx(t)+t_rlx_upd(t/s)trdmp.
20	!!
21	!! The total bias in T/S is partion between the correction to T/S only
22	!! (fb_t) and the correction applied to the pressure gradient (fb_p).
23	!! We impose that (fb_p=1.-fb_t). These factors can be different for the
24	!! different terms(fb_t_rxl,fb_t_asm,fb_t_ofl)
25	!!
26	!! (t/s)bias = fb_t_ofl * (t/s)bias_ofl +
27	!! fb_t_rlx * (t/s)bias_rlx +
28	!! fb_t_asm * (t/s)bias_asm
29	!! (t/s)bias_p =fb_p_ofl * (t/s)bias_ofl+
30	!! fb_p_rlx * (t/s)bias_rlx_p +
31	!! fb_p_asm * (t/s)bias_asm_p
32	!! (t/s)bias is applied directely to correct T and S
33	!! (t/s)bias_p is used to compute rhd_pc and gru/v_pc
34	!!
35	!! Note: the above is an adhoc /simple way of making the partition
36	!! between bias in pressure and bias in T/S. It would be better
37	!! if the partition is done at the time of estimating the time
38	!! evolution of the bias. That would mean doubling the number of
39	!! 3D arrays.
40	!!
41	!! New addtion: (not in Balmaseda et al 2007):
42	!! A more physical alternative to the partition of the bias can be
43	!! done in terms of the inertial frequency: when the time scales of
44	!! adjustment are shorter that >1/f (Equator), the bias correction should
45	!! be in the the pressure term. Otherwise, it can act directly in T/S.
46	!! NOTE that the bias correction in the pressure term here (following
47	!! (Bell et al 2007) is just the "opposite" as the semi-prognostic method
48	!! in Greatbatch et al 2004.
49	!! The use of this partition is controlled by ln_inertial=.true.
50	!!
51	!!
52	!! 2009-03 (M.A. Balmaseda ECMWF)
53	!!======================================================================
54
55	!!----------------------------------------------------------------------
56	!! bias_init : Read in the bias namelist and the bias arrays
57	!! tra_bias : Apply the bias fields on T/S directly
58	!! dyn_bias : Compute density correction for dynamic hpg
59	!! bias_opn : open bias files for restart capabilities
60	!! bias_wrt : write bias fies " " "
61	!!----------------------------------------------------------------------
62	!! * Modules used
63	USE par_kind, ONLY: &
64	& wp
65	USE par_oce, ONLY: &
66	& jpi, &
67	& jpj, &
68	& jpk
69	USE dom_oce, ONLY: &
70	& rdt, &
71	& ln_zps, &
72	& gphit
73	USE phycst, ONLY: &
74	& rday, &
75	& rad
76	USE oce, ONLY: &
77	& tsb, tsn, tsa, &
78	& rhop, &
79	& gtsu, gtsv
80	USE dynhpg, ONLY: &
81	& ln_dynhpg_imp
82	USE tradmp
83	USE dtatsd, ONLY: &
84	& ln_tsd_tradmp
85	USE in_out_manager, ONLY: &
86	& lwp, &
87	& numnam_ref, &
88	& numnam_cfg, &
89	& numond, &
90	& numout, &
91	& lrst_oce, &
92	& nit000
93	USE iom
94	USE eosbn2
95	USE zpshde ! partial step: hor. derivative (zps_hde routine)
96	USE biaspar
97	USE fldread ! read input fields
98	USE lbclnk ! lateral boundary conditions (or mpp link)
99	USE asmpar
100	USE asminc
101	USE lib_mpp, ONLY: &
102	& ctl_stop, &
103	& ctl_nam
104
105	IMPLICIT NONE
106
107	!! * Routine accessibility
108	PRIVATE
109	PUBLIC bias_init, & !: Read in the bias namelist and the bias arrays
110	& tra_bias, & !: Estimate/apply bias on traces
111	& dyn_bias, & !: " density correction for pressure gradient.
112	& bias_opn, &
113	& bias_wrt
114
115
116	!! * Shared variables
117	!! * Private module variables
118	REAL(wp), PRIVATE :: &
119	& bias_time_unit_asm, & !: bias_asm units in s ( per day = 86400 s)
120	& bias_time_unit_rlx, & !: bias_rlx units in s ( 1 second)
121	& bias_time_unit_ofl, & !: bias_ofl units in s ( 1 second)
122	& t_rlx_mem, & !: time param for mem in bias_rlx model
123	& t_rlx_upd, & !: time param for update in bias_rlx model
124	!: (pct of incr for computation of bias)
125	& t_asm_mem, & !: time param for mem in bias_asm model
126	& t_asm_upd, & !: time param for update in bias_asm model
127	!: (pct of incr for computation of bias)
128	& fb_t_rlx, & !: parition of bias in T for rlx bias term
129	& fb_t_asm, & !: parition of bias in T for asm bias term
130	& fb_t_ofl, & !: parition of bias in T for ofl bias term
131	& fb_p_rlx, & !: parition of bias in P for rlx bias term
132	& fb_p_asm, & !: parition of bias in P for asm bias term
133	& fb_p_ofl, & !: parition of bias in P for ofl bias term
134	& fctamp, & !: amplification factor for T if inertial
135	& rn_maxlat_bias, & !: Max lat for latitudinal ramp
136	& rn_minlat_bias !: Min lat for latitudinal ramp
137
138	LOGICAL, PRIVATE :: lalloc
139	REAL(wp), PRIVATE, DIMENSION(:,:,:), ALLOCATABLE :: &
140	& tbias_asm, & !: Temperature bias field
141	& sbias_asm, & !: Salinity bias field
142	& tbias_rlx, & !: Temperature bias field
143	& sbias_rlx, & !: Salinity bias field
144	& tbias_asm_out, & !: Output temperature bias field
145	& sbias_asm_out, & !: Output salinity bias field
146	& tbias_rlx_out, & !: Output temperature bias field
147	& sbias_rlx_out, & !: Output salinity bias field
148	& tbias_p_out, & !: Output temperature bias field for P correction
149	& sbias_p_out, & !: Output salinity bias field for P correction
150
151	INTEGER, PRIVATE :: nn_lat_ramp ! choice of latitude dependent ramp
152	! for the pressure correction.
153	! 1:inertial ramp, 2:linear ramp, else:no ramp
154	LOGICAL, PRIVATE :: ln_bsyncro ! syncronous or assincrous bias correction
155	LOGICAL, PRIVATE :: ln_itdecay ! evolve bias correction at every time step.
156	LOGICAL, PRIVATE :: ln_incpc ! incremental pressure correction
157
158	REAL(wp), PRIVATE, DIMENSION(:,:), ALLOCATABLE :: fbcoef
159	REAL(wp), PRIVATE, DIMENSION(:,:), ALLOCATABLE :: fbcoef_stscale
160
161	INTEGER, PRIVATE :: &
162	& numbias_asm, & ! netcdf id of bias file from assim
163	& numbias_tot, & ! netcdf id of bias file with total bias
164	& nn_bias_itwrt ! time step for outputting bias pressure corr
165
166	CHARACTER(LEN=128), PRIVATE :: &
167	& cn_bias_asm, & ! name of bias file from assim
168	& cn_bias_tot ! name of bias with total/rlx bias
169
170	! Structure of input T and S bias offline (file informations, fields read)
171	TYPE(FLD), PRIVATE, ALLOCATABLE, DIMENSION(:) :: sf_tbias_ofl
172	TYPE(FLD), PRIVATE, ALLOCATABLE, DIMENSION(:) :: sf_sbias_ofl
173
174	TYPE(FLD_N), PRIVATE ::& ! information about the fields to be read
175	& sn_tbias_ofl, sn_sbias_ofl
176
177	CONTAINS
178
179	SUBROUTINE bias_init
180	!!----------------------------------------------------------------------
181	!! * ROUTINE bias_init *
182	!!
183	!! ** Purpose : Read in the bias namelist and read in the bias fields.
184	!!
185	!! ** Method : Read in the bias namelist and read in the bias fields.
186	!!
187	!! ** Action :
188	!!
189	!! History :
190	!! ! 08-05 (D. Lea) Initial version
191	!! ! 08-10 (M. Martin) Tidy
192	!! ! 09-03 (M. Balmaseda). Generalize to estimate the bias
193	!! from relax and offline bias term.
194	!! Introduce parameters to control the
195	!! model for the bias
196	!! (variables and time evolution)
197	!!----------------------------------------------------------------------
198
199	IMPLICIT NONE
200
201	!! * Local declarations
202
203	CHARACTER(len=100) :: cn_dir ! dir for location ofline bias
204	INTEGER :: ierror
205	INTEGER :: ios ! Local integer output status for namelist read
206	REAL(wp) :: eft_rlx, & ! efolding time (bias memory) in days
207	& eft_asm, & ! " "
208	& log2, &
209	& lenscl_bias !lengthscale of the pressure bias decay between minlat and maxlat.
210
211	NAMELIST/nambias/ ln_bias, ln_bias_asm, ln_bias_rlx, ln_bias_ofl, &
212	& ln_bias_ts_app, ln_bias_pc_app, &
213	& fb_t_asm, fb_t_rlx, fb_t_ofl, fb_p_asm, fb_p_rlx, fb_p_ofl, &
214	& eft_rlx, t_rlx_upd, eft_asm, t_asm_upd, nn_lat_ramp, &
215	& bias_time_unit_asm, bias_time_unit_rlx, bias_time_unit_ofl, &
216	& cn_bias_tot, cn_bias_asm, cn_dir, sn_tbias_ofl, sn_sbias_ofl, &
217	& ln_bsyncro, fctamp, rn_maxlat_bias, rn_minlat_bias, &
218	& nn_bias_itwrt, ln_itdecay, ln_incpc
219
220
221	!-----------------------------------------------------------------------
222	! Read Namelist : bias interface
223	!-----------------------------------------------------------------------
224
225
226	REWIND( numnam_ref ) ! Namelist nambias in reference namelist : Bias pressure correction
227	READ ( numnam_ref, nambias, IOSTAT = ios, ERR = 901)
228	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambias in reference namelist', lwp )
229
230
231	! Set additional default values (note that most values are set in the reference namelist)
232
233	IF ( ln_asmiau ) nn_bias_itwrt = nitiaufin
234
235	! ... default values (NB: frequency positive => hours, negative => months)
236	! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or !
237	! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' !
238	sn_tbias_ofl = FLD_N( 'tbias_ofl' , -1. , 'tbias' , .TRUE. , .FALSE. , 'yearly', '', '', '' )
239	sn_sbias_ofl = FLD_N( 'sbias_ofl' , -1. , 'sbias' , .TRUE. , .FALSE. , 'yearly', '', '', '' )
240
241
242	REWIND( numnam_cfg ) ! Namelist nambias in configuration namelist : Bias pressure correction
243	READ ( numnam_cfg, nambias, IOSTAT = ios, ERR = 902 )
244	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambias in configuration namelist', lwp )
245	IF(lwm) WRITE ( numond, nambias )
246
247
248	IF ( ( .NOT. ln_bias_asm ) .AND. ( .NOT. ln_bias_ofl ) .AND. ( .NOT. ln_bias_rlx ) ) THEN
249	ln_bias_ts_app = .FALSE.
250	ln_bias_pc_app = .FALSE.
251	ln_bias = .FALSE.
252	ENDIF
253
254	! set up decay scales
255	log2 = LOG( 2.0_wp )
256	t_rlx_mem = EXP( - log2 * rdt / ( eft_rlx * rday ) )
257	t_asm_mem = EXP( - log2 * bias_time_unit_asm/ ( eft_asm * rday ) )
258
259	! Control print
260	IF(lwp) THEN
261	WRITE(numout,*)
262	WRITE(numout,*) 'bias_init : '
263	WRITE(numout,*) '~~~~~~~~~~~ '
264	WRITE(numout,*) ' Namelist nambias : '
265
266	WRITE(numout,*) ' Bias switches/options/variables '
267	WRITE(numout,*) ' bias main switch ln_bias = ',ln_bias
268	WRITE(numout,*) ' bias from assim ln_bias_asm = ',ln_bias_asm
269	WRITE(numout,*) ' bias from relax ln_bias_rlx = ',ln_bias_rlx
270	WRITE(numout,*) ' bias from offln ln_bias_ofl = ',ln_bias_ofl
271	WRITE(numout,*) ' bias T and S apply ln_bias_ts_app = ',ln_bias_ts_app
272	WRITE(numout,*) ' bias pressure correctn apply ln_bias_pc_app = ',ln_bias_pc_app
273	WRITE(numout,*) ' bias pressure correctn apply ln_bias_pc_app = ',ln_bias_pc_app
274	WRITE(numout,*) ' lat ramp for bias correction nn_lat_ramp = ',nn_lat_ramp
275	WRITE(numout,*) ' time step for writing bias fld nn_bias_itwrt = ',nn_bias_itwrt
276	WRITE(numout,*) ' evolve pcbias at each timestep ln_itdecay = ',ln_itdecay
277	WRITE(numout,*) ' incremental press. correction ln_incpc = ',ln_incpc
278	WRITE(numout,*) ' Parameters for parition of bias term '
279	WRITE(numout,*) ' fb_t_rlx = ',fb_t_rlx
280	WRITE(numout,*) ' fb_t_asm = ',fb_t_asm
281	WRITE(numout,*) ' fb_t_ofl = ',fb_t_ofl
282	WRITE(numout,*) ' fb_p_rlx = ',fb_p_rlx
283	WRITE(numout,*) ' fb_p_asm = ',fb_p_asm
284	WRITE(numout,*) ' fb_p_ofl = ',fb_p_ofl
285	WRITE(numout,*) ' Parameters for time evolution of bias '
286	WRITE(numout,) ' Rlx efolding time (mem) eft_rlx,t_rlx_mem = ', eft_rlx, t_rlx_mem, 1. - log2 rdt / (eft_rlx * rday)
287	WRITE(numout,*) ' uptdate factor t_rlx_upd = ',t_rlx_upd
288	WRITE(numout,) ' Asm efolding time (mem) eft_asm,t_asm_mem = ', eft_asm, t_asm_mem, 1. - log2 rdt / (eft_asm * rday)
289	WRITE(numout,*) ' uptdate factor t_asm_upd = ',t_asm_upd
290	WRITE(numout,*) ' Filenames and input structures'
291	WRITE(numout,*) ' bias_tot filename cn_bias_to = ',cn_bias_tot
292	WRITE(numout,*) ' bias_asm filename cn_bias_asm = ',cn_bias_asm
293	WRITE(numout,*) ' bias_asm time unit (secs) bias_time_unit_asm = ',bias_time_unit_asm
294	WRITE(numout,*) ' structure Tem bias ofl sn_tbias_ofl = ',sn_tbias_ofl
295	WRITE(numout,*) ' structure Sal bias ofl sn_sbias_ofl = ',sn_sbias_ofl
296
297	IF ( ( (.NOT. ln_tsd_tradmp) .OR. (.NOT. ln_tradmp) ) .AND. ln_bias_rlx ) &
298	& CALL ctl_stop (' lk_dtatem, lk_dtasal and lk_tradmp need to be true with ln_bias_rlx' )
299
300	ENDIF
301	IF( .NOT. ln_bias ) RETURN
302
303	IF( .NOT. lalloc ) THEN
304
305	ALLOCATE( tbias(jpi,jpj,jpk) , &
306	& sbias(jpi,jpj,jpk) , &
307	& tbias_p(jpi,jpj,jpk), &
308	& sbias_p(jpi,jpj,jpk), &
309	& rhd_pc(jpi,jpj,jpk) , &
310	& gru_pc(jpi,jpj) , &
311	& grv_pc(jpi,jpj) )
312
313	ALLOCATE( fbcoef(jpi,jpj), fbcoef_stscale(jpi,jpj) )
314
315	IF( ln_bias_asm ) ALLOCATE( tbias_asm(jpi,jpj,jpk), &
316	& sbias_asm(jpi,jpj,jpk), &
317	tbias_asm_out(jpi,jpj,jpk), &
318	sbias_asm_out(jpi,jpj,jpk), &
319	tbias_p_out(jpi,jpj,jpk), &
320	sbias_p_out(jpi,jpj,jpk) )
321
322	IF( ln_bias_rlx ) ALLOCATE( tbias_rlx(jpi,jpj,jpk), &
323	& sbias_rlx(jpi,jpj,jpk), &
324	tbias_rlx_out(jpi,jpj,jpk), &
325	sbias_rlx_out(jpi,jpj,jpk) )
326
327	IF( ln_incpc ) ALLOCATE( tbias_asm_stscale(jpi,jpj,jpk), &
328	& sbias_asm_stscale(jpi,jpj,jpk))
329
330	lalloc = .TRUE.
331
332	ENDIF
333
334	IF( ln_bias_ofl ) THEN ! set sf_tbias_ofl and sf_sbias_ofl strctrs
335	!
336	! tbias
337	!
338	ALLOCATE( sf_tbias_ofl(1), STAT=ierror )
339	IF( ierror > 0 ) THEN
340	CALL ctl_stop( 'bias_init: unable to allocate sf_tbias_ofl structure' ) ; RETURN
341	ENDIF
342	ALLOCATE( sf_tbias_ofl(1)%fnow(jpi,jpj,jpk) )
343	ALLOCATE( sf_tbias_ofl(1)%fdta(jpi,jpj,jpk,2) )
344
345	! fill structure with values and control print
346	CALL fld_fill( sf_tbias_ofl, (/ sn_tbias_ofl /), cn_dir, 'bias_init', 'Offline T bias term ', 'nam_tbias_ofl' )
347	!
348	! salinity bias
349	!
350	ALLOCATE( sf_sbias_ofl(1), STAT=ierror )
351	IF( ierror > 0 ) THEN
352	CALL ctl_stop( 'bias_init: unable to allocate sf_sbias_ofl structure' ) ; RETURN
353	ENDIF
354	ALLOCATE( sf_sbias_ofl(1)%fnow(jpi,jpj,jpk) )
355	ALLOCATE( sf_sbias_ofl(1)%fdta(jpi,jpj,jpk,2) )
356
357	! fill structure with values and control print
358	CALL fld_fill( sf_sbias_ofl, (/ sn_sbias_ofl /), cn_dir, 'bias_init', 'Offline S bias term ', 'nam_sbias_ofl' )
359
360	ENDIF
361
362	! Read total bias
363	IF ( ln_bias ) THEN
364	tbias(:,:,:) = 0.0_wp
365	sbias(:,:,:) = 0.0_wp
366	tbias_p(:,:,:) = 0.0_wp
367	sbias_p(:,:,:) = 0.0_wp
368	gru_pc(:,:) = 0.0_wp
369	grv_pc(:,:) = 0.0_wp
370	IF ( ln_bias_rlx ) THEN
371	tbias_rlx(:,:,:) = 0.0_wp
372	sbias_rlx(:,:,:) = 0.0_wp
373	ENDIF
374	IF ( ln_bias_asm ) THEN !now rlx and asm bias in same file
375	tbias_asm(:,:,:) = 0.0_wp
376	sbias_asm(:,:,:) = 0.0_wp
377	ENDIF
378
379	IF ( ln_incpc ) THEN !incr pressure correction
380	tbias_asm_stscale(:,:,:) = 0.0_wp
381	sbias_asm_stscale(:,:,:) = 0.0_wp
382	ENDIF
383
384
385	numbias_tot = 0
386	! Get bias from file and prevent fail if the file does not exist
387	IF(lwp) WRITE(numout,*) 'Opening ',TRIM( cn_bias_tot )
388	CALL iom_open( cn_bias_tot, numbias_tot, ldstop=.FALSE. )
389
390	IF ( numbias_tot > 0 ) THEN
391	! Could check validity time of bias fields here...
392	! Get the T and S bias data
393	IF(lwp) WRITE(numout,*) 'Reading bias fields from tot...'
394
395	!Search for bias from relaxation term if needed. Use same file
396	IF ( ln_bias_rlx ) THEN
397	IF(lwp) WRITE(numout,*) 'Reading bias fields for bias rlx from file ',cn_bias_tot
398	IF( iom_varid( numbias_tot, 'tbias_rlx' ) > 0 ) THEN
399	! Get the T and S bias data
400	CALL iom_get( numbias_tot, jpdom_autoglo, 'tbias_rlx', tbias_rlx )
401	CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_rlx', sbias_rlx )
402	ELSE
403	CALL ctl_stop( 'Bias relaxation variables not found in ',cn_bias_tot )
404	ENDIF
405	ENDIF
406
407
408	!Search for bias from assim term if needed. Use same file
409	IF ( ln_bias_asm .and. .not.ln_bsyncro ) THEN
410	IF(lwp) WRITE(numout,*) 'Reading a-syncro bias fields for bias asm from file ',cn_bias_tot
411	IF( iom_varid( numbias_tot, 'tbias_asm' ) > 0 ) THEN
412	! Get the T and S bias data
413	CALL iom_get( numbias_tot, jpdom_autoglo, 'tbias_asm', tbias_asm )
414	CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_asm', sbias_asm )
415	ELSE
416	CALL ctl_stop( 'Bias assim variables not found in ',cn_bias_tot )
417	ENDIF
418
419
420	IF ( ln_incpc ) THEN
421	IF(lwp) WRITE(numout,*) 'Reading short time scale bias correction fields for bias asm from file ',cn_bias_tot
422	IF( iom_varid( numbias_tot, 'tbias_asm_stscale' ) > 0 ) THEN
423	! Get the T and S bias data
424	CALL iom_get( numbias_tot, jpdom_autoglo, 'tbias_asm_stscale', tbias_asm_stscale )
425	CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_asm_stscale', sbias_asm_stscale )
426	ELSE
427	CALL ctl_stop( 'Short time scale bias assim variables not found in ',cn_bias_tot )
428	ENDIF
429	ENDIF
430
431	ENDIF
432
433
434
435	! Close the file
436	CALL iom_close(numbias_tot)
437
438	ELSE
439	IF(lwp) WRITE(numout,*) 'No bias file found so T and S bias fields are set to zero'
440	ENDIF
441
442	ENDIF
443
444	! for the time being, the bias_asm is read in the same file as
445	! bias_rlx
446	! Implications: bias_asm is estimated/evolved in time in the second outer
447	! loop only, when the assimilation increments are ready.
448	! bias_asm is kept cte during the first outer loop.
449	! => Assyncronous bias correction.
450	! Alternative: Syncronous bias correction:
451	! bias_asm estimated/evolved in the first outer loop
452	! with the asm incrments of the previous cycle.
453	! bias_asm kept cte during the second outer loop.
454	! Implication: bias_asm should be estimated really in the
455	! inner loop.
456	IF ( ln_bsyncro ) THEN
457	! Read bias from assimilation from a separate file
458	IF ( ln_bias_asm ) THEN
459	tbias_asm(:,:,:) = 0.0_wp
460	sbias_asm(:,:,:) = 0.0_wp
461	numbias_asm = 0
462	! Get bias from file and prevent fail if the file does not exist
463	IF(lwp) WRITE(numout,*) 'Opening file for syncro assim bias ',TRIM( cn_bias_asm )
464	CALL iom_open( cn_bias_asm, numbias_asm, ldstop=.FALSE. )
465
466	IF ( numbias_asm > 0 ) THEN
467	! Could check validity time of bias fields here...
468
469	! Get the T and S bias data
470	IF(lwp) WRITE(numout,*) 'Reading syncro bias fields from asm from file ',cn_bias_asm
471	CALL iom_get( numbias_asm, jpdom_autoglo, 'tbias_asm', tbias_asm )
472	CALL iom_get( numbias_asm, jpdom_autoglo, 'sbias_asm', sbias_asm )
473
474	! this is only applicable if tbias_asm were to be calculated in the inner loop
475	tbias_asm(:,:,:) = tbias_asm(:,:,:) * rdt / bias_time_unit_asm
476	sbias_asm(:,:,:) = sbias_asm(:,:,:) * rdt / bias_time_unit_asm
477
478	! Close the file
479	CALL iom_close(numbias_asm)
480
481	ELSE
482	IF(lwp) WRITE(numout,*) 'No bias file found from asm so T and S bias fields are set to zero'
483	ENDIF
484
485	ENDIF
486
487	ENDIF
488
489	!latitudinal dependence of partition coeficients. Adhoc
490	IF ( nn_lat_ramp == 1 ) THEN
491	! Use the inertial ramp.
492	lenscl_bias = ( rn_maxlat_bias - rn_minlat_bias )*2._wp
493	WHERE ( abs( gphit(:,:) ) <= rn_minlat_bias )
494	fbcoef(:,:) = 0._wp
495	ELSEWHERE ( abs( gphit(:,:) ) >= rn_maxlat_bias )
496	fbcoef(:,:) = 1._wp
497	ELSEWHERE
498	fbcoef(:,:) = 1._wp - exp( -( abs( gphit(:,:) ) - rn_minlat_bias ) &
499	* ( abs( gphit(:,:) ) - rn_minlat_bias ) / lenscl_bias )
500	ENDWHERE
501	ELSEIF ( nn_lat_ramp == 2 ) THEN
502	! Use a linear ramp consist with the geostrophic velocity balance ramp in NEMOVAR
503
504	WHERE ( abs( gphit(:,:) ) <= rn_minlat_bias )
505	fbcoef(:,:) = 0._wp
506	ELSEWHERE ( abs( gphit(:,:) ) >= rn_maxlat_bias )
507	fbcoef(:,:) = 1._wp
508	ELSEWHERE
509	fbcoef(:,:) = 1._wp - ((rn_maxlat_bias - abs( gphit(:,:)))/(rn_maxlat_bias - rn_minlat_bias))
510	ENDWHERE
511	ELSE
512	fbcoef(:,:) = 0.0_wp
513	fctamp = 0.0_wp
514	fbcoef_stscale(:,:) = 0.0_wp
515	ENDIF
516
517
518	IF ( ln_incpc) THEN
519	! not sure if this should be a special case of nn_lat_ramp == 2
520	minlat_bias = 3.0_wp
521	maxlat_bias = 8.0_wp
522	WHERE ( abs( gphit(:,:) ) <= minlat_bias )
523	fbcoef_stscale(:,:)=0._wp
524	ELSEWHERE ( abs( gphit(:,:) ) >= maxlat_bias )
525	fbcoef_stscale(:,:)=1._wp
526	ELSEWHERE
527	fbcoef_stscale(:,:)=1._wp - ((maxlat_bias - abs( gphit(:,:)))/(maxlat_bias-minlat_bias))
528	ENDWHERE
529	ENDIF
530
531
532	END SUBROUTINE bias_init
533
534	SUBROUTINE tra_bias ( kt )
535	!!----------------------------------------------------------------------
536	!! * ROUTINE tra_bias *
537	!!
538	!! ** Purpose : Update bias field every time step
539	!!
540	!! ** Method : add contributions to bias from 3 terms
541	!!
542	!! ** Action : Bias from assimilation (read in bias_init)
543	!! Bias from relaxation term is estimated according to
544	!! the prescribed time evolution of the bias
545	!! Bias from ofl term is read from external file
546	!! The difference contributions are added and the partition
547	!! into direct bias in T/S and pressure perfomed.
548	!!
549	!! History : 09-03 (M. Balmaseda)
550	!!----------------------------------------------------------------------
551	!! called every timestep after dta_sst if ln_bias true.
552
553	IMPLICIT NONE
554
555	!! * Arguments
556	INTEGER, INTENT(IN) :: kt ! ocean time-step index
557	!! * Local variables
558	INTEGER :: ji,jj,jk, it ! local loop index
559	REAL(wp) :: tsclf ! time saling factor
560	REAL(wp) :: fb_t_asm_max, fb_t_rlx_max, fb_t_ofl_max
561	REAL(wp) :: ztfrac, ztsday
562	REAL(wp) :: zfrac, zfrac1 ! temporal weights for inst pcbias (names could be changed)
563	REAL(wp) :: ztscale ! reduce the inst pcbias by this amount per 24 hours
564	REAL(wp) :: zwgt ! Weight for the inst pcorr term
565	REAL(wp) :: zdecay ! used in inst pcorr
566	REAL(wp), DIMENSION(jpi,jpj) :: zcof1, zcof2
567
568	IF ( .NOT. ln_bias ) RETURN
569	fb_t_rlx_max = MIN(fb_t_rlx*fctamp,1.0_wp)
570	fb_t_asm_max = MIN(fb_t_asm*fctamp,1.0_wp)
571	fb_t_ofl_max = MIN(fb_t_ofl*fctamp,1.0_wp)
572
573	tbias(:,:,:) = 0.0_wp
574	sbias(:,:,:) = 0.0_wp
575	tbias_p(:,:,:) = 0.0_wp
576	sbias_p(:,:,:) = 0.0_wp
577
578	ztscale = 0.1_wp
579	zwgt = 1.0_wp
580
581	IF ( ln_bias_asm ) THEN
582
583	tsclf = 1
584	IF ( .NOT.ln_bsyncro ) tsclf = rdt / bias_time_unit_asm
585	zcof1(:,:) = tsclf * ( ( 1.0_wp - fbcoef(:,:) ) * fb_t_asm + &
586	& fbcoef(:,:) * fb_t_asm_max )
587	zcof2(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_p_asm
588
589	IF ( ln_itdecay ) THEN !decay the pressure correction at each time step
590
591	ztsday = rday / real(rdt)
592
593	zdecay = (1-ztscale)**(1/real(ztsday)) ! used in ipc
594	zfrac1 = max(0.0_wp, zdecay**real(kt)) ! used in ipc
595
596
597	IF( ln_asmiau .and. ln_trainc ) THEN ! nudge in increments and decay historical contribution
598
599
600	IF ( kt <= nitiaufin ) THEN ! During IAU calculate the fraction of increments to apply at each time step
601
602	ztfrac = real(kt) / real(nitiaufin) ! nudging factor during the IAU
603
604
605	IF (lwp) THEN
606	WRITE(numout,*) 'tra_bias : bias weights'
607	WRITE(numout,*) '~~~~~~~~~~~~'
608	WRITE(numout,* ) "proportion of increment applied in pcbias ",ztfrac
609	WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(real(kt)/ztsday)
610	ENDIF
611
612	DO jk = 1, jpkm1
613	tbias(:,:,jk) = tbias(:,:,jk) + &
614	& ( t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) + &
615	& ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:)
616	sbias(:,:,jk) = sbias(:,:,jk) + &
617	& ( t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) + &
618	& ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:)
619	tbias_p(:,:,jk) = tbias_p(:,:,jk) + &
620	& ( t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) + &
621	& ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:)
622	sbias_p(:,:,jk) = sbias_p(:,:,jk) + &
623	& ( t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) + &
624	& ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:)
625	ENDDO
626
627
628	IF (ln_incpc) THEN
629
630	DO jk = 1, jpkm1
631	tbias_p(:,:,jk) = tbias_p(:,:,jk) + ( t_bkginc(:,:,jk) * zwgt * ztfrac * (1.0 - fbcoef_stscale(:,:)) ) &
632	& - ( tbias_asm_stscale(:,:,jk) * zfrac1 * (1.0 - fbcoef_stscale(:,:)) )
633	sbias_p(:,:,jk) = sbias_p(:,:,jk) + ( s_bkginc(:,:,jk) * zwgt * ztfrac * (1.0 - fbcoef_stscale(:,:)) ) &
634	& - ( sbias_asm_stscale(:,:,jk) * zfrac1 * (1.0 - fbcoef_stscale(:,:)) )
635	ENDDO
636
637	ENDIF
638
639	IF ( .not.ln_bsyncro ) THEN
640	IF ( kt == nn_bias_itwrt ) THEN
641	DO jk = 1, jpkm1
642	tbias_asm_out(:,:,jk) = t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) + &
643	& ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk)
644	sbias_asm_out(:,:,jk) = t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) + &
645	& ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk)
646	END DO
647	ENDIF
648	ENDIF
649
650	! update the historical component with the increments at the end of the IAU
651	IF ( kt == nitiaufin ) THEN
652	DO jk = 1, jpkm1
653	tbias_asm(:,:,jk) = t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) + &
654	& ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk)
655	sbias_asm(:,:,jk) = t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) + &
656	& ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk)
657	END DO
658	ENDIF
659
660	ELSE ! decay pressure correction from combined historical component and increments after IAU
661
662	ztfrac=t_asm_mem**(real(kt-nitiaufin)/real(nitiaufin)) ! decay from end of IAU
663
664
665	DO jk = 1, jpkm1
666	tbias(:,:,jk) = tbias(:,:,jk) + &
667	& ( ztfrac * tbias_asm(:,:,jk) ) * zcof1(:,:)
668	sbias(:,:,jk) = sbias(:,:,jk) + &
669	& ( ztfrac * sbias_asm(:,:,jk) ) * zcof1(:,:)
670	tbias_p(:,:,jk) = tbias_p(:,:,jk) + &
671	& ( ztfrac * tbias_asm(:,:,jk) ) * zcof2(:,:)
672	sbias_p(:,:,jk) = sbias_p(:,:,jk) + &
673	& ( ztfrac * sbias_asm(:,:,jk) ) * zcof2(:,:)
674
675	ENDDO
676
677	IF (ln_incpc) THEN
678
679	zfrac = max(0.0_wp, zdecay**real(kt-nitiaufin) )
680
681	DO jk = 1, jpkm1
682	tbias_p(:,:,jk) = tbias_p(:,:,jk) + ( t_bkginc(:,:,jk) * zwgt * zfrac * (1.0 - fbcoef_stscale(:,:)) ) &
683	& - ( tbias_asm_stscale(:,:,jk) * zfrac1 * (1.0 - fbcoef_stscale(:,:)) )
684	sbias_p(:,:,jk) = sbias_p(:,:,jk) + ( s_bkginc(:,:,jk) * zwgt * zfrac * (1.0 - fbcoef_stscale(:,:)) ) &
685	& - ( sbias_asm_stscale(:,:,jk) * zfrac1 * (1.0 - fbcoef_stscale(:,:)) )
686	ENDDO
687
688	ENDIF
689
690	IF (lwp) THEN
691	WRITE(numout,*) 'tra_bias : bias weights'
692	WRITE(numout,*) '~~~~~~~~~~~~'
693	WRITE(numout,* ) "proportion of increments + historical pcbias applied ",ztfrac
694	ENDIF
695
696	IF ( ln_trainc .and. .not.ln_bsyncro ) THEN
697	IF ( kt == nn_bias_itwrt ) THEN
698	DO jk = 1, jpkm1
699	tbias_asm_out(:,:,jk) = ztfrac * tbias_asm(:,:,jk)
700	sbias_asm_out(:,:,jk) = ztfrac * sbias_asm(:,:,jk)
701	END DO
702	ENDIF
703	ENDIF
704
705	ENDIF
706
707
708	ELSE ! no assimilation increments, simply decay pressure correction (e.g for forecasts)
709
710	DO jk = 1, jpkm1
711	tbias(:,:,jk) = tbias(:,:,jk) + &
712	& ( t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) ) * zcof1(:,:)
713	sbias(:,:,jk) = sbias(:,:,jk) + &
714	& ( t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) ) * zcof1(:,:)
715	tbias_p(:,:,jk) = tbias_p(:,:,jk) + &
716	& ( t_asm_mem*(real(kt)/ztsday) tbias_asm(:,:,jk) ) * zcof2(:,:)
717	sbias_p(:,:,jk) = sbias_p(:,:,jk) + &
718	& ( t_asm_mem*(real(kt)/ztsday) sbias_asm(:,:,jk) ) * zcof2(:,:)
719	ENDDO
720
721	IF (lwp) THEN
722	WRITE(numout,*) 'tra_bias : bias weights'
723	WRITE(numout,*) '~~~~~~~~~~~~'
724	WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(real(kt)/ztsday)
725	ENDIF
726
727	ENDIF
728
729	ELSE ! maintain a constant pressure correction
730
731	DO jk = 1, jpkm1
732	tbias(:,:,jk) = tbias(:,:,jk) + tbias_asm(:,:,jk) * zcof1(:,:)
733	sbias(:,:,jk) = sbias(:,:,jk) + sbias_asm(:,:,jk) * zcof1(:,:)
734	tbias_p(:,:,jk) = tbias_p(:,:,jk) + tbias_asm(:,:,jk) * zcof2(:,:)
735	sbias_p(:,:,jk) = sbias_p(:,:,jk) + sbias_asm(:,:,jk) * zcof2(:,:)
736	END DO
737
738	IF( ln_asmiau .and. ln_trainc .and. .not.ln_bsyncro ) THEN
739	! if last outer loop (ln_asmiau=true and ln_trainc=true). t/sbias_asm
740	! is updated, only once (end of run) taking into account units.
741	IF ( kt == nn_bias_itwrt ) THEN
742	IF(lwp) WRITE(numout,*)' estimating asm bias at timestep: ',kt
743	DO jk = 1, jpkm1
744	tbias_asm_out(:,:,jk) = t_asm_mem * tbias_asm(:,:,jk) + &
745	& t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk)
746	sbias_asm_out(:,:,jk) = t_asm_mem * sbias_asm(:,:,jk) + &
747	& t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk)
748	END DO
749	ENDIF
750	ENDIF
751
752	ENDIF
753
754	ENDIF
755
756
757	#if defined key_tradmp
758	! Time evolution of bias from relaxation
759	IF ( ln_bias_rlx ) THEN
760	tbias_rlx(:,:,:) = t_rlx_mem * tbias_rlx(:,:,:) + &
761	& t_rlx_upd * ttrdmp(:,:,:) * rdt
762	sbias_rlx(:,:,:) = t_rlx_mem * sbias_rlx(:,:,:) + &
763	& t_rlx_upd * strdmp(:,:,:) * rdt
764	zcof1(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_t_rlx + &
765	& fbcoef(:,:) * fb_t_rlx_max
766	zcof2(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_p_rlx
767	DO jk = 1, jpkm1
768	tbias(:,:,jk) = tbias(:,:,jk) + tbias_rlx(:,:,jk) * zcof1(:,:)
769	sbias(:,:,jk) = sbias(:,:,jk) + sbias_rlx(:,:,jk) * zcof1(:,:)
770	tbias_p(:,:,jk) = tbias_p(:,:,jk) + tbias_rlx(:,:,jk) * zcof2(:,:)
771	sbias_p(:,:,jk) = sbias_p(:,:,jk) + sbias_rlx(:,:,jk) * zcof2(:,:)
772	ENDDO
773	IF ( kt == nn_bias_itwrt ) THEN
774	tbias_rlx_out(:,:,:) = tbias_rlx(:,:,:)
775	sbias_rlx_out(:,:,:) = sbias_rlx(:,:,:)
776	ENDIF
777	ENDIF
778	#endif
779	! ofline bias
780	IF ( kt == nit000 ) THEN
781	IF(lwp) WRITE(numout,*) ' tra_bias: ln_bias_ofl = ',ln_bias_ofl
782	IF(lwp) WRITE(numout,*) ' ~~~~~~~~~'
783	ENDIF
784	IF ( ln_bias_ofl ) THEN
785	IF(lwp) WRITE(numout,*) 'reading offline bias'
786	CALL fld_read( kt, 1, sf_tbias_ofl )
787	CALL fld_read( kt, 1, sf_sbias_ofl )
788
789	zcof1(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_t_ofl + &
790	& fbcoef(:,:) * fb_t_ofl_max
791	zcof2(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_p_ofl
792	DO jk = 1, jpkm1
793	tbias(:,:,jk) = tbias(:,:,jk) + sf_tbias_ofl(1)%fnow(:,:,jk) * zcof1(:,:)
794	sbias(:,:,jk) = sbias(:,:,jk) + sf_sbias_ofl(1)%fnow(:,:,jk) * zcof1(:,:)
795	tbias_p(:,:,jk) = tbias_p(:,:,jk) + sf_tbias_ofl(1)%fnow(:,:,jk) * zcof2(:,:)
796	sbias_p(:,:,jk) = sbias_p(:,:,jk) + sf_sbias_ofl(1)%fnow(:,:,jk) * zcof2(:,:)
797	ENDDO
798	ENDIF
799
800
801	!apply bias on tracers if needed
802	IF ( ln_bias_ts_app ) THEN
803
804	! Add the bias directely to T/s
805	tsa(:,:,:,jp_tem) = tsa(:,:,:,jp_tem) + tmask(:,:,:) * tbias(:,:,:) / rdt
806	tsa(:,:,:,jp_sal) = tsa(:,:,:,jp_sal) + tmask(:,:,:) * sbias(:,:,:) / rdt
807
808	! lateral boundary conditions (is this needed?)
809	CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1.0_wp )
810	CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1.0_wp )
811
812	ENDIF
813
814	IF ( kt == nn_bias_itwrt ) THEN
815	tbias_p_out(:,:,:) = tbias_p(:,:,:)
816	sbias_p_out(:,:,:) = sbias_p(:,:,:)
817	ENDIF
818
819	END SUBROUTINE tra_bias
820
821	SUBROUTINE dyn_bias( kt )
822	!!----------------------------------------------------------------------
823	!! * ROUTINE dyn_bias *
824	!!
825	!! ** Purpose : Computes rhd_pc, gru/v_pc bias corrected
826	!! for hydrostatic pressure gradient
827	!! depending on time step (semi-implicit/centered)
828	!! If partial steps computes bottom pressure gradient.
829	!! These correction terms will affect only the dynamical
830	!! component (pressure gradient calculation), but not
831	!! the isopycnal calculation for the lateral mixing.
832	!!
833	!! ** Method : At this stage of the computation, ta and sa are the
834	!! after temperature and salinity. If semi-implicit, these
835	!! are used to compute rho and bottom pressure gradient.
836	!! If centered, tb,sb are used instead.
837	!! If bias key is activated, the temperature,salinity are
838	!! bias corrected in the calculation of the density fields
839	!! used in the pressure gradient calculation.
840	!!
841	!!
842	!! ** Action : - rhd_pc ready. rhop will be overwriten later
843	!! - if ln_zps, bottom density gradients gru/v_pc ready.
844	!!----------------------------------------------------------------------
845	!!
846	!! * Arguments
847	INTEGER, INTENT(IN) :: kt ! ocean time-step index
848	!! * Local variables
849	REAL(wp) :: tsw(jpi,jpj,jpk,jpts)
850	!!
851	!!----------------------------------------------------------------------
852	!
853	! gtu,gsu,gtv,gsv rhop will be overwritten later in step.
854	!
855	IF( ln_dynhpg_imp ) THEN ! semi-implicit hpg
856	tsw(:,:,:,jp_tem) = tsa(:,:,:,jp_tem) - tbias_p(:,:,:)
857	tsw(:,:,:,jp_sal) = tsa(:,:,:,jp_sal) - sbias_p(:,:,:)
858	ELSE
859	tsw(:,:,:,jp_tem) = tsb(:,:,:,jp_tem) - tbias_p(:,:,:)
860	tsw(:,:,:,jp_sal) = tsb(:,:,:,jp_sal) - sbias_p(:,:,:)
861	ENDIF
862
863	CALL eos( tsw, rhd_pc, rhop )
864
865	! is this needed?
866	CALL lbc_lnk( rhd_pc, 'T', 1.0_wp )
867
868	! Partial steps: now horizontal gradient of t,s,rd
869	! at the bottom ocean level
870	IF( ln_zps ) THEN
871	CALL zps_hde( kt, jpts, tsw, gtsu, gtsv, &
872	& rhd_pc, gru_pc , grv_pc )
873	ENDIF
874
875	END SUBROUTINE dyn_bias
876
877	SUBROUTINE bias_opn( kt )
878	!!---------------------------------------------------------------------
879	!! * ROUTINE bias_opn *
880	!!
881	!! ** Purpose : open bias restart file following the same logic as the
882	!! standard restarts.
883	!!----------------------------------------------------------------------
884	!! * Arguments
885	INTEGER, INTENT(IN) :: kt ! ocean time-step
886	!! * Local variables
887	CHARACTER(LEN=20) :: clbkt ! ocean time-step deine as a character
888	CHARACTER(LEN=50) :: clbias_tot ! total bias restart file name
889	!!----------------------------------------------------------------------
890	!
891	IF( lrst_oce .AND. .NOT.lrst_bias ) THEN ! create bias file
892	IF( nitend > 999999999 ) THEN ; WRITE(clbkt, * ) nitend
893	ELSE ; WRITE(clbkt, '(i8.8)') nitend
894	ENDIF
895	clbias_tot = TRIM(cexper)//"_"//TRIM(ADJUSTL(clbkt))//"_"//TRIM(cn_bias_tot)
896	IF(lwp) THEN
897	WRITE(numout,*)
898	SELECT CASE ( jprstlib )
899	CASE ( jprstdimg ) ; WRITE(numout,*) ' open tot bias restart binary file: '//clbias_tot
900	CASE DEFAULT ; WRITE(numout,*) ' open tot bias restart NetCDF file: '//clbias_tot
901	END SELECT
902	ENDIF
903	CALL iom_open( clbias_tot, numbias_tot , ldwrt = .TRUE., kiolib = jprstlib )
904	lrst_bias=.TRUE.
905
906	ENDIF
907	!
908	END SUBROUTINE bias_opn
909
910	SUBROUTINE bias_wrt( kt )
911	!!---------------------------------------------------------------------
912	!! * ROUTINE bias_wrt *
913	!!
914	!! ** Purpose : Write bias restart fields in the format corresponding to jprstlib
915	!!
916	!! ** Method : Write in numbias_tot when kt == nitend in output
917	!! file, save fields which are necessary for restart.
918	!!
919	!! Changed the timestep for writing to nitend rather than nitrst as this causes a
920	!! problem if the nstock list is used to determine the restart output times and
921	!! means that the bias is not output at all. M. Martin. 08/2011.
922	!! Need to check with Magdalena that this is ok for ECMWF.
923	!!----------------------------------------------------------------------
924	!! * Arguments
925	INTEGER, INTENT(IN) :: kt ! ocean time-step
926	!!----------------------------------------------------------------------
927	! ! the begining of the run [s]
928
929	IF ( ln_bias_rlx ) THEN
930	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'tbias_rlx' , tbias_rlx_out )
931	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'sbias_rlx' , sbias_rlx_out )
932	ENDIF
933
934	IF ( ln_bias_asm ) THEN
935	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'tbias_asm' , tbias_asm_out )
936	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'sbias_asm' , sbias_asm_out )
937	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'tbias_p' , tbias_p_out )
938	CALL iom_rstput( nn_bias_itwrt, nn_bias_itwrt, numbias_tot, 'sbias_p' , sbias_p_out )
939	ENDIF
940
941	IF( kt == nitend ) THEN
942	CALL iom_close( numbias_tot ) ! close the restart file (only at last time step)
943	lrst_bias = .FALSE.
944	ENDIF
945	!
946	END SUBROUTINE bias_wrt
947
948	END MODULE bias

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

Download in other formats: