diff --git a/src/rashmi's code/CellSegmentation.m b/src/rashmi's code/CellSegmentation.m new file mode 100644 index 0000000..9e51dda --- /dev/null +++ b/src/rashmi's code/CellSegmentation.m @@ -0,0 +1,181 @@ +clc ; +close all; +clear all; +%mkdir Live; +mkdir slice; + + +%% Code to write images + + + +%open(NewObj1); +%open(NewObj2); +updating_array=0; +dead_cells=0; +total_dead_cells=0; +total_live_cells=0; + +live_cells=0; +for videos=91:91 + + +z=0; + +VidObj = VideoReader(strcat('D:\All_videos\Yeast_1\yeast00',num2str(videos),'.avi')); + + + +VidObjDat = read(VidObj); +Nframes = VidObj.NumberOfFrames; +vidHeight = VidObj.Height; +vidWidth = VidObj.Width; +mean_count=1; +pixel_count=1; +%% Generate Background +a = read(VidObj,1); +%%% average ten frames to generate background. +bg = 0; +count = 30 ;% set number of frame to be averaged +N = 1; % start frame number +for i = N:N+count + bg = bg +double(rgb2gray(read(VidObj,i))); +end +bg = uint8(bg /count); % Final Background generated. +%imshow(bg); +% set appropriate cropping for the background +% RS = 135;% start value of the row + %RE = 280; % end value of the row + RS = 76;% start value of the row + RE = 280;% start value of the row + +% RS = 493;% start value of the row +% RE = 638; % end value of the row +CS = 1;% start value of the column +CE = 216; % end value of the column +% crop as per the following imshow(CS:CE,RS:RE) +%imshow(bg(:,RS:RE)); + + + +%% subtract all subsequent frames and write into a video. + ccount = 0; + gccount = 0; + framecount = 0; + + + ccount = 0; + MAVector = 0; + MVector = 0; + framenum = 0; + pframenum = 0; + + for k =1:Nframes + Copy=rgb2gray(read(VidObj,k)); + CurrFrame = double(rgb2gray(read(VidObj,k))); + CurrFrame1=read(VidObj,k); + CurrFrame1=CurrFrame1(:,RS:RE); + Sub = double(CurrFrame(:,RS:RE)-double(bg(:,RS:RE))); + Submin = min(Sub(:)); + Submax = max(Sub(:)); + AdjBGSub = uint8( (Sub - Submin)/(Submax-Submin) * 255); + I = AdjBGSub; + [A, B]=size(I); + +% a = uint8(min(Sub(:))*-1+Sub); +[~, threshold] = edge(I, 'sobel'); +fudgeFactor = 1; +BWs = edge(I,'sobel', threshold * fudgeFactor); +se90 = strel('line', 3, 90); +se0 = strel('line', 3, 0); +BWsdil = imdilate(BWs, [se90 se0]); +BWdfill = imfill(BWsdil, 'holes'); +BWdfillopen = imopen(BWdfill,strel('disk',5,4)); +BWnobord = imclearborder(BWdfillopen, 4); +seD = strel('diamond',2); +BWfinal = imerode(BWnobord,seD); +BWfinal = imerode(BWfinal,seD); +CC =bwconncomp(BWfinal); +Areas = regionprops(CC,'Area'); +MALength = regionprops(CC,'MajorAxisLength'); +MinLength = regionprops(CC,'MinorAxisLength'); +Centroids=regionprops(CC,'Centroid'); + +BWthresh = imcomplement(im2bw(I,0.15)); + % imshow(BWthresh); + + % figure; + % imshow(uint8(BWthresh).*I); + % CC_par =bwconncomp(BWthresh); + % Centroids=regionprops(CC_par,'Centroid'); + % Areas_par = regionprops(CC_par,'Area'); + +for il = 1:CC.NumObjects + + if(Areas(il).Area > 40 && Areas(il).Area < 45840 ) + ccount = ccount +1; + gccount = gccount +1; + + CR1=round(Centroids(il).Centroid(1))-15; + CC1=round(Centroids(il).Centroid(2))-15; + CR2=round(Centroids(il).Centroid(1))+15; + CC2=round(Centroids(il).Centroid(2))+15; + if CR1 <= 0 + CR1=1; + end + if CC1 <= 0 + CC1=1; + end + if CR2 > B + CR2=B; + end + if CC2 > A + CC2=A; + end + + sliceOfImage = CurrFrame1(CC1:CC2,CR1:CR2); + + %imwrite(sliceOfImage,strcat('C:\Users\RBCCPS\Documents\MATLAB\slice/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + %mean(mean_count)=mean2(I(CC1:CC2,CR1:CR2)); + %mean_count=mean_count+1; + %mean_slice=mean2(sliceOfImage); + [p, q]=size(sliceOfImage); + + % BWthresh = im2bw(sliceOfImage,0.15); + % imshow(BWthresh); + + % figure; + % imshow(uint8(BWthresh).*I); + %CC1 =bwconncomp(imcomplement(BWthresh)); + % Areas_stained = regionprops(CC1,'Area'); + % C_Percentstained = sum(struct2array(Areas_stained))/max(struct2array(Areas))*100; +% Percentstained = horzcat(Percentstained, C_Percentstained ); + + z=z+1; %if (pixel_count ~=0) + pixel_count_array(z)=pixel_count; + + + + imwrite(sliceOfImage,strcat('C:\Users\RBCCPS\Documents\MATLAB\slice/FrameNumber','-',num2str(videos),'-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + + + end + + end + end + +updating_array=vertcat(updating_array(:),pixel_count_array(:)); +pixel_count_array=0; +z=0; +videos +end + %end +%end + + + + %end + +%close(NewObj1); +%close(NewObj2);-- \ No newline at end of file diff --git a/src/rashmi's code/CellSegmentationRefined.m b/src/rashmi's code/CellSegmentationRefined.m new file mode 100644 index 0000000..daf47cd --- /dev/null +++ b/src/rashmi's code/CellSegmentationRefined.m @@ -0,0 +1,302 @@ +%% +%The code creates a feature table for the different types of cells and +%is used to create a initial training dataset for the classifier + +%% +clc ; +close all; +clear all; +initialFrame=1646; +finalFrame=1646; +folder='D:\videos\May10\n3\'; +mkdir(strcat(folder,'RBCs')); +mkdir(strcat(folder,'Parasites')); +mkdir(strcat(folder,'Clusters')); +mkdir(strcat(folder,'Mask')); +mkdir(strcat(folder,'choppedCells')); +mkdir(strcat(folder,'Platelets')); +mkdir(strcat(folder,'FalseNegatives')); +mkdir(strcat(folder,'Falsepositive')); +%% +RBCCount=0; parCount=0; clusterCount=0; Fcount=0; choppedCells=0;platelets=0; +Percentstained = 0; C_Percentstained = 0; +count = 0; mainA = 0; StainA = 0;circularity=[];circularity_par=[];circularity_rbc=[]; +P_axisRatio=[]; P_majAxis=[]; P_minAxis=[]; P_area=[]; +C_axisRatio=[]; C_majAxis=[]; C_minAxis=[]; C_area=[]; +stainedpercent = []; +area_par=[]; +MajAxis_par=[]; MinAxis_par=[]; AxisRatio_par=[]; +%% +% To create a circular mask to eliminate the spokes channels and perform +% processing of cells in the central imaging region +% Create a logical image of a circle with specified +% diameter, center, and image size. +fcount=0; + +% First create the image. +imageSizeX = 640; %X pixel size of the image +imageSizeY = 480; %Y pixel size of the image +[columnsInImage, rowsInImage] = meshgrid(1:imageSizeX, 1:imageSizeY); +% Next create the circle in the image. +centerX = 330; %center X pixel of the circle +centerY = 260; %center y pixel of the circle +radius = 160; %Radius of the circle +circlePixels = (rowsInImage - centerY).^2 ... + + (columnsInImage - centerX).^2 <= radius.^2; % creates the circular mask + +%% Generate Background +%a = rgb2gray(imread(strcat('E:\Experiments\odroid\step1\1 (',num2str(V),').jpg'))); +%imshow(a); +%%% average ten frames to generate background. +bg = 0; +count = 300 ;% set number of frame to be averaged +N = 2; % start frame number +for i = N:N+count + % bg = bg +double(rgb2gray((imread(strcat('D:\videos\S\S1\1 (',num2str(V+i),').jpg'))))); + + imagename=char(sprintf('%08d',i)); + bg = bg +double(rgb2gray(imread(strcat(folder,imagename,'.jpg')))); +end +bg = uint8(bg /count); % Final Background generated. +imshow(uint8(circlePixels).*bg); +%% +tic +r =0;b=0; +count = 0; +maxglcm = 0; +Cont =0;Corr=0;Homo=0;I_mean=0;I_std=0;I_cir=0;I_area =0;stain=0; + ccount = 0; + C_Percentstained =0;StainMaxlength=0;StainMinlength=0;StainSolid=0; StainNumobj=0; + Features_table = [];F=[]; +CellFeature_table = [];cell=[]; +% Feature_table = table('VariableNames',{'Area', 'ConvexArea', 'Eccentricity', 'EquivDiameter', 'EulerNumber', 'Extent', 'FilledArea', 'MajorAxisLength', 'MinorAxisLength', 'Orientation', 'Perimeter', 'Solidity'}); + for k =initialFrame:finalFrame + k +Shape_Features =[]; Texture_Features1 = []; +% The if condition is used to refresh the background for every 2000 frames +% This helps eliminate any debris/ struck cells in the ROI +if mod(k,2000) == 0 + %if V <19999-30 +% bg1 = bg; +bg =0; +count = 300 ;% set number of frame to be averaged +N = k; % start frame number +for i = N:N+count + % bg = bg +double(rgb2gray((imread(strcat('D:\videos\S\S1\1 (',num2str(V+i),').jpg'))))); + + imagename=char(sprintf('%08d',i)); + bg = bg +double(rgb2gray(imread(strcat(folder,imagename,'.jpg')))); +end +bg = uint8(bg /count); % Final Background generated. +% bg = (bg+bg1)/2.0 +end + % The following section performs the segmentation based on histogram values +imagename=char(sprintf('%08d',k)); + + CurrFrame = (rgb2gray((imread(strcat(folder,imagename,'.jpg'))))); + Sub =double(double(CurrFrame.*uint8(circlePixels))-double(bg.*uint8(circlePixels))); + %imshow(Sub); + Submin = min(Sub(:)); + Submax = max(Sub(:)); + AdjBGSub = uint8( (Sub - Submin)/(Submax-Submin) * 255); + Ia = AdjBGSub; + [A, B]=size(Ia); + I = medfilt2(Ia); % To smoothen the image + I = adapthisteq(I); +I1 = medfilt2(I); +Irgb = cat(3, I1, I1, I1); +meanIntensityValue(k) = mean2(I1); % Finds the mean of the intensities of the image pixels +stdIntensityValue = std2(I1); % Finds the standard deviation of the intensities of the image pixels +% Selection of max and minimum of intensities for the thresholding +%This multipication factor can be varied when you are optimizing the +%thresholding +Imax = meanIntensityValue(k)+stdIntensityValue*4; % mean+ 4*standard deviation +Imin = meanIntensityValue(k)-stdIntensityValue*4; % mean- 4*standard deviation +% Further processing of the thresholded image +Mask = createMask(Irgb,Imax,Imin);%figure;imshow(Mask); +MaskInv = ~Mask; +Maskopen =bwareaopen(MaskInv,150);%figure;imshow(Maskfinal); +% Maskdil = imdilate(Maskopen, [se90 se0]);%figure;imshow(Maskdil); +Maskdil = imclose(Maskopen, strel('disk',5)); +Maskfill = imfill(Maskdil, 'holes');%figure;imshow(Maskfill); +Maskclose = Maskfill; +% Maskclose = imclose(Maskfill, strel('disk',5)); +Maskclear = imclearborder(Maskclose);%figure;imshow(Maskclear); +Maskfinal = Maskclear; +% % Filter image based on image properties. +% Maskfinal = bwpropfilt(Maskfinal, 'Area', [350 + eps(350), Inf]); % Area greater than 350 pixels +% Maskfinal = bwpropfilt(Maskfinal, 'Solidity', [0.6 + eps(0.6), Inf]); % 1 is completely solid region +% Maskfinal = bwpropfilt(Maskfinal, 'EulerNumber', [4.94065646e-324 + eps(4.94065646e-324), Inf]);k + CC =bwconncomp(Maskfinal); +AxisRatio=[]; + +Centroids = regionprops(CC,'Centroid'); +Areas = regionprops(CC,'Area'); + +Perimeters=regionprops(CC,'Perimeter'); +%mainA = horzcat(mainA, max(struct2array(Areas))); + for il = 1:CC.NumObjects + if (CC.NumObjects>0) + if (Areas(il).Area > 500 && Areas(il).Area < 5000 ) + ccount=ccount+1; + X = Centroids(il).Centroid(1); + Y = Centroids(il).Centroid(2); + lowX = round(X)-20; + lowY = round(Y)-20; + HighX = round(X)+20; + HighY = round(Y)+20; + if lowX <= 0 + lowX=1; + end + if lowY <= 0 + lowY=1; + end + if HighX > B + HighX=B; + end + if HighY > A + HighY=A; + end + I1=Ia(lowY:HighY,lowX:HighX); + circularity(ccount) = Perimeters(il).Perimeter ^ 2 / (4 * pi * Areas(il).Area); + + BWthresh = im2bw(I1,0.2); + %imtool(BWthresh); +% Par_dil = imdilate(BWthresh, [se90 se0]); +% Par_fill = imfill(BWthresh, 'holes'); + imwrite(Maskfinal,strcat(folder,'Mask/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + %imtool(BWthresh) +% imshow(BWthresh); + +% figure; +% imshow(uint8(BWthresh).*Ia); + + CC1 =bwconncomp(imcomplement(BWthresh)); + Areas_stained = regionprops(CC1,'Area'); + MajAxis= regionprops(CC1,'MajorAxisLength'); + MinAxis= regionprops(CC1,'MinorAxisLength'); + AxisRatio=struct2array(MajAxis)/struct2array(MinAxis); + struct2array(MajAxis); + struct2array(MinAxis); + Areas_stained.Area; + + C_Percentstained = sum(struct2array(Areas_stained))/max(struct2array(Areas))*100; + Percentstained = horzcat(Percentstained, C_Percentstained ); + StainA = horzcat(StainA, sum(struct2array(Areas_stained))); + MajAxisMat=struct2array(MajAxis); + numMajaxis=size(MajAxisMat,2); + no_majAxis=isempty(MajAxis); + MajAxis_par=[];MinAxis_par=[]; + MajAxis_par=(struct2array(MajAxis)); + MinAxis_par=(struct2array(MinAxis)); + Area_par=struct2array(Areas_stained); + + + i=1; + f=numMajaxis; + % For loop to compute ratio of major and minor axis + % to check if the object in the BWthresh is long or + % round shaped + AxisRatio_par=[]; + + for i=1:f + AxisRatio_par(i)=MajAxis_par(i)/MinAxis_par(i) + end + AxisRatio_par + countgreater2=0; + %Another for loop to check if the object has the + %ratio greater than 2.7, it should be classified as + %RBC and not as a parasite. + + for i=1:f + if (AxisRatio_par(i) > 2.7) + countgreater2=countgreater2+1; + end + end + + + if (Areas(il).Area > 1500 && Areas(il).Area < 5000 ) + clusterCount=clusterCount+2; + imwrite(I1,strcat(folder,'Clusters/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + %If there are major Axis + elseif (no_majAxis == 0) + if (C_Percentstained<40 && C_Percentstained>3 && AxisRatio < 3 && max(struct2array(MajAxis)) < 18 && numMajaxis < 4 ) + stainedpercent(ccount) = C_Percentstained; + area_par(ccount)=Areas(il).Area; + + + + %%If there are no cells having count greater than + %%2.7 i.e countgreater2==0 then the cell is + %%classified as a parasite + if (countgreater2 == 0) + 'Parasite1' + parCount=parCount+1; + imwrite(I1,strcat(folder,'Parasites/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'Falsepositive/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + +% 'default' +% if (max(Area_par) > 20 && max(AxisRatio_par) <4 && max(struct2array(MajAxis)) < 25) +% 'Parasite3' +% parCount=parCount+1; +% imwrite(I1,strcat(folder,'Parasites/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); +% imwrite(BWthresh,strcat(folder,'Falsepositive/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + elseif (max(struct2array(MajAxis))< 10 && C_Percentstained > 5 && max(AxisRatio_par) < 3.5) + + 'Parasite2' + parCount=parCount+1; + imwrite(I1,strcat(folder,'Parasites/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'Falsepositive/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + else + 'RBC1' + RBCCount=RBCCount+1; + circularity_rbc(ccount)=Perimeters(il).Perimeter ^ 2 / (4 * pi * Areas(il).Area); + imwrite(I1,strcat(folder,'RBCs/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'FalseNegatives/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + + end + %RBC + elseif ( max(struct2array(MajAxis)) > 15 || numMajaxis>1 || (C_Percentstained<5 && C_Percentstained>1) || AxisRatio > 2) + + if (max(struct2array(MajAxis))< 10 && C_Percentstained > 10 && max(AxisRatio_par) < 2 && max(Area_par) < 200) + + 'Parasite3' + parCount=parCount+1; + imwrite(I1,strcat(folder,'Parasites/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'Falsepositive/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + else + 'RBC2' + RBCCount=RBCCount+1; + circularity_rbc(ccount)=Perimeters(il).Perimeter ^ 2 / (4 * pi * Areas(il).Area); + imwrite(I1,strcat(folder,'RBCs/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'FalseNegatives/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + + end + end + + elseif (~isempty(Areas_stained)) + if(Areas_stained.Area >25) + parCount=parCount+1; + imwrite(I1,strcat(folder,'Parasites/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'Falsepositive/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + end + else + 'RBC3' + RBCCount=RBCCount+1; + imwrite(I1,strcat(folder,'RBCs/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + imwrite(BWthresh,strcat(folder,'FalseNegatives/FrameNumber','-',num2str(k),'-',num2str(il),'.jpg'),'jpeg'); + + end + end + end + end + end + Fcount=clusterCount+RBCCount+parCount+choppedCells; + +