PROGRAM bin2gem,32
!
!##################################################################
!##################################################################
!###### ######
!###### PROGRAM BIN2GEM ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
!
!-----------------------------------------------------------------------
! PURPOSE:
! Read in binary dump produced by arpspost and write out in GEMPAK
!
! Author : Fanyou Kong
! History: March 31, 2007: First code implementation
!
!-----------------------------------------------------------------------
!
! MODIFIED:
!
!-----------------------------------------------------------------------
! IMPLICIT NONE
INCLUDE 'GEMPRM.PRM'
INCLUDE 'mp.inc'
INTEGER :: nx,ny,nz
INTEGER :: ireturn, iret, ierr, istatus
CHARACTER (LEN=256) :: filename,filnam
CHARACTER (LEN=256) :: binheader,gemoutheader
INTEGER :: lengbf,lenfil,lenbin,lengem
INTEGER :: ihr
CHARACTER (LEN=4) :: ihrc
CHARACTER (LEN=6) :: tlevel(100)
CHARACTER (LEN=6) :: varname,varunit
INTEGER :: i,j,k, itags,itagr
INTEGER :: nf
NAMELIST /bin2gem_data/nf,tlevel,binheader,gemoutheader
REAL :: lat1,lon1,lat2,lon2
REAL :: ctime
INTEGER :: year,month,day,hour,minute
REAL, ALLOCATABLE :: hgtp(:,:,:),uwndp(:,:,:),vwndp(:,:,:),wwndp(:,:,:)
REAL, ALLOCATABLE :: tmpp(:,:,:),sphp(:,:,:)
REAL, ALLOCATABLE :: uwndh(:,:,:),vwndh(:,:,:),wwndh(:,:,:)
REAL, ALLOCATABLE :: psf(:,:),mslp(:,:),vort500(:,:)
REAL, ALLOCATABLE :: temp2m(:,:),dewp2m(:,:),qv2m(:,:)
REAL, ALLOCATABLE :: u10m(:,:),v10m(:,:),hgtsfc(:,:)
REAL, ALLOCATABLE :: refl1km(:,:),refl4km(:,:),cmpref(:,:)
REAL, ALLOCATABLE :: accppt(:,:),convppt(:,:)
REAL, ALLOCATABLE :: cape(:,:),mcape(:,:),cin(:,:),mcin(:,:),lcl(:,:)
REAL, ALLOCATABLE :: srh01(:,:),srh03(:,:),uh25(:,:),sh01(:,:),sh06(:,:)
REAL, ALLOCATABLE :: thck(:,:),li(:,:),brn(:,:),pw(:,:)
!-----------------------------------------------------------------------
! nprgem: number of pressure levels for GEMPAK file.
!
! iprgem: actual pressure levels defined by the user (mb)
!
INTEGER :: nprgem
!
PARAMETER (nprgem = 5)
!
INTEGER :: iprgem(nprgem)
!
DATA iprgem / 850, 700, 600, 500, 250/
!
!-----------------------------------------------------------------------
! nhtgem: number of height levels for GEMPAK file.
!
! ihtgem: actual height levels defined by the user (km)
!
INTEGER :: nhtgem
!
PARAMETER (nhtgem = 6)
!
INTEGER :: ihtgem(nhtgem)
!
DATA ihtgem / 1, 2, 3, 4, 5, 6/
CHARACTER (LEN=6) :: hgtptag(5),uwndptag(5)
CHARACTER (LEN=6) :: vwndptag(5),wwndptag(5)
CHARACTER (LEN=6) :: tmpptag(5),qvptag(5)
CHARACTER (LEN=6) :: uwndhtag(6),vwndhtag(6),wwndhtag(6)
DATA hgtptag / 'hgt850','hgt700','hgt600','hgt500','hgt250'/
DATA uwndptag / 'u850__','u700__','u600__','u500__','u250__'/
DATA vwndptag / 'v850__','v700__','v600__','v500__','v250__'/
DATA wwndptag / 'w850__','w700__','w600__','w500__','w250__'/
DATA tmpptag / 'tmp850','tmp700','tmp600','tmp500','tmp250'/
DATA qvptag / 'sph850','sph700','sph600','sph500','sph250'/
!
DATA uwndhtag /'u1km__','u2km__','u3km__','u4km__','u5km__','u6km__'/
DATA vwndhtag /'v1km__','v2km__','v3km__','v4km__','v5km__','v6km__'/
DATA wwndhtag /'w1km__','w2km__','w3km__','w4km__','w5km__','w6km__'/
!-----------------------------------------------------------------------
!
! GEMPAK variables
!
!-----------------------------------------------------------------------
!
CHARACTER (LEN=72) :: gdfile,gdarea
REAL :: rnvblk (llnnav), anlblk (llnanl)
CHARACTER (LEN=3) :: cproj(0:3)
DATA cproj /'CED','STR','LCC','MER'/
CHARACTER (LEN=72) :: chproj
REAL :: deltan
REAL :: deltax,deltay
REAL :: angle1,angle2,angle3
LOGICAL :: angflg
PARAMETER ( deltan= 1., &
deltax= -9999., &
deltay= -9999., &
angflg=.true.)
REAL :: gbnds(4)
LOGICAL :: wrtflg,rewrite,notopn,gemExist,tppExist,cppExist
COMMON /gemblk/rnvblk,anlblk,notopn
nproc_x = 1
nproc_y = 1
readsplit = 0
joindmp = 1
nprocs = 1
max_fopen = 1
READ(5,bin2gem_data,END=200)
write(6,bin2gem_data)
GO TO 20
200 WRITE(6,'(a)') &
'Error reading NAMELIST block bin2gem_data. ', &
'Default values used.'
20 CONTINUE
!-----------------------------------------------------------------------
!
! Get dimension nx,ny
!
!-----------------------------------------------------------------------
OPEN(4,FILE=trim(binheader)//'.sfpres'//tlevel(1),STATUS='old', &
FORM='unformatted',ERR=9000, IOSTAT=istat)
READ(4, ERR=9000, END=9000, IOSTAT=istat) nx,ny,nz
CLOSE(4)
print *,'nx,ny,nz:',nx,ny,nz
!
!-----------------------------------------------------------------------
!
! Allocate arrays
!
!-----------------------------------------------------------------------
!
ALLOCATE(hgtp(nx,ny,nprgem),STAT=istatus)
hgtp=0
ALLOCATE(uwndp(nx,ny,nprgem),STAT=istatus)
uwndp=0
ALLOCATE(vwndp(nx,ny,nprgem),STAT=istatus)
vwndp=0
ALLOCATE(wwndp(nx,ny,nprgem),STAT=istatus)
wwndp=0
ALLOCATE(tmpp(nx,ny,nprgem),STAT=istatus)
tmpp=0
ALLOCATE(sphp(nx,ny,nprgem),STAT=istatus)
sphp=0
ALLOCATE(uwndh(nx,ny,nhtgem),STAT=istatus)
uwndh=0
ALLOCATE(vwndh(nx,ny,nhtgem),STAT=istatus)
vwndh=0
ALLOCATE(wwndh(nx,ny,nhtgem),STAT=istatus)
wwndh=0
ALLOCATE(psf(nx,ny),STAT=istatus)
psf=0
ALLOCATE(mslp(nx,ny),STAT=istatus)
mslp=0
ALLOCATE(vort500(nx,ny),STAT=istatus)
vort500=0
ALLOCATE(temp2m(nx,ny),STAT=istatus)
temp2m=0
ALLOCATE(dewp2m(nx,ny),STAT=istatus)
dewp2m=0
ALLOCATE(qv2m(nx,ny),STAT=istatus)
qv2m=0
ALLOCATE(u10m(nx,ny),STAT=istatus)
u10m=0
ALLOCATE(v10m(nx,ny),STAT=istatus)
v10m=0
ALLOCATE(hgtsfc(nx,ny),STAT=istatus)
hgtsfc=0
ALLOCATE(refl1km(nx,ny),STAT=istatus)
refl1km=0
ALLOCATE(refl4km(nx,ny),STAT=istatus)
refl4km=0
ALLOCATE(cmpref(nx,ny),STAT=istatus)
cmpref=0
ALLOCATE(accppt(nx,ny),STAT=istatus)
accppt=0
ALLOCATE(convppt(nx,ny),STAT=istatus)
convppt=0
ALLOCATE(cape(nx,ny),STAT=istatus)
cape=0
ALLOCATE(mcape(nx,ny),STAT=istatus)
mcape=0
ALLOCATE(cin(nx,ny),STAT=istatus)
cin=0
ALLOCATE(mcin(nx,ny),STAT=istatus)
mcin=0
ALLOCATE(lcl(nx,ny),STAT=istatus)
lcl=0
ALLOCATE(srh01(nx,ny),STAT=istatus)
srh01=0
ALLOCATE(srh03(nx,ny),STAT=istatus)
srh03=0
ALLOCATE(uh25(nx,ny),STAT=istatus)
uh25=0
ALLOCATE(sh01(nx,ny),STAT=istatus)
sh01=0
ALLOCATE(sh06(nx,ny),STAT=istatus)
sh06=0
ALLOCATE(thck(nx,ny),STAT=istatus)
thck=0
ALLOCATE(li(nx,ny),STAT=istatus)
li=0
ALLOCATE(brn(nx,ny),STAT=istatus)
brn=0
ALLOCATE(pw(nx,ny),STAT=istatus)
pw=0
! Read in map parameters
OPEN(4,file=trim(binheader)//'.map')
READ(4,*) mapproj,trulon,trulat1,trulat2,lat1,lon1,lat2,lon2
CLOSE(4)
WRITE(6,'(/a,4(a,F7.2),a/)') ' ARPS domain ', &
'SW: (',lat1,',',lon1, ') NE: (',lat2,',',lon2,')'
!
!-----------------------------------------------------------------------
!
! Initializing GEMPAK
!
!-----------------------------------------------------------------------
!
! notopn = .true.
DO ntime=1,nf
CALL in_bdta(iret) ! initializing GEMPAK sans TAE
chproj = cproj(mapproj)
angle2 = trulon
!
! modify projection if polar stereographic
!
IF (chproj == 'STR') THEN
angle1=90.0
angle3=0.0
WRITE (*,*) '**********************************'
WRITE (*,*) '**********************************'
WRITE (*,*) ' '
WRITE (*,*) 'SETTING ANGLE1=90. FOR GEMPAK FILE'
WRITE (*,*) 'SETTING ANGLE3=0.0 FOR GEMPAK FILE'
WRITE (*,*) '(GEMPAK origin is at North Pole)'
WRITE (*,*) ' '
WRITE (*,*) '**********************************'
WRITE (*,*) '**********************************'
ELSE
angle1=trulat1
angle3=trulat2
WRITE (*,*) '**********************************'
WRITE (*,*) '**********************************'
WRITE (*,*) ' '
WRITE (*,*) ' SETTING ANGLE1=TRULAT1 '
WRITE (*,*) ' SETTING ANGLE3=TRULAT2 '
WRITE (*,*) ' '
WRITE (*,*) '**********************************'
WRITE (*,*) '**********************************'
END IF
gbnds(1)=lat1
gbnds(2)=lon1
gbnds(3)=lat2
gbnds(4)=lon2
WRITE(gdarea,'(f8.4,a1,f9.4,a1,f8.4,a1,f9.4)') &
lat1,';',lon1,';',lat2,';',lon2
call gr_mnav(chproj, nx,ny,lat1,lon1,lat2,lon2, &
angle1,angle2,angle3,angflg, &
rnvblk, iret)
IF( iret /= 0 ) GO TO 950
!
!-----------------------------------------------------------------------
!
! Build analysis block
!
!-----------------------------------------------------------------------
!
! CALL gr_mban ( deltan, deltax, deltay, &
! gbnds, gbnds, gbnds, anlblk, iret )
!
! IF( iret /= 0 ) GO TO 950
!
! Read 2D Fields ...
filnam = trim(binheader)//'.sfpres'//tlevel(ntime)
CALL binread2d
(nx,ny,trim(filnam),psf)
print *,trim(filnam)
filnam = trim(binheader)//'.mspres'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),mslp)
print *,trim(filnam)
filnam = trim(binheader)//'.accppt'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),accppt)
print *,trim(filnam)
filnam = trim(binheader)//'.pwat__'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),pw)
print *,trim(filnam)
CALL a3dmax0(pw,1,nx,1,nx-1,1,ny,1,ny,1,1,1,1, amax,amin)
print *,'PWAT -- amax, amin: ',amax, amin
filnam = trim(binheader)//'.temp2m'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),temp2m)
print *,trim(filnam)
filnam = trim(binheader)//'.dewp2m'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),dewp2m)
print *,trim(filnam)
filnam = trim(binheader)//'.qv2m__'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),qv2m)
print *,trim(filnam)
filnam = trim(binheader)//'.u10m__'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),u10m)
print *,trim(filnam)
filnam = trim(binheader)//'.v10m__'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),v10m)
print *,trim(filnam)
filnam = trim(binheader)//'.hgtsfc'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),hgtsfc)
print *,trim(filnam)
! filnam = trim(binheader)//'.vrt500'//tlevel(ntime)
! CALL binread2d(nx,ny,trim(filnam),vort500)
! print *,trim(filnam)
filnam = trim(binheader)//'.ref1km'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),refl1km)
print *,trim(filnam)
filnam = trim(binheader)//'.ref4km'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),refl4km)
print *,trim(filnam)
filnam = trim(binheader)//'.cmpref'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),cmpref)
print *,trim(filnam)
filnam = trim(binheader)//'.sbcape'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),cape)
print *,trim(filnam)
filnam = trim(binheader)//'.mucape'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),mcape)
print *,trim(filnam)
filnam = trim(binheader)//'.sbcins'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),cin)
print *,trim(filnam)
filnam = trim(binheader)//'.mucins'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),mcin)
print *,trim(filnam)
filnam = trim(binheader)//'.sblcl_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),lcl)
print *,trim(filnam)
filnam = trim(binheader)//'.srh01_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),srh01)
print *,trim(filnam)
filnam = trim(binheader)//'.srh03_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),srh03)
print *,trim(filnam)
filnam = trim(binheader)//'.vhel__'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),uh25)
print *,trim(filnam)
filnam = trim(binheader)//'.shr01_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),sh01)
print *,trim(filnam)
filnam = trim(binheader)//'.shr06_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),sh06)
print *,trim(filnam)
filnam = trim(binheader)//'.shr06_'//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),sh06)
print *,trim(filnam)
IF( nprgem > 0 ) THEN
DO klev=1,nprgem
filnam = trim(binheader)//'.'//hgtptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),hgtp(1,1,klev))
print *,trim(filnam)
filnam = trim(binheader)//'.'//uwndptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),uwndp(1,1,klev))
print *,trim(filnam)
filnam = trim(binheader)//'.'//vwndptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),vwndp(1,1,klev))
print *,trim(filnam)
filnam = trim(binheader)//'.'//wwndptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),wwndp(1,1,klev))
print *,trim(filnam)
filnam = trim(binheader)//'.'//tmpptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),tmpp(1,1,klev))
print *,trim(filnam)
filnam = trim(binheader)//'.'//qvptag(klev)//tlevel(ntime)
CALL binread2d(nx,ny,trim(filnam),sphp(1,1,klev))
print *,trim(filnam)
END DO ! nprgem-loop
END IF
!CALL a3dmax0(tmpp(1,1,1),1,nx,1,nx-1,1,ny,1,ny,1,1,1,1, amax,amin)
!print *,'T 850 -- amax, amin: ',amax, amin
! Read in time parameters
OPEN(4,file=trim(binheader)//'.time'//tlevel(ntime))
READ(4,*) year,month,day,hour,minute,ctime
CLOSE(4)
print *,'year,month,day,hour,minute,ctime:'
print *,year,month,day,hour,minute,ctime
! Find forecast hour for name construction
ihr = int(ctime/3600.)
WRITE(ihrc,'(a1,i3.3)') 'f',ihr
filename=trim(gemoutheader)//ihrc
lenfil = len_trim(filename)
WRITE(6,'(a,a/)') &
' GEMPAK output file: ', filename(1:lenfil)
!if(0>1) then
CALL enswrtgem (nx,ny,filename(1:lenfil), lenfil, &
ctime,year,month,day,hour,minute, &
iprgem,nprgem,ihtgem,nhtgem, &
hgtp,uwndp,vwndp,wwndp,tmpp,sphp, &
uwndh,vwndh,wwndh, &
psf,mslp,vort500,temp2m,dewp2m,qv2m, &
u10m,v10m,hgtsfc,refl1km,refl4km,cmpref, &
accppt,convppt,cape,mcape,cin,mcin,lcl, &
srh01,srh03,uh25,sh01,sh06,thck,li,brn,pw)
!end if
! CALL getunit( nchout )
OPEN (UNIT=4,FILE=trim(filename)//"_ready")
WRITE (4,'(a)') trim(filename)//"_ready"
CLOSE (4)
! CALL retunit (nchout )
ENDDO ! ntime loop
STOP
9000 CONTINUE ! I/O errors
CLOSE(4)
PRINT *, 'IO ERRORS'
STOP
950 CONTINUE
WRITE(6,'(a)') ' Error setting up GEMPAK file'
STOP
END PROGRAM bin2gem