From b71fddc1f4e67d75860df1716e05d762acb8f142 Mon Sep 17 00:00:00 2001 From: James Kerns Date: Thu, 28 Mar 2024 09:42:12 -0500 Subject: [PATCH 1/5] fix Low equation for couch-kick shift vector --- docs/source/changelog.rst | 3 + docs/source/winston_lutz.rst | 163 ++++++++++++++++----- pylinac/winston_lutz.py | 12 +- tests_basic/test_winstonlutz.py | 244 +++++++++++++++++++++++++++++--- 4 files changed, 364 insertions(+), 58 deletions(-) diff --git a/docs/source/changelog.rst b/docs/source/changelog.rst index c827c884..cf22f6fa 100644 --- a/docs/source/changelog.rst +++ b/docs/source/changelog.rst @@ -47,6 +47,9 @@ Winston Lutz This can be turned off by passing ``zoom=False`` to the ``plot_images`` method. * When using custom BB arrangements, use the new :class:`~pylinac.winston_lutz.BBConfig` class instead of a dictionary. See the updated :ref:``custom-bb-arrangements`` section for more. +* A bug was fixed for the BB shift vector/instructions when analyzing images with couch kicks. + The Low paper which contains the mathematical transforms appears to have missing negative signs. This + has been fixed and validated using the new image generator ability to create images with couch kicks. Image Generator ^^^^^^^^^^^^^^^ diff --git a/docs/source/winston_lutz.rst b/docs/source/winston_lutz.rst index 2a7e5af8..58eacbe3 100644 --- a/docs/source/winston_lutz.rst +++ b/docs/source/winston_lutz.rst @@ -567,6 +567,8 @@ if the algorithm doesn't work on a given image and when commissioning pylinac. I to question the accuracy of the algorithm. Since no linac is perfect and the results are sub-millimeter, discerning what is true error vs algorithmic error can be difficult. The image generator module is a perfect solution since it can remove or reproduce the former error. +See the :ref:`image_generator_module` module for more information and specifically the :func:`~pylinac.core.image_generator.utils.generate_winstonlutz` function. + Perfect Delivery ^^^^^^^^^^^^^^^^ @@ -577,6 +579,7 @@ and a field size of 4x4cm: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -587,17 +590,17 @@ and a field size of 4x4cm: wl_dir = 'wl_dir' generate_winstonlutz( - AS1200Image(), + AS1200Image(1000), FilteredFieldLayer, dir_out=wl_dir, final_layers=[GaussianFilterLayer(),], bb_size_mm=4, - field_size_mm=(40, 40), + field_size_mm=(25, 25), ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4) - wl.plot_images() + wl.plot_images(axis=Axis.GBP_COMBO) which has an output of:: @@ -606,7 +609,8 @@ which has an output of:: Number of images: 4 Maximum 2D CAX->BB distance: 0.00mm Median 2D CAX->BB distance: 0.00mm - Shift to iso: facing gantry, move BB: RIGHT 0.00mm; IN 0.00mm; UP 0.00mm + Mean 2D CAX->BB distance: 0.00mm + Shift to iso: facing gantry, move BB: RIGHT 0.00mm; OUT 0.00mm; DOWN 0.00mm Gantry 3D isocenter diameter: 0.00mm (4/4 images considered) Maximum Gantry RMS deviation (mm): 0.00mm Maximum EPID RMS deviation (mm): 0.00mm @@ -626,6 +630,7 @@ Let's now offset the BB by 1mm to the left: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -636,18 +641,18 @@ Let's now offset the BB by 1mm to the left: wl_dir = 'wl_dir' generate_winstonlutz( - AS1200Image(), + AS1200Image(1000), FilteredFieldLayer, dir_out=wl_dir, final_layers=[GaussianFilterLayer(),], bb_size_mm=4, - field_size_mm=(40, 40), + field_size_mm=(25, 25), offset_mm_left=1, ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4) - wl.plot_images() + wl.plot_images(axis=Axis.GBP_COMBO) with an output of:: @@ -656,7 +661,8 @@ with an output of:: Number of images: 4 Maximum 2D CAX->BB distance: 1.01mm Median 2D CAX->BB distance: 0.50mm - Shift to iso: facing gantry, move BB: RIGHT 1.01mm; IN 0.00mm; UP 0.00mm + Mean 2D CAX->BB distance: 0.50mm + Shift to iso: facing gantry, move BB: RIGHT 1.01mm; OUT 0.00mm; DOWN 0.00mm Gantry 3D isocenter diameter: 0.00mm (4/4 images considered) Maximum Gantry RMS deviation (mm): 1.01mm Maximum EPID RMS deviation (mm): 0.00mm @@ -671,12 +677,17 @@ We have correctly found that the max distance is 1mm and the required shift to i Gantry Tilt ^^^^^^^^^^^ +.. note:: + + The term tilt (vs sag) here represents movement in the Y-axis :ref:`coordinate space `. + We can simulate gantry tilt, where at 0 and 180 the gantry tilts forward and backward respectively. We use a realistic value of 1mm. Note that everything else is perfect: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -687,31 +698,32 @@ We can simulate gantry tilt, where at 0 and 180 the gantry tilts forward and bac wl_dir = 'wl_dir' generate_winstonlutz( - AS1200Image(), + AS1200Image(1000), FilteredFieldLayer, dir_out=wl_dir, final_layers=[GaussianFilterLayer(),], bb_size_mm=4, - field_size_mm=(40, 40), + field_size_mm=(25, 25), gantry_tilt=1, ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4) - wl.plot_images() + wl.plot_images(axis=Axis.GBP_COMBO) with output of:: Winston-Lutz Analysis ================================= Number of images: 4 - Maximum 2D CAX->BB distance: 0.90mm - Median 2D CAX->BB distance: 0.45mm - Shift to iso: facing gantry, move BB: LEFT 0.00mm; IN 0.00mm; UP 0.00mm - Gantry 3D isocenter diameter: 1.79mm (4/4 images considered) - Maximum Gantry RMS deviation (mm): 0.90mm - Maximum EPID RMS deviation (mm): 0.90mm - Gantry+Collimator 3D isocenter diameter: 1.79mm (4/4 images considered) + Maximum 2D CAX->BB distance: 1.01mm + Median 2D CAX->BB distance: 0.50mm + Mean 2D CAX->BB distance: 0.50mm + Shift to iso: facing gantry, move BB: RIGHT 0.00mm; OUT 0.00mm; DOWN 0.00mm + Gantry 3D isocenter diameter: 2.02mm (4/4 images considered) + Maximum Gantry RMS deviation (mm): 1.01mm + Maximum EPID RMS deviation (mm): 1.01mm + Gantry+Collimator 3D isocenter diameter: 2.02mm (4/4 images considered) Collimator 2D isocenter diameter: 0.00mm (1/4 images considered) Maximum Collimator RMS deviation (mm): 0.00 Couch 2D isocenter diameter: 0.00mm (1/4 images considered) @@ -720,16 +732,74 @@ with output of:: Note that since the tilt is symmetric the shift to iso is 0 despite our non-zero median distance. I.e. we are at iso, the iso just isn't perfect and we are thus at the best possible position. +Gantry Sag +^^^^^^^^^^ + +.. note:: + + The term sag (vs tilt) here represents movement in the Z-axis at 90/270 :ref:`coordinate space `. + +We can simulate gantry sag, where at 90 and 270 the gantry tilts towards the floor. We use a realistic value of +1mm. Note that everything else is perfect: + +.. plot:: + + import pylinac + from pylinac.winston_lutz import Axis + from pylinac.core.image_generator import ( + GaussianFilterLayer, + FilteredFieldLayer, + AS1200Image, + RandomNoiseLayer, + generate_winstonlutz, + ) + + wl_dir = 'wl_dir' + generate_winstonlutz( + AS1200Image(1000), + FilteredFieldLayer, + dir_out=wl_dir, + final_layers=[GaussianFilterLayer(),], + bb_size_mm=4, + field_size_mm=(25, 25), + gantry_sag=1, + ) + + wl = pylinac.WinstonLutz(wl_dir) + wl.analyze(bb_size_mm=4) + wl.plot_images(axis=Axis.GBP_COMBO) -Perfect Multi-Axis +with output of:: + + Winston-Lutz Analysis + ================================= + Number of images: 4 + Maximum 2D CAX->BB distance: 1.01mm + Median 2D CAX->BB distance: 0.50mm + Mean 2D CAX->BB distance: 0.50mm + Shift to iso: facing gantry, move BB: LEFT 0.00mm; IN 0.00mm; DOWN 1.01mm + Gantry 3D isocenter diameter: 0.00mm (4/4 images considered) + Maximum Gantry RMS deviation (mm): 1.01mm + Maximum EPID RMS deviation (mm): 1.01mm + Gantry+Collimator 3D isocenter diameter: 0.00mm (4/4 images considered) + Collimator 2D isocenter diameter: 0.00mm (1/4 images considered) + Maximum Collimator RMS deviation (mm): 0.00 + Couch 2D isocenter diameter: 0.00mm (1/4 images considered) + Maximum Couch RMS deviation (mm): 0.00 + +Sag will not realistically vary smoothly with gantry angle but for the purposes of the test it is a good approximation. +This appears as an offset of the BB in the Z-axis. + +Offset Multi-Axis ^^^^^^^^^^^^^^^^^^ We can also vary the axis data for the images produced. Below we create a typical multi-axis WL with varying gantry, collimator, and couch -(cardinal values for axis of interest with all other axes at 0): +values. We offset the BB to the left by 2mm for visualization purposes: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -740,20 +810,43 @@ We can also vary the axis data for the images produced. Below we create a typica wl_dir = 'wl_dir' generate_winstonlutz( - AS1200Image(), + AS1200Image(1000), FilteredFieldLayer, dir_out=wl_dir, - final_layers=[GaussianFilterLayer(),], + final_layers=[GaussianFilterLayer(), ], bb_size_mm=4, - field_size_mm=(40, 40), - image_axes=[(0, 0, 0), (0, 90, 0), (0, 270, 0), + field_size_mm=(20, 20), + offset_mm_left=2, + image_axes=[(0, 0, 0), (0, 45, 0), (0, 270, 0), (90, 0, 0), (180, 0, 0), (270, 0, 0), - (0, 0, 90), (0, 0, 270)] + (0, 0, 90), (0, 0, 315)] ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4) - wl.plot_images() + print(wl.results()) + wl.plot_images(axis=Axis.GBP_COMBO) + +with output of:: + + Winston-Lutz Analysis + ================================= + Number of images: 8 + Maximum 2D CAX->BB distance: 2.01mm + Median 2D CAX->BB distance: 2.01mm + Mean 2D CAX->BB distance: 1.51mm + Shift to iso: facing gantry, move BB: RIGHT 2.01mm; OUT 0.00mm; DOWN 0.00mm + Gantry 3D isocenter diameter: 0.00mm (4/8 images considered) + Maximum Gantry RMS deviation (mm): 2.01mm + Maximum EPID RMS deviation (mm): 0.00mm + Gantry+Collimator 3D isocenter diameter: 0.00mm (6/8 images considered) + Collimator 2D isocenter diameter: 0.00mm (3/8 images considered) + Maximum Collimator RMS deviation (mm): 2.01 + Couch 2D isocenter diameter: 3.72mm (3/8 images considered) + Maximum Couch RMS deviation (mm): 2.01 + +Note the shift to iso is 2mm to the right, as we offset the BB to the left by 2mm, even with our +couch kicks and collimator rotations. Perfect Cone ^^^^^^^^^^^^ @@ -763,6 +856,7 @@ We can also look at simulated cone WL images. Here we use the 17.5mm cone: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -773,17 +867,18 @@ We can also look at simulated cone WL images. Here we use the 17.5mm cone: wl_dir = 'wl_dir' generate_winstonlutz_cone( - AS1200Image(), + AS1200Image(1000), FilterFreeConeLayer, dir_out=wl_dir, - final_layers=[GaussianFilterLayer(),], + final_layers=[GaussianFilterLayer(), ], bb_size_mm=4, cone_size_mm=17.5, ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4) - wl.plot_images() + print(wl.results()) + wl.plot_images(axis=Axis.GBP_COMBO) Low-density BB ^^^^^^^^^^^^^^ @@ -793,6 +888,7 @@ Simulate a low-density BB surrounded by higher-density material: .. plot:: import pylinac + from pylinac.winston_lutz import Axis from pylinac.core.image_generator import ( GaussianFilterLayer, FilteredFieldLayer, @@ -803,19 +899,20 @@ Simulate a low-density BB surrounded by higher-density material: wl_dir = 'wl_dir' generate_winstonlutz( - AS1200Image(), + AS1200Image(1000), FilteredFieldLayer, dir_out=wl_dir, - final_layers=[GaussianFilterLayer(),], + final_layers=[GaussianFilterLayer(), ], bb_size_mm=4, - field_size_mm=(40, 40), + field_size_mm=(25, 25), field_alpha=0.6, # set the field to not max out bb_alpha=0.3 # normally this is negative to attenuate the beam, but here we increase the signal ) wl = pylinac.WinstonLutz(wl_dir) wl.analyze(bb_size_mm=4, low_density_bb=True) - wl.plot_images() + print(wl.results()) + wl.plot_images(axis=Axis.GBP_COMBO) API Documentation diff --git a/pylinac/winston_lutz.py b/pylinac/winston_lutz.py index cfc3752f..d21a1fcd 100644 --- a/pylinac/winston_lutz.py +++ b/pylinac/winston_lutz.py @@ -1227,8 +1227,12 @@ def bb_shift_vector(self) -> Vector: ) A[2 * idx : 2 * idx + 2, :] = np.array( [ - [-cos(couch), -sin(couch), 0], - [-cos(gantry) * sin(couch), cos(gantry) * cos(couch), -sin(gantry)], + [-cos(couch), sin(couch), 0], + [ + -cos(gantry) * -sin(couch), + cos(gantry) * cos(couch), + -sin(gantry), + ], ] ) # equation 6 (minus delta) epsilon[2 * idx : 2 * idx + 2] = np.array( @@ -2341,10 +2345,10 @@ def bb_projection_gantry_plane( gantry_offset = ( offset_up * -sin(gantry) + addtl_up_shift - + offset_left * -cos(gantry) * -cos(couch) + + offset_left * -cos(gantry) * cos(couch) + addtl_left_shift ) - return gantry_offset + couch_long_aspect + return gantry_offset - couch_long_aspect def _bb_projection_with_rotation( diff --git a/tests_basic/test_winstonlutz.py b/tests_basic/test_winstonlutz.py index 4ac0bd5d..560288f4 100644 --- a/tests_basic/test_winstonlutz.py +++ b/tests_basic/test_winstonlutz.py @@ -13,6 +13,12 @@ from pylinac import WinstonLutz from pylinac.core.array_utils import create_dicom_files_from_3d_array from pylinac.core.geometry import Vector, vector_is_close +from pylinac.core.image_generator import ( + AS1200Image, + GaussianFilterLayer, + PerfectFieldLayer, + generate_winstonlutz, +) from pylinac.core.io import TemporaryZipDirectory from pylinac.core.scale import MachineScale from pylinac.metrics.features import is_round @@ -331,14 +337,33 @@ def test_gantry_plane_projection_with_couch(self): bb_projection_gantry_plane( offset_up=0, offset_left=0, sad=1000, gantry=0, couch=90, offset_in=10 ), - 10, + -10, + ) + + # couch 315 should be 0.707 of the offset positively (right) + self.assertAlmostEqual( + bb_projection_gantry_plane( + offset_up=0, offset_left=0, sad=1000, gantry=0, couch=315, offset_in=10 + ), + 7.07, + places=2, + ) + + # couch 45 should be 0.707 of the offset negatively (left) + self.assertAlmostEqual( + bb_projection_gantry_plane( + offset_up=0, offset_left=0, sad=1000, gantry=0, couch=45, offset_in=10 + ), + -7.07, + places=2, ) + # couch 270 with in offset will become fully right offset self.assertAlmostEqual( bb_projection_gantry_plane( offset_up=0, offset_left=0, sad=1000, gantry=0, couch=270, offset_in=10 ), - -10, + 10, ) # couch 45 w/ left offset should be 0.707 of the offset @@ -346,7 +371,7 @@ def test_gantry_plane_projection_with_couch(self): bb_projection_gantry_plane( offset_up=0, offset_left=10, sad=1000, gantry=0, couch=45, offset_in=0 ), - 7.07, + -7.07, places=2, ) @@ -375,21 +400,21 @@ def test_gantry_plane_projection_with_couch(self): 0, ) - # couch 45 will be 0.707 of the offset in + # couch 45 will be 0.707 of the offset in (left, negative) self.assertAlmostEqual( bb_projection_gantry_plane( offset_up=0, offset_left=0, sad=1000, gantry=0, couch=45, offset_in=10 ), - 7.07, + -7.07, places=2, ) - # couch 45 and gantry 45 will be 0.707 * 0.707 of the offset in + # couch 45 and gantry 45 will be 0.707 * 0.707 of the offset in (left, negative) self.assertAlmostEqual( bb_projection_gantry_plane( offset_up=0, offset_left=0, sad=1000, gantry=45, couch=45, offset_in=10 ), - 5, + -5, places=2, ) @@ -616,10 +641,10 @@ def test_str_or_enum(self): def test_bb_shift_instructions(self): move = self.wl.bb_shift_instructions() - self.assertTrue("RIGHT" in move) + self.assertTrue("LEFT" in move) move = self.wl.bb_shift_instructions(couch_vrt=-2, couch_lat=1, couch_lng=100) - self.assertTrue("RIGHT" in move) + self.assertTrue("LEFT" in move) self.assertTrue("VRT" in move) def test_results_data(self): @@ -897,6 +922,183 @@ def test_bb_size_doesnt_change_result(self): ) +class SyntheticWLMixin(WinstonLutzMixin): + """This mixin generates WL images on the fly for testing purposes rather than pulling from the cloud""" + + zip = False + field_size = 20 + penumbra_mm = 1.5 + offset_mm_in = 0 + offset_mm_left = 0 + offset_mm_up = 0 + bb_alpha = -0.8 + gantry_sag = 0 + gantry_tilt = 0 + images_axes: tuple[tuple[int, int, int]] = ( + (0, 0, 0), + (90, 0, 0), + (180, 0, 0), + (270, 0, 0), + (0, 0, 45), + (0, 0, 90), + (0, 0, 270), + (0, 0, 315), + ) + + @classmethod + def get_filename(cls) -> str: + """We generate the files and return a local temp path""" + tmp_path = tempfile.mkdtemp() + generate_winstonlutz( + simulator=AS1200Image(1000), + field_layer=PerfectFieldLayer, + dir_out=tmp_path, + field_size_mm=(cls.field_size, cls.field_size), + final_layers=[GaussianFilterLayer(sigma_mm=cls.penumbra_mm)], + bb_alpha=cls.bb_alpha, + offset_mm_in=cls.offset_mm_in, + offset_mm_up=cls.offset_mm_up, + offset_mm_left=cls.offset_mm_left, + bb_size_mm=cls.bb_size, + gantry_sag=cls.gantry_sag, + gantry_tilt=cls.gantry_tilt, + machine_scale=MachineScale.IEC61217, + image_axes=cls.images_axes, + ) + return tmp_path + + @property + def num_images(self): + """Shortcut the num of images check since we are creating them. No need to check.""" + return len(self.images_axes) + + def test_bb_shift_vector(self): + """The vector should be opposite the set offsets""" + self.assertAlmostEqual( + self.wl.bb_shift_vector.x, self.offset_mm_left, delta=0.05 + ) # no negative; left is negative + self.assertAlmostEqual( + self.wl.bb_shift_vector.y, -self.offset_mm_in, delta=0.05 + ) + self.assertAlmostEqual( + self.wl.bb_shift_vector.z, -self.offset_mm_up, delta=0.05 + ) + + +class Synthetic1mmLeftNoCouch(SyntheticWLMixin, TestCase): + images_axes = ( + (0, 0, 0), + (90, 0, 0), + (180, 0, 0), + (270, 0, 0), + ) + offset_mm_left = 1 + cax2bb_max_distance = 1 + cax2bb_mean_distance = 0.5 + cax2bb_median_distance = 0.5 + + +class Synthetic1mmLeft(SyntheticWLMixin, TestCase): + """This should give the same result as above but including couch images""" + + offset_mm_left = 1 + cax2bb_max_distance = 1 + cax2bb_mean_distance = 0.67 + cax2bb_median_distance = 1 + couch_iso_size = 2 + + +class Synthetic1mmRight(SyntheticWLMixin, TestCase): + offset_mm_left = -1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 0.75 + cax2bb_median_distance = 1 + couch_iso_size = 2 + + +class Synthetic1mmUp(SyntheticWLMixin, TestCase): + offset_mm_up = 1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 0.25 + cax2bb_median_distance = 0 + + +class Synthetic1mmDown(SyntheticWLMixin, TestCase): + offset_mm_up = -1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 0.25 + cax2bb_median_distance = 0 + couch_iso_size = 0 + + +class Synthetic1mmIn(SyntheticWLMixin, TestCase): + offset_mm_in = 1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 1.0 + cax2bb_median_distance = 1 + couch_iso_size = 2.0 + + +class Synthetic1mmOut(SyntheticWLMixin, TestCase): + offset_mm_in = -1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 1.0 + cax2bb_median_distance = 1 + couch_iso_size = 2.0 + + +class Synthetic1mmIn1mmLeft(SyntheticWLMixin, TestCase): + offset_mm_in = 1 + offset_mm_left = 1 + cax2bb_max_distance = 1.41 + cax2bb_mean_distance = 1.3 + cax2bb_median_distance = 1.4 + couch_iso_size = 2.8 + + +class Synthetic1mmOut1mmRight(SyntheticWLMixin, TestCase): + offset_mm_in = -1 + offset_mm_left = -1 + cax2bb_max_distance = 1.41 + cax2bb_mean_distance = 1.3 + cax2bb_median_distance = 1.4 + couch_iso_size = 2.8 + + +class Synthetic2mmUp1mmLeft(SyntheticWLMixin, TestCase): + offset_mm_up = 2 + offset_mm_left = 1 + cax2bb_max_distance = 2.0 + cax2bb_mean_distance = 1.25 + cax2bb_median_distance = 1.0 + couch_iso_size = 2.0 + + +class Synthetic2mmRight1mmDown(SyntheticWLMixin, TestCase): + offset_mm_up = -1 + offset_mm_left = -2 + cax2bb_max_distance = 2.0 + cax2bb_mean_distance = 1.75 + cax2bb_median_distance = 2.0 + couch_iso_size = 4.0 + + +class Synthetic1mmOut1SidedCouch(SyntheticWLMixin, TestCase): + offset_mm_in = -1 + cax2bb_max_distance = 1.0 + cax2bb_mean_distance = 1.0 + cax2bb_median_distance = 1 + couch_iso_size = 1.42 + images_axes = ( + (0, 0, 0), + (90, 0, 0), + (180, 0, 0), + (270, 0, 0), + (0, 0, 45), + (0, 0, 90), # only shift couch one way; should still give same shift result + ) + + class WLDemo(WinstonLutzMixin, TestCase): num_images = 17 gantry_iso_size = 1 @@ -908,7 +1110,7 @@ class WLDemo(WinstonLutzMixin, TestCase): machine_scale = MachineScale.VARIAN_IEC epid_deviation = 1.3 axis_of_rotation = {0: Axis.REFERENCE} - bb_shift_vector = Vector(x=0.4, y=-0.4, z=-0.2) + bb_shift_vector = Vector(x=0, y=-0.25, z=-0.2) delete_file = False @classmethod @@ -921,12 +1123,12 @@ def new_instance(cls) -> WinstonLutz: return WinstonLutz.from_demo_images() def test_different_scale_has_different_shift(self): - assert "RIGHT" in self.wl.bb_shift_instructions() - assert self.wl.bb_shift_vector.x > 0.1 + assert "LEFT" in self.wl.bb_shift_instructions() + assert self.wl.bb_shift_vector.x < 0.0 new_wl = WinstonLutz.from_demo_images() new_wl.analyze(machine_scale=MachineScale.IEC61217) - assert new_wl.bb_shift_vector.x < 0.1 - assert "LEFT" in new_wl.bb_shift_instructions() + assert new_wl.bb_shift_vector.x > 0.1 + assert "RIGHT" in new_wl.bb_shift_instructions() def test_multiple_analyses_gives_same_result(self): original_vector = copy.copy(self.wl.bb_shift_vector) @@ -1098,7 +1300,7 @@ class KatyiX0(WinstonLutzMixin, TestCase): cax2bb_median_distance = 0.8 cax2bb_mean_distance = 0.7 machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=-0.5, y=0.4, z=-0.5) + bb_shift_vector = Vector(x=-0.4, y=0.15, z=-0.5) print_results = True @@ -1124,7 +1326,7 @@ class KatyiX2(WinstonLutzMixin, TestCase): cax2bb_median_distance = 0.5 cax2bb_mean_distance = 0.6 machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=0.4, y=-0.1, z=0.1) + bb_shift_vector = Vector(x=0.15, y=-0.15, z=0.1) class KatyiX3(WinstonLutzMixin, TestCase): @@ -1137,7 +1339,7 @@ class KatyiX3(WinstonLutzMixin, TestCase): cax2bb_median_distance = 0.8 cax2bb_mean_distance = 0.75 machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=-0.3, y=0.4, z=-0.5) + bb_shift_vector = Vector(x=-0.1, y=0.2, z=-0.5) class KatyTB0(WinstonLutzMixin, TestCase): @@ -1151,7 +1353,7 @@ class KatyTB0(WinstonLutzMixin, TestCase): cax2bb_mean_distance = 0.8 machine_scale = MachineScale.VARIAN_IEC axis_of_rotation = {-1: Axis.REFERENCE} - bb_shift_vector = Vector(x=-0.7, y=-0.1, z=-0.2) + bb_shift_vector = Vector(x=-0.4, y=-0.1, z=-0.25) class KatyTB1(WinstonLutzMixin, TestCase): @@ -1165,7 +1367,7 @@ class KatyTB1(WinstonLutzMixin, TestCase): cax2bb_mean_distance = 0.6 axis_of_rotation = {0: Axis.GANTRY} machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=-0.6, y=-0.2) + bb_shift_vector = Vector(x=-0.3, y=-0.2) class KatyTB2(WinstonLutzMixin, TestCase): @@ -1229,7 +1431,7 @@ class SugarLandiX2015(WinstonLutzMixin, TestCase): cax2bb_median_distance = 1.05 cax2bb_mean_distance = 1.1 machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=0.4, y=-0.7, z=0.1) + bb_shift_vector = Vector(x=0.6, y=-0.5, z=0.1) class BayAreaiX0(WinstonLutzMixin, TestCase): @@ -1243,7 +1445,7 @@ class BayAreaiX0(WinstonLutzMixin, TestCase): cax2bb_median_distance = 0.6 cax2bb_mean_distance = 0.6 machine_scale = MachineScale.VARIAN_IEC - bb_shift_vector = Vector(x=0.3, y=-0.4, z=-0.2) + bb_shift_vector = Vector(x=0, y=-0.3, z=-0.2) class DAmoursElektaOffset(WinstonLutzMixin, TestCase): From 2730c4df6161e25420fc632a1a442b7d32c8c8b8 Mon Sep 17 00:00:00 2001 From: James Kerns Date: Thu, 28 Mar 2024 15:46:55 -0500 Subject: [PATCH 2/5] write up docs on couch shift. --- docs/source/changelog.rst | 13 ++++++---- docs/source/winston_lutz.rst | 49 ++++++++++++++++++++++++++++++++++-- pylinac/winston_lutz.py | 2 +- 3 files changed, 56 insertions(+), 8 deletions(-) diff --git a/docs/source/changelog.rst b/docs/source/changelog.rst index f863957e..b4dcaabf 100644 --- a/docs/source/changelog.rst +++ b/docs/source/changelog.rst @@ -10,7 +10,7 @@ Image Metrics * The ``GlobalSizedDiskLocator`` class has added an ``invert`` parameter. This parameter existed for the other locators, but was missing for the global disk locator. Previously, the locator was always inverting the image (assuming images like EPID). Now, the parameter can be used to control this behavior. By - default, the paramter is true for backwards-compatibility. + default, the parameter is true for backwards-compatibility. Image ^^^^^ @@ -35,18 +35,19 @@ Winston Lutz cases of N targets and M fields. For multi-target/multi-field analyses, the algorithm was very memory-intensive because it was creating X*Y analysis objects where X is the number of images and Y is the number of targets. Memory usage has been reduced from this refactor. -* The class ``WinstonLutz2DMultTarget`` has changed to ``WinstonLutzMultiTargetMultiFieldImage``. +* The class ``WinstonLutz2DMultiTarget`` has changed to :class:`~pylinac.winston_lutz.WinstonLutzMultiTargetMultiFieldImage`. Unless you are using the class directly, this change should not affect you. * The :meth:`~pylinac.winston_lutz.WinstonLutzMultiTargetMultiField.plot_images` method has changed. - Instead of returning N figures where each figure is a set of plots is for a single BB, M figures are returned where + Instead of returning N figures where N is the number of BBs where each figure is a set of plots for each BB, M figures are returned where M is the number of images. Each plot will show the image and all detected BBs and fields. This gives - better context about which BB was dected where as it relates to the image as a whole. + better context about which BB was detected where as it relates to the image as a whole. Images within PDFs will also be generated in the same way. * For MultiField analyses, the ``cax2bb_distance()`` and ``cax2epid_distance()`` metrics were giving artificially high values when the metric was ``median`` or ``mean``. This was because the metric was first calculating the maximum distance for a given image, and then taking the median or mean of those values. This was not the intended behavior. The metric now calculates the median or mean of all the distances for all BBs together. I.e. it was doing ``median(max(a1, a2, a3), max(b2, ...), ...)`` instead of ``median(a1, a2, b1, b2, ...)``. + This will result in lower values for the metric compared to previously. * Plots now show a legend of the EPID, BB, and field CAX. The legend can be turned off by passing ``legend=False`` to the ``plot_images`` method. * Plots are now zoomed to fit all the BBs/fields detected. In the simple case of a single BB at isocenter, this hasn't changed. For multi-target/multi-field WL, the plots will now be zoomed to fit all the detected BBs and fields. @@ -54,8 +55,10 @@ Winston Lutz * When using custom BB arrangements, use the new :class:`~pylinac.winston_lutz.BBConfig` class instead of a dictionary. See the updated :ref:``custom-bb-arrangements`` section for more. * A bug was fixed for the BB shift vector/instructions when analyzing images with couch kicks. - The Low paper which contains the mathematical transforms appears to have missing negative signs. This + The Low paper which contains the mathematical transforms appears to have incorrect signs in equation 6. This has been fixed and validated using the new image generator ability to create images with couch kicks. + The bug was causing the BB shift vector to be incorrect when analyzing images with couch kicks. The shift errors + were always in the LAT/LONG plane and for the most part underestimated the shift that would be needed. Image Generator ^^^^^^^^^^^^^^^ diff --git a/docs/source/winston_lutz.rst b/docs/source/winston_lutz.rst index 87781e2c..2a09617d 100644 --- a/docs/source/winston_lutz.rst +++ b/docs/source/winston_lutz.rst @@ -496,8 +496,13 @@ of the field CAX points. They also found that the gantry isocenter could by foun the field CAX as a line in 3D coordinate space, with the BB being the reference point. This method is used to find the gantry isocenter size. -Low determined the geometric transformations to apply to 2D planar images to calculate the shift to apply to the BB. -This method is used to determine the shift instructions. Specifically, equations 6 and 9. +Couch shift +^^^^^^^^^^^ + +`Low et al`_ determined the geometric transformations to apply to 2D planar images to calculate the shift to apply to the BB. +This method is used to determine the shift instructions. Specifically, equations 6, 7, and 9. + + .. note:: @@ -509,6 +514,46 @@ This method is used to determine the shift instructions. Specifically, equations To use a different scale use the ``machine_scale`` parameter, shown here :ref:`passing-a-coordinate-system`. Also see :ref:`scale`. +Implementation for the couch shift is as follows: + +For each image we determine (equation 6a): + +.. math:: + + \textbf{A}(\phi, \theta) = \begin{pmatrix} -\cos(\phi) & \sin(\phi) & 0 \\ \cos(\theta)\sin(\phi) & \cos(\theta)\cos(\phi) & -\sin(\theta) \end{pmatrix} + +where :math:`\theta` is the gantry angle and :math:`\phi` is the couch angle. + +.. warning:: + + The Low paper appears to have incorrect signs for some of the matrix. Using synthetic images with known shifts + can prove the correctness of the algorithm. + +The :math:`\xi` matrix is then calculated (equation 7): + +.. math:: + + \xi = (x_{1},y_{1}, ..., x_{i},y_{i},..., x_{n},y_{n})^{T} + +where :math:`x_{i}` and :math:`y_{i}` are the x and y scalar shifts from the field CAX to the BB for :math:`n` images. + +From equation 9 we can calculate :math:`\textbf{B}(\phi_{1},\theta_{1},..., \phi_{n},\theta_{n})`: + +.. math:: + + \textbf{B} = \begin{pmatrix} \textbf{A}(\phi_{1},\theta_{1}) \\ \vdots \\ \textbf{A}(\phi_{i},\theta_{i}) \\ \vdots \\ \textbf{A}(\phi_{n},\theta_{n}) \end{pmatrix} + +using the above definitions. + +We can then solve for the shift vector :math:`\boldsymbol{\Delta}`: + +.. math:: + + \boldsymbol{\Delta} = \textbf{B}(\phi_{1},\theta_{1},..., \phi_{n},\theta_{n}) \cdot \xi + +Implementation +^^^^^^^^^^^^^^ + The algorithm works like such: **Allowances** diff --git a/pylinac/winston_lutz.py b/pylinac/winston_lutz.py index 60c77f87..7bf209cd 100644 --- a/pylinac/winston_lutz.py +++ b/pylinac/winston_lutz.py @@ -1252,7 +1252,7 @@ def bb_shift_vector(self) -> Vector: [ [-cos(couch), sin(couch), 0], [ - -cos(gantry) * -sin(couch), + cos(gantry) * sin(couch), cos(gantry) * cos(couch), -sin(gantry), ], From 2036f742fe2ec324370eb3b03417db40afa5335c Mon Sep 17 00:00:00 2001 From: James Kerns Date: Thu, 28 Mar 2024 16:48:34 -0500 Subject: [PATCH 3/5] cache and clean --- tests_basic/test_winstonlutz.py | 48 ++++++++++++++++++++------------- 1 file changed, 29 insertions(+), 19 deletions(-) diff --git a/tests_basic/test_winstonlutz.py b/tests_basic/test_winstonlutz.py index 88c221b4..9dc4ce93 100644 --- a/tests_basic/test_winstonlutz.py +++ b/tests_basic/test_winstonlutz.py @@ -3,6 +3,7 @@ import copy import io import math +import shutil import tempfile from pathlib import Path from unittest import TestCase @@ -929,6 +930,7 @@ def test_bb_size_doesnt_change_result(self): class SyntheticWLMixin(WinstonLutzMixin): """This mixin generates WL images on the fly for testing purposes rather than pulling from the cloud""" + tmp_path: str = "" zip = False field_size = 20 penumbra_mm = 1.5 @@ -949,27 +951,35 @@ class SyntheticWLMixin(WinstonLutzMixin): (0, 0, 315), ) + @classmethod + def tearDownClass(cls): + # clean up the folder we created; + # in BB space can be at a premium. + shutil.rmtree(cls.tmp_path, ignore_errors=False) + @classmethod def get_filename(cls) -> str: - """We generate the files and return a local temp path""" - tmp_path = tempfile.mkdtemp() - generate_winstonlutz( - simulator=AS1200Image(1000), - field_layer=PerfectFieldLayer, - dir_out=tmp_path, - field_size_mm=(cls.field_size, cls.field_size), - final_layers=[GaussianFilterLayer(sigma_mm=cls.penumbra_mm)], - bb_alpha=cls.bb_alpha, - offset_mm_in=cls.offset_mm_in, - offset_mm_up=cls.offset_mm_up, - offset_mm_left=cls.offset_mm_left, - bb_size_mm=cls.bb_size, - gantry_sag=cls.gantry_sag, - gantry_tilt=cls.gantry_tilt, - machine_scale=MachineScale.IEC61217, - image_axes=cls.images_axes, - ) - return tmp_path + """We generate the files and return a local temp path. + This may get called multiple times so we do a poor-man's caching""" + if not cls.tmp_path: + cls.tmp_path = tempfile.mkdtemp() + generate_winstonlutz( + simulator=AS1200Image(1000), + field_layer=PerfectFieldLayer, + dir_out=cls.tmp_path, + field_size_mm=(cls.field_size, cls.field_size), + final_layers=[GaussianFilterLayer(sigma_mm=cls.penumbra_mm)], + bb_alpha=cls.bb_alpha, + offset_mm_in=cls.offset_mm_in, + offset_mm_up=cls.offset_mm_up, + offset_mm_left=cls.offset_mm_left, + bb_size_mm=cls.bb_size, + gantry_sag=cls.gantry_sag, + gantry_tilt=cls.gantry_tilt, + machine_scale=MachineScale.IEC61217, + image_axes=cls.images_axes, + ) + return cls.tmp_path @property def num_images(self): From 04456665467ff3d274d2c3a366949e0077b07a8e Mon Sep 17 00:00:00 2001 From: James Kerns Date: Fri, 29 Mar 2024 09:03:28 -0500 Subject: [PATCH 4/5] build bug --- docs/source/winston_lutz.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/winston_lutz.rst b/docs/source/winston_lutz.rst index 2a09617d..4d069185 100644 --- a/docs/source/winston_lutz.rst +++ b/docs/source/winston_lutz.rst @@ -799,7 +799,7 @@ We can simulate gantry sag, where at 90 and 270 the gantry tilts towards the flo generate_winstonlutz, ) - wl_dir = 'wl_dir' + wl_dir = 'wldir' generate_winstonlutz( AS1200Image(1000), FilteredFieldLayer, From b39b820a3efd5efb8d20531c1960e9087c6540c2 Mon Sep 17 00:00:00 2001 From: James Kerns Date: Mon, 1 Apr 2024 13:25:13 -0500 Subject: [PATCH 5/5] add comment --- pylinac/winston_lutz.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/pylinac/winston_lutz.py b/pylinac/winston_lutz.py index 7bf209cd..00694188 100644 --- a/pylinac/winston_lutz.py +++ b/pylinac/winston_lutz.py @@ -1249,6 +1249,8 @@ def bb_shift_vector(self) -> Vector: rotation=img.couch_angle, ) A[2 * idx : 2 * idx + 2, :] = np.array( + # note the signs are different than the paper; based on + # synthetic data we can prove this. See docs [ [-cos(couch), sin(couch), 0], [