JP2015192238A - 画像データ生成装置および画像データ生成方法 - Google Patents

画像データ生成装置および画像データ生成方法 Download PDF

Info

Publication number
JP2015192238A
JP2015192238A JP2014067033A JP2014067033A JP2015192238A JP 2015192238 A JP2015192238 A JP 2015192238A JP 2014067033 A JP2014067033 A JP 2014067033A JP 2014067033 A JP2014067033 A JP 2014067033A JP 2015192238 A JP2015192238 A JP 2015192238A
Authority
JP
Japan
Prior art keywords
image data
viewpoint
scattered light
scattered
filter
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
JP2014067033A
Other languages
English (en)
Inventor
村上 友近
Tomochika Murakami
友近 村上
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Canon Inc
Original Assignee
Canon Inc
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Canon Inc filed Critical Canon Inc
Priority to JP2014067033A priority Critical patent/JP2015192238A/ja
Priority to US14/659,957 priority patent/US20150279033A1/en
Publication of JP2015192238A publication Critical patent/JP2015192238A/ja
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • G06T7/0012Biomedical image inspection
    • G06T7/0014Biomedical image inspection using an image reference approach
    • GPHYSICS
    • G02OPTICS
    • G02BOPTICAL ELEMENTS, SYSTEMS OR APPARATUS
    • G02B21/00Microscopes
    • G02B21/0004Microscopes specially adapted for specific applications
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/50Image enhancement or restoration by the use of more than one image, e.g. averaging, subtraction
    • G06T5/73
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/136Segmentation; Edge detection involving thresholding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/60Analysis of geometric attributes
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10024Color image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10052Images from lightfield camera
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10056Microscopic image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10141Special mode during image acquisition
    • G06T2207/10148Varying focus
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20048Transform domain processing
    • G06T2207/20056Discrete and fast Fourier transform, [DFT, FFT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20076Probabilistic image processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing
    • G06T2207/30024Cell structures in vitro; Tissue sections in vitro
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing
    • G06T2207/30056Liver; Hepatic
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing
    • G06T2207/30096Tumor; Lesion

Abstract

【課題】被写体を撮影して得られた元画像データから、画像処理によって、被写体の観察や診断に適した画像データを生成するための新規な技術を提供する。
【解決手段】被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを、周波数空間のデータに変換する変換手段と、前記変換された周波数空間のデータに対し、フィルタを適用するフィルタ適用手段と、前記フィルタが適用されたデータを実空間の画像データに逆変換する逆変換手段と、を有する。前記フィルタは、前記逆変換された画像データに含まれるレイヤー画像データのいずれかが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されている。
【選択図】図23

Description

本発明は、被写体を撮影して得られた画像データから観察や診断に適した画像データを生成する画像データ生成装置および画像データ生成方法に関する。
病理分野において、病理診断のツールである光学顕微鏡の代替として、プレパラートに載置された被検試料を撮像しデジタル化してディスプレイ上での病理診断を可能とするバーチャル・スライド・システムがある。バーチャル・スライド・システムによる病理診断画像のデジタル化により、従来の被検試料の光学顕微鏡像をデジタルデータとして取り扱える。それによって、遠隔診断の迅速化、デジタル画像データを使った患者への説明、希少症例の共有化、教育・実習の効率化、などのメリットが得られる。
またデジタルデータに対しては様々な画像処理が可能で、バーチャル・スライド・システムで撮影した画像データに対し、病理医の診断を支援する種々の診断支援機能が提案されている。
従来、診断支援機能の一例として、以下の提案がされている。
非特許文献1は、がんを診断する上で重要な所見であるN/C比(細胞質に対して核が占める比率)を算出することを目標とし、デジタル画像処理技術を用いて肝臓の病理組織標本画像データから細胞膜を抽出する方法を開示している。非特許文献1では明視野、暗視野、位相差の3種類の観察像の色情報を組み合わせることで、明視野観察像単独の場合に比べて細胞膜の抽出正解率を向上させている。
また、細胞膜に限らず、細胞境界(細胞と細胞の間の細胞境界には細胞膜以外にも細胞間物質などが存在)や細胞と管や腔との境界を明瞭にすることは、診断を行う上で大きな意味がある。明瞭な境界は、医師が標本から複雑な肝臓の3次元構造を推測することを容易にするので、限られた情報からより精度の高い診断が実現できる。
また、細胞と管や腔との境界はN/C比を精度良く算出する上でも有用な情報である。例えば、肝臓の病理組織標本には、大別して核と細胞質からなる細胞の領域、肝細胞へ物質を供給する血管である類洞の領域があり、正しいN/C比を算出するには細胞が存在しない類洞の領域を正しく除外する必要がある。
特開2007−128009号公報
鳥澤奈美子,高橋正信,中野雅行,"肝病理組織標本画像中の細胞膜抽出におけるマルチイメージング利用の検討",電子情報通信学会総合大会,D−16−9,2009/3 児玉和也,久保田彰,"単一のレンズ系からの多様なボケ味の生成",映像情報メディア学会誌65(3),pp.372−381,2011年3月 児玉和也,久保田彰,"周波数領域上での線型結合に基づくScene Refocusing",映像メディア処理シンポジウム(IMPS2012),I−3.02,pp.45−46,2012年10月 Kazuya KODAMA, Akira KUBOTA, "Efficient Reconstruction of All−in−Focus Images Through Shifted Pinholes from Multi−Focus Images for Dense Light Field Synthesis and Rendering, IEEE Trans. Image Processing, Vol.22, Issue 11, 15pages (2013−11)
しかしながら上述した従来の技術においては、以下のような問題があった。
非特許文献1では、明視野、暗視野、位相差観察像を取得するために、明視野顕微鏡に対して位相差用対物レンズや共用コンデンサを装備し、それらを切り替えて撮影している。その為、明視野観察用の光学顕微鏡に追加の部品が必要になるというコスト的な課題、撮影時に光学系および露出条件の変更の手間が発生し、撮影時間が増加するという課題があった。また前記の手間は、広視野領域を局所的な視野毎に分割撮影して繋ぎ合せるバーチャル・スライド・システムにおいて新たな課題を生む。局所的な視野の撮影ごとにメカ機構を切り替えて広視野領域を撮影する場合、撮影時間の増加だけでなく、光学系を切り替えて撮影するためにその保持機構の耐久性の課題が生じる。一方、明視野、暗視野、位相差観察を変更することなく各々の広視野領域を撮影する場合には、各々の画像での焦点合わせ精度の差やステージ制御などを原因とする視野移動時の撮影位置の累積誤差が発生しやすい。そのため互いの画像データ間に位置ずれが生じやすく、同一の画素位置で各々の画像データの比較が難しくなるという課題が生じる。
本発明は、このような問題点に鑑みてなされたものであり、被写体を撮影して得られた元画像データから、画像処理によって、被写体の観察や診断に適した画像データを生成するための新規な技術を提供することを目的とする。
本発明の第1態様は、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを、周波数空間のデータに変換する変換手段と、前記変換された周波数空間のデータに対し、フィルタを適用するフィルタ適用手段と、前記フィルタが適用されたデータを実空間の画像データに逆変換する逆変換手段と、を有し、前記フィルタは、前記逆変換された画像データに含まれるレイヤー画像データのいずれかが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されていることを特徴とする画像データ生成装置である。
本発明の第2態様は、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを用いて、前記複数のレイヤー画像データのそれぞれを周波数空間のデータに変換する変換手段と、前記変換された周波数空間の複数のデータに対し、複数のフィルタをそれぞれ適用するフィルタ適用手段と、前記フィルタが適用された複数のデータを統合する統合手段と、前記統合されたデータを実空間の画像データに逆変換する逆変換手段と、を有し、前記複数のフィルタは、前記逆変換された画像データが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されていることを特徴とする画像データ生成装置である。
本発明の第3態様は、コンピュータが、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを、周波数空間のデータに変換するステップと、コンピュータが、前記変換された周波数空間のデータに対し、フィルタを適用するステップと、コンピュータが、前記フィルタが適用されたデータを実空間の画像データに逆変換するステップと、を有し、前記フィルタは、前記逆変換され
た画像データに含まれるレイヤー画像データのいずれかが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されていることを特徴とする画像データ生成方法である。
本発明の第4態様は、コンピュータが、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを用いて、前記複数のレイヤー画像データのそれぞれを周波数空間のデータに変換するステップと、コンピュータが、前記変換された周波数空間の複数のデータに対し、複数のフィルタをそれぞれ適用するステップと、コンピュータが、前記フィルタが適用された複数のデータを統合するステップと、コンピュータが、前記統合されたデータを実空間の画像データに逆変換するステップと、を有し、前記複数のフィルタは、前記逆変換された画像データが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されていることを特徴とする画像データ生成方法である。
本発明の第5態様は、本発明に係る画像データ生成方法の各ステップをコンピュータに実行させることを特徴とするプログラムである。
本発明によれば、被写体を撮影して得られた元画像データから、画像処理によって、被写体の観察や診断に適した画像データを生成することができる。
本発明の実施形態の画像データ生成および表示システムの構成図 画像表示アプリケーションの機能を説明する為の表示例 画像データ生成装置の内部構成を示す図 被写体の一例であるプレパラートを示す図 被写体を撮影する撮像装置の構成を模式的に示す図 視点画像データでコントラストが強調される理由を説明する図 視点の偏角および視線方向と光軸のなす角(観察角)の関係を示す図 プレパラート内の病理標本の表面に存在する凹凸を示す図 図8(a)の各面での観察角φにおける散乱光の強度を示す図 物体のZ位置が異なる場合のZスタック画像データを示す図 物体のZ位置が異なる場合の視点画像データと焦点ぼけを示す図 散乱光抽出画像データと散乱光強調画像データの焦点ぼけを示す図 実施例1のGUIの例を示す図 実施例1の散乱画像データ生成処理のフローチャート 実施例1の透過光成分抑制マスクとその生成処理を説明する図 実施例1の散乱光抽出画像データ生成処理のフローチャート 実施例1の透過光成分抑制処理を説明する図 実施例1のN/C比算出処理を示すフローチャート 実施例2の散乱画像データ計算処理を説明する図 実施例2の散乱光強調画像データ生成処理のフローチャート 実施例2の散乱光強調画像データ生成処理を説明する図 実施例3の散乱光強調画像データ生成処理を説明する図 散乱光抽出画像データ生成処理のフローチャート 散乱光抽出画像データ生成処理のフローチャート 視点重み関数および散乱光情報抽出用の視点重み関数の例 標本の表面凹凸を偏角が180度異なる視点から観察する様子を示す図 実施例6の視点の位置、偏角選択関数、偏角選択視点重み関数の例 実施例6の処理のフローチャート
(全体構成)
図1は、本発明の実施形態に係る画像データ生成および表示システムの構成を示している。
画像データ生成装置(ホストコンピュータ)100には、ユーザからの入力を受け付ける操作入力デバイス110と、画像データ生成装置100から出力される画像データなどをユーザに提示するためのディスプレイ120が接続される。操作入力デバイス110としては、キーボード111、マウス112、ユーザの操作性を高めるための専用コントローラ113(例えばトラックボール、タッチパッド)などを利用できる。また画像データ生成装置100には、ハードディスクや光学ドライブ、フラッシュメモリ等の記憶装置130、ネットワークI/Fを通じてアクセスできる他のコンピュータシステム140が接続されている。なお、図1では記憶装置130は画像データ生成装置100の外部に存在するが、画像データ生成装置100に内蔵しても良い。
画像データ生成装置100は、操作入力デバイス110から入力されたユーザの制御信号に従い、記憶装置130から画像データを取得し、画像処理を適用することによって観察に適した観察用画像データを生成したり、診断に必要な情報を抽出したりする。
画像表示アプリケーションおよび画像データ生成プログラム(いずれも不図示)は画像データ生成装置100で実行されるコンピュータプログラムである。これらのプログラムは画像データ生成装置100内の内部記憶装置(不図示)または記憶装置130に格納されている。後述する画像データ生成に関わる機能は画像データ生成プログラムによって提供されるものであり、画像データ生成プログラムの各機能は画像表示アプリケーションを介して呼び出す(利用する)ことができる。画像生成プログラムの処理結果(例えば生成された観察用画像データ)は、画像表示アプリケーションを介して、ユーザに提示される。
(表示画面)
図2は、予め撮影した病理標本の画像データを、画像表示アプリケーションを通じて、ディスプレイ120に表示した場合の一例である。
図2は画像表示アプリケーションの画面レイアウトの基本構成である。表示画面の全体ウィンドウ201内に、表示や操作のステータスと各種画像データの情報を示す情報エリア202、観察対象の病理標本のサムネイル画像203、病理標本画像の詳細観察用の表示領域205、表示領域205の表示倍率206、が配置されている。サムネイル画像203上に描画された枠線204は、詳細観察用の表示領域205に拡大表示している領域の位置および大きさを示している。このサムネイル画像203と枠線204によって、ユーザは病理標本画像全体中のどの部分を観察しているのかを容易に把握できる。
詳細観察用の表示領域205に表示する画像は、入力操作デバイス110による移動操作や拡大・縮小操作によって設定、更新できる。例えば、移動は画面上でのマウスのドラッグ操作により、拡大縮小はマウスホイールの回転等によって実現できる(例えば、ホイールの前方回転を拡大、後方回転を縮小に割り当てる)。また、焦点位置の異なる画像、すなわち奥行き位置の異なる画像への切り替えは、所定のキー(例えばCtrlキー)を押しながらのマウスホイールの回転等で実現できる。例えば、ホイールの前方回転を奥行きが深い画像への移動に、後方回転を奥行きが浅い画像への移動に割り当てる。上記のようなユーザの表示画像の変更操作に伴い、表示領域205、表示倍率206、サムネイル
画像203内の枠線204が更新される。このようにして、ユーザは所望する面内位置、奥行き位置、倍率の画像を観察できる。
(画像データ生成装置)
図3は画像データ生成装置100の内部構成を示す図である。
CPU301はメインメモリ302に格納されているプログラムやデータを用いて画像データ生成装置全体の制御を行う。またCPU301は、以降の実施例で説明する各種演算処理、データ処理、例えば、散乱画像データ生成処理や透過光成分抑制マスク生成処理、散乱光抽出画像データ生成処理、透過光成分抑制処理、散乱光強調画像データ生成処理等を行う。
メインメモリ302は記憶装置130からロードされたプログラムやデータ、他のコンピュータシステム140からネットワークI/F(インターフェース)304を介してダウンロードしたプログラムやデータを一時的に記憶するエリアを備える。またメインメモリ302は、CPU301が各種の処理を行うために必要とするワークエリアを備える。
操作入力デバイス110はキーボード102、マウス103、専用コントローラ113などCPU301に各種の指示を入力することのできるデバイスにより構成される。ユーザは画像データ生成装置100の動作を制御する情報を操作入力デバイス110により入力する。305は操作入力デバイス110を介して入力された各種の指示等をCPU301に通知するためのI/O(入出力)である。
記憶装置130はハードディスクなどの大容量情報記憶装置であり、OS(オペレーティングシステム)や以降の実施例で説明する処理をCPU301に実行させるためのプログラムや画像データなどを記憶する。記憶装置130への情報の書き込みや記憶装置130からの情報の読み出しはI/O306を介して行われる。
ディスプレイ制御装置307は画像データや文字等をディスプレイ120に表示させるための制御処理を行う。ディスプレイ120はユーザに入力を求めるための画面表示を行うとともに、記憶装置130や他のコンピュータシステム140から取得しCPU301で処理した画像データを表示する。
演算処理ボード303は、画像処理など特定の演算機能が強化されたプロセッサおよびバッファメモリ(不図示)を備えている。以降の説明では各種演算処理、データ処理にはCPU301を、メモリ領域としてメインメモリ302を用いるとして説明するが、演算処理ボード内のプロセッサやバッファメモリを用いることも可能であり、本発明の範疇とする。
(被写体)
図4は被写体の一例である病理標本のプレパラート(スライドとも呼ぶ)を表す。病理標本のプレパラートでは、スライドグラス410上に載置した病理標本400が封入剤(不図示)とその上に載せるカバーグラス411によって封入されている。病理標本400の大きさや厚みは病理標本によって異なっている。更にスライドグラス410上には病理標本に関する情報が記録されたラベルエリア412が存在する。ラベルエリア412への情報の記録は、ペンによる記入でもよいし、バーコードや2次元コードの印刷でもよい。また電気的、磁気的、または光学的な方法により情報を記憶可能な記憶媒体をラベルエリア412に設けてもよい。以降の実施形態では、被写体として図4に示す病理標本のプレパラートを例に説明する。
(撮像装置)
図5は被写体を撮影しデジタル画像を取得する撮像装置の構成の一部を模式的に表す。図5に示すように、本実施形態では、病理標本400の表面に平行にx軸とy軸をとり、病理標本400の奥行き方向(光学系の光軸方向)にz軸をとる。
プレパラート(病理標本400)をステージ502上に置き、照明ユニット501から光を照射する。病理標本400を透過した光は、撮像光学系503によって拡大され、撮像センサ504の受光面に結像する。撮像センサ504は複数の光電変換素子を有するラインセンサ(一次元センサ)またはエリアセンサ(二次元センサ)である。病理標本400の光像は撮像センサ504により電気信号に変換され、デジタルデータとして出力される。
一回の撮影で病理標本全体の画像データを取得できない場合には、ステージ502をx方向および/またはy方向に移動しながら複数回の分割撮影を行い、得られた複数の分割画像データを合成(繋ぎ合わせ)して病理標本全体の画像データを生成する。また、ステージ502をz方向に移動しつつ複数回の撮影を行うことで、光軸方向(奥行き方向)の焦点位置が異なる複数枚の画像データ(レイヤー画像データと呼ぶ)を取得する。本明細書では、光軸方向(奥行き方向)の焦点位置が異なる複数枚のレイヤー画像データからなる画像データ群を「Zスタック画像データ」と呼ぶ。また、被写体を撮影することによって取得された画像データであるレイヤー画像データやZスタック画像データを「元画像データ」と呼ぶ。
図2の表示倍率206に表示される倍率の値は、撮像光学系503の倍率に画像表示アプリケーション上での拡大/縮小率を掛けた値である。なお、撮像光学系503の倍率は固定でも良いし、対物レンズの交換によって可変であっても良い。
(視点画像データを生成する技術の説明)
画像データ生成装置100では、暗視野観察や位相差観察など光学系に変更を加える観察・撮像方法を必要としない。その代わりに、Zスタック画像データから画像処理によって中間画像データ(視点画像データ)を生成し、その中間画像データを用いて観察や診断に適した観察用画像データを生成する。まずは、Zスタック画像データから中間画像データとしての視点画像データを生成する処理に利用可能な技術について説明する。
光軸方向の焦点位置を変えて撮影した複数枚の画像データ(Zスタック画像データ)を元に、任意の方向から観察した視点画像データ(任意視点画像データ)を生成できることが知られている。ここで、視点画像データとは所定の観察方向(即ち視点)から被写体を観察した画像データを表す。
例えば特開2007−128009号公報(以降、特許文献1と呼ぶ)には、焦点位置を変えて撮影した焦点ぼけ画像データ群から、任意の視点や任意のぼけの画像データを生成する方法が開示されている。この方法では、焦点ぼけ画像データ群に対し、焦点位置の前後で発生する2次元の焦点ぼけを集めた3次元的な焦点ぼけ(以降、3次元焦点ぼけと呼ぶ)がXYZ位置で不変になるように座標変換処理を施す。そして、得られた直交座標系(XYZ)において3次元のフィルタ処理を適用することで視点やぼけを変更した画像データを得る。
また、非特許文献2には、特許文献1の方法の改良が開示されている。非特許文献2では、視点から視線方向を求め、Zスタック画像データを視線方向に積算することで積算画像データを生成するとともに、同じように3次元焦点ぼけの視線方向の積算画像データも生成する。その後、Zスタック画像データの積算画像データに対し3次元焦点ぼけの積算画像データを逆フィルタ処理することにより、Z方向(レイヤー画像データの枚数)の制約による影響を抑制し、高画質な視点画像データを生成することができる。
また、非特許文献3には、非特許文献2の計算を高速化する方法が開示されている。非
特許文献3の方法では、被写体(シーン)に依存せずに予め定まるフィルタと焦点ぼけ画像データ群の各Z位置のフーリエ変換画像データの線型結合によって、周波数領域上での任意視点画像データや任意焦点ぼけ画像データが効率的に計算できる。
非特許文献4は、非特許文献3の方法をより詳細に記した文献である。
本実施形態では、光軸方向の焦点位置を変えて撮影した複数枚の画像データ(Zスタック画像データ)を元に、任意の方向から観察した視点画像データ(任意視点画像データ)を生成したり、任意の焦点ぼけを持つ画像データを生成する。以降の説明では、これらの手法を総称して、MFI(マルチフォーカスイメージング)任意視点/焦点ぼけ画像生成法と呼ぶこととする。
なお、両側テレセントリックな光学系を持つ顕微鏡で焦点位置を変えて撮影したZスタック画像データは3次元焦点ぼけがxyz位置で不変である。よって、両側テレセントリックな光学系で撮影したZスタック画像データにMFI任意視点/焦点ぼけ画像生成法を適用する場合は、座標変換処理および座標変換処理に伴う画像データの拡縮処理は不要である。
ライトフィールドと呼ばれる4次元の情報(XYの2次元画像データに視点位置の自由度を加えた情報)が記録された画像データを1回の撮影で取得可能な撮像装置が知られている。このような撮像装置はライトフィールドカメラやライトフィールド顕微鏡と呼ばれる。これらの装置では、本来結像面となる位置にレンズアレイが配置されており、それよりも後方のイメージセンサでライトフィールドを撮影する。ライトフィールドが記録された元画像データからも、公知の技術を用いて、任意の焦点位置の画像データや任意の方向から観察した視点画像データ(任意視点画像データ)を生成できる。
本実施例では、撮像装置の被写体に対する方向を物理的に変えることなく、Zスタック画像データ又はライトフィールドなどの撮像画像データを元に、デジタル画像処理によって生成される任意の観察方向の画像データを「視点画像データ」とよぶ。この視点画像データは、被写体の撮影に用いる撮像光学系を通る任意の光線を主光線とし、その主光線を中心とする光束によって撮像面に形成される光学像を模擬した画像データである。主光線の方向が観察方向に対応する。主光線の方向は任意に設定できる。また光束の大きさ(NA)も任意に設定できる。画像診断等が目的の場合には、視点画像データの被写界深度は深いことが望ましいので、視点画像データに対応する光束のNAは0.1以下が望ましい。
なお、デジタル画像処理によって生成(計算)した視点画像データは、撮像光学系の撮影条件(絞りの位置・大きさ)や光軸方向やレンズ等を物理的に変えて撮影した画像データとは、必ずしも一致しない。しかし、現実に撮影した画像データと一致していなくても、視点を変えて被写体を観察したのと同様の特徴をもつ画像データであれば(つまり、観察方向を変えるのと同様の効果をデジタル画像処理によって付与できれば)、画像観察や画像診断等には有用である。従って、現実に光軸方向等を変えて撮影した画像データとは厳密には一致しないが、そのような画像データと同様の特徴が現れるようにデジタル画像処理された画像データも、本実施例の視点画像データに含まれる。
座標変換を施した焦点ぼけ画像群から、実空間上のレンズ面(瞳面に相当)上にある原点O(x,y,z)=(0,0,0)から視点(x,y,z)=(s,t,0)だけズレた位置のピンホールを通して観察した視点画像データを生成できる。MFI任意視点/焦点ぼけ画像生成法では、レンズ面上の視点の位置を変えることで被写体を観察する観察方向、すなわち視線方向を変化させることができる。
視線方向は、結像した像に対応する被写体の所定の位置から発する光束の中で、レンズ
面上の視点の位置(x,y,z)=(s,t,0)を通る直線の傾きである、と定義できる。視線方向は様々な方法で表現可能である。例えば、直線の進行方向を示す3次元的なベクトルによる表現でも良いし、前述の3次元的なベクトルが光軸となす角(観察角)と光軸に垂直な平面に射影したときのベクトルがX軸となす角(偏角)による表現でも良い。
撮像光学系が両側テレセントリックではない場合、撮像面における3次元焦点ぼけは、合焦した被写体の空間的な位置(xyz座標内の位置)によって変化し、レンズ面上の視点の位置(x,y,z)=(s,t,0)を通る直線の傾きは一定とはならない。その場合、特許文献1に記載の座標変換後の直交座標系(XYZ)の上で視線方向を定義すると良く、視線方向は(X,Y,Z)=(−s,−t,1)のベクトルで表すことができる。以下、座標変換後の視線方向の求め方を説明する。
特許文献1には、撮像光学系の焦点が合った任意の位置と撮像装置のレンズ面(瞳面に相当)上の同一の視点の位置(x,y,z)=(s,t,0)を結ぶ光線は、座標変換後の直交座標系(XYZ)では全て互いに平行な光線となることが記載されている。(特許文献1の図1〜3およびその説明を参照)
透視座標系(座標変換前の実空間)の被写体のある点から出た光は、(p+s,q+t,f)(fは焦点距離)を通過し、視点の位置(x,y,z)=(s,t,0)で屈折する。この直線は以下の式で表される。
Figure 2015192238
数1の直線は座標変換後の直交座標系(XYZ)では以下の式で表される。
Figure 2015192238
数2にZ=0(z=f)、Z=1(z=∞)を代入し、3次元座標を求めると、それぞれ(X,Y,Z)=(p+s,q+t,0)、(X,Y,Z)=(p,q,1)となることから、直交座標系(X,Y,Z)での直線の傾きは(−s,−t,1)で表される。
従って、座標変換後の直交座標系での視線方向を表すベクトルは(X,Y,Z)=(−s,−t,1)となる。
なお、撮像光学系が両側テレセントリックである場合、奥行き方向に焦点を変えて撮影した複数枚の画像データ(Zスタック画像データ)中の3次元焦点ぼけはZの位置によらず不変となる。
従って、空間的な位置によらず3次元焦点ぼけを不変とするための座標変換は必要ない。実空間においてピントが合った被写体の所定の位置(x,y,z)=(0,0,za)とレンズ面上の視点の位置(x,y,z)=(s,t,0)を結ぶ直線の傾き(−s,−t,za)をそのまま視線方向と見なしても良い。
(視点と実際に標本を観察したときの偏角θ及び観察角φの対応関係)
図7(a)は実空間上の視点の位置(x,y,z)=(s,t,0)を表す模式図であり、図7(b)は直交座標系(XYZ)において視点の位置(x,y,z)=(s,t,0)を通る光線を表す模式図である。
図7(a)で示す点線の円はレンズ面上(z=0)で光線が通過可能な範囲を表している。偏角θをレンズ面上の視点の位置(x,y,z)=(s,t,0)がレンズ面上(z=0)のx軸となす角、あるいは視線(−s,−t,1)をxy平面に射影したときの直
線がx軸となす角と定義すると、偏角θは下記の式で求められる。
Figure 2015192238
ただし、θはt,sの符号に応じて−180〜+180度の範囲に収まるように調整する。
続いて、図7(b)を用いて視線と変換座標上の観察角φの関係について説明する。
図7(b)では、数2で示す直線と数2に光軸上の点p=0、q=0を代入した場合の直線を太字の矢印で示している。
特許文献1によれば、直交座標系(XYZ)のZ=0は透視座標系(xyz)のz=f(またはz=−∞)に対応し、Z=1はz=∞(またはz=−f)に対応している。そのため、図7(b)は直交座標系XYZにおいて無限遠(Z=1)からの光束が、レンズ面上の手前の焦点面(Z=0)で広がりを持つことを示している。(特許文献1の図3およびその説明を参照)
ここで、変換座標上の観察角φを、視線(−s,−t,1)と光軸(Z軸)がなす角と定義すると、図7(b)からも明らかなように視線は、被写体の位置に依存しないので、観察角φは以下の式で求められる。
Figure 2015192238
なお、図7(b)の2本の点線はレンズ面上の最も端を通る光線を示している。座標変換前の透視座標系(xyz)でのレンズの絞り半径をrとすると、視点の位置(x,y,z)=(s,t,0)が半径rの内部にある場合のみ視点画像データは計算できる。
次に、変換座標上ではなく、実際に標本を観察したときの視線に対応する偏角θ、観察角φについて述べる。スネルの法則では屈折率の異なる境界に光線が入射したとき、光線の入射角と入射側の媒質の屈折率の積は、光線の屈折角と屈折側の媒体の屈折率の積に等しい。標本の屈折率は空気の屈折率よりも大きいことから、空気中の観察角に比べ標本中の観察角は小さくなっている。そのため、屈折した光線で構成される、標本中の3次元焦点ぼけは、空気中での3次元焦点ぼけよりも小さくなっている。しかし、本実施例では標本中での3次元焦点ぼけによる3次元的な結像関係に基づいて視点の位置を計算しているため、標本の屈折率の影響は考える必要はなく、偏角θおよび観察角φはそのまま標本中での観察方向を表している。
撮像光学系が両側テレセントリックである場合には、座標変換を必要としないため、x方向とy方向のセンサ画素ピッチが等しいとすれば、観察角φは、x方向のセンサ画素ピッチΔxと、z方向の移動間隔Δz(単位はμm)を用いて、下記式で表現できる。
Figure 2015192238
なお、撮像光学系が両側テレセントリックでない場合には、数5でΔxとΔzの代わりに、直交座標系(XYZ)におけるX方向のセンサ画素ピッチΔXとZ方向の移動間隔ΔZを用いれば観察角φが求められる。
以上で、実際に標本を観察したときの視線に対応する偏角θ、観察角φについて説明した。
以降の説明では、レンズ面上の視点の位置(x,y,z)=(s,t,0)を視点(s
,t)と略して記載する。また、直交座標系(XYZ)上での画像処理を前提として説明するため、視点の位置(s,t)に言及する場合のみ透視座標系(座標変換前の実空間)での位置を表すものとし、その他は特に断りが無い限りは、直交座標系(XYZ)での位置を表すとする。
図5の撮像装置で取得したZスタック画像データに対しMFI任意視点/焦点ぼけ画像生成法を適用すると、視点の位置、即ち観察方向を変えた視点画像データを生成することが出来る。
MFI任意視点/焦点ぼけ画像生成法で計算する視点画像データには、主に2つの特徴がある。1つは視点画像データが非常に深い(無限大の)被写界深度を持ち、透過率が異なる標本内の物質の境界が明瞭に見える点である。もう1つは、視点画像データは、照明の一部領域から標本を照らして得る偏斜照明の観察像に近く、XY平面での視線方向に沿って変化する凹凸が強調され、標本が立体的に見える点である。視点画像データでは偏斜照明の像と同様、光軸に対して視線方向が傾くほど、即ち、視線の観察角φが大きくなるほど、標本表面の凹凸のコントラストが高くなり、標本表面が立体的に見える。
(ただし、物理的には偏斜照明の像と視点画像データは異なっている。偏斜照明の像は焦点位置の変更に伴って光学的なぼけが発生するが、視点画像データは焦点位置の変更によらず被写界深度は非常に深いままであるという違いがある。なお、視点画像データはピントを合わせるZスタック画像データのZ位置Zfによって変化するが、その変化はXY方向の平行移動で表される。)
次に、1つ目の特徴である、透過率が異なる物質の境界が明瞭に見える理由について説明する。
図6(a)は、直交座標系(XYZ)での光学系の3次元焦点ぼけを表す図である。600は3次元焦点ぼけ形状を示し、焦点位置(2つの円錐の頂点)では焦点ぼけは僅かだが、Z位置が焦点位置から離れるに従い、焦点ぼけが広がる様子を示している。MFI任意視点/焦点ぼけ画像生成法を用いれば、Zスタック画像データから円錐600の内部を通る任意の視線方向(例えば直線610)の光線で構成される視点画像データを生成することが出来る。
図6(b)は、直交座標系(XYZ)での病理標本を異なる方向から見た様子を示す。図6(b)の標本620の内部には斜め方向の空洞630が存在している。
方向631から観察すると空洞630以外の部位が透けて見える為、空洞630の壁面のコントラストは明瞭ではない。また方向632から観察した場合も同様で空洞630のコントラストは不明瞭である。しかし空洞630の壁面に沿った方向633から観察する場合、他の部位の影響を受けないので、空洞630の壁面のコントラストは明瞭となる。なお、視線方向が多少空洞の壁面の方向と異なってもコントラストは比較的高い状態を維持できる。
一方、標本620のZスタック画像データの3次元焦点ぼけは様々な視線方向の光線を統合したものである。そのため、いずれのZ位置(焦点位置)のレイヤー画像データにおいても、方向631〜633の光線を含む多方向の光束の影響で像がぼやけ、空洞の壁面のコントラストは方向633からの観察像に比べて明瞭にはならない。この現象は空洞に限らず、核や細胞膜、線維等においても同様である。
以上、視点画像データにおいて標本内の透過率が異なる物質の境界が明瞭に見える現象について説明した。
病理標本は半透明の物体であり、透過光以外にも散乱光が存在している。その散乱光の存在が視点画像データの2つ目の特徴を生んでいる。
次に、標本での散乱光により、視線の観察角φが大きくなるほど、標本表面の凹凸のコ
ントラストが高くなる理由について説明する。
図8(a)の800はプレパラート内の病理標本の表面に存在する凹凸を示す模式図である。図8(a)に示すxz平面の凹凸は奥行き方向であるy方向にも続いているとする。
組織診用の病理標本はパラフィンで固定された後、マイクロトームで均一な厚みにスライスされ、その後、染色が施されている。しかし、病理標本は完全に均一ではなく、細胞と管や腔の境界、核と細胞質の間の境界などでは組織の構造や物質の成分に起因する凹凸が存在し、病理標本の表面には図8(a)に示すような起伏のある構造が存在している。
図8(a)は簡易的なモデルであり、実際の標本の凹凸には図8(a)のようにとがった部分は少ない。また図8(a)のような凸の構造だけでなく標本内部に凹んだ構造も存在する。また表面が平らでも内部に屈折率が異なる物質が存在する場合は光学的な距離が変わるため、標本内部の屈折率の不連続は表面凹凸と見なすことができる。
なお、実際のプレパラートでは、カバーグラスと標本の間に透明な封入剤が存在している。しかし、封入剤の屈折率と標本の屈折率の差は僅かで影響は少ないため、以降、両者の屈折率は同一として説明する。
図8(a)の811は凹凸の無い面、812は右上がりの斜面、813は右下がりの斜面を示す。斜面812、斜面813がx軸となす傾斜角はそれぞれα(α>0)である。
図9(a)〜図9(c)は図8(a)の811〜813の面での観察角φにおける散乱光の強度を示す模式図である。図9(a)、図9(b)、図9(c)はそれぞれ平面811、斜面812、813での光の散乱を示している。それぞれの面に接する円は、標本表面での光の拡散特性を完全拡散透過面と仮定した場合の散乱方向による散乱光の強度を示す。円内の太い矢印の線の長さは光軸(z軸)からφだけ傾けた角度から観察した場合の散乱光の強度を示している。(実際には標本表面は完全拡散透過面ではなく、光の入射方向・観察方向による強度依存性があるが、説明を簡略化するため、ここでは完全拡散透過面と仮定して説明する。)
完全拡散透過面では、面と直交する法線方向の光の強度をI、観察方向と面の法線のなす角をδとすると、δ方向の散乱光の強度I(δ)はI(δ)=Icosδで表せる。
図9(a)、図9(b)、図9(c)において、観察方向と面の法線のなす角δはそれぞれφ、φ+α、φ−αで表せるため、それぞれの散乱光の強度は、
cosφ、Icos(φ+α)、Icos(φ−α)
となる。
なお、観察方向から見てzの値が増加する斜面(上りの斜面)では傾斜角αを正とし、zの値が減少する斜面(下りの斜面)では傾斜角αを負とすれば、いずれの面でも散乱光の強度は、Icos(φ−α)で表せる。
斜面812、813の観察角φ方向の散乱光の強度を平面811の観察角φ方向の散乱光の強度で割った値をコントラストC(φ,α)として定義すると、コントラストは下記の式となる。
Figure 2015192238
表1にφおよびαを変えたコントラストC(φ,α)の値を示す。
Figure 2015192238
表1より、観察角φが小さいときは、傾斜角αが大きくても斜面812、813の間のコントラストは低いため観察しづらく、観察角φが大きくなるに従い、傾斜角αが小さくてもコントラストは大きくなり、観察しやすくなることが分かる。
以上、標本での散乱光により、視線の観察角φが大きくなるほど、標本表面の凹凸のコントラストが高くなる理由について説明した。
(視点を変えたときの散乱光強度の変化)
次に視点を変えたときの標本表面の散乱光強度の変化について説明する。
図8(a)は、xy面内において表面凹凸のエッジ方向がx軸に直交しており、表面凹凸の高低変化の方向(エッジに垂直な方向)が視点の偏角θ(θ=0)と一致した場合の図である。一方、図8(b)は表面凹凸の高低変化の方向がx軸と一致しない場合の図である。図8(b)に示すように表面凹凸(820)のエッジに垂直な方向(821)とx軸とのなす角を凹凸方向角βとすると、凹凸方向角βと偏角θが一致しない場合には、表面凹凸820を斜め方向から観察することになる。このとき、偏角θ−βの角度を持つ観察方向から見た、斜面813とその反対方向の斜面812の見かけ上の傾斜角α’は
Figure 2015192238
で求まる。数7から、見かけ上の傾斜角α’はαより小さくなり、凹凸方向角βと偏角θの差|θ−β|によって、コントラストCが低下することが分かる。
ここで、|θ−β|が90度を境に、傾斜角αおよびα’の符号は正から負、あるいは負から正に切り替わる点に注意する。これは、視点の偏角を変えた観察方向によって上りの斜面と下りの斜面が変わることに対応する。傾斜角α’は−α〜+αの範囲を取り、例えば、斜面813では視点の偏角が|θ−β|=0のときα’=α(観察方向から見て上りの斜面)、|θ−β|=πのときα’=−α(下りの斜面)となる。
次に、視点の偏角がθ−βの方向から観察した斜面813の散乱光強度を平面811で観察される光の強度で正規化した、散乱光正規化強度V(φ,α)について考える。平面
811で観察される光の強度をIts(φ)で表し、透過光の強度および散乱光の強度の和で構成される観察角φの関数とする。前記より、散乱光正規化強度V(φ,α)は以下の式で表せる。
Figure 2015192238
病理標本では表面凹凸の傾斜角αは十分に小さいと推定できるため、cosα=1、sinα=αの近似式を用いて、数8は以下のように近似できる。
Figure 2015192238
続いて、A(φ)およびB(φ)のとる値について述べる。
光の強度Its(φ)は照明光学系の特性から一般に観察角φに従って減少する傾向があるため、Iを観察角φ=0のときの透過光強度とし、Its(φ)=Icosφとして近似した場合、A(φ)=I/Iとなる。散乱光が透過光に比べて小さい場合、A(φ)は比較的小さな定数とみなすことができる。一方、B(φ)は、φ=0のとき、sinφ=0よりB(φ)=0となる。観察角φの増加に従ってsinφは増加し、また光の強度Its(φ)は減少すると見なせることから、B(φ)は増加関数となる。Its(φ)=Icosφと仮定すれば、B(φ)=I/I×tanφとなる。
MFI任意視点/焦点ぼけ画像生成法で計算する視点画像データでは、平均を維持する周波数フィルタ処理により、視点にかかわらず画像データの平均輝度は変化しない。従って、視点画像データでの斜面813の輝度から染色部位の透過率の影響を受けた透過光強度を除いた値は、正規化強度V(φ,α)にほぼ等しいと考えることができる。
(視点を変えたときの透過光強度の変化)
次に、視点が異なる視点画像データでの透過光強度の変化について考える。
観察する標本において、隣接する物質との輝度差が少ない領域(細胞質や細胞境界)では視線方向が違っても輝度差は少なく、透過光の強度には影響が少ない。また、組織診病理標本のように厚みが4μm程度の比較的薄い標本では、標本表面にピントを合わせて計算した視点画像データでは観察される標本内部の物体の位置は大きくは変わらない。(MFI任意視点/焦点ぼけ画像生成法で計算される視点画像データではピントを合わせたZスタック画像データのZ位置Zfの近傍に存在する物体は、視点によらずほぼ同一の視点画像データのXY位置に現れる。それは複数視点の視点画像データを統合した画像データにおいてXY位置に差がないため、ぼけが少なくなることからも理解できる。)
従って、視点の位置を変えて求めた視点画像データ間での透過光強度の差は僅かと見なすことができ、散乱光正規化強度間の差分と視点画像データ間の差分はほぼ等しいと考えることができる。
(視点の偏角θを変えた視点画像データからの標本表面の散乱光の情報の抽出)
続いて、偏角θを変えた視点画像データ間の演算により、標本表面の散乱光の情報を抽出することを考える。
図8の813に示す傾斜角αの所定の表面凹凸において、数7より、凹凸方向角βと偏角θの差|θ−β|が0度の場合、傾斜角はαであり、散乱光正規化強度V(φ,α)は最大となる。一方、|θ−β|が180度の場合、見かけ上の傾斜角α’は−αとなり、
散乱光正規化強度V(φ,α’)は最小となる。
従って、|θ−β|が0になる視点画像データと|θ−β|が180度となる視点画像データの間で減算をすれば、効率的に標本表面の散乱光の情報(の一部)が抽出できることが分かる。以下に式を示す。
Figure 2015192238
(ただし、|θ−β|=πのときα’=−α)
標本の表面凹凸の傾斜角α、凹凸方向角βは様々であるため、偏角θを変えた様々な視点において偏角θが180度異なる視点の視点画像データ間で差分を求め、その結果を集めれば、標本表面の散乱光の情報が抽出できることが分かる。
なお、数10より、偏角θが180度異なる視点の視点画像データ間以外の演算であっても、散乱光の情報の抽出は可能であることが分かる。例えば、|θ−β|=π/2(90度)のときα’=0であり、視点画像データ間の差分からα×B(φ)が抽出できる。
また、減算だけでなく除算でも標本表面の散乱光の情報が抽出できる。すなわち、V(φ,α)/V(φ,α’)の値も傾斜角αに応じて変化するので、散乱光の情報を表す指標といえる。
(視点の観察角φを変えた視点画像データからの標本表面の散乱光の情報の抽出)
続いて、前述の偏角θの場合と同様に、観察角φを変えた視点画像データ間の演算により、標本表面の散乱光の情報を抽出することを考える。
観察角がφとφ’の2つの視点画像データ間の減算は下記の式で表現できる。(ただし、A(φ)=A(φ’)と近似)
Figure 2015192238
既に述べたように、B(φ)はφの増加に伴って増加する関数であることから、観察角φを変えた2枚の視点画像データ間の演算により、標本表面の散乱光の情報(の一部)が抽出できることが分かる。B(φ)はφ=0のとき0となるため、観察角φの視点画像データと観察角φ’が0になる視点の視点画像データとの減算が最も効果的である。
前述の偏角θの場合と同様、標本中の凹凸方向角βは様々であるため、観察角φだけでなく、偏角θを変えた様々な視点において視点画像データ間で差分を求め、その結果を集めれば、標本表面の散乱光の情報が抽出できることが分かる。
また、減算だけでなく除算でも標本表面の散乱光の情報が抽出できる。これまでの説明から、観察角φにおける視点画像データの輝度は、α×B(φ)+D(Dは定数)で近似できるため、観察角φを変えた視点画像データ間の輝度の除算は、
Figure 2015192238
で表現できる。
DIV(φ,φ’)は、αが0のとき1となり、αが0でないときはαの大きさに応じて強度が変わる。よって、除算DIV(φ,φ’)によっても、傾斜角αの情報が抽出できることになる。
減算の場合と同様に、観察角φだけでなく、偏角θを変えた様々な視点において視点画像データ間で除算をし、その結果を集めれば、標本表面の散乱光の情報が抽出できる。
本実施例のように透過型の顕微鏡により撮影した画像データ(顕微鏡画像データ)には、標本を透過した光による透過光成分と、標本の表面等で散乱した光による散乱光成分とが含まれる。透過光成分の強度は標本の光の透過率に依存するので、標本内の色の違いや屈折率の違いなどを表している。一方、散乱光成分の強度は、前述のように、主に標本表面の凹凸(表面形状)に依存する。元画像データとしての顕微鏡画像データでは透過光成分が支配的なため散乱光成分は殆ど認識できないが、上記のように観察方向の異なる複数の視点画像データの間の差異を抽出ないし強調することで、散乱光成分を抽出ないし強調できる。言い換えれば、観察方向の異なる2つの視点画像データの「差(差分)」および「比」は、元画像データに含まれる散乱光成分(散乱光の情報)を表している特徴量といえる。
なお、「強調」とは、ある部分を(元の状態よりも)目立たせる操作であり、「抽出」とは、ある部分のみを取り出す操作であるが、ある部分に注目するという点では同じ操作であるため、本明細書では二つの用語を特に区別せずに用いる場合もある。また、散乱光成分の強調という操作には、画像データ中の散乱光成分の強度を上げる操作だけでなく、画像データ中の透過光成分の強度を下げることで散乱光成分の強度を相対的に上げる操作も該当する。以降、画像データに含まれる散乱光成分(散乱光の情報)を抽出ないし強調して得られる画像データを散乱画像データ又は特徴画像データとよぶ。また、元画像データから散乱画像データを生成する操作を散乱画像データ生成(特徴画像データ生成)とよぶ。
また本発明では、散乱画像データのうち、被写体から散乱光を抽出した画像データを散乱光抽出画像データ、被写体の散乱光を強調した画像データを散乱光強調画像データとよぶこととする。
なお、視点画像データの差異を抽出ないし強調する演算として減算と除算を例示したが、画像データの差異を強調できればどのような演算を用いてもよい。
散乱画像データは、標本表面の凹凸(表面形状)の観察及びそれに基づく診断に適した観察用画像データとして利用価値があるものと期待される。例えば、標本内のある領域の表面に凹凸があったとしても、その領域内の透過率がほぼ均一な場合には、元画像データではその領域の輝度又は色は一様になるため、目視では凹凸を認識できない。散乱画像データはこのような表面凹凸(輝度や色の変化として現れない表面凹凸)を可視化するのに特に有効である。また、細胞境界や細胞と類洞の境目なども凹凸があるので、散乱画像データはこれらの境界を明瞭化する効果もある。なお、画像特徴を抽出ないし強調する処理の一つにエッジ抽出(強調)処理があるが、画像の輝度や色の変化として現れない表面凹凸については、エッジ抽出(強調)では可視化できない。よって、標本のどのような構造ないし特徴を観察したいかに応じて、散乱画像データ抽出とエッジ抽出を使い分けるとよい。
なお、これまで主に散乱光が標本表面で発生する場合について説明したが、観察対象は必ずしも標本表面に限らない。フォーカスを合わせた位置(Z=Zf)に散乱の原因となる、周囲と屈折率が異なる物質などが存在する場合でも標本表面の散乱現象と同様に考えることができ、散乱画像データを用いて観察できる。
(MFI任意視点/焦点ぼけ画像生成法における3次元的な結像関係)
ここまで被写体の表面(フォーカスを合わせた面)での散乱光が持つ特徴について述べた。しかし、被写体には散乱光だけでなく、透過光も存在する。視点画像データにおける透過光の影響について述べるにあたり、まず、MFI任意視点/焦点ぼけ画像生成法における3次元的な結像関係について説明する。
MFI任意視点/焦点ぼけ画像生成法では、座標変換を施した焦点ぼけ画像群(Zスタック画像データ)が3次元的な被写体と3次元焦点ぼけのコンボリューションによって表
される関係が成立している。3次元被写体をf(X,Y,Z)、3次元焦点ぼけをh(X,Y,Z)、座標変換を施した焦点ぼけ画像群(Zスタック画像データ)をg(X,Y,Z)とすると以下の式が成立する。
Figure 2015192238
(ただし、***は3次元コンボリューションを表す。)
図10(a)〜図10(f)は3次元被写体の異なるZ位置に物体が存在する場合のZスタック画像データの違いを示す模式図である。3次元被写体1001はZ=Zfの位置に物体が存在し、3次元被写体1011ではZ=Zoの位置に物体が存在する。
3次元焦点ぼけ1002を持つ光学系を用いてZ方向の焦点位置を変えながら、3次元被写体1001、1011を撮影することで、それぞれZスタック画像データ1003、1013が得られる。(なお、3次元焦点ぼけh(X,Y,Z)は座標変換によってシフト不変に変換されるため、同一である。)自明だが、Zスタック画像データ1003ではZ=Zfの位置で、Zスタック画像データ1013ではZ=Zoの位置で最も焦点ぼけが少ない画像データが得られる。
(MFI任意視点/焦点ぼけ画像生成法での任意視点画像データの計算)
MFI任意視点/焦点ぼけ画像生成法として特許文献1、非特許文献2,3,4の方法が知られている。一例として非特許文献2の方法による任意視点画像データ生成方法を説明する。
視点(s,t)に対応する視線方向から観察した視点画像データをas,t(X,Y,Zf)とする。as,t(X,Y,Zf)は、同一の視線方向のZスタック画像データの積算画像データb(X,Y,Zf)を、同一の視線方向の撮像光学系の3次元焦点ぼけの積分値cs,t(X,Y,Zf)でデコンボリューションして得られる。(図10で説明される関係、即ち、f(X,Y,Z)とh(X,Y,Z)をコンボリューションすることでg(X,Y,Z)が得られる関係は視線方向への積分にかかわらず成立する。)
式で表現すると以下の関係が成立する。
Figure 2015192238
ただし、As,t(u,v),Bs,t(u,v),Cs,t(u,v)はそれぞれas,t(X,Y,Zf)、bs,t(X,Y,Zf)、cs,t(X,Y,Zf)のフーリエ変換である。u,vはそれぞれX,Y方向の変化に対応する周波数座標である。
s,t(X,Y,Zf)は以下の式で求められる。
Figure 2015192238
ただし、F−1{}はフーリエ逆変換を表す。
次に、MFI任意視点/焦点ぼけ画像生成法における任意視点画像データと任意焦点ぼけ画像データの関係について説明する。
図11(a)〜図11(h)は、3次元被写体の異なるZ位置に物体が存在する場合のMFI任意視点/焦点ぼけ画像生成法で求める任意視点画像データと任意焦点ぼけ画像データの関係を示す模式図である。
図11(a)の1100、図11(e)の1110は3次元被写体であり、それぞれ図10(a)の1001、図10(d)の1011と等しい。3次元被写体1100および1110中の実線で表す光線は、被写体中の物体を通り、かつ光学系のレンズ面上のそれ
ぞれに対応した所定の視点を通る光線を表す。なお、図11(e)中の破線は図11(a)中の実線に対応している。
光線1100aおよび1110aはレンズ面上のある視点aを通る光線であり、視点に対応する視線方向からZ=Zfの位置が焦点位置となるように観察した視点画像データはそれぞれ、図11(b)の1101、図11(f)の1111となる。また同様に光線1100bおよび1110bはレンズ面上のある視点bを通る光線であり、それぞれのZ=Zfの位置が焦点位置となるように観察した視点画像データは図11(c)の1102、図11(g)の1112となる。
視点画像データ1101、1102では物体の像は同一の位置に現れる。これに対し、視点画像データ1111、1112では物体の像の位置は同一の位置に現れず、それぞれ、ZfとZoの間の距離dZ(=Zo−Zf)と視線方向で決まる量だけシフトして現れる。
MFI任意視点/焦点ぼけ画像生成法では、Zスタック画像データのZ=Zfの位置のレイヤー画像データg(X,Y,Zf)は以下の式で表現できる。
Figure 2015192238
ここで、k(s,t)は撮像光学系の3次元焦点ぼけを間接的に表す関数で、レンズ面上の各視点(s,t)を通る光線の相対強度分布を表している。k(s,t)とh(X,Y,Z)の間には以下の関係が成立する。
Figure 2015192238
k(s,t)は以下の式に示すように、レンズ面上に存在する全ての視点(s,t)での総和が1になるように規格化されているものとする。
Figure 2015192238
また、数16のas,t(X,Y,Zf)はZ=Zfにおける視点画像データを表す。
数16から、あるZ位置(Z=Zf)における複数の視点画像データを重み付け合成することによって、そのZ位置(Z=Zf)における任意の焦点ぼけをもつ画像データ、すなわち、任意焦点ぼけ画像データを再構成できることが分かる。
以降の説明では、k(s,t)のように、複数の画像データを統合する際の視点(s,t)ごと(視線方向ごと)の重みを定義する関数を「視点重み関数」とよぶ。視点重み関数は、被写体上のある点から発し視点(s,t)を通る光線の相対強度分布を表しているともいえる。上記のk(s,t)は、被写体の撮影に用いた撮像光学系の3次元焦点ぼけに対応する特性をもつ視点重み関数である。なお、視点重み関数としては、どのような関数を用いてもよい。例えば、MFI任意視点/焦点ぼけ画像生成法で任意焦点ぼけ画像データを生成するための、任意の3次元焦点ぼけに対応するk(s,t)も視点重み関数として用いることができる。また、後述の実施例で例示する円柱形のような任意の形状の関数を視点重み関数として用いることもできる。
なお、数16の式は連続的な値の積分を意味するインテグラル∫で表記されているが、実際の画像処理では離散的な(有限個の)複数の視点(s,t)に対する演算であるためシグマΣで表記するのが正しい。しかし、インテグラル∫の方が一般化できるため、以降の説明でもインテグラル∫を用いて説明する。なお、視点が離散的(有限個)である場合では、数18は演算で用いる全ての視点重み関数k(s,t)あるいはk(s,t)の
和が1となるように規格化することを表す。
(Z=Zfにおけるぼけ像)
次に、MFI任意視点/焦点ぼけ画像生成法を用いて、Z=Zfにおける視点画像データからZスタック画像データのZ=Zfの位置のレイヤー画像データを求める。
ここで、3次元被写体1001および1011において、点物体が光学系の光軸上に存在し、光軸を通る視点から見た視点画像データ(全焦点画像データ)が以下の式で表されるとする。
Figure 2015192238
ただし、関数δはディラックのデルタ関数である。aは点物体の強度(画素値)であり、bは背景の強度である。
以下、数16を用いて、3次元被写体1001および1011を撮影したときの、Zスタック画像データのZ=Zfの位置のレイヤー画像データを求める。
(3次元被写体1001のZ=Zfにおけるレイヤー画像データ)
物体の厚みが無視できるとすると、数16の視点画像データas,t(X,Y,Zf)は視点の位置(s,t)にかかわらず、I(X,Y,Zf)で表せる。従って、数19を数16に代入し、変形すると以下の式で表されるレイヤー画像データが得られる。
Figure 2015192238
(3次元被写体1011のZ=Zfにおけるレイヤー画像データ)
同様に、物体の厚みが無視できるとすると、数16の視点画像データas,t(X,Y,Zf)はI(X,Y,Zf)の平行移動で表せる。従って、視点(s,t)と被写体と焦点位置の距離dZに応じて数19を平行移動した上で数16に代入し、数17、数18を用いて変形すると以下の式となる。
Figure 2015192238
Zスタック画像データg(X,Y,Zf)には被写体と焦点位置のずれdZによる焦点ぼけが含まれることがわかる。
図11(d)および図11(h)の焦点ぼけ画像データ(レイヤー画像データ)1103、1113は、レンズ面上に設定した多数の視点に対応する複数の視点画像データから求めた、Z=Zfにおけるレイヤー画像データの模式図である。レイヤー画像データ1103は数20に対応する画像データを表し、1113は数21に対応する画像データを表している。
レイヤー画像データ1103では、物体の像の位置がシフトしていない視点画像データに視点重み関数を掛け算し、視点数分だけ統合されるため、ぼけの無い像が得られる。一方、レイヤー画像データ1113では、物体の像の位置がシフトした視点画像データに視点重み関数を掛け算し、視点数分だけ統合されるため、光学系のぼけを含んだ像が得られている。
被写体の焦点位置からのずれdZが存在する場合、各々の視点画像データには物体の像の位置のシフトが発生し、物体の周囲との輝度差は(a−b)となるため、dZと輝度差(a−b)が大きいほど、各視点画像データ間の差が大きくなることが分かる。
(視点画像データを用いた散乱画像データの抽出方法)
続いて、視点画像データを用いた散乱画像データの抽出方法について説明する。
被写体の焦点位置からのずれdZが0であるとすると、既に数10で説明したように、所定の視点(s,t)の視点画像データと、同一の観察角φで偏角θが180度異なる視点の視点画像データとの間で減算すると、散乱画像データを抽出することが出来る。
以下、式を用いて説明する。
(偏角θの異なる視点画像データ間の演算)
図7(a)に示す座標系において、処理対象の視点P0の座標を(x,y)=(s,t)とすると、180度回転した視点P1の座標は(x,y)=(−s,−t)となる。(偏角の回転角が180度以外の角度でも散乱画像データの抽出は可能だが、ここでは最も効果の大きな180度の回転の場合について説明する。)
視点P0の視線方向からZ=Zfの位置が焦点位置となるように観察した視点画像データをIP0(X,Y,Zf)、偏角回転視点P1の視線方向からZ=Zfの位置が焦点位置となるように観察した視点画像データをIP1(X,Y,Zf)とする。また、視点P0の視点画像データから求めた散乱画像データである視点散乱画像データをSP0(X,Y,Zf)とする。
その場合、下記の数22または数23に示す演算により、標本表面の散乱光の情報を抽出することができる。
Figure 2015192238
Figure 2015192238
数22において1/2の係数を掛ける理由は、視点P0および視点P1のそれぞれに対応する視線方向から観察した場合に大きくなる散乱光が絶対値によって同時に抽出され、倍の強度となるのを防ぐためである。数23では視点P1と比較して視点P0の視点画像データの強度の大きな画素を散乱光が強く発生している画素とみなし、差分が正になる画素の値だけを用いる。
以下、数22、数23の演算式により、標本表面の散乱光の情報が抽出できる理由を説明する。既に述べたように、視点を変えたときの視点画像データの透過光の強度変化は少ない一方、標本表面の散乱光による輝度変化は偏角θに従って大きく変化する。従って、視点を変えた視点画像データ間の減算により、画像データ中の透過光成分はキャンセルされ、散乱光成分(散乱光の情報)が抽出できる。透過光成分は完全にキャンセルできなくてもよい。画像データ中の透過光成分の強度が弱められ、相対的に画像データ中の散乱光成分が強調された画像データを得ることができれば、そのような操作も散乱画像データの生成ということができる。
偏角回転角が180度である場合、視点P0の偏角で散乱光が最大となる凹凸方向角を持つ表面凹凸の散乱光は、視点P1の偏角で最小となる。同様に、視点P1の偏角で散乱光が最大となる凹凸方向角を持つ表面凹凸の散乱光は視点P0の偏角で最小となる。
従って、視点P0の視点画像データと視点P1の視点画像データの差分D(X,Y,Zf)を求めれば、視点P0およびP1の偏角で散乱光が最大となる凹凸方向角を持つ表面凹凸の散乱光の情報を画像データとして抽出できる。
数22、数23において、視点散乱画像データSP0(X,Y,Zf)の計算に差分D(X,Y,Zf)の絶対値、あるいは正の値を用いる理由は、視点画像データから求めた視点散乱画像データを複数統合する過程で非線型性が必要となるためである。絶対値をとらない場合、様々な視点で視点画像データ間の差を求めた後に、その差分の和を取る操作は、様々な視点の視点画像データの和を取った後に差を取る操作と等しくなる。その場合、視点画像データ同士が互いに打ち消し合って、得られる統合像は実質的に0となる。従って、視点散乱画像データSP0(X,Y,Zf)では前記の打ち消し合いを防ぐため、非線型な関数を用いている。
なお、非線型な関数とは以下の式に示すように、所定の画像データ(視点画像データ)に関数を適用した後で複数加算した結果と、所定の画像データ(視点画像データ)を複数加算した後で関数を適用した結果とが、一致しない性質をもつ関数であるとする。
Figure 2015192238
ただし、Nは視点の数、D(i)は視点iにおける視点画像データと偏角回転視点画像データの差分画像データ、S()は関数であり、数22、数23のSP0(X,Y,Zf)を求める式に対応する。
続いて、様々な視点に対して、数22、数23の視点散乱画像データSP0(X,Y,Zf)を積分した画像データ、即ち、散乱画像データDS(X,Y,Zf)について考える。散乱画像データDS(X,Y,Zf)を求める式は以下で表現される。
Figure 2015192238
ただし、k(s,t)は任意の視点重み関数であるとする。
様々な視点(s,t)に対して数25の積分を行えば、様々な方向の散乱光成分を集めることができる。
なお、数22の演算では、視点P0の視点散乱画像データIP0(X,Y,Zf)と視点P1の視点散乱画像データIP1(X,Y,Zf)が等しくなる。よって、積分すべき全視点(s,t)のうちの半分の視点について視点散乱画像データを求め、それらの重み付き積分を2倍することで、数25と同じ結果がえられる。つまり、数22の場合は積分の計算量が約1/2になるというメリットがある。一方、数23は散乱光が抽出される位置に観察方向の情報も含むため、観察方向と散乱光強度の関係を把握する場合に有効である。ただし、レンズ面上の可能な範囲で偏角θ、観察角φを変えた多数の視点に対して数23の演算を行う場合には、数22、数23のいずれを用いても結果は同一となる。
(観察角φの異なる視点画像データ間の演算)
既に数11で説明したように、ある視点(s,t)の視点画像データと、同一の偏角θで観察角φが異なる視点あるいは視点(0,0)の視点画像データとの間の減算によっても散乱画像データを抽出することができる。
観察角φが異なる視点画像データの間でも、視点P0から観察角変更視点P1を求め、数22または数23に示す演算を行った後、数25を実行することで散乱画像データが生成できる。しかし、観察角φを変更する場合、別の演算でも散乱画像データを生成するこ
とができる。以下に説明する。
図7(a)に示す座標系において、処理対象の視点P0の座標を(x,y)=(s,t)とすると、観察角φが0となる視点P1の座標は(x,y)=(0,0)となる。(観察角を変更した角度は0度以外でも散乱画像データの抽出は可能だが、ここでは最も効果の大きな観察角が0度となる視点P1を用いる場合について説明する。)
視点P0と視点P1が前記の関係にあるとき、下記の式で散乱光が抽出できる。SP0(X,Y,Zf)は視点P0の視点散乱画像データである。
Figure 2015192238
偏角θを回転する場合とは異なり、必ずしも非線型な演算を用いる必要はない。なぜなら、視点P0と視点P1の観察角φが異なる場合には、偏角θが異なる視点(すなわち、図7(a)で動径rが等しくなる視点)で積分を計算した結果は等しくならないため、数26の計算で0になることはないためである。
視点P0の動径rをr1、視点P1の動径rをr2とすると、それぞれ偏角θだけを変えた積分を実行して得られる散乱画像データは以下の式で表わされる。
Figure 2015192238
動径がrとなる視点位置の積分は、動径がrの視点にのみ値を持つ視点重み関数で表現できる。そのため、数27を変形すると以下の式で変形できる。
Figure 2015192238
ただし、kex(s,t)は以下となる。
Figure 2015192238
なお、視点重み関数ka1(s,t)、ka2(s,t)の全視点に対する積分はそれぞれ1となることから、kex(s,t)の全視点に対する積分は0となる。
Figure 2015192238
既に数11で述べたように、差分で散乱光の情報を抽出するためには、IP0(X,Y,Zf)の観察角φがIP1(X,Y,Zf)の観察角φよりも大きい必要がある。
従って、一般化すると、散乱光の情報を抽出するためのkex(s,t)の条件は、所定の動径rthを閾値として、rthよりも大きな動径を持つ領域(すなわち観察角φが大きな領域)では正の値を持ち、rth以下の領域では負の値を持つこととなる。ここで閾値rthが小さくなるほど生成される散乱画像データの強度は大きくなる。効果的な条件の一例としてrth=0が挙げられる。kex(s,t)の条件を以下の式で表現できる。
Figure 2015192238
ただし、rthは、rth≦rを満たす所定の値とする。rは、散乱画像データの計算に用いる視点のうちレンズ面上で最も外側の視点、の原点(0,0)からの距離である。
数30および数31の条件を満たせば、kex(s,t)は自由に設計できる。
数28から、散乱画像データはka1(s,t)およびka2(s,t)を用いて生成される異なる焦点ぼけを持つ2つの任意焦点ぼけ画像データの差分で計算できることが分かる。また、数30および数31の条件を満たすkex(s,t)を用いて生成される任意焦点ぼけ画像データによっても計算できることが分かる。なお、以降、数29によって計算される視点重み関数、または数30および数31の条件を満たす視点重み関数kex(s,t)を、散乱光情報抽出用の視点重み関数と呼ぶこととする。
(散乱画像データ中の散乱光の強調および焦点ぼけの低減)
既に述べたように図10(d)に示すように被写体が焦点位置から離れた位置にある場合、視点画像データ中の被写体の位置が変化する。そのため、数25または数28の計算によって得られる散乱画像データ中にも焦点位置から離れた被写体の焦点ぼけの影響が含まれている。
散乱画像データ中の散乱光を強調しながら、焦点ぼけの影響を低減する方法として、下記に述べる式がある。
Figure 2015192238
ただし、CP0(X,Y,Zf)は以下の式で表わされる画像データ、つまり、視点P0の視点画像データに同じ視点P0の散乱画像データを加算した画像データである。以降、CP0(X,Y,Zf)を視点P0の散乱光強調画像データと呼ぶ。
Figure 2015192238
なお、数33のSP0(X,Y,Zf)は数22または数23で表される視点散乱画像データである。
数32は変形すると以下の式で表現できる。
Figure 2015192238
即ち、数32は散乱光強調画像データの間の差分に相当する。
続いて、図12を用いて、数34の散乱光強調画像データの間の差分によって散乱光の強調と焦点ぼけの低減ができる理由を説明する。
図12の(a)の段は、被写体が焦点位置から離れた位置にあり、周囲より輝度値が低い点像である場合の、(1)全焦点画像データ、(2)任意焦点ぼけ画像データ、(2)散乱光抽出画像データ、(4)散乱光強調画像データを示す。また、(b)の段は、被写体が焦点位置から離れた位置にあり、輝度差のあるエッジである場合の、(1)全焦点画像データ、(2)任意焦点ぼけ画像データ、(3)散乱光抽出画像データ、(4)散乱光強調画像データを示す。
なお、散乱光抽出画像データは数20を用いて計算したものであり、散乱光強調画像データは以下の式で表現する任意焦点ぼけ画像データと散乱光抽出画像データの加算値とする。ただし、図12においてb>a>0が成立するとする。
Figure 2015192238
図12の(a)の段にある画像データ(全焦点画像データ1201、任意焦点ぼけ画像データ1202、散乱光抽出画像データ1203、散乱光強調画像データ1204)を順に数式で表わすと下記のようになる。
Figure 2015192238
Figure 2015192238
ただし、h(X,dZ)はX方向とZ方向についてのみ考慮した2次元の焦点ぼけとする。
Figure 2015192238
Figure 2015192238
散乱光強調画像データ1204では、全焦点画像データ1201で高輝度領域(輝度がbの領域)において焦点ぼけが抑制できることが分かる。
一方、図12の(b)の段にある画像データ(全焦点画像データ1211、任意焦点ぼけ画像データ1212、散乱光抽出画像データ1213、散乱光強調画像データ1214)を順に数式で表わすと下記のようになる。
Figure 2015192238
Figure 2015192238
Figure 2015192238
Figure 2015192238
同様に、散乱光強調画像データ1214、全焦点画像データ1211で輝度の高い領域(輝度がbの領域)において焦点ぼけが抑制できるが、輝度の低い領域では焦点ぼけが打ち消せていないことが分かる。しかし、任意焦点ぼけ画像データ1212中に存在する散乱光成分と散乱光抽出画像データ1213中に存在する散乱光成分は互いに加算されることで、より一層散乱光の情報を強めることが可能となる。なお図12には透過光成分だけを示しているため、散乱光成分は示されていない。
既に述べたように数34は、異なる視点重み関数ka1(s,t)、ka2(s,t)を持つ散乱光強調画像データの間の差を表しており、視点重み関数を変えて計算した散乱光強調画像データ間の差分によって散乱光を抽出できることを意味する。従って、数28と同様、ka1(s,t)の強度が観察角φの大きな位置で相対的に大きく、ka2(s,t)の強度が観察角φの小さな位置で相対的に大きい場合には、散乱光成分が効率よく抽出できることが分かる。
続いて、偏角の変更によって得られる散乱画像データ(数25)と散乱光を強調した散乱画像データ(数32)の関係を分析する。数25の視点重み関数k(s,t)を散乱光情報抽出用の視点重み関数kex(s,t)に変更した以下の式について考える。
Figure 2015192238
ただし、SP0(X,Y,Zf)は数22または数23で表わされる視点散乱画像データである。
数10で説明したように観察角φが小さい場合には、数22または数23で計算される視点散乱画像データSP0(X,Y,Zf)に含まれる散乱光の情報は小さい。従って、kex(s,t)が小さな観察角φのみで負の強度を持つならば、数25で計算する散乱光抽出画像データと数44で計算する散乱光抽出画像データの差は小さい。観察角φが0のときのみ負の強度を持つ場合、数25と数44で計算する散乱光抽出画像データの差は0となる。
なお、数32は、観察角φを変えた散乱光抽出画像データを表す数28と偏角θを変えた場合の視点散乱画像データを表す数44の和となっていることが分かる。即ち、数32は偏角θおよび観察角φの変更によって得られる複合的な散乱光画像抽出データを意味していることが分かる。以上、偏角θと観察角φが異なる視点画像データを用いて散乱画像データが生成できることを述べた。
(病理標本の画像データにおける被写体の焦点ぼけ影響)
次に、病理標本から生成した散乱画像データでの被写体の焦点ぼけの影響について述べる。
図12(b)の1213に示すように、被写体が焦点位置から離れた位置にあり、かつ、被写体に大きな輝度差が存在する場合(例えばエッジ)、散乱光抽出画像データには周囲よりも輝度の低い領域に焦点ぼけの滲みだしたような、散乱光成分とは異なるアーティファクト(焦点ぼけ成分)が発生する。数32を用いて焦点ぼけの異なる2つの散乱光強調画像データの差分から散乱光抽出画像データを求める場合にも、図12(b)の1214の低輝度領域(Xが負の領域)に現れる焦点ぼけが互いに異なるため、その差分には同様の現象が発生する。数28で焦点ぼけの異なる2つの焦点ぼけ画像データの間の差分を求めた場合でも同様に焦点ぼけの差分が残る。
上記の現象は特に輝度差が大きなエッジの低輝度領域で顕著になる傾向がある。例えば、HE染色された病理標本の顕微鏡画像では、核や核の内部が濃青で染色され透過率が低く、細胞質の部分は薄いピンクで染色され透過率が高い。
そのため、細胞質の内部では透過率の差が比較的少ない(図12の場合、aとbの輝度差が小さい)ため、焦点ぼけの影響は僅かであり、低輝度領域の焦点ぼけの影響はほとんど問題にはならない。
しかし、核と細胞質の境界では透過率の差が大きい(図12の場合、aとbの輝度差が大きい)ため、被写体の焦点ぼけが低輝度領域(1214の場合、Xが負の領域)に滲みだし、アーティファクトが観察される場合がある。この現象は散乱光抽出画像データ、散乱光強調画像データの双方で観察される。
以下に述べる実施例では、低輝度領域では散乱光が少なく、高輝度領域では散乱光が多い傾向にある現象に注目することで、低輝度領域に現れる被写体の焦点ぼけの影響を抑え、散乱画像データ中の散乱光の情報をより信頼できるものとする方法を提案する。
[実施例1]
以下、画像データ生成装置100の具体的な実施例について説明する。
(散乱画像データ生成設定画面)
図13(a)、図13(b)は実施例1における散乱画像データ生成機能の設定画面の一例である。図2の画像表示アプリケーションで表示画像中の領域207をマウスで選択した後、マウスの右クリックで表示される機能拡張メニュー208から「散乱画像データ生成」(不図示)の項目を選択する。それに対応して散乱画像データ生成処理前後の画像を示す新規ウィンドウ1300(図13(a))および散乱画像データ生成機能の設定画面1303(図13(b))が表示される。ウィンドウ1300の左側領域1301には領域207内の画像データが表示され、右側領域1302には散乱画像データ生成結果の画像データが表示される。
散乱画像データ生成機能の設定を変更する場合には、設定画面1303を操作する。ユーザがマウスにより視点設定ボタン1304を押下すると、散乱画像データ生成に用いる視線方向(3次元的な観察方向)を決める為の設定画面が表示される。なお、視点の数は1つでも複数でも良い。詳細は後述する。ユーザが散乱画像データ計算設定ボタン1305を押下すると、散乱画像データを計算する方法やパラメータを設定する散乱画像データ計算設定画面が表示される。散乱画像データを生成する方法としては既に述べた種々の方式が選択可能であり、詳細は後述する。また必要ならばオプションで散乱画像データに対するノイズ除去パラメータ等の設定を行うことも可能である。詳細は後述する。ユーザが透過光成分抑制設定ボタン1306を押下すると、散乱画像データ中の低輝度領域のアーティファクト(焦点ぼけ成分)を抑制するための設定画面が表示される。詳細は後述する。
オーバーレイ表示1307はチェックボックスであり、この設定を有効にすると右側領域1302には選択領域207内の画像データと散乱画像データが重ね合わせて表示される。ユーザが、必要に応じて上記設定を行い、実行ボタン1308を押すと、散乱画像データの生成が行われ、処理結果が表示される。詳細は後述する。
1310はウィンドウ1300内で右クリックすることで呼び出し可能な機能拡張メニューである。機能拡張メニュー1310の中にはN/C比算出(不図示)等の画像解析の項目が並んでいる。項目を選択すると、画像解析処理の設定画面(不図示)が表示され、ウィンドウ内の選択領域あるいは全体に対して解析処理を実行し、処理結果を表示する。詳細は後述する。
(散乱画像データ生成処理)
図14(a)は前述の実行ボタン1308を押下した際に実行される散乱画像データ生成処理のフローを示す。この処理は、画像表示アプリケーションとそこから呼び出される画像データ生成プログラムによって実現されるものである。
Zスタック画像データ取得ステップS1501では、画像表示アプリケーションで表示中の画像の選択領域207の座標を元に、メインメモリ302または記憶装置130に格納されたZスタック画像データから必要な範囲のデータを取得する。なお、Zスタック画像データが他のコンピュータシステム140に存在する場合にはネットワークI/F304を通じてデータを取得し、メインメモリ302に格納する。
続いて、散乱画像データ計算処理ステップS1502では、被写体に対する視線方向(観察方向)を決める視点の情報に基づき、Zスタック画像データから複数の視点に対する視点画像データを生成する。そして、各視点画像データを用いて散乱画像データを生成する。詳細は後述する。
続いて、輪郭抽出処理ステップS1503では、散乱画像データから輪郭を抽出した輪郭抽出画像データを生成する。なお、ステップS1503の処理は必須ではなく、不図示の設定に従って適用/不適用が変更できる。詳細は後述する。
画像表示処理ステップS1504では、輪郭抽出画像データ、視点散乱画像データ、又は散乱画像データ(散乱光強調画像データ又は散乱光抽出画像データ)を画像表示アプリケーションの表示倍率に合わせて拡大/縮小し、右側領域1302に表示する。オーバーレイ表示1307が有効である場合、選択領域207内の画像データに輪郭抽出画像データ、視点散乱画像データ、あるいは散乱画像データを重ね合わせて表示する。この際、選択領域207内の画像データに対し、対応する位置の視点散乱画像データあるいは散乱画
像データを合成(加算や減算など)した画像を右側領域1302に表示しても良い。またさらに合成した画像データを、選択領域207内の画像データと明るさが近くなるように諧調補正した画像データを右側領域1302に表示しても良い。複数の視点散乱画像データを一定の時間間隔で切り替えるアニメーション表示を行ってもよい。このとき、輪郭抽出画像データ、視点散乱画像データ、あるいは散乱画像データはチャネル(RGB)毎に色を付けて表示しても良いし、標本の色と重ならない別の色に変更しても良い。ここで表示に用いられる画像データ(輪郭抽出画像データ、視点散乱画像データ、散乱画像データ、又は、それらを元画像データに合成した画像データ)は、いずれも画像観察や画像診断(自動診断を含む)に適した観察用画像データである。
(散乱画像データ計算処理ステップS1502)
次に、実行ボタン1308を押下して実行される散乱画像データ生成処理のメインの処理にあたる散乱画像データ計算処理ステップS1502の内部処理について説明する。
図14(b)は本実施例の散乱画像データ計算処理ステップS1502の内部処理を表すフローチャートであり、図14(c)は散乱画像データ計算処理の機能を実現するブロック図である。図14(c)の透過光成分抑制マスク生成部1701、散乱光抽出画像データ生成部1702、透過光成分抑制処理部1703は、それぞれ図14(b)のステップS1601、S1602、S1603の機能を実現するブロックである。
以下、図14(b)および図14(c)を用いて、散乱画像データの生成方法を説明する。
まず、透過光成分抑制マスク生成ステップS1601では、透過光成分抑制マスク生成部1701が、入力されるZスタック画像データを用いて、透過光成分抑制マスクを生成し、後段の透過光成分抑制処理部1703に出力する。
透過光成分抑制マスクは焦点位置Z=Zfの画像データと同じXY方向のサイズを持ち、各画素に係数が割り当てられたマスクデータである。詳細は後述する。
続いて、散乱光抽出画像データ生成ステップS1602では、散乱光抽出画像データ生成部1702が、入力されるZスタック画像データを用いて、散乱光抽出画像データを生成し、後段の透過光成分抑制処理部1703に出力する。詳細は後述する。
最後に、透過光成分抑制処理ステップS1603では、透過光成分抑制処理部1703が、透過光成分抑制マスクと散乱光抽出画像データを用いて、透過光成分が抑制された補正散乱光抽出画像データを生成する。詳細は後述する。
(透過光成分抑制マスク生成ステップS1601)
続いて、透過光成分抑制マスク生成ステップS1601の内部処理について説明する。図15(a)はステップS1601の内部処理を示すフローチャートである。
透過光成分抑制マスク生成ステップS1601では、透過光成分抑制設定(1403)に基づき、透過光成分抑制マスクを生成する。設定も交えながら処理を説明する。
図13(e)の透過光成分抑制設定画面1403は、透過光成分抑制設定ボタン1306の押下時に表示される設定画面の一例である。ここでは前述の低輝度領域における被写体の焦点ぼけ成分を抑制するための情報を設定する。
設定画面1403では、透過光成分の抑制方法を指定する「方式」と透過光成分抑制マスクの生成時に参照する参照画像データを指定する「参照画像」、参照画像データに対して適用する諧調補正の種類を選ぶ「補正特性」がドロップダウンリストで選択できる。
ステップS1601ではこれらの設定情報に基づき、透過光成分抑制マスクを生成する。
まず参照画像データ生成ステップS1801では、入力されるZスタック画像データから、透過光成分抑制マスク生成用の参照画像データを生成する。参照画像データとしては、元画像データであるZスタック画像データの輝度特徴を表す画像データであれば、いろ
いろな種類又は特性の画像データを利用できる。設定画面1403の「参照画像」において所望の参照画像データの種類又は特性を設定できる。典型的には、Zスタック画像データから生成される、光軸を通る視点から見た視点画像データ、即ち全焦点画像データを参照画像データとして用いるとよい。あるいは、Zスタック画像データから生成される任意焦点ぼけ画像データ(つまり、元のレイヤー画像データとは異なる焦点ぼけをもつ画像データ)でもよい。例えば、元のレイヤー画像データよりも被写界深度を拡大した画像データなどを用いるとよい。あるいは、Zスタック画像データの中から選ばれたいずれかのレイヤー画像データを参照画像データとしてもよい。また、両側テレセントリックな光学系を持つ顕微鏡で取得したZスタック画像データの場合、全焦点画像データにセンサの固定パターンノイズが表れやすい。そのような場合には、全焦点画像データに対しメディアンフィルタ等のノイズ低減処理を行っても良い。また、全焦点画像データは、Zスタック画像データ内の複数のZ位置のレイヤー画像データの各々からXY方向のコントラストの高い領域を抽出し、それらを合成することで生成しても良い。なお、任意視点画像データや任意焦点ぼけ画像データの生成方法については既に述べているため、説明は省略する。
次に輝度変換ステップS1802では、ステップ1801で生成した全焦点画像データを輝度に変換し、マスク参照用の参照画像データを生成する。なお、ステップS1801とステップS1802の順序を入れ替え、Zスタック画像データに対し先に輝度変換を実行してから、全焦点画像データを生成しても良い。ここで輝度の変換処理について述べる。Zスタック画像データがモノクロ画像データである場合には輝度値は画素値そのままとする。Zスタック画像データがRGBカラー画像である場合には、輝度値(Y)は各色の画素値(それぞれR,G,B)から数式によって求める。例えば、アナログテレビ放送で用いられるNTSC信号の場合、RGB値から輝度値(Y)への変換式は、
Y=0.299×R+0.587×G+0.114×B
となる。なお、RGB値から輝度値への変換式は信号を表現する色空間(sRGBやAdobeRGB等)によって異なるため、それぞれの色空間に対応した変換式を用いると良い。
次に諧調補正ステップS1803では、参照画像データに対して諧調補正を行う。ここでは設定画面1403の「補正特性」で設定された諧調補正が適用される。諧調補正の例としては「トーンカーブを用いた補正」や「2値化」や「適応的2値化」などがある。2値化は所定(固定)の閾値による2値化であり、適応的2値化は画像データの輝度分布に応じて動的に適切な閾値を選ぶ2値化である。2値化の閾値はユーザが直接設定しても良い。HE染色された病理標本に適用する場合、核と細胞質の領域を分離できる閾値を設定することが好ましい。色成分を用いて核と細胞質の領域を識別する閾値を決定する場合には、ステップS1802で輝度変換を行わず、ステップS1803ではRやBなどの特定の色成分を入力しても良い。2値化では閾値よりも輝度の高い領域は最大輝度(255)に、輝度の低い領域は最低輝度(0)に設定する。なお、諧調補正ステップS1803は必須の処理ではない(設定画面1403で「諧調補正なし」が選択された場合、ステップS1803はスキップされる)。
最後にマスク係数算出処理ステップS1804では、諧調補正した参照画像データから以下の式を用いて透過光成分抑制マスクM(X,Y,Zf)を生成する。
Figure 2015192238
なお、Iref(X,Y,Zf)は参照画像データを表し、Imaxは最大輝度値(8ビット画像データの場合は255)を表す。すなわち、透過光成分抑制マスクは、参照画像データ(又は諧調補正した参照画像データ)の各画素の輝度値を最大輝度値で正規化した画像データである。
図15(b)は参照画像データの入力輝度(横軸)に対する透過光成分抑制マスクの出力係数(縦軸)の関係を表すグラフである。1901は諧調補正なしの場合、1902はトーンカーブで諧調補正した場合、1903は2値化処理を施した場合の関係を示す。いずれも透過光成分抑制マスクの各画素の係数は0〜1の値域を持ち、低輝度領域では低い係数を持ち、高輝度領域では高い係数を持つ。すなわち、透過光成分抑制マスクの各画素の係数の値は、参照画像データの各画素の輝度と正の相関をもつように設定されている。
図15(b)に示すような透過光成分抑制マスクを用いることで、ステップS1602で生成する散乱光抽出画像データの低輝度領域における被写体の焦点ぼけ成分を抑制することが可能となる。詳細は後述する。
(散乱光抽出画像データ生成ステップS1602)
図16(a)は散乱光抽出画像データ生成ステップS1602の内部処理を表すフローチャートである。
まず散乱光抽出画像データ生成方式決定処理ステップS2001では、設定画面1402で選択した「方式」に基づき、散乱光抽出画像データ生成方式を決定する。前述のように方式が異なると後段の散乱画像データ生成実行処理ステップS2002での計算方法が異なる。なお、設定画面1402では方式として「偏角変更」、「観察角変更」、「複合変更」などが選べる。続いて、散乱画像データ生成実行処理ステップS2002では、設定画面1401および1402で設定されたパラメータに基づいて散乱画像データを生成する。
(散乱画像データ生成実行処理ステップS2002:偏角変更又は複合変更)
図16(b)は設定画面1402で「方式」として「偏角変更」あるいは「複合変更」を選択した場合の散乱画像データ生成実行処理ステップS2002の内部処理を示すフローチャートである。「偏角変更」を選択した場合には数44に相当する演算を行い、「複合変更」を選択した場合には数32に相当する演算を行う。
まず、視点取得処理ステップS2101では、後段のステップS2102において視点画像データの生成に必要な視点の位置情報を取得する。ステップS2101では予め定められた視点の位置情報を、メインメモリ302や記憶装置130、他のコンピュータシステム140から取得しても良い。またステップS2101では画像表示アプリケーション上で設定した情報に基づき視点の位置情報を計算して求めても良い。詳細は後述する。
続いて、視点画像データ生成ステップS2102では、選択領域207のZスタック画像データを元に、ステップS2101で求めた視点に対応する視点画像データを生成する。画像データ生成装置100(CPU301)が実行するこの機能を視点画像データ生成手段とよぶ。なお、Zスタック画像データから任意の視点画像データを生成する手法(MFI任意視点/焦点ぼけ画像生成法)には、前述の特許文献1、非特許文献2、3、4等の方法をはじめとして、如何なる方法を用いてもよい。
続いて、視点散乱画像データ生成処理ステップS2103では、生成した視点画像データに対し、散乱画像データ計算設定(1305)に基づき、散乱画像データ生成処理を行う。視点が複数存在する場合は、視点数分の視点散乱画像データ生成処理を実行する。続いて、視点散乱画像データ統合処理ステップS2104では、散乱画像データ計算設定(1305)に基づき、ステップS2103で生成した複数の視点散乱画像データを統合し、散乱画像データを生成する。画像生成装置100(CPU301)が実行するステップS2103−S2104の機能を散乱画像データ生成手段(特徴画像データ生成手段)とよぶ。ただし、被写体の散乱光を抽出する場合には散乱光抽出画像データ生成手段、被写体の散乱光を強調する場合には散乱光強調画像データ生成手段とよび、両者を区別しても良い。詳細は後述する。
以下、視点取得処理ステップS2101、視点散乱画像データ生成処理ステップS2103、視点散乱画像データ統合処理ステップS2104の詳細を説明する。
(視点取得処理ステップS2101)
以下、視点設定(1304)に基づき、視点取得処理ステップS2101で視点の位置情報を計算する処理について説明する。
図13(c)の視点設定画面1401は、図13(b)の視点設定ボタン1304の押下時に表示される設定画面の一例である。ここでは散乱画像データの生成に用いる視点画像データの視点位置を設定する。
視点設定画面1401では視点の設定方法として「直接設定」と「メッシュ設定」の2つが選択できる。直接設定では、視点数と視点の位置(s,t)をユーザに直接指定させる。一方、メッシュ設定では、外径、内径(中心遮蔽)、離散化ステップをユーザに指定させ、これらの指定値から各視点の位置を計算する。
「外径」には、計算する視点の最大のズレ量(光軸とレンズ面が交わる位置を原点とした場合の視点の原点からの距離)を指定し、「内径(中心遮蔽)」には、計算する視点の最小のズレ量(つまり計算しない視点の最大のズレ量)を指定する。ここでは、レンズ面上の原点を中心とする距離(半径)によって、外径および内径(中心遮蔽)の値を設定する。なお、外径にはレンズ面上での光学系の半径rを越える値は設定できない。「離散化ステップ」は、「外径」で規定される円から「内径」で規定される円を除いたドーナツ状の領域内に、視点画像データを生成する視点の位置を離散的に設定するための刻み間隔である。離散化ステップが細かい程、計算する視点数は増加する。
なお上述の円以外にも様々な形状の設定が可能である。例えば、半径の異なる複数の同心円や中心から放射線上に伸びる直線が設定できる。同心円設定の場合、各々の円の半径、各々の円上の視点の密度を決める離散化ステップ(例えば、角度間隔の設定)が設定できる。また中心から放射線上に伸びる直線の場合、線の間隔(例えば、角度間隔の設定)や放射線上の視点の密度を決める離散化ステップが設定できる。
(視点散乱画像データ生成処理ステップS2103)
図13(d)の散乱画像データ計算設定画面1402は、図13(b)の散乱画像データ計算設定ボタン1305押下時に表示する設定画面の一例である。ここでは、散乱光強調画像データまたは散乱光抽出画像データの生成の種別、その計算方式や視点重み関数などのパラメータを設定する。
設定画面1402の「種別」では、ユーザの使用目的に応じて「種別」の下にあるグループボックスから「散乱光強調画像」または「散乱光抽出画像」のいずれかを選択する。
「方式」の下のドロップダウンリストからは視点散乱画像データ生成処理ステップS2103で用いる散乱画像データの計算方法が選択できる。前述のように本実施例では「偏角変更」、「観察角変更」、「複合変更」の中から計算方法が選択できる。
「視点重み関数1」および「視点重み関数2」の下にあるドロップダウンリストからは、視点散乱画像データ統合処理ステップS2104で用いる視点重み関数を選択することができる。ドロップダウンリストには、ガウス関数で表わされるガウスぼけや円柱形で表わされる円柱ぼけ等の焦点ぼけの形状とその半径が組み合わされたパラメータ設定が複数表示され、ユーザはその中から選ぶことができる。なお、視点散乱画像データ生成処理ステップS2103では「視点重み関数1」と「視点重み関数2」のドロップダウンリストで選択した視点重み関数を用いて、散乱光情報抽出用の視点重み関数を生成する。(「視点重み関数1」は数29のka1(s,t)に対応し、「視点重み関数2」は数29のka2(s,t)に対応する。そのため、散乱光情報抽出用の視点重み関数kex(s,t
)は両者の差分で計算できる。)
また、視点散乱画像データ計算設定画面1402ではノイズ除去設定が可能である。ノイズ除去設定としては、閾値による2値化や、メディアンフィルタ、エッジを保ったノイズ除去が可能なバイラテラルフィルタなどが適用できる。この処理により明瞭なコントラストを持つ散乱画像データを生成し、よりN/C比を検出しやすくすることができる。
(設定できる視点重み関数の例)
図25(a)および図25(c)は視点重み関数の例、図25(c)は散乱光情報抽出用の視点重み関数の例を示す模式図である。図25(a)は円柱形で表される視点重み関数、図25(b)はガウス関数で表される視点重み関数である。rは任意焦点ぼけ画像の生成に用いる視点のうち、レンズ面上で最も外側の視点の原点(0,0)からの距離を表す。
図25(a)の視点重み関数k(s,t)は以下の式で表現できる。
Figure 2015192238
また、図25(b)の視点重み関数k(s,t)は以下の式で表現できる。
Figure 2015192238
σはガウス関数の標準偏差を表す。σによりぼけの広がりを制御できる。
数46および数47に示すk(s,t)の積分は1であり、数18の条件を満たす。
撮像光学系の3次元焦点ぼけは、波動光学的なぼけや各種収差の影響により、厳密に視点重み関数で表現できる訳ではないが、数47に示すガウス関数の視点重み関数により比較的良く近似できる。図25(c)の視点重み関数は、図25(a)の視点重み関数から図25(b)の視点重み関数を減算したものであり、数30および数31の条件を満たしている。即ち、図25(c)の視点重み関数は、すべての視点の積分は0となるとともに、視点の動径rが所定の閾値rth以下の領域の積分は負となり、所定の閾値rthより大きな領域の積分は正となっている。
次に、視点散乱画像データ生成処理ステップS2103における視点散乱画像データの生成方法を説明する。
既に数10で説明したように、所定の視点(s,t)の視点画像データと、同一の観察角φで偏角θが180度異なる視点の視点画像データとの間で減算すると、標本の透過率の違いによる濃淡の情報は打ち消され、効率的に散乱画像データを生成できる。
図16(c)は本実施例における視点散乱画像データ生成処理S2103の内部処理を示すフローチャートである。以下、図16(c)を用いて視点毎の散乱画像データの生成方法を説明する。
まず偏角回転視点算出ステップS2301では、処理対象の視点に対し、偏角θが所定の角度だけ回転した位置にある、偏角回転視点を算出する。なお、回転角度としては180度が最も好ましいため、以降は180度の場合について説明する。
図7(a)に示す座標系において、処理対象の視点P0の座標を(x,y)=(s,t)とすると、180度回転した偏角回転視点P1の座標は(x,y)=(−s,−t)となる。
次に、偏角回転視点画像データ生成ステップS2302ではステップS2301で計算した偏角回転視点P1から観察した視点画像データを生成する。視点画像データの生成方法は視点画像データ生成ステップS2102で既に説明したため、省略する。なお、視点画像データ生成ステップS2102で既に、偏角回転視点P1の視点画像データを計算済みの場合には、改めて計算せずに、画像データ生成装置100の記憶装置130やメインメモリ303に存在するデータを読み出して利用する。
次に視点散乱画像データ生成演算ステップS2303では、処理対象の視点P0の視点画像データIP0(X,Y,Zf)と偏角回転視点P1の視点画像データIP1(X,Y,Zf)とから、視点散乱画像データを生成し、出力する。具体的には、設定画面1402で方式として「偏角変更」が選択された場合は数22または数23に示す演算が行われ、「複合変更」が選択された場合は数33の演算が行われる。
(視点散乱画像データ統合処理ステップS2104)
視点散乱画像データ統合処理ステップS2104では、複数の視点散乱画像データを統合し、散乱画像データ(散乱光強調画像データあるいは散乱光抽出画像データ)を生成する。
具体的には、設定画面1402で方式として「偏角変更」が選択された場合には数44の演算が行われ、「複合変更」が選択された場合には数32の演算が行われる。
このとき、数44、数32のkex(s,t)は、設定画面1402の「視点重み関数1」および「視点重み関数2」で設定されたパラメータからka1(s,t)、ka2(s,t)を計算した後、数29を用いて計算する。
なお、図16(b)のフローチャートでは、視点ループの処理に入る前に散乱画像データの格納用バッファを0で初期化し、ステップS2103で数44または数32の積分内の演算を実行し、その結果を散乱画像データの格納用バッファに加算しても良い。その場合には、ステップS2104の演算は不要となる。
また、数44の代わりに数25を用い、散乱光情報抽出用の視点重み関数kex(s,t)の代わりに任意の視点重み関数k(s,t)を用いて散乱光抽出画像データを生成することも可能であり、本発明の範疇とする。
また、視点散乱画像データ統合処理ステップS2104では、視点散乱画像データ生成処理ステップS2103同様に散乱画像データに含まれるノイズを除去する処理を行っても良い。その場合、散乱画像データ計算設定画面1402でノイズ除去設定を行う。生成した散乱画像データの0未満の数値は0とする。0未満の領域には散乱光成分はほとんどなく、焦点ぼけ成分が多いためである(この処理もノイズ除去の効果がある)。
(散乱画像データ生成実行処理ステップS2002:観察角変更)
図16(d)は設定画面1402で「方式」として「観察角変更」を選択した場合の散乱画像データ生成実行処理ステップS2002の内部処理を示すフローチャートである。
視点取得処理ステップS2201の処理は図16(b)のステップS2101と同一であるため説明は省略する。
次に視点画像データ生成ステップS2202では、Zスタック画像データを元に、ステップS2201で求めた視点に対応する視点画像データを生成する。画像データ生成装置100(CPU301)が実行するこの機能を視点画像データ生成手段とよぶ。Zスタック画像データから任意の視点画像データを生成する手法(MFI任意視点/焦点ぼけ画像生成法)には、前述の特許文献1、非特許文献2、3、4等の方法をはじめとして、如何なる方法を用いてもよい。
ステップS2201で求めた全ての視点に対してステップS2202の処理を終えたら、視点画像データ統合処理ステップS2203に進む。
視点画像データ統合処理では、設定画面1402で選択した「視点重み関数1」と「視点重み関数2」に基づき、数29を用いて計算した散乱光情報抽出用の視点重み関数を用いてステップS2202で求めた視点画像データを統合する。視点画像データの統合においては数28に相当する計算を行う。これにより、散乱画像データが得られる。
視点画像データ統合処理ステップS2203においても、図16(b)のステップS2104と同様に設定画面1402のノイズ除去設定に従って、ノイズ除去を行っても良い。また生成した散乱画像データの0未満の数値は0とする。
(透過光成分抑制処理ステップS1603)
続いて、図14(b)の透過光成分抑制処理ステップS1603の処理について説明する。
ステップS1603では、ステップS1602で生成した散乱画像データの中に含まれるアーティファクト(焦点ぼけ成分)を抑制する処理を行う。前述のように、ここで問題とするアーティファクトは、散乱画像データ(特徴画像データ)のうち、Zスタック画像データ(撮影して得られた元画像データ)において低輝度な領域に顕著に現れる。そこで本実施例では、散乱画像データのうち少なくとも元画像データにおける低輝度領域に対応する画素に対して、輝度を低下させる補正処理を行う、という操作により、上記アーティファクトを低減する。輝度補正後の散乱画像データ(特徴画像データ)を補正画像データと呼ぶ。図17(a)はステップS1603の内部処理を表すフローチャートである。設定も交えながら処理を説明する。
まず、抑制方式決定処理ステップS2401では、図13(e)の透過光成分抑制設定画面1403の「方式」で選択した透過光成分抑制処理の方式(例えば「掛け算」など)に基づき、後段のステップS2402で実行する方式を決定する。
次に、抑制実行処理ステップS2402では、ステップS2401で決定した方式に基づき、散乱光抽出画像データ内のアーティファクト(焦点ぼけ成分)を抑制する。
本実施例の透過光成分抑制処理部1703では、設定画面1403で選択された「方式」である「掛け算」に基づき、透過光成分抑制マスクM(X,Y,Zf)と散乱光抽出画像データDS(X,Y,Zf)の間で掛け算を実行する。これにより、低輝度領域のアーティファクト(焦点ぼけ成分)を抑制した補正散乱光抽出画像データ(補正画像データ)MDS(X,Y,Zf)を生成する。透過光成分抑制マスクM(X,Y,Zf)と散乱光抽出画像データDS(X,Y,Zf)の間の掛け算とは、散乱光抽出画像データの各画素の画素値に対し、透過光成分抑制マスクの対応画素の係数を乗じる操作である。式では下記で表わされる。
Figure 2015192238
図17(b)は方式として「掛け算」を選択した場合の、抑制実行処理ステップS2402内の処理および効果を示す1次元の模式図である。図17(b)を用いて詳細を説明する。
図17(b)の2501は、輝度差の大きなエッジを持つ、焦点位置から離れた位置に存在する被写体の参照画像データを表し、2502は透過光成分抑制マスク生成部1701で計算した透過光成分抑制マスクM(X,Y,Zf)を表す。マスク2502では低輝度領域(Xが負の領域)での係数は小さくなっている。マスク2502の係数は0〜1の値をもつ。なお、参照画像データ2501としては、前述のとおり、全焦点画像データや、任意の焦点ぼけを持つ任意焦点ぼけ画像データなどを用いることができる。
2503は散乱光抽出画像データ生成部1702で計算した散乱光抽出画像データDS(X,Y,Zf)を表す。エッジの低輝度領域(Xが負の領域)ではアーティファクト(
焦点ぼけ成分)が発生している。なお、2503には散乱光成分も同時に含まれているが、図には表示していない。
2504は数48の演算結果である。透過光成分抑制マスク2502は低輝度領域(Xが負の領域)の係数は、高輝度領域(Xが正の領域)の係数よりも相対的に小さい値になっている。そのため、散乱光抽出画像データ2503に係数を掛けた結果の画像データ2504では、低輝度領域(Xが負の領域)のアーティファクト(焦点ぼけ成分)が低減されている。
なお、図13(c)〜図13(e)に示す設定画面はあくまで一例である。ユーザである病理医が設定に煩わされずに素早く観察・診断できるよう、デフォルトの設定や自動的に最適値が設定される機能を備えることが望ましい。
以上、本実施例における散乱画像データ計算処理(図14(a)のS1502)について述べた。
(輪郭抽出処理)
続いて輪郭抽出処理(図14(a)のS1503)の一例について述べる。
散乱画像データでは被写体中の散乱光成分が強調されているものの、ノイズや信号の強弱が存在する。そこで、より輪郭を見やすくするために、輪郭抽出処理を行う。例えば、散乱画像データを2値化し(2値化の閾値は予め決められた値を用いても良いし、動的に決めても良い)、その後、膨張・縮小処理を繰り返すことにより、輪郭を抽出できる。他にも輪郭抽出方法には様々な公知の技術があり、ここではいずれの方法も適用できる。またさらに細線化処理を追加することで輪郭が存在する位置精度を高めることができる。処理の結果、散乱画像データから輪郭抽出像データが得られる。
(画像データの表示・解析)
続いて画像表示処理S1504を経て、画像表示アプリケーション上に視点散乱画像データ、散乱画像データ、あるいは輪郭抽出像データを示すことで、細胞と細胞の間の細胞境界、細胞と類洞の境目などを分かりやすくできる。それによって病理医は患部組織の3次元構造をイメージしやすくなる。
さらに、ウィンドウ1300内でマウスを右クリックすることにより機能拡張メニュー1310を呼び出し、N/C比(核/細胞質比)算出等の項目を選択することで、画像解析を行うことが出来る。
図18はN/C比算出の処理フローの一例である。
N/C比算出にあたっては左側領域1301内の選択領域207の画像データと輪郭抽出像データの2枚の画像データを用いることを前提とする。以下、画像データ中の核の部分を核領域、核を取り巻く細胞質の部分を細胞質領域、核領域と細胞質領域を合わせた全体を細胞領域と呼ぶ。
まず、核領域決定処理ステップS2601では核領域の決定を行う。例としては以下の方法がある。HE染色では核内は濃青に染色される為、輪郭抽出像内の該当閉領域内に位置する選択領域207内の画素が一定以上の比率で所定の範囲の色域に属すか否かで、核領域か否かを判別できる。判別に用いる比率及び色域は、予め複数のサンプルを用いて学習すれば良い。
続いて、細胞質領域決定処理ステップS2602では細胞質領域の決定を行う。HE染色では細胞質はピンク色に染色される。従って、核領域決定処理と同様、輪郭抽出像データ内の該当閉領域内に位置する選択領域207内の画素が一定以上の比率で所定の範囲の色域に属すか否かで、細胞領域か否かを判別できる。その後、細胞領域から、ステップS2601で核領域と見なされた閉領域を除外することにより、細胞質領域を特定する。こ
こでの判別に用いる比率及び色域も、予め複数のサンプルを用いて学習すれば良い。
自動処理では十分な精度が出ない場合には、ユーザが介在(補助)して領域決定を行っても良い。その場合、ステップS2602の後、GUIでユーザに輪郭、核領域、細胞領域を修正できる設定画面を表示する。
最後にN/C比算出処理ステップS2603では、上記で求めた核領域の面積を細胞質領域の面積で割り、N/C比を求める。
上記で述べたN/C比算出フローはあくまで一例であり、その他、様々な変形や改良が可能である。
(本実施例の利点)
以上述べたように、本実施例では、低輝度領域では散乱光が少ない性質を利用し、散乱光抽出画像データに被写体から計算したマスクを掛け算することで低輝度領域のアーティファクト(焦点ぼけ成分)を低減させることが可能となる。これにより、病理標本の核などに存在する低輝度領域のアーティファクト(焦点ぼけ成分)の影響を気にすることなく、細胞質や細胞境界での表面凹凸(散乱光)を観察したいというユーザの要望にも応えることができる。
なお、本実施例で述べる手法はいずれも、Zスタック画像データから標本の散乱画像データを抽出することができる。そのため、画像処理(コンピュータによる後処理)だけで(つまり、撮影時の光学系や露出条件の変更を行うことなく)、標本を観察する上で有用な細胞膜や細胞境界、細胞と管や腔との境界を明瞭にすることができる。さらに、輝度や色の変化としてほとんど現れない表面凹凸についても、コントラストを高めて可視化することができるという効果もある。
それにより診断に有用な画像の提示、N/C比の算出等の診断支援機能が実現できる。
なお、本実施例では実行ボタン1308が押されたときに散乱画像データ生成処理を実行するようにしたが、図13(b)、図13(c)〜図13(e)に示す設定パラメータが変更される度に散乱画像データ生成処理を実行するようにしても良い。そうすると、設定パラメータの変更に同期してリアルタイムに処理結果が表示されることとなる。この構成の場合には、図13(b)、図13(c)〜図13(e)に示す設定項目を1つの設定画面内に展開して配置すると良い。このような実装形態も本発明の範疇となる。
[実施例2]
上述した実施例1では、散乱画像データの一例として散乱光抽出画像データに対する処理について述べた。それに対し、実施例2では、散乱光強調画像データの低輝度領域のアーティファクト(焦点ぼけ成分)を抑える方法について述べる。なお、既に述べたが、散乱光強調画像データとは、任意焦点ぼけ画像データに散乱光抽出画像データを加算した画像データである。
(散乱画像データ計算設定画面1402)
本実施例では、散乱画像データ計算設定画面1402では「種別」の設定として「散乱光強調画像」を選択する。「方式」は実施例1と同様に、「偏角変更」「観察角変更」「複合変更」の中から選択できる。また「視点重み関数1」および「視点重み関数2」としは実施例1と同様、ドロップダウンリストからガウス関数で表わされるガウスぼけや円柱形で表わされる円柱ぼけ等の焦点ぼけの形状とその半径が組み合わされたパラメータ設定が選択できる。なお、視点設定画面1401および透過光成分抑制設定画面1403に関しては実施例1と同じであるため、説明を省略する。
図19(a)は本実施例における散乱画像データ計算処理ステップS1502の内部処理を示すフローチャートである。また図19(b)は散乱画像データ計算処理の機能を実
現するブロック図である。以下、図19(a)、図19(b)を用いて、散乱画像データ(散乱光強調画像データ)の生成方法を説明する。
本実施例のフローチャート(図19(a))と実施例1のフローチャート(図14(b))の違いは、散乱光強調画像データ生成処理ステップS2704が存在する点である。また本実施例の機能ブロック(図19(b))と実施例1の機能ブロック(図14(c))の違いは、散乱光強調画像データ生成部2804が存在する点である。図19(a)のステップS2701〜S2703は図14(b)のステップS1601〜S1603と同じ処理のため、詳細説明は省略する。また図19(b)のブロック2801〜2803は図14(c)のブロック1701〜1703と同じ機能を有するため、詳細説明は省略する。以降、散乱光強調画像データ生成処理ステップS2704および散乱光強調画像データ生成部2804について詳しく説明する。
(散乱光強調画像データ生成部2804)
図19(b)の散乱光強調画像データ生成部2804の入出力について説明する。散乱光強調画像データ生成部2804には、Zスタック画像データおよび、設定画面1401〜1403の情報から計算された補正散乱光抽出画像データMDS(X,Y,Zf)が入力される。補正散乱光抽出画像データMDS(X,Y,Zf)は例えば数48の演算により得られる。
散乱光強調画像データ生成部2804は、散乱光強調画像データ生成処理ステップS2704の処理を実行する。Zスタック画像データと補正散乱光抽出画像データMDS(X,Y,Zf)とから、低輝度領域のアーティファクト(焦点ぼけ成分)を抑制した補正散乱光強調画像データMComp(X,Y,Zf)が生成され出力される。
図20は散乱光強調画像データ生成処理ステップS2704の内部処理を示すフローチャートである。以下、図20の処理を説明する。
任意焦点ぼけ画像データ生成処理ステップS2901では、まずZスタック画像データから任意焦点ぼけ画像データa(X,Y,Zf)を生成する。Zスタック画像データからの任意焦点ぼけ画像データa(X,Y,Zf)の生成には前述の特許文献1、非特許文献2、3、4等の方法をはじめとして、如何なる方法を用いてもよい。任意焦点ぼけ画像データには、Zスタック画像データを構成するレイヤー画像データと異なる焦点ぼけをもつ画像データと、レイヤー画像データと同じ焦点ぼけをもつ画像データとが含まれる。いずれの画像データも特許文献1、非特許文献2、3、4等の方法で生成できるが、後者の画像データについてはZスタック画像データの中から選択したいずれかのレイヤー画像データをそのまま用いてもよい。
次に、補正散乱光強調画像データ生成処理ステップS2902では、ステップS2703で生成した補正散乱光抽出画像データMDS(X,Y,Zf)を係数αで変倍し、任意焦点ぼけ画像データa(X,Y,Zf)に加算する。これにより、補正散乱光強調画像データMComp(X,Y,Zf)が生成される。式で表現すると下記となる。
Figure 2015192238
なお、αは散乱光の強調の度合いを決める合成係数であり、正の値を持つ。αの設定は設定画面1403でユーザにより変更可能である(不図示)。以降の説明ではα=1の場合について説明する。
設定画面1402の「方式」で「偏角変更」を設定した場合、設定画面1402で設定した「視点重み関数1」が任意焦点ぼけ画像データおよび散乱光強調画像データの生成で用いられる。ステップS2902で生成される補正散乱光強調画像データは以下で表現される。
Figure 2015192238
続いて、図21(a)、図21(b)を用いて散乱光強調画像データ生成処理ステップS2704の処理および効果について説明する。
図21(a)は透過光成分抑制処理部2803で行う処理を1次元で表した模式図である。図21(b)は散乱光強調画像データ生成部2804で行う処理を1次元で表した模式図である。
図21(a)の3001は散乱光抽出画像データDS(X,Y,Zf)を表し、3002は透過光成分抑制マスクM(X,Y,Zf)を表し、3003は補正散乱光抽出画像データMDS(X,Y,Zf)を表す。既に実施例1で説明したように、補正散乱光抽出画像データMDS(X,Y,Zf)では、エッジの低輝度領域(Xが負の領域)のアーティファクト(焦点ぼけ成分)は抑制されている。
図21(b)の3004はZスタック画像データから計算した任意焦点ぼけ画像データa(X,Y,Zf)を表す。任意焦点ぼけ画像データ3004に補正散乱光抽出画像データ3003を加算すると、3005に示す補正散乱光強調画像データMComp(X,Y,Zf)が生成できる。
補正散乱光強調画像データ3005では、図12に示した散乱光強調画像データ1214と比べ、低輝度領域(Xが負の領域)のアーティファクト(焦点ぼけ成分)が抑制される。一方で、高輝度領域(Xが正の領域)では散乱光強調画像データ1214と同様に焦点ぼけが抑制される。また散乱光抽出画像データ3003の高輝度領域には散乱光成分(不図示)が多く含まれるため、補正散乱光強調画像データ3005では散乱光成分が強調されている。
なお、図19(a)に示すフローチャートではステップS2702〜S2704のそれぞれで複数視点の積分を実行している。しかし、数50に示すように、視点毎にステップS2702〜S2704の計算を行った後に積分を行うことも可能であり、その方法も本発明の範疇とする。
なお、設定画面1402の「方式」で「観察角変更」や「複合変更」を設定する場合には、任意焦点ぼけ画像データの生成では設定画面1402で設定した「視点重み関数1」を用いる。また、散乱光強調画像データの生成では設定画面1402で設定した「視点重み関数1」および「視点重み関数2」の設定情報から計算した散乱光情報抽出用の視点重み関数を用いる。その場合でも補正散乱光抽出画像データMDS(X,Y,Zf)を用いるため、低輝度領域のアーティファクト(焦点ぼけ成分)を抑えつつ、散乱光を強調することができる。
(本実施例の利点)
本実施例の方法によれば、散乱光強調画像データにおいて低輝度領域の焦点ぼけの影響を低減させることが可能となる。これにより、病理標本の核などの低輝度領域におけるアーティファクト(焦点ぼけ成分)の影響を気にすることなく、細胞質や細胞境界における表面凹凸や標本内の散乱光を観察したいというユーザの要望にも応えることができる。
[実施例3]
本実施例では透過光成分抑制マスクM(X,Y,Zf)の値域を変更することにより、
低輝度領域のアーティファクト(焦点ぼけ成分)をより抑える方法について説明する。
本実施例では図15(a)のステップS1803において閾値処理を行い、輝度の高い領域は最大輝度(255)に、輝度の低い領域は0ではなく負の最大輝度(−255)に設定する。それによりステップS1804では−1〜+1の値域を持つ透過光成分抑制マスクM(X,Y,Zf)が生成できる。
このマスクを用いることで、低輝度領域では焦点ぼけをより一層低減することができ、高輝度領域では散乱光を強調することができる。なお、本実施例では散乱画像データ計算設定1402の「種別」では「散乱光画像強調」を選択し、「方式」では「観察角変更」を選択する。また透過光成分抑制設定画面1403の「方式」において「掛け算」を選択する。
次に上記の透過光成分抑制マスクM(X,Y,Zf)を用いる場合の効果について説明する。散乱画像データの生成処理の流れは実施例2と同様のため、以降の説明では、図19(b)のブロック図も参照する。
図22(a)〜図22(d)はそれぞれ、透過光成分抑制マスク生成部2801、散乱光抽出画像データ生成部2802、透過光成分抑制処理部2803、散乱光強調画像データ生成部2804の内部処理を1次元で説明した模式図である。
図22(a)では、輝度差のある被写体のエッジの全焦点画像データ(参照画像データ)3101から、−1〜+1の値域を持つ透過光成分抑制マスク3102が生成される。
図22(b)では、設定画面1402で選択した「視点重み関数1」および「視点重み関数2」のパラメータを用いて散乱光情報抽出用の視点重み関数を求め、数28に相当する処理を実行する。3111は「視点重み関数1」で選択したパラメータを設定した場合の任意焦点ぼけ画像データを表す。3101は「視点重み関数2」で選択したパラメータから生成される全焦点画像データを表す。3112は、任意焦点ぼけ画像データ3111と全焦点画像データ3101の差をとり、正の値のみを取った散乱光抽出画像データを表す。
図22(c)は散乱光抽出画像データ3112と透過光成分抑制マスク3102の掛け算によって低輝度領域で符号が反転した補正散乱光抽出画像データ3121が生成される様子を表している。
図22(d)は、任意焦点ぼけ画像データ3111と補正散乱光抽出画像データ3121の加算によって、低輝度領域のアーティファクト(焦点ぼけ成分)が抑えられた補正散乱光強調画像データ3131が生成される様子を表している。
補正散乱光強調画像データ3131では、実施例2の補正散乱光強調画像データ3005(図21参照)と比べて、低輝度領域(Xが負の領域)のアーティファクトがより一層抑制できている。補正散乱光強調画像データ3131の高輝度領域(Xが正の領域)では焦点ぼけは低減されてないが、加算によって散乱光成分(不図示)は強調されている。
本実施例の方法をHE染色された病理標本の撮影画像データに適用する場合には、核が存在する低輝度領域ではアーティファクト(焦点ぼけ成分)の影響を抑制することができ、細胞質などが存在する高輝度領域では散乱光を強調することが可能となる。
(本実施例の利点)
本実施例の方法によれば、実施例2に比べて低輝度領域のアーティファクト(焦点ぼけ成分)をより低減させた散乱光強調画像データが生成できる。これにより、病理標本の核などの低輝度領域におけるアーティファクト(焦点ぼけ成分)の影響を気にすることなく、細胞質や細胞境界における表面凹凸や標本内の散乱光を観察したいというユーザの要望にも応えることができる。
[実施例4]
実施例4では、実施例1の設定画面1402で方式として「観察角変更」を選択した場合の散乱光抽出画像データ生成ステップS1602の処理を高速化する手法について述べる。
実施例1の図16(d)のステップS2203に示すフローチャートで計算する数28の演算式は、Zスタック画像データから、視点重み関数kex(s,t)に対応する焦点ぼけを持った任意焦点ぼけ画像データを生成する処理とみなすことができる。
非特許文献2に基づく方法では、視点画像データの計算ごとに数14に示すフーリエ変換上での掛け算と数15に示すフーリエ逆変換が必要になる。そのため、計算精度を上げるために視点数を増やすと、視点数に比例して計算時間がかかる課題がある。
そのため、視点毎のフーリエ変換、フーリエ逆変換を必要としない、特許文献1あるいは非特許文献3,4の手法を用いることで高速化が可能となる。
(特許文献1の手法による計算)
まず、MFI任意視点/焦点ぼけ画像生成法の一つである特許文献1に記載の方法の計算の概要について説明する。
3次元被写体f(X,Y,Z)、3次元焦点ぼけh(X,Y,Z)、Zスタック画像データg(X,Y,Z)の間に数13の関係が成立するとき、空間上でのコンボリューションは周波数上での積で表わされるため、以下の関係が成立する。(撮像光学系が両側テレセントリックではない場合にはZスタック画像データに対して座標変換処理が必要になるが、既に述べたため説明は省略する。)
Figure 2015192238
F(u,v,w)、H(u,v,w)、G(u,v,w)はそれぞれ3次元被写体f(X,Y,Z)、撮像光学系の3次元焦点ぼけh(X,Y,Z)、Zスタック画像データg(X,Y,Z)の3次元フーリエ変換を表す。
所望の3次元焦点ぼけh(X,Y,Z)の3次元フーリエ変換をH(u,v,w)とする。h(X,Y,Z)で表わされる焦点ぼけを持つZスタック画像データg(X,Y,Z)の3次元フーリエ変換G(u,v,w)は以下で求められる。
Figure 2015192238
ただし、C(u,v,w)は3次元周波数フィルタ(単に「3次元フィルタ」とも呼ぶ)であり、以下で表わされる。
Figure 2015192238
数52より所望の3次元焦点ぼけh(X,Y,Z)を持つZスタック画像データは以下の式で計算できる。
Figure 2015192238
即ち、3次元被写体f(X,Y,Z)の情報を用いずに、Zスタック画像データg(X,Y,Z)と撮像光学系の3次元焦点ぼけh(X,Y,Z)と所望の3次元焦点ぼけh(X,Y,Z)の情報から、所望の焦点ぼけ画像データを生成できる。
(X,Y,Z)を安定的に求めるためには、数53で表わされる3次元周波数フィルタが安定である必要があり、3次元焦点ぼけh(X,Y,Z)は撮像光学系内を通る光束の範囲内の光線によって表現される焦点ぼけである必要がある。
既に数17で説明したように、撮像光学系の3次元焦点ぼけh(X,Y,Z)は撮像光学系に対応する視点重み関数k(s,t)から生成される。図10(b)や図10(e)に示すように理想的な3次元焦点ぼけは焦点位置からの距離が離れるほど拡大される。従って、視点重み関数k(s,t)が関数や画像データで表わされる場合、3次元焦点ぼけh(X,Y,Z)の各Z位置における2次元の焦点ぼけは、関数または画像データを焦点位置からの距離に応じて拡縮することで求められる。
撮像光学系の3次元焦点ぼけh(X,Y,Z)および任意の3次元焦点ぼけh(X,Y,Z)は被写体の撮影像に依存しないため、予め計算することが可能である。本実施例では「視点重み関数1」および「視点重み関数2」で選択可能な視点重み関数に対応するha1(X,Y,Z)およびha2(X,Y,Z)のそれぞれのフーリエ変換であるHa1(u,v,w)、Ha2(u,v,w)を事前に計算する。そして、それらのフーリエ変換を、画像データ生成装置100内の記憶装置130やメインメモリ302に予め格納しておく。また撮像光学系の3次元焦点ぼけh(X,Y,Z)のフーリエ変換H(u,v,w)も同様に事前に計算し、記憶装置130やメインメモリ302に格納しておく。
(特許文献1の方法を用いた散乱光抽出画像データ生成ステップS1602)
図23は特許文献1に記載の方式を用いる場合の、散乱光抽出画像データ生成ステップS1602の内部処理を表すフローチャートである。以下、処理を順に説明する。
まず、3次元周波数フィルタ生成処理ステップS3201では3次元周波数フィルタC(u,v,w)を生成する。以下に3次元周波数フィルタC(u,v,w)の生成方法を説明する。
撮像系の3次元焦点ぼけのフーリエ変換H(u,v,w)を記憶装置130等から読み出す。また、設定画面1402での設定情報をインデックスとして「視点重み関数1」と「視点重み関数2」で選択した視点重み関数に対応する3次元焦点ぼけのフーリエ変換Ha1(u,v,w)、Ha2(u,v,w)を記憶装置130等から読み出す。
「視点重み関数1」および「視点重み関数2」から求まる散乱光情報抽出用の視点重み関数kex(s,t)に対応する3次元焦点ぼけのフーリエ変換は、フーリエ変換の線型性の性質を用いて、以下の式で求められる。
Figure 2015192238
H(u,v,w)と数55で得られるH(u,v,w)を数53に代入することで、散乱光情報抽出用の視点重み関数kex(s,t)に対応する3次元周波数フィルタC(u,v,w)が生成できる。
なお、数29の演算でkex(s,t)を求めてから、焦点位置からの距離に応じた関数または画像データの拡縮によりh(X,Y,Z)を求め、最後にh(X,Y,Z)をフーリエ変換することでH(u,v,w)を求めても良い。
なお、実施例1と同様、kex(s,t)は数30および数31の条件を満たす視点重み関数である。
次に、3次元フーリエ変換処理ステップS3202ではZスタック画像データg(X,Y,Z)を3次元フーリエ変換し、Zスタック画像データの3次元フーリエ変換G(u,v,w)を生成する。
次に、3次元フィルタ適用処理ステップS3203では、ステップS3201で計算し
たC(u,v,w)とステップS3202で計算したG(u,v,w)を数52に代入し、G(u,v,w)を求める。
次に、3次元フーリエ逆変換処理ステップS3204では、G(u,v,w)に対して3次元フーリエ逆変換を適用し、所望の3次元焦点ぼけh(X,Y,Z)を持つZスタック画像データg(X,Y,Z)を生成する。
最後に、レイヤー画像データ取得処理ステップS3205では、計算したg(X,Y,Z)から、焦点位置(Z=Zf)のレイヤー画像データg(X,Y,Zf)を取得し、散乱光抽出画像データDS(X,Y,Zf)として出力する。
以上の処理により、視点画像データを個別に求めることなく、散乱光情報抽出用の視点重み関数kex(s,t)を用いて元のZスタック画像データから散乱光抽出画像データを生成することができる。
(非特許文献3,4の方法による計算)
続いて、MFI任意視点/焦点ぼけ画像生成法の一つである非特許文献3,4に記載の方法の計算の概要について説明する。
非特許文献3,4に記載の方法では、以下の式により、Z位置Z=Zfにおける任意焦点ぼけ画像データa(X,Y,Zf)のフーリエ変換A(u,v)が生成できる。
Figure 2015192238
数56のnはレイヤー画像データの番号であり0から始まりN−1で終わる。nはZ=Zfに対応する位置のレイヤー画像データの番号を表す。またG(n)(u,v)はZスタック画像データのn枚目のレイヤー画像データg(n)(X,Y)のフーリエ変換であり、以下の式で表わされる。
Figure 2015192238
なお、F{}はフーリエ変換を表す。
また、数56のH(n)(u,v)はレイヤー画像データごとの2次元周波数フィルタ(単に「2次元フィルタ」とも呼ぶ)を表し、以下の式で表わされる。
Figure 2015192238
なお、k(s,t)は所望の3次元焦点ぼけを生成するための視点重み関数を表わす。また、Cs,t(u,v)は撮像光学系の3次元焦点ぼけh(X,Y,Z)の視線方向への積分値cs,t(X,Y,Zf)のフーリエ変換であり、Cs,t(u,v)−1はその逆数である。またe−2πi(su+tv)nは周波数領域で平行移動を行うフィルタである。
数58は被写体の撮影像G(n)(u,v)に依存しない。そのため、H(n)(u,v)を事前に計算し、画像データ生成装置100内の記憶装置130やメインメモリ302に予め格納しておくことで、数56の演算を大幅に高速化することが可能になる。
本実施例の場合、設定画面1402の「視点重み関数1」および「視点重み関数2」で選択可能な視点重み関数に対しては、事前に数58を用いてH(n)(u,v)を求めておき、記憶装置130又はメインメモリ302に予め格納しておく。
(非特許文献3,4の方法を用いた散乱光抽出画像データ生成ステップS1602)
図24は非特許文献3,4に記載の方法を用いる場合の、散乱光抽出画像データ生成ステップS1602の内部処理を表すフローチャートである。以下、処理を順に説明する。
ステップS3301〜ステップS3303はZスタック画像データのレイヤーの枚数分のループの中に存在し、レイヤー画像データごとに処理を行う。なお、レイヤー番号ループでは0からN−1まで順に変更していく。
レイヤー画像データフーリエ変換処理ステップS3301では、Zスタック画像データから処理対象であるレイヤー番号nのレイヤー画像データg(n)(X,Y)を取得し、数57に示すようにフーリエ変換を実行し、G(n)(u,v)を生成する。
レイヤーフィルタ読み出し処理ステップS3302では、「視点重み関数1」と「視点重み関数2」の選択情報に基づき、記憶装置130又はメインメモリ302からそれぞれの視点重み関数に対応するレイヤー画像データごとの周波数フィルタを読み出す。
続いて、視点重み関数1および視点重み関数2のそれぞれに対応するレイヤー画像データごとの周波数フィルタを用いて散乱光情報抽出用の視点重み関数に相当するレイヤー画像データごとの周波数フィルタを生成する。
例えば、視点重み関数1をka1(s,t)、視点重み関数2をka2(s,t)とし、それぞれのレイヤー画像データごとの周波数フィルタをHa1 (n)(u,v)、Ha2 (n)(u,v)とする。その場合、散乱光情報抽出用の視点重み関数kex(s,t)(=ka1(s,t)−ka2(s,t))に対応する周波数フィルタH(n)(u,v)は以下の式で計算できる。
Figure 2015192238
散乱光情報抽出用の視点重み関数kex(s,t)は実施例1と同様、数30および数31の条件を満たす視点重み関数である。
なお、予め、視点重み関数1および視点重み関数2が選択された場合の散乱光情報抽出用の視点重み関数に相当するレイヤー画像データごとの周波数フィルタを生成し、記憶装置130やメインメモリ302に格納しても良い。その場合、レイヤーフィルタ読み出し処理ステップS3302では視点重み関数1および視点重み関数2の情報に基づき、散乱光情報抽出用の視点重み関数kex(s,t)に対応するH(n)(u,v)を読み出す。
次にレイヤーフィルタ適用処理ステップS3303では、読みだしたレイヤー画像データ毎の周波数フィルタH(n−nf)(u,v)とG(n)(u,v)の間で掛け算を実行する。これにより、レイヤー画像データ毎の周波数フィルタ適用結果H(n−nf)(u,v)×G(n)(u,v)を生成する。
未処理のレイヤー画像データがある場合には、レイヤー番号を増加させてステップS3301に戻る。一方、0〜N−1までの全てのレイヤー画像データに対し、上記のステップS3301〜S3303の処理が終了した場合には、ステップS3304に進む。
フィルタ結果統合処理ステップS3304では、すべてのレイヤー画像データに対する周波数フィルタ適用結果を統合する。具体的には数56に相当する処理を行い、A(u,v)を求める。
最後にステップS3305では、A(u,v)にフーリエ逆変換を行い、任意焦点ぼけ画像データa(X,Y,Zf)を生成する。数28で説明したように散乱光情報抽出用の視点重み関数kex(s,t)を用いて生成した任意焦点ぼけ画像データa(X,Y,Zf)は散乱画像データ(散乱光抽出画像データ)DS(X,Y,Zf)となる。そのため、ステップS3305では計算した画像データを散乱光抽出画像データDS(X,Y,Z
f)として出力する。
以上のように、非特許文献3,4の方法を用いる場合でも、視点毎の視点画像データを求めることなく、散乱光情報抽出用の視点重み関数kex(s,t)を用いて元のZスタック画像データから散乱光抽出画像データを生成することができる。
(本実施例の利点)
本実施例の方法によれば、実施例1の図16(d)で説明した処理フローに比べ、計算負荷の大きな視点数分の視点画像データの生成を必要としないため、高速に散乱画像データ(散乱光抽出画像データ)の生成が可能となる。これにより、迅速に細胞質や細胞境界における表面凹凸や標本内の散乱光を観察したいというユーザの要望にも応えることができる。
[実施例5]
実施例5では、散乱画像データ計算設定1402の「方式」において「観察角変更」を選択する場合における、より散乱光を強調した散乱光強調画像データの生成方法について述べる。
本実施例では、散乱光強調画像データも散乱光抽出画像データと同様、図16(d)、図23または図24に示すフローチャートを用いて計算することができる。ただし、適用する周波数フィルタが前述の実施例とは異なる。
散乱光強調画像データに対応する視点重み関数は、所望の任意焦点ぼけ画像データの視点重み関数k(s,t)と散乱光情報抽出用の視点重み関数kex(s,t)の和で求められる。従って、図23のステップS3201または図24のステップ3302において、散乱光情報抽出用の視点重み関数kex(s,t)の代わりに、下記式の関数k(s,t)を用いて周波数フィルタを生成又は読み出せばよい。k(s,t)を、散乱光情報強調用の視点重み関数とよぶ。
Figure 2015192238
なお、k(s,t)の積分値は数18より1であり、kex(s,t)の積分値は数30より0であることから、k(s,t)の積分値は1となる。
Figure 2015192238
また、k(s,t)が以下の関係を満たす場合に、より散乱光の情報が強調される散乱強調画像データが生成できる。
Figure 2015192238
数62において所定の閾値rth以下の視点重み関数の積分が負の値となることは、kex(s,t)を用いた散乱光情報の抽出効果を打ち消さないことを意味する。また閾値rthよりも大きな視点重み関数の積分が1よりも大きいことは、所定の閾値よりも大きな観察角φを持つ、散乱光成分のコントラストが大きな視点画像データを利用することを意味する。
数61および数62の条件を満たす視点重み関数で生成した任意焦点ぼけ画像データは、数61の条件のみを満たす視点重み関数で生成した任意焦点ぼけ画像データに比べ、より散乱光の情報が強調されることとなる。
本実施例において、散乱光強調画像データを生成する場合には、散乱画像データ計算設定1402で「種別」として「散乱光抽出画像」を選択し、「方式」として「観察角変更」を選択する。また、透過光成分抑制設定1403の「方式」で「抑制処理なし」を選択する。「抑制処理なし」を選択する場合には、実施例2の図19(a)ではステップS2701、ステップS2703およびステップS2704の処理は実行されずに、ステップS2702の処理のみ実行されるとする。また図19(b)では2801、2803および2804の各ブロックで処理は実行されず、散乱光抽出画像データ生成部2802のみ実行し結果を出力するとする。
本実施例では、散乱光抽出画像データ生成処理ステップS2702では、以下の式を用いて、散乱光情報強調用の視点重み関数k(s,t)を生成する。
Figure 2015192238
より散乱光の情報を強調するには、散乱光情報強調用の視点重み関数が数61と数62の条件を満たすことが好ましい。そのため、数61と数62の条件を満たす散乱光情報強調用の視点重み関数を生成可能な「視点重み関数1」と「視点重み関数2」の組み合わせを表示するユーザ補助設定が存在しても良い。例えば、「視点重み関数1」として動径rm以下の全ての視点で等しい重みを持つ「半径rmの円柱ぼけ」を選択し、「視点重み関数2」として視点(0,0)にのみ重みを持つピンホールを表す、「ピンホール(0,0)」を選択する。その場合、数63で計算する散乱光情報強調用の視点重み関数k(s,t)は数61および数62の条件を満たすため、ドロップダウンリストの「ピンホール(0,0)」では色を変えて表示しても良い。
ステップS2702の内部処理は実施例1で述べた図16(d)、図23および図24のいずれを用いても良い。
図16(d)の処理を用いる場合には、ステップS2203において、「視点重み関数1」と「視点重み関数2」に基づき数63の式を用いて計算した散乱光情報強調用の視点重み関数を用いて、ステップS2202で求めた視点画像データを統合する。この処理は、非特許文献2の方式に対応する。
次に、特許文献1の方式に対応する図23のフローチャートを用いる場合について説明する。
散乱光情報強調用の視点重み関数に対応する3次元焦点ぼけのフーリエ変換H(u,v)は、ステップS3201において、以下の式で計算できる。Ha1(u,v,w)およびHa2(u,v,w)は、設定画面1402で選択した「視点重み関数1」と「視点重み関数2」の情報に基づいて読み出したそれぞれの3次元焦点ぼけのフーリエ変換である。
Figure 2015192238
そして、数53を用いて、散乱光情報強調用の視点重み関数k(s,t)に対応する3次元周波数フィルタC(u,v,w)が生成できる。
次に、非特許文献3,4の方式に対応する図24のフローチャートを用いる場合について説明する。
散乱光情報強調用の視点重み関数に対応するレイヤー画像ごとの周波数フィルタH(n)(u,v)は、以下の式で計算できる。Ha1 (n)(u,v)およびHa2 (n)(u,v)は、設定画面1402で選択した「視点重み関数1」と「視点重み関数2」の情報に基づいて読み出したそれぞれのレイヤー画像ごとの周波数フィルタである。
Figure 2015192238
以上、散乱画像データ計算設定1402の「方式」において「観察角変更」を選択する場合における、より散乱光を強調するための条件(数61および数62)およびその散乱光強調画像データの生成方法について述べた。
なお、本実施例では、図19(a)では透過光成分抑制マスク生成処理ステップS2701、透過光成分抑制処理ステップS2703および散乱光強調画像データ生成処理ステップS2704の処理を実行しない場合について説明した。しかし、実施例2と同様に、実行することも可能であり、その場合も本発明の範疇とする。
(本実施例の利点)
本実施例の方法を用いることにより、散乱光の情報をより強調した、散乱光強調画像データを生成することができる。また、図23および図24の処理フローを用いれば、高速に散乱画像データ(散乱光強調画像データ)の生成が可能となる。これにより、迅速に細胞質や細胞境界における表面凹凸や標本内の散乱光を観察したいというユーザの要望にも応えることができる。
[実施例6]
実施例6では、散乱画像データ計算設定1402の「方式」で「偏角変更」または「複合変更」を選択する場合の散乱光抽出画像データまたは散乱光強調画像データの高速な生成方法について述べる。
実施例4では観察角φが異なる視点画像データを計算してから散乱画像データを求めるのではなく、フィルタを用いてまとめて計算する方法を述べた。本実施例では偏角θが異なる場合でもフィルタを用いてまとめて計算する方法を述べる。
図26は図8に示す標本の表面凹凸の斜面813を、同一の観察角φで偏角θが180度異なる視点から観察する様子を示す模式図である。直線3501は光軸と平行な直線であり、観察角φは0となる。直線3511と直線3512の偏角θはそれぞれθおよびθ+π([rad])であり、観察角φはいずれもφである。直線3521と直線3522の偏角θはそれぞれθおよびθ+π([rad])であり、観察角φはいずれもφである。
既に数10で述べたように、同一の観察角で偏角が異なる視点画像データ間の差分により、斜面が存在する所定の位置(X,Y,Zf)での散乱光の情報が抽出できる。数10で述べたように、視線方向が直線3511と平行な視点画像データと視線方向が直線3512と平行な視点画像データの差を求めると、斜面813の傾斜角αとB(φ)に比例
する2α×B(φ)の散乱光の情報が抽出できる。同様に、視線方向が直線3521と平行な視点画像データと視線方向が直線3522と平行な視点画像データの差を求めると、2α×B(φ)の散乱光の情報が抽出できる。
本実施例では以降、偏角θが互いに180度(π[rad])異なる視点に対応する視点画像データ間の差分で散乱光の情報を抽出する場合について説明するが、2つの視点の偏角の差は180度に限定しない。例えば、偏角の差が45度の場合でも、数7より見かけ上の傾斜角α’の変化は約1/2であり、数10ではα−α’=α/2となり、偏角の差が180度の場合の約1/4の強度で散乱画像が抽出できる。観察角φの大きさにも依存するが、2つの視点の偏角の差が45度以上ならば、散乱光の情報の抽出に十分な強度が得られるため、本発明の範疇とする。
撮像光学系が両側テレセントリックの場合、視線の観察角φと視点の動径r(=√(s+t))の間には数5の関係がある。(両側テレセントリックでない場合でも変換座標上で考えれば観察角φと動径rの間には数4の関係がある。)そのため、視線の観察角φの代わりに視点の動径rを用いて表現することもできる。
視点の動径がr(視線の観察角がφ)で偏角がθの視点画像データをIr,θ(X,Y,Zf)と表わす。同一の視点の動径r(視線の観察角φ)を持ち、それぞれ偏角がθとθ+π([rad])の視点画像データ間の差の絶対値を、M個の異なる動径rで加算する操作は以下の式で表わされる。
Figure 2015192238
偏角θの方向から観察した標本のそれぞれの位置の見かけ上の傾斜角をαθ(X,Y,Zf)で表わすと、数66は以下のように変形できる。
Figure 2015192238
偏角θを変えない限り、見かけ上の傾斜角は変化しないため、αθi(X,Y,Zf)は一定である。またB(φ)は、φ=0のとき0で、φの増加につれて増加する増加関数である。(数5よりB(φ)はrの関数とみなすことができるので、B(φ)は、r=0のとき0で、rの増加につれて増加する増加関数ともいえる。)
従って、偏角θを変更しない限り、動径rを変更しても数67の絶対値内の正と負の符号は変わることはない。そのため、数67の右辺の絶対値とΣの順序は入れ替えることができ、以下のように変形できる。
Figure 2015192238
よって数66の右辺は以下のように変形できる。
Figure 2015192238
即ち、同一の動径rで偏角θが異なる2つの視点画像データの差分の絶対値を、複数の動径について集めた散乱光の情報は、同一の動径rで偏角θが異なる2つの視点画像データの差分を複数の動径について集めた後で絶対値を取った結果と等しい。
なお、数66の右辺に視点の動径r,偏角θの視点画像データに対応する、正の値を持つ視点重み関数k(r,θ)を掛け算した場合でも絶対値の内部の符号は変化しない。そのため、視点重み関数k(r,θ)も絶対値の中に入れることができ、以下の関係も成立する。
Figure 2015192238
(ただし、数70の視点重み関数は光軸対称であり、動径r,偏角θ+πの視点に対応する視点画像データの視点重み関数の値はk(r,θ)に等しい。)
続いて、数25のSP0(X,Y,Zf)の非線型性を持つ演算を変形し、高速化するために、数70と数25を比較する。以下に数25のSP0(X,Y,Zf)の計算に数22を用いた場合の式を表す。
Figure 2015192238
sおよびtの積分が、視点P0の偏角がθ、視点P1の偏角がθ+πを満たし、視点P0とP1の動径r(観察角)が同一、という条件下でのみ行われるならば、数70と同様、数71の右辺の絶対値と積分の順序は変更できる。
即ち、数71における積分は以下の式で表現できる。
Figure 2015192238
なお、視点P0および視点P1の座標は極座標表現でP0(r,θ)、P1(r,θ+π)とする。
(絶対値と積分の順序を変更する効果)
続いて、絶対値と積分の順序を変更する効果について述べる。
数72のように積分順序を変更できる場合、実施例4と同様に、数72の絶対値の中の式はフィルタ処理で表わすことができる。そのため、多数の視点に対する視点画像データを求める必要がなく、計算を高速化することができる。詳細は後述する。
また、数72の左辺の式を用いて散乱光の情報を抽出する場合、各々の視点画像データにノイズがある場合、差の絶対値を取ることでノイズは正の値に変換され、積分で累積さ
れる。一方、数72の右辺の式を用いて散乱光の情報を抽出する場合、視点画像データにノイズがあっても複数の視点画像データの差の加算で平均化され、その後に絶対値をとるため、ノイズが累積されにくい。そのため、数72の右辺の式を用いる方が、ノイズが少なく高画質な散乱光の抽出情報が得られる。
(極座標で表わされる視点を用いた高速な散乱画像データの生成方法)
図27(a)は本実施例における視点の位置の例を表す図である。各々の視点の位置は以下に示すようにM個の動径rおよび2N個の偏角θの要素を用いて極座標で表現される。
Figure 2015192238
なお、動径rおよび偏角θの存在する範囲はそれぞれ0≦r≦rおよび0≦θ<2πである。図27(a)では動径方向のサンプリング数はM=5、偏角方向のサンプリング数は2N=16(N=8)であり、動径rおよび偏角θは以下で表わされている。
Figure 2015192238
また図27(a)では、以下に示すように、全ての視点において偏角θに対して偏角を180度(π[rad])回転した視点が存在している。
Figure 2015192238
なお、図27(a)に示すような動径方向と偏角方向による視点位置の指定は視点設定1401の設定(不図示)で行うことができる。
(散乱光情報の抽出)
続いて、所定の偏角θとθ+πの2つの視点画像データの間の差の絶対値を複数の動径rについて加算し、その結果を様々な偏角θについて加算することで、図27(a)に示す多数の視点の位置での散乱光の情報を抽出することを考える。数式で表わすと、以下の数76、数77の計算となる。
Figure 2015192238
Figure 2015192238
視点P0(r,θ)と視点P1(r,θ+π)を入れ替えたとき、視点重み関数の値、絶対値内の結果は等しくなるため、数76および数77はそれぞれ以下の式に変
形できる。数79では加算(Σ)の数を2N個からN個に削減することができる。
Figure 2015192238
Figure 2015192238
以降、数76または数78で表わされるDSθi(X,Y,Zf)を偏角選択散乱画像データと呼ぶ。なお、数77または数79で表わされるDS(X,Y,Zf)はこれまでと同様に散乱画像データまたは散乱光抽出画像データと呼ぶ。
ここで数78の右辺は数70を用いて以下の式に変形できる。
Figure 2015192238
続いて、数80の視点の位置(r,θ)を(s,t)で表わし、シグマ(Σ)をインテグラル(∬)で表現すると、以下の数81〜数84のように表現できる。
Figure 2015192238
Figure 2015192238
ここでLθi(s,t)は以下の式で表わされる。
Figure 2015192238
また、lθi(s,t)は以下の式で表わされる。
Figure 2015192238
即ち、lθi(s,t)は、視点(s,t)の偏角θがθであるときは1を取り、それ以外では0となる関数である。ただし、視点(s,t)=(0,0)の場合は、lθi(,0)=0、Lθi(0,0)=0とする。
以降、k(s,t)Lθi(s,t)を偏角選択視点重み関数と呼び、Lθi(s,t)を偏角選択関数と呼ぶ。
数82のSθi(X,Y,Zf)を求める計算は、視点重み関数が偏角選択視点重み関数k(s,t)Lθi(s,t)である任意焦点ぼけ画像データを求める計算と等しいことを表す。
(偏角選択視点重み関数について)
続いて、図を用いて偏角選択視点重み関数k(s,t)Lθi(s,t)の生成方法を説明する。
図27(b)は偏角選択関数Lθi(s,t)の一例を示す模式図であり、図27(c)は偏角選択視点重み関数k(s,t)Lθi(s,t)の一例を示す模式図である。
図27(a)は視点重み関数k(s,t)の模式図であり、既に述べたように複数の動径rと偏角θの組合せでサンプリング位置が決められている。図27(b)はθ=θの偏角選択関数Lθ1(s,t)を表す模式図である。偏角θの線上(3611)では1の値を持ち、偏角θ+πの線上(3612)では−1の値を持ち、それ以外の領域では0となる。なお、原点(3613)では0となる。
図27(c)は、図27(a)に示す視点重み関数k(s,t)と図27(b)に示す偏角選択関数Lθ1(s,t)を掛け算して得られる、偏角選択視点重み関数k(s,t)Lθ1(s,t)を示す模式図である。偏角がθとなる視点(原点を除く)では正の値を、偏角がθ+πとなる視点(原点を除く)では負の値を持ち、それ以外では0となる。すなわち、数82の右辺の式は偏角がθとなる複数の視点に対応する視点画像データの加算、偏角がθ+πとなる複数の視点に対応する視点画像データの減算を行うことを意味する。図27(c)では偏角θ、偏角θ+πとなる視点はそれぞれ4つあるため、偏角選択視点重み関数に対応するフィルタで処理すれば、合計8つの視点画像データの生成と対応する視点画像データの差分をまとめて処理できることになる。
(散乱光抽出画像データの生成フロー)
図28(a)は散乱画像データ計算設定画面1402において「種別」で「散乱光抽出画像」を選択し、「方式」で「偏角変更」を選択する場合の、本実施例における散乱光抽出画像データ生成ステップS1602の内部処理を示すフローチャートである。以下、図28(a)を用いて処理を説明する。(なお、図28(a)に示す散乱光抽出画像データの生成方法は数25の計算に対応する。)
図28(a)の偏角選択散乱画像データ生成処理ステップS3701は、視点の偏角θのループ内にあり、ループで設定された偏角θに対応した偏角選択散乱画像データを生成する。例えば、図27(a)に示す視点位置の場合、θはθからθN−1まで変更され、N回のループ処理が行われる。
なお、偏角選択散乱画像データDSθi(X,Y,Zf)の計算には数81および数82を用いる。偏角選択視点重み関数k(s,t)Lθi(s,t)に対応する任意焦点ぼけ画像データを生成することでSθi(X,Y,Zf)を求めた後、数81によりSθi(X,Y,Zf)の絶対値を求める。詳細は後で説明する。
予め設定した全ての偏角θに対してステップS3701の処理が終了したら、偏角選択散乱画像データ統合処理ステップS3702に進む。
偏角選択散乱画像データ統合処理ステップS3702では、数79に示す演算を行い、ステップS3701で計算した偏角選択散乱画像データDSθi(X,Y,Zf)から散乱光抽出画像データDS(X,Y,Zf)を生成する。
続いて、偏角選択散乱画像データ生成処理ステップS3701の詳細を説明する。図28(b)は偏角選択散乱画像データ生成処理ステップS3701の内部処理を示すフローチャートである。
偏角選択視点重み関数生成処理ステップS3801では、偏角選択視点重み関数k(s,t)Lθi(s,t)を生成する。具体的には、まず、散乱画像データ計算設定画面1402における「視点重み関数1」の選択情報に基づき、ka1(s,t)(極座標表示ではka1(r,θ))を生成する。次に、処理対象である偏角θに対応する偏角選択関数Lθi(s,t)を、画像データ生成装置100内の記憶装置130又はメインメモリ302から読み出す。そして、偏角選択関数Lθi(s,t)と視点重み関数1に対応するka1(s,t)とを掛け算することで、偏角選択視点重み関数k(s,t)Lθi
(s,t)を生成する。
なお、ステップS3801では、偏角選択視点重み関数k(s,t)Lθi(s,t)を生成せず、偏角選択視点重み関数に対応する複数の2次元フィルタや3次元フィルタを、メインメモリ302等から読み出すためのインデックスを生成しても良い。例えばインデックスとして視点重み関数の番号d1と偏角θの番号d2を定める。
次に任意焦点ぼけ画像生成処理ステップS3802では、偏角選択視点重み関数k(s,t)Lθi(s,t)に基づいて、任意焦点ぼけ画像データを生成し、Sθi(X,Y,Zf)として出力する。なお、ステップS3801で偏角選択視点重み関数を読み出すためのインデックスを作成した場合、番号d1およびd2から、偏角選択視点重み関数に対応する複数の2次元フィルタや3次元フィルタを読み出し、任意焦点ぼけ画像データを生成する。
任意焦点ぼけ画像データの生成方法として、図23の処理フロー(特許文献1)、図24の処理フロー(非特許文献3,4)に示すいずれの方法を用いても良い。詳細は後述する。なお、高速化ではなく、ノイズを低減する目的で図16(d)の処理フロー(非特許文献2)に示す方法を用いても良い。
最後に偏角選択散乱画像データ抽出処理ステップS3803では、ステップS3802が出力するSθi(X,Y,Zf)に対し、数81に示す絶対値を求める演算を行い、偏角選択散乱画像データDSθi(X,Y,Zf)を生成する。
続いて、任意焦点ぼけ画像生成処理ステップS3802の詳細について述べる。
(任意焦点ぼけ画像データの生成に図23の処理フローを用いる場合)
図23の処理フロー(特許文献1)を用いる場合の各ステップでの処理について述べる。
ステップS3201では、偏角選択視点重み関数k(s,t)Lθi(s,t)に対応する3次元焦点ぼけのフーリエ変換H(u,v,w)と撮像光学系の3次元焦点ぼけのフーリエ変換H(u,v,w)をメインメモリ302等から読み出す。これらのフーリエ変換は予め計算され、メインメモリ302や記憶装置130に格納されているものとする。偏角選択視点重み関数に対応する3次元焦点ぼけのフーリエ変換H(u,v,w)はステップS3801で定めたインデックスを用いて読み出すことが可能である。
そして、数53に示す演算により、3次元周波数フィルタC(u,v,w)を生成する。
3次元フーリエ変換処理ステップS3202では、処理対象の偏角θ毎にZスタック画像データの3次元フーリエ変換を実行せず、一度計算したフーリエ変換結果をメインメモリ302や記憶装置130に格納し、それ以降は読み出して処理する。
以降の3次元フィルタ適用処理ステップS3203、3次元フーリエ逆変換処理S3204、レイヤー画像データ取得処理S3205の各々の処理は実施例4の説明と同一であるため、説明は省略する。
(任意焦点ぼけ画像データの生成に図24の処理フローを用いる場合)
次に、図24の処理フロー(非特許文献3,4)を用いる場合の各ステップでの処理について述べる。
ステップS3301では、処理対象の偏角θ毎にレイヤー画像データのフーリエ変換を実行せず、一度計算したフーリエ変換結果をメインメモリ302や記憶装置130に格納し、それ以降は読み出して処理する。
続いて、ステップS3302では、偏角選択視点重み関数k(s,t)Lθi(s,t)に対応するレイヤー画像データ毎の周波数フィルタデータH(n−nf)(u,v)を
メインメモリ302等から読み出す。これらの周波数フィルタデータは予め計算され、メインメモリ302や記憶装置130に格納されているものとする。
レイヤー画像データ毎の周波数フィルタデータH(n−nf)(u,v)はステップS3801で定めたインデックスを用いて読み出すことが可能である。又は、ステップS3801で計算した偏角選択視点重み関数k(s,t)Lθi(s,t)を用いて、数58からレイヤー画像データ毎の周波数フィルタデータH(n−nf)(u,v)を計算しても良い。この場合に、Cs,t(u,v)−1やe−2πi(su+tv)nが予め計算され、メインメモリ302や記憶装置130に格納されているとよい。それらを読み出して計算に用いれば、数58は高速に計算できる。
以降のレイヤーフィルタ適用処理ステップS3303、フィルタ結果統合処理ステップS3304、逆フーリエ変換ステップS3305の各々の処理は実施例4の説明と同一であるため、説明は省略する。
以上、散乱画像データ計算設定画面1402の「方式」として「偏角変更」を選択する場合の散乱光抽出画像データの生成方法について述べた。
散乱光抽出画像データ生成ステップS1602の内部で図28(a)および図28(b)の処理フローを用いることで散乱光抽出画像データの計算を高速化することができる。また散乱光抽出画像データではノイズが平均化されて抑えられる効果があり、生成した散乱光抽出画像データのノイズを少なくすることができる。
(「偏角変更」を選択する場合の散乱光強調画像データの生成方法)
次に、本実施例で述べた散乱光抽出画像データから散乱光強調画像データを求める場合、すなわち、散乱画像データ計算設定画面1402において「種別」で「散乱光強調画像」を選択し、「方式」で「偏角変更」を選択する場合について説明する。実施例2と同様、図19(a)を用いて説明する。
なお、説明を簡略にするため、本実施例でも実施例5と同様、透過光成分抑制設定画面1403の「方式」で「抑制処理なし」を選択する場合について説明する。そのため、実施例2の図19(a)ではステップS2701、ステップS2703の処理は実行されずに、ステップS2702、ステップS2704の処理が実行されるとする。また図19(b)では2801、2803および2804の各ブロックで処理は実行されず、2802のみ実行し結果を出力する。なお、本実施例では説明しないが、実施例2と同様、透過光成分抑制設定画面1403の「方式」で「掛け算」を選択した場合も本発明の範疇とする。
まず、散乱光抽出画像データ生成処理ステップS2702では、「視点重み関数1」の設定情報に対応した視点重み関数ka1(s,t)を生成し、散乱光抽出画像データDS(X,Y,Zf)を生成する。本実施例のステップS2702の内部では図28(a)および図28(b)に示す処理フローを行う。既に説明したため、詳細は省略する。
次に散乱光強調画像データ生成処理ステップS2704ではZスタック画像データから任意焦点ぼけ画像データa(X,Y,Zf)を生成し、散乱光抽出画像データDS(X,Y,Zf)と合成する処理を行う。図20を用いてステップS2704の内部処理を説明する。
任意焦点ぼけ画像データ生成処理ステップS2901では、まずZスタック画像データから任意焦点ぼけ画像データa(X,Y,Zf)を生成する。Zスタック画像データからの任意焦点ぼけ画像データの生成方法として、図23の処理フロー(特許文献1)、図24の処理フロー(非特許文献3,4)に示すいずれの方法を用いても良い。
本実施例では、設定画面1403の「方式」で「抑制処理なし」を選択している。それゆえ、補正散乱光強調画像データ生成処理ステップS2902では、補正散乱光抽出画像
データMDS(X,Y,Zf)としてステップS2702で生成した散乱光抽出画像データDS(X,Y,Zf)を用いる。そして数49に示すように、右辺のMDS(X,Y,Zf)を変倍し、任意焦点ぼけ画像データa(X,Y,Zf)に加算することで散乱光強調画像データを求める。
以上、散乱画像データ計算設定画面1402の「方式」として「偏角変更」を選択する場合の散乱光強調画像データの生成方法について述べた。
ステップS2702で図28(a)および図28(b)の処理フローを用いることで散乱光抽出画像データの計算を高速化し、散乱光強調画像データを高速に求めることができる。また、ステップS2702で求めた散乱光抽出画像データではノイズが平均化されて抑えられる効果があるため、散乱光強調画像データもノイズを少なくすることができる。
(「複合変更」を選択する場合の散乱光抽出画像データの生成方法)
次に、設定画面1402の「方式」で「複合変更」を選択した場合の散乱光抽出画像データの生成方法について述べる。前記の設定では、数32に示すように、散乱光情報抽出用の視点重み関数を用いて、「観察角変更」を選択した場合の散乱光抽出画像データと「偏角変更」を選択した場合の散乱光抽出画像データの加算を求める。
数44を用いて「偏角変更」を選択した場合の散乱光抽出画像データを求める場合、「視点重み関数1」および「視点重み関数2」の設定情報から求めた散乱光情報抽出用の視点重み関数kex(s,t)を用いることで計算できた。しかし、本実施例でこれまで述べた数79および数80に示す計算で散乱光抽出画像データを求める場合には、視点重み関数は絶対値の内部に存在するため、線型性が成立しない。そのため、数80の視点重み関数k(r,θ)をkex(r,θ)に置き換えることで散乱光抽出画像データを計算することはできない。従って、数34に示すように視点重み関数ka1(s,t)および視点重み関数ka2(s,t)のそれぞれの散乱光強調画像データを求め、その差をとることで散乱光抽出画像データを求める。
以下に数34を変形して求めた、本実施例における設定画面1402の「方式」で「複合変更」を選択した場合の散乱光抽出画像データを求める式を示す。
Figure 2015192238
即ち、実施例5で述べた散乱光抽出画像データに、視点重み関数ka1(s,t)を用いて求めた散乱光強調画像データと視点重み関数ka2(s,t)を用いて求めた散乱光強調画像データの差を加算することと等しい。既に述べたように実施例4および本実施例でこれまで述べた散乱光抽出画像データの生成フローは、複数の視点画像データを求めることなく、任意焦点ぼけ画像データの生成として処理するため、高速な計算が可能である。
図28(c)は設定画面1402において「種別」で「散乱光抽出画像」を選択し、「方式」で「複合変更」を選択する場合の、本実施例における散乱光抽出画像データ生成ステップS1602の内部処理を示すフローチャートである。
観察角変更散乱光抽出画像データ生成ステップS3901では、実施例4で説明したように設定画面1402で設定した「視点重み関数1」と「視点重み関数2」の設定情報から散乱光情報抽出用の視点重み関数kex(s,t)を生成する。そして散乱光情報抽出用の視点重み関数に対応する任意焦点ぼけ画像データを生成する。なお、任意焦点ぼけ画像データの生成方法として、図23の処理フロー(特許文献1)、図24の処理フロー(
非特許文献3,4)のいずれの方法を用いても良い。ステップS3901で計算した散乱光抽出画像データを観察角変更散乱光抽出画像データと呼ぶ。
次に、第1の偏角変更散乱光抽出画像データ生成ステップS3902では、「視点重み関数1」の設定情報に基づき、ka1(s,t)又はそのインデックス番号を求め、数81と数79を用いて散乱光抽出画像データDS(X,Y,Zf)を生成する。処理フローは図28(a)および図28(b)で説明したため省略する。ステップS3902で計算した散乱光抽出画像データを第1の偏角変更散乱光抽出画像データと呼ぶ。
第2の偏角変更散乱光抽出画像データ生成ステップS3903でも同様に、「視点重み関数2」の設定情報に基づき、ka2(s,t)又はそのインデックス番号を求め、数81と数79を用いて散乱光抽出画像データDS(X,Y,Zf)を生成する。視点重み関数が異なる点を除き、処理内容はステップS3902と同様である。ステップS3903で計算した散乱光抽出画像データを第2の偏角変更散乱光抽出画像データと呼ぶ。
最後に複合変更散乱光抽出画像データ生成ステップS3904では、数85に示す計算を行う。即ち、まず第1の偏角変更散乱光抽出画像データから第2の偏角変更散乱光抽出画像データを減算し、偏角変更散乱光抽出画像データを計算する。次に、求めた偏角変更散乱光抽出画像データを観察角変更散乱光抽出画像データに加算し、複合変更散乱光抽出画像データを生成する。
なお、複合変更散乱光抽出画像データを画像として表示する場合、0未満のデータには散乱光の情報は少なくノイズとみなせるため、0未満は0として表示すると良い。
以上、設定画面1402の「方式」として「複合変更」を選択する場合の散乱光抽出画像データの生成方法について述べた。実施例4および本実施例で述べた計算方法を用いることで散乱光抽出画像データを高速に計算できる。また、既に述べたように偏角変更散乱光抽出画像データではノイズが抑えられるため、生成する散乱光抽出画像データも低ノイズとなる。
なお、設定画面1402において「種別」で「散乱光抽出画像」を選択し、「方式」で「偏角変更」を選択した場合に、散乱光抽出画像データとして、第1の偏角変更散乱光抽出画像データと第2の偏角変更散乱光抽出画像データの差を出力しても良い。この場合の散乱光抽出画像データは数44の計算に相当する。
(本実施例の利点)
本実施例の方法を用いれば、同一の観察角φで異なる視点の偏角θを持つ視点画像データの差分を計算し、それらを加算によって集める場合でも、対応する視点画像データを計算することなく、フィルタ処理によってまとめて計算することができる。そのため、散乱光抽出画像データおよび散乱光強調画像データの生成を高速化することができる。これにより、迅速に細胞質や細胞境界における表面凹凸や標本内の散乱光を観察したいというユーザの要望にも応えることができる。また、本実施例の方法で求めた散乱光抽出画像データおよび散乱光強調画像データでは、複数の視点画像データで加算を行った後で減算および絶対値を取るため、ノイズ成分は抑制される。そのため、弱い散乱光成分であってもノイズ成分に埋もれることなく観察できる。
以上、本発明の好適な実施例を説明したが、本発明の構成はこれらの実施例に限られない。
例えば、上記実施例では明視野顕微鏡で撮影されたZスタック画像データを元画像データとして用いた場合について説明した。しかし、本発明は、落射照明型顕微鏡、ライトフィールドカメラ、ライトフィールド顕微鏡等で撮影された元画像データに対しても適用可能である。
また上記実施例では被写体として病理標本を例に説明してきたが、被写体はそれに限定されない。落射照明型顕微鏡の観察対象である金属等の反射物体でも構わない。また透過観察型顕微鏡の観察対象である透明な生物標本でも良い。いずれの場合においても特許文献1等で開示される技術を用いれば、被写体の奥行き方向の焦点位置を変えて撮影した複数枚のレイヤー画像データ群から任意の視点画像データが生成でき、本発明が適用できる。反射物体を撮影した元画像データを用いる場合、元画像データには反射光(鏡面反射)成分の像と散乱光成分の像とが含まれるが、紙のように光沢が少ない被写体では散乱光が支配的となる。その場合、上記実施例と同様の処理を行うことで、散乱光成分を抽出ないし強調することができる。
また各実施例で説明した構成を互いに組み合わせてもよい。例えば、ステップS2103で視点散乱画像データを生成する際、処理対象の視点と偏角回転視点と観察角変更視点を用いて、それぞれの視点散乱画像データを求め、両者の間で加算等を行って、視点散乱画像データの強度や信頼度を高めても良い。
また上記実施例で述べた各種の画像データの生成処理については、同じ処理を実空間で演算できる場合と周波数空間で演算できる場合とがある。その場合には実空間で演算してもよいし、周波数空間で演算してもよい。すなわち、本明細書において、画像データという用語は、実空間の画像データと周波数空間の画像データのいずれも含む概念である。
また上記各実施例では、画像データの演算を数式で表現しているが、実際の処理では数式どおりの計算を必ずしも行う必要はない。数式で表現された演算結果に相当する画像データが得られるのであれば、具体的な処理やアルゴリズムはどのように設計してもよい。
100:画像生成装置

Claims (17)

  1. 被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを、周波数空間のデータに変換する変換手段と、
    前記変換された周波数空間のデータに対し、フィルタを適用するフィルタ適用手段と、
    前記フィルタが適用されたデータを実空間の画像データに逆変換する逆変換手段と、を有し、
    前記フィルタは、前記逆変換された画像データに含まれるレイヤー画像データのいずれかが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されている
    ことを特徴とする画像データ生成装置。
  2. 被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを用いて、前記複数のレイヤー画像データのそれぞれを周波数空間のデータに変換する変換手段と、
    前記変換された周波数空間の複数のデータに対し、複数のフィルタをそれぞれ適用するフィルタ適用手段と、
    前記フィルタが適用された複数のデータを統合する統合手段と、
    前記統合されたデータを実空間の画像データに逆変換する逆変換手段と、を有し、
    前記複数のフィルタは、前記逆変換された画像データが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されている
    ことを特徴とする画像データ生成装置。
  3. 前記フィルタは、視点の位置ごとの重みを定義する関数である視点重み関数を用いて生成されるデータである
    ことを特徴とする請求項1又は2に記載の画像データ生成装置。
  4. 視点の位置を(s,t)、視点重み関数をk(s,t)と表す場合に、
    前記視点重み関数k(s,t)は、
    Figure 2015192238
    (ただし、rthは予め定められた閾値)
    という条件を満たす関数である
    ことを特徴とする請求項3に記載の画像データ生成装置。
  5. 視点の位置を(s,t)、視点重み関数をk(s,t)と表す場合に、
    前記視点重み関数k(s,t)は、
    Figure 2015192238
    (ただし、rthは予め定められた閾値)
    という条件を満たす関数である
    ことを特徴とする請求項3に記載の画像データ生成装置。
  6. 前記光軸方向に平行な軸まわりの角を偏角と呼ぶときに、
    前記視点重み関数は、第一の偏角を持つ視点の重みが正の値をとり、前記第一の偏角から所定の角度だけ回転した第二の偏角を持つ視点の重みが負の値をとる関数である
    ことを特徴とする請求項3に記載の画像データ生成装置。
  7. 前記フィルタは、
    Figure 2015192238
    (ただし、(s,t)は視点の位置、k(s,t)は視点重み関数、δはディラックのデルタ関数、XYZはZ軸が光軸方向に平行となる直交座標系である。)
    で表される3次元的な焦点ぼけを周波数空間に変換することで得られる、3次元フィルタである
    ことを特徴とする請求項3〜6のうちいずれか1項に記載の画像データ生成装置。
  8. 前記フィルタは、
    Figure 2015192238
    (ただし、(s,t)は視点の位置、k(s,t)は視点重み関数、Aは周波数空間での平行移動を表す関数、Bは被写体を撮影した光学系の3次元的な焦点ぼけを視点(s,t)を通る視線方向へ積分した積分値を周波数空間に変換した値の逆数である。)
    で表される2次元フィルタである
    ことを特徴とする請求項3〜6のうちいずれか1項に記載の画像データ生成装置。
  9. 予め計算された前記フィルタのデータを記憶する記憶装置をさらに有し、
    前記フィルタ適用手段は、前記変換された周波数空間のデータに適用すべきフィルタのデータを前記記憶装置から読み込む
    ことを特徴とする請求項1〜8のうちいずれか1項に記載の画像データ生成装置。
  10. 前記フィルタの特性を決定するパラメータが、ユーザにより変更可能である
    ことを特徴とする請求項1〜9のうちいずれか1項に記載の画像データ生成装置。
  11. 前記特徴画像データの画素のうち、少なくとも前記元画像データにおける低輝度領域に対応する画素に対して、輝度を低下させる補正処理を行うことにより、補正画像データを生成する補正手段をさらに有する
    ことを特徴とする請求項1〜10のうちいずれか1項に記載の画像データ生成装置。
  12. 前記元画像データは、前記被写体を顕微鏡により撮影することにより得られる顕微鏡画像データである
    ことを特徴とする請求項1〜11のうちいずれか1項に記載の画像データ生成装置。
  13. 前記特徴画像データは、前記元画像データに含まれる散乱光成分の特徴を抽出ないし強調した画像データである
    ことを特徴とする請求項1〜12のうちいずれか1項に記載の画像データ生成装置。
  14. 前記特徴画像データは、前記元画像データに比べて、前記被写体の表面の凹凸のコントラストが高められた画像データである
    ことを特徴とする請求項1〜12のうちいずれか1項に記載の画像データ生成装置。
  15. コンピュータが、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを、周波数空間のデータに変換するステップと、
    コンピュータが、前記変換された周波数空間のデータに対し、フィルタを適用するステップと、
    コンピュータが、前記フィルタが適用されたデータを実空間の画像データに逆変換するステップと、を有し、
    前記フィルタは、前記逆変換された画像データに含まれるレイヤー画像データのいずれかが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されている
    ことを特徴とする画像データ生成方法。
  16. コンピュータが、被写体における光軸方向の位置が異なる複数の層を撮影して得られた複数のレイヤー画像データを含む元画像データを用いて、前記複数のレイヤー画像データのそれぞれを周波数空間のデータに変換するステップと、
    コンピュータが、前記変換された周波数空間の複数のデータに対し、複数のフィルタをそれぞれ適用するステップと、
    コンピュータが、前記フィルタが適用された複数のデータを統合するステップと、
    コンピュータが、前記統合されたデータを実空間の画像データに逆変換するステップと、を有し、
    前記複数のフィルタは、前記逆変換された画像データが、前記被写体に対する視線方向が互いに異なる複数の視点画像データのあいだの差異を抽出または強調した画像データに相当する、特徴画像データとなるように、設計されている
    ことを特徴とする画像データ生成方法。
  17. 請求項15又は16に記載の画像データ生成方法の各ステップをコンピュータに実行させることを特徴とするプログラム。
JP2014067033A 2014-03-27 2014-03-27 画像データ生成装置および画像データ生成方法 Pending JP2015192238A (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2014067033A JP2015192238A (ja) 2014-03-27 2014-03-27 画像データ生成装置および画像データ生成方法
US14/659,957 US20150279033A1 (en) 2014-03-27 2015-03-17 Image data generating apparatus and image data generating method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2014067033A JP2015192238A (ja) 2014-03-27 2014-03-27 画像データ生成装置および画像データ生成方法

Publications (1)

Publication Number Publication Date
JP2015192238A true JP2015192238A (ja) 2015-11-02

Family

ID=54191110

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2014067033A Pending JP2015192238A (ja) 2014-03-27 2014-03-27 画像データ生成装置および画像データ生成方法

Country Status (2)

Country Link
US (1) US20150279033A1 (ja)
JP (1) JP2015192238A (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114283128A (zh) * 2021-12-16 2022-04-05 中国人民解放军海军潜艇学院 一种基于多参数阈值的海洋中尺度涡边缘检测方法及系统
WO2023037657A1 (ja) * 2021-09-09 2023-03-16 株式会社島津製作所 データ解析装置およびデータ解析方法

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6455829B2 (ja) * 2013-04-01 2019-01-23 キヤノン株式会社 画像処理装置、画像処理方法、およびプログラム
JP6179551B2 (ja) * 2015-05-12 2017-08-16 コニカミノルタ株式会社 画像検査装置及び画像形成装置
JP6598660B2 (ja) 2015-12-01 2019-10-30 キヤノン株式会社 画像処理装置および画像処理方法
JP2017099616A (ja) * 2015-12-01 2017-06-08 ソニー株式会社 手術用制御装置、手術用制御方法、およびプログラム、並びに手術システム
JP6742724B2 (ja) * 2015-12-28 2020-08-19 シスメックス株式会社 細胞領域決定方法、細胞撮像システム、細胞画像の処理装置及びコンピュータプログラム
CN107645625B (zh) * 2016-07-22 2020-10-09 松下知识产权经营株式会社 图像生成装置以及图像生成方法
CN106991696B (zh) * 2017-03-09 2020-01-24 Oppo广东移动通信有限公司 逆光图像处理方法、逆光图像处理装置及电子装置
CN110782502B (zh) * 2018-07-31 2023-11-03 通用电气公司 基于深度学习的pet散射估计系统和使用感知神经网络模型的方法
US11893668B2 (en) 2021-03-31 2024-02-06 Leica Camera Ag Imaging system and method for generating a final digital image via applying a profile to image information

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4437228B2 (ja) * 2005-11-07 2010-03-24 大学共同利用機関法人情報・システム研究機構 焦点ぼけ構造を用いたイメージング装置及びイメージング方法
EP2929508A1 (en) * 2012-12-07 2015-10-14 Canon Kabushiki Kaisha Image generating apparatus and image generating method

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2023037657A1 (ja) * 2021-09-09 2023-03-16 株式会社島津製作所 データ解析装置およびデータ解析方法
CN114283128A (zh) * 2021-12-16 2022-04-05 中国人民解放军海军潜艇学院 一种基于多参数阈值的海洋中尺度涡边缘检测方法及系统
CN114283128B (zh) * 2021-12-16 2024-04-26 中国人民解放军海军潜艇学院 一种基于多参数阈值的海洋中尺度涡边缘检测方法及系统

Also Published As

Publication number Publication date
US20150279033A1 (en) 2015-10-01

Similar Documents

Publication Publication Date Title
JP2015192238A (ja) 画像データ生成装置および画像データ生成方法
Broxton et al. Wave optics theory and 3-D deconvolution for the light field microscope
US9881373B2 (en) Image generating apparatus and image generating method
JP5705096B2 (ja) 画像処理装置及び画像処理方法
Kandel et al. Label-free tissue scanner for colorectal cancer screening
Stefanoiu et al. Artifact-free deconvolution in light field microscopy
Cannell et al. Image enhancement by deconvolution
JP6305175B2 (ja) 画像処理装置、画像処理方法、画像処理システム
JP5996334B2 (ja) 顕微鏡システム、標本画像生成方法及びプログラム
JP2002531840A (ja) コンピュータを使用した適応画像形成装置及び方法
JP6362062B2 (ja) 画像生成装置および画像生成方法
JP2015108837A (ja) 画像処理装置及び画像処理方法
JP2015057682A (ja) 画像生成装置および画像生成方法
Bueno et al. An automated system for whole microscopic image acquisition and analysis
Guo et al. Deep learning-enabled whole slide imaging (DeepWSI): oil-immersion quality using dry objectives, longer depth of field, higher system throughput, and better functionality
Zhang et al. Virtual staining of defocused autofluorescence images of unlabeled tissue using deep neural networks
Goodwin Quantitative deconvolution microscopy
JP2015191362A (ja) 画像データ生成装置および画像データ生成方法
Jiang et al. Blind deblurring for microscopic pathology images using deep learning networks
Teranikar et al. Correcting anisotropic intensity in light sheet images using dehazing and image morphology
Prigent et al. SPITFIR (e): a supermaneuverable algorithm for fast denoising and deconvolution of 3D fluorescence microscopy images and videos
Conchello et al. Extended depth-of-focus microscopy via constrained deconvolution
Heintzmann et al. Reconstruction of axial tomographic high resolution data from confocal fluorescence microscopy: a method for improving 3D FISH images
WO2019246478A1 (en) Systems and methods for interferometric multifocus microscopy
Garud et al. Volume visualization approach for depth-of-field extension in digital pathology