wrf-phdf5attr.F90
References to this file elsewhere.
1 !***************************************************************************
2 !* The HDF5 WRF IO module was written by the the HDF Group at NCSA, the *
3 !* National Center for Supercomputing Applications. *
4 !* HDF Group *
5 !* National Center for Supercomputing Applications *
6 !* University of Illinois at Urbana-Champaign *
7 !* 605 E. Springfield, Champaign IL 61820 *
8 !* http://hdf.ncsa.uiuc.edu/ *
9 !* *
10 !* Copyright 2004 by the Board of Trustees, University of Illinois, *
11 !* *
12 !* Redistribution or use of this IO module, with or without modification, *
13 !* is permitted for any purpose, including commercial purposes. *
14 !* *
15 !* This software is an unsupported prototype. Use at your own risk. *
16 !* http://hdf.ncsa.uiuc.edu/apps/WRF-ROMS *
17 !* *
18 !* This work was funded by the MEAD expedition at the National Center *
19 !* for Supercomputing Applications, NCSA. For more information see: *
20 !* http://www.ncsa.uiuc.edu/expeditions/MEAD *
21 !* *
22 !* *
23 !****************************************************************************/
24
25 module get_attrid_routine
26
27 Interface get_attrid
28 module procedure get_attrid
29 end interface
30
31 contains
32
33 subroutine get_attrid(DataHandle,Element,h5_attrid,Status,VAR)
34
35 use wrf_phdf5_data
36 use ext_phdf5_support_routines
37 USE HDF5 ! This module contains all necessary modules
38 implicit none
39 include 'wrf_status_codes.h'
40
41 character*(*) ,intent(in) :: Element
42 integer ,intent(in) :: DataHandle
43 integer(hid_t) ,intent(out) :: h5_attrid
44 integer(hid_t) :: dset_id
45 integer ,intent(out) :: Status
46 character*(*) ,intent(in),optional :: VAR
47 integer(hid_t) :: hdf5err
48 type(wrf_phdf5_data_handle),pointer :: DH
49
50 character(Len = MaxTimeSLen) :: tname
51 character(Len = 512) :: tgroupname
52 integer(hid_t) :: tgroup_id
53
54 call GetDH(DataHandle,DH,Status)
55 if(Status /= WRF_NO_ERR) then
56 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
57 call wrf_debug ( WARN , msg)
58 return
59 endif
60
61 if(present(VAR)) then
62 call numtochar(1,tname)
63 tgroupname = 'TIME_STAMP_'//tname
64 call h5gopen_f(DH%GroupID,tgroupname,tgroup_id,hdf5err)
65 if(hdf5err.lt.0) then
66 Status = WRF_HDF5_ERR_GROUP
67 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
68 call wrf_debug ( WARN , msg)
69 return
70 endif
71 call h5dopen_f(tgroup_id,VAR,dset_id,hdf5err)
72 if(hdf5err.lt.0) then
73 Status = WRF_HDF5_ERR_DATASET_OPEN
74 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
75 call wrf_debug ( WARN , msg)
76 return
77 endif
78
79 call h5aopen_name_f(dset_id,Element,h5_attrid,hdf5err)
80 if(hdf5err.lt.0) then
81 Status = WRF_HDF5_ERR_ATTRIBUTE_OPEN
82 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
83 call wrf_debug ( WARN , msg)
84 return
85 endif
86 call h5dclose_f(dset_id,hdf5err)
87 call h5gclose_f(tgroup_id,hdf5err)
88 else
89 call h5aopen_name_f(DH%GroupID,Element,h5_attrid,hdf5err)
90 write(*,*) "Element ",Element
91 if(hdf5err.lt.0) then
92 Status = WRF_HDF5_ERR_ATTRIBUTE_OPEN
93 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
94 call wrf_debug ( WARN , msg)
95 return
96 endif
97 endif
98 return
99 end subroutine get_attrid
100 end module get_attrid_routine
101
102 subroutine create_phdf5_objid(DataHandle,obj_id,routine_type,var,Status)
103
104 use wrf_phdf5_data
105 use ext_phdf5_support_routines
106 use HDF5
107 implicit none
108 include 'wrf_status_codes.h'
109
110 integer :: i
111 integer ,intent(in) :: DataHandle
112 integer(hid_t) ,intent(out) :: obj_id
113 character*3 ,intent(in) :: routine_type
114 character*(*) ,intent(in) :: var
115 integer ,intent(out) :: Status
116 integer(hid_t) :: hdf5err
117 type(wrf_phdf5_data_handle),pointer :: DH
118
119
120 call GetDH(DataHandle,DH,Status)
121 if(Status /= WRF_NO_ERR) then
122 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,' ',ROUTINE_TYPE,', line', __LINE__
123 call wrf_debug ( WARN , msg)
124 return
125 endif
126
127 if(routine_type == 'DOM') then
128
129 if(DH%FileStatus == WRF_FILE_OPENED_AND_COMMITTED) then
130 obj_id = DH%GroupID
131 endif
132
133 else if(routine_type == 'VAR') then
134
135 if(DH%FileStatus == WRF_FILE_OPENED_AND_COMMITTED) then
136 do i = 1, MaxVars
137 if(DH%VarNames(i) == var) then
138 obj_id = DH%dsetids(i)
139 exit
140 endif
141 enddo
142 endif
143
144 else
145 Status = WRF_HDF5_ERR_DATA_ID_NOTFOUND
146 write(msg,*) 'CANNOT FIND DATASET ID of the attribute in',__FILE__,&
147 ', line',__LINE__
148 endif
149
150 return
151 end subroutine create_phdf5_objid
152
153
154 subroutine create_phdf5_adtypeid(h5_atypeid,routine_datatype,Count,Status,DataHandle)
155
156 use wrf_phdf5_data
157 use ext_phdf5_support_routines
158 use HDF5
159 implicit none
160 include 'wrf_status_codes.h'
161
162 integer :: i
163 integer(hid_t) ,intent(out) :: h5_atypeid
164 integer ,intent(in) :: routine_datatype
165 integer ,intent(in) :: Count
166 integer ,intent(out) :: Status
167 integer(hid_t) :: hdf5err
168 integer, intent(in), optional :: DataHandle
169 integer(size_t) :: count_size
170
171 type(wrf_phdf5_data_handle),pointer :: DH
172
173 if(routine_datatype == WRF_LOGICAL)then
174 call GetDH(DataHandle,DH,Status)
175 if(Status /= WRF_NO_ERR) then
176 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
177 call wrf_debug ( WARN , msg)
178 return
179 endif
180
181 endif
182
183 select case(routine_datatype)
184 case (WRF_REAL)
185 h5_atypeid = H5T_NATIVE_REAL
186 case (WRF_DOUBLE)
187 h5_atypeid = H5T_NATIVE_DOUBLE
188 case (WRF_INTEGER)
189 h5_atypeid = H5T_NATIVE_INTEGER
190 case (WRF_LOGICAL)
191 h5_atypeid = DH%EnumID
192 case (WRF_CHARACTER)
193
194 call h5tcopy_f(H5T_NATIVE_CHARACTER,h5_atypeid,hdf5err)
195 if(hdf5err.lt.0) then
196 Status = WRF_HDF5_ERR_DATATYPE
197 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
198 call wrf_debug ( WARN , msg)
199 return
200 endif
201
202 count_size = count
203 call h5tset_size_f(h5_atypeid,count_size,hdf5err)
204 if(hdf5err.lt.0) then
205 Status = WRF_HDF5_ERR_DATATYPE
206 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
207 call wrf_debug ( WARN , msg)
208 return
209 endif
210
211 call h5tset_strpad_f(h5_atypeid,H5T_STR_SPACEPAD_F,hdf5err)
212 if(hdf5err.lt.0) then
213 Status = WRF_HDF5_ERR_DATATYPE
214 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
215 call wrf_debug ( WARN , msg)
216 return
217 endif
218
219 case default
220 Status = WRF_HDF5_ERR_DATATYPE
221 return
222 end select
223
224 Status = WRF_NO_ERR
225
226 return
227 end subroutine create_phdf5_adtypeid
228
229 subroutine create_phdf5_adspaceid(Count,str_flag,h5_aspaceid,Status)
230
231 use wrf_phdf5_data
232 use HDF5
233 implicit none
234 include 'wrf_status_codes.h'
235
236 integer :: i
237 integer ,intent(in) :: Count
238 integer ,intent(in) :: str_flag
239 integer ,intent(out) :: Status
240
241 integer(hsize_t) , dimension(1) :: adims
242 integer(hid_t) :: hdf5err
243 integer(hid_t) ,intent(out) :: h5_aspaceid
244 integer :: arank = 1
245
246 ! if string, count is always 1
247 if(str_flag == 1) then
248 adims(1) = 1
249 call h5screate_simple_f(arank,adims,h5_aspaceid,hdf5err)
250 if(hdf5err.lt.0) then
251 Status = WRF_HDF5_ERR_DATASPACE
252 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
253 call wrf_debug ( WARN , msg)
254 return
255 endif
256
257 else
258 adims(1) = Count
259 call h5screate_simple_f(arank,adims,h5_aspaceid,hdf5err)
260 if(hdf5err.lt.0) then
261 Status = WRF_HDF5_ERR_DATASPACE
262 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
263 call wrf_debug ( WARN , msg)
264 return
265 endif
266
267 endif
268
269 Status = WRF_NO_ERR
270
271 return
272 end subroutine create_phdf5_adspaceid
273
274
275 subroutine clean_phdf5_attrids(h5_attr_typeid,h5_space_typeid, &
276 h5_attrid,str_flag,Status)
277
278 use wrf_phdf5_data
279 use HDF5
280 implicit none
281 include 'wrf_status_codes.h'
282 integer ,intent(out) :: Status
283 integer(hid_t) ,intent(in) :: h5_attr_typeid
284 integer(hid_t) ,intent(in) :: h5_space_typeid
285 integer(hid_t) ,intent(in) :: h5_attrid
286 integer ,intent(in) :: str_flag
287 integer :: hdf5err
288
289 if(str_flag == 1) then
290 call h5tclose_f(h5_attr_typeid,hdf5err)
291 if(hdf5err.lt.0) then
292 Status = WRF_HDF5_ERR_CLOSE_GENERAL
293 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
294 call wrf_debug ( WARN , msg)
295 return
296 endif
297 endif
298
299 call h5sclose_f(h5_space_typeid,hdf5err)
300 if(hdf5err.lt.0) then
301 Status = WRF_HDF5_ERR_CLOSE_GENERAL
302 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
303 call wrf_debug ( WARN , msg)
304 return
305 endif
306 call h5aclose_f(h5_attrid,hdf5err)
307 if(hdf5err.lt.0) then
308 Status = WRF_HDF5_ERR_ATTRIBUTE_CLOSE
309 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
310 call wrf_debug ( WARN , msg)
311 return
312 endif
313
314 Status = WRF_NO_ERR
315
316 return
317
318 end subroutine clean_phdf5_attrids
319
320
321 subroutine create_h5filetype(dtype_id,Status)
322
323 use wrf_phdf5_data
324 use ext_phdf5_support_routines
325 use hdf5
326 implicit none
327 include 'wrf_status_codes.h'
328
329 integer(hid_t),intent(out) :: dtype_id
330 integer(hid_t) :: dtstr_id
331 integer(size_t) :: type_size
332 integer(size_t) :: type_sizes
333 integer(size_t) :: type_sizei
334 integer(size_t) :: offset
335 integer, intent(out) :: Status
336 integer(hid_t) :: hdf5err
337
338 call h5tcopy_f(H5T_NATIVE_CHARACTER,dtstr_id,hdf5err)
339 if(hdf5err.lt.0) then
340 Status = WRF_HDF5_ERR_DATATYPE
341 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
342 call wrf_debug ( WARN , msg)
343 return
344 endif
345
346 type_size = 256
347
348 call h5tset_size_f(dtstr_id,type_size,hdf5err)
349 if(hdf5err.lt.0) then
350 Status = WRF_HDF5_ERR_DATATYPE
351 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
352 call wrf_debug ( WARN , msg)
353 return
354 endif
355
356 call h5tget_size_f(dtstr_id,type_sizes,hdf5err)
357 if(hdf5err.lt.0) then
358 Status = WRF_HDF5_ERR_DATATYPE
359 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
360 call wrf_debug ( WARN , msg)
361 return
362 endif
363
364 call h5tget_size_f(H5T_NATIVE_INTEGER,type_sizei,hdf5err)
365 if(hdf5err.lt.0) then
366 Status = WRF_HDF5_ERR_DATATYPE
367 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
368 call wrf_debug ( WARN , msg)
369 return
370 endif
371
372 type_size = type_sizes + 2*type_sizei
373
374 call h5tcreate_f(H5T_COMPOUND_F,type_size,dtype_id,hdf5err)
375 if(hdf5err.lt.0) then
376 Status = WRF_HDF5_ERR_DATATYPE
377 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
378 call wrf_debug ( WARN , msg)
379 return
380 endif
381
382
383 offset = 0
384
385 call h5tinsert_f(dtype_id,"dim_name",offset,dtstr_id,hdf5err)
386 if(hdf5err.lt.0) then
387 Status = WRF_HDF5_ERR_DATATYPE
388 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
389 call wrf_debug ( WARN , msg)
390 return
391 endif
392
393 offset = offset + type_sizes
394 call h5tinsert_f(dtype_id,"dim_length",offset,H5T_NATIVE_INTEGER,&
395 hdf5err)
396 if(hdf5err.lt.0) then
397 Status = WRF_HDF5_ERR_DATATYPE
398 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
399 call wrf_debug ( WARN , msg)
400 return
401 endif
402
403 offset = offset + type_sizei
404
405 call h5tinsert_f(dtype_id,"dim_unlimited",offset,H5T_NATIVE_INTEGER,&
406 hdf5err)
407 if(hdf5err.lt.0) then
408 Status = WRF_HDF5_ERR_DATATYPE
409 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
410 call wrf_debug ( WARN , msg)
411 return
412 endif
413
414
415 call h5tclose_f(dtstr_id,hdf5err)
416 if(hdf5err.lt.0) then
417 Status = WRF_HDF5_ERR_CLOSE_GENERAL
418 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
419 call wrf_debug ( WARN , msg)
420 return
421 endif
422
423 Status = WRF_NO_ERR
424 return
425 end subroutine create_h5filetype
426
427 ! check whether two types are equal, attr_type and h5_attrid
428 subroutine check_type(DataHandle,attr_type,h5_attrid,Status)
429
430 use wrf_phdf5_data
431 use ext_phdf5_support_routines
432 USE HDF5 ! This module contains all necessary modules
433 implicit none
434 include 'wrf_status_codes.h'
435
436 integer ,intent(in) :: DataHandle
437 integer(hid_t) ,intent(in) :: attr_type
438 integer(hid_t) ,intent(in) :: h5_attrid
439 integer ,intent(out) :: Status
440 integer(hid_t) :: h5_atypeid
441 integer(hid_t) :: h5_classid
442 integer(hid_t) :: h5_wrfclassid
443 logical :: flag
444 integer :: hdf5err
445 type(wrf_phdf5_data_handle),pointer :: DH
446
447 call GetDH(DataHandle,DH,Status)
448 if(Status /= WRF_NO_ERR) then
449 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
450 call wrf_debug ( WARN , msg)
451 return
452 endif
453
454 call h5aget_type_f(h5_attrid,h5_atypeid,hdf5err)
455 if(hdf5err.lt.0) then
456 Status = WRF_HDF5_ERR_DATATYPE
457 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
458 call wrf_debug ( WARN , msg)
459 return
460 endif
461
462 call h5tget_class_f(h5_atypeid,h5_classid,hdf5err)
463 if(hdf5err.lt.0) then
464 Status = WRF_HDF5_ERR_DATATYPE
465 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
466 call wrf_debug ( WARN , msg)
467 return
468 endif
469
470 call h5tget_class_f(attr_type,h5_wrfclassid,hdf5err)
471 if(hdf5err.lt.0) then
472 Status = WRF_HDF5_ERR_DATATYPE
473 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
474 call wrf_debug ( WARN , msg)
475 return
476 endif
477
478 if((h5_classid==H5T_STRING_F).AND.&
479 (attr_type==H5T_NATIVE_CHARACTER)) then
480 flag = .TRUE.
481 else
482 if(h5_classid .NE. h5_wrfclassid) then
483 flag = .FALSE.
484 else
485 flag = .TRUE.
486 endif
487 endif
488
489 if(flag .EQV. .FALSE.) then
490 Status = WRF_HDF5_ERR_TYPE_MISMATCH
491 return
492 endif
493
494 Status = WRF_NO_ERR
495 return
496 end subroutine check_type
497
498
499 subroutine retrieve_ti_info(DataHandle,h5_attrid,h5_atypeid,Count,OutCount,Status)
500
501 use wrf_phdf5_data
502 use ext_phdf5_support_routines
503 USE HDF5 ! This module contains all necessary modules
504 implicit none
505 include 'wrf_status_codes.h'
506
507 integer ,intent(in) :: DataHandle
508 integer ,intent(in) :: h5_attrid
509 integer(hid_t) ,intent(out) :: h5_atypeid
510 integer ,intent(in) :: Count
511 integer ,intent(out) :: OutCount
512 integer ,intent(out) :: Status
513 integer(hid_t) :: h5_aspaceid
514 integer :: typeclass
515 integer :: hdf5err
516 integer :: rank
517 integer(hsize_t) :: npoints
518
519 type(wrf_phdf5_data_handle),pointer :: DH
520
521
522 call GetDH(DataHandle,DH,Status)
523 if(Status /= WRF_NO_ERR) then
524 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
525 call wrf_debug ( WARN , msg)
526 return
527 endif
528
529 if(DH%FileStatus == WRF_FILE_OPENED_FOR_READ) then
530
531 call h5aget_type_f(h5_attrid,h5_atypeid,hdf5err)
532 if(hdf5err.lt.0) then
533 Status = WRF_HDF5_ERR_DATATYPE
534 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
535 call wrf_debug ( WARN , msg)
536 return
537 endif
538
539 call h5aget_space_f(h5_attrid,h5_aspaceid,hdf5err)
540 if(hdf5err.lt.0) then
541 Status = WRF_HDF5_ERR_DATASPACE
542 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
543 call wrf_debug ( WARN , msg)
544 return
545 endif
546
547 call h5sget_simple_extent_ndims_f(h5_aspaceid,rank,hdf5err)
548 if(hdf5err.lt.0) then
549 Status = WRF_HDF5_ERR_DATASPACE
550 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
551 call wrf_debug ( WARN , msg)
552 return
553 endif
554
555 if(rank > 1) then
556 ! The rank can be either 0 or 1
557 Status = WRF_HDF5_ERR_OTHERS
558 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
559 call wrf_debug ( WARN , msg)
560 return
561 endif
562
563 call h5sget_simple_extent_npoints_f(h5_aspaceid,npoints,hdf5err)
564 if(hdf5err.lt.0) then
565 Status = WRF_HDF5_ERR_DATASPACE
566 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
567 call wrf_debug ( WARN , msg)
568 return
569 endif
570
571 OutCount = npoints
572 call h5tget_class_f(h5_atypeid,typeclass,hdf5err)
573 if(hdf5err.lt.0) then
574 Status = WRF_HDF5_ERR_DATATYPE
575 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
576 call wrf_debug ( WARN , msg)
577 return
578 endif
579 if((npoints > Count).and.(typeclass.ne.H5T_STRING_F)) then
580 OutCount = Count
581 Status = WRF_HDF5_ERR_MORE_DATA_IN_FILE
582 else
583 OutCount = npoints
584 endif
585 endif
586 return
587 end subroutine retrieve_ti_info
588
589 subroutine setup_wrtd_dataset(DataHandle,DataSetName,dtypeid,countmd,&
590 dsetid,dspace_id,fspace_id,tgroupid, &
591 TimeIndex,Status)
592
593 use wrf_phdf5_data
594 use ext_phdf5_support_routines
595 USE HDF5 ! This module contains all necessary modules
596 implicit none
597 include 'wrf_status_codes.h'
598
599 integer ,intent(in) :: DataHandle
600 character*(*) ,intent(in) :: DataSetName
601 integer(hid_t) ,intent(in) :: dtypeid
602 integer ,intent(in) :: countmd
603 integer ,intent(in) :: TimeIndex
604
605 integer(hid_t) ,intent(out) :: dsetid
606 integer(hid_t) ,intent(out) :: dspace_id
607 integer(hid_t) ,intent(out) :: fspace_id
608 integer(hid_t) ,intent(out) :: tgroupid
609 integer(hid_t) :: crp_list
610 integer ,intent(out) :: Status
611
612 integer(hsize_t) ,dimension(1) :: sizes
613 integer(hsize_t) ,dimension(1) :: chunk_dims
614 integer(hsize_t) ,dimension(1) :: dims
615 integer(hsize_t) ,dimension(1) :: hdf5_maxdims
616 integer(hsize_t) ,dimension(1) :: offset
617 integer(hsize_t) ,dimension(1) :: count
618 type(wrf_phdf5_data_handle),pointer :: DH
619
620 character(Len = MaxTimeSLen) :: tname
621 character(Len = 512) :: tgroupname
622 integer :: hdf5err
623
624
625 ! get datahandle
626 call GetDH(DataHandle,DH,Status)
627 if(Status /= WRF_NO_ERR) then
628 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
629 call wrf_debug ( WARN , msg)
630 return
631 endif
632
633
634 chunk_dims(1) = countmd
635
636 dims(1) = countmd
637
638 count(1) = countmd
639
640 offset(1) = 0
641
642 sizes(1) = countmd
643
644 hdf5_maxdims(1) = countmd
645
646 ! create the memory space id
647 call h5screate_simple_f(1,dims,dspace_id,hdf5err,dims)
648 if(hdf5err.lt.0) then
649 Status = WRF_HDF5_ERR_DATASPACE
650 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
651 call wrf_debug ( WARN , msg)
652 return
653 endif
654
655 ! create file space(for parallel module, each dataset per time step)
656 call h5screate_simple_f(1,dims,fspace_id,hdf5err,hdf5_maxdims)
657 if(hdf5err.lt.0) then
658 Status = WRF_HDF5_ERR_DATASPACE
659 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
660 call wrf_debug ( WARN , msg)
661 return
662 endif
663
664 ! obtain the absolute name of the group where the dataset is located
665 call numtochar(TimeIndex,tname)
666 tgroupname = 'TIME_STAMP_'//tname
667 if(DH%TgroupIDs(TimeIndex) /= -1) then
668 tgroupid = DH%TgroupIDs(TimeIndex)
669 ! call h5gopen_f(DH%GroupID,tgroupname,tgroupid,hdf5err)
670 else
671 call h5gcreate_f(DH%GroupID,tgroupname,tgroupid,hdf5err)
672 if(hdf5err.lt.0) then
673 Status = WRF_HDF5_ERR_GROUP
674 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
675 call wrf_debug(WARN,msg)
676 return
677 endif
678 DH%TgroupIDs(TimeIndex) = tgroupid
679 endif
680
681 ! create dataset
682 call h5dcreate_f(tgroupid,DatasetName,dtypeid,fspace_id,&
683 dsetid,hdf5err)
684 if(hdf5err.lt.0) then
685 Status = WRF_HDF5_ERR_DATASET_CREATE
686 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
687 call wrf_debug ( WARN , msg)
688 return
689 endif
690
691 return
692 end subroutine setup_wrtd_dataset
693
694 subroutine extend_wrtd_dataset(DataHandle,TimeIndex,countmd,dsetid,dspaceid,&
695 fspaceid,Status)
696
697 use wrf_phdf5_data
698 use ext_phdf5_support_routines
699 USE HDF5 ! This module contains all necessary modules
700 implicit none
701 include 'wrf_status_codes.h'
702
703 integer ,intent(in) :: DataHandle
704 integer ,intent(in) :: countmd
705 integer ,intent(in) :: TimeIndex
706
707 integer(hid_t) ,intent(out) :: dsetid
708 integer(hid_t) ,intent(out) :: dspaceid
709 integer(hid_t) ,intent(out) :: fspaceid
710 integer ,intent(out) :: Status
711
712 integer(hsize_t) ,dimension(2) :: sizes
713 integer(hsize_t) ,dimension(2) :: chunk_dims
714 integer(hsize_t) ,dimension(2) :: dims
715 integer(hsize_t) ,dimension(2) :: hdf5_maxdims
716 integer(hsize_t) ,dimension(2) :: offset
717 integer(hsize_t) ,dimension(2) :: count
718
719 integer :: hdf5err
720
721 sizes(1) = countmd
722 sizes(2) = TimeIndex
723 offset(1) = 0
724 offset(2) = TimeIndex - 1
725 count(1) = countmd
726 count(2) = 1
727 dims(1) = countmd
728 dims(2) = 1
729
730 ! extend the dataset
731 CALL h5dextend_f(dsetid,sizes,hdf5err)
732 if(hdf5err.lt.0) then
733 Status = WRF_HDF5_ERR_DATASET_GENERAL
734 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
735 call wrf_debug ( WARN , msg)
736 return
737 endif
738
739 ! obtain file space id
740 CALL h5dget_space_f(dsetid,fspaceid,hdf5err)
741 if(hdf5err.lt.0) then
742 Status = WRF_HDF5_ERR_DATASPACE
743 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
744 call wrf_debug ( WARN , msg)
745 return
746 endif
747
748
749 CALL h5sselect_hyperslab_f(fspaceid, H5S_SELECT_SET_F, &
750 offset, count, hdf5err)
751 if(hdf5err.lt.0) then
752 Status = WRF_HDF5_ERR_DATASET_GENERAL
753 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
754 call wrf_debug ( WARN , msg)
755 return
756 endif
757
758 ! create the memory space id
759 call h5screate_simple_f(2,dims,dspaceid,hdf5err,dims)
760 if(hdf5err.lt.0) then
761 Status = WRF_HDF5_ERR_DATASPACE
762 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
763 call wrf_debug ( WARN , msg)
764 return
765 endif
766
767 return
768 end subroutine extend_wrtd_dataset
769
770 subroutine setup_rdtd_dataset(DataHandle,DataSetName,mtypeid,TimeIndex,&
771 countmd,outcountmd,dset_id,memspaceid, &
772 dspace_id,tgroupid,Status)
773
774 use wrf_phdf5_data
775 use ext_phdf5_support_routines
776 USE HDF5 ! This module contains all necessary modules
777 implicit none
778 include 'wrf_status_codes.h'
779
780 integer ,intent(in) :: DataHandle
781 character*(*) ,intent(in) :: DataSetName
782 integer ,intent(in) :: countmd
783 integer ,intent(out) :: outcountmd
784 integer ,intent(inout) :: mtypeid
785 integer ,intent(in) :: TimeIndex
786
787 integer(hid_t) ,intent(out) :: dset_id
788 integer(hid_t) ,intent(out) :: dspace_id
789 integer(hid_t) ,intent(out) :: memspaceid
790 integer(hid_t) ,intent(out) :: tgroupid
791 integer ,intent(out) :: Status
792
793 integer(hid_t) :: dtype_id
794 integer(hid_t) :: class_type
795 integer(hsize_t) ,dimension(1) :: sizes
796 integer(hsize_t) ,dimension(1) :: dims
797 integer(hsize_t) ,dimension(1) :: h5_dims
798 integer(hsize_t) ,dimension(1) :: hdf5_maxdims
799 integer(hsize_t) ,dimension(1) :: offset
800 integer(hsize_t) ,dimension(1) :: count
801 integer :: StoredDim
802 type(wrf_phdf5_data_handle),pointer :: DH
803
804 logical :: flag
805 integer :: hdf5err
806
807 character(Len = MaxTimeSLen) :: tname
808 character(Len = 512) :: tgroupname
809 ! get datahandle
810 call GetDH(DataHandle,DH,Status)
811 if(Status /= WRF_NO_ERR) then
812 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
813 call wrf_debug ( WARN , msg)
814 return
815 endif
816
817 ! obtain the absolute name of the group where the dataset is located
818 call numtochar(TimeIndex,tname)
819 tgroupname = 'TIME_STAMP_'//tname
820 call h5gopen_f(DH%GroupID,tgroupname,tgroupid,hdf5err)
821
822 ! Obtain HDF5 dataset id
823 call h5dopen_f(tgroupid,DataSetName,dset_id,hdf5err)
824 if(hdf5err.lt.0) then
825 Status = WRF_HDF5_ERR_DATASET_OPEN
826 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
827 call wrf_debug ( WARN , msg)
828 return
829 endif
830
831 ! Obtain the datatype
832 call h5dget_type_f(dset_id,dtype_id,hdf5err)
833 if(hdf5err.lt.0) then
834 Status = WRF_HDF5_ERR_DATASET_GENERAL
835 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
836 call wrf_debug ( WARN , msg)
837 return
838 endif
839
840 call h5tget_class_f(dtype_id,class_type,hdf5err)
841 if(hdf5err.lt.0) then
842 Status = WRF_HDF5_ERR_DATATYPE
843 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
844 call wrf_debug ( WARN , msg)
845 return
846 endif
847
848
849 if(mtypeid == H5T_NATIVE_REAL .or. mtypeid == H5T_NATIVE_DOUBLE) then
850 if( class_type /= H5T_FLOAT_F) then
851 Status = WRF_HDF5_ERR_TYPE_MISMATCH
852 write(msg,*) 'Warning TYPE MISMATCH in ',__FILE__,', line', __LINE__
853 call wrf_debug ( WARN , msg)
854 return
855 endif
856 else if(mtypeid ==H5T_NATIVE_CHARACTER) then
857 if(class_type /= H5T_STRING_F) then
858 Status = WRF_HDF5_ERR_TYPE_MISMATCH
859 write(msg,*) 'Warning TYPE MISMATCH in ',__FILE__,', line', __LINE__
860 call wrf_debug ( WARN , msg)
861 return
862 endif
863 else if(mtypeid == H5T_NATIVE_INTEGER) then
864 if(class_type /= H5T_INTEGER_F) then
865 Status = WRF_HDF5_ERR_TYPE_MISMATCH
866 write(msg,*) 'Warning TYPE MISMATCH in ',__FILE__,', line', __LINE__
867 call wrf_debug ( WARN , msg)
868 return
869 endif
870 else if(mtypeid == DH%EnumID) then
871 if(class_type /= H5T_ENUM_F) then
872 Status = WRF_HDF5_ERR_TYPE_MISMATCH
873 write(msg,*) 'Warning TYPE MISMATCH in ',__FILE__,', line', __LINE__
874 call wrf_debug ( WARN , msg)
875 return
876 endif
877 call h5tequal_f(dtype_id,DH%EnumID,flag,hdf5err)
878 if(hdf5err.lt.0) then
879 Status = WRF_HDF5_ERR_DATASET_OPEN
880 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
881 call wrf_debug ( WARN , msg)
882 return
883 endif
884 if(flag .EQV. .FALSE.) then
885 Status = WRF_HDF5_ERR_TYPE_MISMATCH
886 write(msg,*) 'Warning TYPE MISMATCH in ',__FILE__,', line', __LINE__
887 call wrf_debug ( WARN , msg)
888 return
889 endif
890 else
891 Status = WRF_HDF5_ERR_BAD_DATA_TYPE
892 write(msg,*)'Fatal Non-WRF supported TYPE in ',__FILE__,', line',__LINE__
893 call wrf_debug(FATAL, msg)
894 return
895 endif
896 ! update string id
897 if(mtypeid == H5T_NATIVE_CHARACTER) then
898 mtypeid = dtype_id
899 endif
900
901 ! Obtain the dataspace
902 call h5dget_space_f(dset_id,dspace_id,hdf5err)
903 if(hdf5err.lt.0) then
904 Status = WRF_HDF5_ERR_DATASPACE
905 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
906 call wrf_debug ( WARN , msg)
907 return
908 endif
909
910 ! Obtain the rank of the dimension
911 call h5sget_simple_extent_ndims_f(dspace_id,StoredDim,hdf5err)
912 if(hdf5err.lt.0) then
913 Status = WRF_HDF5_ERR_DATASPACE
914 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
915 call wrf_debug ( WARN , msg)
916 return
917 endif
918
919
920 if(StoredDim /=1) then
921 Status = WRF_HDF5_ERR_DATASPACE
922 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
923 call wrf_debug ( WARN , msg)
924 return
925 endif
926
927
928 call h5sget_simple_extent_dims_f(dspace_id,h5_dims,hdf5_maxdims,hdf5err)
929 if(hdf5err.lt.0) then
930 Status = WRF_HDF5_ERR_DATASPACE
931 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
932 call wrf_debug ( WARN , msg)
933 return
934 endif
935
936
937 if(countmd <= 0) then
938 Status = WRF_HDF5_ERR_ZERO_LENGTH_READ
939 write(msg,*) 'Warning ZERO LENGTH READ in ',__FILE__,', line', __LINE__
940 call wrf_debug ( WARN , msg)
941 return
942 endif
943
944 if(countmd .lt. h5_dims(1)) then
945 outcountmd = countmd
946 else
947 outcountmd = h5_dims(1)
948 endif
949
950 ! create memspace_id
951 dims(1) = outcountmd
952
953 call h5screate_simple_f(1,dims,memspaceid,hdf5err)
954 if(hdf5err.lt.0) then
955 Status = WRF_HDF5_ERR_DATASPACE
956 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
957 call wrf_debug ( WARN , msg)
958 return
959 endif
960
961 offset(1) = 0
962 count(1) = outcountmd
963
964 return
965 end subroutine setup_rdtd_dataset
966
967 subroutine make_strid(str_len,str_id,Status)
968
969 use wrf_phdf5_data
970 USE HDF5 ! This module contains all necessary modules
971 implicit none
972 include 'wrf_status_codes.h'
973
974 integer ,intent(in) :: str_len
975 integer(hid_t),intent(out) :: str_id
976 integer ,intent(out) :: Status
977 integer(size_t) :: str_lensize
978 integer :: hdf5err
979
980 Status = WRF_NO_ERR
981 if(str_len <= 0) then
982 Status = WRF_HDF5_ERR_ATTRIBUTE_GENERAL
983 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
984 call wrf_debug ( WARN , msg)
985 return
986 endif
987
988 call h5tcopy_f(H5T_NATIVE_CHARACTER,str_id,hdf5err)
989 if(hdf5err.lt.0) then
990 Status = WRF_HDF5_ERR_DATATYPE
991 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
992 call wrf_debug ( WARN , msg)
993 return
994 endif
995
996 str_lensize = str_len
997 call h5tset_size_f(str_id,str_lensize,hdf5err)
998 if(hdf5err.lt.0) then
999 Status = WRF_HDF5_ERR_DATATYPE
1000 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
1001 call wrf_debug ( WARN , msg)
1002 return
1003 endif
1004
1005 call h5tset_strpad_f(str_id,H5T_STR_SPACEPAD_F,hdf5err)
1006 if(hdf5err.lt.0) then
1007 Status = WRF_HDF5_ERR_DATATYPE
1008 write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
1009 call wrf_debug ( WARN , msg)
1010 return
1011 endif
1012
1013 end subroutine make_strid