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_sla.F90 in branches/UKMO/dev_r4650_general_vert_coord_obsoper_removetides/NEMOGCM/NEMO/OPA_SRC/OBS – NEMO

source: branches/UKMO/dev_r4650_general_vert_coord_obsoper_removetides/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_sla.F90 @ 7375

Last change on this file since 7375 was 7375, checked in by mattmartin, 5 years ago

Added code to branch to allow sla model counterpart to be a time-average.

File size: 21.5 KB
Line 
1MODULE obs_read_sla
2   !!======================================================================
3   !!                       ***  MODULE obs_read_sla  ***
4   !! Observation diagnostics: Read the along track SLA data from
5   !!                          the AVISO/CLS database or any SLA data
6   !!                          from feedback files
7   !!======================================================================
8
9   !!----------------------------------------------------------------------
10   !!   obs_rea_sla : Driver for reading SLA data from the AVISO/CLS/feedback
11   !!----------------------------------------------------------------------
12
13   !! * Modules used   
14   USE par_kind                 ! Precision variables
15   USE in_out_manager           ! I/O manager
16   USE dom_oce                  ! Ocean space and time domain variables
17   USE obs_mpp                  ! MPP support routines for observation diagnostics
18   USE julian                   ! Julian date routines
19   USE obs_utils                ! Observation operator utility functions
20   USE obs_grid                 ! Grid search
21   USE obs_sort                 ! Sorting observation arrays
22   USE obs_surf_def             ! Surface observation definitions
23   USE obs_types                ! Observation type definitions
24   USE obs_sla_io               ! I/O for sla files
25   USE netcdf                   ! NetCDF library
26
27   IMPLICIT NONE
28
29   !! * Routine accessibility
30   PRIVATE
31
32   PUBLIC obs_rea_sla  ! Read the SLA observations from the AVISO/SLA database
33
34   !!----------------------------------------------------------------------
35   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
36   !! $Id$
37   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
38   !!----------------------------------------------------------------------
39
40CONTAINS
41
42   SUBROUTINE obs_rea_sla( kformat, &
43      &                    sladata, knumfiles, cfilenames, &
44      &                    kvars, kextr, kstp, ddobsini, ddobsend, &
45      &                    ldignmis, ldmod, MeanPeriodHours, ldobstd )
46      !!---------------------------------------------------------------------
47      !!
48      !!                   *** ROUTINE obs_rea_sla ***
49      !!
50      !! ** Purpose : Read from file the SLA data
51      !!
52      !! ** Method  : Depending on kformat either AVISO or
53      !!              feedback data files are read
54      !!
55      !! ** Action  :
56      !!
57      !!
58      !! History : 
59      !!      ! :  2009-01 (K. Mogensen) Initial version based on old versions
60      !!----------------------------------------------------------------------
61      !! * Modules used
62
63      !! * Arguments
64      INTEGER :: kformat    ! Format of input data
65      !                     ! 0: Feedback
66      !                     ! 1: AVISO
67      TYPE(obs_surf), INTENT(INOUT) :: sladata    ! SLA data to be read
68      INTEGER, INTENT(IN) :: knumfiles      ! Number of files to read in
69      CHARACTER(LEN=128), INTENT(IN) ::  &
70         & cfilenames(knumfiles) ! File names to read in
71      INTEGER, INTENT(IN) :: kvars     ! Number of variables in sladata
72      INTEGER, INTENT(IN) :: kextr     ! Number of extra fields for each var in sladata
73      INTEGER, INTENT(IN) :: kstp      ! Ocean time-step index
74      LOGICAL, INTENT(IN) :: ldignmis    ! Ignore missing files
75      LOGICAL, INTENT(IN) :: ldmod       ! Initialize model from input data
76      LOGICAL, INTENT(INOUT), optional :: &
77         & ldobstd        ! Read observation standard deviation from fb. file
78      REAL(KIND=dp), INTENT(IN) :: ddobsini  ! Obs. ini time in YYYYMMDD.HHMMSS
79      REAL(KIND=dp), INTENT(IN) :: ddobsend  ! Obs. end time in YYYYMMDD.HHMMSS
80      REAL(KIND=wp), INTENT(IN) :: MeanPeriodHours ! Averaging period in hours
81         
82      !! * Local declarations
83      CHARACTER(LEN=11), PARAMETER :: cpname='obs_rea_sla'
84      INTEGER :: ji
85      INTEGER :: jj
86      INTEGER :: jk
87      INTEGER :: iflag
88      INTEGER :: inobf
89      INTEGER :: i_file_id
90      INTEGER :: inowin
91      INTEGER :: iyea
92      INTEGER :: imon
93      INTEGER :: iday
94      INTEGER :: ihou
95      INTEGER :: imin
96      INTEGER :: isec
97      INTEGER, DIMENSION(knumfiles) :: &
98         & irefdate
99      INTEGER :: iobsmpp
100      INTEGER, DIMENSION(imaxmissions+1) :: &
101         & ityp, &
102         & itypmpp
103      INTEGER, DIMENSION(:), ALLOCATABLE :: &
104         & iobsi,    &
105         & iobsj,    &
106         & iproc,    &
107         & iindx,    &
108         & ifileidx, &
109         & islaidx
110      INTEGER :: itype 
111      REAL(wp), DIMENSION(:), ALLOCATABLE :: &
112         & zphi, &
113         & zlam
114      real(wp), DIMENSION(:), ALLOCATABLE :: &
115         & zdat
116      LOGICAL :: llvalprof
117      LOGICAL :: llobstd
118      TYPE(obfbdata), POINTER, DIMENSION(:) :: &
119         & inpfiles
120      real(wp), DIMENSION(knumfiles) :: &
121         & djulini, &
122         & djulend
123      INTEGER, DIMENSION(knumfiles) :: &
124         & iobspos, &
125         & iobcpos
126      INTEGER :: iobs
127      INTEGER :: iobstot
128      INTEGER :: ios
129      INTEGER :: ioserrcount
130      CHARACTER(len=8) :: cl_refdate
131   
132      ! Local initialization
133      iobs = 0
134      IF ( PRESENT(ldobstd) ) THEN
135         IF (.NOT.ldmod) THEN
136            llobstd = .false.
137         ELSE
138            llobstd = ldobstd
139         ENDIF
140      ELSE
141         llobstd = .FALSE.
142      ENDIF
143 
144      !-----------------------------------------------------------------------
145      ! Check that the model part is just with feedback data files
146      !-----------------------------------------------------------------------
147      IF ( ldmod .AND. ( kformat /= 0 ) ) THEN
148         CALL ctl_stop( 'Model can only be read from feedback data' )
149         RETURN
150      ENDIF
151
152      !-----------------------------------------------------------------------
153      ! Check that the prescribed obs err is just with feedback data files
154      !-----------------------------------------------------------------------
155      IF ( llobstd .AND. ( kformat /= 0 ) ) THEN
156         CALL ctl_stop( 'Observation error can only be read from feedback files' )
157         RETURN
158      ENDIF
159
160      !-----------------------------------------------------------------------
161      ! Count the number of files needed and allocate the obfbdata type
162      !-----------------------------------------------------------------------
163     
164      inobf = knumfiles
165     
166      ALLOCATE( inpfiles(inobf) )
167
168      sla_files : DO jj = 1, inobf
169         
170         !---------------------------------------------------------------------
171         ! Prints
172         !---------------------------------------------------------------------
173         IF(lwp) THEN
174            WRITE(numout,*)
175            WRITE(numout,*) ' obs_rea_sla : Reading from file = ', &
176               & TRIM( TRIM( cfilenames(jj) ) )
177            WRITE(numout,*) ' ~~~~~~~~~~~'
178            WRITE(numout,*)
179         ENDIF
180
181         !---------------------------------------------------------------------
182         !  Initialization: Open file and get dimensions only
183         !---------------------------------------------------------------------
184         
185         iflag = nf90_open( TRIM( TRIM( cfilenames(jj) ) ), nf90_nowrite, &
186            &                      i_file_id )
187         
188         IF ( iflag /= nf90_noerr ) THEN
189
190            IF ( ldignmis ) THEN
191               inpfiles(jj)%nobs = 0
192               CALL ctl_warn( 'File ' // TRIM( TRIM( cfilenames(jj) ) ) // &
193                  &           ' not found' )
194            ELSE
195               CALL ctl_stop( 'File ' // TRIM( TRIM( cfilenames(jj) ) ) // &
196                  &           ' not found' )
197            ENDIF
198
199         ELSE 
200           
201            !------------------------------------------------------------------
202            !  Close the file since it is opened in read_proffile
203            !------------------------------------------------------------------
204           
205            iflag = nf90_close( i_file_id )
206
207            !------------------------------------------------------------------
208            !  Read the profile file into inpfiles
209            !------------------------------------------------------------------
210            IF ( kformat == 0 ) THEN
211               CALL init_obfbdata( inpfiles(jj) )
212               IF(lwp) THEN
213                  WRITE(numout,*)
214                  WRITE(numout,*)'Reading from feedback file :', &
215                     &           TRIM( cfilenames(jj) )
216               ENDIF
217               CALL read_obfbdata( TRIM( cfilenames(jj) ), inpfiles(jj), &
218                  &                ldgrid = .TRUE. )
219               IF ( inpfiles(jj)%nvar < 1 ) THEN
220                  CALL ctl_stop( 'Feedback format error' )
221                  RETURN
222               ENDIF
223               IF ( TRIM(inpfiles(jj)%cname(1)) /= 'SLA' ) THEN
224                  CALL ctl_stop( 'Feedback format error' )
225                  RETURN
226               ENDIF
227               IF ( ldmod .AND. ( ( inpfiles(jj)%nadd == 0 ) .OR.&
228                  &               ( inpfiles(jj)%next == 0 ) ) ) THEN
229                  CALL ctl_stop( 'Model not in input data' )
230                  RETURN
231               ENDIF
232            ELSEIF ( kformat == 1) THEN
233               CALL read_avisofile( TRIM( cfilenames(jj) ), inpfiles(jj), &
234               &                    numout, lwp, .TRUE. )
235            ELSE
236               CALL ctl_stop( 'File format unknown' )
237            ENDIF
238
239            !------------------------------------------------------------------
240            !  Change longitude (-180,180)
241            !------------------------------------------------------------------
242
243            DO ji = 1, inpfiles(jj)%nobs 
244
245               IF ( inpfiles(jj)%plam(ji) < -180. ) &
246                  &   inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) + 360.
247
248               IF ( inpfiles(jj)%plam(ji) >  180. ) &
249                  &   inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) - 360.
250
251            END DO
252
253            !------------------------------------------------------------------
254            !  Calculate the date  (change eventually)
255            !------------------------------------------------------------------
256            cl_refdate=inpfiles(jj)%cdjuldref(1:8)
257            READ(cl_refdate,'(I8)') irefdate(jj)
258           
259            CALL ddatetoymdhms( ddobsini, iyea, imon, iday, ihou, imin, isec )
260            CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulini(jj), &
261               &           krefdate = irefdate(jj) )
262            CALL ddatetoymdhms( ddobsend, iyea, imon, iday, ihou, imin, isec )
263            CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulend(jj), &
264               &           krefdate = irefdate(jj) )
265            IF ( inpfiles(jj)%nobs > 0 ) THEN
266               inpfiles(jj)%iproc = -1
267               inpfiles(jj)%iobsi = -1
268               inpfiles(jj)%iobsj = -1
269            ENDIF
270           
271            !If the observations are representing a time mean then set the time
272            !of the obs to the end of that meaning period relative to the start of the run
273            IF ( MeanPeriodHours > 0._wp ) THEN
274               IF (lwp) WRITE(numout,*)'Adjusting time of sla obs to end of meaning period: ', djulini(jj) + (MeanPeriodHours/24.)
275               DO ji = 1, inpfiles(jj)%nobs
276                  inpfiles(jj)%ptim(ji) = &
277                        & djulini(jj) + (MeanPeriodHours/24.)
278               END DO
279            ENDIF
280
281            inowin = 0
282            DO ji = 1, inpfiles(jj)%nobs
283               IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
284               IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
285               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
286                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
287                  inowin = inowin + 1
288               ENDIF
289            END DO
290            ALLOCATE( zlam(inowin)  )
291            ALLOCATE( zphi(inowin)  )
292            ALLOCATE( iobsi(inowin) )
293            ALLOCATE( iobsj(inowin) )
294            ALLOCATE( iproc(inowin) )
295            inowin = 0
296            DO ji = 1, inpfiles(jj)%nobs
297               IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
298               IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
299               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
300                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
301                  inowin = inowin + 1
302                  zlam(inowin) = inpfiles(jj)%plam(ji)
303                  zphi(inowin) = inpfiles(jj)%pphi(ji)
304               ENDIF
305            END DO
306
307            CALL obs_grid_search( inowin, zlam, zphi, iobsi, iobsj, iproc, 'T' )
308
309            inowin = 0
310            DO ji = 1, inpfiles(jj)%nobs
311               IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
312               IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
313               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
314                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
315                  inowin = inowin + 1
316                  inpfiles(jj)%iproc(ji,1) = iproc(inowin)
317                  inpfiles(jj)%iobsi(ji,1) = iobsi(inowin)
318                  inpfiles(jj)%iobsj(ji,1) = iobsj(inowin)
319               ENDIF
320            END DO
321            DEALLOCATE( zlam, zphi, iobsi, iobsj, iproc )
322
323            DO ji = 1, inpfiles(jj)%nobs
324               IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
325               IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
326               IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
327                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
328                  IF ( nproc == 0 ) THEN
329                     IF ( inpfiles(jj)%iproc(ji,1) >  nproc ) CYCLE
330                  ELSE
331                     IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
332                  ENDIF
333                  llvalprof = .FALSE.
334                  IF ( ( inpfiles(jj)%ivlqc(1,ji,1) == 1 ) .OR. &
335                     & ( inpfiles(jj)%ivlqc(1,ji,1) == 2 ) ) THEN
336                     iobs = iobs + 1
337                  ENDIF
338               ENDIF
339            END DO
340
341         ENDIF
342     
343      END DO sla_files
344
345      IF (llobstd) THEN
346         
347         DO jj = 1, inobf
348            iobspos(jj) = -1
349            iobcpos(jj) = -1
350            DO ji = 1,inpfiles(jj)%nadd
351               IF ( TRIM(inpfiles(jj)%caddname(ji)) == 'OSTD' ) THEN
352                  iobspos(jj)=ji
353               ENDIF
354               IF ( TRIM(inpfiles(jj)%caddname(ji)) == 'OCNT' ) THEN
355                  iobcpos(jj)=ji
356               ENDIF
357            END DO
358         END DO
359         llobstd = ( ( MINVAL(iobspos) > 0 ) .AND. ( MINVAL(iobcpos) > 0 ) )
360         IF (llobstd) THEN
361            IF (lwp) WRITE(numout,*)'SLA superobs information present.'
362         ELSE
363            IF (lwp) WRITE(numout,*)'SLA superobs information not present.'
364         ENDIF
365
366      ENDIF
367     
368      !-----------------------------------------------------------------------
369      ! Get the time ordered indices of the input data
370      !-----------------------------------------------------------------------
371
372      !---------------------------------------------------------------------
373      !  Loop over input data files to count total number of profiles
374      !---------------------------------------------------------------------
375      iobstot = 0
376      DO jj = 1, inobf
377         DO ji = 1, inpfiles(jj)%nobs
378            IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
379            IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
380            IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
381               & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
382               iobstot = iobstot + 1
383            ENDIF
384         END DO
385      END DO
386
387      ALLOCATE( iindx(iobstot), ifileidx(iobstot), &
388         &      islaidx(iobstot), zdat(iobstot) )
389      jk = 0
390      DO jj = 1, inobf
391         DO ji = 1, inpfiles(jj)%nobs
392            IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
393            IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
394            IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND. &
395               & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN
396               jk = jk + 1
397               ifileidx(jk) = jj
398               islaidx(jk) = ji
399               zdat(jk)     = inpfiles(jj)%ptim(ji)
400            ENDIF
401         END DO
402      END DO
403      CALL sort_dp_indx( iobstot, &
404         &               zdat,     &
405         &               iindx   )
406     
407      CALL obs_surf_alloc( sladata, iobs, kvars, kextr, kstp, jpi, jpj )
408     
409      ! * Read obs/positions, QC, all variable and assign to sladata
410 
411      iobs = 0
412
413      ityp   (:) = 0
414      itypmpp(:) = 0
415
416      ioserrcount = 0
417     
418      DO jk = 1, iobstot
419         
420         jj = ifileidx(iindx(jk))
421         ji = islaidx(iindx(jk))
422
423         IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
424         IF ( inpfiles(jj)%ivqc(ji,1) > 2 ) CYCLE
425
426         IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND.  &
427            & ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
428           
429            IF ( nproc == 0 ) THEN
430               IF ( inpfiles(jj)%iproc(ji,1) >  nproc ) CYCLE
431            ELSE
432               IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
433            ENDIF
434           
435            ! Set observation information
436           
437            IF ( ( inpfiles(jj)%ivlqc(1,ji,1) == 1 ) .OR. &
438               & ( inpfiles(jj)%ivlqc(1,ji,1) == 2 ) ) THEN
439
440               iobs = iobs + 1
441
442               CALL jul2greg( isec,                   &
443                  &           imin,                   &
444                  &           ihou,                   &
445                  &           iday,                   &
446                  &           imon,                   &
447                  &           iyea,                   &
448                  &           inpfiles(jj)%ptim(ji), &
449                  &           irefdate(jj) )
450
451
452               ! SLA time coordinates
453               sladata%nyea(iobs) = iyea
454               sladata%nmon(iobs) = imon
455               sladata%nday(iobs) = iday
456               sladata%nhou(iobs) = ihou
457               sladata%nmin(iobs) = imin
458               
459               ! SLA space coordinates
460               sladata%rlam(iobs) = inpfiles(jj)%plam(ji)
461               sladata%rphi(iobs) = inpfiles(jj)%pphi(ji)
462
463               ! Coordinate search parameters
464               sladata%mi  (iobs) = inpfiles(jj)%iobsi(ji,1)
465               sladata%mj  (iobs) = inpfiles(jj)%iobsj(ji,1)
466               
467               ! Instrument type
468               READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 901 ) itype
469901            IF ( ios /= 0 ) THEN
470                  IF (ioserrcount == 0) CALL ctl_warn ( 'Problem converting an instrument type to integer. Setting type to zero' ) 
471                  ioserrcount = ioserrcount + 1
472                  itype = 0
473               ENDIF
474               sladata%ntyp(iobs) = itype
475               ityp(itype+1) = ityp(itype+1) + 1
476
477               ! Identifier
478               sladata%cwmo(iobs) = inpfiles(jj)%cdwmo(ji)
479
480               ! Bookkeeping data to match observations
481               sladata%nsidx(iobs) = iobs
482               sladata%nsfil(iobs) = iindx(jk)
483
484               ! QC flags
485               sladata%nqc(iobs) = inpfiles(jj)%ivqc(ji,1)
486
487               ! Observed value
488               sladata%robs(iobs,1) = inpfiles(jj)%pob(1,ji,1)
489
490
491               ! Model and MDT is set to fbrmdi unless read from file
492               IF ( ldmod ) THEN
493                  sladata%rmod(iobs,1) = inpfiles(jj)%padd(1,ji,1,1)
494                  sladata%rext(iobs,1) = inpfiles(jj)%padd(1,ji,2,1)
495                  sladata%rext(iobs,2) = inpfiles(jj)%pext(1,ji,1)
496                  IF (llobstd) THEN
497                     sladata%rext(iobs,3) = &
498                        & inpfiles(jj)%padd(1,ji,iobspos(jj),1)
499                     sladata%rext(iobs,4) = &
500                        & inpfiles(jj)%padd(1,ji,iobcpos(jj),1)
501                  ENDIF
502               ELSE
503                  sladata%rmod(iobs,1) = fbrmdi
504                  sladata%rext(iobs,:) = fbrmdi
505               ENDIF
506
507            ENDIF
508         ENDIF
509
510      END DO
511
512      !-----------------------------------------------------------------------
513      ! Sum up over processors
514      !-----------------------------------------------------------------------
515     
516      CALL obs_mpp_sum_integer( iobs, iobsmpp )
517      CALL obs_mpp_sum_integers( ityp, itypmpp, imaxmissions + 1 )
518     
519      !-----------------------------------------------------------------------
520      ! Output number of observations.
521      !-----------------------------------------------------------------------
522      IF (lwp) THEN
523
524         WRITE(numout,*)
525         WRITE(numout,'(1X,A)')'Altimeter satellites'
526         WRITE(numout,'(1X,A)')'--------------------'
527         DO jj = 1,8
528            IF ( itypmpp(jj) > 0 ) THEN
529               WRITE(numout,'(1X,A38,A2,I10)')calttyp(jj-1),'= ',itypmpp(jj)
530            ENDIF
531         END DO
532         WRITE(numout,'(1X,A50)')'--------------------------------------------------'
533         WRITE(numout,'(1X,A40,I10)')'Total                                 = ',iobsmpp
534         WRITE(numout,*)
535
536      ENDIF
537
538      !-----------------------------------------------------------------------
539      ! Deallocate temporary data
540      !-----------------------------------------------------------------------
541      DEALLOCATE( ifileidx, islaidx, zdat )
542
543      !-----------------------------------------------------------------------
544      ! Deallocate input data
545      !-----------------------------------------------------------------------
546      DO jj = 1, inobf
547         CALL dealloc_obfbdata( inpfiles(jj) )
548      END DO
549      DEALLOCATE( inpfiles )
550
551      !-----------------------------------------------------------------------
552      ! Reset ldobstd if the data is present
553      !-----------------------------------------------------------------------
554      IF ( PRESENT(ldobstd) ) THEN
555         ldobstd = llobstd
556      ENDIF
557
558   END SUBROUTINE obs_rea_sla
559
560END MODULE obs_read_sla
Note: See TracBrowser for help on using the repository browser.