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_read_surf.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_read_surf.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: 26.5 KB
Line 
1MODULE obs_read_surf
2   !!======================================================================
3   !!                       ***  MODULE obs_read_surf  ***
4   !! Observation diagnostics: Read the surface data from feedback files
5   !!======================================================================
6
7   !!----------------------------------------------------------------------
8   !!   obs_rea_surf : Driver for reading surface data from feedback files
9   !!----------------------------------------------------------------------
10
11   !! * Modules used
12   USE par_kind                 ! Precision variables
13   USE in_out_manager           ! I/O manager
14   USE dom_oce                  ! Ocean space and time domain variables
15   USE obs_mpp                  ! MPP support routines for observation diagnostics
16   USE julian                   ! Julian date routines
17   USE obs_utils                ! Observation operator utility functions
18   USE obs_grid                 ! Grid search
19   USE obs_sort                 ! Sorting observation arrays
20   USE obs_surf_def             ! Surface observation definitions
21   USE obs_types                ! Observation type definitions
22   USE obs_fbm                  ! Feedback routines
23   USE netcdf                   ! NetCDF library
24
25   IMPLICIT NONE
26
27   !! * Routine accessibility
28   PRIVATE
29
30   PUBLIC obs_rea_surf      ! Read the surface observations from the point data
31
32   !!----------------------------------------------------------------------
33   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
34   !! $Id$
35   !! Software governed by the CeCILL license (see ./LICENSE)
36   !!----------------------------------------------------------------------
37
38CONTAINS
39
40   SUBROUTINE obs_rea_surf( surfdata, knumfiles, cdfilenames, &
41      &                     kvars, kadd, kextr, kstp, ddobsini, ddobsend, &
42      &                     ldignmis, ldmod, ldnightav, cdvars )
43      !!---------------------------------------------------------------------
44      !!
45      !!                   *** ROUTINE obs_rea_surf ***
46      !!
47      !! ** Purpose : Read from file the surface data
48      !!
49      !! ** Method  : Read in the data from feedback format files and
50      !!              put into the NEMO internal surface data structure
51      !!
52      !! ** Action  :
53      !!
54      !!
55      !! History : 
56      !!      ! :  2009-01 (K. Mogensen) Initial version based on old versions
57      !!      ! :  2015-02 (M. Martin)   Unify the different surface data type reading.
58      !!----------------------------------------------------------------------
59      !! * Modules used
60
61      !! * Arguments
62      TYPE(obs_surf), INTENT(INOUT) :: &
63         & surfdata                     ! Surface data to be read
64      INTEGER, INTENT(IN) :: knumfiles  ! Number of corio format files to read
65      CHARACTER(LEN=128), INTENT(IN) :: &
66         & cdfilenames(knumfiles)       ! File names to read in
67      INTEGER, INTENT(IN) :: kvars      ! Number of variables in surfdata
68      INTEGER, INTENT(IN) :: kadd       ! Number of additional fields
69                                        !   in addition to those in the input file(s)
70      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields
71                                        !   in addition to those in the input file(s)
72      INTEGER, INTENT(IN) :: kstp       ! Ocean time-step index
73      LOGICAL, INTENT(IN) :: ldignmis   ! Ignore missing files
74      LOGICAL, INTENT(IN) :: ldmod      ! Initialize model from input data
75      LOGICAL, INTENT(IN) :: ldnightav  ! Observations represent a night-time average
76      REAL(dp), INTENT(IN) :: ddobsini   ! Obs. ini time in YYYYMMDD.HHMMSS
77      REAL(dp), INTENT(IN) :: ddobsend   ! Obs. end time in YYYYMMDD.HHMMSS
78      CHARACTER(len=8), DIMENSION(kvars), INTENT(IN) :: cdvars
79
80      !! * Local declarations
81      CHARACTER(LEN=11), PARAMETER :: cpname='obs_rea_surf'
82      CHARACTER(len=8) :: clrefdate
83      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clvarsin
84      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: cllongin
85      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clunitin
86      CHARACTER(len=ilengrid), DIMENSION(:),   ALLOCATABLE :: clgridin
87      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: claddvarsin
88      CHARACTER(len=ilenlong), DIMENSION(:,:), ALLOCATABLE :: claddlongin
89      CHARACTER(len=ilenunit), DIMENSION(:,:), ALLOCATABLE :: claddunitin
90      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clextvarsin
91      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: clextlongin
92      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clextunitin
93      INTEGER :: ji
94      INTEGER :: jj
95      INTEGER :: jk
96      INTEGER :: jvar
97      INTEGER :: jext
98      INTEGER :: jadd
99      INTEGER :: jadd2
100      INTEGER :: iadd
101      INTEGER :: iaddin
102      INTEGER :: iextr
103      INTEGER :: iflag
104      INTEGER :: inobf
105      INTEGER :: i_file_id
106      INTEGER :: inowin
107      INTEGER :: iyea
108      INTEGER :: imon
109      INTEGER :: iday
110      INTEGER :: ihou
111      INTEGER :: imin
112      INTEGER :: isec
113      INTEGER :: itype
114      INTEGER :: iobsmpp
115      INTEGER :: iobs
116      INTEGER :: iobstot
117      INTEGER :: ios
118      INTEGER :: ioserrcount
119      INTEGER, PARAMETER :: jpsurfmaxtype = 1024
120      INTEGER, DIMENSION(knumfiles) :: irefdate
121      INTEGER, DIMENSION(jpsurfmaxtype+1) :: &
122         & ityp, &
123         & itypmpp
124      INTEGER, DIMENSION(:), ALLOCATABLE :: &
125         & iobsi,    &
126         & iobsj,    &
127         & iproc,    &
128         & iindx,    &
129         & ifileidx, &
130         & isurfidx
131      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
132         & zphi, &
133         & zlam
134      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
135         & zdat
136      REAL(wp), DIMENSION(knumfiles) :: &
137         & djulini, &
138         & djulend
139      LOGICAL :: llvalprof
140      TYPE(obfbdata), POINTER, DIMENSION(:) :: &
141         & inpfiles
142      CHARACTER(LEN=256) :: cout1  ! Diagnostic output line
143
144      ! Local initialization
145      iobs = 0
146
147      !-----------------------------------------------------------------------
148      ! Count the number of files needed and allocate the obfbdata type
149      !-----------------------------------------------------------------------
150
151      inobf = knumfiles
152
153      ALLOCATE( inpfiles(inobf) )
154
155      iadd  = 0
156      iextr = 0
157
158      surf_files : DO jj = 1, inobf
159
160         !---------------------------------------------------------------------
161         ! Prints
162         !---------------------------------------------------------------------
163         IF(lwp) THEN
164            WRITE(numout,*)
165            WRITE(numout,*) ' obs_rea_surf : Reading from file = ', &
166               & TRIM( TRIM( cdfilenames(jj) ) )
167            WRITE(numout,*) ' ~~~~~~~~~~~'
168            WRITE(numout,*)
169         ENDIF
170
171         !---------------------------------------------------------------------
172         !  Initialization: Open file and get dimensions only
173         !---------------------------------------------------------------------
174
175         iflag = nf90_open( TRIM( TRIM( cdfilenames(jj) ) ), nf90_nowrite, &
176            &                      i_file_id )
177
178         IF ( iflag /= nf90_noerr ) THEN
179
180            IF ( ldignmis ) THEN
181               inpfiles(jj)%nobs = 0
182               CALL ctl_warn( 'File ' // TRIM( TRIM( cdfilenames(jj) ) ) // &
183                  &           ' not found' )
184            ELSE
185               CALL ctl_stop( 'File ' // TRIM( TRIM( cdfilenames(jj) ) ) // &
186                  &           ' not found' )
187            ENDIF
188
189         ELSE 
190
191            !------------------------------------------------------------------
192            !  Close the file since it is opened in read_obfbdata
193            !------------------------------------------------------------------
194
195            iflag = nf90_close( i_file_id )
196
197            !------------------------------------------------------------------
198            !  Read the surface file into inpfiles
199            !------------------------------------------------------------------
200            CALL init_obfbdata( inpfiles(jj) )
201            CALL read_obfbdata( TRIM( cdfilenames(jj) ), inpfiles(jj), &
202               &                ldgrid = .TRUE. )
203
204            IF ( inpfiles(jj)%nvar /= kvars ) THEN
205               CALL ctl_stop( 'Feedback format error: ', &
206                  &           ' unexpected number of vars in feedback file', &
207                  &           TRIM(cdfilenames(jj)) )
208            ENDIF
209
210            IF ( ldmod .AND. ( inpfiles(jj)%nadd == 0 ) ) THEN
211               CALL ctl_stop( 'Model not in input data in', &
212                  &           TRIM(cdfilenames(jj)) )
213               RETURN
214            ENDIF
215
216            IF ( (iextr > 0) .AND. (inpfiles(jj)%next /= iextr) ) THEN
217               CALL ctl_stop( 'Number of extra variables not consistent', &
218                  &           ' with previous files for this type in', &
219                  &           TRIM(cdfilenames(jj)) )
220            ELSE
221               iextr = inpfiles(jj)%next
222            ENDIF
223
224            ! Ignore model counterpart
225            iaddin = inpfiles(jj)%nadd
226            DO ji = 1, iaddin
227               IF ( TRIM(inpfiles(jj)%caddname(ji)) == 'Hx' ) THEN
228                  iaddin = iaddin - 1
229                  EXIT
230               ENDIF
231            END DO
232            IF ( ldmod .AND. ( inpfiles(jj)%nadd == iaddin ) ) THEN
233               CALL ctl_stop( 'Model not in input data', &
234                  &           TRIM(cdfilenames(jj)) )
235            ENDIF
236
237            IF ( (iadd > 0) .AND. (iaddin /= iadd) ) THEN
238               CALL ctl_stop( 'Number of additional variables not consistent', &
239                  &           ' with previous files for this type in', &
240                  &           TRIM(cdfilenames(jj)) )
241            ELSE
242               iadd = iaddin
243            ENDIF
244
245            IF ( jj == 1 ) THEN
246               ALLOCATE( clvarsin( inpfiles(jj)%nvar ) )
247               ALLOCATE( cllongin( inpfiles(jj)%nvar ) )
248               ALLOCATE( clunitin( inpfiles(jj)%nvar ) )
249               ALLOCATE( clgridin( inpfiles(jj)%nvar ) )
250               DO ji = 1, inpfiles(jj)%nvar
251                 clvarsin(ji) = inpfiles(jj)%cname(ji)
252                 cllongin(ji) = inpfiles(jj)%coblong(ji)
253                 clunitin(ji) = inpfiles(jj)%cobunit(ji)
254                 clgridin(ji) = inpfiles(jj)%cgrid(ji)
255                 IF ( clvarsin(ji) /= cdvars(ji) ) THEN
256                    CALL ctl_stop( 'Feedback file variables do not match', &
257                        &           ' expected variable names for this type' )
258                 ENDIF
259               END DO
260               IF ( iadd > 0 ) THEN
261                  ALLOCATE( claddvarsin( iadd ) )
262                  ALLOCATE( claddlongin( iadd, inpfiles(jj)%nvar ) )
263                  ALLOCATE( claddunitin( iadd, inpfiles(jj)%nvar ) )
264                  jadd = 0
265                  DO ji = 1, inpfiles(jj)%nadd
266                    IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
267                       jadd = jadd + 1
268                       claddvarsin(jadd) = inpfiles(jj)%caddname(ji)
269                       DO jk = 1, inpfiles(jj)%nvar
270                          claddlongin(jadd,jk) = inpfiles(jj)%caddlong(ji,jk)
271                          claddunitin(jadd,jk) = inpfiles(jj)%caddunit(ji,jk)
272                       END DO
273                    ENDIF
274                  END DO
275               ENDIF
276               IF ( iextr > 0 ) THEN
277                  ALLOCATE( clextvarsin( iextr ) )
278                  ALLOCATE( clextlongin( iextr ) )
279                  ALLOCATE( clextunitin( iextr ) )
280                  DO ji = 1, iextr
281                    clextvarsin(ji) = inpfiles(jj)%cextname(ji)
282                    clextlongin(ji) = inpfiles(jj)%cextlong(ji)
283                    clextunitin(ji) = inpfiles(jj)%cextunit(ji)
284                  END DO
285               ENDIF
286            ELSE
287               DO ji = 1, inpfiles(jj)%nvar
288                  IF ( inpfiles(jj)%cname(ji) /= clvarsin(ji) ) THEN
289                     CALL ctl_stop( 'Feedback file variables not consistent', &
290                        &           ' with previous files for this type in', &
291                        &           TRIM(cdfilenames(jj)) )
292                  ENDIF
293               END DO
294               IF ( iadd > 0 ) THEN
295                  jadd = 0
296                  DO ji = 1, inpfiles(jj)%nadd
297                     IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
298                        jadd = jadd + 1
299                        IF ( inpfiles(jj)%caddname(ji) /= claddvarsin(jadd) ) THEN
300                           CALL ctl_stop( 'Feedback file additional variables not consistent', &
301                              &           ' with previous files for this type in', &
302                              &           TRIM(cdfilenames(jj)) )
303                        ENDIF
304                     ENDIF
305                  END DO
306               ENDIF
307               IF ( iextr > 0 ) THEN
308                  DO ji = 1, iextr
309                     IF ( inpfiles(jj)%cextname(ji) /= clextvarsin(ji) ) THEN
310                        CALL ctl_stop( 'Feedback file extra variables not consistent', &
311                           &           ' with previous files for this type in', &
312                           &           TRIM(cdfilenames(jj)) )
313                     ENDIF
314                  END DO
315               ENDIF
316
317            ENDIF
318
319            IF (lwp) WRITE(numout,*)'Observation file contains ',inpfiles(jj)%nobs,' observations'
320
321            !------------------------------------------------------------------
322            !  Change longitude (-180,180)
323            !------------------------------------------------------------------
324
325            DO ji = 1, inpfiles(jj)%nobs
326
327               IF ( inpfiles(jj)%plam(ji) < -180. ) &
328                  &   inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) + 360.
329
330               IF ( inpfiles(jj)%plam(ji) >  180. ) &
331                  &   inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) - 360.
332
333            END DO
334
335            !------------------------------------------------------------------
336            !  Calculate the date  (change eventually)
337            !------------------------------------------------------------------
338            clrefdate=inpfiles(jj)%cdjuldref(1:8)
339            READ(clrefdate,'(I8)') irefdate(jj)
340
341            CALL ddatetoymdhms( ddobsini, iyea, imon, iday, ihou, imin, isec )
342            CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulini(jj), &
343               &           krefdate = irefdate(jj) )
344            CALL ddatetoymdhms( ddobsend, iyea, imon, iday, ihou, imin, isec )
345            CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulend(jj), &
346               &           krefdate = irefdate(jj) )
347
348            IF ( ldnightav ) THEN
349
350               IF ( lwp ) THEN
351                  WRITE(numout,*)'Resetting time of night-time averaged observations', &
352                     &             ' to the end of the day'
353               ENDIF
354
355               DO ji = 1, inpfiles(jj)%nobs
356                  !  for night-time averaged data force the time
357                  !  to be the last time-step of the day, but still within the day.
358                  IF ( inpfiles(jj)%ptim(ji) >= 0. ) THEN
359                     inpfiles(jj)%ptim(ji) = &
360                        & INT(inpfiles(jj)%ptim(ji)) + 0.9999
361                  ELSE
362                     inpfiles(jj)%ptim(ji) = &
363                        & INT(inpfiles(jj)%ptim(ji)) - 0.0001
364                  ENDIF
365               END DO
366            ENDIF
367
368            IF ( inpfiles(jj)%nobs > 0 ) THEN
369               inpfiles(jj)%iproc = -1
370               inpfiles(jj)%iobsi = -1
371               inpfiles(jj)%iobsj = -1
372            ENDIF
373            inowin = 0
374            DO ji = 1, inpfiles(jj)%nobs
375               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
376                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
377                  inowin = inowin + 1
378               ENDIF
379            END DO
380            ALLOCATE( zlam(inowin)  )
381            ALLOCATE( zphi(inowin)  )
382            ALLOCATE( iobsi(inowin) )
383            ALLOCATE( iobsj(inowin) )
384            ALLOCATE( iproc(inowin) )
385            inowin = 0
386            DO ji = 1, inpfiles(jj)%nobs
387               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
388                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
389                  inowin = inowin + 1
390                  zlam(inowin) = inpfiles(jj)%plam(ji)
391                  zphi(inowin) = inpfiles(jj)%pphi(ji)
392               ENDIF
393            END DO
394
395            CALL obs_grid_search( inowin, zlam, zphi, iobsi, iobsj, iproc, 'T' )
396
397            inowin = 0
398            DO ji = 1, inpfiles(jj)%nobs
399               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
400                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
401                  inowin = inowin + 1
402                  inpfiles(jj)%iproc(ji,1) = iproc(inowin)
403                  inpfiles(jj)%iobsi(ji,1) = iobsi(inowin)
404                  inpfiles(jj)%iobsj(ji,1) = iobsj(inowin)
405               ENDIF
406            END DO
407            DEALLOCATE( zlam, zphi, iobsi, iobsj, iproc )
408
409            DO ji = 1, inpfiles(jj)%nobs
410               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
411                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
412                  IF ( nproc == 0 ) THEN
413                     IF ( inpfiles(jj)%iproc(ji,1) >  nproc ) CYCLE
414                  ELSE
415                     IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
416                  ENDIF
417                  llvalprof = .FALSE.
418                  IF ( .NOT. BTEST(inpfiles(jj)%ivlqc(1,ji,1),2) ) THEN
419                     iobs = iobs + 1
420                  ENDIF
421               ENDIF
422            END DO
423
424         ENDIF
425
426      END DO surf_files
427
428      !-----------------------------------------------------------------------
429      ! Get the time ordered indices of the input data
430      !-----------------------------------------------------------------------
431
432      !---------------------------------------------------------------------
433      !  Loop over input data files to count total number of profiles
434      !---------------------------------------------------------------------
435      iobstot = 0
436      DO jj = 1, inobf
437         DO ji = 1, inpfiles(jj)%nobs
438            IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
439               & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
440               iobstot = iobstot + 1
441            ENDIF
442         END DO
443      END DO
444
445      ALLOCATE( iindx(iobstot), ifileidx(iobstot), &
446         &      isurfidx(iobstot), zdat(iobstot) )
447      jk = 0
448      DO jj = 1, inobf
449         DO ji = 1, inpfiles(jj)%nobs
450            IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
451               & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
452               jk = jk + 1
453               ifileidx(jk) = jj
454               isurfidx(jk) = ji
455               zdat(jk)     = inpfiles(jj)%ptim(ji)
456            ENDIF
457         END DO
458      END DO
459      CALL sort_dp_indx( iobstot, &
460         &               zdat,     &
461         &               iindx   )
462
463      CALL obs_surf_alloc( surfdata, iobs, kvars, kadd+iadd, kextr+iextr, kstp, jpi, jpj )
464
465      ! Read obs/positions, QC, all variable and assign to surfdata
466
467      iobs = 0
468
469      surfdata%cvars(:)  = clvarsin(:)
470      surfdata%clong(:)  = cllongin(:)
471      surfdata%cunit(:)  = clunitin(:)
472      surfdata%cgrid(:)  = clgridin(:)
473      IF ( iadd > 0 ) THEN
474         surfdata%caddvars(kadd+1:)   = claddvarsin(:)
475         surfdata%caddlong(kadd+1:,:) = claddlongin(:,:)
476         surfdata%caddunit(kadd+1:,:) = claddunitin(:,:)
477      ENDIF
478      IF ( iextr > 0 ) THEN
479         surfdata%cextvars(kextr+1:) = clextvarsin(:)
480         surfdata%cextlong(kextr+1:) = clextlongin(:)
481         surfdata%cextunit(kextr+1:) = clextunitin(:)
482      ENDIF
483
484      ityp   (:) = 0
485      itypmpp(:) = 0
486
487      ioserrcount = 0
488
489      DO jk = 1, iobstot
490
491         jj = ifileidx(iindx(jk))
492         ji = isurfidx(iindx(jk))
493         IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND.  &
494            & ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
495
496            IF ( nproc == 0 ) THEN
497               IF ( inpfiles(jj)%iproc(ji,1) >  nproc ) CYCLE
498            ELSE
499               IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
500            ENDIF
501
502            ! Set observation information
503
504            IF ( .NOT. BTEST(inpfiles(jj)%ivlqc(1,ji,1),2) ) THEN
505
506               iobs = iobs + 1
507
508               CALL jul2greg( isec,                   &
509                  &           imin,                   &
510                  &           ihou,                   &
511                  &           iday,                   &
512                  &           imon,                   &
513                  &           iyea,                   &
514                  &           inpfiles(jj)%ptim(ji), &
515                  &           irefdate(jj) )
516
517
518               ! Surface time coordinates
519               surfdata%nyea(iobs) = iyea
520               surfdata%nmon(iobs) = imon
521               surfdata%nday(iobs) = iday
522               surfdata%nhou(iobs) = ihou
523               surfdata%nmin(iobs) = imin
524
525               ! Surface space coordinates
526               surfdata%rlam(iobs) = inpfiles(jj)%plam(ji)
527               surfdata%rphi(iobs) = inpfiles(jj)%pphi(ji)
528
529               ! Coordinate search parameters
530               surfdata%mi  (iobs) = inpfiles(jj)%iobsi(ji,1)
531               surfdata%mj  (iobs) = inpfiles(jj)%iobsj(ji,1)
532
533               ! WMO number
534               surfdata%cwmo(iobs) = inpfiles(jj)%cdwmo(ji)
535
536               ! Instrument type
537               READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 901 ) itype
538901            IF ( ios /= 0 ) THEN
539                  IF (ioserrcount == 0) THEN
540                     CALL ctl_warn ( 'Problem converting an instrument type ', &
541                        &            'to integer. Setting type to zero' )
542                  ENDIF
543                  ioserrcount = ioserrcount + 1
544                  itype = 0
545               ENDIF
546               surfdata%ntyp(iobs) = itype
547               IF ( itype < jpsurfmaxtype + 1 ) THEN
548                  ityp(itype+1) = ityp(itype+1) + 1
549               ELSE
550                  IF(lwp)WRITE(numout,*)'WARNING:Increase jpsurfmaxtype in ',&
551                     &                  cpname
552               ENDIF
553
554               ! Bookkeeping data to match observations
555               surfdata%nsidx(iobs) = iobs
556               surfdata%nsfil(iobs) = iindx(jk)
557
558               DO jvar = 1, kvars
559
560                  ! QC flags
561! WHY IS THIS NOT A FUNCTION OF NUM VARS?
562                  surfdata%nqc(iobs) = inpfiles(jj)%ivqc(ji,jvar)
563
564                  ! Observed value
565                  surfdata%robs(iobs,jvar) = inpfiles(jj)%pob(1,ji,jvar)
566
567! THIS NEEDS SORTING
568!                  ! Model and MDT is set to fbrmdi unless read from file
569!                  IF ( ldmod ) THEN
570!                     surfdata%rmod(iobs,jvar) = inpfiles(jj)%padd(1,ji,1,1)
571!                     IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) THEN
572!                        surfdata%rext(iobs,1) = inpfiles(jj)%padd(1,ji,2,1)
573!                        surfdata%rext(iobs,2) = inpfiles(jj)%pext(1,ji,1)
574!                     ENDIF
575!                   ELSE
576!                     surfdata%rmod(iobs,jvar) = fbrmdi
577!                     IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) surfdata%rext(iobs,:) = fbrmdi
578!                  ENDIF
579
580                  ! Additional variables
581                  surfdata%rmod(iobs,jvar) = fbrmdi
582                  IF ( iadd > 0 ) THEN
583                     jadd2 = 0
584                     DO jadd = 1, inpfiles(jj)%nadd
585                        IF ( TRIM(inpfiles(jj)%caddname(jadd)) == 'Hx' ) THEN
586                           IF ( ldmod ) THEN
587                              surfdata%rmod(iobs,jvar) = inpfiles(jj)%padd(1,ji,jadd,jvar)
588                           ENDIF
589                        ELSE
590                           jadd2 = jadd2 + 1
591                           surfdata%radd(iobs,kadd+jadd2,jvar) = &
592                              &                inpfiles(jj)%padd(1,ji,jadd,jvar)
593                        ENDIF
594                     END DO
595                  ENDIF
596
597               END DO
598                 
599               ! Extra variables
600               IF ( iextr > 0 ) THEN
601                  DO jext = 1, iextr
602                     surfdata%rext(iobs,kextr+jext) = inpfiles(jj)%pext(1,ji,jext)
603                  END DO
604               ENDIF
605            ENDIF
606         ENDIF
607
608      END DO
609
610      !-----------------------------------------------------------------------
611      ! Sum up over processors
612      !-----------------------------------------------------------------------
613
614      CALL obs_mpp_sum_integer( iobs, iobsmpp )
615      CALL obs_mpp_sum_integers( ityp, itypmpp, jpsurfmaxtype + 1 )
616
617      !-----------------------------------------------------------------------
618      ! Output number of observations.
619      !-----------------------------------------------------------------------
620      IF (lwp) THEN
621         DO jvar = 1, surfdata%nvar       
622            IF ( jvar == 1 ) THEN
623               cout1=TRIM(surfdata%cvars(1))                 
624            ELSE
625               WRITE(cout1,'(A,A1,A)') TRIM(cout1), '/', TRIM(surfdata%cvars(jvar))           
626            ENDIF
627         END DO
628 
629         WRITE(numout,*)
630         WRITE(numout,'(1X,A)')TRIM( cout1 )//' data'
631         WRITE(numout,'(1X,A)')'--------------'
632         DO jj = 1,8
633            IF ( itypmpp(jj) > 0 ) THEN
634               WRITE(numout,'(1X,A4,I4,A3,I10)')'Type ', jj,' = ',itypmpp(jj)
635            ENDIF
636         END DO
637         WRITE(numout,'(1X,A)') &
638            & '---------------------------------------------------------------'
639         WRITE(numout,'(1X,A,I8)') &
640            & 'Total data for variable '//TRIM( cout1 )// &
641            & '           = ', iobsmpp
642         WRITE(numout,'(1X,A)') &
643            & '---------------------------------------------------------------'
644         WRITE(numout,*)
645
646      ENDIF
647
648      !-----------------------------------------------------------------------
649      ! Deallocate temporary data
650      !-----------------------------------------------------------------------
651      DEALLOCATE( ifileidx, isurfidx, zdat, clvarsin, &
652         &        cllongin, clunitin, clgridin )
653      IF ( iadd > 0 ) THEN
654         DEALLOCATE( claddvarsin, claddlongin, claddunitin)
655      ENDIF
656      IF ( iextr > 0 ) THEN
657         DEALLOCATE( clextvarsin, clextlongin, clextunitin )
658      ENDIF
659
660      !-----------------------------------------------------------------------
661      ! Deallocate input data
662      !-----------------------------------------------------------------------
663      DO jj = 1, inobf
664         IF ( inpfiles(jj)%lalloc ) THEN
665            CALL dealloc_obfbdata( inpfiles(jj) )
666         ENDIF
667      END DO
668      DEALLOCATE( inpfiles )
669
670   END SUBROUTINE obs_rea_surf
671
672END MODULE obs_read_surf
Note: See TracBrowser for help on using the repository browser.