JP2007014706A - 画像融合処理方法、画像融合処理プログラム、画像融合処理装置 - Google Patents

画像融合処理方法、画像融合処理プログラム、画像融合処理装置 Download PDF

Info

Publication number
JP2007014706A
JP2007014706A JP2005202313A JP2005202313A JP2007014706A JP 2007014706 A JP2007014706 A JP 2007014706A JP 2005202313 A JP2005202313 A JP 2005202313A JP 2005202313 A JP2005202313 A JP 2005202313A JP 2007014706 A JP2007014706 A JP 2007014706A
Authority
JP
Japan
Prior art keywords
image
fusion processing
calculating
image fusion
voxel
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.)
Granted
Application number
JP2005202313A
Other languages
English (en)
Other versions
JP4267598B2 (ja
Inventor
Kazuhiko Matsumoto
和彦 松本
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.)
Ziosoft Inc
Original Assignee
Ziosoft 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 Ziosoft Inc filed Critical Ziosoft Inc
Priority to JP2005202313A priority Critical patent/JP4267598B2/ja
Priority to US11/484,059 priority patent/US7817877B2/en
Publication of JP2007014706A publication Critical patent/JP2007014706A/ja
Application granted granted Critical
Publication of JP4267598B2 publication Critical patent/JP4267598B2/ja
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T15/003D [Three Dimensional] image rendering
    • G06T15/08Volume rendering

Abstract

【課題】 情報量が縮退することなく、同一の観察対象から得た複数の3次元以上の画像を融合することのできる画像融合処理方法、画像融合処理プログラム、画像融合処理装置を提供する。
【解決手段】 CPUは第1ボクセル値D1n に対応する光パラメータを算出し(ステップS2−30)、第2ボクセル値D2n に対応する光パラメータを算出する(ステップS2−35)。次に、CPUは現在位置X1n,X2n毎に合成比率決定処理(ステップS2−40)を行って得た合成比率に基づいて合成したそれぞれの光パラメータを算出し(ステップS2−45)、それら光パラメータを用いて残存光In 、反射光En を更新する(ステップS2−50)。現在位置X1n 又は現在位置X2n が終了位置であるとき(ステップS2−60でYES)、CPUは反射光En をピクセル値とし(ステップS2−65)、フレームを構成する1ピクセルについての融合画像データ生成処理を終了する。
【選択図】 図8

Description

本発明は、画像融合処理方法、画像融合処理プログラム、画像融合処理装置に関する。
従来、人体の臓器等の同一の観察対象をCT(Computerized Tomography )画像撮影装置やPET(Positron Emission Tomography)画像撮影装置等の複数の異なるモダリティで撮影し、得られた画像を融合して診察等に利用している。例えば、異常の有無等の臓器の機能情報を示すPET画像と、解像度が高いCT画像とを融合することで、異常部の大きさ、位置関係を一層正確に把握できる。
このような画像融合において、図21に示すように、2つのボリュームデータ(CT画像101、PET画像103)を別々にレンダリングした後、それぞれの投影画像105,107を一定の比率で融合して融合画像109を得る方法がある。しかし、この方法では、2つの画像を2次元的に合成するので色が不鮮明になるとともに、CT画像で精細に表現されていた表面の凹凸が見えにくくなる。
これに対し、図22に示すように、所定のサンプリング点毎に、ボクセル値を合成し、その合成したボクセル値111をレンダリングして融合画像113を得る画像融合処理方法がある(例えば、非特許文献1)。これにより、前後関係の表現が可能となり、得られる画像が鮮明となっている。
マリア・フェレ(Maria Ferre)、アンナ・ピュイグ(Anna Puig)、ダニー・トスト(Dani Tost)、「多様なボリュームデータの合成方法及びレンダリング技術(A framework for fusion methods and rendering techniques of multimodal volume data )」、ジャーナルオブヴィジュアライゼーションアンドコンピュータアニメーション(Journal of Visualization and Computer Animation)、(米国)、2004年5月、15巻(Volume 15)、第2号(Number 2)、p.63-67
しかしながら、この画像融合処理方法では、レンダリング前にボクセル値を合成してしまうので、不透明度、シェーディング係数、色等のボリュームデータから得られる光パラメータが全て合成後のボクセル値に従属することから、複数のボリュームデータで多次元の情報量を有していたボクセル値の情報が1次元情報となり、情報量が縮退していた。すなわち、CT画像、PET画像のいずれかの持つ光パラメータでレンダリングされることとなるため、融合後の画像の表現が限定されていた。特にCT画像の解像度とPET画像の解像度とが異なるとき、得られる画像は不自然になっていた。
本発明は、上記問題点を解消するためになされたものであって、その目的は、情報量が縮退することなく、同一の観察対象から得た複数の3次元以上の画像を融合することのできる画像融合処理方法、画像融合処理プログラム、画像融合処理装置を提供することにある。
上記問題点を解決するために、請求項1に記載の発明は、1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の
3次元以上の画像データを2次元画面上に融合した融合画像を得る画像融合処理方法において、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定段階と、前記合成比率に基づいて合成光パラメータを算出する算出段階と、前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出段階と、前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出段階とを備えた。
この発明によれば、複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを融合画像を生成するために合成する合成比率を仮想光線上で少なくとも二度決定し、合成比率に基づいて合成光パラメータを算出し、合成光パラメータを基に仮想光線に対する部分反射光をサンプリング位置毎に算出し、部分反射光を積算して融合画像の画素値を算出する。この結果、複数の光パラメータを合成する合成比率を仮想光線上で少なくとも二度決定するため、情報量が縮退することなく、同一の観察対象から得た複数の3次元以上の画像を融合することができる。
請求項2に記載の発明は、請求項1に記載の画像融合処理方法において、前記光パラメータは色情報および不透明度情報であって、前記合成比率決定段階は、前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定段階と、前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定段階とを備え、前記算出段階は、前記色合成比率に基づいて合成色情報を算出する第1算出段階と、前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出段階とを備えた。
この発明によれば、光パラメータは色情報および不透明度情報であって、各ボクセル値に対応付けられた色情報を合成し、仮想光線上で少なくとも二度決定する色合成比率に基づいて合成色情報を算出し、各ボクセル値に対応付けられた不透明度情報を合成し、仮想光線上で少なくとも二度決定する不透明度合成比率に基づいて合成不透明度情報を算出する。この結果、仮想光線上で少なくとも二度、色情報によって観察対象の状態を、不透明度情報によって観察対象の形状をそれぞれ出力することができることから、一層正確に観察対象の状態及び形状を把握可能な融合画像を生成することができる。
請求項3に記載の発明は、請求項2に記載の画像融合処理方法において、前記合成比率決定段階は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定段階をさらに備え、前記算出段階は、前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出段階をさらに備え、前記部分反射光算出段階は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する。
この発明によれば、ボクセル値に対応付けられたシェーディング係数を合成し、前記仮想光線上で少なくとも二度決定する陰影合成比率に基づいて合成シェーディング係数を算出し、合成色情報及び合成不透明度情報及び合成シェーディング係数を基に仮想光線に対する部分反射光を算出する。この結果、シェーディング係数を動的に設定することができるので、観察対象の表面の凹凸を出力することができることから、一層、正確に観察対象の凹凸を把握可能な融合画像を生成することができる。
請求項4に記載の発明は、請求項1乃至3のいずれか1項に記載の画像融合処理方法において、前記合成比率決定段階は、前記複数の3次元以上の画像データのうち少なくとも
1つの3次元以上の画像データに投射した仮想光線上のサンプリング位置のボクセル値と第1基準値とを比較する第1比較段階をさらに備え、前記第1比較段階での比較結果に応じて前記合成比率を決定する。
この発明によれば、複数の3次元以上の画像データのうち少なくとも1つの3次元以上の画像データに投射した仮想光線上のサンプリング位置のボクセル値と第1基準値とを比較し、その比較結果に応じて合成比率を決定する。この結果、例えば、ボクセル値が第1基準値より大きい部分を強調表示することができる。従って、一層、観察しやすい融合画像を生成することができる。
請求項5に記載の発明は、請求項1乃至3のいずれか1項に記載の画像融合処理方法において、前記合成比率決定段階は、前記合成比率を多段階で決定する。この発明によれば、合成比率決定段階は合成比率を多段階で決定するので、一層、各3次元以上の画像データの特徴を生かした融合画像を生成することができる。
請求項6に記載の発明は、請求項1乃至5のいずれか1項に記載の画像融合処理方法において、前記合成比率決定段階は、前記合成比率を前記ボクセル値の勾配に応じて設定する。この発明によれば、合成比率をボクセル値の勾配に応じて設定するので、例えば、ボクセル値の勾配(陰影)と組織情報とを区別することができることから、観察対象の輪郭を出力しながらも、観察対象内に3次元的に含まれている組織情報を出力することができる。従って、観察対象上のみならず観察対象内部の状態も観察することができる。
請求項7に記載の発明は、請求項1乃至6のいずれか1項に記載の画像融合処理方法において、前記ボクセル値のうちの1つと第2基準値とを比較する第2比較段階と、前記第2比較段階での比較結果に応じてマスク領域を生成するマスク領域生成段階とをさらに備えた。この発明によれば、ボクセル値のうちの1つと第2基準値とを比較し、その比較結果に応じてマスク領域を生成する。この結果、例えば、ボクセル値が異常である部分のみマスク領域として画面に出力することができる。従って、所望の箇所のみ観察することができるので、一層、観察しやすい融合画像を得ることができる。
請求項8に記載の発明は、請求項1または2に記載の画像融合処理方法において、同一の複数の3次元以上の画像データのサンプリング位置のボクセル値の前記光パラメータをそれぞれ異なる合成比率で生成した複数の融合画像を用いてアニメーションを作成する。この発明によれば、同一の複数の3次元以上の画像データのサンプリング位置のボクセル値の光パラメータをそれぞれ異なる合成比率で生成した複数の融合画像を用いてアニメーションを作成するので、動的に変更した合成比率に対応した融合画像を動的に観察することができることから、一層、正確に融合画像を観察することができる。
請求項9に記載の発明は、1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の3次元以上の画像データを2次元画面上に融合した融合画像を得る画像融合処理プログラムにおいて、前記1つのコンピュータまたは複数のコンピュータを、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定手段と、前記合成比率に基づいて合成光パラメータを算出する算出手段と、前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出手段と、前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出手段として機能させる。この発明によれば、請求項1に記載の発明と同様
の作用効果が得られる。
請求項10に記載の発明は、請求項9に記載の画像融合処理プログラムにおいて、前記光パラメータは色情報および不透明度情報であって、前記合成比率決定手段は、前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定手段と、前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定手段とを備え、前記算出手段は、前記色合成比率に基づいて合成色情報を算出する第1算出手段と、前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出手段とを備えた。この発明によれば、請求項2に記載の発明と同様の作用効果が得られる。
請求項11に記載の発明は、請求項10に記載の画像融合処理プログラムにおいて、前記合成比率決定手段は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定手段をさらに備え、前記算出手段は、前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出手段をさらに備え、前記部分反射光算出手段は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する。この発明によれば、請求項3に記載の発明と同様の作用効果が得られる。
請求項12に記載の発明は、請求項9乃至11のいずれか1項に記載の画像融合処理プログラムにおいて、前記合成比率決定手段は、前記複数の3次元以上の画像データのうち少なくとも1つの3次元以上の画像データに投射した仮想光線上のサンプリング位置のボクセル値と第1基準値とを比較する第1比較手段をさらに備え、前記第1比較手段での比較結果に応じて前記合成比率を決定する。この発明によれば、請求項4に記載の発明と同様の作用効果が得られる。
請求項13に記載の発明は、請求項9乃至11のいずれか1項に記載の画像融合処理プログラムにおいて、前記合成比率決定手段は、前記合成比率を多段階で決定する。この発明によれば、請求項5に記載の発明と同様の作用効果が得られる。
請求項14に記載の発明は、請求項9乃至13のいずれか1項に記載の画像融合処理プログラムにおいて、前記合成比率決定手段は、前記合成比率を前記ボクセル値の勾配に応じて設定する。この発明によれば、請求項6に記載の発明と同様の作用効果が得られる。
請求項15に記載の発明は、請求項9乃至14のいずれか1項に記載の画像融合処理プログラムにおいて、前記ボクセル値のうちの1つと第2基準値とを比較する第2比較手段と、前記第2比較手段での比較結果に応じてマスク領域を生成するマスク領域生成手段とをさらに備えた。この発明によれば、請求項7に記載の発明と同様の作用効果が得られる。
請求項16に記載の発明は、請求項9または10に記載の画像融合処理プログラムにおいて、同一の複数の3次元以上の画像データのサンプリング位置のボクセル値の前記光パラメータをそれぞれ異なる合成比率で生成した複数の融合画像を用いてアニメーションを作成する。この発明によれば、請求項8に記載の発明と同様の作用効果が得られる。
請求項17に記載の発明は、1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の3次元以上の画像データを2
次元画面上に融合した融合画像を得る画像融合処理装置であって、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定手段と、前記合成比率に基づいて合成光パラメータを算出する算出手段と、前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出手段と、前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出手段とを備えた。この発明によれば、請求項1,9に記載の発明と同様の作用効果が得られる。
請求項18に記載の発明は、請求項17に記載の画像融合処理装置において、前記光パラメータは色情報および不透明度情報であって、前記合成比率決定手段は、前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定手段と、前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定手段とを備え、前記算出手段は、前記色合成比率に基づいて合成色情報を算出する第1算出手段と、前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出手段とを備えた。この発明によれば、請求項2,10に記載の発明と同様の作用効果が得られる。
請求項19に記載の発明は、請求項18に記載の画像融合処理装置において、前記合成比率決定手段は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定手段をさらに備え、前記算出手段は、前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出手段をさらに備え、前記部分反射光算出手段は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する。この発明によれば、請求項3,11に記載の発明と同様の作用効果が得られる。
請求項20に記載の発明は、請求項17乃至19のいずれか1項に記載の画像融合処理装置において、前記合成比率決定手段、前記算出手段、前記部分反射光算出手段、前記画素値算出手段の一部または全部は、GPUである。この発明によれば、合成比率決定手段、算出手段、部分反射光算出手段、画素値算出手段の一部または全部は、GPUであるので、合成比率の決定、光パラメータの算出、部分反射光の算出、画素値の算出に要する時間を短縮することができる。
請求項21に記載の発明は、請求項17乃至20のいずれか1項に記載の画像融合処理装置において、前記画素値算出手段が算出した画素値を出力部に出力可能な形式に変換する後処理は、GPUで行う。この発明によれば、画素値算出手段が算出した画素値を出力部に出力可能な形式に変換する後処理は、GPUで行うので、後処理に要する時間を短縮することができる。
本発明によれば、情報量が縮退することなく、同一の観察対象から得た複数の3次元以上の画像を融合することができる。
(第1実施形態)
以下、本発明を具体化した第1実施形態を図1〜図10に従って説明する。
図1に示すように、画像表示装置1は、データベース2から例えば、CT画像撮影装置により撮影されたCT画像データ及びPET画像撮影装置により撮影されたPET画像データを読み取って、医療診断用の画像を生成し画面に表示する。本実施形態では、CT画
像データ及びPET画像データを例に説明するが、これに限定されない。すなわち、使用される画像データは、CT画像データ、PET画像データに限らず、MRI(Magnetic Resonance Imaging)、MRA(Magnetic Resonance Angiography)等の医用画像処理装置より得られる画像データ及びそれらを組み合わせたり加工したりしたものである。
画像表示装置1は、計算機(コンピュータ、ワークステーション、パーソナルコンピュータ)3と、モニタ4と、キーボード5及びマウス6などの入力装置とを備えている。計算機3はデータベース2と接続されている。
図2は、画像表示装置1の概略構成を示す。計算機(コンピュータ)3にはCPU(中央処理装置)7、ハードディスク等からなるメモリ8が備えられており、メモリ8には、画像融合処理を実行するためのプログラム(アプリケーションソフト)9が記憶されている。メモリ8はデータベース2又はハードディスクから読み取ったCT画像データから得られた第1ボクセルデータVD1及び、PET画像データから得られた第2ボクセルデータVD2を一時記憶するメモリ部8aを備えている。また、メモリ部8aは、投影開始点O1,O2、サンプリング間隔S1,S2、サンプリング位置としての現在位置X1n ,X2n 、終了位置を一時記憶するが、詳細は後述する。
メモリ8は、第1及び第2ボクセルデータVD1,VD2に属するボクセル値からそれぞれ得られる情報であって、観察対象の形状、組織情報、表面の凹凸等を表す光パラメータPを記憶する光パラメータ記憶部MPを備えている。また、メモリ8は、第1及び第2ボクセルデータVD1,VD2のいずれに属するボクセル値の光パラメータPを用いるかを判断するための第1基準値を記憶する基準値記憶部SVを備えている。そして、メモリ8は、第1及び第2ボクセルデータVD1,VD2に属する各ボクセルのボクセル値の光パラメータPを合成する比率である合成比率を記憶する合成比率記憶部CVを備えている。さらに、メモリ8は、第1及び第2ボクセルデータVD1,VD2にそれぞれ照射された仮想光線が第1及び第2ボクセルデータVD1,VD2に属する各ボクセルにてそれぞれ反射する反射光を記憶する反射光記憶部RVを備えている。
そして、CPU7は、プログラム9を実行することにより、データベース2から取得したCT画像データ、PET画像データから得られた第1及び第2ボクセルデータVD1,VD2を用いて複数のモダリティで撮影して得られた画像を融合する画像融合処理を実行する。すなわち、本実施形態では、CPU7(計算機3)が、画像融合処理(合成比率決定段階、算出段階、部分反射光算出段階、画素値算出段階、色合成比率決定段階、不透明度合成比率決定段階、第1算出段階、第2算出段階、陰影合成比率決定段階、第3算出段階、第1比較段階、第2比較段階等)の画像融合処理プログラムを実行する。これにより、計算機3は、合成比率決定手段、算出手段、部分反射光算出手段、画素値算出手段、色合成比率決定手段、不透明度合成比率決定手段、第1算出手段、第2算出手段、陰影合成比率決定手段、第3算出手段、第1比較手段、第2比較手段等として機能する。そして、モニタ4(画面4a)には、図10に示すように、画像融合処理を実行した後の画像である融合画像P1が表示されるようになっている。
ここで、3次元以上の画像データとしての第1ボクセルデータVD1、第2ボクセルデータVD2とは、図3に示すように、ボクセル(Voxel )の集合であり、3次元の格子点にボクセル値として濃度値が割り当てられている。本実施形態では、例えば、第1ボクセルデータVD1と第2ボクセルデータVD2とは、同じ患者の同じ部位(肺)を異なるモダリティで撮影したものであるため、第2ボクセルデータVD2は、第1ボクセルデータVD1と符号を同一にし、その詳細な説明を省略する。そして、本実施形態では、CT画像データの画素値、すなわちCT値をそのまま第1ボクセルデータVD1の濃度値としている。また、同様に、PET画像データの画素値、すなわちPET値をそのまま第2ボク
セルデータVD2の濃度値としている。
CT画像データは、患者等の人体を断層撮影したもので、1枚については骨、血管、臓器等の観察対象の2次元断層画像であるが、多数の隣接するスライス(断層)について得られていることから、これら全体では3次元の画像データと言える。従って、以下、CT画像データは、複数のスライスを含んだ3次元の画像データを意味する。
また、CT画像データは、被写体としての組織(骨、血管、臓器等)毎に異なるCT値を持っている。CT値は、水を基準として表現した組織のX線減弱係数であり、CT値により組織や病変の種類等が判断できるようになっている。また、CT画像データには、CT撮影装置によりCTスキャンされた人体の断層画像(スライス画像)の座標データもすべてあり、視線方向(奥行き方向)における異なる組織間の位置関係は、座標データから判別できるようになっている。すなわち、第1ボクセルデータVD1は、CT値(以下、ボクセル値としての第1ボクセル値D1n という。)及び座標データを備えている。
PET画像データは、ポジトロンを放出する放射性同位元素で標識された放射性薬剤を投与された患者の人体を断層撮影したもので、1枚については臓器等の2次元断層画像であるが、多数の隣接するスライス(断層)について得られていることから、これら全体では3次元の画像データと言える。従って、以下、PET画像データは、複数のスライスを含んだ3次元の画像データを意味する。
また、PET画像データは、臓器の機能情報(異常の有無)毎に異なるPET値を持っている。PET値は放射性薬剤の集積度であって、PET値により臓器局所の血流、代謝等の生理的、生化学的機能等が判断できるようになっている。また、PET画像データには、PET撮影装置によりPETスキャンされた人体の断層画像(スライス画像)の座標データもすべてあり、視線方向(奥行き方向)における異なる組織間の位置関係は、座標データから判別できるようになっている。すなわち、第2ボクセルデータVD2は、PET値(以下、ボクセル値としての第2ボクセル値D2n という。)及び座標データを備えている。
融合画像P1とは、第1ボクセルデータVD1と第2ボクセルデータVD2がそれぞれ持つ光パラメータPを予め定めたサンプリング間隔毎に決定される合成比率で合成したときの反射光を2次元平面上に投影した画像であって、本実施形態では、ボリュームレンダリング処理により生成される。
まず、ボリュームレンダリング処理について説明する。ボリュームレンダリングには一般的にレイキャスティング法が用いられる。レイキャスティング法とは、図3に示すように、観察する側(フレームFR側)から光の経路を考えるもので、フレームFR側のピクセルPXから光線(仮想光線R)を飛ばし、一定距離を進むごとにその位置での反射光を計算する(図3では「…、V1、V2、V3、…」の符号が各到達位置のボクセルに対応している)。
視線方向からボクセルデータに1本の仮想光線Rを照射すると、仮想光線Rは第1ボクセルデータVD1に当たって一部反射しながら第1ボクセルデータVD1を構成するボクセルを順次透過して進むことになる。そして、ボクセル毎に光の吸収・反射の離散的な計算をし、その反射光を演算することでフレームFRに投影される画像のピクセル値(画素値)を求めて2次元画像を生成する。
図4は、レイキャスティング法の計算方法を説明するもので、単一画像のボリュームレンダリングに対応した処理である。図4において、投影位置Xn(仮想光線到達位置の現
在位置)が格子上になかったとすると、まず、その投影位置Xn(仮想光線到達位置の現在位置)の周りのボクセルのボクセル値から補間処理を行ってその位置でのボクセル値Dnを計算する。次に、光に対する特性パラメータ(以下、光パラメータPという)を決定する。
光パラメータPとは、不透明度情報としての不透明度(オパシティ(opacity )値)αn およびシェーディング係数βn 、色情報としての色γn であって、それぞれ独立した光学的特性を表す情報である。ここで、不透明度αn は、ボクセルデータの形状を表すパラメータであって、0≦αn ≦1を満たす数値で表され、値(1−αn )は透明度(transparency)を示す。不透明度αn =1は不透明、αn =0は透明、0<αn <1は半透明にそれぞれ対応する。なお、各ボクセル値に対して不透明度αn との関係付けが予めなされており、その関係付け情報に基づきボクセル値から不透明度αn を得ている。例えば、骨のボリュームレンダリング画像を得たい場合、骨に対応したボクセル値には不透明度「1」を関係付け、他のボクセル値には不透明度「0」を関係付けることで、骨を表示することができる。このようにある値に対応する別の値を柔軟に導く関数としてLUT(Look Up
Table)関数がある。シェーディング係数βn は、ボクセルデータの表面の凹凸(陰影)
を表すパラメータである。そして、色γn は、ボクセルデータの組織情報、すなわち、ボクセルデータが骨、血液、内臓、あるいは腫瘍であるかの情報を表している。
図5に示すように、初期入射光(光線)I1は、各ボクセルを順次透過してゆくとともに各ボクセルで一部反射および吸収されることによりその残存光(透過光)は徐々に減衰する。なお、図5は、レイキャスティング法の計算方法を説明するもので、図3における1本の仮想光線に対応した処理である。そして、各ボクセルにおける部分反射光Fn (n=1,2,…)の積算値(積算反射光)が、フレームFR側におけるピクセルPXの輝度に相当する。ここで、減衰光An (n=1,2,…)は、n番目のボクセルの入射光In を用いて、式 An =αn ×In で表されるため、部分反射光Fn は、式 Fn =βn ×γn ×An =βn ×γn ×αn ×In で表される。
よって残存光In 、積算反射光である反射光En は、次式により表される。
In =In-1 −αn-1 ×In-1 =In-1 −An-1
En =En-1 +αn ×βn ×γn ×In =En-1 +Fn
そして、図4に示すように、この反射光Enがサンプリング位置毎に加算され、ピクセル値として画像が生成される。
図6は、本実施形態における画像融合処理方法を説明するものである。なお、本実施形態では、第1ボクセルデータVD1において投影を開始する投影開始点O1が、ユーザのマウス6等の入力装置を介した操作によって入力され、メモリ部8aに予め記憶されているとする。また、第2ボクセルデータVD2において投影開始点O1と同じ座標に相当する投影開始点O2が、CPU7によって算出され、メモリ部8aに予め記憶されているとする。同様に、第1ボクセルデータVD1をサンプリングする間隔であるサンプリング間隔S1と第2ボクセルデータVD2をサンプリングする間隔であるサンプリング間隔S2とがマウス6によって入力され、メモリ部8aに予め記憶されているとする。また、投影位置(仮想光線到達位置の現在位置)X1n は、第1ボクセルデータVD1において、投影開始点O1にサンプリング間隔S1を積算することによって算出され、メモリ部8aに記憶される。同様に、投影位置(仮想光線到達位置の現在位置)X2n は、第2ボクセルデータVD2において、投影開始点O2にサンプリング間隔S2を積算することによって算出され、メモリ部8aに記憶される。
なお、本実施形態では、第1ボクセルデータVD1、第2ボクセルデータVD2がそれぞれ光パラメータPを備えている。すなわち、第1ボクセルデータVD1の光パラメータ
Pは、CT不透明度α1n 、CTシェーディング係数β1n 、CT色γ1n を備えており、第2ボクセルデータVD2の光パラメータPは、PET不透明度α2n 、PETシェーディング係数β2n 、PET色γ2n を備えている。さらに、メモリ部8aには終了位置が予め記憶されている。そして、CPU7は、第1ボクセルデータVD1の光パラメータPと第2ボクセルデータVD2の光パラメータPとをそれぞれ現在位置X1n ,X2n 毎に定めた合成比率CRで投影開始点O1,O2から終了位置まで合成する。
詳述すると、合成比率CRは、不透明度合成比率としての第1合成比率a、陰影合成比率としての第2合成比率b、色合成比率としての第3合成比率cを備えている。第1合成比率aは、CT不透明度α1n とPET不透明度α2n とを合成する比率(合成比率)を表しており、第2合成比率bは、CTシェーディング係数β1n とPETシェーディング係数β2n との合成比率を表している。また、第3合成比率cは、CT色γ1n とPET色γ2n との合成比率を表している。そして、合成後の合成光パラメータ(合成不透明度情報としての合成不透明度α3n 、合成シェーディング係数β3n 、合成色情報としての合成色γ3n )は、この合成比率CR(第1合成比率a、第2合成比率b、第3合成比率c)を用いて次式により表される。
α3n =aα1n +(1−a)α2n
β3n =bβ1n +(1−b)β2n
γ3n =cγ1n +(1−c)γ2n
そして、本実施形態では、図6に示すように、第1合成比率a、第2合成比率b、第3合成比率cのそれぞれをPET画像データDT2の投影位置X2n (仮想光線到達位置の現在位置)の第2ボクセル値D2n に応じて仮想光線到達位置の現在位置毎に決定する。詳述すると、第2ボクセル値D2n が第1基準値T1以下であったとき、第1合成比率a=「1」、すなわち合成不透明度α3n はCT画像データDT1の投影位置(仮想光線到達位置の現在位置)X1n におけるCT不透明度α1n とする。また、第3合成比率c=「1」、すなわち合成色γ3n はCT画像データDT1の現在位置X1n におけるCT色γ1n とする。一方、第2ボクセル値D2n が第1基準値T1以上であったとき、第1合成比率aは同様に「1」とし、第3合成比率c=「0」、すなわち合成色γ3n はPET画像データDT2の投影位置(仮想光線到達位置の現在位置)X2n におけるPET色γ2n とする。なお、本実施形態では、第2ボクセル値D2n の値によらず、第2合成比率b=「1」、すなわち合成シェーディング係数β3n はCT画像データDT1の現在位置X1n におけるCTシェーディング係数β1n とする。これは、PET画像データDT2は異常部分の有無を表現するのみであるが、CT画像データDT1はボクセルデータの表面の凹凸の精細な表現が可能なためである。
従って、第2ボクセル値D2n (PET値)が正常な部分は第1ボクセル値D1n (CT値)に基づいて通常の画像を表示し、第2ボクセル値D2n が異常な部分は第2ボクセル値D2n に従属した合成色γ3n を表示することにより、異常部分のみを強調表示した融合画像データDT3(融合画像P1)を得ることができる。
また、図2に示すように、計算機(コンピュータ)3は、GPU(Graphics Processing Unit)10を備えている。GPU10は、主に3次元の高度なグラフィックス機能をサポートしたグラフィックス用のコントローラ・チップであり、ユーザから与えられたプログラム等に基づいて2次元や3次元のグラフィックス描画を高速化する機能を持つ。本実施形態では、GPU10により融合画像データDT3に対して後処理が実行され、融合画像P1が生成される。これにより、融合画像P1の表示に要する時間を短縮することができる。
後処理は、算出した融合画像データDT3をモニタ4等の出力装置に表示するために色
、コントラスト及び明るさ等を補正する処理である。詳しくは、多くの医用画像装置の出力(CT画像、PET画像等)は12bit階調データであるため、画像融合処理で生成された融合画像P1も12bit階調データであるが、コンピュータ3等のモニタ4はRGB各色を8bitで表現する画像を表示することが多い。そのため、WL変換(Window
Level Transformation )やLUT変換(Color Look-Up Table Transformation)等を行い、色、コントラスト及び明るさ等をモニタ4に表示できるように変換する。また、画面の大きさ等に合わせてアフィン変換(Affine Transformation )を行い、それぞれのモニタ4に表示できるように変換する。
次に、このように構成された画像融合処理の作用について説明する。
図7は、画像融合処理のフローチャートを示している。
まず、CPU7は、融合画像データ生成処理を実行する(ステップS1−10)。融合画像データ生成処理は図8に示す処理が実行される。まず、CPU7は、第1ボクセルデータVD1における投影開始点O1及びサンプリング間隔S1、第2ボクセルデータVD2における投影開始点O2及びサンプリング間隔S2を設定する(ステップS2−10)。このとき、CPU7は、投影開始点O1,O2及びサンプリング間隔S1,S2をメモリ部8aに記憶する。次に、CPU7は、反射光En 、残存光In 、現在位置X1n ,X2n の初期化をする(ステップS2−15)。すなわち、反射光En =0、残存光In =1、第1ボクセルデータVD1における現在位置X1n =0、第2ボクセルデータVD2における現在位置X2n =0とする。そして、CPU7は、現在位置X1n における第1ボクセル値D1n 、勾配g1n を算出する(ステップS2−20)。本実施形態では、現在位置X1n が格子上になかったとすると、現在位置X1n の周りのボクセルのボクセル値から補間処理を行い、その位置での第1ボクセル値D1n を計算する。また、勾配g1n を公知の方法で算出する。同様に、CPU7は、現在位置X2n における第2ボクセル値D2n 、勾配g2n を算出する(ステップS2−25)。
次に、CPU7は、第1ボクセル値D1n に対応するCT不透明度α1n とCT色γ1n 、勾配g1n に対応するCTシェーディング係数β1n を算出する(ステップS2−30)。すなわち、光パラメータ記憶部MPにルックアップテーブル(Look Up Table )として記憶されている現在位置X1n に対応するCT不透明度α1n 、CT色γ1n を読み出す。また、勾配g1n を用いて公知の方法でCTシェーディング係数β1n を算出する。同様に、CPU7は、光パラメータ記憶部MPにルックアップテーブル(Look Up Table )として記憶されている第2ボクセル値D2n に対応するPET不透明度α2n とPET色γ2n を読み出し、また、勾配g2n に対応するPETシェーディング係数β2n を算出する(ステップS2−35)。
そして、CPU7は、合成比率決定処理を実行する(ステップS2−40)。合成比率決定処理は、図9に示すような、第1合成比率a、第2合成比率b、第3合成比率cのそれぞれを現在位置毎に決定する処理が実行される。まず、CPU7は、第1合成比率a=「1」とする(ステップS3−10)。これにより、合成不透明度α3n は、CT画像データDT1のCT不透明度α1n を用いる。そして、CPU7は、第2合成比率b=「1」とする(ステップS3−15)。これにより、合成シェーディング係数β3n は、CT画像データDT1のCTシェーディング係数β1n を用いる。次に、CPU7は、第2ボクセル値D2n が第1基準値T1より大きいか否かを判断する(ステップS3−20)。CPU7は、第1基準値T1を基準値記憶部SVから読み出して、第2ボクセル値D2n と比較する。そして、第2ボクセル値D2n が第1基準値T1より大きいとき(ステップS3−20でYES)、第3合成比率c=「0」とする(ステップS3−25)。すなわち、PET値が異常であるため、合成色γ3n は、PET画像データDT2のPET色γ2n を用いる。一方、第2ボクセル値D2n が第1基準値T1より小さいとき(ステップS3−20でNO)、第3合成比率c=「1」とする(ステップS3−30)。すなわち
、PET値が正常であるため、合成色γ3n は、CT画像データDT1のCT色γ1n を用いる。そして、CPU7は、それら合成比率CRを合成比率記憶部CVに記憶し、この合成比率決定処理を終了する。
合成比率決定処理が終了すると、それら第1合成比率a、第2合成比率b、第3合成比率cに基づいて合成不透明度α3n 、合成シェーディング係数β3n 、合成色γ3n を算出する(図8のステップS2−45)。すなわち、CPU7は、合成比率CRを合成比率記憶部CVから、CT不透明度α1n 、CTシェーディング係数β1n 、CT色γ1n (PET色γ2n )を光パラメータ記憶部MPからそれぞれ読み出し、合成不透明度α3n 、合成シェーディング係数β3n 、合成色γ3n を算出する。そして、CPU7は、それら光パラメータPを用いて残存光In 、反射光En を更新する(ステップS2−50)。詳述すると、CPU7は、まず、減衰光An (=α3n An-1 )、部分反射光Fn (=β3γ3An )を算出する。そして、それら算出した減衰光An 、部分反射光Fn を用いて残存光In 、反射光En を更新する。なお、このステップS2―10〜S2−50の処理を、以下、融合処理という。
次に、CPU7は、現在位置X1n ,X2n をそれぞれ進行させる(ステップS2−55)。すなわち、CPU7は、サンプリング間隔S1,S2をそれぞれメモリ部8aから読み出し、現在位置X1n をサンプリング間隔S1、現在位置X2n をサンプリング間隔S2、それぞれ進行させる。そして、CPU7は、現在位置X1n または現在位置X2n が予め定めた終了位置であるか否かを判断する(ステップS2−60)。すなわち、CPU7は、終了位置をメモリ部8aから読み出し、現在位置X1n ,X2n とそれぞれ比較する。そして、現在位置X1n または現在位置X2n が終了位置でないとき(ステップS2−60でNO)、CPU7は、再び融合処理を実行する。一方、現在位置X1n または現在位置X2n が終了位置であるとき(ステップS2−60でYES)、CPU7は、反射光En をピクセル値とし(ステップS2−65)、フレームFRを構成する1ピクセルについての融合画像データ生成処理を終了する。
フレームFR内の全てのピクセルについて融合画像データ生成処理を終了して融合画像データDT3を得ると、GPU10は、その融合画像データDT3について後処理を行い、融合画像P1を生成する(図7のステップS1−15)。後処理が終了すると、図10に示すように、融合画像P1をモニタ4の画面4aに出力する(図7のステップS1−20)。このとき、融合画像P1は、その光パラメータP(合成不透明度α3n 、合成シェーディング係数β3n 、合成色γ3n )のそれぞれについて、それぞれの現在位置X1n,X2nのPET値が正常である正常部分ではCT画像データDT1のCT不透明度α1n
とCT色γ1n を用いる。そして、PET値が異常である異常部分U1ではPET画像
データDT2のPET色γ2n を用いている。これにより、正常部分ではCT画像データDT1に基づいて通常の画像を表示しながらも、異常部分U1のみをPET画像データDT2に基づいて強調表示した融合画像P1を得ることができる。
すなわち、CT画像データDT1に基づくことによって観察対象の形状を把握し、PET画像データDT2に基づくことによって観察対象の状態を把握することができる。従って、観察対象の形状の観察を妨げることなく状態を把握しながらも、異常部分U1の位置を把握することができるため、一層正確に画像診断をすることができる。
上記第1実施形態によれば、以下のような効果を得ることができる。
(1) 第1実施形態によれば、第1合成比率a、第2合成比率b、第3合成比率cはそれぞれサンプリング位置毎に計算される。この結果、それぞれの光パラメータP(合成不透明度α3n 、合成シェーディング係数β3n 、合成色γ3n )をそれぞれ現在位置X1n ,X2n 毎に第2ボクセル値D2n に応じて決定した合成比率CRで合成することに
よって、融合画像P1を生成することができる。従って、現在位置X1n ,X2n 毎にそれぞれ独立した光学的特性を表す情報である光パラメータPを合成することができるので、情報量を縮退することなく、融合画像P1を生成することができる。
(2) 第1実施形態によれば、第1合成比率a、第2合成比率b、第3合成比率cはそれぞれサンプリング位置毎に計算される。この結果、それぞれの臓器や臓器の状態に適した合成比率CRが動的に計算され、異常部分U1の形状や正常部分の観察を妨げることがなく、異常部分U1を強調表示することができる。従って、一層正確に観察することができる。
(3) 第1実施形態によれば、CTスキャンによって得られたCT画像データDT1に従属している光パラメータPとPETスキャンによって得られたPET画像データDT2に従属している光パラメータPをそれぞれ現在位置X1n ,X2n 毎に第2ボクセル値D2n に応じて決定した合成比率CRで合成することによって、融合画像P1を生成した。この結果、現在位置X1n ,X2n 毎にそれぞれ独立した光学的特性を表す情報である光パラメータPを合成することができるので、情報量を縮退することなく、それぞれの画像データの持つ特性を生かした融合画像P1を生成することができる。
(4) 第1実施形態によれば、合成比率CRは、現在位置X1n ,X2n 毎に決定したので、現在位置X1n ,X2n 毎に第2ボクセル値D2n に応じてCT画像データDT1、PET画像データDT2、それぞれの画像データの持つ特性を生かした画像を生成することができる。従って、一層正確に画像診断をすることができる。
(5) 第1実施形態によれば、CT画像データDT1とPET画像データDT2とを融合した1つの画像である融合画像P1をモニタ4の画面4a上に表示する。従って、例えば、CT画像とPET画像とを並べて観察するときと比較して、簡単でありながらも正確に観察することができる。
(6) 第1実施形態によれば、合成比率決定処理は複雑な計算を行わないため、リアルタイムに複数種類の画像データを融合して融合画像P1を生成し、モニタ4の画面4aに表示することができる。
(7) 第1実施形態によれば、融合画像P1は、CT画像データDT1とPET画像データDT2にそれぞれ従属している光パラメータPを合成するので、例えば、CT画像データDT1とPET画像データDT2とを2次元的に合成するときと比較して、視線方向を変えた観察が可能となる。従って、一層正確に観察することができる。
(8) 第1実施形態によれば、融合画像P1は、CT画像データDT1とPET画像データDT2にそれぞれ従属している光パラメータPを合成するので、観察者が前後関係を把握することができる。従って、一層、正確に観察することができる。
(9) 第1実施形態によれば、第2ボクセル値D2n が第1基準値T1より大きい、すなわち異常部分U1であるとき、合成色γ3n はPET画像データDT2に従属しているPET色γ2n を用いたため、異常部分U1のみ強調表示することができる。
(10) 第1実施形態によれば、光パラメータPのうち、合成不透明度α3n 、合成シェーディング係数β3n はCT画像データDT1のCT不透明度α1n 、CTシェーディング係数β1n を用いたため、異常部分U1の形状や正常部分の観察を妨げることがなく、異常部分U1を強調表示することができる。従って、一層正確に観察することができる。
(第2実施形態)
次に、本発明を具体化した第2実施形態を図10〜図13に従って説明する。なお、第2実施形態は、合成比率決定処理において多段階で合成比率CRを決定し、さらに、例えば骨に転移した腫瘍のように、組織内部に3次元的に含まれた異常部分を表示することに特徴がある。そのため、第1実施形態と同様の部分については同一の符号を付し、その詳細な説明を省略する。
本実施形態では、図11に示すように、合成比率CRをCT画像データDT1の勾配g1n とPET画像データDT2の第2ボクセル値D2n とに応じて決定する。詳述すると、図10に示すように、融合画像上に表示したいPET値(第2ボクセル値D2n )、すなわち異常部分U2が観察対象の組織内部にあったとしても、その組織のCT値(第1ボクセル値D1n )の高い部分や観察者の視線方向に対して影になる部分等は黒く表示されてしまうため、観察者がその異常部分U2を視認することができない。そこで、CT値が高く且つPET値が高いとき、合成不透明度α3n 及び合成色γ3n をPET値に応じた値にすることによって、PET画像データDT2から得られる情報を優先的に表示する。
本実施形態では、例えば、骨への腫瘍(異常部分U2)の転移を強調表示したいとする。そして、CT画像データDT1は骨を表示するデータであって、PET画像データDT2は腫瘍を表示するデータであるとする。このとき、骨及び腫瘍の両方を表示可能な合成比率CRを決定するため、第1中間変数d、第2中間変数e、第3中間変数fを算出してメモリ部8aに記憶する。第1中間変数dは、骨の領域とその周辺を表示するための変数である。そして、CT値(第1ボクセル値D1n )は、骨の領域では1000以上の値を示すことから、第1中間変数dは次式により表される。
第1中間変数d=D1n /1000
このようにすることによって、骨の領域及びその周辺を多少含んだ領域を柔軟に表現することができる。なお、第1中間変数dの最大値を「1」、最小値を「0」とすることによって、範囲外の値によって不適切な結果が得られることを低減することができる。
第2中間変数eは、異常部分U2を表示するための変数であって、異常部分U2ではPET値(第2ボクセル値D2n )が100以上の値であることから、次式により表される。
第2中間変数e=D2n /100
このようにすることによって、異常の度合いをより柔軟に表現することができる。なお、第2中間変数eの最大値を「1」、最小値を「0」とすることによって、範囲外の値によって不適切な結果が得られることを低減することができる。
第3中間変数fは、骨の領域と腫瘍の領域の両方の条件を兼ねる箇所を特にPET画像データDT2の情報(腫瘍の情報)を用いて表示するための変数であって、次式により表される。
第3中間変数f=1−d×e
次に、このように構成された合成比率決定処理の作用について説明する。
図12は、ピクセル毎の合成比率決定処理のフローチャートを示している。
まず、CPU7は、第1中間変数dを算出する(ステップS4−10)。すなわち、第1ボクセル値D1n を用いて第1中間変数dを算出し、メモリ部8aに記憶する。次に、CPU7は、第2中間変数eを算出する(ステップS4−15)。すなわち、第2ボクセ
ル値D2n を用いて第2中間変数eを算出し、メモリ部8aに記憶する。次に、CPU7は、第3中間変数fを算出する(ステップS4−20)。すなわち、メモリ部8aから、第1中間変数d及び第2中間変数eを読み出して第3中間変数fを算出し、メモリ部8aに記憶する。そして、CPU7は、第1合成比率aにその第3中間変数fを代入し、合成比率記憶部CVに記憶する(ステップS4−25)。また、CPU7は、第2合成比率bに「1」を代入し、合成比率記憶部CVに記憶する(ステップS4−30)。さらに、CPU7は、第3合成比率cにその第3中間変数fを代入し、合成比率記憶部CVに記憶し(ステップS4−35)、この合成比率決定処理を終了する。その後、CPU7は、第1実施形態と同様に画像融合処理を実行し、図13に示すように、骨(第1ボクセル値D1n )に転移した腫瘍(第2ボクセル値D2n )のような異常部分U2の観察が可能な融合画像P2を得る。
これにより、PET値(第2ボクセル値D2n )が高く、融合画像P2上に表示したい異常部分U2がCT値(第1ボクセル値D1n )の高い部分の内部に3次元的に含まれているときであっても、融合画像P2上に表示することができる。また、第1中間変数d、第2中間変数e、第3中間変数fによって、現在位置X1n,X2n毎にそれぞれの合成比率CR(第1合成比率a、第2合成比率b、第3合成比率c)を多段階で設定することができるため、一層精細な融合画像P2を生成することができる。
上記第2実施形態によれば、第1実施形態に記載の(1)〜(8)の効果に加えて以下のような効果を得ることができる。
(11) 第2実施形態によれば、融合画像P2上に黒く表示され得る箇所が、CT値(第1ボクセル値D1n )によるものか、第1ボクセル値D1n の勾配g1n (影)によるものかを判断し、CT値によるものであったとき、CT値よりもPET値(第2ボクセル値D2n )を優先して表示した。これにより、例えば、骨に転移した腫瘍等の組織内部の異常部分U2も観察することができる。
(12) 第2実施形態によれば、第1中間変数d、第2中間変数e、第3中間変数fによって多段階で合成比率CRを設定し、その合成比率CRを用いて融合画像P2を生成したので、一層、CT画像データDT1とPET画像データDT2の特性を生かした融合画像P2を生成することができる。
(13) 第2実施形態によれば、第1中間変数dは第1ボクセル値D1n と所定のCT値(「1000」)によって算出し、第2中間変数eは第2ボクセル値D2n と所定のPET値(「100」)によって算出した。この結果、所定のCT値、所定のPET値を変更することによって、合成比率CRの算出基準を変更することができるので、多様な観察対象に対応した融合画像P1を生成することができる。
(14) 第2実施形態によれば、融合画像P2上に黒く表示され得る原因が、CT値(第1ボクセル値D1n )によるものか、勾配g1n (影)によるものかを切り分けた。この結果、観察者の視線の変化に応じて表示される影の部分を生かしながらも、例えば、CT値(第1ボクセル値D1n )の高い部分の内部に3次元的に含まれたPET値(第2ボクセル値D2n )を表示することができる。
(第3実施形態)
次に、本発明を具体化した第3実施形態を図14〜図17に従って説明する。なお、第3実施形態は、マスク領域を生成し、観察対象のうちそのマスク領域のみ融合画像P3として表示することに特徴がある。そのため、第1及び第2実施形態と同様の部分については同一の符号を付し、その詳細な説明を省略する。
本実施形態では、図14及び図15に示すように、第2ボクセル値D2n に応じたマスク領域M1を生成する。すなわち、第2ボクセル値D2n (PET値)が第2基準値T2を超えた箇所は異常部分であり、マスク領域M1内であるとする。一方、第2ボクセル値D2n (PET値)が第2基準値T2以下の箇所は正常部分であり、マスク領域M1外であるとする。そして、このマスク領域M1内のみ、第1合成比率a、第2合成比率b、第3合成比率cのそれぞれを現在位置毎に決定する合成比率決定処理を行い、融合画像P3を生成する。なお、本実施形態では、第2基準値T2は、基準値記憶部SVに予め記憶されている。また、CPU7(計算機3)が、さらに、画像融合処理(マスク領域生成段階等)の画像融合処理プログラムを実行する。これにより、計算機3は、さらにマスク領域生成手段等として機能する。
次に、このように構成された融合画像データ生成処理の作用について説明する。
図16は、ピクセル毎の融合画像データ生成処理のフローチャートを示している。
まず、CPU7は、投影開始点O1及びサンプリング間隔S1、投影開始点O2及びサンプリング間隔S2を設定し(ステップS2−10)、反射光En 、残存光In 、現在位置X1n ,X2n の初期化をする(ステップS2−15)。次に、マスク生成処理を実行する(ステップS5−10)。マスク生成処理は、図17に示すような処理が実行される。まず、CPU7は、現在位置X2n における第2ボクセル値D2n を算出し(ステップS6−10)、その第2ボクセル値D2n が第2基準値T2を超えているか否かを判断する(ステップS6−15)。すなわち、CPU7は、基準値記憶部SVから第2基準値T2を読み出し、第2ボクセル値D2n とその読み出した第2基準値T2とを比較する。そして、第2ボクセル値D2n が第2基準値T2を超えているとき(ステップS6−15でYES)、CPU7は、現在位置X2n はマスク領域M1内であると判断し(ステップS6―20)、融合処理を実行する。一方、第2ボクセル値D2n が第2基準値T2以下であるとき(ステップS6−15でNO)、CPU7は、現在位置X2n はマスク領域M1外であると判断し(ステップS6―25)、融合処理を実行しない。
そして、いずれの場合も、CPU7は現在位置X1n ,X2n を進行し(図16のステップS2−55)、現在位置X1n または現在位置X2n が予め定めた終了位置であるか否かを判断する(ステップS2−60)。そして、現在位置X1n または現在位置X2n が終了位置でないとき(ステップS2−60でNO)、CPU7は、再度マスク生成処理及び融合処理を実行する。一方、現在位置X1n または現在位置X2n が終了位置であるとき(ステップS2−60でYES)、反射光En をピクセル値として(ステップS2−65)、フレームFRを構成する1ピクセルについての融合画像データ生成処理を終了する。
その後、CPU7は、第1及び第2実施形態と同様に画像融合処理を実行し、図15に示すように、腫瘍が疑われる箇所である、PET値(第2ボクセル値D2n )の高い箇所のみを表示した融合画像P3を得る。これにより、観察したい箇所のみ表示した融合画像P3を観察することができる。
上記第3実施形態によれば、第1実施形態に記載の(1)〜(8)の効果に加えて以下のような効果を得ることができる。
(15) 第3実施形態によれば、第2ボクセル値D2n に応じてマスク領域M1を動的に生成したので、異常部分とその周辺のみを融合画像P3上に表示することができる。この結果、観察したい箇所のみ表示するため、容易でありながらも正確に観察対象を観察することができる。
(16) 第3実施形態によれば、第2ボクセル値D2n に応じてマスク領域M1を動的に生成したので、異常部分とその周辺のみを融合画像P3上に表示することができるこ
とから、融合画像P3の表示に要する時間を低減することができる。従って、一層、リアルタイムに融合画像P3を観察することができる。
(17) 第3実施形態によれば、融合処理前にマスク生成処理を実行したので、マスク領域M1外の部分について融合処理を行う必要がないことから、融合処理に要する時間を低減することができる。従って、一層、リアルタイムに融合画像P3を観察することができる。
(第4実施形態)
前記第1実施形態では、1台のワークステーションなどの計算機(コンピュータ)3が単独で画像融合処理を行ったが、本実施形態では、画像融合処理を構成する各処理のうち少なくとも1つの手段を複数のコンピュータが分散処理で行う。以下の実施形態において、前記第1〜第3実施形態と同様の部分については、同一の符号を付し、その詳細な説明は省略する。
例えば、複数台のワークステーションが接続された病院内ネットワークにおいて、少なくとも1つの処理を複数台のワークステーションが分散処理で行う構成を採用できる。以下、画像融合処理を分散処理で行う場合の例として、仮想光線Rの本数を分割する場合と、後処理だけ分割する場合の2つを示す。なお、説明の便宜上、図18及び図19に示すように、2台のワークステーションWS1,WS2を用いて、サイズ512×512の画像をつくる場合を例とするが、これを複数台のワークステーションで分散処理を行ってもよい。また、本実施形態では、ワークステーションWS2にのみGPU10を搭載しているとする。
(例1)図18に示すように、例1は、仮想光線Rを仮想光線RA(仮想光線R1〜仮想光線Rk)と仮想光線RB(仮想光線Rk+1〜仮想光線Rn)とに分割する場合である。仮想光線RAが通過するボクセルV1〜Vkと、仮想光線RBが通過するVk+1〜Vnとに分割する。また、仮想光線RAと仮想光線RBの分割は、第1ボクセルデータVD1(第2ボクセルデータVD2)の一部データが重複するように分割し、冗長性を持たせている。これは、仮想光線RA,RBが格子点上に到達しなかった場合に、第1ボクセル値D1n (第2ボクセル値D2n )を求める際の補間計算等に周囲のボクセルが必要になるからである。この場合、それぞれのワークステーションWS1,WS2で画像融合処理を行う。このようにすれば、ワークステーション毎の合成比率記憶部CV、光パラメータ記憶部MP及び反射光記憶部RVのメモリ量、転送量が全体の融合画像サイズの半分で済む。処理の手順は以下のようになる。
(1−1)ワークステーションWS1は、仮想光線RA上の第1ボクセルデータVD1(第2ボクセルデータVD2)(ボクセルV1〜Vk)について、融合処理を行う。そして、反射光EAn を算出し、その反射光EAn を反射光記憶部RVAに記憶させる。一方、ワークステーションWS2は、仮想光線RB上の第1ボクセルデータVD1(第2ボクセルデータVD2)(ボクセルVk+1〜Vn)について、融合処理を行う。そして、反射光EBn を算出し、その反射光EBn を反射光記憶部RVBに記憶させる。
(1−2)ワークステーションWS2の反射光記憶部RVBに記憶した反射光EBn をワークステーションWS1に転送する。このときの転送サイズは512×256で済む。
(1−3)ワークステーションWS1は、反射光記憶部RVAの反射光EAn と、反射光記憶部RVBの反射光EBn とを合成した反射光En を記憶した反射光記憶部RVに対して後処理を行い、CT画像データDT1とPET画像データDT2とが融合された融合画像P1(P2,P3)を得る。
(例2)図19に示すように、例2は、後処理だけ分割する場合である。この場合、ワークステーションWS1で第1ボクセルデータVD1、第2ボクセルデータVD2全体に
対して融合画像データ生成処理を行う。そして、高速画像処理に適しているGPU10を搭載したワークステーションWS2で後処理を行うようにすれば、後処理に要する時間を短縮できる。処理の手順は以下のようになる。
(2−1)ワークステーションWS1は、第1ボクセルデータVD1及び第2ボクセルデータVD2について、融合画像データ生成処理を行い、その算出した融合画像データDT3をメモリ部8aに記憶する。
(2−2)ワークステーションWS1のメモリ部8aに記憶した融合画像データDT3をワークステーションWS2に転送する。このときの転送サイズは512×512である。(2−3)ワークステーションWS2は、融合画像データDT3に対して後処理を行い、CT画像データDT1とPET画像データDT2とが融合された融合画像P1(P2,P3)を得る。
上記第4実施形態によれば、第1〜第3実施形態の(1)〜(17)の効果に加えて、以下のような効果を得ることができる。
(18)第4実施形態によれば、複数の計算機(コンピュータ)3による分散処理を採用するため画像融合処理の速度向上を図ることができるので、例えば、モニタ4の画面4aに表示される融合画像P1(P2,P3)のリアルタイム性を確保し易くなる。
(19)第4実施形態によれば、複数の計算機(コンピュータ)3による分散処理を採用するため、反射光記憶部RVのために使用する一台あたりのメモリ量を低減することができる。
なお、上記各実施形態は以下のように変更してもよい。
○上記第1実施形態では、合成比率決定処理において、第2ボクセル値D2n が予め定めた第1基準値T1より大きいとき、合成色γ3n はPET画像データDT2に従属しているPET色γ2n を用いたが、第1基準値T1を可変にしてもよい。すなわち、キーボード5等の入力装置を用いてその都度入力してもよいし、観察対象の種類に応じた第1基準値T1を複数種類予め記憶しておき、観察者がその複数種類の第2基準値T2からマウス6等の入力装置を用いて選択してもよい。これにより、一層動的に合成比率を決定することができる。
○上記第1実施形態では、まず、現在位置X1n における第1ボクセル値D1n 、勾配g1n を算出し(図8のステップS2−20)、現在位置X2n における第2ボクセル値D2n 、勾配g2n を算出した(ステップS2−25)。その後、CT不透明度α1n 、CTシェーディング係数β1n 、CT色γ1n を算出し(ステップS2−30)、PET不透明度α2n 、PETシェーディング係数β2n 、PET色γ2n を算出した(ステップS2−35)。これを、まず、現在位置X1n における第1ボクセル値D1n 、勾配g1n を算出し(ステップS2−20)、CT不透明度α1n 、CTシェーディング係数β1n 、CT色γ1n を算出する(ステップS2−30)。その後、現在位置X2n における第2ボクセル値D2n 、勾配g2n を算出し(ステップS2−25)、PET不透明度α2n 、PETシェーディング係数β2n 、PET色γ2n を算出(ステップS2−35)してもよい。
○上記第1実施形態では、PET画像データDT2に従属した光パラメータPを用いることによって、異常部分U1を強調表示した。これを、例えば、異常部分U1を強調していない融合画像P1を生成し、異常部分U1を表示した融合画像P1と交互に表示することによってその異常部分U1を点滅させるアニメーションを作成し、一層、異常部分U1を強調表示してもよい。これにより、観察者が異常部分U1を直感的に把握することができる。
○上記第2実施形態では、第1ボクセル値D1n 及び第2ボクセル値D2n に応じて、多段階で第1合成比率a及び第3合成比率cを決定し、第2合成比率bはCTシェーディング係数β1n を用いた。これを、まず、勾配g1n が所定の基準値以上の値であって、画面上に黒く表示されるか否かを判断してもよい。すなわち、勾配g1n が所定の基準値以上の値であるとき、その現在位置X1n ,X2n は観察時に影になる部分であって、CT値及びPET値に関わらず画面上に黒く表示されるので、第1合成比率a及び第3合成比率cを算出せず、初期値のままとしてもよい。これにより、合成比率決定処理の計算量を低減することができるため、一層、リアルタイムに融合画像P2を表示することができる。
○上記第2実施形態では、第1中間変数d、第2中間変数eは予め定めた定数を用いて算出したが、変数にしてもよい。すなわち、キーボード5等の入力装置を用いてその都度入力してもよいし、観察対象の種類に応じた変数を複数種類、予め記憶しておき、その複数種類の変数からマウス6等の入力装置を用いて選択してもよい。また、CPU7が計算によって求めてもよく、例えば、ボリュームデータのボクセル値の平均値や分散値、あるいは領域抽出処理の結果などを利用して計算することができる。これにより、一層動的に合成比率を決定することができる。
○上記第3実施形態では、マスク生成処理は、第2ボクセル値D2n が予め定めた第2基準値T2を超えたか否かによってマスク領域M1を設定したが、第2基準値T2を可変にしてもよい。すなわち、キーボード5等の入力装置を用いてその都度入力してもよいし、観察対象の種類に応じた第2基準値T2を複数種類、予め記憶しておき、その複数種類の第2基準値T2からマウス6等の入力装置を用いて選択してもよい。
○上記第3実施形態では、マスク生成処理後に融合処理を実行したが、融合処理後にマスク生成処理を実行してもよい。これにより、例えば、一度融合画像P3を観察後に所望の領域をマスク領域M1として設定したり、マスク領域M1を変更したりすることができる。
○上記第3実施形態では、マスク生成処理において、第2ボクセル値D2n が予め定めた第2基準値T2を超えたか否かによってマスク領域M1を設定したが、これを第1ボクセル値D1n 及び第2ボクセル値D2n がそれぞれ予め定めた基準値を超えたか否かによってマスク領域M1を設定してもよい。
詳述すると、マスク生成処理において、図20に示すように、CPU7は、現在位置X1n から半径1cm以内のボクセル値の最小値D1nminを算出し、メモリ部8aに記憶する(ステップS7−10)。次に、CPU7は、現在位置X2n における第2ボクセル値D2n を算出し、メモリ部8aに記憶する(ステップS7−15)。そして、CPU7は、最小値D1nminと第2ボクセル値D2n をそれぞれメモリ部8aから読み出し、最小値D1nminが「0」より小さいか、または第2ボクセル値D2n が「100」より大きいか否かを判断する(ステップS7−20)。いずれかの条件を満たすとき(ステップS7−20でYES)、CPU7は、現在位置X1n ,X2n はマスク領域M1内であると設定する(ステップS7−25)。一方、いずれの条件も満たさないとき(ステップS7−20でNO)、CPU7は、現在位置X1n ,X2n はマスク領域M1外であると設定する(ステップS7−30)。
これにより、例えば、肺に発生した腫瘍の観察等のように、肺には空気が多く含まれるが、実際に病変が発生するのは肺の実質的な組織であるとき、肺の組織及び腫瘍の両方を観察することができる。
○上記第4実施形態では、ネットワークを介して接続されたワークステーションWS1,WS2によってネットワーク分散処理を行った。これを、1台のコンピュータに多数のプロセッサを搭載して分散処理を行ってもよい。
○上記第4実施形態では、画像融合処理を分散処理で行う場合の例として、仮想光線Rの本数を分割する場合と、後処理だけ分割する場合の2つを示した。これを、融合画像P1(P2,P3)をライン(モニタ4の走査線)毎にレンダリングするときのそのライン毎に画像融合処理を分割するようにしてもよい。
○上記各実施形態では、部分反射光Fnは、合成不透明度α3n、合成色γ3n及び合成シェーディング係数β3n(または所定のシェーディング係数)を用いて算出した。これを、合成シェーディング係数β3n(または所定のシェーディング係数)を用いずに、合成不透明度α3n及び合成色γ3nのみを用いて部分反射光Fnを算出してもよい。
○上記各実施形態では、合成比率CR(第1合成比率a、第2合成比率b、第3合成比率c)は、サンプリング位置(現在位置X1n ,X2n )毎に決定した。これを、所定の数のサンプリング位置毎や、ボクセル値等に重要な変化があった場合のみ合成比率CRを決定してもよい。これにより、融合処理に要する計算量を低減し、ひいては融合画像P1(P2,P3)を表示するのに要する時間を短縮することができるため、融合画像P1(P2,P3)の表示のリアルタイム性を確保し易くなる。
○上記各実施形態では、現在位置X1n または現在位置X2n が終了位置に到達したとき融合処理を終了したが、現在位置X1n または現在位置X2n が終了位置に到達していないときであっても、減衰光An が「0」になったとき融合処理を終了してもよい。これにより、減衰光An が「0」であって反射光En が自明であるにも関わらず、終了位置に到達するまでサンプリングすることがないため、融合処理に要する計算量を低減し、ひいては融合画像P1(P2,P3)を表示するのに要する時間を短縮することができる。
○上記各実施形態では、ボクセル値から光パラメータを導出し光パラメータを用いて画像処理を行ったが、光パラメータはボクセル値から一意的に導ける値であれば何でもよく、物理的な光学特性とは何ら関係がない。例えば、ボクセル値が500〜1000の値の範囲内の時には光パラメータを1、範囲外の時には光パラメータは0と言った抽象的な値であってもよい。更に、光パラメータはボクセル値そのものであってもよい。
○上記各実施形態では、CT画像データDT1とPET画像データDT2に対して画像融合処理を行って融合画像を生成した。これをCT画像データDT1、PET画像データDT2以外のMRI画像データ、MRA画像データ等を用いて画像融合処理を行ってもよい。また、画像融合処理を行う画像データは2種類に限らず、3種類以上のモダリティによって撮影された画像データに対して画像融合処理を行って、融合画像を生成してもよい。
○上記各実施形態では、CT画像データDT1とPET画像データDT2に対して画像融合処理を行って融合画像P1(P2,P3)を生成した。これをどちらの画像データも同じ種類のモダリティで撮影した画像であってもよい。例えば、撮影時刻が異なったり、撮影条件が異なったり、造影条件が異なる画像データを融合することにより、診断に有益な情報を得ることができる。また、同一の画像データを2つ以上に複製し、それぞれに異なるフィルターを加えた生成した画像データを融合させてもよい。
○上記各実施形態では、レイキャスト法によってボリュームレンダリングを行った。これを最大値投影法(Maximum Intensity Projection)、最小値投影法(Minimum Intensity P
rojection)、平均値法、合計値法(Raysum)等の他のボリュームレンダリング方法を用いてボリュームレンダリングを行ってもよい。
○上記各実施形態では、3次元画像データに対して画像融合処理を行ったが、これを4次元以上の画像データに対して行ってもよい。
○上記各実施形態では、骨や臓器等の人体の部分について撮影されたCT画像、PET画像に対して画像融合処理を行ったが、CT撮影等が可能であれば、特に人体や動物、植物等の生物の組織に限らず、地質調査、鉱物探査、機械や装置類の構造材、電気回路のパターンを見る画像処理、LSIの故障診断等にも適用することができる。
○上記各実施形態では、骨や臓器等の人体の部分について撮影されたCT画像、PET画像に対して画像融合処理を行ったが、これを、計算機支援エンジニアリングシステムや科学技術計算の結果に適用してもよい。
第1実施形態の画像表示装置の概略構成図。 同じく、画像表示装置の概略構成を示すブロック図。 同じく、ボリュームレンダリングを説明する説明図。 同じく、ボリュームレンダリングを説明する説明図。 同じく、ボリュームレンダリングを説明するブロック図。 同じく、画像融合処理を説明するブロック図。 同じく、画像融合処理を説明するフローチャート。 同じく、融合画像データ生成処理を説明するフローチャート。 同じく、合成比率決定処理を説明するフローチャート。 同じく、融合画像について説明する模式図。 第2実施形態の画像融合処理を説明するブロック図。 同じく、合成比率決定処理を説明するフローチャート。 同じく、融合画像について説明する模式図。 第3実施形態の画像融合処理を説明するブロック図。 同じく、融合画像について説明する模式図。 同じく、融合画像データ生成処理を説明するフローチャート。 同じく、マスク生成処理を説明するフローチャート。 第4実施形態の画像融合処理の分散処理を示すブロック図。 同じく、画像融合処理の分散処理を示すブロック図。 別例のマスク生成処理を説明するフローチャート。 従来の画像融合処理を説明するブロック図。 同じく、画像融合処理を説明するブロック図。
符号の説明
α1n …CT不透明度、α2n …PET不透明度、α3n …合成不透明度、αn …不透明度、β1n …CTシェーディング係数、β2n …PETシェーディング係数、β3n …合成シェーディング係数、βn …シェーディング係数、γ1n …CT色、γ2n …PET色、γ3n …合成色、γn …色、a…第1合成比率、b…第2合成比率、c…第3合成比率、g1n ,g2n …勾配、An …減衰光、CR…合成比率、CV…合成比率記憶部、D1n …第1ボクセル値、D2n …第2ボクセル値、En ,EAn ,EBn …反射光、Fn …部分反射光、M1…マスク領域、MP…光パラメータ記憶部、P…光パラメータ、P1,P2,P3…融合画像、PX…ピクセル、R,RA,RB…仮想光線、RV…反射光記憶部、SV…基準値記憶部、T1…第1基準値、T2…第2基準値、U1,U2…異常部分、V1〜Vn…ボクセル、VD1…第1ボクセルデータ、VD2…第2ボクセルデータ、WS1,WS2…ワークステーション、1…画像表示装置、3…計算機(コンピュータ
)、4…出力部としてのモニタ、7…CPU、10…GPU。

Claims (21)

  1. 1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の3次元以上の画像データを2次元画面上に融合した融合画像を得る画像融合処理方法において、
    前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定段階と、
    前記合成比率に基づいて合成光パラメータを算出する算出段階と、
    前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出段階と、
    前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出段階と
    を備えたことを特徴とする画像融合処理方法。
  2. 請求項1に記載の画像融合処理方法において、
    前記光パラメータは色情報および不透明度情報であって、
    前記合成比率決定段階は、
    前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定段階と、
    前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定段階と
    を備え、
    前記算出段階は、
    前記色合成比率に基づいて合成色情報を算出する第1算出段階と、
    前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出段階と
    を備えたことを特徴とする画像融合処理方法。
  3. 請求項2に記載の画像融合処理方法において、
    前記合成比率決定段階は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定段階をさらに備え、
    前記算出段階は、
    前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出段階をさらに備え、
    前記部分反射光算出段階は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出することを特徴とする画像融合処理方法。
  4. 請求項1乃至3のいずれか1項に記載の画像融合処理方法において、
    前記合成比率決定段階は、
    前記複数の3次元以上の画像データのうち少なくとも1つの3次元以上の画像データに投射した仮想光線上のサンプリング位置のボクセル値と第1基準値とを比較する第1比較段階をさらに備え、
    前記第1比較段階での比較結果に応じて前記合成比率を決定することを特徴とする画像融合処理方法。
  5. 請求項1乃至3のいずれか1項に記載の画像融合処理方法において、
    前記合成比率決定段階は、
    前記合成比率を多段階で決定することを特徴とする画像融合処理方法。
  6. 請求項1乃至5のいずれか1項に記載の画像融合処理方法において、
    前記合成比率決定段階は、前記合成比率を前記ボクセル値の勾配に応じて設定することを特徴とする画像融合処理方法。
  7. 請求項1乃至6のいずれか1項に記載の画像融合処理方法において、
    前記ボクセル値のうちの1つと第2基準値とを比較する第2比較段階と、
    前記第2比較段階での比較結果に応じてマスク領域を生成するマスク領域生成段階と
    をさらに備えたことを特徴とする画像融合処理方法。
  8. 請求項1または2に記載の画像融合処理方法において、
    同一の複数の3次元以上の画像データのサンプリング位置のボクセル値の前記光パラメータをそれぞれ異なる合成比率で生成した複数の融合画像を用いてアニメーションを作成することを特徴とする画像融合処理方法。
  9. 1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の3次元以上の画像データを2次元画面上に融合した融合画像を得る画像融合処理プログラムにおいて、
    前記1つのコンピュータまたは複数のコンピュータを、
    前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定手段と、
    前記合成比率に基づいて合成光パラメータを算出する算出手段と、
    前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出手段と、
    前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出手段と
    して機能させることを特徴とする画像融合処理プログラム。
  10. 請求項9に記載の画像融合処理プログラムにおいて、
    前記光パラメータは色情報および不透明度情報であって、
    前記合成比率決定手段は、
    前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定手段と、
    前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定手段と
    を備え、
    前記算出手段は、
    前記色合成比率に基づいて合成色情報を算出する第1算出手段と、
    前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出手段と
    を備えたことを特徴とする画像融合処理プログラム。
  11. 請求項10に記載の画像融合処理プログラムにおいて、
    前記合成比率決定手段は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定手段をさらに備え、
    前記算出手段は、
    前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出手段をさらに
    備え、
    前記部分反射光算出手段は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出することを特徴とする画像融合処理プログラム。
  12. 請求項9乃至11のいずれか1項に記載の画像融合処理プログラムにおいて、
    前記合成比率決定手段は、
    前記複数の3次元以上の画像データのうち少なくとも1つの3次元以上の画像データに投射した仮想光線上のサンプリング位置のボクセル値と第1基準値とを比較する第1比較手段をさらに備え、
    前記第1比較手段での比較結果に応じて前記合成比率を決定することを特徴とする画像融合処理プログラム。
  13. 請求項9乃至11のいずれか1項に記載の画像融合処理プログラムにおいて、
    前記合成比率決定手段は、
    前記合成比率を多段階で決定することを特徴とする画像融合処理プログラム。
  14. 請求項9乃至13のいずれか1項に記載の画像融合処理プログラムにおいて、
    前記合成比率決定手段は、前記合成比率を前記ボクセル値の勾配に応じて設定することを特徴とする画像融合処理プログラム。
  15. 請求項9乃至14のいずれか1項に記載の画像融合処理プログラムにおいて、
    前記ボクセル値のうちの1つと第2基準値とを比較する第2比較手段と、
    前記第2比較手段での比較結果に応じてマスク領域を生成するマスク領域生成手段と
    をさらに備えたことを特徴とする画像融合処理プログラム。
  16. 請求項9または10に記載の画像融合処理プログラムにおいて、
    同一の複数の3次元以上の画像データのサンプリング位置のボクセル値の前記光パラメータをそれぞれ異なる合成比率で生成した複数の融合画像を用いてアニメーションを作成することを特徴とする画像融合処理プログラム。
  17. 1つのコンピュータが単独処理でまたは複数のコンピュータが単独処理と分散処理とのうち少なくとも1つの処理で、同一の観察対象から得た複数の3次元以上の画像データの位置関係を互いに対応させ、前記複数の3次元以上の画像データにそれぞれ対応する複数の仮想光線を投射し、同複数の3次元以上の画像データを2次元画面上に融合した融合画像を得る画像融合処理装置であって、
    前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値にそれぞれ対応付けられた複数の光パラメータを前記融合画像を生成するために合成する比率を合成比率として前記仮想光線上で少なくとも二度決定する合成比率決定手段と、
    前記合成比率に基づいて合成光パラメータを算出する算出手段と、
    前記合成光パラメータを基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出する部分反射光算出手段と、
    前記部分反射光を積算して前記融合画像の画素値を算出する画素値算出手段と
    を備えたことを特徴とする画像融合処理装置。
  18. 請求項17に記載の画像融合処理装置において、
    前記光パラメータは色情報および不透明度情報であって、
    前記合成比率決定手段は、
    前記ボクセル値に対応付けられた色情報を合成する色合成比率を前記仮想光線上で少なくとも二度決定する色合成比率決定手段と、
    前記ボクセル値に対応付けられた不透明度情報を合成する不透明度合成比率を前記仮想光線上で少なくとも二度決定する不透明度合成比率決定手段と
    を備え、
    前記算出手段は、
    前記色合成比率に基づいて合成色情報を算出する第1算出手段と、
    前記不透明度合成比率に基づいて合成不透明度情報を算出する第2算出手段と
    を備えたことを特徴とする画像融合処理装置。
  19. 請求項18に記載の画像融合処理装置において、
    前記合成比率決定手段は、前記複数の仮想光線上のそれぞれ対応するサンプリング位置のボクセル値に対応付けられたシェーディング係数を合成する陰影合成比率を前記仮想光線上で少なくとも二度決定する陰影合成比率決定手段をさらに備え、
    前記算出手段は、
    前記陰影合成比率に基づいて合成シェーディング係数を算出する第3算出手段をさらに備え、
    前記部分反射光算出手段は、前記合成色情報及び前記合成不透明度情報及び前記合成シェーディング係数を基に前記仮想光線に対する部分反射光を前記サンプリング位置毎に算出することを特徴とする画像融合処理装置。
  20. 請求項17乃至19のいずれか1項に記載の画像融合処理装置において、
    前記合成比率決定手段、前記算出手段、前記部分反射光算出手段、前記画素値算出手段の一部または全部は、GPUであることを特徴とする画像融合処理装置。
  21. 請求項17乃至20のいずれか1項に記載の画像融合処理装置において、
    前記画素値算出手段が算出した画素値を出力部に出力可能な形式に変換する後処理は、GPUで行うことを特徴とする画像融合処理装置。
JP2005202313A 2005-07-11 2005-07-11 画像融合処理方法、画像融合処理プログラム、画像融合処理装置 Expired - Fee Related JP4267598B2 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2005202313A JP4267598B2 (ja) 2005-07-11 2005-07-11 画像融合処理方法、画像融合処理プログラム、画像融合処理装置
US11/484,059 US7817877B2 (en) 2005-07-11 2006-07-10 Image fusion processing method, processing program, and processing device

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2005202313A JP4267598B2 (ja) 2005-07-11 2005-07-11 画像融合処理方法、画像融合処理プログラム、画像融合処理装置

Publications (2)

Publication Number Publication Date
JP2007014706A true JP2007014706A (ja) 2007-01-25
JP4267598B2 JP4267598B2 (ja) 2009-05-27

Family

ID=37752359

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2005202313A Expired - Fee Related JP4267598B2 (ja) 2005-07-11 2005-07-11 画像融合処理方法、画像融合処理プログラム、画像融合処理装置

Country Status (2)

Country Link
US (1) US7817877B2 (ja)
JP (1) JP4267598B2 (ja)

Cited By (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR100766070B1 (ko) 2007-06-28 2007-10-12 (주)한성유아이엔지니어링 지피에스를 이용한 지형변화 모니터링 시스템
EP1988511A2 (en) * 2007-04-12 2008-11-05 FUJIFILM Corporation Method and apparatus for volume rendering
JP2008272470A (ja) * 2007-04-26 2008-11-13 General Electric Co <Ge> 被撮像体内の目的物の視覚性を改善するためのシステム及び方法
JP2009160306A (ja) * 2008-01-09 2009-07-23 Ziosoft Inc 画像表示装置、画像表示装置の制御方法、および画像表示装置の制御プログラム
JP2011525271A (ja) * 2008-05-23 2011-09-15 ザ オーストラリアン ナショナル ユニヴァーシティ 画像データ処理
CN102231844A (zh) * 2011-07-21 2011-11-02 西安电子科技大学 基于结构相似度和人眼视觉的视频图像融合性能评价方法
WO2012165044A1 (ja) * 2011-06-01 2012-12-06 株式会社日立メディコ 画像表示装置、画像表示システムおよび画像表示方法
JP2013040922A (ja) * 2011-07-21 2013-02-28 Ngk Spark Plug Co Ltd ガスセンサ
JP2014000182A (ja) * 2012-06-18 2014-01-09 Aze Ltd 医用画像生成装置及びプログラム
WO2014192829A1 (ja) * 2013-05-28 2014-12-04 株式会社東芝 医用画像処理装置
JP2016142664A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像解析技術
JP2016142665A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像中の腫瘍領域を抽出する技術
JP2016142666A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像中の腫瘍輪郭を抽出する技術
WO2016125349A1 (ja) * 2015-02-04 2016-08-11 日本メジフィジックス株式会社 核医学画像解析技術
JP2017174039A (ja) * 2016-03-23 2017-09-28 富士フイルム株式会社 画像分類装置、方法およびプログラム
US10346974B2 (en) 2017-05-18 2019-07-09 Toshiba Medical Systems Corporation Apparatus and method for medical image processing
CN112184554A (zh) * 2020-10-13 2021-01-05 重庆邮电大学 一种基于残差混合膨胀卷积的遥感图像融合方法

Families Citing this family (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4233547B2 (ja) * 2005-07-06 2009-03-04 ザイオソフト株式会社 画像表示処理方法
US8228400B2 (en) * 2009-04-17 2012-07-24 Sony Corporation Generation of simulated long exposure images in response to multiple short exposures
US9959594B2 (en) 2010-07-22 2018-05-01 Koninklijke Philips N.V. Fusion of multiple images
WO2012048070A1 (en) * 2010-10-07 2012-04-12 Siemens Corporation Non-rigid composition of multiple overlapping medical imaging volumes
EP2647976A4 (en) * 2010-12-03 2017-08-09 Sony Corporation 3d data analysis device, 3d data analysis method, and 3d data analysis program
US9824302B2 (en) * 2011-03-09 2017-11-21 Siemens Healthcare Gmbh Method and system for model-based fusion of multi-modal volumetric images
US9314160B2 (en) * 2011-12-01 2016-04-19 Varian Medical Systems, Inc. Systems and methods for real-time target validation for image-guided radiation therapy
CN104103083A (zh) * 2013-04-03 2014-10-15 株式会社东芝 图像处理装置和方法以及医学成像设备
CN103345746B (zh) * 2013-06-25 2016-12-28 上海交通大学 对ct图片重建出三维图形的方法
US9842424B2 (en) * 2014-02-10 2017-12-12 Pixar Volume rendering using adaptive buckets
US10169909B2 (en) * 2014-08-07 2019-01-01 Pixar Generating a volumetric projection for an object
KR102412122B1 (ko) * 2015-05-27 2022-06-23 삼성전자주식회사 의료 영상 디스플레이 방법 및 장치
JP2017099616A (ja) * 2015-12-01 2017-06-08 ソニー株式会社 手術用制御装置、手術用制御方法、およびプログラム、並びに手術システム
WO2019045144A1 (ko) * 2017-08-31 2019-03-07 (주)레벨소프트 의료용 항법 장치를 위한 의료 영상 처리 장치 및 의료 영상 처리 방법
CN108876783B (zh) 2018-06-27 2021-02-05 上海联影医疗科技股份有限公司 图像融合方法及系统、医疗设备和图像融合终端
CN113421323B (zh) * 2021-06-21 2022-07-12 北京理工大学 一种基于光线投射的双模态体数据感兴趣区域融合方法
CN116363031B (zh) * 2023-02-28 2023-11-17 锋睿领创(珠海)科技有限公司 基于多维光学信息融合的成像方法、装置、设备及介质

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH04183446A (ja) * 1990-11-19 1992-06-30 Res Dev Corp Of Japan 画像合成による手術支援システム
JPH096986A (ja) * 1995-06-19 1997-01-10 Hitachi Ltd マルチモダリティに対応した3次元画像合成表示装置
JP2001291091A (ja) * 2000-01-31 2001-10-19 Mitsubishi Electric Corp 画像処理装置および方法
JP2003036448A (ja) * 2001-07-25 2003-02-07 Ziosoft Inc 3次元画像表示方法
JP2003157444A (ja) * 2001-11-21 2003-05-30 Mitsubishi Electric Corp 3次元情報合成装置、3次元情報表示装置及び3次元情報合成方法
JP2003233600A (ja) * 2001-12-03 2003-08-22 Ziosoft Inc ボリュームレンダリング処理方法、ボリュームレンダリング処理システム、計算機及びプログラム
WO2003077202A1 (en) * 2002-03-06 2003-09-18 Siemens Corporate Research, Inc. Visualization of volume-volume fusion
JP2004215846A (ja) * 2003-01-14 2004-08-05 Ziosoft Inc ボリュームレンダリング画像処理方法、ボリュームレンダリング画像処理装置及びプログラム

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR2597227B1 (fr) * 1986-04-14 1992-09-11 Pixar Procede pour produire un affichage bidimensionnel representant un ensemble de donnees tridimensionnelles
US5839440A (en) * 1994-06-17 1998-11-24 Siemens Corporate Research, Inc. Three-dimensional image registration method for spiral CT angiography
US6639595B1 (en) * 2000-08-23 2003-10-28 Nintendo Co., Ltd. Achromatic lighting in a graphics system and method
JP2005142680A (ja) * 2003-11-04 2005-06-02 Olympus Corp 画像処理装置
US20050195206A1 (en) * 2004-03-04 2005-09-08 Eric Wogsberg Compositing multiple full-motion video streams for display on a video monitor

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH04183446A (ja) * 1990-11-19 1992-06-30 Res Dev Corp Of Japan 画像合成による手術支援システム
JPH096986A (ja) * 1995-06-19 1997-01-10 Hitachi Ltd マルチモダリティに対応した3次元画像合成表示装置
JP2001291091A (ja) * 2000-01-31 2001-10-19 Mitsubishi Electric Corp 画像処理装置および方法
JP2003036448A (ja) * 2001-07-25 2003-02-07 Ziosoft Inc 3次元画像表示方法
JP2003157444A (ja) * 2001-11-21 2003-05-30 Mitsubishi Electric Corp 3次元情報合成装置、3次元情報表示装置及び3次元情報合成方法
JP2003233600A (ja) * 2001-12-03 2003-08-22 Ziosoft Inc ボリュームレンダリング処理方法、ボリュームレンダリング処理システム、計算機及びプログラム
WO2003077202A1 (en) * 2002-03-06 2003-09-18 Siemens Corporate Research, Inc. Visualization of volume-volume fusion
JP2004215846A (ja) * 2003-01-14 2004-08-05 Ziosoft Inc ボリュームレンダリング画像処理方法、ボリュームレンダリング画像処理装置及びプログラム

Cited By (22)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP1988511A2 (en) * 2007-04-12 2008-11-05 FUJIFILM Corporation Method and apparatus for volume rendering
JP2008272470A (ja) * 2007-04-26 2008-11-13 General Electric Co <Ge> 被撮像体内の目的物の視覚性を改善するためのシステム及び方法
KR100766070B1 (ko) 2007-06-28 2007-10-12 (주)한성유아이엔지니어링 지피에스를 이용한 지형변화 모니터링 시스템
JP2009160306A (ja) * 2008-01-09 2009-07-23 Ziosoft Inc 画像表示装置、画像表示装置の制御方法、および画像表示装置の制御プログラム
JP2011525271A (ja) * 2008-05-23 2011-09-15 ザ オーストラリアン ナショナル ユニヴァーシティ 画像データ処理
US9478021B2 (en) 2011-06-01 2016-10-25 Hitachi, Ltd. Image display device, image display system, and image display method
WO2012165044A1 (ja) * 2011-06-01 2012-12-06 株式会社日立メディコ 画像表示装置、画像表示システムおよび画像表示方法
JP2012252435A (ja) * 2011-06-01 2012-12-20 Hitachi Medical Corp 画像表示装置、画像表示システムおよび画像表示方法
CN102231844A (zh) * 2011-07-21 2011-11-02 西安电子科技大学 基于结构相似度和人眼视觉的视频图像融合性能评价方法
JP2013040922A (ja) * 2011-07-21 2013-02-28 Ngk Spark Plug Co Ltd ガスセンサ
JP2014000182A (ja) * 2012-06-18 2014-01-09 Aze Ltd 医用画像生成装置及びプログラム
WO2014192829A1 (ja) * 2013-05-28 2014-12-04 株式会社東芝 医用画像処理装置
US10110874B2 (en) 2013-05-28 2018-10-23 Toshiba Medical Systems Corporation Medical-image processing apparatus generating plural parallax images with different viewpoint positions based on adjusting parallactic angles
JP2016142665A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像中の腫瘍領域を抽出する技術
WO2016125349A1 (ja) * 2015-02-04 2016-08-11 日本メジフィジックス株式会社 核医学画像解析技術
JP2016142666A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像中の腫瘍輪郭を抽出する技術
JP2016142664A (ja) * 2015-02-04 2016-08-08 日本メジフィジックス株式会社 核医学画像解析技術
US10395364B2 (en) 2015-02-04 2019-08-27 Nihon Medi-Physics Co., Ltd. Nuclear medical image analysis technique
JP2017174039A (ja) * 2016-03-23 2017-09-28 富士フイルム株式会社 画像分類装置、方法およびプログラム
US10346974B2 (en) 2017-05-18 2019-07-09 Toshiba Medical Systems Corporation Apparatus and method for medical image processing
CN112184554A (zh) * 2020-10-13 2021-01-05 重庆邮电大学 一种基于残差混合膨胀卷积的遥感图像融合方法
CN112184554B (zh) * 2020-10-13 2022-08-23 重庆邮电大学 一种基于残差混合膨胀卷积的遥感图像融合方法

Also Published As

Publication number Publication date
JP4267598B2 (ja) 2009-05-27
US7817877B2 (en) 2010-10-19
US20070098299A1 (en) 2007-05-03

Similar Documents

Publication Publication Date Title
JP4267598B2 (ja) 画像融合処理方法、画像融合処理プログラム、画像融合処理装置
JP4213135B2 (ja) 展開画像投影方法、展開画像投影プログラム、展開画像投影装置
CN106981098B (zh) 虚拟场景组分的视角表示
JP4335817B2 (ja) 関心領域指定方法、関心領域指定プログラム、関心領域指定装置
JP3570576B2 (ja) マルチモダリティに対応した3次元画像合成表示装置
US7889194B2 (en) System and method for in-context MPR visualization using virtual incision volume visualization
US7424140B2 (en) Method, computer program product, and apparatus for performing rendering
US20090174729A1 (en) Image display device and control method thereof
JP4233547B2 (ja) 画像表示処理方法
JP4109224B2 (ja) 展開画像投影方法、展開画像投影プログラム、展開画像投影装置
JP2009034503A (ja) トモシンセシス画像を表示するための方法及びシステム
US20110082667A1 (en) System and method for view-dependent anatomic surface visualization
US20090003668A1 (en) Image processing method, image processing program, and image processing device
US20080259080A1 (en) Image processing method, apparatus, and program
JP2008543477A (ja) カーブした細長い構造の切断面を視覚化する方法
JP4634437B2 (ja) 展開画像投影方法、展開画像投影プログラム、展開画像投影装置
JP4563326B2 (ja) 画像処理方法および画像処理プログラム
JP4122314B2 (ja) 投影画像処理方法、投影画像処理プログラム、投影画像処理装置
JP4896470B2 (ja) 画像処理装置、医用画像診断装置及び画像処理方法
JP2006136619A (ja) 画像処理方法
JP2010000306A (ja) 医用画像診断装置、画像処理装置、及びプログラム
JP4709603B2 (ja) 医用画像処理装置
JP6106260B2 (ja) ボリュームレンダリング
EP3809376A2 (en) Systems and methods for visualizing anatomical structures
JP5245811B2 (ja) ボクセル配列可視化装置

Legal Events

Date Code Title Description
A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20080715

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20080916

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20081104

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20081225

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20090127

A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20090218

R150 Certificate of patent or registration of utility model

Free format text: JAPANESE INTERMEDIATE CODE: R150

FPAY Renewal fee payment (event date is renewal date of database)

Free format text: PAYMENT UNTIL: 20120227

Year of fee payment: 3

LAPS Cancellation because of no payment of annual fees