(57)【特許請求の範囲】
【請求項１】
撮像部により取得された入力画像データに対してノイズ低減処理をする画像処理装置であって、
前記入力画像データに基づいて、前記入力画像データのヒストグラムを生成する生成手段と、
前記ヒストグラムを先鋭化することにより、回復ヒストグラムを生成する回復手段と、
前記撮像部により得られる画像データを構成する画素が取り得る画素値毎のノイズ特性を示すノイズ特性データを保持する保持手段と、
前記入力画像データにおける注目画素について、ノイズ特性データから得られる前記注目画素の画素値に対応するノイズ特性と前記回復ヒストグラムとを乗算することにより算出される確率分布の重心に相当する値を、前記注目画素の出力値とすることにより、前記入力画像データを補正する補正手段とを有することを特徴とする画像処理装置。
【請求項２】
撮像部により取得された入力画像データに対してノイズ低減処理をする画像処理装置であって、
前記入力画像データに基づいて、前記入力画像データにおける画素値をノイズ混入前の画素値の推定値と対応づけることにより、前記入力画像データにおける画素値を変換するためのルックアップテーブルを作成する作成手段と、
前記ルックアップテーブルを用いて前記入力画像データを構成する各画素の画素値を、該画素値に対応する推定値に変換することにより、前記入力画像データを補正する補正手段とを有し、
前記作成手段は、
前記入力画像データに基づいて、前記入力画像データのヒストグラムを生成する生成手段と、
前記ヒストグラムを先鋭化することにより、回復ヒストグラムを生成する回復手段と、
前記撮像部により得られる画像データを構成する画素が取り得る画素値毎のノイズ特性を示すノイズ特性データを保持する保持手段と、
前記回復ヒストグラムと画素値毎のノイズ特性データとを乗算することにより算出される確率分布の重心に相当する値に基づいて、画素値毎の推定値を決定し、前記ルックアップテーブルを作成するルックアップテーブル作成手段とを有することを特徴とする画像処理装置。
【請求項３】
前記ノイズ特性データは、理想的な画素値がノイズの影響によりどのような実際の画素値に変換されるかという確率を示すデータであることを特徴とする請求項１または２に記載の画像処理装置。
【請求項４】
前記ノイズ特性データは、前記撮像部により略一様な輝度を有する被写体を、露光量を変えて撮影し、各露光量に対応する複数の参照画像を取得し、
前記参照画像それぞれの平均値を前記理想的な画素値とみなし、前記参照画像ごとに、該参照画像を構成する画素の画素値の出現頻度を取得することにより、前記ノイズ特性データは生成されることを特徴とする請求項３に記載の画像処理装置。
【請求項５】
前記回復手段は、前記入力画像データの総画素数、前記ノイズ特性データ、及び前記ヒストグラムに基づいて生成される回復行列を用いて、前記ヒストグラムを先鋭化することを特徴とする請求項１乃至４のいずれか一項に記載の画像処理装置。
【請求項６】
前記重心に相当する値は、前記入力画像データにおける注目画素について、ノイズ特性データから得られる前記注目画素の画素値に対応するノイズ特性と前記回復ヒストグラムとを乗算することにより算出される確率分布の中央値であることを特徴とする請求項１乃至５のいずれか一項に記載の画像処理装置。
【請求項７】
前記重心に相当する値は、前記入力画像データにおける注目画素について、ノイズ特性データから得られる前記注目画素の画素値に対応するノイズ特性と前記回復ヒストグラムとを乗算することにより算出される確率分布の最頻値であることを特徴とする請求項１乃至５のいずれか一項に記載の画像処理装置。
【請求項８】
前記ルックアップテーブル作成手段は、前記撮像部により得られる画像データを構成する画素が取り得る画素値のうち代表的な画素値のみに対して前記推定値を算出し、前記代表的な画素値以外の画素値については、補間処理により推定値を取得することを特徴とする請求項２に記載の画像処理装置。
【請求項９】
撮像部により取得された入力画像データに対してノイズ低減処理をする画像処理方法であって、
前記入力画像データに基づいて、前記入力画像データのヒストグラムを生成し、前記ヒストグラムを先鋭化することにより、回復ヒストグラムを生成し、
前記撮像部により得られる画像データを構成する画素が取り得る画素値毎のノイズ特性を示すノイズ特性データを保持し、
前記入力画像データにおける注目画素について、ノイズ特性データから得られる前記注目画素の画素値に対応するノイズ特性と前記回復ヒストグラムとを乗算することにより算出される確率分布の重心に相当する値を、前記注目画素の出力値とすることにより、前記入力画像データを補正することを特徴とする画像処理方法。
【請求項１０】
撮像部により取得された入力画像データに対してノイズ低減処理をする画像処理方法であって、
前記入力画像データに基づいて、前記入力画像データにおける画素値をノイズ混入前の画素値の推定値と対応づけることにより、前記入力画像データにおける画素値を変換するためのルックアップテーブルを作成し、
前記ルックアップテーブルを用いて前記入力画像データを構成する各画素の画素値を、該画素値に対応する推定値に変換することにより、前記入力画像データを補正する補正手段とを有し、
前記ルックアップテーブルの作成は、
前記入力画像データに基づいて、前記入力画像データのヒストグラムを生成する生成工程と、
前記ヒストグラムを先鋭化することにより、回復ヒストグラムを生成する回復工程と、前記撮像部により得られる画像データを構成する画素が取り得る画素値毎のノイズ特性を示すノイズ特性データを保持する保持工程と、
前記回復ヒストグラムと画素値毎のノイズ特性データとを乗算することにより算出される確率分布の重心に相当する値に基づいて、画素値毎の推定値を決定する決定工程を含むことを特徴とする画像処理方法。
【請求項１１】
コンピュータに読み込ませ実行させることで、前記コンピュータを請求項１乃至８のいずれか一項に記載されている画像処理装置として機能させることを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明はノイズ低減処理に関する。
【背景技術】
【０００２】
通常、デジタルカメラやスキャナ等の撮像装置により撮影された画像にはノイズが含まれている。被写体像にノイズが混入すると、撮影画像の画質を低下させるという問題がある。
【０００３】
従来、様々なノイズ低減処理方法が提案されているが、それらは主として、平滑化フィルタによる処理を基礎にしている。
【０００４】
通常、被写体像は低周波成分に比べて高周波成分の強度が弱いため、画像の高周波領域ではノイズが支配的となる。そこで高周波成分を抑圧する平滑化フィルタを画像に適用すれば、被写体像の低〜中周波成分を保存しつつ、ノイズを低減し得る。
【０００５】
また、特許文献１では撮像画像の高周波成分をノイズらしくなるように補正して推定ノイズ画像を取得し、取得された推定ノイズ画像を撮像画像から減算することでノイズ低減処理を行なっている。
【先行技術文献】
【特許文献】
【０００６】
【特許文献１】特開２００６−３１０９９９
【発明の概要】
【発明が解決しようとする課題】
【０００７】
上述の平滑化フィルタをベースとするノイズ低減技術では、ノイズを低減するとともに被写体像の高周波成分も抑圧してしまう。その結果、被写体像の細かなテクスチャ成分が除去されてしまうという課題があった。
【０００８】
また、撮像画像中の個々の画素値のみから、注目画素が被写体のエッジであるかノイズによるものであるかを区別することができない。よって、特許文献１では、推定ノイズ画像の中に被写体像のエッジ成分が含まれることとなる。その結果、特許文献１では、撮像画像からエッジ成分が含まれる推定ノイズ画像を減算することとなり、ノイズ低減処理後の画像のエッジがぼけてしまうという課題があった。
【課題を解決するための手段】
【０００９】
上記課題を解決するために、本発明の画像処理装置は、撮像部により取得された入力画像データを入力する入力手段と、前記入力画像データに基づいて、前記入力画像データの入力ヒストグラムを生成する第一生成手段と、前記入力画像データの示す画像中の画素の特徴を決定する第一決定手段と、ノイズ特性データと前記入力ヒストグラムと前記画素の特徴とに基づいて、前記入力画像データを補正する補正手段と、を有する画像処理装置であって、前記ノイズ特性データは、注目画素の特性値の理想値が、複数の特性値の各々である確率に応じたデータであることを特徴とする。
【発明の効果】
【００１０】
本発明によれば、被写体像の細かなテクスチャ成分とエッジ成分を保存しつつ良好なノイズ低減を実行することができる。
【図面の簡単な説明】
【００１１】
【図１】実施例１における装置構成を表す図
【図２】実施例１における信号処理部１０３を表す図
【図３】実施例１におけるノイズ低減処理にかかる構成を表す図
【図４】実施例１におけるノイズ特性テーブルを表す図
【図５】ノイズの有無と画像の一次元ヒストグラムの関係を表す概念図
【図６】実施例１におけるヒストグラム回復の概念図
【図７】実施例１におけるＬＵＴ作成方法、及びＬＵＴを説明に用いる図
【図８】図４のノイズ特性テーブルにおいて画素値９０に対応する行を抜き出した図
【図９】グラデーション画像及びノイズを含めた該画像を表す図
【図１０】図９の画像データに対応する一次元ヒストグラムを表す図
【図１１】グラデーション画像が入力である場合のＬＵＴ作成方法の説明する図
【図１２】ノイズの有無とかかる画像の２次元ヒストグラムの関係を表す概念図
【図１３】実施例１におけるＬＵＴ作成方法を説明する図
【図１４】グラデーション中にエッジ部を含む画像を表す図
【図１５】図１４にノイズがのった画像が入力である場合のＬＵＴ作成方法の説明する図
【図１６】実施例３におけるノイズ低減処理にかかる構成を表す図
【図１７】実施例３における領域判定部のブロック図
【図１８】自然画像の多次元ヒストグラムの例
【図１９】自然画像の多次元ヒストグラムの例
【図２０】ヒストグラム生成のブロック図
【図２１】多次元ヒストグラム作成時の位置関係を表す図
【図２２】画素値変換部により画素値を割り当てる際の組合せを模式的に示した図
【図２３】統計ノイズを説明する図
【図２４】ノイズ特性テーブルの作成方法のフローを示す図
【発明を実施するための形態】
【００１２】
［実施例１］
（画像処理装置の概要）
図１は、本実施例における画像処理装置を示す図である。
撮像部１０１は、ズーム／フォーカス／ぶれ補正／レンズ、絞り、シャッター、光学ＬＰＦ、ＩＲカットフィルタ、カラーフィルタ、及びＣＭＯＳやＣＣＤなどのセンサなどから構成され、被写体から入射される光量を取得（撮像）する。
Ａ／Ｄ変換部１０２は、センサが取得した被写体から入射される光量をデジタル値に変換する。
【００１３】
信号処理部１０３は、デジタル値に対して、ノイズ低減処理、ホワイトバランス処理、エッジ強調処理、色変換処理、ガンマ処理を行い、デジタル画像データを生成する。各処理の順序は上述の順序に限らない。
Ｄ／Ａ変換部１０４は、デジタル画像データに対しアナログ変換を行い、表示部１１３に出力する。
【００１４】
エンコーダ部１０５は、デジタル画像データをＪＰＥＧやＭＰＥＧなどのファイルフォーマットの画像データに変換する処理を行う。
メディアインターフェース１０６は、ＰＣその他メディア（例えば、ハードディスク、メモリーカード、ＣＦカード、ＳＤカード、ＵＳＢメモリ）につなぐためのインターフェースである。
【００１５】
ＣＰＵ１０７は、画像処理装置内の各構成の処理に関わり、ＲＯＭ１０８やＲＡＭ１０９に格納された命令を順に読み込み、本実施例の各種処理を実行する。また、ＲＯＭ１０８とＲＡＭ１０９は、その処理に必要なプログラム、データ、作業領域などをＣＰＵ１０７に提供する。
【００１６】
撮像系制御部１１０は、フォーカスの合焦、シャッターの開閉、絞り調節などの、ＣＰＵ１０７から指示された撮像部１０１の制御を行う。
操作部１１１は、タッチパネルやボタン、モードダイヤルなどにより構成され、これらを介して入力されたユーザ指示を受け取る。
【００１７】
キャラクタージェネレーション部１１２は、表示などのための文字やグラフィックなどを生成する。
【００１８】
表示部１１３は、一般的には液晶ディスプレイが広く用いられており、キャラクタージェネレーション部１１２やＤ／Ａ変換部１０４から受け取った画像データで示される画像や文字、グラフィックの表示を行う。また、タッチスクリーン機能を有していても良く、その場合は、ユーザ指示を操作部１１１の入力として扱うことも可能である。
【００１９】
なお、本実施例における装置の構成要素は上記の画像処理装置以外にも考えられる。例えば、デジタルカメラ等の撮像装置で撮影されたデジタル画像データを、ネットワークやメディアを介してコンピュータが取得し、コンピュータにあるアプリケーションが本実施例の一連の処理を実施しても良い。
【００２０】
本実施例のノイズ低減処理は、信号処理部１０３において実行されることが好適である。信号処理部１０３の詳細図を図２に示す。図２に示す信号処理部１０３は入力されたデジタル画像データに対して画質向上のための様々な処理を行う。まず、ノイズ低減処理部２０１は、以下に詳述するノイズ低減処理を実行する。ノイズ低減処理された画像データは、ホワイトバランス制御部２０２により画像のホワイトバランスが調整される。その後、エッジ強調部２０３により画像データにより示される画像の鮮鋭性を向上させるためのエッジ強調処理がなされる。そして、色変換部２０４による色再現性を向上させるための色変換処理や、ガンマ処理部２０５による画像データに対するガンマ処理がなされる。これら一連の信号処理部１０３によって処理された結果はＤ／Ａ変換部１０４やエンコーダ部１０５に送られる。なお、ノイズ低減処理部２０１によるノイズ低減処理は、必ずしもホワイトバランス制御部２０２の前段で行なう必要はない。例えば、色変換部２０４による色変換処理やガンマ処理部２０５によるガンマ処理の後段で、ノイズ低減処理を行なってもよい。ただし、色変換やガンマ処理によって実質的なノイズ特性が変わるため、ノイズ低減処理の前段の処理に応じて後述するノイズ特性テーブルも補正する必要があることは言うまでもない。例えば、入力画像の画素値を２倍する処理がノイズ低減処理の前に存在する場合、ノイズ低減処理の処理対象画像のノイズ量も２倍に増えるため、ノイズ特性テーブルのノイズ量を２倍する補正などである。
【００２１】
（一次元ヒストグラムを用いた場合のノイズ低減方法）
多次元ヒストグラムを用いたノイズ低減方法の原理の説明のため、まず一次元ヒストグラムを用いたノイズ低減方法について以下に述べる。
図３（ａ）に一次元ヒストグラムを用いたノイズ低減処理に必要なブロック構成図を示す。図３（ａ）の各ブロックの動作を説明する。まず、入力部３０１は、入力画像データをＡ／Ｄ変換部１０２から取得する。また、入力部３０１は、入力画像データが示す入力画像の総画素数およびノイズ特性データをＲＯＭ１０８あるいはＲＡＭ１０９から取得する。なお、入力画像の総画素数は、入力された入力画像データを解析することにより取得しても良い。
【００２２】
図４（ａ）は、ＲＯＭ１０８あるいはＲＡＭ１０９に保持されているノイズ特性データ（又はノイズ特性テーブル）の一例である。
図４（ａ）のノイズ特性テーブルは、ある画素値（理想値）が画像処理装置に入力された際に、その画素値がノイズの影響によりどのようなデジタル値（実際値）に変換されるかという確率ヒストグラムを示すテーブルである。
図４（ａ）の例では、入力された画素値５８に対して、画像中の画素の１０％（０．１０）は画素値が５８となる。一方で、５％の画素は画素値が５６となり、９％の画素は画素値が５７となり、９％の画素は画素値５９となる。
【００２３】
図４（ａ）に示すノイズ特性テーブルの作成方法を説明する。図２４はノイズ特性テーブルの作成フローを示している。ステップＳ４０１にて、略一様な輝度を有する被写体を、露光量を変えて撮影し、各露光量に対応した複数枚の画像（参照画像）を取得する。ステップＳ４０２にて、各参照画像の平均画素値を取得する。この平均画素値を前述の中央画素値とみなす。なお、平均画素値（中央画素値）は、各参照画像の画素値の最頻値であっても良い。ステップＳ４０３にて、参照画像ごとに参照画像中の画素値の出現頻度と平均画素値との関係を示すヒストグラムを生成する。そして、ステップＳ４０４にて、参照画像の総画素数で、各ノイズヒストグラムを除算することにより正規化ノイズヒストグラムを生成する。中央画素値がＭである参照画像に対応する正規化ノイズヒストグラムをＨＭとする。ステップＳ４０５にて、正規化ノイズヒストグラムＨＭを列ベクトルとして表し、列ベクトルＨ０から列ベクトルＨＭｍａｘ（Ｍｍａｘは画像処理装置がとりうる画素値の最大値）までを横に並べて図４（ａ）のノイズ特性テーブルを得る。ノイズ特性テーブル中の各要素（画像中の各画素値の画素数を総画素数Ｎで除算したもの）を正規化画素頻度と呼ぶ。なお、略一様な輝度を有する被写体を撮影することにより、撮像光学系の解像度特性の影響を低減した状態でノイズ特性が得られる。また、周辺光量落ちのあるレンズでは、画像中の小領域でノイズ特性を計算するなどが考えられる。
【００２４】
なおノイズ特性テーブル中の各要素は、正規化画素頻度ではなく、画像中の各画素値の画素数であってもよい。
また、代表的な中央画素値Ｍに対する正規化ノイズヒストグラムＨＭを複数取得し、取得された正規化ノイズヒストグラムＨＭに対して補間処理を行なうことによって残りの中央画素値Ｍに対する正規化ノイズヒストグラムＨＭを計算しても良い。
【００２５】
更に、ノイズ特性データは、図４（ａ）のようなテーブル形式としたが、画像処理装置における入力された理想値に対する実際値の確率ヒストグラムが生成できるものであればこの限りではない。例えば、ノイズ特性データを関数で作成し、入力された理想値に対して逐次実際値の確率ヒストグラムを生成するとしても良い。
【００２６】
また、複数の反射率の異なる領域を有するチャートを撮影した画像からノイズ特性データを生成しても良い。この場合は、領域毎に中央画素値Ｍに対する正規化ノイズヒストグラムＨＭを取得することにより、少ない撮影回数でノイズ特性データを生成することができる。
【００２７】
なお、ノイズ特性データは、テーブル形式で記述されるとしたが、他の形式で記述しても良い。例えば、中央画素値（理想値）ごとに標準偏差等の分布を示す統計値を保持しておくことにより、画素値（実際値）を計算に求めるとしても良い。
【００２８】
図５は、画像処理装置におけるある画素値の正規化ノイズヒストグラムの分布を示している。
輝度が略一様な被写体を撮影した場合、ノイズが少なければ図５（ａ）のような分散の少ない画像値の正規化ノイズヒストグラムの分布が得られる。一方で、ノイズが画像に混入した場合、図５（ｂ）のような分散の大きい画素値の正規化ノイズヒストグラムの分布が得られる。この広がりがノイズの特性を表している。実際に得られる画像は図５（ｂ）のようになることが多い。
【００２９】
そこで、画像処理装置に入力された中央画素値６０１に対応する画素にノイズが混入することにより、中央画素値６０１の周囲に画素値（実際値）が広がって分布されると仮定する。前述のとおり、輝度によってノイズ特性が異なるため、様々な輝度で撮影を行い、複数の中央画素値６０１に対する正規化ノイズヒストグラムを並べて図４（ａ）のノイズ特性テーブルを作成する。
【００３０】
また、図４（ａ）のノイズ特性テーブルの他の作成方法として、まずノイズの少ないカメラ（参照カメラ）を用いて撮影した画像をノイズがない画像と見なす。この場合、上述のような輝度が略一様な被写体を撮影する必要はない。例えば、屋内の被写体が固定されているものであれば良い。そして前記参照カメラで撮影した被写体を、ノイズ低減処理を行いたいカメラ（対象カメラ）で撮影する。参照カメラで取得した画像データの平均値を中央画素値に設定し、図２４のＳ４０２の平均画素値（中央値）とみなす。そして、対象カメラにより取得した画像データを図２４のＳ４０１の一つの参照画像とみなす。この処理を参照カメラと対象カメラの露光量を変更しながら複数回行い、複数の中央画素値および参照画像を取得する。あとの処理は、前述のノイズ特性テーブルの作成処理と同一である。
【００３１】
また、同一のカメラにより図４（ａ）のノイズ特性テーブルを作成しても良い。撮影感度を調整し、ノイズが少ない画像を取得した後、同一の被写体についてノイズ低減処理が必要な撮影感度で画像を取得する。これら二つの差分画像を参照画像とみなし、同様にノイズ特性テーブルを作成しても良い。
ノイズ特性は、露光量のみならず、温度に依存する。よって、ノイズ低減処理の精度を高めるために、図４（ａ）に示すノイズ特性テーブルを温度に応じて作成していても良い。
【００３２】
次に、一次元ヒストグラム生成部３０４は、入力画像データの画素値に対する頻度を総画素数Ｎで除算して一次元ヒストグラムを生成する。なお除算は、ビットシフトなどの演算コストの低い除算に相当する処理でもよい。生成された一次元ヒストグラムは、入力画像データ中の画素値がＸである画素の画素数を第Ｘ要素とするベクトルにより表す。
【００３３】
一次元ヒストグラム生成部３０４の詳細を図２０（ａ）に示す。頻度計算部はＡ／Ｄ変換部からの入力画像データを得て、入力画像データに含まれる各画素値の頻度をカウントし、カウントされた頻度を入力画像データの総画素数で除算して一次元ヒストグラムを出力する。なお、出力される一次元ヒストグラムは、総画素数で除算することなく各画素値の画素数の頻度からなるヒストグラムであってもよい。
【００３４】
次に、回復行列生成部３０６により、入力画像データの総画素数Ｎ、ノイズ特性テーブル、及び一次元ヒストグラムから回復行列を生成し、回復行列保持バッファ３０７に保持する。回復行列の生成の詳細は後述する。回復行列を一次元ヒストグラムに適用することにより、ノイズ混入により分散したヒストグラムを、ノイズ混入前の頻度が鮮鋭化した状態に近づける。「鮮鋭化した状態」とは、ヒストグラムをフーリエ変換（周波数変換）した際に、高周波成分が大きくなった状態を言う。以降では、この作用を回復（又はヒストグラム回復）と呼ぶ。
【００３５】
次にヒストグラム回復部３０８において、一次元ヒストグラムに回復行列を適用する事により、回復ヒストグラムを算出する。図６（ａ）にヒストグラム回復の概念図を示す。ヒストグラム回復により、入力画像データの一次元ヒストグラム８０７が、より鮮鋭化し、回復ヒストグラムへと変換される。鮮鋭化された回復ヒストグラムの高周波成分は、入力画像データの一次元ヒストグラムの高周波成分よりも大きい。図６（ｂ）は、入力画像データの一次元ヒストグラムの周波数成分、および先鋭化された回復ヒストグラムの周波数成分を示している。図中、高周波成分に注目すれば明らかなように、一次元ヒストグラムの高周波成分は、回復ヒストグラムの高周波成分よりも大きい。
次に、ＬＵＴ作成部３１０によりノイズ特性データ、及び回復ヒストグラムからＬＵＴ（ルックアップテーブル）を作成し、ＬＵＴ保持バッファ３１１にＬＵＴを保存する。ＬＵＴは入力画像データの画素値を、ノイズ混入前の画素値の推定値に変換するためのテーブルである。
【００３６】
ＬＵＴ作成部３１０の説明のため、例として、画素値が９０の入力画像データに対する、ノイズが混入する前の画素値を導出する方法を以下に説明する。
ノイズの混入前のなんらかの画素値がノイズの混入により、入力画像データの画素値が９０になると推定する。そこで、どのような画素値にノイズが混入することで入力画像データの画素値が９０となったかを考える。まず図４（ａ）のノイズ特性テーブルにおいて、ノイズ混入後に画素値が９０となる行（実際値が９０になる行）を抜き出したものを図８（ａ）に示す。図８（ａ）によれば、画素値が９０になる確率が最も高い元の画素値は９０である。また、元の画素値が８９や９１など画素値が９０に近い画素値は、ノイズが混入すると画素値９０となる確率を有する。しかし、元の画素値が０である画素が、ノイズが混入することにより画素値が９０になる確率は０となっている。このように図４（ａ）のノイズ特性テーブルを横方向に読むことによって、ノイズ混入により画素値９０となり得る画素値の確率が特定され得るのである。
【００３７】
図４（ａ）のノイズ特性テーブルは縦方向には１で正規化されているが、横方向には正規化されていない。そこで、図４（ａ）のテーブルを横方向に正規化することで、入力画像データの画素値Ｘ’として９０を得た場合に、元の画素値の確率を表すテーブルとして利用できる。
【００３８】
図８（ｂ）は、図８（ａ）のノイズ特性テーブルの画素値が９０の行を正規化したものを示している。図８（ｂ）において、画素値が９０の行の正規化画素値頻度を合計すると１となる。
この図８（ｂ）のような横方向の正規化画素頻度の分布をノイズ変動特性分布と呼ぶ。
【００３９】
また同様に、回復ヒストグラムも合計値が１になるように正規化することで正規化回復ヒストグラム８０６を取得する。
【００４０】
図７（ａ）は、正規化回復ヒストグラム８０６及びノイズ変動特性分布８０１を同時に図示している。ノイズ混入前の元の画素値Ｘは、正規化回復ヒストグラムとノイズ変動特性分布の重複する領域に存在する可能性が非常に高い。そこで、正規化回復ヒストグラム８０６とノイズ変動特性分布８０１を画素値ごとに乗算した結果をノイズ混入前画素値の確率分布８０２と見なす。そして、ノイズ混入前画素値の確率分布８０２の重心に相当する画素値をもって、ノイズ混入前の画素値の推定値Ｘとする。なお、ノイズ混入前の画素値の推定値Ｘは、ノイズ混入前画素値の確率分布８０２の中央値や最頻値等の他の統計値であっても良い。
【００４１】
ＬＵＴ作成部３１０では上記のような計算を、画素値が９０以外の他の画素値に対しても行い、入力画像データの画素値Ｘ’から、ノイズ混入前の画素値の推定値Ｘを対応づけるＬＵＴを作成する。
【００４２】
図７（ｂ）は、入力画像データの画素値Ｘ’とノイズ混入前の画素値の推定値Ｘ（理想値Ｘ）との対応関係を示すＬＵＴの一例である。作成されたＬＵＴはＬＵＴ保持バッファ３１１に保持される。なお、代表的な画素値のみに対してノイズ混入前の画素値を推定し、代表的な画素値以外の画素値については、補間処理により推定値を取得することによりＬＵＴを作成しても良い。
【００４３】
次に、画素値変換部３１２では入力画像データにＬＵＴ保持バッファ３１１にあるＬＵＴを適用して処理画像（補正画像）を得る。この処理により、入力画像データの画素にノイズ混入前の画素値の推定値を割り当てた結果、処理画像データはノイズ混入前画像の推定画像データとなる。処理画像は入力画像データが示す画像に比べノイズが低減されている。
【００４４】
最後に、画素値変換部３１２では信号処理部１０３における他処理に処理画像データを出力する。なお、画素値変換部３１２は、ＬＵＴを参照することなく、ノイズ混入前の画素値の推定値を入力画像データの画素毎に算出し、その推定値を出力してもよい。
【００４５】
以上の動作により、ノイズの低減された画像を得る事が出来る。画素値ごとにＬＵＴを適用することでノイズを低減しているため、入力画像の空間周波数に応じて処理内容が変化することはない。従って入力画像の空間周波数に関わらずノイズを低減できる。
【００４６】
（回復行列）
以下、回復行列の生成及び適用方法について説明する。
図４（ａ）に示されたノイズ特性を表すノイズ特性テーブル中の正規化画素頻度を要素として持つ行列をＤと表記する。行列Ｄの要素Ｄｉｊはｊ番目のノイズ混入前の画素値（図４（ａ）の中央画素値）がノイズ混入によってｉ番目の画素値（図４（ａ）の画素値）になる確率を表している。入力画像データの一次元ヒストグラムをベクトル表記でｐ’（ｘ）とする。ｘは画素値である。ノイズ混入前の画像のヒストグラムを同じくベクトル表記でｐ（ｘ）とする。行列Ｄ及びベクトルｐ’（ｘ）とｐ（ｘ）との間には次式の関係が成り立つ。
【００４７】
【数１】
0005832095.tif 000002
【００４８】
式（１）は、ノイズ特性を示す行列Ｄとノイズ混入前の画像のヒストグラムｐ（ｘ）との積により、ノイズ混入後の画像のヒストグラムｐ’（ｘ）（入力画像の一次元ヒストグラム）が表されることを示している。ノイズ混入前の画像のヒストグラムｐ（ｘ）を導出するには、式（２）のように行列Ｄの逆行列Ｄ−１を式（１）の両辺に乗ずれば良い。
【００４９】
【数２】
0005832095.tif 000003
【００５０】
なお、行列Ｄの特性によっては正確にベクトルｐを計算することができない場合がある。例えば、行列Ｄが非正方行列や非正則行列の場合である。この際、Ｄの疑似逆行列Ｑを計算し、Ｑｐ’（ｘ）をもってノイズ混入前の画像のヒストグラムのベクトルｐ（ｘ）の近似されたヒストグラム
【００５１】
【数３】
0005832095.tif 000004
【００５２】
とする。式（３）は、擬似逆行列Ｑによってノイズ混入前の画像のヒストグラムを取得する際の式である。
【００５３】
【数４】
0005832095.tif 000005
【００５４】
この疑似逆行列Ｑを改良し、後述する統計ノイズに対してもロバストな回復行列Ｒを計算してＲｐ’によりノイズ混入前のヒストグラムｐの近似とすることも可能である。
【００５５】
図２３を用いて統計ノイズを説明する。図２３の実線は計算機で確率分布が正規分布に従う乱数を発生させ、そのヒストグラムをプロットしたものである。発生させる乱数が有限であるかぎり、ヒストグラムは正規分布にはならない。発生させる乱数の数を増加させていくと、ヒストグラムは正規分布に近づいていく。このように、サンプルが有限であるが故に発生する、理想ヒストグラムと実際のヒストグラムとの誤差を統計ノイズと称する。この統計ノイズの存在を考慮して疑似逆行列を改良したものが回復行列Ｒである。疑似逆行列Ｑはノイズによるヒストグラムの分散を低減するため、ヒストグラムを鮮鋭にする作用がある。しかし、実際のヒストグラムには統計ノイズが存在するため、単純に疑似逆行列Ｑを入力画像データのヒストグラムｐ’に適用すると、統計ノイズまでも鮮鋭化してしまい、ノイズ混入前のヒストグラムの推定誤差が大きくなる。そこで、統計ノイズが大きい場合には疑似逆行列Ｑによる鮮鋭化の程度を抑えることが望ましい。統計ノイズの大小はヒストグラム作成に用いたサンプルの個数、すなわち総画素数に依存するため、総画素数に応じて疑似逆行列Ｑによる鮮鋭化の程度を調整した回復行列Ｒを用いることが望ましい。回復行列Ｒの計算方法の一例を次式に示す。
【００５６】
【数５】
0005832095.tif 000006
【００５７】
上式の回復行列Ｒは単位行列Ｉと疑似逆行列Ｑの加重平均であることを意味している。係数ａ，ｂはヒストグラムｐ’（ｘ）、及び入力画像データの総画素数Ｎに応じて変化する係数である。総画素数Ｎが多い場合はａが大きくｂが小さくなることによってＲを疑似逆行列に近づけ、回復処理を行う。一方で、総画素数Ｎが少ない場合は、ａを小さくｂを大きくすることにより統計ノイズによる推定誤差を抑える。
【００５８】
係数ａ，ｂが総画素数Ｎ以外にもヒストグラムｐ’（ｘ）に依存する理由を説明する。例えば、ヒストグラムにおいて画素値が多く分布する範囲と、画素値が少なく分布する範囲がある場合、前者の範囲は後者の範囲よりもサンプル数が多いため、より統計ノイズの少ないヒストグラムの範囲である。そのため、画素値の多い範囲ではａの値を上げ、ｂの値を下げて、回復行列Ｒを疑似逆行列Ｑに近づけても統計ノイズによる誤差は少ないと考えられる。一方で、サンプル数が少ない範囲では統計ノイズによる推定誤差が大きくなると考えられるため、係数ｂの値を大きく、係数ａの値を小さく設定することにより統計ノイズによる影響を相対的に小さくすることが望ましい。
【００５９】
式（４）により導出される回復行列Ｒを用いた回復処理を以下の（５）式に示す。
【００６０】
【数６】
0005832095.tif 000007
【００６１】
（一次元ヒストグラムを用いた場合のノイズ低減方法の制限）
上記の一次元ヒストグラムを用いたノイズ低減方法では、実際のノイズ低減には不充分である場合があることを示す。図９（ａ），（ｂ）はグラデーション画像と、グラデーション画像にノイズが重畳された画像をそれぞれ示す。図９（ａ），（ｂ）の一次元ヒストグラムは図１０（ａ），（ｂ）のように平坦なものとなる。これらの一次元ヒストグラムの回復ヒストグラムはやはり平坦であり、回復前後でヒストグラムに差が生じにくい。
【００６２】
図１１は、入力画像データが図１０（ａ），（ｂ）のようなグラデーション画像の場合の回復処理を示している。図１１から明らかなように、（正規化）回復ヒストグラムとノイズ変動特性分布の積はノイズ変動特性分布と同一形状となることが分かる。すると、ノイズ混入前の画素値の確率分布の重心により示される画素値は入力画像データの画素値と同一となる。そのため、実質的にノイズ低減処理を行う事ができない。この問題はノイズ混入という現象が一次元ヒストグラムの変化に現れていないことに起因する。
【００６３】
（多次元ヒストグラムを用いた場合のノイズ低減方法）
そこで本実施例では、多次元ヒストグラムを用いたノイズ低減処理を示す。
図３（ｂ）は、多次元ヒストグラムを用いたノイズ低減処理に必要なブロック構成図を示している。なお、図３（ｂ）の処理の一部は、図３（ａ）に記載の一次元ヒストグラムを用いたノイズ低減処理と同様の処理を行なうため重複する部分の詳細は省略する。
【００６４】
まず、入力部３５１は、入力画像データ、総画素数、及びノイズ特性データを取得する。
そして、多次元ヒストグラム生成部３５４は、入力画像データ、及び総画素数により多次元ヒストグラムを生成する。多次元ヒストグラム生成部の詳細を図２０（ｂ）に示す。周辺画素選択部２０１１は入力画像データの画素値をＡ／Ｄ変換部から逐次取得し、注目画素の画素値と、所定の位置関係にある画素の画素値の組合せを頻度計算部２０１２に出力する。頻度計算部２０１２では組合せの頻度をカウントし、入力部３５１により入力される画素値で除算することにより多次元ヒストグラムを生成する。本実施例において、多次元ヒストグラムとは、複数画素を一組とする頻度分布を言う。即ち、多次元ヒストグラムとは、入力画像データが示す入力画像中の各注目画素の画素値と所定の位置関係にある画素の画素値との組合せの出現頻度に基づいて生成される。本実施例では、一次元ヒストグラムを用いたノイズ低減処理における画素値、および、或いは多次元ヒストグラムを用いたノイズ低減処理における画素値と画素値の組合せを特性値と呼ぶ。
【００６５】
本実施例では、多次元ヒストグラムは二次元ヒストグラムであり、位置関係が「注目画素と右隣の画素」とする。入力画像データの二次元ヒストグラムは、入力画像データを走査することにより、注目画素の画素値がＸであり、かかる注目画素の右隣の画素の画素値がＹである組合せを画像全体に渡ってカウントする。そして、一次元ヒストグラムの場合と同様に、かかる頻度を画像全体の総画素数で除算することで二次元ヒストグラムを取得する。図１２（ａ）は、図９（ａ）のグラデーション画像に対する二次元ヒストグラム、図１２（ｂ）は、図９（ｂ）のノイズを含むグラデーション画像に対する二次元ヒストグラムを示している。Ｘは注目画素の画素値、Ｙは注目画素と所定の位置関係ある画素の画素値を示している。
【００６６】
図１２（ａ）（ｂ）では、実線で囲った範囲内に０でない頻度（出現確率）が分布しているものとする。
図１２（ａ）は、対角方向に０でない頻度が分布している。これはグラデーション画像では隣接する画素の画素値が近いため、Ｘ＝Ｙとなる直線の付近に分布が発生するためである。一方で、Ｘ＝Ｙとなる直線から外れた領域での分布は発生しない。
【００６７】
図１２（ｂ）は、同様に対角方向に０でない頻度が分布している。加えて、ノイズの影響によりＸ＝Ｙとなる直線から離れた領域においても分布が存在する。このような広がった分布は、図９（ｂ）からも明らかなようにノイズを含む画素とノイズを含まない画素との組合せにより生じる分布である。
【００６８】
以上のとおり、二次元ヒストグラムではノイズ混入によるヒストグラムの形状の変化が顕著に見られる。これは、図１０（ａ）（ｂ）に示した、ノイズ混入による一次元ヒストグラムの形状の変化とは大きく異なる点である。
【００６９】
以下、回復テンソル生成部３５６に入力されるノイズ特性テーブルについて説明する。図４（ｂ）は、中央画素値が（３０，４５）の組合せのノイズ特性テーブルを示している。図４（ｂ）は、中央画素値（３０，４５）が入力された場合、画素値が（３０，４５）の組合せとなる確率は０．０４０であることを示している。二次元ヒストグラムを用いてノイズ低減処理を行なう場合は、図４（ｂ）のようなテーブルを全ての画素値の組合せに対して用意することとなる。前述の通り、全ての画素値の組合せについて撮影して取得することが困難な場合は、補間演算などをつかってノイズ特性テーブルを生成しても良い。
回復テンソル生成部３５６にて行われる処理を説明する。まず、二次元ヒストグラムを用いる場合に上述の（１）式に対応する、ヒストグラムとノイズの関係を表すモデル式は次式となる。
【００７０】
【数７】
0005832095.tif 000008
【００７１】
（６）式において、ｐ’（ｘ，ｙ）は多次元ヒストグラム生成部３５４で生成される入力画像データの二次元ヒストグラムである。ｘ、ｉは注目画素の画素値、ｙ、ｊは注目画素の右に隣接する画素（第一の位置関係にある画素）の画素値、を示している。ｐ（ｉ，ｊ）はノイズ混入前画像の二次元ヒストグラムである。また、Ｄｘｉはノイズ混入によってｉ番目の画素値がｘ番目の画素値となる確率である。Ｄｙｊも同様である。（６）式はノイズ混入によって多次元ヒストグラムがＸ方向，Ｙ方向に広がることを意味している。一次元ヒストグラムを用いる場合と同様に行列Ｄの擬似逆行列Ｑを用いて（６）式を以下のように変形することができる。
【００７２】
【数８】
0005832095.tif 000009
【００７３】
Ｑｘｉ、Ｑｙｊそれぞれは、Ｄｘｉ、Ｄｙｊの擬似逆行列である。
【００７４】
【数９】
0005832095.tif 000010
【００７５】
はノイズ混入前の二次元ヒストグラムｐ（ｘ，ｙ）の推定値である。この場合においても統計ノイズを考慮して、（４）式に対応する（８）式を用いて回復テンソルＲｘｉＲｙｊ（変換データ）が求められる。
【００７６】
【数１０】
0005832095.tif 000011
【００７７】
ここで、δｘｉ、δｙｊは単位行列を示す。一次元ヒストグラムを用いたノイズ低減処理の変換データＲは行列形式であるが、二次元以上のヒストグラムを用いたノイズ低減処理の変換データは、テンソル形式となる。一次元ヒストグラムの場合と同様に、係数ａ，ｂはヒストグラムｐ’および入力画像データの総画素数Ｎに応じて変化する係数である。二次元ヒストグラムにおいて画素値が多く分布する範囲と画素値が少なく分布する範囲が存在する場合、以下のような処理をする。即ち、画素値の多い範囲では係数ａの値を上げ、ｂの値を下げて、回復テンソルＲを擬似逆行列Ｑに近づける。一方で、総画素数（サンプル数）が少ない範囲では、統計ノイズによる推定誤差が大きくなると考えられるため、係数ｂの値を大きく、係数ａの値を小さく設定する。以上により、統計ノイズの影響を相対的に小さくすることができる。
【００７８】
以上のように求められた回復テンソルは、回復テンソル保持バッファ３５７に保持される。
そして、ヒストグラム回復部３５８により、回復テンソル保持バッファ３５７に保持されている回復テンソルＲｘｉＲｙｊを用いて、回復ヒストグラム（回復多次元ヒストグラム）
【００７９】
【数１１】
0005832095.tif 000012
【００８０】
を以下の（９）式により求める。
【００８１】
【数１２】
0005832095.tif 000013
【００８２】
多次元ＬＵＴ作成部３６０で作成される多次元ＬＵＴは２つの画素値（Ｘ’、Ｙ’）対２つの画素値（Ｘ、Ｙ）の変換テーブルである。一次元ヒストグラムによるノイズ低減処理の場合と同様に、ノイズ変動特性分布と回復ヒストグラムとの積を画素値の組合せごとに計算し、ノイズ混入前の画素値の確率分布を求める。そして、ノイズ混入前イの画素値の確率分布の重心をノイズ混入前の画素値の推定値（Ｘ，Ｙ）とする。
【００８３】
図１３は二次元ヒストグラムを用いたノイズ低減処理を示す図であり、横軸Ｘが注目画素の画素値を、縦軸Ｙが注目画素の右隣の画素の画素値を表す。また、図示していないが、注目画素の画素値Ｘと注目画素の右隣の画素の画素値Ｙとの組合せの頻度（又は出現確率）を示す軸が紙面と垂直な方向に存在する。なお、実線（１４０３、１４０４）ないし破線（１４０２）で囲った範囲内に０でない頻度が分布しているとする。ヒストグラム回復部３５８によるヒストグラム回復によって、入力画像データの二次元ヒストグラム１４０２に対して回復ヒストグラム１４０３が計算される。回復ヒストグラム１４０３は入力画像の二次元ヒストグラム１４０２よりも頻度が鮮鋭化している。入力画素値１４０１を得た場合、まずノイズ変動特性分布１４０４が行列Ｄにより計算される。具体的には入力画素値１４０１が（Ｘ’，Ｙ’）であった場合、ノイズが加わる前の状態が（Ｘ，Ｙ）である確率はＤＸＸ’ＤＹＹ’で表される。この確率分布がノイズ変動特性分布１４０４である。次にノイズ変動特性分布１４０４と回復ヒストグラム１４０３の積を計算し、ノイズ混入前画素値の確率分布１４０５を得る。そして、ノイズ混入前の画素値の確率分布１４０５の重心１４０６を計算し出力画素値とする。以上の計算を入力画像データ中の各（Ｘ’，Ｙ’）に対して行い、入力画素値に対応する出力画素値（Ｘ，Ｙ）を計算して、二次元ＬＵＴを作成する。作成された二次元ＬＵＴは、多次元ＬＵＴ保持バッファ３６１に保持される。
【００８４】
図７（ｃ）は、二次元ＬＵＴの一例を示している。二次元ＬＵＴ中、左列は入力画像データの画素値の組合せ（Ｘ’，Ｙ’）を表し、右列に推定された画素値の組合せ（出力画素値（Ｘ，Ｙ））を表している。
【００８５】
画素値変換部３６２は、入力画像データに多次元ＬＵＴ保持バッファ３６１にある二次元ＬＵＴを適用して処理画像（補正画像）を得る。具体的には、入力画像データを走査し、順次、注目画素と右隣にある画素の画素値の組合せ（Ｘ’，Ｙ’）を取得する。そして、二次元ＬＵＴを参照することにより対応する出力画素値（Ｘ，Ｙ）を取得し、取得された出力画素値Ｘを注目画素の画素値とする。
【００８６】
ここで、入力画像データが高周波成分を含んでいたとしても、上記の処理によってノイズ低減が良好に行えることを説明する。高周波成分を有する画像として図１４のようなグラデーション画像の中に画素値Ｅ，画素値Ｆの画素から構成されるエッジを有する画像にノイズがのった入力画像データを考える。平滑化フィルタでは図１４に示すエッジ部の周辺で、エッジをぼけを抑えつつノイズを低減することが困難であった。
【００８７】
図１５に入力画像データの図１４のグラデーション画像の二次元ヒストグラムを示す。横軸は注目画素の画素値Ｘ’、縦軸はその右隣の画素値Ｙ’である。図１５において破線で囲った領域は、入力画像データの二次元ヒストグラムにおいて頻度（又は頻度確率）が０でない領域を表す。図１５に示すようにグラデーションとエッジが混ざった画像にノイズがのった画像のヒストグラムは２つの領域に集中して分布する。すなわち、グラデーション領域、平坦領域に属する画素値の組合せは領域１６０１に対応し、エッジ部の画素値の組合せは領域１６０２に対応する。エッジ部近傍にノイズが重畳された画素値の組合せ１６０５を得たとする。回復ヒストグラム処理により、二次元ヒストグラム１６０１に対応する回復ヒストグラム１６０３が生成される。同様に、回復ヒストグラムにより、二次元ヒストグラム１６０２に対応する回復ヒストグラム１６０４が生成される。回復ヒストグラム１６０４、１６０５の実線で囲まれた領域に、０でない頻度（又は頻度確率）が分布する。画素値の組合せ１６０５に対応するノイズ変動特性分布１６０６と回復ヒストグラムとの積を計算しノイズ混入前の画素値の確率分布を取得する。そして、ノイズ混入前の画素値の確率分布の重心を計算することで出力画素値１６０８を得る。さて、入力画素値１６０５はノイズがのったエッジ部近傍の画素の組合せである。この組合せが本実施例のノイズ低減により出力画素値１６０８へと変換される。多次元ヒストグラムにおいてコントラストが高い画素値の組合せは図１５では左上ないし右下に配置される。すると、本実施例のノイズ低減処理によってエッジ部近傍の入力画素値の組合せ１６０５は左上に移動したのであるから、コントラストは下がるどころか、むしろ向上している。また本実施例のノイズ低減により、入力画素値１６０５がノイズが重畳される前の状態の画素値１６０７に近づいているため、ノイズも低減されている事が分かる。以上により、本実施例では入力画像データが高周波成分を含んでいても、良好にノイズ低減が行えることが分かる。
【００８８】
その一方で、低〜中周波成分が主要なグラデーション部分に属する画素値の組合せ１６０９は、回復ヒストグラムと画素値の組合せ１６０９のノイズ変動特性分布１６１０の積の重心である１６１１に変換される。グラデーション部分の注目画素と右隣の画素値の差をＣとして、図１５の横軸をＸ、縦軸をＹとおくと、多次元ヒストグラムの一部領域である１６０１は線分Ｙ＝Ｘ＋Ｃで表される中心軸１６１２付近に０でない値が分布する。ヒストグラム回復により、領域１６０１は中心軸１６１２により近い領域１６０３に変換される。そのため、移動後の画素値の組合せである１６１１は１６０９より中心軸１６１２に近い座標であり、ノイズ混入前の状態であるＹ＝Ｘ＋Ｃ線分上に近づいている。つまり１６１１は入力画素値の組合せ１６０９よりノイズの低減された状態であるから、グラデーション部に対しても良好にノイズが低減されたことになる。
【００８９】
以上により、入力画像の高周波成分であるエッジ部においても、低〜中周波成分であるグラデーション部においても本実施例によりノイズが低減されることが示された。言い換えれば、本件発明では入力画像の空間周波数に関わらず、良好にノイズを低減できることが分かる。
【００９０】
以上、多次元ヒストグラムの例として二次元ヒストグラムを用いてノイズ低減処理を説明した。本実施例においては、画像の右端に位置する注目画素を補正する際、注目画素の右隣の画素値を得る事ができない。その場合には、端に位置する注目画素の右隣に注目画素と同じ画素値を有する画素があるとしてヒストグラムを作成したり、端に位置する注目画素に対して補正処理をしないことが考えられる。入力画像がベイヤ配列を採用したセンサで撮影された画像である場合は、注目画素の右隣が注目画素と同じ色に対応する画素であるとは限らない。そこで入力画像中の画素が、複数の色のいづれかに対応している場合、色毎に本実施例を適用しても良い。また、入力画像を間引いて縮小した画像のノイズ低減に本実施例を適用する場合には、間引かれない画素を対象として本実施例を適用する事が望ましい。また、入力画像の奇数列と偶数列を別々の画像処理装置で処理する場合には、それぞれの演算装置で奇数列のみ、あるいは偶数列のみからなる画像に本実施例を適用しても良い。
【００９１】
あるいは画像を領域分割し、領域ごとに本実施例のノイズ低減処理を適用してもよい。かかる領域とは色、入力画素の位置、入力画素値などで区分された画素の集合のことである。色により画素を区分する場合、例えば、Ｒ値／Ｂ値／Ｇ値で区分したり、あるいは輝度値／色差値で区分してもよい。入力画素の位置により画素を区分する場合、所定のサイズのブロック（メッシュ）で区分したり、画像判定技術（例えば、顔認識）により区分してもよい。領域毎のノイズ低減処理の結果は合成され、最終的な補正画像が取得される。また、ノイズ低減処理を行なった領域と、ノイズ低減処理を行なわなかった領域を合成して最終的な補正画像を取得してもよい。
【００９２】
二次元ヒストグラムを生成する際に注目画素の画素値と右隣の画素の画素値を用いた。しかし、この位置関係は右隣の画素である必要はなく、例えば、左隣の画素、或いは注目画素と一画素を挟んで右の画素などとしてもよい。また、三次元以上のヒストグラムに基づいてノイズ低減処理を行なっても良い。例えば、入力画像データの注目画素の画素値と右隣にある画素の画素値と下隣にある画素の画素値との組合せにより多次元ヒストグラムを生成しても良い。
入力画像と、ノイズ低減処理後の出力画像との関係からみた本実施例の特徴は、入力画像のある画素値とその隣の画素値の組合せによって出力画素の値が決定する事である。さらにかかる決定の方法は、入力画素値の組合せを、入力画像の二次元ヒストグラムにおいて、頻度の高い画素値の組合せに近づけるようになされることである。これはヒストグラム回復により、入力ヒストグラムの山はより急峻になり、山の付近の画素値の組合せはより山に近づくように変換されるからである。
【００９３】
［実施例２］
実施例１のノイズ低減処理では、入力画像を走査しつつ、図７（ｃ）のＬＵＴを用いて入力画像データの画素値の組合せ（Ｘ’、Ｙ’）に対して、ノイズ混入前の画素値の推定値の組合せ（Ｘ，Ｙ）を取得する。そして、画素値変換部３６２は、ＬＵＴを参照することにより注目画像の画素値Ｘ’を推定値Ｘに変換する。
【００９４】
本実施例は、画素値変換部３６２で一つの画素に対して画素を割り当てる際に少なくとも二つの組合せの画素値を使用する。
図２２（ａ）は、画素値変換部３６２により画素値Ｗを割り当てる際の組合せを模式的に示している。ＡからＩは入力画像データの画素を示している。
注目画素が画素Ａの場合、画素Ａの画素値ＸＡ’と画素Ｂの画素値ＹＢ’と基づいて実施例１に記載されたノイズ低減処理によりノイズ混入前の画素値の推定値ＸＡとＹＢとを得る。同様に注目画素が画素Ｂの場合、画素Ｂの画素値ＸＢ’と画素Ｃの画素値ＹＣ’と基づいてノイズ混入前の画素値の推定値ＸＢとＹＣとを得る。
【００９５】
本実施例の画素値変換部３６２は、画素Ｂに対するノイズ低減処理後の画素値の推定値ＷＢを以下の式（１０）により割り当てる。
ＷＢ＝（ＹＢ＋ＸＢ）／２    （１０）
即ち、実施例１では、画素Ｂに対して推定値ＸＢを割り当てられる。一方で、本実施例では、画素Ｂに対して二つの組合せ（１７０１および１７０２）それぞれについて推定された画素Ｂの推定値（ＹＢおよびＸＢ）の平均値ＷＢが割り当てられる。
【００９６】
画素値変換部３６２は、他の画素に対しても式（１０）と同様に、ノイズ低減処理後の画素値の推定値Ｗを割り当てる。
【００９７】
図２２（ｂ）は、本実施例における三次元ヒストグラムを用いた場合のノイズ低減処理を実行した場合における割り当て処理の概要を示した図である。画素Ｅに対して、ノイズ低減処理後の画素値の推定値ＷＥを割り当てる際に、以下の組合せによる画素値の推定値を使用する。
組合せＢ，Ｃ，Ｅにおける画素Ｅの推定画素値Ｖ１
組合せＤ，Ｅ，Ｇにおける画素Ｅの推定画素値Ｖ２
組合せＥ，Ｆ，Ｈにおける画素Ｅの推定画素値Ｖ３
Ｖ１、Ｖ２、Ｖ３を使用してＷＥは以下の式により求める。
ＷＥ＝（Ｖ１＋Ｖ２＋Ｖ３）／３      （１１）
以上のように、本実施例では、より多くの情報を使うためノイズ混入前画像の推定精度が高まることが期待される。
【００９８】
また、入力画像データを大きく変更したくない場合は、ノイズ低減の効果を抑えるために、Ｖ１、Ｖ２、Ｖ３の中で入力画像の画素Ｅの画素値に最も近い、すなわち補正による変化の少ない推定画素値を割り当てても良い。このように割り当て方を様々変える事によってノイズ低減後の画像に所望の特性を持たせることが出来る。
【００９９】
また本実施例では、画像のＤＣが必ずしも保存されない。入力画像の画素値の組合せを（Ｘ’，Ｙ’），ノイズ混入前の画素値の推定値の組合せ（Ｘ，Ｙ）とする。すると（Ｘ’＋Ｙ’）／２と（Ｘ＋Ｙ）／２は必ずしも等しくなるとは限らない。そこで、平均値が保存されるようにＬＵＴを補正してもよい。ＤＣが保存されるよう補正されたノイズ混入前の推定値を（Ｘ２，Ｙ２）と定義する。すると（Ｘ２，Ｙ２）は
Ｘ２＝Ｘ−（Ｘ−Ｘ’＋Ｙ−Ｙ’）／２            （１２）
Ｙ２＝Ｙ−（Ｘ−Ｘ’＋Ｙ−Ｙ’）／２            （１３）
と計算できる。上式によれば
（Ｘ２＋Ｙ２）／２＝（Ｘ＋Ｙ）／２              （１４）
となりＤＣが保存されていることが分かる。
【０１００】
［実施例３］
（多次元ヒストグラムを用いた場合のノイズ低減方法の拡張）
以下、画像を複数の領域に分割して、よりエッジをぼけにくくしながら、ノイズ低減処理を実現する実施例を示す。
【０１０１】
図１８は、本実施例の自然画像の２次元ヒストグラムの一例を示している。図１８のＸ軸、Ｙ軸は、図１３や図１５のＸ軸、Ｙ軸と同様である。自然画像では隣接画素が近い画素値をとる傾向があるため、２次元ヒストグラムでは直線Ｙ＝Ｘの付近の平坦領域頻度分布４３０２の頻度が非常に高い。平坦領域頻度分布４３０２は、画像中の平坦領域に対応する。一方、画像中のエッジ領域に対応する２次元ヒストグラム上のエッジ領域頻度分布４３０１の頻度は平坦領域頻度分布４３０２の頻度に比べて非常に小さい。図１９は線分４３０４に沿った頻度を縦軸にプロットした図である。平坦領域頻度分布４４０２の頻度はエッジ領域頻度分布４４０１の頻度に比べて非常に大きい。今、図１８の４３０３、図１９の４４０３で表される入力画素値の組合せを得たとする。入力画素値の組合せは、例えば、注目画素の画素値とその右隣の画素の画素値との組み合わせである。入力画素値の組合せ４４０３とエッジ領域頻度分布の頻度の中央値との距離４４０６は、入力画素値の組合せ４４０３と平坦領域頻度分布４４０２の頻度の中央値との距離４４０５に比べて短い。実施例１と同様に、図１９において、ノイズ変動特性分布４４０４と、２次元ヒストグラム（４４０１と４４０２）との積をとり補正量を決定する。ここで、平坦領域４４０２の頻度が非常に大きいために距離的にはエッジ領域頻度分布４４０１に近くても、平坦領域頻度分布４４０２に近づくように入力画素値の組合せ４４０３が補正されてしまう。なお、ノイズ変動特性分布４４０４は、図７（ａ）のノイズ変動特性分布８０１に対応する。この場合、エッジ領域頻度分布４４０１とノイズ変動特性分布４４０４との間に重複部分があるため、急激に平坦領域頻度分布４４０２に近づくような補正をされるようなことはない。エッジ領域頻度分布４４０１の頻度は入力画素値の組合せ４４０３がぼけすぎないよう、ブレーキ（抑制）の役割は果たしている。仮に入力画素値の組合せ４４０３がノイズ混入前には平坦領域であった場合は、エッジ領域の頻度は無視して、直線Ｙ＝Ｘに近づけたほうが良い。また、仮に入力画素値の組合せ４４０３がノイズ混入前にはエッジ領域であった場合には平坦領域の頻度は無視してエッジ領域頻度分布４４０１に近づけた方が良い。しかし、事前に入力画素値の組合せ４４０３がどちらの領域に属するかという画素の特徴が分からないため、上記のような中途半端な補正をすることとなる。以下、これらの課題を解決するための本実施例について説明する。
【０１０２】
まず、画像中の全ての画素について予め領域判定を行い、各画素の特徴を決定する。本実施例においては、注目画素の特徴とは、画像中の注目画素とその周辺画素とから、注目画素が平坦領域に属する確率、及び注目画素がエッジ部に属する確率を示す情報である。そして、それら確率に基づいて画素値の補正を行う。図１９の例を用いて説明すると、仮に入力画素値の組合せ４４０３が１００％の確率で平坦領域に属すると判定されていれば、エッジ領域頻度分布４４０１の頻度を無視してノイズ低減処理を行う。逆に入力画素値の組合せ４４０３が１００％の確率でエッジ領域に属すると判定されたならば、平坦領域頻度分布４４０２の頻度を無視して、ぼかすことなくノイズ低減処理がなされる。どちらも１００％でない場合には、平坦領域頻度分布４４０２、エッジ領域頻度分布４４０１の頻度を適正に考慮して補正がなされる。
【０１０３】
図１６は本実施例におけるノイズ低減処理にかかる構成を示すブロック図である。実施例１に記載されている図３（ｂ）との差異の一つは領域判定部４０２０を設けたことである。領域判定部４０２０ではＡ／Ｄ変換部１０２から得た入力画像データの各画素について領域判定データを生成する。
【０１０４】
図１７は、領域判定部４０２０のブロック図を示している。周辺画素信号選択部４２２０は、注目画素の画素値と注目画素の周辺画素の画素値とを入力画像データから取り出す。例えば注目画素を囲む１１ｘ１１ピクセルほどの画像領域（以降、注目領域）の画素値を出力する。
平坦度計算部４２２１は、まず注目領域の平均画素値μ、及び分散σを計算する。仮に注目領域が平坦領域であれば、注目領域内の画素のばらつきはノイズによるものとみなせる。ノイズは正規分布に従うため、注目領域内の画素群は正規分布に従って分布しているはずである。統計学では、平均μ、分散σの正規分布に従うｋ個の独立な変数をｘｉとすると、
【０１０５】
【数１３】
0005832095.tif 000014
【０１０６】
は自由度ｋのカイ二乗分布という確率分布に従うことが知られている。カイ二乗分布を用いれば、注目領域内の画素群から計算したＺがどれくらいの確率で得られるものかが計算できる。例えばＺが得られる確率が０．００１％であるとすると、非常に希な事象が発生したと考えられる。本実施例ではこの０．００１％という確率を注目領域が平坦領域である確率と読み替える。以下、具体的には注目領域が平坦領域である確率Ｐ-ｆｌａｔの計算方法を述べる。カイ二乗分布の累積分布関数をＦ（ｘ，ｋ）とすると
【０１０７】
【数１４】
0005832095.tif 000015
【０１０８】
と表される。ここでγは不完全ガンマ関数、Γはガンマ関数である。Ｐｆｌａｔを以下の式で計算する。
【０１０９】
【数１５】
0005832095.tif 000016
【０１１０】
このように計算されたＰｆｌａｔは平均μ、分散σの正規分布に従うｋ個の独立な変数から計算したＺがＺ以上の値を取る確率に対応している。Ｚが大きいほど希な事象であるため、Ｚが大きくなるに従ってＰｆｌａｔは小さい値をとるようになる。
【０１１１】
以上のとおり平坦度計算部４２２１は、式（１５）乃至（１７）に基づいて、注目領域が平坦領域である確率Ｐ-ｆｌａｔを出力する。注目領域が平坦である確率はＰ-ｆｌａｔであるため、注目領域がエッジ領域である確率は１−Ｐ-ｆｌａｔとなる。本実施例ではさらにエッジ領域を４つの領域に分ける。すなわち、横方向、縦方向、右斜め、左斜めに相関がある領域である。各領域はそれぞれ横のエッジ、縦のエッジ、右斜めのエッジ、左斜めのエッジに相当する領域である。
【０１１２】
図２１（ａ）の位置関係にあるＸＹの頻度に基づいて生成される２次元ヒストグラムを考える。図３（ｂ）に示される多次元ヒストグラムを用いたノイズ低減処理において、複数の画素を一組として頻度を計算する理由は、ノイズ混入前の被写体画像における、かかる複数の画素に相関があるからである。つまり、複数の画素に相関があるために図１３（ａ），（ｂ）のように多次元ヒストグラムにおいて０でない頻度が一部に集中し、その集中の度合いを高めることでヒストグラム回復が可能なのである。ここで集中の度合いを高めるとは、多次元ヒストグラムのピーク周辺の画素をピークに近づける事を意味する。例えば、多次元ヒストグラムを生成する場合に、図２１（ａ）に示されるＸＹの位置関係を採用することが考えられる。しかし仮に注目画素が縦方向に相関の高い領域であることが予め分かっていたとすると、図２１（ｂ）の位置関係にある複数の画素を用いて多次元ヒストグラムを構成した方が、図２１（ａ）の位置関係の場合と比べて密集度は高まるはずである。そこで、領域判定部４０２０では平坦領域に加え、横方向相関、縦方向相関、右斜め相関、左斜め相関の計５つの領域を考え、注目画素がそれぞれの領域に属する度合いを確率として出力する。横相関計算部４２２２，縦相関計算部４２２３，左斜め相関計算部４２２５，右斜め計算部４２２６の出力値をそれぞれＳｈｏｒｉｚｏｎｔａｌ，Ｓｖｅｒｔｉｃａｌ，Ｓｌｅｆｔ＿ｄｉａｇｏｎａｌ，Ｓｒｉｇｈｔ＿ｄｉａｇｏｎａｌとする。これら４つの出力値は、注目画素近傍領域の、各方向への相関の強さを示している。本実施例では、これら４つの相関の強さを示す出力値を共相関値と呼ぶ。共相関値を計算するには、例えば次式により計算を行う方法が挙げられる。
【０１１３】
【数１６】
0005832095.tif 000017
【０１１４】
ただし、Ｉ（ｘ，ｙ）は注目画素の近傍領域の画像データであり、ｘ，ｙの総和する範囲は注目画素の近傍領域である。ｘは画素の縦位置，ｙは画像の横位置を示す。上の式において（ｎ，ｍ）＝（０，１），（１，０），（１，１），（１，−１）の場合のＳがＳｈｏｒｉｚｏｎｔａｌ，Ｓｖｅｒｔｉｃａｌ，Ｓｌｅｆｔ＿ｄｉａｇｏｎａｌ，Ｓｒｉｇｈｔ＿ｄｉａｇｏｎａｌに対応する。
【０１１５】
領域信号生成部４２２４は平坦度計算部４２２１の出力Ｐｆｌａｔと、横相関計算部４２２２〜４２２６の出力値Ｓｈｏｒｉｚｏｎｔａｌ，Ｓｖｅｒｔｉｃａｌ，Ｓｌｅｆｔ＿ｄｉａｇｏｎａｌ，Ｓｒｉｇｈｔ＿ｄｉａｇｏｎａｌから、領域判定データを計算する。領域判定データは注目画素の近傍が、平坦、横相関、縦相関、左斜め相関、右斜め相関の５つの領域に属する確率である。これら確率は、Ｐｆｌａｔ，Ｐｈｏｒｉｚｏｎｔａｌ，Ｐｖｅｒｔｉｃａｌ，Ｐｌｅｆｔ＿ｄｉａｇｏｎａｌ，Ｐｒｉｇｈｔ＿ｄｉａｇｏｎａｌにより表される。Ｐｈｏｒｉｚｏｎｔａｌ，Ｐｖｅｒｔｉｃａｌ，Ｐｌｅｆｔ＿ｄｉａｇｏｎａｌ，Ｐｒｉｇｈｔ＿ｄｉａｇｏｎａｌは次式により計算する。
【０１１６】
【数１７】
0005832095.tif 000018
【０１１７】
以上のとおり、領域判定部４０２０は入力画像データの各画素について、平坦、横相関、縦相関、左斜め相関、右斜め相関の５つの領域夫々に属する確率（領域判定データ）を出力する。なお、領域判定データを算出する手法は、上述の手法に限らない。例えば、各方向に応じたエッジ検出フィルタを適用することにより領域判定データを算出しても良い。
【０１１８】
次に多次元ヒストグラム生成部４００４では入力画像データと、領域判定データから、平坦、横相関、縦相関、左斜め相関、右斜め相関の５つの領域それぞれに対応する５つの多次元ヒストグラムを計算して、出力する。
【０１１９】
多次元ヒストグラム生成部４００４の動作を説明する。例として横相関領域の多次元ヒストグラムの計算方法を述べる。横相関領域の多次元ヒストグラムをＨｈｏｒｉｚｏｎｔａｌ（Ｘ，Ｙ）と表記する。Ｘは注目画素の画素値、Ｙは注目画素の右隣の画素の画素値を示す。まず全てのＸ，ＹについてＨｈｏｒｉｚｏｎｔａｌ（Ｘ，Ｙ）を０で初期化する。次に画像を走査して、Ｘ，Ｙの組合せの頻度を示すＨｈｏｒｉｚｏｎｔａｌ（Ｘ，Ｙ）の値を更新する。更新の方法は以下のようである。走査中の注目画素ｉの画素値をＸｉ、注目画素の横相関領域判定データＰｈｏｒｉｚｏｎｔａｌ，ｉとする。注目画素ｉに対し図２１（ａ）の位置関係にある右隣りの画素の画素値ＹをＹｉとおく。ヒストグラムの更新はＨｈｏｒｉｚｏｎｔａｌ（Ｘｉ，Ｙｉ）の値をＰｈｏｒｉｚｏｎｔａｌ，ｉだけ増やす事で行われる。以上の処理を画像中の全画素に対して行うことによって、横相関領域に属する確率で重み付けされた横相関領域の多次元ヒストグラムが計算される。他の４つの領域についても同様である。但し多次元ヒストグラム生成時のＸＹの位置関係は縦相関領域では図２１（ｂ），左斜め相関領域では図２１（ｃ），右斜め相関領域では図２１（ｄ）を用いる。なお、平坦領域における多次元ヒストグラムの生成に当たっては、ＸＹの位置関係は図２１（ａ）乃至（ｄ）のいずれでも良い。その理由は、そもそも平坦領域では領域内の画素が似たような値をとるため、全方向に相関が高い。そのため、どの位置関係を用いても２次元ヒストグラムは頻度が密集する傾向があるからである。本実施例では図２１（ａ）の位置関係を用いることとする。
【０１２０】
回復テンソル生成部４００６の動作は式（７）、式（９）で示される回復テンソルの生成と同様である。本実施例では多次元ヒストグラムが５つあるので、回復テンソルも各多次元ヒストグラムに応じて５つ作成し、回復テンソル保持バッファ４００７に保持する。個々の回復テンソルの作成方法は、二つの画素値を取得する際の画素の位置関係が図２１（ａ）乃至（ｄ）のそれぞれであるということ以外は、実施例１の通りである。平坦領域における回復テンソルの作成については、本実施例では図２１（ａ）の位置関係を用いることとする。
【０１２１】
ヒストグラム回復部４００８では５つの領域の多次元ヒストグラムを対応する５つの回復テンソルでそれぞれ、回復処理を行う。個々の回復処理の内容は実施例１の通りである。
【０１２２】
ルックアップテーブル作成部４０１０では５つの回復多次元ヒストグラム、ノイズ特性データから５つの領域に対応するルックアップテーブルを作成してルックアップテーブル保持バッファ４０１１に保持する。各ルックアップテーブル作成方法は実施例１の通りである。以降、平坦領域のルックアップテーブルにおける入力値（Ｘ，Ｙ）に対する出力値をＬＵＴＸｆｌａｔ（Ｘ，Ｙ）、横相関領域の出力値をＬＵＴＸｈｏｒｉｚｏｎｔａｌ（Ｘ，Ｙ）、縦相関領域の出力値をＬＵＴＸｖｅｒｔｉｃａｌとする。更に、左斜め相関領域のルックアップテーブルにおける入力値（Ｘ，Ｙ）に対する出力値をＬＵＴＸｌｅｆｔ＿ｄｉａｇｏｎａｌ（Ｘ，Ｙ）、右斜め相関領域の出力値をＬＵＴＸｒｉｇｈｔ＿ｄｉａｇｏｎａｌ（Ｘ，Ｙ）とする。
【０１２３】
画素値変換部４０１２では、入力画像データを走査し、５つの領域に対応するルックアップテーブルと領域判定データとを用いて画素値を変換する。走査中の注目画素ｉの画素値をＸｉとする。注目画素ｉに対して図２１（ａ）の位置関係にある画素の画素値ＹをＹｈｏｒｉｚｏｎｔａｌ＿ｉ、及びＹｆｌａｔ＿ｉとする。同様にして図２１（ｂ）の位置関係にある画素の画素値ＹをＹｖｅｒｔｉｃａｌ＿ｉ図２１（ｃ）の位置関係にある画素の画素値ＹをＹｌｅｆｔ＿ｄｉａｇｏｎａｌ＿ｉ、図２１（ｄ）の位置関係にある画素の画素値をＹｒｉｇｈｔ＿ｄｉａｇｏｎａｌ＿ｉとする。さらに注目画素の領域判定データをＰｆｌａｔ＿ｉ，Ｐｈｏｒｉｚｏｎｔａｌ＿ｉ，Ｐｖｅｒｔｉｃａｌ＿ｉ，Ｐｌｅｆｔ＿ｄｉａｇｏｎａｌ＿ｉ，Ｐｒｉｇｈｔ＿ｄｉａｇｏｎａｌ＿Ｉとする。画素値変換部４０１２は、Ｘｉは次式によりＸｉ’に変換される。
【０１２４】
【数１８】
0005832095.tif 000019
【０１２５】
上式は注目画素Ｘｉの各領域に属する確率で、各領域のルックアップテーブルの出力を重み付けしたものの総和が、注目画素の補正後の画素値Ｘｉ’であることを示している。
【０１２６】
最後に画素値変換部４０１２は変換後のＸｉ’をホワイトバランス制御部へ出力する。
【０１２７】
本実施例を総括する。実施例１及び２ではエッジ領域の画素及び平坦領域の画素を一まとめにして多次元ヒストグラムを構成していた。かかる方式による課題は、自然画像では平坦領域が多いため、回復ヒストグラムも平坦領域の頻度が大きくなり、ノイズ低減処理が平坦領域の影響を多く受け、処理結果が平坦になる傾向があることである。そこで本実施例では、予め平坦領域、及びエッジ領域の画素の特徴を決定し、領域毎に補正画素値を算出する。そして、注目画素が各領域に属する確率をもって前記補正画素値を加重平均して最終的な補正画素値とする。こうすることで、画像中の各領域の特徴に応じた適応的なノイズ低減処理を可能にしている。即ち、例えば、注目画素が縦方向のエッジに属する場合には縦相関領域の補正結果を重視する最終的な補正結果が得られる。一方で、縦相関領域の補正結果は平坦領域の影響が少ないので、入力画素は平坦化の度合いを抑えて補正される。つまり、本実施例のノイズ低減方法は実施例１及び２のノイズ低減方法に比べ、エッジ領域においてエッジがぼけにくいという利点がある。
【０１２８】
なお、本実施例では２次元ヒストグラムで議論を進めたが、３次元ヒストグラムや、その他次元のヒストグラムでも同様の手法が実施できる。一例として３次元ヒストグラム作成時の各領域の位置関係を図２２（ｅ）−（ｉ）に示す。図２２（ｅ）〜（ｉ）はそれぞれ、平坦領域、横相関領域、縦相関領域、左斜め相関領域、右斜め相関領域の２次元ヒストグラム作成に関するＸＹＺの位置関係である。平坦領域の位置関係は、図２２（ｅ）に示されている通り、横方向及び縦方向に関して密集度が高くなるような多次元ヒストグラムを生成することができる。
【０１２９】
＜他の実施例＞
また、本発明は、以下の処理を実行することによっても実現される。即ち、上述した実施形態の機能を実現するソフトウェア（プログラム）を、ネットワーク又は各種記憶媒体を介してシステム或いは装置に供給し、そのシステム或いは装置のコンピュータ（またはＣＰＵやＭＰＵ等）がプログラムを読み出して実行する処理である。
【図１】
0005832095.tif 000020
【図２】
0005832095.tif 000021
【図３】
0005832095.tif 000022
【図４】
0005832095.tif 000023
【図５】
0005832095.tif 000024
【図６】
0005832095.tif 000025
【図７】
0005832095.tif 000026
【図８】
0005832095.tif 000027
【図９】
0005832095.tif 000028
【図１０】
0005832095.tif 000029
【図１１】
0005832095.tif 000030
【図１２】
0005832095.tif 000031
【図１３】
0005832095.tif 000032
【図１４】
0005832095.tif 000033
【図１５】
0005832095.tif 000034
【図１６】
0005832095.tif 000035
【図１７】
0005832095.tif 000036
【図１８】
0005832095.tif 000037
【図１９】
0005832095.tif 000038
【図２０】
0005832095.tif 000039
【図２１】
0005832095.tif 000040
【図２２】
0005832095.tif 000041
【図２３】
0005832095.tif 000042
【図２４】
0005832095.tif 000043