-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathtileQualityControl.m
131 lines (113 loc) · 4.6 KB
/
tileQualityControl.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
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
function [ outputCode,outputMsg,varargout ] = tileQualityControl( tileFile,settingsFile, ~, paramA , paramB, varargin )
%tileQualityControl Monitors quality of aquired image tiles.
%Function reads in provided tiff files and analyzes different aspect of the imaging quality for
%instance the obstruction of the objective.
tileInfo = [];
varargout{1} = [];
%% Load settings file
try
jsonText = fileread(fullfile(settingsFile));
Settings = jsondecode(jsonText);
catch
error('Could not read %s',fullfile(settingsFile));
end
outParam = Settings.outParam;
% Set default output code (Okay).
outputCode = outParam.default.code; outputMsg = outParam.default.msg;
% Process tilefile.
tileFile = strrep(tileFile,'.0.',sprintf('.%i.',Settings.Channel));
tileFile = fullfile(tileFile);
%% Store usefull info on sample.
sampleID = regexp(tileFile,'\d{4}-\d{2}-\d{2}','match');
sampleID = sampleID{1};
tileInfo.sampleID = char(sampleID);
tileInfo.mainFolder = char(regexp(tileFile,'.*(?=\d{4}-\d{2}-\d{2}\\\d{2}\\\d{5}\\)','match'));
tileInfo.folder = char(regexp(tileFile,'\d{4}-\d{2}-\d{2}\\\d{2}\\\d{5}','match'));
qcFolder = fullfile(Settings.storeFolder,tileInfo.sampleID);
if isempty(dir(qcFolder)), mkdir(qcFolder);end
%% Open/create logging file.
logFile = fullfile(qcFolder,'log.txt');
fid = fopen(logFile,'a');
logMessage(fid,sprintf('%s Tile: %s',Settings.Name, tileInfo.folder),true);
c = onCleanup(@()fclose(fid)); % Close log file on cleanup.
%% Try to load previous store or create new.
try
qcFile = fullfile(qcFolder,'QC.mat');
if isempty(dir(qcFile))
QC = [];
save(qcFile,'QC');
else
load(qcFile);
end
catch
[ outputCode,outputMsg ] = processError( outParam.storage,fid,tileInfo,outputCode,outputMsg,Settings );
return
end
%% Read image data.
try
I = Tiff.readTifFast(tileFile, Settings.ImageSize,Settings.FrameAvg(1):Settings.FrameAvg(2), 'uint16');
Iavg = uint16(mean(I,3));
catch
%% Process fault.
[ outputCode,outputMsg ] = processError( outParam.imageFile,fid,tileInfo,outputCode,outputMsg,Settings );
return
end
%% Get info .microscope file (required).
protoLoc = [tileFile(1:end-6),'.microscope'];
try
tileInfo.FOV = ProtoBuf.fastProtoBuf( protoLoc, {'x_size_um','y_size_um','x_overlap_um','y_overlap_um'} );
tileInfo.autotile = ProtoBuf.fastProtoBuf( protoLoc, {{'autotile','intensity_threshold'},{'autotile','area_threshold'}} );
tileInfo.trip_detect = ProtoBuf.fastProtoBuf( protoLoc, {{'trip_detect','intensity_threshold'}} );
tileInfo.pos_mm = ProtoBuf.fastProtoBuf( protoLoc,{{'last_target_mm','x'},{'last_target_mm','y'},{'last_target_mm','z'}});
tileInfo.x_pix_size = tileInfo.FOV.x_size_um/Settings.ImageSize(2); tileInfo.y_pix_size = tileInfo.FOV.y_size_um/Settings.ImageSize(1);
catch
%% Process fault.
[ outputCode,outputMsg ] = processError( outParam.microMissing,fid,tileInfo,outputCode,outputMsg,Settings );
return
end
%% Get Latice info (optional).
protoLoc = [tileFile(1:end-6),'.acquisition'];
try
tileInfo.pos_lat = ProtoBuf.fastProtoBuf( protoLoc, {{'current_lattice_position','x'},{'current_lattice_position','y'},{'current_lattice_position','z'}} );
catch
tileInfo.pos_lat = struct('x',[],'y',[],'z',[]);
end
%% Calls to quality checks %%%%
%% Slice Thickness check.
if ~isempty(QC) && ~isempty(tileInfo.pos_lat.z)
[code,~,~] = sliceCheck( tileInfo, QC, outParam.sliceThick.threshold);
if code ~= 100
[ outputCode,outputMsg ] = processError( outParam.sliceThick,fid,tileInfo,outputCode,outputMsg,Settings );
if outputCode~=100
return
end
end
end
%% Line offset check.
[code, ~] = lineOffsetCheck( I, tileInfo, outParam.lineOff.threshold );
if code ~= 100
[ outputCode,outputMsg ] = processError( outParam.lineOff,fid,tileInfo,outputCode,outputMsg,Settings );
if outputCode~=100
return
end
end
%% Blocked Objective detection.
[ code, ~, tileInfo ] = blockDetection( Iavg, tileInfo,QC, outParam.block.threshold, fid );
%% Store data.
QC = [QC;tileInfo];
try
save(qcFile,'QC');
catch
[ outputCode,outputMsg ] = processError( outParam.storage,fid,tileInfo,outputCode,outputMsg,Settings );
end
%% throwing error code of blocked objective after save for counter persistence.
if code ~=100
[ outputCode,outputMsg ] = processError( outParam.block,fid,tileInfo,outputCode,outputMsg,Settings );
if outputCode~=100
return
end
end
%% Send default output if no error.
if outputCode == 100, varargout{1} = Iavg; end
fprintf('\n');
end