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 branches/UKMO/dev_r4650_general_vert_coord_obsoper/NEMOGCM/NEMO/OPA_SRC/OBS – NEMO

source: branches/UKMO/dev_r4650_general_vert_coord_obsoper/NEMOGCM/NEMO/OPA_SRC/OBS/obs_write.F90 @ 7713

Last change on this file since 7713 was 7713, checked in by dford, 7 years ago

Add observation operator code for surface log10(chlorophyll), SPM, pCO2 and fCO2, for use with FABM-ERSEM, HadOCC and MEDUSA. See internal Met Office NEMO ticket 660.

File size: 55.2 KB
Line 
1MODULE obs_write
2   !!======================================================================
3   !!                       ***  MODULE obs_write   ***
4   !! Observation diagnosticss: Write observation related diagnostics
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   obs_wri_p3d   : Write profile observation diagnostics in NetCDF format
9   !!   obs_wri_sla   : Write SLA observation related diagnostics
10   !!   obs_wri_sst   : Write SST observation related diagnostics
11   !!   obs_wri_seaice: Write seaice observation related diagnostics
12   !!   obs_wri_vel   : Write velocity observation diagnostics in NetCDF format
13   !!   obs_wri_logchl: Write logchl observation related diagnostics
14   !!   obs_wri_spm   : Write spm observation related diagnostics
15   !!   obs_wri_fco2  : Write fco2 observation related diagnostics
16   !!   obs_wri_pco2  : Write fco2 observation related diagnostics
17   !!   obs_wri_stats : Print basic statistics on the data being written out
18   !!----------------------------------------------------------------------
19
20   !! * Modules used
21   USE par_kind, ONLY : &   ! Precision variables
22      & wp
23   USE in_out_manager       ! I/O manager
24   USE dom_oce              ! Ocean space and time domain variables
25   USE obs_types            ! Observation type integer to character translation
26   USE julian, ONLY : &         ! Julian date routines
27      & greg2jul
28   USE obs_utils, ONLY : &  ! Observation operator utility functions
29      & chkerr
30   USE obs_profiles_def     ! Type definitions for profiles
31   USE obs_surf_def         ! Type defintions for surface observations
32   USE obs_fbm              ! Observation feedback I/O
33   USE obs_grid             ! Grid tools
34   USE obs_conv             ! Conversion between units
35   USE obs_const
36   USE obs_sla_types
37   USE obs_rot_vel          ! Rotation of velocities
38   USE obs_mpp              ! MPP support routines for observation diagnostics
39   USE lib_mpp        ! MPP routines
40
41   IMPLICIT NONE
42
43   !! * Routine accessibility
44   PRIVATE
45   PUBLIC obs_wri_p3d, &    ! Write profile observation related diagnostics
46      &   obs_wri_sla, &    ! Write SLA observation related diagnostics
47      &   obs_wri_sst, &    ! Write SST observation related diagnostics
48      &   obs_wri_sss, &    ! Write SSS observation related diagnostics
49      &   obs_wri_seaice, & ! Write seaice observation related diagnostics
50      &   obs_wri_vel, &    ! Write velocity observation related diagnostics
51      &   obs_wri_logchl, & ! Write logchl observation related diagnostics
52      &   obs_wri_spm, &    ! Write spm observation related diagnostics
53      &   obs_wri_fco2, &   ! Write fco2 observation related diagnostics
54      &   obs_wri_pco2, &   ! Write pco2 observation related diagnostics
55      &   obswriinfo
56   
57   TYPE obswriinfo
58      INTEGER :: inum
59      INTEGER, POINTER, DIMENSION(:) :: ipoint
60      CHARACTER(len=ilenname), POINTER, DIMENSION(:) :: cdname
61      CHARACTER(len=ilenlong), POINTER, DIMENSION(:,:) :: cdlong
62      CHARACTER(len=ilenunit), POINTER, DIMENSION(:,:) :: cdunit
63   END TYPE obswriinfo
64
65   !!----------------------------------------------------------------------
66   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
67   !! $Id$
68   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
69   !!----------------------------------------------------------------------
70
71CONTAINS
72
73   SUBROUTINE obs_wri_p3d( cprefix, profdata, padd, pext )
74      !!-----------------------------------------------------------------------
75      !!
76      !!                     *** ROUTINE obs_wri_p3d  ***
77      !!
78      !! ** Purpose : Write temperature and salinity (profile) observation
79      !!              related diagnostics
80      !!
81      !! ** Method  : NetCDF
82      !!
83      !! ** Action  :
84      !!
85      !! History :
86      !!      ! 06-04  (A. Vidard) Original
87      !!      ! 06-04  (A. Vidard) Reformatted
88      !!      ! 06-10  (A. Weaver) Cleanup
89      !!      ! 07-01  (K. Mogensen) Use profile data types
90      !!      ! 07-03  (K. Mogensen) General handling of profiles
91      !!      ! 09-01  (K. Mogensen) New feedback format
92      !!-----------------------------------------------------------------------
93
94      !! * Modules used
95
96      !! * Arguments
97      CHARACTER(LEN=*), INTENT(IN) :: cprefix        ! Prefix for output files
98      TYPE(obs_prof), INTENT(INOUT) :: profdata      ! Full set of profile data
99      TYPE(obswriinfo), OPTIONAL :: padd             ! Additional info for each variable
100      TYPE(obswriinfo), OPTIONAL :: pext             ! Extra info
101     
102      !! * Local declarations
103      TYPE(obfbdata) :: fbdata
104      CHARACTER(LEN=40) :: cfname
105      INTEGER :: ilevel
106      INTEGER :: jvar
107      INTEGER :: jo
108      INTEGER :: jk
109      INTEGER :: ik
110      INTEGER :: ja
111      INTEGER :: je
112      REAL(wp) :: zpres
113      INTEGER :: nadd
114      INTEGER :: next
115
116      IF ( PRESENT( padd ) ) THEN
117         nadd = padd%inum
118      ELSE
119         nadd = 0
120      ENDIF
121
122      IF ( PRESENT( pext ) ) THEN
123         next = pext%inum
124      ELSE
125         next = 0
126      ENDIF
127     
128      CALL init_obfbdata( fbdata )
129
130      ! Find maximum level
131      ilevel = 0
132      DO jvar = 1, 2
133         ilevel = MAX( ilevel, MAXVAL( profdata%var(jvar)%nvlidx(:) ) )
134      END DO
135      CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, &
136         &                 1 + nadd, 1 + next, .TRUE. )
137
138      fbdata%cname(1)      = 'POTM'
139      fbdata%cname(2)      = 'PSAL'
140      fbdata%coblong(1)    = 'Potential temperature'
141      fbdata%coblong(2)    = 'Practical salinity'
142      fbdata%cobunit(1)    = 'Degrees centigrade'
143      fbdata%cobunit(2)    = 'PSU'
144      fbdata%cextname(1)   = 'TEMP'
145      fbdata%cextlong(1)   = 'Insitu temperature'
146      fbdata%cextunit(1)   = 'Degrees centigrade'
147      DO je = 1, next
148         fbdata%cextname(1+je) = pext%cdname(je)
149         fbdata%cextlong(1+je) = pext%cdlong(je,1)
150         fbdata%cextunit(1+je) = pext%cdunit(je,1)
151      END DO
152      fbdata%caddname(1)   = 'Hx'
153      fbdata%caddlong(1,1) = 'Model interpolated potential temperature'
154      fbdata%caddlong(1,2) = 'Model interpolated practical salinity'
155      fbdata%caddunit(1,1) = 'Degrees centigrade'
156      fbdata%caddunit(1,2) = 'PSU'
157      fbdata%cgrid(:)      = 'T'
158      DO ja = 1, nadd
159         fbdata%caddname(1+ja) = padd%cdname(ja)
160         DO jvar = 1, 2
161            fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
162            fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
163         END DO
164      END DO
165         
166      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
167
168      IF(lwp) THEN
169         WRITE(numout,*)
170         WRITE(numout,*)'obs_wri_p3d :'
171         WRITE(numout,*)'~~~~~~~~~~~~~'
172         WRITE(numout,*)'Writing profile feedback file : ',TRIM(cfname)
173      ENDIF
174
175      ! Transform obs_prof data structure into obfbdata structure
176      fbdata%cdjuldref = '19500101000000'
177      DO jo = 1, profdata%nprof
178         fbdata%plam(jo)      = profdata%rlam(jo)
179         fbdata%pphi(jo)      = profdata%rphi(jo)
180         WRITE(fbdata%cdtyp(jo),'(I4)') profdata%ntyp(jo)
181         fbdata%ivqc(jo,:)    = profdata%ivqc(jo,:)
182         fbdata%ivqcf(:,jo,:) = profdata%ivqcf(:,jo,:)
183         IF ( profdata%nqc(jo) > 10 ) THEN
184            fbdata%ioqc(jo)    = 4
185            fbdata%ioqcf(1,jo) = profdata%nqcf(1,jo)
186            fbdata%ioqcf(2,jo) = profdata%nqc(jo) - 10
187         ELSE
188            fbdata%ioqc(jo)    = profdata%nqc(jo)
189            fbdata%ioqcf(:,jo) = profdata%nqcf(:,jo)
190         ENDIF
191         fbdata%ipqc(jo)      = profdata%ipqc(jo)
192         fbdata%ipqcf(:,jo)   = profdata%ipqcf(:,jo)
193         fbdata%itqc(jo)      = profdata%itqc(jo)
194         fbdata%itqcf(:,jo)   = profdata%itqcf(:,jo)
195         fbdata%cdwmo(jo)     = profdata%cwmo(jo)
196         fbdata%kindex(jo)    = profdata%npfil(jo)
197         DO jvar = 1, profdata%nvar
198            IF (ln_grid_global) THEN
199               fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar)
200               fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar)
201            ELSE
202               fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar))
203               fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar))
204            ENDIF
205         END DO
206         CALL greg2jul( 0, &
207            &           profdata%nmin(jo), &
208            &           profdata%nhou(jo), &
209            &           profdata%nday(jo), &
210            &           profdata%nmon(jo), &
211            &           profdata%nyea(jo), &
212            &           fbdata%ptim(jo),   &
213            &           krefdate = 19500101 )
214         ! Reform the profiles arrays for output
215         DO jvar = 1, 2
216            DO jk = profdata%npvsta(jo,jvar), profdata%npvend(jo,jvar)
217               ik = profdata%var(jvar)%nvlidx(jk)
218               fbdata%padd(ik,jo,1,jvar) = profdata%var(jvar)%vmod(jk)
219               fbdata%pob(ik,jo,jvar)    = profdata%var(jvar)%vobs(jk)
220               fbdata%pdep(ik,jo)        = profdata%var(jvar)%vdep(jk)
221               fbdata%idqc(ik,jo)        = profdata%var(jvar)%idqc(jk)
222               fbdata%idqcf(:,ik,jo)     = profdata%var(jvar)%idqcf(:,jk)
223               IF ( profdata%var(jvar)%nvqc(jk) > 10 ) THEN
224                  fbdata%ivlqc(ik,jo,jvar) = 4
225                  fbdata%ivlqcf(1,ik,jo,jvar) = profdata%var(jvar)%nvqcf(1,jk)
226                  fbdata%ivlqcf(2,ik,jo,jvar) = profdata%var(jvar)%nvqc(jk) - 10
227               ELSE
228                  fbdata%ivlqc(ik,jo,jvar) = profdata%var(jvar)%nvqc(jk)
229                  fbdata%ivlqcf(:,ik,jo,jvar) = profdata%var(jvar)%nvqcf(:,jk)
230               ENDIF
231               fbdata%iobsk(ik,jo,jvar)  = profdata%var(jvar)%mvk(jk)
232               DO ja = 1, nadd
233                  fbdata%padd(ik,jo,1+ja,jvar) = &
234                     & profdata%var(jvar)%vext(jk,padd%ipoint(ja))
235               END DO
236               DO je = 1, next
237                  fbdata%pext(ik,jo,1+je) = &
238                     & profdata%var(jvar)%vext(jk,pext%ipoint(je))
239               END DO
240               IF ( jvar == 1 ) THEN
241                  fbdata%pext(ik,jo,1) = profdata%var(jvar)%vext(jk,1)
242               ENDIF
243            END DO
244         END DO
245      END DO
246
247      ! Convert insitu temperature to potential temperature using the model
248      ! salinity if no potential temperature
249      DO jo = 1, fbdata%nobs
250         IF ( fbdata%pphi(jo) < 9999.0 ) THEN
251            DO jk = 1, fbdata%nlev
252               IF ( ( fbdata%pob(jk,jo,1) >= 9999.0 ) .AND. &
253                  & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
254                  & ( fbdata%padd(jk,jo,1,2) < 9999.0 ) .AND. &
255                  & ( fbdata%pext(jk,jo,1) < 9999.0 ) ) THEN
256                  zpres = dep_to_p( REAL(fbdata%pdep(jk,jo),wp), &
257                     &              REAL(fbdata%pphi(jo),wp) )
258                  fbdata%pob(jk,jo,1) = potemp( &
259                     &                     REAL(fbdata%padd(jk,jo,1,2), wp),  &
260                     &                     REAL(fbdata%pext(jk,jo,1), wp), &
261                     &                     zpres, 0.0_wp )
262               ENDIF
263            END DO
264         ENDIF
265      END DO
266     
267      ! Write the obfbdata structure
268      CALL write_obfbdata( cfname, fbdata )
269
270      ! Output some basic statistics
271      CALL obs_wri_stats( fbdata )
272
273      CALL dealloc_obfbdata( fbdata )
274     
275   END SUBROUTINE obs_wri_p3d
276
277   SUBROUTINE obs_wri_sla( cprefix, sladata, padd, pext )
278      !!-----------------------------------------------------------------------
279      !!
280      !!                     *** ROUTINE obs_wri_sla  ***
281      !!
282      !! ** Purpose : Write SLA observation diagnostics
283      !!              related
284      !!
285      !! ** Method  : NetCDF
286      !!
287      !! ** Action  :
288      !!
289      !!      ! 07-03  (K. Mogensen) Original
290      !!      ! 09-01  (K. Mogensen) New feedback format.
291      !!-----------------------------------------------------------------------
292
293      !! * Modules used
294      IMPLICIT NONE
295
296      !! * Arguments
297      CHARACTER(LEN=*), INTENT(IN) :: cprefix          ! Prefix for output files
298      TYPE(obs_surf), INTENT(INOUT) :: sladata         ! Full set of SLAa
299      TYPE(obswriinfo), OPTIONAL :: padd               ! Additional info for each variable
300      TYPE(obswriinfo), OPTIONAL :: pext               ! Extra info
301
302      !! * Local declarations
303      TYPE(obfbdata) :: fbdata
304      CHARACTER(LEN=40) :: cfname         ! netCDF filename
305      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_sla'
306      INTEGER :: jo
307      INTEGER :: ja
308      INTEGER :: je
309      INTEGER :: nadd
310      INTEGER :: next
311
312      IF ( PRESENT( padd ) ) THEN
313         nadd = padd%inum
314      ELSE
315         nadd = 0
316      ENDIF
317
318      IF ( PRESENT( pext ) ) THEN
319         next = pext%inum
320      ELSE
321         next = 0
322      ENDIF
323
324      CALL init_obfbdata( fbdata )
325
326      CALL alloc_obfbdata( fbdata, 1, sladata%nsurf, 1, &
327         &                 2 + nadd, 1 + next, .TRUE. )
328
329      fbdata%cname(1)      = 'SLA'
330      fbdata%coblong(1)    = 'Sea level anomaly'
331      fbdata%cobunit(1)    = 'Metres'
332      fbdata%cextname(1)   = 'MDT'
333      fbdata%cextlong(1)   = 'Mean dynamic topography'
334      fbdata%cextunit(1)   = 'Metres'
335      DO je = 1, next
336         fbdata%cextname(1+je) = pext%cdname(je)
337         fbdata%cextlong(1+je) = pext%cdlong(je,1)
338         fbdata%cextunit(1+je) = pext%cdunit(je,1)
339      END DO
340      fbdata%caddname(1)   = 'Hx'
341      fbdata%caddlong(1,1) = 'Model interpolated SSH - MDT'
342      fbdata%caddunit(1,1) = 'Metres' 
343      fbdata%caddname(2)   = 'SSH'
344      fbdata%caddlong(2,1) = 'Model Sea surface height'
345      fbdata%caddunit(2,1) = 'Metres'
346      fbdata%cgrid(1)      = 'T'
347      DO ja = 1, nadd
348         fbdata%caddname(2+ja) = padd%cdname(ja)
349         fbdata%caddlong(2+ja,1) = padd%cdlong(ja,1)
350         fbdata%caddunit(2+ja,1) = padd%cdunit(ja,1)
351      END DO
352
353      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
354
355      IF(lwp) THEN
356         WRITE(numout,*)
357         WRITE(numout,*)'obs_wri_sla :'
358         WRITE(numout,*)'~~~~~~~~~~~~~'
359         WRITE(numout,*)'Writing SLA feedback file : ',TRIM(cfname)
360      ENDIF
361
362      ! Transform obs_prof data structure into obfbdata structure
363      fbdata%cdjuldref = '19500101000000'
364      DO jo = 1, sladata%nsurf
365         fbdata%plam(jo)      = sladata%rlam(jo)
366         fbdata%pphi(jo)      = sladata%rphi(jo)
367         WRITE(fbdata%cdtyp(jo),'(I4)') sladata%ntyp(jo)
368         fbdata%ivqc(jo,:)    = 0
369         fbdata%ivqcf(:,jo,:) = 0
370         IF ( sladata%nqc(jo) > 10 ) THEN
371            fbdata%ioqc(jo)    = 4
372            fbdata%ioqcf(1,jo) = 0
373            fbdata%ioqcf(2,jo) = sladata%nqc(jo) - 10
374         ELSE
375            fbdata%ioqc(jo)    = sladata%nqc(jo)
376            fbdata%ioqcf(:,jo) = 0
377         ENDIF
378         fbdata%ipqc(jo)      = 0
379         fbdata%ipqcf(:,jo)   = 0
380         fbdata%itqc(jo)      = 0
381         fbdata%itqcf(:,jo)   = 0
382         fbdata%cdwmo(jo)     = sladata%cwmo(jo)
383         fbdata%kindex(jo)    = sladata%nsfil(jo)
384         IF (ln_grid_global) THEN
385            fbdata%iobsi(jo,1) = sladata%mi(jo)
386            fbdata%iobsj(jo,1) = sladata%mj(jo)
387         ELSE
388            fbdata%iobsi(jo,1) = mig(sladata%mi(jo))
389            fbdata%iobsj(jo,1) = mjg(sladata%mj(jo))
390         ENDIF
391         CALL greg2jul( 0, &
392            &           sladata%nmin(jo), &
393            &           sladata%nhou(jo), &
394            &           sladata%nday(jo), &
395            &           sladata%nmon(jo), &
396            &           sladata%nyea(jo), &
397            &           fbdata%ptim(jo),   &
398            &           krefdate = 19500101 )
399         fbdata%padd(1,jo,1,1) = sladata%rmod(jo,1)
400         fbdata%padd(1,jo,2,1) = sladata%rext(jo,1)
401         fbdata%pob(1,jo,1)    = sladata%robs(jo,1) 
402         fbdata%pdep(1,jo)     = 0.0
403         fbdata%idqc(1,jo)     = 0
404         fbdata%idqcf(:,1,jo)  = 0
405         IF ( sladata%nqc(jo) > 10 ) THEN
406            fbdata%ivqc(jo,1)       = 4
407            fbdata%ivlqc(1,jo,1)    = 4
408            fbdata%ivlqcf(1,1,jo,1) = 0
409            fbdata%ivlqcf(2,1,jo,1) = sladata%nqc(jo) - 10
410         ELSE
411            fbdata%ivqc(jo,1)       = sladata%nqc(jo)
412            fbdata%ivlqc(1,jo,1)    = sladata%nqc(jo)
413            fbdata%ivlqcf(:,1,jo,1) = 0
414         ENDIF
415         fbdata%iobsk(1,jo,1)  = 0
416         fbdata%pext(1,jo,1) = sladata%rext(jo,2)
417         DO ja = 1, nadd
418            fbdata%padd(1,jo,2+ja,1) = &
419               & sladata%rext(jo,padd%ipoint(ja))
420         END DO
421         DO je = 1, next
422            fbdata%pext(1,jo,1+je) = &
423               & sladata%rext(jo,pext%ipoint(je))
424         END DO
425      END DO
426
427      ! Write the obfbdata structure
428      CALL write_obfbdata( cfname, fbdata )
429
430      ! Output some basic statistics
431      CALL obs_wri_stats( fbdata )
432
433      CALL dealloc_obfbdata( fbdata )
434
435   END SUBROUTINE obs_wri_sla
436
437   SUBROUTINE obs_wri_sst( cprefix, sstdata, padd, pext )
438      !!-----------------------------------------------------------------------
439      !!
440      !!                     *** ROUTINE obs_wri_sst  ***
441      !!
442      !! ** Purpose : Write SST observation diagnostics
443      !!              related
444      !!
445      !! ** Method  : NetCDF
446      !!
447      !! ** Action  :
448      !!
449      !!      ! 07-07  (S. Ricci) Original
450      !!      ! 09-01  (K. Mogensen) New feedback format.
451      !!-----------------------------------------------------------------------
452
453      !! * Modules used
454      IMPLICIT NONE
455
456      !! * Arguments
457      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
458      TYPE(obs_surf), INTENT(INOUT) :: sstdata      ! Full set of SST
459      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
460      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
461
462      !! * Local declarations
463      TYPE(obfbdata) :: fbdata
464      CHARACTER(LEN=40) ::  cfname             ! netCDF filename
465      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_sst'
466      INTEGER :: jo
467      INTEGER :: ja
468      INTEGER :: je
469      INTEGER :: nadd
470      INTEGER :: next
471
472      IF ( PRESENT( padd ) ) THEN
473         nadd = padd%inum
474      ELSE
475         nadd = 0
476      ENDIF
477
478      IF ( PRESENT( pext ) ) THEN
479         next = pext%inum
480      ELSE
481         next = 0
482      ENDIF
483
484      CALL init_obfbdata( fbdata )
485
486      CALL alloc_obfbdata( fbdata, 1, sstdata%nsurf, 1, &
487         &                 1 + nadd, next, .TRUE. )
488
489      fbdata%cname(1)      = 'SST'
490      fbdata%coblong(1)    = 'Sea surface temperature'
491      fbdata%cobunit(1)    = 'Degree centigrade'
492      DO je = 1, next
493         fbdata%cextname(je) = pext%cdname(je)
494         fbdata%cextlong(je) = pext%cdlong(je,1)
495         fbdata%cextunit(je) = pext%cdunit(je,1)
496      END DO
497      fbdata%caddname(1)   = 'Hx'
498      fbdata%caddlong(1,1) = 'Model interpolated SST'
499      fbdata%caddunit(1,1) = 'Degree centigrade'
500      fbdata%cgrid(1)      = 'T'
501      DO ja = 1, nadd
502         fbdata%caddname(1+ja) = padd%cdname(ja)
503         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
504         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
505      END DO
506
507      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
508
509      IF(lwp) THEN
510         WRITE(numout,*)
511         WRITE(numout,*)'obs_wri_sst :'
512         WRITE(numout,*)'~~~~~~~~~~~~~'
513         WRITE(numout,*)'Writing SST feedback file : ',TRIM(cfname)
514      ENDIF
515
516      ! Transform obs_prof data structure into obfbdata structure
517      fbdata%cdjuldref = '19500101000000'
518      DO jo = 1, sstdata%nsurf
519         fbdata%plam(jo)      = sstdata%rlam(jo)
520         fbdata%pphi(jo)      = sstdata%rphi(jo)
521         WRITE(fbdata%cdtyp(jo),'(I4)') sstdata%ntyp(jo)
522         fbdata%ivqc(jo,:)    = 0
523         fbdata%ivqcf(:,jo,:) = 0
524         IF ( sstdata%nqc(jo) > 10 ) THEN
525            fbdata%ioqc(jo)    = 4
526            fbdata%ioqcf(1,jo) = 0
527            fbdata%ioqcf(2,jo) = sstdata%nqc(jo) - 10
528         ELSE
529            fbdata%ioqc(jo)    = MAX(sstdata%nqc(jo),1)
530            fbdata%ioqcf(:,jo) = 0
531         ENDIF
532         fbdata%ipqc(jo)      = 0
533         fbdata%ipqcf(:,jo)   = 0
534         fbdata%itqc(jo)      = 0
535         fbdata%itqcf(:,jo)   = 0
536         fbdata%cdwmo(jo)     = ''
537         fbdata%kindex(jo)    = sstdata%nsfil(jo)
538         IF (ln_grid_global) THEN
539            fbdata%iobsi(jo,1) = sstdata%mi(jo)
540            fbdata%iobsj(jo,1) = sstdata%mj(jo)
541         ELSE
542            fbdata%iobsi(jo,1) = mig(sstdata%mi(jo))
543            fbdata%iobsj(jo,1) = mjg(sstdata%mj(jo))
544         ENDIF
545         CALL greg2jul( 0, &
546            &           sstdata%nmin(jo), &
547            &           sstdata%nhou(jo), &
548            &           sstdata%nday(jo), &
549            &           sstdata%nmon(jo), &
550            &           sstdata%nyea(jo), &
551            &           fbdata%ptim(jo),   &
552            &           krefdate = 19500101 )
553         fbdata%padd(1,jo,1,1) = sstdata%rmod(jo,1)
554         fbdata%pob(1,jo,1)    = sstdata%robs(jo,1)
555         fbdata%pdep(1,jo)     = 0.0
556         fbdata%idqc(1,jo)     = 0
557         fbdata%idqcf(:,1,jo)  = 0
558         IF ( sstdata%nqc(jo) > 10 ) THEN
559            fbdata%ivqc(jo,1)       = 4
560            fbdata%ivlqc(1,jo,1)    = 4
561            fbdata%ivlqcf(1,1,jo,1) = 0
562            fbdata%ivlqcf(2,1,jo,1) = sstdata%nqc(jo) - 10
563         ELSE
564            fbdata%ivqc(jo,1)       = MAX(sstdata%nqc(jo),1)
565            fbdata%ivlqc(1,jo,1)    = MAX(sstdata%nqc(jo),1)
566            fbdata%ivlqcf(:,1,jo,1) = 0
567         ENDIF
568         fbdata%iobsk(1,jo,1)  = 0
569         DO ja = 1, nadd
570            fbdata%padd(1,jo,1+ja,1) = &
571               & sstdata%rext(jo,padd%ipoint(ja))
572         END DO
573         DO je = 1, next
574            fbdata%pext(1,jo,je) = &
575               & sstdata%rext(jo,pext%ipoint(je))
576         END DO
577
578      END DO
579
580      ! Write the obfbdata structure
581
582      CALL write_obfbdata( cfname, fbdata )
583
584      ! Output some basic statistics
585      CALL obs_wri_stats( fbdata )
586
587      CALL dealloc_obfbdata( fbdata )
588
589   END SUBROUTINE obs_wri_sst
590
591   SUBROUTINE obs_wri_sss
592   END SUBROUTINE obs_wri_sss
593
594   SUBROUTINE obs_wri_seaice( cprefix, seaicedata, padd, pext )
595      !!-----------------------------------------------------------------------
596      !!
597      !!                     *** ROUTINE obs_wri_seaice  ***
598      !!
599      !! ** Purpose : Write sea ice observation diagnostics
600      !!              related
601      !!
602      !! ** Method  : NetCDF
603      !!
604      !! ** Action  :
605      !!
606      !!      ! 07-07  (S. Ricci) Original
607      !!      ! 09-01  (K. Mogensen) New feedback format.
608      !!-----------------------------------------------------------------------
609
610      !! * Modules used
611      IMPLICIT NONE
612
613      !! * Arguments
614      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
615      TYPE(obs_surf), INTENT(INOUT) :: seaicedata   ! Full set of sea ice
616      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
617      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
618
619      !! * Local declarations
620      TYPE(obfbdata) :: fbdata
621      CHARACTER(LEN=40) :: cfname             ! netCDF filename
622      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_seaice'
623      INTEGER :: jo
624      INTEGER :: ja
625      INTEGER :: je
626      INTEGER :: nadd
627      INTEGER :: next
628
629      IF ( PRESENT( padd ) ) THEN
630         nadd = padd%inum
631      ELSE
632         nadd = 0
633      ENDIF
634
635      IF ( PRESENT( pext ) ) THEN
636         next = pext%inum
637      ELSE
638         next = 0
639      ENDIF
640
641      CALL init_obfbdata( fbdata )
642
643      CALL alloc_obfbdata( fbdata, 1, seaicedata%nsurf, 1, 1, 0, .TRUE. )
644
645      fbdata%cname(1)      = 'SEAICE'
646      fbdata%coblong(1)    = 'Sea ice'
647      fbdata%cobunit(1)    = 'Fraction'
648      DO je = 1, next
649         fbdata%cextname(je) = pext%cdname(je)
650         fbdata%cextlong(je) = pext%cdlong(je,1)
651         fbdata%cextunit(je) = pext%cdunit(je,1)
652      END DO
653      fbdata%caddname(1)   = 'Hx'
654      fbdata%caddlong(1,1) = 'Model interpolated ICE'
655      fbdata%caddunit(1,1) = 'Fraction'
656      fbdata%cgrid(1)      = 'T'
657      DO ja = 1, nadd
658         fbdata%caddname(1+ja) = padd%cdname(ja)
659         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
660         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
661      END DO
662
663      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
664
665      IF(lwp) THEN
666         WRITE(numout,*)
667         WRITE(numout,*)'obs_wri_seaice :'
668         WRITE(numout,*)'~~~~~~~~~~~~~~~~'
669         WRITE(numout,*)'Writing SEAICE feedback file : ',TRIM(cfname)
670      ENDIF
671
672      ! Transform obs_prof data structure into obfbdata structure
673      fbdata%cdjuldref = '19500101000000'
674      DO jo = 1, seaicedata%nsurf
675         fbdata%plam(jo)      = seaicedata%rlam(jo)
676         fbdata%pphi(jo)      = seaicedata%rphi(jo)
677         WRITE(fbdata%cdtyp(jo),'(I4)') seaicedata%ntyp(jo)
678         fbdata%ivqc(jo,:)    = 0
679         fbdata%ivqcf(:,jo,:) = 0
680         IF ( seaicedata%nqc(jo) > 10 ) THEN
681            fbdata%ioqc(jo)    = 4
682            fbdata%ioqcf(1,jo) = 0
683            fbdata%ioqcf(2,jo) = seaicedata%nqc(jo) - 10
684         ELSE
685            fbdata%ioqc(jo)    = MAX(seaicedata%nqc(jo),1)
686            fbdata%ioqcf(:,jo) = 0
687         ENDIF
688         fbdata%ipqc(jo)      = 0
689         fbdata%ipqcf(:,jo)   = 0
690         fbdata%itqc(jo)      = 0
691         fbdata%itqcf(:,jo)   = 0
692         fbdata%cdwmo(jo)     = ''
693         fbdata%kindex(jo)    = seaicedata%nsfil(jo)
694         IF (ln_grid_global) THEN
695            fbdata%iobsi(jo,1) = seaicedata%mi(jo)
696            fbdata%iobsj(jo,1) = seaicedata%mj(jo)
697         ELSE
698            fbdata%iobsi(jo,1) = mig(seaicedata%mi(jo))
699            fbdata%iobsj(jo,1) = mjg(seaicedata%mj(jo))
700         ENDIF
701         CALL greg2jul( 0, &
702            &           seaicedata%nmin(jo), &
703            &           seaicedata%nhou(jo), &
704            &           seaicedata%nday(jo), &
705            &           seaicedata%nmon(jo), &
706            &           seaicedata%nyea(jo), &
707            &           fbdata%ptim(jo),   &
708            &           krefdate = 19500101 )
709         fbdata%padd(1,jo,1,1) = seaicedata%rmod(jo,1)
710         fbdata%pob(1,jo,1)    = seaicedata%robs(jo,1)
711         fbdata%pdep(1,jo)     = 0.0
712         fbdata%idqc(1,jo)     = 0
713         fbdata%idqcf(:,1,jo)  = 0
714         IF ( seaicedata%nqc(jo) > 10 ) THEN
715            fbdata%ivlqc(1,jo,1) = 4
716            fbdata%ivlqcf(1,1,jo,1) = 0
717            fbdata%ivlqcf(2,1,jo,1) = seaicedata%nqc(jo) - 10
718         ELSE
719            fbdata%ivlqc(1,jo,1) = MAX(seaicedata%nqc(jo),1)
720            fbdata%ivlqcf(:,1,jo,1) = 0
721         ENDIF
722         fbdata%iobsk(1,jo,1)  = 0
723         DO ja = 1, nadd
724            fbdata%padd(1,jo,1+ja,1) = &
725               & seaicedata%rext(jo,padd%ipoint(ja))
726         END DO
727         DO je = 1, next
728            fbdata%pext(1,jo,je) = &
729               & seaicedata%rext(jo,pext%ipoint(je))
730         END DO
731
732      END DO
733
734      ! Write the obfbdata structure
735      CALL write_obfbdata( cfname, fbdata )
736
737      ! Output some basic statistics
738      CALL obs_wri_stats( fbdata )
739
740      CALL dealloc_obfbdata( fbdata )
741
742   END SUBROUTINE obs_wri_seaice
743
744   SUBROUTINE obs_wri_vel( cprefix, profdata, k2dint, padd, pext )
745      !!-----------------------------------------------------------------------
746      !!
747      !!                     *** ROUTINE obs_wri_vel  ***
748      !!
749      !! ** Purpose : Write current (profile) observation
750      !!              related diagnostics
751      !!
752      !! ** Method  : NetCDF
753      !!
754      !! ** Action  :
755      !!
756      !! History :
757      !!      ! 09-01  (K. Mogensen) New feedback format routine
758      !!-----------------------------------------------------------------------
759
760      !! * Modules used
761
762      !! * Arguments
763      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
764      TYPE(obs_prof), INTENT(INOUT) :: profdata     ! Full set of profile data
765      INTEGER, INTENT(IN) :: k2dint                 ! Horizontal interpolation method
766      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
767      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
768
769      !! * Local declarations
770      TYPE(obfbdata) :: fbdata
771      CHARACTER(LEN=40) :: cfname
772      INTEGER :: ilevel
773      INTEGER :: jvar
774      INTEGER :: jk
775      INTEGER :: ik
776      INTEGER :: jo
777      INTEGER :: ja
778      INTEGER :: je
779      INTEGER :: nadd
780      INTEGER :: next
781      REAL(wp) :: zpres
782      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
783         & zu, &
784         & zv
785
786      IF ( PRESENT( padd ) ) THEN
787         nadd = padd%inum
788      ELSE
789         nadd = 0
790      ENDIF
791
792      IF ( PRESENT( pext ) ) THEN
793         next = pext%inum
794      ELSE
795         next = 0
796      ENDIF
797
798      CALL init_obfbdata( fbdata )
799
800      ! Find maximum level
801      ilevel = 0
802      DO jvar = 1, 2
803         ilevel = MAX( ilevel, MAXVAL( profdata%var(jvar)%nvlidx(:) ) )
804      END DO
805      CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, 2, 0, .TRUE. )
806
807      fbdata%cname(1)      = 'UVEL'
808      fbdata%cname(2)      = 'VVEL'
809      fbdata%coblong(1)    = 'Zonal velocity'
810      fbdata%coblong(2)    = 'Meridional velocity'
811      fbdata%cobunit(1)    = 'm/s'
812      fbdata%cobunit(2)    = 'm/s'
813      DO je = 1, next
814         fbdata%cextname(je) = pext%cdname(je)
815         fbdata%cextlong(je) = pext%cdlong(je,1)
816         fbdata%cextunit(je) = pext%cdunit(je,1)
817      END DO
818      fbdata%caddname(1)   = 'Hx'
819      fbdata%caddlong(1,1) = 'Model interpolated zonal velocity'
820      fbdata%caddlong(1,2) = 'Model interpolated meridional velocity'
821      fbdata%caddunit(1,1) = 'm/s'
822      fbdata%caddunit(1,2) = 'm/s'
823      fbdata%caddname(2)   = 'HxG'
824      fbdata%caddlong(2,1) = 'Model interpolated zonal velocity (model grid)'
825      fbdata%caddlong(2,2) = 'Model interpolated meridional velocity (model grid)'
826      fbdata%caddunit(2,1) = 'm/s'
827      fbdata%caddunit(2,2) = 'm/s' 
828      fbdata%cgrid(1)      = 'U' 
829      fbdata%cgrid(2)      = 'V'
830      DO ja = 1, nadd
831         fbdata%caddname(2+ja) = padd%cdname(ja)
832         fbdata%caddlong(2+ja,1) = padd%cdlong(ja,1)
833         fbdata%caddunit(2+ja,1) = padd%cdunit(ja,1)
834      END DO
835
836      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
837
838      IF(lwp) THEN
839         WRITE(numout,*)
840         WRITE(numout,*)'obs_wri_vel :'
841         WRITE(numout,*)'~~~~~~~~~~~~~'
842         WRITE(numout,*)'Writing velocuty feedback file : ',TRIM(cfname)
843      ENDIF
844
845      ALLOCATE( &
846         & zu(profdata%nvprot(1)), &
847         & zv(profdata%nvprot(2))  &
848         & )
849      CALL obs_rotvel( profdata, k2dint, zu, zv )
850
851      ! Transform obs_prof data structure into obfbdata structure
852      fbdata%cdjuldref = '19500101000000'
853      DO jo = 1, profdata%nprof
854         fbdata%plam(jo)      = profdata%rlam(jo)
855         fbdata%pphi(jo)      = profdata%rphi(jo)
856         WRITE(fbdata%cdtyp(jo),'(I4)') profdata%ntyp(jo)
857         fbdata%ivqc(jo,:)    = profdata%ivqc(jo,:)
858         fbdata%ivqcf(:,jo,:) = profdata%ivqcf(:,jo,:)
859         IF ( profdata%nqc(jo) > 10 ) THEN
860            fbdata%ioqc(jo)    = 4
861            fbdata%ioqcf(1,jo) = profdata%nqcf(1,jo)
862            fbdata%ioqcf(2,jo) = profdata%nqc(jo) - 10
863         ELSE
864            fbdata%ioqc(jo)    = profdata%nqc(jo)
865            fbdata%ioqcf(:,jo) = profdata%nqcf(:,jo)
866         ENDIF
867         fbdata%ipqc(jo)      = profdata%ipqc(jo)
868         fbdata%ipqcf(:,jo)   = profdata%ipqcf(:,jo)
869         fbdata%itqc(jo)      = profdata%itqc(jo)
870         fbdata%itqcf(:,jo)   = profdata%itqcf(:,jo)
871         fbdata%cdwmo(jo)     = profdata%cwmo(jo)
872         fbdata%kindex(jo)    = profdata%npfil(jo)
873         DO jvar = 1, profdata%nvar
874            IF (ln_grid_global) THEN
875               fbdata%iobsi(jo,jvar) = profdata%mi(jo,jvar)
876               fbdata%iobsj(jo,jvar) = profdata%mj(jo,jvar)
877            ELSE
878               fbdata%iobsi(jo,jvar) = mig(profdata%mi(jo,jvar))
879               fbdata%iobsj(jo,jvar) = mjg(profdata%mj(jo,jvar))
880            ENDIF
881         END DO
882         CALL greg2jul( 0, &
883            &           profdata%nmin(jo), &
884            &           profdata%nhou(jo), &
885            &           profdata%nday(jo), &
886            &           profdata%nmon(jo), &
887            &           profdata%nyea(jo), &
888            &           fbdata%ptim(jo),   &
889            &           krefdate = 19500101 )
890         ! Reform the profiles arrays for output
891         DO jvar = 1, 2
892            DO jk = profdata%npvsta(jo,jvar), profdata%npvend(jo,jvar)
893               ik = profdata%var(jvar)%nvlidx(jk)
894               IF ( jvar == 1 ) THEN
895                  fbdata%padd(ik,jo,1,jvar) = zu(jk)
896               ELSE
897                  fbdata%padd(ik,jo,1,jvar) = zv(jk)
898               ENDIF
899               fbdata%padd(ik,jo,2,jvar) = profdata%var(jvar)%vmod(jk)
900               fbdata%pob(ik,jo,jvar)    = profdata%var(jvar)%vobs(jk)
901               fbdata%pdep(ik,jo)        = profdata%var(jvar)%vdep(jk)
902               fbdata%idqc(ik,jo)        = profdata%var(jvar)%idqc(jk)
903               fbdata%idqcf(:,ik,jo)     = profdata%var(jvar)%idqcf(:,jk)
904               IF ( profdata%var(jvar)%nvqc(jk) > 10 ) THEN
905                  fbdata%ivlqc(ik,jo,jvar) = 4
906                  fbdata%ivlqcf(1,ik,jo,jvar) = profdata%var(jvar)%nvqcf(1,jk)
907                  fbdata%ivlqcf(2,ik,jo,jvar) = profdata%var(jvar)%nvqc(jk) - 10
908               ELSE
909                  fbdata%ivlqc(ik,jo,jvar) = profdata%var(jvar)%nvqc(jk)
910                  fbdata%ivlqcf(:,ik,jo,jvar) = profdata%var(jvar)%nvqcf(:,jk)
911               ENDIF
912               fbdata%iobsk(ik,jo,jvar)  = profdata%var(jvar)%mvk(jk)
913               DO ja = 1, nadd
914                  fbdata%padd(ik,jo,2+ja,jvar) = &
915                     & profdata%var(jvar)%vext(jk,padd%ipoint(ja))
916               END DO
917               DO je = 1, next
918                  fbdata%pext(ik,jo,je) = &
919                     & profdata%var(jvar)%vext(jk,pext%ipoint(je))
920               END DO
921            END DO
922         END DO
923      END DO
924
925      ! Write the obfbdata structure
926      CALL write_obfbdata( cfname, fbdata )
927     
928      ! Output some basic statistics
929      CALL obs_wri_stats( fbdata )
930
931      CALL dealloc_obfbdata( fbdata )
932     
933      DEALLOCATE( &
934         & zu, &
935         & zv  &
936         & )
937
938   END SUBROUTINE obs_wri_vel
939
940   SUBROUTINE obs_wri_logchl( cprefix, logchldata, padd, pext )
941      !!-----------------------------------------------------------------------
942      !!
943      !!                     *** ROUTINE obs_wri_logchl  ***
944      !!
945      !! ** Purpose : Write logchl observation diagnostics
946      !!              related
947      !!
948      !! ** Method  : NetCDF
949      !!
950      !! ** Action  :
951      !!
952      !!-----------------------------------------------------------------------
953
954      !! * Modules used
955      IMPLICIT NONE
956
957      !! * Arguments
958      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
959      TYPE(obs_surf), INTENT(INOUT) :: logchldata   ! Full set of logchl
960      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
961      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
962
963      !! * Local declarations
964      TYPE(obfbdata) :: fbdata
965      CHARACTER(LEN=40) :: cfname             ! netCDF filename
966      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_logchl'
967      INTEGER :: jo
968      INTEGER :: ja
969      INTEGER :: je
970      INTEGER :: nadd
971      INTEGER :: next
972
973      IF ( PRESENT( padd ) ) THEN
974         nadd = padd%inum
975      ELSE
976         nadd = 0
977      ENDIF
978
979      IF ( PRESENT( pext ) ) THEN
980         next = pext%inum
981      ELSE
982         next = 0
983      ENDIF
984
985      CALL init_obfbdata( fbdata )
986
987      CALL alloc_obfbdata( fbdata, 1, logchldata%nsurf, 1, &
988         &                 1 + nadd, next, .TRUE. )
989
990      fbdata%cname(1)      = 'LOGCHL'
991      fbdata%coblong(1)    = 'logchl concentration'
992      fbdata%cobunit(1)    = 'mg/m3'
993      DO je = 1, next
994         fbdata%cextname(je) = pext%cdname(je)
995         fbdata%cextlong(je) = pext%cdlong(je,1)
996         fbdata%cextunit(je) = pext%cdunit(je,1)
997      END DO
998      fbdata%caddname(1)   = 'Hx'
999      fbdata%caddlong(1,1) = 'Model interpolated LOGCHL'
1000      fbdata%caddunit(1,1) = 'mg/m3'
1001      fbdata%cgrid(1)      = 'T'
1002      DO ja = 1, nadd
1003         fbdata%caddname(1+ja) = padd%cdname(ja)
1004         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
1005         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
1006      END DO
1007
1008      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
1009
1010      IF(lwp) THEN
1011         WRITE(numout,*)
1012         WRITE(numout,*)'obs_wri_logchl :'
1013         WRITE(numout,*)'~~~~~~~~~~~~~~~~'
1014         WRITE(numout,*)'Writing logchl feedback file : ',TRIM(cfname)
1015      ENDIF
1016
1017      ! Transform obs_prof data structure into obfbdata structure
1018      fbdata%cdjuldref = '19500101000000'
1019      DO jo = 1, logchldata%nsurf
1020         fbdata%plam(jo)      = logchldata%rlam(jo)
1021         fbdata%pphi(jo)      = logchldata%rphi(jo)
1022         WRITE(fbdata%cdtyp(jo),'(I4)') logchldata%ntyp(jo)
1023         fbdata%ivqc(jo,:)    = 0
1024         fbdata%ivqcf(:,jo,:) = 0
1025         IF ( logchldata%nqc(jo) > 10 ) THEN
1026            fbdata%ioqc(jo)    = 4
1027            fbdata%ioqcf(1,jo) = 0
1028            fbdata%ioqcf(2,jo) = logchldata%nqc(jo) - 10
1029         ELSE
1030            fbdata%ioqc(jo)    = MAX(logchldata%nqc(jo),1)
1031            fbdata%ioqcf(:,jo) = 0
1032         ENDIF
1033         fbdata%ipqc(jo)      = 0
1034         fbdata%ipqcf(:,jo)   = 0
1035         fbdata%itqc(jo)      = 0
1036         fbdata%itqcf(:,jo)   = 0
1037         fbdata%cdwmo(jo)     = ''
1038         fbdata%kindex(jo)    = logchldata%nsfil(jo)
1039         IF (ln_grid_global) THEN
1040            fbdata%iobsi(jo,1) = logchldata%mi(jo)
1041            fbdata%iobsj(jo,1) = logchldata%mj(jo)
1042         ELSE
1043            fbdata%iobsi(jo,1) = mig(logchldata%mi(jo))
1044            fbdata%iobsj(jo,1) = mjg(logchldata%mj(jo))
1045         ENDIF
1046         CALL greg2jul( 0, &
1047            &           logchldata%nmin(jo), &
1048            &           logchldata%nhou(jo), &
1049            &           logchldata%nday(jo), &
1050            &           logchldata%nmon(jo), &
1051            &           logchldata%nyea(jo), &
1052            &           fbdata%ptim(jo),   &
1053            &           krefdate = 19500101 )
1054         fbdata%padd(1,jo,1,1) = logchldata%rmod(jo,1)
1055         fbdata%pob(1,jo,1)    = logchldata%robs(jo,1)
1056         fbdata%pdep(1,jo)     = 0.0
1057         fbdata%idqc(1,jo)     = 0
1058         fbdata%idqcf(:,1,jo)  = 0
1059         IF ( logchldata%nqc(jo) > 10 ) THEN
1060            fbdata%ivlqc(1,jo,1) = 4
1061            fbdata%ivlqcf(1,1,jo,1) = 0
1062            fbdata%ivlqcf(2,1,jo,1) = logchldata%nqc(jo) - 10
1063         ELSE
1064            fbdata%ivlqc(1,jo,1) = MAX(logchldata%nqc(jo),1)
1065            fbdata%ivlqcf(:,1,jo,1) = 0
1066         ENDIF
1067         fbdata%iobsk(1,jo,1)  = 0
1068         DO ja = 1, nadd
1069            fbdata%padd(1,jo,1+ja,1) = &
1070               & logchldata%rext(jo,padd%ipoint(ja))
1071         END DO
1072         DO je = 1, next
1073            fbdata%pext(1,jo,je) = &
1074               & logchldata%rext(jo,pext%ipoint(je))
1075         END DO
1076
1077      END DO
1078
1079      ! Write the obfbdata structure
1080      CALL write_obfbdata( cfname, fbdata )
1081     
1082      ! Output some basic statistics
1083      CALL obs_wri_stats( fbdata )
1084
1085      CALL dealloc_obfbdata( fbdata )
1086
1087   END SUBROUTINE obs_wri_logchl
1088
1089   SUBROUTINE obs_wri_spm( cprefix, spmdata, padd, pext )
1090      !!-----------------------------------------------------------------------
1091      !!
1092      !!                     *** ROUTINE obs_wri_spm  ***
1093      !!
1094      !! ** Purpose : Write spm observation diagnostics
1095      !!              related
1096      !!
1097      !! ** Method  : NetCDF
1098      !!
1099      !! ** Action  :
1100      !!
1101      !!-----------------------------------------------------------------------
1102
1103      !! * Modules used
1104      IMPLICIT NONE
1105
1106      !! * Arguments
1107      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
1108      TYPE(obs_surf), INTENT(INOUT) :: spmdata   ! Full set of spm
1109      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
1110      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
1111
1112      !! * Local declarations
1113      TYPE(obfbdata) :: fbdata
1114      CHARACTER(LEN=40) :: cfname             ! netCDF filename
1115      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_spm'
1116      INTEGER :: jo
1117      INTEGER :: ja
1118      INTEGER :: je
1119      INTEGER :: nadd
1120      INTEGER :: next
1121
1122      IF ( PRESENT( padd ) ) THEN
1123         nadd = padd%inum
1124      ELSE
1125         nadd = 0
1126      ENDIF
1127
1128      IF ( PRESENT( pext ) ) THEN
1129         next = pext%inum
1130      ELSE
1131         next = 0
1132      ENDIF
1133
1134      CALL init_obfbdata( fbdata )
1135
1136      CALL alloc_obfbdata( fbdata, 1, spmdata%nsurf, 1, &
1137         &                 1 + nadd, next, .TRUE. )
1138
1139      fbdata%cname(1)      = 'spm'
1140      fbdata%coblong(1)    = 'spm'
1141      fbdata%cobunit(1)    = 'g/m3'
1142      DO je = 1, next
1143         fbdata%cextname(je) = pext%cdname(je)
1144         fbdata%cextlong(je) = pext%cdlong(je,1)
1145         fbdata%cextunit(je) = pext%cdunit(je,1)
1146      END DO
1147      fbdata%caddname(1)   = 'Hx'
1148      fbdata%caddlong(1,1) = 'Model interpolated spm'
1149      fbdata%caddunit(1,1) = 'g/m3'
1150      fbdata%cgrid(1)      = 'T'
1151      DO ja = 1, nadd
1152         fbdata%caddname(1+ja) = padd%cdname(ja)
1153         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
1154         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
1155      END DO
1156
1157      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
1158
1159      IF(lwp) THEN
1160         WRITE(numout,*)
1161         WRITE(numout,*)'obs_wri_spm :'
1162         WRITE(numout,*)'~~~~~~~~~~~~~~~~'
1163         WRITE(numout,*)'Writing spm feedback file : ',TRIM(cfname)
1164      ENDIF
1165
1166      ! Transform obs_prof data structure into obfbdata structure
1167      fbdata%cdjuldref = '19500101000000'
1168      DO jo = 1, spmdata%nsurf
1169         fbdata%plam(jo)      = spmdata%rlam(jo)
1170         fbdata%pphi(jo)      = spmdata%rphi(jo)
1171         WRITE(fbdata%cdtyp(jo),'(I4)') spmdata%ntyp(jo)
1172         fbdata%ivqc(jo,:)    = 0
1173         fbdata%ivqcf(:,jo,:) = 0
1174         IF ( spmdata%nqc(jo) > 10 ) THEN
1175            fbdata%ioqc(jo)    = 4
1176            fbdata%ioqcf(1,jo) = 0
1177            fbdata%ioqcf(2,jo) = spmdata%nqc(jo) - 10
1178         ELSE
1179            fbdata%ioqc(jo)    = MAX(spmdata%nqc(jo),1)
1180            fbdata%ioqcf(:,jo) = 0
1181         ENDIF
1182         fbdata%ipqc(jo)      = 0
1183         fbdata%ipqcf(:,jo)   = 0
1184         fbdata%itqc(jo)      = 0
1185         fbdata%itqcf(:,jo)   = 0
1186         fbdata%cdwmo(jo)     = ''
1187         fbdata%kindex(jo)    = spmdata%nsfil(jo)
1188         IF (ln_grid_global) THEN
1189            fbdata%iobsi(jo,1) = spmdata%mi(jo)
1190            fbdata%iobsj(jo,1) = spmdata%mj(jo)
1191         ELSE
1192            fbdata%iobsi(jo,1) = mig(spmdata%mi(jo))
1193            fbdata%iobsj(jo,1) = mjg(spmdata%mj(jo))
1194         ENDIF
1195         CALL greg2jul( 0, &
1196            &           spmdata%nmin(jo), &
1197            &           spmdata%nhou(jo), &
1198            &           spmdata%nday(jo), &
1199            &           spmdata%nmon(jo), &
1200            &           spmdata%nyea(jo), &
1201            &           fbdata%ptim(jo),   &
1202            &           krefdate = 19500101 )
1203         fbdata%padd(1,jo,1,1) = spmdata%rmod(jo,1)
1204         fbdata%pob(1,jo,1)    = spmdata%robs(jo,1)
1205         fbdata%pdep(1,jo)     = 0.0
1206         fbdata%idqc(1,jo)     = 0
1207         fbdata%idqcf(:,1,jo)  = 0
1208         IF ( spmdata%nqc(jo) > 10 ) THEN
1209            fbdata%ivlqc(1,jo,1) = 4
1210            fbdata%ivlqcf(1,1,jo,1) = 0
1211            fbdata%ivlqcf(2,1,jo,1) = spmdata%nqc(jo) - 10
1212         ELSE
1213            fbdata%ivlqc(1,jo,1) = MAX(spmdata%nqc(jo),1)
1214            fbdata%ivlqcf(:,1,jo,1) = 0
1215         ENDIF
1216         fbdata%iobsk(1,jo,1)  = 0
1217         DO ja = 1, nadd
1218            fbdata%padd(1,jo,1+ja,1) = &
1219               & spmdata%rext(jo,padd%ipoint(ja))
1220         END DO
1221         DO je = 1, next
1222            fbdata%pext(1,jo,je) = &
1223               & spmdata%rext(jo,pext%ipoint(je))
1224         END DO
1225
1226      END DO
1227
1228      ! Write the obfbdata structure
1229      CALL write_obfbdata( cfname, fbdata )
1230     
1231      ! Output some basic statistics
1232      CALL obs_wri_stats( fbdata )
1233
1234      CALL dealloc_obfbdata( fbdata )
1235
1236   END SUBROUTINE obs_wri_spm
1237
1238   SUBROUTINE obs_wri_fco2( cprefix, fco2data, padd, pext )
1239      !!-----------------------------------------------------------------------
1240      !!
1241      !!                     *** ROUTINE obs_wri_fco2  ***
1242      !!
1243      !! ** Purpose : Write fco2 observation diagnostics
1244      !!              related
1245      !!
1246      !! ** Method  : NetCDF
1247      !!
1248      !! ** Action  :
1249      !!
1250      !!-----------------------------------------------------------------------
1251
1252      !! * Modules used
1253      IMPLICIT NONE
1254
1255      !! * Arguments
1256      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
1257      TYPE(obs_surf), INTENT(INOUT) :: fco2data   ! Full set of fco2
1258      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
1259      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
1260
1261      !! * Local declarations
1262      TYPE(obfbdata) :: fbdata
1263      CHARACTER(LEN=40) :: cfname             ! netCDF filename
1264      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_fco2'
1265      INTEGER :: jo
1266      INTEGER :: ja
1267      INTEGER :: je
1268      INTEGER :: nadd
1269      INTEGER :: next
1270
1271      IF ( PRESENT( padd ) ) THEN
1272         nadd = padd%inum
1273      ELSE
1274         nadd = 0
1275      ENDIF
1276
1277      IF ( PRESENT( pext ) ) THEN
1278         next = pext%inum
1279      ELSE
1280         next = 0
1281      ENDIF
1282
1283      CALL init_obfbdata( fbdata )
1284
1285      CALL alloc_obfbdata( fbdata, 1, fco2data%nsurf, 1, &
1286         &                 1 + nadd, next, .TRUE. )
1287
1288      fbdata%cname(1)      = 'fco2'
1289      fbdata%coblong(1)    = 'fco2'
1290      fbdata%cobunit(1)    = 'uatm'
1291      DO je = 1, next
1292         fbdata%cextname(je) = pext%cdname(je)
1293         fbdata%cextlong(je) = pext%cdlong(je,1)
1294         fbdata%cextunit(je) = pext%cdunit(je,1)
1295      END DO
1296      fbdata%caddname(1)   = 'Hx'
1297      fbdata%caddlong(1,1) = 'Model interpolated fco2'
1298      fbdata%caddunit(1,1) = 'uatm'
1299      fbdata%cgrid(1)      = 'T'
1300      DO ja = 1, nadd
1301         fbdata%caddname(1+ja) = padd%cdname(ja)
1302         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
1303         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
1304      END DO
1305
1306      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
1307
1308      IF(lwp) THEN
1309         WRITE(numout,*)
1310         WRITE(numout,*)'obs_wri_fco2 :'
1311         WRITE(numout,*)'~~~~~~~~~~~~~~~~'
1312         WRITE(numout,*)'Writing fco2 feedback file : ',TRIM(cfname)
1313      ENDIF
1314
1315      ! Transform obs_prof data structure into obfbdata structure
1316      fbdata%cdjuldref = '19500101000000'
1317      DO jo = 1, fco2data%nsurf
1318         fbdata%plam(jo)      = fco2data%rlam(jo)
1319         fbdata%pphi(jo)      = fco2data%rphi(jo)
1320         WRITE(fbdata%cdtyp(jo),'(I4)') fco2data%ntyp(jo)
1321         fbdata%ivqc(jo,:)    = 0
1322         fbdata%ivqcf(:,jo,:) = 0
1323         IF ( fco2data%nqc(jo) > 10 ) THEN
1324            fbdata%ioqc(jo)    = 4
1325            fbdata%ioqcf(1,jo) = 0
1326            fbdata%ioqcf(2,jo) = fco2data%nqc(jo) - 10
1327         ELSE
1328            fbdata%ioqc(jo)    = MAX(fco2data%nqc(jo),1)
1329            fbdata%ioqcf(:,jo) = 0
1330         ENDIF
1331         fbdata%ipqc(jo)      = 0
1332         fbdata%ipqcf(:,jo)   = 0
1333         fbdata%itqc(jo)      = 0
1334         fbdata%itqcf(:,jo)   = 0
1335         fbdata%cdwmo(jo)     = ''
1336         fbdata%kindex(jo)    = fco2data%nsfil(jo)
1337         IF (ln_grid_global) THEN
1338            fbdata%iobsi(jo,1) = fco2data%mi(jo)
1339            fbdata%iobsj(jo,1) = fco2data%mj(jo)
1340         ELSE
1341            fbdata%iobsi(jo,1) = mig(fco2data%mi(jo))
1342            fbdata%iobsj(jo,1) = mjg(fco2data%mj(jo))
1343         ENDIF
1344         CALL greg2jul( 0, &
1345            &           fco2data%nmin(jo), &
1346            &           fco2data%nhou(jo), &
1347            &           fco2data%nday(jo), &
1348            &           fco2data%nmon(jo), &
1349            &           fco2data%nyea(jo), &
1350            &           fbdata%ptim(jo),   &
1351            &           krefdate = 19500101 )
1352         fbdata%padd(1,jo,1,1) = fco2data%rmod(jo,1)
1353         fbdata%pob(1,jo,1)    = fco2data%robs(jo,1)
1354         fbdata%pdep(1,jo)     = 0.0
1355         fbdata%idqc(1,jo)     = 0
1356         fbdata%idqcf(:,1,jo)  = 0
1357         IF ( fco2data%nqc(jo) > 10 ) THEN
1358            fbdata%ivlqc(1,jo,1) = 4
1359            fbdata%ivlqcf(1,1,jo,1) = 0
1360            fbdata%ivlqcf(2,1,jo,1) = fco2data%nqc(jo) - 10
1361         ELSE
1362            fbdata%ivlqc(1,jo,1) = MAX(fco2data%nqc(jo),1)
1363            fbdata%ivlqcf(:,1,jo,1) = 0
1364         ENDIF
1365         fbdata%iobsk(1,jo,1)  = 0
1366         DO ja = 1, nadd
1367            fbdata%padd(1,jo,1+ja,1) = &
1368               & fco2data%rext(jo,padd%ipoint(ja))
1369         END DO
1370         DO je = 1, next
1371            fbdata%pext(1,jo,je) = &
1372               & fco2data%rext(jo,pext%ipoint(je))
1373         END DO
1374
1375      END DO
1376
1377      ! Write the obfbdata structure
1378      CALL write_obfbdata( cfname, fbdata )
1379     
1380      ! Output some basic statistics
1381      CALL obs_wri_stats( fbdata )
1382
1383      CALL dealloc_obfbdata( fbdata )
1384
1385   END SUBROUTINE obs_wri_fco2
1386
1387   SUBROUTINE obs_wri_pco2( cprefix, pco2data, padd, pext )
1388      !!-----------------------------------------------------------------------
1389      !!
1390      !!                     *** ROUTINE obs_wri_pco2  ***
1391      !!
1392      !! ** Purpose : Write pco2 observation diagnostics
1393      !!              related
1394      !!
1395      !! ** Method  : NetCDF
1396      !!
1397      !! ** Action  :
1398      !!
1399      !!-----------------------------------------------------------------------
1400
1401      !! * Modules used
1402      IMPLICIT NONE
1403
1404      !! * Arguments
1405      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files
1406      TYPE(obs_surf), INTENT(INOUT) :: pco2data   ! Full set of pco2
1407      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable
1408      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info
1409
1410      !! * Local declarations
1411      TYPE(obfbdata) :: fbdata
1412      CHARACTER(LEN=40) :: cfname             ! netCDF filename
1413      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_pco2'
1414      INTEGER :: jo
1415      INTEGER :: ja
1416      INTEGER :: je
1417      INTEGER :: nadd
1418      INTEGER :: next
1419
1420      IF ( PRESENT( padd ) ) THEN
1421         nadd = padd%inum
1422      ELSE
1423         nadd = 0
1424      ENDIF
1425
1426      IF ( PRESENT( pext ) ) THEN
1427         next = pext%inum
1428      ELSE
1429         next = 0
1430      ENDIF
1431
1432      CALL init_obfbdata( fbdata )
1433
1434      CALL alloc_obfbdata( fbdata, 1, pco2data%nsurf, 1, &
1435         &                 1 + nadd, next, .TRUE. )
1436
1437      fbdata%cname(1)      = 'pco2'
1438      fbdata%coblong(1)    = 'pco2'
1439      fbdata%cobunit(1)    = 'uatm'
1440      DO je = 1, next
1441         fbdata%cextname(je) = pext%cdname(je)
1442         fbdata%cextlong(je) = pext%cdlong(je,1)
1443         fbdata%cextunit(je) = pext%cdunit(je,1)
1444      END DO
1445      fbdata%caddname(1)   = 'Hx'
1446      fbdata%caddlong(1,1) = 'Model interpolated pco2'
1447      fbdata%caddunit(1,1) = 'uatm'
1448      fbdata%cgrid(1)      = 'T'
1449      DO ja = 1, nadd
1450         fbdata%caddname(1+ja) = padd%cdname(ja)
1451         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
1452         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
1453      END DO
1454
1455      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc
1456
1457      IF(lwp) THEN
1458         WRITE(numout,*)
1459         WRITE(numout,*)'obs_wri_pco2 :'
1460         WRITE(numout,*)'~~~~~~~~~~~~~~~~'
1461         WRITE(numout,*)'Writing pco2 feedback file : ',TRIM(cfname)
1462      ENDIF
1463
1464      ! Transform obs_prof data structure into obfbdata structure
1465      fbdata%cdjuldref = '19500101000000'
1466      DO jo = 1, pco2data%nsurf
1467         fbdata%plam(jo)      = pco2data%rlam(jo)
1468         fbdata%pphi(jo)      = pco2data%rphi(jo)
1469         WRITE(fbdata%cdtyp(jo),'(I4)') pco2data%ntyp(jo)
1470         fbdata%ivqc(jo,:)    = 0
1471         fbdata%ivqcf(:,jo,:) = 0
1472         IF ( pco2data%nqc(jo) > 10 ) THEN
1473            fbdata%ioqc(jo)    = 4
1474            fbdata%ioqcf(1,jo) = 0
1475            fbdata%ioqcf(2,jo) = pco2data%nqc(jo) - 10
1476         ELSE
1477            fbdata%ioqc(jo)    = MAX(pco2data%nqc(jo),1)
1478            fbdata%ioqcf(:,jo) = 0
1479         ENDIF
1480         fbdata%ipqc(jo)      = 0
1481         fbdata%ipqcf(:,jo)   = 0
1482         fbdata%itqc(jo)      = 0
1483         fbdata%itqcf(:,jo)   = 0
1484         fbdata%cdwmo(jo)     = ''
1485         fbdata%kindex(jo)    = pco2data%nsfil(jo)
1486         IF (ln_grid_global) THEN
1487            fbdata%iobsi(jo,1) = pco2data%mi(jo)
1488            fbdata%iobsj(jo,1) = pco2data%mj(jo)
1489         ELSE
1490            fbdata%iobsi(jo,1) = mig(pco2data%mi(jo))
1491            fbdata%iobsj(jo,1) = mjg(pco2data%mj(jo))
1492         ENDIF
1493         CALL greg2jul( 0, &
1494            &           pco2data%nmin(jo), &
1495            &           pco2data%nhou(jo), &
1496            &           pco2data%nday(jo), &
1497            &           pco2data%nmon(jo), &
1498            &           pco2data%nyea(jo), &
1499            &           fbdata%ptim(jo),   &
1500            &           krefdate = 19500101 )
1501         fbdata%padd(1,jo,1,1) = pco2data%rmod(jo,1)
1502         fbdata%pob(1,jo,1)    = pco2data%robs(jo,1)
1503         fbdata%pdep(1,jo)     = 0.0
1504         fbdata%idqc(1,jo)     = 0
1505         fbdata%idqcf(:,1,jo)  = 0
1506         IF ( pco2data%nqc(jo) > 10 ) THEN
1507            fbdata%ivlqc(1,jo,1) = 4
1508            fbdata%ivlqcf(1,1,jo,1) = 0
1509            fbdata%ivlqcf(2,1,jo,1) = pco2data%nqc(jo) - 10
1510         ELSE
1511            fbdata%ivlqc(1,jo,1) = MAX(pco2data%nqc(jo),1)
1512            fbdata%ivlqcf(:,1,jo,1) = 0
1513         ENDIF
1514         fbdata%iobsk(1,jo,1)  = 0
1515         DO ja = 1, nadd
1516            fbdata%padd(1,jo,1+ja,1) = &
1517               & pco2data%rext(jo,padd%ipoint(ja))
1518         END DO
1519         DO je = 1, next
1520            fbdata%pext(1,jo,je) = &
1521               & pco2data%rext(jo,pext%ipoint(je))
1522         END DO
1523
1524      END DO
1525
1526      ! Write the obfbdata structure
1527      CALL write_obfbdata( cfname, fbdata )
1528     
1529      ! Output some basic statistics
1530      CALL obs_wri_stats( fbdata )
1531
1532      CALL dealloc_obfbdata( fbdata )
1533
1534   END SUBROUTINE obs_wri_pco2
1535
1536   SUBROUTINE obs_wri_stats( fbdata )
1537      !!-----------------------------------------------------------------------
1538      !!
1539      !!                     *** ROUTINE obs_wri_stats  ***
1540      !!
1541      !! ** Purpose : Output some basic statistics of the data being written out
1542      !!
1543      !! ** Method  :
1544      !!
1545      !! ** Action  :
1546      !!
1547      !!      ! 2014-08  (D. Lea) Initial version
1548      !!-----------------------------------------------------------------------
1549
1550      !! * Arguments
1551      TYPE(obfbdata) :: fbdata
1552
1553      !! * Local declarations
1554      INTEGER :: jvar
1555      INTEGER :: jo
1556      INTEGER :: jk
1557
1558!      INTEGER :: nlev
1559!      INTEGER :: nlevmpp
1560!      INTEGER :: nobsmpp
1561      INTEGER :: numgoodobs
1562      INTEGER :: numgoodobsmpp
1563      REAL(wp) :: zsumx
1564      REAL(wp) :: zsumx2
1565      REAL(wp) :: zomb
1566
1567      IF (lwp) THEN
1568         WRITE(numout,*) ''
1569         WRITE(numout,*) 'obs_wri_stats :'
1570         WRITE(numout,*) '~~~~~~~~~~~~~~~' 
1571      ENDIF
1572
1573      DO jvar = 1, fbdata%nvar
1574         zsumx=0.0_wp
1575         zsumx2=0.0_wp
1576         numgoodobs=0
1577         DO jo = 1, fbdata%nobs
1578            DO jk = 1, fbdata%nlev
1579               IF ( ( fbdata%pob(jk,jo,jvar) < 9999.0 ) .AND. &
1580                  & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
1581                  & ( fbdata%padd(jk,jo,1,jvar) < 9999.0 ) ) THEN
1582       
1583             zomb=fbdata%pob(jk, jo, jvar)-fbdata%padd(jk, jo, 1, jvar)
1584                  zsumx=zsumx+zomb
1585                  zsumx2=zsumx2+zomb**2
1586                  numgoodobs=numgoodobs+1
1587          ENDIF
1588            ENDDO
1589         ENDDO
1590
1591         CALL obs_mpp_sum_integer( numgoodobs, numgoodobsmpp )
1592         CALL mpp_sum(zsumx)
1593         CALL mpp_sum(zsumx2)
1594
1595         IF (lwp) THEN
1596       WRITE(numout,*) 'Type: ',fbdata%cname(jvar),'  Total number of good observations: ',numgoodobsmpp 
1597       WRITE(numout,*) 'Overall mean obs minus model of the good observations: ',zsumx/numgoodobsmpp
1598            WRITE(numout,*) 'Overall RMS obs minus model of the good observations: ',sqrt( zsumx2/numgoodobsmpp )
1599       WRITE(numout,*) ''
1600         ENDIF
1601 
1602      ENDDO
1603
1604   END SUBROUTINE obs_wri_stats
1605
1606END MODULE obs_write
Note: See TracBrowser for help on using the repository browser.