diff --git a/examples/scan_shepp_logan.m b/examples/scan_shepp_logan.m index c7440ee..98bb059 100644 --- a/examples/scan_shepp_logan.m +++ b/examples/scan_shepp_logan.m @@ -2,7 +2,7 @@ vox_arr_center = zeros(3, 1); phantom_size = 500; voxel_size = 1; -water = material("water"); +water = get_material("water"); % Create voxel array voxel_generator = voxel_shepp_logan(vox_arr_center, phantom_size, voxel_size); diff --git a/src/detectors/detector.m b/src/detectors/detector.m index 645cf3d..85db422 100644 --- a/src/detectors/detector.m +++ b/src/detectors/detector.m @@ -96,7 +96,7 @@ function scan = air_scan(self) % air_scan Generate a scan of air dtd = self.dist_to_detector; - air = material("air"); + air = get_material("air"); air_cylinder = @(i,j,k,e) air.get_mu(e); array = voxel_array(zeros(3, 1), zeros(3,1)+dtd*2, dtd/100, air_cylinder); scan = zeros(self.ny_pixels, self.nz_pixels, self.num_rotations); diff --git a/src/get_material.m b/src/get_material.m index e8d2189..719462a 100644 --- a/src/get_material.m +++ b/src/get_material.m @@ -1,15 +1,18 @@ -classdef material +classdef get_material %MATERIAL A class to represent a material and its properties necessary % for ray tracing and scattering properties + get_mu + end + + properties (Access=protected) atomic_numbers mass_fractions - get_mu density end - properties (Access=private) + properties (Access=private, Constant) known_materials = {'air','blood','bone','lung','muscle','water'}; known_densitys = [1.205E-03, 1.06, 1.92, 1.05, 1.05, 1.00]; known_atomic_numbers = {... @@ -47,7 +50,7 @@ end methods - function self = material(varargin) + function self = get_material(varargin) %MATERIAL Construct a material object % material("material_name") creates a material object with the properties of the given material, if it is available % material(atomic_numbers, mass_fractions, density) creates a material object with the given properties, using the PhotonAttenuation package @@ -56,30 +59,28 @@ error(aac, 'MATLAB:notEnoughInputs', 'Not enough input arguments.') elseif nargin == 1 material_name = varargin{1}; - assert(isstring(material_name), 'When only one input is given, it should be a string with the name of the material'); + assert(isstring(material_name), 'assert:failure', 'When only one input is given, it should be a string with the name of the material.'); material_index = find(self.known_materials == lower(material_name)); if any(material_index) self.atomic_numbers = self.known_atomic_numbers{material_index}; self.mass_fractions = self.known_mass_fractions{material_index}; self.density = self.known_densitys(material_index); else - aac = matlab.lang.correction.AppendArgumentsCorrection('"blood"'); - list_available_materials = strjoin(self.known_materials, ', '); - error(aac, 'MATLAB:invalidMaterial', ... + error('MATLAB:invalidMaterial', ... 'The material %s is not available. Available materials are: %s', ... - material_name, list_available_materials); + material_name, strjoin(self.known_materials, ', ')); end elseif nargin == 3 self.atomic_numbers = varargin{1}; self.mass_fractions = varargin{2}; - self.mass_fractions = self.mass_fractions / sum(self.mass_fractions); self.density = varargin{3}; - assert(isnumeric(self.density) && isscalar(self.density), 'The density should be a scalar number'); - assert(isvector(self.atomic_numbers), 'The atomic numbers should be a vector'); - assert(isvector(self.mass_fractions), 'The mass fractions should be a vector'); - assert(length(self.atomic_numbers) == length(self.mass_fractions), 'The atomic numbers and mass fractions should have the same length'); + assert(isvector(self.atomic_numbers), 'assert:failure', 'The atomic numbers should be a vector.'); + assert(isvector(self.mass_fractions), 'assert:failure', 'The mass fractions should be a vector.'); + assert(length(self.atomic_numbers) == length(self.mass_fractions), 'assert:failure', 'The atomic numbers and mass fractions should have the same length.'); + assert(isnumeric(self.density) && isscalar(self.density), 'assert:failure', 'The density should be a scalar number.'); + self.mass_fractions = self.mass_fractions / sum(self.mass_fractions); else - error('MATLAB:invalidInput', 'Invalid number of input arguments'); + error('MATLAB:invalidInput', 'Invalid number of input arguments.'); end self.get_mu = @(E) sum(PhotonAttenuationQ(self.atomic_numbers, E, 1) .* self.mass_fractions) * self.density; end diff --git a/tests/attenuation_tests.m b/tests/attenuation_tests.m index 3e5653f..7ef1046 100644 --- a/tests/attenuation_tests.m +++ b/tests/attenuation_tests.m @@ -1,9 +1,32 @@ classdef attenuation_tests < matlab.unittest.TestCase methods(Test) % Test methods + + function test_constructor(tc) + % Test the constructor + tc.verifyError(@() get_material(), 'MATLAB:notEnoughInputs', 'Not enough input arguments.'); + tc.verifyError(@() get_material("water", 2), 'MATLAB:invalidInput', 'Invalid number of input arguments.'); + tc.verifyError(@() get_material(1, 2, 3, 4), 'MATLAB:invalidInput', 'Invalid number of input arguments.'); + tc.verifyError(@() get_material(1), 'assert:failure', 'When only one input is given, it should be a string with the name of the material.'); + tc.verifyError(@() get_material('water'), 'assert:failure', 'When only one input is given, it should be a string with the name of the material.'); + tc.verifyError(@() get_material("wat"), 'MATLAB:invalidMaterial', 'The material wat is not available. Available materials are: air, blood, bone, lung, muscle, water.'); + tc.verifyError(@() get_material([1, 8], [0.111898, 0.888102], [0, 1]), 'assert:failure', 'The density should be a scalar number.'); + tc.verifyError(@() get_material([1, 8], [0.1, 0.2], "1"), 'assert:failure', 'The density should be a scalar number.'); + tc.verifyError(@() get_material([1, 2;1, 2], [0.111898, 0.888102], 1), 'assert:failure', 'The atomic numbers should be a vector.'); + tc.verifyError(@() get_material([1, 8], [0.1, 0.2; 0.3, 0.4], 1), 'assert:failure', 'The mass fractions should be a vector.'); + tc.verifyError(@() get_material(1, [0.1, 0.2], 1), 'assert:failure', 'The atomic numbers and mass fractions should have the same length.'); + end function test_water_by_z(tc) - water_mat = material("water"); + water_mat = get_material("water"); + water = @water_mat.get_mu; + tc.verifyEqual(water(1e-3 ), 4.078E+03, "RelTol", 1e-3) + tc.verifyEqual(water(1.5e-2), 1.673E+00, "RelTol", 1e-3) + tc.verifyEqual(water(8e-2 ), 1.837E-01, "RelTol", 1e-3) + tc.verifyEqual(water(1.25 ), 6.323E-02, "RelTol", 1e-3) + tc.verifyEqual(water(2e1 ), 1.813E-02, "RelTol", 1e-3) + + water_mat = get_material([1, 8], [0.111898, 0.888102], 1); water = @water_mat.get_mu; tc.verifyEqual(water(1e-3 ), 4.078E+03, "RelTol", 1e-3) tc.verifyEqual(water(1.5e-2), 1.673E+00, "RelTol", 1e-3) @@ -13,7 +36,7 @@ function test_water_by_z(tc) end function test_bone_by_z(tc) - bone_mat = material("bone"); + bone_mat = get_material("bone"); bone = @bone_mat.get_mu; tc.verifyEqual(bone(1e-3 ), 3.781E+03 * 1.92, "RelTol", 1e-3) tc.verifyEqual(bone(1.305E-03), 1.873E+03 * 1.92, "RelTol", 1e-2)% Check this is because of the absorption edge @@ -27,7 +50,7 @@ function test_bone_by_z(tc) end function test_blood_by_z(tc) - blood_mat = material("blood"); + blood_mat = get_material("blood"); blood = @blood_mat.get_mu; tc.verifyEqual(blood(1e-3 ), 3.806E+03 * 1.06, "RelTol", 1e-3) tc.verifyEqual(blood(1.5e-3 ), 1.282E+03 * 1.06, "RelTol", 1e-3) @@ -42,7 +65,7 @@ function test_blood_by_z(tc) end function test_lung_by_z(tc) - lung_mat = material("lung"); + lung_mat = get_material("lung"); lung = @lung_mat.get_mu; tc.verifyEqual(lung(1e-3 ), 3.803E+03 * 1.050, "RelTol", 1e-3) tc.verifyEqual(lung(1.5e-3 ), 1.283E+03 * 1.050, "RelTol", 1e-3) @@ -56,7 +79,7 @@ function test_lung_by_z(tc) end function test_muscle_by_z(tc) - muscle_mat = material("muscle"); + muscle_mat = get_material("muscle"); muscle = @muscle_mat.get_mu; tc.verifyEqual(muscle(1e-3 ), 3.719E+03 * 1.050, "RelTol", 1e-3) tc.verifyEqual(muscle(1.5e-3 ), 1.251E+03 * 1.050, "RelTol", 1e-3) diff --git a/tests/detector_tests.m b/tests/detector_tests.m index 237e9c2..25b752a 100644 --- a/tests/detector_tests.m +++ b/tests/detector_tests.m @@ -117,7 +117,7 @@ function test_get_scan_angles(tc) function test_generate_image(tc) detector = parallel_detector(10, [1, 1], [5, 1], 4); - mat = material("water"); mat.get_mu = @(e) 1; + mat = get_material("water"); mat.get_mu = @(e) 1; my_box = voxel_box([0;0;0], [3;3;3], mat); array = voxel_array(zeros(3, 1), [5; 5; 5], 1, my_box); sq2 = sqrt(2); @@ -137,7 +137,7 @@ function test_generate_image(tc) function test_generate_image_p(tc) detector = parallel_detector(10, [1, 1], [5, 1], 4); - my_mat = material("water"); my_mat.get_mu = @(e) 1; + my_mat = get_material("water"); my_mat.get_mu = @(e) 1; my_box = voxel_box([0;0;0], [3;3;3], my_mat); array = voxel_array(zeros(3, 1), [5; 5; 5], 1, my_box); sq2 = sqrt(2); @@ -154,7 +154,7 @@ function test_curved_generate_image_p(tc) dist_to_detector = 10; num_rotations = 16; detector = curved_detector(dist_to_detector, [pi/500, 0.01], [25, 10], num_rotations); - my_mat = material("water"); my_mat.get_mu = @(e) 1; + my_mat = get_material("water"); my_mat.get_mu = @(e) 1; my_box = voxel_box([0;0;0], [3;3;3], my_mat); array = voxel_array(zeros(3, 1), [5; 5; 5], 1, my_box); image = detector.generate_image(array); diff --git a/tests/ray_tests.m b/tests/ray_tests.m index fddf916..f412e65 100644 --- a/tests/ray_tests.m +++ b/tests/ray_tests.m @@ -84,7 +84,7 @@ function test_siddon_ray(tc) end function test_calculate_mu(tc) - water = material("water"); + water = get_material("water"); my_box = voxel_box([0;0;0], [3;3;3], water); array = voxel_array(zeros(3, 1), [5; 5; 5], 1, my_box); voxel_attenuation = water.get_mu(100/1000); % Default value of ray diff --git a/tests/voxel_shapes_tests.m b/tests/voxel_shapes_tests.m index 796e077..8c4141e 100644 --- a/tests/voxel_shapes_tests.m +++ b/tests/voxel_shapes_tests.m @@ -14,7 +14,7 @@ function test_cylinder(tc) radius = 1; width = 2; - water= material("water"); + water= get_material("water"); my_cyl = voxel_cylinder([0, 0, 0], radius, width, water); for x = -2:0.25:2 if x ^ 2 <= radius ^ 2 @@ -59,7 +59,7 @@ function test_cylinder(tc) end function box(tc) - water = material("water"); + water = get_material("water"); my_box = voxel_box([0,0,0], 100, water); res = my_box(-100:0.5:100, -100:0.5:100, -100:0.5:100, 10); exp = cat(2,zeros(1, 100), zeros(1, 201)+water.get_mu(10), zeros(1, 100)); %101 includes 0 @@ -85,11 +85,11 @@ function shepp_logan(tc) end function test_collection(tc) - mat1 = material("water"); mat1.get_mu = @(e) 10; + mat1 = get_material("water"); mat1.get_mu = @(e) 10; big_box = voxel_box([0,0,0], 10, mat1); - mat2 = material("water"); mat2.get_mu = @(e) 5; + mat2 = get_material("water"); mat2.get_mu = @(e) 5; med_box = voxel_box([0,0,0], 6, mat2); - mat3 = material("water"); mat3.get_mu = @(e) 1; + mat3 = get_material("water"); mat3.get_mu = @(e) 1; small_box = voxel_box([0,0,0], 2, mat3); my_collection = voxel_collection(big_box, med_box, small_box); for x = -5:5