How to find the coordinates of all internal edges of a 3D mesh?

Question

I have a 3D image with dimensions specified, and with the code as follows:

% input image
imageSizeY = 30;
imageSizeX = 20;
imageSizeZ = 5;

% obtain size of internal vertices
ny = imageSizeY - 1;
nx = imageSizeX - 1;
nz = imageSizeZ - 1;

% generate 3D mesh (image coordinate system)
[columnsInImage, rowsInImage, pagesInImage] = meshgrid(1:nx, 1:ny, 1:nz);

% obtain coordinates of vertices
vertexCoords = [rowsInImage(:), columnsInImage(:), pagesInImage(:)];

% obtain count of vertices (from top to bottom)
countVertex= ny*(vertexCoords(:,2) - 1) + vertexCoords(:,1) + ny*nx*(vertexCoords(:,3)-1);

In the attached image, the internal vertices are marked dots while the edges are marked x

I have been able to find the coordinates of all the internal vertices (from top to bottom) of the image as shown in the code. I however will need to find the coordinates of all edges but I am not sure I know exactly how to go about achieving this.

Recent figure

Answer 1

Because it's a regular grid, I used logical operations on the coordinates to find their min's and max's:

% deal vertexes into xyz coordinates
[y,x,z] = deal(vertexCoords(:,1),vertexCoords(:,2),vertexCoords(:,3));
nv = numel(x); % number of vertexes
allV = (1:nv).'; % indexes of all vertexes
% compute cube edges indexes => edge = [point1,point2]
edgeIdxs1 = [allV(y < ny),allV(y < ny) + 1];
edgeIdxs2 = [allV(x < nx),allV(x < nx) + ny];
edgeIdxs3 = [allV(z < nz),allV(z < nz) + ny*nx];
edgeIdxs = [edgeIdxs1;edgeIdxs2;edgeIdxs3];
% xyz coordinates of the edges
xEdges = [x(edgeIdxs(:,1)),x(edgeIdxs(:,2))];
yEdges = [y(edgeIdxs(:,1)),y(edgeIdxs(:,2))];
zEdges = [z(edgeIdxs(:,1)),z(edgeIdxs(:,2))];
% "outer" edges (as defined by you) 
outerEdges = all(xEdges == nx,2) | all(xEdges == 1,2) ...
    | all(yEdges == ny,2) | all(yEdges == 1,2) ...
    | all(zEdges == nz,2) | all(zEdges == 1,2);
% "inner" edges indexes
innerEdges = ~outerEdges;
innerEdgeIdxs = edgeIdxs(innerEdges,:);
% xyz coordinates of the inner edges (to plot x on their middle)
xInnerEdges = [x(innerEdgeIdxs(:,1)),x(innerEdgeIdxs(:,2))];
yInnerEdges = [y(innerEdgeIdxs(:,1)),y(innerEdgeIdxs(:,2))];
zInnerEdges = [z(innerEdgeIdxs(:,1)),z(innerEdgeIdxs(:,2))];
% "outer" vertexes indexes
outerVertexes = (x == nx | x == 1) | (y == ny | y == 1) | (z == nz | z == 1);
% "inner" vertexes indexes (to plot o on them)
innerVertexes = ~outerVertexes;
% plot all edges
plot3([x(edgeIdxs(:,1)),x(edgeIdxs(:,2))]',...
    [y(edgeIdxs(:,1)),y(edgeIdxs(:,2))]',...
    [z(edgeIdxs(:,1)),z(edgeIdxs(:,2))]','k');
hold on
% plot x's
plot3(mean(xInnerEdges,2),mean(yInnerEdges,2),mean(zInnerEdges,2),'xr');
% plot o's
plot3(x(innerVertexes),y(innerVertexes),z(innerVertexes),'oc');
axis equal
xlabel('x');ylabel('y');zlabel('z');

this is an example with:

% input image
imageSizeY = 5;
imageSizeX = 5;
imageSizeZ = 5;