function NCP1BasFxnVal = NC_TriP1_BasFxnVal(TriMesh,pt)
%% Nonconforming FEs: Triangular mesh: Nonconforming P1 bas.fxns.: Values 
% Assume pt has size NumEms-by-2 (one point inside each triangular element) 
% James Liu, ColoState; 2012/07--2017/02 

%% Mesh info 
NumEms = TriMesh.NumEms;
k1 = TriMesh.elem(:,1);  k2 = TriMesh.elem(:,2);  k3 = TriMesh.elem(:,3);
x1 = TriMesh.node(k1,1);  y1 = TriMesh.node(k1,2); 
x2 = TriMesh.node(k2,1);  y2 = TriMesh.node(k2,2);
x3 = TriMesh.node(k3,1);  y3 = TriMesh.node(k3,2);
xm1 = 0.5*(x2+x3);  ym1 = 0.5*(y2+y3);
xm2 = 0.5*(x3+x1);  ym2 = 0.5*(y3+y1);
xm3 = 0.5*(x1+x2);  ym3 = 0.5*(y1+y2);
marea = 0.5*((xm2-xm1).*(ym3-ym1)-(xm3-xm1).*(ym2-ym1));

%% 
NCP1BasFxnVal = zeros(NumEms,3);
x = pt(:,1);  y = pt(:,2);
NCP1BasFxnVal(:,1) = ( (xm2.*ym3-xm3.*ym2)-(ym3-ym2).*x+(xm3-xm2).*y)./(2*marea);
NCP1BasFxnVal(:,2) = (-(xm1.*ym3-xm3.*ym1)+(ym3-ym1).*x-(xm3-xm1).*y)./(2*marea);
NCP1BasFxnVal(:,3) = ( (xm1.*ym2-xm2.*ym1)-(ym2-ym1).*x+(xm2-xm1).*y)./(2*marea);

return;