下のコードのように考えてみました。
モルフォロジー処理を利用していますが、切削形状によっては上手くいかないこともあるかもしれません。
%% ダミーデータ
test = zeros(100,100,100);
test(10:90,10:90,10:90) = 1;
volshow(test)
%% 表面ピクセルの取得
test_erode = imerode(test,strel('cube',3)); % モルフォロジー処理
test_surf = test & ~test_erode; % 元のボリュームからの差分が表面に相当
volshow(test_surf);
%% 復元
test_reconstruct = imfill(test_surf,'holes');
val = test == test_reconstruct;
any(~val(:)) % 復元を確認
%% 削る
test2 = test;
test2(80:90,80:90,80:90) = 0;
volshow(test2);
%% 削った部分
part = test & ~test2;
volshow(part)
%% 削った部分に接する元のボリュームの面
surf_on_part = test2 & imdilate(part,strel('cube',3));
test_label = test2;
test_label(surf_on_part) = 2;
labelvolshow(test_label,test2);
