(57)【要約】解剖学的オブジェクトのコンピュータ支援された四次元（4D）モデリングの方法が、そのオブジェクトの、サイクルを通じた複数の静的状態を表す一組の三次元（3D）モデルを取得することを含む。その一組の3Dモデルに対して4D対応推定が実行され、3Dモデルのどの点が互いに対応することが最も確からしいかを判別する。ここで、4D対応推定は、（ｉ）基準位相を定義する、（ｉｉ）血管指向の対応推定を実行するおよび（ｉｉｉ）4D動きデータの後処理のうちの一つまたは複数を含む。本方法はさらに、波面伝搬アルゴリズムによる自動的な3Dモデリングを含む。
2009504297.tif 000002
【特許請求の範囲】
【請求項１】
解剖学的オブジェクトのコンピュータ支援されたモデリングの方法であって：
解剖学的オブジェクトのゲートされた回転X線投影を取得する段階と；
該ゲートされた回転X線投影から波面伝搬法を使って三次元（3D）血管中心線を自動的に抽出する段階とを有しており、ここで、波面伝搬法は、単一位相の波面伝搬の種々のものにおいて自動的に点を見出すことを含む、方法。
【請求項２】
単一位相の波面伝搬の前記種々のものにおける対応する点を見出すのに反応して、四次元（4D）冠動脈動き場が、前記対応する点の関数として生成されることができる、請求項１記載の方法。
【請求項３】
3D血管中心線を自動的に抽出する前記段階が：
（ｉ）前記ゲートされた回転X線投影を事前フィルタ処理する段階であって、前記ゲートされた投影を諸データセットに分類することを含み、ここで前記ゲートされた投影のデータセットは、すべての心臓サイクルからの、ある所与のゲート点への最近接点の投影を含む、段階と；
（ｉｉ）シード点を発見する段階であって、該シード点は所与の部分体積内で最大の3D血管応答を有するボクセルである、段階と；
（ｉｉｉ）波面伝搬を実行する段階であって、波面伝搬の実行される反復の回数は、（ａ）波面伝搬体積のボクセル解像度からあるいは（ｂ）抽出されたボクセル候補に沿った三次元（3D）応答の減少の解析によって導出される、段階と；
（ｉｖ）抽出された血管候補および対応する下位血管（sub-vessels）について：（ａ）血管端点を発見し、（ｂ）最急勾配をもつ経路に沿って血管中心線を前記シード点まで逆追跡し、（ｃ）クロッピングおよび構造化を行い、ここで、クロッピングおよび構造化は、血管を種々のセグメントに分割し、さらに抽出された中心線の、均一な3D血管応答をもつ諸セクションを決定する、段階と；
（ｖ）弧根を発見する段階であって、該弧根は冠動脈樹の流入ノードに対応する、段階と；
（ｖｉ）関係する血管セグメントを互いにリンクさせる段階であって、対応する後継血管セグメントは、所与の血管セグメントの末端点に幾何学的に最も近い点を選ぶことによって決定され；
（ｖｉｉ）血管セグメントを重み付けする段階であって、各血管セグメントの重み付けは、（ａ）血管セグメントの長さ、（ｂ）3D血管応答および（ｃ）中心線の形および位置を含む一つまたは複数の異なる基準に従って実行される、方法。
【請求項４】
さらに、前記投影データセットは、ECG信号のRピークに対して同じ遅延をもつものである、請求項３記載の方法。
【請求項５】
事前フィルタ処理する前記段階がさらに、前記ゲートされた回転X線投影を、マルチスケール血管度フィルタを使ってフィルタ処理することをさらに含み、該マルチスケール血管度フィルタとは、すべてのスケールのヘシアン行列の固有値の最大値として定義される、請求項３記載の方法。
【請求項６】
事前フィルタ処理する前記段階がさらに、前記投影データセットの幅の約98%の半径をもつ円形マスクを用いて前記投影データセットをクロッピングすることを含む、請求項３記載の方法。
【請求項７】
前記ゲートされた回転X線投影のゲーティングは、同時に記録される心電図（ECG）信号によって実行される、請求項１記載の方法。
【請求項８】
前記ゲートされた回転X線投影を事前フィルタ処理して、該投影を、ECG信号のRピークに対して同じ遅延をもつグループに分類する段階をさらに有する、請求項１記載の方法。
【請求項９】
残留呼吸動きのある前記ゲートされた回転X線投影から最適な心臓位相を判別する段階と；
該ゲートされた回転X線投影から波面伝搬法を使って三次元（3D）血管中心線を自動的に抽出するのが、さらに前記最適心臓位相に応じてである、請求項１記載の方法。
【請求項１０】
3D血管度確率を使って波面伝搬法の速度を制御することをさらに含む、請求項１記載の方法。
【請求項１１】
前記3D血管度確率が、同じ心臓位相のあらゆる血管度フィルタ処理された投影に、考えられているボクセルを前方投影し、二次元（2D）応答ピクセル値を選択し、該2D応答ピクセル値を組み合わせて3D血管度確率にすることによって定義される、請求項１０記載の方法。
【請求項１２】
波面伝搬が、冠動脈に属するボクセルを選択する、請求項１記載の方法。
【請求項１３】
波面伝搬モデルが、二つ以上の単一位相波面伝搬を使って組み合わされた多位相波面伝搬を構築する、請求項１記載の方法。
【請求項１４】
単一位相の波面伝搬の種々のものにおける対応する点を見出す段階と；
四次元（4D）冠動脈動き場を、前記種々の単一位相波面伝搬における対応する点の関数として生成する段階とをさらに有する、請求項１記載の方法。
【請求項１５】
患者の身体部分の、異なる投影方向から記録された複数の回転X線投影を含む投影データセットを生成する手段と、該投影データセットから三次元オブジェクトを再構成するコンピュータ手段とを有しており、該コンピュータ手段は、請求項１記載の方法に従って前記オブジェクトのコンピュータ支援されたモデリングを実行するよう動作するコンピュータ制御を有する、撮像装置。
【請求項１６】
回転X線投影の記録を患者の心臓サイクルに基づいて制御できるECG制御をさらに有する、請求項１５記載の撮像装置。
【請求項１７】
請求項１記載の方法に従ってオブジェクトのコンピュータ支援されたモデリングを実行するための、コンピュータによって実行可能な命令の組を有するコンピュータ可読媒体。
【請求項１８】
解剖学的オブジェクトのコンピュータ支援された四次元（4D）モデリングの方法であって：
そのオブジェクトの、サイクルを通じた複数の静的状態を表す一組の三次元（3D）モデルを取得する段階と；
その一組の3Dモデルに対して4D対応推定を実行し、3Dモデルのどの点が互いに対応することが最も確からしいかを判別する段階とを有しており、ここで、4D対応推定は、（ｉ）基準位相を定義する、（ｉｉ）血管指向の対応推定を実行するおよび（ｉｉｉ）4D動きデータの後処理のうちの一つまたは複数を含む、方法。
【請求項１９】
前記取得する段階が、心臓サイクル全体を通じたすべての静的状態を表す一組の三次元（3D）モデルを取得することを含む、請求項１８記載の方法。
【請求項２０】
前記サイクルが心臓サイクルであり、前記一組の3Dモデルを取得する段階がさらに、3Dモデリング手順を心臓サイクルのいくつかの識別可能な心臓位相について繰り返すことによって取得することをさらに含む、請求項１８記載の方法。
【請求項２１】
識別可能な心臓位相の数が、回転ランの間の最小心拍率および取得フレームレートに依存する、請求項２０記載の方法。
【請求項２２】
前記4D対応推定が、心臓サイクルを通じた血管樹のある部分の動きの推定を可能にする、請求項１８記載の方法。
【請求項２３】
前記基準位相が、前記血管指向の対応推定に先立って定義される所定の安定位相である、請求項１８記載の方法。
【請求項２４】
前記基準位相を定義することが、自動定義または手動定義のうちの一方である、請求項１８記載の方法。
【請求項２５】
前記自動定義が、（ｉ）良好な抽出品質の位相に対応する、低い動きの所望の位相を表す、所与の割合のRRに最も近い3Dモデルまたは（ｉｉ）最長の３本の血管を含む3Dモデルの一方を選択する、請求項２４記載の方法。
【請求項２６】
前記手動定義が、（ｉ）抽出された3Dモデルを目視検査すること、（ｉｉ）該目視検査された3Dモデルから最も好適な心臓位相を手動で定義することおよび（ｉｉｉ）基準位相の該手動定義を用いて4D対応推定を開始することを含む、請求項２４記載の方法。
【請求項２７】
血管指向の対応推定が、各3Dモデルにおける安定点を使って、基準位相におけるあらゆる抽出された血管について独立して実行される、請求項１８記載の方法。
【請求項２８】
最初の血管指向の対応推定については、前記安定点は主分岐であり、血管指向の対応推定のその後の一つまたは複数の反復については、前記安定点は下位分岐点である、請求項２７記載の方法。
【請求項２９】
前記血管指向の対応推定は、（ｉ）任意の血管点の3D座標を血管の弧長λによってパラメータ化し、ここでλは、考えられている位相番号p、考えられている血管番号vおよび血管経路に沿ったボクセル番号iに依存し、（ｉｉ）現在考えられている基準位相血管および現在の標的血管の両方の、所定の間隔を維持する等間隔バージョンを生成し、（ｉｉｉ）安定な弧長基準を与えるよう血管点の座標の低域フィルタ処理を実行し、（ｉｖ）二つの血管を点ごとに比較し、（ｖ）二つの血管の点ごとの比較に基づいて全体的な類似性基準を計算する、請求項２７記載の方法。
【請求項３０】
前記血管指向の対応推定はさらに、源血管および標的位相血管ならびに基準位相以外のあらゆる可能な標的位相について、該血管指向の対応推定のステップ（ｉ）ないし（ｖ）を繰り返すことを含み、さらに、対応する血管のすべての対応する座標を、位相についての添え字および対応する3D点についての添え字をもつ動的動き場配列に保存することを含む、請求項２９記載の方法。
【請求項３１】
4D動きデータの後処理が、諸心臓サイクルを通じて諸点を、逸脱物がないかどうか検査することを含み、特定の位相における弧根の、メジアン位置までの距離が所与の閾値を上回っているのを見出すのに反応して、4D動きデータの後処理はさらに、その心臓位相を4Dモデリングから除外することを含む、請求項１８記載の方法。
【請求項３２】
4D動きデータの後処理が、ある位相に属する各組み合わせの点間のユークリッド距離dを計算し、該距離がある閾値を下回ったら一方の点を破棄することを含む、請求項１８記載の方法。
【請求項３３】
患者の身体部分の、異なる投影方向から記録された複数の二次元投影を含む投影データセットを生成する手段と、該投影データセットから三次元オブジェクトを再構成するコンピュータ手段とを有しており、該コンピュータ手段は、請求項１８記載の方法に従って前記オブジェクトのコンピュータ支援された四次元モデリングおよび動き補償された再構成を実行するよう動作するコンピュータ制御を有する、撮像装置。
【請求項３４】
二次元投影の記録を患者の心臓サイクルに基づいて制御できるECG制御をさらに有する、請求項３３記載の撮像装置。
【請求項３５】
請求項１８記載の方法に従ってオブジェクトのコンピュータ支援された四次元モデリングおよび動き補償された再構成を実行するための、コンピュータによって実行可能な命令の組を有するコンピュータ可読媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本願の実施形態は、概括的には診断画像データからの三次元解剖学的オブジェクトのコンピュータ支援再構成に、より詳細には、自動4D冠動脈モデリングおよび動きベクトル場推定のための方法および装置に関する。
【背景技術】
【０００２】
冠動脈は、造影剤注入後の介入的なX線システムによって撮像できる。冠動脈の動きのため、一組の二次元（2D）投影からの三次元（3D）再構成の生成は、同じ心臓位相に属する限られた数の投影を使ってできるだけであり、これは非常に貧弱な画質につながる。したがって、二つ以上の投影から冠動脈樹（coronary  tree）の3Dモデルを導出するための諸方法が開発されてきた。そうした方法のいくつかは、X線血管造影図の一つにおける初期2D中心線に基づき、同じ心臓位相の他の血管造影図における対応する中心線の点をエピポーラー拘束を利用して検索することに基づいている。結果として、アルゴリズムは呼吸およびその他の残留する非周期的な運きに非常に敏感である。
【０００３】
もう一つの方法は、3Dでの波面伝搬（front  propagation）アルゴリズムに基づくものである。この後者の方法では、波面伝搬を制御するための速度関数が、波面の境界ボクセルが、ある血管に属する確率によって定義される。該確率は、当該ボクセルを、同じ心臓位相の、血管度フィルタ処理された（vesselness-filtered）全投影に前方投影し、応答値をかけ合わせることによって評価される。そのようなアルゴリズムは、異なる血管造影図の間の残留運きの不均一性にそれほど敏感ではないことが注意される。しかしながら、3Dでのそのような波面伝搬アルゴリズムは半自動でしかない。
【０００４】
たとえば、波面伝搬の出発点となる3Dシード点は手動で定義してやる必要がある。各血管についての3D端点を手動で定義してやる必要がある。端点からシード点まで、3D波面伝搬アルゴリズムは、速度関数に関して自動的に最高速の接続経路を探索する。3D波面伝搬アルゴリズムのある側面では、端点は再構成体積の考えられる大きさから導かれる。しかしながら、この基準は特定的でなく、小さく設定されすぎるとアルゴリズムが血管分枝を見逃すことになり、あるいはその値が大きく設定されすぎると波面が血管樹（vessel  tree）体積の境界を越えて伝搬する。たいていの場合、血管樹全体についての上述したアーチファクトを回避する前記基準の単一の値はないというのが普通である。各血管について最適化された、ずっと特定的な基準が必要とされている。
【０００５】
さらに、3D波面伝搬アルゴリズムに関し、種々の血管および血管セグメントの重要度に基づいた検索およびランク付けは、「構造化（structuring）」と称される。3D波面伝搬アルゴリズムの作業フローでは、ユーザーは、手動で特定の血管を選択し、あらゆる血管について手動でシード点および端点を定義することによってランク付けを実行し、そうして手動で「構造化」を達成する。
【０００６】
さらに、3D波面伝搬アルゴリズムは、単一の心臓位相について冠動脈モデルおよび中心線を抽出するのみである。異なる心臓位相からの一組のモデルまたは中心線から四次元（4D）動き場を導出するために、3D中心線上での対応する点を導出するための方法が与えられる必要がある。
【０００７】
図１は、同じ心臓位相においてX線蛍光法によって取得された二枚の二次元（2D）投影１および２からなる診断投影データセットを概略的に示している。心臓位相モニタリングのいかなる好適な型を使うこともでき、たとえばX線投影の取得と並行した心電図（ECG）の記録を使うこともできることを注意しておく。異なる投影角で記録された投影１および２のそれぞれは、患者の分岐した血管３を示す。したがって、投影画像１および２は、異なる視点からの同じ血管３を示している。投影データセットを取得するためには、血管３が投影において黒く見えるよう、造影剤が患者に投与された。
【０００８】
3D波面伝搬法に従って血管３の三次元構造を再構成するためには、初期に再構成体積４内にシード点５が設定される。次いで、それぞれのケースで伝搬基準に従って血管３に属する体積４内の隣接する点を位置特定していくことによって、血管３が体積４内に再構成される。この目的のため、二次元投影１および２それぞれの内部で個別の点５に属する局所領域６および７がそれぞれのケースで個別に数学的解析にかけられる。シード点５に隣接するある点を位置特定したのち、手順は、この点に隣接する諸点について順に反復され、しまいには血管３の全体構造が体積４内に再構成される。
【０００９】
それぞれのケースでそれぞれの伝搬ステップをもって調査された点は、局所領域６および７の数学的解析が、その投影データセットに属する投影（すなわち、この例では投影１および２のそれぞれ）の全部または大半について、陽の結果を与える場合に、血管に属すると識別される。局所領域６および７は、点５を、二つの投影１および２が記録された投影方向に従ってこれら二つの投影の対応する平面に投影することによって決定される。これは図１ではそれぞれ矢印８および９によって示されている。この既知の3D波面伝搬法は、同じ心臓位相の二つの投影に関して記載されてきたが、二つの投影に限られないことを注意しておく。
【発明の開示】
【発明が解決しようとする課題】
【００１０】
したがって、当技術分野における問題を克服するための改良された方法およびシステムが望まれる。
【課題を解決するための手段】
【００１１】
本開示のある実施形態によれば、解剖学的オブジェクトのコンピュータ支援された自動的な四次元（4D）モデリングの方法は、そのオブジェクトのサイクルを通じた複数の静的状態を表す一組の三次元（3D）モデルを自動的に取得することを含む。その一組の3Dモデルに対して4D対応推定が実行され、3Dモデルのどの点が互いに対応することが最も確からしいかを判別する。ここで、4D対応推定は、（ｉ）基準位相を定義する、（ｉｉ）血管指向の対応推定を実行するおよび（ｉｉｉ）4D動きデータの後処理のうちの一つまたは複数を含む。本方法は、撮像システムによっても、コンピュータ・プログラム・プロダクトの形でも実装できる。さらに、本開示のある実施形態に基づく方法は、波面伝搬アルゴリズムを用いた自動3Dモデリングを可能にすることをも含む。
【００１２】
図面において、同様の参照符号は同様の要素を指す。さらに、図面が縮尺通りに描かれていないこともあることを注意しておく。
【発明を実施するための最良の形態】
【００１３】
本開示のある実施形態によれば、波面伝搬法を使ったゲートされた回転血管造影X線投影からの自動3D血管中心線抽出を含む。具体的には、該方法は、ゲートされた3D回転X線投影からの冠動脈中心線樹の自動抽出のための対話的でない、すなわち人間との相互作用がないアルゴリズムを含む。本方法は、冠動脈に属するボクセルを選択するために波面伝搬アプローチを利用する。波面伝搬速度は、3D血管度確率によって制御される。3D血管度確率とは、考慮しているボクセルを、同じ心臓位相の、血管度フィルタ処理された全投影に前方投影し、2D応答ピクセル値を拾い、該値を組み合わせることによって定義される。本方法はさらに、2D応答値を組み合わせて3D血管度確率にする種々の方法を含む。本方法はさらに、いくつかの単一位相モデルを利用して組み合わされた多位相モデルを構築することをも含む。
【００１４】
換言すれば、本方法は、ゲートされた3D回転X線投影からの冠動脈中心線樹を抽出するための完全に自動的なアルゴリズムを含む。該アルゴリズムは、拡張終期の心臓位相における良好な品質の投影を使うときに現実的である。収縮期位相におけるほとんどキスしている血管（almost  kissing  vessels）からの短絡アーチファクトおよびゴースト血管アーチファクトは、波面伝搬アルゴリズムの代替バージョンの使用によって著しく軽減できる。すべてのアルゴリズム・バージョンは、限られた動き補償能力をもつので、最適な心臓位相を見出したのち、残留呼吸運動をもつ投影の中心線抽出が可能である。さらに、最善の心臓位相を決定し、誤って追跡された血管の確率を減らすために単一位相モデルも組み合わされることもできる。さらに、このアプローチを用いた完全な4D冠動脈動き場を生成するために、異なる単一位相モデルにおける対応する点を見出すことができる。
【００１５】
したがって、ここで論じる波面伝搬法により、人間との相互作用なしに冠動脈中心線樹の自動抽出が可能になる。さらに、上記のように、波面伝搬モデルは、特に呼吸によって引き起こされる残留運動に比較的不感である。ある実施形態によれば、一組のECGゲートされたモデルから最小運動の心臓位相における冠動脈形状を表すモデルを決定する必要がある。中心線抽出アルゴリズムでは、アルゴリズムは、波面伝搬アプローチに基づいて完全に自動的な冠動脈中心線抽出を可能にする。
【００１６】
ここで論じるように、自動3D波面伝搬アルゴリズムは、ゲートされた投影を入力として使う。ゲートは、同時に記録される心電図（ECG）信号に従って実行される。該アルゴリズムは複数の準備および解析ステップからなる。それには、（ｉ）ゲートされた投影の事前フィルタ処理、（ｉｉ）シード点の発見、（ｉｉｉ）波面伝搬、（ｉｖ）すべての血管候補について：（ａ）端点の発見、（ｂ）逆追跡および（ｃ）クロッピングおよび構造化、（ｖ）「弧根（root  arc）」の発見、（ｖｉ）リンク、（ｖｉｉ）重み付けならびに（ｖｉｉｉ）出力および出力のためのリンクが含まれる。
【００１７】
〈ゲートされた投影の事前フィルタ処理〉
第一の段階では、投影は、ECG信号のRピークに関して同じ遅延のグループに分類される。ゲートされた投影のデータセットは、すべての心臓サイクルからの、ある所与のゲート点への最近接点投影からなる。アルゴリズムのその後のすべてのステップは、ゲートされた投影セット上で実行される。次のステップでは、投影は、フィルタ幅1ないし7ピクセルのマルチスケールの血管度フィルタを使ってフィルタ処理される。結果は、各ピクセルについて血管に属するかどうかの確率を与える2D応答マトリクスR2Dの組である。マルチスケール血管度フィルタは、あらゆるスケールのヘシアン行列の固有値の最大値として定義される。境界アーチファクトを避けるため、血管フィルタ処理された投影は約(0.98×投影幅)の半径をもつ円形マスクによってクロッピングされることができる。
【００１８】
〈シード点の発見〉
各ボクセル
（外１）
2009504297.tif 000003
について、各投影上の対応するピクセルは、円錐ビーム前方投影を使って計算できる。円錐ビーム前方投影は、nを現在の投影を表し、
（外２）
2009504297.tif 000004
は検出器平面の法線ベクトル、
（外３）
2009504297.tif 000005
は検出器原点、
（外４）
2009504297.tif 000006
は焦点として、各投影についての軌跡データを定義して特徴付けられることができる。
【００１９】
（外５）
2009504297.tif 000007
はボクセル、
（外６）
2009504297.tif 000008
はその投影と考えられる。検出器平面の大きさは、wxおよびwy（mm単位の幅および高さ）ならびにpxおよびpy（ピクセル単位の幅および高さ）によって決定される。
【００２０】
3Dでの検出器平面上の投影されたピクセルは次のように計算される：
【００２１】
【数１】
2009504297.tif 000009
すると、投影上での対応するx,y座標は：
【００２２】
【数２】
2009504297.tif 000010
となる。システム幾何データは各投影に固有なので、ピクセル座標vは現在の投影nにも依存している。
【００２３】
異なる投影の間に動きがないとすると、ボクセル
（外７）
2009504297.tif 000011
が血管内に位置する確率R3Dは、対応する全ピクセルについて2D血管度の結果値R2Dを乗算することによって得ることができる：
【００２４】
【数３】
2009504297.tif 000012
その結果、ある部分体積内で最大の応答をもつボクセルを選ぶことによってシード点が見出される。
【００２５】
現在、このようにして検査される部分体積は、全体積の約11%である。というのも、主血管（理想的には根弧）は、体積の頭蓋側の半分（cranial  half）内かつ中央に位置されるので、部分体積は次のように決定される：
0.25xmax≦x≦0.75xmax
0.25zmax≦z≦0.75zmax        (式4)
0.5ymax≦y≦0.95ymax
ここで、y軸が尾頭方向である。yの値は最大でもymaxに達するべきではない。血管フィルタ処理された投影の残留境界アーチファクトが、適切なシード点の検索に影響することがありうるからである。
【００２６】
さらなる加速のため、各ボクセルについての3D応答値は、N個すべての投影を使って完全に計算されはしない。n個の投影の積を計算したのちに中間的な値が現在最高の応答値を下回れば、残りのN−n個の投影は計算する必要がない。というのも、追加的な乗算のたびに中間的な応答値はさらに減少する一方だからである。これは、ソースデータにもよるが、2ないし5の追加的な加速因子につながる。
【００２７】
〈波面伝搬〉
適切なシード点がみつかったのち、波面伝搬が開始できる。以前に調査された各ボクセルについて、ある特徴的な値が記憶されることになる。その値は、シード点に発して波面がどのくらい「素早く」このボクセルに伝搬したかを示すものである。結果として、この値は時間値と呼ばれ、シード点において0に設定される。したがって、ある任意の経路に沿ったこれらの時間値の増加は、おそらく良好な血管についてはより低く、「悪い」血管およびアーチファクトについてはより高い（より急峻）。
【００２８】
各反復ステップにおいて、現在最低の時間値をもつ波面上のボクセルから始まって、あらゆる隣接ボクセルの3D血管応答値が計算され、その逆数が考えられている開始ボクセルの時間値に加えられる。隣接ボクセルが以前に考慮されていれば、その値は再び計算し直されることはしない。こうして、λ0ステップ後に到達される、ボクセル
（外８）
2009504297.tif 000013
についての時間値
【００２９】
【数４】
2009504297.tif 000014
は、あらゆる先行するボクセルの応答値を含んでいるので、シード点から始まる可能な最良の経路の履歴を表す：
【００３０】
【数５】
2009504297.tif 000015
各ボクセルについての適切な応答値R3Dを計算するためにはいくつかの方法がある。アルゴリズムの全体的な品質は主としてここで使われる手法の品質に依存する。よって、種々の手法が試されてきたが、現実的であると立証されたのは三つだけである。
【００３１】
〈第一の波面伝搬手法（FP1）〉
単純かつ安定な一つの方法は、フィルタ処理された投影それぞれについて対応するピクセルのすべての応答値をかけ合わせることである：
【００３２】
【数６】
2009504297.tif 000016
ここで、nはゲートされた投影をカバーする添え字であり、R2Dは現在のフィルタ処理された投影上での対応するピクセル値であり、該投影の座標は上述したように
（外９）
2009504297.tif 000017
で与えられる。こうして、R3Dはよりよい応答についてはより高くなり、より悪い応答についてはより低くなる。非常に低いR2D応答については乗算は実際上問題にならない。というのも、血管構造を別としても、R2D応答は実際に0に達しはしないからである。
【００３３】
この手法は、セットのほとんど全投影上の血管が同様で、比較的高い品質である場合にはそこそこの結果を与える。弱く、細い血管を追跡するには問題がある。よって、より大きな血管でさえも、端は細くなるので、その実際の末端まで追跡できないことがありうる。波面は「良好な」血管に向かって素早く伝搬するが、次第に弱くなるので、波面の前進はだんだん中立的になっていき、血管の境界（border）に向かって伝搬する傾向がある。したがって、比較的貧弱な品質の投影を使った血管樹全体の納得のいく追跡は、多くの反復を行うことによって多大な計算パワーを消費する（たとえば5123解像度については300ないし500万）。それにもかかわらず、血管の外端（outer  ends）は完全に追跡されないことがありうる。
【００３４】
〈第二の波面伝搬手法（FP2）〉
前節で述べたような細い血管を追跡する問題の解決策は、低応答のボクセルを全然血管上にないことが明らかなボクセルより優先することでありうる。したがって、この第二の波面伝搬手法は、全投影で比較的均一な応答をもつボクセルを、逆投影されたピクセルの応答値の相違がより大きなボクセルに比べて強調することを試みる。「正しい」ボクセルでも、投影によっては、動きまたは悪い投影／前フィルタ処理品質のために悪い応答値をもつことがありうるので、この判定は間違っていることもありうる。あらゆるフィルタ処理された投影は1に規格化されるので、結果は、1より小さな冪乗にすることによって強調でき、1より大きな冪乗にすることによって抑制できる。あるボクセル
（外１０）
2009504297.tif 000018
の2D応答値がどのくらい一様に分散しているかを記述するために、今、指数
（外１１）
2009504297.tif 000019
を、規格化された分散として計算する：
【００３５】
【数７】
2009504297.tif 000020
ここで、
【００３６】
【数８】
2009504297.tif 000021
である。上記の指数は次のように使われる：
【００３７】
【数９】
2009504297.tif 000022
この手法は、弱い血管を優先するが、動き補償能力を低下させる。場合によっては不安定になる傾向がある。
【００３８】
〈第三の波面伝搬手法（FP1）〉
第三の波面伝搬手法は、垂直なビューどうしから抽出された情報を近い角のビューどうしから取られた情報より優先するために、二つの投影mおよびnの間の投影角差αm−αnを考慮に入れるというものである。これは、二つの投影内の奥行き情報の誤解釈を最小化するはずである。三つ以上の投影が利用可能なので、すべての投影(1…n0)は対で考えられ、それぞれの結果は乗算によって組み合わされる。投影の各対についての応答値は、それらの対応する2D応答値をかけ合わせ、それらを投影差角の正弦によって重み付けすることによって計算される：
【００３９】
【数１０】
2009504297.tif 000023
正弦は、体積中心Mから検出器Dに向かう両ベクトルをそれぞれの長さで割ったもののクロス積を計算することによって得られる：
【００４０】
【数１１】
2009504297.tif 000024
この第三の波面伝搬手法は、細い血管を追跡するときによい性能を示し、残留動きを補償する。さらに、第三の波面伝搬手法は、第二の波面伝搬手法より安定でありうる。
【００４１】
〈波面伝搬の終了〉
体積解像度および投影品質に依存して、納得のいく反復回数のおおまかな値がある：
i0,FP1≒0.03×ボクセル数        (式12)
第一の波面伝搬については、2563ボクセルについて、約500k反復が十分であるのに対し、5123なら波面に同様の諸領域に伝搬させるのに約4000k反復を必要とする。しかしながら、この後者の反復回数は約8倍のメモリおよび計算時間を消費する。第二および第三のFP手法は、同様の結果を得るために約半分の反復しか必要としない。
【００４２】
〈血管セグメントの発見〉
端点をみつけたのち、血管中心線が追跡され、クロッピングされ、その諸部分が別個に記憶される。相続く血管が同じ仕方で扱われる。次の三つの段階（１）端点の発見、（２）逆追跡および（３）クロッピングおよび構造化は、よって、各血管候補およびその部分血管それぞれについて行われる。
【００４３】
（１）端点の発見
波面伝搬が終了したのち、あらゆる血管について、適切な端点をみつける必要がある。これは、全体積をn3個の部分体積に分割することによって達成される。ここで、この段階ではn＝50である。各体積内で、最高の時間値をもつボクセルが選ばれる。このボクセルは、血管の外端上に位置している。というのも、波面は各血管の中心で素早く伝搬し、次いでその境界付近でゆっくりと広がる（高い時間値を引き起こす）からである。
【００４４】
（２）逆追跡
逆追跡は、最急勾配法を使って実行される。端点が与えられると、逆追跡は、現在のボクセルに対して時間値の減少が最も大きいボクセルのほうに向かう。ステップ毎に最大の減少をたどっていくことにより、シード点に戻る最適な経路が計算される。波面伝搬の表面で始まり、直接、血管中心に、次いで中心線に沿ってシード点へとつながる。経路が前に以前の反復によってすでに追跡されている場合、その経路は再び追跡されはしない。これは、追跡されるボクセルが、各側に２ボクセルの追加的な安全領域とともにマークされる3Dビットマップによって管理される。これは同様の（平行な）経路をダブって追跡することを防止する。
【００４５】
（３）クロッピングおよび構造化
血管の境界に位置するボクセルは中心線に属さず、よってそのようなボクセルはクロッピングされる必要があることを注意しておく。クロッピングは、再帰的なアルゴリズムによって行われる。ここで、再帰的アルゴリズムの課題は、追跡された中心線を異なる品質の諸セグメントに分割することである。逆追跡が始まった点でのセグメントは、最悪の品質をもち、よって消去される。
【００４６】
再帰的クロッピング・アルゴリズムは、あらゆる血管の品質が、シード点付近で最良で、開始点の逆追跡につれて低下すると想定している。現在の血管ボクセルの最初の四分の一の平均値が計算され、ここで計算された値がその後、追跡開始点に向かって走査する際に閾値として使われる。該閾値は、時折何度か超過されてもよいが、超過ボクセルの数が許容値（たとえば、最大で連続10回まで）を超えると、その特定の点は著しい品質違反であると考えられ、血管は二つの部分に分割される。これは、最悪品質のセグメントは、よりよい品質の血管セグメントから切り去られ、独立した血管として保存されるということを意味している。この第二の血管はその後、同じようにして扱われ、独立した血管についてのセグメントが分離される、などとなる。再帰アルゴリズムは、残りの部分がある最小長（たとえば10ボクセルのオーダー）より短くなると中止される。追跡開始点に位置する境界ボクセルは、最小長基準によって切り去られるか、その長さが10ボクセルを超えていれば、のちに論じる重み付けアルゴリズムによって無視できると評価される。
【００４７】
〈「弧根」の発見〉
先述したように、波面伝搬のためのシード点は必ずしも弧根に対応しない。弧根とは、冠動脈樹の流入ノードである。結果として、あらゆる血管はこの「間違った」出発点に逆追跡される。弧根の真の位置を推定するために、最長の３本の単一の血管セグメントの最も頭寄りの（cranial）点が使われる。シード点とこの新しい最上点との間のリンク血管セグメントは、その後、必要なら他の血管を拡張するのに使われる。
【００４８】
〈リンク〉
ここまでは、血管は互いに関係をもたなかった。各血管の終わりは、次の三つの理由のうちの一つによって引き起こされる：ｉ）弧根に到達した、よってリンクは必要とされない；ｉｉ）その血管が以前により長い血管の一部であり、上記のクロッピングおよび構造化アルゴリズムによって分離された；およびｉｉｉ）分岐がある。この分岐は、逆追跡段階で検出されたもう一つの血管交差があるということを意味している。ここまでは、経路が以前にたどられたかどうかだけが知られており、どの血管がそれを使うかは知られていなかった。正しい後継血管は、あらゆる血管セグメントの端点に幾何学的に最も近い点を選ぶことによって決定される。逆追跡段階ではすべての血管が昇順にインデックス付けされたので、考えられている血管よりも小さなインデックスの血管上の点を検索することが必要なだけである。リンク後、あらゆる血管の全長（端点から弧根まで）は、リンク経路に沿ったすべての血管セグメントの長さを足すことによって簡単に計算できる。
【００４９】
〈重み付け〉
上記のステップにおいて、多数の経路が抽出されたが、そのうちほんの若干のみが本当に存在している血管を表現する。大半の経路は、投影品質の不足、残留動き、遠近法効果などといったアーチファクトによって引き起こされる。したがって、それらのうちのどれが本物の血管を表していることが最も確からしいかを判別しなければならない。抽出された経路候補の全体的な有意性のための尺度Sはいくつかの因子から構成されることができる：ｉ）血管セグメントの長さまたは全長、ｉｉ）時間値によって決定される品質、ｉｉｉ）3D位置（おそらくは所定のモデルの支援とともに）およびｉｖ）形である。有意性値Sに従って、すべての経路候補をソートできる。これにより、最も有意な経路を出力のために選ぶことができる。出力する経路の最大数はシステム・ユーザーが設定できる。有意性値Sの計算は、いまだ改善の余地がある。というのも、ここでの誤判断は誤った（「ゴースト」）血管の出力につながりうるからである。ある実施形態では、Sは次のように計算される：
【００５０】
【数１２】
2009504297.tif 000025
ここで、Yendおよびyroot_arcは、上記のようにして決定された、それぞれ現在の血管セグメント端点および弧根のy座標（尾頭回転軸に沿った）である。量lpartは、ボクセル単位での血管セグメントの長さであり、
（外１２）
2009504297.tif 000026
は血管セグメントの端点の時間値である。たとえば勾配基準を使って、そこそこの数の抽出可能な血管中心線を自動的に推定することが可能となりうる。
【００５１】
〈出力および出力のためのリンク〉
中心線データをファイルに保存するとき、リンクを検査して、血管のいくつかの部分をリンクし直すことが必要となることがありうる。というのも、リンクされた経路の一つまたは複数のセグメントが出力のために選択されていないことがありうるからである。
【００５２】
本開示のある実施形態によれば、改良された波面伝搬アルゴリズムが、半自動的な3Dアルゴリズムの従来既知の方法を、完全に自動的な4Dアルゴリズムに変換する。本方法は、これまでに論じたさまざまな問題に対処し、次のような解決策を提供する。
【００５３】
１．シード点：ある実施形態によれば、シード点は、あらゆる血管造影図において観察可能な3D体積の中央の置かれた頭部（cranial）部分体積における上述の3D血管応答を評価し、最大の3D応答をもつ点を選択することによって自動的に定義される。対応する3D応答のX線投影の取得と並行して、任意の好適な型の心臓位相モニタリングを使うことができる。たとえば、心臓位相モニタリングは、心電図（ECG）の記録を含んでいてもよい。最大3D応答点は、血管樹の上に位置されるが、必ずしも主たる分岐の流入ノードのところではない。代替的な方法は、上述の体積の表面の頭部（cranial  part）上で最大の3D応答をもつ点を選択することである。この後者の例では、これは、大動脈経由で頭側（cranial  side）からはいってくる、造影剤で満たされたカテーテル上に位置されたシード点を与える。
【００５４】
２．波面伝搬の停止：波面伝搬の実行される反復の回数は、（ｉ）波面伝搬体積のボクセル解像度から、または（ｉｉ）抽出された血管に沿っての3D応答値の減少を解析することによって導出される。
【００５５】
３．端点：血管の潜在的な端点は、一つまたは複数の異なる方法によって自動的に決定できる。第一の実施例では、波面伝搬体積は多数の部分体積に分割される（たとえば503すなわち50×50×50）。あらゆる部分体積内で、最も遅い波面到着をもつ点が、逆追跡アルゴリズムのための開始点として選択される。逆追跡アルゴリズムは、最急勾配をもつ経路に沿って速度場を逆向きにたどってシード点に向かう。第二の実施例では、波面伝搬の間に、アルゴリズムは最急勾配に沿った経路を追跡し、3D血管応答の顕著な低下が検出されたら停止する。いずれにせよ、潜在的な血管端点の精密な推定はきわめて決定的なわけではない。というのも、続く構造化のステップにおいて、血管セグメントが解析され、その有意性に基づいて重み付けされるからである。
【００５６】
４．構造化：血管は、動的な構造化アルゴリズムによって種々のセグメントに分割される。動的な構造化アルゴリズムは、抽出された中心線の、均一な3D血管応答をもつセクションを決定する。各血管セグメントの重み付けは、種々の基準に従って実行される：（ｉ）長さ、（ｉｉ）3D血管応答（品質に対応）、（ｉｉｉ）中心線の形および位置（または任意的に、先験的な冠動脈モデルに基づく）。最も有意な重み付けされた諸血管が自動的に選択され、3Dアルゴリズムの出力をなす。図２は、本開示のある実施形態に基づくモデリング方法を用いて得られた、ある根底にある心臓位相の二つの投影（２２および２４）に逆投影された、完全に自動的に抽出された3D中心線の例（２０）を含んでいる。
【００５７】
〈4Dアルゴリズム〉
本開示のある実施形態によれば、自動的な4D冠動脈モデリングおよび動きベクトル場推定法は、上記の手順をあらゆる識別可能な心臓位相について繰り返すことによって心臓サイクル全体を通じたあらゆる静的状態を表す一組の3Dモデルを入力において必要とする。本方法は、異なるモデルどうしの分岐およびその他の形状属性を照合することによって、異なるモデルの対応する点を決定する。4D情報を活用する可能性のある用途は、ゲートされた、あるいは動き補償された3D再構成のための最適な心臓位相を導出することである。
【００５８】
本開示の実施形態に基づく方法は、冠動脈中心線抽出およびモデリングのための、完全に自動的で堅牢な4Dアルゴリズムを提供する。本方法は、同じ心臓位相の血管造影図どうしにおける、残留動きに起因する不一致に対処することができる。さらに、本開示の実施形態に基づく方法は、従来の既知の3D波面伝搬アルゴリズムに対する改善を提供し、そのような改善は、4D動き補償された再構成およびモデリングのような新たな用途を可能にする。
【００５９】
心臓サイクル全体を通じたあらゆる静的状態を表す一組の3Dモデルは、あらゆる識別可能な心臓位相について3Dモデリング手順を繰り返すことによって得られる。回転ラン（rotational  run）の間の最小心拍率rh,min（拍毎分(beats  per  minute)、bmp）および取得（acquisition）フレームレートfa（1/s単位）に依存して、識別可能な心臓位相の数pNは次に等しくなる：
pN＝（60s/min）（fa/fh,min）
これは、pN個の独立な3Dモデルが生成されたということを意味する。この値は、取得フレームレートfaが25fps（フレーム毎秒）で心拍率fhが100bpm（拍毎分）のときの約15からfaが30fpsでfhが45bpm（拍毎分）のときの約40まで範囲がある。4D対応推定の課題は、モデルのどの点が最も互いに対応しあいそうかを決定することである。これは、心臓サイクルを通じて血管樹のある部分の動きを推定することを可能にする。血管の長さ方向の動きおよび3Dモデリングプロセスの間に引き起こされた曖昧さのような問題は、4D対応推定を難しくするものであり、考慮に入れる必要がある。対応推定は、次のステップを実行することによって遂行される：
１．基準位相（安定位相）の定義
２．血管指向の対応推定
３．4D動きデータの後処理。
【００６０】
１．基準位相の定義
安定な4D対応を推定するために、諸ステップの間に抽出される多くの潜在的な血管構造のうち、どれが心臓サイクル全体の間で最も重要かを決定することが必要である。3Dアルゴリズムの間に、血管セグメントは、その推定される有意性に従って重み付けされるが、これはあらゆる単一の3Dモデルについて独立になされる。そのことは、異なる心臓位相で抽出された血管の揺動につながる。したがって、あらゆる所望の血管が抽出される基準位相pr（安定位相）が、対応推定に先立って定義されねばならない。これは、自動的または手動でできる。
【００６１】
自動定義：35%RRに最も近い位相、すなわち動きの低い位相、よって良好な抽出品質の位相である可能性の高い位相を表す3Dモデルが選ばれるか、３つの最長の血管を含むモデルが選ばれるかする。RRとは、ECGの２つの相続くRピークによって定義される時間期間を表すことを注意しておく。ここで、ECGではRピークが優勢であり、各Rピークは心臓の収縮に先立つ電気インパルスを表す。
【００６２】
手動定義：あらゆる抽出された3Dモデルの目視検査（たとえば、図３に示されるような、すべてのモデルの投影をもつ概観プロット３０を使う）に従って、手動で最も好適な心臓位相を定義し、アルゴリズムを再スタートさせることができる。図３は、異なる心臓位相において抽出された血管の二つの投影３０の例を示す。上の列３２は43.5%RRの心臓位相を表しており、三つの正しく抽出された血管を示している。これは、その位相に、潜在的な基準位相としての資格を与える。一方、下の列３４（5%RR）に示される血管の品質はより悪くなっている。
【００６３】
２．血管指向の対応推定
対応推定は、各モデルにおける一つの安定点を使って、基準位相prにおけるあらゆる抽出された血管について独立して実行される。安定点の役をするのは、初めてこのステップを実行するときには主分岐（「弧根」）であるのに対し、のちの反復では、おそらくはより高い精度の、下位分岐点が使われる。アルゴリズムは、心臓サイクルの間、血管の弧長λが著しく変化しない（全体で2%未満）という事実を活用する。任意の血管点の3D座標
【００６４】
【数１３】
2009504297.tif 000027
は、血管の弧長λによってパラメータ化される。このλは、考えられている位相番号p、考えられている血管番号vおよび血管経路に沿ったボクセル番号iに依存する：λ＝λ(p,v,i)。以下の文章で血管全体をいう場合には、ボクセル番号iは省略される。
【００６５】
現在考えられている基準位相血管λ(pr,vr)および現在の標的位相血管λ(p,v)の両方の、所定の間隔s（現在、2mmに設定されている）を維持する等間隔バージョンが生成される。等間隔バージョンを生成するのは、元の3Dモデルの点から点までの距離は、対角方向のボクセル距離およびリンク・ギャップのため、因子√3またはそれ以上変化があるからである。それらの血管は、安定点から血管の末端までの経路全体を表している。波面伝搬のボクセル表現に由来する量子化効果を解消し、安定な弧長基準を提供するために、等間隔化に先立って、血管点の座標が低域フィルタ処理される。血管λ(p,v)の低域バージョンはλ′(p,v)によって表される。二つの血管は点ごとに比較され、全体的な類似性基準Cが計算される：
【００６６】
【数１４】
2009504297.tif 000028
類似性基準Cが小さいほど、二つの現在の血管がよりよく対応していることを示す。よって、Cが最小になる組み合わせの血管どうしが等価と考えられる。この手順は、源血管vrおよび標的位相血管vならびにあらゆる可能な標的位相p≠prについて繰り返される。対応する血管のすべての対応する座標が最終的に、添え字[0..pN-1]（位相）および[0..imax-1]（対応する3D点）をもつ動的配列A(p,i)に保存される（動き場と呼ばれる）。
【００６７】
３．4D動きデータの後処理
対応推定手順の間、あらゆる対応血管の表現は基準点（通常は弧根）から始まる。これは、血管樹のいくつかの部分が複数回表現されることになる。これは、高い局所点密度につながる。これは特異点およびその他の曖昧さを避けるために薄めてやる必要がある。削減は、ある位相に属する各組み合わせの点間のユークリッド距離d
【００６８】
【数１５】
2009504297.tif 000029
を計算し、該距離がある閾値を下回ったら一方の点を消去することによって達成される。閾値は、t＝0.5s＝1mmとして定義される。
【００６９】
すべての心臓サイクルを通じての、結果として得られる対応する「弧根」点は、逸脱物（outliers）がないかどうか検査されることができる。特定の位相における弧根の、メジアン（または平均）位置までの距離が所与の閾値を上回っていれば、この心臓位相はモデルから除外される。同様にして、他のすべての分岐および単独点が扱われることができる。
【００７０】
ここで図４に目を転じると、図示されている撮像装置はCアームX線装置であり、保持器１１によってたとえば天井（図示せず）から懸架されるCアーム１０を有している。X線源１２およびX線画像変換器１３はCアーム１０上で可動的に案内される。Cアーム１０の中心にあるテーブル１４上に横たわる患者１５の複数の二次元投影X線画像が、異なる投影角にて記録されうる。X線源１２およびX線画像変換器１３の同期的な動きは、制御ユニット１６によって制御される。画像記録の間、X線源１２およびX線画像変換器１３は、患者１５のまわりを同期的に移動する。X線画像変換器１３によって生成される画像信号は、制御された画像処理ユニット１７に送信される。患者１５の心拍は、ECG装置１８を使ってモニタリングされる。ECG装置１８は、画像処理ユニット１７に制御信号を送信する。それにより、画像処理ユニット１７は、冠動脈の血管造影調査を実行するために、それぞれ心拍サイクルの同じ位相にある複数の二次元投影を保存する位置にくる。画像処理ユニット１７は、プログラム制御を有しており、それにより、3D波面伝搬法に従って、こうして取得された投影データセットで検出された血管樹の三次元モデルが実行されることができる。さらに、画像処理ユニット１７は、さらなるプログラム制御を有しており、それにより、4Dモデリングが、本開示の諸実施形態に従って実行できる。次いで4Dモデリングは、一つまたは複数の再構成された血管とともに、画像処理ユニット１７に接続されたモニタ１９上に任意の好適な仕方で視覚化されうる。
【００７１】
上記では若干の例示的な実施形態しか詳細に記載されていないが、当業者は、本開示の実施形態の新規な教示および利点から実質的に外れることなく例示的な実施形態には数多くの修正が可能であることをすぐ認識するであろう。たとえば、本開示の実施形態は、心臓血管のような他の周期的に動く構造に、あるいはより一般に樹状構造に適用されることができる。したがって、そのようなすべての修正は、付属の請求項において定義される本開示の実施形態の範囲内に含められることが意図されている。請求項においては、手段プラス機能（means-plus-function）節は、挙げられている機能を実行するものとして、ここに記載される構造ならびに構造的等価物のみならず等価な構造をもカバーすることが意図されている。
【００７２】
さらに、一つまたは複数の請求項において括弧内に参照符号があったとしても、それは請求項を限定するものと解釈してはならない。「有する」および「含む」などの語は、請求項または明細書全体に挙げられている以外の要素やステップの存在を排除しない。要素の単数形での言及は、そのような要素の複数の言及を排除せず、要素の複数形での言及は、そのような要素の単数の言及を排除しない。実施形態の一つまたは複数は、いくつかの相異なる要素を有するハードウェアによって、および／または好適にプログラムされたコンピュータによって実装されてもよい。いくつかの手段を列挙する装置請求項では、それらの手段のいくつかはハードウェアの同一の項目によって具現されてもよい。ある種の施策が互いに異なる従属請求項において挙げられているというだけの事実が、それらの施策の組み合わせが有利に使うことができないことを示すことはない。
【図面の簡単な説明】
【００７３】
【図１】二つの二次元（2D）投影画像からなる診断投影データセットを概略的に示す図である。
【図２】本開示のある実施形態に基づくモデリング方法を用いて得られた、ある根底にある心臓位相の二つの投影画像に逆投影された完全に自動的に抽出された3D中心線の例を示す図である。
【図３】本開示のある実施形態に基づくモデリング方法を用いて得られた、二つの異なる心臓位相における抽出された血管の三つの直交軸に沿った投影の例を示す例示的な図である。
【図４】本開示のもう一つの実施形態に基づく撮像装置の部分ブロック図である。
【図１】
2009504297.tif 000030
【図２】
2009504297.tif 000031
【図３】
2009504297.tif 000032
【図４】
2009504297.tif 000033
【国際調査報告】
2009504297.tif 000034
2009504297.tif 000035