      program read_ssgmag_datafile
* read data for (3+d)-dimensional magnetic superspace groups (d=1,2,3)
      implicit none
      integer i,j,k,m,n,n4

* number of magnetic superspace groups
      integer nmgroups
* for each magnetic superspace group
*   associated nonmagnetic superspace group, same numbering as in
*   ssg_datafile.txt
      integer mgroup_ssg(325127)
*   basic nonmagnetic space group
      integer mgroup_spacegroup(325127)
*   number of modulations
      integer mgroup_nmod(325127)
*   basic magnetic space group number
*   number is given by mgroup_spacegroup(i).mgroup_mag(i): 1.1, 2.1, 2.2, etc.
      integer mgroup_mag(325127)
*   group number label: 1.1.1.1.m1.1, 2,1,1,1.m2.1, etc.
      character mgroup_nlabel(325127)*21
*   group label
      character mgroup_label(325127)*63
*   number of operators
      integer mgroup_nops(325127)
*   time reversal (1=includes time reversal, -1=does not include time reversal)
      integer mgroup_ops_r(96,325127)
*   (d+1)x(d+1) augmented matrix for each operator
*   common denominator in element (d+1,d+1)
      integer*1 mgroup_ops(7,7,96,325127)
*   augmented transformation matrix from the standard setting of the
*   magnetic basic space group to the standard setting of the magnetic
*   superspace group.  This matrix is used to transform Wyckoff positions
*   in the magnetic basic space group to the setting of the magnetic
*   superspace group.
*   common denominator in element (4,4)
      integer*1 mgroup_transmag(4,4,325127)
* open data file
      open(30,file='ssgmag_datafile.txt')
* skip heading
      read(30,*)
      read(30,*)
      read(30,*)
* read number of magnetic superspace groups
      read(30,*)nmgroups
* read each magnetic superspace group
      do m=1,nmgroups
        read(30,*)n,mgroup_ssg(m),mgroup_spacegroup(m),mgroup_nmod(m),
     $       mgroup_mag(m)
        if(n.ne.m)then
          stop 'Read error in read_ssgmag_datafile'
        endif
        read(30,*)mgroup_nlabel(m),mgroup_label(m)
        read(30,*)mgroup_nops(m)
        n4=mgroup_nmod(m)+4
        read(30,*)(mgroup_ops_r(k,m),
     $       ((mgroup_ops(i,j,k,m),i=1,n4),j=1,n4),
     $       k=1,mgroup_nops(m))
        read(30,*)((mgroup_transmag(i,j,m),i=1,4),j=1,4)
      enddo
      close(30)
      end
