(57)【要約】【課題】システムおよび方法は、反射率および照明光Ｌを求める推定問題を解く。【解決手段】  本システムおよび方法は、汎関数の最小値が照明光Ｌの良好な推定値を与えると仮定されるような、未知の照明光Ｌの汎関数を解く。照明光Ｌの良好な推定値が得られることは、反射率Ｒの良好な推定値が得られることを意味する。汎関数は、変分フレームワークを用いて、最適解に対する要件を表現する。この要件は以下のことを含む。１）照明光Ｌが空間的に滑らかであること、２）反射率値が区間［０，１］にあり、したがって、画像Ｓを分解する場合、解は制約Ｌ＞を満たすべきであること、３）すべての可能な解のうち、得られるＲのコントラストが最大になるように、照明光Ｌの推定値が画像Ｓにできるだけ近くなるべきであること、そして４）反射率Ｒが一般的な自然画像挙動（例えば、反射率は区分的に滑らかである）に従うこと。
2005515515.tif 000002
【特許請求の範囲】
【請求項１】
画像向上方法であって、
画像をキャプチャすること（１２１）、
１つまたは複数の解像度レイヤを備える多重解像度構造を構成すること（１２３）、
少なくとも１回の反復を有する反復アルゴリズムを用いて各解像度レイヤを処理すること（１２５）、
それぞれの処理された解像度レイヤを次の解像度レイヤに射影すること（１４９）、
それぞれの射影された解像度レイヤを前記次の解像度レイヤに呼び出すこと（１２９）、および
前記射影された解像度レイヤを使用することであって、それによって、照明画像を推定する（１３０）使用すること
を含む方法。
【請求項２】
１回以上の反復のそれぞれについて、
ペナルティ汎関数の勾配を算出すること（１４３）、
最適直線探索ステップサイズを計算すること（１３２）、
をさらに含む請求項１に記載の画像向上方法。
【請求項３】
前記ペナルティ汎関数は、∂Ω上でｌ≧ｓおよび
【数１】
2005515515.tif 000020
という条件の下で次式：
【数２】
2005515515.tif 000021
で与えられ、
ここでΩは前記画像の台であり、∂Ωは画像境界であり、
【数３】
2005515515.tif 000022
は該画像境界の法線であり、αおよびβは自由な非負の実数である
請求項２に記載の画像向上方法。
【請求項４】
前記ペナルティ汎関数は次式：
【数４】
2005515515.tif 000023
で与えられ、ここでｗ1およびｗ2は前記勾配の非線形関数である
請求項２に記載の画像向上方法。
【請求項５】
前記反復アルゴリズムは、射影正規化最急降下アルゴリズムである
請求項１に記載の画像向上方法。
【請求項６】
前記反復アルゴリズムは、最急降下アルゴリズムである
請求項１に記載の画像向上方法。
【請求項７】
制約のセットが、前記照明光が前記画像の強度より大きい、Ｌ＞Ｓ、であること
を含む請求項１に記載の画像向上方法。
【請求項８】
ペナルティ項を適用することをさらに含み、該ペナルティ項は、
前記照明光が空間的に滑らかであること、
反射率が最大化されること、
前記反射率が区分的に滑らかであること
を含む請求項１に記載の画像向上方法。
【請求項９】
前記キャプチャされた画像および前記推定された照明画像に基づいて反射光画像を計算すること、
ガンマ補正係数を計算すること、
前記ガンマ補正係数を前記推定された照明画像に適用すること、および
前記ガンマ補正された照明画像と前記反射光画像を乗算することであって、それによって、補正画像を生成する、乗算すること
をさらに含む請求項１に記載の画像向上方法。
【請求項１０】
コンピュータ可読媒体に具現化された、ディジタル画像を向上させるシステム（１１２）であって、
入力ディジタル画像Ｓ（１０１）を受け取り該入力ディジタル画像の対数ｓ（１０３）を計算する対数モジュール（１０２）と、
前記入力ディジタル画像Ｓの照明光成分Ｌの推定値ｌ*（１０７）を生成する照明光推定器モジュール（１０４）であって、該推定器モジュール（１０４）は１つまたは複数の解像度レイヤを備える構成と、該１つまたは複数の解像度レイヤのそれぞれを処理する反復アルゴリズムとを使用する、照明光推定器モジュール（１０４）と、
前記対数ｓ（１０３）と、前記推定値ｌ*（１０７）に負号を付けたものとの和をとる加算ノード（１０８）であって、それによって、前記入力ディジタル画像Ｓの反射率成分Ｒの対数の推定値ｒ*（１０９）を生成する、加算ノード（１０８）と
を備え、
処理された解像度レイヤが使用され、それによって、次の解像度レイヤをスケールアップする
ディジタル画像を向上させるシステム。
【発明の詳細な説明】
【技術分野】
【０００１】
技術分野は、ディジタル画像の向上である。
【背景技術】
【０００２】
人間の視覚系は通常、場面間の明暗変化に適応し、非常に異なるダイナミックレンジを有する領域において細部を知覚することができる。
しかし、ほとんどの画像記録システムは、このダイナミックレンジ圧縮作業ができない。
その結果、これらの画像記録システムによって生成される画像は、人間の知覚によって生成される画像に比べて品質が悪いことが多い。
画像記録システムによる実行がうまく行かないことが多いもう１つの作業は、色恒常性に関するものである。
人間は場面照明光とはかなり独立に色を知覚するのに対して、画像記録システムはスペクトルシフトによって強く影響される。
【０００３】
上記の問題は、獲得された画像Ｓ、画像の物体の反射率Ｒ、および照明光Ｌの間の関係を、次式の画素ごとの乗算で記述することによって、数学的に述べることができる：
Ｓ＝Ｌ・Ｒ
この式は、画像Ｓ内の各点で、色値は反射率値と照明光値の積であることを意味している。
画像Ｓが与えられた場合、解くべき問題は、照明光の効果の除去および反射率画像の回復である。
すなわち、Ｓが与えられた場合、ＲおよびＬの両方を求めることである。
しかし、ＲおよびＬの両方の計算は通常不可能である。
そこで解決策は、ＲおよびＬを推定することができる方法を作り出すこととなる。
【０００４】
レティネックス(Retinex)理論は、画像における照明光の補償を扱う。
レティネックス理論の目標は、所与の画像Ｓを反射率画像Ｒと照明画像Ｌに分解し、画像ドメイン内の各点（ｘ，ｙ）でＳ（ｘ，ｙ）がＲ（ｘ，ｙ）とＬ（ｘ，ｙ）の積に等しくなるようにすることである。
このような分解の利益には、逆光／順光の照明光効果を除去し、屋内および屋外区域を含む画像のような空間変化する照明光を含む写真や他の画像のキャプチャ方法を向上させ、照明光により生じた色シフトを除去することにより画像の色を補正するという可能性が含まれる。
【０００５】
上記のように、所与の画像Ｓから照明光Ｌを回復することは数学的には不良設定問題であることが知られており、既知のアルゴリズムは、この限界を克服する方法および有効性においてさまざまである。
レティネックス手法は、１つのそのような方法のためのフレームワークを提供する。
レティネックス法は、人間視覚系についてのEdward Landの研究によって動機づけられたものである。
その研究については、非特許文献１（R. H. Land著「The Retinex Theory of Color Vision」Sci. Amer., Vol.237, pp.108-128 (1977)）に記載されている。
【０００６】
第１のレティネックスアルゴリズムは、ランダムウォーク型のものであった。
その後のレティネックスアルゴリズムは、ホモモルフィックフィルタリングを用いていた。
さらに別のレティネックスアルゴリズムのグループは、ポアソン方程式を解くことに基づいている。
【０００７】
図１は、従来技術のレティネックスアルゴリズムの一般的なクラスを表すアルゴリズム１０のブロック図である。
図１において、入力画像Ｓ１１が対数モジュール１２に入力され、入力画像Ｓおよびその２つの成分、照明光Ｌおよび反射率Ｒ、の対数バージョン１３が生成される。
すなわち、ｓ＝ｌｏｇＳ、ｌ＝ｌｏｇＬ、およびｒ＝ｌｏｇＲであり、したがってｓ＝ｌ＋ｒである。
入力画像の対数バージョン１３を用いて、推定器モジュール１４が、照明光の推定値１７（図１ではｌ*で示す）を計算する。
次に、推定値１７（ｌ*）は加算器１８で対数バージョン１３（ｓ）と組み合わされ、反射率の推定値１９（ｒ*で示す）が生成される。
最後に、推定値１９（ｒ*）は拡大器２０で対数からその対数に対応する数（真数）に変換され、反射率の推定値２１（Ｒ*で示す）として実数値が生成される。
従来技術のレティネックスアルゴリズムは通常、図１に示したのと同じプロセスを使用する。
【０００８】
レティネックス型のアルゴリズムはいくつかの異なる形態をとる。
１つのそのような形態はランダムウォークアルゴリズムである。
これは離散時間ランダム過程であり、「次の画素位置」が現在の画素位置の隣接点からランダムに選択されるものである。
ランダムウォーク型のレティネックスアルゴリズムは、以下の基本的定式化の変種である：多数のウォーカーが対数バージョン１３（ｓ）のランダム位置から出発する。
初期位置としてあるグレイ値を採用する。
ｓと同サイズのアキュムレータ画像Ａがゼロに初期化される。
ウォーカーが歩き回るにつれて、ウォーカーが訪れる各位置（ｘ，ｙ）にウォーカーの値を加えることによって、ウォーカーはＡを更新する。
照明画像は、アキュムレータ画像Ａを正規化することによって、すなわち、アキュムレータ画像Ａの各位置での値をその位置を訪れたウォーカーの個数で割ることによって得られる。
【０００９】
長いパスを有する多数のウォーカーを使用することによって、各位置でのアキュムレータ値がその位置の隣接点のガウス平均に収束することを容易に確かめることができる。
【００１０】
もう１つのタイプのレティネックスアルゴリズムは、ホモモルフィックフィルタリングを使用するものである。
このアルゴリズムでは、ローパスフィルタを用いてｓからｌを再構成する。
ホロモルフィック型レティネックスアルゴリズムは、反射率画像Ｒが画像内の先鋭な細部（すなわちエッジ）に対応するのに対して、照明光Ｌは空間的に滑らかであると予想されるという事実に基づいている。
そこでｌの妥当な推測は、ｌ*＝ＬＰ｛ｓ｝となる。
ここで、ＬＰはある広いガウシアン核との畳込みである。
したがって、ホモモルフィックフィルタリングを用いたレティネックスアルゴリズムは、実際には、単一の直接畳込みによるランダムウォークアルゴリズムと同じプロセスを適用する。
【００１１】
反射光Ｌは空間的に滑らかであると予想されるため、反射光Ｌの微分はゼロに近いはずである。
しかし、反射率Ｒは区分的に一定であるという仮定を考慮すると、反射率Ｒの微分はほとんど至る所で消えるが、画像のエッジに沿って高い値を有するはずである。
和ｓ＝ｌ＋ｒの微分をとり、高い微分ピークをクリッピングして落とすことは、クリッピング後の微分信号が照明光Ｌにのみ対応することを意味する。
モンドリアンワールドモデルに基づくポアソン方程式型レティネックスアルゴリズムは、区分的一定画像であるという反射率Ｒに対する上記の仮定を用いる。
ラプラシアンを作用させ、次のクリッピング操作：
τ（Δs）＝Δs  ここで｜Δs｜＜τ
その他の場合  τ（Δs）＝０
により、次のポアソン方程式が得られる：
Δｌ*＝τ（Δs）
このポアソン方程式の解法としては、実質的にラプラス演算子の逆をとる反復手続きを含むものがあり得る。
この方法への改良としては、ラプラシアンの代わりに画像勾配絶対値から不連続点を抜き出して、より良好な境界条件を提供するものがある。
【００１２】
ある解法では、ｌ≧ｓという制約を強制する「リセット」非線形性を有する反復アルゴリズムを使用する。
このアルゴリズムは、次の反復手順を実行する。
ｌ*n+1＝ｍａｘ｛（ｌ*n＋ｓ）／２，（ｌ*n＋Ｄn［ｌ*n］）／２｝
ここでＤnは並進演算子であり、画像ｓを、螺旋状に減衰する並進ベクトル列のｎ番目の要素分だけシフトする。
ｍａｘ演算子を取り除くと、次の簡略化バージョンが得られる。
ｌ*n+1＝（ｌ*n＋Ｄn［ｌ*n］）／２
上記の方程式は、画像を平滑化する簡単な平均化演算である。
非線形（ｍａｘ）演算は、照明画像ｌが制約ｌ*≧ｓを満たすように強制する。
【非特許文献１】R. H. Land著「The Retinex Theory of Color Vision」Sci. Amer., Vol.237, pp.108-128 (1977)
【発明の開示】
【発明が解決しようとする課題】
【００１３】
これらのアルゴリズムの使用にもかかわらず、現在の画像記録および画像向上のシステムおよび方法は、人間視覚系によって知覚される画像に匹敵するほどの十分な品質の画像を生成することができない。
したがって、より良好に照明光効果を除去し反射率画像を回復する改良された方法およびシステムが必要とされている。
【課題を解決するための手段】
【００１４】
［概要］
システムおよび方法は、反射率Ｒおよび照明光Ｌを求める推定問題を解く。
本システムおよび方法は、汎関数の最小値が照明光Ｌの良好な推定値を与えると仮定されるような、未知の照明光Ｌの汎関数を使用する。
照明光Ｌの良好な推定値が得られることは、反射率Ｒの良好な推定値が得られることを意味する。
【００１５】
汎関数は、変分フレームワークを用いて、最適解に対する要件を表現する。
この要件は以下のことを含む：１）照明光Ｌが空間的に滑らかであること、２）反射率値が区間［０，１］にあること、したがって、画像Ｓを分解する場合、解は制約Ｌ＞Ｓを満たすべきであること、３）すべての可能な解のうち、得られるＲのコントラストが極大になるように、照明光Ｌの推定値が画像Ｓにできるだけ近くなるべきであること、そして４）反射率Ｒが典型的な自然画像挙動（例えば、反射率は区分的に滑らかである）に従うこと。
【００１６】
４つの要件から、凸であり単一の解をもつ２次計画法の、振る舞いのよい最適化問題が得られる。
本方法およびシステムは、この解に対する数値アルゴリズムを使用する。
ここでアルゴリズムは、画像Ｓの多重解像度分解に基づく。
【００１７】
本方法およびシステムは、種々の画像状況に適応可能である。
第１の状況では、ＲＧＢ色空間における色画像が与えられた場合、アルゴリズムを３つの色チャネルに別々に適用することにより色補正が得られる。
代替実施形態では、入力画像Ｓをクロマ−ルミナンス色空間に変換することを含む。
この場合、アルゴリズムはルミナンスのみに適用されることにより、入力色を保存しつつ局所コントラストを改善した画像が得られる。
【００１８】
第２の状況は、照明光Ｌは大域的レベルで空間的に滑らかであると仮定されるが、この仮定は誤っている場合があることの原因である。
この状況は図９に劇的に示されている。
図９では、照明光Ｌは窓境界上に先鋭なエッジを有すると予想される。
図９の例によって示されるような状況では、照明光Ｌが先鋭なエッジを含み得ることは、大域的滑らかさを仮定するよりもむしろ、照明光Ｌに対して区分的滑らかさを強制することによって可能となるであろう。
【００１９】
第３の状況は、反射率Ｒに対する区分的滑らかさの仮定を使用しない。
その代わりに、より良好な品質の事前画像Ｓが使用される。
より良好な事前画像Ｓは、多重解像度および非等方性を考慮することを可能にする。
【００２０】
最後に第４の状況は、記録画像の自然な感じを人間のユーザに与えるためにはある程度の照明光効果が必要であることの原因である。
したがって、本方法およびシステムは、最終出力画像が、ある程度の照明光を含むがコントラストは改善されるように、照明光を部分的に戻すプロセスを含んでもよい。
【００２１】
上に列挙した４つの要件を用いて、変分フレームワークおよびアルゴリズムは、照明光を推定する問題に対する新規な解法を提供する。
変分フレームワークおよびアルゴリズムは、ディジタルカメラにおいて自動照明光補償モジュールとして画像再生品質を改善することを目的とする場合や、スキャナおよびプリンタにおいてスキャン／印刷画像の視覚品質を改善する特殊効果としての場合を含む多くの画像記録状況に適用可能である。
【００２２】
詳細な説明では以下の図面を参照する。
図面において、同じ番号は同じ対象を指す。
【発明を実施するための最良の形態】
【００２３】
［詳細な説明］
人間の視覚系は通常、場面（例えば、影、強い照明光源）間の明暗変化に適応し、非常に異なるダイナミックレンジを有する領域において細部を知覚することができる。
しかし、ほとんどの画像記録システム（例えばカメラ）は、このダイナミックレンジ圧縮作業ができない。
その結果、これらの画像記録システムによって生成される画像は、人間の知覚によって生成される画像に比べて品質が悪いことが多い。
画像記録システムによる実行がうまく行かないことが多いもう１つの作業は、色恒常性に関するものである。
人間は場面照明光とは関係なく色を知覚するのに対して、画像記録システムはスペクトルシフトによって強く影響される。
【００２４】
上記の問題は、獲得された画像Ｓ、画像の物体の反射率Ｒ、および照明光Ｌの間の関係を、次式の画素ごとの乗算で記載することによって、数学的に述べることができる。
すなわち、
Ｓ＝Ｌ・Ｒ
この式は、画像Ｓ内の各点で、色値は反射率値と照明光値の積であることを意味している。
この関係の効果は、簡単な実験によって見ることができる。
室内で、人物または他の物体が昼間に窓際（または他の光源の近く）に配置され、記録されるべき画像Ｓの一部を強い光源（すなわち、明るい日光）が照明するようにする。
したがって、照明光は窓越しに強く、室内で弱い。
次に写真を撮り、記録画像Ｓを反射率Ｒと照明光Ｌに分解する。
人物は、画像Ｓおよび照明光Ｌではほとんど見えないが、反射率Ｒでは見える。
【００２５】
画像Ｓが与えられた場合、解くべき問題は、照明光の効果の除去および反射率画像の回復である。
すなわち、Ｓが与えられた場合、ＲおよびＬの両方を求めることである。
しかし、ＲまたはＬのいずれかの計算は通常不可能である。
次に解決策は、ＲおよびＬを推定することができる方法を作り出すことである。
【００２６】
図２は、照明光Ｌおよび反射率Ｒを推定するために使用可能な改良されたアルゴリズム１００のブロック図である。
図２において、入力画像Ｓ１０１が対数モジュール１０２に入力され、入力画像Ｓ１０１およびその２つの成分、照明光Ｌおよび反射率Ｒ、の対数ｓ１０３が生成される。
すなわち、ｓ＝ｌｏｇＳ、ｌ＝ｌｏｇＬおよびｒ＝ｌｏｇＲであり、したがってｓ＝ｌ＋ｒである。
対数ｓ１０３を用いて、画像処理モジュール１０４は、多重解像度処理とともに射影正規最急降下(Projected Normal Steepest Descent)または類似のアルゴリズムを用いて、照明光の推定値１０７（図２ではｌ*で示す）を計算する。
次に、推定値１０７（ｌ*）は加算器１０８で出力１０３（ｓ）と組み合わされ、反射率の推定値１０９（ｒ*で示す）が生成される。
最後に、推定値１０９（ｒ*）は、指数モジュール１１０で対数からそれに対応する基本数値に変換され、反射率の推定値１１１（Ｒ*で示す）としての数値が生成される。
【００２７】
図３は、本明細書で説明するシステムおよび方法を使用する代表的なコンポーネント（デジタルカメラ）のブロック図である。
コンポーネント１１２は、画像をキャプチャする画像キャプチャデバイス１１４を含む。
画像キャプチャデバイスは、レンズ、ＣＣＤアレイおよび他のサブコンポーネントを含み得る。
画像コンバータ１１５は、キャプチャ画像の前処理に用いられる任意のサブコンポーネントであり得る。
画像プロセッサ１１６は、本明細書で説明する方法およびシステムに従って、キャプチャ画像を向上させるために必要な処理ソフトウェアおよびアルゴリズムを含み得る。
コンポーネント１１２は、画像プロセッサ１１６による処理の前後両方で、キャプチャ画像を格納する画像メモリ１１７を含み得る。
最後に、ディスプレイデバイス１１８は液晶ディスプレイデバイスでよく、例えばコンポーネント１１２のユーザに対して視覚ディスプレイを提供する。
【００２８】
ＲおよびＬを推定する問題に対する解に到達するため、方法およびその方法を適用するシステムは、ＳからＲおよびＬを求める推定問題を最適化問題として定式化することから始める。
すなわち、一方の未知量（ＬまたはＲのいずれか）の汎関数が提供され、その汎関数の最小値が所望の結果を与えるようにされる。
この定式化は、変分フレームワークを用いて、最適解に対する要件を表す。
フレームワークは以下の要件を具現する：
１．照明光Ｌは空間的に滑らかである。
２．反射率Ｒは単位区間に制限される。
この単位区間は制約Ｌ≧Ｓを加える。
対数関数は単調であるため、ｌ≧ｓである。
３．ｌ＝Ｃｏｎｓｔと設定する。
ここで、Ｃｏｎｓｔはｓの最大値より大きい任意の定数であり、前の２つの要件を満たす自明な解が得られる。
したがって、照明画像ｌが強度画像ｓに近くなるような要件が追加される。
すなわち、ｌの値は、ｄｉｓｔ（ｌ，ｓ）、例えばＬ２ノーム（ｌ−ｓ）2の形のペナルティ項を最小化する。
４．反射率画像ｓ＝ｌ−ｒは、高い事前確率を有すると仮定することができる。
自然画像に用いられる最も簡単な事前関数の１つは、空間的に滑らかな画像に高い確率を割り当てるものである。
５．照明光は、画像境界を越えて定数として滑らかに連続する。
これは、境界条件に対して要求される人為的仮定であるが、最終結果に対する効果は小さいであろう。
【００２９】
上記のすべての要件を１つの式にまとめると、次のペナルティ汎関数が得られる：
【００３０】
【数５】
2005515515.tif 000003
【００３１】
ここで、Ωは画像の台であり、∂Ωは画像境界であり、かつ
【数６】
2005515515.tif 000004
はその境界の法線である。
αおよびβは、自由な非負の実パラメータである。
汎関数Ｆ［ｌ］において、第１ペナルティ項（｜∇ｌ｜2）は、照明画像に対して空間的滑らかさを強制する。
∫（｜∇ｌ｜2）ｄｘｄｙを最小化することはオイラー・ラグランジュ（ＥＬ）方程式Δｌ＝０に帰着することに留意すれば、この滑らかさペナルティの選択は自然である。
第１ペナルティ項に対する最急降下解は、初期条件の変動が増大するとともにガウシアン平滑化演算となる。
前のセクションで言及したように、何人かの著者が照明光再構成のためのｓのガウシアン平滑化を提案している。
【００３２】
第２ペナルティ項（ｌ−ｓ）2は、ｌとｓの間の近接性を強制する。
これらの画像間の差は厳密にはｒであり、これはｒのノームが小さくなければならない（すなわち、Ｒが白くなる傾向がある）ことを意味する。
この項は自由パラメータαで重み付けられている。
この項の主な目的は、問題を正則化して、第２ペナルティ項の条件をより良くすることである。
なお、解ｌはｌ≧ｓであるように強制される。
【００３３】
第３ペナルティ項は、ベイジアン事前確率理論によって動機づけられたものであり、反射率画像ｒが「見た目の良い画像」であることを促進する。
実際、この事前確率式は、ｒを空間的に滑らかであるように強制する。
第３ペナルティ項は自由パラメータβで重み付けられている。
先鋭なエッジ、テクスチャ、１／ｆ挙動等を可能にする、より複雑なベイジアン式を用いてもよい。
この式が純粋に２次である限り、上記の最小化問題はかなり簡単なままに留まる。
【００３４】
上で定義された問題は２次計画法（ＱＰ）の形を有する。
最小化問題に対する必要十分条件は、オイラー・ラグランジュ方程式により次のように得られる：
【００３５】
【数７】
2005515515.tif 000005
【００３６】
この微分方程式は、ｌ＝ｓの場合は成り立つ必要はない。
【００３７】
最小化問題は、未知画像ｌに関するＱＰである。
多重分解能技法で加速化した射影正規化最急降下（ＰＮＳＤ）アルゴリズムを用いてこの最小化問題を解く。
【００３８】
ＰＮＳＤアルゴリズムは、汎関数Ｆ［ｌ］を最小化する正規化最急降下（ＮＳＤ）反復を適用し、続いて制約への射影をすることを必要とする。
１回のＮＳＤ反復は次の形式を有する。
すなわち、
ｌj＝ｌj-1−μNSDＧ
ここで、ｌjおよびｌj-1はそれぞれ、ステップｊおよびｊ−１における照明画像である。
ＧはＦ［ｌ］の勾配であり、またμNSDは最適直線探索ステップサイズである。
式（２）に対して、Ｆ［ｌ］の勾配は次式で与えられる。
すなわち、
Ｇ＝∇ｌj-1＋（α−βΔ）（ｌj-1−ｓ）
また、μNSDは次式で与えられる。
すなわち、
μNSD＝∫｜Ｇ｜2／∫（α｜Ｇ｜2＋（１＋β）｜∇Ｇ｜2）
部分積分を用いると、境界条件の自由度を除いて∫｜∇Ｇ｜2＝−∫ＧΔＧである。
【００３９】
代替手法として最急降下（ＳＤ）アルゴリズムがある。
この場合、μNSDは次のような一定値μSDで置き換えられる。
すなわち、
μSD∈（０，２／（λmax｛−（１＋β）Δ＋αＩ｝））
ここで、λmax｛Ａ｝は線形演算子Ａの最大固有値を表す。
この代替方法は、収束がわずかに遅くなることを代償にして計算量を節減する。
【００４０】
最後に、ｌj＝ｍａｘ（ｌj，ｓ）によって、制約ｌ≧ｓへの射影を行う。
実際には、Ｇは次のように計算することができる。
すなわち、
Ｇ＝ＧA＋α（ｌj-1−ｓ）−β（ＧA−ＧB）
ここで、
【００４１】
【数８】
2005515515.tif 000006
【００４２】
同様に、μNSDは次式で与えられる。
すなわち、
μNSD＝μA／（αμA＋（１＋β）μB）
ここで、
【００４３】
【数９】
2005515515.tif 000007
【００４４】
ラプラシアンは、次の核ｋLAPによる線形畳込みによって近似することができる。
【００４５】
【数１０】
2005515515.tif 000008
【００４６】
そして、積分は次のように和で近似される：
【００４７】
【数１１】
2005515515.tif 000009
【００４８】
ここで、Ｇ［ｍ，ｎ］＝Ｇ（ｍΔｘ，ｎΔｙ）である。
式（１）で与えられるような境界条件を取り入れるため、上記の畳込みは、画像Ｇの拡張バージョンに適用される。
拡張は、最初と最後の列および行を複製することによってなされる。
追加した行および列は畳込み後に除去される。
【００４９】
ＰＮＳＤアルゴリズムは、簡単ではあるが、通常は収束が遅い。
一般的な加速方式の代わりに、未知画像ｌは滑らかであると仮定されているという事実を用いてもよい。
具体的には、本方法は、まず粗い解像度画像ｌを推定し、補間により解像度画像ｌを拡張し、かつその結果を次の解像度レイヤの初期化として使用するという多重解像度アルゴリズムを適用する。
このようにすれば、収束のためには各解像度ごとにほんの数回の反復だけが必要である。
【００５０】
上記をまとめて、式（１）の解を求めるアルゴリズム１２０のブロック図を図４に示す。
アルゴリズム１２０は入力ブロック１２１から開始される。
ここでアルゴリズム１２０への入力は、サイズ［Ｎ，Ｍ］の画像ｓと、２個のパラメータαおよびβである。
【００５１】
初期化ブロック１２３で、画像ｓのガウシアンピラミッドが計算される。
こうして構成されたピラミッドは、ｐ個の解像度レイヤ（精細（１）から粗（ｐ）まで）を含み、現在の解像度レイヤ（ｋ）がｐに設定される。
ブロック１２５で、ＰＮＳＤアルゴリズムのＴk回の反復がｋ番目の解像度レイヤに適用され、これはすべての解像度レイヤがチェックされるまで（ブロック１２７）行われる。
ブロック１２９で、次の解像度レイヤが更新される。
すべての解像度レイヤが処理されると、その結果がアルゴリズム１２０の最終出力となる。
【００５２】
図５は、初期化ルーチン１２３をさらに詳細に示すブロック図である。
ブロック１３１で、カウンタがｌ＝１に初期化される。
ブロック１３２で、画像を次の核ｋPYRのような核で平滑化することによって、ガウシアンピラミッドが構成される：
【００５３】
【数１２】
2005515515.tif 000010
【００５４】
ブロック１３３で、画像は、例えば２：１の比でデシメーションされる。
プロセスは反復され（ブロック１３４）、カウンタがｐ回インクリメントされる（ブロック１３５）。
ここで、このセクションではｐ＜ｌg2（ｍｉｎ（Ｍ，Ｎ））である。
このプロセスは、従来「ｓのガウシアンピラミッド」と呼ばれている画像列｛ｓk｝pK=1を生成する。
画像ｓ1は原画像ｓであり、ｓpは、このガウシアンピラミッドに対する最も粗い解像度の画像である。
【００５５】
ブロック１３７で、プロセスはｋ＝ｐと設定し（すなわち、最も粗い解像度レイヤｋから出発する）、初期条件をｌ0＝ｍａｘ｛ｓp｝と設定する。
【００５６】
図６は、ＰＮＳＤ処理メインルーチン１２５のブロック図である。
ルーチン１２５はブロック１４１から開始される。
ここで、ｋ番目の解像度レイヤに対して、ルーチン１２５は次式を計算する。
すなわち、
ＧB＝Δkｓk
ここで、
【００５７】
【数１３】
2005515515.tif 000011
【００５８】
すなわち、ζと上記の核ｋLAPとの畳込みを因子２(2(k-l))で正規化したものである。
【００５９】
次に、ｊ＝１，...，Ｔkに対して、以下を行う。
すなわち、
【００６０】
ブロック１４３で、ルーチン１２５は勾配を計算する：
ＧA＝Δkｌj-1
Ｇ＝ＧA＋α（ｌｊ−１−ｓｋ）−β（ＧA−ＧB）
【００６１】
ブロック１４５で、ルーチン１２５はμNSDを計算する：
μA＝〈Ｇ，Ｇ〉
μB＝−〈Ｇ，ΔkＧ〉
μNSD＝μA／（αμA＋（１＋β）μB）
【００６２】
【数１４】
2005515515.tif 000012
【００６３】
ブロック１４７で、ルーチン１２５は次のようにＮＳＤ反復を完了する：
ｌj＝ｌj-1−μNSD・Ｇ
【００６４】
ブロック１４９で、その結果が制約に射影される：
ｌj＝ｍａｘ｛ｌj，ｓk｝
【００６５】
ループ処理は、次の中間問題を解く：
【００６６】
【数１５】
2005515515.tif 000013
【００６７】
図４に戻り、ブロック１２７で、ｋ＞１の場合、結果ｌTkは画素複製により（２：１の比で）新しいｌ0へとスケールアップされ、次の解像度ｋ−１レイヤに対する初期化がなされる。
解像度レイヤはｋ＝ｋ−１と更新され、アルゴリズムは再びブロック１２５に進むことにより続行される。
ｋ＝１の場合、結果ｌT1はアルゴリズム１２０の最終出力である。
【００６８】
パラメータをα＝β＝０と設定し、制約ｌ≧ｓを取り除くことにより、アルゴリズム１２０はホモモルフィックフィルタリングプロセスに帰着する。
これはランダムウォークアルゴリズムの基本的定式化と等価であることが示されている。
【００６９】
上記の方法およびシステムを用いると、式（１）に対する解、および数値アルゴリズム１２０の収束は一意であることを示すことができる。
すなわち、次式で与えられる変分最適化問題Ｐ、
【００７０】
【数１６】
2005515515.tif 000014
【００７１】
は、α≧０かつβ≧０の場合、一意解を有する。
【００７２】
ここまで、アルゴリズム１２０を単一チャネルに関して説明した。
別の実施形態では、本システムおよび方法は色画像に適用される。
色画像を処理する場合、１つの選択肢は、各色チャネルを別々に扱うことである。
このようなチャネルごとの処理は、「ＲＧＢレティネックス」処理と呼ぶことができる。
Ｒ、ＧおよびＢチャネルを別々に扱うと、通常色補正効果が得られる。
例えば、ＲＧＢレティネックス処理を赤みがかった画像に適用すると、赤の色相が除去され、その結果として得られる画像が明るく補正されるように、照明光が修正されるはずである。
したがって、一部の画像タイプに対しては、ＲＧＢレティネックス処理は色を改善する。
他の画像タイプに対しては、このような色補正は、色シフトを誇張するか、あるいは色飽和度を低下させる色アーティファクトを生じることがある。
【００７３】
さらに別の実施形態では、色は、ＨＳＶのようなルミナンス／クロミナンス色空間にマッピングされる。
次に、レティネックス補正が適用されるが、ルミナンス、すなわちＶレイヤにのみ適用される。
その後、色はＲＧＢドメインにマッピングされて戻される。
この方法は「ＨＳＶレティネックス」と呼ぶことができる。
このような場合のほうが色シフトは起こりにくい。
ＨＳＶレティネックス処理の利点は、単一チャネルのみが処理されることである。
しかし、この実施形態を用いると、色は色相（Ｈ）および飽和度
【数１７】
2005515515.tif 000015
に関しては補正されない。
【００７４】
レティネックスプロセスによって得られる反射率画像Ｒは、過度に向上された画像となることがある。
これは、
ｉ）人間の視覚系は通常、画像において照明光を好むこと、および
ｉｉ）すべての照明光の除去は、原画像の暗領域に存在するかもしれないノイズを露呈すること
という２つの事実によって説明することができる。
【００７５】
さらに別の実施形態では、再構成された照明光の補正バージョンが、構成された反射率画像に加え戻される。
図７はこの動作を記載している。
図７において、アルゴリズム２００は、強度画像Ｓ＝ｅｘｐ（ｓ）および反射率画像Ｒ＝Ｓ／Ｌから照明画像Ｌ＝ｅｘｐ（ｌ）を計算する。
【００７６】
これは、入力画像Ｓ２０１がレティネックスモジュール２０２に加えられ、照明画像Ｌ２０３が生成されることに示されている。
照明画像Ｌ２０３は、自由パラメータγによるガンマ補正動作２０４によって調整され、新たな照明画像Ｌ′２０９が得られる。
【００７７】
反射率モジュール２０６は、入力画像Ｓ２０１および照明画像Ｌ２０３から反射率Ｒ２０７を計算する。
新たな照明画像Ｌ′２０９が２１０で反射率Ｒ２０７に乗じられ、出力画像Ｓ′２１１が生成される。
ここで、Ｓ′＝Ｌ′Ｒである。ガンマ補正は次式により実行される。
Ｌ′＝Ｗ・［Ｌ／Ｗ］1/γ
ここで、Ｗは白色値（８ビット画像では２５５に等しい）である。
【００７８】
したがって、最終結果Ｓ′は次式により与えられる：
Ｓ′＝Ｌ′・Ｒ＝（Ｌ′／Ｌ）Ｓ
＝Ｗ（Ｌ／Ｗ）・Ｓ／Ｌ＝Ｓ／（Ｌ／Ｗ1-1/γ）
γ＝１の場合、照明光全体が加え戻され、したがってＳ′＝Ｓである。
γ＝∞の場合、照明光は全く戻されず、Ｓ′＝ＲＷであり、これは、アルゴリズム１００を区間［０，Ｗ］に引き伸ばした場合に得られるものと同じ反射率画像Ｒである。
この後者の場合は、反射率画像Ｒに最大一定照明光Ｗを乗じたものとみなすこともできる。
【００７９】
最終画像に照明光を加えることの類似物は、ホモモルフィックフィルタリング手法に見出すことができる。
対数ドメインでの照明光計算のための線形フィルタが、ｓのハイパス空間成分を除去するが、ローパス成分もγi分の１に減衰させる（ここでｉは照明光(illumination)を表す）。
Ｓ′に対する式は次の形に書くことができるため、ローパス成分は、γ＝γiとした照明光のガンマ補正の類似物である。
Ｓ′／Ｗ＝（Ｌ／Ｗ）1/r・Ｒ
したがって、
ｓ′−ｗ＝（１／γ）（ｌ−ｗ）＋ｒ
＝（１／γ）（ローパス成分）＋（ハイパス成分）
【００８０】
理想的環境では、照明光Ｌは区分的に滑らかであると仮定される。
すなわち、Ｌは、入力画像Ｓの強不連続点を除いて至る所で滑らかである。
この理想には、いくつかの方法で接近することができる。
第１に、再び図（１）を参照すると、照明光Ｌおよび反射率Ｒが両方とも滑らかであることが要求されている。
しかし、ｌとｒの和はｓに等しくなければならないため、これらの要件は互いに矛盾する。
したがって、ｌとｒは不連続点を共有する。
パラメータβはこの矛盾を調整し、各不連続点の適当な部分がｌとｒに分配されるようにする。
【００８１】
第２の手法は、照明光ｌの滑らかさ要件がｓの不連続点を除く至る所で成り立ち、反射光の滑らかさはｓの不連続点で成り立つようにすることである。
【００８２】
この第２の手法は、画像ｓのほとんどの部分で１に近い値をとり、強不連続点で０に近い値をとる重み関数ｗ（∇s）を計算することを必要とする。
その場合、式（１）は次のようになる：
【００８３】
【数１８】
2005515515.tif 000016
【００８４】
一連の勾配ステップから構成される式（３）に対する射影最急降下解は次式のようになる： 
【００８５】
【数１９】
2005515515.tif 000017
【００８６】
そして射影は、
ｌi＝ｍａｘ｛ｌi，ｓ｝
である。
ここで：
μは定数である。
Ｄxf（・）はｘ方向の前方微分である。
Ｄxb（・）はｘ方向の後方微分である。
Ｄyf（・）はｙ方向の前方微分である。
Ｄyb（・）はｙ方向の後方微分である。
【００８７】
図８は、汎関数（式（３））を解くために使用可能な、アルゴリズム１２０内のアルゴリズム１２５の代替形態であるアルゴリズム１２５ａを示している。
【００８８】
ブロック２２１で、アルゴリズム１２５ａはｗ（∇s）、αkおよびＴkを計算する：
【００８９】
【数２０】
2005515515.tif 000018
【００９０】
ここでＴｈは一定閾値である。
また、αk＝α0・２2k、Ｔk＝Ｔ0・２k-1である。
なお、αkは各解像度レベルごとに変化する。
大きいαkは粗い解像度レベルに対応する。
【００９１】
ブロック２２３で、アルゴリズム１２５ａは次式を計算する：
【００９２】
【数２１】
2005515515.tif 000019
【００９３】
ブロック２２５で、アルゴリズム１２５ａは次のようにＳＤ反復を完了する：
ｌi+1＝ｌi＋μ・Ｇ
【００９４】
ブロック２２７で、その結果が制約に射影される：
ｌj＝ｍａｘ｛ｌj，ｓｋ｝
【００９５】
ブロック２２９で、アルゴリズム１２５ａはｉ＝Ｔkかどうかをチェックし、そうである場合、アルゴリズム２２０による処理は終了する。
【図面の簡単な説明】
【００９６】
【図１】照明光を推定する従来技術のアルゴリズムの図である。
【図２】画像を分解する改良されたアルゴリズムのブロック図である。
【図３】図２の改良されたアルゴリズムを使用するシステムの図である。
【図４】図２のアルゴリズムによるプロセスを示すフローチャートである。
【図５】図４のプロセスのサブルーチンのフローチャートである。
【図６】図４のもう１つのサブルーチンのフローチャートである。
【図７】改良された代替アルゴリズムのブロック図である。
【図８】さらに別の改良された代替アルゴリズムのブロック図である。
【図９】画像に対する照明光の効果を示す画像の３つ組である。
【符号の説明】
【００９７】
１００・・・アルゴリズム、
１０２・・・対数モジュール、
１０４・・・画像処理モジュール、
１０８・・・加算器、
１１０・・・指数モジュール、
１１２・・・コンポーネント、
１１４・・・画像キャプチャデバイス、
１１５・・・画像コンバータ、
１１６・・・画像プロセッサ、
１１７・・・画像メモリ、
２０２・・・レティネックスモジュール、
２０５・・・ガンマ補正動作、
２０６・・・反射率モジュール、
【図１】
2005515515.tif 000024
【図２】
2005515515.tif 000025
【図３】
2005515515.tif 000026
【図４】
2005515515.tif 000027
【図５】
2005515515.tif 000028
【図６】
2005515515.tif 000029
【図７】
2005515515.tif 000030
【図８】
2005515515.tif 000031
【図９】
2005515515.tif 000032
【国際調査報告】
2005515515.tif 000033
2005515515.tif 000034
2005515515.tif 000035