JP6905323B2 - 画像処理装置、画像処理方法、およびプログラム - Google Patents

Description

本発明は、核磁気共鳴映像装置（ＭＲＩ）、Ｘ線コンピュータ断層撮影装置（Ｘ線ＣＴ）、超音波画像診断装置（ＵＳ）など、種々の撮像装置（モダリティ）で撮像した３次元画像を処理する画像処理装置、画像処理方法、およびプログラムに関する。

３次元画像（被検体内部の情報を表す３次元断層画像）を用いた画像診断において、医師は、複数の撮像装置（モダリティ）、異なる体位、時刻、撮像パラメータ等で撮像した画像を対比しながら診断を行う。しかし、画像間で被検体の形状が異なるため、病変部の同定や対比を行うことが困難である。そこで、複数画像間の変形位置合わせ（変形推定）を行うことが試みられている。これにより、一方の画像に変形を施して、画像中に描出されている被検体の位置や形状を他方の画像に略一致させた変形画像を生成することが可能となる。また、一方の画像上で注目した点の、他方の画像上における対応点の位置を算出して提示することが可能となる。その結果、医師は、複数画像間における病変部の同定や対比を容易に行うことが可能となる。医療以外の分野においても、３次元画像を用いて物体の内部状態を検査する目的において、同様の作業が実施される場合がある。

このとき、画像には、様々な臓器や体組織が含まれており、種類によって組織の硬さが異なる。例えば、骨は非常に硬いため、比較する画像間で被検体の姿勢や形状が変わっても、殆ど形状が変形しない。また、腫瘍などの病変部も種類によっては周囲の組織よりも硬く、形状が変形しづらいものが存在する。従って、このような硬い組織を周囲の柔らかい組織と同様に扱って変形位置合わせ（変形推定）を行うと、実際には変形していない硬い領域が変形していると推定されて、誤った位置合わせ結果が得られてしまう場合がある。

この問題に対する解決手段として、特許文献１では、非剛体位置合わせによって得られた変位場において、剛体であるべき硬い関心領域の変位場のみを剛体変換に近付けることで、関心領域が変形しているという誤った推定を回避する技術が提案されている。より具体的には、特許文献１には、関心領域の変形のみを剛体変換に近似し、関心領域に関しては近似した剛体変換の変位場を生成し、非関心領域に関しては元の変形位置合わせの変位場を生成し、それらを空間的に結合した変位場を生成する技術が記載されている。ここで、変位場とは２画像間の各位置の間の変換（変位）を保持した場のことである。このように、特許文献１では、画像内の領域によって性質の異なる変換モデルを用いて、画像間の変形位置合わせ（変形推定）を行う技術が開示されている。ここで、変換モデルとは、位置合わせにおいて画像上の各位置の座標を変換させるためのモデルを表す。

特開２０１３−１４１６０３号公報

Ｋ．Ｒｏｈｒ，Ｈ．Ｓ．Ｓｔｉｅｈｌ，Ｒ．Ｓｐｒｅｎｇｅｌ，Ｔ．Ｍ．Ｂｕｚｕｇ，Ｊ．Ｗｅｅｓｅ，ａｎｄ Ｍ．Ｈ．Ｋｕｈｎ，"Ｌａｎｄｍａｒｋ−Ｂａｓｅｄ Ｅｌａｓｔｉｃ Ｒｅｇｉｓｔｒａｔｉｏｎ Ｕｓｉｎｇ Ａｐｐｒｏｘｉｍａｔｉｎｇ Ｔｈｉｎ−Ｐｌａｔｅ Ｓｐｌｉｎｅｓ"，ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｍｅｄｉｃａｌ Ｉｍａｇｉｎｇ，Ｖｏｌ．２０，Ｎｏ．６ Ｊｕｎｅ２００１．

しかし、特許文献１の技術では、関心領域と非関心領域とで性質の異なる変換モデルで座標変換を表し、その間を補間してつなぎ合わせることで全体の変形位置合わせを行っている。そのため、特にそのつなぎ目で整合の取れた変形が得られない場合があるという課題があった。

本発明は、上記課題に鑑みてなされたものであり、画像全体として整合の取れた変形を推定することを目的とする。

上記目的を達成するための一手段として、本発明の画像処理装置は以下の構成を備える。すなわち、画像処理装置は、被検体を異なる条件で撮像することにより得られた第１の画像と第２の画像を取得するデータ取得手段と、前記第１の画像内の関心領域を取得する関心領域取得手段と、前記第１の画像と前記第２の画像の間の第１の変形変位場を取得する第１の変形取得手段と、前記第１の変形変位場を、前記第１の変形変位場より自由度が少ない近似変換モデルを用いて近似して近似変位場を生成する変形近似手段と、前記関心領域に関して、前記第１の画像と前記第１の画像が前記近似変位場に基づいて変位される画像との間の対応情報を生成する対応情報生成手段と、前記対応情報を用いて、前記第１の画像と前記第２の画像の間の第２の変形変位場を取得する第２の変形取得手段と、を備えることを特徴とする。

本発明によれば、関心領域と非関心領域とをひとつの変換モデルでシームレスに記述できるので、画像全体として整合の取れた変形を推定できる。

第１の実施形態における画像処理システムの構成を示す図。 第１の実施形態における全体の処理手順を示すフローチャート。 同一の被検体に関する第１の画像と第２の画像を示す図。 第１の画像と代表点に対して第１の変形を適用する様子を説明する図。 関心領域の第１の変形を剛体変換で近似する様子を説明する図。 関心領域上に近似仮想対応点を生成する様子を説明する図。 第１の画像に対して第２の変形を適用する様子を説明する図。 第２の実施形態における画像処理装置の構成を示す図。 第２の実施形態における全体の処理手順を示すフローチャート。 第３の実施形態における画像処理装置の構成を示す図。 第３の実施形態における全体の処理手順を示すフローチャート。 中間代表点を示す図。 第４の実施形態における画像処理装置の構成を示す図。 第４の実施形態における全体の処理手順を示すフローチャート。 第５の実施形態における画像処理装置の構成を示す図。 第５の実施形態における全体の処理手順を示すフローチャート。 第５の実施形態における第１の変形を取得する処理の詳細な手順を示すフローチャート。

以下、添付図面に従って本発明における画像処理装置の好ましい実施形態について詳説する。ただし、発明の範囲は図示例に限定されるものではない。

＜第１の実施形態＞
本実施形態における画像処理装置は、３次元画像である第１の画像と第２の画像の変形位置合わせを行う。このとき、硬い病変や骨などの硬くあるべき領域が関心領域として第１の画像上で与えられた場合には、画像処理装置は、関心領域を剛体として維持したまま、画像全体として整合性の取れた変形位置合わせを行う。また、画像処理装置は、剛体に限らず、全体よりも自由度が小さい変形を関心領域に適用させた変形位置合わせを行うことができる。ここで、変形位置合わせとは、第１の画像中の披検体の形状が第２の画像中の披検体の形状に合致するように第１の画像に施す変形を推定することを意味する。また、そのように第１の画像を変形させた変形画像を生成することを意味する。

本実施形態における画像処理装置は、具体的には、まず、位置合わせの適切さを評価する所定の評価関数に基づき、２つの画像間の第１の変形位置合わせを行う。次に、画像処理装置は、第１の変形位置合わせで得た第１の画像内の関心領域の局所的な変形を、剛体変換のような自由度の小さい変換で近似する。そして、画像処理装置は、関心領域上の複数の代表点の夫々について、近似した変換を用いて代表点を第１の画像から第２の画像に変位させた点（変位点）を得て、代表点と変位点とからなる２画像間で仮想的に対応する対応点（対応する点の対）を生成する。最後に、画像処理装置は、生成した複数の仮想的な対応点を２つの画像間の位置合わせの拘束条件として位置合わせの評価関数に追加し、再度、２つの画像間の第２の変形位置合わせを行う。これにより、硬い病変や骨などの硬くあるべき領域を関心領域（変形を抑制したい領域）として設定することで、評価関数に基づく適切な変形を画像全体に適用しつつ、関心領域の形状を維持した変形を推定できる。つまり、画像全体として整合が取れた変形を推定できる。

以下、図１から図７を用いて、本実施形態の構成及び処理を説明する。図１は、本実施形態における画像処理装置１００の構成を示す図である。同図に示すように、画像処理装置１００は、データサーバ１１０および表示部１２０と接続されている。なお、画像処理装置１００は、データサーバ１１０と表示部１２０を含む構成としてもよい。

データサーバ１１０が保持する第１および第２の画像は、異なる条件（異なるモダリティ、撮像モード、日時、体位等）で被検体を予め撮像して得られた３次元断層画像（ボリュームデータ）である。３次元断層画像を撮像するモダリティは、ＭＲＩ（Magnetic Resonance Imaging）装置、Ｘ線ＣＴ（Computed Tomography）装置、３次元超音波撮影装置、光音響トモグラフィ装置、ＰＥＴ（Positron Emission Tomography）／ＳＰＥＣＴ（Single Photon Emission Computed Tomography）、ＯＣＴ（Optical Coherence Tomography）装置などであってもよい。第１および第２の画像は、例えば、異なるモダリティや異なる撮像モードで同時期に撮像されたものであってもよい。また、経過観察のために同一患者を同一モダリティ、同一体位で異なる日時に撮像した画像であってもよい。第１および第２の画像は、データ取得部１０１０を介して画像処理装置１００に入力される。

表示部１２０は液晶モニタ等であり、画像処理装置１００が生成する表示画像等の各種の情報を表示する。また、表示部１２０には、ユーザからの指示を取得するためのＧＵＩ（Graphical User Interface）も配置されている。このＧＵＩは、以下の説明における操作部としても機能する。

画像処理装置１００は、データ取得部１０１０、関心領域取得部１０２０、第１の変形取得部１０３０、変形近似部１０４０、対応情報生成部１０５０、第２の変形取得部１０６０、表示制御部１０７０から構成される。

データ取得部１０１０は、画像処理装置１００へと入力される第１および第２の画像をデータサーバ１１０から取得する。関心領域取得部１０２０は、第１の画像上の変形を抑制したい領域（以下、関心領域）の情報を取得する。第１の変形取得部１０３０は、第１の画像と第２の画像との間の第１の変形位置合わせを行い、第１の画像を変形させるための第１の変形変位場を取得する。変形近似部１０４０は、関心領域における第１の変形変位場を、この変位場よりも自由度の小さい近似変換モデルを用いて近似することで、近似変位場を生成する。対応情報生成部１０５０は、画像間における仮想的な対応情報（仮想対応情報）として、近似変位場に基づき、関心領域上における仮想的な対応点（近似仮想対応点）の情報を生成する。第２の変形取得部１０６０は、画像間の一致度を評価する評価関数に仮想対応情報に基づく拘束条件を追加した評価関数に基づき、第１の画像と第２の画像との間で第２の変形位置合わせを行い、第１の画像を変形させるための第２の変形変位場を取得する。表示制御部１０７０は、第１の画像、第２の画像、位置合わせ結果に基づく変形画像の断層画像などの情報を表示部１２０に表示させる制御を行う。

図２は、画像処理装置１００が行う全体の処理手順を示すフローチャートである。

（Ｓ２０００）（データの取得）
ステップＳ２０００において、データ取得部１０１０は、データサーバ１１０から、第１の画像と第２の画像を取得する。そして、データ取得部１０１０は、第１および第２の画像を、関心領域取得部１０２０、第１の変形取得部１０３０、第２の変形取得部１０６０、および表示制御部１０７０へと出力する。

図３は、同一の被検体に関する第１の画像３０００と第２の画像３０３０を示す図である。図３の例では、第１の画像と第２の画像は異なる種類のモダリティの撮像画像であり、被写体が人体における乳房である場合を示している。それぞれの画像は、例えばＣＴやＭＲＩ、３Ｄ超音波画像などの３次元のボリュームデータを表すが、紙面が２次元である都合上、ボリュームデータからｘ−ｙ平面に平行に切り出された断面画像の形式で示されている。

第１の画像３０００には、被写体の領域３０１０と病変領域３０２０が示されている。同様に、第２の画像３０３０には、被写体の領域３０４０と病変領域３０５０が示されている。病変領域３０２０と病変領域３０５０は、それぞれ解剖学的に対応する領域である。図３の例では、第１の画像と第２の画像の間で、乳房である被検体の領域３０１０と３０４０は大きく変形しているが、病変領域３０２０と３０５０は、解剖学的に硬い病変であるため、あまり変形していない状況を表している。

（Ｓ２０１０）（関心領域の取得）
ステップＳ２０１０において、関心領域取得部１０２０は、データ取得部１０１０から取得した画像上において、関心領域の情報（変形を抑制したい病変部等の情報）を取得する。そして、関心領域取得部１０２０は、取得した関心領域の情報を、対応情報生成部１０５０へと出力する。

ここで、関心領域は、例えば、硬い病変部等の領域や骨の領域、胸壁、大胸筋など、撮像画像間で変形しづらい領域であってもよい。さらには、第１の画像と第２の画像が異なる日時に撮像された経時画像である場合には、関心領域は、経過観察の対象となる病変領域であってもよい。これは、治療効果の判定等を目的として、経時画像間を変形位置合わせした上で病変領域の経過観察をしたい場合に、病変領域まで変形してしまっては、経時変化による病変の大きさ等の変化が観察できないためである。そこで、実際には剛体ではない病変であっても、画像間で病変領域を剛体（形状を変化させない領域）として扱うことで、大きさ等の変化の観察が容易となる。

本実施形態では、例えば、関心領域取得部１０２０は、図３の第１の画像３０００上で輪郭が明瞭に描出されている病変領域３０２０を関心領域として取得する。本実施形態では、ユーザが不図示の操作部を操作して画像上の病変部の輪郭形状を手動で抽出することにより、関心領域取得部１０２０は、関心領域の情報を取得する。取得された関心領域としての病変領域３０２０は例えば、関心領域内外が２値化されたマスク画像Ｉ_ＲＯＩで表現される。なお、関心領域は、関心領域の内外が区別できる形式であれば、どのような形式で表現されてもよい。また、関心領域は必ずしも１箇所である必要はなく、複数の関心領域が取得されてもよい。

なお、関心領域取得部１０２０は、関心領域を取得する処理を、一般的に用いられている画像解析技術により行ってもよい。例えば、変形を抑制する病変領域内の座標をユーザがシード点として与えた後に、関心領域取得部１０２０は、領域拡張法によって領域を抽出してもよい。なお、病変領域の抽出方法はこれに限られるものではなく、ＳＮＡＫＥＳやＬｅｖｅｌＳｅｔ法などの公知の領域分割手法で行ってもよい。また、関心領域取得部１０２０が複数の病変候補領域を検出する一般的な画像処理を第１または第２の画像に施した後に、変形を制御する領域をユーザが不図示の操作部を操作して指定（選択）することにより、関心領域取得部１０２０は関心領域を取得してもよい。

（Ｓ２０２０）（第１の変形の取得）
ステップＳ２０２０において、第１の変形取得部１０３０は、データ取得部１０１０から取得した第１の画像と第２の画像との間の第１の変形位置合わせを実行し、第１の変形変位場を取得する。そして、第１の変形取得部１０３０は、取得した第１の変形変位場を対応情報生成部１０５０へと出力する。さらに、第１の変形取得部１０３０は、第１の変形変位場に基づいて第１の画像を変形させた第１の変形画像を生成して、表示制御部１０７０へと出力する。

このとき、第１の変形取得部１０３０は、変形Φによる位置合わせの適切さを評価する評価関数（コスト関数）を定義し、該評価関数を最小化するような変形Φを推定する。このとき推定された変形Φを第１の変形と称する。変形Φは変位場（変形変位場）によって表現される。ここで、変位場の形態は、画像上の各位置を変位させる変位ベクトルからなる変位ベクトル場であってもよいし、空間的な場で表現されておらずとも、画像上の任意の位置を変位可能な所定の変換モデルが有する変換パラメータの形で表わされていてもよい。

ここで、本ステップで定義する評価関数を第１の評価関数Ｅ１（Φ）と称する。本実施形態では、第１の評価関数を構成する指標として、変形による画像間の一致度を適用する。例えば、第１の画像と第２の画像の間における解剖学的に対応する特徴点（対応点）の位置の誤差を適用する。後述する仮想対応点と区別するために、このような画像間で実際に対応する対応点を実対応点と称する。特に、上述のように解剖学的に画像間で対応する実対応点をランドマーク対応点と称する。すなわち、第１の変形取得部１０３０は、ランドマーク対応点の位置を極力一致させるような変形を推定する。

第１の変形取得部１０３０は、対応点を取得するために、例えば、まず第１の画像、第２の画像に、それぞれインタレストオペレータ等の公知の画像特徴点抽出処理を施して、それぞれの画像から血管分岐等の解剖学的な特徴点を抽出する。その後、第１の変形取得部１０３０は、抽出されたそれぞれの特徴点の間で、１対１に対応付けすることによってランドマーク対応点の組を生成する処理を行う。この処理は、各画像上の特徴点ごとに特徴点近傍に関心領域を設定し、画像間で関心領域の画像類似度が高くなる特徴点同士を対応付けすることで行われる。

なお、ランドマーク対応点の取得方法は、上記の方法に限られるものではない。例えば、表示部１２０上に表示された第１の画像と第２の画像に対して、ユーザが不図示の操作部を用いて手動で入力した対応点により、ランドマーク対応点が取得されてもよい。また、第１の評価関数Ｅ１（Φ）は、上記のように実対応点位置の誤差に限られるものではなく、例えば、変形させた第１の画像と第２の画像との間の画像類似度を指標としてもよい。また、評価関数Ｅ１（Φ）として、上述の画像全体の類似度や実対応点位置の誤差の何れか一方ではなく、両方を併用してもよい。

第１の変形取得部１０３０は、第１の評価関数Ｅ１を最小化するように、変形Φを表現する所定の変換モデルの変換パラメータを最適化する。これにより、画像内の被検体全体の位置合わせとして望ましい最適な変形が推定される。ここで、所定の変換モデルで表現された変換パラメータをｐと定義する。本実施形態では、所定の変換モデルとして、変形の基底関数がＢスプライン関数で表されるＦＦＤ（Free Form Deformation）を適用する。ＦＦＤにおける変換パラメータｐは、画像内に規則的に配置された制御点の制御量で表現される。なお、本ステップで推定された変換パラメータｐを第１の変換パラメータｐ１と表記する。なお、適用する変換モデルは、ＦＦＤに限定されるものではない。例えば、ＴＰＳ（Thin Plate Spline）などの放射基底関数や、ＬＤＤＭＭ（Large Deformation Diffeomorphic Metric Mapping）など、変形を表現できる公知の何れの変換モデルであってもよい。

このように、第１の変形取得部１０３０は、第１の変形変位場として、変換パラメータｐ１を生成する。そして、第１の変形取得部１０３０は、生成した第１の変形変位場に基づき、第２の画像に形状が一致するように第１の画像を変形させた第１の変形画像を生成する。また、必要に応じて、第１の画像の形状に一致するように第２の画像を変形させた変形画像（逆変形画像）を、第１の変形変位場に基づいて生成する。

（Ｓ２０３０）（関心領域における第１の変形を近似変換で近似）
ステップＳ２０３０において、対応情報生成部１０５０および変形近似部１０４０は、関心領域における第１の変形変位場よりも自由度の小さい変換モデルで近似した近似変位場を生成する。そして、対応情報生成部１０５０は、生成した近似変位場を、第２の変形取得部１０６０へと出力する。

以下では、近似変位場を生成するための変換モデル（すなわち、第１の変形変位場よりも自由度の小さい変換モデル）を、近似変換モデルと称する。例えば、近似変換モデルとして、剛体変換モデルを適用することができる。第１の変形変位場は非剛体の変形位置合わせの結果であるため、剛体変換モデルは、第１の変形変位場よりも自由度が小さいモデルといえる。より具体的な処理を以下に説明する。

まず、対応情報生成部１０５０は、関心領域上の複数の代表点の夫々について、代表点と、それを第１の変形変位場を用いて変位させた変位点とからなる、画像間の仮想的な対応点（変形仮想対応点）を生成する。そして、対応情報生成部１０５０は、生成した変形仮想対応点を変形近似部１０４０へと出力する。

この、対応情報生成部１０５０の処理について、具体的な処理を説明する。まず、対応情報生成部１０５０は、関心領域上に代表点を設定する。ここで、代表点を設定する領域は、関心領域の内部領域のみであってもよいし、関心領域を囲む所定範囲内の領域であってもよい。例えば、代表点を設定する領域を、関心領域の輪郭から周囲５ｍｍ以内の領域に設定できる。また、対応情報生成部１０５０は、代表点を、領域内に均等な間隔（０．２５ｍｍ等）で配置してもよいし、ランダムに配置してもよい。例えば、対応情報生成部１０５０は、関心領域上に等間隔に複数の代表点を設定する。次に、対応情報生成部１０５０は、第１の変形変位場を用いて夫々の代表点の位置を変位させて、夫々の代表点に対応する変位点の位置を算出する。算出した変位点を変形代表点と称する。本実施形態では、代表点とそれに対応する変形代表点とを、第１の変形による画像間の仮想的な対応点と定義し、以下ではこれを、変形仮想対応点と称する。

図４は、第１の変形変位場に基づいて変形仮想対応点を生成する様子を示す図である。図４には、第１の変形変位場に基づいて第１の画像３０００を変形させた第１の変形画像４０００が表されている。第１の変形画像４０００には、変形後の被写体の領域４０１０と変形後の関心領域４０２０が示されている。また、代表点４０３０は、関心領域としての病変領域３０２０上に設定された代表点、変形代表点４０４０は、代表点４０３０を第１の変形場に基づいて変形させた代表点を表す。図４の例では、第１の変形位置合わせによって、関心領域４０２０が被検体の領域４０１０と同様に押しつぶされると同時に、並進移動と回転が生じているように推定された様子が示されている。このとき、変形代表点４０４０の分布も、関心領域４０２０と同様に、代表点４０３０が全体的に押しつぶされて、さらに並進移動と回転が生じたような分布として算出される。

次に、変形近似部１０４０は、生成した変形仮想対応点に基づき、変形仮想対応点の変位を近似変換モデルで近似した近似変位場を生成する。より具体的には、変形近似部１０４０は、式（１）に示す誤差の和ｅを最小にする近似変換モデルΦ’のパラメータ（近似変換パラメータ）ｑを算出する。近似変換モデルとして剛体変換モデルを用いる場合には、物体の移動と回転を表す３×４の剛体変換行列Ｔを規定する位置と姿勢の６パラメータが、近似変換パラメータｑに相当する。すなわち、式（２）のように、変形近似部１０４０は、剛体変換行列及び代表点の積と、変形代表点との差のノルムの総和ｅを算出し、総和ｅが最小となる剛体変換行列Ｔを剛体変換の制約下で算出する。

ここで、Φ’（ｘ｜ｑ）は、パラメータｑで表される近似変換モデルΦ’によって座標ｘを変位させた座標を表す。また、（ｘ_Ｄ１ｉ，ｘ_Ｄ２ｉ）は、第１の画像と第２の画像上におけるｉ番目の変形仮想対応点の座標の組を表す。すなわち、ｘ_Ｄ１ｉがｉ番目の代表点の座標、ｘ_Ｄ２ｉがｉ番目の変形代表点の座標を表す。また、Ｎ_Ｄは、変形仮想対応点の数を表す。なお、行列Ｔは特異値分解などを用いる公知の方法で算出することができるので、算出方法の説明は省略する。

このようにして、変形近似部１０４０は、関心領域における第１の変形変位場の近似変位場を取得する。最後に、変形近似部１０４０は、取得した近似変換パラメータｑを、近似変位場の情報として対応情報生成部１０５０へと出力する。

図５は、関心領域としての病変領域３０２０における第１の変形を剛体変換で近似した様子を説明する図である。図５において、座標軸５０００は、関心領域としての病変領域３０２０上の代表点（図４の代表点４０３０）の重心位置を原点とし、第１の画像３０００のｘ‐ｙ軸に平行な座標軸を表している。座標軸５０１０は、剛体変換行列Ｔを用いて座標軸５０００を変換した後の座標軸を表している。但し、紙面が平面である都合上、説明を分かりやすくするためにｚ軸に平行な軸に関しては省略している。また、座標軸５０１０のｘ‐ｙ軸がなす面は第１の変形画像４０００のｘ‐ｙ軸がなす面と必ずしも平行になるものではないが、簡略化して図示している。図に示すように、座標軸５０１０が、図４の代表点４０３０から変形代表点４０４０への変換（第１の変形変位場）から、領域を押しつぶすような変形成分のみを除いた、並進移動・回転成分を抽出した変換に一致することが分かる。

（Ｓ２０４０）（近似変換に基づく近似仮想対応情報の生成）
ステップＳ２０４０において、対応情報生成部１０５０は、画像間における仮想的な対応情報として、関心領域上に近似変位場に基づく複数の仮想的な対応点（仮想対応情報）を生成する。本実施形態において、以下ではこれを、近似仮想対応点と称する。そして、生成した近似仮想対応点を第２の変形取得部１０６０へと出力する。

以下に、具体的な処理を説明する。まず、対応情報生成部１０５０は、ステップＳ２０３０の処理と同様に、関心領域上に設定した複数の代表点の夫々の位置を、近似変位場Φ’を用いて変位させて、夫々の代表点に対応する変位点の位置を算出する。算出した変位点を近似代表点と称する。そして、代表点とそれに対応する近似代表点とを、近似変換モデルによる画像間の仮想的な対応点と定義し、以下ではこれを、近似仮想対応点と称する。

図６は、近似変位場に基づいて近似仮想対応点を生成する様子を示す図である。近似代表点６０００は、代表点４０３０を近似変形場に基づいて変換させた近似代表点を表す。図に示すように、近似代表点６０００の分布が、変形は生じずに、図５の座標軸５０００から座標軸５０１０への変換が示す並進移動と回転のみが生じた分布になっていることが分かる。

（Ｓ２０５０）（近似仮想対応点に基づく第２の変形の取得）
ステップＳ２０５０において、第２の変形取得部１０６０は、ステップＳ２０４０で取得した近似仮想対応点（仮想対応情報）を拘束条件として用いて、第１の画像と第２の画像との間の第２の変形位置合わせを実行する。そして、第２の変形取得部１０６０は、第１の画像と第２の画像との間の変形推定の結果として、生成した第２の変形変位場をデータサーバ１１０へと出力する。データサーバ１１０への出力は、不図示の出力部を用いて行われ得る。また、第２の変形取得部１０６０は、第２の変形変位場に基づいて第１の画像を変形させた第２の変形画像を生成して、表示制御部１０７０へと出力する。

このとき、第２の変形取得部１０６０は、ステップＳ２０２０と同様に、変形Φによる位置合わせの適切さを評価する評価関数を最小化するような変形Φを推定する。このとき推定された変形Φを第２の変形と称する。また、本ステップで定義する評価関数を第２の評価関数Ｅ２（Φ）と称する。本実施形態では、第２の評価関数として、ステップＳ２０２０で用いた第１の評価関数Ｅ１（Φ）に対して近似仮想対応点の対応点位置の誤差項を追加した評価関数を適用する。

ここで、変形Φによる近似仮想対応点の対応点位置の誤差項をＲ（Φ）と定義すると、第２の評価関数Ｅ２（Φ）は以下の式で表される。

ここで、Φ（ｘ｜ｐ２）は、パラメータをｐ２とした所定の変換モデルで表される第２の変形Φによって座標ｘを変位させた座標を表す。また（ｘ_Ｒ１ｉ，ｘ_Ｒ２ｉ）は、第１の画像と第２の画像上におけるｉ番目の近似仮想対応点の座標の組を表す。すなわち、ｘ_Ｒ１ｉがｉ番目の代表点の座標、ｘ_Ｒ２ｉがｉ番目の近似代表点の座標を表す。また、Ｎ_Ｒは、近似仮想対応点の総数を表す。すなわち、Ｒ（Φ）は、第２の変形Φによる変位後の代表点を、対応する近似代表点になるべく一致させるための項である。

第２の変形取得部１０６０は、第２の評価関数Ｅ２（Φ）を最小化するように、変形Φを表現する所定の変換モデルの変換パラメータｐ２を最適化する。これにより、画像内の被検体全体の位置合わせとして望ましい最適な変形が推定される。本実施形態では、所定の変換モデルとして、ステップＳ２０２０で用いたものと同じ変換モデルを適用する。すなわち、ステップＳ２０２０でＦＦＤを適用した場合は、本ステップにおいてもＦＦＤを適用する。ここで、変換モデルがＦＦＤであって、第２の評価関数Ｅ２（Φ）が、実対応点（ランドマーク対応点）位置の誤差項と近似仮想対応点位置の誤差項のみによって表現される場合を考える。この場合、ＦＦＤのパラメータと対応点間の位置の誤差との関係は、線形方程式で表すことができる。したがって、第２の変形取得部１０６０は、対応点間の位置の誤差を最小化する変換パラメータｐ２を公知の線形最適化手法によって推定する。これにより、高速に変換パラメータｐ２を算出可能である。また、変換モデルがＴＰＳの場合も同様に、第２の変形取得部１０６０は、公知の線形最適化手法によってパラメータｐ２が算出できる。なお、評価関数Ｅ２（Φ）として画像間の類似度や非線形な正則化項を考慮する場合や、変換モデルとしてＬＤＤＭＭを用いる場合には、第２の変形取得部１０６０は、公知の非線形最適化によってパラメータｐ２を推定する。このとき、第２の変形取得部１０６０は、ステップＳ２０２０で第１の変形取得部１０３０が算出したパラメータｐ１を、パラメータｐ２の初期値として用いることができる。これによると、パラメータｐ２の算出をより少ない反復処理で行うことができる。

このように、第２の変形取得部１０６０は、第２の変形変位場として、変換パラメータｐ２を生成する。そして、第２の変形取得部１０６０は、生成した第２の変形変位場に基づき、第２の画像に形状が一致するように第１の画像を変形させた第２の変形画像を生成する。また、必要に応じて、第２の変形取得部１０６０は、第１の画像の形状に一致するように第２の画像を変形させた変形画像（逆変形画像）を、第２の変形変位場に基づいて生成する。

このように、本実施形態では、第１の評価関数Ｅ１（Φ）に基づくことで被検体全体として適切な位置合わせを実現しつつ、近似仮想対応点位置の誤差項Ｒ（Φ）に基づくことで関心領域をなるべく近似変換モデルが表す形状（例えば剛体）に保つような変形を推定できる。さらに、画像全体の変形が単一の変換モデル（ＦＦＤ）で表現されるので、全体として整合性のとれた変形が推定できる。

また、本ステップで得られた第２の変形は、関心領域以外の領域を評価する評価関数内の指標、および変換モデルともに、第１の変形と同条件となる。従って、関心領域以外の領域の変形を第１の変形に類似した変形にすることができる。これにより、第２の変形において関心領域の周囲の近似変換成分も第１の変形に類似するものとなる。従って、関心領域の変形に対して、近似仮想対応点の拘束により第１の変形の関心領域の近似変換成分がそのまま適用される結果となっても、近似変換成分が関心領域の周囲と略合致するため、違和感のない整合の取れた位置合わせ結果を取得することができる。

図７は、第１の画像を第２の変形変位場に基づいて変形させた第２の変形画像を生成する様子を示す図である。第２の変形画像７０００には、変形後の被写体の領域７０１０、変形後の関心領域７０２０が示されている。図に示すように、変形後の被検体の領域７０１０は、図４における変形後の被検体の領域４０１０とほぼ同じような形状となる。一方で、変形後の関心領域７０２０は、図４における変形後の関心領域４０２０のような押しつぶされた変形は生じずに、関心領域としての病変領域３０２０の形状が剛体に保たれたまま、並進移動と回転のみ関心領域４０２０と同様の変換が行われた結果になることを示している。これは、第２の位置合わせにおいて、代表点４０３０と近似代表点６０００からなる近似仮想対応点位置の誤差項による、関心領域の変換を近似変換モデル（ここでは剛体）に拘束する拘束条件の効果が発揮されていることを示している。

（Ｓ２０６０）（変形画像の表示）
ステップＳ２０６０において、表示制御部１０７０は、ユーザによる不図示の操作部への操作に応じて、取得した変形画像と第２の画像の断面画像を、表示部１２０へと表示する制御を行う。このとき、表示制御部１０７０は、変形画像として、第１の変形変位場に基づく第１の変形画像と、第２の変形変位場に基づく第２の変形画像と、第２の画像の同一断面を、並べて表示するように構成することができる。また、表示制御部１０７０は、第１の変形画像と第２の変形画像を、ユーザの操作に応じて、１つの表示画面領域内で切り替えて表示するように制御を行ってもよい。なお、表示制御部１０７０は、第１の変形画像の生成や表示を行わなくてもよい。以上が、本実施形態における画像処理装置１００の処理の説明である。

本実施形態によれば、画像処理装置１００は、単一の変換モデルを用いて、画像全体を適切に一致させる指標と、関心領域の変形をより小さい自由度（例えば剛体）に留めるための指標とに基づいて変形を推定する。これにより、画像処理装置１００は、被検体内の関心領域以外の領域を適切に変形させつつ、関心領域に対してはより自由度の小さい変換（例えば剛体変換）を施した上で、画像全体として整合の取れた変形を推定することができる。

（変形例１−１）（近似変換モデルは剛体変換でなくともよい）
第１の実施形態では、ステップＳ２０３０で関心領域の変形を近似する近似変換モデルとして、剛体変換を適用していたが、近似変換モデルは、第１の変形よりも自由度が小さければどのようなモデルであってもよい。例えば、変形近似部１０４０は、アフィン変換によって第１の変形を近似してもよい。

具体的には、まず、対応情報生成部１０５０は、第１の実施形態と同様に、変形仮想対応点を生成する。次に、変形近似部１０４０は、生成した変形仮想対応点の変位をアフィン変換で近似したアフィン変位場を生成する。より具体的には、変形近似部１０４０は、式（２）に示す誤差の和ｅを最小にする３×４の行列Ｔを、アフィン変換の制約化で算出する。この式において、アフィン変換行列Ｔと上記のノルムの間の関係を線形方程式で表せることが一般的に知られている。従って、行列Ｔは、ＬＵ分解やＱＲ分解などの公知の方法で算出することができる。この方法の説明は公知の技術であるため、説明は省略する。これにより、例えば、関心領域が完全に剛体ではないが、複雑な変形はしないような場合に、関心領域の変形を実際に近い形で推定できる。

また、変形近似部１０４０は、関心領域の変形を単なるアフィン変換ではなく、体積保存の拘束条件が加わったアフィン変換で近似するようにしてもよい。このとき、体積保存の拘束条件は、アフィン変換行列Ｔにおける、ｘ、ｙ、ｚ方向のスケールファクター成分をそれぞれ、ｓ_ｘ、ｓ_ｙ、ｓ_ｚとすると、以下の式で表される。

すなわち、変形近似部１０４０は、ステップＳ２０３０におけるパラメータｑを算出するという問題を、式（５）を満たすという条件下で式（２）の評価関数ｅを最小化するようなアフィン変換行列Ｔのパラメータを推定する問題とみなすことができる。この行列Ｔは、公知の非線形最適化手法を用いて算出できる。これにより、例えば、関心領域が複雑な変形はせず、かつ変形の過程で圧縮・膨張が発生しないような場合に、関心領域の変形を実際に近い形で推定できる。

（変形例１−２）（拘束条件は近似仮想対応点を面状に拘束するものであってもよい）
第１の実施形態では、ステップＳ２０５０において、近似仮想対応点を構成する各対応点間の位置関係が一致する拘束条件を、第２の評価関数Ｅ２（Φ）に適用した。しかし、第２の評価関数Ｅ２（Φ）に追加する拘束条件は、必ずしもこれに限定されない。

例えば、関心領域が胸壁や大胸筋面などの面状の領域である場合に、近似仮想対応点を関心領域の面上に拘束する（面上にあれば位置が一致せずともよい）拘束条件を、第２の評価関数Ｅ２（Φ）に適用してもよい。これにより、以下のような場合に効果を期待できる。

例えば、乳房を圧迫するような変形において、乳房に外力が与えられると、乳房内の乳腺や脂肪の組織はその奥にある大胸筋面上を滑って変形していくことが知られている。このとき、大胸筋面は殆ど変形しないため、大胸筋面を関心領域として、大胸筋面を剛体に保つような変形を行うことを想定する。このとき、ステップＳ２０５０において、近似仮想対応点の位置そのものが一致するような拘束条件を第２の評価関数Ｅ２（Φ）に適用すると、大胸筋近傍の組織は、そのまま大胸筋面上の位置に拘束されてしまうため、大胸筋面上を滑るような変形にはならない。そこで、本変形例では、近似仮想対応点を関心領域（大胸筋面）の面上のみに拘束するという条件を第２の評価関数Ｅ２（Φ）に追加する。これにより、近似仮想対応点が関心領域の面上を滑る場合を許容できる。すなわち、関心領域の外面を含む関心領域上に設定された代表点は、近似変位場による変位後も、関心領域の面上を沿って動くことを許容される。これにより、大胸筋近傍の組織が大胸筋面上を滑ることを許容できるようにすることが可能となる。

以下に具体的な方法を説明する。まず、近似仮想対応点を構成する代表点は、面領域である関心領域上に均等に配置されているものとする。このとき、ステップＳ２０５０において、近似仮想対応点を関心領域の面上に拘束する拘束条件は、評価関数Ｅ２（Φ）を構成する関数Ｒ（Φ）として、以下の式で表わされる。

ここで、（６）式では、（４）式と同じ記号は、（４）式と同様の意味を表す。また、Ｃｏｖ（ｘ）は位置ｘにおける３次元の誤差に関する共分散行列を返す関数を表す。以下、共分散Ｃｏｖ（ｘ）をＣｏｖと略して表記する。

本変形例の場合、共分散Ｃｏｖは、例えば、ｉ番目の近似代表点ｘ_Ｒ２ｉ近傍の法線方向には狭く、それに直交する平面の方向には非常に広い分布を有するガウス分布を表す共分散を返す関数にできる。ここで、近似代表点ｘ_Ｒ２ｉ近傍の法線方向とは、面状に分布している近似代表点の中で、近似代表点ｘ_Ｒ２ｉ近傍の近似代表点の集合によって張られる平面の法線ベクトルの方向を意味する。この法線ベクトルは、例えば、公知の手法である、近似代表点ｘ_Ｒ２ｉ近傍の近似代表点の集合の位置を主成分分析し、その第３主成分ベクトルとして算出できる。

なお、法線ベクトルの算出方法は、これに限られるものではなく、公知の何れの方法を採ってもよい。ここで、（６）式において、第２の変形Φによる変位後の代表点と、対応する近似代表点の間の誤差を表す３次元ベクトルをｅと定義する。（６）式により、入力ベクトルにおける法線成分が大きい場合には、Ｃｏｖ^−１と誤差ベクトルｅの積は大きい値を返し、それ以外の成分が大きい場合には非常に小さい値を返す。従って、（６）式により、Ｒ（Φ）は、第２の変形Φによる変位後の代表点と、対応する近似代表点の間の誤差ｅのうち、近似代表点の法線方向の成分にのみペナルティをかけることができる。すなわち、Ｒ（Φ）は、変位後の代表点が、対応する近似代表点から法線方向にのみずれないようにするための項として機能する。

このように定義された評価関数Ｅ２（Φ）を最小化するような変換パラメータｐ２の算出方法としては、ステップＳ２０５０と同様に公知の非線形最適化手法で算出することができる。また、非特許文献１に記載された、上記共分散Ｃｏｖに基づく線形最適化手法によっても算出することができる。

（変形例１−３）
第１の実施形態では、ステップＳ２０４０およびステップＳ２０５０で利用する画像間における仮想対応情報として、仮想的な対応点を適用していた。しかし、利用する仮想対応情報は、必ずしも対応点ではなくてもよい。

例えば、仮想対応情報として、関心領域の形状に関する情報の、画像間における仮想的な対応情報を用いてもよい。具体的には、ステップＳ２０４０において、対応情報生成部１０５０は、関心領域そのものを、近似変形場Φ’を用いて変位させて、対応する変位後の関心領域（近似関心領域）を取得する。より具体的には、対応情報生成部１０５０は、関心領域を表すマスク画像Ｉ_ＲＯＩを、近似変形場Φ’を用いて変位させて、近似関心領域のマスク画像Ｉ_{Ａ＿ＲＯＩ}を取得する。そして、対応情報生成部１０５０は、関心領域と近似関心領域とを、近似変換モデルによる画像間の仮想的な対応領域として生成する。仮想的な対応領域の情報は、関心領域の輪郭形状の情報を保持しているため、本変形例における画像間の仮想対応情報として、これを近似仮想対応形状と称する。

そして、対応情報生成部１０５０は、生成した近似仮想対応形状を第２の変形取得部１０６０へ出力する。そして、ステップＳ２０５０において、第２の変形取得部１０６０は、ステップＳ２０４０で取得した近似仮想対応形状を拘束条件として用いて、第１の画像と第２の画像との間の第２の変形位置合わせを実行する。このとき、本ステップで最小化する第２の評価関数Ｅ２（Φ）において、本変形例では、第１の実施形態における近似仮想対応点の対応点位置の誤差項Ｒ（Φ）の代わりに、近似仮想対応形状の誤差項Ｒ’（Φ）を適用する。近似仮想対応形状の誤差項Ｒ’（Φ）は以下の式で表される。

ここで、関数Ｓ（Ｉ_１，Ｉ_２）は、画像Ｉ_１と画像Ｉ_２の間の画像類似度関数を表す。画像類似度を計算する尺度としては、一般的に用いられているＳＳＤ（Sum of Squared Difference）や相互情報量、相互相関係数などの公知の手法を用いるものとする。これにより、近似仮想対応形状である関心領域のマスク画像Ｉ_ＲＯＩが、近似関心領域のマスク画像Ｉ_{Ａ＿ＲＯＩ}に極力一致するように変形Φを推定することができる。つまり、
関心領域の輪郭形状を近似変換モデルが表す形状に保つような変形を推定できる。

また、仮想対応情報として、次の情報を適用してもよい。まず、対応情報生成部１０５０は、関心領域のマスク画像Ｉ_ＲＯＩと近似関心領域のマスク画像Ｉ_{Ａ＿ＲＯＩ}の夫々から、画像中の関心領域の輪郭位置からの距離場を生成する。続いて、対応情報生成部１０５０は、夫々の距離場内の各距離値を輝度値に変換して画像化した距離場画像を夫々生成する。ここで、Ｉ_ＲＯＩに対応する距離場画像をＩ_ＤＩＳＴ、Ｉ_{Ａ＿ＲＯＩ}に対応する距離場画像をＩ_{Ａ＿ＤＩＳＴ}と定義し、それらを、近似変換モデルによる画像間の仮想対応情報とする。この場合、（６）式において、Ｉ_ＲＯＩ、Ｉ_{Ａ＿ＲＯＩ}の代わりに、Ｉ_ＤＩＳＴ、Ｉ_{Ａ＿ＤＩＳＴ}を適用することで、関心領域の輪郭だけでなく、内部まで輝度情報を類似させることができる。これにより、関心領域内部も含めて関心領域Ｉ_ＲＯＩの形状を近似変換モデルが表す形状に保つような変形を推定できる。

（変形例１−４）
第１の実施形態では、ステップＳ２０３０で生成する変形仮想対応点と、ステップＳ２０４０で生成する近似仮想対応点は、第１の画像上の同じ代表点から生成するものであった。しかし、それらは同じ代表点から生成されるものである必要はなく、変形仮想対応点と近似仮想対応点とが別々の代表点から生成されるものであってもよい。例えば、変形仮想対応点に関しては、関心領域を囲む所定範囲内（周囲５ｍｍ以内）の領域上に設定された代表点から生成し、近似仮想対応点に関しては、厳密に関心領域上に設定された代表点から生成する方法を採ってもよい。これにより、第２の変形において、関心領域の周囲の領域により整合した近似変換成分を抽出でき、かつ厳密に関心領域に対してのみ形状を維持させる拘束条件を与えることができる。

（変形例１−５）（第１の変形推定と第２の変形推定でＥ１（Φ）は同じでなくてもよい）
第１の実施形態では、ステップＳ２０５０における第２の評価関数Ｅ２（Φ）は、（３）式に示すように第１の評価関数Ｅ１（Φ）を用いていた。しかし、ステップＳ２０５０で、第１の評価関数Ｅ１（Φ）を必ずしも用いる必要はない。例えば、ステップＳ２０５０では、より多くの計算コストを要して高精度な評価を行う評価関数を、Ｅ１（Φ）として用いるようにしてもよい。例えば、第１の変形推定ではランドマーク対応点位置の誤差をＥ１（Φ）として用いて、第２の変形推定では画像全体の画像類似度をＥ１（Φ）として用いるようにできる。これは、第１の変形推定は近似変形場を求めるための中間処理に過ぎないが、第２の変形推定は位置合わせの結果を得ることを目的としているからである。

（変形例１−６）（第２の変形における変換モデルは第１の変形と同じでなくてもよい）
第１の実施形態では、ステップＳ２０５０の第２の変形における変換モデルは、ステップＳ２０２０における第１の変形と同じ変形モデルを用いていた。しかし、関心領域の周囲の近似変換成分が第１の変形と類似する結果になるのであれば、第２の変形で利用する変換モデルは、必ずしもステップＳ２０２０と一致する必要はない。

例えば、本実施形態では、ステップＳ２０２０においてＦＦＤを変換モデルとして採用しているが、ＦＦＤのモデルを構成するパラメータの一つである、制御点の間隔を変えたＦＦＤをステップＳ２０５０で採用してもよい。具体的には、ステップＳ２０２０において、制御点を２０ｍｍ間隔で配置している場合に、ステップＳ２０５０において、その半分の１０ｍｍ間隔で配置するようにしてもよい。これにより、Ｂスプライン基底関数により滑らかな変形が得られるという特性は変えずに、変形の自由度を増やすことができる。つまり、より詳細な変形を表現できる。さらに、画像全体の概略の変形を第１の変形と第２の変形の間で類似させるために、例えば、第２の変形では段階的に変形の自由度を上げていく変換モデルを適用する。具体的には、第１の変形と同様の変形の自由度の低い「荒い制御点（２０ｍｍ間隔）」で表現されたＦＦＤから変形の自由度の高い「細かい制御点（１０ｍｍ間隔）」で表現されるＦＦＤの順に適用する、公知の技術である多重解像度ＦＦＤを、変換モデルとして適用する。

これにより、第２の変形において、第１の変形よりもより詳細な位置合わせを実現しながら、概略の変形としては第１の変形の状態が保持されるため、関心領域の周囲の近似変換成分が第１の変形と類似する結果になる。また、関心領域の周囲の近似変換成分が大きく変わらないのであれば、第１の変形と第２の変形の組み合わせとして、ＦＦＤとＴＰＳ、ＦＦＤとＬＤＤＭＭなど、性質の異なる変換モデルを適用してもよい。また、仮想対応点の対応点位置の誤差項Ｒ（Φ）を第２の変形の評価関数Ｅ２（Φ）の一部として用いるものであれば、何れの変換モデルや評価関数を用いてもよい。

＜第２の実施形態＞
第１の実施形態では、予め定められた近似変換モデルを用いて関心領域の第１の変形を近似していた。しかし、関心領域の性質に応じて適応的に近似変換モデルを選択できるようにしてもよい。本実施形態における画像処理装置は、関心領域に対応する部位（臓器や病変等）の種類に応じて適切な近似変換モデルを選択し、選択された近似変換モデルを用いて関心領域の第１の変形を近似する。これにより、予め定められた近似変換モデルが使用されることで、本来関心領域の変形とは異なる性質の近似変換モデルが適用されてしまい、位置合わせ精度が低下することを防ぐことができる。以下、本実施形態における画像処理装置について、第１の実施形態との相違部分について説明する。

図８は、本実施形態における画像処理装置８００の構成を示す図であり、新たに近似変換モデル選択部８０００が追加された点以外は、第１の実施形態における画像処理装置１００（図１）と同様である。また、関心領域取得部１０２０と変形近似部１０４０に関しては、第１の実施形態と機能が異なるため以下でその機能を説明する。その他の構成については第１の実施形態と機能が同じであるため、同じ機能に関しては説明を省略する。

関心領域取得部１０２０は、第１の画像から関心領域の取得と、関心領域に対応する部位の種類の分類（識別）を行う。近似変換モデル選択部８０００は、分類（識別）された部位の種類に応じて、当該部位の変形表現に適切な近似変換モデルを選択する。変形近似部１０４０は、関心領域における第１の変形変位場を、選択された近似変換モデルを用いて近似することで、近似変位場を生成する。

図９は、画像処理装置８００が行う全体の処理手順を示すフローチャートである。但し、本フローチャートにおいて、ステップＳ９０００、ステップＳ９０２０、およびステップＳ９０４０からＳ９０６０までの処理はそれぞれ、図２におけるステップＳ２０００、ステップＳ２０２０、およびステップＳ２０４０からＳ２０６０までの処理と同じであるため説明は省略する。以下、図２のフローチャートとの相違部分についてのみ説明する。

（Ｓ９０１０）（関心領域の取得と対応する部位の種類の分類）
ステップＳ９０１０において、関心領域取得部１０２０は、データ取得部１０１０から取得した画像上において、変形を抑制したい関心領域を表す情報を取得し、該関心領域に対応する部位の種類を分類（識別）する。そして、取得された関心領域とその分類情報を近似変換モデル選択部８０００へ出力する。

本ステップの具体的な処理を説明する。例えば、関心領域が臓器である場合を考える。関心領域取得部１０２０は、臓器の種類ごとの統計アトラスの当てはめ等の公知の画像解析処理によって、第１の画像から、臓器領域の抽出と、対応する臓器の種類の分類を行う。そして、関心領域取得部１０２０は、抽出・分類された臓器領域に基づいて、関心領域を取得する。関心領域として、抽出された臓器領域の少なくとも一つの領域が設定される。これは、ユーザによる画像上の領域の指定などの入力により設定されてもよいし、予め定められた変形を近似すべき臓器の種類（骨や肝臓など）に対応する領域が、抽出された臓器領域の中から設定されてもよい。なお、関心領域が複数存在する場合には、関心領域取得部１０２０は、これらの処理を夫々の関心領域に対して実行する。

なお、関心領域取得部１０２０は、データサーバに格納された被検体の臨床情報に記載された情報に基づいて、臓器領域の分類を行ってもよい。例えば、撮影された画像が腹部造影ＣＴ画像である場合は、画像内に写る臓器は、肝臓や脾臓の領域等であるため、関心領域取得部１０２０は、臓器の種類を絞りこむことができる。そして、例えば、関心領域取得部１０２０は、抽出された領域の体積が大きければ肝臓、小さければ脾臓のように分類できる。

また、関心領域が病変部である場合に、第１の画像に付随する情報として所見情報が記載されている場合に、関心領域取得部１０２０は、該所見情報に応じて病変部の分類情報を取得するようにしてもよい。例えば、「硬癌」「良性腫瘍」「ＤＣＩＳ」やそれに類する記載がある場合には、関心領域取得部１０２０は、その情報を部位の分類情報と設定することができる。また、それらの記載がない場合には、関心領域取得部１０２０は、分類情報を「不明」とすることができる。

（Ｓ９０１５）（近似変換モデルの選択）
ステップＳ９０１５において、近似変換モデル選択部８０００は、関心領域に対応する臓器の種類や性質に応じて、近似変換モデルの種類を選択する。例えば、近似変換モデル選択部８０００は、骨に分類された領域には剛体変換モデルを適用し、肝臓に分類された領域には体積保存のアフィン変換モデルを適用することができる。なお、関心領域が複数存在する場合には、近似変換モデル選択部８０００は、これらの処理を夫々の関心領域に対して実行する。また、近似変換モデル選択部８０００は、「硬癌」には剛体変換モデル、「良性腫瘍」には体積保存のアフィン変換モデル、「ＤＣＩＳ」や「不明」には「近似なし」を適用することができる。なお、近似変換モデル選択部８０００は、「近似なし」を適用した関心領域については、以降の近似処理を行う処理の対象から除外することができる。そして、近似変換モデル選択部８０００は、選択した近似変換モデルの種類を、対応情報生成部１０５０へ出力する。

（Ｓ９０３０）（関心領域における第１の変形を近似変換で近似）
ステップＳ９０３０において、変形近似部１０４０は、関心領域における第１の変形変位場を、この変位場よりも自由度の小さい、近似変換モデル選択部８０００から取得した近似変換モデルで近似した近似変位場を生成する。そして、変形近似部１０４０は、生成した近似変位場を、第２の対応情報生成部１０５０へと出力する。

本ステップの処理は、以下の点についてのみ第１の実施形態のステップＳ２０３０と異なる。すなわち、第１の実施形態のステップＳ２０３０では、変形近似部１０４０は、予め定められた近似変換モデルを用いていた。それに対し、本ステップでは、変形近似部１０４０は、近似変換モデル選択部８０００から取得した近似変換モデルを用いるようにしている。それ以外は、ステップＳ２０３０の処理と同様であるため、詳細な説明は省略する。なお、関心領域が複数存在する場合には、変形近似部１０４０は、これらの処理を夫々の関心領域に対して実行する。以上によって、画像処理装置８００の処理が実施される。

このように、本実施形態における画像処理装置は、関心領域に対応する臓器の種類に応じて適切な近似変換モデルを選択し、選択された近似変換モデルを用いて関心領域の第１の変形を近似する。これにより、予め定められた近似変換モデルが使用されることで、本来関心領域の変形とは異なる性質の近似変換モデルが適用されてしまい、位置合わせ精度が低下することを防ぐことができる。

＜第３の実施形態＞
第１の実施形態では、関心領域の第１の変形を自由度が小さい近似変換で近似して得た仮想的な対応情報を、そのまま第２の変形位置合わせの位置合わせ指標として利用していた。しかし、近似変換から得た仮想的な対応情報を直接、第２の変形位置合わせの位置合わせ指標として利用しなくてもよい。本実施形態では、第１の変形から得た仮想的な対応情報と、近似変換から得た仮想的な対応情報と、の間に位置する中間的な対応情報を、第２の変形位置合わせの位置合わせ指標として利用する。以下、本実施形態における画像処理装置について、第１の実施形態との相違部分について説明する。

図１０は、本実施形態における画像処理装置１０００の構成を示す図であり、新たに形状維持度取得部１００１０が追加された点以外は、第１の実施形態における画像処理装置１００と同様である。しかし、対応情報生成部１０５０、第２の変形取得部１０６０、表示制御部１０７０に関しては、第１の実施形態と機能が異なるため、以下でその機能を説明する。その他の構成については第１の実施形態と機能が同じであるため、同じ機能に関しては説明を省略する。

形状維持度取得部１００１０は、ユーザによる不図示の操作部からの入力により、変形位置合わせにおいて関心領域にどの程度形状を維持させるかを表す形状維持度の値を取得する。言い換えると、形状維持度とは、関心領域に関しての前記第１の変形変位場から前記第２の変形変位場への変化の程度を示す指標である。対応情報生成部１０５０は、第１実施形態と同様に変形仮想対応点と近似仮想対応点を生成する機能に加え、形状維持度取得部１００１０から取得した形状維持度に基づき変形仮想対応点と近似仮想対応点との間に位置する仮想対応点を生成する。これを中間仮想対応点と称する。第２の変形取得部１０６０は、第１実施形態において近似仮想対応点に基づいて行う変形位置合わせの処理を、近似仮想対応点に限定せずに、近似仮想対応点または中間仮想対応点の何れかの仮想対応点に基づくように置き換えた処理を行う。それ以外の機能は第１実施形態と同様であるため、詳細な説明は省略する。表示制御部１０７０は、第１実施形態と同様の機能に加え、ユーザによる不図示の操作部から処理を終了するか否かの命令を取得し、命令“終了する”を取得した場合には、処理は終了する。命令“終了しない”を取得した場合には、ユーザによる入力待ちの状態を維持する。

図１１は、画像処理装置１０００が行う全体の処理手順を示すフローチャートである。但し、本フローチャートにおいて、ステップＳ１１０００からＳ１１０２０、Ｓ１１０５０、Ｓ１１０６０、およびステップＳ１１１００の処理はそれぞれ、図２のフローチャートにおけるステップＳ２０００からＳ２０２０、Ｓ２０３０、Ｓ２０４０、およびステップＳ２０６０の処理と同じであるため説明は省略する。なお、本実施形態において、仮想対応情報は、第１の実施形態において述べたような近似仮想対応点であるとする。以下、図２のフローチャートとの相違部分について説明する。

（Ｓ１１０２５）（第１の変形画像の表示）
ステップＳ１１０２５において、表示制御部１０７０は、ユーザによる不図示の操作部への操作に応じて、取得した第１の変形画像と第２の画像の断面画像を、表示部１２０へと表示する制御を行う。このとき、表示制御部１０７０は、第１の変形画像と第２の画像の同一断面を、並べて表示部１２０へ表示するように制御することができる。

（Ｓ１１０３０）（形状維持度の取得）
ステップＳ１１０３０において、形状維持度取得部１００１０は、ユーザが不図示の操作部から入力した形状維持度を取得する。そして、取得した形状維持度の値を、対応情報生成部１０５０および表示制御部１０７０へと出力する。

以下では、形状維持度の値をλ_Ｒと表記する。本ステップで取得された形状維持度の値λ_Ｒに基づいて、ステップＳ１１０９０の第２の変形位置合わせが実施される。形状維持度の値λ_Ｒは、０≦λ_Ｒ≦１の範囲を満たす。λ_Ｒが０に近いほど、中間仮想対応点が変形仮想対応点の位置に近付き、関心領域の変形は周囲と変わらない変形に近付く。すなわち、関心領域は、周囲と変わらず特に拘束のない変形しやすい状態に近付く。一方、λ_Ｒが１に近いほど、中間仮想対応点の位置が近似仮想対応点の位置に近付き、関心領域の変形は近似変換に近付く。すなわち、関心領域は、周囲とは異なり変形しづらい状態に近付く。

形状維持度取得部１００１０は、λ_Ｒを、例えば、予め定められた複数の値の中からユーザが選択することで取得できる。例えば、予め設定された複数のλ_Ｒのリスト（例えば、｛０、０．２５、０．５、０．７５、１．０｝）が、表示部１２０上に表示されている場合を考える。この場合、ユーザは不図示の操作部を通じて、リストの中の何れかの値（例えば、０．５）を選択すると、形状維持度取得部１００１０は、その選択された値をλ_Ｒとして取得する。

なお、形状維持度の取得方法はこれに限られるものではなく、どのような方法であってもよい。例えば、ユーザが不図示の操作部を通じて直接値を入力し、形状維持度取得部１００１０は、その値を取得するようにしてもよい。この場合、例えばユーザが値１．０と入力すれば、形状維持度取得部１００１０は、λ_Ｒ＝１．０を取得する。

また、ユーザは、ステップＳ１１０３０において表示部１２０上に表示された第１の変形画像と第２の画像とを観察し、その結果に基づいて判断した形状維持度の値を入力することができる。例えば、表示部１２０上に表示された第１の変形画像上の関心領域が第２の画像上の対応する領域にある程度近い場合には、ユーザは、関心領域の変形が第１の変形の変換モデルである程度正確に表現できていると判断して、形状維持度を小さく設定できる。一方、形状が大きくかけ離れている場合には、ユーザは、関心領域の変形が第１の変形の変換モデルからかけ離れていると判断して、形状維持度を大きく設定できる。

前述したように、本実施形態における画像処理装置１０００は、形状維持度に応じて、関心領域の変形のしにくさを調整できる。例えば、近似変換モデルとして剛体変換を適用する場合には、画像処理装置１０００は、形状維持度によって関心領域の大よその硬さを調整することができる。そこで、例えば、予め作成された、定性的な硬さの「度合い」と形状維持度の値「λ_Ｒ」との対応テーブルにおいて、ユーザが定性的な硬さの度合いを選択してもよい。形状維持度取得部１００１０は、選択された「度合い」それに対応するλ_Ｒを取得するようにしてもよい。より具体的には、例えば、“柔らかい”＝０．０、“やや柔らかい”＝０．２５、“中間”＝０．２５、“やや硬い”＝０．７５、“硬い”＝１．０のようなマッピングを行う対応テーブルが予め画像処理装置１０００に保持されているとする。そして、表示制御部１０７０が、表示部１２０上に定性的な硬さの度合いのリスト（例えば、｛“柔らかい”、“やや柔らかい”、“中間”、“やや硬い”、“硬い”｝）を表示する。そして、ユーザが不図示の操作部を通じていずれかの度合い（例えば、“硬い”）を選択した場合に、形状維持度取得部１００１０は、対応テーブルを参照して、該度合いに対応する値（例えば、１．０）をλ_Ｒとして取得することができる。

なお、形状維持度λ_Ｒの取得方法は、ユーザによる手動入力による方法に限定されない。例えば、形状維持度取得部１００１０は、事前に画像処理装置８００内の不図示の記憶部に記憶された所定の値をλ_Ｒとして取得するようにしてもよい。また、形状維持度取得部１００１０は、データサーバ１１０に格納された被検体の臨床情報に基づいてλ_Ｒを取得するようにしてもよい。

例えば、読影レポート等の臨床情報に“病変部が硬い”あるいは“硬癌”、”石灰化”といった意味を示す所見や診断名が記載されている場合は、形状維持度取得部１００１０は、λ_Ｒ＝１．０を設定する。同様に、“病変部が周囲の組織と変わらず柔らかい”あるいは“非浸潤性癌”といった意味を示す所見や診断名が記載されている場合は、形状維持度取得部１００１０は、λ_Ｒ＝０．０を設定する。また、不図示の解析部が、第１の画像における関心領域に相当する領域の画像特徴を解析することで関心領域の定性的な硬さの度合いを判定し、それに基づいて、形状維持度取得部１００１０が、上述したマッピングによりλ_Ｒを取得するようにしてもよい。また、第２の実施形態において説明した図９のステップＳ９０１０と同様の方法で、関心領域取得部１０２０が第１の画像から臓器領域の抽出を行い、形状維持度取得部１００１０が、夫々の臓器に予め定めた形状維持度の値を夫々の領域に設定するようにしてもよい。例えば、形状維持度取得部１００１０は、骨と分類された領域にはλ_Ｒ＝１．０を、肝臓や肺と分類された領域にはλ_Ｒ＝０．７５を設定するといった処理を行うことができる。

（Ｓ１１０４０）（形状維持度が０か否かの判定）
ステップＳ１１０４０において、形状維持度取得部１００１０は、形状維持度の値λ_Ｒが０か否かを判定する。そして、λ_Ｒ≠０の場合は、処理はステップＳ１１０５０へ進む。一方、λ_Ｒ＝０の場合は、処理はステップＳ１１１００へ進む。これにより、形状維持度λ_Ｒが０以外の場合は、画像処理装置１０００は、ステップＳ１１０５０からＳ１１０９０の処理により、形状維持度の値λ_Ｒに基づいて関心領域の形状が維持された第２の変形画像を最終的な位置合わせ結果として取得できる。一方、λ_Ｒが０の場合は関心領域の形状を維持する必要がないため、画像処理装置１０００は、ステップＳ１１０２０で生成された第１の変形画像をそのまま最終的な位置合わせ結果として取得できる。これにより、不要に第２の画像を生成する処理を削減でき、全体の処理フローを効率化することができる。

（Ｓ１１０７０）（形状維持度が１か否かの判定）
ステップＳ１１０７０において、形状維持度取得部１００１０は、形状維持度の値λ_Ｒが１か否かを判定する。そして、λ_Ｒ≠１の場合は、処理はステップＳ１１０８０へ進む。一方、λ_Ｒ＝１の場合は、処理はステップＳ１１０９０へ進む。これにより、λ_Ｒが１以外の場合は、画像処理装置１０００は、ステップＳ１１０８０の処理によって中間仮想対応点を生成し、ステップＳ１１０９０においてそれを仮想対応点位置の誤差項に適用して位置合わせした第２の変形画像を最終的な位置合わせ結果として取得できる。一方、λ_Ｒが１の場合は、画像処理装置１０００は、ステップＳ１１０６０で生成した近似仮想対応点を、そのままステップＳ１１０９０における仮想対応点位置の誤差項に適用して位置合わせした第２の変形画像を最終的な位置合わせ結果として取得できる。従って、不要に中間仮想対応点を生成する処理を削減でき、全体の処理フローを効率化することができる。

（Ｓ１１０８０）（形状維持度に基づき中間仮想対応点を生成）
ステップＳ１１０８０において、対応情報生成部１０５０は、形状維持度取得部１００１０から取得した形状維持度の値λ_Ｒに基づき、関心領域Ｄ_ＲＯＩ上の変形仮想対応点と近似仮想対応点の間に位置する中間仮想対応点を生成する。そして、対応情報生成部１０５０は、生成した中間仮想対応点を第２の変形取得部１０６０へと出力する。

中間仮想対応点の位置は、λ_Ｒが０に近い場合は変形仮想対応点に近付き、λ_Ｒが１に近い場合は近似仮想対応点に近付く。近似仮想対応点は、関心領域の変形を近似変換モデルに従うように拘束するの。それに対し、中間仮想対応点は、関心領域の変形を周囲の変形と変わらない状態と近似変換モデルに従う状態の間の状態に拘束する特徴を持つ。以下に中間仮想対応点の具体的な生成方法を説明する。

まず、対応情報生成部１０５０は、変形仮想対応点と近似仮想対応点とを夫々構成する、第１の画像上の代表点を夫々変換した先の対応する代表点、つまり変形代表点と近似代表点との間に位置する代表点を生成する。これを、中間代表点と称する。次に、対応情報生成部１０５０は、代表点とそれに対応する中間代表点から、中間仮想対応点として、画像間の仮想的な対応点を生成する。

中間代表点の具体的な生成方法として、本実施形態では以下の方法を採用する。変形代表点と近似代表点におけるｉ番目の点を夫々、Ｐ_Ｄｉ、Ｐ_Ｒｉと定義する。点Ｐ_Ｄｉ、Ｐ_Ｒｉの座標は、前述の通り夫々、ｘ_Ｄ２ｉ、ｘ_Ｒ２ｉで表される。ｉは、１≦ｉ≦Ｎ_Ｄ＝Ｎ_Ｒを満たす。このとき、中間代表点のｉ番目の点Ｐ_Ｍｉを、Ｐ_ＤｉとＰ_Ｒｉの間を結んだ線分Ｐ_ＤｉＰ_Ｒｉ上に位置し、かつ点Ｐ_Ｄｉからの距離が線分Ｐ_ＤｉＰ_Ｒｉの距離に対して比率がλ_Ｒとなる位置として算出する。つまり、点Ｐ_Ｍｉの座標をｘ_Ｍ２ｉとすると、ｘ_Ｍ２ｉは以下の式で表される。

これにより、中間代表点として、変形代表点と近似代表点の間に位置する点を取得することができる。

図１２は、中間代表点１２０００を示す図である。図に示すように、中間代表点１２０００は、元の代表点の分布が斜めにつぶれたような分布の変形代表点４０４０と、つぶれていない分布の剛体代表点（近似代表点）６０００の間の分布となっていることが分かる。

（Ｓ１１０９０）（取得した仮想対応点に基づき第２の変形を取得）
ステップＳ１１０９０において、第２の変形取得部１０６０は、取得した仮想対応点に基づき、データ取得部１０１０から取得した第１の画像を第２の画像に対して変形位置合わせする。第１実施形態のステップＳ２０５０で用いる仮想対応点が近似仮想対応点に限定されていたのに対し、本実施形態では、近似仮想対応点と中間仮想対応点の何れかの仮想対応点を用いる点が異なる。すなわち、λ_Ｒ＝１の場合はステップＳ１１０６０で生成した近似仮想対応点を用い、λ_Ｒ≠１の場合はステップＳ１１０８０で生成した近似仮想対応点を用いる。但し、取得された仮想対応点に基づく変形位置合わせ処理に関しては、ステップＳ２０５０と同様であるため、詳細な説明は省略する。

（Ｓ１１１１０）（終了判定）
ステップＳ１１１１０において、表示制御部１０７０は、ユーザによる不図示の操作部から処理を終了するか否かの命令を取得する。そして、命令“終了する”を取得した場合には、処理は終了、命令“終了しない”を取得した場合には、処理はステップＳ１１０３０へ進み、ユーザによる入力待ちの状態が維持される。これにより、ユーザがステップＳ１１１００において位置合わせ結果として第１の変形画像または第２の変形画像を観察した結果、位置合わせを完了したと判断した場合には、不図示の操作部を通じて命令“終了する”を与えることで、画像処理装置１０００は処理を終了できる。一方、ユーザが位置合わせ結果の画像を観察した結果、形状維持度を変更して再度結果を観察したいと考えた場合には、不図示の操作部を通じて命令“終了しない”を与えることで、画像処理装置１０００は、再度、新たな形状維持度を与えて位置合わせ結果を生成し表示させることができる。

以上によって、画像処理装置８００の処理が実施される。本処理フローにより、例えば、ステップＳ１１０２０およびステップＳ１１０９０において、線形最適化で変換パラメータが算出できる場合は、高速に夫々の位置合わせが可能である。そのため、ユーザがインタラクティブに形状維持度を画像処理装置１０００に入力した場合に、該入力に対応する第２の変形画像の表示が可能となる。より具体的には、不図示の操作部が、形状維持度を０から１まで調整可能なスライダーであるとする。このとき、ユーザがスライダーを動かした場合は、命令“終了しない”を画像処理装置１０００に与えるようにすることで、画像処理装置１０００は、スライダー位置の形状維持度に基づく変形画像を生成し表示できる。これにより、例えば、ユーザがスライダー位置を０から１に向けて徐々に動かしていくと、画像処理装置１０００は、それに連動して徐々に形状が維持された第２の変形画像をほぼリアルタイムに表示部１２０に表示できる。従って、ユーザは表示部１２０に表示された第２の画像と第２の変形画像を観察しながら、適切な形状維持度に基づく第２の変形画像を容易に決めることができる。なお、ユーザがインタラクティブに形状維持度の入力を行うための不図示の操作部のインターフェースは、スライダーに限られるものではなく、例えば直接数値を入力するようなものであってもよい。

本実施形態によれば、画像処理装置１０００は、関心領域の変形仮想的対応点と近似仮想対応点との間に位置する中間仮想対応点を、画像間の位置合わせ指標に利用する。これにより、関心領域の変形に関して、第１の変形と近似変換の間の変形状態となるような位置合わせ結果が得られる。そのため、周囲より硬いが変形は生じるような病変部を含む画像であっても、被検体全体を適切に位置合わせできる。

（変形例３−１）（中間仮想対応点を非線形の補間で生成）
第３の実施形態では、ステップＳ１１０８０において、中間仮想対応点を構成する中間代表点を、変形代表点と近似代表点とを結ぶ線分上の点としていた。言い換えると、変形代表点と近似代表点の座標間を線形補間した座標としていた。しかし、中間代表点の生成方法はこれに限られるものではなく、変形代表点と近似代表点の座標間を非線形に補間した座標としてもよい。すなわち、関心領域に関しての第１の変形変位場から第２の変形変位場への変化は線形であっても、非線形であってもよい。

例えば、以下の方法を採ってもよい。まず、対応情報生成部１０５０は、関心領域の変形代表点と近似体表点の夫々から、関心領域の表面上に設定されたもののみを抽出する。それらを、それぞれ変形表面代表点、近似表面代表点と称する。次に、対応情報生成部１０５０は、変形表面代表点と近似表面代表点の間に位置する中間的な表面代表点を、ステップＳ１１０８０と同様に形状維持度の値λ_Ｒに基づいて線形補間によって生成する。これを中間表面代表点と称する。そして、対応情報生成部１０５０は、近似表面代表点と中間表面代表点の位置を境界条件とし、近似表面代表点から中間表面代表点へ向かう方向へ、関心領域内部の変形を物理変形シミュレーションによって推定する。

ここで、近似表面代表点を基準としてシミュレーションを行っているのは、関心領域が近似変換モデルに従っている状態を、仮想的に外力が生じていない状態であるとしたためである。そして、変形表面代表点に対応する関心領域が最も外力が生じている状態であり、中間表面代表点に対応する関心領域がその中間的な外力が生じている状態であると仮定している。また、シミュレーションに必要な関心領域のヤング率・ポアソン比の値は、予め所定の値を設定する。例えば、画像処理装置１０００は、関心領域の種類ごとに予め対応するヤング率・ポアソン比の値のテーブルを用意しておき、対応情報生成部１０５０は、当該症例の関心領域の種類を臨床データ等から取得することで、テーブルから対応するヤング率・ポアソン比の値を取得し、設定する。関心領域の種類の情報は、データサーバ１１０から、データ取得部１０１０、関心領域取得部１０２０を通じて、対応情報生成部１０５０が取得するものとする。そして、対応情報生成部１０５０は、シミュレーションにより得られた変位場を用いて、近似代表点を変換することで、中間表面代表点に対応する、関心領域内部の代表点が取得できる。これを中間内部代表点と称する。そして、対応情報生成部１０５０は、中間表面代表点と中間内部代表点を合わせると、関心領域全体に関する中間的な変換後の代表点が得られるため、これを中間代表点として取得する。

これにより、関心領域の表面に関しては近似代表点と変形代表点との座標を線形補間した座標となる。一方、内部領域に関しては、近似代表点と変形代表点との座標を非線形に（より、物理的に整合性のとれた位置に）補間した座標として取得できる。なお、例えば、関心領域内部のヤング率を領域によって変えることで、領域によって異なった非線形な補間を行うことができる。例えば、関心領域の重心に近いほどヤング率を大きくし、関心領域の表面に近いほどヤング率を小さくすることで、関心領域の重心近くは硬くなり変位量が小さくなる。一方、関心領域の表面近くは柔らかくなり変位量が大きくなる。これにより、中心部が硬く周辺部が柔らかい病変などに関して、現実に近い変形を取得することができる。

（変形例３−２）（一度に複数の形状維持度を取得し複数の第２の画像を生成）
第３の実施形態では、ステップＳ１１０３０において、形状維持度取得部１００１０は、形状維持度を１つだけ取得した。そして、ステップＳ１１０４０からＳ１１０９０を経て、その値に対応する変形画像を１つだけ生成されていた。また、必要であれば、形状維持度を再取得しそれに対応する変形画像を取得する、という処理を繰り返していた。しかしながら、形状維持度取得部１００１０は、形状維持度の値を１度に複数取得し、それに対応する変形画像を１度に複数生成されるようにしてもよい。

例えば、ステップＳ１１０００からＳ１１０２５の後のステップＳ１１０３０において、形状維持度取得部１００１０は、予め定めた複数の形状維持度の値λ_Ｒを取得する。例えば、形状維持度取得部１００１０は、｛０、０．２５、０．５、０．７５、１．０｝の値すべてをλ_Ｒとして取得する。このとき、画像処理装置１０００は、ステップＳ１１０４０及びステップＳ１１０７０の処理は行わない。次に、変形近似部１０４０と対応情報生成部１０５０は、ステップＳ１１０５０およびＳ１１０６０の処理を行い、近似仮想対応点を生成する。次に、ステップＳ１１０８０において、対応情報生成部１０５０は、λ_Ｒ≠０かつλ_Ｒ≠１である全ての形状維持度の夫々に対応する中間仮想対応点を生成する。前記の例では、対応情報生成部１０５０は、｛０．２５、０．５、０．７５｝の夫々に対応する中間仮想対応点３種類を生成する。そして、ステップＳ１１０９０において、第２の変形取得部１０６０は、ステップＳ１１０６０で生成した近似仮想対応点と、ステップＳ１１０８０で生成した全ての中間仮想対応点の夫々に基づいて、対応する第２の変形画像を生成する。前記の例では、中間仮想対応点３種類（λ_Ｒ＝０．２５、０．５、０．７５に対応）と、近似仮想対応点（λ_Ｒ＝１に対応）の、合計４種類に対応する第２の変形画像が生成される。

そして、ステップＳ１１１００において、表示制御部１０７０は、生成された全ての変形画像を表示する。前記の例では、第１の変形画像（λ_Ｒ＝０に対応）と、第２の変形画像４種類（λ_Ｒ＝０．２５、０．５、０．７５、１．０に対応）の、合計５種類の変形画像が表示される。表示方法として、例えば、全ての変形画像を並べて表示してもよいし、ユーザの不図示の操作部からの入力に応じて、同じ表示領域上で切り替えて表示するようにしてもよい。

本変形例によると、画像処理装置１０００は、予めステップＳ１１０００からＳ１１０９０までの処理を、ユーザの入力を介さずオフラインで処理することができる。この場合、画像処理装置１０００は、生成した全ての変形画像を不図示の記憶部に保存しておく。そして、ステップＳ１１１００において、記憶部から変形画像を読み込んで表示するようにする。第３の実施形態では、異なる形状維持度の結果を表示するためには、その度にステップＳ１１０４０からＳ１１０９０までの処理を行う必要があり、ユーザに待ち時間を与えてしまっていた。しかし、本変形例の処理は、ステップＳ１１１００では予め生成された複数の変形画像を読み込んで表示するだけである。そのため、ユーザは待ち時間なく変形画像を観察でき、効率的に作業を行うことができる。

＜第４の実施形態＞
第１の実施形態では、位置合わせの評価関数の例として、ランドマーク対応点位置の誤差を採用していた。しかし、評価関数には、予め１対１に対応する点として取得されない情報の誤差を含むようにしてもよい。本実施形態における画像処理装置は、夫々の画像から抽出した被検体の表面形状の情報に基づいて、ある変換パラメータに対する画像間での表面形状の一致度を求め、これを評価関数の誤差項として併用する。これにより、表面形状が一致するような変換パラメータを算出することができる。以下、本実施形態における画像処理装置について、第１の実施形態との相違部分について説明する。

図１３は、本実施形態における画像処理装置１３００の構成を示す図であり、新たに表面形状取得部１３０１０が追加された点以外は、第１の実施形態における画像処理装置１００（図１）と同様である。また、第１の変形取得部１０３０及び第２の変形取得部１０６０に関しては、第１の実施形態と機能が異なるため以下でその機能を説明する。その他の構成については第１の実施形態と機能が同じであるため、同じ機能に関しては説明を省略する。

表面形状取得部１３０１０は、第１の画像および第２の画像から被検体の表面形状を夫々取得する。第１の変形取得部１０３０及び第２の変形取得部１０６０の機能は第１の実施形態とほぼ同じであるが、取得された夫々の表面形状にさらに基づいて変形位置合わせを行う点のみが、第１の実施形態と異なっている。

図１４は、画像処理装置１３００が行う全体の処理手順を示すフローチャートである。但し、本フローチャートにおいて、ステップＳ１４０００、Ｓ１４０１０、Ｓ１４０３０、Ｓ１４０４０、Ｓ１４０６０の処理はそれぞれ、図２におけるステップＳ２０００、Ｓ２０１０、Ｓ２０３０、Ｓ２０４０、Ｓ２０６０の処理と同じであるため説明は省略する。以下、図２のフローチャートとの相違部分についてのみ説明する。

（Ｓ１４０１５）（表面形状の取得）
ステップＳ１４０１５において、表面形状取得部１３０１０は、データ取得部１０１０から取得した第１の画像と第２の画像の夫々において、被検体の表面形状を表す情報を取得する。ここで、被検体の表面とは、例えば、被検体が乳房の場合には、体表や大胸筋面である。また、被検体が肝臓や腎臓などの臓器の場合には、表面とは当該臓器の表面のことである。ここで、表面形状は、例えば画像に対して面強調フィルタやエッジ検出等の画像処理を行うことで自動的に取得するようにできる。また、画像をユーザが観察できるようにして、ユーザによる入力操作等に基づいて表面形状を取得するようにしてもよい。取得した第１の画像および第２の画像の表面形状を夫々、第１の表面形状および第２の表面形状と表記する。本実施例では表面形状は点群により構成される。但し、表面形状は必ずしも点群である必要はなく、形状を表現できる形式であれば何であってもよい。例えば、表面形状は、多項式などの曲面を表現可能な数式で近似された関数（以降、曲面関数と称する）であってもよい。そして、取得された第１の表面形状および第２の表面形状を第１の変形取得部１０３０及び第２の変形取得部１０６０へと出力する。

（Ｓ１４０２０）（表面形状に基づく第１の変形の取得）
ステップＳ１４０２０において、第１の変形取得部１０３０は、取得した第１の表面形状および第２の表面形状に基づき、取得した第１の画像と第２の画像との間の第１の変形位置合わせを実行し、第１の変形変位場を取得する。そして、第１の変形取得部１０３０は、取得した第１の変形変位場を対応情報生成部１０５０へと出力する。さらに、第１の変形取得部１０３０は、第１の変形変位場に基づいて第１の画像を変形させた第１の変形画像を生成して、表示制御部１０７０へと出力する。

このとき、第１の変形取得部１０３０は、第１の実施形態のステップＳ２０２０と同様、変形Φによる位置合わせの適切さを評価する評価関数Ｅ１（Φ）を定義し、該評価関数を最小化するような変形Φを推定する。本ステップとステップＳ２０２０との処理の違いは、評価関数Ｅ１（Φ）の計算方法のみであり、それ以外の処理はステップＳ２０２０と同様であるため、説明を省略する。ステップＳ２０２０では、評価関数Ｅ１（Φ）において、実対応点位置の誤差はランドマーク対応点位置の誤差で表わされていたのに対し、本ステップでは実対応点位置の誤差は、ランドマーク対応点位置の誤差に加え、画像間の表面形状の誤差も含めた関数で表わされる点で異なる。ここで、表面形状の誤差は、画像間の表面形状を表す点群同士が対応付けられた対応点（表面対応点）の位置の誤差として表現される。

以下に、具体的な算出方法を説明する。表面対応点位置の誤差項は、夫々点群によって表現されている第１の表面形状と第２の表面形状の間で、１対１に対応付けられた表面対応点を生成し、生成された表面対応点の誤差として算出する。しかし、第１の表面形状と第２の表面形状は夫々を構成する点群同士が最初から対応付けられているわけではない。そのため、第１の表面形状を構成する点群により表現される表面形状の曲面と、第２の表面形状を構成する点群により表現される表面形状の曲面との間の位置合わせを行い、点群間の対応付けを行う必要がある。この処理は、例えば、Iterative Closest Point(ICP)法などを公知の手法を用いて実行することができる。また、画像間で元々の表面形状が大きく異なる場合には、まずステップＳ２０２０と同様の処理によってランドマーク対応点の誤差のみを用いた変形変位場を算出し、この変形変位場を第１の表面形状に適用することで、一旦、第２の表面形状に形状を近付けた変形後の第１の表面形状を取得する。そして、変形後の第１の表面形状を表す点群と第２の表面形状を表す点群の間でICP法などで対応付けるようにしてもよい。これにより、画像間で元々の表面形状が大きく異なる場合でも、第１の表面形状と第２の表面形状との間で適切に点群を対応付けることができる。なお、表面形状の対応付け方法はこれに限られるものではない。例えば、一方の表面形状を点群として、他方の表面形状を曲面関数として取得し、これらの点群と曲面関数との間で、一方の点群に対応する他方の曲面関数上の最近傍点を探索しながらICP法などで点群の対応付けを行ってもよい。これにより、一方の表面形状を表す点群は、他方の連続的な形式で表現された表面形状上の点と対応付けられることになる。そのため、両方の画像の表面形状を点群として離散化した場合に生じる、一方の表面形状上の点に厳密に対応する他方の表面形状上の点が存在せず、誤差が発生してしまうという問題が起こらない。

次に、表面対応点の誤差項Ｓｕｒｆ（Φ）は次式で表わされる。

ここで、（９）式では、Φ（ｘ｜ｐ１）は、パラメータをｐ１とした所定の変換モデルで表される第１の変形Φによって座標ｘを変位させた座標を表す。また（ｘ_Ｓ１ｉ，ｘ_Ｓ２ｉ）は、第１の画像と第２の画像上におけるｉ番目の表面対応点の座標の組を表す。また、ＮＳは、表面対応点の総数を表す。また、Ｃｏｖ（ｘ）は、第１の実施形態の変形例１−２と同様、位置ｘにおける３次元の誤差に関する共分散行列を返す関数を表す。（９）式は、（６）式における近似仮想対応点（ｘ_Ｒ１ｉ，ｘ_Ｒ２ｉ）を、表面対応点（ｘ_Ｓ１ｉ，ｘ_Ｓ２ｉ）に、変換パラメータｐ２を変換パラメータｐ１に置き換えただけの式である。従って、変形例１−２においてこのようにパラメータを置き換えると、（９）式は以下のように解釈できる。すなわち、Ｓｕｒｆ（Φ）は、変換パラメータｐ１による変位後の第１の表面形状の表面対応点が、対応する第２の表面形状の表面対応点から法線方向にのみずれないようにするための項として機能する。つまり、表面形状上に設定された表面対応点は、第１の変位場の算出時に、表面形状の面上に沿って動くことを許容される。なお、評価関数Ｅ１（Φ）、つまり実対応点の誤差は、上述の通り、ランドマーク対応点の誤差項に表面対応点の誤差項Ｓｕｒｆ（Φ）を加えることで算出できる。実対応点の誤差は表面対応点の誤差項Ｓｕｒｆ（Φ）のみで構成されてもよい。

このように、評価関数Ｅ１（Φ）は対応点位置の誤差項のみによって表されるため、第１の実施形態のステップＳ２０５０に記載した通り、変換パラメータと対応点間の誤差（誤差）との関係は、線形方程式で表すことができる。従って、対応点間の誤差を最小化する変換パラメータｐ１を線形最適化手法で推定できる。

（Ｓ１４０５０）（表面形状に基づく第２の変形の取得）
ステップＳ１４０５０において、第２の変形取得部１０６０は、ステップＳ２０４０で取得した近似仮想対応点（仮想対応情報）を拘束条件として用いて第２の変形変位場を取得する。このとき、第２の評価関数Ｅ２（Φ）の算出において、ステップＳ１４０２０の処理と同様に、ステップＳ１４０１５で取得した第１の表面形状および第２の表面形状の一致度を考慮する点が第１の実施形態とは異なっている。ここで、第２の評価関数Ｅ２（Φ）の算出において、第１の評価関数Ｅ１（Φ）に追加される近似仮想対応点の誤差項Ｒ（Φ）も対応点位置の誤差項である。従って、最終的な位置合わせの評価関数Ｅ２（Φ）も対応点位置の誤差項のみで表わされるため、変換パラメータと対応点間の誤差との関係を線形方程式で表わすことができる。従って、対応点間の誤差を最小化する変換パラメータｐ２を線形最適化手法によって推定できる。

以上によって、画像処理装置１３００の処理が実施される。本実施形態によれば、被検体の表面形状が一致するような位置合わせを行うことができる。このとき、予め対応点として取得されていない画像間の表面形状の一致度（誤差）を評価関数に含む場合でも、線形最適化によって変換パラメータを算出可能となり、高速に位置合わせの結果を取得することができる。さらに、本実施形態によれば、位置合わせにおいて、被検体の表面形状上に設定された表面対応点は、表面形状の面上に沿って動くことを許容される。これにより、以下のような効果が得られる。本実施形態では、表面形状間をICP法などの方法で対応付けを行っているが、アルゴリズムの不安定さや精度面での不十分さなどにより各画像の表面形状の面上で解剖学的に厳密に正しい位置が常に対応付けられるわけではない。そのため、多少なりとも対応付けの誤りが生じる。従って、表面形状の面上の位置の対応の誤りにより、表面形状上の近傍に存在する対応点間の距離が画像間で大きく食い違ってしまうことがある。これにより、単純に対応点の位置そのものをなるべく一致させるような評価関数で位置合わせすると、表面形状が局所的に不自然に伸び縮みしてしまう。これに対して、表面形状の面上に動くことを許容する評価関数で位置合わせすることで、変形を表現する変換モデルの性質に基づいて、表面形状の面上に沿って誤差を許した状態で位置合わせが行われる。これにより、本実施形態のように滑らかな変形を表現することが知られているＦＦＤやＴＰＳのような変換モデルを適用することで、被検体の表面形状が局所的に不自然な変形にならないような位置合わせ結果を得ることができる。

（変形例４−１）
上記の実施形態では、第１の実施形態から第３の実施形態と同様に、関心領域の変形を制御するための第２の変形を取得していた。しかし、関心領域の変形を制御しない（制御する必要がない）場合であっても、被検体の表面形状を一致させるための処理を位置合わせに用いる上記処理の効果は損なわれるものではない。すなわち、ステップＳ１４０１０、および、Ｓ１４０３０からステップＳ１４０５０の処理は実行せずに、第１の変形を変形結果として用いるようにしてもよい。

＜第５の実施形態＞
上述の第４の実施形態の変形例４−１では、被検体の表面形状に基づき、画像間の表面形状の一致度を、位置合わせの評価関数の誤差項として用いる場合について説明した。しかし、これらの実施形態は本発明の実施のおける一例に過ぎない。本実施形態では変形例４−１とは異なる形態の例について説明する。

図１５は、本実施形態における画像処理装置１５００の構成を示す図である。本実施形態における画像処理装置１５００は、第４の実施形態で説明した画像処理装置１３００の構成要素の一部と同じ構成要素で構成されている。同様の機能を有する構成要素には同一の番号を付し、詳細な説明は省略する。

図１６は、本実施形態における画像処理装置１５００が行う全体の処理手順を示すフローチャートである。

（Ｓ１５０００）データの取得
ステップＳ１５０００において、データ取得部１０１０は、データサーバ１１０から、第１の画像と第２の画像を取得する。この処理は第４の実施形態におけるステップＳ１４０００と同様であるため説明を省略する。

（Ｓ１５０１５）表面形状の取得
ステップＳ１５０１５において、表面形状取得部１３０１０は、データ取得部１０１０から取得した第１の画像と第２の画像の夫々において、被検体の表面形状を表す情報を取得する。この処理は第４の実施形態におけるステップＳ１４０１５と同様であるため説明を省略する。

（Ｓ１５０２０）表面形状に基づく第１の変形の取得
ステップＳ１５０２０において、第１の変形取得部１０３０は、取得した第１の表面形状および第２の表面形状に基づき、取得した第１の画像と第２の画像との間の第１の変形位置合わせを実行し、第１の変形変位場を取得する。さらに、第１の変形取得部１０３０は、第１の変形変位場に基づいて第１の画像を変形させた第１の変形画像を生成して、表示制御部１０７０へと出力する。

図１７は、本処理ステップにおいて第１の変形取得部１０３０が実行する処理をより詳細に説明するためのフローチャートである。

（Ｓ１５０２２）対応関係取得
ステップＳ１５０２２において、第１の変形取得部１０３０は、第１の表面形状および第２の表面形状を構成する点群の間の対応付け処理を実行する。具体的な処理方法は、第４の実施形態におけるステップＳ１４０２０に一例として記載したICP法などを用いることができる。この処理により、表面対応点（ｘ_Ｓ１ｉ，ｘ_Ｓ２ｉ）、1≦ｉ≦Ｎ_Ｓを取得する。ここで、Ｎ_Ｓは表面対応点の総数である。

（Ｓ１５０２４）法線方向取得
ステップＳ１５０２４において、第１の変形取得部１０３０は、ステップＳ１５０２２で取得した表面対応点の夫々について、当該位置における表面形状に関する法線方向の算出処理を実行する。ここで、表面対応点（ｘ_Ｓ１ｉ，ｘ_Ｓ２ｉ）を構成するｘ_Ｓ１ｉおよびｘ_Ｓ２ｉは、夫々、第１の表面形状および第２の表面形状を構成する点である。本処理ステップでは、位置ｘ_Ｓ１ｉにおける第１の表面形状の法線方向を表すベクトル（第１の法線ベクトル）ｎ_１ｉおよび、位置ｘ_Ｓ２ｉにおける第２の表面形状の法線方向を表すベクトル（第２の法線ベクトル）ｎ_２ｉを夫々算出して取得する。表面形状から法線ベクトルを算出する方法は、変形例１−２に記載した方法（当該位置の近傍の点群の位置を主成分分析して得る方法）で算出できる。法線ベクトルの算出方法は、他にも表面形状の距離場を算出し、その距離場の空間勾配に基づいて法線ベクトルを算出するようにしても良い。以上の方法により得た第１の表面形状の位置ｘ_Ｓ１ｉにおける法線ベクトルｎ_１ｉおよび、第２の表面形状の位置ｘ_Ｓ２ｉにおける法線ベクトルｎ_２ｉを得る。なお、ｎ_１ｉおよびｎ_２ｉは夫々、ノルムが1.0に正規化された３次元ベクトルである。

（Ｓ１５０２６）変形算出
ステップＳ１５０２６において第１の変形取得部１０３０は、ステップＳ１５０２２およびステップＳ１５０２４の処理結果に基づいて、第１の変形変位場を取得する。

このとき、第１の変形取得部１０３０は、第４の実施形態のステップＳ１４０２０と同様、変形Φによる位置合わせの適切さを評価する評価関数Ｅ１（Φ）を定義し、該評価関数を最小化するような変形Φを推定する。

（１０）式の、Φ（ｘ｜ｐ１）は、パラメータをｐ１とした所定の変換モデルで表される第１の変形Φによって座標ｘを変位させた座標を表す。本実施形態において、Φ（ｘ｜ｐ１）は、第１の画像の座標系から第２の画像の座標系への変換を意味する。また、Φ^−１（ｘ｜ｐ１）はΦ（ｘ｜ｐ１）の逆関数であり、第２の画像の座標系から第１の画像の座標系への変換を意味する。また、Ｃｏｖ_１ｉ、Ｃｏｖ_２ｉは、ｉ番目の表面対応点の位置における３次元の誤差に関する共分散行列である。Ｃｏｖ_１ｉは、第１の表面形状の位置ｘ_Ｓ１ｉにおける法線ベクトルｎ_１ｉに基づいて算出される。また、Ｃｏｖ_２ｉは、第２の表面形状の位置ｘ_Ｓ２ｉにおける法線ベクトルｎ_２ｉに基づいて算出される。

（Ｓ１５０２８）変形画像生成
ステップＳ１５０２８において第１の変形取得部１０３０は、ステップＳ１５０２６で取得した第１の変形変位場に基づいて第１の画像を変形させた第１の変形画像を生成して、表示制御部１０７０へと出力する。

以上に説明したステップＳ１５０２２からステップＳ１５０２８により、本実施形態におけるステップＳ１５０２０の処理が実行される。

（Ｓ１５０６０）変形画像の表示
ステップＳ１５０６０において、表示制御部１０８７は、ユーザによる不図示の操作部への操作に応じて、取得した変形画像と第２の画像の断面画像を、表示部１２０へと表示する制御を行う。この処理は第４の実施形態におけるステップＳ１４０６０と同様であるため詳細な説明を省略する。

以上によって、第６の実施形態における画像処理装置１５００の処理が実施される。本実施形態によれば、第１の表面形状に関する法線ベクトルおよび、第２の表面形状に関する法線ベクトルの双方に基づいて第１の変形変位場を取得することができる。これにより、第１の表面形状と第２の表面形状の違い等に対して安定的に動作させることができる効果がある。

（変形例５−１）：法線方向の算出方法（画像の輝度勾配も用いる）
上記の実施形態では、第１の表面形状に基づいて第１の法線ベクトルを算出する場合を例として説明した。しかし、本発明の実施はこれに限らない。例えば、第１の法線ベクトルｎ_１ｉは、第１の画像の表面対応点位置ｘ_Ｓ１ｉの近傍の輝度勾配に基づいて算出するようにしても良い。同様に、第２の法線ベクトルｎ_２ｉは第２の画像の表面対応点位置ｘ_Ｓ２ｉの近傍の輝度勾配に基づいて算出するようにしても良い。表面形状を表す点群の分布が空間的に疎な場合など、表面形状に関する情報から法線ベクトルを高精度に算出するのが困難な場合であっても以降の処理を実行することができる。

また、法線ベクトルの算出は、表面形状と画像の輝度勾配の両方に基づいて算出するようにしても良い。例えば、表面形状に基づいて算出した法線ベクトルと、画像の輝度勾配に基づいて算出した法線ベクトルとの中間ベクトルを算出し、そのベクトルを用いて以後の処理を実行するようにしても良い。また、表面形状の粗密度合いに基づいて、表面形状から法線ベクトルを算出するか、画像の輝度勾配から法線ベクトルを算出するかを切り替えるようにしても良い。例えば、表面形状を表す点群の空間密度が所定の閾値よりも大きい場合には表面形状に基づいて法線ベクトルを算出し、そうでない場合には画像の輝度勾配に基づいて法線ベクトルを算出するようにできる。また、画像の輝度勾配の大きさによって上記の切り替えを行うようにしても良い。これらの方法によれば、表面形状を表す点群データや画像の輝度分布に対して、より安定的に法線ベクトルを算出できる効果がある。

（変形例５−２）：２つの法線方向を統合して用いる場合
上記の実施形態では、第１の法線ベクトルおよび第２の法線ベクトルを算出し、それらに基づいて式（１０）に記載した評価関数を用いる場合を例として説明した。しかし、本発明の実施はこれに限らない。例えば、第１の法線ベクトルおよび第２の法線ベクトルの双方に基づいて、これらを統合した法線ベクトルを算出するようにしても良い。そして、統合した法線ベクトルに基づいて第４の実施形態で説明した式（９）に記載した評価関数を用いて第１の変形変換を取得するようにしても良い。具体的には、ステップＳ１５０１５で取得した第１の画像の座標系で定義される第１の法線ベクトルを座標変換Φ（ｘ|ｐ１）に基づいて、第２の画像座標系の法線ベクトルへと変換する。より具体的にはｎ_１ｉ’＝Φ（ｘ_１ｉ＋ｎ_１ｉ|ｐ１）−Φ（ｘ_１ｉ|ｐ１）により法線ベクトルを計算する。そして、ｎ_１ｉ’とｎ_２ｉから中間ベクトルｎ_ｍｉを算出する。そして、第４の実施形態のステップＳ１４０２０で説明した方法により共分散行列Ｃｏｖ（ｘ）を算出し、同実施形態の式（９）に記載し評価関数を用いて第１の変形変位場を取得するようにできる。以上に説明した方法によれば、式（１０）を用いる第６の実施形態に記載の方法に比べて、より簡易な計算により本発明を実施することができる。

また、第２の画像の座標系に変換した第１の法線ベクトルｎ_１ｉ’と第２の法線ベクトルｎ_２ｉとの中間ベクトルを算出する際に、これらのベクトルに対して異なる重み付けをしても良い。例えば、第１の画像および第２の画像の撮影装置（モダリティ）の種別に基づいて、この重みを設定するようにできる。例えば、第１の画像が第２の画像と比べて被検体の表面形状を高い信頼度で算出できるモダリティで撮影されている場合には、第１の法線ベクトルに対する重みを大きくするようにできる。これによれば、より信頼できる表面形状に関する情報を重視して法線ベクトルを算出できるため、第１の変形変位場を高い精度で取得できる効果がある。また前記の重みづけは、モダリティの種別に基づく場合に限らず、例えば第１の表面形状および第２の表面形状の信頼度を取得し、それに基づいて重みを設定するようにしても良い。この信頼度は、例えばユーザによる入力操作によって取得するようにしても良いし、これ以外にも第１および第２の体表形状を表す点群の空間密度や、第１および第２の画像のノイズレベル等に基づいて取得するようにしても良い。この場合、体表対応点の全てに対して同一の重みを設定しても良いし、体表対応点の各々に関して個別に重みを設定するようにしても良い。また設定する重みは０から１までの連続値であっても良いし、０または１といった二値の重みであっても良い。

（変形例５−３）：夫々の法線方向で評価値を統合して評価値を算出する
第６の実施形態では、式（１０）に示す通り、第１の法線ベクトルに基づいて算出される評価値と、第２の法線ベクトルに基づいて算出される評価値との合算により全体の評価値を計算する場合を例として説明した。しかし本発明の実施はこれに限らない。例えば、式（１１）に示すようにして評価関数を構成するようにしても良い。

ここで関数ｍｉｎ（ａ，ｂ）は最小値選択関数であり、関数ｍｉｎの１番目の引数である評価値と、２番目の引数である評価値のうち、評価が良い方（値が小さい方）の評価値が、全体の評価として採用される。これによれば、第１の法線ベクトルと第２の法線ベクトルの向きが大きく異なる場合に、より適切な変形変位場を取得できる効果がある。

＜第６の実施形態＞
本発明は、上述の実施形態の１以上の機能を実現するプログラムを、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサーがプログラムを読出し実行する処理でも実現可能である。また、１以上の機能を実現する回路（例えば、ＡＳＩＣ）によっても実現可能である。

１００ 画像処理装置 １１０ データサーバ １２０ 表示部 １０１０ データ取得部 １０２０ 関心領域取得部 １０３０ 第１の変形取得部 １０４０ 変形近似部 １０５０ 対応情報生成部 １０６０ 第２の変形取得部 １０７０ 表示制御部

Claims (19)

1. 被検体を異なる条件で撮像することにより得られた第１の画像と第２の画像を取得するデータ取得手段と、
   前記第１の画像内の関心領域を取得する関心領域取得手段と、
   前記第１の画像と前記第２の画像の間の第１の変形変位場を取得する第１の変形取得手段と、
   前記関心領域における前記第１の変形変位場を、前記第１の変形変位場より自由度が少ない近似変換モデルを用いて近似して近似変位場を生成する変形近似手段と、
   前記関心領域に関して、前記第１の画像と、前記第１の画像が前記近似変位場に基づいて変位される画像との間の対応情報を生成する対応情報生成手段と、
   前記対応情報を用いて、前記第１の画像と前記第２の画像の間の第２の変形変位場を取得する第２の変形取得手段と、
   を備えることを特徴とする画像処理装置。
2. 前記第１の変形取得手段は、前記第１の変形変位場に基づいて、前記第１の画像を変位させた第１の変形画像を生成し、
   前記第２の変形取得手段は、前記第２の変形変位場に基づいて、前記第１の画像を変位させた第２の変形画像を生成することを特徴とする請求項１に記載の画像処理装置。
3. 前記対応情報は、前記関心領域上に設定された複数の代表点と、該代表点の夫々を前記近似変位場を用いて変位させた複数の代表点とに基づいて取得される、仮想的な対応点であることを特徴とする請求項２に記載の画像処理装置。
4. 前記第２の変形取得手段は、前記設定された複数の代表点と、前記変位させた複数の代表点の何れか一方を変位させて、前記設定された複数の代表点と前記変位させた複数の代表点との位置関係を一致させるように前記第２の変形変位場を取得することを特徴とする請求項３に記載の画像処理装置。
5. 前記変位させた複数の代表点は、前記関心領域の面上に沿って動くことが許容されることを特徴とする請求項３に記載の画像処理装置。
6. 前記対応情報は、前記関心領域上に設定された領域と、該設定された領域を前記近似変位場を用いて変位させた領域とに基づいて取得される、仮想的な形状であることを特徴とする請求項２に記載の画像処理装置。
7. 前記関心領域の性質に応じて、前記近似変位場を生成するために用いられる近似変換モデルを選択する近似変換モデル選択手段を更に備えることを特徴とする請求項１から６のいずれか１項に記載の画像処理装置。
8. 前記関心領域が複数存在する場合は、前記近似変換モデル選択手段は、それぞれの前記関心領域の性質に応じて、前記近似変位場を生成するために用いられる近似変換モデルを選択することを特徴とする請求項７に記載の画像処理装置。
9. 前記関心領域に関しての前記第１の変形変位場から前記第２の変形変位場への変化の程度を示す形状維持度を取得する形状維持度取得手段を更に備え、
   前記対応情報生成手段は、前記第１の画像と、前記第１の画像が前記近似変位場と前記形状維持度とに基づいて変位される画像との間の対応情報を生成することを特徴とする請求項１から８のいずれか１項に記載の画像処理装置。
10. 前記関心領域に関しての前記第１の変形変位場から前記第２の変形変位場への変化は線形であることを特徴とする請求項９に記載の画像処理装置。
11. 前記関心領域に関しての前記第１の変形変位場から前記第２の変形変位場への変化は非線形であることを特徴とする請求項９に記載の画像処理装置。
12. 前記第１の変形変位場は、前記第１の画像と前記第２の画像との間で実際に対応する実対応点の位置の誤差の評価に基づいて算出され、
    前記第２の変形変位場は、前記実対応点の位置の誤差と前記仮想的な対応点の位置の誤差の評価とに基づいて算出されることを特徴とする請求項３に記載の画像処理装置。
13. 前記第１の画像に描出された前記被検体の第１の表面形状を取得し、前記第２の画像に描出された該被検体の第２の表面形状を取得する表面形状取得手段をさらに備え、
    前記実対応点は、前記第１の表面形状と前記第２の表面形状との間で対応する表面対応点であることを特徴とする請求項１２に記載の画像処理装置。
14. 前記実対応点の位置の誤差の評価は、前記第１の表面形状と前記第２の表面形状の何れか一方の曲面の法線方向の位置の誤差は大きく評価し、該曲面に沿った方向の位置の誤差は小さく評価することを特徴とする請求項１３に記載の画像処理装置。
15. 前記近似変換モデルは剛体変換モデルであることを特徴とする請求項１から１４のいずれか１項に記載の画像処理装置。
16. 前記関心領域取得手段は、ユーザによる操作に基づいて前記関心領域を取得することを特徴とする請求項１から１５のいずれか１項に記載の画像処理装置。
17. 前記第１の変形画像と前記第２の変形画像を表示部に表示させるための制御を行う表示制御手段を更に備えることを特徴とする請求項２から６のいずれか１項に記載の画像処理装置。
18. 被検体を異なる条件で撮像することにより得られた第１の画像と第２の画像を取得するデータ取得工程と、
    前記第１の画像内の関心領域を取得する関心領域取得工程と、
    前記第１の画像と前記第２の画像の間の第１の変形変位場を取得する第１の変形取得工程と、
    前記関心領域における前記第１の変形変位場を、前記第１の変形変位場より自由度が少ない近似変換モデルを用いて近似して近似変位場を生成する変形近似工程と、
    前記関心領域に関して、前記第１の画像と、前記第１の画像が前記近似変位場に基づいて変位される画像との間の対応情報を生成する対応情報生成工程と、
    前記対応情報を用いて、前記第１の画像と前記第２の画像の間の第２の変形変位場を取得する第２の変形取得工程と、
    を備えることを特徴とする画像処理方法。
19. コンピュータを、請求項１から１７のいずれか１項に記載の画像処理装置として機能させるためのプログラム。

Images