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.

diaobs.F90 in branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/OBS – NEMO

Context Navigation

source: branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90 @ 5659

Last change on this file since 5659 was 5659, checked in by mattmartin, 9 years ago
Updated OBS simplification branch to the head of the trunk.
Property svn:keywords set to `Id`
File size: 33.3 KB

Line
1	MODULE diaobs
2	!!======================================================================
3	!! * MODULE diaobs *
4	!! Observation diagnostics: Computation of the misfit between data and
5	!! their model equivalent
6	!!======================================================================
7
8	!!----------------------------------------------------------------------
9	!! 'key_diaobs' : Switch on the observation diagnostic computation
10	!!----------------------------------------------------------------------
11	!! dia_obs_init : Reading and prepare observations
12	!! dia_obs : Compute model equivalent to observations
13	!! dia_obs_wri : Write observational diagnostics
14	!! ini_date : Compute the initial date YYYYMMDD.HHMMSS
15	!! fin_date : Compute the final date YYYYMMDD.HHMMSS
16	!!----------------------------------------------------------------------
17	!! * Modules used
18	USE wrk_nemo ! Memory Allocation
19	USE par_kind ! Precision variables
20	USE in_out_manager ! I/O manager
21	USE par_oce
22	USE dom_oce ! Ocean space and time domain variables
23	USE obs_fbm, ONLY: ln_cl4 ! Class 4 diagnostic switch
24	USE obs_read_prof ! Reading and allocation of observations (Coriolis)
25	USE obs_read_surf ! Reading and allocation of SLA observations
26	USE obs_readmdt ! Reading and allocation of MDT for SLA.
27	USE obs_prep ! Preparation of obs. (grid search etc).
28	USE obs_oper ! Observation operators
29	USE obs_write ! Writing of observation related diagnostics
30	USE obs_grid ! Grid searching
31	USE obs_read_altbias ! Bias treatment for altimeter
32	USE obs_profiles_def ! Profile data definitions
33	USE obs_surf_def ! Surface data definitions
34	USE obs_types ! Definitions for observation types
35	USE mpp_map ! MPP mapping
36	USE lib_mpp ! For ctl_warn/stop
37
38	IMPLICIT NONE
39
40	!! * Routine accessibility
41	PRIVATE
42	PUBLIC dia_obs_init, & ! Initialize and read observations
43	& dia_obs, & ! Compute model equivalent to observations
44	& dia_obs_wri, & ! Write model equivalent to observations
45	& dia_obs_dealloc ! Deallocate dia_obs data
46
47	!! * Shared Module variables
48	LOGICAL, PUBLIC, PARAMETER :: &
49	#if defined key_diaobs
50	& lk_diaobs = .TRUE. !: Logical switch for observation diangostics
51	#else
52	& lk_diaobs = .FALSE. !: Logical switch for observation diangostics
53	#endif
54
55	!! * Module variables
56	LOGICAL, PUBLIC :: ln_t3d !: Logical switch for temperature profiles
57	LOGICAL, PUBLIC :: ln_s3d !: Logical switch for salinity profiles
58	LOGICAL, PUBLIC :: ln_sla !: Logical switch for sea level anomalies
59	LOGICAL, PUBLIC :: ln_sst !: Logical switch for sea surface temperature
60	LOGICAL, PUBLIC :: ln_seaice !: Logical switch for sea ice concentration
61	LOGICAL, PUBLIC :: ln_vel3d !: Logical switch for velocity component (u,v) observations
62	LOGICAL, PUBLIC :: ln_ssh !: Logical switch for sea surface height
63	LOGICAL, PUBLIC :: ln_sss !: Logical switch for sea surface salinity
64	LOGICAL, PUBLIC :: ln_sstnight !: Logical switch for night mean SST observations
65	LOGICAL, PUBLIC :: ln_nea !: Remove observations near land
66	LOGICAL, PUBLIC :: ln_altbias !: Logical switch for altimeter bias
67	LOGICAL, PUBLIC :: ln_ignmis !: Logical switch for ignoring missing files
68	LOGICAL, PUBLIC :: ln_s_at_t !: Logical switch to compute model S at T observations
69
70	REAL(KIND=dp), PUBLIC :: dobsini !: Observation window start date YYYYMMDD.HHMMSS
71	REAL(KIND=dp), PUBLIC :: dobsend !: Observation window end date YYYYMMDD.HHMMSS
72
73	INTEGER, PUBLIC :: numobtypes !: Number of observation types to read in.
74	INTEGER, PUBLIC :: n1dint !: Vertical interpolation method
75	INTEGER, PUBLIC :: n2dint !: Horizontal interpolation method
76	INTEGER, DIMENSION(:), ALLOCATABLE :: nvarsprof !Number of profile variables
77	INTEGER, DIMENSION(:), ALLOCATABLE :: nextrprof !Number of profile extra variables
78	INTEGER, DIMENSION(:), ALLOCATABLE :: nvarssurf !Number of surface variables
79	INTEGER, DIMENSION(:), ALLOCATABLE :: nextrsurf !Number of surface extra variables
80	INTEGER, DIMENSION(imaxavtypes) :: &
81	& dailyavtypes !: Data types which are daily average
82
83	TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) :: surfdata ! Initial surface data
84	TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) :: surfdataqc ! Surface data after quality control
85	TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) :: profdata ! Initial profile data
86	TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) :: profdataqc ! Profile data after quality control
87
88	CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE :: obstypesprof
89	CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE :: obstypessurf
90
91
92
93	INTEGER, PARAMETER :: MaxNumFiles = 1000
94
95	LOGICAL, DIMENSION(MaxNumFiles) :: &
96	& ln_profb_ena, & !: Is the feedback files from ENACT data ?
97	! !: If so use dailyavtypes
98	& ln_profb_enatim !: Change tim for 820 enact data set.
99
100	LOGICAL, DIMENSION(MaxNumFiles) :: &
101	& ln_velfb_av !: Is the velocity feedback files daily average?
102	LOGICAL, DIMENSION(:), ALLOCATABLE :: &
103	& ld_enact !: Profile data is ENACT so use dailyavtypes
104	LOGICAL, DIMENSION(:), ALLOCATABLE :: &
105	& ld_velav !: Velocity data is daily averaged
106	LOGICAL, DIMENSION(:), ALLOCATABLE :: &
107	& ld_sstnight !: SST observation corresponds to night mean
108
109	!!----------------------------------------------------------------------
110	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
111	!! $Id$
112	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
113	!!----------------------------------------------------------------------
114
115	CONTAINS
116
117	SUBROUTINE dia_obs_init
118	!!----------------------------------------------------------------------
119	!! * ROUTINE dia_obs_init *
120	!!
121	!! ** Purpose : Initialize and read observations
122	!!
123	!! ** Method : Read the namelist and call reading routines
124	!!
125	!! ** Action : Read the namelist and call reading routines
126	!!
127	!! History :
128	!! ! 06-03 (K. Mogensen) Original code
129	!! ! 06-05 (A. Weaver) Reformatted
130	!! ! 06-10 (A. Weaver) Cleaning and add controls
131	!! ! 07-03 (K. Mogensen) General handling of profiles
132	!! ! 15-02 (M. Martin) Simplification of namelist and code
133	!!----------------------------------------------------------------------
134
135	IMPLICIT NONE
136
137	!! * Local declarations
138	CHARACTER(len=128) :: profbfiles(MaxNumFiles)
139	CHARACTER(len=128) :: sstfbfiles(MaxNumFiles)
140	CHARACTER(len=128) :: slafbfiles(MaxNumFiles)
141	CHARACTER(len=128) :: seaicefbfiles(MaxNumFiles)
142	CHARACTER(len=128) :: velfbfiles(MaxNumFiles)
143	CHARACTER(LEN=128) :: bias_file
144	CHARACTER(LEN=20) :: datestr=" ", timestr=" "
145
146	NAMELIST/namobs/ln_t3d, ln_s3d, ln_sla, ln_sss, ln_ssh, &
147	& ln_sst, ln_seaice, ln_vel3d, &
148	& ln_altbias, ln_nea, ln_grid_global, &
149	& ln_grid_search_lookup, ln_cl4, &
150	& ln_ignmis, ln_s_at_t, ln_sstnight, &
151	& ln_profb_ena, ln_profb_enatim, &
152	& profbfiles, slafbfiles, sssfbfiles, &
153	& sshfbfiles, sstfbfiles, seaicefbfiles, &
154	& velfbfiles, bias_file, grid_search_file, &
155	& dobsini, dobsend, n1dint, n2dint, &
156	& nmsshc, mdtcorr, mdtcutoff, &
157	& grid_search_res, dailyavtypes
158
159	INTEGER :: jtype
160	INTEGER :: ios ! Local integer output status for namelist read
161	INTEGER, DIMENSION(:), ALLOCATABLE :: jnumfilesprof
162	INTEGER, DIMENSION(:), ALLOCATABLE :: jnumfilessurf
163	CHARACTER(len=128), DIMENSION(:,:), ALLOCATABLE :: obsfilesprof
164	CHARACTER(len=128), DIMENSION(:,:), ALLOCATABLE :: obsfilessurf
165	LOGICAL :: lmask(MaxNumFiles)
166	!-----------------------------------------------------------------------
167	! Read namelist parameters
168	!-----------------------------------------------------------------------
169
170	profbfiles(:) = ''
171	slafbfiles(:) = ''
172	sstfbfiles(:) = ''
173	seaicefbfiles(:) = ''
174	velfbfiles(:) = ''
175	dailyavtypes(:) = -1
176	dailyavtypes(1) = 820
177	ln_profb_ena(:) = .FALSE.
178	ln_profb_enatim(:) = .TRUE.
179	ln_velfb_av(:) = .FALSE.
180	ln_ignmis = .FALSE.
181
182	CALL ini_date( dobsini )
183	CALL fin_date( dobsend )
184
185	! Read Namelist namobs : control observation diagnostics
186	REWIND( numnam_ref ) ! Namelist namobs in reference namelist : Diagnostic: control observation
187	READ ( numnam_ref, namobs, IOSTAT = ios, ERR = 901)
188	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namobs in reference namelist', lwp )
189
190	REWIND( numnam_cfg ) ! Namelist namobs in configuration namelist : Diagnostic: control observation
191	READ ( numnam_cfg, namobs, IOSTAT = ios, ERR = 902 )
192	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namobs in configuration namelist', lwp )
193	IF(lwm) WRITE ( numond, namobs )
194
195	!Set up list of observation types to be used
196	numproftypes = COUNT( (/ln_t3d .OR. ln_s3d, ln_vel3d /) )
197	numsurftypes = COUNT( (/ln_sla, ln_sss, ln_sst, ln_seaice /) )
198	IF ( numproftypes > 0 ) THEN
199
200	ALLOCATE( obstypesprof(numproftypes) )
201	ALLOCATE( jnumfilesprof(numproftypes) )
202	ALLOCATE( obsfilesprof(numproftypes, MaxNumFiles) )
203
204	DO jtype = 1, numproftypes
205	IF (ln_t3d .OR. ln_s3d) THEN
206	obsfilesprof(:,jtype) = profbfiles(:)
207	obstypesprof(jtype) = 'prof'
208	ENDIF
209	IF (ln_vel3d) THEN
210	obsfilesprof(:,jtype) = velfbfiles(:)
211	obstypesprof(jtype) = 'vel'
212	ENDIF
213
214	lmask(:) = .FALSE.
215	WHERE (obsfilesprof(jtype,:) /= '') lmask(:) = .TRUE.
216	jnumfilesprof(jtype) = COUNT(lmask)
217	END DO
218
219	ENDIF
220
221	IF ( numsurftypes > 0 ) THEN
222
223	ALLOCATE( obstypessurf(numsurftypes) )
224	ALLOCATE( jnumfilessurf(numproftypes) )
225	ALLOCATE( obsfilessurf(numsurftypes, MaxNumFiles) )
226
227	DO jtype = 1, numsurftypes
228	IF (ln_sla) THEN
229	obsfilessurf(:,jtype) = slafbfiles(:)
230	obstypessurf(jtype) = 'sla'
231	ENDIF
232	IF (ln_sss) THEN
233	obsfilessurf(:,jtype) = sssfbfiles(:)
234	obstypessurf(jtype) = 'sss'
235	ENDIF
236	IF (ln_sst) THEN
237	obsfilessurf(:,jtype) = sstfbfiles(:)
238	obstypessurf(jtype) = 'sst'
239	ENDIF
240	#if defined key_lim2 \|\| defined key_lim3
241	IF (ln_seaice) THEN
242	obsfilessurf(:,jtype) = seaicefbfiles(:)
243	obstypessurf(jtype) = 'seaice'
244	ENDIF
245	#endif
246
247	lmask(:) = .FALSE.
248	WHERE (obsfilessurf(jtype,:) /= '') lmask(:) = .TRUE.
249	jnumfilessurf(jtype) = COUNT(lmask)
250
251	END DO
252
253	ENDIF
254
255	!Write namelist settings to stdout
256	IF(lwp) THEN
257	WRITE(numout,*)
258	WRITE(numout,*) 'dia_obs_init : Observation diagnostic initialization'
259	WRITE(numout,*) '~~~~~~~~~~~~'
260	WRITE(numout,*) ' Namelist namobs : set observation diagnostic parameters'
261	WRITE(numout,*) ' Logical switch for T profile observations ln_t3d = ', ln_t3d
262	WRITE(numout,*) ' Logical switch for S profile observations ln_s3d = ', ln_s3d
263	WRITE(numout,*) ' Logical switch for SLA observations ln_sla = ', ln_sla
264	WRITE(numout,*) ' Logical switch for SSH observations ln_ssh = ', ln_ssh
265	WRITE(numout,*) ' Logical switch for SST observations ln_sst = ', ln_sst
266	WRITE(numout,*) ' Logical switch for night-time SST obs ln_sstnight = ', ln_sstnight
267	WRITE(numout,*) ' Logical switch for SSS observations ln_sss = ', ln_sss
268	WRITE(numout,*) ' Logical switch for Sea Ice observations ln_seaice = ', ln_seaice
269	WRITE(numout,*) ' Logical switch for velocity observations ln_vel3d = ', ln_vel3d
270	WRITE(numout,*) ' Global distribution of observations ln_grid_global = ',ln_grid_global
271	WRITE(numout,*) &
272	' Logical switch for obs grid search w/lookup table ln_grid_search_lookup = ',ln_grid_search_lookup
273	IF (ln_grid_search_lookup) &
274	WRITE(numout,*) ' Grid search lookup file header grid_search_file = ', grid_search_file
275	WRITE(numout,*) ' Initial date in window YYYYMMDD.HHMMSS dobsini = ', dobsin
276	WRITE(numout,*) ' Final date in window YYYYMMDD.HHMMSS dobsend = ', dobsend
277	WRITE(numout,*) ' Type of vertical interpolation method n1dint = ', n1dint
278	WRITE(numout,*) ' Type of horizontal interpolation method n2dint = ', n2dint
279	WRITE(numout,*) ' Rejection of observations near land swithch ln_nea = ', ln_nea
280	WRITE(numout,*) ' MSSH correction scheme nmsshc = ', nmsshc
281	WRITE(numout,*) ' MDT correction mdtcorr = ', mdtcorr
282	WRITE(numout,*) ' MDT cutoff for computed correction mdtcutoff = ', mdtcutoff
283	WRITE(numout,*) ' Logical switch for alt bias ln_altbias = ', ln_altbias
284	WRITE(numout,*) ' Logical switch for ignoring missing files ln_ignmis = ', ln_ignmis
285	WRITE(numout,*) ' Daily average types = ', dailyavtypes
286
287	IF ( numproftypes > 0 ) THEN
288	DO jtype = 1, numproftypes
289	DO ji = 1, jnumfilesprof(jtype)
290	WRITE(numout,'(1X,2A)') ' '//obstypesprof(jtype)//' input observation file names = ', &
291	TRIM(obsfilesprof(jtype,ji))
292	IF ( TRIM(obstypesprof(jtype)) == 'prof' ) &
293	WRITE(numout,'(1X,2A)') ' Enact feedback input time setting switch ln_profb_enatim = ', ln_profb_enatim(ji)
294	END DO
295	END DO
296	ENDIF
297
298	IF ( numsurftypes > 0 ) THEN
299	DO jtype = 1, numsurftypes
300	DO ji = 1, jnumfilessurf(jtype)
301	WRITE(numout,'(1X,2A)') ' '//obstypessurf(jtype)//' input observation file names = ', &
302	TRIM(obsfilessurf(jtype,ji))
303	END DO
304	END DO
305	ENDIF
306
307	ENDIF
308
309	IF ( ln_vel3d .AND. ( .NOT. ln_grid_global ) ) THEN
310	CALL ctl_stop( 'Velocity data only works with ln_grid_global=.true.' )
311	RETURN
312	ENDIF
313
314	CALL obs_typ_init
315
316	CALL mppmap_init
317
318	! Parameter control
319	#if defined key_diaobs
320	IF ( numobtypes == 0 ) THEN
321	IF(lwp) WRITE(numout,cform_war)
322	IF(lwp) WRITE(numout,*) ' key_diaobs is activated but logical flags', &
323	& ' ln_t3d, ln_s3d, ln_sla, ln_ssh, ln_sst, ln_sss, ln_seaice, ln_vel3d are all set to .FALSE.'
324	nwarn = nwarn + 1
325	ENDIF
326	#endif
327
328	CALL obs_grid_setup( )
329	IF ( ( n1dint < 0 ).OR.( n1dint > 1 ) ) THEN
330	CALL ctl_stop(' Choice of vertical (1D) interpolation method', &
331	& ' is not available')
332	ENDIF
333	IF ( ( n2dint < 0 ).OR.( n2dint > 4 ) ) THEN
334	CALL ctl_stop(' Choice of horizontal (2D) interpolation method', &
335	& ' is not available')
336	ENDIF
337
338	!-----------------------------------------------------------------------
339	! Depending on switches read the various observation types
340	!-----------------------------------------------------------------------
341
342	IF ( numproftypes > 0 ) THEN
343
344	ALLOCATE(profdata(numproftypes))
345	ALLOCATE(profdataqc(numproftypes))
346	ALLOCATE(nvarsprof(numproftypes))
347	ALLOCATE(nextrprof(numproftypes))
348
349	DO jtype = 1, numproftypes
350
351	nvarsprof(jtype) = 2
352	IF ( TRIM(obstypesprof(jtype)) == 'prof' ) nextrprof(jtype) = 1
353	IF ( TRIM(obstypesprof(jtype)) == 'vel' ) nextrprof(jtype) = 2
354
355	!Read in profile or velocity obs types
356	CALL obs_rea_prof( profdata(jtype), &
357	& jnumfilesprof(jtype), &
358	& obsfilesprof(jtype,1:jnumfilesprof(jtype)), &
359	& nvarsprof(jtype), nextrprof(jtype), nitend-nit000+2, &
360	& dobsini, dobsend, ln_t3d, ln_s3d, &
361	& ln_ignmis, ln_s_at_t, .TRUE., .FALSE., &
362	& kdailyavtypes = dailyavtypes )
363
364	DO jvar = 1, nvars
365	CALL obs_prof_staend( profdata(jtype), jvar )
366	END DO
367
368	CALL obs_pre_prof( profdata(jtype), profdataqc(jtype), &
369	& ln_t3d, ln_s3d, ln_nea, &
370	& kdailyavtypes = dailyavtypes )
371
372	END DO
373
374	DEALLOCATE( jnumfilesprof, obsfilesprof )
375
376	ENDIF
377
378	IF ( numsurftypes > 0 ) THEN
379
380	ALLOCATE(surfdata(numsurftypes))
381	ALLOCATE(surfdatatqc(numsurftypes))
382	ALLOCATE(nvarssurf(numsurftypes))
383	ALLOCATE(nextrsurf(numsurftypes))
384
385	DO jtype = 1, numsurftypes
386
387	nvarssurf(jtype) = 1
388	nextrsurf(jtype) = 0
389	IF ( TRIM(obstypessurf(jtype)) == 'sla' ) nextrsurf(jtype) = 2
390
391	!Read in surface obs types
392	CALL obs_rea_surf( surfdata(jtype), jnumfilessurf(jtype), &
393	& obsfilessurf(jtype,1:jnumfilessurf(jtype)), &
394	& nvarssurf(jtype), nextrsurf(jtype), nitend-nit000+2, &
395	& dobsini, dobsend, ln_ignmis, .FALSE. )
396
397	CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), ln_nea )
398
399	IF ( TRIM(obstypessurf(jtype)) == 'sla' ) THEN
400	CALL obs_rea_mdt( surfdataqc(jtype), n2dint )
401	IF ( ln_altbias ) CALL obs_rea_altbias ( surfdataqc(jtype), n2dint, bias_file )
402	ENDIF
403
404	DEALLOCATE( jnumfilessurf, obsfilessurf )
405
406	END DO
407
408	END SUBROUTINE dia_obs_init
409
410	SUBROUTINE dia_obs( kstp )
411	!!----------------------------------------------------------------------
412	!! * ROUTINE dia_obs *
413	!!
414	!! ** Purpose : Call the observation operators on each time step
415	!!
416	!! ** Method : Call the observation operators on each time step to
417	!! compute the model equivalent of the following date:
418	!! - T profiles
419	!! - S profiles
420	!! - Sea surface height (referenced to a mean)
421	!! - Sea surface temperature
422	!! - Sea surface salinity
423	!! - Velocity component (U,V) profiles
424	!!
425	!! ** Action :
426	!!
427	!! History :
428	!! ! 06-03 (K. Mogensen) Original code
429	!! ! 06-05 (K. Mogensen) Reformatted
430	!! ! 06-10 (A. Weaver) Cleaning
431	!! ! 07-03 (K. Mogensen) General handling of profiles
432	!! ! 07-04 (G. Smith) Generalized surface operators
433	!! ! 08-10 (M. Valdivieso) obs operator for velocity profiles
434	!!----------------------------------------------------------------------
435	!! * Modules used
436	USE dom_oce, ONLY : & ! Ocean space and time domain variables
437	& rdt, &
438	& gdept_1d, &
439	& tmask, umask, vmask
440	USE phycst, ONLY : & ! Physical constants
441	& rday
442	USE oce, ONLY : & ! Ocean dynamics and tracers variables
443	& tsn, &
444	& un, vn, &
445	& sshn
446	#if defined key_lim3
447	USE ice, ONLY : & ! LIM Ice model variables
448	& frld
449	#endif
450	#if defined key_lim2
451	USE ice_2, ONLY : & ! LIM Ice model variables
452	& frld
453	#endif
454	IMPLICIT NONE
455
456	!! * Arguments
457	INTEGER, INTENT(IN) :: kstp ! Current timestep
458	!! * Local declarations
459	INTEGER :: idaystp ! Number of timesteps per day
460	INTEGER :: jtype ! data loop variable
461	INTEGER :: jvar ! Variable number
462	#if ! defined key_lim2 && ! defined key_lim3
463	REAL(wp), POINTER, DIMENSION(:,:) :: frld
464	#endif
465	CHARACTER(LEN=20) :: datestr=" ",timestr=" "
466
467	#if ! defined key_lim2 && ! defined key_lim3
468	CALL wrk_alloc(jpi,jpj,frld)
469	#endif
470
471	IF(lwp) THEN
472	WRITE(numout,*)
473	WRITE(numout,*) 'dia_obs : Call the observation operators', kstp
474	WRITE(numout,*) '~~~~~~~'
475	ENDIF
476
477	idaystp = NINT( rday / rdt )
478
479	!-----------------------------------------------------------------------
480	! No LIM => frld == 0.0_wp
481	!-----------------------------------------------------------------------
482	#if ! defined key_lim2 && ! defined key_lim3
483	frld(:,:) = 0.0_wp
484	#endif
485	!-----------------------------------------------------------------------
486	! Depending on switches call various observation operators
487	!-----------------------------------------------------------------------
488
489	IF ( numproftypes > 0 ) THEN
490	DO jtype = 1, numproftypes
491
492	SELECT CASE ( TRIM(obstypesprof(jtype)) )
493	CASE('prof')
494	profvar1(:,:,:) = tsn(:,:,:,jp_tem)
495	profvar2(:,:,:) = tsn(:,:,:,jp_sal)
496	profmask1(:,:,:) = tmask(:,:,:)
497	profmask2(:,:,:) = tmask(:,:,:)
498	CASE('vel')
499	profvar1(:,:,:) = un(:,:,:)
500	profvar2(:,:,:) = vn(:,:,:)
501	profmask1(:,:,:) = umask(:,:,:)
502	profmask2(:,:,:) = vmask(:,:,:)
503	END SELECT
504
505	CALL obs_prof_opt( profdataqc(jtype), &
506	& kstp, jpi, jpj, jpk, nit000, idaystp, &
507	& profvar1, profvar2, &
508	& gdept_1d, profmask1, profmask2, n1dint, n2dint, &
509	& kdailyavtypes = dailyavtypes )
510
511	END DO
512
513	ENDIF
514
515	IF ( numsurftypes > 0 ) THEN
516	DO jtype = 1, numsurftypes
517
518	SELECT CASE ( TRIM(obstypessurf(jtype)) )
519	CASE('sst')
520	surfvar(:,:) = tsn(:,:,1,jp_tem)
521	CASE('sla')
522	surfvar(:,:) = sshn(:,:)
523	CASE('sss')
524	surfvar(:,:) = tsn(:,:,1,jp_sal)
525	#if defined key_lim2 \|\| defined key_lim3
526	CASE('seaice')
527	surfvar(:,:) = 1._wp - frld(:,:)
528	#endif
529	END SELECT
530
531	CALL obs_surf_opt( surfdatqc(jtype), &
532	& kstp, jpi, jpj, nit000, surfvar, &
533	& tmask(:,:,1), n2dint, ld_sstnight )
534
535	END DO
536
537	ENDIF
538
539	#if ! defined key_lim2 && ! defined key_lim3
540	CALL wrk_dealloc(jpi,jpj,frld)
541	#endif
542
543	END SUBROUTINE dia_obs
544
545	SUBROUTINE dia_obs_wri
546	!!----------------------------------------------------------------------
547	!! * ROUTINE dia_obs_wri *
548	!!
549	!! ** Purpose : Call observation diagnostic output routines
550	!!
551	!! ** Method : Call observation diagnostic output routines
552	!!
553	!! ** Action :
554	!!
555	!! History :
556	!! ! 06-03 (K. Mogensen) Original code
557	!! ! 06-05 (K. Mogensen) Reformatted
558	!! ! 06-10 (A. Weaver) Cleaning
559	!! ! 07-03 (K. Mogensen) General handling of profiles
560	!! ! 08-09 (M. Valdivieso) Velocity component (U,V) profiles
561	!!----------------------------------------------------------------------
562	IMPLICIT NONE
563
564	!! * Local declarations
565
566	INTEGER :: jtype ! Data set loop variable
567	!-----------------------------------------------------------------------
568	! Depending on switches call various observation output routines
569	!-----------------------------------------------------------------------
570
571	IF ( numproftypes > 0 ) THEN
572	DO jtype = 1, numproftypes
573
574	CALL obs_prof_decompress( profdataqc(jtype), &
575	& profdata(jtype), .TRUE., numout )
576
577	CALL obs_wri_prof( obstypesprof(jtype), profdata(jtype), n2dint )
578
579	END DO
580
581	ENDIF
582
583	IF ( numsurftypes > 0 ) THEN
584	DO jtype = 1, numsurftypes
585
586	CALL obs_surf_decompress( surfdatqc(jtype), &
587	& surfdata(jtype), .TRUE., numout )
588
589	CALL obs_wri_surf( obstypessurf(jtype), surfdata(jtype), n2dint )
590
591	END DO
592
593	ENDIF
594
595
596	END SUBROUTINE dia_obs_wri
597
598	SUBROUTINE dia_obs_dealloc
599	IMPLICIT NONE
600	!!----------------------------------------------------------------------
601	!! * ROUTINE dia_obs_dealloc *
602	!!
603	!! ** Purpose : To deallocate data to enable the obs_oper online loop.
604	!! Specifically: dia_obs_init --> dia_obs --> dia_obs_wri
605	!!
606	!! ** Method : Clean up various arrays left behind by the obs_oper.
607	!!
608	!! ** Action :
609	!!
610	!!----------------------------------------------------------------------
611	!! obs_grid deallocation
612	CALL obs_grid_deallocate
613
614	!! diaobs deallocation
615	IF ( numproftypes > 0 ) DEALLOCATE(profdata, profdataqc, nvarsprof, nextrprof)
616	IF ( numsurftypes > 0 ) DEALLOCATE(surfdata, surfdataqc, nvarssurf, nextrsurf)
617
618	END SUBROUTINE dia_obs_dealloc
619
620	SUBROUTINE ini_date( ddobsini )
621	!!----------------------------------------------------------------------
622	!! * ROUTINE ini_date *
623	!!
624	!! ** Purpose : Get initial data in double precision YYYYMMDD.HHMMSS format
625	!!
626	!! ** Method : Get initial data in double precision YYYYMMDD.HHMMSS format
627	!!
628	!! ** Action : Get initial data in double precision YYYYMMDD.HHMMSS format
629	!!
630	!! History :
631	!! ! 06-03 (K. Mogensen) Original code
632	!! ! 06-05 (K. Mogensen) Reformatted
633	!! ! 06-10 (A. Weaver) Cleaning
634	!! ! 06-10 (G. Smith) Calculates initial date the same as method for final date
635	!! ! 10-05 (D. Lea) Update to month length calculation for NEMO vn3.2
636	!!----------------------------------------------------------------------
637	USE phycst, ONLY : & ! Physical constants
638	& rday
639	! USE daymod, ONLY : & ! Time variables
640	! & nmonth_len
641	USE dom_oce, ONLY : & ! Ocean space and time domain variables
642	& rdt
643
644	IMPLICIT NONE
645
646	!! * Arguments
647	REAL(KIND=dp), INTENT(OUT) :: ddobsini ! Initial date in YYYYMMDD.HHMMSS
648
649	!! * Local declarations
650	INTEGER :: iyea ! date - (year, month, day, hour, minute)
651	INTEGER :: imon
652	INTEGER :: iday
653	INTEGER :: ihou
654	INTEGER :: imin
655	INTEGER :: imday ! Number of days in month.
656	REAL(KIND=wp) :: zdayfrc ! Fraction of day
657
658	INTEGER, DIMENSION(12) :: imonth_len !: length in days of the months of the current year
659
660	!!----------------------------------------------------------------------
661	!! Initial date initialization (year, month, day, hour, minute)
662	!! (This assumes that the initial date is for 00z))
663	!!----------------------------------------------------------------------
664	iyea = ndate0 / 10000
665	imon = ( ndate0 - iyea * 10000 ) / 100
666	iday = ndate0 - iyea * 10000 - imon * 100
667	ihou = 0
668	imin = 0
669
670	!!----------------------------------------------------------------------
671	!! Compute number of days + number of hours + min since initial time
672	!!----------------------------------------------------------------------
673	iday = iday + ( nit000 -1 ) * rdt / rday
674	zdayfrc = ( nit000 -1 ) * rdt / rday
675	zdayfrc = zdayfrc - aint(zdayfrc)
676	ihou = int( zdayfrc * 24 )
677	imin = int( (zdayfrc * 24 - ihou) * 60 )
678
679	!!-----------------------------------------------------------------------
680	!! Convert number of days (iday) into a real date
681	!!----------------------------------------------------------------------
682
683	CALL calc_month_len( iyea, imonth_len )
684
685	DO WHILE ( iday > imonth_len(imon) )
686	iday = iday - imonth_len(imon)
687	imon = imon + 1
688	IF ( imon > 12 ) THEN
689	imon = 1
690	iyea = iyea + 1
691	CALL calc_month_len( iyea, imonth_len ) ! update month lengths
692	ENDIF
693	END DO
694
695	!!----------------------------------------------------------------------
696	!! Convert it into YYYYMMDD.HHMMSS format.
697	!!----------------------------------------------------------------------
698	ddobsini = iyea * 10000_dp + imon * 100_dp + &
699	& iday + ihou * 0.01_dp + imin * 0.0001_dp
700
701
702	END SUBROUTINE ini_date
703
704	SUBROUTINE fin_date( ddobsfin )
705	!!----------------------------------------------------------------------
706	!! * ROUTINE fin_date *
707	!!
708	!! ** Purpose : Get final data in double precision YYYYMMDD.HHMMSS format
709	!!
710	!! ** Method : Get final data in double precision YYYYMMDD.HHMMSS format
711	!!
712	!! ** Action : Get final data in double precision YYYYMMDD.HHMMSS format
713	!!
714	!! History :
715	!! ! 06-03 (K. Mogensen) Original code
716	!! ! 06-05 (K. Mogensen) Reformatted
717	!! ! 06-10 (A. Weaver) Cleaning
718	!! ! 10-05 (D. Lea) Update to month length calculation for NEMO vn3.2
719	!!----------------------------------------------------------------------
720	USE phycst, ONLY : & ! Physical constants
721	& rday
722	! USE daymod, ONLY : & ! Time variables
723	! & nmonth_len
724	USE dom_oce, ONLY : & ! Ocean space and time domain variables
725	& rdt
726
727	IMPLICIT NONE
728
729	!! * Arguments
730	REAL(KIND=dp), INTENT(OUT) :: ddobsfin ! Final date in YYYYMMDD.HHMMSS
731
732	!! * Local declarations
733	INTEGER :: iyea ! date - (year, month, day, hour, minute)
734	INTEGER :: imon
735	INTEGER :: iday
736	INTEGER :: ihou
737	INTEGER :: imin
738	INTEGER :: imday ! Number of days in month.
739	REAL(KIND=wp) :: zdayfrc ! Fraction of day
740
741	INTEGER, DIMENSION(12) :: imonth_len !: length in days of the months of the current year
742
743	!-----------------------------------------------------------------------
744	! Initial date initialization (year, month, day, hour, minute)
745	! (This assumes that the initial date is for 00z)
746	!-----------------------------------------------------------------------
747	iyea = ndate0 / 10000
748	imon = ( ndate0 - iyea * 10000 ) / 100
749	iday = ndate0 - iyea * 10000 - imon * 100
750	ihou = 0
751	imin = 0
752
753	!-----------------------------------------------------------------------
754	! Compute number of days + number of hours + min since initial time
755	!-----------------------------------------------------------------------
756	iday = iday + nitend * rdt / rday
757	zdayfrc = nitend * rdt / rday
758	zdayfrc = zdayfrc - AINT( zdayfrc )
759	ihou = INT( zdayfrc * 24 )
760	imin = INT( ( zdayfrc * 24 - ihou ) * 60 )
761
762	!-----------------------------------------------------------------------
763	! Convert number of days (iday) into a real date
764	!----------------------------------------------------------------------
765
766	CALL calc_month_len( iyea, imonth_len )
767
768	DO WHILE ( iday > imonth_len(imon) )
769	iday = iday - imonth_len(imon)
770	imon = imon + 1
771	IF ( imon > 12 ) THEN
772	imon = 1
773	iyea = iyea + 1
774	CALL calc_month_len( iyea, imonth_len ) ! update month lengths
775	ENDIF
776	END DO
777
778	!-----------------------------------------------------------------------
779	! Convert it into YYYYMMDD.HHMMSS format
780	!-----------------------------------------------------------------------
781	ddobsfin = iyea * 10000_dp + imon * 100_dp + iday &
782	& + ihou * 0.01_dp + imin * 0.0001_dp
783
784	END SUBROUTINE fin_date
785
786	END MODULE diaobs

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

Download in other formats: