%% reads in cell images (file) and parameters (params) and outputs % "blebdata" which is a structure array containing the following fields: % blebdata.bleb(:,:,j,i) = binary matrix for bleb j at time frame i % blebdata.blebsize(:,i) = array of blebsizes % blebdata.alpha(:,i) = array for the "roundness" of cell % blebdata.showAllBlebs(:,:,i) = single binary matrix for all the blebs at time frame i % blebdata.numberofblebs(i) = number of blebs at time frame i % blebdata.numberofimages = total number of images % blebdata.images.intensities(:,:,i) = grayscale image matrix for time frame i % blebdata.images.BW(:,:,i) = black and white image matrix for time frame i % blebdata.images.BW2(:,:,i) = black and white image matrix after filling in the cell and removing isolated pixels % blebdata.indivBlebEdge(:,:,j,i) = matrix of 1s 2s and 3s representing the blebbing intensities of bleb j at time frame i % blebdata.outline(:,:,i) = binary matrix representing the outline of the cell function [ blebdata ] = findBlebs( parameters , filedata ) tic file = filedata.filename; %'agar_2per_5_C1.tif'; x = filedata.numberofimages; % 89; % number of images outputFilename = filedata.outputfilename; %'BLEBDATA_.mat'; outputToFile = filedata.outputtofile; % 0; lowerBlebSizeThreshold = parameters.lowerblebsizethreshold; %100; upperBlebSizeThreshold = parameters.upperblebsizethreshold; %500; alphaThreshold = parameters.alphathreshold; %0.1; % clear unnecessary variables clearvars -except outputToFile file x outputFilename lowerBlebSizeThreshold upperBlebSizeThreshold alphaThreshold a = 0; % first image minIntensity = 25; % intensities of less minIntensity are zeroed out stepSize = 1; %% image input for i = 1 : x I = imread( file ,i+a); % read in image i+a Intensities(:,:,i) = I; I(I < minIntensity ) = 0; I(I > minIntensity - 1 ) = 255; % fill in the cell J = bwareaopen( imfill(I,'holes') , 5000 ).*255 ; % filled and removed "dust" [ blebdata.cellsize(i) p ] = size( find( J ) ); BW(:,:,i) = I; BW2(:,:,i) = J; % get the outline of the cell outline(:,:,i) = imdilate(double(bwperim(BW2(:,:,i),8)).*255, strel('sphere',1)); end %% Looping through to find blebs for i = 1 : x - stepSize % zero of the correct size z = zeros(size(BW(:,:,i))); diff = BW2(:,:,i+stepSize) - BW2(:,:,i); diff(diff < 0) = 0; % set negative pixels to zero diff = imfill(diff,'holes') ; % filled in image % diff = rsg_expand(diff, 1); % expand by a pixel in each direction % find the connected components of the difference matrix [L, numberconnected] = bwlabel( diff , 8 ); numberofblebs(i) = 0; showAllBlebs(:,:,i)= z; count = 0; for j = 1 : numberconnected B = double(L == j); % B is connected component of the difference matrix numberpixels = sum(sum(B)); % number of pixels in the connected component % get indices of non-zero pixels [r c] = find( B ); diam = 0; % calculate the diameter of the connected component for l = 1 : length(r) for k = 1 : length(c) b = sqrt(( r(l) - r(k) ).^2 + ( c(l) - c(k) ).^2); if (b > diam ) diam = b; end end end % calculate the alpha sizetodiam2 = (4 / pi) * numberpixels / ( diam * diam ); if or(numberpixels > upperBlebSizeThreshold - 1 , ... % bleb and( numberpixels > lowerBlebSizeThreshold - 1, sizetodiam2 > alphaThreshold ) ) numberofblebs(i) = numberofblebs(i) + 1; blebsize(numberofblebs(i),i) = numberpixels; alpha(numberofblebs(i),i) = sizetodiam2; bleb(:,:,numberofblebs(i),i) = B; showAllBlebs(:,:,i) = max(showAllBlebs(:,:,i), bleb(:,:,numberofblebs(i),i)); else % not a bleb count = count + 1; end end if (count == numberconnected) % no blebs bleb(:,:,1,i) = z; blebsize(1,i) = 0; alpha(1,i) = 0; end % find the centre of the cell [rowo,colo] = find(outline(:,:,i)); valo = max(size(rowo)); origin(:,i) = [mean(colo) mean(rowo)]; % if no blebs then no bleb edges if numberofblebs(i) == 0 indivEdgeReal(:,:,1,i) = z; else % for each bleb, find the bleb edge and assign the bleb intensities for j = 1 : numberofblebs(i) indivEdgeReal(:,:,j,i) = double(bleb(:,:,j,i)); % get the intersection of the bleb and the cell outline A = double( and(indivEdgeReal(:,:,j,i).*255 == outline(:,:,i), outline(:,:,i) ~= 0) ); indivBlebEdge(:,:,j,i) = double(bwmorph(A,'thin')).*255; % find the endpoints of the bleb - since the bleb edge is not nec. connected, there may be more than 2 of these indivBlebEdgeEndpoints(:,:,j,i) = double(bwmorph(indivBlebEdge(:,:,j,i),'endpoints')); [rrefer,crefer] = find(indivBlebEdgeEndpoints(:,:,j,i)); % get non-zero rows and columns of the endpoints numberofedgeendpoints = max(size(rrefer)); % get the number of endpoints if numberofedgeendpoints < 2 else % calculate the angle between each combination of edge points theta=zeros(numberofedgeendpoints, numberofedgeendpoints); for p= 1 : numberofedgeendpoints pointp = [rrefer(p), crefer(p)]; for q = p + 1 : numberofedgeendpoints pointq = [rrefer(q), crefer(q)]; theta(p,q) = angle_between(pointp,pointq,origin(:,i).'); end end % the actual endpoints are the two points that maximise the angle between two points on bleb edge I_row=[]; I_col=[]; [thetamax,thetaindex] = max(theta(:)); [I_row, I_col] = ind2sub(size(theta),thetaindex); endpoint1(:,j,i) = [rrefer(I_row) crefer(I_row)]; endpoint2(:,j,i) = [rrefer(I_col) crefer(I_col)]; [re,ce] = find(indivBlebEdge(:,:,j,i)); vale = max(size(re)); % for each point in the bleb edge assign a blebbing intensity based on the angle between the point and one of the endpoints for p=1:vale point = [re(p), ce(p)]; if or(angle_between(endpoint1(:,j,i).',point,origin(:,i).') < thetamax/6 , angle_between(endpoint1(:,j,i).',point,origin(:,i).') >= 5*thetamax/6) indivBlebEdge(re(p),ce(p),j,i) = 1; elseif or(angle_between(endpoint1(:,j,i).',point,origin(:,i).') < thetamax/3 , angle_between(endpoint1(:,j,i).',point,origin(:,i).') >= 2*thetamax/3) indivBlebEdge(re(p),ce(p),j,i) = 2; else indivBlebEdge(re(p),ce(p),j,i) = 3; end end end end end end toc parameters.inputFilename = file; parameters.numberOfImages = x; parameters.minimumIntensityThreshold = minIntensity; parameters.lowerBlebSizeThreshold = lowerBlebSizeThreshold; parameters.upperBlebSizeThreshold = upperBlebSizeThreshold; parameters.alphaThreshold = alphaThreshold; %clearvars -except outputToFile outputFilename parameters bleb blebsize alpha showAllBlebs numberofblebs blebdata.bleb = bleb; blebdata.blebsize = blebsize; blebdata.alpha = alpha; blebdata.showAllBlebs = showAllBlebs; blebdata.numberofblebs = numberofblebs; blebdata.numberofimages = x; blebdata.images.intensities = Intensities; blebdata.images.BW = BW; blebdata.images.BW2 = BW2; blebdata.indivBlebEdge = indivBlebEdge; blebdata.outline = outline; if (outputToFile == 1) clearvars -except blebdata parameters outputFilename save( outputFilename ); end clearvars -except blebdata end function angle = angle_between(a,b,o) angle = acos((dot(a-o,b-o)) / (norm(a-o)*norm(b-o))); end