|
1 | | -function results_dict = measure_LCD(signal_present_array, signal_absent_array, ground_truth, observers, n_reader, pct_split, seed_split) |
| 1 | +function results = measure_LCD(varargin) |
2 | 2 | % given a dataset calculate low contrast detectability as auc curves and return as a table ready for saving or plotting |
3 | 3 | % |
| 4 | +% Usage: |
| 5 | +% 1. measure_LCD(signal_present_array, signal_absent_array, ground_truth, observers, n_reader, pct_split, seed_split) |
| 6 | +% 2. measure_LCD(base_directory, observers, ground_truth, offset, n_reader, pct_split, seed_split) |
| 7 | +% |
4 | 8 | % :param signal_present_array: image stack of signal present images |
5 | 9 | % :param signal_absent_array: corresponding image stack of signal absent images |
| 10 | +% :param base_directory: path to directory containing 'signal_present' and 'signal_absent' subfolders |
6 | 11 | % :param observers: list of observer objects or strings of name of observers. Options: LG_CHO_2D, DOG_CHO_2D, GABOR_CHO_2D |
7 | | -% :param ground_truth: image or filename of image with no noise of MITA LCD phantom, see `approximate_groundtruth` for details on how to turn repeat scans into a ground truth image |
| 12 | +% :param ground_truth: image or filename of image with no noise of MITA LCD phantom |
| 13 | +% :param offset: offset value to subtract from loaded images (only used with directory input) |
8 | 14 | % :param n_reader: number of readers (default is 10) |
9 | | -% :param pct_split: percent of images to be used for training, remainder (1 - split_pct) to be used for testing (default is 0.5) |
10 | | -% :param seed_split: 1d vector containing 'nreader' of random seed values. (defaults to randomly selected seed) |
| 15 | +% :param pct_split: percent of images to be used for training (default is 0.5) |
| 16 | +% :param seed_split: 1d vector containing 'nreader' of random seed values |
11 | 17 | % |
12 | | -% :return res_table: table ready for saving or plotting |
| 18 | +% :return results: table (or struct in Octave) ready for saving or plotting |
13 | 19 |
|
14 | 20 | if is_octave |
15 | 21 | if length(pkg('list', 'image')) < 1 |
|
18 | 24 | pkg load image |
19 | 25 | end |
20 | 26 |
|
21 | | -observer_idx = 1; |
22 | | -for i=1:length(observers) |
23 | | - switch upper(observers{i}) |
24 | | - case 'LG_CHO_2D' |
25 | | - observers{observer_idx} = LG_CHO_2D(); |
26 | | - observer_idx = observer_idx + 1; |
27 | | - case 'NPWE_2D' |
28 | | - observers{observer_idx} = NPWE_2D(); |
29 | | - observer_idx = observer_idx + 1; |
30 | | - case 'DOG_CHO_2D' |
31 | | - observers{observer_idx} = DOG_CHO_2D(); |
32 | | - observer_idx = observer_idx + 1; |
33 | | - case 'GABOR_CHO_2D' |
34 | | - observers{observer_idx} = GABOR_CHO_2D(); |
35 | | - observer_idx = observer_idx + 1; |
36 | | - otherwise |
37 | | - error([observers{i} ' not in LG_CHO_2D, DOG_CHO_2D, GABOR_CHO_2D, NPWE_2D']) |
| 27 | +% Parse inputs |
| 28 | +arg1 = varargin{1}; |
| 29 | +dose_level_extracted = NaN; |
| 30 | + |
| 31 | +if ischar(arg1) || (isstring(arg1) && isscalar(arg1)) |
| 32 | + % Signature 2: Directory input |
| 33 | + base_dir = char(arg1); |
| 34 | + observers = varargin{2}; |
| 35 | + ground_truth = varargin{3}; |
| 36 | + if length(varargin) >= 4 |
| 37 | + offset = varargin{4}; |
| 38 | + else |
| 39 | + offset = 0; |
38 | 40 | end |
39 | | -end |
| 41 | + if length(varargin) >= 5; n_reader = varargin{5}; end |
| 42 | + if length(varargin) >= 6; pct_split = varargin{6}; end |
| 43 | + if length(varargin) >= 7; seed_split = varargin{7}; end |
40 | 44 |
|
41 | | -if ~exist('observers', 'var') |
42 | | - observers = {LG_CHO_2D()}; |
43 | | -% observers = {LG_CHO_2D(),... % will need to modify .channel_width attribute later when insert width known |
44 | | -% DOG_CHO_2D(),... |
45 | | -% GABOR_CHO_2D(),... |
46 | | -% }; |
47 | | -end |
48 | | -%% check for ground truth |
49 | | -if ~exist('ground_truth', 'var') |
50 | | - ground_truth_fname = fullfile(base_dir, 'ground_truth.mhd'); |
51 | | - ground_truth = mhd_read_image(ground_truth_fname); %need to build in offset to dataset |
52 | | -end |
| 45 | + % Load images |
| 46 | + sp_dir = fullfile(base_dir, 'signal_present'); |
| 47 | + sa_dir = fullfile(base_dir, 'signal_absent'); |
53 | 48 |
|
54 | | -if ~exist('n_reader','var') |
55 | | - n_reader = 10; |
56 | | -end |
57 | | -if ~exist('pct_split','var') |
58 | | - pct_split = 0.5; |
| 49 | + sp_mhd = fullfile(sp_dir, 'signal_present.mhd'); |
| 50 | + if exist(sp_mhd, 'file') |
| 51 | + signal_present_array = mhd_read_image(sp_mhd); |
| 52 | + else |
| 53 | + error(['Could not find ' sp_mhd]); |
| 54 | + end |
| 55 | + |
| 56 | + sa_mhd = fullfile(sa_dir, 'signal_absent.mhd'); |
| 57 | + if exist(sa_mhd, 'file') |
| 58 | + signal_absent_array = mhd_read_image(sa_mhd); |
| 59 | + else |
| 60 | + error(['Could not find ' sa_mhd]); |
| 61 | + end |
| 62 | + |
| 63 | + if offset ~= 0 |
| 64 | + signal_present_array = signal_present_array - offset; |
| 65 | + signal_absent_array = signal_absent_array - offset; |
| 66 | + end |
| 67 | + |
| 68 | + % Try to extract dose from base_dir path |
| 69 | + [~, name, ~] = fileparts(base_dir); |
| 70 | + % Check for dose_XXX format |
| 71 | + tokens = regexp(name, 'dose_(\d+)', 'tokens'); |
| 72 | + if ~isempty(tokens) |
| 73 | + dose_level_extracted = str2double(tokens{1}{1}); |
| 74 | + end |
| 75 | + |
| 76 | +else |
| 77 | + % Signature 1: Array input |
| 78 | + signal_present_array = varargin{1}; |
| 79 | + signal_absent_array = varargin{2}; |
| 80 | + ground_truth = varargin{3}; |
| 81 | + observers = varargin{4}; |
| 82 | + if length(varargin) >= 5; n_reader = varargin{5}; end |
| 83 | + if length(varargin) >= 6; pct_split = varargin{6}; end |
| 84 | + if length(varargin) >= 7; seed_split = varargin{7}; end |
59 | 85 | end |
60 | | -if ~exist('seed_val','var') |
61 | | - seed_split = randi(10000, n_reader, 1); |
| 86 | + |
| 87 | +% Default values |
| 88 | +if ~exist('n_reader','var'); n_reader = 10; end |
| 89 | +if ~exist('pct_split','var'); pct_split = 0.5; end |
| 90 | +if ~exist('seed_split','var'); seed_split = randi(10000, n_reader, 1); end |
| 91 | + |
| 92 | +% Instantiate observers if strings are passed |
| 93 | +observer_objs = cell(1, length(observers)); |
| 94 | +for i=1:length(observers) |
| 95 | + obs = observers{i}; |
| 96 | + if ischar(obs) || (isstring(obs) && isscalar(obs)) |
| 97 | + switch upper(obs) |
| 98 | + case 'LG_CHO_2D' |
| 99 | + observer_objs{i} = LG_CHO_2D(); |
| 100 | + case 'NPWE_2D' |
| 101 | + observer_objs{i} = NPWE_2D(); |
| 102 | + case 'DOG_CHO_2D' |
| 103 | + observer_objs{i} = DOG_CHO_2D(); |
| 104 | + case 'GABOR_CHO_2D' |
| 105 | + observer_objs{i} = GABOR_CHO_2D(); |
| 106 | + otherwise |
| 107 | + error([obs ' not in LG_CHO_2D, DOG_CHO_2D, GABOR_CHO_2D, NPWE_2D']) |
| 108 | + end |
| 109 | + else |
| 110 | + % Assuming it is already an object |
| 111 | + observer_objs{i} = obs; |
| 112 | + end |
62 | 113 | end |
63 | | -% input is a binary mask specifying signal known exactly (SKE) |
| 114 | +observers = observer_objs; |
64 | 115 |
|
| 116 | +% input is a binary mask specifying signal known exactly (SKE) |
65 | 117 | truth_masks = get_demo_truth_masks(ground_truth); |
66 | 118 |
|
67 | 119 | insert_rs = []; |
|
128 | 180 | end |
129 | 181 | end |
130 | 182 |
|
131 | | -results_dict.observer=observer; |
132 | | -results_dict.insert_HU=insert_HU'; |
133 | | -results_dict.snr=snr'; |
134 | | -results_dict.auc=auc'; |
135 | | -results_dict.reader=reader'; |
136 | | -results_dict; |
| 183 | +results_dict.observer = observer; |
| 184 | +results_dict.insert_HU = insert_HU; |
| 185 | +results_dict.snr = snr; |
| 186 | +results_dict.auc = auc; |
| 187 | +results_dict.reader = reader; |
| 188 | +results_dict.diameter = insert_diameter_pix; |
| 189 | +results_dict.dose_level = repmat(dose_level_extracted, length(auc), 1); |
| 190 | + |
| 191 | +if is_octave |
| 192 | + pkg load tablicious |
| 193 | +end |
| 194 | + |
| 195 | +results = struct2table(results_dict); |
137 | 196 |
|
138 | 197 | end |
0 commit comments