JP6759032B2 - 光音響装置、情報処理方法、及びプログラム - Google Patents

光音響装置、情報処理方法、及びプログラム Download PDF

Info

Publication number
JP6759032B2
JP6759032B2 JP2016188408A JP2016188408A JP6759032B2 JP 6759032 B2 JP6759032 B2 JP 6759032B2 JP 2016188408 A JP2016188408 A JP 2016188408A JP 2016188408 A JP2016188408 A JP 2016188408A JP 6759032 B2 JP6759032 B2 JP 6759032B2
Authority
JP
Japan
Prior art keywords
wavelength
light
information
misalignment
image data
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.)
Active
Application number
JP2016188408A
Other languages
English (en)
Other versions
JP2018050774A (ja
Inventor
孝太郎 梅澤
孝太郎 梅澤
隆一 七海
隆一 七海
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 JP2016188408A priority Critical patent/JP6759032B2/ja
Priority to CN201710881220.3A priority patent/CN107865639B/zh
Priority to KR1020170123862A priority patent/KR20180034273A/ko
Priority to EP17193277.5A priority patent/EP3298958A1/en
Priority to US15/716,150 priority patent/US11529057B2/en
Publication of JP2018050774A publication Critical patent/JP2018050774A/ja
Application granted granted Critical
Publication of JP6759032B2 publication Critical patent/JP6759032B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/0093Detecting, measuring or recording by applying one single type of energy and measuring its conversion into another type of energy
    • A61B5/0095Detecting, measuring or recording by applying one single type of energy and measuring its conversion into another type of energy by applying light and detecting acoustic waves, i.e. photoacoustic measurements
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/145Measuring characteristics of blood in vivo, e.g. gas concentration, pH value; Measuring characteristics of body fluids or tissues, e.g. interstitial fluid, cerebral tissue
    • A61B5/14542Measuring characteristics of blood in vivo, e.g. gas concentration, pH value; Measuring characteristics of body fluids or tissues, e.g. interstitial fluid, cerebral tissue for measuring blood gases
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7271Specific aspects of physiological measurement analysis
    • A61B5/7278Artificial waveform generation or derivation, e.g. synthesising signals from measured signals
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/74Details of notification to user or communication with user or patient ; user input means
    • A61B5/742Details of notification to user or communication with user or patient ; user input means using visual displays
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/74Details of notification to user or communication with user or patient ; user input means
    • A61B5/7475User input or interface means, e.g. keyboard, pointing device, joystick
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/04Analysing solids
    • G01N29/06Visualisation of the interior, e.g. acoustic microscopy
    • G01N29/0654Imaging
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/50Image enhancement or restoration using two or more images, e.g. averaging or subtraction
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/80Geometric correction
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B2576/00Medical imaging apparatus involving image processing or analysis

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Pathology (AREA)
  • Public Health (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Biophysics (AREA)
  • Veterinary Medicine (AREA)
  • Biomedical Technology (AREA)
  • Acoustics & Sound (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Artificial Intelligence (AREA)
  • Physiology (AREA)
  • Psychiatry (AREA)
  • Signal Processing (AREA)
  • Theoretical Computer Science (AREA)
  • Optics & Photonics (AREA)
  • Biochemistry (AREA)
  • Analytical Chemistry (AREA)
  • Immunology (AREA)
  • Chemical & Material Sciences (AREA)
  • Ultra Sonic Daignosis Equipment (AREA)
  • Geometry (AREA)

Description

本発明は、光音響装置に関する。
医療分野において生体の生理的情報である機能情報のイメージングの研究が近年行われている。機能情報のイメージング技術の一つとして光音響イメージング(Photoacoustic Imaging:PAI)がある。
光音響イメージングでは、光が被検体に照射される。照射光は被検体を伝播・拡散し、そのエネルギーは被検体内で吸収される。その結果、光音響効果により音響波(以降、光音響波と呼ぶ)が発生する。この光音響波の受信信号が解析処理されることにより、被検体内部の光学特性値の空間分布が画像データとして取得される。
複数回光を被検体に照射し、複数回の光照射により発生する光音響波に基づいて画像データを生成する場合、複数回の光照射の間に被検体とプローブとの相対位置が変動すると、画像データの解像度が低下する。
特許文献1には、複数回の光照射によって得られる複数の画像データ間の位置ずれ量を算出することにより、被検体とプローブとの相対位置の変動量を推定する技術が開示されている。また、特許文献1には、推定された相対位置の変動量に基づいて、複数の画像データを位置合わせする技術が開示されている。
特開2014−140716号公報
ところで、光音響イメージングにおいては、互いに異なる複数の波長の光を被検体に照射することにより機能情報を取得することができる。
しかしながら、複数の波長のそれぞれの光の照射により発生する光音響波に基づいて複数の画像データを生成する場合、複数の画像データ間の位置ずれ量の推定精度が低下してしまう場合がある。
そこで、本発明は、複数の波長の光照射により発生する光音響波に基づいて複数の画像データを生成する場合にも、複数の画像データ間の位置ずれ情報を精度よく推定することのできる装置を提供することを目的とする。すなわち、本発明は、複数の波長の光照射により発生する光音響波に基づいて複数の画像データを生成する場合にも、被検体とプローブとの相対位置の変動量を精度よく推定することのできる装置を提供することを目的とする。
本発明の一実施形態に係る光音響装置は、第1の波長の光、及び、第1の波長とは異なる第2の波長の光のそれぞれを複数回被検体に照射する光照射手段と、第1の波長の光及び第2の波長の光のそれぞれの複数回の光照射により発生する光音響波を受信することにより、信号群を出力する受信手段と、信号群に基づいて、第1の波長に対応する第1の複数の画像データを生成し、第1の波長に対応する第1の複数の画像データに基づいて、第1の波長の光の照射タイミングにおける第1の位置ずれ情報を取得し、第1の位置ずれ情報に基づいて、第2の波長の光の照射タイミングにおける第2の位置ずれ情報を取得する処理手段と、を有する。
本発明によれば、被検体とプローブとの相対位置の変動量を精度よく推定することができる。
第1の実施形態に係る光音響装置の構成を示す模式図 第1の実施形態に係るプローブの模式図 第1の実施形態に係るコンピュータ及びその周辺機器の構成図 第1の実施形態に係る光音響装置の作動のフローチャート 第1の実施形態に係る光音響装置の測定位置を示す模式図 第1の実施形態に係るパルスボリュームデータを示す模式図 第1の実施形態に係る画像間の位置ずれ量を示す模式図 第1の実施形態に係る位置ずれ量の時間的な補間を説明するためのシーケンス図 第1の実施形態に係る光音響装置の測定位置を示す別の模式図 第1の実施形態に係る光音響装置の測定位置を示す別の模式図 第1の実施形態に係る位置ずれ補正を説明するための模式図 第2の実施形態に係る光音響装置の作動のフローチャート 第2の実施形態に係るGUIを示す模式図 第3の実施形態に係る光音響装置の作動のフローチャート 第4の実施形態に係るGUIを示す模式図
以下、本発明の実施の形態について、図面を用いて説明する。同一の構成要素には原則として同一の符号を付して、説明を省略する。
[第1の実施形態]
本明細書においては、撮像対象とプローブとの相対位置が変動する事象を「位置ずれ」と呼ぶ。また、撮像対象とプローブとの相対位置の変動量を「位置ずれ量」と呼ぶ。また、位置ずれに起因する、複数の画像データに現れる撮像対象の位置の変化量を、「複数の画像データ間の位置ずれ量」と呼ぶ。並進量、回転量、または変形量などの位置ずれを表現するパラメータを「位置ずれ情報」と呼ぶ。
位置ずれは、複数回の光照射の間に被検体が動くことやプローブが動くことにより発生する。例えば、ユーザがハンドヘルド型のプローブを把持して走査することや走査手段がプローブを機械的に走査することなどによりプローブは動く。これらの位置ずれに起因して、複数回の光照射により発生する光音響波に基づいて画像データを生成するときに解像度の低下が生じてしまう。
ところで、光照射により発生する光音響波の強度は、撮像対象の吸収係数に比例して決定される。また、吸収係数は波長依存性を有しているため、波長毎に光フルエンスが同じであっても、波長によって発生する光音響波の強度は変化する。
一方、特許文献1に記載されたように複数の画像データ間の位置ずれを推定する場合、複数の画像データ間で同じ撮像対象の画像強度が異なると、位置ずれの推定精度が低下してしまう場合がある。そのため、互いに異なる複数の波長のそれぞれの光照射によって得られた複数の画像データを用いて位置ずれ量を推定しようとすると、波長間で画像強度や解像度が異なる等の理由から、位置ずれの推定精度が低下してしまう場合がある。
そこで、本発明の一実施形態に係る光音響装置は、互いに異なる複数の波長の一部の波長に対応する複数の画像データ間の位置ずれ情報を取得し、当該位置ずれ情報をその他の波長に対応する位置ずれ情報の取得に用いる。これにより、特定の波長に対応する位置ずれ情報を基準に複数の波長に対応する位置ずれ情報を推定するため、波長間の画像強度の違いによる位置ずれ情報の取得精度に与える影響を抑制することができる。また、特定の波長に対応する位置ずれ情報を取得し、その他の波長に対応する位置ずれ情報を補間により取得することにより、その他の波長に対応する位置ずれ情報の取得に要する処理量を低減することができる。
すなわち、本発明の一実施形態に係る光音響装置においては、互いに異なる第1の波長及び第2の波長のそれぞれの光が複数回被検体に照射される。そして、第1の波長に対応する第1の複数の画像データが生成され、第1の複数の画像データ間の位置ずれ情報が取得される。ここで取得された位置ずれ情報が、第1の波長の光の照射タイミングにおける被検体とプローブとの相対位置の変動量(位置ずれ量)に相当する。そして、第1の位置ずれ情報に基づいて、第2の波長の光照射タイミングにおける被検体とプローブとの変動量(位置ずれ量)が推定される。
また、本発明の一実施形態に係る光音響装置は、上記のように推定された複数波長の光照射に対応する位置ずれ情報と、複数波長の光照射により発生した光音響波に起因する受信データ(信号群または画像データ)とに基づいて、機能情報が取得される。機能情報としては、オキシヘモグロビン濃度、デオキシヘモグロビン濃度、総ヘモグロビン量、または酸素飽和度等の物質の濃度を示す情報が挙げられる。総ヘモグロビン量とは、オキシヘモグロビンおよびデオキシヘモグロビンの和の総量のことである。酸素飽和度とは、全ヘモグロビンに対するオキシヘモグロビンの割合のことである。
なお、機能情報は、空間分布を表す画像データに限らず、数値や文字を表す情報であってもよい。例えば、機能情報は、被検体を構成する物質の平均的な濃度、物質濃度の空間分布の特定位置の値、物質濃度の空間分布の画素値の統計値(平均値や中央値など)などの情報を含む概念である。例えば、機能情報を示す画像として、被検体を構成する物質の平均的な濃度の数値が、表示部160に表示されてもよい。
本発明の一実施形態においては、オキシヘモグロビンのモル吸収係数とデオキシヘモグロビンのモル吸収係数が等しい波長の光を利用して得られる複数の画像データに基づいて、位置ずれ情報を取得してもよい。血管を測定対象とする場合、このような波長の光を用いて血管から発生する光音響波の強度は、酸素飽和度に依存しないため、血管毎の画像強度のばらつきが小さくなる。そのため、このような波長の光を用いることにより、位置ずれ情報の推定精度が高くなる傾向がある。なお、オキシヘモグロビンのモル吸収係数とデオキシヘモグロビンのモル吸収係数が完全に等しい波長(等モル吸収係数の波長)に限らず、略等しい波長であっても位置ずれ情報の推定精度は高くなる傾向がある。例えば、等モル吸収係数の波長の±10nm以内の波長を、略等しい波長として採用してもよい。また、例えば、オキシヘモグロビンのモル吸収係数とデオキシヘモグロビンのモル吸収係数が、等モル吸収係数から±10%以内に相当する波長を略等しい波長として採用してもよい。すなわち、位置ずれ情報の取得に好適な波長の光を利用して得られた位置ずれ情報に基づいて、その他の波長に対応する位置ずれ情報を取得することが好ましい。
以下、第1の実施形態の光音響装置の構成の構成及び処理について説明する。
本実施形態において、光音響装置を用いた例を説明する。図1を用いて本実施形態に係る光音響装置の構成を説明する。図1は、光音響装置全体の概略ブロック図である。本実施形態に係る光音響装置は、駆動部130、信号収集部140、コンピュータ150、表示部160、入力部170、及びプローブ180を有する。
図2は、本実施形態に係るプローブ180の模式図である。プローブ180は、光照射部110、及び受信部120を有する。測定対象は、被検体100である。
駆動部130は、光照射部110と受信部120を駆動し、機械的な走査を行う。光照射部110がパルス光113を被検体100に照射し、被検体100内で音響波が発生する。光に起因して光音響効果により発生する音響波を光音響波とも呼ぶ。受信部120は、光音響波を受信することによりアナログ信号としての電気信号(光音響信号)を出力する。
信号収集部140は、受信部120から出力されたアナログ信号をデジタル信号に変換し、コンピュータ150に出力する。コンピュータ150は、信号収集部140から出力されたデジタル信号を、超音波または光音響波に由来する信号データとして記憶する。
コンピュータ150は、記憶されたデジタル信号に対して信号処理を行うことにより、光音響画像を表す画像データを生成する。また、コンピュータ150は、得られた画像データに対して画像処理を施した後に、画像データを表示部160に出力する。表示部160は、光音響画像を表示する。操作者としての医師や技師等は、表示部160に表示された光音響画像を確認することにより、診断を行うことができる。表示画像は、操作者やコンピュータ150からの保存指示に基づいて、コンピュータ150内のメモリや、モダリティとネットワークで接続されたデータ管理システムなどに保存される。
また、コンピュータ150は、光音響装置に含まれる構成の駆動制御も行う。また、表示部160は、コンピュータ150で生成された画像の他にGUIなどを表示してもよい。入力部170は、操作者が情報を入力できるように構成されている。操作者は、入力部170を用いて測定開始や終了、作成画像の保存指示などの操作を行うことができる。
本実施形態に係る光音響装置により得られる光音響画像は、光照射により発生した光音響波に由来するあらゆる画像を含む概念である。光音響画像は、光音響波の発生音圧(初期音圧)、光吸収エネルギー密度、及び光吸収係数、などの少なくとも1つの空間分布を表す画像データである。
以下、本実施形態に係る光音響装置の各構成の詳細を説明する。
(光照射部110)
光照射部110は、パルス光113を発する光源と、光源から射出されたパルス光113を被検体100へ導く光学系とを含む。なお、パルス光は、いわゆる矩形波、三角波などの光を含む。
光源が発する光のパルス幅としては、1ns以上、100ns以下のパルス幅であってもよい。また、光の波長として400nmから1600nm程度の範囲の波長であってもよい。血管を高解像度でイメージングする場合は、血管での吸収が大きい波長(400nm以上、700nm以下)を用いてもよい。生体の深部をイメージングする場合には、生体の背景組織(水や脂肪など)において典型的に吸収が少ない波長(700nm以上、1100nm以下)の光を用いてもよい。
光源としては、レーザーや発光ダイオードを用いることができる。また、複数波長の光を用いて測定する際には、波長の変換が可能な光源であってもよい。なお、複数波長を被検体に照射する場合、互いに異なる波長の光を発生する複数台の光源を用意し、それぞれの光源から交互に照射することも可能である。複数台の光源を用いた場合もそれらをまとめて光源として表現する。レーザーとしては、固体レーザー、ガスレーザー、色素レーザー、半導体レーザーなど様々なレーザーを使用することができる。例えば、Nd:YAGレーザーやアレキサンドライトレーザーなどのパルスレーザーを光源として用いてもよい。また、Nd:YAGレーザー光を励起光とするTi:saレーザーやOPO(Optical Parametric Oscillators)レーザーを光源として用いてもよい。また、光源としてマイクロウェーブ源を用いてもよい。
光学系には、レンズ、ミラー、光ファイバ等の光学素子を用いることができる。***等を被検体100とする場合、パルス光のビーム径を広げて照射するために、光学系の光出射部は光を拡散させる拡散板等で構成されていてもよい。一方、光音響顕微鏡においては、解像度を上げるために、光学系の光出射部はレンズ等で構成し、ビームをフォーカスして照射してもよい。
なお、光照射部110が光学系を備えずに、光源から直接被検体100にパルス光113を照射してもよい。
(受信部120)
受信部120は、音響波を受信することにより電気信号を出力するトランスデューサ121と、トランスデューサ121を支持する支持体122とを含む。また、トランスデューサ121は、音響波を送信することもできる。便宜上、図2においてはトランスデューサ121を1つしか図示していないが、受信部120は複数のトランスデューサを含んでいてもよい。
トランスデューサ121を構成する部材としては、PZT(チタン酸ジルコン酸鉛)に代表される圧電セラミック材料や、PVDF(ポリフッ化ビニリデン)に代表される高分子圧電膜材料などを用いることができる。また、圧電素子以外の素子を用いてもよい。例えば、静電容量型トランスデューサ(CMUT:Capacitive Micro−machined Ultrasonic Transducers)、ファブリペロー干渉計を用いたトランスデューサなどを用いることができる。なお、音響波を受信することにより電気信号を出力できる限り、いかなるトランスデューサを採用してもよい。また、トランスデューサにより得られる信号は時間分解信号である。つまり、トランスデューサにより得られる信号の振幅は、各時刻にトランスデューサで受信される音圧に基づく値(例えば、音圧に比例した値)を表したものである。
光音響波を構成する周波数成分は、典型的には100KHzから100MHzであり、トランスデューサ121として、これらの周波数を検出することのできるものを採用することができる。
支持体122は、機械的強度が高い金属材料などから構成されていてもよい。照射光を被検体に多く入射させるために、支持体122の被検体100側の表面に鏡面もしくは光散乱させる加工が行われていてもよい。本実施形態において支持体122は半球殻形状であり、半球殻上に複数のトランスデューサ121を支持できるように構成されている。この場合、支持体122に配置されたトランスデューサ121の指向軸は半球の曲率中心付近に集まる。そして、複数のトランスデューサ121から出力された信号を用いて画像化したときに曲率中心付近の画質が高くなる。なお、支持体122はトランスデューサ121を支持できる限り、いかなる構成であってもよい。支持体122は、1Dアレイ、1.5Dアレイ、1.75Dアレイ、2Dアレイと呼ばれるような平面又は曲面内に、複数のトランスデューサを並べて配置してもよい。
また、支持体122は音響マッチング材を貯留する容器として機能してもよい。すなわち、支持体122をトランスデューサ121と被検体100との間に音響マッチング材を配置するための容器としてもよい。
また、受信部120が、トランスデューサ121から出力される時系列のアナログ信号を増幅する増幅器を備えてもよい。また、受信部120が、トランスデューサ121から出力される時系列のアナログ信号を時系列のデジタル信号に変換するA/D変換器を備えてもよい。すなわち、受信部120が後述する信号収集部140を備えてもよい。
なお、音響波を様々な角度で検出できるようにするために、理想的には被検体100を全周囲から囲むようにトランスデューサ121を配置してもよい。ただし、被検体100が大きく全周囲を囲むようにトランスデューサを配置できない場合は、図2に示したように半球状の支持体122上にトランスデューサを配置して全周囲を囲む状態に近づけてもよい。図2においては、トランスデューサ121を1つしか図示していないが、複数のトランスデューサを半球上の支持体122上に配置してもよい。
なお、トランスデューサの配置や数及び支持体の形状は被検体に応じて最適化すればよく、本発明に関してはあらゆる受信部120を採用することができる。
受信部120と被検体100との間の空間は、光音響波が伝播することができる媒質で満たす。この媒質には、音響波が伝搬でき、被検体100やトランスデューサ121との界面において音響特性が整合し、できるだけ光音響波の透過率が高い材料を採用する。例えば、この媒質には、水、超音波ジェルなどを採用することができる。
(駆動部130)
駆動部130は、被検体100と受信部120との相対位置を変更する部分である。本実施形態では、駆動部130は、支持体122をXY方向に移動させる装置であり、ステッピングモータを搭載した電動のXYステージある。駆動部130は、駆動力を発生させるステッピングモータなどのモータと、駆動力を伝達させる駆動機構と、受信部120の位置情報を検出する位置センサとを含む。駆動機構としては、リードスクリュー機構、リンク機構、ギア機構、油圧機構、などを用いることができる。また、位置センサとしては、エンコーダー、可変抵抗器、などを用いたポテンショメータなどを用いることができる。
なお、駆動部130は被検体100と受信部120との相対位置をXY方向(二次元)に変更させるものに限らず、一次元または三次元に変更させてもよい。
なお、駆動部130は、被検体100と受信部120との相対的な位置を変更できれば、受信部120を固定し、被検体100を移動させてもよい。被検体100を移動させる場合は、被検体100を指示する被検体支持部(不図示)を動かすことで被検体100を移動させる構成などが考えられる。また、被検体100と受信部120の両方を移動させてもよい。
また駆動部130は、相対位置を連続的に移動させてもよいし、ステップアンドリピートによって移動させてもよい。駆動部130は、電動ステージであってもよいし、手動ステージであってもよい。
また、本実施形態では、駆動部130は光照射部110と受信部120を同時に駆動して走査を行っているが、光照射部110だけを駆動したり、受信部120だけを駆動したりしてもよい。
(信号収集部140)
信号収集部140は、トランスデューサ121から出力されたアナログ信号である電気信号を増幅するアンプと、アンプから出力されたアナログ信号をデジタル信号に変換するA/D変換器とを含む。信号収集部140は、FPGA(Field Programmable Gate Array)チップなどで構成されてもよい。信号収集部140から出力されるデジタル信号は、コンピュータ150内の記憶部152に記憶される。信号収集部140は、Data Acquisition System(DAS)とも呼ばれる。本明細書において電気信号は、アナログ信号もデジタル信号も含む概念である。なお、信号収集部140は、光照射部110の光射出部に取り付けられた光検出センサと接続されており、パルス光113が光照射部110から射出されたことをトリガーに、同期して処理を開始してもよい。また、信号収集部140は、フリーズボタンなどを用いてなされる指示をトリガーに同期して、当該処理を開始してもよい。
(コンピュータ150)
コンピュータ150は、演算部151、記憶部152、制御部153を含む。各構成の機能については処理フローの説明の際に説明する。
演算部151としての演算機能を担うユニットは、CPUやGPU(Graphics Processing Unit)等のプロセッサ、FPGA(Field Programmable Gate Array)チップ等の演算回路で構成されることができる。これらのユニットは、単一のプロセッサや演算回路から構成されるだけでなく、複数のプロセッサや演算回路から構成されていてもよい。演算部151は、入力部170から、被検体音速や保持カップの構成などの各種パラメータを受けて、受信信号を処理してもよい。
記憶部152は、ROM(Read only memory)、磁気ディスクやフラッシュメモリなどの非一時記憶媒体で構成することができる。また、記憶部152は、RAM(Random Access Memory)などの揮発性の媒体であってもよい。なお、プログラムが格納される記憶媒体は、非一時記憶媒体である。なお、記憶部152は、1つの記憶媒体から構成されるだけでなく、複数の記憶媒体から構成されていてもよい。
記憶部152は、後述する方法で演算部151により生成される光音響画像を示す画像データを保存することができる。
制御部153は、CPUなどの演算素子で構成される。制御部153は、光音響装置の各構成の動作を制御する。制御部153は、入力部170からの測定開始などの各種操作による指示信号を受けて、光音響装置の各構成を制御してもよい。また、制御部153は、記憶部152に格納されたプログラムコードを読み出し、光音響装置の各構成の作動を制御する。
コンピュータ150は専用に設計されたワークステーションであってもよい。また、コンピュータ150の各構成は異なるハードウェアによって構成されてもよい。また、コンピュータ150の少なくとも一部の構成は単一のハードウェアで構成されてもよい。
図3は、本実施形態に係るコンピュータ150の具体的な構成例を示す。本実施形態に係るコンピュータ150は、CPU154、GPU155、RAM156、ROM157、外部記憶装置158から構成される。また、コンピュータ150には、表示部160としての液晶ディスプレイ161、入力部170としてのマウス171、キーボード172が接続されている。
また、コンピュータ150および複数のトランスデューサ121は、共通の筺体に収められた構成で提供されてもよい。ただし、筺体に収められたコンピュータで一部の信号処理を行い、残りの信号処理を筺体の外部に設けられたコンピュータで行ってもよい。この場合、筺体の内部および外部に設けられたコンピュータを総称して、本実施形態に係るコンピュータとすることができる。すなわち、コンピュータを構成するハードウェアが一つの筺体に収められていなくてもよい。
(表示部160)
表示部160は、液晶ディスプレイや有機EL(Electro Luminescence)などのディスプレイである。コンピュータ150により得られた被検体情報等に基づく画像や特定位置の数値等を表示する装置である。表示部160は、画像や装置を操作するためのGUIを表示してもよい。なお、被検体情報の表示にあたっては、表示部160またはコンピュータ150において画像処理(輝度値の調整等)を行った上で表示することもできる。
(入力部170)
入力部170としては、ユーザーが操作可能な、マウスやキーボードなどで構成される操作コンソールを採用することができる。入力できる内容としては、画像再構成の条件の選択であったり、位置ずれ補正の方法の選択であったり、補間方法の選択であったりしてもよい。また、スライダーバーのような形で、合成された画像を確認しながら、位置ずれ量の足し合わせの重みなどを変更できるようにしてもよい。また、表示部160をタッチパネルで構成し、表示部160を入力部170として利用してもよい。
なお、光音響装置の各構成はそれぞれ別の装置として構成されてもよいし、一体となった1つの装置として構成されてもよい。また、光音響装置の少なくとも一部の構成が一体となった1つの装置として構成されてもよい。
(被検体100)
被検体100は光音響装置を構成するものではないが、以下に説明する。本実施形態に係る光音響装置は、人や動物の悪性腫瘍や血管疾患などの診断や化学治療の経過観察などを目的として使用できる。よって、被検体100としては、生体、具体的には人体や動物の***や各臓器、血管網、頭部、頸部、腹部、手指および足指を含む四肢などの診断の対象部位が想定される。例えば、人体が測定対象であれば、オキシヘモグロビンあるいはデオキシヘモグロビンやそれらを含む多く含む血管あるいは腫瘍の近傍に形成される新生血管などを光吸収体の対象としてもよい。また、頸動脈壁のプラークなどを光吸収体の対象としてもよい。また、メチレンブルー(MB)、インドシニアングリーン(ICG)などの色素、金微粒子、またはそれらを集積あるいは化学的に修飾した外部から導入した物質を光吸収体としてもよい。
以下、図4に示すフローチャートに沿って、本実施形態に係る情報処理を含む光音響装置の作動を説明する。
(S100:波長λ1の光及び波長λ2の光をそれぞれ複数回照射し、光音響波を受信する工程)
光照射部110は、互いに異なる波長λ1及び波長λ2の光のそれぞれを複数回照射し、受信部120は、それらの光照射により発生する光音響波を受信する。 制御部153は、走査情報と光照射を示す情報(制御信号)をプローブ180に送信する。駆動部130が受信部120を移動させながら、光照射部110は複数回にわたって複数の波長のパルス光を被検体100に照射する。すなわち、複数回の光照射が行われる期間に、駆動部130は受信部120を移動させる。その結果、駆動部130は、受信部120が各光照射時に互いに異なる位置に位置するように受信部120を移動させることができる。トランスデューサ121は、光照射部110が複数回パルス光を照射することにより発生した光音響波を受信することにより、光照射回数分だけ信号を出力する。以下、複数波長による複数回の光照射回数分出力された信号を総称して、複数波長に対応する信号群と呼ぶ。
以下、N回の光照射を行う場合を説明する。なお、波長λ1で、i回目の光照射により得られた信号を
Figure 0006759032
と表記する。添字iが付いたアイテムは、i回目の光照射に対応するアイテムであることを表す。iは正の整数であり、パルスインデックスとも呼ぶ。
例えば、図5に示すように、プローブ180を移動させながら光照射及び光音響波の受信を行うことができる。符号501で示された丸は、波長λ1の光を照射したときのプローブ180の位置(測定位置)を表す。符号502で示された三角は、波長λ2の光を照射したときのプローブ180の位置(測定位置)を表す。符号503で示された実線は、プローブ180の軌跡を表す。図5に示されたように、本実施形態では、光照射部110が複数の波長の光を交互に照射しながら、駆動部130がプローブ180を動かしてもよい。図5に示す例では最外周の軌跡から内側に向かって渦巻状にプローブ180を走査している。図5では便宜上XY平面で描写しているが、平面にかぎらず立体的にプローブを走査してもよい。
信号収集部140は、トランスデューサ121から出力されたアナログ信号群である複数波長に対応する信号群に対してAD変換処理等を行い、処理後の光音響信号をコンピュータ150に送信する。デジタル信号群としての光音響信号は、記憶部152に記憶される。
(S200:波長λ1に対応する複数の画像データを取得する工程)
演算部151は、S100で取得された信号群に基づいて、第1の波長の光照射に対応する複数の画像データを取得する。演算部151は、各光照射の信号のそれぞれから画像データを生成し、その中から波長λ1の光照射により得られた画像データを取得してもよい。また、演算部151は、信号群の中から波長λ1の光照射に対応する信号を選択的に用いて、波長λ1の光照射により得られた画像データを取得してもよい。
演算部151は、光音響信号に対してUniversal Back−Projection(UBP)等の再構成処理を行うことにより、光音響画像を生成する。なお、光音響画像を生成したところで、記憶部152に保存された光音響信号を削除してもよい。1回のパルス光の照射によって得られる画像データをパルスボリュームデータとも呼ぶ。パルスボリュームデータは、2次元または3次元に配列したボクセル(2次元の場合はピクセルとも呼ぶ)のそれぞれに当該位置における値を収めたボリュームデータの形式で取得される。ボリュームデータは、2次元または3次元ボリューム、2次元または3次元画像、2次元または3次元断層像とも呼べる。
再構成手法については、タイムドメイン再構成手法、フーリエドメイン再構成手法、モデルベース再構成手法(繰り返し再構成手法)などの高知の再構成手法を採用することができる。例えば、PHYSICAL REVIEW E71,016706(2005)に記載されたようなUniversal Back−Projection(UBP)と呼ばれるタイムドメイン再構成手法を採用することができる。
演算部151は、画像データとして初期音圧分布データ
Figure 0006759032
を取得してもよい。aは、波長λaの光照射に対応するアイテムであることを示す波長インデックスである。演算部151は、信号群に加え、光照射時のトランスデューサの位置情報に基づいて、初期音圧分布データを取得してもよい。演算部151は、記憶部152に予め記憶された各光照射時のトランスデューサの位置情報を読み出すことにより位置情報を取得することができる。また、演算部151は、光照射をトリガーとして、駆動部130に備えられた位置センサからの受信部120の位置情報を受け取ることにより、トランスデューサの位置情報を取得してもよい。
図6は、本実施形態に係る波長λ1でのパルスボリュームデータの一部(Pλ1,1〜Pλ1,10)を示す。本実施形態におけるパルスボリュームデータは3次元空間中のボリュームデータであるが、紙面の説明の都合からパルスボリュームデータをXY面で表す。本実施形態においては、時間的に隣接する初期音圧分布データの少なくとも一部の領域が重畳するように、再構成領域が設定される。本実施形態では、半球状の支持体122の曲率中心を中心とした60mm角の立方体領域を、1回の光照射、すなわち1つの電気信号群に基づいて再構成される再構成領域とする。この場合、1回の光照射による再構成領域の大きさ(60mm)は、光照射間の受信部120の移動量と比べて大きい。そのため、図6に示すように、時間的に一続きの光照射に対応する2つ以上のパルスボリュームデータが重畳することとなる。再構成領域の大きさや形状は、予め設定されていてもよい。また、ユーザが入力部170を用いて再構成領域の大きさや形状を指定してもよい。基準位置に対する各パルスボリュームデータの中心位置を、各パルスボリュームデータの位置とする。図6では、一例としてパルスボリュームデータPλ1,10の位置PosPλ1,10を示した。この図の例では光照射ごとに受信部120の位置が異なるため、図6に示された本実施形態で得られる各パルスボリュームデータは、基準位置に対して互いに異なる位置に位置する。
なお、本工程において、演算部151は、被検体内での光フルエンス分布データΦ[Pa・m3/J]を取得してもよい。そして、演算部151は、初期音圧分布データを、光フルエンス分布データとグリュナイゼン係数分布データとで除算することにより、被検体内の光吸収係数分布データμa[1/m]を取得してもよい。この場合、光吸収係数分布データをパルスボリュームデータとしてもよい。
例えば、演算部151は、Proc.of SPIE Vol.7561 756117−1に記載されたように、光拡散方程式を解くことにより光フルエンス分布データを取得してもよい。
また、例えば、グリュナイゼン係数は被検体の種類が決定するとほぼ一意に値が決定されることが知られているため、被検体に対応するグリュナイゼン係数分布データΓを予め記憶部152に記憶しておくことができる。そして、演算部151は、予め記憶部152に記憶されたグリュナイゼン係数分布データΓを読み出すことにより、取得してもよい。
なお、ユーザが、把持部を有するプローブ180を把持し、プローブ180を移動させてもよい。また、複数回の光照射を行っている期間にプローブ180を移動させなくてもよい。また、演算部151は、1回の光照射により得られる電気信号群に基づいて全画像化領域の画像データを取得し、それを複数回の光照射に繰り返してもよい。
(S300:波長λ1の光の照射タイミングにおける位置ずれ情報を取得する工程)
演算部151は、S200で取得された波長λ1に対応する複数の画像データに基づいて、波長λ1の光の照射タイミングにおける位置ずれ情報を取得する。例えば、以下説明するように、パルスボリュームデータ間の位置ずれ量を算出することや、パルスボリュームデータが合成されたボリュームデータの位置ずれ量を算出することにより、当該位置ずれ情報を取得することができる。
まず、演算部151は、波長λ1の光照射間における、対象物と受信部との相対的な位置関係の変動による、各パルスボリュームデータの位置ずれ量を推定する。その際、得られたパルスボリュームデータから、任意のパルスボリュームデータのペアを選択する。k個目のペアをR_kと表記する。また、ペアR_kを構成する波長λ1のパルスボリュームデータの一方をPλ1,k1、もう一方をPλ1,k2と表記する。以下、本実施形態では、K個のペアが選択された場合について説明する。なお、オーバーラップ領域を有する2つのパルスボリュームデータをペアとすることが好ましい。これにより、共通の特徴を有さないパルスボリュームデータ同士を比較することを避けることができるため、冗長な計算を減らすことができる。さらに、オーバーラップ領域の大きいパルスボリュームデータ同士をペアとすることが好ましい。そこで、例えば、演算部151は、パルスボリュームデータ間のオーバーラップ領域の体積の割合が所定の値以上のペアを選択してもよい。また、後述するように、位置ずれ補正に合成ボリュームデータを用いる場合には、パルスボリュームデータの重畳数が多い領域同士が重なるようなペアを選択してもよい。
また、あるパルスボリュームデータに対して、当該パルスボリュームデータのインデックスから、インデックスが所定の範囲に含まれるパルスボリュームデータをペアの対象として選択してもよい。また、インデックスが連続する、すなわち時間的に連続するパルスボリュームデータをペアの対象として選択してもよい。例えば、本実施形態では、演算部151は、オーバーラップ領域が50%以上となるパルスボリュームデータを対象としてペアを選択する。
以下、各パルスボリュームデータの位置ずれ量を推定する方法の例を説明する。
演算部151は、式1に示すように、Pλ1,k1とPλ1,k2との間の類似度関数F_kを取得する。
F_k(x,y,z)=fsimil(P_k,x,y,z)・・・式1
ここで、類似度関数F_kは、ペアR_kを構成する片方のパルスボリュームデータPλ1,k1に対する、もう一方のパルスボリュームデータPλ1,k2の相対位置を(x,y,z)だけ並進させた場合の類似度を算出する関数である。ここで関数fsimilは、画像間の類似度が高い場合には関数値として高い値を返すものとする。類似度関数Fの取得とは、各関数の引数である並進量(x,y,z)、すなわち画像データ間の相対位置を所定の範囲内で離散的に変化させた場合の関数値の取得を意味する。例えば、x,y,zのそれぞれの値を−Lから+Lまでの整数値として変化させた場合のそれぞれについてFが返す(2L+1)×(2L+1)×(2L+1)個の値の集合の取得を意味する。より発展的には、(2L+1)×(2L+1)×(2L+1)個の値の集合をさらにバイリニア法やバイキュービック法などを用いて、より連続関数に近い情報として類似度関数Fを導出し、これを取得するものとしてよい。
なお、Pλ1,k1に対するPλ1,k2の相対位置(光照射間の受信部120の移動量)だけ並進させた位置を基準として、Pλ1,k2の位置を所定の範囲内で離散的に変化させた場合の関数値を取得してもよい。
例えば、類似度を算出する関数としては、SSD(Sum of Squared Difference)やSAD(Sum of Absolute Difference)、相互情報量、相互相関など、任意の類似度尺度が適用できる。また、例えば、パルスボリュームデータから特徴的な形態を抽出し、それらの位置との一致度を測ることによって類似度関数を取得してもよい。抽出する特徴としては、血管などの解剖学的な特徴、エッジ検出やコーナー検出などの画像処理分野で一般的に用いられる公知の技術によって抽出された特徴を用いてもよい。また、抽出する特徴としては、コンピュータビジョン等の技術分野で一般的に使われるSIFT特徴やSURF特徴などのより高次の局所画像特徴等を用いてもよい。これらの方法によれば、パルスボリュームデータ間の輝度分布の相違やノイズの混入によって、より頑健な類似度関数を取得することができると考えられる。
なお、演算部151は、類似度の計算の結果に、重みを乗ずる処理を施すことにより類似度関数を取得してもよい。
また、類似度算出の対象となるパルスボリュームデータ間で正しく類似度を算出することができなかった場合、その結果を以降の処理には使用しなくてもよい。正しく類似度を算出することができない場合としては、いずれに並進させても類似度が十分小さい、または変わらない場合などが考えられる。この処理によれば、同一の特徴が十分に現れているパルスボリュームデータ同士の比較結果(類似度関数)を選択的に、以降の処理に使用することができる。
続いて、演算部151は、式2に示すように、類似度関数F_kの関数値が最大となる並進量M_kの取得を行う。
M_k=arg max(F_k(x,y,z))・・・式2
演算部151は、各ペアについて類似度関数F_kの関数値が最大となる並進量M_kの取得を行う。
パルスボリュームデータの位置を推定する場合に、ペアR_kに対する個別最適値である並進量M_kをなるべく保つような評価関数を定義する。すなわち、Pλ1,k1に対するPλ1,k2の位置が並進量M_kから離れるにつれて、値が低下する評価関数を定義する。式3は、この場合の評価関数E_kの例を表す。
E_k=(M_k−(PosPλ1,k1−PosPλ1,k2
=(M_k(x)−(PosPλ1,k1(x)−PosPλ1,k2(x)
+(M_k(y)−(PosPλ1,k1(y)−PosPλ1,k2(y)
+(M_k(z)−(PosPλ1,k1(z)−PosPλ1,k2(z)
・・・式3
PosPλ1,k1は、基準位置に対するPλ1,k1の位置を表す。PosPλ1,k2は、基準位置に対するPλ1,k2の位置を表す。なお、評価関数を定義する際に、類似度関数F_kを、当該類似度関数F_kにフィットするような二次関数に近似してもよい。また、類似度関数F_kが、並進量M_kの周辺において、二次関数に従って低下すると禁じできる場合には、式3はPλ1,k1とPλ1,k2の位置関係から類似度関数F_kの値を並進量M_kの周辺で近似する関数となる。
続いて、演算部151は、式4のように定義されたコスト関数Eが最小化した時の、基準位置に対する全てのパルスボリュームデータの位置PosP′λ1,jを取得する。ここで、jはあるパルスに対するパルスインデックスである。
Figure 0006759032
コスト関数が最小化したときの、基準位置に対するパルスボリュームデータの位置は、被検体100と受信部120との相対的な位置関係の変動による位置ずれ後のパルスボリュームデータの位置情報を表す。
例えば、演算部151は、式4に示すコスト関数Eを最小化する(0に最近接する)解を線形最小二乗法を解くことにより得る。これにより、一意に各パルスボリュームデータの位置PosP′λ1,jを求めることができるため、計算コストが小さい。
なお、上記で説明した線形最適化によるコスト関数の最適化に限らず、コスト関数の最適化は公知のいかなる方法を用いてもよい。例えば、最急降下法、ニュートン法のような繰り返し計算による非線形最適化の方法などによって最適化してもよい。すなわち、演算部151は、コスト関数が最小化するような各パルスボリュームデータの位置を探索することにより、基準位置に対するパルスボリュームデータの位置ずれ後の位置情報を取得する。
なお、コスト関数は、想定される各パルスボリュームデータの位置の光照射間の変動(動き)に対して正則化をかけるように定義してもよい。被検体として***を考えた場合、呼吸による動きが支配的であると想定される。この場合、被検体の動きは最大で数mm程度の動きであり、その動きは時間的に連続で滑らかなものであることが想定される。また、その動きは周期的な動きとなることが想定される。上記のように想定される被検体の動きから逸脱するような動きが算出されることに対して抑制を働かせるような正則化をかけることができる。
正則化のかけ方はいかなる方法であってもよい。例えば、導出過程の被検体の変動量(移動距離)の総和に所定の重み係数をかけてコスト関数に加算することで、正則化することができる。また、被検体の変動の周波数成分値に基づいて算出された値をコスト関数に加算してもよい。また、被検体の典型的な変動のしかたをモデルとして用意し、そのモデルにおける変動との相違をコストとしてコスト関数に加算するようにしてもよい。
また、「コスト関数を最小化させる」とは、コスト関数が厳密に最小となる場合だけではなく、解の候補を変化させたときにコスト関数の値が所定の値以下となる場合やコスト関数の変化量が所定の値以下となる場合も含む。すなわち、演算部151は、コスト関数が所定の条件を満たすことをもって、コスト関数が最小化したと判断してもよい。また、ユーザが入力部170を用いて、コスト関数が最小化したことを指示してもよい。この場合、演算部151は、入力部170からの指示を受けてコスト関数が最小化したと判断する。
続いて、演算部151は、各パルスボリュームデータに対する、コスト関数が最小化したときの位置ずれ量Moptλ1,jを取得する。この位置ずれ量Moptλ1,jは、被検体100と受信部120との相対的な位置関係の変動による、各パルスボリュームデータの位置ずれ量を表す。
図7は、パルスボリュームデータPλ1,2の位置PosPλ1,2と、前述の方法でコスト関数が最小化したときのパルスボリュームデータP′λ1,2の位置PosP′λ1,2とを示す。図7において、パルスボリュームデータPλ1,2を実線で表し、コスト関数が最小化したときのパルスボリュームデータP′λ1,2を破線で表す。
なお、本工程においては、被検体100と受信部120との相対位置の変動による、パルスボリュームデータの位置ずれ情報を取得することができる限り、いかなる手法を用いてもよい。
なお、演算部151は、パルスボリュームデータを2つ以上合成することで得られる波長λ1でのk番目の合成ボリュームデータGλ1,kを用いて位置ずれ情報を取得してもよい。その際、合成したパルスボリュームデータをすべて包含し、かつ、最小となる矩形の領域を合成ボリュームデータとしたり、少なくとも2つ以上のパルスボリュームデータが重畳した領域を包含する任意の領域を合成ボリュームデータとしてもよい。すなわち、合成ボリュームデータの領域は、合成対象のパルスボリュームデータを全て包含しなくてもよい。
合成ボリュームデータは、選択されたパルスボリュームデータを位置に従って加算することで取得する。その他、選択されたパルスボリュームデータを加算し、重なるパルスボリュームの数で除算することで平均化してもよいし、被検体の特徴をより正確に再現したボリュームデータを取得できるかぎり、あらゆる手法を用いることができる。ただし、光照射間の被検体100と受信部120との相対位置の変動を補正する処理(例えば、パルスボリュームの位置を変更する処理)については、本明細書に係る「合成」に含まれない。
例えば、演算部151は、合成対象のパルスボリュームデータのそれぞれを重み付けした後に加算することにより合成してもよい。また、演算部151は、外れ値除去法などによってノイズを多く含む値を除外したパルスボリュームデータに対して加算値や平均値を算出するようにしてもよい。
これらの合成処理により、各パルスボリュームデータに含まれるノイズが低減され、被検体の特徴をより正確に再現した合成ボリュームデータを取得することができる。さらに、パルスボリュームデータの数よりも合成ボリュームデータの数のほうが少ない場合は、パルスボリュームデータ同士を比較することにより全パルスボリュームデータの位置ずれを推定する手法に比べて、計算量と計算コストを少なくすることもできる。
ただし、このように合成して得られる合成ボリュームデータには、複数回の光照射間の被検体と光音響波の受信部との相対的な位置関係の変動による影響が含まれている。そのため、合成ボリュームデータには当該変動による品質の低下が生じている可能性がある。以下、この品質の低下を抑制するために、合成ボリュームデータの推定位置からパルスボリュームデータの位置を推定する処理を説明する。
演算部151は、波長λ1での各合成ボリュームデータの位置ずれ量MGoptλ1,jに基づいて、各パルスボリュームデータの位置ずれ量Moptλ1,iを推定する。
演算部151は、合成ボリュームデータに対応付けられたパルスボリュームデータの位置ずれ量については、推定された合成ボリュームデータの位置ずれ量割り当てることができる。その他のパルスボリュームデータの位置ずれ量については、演算部151が、割り当てられたパルスボリュームデータの位置ずれ量に対して補間処理を行うことにより推定することができる。補間処理の手法について、線形補間やスプライン補間などの公知の手法を採用することができる。また、想定される被検体の動きから逸脱するような位置を算出しないような制約をかけて補間処理を行ってもよい。
合成の対象となったパルスボリュームデータのうち、任意のパルスボリュームデータを、合成ボリュームデータに対応付けられたパルスボリュームデータとしてもよい。例えば、合成対象のパルスボリュームデータが奇数個である場合、時間的に中心に位置するパルスボリュームデータを合成ボリュームデータに対応付けてもよい。
また、例えば、本実施形態のように合成対象のパルスボリュームデータが偶数個である場合、時間的に中心付近に位置するいずれかのパルスボリュームデータを合成ボリュームデータに対応付けてもよい。例えば、本実施形態のように10個のパルスボリュームデータを合成対象とする場合、波長λ1での合成ボリュームデータGλ1,jの位置ずれ量MGoptλ1,jを、パルスボリュームデータPλ1,5jの位置ずれ量MGoptλ1,5jとして割り当ててもよい。
また、合成対象のパルスボリュームデータが偶数個である場合、時間的に中心に位置する仮想のパルスボリュームデータを合成ボリュームデータに対応付けてもよい。例えば、本実施形態のように10個のパルスボリュームデータを合成対象とする場合、波長λ1での合成ボリュームデータGλ1,jの位置ずれ量を、パルスインデックスが5.5jの仮想パルスボリュームデータの位置ずれ量に割り当ててもよい。
また、重み付けを伴って合成される場合、合成対象のパルスボリュームデータのうち、最も高い重み係数で重み付けられたパルスボリュームデータを、合成ボリュームデータに対応付けてもよい。また、合成対象のパルスボリュームデータのうち、重み係数が中央値となるパルスボリュームデータを、合成ボリュームデータに対応付けてもよい。
以上の処理により、波長λ1に対応する合成ボリュームデータの位置ずれ情報に基づいて、波長λ1に対応するパルスボリュームデータの位置ずれ情報を取得することができる。このようにして、波長λ1の光の照射タイミングにおける位置ずれ情報を取得してもよい。
なお、パルスボリュームデータや合成ボリュームデータの最大値投影画像(Maximum Intensity Projection:MIP)等の2次元の投影データを用いて上記と同様に位置ずれ情報を取得してもよい。以下、その処理の一例を説明する。
演算部151は、波長λ1での合成された初期音圧分布データGλ1,jについて、X方向、Y方向、Z方向のそれぞれの方向に投影した投影データとしてMIPデータを取得する。X方向に投影したMIPデータはY軸とZ軸によって表される2次元の空間情報であり、IGλ1,j(y,z)と表記される。Y方向に投影したMIPデータはZ軸とX軸によって表される2次元の空間分布情報であり、IGλ1,j(z,x)と表記される。Z方向に投影したMIPデータはX軸とY軸によって表される2次元の空間分布情報であり、IGλ1,j(x,y)と表記する。
なお、3次元の画像データを2次元の画像データへと変換できるかぎり、MIP画像以外の投影手法を採用してもよい。例えば、MIP画像にかえて最小値投影(Minimum Intensity Projection:MINP)画像を生成して用いるようにしてもよい。また、投影方向の複数のスライドを加算することにより、投影データを取得してもよい。
続いて、XY面、YZ面、ZX面のそれぞれについて、ペアを構成する片方のMIPデータに対するもう一方のMIPデータの相対位置を並進させ、類似度を算出する。類似度算出の手法については、前述した手法を用いることができる。そして、XY面、YZ面、ZX面のそれぞれについて、Pλ1,k1に対するPλ1,k2の類似度が最大となる並進量MX_k、MY_k、MZ_kと、それぞれの並進量の各座標軸の成分の平均値を、類似度が最大となる三次元の並進量Mの各成分値として算出する。
Figure 0006759032
続いて、演算部151は、式5に示す並進量M_kを用いて、式4に示すコスト関数が最小化したときの各合成ボリュームデータの位置を推定することができる。
以上の処理により、3次元の画像データから変換された2次元の画像データに基づいて、基準位置に対する合成ボリュームデータの位置を取得することができる。3次元画像データから2次元画像データに変換することにより、3次元画像データのまま処理を行う場合と比べて少ない計算コストで、位置ずれ後のボリュームデータの位置を取得することができる。
なお、パルスボリュームとしての3次元画像データを2次元の画像データに変換することにより、上記の方法で波長λ1に対応するパルスボリュームデータの位置ずれ情報を取得してもよい。このようにして、波長λ1の光の照射タイミングにおける位置ずれ情報を取得してもよい。
ここまでは、被検体と受信部との相対的な位置関係の変動として並進が生じる場合の例を説明した。ただし、当該変動として回転や変形が生じる場合についても同様にそれらによる位置ずれ量を推定することができる。
例えば、回転を考慮する場合、演算部151は、並進量に加えて回転量を引数として、各パルスボリュームデータの位置および回転量(位置ずれ量)を推定することができる。続いて、演算部151は、推定された位置および回転量に基づいて各パルスボリュームデータを剛体変換処理(位置ずれ補正処理)した後に合成することにより、合成ボリュームデータを取得することができる。なお、回転量のみを位置ずれ量としてもよいし、算出した2次元もしくは3次元の並進・回転行列等の変換行列や変換のための各種パラメータを位置ずれ量としてもよい。
また、例えば、変形を考慮する際、演算部151は、パルスボリュームデータに設定された各点での変位量(並進および回転量の少なくとも一つ)を引数として変位量を推定することができる。続いて、演算部151は、推定された変位量に基いて各パルスボリュームデータを変形処理(位置ずれ補正処理)した後に合成することにより、合成ボリュームデータを取得することができる。例えば、Free Form Deformation(FFD)やThin Plate Splineなどの変形を表現する手法によって、パルスボリュームデータ間の変位量を算出することができる。これらの処理により、変形を含む高次な変動を考慮して、品質の高い合成ボリュームデータを取得することができる。
本実施形態では全光照射による光音響波の測定が完了した後にパルスボリュームデータの取得を開始する例を説明したが、光照射のたびに逐次パルスボリュームデータを取得してもよい。後者の場合、取得されたパルスボリュームデータを逐次表示部160に表示させてもよい。これにより、ユーザは全測定が完了する前に取得済のパルスボリュームデータを確認することができる。このとき、パルスボリュームデータが重畳した領域については、前述の方法で合成してもよい。
また、本実施形態では全てのパルスボリュームデータを取得した後に位置ずれ量を算出したが、光照射のたびに得られたパルスボリュームデータを用いて逐次位置ずれ量を算出してもよい。また、逐次得られるパルスボリュームデータから一定の数のパルスボリュームデータを合成し合成ボリュームデータを作成して逐次位置ずれ量を算出してもよい。
また、位置ずれ量を算出するためのボリュームデータは、負の値を除去したり、画像強度を正規化したり等の事前処理を行ってから使用してもよい。
(S400:波長λ1の光の照射タイミングにおける位置ずれ情報に基づいて、波長λ2の光の照射タイミングにおける位置ずれ情報を取得する工程)
演算部151は、S300で取得された波長λ1の光の照射タイミングにおける位置ずれ情報に基づいて、波長λ2の光の照射タイミングにおける位置ずれ情報を取得する。すなわち、S300に記載されたような方法で算出された特定の波長に対応する位置ずれ情報から、その他の波長に対応する位置ずれ情報を算出する。例えば、光の照射タイミングが既知である場合、λ1の光の照射タイミングにおける位置ずれ情報を時間的に補間してもよい。また、光の照射タイミングにおけるプローブ180の位置、すなわち受信部120の位置が既知である場合、波長λ1の光の照射タイミングにおける位置ずれ情報を空間的に補間してもよい。演算部151は、これらの補間により、波長λ2の光の照射タイミングにおける位置ずれ情報を取得してもよい。
時間的に補間して算出する場合に関して、図8を用いて説明する。図8は、時間tにおける複数の波長のパルス光の照射タイミングを表すシーケンス図を含む。また、図8は、波長λ1での各パルス光に対応する位置ずれ量(並進量)を表すグラフを含む。符号801の逆三角形は、パルス光の照射タイミングを模式的に表したものである。波長λ1のパルス光は、T11、T21、T31、T41、及びT51の時間にそれぞれ照射されていることが理解される。また、波長λ2のパルス光は、T12、T22、T32、T42、及びT52の時間にそれぞれ照射されることが理解される。また、S300において波長λ1のパルス光の照射タイミングにおける位置ずれ量が算出されているため、その位置ずれ量が図8に丸印802でプロットされている。そこで、演算部151は、波長λ2に対応する位置ずれ量は波長λ1での位置ずれ量から、時間的に補間して算出する。図8においてバツ印803は、時間的に補間して算出された波長λ1のパルス光の照射タイミングにおける位置ずれ量を示す。補間方法はランツォシュ補間を用いるのが理想的に正しいが、その他線形補間、キュービック補間、スプライン補間等、どのような補間方法でもよい。ここでは例として、線形補間、ランツォシュ補間を用いた場合で説明する。例えば、波長λ2でのT12のタイミングでの位置ずれ量は、観測している被検体の動きが波長λ1のパルス照射時間(サンプリング時間)よりも遅い場合、波長λ1に対応するT11とT21での位置ずれ量から算出できる。T11、T21での位置ずれ量をそれぞれM11、M21とした場合、T12での位置ずれ量M12は、下記式6で算出してもよい。
Figure 0006759032
ここで波長λ2の時刻T12のパルスの位置ずれ量M12を算出するのに、M11とM21を用いて算出を行った。しかし、補間方法によっては、さらに時間的に離れたパルスの位置ずれ量を用いてもよい。例えば、ランツォシュ補間は多変数補間の一種で、デジタル信号のサンプリングレートを向上するのに用いられ、最もよく補間を行えることが知られている。1次元の信号に対するランツォシュ補間は、式2のランツォシュカーネルL(t)を用い、式7によって行う。
Figure 0006759032
Figure 0006759032
ここでtは時刻、aはカーネルサイズを決定する正の整数、siは正の整数iに対する1次元信号のサンプル、S(t)は補間された値、[ ]は床関数である。線形補間の場合と異なり、最近傍の2つのパルスの位置ずれ量のみならず、近傍にある複数のパルスの位置ずれ量を用いて補間を行っていることがわかる。
このようにすることで、1つの波長の時間的に近傍にあるパルスの位置ずれ量から、他の波長のパルスの位置ずれ量を算出することができる。また、その過程で、複数の波長の位置ずれ補正を行う必要がなくなるため、計算時間も短縮することができる。
時間的に補間を行う方法の例として、2つの2次元のAffine変換行列から補間して1つの2次元のAffine変換行列を算出する方法を説明する。
取得した複数波長のパルス光による信号群のうち、波長λ1のパルス光を照射して取得した信号から、対応するパルスボリュームデータを取得する。得られた波長λ1のパルスボリュームデータから、Affine変換をベースにした位置ずれ補正を行い、波長λ2での位置ずれ量を算出する。波長λ2での位置ずれ量は、パルスボリュームデータに対してAffine変換行列として算出される。
まず、2つの2次元のAffine変換行列A1とA2があるとする。ここで、A1は、t1の時刻に波長λ1のパルス光を照射して得られるパルスボリュームデータPλ1,0と、t1+Δの時刻に波長λ1のパルス光を照射して得られるパルスボリュームデータPλ1,1間での位置ずれ量(Affine変換行列)である。また、A2は、t1+Δの時刻に波長λ2のパルス光を照射して得られるパルスボリュームデータPλ2,1と、t1+2Δの時刻に波長λ2のパルス光を照射して得られるパルスボリュームデータPλ2,2間での位置ずれ量である。
このとき、求めたいAffine変換行列Bは、t1+Δ/2の時刻に波長λ1のパルスを照射して得られるパルスボリュームデータPλ1,1/2と、t1+3Δ/2の時刻に波長のパルスを照射したパルスボリュームデータPλ1,3/2間での位置ずれ量であるとする。
A1はA1Rという回転成分とA1Sという拡大縮小成分に分解することができ、A2も同様にA2Rという回転成分とA2Sという拡大縮小成分に分解できる。よって、Bの回転成分BRはA1RとA2Rの行列の各要素を時間的に補間することで、Bの拡大縮小成分BSはA1SとA2Sの行列の各要素を時間的に補間することで得られる。このようにして得られたBRとBSを積算することで、位置ずれ量Bが算出できる。
また、空間的に補間を行い、空間的に近傍にあるパルスの位置ずれ量を算出してもよい。図9は、図5と同様に空間的に走査しながら複数の波長のパルスを照射するときのプローブ180の位置(測定位置)を表す模式図である。波長λ1に対応する位置ずれ情報を空間的に補間することにより、波長λ2の光照射のタイミングにおける測定位置502Aに対応する位置ずれ情報を取得する場合を説明する。演算部151は、測定位置502Aの近傍の測定位置に対応する波長λ1の測定位置を決定する。例えば、波長λ2の測定位置502Aから所定の距離に含まれる波長λ1の測定位置を決定する。符号901の一点鎖線は、測定位置502Aから所定の距離の範囲を示す。そして、波長λ1の測定位置501A、501B、501C、501D、及び501Eは、波長λ2の測定位置502Aから所定の距離に含まれる波長λ1の測定位置を示す。よって、演算部151は、波長λ1の測定位置501A、501B、501C、501D、及び501Eに対応する位置ずれ情報を空間的に補間することにより、波長λ2の測定位置502Aに対応する位置ずれ情報を取得する。なお、空間的に補間する方法としては、すでに算出されている位置ずれ情報に対して関数のフィッティングを行い、求めたい位置での位置ずれ情報を関数から求めたてもよい。また、空間的に補間する方法としては、距離に応じた重みをつけ、求めたい位置での位置ずれ情報を算出してもよい。以下では、位置ずれ情報を取得したい測定位置から所定の距離内にあるパルスボリュームデータの位置ずれ情報に対して距離に応じた重みをつけることで、算出したいパルスボリュームデータの位置ずれ情報を算出する方法について説明する。
ここではパルスを照射した位置をパルスボリュームデータの中心の座標とする。また、位置ずれ情報を算出したい波長λ2のパルスボリュームデータPresultの中心座標POSresultから、半径Rmmの範囲に含まれる波長λ1のパルスボリュームデータの位置ずれ情報のパルス数をNとする。また、位置ずれ情報としての位置ずれ量(並進量)を
Figure 0006759032
、中心座標を
Figure 0006759032
とすると、Presultの位置ずれ量Mresultは、以下の式で算出できる。
Figure 0006759032
例えば、図9においては、符号502Aが、位置ずれ情報を算出したい波長λ2のパルスボリュームデータPresultの中心座標POSresultに対応する測定位置を表す。また、符号901が、測定位置502Aから半径Rmmの範囲を表す。また、符号501A〜Eが、範囲901に含まれるパルスパルスボリュームデータの中心座標に対応する測定位置を表す。
ここでωは、距離から計算される重みであり、POSとPOSresultの間の空間的な距離
Figure 0006759032
を用いて次の式で計算できる。
Figure 0006759032
なお、空間的に補間する場合、上記方法に限らず、あらゆるパルスボリュームを用いて空間的に補間してもよい。例えば、再近傍のパルスを用いたり、渦巻状の走査が複数周回している場合は各周回で近傍2つのパルスを用いたり、全パルスを用いたりしてもよい。
このようにすることで、特定の波長の測定位置と空間的に近傍にある測定位置に対応する位置ずれ情報から、他の波長のパルスの位置ずれ情報を算出することができる。
3つ以上の波長の光を照射する場合に本実施形態に係る位置ずれ情報の取得方法を適用してもよい。この場合、特定の波長で得られた位置ずれ情報を、その他の複数波長の位置ずれ情報の取得に用いてもよい。また、位置ずれ情報の取得に好適な波長を含む複数の波長で得られた位置ずれ情報を用いて、その他の波長の位置ずれ情報を取得してもよい。すなわち、少なくとも一つ以上の波長により得られた位置ずれ情報を用いて、その他の波長の位置ずれ情報を求めてもよい。
また、FFD等の変形位置合わせ技術を用いて、変形を加味した位置ずれ情報を取得する場合、その位置ずれ情報は変形場でもよい。変形場とは、変形用画像から基準画像への変形の場である。
位置ずれ情報が変形場の場合の例を、図10を用いて説明する。ここでは、簡便のために3次元画像データ領域を、XYの2次元平面にして説明する。符号1011、符号1031は波長λ1でのパルス照射位置、符号1012、符号1032は前述のパルスの画像データ領域、符号1050の点線の領域は符号1012と符号1032の画像の重なり領域である。また、符号1021は波長λ2でのパルス照射位置、符号1022は符号1021のパルスに対応する画像データ領域である。ここで、符号1011、符号1031、符号1021はそれぞれ、時刻t1、t2、(t1+t2)/2の時刻に照射されたとする。
ここで、符号1012と符号1032の重なり領域1050において、符号1032の画像を基準として符号1012の画像を変形位置合わせを行うことを考える。すると、重なり領域符号1050のサイズの、3次元的な変形場を算出することができる。
そして、符号1022の画像データ領域のうち、符号1050の領域の変形場を、符号1032の画像を基準とした符号1012の画像の変形場から算出することができる。算出方法としては、符号1032の画像を基準として符号1012の画像の変形が、t2−t1の時間で時間の経過に対して線形に発生したものと想定して、変形場を時間で補間することで得られる。より具体的には、(t2−(t1+t2)/2)/(t2−t1)を積算すればよい。
この例では、2画像間の変形場から、その中間の時刻での画像の片方の画像からの変形場を算出した。これによると、複数波長で交互に照射されたパルスの画像データそれぞれに対して、1つの波長の中で重なり領域がある2つの画像の撮像時刻の間に撮像した、その他の波長の画像の変形場を算出できる。これを複数波長の全撮像パルスに適用することで、全ての波長でいずれかのパルスの画像データを基準とした各パルスへの変形場を算出できる。
その他の例としては、複数のパルスの画像データ間で複数の変形場が求められている際に、それら複数の変形場から時間的に補間することで、その近傍の時刻で照射したパルスの変形場を算出することができる。変形場の算出のためには、画像データ中のあるボクセルにおける変形場での変位量に着目し、離散的な時刻で撮像した画像データ間の複数の変形場において、変形場中の各ボクセルがどのように変位していくかを追跡する。そしてその経時的な変位量から、とある時刻での変形場を算出してもよい。
位置ずれ情報を取得するための波長λ1としての基準は、位置ずれ情報の取得精度の高いものを選択するのが好ましい。そこで、例えば、演算部151は、複数波長のそれぞれを用いて生成されたパルスボリュームデータや合成ボリュームデータ、MIP画像等の解像度やSN(Signal to Noise Ratio)などの画質を示す評価指標を算出する。そして、演算部151が、算出された画質を示す評価指標が良好であった波長を波長λ1として、位置ずれ情報を取得してもよい。
また、時間的に補間して位置ずれ情報を取得する場合、演算部151は、所定の期間内に多くの回数照射された光の波長を判定し、その波長を波長λ1として設定してもよい。例えば、演算部151は、所定の期間内に最も多くの回数が照射された光の波長を、当該所定の期間内においては波長λ1として設定してもよい。また、演算部151は、所定の期間内に所定の回数より多く照射された光の波長を、当該所定の期間内においては波長λ1として設定してもよい。なお、波長λ1として複数の波長が設定されてもよい。
また、波長間の照射タイミング実質的に同じと見なせる場合には、波長λ1の光の照射タイミングにおける位置ずれ情報を、波長λ2の光の照射タイミングにおける位置ずれ情報として取得してもよい。例えば、S300で得られた位置ずれ量をその他の波長のパルスにそのまま用いることができる。
例えば、複数波長間でのパルス光の照射間隔が被検体の動きに対して無視できるほど短い場合などには、波長間の照射タイミングが実質的に同じとみなせる。例えば、生体の呼吸による体動を考えた場合に、波長間の照射タイミングを実質的に同じとみなせる照射間隔を説明する。周期が3秒、最大変位量が3mmの呼吸により生体が線形に変位すると仮定し、音響画像の解像度の許容誤差を0.25mmと設定した場合、照射間隔が125ms以内であれば、波長間の照射タイミングが実質的に同じであるとみなしてもよい。なお、このとき、生体の変位量が光音響画像の解像度の誤差にそのまま反映されると仮定している。また、光音響画像の解像度の許容誤差を0.1mmと設定すると、照射間隔が50ms以内であれば、波長間の照射タイミングが実質的に同じであるとみなしてもよい。
また、演算部151は、S300またはS400で得られた位置ずれ情報を、時間的または空間的に補間することにより、その他のモダリティで得られた画像の位置ずれ情報として使用してもよい。
(S500:波長λ1及び波長λ2の光照射タイミングにおける位置ずれ情報に基づいて、位置ずれ補正を行う工程)
演算部151は、S300で得られた波長λ1に対応する位置ずれ情報、及び、S400で得られた波長λ2に対応する位置ずれ情報に基づいて位置ずれ補正を行う。
例えば、演算部151は、S300で得られた波長λ1に対応する位置ずれ情報に基づいて、S200で取得された波長λ1のパルスボリュームデータの位置情報を位置ずれ補正する処理を行ってもよい。また、演算部151は、S100で得られた信号群に基づいて、波長λ2に対応するパルスボリュームデータを生成してもよい。そして、演算部151は、S400で得られた波長λ2に対応する位置ずれ情報に基づいて、波長λ2のパルスボリュームデータの位置を位置ずれ量だけ補正する処理を行ってもよい。さらに、演算部151は、位置ずれが補正された波長λ1のパルスボリュームデータ及び波長λ2のパルスボリュームデータを合成することにより位置合わせし、合成ボリュームデータを生成してもよい。
以下、波長λ1のパルスボリュームの位置ずれを補正する例を説明する。図11は、本工程における位置ずれ補正処理(並列処理)の一例を示す。
図11Aは、本実施形態における並列処理前のパルスボリュームデータの一部(Pλ1,23〜Pλ1,25)を示す。破線は、Pλ1,25の外周およびPλ1,25内の特徴601を表す。実線は、Pλ1,24の外周および内部の特徴602表す。点線は、Pλ1,23の外周および内部の特徴603を表す。なお、特徴601、特徴602、特徴603はいずれも同一の特徴を表している。図11Aの状態では、各パルスボリュームデータ内の特徴は異なる位置に位置している。
図11Bは、並進前のパルスボリュームデータPλ1,25を、前述の方法で推定された位置ずれ量Moptλ1,25だけ並進させた後のパルスボリュームデータP′λ1,25を示す。図11Cは、並進前のパルスボリュームデータPλ1,24を、前述の方法で推定された位置ずれ量Moptλ1,24だけ並進させた後のパルスボリュームデータP′λ1,24を示す。
図11Dは、並進前のパルスボリュームデータPλ1,23を、前述の方法で推定された位置ずれ量Moptλ1,23だけ並進させた後のパルスボリュームデータP′λ1,23を示す。
図11Eは、並進後のパルスボリュームデータP′λ1,23、Pλ1,24、およびP′λ1,25を重ねあわせた様子を示す。図11Eにおいては、各パルスボリュームデータ内の特徴601、602、および603はほぼ同じ位置で重なっている。演算部151は、図11Eに示すように並進処理された各パルスボリュームデータを合成することにより、位置合わせされたボリュームデータを取得することができる。位置合わせされたボリュームデータは、合成されたボリュームデータに相当する。「位置合わせ」とは、位置ずれ補正処理と合成処理の両処理を行うことを指す。
また、例えば、演算部151は、再構成処理において、S300で得られた波長λ1に対応する位置ずれ情報に基づいて、波長λ1の照射タイミングにおける受信部120の位置情報を位置ずれ量だけ補正する処理を行ってもよい。また、演算部151は、再構成処理において、S400で得られた波長λ2に対応する位置ずれ情報に基づいて、波長λ2の照射タイミングにおける受信部120の位置情報を位置ずれ量だけ補正する処理を行ってもよい。そして、演算部151は、S100で得られた信号群、及び、波長λ1及び波長λ2の光の照射タイミングにおける位置ずれが補正された受信部120の位置情報に基づいて再構成処理を行い、合成ボリュームデータを生成してもよい。
ここで、受信部120の位置情報の位置ずれ補正を行う方法の例を説明する。まず、演算部151は、位置ずれを考慮しないときの受信部120の位置情報を取得する。例えば、演算部151は、記憶部152に予め記憶された光照射時の受信部120の位置情報を読み出すことにより、位置ずれを考慮しないときの位置情報を取得してもよい。また、演算部151は、光照射をトリガーとして、駆動部130に備えられた位置センサから受信部120の位置情報を受け取ることにより、位置ずれを考慮しないときの受信部120の位置情報を取得してもよい。
続いて、演算部151は、光照射時における位置ずれを考慮しないときの受信部120の位置情報を、S300またはS400で得られた位置ずれ情報が示す位置ずれ量だけ補正(例えば、並進処理)する。これにより、演算部151は、光照射時における位置ずれ補正された受信部120の位置情報を取得することができる。すなわち、演算部151は、S300またはS400で得られた位置ずれ情報に基づいて、位置ずれ補正された受信部120の位置情報を取得することができる。
演算部151は、S100で得られた信号群と、位置ずれ補正された受信部120の位置情報とに基づいて、合成画像データを取得する。本工程では、演算部151が、1回の光照射に対応する信号からは全画像化領域よりも小さい領域を再構成し、これを複数の光照射に繰り返すことにより1つのボリュームデータを生成してもよい。この場合、本工程において、演算部151は、複数の光照射に対応する複数のパルスボリュームデータを取得し、当該複数の画像データを合成する。また、演算部151は、複数の光照射に対応する信号群から全画像化領域を再構成することにより、1つのボリュームデータを生成してもよい。
なお、複数の波長に対応するデータの合成処理を行う場合にも、上記と同様の処理を行うことができる。また、複数波長に対応するデータの合成処理には、複数波長に対応するデータの加算処理や平均化処理の他に、複数波長に対応するデータを比較演算することにより酸素飽和度などの機能情報を取得する処理も含まれる。以下、複数波長に対応するデータを用いて機能情報としての酸素飽和度を取得する方法を説明する。
波長λ1と波長λ2では、ヘモグロビン以外の光吸収が無視できるほど低いと仮定すると、波長λ1と波長λ2の吸収係数はそれぞれ、オキシヘモグロビンのモル吸光係数とデオキシヘモグロビンのモル吸収係数を用いて式11、式12のように表わされる。
μ(λ)=εox(λ)Cox+εde(λ)Cde・・・式11
μ(λ)=εox(λ)Cox+εde(λ)Cde・・・式12
ここで、μ(λ)は、位置(i、j、k)における波長λ1の光の吸収係数、μ(λ)は位置(i、j、k)における波長λ2の光の吸収係数を示し、単位は[mm−1]で示すことができる。Coxはオキシヘモグロビンの量[mol]、Cdeはデオキシヘモグロビンの量[mol]である。いずれも位置(i、j、k)における値を示すものとする。
εox(λ)とεde(λ)はそれぞれ波長λ1におけるオキシヘモグロビン、デオキシヘモグロビンのモル吸収係数[mm−1mol−1]を示す。εox(λ)とεde(λ)はそれぞれ波長λ2におけるオキシヘモグロビン、デオキシヘモグロビンのモル吸収係数[mm−1mol−1]を示す。εox(λ)、εde(λ)、εox(λ)、εde(λ)はあらかじめ測定や文献値によって得ることができる。
よって、CoxとCdeはそれぞれ、モル吸光係数と、μ(λ)及びμ(λ)と、を用いて式11、12の連立方程式を解くことにより求められる。用いる波長の数が多い場合は、最小二乗法を用いるとよい。また、酸素飽和度SOは式13に示すように、全ヘモグロビン中のオキシヘモグロビンの割合で定義される。よって、酸素飽和度SOは、式11、12、13に基づき、式14で示すことができる。
よって、演算部151は、式14にしたがって、モル吸光係数と、μ(λ)及びμ(λ)とに基づき、位置(i、j、k)における酸素飽和度SOを得ることができる。
Figure 0006759032
Figure 0006759032
このような処理を複数位置に対して行うことにより、複数位置での酸素飽和度が求められ、酸素飽和度分布を取得できる。酸素飽和度分布は、吸収係数分布の比較演算(例えば、非を算出する処理)で得られるものであり、複数の波長での吸収係数の値が相対的に正しければ、酸素飽和度分布は適切に求められる。よって、吸収係数分布が絶対値として正確に求まっていなくてもよい。
ここでは、位置ずれの影響が低減した初期音圧画像を算出し、そこから吸収係数画像、酸素飽和度画像を算出した。なお、その方法に限らず、複数の波長でのパルスボリュームデータを吸収係数分布画像として求め、そのうちの一つの波長での吸収係数分布画像で位置ずれ量を算出してもよい。また、そのようにして算出された位置ずれ量をその他の波長に適用することで、位置ずれの影響が低減した吸収係数画像を算出し、それらの吸収係数画像から酸素飽和度画像を算出してもよい。
このようにすることで、複数の波長での画像間で、位置ずれの影響が低減した吸収係数画像を取得し、さらに位置ずれの影響が低減された吸収係数画像を合成することにより、酸素飽和度画像を算出することが可能となる。
本工程では、光照射間の被検体と光音響波の受信部との相対的な位置関係の変動の影響が抑制された画像データ(位置合わせされたボリュームデータ)を取得することができる。
なお、光照射間で受信部120が移動しない場合であっても上記処理を適用することができる。すなわち、光音響装置が駆動部130を備えない場合であっても上記処理を適用することができる。この場合も、複数回の光照射間の被検体と光音響波の受信部との相対的な位置関係の変動の影響が抑制された画像データを取得することができる。
また、光音響装置とは異なるモダリティ(例えば、超音波診断装置など)で、光音響装置を用いた撮像と併せて撮像を行う場合、光音響装置で得られた位置ずれ情報を他のモダリティにおける位置ずれ情報の取得に用いてもよい。例えば、演算部151は、S300、S400、またはS500で得られた位置ずれ情報を時間的または空間的に補間することにより、他のモダリティで得られた画像データの位置ずれ情報を取得してもよい。また、光音響装置と他のモダリティとが実質的に同じ時間に実質的に同じ位置の画像データを生成する場合、光音響装置で得られた位置ずれ情報を他のモダリティでの位置ずれ情報としてもよい。
例えば、他のモダリティとして超音波診断装置を想定する場合、光音響装置を含む検査システムは、超音波の送受信部を備える。送受信部は、被検体に対して超音波の送信し、送信された超音波のエコー波を受信することにより超音波信号を出力する。送受信部は、音響波を受信することにより電気信号を出力するトランスデューサを含む。送受信部は複数のトランスデューサを含んでいてもよい。なお、送受信部は、超音波を送信するトランスデューサと、音響波を受信するためのトランスデューサとを別に用意してもよい。また、超音波を送信するトランスデューサと、音響波を受信するためのトランスデューサとが、同じトランスデューサで構成されていてもよい。また、超音波を送受信するためのトランスデューサと、光音響波を受信するためのトランスデューサとを別に用意してもよい。また、超音波を送受信するトランスデューサと光音響波を受信するトランスデューサとが、同じトランスデューサで構成されていてもよい。
[第2の実施形態]
以下、第2の実施形態の光音響装置の構成の構成及び処理について説明する。第2の実施形態では、第1の実施形態の光音響装置と同様の装置を用いる。第2の実施形態においては、第1の実施形態の光音響装置と同様の構成については同じ符号を付し、詳細な説明を省略する。
前述したように、複数の波長の光を用いる場合、複数の画像データ間の位置ずれの推定精度が低い波長が含まれている可能性がある。この場合、その波長により得られた複数の画像データのみを用いると、その波長の光の照射タイミングにおける位置ずれ情報を精度よく取得することができない。すなわち、複数の波長の光を用いる場合、波長間での位置ずれの推定精度にバラつきが生じてしまう。
そこで、本発明の一実施形態に係る光音響装置は、互いに異なる複数の波長のそれぞれについて、複数の画像データ間の位置ずれ情報を取得する。そして、複数の波長のそれぞれで得られた位置ずれ情報を組み合わせることにより、先に得られた位置ずれ情報を更新する。これにより、推定精度が低くなってしまう波長が含まれている場合であっても、位置ずれの推定精度が比較的高い波長により得られた位置ずれ情報も用いて、得られた位置ずれ情報を更新するため、位置ずれ情報を精度良く取得することができる。
すなわち、本発明の一実施形態に係る光音響装置においては、互いに異なる第1の波長及び第2の波長のそれぞれの光が複数回被検体に照射される。そして、第1の波長に対応する第1の複数の画像データが生成され、第1の複数の画像データ間の位置ずれ情報が取得される。ここで取得された位置ずれ情報が、第1の波長の光の照射タイミングにおける被検体とプローブとの相対位置の変動量(位置ずれ量)に相当する。また、第2の波長に対応する第2の複数の画像データが生成され、第2の複数の画像データ間の位置ずれ情報が取得される。そして、第1の位置ずれ情報と第2の位置ずれ情報とに基づいて、第1及び第2の波長の光照射タイミングにおける位置ずれ情報が更新される。これにより、位置ずれの推定精度の低い波長の位置ずれ情報の取得の際に、位置ずれの推定精度の高い波長により得られた位置ずれ情報も用いるため、波長間での位置ずれの推定精度のバラつきが小さくなる。
さらに、それぞれの波長の位置ずれ情報を組み合わせるときの位置ずれ情報に対する重みを変更することにより、推定精度の高い波長の推定結果を優先させることもできる。
なお、信号群の取得や各波長の光照射タイミングにおける位置ずれ情報の取得等の方法については、第1の実施形態に記載の方法を適用することができる。
以下、図12に示すフローチャートに沿って、本実施形態に係る情報処理を含む光音響装置の作動を説明する。図4に示すステップと同様のステップには同一の符号を付し、詳細な説明を省略する。
(S210:波長λ2に対応する複数の画像データを生成する工程)
演算部151は、S200と同様に、S100で取得された信号群に基づいて、第2の波長の光照射に対応する複数の画像データを取得する。
(S310:波長λ2の光の照射タイミングにおける位置ずれ情報を取得する工程)
演算部151は、S300と同様に、S210で取得された波長λ2に対応する複数の画像データに基づいて、波長λ2の光の照射タイミングにおける位置ずれ情報を取得する。
(S600:波長λ1及び波長λ2の光の照射タイミングにおける位置ずれ情報に基づいて、位置ずれ情報を更新する)
演算部151は、S300で得られた波長λ1の光の照射タイミングにおける位置ずれ情報と、波長λ2の光の照射タイミングにおける位置ずれ情報とに基づいて、それぞれの位置ずれ情報を更新する。
例えば、演算部151が、S300において波長λ1の光の照射タイミングにおける位置ずれ情報として並進量Mopt_stλ1を算出し、S310において波長λ2の光の照射タイミングにおける位置ずれ情報として並進量Mopt_stλ2を算出した場合を考える。このとき、本工程において、演算部151は、式15に示すように位置ずれ情報を平均化することにより、各波長の位置ずれ情報Mopt_stλ1及びMopt_stλ2を更新してもよい。
Mopt_stλ1=(Moptλ1+Moptλ2)/2
Mopt_stλ2=(Moptλ1+Moptλ2)/2・・・式15
その他、演算部151は、各波長の位置ずれ情報を用いて、第1の実施形態で説明したように、時間的または空間的に補間することにより、各波長の位置ずれ情報を更新してもよい。
また、演算部151は、各波長の位置ずれ情報に重みづけを行うことにより、各波長の位置ずれ情報を更新してもよい。演算部151は、所定の重みを用いて各波長の位置ずれ情報を重みづけすることにより、各波長の位置ずれ情報を更新してもよい。また、演算部151は、ユーザが入力部170を用いて行った指示によって決定された重みを用いて各波長の位置ずれ情報を重みづけしてもよい。
(S700:更新された位置ずれ情報に基づいて、位置ずれ補正を行う工程)
演算部151は、S600で更新された、波長λ1及び波長λ2の光の照射タイミングにおける位置ずれ情報に基づいて、S500で説明した同様の方法で、位置ずれ補正を行う。
本実施形態では、位置ずれの推定精度の低い波長により得られた位置ずれ情報だけでなく、位置ずれの推定精度の高い波長により得られた位置ずれ情報も用いて組み合わせられた位置ずれ情報を利用するため、波長間での位置ずれの推定精度のバラつきが小さくなる。
図13は、各波長の位置ずれ情報に対する重みを決定するためのスライダーを有するGUIの模式図を示す。例えば、表示部160に表示される図13に示すようなGUIを用いて重みを決定することができる。
GUI1001には、光音響画像の表示領域1002、波長λ1の位置ずれ量のグラフ1003、及び波長λ2の位置ずれ量のグラフ1004が表示されている。グラフ1003及び1004の縦軸は位置ずれ量、横軸はパルスインデックスを示している。すなわち、グラフ1003及び1004は、波長λ1及び波長λ2でのパルスボリュームデータに対する位置ずれ量をプロットしたものである。スライダーバー1007は、複数波長の位置ずれ量の重みを決定するためのアイテムである。スライダーバー1007を左に操作すると波長λ1の位置ずれ量が優位に、右に操作すると波長λ2の位置ずれ量が優位に足しあわされる。すなわち、スライダーバー1007を左に操作すると波長λ1の重みが大きくなり、右に操作すると波長λ2の重みが大きくなり、重みづけされた最終的な位置ずれ量が決定される。そして、重みづけされた最終的な位置ずれ量を用いて得られる光音響画像が表示領域1002に表示される。さらに、ユーザがスライダーバー1007を操作すると、最終的な位置ずれ量が再計算され、更新される。そして、更新された位置ずれ量を用いて再度光音響画像が取得され、表示領域1002に表示される光音響画像が更新される。図13においては、表示領域1002に足の光音響画像1008が表示されており、光音響画像1008には血管像1009が含まれている。
このようなGUIを用いることにより、複数の波長での位置ずれ量に対する重みをスライダーバーで変化させて決定することで、表示領域に表示される光音響画像の画質(解像度等)の変化を確認しながら最終的な位置ずれ量を決定することができる。
[第3の実施形態]
以下、第3の実施形態の光音響装置の構成の構成及び処理について説明する。第3の実施形態では、第1また第2の実施形態の光音響装置と同様の装置を用いる。第3の実施形態においては、第1または第2の実施形態の光音響装置と同様の構成については同じ符号を付し、詳細な説明を省略する。
前述したように、互いに異なる複数の波長の光を用いる場合、波長間で画像強度等の画像特性が異なる場合がある。その場合、波長間で画像強度が異なるために位置ずれの推定精度が低下してしまう場合がある。
そこで、本発明の一実施形態に係る光音響装置は、互いに異なる複数の波長のそれぞれの光照射により発生する光音響波に基づいて、各波長の光による画像データを生成する。そして、波長間の画像特性の差が小さくなるような処理を行った後に、複数波長の複数の画像データを用いて位置ずれ情報を取得する。これにより、波長間の画像特性の違いにより生じる位置ずれの推定精度の低下を抑制することができる。
すなわち、本発明の一実施形態に係る光音響装置においては、互いに異なる第1の波長及び第2の波長のそれぞれの光が複数回被検体に照射される。そして、第1の波長に対応する第1の複数の画像データが生成される。また、第2の波長に対応する第2の複数の画像データが生成される。そして、第1の複数の画像データと第2の複数の画像データとの画像特性の差が小さくなるように、第1の複数の画像データ及び第2の複数の画像データの少なくとも一方に画像処理が行われる。そして、画像特性の差が小さくなる処理が行われた後の、第1の複数の画像データ及び第2の複数の画像データを用いて、第1及び第2の波長の光照射タイミングにおける位置ずれ情報が取得される。
なお、上記では、画像データに対する画像処理により、波長間の画像特性の差を小さくする例を説明したが、画像データとする前の信号に対する信号処理により結果的に波長間の画像特性の差が小さくしてもよい。
以下、図14に示すフローチャートに沿って、本実施形態に係る情報処理を含む光音響装置の作動を説明する。図4に示すステップと同様のステップには同一の符号を付し、詳細な説明を省略する。
(S800:波長λ1に対応する複数の画像データと波長λ2に対応する複数の画像データとの画像特性の差を小さくする処理を行う)
演算部151は、S200で得られた波長λ1の複数の画像データと、S210で得られた波長λ2の複数の画像データとの画像特性の差が小さくなるように、波長λ1の複数の画像データ及び波長λ2の複数の画像データに対して画像処理を行う。例えば、波長間の画像特性の差を小さくする処理としては、画像強度の最大値または最小値を揃える処理や複数波長の画像強度の画像平均や分散が同程度となるような処理を行うことが挙げられる。なお、本明細書においては、画像データの画像特性の差を小さくする処理のことを正規化処理と呼ぶ。
ここでは、画像データの画像強度を正規化する例を説明する。まず、S200及びS300で得られた全パルスボリュームデータの画像強度の最大値valmax(全パルスに対して1つの値)を求める。そして、パルスボリュームデータの画像強度のうち0以下の値を0に丸め、その後、パルスボリュームデータのボクセルの画像強度Px,y,zを、求められた画像強度の最大値valmaxが所定値valになるように正規化する。すなわち、演算部151は、式16に示すようにパルスボリュームデータを正規化する。
P′x,y,z=Px,y,z*val/valmax・・・式16
ここでP′x,y,zは、各ボクセルの正規化後の画像強度の値である。ここでは複数波長のパルスボリュームデータの最大値を所定値とするように正規化したが、各波長について最大値valmaxを求め、最大値valmaxが所定値valになるように正規化してもよい。また、0以下の値を0に丸めず、パルスボリュームデータの最小値を取得し、最小値と最大値の間の強度を0からvalに変換したり、最小値と最大値の間の強度をval’からvalに変換したりしてもよい。すなわち、画像強度が所望の数値範囲に含まれるように正規化してもよい。なお、波長間の画像強度の差が小さくなる限り、どのような方法で正規化してもよい。
また、その他の手法としては、画像の強度の平均が0、分散が1などの、特定の値になるように画像を正規化してもよいし、それぞれの値を0、1ではなく、特定の値になるように正規化してもよい。
(S900:画像特性の差を小さくする処理後の画像データに基づいて位置ずれ情報を取得する工程)
演算部151は、S800で波長間の画像特性の差を小さくする処理を行った後の、波長λ1及び波長λ2の複数の画像データに基づいて、第1及び第2の波長の光照射タイミングにおける位置ずれ情報を取得する。画像データを用いた位置ずれ情報の取得方法については、S300で説明した方法と同様の方法を採用することができる。
このように、波長間の画像特性の差が低減された複数の画像データを用いることにより、波長間の画像特性の差に起因する位置ずれの推定精度の低下を抑制することができる。そして、このようにして得られた位置ずれ情報を用いて、精度良く位置ずれ補正を行うことができる(S500)。
なお、本実施形態では、各波長において複数の画像データが生成される例を説明したが、1つの波長では1つの画像データが生成されない場合であっても、複数の画像データが生成される限り、本実施形態で説明した位置ずれ情報の取得方法を適用することができる。
[第4の実施形態]
以下、第4の実施形態の光音響装置の構成の構成及び処理について説明する。第4の実施形態では、第1、第2、または第3の実施形態の光音響装置と同様の装置を用いる。第4の実施形態においては、第1、2、または第3の実施形態の光音響装置と同様の構成については同じ符号を付し、詳細な説明を省略する。
本実施形態では、入力部170を用いたユーザの指示に基づいて、第1、第2、及び第3で説明した位置ずれ情報の取得方法の少なくとも一つの方法を実行する例を説明する。
図15は、本実施形態の表示部160に表示されるGUIを示す。GUI1501には、光音響画像の表示領域1502、波長λ1の位置ずれ量のグラフ1503、及び波長λ2の位置ずれ量のグラフ1504が表示されている。グラフ1503及び1504の縦軸は位置ずれ量、横軸はパルスインデックスを示している。すなわち、グラフ1503及び1504は、波長λ1及び波長λ2でのパルスボリュームデータに対する位置ずれ量をプロットしたものである。図15においては、表示領域1502に足の光音響画像1508が表示されており、光音響画像1508には血管像1509が含まれている。
ユーザは、位置ずれ情報の取得手法の選択ボタン1511、1512、及び1513の中から入力部170を用いて所望の取得方法を選択する。
選択ボタン1511は、第1の実施形態で説明した特定の波長に対応する位置ずれ情報を、その他に対応する位置ずれ情報の取得に用いる方法に対応するラジオボタンである。本実施形態においては、選択ボタン1511が選択されると、波長λ1に対応する位置ずれ情報を補間することにより、波長λ2に対応する位置ずれ情報を取得するモードとなる。さらに、補間方法についても、ユーザが選択ボタン1514を選択することにより、演算部151が時間的に補間するのか、空間的に補間するのかを決定することができる。
また、選択ボタン1512は、第2の実施形態で説明した複数波長のそれぞれに対応する位置ずれ情報を合成し、各位置ずれ情報を更新する方法に対応するラジオボタンである。本実施形態においては、選択ボタン1512が選択されると、波長λ1に対応する位置ずれ情報と波長λ2による位置ずれ情報とを重みづけした後に合成するモードとなる。さらに、ユーザがスライダーバー1515を操作することにより、各位置ずれ情報の重みを変更することができる。なお、スライダーバー1515の機能は、スライダーバー1007と同様である。
また、選択ボタン1513は、第3の実施形態で説明した波長間の画像特性の差が小さくなるように正規化処理された後に位置ずれ情報を取得する方法に対応するラジオボタンである。本実施形態においては、選択ボタン1513が選択されると、波長λ1の光に基づいて得られたデータ及び波長λ2の光に基づいて得られたデータに対して正規化処理を行うモードとなる。さらに、正規化処理の方法についても、ユーザが選択ボタン1516を選択することにより、演算部151が画像強度の差を小さくする正規化処理を行うのか、解像度の差を小さくする正規化処理を行うのかを決定することができる。
このように、ユーザが所望の位置ずれ情報の取得方法を選択し、選択された位置ずれ情報が適用された光音響画像を確認することにより、画像データの特性に適した位置ずれ情報の取得方法で得られた合成画像データを確認することができる。
[その他の実施形態]
また、本発明は、以下の処理を実行することによっても実現される。即ち、上述した実施形態の機能を実現するソフトウェア(プログラム)を、ネットワーク又は各種記憶媒体を介してシステム或いは装置に供給し、そのシステム或いは装置のコンピュータ(またはCPUやMPU等)がプログラムを読み出して実行する処理である。
110 光照射部
120 受信部
150 コンピュータ

Claims (15)

  1. 第1の波長の光、及び、前記第1の波長とは異なる第2の波長の光のそれぞれを複数回被検体に照射する光照射手段と、
    前記第1の波長の光及び前記第2の波長の光のそれぞれの複数回の光照射により発生する光音響波を受信することにより、信号群を出力する受信手段と、
    前記信号群に基づいて、前記第1の波長に対応する第1の複数の画像データを生成し、
    前記第1の波長に対応する第1の複数の画像データに基づいて、前記第1の波長の光の照射タイミングにおける第1の位置ずれ情報を取得し、
    前記第1の位置ずれ情報に基づいて、前記第2の波長の光の照射タイミングにおける第2の位置ずれ情報を取得する処理手段と、
    を有する
    ことを特徴とする光音響装置。
  2. 前記処理手段は、
    前記第1の位置ずれ情報、前記第2の位置ずれ情報、及び前記信号群に基づいて、合成画像データを生成する
    ことを特徴とする請求項1に記載の光音響装置。
  3. 前記処理手段は、
    前記第1の位置ずれ情報に基づいて、前記第1の複数の画像データを位置ずれ補正し、
    前記信号群に基づいて、前記第2の波長に対応する第2の複数の画像データを生成し、
    前記第2の位置ずれ情報に基づいて、前記第2の複数の画像データを位置ずれ補正し、
    位置ずれ補正された前記第1及び第2の複数の画像データに基づいて、合成画像データを生成する
    ことを特徴とする請求項1に記載の光音響装置。
  4. 前記処理手段は、
    前記第1の波長の光の照射タイミングにおける前記受信手段の第1の位置情報を取得し、
    前記第1の位置情報及び前記信号群を用いて、前記第1の複数の画像データを生成する
    ことを特徴とする請求項1から3のいずれか1項に記載の光音響装置。
  5. 前記処理手段は、
    前記第1の波長の光の照射タイミングにおける前記受信手段の第1の位置情報を取得し、
    前記第1の位置ずれ情報に基づいて、前記第1の位置情報を位置ずれ補正し、
    前記第2の波長の光の照射タイミングにおける前記受信手段の第2の位置情報を取得し、
    前記第2の位置ずれ情報に基づいて、前記第2の位置情報を位置ずれ補正し、
    位置ずれ補正された前記第1及び第2の位置情報及び前記信号群を用いて、合成画像データを取得する
    ことを特徴とする請求項1に記載の光音響装置。
  6. 前記処理手段は、
    前記第1の位置ずれ情報を時間的または空間的に補間することにより、前記第2の位置ずれ情報を取得する
    ことを特徴とする請求項1から5のいずれか1項に記載の光音響装置。
  7. 前記処理手段は、
    時間的に近傍する、前記第1及び第2の波長の光の照射タイミングにおける前記第1の位置ずれ情報を、前記第2の位置ずれ情報として取得する
    ことを特徴とする請求項1から5のいずれか1項に記載の光音響装置。
  8. 前記処理手段は、
    空間的に近傍する、前記第1及び第2の波長の光の照射タイミングにおける前記受信手段の位置に対応する前記第1の位置ずれ情報を、前記第2の位置ずれ情報として取得する
    ことを特徴とする請求項1から5のいずれか1項に記載の光音響装置。
  9. 前記処理手段は、
    酸素飽和度の空間分布を示すデータを、前記合成画像データとして生成する
    ことを特徴とする請求項1から8のいずれか1項に記載の光音響装置。
  10. 前記第1の波長は、オキシヘモグロビンのモル吸収係数とデオキシヘモグロビンのモル吸収係数が略等しい波長である
    ことを特徴とする請求項1から9のいずれか1項に記載の光音響装置。
  11. 前記処理手段は、
    前記第1の波長は、前記複数回の光照射のうち、所定の回数より多く照射された光の波長を前記第1の波長として決定する
    ことを特徴とする請求項1から9のいずれか1項に記載の光音響装置。
  12. 前記処理手段は、
    前記複数回の光照射のうち、最も多くの照射された光の波長を判定し、当該波長を前記第1の波長として決定する
    ことを特徴とする請求項1から9のいずれか1項に記載の光音響装置。
  13. 前記処理手段は、
    前記第1の複数の画像データ間の位置ずれ情報を前記第1の位置ずれ情報として取得する
    ことを特徴とする請求項1から11のいずれか1項に記載の光音響装置。
  14. 第1の波長の光、及び、前記第1の波長とは異なる第2の波長の光のそれぞれが複数回被検体に照射されることにより発生する音響波に基づいて生成された、前記第1の波長に対応する第1の複数の画像データを取得し、
    前記第1の複数の画像データに基づいて、前記第1の波長の光の照射タイミングにおける第1の位置ずれ情報を取得し、
    前記第1の位置ずれ情報に基づいて、前記第2の波長の光の照射タイミングにおける第2の位置ずれ情報を取得する
    ことを特徴とする情報処理方法。
  15. 請求項14に記載の情報処理方法をコンピュータに実行させるためのプログラム。
JP2016188408A 2016-09-27 2016-09-27 光音響装置、情報処理方法、及びプログラム Active JP6759032B2 (ja)

Priority Applications (5)

Application Number Priority Date Filing Date Title
JP2016188408A JP6759032B2 (ja) 2016-09-27 2016-09-27 光音響装置、情報処理方法、及びプログラム
CN201710881220.3A CN107865639B (zh) 2016-09-27 2017-09-26 光声装置、信息处理方法和程序
KR1020170123862A KR20180034273A (ko) 2016-09-27 2017-09-26 광음향 장치, 정보처리방법, 및 프로그램
EP17193277.5A EP3298958A1 (en) 2016-09-27 2017-09-26 Photoacoustic apparatus, information processing method, and program
US15/716,150 US11529057B2 (en) 2016-09-27 2017-09-26 Photoacoustic apparatus, information processing method, and program

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2016188408A JP6759032B2 (ja) 2016-09-27 2016-09-27 光音響装置、情報処理方法、及びプログラム

Publications (2)

Publication Number Publication Date
JP2018050774A JP2018050774A (ja) 2018-04-05
JP6759032B2 true JP6759032B2 (ja) 2020-09-23

Family

ID=59997113

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2016188408A Active JP6759032B2 (ja) 2016-09-27 2016-09-27 光音響装置、情報処理方法、及びプログラム

Country Status (5)

Country Link
US (1) US11529057B2 (ja)
EP (1) EP3298958A1 (ja)
JP (1) JP6759032B2 (ja)
KR (1) KR20180034273A (ja)
CN (1) CN107865639B (ja)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2018050775A (ja) * 2016-09-27 2018-04-05 キヤノン株式会社 光音響装置、情報処理方法、及びプログラム
JP6759032B2 (ja) 2016-09-27 2020-09-23 キヤノン株式会社 光音響装置、情報処理方法、及びプログラム
CN110891481B (zh) * 2017-08-24 2023-07-04 松下知识产权经营株式会社 生物体计测装置及头戴式显示器装置
WO2019092950A1 (ja) * 2017-11-13 2019-05-16 ソニー株式会社 画像処理装置、画像処理方法および画像処理システム
JP7195759B2 (ja) * 2018-04-20 2022-12-26 キヤノン株式会社 光音響装置および被検体情報取得方法

Family Cites Families (23)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5818444B2 (ja) 2010-02-04 2015-11-18 キヤノン株式会社 機能情報取得装置、機能情報取得方法、及びプログラム
JP6173402B2 (ja) 2010-02-04 2017-08-02 キヤノン株式会社 機能情報取得装置、機能情報取得方法、及びプログラム
JP5627328B2 (ja) 2010-07-28 2014-11-19 キヤノン株式会社 光音響診断装置
JP5281624B2 (ja) 2010-09-29 2013-09-04 日本電信電話株式会社 画像符号化方法,画像復号方法,画像符号化装置,画像復号装置およびそれらのプログラム
JP2012196308A (ja) * 2011-03-22 2012-10-18 Fujifilm Corp 光音響画像生成装置及び方法
JP5704998B2 (ja) * 2011-04-06 2015-04-22 キヤノン株式会社 光音響装置およびその制御方法
JP5730253B2 (ja) 2011-09-27 2015-06-03 富士フイルム株式会社 レーザ光源ユニット及び光音響画像生成装置
US9814394B2 (en) * 2011-11-02 2017-11-14 Seno Medical Instruments, Inc. Noise suppression in an optoacoustic system
JP5840069B2 (ja) 2012-05-08 2016-01-06 富士フイルム株式会社 光音響画像生成装置、システム、及び方法
KR102170350B1 (ko) 2012-06-13 2020-10-27 세노 메디컬 인스투르먼츠 인코포레이티드 광음향 데이터의 파라메트릭 맵들을 생성하기 위한 시스템 및 방법
US20150351639A1 (en) * 2012-12-28 2015-12-10 Canon Kabushiki Kaisha Subject information obtaining apparatus, method for controlling subject information obtaining apparatus, and program
EP2749209A1 (en) 2012-12-28 2014-07-02 Canon Kabushiki Kaisha Object information acquisition apparatus, display method, and program
JP6192297B2 (ja) 2013-01-16 2017-09-06 キヤノン株式会社 被検体情報取得装置、表示制御方法、およびプログラム
JP6238550B2 (ja) * 2013-04-17 2017-11-29 キヤノン株式会社 被検体情報取得装置、被検体情報取得装置の制御方法
US10143381B2 (en) * 2013-04-19 2018-12-04 Canon Kabushiki Kaisha Object information acquiring apparatus and control method therefor
KR101581689B1 (ko) 2014-04-22 2015-12-31 서강대학교산학협력단 움직임 보정을 이용한 광음향 영상 획득 장치 및 방법
JP6498036B2 (ja) * 2014-06-13 2019-04-10 キヤノン株式会社 光音響装置、信号処理方法、及びプログラム
US9880053B2 (en) * 2014-10-29 2018-01-30 Panasonic Intellectual Property Management Co., Ltd. Image pickup apparatus, spectroscopic system, and spectroscopic method
US20170332913A1 (en) * 2014-11-07 2017-11-23 Canon Kabushiki Kaisha Object information acquiring apparatus
JP2016101393A (ja) * 2014-11-28 2016-06-02 キヤノン株式会社 被検体情報取得装置およびその制御方法
JP2018050776A (ja) 2016-09-27 2018-04-05 キヤノン株式会社 光音響装置、情報処理方法、及びプログラム
JP6759032B2 (ja) 2016-09-27 2020-09-23 キヤノン株式会社 光音響装置、情報処理方法、及びプログラム
JP2018050775A (ja) 2016-09-27 2018-04-05 キヤノン株式会社 光音響装置、情報処理方法、及びプログラム

Also Published As

Publication number Publication date
JP2018050774A (ja) 2018-04-05
US20180085005A1 (en) 2018-03-29
CN107865639A (zh) 2018-04-03
KR20180034273A (ko) 2018-04-04
EP3298958A1 (en) 2018-03-28
CN107865639B (zh) 2021-01-19
US11529057B2 (en) 2022-12-20

Similar Documents

Publication Publication Date Title
JP6759032B2 (ja) 光音響装置、情報処理方法、及びプログラム
JP5528083B2 (ja) 画像生成装置、画像生成方法、及び、プログラム
JP6440140B2 (ja) 被検体情報取得装置、処理装置、および信号処理方法
JP6366272B2 (ja) 被検体情報取得装置、被検体情報取得装置の制御方法、およびプログラム
JP6498036B2 (ja) 光音響装置、信号処理方法、及びプログラム
JP6525565B2 (ja) 被検体情報取得装置および被検体情報取得方法
US10548477B2 (en) Photoacoustic apparatus, information processing method, and storage medium
JP2017119094A (ja) 情報取得装置、情報取得方法、及びプログラム
US10578588B2 (en) Photoacoustic apparatus, information processing method, and storage medium
JP4879872B2 (ja) 画像処理装置、画像処理プログラム、記憶媒体及び超音波診断装置
US10436706B2 (en) Information processing apparatus, information processing method, and storage medium
US10445897B2 (en) Device for acquiring information relating to position displacement of multiple image data sets, method, and program
JP2018089371A (ja) 情報処理装置、情報処理方法、及びプログラム
JP2019092929A (ja) 情報処理装置、情報処理方法、プログラム
US20200113541A1 (en) Information processing apparatus, information processing method, and storage medium
US20200305727A1 (en) Image processing device, image processing method, and program
JP2018143764A (ja) 画像生成装置、画像生成方法、プログラム
JP6682207B2 (ja) 光音響装置、画像処理方法、プログラム
JP2017042603A (ja) 被検体情報取得装置
JP2019122621A (ja) 被検体情報取得装置および被検体情報取得方法
JP2017018284A (ja) 光音響装置

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20190924

TRDD Decision of grant or rejection written
A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20200731

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

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20200804

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20200902

R151 Written notification of patent or utility model registration

Ref document number: 6759032

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R151