generate 3D image stack from pixellist

Hello, I have a list of objects and their pixel list generated in pixlistC by the code below. I now want to make a 3D image stack with the same dimensions as C containing all pixels from pixlistC == 1 and all other pixels ==0. Can anyone suggest how to do this please? Thank you!

[Clabel1,Cnumobj1]=bwlabeln(C);
pixlistC = regionprops(Clabel1,'PixelList');
nobj = numel(pixlistC);
for ii = nobj:-1:1
        multipageC(ii,1) = numel(unique(pixlistC(ii).PixelList(:,3)))>1;    
end
for jj = nobj:-1:1
    if multipageC(jj)==0
        pixlistC(jj) = [];
        jj=jj-1
    end 
end

Answer 1

[Walter Roberson]

Cs = size(C);
RebuiltC = zeros(Cs);

allids = vertcat(pixlistC.PixelList);
ind = sub2ind(Cs, allids(:,1), allids(:,2), allids(:,3));
RebuiltC(ind) = 1;

Answer 2

[Image Analyst]

Maybe I don't understand what you want. If you take all the pixel list coordinates from all the blobs, and create a new binary image, then you'll end up with just your original C. You don't need to do anything - you've already got it. Why do you think you need to do something?

[Tara]

Thanks for getting back to me. I am not getting back the original C as I've removed objects that were only present in one section of the stack. That's the bottom "for" loop. So pixlistC is now a bit different and I want an image of those objects.

[Image Analyst]

Let's examine, step by step, just exactly what this code

for ii = nobj:-1:1
        multipageC(ii,1) = numel(unique(pixlistC(ii).PixelList(:,3)))>1;    
end

is doing. So you're looping over blobs, from the last blob to the first blob. I don't see any reason to do it in that order, but whatever - it will work. Now pixlistC(ii).PixelList(:,3) is a list of all the z coordinates of every voxel in the 3D blob. So maybe the blob has 10,000 voxels and you'll get 10,000 z values. Now, you pass that into unique() which will give you only the slices (planes) where the ii'th blob resides. for example, maybe that blob has 10,000 voxels in it, but the z values only range from 31 to 39, meaning the blob lives in slices 31-39 inclusive. Now you're taking the number of elements from that with the numel() function. So in this example if you go from 31-39 inclusive, you have 9 slices and numel(unique(pixlistC(ii).PixelList(:,3))) will be 9. But now you're comparing it to 1, so any blob that spans 2 or more slices will have numel()>1 be true. multipageC is a 1D vector, so now you're setting the ii'th element of multipageC to true. At this point we're basically making a list of which blob numbers are 2 or more slices tall in the z direction. Okay, kinda weird and useless, but let's continue.

Now our next loop

for jj = nobj:-1:1
    if multipageC(jj)==0
        pixlistC(jj) = [];
        jj=jj-1
    end 
end

looks at every blob, and if that blob's voxels resided in one and only one plane then multipageC will be false. (Recall that if there were 2 or more slices multipageC was true.) If that happens you're decrementing jj. Now the code has completely lost all credibility. Decrementing jj has no effect. When the loop resumes on the next iteration jj will be exactly what it would have been regardless of whether you decremented it or not in the loop! If you don't believe me - just try it. Even if you were able to decrement it, then you'd be skipping checking the blob after the current one, and that makes no sense whatsoever, but like I said that won't happen.

This code makes no sense at all and is so badly constructed that I can't even begin to figure out what is intended. All I can try to deduce is that you said "I now want to make a 3D image stack" but as I said before, you already have the 3D image stack - it's just the original C.

[Image Analyst]

Ok, after reading your comment, what I think you're trying to do is to remove blobs where the z-height is only 1 slice. To do that you'd do this (untested)

[Clabel1, numObjects] = bwlabeln(C);
measurements = regionprops(Clabel1,'BoundingBox');
allBB = [measurements.BoundingBox];
% Get just the z heights of the blobs
zHeights = allBB(:, 6);
% See which are more than one slice
tallBlobIndexes = find(zHeights > 1); % A logical vector.
% Extract only the tall blobs from Clabel1
C_tall = ismember(Clabel1, tallBlobIndexes);
% Convert from a labeled image into a binary image.
C_tall = C_tall > 0;  % Now it's a logical/binary image.

This is a lot simpler with fewer steps and more efficient than what you were trying to do.