      subroutine ahopen(luin,filename,mode,istat)
c---- ahopen opens an ah file for reading
c
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      character*(*) filename
      logical exists
      data luopen/nunits*-1/
      data iaddhead/nunits*0/
      data iadddata/nunits*0/
c
      istat=0
      if(luopen(luin).ne.-1) then
        istat=5
      else
        if(mode.eq.0) then
          inquire(file=filename,exist=exists)
          if(exists) then
            lassign=iopen(filename,0)
            luopen(luin)=lassign
            write(6,"('opening ',i3,' to file ',a40)")lassign,filename
          else
            write(6,"('no such file ',a40)")filename
            istat=4
          endif
        else if(mode.eq.1) then
          ierr=icreat(filename,420)
          lassign=iopen(filename,1)
          luopen(luin)=lassign
          write(6,"('opening ',i3,' to file ',a40)")lassign,filename
        else if(mode.eq.2) then
          inquire(file=filename,exist=exists)
          if(.not.exists) then
            ierr=icreat(filename,0644)
          endif
          lassign=iopen(filename,2)
          luopen(luin)=lassign
          write(6,"('opening ',i3,' to file ',a40)")lassign,filename
        else
          write(6,"('cannot open file with r/w mode',i4)") mode
          istat=4
        endif
      endif
      return
      end
c
      subroutine ahclose(luin,istat)
c---- ahopen opens an ah file for reading
c
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      istat=0
      if(luopen(luin).eq.-1) then
        istat=5
      else
        luopen(luin)=-1
        iaddhead(luin)=0
        iadddata(luin)=0
        close(luin)
      endif
      return
      end
c
      subroutine ahrewind(luin,istat)
c---- ahrewind rewinds an ah file
c
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      istat=0
      if(luopen(luin).eq.-1) then
        istat=5
      else
        iaddhead(luin)=0
        iadddata(luin)=0
      endif
      return
      end
c      
      subroutine rdahhead(luin,cod,chan,astla,astlo,astel,
     # jyr,jmo,jda,jho,jmi,asec,kyr,kmo,kda,kho,kmi,esec,dt
     #,nda,istat)
c---- rdahhead - routine to read an ah-header 
      include 'ahhead_common'
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      character*6 cod,chan
      istat=0
      ifd=luopen(luin)
      if(ifd.eq.-1) then
        write(6,"('logical unit',i3,' not open')")luin
        istat=5    
      else
        i1=iaddhead(luin)
        ierr=iread(ifd,ahheader,i1,1024)
c        if(ierr.ne.1024) go to 102
        if(ierr.ne.1024) go to 102
        read(codeah,"(a6)")cod
        read(chanah,"(a6)")chan
        call c2fstr(cod,6)
        call c2fstr(chan,6)
        iadddata(luin)=iaddhead(luin)+RECW4*4
        iaddhead(luin)=iaddhead(luin)+4*ndatah+RECW4*4
c
        dt=deltah
        nda=ndatah
        jyr=isyrah
        jmo=ismoah
        jda=isdaah
        jho=ishrah
        jmi=ismnah
        asec=ssecah
        astla=slatah
        astlo=slonah
        astel=elevah
        kyr=ieyrah
        kmo=iemoah
        kda=iedaah
        kho=iehrah
        kmi=iemnah
        esec=esecah
        return
  102   continue
        write(6,"('end of file')")
        istat=3
        return
  103   continue
        write(6,"('error')")
        istat=5
      endif
      return
      end
      subroutine rdahhead1(luin,cod,chan,astla,astlo,astel,
     # jyr,jmo,jda,jho,jmi,asec,elate,elone,kyr,kmo,kda,kho,kmi,esec,dt
     #,nda,istat)
c---- rdahhead - routine to read an ah-header 
      include 'ahhead_common'
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      character*6 cod,chan
      istat=0
      ifd=luopen(luin)
      if(ifd.eq.-1) then
        write(6,"('logical unit',i3,' not open')")luin
        istat=5    
      else
        i1=iaddhead(luin)
        ierr=iread(ifd,ahheader,i1,1024)
c        if(ierr.ne.1024) go to 102
        if(ierr.ne.1024) go to 102
        read(codeah,"(a6)")cod
        read(chanah,"(a6)")chan
        call c2fstr(cod,6)
        call c2fstr(chan,6)
        iadddata(luin)=iaddhead(luin)+RECW4*4
        iaddhead(luin)=iaddhead(luin)+4*ndatah+RECW4*4
c
        dt=deltah
        nda=ndatah
        jyr=isyrah
        jmo=ismoah
        jda=isdaah
        jho=ishrah
        jmi=ismnah
        asec=ssecah
        astla=slatah
        astlo=slonah
        astel=elevah
	print*,'astla astlo astel',astla,astlo,astel
        kyr=ieyrah
        kmo=iemoah
        kda=iedaah
        kho=iehrah
        kmi=iemnah
        esec=esecah
	elate=elatah
	elone=elonah
        return
  102   continue
        write(6,"('end of file')")
        istat=3
        return
  103   continue
        write(6,"('error')")
        istat=5
      endif
      return
      end
      subroutine rdahdata(luin,data,nread,istat)
c---- rdahhead - routine to read ah-data 
      include 'ahhead_common'
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      dimension data(1)
      istat=0
      ifd=luopen(luin)
      if(ifd.eq.-1) then
        write(6,"('logical unit',i3,' not open')")luin
        istat=5    
      else if(iadddata(luin).eq.0) then
        write(6,"('header must be read before data')")
        istat=5
      else
        i1=iadddata(luin)
        ntoread=nread*4
c        if(ntoread.gt.16384) then
c          ierr=0
c          iloop=ntoread/16384
c          nremain=ntoread
c          do 10 i=1,iloop+1
c            write(6,"('loop',i4,' bytes remaining',i6)") i,nremain
c            iindex=1+(i-1)*4096
c            if(nremain.gt.16384) then
c              ierr1=iread(ifd,data(iindex),i1,16384)
c              if(ierr1.ne.16384) go to 102
c              i1=i1+16384
c            else if(nremain.le.0) then
c               write(6,"('no more to read')")
c            else
c              ierr1=iread(ifd,data(iindex),i1,nremain)
c              if(ierr1.ne.nremain) go to 102
c              i1=i1+nremain
cc            endif
c            ierr=ierr+ierr1
c            nremain=ntoread-16384*i
c   10     continue
c        else
          ierr=iread(ifd,data,i1,ntoread)
c       print*,'data',(data(kk),kk=1,nread)
c        endif
        if(ierr.ne.ntoread) go to 102
        return
  102   continue
        istat=3
        write(6,"('end of file')")
        return
  103   continue
        write(6,"('error')")
        istat=5
      endif
      return
      end
c
      subroutine wrahhead(lu,istat,dt,nda,ida,ih0,im0,s0,alat,alon,stel)
      include 'ahhead_common'
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
	slatah=alat
	slonah=alon
	elevah=stel
	isdaah=ida
	ishrah=ih0
	ismnah=im0
	ssecah=s0
	ndatah=nda
	deltah=dt
      istat=0
      ifd=luopen(lu)
      if(ifd.eq.-1) then
        write(6,"('logical unit',i3,' not open')")lu
        istat=5    
      else
        i1=iaddhead(lu)
        i2=iadddata(lu)
        if(i1.ge.i2) then
          ierr=iwrite(ifd,ahheader,i1,1024)
          if(ierr.ne.1024) then
            istat=4
          else
c            iadddata(lu)=iaddhead(lu)+1024
            iadddata(lu)=iaddhead(lu)+RECW4*4
          endif
        else
          write(6,"('i1,i2',2I4)")i1,i2
          write(6,"('cannot write another header before data')")
          istat=4
        endif
      endif
      return
      end
      subroutine wrahdata(lu,data,nwrit,istat)
      include 'ahhead_common'
      parameter (nunits=255)
      common/ahinou/ luopen(nunits),iaddhead(nunits),
     #           iadddata(nunits)
      save /ahinou/
      dimension data(1)
      istat=0
      ifd=luopen(lu)
      if(ifd.eq.-1) then
        write(6,"('logical unit',i3,' not open')")lu
        istat=5    
      else
        i1=iaddhead(lu)
        i2=iadddata(lu)
        if(i2.ge.i1) then
          ntowrite=nwrit*4
c          if(ntowrite.gt.1116384) then
c            ierr=0
c            iloop=ntowrite/16384
c            nremain=ntowrite
c            do 10 i=1,iloop+1
c              write(6,"('loop',i4,' bytes remaining',i6)") i,nremain
c              iindex=1+(i-1)*4096
c              if(nremain.gt.16384) then
c                ierr1=iread(ifd,data(iindex),i1,16384)
c                if(ierr1.ne.16384) go to 102
c                i1=i1+16384
c              else if(nremain.eq.0) then
c              else
c                ierr1=iread(ifd,data(iindex),i1,nremain)
c                if(ierr1.ne.nremain) go to 102
c                i1=i1+nremain
c              endif
c              ierr=ierr+ierr1
c              nremain=ntoread-16384*i
c   10       continue
c          else
c
c---- this should not work since ntowrite is > int*2 --- but it works
            ierr=iwrite(ifd,data,i2,ntowrite)
c          endif
          if(ierr.ne.ntowrite) then
            istat=4
          else
            iaddhead(lu)=iadddata(lu)+ntowrite
          endif
        else
          write(6,"('cannot write more data before next header')")
          istat=4
        endif
      endif
      return
c  102 continue
c      istat=4
c      return
      end
c
      subroutine getahpz(ds,a0,np,nz,pole,zero)
      include 'ahhead_common'
      dimension pole(2,1),zero(2,1)
      ds=dsah
      a0=a0ah
      np=int(xnpah+0.5)
      nz=int(xnzah+0.5)
      do 10 i=1,np
        pole(1,i)=pz(1,1,i)
        pole(2,i)=pz(2,1,i)
   10 continue
      do 20 i=1,nz
        zero(1,i)=pz(1,2,i)
        zero(2,i)=pz(2,2,i)
   20 continue
      return
      end
c
      subroutine putahpz(ds,a0,np,nz,pole,zero)
      include 'ahhead_common'
      dimension pole(2,1),zero(2,1)
      dsah=ds
      a0ah=a0
      xnpah=float(np)
      xnzah=float(nz)
      do 10 i=1,np
        pz(1,1,i)=pole(1,i)
        pz(2,1,i)=pole(2,i)
   10 continue
      do 20 i=1,nz
        pz(1,2,i)=zero(1,i)
        pz(2,2,i)=zero(2,i)
   20 continue
      return
      end
c
      subroutine c2fstr(string,len)
      character*(*) string
      do 1 i=1,len
        if(string(i:i).eq.'\0') then
          string(i:i)=' '
        endif
    1 continue
      return
      end
      
