surfbrick(p,t)

function [qu] = surfbrick(p,t)
%SURFbrick: find surface faces from quad 3d mesh
% inputs:
%   p       node coordinates
%   t       node connectivity per each element

% Form all faces, non-duplicates are surface triangles
faces=[t(:,[1,2,6,5]);
       t(:,[1,5,8,4]);
       t(:,[5,6,7,8]);
       t(:,[2,6,7,3]);
       t(:,[1,2,3,4]);
       t(:,[4,3,7,8])];
%node8 = [t(:,8);t(:,7);t(:,6);t(:,5);t(:,4);t(:,3);t(:,2);t(:,1)];
faces2 = sort(faces,2);
[foo,ix,jx] = unique(faces2,'rows');
vec   = histc(jx,1:max(jx));
qx    = find(vec==1);
qu    = faces(ix(qx),:);