特開 2023-70808 被写体撮影画像シミュレーション方法、画像生成装置およびプログラム。

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】公開特許公報(A)

(11)【公開番号】P2023070808

(43)【公開日】2023-05-22

(54)【発明の名称】被写体撮影画像シミュレーション方法、画像生成装置およびプログラム。

(51)【国際特許分類】

   G06T 15/06 20110101AFI20230515BHJP

   G06T 15/50 20110101ALI20230515BHJP

   A61B 1/045 20060101ALI20230515BHJP

【ＦＩ】

G06T15/06

G06T15/50

A61B1/045 610

【審査請求】未請求

【請求項の数】9

【出願形態】ＯＬ

(21)【出願番号】P 2021183156

(22)【出願日】2021-11-10

(71)【出願人】

【識別番号】000000376

【氏名又は名称】オリンパス株式会社

(74)【代理人】

【識別番号】110002147

【氏名又は名称】弁理士法人酒井国際特許事務所

(72)【発明者】

【氏名】関山 健太郎

(72)【発明者】

【氏名】中嶌 奨太

【テーマコード（参考）】

4C161

5B080

【Ｆターム（参考）】

4C161GG11

4C161WW20

5B080AA19

5B080DA06

5B080GA06

5B080GA11

(57)【要約】

【課題】表面散乱の弱い被写体の場合であっても、低ノイズで信頼性の高い被写体の撮影画像のシミュレーションを短時間で正確に生成することができる被写体撮影画像シミュレーション方法、画像生成装置およびプログラムを提供する。
【解決手段】被写体撮影画像シミュレーション方法は、被写体に入射し内部で散乱した後に被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算工程と、光線情報に基づいて、写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定工程と、光線情報と、輝度推定関数と、基づいて、被写体の輝度分布を算出する輝度分布算出工程と、被写体の輝度分布に基づいて、被写体の撮影画像を生成する生成工程と、を含む。
【選択図】図７

【特許請求の範囲】

【請求項1】

照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする被写体撮影画像シミュレーション方法であって、
前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算工程と、
前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定工程と、
前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出工程と、
前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成工程と、
を含む、
被写体撮影画像シミュレーション方法。

【請求項2】

請求項１に記載の被写体撮影画像シミュレーション方法であって、
前記輝度推定関数設定工程は、
前記輝度推定関数としてカーネル関数を設定する、
被写体撮影画像シミュレーション方法。

【請求項3】

請求項２に記載の被写体撮影画像シミュレーション方法であって、
前記カーネル関数は、
矩形関数である、
被写体撮影画像シミュレーション方法。

【請求項4】

請求項２に記載の被写体撮影画像シミュレーション方法であって、
前記カーネル関数は、
ガウス関数である、
被写体撮影画像シミュレーション方法。

【請求項5】

請求項１～４のいずれか一つに記載の被写体撮影画像シミュレーション方法であって、
前記光線追跡計算工程は、
前記照明系から発生する光線を前記被写体に照射し、前記被写体内に入射した光線が前記被写体の内部で散乱した後に前記被写体から出射する光線を計算して前記光線情報を生成し、
前記輝度推定関数設定工程は、
前記被写体表面から前記撮像系に届く光線の輝度を推定する、
被写体撮影画像シミュレーション方法。

【請求項6】

請求項１～４のいずれか一つに記載の被写体撮影画像シミュレーション方法であって、
前記光線追跡計算工程は、
前記撮像系から前記被写体に向かって前記被写体内部で散乱したのちに前記被写体から射出した光線を追跡し、
前記輝度推定関数設定工程は、
前記被写体の表面から前記撮像系に届く光線の輝度を推定し、
輝度分布算出工程は、前記被写体の内部で散乱した後に前記被写体から射出した光線を前記照明系まで追跡して求められた前記照明系の輝度分布から前記撮像系への寄与を算出する、
被写体撮影画像シミュレーション方法。

【請求項7】

請求項５または６に記載の被写体撮影画像シミュレーション方法であって、
前記光線追跡計算工程は、
前記撮像系が備える撮像レンズの入射瞳径と前記照明系が備える照明レンズの射出瞳径および前記撮像系が備える撮像素子の全画素のうちシミュレーションする画素数の比率を決定し、
前記比率に前記入射瞳径を乗じた値が前記射出瞳径よりも大きい場合、前記照明系から光線追跡をおこなう順方向追跡処理を実行する一方、
前記比率に前記入射瞳径を乗じた値が前記射出瞳径以下の場合、前記撮像系から光線追跡をおこなう逆方向追跡処理を実行する、
被写体撮影画像シミュレーション方法。

【請求項8】

照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする画像生成装置であって、
前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算部と、
前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定部と、
前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出部と、
前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成部と、
を備える、
画像生成装置。

【請求項9】

照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする画像生成装置が実行するプログラムであって、
前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算ステップと、
前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定ステップと、
前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出ステップと、
前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成ステップと、
を含む、
プログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本開示は、被写体撮影画像シミュレーション方法、画像生成装置およびプログラムに関する。

【背景技術】

【0002】

従来、内視鏡の設計において、内視鏡が撮影する撮影画像のシミュレーションが試みられている。また、内視鏡で観察する被写体の臓器は、表面のみで散乱が起こるのではなく、体積散乱体であるため、シミュレーションの撮影画像を正確に予測するには、被写体内部の光散乱を計算してシミュレーションの撮影画像に反映させる必要がある。この方法の一つとして、被写体の内部の光散乱を光線追跡により求めるモンテカルロ法により被写体内部の光散乱を計算することが知られているが、計算コストが膨大なため、撮影画像のシミュレーションにおける画像生成に用いることが難しかった。

【0003】

一方、撮影画像のシミュレーションにおける画像生成の効率化手段として表面フォトンマッピングの技術が知れている（例えば特許文献１参照）。この技術では、図２１のように、体積散乱体である皮膚の画像生成をする手法が開示されている。

【先行技術文献】

【特許文献】

【0004】

【特許文献1】特許第６０９２５６５号公報

【発明の概要】

【発明が解決しようとする課題】

【0005】

ところで、上述した特許文献１のような被写体が皮膚の場合、散乱が強いため、表面散乱特性がブロードとなるため、撮像系で高輝度の構成を大量に取得して画像計算を効率化することができる。

【0006】

しかしながら、内視鏡において観察する臓器等の被写体の表面は、粘膜等の液体で覆われていることが多く、表面での散乱が弱くなる。この表面散乱が弱く、出射直前の角度方向以外の光線の輝度が限りなくゼロに近い場合（表面散乱特性がデルタ関数の場合）、低輝度の光線ばかりが撮像系へ到達する（図２２を参照）。このため、上述した特許文献１の技術を用いた場合であっても、低輝度の光線ばかりの撮影画像は、ノイズが多く信頼性が低いことで、ノイズを減らすために膨大な試行回数の計算が必要となり、膨大な時間を要するという問題点があった。

【0007】

本開示は、上記に鑑みてなされたものであって、表面散乱の弱い被写体の場合であっても、低ノイズで信頼性の高い被写体の撮影画像のシミュレーションを短時間で正確に生成することができる被写体撮影画像シミュレーション方法、画像生成装置およびプログラムを提供することを目的とする。

【課題を解決するための手段】

【0008】

上述した課題を解決し、目的を達成するために、本開示に係る被写体撮影画像シミュレーション方法は、照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする被写体撮影画像シミュレーション方法であって、前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算工程と、前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定工程と、前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出工程と、前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成工程と、を含む。

【0009】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記輝度推定関数設定工程は、前記輝度推定関数としてカーネル関数を設定する。

【0010】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記カーネル関数は、矩形関数である。

【0011】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記カーネル関数は、ガウス関数である。

【0012】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記光線追跡計算工程は、前記照明系から発生する光線を前記被写体に照射し、前記被写体内に入射した光線が前記被写体の内部で散乱した後に前記被写体から出射する光線を計算して前記光線情報を生成し、前記輝度推定関数設定工程は、前記被写体表面から前記撮像系に届く光線の輝度を推定する。

【0013】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記光線追跡計算工程は、前記撮像系から前記被写体に向かって前記被写体内部で散乱したのちに前記被写体から射出した光線を追跡し、前記輝度推定関数設定工程は、前記被写体の表面から前記撮像系に届く光線の輝度を推定し、輝度分布算出工程は、前記被写体の内部で散乱した後に前記被写体から射出した光線を前記照明系まで追跡して求められた前記照明系の輝度分布から前記撮像系への寄与を算出する。

【0014】

また、本開示に係る被写体撮影画像シミュレーション方法は、上記開示において、前記光線追跡計算工程は、前記撮像系が備える撮像レンズの入射瞳径と前記照明系が備える照明レンズの射出瞳径および前記撮像系が備える撮像素子の全画素のうちシミュレーションする画素数の比率を決定し、前記比率に前記入射瞳径を乗じた値が前記射出瞳径よりも大きい場合、前記照明系から光線追跡をおこなう順方向追跡処理を実行する一方、前記比率に前記入射瞳径を乗じた値が前記射出瞳径以下の場合、前記撮像系から光線追跡をおこなう逆方向追跡処理を実行する。

【0015】

また、本開示に係る画像生成装置は、照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする画像生成装置であって、前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算部と、前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定部と、前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出部と、前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成部と、を備える。

【0016】

また、本開示に係るプログラムは、照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションする画像生成装置が実行するプログラムであって、前記被写体に入射し内部で散乱した後に前記被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する光線追跡計算ステップと、前記光線情報に基づいて、前記被写体表面から前記撮像系または前記照明系に届く光線の輝度および輝度推定関数を推定する輝度推定関数設定ステップと、前記光線情報と、前記輝度推定関数と、基づいて、前記被写体の輝度分布を算出する輝度分布算出ステップと、前記被写体の輝度分布に基づいて、前記被写体の撮影画像を生成する生成ステップと、を含む。

【発明の効果】

【0017】

本発明によれば、表面散乱の弱い被写体の場合であっても、低ノイズで信頼性の高い被写体の撮影画像のシミュレーションを短時間で正確に生成することができるという効果を奏する。

【図面の簡単な説明】

【0018】

【図1】図１は、一実施の形態に係る画像生成システムの概要を示す模式図である。

【図2】図２は、一実施の形態に係る被写体であるファントムの構造と光学特性との関係の一例を示す図である。

【図3】図３は、一実施の形態に係る照明系の条件の一例を示す図である。

【図4】図４は、一実施の形態に係る撮像系の条件の一例を示す図である。

【図5】図５は、一実施の形態に係る画像生成システムが光線情報からカーネル関数を決定する概要を模式的に示す図である。

【図6】図６は、一実施の形態に係る画像生成システムが光線にカーネル関数を適用した際に分散された構成の概要を模式的に示す図である。

【図7】図７は、一実施の形態に係る画像生成システムが実行する被写体撮影画像シミュレーション方法の処理の概要を示すフローチャートである。

【図8】図８は、一実施の形態に係る光線追跡計算部による被写体内部の散乱の計算を模式的に示す図である。

【図9】図９は、一実施の形態に係る光線情報記録部が記録する光線情報の一部を示す図である。

【図10】図１０は、一実施の形態に係るカーネル関数Ｋ（ｔ）の一例を示す図である。

【図11】図１１は、一実施の形態に係る輝度推定関数設定部が照明系に到達する光線から撮像系の撮像素子における１画素の照度を算出する際の算出方法を模式的に示す図である。

【図12】図１２は、一実施の形態に係る輝度推定関数設定部が撮像系の撮像素子における１画素の照度を算出する際の算出方法を模式的に示す図である。

【図13】図１３は、一実施の形態に係る照度Ｂｉａｓの誤差とカーネル関数角度幅との関係を示す図である。

【図14】図１４は、一実施の形態に係る順方向追跡処理の概要を示すフローチャートである。

【図15】図１５は、一実施の形態に係る順方向追跡処理の概要を模式的に示す図である。

【図16】図１６は、一実施の形態に係る生成部が生成した撮影画像の一例を示す図である。

【図17】図１７は、一実施の形態に係る逆方向追跡処理の概要を模式的に示すフローチャートである。

【図18】図１８は、一実施の形態に係る逆方向追跡処理の概要を模式的に示す図である。

【図19】図１９は、一実施の形態に係る実施例（本計算方法）および比較例（従来計算方法）を用いて各画素の推定輝度から撮像素子における照度を計算したラインプロファイルの比較結果を示す図である。

【図20】図２０は、一実施の形態の変形例１に係るカーネル関数の別の一例を示す図である。

【図21】図２１は、従来技術の被写体の表面散乱が強い場合の一例を模式的に示す図である。

【図22】図２２は、従来技術の被写体の表面散乱が弱い場合の一例を模式的に示す図である。

【発明を実施するための形態】

【0019】

以下、本開示を実施するための形態を図面とともに詳細に説明する。なお、以下の実施の形態により本開示が限定されるものでない。また、以下の説明において参照する各図は、本開示の内容を理解でき得る程度に形状、大きさ、および位置関係を概略的に示してあるに過ぎない。即ち、本開示は、各図で例示された形状、大きさおよび位置関係のみに限定されるものではない。さらに、以下の説明では、被写体撮影画像シミュレーション方法を実行可能な画像生成装置を備える画像生成システムについて説明する。

【0020】

〔画像生成システムの概要〕
図１は、一実施の形態に係る画像生成システムの概要を示す模式図である。図１に示す画像生成システム１は、照明系、撮像系および体積散乱体の被写体で構成されたモデルにおける撮影画像をシミュレーションするものである。画像生成システム１は、体積散乱体の模擬生体（ファントム）である被写体２と、照明系３と、撮像系４と、画像生成装置５と、を備える。

【0021】

被写体２（ファントム）は、均一で吸収のない散乱体である組織の中に、完全吸収体である１本の血管を配置したものである。

【0022】

図２は、被写体２であるファントムの構造と光学特性との関係の一例を示す図である。図２に示す表Ｔ１には、被写体２であるファントムにおける組織および血管の構造と、被写体２であるファントムにおける組織および血管の各々の光学特性と、を対応付けて記載されている。

【0023】

照明系３は、被写体２から所定の距離に配置され、被写体２に向けて光を照射する。照明系３は、少なくとも、光を発する光源と、光源が発した光を集光して被写体２に向けて照射する照明レンズと、を有する。

【0024】

図３は、照明系３の条件の一例を示す図である。図３に示す表Ｔ２には、照明系３における各種の照明条件と、観察距離と、を対応付けて記載されている。具体的には、表Ｔ２には、照明系３の形状と、サイズと、配光と、観察距離と、を対応付けて記載されている。

【0025】

撮像系４は、被写体２から所定の距離に配置され、被写体２からの光を受光して撮像する。撮像系４は、少なくとも、被写体２からの光を集光して被写体像を結像する撮像レンズと、撮像レンズによって結像された被写体像を受光する撮像素子（ＣＭＯＳ（Complementary Metal Oxide Semiconductor）センサまたはＣＣＤ（Charge Coupled Device）センサ）と、を有する。

【0026】

図４は、撮像系４の条件の一例を示す図である。図４に示す表Ｔ３には、撮像系４における各種の撮像条件と、観察距離と、を対応付けて記載されている。具体的には、表Ｔ３には、撮像素子のサイズと、画素数と、観察距離と、を対応付けて記載されている。

【0027】

画像生成装置５は、撮影画像をシミュレーションする。画像生成装置５は、通信部５１と、表示部５２と、入力部５３と、記録部５４と、制御部５５と、を備える。

【0028】

通信部５１は、制御部５５から入力された各種データを照明系３および撮像系４へ出力し、撮像系４から入力された被写体撮影画像データを制御部５５へ出力する。通信部５１は、通信インターフェース等を用いて構成される。

【0029】

表示部５２は、制御部５５の制御のもと、画像生成システム１に関する各種情報および被写体撮影画像を表示する。表示部５２は、有機ＥＬディスプレイ（Organic Electroluminescent Display）および液晶ディスプレイ等を用いて構成される。

【0030】

入力部５３は、操作者の各種の操作の入力を受け付け、受け付けた操作に応じた信号を制御部５５へ出力する。入力部５３は、キーボード、マウスおよびタッチパネル等を用いて構成される。

【0031】

記録部５４は、画像生成システム１に関する各種の情報を記録する。記録部５４は、ＲＡＭ（Random Access Memory）、ＲＯＭ（Read Only Memory）、ＨＤＤ（Hard Disk Drive）およびＳＳＤ（Solid State Drive）等を用いて構成される。記録部５４は、画像生成装置５が実行する各種のプログラムを記録するプログラム記録部５４１と、後述する光線情報を記録する光線情報記録部５４２と、誤差発生量テーブルデータを記録する誤差発生量テーブルデータ記録部５４３と、を有する。

【0032】

制御部５５は、画像生成システム１を構成する各部を制御する。制御部５５は、メモリと、ＣＰＵ（Central Processing Unit）およびＧＰＵ（Graphics Processing Unit）等のハードウェアを有するプロセッサを用いて構成される。制御部５５は、光線追跡計算部５５１と、輝度推定関数設定部５５２と、輝度分布算出部５５３と、生成部５５４と、を有する。

【0033】

光線追跡計算部５５１は、被写体に入射し内部で散乱した後に被写体表面から出射する光線を計算して少なくとも光線の角度分布を含む光線情報を生成する。ここで、光線の角度分布には、散乱体である被写体の表面から出射する出射光線の位置、出射光線の方向余弦および出射光線の光束が含まれる。

【0034】

輝度推定関数設定部５５２は、光線追跡計算部５５１が生成した光線情報に基づいて、被写体表面から撮像系４または照明系３に届く光線の輝度および輝度推定関数を推定する。具体的には、輝度推定関数設定部５５２は、輝度推定関数としてカーネル関数を推定する。

【0035】

輝度分布算出部５５３は、光線追跡計算部５５１が生成した光線情報と、輝度推定関数設定部５５２が推定した輝度推定関数と、基づいて、被写体の輝度分布を算出する。

【0036】

生成部５５４は、輝度分布算出部５５３が算出した被写体の輝度分布基づいて、被写体の撮影画像を生成する。

【0037】

〔画像生成システム１による処理の概要〕
次に、画像生成システム１が実行する被写体撮影画像シミュレーション方法の概要について説明する。図５は、画像生成システム１が光線情報からカーネル関数を決定する概要を模式的に示す図である。図６は、画像生成システム１が光線にカーネル関数を適用した際に分散された構成の概要を模式的に示す図である。

【0038】

図５に示すように、画像生成システム１は、被写体２表面の散乱特性ではなく、被写体２内部で散乱後に表面から出射する光線Ｒ１（以下、「出射光線Ｒ１」という）の光線群Ｄ１である光線情報(少なくとも、出射光線Ｒ１の位置、出射光線Ｒ１の方向、出射光線Ｒ１の光束を含む出射光線Ｒ１の角度分布を含む)に基づいて、カーネル関数Ｋ（ｔ）を決定する。具体的には、図６に示すように、被写体２内部で多重散乱した出射光線Ｒ１は、散乱ごとに光線方向が一様になり、被写体２表面からの射出時に適度に分散した出射光線Ｒ２となる。

【0039】

従って、画像生成システム１は、被写体２内部で散乱後に表面から出射する出射光線Ｒ１の光線群Ｄ１である光線情報に基づいて決定したカーネル関数Ｋ（ｔ）を用いて出射光線Ｒ１を分散して出射し、撮像系４へ到達させることで、出射光線Ｒ２を大量に撮影画像へ寄与させ、計算の効率化を実現する。

【0040】

これにより、画像生成システム１は、被写体２表面の散乱特性(ＢＴＤＦ)ではなく、被写体２の内部で散乱した後に被写体２表面から射出する出射光線Ｒ１の光線情報から求めたカーネル関数Ｋ（ｔ）を用いて被写体２表面から撮像系４または照明系３に届く光線の輝度の推定を行う。この結果、画像生成システム１は、消化管の表面等で液面に覆われており、表面の散乱特性が非常に弱い被写体２の場合であっても、光線情報が適度に分散していることによって、撮像系４に入射する成分が多くなり、輝度推定関数としてのカーネル関数Ｋ（ｔ）の関数幅が推定に適した値となる。

【0041】

一般的な内視鏡では、被写体２を照明する場合、光源が小さいため、被写体２表面における照明光のＮＡも小さくなり、入射光の角度が局在している。しかしながら、被写体２は、体積散乱体なので、光が散乱するごとに光線の方向が一様になっていく。従って、画像生成システム１によれば、入射光線の方向の指向性が非常に強い場合でも、内部で散乱を繰り返すことにより被写体２で散乱した後に被写体２表面での光線情報が適度に分散するので、カーネル関数Ｋ（ｔ）の関数幅が分布をもった推定に適した値となる。この結果、撮像系４に入射する出射光線Ｒ２の輝度が高くなり、計算の効率を向上することができる。

【0042】

このように画像生成システム１によれば、被写体２の表面散乱が非常に弱い場合であっても、上述した輝度推定関数設定方法を用いた画像生成を行うことで、臓器等の表面が液体に覆われ被写体２表面の散乱特性が非常に弱い生体の撮影画像のシミュレーションを正確かつ効率的に行うことができる。

【0043】

〔画像生成システムが実行する処理〕
次に、画像生成システム１が実行する被写体撮影画像シミュレーション方法の処理について説明する。図７は、画像生成システム１が実行する被写体撮影画像シミュレーション方法の処理の概要を示すフローチャートである。

【0044】

図７に示すように、まず、光線追跡計算部５５１は、照明系３、撮像系４および被写体モデルを決定する（ステップＳ１０１）。具体的には、光線追跡計算部５５１は、上述した図２～図４における表Ｔ１～Ｔ３の照明系３、撮像系４および被写体モデルを決定する。

【0045】

続いて、光線追跡計算部５５１は、撮像系４における撮像素子の軸上１画素から射出瞳方向へ光線追跡を行い（ステップＳ１０２）、被写体２内部の散乱を計算し（ステップＳ１０３）、被写体２である散乱体表面から出射する出射光線を光線情報として光線情報記録部５４２に記録する（ステップＳ１０４）。具体的には、光線追跡計算部５５１は、被写体２である散乱体表面から出射する出射光線の角度分布を含む光線情報を光線情報記録部５４２に記録する。ここで、出射光線の角度分布には、散乱体である被写体２の表面２１から出射する出射光線の位置（ｘ、ｙ、ｚ）、出射光線の方向余弦α、β、γおよび出射光線の光束Φが含まれる。

【0046】

図８は、光線追跡計算部５５１による被写体２内部の散乱の計算を模式的に示す図である。図８に示すように、光線追跡計算部５５１は、撮像系４における撮像素子の軸上１画素から射出瞳方向へ光線Ｒ３の追跡し、被写体２内部の散乱を計算後、散乱体である被写体２の表面２１から出射する出射光線Ｒ１の位置Ｐ１（ｘ、ｙ、ｚ）、出射光線Ｒ１の方向余弦α、β、γおよび出射光線Ｒ１の光束Φを光線情報として生成して光線情報記録部５４２に記録する。

【0047】

図９は、光線情報記録部５４２が記録する光線情報の一部を示す図である。図９に示す光線情報Ｔ４には、被写体２の表面２１における出射光線Ｒ１の位置Ｐ１であるｘ、ｙ、ｚ、出射光線Ｒ１の方向余弦α、β、γおよび出射光線Ｒ１の光束Φが含まれる。また、光線追跡計算部５５１が追跡する光線本数は、目標物（ゴール）である照明系３へ１０，０００本到達するまで（統計誤差１％）とする。このため、一実施の形態では、光線追跡計算部５５１が撮像系４から３００，０００本の光線を追跡するものとする。さらに、一実施の形態では、光線追跡計算部５５１による光線追跡をＬｉｇｈｔＴｏｏｌｓ（登録商標）（Ｓｙｎｏｐｓｙｓ，Ｉｎｃ.（登録商標））のＶｅｒ８.５を用いて行った。

【0048】

ステップＳ１０４の後、輝度推定関数設定部５５２は、光線情報記録部５４２が記録する光線情報をそのままに照明系３へ到達する光線と、照明系３の輝度と、に基づいて、撮像系４の撮像素子における１画素の照度を算出する（ステップＳ１０５）。

【0049】

続いて、輝度推定関数設定部５５２は、カーネル関数を誤差発生量テーブルデータ記録部５４３が記録する人臓器散乱体の初期値に設定する（ステップＳ１０６）。ここで、カーネル関数角度幅の初期値の決定方法は、過去のヒト臓器の計算実施例を用いて、カーネル角度幅に対する誤差発生量を、テーブルデータ化して誤差発生量テーブルデータ記録部５４３に記録したものである。即ち、輝度推定関数設定部５５２は、誤差発生量テーブルデータ記録部５４３が記録する誤差発生量テーブルデータを参照し、目標とする誤差量に対して適当な関数幅を初期値として設定する。これにより、目標値から大きく外れていない範囲で角度幅の最適化を行うことができる。

【0050】

図１０は、カーネル関数Ｋ（ｔ）の一例を示す図である。図１０に示すように、カーネル関数Ｋ（ｔ）の形状は、矩形関数を使用し、以下の式（１）を満たす。また、カーネル関数Ｋ（ｔ）は、矩形関数のため、関数の幅Δが設定すべき変数となる。

【数1】

【0051】

ステップＳ１０６の後、輝度推定関数設定部５５２は、光線情報記録部５４２が記録する光線情報に基づいて、光線情報そのままに照明系３へ到達する光線を追跡し、照明系３に到達する光線から撮像系４の撮像素子における１画素の照度を算出する（ステップＳ１０７）。

【0052】

図１１は、輝度推定関数設定部５５２が照明系３に到達する光線から撮像系４の撮像素子における１画素の照度を算出する際の算出方法を模式的に示す図である。図１１に示すように、輝度推定関数設定部５５２は、光線情報記録部５４２が記録する光線情報をそのままに照明系３へ到達する光線Ｒ４を追跡し、照明系３に到達する光線Ｒ４から撮像系４の撮像素子における１画素の照度を算出する。

【0053】

ステップＳ１０７の後、輝度推定関数設定部５５２は、あるカーネル関数を仮定し、それぞれの光線をカーネル関数に応じて分散させ（ステップＳ１０８）、照明系３に到達する分散光線と照明系３の輝度とから撮像系４の撮像素子における１画素の照度を算出する（ステップＳ１０９）。

【0054】

図１２は、輝度推定関数設定部５５２が撮像系４の撮像素子における１画素の照度を算出する際の算出方法を模式的に示す図である。図１２に示すように、輝度推定関数設定部５５２は、あるカーネル関数Ｋ（ｔ）を仮定し、それぞれの光線をカーネル関数Ｋ（ｔ）に応じて分散させた出射光線Ｒ２（分散光線）と、照明系３に届く光線の輝度と、に基づいて、撮像系４の撮像素子における１画素の照度を算出する。なお、出射光線Ｒ２（分散光線）については、本実施例においてはカーネル関数幅の範囲でランダムに光線を１０００本生成し、光線１本のエネルギーを１／１０００とした。エネルギー保存の観点から、生成する光線本数と、エネルギーの減少率は逆数の関係になっている必要がある。

【0055】

ステップＳ１０９の後、輝度推定関数設定部５５２は、上述したステップＳ１０７において算出した撮像系４の撮像素子における１画素の照度と、上述したステップＳ１０９において算出した撮像系４の撮像素子における１画素の照度との誤差（以下、「照度Ｂｉａｓの誤差」と表記する）を算出する（ステップＳ１１０）。この場合、輝度推定関数設定部５５２は、各光線が内部で複数回散乱し、出射角度が均一に近くなるため、ステップＳ１０８においてカーネル関数Ｋ（ｔ）を適用して光線をある角度範囲で分散しても、照度Ｂｉａｓの誤差が抑えられたまま、高輝度光線を効率よく取得することができる。

【0056】

続いて、輝度推定関数設定部５５２は、照度Ｂｉａｓの誤差が予め設定した許容量以下であるか否かを判定する（ステップＳ１１１）。具体的には、輝度推定関数設定部５５２は、照度Ｂｉａｓの誤差が予め設定した許容量±５％であるか否かを判定する。この許容量は、ユーザが適宜設定することができる。例えば、ユーザは、照度Ｂｉａｓの誤差Ｓｉｍと比較する際に、予め設定した許容量が実際に用いる計測器の誤差量内に抑えることができる目標を設定する。輝度推定関数設定部５５２は、照度Ｂｉａｓの誤差が予め設定した許容量以下であると判定した場合（ステップＳ１１１：Ｙｅｓ）、カーネル関数Ｋ（ｔ）を最終値（以下、「カーネル関数Ｋ（θ）」と表記する）に決定する（ステップＳ１１２）。ステップＳ１１２の後、画像生成システム１は、後述するステップＳ１１５へ移行する。

【0057】

ステップＳ１１１において、輝度推定関数設定部５５２は、照度Ｂｉａｓの誤差が予め設定した許容量を超えると判定した場合（ステップＳ１１１：Ｎｏ）、画像生成システム１は、ステップＳ１１２へ移行する。

【0058】

続いて、輝度推定関数設定部５５２は、誤差発生量テーブルデータ記録部５４３が記録する誤差発生量テーブルデータを参照し、現状の照度Ｂｉａｓの誤差量からカーネル関数角度幅の変更量を決定し（ステップＳ１１３）、カーネル関数を変更して設定する（ステップＳ１１４）。ステップＳ１１４の後、画像生成システム１は、上述したステップＳ１０７へ戻り、照度Ｂｉａｓの誤差が予め設定した許容量以下となるまで、ステップＳ１０７～ステップＳ１０９を繰り返す。

【0059】

図１３は、照度Ｂｉａｓの誤差とカーネル関数角度幅との関係を示す図である。図１３において、縦軸が照度Ｂｉａｓの誤差を示し、横軸がカーネル関数角度幅［ｄｅｇ］を示す。また、図１３において、折れ線Ｌ１が照度Ｂｉａｓの誤差とカーネル関数角度幅との関係を示す。

【0060】

図１３の折れ線Ｌ１に示すように、輝度推定関数設定部５５２は、カーネル関数角度幅の初期値を６３[ｄｅｇ]とし、許容量以上の照度Ｂｉａｓの誤差が発生する毎に、角度幅を９[ｄｅｇ]減らして再計算を行う。図１３の折れ線Ｌ１に示すように、輝度推定関数設定部５５２は、角度幅が４５[ｄｅｇ]において、照度Ｂｉａｓの誤差が許容量の５％を達成したため、これをカーネル関数角度幅として採用する。このように、輝度推定関数設定部５５２は、カーネル関数幅をより小さく抑え、再設定したカーネル関数に応じて光線を分散させ、照度を算出し、再度、照度Ｂｉａｓの誤差を評価し、許容量との比較を行う。そして、輝度推定関数設定部５５２は、最終的に、照度Ｂｉａｓの誤差が許容量に収まるカーネル関数幅を決定する。

【0061】

ステップＳ１１１の後、光線追跡計算部５５１は、撮像系４の撮像レンズの入射瞳径Ｐ_ｉｎと照明系３の照明レンズの出射瞳径Ｐ_ｏｕｔとの比を算出する（ステップＳ１１５）。

【0062】

続いて、光線追跡計算部５５１は、撮像系４における撮像素子の全画素のうち、シミュレーションする画素数の比率αを決定する（０＜α≦１）（ステップＳ１１６）。

【0063】

その後、光線追跡計算部５５１は、比率αに撮像系４の撮像の入射瞳径Ｐ_ｉｎを乗じた値（α×Ｐ_ｉｎ）が照明系３の照明レンズの出射瞳径Ｐ_ｏｕｔより大きいか否かを判定する（ステップＳ１１７）。光線追跡計算部５５１が比率αに撮像系４の撮像レンズの入射瞳径Ｐ_ｉｎを乗じた値（α×Ｐ_ｉｎ）が照明系３の照明レンズの出射瞳径Ｐ_ｏｕｔより大きいと判定した場合（ステップＳ１１７：Ｙｅｓ）、画像生成システム１は、後述するステップＳ１１８へ移行する。これに対して、光線追跡計算部５５１が比率αに撮像系４の撮像レンズの入射瞳径Ｐ_ｉｎを乗じた値（α×Ｐ_ｉｎ）が照明系３の照明レンズの出射瞳径Ｐ_ｏｕｔより大きくないと判定した場合（ステップＳ１１７：Ｎｏ）、画像生成システム１は、後述するステップＳ１１９へ移行する。

【0064】

ステップＳ１１８において、照明系３から光線追跡を開始する順方向追跡処理を実行する。なお、順方向追跡処理の詳細は、後述する。ステップＳ１１８の後、画像生成システム１は、本処理を終了する。

【0065】

ステップＳ１１９において、撮像系４から光線を飛ばし照明系３まで光線追跡を開始する逆方向追跡処理を実行する。なお、逆方向追跡処理の詳細は、後述する。ステップＳ１１９の後、画像生成システム１は、本処理を終了する。

【0066】

〔順方向追跡処理の詳細〕
次に、上述したステップＳ１１８で説明した順方向追跡処理の詳細について説明する。図１４は、順方向追跡処理の概要を示すフローチャートである。図１５は、順方向追跡処理の概要を模式的に示す図である。ここで、順方向追跡処理とは、照明系３から光線追跡を行う手法の処理である。

【0067】

図１４に示すように、まず、輝度推定関数設定部５５２は、光線情報記録部５４２が記録する光線情報に基づいて、照明系３から光線追跡し（ステップＳ２０１）、内部散乱を計算する（ステップＳ２０２）。

【0068】

続いて、輝度推定関数設定部５５２は、散乱体である被写体２の生体内部で散乱した後に、生体である被写体２の表面から出射する際に、上述したステップＳ１１２で決定したカーネル関数Ｋ（θ）に応じて光線を分散する（ステップＳ２０３）。具体的には、図１５に示すように、輝度推定関数設定部５５２は、散乱体である被写体２内部で散乱した後に、被写体２の表面から出射する光線に対してカーネル関数Ｋ（θ）によって分散させる。なお、光線の分散方法は、上述したステップＳ１０８のカーネル関数設定時と同様のため（図１２を参照）、詳細な説明は省略する。

【0069】

その後、輝度分布算出部５５３は、輝度推定関数設定部５５２がカーネル関数Ｋ（θ）を用いて分散した光線を撮像系４まで追跡し、撮像系４の撮像素子の各画素の照度を算出する（ステップＳ２０４）。具体的には、輝度分布算出部５５３は、分散した光線を撮像系４まで追跡し、撮像系４の撮像素子の各画素の照度を算出することによって、被写体２の輝度分布を算出する。これにより、例えば、内視鏡では、撮像レンズにより撮像素子へと光線が結像され、各画素の照度（照度値）の輝度分布をシミュレーションすることができる。

【0070】

続いて、生成部５５４は、撮像系４の撮像素子に到達した光線に基づいて、被写体２の撮影画像を生成する（ステップＳ２０５）。具体的には、生成部５５４は、輝度分布算出部５５３が算出した被写体の輝度分布に対応する撮像系４の撮像素子に到達した光線に基づいて、被写体２の撮影画像を生成する。ステップＳ２０５の後、画像生成システム１は、図７のメインルーチンへ戻り、本処理を終了する。

【0071】

図１６は、生成部５５４が生成した撮影画像の一例を示す図である。図１６に示すように、生成部５５４は、輝度分布算出部５５３が算出した被写体の輝度分布に対応する撮像系４の撮像素子に到達した光線に基づいて、被写体撮影画像をシミュレーションした撮影画像ｉ１を生成する。なお、図１６における撮影画像ｉ１では、画素数が７００×５００（Pixel）の撮像素子を用いた一例を示す。さらに、一実施の形態では、ＣＰＵのコア数が１２コアの制御部５５によって、ピーク統計誤差３％の撮影画像ｉ１を約８時間で生成することができた。これに対して、従来技術のモンテカルロ計算では、同様の計算に約２０００日が必要であり、約６，０００倍の効率化を達成することができた。

【0072】

〔逆方向追跡処理の詳細〕
次に、上述したステップＳ１１９で説明した逆方向追跡処理の詳細について説明する。図１７は、逆方向追跡処理の概要を模式的に示すフローチャートである。図１８は、逆方向追跡処理の概要を模式的に示す図である。ここで、逆方向追跡処理とは、撮像系４側から光線を飛ばし、照明系３まで追跡する手法の処理である。即ち、逆方向追跡処理は、光線の通過したパスの情報をそのまま活用し、照明系３へ到達した際の位置および角度と、照明系３の配光特性を照合して、撮像系４の撮像素子における１画素の照度を算出する。

【0073】

図１７に示すように、まず、輝度推定関数設定部５５２は、光線情報記録部５４２が記録する光線情報に基づいて、撮像系４の撮像素子における各画素から光線追跡し（ステップＳ３０１）、内部散乱を計算する（ステップＳ３０２）。

【0074】

続いて、輝度推定関数設定部５５２は、散乱体である被写体２の生体内部で散乱した後に、散乱体である被写体２表面から出射する際に、上述したステップＳ１１２で設定したカーネル関数Ｋ（θ）に応じて光線を分散する（ステップＳ３０３）。具体的には、図１８に示すように、輝度推定関数設定部５５２は、散乱体である被写体２内部で散乱した後に、被写体２表面から出射する光線に対してカーネル関数Ｋ（θ）によって分散させる。なお、光線の分散方法は、上述したステップＳ１０８のカーネル関数設定時と同様のため（図１２を参照）、詳細な説明は省略する。

【0075】

その後、輝度分布算出部５５３は、撮像系４側から光線を飛ばし、照明系３まで追跡し、光線の通過したパスの情報をそのまま活用し、照明系３へ到達した際の位置および角度と、照明系３の配光特性と、を照合して、撮像系４の照度を逆算する（ステップＳ３０４）。具体的には、輝度分布算出部５５３は、撮像系４側から光線を飛ばし、照明系３まで追跡し、照明系３へ到達した際の位置および角度と、照明系３の配光特性と、を照合して、撮像系４の照度を逆算することによって、被写体２の輝度分布を算出する。この場合、輝度分布算出部５５３は、照明系３が持つ配光分布と面内分布を考慮するため、照明系３への光線の入射位置および入射角度の情報を予め把握しておくことが好ましい。

【0076】

続いて、生成部５５４は、撮像系４の撮像素子に到達した光線に基づいて、画像を生成する（ステップＳ３０５）。具体的には、生成部５５４は、輝度分布算出部５５３が算出した被写体の輝度分布に対応する撮像系４の撮像素子に到達した光線に基づいて、被写体２の撮影画像を生成する。ステップＳ３０５の後、画像生成システム１は、図７のメインルーチンへ戻り、本処理を終了する。

【0077】

このように、逆方向追跡処理によれば、内視鏡のように受光器（撮像レンスおよび撮像素子）の面積およびＮＡが極端に小さい場合、受光器よりも相対的に大きい光源を光線追跡の目標物（ゴール）となり、撮影画像に寄与する光線を取り込みやすくなる。

【0078】

さらに、逆方向追跡処理によれば、撮像系４における撮像素子の特定の画素が受光する照度を選択的に計算することができ、必要な画素数に応じて計算時間を短縮できる。

【0079】

（実施例）
次に、実施例と比較例とを参照して、上記効果についてより具体的に説明する。図１９は、実施例（本計算方法）および比較例（従来計算方法）を用いて各画素の推定輝度から撮像素子における照度を計算したラインプロファイルの比較結果を示す図である。なお、図１９では、２００画素１ラインの撮像素子を用いた場合について説明する。さらに、図１９において、曲線Ｌ１１が実施例（本計算方法）を示し、曲線Ｌ１２が比較例（従来計算方法）を示す。なお、図１９では、制御部５５のＣＰＵのコア数が１２個の場合で、約２０分間計算した結果を示す。

【0080】

図１９に示すように、比較例（従来来法）は、ピーク統計誤差約９.４％である（曲線Ｌ１２を参照）。これに対して、実施例は、ピーク統計誤差約０.５％であり（曲線Ｌ１１を参照）、比較例に対して約３５０倍の効率化を達成している。

【0081】

以上説明した一実施の形態によれば、表面散乱の弱い被写体２の場合であっても、低輝度の光線ばかりが撮像系４へ到達することなく、高輝度の光線を大量に取得することがで、低ノイズで信頼性の高い被写体２の撮影画像ｉ１のシミュレーションを短時間で正確に生成することができる。

【0082】

（変形例１）
また、一実施の形態では、カーネル関数が矩形状をなしていたが、これに限定されることなく、種々の形状を適用することができる。図２０は、一実施の形態の変形例１に係るカーネル関数の別の一例を示す図である。図２０に示すように、角度関数としてガウス関数であっても適用することができる。図２０のガウス関数は、σ＝３０ｄｅｇの場合のものである。

【0083】

（その他の実施の形態）
上述した本開示の一実施の形態に係る画像生成システムに開示されている複数の構成要素を適宜組み合わせることによって、種々の発明を形成することができる。例えば、上述した本開示の実施の形態に係る画像生成システムに記載した全構成要素からいくつかの構成要素を削除してもよい。さらに、上述した本開示の実施の形態に係る画像生成システムで説明した構成要素を適宜組み合わせてもよい。

【0084】

また、本開示の一実施の形態に係る画像生成システムでは、上述してきた「部」は、「手段」や「回路」などに読み替えることができる。例えば、制御部は、制御手段や制御回路に読み替えることができる。

【0085】

また、本開示の一実施の形態に係る画像生成システムに実行させるプログラムは、インストール可能な形式または実行可能な形式のファイルデータでＣＤ－ＲＯＭ、フレキシブルディスク（ＦＤ）、ＣＤ－Ｒ、ＤＶＤ（Digital Versatile Disk）、ＵＳＢ媒体、フラッシュメモリ等のコンピュータで読み取り可能な記録媒体に記録されて提供される。

【0086】

また、本開示の一実施の形態に係る画像生成システムに実行させるプログラムは、インターネット等のネットワークに接続されたコンピュータ上に格納し、ネットワーク経由でダウンロードさせることにより提供するように構成してもよい。

【0087】

なお、本明細書におけるフローチャートの説明では、「まず」、「その後」、「続いて」等の表現を用いてステップ間の処理の前後関係を明示していたが、本発明を実施するために必要な処理の順序は、それらの表現によって一意的に定められるわけではない。即ち、本明細書で記載したフローチャートにおける処理の順序は、矛盾のない範囲で変更することができる。また、こうした、単純な分岐処理からなるプログラムに限らず、より多くの判定項目を総合的に判定して分岐させてもよい。

【0088】

以上、本願の実施の形態のいくつかを図面に基づいて詳細に説明したが、これらは例示であり、本発明の開示の欄に記載の態様を始めとして、当業者の知識に基づいて種々の変形、改良を施した他の形態で本発明を実施することが可能である。

【符号の説明】

【0089】

１ 画像生成システム
２ 被写体
３ 照明系
４ 撮像系
５ 画像生成装置
５１ 通信部
５２ 表示部
５３ 入力部
５４ 記録部
５５ 制御部
５６ ：輝度推定関数設定部
５４１ プログラム記録部
５４２ 光線情報記録部
５４３ 差発生量テーブルデータ記録部
５５１ 光線追跡計算部
５５２ 輝度推定関数設定部
５５３ 輝度分布算出部
５５４ 生成部

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】

【図17】

【図18】

【図19】

【図20】

【図21】

【図22】