SUBROUTINE join_hdf (fileheader,nxsm,nysm,nz,nzsoil,nstyps,nxlg,nylg, & 3,33 buf_r,buf_rsm,buf_rsoil,buf_rsmsoil, & buf_i,buf_ism,buf_i16,buf_i16sm, & buf_i16soil,buf_i16soilsm, & buf_r1,buf_r2,sstat) IMPLICIT NONE INCLUDE 'mp.inc' INCLUDE 'hdf.f90' ! HDF4 library include file CHARACTER(LEN=*), INTENT(IN) :: fileheader INTEGER, INTENT(IN) :: nxsm,nysm,nz,nzsoil,nstyps,nxlg,nylg REAL, INTENT(OUT) :: buf_r(nxlg,nylg,nz), buf_rsm(nxsm,nysm,nz) REAL, INTENT(OUT) :: buf_rsoil(nxlg,nylg,nzsoil,0:nstyps), & buf_rsmsoil(nxsm,nysm,nzsoil,0:nstyps) INTEGER, INTENT(OUT) :: buf_i(nxlg,nylg,nz), buf_ism(nxsm,nysm,nz) INTEGER(KIND=selected_int_kind(4)), INTENT(OUT) :: & buf_i16(nxlg,nylg,nz), & buf_i16sm(nxsm,nysm,nz) INTEGER(KIND=selected_int_kind(4)), INTENT(OUT) :: & buf_i16soil(nxlg,nylg,nzsoil,0:nstyps), & buf_i16soilsm(nxsm,nysm,nzsoil,0:nstyps) REAL, INTENT(OUT) :: buf_r1(nxsm+nysm+nz), buf_r2(nxlg+nylg+nz) ! 1D arrays or max and min INTEGER, INTENT(OUT) :: sstat ! join status: sstat=1 if error encountered !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- CHARACTER (LEN=256) :: filename INTEGER :: fi, fj, i, j, k, l INTEGER :: kk INTEGER :: nxin, nyin, nzin, nzsoilin REAL :: amax, amin, scale INTEGER :: ierr INTEGER :: sd_id,sd_id1,sd_id2,ndata,nattr,istat,aindex,dindex INTEGER :: sds_id,sds_id1,sds_id2,rank,dims(6),dtype,ndattr,adtype INTEGER :: sds_index CHARACTER (LEN=256 ) :: name, aname CHARACTER (LEN=1024) :: char_attr CHARACTER (LEN=1), ALLOCATABLE :: lchar_attr(:) INTEGER :: nvalues INTEGER :: size(4),start(4),stride(4),strideout(4),sizeout(4) INTEGER :: x_off,y_off,i0,j0 INTEGER :: comp_code, comp_prm(1) INTEGER :: itmp ! fmtver??: to label each data a version. ! intver??: an integer to allow faster comparison than fmtver??, ! which are strings. ! ! Version 5.00: significant change in soil variables since version 4.10. ! CHARACTER (LEN=40) :: fmtver,fmtver410a,fmtver500a CHARACTER (LEN=40) :: fmtver410b,fmtver500b INTEGER :: intver,intver410,intver500 PARAMETER (fmtver410a='* 004.10 GrADS Soilvar Data',intver410=410) PARAMETER (fmtver500a='* 005.00 GrADS Soilvar Data',intver500=500) PARAMETER (fmtver410b='004.10 HDF4 Coded Data') PARAMETER (fmtver500b='005.00 HDF4 Coded Data') CHARACTER (LEN=40) :: fmtverin !----------------------------------------------------------------------- ! ! Functions ! !----------------------------------------------------------------------- INTEGER :: sfcreate, sfrdata, sfrnatt, sfscompress, sfselect, & sfsnatt, sfwdata, sffattr, sfendacc, sfn2index !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! nxlg = (nxsm-3)*nproc_x+3 ! nylg = (nysm-3)*nproc_y+3 stride = 1 start = 0 sstat = 0 !----------------------------------------------------------------------- ! ! Open file 0101, read in its attributes and check the ! dimensions in the file against the dimensions passed in. ! Open the joined file as well. ! !----------------------------------------------------------------------- CALL hdfopen(trim(fileheader)//'_0101',1,sd_id1) IF (sd_id1 < 0) THEN WRITE (6,*) "JOIN_HDF: ERROR opening ", & trim(fileheader)//'_0101'," for reading." sstat = 1 RETURN END IF CALL hdfopen(trim(fileheader),2,sd_id2) IF (sd_id2 < 0) THEN WRITE (6,*) "JOIN_HDF: ERROR creating HDF4 file: ", & trim(fileheader) sstat = 1 GO TO 600 END IF CALL hdfinfo(sd_id1,ndata,nattr,istat) fmtverin = "" CALL hdfrdc(sd_id1,40,"fmtver",fmtverin,istat) ! The following code is a dangerous practice ! but it may be the only way to distinguish ! versions prior to 500. ! IF (fmtverin == fmtver500a .OR. fmtverin == fmtver500b) THEN intver=intver500 ELSE intver=intver410 ! prior to 500, there is no fmtver variable istat=0 WRITE(6,'(/1x,a/,a/)') 'WARNING: Incoming data format are ' // & 'older than version 5.00!!! ', & 'It is to be read as if it version 4.10!!! ' END IF CALL hdfrdi(sd_id1,"nx",nxin,istat) CALL hdfrdi(sd_id1,"ny",nyin,istat) IF (nz > 1) THEN CALL hdfrdi(sd_id1,"nz",nzin,istat) ELSE nzin = nz ! the file is 2-D so it won't have nz ENDIF IF (nzsoil > 1) THEN ! Soil levels expected. IF (intver >= intver500) THEN ! New version with OUSoil model. CALL hdfrdi(sd_id1,"nzsoil",nzsoilin,istat) ELSE nzsoilin = 2 ! prior to version 500, it is 2 layer soil END IF ELSE nzsoilin = nzsoil END IF IF ((nxin /= nxsm).OR.(nyin /= nysm).OR.(nzin /= nz).OR. & (nzsoilin /= nzsoil)) THEN WRITE (*,*) "ERROR: mismatch in sizes." WRITE (*,*) "nxin,nyin,nzin,nzsoilin: ",nxin,nyin,nzin,nzsoilin WRITE (*,*) "nxsm,nysm,nz,nzsoil: ",nxsm,nysm,nz,nzsoil sstat = 1 RETURN END IF !----------------------------------------------------------------------- ! ! Read/write header info. ! !----------------------------------------------------------------------- DO aindex = 0,nattr-1 CALL hdfainfo(sd_id1,aindex,name,dtype,nvalues,istat) IF (dtype == dfnt_char8) THEN IF (nvalues > 1024) THEN ALLOCATE (lchar_attr(nvalues)) CALL hdfrdc(sd_id1,nvalues,name,lchar_attr,istat) CALL hdfwrtc(sd_id2,nvalues,name,lchar_attr,istat) DEALLOCATE(lchar_attr) ELSE CALL hdfrdc(sd_id1,nvalues,name,char_attr,istat) CALL hdfwrtc(sd_id2,nvalues,name,char_attr,istat) ENDIF ELSE IF (dtype == dfnt_float32) THEN istat = sfrnatt(sd_id1, aindex, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading attribute ",trim(name) sstat = 1 GOTO 600 ENDIF istat = sfsnatt(sd_id2, trim(name), dfnt_float32, nvalues, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR, writing attribute ",trim(name) sstat = 1 GOTO 600 ENDIF ELSE IF (dtype == dfnt_int32) THEN IF (trim(name) == 'nx') THEN CALL hdfwrti(sd_id2, 'nx', nxlg, istat) ELSE IF (trim(name) == 'ny') THEN CALL hdfwrti(sd_id2, 'ny', nylg, istat) ELSE istat = sfrnatt(sd_id1, aindex, buf_i) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading attribute ",trim(name) sstat = 1 GOTO 600 ENDIF istat = sfsnatt(sd_id2, trim(name), dfnt_int32, nvalues, buf_i) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR, writing attribute ",trim(name) sstat = 1 GOTO 600 ENDIF ENDIF ELSE WRITE (6,*) "JOIN_HDF: ERROR, unknown data type for ", & "attribute ", trim(name) sstat = 1 GOTO 600 ENDIF END DO !----------------------------------------------------------------------- ! ! Read/write each data set. ! !----------------------------------------------------------------------- DO dindex = 0,ndata-1 sds_id1 = sfselect(sd_id1,dindex) ! data set CALL hdfdinfo(sds_id1,name,rank,dims,dtype,ndattr,istat) IF (rank /= 1) THEN ! x, y, z etc. 1d arrays do not have comp_prm ! and hdf_comp_code - WYH - why? CALL hdfrdi(sds_id1,"hdf_comp_prm",comp_prm,istat) IF (comp_prm(1) /= 0) & CALL hdfrdi(sds_id1,"hdf_comp_code",comp_code,istat) END IF DO fj=1,nproc_y DO fi=1,nproc_x WRITE (filename, '(a,a,2i2.2)') trim(fileheader),'_',fi,fj CALL hdfopen(filename,1,sd_id) IF (sd_id < 0) THEN WRITE (6,*) "JOIN_HDF: ERROR opening HDF4 file: ", & trim(filename) sstat = 1 GO TO 600 END IF sds_index = sfn2index(sd_id, trim(name)) IF (sds_index == -1) THEN WRITE (6,*) "JOIN_HDF: ERROR, variable ", & trim(name)," not found in file",trim(filename),"." sstat = 1 GO TO 600 END IF sds_id = sfselect(sd_id, sds_index) x_off = (fi-1) * (nxsm-3) y_off = (fj-1) * (nysm-3) i0 = min(2,fi) j0 = min(2,fj) size(1) = nxsm size(2) = nysm size(3) = dims(3) size(4) = dims(4) ! EMK Test sizeout(1) = nxlg sizeout(2) = nylg sizeout(3) = dims(3) ! sizeout(4) = dims(4) ! EMK Test sizeout(4) = nstyps+1 IF (rank == 1) THEN IF (dtype /= dfnt_float32) THEN WRITE (6,*) "JOIN_HDF: ERROR, unsupported data type ", & "for 1-d variable ",trim(name) sstat = 1 GOTO 600 ENDIF IF (trim(name) == 'x') THEN ! x IF (fj == 1) THEN size(1) = nxsm istat = sfrdata(sds_id, start, stride, size, buf_r1) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),"." sstat = 1 GOTO 600 ENDIF DO i=i0,nxsm buf_r2(i+x_off) = buf_r1(i) END DO ENDIF ELSE IF (trim(name) == 'y') THEN ! y IF (fi == 1) THEN size(1) = nysm istat = sfrdata(sds_id, start, stride, size, buf_r1) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),"." sstat = 1 GOTO 600 ENDIF DO j=j0,nysm buf_r2(j+y_off) = buf_r1(j) END DO ENDIF ELSE IF (trim(name) == 'z') THEN ! z IF (fi == 1 .and. fj == 1) THEN size(1) = nz istat = sfrdata(sds_id, start, stride, dims, buf_r1) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),"." sstat = 1 GOTO 600 ENDIF ENDIF ELSE IF (fi == 1 .and. fj == 1) THEN istat = sfrdata(sds_id, start, stride, dims, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),"." sstat = 1 GOTO 600 ENDIF ENDIF ENDIF ELSE IF (dtype == dfnt_float32) THEN IF(intver <= intver410 .AND. & (TRIM(name) == "tsfc" .OR. TRIM(name) == "tsoil" .OR. & TRIM(name) == "wetsfc" .OR. TRIM(name) == "wetdp")) THEN istat = sfrdata(sds_id, start, stride, size, buf_rsm) IF (istat /= 0) THEN WRITE (6,*) "SPLIT_HDF: ERROR reading ",trim(name),"." sstat = 1 GOTO 600 ENDIF DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,1) = buf_rsm(i,j,1) END DO END DO ELSE IF(intver >= intver500 .AND. & (TRIM(name) == "tsoil" .OR. TRIM(name) == "qsoil")) THEN istat = sfrdata(sds_id, start, stride, size, buf_rsmsoil) IF (istat /= 0) THEN WRITE (6,*) "SPLIT_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF DO kk=0,nstyps DO k=1,nzsoil DO j=j0,nysm DO i=i0,nxsm buf_rsoil(i+x_off,j+y_off,k,kk) = & buf_rsmsoil(i,j,k,kk) END DO END DO END DO END DO ELSE IF (intver >= intver500 .AND. & (TRIM(name) == "zpsoil")) THEN istat = sfrdata(sds_id, start, stride,size,buf_rsmsoil(1,1,1,0)) IF (istat /= 0) THEN WRITE (6,*) "SPLIT_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF DO k=1,nzsoil DO j=j0,nysm DO i=i0,nxsm buf_rsoil(i+x_off,j+y_off,k,0) = & buf_rsmsoil(i,j,k,0) END DO END DO END DO ELSE ! Not a soil variable istat = sfrdata(sds_id, start, stride, size, buf_rsm) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF DO k=1,nz DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,k) = buf_rsm(i,j,k) END DO END DO END DO END IF ELSE IF (dtype == dfnt_int32) THEN istat = sfrdata(sds_id, start, stride, size, buf_ism) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF DO k=1,nz DO j=j0,nysm DO i=i0,nxsm buf_i(i+x_off,j+y_off,k) = buf_ism(i,j,k) END DO END DO END DO ELSE IF (dtype == dfnt_int16) THEN IF (rank == 2) THEN aindex = sffattr(sds_id, "max") istat = sfrnatt(sds_id, aindex, amax) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading max for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF aindex = sffattr(sds_id, "min") istat = sfrnatt(sds_id, aindex, amin) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading min for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF istat = sfrdata(sds_id, start, stride, size, buf_i16sm) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF scale = (amax - amin) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,1) = & scale * (buf_i16sm(i,j,1) + 32767) + amin END DO END DO ELSE IF (rank == 3) THEN aindex = sffattr(sds_id, "max") istat = sfrnatt(sds_id, aindex, buf_r2) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading max for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF aindex = sffattr(sds_id, "min") istat = sfrnatt(sds_id, aindex, buf_r1) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading min for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF IF( TRIM(name) == "zpsoil") THEN istat = sfrdata(sds_id, start, stride, size, buf_i16soilsm(:,:,:,0)) ELSE istat = sfrdata(sds_id, start, stride, size, buf_i16sm) END IF IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF IF( TRIM(name) == "zpsoil" ) THEN DO k = 1, nzsoil scale = (buf_r2(k) - buf_r1(k)) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_rsoil(i+x_off,j+y_off,k,0) = scale * & (buf_i16soilsm(i,j,k,0) + 32767) + buf_r1(k) END DO END DO END DO ELSE IF( TRIM(name) == "stypfrct" ) THEN DO k=1, nstyps scale = (buf_r2(k) - buf_r1(k)) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,k) = & scale * (buf_i16sm(i,j,k) + 32767) + buf_r1(k) END DO END DO END DO ELSE IF( TRIM(name) == "wetcanp" ) THEN DO k=1, nstyps+1 scale = (buf_r2(k) - buf_r1(k)) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,k) = & scale * (buf_i16sm(i,j,k) + 32767) + buf_r1(k) END DO END DO END DO ELSE DO k=1,nz scale = (buf_r2(k) - buf_r1(k)) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_r(i+x_off,j+y_off,k) = & scale * (buf_i16sm(i,j,k) + 32767) + buf_r1(k) END DO END DO END DO END IF ELSE IF(rank == 4) THEN aindex = sffattr(sds_id, "max") istat = sfrnatt(sds_id, aindex, buf_r2) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading max for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF aindex = sffattr(sds_id, "min") istat = sfrnatt(sds_id, aindex, buf_r1) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading min for ", & trim(name),", exiting" sstat = 1 GOTO 600 ENDIF istat = sfrdata(sds_id, start, stride, size, buf_i16soilsm) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading ",trim(name),", exiting" sstat = 1 GOTO 600 ENDIF DO l = 0, nstyps DO k = 1, nzsoil scale = (buf_r2(k) - buf_r1(k)) / 65534.0 DO j=j0,nysm DO i=i0,nxsm buf_rsoil(i+x_off,j+y_off,k,l) = scale * & (buf_i16soilsm(i,j,k,l) + 32767) + buf_r1(k) END DO END DO END DO END DO ELSE WRITE (6,*) "JOIN_HDF: ERROR, unsupported rank,",rank, & " for 16 bit remapped data" GOTO 600 ENDIF ELSE WRITE (6,*) "JOIN_HDF: ERROR, unknown data type for ", & "attribute ", trim(aname), " of sds ",trim(name) sstat = 1 GOTO 600 ENDIF ENDIF CALL hdfclose(sd_id,istat) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR on close of file ",trim(filename) sstat = 1 GOTO 600 ENDIF END DO ! fi END DO ! fj IF (rank == 1) THEN IF (trim(name) == 'x') THEN ! x sizeout(1) = nxlg sds_id2 = sfcreate(sd_id2, trim(name), dfnt_float32, 1, sizeout) istat = sfwdata(sds_id2, start, stride, sizeout, buf_r2) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE IF (trim(name) == 'y') THEN ! y sizeout(1) = nylg sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,1,sizeout) istat = sfwdata(sds_id2, start, stride, sizeout, buf_r2) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE IF (trim(name) == 'z') THEN ! z sizeout(1) = nz sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,1,sizeout) istat = sfwdata(sds_id2, start, stride, sizeout, buf_r1) ! note, buf_r1 not buf_r2 IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,1,dims) istat = sfwdata(sds_id2, start, stride, dims, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ENDIF ELSE IF (dtype == dfnt_float32) THEN IF (intver <= intver410 .AND. & (TRIM(name) == "tsfc" .OR. TRIM(name) == "tsoil" .OR. & TRIM(name) == "wetsfc" .OR. TRIM(name) == "wetdp")) THEN sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF istat = sfwdata(sds_id2, start, stride, sizeout, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE IF (intver >= intver500 .AND. & (TRIM(name) == "tsoil" .OR. TRIM(name) == "qsoil")) THEN sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF istat = sfwdata(sds_id2, start, stride, sizeout, buf_rsoil) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE IF (intver >= intver500 .AND. (TRIM(name) == "zpsoil")) THEN sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF istat = sfwdata(sds_id2, start, stride, sizeout, buf_rsoil(1,1,1,0)) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE sds_id2 = sfcreate(sd_id2,trim(name),dfnt_float32,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF istat = sfwdata(sds_id2, start, stride, sizeout, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF END IF ELSE IF (dtype == dfnt_int32) THEN sds_id2 = sfcreate(sd_id2,trim(name),dfnt_int32,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF istat = sfwdata(sds_id2, start, stride, sizeout, buf_i) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ELSE IF (dtype == dfnt_int16) THEN sds_id2 = sfcreate(sd_id2,trim(name),dfnt_int16,rank,sizeout) IF (comp_prm(1) > 0) THEN istat = sfscompress(sds_id2, comp_code, comp_prm) ENDIF IF (rank == 2) THEN CALL a3dmax0lcl(buf_r,1,nxlg,1,nxlg,1,nylg,1,nylg,1,1,1,1,amax,amin) IF (ABS(amax-amin) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (amax - amin) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_r(i,j,1) - amin)) - 32767 buf_i16(i,j,1) = itmp END DO END DO buf_r1(1) = amin buf_r2(1) = amax kk = 1 ELSE IF (rank == 3) THEN IF (TRIM(name) == "zpsoil" ) THEN DO k = 1, nzsoil CALL a3dmax0lcl(buf_rsoil(1,1,k,0),1,nxlg,1,nxlg,1,nylg,1,nylg, & 1,nzsoil,1,1,amax,amin) buf_r2(k) = amax buf_r1(k) = amin IF (ABS(buf_r2(k)-buf_r1(k)) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (buf_r2(k) - buf_r1(k)) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_rsoil(i,j,k,0) - buf_r1(k))) - 32767 buf_i16soil(i,j,k,0) = itmp END DO END DO END DO kk = nzsoil ELSE IF (TRIM(name) == "stypfrct" ) THEN DO k = 1, nstyps CALL a3dmax0lcl(buf_r(1,1,k),1,nxlg,1,nxlg,1,nylg,1,nylg, & 1,nstyps,1,1,amax,amin) buf_r2(k) = amax buf_r1(k) = amin IF (ABS(buf_r2(k)-buf_r1(k)) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (buf_r2(k) - buf_r1(k)) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_r(i,j,k) - buf_r1(k))) - 32767 buf_i16(i,j,k) = itmp END DO END DO END DO kk = nstyps ELSE IF (TRIM(name) == "wetcanp" ) THEN DO k = 1, nstyps+1 CALL a3dmax0lcl(buf_r(1,1,k),1,nxlg,1,nxlg,1,nylg,1,nylg, & 1,nstyps+1,1,1,amax,amin) buf_r2(k) = amax buf_r1(k) = amin IF (ABS(buf_r2(k)-buf_r1(k)) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (buf_r2(k) - buf_r1(k)) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_r(i,j,k) - buf_r1(k))) - 32767 buf_i16(i,j,k) = itmp END DO END DO END DO kk = nstyps + 1 ELSE DO k=1,nz CALL a3dmax0lcl(buf_r(1,1,k),1,nxlg,1,nxlg,1,nylg,1,nylg, & 1,nz,1,1,buf_r2(k),buf_r1(k)) IF (ABS(buf_r2(k)-buf_r1(k)) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (buf_r2(k) - buf_r1(k)) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_r(i,j,k) - buf_r1(k))) - 32767 buf_i16(i,j,k) = itmp END DO END DO END DO kk = nz END IF ELSE IF (rank == 4) THEN DO l = 0, nstyps DO k = 1, nzsoil CALL a3dmax0lcl(buf_rsoil(1,1,k,l),1,nxlg,1,nxlg,1,nylg,1,nylg, & 1,nzsoil,1,1,amax,amin) IF(l == 0) THEN buf_r2(k) = amax buf_r1(k) = amin END IF buf_r2(k) = MAX(buf_r2(k),amax) buf_r1(k) = MIN(buf_r1(k),amin) END DO END DO DO l = 0, nstyps DO k = 1, nzsoil IF (ABS(buf_r2(k)-buf_r1(k)) < 1.0E-10) THEN scale = 65534.0 ELSE scale = 65534.0 / (buf_r2(k) - buf_r1(k)) END IF DO j=1,nylg DO i=1,nxlg itmp = nint(scale * (buf_rsoil(i,j,k,l) - buf_r1(k))) - 32767 buf_i16soil(i,j,k,l) = itmp END DO END DO END DO END DO kk = nzsoil ENDIF IF(rank < 4 .AND. TRIM(name) /= "zpsoil" ) THEN istat = sfwdata(sds_id2, start, stride, sizeout, buf_i16) ELSE istat = sfwdata(sds_id2, start, stride, sizeout, buf_i16soil) END IF IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing ",trim(name), & " to file ",trim(fileheader)," , exiting" sstat = 1 GOTO 600 ENDIF ENDIF ENDIF ! data set attributes DO aindex = 0,ndattr-1 CALL hdfainfo(sds_id1,aindex,aname,adtype,nvalues,istat) IF ( dtype == dfnt_int16 .AND. TRIM(aname) == "packed16" ) THEN istat = sfsnatt(sds_id2, 'packed16', dfnt_int32, 1, 1) CYCLE ELSE IF ( dtype == dfnt_int16 .AND. TRIM(aname) == "min" ) THEN istat = sfsnatt(sds_id2, 'min', dfnt_float32, kk, buf_r1) CYCLE ELSE IF ( dtype == dfnt_int16 .AND. TRIM(aname) == "max" ) THEN istat = sfsnatt(sds_id2, 'max', dfnt_float32, kk, buf_r2) CYCLE END IF IF (adtype == dfnt_char8) THEN IF (nvalues > 1024) THEN ALLOCATE (lchar_attr(nvalues)) CALL hdfrdc(sds_id1,nvalues,aname,lchar_attr,istat) CALL hdfwrtc(sds_id2,nvalues,aname,lchar_attr,istat) DEALLOCATE(lchar_attr) ELSE CALL hdfrdc(sds_id1,nvalues,aname,char_attr,istat) CALL hdfwrtc(sds_id2,nvalues,aname,char_attr,istat) ENDIF ELSE IF (adtype == dfnt_float32) THEN IF (size(3) == nzsoil) THEN istat = sfrnatt(sds_id1, aindex, buf_rsoil) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF istat = sfsnatt(sds_id2, trim(aname), dfnt_float32, nvalues, & buf_rsoil) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR, writing attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF ELSE istat = sfrnatt(sds_id1, aindex, buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF istat = sfsnatt(sds_id2, trim(aname), dfnt_float32, nvalues, & buf_r) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR, writing attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF END IF ELSE IF (adtype == dfnt_int32) THEN istat = sfrnatt(sds_id1, aindex, buf_i) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR reading attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF istat = sfsnatt(sds_id2, trim(aname), dfnt_int32, nvalues, buf_i) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR, writing attribute ",trim(aname) sstat = 1 GOTO 600 ENDIF ELSE WRITE (6,*) "JOIN_HDF: ERROR, unknown data type for ", & "attribute ", trim(aname) sstat = 1 GOTO 600 ENDIF ! adtype == END DO istat = sfendacc(sds_id2) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR writing variable ",trim(name) END IF END DO ! dindex !----------------------------------------------------------------------- ! ! Close I/O and deallocate. ! !---------------------------------------------------------------------- 600 CONTINUE CALL hdfclose(sd_id1,istat) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR on close of file ",trim(fileheader)//'_0101' sstat = 1 ENDIF CALL hdfclose(sd_id2,istat) IF (istat /= 0) THEN WRITE (6,*) "JOIN_HDF: ERROR on close of file ",trim(fileheader) sstat = 1 ENDIF RETURN END SUBROUTINE join_hdf