diff --git a/auto_examples/DICOMExample/create3DSeqFromDicom.ipynb b/auto_examples/DICOMExample/create3DSeqFromDicom.ipynb new file mode 100644 index 0000000..cacb09d --- /dev/null +++ b/auto_examples/DICOMExample/create3DSeqFromDicom.ipynb @@ -0,0 +1,158 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Create 3D Sequence from DICOM\nauthor: OpenTPS team\n\nThis example shows how to read data from a 4DCT folder, create a dynamic 3D sequence with the 4DCT data, and save this sequence in serialized format on the drive.\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nfrom pathlib import Path\nimport sys" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dataLoader import readData\nfrom opentps.core.data.dynamicData._dynamic3DSequence import Dynamic3DSequence\nfrom opentps.core.io.serializedObjectIO import saveSerializedObjects" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Get the current working directory, its parent, then add the testData folder at the end of it\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "testDataPath = os.path.join(Path(os.getcwd()).parent.absolute(), 'testData/')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# read a serialized dynamic sequence\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dataPath = 'Path_to_your_4DCT_data/' # replace with the path to your 4DCT data folder\n\nprint('Datas present in ' + dataPath + 'are loaded.')\ndataList = readData(dataPath)\nprint(len(dataList), 'images found in the folder')\nprint('Image type =', type(dataList[0]))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "create a Dynamic3DSequence and change its name\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dynseq = Dynamic3DSequence(dyn3DImageList=dataList)\nprint('Type of the created object =', type(dynseq))\nprint('Sequence name =', dynseq.name)\ndynseq.name = 'Light4DCT'\nprint('Sequence name = ', dynseq.name)\nprint('Sequence lenght =', len(dynseq.dyn3DImageList))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "save it as a serialized object\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "savingPath = testDataPath + 'lightDynSeq'\nsaveSerializedObjects(dynseq, savingPath)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/create3DSeqFromImages.ipynb b/auto_examples/DICOMExample/create3DSeqFromImages.ipynb new file mode 100644 index 0000000..8126ebd --- /dev/null +++ b/auto_examples/DICOMExample/create3DSeqFromImages.ipynb @@ -0,0 +1,158 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Create 3D Sequence from Images\nauthor: OpenTPS team\n\nThis example shows how to read data from a 4DCT folder, create a dynamic 3D sequence with the 4DCT data, and save this sequence in serialized format on the drive.\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nfrom pathlib import Path\nimport sys" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dataLoader import readData\nfrom opentps.core.data.dynamicData._dynamic3DSequence import Dynamic3DSequence\nfrom opentps.core.io.serializedObjectIO import saveSerializedObjects" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Get the current working directory, its parent, then add the testData folder at the end of it\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "testDataPath = os.path.join(Path(os.getcwd()).parent.absolute(), 'testData/')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "read a serialized dynamic sequence\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dataPath = 'Path_to_your_4DCT_data/' # replace with the path to your 4DCT data folder\n\nprint('Datas present in ' + dataPath + 'are loaded.')\ndataList = readData(dataPath)\nprint(len(dataList), 'images found in the folder')\nprint('Image type =', type(dataList[0]))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "create a Dynamic3DSequence and change its name\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dynseq = Dynamic3DSequence(dyn3DImageList=dataList)\nprint('Type of the created object =', type(dynseq))\nprint('Sequence name =', dynseq.name)\ndynseq.name = 'Light4DCT'\nprint('Sequence name = ', dynseq.name)\nprint('Sequence lenght =', len(dynseq.dyn3DImageList))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "save it as a serialized object\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "savingPath = testDataPath + 'lightDynSeq'\nsaveSerializedObjects(dynseq, savingPath)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/createDynamic3DModelFromDicomFields.ipynb b/auto_examples/DICOMExample/createDynamic3DModelFromDicomFields.ipynb new file mode 100644 index 0000000..9aa5e55 --- /dev/null +++ b/auto_examples/DICOMExample/createDynamic3DModelFromDicomFields.ipynb @@ -0,0 +1,158 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Create Dynamic 3D Model from DICOM Fields\nauthor: OpenTPS team\n\nThis example shows how to read a DICOM CT and deformation fields, create a dynamic 3D model with the mid-position CT and the deformation fields, and print the model information.\n\nrunning time: ~ 5 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport sys" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dataLoader import readData\nfrom opentps.core.data.images._deformation3D import Deformation3D\nfrom opentps.core.data.dynamicData._dynamic3DModel import Dynamic3DModel" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Load DICOM CT\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "inputPaths = 'Path_to_your_CT_data/' # replace with the path to your CT data folder\ndataList = readData(inputPaths, maxDepth=0)\nmidP = dataList[0]\nprint(type(midP))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Load DICOM Deformation Fields\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "inputPaths = 'Path_to_your_deformation_fields/' # replace with the path to your deformation fields folder\ndefList = readData(inputPaths, maxDepth=0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Transform VectorField3D to deformation3D\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "deformationList = []\nfor df in defList:\n df2 = Deformation3D()\n df2.initFromVelocityField(df)\n deformationList.append(df2)\ndel defList\nprint(deformationList)\n\npatient_name = 'OpenTPS_Patient'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Create Dynamic 3D Model\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "model3D = Dynamic3DModel(name=patient_name, midp=midP, deformationList=deformationList)\nprint(model3D)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/createModelWithROIsFromDicomImages.ipynb b/auto_examples/DICOMExample/createModelWithROIsFromDicomImages.ipynb new file mode 100644 index 0000000..9c1c990 --- /dev/null +++ b/auto_examples/DICOMExample/createModelWithROIsFromDicomImages.ipynb @@ -0,0 +1,230 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Create Model with ROIs from DICOM Images\nauthor: OpenTPS team\n\nThis example shows how to: \n- read dicom data from a 4DCT folder\n- create a dynamic 3D sequence with the 4DCT data\n- read an rtStruct dicom file\n- create a dynamic 3D model and compute the midP image with the dynamic 3D sequence\n- create a patient, give him the model and rtStruct and save it as serialized data \n!!! does not work with public data for now since there is no struct in the public data !!!\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport sys\nimport time\nimport numpy as np" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from pydicom.uid import generate_uid\nfrom opentps.core.io.dataLoader import readData\nfrom opentps.core.data.dynamicData._dynamic3DSequence import Dynamic3DSequence\nfrom opentps.core.io.serializedObjectIO import saveSerializedObjects\nfrom opentps.core.data.dynamicData._dynamic3DModel import Dynamic3DModel\nfrom opentps.core.data._patient import Patient" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "chose the patient folder, which will be used as the patient name\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patientName = 'Patient_0'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "chose the 4DCT data folder\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "data4DPath = 'Path_to_your_4DCT_data/' # replace with the path to your 4DCT data folder\n# chose the dicom rtStruct file\ndataStructPath = 'Path_to_your_rtStruct_data/' # replace with the path to your rtStruct data folder\n# chose a path to save the results\nsavingPath = 'Path_to_your_saving_path/' # replace with the path where you want to save the results" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "load the 4DCT data\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "data4DList = readData(data4DPath)\nprint(len(data4DList), 'images found in the folder')\nprint('Image type =', type(data4DList[0]))\nprint('Image 0 shape =', data4DList[0].gridSize)\n\n# create a Dynamic3DSequence and change its name\ndynSeq = Dynamic3DSequence(dyn3DImageList=data4DList)\ndynSeq.name = '4DCT'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "load the rtStruct data and print its content\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "structData = readData(dataStructPath)[0]\nprint('Available ROIs')\nstructData.print_ROINames()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "create Dynamic3DModel\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "model3D = Dynamic3DModel()\n\n# change its name\nmodel3D.name = 'MidP'\n\n# give it an seriesInstanceUID\nmodel3D.seriesInstanceUID = generate_uid()\n\n# give it an seriesInstanceUID\nmodel3D.seriesInstanceUID = generate_uid()\n\n# generate the midP image and deformation fields from the dynamic 3D sequence\nstartTime = time.time()\nmodel3D.computeMidPositionImage(dynSeq, tryGPU=True)\nstopTime = time.time()\n\nprint(model3D.midp.name)\nprint('MidP computed in ', np.round(stopTime-startTime))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Create a patient and give it the patient name\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = patientName" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Add the model and rtStruct to the patient\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient.appendPatientData(model3D)\npatient.appendPatientData(structData)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Save it as a serialized object\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "saveSerializedObjects(patient, savingPath)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/crop3DDataAroundROI.ipynb b/auto_examples/DICOMExample/crop3DDataAroundROI.ipynb new file mode 100644 index 0000000..9174893 --- /dev/null +++ b/auto_examples/DICOMExample/crop3DDataAroundROI.ipynb @@ -0,0 +1,158 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Crop 3D Data Around ROI\nauthor: OpenTPS team\n\nThis example shows how to: \n- Read a serialized patient with a Dynamic3DSequence, a Dynamic3DModel and an RTStruct\n!! The data is not given in the test data folder of the project !!\n- Select an ROI from the RTStruct object\n- Get the ROI as an ROIMask\n- Get the box around the ROI in scanner coordinates\n- Crop the dynamic sequence and the dynamic model around the box\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport sys" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.processing.imageProcessing.resampler3D import crop3DDataAroundBox\nfrom opentps.core.processing.segmentation.segmentation3D import getBoxAroundROI\nfrom opentps.core.io.serializedObjectIO import loadDataStructure" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Load the serialized patient data\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dataPath = 'Path_to_your_serialized_patient_data/' # replace with the path to your serialized patient data\npatient = loadDataStructure(dataPath)[0]\n\ndynSeq = patient.getPatientDataOfType(\"Dynamic3DSequence\")[0]\ndynMod = patient.getPatientDataOfType(\"Dynamic3DModel\")[0]\nrtStruct = patient.getPatientDataOfType(\"RTStruct\")[0]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "get the ROI and mask on which we want to apply the motion signal\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "print('Available ROIs')\nrtStruct.print_ROINames()\nbodyContour = rtStruct.getContourByName('body')\nROIMask = bodyContour.getBinaryMask(origin=dynMod.midp.origin, gridSize=dynMod.midp.gridSize, spacing=dynMod.midp.spacing)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "get the box around the ROI\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "box = getBoxAroundROI(ROIMask)\nmarginInMM = [10, 10, 10]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "crop the dynamic sequence and the dynamic model around the box\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "crop3DDataAroundBox(dynSeq, box, marginInMM=marginInMM)\nprint('-'*50)\ncrop3DDataAroundBox(dynMod, box, marginInMM=marginInMM)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/exampleImageResampling.ipynb b/auto_examples/DICOMExample/exampleImageResampling.ipynb new file mode 100644 index 0000000..70c593f --- /dev/null +++ b/auto_examples/DICOMExample/exampleImageResampling.ipynb @@ -0,0 +1,302 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Image Resampling \nauthor: OpenTPS team\n\nThis example shows how the resampling function can be used on image3D objects.\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nimport matplotlib.pyplot as plt\nimport os\nfrom pathlib import Path" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dataLoader import readData\nfrom opentps.core.processing.imageProcessing.resampler3D import resample, resampleOnImage3D" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Get the current working directory, its parent, then add the testData folder at the end of it\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "testDataPath = os.path.join(Path(os.getcwd()).parent.absolute(), 'testData/')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "load an image to use as example\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dataPath = \"Path_to_your_image_data/\" # replace with the path to your image data folder\nimg = readData(dataPath)[0]\nprint('Image type =', type(img))\nzSlice = int(img.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "let's resample the image using a specific spacing (upsampling or downsampling)\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "spacingResampledDown = resample(img, img.spacing * 1.5)\nspacingResampledDownZSlice = int(spacingResampledDown.gridSize[2] / 2)\n\"\"\"\nNote that as the spacing is the second argument in the resample function, it can be use without specifying the argument name if put in the second position (here above)\nIf you prefer to be sure to specify the correct argument, you can use the name as the example here under\n\"\"\"\nspacingResampledUp = resample(img, spacing=img.spacing * 0.5)\nspacingResampledUpZSlice = int(spacingResampledUp.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "display results with spacing, gridSize and origin\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 3)\n\nax[0].imshow(img.imageArray[:, :, zSlice])\nax[0].set_title(f\"Original image\")\nax[0].set_xlabel(f\"Spacing, {img.spacing} \\n Grid Size {img.gridSize} \\n Origin {img.origin}\")\n\nax[1].imshow(spacingResampledDown.imageArray[:, :, spacingResampledDownZSlice])\nax[1].set_title(f\"Downsampled using spacing\")\nax[1].set_xlabel(f\"Spacing, {spacingResampledDown.spacing} \\n Grid Size {spacingResampledDown.gridSize} \\n Origin {spacingResampledDown.origin}\")\n\nax[2].imshow(spacingResampledUp.imageArray[:, :, spacingResampledUpZSlice])\nax[2].set_title(f\"Upsampled using spacing\")\nax[2].set_xlabel(f\"Spacing, {spacingResampledUp.spacing} \\n Grid Size {spacingResampledUp.gridSize} \\n Origin {spacingResampledUp.origin}\")\n\nplt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "now let's resample the image using a specific gridSize\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "gridSizeResampled = resample(img, gridSize=(200, 200, 200))\ngridSizeResampledZSlice = int(gridSizeResampled.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "and using BOTH a specific spacing AND gridSize\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "\"\"\"\nNote that this is not recomanded as it can push parts of the image outside the new array grid and have the same effect as cropping the data\n\"\"\"\nspacingAndGSResampled = resample(img, gridSize=(100, 100, 100), spacing=(2, 2, 2))\nspacingAndGSResampledZSlice = int(spacingAndGSResampled.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "display results with spacing, gridSize and origin\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 3)\n\nax[0].imshow(img.imageArray[:, :, zSlice])\nax[0].set_title(f\"Original image\")\nax[0].set_xlabel(f\"Spacing, {img.spacing} \\n Grid Size {img.gridSize} \\n Origin {img.origin}\")\n\nax[1].imshow(gridSizeResampled.imageArray[:, :, gridSizeResampledZSlice])\nax[1].set_title(f\"Resampled using gridSize\")\nax[1].set_xlabel(f\"Spacing, {gridSizeResampled.spacing} \\n Grid Size {gridSizeResampled.gridSize} \\n Origin {gridSizeResampled.origin}\")\n\nax[2].imshow(spacingAndGSResampled.imageArray[:, :, spacingAndGSResampledZSlice])\nax[2].set_title(f\"Resampled using gridSize AND spacing\")\nax[2].set_xlabel(f\"Spacing, {spacingAndGSResampled.spacing} \\n Grid Size {spacingAndGSResampled.gridSize} \\n Origin {spacingAndGSResampled.origin}\")\n\nplt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "now let's try using the origin, which corresponds to a translation of the image\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "originResampled = resample(img, origin=(-220, -200, -200))\noriginResampledZSlice = int(originResampled.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "and using BOTH the origin AND gridSize\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "\"\"\"\nNote that this is not recomanded as it can push parts of the image outside the new array grid and have the same effect as cropping the data\n\"\"\"\noriginAndGSResampled = resample(img, gridSize=(50, 50, 50), origin=(-220, -200, -200))\noriginAndGSResampledZSlice = int(originAndGSResampled.gridSize[2] / 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "display results with spacing, gridSize and origin\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 3)\n\nax[0].imshow(img.imageArray[:, :, zSlice])\nax[0].set_title(f\"Original image\")\nax[0].set_xlabel(f\"Spacing, {img.spacing} \\n Grid Size {img.gridSize} \\n Origin {img.origin}\")\n\nax[1].imshow(originResampled.imageArray[:, :, originResampledZSlice])\nax[1].set_title(f\"Resampled using origin\")\nax[1].set_xlabel(f\"Spacing, {originResampled.spacing} \\n Grid Size {originResampled.gridSize} \\n Origin {originResampled.origin}\")\n\nax[2].imshow(originAndGSResampled.imageArray[:, :, originAndGSResampledZSlice])\nax[2].set_title(f\"Resampled using gridSize AND origin\")\nax[2].set_xlabel(f\"Spacing, {originAndGSResampled.spacing} \\n Grid Size {originAndGSResampled.gridSize} \\n Origin {originAndGSResampled.origin}\")\n\nplt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now you can also use the following function if you need to resample an image on the grid of another image\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "resampledOnGrid = resampleOnImage3D(gridSizeResampled, fixedImage=img)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "where the spacing, gridSize and origin of the fixedImage is used to resample the data (first argument)\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "resampledOnGridZSlice = int(resampledOnGrid.gridSize[2] / 2)\n\nfig, ax = plt.subplots(1, 2)\n\nax[0].imshow(gridSizeResampled.imageArray[:, :, gridSizeResampledZSlice])\nax[0].set_title(f\"Before resampling\")\nax[0].set_xlabel(f\"Spacing, {gridSizeResampled.spacing} \\n Grid Size {gridSizeResampled.gridSize} \\n Origin {gridSizeResampled.origin}\")\n\nax[1].imshow(resampledOnGrid.imageArray[:, :, resampledOnGridZSlice])\nax[1].set_title(f\"After resampling\")\nax[1].set_xlabel(f\"Spacing, {resampledOnGrid.spacing} \\n Grid Size {resampledOnGrid.gridSize} \\n Origin {resampledOnGrid.origin}\")\n\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/generateDRRAndGTVMasks.ipynb b/auto_examples/DICOMExample/generateDRRAndGTVMasks.ipynb new file mode 100644 index 0000000..c52c244 --- /dev/null +++ b/auto_examples/DICOMExample/generateDRRAndGTVMasks.ipynb @@ -0,0 +1,302 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Generate DRR and GTV Masks\nauthor: OpenTPS team\n\nThis example shows how to:\n- read model + ROI data from a serialized file\n- create a breathing signal using the motion amplitude present in the model\n- chose an ROI to apply the breathing signal to its center of mass\n-\n\n!!! does not work with public data for now since there is no struct in the public data !!!\n\nrunning time: ~ 10 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\nfrom scipy.ndimage import zoom\nimport math\nimport time\nimport concurrent\nfrom itertools import repeat\nimport os\nimport sys\nimport numpy as np" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.processing.imageProcessing.resampler3D import crop3DDataAroundBox\nfrom opentps.core.io.serializedObjectIO import saveSerializedObjects, loadDataStructure\nfrom opentps.core.data.dynamicData._breathingSignals import SyntheticBreathingSignal\nfrom opentps.core.processing.deformableDataAugmentationToolBox.generateDynamicSequencesFromModel import generateDeformationListFromBreathingSignalsAndModel\nfrom opentps.core.processing.imageSimulation.DRRToolBox import forwardProjection\nfrom opentps.core.processing.imageProcessing.image2DManip import getBinaryMaskFromROIDRR, get2DMaskCenterOfMass\nfrom opentps.core.processing.imageProcessing.crop2D import getBoxAroundROI\nfrom opentps.core.processing.imageProcessing.imageTransform3D import getVoxelIndexFromPosition\nfrom opentps.core.processing.deformableDataAugmentationToolBox.modelManipFunctions import getAverageModelValuesAroundPosition" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "paths selection\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patientFolder = 'Patient_0'\npatientFolderComplement = ''\norgan = 'liver'\nbasePath = 'D:/ImageData/'\n\ndataSetFolder = '/test/'\ndataSetDataFolder = 'data/'\n\ndataPath = \"Path_to_your_serialized_patient_data/\" # replace with the path to your serialized patient data\nsavingPath = \"Path_to_your_saving_path/\" # replace with the path where you want to save the results\nif not os.path.exists(savingPath):\n os.umask(0)\n os.makedirs(savingPath) # Create a new directory because it does not exist\n os.makedirs(savingPath + dataSetDataFolder) # Create a new directory because it does not exist\n print(\"New directory created to save the data: \", savingPath)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "parameters selection\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "sequenceDurationInSecs = 10\nsamplingFrequency = 4\nsubSequenceSize = 50\noutputSize = [64, 64]\nbodyContourToUse = 'body'\notherContourToUse = 'MidP CT GTV'\nmarginInMM = [50, 0, 100]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "breathing signal parameters\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "amplitude = 'model'\nvariationAmplitude = 2\nbreathingPeriod = 4\nvariationFrequency = 0.1\nshift = 2\nmeanNoise = 0\nvarianceNoise = 0.5\nsamplingPeriod = 1 / samplingFrequency\nsimulationTime = sequenceDurationInSecs\nmeanEvent = 2 / 30\n# use Z - 0 for Coronal and Z - 90 for sagittal\nprojAngle = 0\nprojAxis = 'Z'\n\nmultiprocessing = False\nmaxMultiProcUse = 4\ntryGPU = True" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This function is specific to this example and used to :\n- deform a CTImage and an ROIMask,\n- create DRR's for both,\n- binarize the DRR of the ROIMask\n- compute its center of mass\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "def deformImageAndMaskAndComputeDRRs(img, ROIMask, deformation, projectionAngle=0, projectionAxis='Z', tryGPU=True, outputSize=[]):\n print('Start deformations and projections for deformation', deformation.name)\n image = deformation.deformImage(img, fillValue='closest', outputType=np.int16, tryGPU=tryGPU)\n # print(image.imageArray.shape, np.min(image.imageArray), np.max(image.imageArray), np.mean(image.imageArray))\n mask = deformation.deformImage(ROIMask, fillValue='closest', outputType=np.int16, tryGPU=tryGPU)\n centerOfMass3D = mask.centerOfMass\n\n DRR = forwardProjection(image, projectionAngle, axis=projectionAxis)\n DRRMask = forwardProjection(mask, projectionAngle, axis=projectionAxis)\n\n halfDiff = int((DRR.shape[1] - image.gridSize[2])/2) ## not sure this will work if orientation is changed\n croppedDRR = DRR[:, halfDiff + 1:DRR.shape[1] - halfDiff - 1] ## not sure this will work if orientation is changed\n croppedDRRMask = DRRMask[:, halfDiff + 1:DRRMask.shape[1] - halfDiff - 1] ## not sure this will work if orientation is changed\n\n if outputSize:\n # print('Before resampling')\n # print(croppedDRR.shape, np.min(croppedDRR), np.max(croppedDRR), np.mean(croppedDRR))\n ratio = [outputSize[0]/croppedDRR.shape[0], outputSize[1]/croppedDRR.shape[1]]\n croppedDRR = zoom(croppedDRR, ratio)\n croppedDRRMask = zoom(croppedDRRMask, ratio)\n # print('After resampling')\n # print(croppedDRR.shape, np.min(croppedDRR), np.max(croppedDRR), np.mean(croppedDRR))\n\n binaryDRRMask = getBinaryMaskFromROIDRR(croppedDRRMask)\n centerOfMass = get2DMaskCenterOfMass(binaryDRRMask)\n # print('CenterOfMass:', centerOfMass)\n\n del image # to release the RAM\n del mask # to release the RAM\n\n print('Deformations and projections finished for deformation', deformation.name)\n\n # plt.figure()\n # plt.subplot(1, 5, 1)\n # plt.imshow(DRR)\n # plt.subplot(1, 5, 2)\n # plt.imshow(croppedDRR)\n # plt.subplot(1, 5, 3)\n # plt.imshow(DRRMask)\n # plt.subplot(1, 5, 4)\n # plt.imshow(croppedDRRMask)\n # plt.subplot(1, 5, 5)\n # plt.imshow(binaryDRRMask)\n # plt.show()\n\n return [croppedDRR, binaryDRRMask, centerOfMass, centerOfMass3D]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "load the patient data structure and get the model and RTStruct\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = loadDataStructure(dataPath)[0]\ndynMod = patient.getPatientDataOfType(\"Dynamic3DModel\")[0]\nrtStruct = patient.getPatientDataOfType(\"RTStruct\")[0]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "get the ROI and mask on which we want to apply the motion signal\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "print('Available ROIs')\nrtStruct.print_ROINames()\n\ngtvContour = rtStruct.getContourByName(otherContourToUse)\nGTVMask = gtvContour.getBinaryMask(origin=dynMod.midp.origin, gridSize=dynMod.midp.gridSize, spacing=dynMod.midp.spacing)\ngtvBox = getBoxAroundROI(GTVMask)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "get the body contour to adjust the crop in the direction of the DRR projection\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "bodyContour = rtStruct.getContourByName(bodyContourToUse)\nbodyMask = bodyContour.getBinaryMask(origin=dynMod.midp.origin, gridSize=dynMod.midp.gridSize, spacing=dynMod.midp.spacing)\nbodyBox = getBoxAroundROI(bodyMask)\n\nif projAngle == 0 and projAxis == 'Z': # coronal\n croppingBox = [gtvBox[0], bodyBox[1], gtvBox[2]] ## create the used box combining the two boxes\nelif projAngle == 90 and projAxis == 'Z': # sagittal\n croppingBox = [bodyBox[0], gtvBox[1], gtvBox[2]]\nelif projAngle == 0 and projAxis == 'X': # coronal\n croppingBox = [gtvBox[0], bodyBox[1], gtvBox[2]]\nelif projAngle == 0 and projAxis == 'Y': # sagittal\n croppingBox = [bodyBox[0], gtvBox[1], gtvBox[2]]\nelse:\n print('Do not know how to handle crop in this axis/angle configuration, so the body is used')\n croppingBox = [bodyBox[0], bodyBox[1], bodyBox[2]]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "crop the model data using the box\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "crop3DDataAroundBox(dynMod, croppingBox, marginInMM=marginInMM)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "if you want to see the crop in the opentps_core you can save the data in cropped version\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "saveSerializedObjects(patient, savingPath + 'croppedModelAndROIs')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "get the 3D center of mass of this ROI\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "gtvCenterOfMass = gtvContour.getCenterOfMass(dynMod.midp.origin, dynMod.midp.gridSize, dynMod.midp.spacing)\ngtvCenterOfMassInVoxels = getVoxelIndexFromPosition(gtvCenterOfMass, dynMod.midp)\nprint('Used ROI name', gtvContour.name)\nprint('Used ROI center of mass :', gtvCenterOfMass)\nprint('Used ROI center of mass in voxels:', gtvCenterOfMassInVoxels)\n\nif amplitude == 'model':\n ## to get amplitude from model !!! it takes some time because 10 displacement fields must be computed just for this\n modelValues = getAverageModelValuesAroundPosition(gtvCenterOfMass, dynMod, dimensionUsed='Z')\n amplitude = np.max(modelValues) - np.min(modelValues)\n print('Amplitude of deformation at ROI center of mass', amplitude)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Signal creation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "newSignal = SyntheticBreathingSignal(amplitude=amplitude,\n variationAmplitude=variationAmplitude,\n breathingPeriod=breathingPeriod,\n variationFrequency=variationFrequency,\n shift=shift,\n meanNoise=meanNoise,\n varianceNoise=varianceNoise,\n samplingPeriod=samplingPeriod,\n simulationTime=sequenceDurationInSecs,\n meanEvent=meanEvent)\n\nnewSignal.generate1DBreathingSignal()\n\npointList = [gtvCenterOfMass]\npointVoxelList = [gtvCenterOfMassInVoxels]\nsignalList = [newSignal]\n\nsaveSerializedObjects([signalList, pointList], savingPath + 'ROIsAndSignalObjects')\nfor signalIndex in range(len(signalList)):\n signalList[signalIndex] = signalList[signalIndex].breathingSignal" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "show signals and ROIs\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "prop_cycle = plt.rcParams['axes.prop_cycle']\ncolors = prop_cycle.by_key()['color']\nplt.figure(figsize=(12, 6))\nsignalAx = plt.subplot(2, 1, 2)\n\nfor pointIndex, point in enumerate(pointList):\n ax = plt.subplot(2, 2 * len(pointList), 2 * pointIndex + 1)\n ax.set_title('Slice Y:' + str(pointVoxelList[pointIndex][1]))\n ax.imshow(np.rot90(dynMod.midp.imageArray[:, pointVoxelList[pointIndex][1], :]))\n ax.scatter([pointVoxelList[pointIndex][0]], [dynMod.midp.imageArray.shape[2] - pointVoxelList[pointIndex][2]],\n c=colors[pointIndex], marker=\"x\", s=100)\n ax2 = plt.subplot(2, 2 * len(pointList), 2 * pointIndex + 2)\n ax2.set_title('Slice Z:' + str(pointVoxelList[pointIndex][2]))\n ax2.imshow(np.rot90(dynMod.midp.imageArray[:, :, pointVoxelList[pointIndex][2]]))\n ax2.scatter([pointVoxelList[pointIndex][0]], [dynMod.midp.imageArray.shape[2] - pointVoxelList[pointIndex][2]],\n c=colors[pointIndex], marker=\"x\", s=100)\n signalAx.plot(newSignal.timestamps / 1000, signalList[pointIndex], c=colors[pointIndex])\n\nsignalAx.set_xlabel('Time (s)')\nsignalAx.set_ylabel('Deformation amplitude in Z direction (mm)')\nplt.savefig(savingPath + 'ROI_And_Signals_fig.pdf', dpi=300)\nplt.show()\n\n## -------------------------------------------------------------\n\nsequenceSize = newSignal.breathingSignal.shape[0]\nprint('Sequence Size =', sequenceSize, 'split by stack of ', subSequenceSize, '. Multiprocessing =', multiprocessing)\n\nsubSequencesIndexes = [subSequenceSize * i for i in range(math.ceil(sequenceSize/subSequenceSize))]\nsubSequencesIndexes.append(sequenceSize)\nprint('Sub sequences indexes', subSequencesIndexes)\n\nstartTime = time.time()\n\nif multiprocessing == False:\n\n resultList = []\n\n for i in range(len(subSequencesIndexes)-1):\n print('Creating deformations for images', subSequencesIndexes[i], 'to', subSequencesIndexes[i + 1] - 1)\n\n deformationList = generateDeformationListFromBreathingSignalsAndModel(dynMod,\n signalList,\n pointList,\n signalIdxUsed=[subSequencesIndexes[i], subSequencesIndexes[i+1]],\n dimensionUsed='Z',\n outputType=np.float32)\n\n for deformationIndex, deformation in enumerate(deformationList):\n resultList.append(deformImageAndMaskAndComputeDRRs(dynMod.midp,\n GTVMask,\n deformation,\n projectionAngle=projAngle,\n projectionAxis=projAxis,\n outputSize=outputSize,\n tryGPU=True))\n\n\n savingPath += dataSetDataFolder + f'Patient_0_{sequenceSize}_DRRMasksAndCOM'\n saveSerializedObjects(resultList, savingPath + str(sequenceSize))\n\n\nelif multiprocessing == True:\n\n resultList = []\n\n if subSequenceSize > maxMultiProcUse: ## re-adjust the subSequenceSize since this will be done in multi processing\n subSequenceSize = maxMultiProcUse\n print('SubSequenceSize put to', maxMultiProcUse, 'for multiprocessing.')\n print('Sequence Size =', sequenceSize, 'split by stack of ', subSequenceSize, '. Multiprocessing =', multiprocessing)\n subSequencesIndexes = [subSequenceSize * i for i in range(math.ceil(sequenceSize / subSequenceSize))]\n subSequencesIndexes.append(sequenceSize)\n\n for i in range(len(subSequencesIndexes)-1):\n print('Creating deformations for images', subSequencesIndexes[i], 'to', subSequencesIndexes[i + 1] - 1)\n\n deformationList = generateDeformationListFromBreathingSignalsAndModel(dynMod,\n signalList,\n pointList,\n signalIdxUsed=[subSequencesIndexes[i], subSequencesIndexes[i+1]],\n dimensionUsed='Z',\n outputType=np.float32)\n\n print('Start multi process deformation with', len(deformationList), 'deformations')\n with concurrent.futures.ProcessPoolExecutor() as executor:\n\n results = executor.map(deformImageAndMaskAndComputeDRRs, repeat(dynMod.midp), repeat(GTVMask), deformationList, repeat(projAngle), repeat(projAxis), repeat(tryGPU), repeat(outputSize))\n resultList += results\n\n print('ResultList lenght', len(resultList))\n\n savingPath += dataSetDataFolder + f'Patient_0_{sequenceSize}_DRRMasksAndCOM_multiProcTest'\n saveSerializedObjects(resultList, savingPath)\n\nstopTime = time.time()\nprint('Test with multiprocessing =', multiprocessing, '. Sub-sequence size:', str(subSequenceSize), 'finished in', np.round(stopTime - startTime, 2) / 60, 'minutes')\nprint(np.round((stopTime - startTime)/len(resultList), 2), 'sec per sample')" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/DICOMExample/run_generateRandomSamplesFromModel.ipynb b/auto_examples/DICOMExample/run_generateRandomSamplesFromModel.ipynb new file mode 100644 index 0000000..f20d67c --- /dev/null +++ b/auto_examples/DICOMExample/run_generateRandomSamplesFromModel.ipynb @@ -0,0 +1,133 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Generate random samples from a model\nauthor: OpenTPS team\n\nThis example demonstrates how to generate random samples from a model using OpenTPS.\n\nrunning time: ~ 5 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport sys\nimport time\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom pathlib import Path\nimport opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.serializedObjectIO import loadDataStructure\nfrom opentps.core.processing.deformableDataAugmentationToolBox.generateRandomSamplesFromModel import generateRandomImagesFromModel, generateRandomDeformationsFromModel" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Data path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "testDataPath = os.path.join(Path(os.getcwd()).parent.absolute(), 'testData/')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "read a serialized dynamic sequence\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "dataPath = testDataPath + \"veryLightDynMod.p\"\ndynMod = loadDataStructure(dataPath)[0]\n\ntryGPU = True\nnumberOfSamples = 50\n\nimageList = []\n\nstartTime = time.time()\n\ndefList = generateRandomDeformationsFromModel(dynMod, numberOfSamples=numberOfSamples, ampDistribution='gaussian')\nfor deformation in defList:\n imageList.append(deformation.deformImage(dynMod.midp, fillValue='closest', tryGPU=tryGPU))\n\nprint(len(imageList))\nprint('first test done in ', np.round(time.time() - startTime, 2))\n\nplt.figure()\nplt.imshow(imageList[0].imageArray[:, 50, :])\nplt.show()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "startTime = time.time()\nimageList = generateRandomImagesFromModel(dynMod, numberOfSamples=numberOfSamples, ampDistribution='gaussian', tryGPU=tryGPU)\nprint('second test done in ', np.round(time.time() - startTime, 2))\n\nplt.figure()\nplt.imshow(imageList[0].imageArray[:, 50, :])\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_4DProtonOptimization.ipynb b/auto_examples/PlanOptimization/run_4DProtonOptimization.ipynb new file mode 100644 index 0000000..0f765e8 --- /dev/null +++ b/auto_examples/PlanOptimization/run_4DProtonOptimization.ipynb @@ -0,0 +1,284 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# 4D Proton Optimization\nauthor: OpenTPS team\n\nThis example shows how to create and optimize a 4D proton plan using OpenTPS.\n\nrunning time: ~ 1 hours\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport logging\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport sys\nimport pydicom\nimport datetime\nsys.path.append('..')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.dataLoader import readData\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import loadRTPlan, saveRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import IntensityModulationOptimizer\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D\nfrom opentps.core.io.dicomIO import writeRTDose, readDicomDose\n\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = os.path.join(os.getcwd(), 'Exemple_Robust4DOptimization')\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT and ROI creation\n Generic example: 4DCT composed of 3 CTs : 2 phases and the MidP. \n The anatomy consists of a square target moving vertically, with an organ at risk and soft tissue (muscle) in front of it. \n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "CT4D = []\nROI4D = []\nfor i in range(0, 3):\n # ++++Don't delete UIDs to build the simple study+++++++++++++++++++\n studyInstanceUID = pydicom.uid.generate_uid()\n ctSeriesInstanceUID = pydicom.uid.generate_uid()\n frameOfReferenceUID = pydicom.uid.generate_uid()\n # structSeriesInstanceUID = pydicom.uid.generate_uid()\n dt = datetime.datetime.now()\n #++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n\n # CT\n patient = Patient()\n patient.name = f'Miro_OpenTPS_4DCT'\n Patient.id = f'12082024'\n Patient.birthDate = dt.strftime('%Y%m%d')\n patient.sex = \"\"\n \n ctSize = 150\n ct = CTImage(seriesInstanceUID=ctSeriesInstanceUID, frameOfReferenceUID=frameOfReferenceUID)\n ct.name = f'CT_Phase_{i}'\n ct.patient = patient\n ct.studyInstanceUID = studyInstanceUID\n\n huWater = 50\n huTarget = 100\n huMuscle = 200\n data = huWater * np.ones((ctSize, ctSize, ctSize))\n\n # Muscle\n data[100:140, 20:130, 55:95] = huMuscle\n # OAR\n data[70:80, 70:80, 65:85] = huTarget\n # TargetVolume\n if i == 0 :\n data[25:45, 70:100, 65:85] = huTarget\n if i == 1 :\n data[25:45, 60:90, 65:85] = huTarget\n if i == 2 :\n data[25:45, 50:80, 65:85] = huTarget\n ct.imageArray = data\n # writeDicomCT(ct, output_path)\n\n #---------------------ROI\n ROI = []\n\n # TargetVolume\n TV = ROIMask()\n TV.patient = patient\n TV.name = 'TV'\n TV.color = (255, 0, 0) # red\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n if i == 0 :\n data[25:45, 70:100, 65:85] = True\n if i == 1 :\n data[25:45, 60:90, 65:85] = True\n if i == 2 :\n data[25:45, 50:80, 65:85] = True\n TV.imageArray = data\n ROI.append(TV)\n \n # Muscle\n Muscle = ROIMask()\n Muscle.patient = patient\n Muscle.name = 'Muscle'\n Muscle.color = (150, 0, 0)\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n data[100:140, 20:130, 55:95] = True\n Muscle.imageArray = data\n ROI.append(Muscle)\n\n # OAR\n OAR = ROIMask()\n OAR.patient = patient\n OAR.name = 'OAR'\n OAR.color = (100, 0, 0)\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n data[70:80, 70:80, 65:85] = True\n OAR.imageArray = data\n ROI.append(OAR)\n\n # Body\n BODY = ROIMask()\n BODY.patient = patient\n BODY.name = 'Body'\n BODY.color = (100, 0, 0)\n data = np.ones((ctSize, ctSize, ctSize)).astype(bool)\n data[np.where(OAR.imageArray)] = False\n data[np.where(Muscle.imageArray)] = False\n data[np.where(TV.imageArray)] = False\n BODY.imageArray = data\n ROI.append(BODY)\n\n CT4D.append(ct)\n ROI4D.append(ROI)\n\nRefCT = CT4D[1]\nRefTV = ROI4D[1][0]\nRefOAR = ROI4D[1][2]\nRefBody = ROI4D[1][3]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Design plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "beamNames = [\"Beam1\"]\ngantryAngles = [90.]\ncouchAngles = [0.]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configure MCsquare\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 1e3\nmc2.ctCalibration = ctCalibration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate new plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan_file = os.path.join(output_path, f\"RobustPlan_4D.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n logger.info('Plan loaded')\n\nelse:\n planDesign = ProtonPlanDesign()\n planDesign.ct = RefCT # Here, it's the MidP\n planDesign.targetMask = RefTV\n planDesign.gantryAngles = gantryAngles\n planDesign.beamNames = beamNames\n planDesign.couchAngles = couchAngles\n planDesign.calibration = ctCalibration\n\n # Robustness settings\n planDesign.robustness.setupSystematicError = [1.6, 1.6, 1.6] # mm (sigma)\n planDesign.robustness.setupRandomError = [0.0, 0.0, 0.0] # mm (sigma)\n planDesign.robustness.rangeSystematicError = 3.0 # %\n\n # 4D Evaluation mode\n planDesign.robustness.Mode4D = planDesign.robustness.Mode4D.MCsquareAccumulation # Or MCsquareSystematic\n planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.REDUCED_SET # RANDOM not available for MCsquareSystematic\n # planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.ALL (includes diagonals on sphere)\n # planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.RANDOM\n planDesign.robustness.numScenarios = 50 # Specify how many random scenarios to simulate, default = 100\n\n # # 4D settings : only for the mode MCsquareAccumulation with the RANDOM strategie\n # planDesign.robustness.Create4DCTfromRef = True\n # planDesign.robustness.SystematicAmplitudeError = 5.0 # % # Only with RANDOM strategie\n # planDesign.robustness.RandomAmplitudeError = 5.0 # %\n # planDesign.robustness.Dynamic_delivery = True\n # planDesign.robustness.SystematicPeriodError = 5.0 # % # Spot timing required. If not, we calculate them with SimpleBeamDeliveryTimings()\n # planDesign.robustness.RandomPeriodError = 5.0 # %\n # planDesign.robustness.Breathing_period = 1 # x100% # default value\n\n planDesign.spotSpacing = 10.0 \n planDesign.layerSpacing = 10.0 \n planDesign.targetMargin = 7 # Enough to encompass target motion\n\n planDesign.defineTargetMaskAndPrescription(target = RefTV, targetPrescription = 60.)\n\n plan = planDesign.buildPlan()\n plan.rtPlanName = f\"RobustPlan_4D\"\n\n # refIndex : \n # ACCUMULATED -> Index of the Image in the 4DCT one wish we will accumulate the dose.\n ## SYSTEMATIC -> Index of the Image in the 4DCT who will be used as the nominal. So the one closer to the MidP. Or the Midp.\n\n nominal, scenarios = mc2.compute4DRobustScenarioBeamlets(CT4D, plan, refIndex=1, roi=ROI4D, storePath=output_path)\n\n plan.planDesign.beamlets = nominal\n plan.planDesign.robustness.scenarios = scenarios\n plan.planDesign.robustness.numScenarios = len(scenarios)\n saveRTPlan(plan, plan_file)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Set objectives\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan.planDesign.objectives.addFidObjective(RefTV, FidObjective.Metrics.DMAX, limitValue = 63.0, weight = 100.0, robust=True)\nplan.planDesign.objectives.addFidObjective(RefTV, FidObjective.Metrics.DMIN, limitValue = 60.0, weight = 100.0, robust=True)\nplan.planDesign.objectives.addFidObjective(RefOAR, FidObjective.Metrics.DMAX, limitValue = 40.0, weight = 80.0)\nplan.planDesign.objectives.addFidObjective(RefBody, FidObjective.Metrics.DMAX, limitValue = 40.0, weight = 80.0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize treatment plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "DoseFile = 'DoseRobustPlan4D'\nDose_file = os.path.join(output_path, DoseFile + '.dcm')\nif os.path.isfile(Dose_file):\n doseImage = readDicomDose(Dose_file)\n print('Dose imported')\nelse :\n plan.planDesign.ROI_cropping = False\n solver = IntensityModulationOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=150)\n # Optimize treatment plan\n doseImage, ps = solver.optimize()\n saveRTPlan(plan, os.path.join(output_path, \"RobustPlan_4D_weighted.tps\"))\n writeRTDose(doseImage, output_path, DoseFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Display results\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(RefTV, doseImage)\nprint('TV -> D95 = ' + str(target_DVH.D95) + ' Gy')\nprint('TV -> D5 = ' + str(target_DVH.D5) + ' Gy')\nprint('TV -> D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))\n\noar = resampleImage3DOnImage3D(RefOAR, ct)\noar_DVH = DVH(oar, doseImage)\nprint('OAR -> D95 = ' + str(oar_DVH.D95) + ' Gy')\nprint('OAR -> DMAX = ' + str(oar_DVH.Dmax) + ' Gy')\n\nBody = resampleImage3DOnImage3D(RefBody, ct)\nBody_DVH = DVH(Body, doseImage)\nprint('Body -> D95 = ' + str(Body_DVH.D95) + ' Gy')\nprint('Body -> DMAX = ' + str(Body_DVH.Dmax) + ' Gy')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## center of mass\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "RefTV = resampleImage3DOnImage3D(RefTV, RefCT)\nCOM_coord = RefTV.centerOfMass\nCOM_index = RefTV.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = RefCT.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = RefTV.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, RefCT)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)\n\ncontourTargetMask0 = ROI4D[0][0].getBinaryContourMask()\nimg_maskP1 = contourTargetMask0.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask2 = ROI4D[2][0].getBinaryContourMask()\nimg_maskP2 = contourTargetMask2.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourOAR = RefOAR.getBinaryContourMask()\nimg_OAR = contourOAR.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Display dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(12, 5))\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\nax[0].imshow(img_maskP1, alpha=.2, cmap='binary')\nax[0].imshow(img_maskP2, alpha=.2, cmap='binary')\nax[0].imshow(img_OAR, alpha=.2, cmap='binary')\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].plot(oar_DVH.histogram[0], oar_DVH.histogram[1], label=oar_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\n\nplt.savefig(f'{output_path}/DoseRobustOptimization_4D.png', format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_4DprotonEvaluation.ipynb b/auto_examples/PlanOptimization/run_4DprotonEvaluation.ipynb new file mode 100644 index 0000000..252956d --- /dev/null +++ b/auto_examples/PlanOptimization/run_4DprotonEvaluation.ipynb @@ -0,0 +1,230 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# 4D Proton Evaluation\nauthor: OpenTPS team\n\nThis example shows how to evaluate a 4D proton plan using OpenTPS.\n\nrunning time: ~ 20 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport datetime\nimport logging\nimport pydicom\nimport datetime\n\nimport numpy as np\nfrom matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data import Patient\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.planEvaluation.robustnessEvaluation import RobustnessEval\nfrom opentps.core.io.dataLoader import readData\nfrom opentps.core.data import DVH\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = os.path.join(os.getcwd(), 'Exemple_Robust4DOptimization')\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT and ROI creation\n Generic example: 4DCT composed of 3 CTs : 2 phases and the MidP. \n The anatomy consists of a square target moving vertically, with an organ at risk and soft tissue (muscle) in front of it. \n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "CT4D = []\nROI4D = []\nfor i in range(0, 3):\n # ++++Don't delete UIDs to build the simple study+++++++++++++++++++\n studyInstanceUID = pydicom.uid.generate_uid()\n ctSeriesInstanceUID = pydicom.uid.generate_uid()\n frameOfReferenceUID = pydicom.uid.generate_uid()\n # structSeriesInstanceUID = pydicom.uid.generate_uid()\n dt = datetime.datetime.now()\n #++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n\n # CT\n patient = Patient()\n patient.name = f'Miro_OpenTPS_4DCT'\n Patient.id = f'12082024'\n Patient.birthDate = dt.strftime('%Y%m%d')\n patient.sex = \"\"\n \n ctSize = 150\n ct = CTImage(seriesInstanceUID=ctSeriesInstanceUID, frameOfReferenceUID=frameOfReferenceUID)\n ct.name = f'CT_Phase_{i}'\n ct.patient = patient\n ct.studyInstanceUID = studyInstanceUID\n\n huWater = 50\n huTarget = 100\n huMuscle = 200\n data = huWater * np.ones((ctSize, ctSize, ctSize))\n\n # Muscle\n data[100:140, 20:130, 55:95] = huMuscle\n # OAR\n data[70:80, 70:80, 65:85] = huTarget\n # TargetVolume\n if i == 0 :\n data[25:45, 70:100, 65:85] = huTarget\n if i == 1 :\n data[25:45, 60:90, 65:85] = huTarget\n if i == 2 :\n data[25:45, 50:80, 65:85] = huTarget\n ct.imageArray = data\n # writeDicomCT(ct, output_path)\n\n #---------------------ROI\n ROI = []\n\n # TargetVolume\n TV = ROIMask()\n TV.patient = patient\n TV.name = 'TV'\n TV.color = (255, 0, 0) # red\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n if i == 0 :\n data[25:45, 70:100, 65:85] = True\n if i == 1 :\n data[25:45, 60:90, 65:85] = True\n if i == 2 :\n data[25:45, 50:80, 65:85] = True\n TV.imageArray = data\n ROI.append(TV)\n \n # Muscle\n Muscle = ROIMask()\n Muscle.patient = patient\n Muscle.name = 'Muscle'\n Muscle.color = (150, 0, 0)\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n data[100:140, 20:130, 55:95] = True\n Muscle.imageArray = data\n ROI.append(Muscle)\n\n # OAR\n OAR = ROIMask()\n OAR.patient = patient\n OAR.name = 'OAR'\n OAR.color = (100, 0, 0)\n data = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\n data[70:80, 70:80, 65:85] = True\n OAR.imageArray = data\n ROI.append(OAR)\n\n # Body\n BODY = ROIMask()\n BODY.patient = patient\n BODY.name = 'Body'\n BODY.color = (100, 0, 0)\n data = np.ones((ctSize, ctSize, ctSize)).astype(bool)\n data[np.where(OAR.imageArray)] = False\n data[np.where(Muscle.imageArray)] = False\n data[np.where(TV.imageArray)] = False\n BODY.imageArray = data\n ROI.append(BODY)\n\n CT4D.append(ct)\n ROI4D.append(ROI)\n\nRefCT = CT4D[1]\nRefTV = ROI4D[1][0]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Design plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "beamNames = [\"Beam1\"]\ngantryAngles = [90.]\ncouchAngles = [0.]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configure MCsquare\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 1e3\nmc2.statUncertainty = 2.\nmc2.ctCalibration = ctCalibration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate new plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan_file = os.path.join(output_path, f\"RobustPlan_4D_weighted.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n logger.info('Plan weighted loaded')\nelse:\n planDesign = ProtonPlanDesign()\n planDesign.ct = RefCT # Here, it's the MidP\n planDesign.targetMask = RefTV\n planDesign.gantryAngles = gantryAngles\n planDesign.beamNames = beamNames\n planDesign.couchAngles = couchAngles\n planDesign.calibration = ctCalibration\n\n planDesign.spotSpacing = 7.0 \n planDesign.layerSpacing = 7.0 \n planDesign.targetMargin = 10 # Enough to encompass target motion\n\n planDesign.defineTargetMaskAndPrescription(target = RefTV, targetPrescription = 60.)\n\n plan = planDesign.buildPlan()\n plan.rtPlanName = f\"RobustPlan_4D\"" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate scenarios\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "scenario_folder = os.path.join(output_path, 'Robustness4D_Test')\nif os.path.isdir(scenario_folder):\n scenarios = RobustnessEval()\n scenarios.selectionStrategy = RobustnessEval.Strategies.REDUCED_SET\n # scenarios.selectionStrategy = RobustnessEval.Strategies.ALL\n # scenarios.selectionStrategy = RobustnessEval.Strategies.RANDOM\n scenarios.setupSystematicError = plan.planDesign.robustnessEval.setupSystematicError\n scenarios.setupRandomError = plan.planDesign.robustnessEval.setupRandomError\n scenarios.rangeSystematicError = plan.planDesign.robustnessEval.rangeSystematicError\n scenarios.load(scenario_folder)\nelse:\n # MCsquare config for scenario dose computation\n mc2.nbPrimaries = 1e6\n plan.planDesign.robustnessEval = RobustnessEval()\n plan.planDesign.robustnessEval.setupSystematicError = [1.6, 1.6, 1.6] # mm (sigma)\n plan.planDesign.robustnessEval.setupRandomError = [0.0, 0.0, 0.0] # mm (sigma)\n plan.planDesign.robustnessEval.rangeSystematicError = 3.0 # %\n \n # 4D Evaluation mode\n plan.planDesign.robustnessEval.Mode4D = plan.planDesign.robustnessEval.Mode4D.MCsquareAccumulation # Or MCsquareSystematic\n\n # # 4D settings : only for the mode MCsquareAccumulation with the RANDOM strategie\n # plan.planDesign.robustnessEval.Create4DCTfromRef = True\n # plan.planDesign.robustnessEval.SystematicAmplitudeError = 5.0 # %\n # plan.planDesign.robustnessEval.RandomAmplitudeError = 5.0 # %\n # plan.planDesign.robustnessEval.Dynamic_delivery = True\n # plan.planDesign.robustnessEval.SystematicPeriodError = 5.0 # % # Spot timing required. If not, we calculate them with SimpleBeamDeliveryTimings()\n # plan.planDesign.robustnessEval.RandomPeriodError = 5.0 # %\n # plan.planDesign.robustnessEval.Breathing_period = 1 # x100% \n\n # Regular scenario sampling\n plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.REDUCED_SET\n\n # All scenarios (includes diagonals on sphere)\n # plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.ALL\n\n # Random scenario sampling \n # plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.RANDOM\n plan.planDesign.robustnessEval.numScenarios = 50 # Specify how many random scenarios to simulate, default = 100\n \n # Run MCsquare simulation\n scenarios = mc2.compute4DRobustScenario(CT4D, plan = plan, refIndex = 1, roi = ROI4D) # 4D method\n output_folder = os.path.join(output_path, 'Robustness4D_Test')\n scenarios.save(output_folder)\n\n\nscenarios.analyzeErrorSpace(RefCT, \"D95\", RefTV, plan.planDesign.objectives.targetPrescription)\nscenarios.printInfo()\nscenarios.recomputeDVH([RefTV])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Show results\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 1, figsize=(5, 5))\nfor i, dvh_band in enumerate(scenarios.dvhBands):\n color = f'C{i % 10}'\n phigh = ax.plot(dvh_band._dose, dvh_band._volumeHigh, alpha=0)\n plow = ax.plot(dvh_band._dose, dvh_band._volumeLow, alpha=0)\n pNominal = ax.plot(dvh_band._nominalDVH._dose, dvh_band._nominalDVH._volume, label=dvh_band._roiName, color = color)\n pfill = ax.fill_between(dvh_band._dose, dvh_band._volumeHigh, dvh_band._volumeLow, alpha=0.2, color=color)\nax.set_xlabel(\"Dose (Gy)\")\nax.set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.savefig(f'{output_path}/Dose4DEvaluation.png', format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_SimpleOptimizationPhoton.ipynb b/auto_examples/PlanOptimization/run_SimpleOptimizationPhoton.ipynb new file mode 100644 index 0000000..f07cb2f --- /dev/null +++ b/auto_examples/PlanOptimization/run_SimpleOptimizationPhoton.ipynb @@ -0,0 +1,230 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Simple IMPT photon plan optimization\nauthor: OpenTPS team\n\nIn this example, we will create and optimize a simple Photons plan.\n\nrunning time: ~ 15 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport logging\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport sys\nimport copy\nfrom scipy.sparse import csc_matrix\nsys.path.append('..')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dicomIO import writeRTDose, writeDicomCT, writeRTPlan, writeRTStruct\nfrom opentps.core.processing.planOptimization.tools import evaluateClinical\nfrom opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import IntensityModulationOptimizer\nfrom opentps.core.processing.doseCalculation.photons.cccDoseCalculator import CCCDoseCalculator\nfrom opentps.core.data.plan import PhotonPlanDesign\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create synthetic CT and ROI\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Simple_Patient'\n\nctSize = 150\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n# ct.origin = -ctSize/2 * ct.spacing\n\nhuAir = -1024.\nhuWater = 0\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\n#writeDicomCT(ct, output_path)\n\n# Struct\nroi = ROIMask()\nroi.patient = patient\n# roi.origin = -ctSize/2 * ct.spacing\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = 'Output'\nif not os.path.exists(output_path):\n os.makedirs(output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plan Creation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Design plan\nbeamNames = [\"Beam1\", \"Beam2\"]\ngantryAngles = [0., 90.]\ncouchAngles = [0.,0]\n\n## Dose computation from plan\nccc = CCCDoseCalculator(batchSize= 30)\nccc.ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\n\n# Load / Generate new plan\nplan_file = os.path.join(output_path, \"PhotonPlan_WaterPhantom_cropped_resampled.tps\")\n\nif os.path.isfile(plan_file): ### Remove the False to load the plan\n plan = loadRTPlan(plan_file, radiationType='photon')\n logger.info('Plan loaded')\nelse:\n planDesign = PhotonPlanDesign()\n planDesign.ct = ct\n planDesign.targetMask = roi\n planDesign.isocenterPosition_mm = None # None take the center of mass of the target\n planDesign.gantryAngles = gantryAngles\n planDesign.couchAngles = couchAngles\n planDesign.beamNames = beamNames\n planDesign.calibration = ctCalibration\n planDesign.xBeamletSpacing_mm = 5\n planDesign.yBeamletSpacing_mm = 5\n planDesign.targetMargin = 5.0\n planDesign.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \n plan = planDesign.buildPlan() \n\n beamlets = ccc.computeBeamlets(ct, plan) \n doseInfluenceMatrix = copy.deepcopy(beamlets)\n \n plan.planDesign.beamlets = beamlets\n beamlets.storeOnFS(os.path.join(output_path, \"BeamletMatrix_\" + plan.seriesInstanceUID + \".blm\"))\n # Save plan with initial spot weights in serialized format (OpenTPS format)\n saveRTPlan(plan, plan_file)\n\n\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.5, 1.0)\n\nplan.numberOfFractionsPlanned = 30\n\nplan.planDesign.ROI_cropping = False # False, not cropping allows you to keep the dose outside the ROIs and then use the 'shift' evaluation method, which simply shifts the beamlets.\nsolver = IntensityModulationOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=1000)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize treatment plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "doseImage, ps = solver.optimize()\nwriteRTDose(doseImage, output_path)\n\n# Save plan with updated spot weights in serialized format (OpenTPS format)\nplan_file_optimized = os.path.join(output_path, \"Plan_WaterPhantom_cropped_resampled_optimized.tps\")\nsaveRTPlan(plan, plan_file_optimized)\n# Save plan with updated spot weights in dicom format\nplan.patient = patient\n# writeRTPlan(plan, output_path, outputFilename = plan.name )\n# writeDicomPhotonRTPlan(plan, output_path )" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Compute DVH on resampled contour\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(roi, doseImage)\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))\nclinROI = [roi.name, roi.name]\nclinMetric = [\"Dmin\", \"Dmax\"]\nclinLimit = [19., 21.]\nclinObj = {'ROI': clinROI, 'Metric': clinMetric, 'Limit': clinLimit}\nprint('Clinical evaluation')\nevaluateClinical(doseImage, [roi], clinObj)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## center of mass\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Display dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 3, figsize=(15, 5))\nax[0].axes.get_xaxis().set_visible(False)\nax[0].axes.get_yaxis().set_visible(False)\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nax[1].grid(True)\nax[1].legend()\n\nconvData = solver.getConvergenceData()\nx_data = np.linspace(0, convData['time'], len(convData['func_0']))\ny_data = convData['func_0']\nax[2].plot(x_data, y_data , 'bo-', lw=2, label='Fidelity')\nax[2].set_xlabel('Time (s)')\nax[2].set_ylabel('Cost')\nax[2].set_yscale('symlog')\nax2 = ax[2].twiny()\nax2.set_xlabel('Iterations')\nax2.set_xlim(0, convData['nIter'])\nax[2].grid(True)\nplt.savefig(os.path.join(output_path, 'Dose_SimpleOptimizationPhotons.png'))\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_SimpleOptimizationProton.ipynb b/auto_examples/PlanOptimization/run_SimpleOptimizationProton.ipynb new file mode 100644 index 0000000..de6acd3 --- /dev/null +++ b/auto_examples/PlanOptimization/run_SimpleOptimizationProton.ipynb @@ -0,0 +1,230 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Simple IMPT proton plan optimization\nauthor: OpenTPS team\n\nIn this example, we will create and optimize a simple Protons plan.\n\nrunning time: ~ 12 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import math\nimport os\nimport sys\n\nimport numpy as np\nfrom matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan import ObjectivesList\nfrom opentps.core.data.plan import ProtonPlanDesign\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D, resampleImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import IntensityModulationOptimizer" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create synthetic CT and ROI\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Patient'\n\nctSize = 150\n\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\n\nhuAir = -1024.\nhuWater = ctCalibration.convertRSP2HU(1.)\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\n\nroi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configure dose engine\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 5e4\nmc2.ctCalibration = ctCalibration\n\nmc2._independentScoringGrid = True\nscoringSpacing = [2, 2, 2]\nmc2._scoringVoxelSpacing = scoringSpacing" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Design plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "beamNames = [\"Beam1\"]\ngantryAngles = [0.]\ncouchAngles = [0.]\n\nplanInit = ProtonPlanDesign()\nplanInit.ct = ct\nplanInit.gantryAngles = gantryAngles\nplanInit.beamNames = beamNames\nplanInit.couchAngles = couchAngles\nplanInit.calibration = ctCalibration\nplanInit.spotSpacing = 6.0\nplanInit.layerSpacing = 6.0\nplanInit.targetMargin = 0.0\nplanInit.setScoringParameters(scoringSpacing=[2, 2, 2], adapt_gridSize_to_new_spacing=True)\n# needs to be called after scoringGrid settings but prior to spot placement\nplanInit.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \nplan = planInit.buildPlan() # Spot placement\nplan.PlanName = \"NewPlan\"\n\nbeamlets = mc2.computeBeamlets(ct, plan, roi=[roi])\nplan.planDesign.beamlets = beamlets\n# doseImageRef = beamlets.toDoseImage()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## objectives\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.5, 1.0)\n# Other examples of objectives\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMEAN, 20, 1.0) \n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DUNIFORM, 20, 1.0)\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DVHMIN, 19, 1.0, volume = 95)\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DVHMAX, 21, 1.0, volume = 5)\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.EUDMIN, 19.5, 1.0, EUDa = 0.2)\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.EUDMAX, 20, 1.0, EUDa = 1)\n# plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.EUDUNIFORM, 20.5, 1.0, EUDa = 0.5)\n# plan.planDesign.objectives.addFidObjective(BODY, FidObjective.Metrics.DFALLOFF, weight=10, fallOffDistance=1, fallOffLowDoseLevel=0, fallOffHighDoseLevel=21)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "solver = IntensityModulationOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=50)\ndoseImage, ps = solver.optimize()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Final dose computation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2.nbPrimaries = 1e7\ndoseImage = mc2.computeDose(ct, plan)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plots\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Compute DVH on resampled contour\nroiResampled = resampleImage3D(roi, origin=ct.origin, spacing=scoringSpacing)\ntarget_DVH = DVH(roiResampled, doseImage)\nprint('D95 = ' + str(target_DVH.D95) + ' Gy')\nprint('D5 = ' + str(target_DVH.D5) + ' Gy')\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))\n\n# center of mass\nroi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)\n\n#Output path\noutput_path = 'Output'\nif not os.path.exists(output_path):\n os.makedirs(output_path)\n\n# Display dose\nfig, ax = plt.subplots(1, 2, figsize=(12, 5))\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.savefig(os.path.join(output_path, 'SimpleOpti1.png'),format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_beamletFreeOpti.ipynb b/auto_examples/PlanOptimization/run_beamletFreeOpti.ipynb new file mode 100644 index 0000000..6468e4f --- /dev/null +++ b/auto_examples/PlanOptimization/run_beamletFreeOpti.ipynb @@ -0,0 +1,248 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Beamlet Free Optimization\nauthor: OpenTPS team\n\nThis example will present the basis of beamlet optimization with openTPS core.\n\nrunning time: ~ 15 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import numpy as np\nimport os\nfrom matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan import ObjectivesList\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = 'Output'\nif not os.path.exists(output_path):\n os.makedirs(output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Generic CT creation\nwe will first create a generic CT of a box fill with water and air\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "#calibratioin of the CT\nctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)\n\n#creation of the patient object\npatient = Patient()\npatient.name = 'Patient'\n\n#size of the 3D box\nctSize = 150\n\n#creation of the CTImage object\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\nhuAir = -1024. #Hounsfield unit of water\nhuWater = ctCalibration.convertRSP2HU(1.) #convert a stopping power of 1. to HU units\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data #the CT generic image is created" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Region of interest\nwe will now create a region of interest which is a small 3D box of size 20*20*20\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data\n\nimage = plt.imshow(ct.imageArray[110,:,:],cmap='Blues')\nplt.colorbar(image)\nplt.contour(roi.imageArray[110,:,:],colors=\"red\")\nplt.title(\"Created CT with ROI\")\nplt.text(5,40,\"Air\",color= 'black')\nplt.text(5,100,\"Water\",color = 'white')\nplt.text(106,111,\"TV\",color ='red')\nplt.savefig(os.path.join(output_path,'beamFree1.png'),format = 'png')\nplt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configuration of Mcsquare\nTo configure the MCsquare calculator we need to calibrate it with the CT calibration obtained above\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.ctCalibration = ctCalibration\nmc2.nbPrimaries = 1e7" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plan Creation\nWe will now create a plan and set objectives for the optimization and set a goal of 20Gy to the target\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Design plan\n\nbeamNames = [\"Beam1\"]\ngantryAngles = [0.]\ncouchAngles = [0.]\n\n# Generate new plan\n\nplanDesign = ProtonPlanDesign() #create a new plan\nplanDesign.ct = ct\nplanDesign.gantryAngles = gantryAngles\nplanDesign.beamNames = beamNames\nplanDesign.couchAngles = couchAngles\nplanDesign.calibration = ctCalibration\nplanDesign.spotSpacing = 5.0\nplanDesign.layerSpacing = 5.0\nplanDesign.targetMargin = 5.0\n# needs to be called prior to spot placement\nplanDesign.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \n\nplan = planDesign.buildPlan() # Spot placement\nplan.PlanName = \"NewPlan\"\n\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.0, 1.0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Mcsquare beamlet free planOptimization\nNow that we have every needed objects we can compute the optimization through MCsquare. :warning: It may take some time to compute.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "doseImage = mc2.optimizeBeamletFree(ct, plan, [roi])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Dose volume histogram\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(roi, doseImage)\nprint('D95 = ' + str(target_DVH.D95) + ' Gy')\nprint('D5 = ' + str(target_DVH.D5) + ' Gy')\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Center of mass\nHere we look at the part of the 3D CT image where \"stuff is happening\" by getting the CoM. We use the function resampleImage3DOnImage3D to the same array size for both images.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plot of the dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(12, 5))\nax[0].axes.get_xaxis().set_visible(False)\nax[0].axes.get_yaxis().set_visible(False)\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nax[0].set_title(\"Computed dose\")\nax[1].set_title(\"DVH\")\nplt.grid(True)\nplt.legend()\nplt.savefig(os.path.join(output_path,'beamFree2.png'),format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_boundConstraintsOpti.ipynb b/auto_examples/PlanOptimization/run_boundConstraintsOpti.ipynb new file mode 100644 index 0000000..1fdcdb0 --- /dev/null +++ b/auto_examples/PlanOptimization/run_boundConstraintsOpti.ipynb @@ -0,0 +1,248 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Bound Constraint Optimization\nauthor: OpenTPS team\n\nIn this example, we optimize an ion plan (Protons) using the BoundConstraintsOptimizer function.\nThis function allows optimization with constraints on the Monitor Unit (MU) values of each spot.\nIt helps to stay as close as possible to reality when certain machines cannot accept MU/spot values that are too high or too low.\n\nrunning time: ~ 15 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import math\nimport os\nimport logging\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport sys\nsys.path.append('..')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan import ObjectivesList\nfrom opentps.core.data.plan import ProtonPlanDesign\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D, resampleImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import BoundConstraintsOptimizer, IntensityModulationOptimizer\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\nWe will create an output folder to store the results of this example\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = 'Output'\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Generic CT creation\nwe will first create a generic CT of a box fill with water and air\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)\n\npatient = Patient()\npatient.name = 'Patient'\n\nctSize = 150\n\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\nhuAir = -1024.\nhuWater = ctCalibration.convertRSP2HU(1.)\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Region of interest\nwe will now create a region of interest wich is a small 3D box of size 20*20*20\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configuration of Mcsquare\nTo configure the MCsquare calculator we need to calibrate it with the CT calibration obtained above\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.ctCalibration = ctCalibration\nmc2.nbPrimaries = 5e4" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plan Creation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Design plan \nbeamNames = [\"Beam1\"]\ngantryAngles = [0.]\ncouchAngles = [0.]\n\n# Load / Generate new plan\nplan_file = os.path.join(output_path,\"Plan_WaterPhantom_cropped_resampled.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n logger.info('Plan loaded')\nelse:\n planInit = ProtonPlanDesign()\n planInit.ct = ct\n planInit.gantryAngles = gantryAngles\n planInit.beamNames = beamNames\n planInit.couchAngles = couchAngles\n planInit.calibration = ctCalibration\n planInit.spotSpacing = 5.0\n planInit.layerSpacing = 5.0\n planInit.targetMargin = 5.0\n planInit.setScoringParameters(scoringSpacing=[2, 2, 2], adapt_gridSize_to_new_spacing=True)\n # needs to be called after scoringGrid settings but prior to spot placement\n planInit.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \n plan = planInit.buildPlan() # Spot placement\n plan.PlanName = \"NewPlan\"\n\n beamlets = mc2.computeBeamlets(ct, plan, roi=[roi])\n plan.planDesign.beamlets = beamlets\n\n beamlets.storeOnFS(os.path.join(output_path, \"BeamletMatrix_\" + plan.seriesInstanceUID + \".blm\"))\n\n saveRTPlan(plan, plan_file)\n\n# Set objectives (attribut is already initialized in planDesign object)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.5, 1.0)\n\nsolver = BoundConstraintsOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=50, bounds=(0.2, 50))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize treatment plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "doseImage, ps = solver.optimize()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Dose volume histogram\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(roi, doseImage)\nprint('D95 = ' + str(target_DVH.D95) + ' Gy')\nprint('D5 = ' + str(target_DVH.D5) + ' Gy')\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Center of mass\nHere we look at the part of the 3D CT image where \"stuff is happening\" by getting the CoM. We use the function resampleImage3DOnImage3D to the same array size for both images.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plot of the dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(12, 5))\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') \ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.show()\nplt.savefig(f'{output_path}/Dose_BoundContraintOpti.png', format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_evaluatePhotonRobustness.ipynb b/auto_examples/PlanOptimization/run_evaluatePhotonRobustness.ipynb new file mode 100644 index 0000000..fb7f880 --- /dev/null +++ b/auto_examples/PlanOptimization/run_evaluatePhotonRobustness.ipynb @@ -0,0 +1,248 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Evaluate Photon Plan Robustness\nauthor: OpenTPS team\n\nThis example shows how to evaluate a photon plan robustness using OpenTPS.\n\nrunning time: ~ 15 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport logging\nimport numpy as np\n\nfrom matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data import Patient\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.planEvaluation.robustnessEvaluation import RobustnessEvalPhoton\nfrom opentps.core.processing.doseCalculation.photons.cccDoseCalculator import CCCDoseCalculator\nfrom opentps.core.data.plan import PhotonPlanDesign\n\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = os.path.join(os.getcwd(), 'Photon_Robust_Output_Example')\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Generic example: box of water with squared target\nctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT and ROI creation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Patient'\n\nctSize = 150\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\nhuAir = -1024.\nhuWater = 0\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\n\nroi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Create output folder\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "if not os.path.isdir(output_path):\n os.mkdir(output_path)\n\n# Design plan\nbeamNames = [\"Beam1\", \"Beam2\"]\ngantryAngles = [0., 90.]\ncouchAngles = [0.,0]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Configure MCsquare\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ccc = CCCDoseCalculator(batchSize= 30)\nccc.ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Load / Generate new plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan_file = os.path.join(output_path, \"Plan_Photon_WaterPhantom_notCropped_optimized.tps\")\n# plan_file = os.path.join(output_path, \"Plan_Photon_WaterPhantom_cropped_optimized.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file, 'photon')\n logger.info('Plan loaded')\nelse:\n planDesign = PhotonPlanDesign()\n planDesign.ct = ct\n planDesign.targetMask = roi\n planDesign.isocenterPosition_mm = None # None take the center of mass of the target\n planDesign.gantryAngles = gantryAngles\n planDesign.couchAngles = couchAngles\n planDesign.beamNames = beamNames\n planDesign.calibration = ctCalibration\n planDesign.xBeamletSpacing_mm = 5\n planDesign.yBeamletSpacing_mm = 5\n planDesign.targetMargin = 5.0\n planDesign.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \n plan = planDesign.buildPlan()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Load / Generate scenarios\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "scenario_folder = os.path.join(output_path, \"RobustnessTest\")\nif os.path.isdir(scenario_folder):\n scenarios = RobustnessEvalPhoton()\n scenarios.selectionStrategy = RobustnessEvalPhoton.Strategies.DEFAULT\n scenarios.setupSystematicError = plan.planDesign.robustnessEval.setupSystematicError\n scenarios.setupRandomError = plan.planDesign.robustnessEval.setupRandomError\n scenarios.load(scenario_folder)\nelse:\n \n # Robust config for scenario dose computation\n plan.planDesign.robustnessEval = RobustnessEvalPhoton()\n plan.planDesign.robustnessEval.setupSystematicError = [1.6, 1.6, 1.6] #sigma (mm)\n plan.planDesign.robustnessEval.setupRandomError = [1.4, 1.4, 1.4] #sigma (mm)\n\n # Strategy selection \n plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.REDUCED_SET\n # plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.ALL\n # plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.RANDOM\n # plan.planDesign.robustnessEval.NumScenarios = 50\n\n plan.planDesign.robustnessEval.doseDistributionType = \"Nominal\"\n\n plan.patient = None\n\n # run MCsquare simulation\n scenarios = ccc.computeRobustScenario(ct, plan, roi = [roi], robustMode = \"Shift\") # 'Simulation' for total recomputation\n output_folder = os.path.join(output_path, \"RobustnessTest\")\n scenarios.save(output_folder)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Robustness analysis\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "scenarios.recomputeDVH([roi])\nscenarios.analyzeErrorSpace(ct, \"D95\", roi, plan.planDesign.objectives.targetPrescription)\nscenarios.printInfo()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Display DVH + DVH-bands\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 1, figsize=(5, 5))\nfor dvh_band in scenarios.dvhBands:\n phigh = ax.plot(dvh_band._dose, dvh_band._volumeHigh, alpha=0)\n plow = ax.plot(dvh_band._dose, dvh_band._volumeLow, alpha=0)\n pNominal = ax.plot(dvh_band._nominalDVH._dose, dvh_band._nominalDVH._volume, label=dvh_band._roiName, color = 'C0')\n pfill = ax.fill_between(dvh_band._dose, dvh_band._volumeHigh, dvh_band._volumeLow, alpha=0.2, color='C0')\nax.set_xlabel(\"Dose (Gy)\")\nax.set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.savefig(os.path.join(output_path, 'Evaluation_RobustOptimizationPhotons.png'))\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_evaluateProtonRobustness.ipynb b/auto_examples/PlanOptimization/run_evaluateProtonRobustness.ipynb new file mode 100644 index 0000000..7cfd1f3 --- /dev/null +++ b/auto_examples/PlanOptimization/run_evaluateProtonRobustness.ipynb @@ -0,0 +1,248 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Evaluate Proton Plan Robustness\nauthor: OpenTPS team\n\nIn this example, we evaluate an optimized ion plan. \nIt is possible to assess range and setup errors and generate DVHs.\n\nrunning time: ~ 45 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport datetime\nimport logging\n\nimport numpy as np\nfrom matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data import Patient\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.planEvaluation.robustnessEvaluation import RobustnessEval\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = os.path.join(os.getcwd(), 'Proton_Robust_Output_Example')\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT and ROI creation\nGeneric example: box of water with squared target\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Patient'\n\nctSize = 150\n\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\n\nhuAir = -1024.\nhuWater = ctCalibration.convertRSP2HU(1.)\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\n\nroi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Design plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "beamNames = [\"Beam1\"]\ngantryAngles = [90.]\ncouchAngles = [0.]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configure MCsquare\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 1e3\nmc2.statUncertainty = 2.\nmc2.ctCalibration = ctCalibration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate new plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan_file = os.path.join(output_path, \"Plan_Proton_WaterPhantom_cropped_optimized.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n print('Plan loaded')\nelse:\n planInit = ProtonPlanDesign()\n planInit.ct = ct\n planInit.gantryAngles = gantryAngles\n planInit.beamNames = beamNames\n planInit.couchAngles = couchAngles\n planInit.calibration = ctCalibration\n planInit.spotSpacing = 10.0\n planInit.layerSpacing = 10.0\n planInit.targetMargin = 5.0\n planInit.setScoringParameters(scoringSpacing=[2, 2, 2], adapt_gridSize_to_new_spacing=True)\n # needs to be called after scoringGrid settings but prior to spot placement\n planInit.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) \n \n plan = planInit.buildPlan() # Spot placement\n plan.PlanName = \"NewPlan\"" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate scenarios\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "scenario_folder = os.path.join(output_path,'RobustnessTest')\nif os.path.isdir(scenario_folder):\n scenarios = RobustnessEval()\n scenarios.selectionStrategy = RobustnessEval.Strategies.ALL\n scenarios.setupSystematicError = plan.planDesign.robustnessEval.setupSystematicError\n scenarios.setupRandomError = plan.planDesign.robustnessEval.setupRandomError\n scenarios.rangeSystematicError = plan.planDesign.robustnessEval.rangeSystematicError\n scenarios.load(scenario_folder)\nelse:\n # MCsquare config for scenario dose computation\n mc2.nbPrimaries = 1e7\n plan.planDesign.robustnessEval = RobustnessEval()\n plan.planDesign.robustnessEval.setupSystematicError = [5.0, 5.0, 5.0] # mm\n plan.planDesign.robustnessEval.setupRandomError = [0.0, 0.0, 0.0] # mm (sigma)\n plan.planDesign.robustnessEval.rangeSystematicError = 3.0 # %\n\n # Regular scenario sampling\n plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.REDUCED_SET\n\n # All scenarios (includes diagonals on sphere)\n # plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.ALL\n\n # Random scenario sampling \n #plan.planDesign.robustnessEval.selectionStrategy = plan.planDesign.robustnessEval.Strategies.RANDOM\n plan.planDesign.robustnessEval.numScenarios = 30 # specify how many random scenarios to simulate, default = 100\n \n plan.patient = None\n # run MCsquare simulation\n scenarios = mc2.computeRobustScenario(ct, plan, [roi])\n output_folder = os.path.join(output_path, \"RobustnessTest\")\n scenarios.save(output_folder)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Robustness analysis\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "scenarios.analyzeErrorSpace(ct, \"D95\", roi, plan.planDesign.objectives.targetPrescription)\nscenarios.printInfo()\nscenarios.recomputeDVH([roi])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Display DVH + DVH-bands\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 1, figsize=(5, 5))\nfor dvh_band in scenarios.dvhBands:\n phigh = ax.plot(dvh_band._dose, dvh_band._volumeHigh, alpha=0)\n plow = ax.plot(dvh_band._dose, dvh_band._volumeLow, alpha=0)\n pNominal = ax.plot(dvh_band._nominalDVH._dose, dvh_band._nominalDVH._volume, label=dvh_band._roiName, color = 'C0')\n pfill = ax.fill_between(dvh_band._dose, dvh_band._volumeHigh, dvh_band._volumeLow, alpha=0.2, color='C0')\nax.set_xlabel(\"Dose (Gy)\")\nax.set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.savefig(f'{output_path}/EvaluateRobustness.png', format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_robustOptimizationProtons.ipynb b/auto_examples/PlanOptimization/run_robustOptimizationProtons.ipynb new file mode 100644 index 0000000..c866998 --- /dev/null +++ b/auto_examples/PlanOptimization/run_robustOptimizationProtons.ipynb @@ -0,0 +1,302 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Robust Proton Plan Optimization\nauthor: OpenTPS team\n\nIn this example, we create and optimize a robust proton plan. \nThe setup and range errors are configurable.\n\nrunning time: ~ 20 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport logging\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport sys\nsys.path.append('..')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.data.images import CTImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data.plan import RobustnessProton\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import loadRTPlan, saveRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import IntensityModulationOptimizer\n\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = os.path.join(os.getcwd(), 'Proton_Robust_Output_Example')\nif not os.path.exists(output_path):\n os.makedirs(output_path)\nlogger.info('Files will be stored in {}'.format(output_path))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create synthetic CT and ROI\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Patient'\n\nctSize = 150\n\nct = CTImage()\nct.name = 'CT'\nct.patient = patient\n\nhuAir = -1024.\nhuWater = ctCalibration.convertRSP2HU(1.)\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\n\nroi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Design plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "beamNames = [\"Beam1\"]\ngantryAngles = [0.]\ncouchAngles = [0.]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create output folder\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "if not os.path.isdir(output_path):\n os.mkdir(output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configure MCsquare\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 1e3\nmc2.ctCalibration = ctCalibration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load / Generate new plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan_file = os.path.join(output_path, \"RobustPlan_notCropped.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n logger.info('Plan loaded')\nelse:\n planDesign = ProtonPlanDesign()\n planDesign.ct = ct\n planDesign.gantryAngles = gantryAngles\n planDesign.beamNames = beamNames\n planDesign.couchAngles = couchAngles\n planDesign.calibration = ctCalibration\n # Robustness settings\n planDesign.robustness = RobustnessProton()\n planDesign.robustness.setupSystematicError = [1.6, 1.6, 1.6] # mm\n planDesign.robustness.setupRandomError = [0.0, 0.0, 0.0] # mm (sigma)\n planDesign.robustness.rangeSystematicError = 5.0 # %\n\n # Regular scenario sampling\n planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.REDUCED_SET\n\n # All scenarios (includes diagonals on sphere)\n # planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.ALL\n\n # Random scenario sampling \n # planDesign.robustness.selectionStrategy = planDesign.robustness.Strategies.RANDOM\n planDesign.robustness.numScenarios = 5 # specify how many random scenarios to simulate, default = 100\n\n planDesign.spotSpacing = 7.0\n planDesign.layerSpacing = 6.0\n planDesign.targetMargin = max(planDesign.spotSpacing, planDesign.layerSpacing) + max(planDesign.robustness.setupSystematicError)\n # scoringGridSize = [int(math.floor(i / j * k)) for i, j, k in zip(ct.gridSize, scoringSpacing, ct.spacing)]\n # planDesign.objectives.setScoringParameters(ct, scoringGridSize, scoringSpacing)\n planDesign.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) # needs to be called prior spot placement\n plan = planDesign.buildPlan() # Spot placement\n plan.PlanName = \"RobustPlan\"\n\n nominal, scenarios = mc2.computeRobustScenarioBeamlets(ct, plan, roi=[roi], storePath=output_path)\n plan.planDesign.beamlets = nominal\n plan.planDesign.robustness.scenarios = scenarios\n plan.planDesign.robustness.numScenarios = len(scenarios)\n \n\n #saveRTPlan(plan, plan_file)\n\n\n\nsaveRTPlan(plan, plan_file)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Set objectives (attribut is already initialized in planDesign object)\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0, robust=True)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.5, 1.0, robust=True)\n\nsolver = IntensityModulationOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=50)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize treatment plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "doseImage, ps = solver.optimize()\n\nplan_file = os.path.join(output_path, \"Plan_Proton_WaterPhantom_cropped_optimized.tps\")\nsaveRTPlan(plan, plan_file, unloadBeamlets=False)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Compute DVH\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(roi, doseImage)\nprint('D95 = ' + str(target_DVH.D95) + ' Gy')\nprint('D5 = ' + str(target_DVH.D5) + ' Gy')\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Center of mass\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Display dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(12, 5))\nax[0].axes.get_xaxis().set_visible(False)\nax[0].axes.get_yaxis().set_visible(False)\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nplt.grid(True)\nplt.legend()\nplt.savefig(os.path.join(output_path, 'Dose_RobustOptimizationProtons.png'))\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/auto_examples/PlanOptimization/run_simpleOptimization_createDicomStudy.ipynb b/auto_examples/PlanOptimization/run_simpleOptimization_createDicomStudy.ipynb new file mode 100644 index 0000000..4e50910 --- /dev/null +++ b/auto_examples/PlanOptimization/run_simpleOptimization_createDicomStudy.ipynb @@ -0,0 +1,302 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n# Simple ion plan optimization and DICOM study creation\nauthor: OpenTPS team\n\nIn this example, we will create and optimize a simple ion (Proton) plan. \nThe generated CT, the plan, and the dose will be saved as DICOM files.\n\nrunning time: ~ 15 minutes\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Setting up the environment in google collab\n First you need to change the type of execution in the bottom left from processor to GPU. Then you can run the example.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import sys\nif \"google.colab\" in sys.modules:\n from IPython import get_ipython\n get_ipython().system('git clone https://gitlab.com/openmcsquare/opentps.git')\n get_ipython().system('pip install ./opentps')\n get_ipython().system('pip install scipy==1.10.1')\n get_ipython().system('pip install cupy-cuda12x')\n import opentps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "imports\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "import os\nimport logging\nimport numpy as np\nfrom matplotlib import pyplot as plt\nimport sys\nimport datetime\nimport pydicom\nsys.path.append('..')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import the needed opentps.core packages\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "from opentps.core.io.dicomIO import writeRTPlan, writeDicomCT, writeRTDose, writeRTStruct\nfrom opentps.core.processing.planOptimization.tools import evaluateClinical\nfrom opentps.core.data.images import CTImage, DoseImage\nfrom opentps.core.data.images import ROIMask\nfrom opentps.core.data.plan import ObjectivesList\nfrom opentps.core.data.plan._protonPlanDesign import ProtonPlanDesign\nfrom opentps.core.data import DVH\nfrom opentps.core.data import Patient\nfrom opentps.core.data import RTStruct\nfrom opentps.core.data.plan import FidObjective\nfrom opentps.core.io import mcsquareIO\nfrom opentps.core.io.scannerReader import readScanner\nfrom opentps.core.io.serializedObjectIO import saveRTPlan, loadRTPlan\nfrom opentps.core.processing.doseCalculation.doseCalculationConfig import DoseCalculationConfig\nfrom opentps.core.processing.doseCalculation.protons.mcsquareDoseCalculator import MCsquareDoseCalculator\nfrom opentps.core.processing.imageProcessing.resampler3D import resampleImage3DOnImage3D, resampleImage3D\nfrom opentps.core.processing.planOptimization.planOptimization import IntensityModulationOptimizer\nfrom opentps.core.data.plan import ProtonPlan\n\nlogger = logging.getLogger(__name__)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## CT calibration and BDL\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "ctCalibration = readScanner(DoseCalculationConfig().scannerFolder)\nbdl = mcsquareIO.readBDL(DoseCalculationConfig().bdlFile)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create the DICOM files\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# ++++Don't delete UIDs to build the simple study+++++++++++++++++++\nstudyInstanceUID = pydicom.uid.generate_uid()\ndoseSeriesInstanceUID = pydicom.uid.generate_uid()\nplanSeriesInstanceUID = pydicom.uid.generate_uid()\nctSeriesInstanceUID = pydicom.uid.generate_uid()\nframeOfReferenceUID = pydicom.uid.generate_uid()\n# structSeriesInstanceUID = pydicom.uid.generate_uid()\ndt = datetime.datetime.now()\n#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output path\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "output_path = 'Output'\nif not os.path.exists(output_path):\n os.makedirs(output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Generic CT creation\nwe will first create a generic CT of a box fill with water and air\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "patient = Patient()\npatient.name = 'Simple_Patient'\nPatient.id = 'Simple_Patient'\nPatient.birthDate = dt.strftime('%Y%m%d')\npatient.sex = \"\"\n\nctSize = 150\nct = CTImage(seriesInstanceUID=ctSeriesInstanceUID, frameOfReferenceUID=frameOfReferenceUID)\nct.name = 'CT'\nct.patient = patient\nct.studyInstanceUID = studyInstanceUID\n\nhuAir = -1024.\nhuWater = ctCalibration.convertRSP2HU(1.)\ndata = huAir * np.ones((ctSize, ctSize, ctSize))\ndata[:, 50:, :] = huWater\nct.imageArray = data\nwriteDicomCT(ct, output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Region of interest\nwe will now create a region of interest wich is a small 3D box of size 20*20*20\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = ROIMask()\nroi.patient = patient\nroi.name = 'TV'\nroi.color = (255, 0, 0) # red\ndata = np.zeros((ctSize, ctSize, ctSize)).astype(bool)\ndata[100:120, 100:120, 100:120] = True\nroi.imageArray = data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Configuration of Mcsquare\nTo configure the MCsquare calculator we need to calibrate it with the CT calibration obtained above\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "mc2 = MCsquareDoseCalculator()\nmc2.beamModel = bdl\nmc2.nbPrimaries = 5e4\nmc2.ctCalibration = ctCalibration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plan Creation\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Design plan\nbeamNames = [\"Beam1\"]\ngantryAngles = [0.]\ncouchAngles = [0.]\n\n# Load / Generate new plan\nplan_file = os.path.join(output_path, \"Plan_WaterPhantom_cropped_resampled.tps\")\n\nif os.path.isfile(plan_file):\n plan = loadRTPlan(plan_file)\n logger.info('Plan loaded')\nelse:\n planDesign = ProtonPlanDesign()\n planDesign.ct = ct\n planDesign.targetMask = roi\n planDesign.gantryAngles = gantryAngles\n planDesign.beamNames = beamNames\n planDesign.couchAngles = couchAngles\n planDesign.calibration = ctCalibration\n planDesign.spotSpacing = 5.0\n planDesign.layerSpacing = 5.0\n planDesign.targetMargin = 5.0\n planDesign.setScoringParameters(scoringSpacing=[2, 2, 2], adapt_gridSize_to_new_spacing=True)\n planDesign.defineTargetMaskAndPrescription(target = roi, targetPrescription = 20.) # needs to be called prior spot placement\n \n plan = planDesign.buildPlan() # Spot placement\n plan.rtPlanName = \"Simple_Patient\"\n\n beamlets = mc2.computeBeamlets(ct, plan, roi=[roi])\n plan.planDesign.beamlets = beamlets\n beamlets.storeOnFS(os.path.join(output_path, \"BeamletMatrix_\" + plan.seriesInstanceUID + \".blm\"))\n # Save plan with initial spot weights in serialized format (OpenTPS format)\n saveRTPlan(plan, plan_file)\n writeRTPlan(plan, output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## objectives\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "# Set objectives (attribut is already initialized in planDesign object)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMAX, 20.0, 1.0)\nplan.planDesign.objectives.addFidObjective(roi, FidObjective.Metrics.DMIN, 20.5, 1.0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Optimize plan\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "plan.seriesInstanceUID = planSeriesInstanceUID\nplan.studyInstanceUID = studyInstanceUID\nplan.frameOfReferenceUID = frameOfReferenceUID\nplan.rtPlanGeometry = \"TREATMENT_DEVICE\"\n\nsolver = IntensityModulationOptimizer(method='Scipy_L-BFGS-B', plan=plan, maxiter=1000)\n# Optimize treatment plan\ndoseImage, ps = solver.optimize()\n\n# Save plan with updated spot weights in serialized format (OpenTPS format)\nplan_file_optimized = os.path.join(output_path, \"Plan_WaterPhantom_cropped_resampled_optimized.tps\")\nsaveRTPlan(plan, plan_file_optimized)\n# Save plan with updated spot weights in dicom format\nplan.patient = patient\nwriteRTPlan(plan, output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Dose volume histogram\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "target_DVH = DVH(roi, doseImage)\nprint('D5 - D95 = {} Gy'.format(target_DVH.D5 - target_DVH.D95))\nclinROI = [roi.name, roi.name]\nclinMetric = [\"Dmin\", \"Dmax\"]\nclinLimit = [19., 21.]\nclinObj = {'ROI': clinROI, 'Metric': clinMetric, 'Limit': clinLimit}\nprint('Clinical evaluation')\nevaluateClinical(doseImage, [roi], clinObj)\n\ndoseImage.referencePlan = plan\ndoseImage.referenceCT = ct\ndoseImage.patient = patient\ndoseImage.studyInstanceUID = studyInstanceUID\ndoseImage.frameOfReferenceUID = frameOfReferenceUID \ndoseImage.sopClassUID = '1.2.840.10008.5.1.4.1.1.481.2'\ndoseImage.mediaStorageSOPClassUID = '1.2.840.10008.5.1.4.1.1.481.2'\ndoseImage.sopInstanceUID = pydicom.uid.generate_uid()\ndoseImage.studyTime = dt.strftime('%H%M%S.%f')\ndoseImage.studyDate = dt.strftime('%Y%m%d')\ndoseImage.SOPInstanceUID = doseImage.sopInstanceUID\n\nif not hasattr(ProtonPlan, \"SOPInstanceUID\"):\n ProtonPlan.SOPInstanceUID = property(lambda self: self.sopInstanceUID)\n\n\nwriteRTDose(doseImage, output_path)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Center of mass\nHere we look at the part of the 3D CT image where \"stuff is happening\" by getting the CoM. We use the function resampleImage3DOnImage3D to the same array size for both images.\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "roi = resampleImage3DOnImage3D(roi, ct)\nCOM_coord = roi.centerOfMass\nCOM_index = roi.getVoxelIndexFromPosition(COM_coord)\nZ_coord = COM_index[2]\n\nimg_ct = ct.imageArray[:, :, Z_coord].transpose(1, 0)\ncontourTargetMask = roi.getBinaryContourMask()\nimg_mask = contourTargetMask.imageArray[:, :, Z_coord].transpose(1, 0)\nimg_dose = resampleImage3DOnImage3D(doseImage, ct)\nimg_dose = img_dose.imageArray[:, :, Z_coord].transpose(1, 0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Plot of the dose\n\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 3, figsize=(15, 5))\nax[0].axes.get_xaxis().set_visible(False)\nax[0].axes.get_yaxis().set_visible(False)\nax[0].imshow(img_ct, cmap='gray')\nax[0].imshow(img_mask, alpha=.2, cmap='binary') # PTV\ndose = ax[0].imshow(img_dose, cmap='jet', alpha=.2)\nplt.colorbar(dose, ax=ax[0])\nax[1].plot(target_DVH.histogram[0], target_DVH.histogram[1], label=target_DVH.name)\nax[1].set_xlabel(\"Dose (Gy)\")\nax[1].set_ylabel(\"Volume (%)\")\nax[1].grid(True)\nax[1].legend()\n\nconvData = solver.getConvergenceData()\nax[2].plot(np.linspace(0, convData['time'], len(convData['func_0'])), convData['func_0'], 'bo-', lw=2,\n label='Fidelity')\nax[2].set_xlabel('Time (s)')\nax[2].set_ylabel('Cost')\nax[2].set_yscale('symlog')\nax2 = ax[2].twiny()\nax2.set_xlabel('Iterations')\nax2.set_xlim(0, convData['nIter'])\nax[2].grid(True)\n\nplt.savefig(f'{output_path}/Dose_SimpleOptimization.png', format = 'png')\nplt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.13" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} \ No newline at end of file diff --git a/sg_execution_times.rst b/sg_execution_times.rst index 4f6f295..ebaa86a 100644 --- a/sg_execution_times.rst +++ b/sg_execution_times.rst @@ -6,7 +6,7 @@ Computation times ================= -**00:14.683** total execution time for 44 files **from all galleries**: +**00:33.871** total execution time for 44 files **from all galleries**: .. container:: @@ -32,38 +32,50 @@ Computation times * - Example - Time - Mem (MB) + * - :ref:`sphx_glr_auto_examples_DICOMExample_run_generateRandomSamplesFromModel.py` (``examples/DICOMExample/run_generateRandomSamplesFromModel.py``) + - 00:11.194 + - 0.0 + * - :ref:`sphx_glr_auto_examples_DoseComputation_run_UseOfRangeShifterForProtonPlan.py` (``examples/DoseComputation/run_UseOfRangeShifterForProtonPlan.py``) + - 00:02.770 + - 0.0 * - :ref:`sphx_glr_auto_examples_DoseComputation_run_SimpleDoseCalculation.py` (``examples/DoseComputation/run_SimpleDoseCalculation.py``) - - 00:02.595 + - 00:02.749 - 0.0 * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_simpleOptimization_createDicomStudy.py` (``examples/PlanOptimization/run_simpleOptimization_createDicomStudy.py``) - - 00:02.053 + - 00:02.321 - 0.0 - * - :ref:`sphx_glr_auto_examples_DoseComputation_run_UseOfRangeShifterForProtonPlan.py` (``examples/DoseComputation/run_UseOfRangeShifterForProtonPlan.py``) - - 00:02.005 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_4DProtonOptimization.py` (``examples/PlanOptimization/run_4DProtonOptimization.py``) + - 00:02.199 + - 0.0 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_SimpleOptimizationProton.py` (``examples/PlanOptimization/run_SimpleOptimizationProton.py``) + - 00:02.082 - 0.0 * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_boundConstraintsOpti.py` (``examples/PlanOptimization/run_boundConstraintsOpti.py``) - - 00:01.844 + - 00:02.079 - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_SimpleOptimizationProton.py` (``examples/PlanOptimization/run_SimpleOptimizationProton.py``) - - 00:01.842 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_robustOptimizationProtons.py` (``examples/PlanOptimization/run_robustOptimizationProtons.py``) + - 00:02.076 + - 0.0 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_beamletFreeOpti.py` (``examples/PlanOptimization/run_beamletFreeOpti.py``) + - 00:01.482 - 0.0 * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_4DprotonEvaluation.py` (``examples/PlanOptimization/run_4DprotonEvaluation.py``) - - 00:01.376 + - 00:01.425 - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_beamletFreeOpti.py` (``examples/PlanOptimization/run_beamletFreeOpti.py``) - - 00:01.193 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_evaluateProtonRobustness.py` (``examples/PlanOptimization/run_evaluateProtonRobustness.py``) + - 00:01.318 - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_SimpleOptimizationPhoton.py` (``examples/PlanOptimization/run_SimpleOptimizationPhoton.py``) - - 00:00.533 + * - :ref:`sphx_glr_auto_examples_DoseComputation_run_protonPlanCreationAndDoseCalculation.py` (``examples/DoseComputation/run_protonPlanCreationAndDoseCalculation.py``) + - 00:00.808 - 0.0 * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_evaluatePhotonRobustness.py` (``examples/PlanOptimization/run_evaluatePhotonRobustness.py``) - - 00:00.528 + - 00:00.543 - 0.0 - * - :ref:`sphx_glr_auto_examples_DoseComputation_run_protonPlanCreationAndDoseCalculation.py` (``examples/DoseComputation/run_protonPlanCreationAndDoseCalculation.py``) - - 00:00.445 + * - :ref:`sphx_glr_auto_examples_PlanOptimization_run_SimpleOptimizationPhoton.py` (``examples/PlanOptimization/run_SimpleOptimizationPhoton.py``) + - 00:00.534 - 0.0 * - :ref:`sphx_glr_auto_examples_DoseComputation_run_photonPlanCreationAndDoseCalculation.py` (``examples/DoseComputation/run_photonPlanCreationAndDoseCalculation.py``) - - 00:00.268 + - 00:00.288 - 0.0 * - :ref:`sphx_glr_auto_community_CommunityExample_SimpleDoseComputationAndOptOnCT.py` (``community/CommunityExample/SimpleDoseComputationAndOptOnCT.py``) - 00:00.000 @@ -71,43 +83,31 @@ Computation times * - :ref:`sphx_glr_auto_community_Template_run_template.py` (``community/Template/run_template.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_DoseDeliverySimulation_PlanDeliverySimulation.py` (``examples/DoseDeliverySimulation/PlanDeliverySimulation.py``) - - 00:00.000 - - 0.0 - * - :ref:`sphx_glr_auto_examples_DoseDeliverySimulation_run_PBSDeliveryTimings.py` (``examples/DoseDeliverySimulation/run_PBSDeliveryTimings.py``) - - 00:00.000 - - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_create3DSeqFromDicom.py` (``examples/OthersExample/create3DSeqFromDicom.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_create3DSeqFromDicom.py` (``examples/DICOMExample/create3DSeqFromDicom.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_create3DSeqFromImages.py` (``examples/OthersExample/create3DSeqFromImages.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_create3DSeqFromImages.py` (``examples/DICOMExample/create3DSeqFromImages.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_createDynamic3DModelFromDicomFields.py` (``examples/OthersExample/createDynamic3DModelFromDicomFields.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_createDynamic3DModelFromDicomFields.py` (``examples/DICOMExample/createDynamic3DModelFromDicomFields.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_createModelWithROIsFromDicomImages.py` (``examples/OthersExample/createModelWithROIsFromDicomImages.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_createModelWithROIsFromDicomImages.py` (``examples/DICOMExample/createModelWithROIsFromDicomImages.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_crop3DDataAroundROI.py` (``examples/OthersExample/crop3DDataAroundROI.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_crop3DDataAroundROI.py` (``examples/DICOMExample/crop3DDataAroundROI.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_exampleImageResampling.py` (``examples/OthersExample/exampleImageResampling.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_exampleImageResampling.py` (``examples/DICOMExample/exampleImageResampling.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_generateDRRAndGTVMasks.py` (``examples/OthersExample/generateDRRAndGTVMasks.py``) + * - :ref:`sphx_glr_auto_examples_DICOMExample_generateDRRAndGTVMasks.py` (``examples/DICOMExample/generateDRRAndGTVMasks.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_OthersExample_run_generateRandomSamplesFromModel.py` (``examples/OthersExample/run_generateRandomSamplesFromModel.py``) - - 00:00.000 - - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_4DProtonOptimization.py` (``examples/PlanOptimization/4DProtonOptimization.py``) - - 00:00.000 - - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_evaluateProtonRobustness.py` (``examples/PlanOptimization/evaluateProtonRobustness.py``) + * - :ref:`sphx_glr_auto_examples_DoseDeliverySimulation_PlanDeliverySimulation.py` (``examples/DoseDeliverySimulation/PlanDeliverySimulation.py``) - 00:00.000 - 0.0 - * - :ref:`sphx_glr_auto_examples_PlanOptimization_robustOptimizationProtons.py` (``examples/PlanOptimization/robustOptimizationProtons.py``) + * - :ref:`sphx_glr_auto_examples_DoseDeliverySimulation_run_PBSDeliveryTimings.py` (``examples/DoseDeliverySimulation/run_PBSDeliveryTimings.py``) - 00:00.000 - 0.0 * - :ref:`sphx_glr_auto_examples_Segmentation_run_Segmentation.py` (``examples/Segmentation/run_Segmentation.py``)