program readsnap
      implicit none

      character(len=1000) :: filename
      character(len=4) :: label
      integer, dimension(6) :: nparttab
      double precision :: redshift,time
      double precision, dimension(6) :: masstab
      real,    dimension(:,:), allocatable :: pos,vel
      integer, dimension(:),   allocatable :: id
      integer :: size,npart,i

      call get_command_argument(1,filename)
      open(1,file=filename,status='old',form='unformatted')
      read(1,end=2) label,size
      write(0,*) label,size-8
      read(1) nparttab,masstab,redshift,time
      npart = nparttab(2)
      write(0,'(i6," particles")') npart
      allocate(id(npart),pos(3,npart),vel(3,npart))
    1 continue
      read(1,end=2) label,size
      write(0,*) label,size-8
      if(label.eq.'ID  ') THEN
        read(1) id
      elseif(label.eq.'POS ') THEN
        read(1) pos
      elseif(label.eq.'VEL ') THEN
        read(1) vel
      else
        read(1)
      endif
      goto 1
    2 continue
      close(1)

      write(*,'("# ",2f7.4)') redshift,time
      write(*,'(i6,6es12.4)') (id(i),pos(:,i),vel(:,i),i=1,npart)

      end