1 |
MODULE flinget_m |
2 |
|
3 |
! From flincom.f90, version 2.2 2006/03/07 09:21:51 |
4 |
|
5 |
IMPLICIT NONE |
6 |
|
7 |
PRIVATE |
8 |
PUBLIC flinget |
9 |
|
10 |
INTERFACE flinget |
11 |
MODULE PROCEDURE flinget_r3d, flinget_r2d |
12 |
END INTERFACE |
13 |
|
14 |
CONTAINS |
15 |
|
16 |
SUBROUTINE flinget_r2d & |
17 |
(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
18 |
!--------------------------------------------------------------------- |
19 |
IMPLICIT NONE |
20 |
|
21 |
INTEGER :: fid_in |
22 |
CHARACTER(LEN=*) :: varname |
23 |
INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin |
24 |
REAL :: var(:,:) |
25 |
|
26 |
INTEGER :: jl, jj, ji |
27 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
28 |
LOGICAL :: check = .FALSE. |
29 |
!--------------------------------------------------------------------- |
30 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
31 |
IF (check) WRITE(*,*) & |
32 |
"flinget_r2d : allocate buff_tmp for buff_sz = ",SIZE(var) |
33 |
ALLOCATE (buff_tmp(SIZE(var))) |
34 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
35 |
IF (check) WRITE(*,*) & |
36 |
"flinget_r2d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
37 |
DEALLOCATE (buff_tmp) |
38 |
ALLOCATE (buff_tmp(SIZE(var))) |
39 |
ENDIF |
40 |
|
41 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
42 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
43 |
|
44 |
jl=0 |
45 |
DO jj=1,SIZE(var,2) |
46 |
DO ji=1,SIZE(var,1) |
47 |
jl=jl+1 |
48 |
var(ji,jj) = buff_tmp(jl) |
49 |
ENDDO |
50 |
ENDDO |
51 |
!------------------------- |
52 |
END SUBROUTINE flinget_r2d |
53 |
|
54 |
!=== |
55 |
|
56 |
SUBROUTINE flinget_r3d(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
57 |
!--------------------------------------------------------------------- |
58 |
IMPLICIT NONE |
59 |
|
60 |
INTEGER, intent(in):: fid_in |
61 |
CHARACTER(LEN=*), intent(in):: varname |
62 |
INTEGER, intent(in):: iim, jjm, llm, ttm, itau_dep, itau_fin |
63 |
REAL, intent(out):: var(:,:,:) |
64 |
|
65 |
INTEGER :: jl, jk, jj, ji |
66 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
67 |
LOGICAL :: check = .FALSE. |
68 |
!--------------------------------------------------------------------- |
69 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
70 |
IF (check) WRITE(*,*) & |
71 |
"flinget_r3d : allocate buff_tmp for buff_sz = ",SIZE(var) |
72 |
ALLOCATE (buff_tmp(SIZE(var))) |
73 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
74 |
IF (check) WRITE(*,*) & |
75 |
"flinget_r3d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
76 |
DEALLOCATE (buff_tmp) |
77 |
ALLOCATE (buff_tmp(SIZE(var))) |
78 |
ENDIF |
79 |
|
80 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
81 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
82 |
|
83 |
jl=0 |
84 |
DO jk=1,SIZE(var,3) |
85 |
DO jj=1,SIZE(var,2) |
86 |
DO ji=1,SIZE(var,1) |
87 |
jl=jl+1 |
88 |
var(ji,jj,jk) = buff_tmp(jl) |
89 |
ENDDO |
90 |
ENDDO |
91 |
ENDDO |
92 |
!------------------------- |
93 |
END SUBROUTINE flinget_r3d |
94 |
|
95 |
!=== |
96 |
|
97 |
SUBROUTINE flinget_mat & |
98 |
(fid_in, varname, iim, jjm, llm, ttm, itau_dep, & |
99 |
itau_fin, iideb, iilen, jjdeb, jjlen, var) |
100 |
!--------------------------------------------------------------------- |
101 |
!- This subroutine will read the variable named varname from |
102 |
!- the file previously opened by flinopen and identified by fid |
103 |
|
104 |
!- It is checked that the dimensions of the variable to be read |
105 |
!- correspond to what the user requested when he specified |
106 |
!- iim, jjm and llm. The only exception which is allowed is |
107 |
!- for compressed data where the horizontal grid is not expected |
108 |
!- to be iim x jjm. |
109 |
|
110 |
!- If variable is of size zero a global attribute is read. |
111 |
!- This global attribute will be of type real |
112 |
|
113 |
!- INPUT |
114 |
|
115 |
!- fid : File ID returned by flinopen |
116 |
!- varname : Name of the variable to be read from the file |
117 |
!- iim : | These three variables give the size of the variables |
118 |
!- jjm : | to be read. It will be verified that the variables |
119 |
!- llm : | fits in there. |
120 |
!- ttm : | |
121 |
!- itau_dep : Time step at which we will start to read |
122 |
!- itau_fin : Time step until which we are going to read |
123 |
!- For the moment this is done on indexes |
124 |
!- but it should be in the physical space. |
125 |
!- If there is no time-axis in the file then use a |
126 |
!- itau_fin < itau_dep, this will tell flinget not to |
127 |
!- expect a time-axis in the file. |
128 |
!- iideb : index i for zoom |
129 |
!- iilen : length of zoom |
130 |
!- jjdeb : index j for zoom |
131 |
!- jjlen : length of zoom |
132 |
|
133 |
!- OUTPUT |
134 |
|
135 |
!- var : array that will contain the data |
136 |
!--------------------------------------------------------------------- |
137 |
USE strlowercase_m, ONLY : strlowercase |
138 |
USE errioipsl, ONLY : histerr |
139 |
USE netcdf, ONLY : nf90_byte, nf90_double, nf90_float, nf90_get_att, & |
140 |
nf90_get_var, nf90_inquire_attribute, nf90_inquire_dimension, & |
141 |
nf90_inquire_variable, nf90_inq_attname, nf90_inq_varid, nf90_int, & |
142 |
nf90_max_var_dims, nf90_noerr, nf90_short, nf90_strerror |
143 |
use flincom, only: ncids |
144 |
|
145 |
IMPLICIT NONE |
146 |
|
147 |
! ARGUMENTS |
148 |
|
149 |
INTEGER, intent(in):: fid_in |
150 |
CHARACTER(LEN=*), intent(in):: varname |
151 |
INTEGER, intent(in):: iim, jjm, llm, ttm, itau_dep, itau_fin |
152 |
INTEGER :: iideb |
153 |
integer, intent(in):: iilen |
154 |
integer jjdeb |
155 |
integer, intent(in):: jjlen |
156 |
REAL :: var(:) |
157 |
|
158 |
! LOCAL |
159 |
|
160 |
INTEGER, SAVE :: cind_vid |
161 |
INTEGER, SAVE :: cind_fid |
162 |
INTEGER, SAVE :: cind_len |
163 |
INTEGER,DIMENSION(:),ALLOCATABLE,SAVE :: cindex |
164 |
INTEGER,DIMENSION(4) :: w_sta, w_len, w_dim |
165 |
INTEGER :: iret, fid |
166 |
INTEGER :: vid, cvid, clen |
167 |
CHARACTER(LEN=70) :: str1 |
168 |
CHARACTER(LEN=250) :: att_n, tmp_n |
169 |
INTEGER :: tmp_i |
170 |
REAL,SAVE :: mis_v=0. |
171 |
REAL :: tmp_r |
172 |
INTEGER :: ndims, x_typ, nb_atts |
173 |
INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: dimids |
174 |
INTEGER :: i, i2d, cnd |
175 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: var_tmp |
176 |
LOGICAL :: uncompress = .FALSE. |
177 |
LOGICAL :: check = .FALSE. |
178 |
!--------------------------------------------------------------------- |
179 |
fid = ncids(fid_in) |
180 |
|
181 |
IF (check) THEN |
182 |
WRITE(*,*) & |
183 |
'flinget_mat : fid_in, fid, varname :', fid_in, fid, TRIM(varname) |
184 |
WRITE(*,*) & |
185 |
'flinget_mat : iim, jjm, llm, ttm, itau_dep, itau_fin :', & |
186 |
iim, jjm, llm, ttm, itau_dep, itau_fin |
187 |
WRITE(*,*) & |
188 |
'flinget_mat : iideb, iilen, jjdeb, jjlen :', & |
189 |
iideb, iilen, jjdeb, jjlen |
190 |
ENDIF |
191 |
|
192 |
uncompress = .FALSE. |
193 |
|
194 |
! 1.0 We get first all the details on this variable from the file |
195 |
|
196 |
vid = -1 |
197 |
iret = NF90_INQ_VARID (fid, varname, vid) |
198 |
|
199 |
IF (vid < 0 .OR. iret /= NF90_NOERR) THEN |
200 |
CALL histerr (3,'flinget', & |
201 |
'Variable '//TRIM(varname)//' not found in file',' ',' ') |
202 |
ENDIF |
203 |
|
204 |
iret = NF90_INQUIRE_VARIABLE (fid, vid, & |
205 |
ndims=ndims, dimids=dimids, nAtts=nb_atts) |
206 |
IF (check) THEN |
207 |
WRITE(*,*) & |
208 |
'flinget_mat : fid, vid :', fid, vid |
209 |
WRITE(*,*) & |
210 |
'flinget_mat : ndims, dimids(1:ndims), nb_atts :', & |
211 |
ndims, dimids(1:ndims), nb_atts |
212 |
ENDIF |
213 |
|
214 |
w_dim(:) = 0 |
215 |
DO i=1,ndims |
216 |
iret = NF90_INQUIRE_DIMENSION (fid, dimids(i), len=w_dim(i)) |
217 |
ENDDO |
218 |
IF (check) WRITE(*,*) & |
219 |
'flinget_mat : w_dim :', w_dim(1:ndims) |
220 |
|
221 |
mis_v = 0.0 |
222 |
|
223 |
IF (nb_atts > 0) THEN |
224 |
IF (check) THEN |
225 |
WRITE(*,*) 'flinget_mat : attributes for variable :' |
226 |
ENDIF |
227 |
ENDIF |
228 |
DO i=1,nb_atts |
229 |
iret = NF90_INQ_ATTNAME (fid, vid, i, att_n) |
230 |
iret = NF90_INQUIRE_ATTRIBUTE (fid, vid, att_n, xtype=x_typ) |
231 |
CALL strlowercase (att_n) |
232 |
IF ( (x_typ == NF90_INT).OR.(x_typ == NF90_SHORT) & |
233 |
.OR.(x_typ == NF90_BYTE) ) THEN |
234 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_i) |
235 |
IF (check) THEN |
236 |
WRITE(*,*) ' ',TRIM(att_n),' : ',tmp_i |
237 |
ENDIF |
238 |
ELSE IF ( (x_typ == NF90_FLOAT).OR.(x_typ == NF90_DOUBLE) ) THEN |
239 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_r) |
240 |
IF (check) THEN |
241 |
WRITE(*,*) ' ',TRIM(att_n),' : ',tmp_r |
242 |
ENDIF |
243 |
IF (index(att_n,'missing_value') > 0) THEN |
244 |
mis_v = tmp_r |
245 |
ENDIF |
246 |
ELSE |
247 |
tmp_n = '' |
248 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_n) |
249 |
IF (check) THEN |
250 |
WRITE(*,*) ' ',TRIM(att_n),' : ',TRIM(tmp_n) |
251 |
ENDIF |
252 |
ENDIF |
253 |
ENDDO |
254 |
!? |
255 |
!!!!!!!!!! We will need a verification on the type of the variable |
256 |
!? |
257 |
|
258 |
! 2.0 The dimensions are analysed to determine what is to be read |
259 |
|
260 |
! 2.1 the longitudes |
261 |
|
262 |
IF ( w_dim(1) /= iim .OR. w_dim(2) /= jjm) THEN |
263 |
!--- |
264 |
!-- There is a possibility that we have to deal with a compressed axis ! |
265 |
!--- |
266 |
iret = NF90_INQUIRE_DIMENSION (fid, dimids(1), & |
267 |
name=tmp_n, len=clen) |
268 |
iret = NF90_INQ_VARID (fid, tmp_n, cvid) |
269 |
!--- |
270 |
IF (check) WRITE(*,*) & |
271 |
'Dimname, iret , NF90_NOERR : ',TRIM(tmp_n),iret,NF90_NOERR |
272 |
!--- |
273 |
!-- If we have an axis which has the same name |
274 |
!-- as the dimension we can see if it is compressed |
275 |
!--- |
276 |
!-- TODO TODO for zoom2d |
277 |
!--- |
278 |
IF (iret == NF90_NOERR) THEN |
279 |
iret = NF90_GET_ATT (fid, cvid, 'compress', str1) |
280 |
!----- |
281 |
IF (iret == NF90_NOERR) THEN |
282 |
iret = NF90_INQUIRE_VARIABLE (fid,cvid,xtype=x_typ,ndims=cnd) |
283 |
!------- |
284 |
IF ( cnd /= 1 .AND. x_typ /= NF90_INT) THEN |
285 |
CALL histerr (3,'flinget', & |
286 |
'Variable '//TRIM(tmp_n)//' can not be a compressed axis', & |
287 |
'Either it has too many dimensions'// & |
288 |
' or it is not of type integer', ' ') |
289 |
ELSE |
290 |
!--------- |
291 |
!-------- Let us see if we already have that index table |
292 |
!--------- |
293 |
IF ( (cind_len /= clen).OR.(cind_vid /= cvid) & |
294 |
.OR.(cind_fid /= fid) ) THEN |
295 |
IF (ALLOCATED(cindex)) DEALLOCATE(cindex) |
296 |
ALLOCATE(cindex(clen)) |
297 |
cind_len = clen |
298 |
cind_vid = cvid |
299 |
cind_fid = fid |
300 |
iret = NF90_GET_VAR (fid, cvid, cindex) |
301 |
ENDIF |
302 |
!--------- |
303 |
!-------- In any case we need to set the slab of data to be read |
304 |
!--------- |
305 |
uncompress = .TRUE. |
306 |
w_sta(1) = 1 |
307 |
w_len(1) = clen |
308 |
i2d = 1 |
309 |
ENDIF |
310 |
ELSE |
311 |
str1 = 'The horizontal dimensions of '//varname |
312 |
CALL histerr (3,'flinget',str1, & |
313 |
'is not compressed and does not'// & |
314 |
' correspond to the requested size',' ') |
315 |
ENDIF |
316 |
ELSE |
317 |
IF (w_dim(1) /= iim) THEN |
318 |
str1 = 'The longitude dimension of '//varname |
319 |
CALL histerr (3,'flinget',str1, & |
320 |
'in the file is not equal to the dimension', & |
321 |
'that should be read') |
322 |
ENDIF |
323 |
IF (w_dim(2) /= jjm) THEN |
324 |
str1 = 'The latitude dimension of '//varname |
325 |
CALL histerr (3,'flinget',str1, & |
326 |
'in the file is not equal to the dimension', & |
327 |
'that should be read') |
328 |
ENDIF |
329 |
ENDIF |
330 |
ELSE |
331 |
w_sta(1:2) = (/ iideb, jjdeb /) |
332 |
w_len(1:2) = (/ iilen, jjlen /) |
333 |
i2d = 2 |
334 |
ENDIF |
335 |
|
336 |
! 2.3 Now the difficult part, the 3rd dimension which can be |
337 |
! time or levels. |
338 |
|
339 |
! Priority is given to the time axis if only three axes are present. |
340 |
|
341 |
IF (ndims > i2d) THEN |
342 |
!--- |
343 |
!-- 2.3.1 We have a vertical axis |
344 |
!--- |
345 |
IF (llm == 1 .AND. ndims == i2d+2 .OR. llm == w_dim(i2d+1)) THEN |
346 |
!----- |
347 |
IF (w_dim(i2d+1) /= llm) THEN |
348 |
CALL histerr (3,'flinget', & |
349 |
'The vertical dimension of '//varname, & |
350 |
'in the file is not equal to the dimension', & |
351 |
'that should be read') |
352 |
ELSE |
353 |
w_sta(i2d+1) = 1 |
354 |
IF (llm > 0) THEN |
355 |
w_len(i2d+1) = llm |
356 |
ELSE |
357 |
w_len(i2d+1) = w_sta(i2d+1) |
358 |
ENDIF |
359 |
ENDIF |
360 |
!----- |
361 |
IF ((itau_fin-itau_dep) >= 0) THEN |
362 |
IF (ndims /= i2d+2) THEN |
363 |
CALL histerr (3,'flinget', & |
364 |
'You attempt to read a time slab', & |
365 |
'but there is no time axis on this variable', varname) |
366 |
ELSE IF ((itau_fin - itau_dep) <= w_dim(i2d+2)) THEN |
367 |
w_sta(i2d+2) = itau_dep |
368 |
w_len(i2d+2) = itau_fin-itau_dep+1 |
369 |
ELSE |
370 |
CALL histerr (3,'flinget', & |
371 |
'The time step you try to read is not', & |
372 |
'in the file (1)', varname) |
373 |
ENDIF |
374 |
ELSE IF (ndims == i2d+2 .AND. w_dim(i2d+2) > 1) THEN |
375 |
CALL histerr (3,'flinget', & |
376 |
'There is a time axis in the file but no', & |
377 |
'time step give in the call', varname) |
378 |
ELSE |
379 |
w_sta(i2d+2) = 1 |
380 |
w_len(i2d+2) = 1 |
381 |
ENDIF |
382 |
ELSE |
383 |
!----- |
384 |
!---- 2.3.2 We do not have any vertical axis |
385 |
!----- |
386 |
IF (ndims == i2d+2) THEN |
387 |
CALL histerr (3,'flinget', & |
388 |
'The file contains 4 dimensions', & |
389 |
'but only 3 are requestes for variable ', varname) |
390 |
ENDIF |
391 |
IF ((itau_fin-itau_dep) >= 0) THEN |
392 |
IF (ndims == i2d+1) THEN |
393 |
IF ((itau_fin-itau_dep) < w_dim(i2d+1) ) THEN |
394 |
w_sta(i2d+1) = itau_dep |
395 |
w_len(i2d+1) = itau_fin-itau_dep+1 |
396 |
ELSE |
397 |
CALL histerr (3,'flinget', & |
398 |
'The time step you try to read is not', & |
399 |
'in the file (2)', varname) |
400 |
ENDIF |
401 |
ELSE |
402 |
CALL histerr (3,'flinget', & |
403 |
'From your input you sould have 3 dimensions', & |
404 |
'in the file but there are 4', varname) |
405 |
ENDIF |
406 |
ELSE |
407 |
IF (ndims == i2d+1 .AND. w_dim(i2d+1) > 1) THEN |
408 |
CALL histerr (3,'flinget', & |
409 |
'There is a time axis in the file but no', & |
410 |
'time step given in the call', varname) |
411 |
ELSE |
412 |
w_sta(i2d+1) = 1 |
413 |
w_len(i2d+1) = 1 |
414 |
ENDIF |
415 |
ENDIF |
416 |
ENDIF |
417 |
ELSE |
418 |
!--- |
419 |
!-- 2.3.3 We do not have any vertical axis |
420 |
!--- |
421 |
w_sta(i2d+1:i2d+2) = (/ 0, 0 /) |
422 |
w_len(i2d+1:i2d+2) = (/ 0, 0 /) |
423 |
ENDIF |
424 |
|
425 |
! 3.0 Reading the data |
426 |
|
427 |
IF (check) WRITE(*,*) & |
428 |
'flinget_mat 3.0 : ', uncompress, w_sta, w_len |
429 |
!--- |
430 |
IF (uncompress) THEN |
431 |
!--- |
432 |
IF (ALLOCATED(var_tmp)) THEN |
433 |
IF (SIZE(var_tmp) < clen) THEN |
434 |
DEALLOCATE(var_tmp) |
435 |
ALLOCATE(var_tmp(clen)) |
436 |
ENDIF |
437 |
ELSE |
438 |
ALLOCATE(var_tmp(clen)) |
439 |
ENDIF |
440 |
!--- |
441 |
iret = NF90_GET_VAR (fid, vid, var_tmp, & |
442 |
start=w_sta(:), count=w_len(:)) |
443 |
!--- |
444 |
var(:) = mis_v |
445 |
var(cindex(:)) = var_tmp(:) |
446 |
!--- |
447 |
ELSE |
448 |
iret = NF90_GET_VAR (fid, vid, var, & |
449 |
start=w_sta(:), count=w_len(:)) |
450 |
ENDIF |
451 |
|
452 |
IF (check) WRITE(*,*) 'flinget_mat 3.1 : ',NF90_STRERROR (iret) |
453 |
!-------------------------- |
454 |
END SUBROUTINE flinget_mat |
455 |
|
456 |
END MODULE flinget_m |