JP4818212B2 - 画像処理装置及び方法 - Google Patents

Description

本発明は、Ｘ線画像等の画像の高周波成分を保持した状態で階調変換を行う画像処理装置及び方法に関するものである。

例えば、Ｘ線胸部画像は、Ｘ線が透過しやすい肺野の画像、及びＸ線が非常に透過しにくい縦隔部の画像より構成されるため、画素値の存在するレンジが非常に広い。このため、肺野及び縦隔部の両方を同時に観察することが可能なＸ線胸部画像を得ることは困難であるとされてきた。

そこで、この問題を回避する方法として、従来より次のような種々の方法が提案されている。
まず、非特許文献１に記載される方法がある。この方法は処理後の画素値Ｓ_D 、オリジナル画素値（入力画素値）Ｓ_ORG、オリジナル画像（入力画像）の低周波画像の画素値Ｓ_US、定数Ａ，Ｂ，Ｃ（例えばＡ＝３、Ｂ＝０．７）をもって、
Ｓ_D＝Ａ〔Ｓ_ORG−Ｓ_US＋Ｂ（Ｓ_US）〕＋Ｃ・・・（１）
なる式（１）で表わされるものである。

この方法は、高周波成分（第１項）、低周波成分（第２項）の重み付けを変えることが可能で、例えばＡ＝３、Ｂ＝０．７では高周波成分を強調し、かつ全体のダイナミックレンジを圧縮する効果が得られるものである。この方法は、処理無し画像と比較して診断に有効であるという評価が５人の放射線医により得られている。

また、特許文献１には、処理後の画素値Ｓ_Dオリジナル画素値（入力画素値）Ｓ_ORG、オリジナル画像（入力画像）のＹ方向プロファイルの平均プロファイルＰｙとＸ方向プロファイルの平均プロファイルＰｘをもって、
Ｓ_D＝Ｓ_ORG＋Ｆ［Ｇ（Ｐｘ，Ｐｙ）・・・（２）
なる式（２）で表される方法が記載されている。

ここで、関数ｆ（ｘ）が有する特性について説明する。まず、「ｘ＞Ｄｔｈ」ではｆ（０）が「０」となり、「０≦ｘ≦Ｄｔｈ」ではｆ（ｘ）が切片を「Ｅ」、傾き「Ｅ／Ｄｔｈ」として単調減少するものであり、（３）式で表されるものである。
Ｆ〔ｘ〕＝Ｅ−（Ｅ／ｔｈ）Ｘ ・・・（３）
Ｐｙ＝（ΣＰｙｉ）／ｎ・・・（４）
Ｐｘ＝（ΣＰｘｉ）／ｎ・・・（５）
但し、（ｉ＝１〜ｎ）、Ｐｙｉ、Ｐｘｉはプロファイル。そして例えば
Ｇ＝（Ｐｘ，Ｐｙ）＝ｍａｘ（ｐｘ，ｐｙ）・・・（６）
で表されるものである。この方法は、低周波画像の画素値Ｄｔｈ以下の濃度レンジが圧縮されるものである。

また、上記公報と同様な方法として、「自己補償デジタルフィルタ」と呼ばれる方法がある（非特許文献２）。この方法は、補償後（処理後）の画素値Ｓ_D、オリジナル画素値（入力画素値）Ｓ_ORG、オリジナル画像（入力画像）をマスクサイズＭ×Ｍ画素で移動平均をとった時の平均画素値Ｓ_US、図１６の単調減少関数ｆ（Ｘ）を用い、
Ｓ_D＝Ｓ_ORG＋ｆ（Ｓ_US）・・・（７）
Ｓ_US＝ΣＳ_ORG／Ｍ²・・・（８）
なる式（７）、（８）で表わされるものである。

ここで、関数ｆ（Ｓ_US）が有する特性について説明すると、まず、図１６に示す特性は、「Ｓ_US＞ＢＡＳＥ」ではｆ（Ｓ_US）が「０」となり、「０≦Ｓ_US≦ＢＡＳＥ」ではｆ（Ｓ_US）が切片を「しきい値ＢＡＳＥ」、傾き「ＳＬＯＰＥ」として単調減少するものである。従って、オリジナルの画素値Ｓ_ORGを濃度相当量として、上記（７）式を実行した際には、画像の平均濃度の低いところで濃度を持ち上げる、という画像に対する効果が得られる。

この方法は、上記（２）式と低周波画像の作成方法が異なり、（２）式では１次元データで低周波画像を作成していたのに対して、２次元データで低周波画像を作成するものである。この方法も低周波画像の画素値でＤｔｈ以下の濃度値を圧縮するものである。

また、特許文献２には、単調増加関数ｆ１（Ｘ）をもって、
Ｓ_D＝Ｓ_ORG＋ｆ１（Ｓ_US）・・・（９）
Ｓ_US＝ΣＳ_ORG／Ｍ²・・・（１０）
で表される方法が記載されている。

ここで、関数ｆ１（ｘ）が有する特性について説明する。まず、「ｘ＜Ｄｔｈ」ではｆ１（ｘ）が「０」となり、「Ｄｔｈ≦ｘ」ではｆ１（ｘ）が切片を「Ｅ」、傾き「Ｅ／Ｄｔｈ」として単調減少するものであり、（１１）式で表されるものである。
ｆ１〔ｘ〕＝Ｅ−（Ｅ／ｔｈ）Ｘ・・・（１１）

また、特許文献３には、一定濃度値以上の画像の高周波成分を強調する鮮鋭化方法が示されている。この方法は、超低周波数成分を強調すると共に、雑音の占める割合が大きい高周波成分を相対的に低減し、視覚的に見やすい画像が得られるようにすること及び偽画像を防止し、また雑音の増大を防止して診断性能を向上させることを目的とするものである。

定数ＢをＳ_ORG又はＳ_USの値の増大に応じて単調増加する変数とし、
Ｓ_D＝Ｓ_ORG＋Ｂ（Ｓ_ORG−Ｓ_US）・・・（１２）
Ｓ_US＝ΣＳ_ORG／Ｍ²・・・（１３）
なる式（１２）（１３）で表わされるものがある。上記（１２）式を実行した際には、画像の高周波成分を強調できるという効果がある。

特許第２５０９５０３号公報 特許第２６６３１８９号公報 特許第１５３０８３２号公報 ＳＰＩＥＶｏｌ．６２６．ＭｅｄｉｃｉｎｅＸＩＶ／ＰＡＣＳＩＶ（１９８６） 「自己補償デジタルフィルタ」（国立がんセンター、阿南氏他、日本放射線技術学会誌 第４５巻第８号１９８９年８月 １０３０頁）

しかしながら、上述したＳＰＩＥＶｏ１．６２６．ＭｅｄｉｃｉｎｅＸＩＶ／ＰＡＣＳＩＶ（１９８６）に記載される方法では、一定濃度範囲のダイナミックレンジを圧縮する思想がないので、画像全体のダイナミックレンジを均等に圧縮する。そのため、一定濃度範囲だけを圧縮することができず、そのため、例えば肺正面画像に本手法を用いた場合、縦隔部のみならず診断に有効な肺部の濃度レンジをも圧縮してしまうという問題があり、縦隔部のみを圧縮した場合よりも診断能が下がるという問題があった。

また、上記「自己補償デジタルフィルタ」による方法では、上記関数ｆ（Ｓ_US）の形状を、例えばＢＡＳＥまでを一定比率で減少する（線形である必要がある）ようにしないと、高周波成分に不自然な歪みが起きる問題がある。従って、高周波成分の振幅の大きさをオリジナル画像（入力画像）の高周波成分の振幅に保ったまま、階調の圧縮を非線形に自在に行えないという問題があった。

また、一般にダイナミックレンジ圧縮された画像は、ＣＲＴ表示、フィルム出力する場合に、再度階調変換する。上記「自己補償デジタルフィルタ」による方法等では、階調変換後の画像の高周波成分の振幅を調整する思想が無いため、ダイナミックレンジ圧縮した画像は、フィルム出力、画像表示する場合にさらに非線形な階調変換が行われる。そのため、階調変換曲線の傾きに依存し、高周波成分の振幅の大きさが変動する。従って、階調変換後の振幅が非線形に歪む問題や、高周波成分の振幅の大きさを保ってダイナミックレンジ圧縮を行っても、階調変換曲線の傾きが低いところでは高周波成分の振幅が小さくなり、有益情報が消失するという問題があった。また、エッジ部分でオーバーシュートやアンダーシュートがおきるという問題があった。

また、上記従来の高周波成分を強調する鮮鋭化方法では、高周波成分の足し込みの強さを自由に調整できるが、ダイナミックレンジを圧縮する思想がなく、濃度分布の広い画像を一枚のフィルムなどで表示できないという問題があった。

本発明は上記のような問題を解決するためになされたもので、画像の濃度分布の幅と高周波成分の振幅を自在に調節できるようにすることを目的とする。

本発明による画像処理装置の一つは、原画像を階調変換特性に基づいて階調変換する階調変換手段と、前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づいて、前記変動を補償した高周波成分を得る手段と、前記変動を補償した高周波成分と前記階調変換した原画像とを加算する加算手段と、を備えることを特徴とする。
本発明による画像処理装置の一つは、原画像を階調変換特性に基づいて階調変換する階調変換手段と、前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づき、前記変動を補償した高周波成分を前記階調変換後の原画像から得る手段と、前記変動を補償した高周波成分と前記階調変換後の原画像とを加算する加算手段と、を備えることを特徴とする。

本発明によれば、階調変換でダイナミックレンジを変更するためオーバシュートなどがでず、変動情報に基づき振幅を変更した高周波成分を加算することで診断情報も保存する効果を有する。

以下、本発明の実施の形態を図面と共に説明する。
まず、第１の実施の形態を説明する。
図２は、第１の実施の形態で用いる階調変換関数Ｆ１（）を示す。ここで、ｆ１（ｘ，ｙ）を２次元の入力原画像の濃度値とし、ｆ０（ｘ，ｙ）を２次元の階調変換後の出力画像（変換画像）の濃度値とし、ｘ，ｙを２次元上の座標とする。また、横軸を入力画像の濃度値ｆ１（ｘ，ｙ）、縦軸を出力画像（変換画像）の濃度値ｆ０（ｘ，ｙ）とする。この階調変換曲線では、入力濃度値２５００以下の傾き（ＳＬＯＰＥ）が０．２、入力濃度値２５００より上の傾きが１となっている。

図３において、実線が入力画像のプロファイル、ｆ１（Ｘ）、点線が入力画像の平滑化（低周波）画像ｆｕｓ（ｘ，ｙ）のプロファイルｆｕｓ（Ｘ）を示す。図４において、実線が入力画像を図２で示す階調変換曲線で階調変換した画像のプロファイルｆ０（Ｘ）、点線が階調変換された画像の平滑化（低周波）画像のプロファイルｆｕｓｏ（Ｘ）を示す。ここでＸは定数とする。
図５は、本実施の形態による画像処理方法の結果を示し、実線が入力画像のプロファイルｆ１（Ｘ）、点線が本実施の形態により処理された処理済み画像のプロファイルｆｄ（Ｘ）を示す。

次に動作について説明する。
まず、入力画像ｆ１（ｘ，ｙ）を図２で示す階調変換関数Ｆ１（）で（１４）式に示すように階調変換を行い、出力画像ｆ０（ｘ，ｙ）を得る。
ｆ０（ｘ，ｙ）＝Ｆ１（ｆ１（ｘ，ｙ））・・・（１４）

処理済み画像の画素値ｆｄ（ｘ，ｙ）を（５）式に従い得る。ここでＦ（ｘ，ｙ）は座標に依存した処理効果を表わす関数で、本実施の形態ではＦ（ｘ，ｙ）＝１とする。
ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋Ｆ（ｘ，ｙ）×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（１５）

ここで、ｆｕｓｏ（ｘ，ｙ）は出力画像（変換画像）ｆ０（ｘ，ｙ）の平滑化（低周波）画像、ｆｕｓ（ｘ，ｙ）は入力画像ｆ１（ｘ，ｙ）の平滑化（低周波）画像であり、例えば後述する（１６）〜（２０）式で求められる。平滑化には平均濃度を用いても、Ｅｒｏｓｉｏｎ、Ｄａｉ１ａｔｉｏｎ、０ｐｅｎｉｎｇ，Ｃｌｏｓｉｎｎｇ等のモルフォロジカルフィルタを用いてもよい。

本実施の形態で用いた階調変換曲線Ｆ１（）（図２）では、出力画像ｆ０（ｘ，ｙ）の濃度値２５００以下の高周波成分の振幅は２０％に圧縮され、濃度値２５００より高ければ、高周波成分の振幅は入力画像の振幅が保存される（図４実線）。

上記ｆｕｓ（ｘ，ｙ）の計算方法を（１６）〜（２０）式とする。ｆ１（ｘ，ｙ）を２次元の入力原画像とすると、
ｆ２（ｘ，ｙ）＝ｍｉｎ｛ｆ１（ｘ＋ｘ１，ｙ＋ｙ１）−Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（１６）
ｆ３（ｘ，ｙ）＝ｍａｘ｛ｆ２（ｘ＋ｘ１，ｙ＋ｙ１）＋Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（１７）
ｆ４（ｘ，ｙ）＝ｍａｘ｛ｆ３（ｘ＋ｘ１，ｙ＋ｙ１）＋Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（１８）
ｆｕｓ（ｘ，ｙ）＝ｍｉｎ｛ｆ４（ｘ＋ｘ１，ｙ＋ｙ１）−Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（１９）
で示される。

ここで、Ｄ（ｘ，ｙ）を円盤状フィルタ、ｒ１を任意の定数とし、入力画像に応じて選択される。
Ｄ（ｘ，ｙ）＝０、ｘ×ｘ＋ｙ×ｙ≦ｒ１×ｒ１
＝一∞、その他・・・（２０）

ここで得られたｆｕｓ（ｘ，ｙ）のプロファイルｆｕｓ（Ｘ）（図３点線）はエッジ構造を保存しているものであり、従来のダイナミックレンジ圧縮の欠点であるオーバーシュート、アンダーシュートが起きないものである。

同様に、ｆｕｓｏ（ｘ，ｙ）の計算方法を（２１）〜（２４）式とする。ある。ｆ０（ｘ，ｙ）を階調変換後の画像とする。
ｆ５（ｘ，ｙ）＝ｍｉｎ｛ｆ０（ｘ＋ｘ１，ｙ＋ｙ１）−Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（２１）
ｆ６（ｘ，ｙ）＝ｍａｘ｛ｆ５（ｘ＋ｘ１，ｙ＋ｙ１）＋Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（２２）
ｆ７（ｘ，ｙ）＝ｍａｘ｛ｆ６（ｘ＋ｘ１，ｙ＋ｙ１）＋Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙｌ×ｙ１≦ｒ１×ｒ１｝・・・（２３）
ｆｕｓｏ（ｘ，ｙ）＝ｍｉｎ｛ｆ８（ｘ＋ｘ１，ｙ＋ｙ１）−Ｄ（ｘ１，ｙ１）｜ｘ１×ｘ１＋ｙ１×ｙ１≦ｒ１×ｒ１｝・・・（２４）

ここで得られたｆｕｓｏ（ｘ，ｙ）のプロファイルｆｕｓｏ（Ｘ）（図４点線）はエッジ構造を保存しているものであり、ｆ０（Ｘ）との交点位置が、ｆｕｓ（Ｘ）とｆ１（Ｘ）との交点位置と同一点となるものである。

そして、図５の点線が得られた処理済み画像ｆｄ（ｘ，ｙ）のプロファイルｆｄ（Ｘ）である。濃度値２５００以下の濃度分布幅は入力画像の２０％に圧縮され、かつ高周波成分の振幅は入力画像の振幅を保持している。

図１は本実施の形態による画像処理装置の構成を示すブロック図であり、前記（１５）式を実現するものである。
図１において、入力画像ｆ１は階調変換手段１０１で図２の関数に基づいて階調変換され、変換画像ｆ０が得られる。この変換画像ｆ０は次に平滑手段１０２で平滑されて平滑化画像ｆｕｓｏが得られ、加算手段１０３に送られる。

一方、上記入力画像ｆ１は他の平滑手段１０４で平滑されて平滑化画像ｆｕｓとなる。次に減算手段１０５において、入力画像ｆ１から上記平滑化画像ｆｕｓが減算されることにより、高周波成分画像が得られる。この高周波成分画像は乗算手段１０６で定数が乗算された後、加算手段１０３において上記平滑化画像ｆｕｓｏと加算されることにより、処理済み画像ｆｄが得られる。

以上のように、第１の実施の形態によれば、入力画像の任意の階調領域の濃度分布幅を圧縮、伸張することができ、かつ階調変換後の高周波成分の振幅を自在に調整できる効果がある。

次に、第２の実施の形態を説明する。
図６は、第２の実施の形態による画像処理方法で用いる階調変換関数Ｆ１（）を示す。ここで、ｆ１（ｘ，ｙ）を２次元の入力原画像の濃度値とし、入力画像の平滑化（低周波）画像をｆｕｓ（ｘ，ｙ）とし、階調変換後の出力画像（変換画像）をｆｕｓ０（ｘ，ｙ）とする。

次に動作について説明する。
まず、入力画像ｆ１（ｘ，ｙ）の平滑画像ｆｕｓ（ｘ，ｙ）を例えば式（１６）〜（２０）で作成し、図６で示す階調変換関数Ｆ１（）で（２５）式で示すように階調変換を行い、出力画像ｆｕｓ０（ｘ，ｙ）を得る。
ｆｕｓ０（ｘ，ｙ）＝Ｆ１（ｆｕｓ（ｘ，ｙ））・・・（２５）

処理済み画像の画素値ｆｄ（ｘ，ｙ）を（２６）式に従い得る。ここでＦ（ｘ，ｙ）は座標に依存した処理効果を表わす関数とする。
ｆｄ（ｘ，ｙ）＝ｆｕｓ０（ｘ，ｙ）＋Ｆ（ｘ，ｙ）×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（２６）

図７は本実施の形態による画像処理装置の構成を示すもので、上記（２６）式を実現するものである。
図７において、入力画像ｆ１は平滑手段２０１で平滑されて平滑化画像ｆｕｓとなり、この平滑化画像ｆｕｓは減算手段２０４に送られると共に、階調変換手段２０２で図６の関数に基づいて階調変換され、変換画像ｆｕｓ０となり、加算手段２０３に送られる。

一方、入力画像ｆ１は減算手段２０４で上記平滑化画像ｆｕｓが減算されて、高周波成分画像が得られる。この高周波成分画像は乗算手段２０５で定数が乗算された後、加算手段２０３で上記出力画像ｆｕｓ０と加算されることにより、処理済み画像ｆｄが得られる。

以上のように、第２の実施の形態によれば、階調変換後に平滑化画像を作成する時間を省略でき、第１の実施の形態に比べて計算時間が早くなる効果がある。また、入力画像の任意の階調領域の濃度分布幅を圧縮、伸張することができ、かつ階調変換後の高周波成分の振幅を自在に調整できる効果がある。

次に、第３の実施の形態を説明する。
図８は、第３の実施の形態による画像処理方法で用いる階調変換関数Ｆ１（）を示す。ここで、ｆ１（ｘ，ｙ）を２次元の入力原画像の濃度値とし、ｆ０（ｘ，ｙ）を２次元の階調変換後の出力画像の濃度値とし、ｘ，ｙを２次元上の座標とする。また、横軸が入力画像の濃度値ｆ１（ｘ，ｙ）、縦軸が出力画像の濃度値ｆ０（ｘ，ｙ）を示す。

図９において、実線が入力画像のプロファイルｆ１（Ｘ）、点線が入力画像の平滑化（低周波）画像ｆｕｓ（ｘ，ｙ）のプロファイルｆｕｓ（Ｘ）、一点鎖線が処理画像のプロファイルｆｄ（Ｘ）である。

次に動作について説明する。
まず、入力画像ｆ１（ｘ，ｙ）を図８で示す階調変換関数Ｆ１（）で（２７）式で示すように階調変換を行い、出力画像（変換画像）ｆ０（ｘ，ｙ）を得る。
ｆ０（ｘ，ｙ）＝Ｆ１（ｆ１（ｘ，ｙ））・・・（２７）

処理済み画像の画素値ｆｄ（ｘ，ｙ）を（２８）、（２９）式に従い得る。ここでｃ（ｘ，ｙ）は階調変換率で（２８）式で定義される。
ｃ（ｘ，ｙ）＝∂Ｆ１（ｆ１（ｘ，ｙ））／∂ｆ１（ｘ，ｙ）・・・（２８）
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）＋（１−ｃ（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（２９）

ここで、ｆｕｓ（ｘ，ｙ）を入力画像ｆ１（ｘ，ｙ）の平滑化（低周波）画像であり、例えば式（３０）で示される。

尚、上記平滑化にはどのような手法を用いてもよく、例えばＥｒｏｓｉｏｎ、Ｄａｉ１ａｔｉｏｎ、０ｐｅｎｉｎｇ，Ｃｌｏｓｉｎｎｇ等のモルフォロジカルフィルタを用いてもよい。

図１０は本実施の形態による画像処理装置の構成示すもので、上記（２９）式を実現するものである。
図１０において、入力画像ｆ１は階調変換手段３０１で図８の階調変換曲線に基づいて階調変換されて変換画像ｆ０が得られ、加算手段３０２に送られる。一方、上記入力画像ｆ１は平滑手段３０３で平滑されて平滑化画像ｆｕｓとなる。減算手段３０４は入力画像ｆ１から平滑化画像ｆｕｓを減算して高周波画像が得られる。この高周波画像は加算手段３０２で上記変換画像ｆ０と加算されることにより、処理済み画像ｆｄが得られる。

以上のように、第３の実施の形態によれば、入力画像の任意の階調領域の濃度分布幅を圧縮、伸張することができ、かつ階調変換後の高周波成分の振幅を入力画像の高周波成分の振幅と同一にできる効果がある。また、平滑化処理を一回しか行う必要がなく、計算時間が短縮できる効果がある。さらに、平均濃度を用いる平滑化方法はモルフォロジカルフィルタ処理を計算時間が短縮できる効果がある。

次に、第４の実施の形態を説明する。
本実施の形態は、処理済み画像の画素値ｆｄ（ｘ，ｙ）、第１の画像（階調変換後の出力画像）ｆ０（ｘ，ｙ）の平滑化画像の画素値ｆｕｓｏ（ｘ，ｙ）、第２の画像（入力画像）の画素値ｆ１（ｘ，ｙ）、第２の画像の平滑化（低周波）画像の画素値ｆｕｓ（ｘ，ｙ）、処理効果を制御する関数Ｆ（）、画像上の座標ｘ，ｙを用い、
ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋Ｆ（ｆ１（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３１）
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）＋Ｆ（ｆ１（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３２）
なる式（２１）もしくは式（３２）なる演算式で表わすものである。

以上のように、第４の実施の形態によれば、処理効果を制御する関数Ｆ（）を第２の画像の濃度値ｆ１（ｘ，ｙ）に依存するようにするため、高周波成分の振幅を第２の画像の濃度値に応じて変更できる効果がある。

次に第５の実施の形態を説明する。
本実施の形態は、処理済み画像の画素値ｆｄ（ｘ，ｙ）、第１の画像ｆ０（ｘ，ｙ）の平滑化画像の画素値ｆｕｓｏ（ｘ，ｙ）、第２の画像の画素値ｆ１（ｘ，ｙ）、第２の画像の平滑化（低周波）画像の画素値ｆｕｓ（ｘ，ｙ）、処理効果を制御する関数Ｆ（）、画像上の座標ｘ，ｙを用い、
ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋Ｆ（ｆ０（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３３）
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）＋Ｆ（ｆ０（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３４）
なる式（３３）もしくは式（３４）なる演算式で表わすものである。

以上のように、第５の実施の形態によれば、処理効果を制御する関数Ｆ（）を第１の画像の濃度値ｆ０（ｘ，ｙ）に依存するようにするため、変換画像の高周波成分の振幅を第１の画像の濃度値に応じて変更できる効果がある。

次に、第６の実施の形態を説明する。
本実施の形態は、処理済み画像の画素値ｆｄ（ｘ，ｙ）、第１の画像ｆ０（ｘ，ｙ）の平滑化画像の画素値ｆｕｓｏ（ｘ，ｙ）、第２の画像の画素値ｆ１（ｘ，ｙ）、第２の画像の平滑化（低周波）画像の画素値ｆｕｓ（ｘ，ｙ）、処理効果を制御する関数Ｆ（）、画像上の座標ｘ，ｙを用い、
ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋Ｆ（ｆｕｓ（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３５）
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）＋Ｆ（ｆｕｓ（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３６）
なる式（３５）もしくは式（３６）なる演算式で表わすものである。

以上のように、第６の実施の形態によれば、処理効果を制御する関数（）を第２の画像の平滑化（低周波）画像の濃度値ｆｕｓ（ｘ，ｙ）に依存するようにするため、変換画像の高周波成分の振幅を第２の画像の平滑化（低周波）画像の濃度値に応じて変更できる効果がある。さらに、平滑化した画像の濃度に依存して高周波成分の振幅を調整するため、第２の画像の高周波成分の振幅に影響を受けない効果がある。

次に、第７の実施の形態を説明する。
本実施の形態は、処理済み画像の画素値ｆｄ（ｘ，ｙ）、第１の画像ｆ０（ｘ，ｙ）の平滑化画像の画素値ｆｕｓｏ（ｘ，ｙ）、第２の画像の画素値ｆ１（ｘ，ｙ）、第２の画像の平滑化（低周波）画像の画素ｆｕｓ（ｘ，ｙ）、処理効果を制御する関数Ｆ（）、画像上の座標ｘ，ｙを用い、
ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋Ｆ（ｆｕｓｏ（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３７）
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）＋Ｆ（ｆｕｓｏ（ｘ，ｙ））×（ｆ１（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ））・・・（３８）
なる式（３７）もしくは式（３８）なる演算式で表わすものである。

以上のように、第７の実施の形態によれば、処理効果を制御する関数Ｆ（）を第１の画像の平滑化（低周波）画像の濃度値ｆｕｓｏ（ｘ，ｙ）に依存するようにするため、変換画像の高周波成分の振幅を第１の画像の平滑化（低周波）画像の濃度値に応じて変更できる効果がある。さらに、平滑化した画像の濃度に依存して高周波成分の振幅を調整するため、第１の画像の高周波成分の振幅に影響を受けない効果がある。

次に、第８の実施の形態を説明する。本実施の形態では、第１の実施の形態で用いた図２、図４、図５を用いる。
図２においては、第１の実施の形態と同様に、ｆ１（ｘ，ｙ）を２次元の入力原画像の濃度値、ｆ０（ｘ，ｙ）を２次元の階調変換後の出力画像の濃度値とし、ｘ，ｙは２次元上の座標を示す。また、横軸が入力画像の濃度値ｆ１（ｘ，ｙ）、縦軸が出力画像の濃度値ｆ０（ｘ，ｙ）を示す。入力濃度値２５００以下の傾きが０．２、入力濃度値２５００より上の傾きが１となっている。

また、図４の実線が階調変換後の出力画像のプロファイルｆ０（Ｘ９）、点線が出力画像の平滑化（低周波）画像のプロファイルｆｕｓｏ（Ｘ）を示す。
また、図５の実線が入力画像のプロファイルｆ１（Ｘ）、点線が本実施の形態による画像処理方法の結果である処理済み画像のプロファイルｆｄ（Ｘ）を示す。

次に動作について説明する。
まず、入力画像ｆ１（ｘ，ｙ）を図２で示す階調変換関数Ｆ１（）で前記（４）式で示すように階調変換を行い、出力画像ｆ０（ｘ，ｙ）を得る。
次に処理済み画像の画素値ｆｄ（ｘ，ｙ）を（３９）式に従い得る。ここでｃ（ｘ，ｙ）は階調変換曲線の傾きを表わす関数で、前記（２８）式で示される。

ｆｄ（ｘ，ｙ）＝ｆｕｓｏ（ｘ，ｙ）＋ａ×（１／ｃ（ｘ，ｙ））×（ｆ０（ｘ，ｙ）−ｆｕｓｏ（ｘ，ｙ））・・・（３９）
ここで、ａは定数、ｆｕｓｏ（ｘ，ｙ）は出力画像ｆ０（ｘ，ｙ）の平滑化（低周波）画像の画素値であり、例えば前記（１６）〜（２０）式で示される。

本実施の形態で用いた階調変換曲線Ｆ１（）（図２）では、出力画像（変換画像）ｆ０（ｘ，ｙ）の濃度値２５００以下の高周波成分の振幅は２０％に圧縮され、濃度値２５００より高ければ、高周波成分の振幅は入力画像の振幅が保存される（図４実線）。
尚、第１の実施の形態と同様に、平滑化画像は例えば前記（３０）式による平均濃度を用いてもよく、上記モルフォロジカルフィルタを用いてもよい。

そして得られた処理画像ｆｄ（ｘ，ｙ）のプロファイルが図５の点線である。濃度値２５００以下の濃度分布幅は入力画像の２０％に圧縮され、かつ高周波成分の振幅は入力画像の振幅を保持している。

図１１は本実施の形態による画像処理装置の構成を示すもので、上記式（３９）式を実現するものである。
図１１において、入力画像ｆ１は階調変換手段８０１で階調変換されて変換画像ｆ０が得られ、減算手段８０５に送られると共に、平滑手段８０２で平滑されて平滑化画像ｆｕｓｏが得られる。この平滑化画像ｆｕｓｏは加算手段８０３に送られる。

一方、減算手段８０５は上記変換画像ｆ０から平滑化画像ｆｕｓｏを減算することにより、高周波成分が得られる。この高周波成分は乗算手段８０６で定数が乗算された後、加算手段８０３で上記平滑化画像ｆｕｓｏと加算されることにより、処理済み画像ｆｄが得られる。

以上のように、第８の実施の形態によれば、入力画像の任意の階調領域の濃度分布幅を圧縮、伸張することができ、かつ階調変換後の高周波成分の振幅を階調変換前の画像の高周波成分の振幅に保持できる効果がある。さらに、オーバーシュートやアンダーシュートが起きない効果がある。また、平滑化画像に濃度平均を用いた場合は計算時間を短縮できる効果がある。

図１２は第９の実施の形態による画像処理装置の構成を示す。
図１２において、９００は原画像としての入力画像を示し、９０１は原画像９００の階調を変換する階調変換手段、９０２は階調変換した変換画像、９０３は階調変換手段９０１で用いる階調変換曲線の微係数を記憶する微係数記憶手段、９０４は上記変換画像９０２の平滑化画像（低周波画像）９０５を作成する平滑化画像作成手段、９０６は変換画像９０２と平滑化画像９０５との差分を計算する高周波成分作成手段、９０７は高周波成分作成手段１０６で作成された高周波成分を、微係数記憶手段９０３に記憶した階調変換曲線の微係数に基づいて変換画像９０２に足し込む高周波成分足し込み手段である。

図１３は本実施の形態の処理の流れを示すフローチャートである。図１４、図１５は、階調変換手段９０１で用いる階調変換曲線を示し、横軸が入力画像の画素値、縦軸が出力画像の画素値を示す。図１４はＳ字型の階調変換曲線、図１５は入力濃度値Ｂ以下の傾きがＡ／Ｂ、入力濃度値Ｂより上の傾きが１となっている。

次に動作について図１３の処理の流れに従い説明する。
階調変換手段９０１は原画像９００を例えば図１４又は図１５に示す階調変換曲線に基づき（４０）式で示すように階調を変換する（ステップＳ２０１）。ここで、ｆ１（ｘ，ｙ）を２次元の入力原画像９００の濃度値とし、ｆ０（ｘ，ｙ）を２次元の階調変換後の変換画像９０２の濃度値とし、Ｆ１（）を階調変換曲線とする。ｘ，ｙは２次元上の座標を示す。
ｆ０（ｘ，ｙ）＝Ｆ１（ｆ１（ｘ，ｙ））・・・（４０）

そして微係数記憶手段９０３は（４１）式で示される、階調変換曲線の微係数を計算し、濃度値をテーブルｃ（ｘ）として記憶する（Ｓ２０２）。
ｃ（Ｆ１（ｘ））＝１−〔∂Ｆ１（ｘ）／∂ｘ〕・・・（４１）

次に、平滑化画像作成手段９０４では、前記（３０）式により画像９０２から平滑化画像９０５を計算する（Ｓ２０３）。
次に、高周波成分作成手段９０６は（４２）式で示すように、階調変換後画像９０２と平滑化画像９０５から高周波画像を計算する（Ｓ２０４）。ここでｆｈ（ｘ，ｙ）を高周波画像の画素値とする。
ｆｈ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）−ｆｕｓ（ｘ，ｙ）・・・（４２）

そして高周波成分足し込み手段９０７は、（４３）式で示すように階調変換後の変換画像９０２に、高周波成分作成手段９０６で計算された高周波成分を、微係数記憶手段１０３に記憶された微係数に基づき足し込んで、処理後画像ｆｄ（ｘ，ｙ）を得る（Ｓ２０５）。
ｆｄ（ｘ，ｙ）＝ｆ０（ｘ，ｙ）十ａ×ｃ（ｆ０（ｘ，ｙ））×ｆｈ（ｘ，ｙ）・・・（４３）
ここで、ａは定数である。

尚、平滑化画像９０５を前記（１６）〜（２０）式により、モルフォジ演算を用いて計算してもよい。
ここで得られたｆｕｓ（ｘ，ｙ）のプロファイルはエッジ構造を保存しているものであり、従来ダイナミックレンジ圧縮の欠点であるオーバーシュート、アンダーシュートが起きないものである。

本実施の形態によれば、入力画像の任意の階調領域の濃度分布幅を圧縮、伸張することができ、かつ階調変換後の高周波成分の振幅を階調変換前の画像の高周波成分の振幅に保持できる効果がある。また、オーバーシュートやアンダーシュートが起きない効果がある。さらに、平滑化画像に濃度平均を用いた場合は計算時間を短縮できる効果がある。

次に本発明による記憶媒体について説明する。
図１、図７、図１０、図１１、図１２等を含む各実施の形態によるシステムは、ハード的に構成してもよく、また、ＣＰＵやメモリ等からなるコンピュータシステムに構成してもよい。コンピュータシステムに構成する場合、上記メモリは本発明による記憶媒体を構成する。この記憶媒体媒体には、図１３を含む前述した各実施の形態の処理を実行するためのプログラムが記憶される。

また、この記憶媒体としては、ＲＯＭ、ＲＡＭ等の半導体メモリ、光ディスク、光磁気ディスク、磁気記憶媒体等を用いてよく、これらをＣＤ−ＲＯＭ、ＦＤ、磁気カード、磁気テープ、不揮発性メモリカード等に構成して用いてよい。

従って、この記憶媒体を上記各図を含む各実施の形態によるシステム以外の他のシステムあるいは装置で用い、そのシステムあるいはコンピュータがこの記憶媒体に格納されたプログラムコードを読み出し、実行することによっても、前述した各実施の形態と同等の機能を実現できると共に、同等の効果を得ることができ、本発明の目的を達成することができる。

また、コンピュータ上で稼働しているＯＳ等が処理の一部又は全部を行う場合、あるいは記憶媒体から読み出されたプログラムコードが、コンピュータに挿入された拡張機能ボードやコンピュータに接続された拡張機能ユニットに備わるメモリに書き込まれた後、そのプログラムコードの指示に基づいて、上記拡張機能ボードや拡張機能ユニットに備わるＣＰＵ等が処理の一部又は全部を行う場合にも、各実施の形態と同等の機能を実現できると共に、同等の効果を得ることができ、本発明の目的を達成することができる。

本発明の第１の実施の形態による画像処理装置のブロック図である。 第１、第８の実施の形態による階調変換関数Ｆ１（）を示す特性図である。 第１の実施の形による入力画像とモルフォロジカルフィルタによる平滑化画像のプロファイルを示す特性図である。 第１、第８の実施の形態による階調変換画像とモルフォロジカルフィルタによる平滑化画像のプロファイルを示す特性図である。 第１、第８の実施の形態による入力画像と処理画像のプロファイルを示す特性図である。 第２の実施の形態による階調変換関数Ｆ１（）を示す特性図である。 第２の実施の形態による画像処理装置のブロック図である。 第３の実施の形態による階調変換関数Ｆ１（）を示す特性図である。 第３の実施の形態による入力画像と入力画像の平滑化画像と処理画像のプロファイルを示す特性図である。 第３の実施の形態による画像処理装置のブロック図である。 第８の実施の形態による画像処理装置のブロック図である。 第９の実施の形態による画像処理装置のブロック図である。 第９の実施の形態による処理を示すフローチャートである。 第９の実施の形態による階調変換関数Ｆ１（）を示す特性図である。 第９の実施の形態による他の階調変換関数Ｆ１（）を示す特性図である。 従来のダイナミックレンジ圧縮に用いる単調減少関数を示す特性図である。

符号の説明

１０１、２０２、３０１、８０１、９０１ 階調変換手段
１０２、１０４、２０１、３０３、８０２、９０４ 平滑手段
１０３、２０３、３０２、８０３ 加算手段
１０５、２０４、３０４、８０５ 減算手段
１０６、２０５、３０５、８０６ 乗算手段
９００、ｆ１ 入力画像
９０２、ｆ０、ｆｕｓ０ 出力画像（変換画像）
９０５、ｆｕｓ、ｆｕｓｏ 平滑化画像
９０６ 高周波成分作成手段
９０７ 高周波成分足し込み手段
ｆｄ 処理済み画像

Claims (8)

1. 原画像を階調変換特性に基づいて階調変換する階調変換手段と、
   前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づいて、前記変動を補償した高周波成分を得る手段と、
   前記変動を補償した高周波成分と前記階調変換した原画像とを加算する加算手段と、を備えることを特徴とする画像処理装置。
2. 前記原画像の高周波成分は原画像のエッジ構造を保存するフィルタリング処理で得られたことを特徴とする請求項１に記載の画像処理装置。
3. 原画像を階調変換特性に基づいて変換する階調変換工程と、
   前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づいて、前記変動を補償した高周波成分と前記階調変換した原画像とを加算する加算工程と、を備えることを特徴とする画像処理方法。
4. 請求項３に記載の画像処理方法をコンピュータで実行するためのプログラム。
5. 原画像を階調変換特性に基づいて階調変換する階調変換手段と、
   前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づき、前記変動を補償した高周波成分を前記階調変換後の原画像から得る手段と、
   前記変動を補償した高周波成分と前記階調変換後の原画像とを加算する加算手段と、を備えることを特徴とする画像処理装置。
6. 前記階調変換後の原画像の高周波成分は前記階調変換後の原画像のエッジ構造を保存するフィルタリング処理で得られたことを特徴とする請求項５に記載の画像処理装置。
7. 原画像を階調変換特性に基づいて階調変換する階調変換工程と、
   前記階調変換によって変動する前記原画像の高周波成分の振幅の変動情報に基づき、前記階調変換後の原画像から得た前記変動を補償した高周波成分と前記階調変換後の原画像とを加算する加算工程と、を備えることを特徴とする画像処理方法。
8. 請求項７に記載の画像処理方法をコンピュータで実行するためのプログラム。

Images