-
Notifications
You must be signed in to change notification settings - Fork 2
/
FociMeasurement.m
87 lines (73 loc) · 15.3 KB
/
FociMeasurement.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
function [fociVols, fociNum, fociNum_norm, volCells]= FociMeasurement(folder, segPlots, timepoints, zSteps, fociThresholdScale, backgroundThreshScale, cellVolThresh, cellHoleThresh, minFociVox, maxFociVox, scrsz, saveSeg, stamp)
%% Adam Tyson | 22/11/2017 | [email protected]
%% get the first file
% bioformats should (usually) pick up all other files in the time series(if exported as OME-TIFF)
cd(folder)
fileDir=dir('*T00*.tif');
file=fileDir.name;
%% Load and convert data
disp('Loading data')
datacell = bfopen(file); % load a single time point (3D), exported from slidebook
arraySize=[size(datacell{1,1}{1,1}) zSteps timepoints];
rawData=zeros(arraySize);
for t=1:timepoints
for z=1:zSteps
% arrayPos=z+(t-1)*timepoints;
arrayPos=z+(t-1)*zSteps;
rawData(:,:,z,t)=double(datacell{1,1}{arrayPos,1}); % make a 4D array, x,y,z,t
end
end
%% Segment foci
disp('Segmenting foci')
[fociSegData, ~]=segFoci4D(rawData,'bi', fociThresholdScale, minFociVox, maxFociVox);
%% Display segmentation
switch segPlots
case 'All'
% Calculate display size
dispScale=(scrsz(4)/arraySize(1))*0.8;
screenSize=[10 10 dispScale*arraySize(2) dispScale*arraySize(1)];
for t=1:timepoints
figTitle=[file '_TimePoint_' num2str(t)];
doubleColorMap(rawData(:,:,:,t), fociSegData(:,:,:,t), screenSize,figTitle)
end
case 'MIP'
% Calculate display size
dispScale=(scrsz(4)/arraySize(1))*0.8;
screenSize=[10 10 dispScale*arraySize(2) dispScale*arraySize(1)];
figTitle=[file '_MIP'];
rawMIP = squeeze(max(rawData, [], 3));
segMIP = squeeze(max(fociSegData, [], 3));
doubleColorMap(rawMIP, segMIP, screenSize,figTitle)
otherwise
end
%% analyse foci
disp('Analysing foci')
for t=1:timepoints
fociProps(t) =analyseFoci(rawData(:,:,:,t), fociSegData(:,:,:,t));
end
%% normalise number of foci to cell volume
disp('Segmenting cells')
fociVols=zeros(timepoints,1);
fociNum_norm=zeros(timepoints,1);
fociNum=zeros(timepoints,1);
[cellSeg, volCells]=cellsegT0(rawData(:,:,:,1),backgroundThreshScale, cellVolThresh, cellHoleThresh);
for t=1:timepoints
fociVols(t)=fociProps(t).meanVol;
fociNum(t)=fociProps(t).numFoci;
fociNum_norm(t)=volCells*fociProps(t).numFoci;
end
%% save segmentation files
[~, file, ~]=fileparts(file); if strcmp(saveSeg, 'Yes')
disp('Saving segmentation')
for t=1:timepoints
for frame=1:zSteps
outfileFoci=['Foci_' file '_t_' num2str(t-1) '_' stamp '.tif'];
imwrite(fociSegData(:,:,frame,t),outfileFoci, 'tif', 'WriteMode', 'append', 'compression', 'none');
end
end
for frame=1:zSteps
outfileCell=['CellSegmentation_' file '_' stamp '.tif'];
imwrite(cellSeg(:,:,frame),outfileCell, 'tif', 'WriteMode', 'append', 'compression', 'none');
end
end
end