Hello.. I am working on detecting lung nodules and I would like to know how to convert the lung nodule extracted from a CT image into a 3D image.

clear
clc
I=imread('lung1.jpg');
%I= imresize(J,[512 512]);
a=graythresh(I);
j2=im2bw(I,a);
figure;imshow(I)
figure;imshow(j2)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
c=bwconncomp(j2);
stats=regionprops(c,'basic');
AA=[stats.Area];
[~,big]=max(AA);
j2(labelmatrix(c)~=big)=0;
figure;imshow(j2)
j2=~j2;
figure;imshow(j2)
%%%%%%%%%%%%%%%%%%
b1=strel('disk',8);
aa=imclose(j2,b1);
figure;imshow(aa)
%%%%%%%%%%%%%%%%%%%%%%%%
aa=~aa;
hh=imfill(aa,'holes');
figure;imshow(hh)
q=~aa;
q(~hh)=0;
figure;imshow(q)
%%%%%%%%%%%%%%%%%%%%%%%%
b2=strel('disk',10);
c2=imopen(q,b2);
c3=imclose(c2,b2);
figure;imshow(c3)
%y=rgb2gray(I);
y=I
y(~c3)=0;
figure;imshow(y)
%%%%%%%%%%%%%%%%%%%%%%%
figure;imhist(y)
axis([0 250 0 2500]);
%z=rgb2gray(I);
%figure;imhist(z)
j1=im2bw(y,150/255);
figure;imshow(j1)
%%%%%%%%%%%%%%%%%%
cc1=bwconncomp(j1);
l=labelmatrix(cc1);
[l2,num]=bwlabel(j1);
stats11=regionprops(cc1,'Area');
Z=[stats11.Area];
for a=1:num
if Z(a)<15
j1(l==a)=0;
end
end
cc2=bwconncomp(j1);
l=labelmatrix(cc2);
[l2,num2]=bwlabel(j1);
figure;imshow(l2)
%%%%%%%%%%%%%%%%%%%%%%%%%%%
j2=j1;
cc=bwconncomp(j2);
l=labelmatrix(cc);
[l2,num]=bwlabel(j2);
stats=regionprops(cc,'Area','Perimeter','MajorAxisLength','MinorAxisLength','Solidity','Extent');
A1=[stats.Area];
A2=[stats.Perimeter];
A3=[stats.MajorAxisLength];
A4=[stats.MinorAxisLength];
s3=[stats.Solidity];
s4=[stats.Extent];
s1=zeros(1,num);
s2=zeros(1,num);
s5=zeros(1,num);
s6=zeros(1,num);
for a=1:num
s1(1,a)=4*pi*A1(a)/A2(a)^2;
s2(1,a)=4*A1(a)/pi*A3(a);
s5(1,a)=s2(1,a)^0.5;
s6(1,a)=A3(a)/A4(a);
end
x=zeros(6,num);
x(1,:)=s1(1,:);
x(2,:)=s2(1,:);
x(3,:)=s3(1,:);
x(4,:)=s4(1,:);
x(5,:)=s5(1,:);
x(6,:)=s6(1,:);
for a=1:num;
if x(1,a)<0.44||x(1,a)>1.251||x(3,a)<0.68||x(4,a)<0.5||x(4,a)>0.86||x(5,a)<13.5||x(6,a)<1||x(6,a)>1.6
j2(l==a)=0;
end
end
cc1=bwconncomp(j2);
l=labelmatrix(cc1);
[l2,num]=bwlabel(j2);
figure;imshow(l2)
%%%%%%%%%%%%%%%%%%%%%5
stats11=regionprops(cc1,'Area');
Z=[stats11.Area];
for a=1:num
    if Z(a)>250
j2(l==a)=0;
end
end
figure;imshow(j2)
%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Lung nodules diagnosis.
cc1=bwconncomp(j2);
l=labelmatrix(cc1);
[l2,num]=bwlabel(j2);
stats22=regionprops(cc1,'MinorAxisLength','Perimeter');
Z1=[stats22.MinorAxisLength];
Z2=[stats22.Perimeter];
T=zeros(1,num);
q1=j2;
q2=j2;
for a=1:num
T(a)=Z2(a)/Z1(a);
if T(a)>3.5
q1(l==a)=0;
end
if T(a)<=3.5
q2(l==a)=0;
end
end
y=strel('disk',2);
q1=imdilate(q1,y);
q2=imdilate(q2,y);
q1=edge(q1);
q2=edge(q2);
y=strel('disk',2);
q1=imdilate(q1,y);
q2=imdilate(q2,y);
[m,n]=size(j2);
z=zeros(m,n,3);
k=zeros(m,n);
l=zeros(m,n);
z(:,:,1)=im2double(j2);
z(:,:,2)=z(:,:,1);
z(:,:,3)=z(:,:,1);
k=z(:,:,1);
l=z(:,:,2);
k(q2)=255;
l(q1)=255;
z(:,:,1)=k;
z(:,:,2)=l;
figure,imshow(z)

[Walter Roberson]

I=imread('lung1.jpg');

In this code, you only have a single slice. You cannot turn a single slice into a 3D matrix -- not unless all you are looking for is to extrude the shape into the third dimension.

[noor]

Yes, I want to convert the shape to 3D.. Can you help me write the rest of the code to display the image in 3D?

[Walter Roberson]

Are you sure that you want to take the regions from the one slice, and create a 3d solid from it?

For example if all you had was the red string on the left, and you wanted to create the ribbon from it??

It seems much more likely to me that you would have additional image files such as lung2.jpg, lung3.jpg, and so on, and that you would want to stack the identified regions so as to create a 3D outline of the nodule based upon different slices of data.

[noor]

Yes, I want to show 3D for each image individually. Is there a way to do that?

Answer 1

[Walter Roberson]

Take the binary image and use bwtraceboundary or bwboundaries to get the coordinates of the edge pixels.

Take each outer edge's row coordinates and call those y. Take each outer edge's column coordinates and call those x.

Construct a list of coordinates, something like this

Zbase = 0;   %where the bottom plate will be in z
Zupper = 2;  %where the top plate will be in z
Z = [Zbase * ones(numel(x), 1); Zupper * ones(numel(x),1)];
XY = [x(:), y(:)];
[uXY, ~, uXYidx] = unique(XY, 'rows');
nuXY = size(uXY,1);

uXYZ_lower = [uXY, Zbase * ones(nuXY), 1)];
uXYZ_upper = [uXY, Zupper * ones(nuXY), 1)];
vertices = [uXYZ_lower; uXYZ_upper];
LL = uXYidx;
LR = circshift(LL,-1);
UL = 2*nuXY - LL + 1;
UR = circshift(UL, 1);
Faces = [LL, LR, UR, UL];

p = patch('Faces', Faces, 'Vertices', vertices);

This is taking the XY coordinates with Z = Zbase, and the same XY coordinates with Z = Zupper. Then, if I have managed to get the boundary conditions right, it draws the edge faces as rectangles.

You will probably also want to draw the top and bottom caps. That will call for two rows in Faces with long vertices list. The columns in Faces that define the quads will need to be nan padded to match. (You might possibly need to add an explicit return to LL when you nan pad.)

The above code is not tested; I might not have constructed the vertex indices for the faces correctly.