(57)【特許請求の範囲】
【請求項１】
エネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップと、ｍ個のエネルギー源密度設定位置ｘi（ｉ＝１，．．．ｍ）に前記最小のエネルギー源密度ｑminを設定する第２のステップと、ｑi＋Δｑ＞ｑmaxとなるｉを除きｉ＝１からｍの各ｉについてそれぞれ、エネルギー源密度設定位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ増加させ、エネルギー分布ベクトルＵiを計算する第３のステップと、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶからパターン一致度γiを計算する第４のステップと、最も大きなパターン一致度を与えるエネルギー源密度設定位置ｘiのエネルギー源密度をｑi＋Δｑに変更する第５のステップと、すべてのエネルギー源密度設定位置ｘiのエネルギー源密度が最大のエネルギー源密度ｑmaxに到達するまで第３のステップから第５のステップまでを繰り返しながら、パターン一致度の最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップと、第６のステップで探索したエネルギー源密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの比ａを計算し、Ｐ／ａをエネルギー源密度分布とする第７のステップとを備えたエネルギー分布を形成する方法。
【請求項２】
エネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップと、ｍ個のエネルギー源密度設定位置ｘi（ｉ＝１，．．．ｍ）に前記最大のエネルギー源密度ｑmaxを設定する第２のステップと、ｑi−Δｑ＜ｑminとなるｉを除きｉ＝１からｍの各ｉについてそれぞれ、エネルギー源密度設定位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ減少させ、エネルギー分布ベクトルＵiを計算する第３のステップと、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶからパターン一致度γiを計算する第４のステップと、最も大きなパターン一致度を与えるエネルギー源密度設定位置ｘiのエネルギー源密度をｑi−Δｑに変更する第５のステップと、すべての位置ｘiのエネルギー源密度が最小のエネルギー源密度ｑminに到達するまで第３のステップから第５のステップまでを繰り返しながら、パターン一致度の最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップと、第６のステップで探索した密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの概略の比ａを計算し、Ｐ／ａをエネルギー源密度分布とする第７のステップとを備えたエネルギー分布を形成する方法。
【請求項３】
エネルギー源分布として、電荷密度分布・粒子線ビーム強度分布・電流密度分布・電圧源分布・電磁界源分布・放射線源分布・熱源分布・光源分布・荷重分布・音源分布・磁化分布のいずれかを選択することを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項４】
エネルギー分布として、電界分布・粒子線線量分布・電位分布・電磁界分布・応力分布・変位分布・温度分布・流速分布・音圧分布・放射線強度分布のいずれかを選択することを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項５】
パターン一致度として、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のコサイン（ｃｏｓ）成分を用いることを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項６】
パターン一致度として、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のサイン（ｓｉｎ）成分を用いることを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項７】
粒子線線量分布として陽子線線量分布、粒子線ビーム強度分布として陽子線ビーム強度分布を選択することを特徴とする請求項４記載のエネルギー分布を形成する方法。
【請求項８】
エネルギー源分布として、粒子線ビームの照射線量分布を選択し、エネルギー分布として、体内の粒子線線量分布を選択することを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項９】
粒子線として、陽子線を選択することを特徴とする請求項１記載のエネルギー分布を形成する方法。
【請求項１０】
ｉ＝１からｍに対する各△ｑの選び方は、△ｑj（ｊ＝１．．．ｎ）に対して
0003784419.tif 000005
となるように決めたことを特徴とする請求項１または請求項２記載のエネルギー分布を形成する方法。
【発明の詳細な説明】
技術分野
この発明は、対象システムで生成されるエネルギー分布が所望の分布に一致するように、前記システムへ付与するエネルギー源を決定し、エネルギー分布を形成する方法に関するものである。
背景技術
例えば、１９９４年発行、電気学会マグネティックス研究会資料ＭＡＧ−９４−２６の６３−６９頁に記載された従来のエネルギー分布を形成する方法は、サンプルド・パターン・マッチング法（Ｓａｍｐｌｅｄ  Ｐａｔｔｅｒｎ  Ｍａｔｃｈｉｎｇ法；以下、ＳＰＭ法と略称する）と呼ばれている。
これは、所望の場の分布（電界・磁界分布、放射線強度分布、温度分布など）が与えられたとき、その場を形成するエネルギー源の分布パターンを決定する手法である。具体的には、空間に格子点を設け、どの格子点に単位のエネルギー源を仮定すると、所望の場の分布に近いパターンが得られるかを逐次的に探索する手法である。パターン一致度が最大に達したときに、探索を終了する。
なお、パターン一致度は測定分布パターンベクトルＶと計算した分布パターンベクトルＵのなす角度θで評価する。具体的には、パターン一致度γは次式（１）で定義される。
γ＝ｃｏｓθ＝Ｕ・Ｖ／（｜Ｕ｜｜Ｖ｜）            ・・・（１）
ここでＵ・Ｖは内積、｜Ｕ｜，｜Ｖ｜はベクトルＵ，Ｖの大きさを表す。
図１（ａ）は目標とする磁界分布１０１を示すものであり、図１（ｂ）はＳＰＭ法で推定した前記目標磁界分布に近い分布を形成するための円弧状コイル１０２である。図１（ｃ）は図１（ｂ）の円弧状コイル１０２が発生する磁界分布１０３である。円の大きさが磁界強度に比例するように表示されている。これは、図１（ａ）の磁界パターンを近似的に再現しており、ＳＰＭ法の有効性を示している。ただし、コイルに流す電流値は一定でなければならない。これは、上記方法では単位の電流源を分布させるだけであり、原理的に電流の振幅を変化できないからである。
従来のエネルギー分布を形成する方法は以下のように構成されているので、強度を一定としたエネルギー源の位置は決定できたが、エネルギー源の強度を可変にできないという課題があった。
この発明は上記のような課題を解決するためになされたもので、エネルギー源の分布の強度を可変にできるエネルギー分布を形成する方法を得ることを目的とする。
発明の開示
この発明は、エネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップと、ｍ個のエネルギー源密度設定位置ｘi（ｉ＝１，．．．ｍ）に前記最小のエネルギー源密度ｑminを設定する第２のステップと、ｑi＋Δｑ＞ｑmaxとなるｉを除きｉ＝１からｍの各ｉについてそれぞれ、位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ増加させ、エネルギー分布ベクトルＵiを計算する第３のステップと、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶからパターン一致度γiを計算する第４のステップと、最も大きなパターン一致度を与える位置ｘiのエネルギー源密度をｑi＋Δｑに変更する第５のステップと、すべての位置ｘiのエネルギー源密度が最大のエネルギー源密度ｑmaxに到達するまで第３のステップから第５のステップまでを繰り返しながら、パターン一致度の最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップとを備え、第７のステップは第６のステップで探索した密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの比ａを計算し、Ｐ／ａをエネルギー源密度分布とするようにしたものである。このことによって、所望のエネルギー分布を与える入力エネルギーを効率よく得ることができる。
また、この発明は、エネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップと、ｍ個のエネルギー源密度設定位置ｘi（ｉ＝１，．．．ｍ）に前記最大のエネルギー源密度ｑmaxを設定する第２のステップと、ｑi−Δｑ＜ｑminとなるｉを除きｉ＝１からｍの各ｉについてそれぞれ、位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ減少させ、エネルギー分布ベクトルＵiを計算する第３のステップと、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶからパターン一致度γiを計算する第４のステップと、最も大きなパターン一致度を与える位置ｘiのエネルギー源密度をｑi−Δｑに変更する第５のステップと、すべての位置ｘiのエネルギー源密度が最小のエネルギー源密度ｑminに到達するまで第３のステップから第５のステップまでを繰り返しながら、パターン一致度の最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップとを備え、第７のステップ、第６のステップで探索した密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの比ａを計算し、Ｐ／ａをエネルギー源密度分布とするようにしたものである。このことによって、所望のエネルギー分布を与える入力エネルギーを効率よく得ることができる。
また、この発明は、エネルギー源分布として、電荷密度分布・粒子線ビーム強度分布・電流密度分布・電圧源分布・電磁界源分布・放射線源分布・熱源分布・光源分布・荷重分布・音源分布・磁化分布のいずれかを選択するようにしたものである。このことによって、エネルギー源の選択が容易である。
また、この発明は、エネルギー分布として、電界分布・粒子線線量分布・電位分布・電磁界分布・応力分布・変位分布・温度分布・流速分布・音圧分布・放射線強度分布のいずれかを選択するようにしたものである。このことによって、エネルギー分布の選択が容易である。
また、この発明は、パターン一致度として、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のコサイン（ｃｏｓ）成分を用いるようにしたものである。このことによって、所望のエネルギー分布を与える入力エネルギーを効率よく得ることができる。
また、この発明は、パターン一致度として、計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のサイン（ｓｉｎ）成分を用いるようにしたものである。このことによって、所望のエネルギー分布を与える入力エネルギーを効率よく得ることができる。
また、この発明は、粒子線線量分布として陽子線線量分布、粒子線ビーム強度分布として陽子線ビーム強度分布を選択するようにしたものである。このことによって、立ち上がり特性および立ち下がり特性が急峻なエネルギー強度分布を容易に得ることができる。
また、この発明は、エネルギー源分布として、粒子線ビームの照射線量分布を選択し、エネルギー分布として、体内の粒子線線量分布を選択するようにしたものである。このことによって、体内の放射線線量分布を望みの分布に設定できる。例えば、腫瘍内部の線量を一様にするとともに、周囲の正常組織に与える放射線線量を最小にできるので、腫瘍全体を効果的に治療でき、さらに副作用も低減できる。
また、この発明は、粒子線として、陽子線を選択するようにしたものである。このことによって、陽子線治療において、体内の放射線線量分布を望みの分布に設定できる。例えば、腫瘍内部の線量を一様にするとともに、周囲の正常組織に与える放射線線量を最小にできるので、腫瘍全体を効果的に治療でき、さらに副作用も低減できる。
【図面の簡単な説明】
第１図は従来のＳＰＭ法の説明図、第２図はこの発明の実施の形態１によるエネルギー分布を形成する方法を示す流れ図、第３図は空間におけるｍ個の電荷密度とｎ個の電界設定点の関係を示す説明図、第４図は電界設定位置対所望の電界強度の分布図、第５図は電界分布を決定する動作過程の説明図、第６図はこの発明の実施の形態１によるエネルギー分布を形成する方法を説明する電荷位置に対する電荷強度の分布図、第７図はこの発明の実施の形態１によるエネルギー分布を形成する方法を説明する位置に対する電界強度の分布図、第８図はこの発明の実施の形態２によるエネルギー分布を形成する方法を示す流れ図、第９図はこの発明の実施の形態２によるエネルギー分布を形成する方法による電荷位置に対する電荷強度の分布図、第１０図はこの発明の実施の形態２によるエネルギー分布を形成する方法による位置に対する電界強度の分布図、第１１図はこの発明の実施の形態３によるエネルギー分布を形成する方法を説明する図、第１２図は横方向の位置に対する放射線線量分布図、第１３図は人体表面の設定点に対する放射線線量分布図、第１４図は横方向の位置に対する陽子線ビーム強度分布図、第１５図は人体表面の設定点に対する放射線線量分布図、第１６図はこの発明の実施の形態４によるエネルギー分布を形成する方法を適用した陽子線治療装置を示す概要図、第１７図は体内線量の深さ方向分布図、第１８図は陽子線照射方法の説明図、第１９図は陽子線の照射線量分布を決定するための流れ図、第２０図はウェッジアブソーバの厚さに対する陽子線の照射線量を決定する過程の説明図、第２１図はウェッジアブソーバの厚さに対する陽子線の照射線量分布図、第２２図は所望の体内線量分布と計算で求められた体内線量分布を示す図、第２３図は計算で求められた体内線量分布と単一エネルギーの陽子線が与えられる体内線量分布を重ねて表示する図、第２４図は体内線量分布図、第２５図は体内線量分布の等高線表示図、第２６図はウェッジアブソーバの厚さに対する陽子線量を示す説明図である。
発明を実施するための最良の形態
以下、この発明をより詳細に説明するために、この発明を実施するための最良の形態について、添付の図面に従ってこれを説明する。
実施の形態１．
図２はこの発明の実施の形態１によるエネルギー分布を形成する方法を示す流れ図である。図において、ＳＴ１は決定するエネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップ、ＳＴ２はｍ個のエネルギー源密度設定位置（以下、位置と略称する）ｘi（ｉ＝１，．．．ｍ）に前記最小のエネルギー源密度ｑminを設定する第２のステップ、ＳＴ３はｉ＝１からｍの各ｉについてそれぞれ、位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ増加させ、エネルギー分布ベクトルＵiを計算する第３のステップである（ただし、ｑi＋Δｑ＞ｑmaxとなるｉを除く）。
ＳＴ４は計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のコサイン（ｃｏｓ）成分を計算し、パターン一致度γiとする第４のステップ、ＳＴ５は最も大きなパターン一致度γiを与える位置ｘiのエネルギー源密度をｑi＋Δｑに変更する第５のステップ、ＳＴ６はすべての位置ｘiのエネルギー源密度が最大のエネルギー源密度ｑmaxに到達するまで第３のステップＳＴ３から第５のステップＳＴ５を繰り返しながら、パターン一致度γiの最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップ、ＳＴ７は第６のステップＳＴ６で探索した密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの概略の比ａを計算し、Ｐ／ａをエネルギー源密度分布の最良推定値とする第７のステップである。
次に動作について図３〜図６を参照しながら説明する。
図３において、１０は電荷ｑl〜ｑmが分布する線状の空間である。ｍ個の電荷ｑl〜ｑmが存在する空間１０の上方の観測地点１１に、所望の電界強度分布を与える電荷分布を決定するｎ個の電界設定点での電界強度ｅl〜ｅnを考える。なお、ｘi〜ｘmはｍ個のエネルギー源密度設定位置（以下、位置と略称する）であり、電荷の数ｍ＝５０とする。
まず第１のステップＳＴ１において、決定する電荷ｑl〜ｑmの最小値ｑminと最大値ｑmaxを決める。例えば、正の電荷だけから所望の電界分布を形成したい場合は、ｑmin＝０とする。ｑmaxは任意の正の数でよいが、説明を簡単にするために１とする。つまり、まず相対的な振幅分布を求めて、後で比例定数を掛ければ振幅の絶対値を求め得る。もし、電荷ｑl〜ｑmの符号の制約がなければ、ｑmin＝−１、ｑmax＝１とすればよい。あるいは、負の電荷だけから形成したい場合は、ｑmax＝０とし、ｑminは任意の正の数でよいが、説明を簡単にするために−１とする。
次に第２のステップＳＴ２において、電荷ｑl〜ｑmの存在する空間１０に等間隔にｍ個の位置ｘi〜ｘmを定義し、各位置の電荷ｑl〜ｑmを仮にｑminとする。実施の形態１では、正の電荷だけで電界分布を形成することを考える。したがって、各位置ｘiの電荷を仮にゼロと設定する。すなわち、ｑmin＝０，ｑmax＝１とする。
然る後、第３のステップＳＴ３において、各位置ｘiの電荷を個別に最小単位Δｑ増加させながら、電界強度分布ベクトルを計算していく。まず、空間１０に電荷ｑi（ｉ＝１，．．．５０）を仮定したときの、各観測地点１１での電界強度ｅj（ｊ＝１，．．．１０）を理論計算して、システム行列Ａを定義する。システム行列のｊ行ｉ列成分Ａjiは次式（２）で与えられる。
Ａji＝１／（４πεｒji2）                        ・・・（２）
ここで、ｒjiは電荷ｑiの位置ｘiと電界強度ｅjを与える位置間の距離、つまり空間１０から観測地点１１までの距離である。この時、電界強度分布ベクトルＵ＝（ｅ1，ｅ2，．．．ｅ10）Tと電荷分布ベクトルＱ＝（ｑ1，ｑ2，．．．．ｑ50）Tは次式（３）で関係づけられる。Ｔはベクトルの転値を表す。
Ｕ＝ＡＱ                                          ・・・（３）
ここで、Ｕは１０行のベクトル、Ａは１０行５０列の行列、Ｑは５０行のベクトルである。
ステップＳＴ２では図５（ａ）に示すように、Ｑ2＝（０，０，．．．０）Tと初期設定されていることに注意する。また、第３のステップＳＴ３では、Ｑ1＝（Δｑ，０，０，．．．０）T，Ｑ＝（０，Δｑ，０，．．．０）T，Ｑ3＝（０，０，Δｑ，０，．．．０）T，．．．Ｑ50＝（０，．．．０，Δｑ）Tの５０通りの可能性のある選択肢を用意して、それぞれの電荷分布Ｑ1ないしＱ50に対して式（３）により電界強度分布ベクトルＵ1，．．．Ｕ50を次式（４）で計算する。
Ｕi＝ＡＱi    （ｉ＝１，．．．５０）            ・・・（４）
ここで、Δｑとしては（ｑmax−ｑmin）を適当な正の整数で割った値を選ぶ。ここでは、Δｑ＝１／５０とする。これはちょうど、強度情報を５０段階に量子化して決定することに相当している。従来のＳＰＭ法ではこの考慮がされていなかったために、強度を可変にできなかった。
次に第４のステップＳＴ４において、所望の電界強度分布ベクトルをＶ＝（ｅ1，ｅ2，．．．ｅ10）Tと定義し、第３のステップＳＴ３で計算した分布パターンベクトルＵiとのなす角度θiを評価する（ｉ＝１，．．．５０）。パターン一致度γiを次式（５）で計算する。
γi＝ｃｏｓθi＝Ｕi・Ｖ／（｜Ｕi｜｜Ｖ｜）    ・・・（５）
ここで，Ｕi・Ｖはベクトルの内積演算、｜Ｕi｜，｜Ｖ｜はベクトルＵi，Ｖの大きさを表す。
第５のステップＳＴ５において、最も大きなパターン一致度γiを与える位置ｘiのエネルギー源密度をΔｑに変更する。例えば、ｉ＝２とすれば、図５（ｂ）に示すように、新しい電荷分布ベクトルは、Ｑ2＝（０，Δｑ，０，．．．０）Tが選ばれるわけである。
第６のステップＳＴ６において、第５のステップＳＴ５で決められた新しい電荷分布ベクトルに対して、再び一カ所だけの電荷をΔｑ増加させた５０個の電荷の候補ベクトルＱ1＝（Δｑ，Δｑ，０，．．．０）T、Ｑ2＝（０，２Δｑ，０，．．．０）T、Ｑ3＝（０，Δｑ，Δｑ，０，．．．０）T．．．．Ｑ50＝（０，Δｑ，０，．．．０，Δｑ）Tを新たに考える。パターン一致度γiを第３，４のステップＳＴ３，ＳＴ４の計算により求める。最も大きなパターン一致度γiを与えるｉについて位置ｘiの電荷をΔｑだけ増加させる。
然る後、今度はｉ＝３とすれば、図５（ｃ）に示すように、新しい電荷分布ベクトルは、Ｑ3＝（０，Δｑ，Δｑ，０，．．．０）Tが選ばれるわけである。すべての電荷が最大値の１になるまでこれを繰り返していく。ただし、途中である位置ｘiにおける電荷が最大設定値１に到達した場合は、それ以上増加させることはモデルの仮定に反するので、その位置ｘiの電荷は固定しておき、他の位置の電荷密度だけに対して、図５（ｄ）、図５（ｅ）に示すように、探索プロセスを続行する。この計算プロセス中、パターン一致度の最大値を与える電荷分布Ｐを探索しておく。
第７のステップＳＴ７において、第６のステップＳＴ６で探索した電荷分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの概略の比ａを計算し、Ｐ／ａを電荷分布と決定する。これは、もともと電荷の最大値を１と仮定して計算したため、電荷分布の絶対値を求めるために必要な手段である。
図６は上記のように決定した電荷位置に対する電荷密度の分布図である。図７はこの電荷位置に対する電界強度分布図であり、この電界強度分布は前記図４に示す電界設定位置に対する所望の電界強度の分布図によく一致している。すなわち、所望の電界強度分布を形成できたことがわかる。
上記の一連のステップはこの発明を実施する一例に過ぎず、数学的に等価な操作に置き換えても同様な効果が得られる。例えば、第４のステップＳＴ４で計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度そのものを評価関数とし、第６のステップＳＴ６でこの角度を最小にするエネルギー源密度分布Ｐを探索してもよい。あるいは、第４のステップＳＴ４で計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のサイン（ｓｉｎ）成分を評価関数とし、第６のステップＳＴ６でこの評価関数を最小にするエネルギー源密度分布Ｐを探索してもよい。
実施の形態２．
図８はこの発明の実施の形態２によるエネルギー分布を形成する方法を示す流れ図であり、ＳＴ１は決定するエネルギー源密度の最小値ｑminと最大値ｑmaxを決める第１のステップ、ＳＴ２’はｍ個のエネルギー源密度設定位置ｘi（ｉ＝１，．．．ｍ）に前記最大のエネルギー源密度ｑminを設定する第２のステップ、ＳＴ３’はｉ＝１からｍの各ｉについてそれぞれ、位置ｘiのエネルギー源密度ｑiだけを所定の値Δｑ減少させ、エネルギー分布ベクトルＵiを計算する第３のステップである（ただし、ｑi−Δｑ＜ｑminとなるｉを除く）。
ＳＴ４は計算したエネルギー分布ベクトルＵiと所望のエネルギー分布ベクトルＶのなす角度のコサイン（ｃｏｓ）成分を計算し、パターン一致度γiとする第４のステップ、ＳＴ５’は最も大きなパターン一致度γiを与える位置ｘiのエネルギー源密度をｑi−Δｑに変更する第５のステップ、ＳＴ６’はすべての位置ｘiのエネルギー源密度が最小のエネルギー源密度ｑminに到達するまで第３のステップＳＴ３’から第５のステップＳＴ５’を繰り返しながら、パターン一致度γiの最大値を与えるエネルギー源密度分布Ｐを探索する第６のステップ、ＳＴ７’は第６のステップＳＴ６’で探索した密度分布Ｐを用いて計算したエネルギー分布ベクトルＵPと所望のエネルギー分布ベクトルＶの概略の比ａを計算し、Ｐ／ａをエネルギー源密度分布の最良推定値とする第７のステップである。
つまり、本実施の形態２が前記図２と機能的に異なるのは、すべての位置の電荷を最大値に設定して、その状態から第３のステップＳＴ３’において逐次Δｑだけ減少させていく探索手順をとることだけである。
したがって、図２に示す実施の形態１の手順を追加法と呼べば、図８に示す実施の形態２の手順は削除法と呼ぶことができる。図９はこの削除法で決定した電荷位置に対する電荷密度の分布図である。図１０は電界設定位置に対する電界強度の分布図であり、この電界強度分布は前記図４に示す電界設定位置に対する所望の電界強度の分布図によく一致している。すなわち、所望の電界強度分布を形成できたことがわかる。
実施の形態３．
図１１はこの発明の実施の形態３によるエネルギー分布を形成する方法を示す説明図であり、本実施の形態では、例えば、ペンシルビーム（陽子線ビーム）３１を被照射対としての人体３２に照射するとき、所望の放射線線量分布（エネルギー分布）を与えるようにペンシルビーム３１の照射強度分布（エネルギー源分布）を決定するものである。
所望の放射線線量分布として、ペンシルビーム３１に垂直な方向（横方向）の分布を図１２のように与える。すなわち、幅Ｗの範囲では一様分布であり、その外側ではゼロになるように設定する。いま、入射エネルギーＥのペンシルビーム３１を考える。ペンシルビーム３１の横方向分布（ｘ方向）をガウス分布とすれば、ペンシルビーム３１を横方向に平行スキャニングしたときの人体表面における線量分布Ｄ（ｘ，０）は次式（６）で与えられる。
0003784419.tif 000002
ここで、Ｉ（ｘ）は座標ｘにおけるペンシルビーム強度、σはペンシルビーム３１を媒質としての人体３２に入射したときの人体表面における１／ｅ散乱半径である。Ｓ（Ｅ，０）は人体表面における入射エネルギーＥのビームの阻止能である。
式（６）でＳ（Ｅ，０）＝１として離散化する。すなわち、線量分布ベクトルＵ、ペンシルビーム強度ベクトルＱ、システム行列に対する行列方程式Ｕ＝ＡＱを考える。このときのシステム行列Ａの成分は式（６）のガウス分布関数となる。たとえば、ペンシルビーム３１を横方向に１００ステップ走査して、人体３２に照射することを考え、人体表面における所望の線量分布ベクトルＶを図１３のように、２０点で設定したとする。
図２の流れ図において、エネルギー源密度として、今回は陽子線ビーム強度を考えればよい。陽子線ビーム強度は常に非負（０または正）であるから、図２のステップ１ではｑmin＝０，ｑmax＝１とする。第２のステップＳＴ２ないし第７のステップＳＴ７においては、式（６）で計算した線量分布Ｕと所望の線量分布Ｖのパターン一致度を評価する。その他は他の実施の形態と同一である。
図１４は本実施の形態３では決定された横方向の位置に対する陽子線ビーム強度分布図であり、このとき人体表面の設定点に形成される放射線線量分布は図１５の太線のようになる。これは図１２に示す所望の分布に近い分布となっている。図１５の細線は陽子線ビーム強度を一様に照射した場合の分布である。なお、陽子線ビーム強度分布として図１４を得たが、実際の放射線線量は陽子線ビーム強度と照射時間の積に比例するので、図１４を陽子線ビーム強度を一定としたときの陽子線ビーム照射時間のパターン図と考えることもできる。
一般に、陽子線治療では、線量の立ち下がりは急峻であることが望まれており、本実施の形態３によって決定された陽子線ビーム強度分布がより優れた急峻な立ち下がり特性をもつことがわかった。
以上の各実施の形態では、対象システムのエネルギー源分布が電荷分布または陽子線ビーム強度分布であったが、これに限るものではなく、陽子線以外の粒子線（炭素など）でもよい。
さらに、そのエネルギー源が対象システム内に形成するエネルギー分布（電界分布・磁界分布・磁束密度分布・応力分布・変位分布・温度分布・熱流速分布・音圧分布・放射線強度分布など）が所望の分布となるように、対応するエネルギー源分布（電流密度分布・磁化分布・熱源分布・音源分布・電圧源分布・電磁界源分布・光源分布・荷重分布・電波発信源分布・放射線源分布など）を本発明の方法で決定することができる。
なお、上記各実施の形態では、わかりやすさのために１次元のエネルギー源分布を決定したが、２次元または３次元のエネルギー源分布を決定することも同様に可能であることは言うまでもない。
実施の形態４．
図１６はこの発明のエネルギー分布を形成する方法を示す説明図であり、ＭｅｄｉｃａｌＰｈｙｓｉｃｓ第１０巻、３４４頁（１９８３年）に記載された従来の陽子線治療装置を改変した装置、つまり人体内部の腫瘍などのターゲットに陽子線を照射して治療する装置を示している。
図１６において、４１は陽子線４２が平行ビームとして入力されるウェッジアブソーバ、４３は矢印方向に移動して陽子線４２の通過範囲を決定する多葉コリメータ、４４は陽子線のエネルギーを変化させるボーラス、４５は照射線量を計測する線量モニタ、４６は線量モニタ４５の計測量を入力する計算機、４７は人体の体表、４８は人体内部の腫瘍などのターゲットである。
次に動作について説明する。
陽子線４２は上方から一様な平行ビームで入射し、ウェッジアブソーバ４１により陽子線４２のエネルギーを減衰させる。多葉コリメータ４３は腫瘍の幅にあわせて移動させる。ボーラス４４により腫瘍形状にあわせて陽子線４２のエネルギーを変化させる。ウェッジアブソーバ４１およびボーラス４４はアクリルなど人体とほぼ同じ陽子線減衰特性を有する材料で構成する。
線量モニタ４５により照射線量を計測して計算機４６に取り込み、ウェッジアブソーバ４１および多葉コリメータ４３を移動するタイミングを制御する。ターゲットの最も深い部分から順に仮想的なスライス１、２、３、・・・を作り、スライスごとに順に照射して各スライス部の線量分布を形成していく。つまり、線量モニタ４５で各スライスに対する照射線量を測定し、あらかじめ計算したパターンに一致するように制御しながら順次照射していくわけである。なお、照射線量とは、照射する放射線強度の時間積分に比例する。陽子線を生成する加速器がサイクロトロンの場合は、照射ビーム強度（放射線強度）が一定なので、スライスごとの照射時間を変化させて、照射線量を制御することが多い。加速器がシンクロトロンの場合は、照射ビーム強度は時間的に一定ではないので、照射線量を線量モニタで測定することにより、照射線量を制御することが多い。陽子線４２を人体に照射するとき、所望の体内放射線線量分布（エネルギー分布）をあたえるように陽子線４２の照射線量分布（エネルギー源分布）を決定することを考える。
体内線量の深さ方向分布は一般に図１７で表される。図１７は前記文献記載の図をわかりやすく追加改変したものである。図において、スライス１、スライス２、スライス３に対する照射線量（エネルギー源分布）はそれぞれｑ１，ｑ２，ｑ３で表現されている。また、この図は人体にほぼ等価な水モデルを用いた実験データに基づいている。点線は従来の方法で最適化された体内線量分布を表すもので、平坦化領域における平坦度が不十分であった。
図１８にさらに詳細な照射方法を示すもので、腫瘍などのターゲット４８の後部形状に沿ったスライスの枚数は図のターゲット４８の代表的長さＬをスライス１枚の厚さで割って求める。例えば、Ｌ＝６ｍｍでスライス厚さが２ｍｍであれば、スライス枚数は３枚である。
図１９は陽子線４２の照射線量分布（エネルギー源分布）を決定するための流れ図である。目標とする体内線量分布（エネルギー分布）を与える陽子線４２の照射線量分布を決定する流れ図となっているところを除いて図２と同じである。
図２０は各ウェッジアブソーバ４１の厚さに対する陽子ビームの照射線量を決定する過程を示す図である。図１９のｉ＝１，ｍがウェッジアブソーバ４１の厚さと対応する番号になっている。ここでは、ウェッジアブソーバ４１の厚さは２ｍｍずつ増加させていく場合を考えるので、ｉ＝１はウェッジアブソーバ４１の厚さが２ｍｍ，ｉ＝２５がウェッジアブソーバ４１の厚さが５０ｍｍに対応することになる。その他は図５と同じである。
図２１は図１９、図２０のアルゴリズムにより計算されたウェッジアブソーバ４１の厚さに対する照射線量分布である．図２１はパターン一致度の最大値を与える解を示しており、図１９のステップＳＴ１６の照射線量分布Ｐに相当する。
図２２は目標とした所望の体内線量分布Ｖと本発明の計算で求められた照射線量分布で照射した場合の体内線量分布Ｕｐを示す。横軸は体内深さ（ただし水で置き換えた水中深さ）である。所望の体内線量分布を与えていることがわかる。
図２３はこの発明の計算で求められた照射線量分布で照射した場合の体内線量分布Ｕｐと単一エネルギーの陽子線４２が与える体内線量分布を重ねて表示したものである。今回計算した照射線量分布によって、目標領域より深い部位の線量の減衰特性は単一エネルギーの陽子線４２の減衰特性と同一になっていることがわかった。これは、腫瘍の後部に放射線を照射してはいけない脊髄などの重要臓器がある場合に必要な特性である。この発明により形成された体内線量分布は所定の範囲（腫瘍の位置）で体内線量分布の平坦度をほぼ２％以下にでき、さらに腫瘍後部の線量減衰も理想的な特性となった。単一エネルギーの陽子線４２の減衰特性を越える減衰速度を得ることは物理的に不可能であるから、この発明の結果は理論限界を与えている。
この発明で得られた照射線量分布を用いて、図１８の腫瘍が水中深く存在する場合の体内線量分布を図２４に示す。図１６において、照射はもっとも深いスライス１から順に行い、ウェッジアブソーバ４１と多葉コリメータ４３をスライスごとに動かす。腫瘍部位では、一様な体内線量が形成されていることがわかる。
図２５は体内線量分布を等高線表示したものであり、腫瘍部位では１００％体内線量を形成していることがわかり、深さ方向と横方向に対する体内線量の減衰特性は非常に急峻な理想的なものになっていることがわかる。
なお、ここでは粒子線として陽子線を例にして説明したが、ヘリウムビーム、炭素ビーム、ネオンビームを初めとする他の粒子線でも同じである。
さらに、この実施の形態では、説明を簡単にするために、図１６で示された水モデルによる体内線量分布を図１８の深さ方向の線量分布Ｕiとして採用したが、実際の治療において患者のＸ線ＣＴ画像データに基づいて体内線量分布を高精度に計算して図１８の深さ方向の線量分布Ｕiとして用いることが望ましい。患者のＸ線ＣＴ画像データに基づいて体内線量分布を高精度に計算する方法にはモンテカルロ法を用いる方法（１９７８年に雑誌Ｒａｄｉａｔｉｏｎ  Ｒｅｓｅａｒｃｈ  第７４巻２１７〜２３０ページにゴイテンらが発表）や無限スラブモデル（１９９５年に雑誌Ｐｈｙｓｉｃｓ  ｉｎ  Ｍｅｄｉｃｉｎｅ  ａｎｄ  Ｂｉｏｌｏｇｙ  第４０巻１０３１〜１０４３ページにラッセルらが発表）を用いる方法などが公知であり、これらにより深さ方向の線量分布Ｕiを計算して、この発明に適用すればよい。
さらに、上記は深さ方向の線量分布Ｕiだけを変化させたが、人体内の腫瘍を細かな体積要素に分割して、各体積要素ごとに照射線量を変化させる治療も行われている（１９８３年に雑誌Ｎｕｃｌｅａｒ  Ｉｎｓｔｒｕｍｅｎｔｓ  ａｎｄ  Ｍｅｔｈｏｄ  第２１４巻４９１〜４９６ページに金井らが発表）。
この場合、図１６のシステムに直径６ｍｍ程度の陽子線を入射させて、図示したｘ方向および紙面に垂直な図示しないｙ方向の２方向に前記直径６ｍｍ程度の陽子線をスキャニングする磁石が設置される。人体に照射される陽子線エネルギーは図１６のウェッジアブソーバまたはこれと同様な機能をもつエネルギーディグレーダにより調節されるか、あるいは陽子線を生成加速する加速器内で調節される。
この結果、前記陽子線は３次元的にスキャニングされる。このような治療装置に図１８で開示されるこの発明を適用すれば、細かな体積要素ごとに体内線量を制御でき、その線量分布が望みの３次元分布に設定できることは言うまでもない。
なお、図１９では陽子ビームの照射線量のステップ幅△ｑを常に一定値の１／５０として説明したが、図２６のように△ｑ1、△ｑ2、．．．、△ｑnのうち少なくとも１つを他と異なる値に設定してもよく、
0003784419.tif 000003
とすれば自由に選ぶことができる。装置制約上、照射線量の最小設定値があまり小さくできない場合に有効である。
図２６では照射線量の最小設定値を０．１として、そのあとは０．０１きざみで最大値１．０まで増加させる例を示している。このように設定すれば、装置の制約上０．１以下が設定できないときにおいても、この制約下で準最適な線量で照射制御ができる。  つまり、
0003784419.tif 000004
の範囲内できざみ幅△ｑは自由に決めることができる。なお、図２６では、ｑmax＝１、ｑmin＝０としている。
産業上の利用可能性
以上のように、この発明にかかるエネルギー分布を形成する方法は、そのエネルギー分布の強度を可変できるようにしたもので、人体内部の腫瘍などのターゲットにエネルギー分布を照射して治療する場合、その腫瘍の後部にある重要臓器にエネルギー分布を照射しないで治療するためのエネルギー分布を形成する方法として有用である。
【図１】
0003784419.tif 000006
【図２】
0003784419.tif 000007
【図３】
0003784419.tif 000008
【図４】
0003784419.tif 000009
【図５】
0003784419.tif 000010
【図６】
0003784419.tif 000011
【図７】
0003784419.tif 000012
【図８】
0003784419.tif 000013
【図９】
0003784419.tif 000014
【図１０】
0003784419.tif 000015
【図１１】
0003784419.tif 000016
【図１２】
0003784419.tif 000017
【図１３】
0003784419.tif 000018
【図１４】
0003784419.tif 000019
【図１５】
0003784419.tif 000020
【図１６】
0003784419.tif 000021
【図１７】
0003784419.tif 000022
【図１８】
0003784419.tif 000023
【図１９】
0003784419.tif 000024
【図２０】
0003784419.tif 000025
【図２１】
0003784419.tif 000026
【図２２】
0003784419.tif 000027
【図２３】
0003784419.tif 000028
【図２４】
0003784419.tif 000029
【図２５】
0003784419.tif 000030
【図２６】
0003784419.tif 000031