Self-answer after topic ping:
Binarization and component analysis are distractions. If every pixel/voxel needs one and only one label, cat+max is the solution. the second output of max() gives the indices of the max values, giving full overage and straightforwardly dealing with any form of overlapping. Basic example use (vol as a cell array storing some number of 3d arrays)
%concatenate in 4th dim, zero vol is so that background is distinct from 1st vol
tmp = cat(4,zeros(size(vol{1}),vol{:}));
[~,atlas] = max(tmp,[],4); %get indices of max values in 4d array
atlas = atlas-1; %makes background 0 rather than 1