c     Program to read *.amc file and put all atomic parameters into arrays.
c     The read_amc subroutine reads unit 4 (a *.amc file) and returns with
c        a character string containing the name of the compound (name)
c        the number of lines containing authors (num_authors)
c        a character array of lines of author names (author)
c        a text string containing the reference to the paper (ref)
c        the number of lines containing the title (num_title)
c        a character array of lines containing the title of the paper (title)
c        a real array containing the cell parameters (cell)
c        a character string containing the space group (spgroup)
c        if tranflag=.true. then transl(3) contains non-standard space group
c           translation
c        if rotflag=.true. then rotate(3,3) contains non-standard space group
c           rotation
c        the number of atoms (numatoms)
c        a character array of atom names (atom)
c        a real 2-dim array of atomic parameters whose contents are described in
c           the data statement below (parm):
c      data f/'x','y','z','occ','Biso','Uiso',
c     +       'B(1,1)','B(2,2)','B(3,3)','B(1,2)','B(1,3)','B(2,3)',
c     +       'U(1,1)','U(2,2)','U(3,3)','U(1,2)','U(1,3)','U(2,3)'/
c        a logical array whose entries are only .true. if a parameter exists for
c           the specified field.

c	Change on May 6, 2003, to read any amc file with references of format
c	journal name xx (1999) xx-yy

	implicit none
	integer maxatoms,num_author,num_title,numatoms,i,j,n,nn
      parameter (maxatoms=200)
      real cell(6) ,parm(maxatoms,18) ,transl(3),rotate(3,3)
      character name*100,author(3)*100,ref*100,title(5)*100,spgroup*100
      character atom(maxatoms)*20 ,line*100 ,temp*20
      logical parmL(18),tranflag,rotflag

      call read_amc(name,num_author,author,ref,num_title,
     &   title,cell,spgroup,numatoms,atom,parm,parmL,
     &   tranflag,rotflag,transl,rotate)

      write(7,'(a)')name(1:len_trim(name))
      do i=1,num_author
         write(7,'(a)') author(i)(1:len_trim(author(i)))
      enddo
      write(7,'(a)') ref(1:len_trim(ref))
      do i=1,num_title
         write(7,'(a)') title(i)(1:len_trim(title(i)))
      enddo
      write(7,'(6f10.5,3x,a)')cell,spgroup(1:len_trim(spgroup))
      if(tranflag)write(7,'(3f10.5)')(transl(i),i=1,3)
      if(rotflag)then
         do i = 1,3
            write(7,'(3f10.5)')(rotate(i,j),j=1,3)
         enddo
      endif
      do i=1,numatoms
         n=len_trim(atom(i))
		line=' '
         line(1:n)=atom(i)(1:n)
         do j=1,18
            if(parmL(j))then
               write(temp,'(f10.5)')parm(i,j)
               nn=len_trim(temp)
               line(n+1:n+nn)=temp(1:nn)
               n=n+nn
            endif
         enddo
         write(7,'(a)')trim(line)
      enddo
      stop
      end
c******************************************************************************
      subroutine read_amc(name,num_author,author,ref,num_title,
     &   title,cell,spgroup,numatoms,atom,parm,parmL,
     &   tranflag,rotflag,transl,rotate)
	implicit none
      integer nf(18,2),maxatoms,num_author,num_title,numatoms,io,n,i
	integer nwidth,npos
      parameter (maxatoms=200)
      real cell(6),parm(maxatoms,18),transl(3),rotate(3,3)
      character name*100 ,author(3)*100 ,line*100 ,ref*100 ,title(5)*100
      character spgroup*100 ,word*100 ,atom(maxatoms)*20 ,form*7
      logical is_cell ,parmL(18)
      logical tranflag,rotflag,CheckIfReference

      spgroup=' '
      tranflag=.false.
      rotflag=.false.
      read(4,'(a)',end=999)name                             ! mineral name
      read(4,'(a)',end=999)author(1)                        ! authors
      num_author=1
100   read(4,'(a)')line
	if(CheckIfReference(line))then				! reference
         ref=line
      else
         num_author=num_author+1
         author(num_author)=line
         goto 100
      endif
      read(4,'(a)')line
      num_title=0
      is_cell=.false.
      do while(.not.is_cell)                                 ! title lines
         num_title=min(num_title+1,5)
         title(num_title)=line
         read(4,'(a)')line
         call check_if_cell(line,is_cell,cell,spgroup)       ! cell parameters
      enddo
      if(spgroup(1:1).eq.'*')then
         call get_nonstandardsg(transl,rotate,tranflag,rotflag,line)
      else
         read(4,'(a)')line                            ! atom parameter header
      endif
c     Now examine atom parameter header
catom  x  y  z
      call get_atom_format(line,nf,parmL)
c     Now read atom parameters till END of file      ! atom parameters
      numatoms=0
      read(4,'(a)',iostat=io)line
      do while(io.eq.0)
         numatoms=numatoms+1
         call get_columns(line,n)
         atom(numatoms)=line(1:n)
         call trim_word(line) ! remove atom name
         do i=1,3                                    ! now get the atom coord
            if(nf(i,1).ne.0)then
               call readcoor(i,nf,parm,line,numatoms)
            endif
         enddo
         do i=4,18                                   ! now get the rest
            if(nf(i,1).ne.0)then
               nwidth=nf(i,2)-nf(i,1)+1
               if(nwidth.lt.10)then
                  write(form,1)'(f',nwidth,'.0)'
1                 format(a2,i1,a3)
               else
                  write(form,2)'(f',nwidth,'.0)'
2                 format(a2,i2,a3)
               endif
               word(1:nwidth)=line(nf(i,1):nf(i,2))
				npos=index(word,'/') ! parse for fractions
				if(npos.ne.0)then
					call convfrac(nwidth,npos,word,parm,numatoms,i)
				else
					read(word,form)parm(numatoms,i)
				endif
               if(i.eq.4.and.abs(parm(numatoms,4)).lt. 0.0001)    ! occ
     +            parm(numatoms,4)=1.0
            endif
         enddo
         read(4,'(a)',iostat=io)line
         if(ichar(line(1:1)).lt.33)io=-1
         if(line(1:3).eq.'END'.or.line(1:3).eq.'end'.or.
     &      line(1:3).eq.'End')io=-1
      enddo
      return
999   stop
      end
c******************************************************************************
c     read non-standard space group info
      subroutine get_nonstandardsg(transl,rotate,tranflag,rotflag,line)
	implicit none
	integer i,j,n,itm,ist
      logical tranflag,rotflag
      real transl(3),rotate(3,3)
      character line*100 ,nonstand(4)*100 ,temp*100 ,form*20

      tranflag=.false.
      rotflag=.false.
      i=0
      read(4,'(a)')line
      do while(line(1:4).ne.'atom')
         i=i+1
         nonstand(i)=line
         read(4,'(a)')line
      enddo
      if(i.eq.1)tranflag=.true.
      if(i.eq.4)tranflag=.true.
      if(i.eq.3)rotflag=.true.
      if(i.eq.4)rotflag=.true.
      if(i.eq.1.or.i.eq.4)then                     ! read translation
         temp=nonstand(1)
         do j=1,3
            call shift_left(temp)
            call get_columns(temp,n)
            if(n.lt.1)return
            if(n.lt.10)then
               write(form,1)'(F',n,'.0)'
1              format(a2,i1,a3)
            else
               write(form,2)'(F',n,'.0)'
2              format(a2,i2,a3)
            endif
            read(temp,form,err=999)transl(j)
            call trim_word(temp)
         enddo
      endif
      if(i.eq.3.or.i.eq.4)then                     ! read rotation
         itm=0
         do ist=i-2,i
            itm=itm+1
            temp=nonstand(ist)
            do j=1,3
               call shift_left(temp)
               call get_columns(temp,n)
               if(n.lt.1)return
               if(n.lt.10)then
                  write(form,1)'(F',n,'.0)'
               else
                  write(form,2)'(F',n,'.0)'
               endif
               read(temp,form,err=999)rotate(itm,j)
               call trim_word(temp)
            enddo
         enddo
      endif

999   return
      end
c******************************************************************************
c     Array nf contains the beginning and ending column for each type of header
      subroutine get_atom_format(line,nf,parmL)
	implicit none
      character f(18)*6 ,line*100
      integer nf(18,2),i,n,j
      logical dummyL ,parmL(18)
      data f/'x','y','z','occ','Biso','Uiso',
     +       'B(1,1)','B(2,2)','B(3,3)','B(1,2)','B(1,3)','B(2,3)',
     +       'U(1,1)','U(2,2)','U(3,3)','U(1,2)','U(1,3)','U(2,3)'/
      do i=1,18
         nf(i,1)=0
         nf(i,2)=0
         n=len_trim(f(i))
         nf(i,2)=index(line,f(i)(1:n))
         if(nf(i,2).ne.0)nf(i,2)=nf(i,2)+n-1
         parmL(i)=.true.
         if(nf(i,2).ne.0)then
            dummyL=.false.
            do j=nf(i,2),1,-1
               if(line(j:j).eq.' '.and. .not.dummyL) dummyL=.true.
               if(line(j:j).ne.' '.and. dummyL) then
                  nf(i,1)=j+1
                  goto 999
               endif
            enddo
         else
            parmL(i)=.false.
         endif
999   enddo
      nf(1,1)=1
c      do i=1,18
c         if(nf(i,1).ne.0)write(*,*)nf(i,1),nf(i,2)
c      enddo
      return
      end
c******************************************************************************
      subroutine check_if_cell(line,is_cell,cell,spgroup)
	implicit none
	integer i,n
      real cell(6)
      character line*100 ,temp*100 ,spgroup*100 ,form*20
      logical is_cell

      temp=line
      do i=1,6
         call shift_left(temp)
         call get_columns(temp,n)
         if(n.lt.1)return
         if(n.lt.10)then
            write(form,1)'(F',n,'.0)'
1           format(a2,i1,a3)
         else
            write(form,2)'(F',n,'.0)'
2           format(a2,i2,a3)
         endif
         read(temp,form,err=999)cell(i)
         call trim_word(temp)
      enddo
      call shift_left(temp)
      call get_columns(temp,n)
      spgroup=temp(1:n)
      call trim_word(temp)
      call shift_left(temp)
      call get_columns(temp,n)
      if(n.eq.0)is_cell=.true.

999   return
      end
c******************************************************************************
      subroutine get_columns(line,n) ! number of character spaces for the word
	implicit none
	integer i,n
      character line*100
      do i=1,100
         if(line(i:i).eq.' ' .or. ichar(line(i:i)).eq.0 .or.
     +      ichar(line(i:i)).eq.255)then
            n=i-1
            return
         endif
      enddo
      n=100
      return
      end
********************************************************************************
      subroutine shift_left(line)
	implicit none
	integer i
      character line*100
100   if(line(1:1).eq.' ')then
         i=len_trim(line)
         if(i.eq.0)return
         line(1:i-1)=line(2:i)
         line(i:i)=' '
         goto 100
      endif
      return
      end
********************************************************************************
      subroutine trim_word(line) ! replace first word of a line with blanks
	implicit none
	integer i
      character line*100
      do i=1,100
         if(line(i:i).eq.' ')return
         line(i:i)=' '
      enddo
      return
      end
********************************************************************************
      subroutine readcoor(i,nf,parm,line,numatoms)
	implicit none
	integer maxatoms,nf(18,2),npos,i,nwidth,numatoms
      parameter (maxatoms=200)
      real parm(maxatoms,18)
      character line*100 ,word*100 ,form*7

      nwidth=nf(i,2)-nf(i,1)+1
      if(nwidth.lt.10)then
         write(form,1)'(f',nwidth,'.0)'
      else
         write(form,2)'(f',nwidth,'.0)'
      endif
      word=' '
      word(1:nwidth)=line(nf(i,1):nf(i,2))
c     check if word contains / in case 1/2,1/3,1/4 etc have been input.
      npos=index(word,'/')
      if(npos.ne.0)then
		call convfrac(nwidth,npos,word,parm,numatoms,i)
      else
         read(word,form)parm(numatoms,i)
      endif
1     format(a2,i1,a3)
2     format(a2,i2,a3)
      return
      end
********************************************************************************
	subroutine convfrac(nwidth,npos,word,parm,numatoms,i)
	implicit none
	integer maxatoms
      parameter (maxatoms=200)
      integer nf(18,2),npos,n,nwidth,numatoms,i
      real parm(maxatoms,18),xu,xl
      character line*100 ,word*100 ,form*7
      n=npos-1
      if(n.lt.10)then                                ! read numerator
         write(form,1)'(f',n,'.0)'
      else
         write(form,2)'(f',n,'.0)'
      endif
      read(word(1:n),form)xu

      n=n+2
      if(n.lt.10)then                                ! read denominator
         write(form,1)'(f',n,'.0)'
      else
         write(form,2)'(f',n,'.0)'
      endif
      read(word(n:nwidth),form)xl

      if(abs(xl).gt.0.0001)then
         parm(numatoms,i)=xu/xl
      else
         parm(numatoms,i)=0.0
      endif
1     format(a2,i1,a3)
2     format(a2,i2,a3)
	return
	end
********************************************************************************
	logical function CheckIfReference(line)
	character line*100
	integer i
	CheckIfReference=.true.
	if(line(1:5).eq.'COD, ')return
	if(line(1:18).eq.'TSL Structure File')return
!  after this point the line must contain typical reference, line Am Min 88 (19xx) 123-124
	CheckIfReference=.false.
	i=index(line,' (19')	! check if string (1920) etc exists
	if(i.le.0)i=index(line,' (20')	! check if string (2002) etc exists
	if(i.le.0)return
	if(line(i+6:i+7).ne.') ')return
	CheckIfReference=.true.
	end function CheckIfReference
!********************************************************************************