program read_incommensurate_irs
* read the data file containing the IRs
      implicit none
      integer i,j,k,m,n,ng
      real sign
* space group number (1-230)
      integer ispacegroup(5508)
* space group symbol: P1, P-1, P2, etc.
      character spacegroupstring(5508)*10
* IR symbol: GM1, DT1, etc.
      character irstring(5508)*8
* k vector symbol: GM, DT, etc.
      character kstring(5508)*2
* k vector string:  (a,b,g), (1/2,0,g), etc.
      character kvectorstring(5508)*12
* k vector components:  
* kvector(i,1,n) is the ith rational component
* kvector(i,2,n) is the ith alpha component
* kvector(i,3,n) is the ith beta component
* kvector(i,4,n) is the ith gamma component
      real kvector(3,4,5508)
* H image string: A1aA1a, A1aA2a, etc.
      character imagestring(5508)*12
* dimension of IR
      integer irdim(5508)
* IR matrix for translation delta along the a4 axis
* 1=0
* 2=cos(2*pi*delta)
* 3=-cos(2*pi*delta)
* 4=sin(2*pi*delta)
* 5=-sin(2*pi*delta)
      integer matrixdelta(8,8,5508)
* number of generators of H group
      integer ngenerators(5508)
* affine transformation for each generator of H group
      real generator(5,5,4,5508)
* IR matrix for each generator of H group
      real genmatrix(8,8,4,5508)


* open data file
      open(30,file='incommensurate_irs_data.txt')
* skip heading
      read(30,*)
      read(30,*)
      read(30,*)
* read each IR
      do m=1,5508
        read(30,*)n,ispacegroup(n),spacegroupstring(n),irstring(n),
     $       kstring(n),kvectorstring(n),
     $       ((kvector(i,j,n),i=1,3),j=1,4),imagestring(n),irdim(n),
     $       ((matrixdelta(i,j,n),j=1,irdim(n)),i=1,irdim(n)),
     $       ngenerators(n),
     $       (((generator(i,j,k,n),j=1,5),i=1,5),
     $       ((genmatrix(i,j,k,n),j=1,irdim(n)),i=1,irdim(n)),
     $       k=1,ngenerators(n))
* fix fractions that have been rounded off
* k vector
        do j=1,4
          do i=1,3
            ng=nint(kvector(i,j,n)*12)
            kvector(i,j,n)=1.0*ng/12
          enddo
        enddo
* generators of H group
        do k=1,ngenerators(n)
          do j=1,5
            do i=1,5
              ng=nint(generator(i,j,k,n)*24)
              generator(i,j,k,n)=1.0*ng/24
            enddo
          enddo
        enddo
* generating matrices of H group
        do k=1,ngenerators(n)
          do j=1,irdim(n)
            do i=1,irdim(n)
              if(abs(genmatrix(i,j,k,n)).lt.0.01)then
                genmatrix(i,j,k,n)=0
              else
                sign=genmatrix(i,j,k,n)/abs(genmatrix(i,j,k,n))
                ng=nint(genmatrix(i,j,k,n)**2*12)
                genmatrix(i,j,k,n)=sign*sqrt(1.0*ng/12)
              endif
            enddo
          enddo
        enddo
      enddo
      end