anchor.diff:              stop
            end if
            !Calculate distance (-cos(angle))
  -         d=-dot_product(mesh%v(kc)%p, mesh%edhx(eds(j))%c%p)
  +         if (var%pos==3) then
  +             d=-dot_product(mesh%v(kc)%p, mesh%edhx(eds(j))%c%p)
  +          else if (var%pos==6) then
  +             d=-dot_product(mesh%v(kc)%p, mesh%ed(eds(j))%c%p)
  +         end if
            dmax=max(d, dmax)
            dsum = dsum + 1. - d ** 2
         end do
  
positive.code:       !   matrix (transposed) as the columns of A, and zero out
        !   the lower triangle (upper triangle of A) with Givens
        !   rotations.
        do  i = 1, 6
           ed = eds (i)
           !Apply the rotation for this point
           call aplyr (mesh%edhx(ed)%c%p(1), mesh%edhx(ed)%c%p(2), &
                mesh%edhx(ed)%c%p(3), cx, sx, cy, sy, xp, yp, zp)
  
           !Define the weight of the edge
           wt = 1. / (1. - zp) - rinv
  
negative.code: 
            !Get triangle
            k=gettr(p, mesh)
         case default
            print*, "vecrecon_lsq error: unknown stencil", trim(stencil)
            stop
         end select
      end if
      !print*, k
  
      call vecrecon_lsqfitpol (k, stencil, var, mesh)
  
      vecrecon_lsq=vecrecon_lsqeval(k, p, var)
  