program tbase_convert

! convert tbase.dat from binary to netcdf
! xlf90 -q64 -qsuffix=f=f90:cpp=F90 -O2 -I/usr/local/include -L/usr/local/lib -lnetcdf tbase_conver.F90


    implicit none

    include 'netcdf.inc'

    integer, parameter :: nlon = 4320, nlat=2161
    integer            :: itopo(nlon,nlat)
    real*8             :: latc(nlat),lonc(nlon),topo(nlon,nlat)
    real*8             :: edges,edgen,edgew,edgee
    real*8             :: dx,dy

    integer            :: i,j,im,jm
    integer            :: ncid,ierr,dimid(2)
    integer            :: nid,sid,eid,wid,latid,lonid,topoid
    character*128      :: str

    open(20,file='tbase.dat',form='unformatted')
    read(20) itopo
    close(20)

    edges =  -90.0
    edgen =   90.0
    edgew =    0.0
    edgee =  360.0

    dx = (edgee - edgew) / float(nlon)
    dy = (edgen - edges) / float(nlat)

! set lats/lons
! change north-south to south-north in itopo data
    do j = 1, nlat
    do i = 1, nlon
      lonc(i) = edgew + (2*i-1)*dx / 2.
      latc(j) = edges + (2*j-1)*dy / 2
      im = i
      jm = nlat - j + 1
      topo(i,j) = float(itopo(im,jm))
    enddo
    enddo

    write(6,*) 'min/max lat',minval(latc),maxval(latc)
    write(6,*) 'min/max lon',minval(lonc),maxval(lonc)
    write(6,*) 'min/max top',minval(topo),maxval(topo)

    ierr = nf_create('tbase.dat.nc',nf_clobber,ncid)
    ierr = nf_def_dim(ncid,'ni',nlon,dimid(1))
    ierr = nf_def_dim(ncid,'nj',nlat,dimid(2))

    str = 'conversion of tbase.dat from binary to netcdf'
    ierr = nf_put_att_text(ncid,NF_GLOBAL,'source',len(trim(str)),trim(str))

    ierr = nf_def_var(ncid,'EDGEN',nf_double,0,0,nid)
    ierr = nf_def_var(ncid,'EDGES',nf_double,0,0,sid)
    ierr = nf_def_var(ncid,'EDGEE',nf_double,0,0,eid)
    ierr = nf_def_var(ncid,'EDGEW',nf_double,0,0,wid)
    ierr = nf_def_var(ncid,'lon',nf_double,1,dimid(1),lonid)
    ierr = nf_def_var(ncid,'lat',nf_double,1,dimid(2),latid)
    ierr = nf_def_var(ncid,'TOPO',nf_double,2,dimid,topoid)

    write(6,*) 'ids ',dimid,nid,sid,eid,wid,latid,lonid,topoid

    str = 'longitude'
    ierr = nf_put_att_text(ncid,lonid,'long_name',len(trim(str)),trim(str))
    str = 'degrees east'
    ierr = nf_put_att_text(ncid,lonid,'units',len(trim(str)),trim(str))

    str = 'latitude'
    ierr = nf_put_att_text(ncid,latid,'long_name',len(trim(str)),trim(str))
    str = 'degrees north'
    ierr = nf_put_att_text(ncid,latid,'units',len(trim(str)),trim(str))

    str = 'topographic elevation'
    ierr = nf_put_att_text(ncid,topoid,'long_name',len(trim(str)),trim(str))
    str = 'm'
    ierr = nf_put_att_text(ncid,topoid,'units',len(trim(str)),trim(str))

    ierr = nf_enddef(ncid)

    ierr = nf_put_var1_double(ncid,nid,1,edgen)
    ierr = nf_put_var1_double(ncid,sid,1,edges)
    ierr = nf_put_var1_double(ncid,eid,1,edgee)
    ierr = nf_put_var1_double(ncid,wid,1,edgew)

    ierr = nf_put_var_double(ncid,lonid,lonc)
    ierr = nf_put_var_double(ncid,latid,latc)
    ierr = nf_put_var_double(ncid,topoid,topo)

    ierr = nf_sync(ncid)
    ierr = nf_close(ncid)

    end