(57)【要約】【課題】混合場のＲＢＥを短時間で精度よく予測して照射パラメータを決定できる照射計画装置、照射計画方法、および荷電粒子照射システムを提供する【解決手段】照射計画装置１０の記憶装置２５に、ドメインの線量平均比エネルギーであるドメイン線量平均比エネルギーと、前記ドメインが多数含まれる細胞核の線量平均比エネルギーである細胞核線量平均比エネルギーと、前記ドメインの飽和線量平均比エネルギーであるドメイン飽和線量平均比エネルギーを記憶し、照射計画装置２０の演算部３３により、前記ペンシルビームから着目位置に与えられる前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーから前記着目位置での生物効果を予測し、前記生物効果に基づいて前記照射パラメータを決定する構成とした。【選択図】図１
2019180908.tif 000002
【特許請求の範囲】
【請求項１】
データを記憶する記憶手段と、演算を行う演算手段とを備え、荷電粒子をペンシルビームとして照射する際の照射パラメータを作成する照射計画装置であって、
前記記憶手段は、
ドメインの線量平均比エネルギーであるドメイン線量平均比エネルギーと、
前記ドメインが多数含まれる細胞核の線量平均比エネルギーである細胞核線量平均比エネルギーと、
前記ドメインの飽和線量平均比エネルギーであるドメイン飽和線量平均比エネルギーを記憶する構成であり、
前記演算手段は、
前記ペンシルビームから着目位置に与えられる前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーから前記着目位置での生物効果を予測し、前記生物効果に基づいて前記照射パラメータを決定する
照射計画装置。
【請求項２】
前記記憶手段は、
前記ペンシルビームの線源データとしての深部線量分布と、
各深さでの前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーを記憶している
請求項１記載の照射計画装置。
【請求項３】
前記演算手段による前記生物効果の予測は、
複数の前記ペンシルビームから着目位置に与えられる前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーを線量加重平均し、前記線量加重平均した値に基づいて前記ペンシルビームの照射による生物効果を予測する構成である
請求項１または２記載の照射計画装置。
【請求項４】
前記生物効果の予測は、次の［数１］を用いる
請求項３記載の照射計画装置。
【数１】
2019180908.tif 000047
【請求項５】
前記照射パラメータの決定は、
前記加重平均した値から前記着目位置での生残率を前記数１により算出してＲＢＥを決定する構成である
請求項４記載の照射計画装置。
【請求項６】
データを記憶する記憶手段と、演算を行う演算手段とを備え、荷電粒子をペンシルビームとして照射する際の照射パラメータを作成する照射計画装置による照射パラメータを決定する照射計画方法であって、
前記記憶手段に、
ドメインの線量平均比エネルギーであるドメイン線量平均比エネルギーと、
前記ドメインが多数含まれる細胞核の線量平均比エネルギーである細胞核線量平均比エネルギーと、
前記ドメインの飽和線量平均比エネルギーであるドメイン飽和線量平均比エネルギーを記憶しておき、
前記演算手段により、
前記ペンシルビームから着目位置に与えられる前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーから前記着目位置での生物効果を予測し、前記生物効果に基づいて前記照射パラメータを決定する。
照射計画方法。
【請求項７】
イオン源で生成した荷電粒子を加速器により加速してペンシルビームとして標的に照射する荷電粒子照射システムであって、
請求項１から５のいずれか１つに記載の照射計画装置により決定された照射パラメータに従って前記ペンシルビームを照射する構成である
荷電粒子照射システム。
【発明の詳細な説明】
【技術分野】
【０００１】
この発明は、例えば、イオン源で生成した荷電粒子を加速器で加速して標的に照射する荷電粒子照射システムの照射パラメータを決定するような照射計画装置、照射計画方法、および荷電粒子照射システムに関する。
【背景技術】
【０００２】
近年、重粒子線治療場（混合場）のＲＢＥを推定する方法として、確率論的マイクロドシメトリックキネティックモデル（ＳＭＫモデル、非特許文献１，２参照）が提案されている。このＳＭＫモデルは、高線量、高ＬＥＴの放射線まで含めて実験値をよく再現することが知られている。
【０００３】
しかし、このモデルは、放射線によりドメインおよび細胞核に与えられる比エネルギーからその放射線の細胞殺傷効果（生物効果）を推定する方法である。そして、このモデルに従って混合場のＲＢＥを予測しようとすると、照射野内の各位置で、重粒子線によってドメインおよび細胞核に与えられる比エネルギーの分布（比エネルギースペクトル）を計算する必要がある。このため、逐次近似繰り返し演算を必要とするスキャニング照射法の治療計画においては、計算時間がかかりすぎ、このモデルを実際の治療の現場で応用するには困難であった。
【先行技術文献】
【特許文献】
【０００４】
【非特許文献１】Tatsuhiko Sato, Yoshiya Furusawa著, “Cell survival fraction estimation based on the probability densities of domain and cell nucleus specific energies using improved microdosimetric kinetic models”, Radiation research誌, vol.178 (2012) pp. 341-356
【非特許文献２】Lorenzo Manganaro, Germano Russo, Roberto Cirio, Federico Dalmasso, Simona Giordanengo, Vincenzo Monaco, Silvia Muraro, Roberto Sacchi, Anna Vignati, Andrea Attili著, “A Monte Carlo approach to the microdosimetric kinetic model to account for dose rate time structure effects in ion beam therapy with application in treatment planning simulations”, Med. Phys.誌、vol. 44 (2017) pp. 1577-1589
【発明の概要】
【発明が解決しようとする課題】
【０００５】
この発明は、上述の問題に鑑みて、混合場のＲＢＥを短時間で精度よく予測して照射パラメータを決定できる照射計画装置、照射計画方法、および荷電粒子照射システムを提供することを目的とする。
【課題を解決するための手段】
【０００６】
この発明は、データを記憶する記憶手段と、演算を行う演算手段とを備え、荷電粒子をペンシルビームとして照射する際の照射パラメータを作成する照射計画装置であって、前記記憶手段は、ドメインの線量平均比エネルギーであるドメイン線量平均比エネルギーと、前記ドメインが多数含まれる細胞核の線量平均比エネルギーである細胞核線量平均比エネルギーと、前記ドメインの飽和線量平均比エネルギーであるドメイン飽和線量平均比エネルギーを記憶する構成であり、前記演算手段は、前記ペンシルビームから着目位置に与えられる前記ドメイン線量平均比エネルギーと前記細胞核線量平均比エネルギーと前記ドメイン飽和線量平均比エネルギーから前記着目位置での生物効果を予測し、前記生物効果に基づいて前記照射パラメータを決定する照射計画装置、照射計画方法、およびこれを用いた荷電粒子照射システムであることを特徴とする。
【発明の効果】
【０００７】
この発明により、混合場のＲＢＥを短時間で精度よく予測して照射パラメータを決定できる照射計画装置、照射計画方法、および荷電粒子照射システムを提供できる。
【図面の簡単な説明】
【０００８】
【図１】荷電粒子照射システムのシステム構成図。
【図２】測定された細胞生存率と推定した細胞生存率を説明する説明図。
【図３】Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、Ｚ−ｎ，Ｄのイベント毎の線量平均比エネルギーを示すグラフの説明図。
【図４】２種類のビームについて深さの関数としてｄ、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、Ｚ−ｎ，Ｄの線量平均比エネルギーを説明する説明図。
【図５】２種類のビームについてある深さのｄ、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、Ｚ−ｎ，Ｄの線量平均比エネルギーを説明する説明図。
【図６】２つのビームが重ねられた場合の線量平均比エネルギーを深さの関数として表すグラフ。
【図７】飽和線量平均比エネルギーを求めてペンシルビーム線源データのビーム軸方向成分を求めるイメージを示す説明図。
【発明を実施するための形態】
【０００９】
以下、この発明の一実施形態について説明する。  
本発明者は、１回の治療で照射する放射線量を増大させ、治療期間を短縮するべく、炭素線やそれよりも重い粒子を用い、あるいは炭素線やそれよりも重い粒子に軽い粒子を組み合わせて用いた放射線治療を実現できるように鋭意研究を行った。そして、従来の炭素線の場合に用いていた演算よりも高精度で、かつ実際の治療の現場で利用できる演算速度で演算できる演算式を発明し、この演算式を利用した照射計画装置と、この照射計画装置により決定した照射パラメータを用いて粒子線を照射する粒子線照射システムを発明した。
【００１０】
以下、まず演算式について説明し、その後に照射計画装置および粒子線照射システムの実施例について説明する。  
＜＜演算式の説明＞＞
＜ＳＭＫモデル＞
ＳＭＫモデル（非特許文献１参照）では、細胞核をドメインと呼ばれる多数の微視的サイトに分けることができる。
【００１１】
細胞集団に巨視的な線量Ｄの電離放射線が照射されると、個々のドメインによる吸収線量である比エネルギーは、細胞集団全体にわたってドメインによって変化するランダム変数となる。電離放射線は、致死損傷および亜致死損傷の２つのタイプの損傷をドメイン内に引き起こすと推定される。
【００１２】
致死損傷は、その損傷が含まれるドメインに対して致死的である。亜致死損傷は生成時に非致死的であり、ドメイン内に生じた他の亜致死損傷と結合して致死損傷になるか、自発的に致死損傷に変化するか、あるいは自然治癒する可能性がある。ドメインの死は、ドメインを含む細胞の不活性化を引き起こす。各損傷の生成数は、［数１］で与えられるドメインの飽和比エネルギーｚ’ｄに比例する。
【数１】
2019180908.tif 000003
ここで、ｚ０は、高ＬＥＴ領域における線量当たりの複合ＤＮＡ損傷の数の減少を表す飽和パラメータである（ＨａｄａおよびＧｅｏｒｇａｋｉｌａｓ  ２００８  ※１）。
※１：Hada M and Georgakilas AG 2008 Formation of clustered DNA damage after high-LET irradiation: a review J. Radial. Res. 49 203-210
損傷の時間変化を考慮すると、ｚｄを受けるドメインの生残率ｓｄは、ドメイン内に致死損傷がない確率として決定され、その自然対数は［数２］により計算される。
【数２】
2019180908.tif 000004
ここで、ＡおよびＢは、放射線によって付与されたエネルギーから独立したパラメータである。細胞核がｐ個のドメインを含むと仮定すると、細胞核の比エネルギーｚｎを受ける細胞の生残率Ｓｎの自然対数は、［数３］により表すことができる。
【数３】
2019180908.tif 000005
ここで、α０およびβ０は、それぞれ、ＬＥＴ＝０の極限におけるＬＱモデルのαおよびβのパラメータを表すｐＡおよびｐＢである。ｆｄ（ｚｄ，ｚｎ）は、細胞核の比エネルギーｚｎを受ける細胞核内のドメインの比エネルギーｚｄの確率密度である。
【００１３】
細胞集団内のｚｎの確率的変動を考慮して、Ｄを照射された細胞の生残率の自然対数は、［数４］により計算される。
【数４】
2019180908.tif 000006
ここで、ｆｎ（ｚｎ，Ｄ）は、巨視的線量Ｄを照射された細胞集団内のｚｎの確率密度であり、次の［数５］および［数６］の関係を有する。
【数５】
2019180908.tif 000007
【数６】
2019180908.tif 000008
【００１４】
Ｚ−ｎ，Ｄは、細胞核に吸収されるイベント毎の線量平均比エネルギーであり、ｚｎおよび単一イベントの確率密度ｆｎ，１（ｚｎ）と共に［数７］により計算される。
【数７】
2019180908.tif 000009
【００１５】
非特許文献１では、与えられた照射条件におけるｆｄ（ｚｄ，ｚｎ）およびｆｎ（ｚｎ，Ｄ）を計算するための計算モデルを開発され、［数３］および［数４］を数値的に解いた。複数イベント照射のｆｄ（ｚｄ，ｚｎ）およびｆｎ（ｚｎ，Ｄ）の計算は、それぞれの照射条件についてモンテカルロ法による巨視的ビーム輸送シミュレーションを必要とし、さらに、ｚｄおよびｚｎの単一イベント確率密度ｆｄ，１（ｚｎ）およびのｆｎ，１（ｚｎ）の畳み込み積分を必要とするので、計算時間が膨大であった。
【００１６】
ｆｄ（ｚｄ，ｚｎ）の計算のために、非特許文献１では、さらに、ドメインへの複数の照射イベントによって引き起こされる飽和効果は無視できるほど小さいと仮定して、ｆｄ，１（ｚｄ）の畳み込み積分を省略した。すると、ｌｎＳｎは、［数８］のように簡略化できる。
【数８】
2019180908.tif 000010
【００１７】
イベント毎の頻度平均比エネルギーは、［数９］で表される。
【数９】
2019180908.tif 000011
【００１８】
イベント毎の飽和頻度平均比エネルギーは、［数１０］で表される。
【数１０】
2019180908.tif 000012
【００１９】
イベント毎の線量平均比エネルギーは、［数１１］で表される。
【数１１】
2019180908.tif 000013
【００２０】
イベント毎の飽和線量平均比エネルギーは、［数１２］で表される。
【数１２】
2019180908.tif 000014
【００２１】
次に、［数８］と［数４］は、細胞生残率を推定するために解かれるＳＭＫモデルの基本方程式である。非特許文献１は、ＳＭＫモデルを用いて細胞生存率を推定し、それらが大部分の荷電粒子線治療で使用される線量範囲（例えば、治療用炭素イオンビームであればＤが１０Ｇｙと同等かそれ以下の吸収線量）について［数３］および［数４］を直接解く場合ほぼ一致するものであることを見出した。
【００２２】
従来のＳＭＫモデルベースの計算では、ｆｎ（ｚｎ，Ｄ）の導出において、各治療照射についてのモンテカルロ法による巨視的ビーム輸送シミュレーションとｆｎ，１（ｚｎ）の畳み込み積分が必要である。したがって、患者内の各位置で生残率を予測する必要がある日々の治療計画にこのモデルを適応するには、依然として膨大な演算が要求される。更に、スキャニング荷電粒子線治療の治療計画にＳＭＫモデルを適合させるためには、ターゲットボリュームにわたるすべてのスポットのビームレットの各位置でのｆｎ，１（ｚｎ）をメモリ空間に格納しなければならず、それには膨大な記憶領域が必要になる。
【００２３】
次に、格納されたｆｎ，１（ｚｎ）は逐次近似演算の各反復においてフィールドの各位置でｆｎ（ｚｎ，Ｄ）を更新するために使用される。したがって、少なくとも商業化された治療計画システムに一般に使用されるコンピュータでは、スキャニング荷電粒子線治療の治療計画へのＳＭＫモデルの適応は、時間および記憶領域の両面から困難である。
【００２４】
前記の従来のSMKモデルでは、スキャニングビーム照射を前提とする粒子線治療装置に適用できないという問題に鑑み、発明者らは、以下説明する〈修正SMKモデル〉として新しいモデルを案出し、課題の解決を図った。
＜修正ＳＭＫモデル＞
荷電粒子療法では、ドメイン比エネルギーｚｄは、一般に、多数の低エネルギー付与イベントからなり、複雑なＤＮＡ損傷の飽和を誘発するイベントは稀である。言い換えれば、ｚｄ＞ｚｏとなるイベントはまれである。
【００２５】
従って、頻度平均比エネルギーについて、［数１３］の近似が成立する。
【数１３】
2019180908.tif 000015
【００２６】
そして、［数８］は［数１５］と［数１６］により［数１４］に変形できる。
【数１４】
2019180908.tif 000016
【数１５】
2019180908.tif 000017
【数１６】
2019180908.tif 000018
【００２７】
生残率の形で［数１４］を記述すると、［数１７］となる。
【数１７】
2019180908.tif 000019
【００２８】
細胞核へのエネルギー付与イベントの数が多いと仮定すると、ｚｎ＝Ｄ周辺のＳｎの二次テイラー展開は、［数１９］と［数２０］と［数２１］を用いて［数１８］により表せる。
【数１８】
2019180908.tif 000020
【数１９】
2019180908.tif 000021
【数２０】
2019180908.tif 000022
【数２１】
2019180908.tif 000023
【００２９】
［数４］に［数１８］を代入することにより、Ｄが照射された細胞の生残率は、［数５］と［数６］の関係を用いて、［数２２］として計算される。
【数２２】
2019180908.tif 000024
【００３０】
次いで、ｌｎＳは、［数１９］〜［数２１］を用いて［数２３］に変形でき、あるいは、細胞生残率の形式で［数２４］に変形することができる。
【数２３】
2019180908.tif 000025
【数２４】
2019180908.tif 000026
【００３１】
＜修正ＳＭＫモデルのパラメータの決定＞
修正ＳＭＫモデルでは、細胞生残率を推定するために必要なパラメータは、それぞれ、［数１５］のα０、［数１６］のβ０、［数１］の飽和パラメータｚ０、ドメインの半径ｒｄ、および細胞核の半径Ｒｎである。ｒｄはＺ−ｄ，ＤおよびＺ−＊ｄ，Ｄを導出するために使用され、ＲｎはＺ−ｎ，Ｄを導出するために使用される。これらのパラメータは、使用される細胞の種類やその状態のみに依存し、放射線の種類に依存しない。これらのパラメータのうち、半径Ｒｎは、光学顕微鏡で直接測定することができる。
【００３２】
なお、これらの「Ｚ−ｄ，Ｄ，Ｚ−ｎ，Ｄ，Ｚ−＊ｄ，Ｄ」について、本来であれば数１５および数１６に示すように平均値を示すバーを「ｚ」の上に表示すべきであるが、明細書本文中ではバーを表示できないため、「Ｚ−」と表示している。以下、本明細書の本文中の「Ｚ−」は全て数式において「ｚ」の上にバーを表示しているものを表す。
【００３３】
Ｒｎは、１６個のヒト細胞株のＲｎを測定した結果として、６．７μｍ〜９．５μｍの範囲であり、平均で８．１μｍであることが報告されている（Ｓｕｚｕｋｉ  ｅｔ  ａｌ  ２０００  ※２）。
※２：Suzuki M, Kase Y, Yamaguchi H, Kanai T and Ando K 2000 Relative biological effectiveness for cell-killing effect on various human cell lines irradiated with heavy-ion medical accelerator in chiba (HIMAC) carbon-ion beams Int. J Radiat. Oneal. Biol. Phys. 48 241-250
【００３４】
推定された細胞生存率は、Ｒｎに対して鈍感であるため（非特許文献１参照）、すべてのヒト細胞株についてＲ０＝８．１μｍと仮定し、ＳＭＫ計算のために決定されるパラメータの数を減らした。
他のパラメータの数値は、様々な照射条件下での推定された細胞生残率と測定された細胞生残率とを最小自乗回帰によって比較することによって決定した。
【００３５】
修正ＳＭＫモデルでは、細胞生存率は、ドメインに与えられたイベント毎の線量平均比エネルギーＺ−ｄ，ＤおよびＺ−＊ｄ，Ｄおよび細胞核に与えられたイベント毎の線量平均比エネルギーＺ−ｎ，Ｄから［数２４］を用いて計算される。比エネルギーは、公知の計算手順（Ｉｎａｎｉｗａ  ｅｔ  ａｌ  ２０１０  ※３）により導出する。
※３：Inaniwa T, Furukawa T, Kase Y, Matsufuji N, Toshito T, Matsumoto Y, Furusawa Y and Noda K 2010 Treatment planning for a scanned carbon ion beam with a modified microdosimetric kinetic model. Phys. Med. Biol. 55 6721-37
【００３６】
ドメインおよび細胞核の有感体積は、それぞれ半径ｒｄおよびＲｎ、長さ２ｒｄおよび２Ｒｎの水の円柱として仮定する。入射イオンは、均一な確率で有感体積の周りを横切り、それらの経路は常に円柱軸に平行である。イオンの軌道周りの放射状線量分布はＫｉｅｆｅｒ−Ｃｈａｔｔｅｒｊｅｅアモルファスイオントラック構造モデル（Ｃｈａｔｔｅｒｊｅｅ  ａｎｄ  Ｓｃｈａｅｆｅｒ  １９７６  ※４，Ｋｉｅｆｅｒ  ａｎｄ  Ｓｔｒａａｔｅｎ  １９８６  ※５，Ｋａｓｅ  ｅｔ  ａｌ  ２００６  ※６）によって記述されている。
※４：Chatterjee A and Schaefer H J 1976 Microdosimetric structure of heavy ion tracks in tissue Radial. Environ. Biophys. 13 215-227
※５：Kiefer J and Straaten H 1986 A model of ion track structure based on classical collision dynamics. Phys. Med. Biol. 311201-1209
※６：Kase Y, Kanai T, Matsumoto Y, Furusawa Y, Okamoto H, Asaba T, Sakama M and Shinoda H 2006 Microdosimetric measurements and estimation of human cell survival for heavy-ion beams Radial. Res. 166 629-638
【００３７】
測定された細胞生残率との比較のために、１２および１３５ＭｅＶ／ｕの入射エネルギーの３Ｈｅ、１２Ｃおよび２０Ｎｅイオンビームが異なる厚さのポリメチルメタクリレート板によって減速された実際の実験とは異なり、報告されたＬＥＴと３Ｈｅ、１２Ｃおよび２０Ｎｅイオンビームの単一エネルギースペクトルを仮定した（Ｆｕｒｕｓａｗａ  ｅｔ  ａｌ  ２０００  ※７）。
※７：Furusawa Y, Fukutsu K, Aoki M, ltsukaichi H, Eguchi-Kasai K, Ohara H, Yatagai F, Kanai T and Ando K 2000 Inactivation of aerobic and hypoxic cells from three different cell lines by accelerated He-, C- and Ne-ion beams Radial. Res. 154 485-96
【００３８】
この仮定は、低エネルギービームでは、プレートによる非弾性散乱やエネルギーストラグリングが、Ｚ−ｄ，ＤおよびＺ−＊ｄ，ＤおよびＺ−ｎ，Ｄに及ぼす影響が小さいため、妥当である。
【００３９】
最小自乗回帰では、Ｒｎ以外のＳＭＫパラメータを段階的に変更し、［数２５］で計算されたｌｏｇ１０Ｓの全二乗偏差が最小になる各パラメータ値を、それらの至適な値として決定する。
【数２５】
2019180908.tif 000027
ここで、Ｓｉ，ｅｘｐ、Ｓｉ，ｃａｌは測定され、ｉ番目の照射条件で推定された生残率を示し、各細胞株について回帰で使用された生残率データの数ｎはおおよそ３００であった。
【００４０】
公知のＭＫモデルのパラメータは、ＨＳＧおよびＶ７９細胞の同じインビトロ実験データを再現するために別個に決定した。ＭＫモデルパラメータの決定において、３Ｈｅ、１２Ｃおよび２０Ｎｅイオンビームの測定された細胞生存データが、データから得られた１０％生存線量ではなく、この実施形態で直接使用されたものであること以外は、Ｉｎａｎｉｗａ  ｅｔ  ａｌ  ２０１０※３の手順に従った。
【００４１】
＜修正ＳＭＫモデルに基づく生物線量の計算＞
非常に高い線量であっても基準放射線、すなわち２００−ｋＶのｘ線の生存曲線がＬＱモデルに従うと仮定すると、位置ｘでの生物線量は［数２６］で表される。
【数２６】
2019180908.tif 000028
【００４２】
ここで、αｘおよびβｘは、基準放射線のＬＱモデルの線形および二次係数である。修正ＳＭＫモデルでは、位置ｘにおける治療荷電粒子ビームの治療照射場（混合場）での生存率の自然対数、ｌｎＳ（Ｄ）、は、［数２７］で表される。
【数２７】
2019180908.tif 000029
ここで、Ｚ−ｄ，Ｄ，ｍｉｘおよびＺ−＊ｄ，Ｄ，ｍｉｘは、ドメインに吸収されるイベント当たりの線量平均比エネルギーおよび飽和線量平均比エネルギーであり、Ｚ−ｎ，Ｄ，ｍｉｘは、混合場の細胞核によって吸収されるイベント当たりの線量平均比エネルギーである。
【００４３】
スキャニング荷電粒子線治療の場合、ｘにおける吸収線量Ｄは、［数２８］で表される。
【数２８】
2019180908.tif 000030
ここで、ｄｊは、ｘに対するｊ番目のスポット（ｊ番目のビームレット）のスキャニングペンシルビームによって付与した線量であり、ｗｊは、ｊ番目のビームレットの入射イオンの数である。
【００４４】
混合場でのｘにおける比エネルギーは、以下のように記述される。
【数２９】
2019180908.tif 000031
【数３０】
2019180908.tif 000032
【数３１】
2019180908.tif 000033
ここで、Ｚ−ｄ，Ｄ，ｊおよびＺ−＊ｄ，Ｄ，ｊは、ｘに対するｊ番目のビームレットによって付与されたドメインのイベント毎の線量平均比エネルギーおよび飽和線量平均比エネルギーであり、Ｚ−ｄ，Ｄ，ｊは、ｘに対するｊ番目のビームレットによる細胞核のイベント毎の線量平均比エネルギーである。
【００４５】
スキャニング荷電粒子治療計画では、患者に対して望まれる生物線量分布を実現するために逐次近似繰り返し演算によってすべてのビームレットの粒子数ｗを決定しなければならない。粒子数ｗｊに対する生物線量の勾配の分析解（後述する［数３６］参照）は、逐次近似繰り返し演算アルゴリズム、例えば準ニュートン法、で用いられる。この詳細は、後述の補足説明に示す。
【００４６】
＜ペンシルビームデータ＞
スキャニング荷電粒子線治療計画において［数２６］に記載の前記生物線量分布を［数２８］〜［数３１］から得られたＤ、Ｚ−ｄ，Ｄ，ｍｉｘ、Ｚ−＊ｄ，Ｄ，ｍｉｘおよびＺ−ｎ，Ｄ，ｍｉｘによって計算するには、ｄｊ、Ｚ−ｄ，Ｄ，ｊ、Ｚ−＊ｄ，Ｄ，ｊおよびＺ−ｎ，Ｄ，ｊの質量を各ビームレットに関して計算しなければならない。
【００４７】
治療計画ソフトウェアでは、水中のペンシルビームのｄ、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，ＤおよびＺ−ｎ，Ｄの分布が予め決定され、ペンシルビームデータとして前記治療計画ソフトゥエアに登録される。前記データは、身体組織の水に対する阻止能比を用いた密度スケーリングにより患者線量計算に適用され、各治療計画の第ｊ番目のビームレットの量であるｄｊ、Ｚ−ｄ，Ｄ，ｊ、Ｚ−＊ｄ，Ｄ，ｊおよびＺ−ｎ，Ｄ，ｊが計算される。
【００４８】
前記水中でのペンシルビームのｄ、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，ＤおよびＺ−ｎ，Ｄの分布は、［数３２］、［数３３］、［数３４］、および［数３５］のそれぞれにおけるペンシルビームのｘにエネルギーを供給するすべてのイオントラックから求めることができる。
【数３２】
2019180908.tif 000034
【数３３】
2019180908.tif 000035
【数３４】
2019180908.tif 000036
【数３５】
2019180908.tif 000037
【００４９】
ここで、ｅｋはペンシルビームのｋ番目のトラックによってｘに与えられるエネルギーである。（Ｚ−ｄ，Ｄ）ｋおよび（Ｚ−＊ｄ，Ｄ）ｋは、それぞれｋ番目のトラックのｘにおけるドメインのイベント当たりの線量平均および飽和線量平均比エネルギーであり、上述したアモルファスイオントラック構造と円柱状のドメインを仮定する方法により求められる。
【００５０】
（Ｚ−ｎ，Ｄ）ｋは、ｋ番目のトラックによるｘにおける細胞核のイベント当たりの線量平均比エネルギーである。前記ｅｋ、（Ｚ−ｄ，Ｄ）ｋ、（Ｚ−＊ｄ，Ｄ）ｋ、および（Ｚ−ｎ，Ｄ）ｋの値は、各施設の照射システムを模したイオンビームのモンテカルロビーム輸送シミュレーションを用いて求めることができる。
【００５１】
＜ビーム輸送シミュレーション＞
荷電粒子療法では治療用ビームとして、４Ｈｅ、１２Ｃおよび２０Ｎｅイオンが、使用されているかまたは使用される計画がある。そこで、この数値計算研究では、治療用ビームのイオン種としてそれら3種が選択された。イオンペンシルビームのｅｋ、（Ｚ−ｄ，Ｄ）ｋ、（Ｚ−＊ｄ，Ｄ）ｋ、および（Ｚ−ｎ，Ｄ）ｋを導出するために、モンテカルロソフトウェアのＰＴＳｉｍを使用した（Ａｓｏ  ｅｔ  ａｌ  ２００７  ※８）。
※８：Aso T, Kimura A, Kameoka S, Murakami K, Sasaki T and Yamashita T 2007 Geant4 based simulation framework for particle therapy system. Nuclear Science Symposium IEEE, Conference Record 4 2564-2567
【００５２】
これはＧｅａｎｔ４ツールキット（バージョン９．２とパッチ０１）で構築された粒子治療のシミュレーションコードである（Ａｇｏｓｔｉｎｅｌｌｉ  ｅｔ  ａｌ  ２００３  ※９）。
※９：Agostinelli Set al. 2003 Geant4 - a simulation toolkit. Nucl. Jnstrum. Methods. Phys. Res. A 506 250-303
物理的相互作用のパッケージ、スキャニング照射システムのモデル、およびモンテカルロシミュレーションに用いられた水ファントムの幾何学的形状は、Ｉｎａｎｉｗａら（２０１７）※１０に記載されたものと同じである。
※１０：lnaniwa T, Kanematsu N, Noda Kand Kamada T 2017 Treatment planning of intensity modulated composite particle therapy with dose and linear energy transfer optimization Phys. Med. Biol. 62 5180-97
【００５３】
シミュレーションでは、１０ｍｍから３００ｍｍまでの範囲の２ｍｍステップの水中飛程を持つような初期エネルギーＥ０、すなわち１４６個のエネルギー、を有する３種のイオン種が、数分布が標準偏差２ｍｍの２Ｄ対称ガウス分布となるように走査磁石の直ぐ上流に生成された。生成されたイオンは、スキャニング照射システムを通過し、２００×２００×４００ｍｍ３の水ファントムに入射する。ファントム容積を、１．０×１．０×０．５ｍｍ３単位（「ボクセル」と呼ばれる）に分割して、シミュレーションしたイオンの様々な量、すなわち質量数Ａｐおよび原子番号Ｚｐで定義されるイオン種、ボクセル位置、イオンの運動エネルギーＴ、およびボクセルに与えられるエネルギーｅの空間分布を記録した。
【００５４】
ペンシルビームの線量分布ｄは、ｅの記録された分布から［数３２］によって単純に導出された。Ｚｐ、Ｔ、ｅの記録量からそれぞれ［数３３］、［数３４］、［数３５］のビームに対するＺ−ｄ，Ｄ、Ｚ−＋ｄ，Ｄ、およびＺ−ｎ，Ｄの分布を効率よく導出するために、Ｚｐ＝１−１０の単一エネルギーイオンのＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄをそれらの運動エネルギーＴの関数として表にした。
【００５５】
＜ビームモデリング分析＞
Ｉｎａｎｉｗａ  ａｎｄ  Ｋａｎｅｍａｔｓｕ（２０１５）※１１に記載されたフィッティング手順を、線量ｄと比エネルギーＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、Ｚ−ｎ，Ｄのシミュレーション分布に適用して、三重ガウストリクロームビームモデルを構築した。
※１１：lnaniwa T and Kanematsu N 2015 A trichrome beam model for biological dose calculation in scanned carbon-ion radiotherapy treatment planning Phys. Med. Biol. 60 437-451
【００５６】
ビームモデルにおいて、ビームの横断線量分布は３つのガウス分布の重ね合わせとして表され、３つのガウス成分に異なる比エネルギーが割り当てられて、ビーム断面にわたって比エネルギーの半径方向の変化を表した。１２Ｃおよび２０Ｎｅイオンビームにでは、一次イオン、原子番号Ｚｐ≧３の重いフラグメント、およびＺｐ≦２の軽いフラグメントのＺ−ｄ，Ｄ、Ｚ−＊ｄ，ＤおよびＺ−ｎ，Ｄは、第１、第２および第３のガウス成分がそれぞれ割り当てられた。４Ｈｅイオンビームに対する、一次イオンのＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、Ｚ−ｎ，Ｄは第１の成分に割り当てられ、他のフラグメントのそれらは第２と第３のガウス成分に割り当てられた（Ｉｎａｎｉｗａ  ｅｔ  ａｌ  ２０１７  ※１０）。
この分析ビームモデルでは、大角度散乱粒子の線量および比エネルギー分布に対する影響を、ペンシルビームアルゴリズムを用いた治療計画において考慮できる。
【００５７】
このようにして、ビームモデルを構築した後、臨床事例に適用したところ、従来法よりも大幅に短時間で演算して照射パラメータを決定でき、かつ、精度も高いことが確認された。よって、下記補足説明の後に、照射計画装置とこれを用いた粒子線照射システムの実施例について説明する。
【００５８】
＜補足説明＞
治療荷電粒子線の混合放射場における位置ｘにおける（∇Ｄｂｉｏ）Ｊの分析解［数３６］は、［数３７］として記載される。
【数３６】
2019180908.tif 000038
【数３７】
2019180908.tif 000039
ｗｊに対するｌｎＳ（Ｄ）の導関数は［数３８］、［数３９］、［数４０］、［数４１］および［数４２］によって与えられる。
【数３８】
2019180908.tif 000040
【数３９】
2019180908.tif 000041
【数４０】
2019180908.tif 000042
【数４１】
2019180908.tif 000043
【数４２】
2019180908.tif 000044
【００５９】
次に、この発明の一実施形態として、上述した演算式を用いた実施例を図面とともに説明する。
【実施例】
【００６０】
図１は、粒子線照射システム１のシステム構成を示すシステム構成図である。  
粒子線照射システム１は、イオン源２から照射された荷電粒子ビーム３を加速して出射する加速器４と、該加速器４から出射された荷電粒子ビーム３を輸送するビーム輸送系５と、該ビーム輸送系５を経た荷電粒子ビーム３を患者７の照射対象であるターゲット部８（例えば、腫瘍部）に照射する照射装置（スキャニング照射装置）６と、粒子線照射システム１を制御する制御装置１０と、粒子線照射システム１の照射パラメータを決定するコンピュータとしての照射計画装置２０とを備えている。なお、この実施例では、イオン源２から照射する荷電粒子ビーム３の核種（イオン種）として炭素とヘリウムの例で説明しているが、これに限らず、ネオンや酸素や陽子などを核種（イオン種）とする様々な荷電粒子ビームを照射する粒子線照射システム１に本発明を適用できる。また、この粒子線照射システム１は、スポットスキャン方式のものを用いているが、ラスタースキャン方式等の他のスキャニング照射方式のものとしてもよい。
【００６１】
加速器４は荷電粒子ビーム３の強度を調整する。  
照射装置６は、荷電粒子ビーム３をビーム進行方向（Ｚ方向）に垂直な平面を形成するＸ−Ｙ方向に偏向させるスキャニングマグネット（図示省略）と、荷電粒子ビーム３の位置を監視する線量モニタ（図示省略）と、Ｚ方向の荷電粒子ビーム３の停止位置を調整するレンジシフタ（図示省略）とを備え、ターゲット部８に対しスキャン軌道沿って荷電粒子ビーム３をスキャンする。
【００６２】
制御装置１０は、加速器４からの荷電粒子ビーム３の強度や、ビーム輸送系５内での荷電粒子ビーム３の位置修正や、照射装置６のスキャニングマグネット（図示省略）によるスキャニングや、レンジシフタ（図示省略）によるビーム停止位置等を制御する。
【００６３】
照射計画装置２０は、キーボードおよびマウス等で構成される入力装置２１、液晶ディスプレイまたはＣＲＴディスプレイ等で構成される表示装置２２、ＣＰＵおよびＲＯＭおよびＲＡＭで構成される制御装置２３、ＣＤ−ＲＯＭおよびＤＶＤ−ＲＯＭ等の記憶媒体２９に対するデータの読み書きを行うディスクドライブ等で構成される媒体処理装置２４、および、ハードディスク等で構成される記憶装置２５（記憶手段）を備えている。
【００６４】
制御装置２３は、記憶装置２５に記憶されている照射計画プログラム３９を読み込み、領域設定処理部３１、処方データ入力処理部３２、演算部３３（演算手段）、出力処理部３４、および３次元ＣＴ値データ取得部３６として機能する。
【００６５】
記憶部２５は、放射線核種毎（イオン種毎）に予め設定されたペンシルビーム線源データ４１を記憶している。ペンシルビーム線源データ４１は、ペンシルビームのビーム軸方向の情報として予め求められた水中での深部線量分布ｄ、ドメインと細胞核の線量平均比エネルギー（Ｚ−ｄ，Ｄ，Ｚ−ｎ，Ｄ）、およびドメインの飽和線量平均比エネルギー（Ｚ−＊ｄ，Ｄ）を有している。
【００６６】
このように構成された照射計画装置２０は、照射計画プログラム３９に従って、各機能部が次のように動作する。
【００６７】
まず、３次元ＣＴ値データ取得部３６は、別途のＣＴ装置から照射対象（患者）の３次元ＣＴ値データを取得する。
【００６８】
領域設定処理部３１は、表示装置２２に３次元ＣＴ値データを画像表示し、計画作成者が入力装置２１で入力する領域指定（ターゲット部８の指定）を受け付ける。
【００６９】
処方データ入力処理部３２は、表示装置２２に処方入力用画面を表示し、計画作成者が入力装置２１で入力する処方データを受け付ける。この処方データは、３次元ＣＴ値データの各座標における粒子線の照射位置とその照射位置で所望される生残率（またはそれと等価な臨床線量）、ビームの照射方向、粒子線の種類（核種）を示すデータである。また、その照射位置周辺に対する影響をできるだけ小さくするなどの種々の設定も処方データとして入力される。
【００７０】
演算部３３は、処方データ、およびペンシルビーム線源データ４１を受け取り、これらに基づいて照射パラメータおよび線量分布を作成する。すなわち、処方データの照射位置の照射対象（例えば癌細胞などの腫瘍）が処方データの生残率となる照射を行うために、粒子線照射システム１から照射すべき粒子線の種類と量（粒子数）をペンシルビーム線源データ４１を用いて逆算し、粒子線照射システム１から照射する粒子線の照射パラメータを算出する。この演算については後述する。
【００７１】
出力処理部３４は、算出した照射パラメータおよび線量分布を表示装置２２に出力して表示する。また、出力処理部３４は、照射パラメータおよび線量分布を、粒子線照射システム１を制御する制御装置１０に送信する。
【００７２】
次に、ペンシルビーム線源データ４１の作成方法について詳述する。  
ペンシルビーム線源データ４１は、モンテカルロシミュレーションで取得したｅ，Ｚｐ，Ｔの空間分布から、上述した［数３２］、［数３３］、［数３４］および［数３５］により決定されている。ただし、［数３３］、［数３４］および［数３５］の簡単化のために、核種毎にテーブルデータが作成され記憶装置２５に記憶されている。
【００７３】
核種毎のテーブルデータは、ドメイン線量平均比エネルギーＺ−ｄ，Ｄ、ドメイン飽和線量平均比エネルギーＺ−＊ｄ，Ｄ、および細胞核線量平均比エネルギーＺ−ｎ，Ｄのそれぞれついて、ある核種が、ある運動エネルギー（あるいは速度）で照射されると、どれだけのエネルギーを与えるかをテーブル化したものである。
このテーブルは、［数２４］を用いて、次のようにして予め作成されたものである。
【００７４】
すなわち、［数２４］を用いて細胞生残率Ｓ（Ｄ）を推定するために必要なパラメータは、それぞれ、［数１５］のα０、［数１６］のβ０、［数１］の飽和パラメータｚ０、ドメインの半径ｒｄ、および細胞核の半径Ｒｎである。ここで、細胞核の半径Ｒｎは光学顕微鏡で直接測定することができるため、測定によって与えられる。このため、決定しなければならない値は、α０、β０、飽和パラメータｚ０、ドメインの半径ｒｄである。なお、α０、β０、ｚ０、ｒｄは、細胞種によって決まるパラメータである。
【００７５】
そこで、粒子線照射システム１にて使用できる複数種類の核種について、核種毎に、当該核種を異なるエネルギーで複数回照射した際のそれぞれの細胞生残率を測定し、その測定値と、［数２４］を用いて計算した計算値の偏差が全ての核種（複数種類の核種）、全てのエネルギー（複数の異なるエネルギー）で最小になるように、最小二乗法等の適宜の最適化手法によって最適なα０、β０、飽和パラメータｚ０、ドメインの半径ｒｄを導出する。
【００７６】
図２は、このように測定値と計算値を修正ＳＭＫモデルにより近似させたグラフの説明図であり、図２（Ａ）〜図２（Ｄ）は、いずれも縦軸を生残率、横軸を照射線量とするグラフである。
【００７７】
図２（Ａ）及び図２（Ｂ）は、ある細胞種について、核種としてヘリウムイオンについての測定値（図中の黒点）と計算値（図中の実線）を示しており、図２（Ｃ）及び図２（Ｄ）は、ある細胞種（図２（Ａ）及び図２（Ｂ）と同じ細胞種）について、核種として炭素イオンについての測定値（図中の黒点）と計算値（図中の実線）を示している。また、図２（Ａ）はＬＥＴ＝１８．６ｋｅＶ／μｍ、図２（Ｂ）はＬＥＴ＝３３．００ｋｅＶ／μｍ、図２（Ｃ）はＬＥＴ＝２２．５ｋｅＶ／μｍ、図２（Ｄ）はＬＥＴ＝１３７ｋｅＶ／μｍで放射線を照射した場合を示している。
【００７８】
このように、ある細胞種に対して［数２４］にて求めた計算値（生残率）が、粒子線照射システム１にて使用する全ての種類の核種でＬＥＴの異なる全ての種類の放射線で線量を異ならせて照射した際の測定値（生残率）を最も良く再現する一組の（α０、β０、ｚ０、ｒｄ）が決定されている。
【００７９】
この図２は、幅広い範囲の線量およびＬＥＴにわたっての３Ｈｅおよび１２Ｃイオンビームに暴露されたＨＳＧおよびＶ７９細胞の測定された細胞生存率および推定した細胞生存率の例をそれぞれ示している。なお、図２（Ａ）〜図２（Ｄ）の点線は、従来のＭＫモデルにより演算した場合を示している。図示では差が小さく見えるが、例えば、炭素イオンを３３３ｋｅＶ／μｍで照射した場合や、ネオンイオンを６５４ｋｅＶ／μｍで照射した場合等には、顕著な差が現れる。
【００８０】
このようにして、細胞種毎に、核種やそのＬＥＴおよび線量を異ならせても測定値と計算値の偏差が最小となるようにα０、β０、ｚ０、ｒｄを求めることで、細胞種毎に［数２４］の変数であったα０、β０、ｚ０、ｒｄが固定パラメータとして定まる。従って、全ての核種について、照射位置での線量Ｄを決定すれば、［数２４］によりその照射位置での生残率Ｓ（Ｄ）を予測できるようになる。
【００８１】
なお、測定値は、生残率が０より大きく１より小さい範囲で線量を異ならせて複数回（４回以上が好ましく、５回以上がより好ましく、６回程度が好適である）測定して図示している。核種別に測定値を得る線量は、最も生残率が高い線量のときに生残率が０．１以上（好ましくは０．３以上）となり、最も生残率が低い線量のときに生残率が０．０５以下（好ましくは０．０３以下）となるようにして、その両端の範囲内でどちらか一端側に全て偏ることが無いようにばらつかせている。
【００８２】
また、図２では２種類のＬＥＴで示しているが、より広範囲のＬＥＴで測定値と計算値が最もよく再現されている一組の（α０、β０、ｚ０、ｒｄ）を決定することが好ましい。
【００８３】
次の［表１］は、ＨＳＧおよびＶ７９細胞の測定された細胞生存率を再現するために決定された修正ＳＭＫモデルでの固定パラメータの例を示している。なお、従来のＭＫモデルのパラメータもここに表示されている。それぞれの細胞株について、Ｒｎを除いて、修正ＳＭＫモデルとＭＫモデルとの間で同様のパラメータ値が決定されている。各細胞株およびモデルについて、データポイントあたりの［数２５］によって得られる平均二乗偏差も表１に示す。この修正ＳＭＫモデルパラメータ（固定パラメータ）として定められたα０、β０、ｒｄ、Ｒｎ、ｚ０は、ペンシルビーム線源データ４１の一部として記憶部２５に記憶されている。なお、表中のＸ２／ｎは、フィッティングの精度を示しており、固定パラメータには含まれていない。
【表１】
2019180908.tif 000045
【００８４】
次に、α０、β０、ｚ０、ｒｄが細胞種毎に定まったことを利用して、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄのそれぞれについて、粒子の運動エネルギー（もしくは粒子が照射される速度）に対してどれだけのエネルギーが与えられるかを演算してグラフ化し、このグラフを用いてテーブルデータとする。
【００８５】
図３は、そのグラフを説明する説明図であり、図３（Ａ）は核種毎のＺ−ｄ，Ｄのグラフ、図３（Ｂ）は核種毎のＺ−＊ｄ，Ｄのグラフ、図３（Ｃ）は各種別のＺ−ｎ，Ｄのグラフを示す。
【００８６】
図示するように、ペンシルビームの速度とペンシルビームから与えられるエネルギーの関係がグラフ化されており、これを用いてどの核種がどの速度で照射されるとどれだけのエネルギーが与えられるかをＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄ毎にテーブルデータとして登録する。なお、このテーブルデータは、表形式としてもよく、あるいは計算式にて都度計算できる方程式としても良い。
【００８７】
このようにして求めた各データを用いて、図４の説明図に示す積分線量分布ｄ（図４（Ａ）参照）、Ｚ−ｄ，Ｄ（図４（Ｂ）参照）、Ｚ−ｎ，Ｄ（図４（Ｃ）参照）およびＺ−＊ｄ，Ｄ（図４（Ｄ）参照）を求める。図４（Ａ）〜図４（Ｄ）は、いずれも、２種類のエネルギーのペンシルビームについて、それぞれ深さの関数として示されている。この図４（Ａ）〜図４（Ｄ）は、横軸が深さ、縦軸が線量平均比エネルギーを示すと共に、点線が第１のビームを示し、実線が第２のビームを示す。
【００８８】
図５は、上述した第１のビームｄ１と第２のビームｄ２について、ある深さ（Ｌ）でのそれぞれの線量平均比エネルギーを示している。すなわち、積分線量分布ｄにおける深さＬでの各ビームの線量ｄ１（Ｌ），ｄ２（Ｌ）（図５（Ａ）参照）、Ｚ−ｄ，Ｄにおける深さＬでの各線量平均比エネルギーＺ−ｄ，Ｄ１（Ｌ），Ｚ−ｄ，Ｄ２（Ｌ）（図５（Ｂ）参照）、Ｚ−ｎ，Ｄにおける深さＬでの各線量平均比エネルギーＺ−ｎ，Ｄ１（Ｌ），Ｚ−ｎ，Ｄ２（Ｌ）（図５（Ｃ）参照）、およびＺ−＊ｄ，Ｄにおける深さＬでの各線量平均比エネルギーＺ−＊ｄ，Ｄ１（Ｌ），Ｚ−＊ｄ，Ｄ２（Ｌ）（図５（Ｄ）参照）を示している。
【００８９】
このように、スキャニング照射法において重ね合わせるペンシルビームのＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄを深さの関数として予め計算し、ペンシルビーム線源データ４１として記憶している。
【００９０】
この第１のビームと第２のビームを重ね合わせると、図６に示すグラフとなる。
図６は、縦軸が線量、横軸が深さを示しており、第１のビームと第２のビームが重ねられた場合のグラフを示している。この図６から、次の［数４３］を求めることができる。
【数４３】
2019180908.tif 000046
【００９１】
すなわち、複数のペンシルビームにより着目位置Lに与えられるＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄを各々線量加重平均することで混合場のＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄを計算する。この混合場のＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄを用いて、上述した［数２４］により混合場の生物効果（生残率）やＲＢＥ（生物学的効果比）を決定する。
【００９２】
ここで、Ｚ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄは測定可能な物理量であるから、混合場の各位置でのＺ−ｄ，Ｄ、Ｚ−＊ｄ，Ｄ、およびＺ−ｎ，Ｄを測定することで、その位置でのＲＢＥ（生物学的効果比）を測定値からＳＭＫモデルの理論を用いて予測できる。
【００９３】
このようにして、混合場のＲＢＥを求めることが可能となる。このグラフデータを用いてペンシルビーム毎に決定した水中でのＺ−ｄ，Ｄ（図４（Ｂ）参照）、Ｚ−ｎ，Ｄ（図４（Ｃ）参照）およびＺ−＊ｄ，Ｄ（図４（Ｄ）参照）をペンシルビーム線源データ４１として登録している。
【００９４】
次に、演算部３３による演算について詳述する。  
まず、どの位置にどれくらいの生残率になるようにビームを照射したいのかを、オペレータに入力される。なお、一般には、生残率を直接入力するのではなく、それに相当する臨床線量を入力するように構成しており、入力された臨床線量を生残率に置換または変換して演算に使用する。このとき、ビームの方向および使用する核種の種類（複数）もオペレータが入力する。
【００９５】
演算部３３は、上述した［数２４］と、照射する核種の固定パラメータ（α０、β０、ｚ０、ｒｄ）を用いて、その核種をその位置にどれくらい照射するとどれくらいの生残率になるかを計算する。
【００９６】
演算部３３は、この計算をビームの線量や核種を変えながら繰り返し、入力された処方データに対してどの核種をどの線量でどの位置に照射すると最も良い結果が得られるかを演算する。そして、指定された核種を複数組み合わせて最適な照射パラメータを決定する。このように繰り返し演算をして最適な照射パラメータを決定する手法は、従来の適宜の手法を用いる。
【００９７】
照射計画装置２０は、このようにして決定した照射パラメータを制御装置１０へ送信し、制御装置１０が、この照射パラメータを用いて粒子線照射システム１による荷電粒子ビームの照射を実施する。
【００９８】
以上に説明した照射計画装置２０により、混合場のＲＢＥを、測定可能な物理量から決定し、細胞照射実験を行うことなく混合場の生物効果を予測できる。これにより、混合場のＲＢＥを短時間で精度よく予測して照射パラメータを決定できる。また、確率論的ＳＭＫモデルに基づいて、高ＬＥＴ・高線量の重粒子線治療場の生物効果や生物線量の計算、照射計画および治療計画を、計算時間の延長を招くことなく短時間に立案することができる。また、重粒子線治療場の生物効果をＳＭＫモデルの理論に基づいて測定値から確認できる。
【００９９】
また、図７の説明図に示すように、図７（Ａ）に示すスキャニング照射において重ね合わせるペンシルビームに関して、図７（Ｂ）に示すように各深さでの過剰殺傷効果を補正した比エネルギースペクトルからその深さでの飽和線量平均比エネルギーを求め、図７（Ｃ）に示すように積分線量分布とともにペンシルビーム線源データ４１のビーム軸方向成分として記憶部２５に登録しておくことができる。
【０１００】
また、この照射計画装置２０を用いた粒子線照射システム１は、高ＬＥＴ・高線量の照射場の生物効果を高精度の予測できるため、炭素だけでなく酸素やネオン等の高ＬＥＴ放射線を用いた大線量による短期照射を実現し、治療期間を短縮することができる。
【０１０１】
また、粒子線として重いもの（炭素や酸素やネオン等）によるペンシルビームだけでなく軽いもの（ヘリウムや陽子等）によるペンシルビームも混合して照射することもできる。この場合も短時間で高精度の予測をして適切な治療計画を立案することができる。
【０１０２】
また、逐次近似繰り返し演算の中で着目位置での各比エネルギースペクトルを導出する必要がないため、計算時間及びメモリ使用領域（制御装置２３のＲＡＭの使用領域）を大幅に軽減できる。従って、逐次近似繰り返し演算による照射計画の立案を実際の現場での使用に耐える時間内に実現することができる。
【０１０３】
なお、この発明は、上述の実施形態の構成のみに限定されるものではなく、多くの実施の形態を得ることができる。
【産業上の利用可能性】
【０１０４】
この発明は、粒子線を加速器で加速して照射する荷電粒子照射システムに利用することができ、特に、スポットスキャン方式およびラスタースキャン方式等のスキャニング照射方式の荷電粒子照射システムに利用することができる。
【符号の説明】
【０１０５】
１…粒子線照射システム
２０…照射計画装置
２５…記憶装置
３３…演算部
４１…ペンシルビーム線源データ
Ｚ−ｄ，Ｄ…ドメイン線量平均比エネルギー
Ｚ−ｎ，Ｄ…細胞核線量平均比エネルギー
Ｚ−＊ｄ，Ｄ…ドメイン飽和線量平均比エネルギー
【図１】
2019180908.tif 000048
【図２】
2019180908.tif 000049
【図３】
2019180908.tif 000050
【図４】
2019180908.tif 000051
【図５】
2019180908.tif 000052
【図６】
2019180908.tif 000053
【図７】
2019180908.tif 000054