anchor.diff:           allocate(automorphism(n2))
            call find_permutation_of_group(g1,g2,automorphism)
            !3) change ilabeling2 into the equivalent ilabeling1
  +         write(*,'("automorphism: ",2x,7i1)') automorphism
            write(*,'("before:  ",14i1)') ilabeling2
            do i = 0, n*nD2-1, n
   !            write(*,'("i  ",i2," vs: ",7i2)')i,(automorphism+i)
               ilabeling2(1+i:n+i) = ilabeling2(automorphism+i)
            enddo
  -         write(*,'("before:  ",14i1)') ilabeling2
  +         write(*,'("after:   ",14i1)') ilabeling2
            !4) proceed as normal
         endif
            foundLab = .false.
  -         nL = size(pLabel,1)
  -         outer: do iL = 1, nL
  +!July 29         nL = size(pLabel,1)
  +         !outer: !do iL = 1, nL
  +            ! GLWH July 2019. As of today, this logic is wrong from earlier attempts. We do not want to loop over all permutations of both the file1 and file2 label. It is sufficient, now that we have the self-complementary equiatomic stuff sorted out, to just check each file2 labeling with all permutations of the file1 labeling (or the other way around.)
  +
  +
               ! We need to loop over all permutations (to get the full orbit) because we do not know
               ! which member of the orbit is the representative.
               ! Also, for 50-50 labelings that are self-complementary (more generally, any
  
positive.code:    iStr2 = 0
     ! You can get a subtle bug here, hard to chase down, if you have another module opening files with the same unit number. Make sure that the unit numbers for files opened by this program are unique.  GLWH 2019 (that's why f2 is set to 14 above)
     compareloop: do
        write(*,'("Compare loop: file1 strN ",i6)') strN
        write(13,'("Compare loop: file1 strN ",i6)') strN
        iStr2 = iStr2 + 1
        !write(*,'("iStr2: ",i4)') iStr2
        read(f2,*,iostat=ioerr) strN2, hnfN2, hdgenfact, labdgen, totdgen, sizeN, n2, pgOps2, diag2, a,b,c,d,e,f, L2(:,:,1), labeling
        write(*,'("read from file 2: strN2:",i5)') strN2
        if(ioerr/=0) exit
        read(labeling,'(500i1)') ilabeling2
        L2(:,:,1) = transpose(L2(:,:,1)) ! Written out column-wise but read in row-wise. So fix it by transposing
        ! Check to see if the two files have different left transform matrices. If they do, they may !  still define equivalent structures so we'll have to do some work.
        if (any(L1/=L2)) then
           !make the qlabel list of permutations
           call get_dvector_permutations(pLV2,dset2,nD2,rot,shift,dRotList2,LatDim2,eps)
           call remove_duplicate_lattices(HNFin,LatDim2,pLV2,rot,shift,dset2,dRotList2,HNFout,fixOp,&
  
negative.code: ! color is exactly the same. If so, we need some extra logic when two labeling are compared
  ! (see email discussion with Rod Forcade on May 16 2019, subject: "Monk?")
     equiatomic = .false.
     if (mod(n*nD1,k)== 0) then ! necessary but not sufficient condition for equiatomic.
        equiatomic = .true.
        equiN = n*nD1/k
        allocate(integerLab(n*nD1))
        do j = 1, n*nD1
           read(labeling(j:j),'(i1)') integerLab(j)
        enddo
        write(*,'("IntLab",14i1)') integerLab
        print*,count(integerLab==0)
        do i = 0,k-1
        if (count(integerLab==i)/=equiN) then ! can't be equiatomic
           print*,"Not equiatomic"
           equiatomic = .false.
           exit
        endif
        if (equiatomic) print*,"Equiatomic"
        enddo
     endif
  
  