dataset_common_hdf5.F90       coverage:  92.31 %func     82.58 %block


     1) module Dataset_Common_HDF5_class
     2)  
     3)   use Dataset_Base_class
     4)   
     5)   use PFLOTRAN_Constants_module
     6) 
     7)   implicit none
     8) 
     9)   private
    10) 
    11) #include "petsc/finclude/petscsys.h"
    12) 
    13)   type, public, extends(dataset_base_type) :: dataset_common_hdf5_type
    14)     character(len=MAXWORDLENGTH) :: hdf5_dataset_name
    15)     PetscBool :: realization_dependent
    16)     PetscInt :: max_buffer_size
    17)     PetscBool :: is_cell_indexed
    18)   end type dataset_common_hdf5_type
    19) 
    20)   public :: DatasetCommonHDF5Create, &
    21)             DatasetCommonHDF5Init, &
    22)             DatasetCommonHDF5Copy, &
    23)             DatasetCommonHDF5Cast, &
    24)             DatasetCommonHDF5Read, &
    25)             DatasetCommonHDF5ReadSelectCase, &
    26)             DatasetCommonHDF5Load, &
    27)             DatasetCommonHDF5IsCellIndexed, &
    28)             DatasetCommonHDF5Print, &
    29)             DatasetCommonHDF5Strip, &
    30)             DatasetCommonHDF5Destroy
    31)   
    32) #if defined(PETSC_HAVE_HDF5)   
    33)   public :: DatasetCommonHDF5ReadTimes
    34) #endif  
    35) 
    36) contains
    37) 
    38) ! ************************************************************************** !
    39) 
    40) function DatasetCommonHDF5Create()
    41)   ! 
    42)   ! Creates members of common hdf5 database class
    43)   ! 
    44)   ! Author: Glenn Hammond
    45)   ! Date: 05/03/13
    46)   ! 
    47)   
    48)   implicit none
    49)   
    50)   class(dataset_common_hdf5_type), pointer :: dataset
    51) 
    52)   class(dataset_common_hdf5_type), pointer :: DatasetCommonHDF5Create
    53)   
    54)   allocate(dataset)
    55)   call DatasetCommonHDF5Init(dataset)
    56) 
    57)   DatasetCommonHDF5Create => dataset
    58)     
    59) end function DatasetCommonHDF5Create
    60) 
    61) ! ************************************************************************** !
    62) 
    63) subroutine DatasetCommonHDF5Init(this)
    64)   ! 
    65)   ! Initializes members of common hdf5 dataset class
    66)   ! 
    67)   ! Author: Glenn Hammond
    68)   ! Date: 05/03/13
    69)   ! 
    70)   
    71)   implicit none
    72)   
    73)   class(dataset_common_hdf5_type) :: this
    74)   
    75)   call DatasetBaseInit(this)
    76)   this%hdf5_dataset_name = ''
    77)   this%realization_dependent = PETSC_FALSE
    78)   this%max_buffer_size = UNINITIALIZED_INTEGER
    79)   this%is_cell_indexed = PETSC_FALSE
    80)   this%data_type = DATASET_REAL
    81)     
    82) end subroutine DatasetCommonHDF5Init
    83) 
    84) ! ************************************************************************** !
    85) 
    86) subroutine DatasetCommonHDF5Copy(this, that)
    87)   ! 
    88)   ! Copies members of common hdf5 dataset class
    89)   ! 
    90)   ! Author: Glenn Hammond
    91)   ! Date: 05/03/13
    92)   ! 
    93)   
    94)   implicit none
    95)   
    96)   class(dataset_common_hdf5_type) :: this
    97)   class(dataset_common_hdf5_type) :: that
    98)   
    99)   call DatasetBaseCopy(this,that)
   100)   that%hdf5_dataset_name = this%hdf5_dataset_name
   101)   that%realization_dependent = this%realization_dependent
   102)   that%max_buffer_size = this%max_buffer_size
   103)   that%is_cell_indexed = this%is_cell_indexed
   104)     
   105) end subroutine DatasetCommonHDF5Copy
   106) 
   107) ! ************************************************************************** !
   108) 
   109) function DatasetCommonHDF5Cast(this)
   110)   ! 
   111)   ! Casts a dataset_base_type to dataset_common_hdf5_type
   112)   ! 
   113)   ! Author: Glenn Hammond
   114)   ! Date: 05/03/13
   115)   ! 
   116)   
   117)   implicit none
   118) 
   119)   class(dataset_base_type), pointer :: this
   120) 
   121)   class(dataset_common_hdf5_type), pointer :: DatasetCommonHDF5Cast
   122)   
   123)   nullify(DatasetCommonHDF5Cast)
   124)   select type (this)
   125)     class is (dataset_common_hdf5_type)
   126)       DatasetCommonHDF5Cast => this
   127)   end select
   128)     
   129) end function DatasetCommonHDF5Cast
   130) 
   131) ! ************************************************************************** !
   132) 
   133) subroutine DatasetCommonHDF5Read(this,input,option)
   134)   ! 
   135)   ! Reads in contents of a dataset card
   136)   ! 
   137)   ! Author: Glenn Hammond
   138)   ! Date: 01/12/11, 06/04/13
   139)   ! 
   140) 
   141)   use Option_module
   142)   use Input_Aux_module
   143)   use String_module
   144) 
   145)   implicit none
   146)   
   147)   class(dataset_common_hdf5_type) :: this
   148)   type(input_type), pointer :: input
   149)   type(option_type) :: option
   150)   
   151)   character(len=MAXWORDLENGTH) :: keyword
   152)   PetscBool :: found
   153) 
   154)   input%ierr = 0
   155)   do
   156)   
   157)     call InputReadPflotranString(input,option)
   158) 
   159)     if (InputCheckExit(input,option)) exit  
   160) 
   161)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   162)     call InputErrorMsg(input,option,'keyword','DATASET')
   163)     call StringToUpper(keyword)   
   164)       
   165)     call DatasetCommonHDF5ReadSelectCase(this,input,keyword,found,option)
   166) 
   167)     if (.not.found) then
   168)       call InputKeywordUnrecognized(keyword,'dataset',option)
   169)     endif
   170)   
   171)   enddo
   172)   
   173)   if (len_trim(this%hdf5_dataset_name) < 1) then
   174)     this%hdf5_dataset_name = this%name
   175)   endif
   176)   
   177) end subroutine DatasetCommonHDF5Read
   178) 
   179) ! ************************************************************************** !
   180) 
   181) subroutine DatasetCommonHDF5ReadSelectCase(this,input,keyword,found,option)
   182)   ! 
   183)   ! Compares keyword against HDF5 common
   184)   ! keywords
   185)   ! 
   186)   ! Author: Glenn Hammond
   187)   ! Date: 06/04/13
   188)   ! 
   189) 
   190)   use Option_module
   191)   use Input_Aux_module
   192)   use String_module
   193) 
   194)   implicit none
   195)   
   196)   class(dataset_common_hdf5_type) :: this
   197)   type(input_type) :: input
   198)   character(len=MAXWORDLENGTH) :: keyword
   199)   PetscBool :: found
   200)   type(option_type) :: option
   201) 
   202)   found = PETSC_TRUE
   203)   select case(trim(keyword))
   204)     case('NAME') 
   205)       call InputReadWord(input,option,this%name,PETSC_TRUE)
   206)       call InputErrorMsg(input,option,'name','DATASET')
   207)     case('HDF5_DATASET_NAME') 
   208)       call InputReadWord(input,option,this%hdf5_dataset_name,PETSC_TRUE)
   209)       call InputErrorMsg(input,option,'hdf5_dataset_name','DATASET')
   210)     case('FILENAME') 
   211)       call InputReadNChars(input,option,this%filename, &
   212)                             MAXSTRINGLENGTH,PETSC_TRUE)
   213)       call InputErrorMsg(input,option,'name','DATASET')
   214)     case('REALIZATION_DEPENDENT')
   215)       this%realization_dependent = PETSC_TRUE
   216)     case('MAX_BUFFER_SIZE') 
   217)       call InputReadInt(input,option,this%max_buffer_size)
   218)       call InputErrorMsg(input,option,'max_buffer_size','DATASET')
   219)     case default
   220)       found = PETSC_FALSE
   221)   end select  
   222)   
   223) end subroutine DatasetCommonHDF5ReadSelectCase
   224) 
   225) #if defined(PETSC_HAVE_HDF5)
   226) 
   227) ! ************************************************************************** !
   228) 
   229) subroutine DatasetCommonHDF5ReadTimes(filename,dataset_name,time_storage, &
   230)                                       option)
   231)   ! 
   232)   ! DatasetGlobalReadTimes: Read dataset times into time storage
   233)   ! 
   234)   ! Author: Glenn Hammond
   235)   ! Date: 01/12/08
   236)   ! 
   237)                          
   238)   use hdf5
   239)   use Time_Storage_module
   240)   use Units_module, only : UnitsConvertToInternal
   241)   use Option_module
   242)   use Logging_module
   243)   use HDF5_Aux_module
   244)   
   245)   implicit none
   246) 
   247)   character(len=MAXSTRINGLENGTH) :: filename
   248)   character(len=MAXWORDLENGTH) :: dataset_name
   249)   type(time_storage_type), pointer :: time_storage
   250)   type(option_type) :: option
   251)   
   252)   integer(HID_T) :: file_id
   253)   integer(HID_T) :: file_space_id
   254)   integer(HID_T) :: memory_space_id
   255)   integer(HID_T) :: dataset_id
   256)   integer(HID_T) :: prop_id
   257)   integer(HID_T) :: grp_id
   258)   integer(HID_T) :: atype_id  
   259)   integer(HID_T) :: attribute_id  
   260)   integer(HSIZE_T) :: num_times  
   261)   integer(HSIZE_T) :: length(1)
   262)   integer(HSIZE_T) :: attribute_dim(1)
   263)   integer(SIZE_T) :: size_t_int  
   264)   PetscMPIInt :: array_rank_mpi
   265)   character(len=MAXSTRINGLENGTH) :: string
   266)   character(len=MAXWORDLENGTH) :: attribute_name, time_units
   267)   character(len=MAXWORDLENGTH) :: internal_units
   268)   PetscMPIInt :: int_mpi
   269)   PetscInt :: temp_int, num_times_read_by_iorank
   270)   PetscMPIInt :: hdf5_err, h5fopen_err
   271)   PetscBool :: attribute_exists, group_exists
   272)   PetscLogDouble :: tstart, tend  
   273)   PetscErrorCode :: ierr
   274)   
   275)   call PetscLogEventBegin(logging%event_read_array_hdf5,ierr);CHKERRQ(ierr)
   276)   
   277) !#define TIME_READING_TIMES
   278) #ifdef TIME_READING_TIMES
   279)   call PetscTime(tstart,ierr);CHKERRQ(ierr)
   280) #endif  
   281)   
   282)   h5fopen_err = 0
   283)   if (option%myrank == option%io_rank) then
   284)     ! open the file
   285)     call h5open_f(hdf5_err)
   286)     option%io_buffer = 'Opening hdf5 file: ' // trim(filename)
   287)     call printMsg(option)
   288)   
   289)     ! set read file access property
   290)     call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   291)     !geh: DO NOT call h5pset_fapl_mpio_f here since the attribute is only being
   292)     !     read by the io_rank process.
   293) !#ifndef SERIAL_HDF5
   294) !    call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   295) !#endif
   296)     call HDF5OpenFileReadOnly(filename,file_id,prop_id,option)
   297)     h5fopen_err = hdf5_err
   298)     call h5pclose_f(prop_id,hdf5_err)
   299) 
   300)     if (h5fopen_err == 0) then
   301)       option%io_buffer = 'Opening hdf5 group: ' // trim(dataset_name)
   302)       call printMsg(option)
   303)       call h5gopen_f(file_id,dataset_name,grp_id,hdf5_err)
   304) 
   305)       attribute_name = "Time Units"
   306)       call H5aexists_f(grp_id,attribute_name,attribute_exists,hdf5_err)
   307)       if (attribute_exists) then
   308)         attribute_dim = 1
   309)         call h5tcopy_f(H5T_NATIVE_CHARACTER,atype_id,hdf5_err)
   310)         size_t_int = MAXWORDLENGTH
   311)         call h5tset_size_f(atype_id,size_t_int,hdf5_err)
   312)         call h5aopen_f(grp_id,attribute_name,attribute_id,hdf5_err)
   313)         call h5aread_f(attribute_id,atype_id,time_units,attribute_dim,hdf5_err)
   314)         call h5aclose_f(attribute_id,hdf5_err)
   315)       else
   316)         option%io_buffer = 'Time Units assumed to be seconds.'
   317)         call printWrnMsg(option)
   318)         time_units = 's'
   319)       endif
   320) 
   321)       ! Check whether a time array actually exists
   322)       string = 'Times'
   323)       call h5lexists_f(grp_id,string,group_exists,hdf5_err)
   324) 
   325)       if (group_exists) then
   326)         !geh: Should check to see if "Times" dataset exists.
   327)         option%io_buffer = 'Opening data set: ' // trim(string)
   328)         call printMsg(option)  
   329)         call h5dopen_f(grp_id,string,dataset_id,hdf5_err)
   330)         call h5dget_space_f(dataset_id,file_space_id,hdf5_err)
   331)         call h5sget_simple_extent_npoints_f(file_space_id,num_times,hdf5_err)
   332)         num_times_read_by_iorank = num_times
   333)       else
   334)         num_times_read_by_iorank = -1
   335)       endif
   336)     endif ! h5fopen_err == 0
   337)   endif
   338) 
   339)   ! Need to catch errors in opening the file.  Since the file is only opened
   340)   ! by the I/O rank, need to broadcast the error flag.
   341)   temp_int = h5fopen_err
   342)   int_mpi = 1
   343)   call MPI_Bcast(temp_int,int_mpi,MPI_INTEGER,option%io_rank, &
   344)                  option%mycomm,ierr)
   345)   if (temp_int < 0) then ! actually h5fopen_err
   346)     option%io_buffer = 'Error opening file: ' // trim(filename)
   347)     call printErrMsg(option)
   348)   endif
   349)   
   350)   int_mpi = 1
   351)   call MPI_Bcast(num_times_read_by_iorank,int_mpi,MPI_INTEGER, &
   352)                  option%io_rank,option%mycomm,ierr)
   353)   num_times = num_times_read_by_iorank
   354)   
   355)   if (num_times == -1) then
   356)     ! no times exist, simply return
   357)     return
   358)   endif
   359)   
   360)   time_storage => TimeStorageCreate()
   361)   time_storage%max_time_index = num_times
   362)   allocate(time_storage%times(num_times))
   363)   time_storage%times = 0.d0
   364)   
   365)   if (option%myrank == option%io_rank) then 
   366)     array_rank_mpi = 1
   367)     length(1) = num_times
   368)     call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
   369) #ifndef SERIAL_HDF5
   370)     call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F,hdf5_err)
   371) #endif  
   372)     call h5screate_simple_f(array_rank_mpi,length,memory_space_id,hdf5_err,length)    
   373)     call h5dread_f(dataset_id,H5T_NATIVE_DOUBLE,time_storage%times, &
   374)                     length,hdf5_err,memory_space_id,file_space_id,prop_id)
   375) 
   376)       
   377)     call h5pclose_f(prop_id,hdf5_err)
   378)     if (memory_space_id > -1) call h5sclose_f(memory_space_id,hdf5_err)
   379)     call h5sclose_f(file_space_id,hdf5_err)
   380)     call h5dclose_f(dataset_id,hdf5_err)  
   381)     option%io_buffer = 'Closing group: ' // trim(dataset_name)
   382)     call printMsg(option)  
   383)     call h5gclose_f(grp_id,hdf5_err)  
   384)     option%io_buffer = 'Closing hdf5 file: ' // trim(filename)
   385)     call printMsg(option)  
   386)     call h5fclose_f(file_id,hdf5_err)
   387)     call h5close_f(hdf5_err) 
   388)     internal_units = 'sec'
   389)     time_storage%times = time_storage%times * &
   390)       UnitsConvertToInternal(time_units,internal_units,option)
   391)   endif
   392) 
   393)   int_mpi = num_times
   394)   call MPI_Bcast(time_storage%times,int_mpi,MPI_DOUBLE_PRECISION, &
   395)                  option%io_rank,option%mycomm,ierr)
   396) 
   397) #ifdef TIME_READING_TIMES
   398)   call MPI_Barrier(option%mycomm,ierr)
   399)   call PetscTime(tend,ierr);CHKERRQ(ierr)
   400)   write(option%io_buffer,'(f6.2," Seconds to read dataset times",a,".")') &
   401)     tend-tstart, trim(dataset_name) // ' (' // trim(option%group_prefix) // &
   402)     ')'
   403)   if (option%myrank == option%io_rank) then
   404)     print *, trim(option%io_buffer)
   405)   endif
   406) #endif
   407) 
   408)   call PetscLogEventEnd(logging%event_read_array_hdf5,ierr);CHKERRQ(ierr)
   409)   
   410) end subroutine DatasetCommonHDF5ReadTimes
   411) #endif
   412) 
   413) ! ************************************************************************** !
   414) 
   415) function DatasetCommonHDF5Load(this,option)
   416)   ! 
   417)   ! Updates indices and returns whether to load new data.
   418)   ! 
   419)   ! Author: Glenn Hammond
   420)   ! Date: 05/03/13
   421)   ! 
   422)   
   423) #if defined(PETSC_HAVE_HDF5)    
   424)   use hdf5, only : H5T_NATIVE_DOUBLE
   425) #endif
   426)   use Option_module
   427)   use Time_Storage_module
   428) 
   429)   implicit none
   430)   
   431)   PetscBool :: DatasetCommonHDF5Load
   432) 
   433)   class(dataset_common_hdf5_type) :: this
   434)   type(option_type) :: option
   435)   
   436)   PetscBool :: read_due_to_time
   437)   
   438)   DatasetCommonHDF5Load = PETSC_FALSE
   439)   
   440)   if (.not.associated(this%time_storage)) then
   441) #if defined(PETSC_HAVE_HDF5)    
   442)     call DatasetCommonHDF5ReadTimes(this%filename,this%hdf5_dataset_name, &
   443)                                     this%time_storage,option)
   444) #endif
   445)     ! if no times are read, this%time_storage will be null coming out of
   446)     ! DatasetCommonHDF5ReadTimes()
   447)   endif
   448)   
   449)   if (associated(this%time_storage)) then
   450)     this%time_storage%cur_time = option%time
   451)     ! sets correct cur_time_index
   452)     call TimeStorageUpdate(this%time_storage)
   453)     read_due_to_time = &
   454)       (this%time_storage%cur_time_index > 0 & 
   455)        .and. &
   456)        this%time_storage%cur_time_index >= &
   457)         ! both of the below will be zero initially
   458)          this%buffer_slice_offset + this%buffer_nslice &
   459)        .and. &
   460)        this%time_storage%cur_time_index < &
   461)          this%time_storage%max_time_index)
   462)   else
   463)     read_due_to_time = PETSC_FALSE
   464)   endif
   465)   
   466)   if (read_due_to_time .or. &
   467)        ! essentially gets the data set read if only one time slice
   468)       .not.associated(this%rarray)) then
   469)     if (associated(this%time_storage)) then
   470)       if (this%time_storage%cur_time_index > 0) then
   471)         this%buffer_slice_offset = this%time_storage%cur_time_index - 1
   472)       endif
   473)     endif
   474)     DatasetCommonHDF5Load = PETSC_TRUE
   475)   endif
   476)     
   477) end function DatasetCommonHDF5Load
   478) 
   479) ! ************************************************************************** !
   480) 
   481) function DatasetCommonHDF5IsCellIndexed(dataset,option)
   482)   ! 
   483)   ! Determine whether a dataset is indexed by
   484)   ! cell ids
   485)   ! 
   486)   ! Author: Glenn Hammond
   487)   ! Date: 05/03/13
   488)   ! 
   489) 
   490)   use Option_module
   491)   use HDF5_Aux_module
   492) 
   493)   implicit none
   494)   
   495)   class(dataset_common_hdf5_type) :: dataset
   496)   type(option_type) :: option
   497)   
   498)   PetscBool :: DatasetCommonHDF5IsCellIndexed
   499)   
   500) #if defined(PETSC_HAVE_HDF5)
   501) 
   502)   DatasetCommonHDF5IsCellIndexed = &
   503)     .not.HDF5GroupExists(dataset%filename,dataset%hdf5_dataset_name,option)
   504) 
   505) #endif
   506)  
   507) end function DatasetCommonHDF5IsCellIndexed
   508) 
   509) ! ************************************************************************** !
   510) 
   511) function DatasetCommonHDF5GetPointer(dataset_list, dataset_name, &
   512)                                      debug_string, option)
   513)   ! 
   514)   ! Returns the pointer to the dataset named "name"
   515)   ! 
   516)   ! Author: Glenn Hammond
   517)   ! Date: 05/03/13
   518)   ! 
   519) 
   520)   use Option_module
   521)   use String_module
   522)   
   523)   class(dataset_base_type), pointer :: dataset_list
   524)   character(len=MAXWORDLENGTH) :: dataset_name
   525)   character(len=MAXSTRINGLENGTH) :: debug_string
   526)   type(option_type) :: option
   527) 
   528)   class(dataset_common_hdf5_type), pointer :: DatasetCommonHDF5GetPointer
   529)   
   530)   class(dataset_base_type), pointer :: dataset
   531) 
   532)   nullify(DatasetCommonHDF5GetPointer)
   533)   dataset => DatasetBaseGetPointer(dataset_list, dataset_name, &
   534)                                    debug_string, option)
   535)   select type(dataset)
   536)     class is (dataset_common_hdf5_type)
   537)       DatasetCommonHDF5GetPointer => dataset
   538)     class default
   539)       option%io_buffer = 'Dataset "' // trim(dataset_name) // '" in "' // &
   540)              trim(debug_string) // '" not of type Common HDF5.'
   541)       call printErrMsg(option)    
   542)   end select
   543) 
   544) end function DatasetCommonHDF5GetPointer
   545) 
   546) ! ************************************************************************** !
   547) 
   548) subroutine DatasetCommonHDF5Print(this,option)
   549)   ! 
   550)   ! Prints dataset info
   551)   ! 
   552)   ! Author: Glenn Hammond
   553)   ! Date: 10/22/13
   554)   ! 
   555) 
   556)   use Option_module
   557) 
   558)   implicit none
   559)   
   560)   class(dataset_common_hdf5_type) :: this
   561)   type(option_type) :: option
   562)     
   563)   if (len_trim(this%hdf5_dataset_name) > 0) then
   564)     write(option%fid_out,'(10x,''HDF5 Dataset Name: '',a)') &
   565)       trim(this%hdf5_dataset_name)
   566)   endif
   567)   if (this%realization_dependent) then
   568)     write(option%fid_out,'(10x,''Realization Dependent: yes'')')
   569)   else
   570)     write(option%fid_out,'(10x,''Realization Dependent: no'')') 
   571)   endif
   572)   if (this%is_cell_indexed) then
   573)     write(option%fid_out,'(10x,''Cell Indexed: yes'')')
   574)   else
   575)     write(option%fid_out,'(10x,''Cell Indexed: no'')')
   576)   endif
   577)   write(option%fid_out,'(10x,''Maximum Buffer Size: '',i3)') &
   578)     this%max_buffer_size
   579)   
   580) end subroutine DatasetCommonHDF5Print
   581) 
   582) ! ************************************************************************** !
   583) 
   584) subroutine DatasetCommonHDF5Strip(this)
   585)   ! 
   586)   ! Strips allocated objects within common hdf5 dataset
   587)   ! object
   588)   ! 
   589)   ! Author: Glenn Hammond
   590)   ! Date: 05/03/13
   591)   ! 
   592) 
   593)   implicit none
   594)   
   595)   class(dataset_common_hdf5_type) :: this
   596)   
   597)   call DatasetBaseStrip(this)
   598)   
   599) end subroutine DatasetCommonHDF5Strip
   600) 
   601) ! ************************************************************************** !
   602) 
   603) subroutine DatasetCommonHDF5Destroy(this)
   604)   ! 
   605)   ! Destroys a dataset
   606)   ! 
   607)   ! Author: Glenn Hammond
   608)   ! Date: 01/12/11
   609)   ! 
   610) 
   611)   implicit none
   612)   
   613)   class(dataset_common_hdf5_type), pointer :: this
   614)   
   615)   if (.not.associated(this)) return
   616)   
   617)   call DatasetCommonHDF5Strip(this)
   618)   
   619)   deallocate(this)
   620)   nullify(this)
   621)   
   622) end subroutine DatasetCommonHDF5Destroy
   623) 
   624) end module Dataset_Common_HDF5_class

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody