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

source: NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_write.F90 @ 15187

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

Update treatment of SLA and POTM additional/extra variables.

File size: 20.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   !!----------------------------------------------------------------------
49   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
50   !! $Id$
51   !! Software governed by the CeCILL license (see ./LICENSE)
52   !!----------------------------------------------------------------------
53
54CONTAINS
55
56   SUBROUTINE obs_wri_prof( profdata, clfiletype, padd, pext )
57      !!-----------------------------------------------------------------------
58      !!
59      !!                     *** ROUTINE obs_wri_prof  ***
60      !!
61      !! ** Purpose : Write profile feedback files
62      !!
63      !! ** Method  : NetCDF
64      !!
65      !! ** Action  :
66      !!
67      !! History :
68      !!      ! 06-04  (A. Vidard) Original
69      !!      ! 06-04  (A. Vidard) Reformatted
70      !!      ! 06-10  (A. Weaver) Cleanup
71      !!      ! 07-01  (K. Mogensen) Use profile data types
72      !!      ! 07-03  (K. Mogensen) General handling of profiles
73      !!      ! 09-01  (K. Mogensen) New feedback format
74      !!      ! 15-02  (M. Martin) Combined routine for writing profiles
75      !!-----------------------------------------------------------------------
76
77      !! * Arguments
78      TYPE(obs_prof), INTENT(INOUT) :: profdata      ! Full set of profile data
79      CHARACTER(LEN=25), INTENT(IN) :: clfiletype    ! Base name for file name
80      TYPE(obswriinfo), OPTIONAL    :: padd          ! Additional info for each variable
81      TYPE(obswriinfo), OPTIONAL    :: pext          ! Extra info
82
83      !! * Local declarations
84      TYPE(obfbdata) :: fbdata
85      CHARACTER(LEN=40) :: clfname
86      CHARACTER(LEN=ilenlong) :: cllongname  ! Long name of variable
87      CHARACTER(LEN=ilenunit) :: clunits     ! Units of variable
88      CHARACTER(LEN=ilengrid) :: clgrid      ! Grid of variable
89      CHARACTER(LEN=12) :: clfmt            ! writing format
90      INTEGER :: idg                        ! number of digits
91      INTEGER :: ilevel
92      INTEGER :: jvar
93      INTEGER :: jvar2
94      INTEGER :: jsal
95      INTEGER :: jo
96      INTEGER :: jk
97      INTEGER :: ik
98      INTEGER :: ja
99      INTEGER :: je
100      INTEGER :: iadd
101      INTEGER :: iext
102      REAL(wp) :: zpres
103
104      IF ( PRESENT( padd ) ) THEN
105         iadd = padd%inum
106      ELSE
107         iadd = 0
108      ENDIF
109
110      IF ( PRESENT( pext ) ) THEN
111         iext = pext%inum
112      ELSE
113         iext = 0
114      ENDIF
115
116      CALL init_obfbdata( fbdata )
117
118      ! Find maximum level
119      ilevel = 0
120      DO jvar = 1, profdata%nvar
121         ilevel = MAX( ilevel, MAXVAL( profdata%var(jvar)%nvlidx(:) ) )
122      END DO
123
124      CALL alloc_obfbdata( fbdata, profdata%nvar, profdata%nprof, ilevel, &
125            &                 1 + iadd, iext, .TRUE. )
126      fbdata%caddname(1)   = 'Hx'
127      DO jvar = 1, profdata%nvar
128         fbdata%cname(jvar)      = profdata%cvars(jvar)
129         fbdata%coblong(jvar)    = profdata%clong(jvar)
130         fbdata%cobunit(jvar)    = profdata%cunit(jvar)
131         fbdata%cgrid(jvar)      = profdata%cgrid(jvar)
132         fbdata%caddlong(1,jvar) = 'Model interpolated ' // TRIM(profdata%clong(jvar))
133         fbdata%caddunit(1,jvar) = profdata%cunit(jvar)
134         IF (iadd > 0) THEN
135            DO ja = 1, iadd
136               fbdata%caddname(1+ja) = padd%cdname(ja)
137               fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
138               fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
139            END DO
140         ENDIF
141         IF (iext > 0) THEN
142            DO je = 1, iext
143               fbdata%cextname(je) = pext%cdname(je)
144               fbdata%cextlong(je) = pext%cdlong(je,1)
145               fbdata%cextunit(je) = pext%cdunit(je,1)
146            END DO
147         ENDIF
148      END DO
149
150      WRITE(clfname, FMT="(A,'fb_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
151
152      IF(lwp) THEN
153         WRITE(numout,*)
154         WRITE(numout,*)'obs_wri_prof :'
155         WRITE(numout,*)'~~~~~~~~~~~~~'
156         WRITE(numout,*)'Writing '//TRIM(clfiletype)//' feedback file : ',TRIM(clfname)
157      ENDIF
158
159      ! Transform obs_prof data structure into obfb data structure
160      fbdata%cdjuldref = '19500101000000'
161      DO jo = 1, profdata%nprof
162         fbdata%plam(jo)      = profdata%rlam(jo)
163         fbdata%pphi(jo)      = profdata%rphi(jo)
164         WRITE(fbdata%cdtyp(jo),'(I4)') profdata%ntyp(jo)
165         fbdata%ivqc(jo,:)    = profdata%ivqc(jo,:)
166         fbdata%ivqcf(:,jo,:) = profdata%ivqcf(:,jo,:)
167         IF ( profdata%nqc(jo) > 255 ) THEN
168            fbdata%ioqc(jo)    = IBSET(profdata%nqc(jo),2)
169            fbdata%ioqcf(1,jo) = profdata%nqcf(1,jo)
170            fbdata%ioqcf(2,jo) = profdata%nqc(jo)
171         ELSE
172            fbdata%ioqc(jo)    = profdata%nqc(jo)
173            fbdata%ioqcf(:,jo) = profdata%nqcf(:,jo)
174         ENDIF
175         fbdata%ipqc(jo)      = profdata%ipqc(jo)
176         fbdata%ipqcf(:,jo)   = profdata%ipqcf(:,jo)
177         fbdata%itqc(jo)      = profdata%itqc(jo)
178         fbdata%itqcf(:,jo)   = profdata%itqcf(:,jo)
179         fbdata%cdwmo(jo)     = profdata%cwmo(jo)
180         fbdata%kindex(jo)    = profdata%npfil(jo)
181         DO jvar = 1, profdata%nvar
182            IF (ln_grid_global) THEN
183               fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar)
184               fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar)
185            ELSE
186               fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar))
187               fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar))
188            ENDIF
189         END DO
190         CALL greg2jul( 0, &
191            &           profdata%nmin(jo), &
192            &           profdata%nhou(jo), &
193            &           profdata%nday(jo), &
194            &           profdata%nmon(jo), &
195            &           profdata%nyea(jo), &
196            &           fbdata%ptim(jo),   &
197            &           krefdate = 19500101 )
198         ! Reform the profiles arrays for output
199         DO jvar = 1, profdata%nvar
200            DO jk = profdata%npvsta(jo,jvar), profdata%npvend(jo,jvar)
201               ik = profdata%var(jvar)%nvlidx(jk)
202               fbdata%padd(ik,jo,1,jvar) = profdata%var(jvar)%vmod(jk)
203               fbdata%pob(ik,jo,jvar)    = profdata%var(jvar)%vobs(jk)
204               fbdata%pdep(ik,jo)        = profdata%var(jvar)%vdep(jk)
205               fbdata%idqc(ik,jo)        = profdata%var(jvar)%idqc(jk)
206               fbdata%idqcf(:,ik,jo)     = profdata%var(jvar)%idqcf(:,jk)
207               IF ( profdata%var(jvar)%nvqc(jk) > 255 ) THEN
208                  fbdata%ivlqc(ik,jo,jvar) = IBSET(profdata%var(jvar)%nvqc(jk),2)
209                  fbdata%ivlqcf(1,ik,jo,jvar) = profdata%var(jvar)%nvqcf(1,jk)
210!$AGRIF_DO_NOT_TREAT
211                  fbdata%ivlqcf(2,ik,jo,jvar) = IAND(profdata%var(jvar)%nvqc(jk),b'0000000011111111')
212!$AGRIF_END_DO_NOT_TREAT
213               ELSE
214                  fbdata%ivlqc(ik,jo,jvar) = profdata%var(jvar)%nvqc(jk)
215                  fbdata%ivlqcf(:,ik,jo,jvar) = profdata%var(jvar)%nvqcf(:,jk)
216               ENDIF
217               fbdata%iobsk(ik,jo,jvar)  = profdata%var(jvar)%mvk(jk)
218               IF (iadd > 0) THEN
219                  DO ja = 1, iadd
220                     fbdata%padd(ik,jo,1+ja,jvar) = &
221                        & profdata%var(jvar)%vadd(jk,padd%ipoint(ja))
222                  END DO
223               ENDIF
224! MOVE OUTSIDE JVAR LOOP?
225               IF (iext > 0) THEN
226                  DO je = 1, iext
227                     fbdata%pext(ik,jo,je) = &
228                        & profdata%vext(jk,pext%ipoint(je))
229                  END DO
230               ENDIF
231            END DO
232         END DO
233      END DO
234
235      ! Convert insitu temperature to potential temperature using the model
236      ! salinity if no potential temperature
237      IF (iext > 0) THEN
238         DO jvar = 1, profdata%nvar
239            IF ( TRIM(profdata%cvars(jvar)) == 'POTM' ) THEN
240               jsal = 0
241               DO jvar2 = 1, profdata%nvar
242                  IF ( TRIM(profdata%cvars(jvar2)) == 'PSAL' ) THEN
243                     jsal = jvar2
244                     EXIT
245                  ENDIF
246               END DO
247               IF (jsal > 0) THEN
248                  DO je = 1, iext
249                     IF ( TRIM(fbdata%cextname(je)) == 'TEMP' ) THEN
250                        DO jo = 1, fbdata%nobs
251                           IF ( fbdata%pphi(jo) < 9999.0 ) THEN
252                              DO jk = 1, fbdata%nlev
253                                 IF ( ( fbdata%pob(jk,jo,jvar)   >= 9999.0 ) .AND. &
254                                    & ( fbdata%pdep(jk,jo)        < 9999.0 ) .AND. &
255                                    & ( fbdata%padd(jk,jo,1,jsal) < 9999.0 ) .AND. &
256                                    & ( fbdata%pext(jk,jo,je)     < 9999.0 ) ) THEN
257                                    zpres = dep_to_p( REAL(fbdata%pdep(jk,jo),wp), &
258                                       &              REAL(fbdata%pphi(jo),wp) )
259                                    fbdata%pob(jk,jo,jvar) = potemp( &
260                                       &                     REAL(fbdata%padd(jk,jo,1,jsal), wp), &
261                                       &                     REAL(fbdata%pext(jk,jo,je), wp),     &
262                                       &                     zpres, 0.0_wp )
263                                 ENDIF
264                              END DO
265                           ENDIF
266                        END DO
267                        EXIT
268                     ENDIF
269                  END DO
270               ENDIF
271               EXIT
272            ENDIF
273         END DO
274      ENDIF
275
276      ! Write the obfbdata structure
277      CALL write_obfbdata( clfname, fbdata )
278
279      ! Output some basic statistics
280      CALL obs_wri_stats( fbdata )
281
282      CALL dealloc_obfbdata( fbdata )
283
284   END SUBROUTINE obs_wri_prof
285
286   SUBROUTINE obs_wri_surf( surfdata, clfiletype, padd, pext )
287      !!-----------------------------------------------------------------------
288      !!
289      !!                     *** ROUTINE obs_wri_surf  ***
290      !!
291      !! ** Purpose : Write surface observation files
292      !!
293      !! ** Method  : NetCDF
294      !!
295      !! ** Action  :
296      !!
297      !!      ! 07-03  (K. Mogensen) Original
298      !!      ! 09-01  (K. Mogensen) New feedback format.
299      !!      ! 15-02  (M. Martin) Combined surface writing routine.
300      !!-----------------------------------------------------------------------
301
302      !! * Modules used
303      IMPLICIT NONE
304
305      !! * Arguments
306      TYPE(obs_surf), INTENT(INOUT) :: surfdata      ! Full set of surface data
307      CHARACTER(LEN=25), INTENT(IN) :: clfiletype    ! Base name for file name
308      TYPE(obswriinfo), OPTIONAL    :: padd          ! Additional info for each variable
309      TYPE(obswriinfo), OPTIONAL    :: pext          ! Extra info
310
311      !! * Local declarations
312      TYPE(obfbdata) :: fbdata
313      CHARACTER(LEN=40) :: clfname         ! netCDF filename
314      CHARACTER(LEN=ilenlong) :: cllongname  ! Long name of variable
315      CHARACTER(LEN=ilenunit) :: clunits     ! Units of variable
316      CHARACTER(LEN=ilengrid) :: clgrid      ! Grid of variable
317      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_surf'
318      INTEGER :: jo
319      INTEGER :: ja
320      INTEGER :: je
321      INTEGER :: jvar
322      INTEGER :: iadd
323      INTEGER :: iext
324
325      IF ( PRESENT( padd ) ) THEN
326         iadd = padd%inum
327      ELSE
328         iadd = 0
329      ENDIF
330
331      IF ( PRESENT( pext ) ) THEN
332         iext = pext%inum
333      ELSE
334         iext = 0
335      ENDIF
336
337      CALL init_obfbdata( fbdata )
338
339      CALL alloc_obfbdata( fbdata, surfdata%nvar, surfdata%nsurf, 1, &
340            &                 1 + iadd, iext, .TRUE. )
341      fbdata%caddname(1)   = 'Hx'
342      DO jvar = 1, surfdata%nvar
343         fbdata%cname(jvar)      = surfdata%cvars(jvar)
344         fbdata%coblong(jvar)    = surfdata%clong(jvar)
345         fbdata%cobunit(jvar)    = surfdata%cunit(jvar)
346         fbdata%cgrid(jvar)      = surfdata%cgrid(jvar)
347         fbdata%caddlong(1,jvar) = 'Model interpolated ' // TRIM(surfdata%clong(jvar))
348         fbdata%caddunit(1,jvar) = surfdata%cunit(jvar)
349         IF (iadd > 0) THEN
350            DO ja = 1, iadd
351               fbdata%caddname(1+ja) = padd%cdname(ja)
352               fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
353               fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
354            END DO
355         ENDIF
356         IF (iext > 0) THEN
357            DO je = 1, iext
358               fbdata%cextname(je) = pext%cdname(je)
359               fbdata%cextlong(je) = pext%cdlong(je,1)
360               fbdata%cextunit(je) = pext%cdunit(je,1)
361            END DO
362         ENDIF
363      END DO
364!fbdata%cname(1)      = surfdata%cvars(1)
365!fbdata%coblong(1)    = 'Sea level anomaly'
366!fbdata%cobunit(1)    = 'Metres'
367!fbdata%cextname(1)   = 'MDT'
368!fbdata%cextlong(1)   = 'Mean dynamic topography'
369!fbdata%cextunit(1)   = 'Metres'
370!fbdata%caddname(2)   = 'SSH'
371!fbdata%caddlong(2,1) = 'Model Sea surface height'
372!fbdata%caddunit(2,1) = 'Metres'
373!IF ( TRIM(surfdata%cvars(1)) == 'ICECONC' ) THEN
374!   fbdata%caddlong(1,1) = 'Model interpolated ICE'
375!ELSE
376!   fbdata%caddlong(1,1) = 'Model interpolated ' // TRIM(surfdata%cvars(1))
377!ENDIF
378
379      WRITE(clfname, FMT="(A,'fb_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
380
381      IF(lwp) THEN
382         WRITE(numout,*)
383         WRITE(numout,*)'obs_wri_surf :'
384         WRITE(numout,*)'~~~~~~~~~~~~~'
385         WRITE(numout,*)'Writing '//TRIM(surfdata%cvars(1))//' feedback file : ',TRIM(clfname)
386      ENDIF
387
388      ! Transform surf data structure into obfbdata structure
389      fbdata%cdjuldref = '19500101000000'
390      DO jo = 1, surfdata%nsurf
391         fbdata%plam(jo)      = surfdata%rlam(jo)
392         fbdata%pphi(jo)      = surfdata%rphi(jo)
393         WRITE(fbdata%cdtyp(jo),'(I4)') surfdata%ntyp(jo)
394         fbdata%ivqc(jo,:)    = 0
395         fbdata%ivqcf(:,jo,:) = 0
396         IF ( surfdata%nqc(jo) > 255 ) THEN
397            fbdata%ioqc(jo)    = 4
398            fbdata%ioqcf(1,jo) = 0
399!$AGRIF_DO_NOT_TREAT
400            fbdata%ioqcf(2,jo) = IAND(surfdata%nqc(jo),b'0000000011111111')
401!$AGRIF_END_DO_NOT_TREAT
402         ELSE
403            fbdata%ioqc(jo)    = surfdata%nqc(jo)
404            fbdata%ioqcf(:,jo) = 0
405         ENDIF
406         fbdata%ipqc(jo)      = 0
407         fbdata%ipqcf(:,jo)   = 0
408         fbdata%itqc(jo)      = 0
409         fbdata%itqcf(:,jo)   = 0
410         fbdata%cdwmo(jo)     = surfdata%cwmo(jo)
411         fbdata%kindex(jo)    = surfdata%nsfil(jo)
412         IF (ln_grid_global) THEN
413            fbdata%iobsi(jo,1) = surfdata%mi(jo)
414            fbdata%iobsj(jo,1) = surfdata%mj(jo)
415         ELSE
416            fbdata%iobsi(jo,1) = mig(surfdata%mi(jo))
417            fbdata%iobsj(jo,1) = mjg(surfdata%mj(jo))
418         ENDIF
419         CALL greg2jul( 0, &
420            &           surfdata%nmin(jo), &
421            &           surfdata%nhou(jo), &
422            &           surfdata%nday(jo), &
423            &           surfdata%nmon(jo), &
424            &           surfdata%nyea(jo), &
425            &           fbdata%ptim(jo),   &
426            &           krefdate = 19500101 ) 
427         fbdata%pdep(1,jo)     = 0.0
428         fbdata%idqc(1,jo)     = 0
429         fbdata%idqcf(:,1,jo)  = 0
430         DO jvar = 1, surfdata%nvar
431            fbdata%padd(1,jo,1,jvar) = surfdata%rmod(jo,jvar)
432!IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%padd(1,jo,2,1) = surfdata%rext(jo,1)
433            fbdata%pob(1,jo,jvar)    = surfdata%robs(jo,jvar)
434            IF ( surfdata%nqc(jo) > 255 ) THEN
435               fbdata%ivqc(jo,jvar)       = 4
436               fbdata%ivlqc(1,jo,jvar)    = 4
437               fbdata%ivlqcf(1,1,jo,jvar) = 0
438!$AGRIF_DO_NOT_TREAT
439               fbdata%ivlqcf(2,1,jo,jvar) = IAND(surfdata%nqc(jo),b'0000000011111111')
440!$AGRIF_END_DO_NOT_TREAT
441            ELSE
442               fbdata%ivqc(jo,jvar)       = surfdata%nqc(jo)
443               fbdata%ivlqc(1,jo,jvar)    = surfdata%nqc(jo)
444               fbdata%ivlqcf(:,1,jo,jvar) = 0
445            ENDIF
446            fbdata%iobsk(1,jo,jvar)  = 0
447!IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%pext(1,jo,1) = surfdata%rext(jo,2)
448            IF (iadd > 0) THEN
449               DO ja = 1, iadd
450                  fbdata%padd(1,jo,1+ja,jvar) = &
451                     & surfdata%radd(jo,padd%ipoint(ja),jvar)
452               END DO
453            ENDIF
454         END DO
455         IF (iext > 0) THEN
456            DO je = 1, iext
457               fbdata%pext(1,jo,je) = &
458                  & surfdata%rext(jo,pext%ipoint(je))
459            END DO
460         ENDIF
461      END DO
462
463      ! Write the obfbdata structure
464      CALL write_obfbdata( clfname, fbdata )
465
466      ! Output some basic statistics
467      CALL obs_wri_stats( fbdata )
468
469      CALL dealloc_obfbdata( fbdata )
470
471   END SUBROUTINE obs_wri_surf
472
473   SUBROUTINE obs_wri_stats( fbdata )
474      !!-----------------------------------------------------------------------
475      !!
476      !!                     *** ROUTINE obs_wri_stats  ***
477      !!
478      !! ** Purpose : Output some basic statistics of the data being written out
479      !!
480      !! ** Method  :
481      !!
482      !! ** Action  :
483      !!
484      !!      ! 2014-08  (D. Lea) Initial version
485      !!-----------------------------------------------------------------------
486
487      !! * Arguments
488      TYPE(obfbdata) :: fbdata
489
490      !! * Local declarations
491      INTEGER :: jvar
492      INTEGER :: jo
493      INTEGER :: jk
494      INTEGER :: inumgoodobs
495      INTEGER :: inumgoodobsmpp
496      REAL(wp) :: zsumx
497      REAL(wp) :: zsumx2
498      REAL(wp) :: zomb
499     
500
501      IF (lwp) THEN
502         WRITE(numout,*) ''
503         WRITE(numout,*) 'obs_wri_stats :'
504         WRITE(numout,*) '~~~~~~~~~~~~~~~'
505      ENDIF
506
507      DO jvar = 1, fbdata%nvar
508         zsumx=0.0_wp
509         zsumx2=0.0_wp
510         inumgoodobs=0
511         DO jo = 1, fbdata%nobs
512            DO jk = 1, fbdata%nlev
513               IF ( ( fbdata%pob(jk,jo,jvar) < 9999.0 ) .AND. &
514                  & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
515                  & ( fbdata%padd(jk,jo,1,jvar) < 9999.0 ) ) THEN
516
517                  zomb=fbdata%pob(jk, jo, jvar)-fbdata%padd(jk, jo, 1, jvar)
518                  zsumx=zsumx+zomb
519                  zsumx2=zsumx2+zomb**2
520                  inumgoodobs=inumgoodobs+1
521               ENDIF
522            ENDDO
523         ENDDO
524
525         CALL obs_mpp_sum_integer( inumgoodobs, inumgoodobsmpp )
526         CALL mpp_sum('obs_write', zsumx)
527         CALL mpp_sum('obs_write', zsumx2)
528
529         IF (lwp) THEN
530            WRITE(numout,*) 'Type: ',fbdata%cname(jvar),'  Total number of good observations: ',inumgoodobsmpp 
531            WRITE(numout,*) 'Overall mean obs minus model of the good observations: ',zsumx/inumgoodobsmpp
532            WRITE(numout,*) 'Overall RMS obs minus model of the good observations: ',sqrt( zsumx2/inumgoodobsmpp )
533            WRITE(numout,*) ''
534         ENDIF
535
536      ENDDO
537
538   END SUBROUTINE obs_wri_stats
539
540END MODULE obs_write
Note: See TracBrowser for help on using the repository browser.