c test to read compact summary files produced for Mark Everett june 96
c
      implicit real*8 (a-h,o-z)
      common/phys/ temp, dens
      common/populations/ pop(348)
      dimension tempa(15),densa(13)
      data tempa/4.5d2,6.d2,7.8d2,
     1          1.d3,1.4d3,2.d3,3.d3,4.5d3,6.d3,7.8d3,
     1          1.d4,1.4d4,2.d4,3.d4,4.5d4/
      data densa/1.d-3,1.d-2,1.d-1,1.d0,1.d1,1.d2,
     $     1.d3,1.d4,1.d5,1.d6,1.d7,1.d8,1.d9/
c
      write(6,*) ' H2 populations from master equation solutions from'
      write(6,*) 
     $     ' Martin, Schwarz, and Mandy 1996, ApJ, 461, 265'
      write(6,*) ' '
      write(6,*) ' you may select temperature from the following'
      write(6,10) (tempa(i),i=1,15)
      write(6,*) ' '
      write(6,*) ' you may select density from the following'
      write(6,10) (densa(i),i=1,13)
 10   format(1x,5(1pd8.1))
c
      write(6,*) ' '
      write(6,*) ' what temperature and density?'
      read(5,*) temp, dens
c
      write(6,*) ' reading populations'
      call readpop
      write(6,*) ' '
c
      write(6,*) ' begin interactive request for population '
      write(6,*) ' '
 1    continue
      write(6,*) ' what (v,j) state? -- -1 0 to exit'
      read(5,*) iv,ij
      if(iv.lt.0) go to 2
      call vjcol(iv,ij,icol)
      if(icol.eq.-1 .or. icol.eq.349) go to 1
      write(6,3) iv,ij,icol,pop(icol)
 3    format(3i8,3x,1pd23.15)
      go to 1
 2    continue
      write(6,*) ' goodbye'
      stop
      end
c
      subroutine readpop
c
c get the populations from the master equation file column format
c
      implicit real*8 (a-h,o-z)
      character*15 name1
      common/phys/  temp, dens
      common/populations/ pop(348)
      pathname='msm'
      write(6,*) ' temperature and density', temp, dens
c
c prepare input filename
c
      call tname(temp,tt,it)
      call tname(dens,td,id)
      if(id.ge.0) then
         write(name1,98) pathname, tt,it, td,id
 98      format(a3,f3.1,'t',i1,'x',f3.1,'np',i1)
      else
         id=-id
         write(name1,99) pathname, tt,it, td,id
 99      format(a3,f3.1,'t',i1,'x',f3.1,'nm',i1)
      endif
      write(6,31) name1,temp,dens
 31   format(1x,a15,'   temp, dens ',f8.0,1pd8.1)
      open(unit=1,status='old',file=name1)
c
c read and check temperature, density
c
      read(1,30) temper,density
      if(dabs((density-dens)/dens).gt.1.d-4) then
         write(6,*) ' densities do not match', dens, density
         stop
      endif
      if(dabs((temper-temp)/temp).gt.1.d-4) then
         write(6,*) ' temperatures do not match', temp, temper
         stop
      endif
c
c read in the populations
c
      read(1,30) (pop(i),i=1,348)
 30   format(3(1pd23.15))
      close(unit=1)
c
      return
      end
c
      subroutine tname(t,tt,it)
      implicit real*8 (a-h,o-z)
c parse the temperature/density for the new dataset same
      at=dlog10(t)
      it=idint(at)
      aat=at-it
      tt=10.d0**aat
      return
      end
c
      subroutine vjcol(n,j,icol)
c  routine using general sindex to find the column index for populations
c
c  input
c  n      state vibrational quantum number 0, 14
c  j      state rotational quantum number 0, 38 or 1, 37
c  n=j=99, ind=1 for dissociation (used in other programs)
c
c  output
c  icol:  column index  1, 349
c     para:             1, 177
c     ortho:            178, 348 (171)
c     dissociation:     349
c
      ioff=177
      if(n.lt.99) then
         call sindex(n,j,icol)
         if(icol.lt.0) go to 1
         if(mod(j,2).eq.1) icol=icol+ioff
      else
         icol=349
      end if
 1    continue
      return
      end
c
      subroutine sindex(n,j,indvj)
c  this subroutine returns indvj [1-178] given 
c   n     vibrational quantum number
c   j     rotational quantum number
      dimension iv(39)
c      data iv/15,15,15,15,15,14,14,14,14,13,13,13,13,12,12,11,11,10,10,
      data iv/15,15,15,15,15,14,14,14,14,13,13,13,12,12,12,11,11,10,10,
     19,9,9,8,8,7,7,6,6,5,5,4,4,3,3,2,2,1,1,1/
      if( n.gt.14 .or. n.lt.0) go to 200
      if( j.gt.38 .or. j.lt.0) go to 200
      if( (n+1).gt.iv(j+1) ) go to 200
      istart=1
      if(mod(j,2).eq.1)istart=2
      icount=1
      istop=j
c avoid the uncertain convention of a single pass or not on negative first step
      if(istart.gt.istop)go to 20
      do 10 i=istart,istop,2
      icount=icount+iv(i)
   10 continue
   20 continue
      icount=icount+n
      indvj=icount
      return
 200  continue
      write(6,9998) n,j
 9998 format(1x,'n, j = ',2i5,' illegal v,j combination')
      indvj=-1
      return
      end
c
