source: utils/tools/OBSTOOLS/src/obsvel_io.h90 @ 10841

Last change on this file since 10841 was 10841, checked in by cbricaud, 19 months ago

fix for ticket #2268

File size: 14.9 KB
Line 
1   !!----------------------------------------------------------------------
2   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
3   !! $Id: obsvel_io.h90 2287 2010-10-18 07:53:52Z smasson $
4   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
5   !!----------------------------------------------------------------------
6
7   SUBROUTINE read_taondbc( cdfilename, inpfile, kunit, ldwp, ldgrid )
8      !!---------------------------------------------------------------------
9      !!
10      !!                     ** ROUTINE read_enactfile **
11      !!
12      !! ** Purpose : Read from file the TAO data fro NDBC.
13      !!
14      !! ** Method  : The data file is a NetCDF file.
15      !!
16      !! ** Action  :
17      !!
18      !! ** Reference : http://tao.noaa.gov/tao/data_deliv/deliv_ndbc.shtml
19      !! History :
20      !!          ! 09-01 (K. Mogensen) Original version.
21      !!----------------------------------------------------------------------
22      !! * Arguments
23      CHARACTER(LEN=*) :: cdfilename ! Input filename
24      TYPE(obfbdata)   :: inpfile    ! Output obfbdata structure
25      INTEGER          :: kunit      ! Unit for output
26      LOGICAL          :: ldwp       ! Print info
27      LOGICAL          :: ldgrid     ! Save grid info in data structure
28      !! * Local declarations
29      INTEGER :: iobs                ! Number of observations
30      INTEGER :: ilev                ! Number of levels
31      INTEGER :: ilat                ! Number of latitudes
32      INTEGER :: ilon                ! Number of longtudes
33      INTEGER :: itim                ! Number of obs. times
34      INTEGER :: i_file_id
35      INTEGER :: i_dimid_id
36      INTEGER :: i_phi_id
37      INTEGER :: i_lam_id
38      INTEGER :: i_depth_id
39      INTEGER :: i_var_id
40      INTEGER :: i_time_id
41      INTEGER :: i_time2_id
42      INTEGER :: i_qc_var_id
43      CHARACTER(LEN=40) :: cl_fld_lam
44      CHARACTER(LEN=40) :: cl_fld_phi
45      CHARACTER(LEN=40) :: cl_fld_depth
46      CHARACTER(LEN=40) :: cl_fld_var_u
47      CHARACTER(LEN=40) :: cl_fld_var_v
48      CHARACTER(LEN=40) :: cl_fld_var_qc_uv1
49      CHARACTER(LEN=40) :: cl_fld_var_qc_uv2
50      CHARACTER(LEN=40) :: cl_fld_time
51      CHARACTER(LEN=40) :: cl_fld_time2
52      INTEGER :: ja
53      INTEGER :: jo
54      INTEGER :: jk
55      INTEGER :: jt
56      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: &
57         & zv,     &
58         & zu,     &
59         & zuv1qc, &
60         & zuv2qc
61      REAL(wp), ALLOCATABLE, DIMENSION(:) :: &
62         & zdep, &
63         & zlat, &
64         & zlon, &
65         & zjuld
66      REAL(wp) :: zl
67      INTEGER, ALLOCATABLE, DIMENSION(:) :: &
68         & itime, &
69         & itime2
70      CHARACTER(LEN=50) :: cdjulref
71      CHARACTER(LEN=12), PARAMETER :: cl_name = 'read_taondbc'
72      CHARACTER(len=1) :: cns, cew
73
74      !-----------------------------------------------------------------------
75      ! Initialization
76      !-----------------------------------------------------------------------
77      cl_fld_lam                 = 'lon'
78      cl_fld_phi                 = 'lat'
79      cl_fld_depth               = 'depth'
80      cl_fld_time                = 'time'
81      cl_fld_time2               = 'time2'
82
83      !-----------------------------------------------------------------------
84      ! Open file
85      !-----------------------------------------------------------------------
86
87      CALL chkerr( nf90_open( TRIM( cdfilename ), nf90_nowrite, &
88            &      i_file_id ),           cl_name, __LINE__ )
89
90      !-----------------------------------------------------------------------
91      ! Read the heading of the file
92      !-----------------------------------------------------------------------
93      IF(ldwp) WRITE(kunit,*)
94      IF(ldwp) WRITE(kunit,*) ' read_taondbc :'
95      IF(ldwp) WRITE(kunit,*) ' ~~~~~~~~~~~~'
96     
97      !---------------------------------------------------------------------
98      ! Read the number of observations and of levels to allocate array
99      !---------------------------------------------------------------------
100      CALL chkerr( nf90_inq_dimid        ( i_file_id, 'time', i_dimid_id ),        &
101         &         cl_name, __LINE__ )
102      CALL chkerr( nf90_inquire_dimension( i_file_id, i_dimid_id, len = itim ),    &
103         &         cl_name, __LINE__ )
104      CALL chkerr( nf90_inq_dimid        ( i_file_id, 'depth', i_dimid_id ),       &
105         &         cl_name, __LINE__ )
106      CALL chkerr( nf90_inquire_dimension( i_file_id, i_dimid_id, len = ilev ),    &
107         &         cl_name, __LINE__ )
108      CALL chkerr( nf90_inq_dimid        ( i_file_id, 'lat', i_dimid_id ),           &
109         &         cl_name, __LINE__ )
110      CALL chkerr( nf90_inquire_dimension( i_file_id, i_dimid_id, len = ilat ),    &
111         &         cl_name, __LINE__ )
112      CALL chkerr( nf90_inq_dimid        ( i_file_id, 'lon', i_dimid_id ),         &
113         &         cl_name, __LINE__ )
114      CALL chkerr( nf90_inquire_dimension( i_file_id, i_dimid_id, len = ilon ),    &
115         &         cl_name, __LINE__ )
116
117      iobs = itim * ilat * ilon
118      IF(ldwp)WRITE(kunit,*) '         No. of data records = ', iobs
119      IF(ldwp)WRITE(kunit,*) '         No. of levels       = ', ilev
120      IF(ldwp)WRITE(kunit,*)
121
122      !---------------------------------------------------------------------
123      ! Allocate arrays
124      !---------------------------------------------------------------------
125
126      CALL init_obfbdata( inpfile )
127      CALL alloc_obfbdata( inpfile, 2, iobs, ilev, 0, 0, ldgrid )
128      inpfile%cname(1) = 'UVEL'
129      inpfile%cname(2) = 'VVEL'
130      inpfile%coblong(1) = 'Zonal current'
131      inpfile%coblong(2) = 'Meridional current'
132      inpfile%cobunit(1) = 'Meters per second'
133      inpfile%cobunit(2) = 'Meters per second'
134
135      ALLOCATE( &
136         & zu(ilon,ilat,ilev,itim),     &
137         & zv(ilon,ilat,ilev,itim),     &
138         & zdep(ilev),                  &
139         & zuv1qc(ilon,ilat,ilev,itim), &
140         & zuv2qc(ilon,ilat,ilev,itim), &
141         & itime(itim),                 &
142         & itime2(itim),                &
143         & zlat(ilat),                  &
144         & zlon(ilon),                  &
145         & zjuld(itim)                  &
146         & )
147
148      !---------------------------------------------------------------------
149      ! Read the time/position variables
150      !---------------------------------------------------------------------
151     
152      CALL chkerr( nf90_inq_varid( i_file_id, cl_fld_time, i_time_id ),                               &
153         &         cl_name, __LINE__ )
154      CALL chkerr( nf90_get_var  ( i_file_id, i_time_id, itime ),                                     &
155         &         cl_name, __LINE__ )
156
157      CALL chkerr( nf90_inq_varid( i_file_id, cl_fld_time2, i_time2_id ),                             &
158         &         cl_name, __LINE__ )
159      CALL chkerr( nf90_get_var  ( i_file_id, i_time2_id, itime2 ),                                   &
160         &         cl_name, __LINE__ )
161     
162      CALL chkerr( nf90_inq_varid( i_file_id, cl_fld_depth, i_depth_id ),                             &
163            &         cl_name, __LINE__ )         
164      CALL chkerr( nf90_get_var  ( i_file_id, i_depth_id, zdep ),                                     &
165         &         cl_name, __LINE__ )
166     
167      CALL chkerr( nf90_inq_varid( i_file_id, cl_fld_phi, i_phi_id ),                                 &
168         &         cl_name, __LINE__ )
169      CALL chkerr( nf90_get_var  ( i_file_id, i_phi_id, zlat ),                                       &
170         &         cl_name, __LINE__ )
171     
172      CALL chkerr( nf90_inq_varid( i_file_id, cl_fld_lam, i_lam_id ),                                 &
173         &         cl_name, __LINE__ )
174      CALL chkerr( nf90_get_var  ( i_file_id, i_lam_id, zlon ),                                       &
175         &         cl_name, __LINE__ )
176     
177      !---------------------------------------------------------------------
178      ! Read the variables
179      !---------------------------------------------------------------------
180
181      ! ADCP format assumed
182      cl_fld_var_u = 'u_1205'
183      IF ( nf90_inq_varid( i_file_id, cl_fld_var_u, i_var_id ) /= nf90_noerr ) THEN
184         ! Try again with current meter format
185         cl_fld_var_u = 'U_320'
186         IF ( nf90_inq_varid( i_file_id, cl_fld_var_u, i_var_id ) /= nf90_noerr ) THEN
187            CALL fatal_error( 'Unknown format in read_taondbc', __LINE__ )
188         ENDIF
189      ENDIF
190      CALL chkerr( nf90_get_var  ( i_file_id, i_var_id, zu ),                                         &
191         &         cl_name, __LINE__ )
192     
193      ! ADCP format assumed
194      cl_fld_var_v = 'v_1206'
195      IF ( nf90_inq_varid( i_file_id, cl_fld_var_v, i_var_id ) /= nf90_noerr ) THEN
196         ! Try again with current meter format
197         cl_fld_var_v = 'V_321'
198         IF ( nf90_inq_varid( i_file_id, cl_fld_var_v, i_var_id ) /= nf90_noerr ) THEN
199            CALL fatal_error( 'Unknown format in read_taondbc', __LINE__ )
200         ENDIF
201      ENDIF
202      CALL chkerr( nf90_get_var  ( i_file_id, i_var_id, zv ),                                         &
203         &         cl_name, __LINE__ )
204
205      !---------------------------------------------------------------------
206      ! Read the QC attributes
207      !---------------------------------------------------------------------
208     
209      ! ADCP format assumed
210      cl_fld_var_qc_uv1 = 'QU_5205'
211      IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv1, i_qc_var_id ) /= nf90_noerr ) THEN
212         ! Try again with current meter format
213         cl_fld_var_qc_uv1 = 'QCS_5300'
214         IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv1, i_qc_var_id ) /= nf90_noerr ) THEN
215            ! Try again with high freq. current meter format
216            cl_fld_var_qc_uv1 = 'QCU_5320'
217            IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv1, i_qc_var_id ) /= nf90_noerr ) THEN
218               CALL fatal_error( 'Unknown format in read_taondbc', __LINE__ )
219            ENDIF
220         ENDIF
221      ENDIF
222      CALL chkerr( nf90_get_var  ( i_file_id, i_qc_var_id, zuv1qc),                                   &
223         &         cl_name, __LINE__ )
224
225      ! ADCP format assumed
226      cl_fld_var_qc_uv2 = 'QV_5206'
227      IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv2, i_qc_var_id ) /= nf90_noerr ) THEN
228         ! Try again with current meter format
229         cl_fld_var_qc_uv2 = 'QCD_5310'
230         IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv2, i_qc_var_id ) /= nf90_noerr ) THEN
231            ! Try again with high freq. current meter format
232            cl_fld_var_qc_uv2 = 'QCV_5321'
233            IF ( nf90_inq_varid( i_file_id, cl_fld_var_qc_uv2, i_qc_var_id ) /= nf90_noerr ) THEN
234               CALL fatal_error( 'Unknown format in read_taondbc', __LINE__ )
235            ENDIF
236         ENDIF
237      ENDIF
238      CALL chkerr( nf90_get_var  ( i_file_id, i_qc_var_id, zuv2qc),                                   &
239         &         cl_name, __LINE__ )
240
241      !---------------------------------------------------------------------
242      ! Close file
243      !---------------------------------------------------------------------
244
245      CALL chkerr( nf90_close( i_file_id ),           cl_name, __LINE__ )
246
247      !---------------------------------------------------------------------
248      ! Convert to to 19500101 based Julian date
249      !---------------------------------------------------------------------
250      DO jt = 1, itim
251         zjuld(jt) = REAL(itime(jt),wp) + REAL(itime2(jt),wp)/86400000.0_wp &
252            &           - 2433283.0_wp
253      END DO
254      inpfile%cdjuldref = '19500101000000'
255
256      !---------------------------------------------------------------------
257      ! Copy info to obfbdata structure
258      !---------------------------------------------------------------------
259
260      iobs = 0
261      DO jt = 1, itim
262         DO ja = 1, ilat
263            DO jo = 1, ilon
264               iobs = iobs + 1
265               zl = zlon(jo)
266               IF ( zl > 180.0_wp ) zl = zl - 360.0_wp
267               IF ( zl < 0 ) THEN
268                  cew = 'w'
269               ELSE
270                  cew = 'e'
271               ENDIF
272               IF ( zlat(jo) < 0 ) THEN
273                  cns = 's'
274               ELSE
275                  cns = 'n'
276               ENDIF
277               WRITE(inpfile%cdwmo(iobs),'(A1,I2.2,A1,I3.3)') &
278                  & cns, ABS(NINT(zlat(ja))), cew, ABS(NINT(zl))
279               DO jk = 1, ilev
280                  inpfile%pob(jk,iobs,1)     = zu(jo,ja,jk,jt)
281                  inpfile%pob(jk,iobs,2)     = zv(jo,ja,jk,jt)
282                  inpfile%pdep(jk,iobs)      = zdep(jk)
283                  inpfile%ivlqc(jk,iobs,1:2) = INT( MAX( zuv1qc(jo,ja,jk,jt), zuv2qc(jo,ja,jk,jt) ) )
284               END DO
285               inpfile%plam(iobs) = zlon(jo)
286               inpfile%pphi(iobs) = zlat(ja)
287               inpfile%ptim(iobs) = zjuld(jt)
288            END DO
289         END DO
290      END DO
291
292      ! No position, time, depth and variable QC in input files
293      DO jo = 1, iobs
294         inpfile%ipqc(jo) = 1
295         inpfile%itqc(jo) = 1
296         inpfile%ivqc(jo,1:2) = 1
297         DO jk = 1, ilev
298            inpfile%idqc(jk,jo) = 1
299         END DO
300      END DO
301
302      !---------------------------------------------------------------------
303      ! Set the platform information
304      !---------------------------------------------------------------------
305      inpfile%cdtyp(:)=' 820'
306
307      !---------------------------------------------------------------------
308      ! Set QC flags for missing data and rescale to m/s
309      !---------------------------------------------------------------------
310
311      DO jo = 1, iobs
312         DO jk = 1, ilev
313            IF ( ( ABS(inpfile%pob(jk,jo,1)) > 10000.0_wp ) .OR. &
314               & ( ABS(inpfile%pob(jk,jo,2)) > 10000.0_wp ) ) THEN
315               inpfile%ivlqc(jk,jo,:) = 4
316               inpfile%pob(jk,jo,1) = fbrmdi
317               inpfile%pob(jk,jo,2) = fbrmdi
318            ELSE
319               inpfile%pob(jk,jo,1) = 0.01 * inpfile%pob(jk,jo,1)
320               inpfile%pob(jk,jo,2) = 0.01 * inpfile%pob(jk,jo,2)
321            ENDIF
322         END DO
323      END DO
324
325      !---------------------------------------------------------------------
326      ! Set file indexes
327      !---------------------------------------------------------------------
328
329      DO jo = 1, inpfile%nobs
330         inpfile%kindex(jo) = jo
331      END DO
332
333      !---------------------------------------------------------------------
334      ! Initialize flags since they are not in the TAO input files
335      !---------------------------------------------------------------------
336
337      inpfile%ioqcf(:,:)      = 0
338      inpfile%ipqcf(:,:)      = 0
339      inpfile%itqcf(:,:)      = 0
340      inpfile%idqcf(:,:,:)    = 0
341      inpfile%ivqcf(:,:,:)    = 0
342      inpfile%ivlqcf(:,:,:,:) = 0
343
344      !---------------------------------------------------------------------
345      ! Deallocate data
346      !---------------------------------------------------------------------
347      DEALLOCATE( &
348         & zu,     &
349         & zv,     &
350         & zdep,   &
351         & zuv1qc, &
352         & zuv2qc, &
353         & itime,  &
354         & itime2, &
355         & zlat,   &
356         & zlon,   &
357         & zjuld   &
358         & )
359
360   END SUBROUTINE read_taondbc
Note: See TracBrowser for help on using the repository browser.