(57)【要約】  画像処理システムＩＰＳ及び関連する方法。画像処理システムＩＰＳは、撮像像置ＩＡによって取得される対象ＯＢの画像ＩＭを受け取る入力インタフェースＩＮを有する。画像処理システムＩＰＳのカーネルプロバイダＫＰは、少なくとも２つの散乱タイプについて個々の散乱カーネルを提供するように構成される。画像処理システムＩＰＳの散乱補正モジュールＳＣＭは、提供された少なくとも２つのカーネルに基づいて画像内で補正を実行するように構成される。
2021525141.tif 000002
【特許請求の範囲】
【請求項１】
撮像装置によって取得される物体の画像を受け取る入力インタフェースと、
少なくとも２つの散乱タイプについてそれぞれの散乱カーネルを提供するカーネルプロバイダと、
提供される少なくとも２つの前記散乱カーネルに基づいて前記画像内で補正を行う散乱補正モジュールであって、前記少なくとも２つの散乱タイプが、それぞれ異なる数の散乱イベントを含み、前記散乱タイプのうち一方の散乱タイプが、単一散乱イベントを含み、他方の散乱タイプが多重散乱イベントを含み、前記カーネルプロバイダが、前記少なくとも２つの散乱カーネルを異なる方法で選択し又は計算するように構成され、前記カーネルプロバイダは、所定の閾値で勾配値をクリッピングすることによって、前記画像の所与のパッチ内のクリッピングされた信号勾配に基づいて、多重散乱イベント用の少なくとも１つの散乱カーネルを提供する、散乱補正モジュールと、
を有する画像処理システム。
【請求項２】
前記カーネルプロバイダは、前記画像内の所与のパッチ内の信号強度及び／又は信号勾配に基づいて、単一散乱イベント用の少なくとも１つのカーネルを選択し又は計算する、請求項１に記載の画像処理システム。
【請求項３】
前記カーネルプロバイダは、前記画像を最初にローパスフィルタリングし、次いで、前記フィルタリングされた画像の所与のパッチにおける信号強度又は勾配を考慮することに基づいて、多重散乱イベント用の前記少なくとも１つのカーネルを提供する、請求項１又は２に記載の画像処理システム。
【請求項４】
前記カーネルプロバイダは、前記画像全体の信号強度及び／又は信号勾配に基づいて、多重散乱イベント用の前記少なくとも１つのカーネルを提供する、請求項１乃至３のいずれか１項に記載の画像処理システム。
【請求項５】
前記カーネルプロバイダは、前記対象の物理的特性に基づいて、前記散乱カーネルのうちの少なくとも１つを提供する、請求項１乃至４のいずれか１項に記載の画像処理システム。
【請求項６】
前記カーネルプロバイダは、前記対象の少なくとも部分形状に基づいて、多重散乱イベント用の前記少なくとも１つのカーネルを提供する、請求項５に記載の画像処理システム。
【請求項７】
前記カーネルプロバイダは、予めに計算されたカーネルのライブラリに基づいて、前記散乱カーネルの少なくとも１つを提供する、請求項１乃至６のいずれか１項に記載の画像処理システム。
【請求項８】
前記画像に基づいて組織タイプを決定する組織決定器を有し、前記カーネルプロバイダは、前記少なくとも２つのカーネルのうちの少なくとも１つを提供する際に、前記組織タイプを考慮に入れる、請求項１乃至７のいずれか１項に記載の画像処理システム。
【請求項９】
撮像装置により取得される対象の画像を受け取るステップと、
カーネルプロバイダにより、少なくとも２つの散乱タイプについてそれぞれの散乱カーネルを提供するステップであって、前記少なくとも２つの散乱タイプは、それぞれ異なる数の散乱イベントを含み、一方の散乱タイプは、単一散乱イベントを含み、他方の散乱タイプは、多重散乱イベントを含む、ステップと、
前記少なくとも２つの散乱カーネルを異なる方法で選択し又は計算するステップと、
予め規定された閾値で勾配値をクリッピングすることによって、前記画像の所定のパッチ内のクリッピングされた信号勾配に基づいて、多重散乱イベント用の少なくとも１つの散乱カーネルを提供するステップと、
提供された前記少なくとも２つの散乱カーネルに基づいて、前記画像内で補正を実施するステップと、
を有する方法。
【請求項１０】
少なくとも１つの処理ユニットによって実行されるとき、前記処理ユニットに請求項９に記載の方法を実行させるコンピュータプログラム。
【請求項１１】
請求項１０に記載のコンピュータプログラムを記憶したコンピュータ可読媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理システム、画像処理方法、撮像装置、コンピュータプログラム、及びコンピュータ可読媒体に関する。
【背景技術】
【０００２】
Ｘ線イメージングは医療分野における診断又は介入目的のための最も重要な撮像モダリティであるが、Ｘ線イメージングは、材料試験などの他の分野でも使用されている。
【０００３】
特に診断医に利用可能な高コントラストＸ線画像を有することは、診断エラーレートを低減することを可能にする。散乱はＸ線撮像の最中に発生する現象であり、画像品質、特にＸ線コントラストを低下させることが知られている。散乱による損傷は、決定的な画像の詳細を隠す可能性があり、従って、例えば微細な微小石灰化が検査されることもあるマンモグラフィにおいて誤診につながることがある。
【発明の概要】
【発明が解決しようとする課題】
【０００４】
散乱の影響に対処するために、さまざまな機構が過去に提案されてきた。そのようなアプローチの１つは、散乱放射線の少なくとも一部を遮るように設計された散乱防止グリッド（ＡＳＧ）を使用することである。他のアプローチは画像処理に基づいており、そのいくつかはカーネルベースである。これらのアプローチでは、散乱の影響を近似するために特定の関数が使用される。これらの方法では、ある基本的な形状、例えば異なるサイズの均質な水球の散乱によるペンシルビームの散乱放射線が予め決定される。次いで、実際のＸ線画像における散乱は、例えば信号強度及び信号勾配に基づいて、局所的にベストマッチするボールを各ビームセクションごとに見つけ、そのボールによって生成される散乱を、一次ビームによってその位置に生成される実散乱の推定値として使用することによって、推定される。水球に基づくカーネル散乱モデリングは、出願人の国際公開第２００７／１４８２６３号パンフレットに記載されている。
【０００５】
当技術分野では、Ｘ線散乱補正を更に洗練し又は改善する必要がありうる。
【課題を解決するための手段】
【０００６】
本発明の目的は、更なる実施形態が従属請求項に組み込まれる独立請求項の主題によって解決される。以下に説明する本発明の態様は、画像処理方法、コンピュータプログラム要素、及びコンピュータ可読媒体に等しく適用されることに留意されたい。
【０００７】
本発明の第１の態様によれば、画像処理システム（ＩＰＳ）であって、撮像装置によって取得された対象の画像を受信するための入力インタフェースと、少なくとも２つの散乱タイプに対してそれぞれの散乱カーネルを提供するように構成されるカーネルプロバイダと、提供された少なくとも２つのカーネルに基づいて画像内で補正を実行するように構成される散乱補正モジュールと、を有する画像処理システムが提供される。
【０００８】
カーネルは、実際の検査対象又は概念的な検査対象（「ファントム」とも呼ばれる）から導かれる明示的又は暗黙的な関数である。カーネルは、ファントムに関連付けられていると言える。カーネルは、点拡がり関数（ＰＳＦ）である。ＰＳＦは一般に、ファントムが検出器のある画像位置の周りにどのようにペンシルビームを拡げるかを記述する２Ｄスカラ場である。カーネルは、ファントムの散乱挙動を記述する。撮像される対象（ファントムとは異なる）が少なくとも局所的にファントムに対応する場合、同様の点拡がり又は散乱挙動を予想することができると仮定される。カーネルは、２Ｄスカラ場であり、画像として考えられることができ、かかる画像は、実際の対象の画像からのものと区別するために、本明細書では「散乱画像」と呼ぶことができる。
【０００９】
少なくとも２つの散乱タイプは、それぞれ異なる数の散乱イベントを含む。好ましい実施形態では、１つの散乱タイプは、単一の散乱イベントを有し、他の散乱タイプは、多重散乱イベントを有する。
【００１０】
単一散乱イベントの推定は、一次ビームに近いところの撮像対象の特性にほぼ専ら依存するが、多重散乱の推定は、更に離れたところの特性にも依存することが分かった。従って、ここでは、散乱カーネルを単一散乱カーネルと多重散乱カーネルとに分けることを提案する。これは、単一散乱カーネルの選択及び多重散乱カーネルの選択のために異なる対象特性を使用することができるという利点を提供する。従って、組み合わされた散乱イベントを最適化するのではなく、散乱のタイプをより良く最適化することができる。言い換えると、異なる散乱タイプ、特に単一散乱と多重散乱とを区別する、より洗練された散乱モデリングが提案され、それにより、そのような区別がなされないアプローチと比べてより良い散乱補正が与えられる。
【００１１】
カーネルプロバイダによる提供は、入力画像内の画像情報（又は勾配、フィルタリングされた又はクリップされた信号又は勾配などの、それに関連する情報）に対して、予め計算されたカーネルの２つの異なる組に関する関連ファントム同士をマッチングすることを含むことができる。２つの異なる組は、単一散乱イベントをシミュレートすることによって得られるカーネルと、多重散乱イベントをシミュレートすることによって得られるカーネルとをそれぞれ含む。各々の組において、それぞれのカーネルは、関連するファントム又は材料のタイプのサイズ／形状に従って、更に階層化されることができる。好ましくは、簡潔にするために、水ベースのファントムが使用される。
【００１２】
カーネルプロバイダは、２つの散乱タイプの異なる特性をより適切に考慮するために、少なくとも２つのカーネルを異なる方法で選択し又は計算するように構成される。従って、少なくとも２つのカーネルのうち少なくとも１つのカーネルが異なる方法で選択されることができる。カーネルプロバイダが少なくとも２つのカーネルのうち少なくとも１つのカーネルを異なる方法で選択し又は計算するように構成される方法の例及び実施形態が、以下でより詳しく説明される。
【００１３】
単一散乱タイプのカーネルの選択又は計算は、好ましくは局所的（local）なアプローチに従い、一方、多重散乱タイプのカーネルの選択は、局所的なアプローチ（ただし、単一散乱の局所的アプローチとは異なる）又は大域的（global）なアプローチのいずれかに従う。
【００１４】
具体的には、実施形態ではカーネルプロバイダは、画像内の所与のパッチ（又はサブセット）内の信号強度及び／又は信号勾配に基づいて、単一散乱イベント用の少なくとも１つのカーネルを選択し又は計算するように構成される。信号強度又は信号勾配に加えて、又はその代わりに、フィルタリングされた勾配、又はより大きなピクセル近傍における信号値又は信号値の関数など、入力画像内の他の画像特徴を使用して、カーネルプロバイダによるカーネル選択をガイドすることが可能である。
【００１５】
実施形態では、カーネルプロバイダが画像の所与のパッチ内のクリップされた信号勾配に基づいて、多重散乱イベント用の少なくとも１つのカーネルを提供するように構成される。
【００１６】
実施形態では、カーネルプロバイダは、最初に画像をフィルタリングし、次いでフィルタリングされた画像の所与のパッチ内の信号強度及び／又は勾配を考慮することに基づいて、多重散乱イベント用の少なくとも１つのカーネルを提供するように構成される。フィルタリングは、ローパスフィルタリングを含む。フィルタリング又は勾配クリッピングの目的は、フィッティング動作が多重散乱よりも単一散乱により関連する高空間周波数構造に対してあまり敏感にならないことを保証することである。
【００１７】
ここで大域的アプローチに目を向けると、実施形態ではカーネルプロバイダは、画像全体の信号強度及び／又は信号勾配に基づいて、多重散乱イベント用の少なくとも１つのカーネルを提供するように構成される。
【００１８】
実施形態では、カーネルプロバイダが、対象の物理的特性に基づいてカーネルの少なくとも１つを提供するように構成される。これは、多重散乱イベントタイプのカーネル選択に特に便利である。具体的には、多重散乱のより大域的な範囲(reach)をより適切に考慮するために、カーネルプロバイダは、対象の少なくとも一部の形状又はサイズに基づいて、多重散乱イベント用の少なくとも１つのカーネルを提供するように構成される。
【００１９】
実施形態では、カーネルプロバイダは、予め計算されたカーネルタイプのライブラリに基づいて、カーネルの少なくとも１つを提供するように構成される。
【００２０】
実施形態では、システムが、画像に基づいて組織タイプを決定するように構成された組織決定器を有し、カーネルプロバイダは、少なくとも２つのカーネルのうちの少なくとも１つを提供するときに前記組織タイプを考慮する、更なる改良が想定される。これは、特に、撮像された解剖学的構造内の骨の存在によって引き起こされうる不均質な密度分布を有する対象について、一層正確な散乱補正を得ることを可能にする。より具体的には、軟組織の光電効果及びコンプトン散乱に対する質量減衰係数は、水の係数と非常に似ているが、骨のＸ線吸収及び散乱特性は大きく異なる。言い換えると、同じ一次吸収に関して、水対象からの散乱は、骨対象からの散乱よりも大きい。組織決定器は、特に骨及び軟組織の異なる散乱特性を考慮するように、局所材料組成に適応された散乱カーネルを形成することを可能にする。
【００２１】
第２の態様によれば、撮像装置によって取得される対象の画像を受け取るステップと、少なくとも２つの散乱タイプについてそれぞれの散乱カーネルを提供するステップと、提供された少なくとも２つのカーネルに基づいて画像内で補正を実行するステップと、を有する画像処理方法が提供される。
【００２２】
第３の態様によれば、少なくとも１つの処理ユニットによって実行されるとき、処理ユニットにこの方法を実行させるように適応されるコンピュータプログラムが提供される。
【００２３】
第４の態様によれば、プログラムを記憶したコンピュータ可読媒体が提供される。
【００２４】
入力画像は投影ドメインから描出されることができるが、提案されるシステム及び方法の再構成画像への適用は、本明細書において除外されない。画像ドメインにおける散乱補正は、より精度の低い散乱補正を用いて、投影画像からの第１の再構成を実行することによって行われることができる。散乱は上述のように推定されるが、この場合は再構成画像に基づいて推定される。補正された再構成は、散乱カーネルのより良い選択のために使用されることができ、散乱カーネルは次いで、第２の再構成のために使用される、より正確な散乱のない投影画像を生成するために使用される。提案される方法は、肺組織のような空気封入部のポケットを有する海綿状組織が骨のような高密度組織に近接している胸部／肺画像において特に有益であることが分かった。
【００２５】
本開示は主に医用Ｘ線画像に関するが、手荷物スクリーニング、又は非破壊材料試験などの非医療用途は本明細書において除外されない。
【００２６】
本開示における「１Ｄ」、「２Ｄ」、又は「３Ｄ」の記述は、それぞれ、空間的な１次元、２次元、又は３次元の省略形である。
【図面の簡単な説明】
【００２７】
【図１】撮像装置の概略ブロック図。
【図２】単一散乱イベント及び多重散乱イベントの図。
【図３】単一散乱イベントＡ）及び他多重散乱イベントＢ）の点拡がり画像を示す図。
【図４】一実施形態において本明細書で想定される画像処理システムの概略ブロック図。
【図５】一実施形態において本明細書で想定されるカーネルモデルを示す図。
【図６】一実施形態による画像処理方法のフローチャートを示す図。
【発明を実施するための形態】
【００２８】
本発明の例示的な実施形態は、一定の縮尺で描かれていない以下の図面を参照して説明される。
【００２９】
図１を参照して、これは、Ｘ線撮像装置ＩＡ及び画像処理システムＩＰＳを備えるＸ線撮像装置の概略図である。
【００３０】
撮像処理システムは以下で更に詳細に説明するように、撮像装置ＩＡによって提供される画像データＩＭを処理するように構成される。
【００３１】
まず、撮像装置ＩＡを見ると、これは図１において側面立面概略図で示されている。撮像装置ＩＡの実施形態は、対象ＯＢのＸ線画像ＩＭを取得するように構成されたＣアーム撮像装置、ＣＢＣＴ（Cone Beam CT）スキャナ、ＣＴスキャナ、マンモグラフィ装置又は放射線撮影装置等を含む。
【００３２】
対象ＯＢは、撮像軸ｚを有する断面で示されている。図１において、撮像軸ｚは、描画平面に対し延びている。本明細書で想定されるＸ線イメージャの主な用途は医療分野であるが、非破壊材料試験又は手荷物スクリーニングなどの非医療の状況も本明細書において除外されない。従って、「対象ＯＢ」という用語は本明細書では人間又は動物の患者、又はその解剖学的部分などの生気のある「対象」を含むように一般的な意味で使用されるが、無生物の対象も含む。
【００３３】
より詳細には、Ｘ線撮像装置ＩＡは、Ｘ線源ＸＳと、Ｘ線感応検出器ＸＤとを有する。使用時には、対象ＯＢが、Ｘ線源ＸＳ及びＸ線検出器ＸＤの範囲内の検査領域内で軸ｚに沿って位置決めされる。Ｘ線源ＸＳは、エックス線ビームＸＢを生成するように通電され、これは焦点スポットＦＳから放射され、検査領域、従って対象ＯＢの少なくとも関心領域を横切る。図１において、ＯＸは光軸であり、Ｘ線ビームＸＢの主伝搬方向であり、軸ｙはＯＸに平行である。具体的には、光軸ＯＸは、焦点ＦＳから検出器ＸＤまで延び、ある点で検出器Ｄと交差する。画像表面は座標ｘ，ｚを有し、ｘは撮像軸ｚに直交し、ｘ，ｚは両方とも光軸ＯＸ又は軸ｙに直交する。
【００３４】
Ｘ線は、対象ＯＢの物質（例えば、組織、骨等）と相互作用する。相互作用後、Ｘ線放射線は、対象ＯＢの向こう側に現れ、Ｘ線検出器ＸＤに衝突する。衝突するＸ線は、検出器ＸＤによって検出され、電気信号に変換される。電気信号は、適当な変換回路（図示せず）によって画像値に変換され、次いで画像値は、画像処理装置ＩＰＳによってＸ線画像に処理されうる。画像値は、Ｘ線画像を形成する行及び／又は列に配置されることができる。適切なイメージングソフトウェアが、１つ又は複数のディスプレイデバイスＤＤ上に画像を視覚化するために使用されることができる。画像は、記憶されるか、又は他の方法で処理されることができる。Ｘ線画像は、撮像対象ＯＢの内部の詳細を明らかにすることができる。これは、撮像された対象ＯＢの診断及び治療又は他の検査に役立つことができる。
【００３５】
Ｘ線撮像装置によって提供されるＸ線画像は、投影画像又は再構成画像のいずれかである。従来のラジオグラフィは、投影画像ベースである。再構成画像は例えば、ＣＴ又はＣアームイメージングの状況において現れる。この場合、信号処理チェーンは、フィルタード逆投影などの再構成アルゴリズムを実行する再構成器ＲＥＣＯＮを有する。再構成器ＲＥＣＯＮは、検出器の平面内に位置する投影画像を処理し、投影画像を、検出器平面の平面に垂直な平面内の画像ドメイン又は検査領域内に位置する１つ又は複数のスライス画像に変換する。従って、再構成器は、投影ドメインからの投影画像を、画像ドメインにおけるスライス画像にマッピングする。一般に、複数の投影画像が、対象の複数の異なる方向から取得されることが理解されるのであろう。ＣＴでは、これは、それぞれ異なる方向に沿って複数の投影画像の組を取得するためにＸ線源及び／又は検出器が対象ＯＢの周りを回転する回転機構によって達成される。この投影画像の組は、１つ又は複数のスライス画像を取得するために再構成器ＲＥＣＯＮに供給され、それらスライス画像は一緒に、対象ＯＢを含む検査領域の画像近似である画像ボリュームを構築すると考えることができる。
【００３６】
撮像装置ＩＡはまた、撮像されるべき関心領域により良く適合するようにＸ線ビームを成形するために、Ｘ線源に配置されるコリメータ（図示せず）を有することができる。検出器ＸＤは、フラットパネル検出器であってもよいが、本明細書では湾曲した実施形態も想定される。検出器ＸＤは、直接変換タイプでも、間接変換タイプでもよい。
【００３７】
差込図１Ａには、Ｘ線検出器ＸＤの画像表面を、軸ｘ、ｚ及び画像表面の４つの外縁で表す概略的なグリッドレイアウトが平面図で示されている。軸ｚは、上記で述べた撮像軸である。図１Ａのグリッドレイアウトにおける小さな正方形は複数の「行」（慣例により、方向ｘに沿う）及び「列」（慣例により、方向ｚに沿う）に配列された検出器セル又はピクセルを表す。代替的に、単一のライン／行の検出器もまた本明細書において想定される。前述のように、画像表面は、図１Ａに示すように平坦であってもよく、又はＸ線源ＸＳの焦点ＦＳを通り軸ｚに平行な長手軸（以下、「曲率軸」と呼ぶ）を有する仮想円筒の部分側面を形成するように湾曲されてもよい。これらの実施形態では、検出器ＸＤの画像表面は、全体としてＸ線源ＸＳの焦点スポットに対しフォーカスされる。
【００３８】
焦点ＦＳから見たとき、Ｘ線検出器の頂部には、散乱放射線への曝露を低減するために、散乱防止グリッドＡＳＧを配置することができる。ＡＳＧは、放射線不透過材料から形成されるストリップを有する。ストリップは、ストリップのそれぞれの長手方向エッジがＸ線源ＸＳに向かうように、１つの組の中で平行に配置される。ストリップは、平行に配置された面を有することができる。代替的に（図１又は図２には示されていない）、それらの面は、一般に平行でなく、ストリップが焦点に向けてフォーカスされるように角度を付けられる。この場合、ストリップのそれぞれの延長部である仮想平面は、仮想線で交差する。仮想線は焦点スポットを通り、ｚ軸に平行に延びる。
【００３９】
図１の側面立面図において、ストリップは、軸ｚに沿って図面の平面内に延びる。この配置は、一般に１Ｄと呼ばれる。代替として、２組のストリップを有する２Ｄレイアウトがある。ストリップは１つの組の中で（例えば、ｘに沿って）平行であり、軸ｚに平行であるストリップの第２の組を横切って延びる。代替として、１Ｄレイアウトが折りたたまれて２Ｄレイアウトを構成してもよい。
【００４０】
本明細書で想定されるように、提案される撮像装置は、好ましくは、撮像結果に対する散乱放射線の悪影響を低減するための２つの機構を有する。１つの機構は、散乱に対するハードウェア防御を形成する上述の散乱防止グリッドＡＳＧである。ＡＳＧの代わりに、好ましくはそれに加えて、上述の処理システムＩＰＳの形態の第２の散乱低減メカニズムも存在する。散乱放射線（本明細書では単に「散乱」と呼ぶ）が考慮されない場合、より低い画像品質、特により低いコントラストにつながることが知られている。
【００４１】
画像処理システムＩＰＳは、画像ＩＭに適用される補正値を計算することによって、画像に対する散乱効果を低減するように構成される。言い換えると、画像処理システムＩＰＳは、散乱の影響を低減するように計算的に作用する。
【００４２】
提案する画像散乱低減システムＩＰＳの動作をより詳細に説明する前に、ここで、散乱現象をより詳細に見る図２を参照する。
【００４３】
散乱は物質とｘ線との間の相互作用の一方のモードであり、他方は吸収である。位相変化など、本明細書では考慮されない他のモードもある。散乱現象自体は、小角（又はレイリー）散乱とコンプトン散乱とに分類されることができる。以下の議論はコンプトン散乱に焦点を当てることがある。
【００４４】
焦点ＦＳから放射されるＸ線ビームＸＢは検出器Ｄ上の異なる検出器ピクセル位置に焦点ＦＳからキャストされうるそれぞれの仮想の幾何学的光線の組に沿って伝搬する個々のペンシルビームの発散束と考えることができる。Ｘ線光子は、幾何学的光線のそれぞれに沿って進み、それらのペンシルビームを形成する。光子は物質に衝突する前に、Ｘ線源ＸＳの焦点ＦＳから出てくるときに初期エネルギーを有する。所与の幾何学的光線に沿って進むいくつかの光子は、吸収され、従って、それぞれの幾何学的光線の端部に位置するそれぞれの検出器ピクセルで見られるように、Ｘ線の強度を低下させる。しかしながら、低下した強度の全てが吸収イベントに起因するわけではない。所与のペンシルビーム内の他の光子は散乱されることができるので、図２Ａに示されるように、幾何学的光線に沿って意図された経路から逸れて発散する。そのような幾何学的光線の１つ、又はペンシルビームＰＢの意図される経路は、所与のピクセルにおいて検出器面と交差する太矢印として示される。散乱イベントＳＥのため、光子は意図された経路から逸れて発散し、その代わりに、意図されたピクセルにかなり近い別の検出器ピクセルに入射する。従って、意図されたピクセルで見られるべきであったある量の強度は代わりに、散乱光子ｓｐｈ＿ｄとして別のピクセル位置で記録される。しかし、散乱イベントが、散乱光子ｓｐｈ＿ｇによって示されるように、光子がある角度を超えてその意図された経路から発散するようなものである場合、光子は、ＡＳＧのストリップによってフィルタ除去される。図２のＡ）の状況は、任意の所与の光子に対して単一の散乱イベントしか存在しないという命題を反映している。この仮定の下で、単一散乱光子は、所与の光子の幾何学的光線が検出器Ｄと交差するところの近傍ピクセルに局所的に記録されるか、又はＡＳＧのフィルタ作用のために全く記録されない。
【００４５】
しかしながら、図２のＢ）に示すように、光子は、ＳＥ１及びＳＥ２のような多重散乱イベントに遭遇することがよくある。多重散乱イベントの場合、光子ｓｐｈ＿Ｄは必ずしも散乱防止グリッドによってフィルタリングされず、その代わりに、その初期経路の幾何学的光線が検出器Ｄと交差する意図されたピクセル位置から比較的遠く離れた検出器ピクセル位置で検出される。
【００４６】
言い換えると、散乱は、単一散乱イベントと多重散乱イベントの２つのタイプに区別されることができる。それらの各々は検出器に異なる影響を与え、従って、検出器によって測定される全体の放射線に対して異なる寄与をする。多重散乱イベントは測定に大域的に寄与する傾向があり、一方、単一の散乱イベントはより局所的に寄与する。以下では、簡潔にするために、単一散乱イベントを「第１タイプ／カテゴリ」の散乱と呼び、多重散乱イベントを「第２タイプ／カテゴリ」の散乱と呼ぶことがある。
【００４７】
２つのタイプの散乱に対する点拡がり関数ＰＳＦの例を図３のＡ）、Ｂ）に示し、図３のＡ）は単一散乱イベントの点拡がり関数を示し、一方、図３のＢ）は多重散乱イベント点拡がりを示す。図３のＡ）、Ｂ）に示される検出器面内の投影像は、焦点ＦＳから、上方から見た平面図である。隣接するそれぞれの図は、それぞれｘ及びｚ方向の強度プロファイルを示す。単一散乱イベントはより集中した点拡がり関数につながり、一方、多重散乱イベントは、より広いリッジを有するより広がった点拡がり関数につながる。
【００４８】
本明細書で提案される画像処理システムＩＰＳは、多重散乱及び単一散乱の２つのタイプの散乱を異なって考慮するように構成される。画像処理システムＩＰＳは、２つのタイプの散乱を考慮に入れるために、補正されるべき画像ＩＭに対して計算上異なるように作用するように構成される。より具体的には、ここでは多重散乱イベントと単一散乱イベントとを区別することができる散乱モデルを洗練させることが提案される。
【００４９】
提案する処理システムＩＰＳの概略ブロック図を図４に示す。画像処理システムは、検出器から受け取った投影画像を処理するように構成されることができる。代替的に、画像処理システムは、再構成された画像に作用することもできる。入力画像ＩＭ（すなわち、投影画像又は再構成画像）は、画像処理システムＩＰＳの入力ポートＩＮで受け取られる。実施形態において、提案される画像処理システムはカーネルベースであり、言い換えると、補正値を計算するために散乱点拡がり関数ＰＳＦの近似を使用する。実施形態では、カーネルは予め計算され、データベースＦＤＢに保持される。
【００５０】
本明細書で想定されるように、カーネルプロバイダＫＰは、入力画像ＩＭに記録された信号にそれぞれが適合する単一散乱及び多重散乱用の専用カーネルをそれぞれ提供する。２つの散乱タイプに対するカーネルが、それぞれ、入力画像に記録された単一散乱寄与及び多重寄与を推定するために使用される。次に、そのように推定された多重散乱量及び単一散乱量が、補正器モジュールＳＣＭによって入力画像ＩＭに適用される。好ましくは、多重及び単一散乱の推定量が入力画像ＩＭから点ごとに減算される。その後、補正された画像ＩＭ＿Ｃが出力ポートＯＵＴに出力される。ＩＭ＿Ｃは、そののち可視化され、記憶され、又は他のやり方で（更に）処理されることができる。
【００５１】
データベースＦＤＢに保持されるカーネルは、単一散乱イベント用と多重散乱イベント用の２タイプである。更に、単一散乱又は多重散乱の各カーネルタイプ又はカテゴリに関して、更なる階層化が、特定のサンプルボディのサイズ及び形状に従って存在し、かかるサイズ又は形状に基づいてそれぞれのカーネルが予め計算される。第１及び第２タイプの散乱の上記の記述と一致して、本明細書では、カーネルが単一散乱又は多重散乱のどちらの特性を捕らえるよう構成されているかに応じて、「第１」又は「第２」タイプのカーネルを参照することができる。
【００５２】
一実施形態では、サンプル本体又はファントムが、図５の側面図及び概略図で示す水球（又は他の材料）のモデルと考えることができる。強度を記録することとは別に、入力画像における画像値は、画像値間の局所的な差の近似である勾配も記録する。勾配は、対象ＯＢにおける組織遷移又は組織ー空気遷移、又は撮像された患者の厚さの変化、特にエッジ又は境界部分に向かう変化を表す。
【００５３】
異なる強度及び勾配をモデル化するために、サンプリングファントムＷＳは、図５に示すように３つのパラメータによってパラメータ化される。１つのパラメータＲは、ファントム本体ＷＳのサイズをパラメータ化する。水球のようなボールの場合、Ｒはその半径である。別のパラメータφは、散乱防止グリッドＡＳＧのストリップに対する仮想ボディＷＳの相対的な配置をパラメータ化する。これは、ボディＷＳを通過する幾何学的光線ＧＲの周りの回転角φで記述できる。第３のパラメータは、ファントム本体ＷＳを通る光線ＧＲに沿った経路長Ｔである。経路長は、所与の位置に対する所与の光線ＧＲからのオフセットＯによって決定される。経路長Ｔは、異なる勾配をモデル化することを可能にする。パラメータの組は２つのコンポーネント、すなわち、ファントムの大きさについての１つ又は複数のコンポーネント（ボール形状ファントムの実施形態におけるパラメータＲなど）と、それぞれの画像位置に対するファントムの相対位置に関連する１つ又は複数のコンポーネント（Ｔ（又はＯ）及びφなど）とに区別されることができる。前者のコンポーネントは、「サイズコンポーネント」と呼ばれることがある。再び、この概念は、任意の形状のファントムに一般化することができ、ボールは一実施形態である。
【００５４】
パラメータ（Ｔ，φ，Ｒ）［又は（Ｏ，φ，Ｒ）］の各組について、点拡がり関数ＰＳＦとしても知られているカーネルは、水のようなファントムボディの物質特性（濃度、Ｚ数）が知られているので、実験的に又は好適には理論的シミュレーションによって、得られることができる。カーネル又は点拡がり関数ＰＳＦは、パラメータ化されたファントム本体ＷＳの「散乱シグネチャ」と考えることができる。以下でより詳細に説明するように、カーネルはファントムが曝露される概念的なペンシルビームによって引き起こされる散乱として現れる。次に、カーネルは、ファントムＷＳ対ペンシルビーム相互作用によって引き起こされる散乱を記述し、この散乱が画像構造に関してどのように現れるかを記述する。点拡がり関数ＰＳＦは２Ｄスカラ場であり、それ自体「散乱画像」として表されることができる。異なるパラメータを有する一連の異なるファントム身体ＷＳのカーネルは、好ましくは予め計算される。カーネルのライブラリ全体を得るために、パラメータ（Ｔ，φ，Ｒ）は、適切なステップでインクリメントされる。
【００５５】
ここで、シミュレーションに基づいてカーネルを得る１つの方法をより詳細に説明する。国際公開第２００７／１４８２６３Ａ１号パンフレットなどの以前のアプローチとは異なり、本明細書で提案されるカーネルライブラリは、単一散乱イベント及び多重散乱イベントそれぞれについて異なるカーネルを計算することによって改良される。これらの異なるカーネルタイプの計算は、パラメータ化されたカーネルＷＳの各々に基づいてモンテカルロ（ＭＣ）シミュレーションを実行することによって行うことができる。ＭＣシミュレーションでは、線源ＦＳから始まるペンシルビームとして個々の光子がシミュレートされる。次に、これらは、ｉ）ファントムＷＳによって吸収されるか、又はｉｉ）検出器（又は散乱防止グリッド）に当たるか、又はｉｉｉ）シミュレーションボリュームを出るまで、シミュレーション全体にわたって追跡される。光子と物質との相互作用は確率論的であり、確率と確率分布は理論から知られており、これらは光子エネルギー、ファントム物質の密度、Ｚ数に関する物質タイプの関数である。好ましくは、水及びその特性が使用される。次いで、光子の経路が現在の位置から次の相互作用までランダムに選択されるが、このランダム選択は正しい確率分布に従って行われる。これは、上で述べたイベントｉ）〜ｉｉｉ）のいずれかが起こるまで、各光子について行われる。このようにして、検出器に当たる各光子について、検出器が以前に何回相互作用したかを正確に知ることができる。この追跡とシミュレーションは、多数のシミュレートされる光子について行われる。
【００５６】
シミュレーションは、最初にN Metropolis et al in "Equation of State Calculations by Fast Computing Machines", Journal of Chemical Physics. 21 (6), pp 1087-1092で報告されたメトロポリス−ヘイスティング(MetropolisーHastings)アルゴリズム（又はその変形）を使用して効率的に実施されることができる。
【００５７】
シミュレートされる光子の個数は、１０５又は１０６以上のオーダーでありうる。次いで、検出器に当たる全てのシミュレートされる光子は、２つのグループ、すなわち、検出器に当たる前に１回相互作用したグループと、複数回相互作用したグループとに分けられることができる。これら２つのグループを使用して、２つのカーネルカテゴリが作成される。具体的には、１回相互作用のグループで検出された光子の信号はそれらが一緒に、ファントムＷＳの特定のパラメータのための、単一散乱イベント用のＰＳＦを形成する。複数回相互作用グループで検出された光子は一緒に、ファントムＷＳの特定のパラメータのために、多重散乱イベント用のＰＳＦを形成する。次に、パラメータインクリメント（例えば、水球モデルＷＳにおける異なるＲ、異なるφ、異なるＴ）を変更することによって、異なるカーネルパラメータ化について上記のシミュレーションを繰り返し、それによって、２層カーネルライブラリを構築することができる。
【００５８】
以上から、各ファントムＷＳ（そのパラメータｉを介して識別可能）は、１つは多重散乱用であり及び１つは単一散乱用である２つのそれぞれのカーネル又は「散乱画像」と関連付けることができる。更に、各ファントムは、それぞれの「理想画像」と呼ぶことができる別のデータコンポーネントとも関連付けられることができることが理解されよう。その結果、各カーネルは、所与のファントムＷＳのそれぞれの理想画像と関連付け可能でなる。この理想画像は、それぞれのファントムＷＳをフルビーム（単にペンシルビームではない）で（理論的に又は実験的に）曝射するとき、及び散乱効果が無視されるとき、導き出されることができる。このような理想画像は、ＭＣシミュレーションにより得ることができる。代替として、本出願人の米国特許出願公開第２０１０／０１８９３７６号（式（２）及び（３））に示されているように、理想的な画像の画像特徴（信号レベル及び／又は信号勾配）を、仮想ＷＳの空間パラメータ（例えば、例示的な水球の実施形態の場合、（Ｔ，φ，Ｒ）又は（Ｏ，φ，Ｒ））に直接的に関連付けることができる。この後者の実施形態では、シミュレーションは不要である。２つのカテゴリのカーネルを予め計算し、記憶することに加えて、各ファントムのそれぞれの理想画像を同様に予め計算し、記憶することが好ましい。
【００５９】
このようにして、各パラメータごとに、単一散乱イベント用及び多重散乱イベント用の、別個のカーネルＰＳＦを得る。従って、カーネルライブラリＦＤＢは、第１及び第２タイプのカーネルを有し、それぞれのタイプが異なったパラメータ（Ｔ，φ，Ｒ）によって階層化される、２層形式で編成されることができる。
【００６０】
なお、上記の水球に基づくことは必須ではなく、他の形状及び／又は他の材料も考えられる。実際に、実施形態においてボール又は球を使用する代わりに、以下でより詳細に説明されるように、カーネルＰＳＦを生成するためのファントムの形状は、撮像されるボディのボディタイプ又は領域に従って選択される。先に述べたように、一実施形態では２つのタイプのカーネル、単一散乱用のカーネル及び多重散乱用のカーネルが予めに計算され、データベースに記憶される。しかしながら、代替的に、カーネルは、オンザフライで、オンデマンドで計算されることができる。これは、例えば、ＡＳＩＣ等のような適切に同調された専用プロセッサを使用することによって達成されることができる。しかしながら、好ましくは、カーネルは、散乱補正画像ＩＭ＿Ｃが展開中により迅速に生成されることができるように、ターンアラウンド時間を減少させるように予め計算される。
【００６１】
実施形態では、カーネルプロバイダＫＰは、単一散乱イベント及び多重散乱イベントにについてそれぞれ別々にカーネルを選択し最適化するように構成される。これは、両方の散乱タイプを考慮するために１つのカーネルを選択するのではなく、散乱イベントタイプについて個別にかつ別々に最適化することを可能にする。先に述べたように、単一散乱イベントは、考慮される画像位置の周りでより局所的な効果を有し、一方、多重散乱イベントは、より大域的な範囲を有する。大域的範囲と局所的範囲との間のこの区別は、好ましくは提案するシステムＩＰＳにおいて、単一散乱及び多重散乱の２つのタイプに関するカーネルの別々の選択によって反映される。従って、カーネルプロバイダは２段階で動作する。
【００６２】
いずれかの段階におけるカーネル選択は一般に、ファントムを介して、より具体的にはファントムの関連する理想画像を介して行われる。理想的な画像は、入力画像内のパッチに対してマッチングされる。ベストマッチが見つかると、ベストマッチする理想的な画像に関連付けられたカーネルが、所定のパッチに対して選択され、プロシージャが、入力画像全体にわたって一部又は全てのパッチに対して繰り返される。次いで、選択されたカーネルは、所与のパッチに関連付けられるペンシルビームの単一散乱イベント又は多重散乱イベント推定を形成する。散乱画像は一般に、考慮されるパッチよりも大きい。
【００６３】
まず、カーネルプロバイダＫＰが単一散乱イベントカーネルを選択するように動作する第１の段階に移ると、入力画像内のいくつかの画像位置又は各画像位置に対して、適切な画像パッチ、例えば１０×１０ピクセルの二乗が規定される。他のパッチサイズ又は形状も想定される。パッチのサイズ及び形状は、一般にデータベースＦＤＢに保持されているカーネルのサイズ及び形状の関数である。
【００６４】
次に、それぞれのパッチにおける記録された画像ＩＭ内の画像情報を用いて、データベースＦＤＢに記憶されている単一イベントカテゴリから、ベストマッチするカーネルを見つける。最適化アルゴリズムを実行して、関連する画像パッチにベストマッチするカーネルを見つけることができる。一実施形態では、これは、第１のカテゴリ内の全てのカーネルにわたって繰り返し、カーネル及び画像パッチを計算するために使用されるファントムＷＳの理想的な画像の各々の間の各ピクセルごとの差を取ることによって形成される個々の差画像を比較することによって、簡単に行うことができる。上記のように、理想的な画像は、検出器上に散乱放射線がないとした場合に、それぞれのファントムＷＳが検出器（又は検出器パッチ）上に生じさせるであろう画像である。次に、そのような差画像を最小にするカーネルが、ライブラリＦＤＢ内の予め計算された第１カテゴリのカーネルの有限集合から選択され、ここで、「最小性」は、適切な尺度ｄによって測定される。適切な尺度ｄは、場合によっては重み付けされた（パッチ又は全領域にわたって合計された）二乗差又はパッチ等の最大差分値を含む。尺度ｄは、画像に記録される強度に単に基づいて規定されることができる。代替として又は付加的に、パッチを使用して、例えばフィルタリングによってペンシルビームが当たるそれぞれのピクセルｊ内の局所勾配を推定し、これをカーネルのその位置ｊにおけるピクセル強度及び／又は勾配にマッチングさせることができる。
【００６５】
しかしながら、大きいカーネルライブラリの場合、予め計算されたカーネル又はファントムの有限集合における探索の上述のアプローチは、遅すぎる可能性があり、共役勾配、勾配降下、ネルダ−ミード(NelderーMead)法などの最適化スキームを使用して、コスト関数が最小化されることができる。コスト関数Ｆは、上述したように差分画像及びメトリックｄに関して定式化されることができ、これは以下のように定式化されることができる：
2021525141.tif 000003
カーネルの最適化に関して、
2021525141.tif 000004
ここで、
2021525141.tif 000005
ここで、ｊは、考慮される画像位置であり、ｉは、第１タイプのＰＳＦ１のカーネルのインデックス／パラメータであり、Ｐ（ＩＭ）ｊは、ｊを中心とする入力画像のパッチであり、Ｋ＿ＩＭｉ１は、検出器（又はそのパッチ）上のファントムの理想的な画像であり、Ｋ１＿ｊは、実画像パッチ（メトリックｄに従う）にベストマッチする理想的な画像をもつファントムに関連付けられた単一散乱カテゴリのカーネルであり、ｄ（．，．）はスカラへのマッピングによって差分を定量化する差分尺度である。
【００６６】
水球の実施形態では、ｉ＝（Ｔ，φ，Ｒ）である。最適化（２）が、カーネルが利用可能でないことを示すパラメータ集合ｉ＊＝（Ｔ＊，φ＊，Ｒ＊）を返す場合、これは、ほとんどの場合に当てはまる可能性があり、予め記憶されたカーネルから所望の最適カーネルｉ＊を補間するために、線形補間がｉにわたって実行されることができる。この場合、予め計算されたカーネルを線形ルックアップテーブルＬＵＴに記憶することが有利である。代替として、最も近い利用可能なカーネルを選択し、例えばペンシルビームの全吸収にマッチするように適切にリスケーリングされてもよい。上記及び以下の全ての実施形態では、位置ｊは、物理的ピクセル又は非物理的ピクセル（ダウンサンプリング後など）を指すことができる。
【００６７】
楕円体又は他のものような、上述のボールモデル以外の形状プリミティブが使用される場合、フィットするべき３より多くの自由パラメータが存在し得る。この場合、画像からのより多くの画像特徴を使用されることができる。実施形態では、所与のピクセル位置ｊの４つの直接近傍の等価水厚さが、５つのパラメータについてフィットさせるために使用されることができる。近傍からの更に多くのピクセルを使用して、カーネルの「形状プリミティブ」を局所的な画像パッチにフィットさせることができる。一般に、ある画像位置の散乱カーネルは、ピクセル値又はそれらのピクセル値の関数（例えば勾配）を使用して、近傍のピクセル群に基づいて選択されることができる。
【００６８】
次いで、単一散乱タイプのカーネルに対する上述のフィッティングプロシージャが、入力画像ＩＭ内の各位置ｊ及び全ての位置ｊについて実行されることができる。しかしながら、代替の実施形態において、これは入力画像ＩＭを粗くするようサブ領域に分け、ピクセルｊあたりではなくサブ領域あたりのみで上述の最適化スキームを実行することによって、単純化されることができる。
【００６９】
上述のフィッティングプロシージャの他の変形が企図され、その場合、サイズコンポーネントなどのパラメータｉのサブセットに基づいて、ベストマッチが見つけられ、次いで、このベストマッチが、各位置ｊについて（１）、（２）に従って残りのパラメータを最適化することによって局所化される。
【００７０】
単一散乱イベントカテゴリに対する全てのカーネルが計算されると、これらを使用して、入力画像ＩＭのピクセルごとの補正を実行し、単一散乱の損傷を補正することが可能である。これは、例えば、以下のように行うことができる：
2021525141.tif 000006
ここで、ｊは、それぞれの画像位置、又はサブ領域などについて選択された全ての最適化された第１タイプのカーネルにわたって実行される。選択されたカーネルは加算され、画像ＩＭから減算される。この段階で得られた画像ＩＭ＿Ｃｐは、これまでに単一散乱の損傷のみが考慮されてきたので、部分的に補正された画像と呼ばれることができる。
【００７１】
（３）、（４）のプロシージャは、所定の数のステップの間、又は収束に達するまで、反復的に繰り返されることができる。
【００７２】
多重散乱は、カーネルプロバイダＫＰの次の動作段階で処理される。第１段階とは別に、第２段階では、カーネルプロバイダＫＰは、一般に、画像内の多重散乱シグニチャとマッチングするための複数の異なるカーネルを選択する。第２段階の処理は、多重散乱の大域的特性にあわせて調整される。具体的には、多重散乱イベントによって引き起こされる画像勾配が、画像ＩＭにおいてあまり目立たなくなることが予想され、高い空間周波数の細部が低減されたより大域的な構造に関連する。この洞察に基づいて、実施形態では、初期画像ＩＭは最初に、画像内のダイナミクスを平坦化することによって前処理を受ける。これを実行するための、実施形態において想定される１つの方法は、画像内の勾配を計算し、次いで、所定の閾値で勾配値をクリッピングすることである。この平坦化された画像では、閾値よりも大きい寄与が除去される。
【００７３】
平坦化された画像を得るための別の方法は、初期画像又は部分的に補正された画像に対して、（局所的又は大域的な）ローパスフィルタＦＩＬを実行することである。ローパスフィルタＦＩＬによるこのフィルタリングは、画像処理チェーン内の更に上流において他の目的のために既に行われている可能性がある任意の他のフィルタリング（ノイズフィルタリングなど）に加えて、又はその代わりに行われることができる。特に、平坦化のためのフィルタモジュールＦＩＬは、より積極的に設定されており、上述の第１段階を実行する前に画像ＩＭに適用されることができる任意の可能なフィルタリングよりも、大きな画像部分をカバーする。例えば、好ましくは、入力画像は、ノイズを除去するためにフィルタリングされる。このノイズフィルタリングは、上述のローパスフィルタＦＩＬに使用されるフィルタカーネルサイズよりも小さいサイズの畳み込みカーネルを用いて局所的に行われる。
【００７４】
この平坦化された画像に基づいて、前述したように画像パッチを規定し及び（１）、（２）に関して上述した最適化プロシージャに基づいて尺度ｄに関してベストマッチを見つけることにより、ライブラリＦＤＢ内の第２のカテゴリカーネルからカーネルが選択されるが、今回は、第１タイプのＰＳＦｉ１の代わりに第２カテゴリの点拡がり関数ＰＳＦｉ２が用いられる。マッチングが勾配ベースである場合、最初にローパスフィルタリングが実行され、そののち勾配が計算される。いずれの実施形態においても、補正モジュールＳＣＭは、次に、画像ＩＭ＿Ｃｐを補正するために、選択された全ての第２タイプのカーネルＫ２＿ｊを使用して、例えば以下のように、補正画像ＩＭ＿Ｃｐを形成する：
2021525141.tif 000007
【００７５】
第２タイプのカーネルを選択するための別の実施形態では、カーネルマッチングが個々の画像パッチに対しては実行されず、画像値又は勾配に基づいて、画像ＩＭ全体に対して大域的に実行される。この場合、本明細書では「大域的ファントム」と呼ぶことができる（単一の）最良マッチングファントムを、そのパラメータが撮像されるべき患者の部分に最良に対応するように、画像全体に対して選択することができる。一実施形態では（全ての実施形態ではないが）、大域的ファントムのフィッティングは、パラメータのサブセットのみを調整することによって行われてもよい。このための１つの例示的な実施形態では、大域的ファントムを見つけるために、サイズ成分のみを調整することができる。より具体的には、水球ファントム（又は他のボール形状ファントム）に基づく実施形態では、水球モデルの半径Ｒ及び／又はオフセットＯの中心位置のみが、頭部画像などの撮像される身体部分とほぼ釣り合うように選択されることができる。
【００７６】
ベストマッチする大域的ファントムが選択されると、カーネルは、このベストマッチするファントムに基づいて各ピクセルごとに選択され、これらは、式（４）のように、画像ＩＭ＿Ｃｐを補正するために使用される。より具体的には大域的ファントムの大きさパラメータがピクセル位置ごとに保持され、一方、ボールの実施形態のφ及び／又はＴなどの残りのパラメータは局所的なピクセルごとの適合を達成するように調整される。例えば、上の図５におけるボールの実施形態では各ピクセルについての残りのパラメータ（Ｔ，φ）が半径Ｒから、及び任意選択的にはボールの中心に対するピクセルの相対位置Ｏから追従する。
【００７７】
代わりに、ボール以外の形状のファントムを使用してもよい。より具体的には、楕円形又は円筒形のファントムが胸部、腹部、又は脚／腕の画像にそれぞれ適合するのにより適しているので、本明細書では特に想定される。この実施形態では、データベースＦＤＢに保持されたカーネルが２つのカテゴリだけでなく、多重散乱イベントのための第２のカテゴリで層化され、楕円形、円筒形、球形、又は他のいずれか１つ又は複数などの異なる形状タイプのための異なるパラメータに関して、更なる階層化がある。この実施形態では、画像処理システムＩＰＳは、ユーザが所定の画像タスクに使用するカーネル形状を指定できるように構成されたユーザインタフェースを有することができる。代替として、カーネル形状とサイズの選択は、ユーザが指定した画像プロトコルに基づいて自動的に行われる。代替として、関心領域の直径、長軸／短軸、周囲長、長さなどの空間測定値が好ましくは形状仕様を有するユーザインタフェースを介してユーザによって提供され、次いで、カーネルプロバイダが最良のカーネルサイズを取り出す。再び、他の実施形態において他のファントム形状又はファントム材料が選択される場合、上記のいずれかが適合されることができる。
【００７８】
ＩＰＳの上述の動作は、本明細書で全て想定される、多くの異なる変形を可能にする。例えば、第２タイプのカーネルＫ２は、第１の段階で最初に選択されることができ、第１タイプのカーネルが、第２の段階で後に選択されることが理解されよう。付加的に又は代替として、部分画像ＩＭ＿Ｃｐは必ずしも形成されないかもしれない。なぜなら、２つの段階におけるカーネル選択は両方とも、最初に完了され、第１及び第２タイプの選択されたカーネルが、入力画像ＩＭ＿Ｃを補正するために、一度に使用されるからである。
【００７９】
しかしながら、性能上の理由から、第１又は第２タイプのカーネルを最初に適用することによって部分的に補正された画像ＩＭ＿Ｃｐを形成し、次に他のタイプのカーネルを後で適用するインターリーブアプローチが好ましいことがある。
【００８０】
部分的に補正された画像ＩＭ＿Ｃｐが形成される実施形態では上述のように、初期入力画像ＩＭ上ではなく、上記の部分的に補正された画像ＩＭ＿Ｃｐ上の他のタイプのカーネルの選択に基づくことが、上記の代替として想定されることができる。これは、好ましくはより良好な補正結果を達成するためでありうる。より具体的には、、部分的に補正された画像ＩＭ＿Ｃｐが、第１（又は第２）カテゴリのカーネル選択に基づいて取得されると、第２（又は第１）カテゴリからのカーネルのマッチングが部分的に補正された画像ＩＭ＿Ｃｐに基づいて行われる。好ましくは、このプロシージャが以下により詳細に説明されるように、反復的に行われ得る。いくつかの実施形態では、カーネルプロバイダＫＰによるカーネル選択が撮像される関心領域に存在する組織タイプを考慮に入れることによって更に洗練される、更なる洗練が想定される。より具体的には、一実施形態では骨組織は、骨割合画像を計算することによって考慮される。骨割合画像は、二重層検出器又は他のスペクトル画像装置からのスキャンの物質分離によって、又は単一エネルギスキャン画像の中間第１パス再構成の前方投影から、計算されることができる。
【００８１】
この実施形態では、画像処理システムは、骨画像などの組織画像を計算するための組織決定器ＴＩＤを有する。組織タイプの洗練を伴うこの実施形態では、カーネル選択が、骨割合画像などの組織画像から導出される係数を用いた重み付けスキームを使用することにより洗練される。より詳細には、骨割合が高いほど、散乱カーネルの重み係数が小さくなる。実施形態では、画像位置当たりの骨割合情報が焦点と画像位置との間の経路における骨の相対量を記述する重み付け係数（「骨割合画像」）にマッピングされる。
【００８２】
次に、図６のフローチャートを参照して、画像処理システムＩＰＳによって実現される画像処理の方法について説明する。しかしながら、以下に説明されるステップは必ずしも上述のような画像処理システムＩＰＳのアーキテクチャに結び付けられるわけではなく、以下のステップはそれ自体が教示するものとして理解されるのであろう。
【００８３】
ステップＳ６１０において、散乱補正されるべき入力Ｘ線画像ＩＭが受信される。入力画像ＩＭはＨＩＳ（病院情報システム）におけるＰＡＣＳ（画像保管及び通信システム）のような画像データベースから取り出すことができ、又はＸ線撮像装置から直接供給することができる。
【００８４】
散乱処理されるべき入力画像ＩＭは、無線で、又は有線接続を介して受信されてもよい。入力Ｘ線画像ＩＭは、投影領域からの投影画像であるか、又は投影画像から再構成アルゴリズムによって再構成された撮像領域内の再構成画像であってもよい。入力画像が再構成されたものである場合、散乱補正は、好ましくは依然として投影画像に適用される。入力画像ＩＭは、２Ｄでも３Ｄでもよい。画像ＩＭは、頭蓋骨、胸部、又は心臓、又は他のものような、人間又は動物の対象の関心領域のような対象の画像ＩＭを取得することができる。場合によっては、軟組織のコントラストを増加させるために、適切な造影剤が予めに注入される。
【００８５】
提案した方法は単一及び多重散乱イベントをそれぞれ説明するために、精密化した信号モデルに基づいている。
【００８６】
これに関して、ステップＳ６２０において、散乱補正のためのそれぞれのカーネルが、多重散乱及び単一散乱をそれぞれ考慮するために提供される。言い換えれば、少なくとも２つの異なるカーネルタイプが提供される。カーネルの提供は、予め計算されたカーネルのライブラリから選択することを含み、又はオンザフライで要求に応じてカーネルを計算することを含み得る。所定のカーネルのライブラリからの選択は、既存の予め計算されたカーネルからの補間を含むことができる。カーネルをそれぞれ異なるプロシージャで別個に提供することは、多重散乱及び単一散乱のそれぞれについて最適化することを可能にする。
【００８７】
多重散乱用のカーネルと、単一散乱用のカーネルである２つの異なるタイプのカーネルを提供するステップＳ６２０は、単一散乱カーネルを選択するための局所的アプローチと、多重散乱イベントカーネルを選択するための大域的アプローチと、を使用することが好ましい。多重散乱カーネルを選択するための局所的アプローチも考えられる。
【００８８】
具体的には、実施形態では第１又は第２タイプのカーネルを提供するステップＳ６２０が、パッチ内の画像位置の局所的な画像特性に対して（尺度ｄに基づく）ベストマッチする既存の又は補間されたカーネルを見つけることを含む。前記特性は、パッチ及び／又は勾配に記録された強度値を有する。実施形態において、多重散乱イベント用カーネルの選択は、個々の画像パッチではなく、画像全体に対するマッチングを含む。
【００８９】
画像パッチ内の局所画像特性とのマッチングが採用される場合、これは、単一の画像位置ごとに行われてもよいが、より粗い粒度のアプローチで行われてもよく、全ての画像位置ではなくサブ選択のみが考慮される。例えば、ランダム又は他の選択（空間方向ｘ，ｚの１つ又は両方において１ピクセルおき又はｋピクセルごと）が行われることができる。ランダムモードでは、適切なランダムジェネレータを実行することにより、解析のための特定の数の画像位置をランダムに選択して、単一散乱又は多重散乱状況のうちの少なくとも一方についてカーネルをマッチングするようにする。これの他の変形例において、画像全体が、サブ部分に分割され、マッチングはそれらのサブ部分においてのみ実行される。例えば、画像は、４、８又はそれ以上の象限のようなｋ個の象限に分割されることができる。次に、マッチングは象限ごとに行われる。カーネル選択プロシージャを高速化するために、他のサブ分割も考えられる。それに加えて、又はその代わりに、入力画像が最初にダウンサンプリングされることができ、次いで、ステップＳ６２０が、ダウンサンプリングされた入力画像に対して実行される。
【００９０】
本明細書で想定されるデータベースは、単一散乱イベントカテゴリ用のカーネルの組と、多重散乱イベントカテゴリ用のカーネルの異なる組とを保持する。各カテゴリにおいて、カーネルは、関連するファントムＷＳのサイズ、形状、又は他のものによって更に区別可能でありうる。具体的には、一実施形態において、カーネルは、出願人による国際公開第２００７／１４８２６３号パンフレットに記述されているような水球モデルに基づく。しかしながら、水以外のモデル及び／又はボール以外の形状に基づく他のプロトタイプファントムもまた、代替実施形態において想定される。特に、円筒形又は楕円形のファントム（水特性又は他の材料に基づく）は、以下で更に詳述されるように多重散乱イベントタイプのカーネルについて特に想定される。
【００９１】
ここでは単一散乱イベントタイプのカーネル選択について、より詳細に説明する。これは、画像パッチを特定の画像位置に対して定義することによって、上記の局所的アプローチのように行うことができる。次に、画像内の強度値及び／又は勾配値に基づいてマッチングが行われる。画像パッチは、第１のカテゴリの各カーネルに関連付けられたそれぞれの理想的な画像に対しマッチングされる。任意の最適化スキームを使用して、ベストマッチを見つけることができる。コスト関数は、適切な尺度ｄに基づいて、局所的パッチと点拡がり画像ＰＳＦとがどのように異なるかを測定するように設定されることができる。ベストマッチするカーネルパラメータｉ、又は少なくとも閾値よりも良いマッチが、最適化処理において決定され、関連するカーネルが選択され、単一散乱イベント推定用のカーネルとして提供される。最適化が実行される基礎となるカーネルパラメータは例えば、上述のような（Ｔ，φ，Ｒ）、又はカーネルの形状に応じた他のものを含む。
【００９２】
ステップＳ６２０において、第２のカテゴリ（すなわち、多重散乱イベントタイプカーネル）からのカーネルに対してマッチングする場合、以下の実施形態が想定される。
【００９３】
一実施形態は、第１カテゴリの選択について上述したように依然として局所的であるが、この実施形態では、画像は最初にローパスフィルタリングを受ける。このローパスフィルタリングを用いて、局所散乱に対してより典型的な任意の構造が除去されるか、又は少なくとも抑制されて、上記で説明した平坦化された画像を形成する。多重散乱イベントタイプカーネルを選択する目的のために、このローパスフィルタリングは、捕捉からの単一処理チェーンの中の更に上流で行われる任意のフィルタリングに加えて行われることができる。特に、このフィルタリングは、第１カテゴリのカーネルを選択する前に行われる他のフィルタリングに加えて行われてもよい。実施形態では、フィルタリングは局所的であり、すなわち、フィルタリングは、それぞれのパッチに対するマッチング処理を実行する前に、個別に調整されたフィルタ係数及びフィルタカーネルを用いて、各パッチごとに個別に適用される。代替として、画像全体が、同じフィルタ係数を用いて全体的にフィルタリングされる。フィルタリングは、ユーザにより調整可能でありうる。正確なフィルタパラメータは、散乱アーチファクトを調べ、アーチファクトを低減するためにフィルタパラメータを調整することにより、実験的に得ることができる。これを効率的に行うために、自動化された試験プロトコルが使用されることができる。
【００９４】
フィルタリングに加えて又はフィルタリングの代わりに、勾配画像を形成し、固定の閾値又はユーザにより調整可能な閾値により勾配値を閾値処理し、高い勾配を除去して高い勾配をクリップすることによって、平坦化された画像が形成されることができる。フィルタリング又はクリッピングのいずれかの結果は、構造的に平坦化された画像であり、平坦化された画像は、上述したように選択プロシージャの中でパッチごとにカーネルに対しマッチングされる。選択は、式（１）、（２）のようにパラメータ空間のパラメータにわたる最適化であり、又はベストマッチするカーネルを見つけるために比較のための適切な尺度ｄ（）を使用して、データベース内のカーネルを介したより単純な比較ループに基づく。
【００９５】
代替の実施形態では、より大域的なアプローチが使用され、それに従って、画像全体が、第２カテゴリタイプのカーネルに対しマッチングされる。この実施形態では、ファントム形状が、撮像される関心領域に合致するように選択されることが好ましい。例えば、頭部画像の場合、人間の頭部の大きさに合わせてボール形状のファントムが選択される。胸部画像の場合、患者の胸部測定値に相応する軸をもつ楕円形状のファントムが使用され、他のＲＯＩの他の幾何学的形状についても同様である。言い換えれば、この大域的なアプローチでは、多重散乱を推定するために、（画像全体について）ファントムパラメータｉの単一の組のみが決定される。この単一の最良のファントムは、サイズコンポーネントのみに基づくなど、全てのパラメータｉのサブセットのみに関連してマッチングすることにより、見つけられることができる。この単一最良のファントムは、それぞれの画像位置ｊに対し局所的に調整されることもできる。実施形態では、この局所的な調整は、式（１）、（２）のように、しかし、第２カテゴリのカーネルに関して、他のパラメータｉ（すなわち、サイズ成分の組補数）についてのみ最適化を行って、実施されることができる。
【００９６】
上記のカーネル提供／選択スキームのいずれかと組み合わせて、ステップＳ６２０は更に、関心領域における組織タイプを決定して、第１カテゴリタイプのカーネル又は第２カテゴリタイプのカーネルのみならず、組織タイプにも基づいて、選択を更に最適化することを含み得る。このために、入力画像ＩＭによる組織タイプの分布が推定される。関心のある、散乱に関連する１つの特定の組織タイプは、骨などの比較的高密度の組織である。骨分布を推定するために、図１〜図５に関連して上述したように、ステップＳ６２０において骨割合画像が決定されることができる。骨の割合に応じて、上述した最適化において見つけられたそれぞれの第１カテゴリのカーネル及び／又は第２カテゴリのカーネルに係数が適用されることができる。この係数は、骨及び水の異なる散乱特性を考慮するために、水ベースのカーネルに適用される。重み係数は、例えばＭＣシミュレーションを使用することによってオフラインで計算されることができ、又は骨及び水の散乱断面積から導出されることができる。代替の実施形態では、骨密度割合は、予め計算されたカーネルセットの追加パラメータであり、すなわち、カーネルセットは、全ての可能な骨割合（例えば、１０％、２０％など）のカーネルを含む。次いで、骨密度割合は、カーネル選択を実行する際に、（空間パラメータに加えて）追加パラメータとして使用される。
【００９７】
単一散乱及び多重散乱イベントのための少なくとも２つのカーネルが見つかると、各カーネルは、組織分布に従う係数で乗算され、これらが、単一散乱及び多重散乱イベント寄与のための推定値をそれぞれ構成する。次いで、カーネルは、補正処理の中で初期画像から減算されてもよく、又はカーネルに関する補正は、別の方法で行われてもよい。このステップは、オリジナル画像の各画像値に符号化された散乱破損が、各多重散乱及び単一散乱からの加法的寄与から構成されるというモデル命題に基づいている。従って、誤り訂正のステップＳ６３０において、単一散乱及び多重散乱のＰＳＦが、入力画像ＩＭから減算されるなどして、入力画像ＩＭに適用され、補正された画像ＩＭ＿Ｃが得られる。このように散乱補正された画像ＩＭ＿Ｃは、その後、ディスプレイデバイス上に表示され、更なる参照のために記憶され、又は他のやり方で処理されることができる。
【００９８】
上述のステップＳ６２０に対する変形例として、第１カテゴリのカーネルが発見されると、これらを初期入力画像から減算して、部分的に補正された画像ＩＭ＿Ｃｐを形成することができる。この段階では、画像ＩＭは、単一散乱イベントに関してのみ補正されているので、この画像は部分的に補正されたものである。第２カテゴリのカーネルを発見するためのマッチングは、この部分的に補正された画像ＩＭ＿Ｃｐに対して行われ、第２カテゴリのカーネルが発見されると、これは、部分的に補正された画像ＩＭ＿Ｃｐから減算されて、（現在完全に）補正された画像ＩＭ＿Ｃを生成する。このアプローチは逆にされることもでき、すなわち、第２カテゴリのカーネルが最初に決定され、初期画像ＩＭから減算されて、部分的に補正された画像ＩＭ＿Ｃｐを形成する。次に、この部分的に補正された画像ＩＭ＿Ｃｐが、上述したように第１カテゴリのカーネルを見つけるために使用される。理解されるように、上述のステップＳ６２０及びＳ６３０は、１つ又は複数の固定数のステップで、又は十分な収束が確立されるまで、反復的に実行されることができる。言い換えると、ステップｋにおいて、散乱補正画像ＩＭ＿Ｃｋが取得されると、これは新しい入力画像としてステップＳ６２０にフィードバックされ、次いで、同じやり方で新しい散乱補正画像ＩＭ＿Ｃｋ＋１＝ＩＭ＿Ｃｋ−ΣｊＫ１＿ｊ＋ΣｊＫ２＿ｊが取得され、ここで、カーネルＫ１＿ｊ、Ｋ２＿ｊは、前の繰り返しからのＩＭ＿Ｃｋに対するマッチングにおいて取得される。
【００９９】
上記の実施形態の全てにおいて、ステップＳ６３０で補正を適用する場合、上記の式（３）、（４）に従って減算スキームを使用する代わりに、代替実施形態において乗算スキームＩＭ＿Ｃ＝（ＩＭ２）／（ＩＭ＋ΣｊＫ１＿ｊ＋ΣｊＫ２＿ｊ）も想定され、これは、出願人の米国特許出願公開第２０１０／０１８９３７６号（式１１）に記載されるものと同様である。反復的に行われるとき、乗算スキームは、ＩＭ＿Ｃｋ＋１＝（ＩＭ＊ＩＭ＿Ｃｋ）／（ＩＭ＿Ｃｋ＋ΣｊＫ１＿ｊ＋ΣｊＫ２＿ｊ）と書かれることがあり、この場合、カーネルＫ１＿ｊ、Ｋ２＿ｊは、以前の繰り返しからのＩＭ＿Ｃｋに対するマッチングにおいて得られる。乗算アプローチは、負の値を回避することを可能にし、従って、プロシージャをよりロバストにする。
【０１００】
上記では単一散乱イベントと多重散乱イベントとを区別したが、最大Ｎ（自然数＞１）のイベントから構成される散乱が１つのカテゴリを形成し、第２カテゴリがＮを超える散乱イベントから構成される散乱によって形成される、より一般的な区別を考えることができる。実際、一実施形態では、これは２つのカテゴリのみを考慮するのではなく、Ｍ個の異なるカテゴリを考慮することによって更に一般化することができ、この場合、各カテゴリは、あるレンジ内の散乱イベントを含む。これは、特定のレンジ内の散乱イベントを含む特定の散乱タイプをよりよく考慮できるように、信号モデルを更に詳細化することを可能にする。
【０１０１】
画像処理システムＩＰＳのコンポーネントは、単一のソフトウェアスイート内のソフトウェアモジュール又はルーチンとして実装され、イメージャＩＭに関連するワークステーション又はイメージャのグループに関連するサーバコンピュータなどの汎用計算ユニットＰＵ上で実行されることができる。代替として、画像処理システムＩＰＳのコンポーネントは、分散アーキテクチャに配置され、適切な通信ネットワークに接続されることができる。ＩＰＳは、Ｘ線撮像装置に完全に一体化されることができる。
【０１０２】
ベストマッチするファントム（従って、カーネル）を見つけるために上記で使用された「最適化」又は「オプティマイザ」という用語は、数値最適化スキーム／アルゴリズムに基づいており、その例は上記で説明されていることを理解されたい。しかしながら、これらのいずれかを用いることは、出力が絶対的な意味で最適であることを必ずしも意味しない。１つの理由で、使用されるアルゴリズムは、反復中、大域的な最小値に向かって必ずしも収束しなくてもよく、局所的な最小値に収束してもよい。第２に、ここで反復は、リソースを節約するために、又は何らかの理由で、数回の反復の後に中止されることができる。例えば、１つの中止条件は、予めに規定された数の反復ステップの後に出力に大きな変化がなければ、反復が中止されることができる。更に、上記の最適化は、コスト関数Ｆの最小化に関して定式化されているが、コンポーネントれは制限的ではなく、効用関数の最大化に関して同等の定式化が他の状況又はセットアップにおいて好ましい場合があり、これらは、実施形態において本明細書で等しく想定される。
【０１０３】
ＩＰＳの幾つか又は全てのコンポーネントは、適切にプログラムされたＦＰＧＡ（フィールドプログラマブルゲートアレイ）又はハードワイヤードＩＣチップのようなハードウェアに配置されることができる。
【０１０４】
本明細書で開示されるＩＰＳの１つ又は複数の特徴は、コンピュータ可読媒体内に符号化された回路として／回路を用いて、及び／又はそれらの組み合わせとして、構成され又は実装されることができる。回路は、ディスクリート及び／又は集積回路、特定用途向け集積回路（ＡＳＩＣ）、システムオンチップ（ＳＯＣ）、及びそれらの組み合わせ、機械、コンピュータシステム、プロセッサ及びメモリ、コンピュータプログラムを含むことができる。
【０１０５】
本発明の別の例示的な実施形態において、適切なシステム上で、前述の実施形態のうちの１つによる方法の方法ステップを実行するように適応されることを特徴とするコンピュータプログラム又はコンピュータプログラム素子が提供される。
【０１０６】
従って、コンピュータプログラム素子は、本発明の実施形態の一部でありうるコンピュータユニットに記憶されることができる。この計算ユニットは、上述の方法のステップを実行し又はその実行を誘導するように適応されることができる。更に、計算ユニットは、上記装置の各コンポーネントを動作させるように適応されることができる。計算ユニットは自動的に動作するように、及び／又はユーザの注文を実行するように、適応されることができる。コンピュータプログラムは、データプロセッサの作業メモリにロードされることができる。従って、データプロセッサは、本発明の方法を実行するように装備されることができる。
【０１０７】
本発明のこの例示的な実施形態は、最初から本発明を使用するコンピュータプログラムと、最新の手段をもつ既存のプログラムを、本発明を使用するプログラムに変えるコンピュータプログラムとの両方を包含する。
【０１０８】
更に、コンピュータプログラム素子は、上述の方法の例示的な実施形態のプロシージャを満たすために必要な全てのステップを提供することができる。
【０１０９】
本発明の更なる例示的な実施形態によれば、ＣＤーＲＯＭなどのコンピュータ可読媒体が提示され、コンピュータ可読媒体はその上に記憶されたコンピュータプログラム素子を有し、そのコンピュータプログラム素子は、前述のセクションによって説明される。
【０１１０】
コンピュータプログラムは、他のハードウェアと一緒に又はその一部として供給される光記憶媒体又はソリッドステート媒体などの適切な媒体（特に、必ずしも非一時的媒体である必要はないが、非一時的媒体）上に記憶され及び／又は配布されることができるが、他の形態で、例えばインターネット又は他の有線もしくは無線電気通信システムなどを介して、配布されることもできる。
【０１１１】
しかしながら、コンピュータプログラムは、ワールドワイドウェブのようなネットワークを通じて提示されることもでき、そのようなネットワークからデータプロセッサの作業メモリにダウンロードされることもできる。本発明の更なる例示的な実施形態によれば、コンピュータプログラム素子をダウンロードのために利用可能にするための媒体が提供され、このコンピュータプログラム素子は、本発明の前述の実施形態のうちの１つによる方法を実行するように構成される。
【０１１２】
本発明の実施形態は、それぞれ異なる主題を参照して説明されることに留意されたい。特に、いくつかの実施形態は、方法タイプの請求項を参照して説明され、他の実施形態は、装置タイプの請求項を参照して説明される。しかしながら、当業者は、上記及び以下の説明から、別段の通知がない限り、１つのタイプの主題に属する特徴の任意の組み合わせに加えて、異なる主題に関する特徴間の任意の組み合わせも、本出願で開示されると考えられることを理解するのであろう。しかしながら、全ての特徴を組み合わせて、特徴の単純な合計よりも高い相乗効果を提供することができる。
【０１１３】
本発明は、図面及び前述の説明において詳しく図示され説明されてきたが、そのような図示及び説明は、例示的又は説明的であり、限定的ではないと考えられるべきである。本発明は、開示された実施形態に限定されない。開示された実施形態に対する他の変形は、図面、開示及び従属請求項の検討から、請求項に記載の発明を実施する際に当業者によって理解され達成されることができる。
【０１１４】
請求項において、単語「有する、含む（comprising）」は、他の構成要素又はステップを排除せず、不定冠詞「a」又は「an」は複数性を排除しない。単一のプロセッサ又は他のユニットは、請求項の中で言及される幾つかの項目の機能を果たすことができる。
【０１１５】
特定の手段が相互に異なる従属請求項において言及されるという単なる事実は、これらの手段の組み合わせが有利に使用されることができないことを示すものではない。
【０１１６】
請求項におけるいかなる参照符号も、請求項の範囲を限定するものとして解釈されるべきではない。
【図１】
2021525141.tif 000008
【図１Ａ】
2021525141.tif 000009
【図２】
2021525141.tif 000010
【図３Ａ）】
2021525141.tif 000011
【図３Ｂ）】
2021525141.tif 000012
【図４】
2021525141.tif 000013
【図５】
2021525141.tif 000014
【図６】
2021525141.tif 000015
【国際調査報告】
2021525141.tif 000016
2021525141.tif 000017
2021525141.tif 000018