source: lmdz_wrf/WRFV3/external/io_phdf5/wrf-phdf5support.F90 @ 1

Last change on this file since 1 was 1, checked in by lfita, 10 years ago
  • -- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

File size: 46.1 KB
Line 
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
25module wrf_phdf5_data
26
27  use HDF5
28  integer                , parameter      :: FATAL            = 1
29  integer                , parameter      :: WARN             = 1
30  integer                , parameter      :: WrfDataHandleMax = 99
31  integer                , parameter      :: MaxDims          = 2000 ! = NF_MAX_VARS
32  integer                , parameter      :: MaxTabDims       = 100  ! temporary,changable
33  integer                , parameter      :: MaxVars          = 2000
34  integer                , parameter      :: MaxTimes         = 9999 ! temporary, changable
35  integer                , parameter      :: MaxTimeSLen      = 6    ! not exceed 1,000,000 timestamp
36  integer                , parameter      :: DateStrLen       = 19
37  integer                , parameter      :: VarNameLen       = 31
38  integer                , parameter      :: NO_DIM           = 0
39  integer                , parameter      :: NVarDims         = 4
40  integer                , parameter      :: NMDVarDims       = 2
41  integer                , parameter      :: CompDsetSize     = 64256 ! set to 63K
42  character (8)          , parameter      :: NO_NAME          = 'NULL'
43  character(4)           , parameter      :: hdf5_true        ='TRUE'
44  character(5)           , parameter      :: hdf5_false       ='FALSE'
45  integer                , parameter      :: MemOrdLen        = 3
46  character (DateStrLen) , parameter      :: ZeroDate = '0000-00-00-00:00:00'
47
48#include "wrf_io_flags.h"
49! This is a hack.  WRF IOAPI no longer supports WRF_CHARACTER.  Rip this out! 
50  integer, parameter  :: WRF_CHARACTER                        = 1080
51
52  character (120)                         :: msg
53
54  ! derived data type for dimensional table
55  type :: dim_scale
56     character (len = 256) :: dim_name
57     integer              :: length
58     integer              :: unlimited
59  end type dim_scale
60
61  type :: wrf_phdf5_data_handle
62     character (256)                        :: FileName
63     integer                               :: FileStatus
64     integer                               :: Comm
65     integer(hid_t)                        :: FileID
66     integer(hid_t)                        :: GroupID
67     integer(hid_t)                        :: DimGroupID
68     integer(hid_t)                        :: EnumID
69     character (256)                       :: GroupName
70     character (256)                       :: DimGroupName
71     logical                               :: Free
72     logical                               :: Write
73     character (5)                         :: TimesName
74     integer                               :: TimeIndex
75     integer                               :: MaxTimeCount
76     integer                               :: CurrentTime  !Only used for read
77     integer                               :: NumberTimes  !Only used for read
78     character (DateStrLen), pointer       :: Times(:)
79     integer(hid_t)                        :: TimesID
80     integer(hid_t)                        :: str_id
81     integer               , pointer       :: DimLengths(:)
82     integer               , pointer       :: DimIDs(:)
83     character (31)        , pointer       :: DimNames(:)
84     integer                               :: DimUnlimID
85     character (9)                         :: DimUnlimName
86     type (dim_scale)      , pointer       :: DIMTABLE(:)
87     integer       , dimension(NVarDims)   :: DimID
88     integer       , dimension(NVarDims)   :: Dimension
89     !     integer               , pointer       :: MDDsetIDs(:)
90     integer               , pointer       :: MDVarDimLens(:)
91     character (256)        , pointer       :: MDVarNames(:)
92     integer(hid_t)        , pointer       :: TgroupIDs(:)
93     integer(hid_t)        , pointer       :: DsetIDs(:)
94     integer(hid_t)        , pointer       :: MDDsetIDs(:)
95     !     integer(hid_t)                        :: DimTableID
96     integer               , pointer       :: VarDimLens(:,:)
97     character (VarNameLen), pointer       :: VarNames(:)
98     integer                               :: CurrentVariable  !Only used for read
99     integer                               :: NumVars
100! first_operation is set to .TRUE. when a new handle is allocated
101! or when open-for-write or open-for-read are committed.  It is set
102! to .FALSE. when the first field is read or written.
103     logical                               :: first_operation
104  end type wrf_phdf5_data_handle
105  type(wrf_phdf5_data_handle),target        :: WrfDataHandles(WrfDataHandleMax)
106
107end module wrf_phdf5_data
108
109
110module ext_phdf5_support_routines
111
112  implicit none
113
114CONTAINS
115
116  subroutine allocHandle(DataHandle,DH,Comm,Status)
117
118    use wrf_phdf5_data
119    use HDF5
120    include 'wrf_status_codes.h'
121
122    integer              ,intent(out) :: DataHandle
123    type(wrf_phdf5_data_handle),pointer:: DH
124    integer              ,intent(IN)  :: Comm
125    integer              ,intent(out) :: Status
126    integer                           :: i
127    integer                           :: j
128    integer                           :: stat
129    integer(hid_t)                    :: enum_type 
130    !    character (256)                    :: NullName
131
132    !    NullName = char(0)
133
134    do i=1,WrfDataHandleMax
135       if(WrfDataHandles(i)%Free) then
136          DH => WrfDataHandles(i)
137          DataHandle = i
138          DH%MaxTimeCount = 1
139
140          DH%FileID = -1
141          DH%GroupID = -1
142          DH%DimGroupID = -1
143
144          call SetUp_EnumID(enum_type,Status)
145          if(Status /= 0) then
146             Status = WRF_HDF5_ERR_ALLOCATION
147             write(msg,*) 'Fatal enum ALLOCATION ERROR in ',__FILE__,', line',__LINE__
148             call wrf_debug ( FATAL , msg)
149             return
150          endif
151          DH%EnumID = enum_type
152
153          allocate(DH%Times(MaxTimes), STAT=stat)
154          if(stat/= 0) then
155             Status = WRF_HDF5_ERR_ALLOCATION
156             write(msg,*) 'Fatal ALLOCATION ERROR in ',__FILE__,', line',__LINE__
157             call wrf_debug ( FATAL , msg)
158             return
159          endif
160          ! wait in the future
161          !          DH%Times(1:MaxTimes) = NullName
162
163          allocate(DH%DimLengths(MaxDims), STAT=stat)
164          if(stat/= 0) then
165             Status = WRF_HDF5_ERR_ALLOCATION
166             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line',__LINE__
167
168             call wrf_debug ( FATAL , msg)
169             return
170          endif
171
172          allocate(DH%DimIDs(MaxDims), STAT=stat)
173          if(stat/= 0) then
174             Status = WRF_HDF5_ERR_ALLOCATION
175             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
176             call wrf_debug ( FATAL , msg)
177             return
178          endif
179
180          allocate(DH%DimNames(MaxDims), STAT=stat)
181          if(stat/= 0) then
182             Status = WRF_HDF5_ERR_ALLOCATION
183             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
184             call wrf_debug ( FATAL , msg)
185             return
186          endif
187
188          allocate(DH%DIMTABLE(MaxTabDims), STAT = stat)
189          if(stat/= 0) then
190             Status = WRF_HDF5_ERR_ALLOCATION
191             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
192             call wrf_debug ( FATAL , msg)
193             return
194          endif
195
196          do j =1,MaxTabDims
197             DH%DIMTABLE(j)%dim_name = NO_NAME
198             DH%DIMTABLE(j)%unlimited = -1
199          enddo
200
201          allocate(DH%MDDsetIDs(MaxVars), STAT=stat)
202          if(stat/= 0) then
203             Status = WRF_HDF5_ERR_ALLOCATION
204             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
205             call wrf_debug ( FATAL , msg)
206             return
207          endif
208
209          allocate(DH%MDVarDimLens(MaxVars), STAT=stat)
210          if(stat/= 0) then
211             Status = WRF_HDF5_ERR_ALLOCATION
212             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
213             call wrf_debug ( FATAL , msg)
214             return
215          endif
216
217          allocate(DH%MDVarNames(MaxVars), STAT=stat)
218          if(stat/= 0) then
219             Status = WRF_HDF5_ERR_ALLOCATION
220             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
221             call wrf_debug ( FATAL , msg)
222             return
223          endif
224
225          allocate(DH%DsetIDs(MaxVars), STAT=stat)
226          if(stat/= 0) then
227             Status = WRF_HDF5_ERR_ALLOCATION
228             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
229             call wrf_debug ( FATAL , msg)
230             return
231          endif
232          DH%DsetIDs = -1
233
234          allocate(DH%TgroupIDs(MaxTimes), STAT=stat)
235          if(stat/= 0) then
236             Status = WRF_HDF5_ERR_ALLOCATION
237             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
238             call wrf_debug ( FATAL , msg)
239             return
240          endif
241          DH%TgroupIDs = -1
242
243          allocate(DH%VarDimLens(NVarDims-1,MaxVars), STAT=stat)
244          if(stat/= 0) then
245             Status = WRF_HDF5_ERR_ALLOCATION
246             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
247             call wrf_debug ( FATAL , msg)
248             return
249          endif
250
251          allocate(DH%VarNames(MaxVars), STAT=stat)
252          if(stat/= 0) then
253             Status = WRF_HDF5_ERR_ALLOCATION
254             write(msg,*) 'Fatal ALLOCATION ERROR in ',"__FILE__",', line', __LINE__
255             call wrf_debug ( FATAL , msg)
256             return
257          endif
258          exit
259       endif
260
261       if(i==WrfDataHandleMax) then
262          Status = WRF_HDF5_ERR_TOO_MANY_FILES
263          write(msg,*) 'Warning TOO MANY FILES in ',"__FILE__",', line', __LINE__
264          call wrf_debug ( WARN , msg)
265          return
266       endif
267    enddo
268
269
270    DH%Free      =.false.
271    DH%Comm      = Comm
272    DH%Write     =.false.
273    DH%first_operation  = .TRUE.
274    Status       = WRF_NO_ERR
275  end subroutine allocHandle
276
277  ! Obtain data handler
278  subroutine GetDH(DataHandle,DH,Status)
279
280    use wrf_phdf5_data
281    include 'wrf_status_codes.h'
282    integer               ,intent(in)          :: DataHandle
283    type(wrf_phdf5_data_handle) ,pointer        :: DH
284    integer               ,intent(out)         :: Status
285
286    if(DataHandle < 1 .or. DataHandle > WrfDataHandleMax) then
287       Status = WRF_HDF5_ERR_BAD_DATA_HANDLE
288       return
289    endif
290    DH => WrfDataHandles(DataHandle)
291    if(DH%Free) then
292       Status = WRF_HDF5_ERR_BAD_DATA_HANDLE
293       return
294    endif
295    Status = WRF_NO_ERR
296    return
297  end subroutine GetDH
298
299  ! Set up eumerate datatype for possible logical type
300  subroutine SetUp_EnumID(enum_type,Status)
301
302    use wrf_phdf5_data
303    use HDF5
304    implicit none
305    include 'wrf_status_codes.h'
306    integer(hid_t)              ,intent(out)  :: enum_type
307    integer                     ,intent(out)  :: Status
308    integer                                   :: hdf5err
309    integer, dimension(2)                     :: data
310
311    data(1) = 1
312    data(2) = 0
313
314    call h5tenum_create_f(H5T_NATIVE_INTEGER,enum_type,hdf5err)
315    if(hdf5err.lt.0) then
316       Status =  WRF_HDF5_ERR_DATATYPE
317       write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
318       call wrf_debug ( WARN , msg)
319       return
320    endif
321
322    call h5tenum_insert_f(enum_type,hdf5_true,data(1),hdf5err)
323    if(hdf5err.lt.0) then
324       Status =  WRF_HDF5_ERR_DATATYPE
325       write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
326       call wrf_debug ( WARN , msg)
327       return
328    endif
329
330    call h5tenum_insert_f(enum_type,hdf5_false,data(2),Status)
331    if(hdf5err.lt.0) then
332       Status =  WRF_HDF5_ERR_DATATYPE
333       write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
334       call wrf_debug ( WARN , msg)
335       return
336    endif
337
338    Status = WRF_NO_ERR
339    return
340  end subroutine SetUp_EnumID
341
342! Returns .TRUE. iff it is OK to write time-independent domain metadata to the
343! file referenced by DataHandle.  If DataHandle is invalid, .FALSE. is
344! returned.
345LOGICAL FUNCTION phdf5_ok_to_put_dom_ti( DataHandle )
346    use wrf_phdf5_data
347    include 'wrf_status_codes.h'
348    INTEGER, INTENT(IN) :: DataHandle
349    CHARACTER*80 :: fname
350    INTEGER :: filestate
351    INTEGER :: Status
352    LOGICAL :: dryrun, first_output, retval
353    call ext_phdf5_inquire_filename( DataHandle, fname, filestate, Status )
354    IF ( Status /= WRF_NO_ERR ) THEN
355      write(msg,*) 'Warning Status = ',Status,' in ',__FILE__, &
356                   ', line', __LINE__
357      call wrf_debug ( WARN , TRIM(msg) )
358      retval = .FALSE.
359    ELSE
360      dryrun       = ( filestate .EQ. WRF_FILE_OPENED_NOT_COMMITTED )
361      first_output = phdf5_is_first_operation( DataHandle )
362      retval = .NOT. dryrun .AND. first_output
363    ENDIF
364    phdf5_ok_to_put_dom_ti = retval
365    RETURN
366END FUNCTION phdf5_ok_to_put_dom_ti
367
368! Returns .TRUE. iff it is OK to read time-independent domain metadata from the
369! file referenced by DataHandle.  If DataHandle is invalid, .FALSE. is
370! returned.
371LOGICAL FUNCTION phdf5_ok_to_get_dom_ti( DataHandle )
372    use wrf_phdf5_data
373    include 'wrf_status_codes.h'
374    INTEGER, INTENT(IN) :: DataHandle
375    CHARACTER*80 :: fname
376    INTEGER :: filestate
377    INTEGER :: Status
378    LOGICAL :: dryrun, retval
379    call ext_phdf5_inquire_filename( DataHandle, fname, filestate, Status )
380    IF ( Status /= WRF_NO_ERR ) THEN
381      write(msg,*) 'Warning Status = ',Status,' in ',__FILE__, &
382                   ', line', __LINE__
383      call wrf_debug ( WARN , TRIM(msg) )
384      retval = .FALSE.
385    ELSE
386      dryrun       = ( filestate .EQ. WRF_FILE_OPENED_NOT_COMMITTED )
387      retval = .NOT. dryrun
388    ENDIF
389    phdf5_ok_to_get_dom_ti = retval
390    RETURN
391END FUNCTION phdf5_ok_to_get_dom_ti
392
393! Returns .TRUE. iff nothing has been read from or written to the file
394! referenced by DataHandle.  If DataHandle is invalid, .FALSE. is returned.
395LOGICAL FUNCTION phdf5_is_first_operation( DataHandle )
396    use wrf_phdf5_data
397    INCLUDE 'wrf_status_codes.h'
398    INTEGER, INTENT(IN) :: DataHandle
399    TYPE(wrf_phdf5_data_handle) ,POINTER :: DH
400    INTEGER :: Status
401    LOGICAL :: retval
402    CALL GetDH( DataHandle, DH, Status )
403    IF ( Status /= WRF_NO_ERR ) THEN
404      write(msg,*) 'Warning Status = ',Status,' in ',__FILE__, &
405                   ', line', __LINE__
406      call wrf_debug ( WARN , TRIM(msg) )
407      retval = .FALSE.
408    ELSE
409      retval = DH%first_operation
410    ENDIF
411    phdf5_is_first_operation = retval
412    RETURN
413END FUNCTION phdf5_is_first_operation
414
415end module ext_phdf5_support_routines
416
417!module wrf_phdf5_opt_data
418!  integer       ,parameter  :: MaxOptVars = 100
419!end module wrf_phdf5_opt_data
420
421!module opt_data_module
422
423!use wrf_phdf5_opt_data
424!   type :: field
425
426!       logical         :: Free
427!       integer,pointer :: darrays(:)
428!       integer         :: index
429!   end type field
430!   type(field),target :: fieldhandle(MaxOptVars)
431!end module opt_data_module
432
433!module opt_support_module
434!   implicit none
435!contains
436!   subroutine alloc_opt_handle(ODH)
437!   use opt_data_module
438!   type(field),pointer::DH
439!   integer            :: i
440
441!   do i =1,MaxOptVars
442!       DH=>fieldhandle(i)
443!        DH%index = 0
444!   enddo
445!end module opt_support_module
446
447! check the date, only use the length
448subroutine DateCheck(Date,Status)
449  use wrf_phdf5_data
450  include 'wrf_status_codes.h'
451  character*(*) ,intent(in)      :: Date
452  integer       ,intent(out)     :: Status
453
454  if(len(Date) /= DateStrLen) then
455     Status = WRF_HDF5_ERR_DATESTR_BAD_LENGTH
456  else 
457     Status = WRF_NO_ERR
458  endif
459  return
460end subroutine DateCheck
461
462! This routine is for meta-data time dependent varible attribute
463subroutine GetName(Element,Var,Name,Status)
464
465  use wrf_phdf5_data
466  include 'wrf_status_codes.h'
467  character*(*) ,intent(in)     :: Element
468  character*(*) ,intent(in)     :: Var
469  character*(*) ,intent(out)    :: Name
470  integer       ,intent(out)    :: Status
471  character (VarNameLen)        :: VarName
472  character (1)                 :: c
473  integer                       :: i
474  integer, parameter            ::  upper_to_lower =IACHAR('a')-IACHAR('A')
475
476  VarName = Var
477  Name = 'MD___'//trim(Element)//VarName
478  do i=1,len(Name)
479     c=Name(i:i)
480     if('A'<=c .and. c <='Z') Name(i:i)=achar(iachar(c)+upper_to_lower)
481     if(c=='-'.or.c==':') Name(i:i)='_'
482  enddo
483  Status = WRF_NO_ERR
484  return
485end subroutine GetName
486
487! Obtain TimeIndex
488subroutine GetDataTimeIndex(IO,DataHandle,DateStr,TimeIndex,Status)
489
490  use HDF5
491  use wrf_phdf5_data
492  use ext_phdf5_support_routines
493
494  implicit none
495  include 'wrf_status_codes.h'
496
497  character (*)         ,intent(in)          :: IO
498  integer               ,intent(in)          :: DataHandle
499  character*(*)         ,intent(in)          :: DateStr
500  character (DateStrLen), pointer            :: TempTimes(:)
501  integer               ,intent(out)         :: TimeIndex
502  integer               ,intent(out)         :: Status
503
504  type(wrf_phdf5_data_handle) ,pointer        :: DH
505  integer                                    :: VStart(2)
506  integer                                    :: VCount(2)
507  integer                                    :: stat
508  integer                                    :: i
509  integer                                    :: PreTimeCount
510
511  integer                                    :: rank
512  integer(hsize_t), dimension(1)             :: chunk_dims =(/1/)
513  integer(hsize_t), dimension(1)             :: dims
514  integer(hsize_t), dimension(1)             :: hdf5_maxdims
515  integer(hsize_t), dimension(1)             :: offset
516  integer(hsize_t), dimension(1)             :: count
517  integer(hsize_t), dimension(1)             :: sizes
518
519  INTEGER(HID_T)                             :: dset_id   ! Dataset ID
520  INTEGER(HID_T)                             :: dspace_id ! Dataspace ID
521  INTEGER(HID_T)                             :: fspace_id ! Dataspace ID
522  INTEGER(HID_T)                             :: crp_list  ! chunk ID
523  integer(hid_t)                             :: str_id    ! string ID
524  integer                                    :: hdf5err
525
526  integer(hid_t)                             :: group_id
527  character(Len = 512)                       :: groupname
528
529  ! for debug
530
531  character(len=100)                        :: buf
532  integer(size_t)                            :: name_size
533  integer(size_t)                            :: datelen_size
534  ! suppose the output will not exceed 100,0000 timesteps.
535  character(Len = MaxTimeSLen)               :: tname
536
537
538  !  DH => WrfDataHandles(DataHandle), don't know why NetCDF doesn't use GetDH
539  call GetDH(DataHandle,DH,Status)
540  if(Status /= WRF_NO_ERR) then
541     write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
542     call wrf_debug ( WARN , msg)
543     return
544  endif
545
546  call DateCheck(DateStr,Status)
547  if(Status /= WRF_NO_ERR) then
548     Status =  WRF_HDF5_ERR_DATESTR_ERROR
549     write(msg,*) 'Warning DATE STRING ERROR in ',"__FILE__",', line', __LINE__
550     call wrf_debug ( WARN , msg)
551     return
552  endif
553
554  if(IO == 'write') then
555     TimeIndex = DH%TimeIndex
556     if(TimeIndex <= 0) then
557        TimeIndex = 1
558     elseif(DateStr < DH%Times(TimeIndex)) then
559        Status = WRF_HDF5_ERR_DATE_LT_LAST_DATE
560        write(msg,*) 'Warning DATE < LAST DATE in ',"__FILE__",', line', __LINE__
561        call wrf_debug ( WARN , msg)
562        return
563     elseif(DateStr == DH%Times(TimeIndex)) then
564        Status = WRF_NO_ERR
565        return
566     else
567        TimeIndex = TimeIndex + 1
568        !       If exceeding the maximum timestep, updating the maximum timestep
569        if(TimeIndex > MaxTimes*(DH%MaxTimeCount)) then
570           PreTimeCount = DH%MaxTimeCount
571           allocate(TempTimes(PreTimeCount*MaxTimes))
572           TempTimes(1:MaxTimes*PreTimeCount)=DH%Times(1:MaxTimes &
573                *PreTimeCount)
574           DH%MaxTimeCount = DH%MaxTimeCount +1
575           deallocate(DH%Times)
576           allocate(DH%Times(DH%MaxTimeCount*MaxTimes))
577           DH%Times(1:MaxTimes*PreTimeCount)=TempTimes(1:MaxTimes &
578                *PreTimeCount)
579           deallocate(TempTimes)
580        endif
581     endif
582     DH%TimeIndex        = TimeIndex
583     DH%Times(TimeIndex) = DateStr
584     !    From NetCDF implementation, keep it in case it can be used.
585     !     VStart(1) = 1
586     !     VStart(2) = TimeIndex
587     !     VCount(1) = DateStrLen
588     !     VCount(2) = 1
589
590     ! create memory dataspace id and file dataspace id
591     dims(1)   = 1
592     count(1)  = 1
593     offset(1) = TimeIndex -1
594     sizes(1)  = TimeIndex
595
596     ! create group id for different time stamp
597     call numtochar(TimeIndex,tname)
598     groupname = 'TIME_STAMP_'//tname
599!     call h5gn_members_f(DH%GroupID,DH%GroupName,nmembers,hdf5err)
600!     do i = 0, nmembers - 1
601!        call h5gget_obj_info_idx_f(DH%GroupID,DH%GroupName,i,ObjName, ObjType, &
602!                                   hdf5err)
603       
604!        if(ObjName(1:17) == groupname) then
605!          call h5gopen_f(DH%GroupID,groupname,tgroupid,hdf5err)
606!          exit
607!        endif
608!     enddo
609           
610     if(DH%Tgroupids(TimeIndex) == -1) then
611       call h5gcreate_f(DH%groupid,groupname,group_id,hdf5err)
612       if(hdf5err .lt. 0) then
613        Status = WRF_HDF5_ERR_GROUP
614        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
615        call wrf_debug ( WARN , msg)
616        return
617       endif
618       DH%Tgroupids(TimeIndex) = group_id
619     else
620!        call h5gopen_f(DH%groupid,groupname,group_id,
621       group_id = DH%Tgroupids(TimeIndex)
622     endif
623
624     call h5screate_simple_f(1,dims,dspace_id,hdf5err,dims)
625     if(hdf5err.lt.0) then
626        Status =  WRF_HDF5_ERR_DATASPACE
627        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
628        call wrf_debug ( WARN , msg)
629        return
630     endif
631
632
633     ! create HDF5 string handler for time
634     if(TimeIndex == 1) then
635        call h5tcopy_f(H5T_NATIVE_CHARACTER, str_id, hdf5err)
636        if(hdf5err.lt.0) then
637           Status =  WRF_HDF5_ERR_DATATYPE
638           write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
639           call wrf_debug ( WARN , msg)
640           return
641        endif
642
643        datelen_size = DateStrLen
644        call h5tset_size_f(str_id,datelen_size,hdf5err)
645        if(hdf5err.lt.0) then
646           Status =  WRF_HDF5_ERR_DATATYPE
647           write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
648           call wrf_debug ( WARN , msg)
649           return
650        endif
651     else
652        str_id = DH%str_id
653     endif
654
655     call h5dcreate_f(group_id,DH%TimesName,str_id,dspace_id,&
656          DH%TimesID, hdf5err)
657     if(hdf5err.lt.0) then
658        Status =  WRF_HDF5_ERR_DATASET_CREATE
659        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
660        call wrf_debug ( WARN , msg)
661        return
662     endif
663
664
665     ! write the data in memory space to file space
666     CALL h5dwrite_f(DH%TimesID,str_id,DateStr,dims,hdf5err,dspace_id,dspace_id)
667     if(hdf5err.lt.0) then
668        Status =  WRF_HDF5_ERR_DATASET_WRITE
669        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
670        call wrf_debug ( WARN , msg)
671        return
672     endif
673
674     if(TimeIndex == 1) then
675        DH%str_id = str_id
676     endif
677
678
679     call h5sclose_f(dspace_id,hdf5err)
680     if(hdf5err.lt.0) then
681        Status =  WRF_HDF5_ERR_DATASPACE
682        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
683        call wrf_debug ( WARN , msg)
684        return
685     endif
686
687     call h5dclose_f(DH%TimesID,hdf5err)
688     if(hdf5err.lt.0) then
689        Status = WRF_HDF5_ERR_DATASET_GENERAL
690        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
691        call wrf_debug ( WARN , msg)
692        return
693     endif
694
695  else
696     ! This is for IO read
697     ! Find the timeIndex(very expensive for large
698     ! time stamp, should use hashing table)
699
700     do i=1,MaxTimes*DH%MaxTimeCount
701
702        !       For handling reading maximum timestamp greater than 9000 in the future       
703        !        if(DH%Times(i) == NullName) then
704        !           Status = WRF_HDF5_ERR_TIME
705        !           write(msg,*) 'Warning TIME ',DateStr,' NOT FOUND in ',"__FILE__",&
706        !               ', line', __LINE__
707        !           call wrf_debug ( WARN , msg)
708        !           return
709        !        endif
710
711        if(DH%Times(i) == DateStr) then
712           Status = WRF_NO_ERR
713           TimeIndex = i
714           exit
715        endif
716
717        ! Need a recursive function to handle this
718        ! This is a potential bug
719        if(i == MaxTimes*DH%MaxTimeCount) then
720           !           PreTimeCount = DH%MaxTimeCount
721           !           allocate(TempTimes(PreTimeCount*MaxTimes))
722           !           TempTimes(1:MaxTimes*PreTimeCount)=DH%Times(1:MaxTimes &
723           !                *PreTimeCount)
724           !           DH%MaxTimeCount = DH%MaxTimeCount +1
725           !           deallocate(DH%Times)
726           !           allocate(DH%Times(DH%MaxTimeCount*MaxTimes))
727           !           DH%Times(1:MaxTimes*PreTimeCount)=TempTimes(1:MaxTimes &
728           !                *PreTimeCount)
729           !           deallocate(TempTimes)
730           Status = WRF_HDF5_ERR_TIME
731           write(msg,*) 'Warning TIME ',DateStr,' NOT FOUND in ',"__FILE__",&
732                ', line', __LINE__
733           call wrf_debug ( WARN , msg)
734           return
735        endif
736     enddo
737
738     ! do the hyperslab selection
739
740  endif
741  return
742end subroutine GetDataTimeIndex
743
744subroutine GetAttrTimeIndex(IO,DataHandle,DateStr,TimeIndex,Status)
745
746  use HDF5
747  use wrf_phdf5_data
748  use ext_phdf5_support_routines
749
750  implicit none
751  include 'wrf_status_codes.h'
752
753  character (*)         ,intent(in)          :: IO
754  integer               ,intent(in)          :: DataHandle
755  character*(*)         ,intent(in)          :: DateStr
756  character (DateStrLen), pointer            :: TempTimes(:)
757  integer               ,intent(out)         :: TimeIndex
758  integer               ,intent(out)         :: Status
759
760  type(wrf_phdf5_data_handle) ,pointer        :: DH
761  integer                                    :: VStart(2)
762  integer                                    :: VCount(2)
763  integer                                    :: stat
764  integer                                    :: i
765  integer                                    :: PreTimeCount
766
767  integer                                    :: rank
768  integer(hsize_t), dimension(1)             :: chunk_dims =(/1/)
769  integer(hsize_t), dimension(1)             :: dims
770  integer(hsize_t), dimension(1)             :: hdf5_maxdims
771  integer(hsize_t), dimension(1)             :: offset
772  integer(hsize_t), dimension(1)             :: count
773  integer(hsize_t), dimension(1)             :: sizes
774
775  INTEGER(HID_T)                             :: dset_id   ! Dataset ID
776  INTEGER(HID_T)                             :: dspace_id ! Dataspace ID
777  INTEGER(HID_T)                             :: fspace_id ! Dataspace ID
778  INTEGER(HID_T)                             :: crp_list  ! chunk ID
779  integer(hid_t)                             :: str_id    ! string ID
780  integer                                    :: hdf5err
781
782  integer(size_t)                            :: datelen_size
783  integer(hid_t)                             :: group_id
784  character(Len = 512)                       :: groupname
785
786  ! suppose the output will not exceed 100,0000 timesteps.
787  character(Len = MaxTimeSLen)               :: tname 
788
789  !  DH => WrfDataHandles(DataHandle), don't know why NetCDF doesn't use GetDH
790  call GetDH(DataHandle,DH,Status)
791  if(Status /= WRF_NO_ERR) then
792     write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
793     call wrf_debug ( WARN , msg)
794     return
795  endif
796
797  call DateCheck(DateStr,Status)
798  if(Status /= WRF_NO_ERR) then
799     Status =  WRF_HDF5_ERR_DATESTR_ERROR
800     write(msg,*) 'Warning DATE STRING ERROR in ',"__FILE__",', line', __LINE__
801     call wrf_debug ( WARN , msg)
802     return
803  endif
804
805  if(IO == 'write') then
806     TimeIndex = DH%TimeIndex
807     if(TimeIndex <= 0) then
808        TimeIndex = 1
809     elseif(DateStr < DH%Times(TimeIndex)) then
810        Status = WRF_HDF5_ERR_DATE_LT_LAST_DATE
811        write(msg,*) 'Warning DATE < LAST DATE in ',"__FILE__",', line', __LINE__
812        call wrf_debug ( WARN , msg)
813        return
814     elseif(DateStr == DH%Times(TimeIndex)) then
815        Status = WRF_NO_ERR
816        return
817     else
818        TimeIndex = TimeIndex + 1
819        !       If exceeding the maximum timestep, updating the maximum timestep
820        if(TimeIndex > MaxTimes*(DH%MaxTimeCount)) then
821           PreTimeCount = DH%MaxTimeCount
822           allocate(TempTimes(PreTimeCount*MaxTimes))
823           TempTimes(1:MaxTimes*PreTimeCount)=DH%Times(1:MaxTimes &
824                *PreTimeCount)
825           DH%MaxTimeCount = DH%MaxTimeCount +1
826           deallocate(DH%Times)
827           allocate(DH%Times(DH%MaxTimeCount*MaxTimes))
828           DH%Times(1:MaxTimes*PreTimeCount)=TempTimes(1:MaxTimes &
829                *PreTimeCount)
830           deallocate(TempTimes)
831        endif
832     endif
833     DH%TimeIndex        = TimeIndex
834     DH%Times(TimeIndex) = DateStr
835
836     !    From NetCDF implementation, keep it in case it can be used.
837     !     VStart(1) = 1
838     !     VStart(2) = TimeIndex
839     !     VCount(1) = DateStrLen
840     !     VCount(2) = 1
841
842     ! create memory dataspace id and file dataspace id
843     dims(1)   = 1
844     count(1)  = 1
845     offset(1) = TimeIndex -1
846     sizes(1)  = TimeIndex
847
848     ! create group id for different time stamp
849     call numtochar(TimeIndex,tname)
850     groupname = 'TIME_STAMP_'//tname
851           
852     if(DH%Tgroupids(TimeIndex) == -1) then
853       call h5gcreate_f(DH%groupid,groupname,group_id,hdf5err)
854       if(hdf5err .lt. 0) then
855        Status = WRF_HDF5_ERR_GROUP
856        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
857        call wrf_debug ( WARN , msg)
858        return
859       endif
860       DH%Tgroupids(TimeIndex) = group_id
861     else
862!        call h5gopen_f(DH%groupid,groupname,group_id,
863       group_id = DH%Tgroupids(TimeIndex)
864     endif
865
866     call h5screate_simple_f(1,dims,dspace_id,hdf5err,dims)
867     if(hdf5err.lt.0) then
868        Status =  WRF_HDF5_ERR_DATASPACE
869        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
870        call wrf_debug ( WARN , msg)
871        return
872     endif
873
874
875     ! create HDF5 string handler for time
876     if(TimeIndex == 1) then
877        call h5tcopy_f(H5T_NATIVE_CHARACTER, str_id, hdf5err)
878        if(hdf5err.lt.0) then
879           Status =  WRF_HDF5_ERR_DATATYPE
880           write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
881           call wrf_debug ( WARN , msg)
882           return
883        endif
884
885        datelen_size = DateStrLen
886        call h5tset_size_f(str_id,datelen_size,hdf5err)
887        if(hdf5err.lt.0) then
888           Status =  WRF_HDF5_ERR_DATATYPE
889           write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
890           call wrf_debug ( WARN , msg)
891           return
892        endif
893     else
894        str_id = DH%str_id
895     endif
896
897     call h5dcreate_f(group_id,DH%TimesName,str_id,dspace_id,&
898          DH%TimesID, hdf5err)
899     if(hdf5err.lt.0) then
900        Status =  WRF_HDF5_ERR_DATASET_CREATE
901        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
902        call wrf_debug ( WARN , msg)
903        return
904     endif
905
906
907     ! write the data in memory space to file space
908     CALL h5dwrite_f(DH%TimesID,str_id,DateStr,dims,hdf5err,dspace_id,dspace_id)
909     if(hdf5err.lt.0) then
910        Status =  WRF_HDF5_ERR_DATASET_WRITE
911        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
912        call wrf_debug ( WARN , msg)
913        return
914     endif
915
916     if(TimeIndex == 1) then
917        DH%str_id = str_id
918     endif
919
920
921     call h5sclose_f(dspace_id,hdf5err)
922     if(hdf5err.lt.0) then
923        Status =  WRF_HDF5_ERR_DATASPACE
924        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
925        call wrf_debug ( WARN , msg)
926        return
927     endif
928
929     call h5dclose_f(DH%TimesID,hdf5err)
930     if(hdf5err.lt.0) then
931        Status = WRF_HDF5_ERR_DATASET_GENERAL
932        write(msg,*) 'Warning Status = ',Status,' in ',__FILE__,', line', __LINE__
933        call wrf_debug ( WARN , msg)
934        return
935     endif
936
937  else
938     ! This is for IO read
939     ! Find the timeIndex(very expensive for large
940     ! time stamp, should use hashing table)
941
942     do i=1,MaxTimes*DH%MaxTimeCount
943
944
945        if(DH%Times(i) == DateStr) then
946           Status = WRF_NO_ERR
947           TimeIndex = i
948           exit
949        endif
950
951        ! Need a recursive function to handle this
952        ! This is a potential bug
953        if(i == MaxTimes*DH%MaxTimeCount) then
954           Status = WRF_HDF5_ERR_TIME
955           write(msg,*) 'Warning TIME ',DateStr,' NOT FOUND in ',"__FILE__",&
956                ', line', __LINE__
957           call wrf_debug ( WARN , msg)
958           return
959        endif
960     enddo
961
962     ! do the hyperslab selection
963
964  endif
965  return
966end subroutine GetAttrTimeIndex
967
968
969! Obtain the rank of the dimension
970subroutine GetDim(MemoryOrder,NDim,Status)
971
972  include 'wrf_status_codes.h'
973  character*(*) ,intent(in)  :: MemoryOrder
974  integer       ,intent(out) :: NDim
975  integer       ,intent(out) :: Status
976  character*3                :: MemOrd
977
978  call LowerCase(MemoryOrder,MemOrd)
979  select case (MemOrd)
980  case ('xyz','xzy','yxz','yzx','zxy','zyx','xsz','xez','ysz','yez')
981     NDim = 3
982  case ('xy','yx','xs','xe','ys','ye')
983     NDim = 2
984  case ('z','c','0')
985     NDim = 1
986  case default
987     Status = WRF_HDF5_ERR_BAD_MEMORYORDER
988     return
989  end select
990  Status = WRF_NO_ERR
991  return
992end subroutine GetDim
993
994! Obtain the index for transposing
995subroutine GetIndices(NDim,Start,End,i1,i2,j1,j2,k1,k2)
996  integer              ,intent(in)  :: NDim
997  integer ,dimension(*),intent(in)  :: Start,End
998  integer              ,intent(out) :: i1,i2,j1,j2,k1,k2
999
1000  i1=1
1001  i2=1
1002  j1=1
1003  j2=1
1004  k1=1
1005  k2=1
1006  i1 = Start(1)
1007  i2 = End  (1)
1008  if(NDim == 1) return
1009  j1 = Start(2)
1010  j2 = End  (2)
1011  if(NDim == 2) return
1012  k1 = Start(3)
1013  k2 = End  (3)
1014  return
1015end subroutine GetIndices
1016
1017! shuffling the memory order to XYZ order
1018subroutine ExtOrder(MemoryOrder,Vector,Status)
1019  use wrf_phdf5_data
1020  include 'wrf_status_codes.h'
1021  character*(*)              ,intent(in)    :: MemoryOrder
1022  integer,dimension(*)       ,intent(inout) :: Vector
1023  integer                    ,intent(out)   :: Status
1024  integer                                   :: NDim
1025  integer,dimension(NVarDims)               :: temp
1026  character*3                               :: MemOrd
1027
1028  call GetDim(MemoryOrder,NDim,Status)
1029  temp(1:NDim) = Vector(1:NDim)
1030  call LowerCase(MemoryOrder,MemOrd)
1031  select case (MemOrd)
1032
1033  case ('xyz','xsz','xez','ysz','yez','xy','xs','xe','ys','ye','z','c')
1034     continue
1035  case ('0')
1036     Vector(1) = 1
1037  case ('xzy')
1038     Vector(2) = temp(3)
1039     Vector(3) = temp(2)
1040  case ('yxz')
1041     Vector(1) = temp(2)
1042     Vector(2) = temp(1)
1043  case ('yzx')
1044     Vector(1) = temp(3)
1045     Vector(2) = temp(1)
1046     Vector(3) = temp(2)
1047  case ('zxy')
1048     Vector(1) = temp(2)
1049     Vector(2) = temp(3)
1050     Vector(3) = temp(1)
1051  case ('zyx')
1052     Vector(1) = temp(3)
1053     Vector(3) = temp(1)
1054  case ('yx')
1055     Vector(1) = temp(2)
1056     Vector(2) = temp(1)
1057  case default
1058     Status = WRF_HDF5_ERR_BAD_MEMORYORDER
1059     return
1060  end select
1061  Status = WRF_NO_ERR
1062  return
1063end subroutine ExtOrder
1064
1065! shuffling the dimensional name order
1066subroutine ExtOrderStr(MemoryOrder,Vector,ROVector,Status)
1067  use wrf_phdf5_data
1068  include 'wrf_status_codes.h'
1069  character*(*)                    ,intent(in)    :: MemoryOrder
1070  character*(*),dimension(*)       ,intent(in)    :: Vector
1071  character(256),dimension(NVarDims),intent(out)   :: ROVector
1072  integer                          ,intent(out)   :: Status
1073  integer                                         :: NDim
1074  character*3                                     :: MemOrd
1075
1076  call GetDim(MemoryOrder,NDim,Status)
1077  ROVector(1:NDim) = Vector(1:NDim)
1078  call LowerCase(MemoryOrder,MemOrd)
1079  select case (MemOrd)
1080
1081  case ('xyz','xsz','xez','ysz','yez','xy','xs','xe','ys','ye','z','c')
1082     continue
1083  case ('0')
1084     ROVector(1) = 'ext_scalar'
1085  case ('xzy')
1086     ROVector(2) = Vector(3)
1087     ROVector(3) = Vector(2)
1088  case ('yxz')
1089     ROVector(1) = Vector(2)
1090     ROVector(2) = Vector(1)
1091  case ('yzx')
1092     ROVector(1) = Vector(3)
1093     ROVector(2) = Vector(1)
1094     ROVector(3) = Vector(2)
1095  case ('zxy')
1096     ROVector(1) = Vector(2)
1097     ROVector(2) = Vector(3)
1098     ROVector(3) = Vector(1)
1099  case ('zyx')
1100     ROVector(1) = Vector(3)
1101     ROVector(3) = Vector(1)
1102  case ('yx')
1103     ROVector(1) = Vector(2)
1104     ROVector(2) = Vector(1)
1105  case default
1106     Status = WRF_HDF5_ERR_BAD_MEMORYORDER
1107     return
1108  end select
1109  Status = WRF_NO_ERR
1110  return
1111end subroutine ExtOrderStr
1112
1113
1114subroutine LowerCase(MemoryOrder,MemOrd)
1115  character*(*) ,intent(in)  :: MemoryOrder
1116  character*(*) ,intent(out) :: MemOrd
1117  character*3                :: c
1118  integer       ,parameter   :: upper_to_lower =IACHAR('a')-IACHAR('A')
1119  integer                    :: i,N
1120
1121  MemOrd = ' '
1122  N = len(MemoryOrder)
1123  MemOrd(1:N) = MemoryOrder(1:N)
1124  do i=1,N
1125     c = MemoryOrder(i:i)
1126     if('A'<=c .and. c <='Z') MemOrd(i:i)=achar(iachar(c)+upper_to_lower)
1127  enddo
1128  return
1129end subroutine LowerCase
1130
1131subroutine UpperCase(MemoryOrder,MemOrd)
1132  character*(*) ,intent(in)  :: MemoryOrder
1133  character*(*) ,intent(out) :: MemOrd
1134  character*3                :: c
1135  integer     ,parameter     :: lower_to_upper =IACHAR('A')-IACHAR('a')
1136  integer                    :: i,N
1137
1138  MemOrd = ' '
1139  N = len(MemoryOrder)
1140  MemOrd(1:N) = MemoryOrder(1:N)
1141  do i=1,N
1142     c = MemoryOrder(i:i)
1143     if('a'<=c .and. c <='z') MemOrd(i:i)=achar(iachar(c)+lower_to_upper)
1144  enddo
1145  return
1146end subroutine UpperCase
1147
1148! subroutine used in transpose routine
1149subroutine reorder (MemoryOrder,MemO)
1150  character*(*)     ,intent(in)    :: MemoryOrder
1151  character*3       ,intent(out)   :: MemO
1152  character*3                      :: MemOrd
1153  integer                          :: N,i,i1,i2,i3
1154
1155  MemO = MemoryOrder
1156  N = len_trim(MemoryOrder)
1157  if(N == 1) return
1158  call lowercase(MemoryOrder,MemOrd)
1159  i1 = 1
1160  i3 = 1
1161  do i=2,N
1162     if(ichar(MemOrd(i:i)) < ichar(MemOrd(i1:i1))) I1 = i
1163     if(ichar(MemOrd(i:i)) > ichar(MemOrd(i3:i3))) I3 = i
1164  enddo
1165  if(N == 2) then
1166     i2=i3
1167  else
1168     i2 = 6-i1-i3
1169  endif
1170  MemO(1:1) = MemoryOrder(i1:i1)
1171  MemO(2:2) = MemoryOrder(i2:i2)
1172  if(N == 3) MemO(3:3) = MemoryOrder(i3:i3)
1173  if(MemOrd(i1:i1) == 's' .or. MemOrd(i1:i1) == 'e') then
1174     MemO(1:N-1) = MemO(2:N)
1175     MemO(N:N  ) = MemoryOrder(i1:i1)
1176  endif
1177  return
1178end subroutine reorder
1179
1180subroutine Transpose_hdf5(IO,MemoryOrder,di, Field,l1,l2,m1,m2,n1,n2 &
1181     ,XField,x1,x2,y1,y2,z1,z2 &
1182     ,i1,i2,j1,j2,k1,k2 )
1183  character*(*)     ,intent(in)    :: IO
1184  character*(*)     ,intent(in)    :: MemoryOrder
1185  integer           ,intent(in)    :: l1,l2,m1,m2,n1,n2
1186  integer           ,intent(in)    :: di
1187  integer           ,intent(in)    :: x1,x2,y1,y2,z1,z2
1188  integer           ,intent(in)    :: i1,i2,j1,j2,k1,k2
1189  integer           ,intent(inout) ::  Field(di,l1:l2,m1:m2,n1:n2)
1190  !jm 010827  integer           ,intent(inout) :: XField(di,x1:x2,y1:y2,z1:z2)
1191  integer           ,intent(inout) :: XField(di,(i2-i1+1)*(j2-j1+1)*(k2-k1+1))
1192  character*3                      :: MemOrd
1193  character*3                      :: MemO
1194  integer           ,parameter     :: MaxUpperCase=IACHAR('Z')
1195  integer                          :: i,j,k,ix,jx,kx
1196
1197  call LowerCase(MemoryOrder,MemOrd)
1198  select case (MemOrd)
1199
1200     !#define XDEX(A,B,C) A-A ## 1+1+(A ## 2-A ## 1+1)*((B-B ## 1)+(C-C ## 1)*(B ## 2-B ## 1+1))
1201
1202
1203  case ('xzy')
1204
1205     ix=0
1206     jx=0
1207     kx=0
1208     call reorder(MemoryOrder,MemO)
1209     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1210     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1211     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1212     do k=k1,k2
1213        do j=j1,j2
1214           do i=i1,i2
1215              if(IO == 'write') then
1216                 XField(1:di,(i-i1+1+(i2-i1+1)*((k-k1)+(j-j1)*(k2-k1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1217              else
1218                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(i-i1+1+(i2-i1+1)*((k-k1)+(j-j1)*(k2-k1+1))))
1219              endif
1220           enddo
1221        enddo
1222     enddo
1223     return
1224
1225  case ('xyz','xsz','xez','ysz','yez','xy','xs','xe','ys','ye','z','c','0')
1226
1227
1228     ix=0
1229     jx=0
1230     kx=0
1231     call reorder(MemoryOrder,MemO)
1232     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1233     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1234     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1235     do k=k1,k2
1236        do j=j1,j2
1237           do i=i1,i2
1238              if(IO == 'write') then
1239                 XField(1:di,(i-i1+1+(i2-i1+1)*((j-j1)+(k-k1)*(j2-j1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1240              else
1241                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(i-i1+1+(i2-i1+1)*((j-j1)+(k-k1)*(j2-j1+1))))
1242              endif
1243           enddo
1244        enddo
1245     enddo
1246     return
1247
1248  case ('yxz')
1249
1250
1251     ix=0
1252     jx=0
1253     kx=0
1254     call reorder(MemoryOrder,MemO)
1255     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1256     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1257     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1258     do k=k1,k2
1259        do j=j1,j2
1260           do i=i1,i2
1261              if(IO == 'write') then
1262                 XField(1:di,(j-j1+1+(j2-j1+1)*((i-i1)+(k-k1)*(i2-i1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1263              else
1264                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(j-j1+1+(j2-j1+1)*((i-i1)+(k-k1)*(i2-i1+1))))
1265              endif
1266           enddo
1267        enddo
1268     enddo
1269     return
1270
1271  case ('zxy')
1272
1273
1274     ix=0
1275     jx=0
1276     kx=0
1277     call reorder(MemoryOrder,MemO)
1278     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1279     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1280     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1281     do k=k1,k2
1282        do j=j1,j2
1283           do i=i1,i2
1284              if(IO == 'write') then
1285                 XField(1:di,(k-k1+1+(k2-k1+1)*((i-i1)+(j-j1)*(i2-i1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1286              else
1287                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(k-k1+1+(k2-k1+1)*((i-i1)+(j-j1)*(i2-i1+1))))
1288              endif
1289           enddo
1290        enddo
1291     enddo
1292     return
1293
1294  case ('yzx')
1295
1296
1297     ix=0
1298     jx=0
1299     kx=0
1300     call reorder(MemoryOrder,MemO)
1301     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1302     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1303     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1304     do k=k1,k2
1305        do j=j1,j2
1306           do i=i1,i2
1307              if(IO == 'write') then
1308                 XField(1:di,(j-j1+1+(j2-j1+1)*((k-k1)+(i-i1)*(k2-k1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1309              else
1310                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(j-j1+1+(j2-j1+1)*((k-k1)+(i-i1)*(k2-k1+1))))
1311              endif
1312           enddo
1313        enddo
1314     enddo
1315     return
1316
1317  case ('zyx')
1318
1319
1320     ix=0
1321     jx=0
1322     kx=0
1323     call reorder(MemoryOrder,MemO)
1324     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1325     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1326     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1327     do k=k1,k2
1328        do j=j1,j2
1329           do i=i1,i2
1330              if(IO == 'write') then
1331                 XField(1:di,(k-k1+1+(k2-k1+1)*((j-j1)+(i-i1)*(j2-j1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1332              else
1333                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(k-k1+1+(k2-k1+1)*((j-j1)+(i-i1)*(j2-j1+1))))
1334              endif
1335           enddo
1336        enddo
1337     enddo
1338     return
1339
1340  case ('yx')
1341
1342
1343     ix=0
1344     jx=0
1345     kx=0
1346     call reorder(MemoryOrder,MemO)
1347     if(IACHAR(MemO(1:1)) > MaxUpperCase) ix=i2+i1
1348     if(IACHAR(MemO(2:2)) > MaxUpperCase) jx=j2+j1
1349     if(IACHAR(MemO(3:3)) > MaxUpperCase) kx=k2+k1
1350     do k=k1,k2
1351        do j=j1,j2
1352           do i=i1,i2
1353              if(IO == 'write') then
1354                 XField(1:di,(j-j1+1+(j2-j1+1)*((i-i1)+(k-k1)*(i2-i1+1)))) = Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k))
1355              else
1356                 Field(1:di,abs(ix-i),abs(jx-j),abs(kx-k)) = XField(1:di,(j-j1+1+(j2-j1+1)*((i-i1)+(k-k1)*(i2-i1+1))))
1357              endif
1358           enddo
1359        enddo
1360     enddo
1361     return
1362
1363  end select
1364  return
1365end subroutine Transpose_hdf5
1366
1367subroutine numtochar(TimeIndex,tname,Status)
1368
1369  use wrf_phdf5_data
1370  integer, intent(in) :: TimeIndex
1371  character(len=MaxTimeSLen),intent(out)::tname
1372  integer                   ,intent(out)::Status
1373  integer             :: i,ten_pow,temp
1374  integer             :: maxtimestep
1375
1376  maxtimestep =1
1377  do i =1,MaxTimeSLen
1378     maxtimestep = maxtimestep * 10
1379  enddo
1380  if(TimeIndex >= maxtimestep) then
1381     Status = WRF_HDF5_ERR_OTHERS
1382     write(msg,*) 'Cannot exceed the maximum timestep',maxtimestep,'in',__FILE__,' line',__LINE__
1383     call wrf_debug(FATAL,msg)
1384     return
1385  endif
1386
1387  ten_pow = 1
1388  temp =10
1389  do i =1,MaxTimeSLen
1390     tname(MaxTimeSLen+1-i:MaxTimeSLen+1-i) = achar(modulo(TimeIndex/ten_pow,temp)+iachar('0'))
1391     ten_pow  = 10* ten_pow
1392  enddo
1393
1394  return
1395end subroutine numtochar
Note: See TracBrowser for help on using the repository browser.