1 |
MODULE flincom |
2 |
|
3 |
! From flincom.f90, version 2.2 2006/03/07 09:21:51 |
4 |
|
5 |
USE netcdf |
6 |
|
7 |
USE calendar, ONLY : ju2ymds, ymds2ju, ioconf_calendar |
8 |
USE errioipsl, ONLY : histerr |
9 |
USE stringop, ONLY : strlowercase |
10 |
|
11 |
IMPLICIT NONE |
12 |
|
13 |
PRIVATE |
14 |
PUBLIC flinput, flincre, flinget, flinget_zoom2d, flinclo |
15 |
public flinopen_nozoom |
16 |
public flininfo, flininspect, flinquery_var |
17 |
|
18 |
INTERFACE flinput |
19 |
!--------------------------------------------------------------------- |
20 |
!- The "flinput" routines will put a variable |
21 |
!- on the netCDF file created by flincre. |
22 |
!- If the sizes of the axis do not match the one of the IDs |
23 |
!- then a new axis is created. |
24 |
!- That is we loose the possibility of writting hyperslabs of data. |
25 |
|
26 |
!- Again here if iim = jjm = llm = ttm = 0 |
27 |
!- then a global attribute is added to the file. |
28 |
|
29 |
!- INPUT |
30 |
|
31 |
!- fid : Identification of the file in which we will write |
32 |
!- varname : Name of variable to be written |
33 |
!- iim : size in x of variable |
34 |
!- nlonid : ID of x axis which could fit for this axis |
35 |
!- jjm : size in y of variable |
36 |
!- nlatid : ID of y axis which could fit for this axis |
37 |
!- llm : size in z of variable |
38 |
!- zdimid : ID of z axis which could fit for this axis |
39 |
!- ttm : size in t of variable |
40 |
!- tdimid : ID of t axis which could fit for this axis |
41 |
|
42 |
!- OUTPUT |
43 |
|
44 |
!- NONE |
45 |
!--------------------------------------------------------------------- |
46 |
MODULE PROCEDURE flinput_r4d, flinput_r3d, flinput_r2d, & |
47 |
flinput_r1d, flinput_scal |
48 |
END INTERFACE |
49 |
|
50 |
INTERFACE flinget |
51 |
MODULE PROCEDURE flinget_r4d, flinget_r3d, flinget_r2d, flinget_r1d, & |
52 |
flinget_scal |
53 |
END INTERFACE |
54 |
INTERFACE flinget_zoom2d |
55 |
MODULE PROCEDURE flinget_r4d_zoom2d, flinget_r3d_zoom2d, & |
56 |
flinget_r2d_zoom2d |
57 |
END INTERFACE |
58 |
|
59 |
! This is the data we keep on each file we open |
60 |
|
61 |
INTEGER, PARAMETER :: nbfile_max = 200 |
62 |
INTEGER, SAVE :: nbfiles = 0 |
63 |
INTEGER, SAVE :: ncids(nbfile_max), ncnbd(nbfile_max), & |
64 |
ncfunli(nbfile_max), ncnba(nbfile_max) |
65 |
INTEGER, SAVE :: ncnbva(nbfile_max), ncdims(nbfile_max,4) |
66 |
LOGICAL, SAVE :: ncfileopen(nbfile_max)=.FALSE. |
67 |
|
68 |
INTEGER, SAVE :: cind_vid, cind_fid, cind_len |
69 |
INTEGER,DIMENSION(:),ALLOCATABLE,SAVE :: cindex |
70 |
|
71 |
INTEGER,DIMENSION(4) :: w_sta, w_len, w_dim |
72 |
|
73 |
CONTAINS |
74 |
|
75 |
SUBROUTINE flincre & |
76 |
(filename, iim1, jjm1, lon1, lat1, llm1, lev1, ttm1, itaus, & |
77 |
time0, dt, fid_out, nlonid1, nlatid1, zdimid1, tdimid1) |
78 |
!--------------------------------------------------------------------- |
79 |
!- This is a "low level" subroutine for opening netCDF files wich |
80 |
!- contain the major coordinate system of the model. |
81 |
!- Other coordinates needed for other variables |
82 |
!- will be added as they are needed. |
83 |
|
84 |
!- INPUT |
85 |
|
86 |
!- filename : Name of the file to be created |
87 |
!- iim1, jjm1 : Horizontal size of the grid |
88 |
!- which will be stored in the file |
89 |
!- lon1, lat1 : Horizontal grids |
90 |
!- llm1 : Size of the vertical grid |
91 |
!- lev1 : Vertical grid |
92 |
!- ttm1 : Size of time axis |
93 |
!- itaus : time steps on the time axis |
94 |
!- time0 : Time in julian days at which itau = 0 |
95 |
!- dt : time step in seconds between itaus |
96 |
!- (one step of itau) |
97 |
|
98 |
!- OUTPUT |
99 |
|
100 |
!- fid : File identification |
101 |
!- nlonid1 : Identification of longitudinal axis |
102 |
!- nlatid1 : Identification of latitudinal axis |
103 |
!- zdimid1 : ID of vertical axis |
104 |
!- tdimid1 : ID of time axis |
105 |
!--------------------------------------------------------------------- |
106 |
IMPLICIT NONE |
107 |
|
108 |
! ARGUMENTS |
109 |
|
110 |
CHARACTER(LEN=*) :: filename |
111 |
INTEGER :: iim1, jjm1, llm1, ttm1 |
112 |
REAL :: lon1(iim1,jjm1) |
113 |
REAL :: lat1(iim1,jjm1) |
114 |
REAL :: lev1(llm1) |
115 |
INTEGER :: itaus(ttm1) |
116 |
REAL :: time0 |
117 |
REAL :: dt |
118 |
INTEGER :: fid_out, zdimid1, nlonid1, nlatid1, tdimid1 |
119 |
|
120 |
! LOCAL |
121 |
|
122 |
INTEGER :: iret, lll, fid |
123 |
INTEGER :: lonid, latid, levid, timeid |
124 |
INTEGER :: year, month, day |
125 |
REAL :: sec |
126 |
CHARACTER(LEN=250):: name |
127 |
|
128 |
LOGICAL :: check = .FALSE. |
129 |
!--------------------------------------------------------------------- |
130 |
lll = LEN_TRIM(filename) |
131 |
IF (filename(lll-2:lll) /= '.nc') THEN |
132 |
name=filename(1:lll)//'.nc' |
133 |
ELSE |
134 |
name=filename(1:lll) |
135 |
ENDIF |
136 |
|
137 |
iret = NF90_CREATE (name, NF90_CLOBBER, fid) |
138 |
|
139 |
iret = NF90_DEF_DIM (fid, 'x', iim1, nlonid1) |
140 |
iret = NF90_DEF_DIM (fid, 'y', jjm1, nlatid1) |
141 |
iret = NF90_DEF_DIM (fid, 'lev', llm1, zdimid1) |
142 |
iret = NF90_DEF_DIM (fid, 'tstep', ttm1, tdimid1) |
143 |
|
144 |
! Vertical axis |
145 |
|
146 |
IF (check) WRITE(*,*) 'flincre Vertical axis' |
147 |
|
148 |
iret = NF90_DEF_VAR (fid, 'lev', NF90_FLOAT, zdimid1, levid) |
149 |
iret = NF90_PUT_ATT (fid, levid, 'units', '-') |
150 |
iret = NF90_PUT_ATT (fid, levid, 'title', 'levels') |
151 |
iret = NF90_PUT_ATT (fid, levid, 'long_name', 'Sigma Levels') |
152 |
|
153 |
! Time axis |
154 |
|
155 |
IF (check) WRITE(*,*) 'flincre time axis' |
156 |
|
157 |
iret = NF90_DEF_VAR (fid, 'tstep', NF90_FLOAT, tdimid1, timeid) |
158 |
iret = NF90_PUT_ATT (fid, timeid, 'units', '-') |
159 |
iret = NF90_PUT_ATT (fid, timeid, 'title', 'time') |
160 |
iret = NF90_PUT_ATT (fid, timeid, 'long_name', 'time steps') |
161 |
|
162 |
! The longitude |
163 |
|
164 |
IF (check) WRITE(*,*) 'flincre Longitude axis' |
165 |
|
166 |
iret = NF90_DEF_VAR (fid, "nav_lon", NF90_FLOAT, & |
167 |
(/ nlonid1, nlatid1 /), lonid) |
168 |
iret = NF90_PUT_ATT (fid, lonid, 'units', "degrees_east") |
169 |
iret = NF90_PUT_ATT (fid, lonid, 'title', "Longitude") |
170 |
iret = NF90_PUT_ATT (fid, lonid, 'nav_model', & |
171 |
"Lambert projection of PROMES") |
172 |
iret = NF90_PUT_ATT (fid, lonid, 'valid_min', & |
173 |
REAL(MINVAL(lon1))) |
174 |
iret = NF90_PUT_ATT (fid, lonid, 'valid_max', & |
175 |
REAL(MAXVAL(lon1))) |
176 |
|
177 |
! The Latitude |
178 |
|
179 |
IF (check) WRITE(*,*) 'flincre Latitude axis' |
180 |
|
181 |
iret = NF90_DEF_VAR (fid, "nav_lat", NF90_FLOAT, & |
182 |
(/ nlonid1, nlatid1 /), latid) |
183 |
iret = NF90_PUT_ATT (fid, latid, 'units', "degrees_north") |
184 |
iret = NF90_PUT_ATT (fid, latid, 'title', "Latitude") |
185 |
iret = NF90_PUT_ATT (fid, latid, 'nav_model', & |
186 |
"Lambert projection of PROMES") |
187 |
iret = NF90_PUT_ATT (fid, latid, 'valid_min', & |
188 |
REAL(MINVAL(lat1))) |
189 |
iret = NF90_PUT_ATT (fid, latid, 'valid_max', & |
190 |
REAL(MAXVAL(lat1))) |
191 |
|
192 |
! The time coordinates |
193 |
|
194 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', & |
195 |
REAL(dt)) |
196 |
|
197 |
CALL ju2ymds (time0, year, month, day, sec) |
198 |
|
199 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'year0', REAL(year)) |
200 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'month0', REAL(month)) |
201 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'day0', REAL(day)) |
202 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'sec0', REAL(sec)) |
203 |
|
204 |
iret = NF90_ENDDEF (fid) |
205 |
|
206 |
IF (check) WRITE(*,*) 'flincre Variable' |
207 |
|
208 |
iret = NF90_PUT_VAR (fid, levid, lev1(1:llm1)) |
209 |
|
210 |
IF (check) WRITE(*,*) 'flincre Time Variable' |
211 |
|
212 |
iret = NF90_PUT_VAR (fid, timeid, REAL(itaus(1:ttm1))) |
213 |
|
214 |
IF (check) WRITE(*,*) 'flincre Longitude' |
215 |
|
216 |
iret = NF90_PUT_VAR (fid, lonid, lon1(1:iim1,1:jjm1)) |
217 |
|
218 |
IF (check) WRITE(*,*) 'flincre Latitude' |
219 |
|
220 |
iret = NF90_PUT_VAR (fid, latid, lat1(1:iim1,1:jjm1)) |
221 |
|
222 |
! Keep all this information |
223 |
|
224 |
nbfiles = nbfiles+1 |
225 |
|
226 |
IF (nbfiles > nbfile_max) THEN |
227 |
CALL histerr (3,'flincre', & |
228 |
'Too many files. Please increase nbfil_max', & |
229 |
'in program flincom.F90.',' ') |
230 |
ENDIF |
231 |
|
232 |
ncids(nbfiles) = fid |
233 |
ncnbd(nbfiles) = 4 |
234 |
|
235 |
ncdims(nbfiles,1:4) = (/ iim1, jjm1, llm1, ttm1 /) |
236 |
|
237 |
ncfunli(nbfiles) = -1 |
238 |
ncnba(nbfiles) = 4 |
239 |
ncnbva(nbfiles) = 0 |
240 |
ncfileopen(nbfiles) = .TRUE. |
241 |
|
242 |
fid_out = nbfiles |
243 |
!--------------------- |
244 |
END SUBROUTINE flincre |
245 |
|
246 |
!=== |
247 |
|
248 |
SUBROUTINE flinopen_nozoom(filename, iim, jjm, llm, lon, lat, lev, & |
249 |
ttm, itaus, date0, dt, fid_out) |
250 |
|
251 |
!- The routine will open an input file |
252 |
!- INPUT |
253 |
!- filename : Name of the netCDF file to be opened |
254 |
|
255 |
!- There is no test of the content of the file against the input |
256 |
! from the model |
257 |
|
258 |
!- iim : size in the x direction in the file (longitude) |
259 |
!- jjm : size in the y direction |
260 |
!- llm : number of levels |
261 |
!- (llm = 0 means no axis to be expected) |
262 |
|
263 |
!- WARNING : |
264 |
!- It is for the user to check |
265 |
!- that the dimensions of lon lat and lev are correct when passed to |
266 |
!- flinopen. This can be done after the call when iim and jjm have |
267 |
!- been retrieved from the netCDF file. In F90 this problem will |
268 |
!- be solved with an internal assign |
269 |
!- IF iim, jjm, llm or ttm are parameters in the calling program |
270 |
!- it will create a segmentation fault |
271 |
|
272 |
!- ttm : size of time axis |
273 |
|
274 |
!- OUTPUT |
275 |
|
276 |
!- lon : array of (iim,jjm), |
277 |
!- that contains the longitude of each point |
278 |
!- lat : same for latitude |
279 |
!- lev : An array of llm for the latitude |
280 |
!- itaus : Time steps within this file |
281 |
!- date0 : Julian date at which itau = 0 |
282 |
!- dt : length of the time steps of the data |
283 |
|
284 |
!--------------------------------------------------------------------- |
285 |
|
286 |
IMPLICIT NONE |
287 |
|
288 |
! ARGUMENTS |
289 |
|
290 |
CHARACTER(LEN=*), intent(in):: filename |
291 |
INTEGER, intent(in) :: iim, jjm, llm, ttm |
292 |
real, intent(out):: lon(iim,jjm), lat(iim,jjm), lev(llm) |
293 |
INTEGER, intent(out):: itaus(ttm) |
294 |
REAL, intent(out):: date0, dt |
295 |
|
296 |
INTEGER, intent(in):: fid_out |
297 |
! (file ID which is later used to read the data) |
298 |
|
299 |
! LOCAL |
300 |
|
301 |
INTEGER :: iret, vid, fid, nbdim, i |
302 |
INTEGER :: gdtt_id, old_id, iv, gdtmaf_id |
303 |
CHARACTER(LEN=250) :: name |
304 |
CHARACTER(LEN=80) :: units, calendar |
305 |
INTEGER :: year, month, day |
306 |
REAL :: r_year, r_month, r_day |
307 |
INTEGER :: year0, month0, day0, hours0, minutes0, seci |
308 |
REAL :: sec, sec0 |
309 |
CHARACTER :: strc |
310 |
|
311 |
REAL,DIMENSION(:),ALLOCATABLE :: vec_tmp |
312 |
|
313 |
!--------------------------------------------------------------------- |
314 |
|
315 |
IF ( (fid_out < 1).OR.(fid_out > nbfile_max) ) THEN |
316 |
! Either the fid_out has not been initialized (0 or very large) |
317 |
! then we have to open anyway. Else we only need to open the file |
318 |
! if it has not been opened before. |
319 |
print *, "Call flinfo before flinopen" |
320 |
stop 1 |
321 |
end IF |
322 |
IF (.NOT.ncfileopen(fid_out)) THEN |
323 |
print *, "Call flinfo before flinopen" |
324 |
stop 1 |
325 |
end IF |
326 |
|
327 |
!-- The user has already opened the file |
328 |
!-- and we trust that he knows the dimensions |
329 |
|
330 |
fid = ncids(fid_out) |
331 |
|
332 |
! 2.0 get the sizes and names of the different coordinates |
333 |
! and do a first set of verification. |
334 |
|
335 |
! 3.0 Check if we are realy talking about the same coodinate system |
336 |
! if not then we get the lon, lat and lev variables from the file |
337 |
|
338 |
!-- 4.0 extracting the coordinates |
339 |
!--- |
340 |
CALL flinfindcood (fid_out, 'lon', vid, nbdim) |
341 |
IF (nbdim == 2) THEN |
342 |
iret = NF90_GET_VAR (fid, vid, lon, & |
343 |
start=(/ 1, 1 /), count=(/ iim, jjm /)) |
344 |
ELSE |
345 |
ALLOCATE(vec_tmp(iim)) |
346 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
347 |
start=(/ 1 /), count=(/ iim /)) |
348 |
DO i=1,jjm |
349 |
lon(:,i) = vec_tmp(:) |
350 |
ENDDO |
351 |
DEALLOCATE(vec_tmp) |
352 |
ENDIF |
353 |
!--- |
354 |
CALL flinfindcood (fid_out, 'lat', vid, nbdim) |
355 |
IF (nbdim == 2) THEN |
356 |
iret = NF90_GET_VAR (fid, vid, lat, & |
357 |
start=(/ 1, 1 /), count=(/ iim, jjm /)) |
358 |
ELSE |
359 |
ALLOCATE(vec_tmp(jjm)) |
360 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
361 |
start=(/ 1 /), count=(/ jjm /)) |
362 |
DO i=1,iim |
363 |
lat(i,:) = vec_tmp(:) |
364 |
ENDDO |
365 |
DEALLOCATE(vec_tmp) |
366 |
ENDIF |
367 |
!--- |
368 |
IF (llm > 0) THEN |
369 |
CALL flinfindcood (fid_out, 'lev', vid, nbdim) |
370 |
IF (nbdim == 1) THEN |
371 |
iret = NF90_GET_VAR (fid, vid, lev, & |
372 |
start=(/ 1 /), count=(/ llm /)) |
373 |
ELSE |
374 |
CALL histerr (3,'flinopen', & |
375 |
'Can not handle vertical coordinates that have more',& |
376 |
'than 1 dimension',' ') |
377 |
ENDIF |
378 |
ENDIF |
379 |
|
380 |
! 5.0 Get all the details for the time if possible needed |
381 |
|
382 |
IF (ttm > 0) THEN |
383 |
!--- |
384 |
!-- 5.1 Find the time axis. Prefered method is the 'timestep since' |
385 |
!--- |
386 |
gdtmaf_id = -1 |
387 |
gdtt_id = -1 |
388 |
old_id = -1 |
389 |
DO iv=1,ncnbva(fid_out) |
390 |
name='' |
391 |
iret = NF90_INQUIRE_VARIABLE (fid, iv, name=name) |
392 |
units='' |
393 |
iret = NF90_GET_ATT (fid, iv, 'units', units) |
394 |
IF (INDEX(units,'seconds since') > 0) gdtmaf_id = iv |
395 |
IF (INDEX(units,'timesteps since') > 0) gdtt_id = iv |
396 |
IF (INDEX(name, 'tstep') > 0) old_id = iv |
397 |
ENDDO |
398 |
!--- |
399 |
IF (gdtt_id > 0) THEN |
400 |
vid = gdtt_id |
401 |
ELSE IF (gdtmaf_id > 0) THEN |
402 |
vid = gdtmaf_id |
403 |
ELSE IF (old_id > 0) THEN |
404 |
vid = old_id |
405 |
ELSE |
406 |
CALL histerr (3, 'flinopen', 'No time axis found',' ',' ') |
407 |
ENDIF |
408 |
!--- |
409 |
ALLOCATE(vec_tmp(ttm)) |
410 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
411 |
start=(/ 1 /), count=(/ ttm /)) |
412 |
itaus(1:ttm) = NINT(vec_tmp(1:ttm)) |
413 |
DEALLOCATE(vec_tmp) |
414 |
!--- |
415 |
!-- Getting all the details for the time axis |
416 |
!--- |
417 |
!-- Find the calendar |
418 |
calendar='XXXX' |
419 |
iret = NF90_GET_ATT (fid, gdtmaf_id, 'calendar', calendar) |
420 |
IF ( INDEX(calendar,'XXXX') < 1 ) THEN |
421 |
CALL ioconf_calendar(calendar) |
422 |
ENDIF |
423 |
!-- |
424 |
units = '' |
425 |
iret = NF90_GET_ATT (fid, vid, 'units', units) |
426 |
IF (gdtt_id > 0) THEN |
427 |
units = units(INDEX(units,'since')+6:LEN_TRIM(units)) |
428 |
READ (units,'(I4.4,5(a,I2.2))') & |
429 |
year0, strc, month0, strc, day0, & |
430 |
strc, hours0, strc, minutes0, strc, seci |
431 |
sec0 = hours0*3600. + minutes0*60. + seci |
432 |
CALL ymds2ju (year0, month0, day0, sec0, date0) |
433 |
iret = NF90_GET_ATT (fid, gdtt_id, 'tstep_sec', dt) |
434 |
ELSE IF (gdtmaf_id > 0) THEN |
435 |
units = units(INDEX(units,'since')+6:LEN_TRIM(units)) |
436 |
READ (units,'(I4.4,5(a,I2.2))') & |
437 |
year0, strc, month0, strc, day0, & |
438 |
strc, hours0, strc, minutes0, strc, seci |
439 |
sec0 = hours0*3600. + minutes0*60. + seci |
440 |
CALL ymds2ju (year0, month0, day0, sec0, date0) |
441 |
!----- |
442 |
ELSE IF (old_id > 0) THEN |
443 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', dt) |
444 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'day0', r_day) |
445 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'sec0', sec) |
446 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'year0', r_year) |
447 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'month0', r_month) |
448 |
!----- |
449 |
day = INT(r_day) |
450 |
month = INT(r_month) |
451 |
year = INT(r_year) |
452 |
!----- |
453 |
CALL ymds2ju (year, month, day, sec, date0) |
454 |
ENDIF |
455 |
ENDIF |
456 |
|
457 |
END SUBROUTINE flinopen_nozoom |
458 |
|
459 |
!=== |
460 |
|
461 |
SUBROUTINE flininfo (filename, iim, jjm, llm, ttm, fid_out) |
462 |
!--------------------------------------------------------------------- |
463 |
!- This subroutine allows to get some information. |
464 |
!- It is usualy done within flinopen but the user may want to call |
465 |
!- it before in order to allocate the space needed to extract the |
466 |
!- data from the file. |
467 |
!--------------------------------------------------------------------- |
468 |
IMPLICIT NONE |
469 |
|
470 |
! ARGUMENTS |
471 |
|
472 |
CHARACTER(LEN=*), intent(in):: filename |
473 |
INTEGER, intent(out):: iim, jjm, llm, ttm, fid_out |
474 |
|
475 |
! LOCAL |
476 |
|
477 |
INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim |
478 |
INTEGER :: iv, lll |
479 |
CHARACTER(LEN=80) :: name |
480 |
CHARACTER(LEN=30) :: axname |
481 |
|
482 |
LOGICAL :: check = .FALSE. |
483 |
!--------------------------------------------------------------------- |
484 |
lll = LEN_TRIM(filename) |
485 |
IF (filename(lll-2:lll) /= '.nc') THEN |
486 |
name = filename(1:lll)//'.nc' |
487 |
ELSE |
488 |
name = filename(1:lll) |
489 |
ENDIF |
490 |
|
491 |
iret = NF90_OPEN (name, NF90_NOWRITE, fid) |
492 |
IF (iret /= NF90_NOERR) THEN |
493 |
CALL histerr(3, 'flininfo','Could not open file :',TRIM(name),' ') |
494 |
ENDIF |
495 |
|
496 |
iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, & |
497 |
nAttributes=nb_atts, unlimitedDimId=id_unlim) |
498 |
|
499 |
iim = 0; |
500 |
jjm = 0; |
501 |
llm = 0; |
502 |
ttm = 0; |
503 |
|
504 |
DO iv=1,ndims |
505 |
!--- |
506 |
iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll) |
507 |
CALL strlowercase (axname) |
508 |
axname = ADJUSTL(axname) |
509 |
!--- |
510 |
IF (check) WRITE(*,*) & |
511 |
'flininfo - getting axname',iv,axname,lll |
512 |
!--- |
513 |
IF ( (INDEX(axname,'x') == 1) & |
514 |
.OR.(INDEX(axname,'lon') == 1) ) THEN |
515 |
iim = lll; |
516 |
ELSE IF ( (INDEX(axname,'y') == 1) & |
517 |
.OR.(INDEX(axname,'lat') == 1) ) THEN |
518 |
jjm = lll; |
519 |
ELSE IF ( (INDEX(axname,'lev') == 1) & |
520 |
.OR.(INDEX(axname,'plev') == 1) & |
521 |
.OR.(INDEX(axname,'z') == 1) & |
522 |
.OR.(INDEX(axname,'depth') == 1) ) THEN |
523 |
llm = lll; |
524 |
ELSE IF ( (INDEX(axname,'tstep') == 1) & |
525 |
.OR.(INDEX(axname,'time_counter') == 1) ) THEN |
526 |
!---- For the time we certainly need to allow for other names |
527 |
ttm = lll; |
528 |
ELSE IF (ndims == 1) THEN |
529 |
!---- Nothing was found and ndims=1 then we have a vector of data |
530 |
iim = lll; |
531 |
ENDIF |
532 |
!--- |
533 |
ENDDO |
534 |
|
535 |
! Keep all this information |
536 |
|
537 |
nbfiles = nbfiles+1 |
538 |
|
539 |
IF (nbfiles > nbfile_max) THEN |
540 |
CALL histerr (3,'flininfo', & |
541 |
'Too many files. Please increase nbfil_max', & |
542 |
'in program flincom.F90.',' ') |
543 |
ENDIF |
544 |
|
545 |
ncids(nbfiles) = fid |
546 |
ncnbd(nbfiles) = ndims |
547 |
|
548 |
ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /) |
549 |
|
550 |
ncfunli(nbfiles) = id_unlim |
551 |
ncnba(nbfiles) = nb_atts |
552 |
ncnbva(nbfiles) = nvars |
553 |
ncfileopen(nbfiles) = .TRUE. |
554 |
|
555 |
fid_out = nbfiles |
556 |
!---------------------- |
557 |
END SUBROUTINE flininfo |
558 |
|
559 |
!=== |
560 |
|
561 |
SUBROUTINE flinput_r1d & |
562 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) |
563 |
!--------------------------------------------------------------------- |
564 |
IMPLICIT NONE |
565 |
|
566 |
INTEGER :: fid_in |
567 |
CHARACTER(LEN=*) :: varname |
568 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
569 |
REAL :: var(:) |
570 |
|
571 |
INTEGER :: fid, ncvarid, ndim, iret |
572 |
LOGICAL :: check = .FALSE. |
573 |
!--------------------------------------------------------------------- |
574 |
IF (check) WRITE(*,*) & |
575 |
"flinput_r1d : SIZE(var) = ",SIZE(var) |
576 |
|
577 |
CALL flinput_mat & |
578 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & |
579 |
fid,ncvarid,ndim) |
580 |
|
581 |
iret = NF90_PUT_VAR (fid, ncvarid, var, & |
582 |
start=w_sta(1:ndim), count=w_len(1:ndim)) |
583 |
!------------------------- |
584 |
END SUBROUTINE flinput_r1d |
585 |
|
586 |
!=== |
587 |
|
588 |
SUBROUTINE flinput_r2d & |
589 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) |
590 |
!--------------------------------------------------------------------- |
591 |
IMPLICIT NONE |
592 |
|
593 |
INTEGER :: fid_in |
594 |
CHARACTER(LEN=*) :: varname |
595 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
596 |
REAL :: var(:,:) |
597 |
|
598 |
INTEGER :: fid, ncvarid, ndim, iret |
599 |
LOGICAL :: check = .FALSE. |
600 |
!--------------------------------------------------------------------- |
601 |
IF (check) WRITE(*,*) & |
602 |
"flinput_r2d : SIZE(var) = ",SIZE(var) |
603 |
|
604 |
CALL flinput_mat & |
605 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & |
606 |
fid,ncvarid,ndim) |
607 |
|
608 |
iret = NF90_PUT_VAR (fid, ncvarid, var, & |
609 |
start=w_sta(1:ndim), count=w_len(1:ndim)) |
610 |
!------------------------- |
611 |
END SUBROUTINE flinput_r2d |
612 |
|
613 |
!=== |
614 |
|
615 |
SUBROUTINE flinput_r3d & |
616 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) |
617 |
!--------------------------------------------------------------------- |
618 |
IMPLICIT NONE |
619 |
|
620 |
INTEGER :: fid_in |
621 |
CHARACTER(LEN=*) :: varname |
622 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
623 |
REAL :: var(:,:,:) |
624 |
|
625 |
INTEGER :: fid, ncvarid, ndim, iret |
626 |
LOGICAL :: check = .FALSE. |
627 |
!--------------------------------------------------------------------- |
628 |
IF (check) WRITE(*,*) & |
629 |
"flinput_r3d : SIZE(var) = ",SIZE(var) |
630 |
|
631 |
CALL flinput_mat & |
632 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & |
633 |
fid,ncvarid,ndim) |
634 |
|
635 |
iret = NF90_PUT_VAR (fid, ncvarid, var, & |
636 |
start=w_sta(1:ndim), count=w_len(1:ndim)) |
637 |
!------------------------- |
638 |
END SUBROUTINE flinput_r3d |
639 |
|
640 |
!=== |
641 |
|
642 |
SUBROUTINE flinput_r4d & |
643 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) |
644 |
!--------------------------------------------------------------------- |
645 |
IMPLICIT NONE |
646 |
|
647 |
INTEGER :: fid_in |
648 |
CHARACTER(LEN=*) :: varname |
649 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
650 |
REAL :: var(:,:,:,:) |
651 |
|
652 |
INTEGER :: fid, ncvarid, ndim, iret |
653 |
LOGICAL :: check = .FALSE. |
654 |
!--------------------------------------------------------------------- |
655 |
IF (check) WRITE(*,*) & |
656 |
"flinput_r4d : SIZE(var) = ",SIZE(var) |
657 |
|
658 |
CALL flinput_mat & |
659 |
(fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & |
660 |
fid,ncvarid,ndim) |
661 |
|
662 |
iret = NF90_PUT_VAR (fid, ncvarid, var, & |
663 |
start=w_sta(1:ndim), count=w_len(1:ndim)) |
664 |
!------------------------- |
665 |
END SUBROUTINE flinput_r4d |
666 |
|
667 |
!=== |
668 |
|
669 |
SUBROUTINE flinput_mat & |
670 |
(fid_in,varname,iim,nlonid,jjm,nlatid, & |
671 |
llm,zdimid,ttm,tdimid,fid,ncvarid,ndim) |
672 |
!--------------------------------------------------------------------- |
673 |
IMPLICIT NONE |
674 |
|
675 |
INTEGER :: fid_in |
676 |
CHARACTER(LEN=*) :: varname |
677 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
678 |
INTEGER :: fid, ncvarid, ndim |
679 |
|
680 |
! LOCAL |
681 |
|
682 |
INTEGER :: iret |
683 |
!--------------------------------------------------------------------- |
684 |
fid = ncids(fid_in) |
685 |
|
686 |
w_sta(1:4) = (/ 1, 1, 1, 1 /) |
687 |
w_len(1:2) = (/ iim, jjm /) |
688 |
w_dim(1:2) = (/ nlonid, nlatid /) |
689 |
|
690 |
IF ( (llm > 0).AND.(ttm > 0) ) THEN |
691 |
ndim = 4 |
692 |
w_len(3:4) = (/ llm, ttm /) |
693 |
w_dim(3:4) = (/ zdimid, tdimid /) |
694 |
ELSE IF (llm > 0) THEN |
695 |
ndim = 3 |
696 |
w_dim(3) = zdimid |
697 |
w_len(3) = llm |
698 |
ELSE IF (ttm > 0) THEN |
699 |
ndim = 3 |
700 |
w_dim(3) = tdimid |
701 |
w_len(3) = ttm |
702 |
ELSE |
703 |
ndim = 2 |
704 |
ENDIF |
705 |
|
706 |
iret = NF90_REDEF (fid) |
707 |
iret = NF90_DEF_VAR (fid,varname,NF90_FLOAT,w_dim(1:ndim),ncvarid) |
708 |
iret = NF90_PUT_ATT (fid,ncvarid,'short_name',TRIM(varname)) |
709 |
iret = NF90_ENDDEF (fid) |
710 |
!-------------------------- |
711 |
END SUBROUTINE flinput_mat |
712 |
|
713 |
!=== |
714 |
|
715 |
SUBROUTINE flinput_scal & |
716 |
(fid_in, varname, iim, nlonid, jjm, nlatid, & |
717 |
llm, zdimid, ttm, tdimid, var) |
718 |
!--------------------------------------------------------------------- |
719 |
IMPLICIT NONE |
720 |
|
721 |
INTEGER :: fid_in |
722 |
CHARACTER(LEN=*) :: varname |
723 |
INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid |
724 |
REAL :: var |
725 |
|
726 |
! LOCAL |
727 |
|
728 |
INTEGER :: fid, iret |
729 |
!--------------------------------------------------------------------- |
730 |
fid = ncids(fid_in) |
731 |
|
732 |
iret = NF90_REDEF (fid) |
733 |
iret = NF90_PUT_ATT (fid, NF90_GLOBAL, varname, REAL(var)) |
734 |
iret = NF90_ENDDEF (fid) |
735 |
!--------------------------- |
736 |
END SUBROUTINE flinput_scal |
737 |
|
738 |
!=== |
739 |
|
740 |
SUBROUTINE flinget_r1d & |
741 |
(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
742 |
!--------------------------------------------------------------------- |
743 |
IMPLICIT NONE |
744 |
|
745 |
INTEGER :: fid_in |
746 |
CHARACTER(LEN=*) :: varname |
747 |
INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin |
748 |
REAL :: var(:) |
749 |
|
750 |
INTEGER :: jl, ji |
751 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
752 |
LOGICAL :: check = .FALSE. |
753 |
!--------------------------------------------------------------------- |
754 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
755 |
IF (check) WRITE(*,*) & |
756 |
"flinget_r1d : allocate buff_tmp for buff_sz = ",SIZE(var) |
757 |
ALLOCATE (buff_tmp(SIZE(var))) |
758 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
759 |
IF (check) WRITE(*,*) & |
760 |
"flinget_r1d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
761 |
DEALLOCATE (buff_tmp) |
762 |
ALLOCATE (buff_tmp(SIZE(var))) |
763 |
ENDIF |
764 |
|
765 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
766 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
767 |
|
768 |
jl=0 |
769 |
DO ji=1,SIZE(var,1) |
770 |
jl=jl+1 |
771 |
var(ji) = buff_tmp(jl) |
772 |
ENDDO |
773 |
!------------------------- |
774 |
END SUBROUTINE flinget_r1d |
775 |
|
776 |
!=== |
777 |
|
778 |
SUBROUTINE flinget_r2d & |
779 |
(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
780 |
!--------------------------------------------------------------------- |
781 |
IMPLICIT NONE |
782 |
|
783 |
INTEGER :: fid_in |
784 |
CHARACTER(LEN=*) :: varname |
785 |
INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin |
786 |
REAL :: var(:,:) |
787 |
|
788 |
INTEGER :: jl, jj, ji |
789 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
790 |
LOGICAL :: check = .FALSE. |
791 |
!--------------------------------------------------------------------- |
792 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
793 |
IF (check) WRITE(*,*) & |
794 |
"flinget_r2d : allocate buff_tmp for buff_sz = ",SIZE(var) |
795 |
ALLOCATE (buff_tmp(SIZE(var))) |
796 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
797 |
IF (check) WRITE(*,*) & |
798 |
"flinget_r2d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
799 |
DEALLOCATE (buff_tmp) |
800 |
ALLOCATE (buff_tmp(SIZE(var))) |
801 |
ENDIF |
802 |
|
803 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
804 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
805 |
|
806 |
jl=0 |
807 |
DO jj=1,SIZE(var,2) |
808 |
DO ji=1,SIZE(var,1) |
809 |
jl=jl+1 |
810 |
var(ji,jj) = buff_tmp(jl) |
811 |
ENDDO |
812 |
ENDDO |
813 |
!------------------------- |
814 |
END SUBROUTINE flinget_r2d |
815 |
|
816 |
!=== |
817 |
|
818 |
SUBROUTINE flinget_r2d_zoom2d & |
819 |
(fid_in,varname,iim,jjm,llm,ttm, & |
820 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) |
821 |
!--------------------------------------------------------------------- |
822 |
IMPLICIT NONE |
823 |
|
824 |
INTEGER :: fid_in |
825 |
CHARACTER(LEN=*) :: varname |
826 |
INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen |
827 |
REAL :: var(:,:) |
828 |
|
829 |
INTEGER :: jl, jj, ji |
830 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
831 |
LOGICAL :: check = .FALSE. |
832 |
!--------------------------------------------------------------------- |
833 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
834 |
IF (check) WRITE(*,*) & |
835 |
"flinget_r2d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) |
836 |
ALLOCATE (buff_tmp(SIZE(var))) |
837 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
838 |
IF (check) WRITE(*,*) & |
839 |
"flinget_r2d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
840 |
DEALLOCATE (buff_tmp) |
841 |
ALLOCATE (buff_tmp(SIZE(var))) |
842 |
ENDIF |
843 |
|
844 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
845 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) |
846 |
|
847 |
jl=0 |
848 |
DO jj=1,SIZE(var,2) |
849 |
DO ji=1,SIZE(var,1) |
850 |
jl=jl+1 |
851 |
var(ji,jj) = buff_tmp(jl) |
852 |
ENDDO |
853 |
ENDDO |
854 |
!-------------------------------- |
855 |
END SUBROUTINE flinget_r2d_zoom2d |
856 |
|
857 |
!=== |
858 |
|
859 |
SUBROUTINE flinget_r3d(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
860 |
!--------------------------------------------------------------------- |
861 |
IMPLICIT NONE |
862 |
|
863 |
INTEGER, intent(in):: fid_in |
864 |
CHARACTER(LEN=*), intent(in):: varname |
865 |
INTEGER, intent(in):: iim, jjm, llm, ttm, itau_dep, itau_fin |
866 |
REAL, intent(out):: var(:,:,:) |
867 |
|
868 |
INTEGER :: jl, jk, jj, ji |
869 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
870 |
LOGICAL :: check = .FALSE. |
871 |
!--------------------------------------------------------------------- |
872 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
873 |
IF (check) WRITE(*,*) & |
874 |
"flinget_r3d : allocate buff_tmp for buff_sz = ",SIZE(var) |
875 |
ALLOCATE (buff_tmp(SIZE(var))) |
876 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
877 |
IF (check) WRITE(*,*) & |
878 |
"flinget_r3d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
879 |
DEALLOCATE (buff_tmp) |
880 |
ALLOCATE (buff_tmp(SIZE(var))) |
881 |
ENDIF |
882 |
|
883 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
884 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
885 |
|
886 |
jl=0 |
887 |
DO jk=1,SIZE(var,3) |
888 |
DO jj=1,SIZE(var,2) |
889 |
DO ji=1,SIZE(var,1) |
890 |
jl=jl+1 |
891 |
var(ji,jj,jk) = buff_tmp(jl) |
892 |
ENDDO |
893 |
ENDDO |
894 |
ENDDO |
895 |
!------------------------- |
896 |
END SUBROUTINE flinget_r3d |
897 |
|
898 |
!=== |
899 |
|
900 |
SUBROUTINE flinget_r3d_zoom2d & |
901 |
(fid_in,varname,iim,jjm,llm,ttm, & |
902 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) |
903 |
!--------------------------------------------------------------------- |
904 |
IMPLICIT NONE |
905 |
|
906 |
INTEGER :: fid_in |
907 |
CHARACTER(LEN=*) :: varname |
908 |
INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen |
909 |
REAL :: var(:,:,:) |
910 |
|
911 |
INTEGER :: jl, jk, jj, ji |
912 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
913 |
LOGICAL :: check = .FALSE. |
914 |
!--------------------------------------------------------------------- |
915 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
916 |
IF (check) WRITE(*,*) & |
917 |
"flinget_r3d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) |
918 |
ALLOCATE (buff_tmp(SIZE(var))) |
919 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
920 |
IF (check) WRITE(*,*) & |
921 |
"flinget_r3d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
922 |
DEALLOCATE (buff_tmp) |
923 |
ALLOCATE (buff_tmp(SIZE(var))) |
924 |
ENDIF |
925 |
|
926 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
927 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) |
928 |
|
929 |
jl=0 |
930 |
DO jk=1,SIZE(var,3) |
931 |
DO jj=1,SIZE(var,2) |
932 |
DO ji=1,SIZE(var,1) |
933 |
jl=jl+1 |
934 |
var(ji,jj,jk) = buff_tmp(jl) |
935 |
ENDDO |
936 |
ENDDO |
937 |
ENDDO |
938 |
!-------------------------------- |
939 |
END SUBROUTINE flinget_r3d_zoom2d |
940 |
|
941 |
!=== |
942 |
|
943 |
SUBROUTINE flinget_r4d & |
944 |
(fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) |
945 |
!--------------------------------------------------------------------- |
946 |
IMPLICIT NONE |
947 |
|
948 |
INTEGER :: fid_in |
949 |
CHARACTER(LEN=*) :: varname |
950 |
INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin |
951 |
REAL :: var(:,:,:,:) |
952 |
|
953 |
INTEGER :: jl, jk, jj, ji, jm |
954 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
955 |
LOGICAL :: check = .FALSE. |
956 |
!--------------------------------------------------------------------- |
957 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
958 |
IF (check) WRITE(*,*) & |
959 |
"flinget_r4d : allocate buff_tmp for buff_sz = ",SIZE(var) |
960 |
ALLOCATE (buff_tmp(SIZE(var))) |
961 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
962 |
IF (check) WRITE(*,*) & |
963 |
"flinget_r4d : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
964 |
DEALLOCATE (buff_tmp) |
965 |
ALLOCATE (buff_tmp(SIZE(var))) |
966 |
ENDIF |
967 |
|
968 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
969 |
itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) |
970 |
|
971 |
jl=0 |
972 |
DO jm=1,SIZE(var,4) |
973 |
DO jk=1,SIZE(var,3) |
974 |
DO jj=1,SIZE(var,2) |
975 |
DO ji=1,SIZE(var,1) |
976 |
jl=jl+1 |
977 |
var(ji,jj,jk,jm) = buff_tmp(jl) |
978 |
ENDDO |
979 |
ENDDO |
980 |
ENDDO |
981 |
ENDDO |
982 |
!------------------------- |
983 |
END SUBROUTINE flinget_r4d |
984 |
|
985 |
!=== |
986 |
|
987 |
SUBROUTINE flinget_r4d_zoom2d & |
988 |
(fid_in,varname,iim,jjm,llm,ttm, & |
989 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) |
990 |
!--------------------------------------------------------------------- |
991 |
IMPLICIT NONE |
992 |
|
993 |
INTEGER :: fid_in |
994 |
CHARACTER(LEN=*) :: varname |
995 |
INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen |
996 |
REAL :: var(:,:,:,:) |
997 |
|
998 |
INTEGER :: jl, jk, jj, ji, jm |
999 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp |
1000 |
LOGICAL :: check = .FALSE. |
1001 |
!--------------------------------------------------------------------- |
1002 |
IF (.NOT.ALLOCATED(buff_tmp)) THEN |
1003 |
IF (check) WRITE(*,*) & |
1004 |
"flinget_r4d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) |
1005 |
ALLOCATE (buff_tmp(SIZE(var))) |
1006 |
ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN |
1007 |
IF (check) WRITE(*,*) & |
1008 |
"flinget_r4d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) |
1009 |
DEALLOCATE (buff_tmp) |
1010 |
ALLOCATE (buff_tmp(SIZE(var))) |
1011 |
ENDIF |
1012 |
|
1013 |
CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & |
1014 |
itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) |
1015 |
|
1016 |
jl=0 |
1017 |
DO jm=1,SIZE(var,4) |
1018 |
DO jk=1,SIZE(var,3) |
1019 |
DO jj=1,SIZE(var,2) |
1020 |
DO ji=1,SIZE(var,1) |
1021 |
jl=jl+1 |
1022 |
var(ji,jj,jk,jm) = buff_tmp(jl) |
1023 |
ENDDO |
1024 |
ENDDO |
1025 |
ENDDO |
1026 |
ENDDO |
1027 |
!-------------------------------- |
1028 |
END SUBROUTINE flinget_r4d_zoom2d |
1029 |
|
1030 |
!=== |
1031 |
|
1032 |
SUBROUTINE flinget_mat & |
1033 |
(fid_in, varname, iim, jjm, llm, ttm, itau_dep, & |
1034 |
itau_fin, iideb, iilen, jjdeb, jjlen, var) |
1035 |
!--------------------------------------------------------------------- |
1036 |
!- This subroutine will read the variable named varname from |
1037 |
!- the file previously opened by flinopen and identified by fid |
1038 |
|
1039 |
!- It is checked that the dimensions of the variable to be read |
1040 |
!- correspond to what the user requested when he specified |
1041 |
!- iim, jjm and llm. The only exception which is allowed is |
1042 |
!- for compressed data where the horizontal grid is not expected |
1043 |
!- to be iim x jjm. |
1044 |
|
1045 |
!- If variable is of size zero a global attribute is read. |
1046 |
!- This global attribute will be of type real |
1047 |
|
1048 |
!- INPUT |
1049 |
|
1050 |
!- fid : File ID returned by flinopen |
1051 |
!- varname : Name of the variable to be read from the file |
1052 |
!- iim : | These three variables give the size of the variables |
1053 |
!- jjm : | to be read. It will be verified that the variables |
1054 |
!- llm : | fits in there. |
1055 |
!- ttm : | |
1056 |
!- itau_dep : Time step at which we will start to read |
1057 |
!- itau_fin : Time step until which we are going to read |
1058 |
!- For the moment this is done on indexes |
1059 |
!- but it should be in the physical space. |
1060 |
!- If there is no time-axis in the file then use a |
1061 |
!- itau_fin < itau_dep, this will tell flinget not to |
1062 |
!- expect a time-axis in the file. |
1063 |
!- iideb : index i for zoom |
1064 |
!- iilen : length of zoom |
1065 |
!- jjdeb : index j for zoom |
1066 |
!- jjlen : length of zoom |
1067 |
|
1068 |
!- OUTPUT |
1069 |
|
1070 |
!- var : array that will contain the data |
1071 |
!--------------------------------------------------------------------- |
1072 |
IMPLICIT NONE |
1073 |
|
1074 |
! ARGUMENTS |
1075 |
|
1076 |
INTEGER, intent(in):: fid_in |
1077 |
CHARACTER(LEN=*), intent(in):: varname |
1078 |
INTEGER, intent(in):: iim, jjm, llm, ttm, itau_dep, itau_fin |
1079 |
INTEGER :: iideb |
1080 |
integer, intent(in):: iilen |
1081 |
integer jjdeb |
1082 |
integer, intent(in):: jjlen |
1083 |
REAL :: var(:) |
1084 |
|
1085 |
! LOCAL |
1086 |
|
1087 |
INTEGER :: iret, fid |
1088 |
INTEGER :: vid, cvid, clen |
1089 |
CHARACTER(LEN=70) :: str1 |
1090 |
CHARACTER(LEN=250) :: att_n, tmp_n |
1091 |
CHARACTER(LEN=5) :: axs_l |
1092 |
INTEGER :: tmp_i |
1093 |
REAL,SAVE :: mis_v=0. |
1094 |
REAL :: tmp_r |
1095 |
INTEGER :: ndims, x_typ, nb_atts |
1096 |
INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: dimids |
1097 |
INTEGER :: i, iv, nvars, i2d, cnd |
1098 |
REAL,DIMENSION(:),ALLOCATABLE,SAVE :: var_tmp |
1099 |
LOGICAL :: uncompress = .FALSE. |
1100 |
LOGICAL :: check = .FALSE. |
1101 |
!--------------------------------------------------------------------- |
1102 |
fid = ncids(fid_in) |
1103 |
|
1104 |
IF (check) THEN |
1105 |
WRITE(*,*) & |
1106 |
'flinget_mat : fid_in, fid, varname :', fid_in, fid, TRIM(varname) |
1107 |
WRITE(*,*) & |
1108 |
'flinget_mat : iim, jjm, llm, ttm, itau_dep, itau_fin :', & |
1109 |
iim, jjm, llm, ttm, itau_dep, itau_fin |
1110 |
WRITE(*,*) & |
1111 |
'flinget_mat : iideb, iilen, jjdeb, jjlen :', & |
1112 |
iideb, iilen, jjdeb, jjlen |
1113 |
ENDIF |
1114 |
|
1115 |
uncompress = .FALSE. |
1116 |
|
1117 |
! 1.0 We get first all the details on this variable from the file |
1118 |
|
1119 |
nvars = ncnbva(fid_in) |
1120 |
|
1121 |
vid = -1 |
1122 |
iret = NF90_INQ_VARID (fid, varname, vid) |
1123 |
|
1124 |
IF (vid < 0 .OR. iret /= NF90_NOERR) THEN |
1125 |
CALL histerr (3,'flinget', & |
1126 |
'Variable '//TRIM(varname)//' not found in file',' ',' ') |
1127 |
ENDIF |
1128 |
|
1129 |
iret = NF90_INQUIRE_VARIABLE (fid, vid, & |
1130 |
ndims=ndims, dimids=dimids, nAtts=nb_atts) |
1131 |
IF (check) THEN |
1132 |
WRITE(*,*) & |
1133 |
'flinget_mat : fid, vid :', fid, vid |
1134 |
WRITE(*,*) & |
1135 |
'flinget_mat : ndims, dimids(1:ndims), nb_atts :', & |
1136 |
ndims, dimids(1:ndims), nb_atts |
1137 |
ENDIF |
1138 |
|
1139 |
w_dim(:) = 0 |
1140 |
DO i=1,ndims |
1141 |
iret = NF90_INQUIRE_DIMENSION (fid, dimids(i), len=w_dim(i)) |
1142 |
ENDDO |
1143 |
IF (check) WRITE(*,*) & |
1144 |
'flinget_mat : w_dim :', w_dim(1:ndims) |
1145 |
|
1146 |
mis_v = 0.0; axs_l = ' '; |
1147 |
|
1148 |
IF (nb_atts > 0) THEN |
1149 |
IF (check) THEN |
1150 |
WRITE(*,*) 'flinget_mat : attributes for variable :' |
1151 |
ENDIF |
1152 |
ENDIF |
1153 |
DO i=1,nb_atts |
1154 |
iret = NF90_INQ_ATTNAME (fid, vid, i, att_n) |
1155 |
iret = NF90_INQUIRE_ATTRIBUTE (fid, vid, att_n, xtype=x_typ) |
1156 |
CALL strlowercase (att_n) |
1157 |
IF ( (x_typ == NF90_INT).OR.(x_typ == NF90_SHORT) & |
1158 |
.OR.(x_typ == NF90_BYTE) ) THEN |
1159 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_i) |
1160 |
IF (check) THEN |
1161 |
WRITE(*,*) ' ',TRIM(att_n),' : ',tmp_i |
1162 |
ENDIF |
1163 |
ELSE IF ( (x_typ == NF90_FLOAT).OR.(x_typ == NF90_DOUBLE) ) THEN |
1164 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_r) |
1165 |
IF (check) THEN |
1166 |
WRITE(*,*) ' ',TRIM(att_n),' : ',tmp_r |
1167 |
ENDIF |
1168 |
IF (index(att_n,'missing_value') > 0) THEN |
1169 |
mis_v = tmp_r |
1170 |
ENDIF |
1171 |
ELSE |
1172 |
tmp_n = '' |
1173 |
iret = NF90_GET_ATT (fid, vid, att_n, tmp_n) |
1174 |
IF (check) THEN |
1175 |
WRITE(*,*) ' ',TRIM(att_n),' : ',TRIM(tmp_n) |
1176 |
ENDIF |
1177 |
IF (index(att_n,'axis') > 0) THEN |
1178 |
axs_l = tmp_n |
1179 |
ENDIF |
1180 |
ENDIF |
1181 |
ENDDO |
1182 |
!? |
1183 |
!!!!!!!!!! We will need a verification on the type of the variable |
1184 |
!? |
1185 |
|
1186 |
! 2.0 The dimensions are analysed to determine what is to be read |
1187 |
|
1188 |
! 2.1 the longitudes |
1189 |
|
1190 |
IF ( w_dim(1) /= iim .OR. w_dim(2) /= jjm) THEN |
1191 |
!--- |
1192 |
!-- There is a possibility that we have to deal with a compressed axis ! |
1193 |
!--- |
1194 |
iret = NF90_INQUIRE_DIMENSION (fid, dimids(1), & |
1195 |
name=tmp_n, len=clen) |
1196 |
iret = NF90_INQ_VARID (fid, tmp_n, cvid) |
1197 |
!--- |
1198 |
IF (check) WRITE(*,*) & |
1199 |
'Dimname, iret , NF90_NOERR : ',TRIM(tmp_n),iret,NF90_NOERR |
1200 |
!--- |
1201 |
!-- If we have an axis which has the same name |
1202 |
!-- as the dimension we can see if it is compressed |
1203 |
!--- |
1204 |
!-- TODO TODO for zoom2d |
1205 |
!--- |
1206 |
IF (iret == NF90_NOERR) THEN |
1207 |
iret = NF90_GET_ATT (fid, cvid, 'compress', str1) |
1208 |
!----- |
1209 |
IF (iret == NF90_NOERR) THEN |
1210 |
iret = NF90_INQUIRE_VARIABLE (fid,cvid,xtype=x_typ,ndims=cnd) |
1211 |
!------- |
1212 |
IF ( cnd /= 1 .AND. x_typ /= NF90_INT) THEN |
1213 |
CALL histerr (3,'flinget', & |
1214 |
'Variable '//TRIM(tmp_n)//' can not be a compressed axis', & |
1215 |
'Either it has too many dimensions'// & |
1216 |
' or it is not of type integer', ' ') |
1217 |
ELSE |
1218 |
!--------- |
1219 |
!-------- Let us see if we already have that index table |
1220 |
!--------- |
1221 |
IF ( (cind_len /= clen).OR.(cind_vid /= cvid) & |
1222 |
.OR.(cind_fid /= fid) ) THEN |
1223 |
IF (ALLOCATED(cindex)) DEALLOCATE(cindex) |
1224 |
ALLOCATE(cindex(clen)) |
1225 |
cind_len = clen |
1226 |
cind_vid = cvid |
1227 |
cind_fid = fid |
1228 |
iret = NF90_GET_VAR (fid, cvid, cindex) |
1229 |
ENDIF |
1230 |
!--------- |
1231 |
!-------- In any case we need to set the slab of data to be read |
1232 |
!--------- |
1233 |
uncompress = .TRUE. |
1234 |
w_sta(1) = 1 |
1235 |
w_len(1) = clen |
1236 |
i2d = 1 |
1237 |
ENDIF |
1238 |
ELSE |
1239 |
str1 = 'The horizontal dimensions of '//varname |
1240 |
CALL histerr (3,'flinget',str1, & |
1241 |
'is not compressed and does not'// & |
1242 |
' correspond to the requested size',' ') |
1243 |
ENDIF |
1244 |
ELSE |
1245 |
IF (w_dim(1) /= iim) THEN |
1246 |
str1 = 'The longitude dimension of '//varname |
1247 |
CALL histerr (3,'flinget',str1, & |
1248 |
'in the file is not equal to the dimension', & |
1249 |
'that should be read') |
1250 |
ENDIF |
1251 |
IF (w_dim(2) /= jjm) THEN |
1252 |
str1 = 'The latitude dimension of '//varname |
1253 |
CALL histerr (3,'flinget',str1, & |
1254 |
'in the file is not equal to the dimension', & |
1255 |
'that should be read') |
1256 |
ENDIF |
1257 |
ENDIF |
1258 |
ELSE |
1259 |
w_sta(1:2) = (/ iideb, jjdeb /) |
1260 |
w_len(1:2) = (/ iilen, jjlen /) |
1261 |
i2d = 2 |
1262 |
ENDIF |
1263 |
|
1264 |
! 2.3 Now the difficult part, the 3rd dimension which can be |
1265 |
! time or levels. |
1266 |
|
1267 |
! Priority is given to the time axis if only three axes are present. |
1268 |
|
1269 |
IF (ndims > i2d) THEN |
1270 |
!--- |
1271 |
!-- 2.3.1 We have a vertical axis |
1272 |
!--- |
1273 |
IF (llm == 1 .AND. ndims == i2d+2 .OR. llm == w_dim(i2d+1)) THEN |
1274 |
!----- |
1275 |
IF (w_dim(i2d+1) /= llm) THEN |
1276 |
CALL histerr (3,'flinget', & |
1277 |
'The vertical dimension of '//varname, & |
1278 |
'in the file is not equal to the dimension', & |
1279 |
'that should be read') |
1280 |
ELSE |
1281 |
w_sta(i2d+1) = 1 |
1282 |
IF (llm > 0) THEN |
1283 |
w_len(i2d+1) = llm |
1284 |
ELSE |
1285 |
w_len(i2d+1) = w_sta(i2d+1) |
1286 |
ENDIF |
1287 |
ENDIF |
1288 |
!----- |
1289 |
IF ((itau_fin-itau_dep) >= 0) THEN |
1290 |
IF (ndims /= i2d+2) THEN |
1291 |
CALL histerr (3,'flinget', & |
1292 |
'You attempt to read a time slab', & |
1293 |
'but there is no time axis on this variable', varname) |
1294 |
ELSE IF ((itau_fin - itau_dep) <= w_dim(i2d+2)) THEN |
1295 |
w_sta(i2d+2) = itau_dep |
1296 |
w_len(i2d+2) = itau_fin-itau_dep+1 |
1297 |
ELSE |
1298 |
CALL histerr (3,'flinget', & |
1299 |
'The time step you try to read is not', & |
1300 |
'in the file (1)', varname) |
1301 |
ENDIF |
1302 |
ELSE IF (ndims == i2d+2 .AND. w_dim(i2d+2) > 1) THEN |
1303 |
CALL histerr (3,'flinget', & |
1304 |
'There is a time axis in the file but no', & |
1305 |
'time step give in the call', varname) |
1306 |
ELSE |
1307 |
w_sta(i2d+2) = 1 |
1308 |
w_len(i2d+2) = 1 |
1309 |
ENDIF |
1310 |
ELSE |
1311 |
!----- |
1312 |
!---- 2.3.2 We do not have any vertical axis |
1313 |
!----- |
1314 |
IF (ndims == i2d+2) THEN |
1315 |
CALL histerr (3,'flinget', & |
1316 |
'The file contains 4 dimensions', & |
1317 |
'but only 3 are requestes for variable ', varname) |
1318 |
ENDIF |
1319 |
IF ((itau_fin-itau_dep) >= 0) THEN |
1320 |
IF (ndims == i2d+1) THEN |
1321 |
IF ((itau_fin-itau_dep) < w_dim(i2d+1) ) THEN |
1322 |
w_sta(i2d+1) = itau_dep |
1323 |
w_len(i2d+1) = itau_fin-itau_dep+1 |
1324 |
ELSE |
1325 |
CALL histerr (3,'flinget', & |
1326 |
'The time step you try to read is not', & |
1327 |
'in the file (2)', varname) |
1328 |
ENDIF |
1329 |
ELSE |
1330 |
CALL histerr (3,'flinget', & |
1331 |
'From your input you sould have 3 dimensions', & |
1332 |
'in the file but there are 4', varname) |
1333 |
ENDIF |
1334 |
ELSE |
1335 |
IF (ndims == i2d+1 .AND. w_dim(i2d+1) > 1) THEN |
1336 |
CALL histerr (3,'flinget', & |
1337 |
'There is a time axis in the file but no', & |
1338 |
'time step given in the call', varname) |
1339 |
ELSE |
1340 |
w_sta(i2d+1) = 1 |
1341 |
w_len(i2d+1) = 1 |
1342 |
ENDIF |
1343 |
ENDIF |
1344 |
ENDIF |
1345 |
ELSE |
1346 |
!--- |
1347 |
!-- 2.3.3 We do not have any vertical axis |
1348 |
!--- |
1349 |
w_sta(i2d+1:i2d+2) = (/ 0, 0 /) |
1350 |
w_len(i2d+1:i2d+2) = (/ 0, 0 /) |
1351 |
ENDIF |
1352 |
|
1353 |
! 3.0 Reading the data |
1354 |
|
1355 |
IF (check) WRITE(*,*) & |
1356 |
'flinget_mat 3.0 : ', uncompress, w_sta, w_len |
1357 |
!--- |
1358 |
IF (uncompress) THEN |
1359 |
!--- |
1360 |
IF (ALLOCATED(var_tmp)) THEN |
1361 |
IF (SIZE(var_tmp) < clen) THEN |
1362 |
DEALLOCATE(var_tmp) |
1363 |
ALLOCATE(var_tmp(clen)) |
1364 |
ENDIF |
1365 |
ELSE |
1366 |
ALLOCATE(var_tmp(clen)) |
1367 |
ENDIF |
1368 |
!--- |
1369 |
iret = NF90_GET_VAR (fid, vid, var_tmp, & |
1370 |
start=w_sta(:), count=w_len(:)) |
1371 |
!--- |
1372 |
var(:) = mis_v |
1373 |
var(cindex(:)) = var_tmp(:) |
1374 |
!--- |
1375 |
ELSE |
1376 |
iret = NF90_GET_VAR (fid, vid, var, & |
1377 |
start=w_sta(:), count=w_len(:)) |
1378 |
ENDIF |
1379 |
|
1380 |
IF (check) WRITE(*,*) 'flinget_mat 3.1 : ',NF90_STRERROR (iret) |
1381 |
!-------------------------- |
1382 |
END SUBROUTINE flinget_mat |
1383 |
|
1384 |
!=== |
1385 |
|
1386 |
SUBROUTINE flinget_scal & |
1387 |
(fid_in, varname, iim, jjm, llm, ttm, itau_dep, itau_fin, var) |
1388 |
!--------------------------------------------------------------------- |
1389 |
!- This subroutine will read the variable named varname from |
1390 |
!- the file previously opened by flinopen and identified by fid |
1391 |
|
1392 |
!- If variable is of size zero a global attribute is read. This |
1393 |
!- global attribute will be of type real |
1394 |
|
1395 |
!- INPUT |
1396 |
|
1397 |
!- fid : File ID returned by flinopen |
1398 |
!- varname : Name of the variable to be read from the file |
1399 |
!- iim : | These three variables give the size of the variables |
1400 |
!- jjm : | to be read. It will be verified that the variables |
1401 |
!- llm : | fits in there. |
1402 |
!- ttm : | |
1403 |
!- itau_dep : Time step at which we will start to read |
1404 |
!- itau_fin : Time step until which we are going to read |
1405 |
!- For the moment this is done on indeces but it should be |
1406 |
!- in the physical space |
1407 |
!- If there is no time-axis in the file then use a |
1408 |
!- itau_fin < itau_dep, this will tell flinget not to |
1409 |
!- expect a time-axis in the file. |
1410 |
|
1411 |
!- OUTPUT |
1412 |
|
1413 |
!- var : scalar that will contain the data |
1414 |
!--------------------------------------------------------------------- |
1415 |
IMPLICIT NONE |
1416 |
|
1417 |
! ARGUMENTS |
1418 |
|
1419 |
INTEGER :: fid_in |
1420 |
CHARACTER(LEN=*) :: varname |
1421 |
INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin |
1422 |
REAL :: var |
1423 |
|
1424 |
! LOCAL |
1425 |
|
1426 |
INTEGER :: iret, fid |
1427 |
|
1428 |
LOGICAL :: check = .FALSE. |
1429 |
!--------------------------------------------------------------------- |
1430 |
fid = ncids(fid_in) |
1431 |
|
1432 |
! 1.0 Reading a global attribute |
1433 |
|
1434 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, varname, var) |
1435 |
!--------------------------- |
1436 |
END SUBROUTINE flinget_scal |
1437 |
|
1438 |
!=== |
1439 |
|
1440 |
SUBROUTINE flinfindcood (fid_in, axtype, vid, ndim) |
1441 |
!--------------------------------------------------------------------- |
1442 |
!- This subroutine explores the file in order to find |
1443 |
!- the coordinate according to a number of rules |
1444 |
!--------------------------------------------------------------------- |
1445 |
IMPLICIT NONE |
1446 |
|
1447 |
! ARGUMENTS |
1448 |
|
1449 |
INTEGER, intent(in):: fid_in |
1450 |
integer vid, ndim |
1451 |
CHARACTER(LEN=3) :: axtype |
1452 |
|
1453 |
! LOCAL |
1454 |
|
1455 |
INTEGER :: iv, iret, dimnb |
1456 |
CHARACTER(LEN=40) :: dimname, dimuni1, dimuni2, dimuni3 |
1457 |
CHARACTER(LEN=30) :: str1 |
1458 |
LOGICAL :: found_rule = .FALSE. |
1459 |
!--------------------------------------------------------------------- |
1460 |
vid = -1 |
1461 |
|
1462 |
! Make sure all strings are invalid |
1463 |
|
1464 |
dimname = '?-?' |
1465 |
dimuni1 = '?-?' |
1466 |
dimuni2 = '?-?' |
1467 |
dimuni3 = '?-?' |
1468 |
|
1469 |
! First rule : we look for the correct units |
1470 |
! lon : east |
1471 |
! lat : north |
1472 |
! We make an exact check as it would be too easy to mistake |
1473 |
! some units by just comparing the substrings. |
1474 |
|
1475 |
SELECTCASE(axtype) |
1476 |
CASE ('lon') |
1477 |
dimuni1 = 'degree_e' |
1478 |
dimuni2 = 'degrees_e' |
1479 |
found_rule = .TRUE. |
1480 |
CASE('lat') |
1481 |
dimuni1 = 'degree_n' |
1482 |
dimuni2 = 'degrees_n' |
1483 |
found_rule = .TRUE. |
1484 |
CASE('lev') |
1485 |
dimuni1 = 'm' |
1486 |
dimuni2 = 'km' |
1487 |
dimuni3 = 'hpa' |
1488 |
found_rule = .TRUE. |
1489 |
CASE DEFAULT |
1490 |
found_rule = .FALSE. |
1491 |
END SELECT |
1492 |
|
1493 |
IF (found_rule) THEN |
1494 |
iv = 0 |
1495 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
1496 |
iv = iv+1 |
1497 |
str1 = '' |
1498 |
iret = NF90_GET_ATT (ncids(fid_in), iv, 'units', str1) |
1499 |
IF (iret == NF90_NOERR) THEN |
1500 |
CALL strlowercase (str1) |
1501 |
IF ( (INDEX(str1, TRIM(dimuni1)) == 1) & |
1502 |
.OR.(INDEX(str1, TRIM(dimuni2)) == 1) & |
1503 |
.OR.(INDEX(str1, TRIM(dimuni3)) == 1) ) THEN |
1504 |
vid = iv |
1505 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, ndims=ndim) |
1506 |
ENDIF |
1507 |
ENDIF |
1508 |
ENDDO |
1509 |
ENDIF |
1510 |
|
1511 |
! Second rule : we find specific names : |
1512 |
! lon : nav_lon |
1513 |
! lat : nav_lat |
1514 |
! Here we can check if we find the substring as the |
1515 |
! names are more specific. |
1516 |
|
1517 |
SELECTCASE(axtype) |
1518 |
CASE ('lon') |
1519 |
dimname = 'nav_lon lon longitude' |
1520 |
found_rule = .TRUE. |
1521 |
CASE('lat') |
1522 |
dimname = 'nav_lat lat latitude' |
1523 |
found_rule = .TRUE. |
1524 |
CASE('lev') |
1525 |
dimname = 'plev level depth deptht' |
1526 |
found_rule = .TRUE. |
1527 |
CASE DEFAULT |
1528 |
found_rule = .FALSE. |
1529 |
END SELECT |
1530 |
|
1531 |
IF (found_rule) THEN |
1532 |
iv = 0 |
1533 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
1534 |
iv = iv+1 |
1535 |
str1='' |
1536 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & |
1537 |
name=str1, ndims=ndim) |
1538 |
IF (INDEX(dimname,TRIM(str1)) >= 1) THEN |
1539 |
vid = iv |
1540 |
ENDIF |
1541 |
ENDDO |
1542 |
ENDIF |
1543 |
|
1544 |
! Third rule : we find a variable with the same name as the dimension |
1545 |
! lon = 1 |
1546 |
! lat = 2 |
1547 |
! lev = 3 |
1548 |
|
1549 |
IF (vid < 0) THEN |
1550 |
SELECTCASE(axtype) |
1551 |
CASE ('lon') |
1552 |
dimnb = 1 |
1553 |
found_rule = .TRUE. |
1554 |
CASE('lat') |
1555 |
dimnb = 2 |
1556 |
found_rule = .TRUE. |
1557 |
CASE('lev') |
1558 |
dimnb = 3 |
1559 |
found_rule = .TRUE. |
1560 |
CASE DEFAULT |
1561 |
found_rule = .FALSE. |
1562 |
END SELECT |
1563 |
!--- |
1564 |
IF (found_rule) THEN |
1565 |
iret = NF90_INQUIRE_DIMENSION (ncids(fid_in), dimnb, name=dimname) |
1566 |
iv = 0 |
1567 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
1568 |
iv = iv+1 |
1569 |
str1='' |
1570 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & |
1571 |
name=str1, ndims=ndim) |
1572 |
IF (INDEX(dimname,TRIM(str1)) == 1) THEN |
1573 |
vid = iv |
1574 |
ENDIF |
1575 |
ENDDO |
1576 |
ENDIF |
1577 |
ENDIF |
1578 |
|
1579 |
! Stop the program if no coordinate was found |
1580 |
|
1581 |
IF (vid < 0) THEN |
1582 |
CALL histerr (3,'flinfindcood', & |
1583 |
'No coordinate axis was found in the file', & |
1584 |
'The data in this file can not be used', axtype) |
1585 |
ENDIF |
1586 |
!-------------------------- |
1587 |
END SUBROUTINE flinfindcood |
1588 |
|
1589 |
!=== |
1590 |
|
1591 |
SUBROUTINE flinclo (fid_in) |
1592 |
!--------------------------------------------------------------------- |
1593 |
IMPLICIT NONE |
1594 |
|
1595 |
INTEGER :: fid_in |
1596 |
|
1597 |
INTEGER :: iret |
1598 |
!--------------------------------------------------------------------- |
1599 |
iret = NF90_CLOSE (ncids(fid_in)) |
1600 |
ncfileopen(fid_in) = .FALSE. |
1601 |
!--------------------- |
1602 |
END SUBROUTINE flinclo |
1603 |
|
1604 |
!=== |
1605 |
|
1606 |
SUBROUTINE flinquery_var(fid_in, varname, exists) |
1607 |
!--------------------------------------------------------------------- |
1608 |
!- Queries the existance of a variable in the file. |
1609 |
!--------------------------------------------------------------------- |
1610 |
IMPLICIT NONE |
1611 |
|
1612 |
INTEGER :: fid_in |
1613 |
CHARACTER(LEN=*) varname |
1614 |
LOGICAL :: exists |
1615 |
|
1616 |
INTEGER :: iret, fid, vid |
1617 |
!--------------------------------------------------------------------- |
1618 |
fid = ncids(fid_in) |
1619 |
vid = -1 |
1620 |
iret = NF90_INQ_VARID (fid, varname, vid) |
1621 |
|
1622 |
exists = ( (vid >= 0).AND.(iret == NF90_NOERR) ) |
1623 |
!--------------------------- |
1624 |
END SUBROUTINE flinquery_var |
1625 |
|
1626 |
!=== |
1627 |
|
1628 |
SUBROUTINE flininspect (fid_in) |
1629 |
!--------------------------------------------------------------------- |
1630 |
IMPLICIT NONE |
1631 |
|
1632 |
! fid : File id to inspect |
1633 |
|
1634 |
INTEGER :: fid_in |
1635 |
|
1636 |
!- LOCAL |
1637 |
|
1638 |
INTEGER :: iim, jjm, llm, ttm |
1639 |
INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim |
1640 |
INTEGER :: iv, in, lll |
1641 |
INTEGER :: xid, yid, zid, tid |
1642 |
INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: idimid |
1643 |
CHARACTER(LEN=80) :: name |
1644 |
CHARACTER(LEN=30) :: axname |
1645 |
!--------------------------------------------------------------------- |
1646 |
fid = ncids(fid_in) |
1647 |
|
1648 |
iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, & |
1649 |
nAttributes=nb_atts, unlimitedDimId=id_unlim) |
1650 |
|
1651 |
WRITE (*,*) 'IOIPSL ID : ',fid_in |
1652 |
WRITE (*,*) 'NetCDF ID : ',fid |
1653 |
WRITE (*,*) 'Number of dimensions : ',ndims |
1654 |
WRITE (*,*) 'Number of variables : ',nvars |
1655 |
WRITE (*,*) 'Number of global attributes : ',nb_atts |
1656 |
WRITE (*,*) 'ID unlimited : ',id_unlim |
1657 |
|
1658 |
xid = -1; iim = 0; |
1659 |
yid = -1; jjm = 0; |
1660 |
zid = -1; llm = 0; |
1661 |
tid = -1; ttm = 0; |
1662 |
|
1663 |
DO iv=1,ndims |
1664 |
!--- |
1665 |
iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll) |
1666 |
CALL strlowercase (axname) |
1667 |
axname = ADJUSTL(axname) |
1668 |
!--- |
1669 |
WRITE (*,*) 'Dimension number : ',iv |
1670 |
WRITE (*,*) 'Dimension name : ',TRIM(axname) |
1671 |
!--- |
1672 |
IF ( (INDEX(axname,'x') == 1) & |
1673 |
.OR.(INDEX(axname,'lon') == 1)) THEN |
1674 |
xid = iv; iim = lll; |
1675 |
WRITE (*,*) 'Dimension X size : ',iim |
1676 |
ELSE IF ( (INDEX(axname,'y') == 1) & |
1677 |
.OR.(INDEX(axname,'lat') == 1)) THEN |
1678 |
yid = iv; jjm = lll; |
1679 |
WRITE (*,*) 'Dimension Y size : ',jjm |
1680 |
ELSE IF ( (INDEX(axname,'lev') == 1) & |
1681 |
.OR.(INDEX(axname,'plev') == 1) & |
1682 |
.OR.(INDEX(axname,'z') == 1) & |
1683 |
.OR.(INDEX(axname,'depth') == 1)) THEN |
1684 |
zid = iv; llm = lll; |
1685 |
WRITE (*,*) 'Dimension Z size : ',llm |
1686 |
ELSE IF ( (INDEX(axname,'tstep') == 1) & |
1687 |
.OR.(INDEX(axname,'time_counter') == 1)) THEN |
1688 |
!---- For the time we certainly need to allow for other names |
1689 |
tid = iv; ttm = lll; |
1690 |
ELSE IF (ndims == 1) THEN |
1691 |
!---- Nothing was found and ndims=1 then we have a vector of data |
1692 |
xid = 1; iim = lll; |
1693 |
ENDIF |
1694 |
!--- |
1695 |
ENDDO |
1696 |
|
1697 |
! Keep all this information |
1698 |
|
1699 |
nbfiles = nbfiles+1 |
1700 |
|
1701 |
IF (nbfiles > nbfile_max) THEN |
1702 |
CALL histerr(3,'flininspect', & |
1703 |
'Too many files. Please increase nbfil_max', & |
1704 |
'in program flincom.F90.',' ') |
1705 |
ENDIF |
1706 |
|
1707 |
ncids(nbfiles) = fid |
1708 |
ncnbd(nbfiles) = ndims |
1709 |
|
1710 |
ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /) |
1711 |
|
1712 |
ncfunli(nbfiles) = id_unlim |
1713 |
ncnba(nbfiles) = nb_atts |
1714 |
ncnbva(nbfiles) = nvars |
1715 |
ncfileopen(nbfiles) = .TRUE. |
1716 |
|
1717 |
DO in=1,nvars |
1718 |
iret = NF90_INQUIRE_VARIABLE (fid, in, & |
1719 |
name=name, ndims=ndims, dimids=idimid, nAtts=nb_atts) |
1720 |
WRITE (*,*) 'Variable number ------------ > ', in |
1721 |
WRITE (*,*) 'Variable name : ', TRIM(name) |
1722 |
WRITE (*,*) 'Number of dimensions : ', ndims |
1723 |
WRITE (*,*) 'Dimensions ID''s : ', idimid(1:ndims) |
1724 |
WRITE (*,*) 'Number of attributes : ', nb_atts |
1725 |
ENDDO |
1726 |
!------------------------- |
1727 |
END SUBROUTINE flininspect |
1728 |
|
1729 |
!=== |
1730 |
|
1731 |
END MODULE flincom |