anchor.diff:               !Save index of edge
                eds(j)=mesh%tr(k)%ed(l)
                !Calculate distance (-cos(angle))
  -             d=-dot_product(pref, mesh%edhx(eds(j))%c%p)
  +             if (var%pos==3) then
  +                d=-dot_product(pref, mesh%edhx(eds(j))%c%p)
  +             else
  +                d=-dot_product(pref, mesh%ed(eds(j))%c%p)
  +             endif
                dmax=max(d, dmax)
                dsum = dsum + 1. - d ** 2
             end do
  
positive.code:           !Save index of edge
            eds(i)=mesh%v(kc)%ed(i)
            !Calculate distance (-cos(angle))
            d=-dot_product(mesh%v(kc)%p, mesh%edhx(eds(i))%c%p)
            dmax=max(d, dmax)
            dsum = dsum + 1. - d ** 2
         end do
  
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)
  