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_sla_io.h90 in branches/dev_1784_OBS/NEMO/OPA_SRC/OBS – NEMO

Context Navigation

source: branches/dev_1784_OBS/NEMO/OPA_SRC/OBS/obs_sla_io.h90 @ 2001

Last change on this file since 2001 was 2001, checked in by djlea, 14 years ago
Adding observation operator code
File size: 13.4 KB

Line
1	SUBROUTINE read_avisofile( cdfilename, inpfile, kunit, ldwp, ldgrid )
2	!!---------------------------------------------------------------------
3	!!
4	!! ROUTINE read_avisofile
5	!!
6	!! ** Purpose : Read from file the AVISO SLA observations.
7	!!
8	!! ** Method : The data file is a NetCDF file.
9	!!
10	!! ** Action :
11	!!
12	!! References : http://www.aviso.oceanobs.com
13	!!
14	!! History :
15	!! ! 09-01 (K. Mogensen) Original based on old versions
16	!!----------------------------------------------------------------------
17	!! * Arguments
18	CHARACTER(LEN=*) :: cdfilename ! Input filename
19	TYPE(obfbdata) :: inpfile ! Output obfbdata structure
20	INTEGER :: kunit ! Unit for output
21	LOGICAL :: ldwp ! Print info
22	LOGICAL :: ldgrid ! Save grid info in data structure
23	!! * Local declarations
24	CHARACTER(LEN=14),PARAMETER :: cpname = 'read_avisofile'
25	INTEGER :: &
26	& i_file_id, & ! netcdf IDS
27	& i_tracks_id, &
28	& i_cycles_id, &
29	& i_data_id, &
30	& i_var_id
31	INTEGER, PARAMETER :: &
32	& imaxdim = 2 ! Assumed maximum for no. dims. in file
33	INTEGER, DIMENSION(2) :: &
34	& idims ! Dimensions in file
35	INTEGER :: &
36	& iilen, & ! Length of netCDF attributes
37	& itype ! Typeof netCDF attributes
38	REAL(fbdp) :: &
39	& zsca, & ! Scale factor
40	& zfill ! Fill value
41	CHARACTER(len=3) :: &
42	& cmission ! Mission global attribute
43	INTEGER :: &
44	& itracks, & ! Maximum number of passes in file
45	& icycles, & ! Maximum number of cycles for each pass
46	& idata ! Number of data per parameter in current file
47	REAL(fbdp) :: &
48	& zdeltat ! Time gap getween two measurements in seconds
49	INTEGER, DIMENSION(:), POINTER :: &
50	& iptracks, & ! List of passes contained in current file
51	& ipnbpoints, & ! Number of points per pass
52	& ipdataindexes ! Index of data in theoretical profile
53	INTEGER, DIMENSION(:,:), POINTER :: &
54	& ipcycles ! List of cycles per pass
55	REAL(fbdp), DIMENSION(:), POINTER :: &
56	& zphi, & ! Latitudes
57	& zlam ! Longitudes
58	REAL(fbdp), DIMENSION(:,:), POINTER :: &
59	& zbegindates ! Date of point with index 0
60	REAL(fbsp), DIMENSION(:,:), POINTER :: &
61	& zsla ! SLA data
62	REAL(fbdp), DIMENSION(:), POINTER :: &
63	& zjuld ! Julian date
64	CHARACTER(len=14) :: &
65	& cdjuldref ! Julian data reference
66	INTEGER :: &
67	& imission ! Mission number converted from Mission global
68	! netCDF atttribute.
69	CHARACTER(len=255) :: &
70	& ctmp
71	INTEGER :: &
72	& iobs
73	INTEGER :: &
74	& jl, &
75	& jm, &
76	& jj, &
77	& ji, &
78	& jk, &
79	& jobs, &
80	& jcycle
81
82	! Open the file
83
84	CALL chkerr( nf90_open( TRIM( cdfilename ), nf90_nowrite, i_file_id ), &
85	& cpname, __LINE__ )
86
87	! Get the netCDF dimensions
88
89	CALL chkerr( nf90_inq_dimid( i_file_id, 'Tracks', i_tracks_id ), &
90	& cpname, __LINE__ )
91	CALL chkerr( nf90_inquire_dimension( i_file_id, i_tracks_id, &
92	& len = itracks ), &
93	& cpname, __LINE__ )
94
95	CALL chkerr( nf90_inq_dimid( i_file_id, 'Cycles', i_cycles_id ), &
96	& cpname, __LINE__ )
97	CALL chkerr( nf90_inquire_dimension( i_file_id, i_cycles_id, &
98	& len = icycles ), &
99	& cpname, __LINE__ )
100
101	CALL chkerr( nf90_inq_dimid( i_file_id, 'Data', i_data_id ), &
102	& cpname, __LINE__ )
103	CALL chkerr( nf90_inquire_dimension( i_file_id, i_data_id, &
104	& len = idata ), &
105	& cpname, __LINE__ )
106
107	! Allocate memory for input data
108
109	ALLOCATE( &
110	& iptracks ( itracks ), &
111	& ipnbpoints ( itracks ), &
112	& ipcycles ( icycles, itracks ), &
113	& zphi ( idata ), &
114	& zlam ( idata ), &
115	& zbegindates ( icycles, itracks ), &
116	& ipdataindexes( idata ), &
117	& zsla ( icycles, idata ), &
118	& zjuld ( idata*icycles ) &
119	& )
120
121	! Get time gap getween two measurements in seconds
122
123	CALL chkerr( nf90_inq_varid( i_file_id, 'DeltaT', i_var_id ), &
124	& cpname, __LINE__ )
125	CALL chkdim( i_file_id, i_var_id, 0, idims, cpname, __LINE__ )
126	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, zdeltat), &
127	& cpname, __LINE__ )
128	zsca = 1.0
129	IF (nf90_inquire_attribute( i_file_id, i_var_id, 'scale_factor') &
130	& == nf90_noerr) THEN
131	CALL chkerr( nf90_get_att( i_file_id, i_var_id, &
132	& 'scale_factor',zsca), cpname, __LINE__)
133	ENDIF
134	zdeltat = zsca * zdeltat
135
136	! Get List of passes contained in current file
137
138	CALL chkerr( nf90_inq_varid( i_file_id, 'Tracks', i_var_id ), &
139	& cpname, __LINE__ )
140	idims(1) = itracks
141	CALL chkdim( i_file_id, i_var_id, 1, idims, cpname, __LINE__ )
142	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, iptracks), &
143	& cpname, __LINE__ )
144
145	! Get Number of points per pass
146
147	CALL chkerr( nf90_inq_varid( i_file_id, 'NbPoints', i_var_id ), &
148	& cpname, __LINE__ )
149	idims(1) = itracks
150	CALL chkdim( i_file_id, i_var_id, 1, idims, cpname, __LINE__ )
151	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, ipnbpoints),&
152	& cpname, __LINE__ )
153
154	! Get list of cycles per pass
155
156	CALL chkerr( nf90_inq_varid( i_file_id, 'Cycles', i_var_id ), &
157	& cpname, __LINE__ )
158	idims(1) = icycles
159	idims(2) = itracks
160	CALL chkdim( i_file_id, i_var_id, 2, idims, cpname, __LINE__ )
161	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, ipcycles), &
162	& cpname, __LINE__ )
163
164	! Get longitudes
165
166	CALL chkerr( nf90_inq_varid( i_file_id, 'Longitudes', i_var_id ), &
167	& cpname, __LINE__ )
168	idims(1) = idata
169	CALL chkdim( i_file_id, i_var_id, 1, idims, cpname, __LINE__ )
170	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, zlam), &
171	& cpname, __LINE__ )
172	zsca = 1.0
173	IF (nf90_inquire_attribute( i_file_id, i_var_id, 'scale_factor') &
174	& == nf90_noerr) THEN
175	CALL chkerr( nf90_get_att( i_file_id, i_var_id, &
176	& 'scale_factor',zsca), cpname, __LINE__)
177	ENDIF
178	zlam(:) = zsca * zlam(:)
179
180	! Get latitudes
181
182	CALL chkerr( nf90_inq_varid( i_file_id, 'Latitudes', i_var_id ), &
183	& cpname, __LINE__ )
184	idims(1) = idata
185	CALL chkdim( i_file_id, i_var_id, 1, idims, cpname, __LINE__ )
186	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, zphi), &
187	& cpname, __LINE__ )
188	zsca = 1.0
189	IF (nf90_inquire_attribute( i_file_id, i_var_id, 'scale_factor') &
190	& == nf90_noerr) THEN
191	CALL chkerr( nf90_get_att( i_file_id, i_var_id, &
192	& 'scale_factor',zsca), cpname, __LINE__)
193	ENDIF
194	zphi(:) = zsca * zphi(:)
195
196	! Get date of point with index 0
197
198	CALL chkerr( nf90_inq_varid( i_file_id, 'BeginDates', i_var_id ), &
199	& cpname, __LINE__ )
200	idims(1) = icycles
201	idims(2) = itracks
202	CALL chkdim( i_file_id, i_var_id, 2, idims, cpname, __LINE__ )
203	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, zbegindates), &
204	& cpname, __LINE__ )
205	CALL chkerr( nf90_inquire_attribute( i_file_id, i_var_id, &
206	& 'units', len = iilen, &
207	& xtype = itype), cpname, __LINE__ )
208	IF (( itype /= NF90_CHAR ) .OR. ( iilen > 255 )) THEN
209	CALL fatal_error('Error decoding BeginDates unit', __LINE__ )
210	ENDIF
211	CALL chkerr( nf90_get_att( i_file_id, i_var_id, 'units', &
212	& ctmp ), cpname, __LINE__ )
213	jl=1
214	DO jm = 1, LEN(TRIM(ctmp))
215	IF ((ctmp(jm:jm)>='0').AND.(ctmp(jm:jm)<='9')) THEN
216	cdjuldref(jl:jl) = ctmp(jm:jm)
217	jl=jl+1
218	ENDIF
219	IF (jl>14) EXIT
220	ENDDO
221
222	! Get indices of data in theoretical profile
223
224	CALL chkerr( nf90_inq_varid( i_file_id, 'DataIndexes', i_var_id ), &
225	& cpname, __LINE__ )
226	idims(1) = idata
227	CALL chkdim( i_file_id, i_var_id, 1, idims, cpname, __LINE__ )
228	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, ipdataindexes), &
229	& cpname, __LINE__ )
230
231	! Get SLA data
232
233	CALL chkerr( nf90_inq_varid( i_file_id, 'SLA', i_var_id ), &
234	& cpname, __LINE__ )
235	idims(1) = icycles
236	idims(2) = idata
237	CALL chkdim( i_file_id, i_var_id, 2, idims, cpname, __LINE__ )
238	CALL chkerr( nf90_get_var ( i_file_id, i_var_id, zsla), &
239	& cpname, __LINE__ )
240	zsca = 1.0
241	IF (nf90_inquire_attribute( i_file_id, i_var_id, 'scale_factor') &
242	& == nf90_noerr) THEN
243	CALL chkerr( nf90_get_att( i_file_id, i_var_id, &
244	& 'scale_factor',zsca), cpname, __LINE__ )
245	ENDIF
246	zfill = 0.0
247	IF (nf90_inquire_attribute( i_file_id, i_var_id, '_FillValue') &
248	& == nf90_noerr) THEN
249	CALL chkerr( nf90_get_att( i_file_id, i_var_id, &
250	& '_FillValue',zfill), cpname, __LINE__ )
251	ENDIF
252	WHERE(zsla(:,:) /= zfill)
253	zsla(:,:) = zsca * zsla(:,:)
254	ELSEWHERE
255	zsla(:,:) = fbrmdi
256	END WHERE
257
258	! Get the global Mission netCDF attribute
259
260	cmission=' '
261	CALL chkerr( nf90_inquire_attribute( i_file_id, nf90_global, &
262	& 'Mission', len = iilen, &
263	& xtype = itype), cpname, __LINE__ )
264	IF (( itype /= NF90_CHAR ) .OR. ( iilen > 3 )) THEN
265	CALL fatal_error('Error decoding Mission global attribute', __LINE__ )
266	ENDIF
267	CALL chkerr( nf90_get_att( i_file_id, nf90_global, 'Mission', &
268	& cmission ), cpname, __LINE__ )
269
270	! Convert it to an integer
271	imission = 0
272	DO jm = 1, imaxmissions
273	IF (cmission==cmissions(jm)) THEN
274	imission = jm
275	EXIT
276	ENDIF
277	ENDDO
278
279	! Close the file
280
281	CALL chkerr( nf90_close( i_file_id ), cpname, __LINE__ )
282
283	! Compute Julian dates for all observations
284
285	jm = 0
286	jl = 0
287	DO jj = 1, itracks
288	DO ji = 1, ipnbpoints(jj)
289	jl = jl + 1
290	DO jk = 1, icycles
291	jm = jm + 1
292	zjuld(jm) = zbegindates(jk,jj) + &
293	& (ipdataindexes(jl) * zdeltat / 86400._wp )
294	END DO
295	END DO
296	END DO
297
298	! Allocate obfbdata
299
300	iobs = idata * icycles
301	CALL init_obfbdata( inpfile )
302	CALL alloc_obfbdata( inpfile, 1, iobs, 1, 0, 0, ldgrid )
303	inpfile%cname(1) = 'SLA'
304
305	! Fill the obfbdata structure from input data
306
307	inpfile%cdjuldref = cdjuldref
308	iobs = 0
309	DO jobs = 1, idata
310	DO jcycle = 1, icycles
311	iobs = iobs + 1
312	! Characters
313	WRITE(inpfile%cdwmo(iobs),'(A3,A5)') cmissions(imission), ' '
314	WRITE(inpfile%cdtyp(iobs),'(I4)') imission
315	! Real values
316	inpfile%plam(iobs) = zlam(jobs)
317	inpfile%pphi(iobs) = zphi(jobs)
318	inpfile%pob(1,iobs,1) = zsla(jcycle,jobs)
319	inpfile%ptim(iobs) = zjuld(iobs)
320	inpfile%pdep(1,iobs) = 0.0
321	! Integers
322	inpfile%kindex(iobs) = iobs
323	IF ( zsla(jcycle,jobs) == fbrmdi ) THEN
324	inpfile%ioqc(iobs) = 4
325	inpfile%ivqc(iobs,1) = 4
326	inpfile%ivlqc(1,iobs,1) = 4
327	ELSE
328	inpfile%ioqc(iobs) = 1
329	inpfile%ivqc(iobs,1) = 1
330	inpfile%ivlqc(1,iobs,1) = 1
331	ENDIF
332	inpfile%ipqc(iobs) = 0
333	inpfile%ipqcf(:,iobs) = 0
334	inpfile%itqc(iobs) = 0
335	inpfile%itqcf(:,iobs) = 0
336	inpfile%ivqcf(:,iobs,1) = 0
337	inpfile%ioqcf(:,iobs) = 0
338	inpfile%idqc(1,iobs) = 0
339	inpfile%idqcf(1,1,iobs) = 0
340	inpfile%ivlqcf(:,1,iobs,1) = 0
341	END DO
342	END DO
343
344
345	! Deallocate memory for input data
346
347	DEALLOCATE( &
348	& iptracks, &
349	& ipnbpoints, &
350	& ipcycles, &
351	& zphi, &
352	& zlam, &
353	& zbegindates, &
354	& ipdataindexes, &
355	& zsla, &
356	& zjuld &
357	& )
358
359	END SUBROUTINE read_avisofile
360

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format