Commit e509701b authored by Antonio Politi's avatar Antonio Politi
Browse files

compatible new bioformat

parent 8b40299a
function [bgVal, imgfinal] = bgCompLiveCell(imgfile, h5segdir, segdir, outdir, upfl, iC)
%%
% bgCompLiveCell(lsmfile, segdir, outdir)
% - imgfile: full path to lsm file or h5 file
......@@ -14,22 +15,22 @@
% Antonio Politi EMBL September 2016, December 2016 (changed to read hdf5
% too)
%%
function [bgVal, imgfinal] = bgCompLiveCell(imgfile, h5segdir, segdir, outdir, upfl, iC)
global h5
%% input data and directories
if nargin < 1
% imgfile can be a tif, lsm or h5
imgfile = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160701_gfpNUP107z26z31\Result\BGComp\h5\160701-HeLa-Nup107-cell-0-t6.h5';
imgfile = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160706_gfpNUP107z26z31\Result\BGComp\h5\160706-HeLa107-4D-cell-5-t6.h5';
%imgfile = 'W:\Shotaro\Live-cell-imaging\HeLa-FCS\160701-HeLa-Nup107\160701-HeLa-Nup107-cell-0-t6.tif';
end
if nargin < 2
h5segdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160701_gfpNUP107z26z31\Result\BGComp\h5_SimpleSegmentation';
h5segdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160706_gfpNUP107z26z31\Result\BGComp\h5_SimpleSegmentation';
end
if nargin < 3
segdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160701_gfpNUP107z26z31\Result\Segmentation\';
segdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160706_gfpNUP107z26z31\Result\Segmentation\';
end
if nargin < 4
outdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160701_gfpNUP107z26z31\Result\BGComp\';
outdir = 'C:\Users\toni\Dropbox\NPCMaturation\expdata\livecellCalibrated\160706_gfpNUP107z26z31\Result\BGComp\';
end
if nargin < 5
upfl = 600;
......@@ -77,10 +78,11 @@ erodFil = zeros(3,3,3);
erodFil(:,:,1) = [0 0 0; 0 1 0; 0 0 0];
erodFil(:,:,2) = [0 1 0; 1 1 1; 0 1 0];
erodFil(:,:,3) = [0 0 0; 0 1 0; 0 0 0];
maxErosion = 5;
maxErosion = 3;
for iT = 1:dim.Nt
try
% the scripts adds a digit code for time point. Default is 001
load(fullfile(segdir, [name sprintf('%03d',iT) '.mat']));
catch
display([name ' has not been segmented yet ']);
......@@ -94,16 +96,26 @@ for iT = 1:dim.Nt
% detect glass surface. This is a slice that has everywhere a higher
% intensity in background
for iZ = 1:dim.Nz
gfpC = imfilter(getSlice(reader, iZ, iT, iC), mF);
gfpC = imfilter(getSlice(reader, iZ, iT, iC, dim), mF);
sV = sort(gfpC(:));
sumGfp(iZ) = mean(sV(1:20));
end
[val, iZglass] = max(sumGfp);
%% for each nuclei
for inuc = 1:numNuc
% check if indeed 2 nuclei have been segmented
levels = unique(nucMaskNonIso);
nucVec = 1:numNuc;
if length(levels) < numNuc +1
display('nuclei have not been segmented properly!')
nucVec = max(levels);
end
for inuc = nucVec
%% nucleus mask (J. Hossain segmentation)
binNuc = nucMaskNonIso==inuc;
%% find slice with largest perimeter avoiding smallest segmented parts
for iZ=1:dim.Nz
props = regionprops(binNuc(:,:,iZ));
......@@ -121,26 +133,38 @@ for iT = 1:dim.Nt
iZ = round(p(1).Centroid(3));
end
%% avoid being too close to bottom of slide
iZn = iZ*(iZ>=iZglass + 2) + (iZglass+2)*(iZ<iZglass + 2);
iZc = (iZ-2)*(iZ-2>=iZglass + 2) + (iZglass+2)*(iZ-2<iZglass + 2);
%iZn = iZ*(iZ >= (iZglass + 2)) + (iZglass+2)*(iZ< (iZglass + 2));
iZn = iZ;
% cytoplasmic slice is twice below nucleus slice
%iZc = (iZ-2)*(iZ-2>=iZglass + 2) + (iZglass+2)*(iZ-2<iZglass + 2);
iZc = iZ-2;
% find slice that contains center of mass (this is more prone to be sometime close to lower or upper part of nucleus
% props = regionprops(binNuc);
% iZ = round(props.Centroid(3))+1;
%% find intensity inside nucleus. Skeletonze nucleus and dilate skeleton
%% find intensity inside 2D nucleus. Skeletonze nucleus and dilate skeleton
% use segmented nuclei instead
% get segmented slice. This is a labelled image
gfpC = getSlice(reader, iZn, iT, iC);
gfpC = getSlice(reader, iZn, iT, iC, dim);
%avoid very high fluorescent spots
gfpCM = (imgaussian(gfpC,2) <= upfl);
nucIntM = skeleton((reshape(ilast(1,:,:,iZn,iT),dim.Ny, dim.Nx) ==2).*binNuc(:,:,iZn))> minSkel; %avoids problems in anaphase
%nucIntM = skeleton((reshape(ilast(1,:,:,iZn,iT),dim.Ny, dim.Nx) ==2).*binNuc(:,:,iZn))> minSkel; %avoids problems in anaphase
nucIntM = skeleton(imerode(binNuc(:,:,iZn), erodFil).*gfpCM)> minSkel; %minSkel; %avoids problems in anaphase
if sum(nucIntM(:)) == 0
nucIntM = skeleton(binNuc(:,:,iZn).*gfpCM)> 1;
end
%nucIntM = skeleton(imfill(binNuc(:,:,iZn), 'holes'))>minSkel;
% dilate once the mask of skeleton
nucIntM = (gfpCM).*imdilate(nucIntM, strel('disk',1));
nucInt = double(nucIntM).*double(gfpC);
bgVal(iT, 1+(inuc-1)*3) = sum(nucInt(:))/sum(nucIntM(:));
%% find intensity of eroded 3D nuclei
nucIntMEr = (reshape(ilast(1,:,:,:,iT),dim.Ny, dim.Nx, dim.Nz) == 2).*binNuc;
%% find intensity of eroded 2D nuclei
gfpC = getSlice(reader, iZn, iT, iC, dim);
gfpCM = (imgaussian(gfpC,2) <= upfl);
nucIntMEr = binNuc(:,:,iZn).*gfpCM;
%nucIntMEr = (reshape(ilast(1,:,:,:,iT),dim.Ny, dim.Nx, dim.Nz) == 2).*binNuc;
for iErode = 1:maxErosion
nucIntMLoc = imerode(nucIntMEr, erodFil);
%avoid erode up to 0
......@@ -150,20 +174,19 @@ for iT = 1:dim.Nt
break
end
end
gfpC = getSlice(reader, [1:dim.Nz], iT, iC);
gfpCM = (imgaussian(gfpC,2) <= upfl);
nucIntMEr = gfpCM.*nucIntMEr;
nucIntEr = double(nucIntMEr).*double(gfpC);
bgVal(iT, 2+(inuc-1)*3) = sum(nucIntEr(:))/sum(nucIntMEr(:));
%% find intensity of cytoplasm
segSlice = reshape(ilast(1,:,:,iZc,iT),dim.Ny, dim.Nx);
gfpC = getSlice(reader, iZn, iT, iC);
gfpC = getSlice(reader, iZc, iT, iC, dim);
gfpCM = (imgaussian(gfpC,2) <= upfl);
% use both nucleus and cytoplasmic mask
bwerode = imerode((segSlice==2) + (segSlice==3),strel('disk',10));
bwerode = bwlabel(bwerode);
bwerode = imdilate(bwerode, strel('disk',9)); %this should separate touching objects
% perform addiotinal marked seed watershed
bwerode = imfill(bwerode);
% perform additional marked seed watershed
D = bwdist(bwerode);
% boundary should also be a minimum
D(:,[1 end]) = -100;
......@@ -179,12 +202,16 @@ for iT = 1:dim.Nt
cytIntM = (gfpCM).*(labimg==lab).*~(imdilate(imfill(binNuc(:,:,iZ), 'holes'), strel('disk',4)));
cytInt = double(cytIntM).*double(gfpC);
bgVal(iT, 3+(inuc-1)*3) = sum(cytInt(:))/sum(cytIntM(:));
if numNuc == 1
bgVal(iT, [4 5 6]) = bgVal(iT, [1 2 3]);
end
% create output image
imgfinal(:,:,iT) = imgfinal(:,:,iT) + binNuc(:,:,iZ)*(1+(inuc-1)*3) + nucIntM*(2+(inuc-1)*3) + (labimg==lab)*(3+(inuc-1)*3) ;
end
if bgVal(iT, 4) == 0
bgVal(iT,[4 5 6]) = bgVal(iT, [1 2 3]);
end
if bgVal(iT, 1) == 0
bgVal(iT, [1 2 3]) = bgVal(iT,[4 5 6]);
end
end
%%
......@@ -210,13 +237,19 @@ end
% iZ : index of Z slice
% iT : index of time point
% iC : index of channel
function slice = getSlice(reader, iZ, iT, iC)
function slice = getSlice(reader, Z, iT, iC, dim)
global h5
if h5
slice = reader(:,:,iZ, iT, iC);
slice = reader(:,:,Z, iT, iC);
else
iPlane = reader.getIndex(iZ-1,iC-1,iT-1) + 1;
slice = bfGetPlane(reader,iPlane);
slice = zeros(dim.Ny, dim.Nx, length(Z));
for iZ = 1:length(Z)
iPlane = reader.getIndex(Z(iZ)-1,iC-1,iT-1) + 1;
slice(:,:,iZ)= bfGetPlane(reader,iPlane);
end
end
end
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment