phantom.m
1.11 KB
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
ps = [128, 128, 128];
n = 50;
I1 = zeros(ps);
for i = 1:n
r1 = randi([4 10]);
r2 = r1;
range1 = -r1:r1;
range2 = -r2:r2;
[cx, cy, cz] = meshgrid(range1, range1, range2);
c1 = zeros(2*r1+1,2*r1+1,2*r2+1);
c1(cx.^2/r1^2 + cy.^2/r1^2 + cz.^2/r2^2 <= 1) = 1;
% c2 = zeros(2*r+1,2*r+1,2*r+1);
% c2(cx.^2 + cy.^2 + cz.^2 < (r-2)^2) = 230;
% c = c1 - c2;
% Get random coordinates in the big image where we will place the upper left corner of the small circle.
dx = ps(1) - 2 * r1;
dy = ps(2) - 2 * r1;
dz = ps(3) - 2 * r2;
sx = dx * rand(n, 1);
sy = dy * rand(n, 1);
sz = dz * rand(n, 1);
% Find the square in the big image where we're going to add a small circle.
x1 = int16(sx(i));
y1 = int16(sy(i));
z1 = int16(sz(i));
x2 = int16(x1 + 2*r1);
y2 = int16(y1 + 2*r1);
z2 = int16(z1 + 2*r2);
% Add in one small circle to the existing big image
I1(x1:x2, y1:y2, z1:z2) = I1(x1:x2, y1:y2, z1:z2) + c1;
I = single(I1);
% I1(I1==0) = 150;
end
I (I>1) = 1;
fid = fopen('128-128-128/s50-r4-10.vol', 'w');
fwrite(fid, I1, '*single');
fclose(fid);
% -------------------------------------------------------------------