JP6491011B2 - 医用画像処理装置および医用画像処理方法 - Google Patents

Description

本発明の実施形態は、医用画像を補正する医用画像処理装置および医用画像処理方法に関する。

脳磁気共鳴イメージング（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ：ＭＲＩ）において、Ｔ１Ｗ、Ｔ２Ｗなどの脳の３次元のＶｏｌｕｍｅ ｄａｔａに対して、ＳＰＭ（Ｓｔａｔｉｓｔｉｃａｌ Ｐａｒａｍｅｔｒｉｃ Ｍａｐｐｉｎｇ）などの統計解析ツールを用いたＶｏｘｅｌ−Ｂａｓｅｄ Ｍｏｒｐｈｏｍｅｔｒｙ （ＶＢＭ）と称する分野がある。ＶＢＭでは、位置合わせ後に個人脳の標準脳への解剖学的標準化が行われるが、一般に強いスムージングを施し微細な変化を消すようにしている。すなわち、ＶＢＭでは、対象画像を意図的にぼかすことで、標準脳への解剖学的標準化（位置合わせ）が実行される。

このＶＢＭという分野では、時間経過に伴う形態の微細な変化を追跡する必要がある認知症や精神神経疾患に応用されてきている。例として、アルツハイマー（Ａｌｚｈｅｉｍｅｒ：ＡＤ）など痴呆症の診断などでは、皮質の厚さをサブｍｍのオーダーで、また海馬などの脳容積を数％のオーダーで、ＭＲＩを用いて、経時的な微小な変化を、１〜数ヶ月単位で数年間〜生涯に亘って随時追跡する必要がある。

また、ＶＢＭと類似の解析手法は、構造データのみならず、拡散画像、機能的（ｆｕｎｃｔｉｏｎａｌ）ＭＲＩ（ｆＭＲＩ）、脳血流（Ａｒｔｅｒｉａｌ Ｓｐｉｎ Ｌａｂｅｌｉｎｇ：ＡＳＬ）などの機能データでも広く応用されてきている。モダリティ別ではＭＲＩ装置に限定されず単光子放出コンピュータ断層撮影（Ｓｉｎｇｌｅ Ｐｈｏｔｏｎ Ｅｍｉｓｓｉｏｎ Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ：ＳＰＥＣＴ）装置、陽電子放出コンピュータ断層撮影（Ｐｏｓｉｔｉｒｏｎ Ｅｍｉｓｓｉｏｎ Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ：ＰＥＴ）装置、Ｘ線ＣＴ（Ｘ線コンピュータ断層撮影）装置などでの応用も多い。

胸部Ｘ線画像では、ある期間毎に撮像した画像を用いて、計時差分をおこなって、以前との変化分を把握しやすくする方法が実用化されている。この場合、骨などの情報を用いて、平面内での２次元的な撮像位置に加え、透過画像特有の奥行き方向の患者撮像体位の補正などが行われている。また、Ｘ線ＣＴ、ＳＰＥＣＴ、ＰＥＴ、ＭＲＩなどの断面像を取得可能な撮像装置においても、以前の断面に切りなおすｒｅ−ｓｌｉｃｉｎｇ程度のことは行われている。

しかしながら、撮像時刻が異なる複数の医用画像において、装置の違いや撮像および画像処理に関するフィルターやマトリクスサイズの違いがあると、患者由来の経時変化因子による正しい経時的変化の追跡が難しい問題がある。また、同一装置条件（撮像条件）であっても、撮像体位の変化を補正する場合、画像の補間を行う必要がある。たとえフィルターが等価であったとしても、例えば、補間回数の相違により、ボケなどの効果が現れてしまうことがあり、現在それらの影響は考慮されていない問題がある。

特許第４０２２５８７号公報 特許第５５３２７３０号公報 特開２０１４−０４２６８４号公報 特開２０１０−０４６１０３号公報

目的は、対象物の解剖学的変化を精度よく追跡可能にするために、医用画像を補正可能な医用画像処理装置と医用画像処理方法とを提供することにある。

本実施形態に係る医用画像処理装置は、第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、を具備し、前記領域設定部は、前記第１画像の画素値に基づいて、前記第１領域を設定し、前記第２画像の画素値に基づいて、前記第２領域を設定することを特徴とする。

図１は、第１の実施形態に係る医用画像処理装置の構成を示す構成図である。 図２は、第１の実施形態に係り、第１補間画像および第２補間画像において、第１領域または第２領域の設定の一例を示す図である。 図３は、第１の実施形態に係り、第１領域または第２領域において、位置（例えば、空気部分）に対する画素値の変化（ノイズプロファイル）の一例を示す図である。 図４は、第１の実施形態に係り、ノイズプロファイルの逆フーリエ変換を用いて発生された第１領域周波数分布に対応するパワースペクトル、第２領域周波数分布に対応するパワースペクトル、およびこれらのパワースペクトルをそれぞれ平滑化した第１平滑スペクトル、第２平滑スペクトルの一例を示す図である。 図５は、第１の実施形態に係り、第１平滑スペクトルを第２平滑スペクトルに変換する変換関数の一例を示す図である。 図６は、第１の実施形態に係り、変換関数を用いて補正された第１周波数分布に対応するパワースペクトルの一例を示す図である。 図７は、第１の実施形態に係り、周波数分布補正処理の手順の一例を示すフローチャートである。 図８は、第２の実施形態に係る医用画像処理装置の構成を示す構成図である。 図９は、第２の実施形態に係り、第１補間画像および第２補間画像において、第１領域または第２領域の設定の一例を示す図である。 図１０は、第２の実施形態に係り、図９における線分ｅｅ’上の画素値の変化（実信号プロファイル）の一例を示す図である。 図１１は、第２の実施形態に係り、図９における線分ｄｄ’または図１０における抽出部分（ｆｆ’）上の画素値の変化（経時不変プロファイル）の一例を示す図である。 図１２は、第２の実施形態に係り、経時不変プロファイルの逆フーリエ変換の絶対値の正規化により発生された第１正規化周波数分布と第２正規化周波数分布との一例を示す図である。 図１３は、第２の実施形態に係り、第１正規化周波数分布を第２正規化周波数分布に変換する変換関数の一例を示す図である。 図１４は、第２の実施形態に係り、変換関数を用いて補正された第１周波数分布の絶対値の一例を示す図である。 図１５は、第２の実施形態に係り、周波数分布補正処理の手順の一例を示すフローチャートである。

以下、本実施形態に係る医用画像処理装置について、図面を参照して説明する。なお、以下の説明において、略同一の機能および構成を有する構成要素については、同一符号を付し、重複説明は必要な場合に行う。

（第１の実施形態）
図１は、第１の実施形態に係る医用画像処理装置１の構成を示すブロック図である。医用画像処理装置１は、補間処理部１１と、領域設定部１３と、周波数分布発生部１５と、補正部１７と、画像発生部１９と、表示部２１と、インターフェース部２３と、入力部２５と、記憶部２７と、制御部２９とを有する。図１におけるＡは、本医用画像処理装置１におけるインターフェース部２３を介して、外部の各種モダリティ、医用画像保管装置などと接続するネットワークを表している。

入力部２５を介した操作者の指示により、第１画像と第２画像とが選択されると、補間処理部１１は、第１画像と第２画像とに対して（補間関数を用いた）補間処理を実行する。補間処理とは、例えば、第１画像と第２画像とにおいて、ボクセル（またはピクセル）サイズの統一、空間位置補正、位置合わせなどである。第１画像は、第２画像に対して、比較対象となる医用画像（最新または現在の医用画像）である。第２画像は、第１画像に対して基準（Ｂａｓｅ）となる医用画像である。具体的には、第２画像は、通常は現在（第１画像）と比較され、第１画像よりも過去に取得された医用画像、または過去に収集されたボリュームデータを基準断面に対して再スライス（ｒｅ−ｓｌｉｃｅ）された医用画像（例えば、ＭＰＲ（Ｍｕｌｔｉ−Ｐｌａｎａｒ Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ）画像）である。

第１画像および第２画像は、被検体に関する同一の解剖学的領域を有する医用画像である。以下、第１画像および第２画像は、脳画像として説明する。なお、第１画像および第２画像は、脳画像に限定されず、任意の解剖学的領域を有する画像であってもよい。第１画像と第２画像とは、例えば、癌に対応する領域を有する医用画像、骨粗鬆症に関する医用画像、癌の進行の評価（浸潤度）として用いられる医用画像、フォローアップ検査に関する医用画像などである。また、第１画像および第２画像は、同一モダリティにより発生された医用画像である。

具体的には、第２画像の選択の後に、補間処理部１１は、例えば、第２画像に基づいて、基準（Ｂａｓｅ）となる画像条件を決定する。画像条件とは、例えば、マトリクスサイズ（ｍａｔｒｉｘ ｓｉｚｅ）、撮像視野（Ｆｉｅｌｄ ｏｆ ｖｉｅｗ：ＦＯＶ）に基づくボクセル（ｖｏｘｅｌ）サイズ、患者の撮像体位などである。なお、基準となる画像条件は、第１画像の付帯情報に基づいて決定されてもよい。また、基準となる画像条件は、第２画像の付帯情報と第１画像の付帯情報とに基づいて決定されてもよい。

例えば、第１画像および第２画像における空間位置の相違に対する補間処理は、第１画像および第２画像が２次元画像なら面内での位置の補正、第１画像および第２画像が３次元画像ならスライス断面の補正に相当する。また、ボクセルサイズおよび空間位置の相違に対する処理は「補間」であるため、第１画像および第２画像に対して同様な処理（同様な補間関数を用いた補間処理）を実行することが望ましい。

なお、基準画像（第2画像）に対して対象画像（第1画像）を位置合わせする場合の補間処理は、極力高精度の周波数劣化の小さな手法で行うことが重要であると同時にその補間処理時のパワースペクトルの劣化を加味して変換フィルター特性を決める必要がある。

領域設定部１３は、補間処理後の第１画像（以下、第１補間画像と呼ぶ）の一部分に対応する第１領域を設定する。領域設定部１３は、補間処理後の第２画像（以下、第２補間画像と呼ぶ）に対応する第２領域を設定する。第１領域は、第１補間画像において、空間的に均一な領域である。第２領域は、第２補間画像において、空間的に均一である。

すなわち、第１部域および第２領域は、例えば、ノイズ部分に関する領域であって、無信号部分（空気領域）および構造がない（非構造であって）空間的に均一な部分に関する領域である。第１領域および第２領域は、例えば、脳画像において、空気に対応する領域および脳室に対応する領域などである。以下、第１領域および第２領域は、非構造であって空間的に均一な部分（無信号部分、ノイズ部分）に関する領域、すなわち空気領域とする。なお、第１領域および第２領域が経時的に変化の小さい領域（経時的に不変な領域）である場合については、第２の実施形態で詳述する。

図２は、第１補間画像および第２補間画像において、第１領域または第２領域の設定の一例を示す図である。図２におけるａは、補間画像における空気の領域を示している。図２におけるｂは、補間画像における骨の領域を示している。図２におけるｃは、補間画像における頭皮の領域を示している。図２におけるＧＭは、補間画像における灰白質（Ｇｒａｙ Ｍａｔｔｅｒ）の領域を示している。図２におけるＷＭは、補間画像における白質（Ｗｈｉｔｅ Ｍａｔｔｅｒ）の領域を示している。図２におけるＣＳＦは、補間画像における脳脊髄液（ＣｅｒｅｂｒｏＳｐｉｎａｌ Ｆｌｕｉｄ：ＣＳＦ）を有する領域を示している。図２において、左右対称であって画像中心部分のＣＳＦは、脳室に対応する。図２におけるＮは、例えば、第１領域および第２領域が設定される領域に対応する。

領域設定部１３は、第1補間画像の画素値に基づいて第１領域を設定し、第２補間画像の画素値に基づいて前記第２領域を設定する。領域設定部１３は、例えば、第１補間画像における複数の画素値に対する閾値により、第１領域を第１補間画像に設定する。領域設定部１３は、第２補間画像における複数の画素値に対する閾値により、第２領域を第２補間画像に設定する。上記閾値とは、例えば、空気を抽出可能な画素値である。この\閾値は、記憶部２７に記憶される。

以下、説明を簡単にするために、第１領域と第２領域との形状は、矩形であるものとする。なお、第１領域と第２領域との形状は、矩形に限定されず、円形、楕円形等の任意の形状であってもよい。また、第１領域および第２領域は、空気の領域またはＣＳＦの領域であってもよい。なお、第１領域または第２領域を空気の領域に設定できない場合、領域設定部１３は、例えば、ＣＳＦの領域に第１領域と第２領域とを設定してもよい。

第１領域の大きさと第２領域の大きさとは、同一の大きさが好ましい。なお、第１領域の大きさと第２領域の大きさとは、異なっていてもよい。第１領域および第２領域の形状および大きさは、予め設定されているものとする。なお、第１領域および第２領域の形状および大きさは、入力部２５を介した操作者の指示により、適宜変更可能である。なお、第１領域と第２領域とは、入力部２５を介した操作者の指示により、設定されてもよい。第１領域の大きさと第２領域の大きさとは、例えば、数百の画素を有する。

なお、第１補間画像および第２補間画像が振幅画像である場合、信号対雑音比（Ｓｉｇｎａｌ−Ｎｏｉｓｅ Ｒａｔｉｏ：ＳＮＲ）は十分大きい（ＳＮＲ＞５）ため、非空気領域（信号部）の複数の信号値は、ガウシアン（Ｇａｕｓｓｉａｎ）分布を有する。一方、空気領域（無信号部）における複数の信号値（ｎｏｉｓｅ）は、例えば、負値がない（非負の）ライシアン（Ｒｉｃｉａｎ）分布となる。このため、第１領域（および第２領域）の設定後の処理において、信号部と無信号部とを別々に扱うか、第１補間画像（第２補間画像）から第１補間画像（第２補間画像）のＤＣ（Ｄｉｒｅｃｔ Ｃｕｒｅｎｔ：直流）成分（平均画素値）を差分した差分画像において（ＤＣ成分除去後）、絶対値化（ｍａｇｎｉｔｕｄｅ化）することにより、第１領域（および第２領域）における無信号部分をライシアンノイズに変換する必要がある。

以下、第１補間画像に設定された第１領域に含まれる複数の画素にそれぞれ対応する複数の画素値を、Ｉ_ｏｒｉｇ（ｘ）とする。また、第２補間画像に設定された第２領域に含まれる複数の画素にそれぞれ対応する複数の画素値を、Ｉ_ｂａｓｅ（ｘ）とする。説明を簡単にするために、第１領域に含まれる複数の画素の位置と、第２領域に含まれる複数の画素の位置とは同じであるものとする。このとき、第１補間画像と第２補間画像とは、補間処理部１１により、予め位置合わせされる。また、第１領域と第２領域とのうち一方が設定されると、他方が自動的に設定されてもよい。Ｉ_ｏｒｉｇ（ｘ）とＩ_ｂａｓｅ（ｘ）とにおける引数ｘは、第１領域および第２領域に含まれる画素の位置を示す。また、引数ｘが含まれる範囲は同一となる。なお、一般的に、第１領域と第２領域とが異なる位置および異なる大きさで設定された場合、引数はそれぞれ異なる文字（例えば、Ｉ_ｏｒｉｇ（ｘ）とＩ_ｂａｓｅ（ｙ）など）により規定され、引数の範囲も異なる。

図３は、第１領域（または第２領域）において、位置（例えば、空気部分）に対する画素値の変化（ノイズプロファイル）の一例を示す図である。図３の縦軸は、画素値を示している。図３の縦軸におけるａｖｅは、第１領域（または第２領域）の平均画素値（ＤＣ成分）に相当する。

第１領域および第２領域が２次元の場合、図３に対応するノイズプロファイルは、３次元的なグラフ（ｘ軸及びｙ軸は画素の位置、ｚ軸はノイズの大きさ）となる。また、第１画像および第２画像がボリュームデータである場合、図３に対応するノイズプロファイルは、４次元的なグラフ（ｘ軸、ｙ軸、ｚ軸は画素の位置、ｗ軸はノイズの大きさ）となる。

周波数分布発生部１５は、第１領域に含まれる複数の画素値に基づいて、第１領域の周波数分布の大きさの平方を示す第１パワースペクトルを発生し、第２領域における複数の画素値に基づいて、第２領域の周波数分布の大きさの平方を示す第２パワースペクトルを発生する。具体的には、周波数分布発生部１５は、第１領域および第２領域各々に含まれる複数の画素値に対して、逆フーリエ変換を実行する。次いで、周波数分布発生部１５は、逆フーリエ変換の結果に対して、絶対値および平方をそれぞれ適用することにより、第１パワースペクトルと第２パワースペクトルとを発生する。

第１パワースペクトルをＰ_ｏｒｉｇ（ｆ）、第２パワースペクトルをＰ_ｂａｓｅ（ｆ）と表記すると、周波数分布発生部１５において実行される計算は、例えば、以下のように表すことができる。ここで、ｆは、周波数であり第１パワースペクトルと第２パワースペクトルとにおける引数である。
Ｐ_ｏｒｉｇ（ｆ）＝｜ＩＦＴ［Ｉ_ｏｒｉｇ（ｘ）］｜^２
Ｐ_ｂａｓｅ（ｆ）＝｜ＩＦＴ［Ｉ_ｂａｓｅ（ｘ）］｜^２ ・・・（１）
上式におけるＩＦＴは、逆フーリエ変換を示している。

すなわち、周波数分布発生部１５は、第１領域および第２領域各々における複数の画素値（画像データ）に対して、第１（補間）画像の次元に応じた２次元または３次元の逆フーリエ変換を実行する。逆フーリエ変換に対する絶対値演算及び平方演算により、第１領域及び第２領域にそれぞれ対応する第１パワースペクトルと第２パワースペクトルとを発生する。

一般に、画像化時のゲインに依存して、第１画像と第２画像とにおける画素値のスケーリング（ｓｃａｌｉｎｇ）は異なる。このため、周波数分布発生部１５は、第１パワースペクトルと第２パワースペクトルとを、基準となるパワースペクトルの最大値（通常はＤＣ（ｆ＝０）成分）で正規化（ｎｏｒｍａｌｉｚｅ）する。説明を簡単にするために、基準となるパワースペクトルは、第２パワースペクトルの最大値Ｐ_ｂａｓｅ（０）であるものとする。なお、基準となるパワースペクトルが第１パワースペクトルである場合、および所定の周波数分布に対応するパワースペクトルである場合については、後ほど変形例で説明する。

具体的には、周波数分布発生部１５は、第１パワースペクトルＰ_ｏｒｉｇ（ｆ）を第２パワースペクトルの最大値Ｐ_ｂａｓｅ（０）で除することにより、第１パワースペクトルを正規化する。また、周波数分布発生部１５は、第２パワースペクトルＰ_ｂａｓｅ（ｆ）を第２パワースペクトルの最大値Ｐ_ｂａｓｅ（０）で除することにより、第２パワースペクトルを正規化する。

すなわち、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）は、Ｐ_ｏｒｉｇ（ｆ）／Ｐ_ｂａｓｅ（０）として計算される。正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）は、Ｐ_ｂａｓｅ（ｆ）／Ｐ_ｂａｓｅ（０）として計算される。周波数分布発生部１５は、正規化（ｎｏｒｍａｌｉｚｅ）係数αを計算する。正規化係数αは、例えば、以下に示す式のように、計算される。
α＝Ｐ_ｂａｓｅ（０）／Ｐ_ｏｒｉｇ（０）・・・（２）
なお、第１領域および第２領域が脳室などのノイズ以外の信号成分を有する場合、第１パワースペクトルと第２パワースペクトルとにおいて、ＤＣ成分が発生する。このため、第１領域と第２領域とが同一部分（同位置）であれば、正規化係数αは、上式（２）と同一の式で計算される。また、第１領域および第２領域において、ｆ＝０における信号成分の有無にかかわらず正規化係数αを上式（２）で計算する場合、Ｐ_ｂａｓｅ（０）とＰ_ｏｒｉｇ（０）とは、ｆ＞０におけるＰ_ｂａｓｅ（ｆ）とＰ_ｏｒｉｇ（ｆ）とからそれぞれ外挿により生成する。

周波数分布発生部１５は、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）と正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）とに対して、強いスムージング（平滑化）を実行する。このスムージングにより、周波数分布発生部１５は、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）を平滑化した第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）と、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）を平滑化した第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）とを発生する。なお、周波数分布発生部１５は、スムージングの代わりに、所定の関数の多項式、ガウシアン（Ｇａｕｓｓｉａｎ）関数、またはライシアン（Ｒｉｃｉａｎ）関数などを用いた曲線近似により、第１平滑スペクトルと第２平滑スペクトルとを発生してもよい。

また、周波数分布発生部１５は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）と、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）との発生において、カットオフ周波数を設けてもよい。第１画像（第１補間画像）および第２画像（第２補間画像）がデジタル画像である場合、一般に高周波成分は、有限なサンプリング（ｓａｍｐｌｉｎｇ）周波数の制約などから、折り返し歪（エイリアジングによる歪み）を含みやすい。このため、周波数分布発生部１５は、ナイキスト（Ｎｙｑｕｉｓｔ）周波数ｆｎ以下にカットオフ周波数ｆｃを設定し（ｆｎ≧ｆｃ）、過剰な高周波（ｆ＞ｆｎ≧ｆｃ）をカットする。

周波数分布発生部１５は、通常、第１画像と第２画像とのうちサンプリング周波数の低い側にカットオフ周波数を合わせる。サンプリング周波数が低いのは過去画像であるｂａｓｅ画像（第２画像）である場合が多いため、周波数分布発生部１５は、第２画像のサンプリング周波数ｆｎ以下に、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）におけるカットオフ周波数を決定する。

図４は、例えば、第１領域および第２領域にそれぞれ対応する図３に示すノイズプロファイルの逆フーリエ変換を用いて発生され、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）、および第１、第２パワースペクトルをそれぞれ平滑化した第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）の一例を示す図である。図４に示すように、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）は、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）を平滑化したパワースペクトルである。第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）は、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）を平滑化したパワースペクトルである。図４におけるｆｃは、カットオフ（ｃｕｔ ｏｆｆ）周波数を示している。

周波数分布発生部１５は、第１画像に対応する第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を発生する。具体的には、周波数分布発生部１５は、第１補間画像または第１画像の画素値Ｉ_ｏｒｉｇ（ｘ）に対して逆フーリエ変換することにより、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を発生する。すなわち、上記記載によれば、周波数分布発生部１５は、例えば以下に示す式で、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を計算する。

Ｓ_ｏｒｉｇ（ｆ）＝ＩＦＴ［Ｉ_ｏｒｉｇ（ｘ）］
周波数分布発生部１５は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）と、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）と、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）とを、補正部１７に出力する。なお、周波数分布発生部１５は、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）と、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）と、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）とを、補正部１７に出力してもよい。このとき、周波数分布発生部１５において、関数近似は実行されない。

補正部１７は、近似関数であるＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）とＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）との比により、第１領域におけるノイズレベル（例えば、ノイズの標準偏差）を第２領域におけるノイズレベルに合わせるための関数ＰＮＲ（ｆ）を決定する。具体的には、補正部１７は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）を、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）に変換する（合わせるまたは一致させる）変換関数ＰＮＲ（ｆ）を決定する。

すなわち、補正部１７は、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）を第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）で除することで、変換関数ＰＮＲ（ｆ）を決定する。変換関数ＰＮＲ（ｆ）は、例えば、ノイズのパワーの比（Ｐｏｗｅｒ Ｎｏｉｓｅ Ｒａｔｉｏ）に相当する。以下に示す式により、変換関数ＰＮＲ（ｆ）は決定される。なお、変換関数ＰＮＲ（ｆ）の高周波側は、カットオフ周波数ｆｃ以下に滑らかに低減される。
ＰＮＲ（ｆ）＝Ｐｎ＿ｅｓｔ_ｂａｓｅ（ｆ）／Ｐｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）・・・（３）
上式（３）において、ゼロによる除算が生じないように工夫する必要がある。例えば、周波数分布発生部１５および図４では、カットオフ周波数ｆｃがナイキスト周波数ｆｎより大きいため、Ｐｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）は、非ゼロとなる。

なお、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）と、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）とに対して関数近似が実行されない場合、補正部１７は、例えば、以下のようにして、変換関数ＰＮＲ（ｆ）を決定することができる。補正部１７は、Ｌ１、Ｌ２ノルム最適化などの最適化問題を解くためのアルゴリズム（最適化アルゴリズムなど）を適用した繰り返し演算（逐次近似などの演算）により、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）と、正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）とを比較しながら、変換関数ＰＮＲ（ｆ）を決定する。

図５は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）を、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）に変換する変換関数の一例を示す図である。図５に示すｆｃは、カットオフ周波数である。図４に対応するＰＮＲ（ｆ）では、図５におけるＰＮＲ（ｆ）の最大値未満とＰＮＲ（ｆ）がαより大きい値との間において、ＰＮＲ（ｆ）＝１の直線が存在する。

補正部１７は、変換関数ＰＮＲ（ｆ）を用いて、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、補正した第１周波数分布Ｓ_ｃｏｒ（ｆ）を画像発生部１９に出力する。具体的には、補正部１７は、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）に変換関数ＰＮＲ（ｆ）を乗ずることにより、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、例えば、以下に示す式に従って第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。

Ｓ_ｃｏｒ（ｆ）＝Ｓ_ｏｒｉｇ（ｆ）×ＰＮＲ（ｆ）・・・（４）
図６は、変換関数を用いて補正された第１周波数分布のパワースペクトルＰｎ_{ｏｒｉｇ＿ｃｏｒ}（ｆ）の一例を示す図である。なお、図６におけるＰｎ_{ｏｒｉｇ＿ｃｏｒ}（ｆ）は、第２パワースペクトルの最大値Ｐ_ｂａｓｅ（０）で正規化されている。このため、図６に示すように、グラフの縦軸とＰｎ_{ｏｒｉｇ＿ｃｏｒ}（ｆ）との交点は、１となる。

画像発生部１９は、補正された第１周波数分布Ｓ_ｃｏｒ（ｆ）に対してフーリエ変換ＦＴを適用することで、第１画像（第１補間画像）を補正した画像Ｉ_ｃｏｒ（ｘ）（以下、補正第１画像と呼ぶ）を発生する。具体的には、画像発生部１９は、以下に示す式に従って、補正第１画像を発生する。

Ｉ_ｃｏｒ（ｘ）＝ＦＴ［Ｓ_ｃｏｒ（ｆ）］・・・（５）
上式（１）乃至（５）は説明を簡単にするため、１次元の式として示したが、補正対象の第１画像（または第２画像）は、２次元画像または３次元画像である。このため、上式（１）乃至（５）は画像の次元に対応する式となる。例えば、補正対象の次元が２次元である場合、本補正は、ｘ、ｙ方向ごとに実行される。また、補正対象の次元が３次元である場合、本補正は、ｘ、ｙ、ｚ方向ごとに実行されるため、高次元の方がより高精度な補正として実行される。

画像発生部１９は、発生した補正第１画像を、表示部２１に出力する。画像発生部１９は、補正第１画像を記憶部２７に出力してもよい。なお、画像発生部１９は、補正第１画像と第２画像とに基づいて、第２画像に対する補正第１画像の経時的変化を示す画像（以下、経時的変化画像と呼ぶ）を発生してもよい。画像発生部１９は、経時的変化画像を、表示部２１および記憶部２７に出力する。

表示部２１は、画像発生部１９により発生された各種医用画像、記憶部２７に記憶された医用画像を表示する。表示部２１は、入力操作に係る各種入力画面を表示する。

インターフェース部２３は、ネットワークＡ、図示していない外部記憶装置などに関するインターフェースである。インターフェース部２３は、図示していない医用画像撮影装置および医用画像保管装置に対して電子的通信回線、典型的にはローカルエリアネットワーク（Ｌｏｃａｌ Ａｒｅａ Ｎｅｔｗｏｒｋ：以下ＬＡＮと呼ぶ）などを介して接続される。医用画像撮影装置で発生されたボリュームデータおよび画像データ、または医用画像保管装置で保管されたボリュームデータおよび画像データは、インターフェース部２３を介して、記憶部２７に記憶される。また、本医用画像処理装置１で補正された医用画像のデータ、ボリュームデータなどは、ネットワークを介して図示していない医用画像保管装置などの他の装置に転送可能である。

入力部２５は、操作者などからの各種指示・命令・情報・選択・設定などを制御部２９に入力する。具体的には、入力部２５は、上記各種指示・命令・情報・選択・設定などを入力するためのトラックボール、スイッチボタン、マウス、マウス・ホイール、キーボード等の入力デバイスを有する。なお、入力デバイスは、表示部２１における表示画面を覆うタッチパネルでもよい。入力部２５は、入力デバイスを介して、補正対象の第１画像と、第１画像より過去に同種のモダリティにより収集され第２画像との選択操作を入力する。なお、入力部２５は、第１領域および第２領域の位置、形状および大きさを入力することも可能である。

記憶部２７は、インターフェース部２３を介して転送された第１画像および第２画像、周波数分布発生部１５により発生された各種周波数分布、パワースペクトル、補正部１７により決定された変換関数、画像発生部１９により発生された各種医用画像などを記憶する。なお、記憶部２７は、領域設定部１３で用いられる第１領域および第２領域の形状および大きさを記憶してもよい。記憶部２７は、制御部２９により用いられる各種制御プログラムを記憶する。また、記憶部２７は、後述する周波数分布補正処理に関する補正プログラム、例えば、補間処理、第１、第２領域の設定、周波数分布および各種パワースペクトルの発生、パワースペクトルの正規化および平滑化、変換関数の決定、周波数分布の補正、補正画像の発生等に関する少なくとも一つ以上のプログラムを記憶する。

制御部２９は、図示していない中央処理装置（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ：以下、ＣＰＵと呼ぶ）、メモリ等を有する。制御部２９は、本医用画像処理装置１の中枢として機能する。具体的には、制御部２９は、入力部２５を介して入力された入力操作に応じて、記憶部２７から各種プログラムなどを読み出して自身のメモリ上に展開する。制御部２９は、展開したプログラムを実行することにより、本医用画像処理装置１の各部を制御する。また、例えば、制御部２９は、記憶部２７から補正プログラムを読み出して実行することにより、補間処理部１１、領域設定部１３、周波数分布発生部１５、補正部１７、画像発生部１９などの各部を制御する。

なお、本医用画像処理装置１は、各種モダリティ、例えば、磁気共鳴イメージング（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ：ＭＲＩ）装置、Ｘ線コンピュータ断層撮影（Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ：ＣＴ）装置、単光子放出コンピュータ断層撮影（Ｓｉｎｇｌｅ Ｐｈｏｔｏｎ Ｅｍｉｓｓｉｏｎ Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ：ＳＰＥＣＴ）装置、陽電子放出コンピュータ断層撮影（Ｐｏｓｉｔｉｒｏｎ Ｅｍｉｓｓｉｏｎ Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ：ＰＥＴ）装置、Ｘ線診断装置などのコンソール（操作卓等）に搭載されてもよい。また、本医用画像処理装置１は、ワークステーションとして実現することも可能である。また、医用画像処理装置１は、例えば、シンクライアント（Ｔｈｉｎ Ｃｌｉｅｎｔ）およびクラウドシステムにおけるサーバとして実現されてもよい。

（周波数分布補正機能）
周波数分布補正機能とは、補正対象の医用画像（第１画像）のノイズレベルと、基準となる医用画像（第２画像）のノイズレベルとを一致させるように、補正対象の医用画像の周波数分布を補正する機能である。以下、周波数分布補正機能に関する処理（以下、周波数分布補正処理と呼ぶ）について説明する。

図７は、周波数分布補正処理に関する手順の一例を示すフローチャートである。
入力部２５を介した操作者の指示により、選択された第１画像と第２画像とが、記憶部２７から読み出される。なお、選択された第１画像および第２画像は、インターフェース部２３とネットワークとを介して、図示していない医用画像保管装置、各種モダリティ等から読み出されてもよい。

読み出された第１画像における複数の画素値に基づいて、空間的に均一な第１領域が設定される（ステップＳａ１）。読み出された第２画像における複数の画素値に基づいて、空間的に均一な第２領域が設定される（ステップＳａ２）。

第１領域における複数の画素値に基づいて、第１領域周波数分布に対応する第１パワースペクトルが発生される（ステップＳａ３）。第２領域における複数の画素値に基づいて、第２領域周波数分布に対応する第２パワースペクトルが発生される（ステップＳａ３）。第２パワースペクトルの最大値を用いて、第１、第２パワースペクトル各々が正規化される（ステップＳａ４）。

正規化された第１パワースペクトルを平滑化することにより、第１平滑スペクトルが発生される（ステップＳａ５）。正規化された第２パワースペクトルを平滑化することにより、第２平滑スペクトルが発生される（ステップＳａ５）。第１平滑スペクトルを第２平滑スペクトルに変換する変換関数が決定される（ステップＳａ６）。

第１画像（第１補間画像）の全域における複数の画素値に基づいて、第１周波数分布が発生される（ステップＳａ７）。変換関数を用いて、第１周波数分布が補正される（ステップＳａ８）。補正された第１周波数分布に基づいて、第１画像を補正した補正第１画像が発生される（ステップＳａ９）。

（第１の変形例）
第１の実施形態との相違は、基準として第１パワースペクトルを用いることにある。
周波数分布発生部１５は、第１パワースペクトルＰ_ｏｒｉｇ（ｆ）を第１パワースペクトルの最大値Ｐ_ｏｒｉｇ（０）で除することにより、第１パワースペクトルを正規化する。また、周波数分布発生部１５は、第２パワースペクトルＰ_ｂａｓｅ（ｆ）を第１パワースペクトルの最大値Ｐ_ｏｒｉｇ（０）で除することにより、第２パワースペクトルを正規化する。

すなわち、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）は、Ｐ_ｏｒｉｇ（ｆ）／Ｐ_ｏｒｉｇ（０）として計算される。正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）は、Ｐ_ｂａｓｅ（ｆ）／Ｐ_ｏｒｉｇ（０）として計算される。周波数分布発生部１５は、正規化（ｎｏｒｍａｌｉｚｅ）係数αを計算する。正規化係数αは、例えば、以下に示す式のように、計算される。

α＝Ｐ_ｏｒｉｇ（０）／Ｐ_ｂａｓｅ（０）
周波数分布発生部１５は、第２画像に対応する第２周波数分布Ｓ_ｂａｓｅ（ｆ）を発生する。具体的には、周波数分布発生部１５は、第２補間画像または第２画像の画素値Ｉ_ｂａｓｅ（ｘ）に対して逆フーリエ変換することにより、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を発生する。すなわち、上記記載によれば、周波数分布発生部１５は、例えば以下に示す式で、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を計算する。

Ｓ_ｂａｓｅ（ｆ）＝ＩＦＴ［Ｉ_ｂａｓｅ（ｘ）］
補正部１７は、近似関数であるＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）とＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）との比により、第２領域におけるノイズレベルを第１領域におけるノイズレベルに合わせるための関数ＰＮＲ（ｆ）を決定する。具体的には、補正部１７は、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）を第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）に変換する（合わせるまたは一致させる）変換関数ＰＮＲ（ｆ）を決定する。

すなわち、補正部１７は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）を第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）で除することで、変換関数ＰＮＲ（ｆ）を決定する。変換関数ＰＮＲ（ｆ）は、例えば、ノイズのパワーの比（Ｐｏｗｅｒ Ｎｏｉｓｅ Ｒａｔｉｏ）に相当する。以下に示す式により、変換関数ＰＮＲ（ｆ）は決定される。
ＰＮＲ（ｆ）＝Ｐｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）／Ｐｎ＿ｅｓｔ_ｂａｓｅ（ｆ）
上式において、ゼロによる除算が生じないように、例えば、周波数分布発生部１５において、カットオフ周波数ｆｃを決定してもよい。なお、変換関数は、第１画像および第２画像の発生時における補間の周波数分布、補間回数、補間アルゴリズム、補間サイズ、適応されたフィルターの種類及び形状等の各種情報が既知であれば、補間の周波数分布、補間回数、適応されたフィルター等を用いて変換関数を決定してもよい。

補正部１７は、変換関数ＰＮＲ（ｆ）を用いて、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。補正部１７は、補正した第２周波数分布Ｓ_ｃｏｒ２（ｆ）を画像発生部１９に出力する。具体的には、補正部１７は、第２周波数分布Ｓ_ｂａｓｅ（ｆ）に変換関数ＰＮＲ（ｆ）を乗ずることにより、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。補正部１７は、例えば、以下に示す式に従って第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。

Ｓ_ｃｏｒ２（ｆ）＝Ｓ_ｂａｓｅ（ｆ）×ＰＮＲ（ｆ）
画像発生部１９は、補正された第２周波数分布Ｓ_ｃｏｒ２（ｆ）に対してフーリエ変換ＦＴを適用することで、第２画像（第２補間画像）を補正した画像Ｉ_ｃｏｒ２（ｘ）（以下、補正第２画像と呼ぶ）を発生する。具体的には、画像発生部１９は、以下に示す式に従って、補正第２画像を発生する。

Ｉ_ｃｏｒ２（ｘ）＝ＦＴ［Ｓ_ｃｏｒ２（ｆ）］
画像発生部１９は、発生した補正第２画像を、表示部２１に出力する。画像発生部１９は、補正第２画像を記憶部２７に出力してもよい。なお、画像発生部１９は、補正第２画像と第１画像とに基づいて、第１画像に対する補正第２画像の経時的変化を示す画像を発生してもよい。

（第２の変形例）
第１の実施形態および第１の変形例との相違は、基準として所定の周波数分布に対応するパワースペクトルを用いることにある。
記憶部２７は、所定の周波数分布に対応するパワースペクトル（以下、所定のパワースペクトルと呼ぶ）Ｐ_{ｐｒｅｄｅｔ}（ｆ）を記憶する。所定の周波数分布とは、例えば、ホワイトノイズに対応する周波数分布、または固定パターンノイズに対応する周波数分布などである。また、記憶部２７は、正規化された所定のパワースペクトル（以下、正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）を記憶する。記憶された所定のパワースペクトルＰ_{ｐｒｅｄｅｔ}（ｆ）は、周波数分布発生部１５により、記憶部２７から周波数分布発生部１５に読み出される。記憶された正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）は、補正部１７により、記憶部２７から補正部１７に読み出される。

周波数分布発生部１５は、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）と、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を発生する。周波数分布発生部１５は、第１パワースペクトルＰ_ｏｒｉｇ（ｆ）を所定のパワースペクトルの最大値Ｐ_{ｐｒｅｄｅｔ}（０）で除することにより、第１パワースペクトルを正規化する。周波数分布発生部１５は、第２パワースペクトルＰ_ｂａｓｅ（ｆ）を所定のパワースペクトルの最大値Ｐ_{ｐｒｅｄｅｔ}（０）で除することにより、第２パワースペクトルを正規化する。正規化された第１パワースペクトルと、正規化された第２パワースペクトルとは、補正部１７に出力される。

具体的には、すなわち、正規化された第１パワースペクトルＰｎ_ｏｒｉｇ（ｆ）は、Ｐ_ｏｒｉｇ（ｆ）／Ｐ_{ｐｒｅｄｅｔ}（０）として計算される。正規化された第２パワースペクトルＰｎ_ｂａｓｅ（ｆ）は、Ｐ_ｂａｓｅ（ｆ）／Ｐ_{ｐｒｅｄｅｔ}（０）として計算される。周波数分布発生部１５は、第１パワースペクトルに対応する正規化（ｎｏｒｍａｌｉｚｅ）係数α１を計算する。周波数分布発生部１５は、第２パワースペクトルに対応する正規化（ｎｏｒｍａｌｉｚｅ）係数α２を計算する。正規化係数α１は、例えば、以下に示す式のように、計算される。

α１＝Ｐ_{ｐｒｅｄｅｔ}（０）／Ｐ_ｏｒｉｇ（０）
また、正規化係数α２は、例えば、以下に示す式のように、計算される。

α２＝Ｐ_{ｐｒｅｄｅｔ}（０）／Ｐ_ｂａｓｅ（０）
補正部１７は、正規化された第１パワースペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）と正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）との比により、第１領域におけるノイズレベルを、所定の周波数分布に関するノイズレベル（以下、所定のノイズレベルと呼ぶ）に合わせるための第１変換関数ＰＮＲ１（ｆ）を決定する。補正部１７は、正規化された第２パワースペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）と正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）との比により、第２領域におけるノイズレベルを、所定のノイズレベルに合わせるための第２変換関数ＰＮＲ２（ｆ）を決定する。

具体的には、補正部１７は、第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）を正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）に変換する（合わせるまたは一致させる）第１変換関数ＰＮＲ（ｆ）を決定する。すなわち、補正部１７は、正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）を第１平滑スペクトルＰｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）で除することで、第１変換関数ＰＮＲ１（ｆ）を決定する。補正部１７は、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）を正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）に変換する（合わせるまたは一致させる）第２変換関数ＰＮＲ（ｆ）を決定する。すなわち、補正部１７は、正規化パワースペクトルＰｎ_{ｐｒｅｄｅｔ}（ｆ）を第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）で除することで、第２変換関数ＰＮＲ（ｆ）を決定する。

第１変換関数ＰＮＲ１（ｆ）および第２変換関数ＰＮＲ２（ｆ）は、例えば、ノイズのパワーの比（Ｐｏｗｅｒ Ｎｏｉｓｅ Ｒａｔｉｏ）に相当する。以下に示す２式により、第１変換関数ＰＮＲ１（ｆ）と第２変換関数ＰＮＲ２（ｆ）とは、それぞれ決定される。

ＰＮＲ１（ｆ）＝Ｐｎ_{ｐｒｅｄｅｔ}（ｆ）／Ｐｎ＿ｅｓｔ_ｏｒｉｇ（ｆ）
ＰＮＲ２（ｆ）＝Ｐｎ_{ｐｒｅｄｅｔ}（ｆ）／Ｐｎ＿ｅｓｔ_ｂａｓｅ（ｆ）
補正部１７は、第１変換関数ＰＮＲ１（ｆ）を用いて、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、第２変換関数ＰＮＲ２（ｆ）を用いて、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。補正部１７は、補正した第１周波数分布Ｓ_ｃｏｒ１（ｆ）と、補正した第２周波数分布Ｓ_ｃｏｒ２（ｆ）とを画像発生部１９に出力する。

具体的には、補正部１７は、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）に第１変換関数ＰＮＲ１（ｆ）を乗ずることにより、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、例えば、以下に示す式に従って第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。

Ｓ_ｃｏｒ１（ｆ）＝Ｓ_ｏｒｉｇ（ｆ）×ＰＮＲ１（ｆ）
また、補正部１７は、第２周波数分布Ｓ_ｂａｓｅ（ｆ）に第２変換関数ＰＮＲ２（ｆ）を乗ずることにより、第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。補正部１７は、例えば、以下に示す式に従って第２周波数分布Ｓ_ｂａｓｅ（ｆ）を補正する。

Ｓ_ｃｏｒ２（ｆ）＝Ｓ_ｂａｓｅ（ｆ）×ＰＮＲ２（ｆ）
画像発生部１９は、補正された第１周波数分布Ｓ_ｃｏｒ１（ｆ）に対してフーリエ変換ＦＴを適用することで、第１画像（第１補間画像）を補正した補正第１画像Ｉ_ｃｏｒ１（ｘ）を発生する。具体的には、画像発生部１９は、以下に示す式に従って、補正第１画像を発生する。

Ｉ_ｃｏｒ１（ｘ）＝ＦＴ［Ｓ_ｃｏｒ１（ｆ）］
画像発生部１９は、補正された第２周波数分布Ｓ_ｃｏｒ２（ｆ）に対してフーリエ変換ＦＴを適用することで、第２画像（第２補間画像）を補正した補正第２画像Ｉ_ｃｏｒ２（ｘ）を発生する。具体的には、画像発生部１９は、以下に示す式に従って、補正第２画像を発生する。

Ｉ_ｃｏｒ２（ｘ）＝ＦＴ［Ｓ_ｃｏｒ２（ｆ）］
画像発生部１９は、補正第１画像と補正第２画像とを、表示部２１に出力する。画像発生部１９は、補正第１画像と補正第２画像とを記憶部２７に出力してもよい。なお、画像発生部１９は、補正第１画像と補正第２画像とに基づいて、補正第２画像に対する補正第１画像の経時的変化を示す画像を発生してもよい。
以上に述べた構成によれば、以下の効果を得ることができる。
本実施形態に係る医用画像処理装置１によれば、非構造的であって空間的に均一な領域である第１領域と第２領域とを、第１画像と第２画像とにそれぞれ設定し、設定された第１領域における第１パワースペクトルと第２パワースペクトルとに基づいて第２画像のノイズレベルに第１画像のノイズレベルを合わせた補正第１画像を発生することができる。特に、本実施形態に係る周波数分布補正処理は、画像再構成手法やフィルターが線形アルゴリズムの場合に有効である。

また、本実施形態に係る第１の変形例によれば、第１パワースペクトルと第２パワースペクトルとに基づいて第１画像のノイズレベルに第２画像のノイズレベルを合わせた補正第２画像を発生することができる。

また、本実施形態に係る第２の変形例によれば、第１パワースペクトルと第２パワースペクトルと所定の周波数分布に関するパワースペクトルに基づいて、所定のノイズレベルに第１画像のノイズレベルと第２画像のノイズレベルとをそれぞれ合わせた補正第１画像と補正第２画像とを発生することができる。

すなわち、本医用画像処理装置１によれば、第１画像および第２画像の収集装置（モダリティ）、および収集条件のパラメータなどの装置由来の画像の変化の影響を除外した補正画像を、画像以外の情報がなくても発生することができる。これにより、意図的に医用画像をぼかすことなく、被検体の経時的な解剖学的・機能的変化、すなわち、本質的な経時変化を表示することができる。

以上のことから、本医用画像処理装置１によれば、ＶＢＭなどの分野を含んで広く、経時的な画像変化を追跡する場合において、撮像装置由来の因子を除外し患者本来の因子による経時変化を追跡することができる。これにより、操作者は、被検体に関する医用画像において、本来の経時的変化を、より詳細に追跡することが可能となる。

（第２の実施形態）
第１の実施形態との相違は、第１領域および第２領域が、経時的に変化の小さい領域または経時的に不変な領域（以下、経時的不変領域と呼ぶ）であって、ノイズ以外の信号成分を有することにある。
図８は、本実施形態に係る医用画像処理装置１の構成を示す構成図である。図１との相違は、エッジ検出部１０を有することにある。

エッジ検出部１０は、第１画像（または第１補間画像）において、所定のエッジ検出処理により、複数の第１エッジを検出する。エッジ検出部１０は、第２画像（または第２補間画像）において、所定のエッジ検出処理により、複数の第２エッジを検出する。

記憶部２７は、所定の閾値を記憶する、所定の閾値とは、第１画像および第２画像において、経時的不変領域を規定するエッジに対応する画素値である。経時的不変領域は、二つのエッジで挟まれた略一定（平坦）の画素値を有する領域であって、例えば、第１画像および第２画像において、被検体とともに撮影された所定のファントムを示す領域、骨の形状を示す領域（例えば、脳画像において頭蓋骨に対応する領域）、主として水を示す領域（例えば、脳画像におけるＣＳＦなど）などである。以下、説明を簡単にするために、経時的不変領域は、ファントムに対応する領域またはＣＳＦに対応する領域であるものとする。

領域設定部１３は、第１領域における対向する２辺が第１エッジに直交し、かつ所定の閾値で挟まれた画素値範囲を、第１領域として、第１補間画像に設定する。領域設定部１３は、第２領域における対向する２辺が第２エッジに直交し、かつ所定の閾値で挟まれた画素値範囲を、第２領域として、第２補間画像に設定する。第１領域と第２領域とは、画素値範囲が同じ領域（例えば、ファントム、骨、ＣＳＦ、水など）が好ましい。なお、第１領域および第２領域としてファントムを用いる場合、安定した周波数分布補正を実行するため、第１画像および第２画像各々の収集時と同時に、またはできるだけ小さい時間差で、ファントム画像を収集した方が好ましい。なお、第１領域および第２領域は、経時的変化が小さい骨、水抽出などにより、設定されてもよい。また、第１領域および第２領域は、入力部２５を介した操作者の指示により設定されてもよい。

図９は、第１補間画像および第２補間画像において、第１領域または第２領域の設定の一例を示す図である。図９におけるｆはファントムを示している。図９における線分ｄｄ’は、ファントムｆに対して設定される第１領域または第２領域の一例を示している。図９における線分ｅｅ’は、頭蓋骨に対応する領域およびＣＳＦに対応する領域に対して設定される第１領域または第２領域の一例を示している。

図１０は、図９における線分ｅｅ’上において一次元に配列した複数の画素にそれぞれ対応する複数の画素値による画素値の変化（実信号プロファイル）の一例を示す図である。図１０の横軸は、画素の位置ｘを示している。図１０の縦軸は、画像Ｉにおける画素の位置ｘの画素値Ｉ（ｘ）を示している。図１０において、第１領域および第２領域は、例えば、線分ｆｆ’で挟まれた領域として設定される。

図１１は、図９における線分ｄｄ’または図１０における抽出部分（ｆｆ’）上の画素値の変化（経時不変プロファイル）の一例を示す図である。図１１の横軸は、画素の位置ｘを示している。図１１の縦軸は、画像Ｉにおける画素の位置ｘの画素値Ｉ（ｘ）を示している。図１１におけるＩ_ｏｒｉｇ（ｘ）は、第１画像において線分ｄｄ’または線分ｆｆ’に位置する複数の画素にそれぞれ対応する複数の画素値のプロファイルを示している。図１１におけるＩ_ｂａｓｅ（ｘ）は、第２画像において線分ｄｄ’または線分ｆｆ’に位置する複数の画素にそれぞれ対応する複数の画素値のプロファイルを示している。図９乃至図１１において、第１領域および第２領域は１次元プロファイルとして説明しているが、２次元または３次元プロファイルであってもよい。

周波数分布発生部１５は、第１領域に対応する第１領域周波数分布Ｓ１_ｏｒｉｇ（ｆ）の絶対値｜Ｓ１_ｏｒｉｇ（ｆ）｜を発生する。具体的には、周波数分布発生部１５は、第１領域における複数の画素値Ｉ１_ｏｒｉｇ（ｘ）に対して逆フーリエ変換を実行して絶対値をとることにより、第１領域周波数分布の絶対値｜Ｓ１_ｏｒｉｇ（ｆ）｜を発生する。すなわち、上記記載によれば、周波数分布発生部１５は、例えば以下に示す式で、第１領域周波数分布の絶対値｜Ｓ１_ｏｒｉｇ（ｆ）｜を計算する。

｜Ｓ１_ｏｒｉｇ（ｆ）｜＝｜ＩＦＴ［Ｉ１_ｏｒｉｇ（ｘ）］｜
周波数分布発生部１５は、第２領域に対応する第２領域周波数分布Ｓ２_ｂａｓｅ（ｆ）の絶対値｜Ｓ２_ｂａｓｅ（ｆ）｜を発生する。具体的には、周波数分布発生部１５は、第２領域における複数の画素値Ｉ２_ｂａｓｅ（ｘ）に対して逆フーリエ変換を実行して絶対値をとることにより、第２領域周波数分布の絶対値｜Ｓ２_ｂａｓｅ（ｆ）｜を発生する。すなわち、上記記載によれば、周波数分布発生部１５は、例えば以下に示す式で、第２領域周波数分布の絶対値の｜Ｓ２_ｂａｓｅ（ｆ）｜を計算する。

｜Ｓ２_ｂａｓｅ（ｆ）｜＝｜ＩＦＴ［Ｉ２_ｂａｓｅ（ｘ）］｜
周波数分布発生部１５は、第１領域周波数分布の絶対値を、第２領域周波数分布の絶対値の最大値｜Ｓ２_ｂａｓｅ（０）｜で除することにより、第１領域周波数分布の絶対値を正規化する。また、周波数分布発生部１５は、第２領域周波数分布の絶対値を第２領域周波数分布の絶対値の最大値｜Ｓ２_ｂａｓｅ（０）｜で除することにより、第２領域周波数分布の絶対値を正規化する。

すなわち、正規化された第１領域周波数分布の絶対値（以下、第１正規化周波数分布と呼ぶ）｜Ｓｎ１_ｏｒｉｇ（ｆ）｜は、｜Ｓ１_ｏｒｉｇ（ｆ）｜／｜Ｓ２_ｂａｓｅ（０）｜として計算される。正規化された第２領域周波数分布の絶対値（以下、第２正規化周波数分布と呼ぶ）｜Ｓｎ２_ｂａｓｅ（ｆ）｜は、｜Ｓ２_ｂａｓｅ（ｆ）｜／｜Ｓ２_ｂａｓｅ（０）｜として計算される。周波数分布発生部１５は、正規化（ｎｏｒｍａｌｉｚｅ）係数αを計算する。正規化係数αは、例えば、以下に示す式のように、計算される。
α＝｜Ｓ２_ｂａｓｅ（０）｜／｜Ｓ１_ｏｒｉｇ（０）｜
周波数分布発生部１５は、第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜と第２正規化周波数分布｜Ｓｎ２_ｂａｓｅ（ｆ）｜とを、補正部１７に出力する。

図１２は、経時不変プロファイルの逆フーリエ変換の絶対値の正規化により発生された第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜と第２正規化周波数分布｜Ｓｎ２_ｂａｓｅ（ｆ）｜との一例を示す図である。図１２における第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜と第２正規化周波数分布｜Ｓｎ２_ｂａｓｅ（ｆ）｜とにおいて、ノイズの影響は実施形態１に比べて比較的小さいため、実施形態１におけるスムージングや関数近似のコストは、実施形態１に比べて小さい。説明を簡単にするため、本実施形態ではスムージングまたは関数近似の処理について言及していないが、｜Ｓ１_ｏｒｉｇ（ｆ）｜と｜Ｓ２_ｂａｓｅ（ｆ）｜とに対して、スムージングまたは関数近似が適宜実行されてもよい。

補正部１７は、第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜と第２正規化周波数分布｜Ｓｎ２_ｂａｓｅ（ｆ）｜との比により、第１領域におけるノイズレベルを第２領域におけるノイズレベルに合わせるための関数Ｆ_ｃｏｒ（ｆ）を決定する。具体的には、補正部１７は、第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜を、第２平滑スペクトルＰｎ＿ｅｓｔ_ｂａｓｅ（ｆ）に変換する（合わせるまたは一致させる）変換関数Ｆ_ｃｏｒ（ｆ）を決定する。

すなわち、補正部１７は、第２正規化周波数分布｜Ｓｎ２_ｂａｓｅ（ｆ）｜を第１正規化周波数分布｜Ｓｎ１_ｏｒｉｇ（ｆ）｜で除することで、変換関数Ｆ_ｃｏｒ（ｆ）を決定する。具体的には、以下に示す式により、変換関数Ｆ_ｃｏｒ（ｆ）は決定される。なお、変換関数Ｆ_ｃｏｒ（ｆ）の高周波側は、カットオフ周波数ｆｃ以下に滑らかに低減される。

Ｆ_ｃｏｒ（ｆ）＝｜Ｓｎ２_ｂａｓｅ（ｆ）｜／｜Ｓｎ１_ｏｒｉｇ（ｆ）｜
上式において、ゼロによる除算が生じないように工夫する必要がある。例えば、周波数分布発生部１５および図１２乃至図１３では、カットオフ周波数ｆｃがナイキスト周波数ｆｎより大きいため、｜Ｓｎ１_ｏｒｉｇ（ｆ）｜は、非ゼロとなる。

図１３は、第１正規化周波数分布を第２正規化周波数分布に変換する変換関数Ｆ_ｃｏｒ（ｆ）の一例を示す図である。図１３に示すｆｃは、カットオフ周波数である。

補正部１７は、変換関数Ｆ_ｃｏｒ（ｆ）を用いて、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、補正した第１周波数分布Ｓ_ｃｏｒ（ｆ）を画像発生部１９に出力する。具体的には、補正部１７は、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）に変換関数Ｆ_ｃｏｒ（ｆ）を乗ずることにより、第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。補正部１７は、例えば、以下に示す式に従って第１周波数分布Ｓ_ｏｒｉｇ（ｆ）を補正する。

Ｓ_ｃｏｒ（ｆ）＝Ｓ_ｏｒｉｇ（ｆ）×Ｆ_ｃｏｒ（ｆ）
図１４は、変換関数Ｆ_ｃｏｒ（ｆ）を用いて補正された第１周波数分布の絶対値｜Ｓ_ｃｏｒ（ｆ）｜の一例を示す図である。図１４に示すｆｃは、カットオフ周波数である。

（周波数分布補正機能）
図１５は、周波数分布補正処理に関する手順の一例を示すフローチャートである。
第１画像において、複数の第１エッジが検出される（ステップＳｂ１）。第２画像において、複数の第２エッジが検出される（ステップＳｂ２）。第１領域の対辺が第１エッジに直交し、経時的に変化の小さい領域が、第１画像に第１領域として設定される（ステップＳｂ３）。第２領域の対辺が第２エッジに直交し、経時的に変化の小さい領域が、第２画像に第２領域として設定される（ステップＳｂ３）。

第１領域における複数の画素値に基づいて、第１領域周波数分布の絶対値が発生される（ステップＳｂ４）。第２領域における複数の画素値に基づいて、第２領域周波数分布の絶対値が発生される（ステップＳｂ４）。第２領域周波数分布の絶対値の最大値を用いて、第１、第２領域周波数分布の絶対値各々が正規化される（ステップＳｂ５）。

第１正規化周波数分布を第２正規化周波数分布に変換する変換関数が決定される（ステップＳｂ７）。第１画像の全域における複数の画素値に基づいて、第１周波数分布が発生される（ステップＳｂ８）。変換関数を用いて、第１周波数分布が補正される（ステップＳｂ９）。補正された第１周波数分布に基づいて、第１画像を補正した補正第１画像が発生される（ステップＳｂ１０）。

以上に述べた構成によれば、以下の効果を得ることができる。
本実施形態に係る医用画像処理装置１によれば、非構造的であって経時的に不変な領域である第１領域と第２領域とを、第１画像と第２画像とにそれぞれ設定し、第１領域における第１領域周波数分布と第２領域における第２領域周波数分布とに基づいて第２画像のノイズレベルに第１画像のノイズレベルを合わせた補正第１画像を発生することができる。

本実施形態による周波数分布補正処理は、ノイズの特性と画像構造部分の特性との類推が困難となる場合、具体的には、例えば構造の複雑な部分を除外して、画素値が平坦な部分だけを強く平滑化するような非線形のアルゴリズムを使用している場合において、第１の実施形態に係る周波数分布補正処理に比べて有効となる。

なお、第１の実施形態に係る周波数分布補正処理と第２の実施形態に係る周波数分布補正処理とは、組み合わせて用いることができる。例えば、通常は、対象画像（第１画像）に対して過去に取得されたｂａｓｅ画像（第２画像）が従来の医用画像発生装置で生成されている場合、線形のアルゴリズムを用いて画像生成が実行されている可能性が大きい。また対象画像となりうる最近の装置による画像では、空気部分をマスクしている場合も多いため、第１の実施形態に係る周波数分布補正処理ではパワースペクトルの取得は困難となる。このような場合、第１の実施形態に係る周波数分布補正処理を用いて、ｂａｓｅ画像（第２画像）からパワースペクトルを取得し、第２の実施形態に係る周波数分布補正処理を用いて、対象画像（第１画像）からパワースペクトルを取得することで、第１周波数分布を補正することが可能である。

以上のことから、第１、第２の実施形態に係る医用画像処理装置１によれば、経時的な変化を過去画像と比較する場合に、装置依存の変動を除外するため、すなわち経時変化における撮像装置に由来する因子を除外するために、周波数分布を補正することができる。これにより、本医用画像処理装置１によれば、操作者は、患者本来の因子による解剖学的な経時変化の追跡を、高精度で実行することができる。

（第３の実施形態）
第１の実施形態と第２実施形態との相違は、第１画像および第２画像のパワースペクトル特性に関係するフィルター関数や補間処理に関する補間関数と、補間回数とに関する各種情報が既知であるならば、この各種情報を用いて、補正部１７で用いられる変換関数Ｆ_ｃｏｒ（ｆ）、第１変換関数ＰＮＲ１（ｆ）、または第２変換関数ＰＮＲ２（ｆ）などの補正関数を生成して、補正することにある。

記憶部２７は、第１画像および第２画像のパワースペクトル特性に関係するフィルター関数や補間処理に関する補間関数と、補間回数とに関する各種情報を記憶する。

補正部１７は、記憶部２７から上記各種情報を読み出して、上記各種情報に基づいて、変換関数Ｆ_ｃｏｒ（ｆ）、第１変換関数ＰＮＲ１（ｆ）、または第２変換関数ＰＮＲ２（ｆ）などの補正関数を決定する。例えば、補正部１７は、第１領域周波数分布および第２領域周波数分布に関するフィルター関数と補間関数と補間回数等とに基づいて補正関数を決定する。次いで、補正部１７は、決定した補正関数を用いて、第１周波数分布と第２周波数分布とのうち少なくとも一方等を補正する。

以上のことから、第３の実施形態に係る医用画像処理装置１によれば、第１、第２の実施形態に係る効果に加えて、以下に示す効果を更に有する。第３の実施形態に係る医用画像処理装置１は、第１、第２の実施形態に比べて、周波数分布発生部１５等における各種パワースペクトルの特性に係る各種計算（パワースペクトル特性測定）は、不要となり、より高精度に、周波数分布、パワースペクトル等を補正することができる。

加えて、実施形態に係る各機能は、当該処理を実行するプログラムをワークステーション等のコンピュータにインストールし、これらをメモリ上で展開することによっても実現することができる。このとき、コンピュータに当該手法を実行させることのできるプログラムは、磁気ディスク（フロッピー（登録商標）ディスク、ハードディスクなど）、光ディスク（ＣＤ−ＲＯＭ、ＤＶＤなど）、半導体メモリなどの記憶媒体に格納して頒布することも可能である。

以上、本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら新規な実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれるとともに、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。

１…医用画像処理装置、１０…エッジ検出部、１１…補間処理部、１３…領域設定部、１５…周波数分布発生部、１７…補正部、１９…画像発生部、２１…表示部、２３…インターフェース部、２５…入力部、２７…記憶部、２９…制御部。

Claims (9)

1. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記領域設定部は、
   前記第１画像の画素値に基づいて、前記第１領域を設定し、
   前記第２画像の画素値に基づいて、前記第２領域を設定すること、
   を特徴とする医用画像処理装置。
2. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記周波数分布発生部は、
   前記第１領域周波数分布に対応するパワースペクトルを平滑化することにより第１平滑スペクトルを発生し、
   前記第２領域周波数分布に対応するパワースペクトルを平滑化することにより第２平滑スペクトルを発生し、
   前記補正部は、
   前記第１平滑スペクトルを前記第２平滑スペクトルに変換する変換関数を決定し、
   前記変換関数を用いて、前記第１周波数分布を補正すること、
   を特徴とする医用画像処理装置。
3. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記周波数分布発生部は、
   前記第１領域周波数分布に対応するパワースペクトルを平滑化することにより第１平滑スペクトルを発生し、
   前記第２領域周波数分布に対応するパワースペクトルを平滑化することにより第２平滑スペクトルを発生し、
   前記補正部は、
   前記第２平滑スペクトルを前記第１平滑スペクトルに変換する変換関数を決定し、
   前記変換関数を用いて、前記第２周波数分布を補正すること、
   を特徴とする医用画像処理装置。
4. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記周波数分布発生部は、
   前記第１領域周波数分布に対応するパワースペクトルを平滑化することにより第１平滑スペクトルを発生し、
   前記第２領域周波数分布に対応するパワースペクトルを平滑化することにより第２平滑スペクトルを発生し、
   前記補正部は、
   前記第１平滑スペクトルを前記所定の周波数分布に対応するパワースペクトルに変換する第１変換関数を決定し、
   前記第２平滑スペクトルを前記所定の周波数分布に対応するパワースペクトルに変換する第２変換関数を決定し、
   前記第１変換関数を用いて前記第１周波数分布を補正し、
   前記第２変換関数を用いて前記第２周波数分布を補正すること、
   を特徴とする医用画像処理装置。
5. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記周波数分布発生部は、
   前記第１領域周波数分布に対応する第１パワースペクトルの最大値と前記第２領域周波数分布に対応する第２パワースペクトルの最大値とを用いて、前記第１パワースペクトルと前記第２パワースペクトルとを正規化し、
   前記補正部は、
   前記第１パワースペクトル、前記第２パワースペクトル、または前記所定の周波数分布に対応するパワースペクトルを基準として、前記第１周波数分布と前記第２周波数分布とのうち少なくとも一方を補正すること、
   を特徴とする医用画像処理装置。
6. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記第１画像において複数の第１エッジを検出し、前記第２画像において複数の第２エッジを検出するエッジ検出部を更に具備し、
   前記領域設定部は、
   前記第１領域における対向する２辺が前記第１エッジに直交し、かつ所定の閾値により規定される画素値範囲に対応する領域を、前記第１領域として設定し、
   前記第２領域における前記２辺が前記第２エッジに直交し、かつ前記画素値範囲に対応する領域を、前記第２領域として設定すること、
   を特徴とする医用画像処理装置。
7. 第１画像の一部分に対応する第１領域と、前記第１画像より過去に収集された第２画像の一部分に対応する第２領域とを設定する領域設定部と、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生する周波数分布発生部と、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正する補正部と、
   を具備し、
   前記補正部は、
   前記第１領域周波数分布および前記第２領域周波数分布に関するフィルター関数と補間関数と補間回数とに基づいて補正関数を決定し、決定した補正関数を用いて、前記第１周波数分布と前記第２周波数分布とのうち少なくとも一方を補正すること、
   を特徴とする医用画像処理装置。
8. 前記領域設定部は、
   前記第１画像における空間的に均一な領域または経時的に不変な領域を、前記第１領域として設定し、
   前記第２画像における空間的に均一な領域、または経時的に不変な領域を、前記第２領域として設定すること、
   を特徴とする請求項１乃至７のうちいずれか一項に記載の医用画像処理装置。
9. 第１画像の一部分に対応する第１領域を前記第１画像の画素値に基づいて設定し、
   前記第１画像より過去に収集された第２画像の一部分に対応する第２領域を前記第２画像の画素値に基づいて設定し、
   前記第１領域における複数の画素値に基づいて第１領域周波数分布を発生し、
   前記第２領域における複数の画素値に基づいて第２領域周波数分布を発生し、
   前記第１領域周波数分布、前記第２領域周波数分布、または所定の周波数分布を基準として、前記第１画像に対応する第１周波数分布と前記第２画像に対応する第２周波数分布とのうち少なくとも一方を補正すること、
   を具備することを特徴とする医用画像処理方法。