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/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/OBS – NEMO

source: NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/OBS/diaobs.F90 @ 11671

Last change on this file since 11671 was 11671, checked in by acc, 5 years ago

Branch 2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles. Final, non-substantive changes to complete this branch. These changes remove all REWIND statements on the old namelist fortran units (now character variables for internal files). These changes have been left until last since they are easily repeated via a script and it may be preferable to use the previous revision for merge purposes and reapply these last changes separately. This branch has been fully SETTE tested.
  • Property svn:keywords set to Id
File size: 46.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 par_oce        ! ocean parameter
31   USE dom_oce        ! Ocean space and time domain variables
32   USE sbc_oce        ! Sea-ice fraction
33   !
34   USE obs_read_prof  ! Reading and allocation of profile obs
35   USE obs_read_surf  ! Reading and allocation of surface obs
36   USE obs_sstbias    ! Bias correction routine for SST
37   USE obs_readmdt    ! Reading and allocation of MDT for SLA.
38   USE obs_prep       ! Preparation of obs. (grid search etc).
39   USE obs_oper       ! Observation operators
40   USE obs_write      ! Writing of observation related diagnostics
41   USE obs_grid       ! Grid searching
42   USE obs_read_altbias ! Bias treatment for altimeter
43   USE obs_profiles_def ! Profile data definitions
44   USE obs_surf_def   ! Surface data definitions
45   USE obs_types      ! Definitions for observation types
46   !
47   USE mpp_map        ! MPP mapping
48   USE lib_mpp        ! For ctl_warn/stop
49
50   IMPLICIT NONE
51   PRIVATE
52
53   PUBLIC dia_obs_init     ! Initialize and read observations
54   PUBLIC dia_obs          ! Compute model equivalent to observations
55   PUBLIC dia_obs_wri      ! Write model equivalent to observations
56   PUBLIC dia_obs_dealloc  ! Deallocate dia_obs data
57   PUBLIC calc_date        ! Compute the date of a timestep
58
59   LOGICAL, PUBLIC :: ln_diaobs          !: Logical switch for the obs operator
60   LOGICAL         :: ln_sstnight        !  Logical switch for night mean SST obs
61   LOGICAL         :: ln_sla_fp_indegs   !  T=> SLA obs footprint size specified in degrees, F=> in metres
62   LOGICAL         :: ln_sst_fp_indegs   !  T=> SST obs footprint size specified in degrees, F=> in metres
63   LOGICAL         :: ln_sss_fp_indegs   !  T=> SSS obs footprint size specified in degrees, F=> in metres
64   LOGICAL         :: ln_sic_fp_indegs   !  T=> sea-ice obs footprint size specified in degrees, F=> in metres
65
66   REAL(wp) ::   rn_sla_avglamscl   ! E/W diameter of SLA observation footprint (metres)
67   REAL(wp) ::   rn_sla_avgphiscl   ! N/S diameter of SLA observation footprint (metres)
68   REAL(wp) ::   rn_sst_avglamscl   ! E/W diameter of SST observation footprint (metres)
69   REAL(wp) ::   rn_sst_avgphiscl   ! N/S diameter of SST observation footprint (metres)
70   REAL(wp) ::   rn_sss_avglamscl   ! E/W diameter of SSS observation footprint (metres)
71   REAL(wp) ::   rn_sss_avgphiscl   ! N/S diameter of SSS observation footprint (metres)
72   REAL(wp) ::   rn_sic_avglamscl   ! E/W diameter of sea-ice observation footprint (metres)
73   REAL(wp) ::   rn_sic_avgphiscl   ! N/S diameter of sea-ice observation footprint (metres)
74
75   INTEGER :: nn_1dint       ! Vertical interpolation method
76   INTEGER :: nn_2dint       ! Default horizontal interpolation method
77   INTEGER :: nn_2dint_sla   ! SLA horizontal interpolation method
78   INTEGER :: nn_2dint_sst   ! SST horizontal interpolation method
79   INTEGER :: nn_2dint_sss   ! SSS horizontal interpolation method
80   INTEGER :: nn_2dint_sic   ! Seaice horizontal interpolation method
81   INTEGER, DIMENSION(imaxavtypes) ::   nn_profdavtypes   ! Profile data types representing a daily average
82   INTEGER :: nproftypes     ! Number of profile obs types
83   INTEGER :: nsurftypes     ! Number of surface obs types
84   INTEGER , DIMENSION(:), ALLOCATABLE ::   nvarsprof, nvarssurf   ! Number of profile & surface variables
85   INTEGER , DIMENSION(:), ALLOCATABLE ::   nextrprof, nextrsurf   ! Number of profile & surface extra variables
86   INTEGER , DIMENSION(:), ALLOCATABLE ::   n2dintsurf             ! Interpolation option for surface variables
87   REAL(wp), DIMENSION(:), ALLOCATABLE ::   zavglamscl, zavgphiscl ! E/W & N/S diameter of averaging footprint for surface variables
88   LOGICAL , DIMENSION(:), ALLOCATABLE ::   lfpindegs              ! T=> surface obs footprint size specified in degrees, F=> in metres
89   LOGICAL , DIMENSION(:), ALLOCATABLE ::   llnightav              ! Logical for calculating night-time averages
90
91   TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) ::   surfdata     !: Initial surface data
92   TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) ::   surfdataqc   !: Surface data after quality control
93   TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) ::   profdata     !: Initial profile data
94   TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) ::   profdataqc   !: Profile data after quality control
95
96   CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types
97
98   !!----------------------------------------------------------------------
99   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
100   !! $Id$
101   !! Software governed by the CeCILL license (see ./LICENSE)
102   !!----------------------------------------------------------------------
103CONTAINS
104
105   SUBROUTINE dia_obs_init
106      !!----------------------------------------------------------------------
107      !!                    ***  ROUTINE dia_obs_init  ***
108      !!         
109      !! ** Purpose : Initialize and read observations
110      !!
111      !! ** Method  : Read the namelist and call reading routines
112      !!
113      !! ** Action  : Read the namelist and call reading routines
114      !!
115      !!----------------------------------------------------------------------
116      INTEGER, PARAMETER ::   jpmaxnfiles = 1000    ! Maximum number of files for each obs type
117      INTEGER, DIMENSION(:), ALLOCATABLE ::   ifilesprof, ifilessurf   ! Number of profile & surface files
118      INTEGER :: ios             ! Local integer output status for namelist read
119      INTEGER :: jtype           ! Counter for obs types
120      INTEGER :: jvar            ! Counter for variables
121      INTEGER :: jfile           ! Counter for files
122      INTEGER :: jnumsstbias
123      !
124      CHARACTER(len=128), DIMENSION(jpmaxnfiles) :: &
125         & cn_profbfiles, &      ! T/S profile input filenames
126         & cn_sstfbfiles, &      ! Sea surface temperature input filenames
127         & cn_sssfbfiles, &      ! Sea surface salinity input filenames
128         & cn_slafbfiles, &      ! Sea level anomaly input filenames
129         & cn_sicfbfiles, &      ! Seaice concentration input filenames
130         & cn_velfbfiles, &      ! Velocity profile input filenames
131         & cn_sstbiasfiles      ! SST bias input filenames
132      CHARACTER(LEN=128) :: &
133         & cn_altbiasfile        ! Altimeter bias input filename
134      CHARACTER(len=128), DIMENSION(:,:), ALLOCATABLE :: &
135         & clproffiles, &        ! Profile filenames
136         & clsurffiles           ! Surface filenames
137         !
138      LOGICAL :: ln_t3d          ! Logical switch for temperature profiles
139      LOGICAL :: ln_s3d          ! Logical switch for salinity profiles
140      LOGICAL :: ln_sla          ! Logical switch for sea level anomalies
141      LOGICAL :: ln_sst          ! Logical switch for sea surface temperature
142      LOGICAL :: ln_sss          ! Logical switch for sea surface salinity
143      LOGICAL :: ln_sic          ! Logical switch for sea ice concentration
144      LOGICAL :: ln_vel3d        ! Logical switch for velocity (u,v) obs
145      LOGICAL :: ln_nea          ! Logical switch to remove obs near land
146      LOGICAL :: ln_altbias      ! Logical switch for altimeter bias
147      LOGICAL :: ln_sstbias      ! Logical switch for bias corection of SST
148      LOGICAL :: ln_ignmis       ! Logical switch for ignoring missing files
149      LOGICAL :: ln_s_at_t       ! Logical switch to compute model S at T obs
150      LOGICAL :: ln_bound_reject ! Logical to remove obs near boundaries in LAMs.
151      LOGICAL :: llvar1          ! Logical for profile variable 1
152      LOGICAL :: llvar2          ! Logical for profile variable 1
153      LOGICAL, DIMENSION(jpmaxnfiles) :: lmask ! Used for finding number of sstbias files
154      !
155      REAL(dp) :: rn_dobsini     ! Obs window start date YYYYMMDD.HHMMSS
156      REAL(dp) :: rn_dobsend     ! Obs window end date   YYYYMMDD.HHMMSS
157      REAL(wp), DIMENSION(jpi,jpj)     ::   zglam1, zglam2   ! Model longitudes for profile variable 1 & 2
158      REAL(wp), DIMENSION(jpi,jpj)     ::   zgphi1, zgphi2   ! Model latitudes  for profile variable 1 & 2
159      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmask1, zmask2   ! Model land/sea mask associated with variable 1 & 2
160      !!
161      NAMELIST/namobs/ln_diaobs, ln_t3d, ln_s3d, ln_sla,              &
162         &            ln_sst, ln_sic, ln_sss, ln_vel3d,               &
163         &            ln_altbias, ln_sstbias, ln_nea,                 &
164         &            ln_grid_global, ln_grid_search_lookup,          &
165         &            ln_ignmis, ln_s_at_t, ln_bound_reject,          &
166         &            ln_sstnight,                                    &
167         &            ln_sla_fp_indegs, ln_sst_fp_indegs,             &
168         &            ln_sss_fp_indegs, ln_sic_fp_indegs,             &
169         &            cn_profbfiles, cn_slafbfiles,                   &
170         &            cn_sstfbfiles, cn_sicfbfiles,                   &
171         &            cn_velfbfiles, cn_sssfbfiles,                   &
172         &            cn_sstbiasfiles, cn_altbiasfile,                &
173         &            cn_gridsearchfile, rn_gridsearchres,            &
174         &            rn_dobsini, rn_dobsend,                         &
175         &            rn_sla_avglamscl, rn_sla_avgphiscl,             &
176         &            rn_sst_avglamscl, rn_sst_avgphiscl,             &
177         &            rn_sss_avglamscl, rn_sss_avgphiscl,             &
178         &            rn_sic_avglamscl, rn_sic_avgphiscl,             &
179         &            nn_1dint, nn_2dint,                             &
180         &            nn_2dint_sla, nn_2dint_sst,                     &
181         &            nn_2dint_sss, nn_2dint_sic,                     &
182         &            nn_msshc, rn_mdtcorr, rn_mdtcutoff,             &
183         &            nn_profdavtypes
184      !-----------------------------------------------------------------------
185
186      !-----------------------------------------------------------------------
187      ! Read namelist parameters
188      !-----------------------------------------------------------------------
189      ! Some namelist arrays need initialising
190      cn_profbfiles  (:) = ''
191      cn_slafbfiles  (:) = ''
192      cn_sstfbfiles  (:) = ''
193      cn_sicfbfiles  (:) = ''
194      cn_velfbfiles  (:) = ''
195      cn_sssfbfiles  (:) = ''
196      cn_sstbiasfiles(:) = ''
197      nn_profdavtypes(:) = -1
198
199      CALL ini_date( rn_dobsini )
200      CALL fin_date( rn_dobsend )
201
202      ! Read namelist namobs : control observation diagnostics
203      READ  ( numnam_ref, namobs, IOSTAT = ios, ERR = 901)
204901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namobs in reference namelist' )
205      READ  ( numnam_cfg, namobs, IOSTAT = ios, ERR = 902 )
206902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namobs in configuration namelist' )
207      IF(lwm) WRITE ( numond, namobs )
208
209      IF( .NOT.ln_diaobs ) THEN
210         IF(lwp) WRITE(numout,*)
211         IF(lwp) WRITE(numout,*) 'dia_obs_init : NO Observation diagnostic used'
212         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
213         RETURN
214      ENDIF
215
216      IF(lwp) THEN
217         WRITE(numout,*)
218         WRITE(numout,*) 'dia_obs_init : Observation diagnostic initialization'
219         WRITE(numout,*) '~~~~~~~~~~~~'
220         WRITE(numout,*) '   Namelist namobs : set observation diagnostic parameters' 
221         WRITE(numout,*) '      Logical switch for T profile observations                ln_t3d = ', ln_t3d
222         WRITE(numout,*) '      Logical switch for S profile observations                ln_s3d = ', ln_s3d
223         WRITE(numout,*) '      Logical switch for SLA observations                      ln_sla = ', ln_sla
224         WRITE(numout,*) '      Logical switch for SST observations                      ln_sst = ', ln_sst
225         WRITE(numout,*) '      Logical switch for Sea Ice observations                  ln_sic = ', ln_sic
226         WRITE(numout,*) '      Logical switch for velocity observations               ln_vel3d = ', ln_vel3d
227         WRITE(numout,*) '      Logical switch for SSS observations                      ln_sss = ', ln_sss
228         WRITE(numout,*) '      Global distribution of observations              ln_grid_global = ', ln_grid_global
229         WRITE(numout,*) '      Logical switch for obs grid search lookup ln_grid_search_lookup = ', ln_grid_search_lookup
230         IF (ln_grid_search_lookup) &
231            WRITE(numout,*) '      Grid search lookup file header                cn_gridsearchfile = ', cn_gridsearchfile
232         WRITE(numout,*) '      Initial date in window YYYYMMDD.HHMMSS               rn_dobsini = ', rn_dobsini
233         WRITE(numout,*) '      Final date in window YYYYMMDD.HHMMSS                 rn_dobsend = ', rn_dobsend
234         WRITE(numout,*) '      Type of vertical interpolation method                  nn_1dint = ', nn_1dint
235         WRITE(numout,*) '      Type of horizontal interpolation method                nn_2dint = ', nn_2dint
236         WRITE(numout,*) '      Rejection of observations near land switch               ln_nea = ', ln_nea
237         WRITE(numout,*) '      Rejection of obs near open bdys                 ln_bound_reject = ', ln_bound_reject
238         WRITE(numout,*) '      MSSH correction scheme                                 nn_msshc = ', nn_msshc
239         WRITE(numout,*) '      MDT  correction                                      rn_mdtcorr = ', rn_mdtcorr
240         WRITE(numout,*) '      MDT cutoff for computed correction                 rn_mdtcutoff = ', rn_mdtcutoff
241         WRITE(numout,*) '      Logical switch for alt bias                          ln_altbias = ', ln_altbias
242         WRITE(numout,*) '      Logical switch for sst bias                          ln_sstbias = ', ln_sstbias
243         WRITE(numout,*) '      Logical switch for ignoring missing files             ln_ignmis = ', ln_ignmis
244         WRITE(numout,*) '      Daily average types                             nn_profdavtypes = ', nn_profdavtypes
245         WRITE(numout,*) '      Logical switch for night-time SST obs               ln_sstnight = ', ln_sstnight
246      ENDIF
247      !-----------------------------------------------------------------------
248      ! Set up list of observation types to be used
249      ! and the files associated with each type
250      !-----------------------------------------------------------------------
251
252      nproftypes = COUNT( (/ln_t3d .OR. ln_s3d, ln_vel3d /) )
253      nsurftypes = COUNT( (/ln_sla, ln_sst, ln_sic, ln_sss /) )
254
255      IF( ln_sstbias ) THEN
256         lmask(:) = .FALSE. 
257         WHERE( cn_sstbiasfiles(:) /= '' )   lmask(:) = .TRUE. 
258         jnumsstbias = COUNT(lmask) 
259         lmask(:) = .FALSE. 
260      ENDIF     
261
262      IF( nproftypes == 0 .AND. nsurftypes == 0 ) THEN
263         CALL ctl_warn( 'dia_obs_init: ln_diaobs is set to true, but all obs operator logical flags',   &
264            &           ' (ln_t3d, ln_s3d, ln_sla, ln_sst, ln_sic, ln_vel3d)',                          &
265            &           ' are set to .FALSE. so turning off calls to dia_obs'  )
266         ln_diaobs = .FALSE.
267         RETURN
268      ENDIF
269
270      IF( nproftypes > 0 ) THEN
271         !
272         ALLOCATE( cobstypesprof(nproftypes)             )
273         ALLOCATE( ifilesprof   (nproftypes)             )
274         ALLOCATE( clproffiles  (nproftypes,jpmaxnfiles) )
275         !
276         jtype = 0
277         IF( ln_t3d .OR. ln_s3d ) THEN
278            jtype = jtype + 1
279            CALL obs_settypefiles( nproftypes, jpmaxnfiles, jtype, 'prof  ', &
280               &                   cn_profbfiles, ifilesprof, cobstypesprof, clproffiles )
281         ENDIF
282         IF( ln_vel3d ) THEN
283            jtype = jtype + 1
284            CALL obs_settypefiles( nproftypes, jpmaxnfiles, jtype, 'vel   ', &
285               &                   cn_velfbfiles, ifilesprof, cobstypesprof, clproffiles )
286         ENDIF
287         !
288      ENDIF
289
290      IF( nsurftypes > 0 ) THEN
291         !
292         ALLOCATE( cobstypessurf(nsurftypes)             )
293         ALLOCATE( ifilessurf   (nsurftypes)             )
294         ALLOCATE( clsurffiles  (nsurftypes,jpmaxnfiles) )
295         ALLOCATE( n2dintsurf   (nsurftypes)             )
296         ALLOCATE( zavglamscl   (nsurftypes)             )
297         ALLOCATE( zavgphiscl   (nsurftypes)             )
298         ALLOCATE( lfpindegs    (nsurftypes)             )
299         ALLOCATE( llnightav    (nsurftypes)             )
300         !
301         jtype = 0
302         IF( ln_sla ) THEN
303            jtype = jtype + 1
304            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sla   ', &
305               &                   cn_slafbfiles, ifilessurf, cobstypessurf, clsurffiles )
306            CALL obs_setinterpopts( nsurftypes, jtype, 'sla   ',      &
307               &                  nn_2dint, nn_2dint_sla,             &
308               &                  rn_sla_avglamscl, rn_sla_avgphiscl, &
309               &                  ln_sla_fp_indegs, .FALSE.,          &
310               &                  n2dintsurf, zavglamscl, zavgphiscl, &
311               &                  lfpindegs, llnightav )
312         ENDIF
313         IF( ln_sst ) THEN
314            jtype = jtype + 1
315            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sst   ', &
316               &                   cn_sstfbfiles, ifilessurf, cobstypessurf, clsurffiles )
317            CALL obs_setinterpopts( nsurftypes, jtype, 'sst   ',      &
318               &                  nn_2dint, nn_2dint_sst,             &
319               &                  rn_sst_avglamscl, rn_sst_avgphiscl, &
320               &                  ln_sst_fp_indegs, ln_sstnight,      &
321               &                  n2dintsurf, zavglamscl, zavgphiscl, &
322               &                  lfpindegs, llnightav )
323         ENDIF
324#if defined key_si3 || defined key_cice
325         IF( ln_sic ) THEN
326            jtype = jtype + 1
327            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sic   ', &
328               &                   cn_sicfbfiles, ifilessurf, cobstypessurf, clsurffiles )
329            CALL obs_setinterpopts( nsurftypes, jtype, 'sic   ',      &
330               &                  nn_2dint, nn_2dint_sic,             &
331               &                  rn_sic_avglamscl, rn_sic_avgphiscl, &
332               &                  ln_sic_fp_indegs, .FALSE.,          &
333               &                  n2dintsurf, zavglamscl, zavgphiscl, &
334               &                  lfpindegs, llnightav )
335         ENDIF
336#endif
337         IF( ln_sss ) THEN
338            jtype = jtype + 1
339            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sss   ', &
340               &                   cn_sssfbfiles, ifilessurf, cobstypessurf, clsurffiles )
341            CALL obs_setinterpopts( nsurftypes, jtype, 'sss   ',      &
342               &                  nn_2dint, nn_2dint_sss,             &
343               &                  rn_sss_avglamscl, rn_sss_avgphiscl, &
344               &                  ln_sss_fp_indegs, .FALSE.,          &
345               &                  n2dintsurf, zavglamscl, zavgphiscl, &
346               &                  lfpindegs, llnightav )
347         ENDIF
348         !
349      ENDIF
350
351
352      !-----------------------------------------------------------------------
353      ! Obs operator parameter checking and initialisations
354      !-----------------------------------------------------------------------
355      !
356      IF( ln_vel3d  .AND.  .NOT.ln_grid_global ) THEN
357         CALL ctl_stop( 'Velocity data only works with ln_grid_global=.true.' )
358         RETURN
359      ENDIF
360      !
361      IF( ln_grid_global ) THEN
362         CALL ctl_warn( 'dia_obs_init: ln_grid_global=T may cause memory issues when used with a large number of processors' )
363      ENDIF
364      !
365      IF( nn_1dint < 0  .OR.  nn_1dint > 1 ) THEN
366         CALL ctl_stop('dia_obs_init: Choice of vertical (1D) interpolation method is not available')
367      ENDIF
368      !
369      IF( nn_2dint < 0  .OR.  nn_2dint > 6  ) THEN
370         CALL ctl_stop('dia_obs_init: Choice of horizontal (2D) interpolation method is not available')
371      ENDIF
372      !
373      CALL obs_typ_init
374      IF( ln_grid_global )   CALL mppmap_init
375      !
376      CALL obs_grid_setup( )
377
378      !-----------------------------------------------------------------------
379      ! Depending on switches read the various observation types
380      !-----------------------------------------------------------------------
381      !
382      IF( nproftypes > 0 ) THEN
383         !
384         ALLOCATE( profdata  (nproftypes) , nvarsprof (nproftypes) )
385         ALLOCATE( profdataqc(nproftypes) , nextrprof (nproftypes) )
386         !
387         DO jtype = 1, nproftypes
388            !
389            nvarsprof(jtype) = 2
390            IF ( TRIM(cobstypesprof(jtype)) == 'prof' ) THEN
391               nextrprof(jtype) = 1
392               llvar1 = ln_t3d
393               llvar2 = ln_s3d
394               zglam1 = glamt
395               zgphi1 = gphit
396               zmask1 = tmask
397               zglam2 = glamt
398               zgphi2 = gphit
399               zmask2 = tmask
400            ENDIF
401            IF ( TRIM(cobstypesprof(jtype)) == 'vel' )  THEN
402               nextrprof(jtype) = 2
403               llvar1 = ln_vel3d
404               llvar2 = ln_vel3d
405               zglam1 = glamu
406               zgphi1 = gphiu
407               zmask1 = umask
408               zglam2 = glamv
409               zgphi2 = gphiv
410               zmask2 = vmask
411            ENDIF
412            !
413            ! Read in profile or profile obs types
414            CALL obs_rea_prof( profdata(jtype), ifilesprof(jtype),       &
415               &               clproffiles(jtype,1:ifilesprof(jtype)), &
416               &               nvarsprof(jtype), nextrprof(jtype), nitend-nit000+2, &
417               &               rn_dobsini, rn_dobsend, llvar1, llvar2, &
418               &               ln_ignmis, ln_s_at_t, .FALSE., &
419               &               kdailyavtypes = nn_profdavtypes )
420               !
421            DO jvar = 1, nvarsprof(jtype)
422               CALL obs_prof_staend( profdata(jtype), jvar )
423            END DO
424            !
425            CALL obs_pre_prof( profdata(jtype), profdataqc(jtype), &
426               &               llvar1, llvar2, &
427               &               jpi, jpj, jpk, &
428               &               zmask1, zglam1, zgphi1, zmask2, zglam2, zgphi2,  &
429               &               ln_nea, ln_bound_reject, &
430               &               kdailyavtypes = nn_profdavtypes )
431         END DO
432         !
433         DEALLOCATE( ifilesprof, clproffiles )
434         !
435      ENDIF
436      !
437      IF( nsurftypes > 0 ) THEN
438         !
439         ALLOCATE( surfdata  (nsurftypes) , nvarssurf(nsurftypes) )
440         ALLOCATE( surfdataqc(nsurftypes) , nextrsurf(nsurftypes) )
441         !
442         DO jtype = 1, nsurftypes
443            !
444            nvarssurf(jtype) = 1
445            nextrsurf(jtype) = 0
446            llnightav(jtype) = .FALSE.
447            IF( TRIM(cobstypessurf(jtype)) == 'sla' )   nextrsurf(jtype) = 2
448            IF( TRIM(cobstypessurf(jtype)) == 'sst' )   llnightav(jtype) = ln_sstnight
449            !
450            ! Read in surface obs types
451            CALL obs_rea_surf( surfdata(jtype), ifilessurf(jtype), &
452               &               clsurffiles(jtype,1:ifilessurf(jtype)), &
453               &               nvarssurf(jtype), nextrsurf(jtype), nitend-nit000+2, &
454               &               rn_dobsini, rn_dobsend, ln_ignmis, .FALSE., llnightav(jtype) )
455               !
456            CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), ln_nea, ln_bound_reject )
457            !
458            IF( TRIM(cobstypessurf(jtype)) == 'sla' ) THEN
459               CALL obs_rea_mdt( surfdataqc(jtype), n2dintsurf(jtype) )
460               IF( ln_altbias )   &
461                  & CALL obs_rea_altbias ( surfdataqc(jtype), n2dintsurf(jtype), cn_altbiasfile )
462            ENDIF
463            !
464            IF( TRIM(cobstypessurf(jtype)) == 'sst' .AND. ln_sstbias ) THEN
465               jnumsstbias = 0
466               DO jfile = 1, jpmaxnfiles
467                  IF( TRIM(cn_sstbiasfiles(jfile)) /= '' )   jnumsstbias = jnumsstbias + 1
468               END DO
469               IF( jnumsstbias == 0 )   CALL ctl_stop("ln_sstbias set but no bias files to read in")   
470               !
471               CALL obs_app_sstbias( surfdataqc(jtype), n2dintsurf(jtype)             ,   & 
472                  &                  jnumsstbias      , cn_sstbiasfiles(1:jnumsstbias) ) 
473            ENDIF
474         END DO
475         !
476         DEALLOCATE( ifilessurf, clsurffiles )
477         !
478      ENDIF
479      !
480   END SUBROUTINE dia_obs_init
481
482
483   SUBROUTINE dia_obs( kstp )
484      !!----------------------------------------------------------------------
485      !!                    ***  ROUTINE dia_obs  ***
486      !!         
487      !! ** Purpose : Call the observation operators on each time step
488      !!
489      !! ** Method  : Call the observation operators on each time step to
490      !!              compute the model equivalent of the following data:
491      !!               - Profile data, currently T/S or U/V
492      !!               - Surface data, currently SST, SLA or sea-ice concentration.
493      !!
494      !! ** Action  :
495      !!----------------------------------------------------------------------
496      USE dom_oce, ONLY : gdept_n, gdept_1d   ! Ocean space and time domain variables
497      USE phycst , ONLY : rday                ! Physical constants
498      USE oce    , ONLY : tsn, un, vn, sshn   ! Ocean dynamics and tracers variables
499      USE phycst , ONLY : rday                ! Physical constants
500#if defined  key_si3
501      USE ice    , ONLY : at_i                ! SI3 Ice model variables
502#endif
503#if defined key_cice
504      USE sbc_oce, ONLY : fr_i     ! ice fraction
505#endif
506
507      IMPLICIT NONE
508
509      !! * Arguments
510      INTEGER, INTENT(IN) :: kstp  ! Current timestep
511      !! * Local declarations
512      INTEGER :: idaystp           ! Number of timesteps per day
513      INTEGER :: jtype             ! Data loop variable
514      INTEGER :: jvar              ! Variable number
515      INTEGER :: ji, jj            ! Loop counters
516      REAL(wp), DIMENSION(jpi,jpj,jpk) :: &
517         & zprofvar1, &            ! Model values for 1st variable in a prof ob
518         & zprofvar2               ! Model values for 2nd variable in a prof ob
519      REAL(wp), DIMENSION(jpi,jpj,jpk) :: &
520         & zprofmask1, &           ! Mask associated with zprofvar1
521         & zprofmask2              ! Mask associated with zprofvar2
522      REAL(wp), DIMENSION(jpi,jpj) :: &
523         & zsurfvar, &             ! Model values equivalent to surface ob.
524         & zsurfmask               ! Mask associated with surface variable
525      REAL(wp), DIMENSION(jpi,jpj) :: &
526         & zglam1,    &            ! Model longitudes for prof variable 1
527         & zglam2,    &            ! Model longitudes for prof variable 2
528         & zgphi1,    &            ! Model latitudes for prof variable 1
529         & zgphi2                  ! Model latitudes for prof variable 2
530
531      !-----------------------------------------------------------------------
532
533      IF(lwp) THEN
534         WRITE(numout,*)
535         WRITE(numout,*) 'dia_obs : Call the observation operators', kstp
536         WRITE(numout,*) '~~~~~~~'
537      ENDIF
538
539      idaystp = NINT( rday / rdt )
540
541      !-----------------------------------------------------------------------
542      ! Call the profile and surface observation operators
543      !-----------------------------------------------------------------------
544
545      IF ( nproftypes > 0 ) THEN
546
547         DO jtype = 1, nproftypes
548
549            SELECT CASE ( TRIM(cobstypesprof(jtype)) )
550            CASE('prof')
551               zprofvar1(:,:,:) = tsn(:,:,:,jp_tem)
552               zprofvar2(:,:,:) = tsn(:,:,:,jp_sal)
553               zprofmask1(:,:,:) = tmask(:,:,:)
554               zprofmask2(:,:,:) = tmask(:,:,:)
555               zglam1(:,:) = glamt(:,:)
556               zglam2(:,:) = glamt(:,:)
557               zgphi1(:,:) = gphit(:,:)
558               zgphi2(:,:) = gphit(:,:)
559            CASE('vel')
560               zprofvar1(:,:,:) = un(:,:,:)
561               zprofvar2(:,:,:) = vn(:,:,:)
562               zprofmask1(:,:,:) = umask(:,:,:)
563               zprofmask2(:,:,:) = vmask(:,:,:)
564               zglam1(:,:) = glamu(:,:)
565               zglam2(:,:) = glamv(:,:)
566               zgphi1(:,:) = gphiu(:,:)
567               zgphi2(:,:) = gphiv(:,:)
568            CASE DEFAULT
569               CALL ctl_stop( 'Unknown profile observation type '//TRIM(cobstypesprof(jtype))//' in dia_obs' )
570            END SELECT
571
572            CALL obs_prof_opt( profdataqc(jtype), kstp, jpi, jpj, jpk,  &
573               &               nit000, idaystp,                         &
574               &               zprofvar1, zprofvar2,                    &
575               &               gdept_n(:,:,:), gdepw_n(:,:,:),            & 
576               &               zprofmask1, zprofmask2,                  &
577               &               zglam1, zglam2, zgphi1, zgphi2,          &
578               &               nn_1dint, nn_2dint,                      &
579               &               kdailyavtypes = nn_profdavtypes )
580
581         END DO
582
583      ENDIF
584
585      IF ( nsurftypes > 0 ) THEN
586
587         DO jtype = 1, nsurftypes
588
589            !Defaults which might be changed
590            zsurfmask(:,:) = tmask(:,:,1)
591
592            SELECT CASE ( TRIM(cobstypessurf(jtype)) )
593            CASE('sst')
594               zsurfvar(:,:) = tsn(:,:,1,jp_tem)
595            CASE('sla')
596               zsurfvar(:,:) = sshn(:,:)
597            CASE('sss')
598               zsurfvar(:,:) = tsn(:,:,1,jp_sal)
599            CASE('sic')
600               IF ( kstp == 0 ) THEN
601                  IF ( lwp .AND. surfdataqc(jtype)%nsstpmpp(1) > 0 ) THEN
602                     CALL ctl_warn( 'Sea-ice not initialised on zeroth '// &
603                        &           'time-step but some obs are valid then.' )
604                     WRITE(numout,*)surfdataqc(jtype)%nsstpmpp(1), &
605                        &           ' sea-ice obs will be missed'
606                  ENDIF
607                  surfdataqc(jtype)%nsurfup = surfdataqc(jtype)%nsurfup + &
608                     &                        surfdataqc(jtype)%nsstp(1)
609                  CYCLE
610               ELSE
611#if defined key_cice || defined key_si3
612                  zsurfvar(:,:) = fr_i(:,:)
613#else
614                  CALL ctl_stop( ' Trying to run sea-ice observation operator', &
615                     &           ' but no sea-ice model appears to have been defined' )
616#endif
617               ENDIF
618
619            END SELECT
620
621            CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj,       &
622               &               nit000, idaystp, zsurfvar, zsurfmask,    &
623               &               n2dintsurf(jtype), llnightav(jtype),     &
624               &               zavglamscl(jtype), zavgphiscl(jtype),     &
625               &               lfpindegs(jtype) )
626
627         END DO
628
629      ENDIF
630
631   END SUBROUTINE dia_obs
632
633   SUBROUTINE dia_obs_wri
634      !!----------------------------------------------------------------------
635      !!                    ***  ROUTINE dia_obs_wri  ***
636      !!         
637      !! ** Purpose : Call observation diagnostic output routines
638      !!
639      !! ** Method  : Call observation diagnostic output routines
640      !!
641      !! ** Action  :
642      !!
643      !! History :
644      !!        !  06-03  (K. Mogensen) Original code
645      !!        !  06-05  (K. Mogensen) Reformatted
646      !!        !  06-10  (A. Weaver) Cleaning
647      !!        !  07-03  (K. Mogensen) General handling of profiles
648      !!        !  08-09  (M. Valdivieso) Velocity component (U,V) profiles
649      !!        !  15-08  (M. Martin) Combined writing for prof and surf types
650      !!----------------------------------------------------------------------
651      !! * Modules used
652      USE obs_rot_vel          ! Rotation of velocities
653
654      IMPLICIT NONE
655
656      !! * Local declarations
657      INTEGER :: jtype                    ! Data set loop variable
658      INTEGER :: jo, jvar, jk
659      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
660         & zu, &
661         & zv
662
663      !-----------------------------------------------------------------------
664      ! Depending on switches call various observation output routines
665      !-----------------------------------------------------------------------
666
667      IF ( nproftypes > 0 ) THEN
668
669         DO jtype = 1, nproftypes
670
671            IF ( TRIM(cobstypesprof(jtype)) == 'vel' ) THEN
672
673               ! For velocity data, rotate the model velocities to N/S, E/W
674               ! using the compressed data structure.
675               ALLOCATE( &
676                  & zu(profdataqc(jtype)%nvprot(1)), &
677                  & zv(profdataqc(jtype)%nvprot(2))  &
678                  & )
679
680               CALL obs_rotvel( profdataqc(jtype), nn_2dint, zu, zv )
681
682               DO jo = 1, profdataqc(jtype)%nprof
683                  DO jvar = 1, 2
684                     DO jk = profdataqc(jtype)%npvsta(jo,jvar), profdataqc(jtype)%npvend(jo,jvar)
685
686                        IF ( jvar == 1 ) THEN
687                           profdataqc(jtype)%var(jvar)%vmod(jk) = zu(jk)
688                        ELSE
689                           profdataqc(jtype)%var(jvar)%vmod(jk) = zv(jk)
690                        ENDIF
691
692                     END DO
693                  END DO
694               END DO
695
696               DEALLOCATE( zu )
697               DEALLOCATE( zv )
698
699            END IF
700
701            CALL obs_prof_decompress( profdataqc(jtype), &
702               &                      profdata(jtype), .TRUE., numout )
703
704            CALL obs_wri_prof( profdata(jtype) )
705
706         END DO
707
708      ENDIF
709
710      IF ( nsurftypes > 0 ) THEN
711
712         DO jtype = 1, nsurftypes
713
714            CALL obs_surf_decompress( surfdataqc(jtype), &
715               &                      surfdata(jtype), .TRUE., numout )
716
717            CALL obs_wri_surf( surfdata(jtype) )
718
719         END DO
720
721      ENDIF
722
723   END SUBROUTINE dia_obs_wri
724
725   SUBROUTINE dia_obs_dealloc
726      IMPLICIT NONE
727      !!----------------------------------------------------------------------
728      !!                    *** ROUTINE dia_obs_dealloc ***
729      !!
730      !!  ** Purpose : To deallocate data to enable the obs_oper online loop.
731      !!               Specifically: dia_obs_init --> dia_obs --> dia_obs_wri
732      !!
733      !!  ** Method : Clean up various arrays left behind by the obs_oper.
734      !!
735      !!  ** Action :
736      !!
737      !!----------------------------------------------------------------------
738      ! obs_grid deallocation
739      CALL obs_grid_deallocate
740
741      ! diaobs deallocation
742      IF ( nproftypes > 0 ) &
743         &   DEALLOCATE( cobstypesprof, profdata, profdataqc, nvarsprof, nextrprof )
744
745      IF ( nsurftypes > 0 ) &
746         &   DEALLOCATE( cobstypessurf, surfdata, surfdataqc, nvarssurf, nextrsurf, &
747         &               n2dintsurf, zavglamscl, zavgphiscl, lfpindegs, llnightav )
748
749   END SUBROUTINE dia_obs_dealloc
750
751   SUBROUTINE calc_date( kstp, ddobs )
752      !!----------------------------------------------------------------------
753      !!                    ***  ROUTINE calc_date  ***
754      !!         
755      !! ** Purpose : Get date in double precision YYYYMMDD.HHMMSS format
756      !!
757      !! ** Method  : Get date in double precision YYYYMMDD.HHMMSS format
758      !!
759      !! ** Action  : Get date in double precision YYYYMMDD.HHMMSS format
760      !!
761      !! ** Action  : Get initial date in double precision YYYYMMDD.HHMMSS format
762      !!
763      !! History :
764      !!        !  06-03  (K. Mogensen)  Original code
765      !!        !  06-05  (K. Mogensen)  Reformatted
766      !!        !  06-10  (A. Weaver) Cleaning
767      !!        !  06-10  (G. Smith) Calculates initial date the same as method for final date
768      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
769      !!        !  2014-09  (D. Lea) New generic routine now deals with arbitrary initial time of day
770      !!----------------------------------------------------------------------
771      USE phycst, ONLY : &            ! Physical constants
772         & rday
773      USE dom_oce, ONLY : &           ! Ocean space and time domain variables
774         & rdt
775
776      IMPLICIT NONE
777
778      !! * Arguments
779      REAL(KIND=dp), INTENT(OUT) :: ddobs                        ! Date in YYYYMMDD.HHMMSS
780      INTEGER :: kstp
781
782      !! * Local declarations
783      INTEGER :: iyea        ! date - (year, month, day, hour, minute)
784      INTEGER :: imon
785      INTEGER :: iday
786      INTEGER :: ihou
787      INTEGER :: imin
788      INTEGER :: imday       ! Number of days in month.
789      REAL(wp) :: zdayfrc    ! Fraction of day
790
791      INTEGER, DIMENSION(12) ::   imonth_len    !: length in days of the months of the current year
792
793      !!----------------------------------------------------------------------
794      !! Initial date initialization (year, month, day, hour, minute)
795      !!----------------------------------------------------------------------
796      iyea =   ndate0 / 10000
797      imon = ( ndate0 - iyea * 10000 ) / 100
798      iday =   ndate0 - iyea * 10000 - imon * 100
799      ihou =   nn_time0 / 100
800      imin = ( nn_time0 - ihou * 100 ) 
801
802      !!----------------------------------------------------------------------
803      !! Compute number of days + number of hours + min since initial time
804      !!----------------------------------------------------------------------
805      zdayfrc = kstp * rdt / rday
806      zdayfrc = zdayfrc - aint(zdayfrc)
807      imin = imin + int( zdayfrc * 24 * 60 ) 
808      DO WHILE (imin >= 60) 
809        imin=imin-60
810        ihou=ihou+1
811      END DO
812      DO WHILE (ihou >= 24)
813        ihou=ihou-24
814        iday=iday+1
815      END DO
816      iday = iday + kstp * rdt / rday 
817
818      !-----------------------------------------------------------------------
819      ! Convert number of days (iday) into a real date
820      !----------------------------------------------------------------------
821
822      CALL calc_month_len( iyea, imonth_len )
823
824      DO WHILE ( iday > imonth_len(imon) )
825         iday = iday - imonth_len(imon)
826         imon = imon + 1 
827         IF ( imon > 12 ) THEN
828            imon = 1
829            iyea = iyea + 1
830            CALL calc_month_len( iyea, imonth_len )  ! update month lengths
831         ENDIF
832      END DO
833
834      !----------------------------------------------------------------------
835      ! Convert it into YYYYMMDD.HHMMSS format.
836      !----------------------------------------------------------------------
837      ddobs = iyea * 10000_dp + imon * 100_dp + &
838         &    iday + ihou * 0.01_dp + imin * 0.0001_dp
839
840   END SUBROUTINE calc_date
841
842   SUBROUTINE ini_date( ddobsini )
843      !!----------------------------------------------------------------------
844      !!                    ***  ROUTINE ini_date  ***
845      !!         
846      !! ** Purpose : Get initial date in double precision YYYYMMDD.HHMMSS format
847      !!
848      !! ** Method  :
849      !!
850      !! ** Action  :
851      !!
852      !! History :
853      !!        !  06-03  (K. Mogensen)  Original code
854      !!        !  06-05  (K. Mogensen)  Reformatted
855      !!        !  06-10  (A. Weaver) Cleaning
856      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
857      !!        !  2014-09  (D. Lea) Change to call generic routine calc_date
858      !!----------------------------------------------------------------------
859
860      IMPLICIT NONE
861
862      !! * Arguments
863      REAL(KIND=dp), INTENT(OUT) :: ddobsini                   ! Initial date in YYYYMMDD.HHMMSS
864
865      CALL calc_date( nit000 - 1, ddobsini )
866
867   END SUBROUTINE ini_date
868
869   SUBROUTINE fin_date( ddobsfin )
870      !!----------------------------------------------------------------------
871      !!                    ***  ROUTINE fin_date  ***
872      !!         
873      !! ** Purpose : Get final date in double precision YYYYMMDD.HHMMSS format
874      !!
875      !! ** Method  :
876      !!
877      !! ** Action  :
878      !!
879      !! History :
880      !!        !  06-03  (K. Mogensen)  Original code
881      !!        !  06-05  (K. Mogensen)  Reformatted
882      !!        !  06-10  (A. Weaver) Cleaning
883      !!        !  10-05  (D. Lea) Update to month length calculation for NEMO vn3.2
884      !!        !  2014-09  (D. Lea) Change to call generic routine calc_date
885      !!----------------------------------------------------------------------
886
887      IMPLICIT NONE
888
889      !! * Arguments
890      REAL(dp), INTENT(OUT) :: ddobsfin ! Final date in YYYYMMDD.HHMMSS
891
892      CALL calc_date( nitend, ddobsfin )
893
894   END SUBROUTINE fin_date
895   
896    SUBROUTINE obs_settypefiles( ntypes, jpmaxnfiles, jtype, ctypein, &
897       &                         cfilestype, ifiles, cobstypes, cfiles )
898
899    INTEGER, INTENT(IN) :: ntypes      ! Total number of obs types
900    INTEGER, INTENT(IN) :: jpmaxnfiles ! Maximum number of files allowed for each type
901    INTEGER, INTENT(IN) :: jtype       ! Index of the current type of obs
902    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: &
903       &                   ifiles      ! Out appended number of files for this type
904
905    CHARACTER(len=6), INTENT(IN) :: ctypein 
906    CHARACTER(len=128), DIMENSION(jpmaxnfiles), INTENT(IN) :: &
907       &                   cfilestype  ! In list of files for this obs type
908    CHARACTER(len=6), DIMENSION(ntypes), INTENT(INOUT) :: &
909       &                   cobstypes   ! Out appended list of obs types
910    CHARACTER(len=128), DIMENSION(ntypes, jpmaxnfiles), INTENT(INOUT) :: &
911       &                   cfiles      ! Out appended list of files for all types
912
913    !Local variables
914    INTEGER :: jfile
915
916    cfiles(jtype,:) = cfilestype(:)
917    cobstypes(jtype) = ctypein
918    ifiles(jtype) = 0
919    DO jfile = 1, jpmaxnfiles
920       IF ( trim(cfiles(jtype,jfile)) /= '' ) &
921                 ifiles(jtype) = ifiles(jtype) + 1
922    END DO
923
924    IF ( ifiles(jtype) == 0 ) THEN
925         CALL ctl_stop( 'Logical for observation type '//TRIM(ctypein)//   &
926            &           ' set to true but no files available to read' )
927    ENDIF
928
929    IF(lwp) THEN   
930       WRITE(numout,*) '             '//cobstypes(jtype)//' input observation file names:'
931       DO jfile = 1, ifiles(jtype)
932          WRITE(numout,*) '                '//TRIM(cfiles(jtype,jfile))
933       END DO
934    ENDIF
935
936    END SUBROUTINE obs_settypefiles
937
938    SUBROUTINE obs_setinterpopts( ntypes, jtype, ctypein,             &
939               &                  n2dint_default, n2dint_type,        &
940               &                  zavglamscl_type, zavgphiscl_type,   &
941               &                  lfp_indegs_type, lavnight_type,     &
942               &                  n2dint, zavglamscl, zavgphiscl,     &
943               &                  lfpindegs, lavnight )
944
945    INTEGER, INTENT(IN)  :: ntypes             ! Total number of obs types
946    INTEGER, INTENT(IN)  :: jtype              ! Index of the current type of obs
947    INTEGER, INTENT(IN)  :: n2dint_default     ! Default option for interpolation type
948    INTEGER, INTENT(IN)  :: n2dint_type        ! Option for interpolation type
949    REAL(wp), INTENT(IN) :: &
950       &                    zavglamscl_type, & !E/W diameter of obs footprint for this type
951       &                    zavgphiscl_type    !N/S diameter of obs footprint for this type
952    LOGICAL, INTENT(IN)  :: lfp_indegs_type    !T=> footprint in degrees, F=> in metres
953    LOGICAL, INTENT(IN)  :: lavnight_type      !T=> obs represent night time average
954    CHARACTER(len=6), INTENT(IN) :: ctypein 
955
956    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: &
957       &                    n2dint 
958    REAL(wp), DIMENSION(ntypes), INTENT(INOUT) :: &
959       &                    zavglamscl, zavgphiscl
960    LOGICAL, DIMENSION(ntypes), INTENT(INOUT) :: &
961       &                    lfpindegs, lavnight
962
963    lavnight(jtype) = lavnight_type
964
965    IF ( (n2dint_type >= 1) .AND. (n2dint_type <= 6) ) THEN
966       n2dint(jtype) = n2dint_type
967    ELSE
968       n2dint(jtype) = n2dint_default
969    ENDIF
970
971    ! For averaging observation footprints set options for size of footprint
972    IF ( (n2dint(jtype) > 4) .AND. (n2dint(jtype) <= 6) ) THEN
973       IF ( zavglamscl_type > 0._wp ) THEN
974          zavglamscl(jtype) = zavglamscl_type
975       ELSE
976          CALL ctl_stop( 'Incorrect value set for averaging footprint '// &
977                         'scale (zavglamscl) for observation type '//TRIM(ctypein) )     
978       ENDIF
979
980       IF ( zavgphiscl_type > 0._wp ) THEN
981          zavgphiscl(jtype) = zavgphiscl_type
982       ELSE
983          CALL ctl_stop( 'Incorrect value set for averaging footprint '// &
984                         'scale (zavgphiscl) for observation type '//TRIM(ctypein) )     
985       ENDIF
986
987       lfpindegs(jtype) = lfp_indegs_type 
988
989    ENDIF
990
991    ! Write out info
992    IF(lwp) THEN
993       IF ( n2dint(jtype) <= 4 ) THEN
994          WRITE(numout,*) '             '//TRIM(ctypein)// &
995             &            ' model counterparts will be interpolated horizontally'
996       ELSE IF ( n2dint(jtype) <= 6 ) THEN
997          WRITE(numout,*) '             '//TRIM(ctypein)// &
998             &            ' model counterparts will be averaged horizontally'
999          WRITE(numout,*) '             '//'    with E/W scale: ',zavglamscl(jtype)
1000          WRITE(numout,*) '             '//'    with N/S scale: ',zavgphiscl(jtype)
1001          IF ( lfpindegs(jtype) ) THEN
1002              WRITE(numout,*) '             '//'    (in degrees)'
1003          ELSE
1004              WRITE(numout,*) '             '//'    (in metres)'
1005          ENDIF
1006       ENDIF
1007    ENDIF
1008
1009    END SUBROUTINE obs_setinterpopts
1010
1011END MODULE diaobs
Note: See TracBrowser for help on using the repository browser.