| 1 |
MODULE flincom |
| 2 |
|
| 3 |
! From flincom.f90, version 2.2 2006/03/07 09:21:51 |
| 4 |
|
| 5 |
IMPLICIT NONE |
| 6 |
|
| 7 |
PRIVATE |
| 8 |
PUBLIC flinclo, flinopen_nozoom, flininfo, ncids |
| 9 |
|
| 10 |
! This is the data we keep on each file we open: |
| 11 |
INTEGER, PARAMETER:: nbfile_max = 200 |
| 12 |
INTEGER, SAVE:: ncids(nbfile_max) |
| 13 |
INTEGER, SAVE:: ncnbva(nbfile_max) |
| 14 |
LOGICAL, SAVE:: ncfileopen(nbfile_max)=.FALSE. |
| 15 |
|
| 16 |
CONTAINS |
| 17 |
|
| 18 |
SUBROUTINE flinopen_nozoom(iim, jjm, llm, lon, lat, lev, & |
| 19 |
ttm, itaus, date0, dt, fid_out) |
| 20 |
|
| 21 |
! The routine will open an input file |
| 22 |
! INPUT |
| 23 |
! There is no test of the content of the file against the input |
| 24 |
! from the model |
| 25 |
|
| 26 |
! iim: size in the x direction in the file (longitude) |
| 27 |
! jjm: size in the y direction |
| 28 |
! llm: number of levels |
| 29 |
! (llm = 0 means no axis to be expected) |
| 30 |
|
| 31 |
! WARNING: |
| 32 |
! It is for the user to check |
| 33 |
! that the dimensions of lon lat and lev are correct when passed to |
| 34 |
! flinopen. This can be done after the call when iim and jjm have |
| 35 |
! been retrieved from the netCDF file. In F90 this problem will |
| 36 |
! be solved with an internal assign |
| 37 |
! IF iim, jjm, llm or ttm are parameters in the calling program |
| 38 |
! it will create a segmentation fault |
| 39 |
|
| 40 |
! ttm: size of time axis |
| 41 |
|
| 42 |
! OUTPUT |
| 43 |
|
| 44 |
! lon: array of (iim, jjm), |
| 45 |
! that contains the longitude of each point |
| 46 |
! lat: same for latitude |
| 47 |
! lev: An array of llm for the latitude |
| 48 |
! itaus: Time steps within this file |
| 49 |
! date0: Julian date at which itau = 0 |
| 50 |
! dt: length of the time steps of the data |
| 51 |
|
| 52 |
!--------------------------------------------------------------------- |
| 53 |
|
| 54 |
USE calendar, ONLY: ymds2ju, ioconf_calendar |
| 55 |
USE errioipsl, ONLY: histerr |
| 56 |
USE netcdf, ONLY: nf90_get_att, nf90_get_var, nf90_global, & |
| 57 |
nf90_inquire_variable |
| 58 |
|
| 59 |
! ARGUMENTS |
| 60 |
|
| 61 |
INTEGER, intent(in):: iim, jjm, llm, ttm |
| 62 |
real, intent(out):: lon(iim, jjm), lat(iim, jjm), lev(llm) |
| 63 |
INTEGER, intent(out):: itaus(ttm) |
| 64 |
REAL, intent(out):: date0, dt |
| 65 |
|
| 66 |
INTEGER, intent(in):: fid_out |
| 67 |
! (file ID which is later used to read the data) |
| 68 |
|
| 69 |
! LOCAL |
| 70 |
|
| 71 |
INTEGER:: iret, vid, fid, nbdim, i |
| 72 |
INTEGER:: gdtt_id, old_id, iv, gdtmaf_id |
| 73 |
CHARACTER(LEN=250):: name |
| 74 |
CHARACTER(LEN=80):: units, my_calendar |
| 75 |
INTEGER:: year, month, day |
| 76 |
REAL:: r_year, r_month, r_day |
| 77 |
INTEGER:: year0, month0, day0, hours0, minutes0, seci |
| 78 |
REAL:: sec, sec0 |
| 79 |
CHARACTER:: strc |
| 80 |
|
| 81 |
REAL, DIMENSION(:), ALLOCATABLE:: vec_tmp |
| 82 |
|
| 83 |
!--------------------------------------------------------------------- |
| 84 |
|
| 85 |
IF ( (fid_out < 1).OR.(fid_out > nbfile_max) ) THEN |
| 86 |
! Either the fid_out has not been initialized (0 or very large) |
| 87 |
! then we have to open anyway. Else we only need to open the file |
| 88 |
! if it has not been opened before. |
| 89 |
print *, "Call flinfo before flinopen" |
| 90 |
stop 1 |
| 91 |
end IF |
| 92 |
IF (.NOT.ncfileopen(fid_out)) THEN |
| 93 |
print *, "Call flinfo before flinopen" |
| 94 |
stop 1 |
| 95 |
end IF |
| 96 |
|
| 97 |
! The user has already opened the file |
| 98 |
! and we trust that he knows the dimensions |
| 99 |
|
| 100 |
fid = ncids(fid_out) |
| 101 |
|
| 102 |
! 2.0 get the sizes and names of the different coordinates |
| 103 |
! and do a first set of verification. |
| 104 |
|
| 105 |
! 3.0 Check if we are realy talking about the same coodinate system |
| 106 |
! if not then we get the lon, lat and lev variables from the file |
| 107 |
|
| 108 |
! 4.0 extracting the coordinates |
| 109 |
|
| 110 |
CALL flinfindcood (fid_out, 'lon', vid, nbdim) |
| 111 |
IF (nbdim == 2) THEN |
| 112 |
iret = NF90_GET_VAR (fid, vid, lon, & |
| 113 |
start=(/ 1, 1 /), count=(/ iim, jjm /)) |
| 114 |
ELSE |
| 115 |
ALLOCATE(vec_tmp(iim)) |
| 116 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
| 117 |
start=(/ 1 /), count=(/ iim /)) |
| 118 |
DO i=1, jjm |
| 119 |
lon(:, i) = vec_tmp(:) |
| 120 |
ENDDO |
| 121 |
DEALLOCATE(vec_tmp) |
| 122 |
ENDIF |
| 123 |
|
| 124 |
CALL flinfindcood (fid_out, 'lat', vid, nbdim) |
| 125 |
IF (nbdim == 2) THEN |
| 126 |
iret = NF90_GET_VAR (fid, vid, lat, & |
| 127 |
start=(/ 1, 1 /), count=(/ iim, jjm /)) |
| 128 |
ELSE |
| 129 |
ALLOCATE(vec_tmp(jjm)) |
| 130 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
| 131 |
start=(/ 1 /), count=(/ jjm /)) |
| 132 |
DO i=1, iim |
| 133 |
lat(i, :) = vec_tmp(:) |
| 134 |
ENDDO |
| 135 |
DEALLOCATE(vec_tmp) |
| 136 |
ENDIF |
| 137 |
|
| 138 |
IF (llm > 0) THEN |
| 139 |
CALL flinfindcood (fid_out, 'lev', vid, nbdim) |
| 140 |
IF (nbdim == 1) THEN |
| 141 |
iret = NF90_GET_VAR (fid, vid, lev, & |
| 142 |
start=(/ 1 /), count=(/ llm /)) |
| 143 |
ELSE |
| 144 |
CALL histerr (3, 'flinopen', & |
| 145 |
'Can not handle vertical coordinates that have more', & |
| 146 |
'than 1 dimension', ' ') |
| 147 |
ENDIF |
| 148 |
ENDIF |
| 149 |
|
| 150 |
! 5.0 Get all the details for the time if possible needed |
| 151 |
|
| 152 |
IF (ttm > 0) THEN |
| 153 |
|
| 154 |
! 5.1 Find the time axis. Prefered method is the 'timestep since' |
| 155 |
|
| 156 |
gdtmaf_id = -1 |
| 157 |
gdtt_id = -1 |
| 158 |
old_id = -1 |
| 159 |
DO iv=1, ncnbva(fid_out) |
| 160 |
name='' |
| 161 |
iret = NF90_INQUIRE_VARIABLE (fid, iv, name=name) |
| 162 |
units='' |
| 163 |
iret = NF90_GET_ATT (fid, iv, 'units', units) |
| 164 |
IF (INDEX(units, 'seconds since') > 0) gdtmaf_id = iv |
| 165 |
IF (INDEX(units, 'timesteps since') > 0) gdtt_id = iv |
| 166 |
IF (INDEX(name, 'tstep') > 0) old_id = iv |
| 167 |
ENDDO |
| 168 |
|
| 169 |
IF (gdtt_id > 0) THEN |
| 170 |
vid = gdtt_id |
| 171 |
ELSE IF (gdtmaf_id > 0) THEN |
| 172 |
vid = gdtmaf_id |
| 173 |
ELSE IF (old_id > 0) THEN |
| 174 |
vid = old_id |
| 175 |
ELSE |
| 176 |
CALL histerr (3, 'flinopen', 'No time axis found', ' ', ' ') |
| 177 |
ENDIF |
| 178 |
|
| 179 |
ALLOCATE(vec_tmp(ttm)) |
| 180 |
iret = NF90_GET_VAR (fid, vid, vec_tmp, & |
| 181 |
start=(/ 1 /), count=(/ ttm /)) |
| 182 |
itaus(1:ttm) = NINT(vec_tmp(1:ttm)) |
| 183 |
DEALLOCATE(vec_tmp) |
| 184 |
|
| 185 |
! Getting all the details for the time axis |
| 186 |
|
| 187 |
! Find the calendar |
| 188 |
my_calendar='XXXX' |
| 189 |
iret = NF90_GET_ATT (fid, gdtmaf_id, 'calendar', my_calendar) |
| 190 |
IF ( INDEX(my_calendar, 'XXXX') < 1 ) THEN |
| 191 |
CALL ioconf_calendar(my_calendar) |
| 192 |
ENDIF |
| 193 |
|
| 194 |
units = '' |
| 195 |
iret = NF90_GET_ATT (fid, vid, 'units', units) |
| 196 |
IF (gdtt_id > 0) THEN |
| 197 |
units = units(INDEX(units, 'since')+6:LEN_TRIM(units)) |
| 198 |
READ (units, '(I4.4, 5(a, I2.2))') & |
| 199 |
year0, strc, month0, strc, day0, & |
| 200 |
strc, hours0, strc, minutes0, strc, seci |
| 201 |
sec0 = hours0*3600. + minutes0*60. + seci |
| 202 |
CALL ymds2ju (year0, month0, day0, sec0, date0) |
| 203 |
iret = NF90_GET_ATT (fid, gdtt_id, 'tstep_sec', dt) |
| 204 |
ELSE IF (gdtmaf_id > 0) THEN |
| 205 |
units = units(INDEX(units, 'since')+6:LEN_TRIM(units)) |
| 206 |
READ (units, '(I4.4, 5(a, I2.2))') & |
| 207 |
year0, strc, month0, strc, day0, & |
| 208 |
strc, hours0, strc, minutes0, strc, seci |
| 209 |
sec0 = hours0*3600. + minutes0*60. + seci |
| 210 |
CALL ymds2ju (year0, month0, day0, sec0, date0) |
| 211 |
|
| 212 |
ELSE IF (old_id > 0) THEN |
| 213 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', dt) |
| 214 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'day0', r_day) |
| 215 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'sec0', sec) |
| 216 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'year0', r_year) |
| 217 |
iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'month0', r_month) |
| 218 |
|
| 219 |
day = INT(r_day) |
| 220 |
month = INT(r_month) |
| 221 |
year = INT(r_year) |
| 222 |
|
| 223 |
CALL ymds2ju (year, month, day, sec, date0) |
| 224 |
ENDIF |
| 225 |
ENDIF |
| 226 |
|
| 227 |
END SUBROUTINE flinopen_nozoom |
| 228 |
|
| 229 |
!*************************************************************** |
| 230 |
|
| 231 |
SUBROUTINE flininfo(filename, iim, jjm, llm, ttm, fid_out) |
| 232 |
|
| 233 |
! This subroutine allows to get some information. |
| 234 |
! It is usualy done within flinopen but the user may want to call |
| 235 |
! it before in order to allocate the space needed to extract the |
| 236 |
! data from the file. |
| 237 |
|
| 238 |
USE strlowercase_m, ONLY: strlowercase |
| 239 |
USE errioipsl, ONLY: histerr |
| 240 |
USE netcdf, ONLY: nf90_inquire, nf90_inquire_dimension, nf90_noerr, & |
| 241 |
nf90_nowrite |
| 242 |
use netcdf95, only: nf95_open |
| 243 |
|
| 244 |
CHARACTER(LEN=*), intent(in):: filename |
| 245 |
INTEGER, intent(out):: iim, jjm, llm, ttm, fid_out |
| 246 |
|
| 247 |
! LOCAL |
| 248 |
|
| 249 |
INTEGER, SAVE:: nbfiles = 0 |
| 250 |
INTEGER, SAVE:: ncdims(nbfile_max, 4) |
| 251 |
INTEGER:: iret, fid, ndims, nvars, nb_atts, id_unlim |
| 252 |
INTEGER:: iv, lll |
| 253 |
CHARACTER(LEN=80):: name |
| 254 |
CHARACTER(LEN=30):: axname |
| 255 |
|
| 256 |
!--------------------------------------------------------------------- |
| 257 |
|
| 258 |
lll = LEN_TRIM(filename) |
| 259 |
IF (filename(lll-2:lll) /= '.nc') THEN |
| 260 |
name = filename(1:lll)//'.nc' |
| 261 |
ELSE |
| 262 |
name = filename(1:lll) |
| 263 |
ENDIF |
| 264 |
|
| 265 |
call NF95_OPEN(name, NF90_NOWRITE, fid) |
| 266 |
iret = NF90_INQUIRE(fid, nDimensions=ndims, nVariables=nvars, & |
| 267 |
nAttributes=nb_atts, unlimitedDimId=id_unlim) |
| 268 |
|
| 269 |
iim = 0 |
| 270 |
jjm = 0 |
| 271 |
llm = 0 |
| 272 |
ttm = 0 |
| 273 |
|
| 274 |
DO iv=1, ndims |
| 275 |
iret = NF90_INQUIRE_DIMENSION(fid, iv, name=axname, len=lll) |
| 276 |
CALL strlowercase(axname) |
| 277 |
axname = ADJUSTL(axname) |
| 278 |
|
| 279 |
IF ((INDEX(axname, 'x') == 1) .OR. (INDEX(axname, 'lon') == 1)) THEN |
| 280 |
iim = lll |
| 281 |
ELSE IF ((INDEX(axname, 'y') == 1) & |
| 282 |
.OR. (INDEX(axname, 'lat') == 1)) THEN |
| 283 |
jjm = lll |
| 284 |
ELSE IF ((INDEX(axname, 'lev') == 1) .OR. (INDEX(axname, 'plev') == 1) & |
| 285 |
.OR. (INDEX(axname, 'z') == 1) & |
| 286 |
.OR. (INDEX(axname, 'depth') == 1)) THEN |
| 287 |
llm = lll |
| 288 |
ELSE IF ((INDEX(axname, 'tstep') == 1) & |
| 289 |
.OR. (INDEX(axname, 'time_counter') == 1)) THEN |
| 290 |
! For the time we certainly need to allow for other names |
| 291 |
ttm = lll |
| 292 |
ELSE IF (ndims == 1) THEN |
| 293 |
! Nothing was found and ndims=1 then we have a vector of data |
| 294 |
iim = lll |
| 295 |
ENDIF |
| 296 |
ENDDO |
| 297 |
|
| 298 |
! Keep all this information |
| 299 |
|
| 300 |
nbfiles = nbfiles+1 |
| 301 |
|
| 302 |
IF (nbfiles > nbfile_max) THEN |
| 303 |
CALL histerr(3, 'flininfo', & |
| 304 |
'Too many files. Please increase nbfil_max', & |
| 305 |
'in program flincom.F90.', ' ') |
| 306 |
ENDIF |
| 307 |
|
| 308 |
ncids(nbfiles) = fid |
| 309 |
ncdims(nbfiles, :) = (/ iim, jjm, llm, ttm /) |
| 310 |
ncnbva(nbfiles) = nvars |
| 311 |
ncfileopen(nbfiles) = .TRUE. |
| 312 |
fid_out = nbfiles |
| 313 |
|
| 314 |
END SUBROUTINE flininfo |
| 315 |
|
| 316 |
!*************************************************************** |
| 317 |
|
| 318 |
SUBROUTINE flinfindcood (fid_in, axtype, vid, ndim) |
| 319 |
|
| 320 |
! This subroutine explores the file in order to find |
| 321 |
! the coordinate according to a number of rules |
| 322 |
|
| 323 |
USE strlowercase_m, ONLY: strlowercase |
| 324 |
USE errioipsl, ONLY: histerr |
| 325 |
USE netcdf, ONLY: nf90_get_att, nf90_inquire_dimension, & |
| 326 |
nf90_inquire_variable, nf90_noerr |
| 327 |
|
| 328 |
! ARGUMENTS |
| 329 |
|
| 330 |
INTEGER, intent(in):: fid_in |
| 331 |
integer vid, ndim |
| 332 |
CHARACTER(LEN=3):: axtype |
| 333 |
|
| 334 |
! LOCAL |
| 335 |
|
| 336 |
INTEGER:: iv, iret, dimnb |
| 337 |
CHARACTER(LEN=40):: dimname, dimuni1, dimuni2, dimuni3 |
| 338 |
CHARACTER(LEN=30):: str1 |
| 339 |
LOGICAL:: found_rule = .FALSE. |
| 340 |
!--------------------------------------------------------------------- |
| 341 |
vid = -1 |
| 342 |
|
| 343 |
! Make sure all strings are invalid |
| 344 |
|
| 345 |
dimname = '?-?' |
| 346 |
dimuni1 = '?-?' |
| 347 |
dimuni2 = '?-?' |
| 348 |
dimuni3 = '?-?' |
| 349 |
|
| 350 |
! First rule: we look for the correct units |
| 351 |
! lon: east |
| 352 |
! lat: north |
| 353 |
! We make an exact check as it would be too easy to mistake |
| 354 |
! some units by just comparing the substrings. |
| 355 |
|
| 356 |
SELECTCASE(axtype) |
| 357 |
CASE ('lon') |
| 358 |
dimuni1 = 'degree_e' |
| 359 |
dimuni2 = 'degrees_e' |
| 360 |
found_rule = .TRUE. |
| 361 |
CASE('lat') |
| 362 |
dimuni1 = 'degree_n' |
| 363 |
dimuni2 = 'degrees_n' |
| 364 |
found_rule = .TRUE. |
| 365 |
CASE('lev') |
| 366 |
dimuni1 = 'm' |
| 367 |
dimuni2 = 'km' |
| 368 |
dimuni3 = 'hpa' |
| 369 |
found_rule = .TRUE. |
| 370 |
CASE DEFAULT |
| 371 |
found_rule = .FALSE. |
| 372 |
END SELECT |
| 373 |
|
| 374 |
IF (found_rule) THEN |
| 375 |
iv = 0 |
| 376 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
| 377 |
iv = iv+1 |
| 378 |
str1 = '' |
| 379 |
iret = NF90_GET_ATT (ncids(fid_in), iv, 'units', str1) |
| 380 |
IF (iret == NF90_NOERR) THEN |
| 381 |
CALL strlowercase (str1) |
| 382 |
IF ( (INDEX(str1, TRIM(dimuni1)) == 1) & |
| 383 |
.OR.(INDEX(str1, TRIM(dimuni2)) == 1) & |
| 384 |
.OR.(INDEX(str1, TRIM(dimuni3)) == 1) ) THEN |
| 385 |
vid = iv |
| 386 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, ndims=ndim) |
| 387 |
ENDIF |
| 388 |
ENDIF |
| 389 |
ENDDO |
| 390 |
ENDIF |
| 391 |
|
| 392 |
! Second rule: we find specific names: |
| 393 |
! lon: nav_lon |
| 394 |
! lat: nav_lat |
| 395 |
! Here we can check if we find the substring as the |
| 396 |
! names are more specific. |
| 397 |
|
| 398 |
SELECTCASE(axtype) |
| 399 |
CASE ('lon') |
| 400 |
dimname = 'nav_lon lon longitude' |
| 401 |
found_rule = .TRUE. |
| 402 |
CASE('lat') |
| 403 |
dimname = 'nav_lat lat latitude' |
| 404 |
found_rule = .TRUE. |
| 405 |
CASE('lev') |
| 406 |
dimname = 'plev level depth deptht' |
| 407 |
found_rule = .TRUE. |
| 408 |
CASE DEFAULT |
| 409 |
found_rule = .FALSE. |
| 410 |
END SELECT |
| 411 |
|
| 412 |
IF (found_rule) THEN |
| 413 |
iv = 0 |
| 414 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
| 415 |
iv = iv+1 |
| 416 |
str1='' |
| 417 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & |
| 418 |
name=str1, ndims=ndim) |
| 419 |
IF (INDEX(dimname, TRIM(str1)) >= 1) THEN |
| 420 |
vid = iv |
| 421 |
ENDIF |
| 422 |
ENDDO |
| 423 |
ENDIF |
| 424 |
|
| 425 |
! Third rule: we find a variable with the same name as the dimension |
| 426 |
! lon = 1 |
| 427 |
! lat = 2 |
| 428 |
! lev = 3 |
| 429 |
|
| 430 |
IF (vid < 0) THEN |
| 431 |
SELECTCASE(axtype) |
| 432 |
CASE ('lon') |
| 433 |
dimnb = 1 |
| 434 |
found_rule = .TRUE. |
| 435 |
CASE('lat') |
| 436 |
dimnb = 2 |
| 437 |
found_rule = .TRUE. |
| 438 |
CASE('lev') |
| 439 |
dimnb = 3 |
| 440 |
found_rule = .TRUE. |
| 441 |
CASE DEFAULT |
| 442 |
found_rule = .FALSE. |
| 443 |
END SELECT |
| 444 |
|
| 445 |
IF (found_rule) THEN |
| 446 |
iret = NF90_INQUIRE_DIMENSION (ncids(fid_in), dimnb, name=dimname) |
| 447 |
iv = 0 |
| 448 |
DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) |
| 449 |
iv = iv+1 |
| 450 |
str1='' |
| 451 |
iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & |
| 452 |
name=str1, ndims=ndim) |
| 453 |
IF (INDEX(dimname, TRIM(str1)) == 1) THEN |
| 454 |
vid = iv |
| 455 |
ENDIF |
| 456 |
ENDDO |
| 457 |
ENDIF |
| 458 |
ENDIF |
| 459 |
|
| 460 |
! Stop the program if no coordinate was found |
| 461 |
|
| 462 |
IF (vid < 0) THEN |
| 463 |
CALL histerr (3, 'flinfindcood', & |
| 464 |
'No coordinate axis was found in the file', & |
| 465 |
'The data in this file can not be used', axtype) |
| 466 |
ENDIF |
| 467 |
|
| 468 |
END SUBROUTINE flinfindcood |
| 469 |
|
| 470 |
!*************************************************************** |
| 471 |
|
| 472 |
SUBROUTINE flinclo (fid_in) |
| 473 |
|
| 474 |
USE netcdf, ONLY: nf90_close |
| 475 |
|
| 476 |
INTEGER:: fid_in |
| 477 |
|
| 478 |
INTEGER:: iret |
| 479 |
|
| 480 |
!--------------------------------------------------------------------- |
| 481 |
|
| 482 |
iret = NF90_CLOSE (ncids(fid_in)) |
| 483 |
ncfileopen(fid_in) = .FALSE. |
| 484 |
|
| 485 |
END SUBROUTINE flinclo |
| 486 |
|
| 487 |
END MODULE flincom |
Parent Directory
| 
Revision Log