!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFREAD ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfread(nx,ny,nz,nzsoil,nstyps, grdbas, filename, time, & 2,121
x,y,z,zp,zpsoil, &
uprt, vprt, wprt, ptprt, pprt, &
qvprt, qc, qr, qi, qs, qh, tke, kmh,kmv, &
ubar, vbar, wbar, ptbar, pbar, rhobar, qvbar, &
soiltyp,stypfrct,vegtyp,lai,roufns,veg, &
tsoil,qsoil,wetcanp,snowdpth, &
raing,rainc,prcrate, &
radfrc,radsw,rnflx,radswnet,radlwin, &
usflx,vsflx,ptsflx,qvsflx, &
ireturn, tem1)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in history data in the NCSA HDF4 format.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/15
!
! MODIFICATION HISTORY:
!
! 05/15/2002 (J. Brotzge)
! Added to allow for multiple soil schemes
!
!-----------------------------------------------------------------------
!
! INPUT :
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
! nzsoil Number of grid points in the soil
!
! grdbas Data read flag.
! =1, only grid and base state arrays will be read
! =0, all arrays will be read based on data
! parameter setting.
! filename Character variable nhming the input HDF file
!
! OUTPUT:
!
! time Time in seconds of data read from "filename"
!
! x x coordinate of grid points in physical/comp. space (m)
! y y coordinate of grid points in physical/comp. space (m)
! z z coordinate of grid points in computational space (m)
! zp z coordinate of grid points in physical space (m)
! zpsoil z coordinate of grid points in the soil (m)
!
! uprt x component of perturbation velocity (m/s)
! vprt y component of perturbation velocity (m/s)
! wprt Vertical component of perturbation velocity in Cartesian
! coordinates (m/s).
! ptprt Perturbation potential temperature (K)
! pprt Perturbation pressure (Pascal)
!
! qvprt Perturbation water vapor mixing ratio (kg/kg)
! qc Cloud water mixing ratio (kg/kg)
! qr Rainwater mixing ratio (kg/kg)
! qi Cloud ice mixing ratio (kg/kg)
! qs Snow mixing ratio (kg/kg)
! qh Hail mixing ratio (kg/kg)
! tke Turbulent Kinetic Energy ((m/s)**2)
! kmh Horizontal turb. mixing coef. for momentum ( m**2/s )
! kmv Vertical turb. mixing coef. for momentum ( m**2/s )
!
! ubar Base state x velocity component (m/s)
! vbar Base state y velocity component (m/s)
! wbar Base state z velocity component (m/s)
! ptbar Base state potential temperature (K)
! pbar Base state pressure (Pascal)
! rhobar Base state air density (kg/m**3)
! qvbar Base state water vapor mixing ratio (kg/kg)
!
! soiltyp Soil type
! vegtyp Vegetation type
! lai Leaf Area Index
! roufns Surface roughness
! veg Vegetation fraction
!
! tsoil Soil temperature (K)
! qsoil Soil moisture (m**3/m**3)
! wetcanp Canopy water amount
!
! raing Grid supersaturation rain
! rainc Cumulus convective rain
! prcrate Precipitation rates
!
! radfrc Radiation forcing (K/s)
! radsw Solar radiation reaching the surface
! rnflx Net radiation flux absorbed by surface
! radswnet Net shortwave radiation
! radlwin Incoming longwave radiation
!
! usflx Surface flux of u-momentum (kg/(m*s**2))
! vsflx Surface flux of v-momentum (kg/(m*s**2))
! ptsflx Surface heat flux (K*kg/(m**2 * s ))
! qvsflx Surface moisture flux of (kg/(m**2 * s))
!
! ireturn Return status indicator
!
! WORK ARRAY
!
! tem1
! tem2 work array.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz
INTEGER :: nzsoil
INTEGER :: grdbas
CHARACTER (LEN=*) :: filename
REAL :: x (nx) ! x coord.
REAL :: y (ny) ! y coord.
REAL :: z (nz) ! z coord.
REAL :: zp (nx,ny,nz) ! physical x coord.
REAL :: zpsoil(nx,ny,nzsoil) ! physical x coord. for soil (m)
REAL :: uprt (nx,ny,nz) ! Perturbation u-velocity (m/s)
REAL :: vprt (nx,ny,nz) ! Perturbation v-velocity (m/s)
REAL :: wprt (nx,ny,nz) ! Perturbation w-velocity (m/s)
REAL :: ptprt (nx,ny,nz) ! Perturbation potential temperature (K)
REAL :: pprt (nx,ny,nz) ! Perturbation pressure (Pascal)
REAL :: qvprt (nx,ny,nz) ! Perturbation water vapor mixing ratio (kg/kg)
REAL :: qc (nx,ny,nz) ! Cloud water mixing ratio (kg/kg)
REAL :: qr (nx,ny,nz) ! Rain water mixing ratio (kg/kg)
REAL :: qi (nx,ny,nz) ! Cloud ice mixing ratio (kg/kg)
REAL :: qs (nx,ny,nz) ! Snow mixing ratio (kg/kg)
REAL :: qh (nx,ny,nz) ! Hail mixing ratio (kg/kg)
REAL :: tke (nx,ny,nz) ! Turbulent Kinetic Energy ((m/s)**2)
REAL :: kmh (nx,ny,nz) ! Horizontal turb. mixing coef. for
! momentum. ( m**2/s )
REAL :: kmv (nx,ny,nz) ! Vertical turb. mixing coef. for
! momentum. ( m**2/s )
REAL :: ubar (nx,ny,nz) ! Base state u-velocity (m/s)
REAL :: vbar (nx,ny,nz) ! Base state v-velocity (m/s)
REAL :: wbar (nx,ny,nz) ! Base state w-velocity (m/s)
REAL :: ptbar (nx,ny,nz) ! Base state potential temperature (K)
REAL :: pbar (nx,ny,nz) ! Base state pressure (Pascal).
REAL :: rhobar(nx,ny,nz) ! Base state air density (kg/m**3)
REAL :: qvbar (nx,ny,nz) ! Base state water vapor mixing ratio
INTEGER :: nstyps ! Number of soil type
INTEGER :: soiltyp (nx,ny,nstyps) ! Soil type
REAL :: stypfrct(nx,ny,nstyps) ! Fraction of soil types
INTEGER :: vegtyp(nx,ny) ! Vegetation type
REAL :: lai (nx,ny) ! Leaf Area Index
REAL :: roufns (nx,ny) ! Surface roughness
REAL :: veg (nx,ny) ! Vegetation fraction
REAL :: tsoil (nx,ny,nzsoil,0:nstyps) ! Soil temperature (K)
REAL :: qsoil (nx,ny,nzsoil,0:nstyps) ! Soil moisture (m**3/m**3)
REAL :: wetcanp(nx,ny,0:nstyps) ! Canopy water amount
REAL :: snowdpth(nx,ny) ! Snow depth (m)
REAL :: raing(nx,ny) ! Grid supersaturation rain
REAL :: rainc(nx,ny) ! Cumulus convective rain
REAL :: prcrate(nx,ny,4) ! precipitation rates (kg/(m**2*s))
! prcrate(1,1,1) = total precip. rate
! prcrate(1,1,2) = grid scale precip. rate
! prcrate(1,1,3) = cumulative precip. rate
! prcrate(1,1,4) = microphysics precip. rate
REAL :: radfrc(nx,ny,nz) ! Radiation forcing (K/s)
REAL :: radsw (nx,ny) ! Solar radiation reaching the surface
REAL :: rnflx (nx,ny) ! Net radiation flux absorbed by surface
REAL :: radswnet(nx,ny) ! Net shortwave radiation
REAL :: radlwin(nx,ny) ! Incoming longwave radiation
REAL :: usflx (nx,ny) ! Surface flux of u-momentum (kg/(m*s**2))
REAL :: vsflx (nx,ny) ! Surface flux of v-momentum (kg/(m*s**2))
REAL :: ptsflx(nx,ny) ! Surface heat flux (K*kg/(m*s**2))
REAL :: qvsflx(nx,ny) ! Surface moisture flux (kg/(m**2*s))
REAL :: tem1 (nx,ny,nz) ! Temporary work array
INTEGER (KIND=selected_int_kind(4)), ALLOCATABLE :: itmp(:,:,:) ! Temporary array
REAL, ALLOCATABLE :: hmax(:), hmin(:) ! Temporary array
REAL, ALLOCATABLE :: hmaxsoil(:), hminsoil(:) ! Temporary array
INTEGER (KIND=selected_int_kind(4)), ALLOCATABLE :: itmpsoil(:,:,:,:)
!
! 06/28/2002 Zuwen He
!
! Create a tem2 which will be used when dump data in previous version.
! It may be better to pass the array from the argument list
! in the future.
!
REAL, allocatable :: tem2(:,:,:) ! Temporary array
INTEGER :: ireturn
REAL :: time
!-----------------------------------------------------------------------
!
! Parameters describing routine that wrote the gridded data
!
!-----------------------------------------------------------------------
!
! 06/28/2002 Zuwen He
!
! fmtver??: to label each data a version.
! intver??: an integer to allow faster comparison than fmtver??,
! which are strings.
!
! Verion 5.00: significant change in soil variables since version 4.10.
!
!-----------------------------------------------------------------------
CHARACTER (LEN=40) :: fmtver410,fmtver500
INTEGER :: intver,intver410,intver500
PARAMETER (fmtver410='004.10 HDF4 Coded Data',intver410=410)
PARAMETER (fmtver500='005.00 HDF4 Coded Data',intver500=500)
CHARACTER (LEN=40) :: fmtverin
CHARACTER (LEN=10) :: tmunit
!-----------------------------------------------------------------------
!
! Misc. local variables
!
!-----------------------------------------------------------------------
INTEGER :: lchanl
PARAMETER (lchanl=6) ! Channel number for formatted printing.
INTEGER :: i,j,k,is
INTEGER :: nxin,nyin,nzin,nzsoilin
INTEGER :: bgrdin,bbasin,bvarin,bicein,btrbin,btkein
INTEGER :: istat, sd_id
INTEGER :: nstyp1,nstypin
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'globcst.inc'
INCLUDE 'grid.inc' ! Grid parameters
INCLUDE 'indtflg.inc'
INCLUDE 'alloc.inc' ! allocation parameters & declarations
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
IF(myproc == 0) &
WRITE(*,*) 'HDFREAD: Reading HDF file: ', trim(filename)
ALLOCATE (itmp(nx,ny,nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating itmp, returning"
RETURN
END IF
ALLOCATE (hmax(nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating hmax, returning"
RETURN
END IF
ALLOCATE (hmin(nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating hmin, returning"
RETURN
END IF
ALLOCATE (itmpsoil(nx,ny,nzsoil,0:nstyps),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating itmp, returning"
RETURN
END IF
ALLOCATE (hmaxsoil(nzsoil),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating hmax, returning"
RETURN
END IF
ALLOCATE (hminsoil(nzsoil),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR allocating hmin, returning"
RETURN
END IF
!-----------------------------------------------------------------------
!
! Read header info
!
!-----------------------------------------------------------------------
CALL hdfopen
(filename,1,sd_id)
IF (sd_id < 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: ERROR opening ", &
trim(filename)," for reading."
GO TO 110
END IF
CALL hdfrdc(sd_id,40,"fmtver",fmtverin,istat)
IF (fmtverin == fmtver410) THEN
intver=intver410
ELSE IF (fmtverin == fmtver500) THEN
intver=intver500
ELSE
IF (myproc == 0) WRITE(6,'(/1x,a,a,a/)') &
'Incoming data format, fmtverin=',fmtverin, &
', not found. The Job stopped.'
CALL arpsstop('arpstop called from HDFREAD. ',1)
END IF
IF (myproc == 0) WRITE(6,'(/1x,a,a/)') &
'Incoming data format, fmtverin=',fmtverin
CALL hdfrdc(sd_id,40,"runname",runname,istat)
CALL hdfrdi(sd_id,"nocmnt",nocmnt,istat)
IF( nocmnt > 0 ) THEN
CALL hdfrdc(sd_id,80*nocmnt,"cmnt",cmnt,istat)
END IF
IF(myproc == 0) &
WRITE(6,'(//'' THE NAME OF THIS RUN IS: '',A//)') trim(runname)
IF(myproc == 0) &
WRITE (6,*) "Comments:"
IF( nocmnt > 0 ) THEN
DO i=1,nocmnt
IF(myproc == 0) &
WRITE(6,'(1x,a)') cmnt(i)
END DO
END IF
IF(myproc == 0) &
WRITE (6,*) " "
CALL hdfrdc(sd_id,10,"tmunit",tmunit,istat)
CALL hdfrdr(sd_id,"time",time,istat)
!-----------------------------------------------------------------------
!
! Get dimensions of data in binary file and check against
! the dimensions passed to HDFREAD
!
!-----------------------------------------------------------------------
CALL hdfrdi(sd_id,"nx",nxin,istat)
CALL hdfrdi(sd_id,"ny",nyin,istat)
CALL hdfrdi(sd_id,"nz",nzin,istat)
IF ( nxin /= nx .OR. nyin /= ny .OR. nzin /= nz ) THEN
WRITE(6,'(1x,a)') ' Dimensions in HDFREAD inconsistent with data.'
WRITE(6,'(1x,a,3I15)') ' Read were: ', nxin, nyin, nzin
WRITE(6,'(1x,a,3I15)') ' Expected: ', nx, ny, nz
WRITE(6,'(1x,a)') ' Program aborted in HDFREAD.'
CALL arpsstop('arpsstop called from HDFREAD due to nxin...',1)
END IF
IF (intver <= intver410) THEN
nzsoilin = 2 ! for versions earlier than 410, it is actually
! a 2 level soil model.
ELSE IF (intver >= intver500) THEN
CALL hdfrdi(sd_id,"nzsoil",nzsoilin,istat)
END IF
IF (nzsoilin /= nzsoil) THEN
IF (intver <= intver410) THEN
IF(myproc == 0) WRITE(6,'(1x,a,a/,2(1x,a/))') &
' The incoming data version is ', fmtverin, &
' In the input file, nzsoil must be set to 2. ', &
' Program aborted in HDFREAD.'
ELSE IF (intver >= intver500) THEN
WRITE(6,'(1x,a)') &
' Dimensions in BINREAD inconsistent with data.'
WRITE(6,'(1x,a,I15)') ' Read were: ', nzsoilin
WRITE(6,'(1x,a,I15)') ' Expected: ', nzsoil
WRITE(6,'(1x,a)') ' Program aborted in HDFREAD.'
WRITE(6,*)' myproc = ',myproc
END IF
CALL arpsstop('arpsstop called from HDFREAD due to nzsoilin...',1)
END IF
!-----------------------------------------------------------------------
!
! Read in flags for different data groups
!
!-----------------------------------------------------------------------
IF ( grdbas == 1 ) THEN ! Read grid and base state arrays
IF(myproc == 0) &
WRITE(lchanl,'(1x,a,f8.1,a,f8.3,a/)') &
'To read grid and base state data at time ', time, &
' secs = ',(time/60.),' mins.'
CALL hdfrdi(sd_id,"grdflg",bgrdin,istat)
CALL hdfrdi(sd_id,"basflg",bbasin,istat)
CALL hdfrdi(sd_id,"varflg",bvarin,istat)
CALL hdfrdi(sd_id,"mstflg",mstin,istat)
CALL hdfrdi(sd_id,"iceflg",bicein,istat)
CALL hdfrdi(sd_id,"trbflg",btrbin,istat)
CALL hdfrdi(sd_id,"landflg",landin,istat)
CALL hdfrdi(sd_id,"totflg",totin,istat)
CALL hdfrdi(sd_id,"tkeflg",btkein,istat)
ELSE ! Normal data reading
IF(myproc == 0) &
WRITE(lchanl,'(1x,a,f8.1,a,f8.3,a/)')'To read data for time:', &
time,' secs = ',(time/60.),' mins.'
CALL hdfrdi(sd_id,"grdflg",grdin,istat)
CALL hdfrdi(sd_id,"basflg",basin,istat)
CALL hdfrdi(sd_id,"varflg",varin,istat)
CALL hdfrdi(sd_id,"mstflg",mstin,istat)
CALL hdfrdi(sd_id,"iceflg",icein,istat)
CALL hdfrdi(sd_id,"trbflg",trbin,istat)
CALL hdfrdi(sd_id,"sfcflg",sfcin,istat)
CALL hdfrdi(sd_id,"rainflg",rainin,istat)
!Don't read landin (land data only written to base file,
!landin flag set there).
!CALL hdfrdi(sd_id,"landflg",landin,istat)
CALL hdfrdi(sd_id,"totflg",totin,istat)
CALL hdfrdi(sd_id,"tkeflg",tkein,istat)
END IF
CALL hdfrdi(sd_id,"nstyp",nstyp1,istat)
IF ( nstyp1 < 1 ) THEN
nstyp1 = 1
END IF
IF (nstyp1 > nstyp) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFREAD: WARNING, nstyp in file (",nstyp1, &
") greater than that specified in input file (",nstyp, &
"), using only",nstyp
nstypin = nstyp
ELSE
nstypin = nstyp1
ENDIF
CALL hdfrdi(sd_id,"prcflg",prcin,istat)
CALL hdfrdi(sd_id,"radflg",radin,istat)
CALL hdfrdi(sd_id,"flxflg",flxin,istat)
CALL hdfrdi(sd_id,"snowflg",snowin,istat)
CALL hdfrdi(sd_id,"month",month,istat)
CALL hdfrdi(sd_id,"day",day,istat)
CALL hdfrdi(sd_id,"year",year,istat)
CALL hdfrdi(sd_id,"hour",hour,istat)
CALL hdfrdi(sd_id,"minute",minute,istat)
CALL hdfrdi(sd_id,"second",second,istat)
CALL hdfrdr(sd_id,"umove",umove,istat)
CALL hdfrdr(sd_id,"vmove",vmove,istat)
CALL hdfrdr(sd_id,"xgrdorg",xgrdorg,istat)
CALL hdfrdr(sd_id,"ygrdorg",ygrdorg,istat)
CALL hdfrdi(sd_id,"mapproj",mapproj,istat)
CALL hdfrdr(sd_id,"trulat1",trulat1,istat)
CALL hdfrdr(sd_id,"trulat2",trulat2,istat)
CALL hdfrdr(sd_id,"trulon",trulon,istat)
CALL hdfrdr(sd_id,"sclfct",sclfct,istat)
CALL hdfrdr(sd_id,"tstop",tstop,istat)
CALL hdfrdr(sd_id,"thisdmp",thisdmp,istat)
CALL hdfrdr(sd_id,"latitud",latitud,istat)
CALL hdfrdr(sd_id,"ctrlat",ctrlat,istat)
CALL hdfrdr(sd_id,"ctrlon",ctrlon,istat)
!-----------------------------------------------------------------------
!
! Read in x,y and z at grid cell centers (scalar points).
!
!-----------------------------------------------------------------------
IF( grdin == 1 .OR. grdbas == 1 ) THEN
CALL hdfrd1d(sd_id,"x",nx,x,istat)
IF (istat /= 0) GO TO 110
CALL hdfrd1d(sd_id,"y",ny,y,istat)
IF (istat /= 0) GO TO 110
CALL hdfrd1d(sd_id,"z",nz,z,istat)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"zp",nx,ny,nz,zp,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
IF (intver <= intver410) THEN
!
! nzsoil must equal to 2, 06/28/2002, Zuwen
! assume zpsoil(,,2) is one meter below the surface.
!
DO j=1,ny
DO i=1,nx
zpsoil(i,j,1)=zp(i,j,2)
zpsoil(i,j,2)=zpsoil(i,j,1)-1.
END DO
END DO
IF (myproc == 0) THEN
WRITE(6,*) ' Assign zpsoil. '
WRITE(6,*) ' Assume zpsoil(,,1) is zp(,,2). '
WRITE(6,*) ' Assume zpsoil(,,2) is zp(,,2)-1. '
END IF
ELSE IF (intver >= intver500) THEN
CALL hdfrd3d(sd_id,"zpsoil",nx,ny,nzsoil,zpsoil,istat, &
itmpsoil,hmaxsoil,hminsoil)
END IF
IF (istat /= 0) GO TO 110
END IF ! grdin
!-----------------------------------------------------------------------
!
! Read in base state fields
!
!-----------------------------------------------------------------------
IF( basin == 1 .OR. grdbas == 1 ) THEN
CALL hdfrd3d(sd_id,"ubar",nx,ny,nz,ubar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"vbar",nx,ny,nz,vbar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"wbar",nx,ny,nz,wbar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"ptbar",nx,ny,nz,ptbar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"pbar",nx,ny,nz,pbar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
IF( mstin == 1) THEN
CALL hdfrd3d(sd_id,"qvbar",nx,ny,nz,qvbar,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
END IF
IF (landin == 1) THEN
CALL hdfrd3di(sd_id,"soiltyp",nx,ny,nstypin,soiltyp(1,1,1),istat)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"stypfrct",nx,ny,nstypin, &
stypfrct(1,1,1),istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL fix_stypfrct_nstyp(nx,ny,nstyp1,nstyp,stypfrct)
CALL hdfrd2di(sd_id,"vegtyp",nx,ny,vegtyp,istat)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"lai",nx,ny,lai,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"roufns",nx,ny,roufns,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"veg",nx,ny,veg,istat,itmp)
IF (istat /= 0) GO TO 110
END IF
END IF
IF( grdbas == 1 ) GO TO 930
IF( varin == 1 ) THEN
IF ( totin == 0 ) THEN
!-----------------------------------------------------------------------
!
! Read in perturbations from history dump
!
!-----------------------------------------------------------------------
CALL hdfrd3d(sd_id,"uprt",nx,ny,nz,uprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"vprt",nx,ny,nz,vprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"wprt",nx,ny,nz,wprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
!-----------------------------------------------------------------------
!
! Read in scalars
!
!-----------------------------------------------------------------------
CALL hdfrd3d(sd_id,"ptprt",nx,ny,nz,ptprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"pprt",nx,ny,nz,pprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
ELSE
!-----------------------------------------------------------------------
!
! Read in total values of variables from history dump
!
!-----------------------------------------------------------------------
CALL hdfrd3d(sd_id,"u",nx,ny,nz,uprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx
uprt(i,j,k) = uprt(i,j,k) - ubar(i,j,k)
END DO
END DO
END DO
CALL hdfrd3d(sd_id,"v",nx,ny,nz,vprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO k=1,nz-1
DO j=1,ny
DO i=1,nx-1
vprt(i,j,k) = vprt(i,j,k) - vbar(i,j,k)
END DO
END DO
END DO
CALL hdfrd3d(sd_id,"w",nx,ny,nz,wprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"pt",nx,ny,nz,ptprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
ptprt(i,j,k) = ptprt(i,j,k) - ptbar(i,j,k)
END DO
END DO
END DO
CALL hdfrd3d(sd_id,"p",nx,ny,nz,pprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
pprt(i,j,k) = pprt(i,j,k) - pbar(i,j,k)
END DO
END DO
END DO
END IF
END IF
!-----------------------------------------------------------------------
!
! Read in moisture variables
!
!-----------------------------------------------------------------------
IF( mstin == 1 ) THEN
IF ( totin == 0 ) THEN
CALL hdfrd3d(sd_id,"qvprt",nx,ny,nz,qvprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
ELSE
CALL hdfrd3d(sd_id,"qv",nx,ny,nz,qvprt,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
qvprt(i,j,k) = qvprt(i,j,k) - qvbar(i,j,k)
END DO
END DO
END DO
END IF
CALL hdfrd3d(sd_id,"qc",nx,ny,nz,qc,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"qr",nx,ny,nz,qr,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
IF( rainin == 1 ) THEN
CALL hdfrd2d(sd_id,"raing",nx,ny,raing,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"rainc",nx,ny,rainc,istat,itmp)
IF (istat /= 0) GO TO 110
END IF
IF( prcin == 1 ) THEN
CALL hdfrd2d(sd_id,"prcrate1",nx,ny,prcrate(1,1,1),istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"prcrate2",nx,ny,prcrate(1,1,2),istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"prcrate3",nx,ny,prcrate(1,1,3),istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"prcrate4",nx,ny,prcrate(1,1,4),istat,itmp)
IF (istat /= 0) GO TO 110
END IF
IF( icein == 1 ) THEN
CALL hdfrd3d(sd_id,"qi",nx,ny,nz,qi,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"qs",nx,ny,nz,qs,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"qh",nx,ny,nz,qh,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
END IF
END IF
IF( tkein == 1 ) THEN
CALL hdfrd3d(sd_id,"tke",nx,ny,nz,tke,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
END IF
IF( trbin == 1 ) THEN
CALL hdfrd3d(sd_id,"kmh",nx,ny,nz,kmh,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd3d(sd_id,"kmv",nx,ny,nz,kmv,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
END IF
IF( sfcin == 1 ) THEN
IF (intver <= intver410) THEN
ALLOCATE (tem2(nx,ny,0:nstyps),stat=istat)
tem2=0.
CALL hdfrd3d(sd_id,"tsfc",nx,ny,nstypin+1,tem2,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO is=0,nstypin
DO j=1,ny
DO i=1,nx
tsoil(i,j,1,is)=tem2(i,j,is)
END DO
END DO
END DO
tem2=0.
CALL hdfrd3d(sd_id,"tsoil",nx,ny,nstypin+1,tem2,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO is=0,nstypin
DO j=1,ny
DO i=1,nx
tsoil(i,j,2,is)=tem2(i,j,is)
END DO
END DO
END DO
tem2=0.
CALL hdfrd3d(sd_id,"wetsfc",nx,ny,nstypin+1,tem2,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO is=0,nstypin
DO j=1,ny
DO i=1,nx
qsoil(i,j,1,is)=tem2(i,j,is)
END DO
END DO
END DO
tem2=0.
CALL hdfrd3d(sd_id,"wetdp",nx,ny,nstypin+1,tem2,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
DO is=0,nstypin
DO j=1,ny
DO i=1,nx
qsoil(i,j,2,is)=tem2(i,j,is)
END DO
END DO
END DO
DEALLOCATE (tem2,stat=istat)
ELSE IF (intver >= intver500) THEN
CALL hdfrd4d(sd_id,"tsoil",nx,ny,nzsoil,nstypin+1,tsoil,istat, &
itmpsoil,hmaxsoil,hminsoil)
IF (istat /= 0) GO TO 110
CALL hdfrd4d(sd_id,"qsoil",nx,ny,nzsoil,nstypin+1,qsoil,istat, &
itmpsoil,hmaxsoil,hminsoil)
IF (istat /= 0) GO TO 110
END IF
CALL hdfrd3d(sd_id,"wetcanp",nx,ny,nstypin+1,wetcanp,istat,itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL fix_soil_nstyp(nx,ny,nzsoil,nstyp1,nstyp,tsoil,qsoil,wetcanp)
IF (snowin == 1) THEN
CALL hdfrd2d(sd_id,"snowdpth",nx,ny,snowdpth,istat,itmp)
IF (istat /= 0) GO TO 110
END IF
END IF
IF( radin == 1 ) THEN
CALL hdfrd3d(sd_id,"radfrc",nx,ny,nz,radfrc,istat, &
itmp,hmax,hmin)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"radsw",nx,ny,radsw,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"rnflx",nx,ny,rnflx,istat,itmp)
IF (istat /= 0) GO TO 110
IF (intver <= intver410) THEN
radswnet=0.
radlwin=0.
ELSE IF (intver >= intver500) THEN
CALL hdfrd2d(sd_id,"radswnet",nx,ny,radswnet,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"radlwin",nx,ny,radlwin,istat,itmp)
IF (istat /= 0) GO TO 110
END IF
END IF
IF( flxin == 1 ) THEN
CALL hdfrd2d(sd_id,"usflx",nx,ny,usflx,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"vsflx",nx,ny,vsflx,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"ptsflx",nx,ny,ptsflx,istat,itmp)
IF (istat /= 0) GO TO 110
CALL hdfrd2d(sd_id,"qvsflx",nx,ny,qvsflx,istat,itmp)
IF (istat /= 0) GO TO 110
END IF
!-----------------------------------------------------------------------
!
! Friendly exit message
!
!-----------------------------------------------------------------------
930 CONTINUE
IF(myproc == 0) &
WRITE(6,'(/a,F8.1,a/)') &
' Data at time=', time/60,' (min) were successfully read.'
ireturn = 0
GO TO 130
!-----------------------------------------------------------------------
!
! Error during read
!
!-----------------------------------------------------------------------
110 CONTINUE
IF(myproc == 0) &
WRITE(6,'(/a/)') ' Error reading data in HDFREAD'
ireturn=1
GO TO 130
!-----------------------------------------------------------------------
!
! End-of-file during read
!
!-----------------------------------------------------------------------
! 120 CONTINUE
! IF(myproc == 0) &
! WRITE(6,'(/a/)') ' End of file reached in HDFREAD'
! ireturn=2
130 CONTINUE
!tmp istat = sfendacc(sd_id) ! is this necessary?
CALL hdfclose(sd_id,istat)
IF (ireturn == 0) THEN
IF (istat == 0) THEN
IF(myproc == 0) &
WRITE(*,*) "HDFREAD: Successfully read ", trim(filename)
ELSE
IF(myproc == 0) &
WRITE(*,*) "HDFREAD: ERROR (status=", istat, ") closing ", trim(filename)
END IF
END IF
DEALLOCATE (itmp,stat=istat)
DEALLOCATE (hmax,stat=istat)
DEALLOCATE (hmin,stat=istat)
DEALLOCATE (itmpsoil,stat=istat)
DEALLOCATE (hmaxsoil,stat=istat)
DEALLOCATE (hminsoil,stat=istat)
RETURN
END SUBROUTINE hdfread
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFDUMP ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfdump(nx,ny,nz,nzsoil,nstyps, filename , grdbas, & 1,169
u,v,w,ptprt,pprt,qv,qc,qr,qi,qs,qh,tke,kmh,kmv, &
ubar,vbar,ptbar,pbar,rhobar,qvbar, &
x,y,z,zp,zpsoil, &
soiltyp,stypfrct,vegtyp,lai,roufns,veg, &
tsoil,qsoil,wetcanp,snowdpth, &
raing,rainc,prcrate, &
radfrc,radsw,rnflx,radswnet,radlwin, &
usflx,vsflx,ptsflx,qvsflx, &
tem1)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Produces a history data file "filename" in the NCSA HDF4 format by
! calling HDF library subroutines.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
! 05/15/2002 (J. Brotzge)
! Added to allow for multiple soil schemes
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
! nzsoil Number of grid points in the vertical
!
! filename File name of history dump data.
! grdbas Flag indicating if this is a call for the data dump
! of grid and base state arrays only. If so, grdbas=1.
!
! u x component of velocity at a given time level (m/s)
! v y component of velocity at a given time level (m/s)
! w Vertical component of Cartesian velocity at a given
! time level (m/s)
! ptprt Perturbation potential temperature at a given time
! level (K)
! pprt Perturbation pressure at a given time level (Pascal)
! qv Water vapor specific humidity at a given time level (kg/kg)
! qc Cloud water mixing ratio at a given time level (kg/kg)
! qr Rainwater mixing ratio at a given time level (kg/kg)
! qi Cloud ice mixing ratio at a given time level (kg/kg)
! qs Snow mixing ratio at a given time level (kg/kg)
! qh Hail mixing ratio at a given time level (kg/kg)
! tke Turbulent Kinetic Energy ((m/s)**2)
!
! kmh Horizontal turb. mixing coef. for momentum ( m**2/s )
! kmv Vertical turb. mixing coef. for momentum ( m**2/s )
!
! ubar Base state zonal velocity component (m/s)
! vbar Base state meridional velocity component (m/s)
! ptbar Base state potential temperature (K)
! pbar Base state pressure (Pascal)
! rhobar Base state density (kg/m**3)
! qvbar Base state water vapor specific humidity (kg/kg)
!
! x x coordinate of grid points in physical/comp. space (m)
! y y coordinate of grid points in physical/comp. space (m)
! z z coordinate of grid points in computational space (m)
! zp Vertical coordinate of grid points in physical space (m)
! zpsoil Vertical coordinate of grid points in the soil (m)
!
! soiltyp Soil type
! vegtyp Vegetation type
! lai Leaf Area Index
! roufns Surface roughness
! veg Vegetation fraction
!
! tsoil Soil temperature (K)
! qsoil Soil moisture (m**3/m**3)
! wetcanp Canopy water amount
!
! raing Grid supersaturation rain
! rainc Cumulus convective rain
! prcrate Precipitation rates
!
! radfrc Radiation forcing (K/s)
! radsw Solar radiation reaching the surface
! rnflx Net radiation flux absorbed by surface
! radswnet Net shortwave radiation
! radlwin Incoming longwave radiation
!
! usflx Surface flux of u-momentum (kg/(m*s**2))
! vsflx Surface flux of v-momentum (kg/(m*s**2))
! ptsflx Surface heat flux (K*kg/(m**2 * s ))
! qvsflx Surface moisture flux of (kg/(m**2 * s))
!
! OUTPUT:
!
! None.
!
! WORK ARRAY:
!
! tem1 work array.
! tem2 work array.
!
!-----------------------------------------------------------------------
!
! The following parameters are passed into this subroutine through
! a common block in globcst.inc, and they determine which
! variables are output.
!
! grdout =0 or 1. If grdout=0, grid variables are not dumped.
! basout =0 or 1. If basout=0, base state variables are not dumped.
! varout =0 or 1. If varout=0, perturbation variables are not dumped.
! mstout =0 or 1. If mstout=0, water variables are not dumped.
! rainout=0 or 1. If rainout=0, rain variables are not dumped.
! prcout =0 or 1. If prcout=0, precipitation rates are not dumped.
! iceout =0 or 1. If iceout=0, qi, qs and qh are not dumped.
! tkeout =0 or 1. If tkeout=0, tke is not dumped.
! trbout =0 or 1. If trbout=0, the eddy viscosity km is not dumped.
! radout =0 or 1. If radout=0, the radiation arrays are not dumped.
! flxout =0 or 1. If flxout=0, the surface fluxes are not dumped.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz ! Number of grid points in 3 directions
INTEGER :: nzsoil ! Number of grid points in the soil
CHARACTER (LEN=*) :: filename
INTEGER :: grdbas ! If this is a grid/base state dump
REAL :: u (nx,ny,nz) ! Total u-velocity (m/s)
REAL :: v (nx,ny,nz) ! Total v-velocity (m/s)
REAL :: w (nx,ny,nz) ! Total w-velocity (m/s)
REAL :: ptprt (nx,ny,nz) ! Perturbation potential temperature (K)
REAL :: pprt (nx,ny,nz) ! Perturbation pressure (Pascal)
REAL :: qv (nx,ny,nz) ! Water vapor specific humidity (kg/kg)
REAL :: qc (nx,ny,nz) ! Cloud water mixing ratio (kg/kg)
REAL :: qr (nx,ny,nz) ! Rain water mixing ratio (kg/kg)
REAL :: qi (nx,ny,nz) ! Cloud ice mixing ratio (kg/kg)
REAL :: qs (nx,ny,nz) ! Snow mixing ratio (kg/kg)
REAL :: qh (nx,ny,nz) ! Hail mixing ratio (kg/kg)
REAL :: tke (nx,ny,nz) ! Turbulent Kinetic Energy ((m/s)**2)
REAL :: kmh (nx,ny,nz) ! Horizontal turb. mixing coef. for
! momentum. ( m**2/s )
REAL :: kmv (nx,ny,nz) ! Vertical turb. mixing coef. for
! momentum. ( m**2/s )
REAL :: ubar (nx,ny,nz) ! Base state u-velocity (m/s)
REAL :: vbar (nx,ny,nz) ! Base state v-velocity (m/s)
REAL :: ptbar (nx,ny,nz) ! Base state potential temperature (K)
REAL :: pbar (nx,ny,nz) ! Base state pressure (Pascal)
REAL :: rhobar(nx,ny,nz) ! Base state air density (kg/m**3)
REAL :: qvbar (nx,ny,nz) ! Base state water vapor specific humidity
! (kg/kg)
REAL :: x (nx) ! The x-coord. of the physical and
! computational grid. Defined at u-point.
REAL :: y (ny) ! The y-coord. of the physical and
! computational grid. Defined at v-point.
REAL :: z (nz) ! The z-coord. of the computational grid.
! Defined at w-point on the staggered grid.
REAL :: zp (nx,ny,nz) ! The physical height coordinate defined at
! w-point of the staggered grid.
REAL :: zpsoil(nx,ny,nzsoil) ! The physical height coordinate defined at
! w-point of the soil.
INTEGER :: nstyps ! Number of soil type
INTEGER :: soiltyp (nx,ny,nstyps) ! Soil type
REAL :: stypfrct(nx,ny,nstyps) ! Fraction of soil types
INTEGER :: vegtyp (nx,ny) ! Vegetation type
REAL :: lai (nx,ny) ! Leaf Area Index
REAL :: roufns (nx,ny) ! Surface roughness
REAL :: veg (nx,ny) ! Vegetation fraction
REAL :: tsoil (nx,ny,nzsoil,0:nstyps) ! Soil temperature (K)
REAL :: qsoil (nx,ny,nzsoil,0:nstyps) ! Soil moisture (m**3/m**3)
REAL :: wetcanp(nx,ny,0:nstyps) ! Canopy water amount
REAL :: snowdpth(nx,ny) ! Snow depth (m)
REAL :: raing(nx,ny) ! Grid supersaturation rain
REAL :: rainc(nx,ny) ! Cumulus convective rain
REAL :: prcrate(nx,ny,4) ! precipitation rates (kg/(m**2*s))
! prcrate(1,1,1) = total precip. rate
! prcrate(1,1,2) = grid scale precip. rate
! prcrate(1,1,3) = cumulative precip. rate
! prcrate(1,1,4) = microphysics precip. rate
REAL :: radfrc(nx,ny,nz) ! Radiation forcing (K/s)
REAL :: radsw (nx,ny) ! Solar radiation reaching the surface
REAL :: rnflx (nx,ny) ! Net radiation flux absorbed by surface
REAL :: radswnet(nx,ny) ! Net shortwave radiation
REAL :: radlwin(nx,ny) ! Incoming longwave radiation
REAL :: usflx (nx,ny) ! Surface flux of u-momentum (kg/(m*s**2))
REAL :: vsflx (nx,ny) ! Surface flux of v-momentum (kg/(m*s**2))
REAL :: ptsflx(nx,ny) ! Surface heat flux (K*kg/(m*s**2))
REAL :: qvsflx(nx,ny) ! Surface moisture flux (kg/(m**2*s))
REAL :: tem1 (nx,ny,nz) ! Temporary work array
INTEGER(2), allocatable :: itmp(:,:,:) ! Temporary array
REAL, allocatable :: hmax(:), hmin(:) ! Temporary array
INTEGER(2), allocatable :: itmpsoil(:,:,:,:) ! Temporary array
REAL, allocatable :: hmaxsoil(:), hminsoil(:)! Temporary array
!
! 06/28/2002 Zuwen He
!
! Create a tem2 which will be used when dump data in previous version.
! It may be better to pass the array from the argument list
! in the future.
!
REAL, allocatable :: tem2(:,:,:) ! Temporary array
REAL :: dx_out,dy_out
!-----------------------------------------------------------------------
!
! Parameters describing this routine
!
!-----------------------------------------------------------------------
! 06/28/2002 Zuwen He
!
! fmtver??: to label each data a version.
! intver??: an integer to allow faster comparison than fmtver??,
! which are strings.
!
! Verion 5.00: significant change in soil variables since version 4.10.
!
CHARACTER (LEN=40) :: fmtver,fmtver410,fmtver500
INTEGER :: intver,intver410,intver500
PARAMETER (fmtver410='004.10 HDF4 Coded Data',intver410=410)
PARAMETER (fmtver500='005.00 HDF4 Coded Data',intver500=500)
CHARACTER (LEN=10) :: tmunit
PARAMETER (tmunit='seconds ')
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: i,j,k,is
INTEGER :: nstypout
INTEGER :: istat, sd_id
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'globcst.inc'
INCLUDE 'grid.inc' ! Grid parameters
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
IF (hdfcompr > 3) THEN
ALLOCATE (itmp(nx,ny,nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating itmp, returning"
RETURN
END IF
ALLOCATE (hmax(nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating hmax, returning"
RETURN
END IF
ALLOCATE (hmin(nz),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating hmin, returning"
RETURN
END IF
ALLOCATE (itmpsoil(nx,ny,nzsoil,0:nstyps),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating itmp, returning"
RETURN
END IF
ALLOCATE (hmaxsoil(nzsoil),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating hmax, returning"
RETURN
END IF
ALLOCATE (hminsoil(nzsoil),stat=istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR allocating hmin, returning"
RETURN
END IF
END IF
IF (myproc == 0) &
WRITE(6,'(1x,a,f13.3,a,a/)') &
'Writing HDF4 data at time=', curtim,' into file ',filename
!-----------------------------------------------------------------------
!
! Create the HDF4 file.
!
!-----------------------------------------------------------------------
CALL hdfopen(filename,2,sd_id)
IF (sd_id < 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR creating HDF4 file: ", &
trim(filename)
GO TO 600
END IF
intver = intver500 ! for the time being, in the future, we will
! allow to dump data in the different version
! intver will be assigned from input file
IF (intver == intver410) THEN
fmtver=fmtver410
ELSE IF (intver == intver500) THEN
fmtver=fmtver500
ELSE
IF (myproc == 0) WRITE(6,'(/1x,a,i10,a/)') &
' Data format, intver=',intver, &
', not found. The Job stopped.'
CALL arpsstop('arpstop called from HDFDUMP. ',1)
END IF
CALL hdfwrtc(sd_id, 40, 'fmtver', fmtver, istat)
CALL hdfwrtc(sd_id, 40, 'runname', runname, istat)
CALL hdfwrti(sd_id, 'nocmnt', nocmnt, istat)
IF( nocmnt > 0 ) THEN
CALL hdfwrtc(sd_id, 80*nocmnt, 'cmnt', cmnt, istat)
END IF
CALL hdfwrtc(sd_id, 7, 'tmunit', 'seconds', istat)
CALL hdfwrtr(sd_id, 'time', curtim, istat)
CALL hdfwrti(sd_id, 'nx', nx, istat)
CALL hdfwrti(sd_id, 'ny', ny, istat)
CALL hdfwrti(sd_id, 'nz', nz, istat)
IF (intver >= intver500) THEN
CALL hdfwrti(sd_id, 'nzsoil', nzsoil, istat)
END IF
IF( grdbas == 1 ) THEN
CALL hdfwrti(sd_id, 'grdflg', 1, istat)
CALL hdfwrti(sd_id, 'basflg', 1, istat)
CALL hdfwrti(sd_id, 'varflg', 0, istat)
CALL hdfwrti(sd_id, 'mstflg', 1, istat)
CALL hdfwrti(sd_id, 'iceflg', 0, istat)
CALL hdfwrti(sd_id, 'trbflg', 0, istat)
CALL hdfwrti(sd_id, 'sfcflg', 0, istat)
CALL hdfwrti(sd_id, 'rainflg', 0, istat)
CALL hdfwrti(sd_id, 'landflg', 1, istat)
CALL hdfwrti(sd_id, 'totflg', 1, istat)
CALL hdfwrti(sd_id, 'tkeflg', 0, istat)
ELSE
CALL hdfwrti(sd_id, 'grdflg', grdout, istat)
CALL hdfwrti(sd_id, 'basflg', basout, istat)
CALL hdfwrti(sd_id, 'varflg', varout, istat)
CALL hdfwrti(sd_id, 'mstflg', mstout, istat)
CALL hdfwrti(sd_id, 'iceflg', iceout, istat)
CALL hdfwrti(sd_id, 'trbflg', trbout, istat)
CALL hdfwrti(sd_id, 'sfcflg', sfcout, istat)
CALL hdfwrti(sd_id, 'rainflg', rainout, istat)
CALL hdfwrti(sd_id, 'landflg', landout*basout, istat)
CALL hdfwrti(sd_id, 'totflg', totout, istat)
CALL hdfwrti(sd_id, 'tkeflg', tkeout, istat)
END IF
nstypout = max(1,nstyp)
CALL hdfwrti(sd_id, 'nstyp', nstypout, istat)
CALL hdfwrti(sd_id, 'prcflg', prcout, istat)
CALL hdfwrti(sd_id, 'radflg', radout, istat)
CALL hdfwrti(sd_id, 'flxflg', flxout, istat)
CALL hdfwrti(sd_id, 'snowflg', snowout, istat)
CALL hdfwrti(sd_id, 'day', day, istat)
CALL hdfwrti(sd_id, 'year', year, istat)
CALL hdfwrti(sd_id, 'month', month, istat)
CALL hdfwrti(sd_id, 'hour', hour, istat)
CALL hdfwrti(sd_id, 'minute', minute, istat)
CALL hdfwrti(sd_id, 'second', second, istat)
CALL hdfwrtr(sd_id, 'umove', umove, istat)
CALL hdfwrtr(sd_id, 'vmove', vmove, istat)
CALL hdfwrtr(sd_id, 'xgrdorg', xgrdorg, istat)
CALL hdfwrtr(sd_id, 'ygrdorg', ygrdorg, istat)
CALL hdfwrti(sd_id, 'mapproj', mapproj, istat)
CALL hdfwrtr(sd_id, 'trulat1', trulat1, istat)
CALL hdfwrtr(sd_id, 'trulat2', trulat2, istat)
CALL hdfwrtr(sd_id, 'trulon', trulon, istat)
CALL hdfwrtr(sd_id, 'sclfct', sclfct, istat)
CALL hdfwrtr(sd_id, 'tstop', tstop, istat)
CALL hdfwrtr(sd_id, 'thisdmp', thisdmp, istat)
CALL hdfwrtr(sd_id, 'latitud', latitud, istat)
CALL hdfwrtr(sd_id, 'ctrlat', ctrlat, istat)
CALL hdfwrtr(sd_id, 'ctrlon', ctrlon, istat)
dx_out = x(2) - x(1)
dy_out = y(2) - y(1)
CALL hdfwrtr(sd_id, 'dx', dx_out, istat)
CALL hdfwrtr(sd_id, 'dy', dy_out, istat)
!-----------------------------------------------------------------------
!
! If grdout=1 or grdbas=1, write out grid variables
!
!-----------------------------------------------------------------------
IF(grdout == 1 .OR. grdbas == 1 ) THEN
CALL hdfwrt1d(x,nx,sd_id,'x','x coordinate','m')
CALL hdfwrt1d(y,ny,sd_id,'y','y coordinate','m')
CALL hdfwrt1d(z,nz,sd_id,'z','z coordinate','m')
CALL hdfwrt3d(zp,nx,ny,nz,sd_id,0,hdfcompr, &
'zp','Physical height coordinate','m', &
itmp,hmax,hmin)
IF (intver >= intver500) THEN
CALL hdfwrt3d(zpsoil,nx,ny,nzsoil,sd_id,0,hdfcompr, &
'zpsoil','Physical height coordinate (soil)','m', &
itmpsoil,hmaxsoil,hminsoil)
END IF
END IF
!-----------------------------------------------------------------------
!
! If basout=1, write out base state variables.
!
!-----------------------------------------------------------------------
IF(basout == 1 .OR. grdbas == 1 ) THEN
CALL edgfill(ubar,1,nx,1,nx, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(ubar,nx,ny,nz,sd_id,1,hdfcompr, &
'ubar','Base state u-velocity','m/s', &
itmp,hmax,hmin)
CALL edgfill(vbar,1,nx,1,nx-1, 1,ny,1,ny, 1,nz,1,nz-1)
CALL hdfwrt3d(vbar,nx,ny,nz,sd_id,2,hdfcompr, &
'vbar','Base state v-velocity','m/s', &
itmp,hmax,hmin)
DO k=1,nz
DO j=1,ny
DO i=1,nx
tem1(i,j,k) = 0.0
END DO
END DO
END DO
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,3,hdfcompr, &
'wbar','Base state w-velocity','m/s', &
itmp,hmax,hmin)
CALL edgfill(ptbar,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(ptbar,nx,ny,nz,sd_id,0,hdfcompr, &
'ptbar','Base state potential temperature','K', &
itmp,hmax,hmin)
CALL edgfill(pbar,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(pbar,nx,ny,nz,sd_id,0,hdfcompr, &
'pbar','Base state pressure','Pascal', &
itmp,hmax,hmin)
IF(mstout == 1) THEN
CALL edgfill(qvbar,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qvbar,nx,ny,nz,sd_id,0,hdfcompr, &
'qvbar','Base state water vapor specific humidity','kg/kg', &
itmp,hmax,hmin)
END IF
IF(landout == 1) THEN
CALL iedgfill(soiltyp(1,1,1),1,nx,1,nx-1, 1,ny,1,ny-1, &
1,nstypout,1,nstypout)
CALL hdfwrt3di(soiltyp,nx,ny,nstypout,sd_id,0,0, &
'soiltyp','Soil type','index')
CALL edgfill(stypfrct(1,1,1),1,nx,1,nx-1, 1,ny,1,ny-1, &
1,nstypout,1,nstypout)
CALL hdfwrt3d(stypfrct(1,1,1),nx,ny,nstypout,sd_id,0,hdfcompr, &
'stypfrct','Soil type fractional coverage','fraction', &
itmp,hmax,hmin)
CALL iedgfill(vegtyp ,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2di(vegtyp,nx,ny,sd_id,0,0, &
'vegtyp','Vegetation type','index')
CALL edgfill(lai ,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(lai,nx,ny,sd_id,0,hdfcompr, &
'lai','Leaf Area Index','index',itmp)
CALL edgfill(roufns ,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(roufns,nx,ny,sd_id,0,hdfcompr, &
'roufns','Surface roughness','0-1',itmp)
CALL edgfill(veg ,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(veg,nx,ny,sd_id,0,hdfcompr, &
'veg','Vegetation fraction','fraction',itmp)
END IF
END IF
IF ( grdbas == 1 ) GO TO 600
!-----------------------------------------------------------------------
!
! If varout = 1, Write out uprt, vprt, wprt, ptprt, pprt.
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Write out u,v and w.
!
!-----------------------------------------------------------------------
IF(varout == 1) THEN
IF ( totout == 0 ) THEN
!-----------------------------------------------------------------------
!
! Write out perturbations to history dump
!
!-----------------------------------------------------------------------
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx
tem1(i,j,k)=u(i,j,k)-ubar(i,j,k)
END DO
END DO
END DO
CALL edgfill(tem1,1,nx,1,nx, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,1,hdfcompr, &
'uprt','Perturbation u-velocity','m/s', &
itmp,hmax,hmin)
DO k=1,nz-1
DO j=1,ny
DO i=1,nx-1
tem1(i,j,k)=v(i,j,k)-vbar(i,j,k)
END DO
END DO
END DO
CALL edgfill(tem1,1,nx,1,nx-1, 1,ny,1,ny, 1,nz,1,nz-1)
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,2,hdfcompr, &
'vprt','Perturbation v-velocity','m/s', &
itmp,hmax,hmin)
CALL edgfill(w,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz)
CALL hdfwrt3d(w,nx,ny,nz,sd_id,3,hdfcompr, &
'wprt','Perturbation w-velocity','m/s', &
itmp,hmax,hmin)
!-----------------------------------------------------------------------
!
! Write out scalars
!
!-----------------------------------------------------------------------
CALL edgfill(ptprt,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(ptprt,nx,ny,nz,sd_id,0,hdfcompr, &
'ptprt','Perturbation potential temperature','K', &
itmp,hmax,hmin)
CALL edgfill(pprt,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(pprt,nx,ny,nz,sd_id,0,hdfcompr, &
'pprt','Perturbation pressure','Pascal', &
itmp,hmax,hmin)
ELSE
!-----------------------------------------------------------------------
!
! Write out total values to history dump
!
!-----------------------------------------------------------------------
CALL edgfill(u,1,nx,1,nx, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(u,nx,ny,nz,sd_id,1,hdfcompr, &
'u','u-velocity','m/s', &
itmp,hmax,hmin)
CALL edgfill(v,1,nx,1,nx-1, 1,ny,1,ny, 1,nz,1,nz-1)
CALL hdfwrt3d(v,nx,ny,nz,sd_id,2,hdfcompr, &
'v','v-velocity','m/s', &
itmp,hmax,hmin)
CALL edgfill(w,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz)
CALL hdfwrt3d(w,nx,ny,nz,sd_id,3,hdfcompr, &
'w','w-velocity','m/s', &
itmp,hmax,hmin)
!-----------------------------------------------------------------------
!
! Write out scalars
!
!-----------------------------------------------------------------------
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
tem1(i,j,k) = ptbar(i,j,k) + ptprt(i,j,k)
END DO
END DO
END DO
CALL edgfill(tem1,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,0,hdfcompr, &
'pt','Potential temperature','K', &
itmp,hmax,hmin)
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
tem1(i,j,k) = pbar(i,j,k) + pprt(i,j,k)
END DO
END DO
END DO
CALL edgfill(tem1,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,0,hdfcompr, &
'p','Pressure','Pascal', &
itmp,hmax,hmin)
END IF
END IF ! varout
!-----------------------------------------------------------------------
!
! If mstout = 1, write out moisture scalars.
!
!-----------------------------------------------------------------------
IF(mstout == 1) THEN
IF( totout == 0 ) THEN
!-----------------------------------------------------------------------
!
! Write out perturbation to history dump
!
!-----------------------------------------------------------------------
DO k=1,nz-1
DO j=1,ny-1
DO i=1,nx-1
tem1(i,j,k)=qv(i,j,k)-qvbar(i,j,k)
END DO
END DO
END DO
CALL edgfill(tem1,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(tem1,nx,ny,nz,sd_id,0,hdfcompr, &
'qvprt','Pert. water vapor specific humidity','kg/kg', &
itmp,hmax,hmin)
ELSE
!-----------------------------------------------------------------------
!
! Write out total values to history dump
!
!-----------------------------------------------------------------------
CALL edgfill(qv,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qv,nx,ny,nz,sd_id,0,hdfcompr, &
'qv','Water vapor specific humidity','kg/kg', &
itmp,hmax,hmin)
END IF
CALL edgfill(qc,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qc,nx,ny,nz,sd_id,0,hdfcompr, &
'qc','Cloud water mixing ratio','kg/kg', &
itmp,hmax,hmin)
CALL edgfill(qr,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qr,nx,ny,nz,sd_id,0,hdfcompr, &
'qr','Rain water mixing ratio','kg/kg', &
itmp,hmax,hmin)
IF(rainout == 1) THEN
CALL edgfill(raing, 1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(raing,nx,ny,sd_id,0,hdfcompr, &
'raing','Grid supersaturation rain','mm',itmp)
CALL edgfill(rainc,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(rainc,nx,ny,sd_id,0,hdfcompr, &
'rainc','Cumulus convective rain','mm',itmp)
END IF !rainout
IF ( prcout == 1 ) THEN
CALL edgfill(prcrate,1,nx,1,nx-1, 1,ny,1,ny-1, 1,4,1,4)
CALL hdfwrt2d(prcrate(1,1,1),nx,ny,sd_id,0,hdfcompr, &
'prcrate1','Total precip. rate','kg/(m**2*s)',itmp)
CALL hdfwrt2d(prcrate(1,1,2),nx,ny,sd_id,0,hdfcompr, &
'prcrate2','Grid scale precip. rate','kg/(m**2*s)',itmp)
CALL hdfwrt2d(prcrate(1,1,3),nx,ny,sd_id,0,hdfcompr, &
'prcrate3','Cumulative precip. rate','kg/(m**2*s)',itmp)
CALL hdfwrt2d(prcrate(1,1,4),nx,ny,sd_id,0,hdfcompr, &
'prcrate4','Microphysics precip. rate','kg/(m**2*s)',itmp)
END IF
IF(iceout == 1) THEN
CALL edgfill(qi,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qi,nx,ny,nz,sd_id,0,hdfcompr, &
'qi','Cloud ice mixing ratio','kg/kg', &
itmp,hmax,hmin)
CALL edgfill(qs,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qs,nx,ny,nz,sd_id,0,hdfcompr, &
'qs','Snow mixing ratio','kg/kg', &
itmp,hmax,hmin)
CALL edgfill(qh,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(qh,nx,ny,nz,sd_id,0,hdfcompr, &
'qh','Hail mixing ratio','kg/kg', &
itmp,hmax,hmin)
END IF !iceout
END IF !mstout
!-----------------------------------------------------------------------
!
! If tkeout = 1, write out tke.
!
!-----------------------------------------------------------------------
IF( tkeout == 1 ) THEN
CALL edgfill(tke,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(tke,nx,ny,nz,sd_id,0,hdfcompr, &
'tke','Turbulent Kinetic Energy','(m/s)**2', &
itmp,hmax,hmin)
END IF
!-----------------------------------------------------------------------
!
! If trbout = 1, write out the turbulence parameter, km.
!
!-----------------------------------------------------------------------
IF( trbout == 1 ) THEN
CALL edgfill(kmh,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(kmh,nx,ny,nz,sd_id,0,hdfcompr, &
'kmh','Hori. turb. mixing coef. for momentum','m**2/s', &
itmp,hmax,hmin)
CALL edgfill(kmv,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(kmv,nx,ny,nz,sd_id,0,hdfcompr, &
'kmv','Vert. turb. mixing coef. for momentum','m**2/s', &
itmp,hmax,hmin)
END IF ! trbout
!-----------------------------------------------------------------------
!
! If sfcout = 1, write out the surface variables,
! tsoil, qsoil, and wetcanp.
!
!-----------------------------------------------------------------------
IF( sfcout == 1) THEN
DO is=0,nstypout
CALL edgfill(tsoil(1,1,1,is), 1,nx,1,nx-1, 1,ny,1,ny-1, &
1,nzsoil,1,nzsoil)
CALL edgfill(qsoil(1,1,1,is), 1,nx,1,nx-1, 1,ny,1,ny-1, &
1,nzsoil,1,nzsoil)
CALL edgfill(wetcanp(1,1,is),1,nx,1,nx-1, 1,ny,1,ny-1, &
1,1,1,1)
END DO
IF (intver == intver410) THEN
!
! 06/28/2002 Zuwen He
!
! In version 410, tsoil is the average of tsoil from 2 to nzsoil in later
! version, and wetdp is similar.
!
ALLOCATE (tem2(nx,ny,0:nstyps),stat=istat)
tem2=0.
DO is=0,nstypout
DO j=1,ny
DO i=1,nx
tem2(i,j,is)=tsoil(i,j,1,is)
END DO
END DO
END DO
CALL hdfwrt3d(tem2,nx,ny,nstypout+1,sd_id,0,hdfcompr, &
'tsfc','Surface ground temperature','K', &
itmp,hmax,hmin)
tem2=0.
DO is=0,nstypout
DO k=2,nzsoil
DO j=1,ny
DO i=1,nx
tem2(i,j,is)=tem2(i,j,is)+tsoil(i,j,k,is)
END DO
END DO
END DO
DO j=1,ny
DO i=1,nx
tem1(i,j,is)=tem1(i,j,is)/float(nzsoil-1)
END DO
END DO
END DO
CALL hdfwrt3d(tem1,nx,ny,nstypout+1,sd_id,0,hdfcompr, &
'tsoil','Deep soil temperature','K', &
itmp,hmax,hmin)
tem2=0.
DO is=0,nstypout
DO j=1,ny
DO i=1,nx
tem2(i,j,is)=qsoil(i,j,1,is)
END DO
END DO
END DO
CALL hdfwrt3d(tem2,nx,ny,nstypout+1,sd_id,0,hdfcompr, &
'wetsfc','Surface soil moisture','fraction', &
itmp,hmax,hmin)
tem2=0.
DO is=0,nstypout
DO k=2,nzsoil
DO j=1,ny
DO i=1,nx
tem2(i,j,is)=tem2(i,j,is)+qsoil(i,j,k,is)
END DO
END DO
END DO
DO j=1,ny
DO i=1,nx
tem1(i,j,is)=tem1(i,j,is)/float(nzsoil-1)
END DO
END DO
END DO
CALL hdfwrt3d(tem1,nx,ny,nstypout+1,sd_id,0,hdfcompr, &
'wetdp','Deep soil moisture','fraction', &
itmp,hmax,hmin)
DEALLOCATE (tem2,stat=istat)
ELSE IF (intver == intver500) THEN
CALL hdfwrt4d(tsoil,nx,ny,nzsoil,nstypout+1,sd_id,0,hdfcompr, &
'tsoil','Soil temperature','K', &
itmpsoil,hmaxsoil,hminsoil)
CALL hdfwrt4d(qsoil,nx,ny,nzsoil,nstypout+1,sd_id,0,hdfcompr, &
'qsoil','Soil moisture','fraction', &
itmpsoil,hmaxsoil,hminsoil)
END IF
CALL hdfwrt3d(wetcanp,nx,ny,nstypout+1,sd_id,0,hdfcompr, &
'wetcanp','Canopy water amount','fraction', &
itmp,hmax,hmin)
IF (snowout == 1) THEN
CALL edgfill(snowdpth,1,nx,1,nx-1, 1,ny,1,ny-1,1,1,1,1)
CALL hdfwrt2d(snowdpth,nx,ny,sd_id,0,hdfcompr, &
'snowdpth','Snow depth','m',itmp)
END IF
END IF
!-----------------------------------------------------------------------
!
! If radout = 1, write out the radiation arrays
!
!-----------------------------------------------------------------------
IF( radout == 1 ) THEN
CALL edgfill(radfrc,1,nx,1,nx-1, 1,ny,1,ny-1, 1,nz,1,nz-1)
CALL hdfwrt3d(radfrc,nx,ny,nz,sd_id,0,hdfcompr, &
'radfrc','Radiation forcing','K/s', &
itmp,hmax,hmin)
CALL edgfill(radsw,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(radsw,nx,ny,sd_id,0,hdfcompr, &
'radsw','Solar radiation reaching the surface','W/m**2',itmp)
CALL edgfill(rnflx,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(rnflx,nx,ny,sd_id,0,hdfcompr, &
'rnflx','Net radiation flux absorbed by surface','W/m**2', &
itmp)
IF (intver >= intver500) THEN
CALL edgfill(radswnet,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(radswnet,nx,ny,sd_id,0,hdfcompr, &
'radswnet','Net solar radiation','W/m**2',itmp)
CALL edgfill(radlwin,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(radlwin,nx,ny,sd_id,0,hdfcompr, &
'radlwin','Incoming longwave radiation','W/m**2', &
itmp)
END IF
END IF ! radout
!-----------------------------------------------------------------------
!
! If flxout = 1, write out the surface fluxes
!
!-----------------------------------------------------------------------
IF( flxout == 1 ) THEN
CALL edgfill(usflx,1,nx,1,nx, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(usflx,nx,ny,sd_id,0,hdfcompr, &
'usflx','Surface flux of u-momentum','kg/(m*s**2)',itmp)
CALL edgfill(vsflx,1,nx,1,nx-1, 1,ny,1,ny, 1,1,1,1)
CALL hdfwrt2d(vsflx,nx,ny,sd_id,0,hdfcompr, &
'vsflx','Surface flux of v-momentum','kg/(m*s**2)',itmp)
CALL edgfill(ptsflx,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(ptsflx,nx,ny,sd_id,0,hdfcompr, &
'ptsflx','Surface heat flux','K*kg/(m*s**2)',itmp)
CALL edgfill(qvsflx,1,nx,1,nx-1, 1,ny,1,ny-1, 1,1,1,1)
CALL hdfwrt2d(qvsflx,nx,ny,sd_id,0,hdfcompr, &
'qvsflx','Surface moisture flux','kg/(m**2*s)',itmp)
END IF ! flxout
600 CONTINUE
CALL hdfclose(sd_id,istat)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFDUMP: ERROR on closing file ",trim(filename), &
" (status",istat,")"
END IF
IF (hdfcompr > 3) THEN
DEALLOCATE (itmp,stat=istat)
DEALLOCATE (hmax,stat=istat)
DEALLOCATE (hmin,stat=istat)
DEALLOCATE (itmpsoil,stat=istat)
DEALLOCATE (hmaxsoil,stat=istat)
DEALLOCATE (hminsoil,stat=istat)
END IF
RETURN
END SUBROUTINE hdfdump
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT3D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
SUBROUTINE hdfwrt3d(var,nx,ny,nz,sd_id,stag_dim,hdfcompr, & 81,1
name,comment,units,itmp,hmax,hmin)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 3-D real array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
!
! sd_id HDF id of the output file
!
! stag_dim Dimension of grid staggering (0-none, 1-x, 2-y, 3-z)
!
! hdfcompr Compression flag (0-none)
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
! hmax Used to store maximum values as a function of z
! hmin Used to store minimum values as a function of z
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz
REAL :: var(nx,ny,nz)
INTEGER :: sd_id, stag_dim, hdfcompr
CHARACTER (LEN=*) :: name, comment, units
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny,nz)
REAL :: hmax(nz), hmin(nz)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(3),start(3),stride(3)
INTEGER :: comp_prm(1)
INTEGER :: sds_id
INTEGER :: i,j,k
REAL :: scale
INTEGER :: itmp1
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
start(1) = 0
start(2) = 0
start(3) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
IF(myproc == 0) &
WRITE (*,*) "HDFWRT3D: Writing variable ",trim(name)
IF (hdfcompr <= 3) THEN
sds_id = sfcreate(sd_id, trim(name), dfnt_float32, 3, dims)
ELSE
sds_id = sfcreate(sd_id, trim(name), dfnt_int16, 3, dims)
END IF
comp_prm(1) = 0
IF (hdfcompr == 1 .OR. hdfcompr == 5) THEN ! quick gzip
comp_prm(1) = 1
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 2 .OR. hdfcompr == 6) THEN ! high gzip
comp_prm(1) = 6
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 3 .OR. hdfcompr == 7) THEN ! huffman
comp_prm(1) = 4 ! this may be a problem on a Cray
istat = sfscompress(sds_id, 3, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 3)
END IF
istat = sfscatt(sds_id, 'hdf_comp_prm', dfnt_int32, 1, comp_prm)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
istat = sfsnatt(sds_id, 'stag_dim', dfnt_int32, 1, stag_dim)
!-----------------------------------------------------------------------
!
! If called for, map reals to 16 bit integers
!
!-----------------------------------------------------------------------
IF (hdfcompr > 3) THEN
DO k=1,nz
CALL a3dmax0lcl(var(1,1,k),1,nx,1,nx,1,ny,1,ny,1,1,1,1, &
hmax(k),hmin(k))
IF (ABS(hmax(k)-hmin(k)) < 1.0E-10) hmax(k) = 1.1 * hmin(k) + 1.
scale = 65534.0 / (hmax(k) - hmin(k))
DO j=1,ny
DO i=1,nx
itmp1 = nint(scale * (var(i,j,k) - hmin(k))) - 32767
itmp(i,j,k) = itmp1
END DO
END DO
END DO
istat = sfsnatt(sds_id, 'packed16', dfnt_int32, 1, 1)
istat = sfsnatt(sds_id, 'max', dfnt_float32, nz, hmax)
istat = sfsnatt(sds_id, 'min', dfnt_float32, nz, hmin)
END IF
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
IF (hdfcompr <= 3) THEN
istat = sfwdata(sds_id, start, stride, dims, var)
ELSE
istat = sfwdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT3D: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
WRITE (6,*) "HDFWRT3D: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt3d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT4D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
SUBROUTINE hdfwrt4d(var,nx,ny,nz,nl,sd_id,stag_dim,hdfcompr,& 6,1
name,comment,units,itmp,hmax,hmin)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 4-D real array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Eric Kemp, 23 May 2002, Based on HDFWRT3D
!
! MODIFICATION HISTORY:
!
! 07/02/2002 Zuwen He
!
! Clean up and make it a general 4d hdfwrt.
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
! nl Number of grid points in the 4th dimension
!
! sd_id HDF id of the output file
!
! stag_dim Dimension of grid staggering (0-none, 1-x, 2-y, 3-z)
!
! hdfcompr Compression flag (0-none)
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
! hmax Used to store maximum values as a function of z
! hmin Used to store minimum values as a function of z
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz,nl
REAL :: var(nx,ny,nz,nl)
INTEGER :: sd_id, stag_dim, hdfcompr
CHARACTER (LEN=*) :: name, comment, units
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny,nz,nl)
REAL :: hmax(nz), hmin(nz)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(4),start(4),stride(4)
INTEGER :: comp_prm(1)
INTEGER :: sds_id
INTEGER :: i,j,k,l
REAL :: scale
INTEGER :: itmp1
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
dims(4) = nl
start(1) = 0
start(2) = 0
start(3) = 0
start(4) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
stride(4) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
! IF (myproc == 0) WRITE (*,*) "HDFWRT4D: Writing variable ",trim(name)
IF (hdfcompr <= 3) THEN
sds_id = sfcreate(sd_id, trim(name), dfnt_float32, 4, dims)
ELSE
sds_id = sfcreate(sd_id, trim(name), dfnt_int16, 4, dims)
END IF
comp_prm(1) = 0
IF (hdfcompr == 1 .OR. hdfcompr == 5) THEN ! quick gzip
comp_prm(1) = 1
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 2 .OR. hdfcompr == 6) THEN ! high gzip
comp_prm(1) = 6
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 3 .OR. hdfcompr == 7) THEN ! huffman
comp_prm(1) = 4 ! this may be a problem on a Cray
istat = sfscompress(sds_id, 3, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 3)
END IF
istat = sfscatt(sds_id, 'hdf_comp_prm', dfnt_int32, 1, comp_prm)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
istat = sfsnatt(sds_id, 'stag_dim', dfnt_int32, 1, stag_dim)
!-----------------------------------------------------------------------
!
! If called for, map reals to 16 bit integers
!
!-----------------------------------------------------------------------
IF (hdfcompr > 3) THEN
DO l=1,nl
DO k=1,nz
CALL a3dmax0lcl(var(1,1,k,l),1,nx,1,nx,1,ny,1,ny,1,1,1,1, &
hmax(k),hmin(k))
IF (ABS(hmax(k)-hmin(k)) < 1.0E-10) hmax(k) = 1.1 * hmin(k) + 1.
scale = 65534.0 / (hmax(k) - hmin(k))
DO j=1,ny
DO i=1,nx
itmp1 = nint(scale * (var(i,j,k,l) - hmin(k))) - 32767
itmp(i,j,k,l) = itmp1
END DO
END DO
END DO
END DO
istat = sfsnatt(sds_id, 'packed16', dfnt_int32, 1, 1)
istat = sfsnatt(sds_id, 'max', dfnt_float32, nz, hmax)
istat = sfsnatt(sds_id, 'min', dfnt_float32, nz, hmin)
END IF
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
IF (hdfcompr <= 3) THEN
istat = sfwdata(sds_id, start, stride, dims, var)
ELSE
istat = sfwdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat /= 0) THEN
WRITE (6,*) "HDFWRT4D: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
WRITE (6,*) "HDFWRT4D: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt4d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT3DI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
SUBROUTINE hdfwrt3di(var,nx,ny,nz,sd_id,stag_dim,hdfcompr, & 4
name,comment,units)
!
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 3-D integer array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
!
! sd_id HDF id of the output file
!
! stag_dim Dimension of grid staggering (0-none, 1-x, 2-y, 3-z)
!
! hdfcompr Compression flag (0-none)
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz
INTEGER :: var(nx,ny,nz)
INTEGER :: sd_id, stag_dim, hdfcompr
CHARACTER (LEN=*) :: name, comment, units
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(3),start(3),stride(3)
INTEGER :: sds_id
INTEGER :: comp_prm(1)
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
start(1) = 0
start(2) = 0
start(3) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
IF (myproc == 0) &
WRITE (*,*) "HDFWRT3DI: Writing variable ",trim(name)
sds_id = sfcreate(sd_id, name, dfnt_int32, 3, dims)
comp_prm(1) = 0
IF (hdfcompr == 1 .OR. hdfcompr == 5) THEN ! quick gzip
comp_prm(1) = 1
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 2 .OR. hdfcompr == 6) THEN ! high gzip
comp_prm(1) = 6
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 3 .OR. hdfcompr == 7) THEN ! huffman
comp_prm(1) = 4 ! this may be a problem on a Cray
istat = sfscompress(sds_id, 3, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 3)
END IF
istat = sfscatt(sds_id, 'hdf_comp_prm', dfnt_int32, 1, comp_prm)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
istat = sfsnatt(sds_id, 'stag_dim', dfnt_int32, 1, stag_dim)
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
istat = sfwdata(sds_id, start, stride, dims, var)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT3DI: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT3DI: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt3di
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT2D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
SUBROUTINE hdfwrt2d(var,nx,ny,sd_id,stag_dim,hdfcompr, & 37,1
name,comment,units,itmp)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 2-D real array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
!
! sd_id HDF id of the output file
!
! stag_dim Dimension of grid staggering (0-none, 1-x, 2-y, 3-z)
!
! hdfcompr Compression flag (0-none)
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny
REAL :: var(nx,ny)
INTEGER :: sd_id, stag_dim, hdfcompr
CHARACTER (LEN=*) :: name, comment, units
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(2),start(2),stride(2)
INTEGER :: sds_id
INTEGER :: comp_prm(1)
REAL :: amax,amin,scale
INTEGER :: i,j
INTEGER :: itmp1
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
IF (myproc == 0) &
WRITE (*,*) "HDFWRT2D: Writing variable ",trim(name)
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
start(1) = 0
start(2) = 0
stride(1) = 1
stride(2) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
IF (hdfcompr <= 3) THEN
sds_id = sfcreate(sd_id, trim(name), dfnt_float32, 2, dims)
ELSE
sds_id = sfcreate(sd_id, trim(name), dfnt_int16, 2, dims)
END IF
comp_prm(1) = 0
IF (hdfcompr == 1 .OR. hdfcompr == 5) THEN ! quick gzip
comp_prm(1) = 1
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 2 .OR. hdfcompr == 6) THEN ! high gzip
comp_prm(1) = 6
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 3 .OR. hdfcompr == 7) THEN ! huffman
comp_prm(1) = 4 ! this may be a problem on a Cray
istat = sfscompress(sds_id, 3, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 3)
END IF
istat = sfscatt(sds_id, 'hdf_comp_prm', dfnt_int32, 1, comp_prm)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
istat = sfsnatt(sds_id, 'stag_dim', dfnt_int32, 1, stag_dim)
!-----------------------------------------------------------------------
!
! If called for, map reals to 16 bit integers
!
!-----------------------------------------------------------------------
IF (hdfcompr > 3) THEN
CALL a3dmax0lcl(var,1,nx,1,nx,1,ny,1,ny,1,1,1,1,amax,amin)
IF (ABS(amax-amin) < 1.0E-10) amax = 1.1 * amin + 1.
scale = 65534.0 / (amax - amin)
DO j=1,ny
DO i=1,nx
itmp1 = nint(scale * (var(i,j) - amin)) - 32767
itmp(i,j) = itmp1
END DO
END DO
istat = sfsnatt(sds_id, 'packed16', dfnt_int32, 1, 1)
istat = sfsnatt(sds_id, 'max', dfnt_float32, 1, amax)
istat = sfsnatt(sds_id, 'min', dfnt_float32, 1, amin)
END IF
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
IF (hdfcompr <= 3) THEN
istat = sfwdata(sds_id, start, stride, dims, var)
ELSE
istat = sfwdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT2D: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT2D: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt2d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT2DI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfwrt2di(var,nx,ny,sd_id,stag_dim,hdfcompr, & 3
name,comment,units)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 2-D integer array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
!
! sd_id HDF id of the output file
!
! stag_dim Dimension of grid staggering (0-none, 1-x, 2-y, 3-z)
!
! hdfcompr Compression flag (0-none)
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
!
IMPLICIT NONE
INTEGER :: nx,ny
INTEGER :: var(nx,ny)
INTEGER :: sd_id, stag_dim, hdfcompr
CHARACTER (LEN=*) :: name, comment, units
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(2),start(2),stride(2)
INTEGER :: sds_id
INTEGER :: comp_prm(1)
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
start(1) = 0
start(2) = 0
stride(1) = 1
stride(2) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
IF (myproc == 0) &
WRITE (*,*) "HDFWRT2DI: Writing variable ",trim(name)
sds_id = sfcreate(sd_id, trim(name), dfnt_int32, 2, dims)
comp_prm(1) = 0
IF (hdfcompr == 1 .OR. hdfcompr == 5) THEN ! quick gzip
comp_prm(1) = 1
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 2 .OR. hdfcompr == 6) THEN ! high gzip
comp_prm(1) = 6
istat = sfscompress(sds_id, 4, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 4)
ELSE IF (hdfcompr == 3 .OR. hdfcompr == 7) THEN ! huffman
comp_prm(1) = 4 ! this may be a problem on a Cray
istat = sfscompress(sds_id, 3, comp_prm)
istat = sfscatt(sds_id, 'hdf_comp_code', dfnt_int32, 1, 3)
END IF
istat = sfscatt(sds_id, 'hdf_comp_prm', dfnt_int32, 1, comp_prm)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
istat = sfsnatt(sds_id, 'stag_dim', dfnt_int32, 1, stag_dim)
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
istat = sfwdata(sds_id, start, stride, dims, var)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT2DI: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT2DI: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt2di
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRT1D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfwrt1d(var,num,sd_id,name,comment,units) 17
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a 1-D real array to an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! var Array to be written to the file
!
! num Number of grid points
!
! sd_id HDF id of the output file
!
! name Variable name
! comment Destriptive string
! units String destribing units of var
!
! OUTPUT:
!
! None.
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: num
REAL :: var(num)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name, comment, units
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: istat
INTEGER :: dims(1),start(1),stride(1)
INTEGER :: sds_id
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfcreate, sfscompress, sfscatt, sfsnatt, &
sfwdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = num
start(1) = 0
stride(1) = 1
!-----------------------------------------------------------------------
!
! Create an entry for the variable, set compression and add attributes.
!
!-----------------------------------------------------------------------
IF (myproc == 0) &
WRITE (*,*) "HDFWRT1D: Writing variable ",trim(name)
sds_id = sfcreate(sd_id, trim(name), dfnt_float32, 1, dims)
IF (len_trim(comment) > 0) THEN
istat = sfscatt(sds_id, 'comment', dfnt_char8, &
len_trim(comment), comment)
ELSE
istat = sfscatt(sds_id, 'comment', dfnt_char8, 1, ' ')
END IF
IF (len_trim(units) > 0) THEN
istat = sfscatt(sds_id, 'units', dfnt_char8, &
len_trim(units), units)
ELSE
istat = sfscatt(sds_id, 'units', dfnt_char8, 1, ' ')
END IF
!-----------------------------------------------------------------------
!
! Write the data in var out to the file
!
!-----------------------------------------------------------------------
istat = sfwdata(sds_id, start, stride, dims, var)
IF (istat /= 0) THEN
IF(myproc == 0) &
WRITE (6,*) "HDFWRT1D: ERROR writing variable ",trim(name)
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
WRITE (6,*) "HDFWRT1D: ERROR writing variable ",trim(name)
END IF
RETURN
END SUBROUTINE hdfwrt1d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRTR ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfwrtr(sd_id,name,val,istat) 109
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a real attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/05
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! val The value of the attribute
! istat Status of the read (0-okay, 1-write error)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
REAL :: val
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc' ! mpi include file
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfsnatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sfsnatt(sd_id, trim(name), dfnt_float32, 1, val)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFWRTR: ERROR writing variable ",trim(name),"."
istat = 1
ELSE
! IF(myproc == 0) &
! WRITE (*,*) "HDFWRTR: Wrote variable ",trim(name)," value:",val
END IF
RETURN
END SUBROUTINE hdfwrtr
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRTI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfwrti(sd_id,name,val,istat) 128
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out an integer attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/05
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! val The value of the attribute
! istat Status of the read (0-okay, 1-write error)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: val
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc' ! mpi parameters
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfsnatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sfsnatt(sd_id, trim(name), dfnt_int32, 1, val)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFWRTI: ERROR writing variable ",trim(name),"."
istat = 1
ELSE
!IF (myproc == 0) &
!WRITE (*,*) "HDFWRTI: Wrote variable ",trim(name)," value:",val
END IF
RETURN
END SUBROUTINE hdfwrti
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFWRTC ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfwrtc(sd_id,strlen,name,string,istat) 23
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Write out a string attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/11
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! strlen Length of string to be written out (set to len(string)
! if strlen passed in as 0)
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! string The string to be written out
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
CHARACTER (LEN=*) :: string
INTEGER :: strlen, istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
! INTEGER :: tmplen
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc' ! mpi parameters
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfscatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
! IF (strlen == 0) THEN
! tmplen = LEN(string)
! ELSE
! tmplen = strlen
! END IF
IF (strlen ==0) strlen = LEN(string)
istat = sfscatt(sd_id, trim(name), dfnt_char8, strlen, string)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFWRTC: ERROR writing variable ",trim(name),"."
istat = 1
ELSE
!IF (myproc == 0) &
!WRITE (*,*) "HDFWRTC: Wrote variable ",trim(name)," value:", trim(string)
END IF
RETURN
END SUBROUTINE hdfwrtc
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD3D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd3d(sd_id,name,nx,ny,nz,var,istat,itmp,hmax,hmin) 74
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 3-D real array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
! hmax Used to store maximum values as a function of z
! hmin Used to store minimum values as a function of z
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz
REAL :: var(nx,ny,nz)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny,nz)
REAL :: hmax(nz),hmin(nz)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: dims(3),start(3),stride(3)
INTEGER :: sds_index,sds_id,attr_index
INTEGER :: packed16
INTEGER :: istat1,istat2,istat3
REAL :: scale
INTEGER :: i,j,k
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sffattr, sfrnatt, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
start(1) = 0
start(2) = 0
start(3) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3D: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
attr_index = sffattr(sds_id, "packed16")
IF (attr_index >= 0) THEN
istat1 = sfrnatt(sds_id, attr_index, packed16)
attr_index = sffattr(sds_id, "max")
istat2 = sfrnatt(sds_id, attr_index, hmax)
attr_index = sffattr(sds_id, "min")
istat3 = sfrnatt(sds_id, attr_index, hmin)
IF (istat1 == -1 .OR. istat2 == -1 .OR. istat3 == -1) THEN
WRITE (6,*) "HDFRD3D: ERROR reading max/min for ",trim(name)
istat = 2
RETURN
END IF
ELSE
packed16 = 0
END IF
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
!IF (myproc == 0) WRITE (*,*) "HDFRD3D: Reading variable ",trim(name)
IF (packed16 == 0) THEN
istat = sfrdata(sds_id, start, stride, dims, var)
ELSE
istat = sfrdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3D: ERROR reading variable ",trim(name),"."
istat = 2
RETURN
END IF
!-----------------------------------------------------------------------
!
! If called for, map 16 bit integers to reals
!
!-----------------------------------------------------------------------
IF (packed16 /= 0) THEN
DO k=1,nz
scale = (hmax(k) - hmin(k)) / 65534.0
DO j=1,ny
DO i=1,nx
var(i,j,k) = scale * (itmp(i,j,k) + 32767) + hmin(k)
END DO
END DO
END DO
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3D: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd3d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD4DSOIL ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd4d(sd_id,name,nx,ny,nz,nl,var,istat, & 4
itmp,hmax,hmin)
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 3-D real array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: J. Brotzge (modified from hdfrd3d)
! 05/23/2002
!
! MODIFICATION HISTORY:
!
! 07/02/2002 Zuwen He
!
! Clean up and make it a general 4d HDF writer
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
! nl Number of grid points in the 4th dimension
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
! hmax Used to store maximum values as a function of z
! hmin Used to store minimum values as a function of z
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz,nl
REAL :: var(nx,ny,nz,nl)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny,nz,nl)
REAL :: hmax(nz),hmin(nz)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
! INTEGER :: dims(3),start(3),stride(3)
INTEGER :: dims(4),start(4),stride(4)
INTEGER :: sds_index,sds_id,attr_index
INTEGER :: packed16
INTEGER :: istat1,istat2,istat3
REAL :: scale
INTEGER :: i,j,k,l
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sffattr, sfrnatt, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
dims(4) = nl
start(1) = 0
start(2) = 0
start(3) = 0
start(4) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
stride(4) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
WRITE (6,*) "HDFRD4DSOIL: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
attr_index = sffattr(sds_id, "packed16")
IF (attr_index >= 0) THEN
istat1 = sfrnatt(sds_id, attr_index, packed16)
attr_index = sffattr(sds_id, "max")
istat2 = sfrnatt(sds_id, attr_index, hmax)
attr_index = sffattr(sds_id, "min")
istat3 = sfrnatt(sds_id, attr_index, hmin)
IF (istat1 == -1 .OR. istat2 == -1 .OR. istat3 == -1) THEN
WRITE (6,*) "HDFRD4DSOIL: ERROR reading max/min for ",trim(name)
istat = 2
RETURN
END IF
ELSE
packed16 = 0
END IF
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
!IF (myproc == 0) WRITE (*,*) "HDFRD4DSOIL: Reading variable ",trim(name)
IF (packed16 == 0) THEN
istat = sfrdata(sds_id, start, stride, dims, var)
ELSE
istat = sfrdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat == -1) THEN
WRITE (6,*) "HDFRD3D: ERROR reading variable ",trim(name),"."
istat = 2
RETURN
END IF
!-----------------------------------------------------------------------
!
! If called for, map 16 bit integers to reals
!
!-----------------------------------------------------------------------
IF (packed16 /= 0) THEN
DO l=1,nl
DO k=1,nz
scale = (hmax(k) - hmin(k)) / 65534.0
DO j=1,ny
DO i=1,nx
var(i,j,k,l) = scale * (itmp(i,j,k,l) + 32767) + hmin(k)
END DO
END DO
END DO
END DO
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
WRITE (6,*) "HDFRD3D: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd4d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD3DI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd3di(sd_id,name,nx,ny,nz,var,istat) 3
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 3-D integer array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
! nz Number of grid points in the vertical
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny,nz
INTEGER :: var(nx,ny,nz)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: dims(3),start(3),stride(3)
INTEGER :: sds_index,sds_id
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
dims(3) = nz
start(1) = 0
start(2) = 0
start(3) = 0
stride(1) = 1
stride(2) = 1
stride(3) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3DI: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
! IF (myproc == 0) &
! WRITE (*,*) "HDFRD3DI: Reading variable ",trim(name)
istat = sfrdata(sds_id, start, stride, dims, var)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3DI: ERROR reading variable ",trim(name),"."
istat = 2
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD3DI: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd3di
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD2D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd2d(sd_id,name,nx,ny,var,istat,itmp) 26
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 2-D real array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
!
! itmp Scratch array for mapping reals to integers (used for
! some values of hdfcompr)
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny
REAL :: var(nx,ny)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
INTEGER (KIND=selected_int_kind(4)) :: itmp(nx,ny)
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: dims(2),start(2),stride(2)
INTEGER :: sds_index,sds_id,attr_index
REAL :: amax,amin, scale
INTEGER :: istat1,istat2,istat3
INTEGER :: i,j
INTEGER :: packed16
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sffattr, sfrnatt, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
start(1) = 0
start(2) = 0
stride(1) = 1
stride(2) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2D: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
attr_index = sffattr(sds_id, "packed16")
IF (attr_index >= 0) THEN
istat1 = sfrnatt(sds_id, attr_index, packed16)
attr_index = sffattr(sds_id, "max")
istat2 = sfrnatt(sds_id, attr_index, amax)
attr_index = sffattr(sds_id, "min")
istat3 = sfrnatt(sds_id, attr_index, amin)
IF (istat1 == -1 .OR. istat2 == -1 .OR. istat3 == -1) THEN
WRITE (6,*) "HDFRD2D: ERROR reading max/min for ",trim(name)
istat = sfendacc(sds_id)
istat = 2
RETURN
END IF
ELSE
packed16 = 0
END IF
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
! IF (myproc == 0) WRITE (*,*) "HDFRD2D: Reading variable ",trim(name)
IF (packed16 == 0) THEN
istat = sfrdata(sds_id, start, stride, dims, var)
ELSE
istat = sfrdata(sds_id, start, stride, dims, itmp)
END IF
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2D: ERROR reading variable ",trim(name),"."
istat = sfendacc(sds_id)
istat = 2
RETURN
END IF
!-----------------------------------------------------------------------
!
! If called for, map 16 bit integers to reals
!
!-----------------------------------------------------------------------
IF (packed16 /= 0) THEN
scale = (amax - amin) / 65534.0
DO j=1,ny
DO i=1,nx
var(i,j) = scale * (itmp(i,j) + 32767) + amin
END DO
END DO
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2D: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd2d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD2DI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd2di(sd_id,name,nx,ny,var,istat) 2
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 2-D integer array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! nx Number of grid points in the x-direction (east/west)
! ny Number of grid points in the y-direction (north/south)
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: nx,ny
INTEGER :: var(nx,ny)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: dims(2),start(2),stride(2)
INTEGER :: sds_index,sds_id
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = nx
dims(2) = ny
start(1) = 0
start(2) = 0
stride(1) = 1
stride(2) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2DI: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
! IF (myproc == 0) WRITE (*,*) "HDFRD2DI: Reading variable ",trim(name)
istat = sfrdata(sds_id, start, stride, dims, var)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2DI: ERROR reading variable ",trim(name),"."
istat = sfendacc(sds_id)
istat = 2
RETURN
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD2DI: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd2di
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRD1D ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrd1d(sd_id,name,num,var,istat) 5
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a 2-D real array from an HDF4 file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the output file
!
! num Number of grid points
!
! OUTPUT:
!
! var Array to be read in
!
! istat Status of read (0-okay, 1-error when reading)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: num
REAL :: var(num)
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: dims(1),start(1),stride(1)
INTEGER :: sds_index,sds_id
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfn2index, sfselect, sfrdata, sfendacc
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!-----------------------------------------------------------------------
!
! Initialize dimension parameters
!
!-----------------------------------------------------------------------
dims(1) = num
start(1) = 0
stride(1) = 1
!-----------------------------------------------------------------------
!
! Get the SDS ID for the variable.
!
!-----------------------------------------------------------------------
sds_index = sfn2index(sd_id, trim(name))
IF (sds_index == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD1D: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
sds_id = sfselect(sd_id, sds_index)
!-----------------------------------------------------------------------
!
! Read data into var.
!
!-----------------------------------------------------------------------
! IF (myproc == 0) WRITE (*,*) "HDFRD1D: Reading variable ",trim(name)
istat = sfrdata(sds_id, start, stride, dims, var)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD1D: ERROR reading variable ",trim(name),"."
istat = sfendacc(sds_id)
istat = 2
RETURN
END IF
istat = sfendacc(sds_id)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRD1D: ERROR reading variable ",trim(name)
istat = 2
END IF
RETURN
END SUBROUTINE hdfrd1d
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRDR ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrdr(sd_id,name,val,istat) 121
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a real attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! val The value of the attribute
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
REAL :: val
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: attr_index, istat
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sffattr, sfrnatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
attr_index = sffattr(sd_id, trim(name))
IF (attr_index < 0) THEN
WRITE (6,*) "HDFRDR: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
istat = sfrnatt(sd_id, attr_index, val)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDR: ERROR reading variable ",trim(name),"."
istat = 2
ELSE
! IF (myproc == 0) &
! WRITE (*,*) "HDFRDR: Read variable ",trim(name)," value:",val
END IF
RETURN
END SUBROUTINE hdfrdr
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRDI ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrdi(sd_id,name,val,istat) 149
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in an integer attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! val The value of the attribute
! istat Status indicator (0-okay, 1-variable not found, 2-read error)
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: val
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: attr_index, istat
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sffattr, sfrnatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
attr_index = sffattr(sd_id, trim(name))
IF (attr_index < 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDI: WARNING, variable ", &
trim(name)," not found in file."
istat = 1
RETURN
END IF
istat = sfrnatt(sd_id, attr_index, val)
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDI: ERROR reading variable ",trim(name),"."
istat = 2
val = -1
ELSE
! IF (myproc == 0) &
! WRITE (*,*) "HDFRDI: Read variable ",trim(name)," value:",val
END IF
RETURN
END SUBROUTINE hdfrdi
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFRDC ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfrdc(sd_id,max_len,name,string,istat) 26
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Read in a string attribute
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/03/15
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! name Variable name
!
! max_len Maximum allowable length of string to be read in
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
! OUTPUT:
!
! string The value of the attribute
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id
CHARACTER (LEN=*) :: name
CHARACTER (LEN=*) :: string
INTEGER :: max_len
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: attr_index, istat
INTEGER :: data_type, n_values
CHARACTER (LEN=128) :: attr_name
INTEGER :: i
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc'
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sffattr, sfgainfo, sfrnatt
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
attr_index = sffattr(sd_id, trim(name))
IF (attr_index < 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDC: WARNING, variable ",trim(name), &
" not located."
istat = 1
RETURN
END IF
istat = sfgainfo(sd_id, attr_index, attr_name, data_type, &
n_values)
IF (istat /= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDC: ERROR in looking up attributes for ", &
trim(name)
istat = 2
RETURN
END IF
IF (n_values <= max_len) THEN
istat = sfrnatt(sd_id, attr_index, string)
ELSE
IF (myproc == 0) &
WRITE (6,*) "HDFRDC: ERROR: string length for variable ", &
trim(name),", ",max_len, &
", is less than string in file:",n_values,","
IF (myproc == 0) &
WRITE (6,*) " value not read in."
END IF
DO i=n_values+1,len(string)
string(i:i) = " "
END DO
IF (istat == -1) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFRDC: ERROR reading variable ",trim(name),"."
istat = 1
ELSE
! IF (myproc == 0) &
! WRITE (*,*) "HDFRDC: Read variable ",trim(name)," value:", trim(string)
END IF
RETURN
END SUBROUTINE hdfrdc
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFOPEN ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfopen(filename,rdwrtflg,sd_id) 28
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Open a HDF file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/11
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! filename File name
!
! rdwrtflg Read/write flag (1-read, 2-write)
!
! OUTPUT:
!
! sd_id HDF id of the file or variable containing the
! named attribute
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id
CHARACTER (LEN=*) :: filename
INTEGER :: rdwrtflg
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
INCLUDE 'mp.inc' ! mpi parameters
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfstart
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
IF (rdwrtflg == 1) THEN
sd_id = sfstart(trim(filename), dfacc_read)
ELSE IF (rdwrtflg == 2) THEN
sd_id = sfstart(trim(filename), dfacc_create)
ELSE
sd_id = -1
IF (myproc == 0) &
WRITE (6,*) "HDFOPEN: ERROR, unsupported rdwrtflg value of", &
rdwrtflg
END IF
IF (sd_id <= 0) THEN
IF (myproc == 0) &
WRITE (6,*) "HDFOPEN: ERROR opening ",trim(filename)
ENDIF
RETURN
END SUBROUTINE hdfopen
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFCLOSE ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfclose(sd_id,istat) 26
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Close a HDF file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/04/11
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! sd_id HDF id of file to close
!
! OUTPUT:
!
! istat Status returned by close
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id
INTEGER :: istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfend
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sfend(sd_id)
RETURN
END SUBROUTINE hdfclose
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFINFO ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfinfo(sd_id,ndata,nattr,istat) 2
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Get the number of data sets and attributes contained in a HDF data file.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/10/18
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! sd_id HDF id of the data file
!
! OUTPUT:
!
! ndata Number of data sets in the file
! nattr Number of attributes in the file
! istat Status returned by sffinfo
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id
INTEGER :: ndata,nattr,istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sffinfo
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sffinfo(sd_id,ndata,nattr)
RETURN
END SUBROUTINE hdfinfo
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFDINFO ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfdinfo(sds_id,name,rank,dims,dtype,nattr,istat) 2
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Get the information about a data set
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/10/18
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! sds_id HDF id of the data set
!
! OUTPUT:
!
! name Name of the data set
! rank Number of dimensions in the data set
! dims Number of points for each dimension
! dtype Data type
! nattr Number of attributes in the data set
! istat Status returned by sfginfo
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sds_id
CHARACTER (LEN=*) :: name
INTEGER :: rank,dims(6),dtype,nattr,istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfginfo
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sfginfo(sds_id,name,rank,dims,dtype,nattr)
RETURN
END SUBROUTINE hdfdinfo
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE HDFAINFO ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
SUBROUTINE hdfainfo(sd_id,aindex,name,dtype,nvalues,istat) 4
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Get the information about an attribute.
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2000/10/25
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! sd_id HDF id of the data set (for HDF file or SDS data set in a file)
! aindex Index number for the attribute
!
! OUTPUT:
!
! name Name of the data set
! dtype Data type
! nvalues Number of values (number of characters if a string)
! istat Status returned by sfgainfo
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
INTEGER :: sd_id,aindex
CHARACTER (LEN=*) :: name
INTEGER :: dtype,nvalues,istat
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
INCLUDE 'hdf.f90' ! HDF4 library include file
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
INTEGER :: sfgainfo
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
istat = sfgainfo(sd_id,aindex,name,dtype,nvalues)
RETURN
END SUBROUTINE hdfainfo
!##################################################################
!##################################################################
!###### ######
!###### SUBROUTINE GET_DIMS_FROM_HDF ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!SUBROUTINE get_dims_from_hdf(filename,nx,ny,nz,nstyps,ireturn)
SUBROUTINE get_dims_from_hdf(filename,nx,ny,nz,nzsoil,nstyps,ireturn) 2,9
!-----------------------------------------------------------------------
!
! PURPOSE:
!
! Get the grid dimensions from an ARPS hdf file. (Similar to
! get_dims_from_data.)
!
!-----------------------------------------------------------------------
!
! AUTHOR: Gene Bassett
! 2001/04/23
!
! MODIFICATION HISTORY:
!
!-----------------------------------------------------------------------
!
! INPUT:
!
! filename File name
!
! OUTPUT:
!
! nx,ny,nz Grid dimensions
! nstyps Number of soil types
! ireturn Return status indicator
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
IMPLICIT NONE
CHARACTER (LEN=*) :: filename
! INTEGER :: nx,ny,nz,nstyps,ireturn
INTEGER :: nx,ny,nz,nzsoil,nstyps,ireturn
! 06/28/2002 Zuwen He
!
! fmtver??: to label each data a version.
! intver??: an integer to allow faster comparison than fmtver??,
! which are strings.
!
! Verion 5.00: significant change in soil variables since version 4.10.
!
CHARACTER (LEN=40) :: fmtver410,fmtver500
INTEGER :: intver,intver410,intver500
PARAMETER (fmtver410='004.10 HDF4 Coded Data',intver410=410)
PARAMETER (fmtver500='005.00 HDF4 Coded Data',intver500=500)
CHARACTER (LEN=40) :: fmtverin
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
INTEGER :: sd_id,istat
!-----------------------------------------------------------------------
!
! Include files:
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Functions
!
!-----------------------------------------------------------------------
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
ireturn = 0
CALL hdfopen(filename,1,sd_id)
IF (sd_id <= 0) THEN
ireturn = 1
RETURN
ENDIF
!EMK 12 July 2002
CALL hdfrdc(sd_id,40,"fmtver",fmtverin,istat)
IF (fmtverin == fmtver410) THEN
intver=intver410
ELSE IF (fmtverin == fmtver500) THEN
intver=intver500
ELSE
WRITE(6,'(/1x,a,a,a/)') &
'Incoming data format, fmtverin=',fmtverin, &
', not found. The Job stopped.'
CALL arpsstop('arpstop called from HDFREAD. ',1)
END IF
WRITE(6,'(/1x,a,a/)')'Incoming data format, fmtverin=',fmtverin
!EMK END 12 July 2002
CALL hdfrdi(sd_id,"nx",nx,istat)
IF (istat > 0) ireturn = 1
CALL hdfrdi(sd_id,"ny",ny,istat)
IF (istat > 0) ireturn = 1
CALL hdfrdi(sd_id,"nz",nz,istat)
IF (istat > 0) ireturn = 1
CALL hdfrdi(sd_id,"nstyp",nstyps,istat)
IF (istat > 0) ireturn = 1
!EMK 12 July 2002
IF (intver == intver500) THEN
CALL hdfrdi(sd_id,"nzsoil",nzsoil,istat)
IF (istat > 0) ireturn = 1
ELSE
nzsoil = 2
END IF
!EMK END 12 July 2002
CALL hdfclose(sd_id,istat)
RETURN
END SUBROUTINE get_dims_from_hdf