       program testmodels
       open(33,file='CH4tit5exo',status='new')
       open(32,file='oT5atm',status='old')
       tdum=0.
       istart=0
       write(*,'("enter species id: 1-5")')
       read(*,*)ispecn
c above is to read ascii file on first call
1      read(32,*,end=99)x,y,z,tdum
       r=sqrt(x**2+y**2+z**2)
       theta=acos(z/r)
       phi=atan2(y,x)
       theta=theta*45./atan(1.)
       phi=phi*45./atan(1.)
c insert model sub call here
c required coord system is TSO (Titan-centered, referenced to Sun)
       call titexo(r,x,ispecn,xnn,istart)   
       write(33,*)tdum,xnn
       write(*,*)tdum,xnn
       tdum=tdum+1.
       go to 1
99     stop
       end
c***********************************************************
       subroutine titexo(r,x,ispecn,xnn,istart)
c required coord system is TSO (Titan-centered, Sun referenced x axis)
c Titan neutral exosphere models from Amsif, Dandouras, and Roelof
c (JGR in press 1995)
c
c Species identifiers:
c 1=CH4, 2=N2, 3=N, 4=H2, 5=H
c altitude range is 1600-8000 km
c
c densities in cm-3
c
       common/modelexo/dnexo(5,642)
       if(istart.gt.0)go to 11
           istart=1
101        format(a)
           open(30,file='titanexo3',status='old')
            do 10 jalt=1,642
              read(30,110)altdum,(dnexo(k,jalt),k=1,5)
110           format(6e16.7)
10          continue
          close(30)
11     continue
c      compute altitude and solar zenith angle for interpolation
       alt=2575.*r-2575.
       call exomod(denexo,alt,ispecn)
       xnn=denexo
       return
       end
C------------------------------------------------------------ 
      subroutine exomod(denexo,alt,ISPECN)
c     this program calculates the neutral densities according to
c the model of Amsif et al. by interpolation of their results.
      common/modelexo/dnexo(5,642)
      denexo=0.
      nalt=(alt-1600.)/10.+1
       if(alt.lt.1600..or.nalt.gt.641)go to 99
       stuff1=dnexo(ISPECN,nalt+1)
       stuff2=dnexo(ISPECN,nalt)
      denexo=stuff2+(alt-1600.-(float(nalt)-1.)*10.)/10.*(stuff1-stuff2)
99     continue
       return
       end
