From 40979a6f4ab678a2f57ccaf8aede1156713e3bf8 Mon Sep 17 00:00:00 2001 From: Tim Date: Fri, 16 Nov 2018 10:14:17 +0100 Subject: add astra_plot_geom command and sample s024 Signed-off-by: Tim --- matlab/algorithms/plot_geom/__main.m | 44 ------ .../plot_geom/astra_create_example_cone.m | 165 +++++++++++++++++++++ .../plot_geom/astra_create_example_fanflat.m | 50 +++++++ .../plot_geom/astra_create_example_parallel3d.m | 30 ++++ matlab/algorithms/plot_geom/create_example_cone.m | 165 --------------------- .../algorithms/plot_geom/create_example_fanflat.m | 50 ------- .../plot_geom/create_example_parallel3d.m | 30 ---- .../algorithms/plot_geom/draw/create_gui_figure.m | 58 -------- .../algorithms/plot_geom/draw/draw_cad_phantom.m | 8 +- .../algorithms/plot_geom/draw/draw_cone_vec_geom.m | 26 ++-- .../plot_geom/draw/draw_fanflat_vec_geom.m | 17 +-- .../plot_geom/draw/draw_parallel3d_vec_geom.m | 14 +- matlab/algorithms/plot_geom/draw/draw_proj_geom.m | 138 +++++++++++++++++ .../algorithms/plot_geom/draw/draw_proj_geometry.m | 142 ------------------ matlab/algorithms/plot_geom/draw/draw_vol_geom.m | 22 ++- matlab/algorithms/plot_geom/stl/bunny.stl | Bin 47484 -> 0 bytes matlab/tools/astra_plot_geom.m | 91 ++++++++++++ 17 files changed, 519 insertions(+), 531 deletions(-) delete mode 100644 matlab/algorithms/plot_geom/__main.m create mode 100644 matlab/algorithms/plot_geom/astra_create_example_cone.m create mode 100644 matlab/algorithms/plot_geom/astra_create_example_fanflat.m create mode 100644 matlab/algorithms/plot_geom/astra_create_example_parallel3d.m delete mode 100644 matlab/algorithms/plot_geom/create_example_cone.m delete mode 100644 matlab/algorithms/plot_geom/create_example_fanflat.m delete mode 100644 matlab/algorithms/plot_geom/create_example_parallel3d.m delete mode 100644 matlab/algorithms/plot_geom/draw/create_gui_figure.m create mode 100644 matlab/algorithms/plot_geom/draw/draw_proj_geom.m delete mode 100644 matlab/algorithms/plot_geom/draw/draw_proj_geometry.m delete mode 100644 matlab/algorithms/plot_geom/stl/bunny.stl create mode 100644 matlab/tools/astra_plot_geom.m (limited to 'matlab') diff --git a/matlab/algorithms/plot_geom/__main.m b/matlab/algorithms/plot_geom/__main.m deleted file mode 100644 index 607a242..0000000 --- a/matlab/algorithms/plot_geom/__main.m +++ /dev/null @@ -1,44 +0,0 @@ -%% main.m -% brief small tool to render astra geometries in matlab -% -% date 20.06.2018 -% author Tim Elberfeld -% imec VisionLab -% University of Antwerp -% -% - last update 09.07.2018 -%% -% close all; - -[h_gui, is_running] = create_gui_figure(); - -%proj_geom = create_example_cone('vec'); -%proj_geom = create_example_cone('normal'); -%proj_geom = create_example_cone('helix'); -proj_geom = create_example_parallel3d('vec'); -%proj_geom = create_example_fanflat('vec'); -%proj_geom = create_example_fanflat(); -%proj_geom = create_example_parallel3d(); - -draw_proj_geometry(proj_geom, h_gui); -hold on; - -vol_magn = 10; -phantom_size = 5; -phantom_px = 1500; -vx_size = phantom_size / phantom_px; % voxel size -vol_geom = astra_create_vol_geom(phantom_px, phantom_px, phantom_px); -line_width = 1; % line width for phantom -draw_vol_geom(vol_geom, vx_size, h_gui, 'Magnification', vol_magn,... - 'LineWidth', line_width, 'Color', 'r'); - -% this magnification is empirically chosen to fit the stl file -cad_magn = 350; -draw_cad_phantom('stl/bunny.stl', cad_magn, h_gui); - -proj_geom = create_example_cone('deform_vec'); -draw_proj_geometry(proj_geom, h_gui, 'VectorIdx', 3, 'Color', 'b',... - 'RotationAxis', [0,0,1]); - -hold off; -axis equal; diff --git a/matlab/algorithms/plot_geom/astra_create_example_cone.m b/matlab/algorithms/plot_geom/astra_create_example_cone.m new file mode 100644 index 0000000..598899f --- /dev/null +++ b/matlab/algorithms/plot_geom/astra_create_example_cone.m @@ -0,0 +1,165 @@ +function [ proj_geom ] = create_example_cone(type) +%% create_example_cone.m +% create an example standard cone beam geometry +% param type type of geometry to create. provided as one of the +% following strings (not case sensitive): +% 'normal' - standard (non vector) geometry +% 'vec' - example vector geometry +% 'helix' - helical trajectory vector geometry +% 'deform_vec' - deformed vector geometry example, +% courtesy of Van Nguyen +% return proj_geom - the geometry that was created +% +% date 21.06.2018 +% author Tim Elberfeld, Alice Presenti, Van Nguyen +% imec VisionLab +% University of Antwerp +% last mod 07.11.2018 +%% + if(nargin < 1) + type = 'normal'; + end + if strcmpi(type, 'normal') + % first, give measurements in mm + det_spacing = [0.035, 0.035]; + detector_px = [1200, 1200]; + angles = linspace2(0, 2*pi, 100); + source_origin = 30; + origin_det = 200; + phantom_size = 5; + + phantom_px = 150; % voxels for the phantom + vx_size = phantom_size / phantom_px; % voxel size + + % now express all measurements in terms of the voxel size + det_spacing = det_spacing ./ vx_size; + origin_det = origin_det ./ vx_size; + source_origin = source_origin ./ vx_size; + + proj_geom = astra_create_proj_geom('cone', det_spacing(1), ... + det_spacing(2), detector_px(1), detector_px(2), angles,... + source_origin, origin_det); + + elseif strcmpi(type, 'vec') + proj_geom = create_standard_vec_geom(); + + elseif strcmpi(type, 'deform_vec') + proj_geom = create_deform_vec_geom(); + + elseif strcmpi(type, 'helix') + detRows = 220; + detCols = 220; + detPitchX = 10^-3 * 11; % microns to mm + detPitchY = 10^-3 * 11; % microns to mm + objSrcDist = 500; % mm + rotStep = 3.6; % deg + numAngles = 200; + zTranslation = 0.5; % mm + zDist = zTranslation * numAngles; + vectors = zeros(numAngles, 12); + translation = -zDist + (zTranslation * (0:(numAngles-1))); + + minAngle = 0; % just assume we start at 0 + maxAngle = (numAngles * rotStep)/180 * pi; % convert deg to rad + angles = linspace(minAngle, maxAngle, 200); + + % source position per angle + vectors(:, 1) = sin(angles) * objSrcDist; + vectors(:, 2) = -cos(angles) * objSrcDist; + vectors(:, 3) = translation; + + % detector position per angle + vectors(:, 4) = 0; + vectors(:, 5) = 0; + vectors(:, 6) = translation; + + % vector from detector pixel (0,0) to (0,1) + vectors(:, 7) = cos(angles); + vectors(:, 8) = sin(angles); + vectors(:, 9) = zeros(numAngles, 1); + + % vector from detector pixel (0,0) to (1, 0) + vectors(:, 10) = zeros(numAngles, 1); + vectors(:, 11) = zeros(numAngles, 1); + vectors(:, 12) = ones(numAngles, 1); + + proj_geom = astra_create_proj_geom('cone_vec', detCols, detRows, vectors); + + else + proj_geom = create_standard_vec_geom(); + end + + function [proj_geom, z_axis] = create_standard_vec_geom() + % geometry settings taken from code from Alice Presenti + settings = struct; + settings.detectorPixelSize = 0.0748; + % detector size and number of projections + settings.projectionSize = [1536 1944 21]; + settings.SOD = 679.238020; %[mm] + settings.SDD = 791.365618; %[mm] + settings.voxelSize = 1; %[mm] + settings.gamma = linspace(0,300,21)*pi/180; + + S0 = zeros(settings.projectionSize(3), 12); + Sorig = [-settings.SOD, 0, 0,... % the ray origin vector (source) + (settings.SDD-settings.SOD), 0, 0,... % detector center + 0, -settings.detectorPixelSize, 0,... % detector u axis + 0, 0, settings.detectorPixelSize]; % detector v axis + + z_axis = [0, 0, 1]; + for i = 1:settings.projectionSize(3) + S0(i,:) = Sorig(:); + S0(i,1:3) = rotate_around3d(S0(i,1:3),... + z_axis, settings.gamma(i)); + S0(i,4:6) = rotate_around3d(S0(i,4:6), z_axis,... + settings.gamma(i)); + S0(i,7:9) = rotate_around3d(S0(i,7:9), z_axis,... + settings.gamma(i)); + S0(i,10:12) = rotate_around3d(S0(i,10:12),... + z_axis, settings.gamma(i)); + end + proj_geom = astra_create_proj_geom('cone_vec', ... + settings.projectionSize(2), settings.projectionSize(1), S0); + end + + function [proj_geom, z_axis] = create_deform_vec_geom() + settings = struct; + settings.detectorPixelSize = 0.0748; + settings.projectionSize = [1536 1944 21]; + settings.SOD = 679.238020; %[mm] + settings.SDD = 791.365618; %[mm] + settings.voxelSize = 1; %[mm] + settings.gamma = linspace(0,300,21)*pi/180; + + S0 = zeros(settings.projectionSize(3), 12); + Sorig = [-settings.SOD, 0, 0,... % the ray origin vector (source) + (settings.SDD-settings.SOD), 0, 0,... % detector center + 0, -settings.detectorPixelSize, 0,... % detector u axis + 0, 0, settings.detectorPixelSize]; % detector v axis + + z_axis = [0, 0, 1]; + for i = 1:settings.projectionSize(3) + S0(i,:) = Sorig(:); + end + + S0(:, 1:3) = translate_3d(S0(:, 1:3), [100, 150, 0]); + S0(:, 4:6) = translate_3d(S0(:, 4:6), [100, 150, 0]); + + S0 = rotate_detector(S0, [0.48,0.32,0]); + S0 = magnify_proj(S0, 100); + for i = 1:settings.projectionSize(3) + S0(i,1:3) = rotate_around3d(S0(i,1:3), z_axis,... + settings.gamma(i)); + S0(i,4:6) = rotate_around3d(S0(i,4:6), z_axis,... + settings.gamma(i)); + S0(i,7:9) = rotate_around3d(S0(i,7:9), z_axis,... + settings.gamma(i)); + S0(i,10:12) = rotate_around3d(S0(i,10:12), z_axis,... + settings.gamma(i)); + end + + proj_geom = astra_create_proj_geom('cone_vec',... + settings.projectionSize(2), settings.projectionSize(1), S0); + end + +end diff --git a/matlab/algorithms/plot_geom/astra_create_example_fanflat.m b/matlab/algorithms/plot_geom/astra_create_example_fanflat.m new file mode 100644 index 0000000..703a618 --- /dev/null +++ b/matlab/algorithms/plot_geom/astra_create_example_fanflat.m @@ -0,0 +1,50 @@ +function [ geom ] = create_example_fanflat( type ) +%% create_example_fanflat.m +% brief create an example geometry of type 'fanflat' +% param type type of geometry to create. provided as one of the +% following strings (not case sensitive): +% 'normal' - standard fanflat geometry +% 'vec' - example vector geometry +% return proj_geom - the geometry that was created +% date 22.06.2018 +% author Tim Elberfeld +% imec VisionLab +% University of Antwerp +% last mod 07.11.2018 +%% + if nargin < 1 + type = 'normal'; + end + + if strcmpi(type, 'normal') + geom = make_normal_geometry(); + elseif strcmpi(type, 'vec') + geom = create_example_fanflat('normal'); + geom = fanflat_to_fanflat_vec(geom); + + else + geom = make_normal_geometry(); + end + + function [geom] = make_normal_geometry() + % first, give measurements in mm + det_spacing = 0.035; + detector_px = 1200; + angles = linspace2(0, 2*pi, 100); + source_origin = 30; + origin_det = 200; + phantom_size = 5; + + phantom_px = 150; % voxels for the phantom + vx_size = phantom_size / phantom_px; % voxel size + + % now express all measurements in terms of the voxel size + det_spacing = det_spacing ./ vx_size; + origin_det = origin_det ./ vx_size; + source_origin = source_origin ./ vx_size; + + geom = astra_create_proj_geom('fanflat', det_spacing, ... + detector_px, angles, source_origin, origin_det); + end +end + diff --git a/matlab/algorithms/plot_geom/astra_create_example_parallel3d.m b/matlab/algorithms/plot_geom/astra_create_example_parallel3d.m new file mode 100644 index 0000000..e6fdf8e --- /dev/null +++ b/matlab/algorithms/plot_geom/astra_create_example_parallel3d.m @@ -0,0 +1,30 @@ +function [ proj_geom ] = create_example_parallel3d( type ) +%% create_example_parallel3d.m +% brief create an example geometry of type 'parallel3d' +% param type type of geometry to create. provided as one of the +% following strings (not case sensitive). +% no arguments will create a standard fanflat geometry. +% 'vec' - example vector geometry + +% return proj_geom - the geometry that was created +% date 22.06.2018 +% author Tim Elberfeld +% imec VisionLab +% University of Antwerp +% last mod 07.11.2018 +%% + if nargin < 1 + type = 'nothing'; + end + + det_spacing = [0.035, 0.035]; + detector_px = [1000, 1000]; + angles = linspace(0, 2*pi, 100); + + proj_geom = astra_create_proj_geom('parallel3d', det_spacing(1),... + det_spacing(2), detector_px(1), detector_px(2), angles); + + if strcmp(type, 'vec') + proj_geom = astra_geom_2vec(proj_geom); + end +end diff --git a/matlab/algorithms/plot_geom/create_example_cone.m b/matlab/algorithms/plot_geom/create_example_cone.m deleted file mode 100644 index 598899f..0000000 --- a/matlab/algorithms/plot_geom/create_example_cone.m +++ /dev/null @@ -1,165 +0,0 @@ -function [ proj_geom ] = create_example_cone(type) -%% create_example_cone.m -% create an example standard cone beam geometry -% param type type of geometry to create. provided as one of the -% following strings (not case sensitive): -% 'normal' - standard (non vector) geometry -% 'vec' - example vector geometry -% 'helix' - helical trajectory vector geometry -% 'deform_vec' - deformed vector geometry example, -% courtesy of Van Nguyen -% return proj_geom - the geometry that was created -% -% date 21.06.2018 -% author Tim Elberfeld, Alice Presenti, Van Nguyen -% imec VisionLab -% University of Antwerp -% last mod 07.11.2018 -%% - if(nargin < 1) - type = 'normal'; - end - if strcmpi(type, 'normal') - % first, give measurements in mm - det_spacing = [0.035, 0.035]; - detector_px = [1200, 1200]; - angles = linspace2(0, 2*pi, 100); - source_origin = 30; - origin_det = 200; - phantom_size = 5; - - phantom_px = 150; % voxels for the phantom - vx_size = phantom_size / phantom_px; % voxel size - - % now express all measurements in terms of the voxel size - det_spacing = det_spacing ./ vx_size; - origin_det = origin_det ./ vx_size; - source_origin = source_origin ./ vx_size; - - proj_geom = astra_create_proj_geom('cone', det_spacing(1), ... - det_spacing(2), detector_px(1), detector_px(2), angles,... - source_origin, origin_det); - - elseif strcmpi(type, 'vec') - proj_geom = create_standard_vec_geom(); - - elseif strcmpi(type, 'deform_vec') - proj_geom = create_deform_vec_geom(); - - elseif strcmpi(type, 'helix') - detRows = 220; - detCols = 220; - detPitchX = 10^-3 * 11; % microns to mm - detPitchY = 10^-3 * 11; % microns to mm - objSrcDist = 500; % mm - rotStep = 3.6; % deg - numAngles = 200; - zTranslation = 0.5; % mm - zDist = zTranslation * numAngles; - vectors = zeros(numAngles, 12); - translation = -zDist + (zTranslation * (0:(numAngles-1))); - - minAngle = 0; % just assume we start at 0 - maxAngle = (numAngles * rotStep)/180 * pi; % convert deg to rad - angles = linspace(minAngle, maxAngle, 200); - - % source position per angle - vectors(:, 1) = sin(angles) * objSrcDist; - vectors(:, 2) = -cos(angles) * objSrcDist; - vectors(:, 3) = translation; - - % detector position per angle - vectors(:, 4) = 0; - vectors(:, 5) = 0; - vectors(:, 6) = translation; - - % vector from detector pixel (0,0) to (0,1) - vectors(:, 7) = cos(angles); - vectors(:, 8) = sin(angles); - vectors(:, 9) = zeros(numAngles, 1); - - % vector from detector pixel (0,0) to (1, 0) - vectors(:, 10) = zeros(numAngles, 1); - vectors(:, 11) = zeros(numAngles, 1); - vectors(:, 12) = ones(numAngles, 1); - - proj_geom = astra_create_proj_geom('cone_vec', detCols, detRows, vectors); - - else - proj_geom = create_standard_vec_geom(); - end - - function [proj_geom, z_axis] = create_standard_vec_geom() - % geometry settings taken from code from Alice Presenti - settings = struct; - settings.detectorPixelSize = 0.0748; - % detector size and number of projections - settings.projectionSize = [1536 1944 21]; - settings.SOD = 679.238020; %[mm] - settings.SDD = 791.365618; %[mm] - settings.voxelSize = 1; %[mm] - settings.gamma = linspace(0,300,21)*pi/180; - - S0 = zeros(settings.projectionSize(3), 12); - Sorig = [-settings.SOD, 0, 0,... % the ray origin vector (source) - (settings.SDD-settings.SOD), 0, 0,... % detector center - 0, -settings.detectorPixelSize, 0,... % detector u axis - 0, 0, settings.detectorPixelSize]; % detector v axis - - z_axis = [0, 0, 1]; - for i = 1:settings.projectionSize(3) - S0(i,:) = Sorig(:); - S0(i,1:3) = rotate_around3d(S0(i,1:3),... - z_axis, settings.gamma(i)); - S0(i,4:6) = rotate_around3d(S0(i,4:6), z_axis,... - settings.gamma(i)); - S0(i,7:9) = rotate_around3d(S0(i,7:9), z_axis,... - settings.gamma(i)); - S0(i,10:12) = rotate_around3d(S0(i,10:12),... - z_axis, settings.gamma(i)); - end - proj_geom = astra_create_proj_geom('cone_vec', ... - settings.projectionSize(2), settings.projectionSize(1), S0); - end - - function [proj_geom, z_axis] = create_deform_vec_geom() - settings = struct; - settings.detectorPixelSize = 0.0748; - settings.projectionSize = [1536 1944 21]; - settings.SOD = 679.238020; %[mm] - settings.SDD = 791.365618; %[mm] - settings.voxelSize = 1; %[mm] - settings.gamma = linspace(0,300,21)*pi/180; - - S0 = zeros(settings.projectionSize(3), 12); - Sorig = [-settings.SOD, 0, 0,... % the ray origin vector (source) - (settings.SDD-settings.SOD), 0, 0,... % detector center - 0, -settings.detectorPixelSize, 0,... % detector u axis - 0, 0, settings.detectorPixelSize]; % detector v axis - - z_axis = [0, 0, 1]; - for i = 1:settings.projectionSize(3) - S0(i,:) = Sorig(:); - end - - S0(:, 1:3) = translate_3d(S0(:, 1:3), [100, 150, 0]); - S0(:, 4:6) = translate_3d(S0(:, 4:6), [100, 150, 0]); - - S0 = rotate_detector(S0, [0.48,0.32,0]); - S0 = magnify_proj(S0, 100); - for i = 1:settings.projectionSize(3) - S0(i,1:3) = rotate_around3d(S0(i,1:3), z_axis,... - settings.gamma(i)); - S0(i,4:6) = rotate_around3d(S0(i,4:6), z_axis,... - settings.gamma(i)); - S0(i,7:9) = rotate_around3d(S0(i,7:9), z_axis,... - settings.gamma(i)); - S0(i,10:12) = rotate_around3d(S0(i,10:12), z_axis,... - settings.gamma(i)); - end - - proj_geom = astra_create_proj_geom('cone_vec',... - settings.projectionSize(2), settings.projectionSize(1), S0); - end - -end diff --git a/matlab/algorithms/plot_geom/create_example_fanflat.m b/matlab/algorithms/plot_geom/create_example_fanflat.m deleted file mode 100644 index 703a618..0000000 --- a/matlab/algorithms/plot_geom/create_example_fanflat.m +++ /dev/null @@ -1,50 +0,0 @@ -function [ geom ] = create_example_fanflat( type ) -%% create_example_fanflat.m -% brief create an example geometry of type 'fanflat' -% param type type of geometry to create. provided as one of the -% following strings (not case sensitive): -% 'normal' - standard fanflat geometry -% 'vec' - example vector geometry -% return proj_geom - the geometry that was created -% date 22.06.2018 -% author Tim Elberfeld -% imec VisionLab -% University of Antwerp -% last mod 07.11.2018 -%% - if nargin < 1 - type = 'normal'; - end - - if strcmpi(type, 'normal') - geom = make_normal_geometry(); - elseif strcmpi(type, 'vec') - geom = create_example_fanflat('normal'); - geom = fanflat_to_fanflat_vec(geom); - - else - geom = make_normal_geometry(); - end - - function [geom] = make_normal_geometry() - % first, give measurements in mm - det_spacing = 0.035; - detector_px = 1200; - angles = linspace2(0, 2*pi, 100); - source_origin = 30; - origin_det = 200; - phantom_size = 5; - - phantom_px = 150; % voxels for the phantom - vx_size = phantom_size / phantom_px; % voxel size - - % now express all measurements in terms of the voxel size - det_spacing = det_spacing ./ vx_size; - origin_det = origin_det ./ vx_size; - source_origin = source_origin ./ vx_size; - - geom = astra_create_proj_geom('fanflat', det_spacing, ... - detector_px, angles, source_origin, origin_det); - end -end - diff --git a/matlab/algorithms/plot_geom/create_example_parallel3d.m b/matlab/algorithms/plot_geom/create_example_parallel3d.m deleted file mode 100644 index e6fdf8e..0000000 --- a/matlab/algorithms/plot_geom/create_example_parallel3d.m +++ /dev/null @@ -1,30 +0,0 @@ -function [ proj_geom ] = create_example_parallel3d( type ) -%% create_example_parallel3d.m -% brief create an example geometry of type 'parallel3d' -% param type type of geometry to create. provided as one of the -% following strings (not case sensitive). -% no arguments will create a standard fanflat geometry. -% 'vec' - example vector geometry - -% return proj_geom - the geometry that was created -% date 22.06.2018 -% author Tim Elberfeld -% imec VisionLab -% University of Antwerp -% last mod 07.11.2018 -%% - if nargin < 1 - type = 'nothing'; - end - - det_spacing = [0.035, 0.035]; - detector_px = [1000, 1000]; - angles = linspace(0, 2*pi, 100); - - proj_geom = astra_create_proj_geom('parallel3d', det_spacing(1),... - det_spacing(2), detector_px(1), detector_px(2), angles); - - if strcmp(type, 'vec') - proj_geom = astra_geom_2vec(proj_geom); - end -end diff --git a/matlab/algorithms/plot_geom/draw/create_gui_figure.m b/matlab/algorithms/plot_geom/draw/create_gui_figure.m deleted file mode 100644 index a9acecf..0000000 --- a/matlab/algorithms/plot_geom/draw/create_gui_figure.m +++ /dev/null @@ -1,58 +0,0 @@ -function [h_ax, running] = create_gui_figure() -%% create_gui_figure.m -% brief gui for the geometry drawing functions -% return h_ax axis to draw into with draw_*() functions -% running state of the rotation, needs to be passed so that -% the callback will work -% date 20.06.2018 -% author Tim Elberfeld -% imec VisionLab -% University of Antwerp -%% - h_figure = figure('name', 'geometry render', 'ButtonDownFcn', @toggle_rotation); - h_ax = axes(h_figure); - - set(h_figure,'CloseRequestFcn', @handle_close_fig) - - xlabel('x axis') - ylabel('y axis') - zlabel('z axis') - - grid on - box off - axis vis3d - axis equal - view(0,0) - - running = false; - do_rotation(); - - function [] = handle_close_fig(h_figure,~) - % this is necessary to stop the rotation before closing the figure - if running - toggle_rotation() - end - - delete(h_figure); - end - - function [] = toggle_rotation(~, ~) - % toggle rotation state - running = ~running; - if running - view(0,0) - do_rotation(); - else - view(45,45); - end - end - - function [] = do_rotation() - % rotate the rendered geometry around the origin - camtarget([0,0,0]); % make origin the camera target and point around which to rotate - while running - camorbit(0.5,0,'camera') - drawnow - end - end -end \ No newline at end of file diff --git a/matlab/algorithms/plot_geom/draw/draw_cad_phantom.m b/matlab/algorithms/plot_geom/draw/draw_cad_phantom.m index 18873a8..92af844 100644 --- a/matlab/algorithms/plot_geom/draw/draw_cad_phantom.m +++ b/matlab/algorithms/plot_geom/draw/draw_cad_phantom.m @@ -1,9 +1,9 @@ -function [] = draw_cad_phantom(filename, magn, h_ax) +function [] = draw_cad_phantom(filename, magn) %% draw_cad_phantom.m % brief render an stl model into a 3d axis object % param vol_geom volume geometry describing the phantom -% param magn magnification multiplier of the phantom. default = 1 % param h_ax handle to axis to plot into +% param magn magnification multiplier of the phantom. default = 1 % % date 02.07.2018 % author Alice Presenti @@ -11,12 +11,10 @@ function [] = draw_cad_phantom(filename, magn, h_ax) % University of Antwerp % Modified by Tim Elberfeld %% + h_ax = gca; if nargin == 1 magn = 1; end - if nargin == 2 - h_ax = axes(gcf); - end [v,f,~,~] = stlRead(filename); m = mean(v); % to center the CAD model! diff --git a/matlab/algorithms/plot_geom/draw/draw_cone_vec_geom.m b/matlab/algorithms/plot_geom/draw/draw_cone_vec_geom.m index c1ed5b0..bee83b7 100644 --- a/matlab/algorithms/plot_geom/draw/draw_cone_vec_geom.m +++ b/matlab/algorithms/plot_geom/draw/draw_cone_vec_geom.m @@ -34,17 +34,15 @@ function [] = draw_cone_vec_geom(h_ax, geom, options) % draw the points and connect with lines hold on; num_angles = size(vectors, 1); - for jj = 1:num_angles - s_source = scatter3(h_ax, xray_source(jj, 1), xray_source(jj, 2),... - xray_source(jj, 3)); - s_source.Marker = options.SourceMarker; - s_source.MarkerEdgeColor = options.SourceMarkerColor; - - s_det = scatter3(h_ax, detector_center(jj, 1),... - detector_center(jj, 2), detector_center(jj, 3)); - s_det.MarkerEdgeColor = options.DetectorMarkerColor; - s_det.Marker = options.DetectorMarker; - end + s_source = scatter3(h_ax, xray_source(:, 1), xray_source(:, 2),... + xray_source(:, 3)); + s_source.Marker = options.SourceMarker; + s_source.MarkerEdgeColor = options.SourceMarkerColor; + + s_det = scatter3(h_ax, detector_center(:, 1),... + detector_center(:, 2), detector_center(:, 3)); + s_det.MarkerEdgeColor = options.DetectorMarkerColor; + s_det.Marker = options.DetectorMarker; detector = struct; detector.u = vectors(options.VectorIdx, 7:9); @@ -56,10 +54,12 @@ function [] = draw_cone_vec_geom(h_ax, geom, options) vertices = draw_detector_vec(h_ax, detector, options); connect_source_detector(h_ax, vertices, detector_center, ... xray_source, options); - distances = eucl_dist3d(detector_center, xray_source); - mean_sdd = mean(distances(:)); % mean source detector distance + % rotation axis will be roughly as long as the source detector distance + distances = eucl_dist3d(detector_center, xray_source); + mean_sdd = mean(distances(:)); % mean source detector distance draw_rotation_axis(h_ax, mean_sdd, options); + text(h_ax, xray_source(options.VectorIdx, 1),... xray_source(options.VectorIdx, 2),... xray_source(options.VectorIdx, 3), 'x-ray source'); diff --git a/matlab/algorithms/plot_geom/draw/draw_fanflat_vec_geom.m b/matlab/algorithms/plot_geom/draw/draw_fanflat_vec_geom.m index 9f5d620..b8261cf 100644 --- a/matlab/algorithms/plot_geom/draw/draw_fanflat_vec_geom.m +++ b/matlab/algorithms/plot_geom/draw/draw_fanflat_vec_geom.m @@ -29,15 +29,13 @@ function [ output_args ] = draw_fanflat_vec_geom( h_ax, geom, options) % draw the points and connect with lines hold on; - for jj = 1:num_angles - scatter3(h_ax, xray_source(jj, 1), xray_source(jj, 2),... - xray_source(jj, 3), options.SourceMarker,... - options.SourceMarkerColor); - scatter3(h_ax, detector_center(jj, 1),... - detector_center(jj, 2), detector_center(jj, 3),... - options.DetectorMarker, options.DetectorMarkerColor); - end - + scatter3(h_ax, xray_source(:, 1), xray_source(:, 2),... + xray_source(:, 3), options.SourceMarker,... + options.SourceMarkerColor); + scatter3(h_ax, detector_center(:, 1),... + detector_center(:, 2), detector_center(:, 3),... + options.DetectorMarker, options.DetectorMarkerColor); + detector = struct; detector.u = [vectors(options.VectorIdx, 5:6), 0]; detector.v = fliplr(detector.u); @@ -49,6 +47,7 @@ function [ output_args ] = draw_fanflat_vec_geom( h_ax, geom, options) connect_source_detector(h_ax, vertices, detector_center,... xray_source, options); + % rotation axis will be roughly as long as the source detector distance distances = eucl_dist3d(detector_center, xray_source); mean_sdd = mean(distances(:)); % mean source detector distance draw_rotation_axis(h_ax, mean_sdd, options); diff --git a/matlab/algorithms/plot_geom/draw/draw_parallel3d_vec_geom.m b/matlab/algorithms/plot_geom/draw/draw_parallel3d_vec_geom.m index 351ec42..ab63e43 100644 --- a/matlab/algorithms/plot_geom/draw/draw_parallel3d_vec_geom.m +++ b/matlab/algorithms/plot_geom/draw/draw_parallel3d_vec_geom.m @@ -23,13 +23,11 @@ function [ ] = draw_parallel3d_vec_geom( h_ax, geom, options) % draw the points and connect with lines hold on; num_angles = size(vectors, 1); - for jj = 1:num_angles - scatter3(h_ax, xray_source(jj, 1), xray_source(jj, 2), xray_source(jj, 3),... - options.SourceMarkerColor,... - options.SourceMarker); - scatter3(h_ax, detector_center(jj, 1), detector_center(jj, 2),... - detector_center(jj, 3), 'k.'); - end + scatter3(h_ax, xray_source(:, 1), xray_source(:, 2),... + xray_source(:, 3), options.SourceMarkerColor,... + options.SourceMarker); + scatter3(h_ax, detector_center(:, 1), detector_center(:, 2),... + detector_center(:, 3), 'k.'); detector = struct; detector.u = vectors(options.VectorIdx, 7:9); @@ -42,9 +40,11 @@ function [ ] = draw_parallel3d_vec_geom( h_ax, geom, options) connect_source_detector(h_ax, detector, detector_center, xray_source,... options); + % rotation axis will be roughly as long as the source detector distance distances = eucl_dist3d(detector_center, xray_source); mean_sdd = mean(distances(:)); % mean source detector distance draw_rotation_axis(h_ax, mean_sdd, options); + text(h_ax, xray_source(options.VectorIdx, 1),... xray_source(options.VectorIdx, 2), ... xray_source(options.VectorIdx, 3), 'x-ray source'); diff --git a/matlab/algorithms/plot_geom/draw/draw_proj_geom.m b/matlab/algorithms/plot_geom/draw/draw_proj_geom.m new file mode 100644 index 0000000..c20ddc8 --- /dev/null +++ b/matlab/algorithms/plot_geom/draw/draw_proj_geom.m @@ -0,0 +1,138 @@ +function [] = draw_proj_geom(geom, varargin) +%% draw_proj_geom.m +% brief rendering function for astra geometries. +% param geom the geometry to plot. If geometry type +% is not supported, throws error +% ------------------------------ +% optional parameters that can be provided as string value pairs: +% +% param RotationAxis if specified, will change the drawn +% rotation axis to provided axis. +% Must be 3-vector. Default value is +% [NaN, NaN, NaN], (meaning do not draw). +% param RotationAxisOffset if specified, will translate the drawn +% rotation axis by the provided vector. +% Default = [0, 0, 0] +% param VectorIdx index of the vector to visualize if geom +% is a vector geometry type. Default = 1 +% param Color Color for all markers and lines if not +% otherwise specified +% param DetectorMarker marker for the detector locations. +% Default = '.' +% param DetectorMarkerColor color specifier for the detector marker. +% Default = 'k' +% param DetectorLineColor color for the lines drawing the detector +% outline +% param DetectorLineWidth line width of detector rectangle +% param SourceMarker marker for the source locations +% param SourceMarkerColor color specifier for the source marker +% param SourceDistance (only for parallel3d and parallel3d_vec) +% distance of source to origin +% param OpticalAxisColor Color for drawing the optical axis +% +% date 20.06.2018 +% author Tim Elberfeld +% imec VisionLab +% University of Antwerp +% +% - last update 07.11.2018 +%% + h_ax = gca; + options = parseoptions(varargin); + + switch geom.type + case 'parallel3d' + disp('type: parallel3d') + disp(['detector spacing: [' num2str(geom.DetectorSpacingX), ', '... + num2str(geom.DetectorSpacingY) ']']); + disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... + num2str(geom.DetectorColCount) ']']); + disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); + disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); + disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); + disp('DistanceOriginDetector inf'); + disp('DistanceOriginSource inf'); + + draw_parallel3d_geom(h_ax, geom, options); + case 'parallel3d_vec' + disp('type: parallel3d_vec') + disp(['detector px: [' num2str(geom.DetectorRowCount), ', '... + num2str(geom.DetectorColCount) ']']); + disp(['# angles: ' num2str(size(geom.Vectors, 1))]); + + draw_parallel3d_vec_geom(h_ax, geom, options); + case 'cone' + disp('type: cone'); + disp(['detector spacing: [' num2str(geom.DetectorSpacingX), ', '... + num2str(geom.DetectorSpacingY) ']']); + disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... + num2str(geom.DetectorColCount) ']']); + disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); + disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); + disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); + disp(['DistanceOriginDetector ' num2str(geom.DistanceOriginDetector)]); + disp(['DistanceOriginSource ' num2str(geom.DistanceOriginSource)]); + + draw_cone_geom(h_ax, geom, options); + case 'cone_vec' + disp('type: cone_vec'); + disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... + num2str(geom.DetectorColCount) ']']); + disp(['# angles: ' num2str(size(geom.Vectors, 1))]); + + draw_cone_vec_geom(h_ax, geom, options); + case 'fanflat' + disp('type: fanflat'); + disp(['detector px: ' num2str(geom.DetectorCount)]); + disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); + disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); + disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); + disp(['DistanceOriginDetector '... + num2str(geom.DistanceOriginDetector)]); + disp(['DistanceOriginSource '... + num2str(geom.DistanceOriginSource)]); + + draw_fanflat_geom(h_ax, geom, options); + case 'fanflat_vec' + disp('type: fanflat_vec'); + disp(['detector px: ' num2str(geom.DetectorCount)]); + disp(['# angles: ' num2str(size(geom.Vectors, 1))]); + + draw_fanflat_vec_geom(h_ax, geom, options); + otherwise + error(['Unknown geometry type ' geom.type]) + end + view(45, 25); % gives nicer default view angle + + function [options] = parseoptions(input_args) + % make an options struct + options = struct; + options.RotationAxis = [NaN, NaN, NaN]; + options.RotationAxisOffset = [0, 0, 0]; + options.VectorIdx = 1; + options.Color = 'k'; + options.DetectorMarker = '.'; + options.DetectorMarkerColor = ''; + options.DetectorLineColor = ''; + options.DetectorLineWidth = 1; + options.SourceMarker = '*'; + options.SourceMarkerColor = ''; + options.SourceDistance = 100; + options.OpticalAxisColor = ''; + options = parseargs(options, input_args{:}); + + % if the color is still empty, replace by global color + if strcmpi(options.DetectorMarkerColor , '') + options.DetectorMarkerColor = options.Color; + end + if strcmpi(options.DetectorLineColor , '') + options.DetectorLineColor = options.Color; + end + if strcmpi(options.SourceMarkerColor , '') + options.SourceMarkerColor = options.Color; + end + if strcmpi(options.OpticalAxisColor , '') + options.OpticalAxisColor = options.Color; + end + end +end diff --git a/matlab/algorithms/plot_geom/draw/draw_proj_geometry.m b/matlab/algorithms/plot_geom/draw/draw_proj_geometry.m deleted file mode 100644 index 0d31fa2..0000000 --- a/matlab/algorithms/plot_geom/draw/draw_proj_geometry.m +++ /dev/null @@ -1,142 +0,0 @@ -function [] = draw_proj_geometry(geom, h_ax, varargin) -%% draw_proj_geometry.m -% brief rendering function for astra geometries. -% param geom the geometry to plot. If geometry type -% is not supported, throws error -% param h_ax handle to axis to plot into -% ------------------------------ -% optional parameters that can be provided as string value pairs: -% -% param RotationAxis if specified, will change the drawn -% rotation axis to provided axis. -% Must be 3-vector. Default value is -% [NaN, NaN, NaN], (meaning do not draw). -% param RotationAxisOffset if specified, will translate the drawn -% rotation axis by the provided vector. -% Default = [0, 0, 0] -% param VectorIdx index of the vector to visualize if geom -% is a vector geometry type. Default = 1 -% param Color Color for all markers and lines if not -% otherwise specified -% param DetectorMarker marker for the detector locations. -% Default = '.' -% param DetectorMarkerColor color specifier for the detector marker. -% Default = 'k' -% param DetectorLineColor color for the lines drawing the detector -% outline -% param DetectorLineWidth line width of detector rectangle -% param SourceMarker marker for the source locations -% param SourceMarkerColor color specifier for the source marker -% param SourceDistance (only for parallel3d and parallel3d_vec) -% distance of source to origin -% param OpticalAxisColor Color for drawing the optical axis -% -% date 20.06.2018 -% author Tim Elberfeld -% imec VisionLab -% University of Antwerp -% -% - last update 07.11.2018 -%% - if nargin == 1 - h_ax = axes(gcf); - end - - options = parseoptions(varargin); - - switch geom.type - case 'parallel3d' - disp('type: parallel3d') - disp(['detector spacing: [' num2str(geom.DetectorSpacingX), ', '... - num2str(geom.DetectorSpacingY) ']']); - disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... - num2str(geom.DetectorColCount) ']']); - disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); - disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); - disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); - disp('DistanceOriginDetector inf'); - disp('DistanceOriginSource inf'); - - draw_parallel3d_geom(h_ax, geom, options); - case 'parallel3d_vec' - disp('type: parallel3d_vec') - disp(['detector px: [' num2str(geom.DetectorRowCount), ', '... - num2str(geom.DetectorColCount) ']']); - disp(['# angles: ' num2str(size(geom.Vectors, 1))]); - - draw_parallel3d_vec_geom(h_ax, geom, options); - case 'cone' - disp('type: cone'); - disp(['detector spacing: [' num2str(geom.DetectorSpacingX), ', '... - num2str(geom.DetectorSpacingY) ']']); - disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... - num2str(geom.DetectorColCount) ']']); - disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); - disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); - disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); - disp(['DistanceOriginDetector ' num2str(geom.DistanceOriginDetector)]); - disp(['DistanceOriginSource ' num2str(geom.DistanceOriginSource)]); - - draw_cone_geom(h_ax, geom, options); - case 'cone_vec' - disp('type: cone_vec'); - disp(['detector px: [' num2str(geom.DetectorRowCount), ', ' ... - num2str(geom.DetectorColCount) ']']); - disp(['# angles: ' num2str(size(geom.Vectors, 1))]); - - draw_cone_vec_geom(h_ax, geom, options); - case 'fanflat' - disp('type: fanflat'); - disp(['detector px: ' num2str(geom.DetectorCount)]); - disp(['angle lo: ' num2str(geom.ProjectionAngles(1))]); - disp(['angle hi: ' num2str(geom.ProjectionAngles(end))]); - disp(['# angles: ' num2str(numel(geom.ProjectionAngles))]); - disp(['DistanceOriginDetector '... - num2str(geom.DistanceOriginDetector)]); - disp(['DistanceOriginSource '... - num2str(geom.DistanceOriginSource)]); - - draw_fanflat_geom(h_ax, geom, options); - case 'fanflat_vec' - disp('type: fanflat_vec'); - disp(['detector px: ' num2str(geom.DetectorCount)]); - disp(['# angles: ' num2str(size(geom.Vectors, 1))]); - - draw_fanflat_vec_geom(h_ax, geom, options); - otherwise - error(['Unknown geometry type ' geom.type]) - end - view(45, 25); % gives nicer default view angle - - function [options] = parseoptions(input_args) - % make an options struct - options = struct; - options.RotationAxis = [NaN, NaN, NaN]; - options.RotationAxisOffset = [0, 0, 0]; - options.VectorIdx = 1; - options.Color = 'k'; - options.DetectorMarker = '.'; - options.DetectorMarkerColor = ''; - options.DetectorLineColor = ''; - options.DetectorLineWidth = 1; - options.SourceMarker = '*'; - options.SourceMarkerColor = ''; - options.SourceDistance = 100; - options.OpticalAxisColor = ''; - options = parseargs(options, input_args{:}); - - % if the color is still empty, replace by global color - if strcmpi(options.DetectorMarkerColor , '') - options.DetectorMarkerColor = options.Color; - end - if strcmpi(options.DetectorLineColor , '') - options.DetectorLineColor = options.Color; - end - if strcmpi(options.SourceMarkerColor , '') - options.SourceMarkerColor = options.Color; - end - if strcmpi(options.OpticalAxisColor , '') - options.OpticalAxisColor = options.Color; - end - end -end diff --git a/matlab/algorithms/plot_geom/draw/draw_vol_geom.m b/matlab/algorithms/plot_geom/draw/draw_vol_geom.m index e6d550d..0b4fe2f 100644 --- a/matlab/algorithms/plot_geom/draw/draw_vol_geom.m +++ b/matlab/algorithms/plot_geom/draw/draw_vol_geom.m @@ -1,17 +1,17 @@ -function [] = draw_vol_geom( vol_geom, vx_size, h_ax, varargin) +function [] = draw_vol_geom( vol_geom, varargin) %% draw_vol_geom.m % brief rendering function for astra volume geometries % describing a phantom. % param vol_geom volume geometry describing the phantom % param vx_size voxel size in unit of preference. must be same unit % that was used to scale the projection geometry -% param h_ax handle to axis to plot into % ------------------------------ % optional parameters that can be provided as string value pairs: % % param Magnification magnification factor for the phantom. for small -% phantoms it might be necessary to scale the render up -% as otherwise it won't show up in the plot. Default = 1 +% phantoms it might be necessary to scale the render +% up as otherwise it won't show up in the plot. +% Default = 1 % param LineWidth line width for the box wireframe. Default = 2 % param Color color of the wireframe. Default = 'r' % @@ -20,12 +20,18 @@ function [] = draw_vol_geom( vol_geom, vx_size, h_ax, varargin) % imec VisionLab % University of Antwerp % -% - last update 09.07.2018 +% - last update 16.11.2018 %% - if nargin < 4 - h_ax = axes(gcf); - end + h_ax = gca; + if mod(size(varargin), 2) ~= 0 + vx_size = varargin{1}; + varargin = varargin(2:end); % consumed vx_size from arg list + else + vx_size = 1; + end + + options = struct; options.Color = 'r'; options.LineWidth = 2; diff --git a/matlab/algorithms/plot_geom/stl/bunny.stl b/matlab/algorithms/plot_geom/stl/bunny.stl deleted file mode 100644 index 0b7fcaa..0000000 Binary files a/matlab/algorithms/plot_geom/stl/bunny.stl and /dev/null differ diff --git a/matlab/tools/astra_plot_geom.m b/matlab/tools/astra_plot_geom.m new file mode 100644 index 0000000..3a654fd --- /dev/null +++ b/matlab/tools/astra_plot_geom.m @@ -0,0 +1,91 @@ +function [] = astra_plot_geom(geometry, varargin) +%-------------------------------------------------------------------------- +% [] = astra_plot_geometry(geometry, varargin) +% +% plot an astra geometry +% +% geometry: any astra geometry, either volume geometry, projection +% geometry or an *.stl file (powered by stlRead). +% varargin: supports a variable number of (ordered and unordered) +% arguments. +% +% the remaining arguments depend on the input: +% if 'geometry' is +% - a volume geometry +% vx_size voxel size in unit of preference. Must +% be same unit that was used to scale the +% projection geometry. +% and as unorderd string-value-pairs +% Magnification magnification factor for the phantom. +% For small phantoms it might be +% necessary to scale the render up as +% otherwise it won't show up in the plot. +% Default = 1 +% LineWidth line width for the box wireframe. +% Default = 2 +% Color color of the wireframe. Default = 'r' +% +% - a projection geometry (as unordered string-value-pairs) +% RotationAxis if specified, will change the drawn +% rotation axis to provided axis. +% Must be 3-vector. Default value is +% [NaN, NaN, NaN], (meaning do not draw). +% RotationAxisOffset if specified, will translate the drawn +% rotation axis by the provided vector. +% Default = [0, 0, 0] +% VectorIdx index of the vector to visualize if +% vector geometry type. Default = 1 +% Color Color for all markers and lines if not +% otherwise specified +% DetectorMarker marker for the detector locations. +% Default = '.' +% DetectorMarkerColor color specifier for the detector marker. +% Default = 'k' +% DetectorLineColor color for the lines drawing the +% detector outline +% DetectorLineWidth line width of detector rectangle +% SourceMarker marker for the source locations +% SourceMarkerColor color specifier for the source marker +% SourceDistance (only for parallel3d and parallel3d_vec) +% distance of source to origin +% OpticalAxisColor Color for drawing the optical axis +% +% - a path to an *.stl file +% magn - magnification factor for vertices in +% CAD file. Default value = 1 +% +%-------------------------------------------------------------------------- +%-------------------------------------------------------------------------- +% This file is part of the ASTRA Toolbox +% +% Copyright: 2010-2018, imec Vision Lab, University of Antwerp +% 2014-2018, CWI, Amsterdam +% License: Open Source under GPLv3 +% Contact: astra@astra-toolbox.com +% Website: http://www.astra-toolbox.com/ +%-------------------------------------------------------------------------- + addpath(genpath('../algorithms/plot_geom')); % add plot_geom tools to matlab path + + if is_vol_geom(geometry) + draw_vol_geom(geometry, varargin{:}); + elseif is_proj_geom(geometry) + draw_proj_geom(geometry, varargin{:}); + elseif ischar(geometry) % assume 'geometry' is a path to a CAD file + draw_cad_phantom(geometry, varargin{:}); + end + + % ---- helper functions ---- + function [ res ] = is_vol_geom(geom) + res = false; + if sum(isfield(geom, {'GridRowCount', 'GridColCount'})) == 2 + res = true; + end + end + + function [ res ] = is_proj_geom(geom) + res = false; + if isfield(geom, 'type') + res = true; + end + end +end -- cgit v1.2.3