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.
obs_write.F90 in NEMO/branches/2020/dev_r12563_ASINTER-06_ABL_improvement/src/OCE/OBS – NEMO

source: NEMO/branches/2020/dev_r12563_ASINTER-06_ABL_improvement/src/OCE/OBS/obs_write.F90 @ 12808

Last change on this file since 12808 was 12377, checked in by acc, 4 years ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

  • Property svn:keywords set to Id
File size: 22.2 KB
Line 
1MODULE obs_write
2   !!======================================================================
3   !!                       ***  MODULE obs_write   ***
4   !! Observation diagnosticss: Write observation related diagnostics
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   obs_wri_prof   : Write profile observations in feedback format
9   !!   obs_wri_surf   : Write surface observations in feedback format
10   !!   obs_wri_stats  : Print basic statistics on the data being written out
11   !!----------------------------------------------------------------------
12
13   !! * Modules used
14   USE par_kind, ONLY : &   ! Precision variables
15      & wp
16   USE in_out_manager       ! I/O manager
17   USE dom_oce              ! Ocean space and time domain variables
18   USE obs_types            ! Observation type integer to character translation
19   USE julian, ONLY : &         ! Julian date routines
20      & greg2jul
21   USE obs_utils, ONLY : &  ! Observation operator utility functions
22      & chkerr
23   USE obs_profiles_def     ! Type definitions for profiles
24   USE obs_surf_def         ! Type defintions for surface observations
25   USE obs_fbm              ! Observation feedback I/O
26   USE obs_grid             ! Grid tools
27   USE obs_conv             ! Conversion between units
28   USE obs_const
29   USE obs_mpp              ! MPP support routines for observation diagnostics
30   USE lib_mpp        ! MPP routines
31
32   IMPLICIT NONE
33
34   !! * Routine accessibility
35   PRIVATE
36   PUBLIC obs_wri_prof, &    ! Write profile observation files
37      &   obs_wri_surf, &    ! Write surface observation files
38      &   obswriinfo
39   
40   TYPE obswriinfo
41      INTEGER :: inum
42      INTEGER, POINTER, DIMENSION(:) :: ipoint
43      CHARACTER(len=ilenname), POINTER, DIMENSION(:) :: cdname
44      CHARACTER(len=ilenlong), POINTER, DIMENSION(:,:) :: cdlong
45      CHARACTER(len=ilenunit), POINTER, DIMENSION(:,:) :: cdunit
46   END TYPE obswriinfo
47
48   !! * Substitutions
49#  include "do_loop_substitute.h90"
50   !!----------------------------------------------------------------------
51   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
52   !! $Id$
53   !! Software governed by the CeCILL license (see ./LICENSE)
54   !!----------------------------------------------------------------------
55
56CONTAINS
57
58   SUBROUTINE obs_wri_prof( profdata, padd, pext )
59      !!-----------------------------------------------------------------------
60      !!
61      !!                     *** ROUTINE obs_wri_prof  ***
62      !!
63      !! ** Purpose : Write profile feedback files
64      !!
65      !! ** Method  : NetCDF
66      !!
67      !! ** Action  :
68      !!
69      !! History :
70      !!      ! 06-04  (A. Vidard) Original
71      !!      ! 06-04  (A. Vidard) Reformatted
72      !!      ! 06-10  (A. Weaver) Cleanup
73      !!      ! 07-01  (K. Mogensen) Use profile data types
74      !!      ! 07-03  (K. Mogensen) General handling of profiles
75      !!      ! 09-01  (K. Mogensen) New feedback format
76      !!      ! 15-02  (M. Martin) Combined routine for writing profiles
77      !!-----------------------------------------------------------------------
78
79      !! * Arguments
80      TYPE(obs_prof), INTENT(INOUT) :: profdata      ! Full set of profile data
81      TYPE(obswriinfo), OPTIONAL :: padd             ! Additional info for each variable
82      TYPE(obswriinfo), OPTIONAL :: pext             ! Extra info
83
84      !! * Local declarations
85      TYPE(obfbdata) :: fbdata
86      CHARACTER(LEN=40) :: clfname
87      CHARACTER(LEN=10) :: clfiletype
88      INTEGER :: ilevel
89      INTEGER :: jvar
90      INTEGER :: jo
91      INTEGER :: jk
92      INTEGER :: ik
93      INTEGER :: ja
94      INTEGER :: je
95      INTEGER :: iadd
96      INTEGER :: iext
97      REAL(wp) :: zpres
98
99      IF ( PRESENT( padd ) ) THEN
100         iadd = padd%inum
101      ELSE
102         iadd = 0
103      ENDIF
104
105      IF ( PRESENT( pext ) ) THEN
106         iext = pext%inum
107      ELSE
108         iext = 0
109      ENDIF
110
111      CALL init_obfbdata( fbdata )
112
113      ! Find maximum level
114      ilevel = 0
115      DO jvar = 1, 2
116         ilevel = MAX( ilevel, MAXVAL( profdata%var(jvar)%nvlidx(:) ) )
117      END DO
118
119      SELECT CASE ( TRIM(profdata%cvars(1)) )
120      CASE('POTM')
121
122         clfiletype='profb'
123         CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, &
124            &                 1 + iadd, 1 + iext, .TRUE. )
125         fbdata%cname(1)      = profdata%cvars(1)
126         fbdata%cname(2)      = profdata%cvars(2)
127         fbdata%coblong(1)    = 'Potential temperature'
128         fbdata%coblong(2)    = 'Practical salinity'
129         fbdata%cobunit(1)    = 'Degrees centigrade'
130         fbdata%cobunit(2)    = 'PSU'
131         fbdata%cextname(1)   = 'TEMP'
132         fbdata%cextlong(1)   = 'Insitu temperature'
133         fbdata%cextunit(1)   = 'Degrees centigrade'
134         fbdata%caddlong(1,1) = 'Model interpolated potential temperature'
135         fbdata%caddlong(1,2) = 'Model interpolated practical salinity'
136         fbdata%caddunit(1,1) = 'Degrees centigrade'
137         fbdata%caddunit(1,2) = 'PSU'
138         fbdata%cgrid(:)      = 'T'
139         DO je = 1, iext
140            fbdata%cextname(1+je) = pext%cdname(je)
141            fbdata%cextlong(1+je) = pext%cdlong(je,1)
142            fbdata%cextunit(1+je) = pext%cdunit(je,1)
143         END DO
144         DO ja = 1, iadd
145            fbdata%caddname(1+ja) = padd%cdname(ja)
146            DO jvar = 1, 2
147               fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
148               fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
149            END DO
150         END DO
151
152      CASE('UVEL')
153
154         clfiletype='velfb'
155         CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, 1, 0, .TRUE. )
156         fbdata%cname(1)      = profdata%cvars(1)
157         fbdata%cname(2)      = profdata%cvars(2)
158         fbdata%coblong(1)    = 'Zonal velocity'
159         fbdata%coblong(2)    = 'Meridional velocity'
160         fbdata%cobunit(1)    = 'm/s'
161         fbdata%cobunit(2)    = 'm/s'
162         DO je = 1, iext
163            fbdata%cextname(je) = pext%cdname(je)
164            fbdata%cextlong(je) = pext%cdlong(je,1)
165            fbdata%cextunit(je) = pext%cdunit(je,1)
166         END DO
167         fbdata%caddlong(1,1) = 'Model interpolated zonal velocity'
168         fbdata%caddlong(1,2) = 'Model interpolated meridional velocity'
169         fbdata%caddunit(1,1) = 'm/s'
170         fbdata%caddunit(1,2) = 'm/s'
171         fbdata%cgrid(1)      = 'U' 
172         fbdata%cgrid(2)      = 'V'
173         DO ja = 1, iadd
174            fbdata%caddname(1+ja) = padd%cdname(ja)
175            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
176            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
177         END DO
178
179      END SELECT
180
181      fbdata%caddname(1)   = 'Hx'
182
183      WRITE(clfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
184
185      IF(lwp) THEN
186         WRITE(numout,*)
187         WRITE(numout,*)'obs_wri_prof :'
188         WRITE(numout,*)'~~~~~~~~~~~~~'
189         WRITE(numout,*)'Writing '//TRIM(clfiletype)//' feedback file : ',TRIM(clfname)
190      ENDIF
191
192      ! Transform obs_prof data structure into obfb data structure
193      fbdata%cdjuldref = '19500101000000'
194      DO jo = 1, profdata%nprof
195         fbdata%plam(jo)      = profdata%rlam(jo)
196         fbdata%pphi(jo)      = profdata%rphi(jo)
197         WRITE(fbdata%cdtyp(jo),'(I4)') profdata%ntyp(jo)
198         fbdata%ivqc(jo,:)    = profdata%ivqc(jo,:)
199         fbdata%ivqcf(:,jo,:) = profdata%ivqcf(:,jo,:)
200         IF ( profdata%nqc(jo) > 255 ) THEN
201            fbdata%ioqc(jo)    = IBSET(profdata%nqc(jo),2)
202            fbdata%ioqcf(1,jo) = profdata%nqcf(1,jo)
203            fbdata%ioqcf(2,jo) = profdata%nqc(jo)
204         ELSE
205            fbdata%ioqc(jo)    = profdata%nqc(jo)
206            fbdata%ioqcf(:,jo) = profdata%nqcf(:,jo)
207         ENDIF
208         fbdata%ipqc(jo)      = profdata%ipqc(jo)
209         fbdata%ipqcf(:,jo)   = profdata%ipqcf(:,jo)
210         fbdata%itqc(jo)      = profdata%itqc(jo)
211         fbdata%itqcf(:,jo)   = profdata%itqcf(:,jo)
212         fbdata%cdwmo(jo)     = profdata%cwmo(jo)
213         fbdata%kindex(jo)    = profdata%npfil(jo)
214         DO jvar = 1, profdata%nvar
215            IF (ln_grid_global) THEN
216               fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar)
217               fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar)
218            ELSE
219               fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar))
220               fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar))
221            ENDIF
222         END DO
223         CALL greg2jul( 0, &
224            &           profdata%nmin(jo), &
225            &           profdata%nhou(jo), &
226            &           profdata%nday(jo), &
227            &           profdata%nmon(jo), &
228            &           profdata%nyea(jo), &
229            &           fbdata%ptim(jo),   &
230            &           krefdate = 19500101 )
231         ! Reform the profiles arrays for output
232         DO jvar = 1, 2
233            DO jk = profdata%npvsta(jo,jvar), profdata%npvend(jo,jvar)
234               ik = profdata%var(jvar)%nvlidx(jk)
235               fbdata%padd(ik,jo,1,jvar) = profdata%var(jvar)%vmod(jk)
236               fbdata%pob(ik,jo,jvar)    = profdata%var(jvar)%vobs(jk)
237               fbdata%pdep(ik,jo)        = profdata%var(jvar)%vdep(jk)
238               fbdata%idqc(ik,jo)        = profdata%var(jvar)%idqc(jk)
239               fbdata%idqcf(:,ik,jo)     = profdata%var(jvar)%idqcf(:,jk)
240               IF ( profdata%var(jvar)%nvqc(jk) > 255 ) THEN
241                  fbdata%ivlqc(ik,jo,jvar) = IBSET(profdata%var(jvar)%nvqc(jk),2)
242                  fbdata%ivlqcf(1,ik,jo,jvar) = profdata%var(jvar)%nvqcf(1,jk)
243!$AGRIF_DO_NOT_TREAT
244                  fbdata%ivlqcf(2,ik,jo,jvar) = IAND(profdata%var(jvar)%nvqc(jk),b'0000000011111111')
245!$AGRIF_END_DO_NOT_TREAT
246               ELSE
247                  fbdata%ivlqc(ik,jo,jvar) = profdata%var(jvar)%nvqc(jk)
248                  fbdata%ivlqcf(:,ik,jo,jvar) = profdata%var(jvar)%nvqcf(:,jk)
249               ENDIF
250               fbdata%iobsk(ik,jo,jvar)  = profdata%var(jvar)%mvk(jk)
251               DO ja = 1, iadd
252                  fbdata%padd(ik,jo,1+ja,jvar) = &
253                     & profdata%var(jvar)%vext(jk,padd%ipoint(ja))
254               END DO
255               DO je = 1, iext
256                  fbdata%pext(ik,jo,1+je) = &
257                     & profdata%var(jvar)%vext(jk,pext%ipoint(je))
258               END DO
259               IF ( ( jvar == 1 ) .AND. &
260                  & ( TRIM(profdata%cvars(1)) == 'POTM' ) ) THEN
261                  fbdata%pext(ik,jo,1) = profdata%var(jvar)%vext(jk,1)
262               ENDIF
263            END DO
264         END DO
265      END DO
266
267      IF ( TRIM(profdata%cvars(1)) == 'POTM' ) THEN
268         ! Convert insitu temperature to potential temperature using the model
269         ! salinity if no potential temperature
270         DO jo = 1, fbdata%nobs
271            IF ( fbdata%pphi(jo) < 9999.0 ) THEN
272               DO jk = 1, fbdata%nlev
273                  IF ( ( fbdata%pob(jk,jo,1) >= 9999.0 ) .AND. &
274                     & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
275                     & ( fbdata%padd(jk,jo,1,2) < 9999.0 ) .AND. &
276                     & ( fbdata%pext(jk,jo,1) < 9999.0 ) ) THEN
277                     zpres = dep_to_p( REAL(fbdata%pdep(jk,jo),wp), &
278                        &              REAL(fbdata%pphi(jo),wp) )
279                     fbdata%pob(jk,jo,1) = potemp( &
280                        &                     REAL(fbdata%padd(jk,jo,1,2), wp),  &
281                        &                     REAL(fbdata%pext(jk,jo,1), wp), &
282                        &                     zpres, 0.0_wp )
283                  ENDIF
284               END DO
285            ENDIF
286         END DO
287      ENDIF
288
289      ! Write the obfbdata structure
290      CALL write_obfbdata( clfname, fbdata )
291
292      ! Output some basic statistics
293      CALL obs_wri_stats( fbdata )
294
295      CALL dealloc_obfbdata( fbdata )
296
297   END SUBROUTINE obs_wri_prof
298
299   SUBROUTINE obs_wri_surf( surfdata, padd, pext )
300      !!-----------------------------------------------------------------------
301      !!
302      !!                     *** ROUTINE obs_wri_surf  ***
303      !!
304      !! ** Purpose : Write surface observation files
305      !!
306      !! ** Method  : NetCDF
307      !!
308      !! ** Action  :
309      !!
310      !!      ! 07-03  (K. Mogensen) Original
311      !!      ! 09-01  (K. Mogensen) New feedback format.
312      !!      ! 15-02  (M. Martin) Combined surface writing routine.
313      !!-----------------------------------------------------------------------
314
315      !! * Modules used
316      IMPLICIT NONE
317
318      !! * Arguments
319      TYPE(obs_surf), INTENT(INOUT) :: surfdata         ! Full set of surface data
320      TYPE(obswriinfo), OPTIONAL :: padd               ! Additional info for each variable
321      TYPE(obswriinfo), OPTIONAL :: pext               ! Extra info
322
323      !! * Local declarations
324      TYPE(obfbdata) :: fbdata
325      CHARACTER(LEN=40) :: clfname         ! netCDF filename
326      CHARACTER(LEN=10) :: clfiletype
327      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_surf'
328      INTEGER :: jo
329      INTEGER :: ja
330      INTEGER :: je
331      INTEGER :: iadd
332      INTEGER :: iext
333
334      IF ( PRESENT( padd ) ) THEN
335         iadd = padd%inum
336      ELSE
337         iadd = 0
338      ENDIF
339
340      IF ( PRESENT( pext ) ) THEN
341         iext = pext%inum
342      ELSE
343         iext = 0
344      ENDIF
345
346      CALL init_obfbdata( fbdata )
347
348      SELECT CASE ( TRIM(surfdata%cvars(1)) )
349      CASE('SLA')
350
351         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
352            &                 2 + iadd, 1 + iext, .TRUE. )
353
354         clfiletype = 'slafb'
355         fbdata%cname(1)      = surfdata%cvars(1)
356         fbdata%coblong(1)    = 'Sea level anomaly'
357         fbdata%cobunit(1)    = 'Metres'
358         fbdata%cextname(1)   = 'MDT'
359         fbdata%cextlong(1)   = 'Mean dynamic topography'
360         fbdata%cextunit(1)   = 'Metres'
361         DO je = 1, iext
362            fbdata%cextname(je) = pext%cdname(je)
363            fbdata%cextlong(je) = pext%cdlong(je,1)
364            fbdata%cextunit(je) = pext%cdunit(je,1)
365         END DO
366         fbdata%caddlong(1,1) = 'Model interpolated SSH - MDT'
367         fbdata%caddunit(1,1) = 'Metres' 
368         fbdata%caddname(2)   = 'SSH'
369         fbdata%caddlong(2,1) = 'Model Sea surface height'
370         fbdata%caddunit(2,1) = 'Metres'
371         fbdata%cgrid(1)      = 'T'
372         DO ja = 1, iadd
373            fbdata%caddname(2+ja) = padd%cdname(ja)
374            fbdata%caddlong(2+ja,1) = padd%cdlong(ja,1)
375            fbdata%caddunit(2+ja,1) = padd%cdunit(ja,1)
376         END DO
377
378      CASE('SST')
379
380         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
381            &                 1 + iadd, iext, .TRUE. )
382
383         clfiletype = 'sstfb'
384         fbdata%cname(1)      = surfdata%cvars(1)
385         fbdata%coblong(1)    = 'Sea surface temperature'
386         fbdata%cobunit(1)    = 'Degree centigrade'
387         DO je = 1, iext
388            fbdata%cextname(je) = pext%cdname(je)
389            fbdata%cextlong(je) = pext%cdlong(je,1)
390            fbdata%cextunit(je) = pext%cdunit(je,1)
391         END DO
392         fbdata%caddlong(1,1) = 'Model interpolated SST'
393         fbdata%caddunit(1,1) = 'Degree centigrade'
394         fbdata%cgrid(1)      = 'T'
395         DO ja = 1, iadd
396            fbdata%caddname(1+ja) = padd%cdname(ja)
397            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
398            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
399         END DO
400
401      CASE('ICECONC')
402
403         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
404            &                 1 + iadd, iext, .TRUE. )
405
406         clfiletype = 'sicfb'
407         fbdata%cname(1)      = surfdata%cvars(1)
408         fbdata%coblong(1)    = 'Sea ice'
409         fbdata%cobunit(1)    = 'Fraction'
410         DO je = 1, iext
411            fbdata%cextname(je) = pext%cdname(je)
412            fbdata%cextlong(je) = pext%cdlong(je,1)
413            fbdata%cextunit(je) = pext%cdunit(je,1)
414         END DO
415         fbdata%caddlong(1,1) = 'Model interpolated ICE'
416         fbdata%caddunit(1,1) = 'Fraction'
417         fbdata%cgrid(1)      = 'T'
418         DO ja = 1, iadd
419            fbdata%caddname(1+ja) = padd%cdname(ja)
420            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
421            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
422         END DO
423
424      CASE('SSS')
425
426         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
427            &                 1 + iadd, iext, .TRUE. )
428
429         clfiletype = 'sssfb'
430         fbdata%cname(1)      = surfdata%cvars(1)
431         fbdata%coblong(1)    = 'Sea surface salinity'
432         fbdata%cobunit(1)    = 'psu'
433         DO je = 1, iext
434            fbdata%cextname(je) = pext%cdname(je)
435            fbdata%cextlong(je) = pext%cdlong(je,1)
436            fbdata%cextunit(je) = pext%cdunit(je,1)
437         END DO
438         fbdata%caddlong(1,1) = 'Model interpolated SSS'
439         fbdata%caddunit(1,1) = 'psu'
440         fbdata%cgrid(1)      = 'T'
441         DO ja = 1, iadd
442            fbdata%caddname(1+ja) = padd%cdname(ja)
443            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
444            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
445         END DO
446
447      CASE DEFAULT
448
449         CALL ctl_stop( 'Unknown observation type '//TRIM(surfdata%cvars(1))//' in obs_wri_surf' )
450
451      END SELECT
452
453      fbdata%caddname(1)   = 'Hx'
454
455      WRITE(clfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
456
457      IF(lwp) THEN
458         WRITE(numout,*)
459         WRITE(numout,*)'obs_wri_surf :'
460         WRITE(numout,*)'~~~~~~~~~~~~~'
461         WRITE(numout,*)'Writing '//TRIM(surfdata%cvars(1))//' feedback file : ',TRIM(clfname)
462      ENDIF
463
464      ! Transform surf data structure into obfbdata structure
465      fbdata%cdjuldref = '19500101000000'
466      DO jo = 1, surfdata%nsurf
467         fbdata%plam(jo)      = surfdata%rlam(jo)
468         fbdata%pphi(jo)      = surfdata%rphi(jo)
469         WRITE(fbdata%cdtyp(jo),'(I4)') surfdata%ntyp(jo)
470         fbdata%ivqc(jo,:)    = 0
471         fbdata%ivqcf(:,jo,:) = 0
472         IF ( surfdata%nqc(jo) > 255 ) THEN
473            fbdata%ioqc(jo)    = 4
474            fbdata%ioqcf(1,jo) = 0
475!$AGRIF_DO_NOT_TREAT
476            fbdata%ioqcf(2,jo) = IAND(surfdata%nqc(jo),b'0000000011111111')
477!$AGRIF_END_DO_NOT_TREAT
478         ELSE
479            fbdata%ioqc(jo)    = surfdata%nqc(jo)
480            fbdata%ioqcf(:,jo) = 0
481         ENDIF
482         fbdata%ipqc(jo)      = 0
483         fbdata%ipqcf(:,jo)   = 0
484         fbdata%itqc(jo)      = 0
485         fbdata%itqcf(:,jo)   = 0
486         fbdata%cdwmo(jo)     = surfdata%cwmo(jo)
487         fbdata%kindex(jo)    = surfdata%nsfil(jo)
488         IF (ln_grid_global) THEN
489            fbdata%iobsi(jo,1) = surfdata%mi(jo)
490            fbdata%iobsj(jo,1) = surfdata%mj(jo)
491         ELSE
492            fbdata%iobsi(jo,1) = mig(surfdata%mi(jo))
493            fbdata%iobsj(jo,1) = mjg(surfdata%mj(jo))
494         ENDIF
495         CALL greg2jul( 0, &
496            &           surfdata%nmin(jo), &
497            &           surfdata%nhou(jo), &
498            &           surfdata%nday(jo), &
499            &           surfdata%nmon(jo), &
500            &           surfdata%nyea(jo), &
501            &           fbdata%ptim(jo),   &
502            &           krefdate = 19500101 )
503         fbdata%padd(1,jo,1,1) = surfdata%rmod(jo,1)
504         IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%padd(1,jo,2,1) = surfdata%rext(jo,1)
505         fbdata%pob(1,jo,1)    = surfdata%robs(jo,1) 
506         fbdata%pdep(1,jo)     = 0.0
507         fbdata%idqc(1,jo)     = 0
508         fbdata%idqcf(:,1,jo)  = 0
509         IF ( surfdata%nqc(jo) > 255 ) THEN
510            fbdata%ivqc(jo,1)       = 4
511            fbdata%ivlqc(1,jo,1)    = 4
512            fbdata%ivlqcf(1,1,jo,1) = 0
513!$AGRIF_DO_NOT_TREAT
514            fbdata%ivlqcf(2,1,jo,1) = IAND(surfdata%nqc(jo),b'0000000011111111')
515!$AGRIF_END_DO_NOT_TREAT
516         ELSE
517            fbdata%ivqc(jo,1)       = surfdata%nqc(jo)
518            fbdata%ivlqc(1,jo,1)    = surfdata%nqc(jo)
519            fbdata%ivlqcf(:,1,jo,1) = 0
520         ENDIF
521         fbdata%iobsk(1,jo,1)  = 0
522         IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%pext(1,jo,1) = surfdata%rext(jo,2)
523         DO ja = 1, iadd
524            fbdata%padd(1,jo,2+ja,1) = &
525               & surfdata%rext(jo,padd%ipoint(ja))
526         END DO
527         DO je = 1, iext
528            fbdata%pext(1,jo,1+je) = &
529               & surfdata%rext(jo,pext%ipoint(je))
530         END DO
531      END DO
532
533      ! Write the obfbdata structure
534      CALL write_obfbdata( clfname, fbdata )
535
536      ! Output some basic statistics
537      CALL obs_wri_stats( fbdata )
538
539      CALL dealloc_obfbdata( fbdata )
540
541   END SUBROUTINE obs_wri_surf
542
543   SUBROUTINE obs_wri_stats( fbdata )
544      !!-----------------------------------------------------------------------
545      !!
546      !!                     *** ROUTINE obs_wri_stats  ***
547      !!
548      !! ** Purpose : Output some basic statistics of the data being written out
549      !!
550      !! ** Method  :
551      !!
552      !! ** Action  :
553      !!
554      !!      ! 2014-08  (D. Lea) Initial version
555      !!-----------------------------------------------------------------------
556
557      !! * Arguments
558      TYPE(obfbdata) :: fbdata
559
560      !! * Local declarations
561      INTEGER :: jvar
562      INTEGER :: jo
563      INTEGER :: jk
564      INTEGER :: inumgoodobs
565      INTEGER :: inumgoodobsmpp
566      REAL(wp) :: zsumx
567      REAL(wp) :: zsumx2
568      REAL(wp) :: zomb
569     
570
571      IF (lwp) THEN
572         WRITE(numout,*) ''
573         WRITE(numout,*) 'obs_wri_stats :'
574         WRITE(numout,*) '~~~~~~~~~~~~~~~'
575      ENDIF
576
577      DO jvar = 1, fbdata%nvar
578         zsumx=0.0_wp
579         zsumx2=0.0_wp
580         inumgoodobs=0
581         DO jo = 1, fbdata%nobs
582            DO jk = 1, fbdata%nlev
583               IF ( ( fbdata%pob(jk,jo,jvar) < 9999.0 ) .AND. &
584                  & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
585                  & ( fbdata%padd(jk,jo,1,jvar) < 9999.0 ) ) THEN
586
587                  zomb=fbdata%pob(jk, jo, jvar)-fbdata%padd(jk, jo, 1, jvar)
588                  zsumx=zsumx+zomb
589                  zsumx2=zsumx2+zomb**2
590                  inumgoodobs=inumgoodobs+1
591               ENDIF
592            ENDDO
593         ENDDO
594
595         CALL obs_mpp_sum_integer( inumgoodobs, inumgoodobsmpp )
596         CALL mpp_sum('obs_write', zsumx)
597         CALL mpp_sum('obs_write', zsumx2)
598
599         IF (lwp) THEN
600            WRITE(numout,*) 'Type: ',fbdata%cname(jvar),'  Total number of good observations: ',inumgoodobsmpp 
601            WRITE(numout,*) 'Overall mean obs minus model of the good observations: ',zsumx/inumgoodobsmpp
602            WRITE(numout,*) 'Overall RMS obs minus model of the good observations: ',sqrt( zsumx2/inumgoodobsmpp )
603            WRITE(numout,*) ''
604         ENDIF
605
606      ENDDO
607
608   END SUBROUTINE obs_wri_stats
609
610END MODULE obs_write
Note: See TracBrowser for help on using the repository browser.