JP4533587B2 - Medical image laminating device - Google Patents

Description

【００２６】
画像の変形は、ずれ変形と拡大・縮小と、さらにＸ線撮影特有の周辺部変形により複合的に決まる。図６には、ずれ変形に関するパラメータを示している。左右方向のずれ変形はｈ_ｘ、上下方向のずれ変形はｈ_ｙで与えられる。
【００２７】
【数２】

【００２８】
図７に拡大・縮小に関するパラメータを示している。左右方向の拡大縮小変形はＳ_ｘ、上下方向の拡大縮小変形はＳ_ｙで与えられる。
【００２９】
【数３】

【００３０】
図８に周辺部変形に関するパラメータを示している。画像の上半分と下半分とにそれぞれ独立した変形量ｐ_ｕ、ｐ_ｌが設定される。
【００３１】
【数４】

逆変換作用素Ｐ⁻ _ｉは、次の式で表される。
【００３２】
【数５】

【００３３】
この周辺部変形方法では、画像の上半分と下半分に独立した変形量ｐ_ｕ、ｐ_ｌが設定される。これにより画像のひずみが複雑で画像の上と下でひずみ方が異なる場合であっても、良好に位置合わせを行うことができるという効果がある。また、この方法は画像の変形が左右方向の拡大・縮小のみであるため、最適化の際、変形量と２つの画像の縦方向の平行移動量とが独立になり、さらに離れたステージ間の変形量も独立であるため、最適化が容易になり高速に位置あわせができるようになるという効果がある。
【００３４】
次に、以下の説明で使用する各処理段階での座標系について説明する。
図９には初期的な画像座標系（２次元）として、撮影したｍ×ｎ画素の画像内の画素の位置を表す座標系を示しており、実際には負の値や非整数値も扱うことができる。図１０には、変形後画像座標系として、画像に周辺部変形を施した後の座標系を示しており、もちろん周辺部変形を施さなければ画像座標系に一致する。
【００３５】
複数のステージの画像を位置合わせする際に、これらの画像を位置を変えながら配置していくための仮想的で全体的座標系を用いる。本実施形態ではこれをグローバル座標系と呼ぶ。図１１は、グローバル座標系に各ステージの画像を配置した様子を示している。配置した各画像を合成した（貼り合わせた）パノラミックビュー画像の領域もグローバル座標系内に設定される。
【００３６】
これら３種の座標系の間では当然にして座標変換処理が必要とされる。ここでは、各２次元座標系の間の変換をアフィン変換にて表記する。従って位置を表すベクトルを、ｘ，ｙ座標に加えて常に値が１である要素を加えた３次元のベクトルとして扱う。
【００３７】
グローバル座標系の座標ｙ_０、ｉ番目(ｉ＝１，２，・・・Ｎ_ｓ)のステージの変形後座標系の座標をｙ_ｉとしたとき、各座標系の関係は、
【数６】

で表される。Ｅ_ｉはステージｉ−１変形後座標系からステージｉ変形後座標系への座標変換行列であり、マニュアル平行移動Ｍ_ｉ、自動位置合わせ平行移動Ｔ_ｉ、ずれ変形Ｈ_ｉ、拡大・縮小Ｓ_ｉ、画像左上隅から画像中心への平行移動Ｃ_ｉを用いて、
【数７】

で表される。変形後画像座標系から変形前の画像座標系への座標変換は、
【数８】

で表される。ｘ_ｉは変形前の画像座標系での座標値である。座標変換の各要素は次の通りである。
【００３８】
【数９】

Ｌは正方形画像の幅である。Ｐ_ｉは非線形作用素であり、
【数１０】

により定義される。これらより、ステージｉ−１変形後座標系からステージｉ変形後座標系への座標変換行列Ｅ_ｉは
【数１１】

となる。隣接するステージ間の座標変換を、変換作用素Ｄ_ｉを用いて表すと
【数１２】

となる。ｘ_ｉはｉ番目のステージの変形前の画像座標系の座標、ｘ_ｉ＋１はｉ＋１番目のステージの変形前の画像座標系の座標である。変換作用素Ｄ_ｉは、次の通り表される。
【００３９】
【数１３】

【００４０】
次に実際に自動位置あわせ処理について説明する。自動位置あわせ処理には、その前処理として「自動ＲＯＩ設定」と、「フィルタ処理」とが実行される。
【００４１】
まず、自動ＲＯＩ設定について説明する。
図１２には画像有効データ領域とＲＯＩ設定有効領域を示している。ｍ×ｎ画素の長方形の画像データの領域のうち撮影したＸ線画像の有効なデータが存在する領域である。ここでは画像有効データ領域は、画像座標系にて中心点（ｍ／２，ｎ／２）を中心とする半径Ｒvalの円内が画像データ有効領域とする。
【００４２】
テストプログラムでは、自動位置合わせ、自動ＲＯＩ設定、フィルタリングなどの全ての処理の前（画像読み込み直後）に、全てのステージの画像は、各画像毎に画像データ有効領域内の画素値の平均値を各画素から引き、平均値が０になるように正規化されており、画像データ有効領域の外の画素は、画素値０に置き換えられる。これは、テストプログラムではマニュアルで位置合わせする際にステージ間の画像を重畳表示しているためである。重畳表示が自然に行われるために、上側のステージはそのまま表示し、下側のステージは画像の平行移動や変形を行なった後、画素値の正負を反転して表示し、上下のステージの画像の重複する領域は下側画像の反転した画素値を上側画像に加えている。このとき、血管以外のバックグラウンドの領域は０に近い値であると重ね合わせが自然に行なわれる。上下のステージの画像を重ねて表示しない場合には、画像読み込み時に画素値の平均値を引く処理は不要である。この場合画像データ有効領域の外の画素は平均値に置き換えた方が良い。
【００４３】
ＲＯＩ設定有効領域は、画像データ有効領域よりも幾分小さい半径Ｒroi内の円状の領域である。この領域は自動ＲＯＩ設定処理にて用いられる。テストプログラムでは入力画像のマトリクスサイズが５１２×５１２であることを仮定しており、Ｒvalには２３２（ｐｉｘｅｌ），Ｒroiには２１５（ｐｉｘｅｌ）という値を用いている。この値は、画像サイズ毎に変える必要がある。
【００４４】
図１３には、しきい値処理による自動ＲＯＩ設定について示している。自動位置合わせアルゴリズムでは隣接するステージ間の画像の対応する点どうしの画素値を比較し、この画素値の差が小さくなるように上下の画像の間で位置の対応関係を求める。画素値の差を計算する評価点は、血管の近辺のみに限定しなければ血管以外の背景部のコントラスト揺らぎの影響を受けてしまうため安定に位置合わせすることが出来ない。また、高速化のためにも血管領域の周辺のみに評価点を限定することが望ましい。自動ＲＯＩ設定処理では位置合わせ処理の前処理として、位置合わせのための関心領域(ＲＯＩ)を設定する処理を行う。
【００４５】
位置合わせ処理では、２つのステージの上側の画像に設定されたＲＯＩ内部の各点に対応する下側画像内の点を求め、両者の画素値の差が小さくなるようにする。その前処理として、自動ＲＯＩ設定では上側画像の血管の周囲にＲＯＩを設定する。
【００４６】
入力は上側画像と初期位置合わせパラメータ、出力は設定されたＲＯＩの形状である。設定されるＲＯＩは、上下のステージの画像を初期位置合わせパラメータの平行移動量を用いて重ね合わせたとき、上下の両方の画像にてＲＯＩ有効領域内となる領域のみに設定される。ＲＯＩ有効領域はＸ線画像の有効データが存在する円より小さめに設定された、円状の領域である。
【００４７】
図１４には自動ＲＯＩ設定の処理手順を示している。２つのステージの間の自動ＲＯＩ設定は次の手順で実行される。
１．上側ステージの画像の画像データ有効領域内の画素値の平均値ａと分散ｖ^２を求める。
２．上側ステージ画像にて画素値が、（ａ−ｎ_ｔｖ）以下の領域を閾値処理により抽出する(領域Ａ)。
ｎ_ｔは閾値を決定する定数でテストプログラムでは２を用いている。従って、閾値は平均値より２ＳＤ（標準偏差の２倍）だけ小さい値となる。
【００４８】
３．上側ステージのＲＯＩ設定有効領域と下側ステージのＲＯＩ設定有効領域との論理積演算を行ってＲＯＩ設定有効領域を決定する（領域Ｅ）。
４．領域ＡとＲＯＩ設定有効領域との論理積演算を行う（領域Ｂ）。
５．領域Ｂをｎ_ｄピクセル分拡大する（領域Ｃ）。
ｎ_ｄは拡大する大きさを決定する定数で、テストプログラムでは５ｐｉｘｅｌを用いている。
５．領域Ｃと上側ステージのＲＯＩ設定有効領域との論理積演算を行う（領域Ｄ）。
【００４９】
テストプログラムでは、このように抽出した領域Ｄを、ランレングスエンコーディングしてメモリに記憶している。自動位置合わせ処理の際は、このランレングスエンコーディングされた領域Ｄに対応するＲＯＩデータを参照して処理が行われている。これにより、自動位置合わせの際に全部の画素をスキャンする必要がなくなるため、自動位置合わせの処理が高速に行われるようになるという効果がある。
【００５０】
次に、もう一つの前処理としてのフィルタ処理について説明する。
図１５にはフィルタ処理について示している。自動位置合わせ処理ではＲＯＩ内の対応する画素値どうしを比較し、画素値の差が小さくなるように位置合わせする。従って、位置合わせ前の初期状態で上下のステージの画像に血管の重なりが無いと、正しく位置合わせされた状態に収束しない。このような結果になる確率を小さくするために、自動位置合わせの前処理として、元画像をぼかすフィルタ処理（ローパスフィルタ処理）を行うことが効果的である。
【００５１】
フィルタのカーネルは、ｎ_ｆ×ｎ_ｆであるが、演算回数の低減のため、縦方向と横方向の畳み込みに分解できるように、Ｋ_ｉｊ＝Ｋ_ｉＫ_ｊで表されるカーネルを採用した。Ｋ_ｉは分解した１次元のフィルタ係数で、次に示す三角形状の係数を用いた。
【００５２】
【数１４】

【００５３】
このカーネルを用いると、ローパスフィルタ処理は２段階の１次元の畳み込みに分解できる。
【００５４】
（１）横方向の畳み込み
【数１５】

（２）縦方向の畳み込み
【数１６】

は、それぞれ、元画像、横方向の畳み込み結果画像、フィルタ処理結果画像の（ｉ，ｊ）の画素の画素値である。
【００５５】
全ステージの画像の前面にフィルタ処理を施すのではなく、位置合わせで参照される可能性のある部分のみにフィルタ処理を限定することで、フィルタリングの所要時間を大幅に削減することが出来る。テストプログラムではフィルタリングを行う領域を図１６に示すように以下の２つの部分に限定している。
【００５６】
・抽出されたＲＯＩ領域（領域Ｄ）
・フィルタリング有効領域と上側画像のフィルタリング有効領域の積領域
フィルタリング有効領域は、画像データ有効領域より一回り大きい円領域であり、その半径はテストプログラムでは画像データ有効領域の半径より１０ｐｉｘｅｌ大きい値（２４２pixel）を用いている。
【００５７】
以上の前処理完了後に、「自動位置合わせ処理」が実行される。全ステージ数Ｎ_ｓの位置合わせは、隣接する２つのステージの間ごとに行われる。２つのステージ（番号ｉとｉ＋１）の画像に加える位置合わせ・変形は、次の４個のパラメータで規定される。
【００５８】
【数１７】

【００５９】
自動位置合わせの際には、上下のステージの画像にこれら６つのパラメータ（ｐ^ｌ、ｔ_ｘ、ｔ_ｙ、ｐ^ｕ、α、β）による変換（位置合わせ、変形、コントラスト変換）が加えられた後、画素値が比較される。そして、画素値の差が小さくなるようなパラメータを求めることで位置合わせ処理が行われる。
【００６０】
ここで画像左上隅から画像中心への平行移動ｃ_ｉｘ、ｃ_ｉｙは、上側ステージの変形に関わるが、この値はステージｉ、ｉ＋１の間の位置合わせだけでなく、ｉ−１、ｉの間の位置合わせにも関わり、これらの値を位置合わせパラメータとして着目している２ステージの位置合わせを行うと、それ以外のステージ間の位置合わせ結果を壊してしまうことになる。本実施形態の位置合わせ方法はそのようなことが無いように隣接２ステージ毎に位置合わせする方法を採用しており、ｃ_ｉｘとｃ_ｉｙの値は一定値０として位置合わせ過程で変更されることが無いようにしている。すなわち、
【数１８】

【００６１】
上述したように、自動位置合わせ処理では設定されたＲＯＩ内の全ての画素を、次のステージの対応する画素値と比較するのではなく、ＲＯＩ内の全ての画素からランダムにピックアップした画素について、次のステージの対応する画素値と比較する。これにより自動位置合わせの高速化がはかられている。ＲＯＩ内の画素から評価点をピックアップする処理を図１７に示す。
【００６２】
図１７の各升目は、ステージｉの画像の各画素を示し、太線で示した枠がＲＯＩ内の領域を示す。ＲＯＩの各画素は、画像の横方向にスキャンされていく。これらの画素をスキャン順に並べたとき、評価点として採用される画素（斜線）の間隔ｓは平均でＮ_ａピクセルになるようにピックアップされる。ｓは平均Ｎ_ａ、幅Ｎ_ａの一様乱数で、ｒが０から１の間を取る一様乱数であるとき、
【数１９】

と書ける。
【００６３】
ここで、ステージｉ（上側ステージ）とステージｉ＋１（下側ステージ）の間での位置合わせを行う場合を考える。目的関数は、上下のステージの画像の画素値の２乗誤差とペナルティー関数からなる。まず、２乗誤差の計算方法について説明する。２乗誤差を計算する際、下側ステージの画像の画素値には平行移動や変形の他に
【数２０】

という画素値変換がなされる。上側ステージには画素値変換が行われない。ステージ数Ｎ_ｓ、ステージｉの画像内の評価点の数Ｃ_ｉ、ステージｉのｊ番目の点のある評価点の画像座標系での座標値を
【数２１】

、上側ステージの画素値
【数２２】

とコントラスト変換が施された下側ステージの画素値
【数２３】

の２乗誤差Ｓ_ｉは
【数２４】

と表される。
【００６４】
画像読み込み時に各ステージの画像の画素値の平均値が０になるように読み込んだ画素値から平均値を引き算してない場合は、２乗誤差Ｓ_ｉを計算する際に平均値を引き算する。
【００６５】
【数２５】

【００６６】
ここで、
【数２６】

はステージｉの画像の画素値の平均値である。目的関数は２乗誤差Ｓ_ｉにペナルティー関数Ｑ_ｉを加えたものを、ステージｉの着目画素（ピックアップされた評価点）の数Ｃ_ｉで割ったものである。
【００６７】
【数２７】

ε_ｉの後ろの括弧は目的関数ε_ｉが６つの変数を持ち、これらに関して最適化されることを明示的に表している。最適化する変数は、上側ステージｉの変換に関するものが、
【数２８】

の１個、下側ステージｉ＋１に関するものが、
【数２９】

の５個である。位置合わせの過程では、これら６つのパラメータを調節して目的関数が最小になるようにする。
【００６８】
【数３０】

はすべて０に固定される。ペナルティー関数は以下の式で表される。
【００６９】
【数３１】

【００７０】
目的関数にペナルティー項が加わっていることは、位置合わせが単に画素値の２乗誤差を小さくするように行われるのではないことを意味する。ペナルティー項はパラメータ
【数３２】

の絶対値がそれぞれの基準値
【数３３】

に比べて大きな値になりすぎないように値を制限する働きをもつ。また、
【数３４】

との差が各々の基準値σ_ｐ、σ_αに比べて大きくなりすぎないように値を制限する働きをもつ。
【００７１】
ペナルティー関数の最後の
【数３５】

の差に関する項は、上下の変形率の差があまり大きくならないように制限する働きがある。この項の存在により、パノラミックビュー画像の作成時に、変形を取り除いて貼り合わせた場合、算出される接合位置が上下の画像の重なり領域の中央から大きく外れることが防止される。
【００７２】
測定データにあるパラメータを決定するだけの情報がない場合、そのパラメータの値が変わっても２乗誤差の変化は小さく、このとき目的関数を最小にするパラメータ値はペナルティー関数が小さな値になるように選ばれる。選ばれる値は狭い範囲に分布することになり、調整範囲が狭いことになる。逆にパラメータを決定する十分な情報がある場合には広い調整範囲の中から適切なパラメータを選ぶことになる。このように、ペナルティー関数の導入によりどのパラメータをどの範囲から選ぶかが与えられた画像から自動的に判断されるようになるという効果がある。たとえば、与えられた２つのステージの間の血管が主に縦方向に走行している場合には縦位置を合わせる情報が画像中にあまり含まれないため、この画像の場合には縦位置の調整範囲が自動的に小さくされることになり、位置合わせ結果がより安定する。
【００７３】
次に、画像補間について説明する。図１８には２×２カーネルの場合の２次元補間を示している。目的関数の計算のためには任意の座標ｘでの画素値を求めることが必要である。一方、元画像はｍ×ｎ個の画素から構成された画像Ｉ（ｊ，ｉ）(i＝0,...,m-1, j＝0,...,n-1)として与えられているので、実数ｘ，ｙに対する画素値Ｉ（ｘ、ｙ）は何らかの補間により求めることが必要である。ここでは関数値と導関数がともに連続であり、かつカーネルサイズが２×２と小さい補間方法を選択した（２×２キュービック補間）。一般に、カーネルサイズが大きいほど補間精度が良くなるので、必要であればよりカーネルサイズが大きい補間方法を選択することも考慮すべきである。
【００７４】
ここで選択した２×２キュービック補間の２次元のカーネルｈ´（ｘ、ｙ）は、ｘおよびｙ方向のみに依存する１次元のカーネルｈ（ｘ）とｈ（ｙ）との積で表される。
【００７５】
【数３６】

【００７６】
このように自動位置合わせ処理は目的関数ε_ｉを最小にする変数
【数３７】

の組を求めることであると定式化された。これを実際に行うアルゴリズムが最適化である。特に今回の目的関数ε_ｉのように非線形の関数の最適化は非線形最適化と呼ばれる。非線形最適化には多くのアルゴリズムが知られている。たとえば、共役法(conjugate direction 法)、準ニュートン法、シンプレックス法やローゼンブロック法(Rosenbrock法)、シミュレーテッドアニーリング法、遺伝的アルゴリズム、タブーサーチ法(tabu search法)などを用いることができる。
以上説明した自動位置合わせ処理のデータフローを図１９に示している。
【００７７】
自動位置合わせ処理が終了した後に、その処理で決定したパラメータを用いて実際にパノラミックビュー画像合成処理が実行される。
パノラミックビュー画像合成処理は、Ｎ_ｓ個のステージの画像を位置合わせした後、ひとつの大きな画像として合成する処理である。各ステージの画像の位置関係は、自動位置合わせまたはマニュアル位置合わせにより得られた位置合わせパラメータにより与えられる。ただし、コントラスト変換のパラメータα_ｉ＋１、β_ｉ＋１は、貼り合わせの際には使用しない。従って、貼り合わせた結果の画像のコントラストを、入力画像と同じコントラストに維持する。
【００７８】
隣接する画像どうしは重なり合う位置関係にあるが、本実施形態で採用した合成方法では、上下の画像とも接合位置で切り取って接合することにより、合成した後の画像には重なりが無い。
【００７９】
図２０にパノラミックビュー画像の合成方法を模式的に示す。図２１にはパノラミックビュー画像合成の手順を示している。順番に説明する。
（１）まず、各ステージの変形後画像の左上隅座標と右下隅座標をグローバル座標系にて求める。
（２）次に、合成画像の左上隅座標と右下隅座標をグローバル座標系にて求める。
（３）そして、合成画像領域を確保し画素値を０で初期化する。
【００８０】
（４）続いて、各ステージの領域毎に、リサンプリングを行う。各ステージの左上隅座標と右下隅座標から、そのステージの画像が表示される画素の四角形領域をグローバル座標系にて求める（図２０参照）。２つのステージの間の重なり範囲の上半分は上側のステージ、下半分は下側のステージの領域に割り振られる。四角形領域を設定したら、この領域の各画素の画素値を該当ステージの画像からリサンプリングにより求める。
【００８１】
（５）画素値を入力画素値の値に戻す。つまり、リサンプリング後、各画素値に対して画像読み込み時に引き算した画素値を再び加えて、元の画素値に戻す。画像読み込み時に平均値を引き算しない場合は、この処理は不要である。画像データ有効領域の外を黒く表示させたい場合は、画像データ有効領域の外の画素値を０またはこれに近い特定の画素値に置き換える処理を行う。
【００８２】
図２２には画像の変形を行わないパノラミックビュー画像の合成方法について示している。画像の貼りあわせ時に、上下の画像の変形率ｐ_ｕ、ｐ_ｌを０とおくことにより、貼り合わせた画像に画像の変形が施されないように出来る。このとき、上下の画像の接合位置を次のように定めると、継ぎ合わせた部分でのずれが生じないようになる。上側画像の中心を（ｘ_１，ｙ₁）、下側画像の中心を（ｘ_２，ｙ_２）とし、接合位置と上側画像の中心との距離をｖとする。上下の画像の間隔は、ｈ＝ｙ_２−ｙ_１となる。
【００８３】
【数３８】

【００８４】
ただし、このように接合位置を決定すると、上下の変形率の値によっては接合位置が上下の画像からはみ出ることがある。このようなことが起こらないように、実際に継ぎ合わせるときは接合位置と上下の画像の端部とを比較し、継ぎ目が画像の外にならないように調整する必要がある。
【００８５】
画像の変形を元に戻すことにより下側画像は横位置を修正することが必要になる。修正後の上下の画像の横ずれ量（ｘ´２−ｘ１）は元の画像の横ずれ（ｘ２−ｘ１）に接合部での拡大率ｇ_ｖを乗じたものになる。
【００８６】
【数３９】

【００８７】
次に、マニュアルで接合位置を決定するパノラミックビュー画像の合成方法について説明する。上述では接合位置ｖを自動決定したが、これをマニュアルで設定する方法も考えられる。まず、隣接する２つのステージの画像を表示する。この際、上下の画像は既定の接合位置でつき合わせられる。既定の接合位置はたとえば上下の画像の切り取り幅が同じになる位置（ｈ／２）や前述の方法で決定した値が用いられる。操作者は表示された画像の接合状況を観察し必要なら接合位置と上下の画像の平行移動量を修正する。修正結果は直ちに画像として再表示される。このようにして最終的に操作者が選んだ接合位置と平行移動量がパノラミックビュー画像の合成に用いられる。
【００８８】
図２３は表示される画像の例である。図２３（ａ）は上下の画像を接合位置にて切断し貼りあわせた表示方法である。図２３（ｂ）は接合位置より外側の部分を反転した後に上下の画像を加算したものである。図２３（ｃ）は片方の画像（この場合上側画像）のみ反転して上下の画像を加算したものである。図２３（ａ）によれば貼りあわせた最終結果がどうなるかを判断することができる。図２３（ｂ）または図２３（ｃ）の表示方法によれば、接合位置や平行移動量をどの方向にどの程度修正すればよいかを容易に判断できるようになる。接合位置をマウスのドラッグ操作で修正する場合、通常は図２３（ａ）の表示を行い、ドラッグ中のみ図２３（ｂ）または図２３（ｃ）の表示を行うことにより、最終結果の判断と修正方向の判断の両方に有用な画像表示を行うことができる。
【００８９】
また、接合位置を決定するために半自動的な方法をとることもできる。上下の画像を貼りあわせた画像が表示された状態で、操作者が上側画像を選択すると、接合位置が暫定的にｖ＝ｈとされて、上側画像の全領域が表示されるようになる。この状態で、操作者はスムースにつなげたい血管をトレースする。２本の血管がある場合には２本の血管をトレースする。トレースする範囲はおおよそ上下両方の画像に含まれると思われる範囲で良い。次に、操作者は下側画像を選択すると接合位置が暫定的にｖ＝０とされて、下側画像の全領域が表示される。同様に操作者は血管をトレースする。上下画像の血管のトレースが終了すると、上下の血管が交差する交差点をすべて算出し、これら交差点の縦位置がなるべく等しくなるような横方向の平行移動量を決定する。このときの交差点の縦位置を接合位置として、貼り合せ画像が再描画される。必要なら操作者は縦方向の平行移動量を修正することができる。縦方向の修正を行うと、上記の処理を再実行し再び横方向の平行移動量を算出しなおし、接合位置を際算出し、貼り合せ画像が描画される。このようにして最終的に決定された平行移動量（縦位置と横位置）と接合位置が、パノラミックビューの画像の合成に用いられる。
【００９０】
画像読み込み時に、各画像の画像データ有効領域内の平均値を求め、各画像の画素値から平均値を引き算する処理を行う。算出した各画像の平均値は一旦記憶され、パノラミックビュー画像合成時に元の画素値に戻すために使われる。マニュアル位置合わせ時に上下の画像を重畳表示しない場合は、平均値を引く処理は不要である。
【００９１】
以上説明したように、本実施形態によれば、各ステージの位置合わせが自動的に行われるため、位置合わせに要する手数を軽減できるという効果がある。自動位置合わせにおいては画像の変形を考慮して画像がスムーズに接続する変形量を求めているが、画像貼り合わせ処理においては、変形せずに画像がスムーズに接続する最適な接合位置を求めてその位置で画像をつなぎ合わせるので、最終的に貼りあわせた画像は変形がなく、かつ、接合位置にて画像がスムーズに接合される。従って、本実施形態によれば画像がスムーズに接続されているにもかかわらず、オリジナル画像の長さや角度が保存されていて診断にも適する画像を提供できるという効果がある。
【００９２】
（変形例）
本発明は、上述した実施形態に限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で種々変形して実施することが可能である。さらに、上記実施形態には種々の段階が含まれており、開示される複数の構成要件における適宜な組み合わせにより種々の発明が抽出され得る。例えば、実施形態に示される全構成要件から幾つかの構成要件が削除されてもよい。
【００９３】
【発明の効果】
本発明によれば、隣り合うステージどうしの画像の貼り合わせ位置を最適化することができる。それにより、貼り合わせ位置調整の手間が軽減され、しかも画像に変形を加えることなく血管が自然につながるように画像を貼り合わせることができる。
【図面の簡単な説明】
【図１】本発明の実施形態に係るＸ線診断装置の構成を示す図。
【図２】本実施形態において、ステッピングＤＳＡの撮影時のフレームの移動を示す図。
【図３】本実施形態において、５つのステージで撮影したＤＳＡ画像と、それらを貼り合せたパノラミックビュー画像の例を示す中間調画像。
【図４】本実施形態において、自動位置合わせとマニュアル位置合わせを含むパノラミックビュー画像作成の全体の処理フローを示す図。
【図５】本実施形態において、平行移動に関するパラメータを示す図。
【図６】本実施形態において、ずれ変形に関するパラメータを示す図。
【図７】本実施形態において、拡大・縮小に関するパラメータを示す図。
【図８】本実施形態において、周辺部変形に関するパラメータを示す図。
【図９】本実施形態において、初期的な画像座標系として、撮影したｍ×ｎ画素の画像内の画素の位置を表す座標系を示す図。
【図１０】本実施形態において、変形後画像座標系として、画像に周辺部変形を施した後の座標系を示す図。
【図１１】本実施形態において、グローバル座標系に各ステージの画像を配置した様子を示す図。
【図１２】本実施形態において、画像有効データ領域とＲＯＩ設定有効領域を示す図。
【図１３】本実施形態において、しきい値処理による自動ＲＯＩ設定について示す図。
【図１４】本実施形態において、自動ＲＯＩ設定の処理手順を示す図。
【図１５】本実施形態において、フィルタ処理について示す図。
【図１６】本実施形態において、テストプログラムでフィルタリングを行う領域として、抽出されたＲＯＩ領域と、フィルタリング有効領域と上側画像のフィルタリング有効領域との積領域とを示す図。
【図１７】本実施形態において、ＲＯＩ内の画素から評価点をピックアップする処理を示す図。
【図１８】本実施形態において、２×２カーネルの場合の２次元補間を示す図。
【図１９】本実施形態において、自動位置合わせ処理のデータフローを示す図。
【図２０】本実施形態において、パノラミックビュー画像の合成処理を示す図。
【図２１】本実施形態において、パノラミックビュー画像合成の手順を示す図。
【図２２】本実施形態において、画像の変形を行わないパノラミックビュー画像の合成処理を示す図。
【図２３】本実施形態において、様々な貼り合わせ方法それぞれの表示画像の例を示す図。
【符号の説明】
１１…Ｘ線管、１２…Ｘ線絞り、１３…天板、１４…検出器、１５…フレーム、１６…アナログディジタル変換器、１７…記憶装置、１８…ディジタル画像処理装置、１９…表示モニタ、２０…操作部、２１…撮影系幾何学的情報発生装置。[0001]
BACKGROUND OF THE INVENTION
The present invention relates to a bonding apparatus for a plurality of medical images having different imaging positions.
[0002]
[Prior art]
When the outer area of the lower extremity, such as an artery of the lower limb, exceeds the size of the detector, contrast images are taken at a plurality of table positions (called stages) while moving the table to follow the movement of the contrast agent. Then, stepping DSA (digital subtraction angiography) for subtracting the pre-contrast mask images corresponding to the stages from these contrast images is performed. A combination of a plurality of DSA images thus captured so that the entire imaging region can be seen is called a long display. However, because the medical images taken are distorted, and because the blood vessels meander and run, the geometry is not always constant. It becomes a problem to become.
[0003]
In order to connect the blood vessels smoothly at the joint, a method of manually pasting the positions of the images on the adjacent stages and then pasting them together can be considered. However, in reality, due to the distortion of the X-ray detector and the difference in geometric magnification between the images, there is a problem that the upper and lower images are generally not connected cleanly simply by translation. If the image is transformed, this problem can be solved and the images can be smoothly connected. However, if the image is transformed, the pasted image is deformed. There is a problem that the length is different and it is not suitable for diagnosis.
[0004]
[Problems to be solved by the invention]
The object of the present invention is to reduce the labor for adjusting the bonding position by optimizing the bonding position of the images between the adjacent stages, and also to connect the images so that the blood vessels are naturally connected without deforming the images. There is to stick together.
[0005]
[Means for Solving the Problems]
  One aspect of the present invention provides a position shift of an image between a joint position between a plurality of medical images partially overlapped on a subject.amount, Image deformationamount, Image contrast conversion ratioBetween the medical images while adjusting the rateOverlapping partPart ofMeans for determining the difference between the pixel values to be minimum, and means for combining the plurality of medical images based on the determined joining position by enlarging and reducing the original image without deformation. A medical image laminating device is provided.
[0006]
DETAILED DESCRIPTION OF THE INVENTION
Hereinafter, an X-ray diagnostic apparatus according to the present invention will be described with reference to the drawings according to preferred embodiments.
FIG. 1 is a diagram showing a configuration of an X-ray diagnostic apparatus according to an embodiment of the present invention. An X-ray tube 11 is attached to one end of the U-shaped frame 15. An X-ray diaphragm 12 for limiting the X-ray irradiation field is attached to the X-ray irradiation window of the X-ray tube 11. The other end of the vertically moving frame 15 is attached with a combination of an image intensifier and a TV camera, or an X-ray detector 14 which is a flat detector that directly converts X-rays into electric signals. The X-ray detector 14 is provided so as to face the X-ray tube 11, and the subject P placed on the couch top 13 is arranged between them. When photographing the stepping DSA, as shown in FIG. 2, the vertically moving frame 15 is withdrawn from the ceiling rail so that it can move over a wide range from the top of the subject to the tip of the foot along the body axis of the subject P. Alternatively, it is provided on the floor rail, or provided on the bed frame rail. Thereby, without moving the subject P, the imaging is repeated while tracking the contrast medium injected into the subject, and a plurality of X-rays continuously and partially overlapping along the body axis of the subject. Images can be acquired.
[0007]
An analog X-ray image generated by the X-ray detector 14 is converted into digital X-ray image data by an analog-digital converter 16 and sent to a digital image processing device 18. The digital image processing apparatus 18 performs a pasting process on a plurality of stages of X-ray image data. The X-ray image data of a plurality of stages and the combined image (panoramic view image) are displayed on the display monitor 19 as a grayscale image. In addition, the digital X-ray image data is stored in a storage device 17 such as a magneto-optical disk device.
[0008]
The digital image processing apparatus 18 is an imaging system comprising a frame 15, an X-ray tube 11, an X-ray diaphragm 12, an X-ray detector 14, a couch top plate 13, and the like in the process of laminating a plurality of stages of X-ray images. Geometric information is supplied from the geometric information generator 21 and used together with information input by the operator via the operation unit 20. The geometric information includes a geometric enlargement ratio determined by a distance from the X-ray tube 11 to the top plate 13 and a distance from the X-ray tube 11 to the instrument 14, a position of the top plate 13, and the like. The geometric information generating device 21 includes, for example, a plurality of rotary encoders installed at appropriate positions required for acquiring geometric information.
[0009]
Here, the outline | summary is demonstrated about the bonding method of the X-ray image which concerns on this embodiment.
[0010]
"Outline of automatic bonding method for stepping DSA images"
FIG. 3 is an example of DSA images photographed at five stages (positions of the top plate 13) and a panoramic view image obtained by pasting them using this method. Typically, an arteriographic image with several regions of the lower limb (distal from the pelvis) as an imaging region and a subtraction image created from a non-contrast image (stepping DSA image) are pasted. In addition, a method for creating a long image of the lower limb artery will be described. Of course, the arterial contrast image is taken at the same position as the non-contrast image, and each photographing position is set so that adjacent photographing regions partially overlap.
[0011]
FIG. 4 shows an overall processing flow for creating a panoramic view image including automatic alignment and manual alignment. As an overall processing procedure, a plurality of images are roughly aligned manually (manually), and then automatic alignment is performed on a plurality of manually aligned images. Region of interest setting processing (automatic ROI setting) and filter processing are automatically executed as preprocessing, and the alignment amount and joining position are automatically determined for a plurality of preprocessed images. (Automatic alignment) Finally, a panoramic view image is created by combining a plurality of images in accordance with the automatically determined detailed positioning amount and the joining position (panoramic view image synthesis).
[0012]
The input of the entire processing is the geometric information (such as the top plate movement amount and the geometric magnification) acquired by the geometric information generating device 21 at the time of stepping DSA imaging, and N_sThis is a subtraction image (hereinafter simply referred to as an image) of each stage. The output is a panoramic view image obtained by pasting the subtraction images of a plurality of stages. FIG. 4 shows a case where the number of stages is three.
[0013]
In the above processing flow, a panoramic view image creation process is performed by the following operation. First, a plurality of input images of a plurality of stages are displayed on the screen in a state where they are pasted together at a default pasting position. The default bonding position is determined based on geometric information at the time of photographing the stepping DSA image. The image of the two adjacent stages is displayed on the screen, and fine adjustment can be performed by translating the position of the image of the lower stage by manual operation.
[0014]
When automatic alignment is executed in this state, alignment between the two stages displayed on the screen is automatically performed. The alignment result is displayed again on the screen, and if the alignment result is not satisfactory, manual alignment can be performed again. When the alignment of the two displayed stages is completed, another two adjacent stages are displayed on the screen, and manual alignment and automatic alignment are performed in the same manner. When the alignment between all the stages is completed, an image (a panoramic view image) in which all the stages are bonded together is created.
[0015]
In manual positioning, the amount of translation (vertical and horizontal directions) between the two displayed stages can be adjusted. The “vertical translation amount in the vertical direction (vertical translation amount)” and “parallel translation amount in the horizontal direction (horizontal translation amount)” are referred to as manual alignment parameters. Furthermore, the joining position of the images on the upper and lower stages can be adjusted.
[0016]
In the automatic alignment, the vertical translation amount and the horizontal translation amount are automatically adjusted with respect to the vertical translation amount and the horizontal translation amount of the manual, and the deformation amount of the image is adjusted. The deformation of the image is represented by two parameters of the “upper half deformation rate” and the “lower half deformation rate” of the image. Between the two adjacent stages, the “lower half deformation rate” of the upper stage image (for example, the head side) and the “upper half deformation rate” of the lower (for example, lower limb) stage image are adjusted. The
[0017]
In addition, only the lower image is subjected to contrast conversion at the time of alignment. Contrast conversion is represented by two parameters of “pixel value conversion ratio” and “offset”.
[0018]
"Vertical translation amount" and "Horizontal translation amount" related to translation, "Upper half deformation rate" and "Lower half deformation rate" related to image deformation, "Pixel value conversion ratio" and "Offset" related to contrast conversion, These 6 parameters are adjusted to determine the position of pasting between images. However, the two parameters for contrast conversion are ignored when the aligned images are finally pasted together. Therefore, the panoramic view image pasted with the original image has the same contrast. In the automatic alignment process, six parameters of the vertical translation amount, the horizontal translation amount, the upper half deformation rate, the lower half deformation rate, the contrast conversion ratio, and the offset of the image are used as control factors. Then, the optimum joining position is automatically determined using the adjusted values of the upper half deformation rate and the lower half deformation rate.
[0019]
When these six parameters are adjusted, the result is displayed each time, and there are a method of deforming the image and a method of not deforming the image. When the image is not deformed, the result of pasting the images using the determined joining position is presented.
[0020]
In the image pasting process, an image is created by pasting all the stages using the parameters adjusted in the automatic positioning process. In that case, a method is used in which an image without deformation is pasted using the determined joining position. When specified by the operator, it is possible to perform image pasting including deformation.
[0021]
The automatic alignment process compares pixel values between the position of the upper image and the lower image corresponding to it, and optimizes the automatic alignment parameters so that the difference between the pixel values is small. To go. The pixel value comparison is actually performed only around the blood vessel. First, an evaluation point is picked up from each point of the region of interest ROI set around the blood vessel in the upper image. Next, the position of the lower image corresponding to each evaluation point is obtained from the alignment parameter. Then, the pixel values of both images at the corresponding points are compared. In the “automatic ROI setting” process, a blood vessel in the upper image is automatically extracted, and a process of setting an ROI around it is executed.
[0022]
Another pre-processing for automatic alignment is “filter processing”. By applying a low-pass filter in advance to an image to be compared with pixel values, the alignment process can be operated more stably.
[0023]
The main point of the present embodiment is to automatically perform image alignment in the “automatic alignment process” and the “image pasting process” in FIG. 4. In addition, an optimum joining position is automatically determined at the time of image bonding, and bonding is performed without deformation.
[0024]
Hereinafter, a method for pasting X-ray images according to the present embodiment will be described in detail. First, the alignment and deformation parameters will be described.
FIG. 5 shows parameters related to the parallel movement. In order to move the point x to x ′, two types of parameters t and m are used. t indicates an automatic parallel movement amount, and m indicates a manual parallel movement amount. t and m are respectively composed of a vertical translation amount and a horizontal translation amount.
[0025]
[Expression 1]

[0026]
The deformation of the image is determined in combination by displacement deformation, enlargement / reduction, and peripheral deformation peculiar to X-ray imaging. FIG. 6 shows parameters regarding displacement deformation. Horizontal displacement is h_xThe vertical displacement is h_yGiven in.
[0027]
[Expression 2]

[0028]
FIG. 7 shows parameters relating to enlargement / reduction. S scaling in the horizontal direction is S_xThe vertical scaling deformation is S_yGiven in.
[0029]
[Equation 3]

[0030]
FIG. 8 shows parameters related to peripheral deformation. Independent deformation amount p for the upper half and the lower half of the image_u, P_lIs set.
[0031]
[Expression 4]

Inverse operator P⁻ _iIs expressed by the following equation.
[0032]
[Equation 5]

[0033]
In this peripheral deformation method, the deformation amount p independent of the upper half and the lower half of the image_u, P_lIs set. As a result, even if the distortion of the image is complicated and the distortion is different between the upper and lower sides of the image, there is an effect that the alignment can be performed satisfactorily. Also, since this method only deforms the image in the horizontal direction, the amount of deformation and the amount of translation in the vertical direction of the two images become independent at the time of optimization. Since the amount of deformation is also independent, there is an effect that optimization becomes easy and positioning can be performed at high speed.
[0034]
Next, a coordinate system at each processing stage used in the following description will be described.
FIG. 9 shows an initial image coordinate system (two-dimensional) that represents a position of a pixel in a photographed m × n pixel image, and actually handles negative values and non-integer values. be able to. FIG. 10 shows a coordinate system after the peripheral deformation is applied to the image as the post-deformation image coordinate system. Of course, if the peripheral deformation is not performed, the image coordinate system is the same.
[0035]
When aligning images of a plurality of stages, a virtual overall coordinate system is used for arranging these images while changing their positions. In the present embodiment, this is called a global coordinate system. FIG. 11 shows a state in which images of each stage are arranged in the global coordinate system. The area of the panoramic view image obtained by combining (pasting) the arranged images is also set in the global coordinate system.
[0036]
Of course, coordinate conversion processing is required between these three coordinate systems. Here, conversion between each two-dimensional coordinate system is expressed by affine transformation. Therefore, the vector representing the position is treated as a three-dimensional vector in which elements having a value of 1 are always added in addition to the x and y coordinates.
[0037]
Coordinate in global coordinate system y₀, I-th (i = 1, 2,... N_s) Stage coordinate after deformation of the stage y_iThe relationship between each coordinate system is
[Formula 6]

It is represented by E_iIs a coordinate transformation matrix from the stage i-1 post-deformation coordinate system to the stage i post-deformation coordinate system._i, Automatic alignment translation T_i, Displacement deformation H_i, Enlargement / reduction S_i, Translation from the upper left corner of the image to the center of the image C_iUsing,
[Expression 7]

It is represented by The coordinate conversion from the image coordinate system after deformation to the image coordinate system before deformation is
[Equation 8]

It is represented by x_iIs a coordinate value in the image coordinate system before deformation. Each element of coordinate transformation is as follows.
[0038]
[Equation 9]

L is the width of the square image. P_iIs a nonlinear operator,
[Expression 10]

Defined by From these, the coordinate transformation matrix E from the stage i-1 post-deformation coordinate system to the stage i post-deformation coordinate system_iIs
## EQU11 ##

It becomes. Coordinate transformation between adjacent stages is transformed into transformation operator D_iWhen expressed using
[Expression 12]

It becomes. x_iIs the coordinate of the image coordinate system of the i-th stage before transformation, x_{i + 1}Is the coordinates of the image coordinate system before the deformation of the i + 1th stage. Transformation operator D_iIs expressed as follows.
[0039]
[Formula 13]

[0040]
Next, actual automatic alignment processing will be described. In the automatic alignment processing, “automatic ROI setting” and “filter processing” are executed as pre-processing.
[0041]
First, automatic ROI setting will be described.
FIG. 12 shows an image effective data area and an ROI setting effective area. This is an area where valid data of a photographed X-ray image exists in an area of rectangular image data of m × n pixels. Here, the image valid data area is defined as an image data valid area in a circle having a radius Rval centered at the center point (m / 2, n / 2) in the image coordinate system.
[0042]
In the test program, before all processing such as automatic alignment, automatic ROI setting, and filtering (immediately after image reading), the images of all stages are averaged for the pixel values in the image data effective area for each image. It is subtracted from each pixel and normalized so that the average value becomes 0, and pixels outside the image data effective area are replaced with the pixel value 0. This is because in the test program, images between stages are superimposed and displayed when manual alignment is performed. In order for the superimposed display to occur naturally, the upper stage is displayed as it is, the lower stage performs parallel movement and deformation of the image, and then displays the pixel values with the positive and negative reversed, and the images of the upper and lower stages. In the overlapping area, the inverted pixel value of the lower image is added to the upper image. At this time, if the background region other than the blood vessel has a value close to 0, the superposition is naturally performed. When the images of the upper and lower stages are not displayed in an overlapping manner, the process of subtracting the average value of the pixel values when reading the image is unnecessary. In this case, it is better to replace pixels outside the image data effective area with average values.
[0043]
The ROI setting effective area is a circular area within a radius Rroi that is somewhat smaller than the image data effective area. This area is used in the automatic ROI setting process. The test program assumes that the matrix size of the input image is 512 × 512, and uses a value of 232 (pixel) for Rval and 215 (pixel) for Rroi. This value needs to be changed for each image size.
[0044]
FIG. 13 shows automatic ROI setting by threshold processing. In the automatic alignment algorithm, the pixel values of corresponding points in the image between adjacent stages are compared, and the positional relationship between the upper and lower images is obtained so that the difference in pixel value is reduced. If the evaluation points for calculating the difference in pixel values are not limited only to the vicinity of the blood vessel, the evaluation points are affected by the contrast fluctuation of the background portion other than the blood vessel, and thus cannot be stably positioned. Also, it is desirable to limit the evaluation points only to the periphery of the blood vessel region for speeding up. In the automatic ROI setting process, a process of setting a region of interest (ROI) for alignment is performed as pre-processing of the alignment process.
[0045]
In the alignment process, a point in the lower image corresponding to each point in the ROI set in the upper image of the two stages is obtained, and the difference between the two pixel values is reduced. As pre-processing, ROI is set around the blood vessel in the upper image in the automatic ROI setting.
[0046]
The input is the upper image and initial alignment parameters, and the output is the set ROI shape. The ROI to be set is set only to an area that falls within the ROI effective area in both the upper and lower images when the images of the upper and lower stages are superimposed using the parallel movement amount of the initial alignment parameter. The ROI effective area is a circular area set to be smaller than the circle in which the effective data of the X-ray image exists.
[0047]
FIG. 14 shows a processing procedure for automatic ROI setting. Automatic ROI setting between the two stages is performed in the following procedure.
1. Average value a and variance v of pixel values in the image data effective area of the upper stage image²Ask for.
2. In the upper stage image, the pixel value is (an_tv) The following areas are extracted by threshold processing (area A).
n_tIs a constant that determines the threshold, and 2 is used in the test program. Accordingly, the threshold value is 2SD (twice the standard deviation) smaller than the average value.
[0048]
3. The ROI setting effective area is determined by performing an AND operation on the ROI setting effective area of the upper stage and the ROI setting effective area of the lower stage (area E).
4). A logical product operation between the area A and the ROI setting effective area is performed (area B).
5. Region B is n_dEnlarge by pixel (area C).
n_dIs a constant that determines the size of enlargement, and the test program uses 5 pixels.
5. A logical product operation is performed between the area C and the ROI setting effective area of the upper stage (area D).
[0049]
In the test program, the area D thus extracted is run-length encoded and stored in the memory. In the automatic alignment process, the process is performed with reference to the ROI data corresponding to the run-length encoded region D. This eliminates the need to scan all the pixels at the time of automatic alignment, so that the automatic alignment process can be performed at high speed.
[0050]
Next, filter processing as another preprocessing will be described.
FIG. 15 shows filter processing. In the automatic alignment process, corresponding pixel values in the ROI are compared with each other so that the difference between the pixel values is reduced. Therefore, if there is no overlapping of blood vessels in the upper and lower stage images in the initial state before alignment, the image will not converge to the correctly aligned state. In order to reduce the probability of such a result, it is effective to perform filter processing (low-pass filter processing) to blur the original image as pre-processing for automatic alignment.
[0051]
The filter kernel is n_f× n_fHowever, in order to reduce the number of operations, K can be broken down into vertical and horizontal convolutions._ij= K_iK_jThe kernel represented by is adopted. K_iIs a decomposed one-dimensional filter coefficient, and the following triangular coefficient is used.
[0052]
[Expression 14]

[0053]
Using this kernel, the low-pass filtering can be decomposed into two stages of one-dimensional convolution.
[0054]
(1) Horizontal convolution
[Expression 15]

(2) Vertical convolution
[Expression 16]

Are pixel values of (i, j) pixels of the original image, the horizontal convolution result image, and the filter processing result image, respectively.
[0055]
The filtering time can be significantly reduced by limiting the filtering process only to the portion that may be referred to in the alignment, instead of performing the filtering process on the front surface of the images of all the stages. In the test program, the area to be filtered is limited to the following two parts as shown in FIG.
[0056]
Extracted ROI region (region D)
-Product area of filtering effective area and filtering effective area of upper image
The filtering effective area is a circular area that is slightly larger than the image data effective area, and the radius of the test program is 10 pixels larger than the radius of the image data effective area (242 pixels).
[0057]
After completion of the above preprocessing, “automatic alignment processing” is executed. Total number of stages N_sAre aligned every two adjacent stages. The alignment and deformation applied to the images of the two stages (numbers i and i + 1) are defined by the following four parameters.
[0058]
[Expression 17]

[0059]
During automatic alignment, these six parameters (p^l, T_x, T_y, P^u, Α, β) are added (positioning, deformation, contrast conversion), and then the pixel values are compared. Then, alignment processing is performed by obtaining a parameter that reduces the difference in pixel values.
[0060]
Here, translation from the upper left corner of the image to the center of the image c_ix, C_iyIs related to the deformation of the upper stage, but this value is related not only to the alignment between stages i and i + 1, but also to the alignment between i-1 and i, and these values are considered as alignment parameters. If the two stages are aligned, the alignment result between the other stages is destroyed. The alignment method of this embodiment employs a method of aligning every two adjacent stages so as not to cause such a situation, and c_ixAnd c_iyIs set to a constant value of 0 so that it is not changed in the alignment process. That is,
[Expression 18]

[0061]
As described above, in the automatic alignment process, instead of comparing all the pixels in the set ROI with the corresponding pixel values in the next stage, for pixels that are randomly picked up from all the pixels in the ROI, Compare with the corresponding pixel value in the next stage. As a result, the speed of the automatic alignment is increased. FIG. 17 shows a process of picking up evaluation points from the pixels in the ROI.
[0062]
Each square in FIG. 17 indicates each pixel of the image of stage i, and a frame indicated by a thick line indicates an area in the ROI. Each pixel of the ROI is scanned in the horizontal direction of the image. When these pixels are arranged in the scan order, the interval s between pixels (hatched lines) used as evaluation points is N on average._aPicked up to be a pixel. s is average N_a, Width N_aAnd when r is a uniform random number between 0 and 1,
[Equation 19]

Can be written.
[0063]
Here, consider a case where alignment is performed between stage i (upper stage) and stage i + 1 (lower stage). The objective function includes a square error of the pixel values of the upper and lower stage images and a penalty function. First, a method for calculating the square error will be described. When calculating the square error, in addition to translation and deformation, the pixel values of the lower stage image
[Expression 20]

Pixel value conversion is performed. Pixel value conversion is not performed on the upper stage. Number of stages N_s, The number of evaluation points C in the image of stage i_i, The coordinate value in the image coordinate system of the evaluation point with the j-th point of stage i
[Expression 21]

, Pixel value of upper stage
[Expression 22]

And lower stage pixel values that have undergone contrast conversion
[Expression 23]

Square error S_iIs
[Expression 24]

It is expressed.
[0064]
When the average value is not subtracted from the read pixel value so that the average value of the pixel values of the images of each stage becomes 0 at the time of image reading, the square error S_iWhen calculating, subtract the average value.
[0065]
[Expression 25]

[0066]
here,
[Equation 26]

Is the average value of the pixel values of the image of stage i. The objective function is the square error S_iPenalty function Q_iIs the number C of the pixel of interest (picked evaluation points) of stage i_iDivided by.
[0067]
[Expression 27]

ε_iIs the objective function ε_iExpresses that it has 6 variables and is optimized for these. The variables to optimize are related to the upper stage i transformation,
[Expression 28]

, One related to the lower stage i + 1
[Expression 29]

It is five. During the alignment process, these six parameters are adjusted to minimize the objective function.
[0068]
[30]

Are all fixed to zero. The penalty function is expressed by the following formula.
[0069]
[31]

[0070]
The addition of a penalty term to the objective function means that the alignment is not simply performed to reduce the square error of the pixel value. The penalty term is a parameter
[Expression 32]

The absolute value of each is the reference value
[Expression 33]

It has the function of limiting the value so that it does not become too large compared to. Also,
[Expression 34]

Is the difference between each reference value σ_p, Σ_αIt works to limit the value so that it is not too large compared to.
[0071]
At the end of the penalty function
[Expression 35]

The term relating to the difference between the two functions serves to limit the difference between the upper and lower deformation ratios from becoming too large. The existence of this term prevents the calculated joining position from greatly deviating from the center of the overlapping region of the upper and lower images when the panoramic view image is created with the deformation removed and pasted together.
[0072]
If there is not enough information to determine a parameter in the measurement data, the change in the square error is small even if the parameter value changes. At this time, the parameter value that minimizes the objective function has a small penalty function. Chosen. The selected value will be distributed in a narrow range, and the adjustment range will be narrow. On the other hand, when there is sufficient information for determining the parameter, an appropriate parameter is selected from a wide adjustment range. As described above, there is an effect that it is automatically determined from the given image which parameter is selected from which range by introducing the penalty function. For example, when the blood vessel between two given stages is traveling mainly in the vertical direction, the information for adjusting the vertical position is not included in the image, so in this image, the adjustment of the vertical position is performed. The range is automatically reduced, and the alignment result is more stable.
[0073]
Next, image interpolation will be described. FIG. 18 shows two-dimensional interpolation in the case of a 2 × 2 kernel. In order to calculate the objective function, it is necessary to obtain a pixel value at an arbitrary coordinate x. On the other hand, the original image is given as an image I (j, i) (i = 0,..., M−1, j = 0,..., N−1) composed of m × n pixels. Therefore, it is necessary to obtain the pixel value I (x, y) for the real numbers x and y by some interpolation. Here, an interpolation method in which both the function value and the derivative are continuous and the kernel size is as small as 2 × 2 is selected (2 × 2 cubic interpolation). In general, the larger the kernel size, the better the interpolation accuracy. If necessary, it should be considered to select an interpolation method having a larger kernel size.
[0074]
The two-dimensional kernel h ′ (x, y) of the 2 × 2 cubic interpolation selected here is expressed by the product of the one-dimensional kernel h (x) and h (y) that depends only on the x and y directions. The
[0075]
[Expression 36]

[0076]
In this way, the automatic alignment process performs the objective function ε_iThat minimizes
[Expression 37]

It was formulated to be a set of The algorithm that actually does this is optimization. In particular, this objective function ε_iThe optimization of a nonlinear function like this is called nonlinear optimization. Many algorithms are known for nonlinear optimization. For example, a conjugate method (conjugate direction method), a quasi-Newton method, a simplex method, a Rosenblock method (Rosenbrock method), a simulated annealing method, a genetic algorithm, a tabu search method (tabu search method), or the like can be used.
FIG. 19 shows a data flow of the automatic alignment processing described above.
[0077]
After the automatic alignment process is completed, the panoramic view image synthesis process is actually executed using the parameters determined by the process.
Panoramic view image composition processing is N_sThis is a process of aligning the images of individual stages and then combining them as one large image. The positional relationship between the images on each stage is given by the alignment parameters obtained by automatic alignment or manual alignment. However, contrast conversion parameter α_{i + 1}, Β_{i + 1}Is not used for pasting. Therefore, the contrast of the images resulting from the pasting is maintained at the same contrast as the input image.
[0078]
Adjacent images are in an overlapping position. However, in the combining method employed in the present embodiment, the upper and lower images are cut and bonded at the bonding position, so that the combined images do not overlap.
[0079]
FIG. 20 schematically shows a method for synthesizing a panoramic view image. FIG. 21 shows a procedure for synthesizing a panoramic view image. We will explain in order.
(1) First, the upper left corner coordinates and lower right corner coordinates of the deformed image of each stage are obtained in the global coordinate system.
(2) Next, the upper left corner coordinate and the lower right corner coordinate of the composite image are obtained in the global coordinate system.
(3) A composite image area is secured and the pixel value is initialized to zero.
[0080]
(4) Subsequently, resampling is performed for each stage area. From the upper left corner coordinates and lower right corner coordinates of each stage, a rectangular area of the pixel on which the image of the stage is displayed is obtained in the global coordinate system (see FIG. 20). The upper half of the overlapping range between the two stages is assigned to the upper stage area, and the lower half is assigned to the lower stage area. When the rectangular area is set, the pixel value of each pixel in this area is obtained by resampling from the image of the corresponding stage.
[0081]
(5) Return the pixel value to the value of the input pixel value. That is, after resampling, the pixel value subtracted at the time of image reading is added again to each pixel value to return to the original pixel value. If the average value is not subtracted at the time of image reading, this process is unnecessary. When it is desired to display the outside of the image data effective area in black, a process of replacing the pixel value outside the image data effective area with a specific pixel value close to 0 or this is performed.
[0082]
FIG. 22 shows a method for synthesizing a panoramic view image without image deformation. Deformation rate p of top and bottom images when pasting images_u, P_lBy setting 0 to 0, it is possible to prevent the image from being deformed on the bonded image. At this time, if the joining positions of the upper and lower images are determined as follows, a shift at the joined portion does not occur. The center of the upper image is (x₁, Y₁), The center of the lower image (x₂, Y₂), And v is the distance between the joint position and the center of the upper image. The interval between the upper and lower images is h = y₂-Y₁It becomes.
[0083]
[Formula 38]

[0084]
However, when the joining position is determined in this way, the joining position may protrude from the upper and lower images depending on the values of the upper and lower deformation rates. In order to prevent this from happening, it is necessary to compare the joining position with the edge portions of the upper and lower images and adjust the seam so that it does not come out of the image.
[0085]
It is necessary to correct the horizontal position of the lower image by restoring the deformation of the image. The lateral shift amount (x′2-x1) of the upper and lower images after correction is the enlargement ratio g at the joint portion to the lateral shift (x2-x1) of the original image_vMultiplied by.
[0086]
[39]

[0087]
Next, a method for synthesizing a panoramic view image for manually determining a joining position will be described. In the above description, the joining position v is automatically determined, but a method of manually setting this is also conceivable. First, images of two adjacent stages are displayed. At this time, the upper and lower images are brought together at a predetermined joining position. As the predetermined joining position, for example, a position (h / 2) where the upper and lower images are cut out at the same width or a value determined by the above-described method is used. The operator observes the joining state of the displayed images and corrects the joining position and the parallel movement amount of the upper and lower images if necessary. The correction result is immediately redisplayed as an image. Thus, the joint position and the amount of parallel movement finally selected by the operator are used for synthesizing the panoramic view image.
[0088]
FIG. 23 is an example of a displayed image. FIG. 23A shows a display method in which the upper and lower images are cut and bonded at the joining position. FIG. 23B is obtained by adding the upper and lower images after inverting the portion outside the joining position. FIG. 23C is obtained by inverting only one image (in this case, the upper image) and adding the upper and lower images. According to FIG. 23A, it can be determined what the final result of pasting will be. According to the display method of FIG. 23B or FIG. 23C, it is possible to easily determine in what direction and how much the joining position and the parallel movement amount should be corrected. When the joint position is corrected by dragging the mouse, the display of FIG. 23A is usually performed, and the display of FIG. 23B or FIG. Image display useful for both determination of the correction direction can be performed.
[0089]
A semi-automatic method can also be used to determine the joining position. When the operator selects the upper image in a state where the upper and lower images are displayed, the joint position is temporarily set to v = h, and the entire area of the upper image is displayed. In this state, the operator traces the blood vessel to be connected smoothly. If there are two blood vessels, trace the two blood vessels. The range to be traced may be a range that seems to be included in both the upper and lower images. Next, when the operator selects the lower image, the joint position is temporarily set to v = 0, and the entire area of the lower image is displayed. Similarly, the operator traces the blood vessel. When the tracing of the blood vessels in the upper and lower images is completed, all intersections where the upper and lower blood vessels intersect are calculated, and a lateral translation amount is determined so that the vertical positions of these intersections are as equal as possible. The combined image is redrawn with the vertical position of the intersection at this time as the joining position. If necessary, the operator can correct the amount of vertical translation. When the correction in the vertical direction is performed, the above processing is re-executed, the parallel movement amount in the horizontal direction is calculated again, the joint position is calculated, and the composite image is drawn. The translation amount (vertical position and horizontal position) and the joint position finally determined in this way are used for synthesizing the panoramic view image.
[0090]
At the time of image reading, an average value in the image data effective area of each image is obtained, and a process of subtracting the average value from the pixel value of each image is performed. The calculated average value of each image is temporarily stored and used to return to the original pixel value when the panoramic view image is synthesized. When the upper and lower images are not superimposed and displayed at the time of manual alignment, the process of subtracting the average value is unnecessary.
[0091]
As described above, according to the present embodiment, since the alignment of each stage is automatically performed, there is an effect that the number of steps required for the alignment can be reduced. In automatic alignment, the amount of deformation that connects images smoothly is determined in consideration of image deformation, but in the image pasting process, the optimal joint position that connects images smoothly without deformation is determined. Since the images are joined at that position, the finally joined images are not deformed, and the images are smoothly joined at the joining position. Therefore, according to the present embodiment, there is an effect that it is possible to provide an image suitable for diagnosis because the length and angle of the original image are preserved even though the images are smoothly connected.
[0092]
(Modification)
The present invention is not limited to the above-described embodiments, and various modifications can be made without departing from the scope of the invention at the stage of implementation. Furthermore, the above embodiment includes various stages, and various inventions can be extracted by appropriately combining a plurality of disclosed constituent elements. For example, some constituent requirements may be deleted from all the constituent requirements shown in the embodiment.
[0093]
【The invention's effect】
According to the present invention, it is possible to optimize the bonding position of images between adjacent stages. Thereby, the trouble of adjusting the bonding position is reduced, and the images can be bonded so that the blood vessels are naturally connected without deforming the images.
[Brief description of the drawings]
FIG. 1 is a diagram showing a configuration of an X-ray diagnostic apparatus according to an embodiment of the present invention.
FIG. 2 is a diagram showing frame movement during stepping DSA shooting in the present embodiment.
FIG. 3 is a halftone image showing an example of DSA images photographed at five stages and a panoramic view image obtained by combining them in the present embodiment.
FIG. 4 is a diagram showing an overall processing flow of creating a panoramic view image including automatic alignment and manual alignment in the present embodiment.
FIG. 5 is a diagram showing parameters related to parallel movement in the present embodiment.
FIG. 6 is a diagram showing parameters related to displacement deformation in the present embodiment.
FIG. 7 is a view showing parameters relating to enlargement / reduction in the present embodiment;
FIG. 8 is a diagram showing parameters related to peripheral portion deformation in the present embodiment.
FIG. 9 is a diagram showing a coordinate system representing the position of a pixel in a photographed m × n pixel image as an initial image coordinate system in the present embodiment.
FIG. 10 is a diagram showing a coordinate system after performing peripheral deformation on an image as a post-deformation image coordinate system in the present embodiment;
FIG. 11 is a diagram showing a state in which images of respective stages are arranged in a global coordinate system in the present embodiment.
FIG. 12 is a diagram showing an image valid data area and an ROI setting valid area in the present embodiment.
FIG. 13 is a diagram showing automatic ROI setting by threshold processing in the present embodiment.
FIG. 14 is a diagram showing a processing procedure for automatic ROI setting in the present embodiment.
FIG. 15 is a diagram showing filter processing in the present embodiment.
FIG. 16 is a diagram showing an extracted ROI region and a product region of the filtering effective region and the filtering effective region of the upper image as regions to be filtered by the test program in the present embodiment.
FIG. 17 is a diagram showing processing for picking up evaluation points from pixels in the ROI in the present embodiment;
FIG. 18 is a diagram showing two-dimensional interpolation in the case of a 2 × 2 kernel in the present embodiment.
FIG. 19 is a diagram showing a data flow of automatic alignment processing in the present embodiment.
FIG. 20 is a view showing a synthesis process of a panoramic view image in the present embodiment.
FIG. 21 is a view showing a procedure for synthesizing a panoramic view image in the present embodiment.
FIG. 22 is a view showing a panoramic view image synthesis process in which the image is not deformed in the present embodiment;
FIG. 23 is a diagram showing examples of display images of various bonding methods in the present embodiment.
[Explanation of symbols]
DESCRIPTION OF SYMBOLS 11 ... X-ray tube, 12 ... X-ray aperture, 13 ... Top plate, 14 ... Detector, 15 ... Frame, 16 ... Analog-digital converter, 17 ... Memory | storage device, 18 ... Digital image processing apparatus, 19 ... Display monitor, 20 ... operation unit, 21 ... photographing system geometric information generator.

Claims (20)

The bonding position between a plurality of medical images partially overlapping about the subject, the position shift amount of the image, the amount of deformation of the image, while adjusting the contrast conversion ratios of images, in overlapping the portion between the medical image Means for determining the difference in pixel values between parts to be minimum;
An apparatus for pasting medical images, comprising: means for pasting the plurality of medical images based on the determined joint positions by enlarging and reducing the original image without deformation.   2. The medical image pasting apparatus according to claim 1, wherein the position shift includes a parallel movement amount of the medical image, and the deformation includes a deformation amount.   The medical image pasting according to claim 2, wherein the deformation amount includes a deformation amount relating to the upper side of the medical image and a deformation amount relating to the lower side of the medical image independent of the deformation amount relating to the upper side. Alignment device.   The medical image pasting apparatus according to claim 1, wherein the deformation amount includes a deformation amount in the horizontal direction of the medical image and does not include a deformation amount in the vertical direction of the medical image.   2. The medical image pasting apparatus according to claim 1, further comprising means for initially pasting and displaying the plurality of medical images based on the position of the top plate at the time of photographing.   The medical image pasting apparatus according to claim 1, wherein the joining position is determined by being limited to a region of interest in the medical image.   Extracting a region of interest from the region where the region of interest setting effective region in the image data effective region in the upper medical image and the region of interest setting effective region in the image data effective region in the lower medical image overlap. 7. The medical image pasting apparatus according to claim 6, wherein the region of interest is determined by:   8. The medical image according to claim 7, wherein a region showing brightness higher or lower than a predetermined threshold is extracted from the overlapping region, and the region of interest is determined by enlarging the extracted region. Bonding equipment.   9. The medical image pasting apparatus according to claim 8, wherein the determined region of interest is run-length encoded.   The medical image pasting apparatus according to claim 1, further comprising means for performing a filtering process as a pre-process on each of the medical images.   11. The medical image pasting apparatus according to claim 10, wherein the filtering process is limited to each region of interest with respect to each of the medical images.   2. The medical image pasting apparatus according to claim 1, wherein the final registration amount is determined so that a difference in pixel values at corresponding positions between adjacent medical images is reduced.   For coordinate conversion of medical images that are vertically adjacent in the plurality of medical images, parameters relating to the deformation amount applied to the upper medical image, parameters relating to the deformation amount applied to the lower medical image, and parallelism between the medical images adjacent to each other vertically. The medical image pasting apparatus according to claim 1, further comprising a parameter relating to a movement amount.   14. The medical image pasting apparatus according to claim 13, wherein the coordinate conversion of the medical images adjacent above and below includes a parameter related to contrast conversion with respect to the lower medical image. The final alignment amount is determined so that a difference in pixel values at corresponding positions between adjacent medical images is small.
The corresponding position is calculated for each pixel in the region of interest by scanning pixels based on run-length encoded region-of-interest data set in the overlapping region of the adjacent medical images. The medical image pasting apparatus according to claim 1, wherein: The joint position is determined so that a difference in pixel value at a corresponding position between adjacent medical images is small,
The corresponding position is scanned based on the run-length encoded region-of-interest data set in the overlapping region of the adjacent medical images, and is randomly picked up from the pixels in the region of interest. The medical image combining device according to claim 1, wherein the device is calculated for a plurality of pixels.   The joint position is determined so as to minimize or maximize an objective function including a square error of a pixel value at a corresponding position between adjacent medical images and a penalty function. Medical image pasting device.   The medical image pasting apparatus according to claim 17, wherein the objective function includes a term related to a difference in deformation rate between adjacent medical images.   2. The medical image pasting apparatus according to claim 1, wherein the joining position is determined using one deformation rate and the other deformation rate of adjacent medical images.   The joint position is obtained by manually tracing a line obtained by manually tracing a blood vessel image included in one medical image of medical images adjacent in the vertical direction, and a blood vessel image included in the other medical image. 2. The medical image pasting apparatus according to claim 1, wherein the distance in the vertical direction between a plurality of intersections with the line is determined to be small.

Images