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 NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS – NEMO

source: NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/diaobs.F90 @ 15240

Last change on this file since 15240 was 15240, checked in by dford, 3 years ago

Add option to use time mean background.

File size: 47.1 KB
Line 
1MODULE diaobs
2   !!======================================================================
3   !!                       ***  MODULE diaobs  ***
4   !! Observation diagnostics: Computation of the misfit between data and
5   !!                          their model equivalent
6   !!======================================================================
7   !! History :  1.0  !  2006-03  (K. Mogensen) Original code
8   !!             -   !  2006-05  (K. Mogensen, A. Weaver) Reformatted
9   !!             -   !  2006-10  (A. Weaver) Cleaning and add controls
10   !!             -   !  2007-03  (K. Mogensen) General handling of profiles
11   !!             -   !  2007-04  (G. Smith) Generalized surface operators
12   !!            2.0  !  2008-10  (M. Valdivieso) obs operator for velocity profiles
13   !!            3.4  !  2014-08  (J. While) observation operator for profiles in all vertical coordinates
14   !!             -   !                      Incorporated SST bias correction 
15   !!            3.6  !  2015-02  (M. Martin) Simplification of namelist and code
16   !!             -   !  2015-08  (M. Martin) Combined surface/profile routines.
17   !!            4.0  !  2017-11  (G. Madec) style only
18   !!----------------------------------------------------------------------
19
20   !!----------------------------------------------------------------------
21   !!   dia_obs_init  : Reading and prepare observations
22   !!   dia_obs       : Compute model equivalent to observations
23   !!   dia_obs_wri   : Write observational diagnostics
24   !!   calc_date     : Compute the date of timestep in YYYYMMDD.HHMMSS format
25   !!   ini_date      : Compute the initial date YYYYMMDD.HHMMSS
26   !!   fin_date      : Compute the final date YYYYMMDD.HHMMSS
27   !!----------------------------------------------------------------------
28   USE par_kind          ! Precision variables
29   USE in_out_manager    ! I/O manager
30   USE timing            ! Timing
31   USE par_oce           ! ocean parameter
32   USE dom_oce           ! Ocean space and time domain variables
33   USE sbc_oce           ! Sea-ice fraction
34   !
35   USE obs_read_prof     ! Reading and allocation of profile obs
36   USE obs_read_surf     ! Reading and allocation of surface obs
37   USE obs_bias          ! Bias correction routine
38   USE obs_readmdt       ! Reading and allocation of MDT for SLA.
39   USE obs_readsnowdepth ! Get model snow depth for conversion of freeboard to ice thickness
40   USE obs_prep          ! Preparation of obs. (grid search etc).
41   USE obs_oper          ! Observation operators
42   USE obs_write         ! Writing of observation related diagnostics
43   USE obs_grid          ! Grid searching
44   USE obs_read_altbias  ! Bias treatment for altimeter
45   USE obs_profiles_def  ! Profile data definitions
46   USE obs_surf_def      ! Surface data definitions
47   USE obs_types         ! Definitions for observation types
48   USE obs_group_def     ! Definitions for observation groups
49   !
50   USE mpp_map           ! MPP mapping
51   USE lib_mpp           ! For ctl_warn/stop
52
53   IMPLICIT NONE
54   PRIVATE
55
56   PUBLIC dia_obs_init     ! Initialize and read observations
57   PUBLIC dia_obs          ! Compute model equivalent to observations
58   PUBLIC dia_obs_wri      ! Write model equivalent to observations
59   PUBLIC calc_date        ! Compute the date of a timestep
60
61   LOGICAL, PUBLIC :: ln_diaobs            !: Logical switch for the obs operator
62   
63   INTEGER :: nn_obsgroups
64
65   TYPE(obs_group), DIMENSION(:), ALLOCATABLE ::   sobsgroups   ! Obs groups
66
67   !!----------------------------------------------------------------------
68   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
69   !! $Id$
70   !! Software governed by the CeCILL license (see ./LICENSE)
71   !!----------------------------------------------------------------------
72CONTAINS
73
74   SUBROUTINE dia_obs_init
75      !!----------------------------------------------------------------------
76      !!                    ***  ROUTINE dia_obs_init  ***
77      !!         
78      !! ** Purpose : Initialize and read observations
79      !!
80      !! ** Method  : Read the namelist and call reading routines
81      !!
82      !! ** Action  : Read the namelist and call reading routines
83      !!
84      !!----------------------------------------------------------------------
85#if defined key_si3
86   USE ice, ONLY : &     ! Sea ice variables
87      & hm_s             ! Snow depth for freeboard conversion
88#elif defined key_cice
89   USE sbc_oce, ONLY : & ! Sea ice variables
90      & thick_s          ! Snow depth for freeboard conversion
91#endif
92
93      IMPLICIT NONE
94
95      INTEGER, PARAMETER ::   jpmaxngroups = 1000    ! Maximum number of obs groups
96      INTEGER :: ios             ! Local integer output status for namelist read
97      INTEGER :: jtype           ! Counter for obs types
98      INTEGER :: jvar            ! Counter for variables
99      INTEGER :: jfile           ! Counter for files
100      INTEGER :: jenabled
101      INTEGER :: jgroup
102      !
103      LOGICAL, DIMENSION(:), ALLOCATABLE :: llvar
104      !
105      REAL(dp) :: rn_dobsini      ! Obs window start date YYYYMMDD.HHMMSS
106      REAL(dp) :: rn_dobsend      ! Obs window end date   YYYYMMDD.HHMMSS
107      !!
108      NAMELIST/namobs/ln_diaobs, nn_obsgroups,               &
109         &            ln_grid_global, ln_grid_search_lookup, &
110         &            cn_gridsearchfile, rn_gridsearchres,   &
111         &            rn_dobsini, rn_dobsend
112      !-----------------------------------------------------------------------
113
114      !-----------------------------------------------------------------------
115      ! Read namelist parameters
116      !-----------------------------------------------------------------------
117      ! Initialise time window
118      CALL ini_date( rn_dobsini )
119      CALL fin_date( rn_dobsend )
120
121      ! Read namelist namobs : control observation diagnostics
122      REWIND( numnam_ref )   ! Namelist namobs in reference namelist
123      READ  ( numnam_ref, namobs, IOSTAT = ios, ERR = 901)
124901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namobs in reference namelist' )
125      REWIND( numnam_cfg )   ! Namelist namobs in configuration namelist
126      READ  ( numnam_cfg, namobs, IOSTAT = ios, ERR = 902 )
127902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namobs in configuration namelist' )
128      IF(lwm) WRITE ( numond, namobs )
129
130      IF( .NOT.ln_diaobs ) THEN
131         IF(lwp) WRITE(numout,*)
132         IF(lwp) WRITE(numout,*) 'dia_obs_init : NO Observation diagnostic used'
133         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
134         RETURN
135      ENDIF
136
137      IF(lwp) THEN
138         WRITE(numout,*)
139         WRITE(numout,*) 'dia_obs_init : Observation diagnostic initialization'
140         WRITE(numout,*) '~~~~~~~~~~~~'
141         WRITE(numout,*) '   Namelist namobs : set observation diagnostic parameters'
142         WRITE(numout,*) '      Number of namobs_dta namelists to read             nn_obsgroups = ', nn_obsgroups
143         WRITE(numout,*) '      Global distribution of observations              ln_grid_global = ', ln_grid_global
144         WRITE(numout,*) '      Logical switch for obs grid search lookup ln_grid_search_lookup = ', ln_grid_search_lookup
145         IF (ln_grid_search_lookup) THEN
146            WRITE(numout,*) '      Grid search lookup file header                cn_gridsearchfile = ', cn_gridsearchfile
147            WRITE(numout,*) '      Grid search resolution                         rn_gridsearchres = ', rn_gridsearchres
148         ENDIF
149      ENDIF
150
151      IF( ln_grid_global ) THEN
152         CALL ctl_warn( 'dia_obs_init: ln_grid_global=T may cause memory issues when used with a large number of processors' )
153      ENDIF
154
155      !-----------------------------------------------------------------------
156      ! Read namobs_dta namelists and set up observation groups
157      !-----------------------------------------------------------------------
158
159      IF( nn_obsgroups == 0 ) THEN
160         CALL ctl_warn( 'dia_obs_init: ln_diaobs is set to true, but nn_obsgroups == 0',   &
161            &           ' so turning off calls to dia_obs'  )
162         ln_diaobs = .FALSE.
163         RETURN
164      ENDIF
165
166      ALLOCATE( sobsgroups(nn_obsgroups) )
167
168      jenabled = 0
169      DO jgroup = 1, nn_obsgroups
170         CALL obs_group_read ( sobsgroups(jgroup) )
171         CALL obs_group_check( sobsgroups(jgroup), jgroup )
172         IF (sobsgroups(jgroup)%lenabled) THEN
173            jenabled = jenabled + 1
174            IF( sobsgroups(jgroup)%lvel  .AND.  .NOT.ln_grid_global ) THEN
175               CALL ctl_stop( 'Velocity data only works with ln_grid_global=.true.' )
176            ENDIF
177         ENDIF
178      END DO
179
180      IF( jenabled == 0 ) THEN
181         CALL ctl_warn( 'dia_obs_init: ln_diaobs is set to true, and nn_obsgroups > 0',   &
182            &           ' but no groups are enabled so turning off calls to dia_obs'  )
183         ln_diaobs = .FALSE.
184         RETURN
185      ENDIF
186
187
188      !-----------------------------------------------------------------------
189      ! Obs operator parameter checking and initialisations
190      !-----------------------------------------------------------------------
191      !
192      CALL obs_typ_init
193      IF( ln_grid_global )   CALL mppmap_init
194      !
195      CALL obs_grid_setup( )
196
197      !-----------------------------------------------------------------------
198      ! Depending on switches read the various observation types
199      !-----------------------------------------------------------------------
200      !
201      DO jgroup = 1, nn_obsgroups
202         IF (sobsgroups(jgroup)%lenabled) THEN
203            IF (sobsgroups(jgroup)%lprof) THEN
204               !
205               ! Read in profile or profile obs types
206               !
207               ALLOCATE( llvar(sobsgroups(jgroup)%nobstypes) )
208               llvar(:) = .TRUE.
209               !
210               CALL obs_rea_prof( sobsgroups(jgroup)%sprofdata,   &
211                  &               sobsgroups(jgroup)%nobsfiles,   &
212                  &               sobsgroups(jgroup)%cobsfiles,   &
213                  &               sobsgroups(jgroup)%nobstypes,   &
214                  &               sobsgroups(jgroup)%naddvars,    &
215                  &               sobsgroups(jgroup)%nextvars,    &
216                  &               nitend-nit000+2,                &
217                  &               rn_dobsini,                     &
218                  &               rn_dobsend,                     &
219                  &               llvar,                          &
220                  &               sobsgroups(jgroup)%lignmis,     &
221                  &               sobsgroups(jgroup)%lall_at_all, &
222                  &               .FALSE.,                        &
223                  &               sobsgroups(jgroup)%cobstypes,   &
224                  &               kdailyavtypes = sobsgroups(jgroup)%nprofdavtypes )
225               !
226               DO jvar = 1, sobsgroups(jgroup)%nobstypes
227                  CALL obs_prof_staend( sobsgroups(jgroup)%sprofdata, jvar )
228               END DO
229               !
230               IF ( sobsgroups(jgroup)%sprofdata%next > 0 ) THEN
231                  CALL obs_prof_staend_ext( sobsgroups(jgroup)%sprofdata )
232               ENDIF
233               !
234               CALL obs_pre_prof( sobsgroups(jgroup)%sprofdata,     &
235                  &               sobsgroups(jgroup)%sprofdataqc,   &
236                  &               llvar,                            &
237                  &               jpi, jpj, jpk,                    &
238                  &               sobsgroups(jgroup)%rmask,         &
239                  &               sobsgroups(jgroup)%rglam,         &
240                  &               sobsgroups(jgroup)%rgphi,         &
241                  &               sobsgroups(jgroup)%lnea,          &
242                  &               sobsgroups(jgroup)%lbound_reject, &
243                  &               kdailyavtypes = sobsgroups(jgroup)%nprofdavtypes )
244               !
245               DEALLOCATE( llvar )
246               !
247            ELSEIF (sobsgroups(jgroup)%lsurf) THEN
248               !
249               ! Read in surface obs types
250               !
251               CALL obs_rea_surf( sobsgroups(jgroup)%ssurfdata,         &
252                  &               sobsgroups(jgroup)%nobsfiles,         &
253                  &               sobsgroups(jgroup)%cobsfiles,         &
254                  &               sobsgroups(jgroup)%nobstypes,         &
255                  &               sobsgroups(jgroup)%naddvars,          &
256                  &               sobsgroups(jgroup)%nextvars,          &
257                  &               nitend-nit000+2,                      &
258                  &               rn_dobsini,                           &
259                  &               rn_dobsend,                           &
260                  &               sobsgroups(jgroup)%rtime_mean_period, &
261                  &               sobsgroups(jgroup)%ltime_mean_bkg,    &
262                  &               sobsgroups(jgroup)%lignmis,           &
263                  &               .FALSE.,                              &
264                  &               sobsgroups(jgroup)%lnight,            &
265                  &               sobsgroups(jgroup)%cobstypes )
266               !
267               IF( sobsgroups(jgroup)%lsla ) THEN
268                  sobsgroups(jgroup)%ssurfdata%cextvars(sobsgroups(jgroup)%next_mdt) = 'MDT'
269                  sobsgroups(jgroup)%ssurfdata%cextlong(sobsgroups(jgroup)%next_mdt) = 'Mean dynamic topography'
270                  sobsgroups(jgroup)%ssurfdata%cextunit(sobsgroups(jgroup)%next_mdt) = 'Metres'
271                  sobsgroups(jgroup)%ssurfdata%caddvars(sobsgroups(jgroup)%nadd_ssh) = 'SSH'
272                  DO jvar = 1, sobsgroups(jgroup)%nobstypes
273                     sobsgroups(jgroup)%ssurfdata%caddlong(sobsgroups(jgroup)%nadd_ssh,jvar) = 'Model Sea surface height'
274                     sobsgroups(jgroup)%ssurfdata%caddunit(sobsgroups(jgroup)%nadd_ssh,jvar) = 'Metres'
275                  END DO
276               ENDIF
277               !
278               IF( sobsgroups(jgroup)%lfbd ) THEN
279                  sobsgroups(jgroup)%ssurfdata%cextvars(sobsgroups(jgroup)%next_snow) = 'SNOW'
280                  sobsgroups(jgroup)%ssurfdata%cextlong(sobsgroups(jgroup)%next_snow) = 'Snow thickness'
281                  sobsgroups(jgroup)%ssurfdata%cextunit(sobsgroups(jgroup)%next_snow) = 'Metres'
282                  sobsgroups(jgroup)%ssurfdata%caddvars(sobsgroups(jgroup)%nadd_fbd)  = 'FBD'
283                  DO jvar = 1, sobsgroups(jgroup)%nobstypes
284                     sobsgroups(jgroup)%ssurfdata%caddlong(sobsgroups(jgroup)%nadd_fbd,jvar) = 'Freeboard'
285                     sobsgroups(jgroup)%ssurfdata%caddunit(sobsgroups(jgroup)%nadd_fbd,jvar) = 'Metres'
286                  END DO
287               ENDIF
288               !
289               CALL obs_pre_surf( sobsgroups(jgroup)%ssurfdata,      &
290                  &               sobsgroups(jgroup)%ssurfdataqc,    &
291                  &               jpi, jpj,                          &
292                  &               sobsgroups(jgroup)%rmask(:,:,1,:), &
293                  &               sobsgroups(jgroup)%rglam,          &
294                  &               sobsgroups(jgroup)%rgphi,          &
295                  &               sobsgroups(jgroup)%lnea,           &
296                  &               sobsgroups(jgroup)%lbound_reject )
297               !
298               IF( sobsgroups(jgroup)%lsla ) THEN
299                  CALL obs_rea_mdt( sobsgroups(jgroup)%ssurfdataqc, &
300                     &              sobsgroups(jgroup)%n2dint,      &
301                     &              sobsgroups(jgroup)%next_mdt )
302                  IF( sobsgroups(jgroup)%laltbias ) THEN
303                     !CALL obs_rea_altbias( sobsgroups(jgroup)%ssurfdataqc, &
304                     !   &                  sobsgroups(jgroup)%n2dint,      &
305                     !   &                  sobsgroups(jgroup)%caltbiasfile )
306                     CALL obs_app_bias( sobsgroups(jgroup)%ssurfdataqc,   &
307                        &               sobsgroups(jgroup)%next_mdt,      & 
308                        &               sobsgroups(jgroup)%n2dint,        & 
309                        &               1,                                &
310                        &               sobsgroups(jgroup)%caltbiasfile,  &
311                        &               'altbias',                        &
312                        &               ld_extvar=.TRUE. ) 
313                  ENDIF
314               ENDIF
315               !
316#if defined key_si3
317               IF( sobsgroups(jgroup)%lfbd ) THEN
318                  CALL obs_rea_snowdepth( sobsgroups(jgroup)%ssurfdataqc, &
319                     &                    sobsgroups(jgroup)%n2dint,      &
320                     &                    sobsgroups(jgroup)%next_snow,   &
321                     &                    hm_s(:,:) )
322               ENDIF
323#elif defined key_cice
324               IF( sobsgroups(jgroup)%lfbd ) THEN
325                  CALL obs_rea_snowdepth( sobsgroups(jgroup)%ssurfdataqc, &
326                     &                    sobsgroups(jgroup)%n2dint,      &
327                     &                    sobsgroups(jgroup)%next_snow,   &
328                     &                    thick_s(:,:) )
329               ENDIF
330#endif
331               !
332               IF( sobsgroups(jgroup)%lobsbias ) THEN
333                  CALL obs_app_bias( sobsgroups(jgroup)%ssurfdataqc,   &
334                     &               sobsgroups(jgroup)%nbiasvar,      & 
335                     &               sobsgroups(jgroup)%n2dint,        & 
336                     &               sobsgroups(jgroup)%nobsbiasfiles, &
337                     &               sobsgroups(jgroup)%cobsbiasfiles, &
338                     &               sobsgroups(jgroup)%cbiasvarname ) 
339               ENDIF
340               !
341            ENDIF
342         ENDIF
343      END DO
344      !
345   END SUBROUTINE dia_obs_init
346
347
348   SUBROUTINE dia_obs( kstp )
349      !!----------------------------------------------------------------------
350      !!                    ***  ROUTINE dia_obs  ***
351      !!         
352      !! ** Purpose : Call the observation operators on each time step
353      !!
354      !! ** Method  : Call the observation operators on each time step to
355      !!              compute the model equivalent of the following data:
356      !!               - Profile data, currently T/S or U/V
357      !!               - Surface data, currently SST, SLA or sea-ice concentration.
358      !!
359      !! ** Action  :
360      !!----------------------------------------------------------------------
361      USE dom_oce, ONLY : gdept_n, gdept_1d   ! Ocean space and time domain variables
362      USE phycst , ONLY : rday                ! Physical constants
363      USE oce    , ONLY : tsn, un, vn, sshn   ! Ocean dynamics and tracers variables
364      USE phycst , ONLY : rday                ! Physical constants
365#if defined key_si3
366      USE ice    , ONLY : at_i, hm_i          ! SI3 Ice model variables
367#elif defined key_cice
368      USE sbc_oce, ONLY : fr_i, thick_i       ! CICE Ice model variables
369#endif
370
371      IMPLICIT NONE
372
373      !! * Arguments
374      INTEGER, INTENT(IN) :: kstp  ! Current timestep
375      !! * Local declarations
376      INTEGER :: idaystp           ! Number of timesteps per day
377      INTEGER :: imeanstp          ! Number of timesteps for time averaging
378      INTEGER :: jtype             ! Data loop variable
379      INTEGER :: jvar              ! Variable number
380      INTEGER :: jgroup
381      INTEGER :: ji, jj, jobs      ! Loop counters
382      LOGICAL :: lstp0             ! Flag special treatment on zeroth time step
383      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: &
384         & zprofvar                ! Model values for variables in a prof ob
385      REAL(wp), DIMENSION(:,:), ALLOCATABLE :: &
386         & zsurfvar                ! Model values for variables in a surf ob
387
388      !-----------------------------------------------------------------------
389
390      IF( ln_timing )   CALL timing_start('dia_obs')
391
392      IF(lwp) THEN
393         WRITE(numout,*)
394         WRITE(numout,*) 'dia_obs : Call the observation operators', kstp
395         WRITE(numout,*) '~~~~~~~'
396      ENDIF
397
398      idaystp = NINT( rday / rdt )
399
400      !-----------------------------------------------------------------------
401      ! Call the profile and surface observation operators
402      !-----------------------------------------------------------------------
403
404      DO jgroup = 1, nn_obsgroups
405         IF (sobsgroups(jgroup)%lenabled) THEN
406
407            IF (sobsgroups(jgroup)%lprof) THEN
408
409               ALLOCATE( zprofvar(jpi,jpj,jpk) )
410
411               DO jvar = 1, sobsgroups(jgroup)%nobstypes
412
413                  SELECT CASE ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) )
414                  CASE('POTM')
415                     zprofvar(:,:,:) = tsn(:,:,:,jp_tem)
416                  CASE('PSAL')
417                     zprofvar(:,:,:) = tsn(:,:,:,jp_sal)
418                  CASE('UVEL')
419                     zprofvar(:,:,:) = un(:,:,:)
420                  CASE('VVEL')
421                     zprofvar(:,:,:) = vn(:,:,:)
422                  CASE DEFAULT
423                     CALL ctl_stop( 'Unknown profile observation type '//TRIM(sobsgroups(jgroup)%cobstypes(jvar))//' in dia_obs' )
424                  END SELECT
425
426                  CALL obs_prof_opt( sobsgroups(jgroup)%sprofdataqc,       &
427                     &               kstp, jpi, jpj, jpk,                  &
428                     &               nit000, idaystp, jvar,                &
429                     &               zprofvar,                             &
430                     &               gdept_n,                              &
431                     &               gdepw_n,                              & 
432                     &               sobsgroups(jgroup)%rmask(:,:,:,jvar), &
433                     &               sobsgroups(jgroup)%rglam(:,:,jvar),   &
434                     &               sobsgroups(jgroup)%rgphi(:,:,jvar),   &
435                     &               sobsgroups(jgroup)%n1dint,            &
436                     &               sobsgroups(jgroup)%n2dint,            &
437                     &               kdailyavtypes = sobsgroups(jgroup)%nprofdavtypes )
438
439               END DO
440
441               DEALLOCATE( zprofvar )
442
443            ELSEIF (sobsgroups(jgroup)%lsurf) THEN
444
445               ALLOCATE( zsurfvar(jpi,jpj) )
446
447               DO jvar = 1, sobsgroups(jgroup)%nobstypes
448
449                  lstp0 = .FALSE.
450                  SELECT CASE ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) )
451                  CASE('SST')
452                     zsurfvar(:,:) = tsn(:,:,1,jp_tem)
453                  CASE('SLA')
454                     zsurfvar(:,:) = sshn(:,:)
455                  CASE('SSS')
456                     zsurfvar(:,:) = tsn(:,:,1,jp_sal)
457                  CASE('UVEL')
458                     zsurfvar(:,:) = un(:,:,1)
459                  CASE('VVEL')
460                     zsurfvar(:,:) = vn(:,:,1)
461                  CASE('ICECONC')
462                     IF ( kstp == 0 ) THEN
463                        lstp0 = .TRUE.
464                     ELSE
465#if defined key_si3
466                        zsurfvar(:,:) = at_i(:,:)
467#elif defined key_cice
468                        zsurfvar(:,:) = fr_i(:,:)
469#else
470                        CALL ctl_stop( ' Trying to run sea-ice observation operator', &
471                           &           ' but no sea-ice model appears to have been defined' )
472#endif
473                     ENDIF
474                  CASE('SIT','FBD')
475                     IF ( kstp == 0 ) THEN
476                        lstp0 = .TRUE.
477                     ELSE
478#if defined key_si3
479                        zsurfvar(:,:) = hm_i(:,:)
480#elif defined key_cice
481                        zsurfvar(:,:) = thick_i(:,:)
482#else
483                        CALL ctl_stop( ' Trying to run sea-ice observation operator', &
484                           &           ' but no sea-ice model appears to have been defined' )
485#endif
486                     ENDIF
487                  END SELECT
488
489                  IF ( lstp0 ) THEN
490                     IF ( sobsgroups(jgroup)%ssurfdataqc%nsstpmpp(1) > 0 ) THEN
491                        DO jobs = sobsgroups(jgroup)%ssurfdataqc%nsurfup + 1, &
492                           &      sobsgroups(jgroup)%ssurfdataqc%nsurfup + sobsgroups(jgroup)%ssurfdataqc%nsstp(1)
493                           sobsgroups(jgroup)%ssurfdata%nqc(jobs) = IBSET(sobsgroups(jgroup)%ssurfdata%nqc(jobs),13)
494                        END DO
495                        IF ( lwp ) THEN
496                           CALL ctl_warn( TRIM(sobsgroups(jgroup)%cobstypes(jvar))// &
497                              &           ' not initialised on zeroth '           // &
498                              &           'time-step but some obs are valid then.' )
499                           WRITE(numout,*)sobsgroups(jgroup)%ssurfdataqc%nsstpmpp(1), &
500                              &           TRIM(sobsgroups(jgroup)%cobstypes(jvar)),   &
501                              &           'observations will be flagged as bad'
502                        ENDIF
503                     ENDIF
504                     IF ( jvar == sobsgroups(jgroup)%ssurfdataqc%nvar ) THEN
505                        sobsgroups(jgroup)%ssurfdataqc%nsurfup = sobsgroups(jgroup)%ssurfdataqc%nsurfup + &
506                           &                                     sobsgroups(jgroup)%ssurfdataqc%nsstp(1)
507                     ENDIF
508                  ELSE
509                     IF ( sobsgroups(jgroup)%ltime_mean_bkg ) THEN
510                        ! Number of time-steps in meaning period
511                        imeanstp = NINT( ( sobsgroups(jgroup)%rtime_mean_period * 60.0 * 60.0 ) / rdt )
512                     ENDIF
513                     CALL obs_surf_opt( sobsgroups(jgroup)%ssurfdataqc,       &
514                        &               kstp, jpi, jpj,                       &
515                        &               nit000, idaystp,                      &
516                        &               jvar, zsurfvar,                       &
517                        &               sobsgroups(jgroup)%rmask(:,:,1,jvar), &
518                        &               sobsgroups(jgroup)%n2dint,            &
519                        &               sobsgroups(jgroup)%lnight,            &
520                        &               sobsgroups(jgroup)%ravglamscl,        &
521                        &               sobsgroups(jgroup)%ravgphiscl,        &
522                        &               sobsgroups(jgroup)%lfp_indegs,        &
523                        &               sobsgroups(jgroup)%ltime_mean_bkg,    &
524                        &               imeanstp,                             &
525                        &               kssh=sobsgroups(jgroup)%nadd_ssh,     &
526                        &               kmdt=sobsgroups(jgroup)%next_mdt,     &
527                        &               kfbd=sobsgroups(jgroup)%nadd_fbd,     &
528                        &               ksnow=sobsgroups(jgroup)%next_snow )
529                  ENDIF
530
531               END DO
532
533               DEALLOCATE( zsurfvar )
534
535            ENDIF
536
537         ENDIF
538      END DO
539
540      IF( ln_timing )   CALL timing_stop('dia_obs')
541
542   END SUBROUTINE dia_obs
543
544   SUBROUTINE dia_obs_wri
545      !!----------------------------------------------------------------------
546      !!                    ***  ROUTINE dia_obs_wri  ***
547      !!         
548      !! ** Purpose : Call observation diagnostic output routines
549      !!
550      !! ** Method  : Call observation diagnostic output routines
551      !!
552      !! ** Action  :
553      !!
554      !! History :
555      !!        !  06-03  (K. Mogensen) Original code
556      !!        !  06-05  (K. Mogensen) Reformatted
557      !!        !  06-10  (A. Weaver) Cleaning
558      !!        !  07-03  (K. Mogensen) General handling of profiles
559      !!        !  08-09  (M. Valdivieso) Velocity component (U,V) profiles
560      !!        !  15-08  (M. Martin) Combined writing for prof and surf types
561      !!----------------------------------------------------------------------
562      !! * Modules used
563      USE obs_rot_vel          ! Rotation of velocities
564
565      IMPLICIT NONE
566
567      !! * Local declarations
568      INTEGER :: jgroup                   ! Data set loop variable
569      INTEGER :: jo, jvar, jk, jadd, jext, jadd2, jext2
570      INTEGER :: iuvar, ivvar
571      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
572         & zu, &
573         & zv
574      LOGICAL, DIMENSION(:), ALLOCATABLE :: ll_write
575      TYPE(obswriinfo) :: sladd, slext
576
577      IF( ln_timing )   CALL timing_start('dia_obs_wri')
578
579      !-----------------------------------------------------------------------
580      ! Depending on switches call various observation output routines
581      !-----------------------------------------------------------------------
582
583      DO jgroup = 1, nn_obsgroups
584         IF (sobsgroups(jgroup)%lenabled) THEN
585
586            IF (sobsgroups(jgroup)%lprof) THEN
587
588               IF (sobsgroups(jgroup)%lvel) THEN
589                  iuvar = 0
590                  ivvar = 0
591                  DO jvar = 1, sobsgroups(jgroup)%nobstypes
592                     IF ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) == cobsname_uvel ) THEN
593                        iuvar = jvar
594                     ELSEIF ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) == cobsname_vvel ) THEN
595                        ivvar = jvar
596                     ENDIF
597                  END DO
598                  IF ( (iuvar > 0) .AND. (ivvar > 0) ) THEN
599
600                     ! For velocity data, rotate the model velocities to N/S, E/W
601                     ! using the compressed data structure.
602                     ALLOCATE( &
603                        & zu(sobsgroups(jgroup)%sprofdataqc%nvprot(iuvar)), &
604                        & zv(sobsgroups(jgroup)%sprofdataqc%nvprot(ivvar))  &
605                        & )
606
607                     CALL obs_rotvel_pro( sobsgroups(jgroup)%sprofdataqc, sobsgroups(jgroup)%n2dint, &
608                        &                 iuvar, ivvar, zu, zv )
609
610                     DO jo = 1, sobsgroups(jgroup)%sprofdataqc%nprof
611                        DO jk = sobsgroups(jgroup)%sprofdataqc%npvsta(jo,iuvar), sobsgroups(jgroup)%sprofdataqc%npvend(jo,iuvar)
612                           sobsgroups(jgroup)%sprofdataqc%var(iuvar)%vmod(jk) = zu(jk)
613                        END DO
614                        DO jk = sobsgroups(jgroup)%sprofdataqc%npvsta(jo,ivvar), sobsgroups(jgroup)%sprofdataqc%npvend(jo,ivvar)
615                           sobsgroups(jgroup)%sprofdataqc%var(ivvar)%vmod(jk) = zv(jk)
616                        END DO
617                     END DO
618
619                     DEALLOCATE( zu )
620                     DEALLOCATE( zv )
621
622                  ELSE
623                     CALL ctl_stop( 'Could not identify velocity observation variables to rotate' )
624                  END IF
625               END IF
626
627               CALL obs_prof_decompress( sobsgroups(jgroup)%sprofdataqc, &
628                  &                      sobsgroups(jgroup)%sprofdata, .TRUE., numout )
629
630               ! Put additional and extra variable information into obswriinfo structure
631               ! used by obs_write.
632               ! add/ext variables generated by the OBS code (1...sobsgroups(jgroup)%naddvars)
633               ! may duplicate ones read in (%naddvars+1...sobsgroups(jgroup)%sprofdata%nadd)
634               ! Check for this, and if so only write out the version generated by the OBS code
635               sladd%inum = sobsgroups(jgroup)%sprofdata%nadd
636               ALLOCATE( ll_write(sobsgroups(jgroup)%sprofdata%nadd) )
637               ll_write(:) = .TRUE.
638               IF ( (sobsgroups(jgroup)%naddvars > 0) .AND. &
639                  & (sobsgroups(jgroup)%sprofdata%nadd > sobsgroups(jgroup)%naddvars) ) THEN
640                  DO jadd = sobsgroups(jgroup)%naddvars + 1, sobsgroups(jgroup)%sprofdata%nadd
641                     DO jadd2 = 1, sobsgroups(jgroup)%naddvars
642                        IF ( TRIM(sobsgroups(jgroup)%sprofdata%caddvars(jadd )) == &
643                           & TRIM(sobsgroups(jgroup)%sprofdata%caddvars(jadd2)) ) THEN
644                           sladd%inum = sladd%inum - 1
645                           ll_write(jadd) = .FALSE.
646                        ENDIF
647                     END DO
648                  END DO
649               ENDIF
650               IF ( sladd%inum > 0 ) THEN
651                  ALLOCATE( sladd%ipoint(sladd%inum),                                   &
652                     &      sladd%cdname(sladd%inum),                                   &
653                     &      sladd%cdlong(sladd%inum,sobsgroups(jgroup)%sprofdata%nvar), &
654                     &      sladd%cdunit(sladd%inum,sobsgroups(jgroup)%sprofdata%nvar) )
655                  jadd2 = 0
656                  DO jadd = 1, sobsgroups(jgroup)%sprofdata%nadd
657                     IF ( ll_write(jadd) ) THEN
658                        jadd2 = jadd2 + 1
659                        sladd%ipoint(jadd2) = jadd
660                        sladd%cdname(jadd2) = sobsgroups(jgroup)%sprofdata%caddvars(jadd)
661                        DO jvar = 1, sobsgroups(jgroup)%sprofdata%nvar
662                           sladd%cdlong(jadd2,jvar) = sobsgroups(jgroup)%sprofdata%caddlong(jadd,jvar)
663                           sladd%cdunit(jadd2,jvar) = sobsgroups(jgroup)%sprofdata%caddunit(jadd,jvar)
664                        END DO
665                     ENDIF
666                  END DO
667               ENDIF
668               DEALLOCATE( ll_write )
669               
670               slext%inum = sobsgroups(jgroup)%sprofdata%next
671               ALLOCATE( ll_write(sobsgroups(jgroup)%sprofdata%next) )
672               ll_write(:) = .TRUE.
673               IF ( (sobsgroups(jgroup)%nextvars > 0) .AND. &
674                  & (sobsgroups(jgroup)%sprofdata%next > sobsgroups(jgroup)%nextvars) ) THEN
675                  DO jext = sobsgroups(jgroup)%nextvars + 1, sobsgroups(jgroup)%sprofdata%next
676                     DO jext2 = 1, sobsgroups(jgroup)%nextvars
677                        IF ( TRIM(sobsgroups(jgroup)%sprofdata%cextvars(jext )) == &
678                           & TRIM(sobsgroups(jgroup)%sprofdata%cextvars(jext2)) ) THEN
679                           slext%inum = slext%inum - 1
680                           ll_write(jext) = .FALSE.
681                        ENDIF
682                     END DO
683                  END DO
684               ENDIF
685               IF ( slext%inum > 0 ) THEN
686                  ALLOCATE( slext%ipoint(slext%inum),   &
687                     &      slext%cdname(slext%inum),   &
688                     &      slext%cdlong(slext%inum,1), &
689                     &      slext%cdunit(slext%inum,1) )
690                  jext2 = 0
691                  DO jext = 1, sobsgroups(jgroup)%sprofdata%next
692                     IF ( ll_write(jext) ) THEN
693                        jext2 = jext2 + 1
694                        slext%ipoint(jext2)   = jext
695                        slext%cdname(jext2)   = sobsgroups(jgroup)%sprofdata%cextvars(jext)
696                        slext%cdlong(jext2,1) = sobsgroups(jgroup)%sprofdata%cextlong(jext)
697                        slext%cdunit(jext2,1) = sobsgroups(jgroup)%sprofdata%cextunit(jext)
698                     ENDIF
699                  END DO
700               ENDIF
701               DEALLOCATE( ll_write )
702
703               CALL obs_wri_prof( sobsgroups(jgroup)%sprofdata, sobsgroups(jgroup)%cgroupname, sladd, slext )
704
705               IF ( sladd%inum > 0 ) THEN
706                  DEALLOCATE( sladd%ipoint, sladd%cdname, sladd%cdlong, sladd%cdunit )
707               ENDIF
708               IF ( slext%inum > 0 ) THEN
709                  DEALLOCATE( slext%ipoint, slext%cdname, slext%cdlong, slext%cdunit )
710               ENDIF
711
712            ELSEIF (sobsgroups(jgroup)%lsurf) THEN
713
714               IF (sobsgroups(jgroup)%lvel) THEN
715                  iuvar = 0
716                  ivvar = 0
717                  DO jvar = 1, sobsgroups(jgroup)%nobstypes
718                     IF ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) == cobsname_uvel ) THEN
719                        iuvar = jvar
720                     ELSEIF ( TRIM(sobsgroups(jgroup)%cobstypes(jvar)) == cobsname_vvel ) THEN
721                        ivvar = jvar
722                     ENDIF
723                  END DO
724                  IF ( (iuvar > 0) .AND. (ivvar > 0) ) THEN
725
726                     ! For velocity data, rotate the model velocities to N/S, E/W
727                     ! using the compressed data structure.
728                     ALLOCATE( &
729                        & zu(sobsgroups(jgroup)%ssurfdataqc%nsurf), &
730                        & zv(sobsgroups(jgroup)%ssurfdataqc%nsurf)  &
731                        & )
732
733                     CALL obs_rotvel_surf( sobsgroups(jgroup)%ssurfdataqc, sobsgroups(jgroup)%n2dint, &
734                        &                  iuvar, ivvar, zu, zv )
735
736                     DO jo = 1, sobsgroups(jgroup)%ssurfdataqc%nsurf
737                        sobsgroups(jgroup)%ssurfdataqc%rmod(jo,iuvar) = zu(jo)
738                        sobsgroups(jgroup)%ssurfdataqc%rmod(jo,ivvar) = zv(jo)
739                     END DO
740
741                     DEALLOCATE( zu )
742                     DEALLOCATE( zv )
743
744                  ELSE
745                     CALL ctl_stop( 'Could not identify velocity observation variables to rotate' )
746                  END IF
747               END IF
748
749               CALL obs_surf_decompress( sobsgroups(jgroup)%ssurfdataqc, &
750                  &                      sobsgroups(jgroup)%ssurfdata, .TRUE., numout )
751
752               IF (sobsgroups(jgroup)%lfbd) THEN
753                  ! Input observations were freeboard, but we're putputting ice thickness
754                  DO jvar = 1, sobsgroups(jgroup)%nobstypes
755                     IF ( sobsgroups(jgroup)%cobstypes(jvar) == cobsname_fbd ) THEN
756                        sobsgroups(jgroup)%ssurfdata%cvars(jvar) = 'SIT'
757                        sobsgroups(jgroup)%ssurfdata%clong(jvar) = 'Sea ice thickness'
758                        sobsgroups(jgroup)%ssurfdata%cunit(jvar) = 'm'
759                        EXIT
760                     ENDIF
761                  END DO
762               ENDIF
763
764               ! Put additional and extra variable information into obswriinfo structure
765               ! used by obs_write.
766               ! add/ext variables generated by the OBS code (1...sobsgroups(jgroup)%naddvars)
767               ! may duplicate ones read in (%naddvars+1...sobsgroups(jgroup)%ssurfdata%nadd)
768               ! Check for this, and if so only write out the version generated by the OBS code
769               sladd%inum = sobsgroups(jgroup)%ssurfdata%nadd
770               ALLOCATE( ll_write(sobsgroups(jgroup)%ssurfdata%nadd) )
771               ll_write(:) = .TRUE.
772               IF ( (sobsgroups(jgroup)%naddvars > 0) .AND. &
773                  & (sobsgroups(jgroup)%ssurfdata%nadd > sobsgroups(jgroup)%naddvars) ) THEN
774                  DO jadd = sobsgroups(jgroup)%naddvars + 1, sobsgroups(jgroup)%ssurfdata%nadd
775                     DO jadd2 = 1, sobsgroups(jgroup)%naddvars
776                        IF ( TRIM(sobsgroups(jgroup)%ssurfdata%caddvars(jadd )) == &
777                           & TRIM(sobsgroups(jgroup)%ssurfdata%caddvars(jadd2)) ) THEN
778                           sladd%inum = sladd%inum - 1
779                           ll_write(jadd) = .FALSE.
780                        ENDIF
781                     END DO
782                  END DO
783               ENDIF
784               IF ( sladd%inum > 0 ) THEN
785                  ALLOCATE( sladd%ipoint(sladd%inum),                                   &
786                     &      sladd%cdname(sladd%inum),                                   &
787                     &      sladd%cdlong(sladd%inum,sobsgroups(jgroup)%ssurfdata%nvar), &
788                     &      sladd%cdunit(sladd%inum,sobsgroups(jgroup)%ssurfdata%nvar) )
789                  jadd2 = 0
790                  DO jadd = 1, sobsgroups(jgroup)%ssurfdata%nadd
791                     IF ( ll_write(jadd) ) THEN
792                        jadd2 = jadd2 + 1
793                        sladd%ipoint(jadd2) = jadd
794                        sladd%cdname(jadd2) = sobsgroups(jgroup)%ssurfdata%caddvars(jadd)
795                        DO jvar = 1, sobsgroups(jgroup)%ssurfdata%nvar
796                           sladd%cdlong(jadd2,jvar) = sobsgroups(jgroup)%ssurfdata%caddlong(jadd,jvar)
797                           sladd%cdunit(jadd2,jvar) = sobsgroups(jgroup)%ssurfdata%caddunit(jadd,jvar)
798                        END DO
799                     ENDIF
800                  END DO
801               ENDIF
802               DEALLOCATE( ll_write )
803               
804               slext%inum = sobsgroups(jgroup)%ssurfdata%nextra
805               ALLOCATE( ll_write(sobsgroups(jgroup)%ssurfdata%nextra) )
806               ll_write(:) = .TRUE.
807               IF ( (sobsgroups(jgroup)%nextvars > 0) .AND. &
808                  & (sobsgroups(jgroup)%ssurfdata%nextra > sobsgroups(jgroup)%nextvars) ) THEN
809                  DO jext = sobsgroups(jgroup)%nextvars + 1, sobsgroups(jgroup)%ssurfdata%nextra
810                     DO jext2 = 1, sobsgroups(jgroup)%nextvars
811                        IF ( TRIM(sobsgroups(jgroup)%ssurfdata%cextvars(jext )) == &
812                           & TRIM(sobsgroups(jgroup)%ssurfdata%cextvars(jext2)) ) THEN
813                           slext%inum = slext%inum - 1
814                           ll_write(jext) = .FALSE.
815                        ENDIF
816                     END DO
817                  END DO
818               ENDIF
819               IF ( slext%inum > 0 ) THEN
820                  ALLOCATE( slext%ipoint(slext%inum),   &
821                     &      slext%cdname(slext%inum),   &
822                     &      slext%cdlong(slext%inum,1), &
823                     &      slext%cdunit(slext%inum,1) )
824                  jext2 = 0
825                  DO jext = 1, sobsgroups(jgroup)%ssurfdata%nextra
826                     IF ( ll_write(jext) ) THEN
827                        jext2 = jext2 + 1
828                        slext%ipoint(jext2)   = jext
829                        slext%cdname(jext2)   = sobsgroups(jgroup)%ssurfdata%cextvars(jext)
830                        slext%cdlong(jext2,1) = sobsgroups(jgroup)%ssurfdata%cextlong(jext)
831                        slext%cdunit(jext2,1) = sobsgroups(jgroup)%ssurfdata%cextunit(jext)
832                     ENDIF
833                  END DO
834               ENDIF
835               DEALLOCATE( ll_write )
836
837               CALL obs_wri_surf( sobsgroups(jgroup)%ssurfdata, sobsgroups(jgroup)%cgroupname, sladd, slext )
838
839               IF ( sladd%inum > 0 ) THEN
840                  DEALLOCATE( sladd%ipoint, sladd%cdname, sladd%cdlong, sladd%cdunit )
841               ENDIF
842               IF ( slext%inum > 0 ) THEN
843                  DEALLOCATE( slext%ipoint, slext%cdname, slext%cdlong, slext%cdunit )
844               ENDIF
845
846            ENDIF
847
848         ENDIF
849
850      END DO
851
852      IF( ln_timing )   CALL timing_stop('dia_obs_wri')
853
854   END SUBROUTINE dia_obs_wri
855
856   SUBROUTINE calc_date( kstp, ddobs )
857      !!----------------------------------------------------------------------
858      !!                    ***  ROUTINE calc_date  ***
859      !!         
860      !! ** Purpose : Get date in double precision YYYYMMDD.HHMMSS format
861      !!
862      !! ** Method  : Get date in double precision YYYYMMDD.HHMMSS format
863      !!
864      !! ** Action  : Get date in double precision YYYYMMDD.HHMMSS format
865      !!
866      !! ** Action  : Get initial date in double precision YYYYMMDD.HHMMSS format
867      !!
868      !! History :
869      !!        !  06-03  (K. Mogensen)  Original code
870      !!        !  06-05  (K. Mogensen)  Reformatted
871      !!        !  06-10  (A. Weaver) Cleaning
872      !!        !  06-10  (G. Smith) Calculates initial date the same as method for final date
873      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
874      !!        !  2014-09  (D. Lea) New generic routine now deals with arbitrary initial time of day
875      !!----------------------------------------------------------------------
876      USE phycst, ONLY : &            ! Physical constants
877         & rday
878      USE dom_oce, ONLY : &           ! Ocean space and time domain variables
879         & rdt
880
881      IMPLICIT NONE
882
883      !! * Arguments
884      REAL(KIND=dp), INTENT(OUT) :: ddobs                        ! Date in YYYYMMDD.HHMMSS
885      INTEGER :: kstp
886
887      !! * Local declarations
888      INTEGER :: iyea        ! date - (year, month, day, hour, minute)
889      INTEGER :: imon
890      INTEGER :: iday
891      INTEGER :: ihou
892      INTEGER :: imin
893      INTEGER :: imday       ! Number of days in month.
894      REAL(wp) :: zdayfrc    ! Fraction of day
895
896      INTEGER, DIMENSION(12) ::   imonth_len    !: length in days of the months of the current year
897
898      !!----------------------------------------------------------------------
899      !! Initial date initialization (year, month, day, hour, minute)
900      !!----------------------------------------------------------------------
901      iyea =   ndate0 / 10000
902      imon = ( ndate0 - iyea * 10000 ) / 100
903      iday =   ndate0 - iyea * 10000 - imon * 100
904      ihou =   nn_time0 / 100
905      imin = ( nn_time0 - ihou * 100 ) 
906
907      !!----------------------------------------------------------------------
908      !! Compute number of days + number of hours + min since initial time
909      !!----------------------------------------------------------------------
910      zdayfrc = kstp * rdt / rday
911      zdayfrc = zdayfrc - aint(zdayfrc)
912      imin = imin + int( zdayfrc * 24 * 60 ) 
913      DO WHILE (imin >= 60) 
914        imin=imin-60
915        ihou=ihou+1
916      END DO
917      DO WHILE (ihou >= 24)
918        ihou=ihou-24
919        iday=iday+1
920      END DO
921      iday = iday + kstp * rdt / rday 
922
923      !-----------------------------------------------------------------------
924      ! Convert number of days (iday) into a real date
925      !----------------------------------------------------------------------
926
927      CALL calc_month_len( iyea, imonth_len )
928
929      DO WHILE ( iday > imonth_len(imon) )
930         iday = iday - imonth_len(imon)
931         imon = imon + 1 
932         IF ( imon > 12 ) THEN
933            imon = 1
934            iyea = iyea + 1
935            CALL calc_month_len( iyea, imonth_len )  ! update month lengths
936         ENDIF
937      END DO
938
939      !----------------------------------------------------------------------
940      ! Convert it into YYYYMMDD.HHMMSS format.
941      !----------------------------------------------------------------------
942      ddobs = iyea * 10000_dp + imon * 100_dp + &
943         &    iday + ihou * 0.01_dp + imin * 0.0001_dp
944
945   END SUBROUTINE calc_date
946
947   SUBROUTINE ini_date( ddobsini )
948      !!----------------------------------------------------------------------
949      !!                    ***  ROUTINE ini_date  ***
950      !!         
951      !! ** Purpose : Get initial date in double precision YYYYMMDD.HHMMSS format
952      !!
953      !! ** Method  :
954      !!
955      !! ** Action  :
956      !!
957      !! History :
958      !!        !  06-03  (K. Mogensen)  Original code
959      !!        !  06-05  (K. Mogensen)  Reformatted
960      !!        !  06-10  (A. Weaver) Cleaning
961      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
962      !!        !  2014-09  (D. Lea) Change to call generic routine calc_date
963      !!----------------------------------------------------------------------
964
965      IMPLICIT NONE
966
967      !! * Arguments
968      REAL(KIND=dp), INTENT(OUT) :: ddobsini                   ! Initial date in YYYYMMDD.HHMMSS
969
970      CALL calc_date( nit000 - 1, ddobsini )
971
972   END SUBROUTINE ini_date
973
974   SUBROUTINE fin_date( ddobsfin )
975      !!----------------------------------------------------------------------
976      !!                    ***  ROUTINE fin_date  ***
977      !!         
978      !! ** Purpose : Get final date in double precision YYYYMMDD.HHMMSS format
979      !!
980      !! ** Method  :
981      !!
982      !! ** Action  :
983      !!
984      !! History :
985      !!        !  06-03  (K. Mogensen)  Original code
986      !!        !  06-05  (K. Mogensen)  Reformatted
987      !!        !  06-10  (A. Weaver) Cleaning
988      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
989      !!        !  2014-09  (D. Lea) Change to call generic routine calc_date
990      !!----------------------------------------------------------------------
991
992      IMPLICIT NONE
993
994      !! * Arguments
995      REAL(dp), INTENT(OUT) :: ddobsfin ! Final date in YYYYMMDD.HHMMSS
996
997      CALL calc_date( nitend, ddobsfin )
998
999   END SUBROUTINE fin_date
1000
1001END MODULE diaobs
Note: See TracBrowser for help on using the repository browser.