Hi, I'd like to register as a developer on the Octave-forge project, SF
username is lxop.
Attached is my version of the Matlab function `phantom' from the image
processing toolbox.
Thanks
Alex
function p = phantom (varargin)
% P = phantom ('Shepp-Logan', n) produces the Shepp-Logan phantom, with
% size n x n. If n is omitted, 256 is used.
%
% P = phantom ('Modified Shepp-Logan', n) produces a modified version of
% the Shepp-Logan phantom which has better contrast than the original,
% with size n x n. If n is omitted, 256 is used.
%
% P = phantom (ellipses, n) produces a custom phantom using the ellipses
% described in `ellipses'. Each row of `ellipses' describes one ellipse,
% and must have 6 columns:
% I, a, b, x0, y0, phi
% where I is the additive intensity of the ellipse
% a is the length of the major axis
% b is the length of the minor axis
% x0 is the horizontal offset of the centre of the ellipse
% y0 is the vercal offset of the centre of the ellipse
% phi is the counterclockwise rotation of the ellipse in degrees,
% measured as the angle between the x axis and the ellipse major axis.
% The image bounding box in the algorithm is {[-1, -1], [1, 1]}, so the
% values of a, b, x0, y0 should all be specified with this in mind.
% If n is omitted, 256 is used.
%
% P = phantom (n) creates a modified Shepp-Logan phantom with size n x n.
%
% P = phantom () creates a modified Shepp-Logan phantom with size 256 x 256.
%
%
% Written by Alex Opie, 2010/03/03
% References:
% Shepp, L. A.; Logan, B. F.; Reconstructing Interior Head Tissue
% from X-Ray Transmissions, IEEE Transactions on Nuclear Science,
% Feb. 1974, p. 232.
%
% Toft, P.; "The Radon Transform - Theory and Implementation", Ph.D. thesis,
% Department of Mathematical Modelling, Technical University
% of Denmark, June 1996.
[n, ellipses] = read_args (varargin {:});
% Blank image
p = zeros (n);
% Create the pixel grid
xvals = (-1 : 2 / (n - 1) : 1);
xgrid = repmat (xvals, n, 1);
for i = 1:size (ellipses, 1)
I = ellipses (i, 1);
a2 = ellipses (i, 2)^2;
b2 = ellipses (i, 3)^2;
x0 = ellipses (i, 4);
y0 = ellipses (i, 5);
phi = ellipses (i, 6) * pi / 180; % Rotation angle in radians
% Create the offset x and y values for the grid
x = xgrid - x0;
y = rot90 (xgrid) - y0;
cos_p = cos (phi);
sin_p = sin (phi);
% Find the pixels within the ellipse
locs = find (((x .* cos_p + y .* sin_p).^2) ./ a2 ...
+ ((y .* cos_p - x .* sin_p).^2) ./ b2 <= 1);
% Add the ellipse intensity to those pixels
p (locs) = p (locs) + I;
end
function [n, ellip] = read_args (varargin)
n = 256;
ellip = mod_shepp_logan ();
if (nargin == 1)
if (ischar (varargin {1}))
ellip = select_phantom (varargin {1});
elseif (numel (varargin {1} == 1))
n = varargin {1};
else
if (size (varargin {1}, 1) != 6)
error ("Wrong number of columns in user
phantom");
end
ellip = varargin {1};
end
elseif (nargin == 2)
n = varargin {2};
if (ischar (varargin {1}))
ellip = select_phantom (varargin {1});
else
if (size (varargin {1}, 1) != 6)
error ("Wrong number of columns in user
phantom");
end
ellip = varargin {1};
end
elseif (nargin > 2)
warning ("Extra arguments passed to phantom were ignored");
end
end
function e = select_phantom (name)
if (strcmpi (name, 'shepp-logan'))
e = shepp_logan ();
elseif (strcmpi (name, 'modified shepp-logan'))
e = mod_shepp_logan ();
else
error ("Unknown phantom type: %s", name);
end
end
function e = shepp_logan
% Standard head phantom, taken from Shepp & Logan
e = [ 2, .69, .92, 0, 0, 0;
-.98, .6624, .8740, 0, -.0184, 0;
-.02, .1100, .3100, .22, 0, -18;
-.02, .1600, .4100, -.22, 0, 18;
.01, .2100, .2500, 0, .35, 0;
.01, .0460, .0460, 0, .1, 0;
.02, .0460, .0460, 0, -.1, 0;
.01, .0460, .0230, -.08, -.605, 0;
.01, .0230, .0230, 0, -.606, 0;
.01, .0230, .0460, .06, -.605, 0];
end
function e = mod_shepp_logan
% Modified version of Shepp & Logan's head phantom,
% adjusted to improve contrast. Taken from Toft.
e = [ 1, .69, .92, 0, 0, 0;
-.80, .6624, .8740, 0, -.0184, 0;
-.20, .1100, .3100, .22, 0, -18;
-.20, .1600, .4100, -.22, 0, 18;
.10, .2100, .2500, 0, .35, 0;
.10, .0460, .0460, 0, .1, 0;
.10, .0460, .0460, 0, -.1, 0;
.10, .0460, .0230, -.08, -.605, 0;
.10, .0230, .0230, 0, -.606, 0;
.10, .0230, .0460, .06, -.605, 0];
end
%function e = ??
% % Add any further phantoms of interest here
% e = [ 0, 0, 0, 0, 0, 0;
% 0, 0, 0, 0, 0, 0];
%end
------------------------------------------------------------------------------
Download Intel® Parallel Studio Eval
Try the new software tools for yourself. Speed compiling, find bugs
proactively, and fine-tune applications for parallel performance.
See why Intel Parallel Studio got high marks during beta.
http://p.sf.net/sfu/intel-sw-dev
_______________________________________________
Octave-dev mailing list
Octave-dev@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/octave-dev
