anchor.diff:            write(6,*) 'read harm transform coords u ',minval(u%data(1)),minval(u%data(2))
             write(6,*) 'read harm transform coords u size ',size(u),hslope
   !          if(BL.eq.3) then
  -! HARD-CODED MMKS OR MKS
  -          uks = ummks2uks(u,dble(x1_arr),dble(x2_arr))
  -          bks = ummks2uks(b,dble(x1_arr),dble(x2_arr))
  -!          else
  -!          uks = umksh2uks(u,r_arr,x2_arr)
  -!          bks = umksh2uks(b,r_arr,x2_arr)
  -!          end if
  +          if(IS_MMKS.eq.1) then
  +             uks = ummks2uks(u,dble(x1_arr),dble(x2_arr))
  +             bks = ummks2uks(b,dble(x1_arr),dble(x2_arr))
  +          else
  +             uks = umksh2uks(u,r_arr,x2_arr)
  +             bks = umksh2uks(b,r_arr,x2_arr)
  +          end if
             write(6,*) 'after uks ',minval(uks%data(1))
             u = uks2ubl(uks,dble(r_arr),dble(asim))
             write(6,*) 'read harm transform coords b', minval(u%data(3)),maxval(u%data(2)),&
  
positive.code: 
        character(len=100) :: dfile, hfile, gfile
        integer :: n, ndumps, nt, indf, dlen, nhead
        integer :: ndim=3
        integer :: IS_MMKS=1
        integer, dimension(:), allocatable :: dumps
        real(kind=8) :: tstep=20d0, &
             toffset=0d0,tcur,dx1,dx2, &
  
negative.code:           poly_alpha=header(10)
            startx1=header(11)
            metric=header(12)
  ! flag for MMKS: metric=0 for MKS, metric=1 for MMKS (=FMKS)
            IS_MMKS=metric
            eHEAT=header(13)
            eRAD=header(14)
            n_prim=header(15)
  