JP2013128759A - 光音響信号処理装置及び方法 - Google Patents

Abstract

【課題】光音響信号処理装置において、光パルス波形に変化が生じるような場合でも、検出信号から正しく吸収分布を求めるようにする。
【解決手段】被検体に対して光が照射されることで、被検体内で光音響信号が生じる。プローブ１１で検出された光音響信号は、受信回路２１を介して、ＡＤ変換手段２２でサンプリングされる。光波形微分信号生成手段２５は、パルス光の光強度の時間波形又はその時間波形を微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、被検体に照射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成する。光波形微分信号逆畳込み手段２６は、光音響信号から光波形微分信号をデコンボリューションする。このデコンボリューションにより、吸収分布が得られる。
【選択図】図４

Description

本発明は、光音響信号処理装置及び方法に関し、更に詳しくは、被検体に照射されたパルス光により被検体内で生じた光音響信号を信号処理する光音響信号処理装置及び方法に関する。

生体内部の状態を非侵襲で検査できる画像検査法の一種として、超音波検査法が知られている。超音波検査では、超音波の送信及び受信が可能な超音波探触子を用いる。超音波探触子から被検体（生体）に超音波を送信させると、その超音波は生体内部を進んでいき、組織界面で反射する。超音波探触子でその反射音波を受信し、反射超音波が超音波探触子に戻ってくるまでの時間に基づいて距離を計算することで、内部の様子を画像化することができる。

また、光音響効果を利用して生体の内部を画像化する光音響イメージングが知られている。一般に光音響イメージングでは、パルスレーザ光を生体内に照射する。生体内部では、例えば生体組織がパルスレーザ光のエネルギーを吸収し、そのエネルギーによる断熱膨張により超音波（光音響信号）が発生する。この光音響信号を超音波プローブなどで検出し、検出信号に基づいて光音響画像を構成することで、光音響信号に基づく生体内の可視化が可能である。

ここで、通常の再構成方法（例えばフーリエドメイン法（ＦＴＡ法）や遅延加算法（Delay&Sum法））は、実質的に圧力分布を画像化しており、生成された光音響画像は光吸収体の分布画像にはなっていない。圧力分布を画像化した光音響画像では、１本の血管が二重に表示されることがある。この現象は、以下のように説明できる。すなわち、微小要素（微小吸収体）が発生する圧力波形をミクロ波形とすると、マクロ構造物としての血管から発生する圧力波形は、ミクロ波形の重ね合わせとして捉えることができる。微小吸収体で発生したミクロ波形を重ね合わせていくと、血管内部で発生したミクロの圧力波形は近接した前後のミクロ波形と重ね合わされてキャンセルされる。血管から出たミクロの圧力波形を全て重ね合わせると、最終的に血管表面で発生した圧力波形の一部（前端）と裏面で発生した圧力波形の一部（後端）とがキャンセルされずに残り、これらが観測されることになる。そのため、血管が二重に表示される。１本の血管が二重に表示されると、画像判定上、血管の位置の確認が困難になるという問題が生じる。また、複数の異なる波長で取得した画像同士を計算する際に、位置がずれやすく適切な結果が得られないなどの問題も生じる。

圧力分布画像に代えて、吸収分布画像を生成する技術がこれまでにいくつか知られている。例えば特許文献１には、試料の熱的インパルス応答から光音響画像の分解能劣化を修復する逆フィルタを求め、得られた光音響画像に逆フィルタを作用させることで、理想的な光音響画像、すなわち試料表面の点光源によって励起され、検出されるその点（無限小）における熱的インピーダンス情報（＝無限小なる点熱源の集合）を得ることが記載されている。

特許文献１では、より詳細には、まず、試料の熱インパルス応答ｈ（ｘ，ｙ）を計算し、次いで光音響画像ｐ（ｘ，ｙ）を構成する。熱的インパルス応答は、無限小なる一点の温度変化が試料表面の微小変位に変換されるまでの伝達関数と定義されている。その後、熱インパルス応答ｈ（ｘ，ｙ）と光音響画像ｐ（ｘ，ｙ）とをそれぞれフーリエ変換し、フーリエ変換像Ｈ（μ，υ）、Ｐ（μ，υ）を得る。１／Ｈ（μ，υ）を逆フィルタとして用い、Ｑ（μ，υ）＝Ｐ（μ，υ）・（１／Ｈ（μ，υ））により、Ｑ（μ，υ）を計算する。そのようにして計算されたＱ（μ，υ）をフーリエ逆変換することで、理想的な光音響画像ｑ（ｘ，ｙ）が得られる。

上記特許文献の他にも、非特許文献１には、論理的には、有限の時間幅を持つ光パルスη（ｔ）をフーリエ変換したものη（ｋ）としたとき、その微分をｉη（ｋ）として考慮することが記載されている。実験的には、超音波検出素子であるＰＺＴ（チタン酸ジルコン酸鉛）での検出帯域以内に励起光パルス波形が入るように光パルス幅を長くしたマイクロ波を被検体に照射し、通常のＰＺＴプローブで光音響信号を検出し、吸収分布を再構成している。

また、非特許文献２には、被写体中の微小要素からの圧力波形として被写体中の微小要素からの圧力波形として光パルス微分関数と装置インパルス応答関数とを合わせたミクロ波形と、吸収分布とを、観測圧力波形に関連付けることが記載されている。吸収像再構成は、光微分とシステム応答とを不可分な状態で含むｐ_ｄ０を測定して、各素子の圧力波形からｐ_ｄ０をデコンボリューションしたのちに、フィルタ補正逆投影法（Filtered Backprojection法）を用いて行う。実験的には、光パルス幅の短いパルスレーザ光で励起し、超音波の検出帯域を通常の超音波診断装置よりも広げて、ハイドロホンとオシロスコープとで光音響信号を検出し、吸収分布を再構成する。

特開平３−１５６３６２号公報

Yuan Xu, et al., IEEE Transactions on Medical Imaging, Volume 21 (2002),p.823-828. Yi Wang, et al., Physics in Medicine and Biology, Volume 49 (2004),p.3117-3124

非特許文献１では、光パルスをフーリエ変換したη（ｋ）を位置依存の関数η（ｒ，ｋとして扱っていない。このため、ｔ＝０に再構成できる（検出した波形を基に時間を遡って光が入射した瞬間（ｔ＝０）に発生した圧力分布を推定する計算を行うことができる）場合には正確な吸収分布を得ることができるものの、ｔ≠０に再構成した場合（光が入射した瞬間（ｔ＝０）の圧力分布を再構成できずに、時刻ｔ＝０からしばらく時間が経った後の圧力分布を推定することになる場合）には、光パルス幅成分が除去できず、圧力分布となる。

また、非特許文献１及び２では、励起レーザと超音波検出装置とのどちらかを実用的な範囲から外してレーザ発光時間と超音波検出時とを合わせて再構成している。このため、非特許文献１及び２では問題は明確には現れないものの、実用的な装置構成を考えた場合、非特許文献１及び２の手法では、ｔ＝０とする再構成が困難である。すなわち、例えば実用的な装置構成として、
・サンプリング周波数１００ＭＨｚ以下で、ＰＺＴなどを利用した狭帯域プローブを用いた超音波検出装置
・強い光音響信号が出る１−１００ｎｓオーダーの光パルス幅を持つ励起レーザ
を用いた場合、レーザ発光が超音波検出時間と比較して短時間の現象のため、ｔ＝０に相当する状態（吸収分布と圧力分布とが比例する時間帯）に正確に再構成できない。

ここで、「時刻ｔ＝０の圧力分布」は「吸収分布」を表すので、時刻ｔ＝０圧力分布が求められれば吸収分布を得ることができる。しかし、一般的な超音波検出装置のサンプリング間隔は２５ｎｓ程度であり、光が当たった瞬間の時刻ｔ＝０のつもりで計算しても、実際はｔ＝±１２．５ｎｓ程度の時間幅でずれが生じる。例えば光パルス幅が１００ｎｓと長い場合には、上記のずれ（±１２．５ｎｓ）は誤差と考えればよいものの、光パルス幅が１０ｎｓであれば、上記のずれは誤差とは呼べなくなり、光が当たった瞬間の圧力分布というよりも圧力波の伝播過程の圧力分布に移行することとなる。その「圧力波の伝播過程の圧力分布」は、「吸収分布」に一致しない。

更に、実験的な生体等のサンプルにおいては、ｔ＝０の圧力分布を定義することが困難である。生体内の音速を例えば１５３０ｍ／ｓと仮定し、検出時刻とレーザ照射時刻との差を伝播時間とすると、伝播時間から伝播距離が求める。生体内の音速が１５３０ｍ／ｓで一定であれば、伝播時間から求めた伝播距離は実際の伝播距離と一致する。しかし、実際には、音速は生体内で一様ではなく、計算上の伝播距離と実際の伝播距離とにずれが生じる。従って、検出信号から伝播距離を推定する場合には、音速差に起因する伝搬距離のあいまいさが残る。生体内の伝播距離のあいまいさを伝播時間のあいまいさと捉えると、時刻ｔ＝０も曖昧性を持つことになり、ｔ＝０の圧力分布というのもあいまいになり、定義が困難となる。ｔ＝０の分布が吸収分布であるのに対し、ｔ＞０の分布は伝播時の圧力分布であり、これらが混ざると吸収分布とは言えなくなる。

上記考察から、本発明者は、観測波形又は再構成後の圧力波形から、光パルス波形の微分波形をデコンボリューションすることが有利であることを見出した。観測波形から光パルス微分波形をデコンボリューションすることで、吸収分布を得ることができる。しかしながら、被検体に照射されるパルスレーザ光の時間波形は、常に一定であるとは限らない。例えば複数波長でパルスレーザ光を照射する場合に、ある波長のパルスレーザ光の時間波形と、別の波長のパルスレーザ光の時間波形とに違いが生じることがある。その違いを考慮せずに、双方の波長に対して同一波形の光微分波形をデコンボリューションすると、観測波形に重畳されている実際の光微分波形とは異なる波形でデコンボリューションが行われ、正しく吸収分布を求めることができなくなる。

本発明は、上記に鑑み、光パルス波形に変化が生じるような場合でも、検出信号から正しく吸収分布を求めることができる光音響信号処理装置及び方法を提供することを目的とする。

上記課題を解決するために、本発明は、被被検体に対するパルス光の出射後に、パルス光より被検体内で生じた光音響信号をサンプリングするサンプリング手段と、パルス光の光強度の時間波形又はそれを微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、被検体に対して出射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成する光波形微分信号生成手段と、サンプリングされた光音響信号から光波形微分信号生成手段で生成された光波形微分信号をデコンボリューションした信号を生成する光波形微分信号逆畳込み手段とを備えたことを特徴とする光音響信号処理装置を提供する。

本発明の光音響信号処理装置では、光波形微分信号生成手段が、パルス光に関連した発光条件と、その発光条件下で被検体に出射されるパルス光のパルス時間幅との対応関係を示す光パルス幅情報を参照し、その光パルス幅情報と被検体に出射されたパルス光の発光条件とに基づいて、被検体に出射されたパルス光のパルス時間幅を推定する光パルス幅推定部と、推定されたパルス時間幅を関数に適用する関数適用部とを含む構成を採用することができる。

光パルス幅情報は、発光条件とパルス時間幅とを対応付けたテーブルを含んでいてもよい。

光パルス幅推定部は、被検体に出射されたパルス光の発光条件に合致する発光条件がテーブルに存在しないときは、テーブルに存在するパルス時間幅に基づいて補間により被検体に出射されたパルス光のパルス時間幅を推定することとしてもよい。

パルス光は、パルスレーザ光であってもよい。

発光条件は、パルスレーザ光の波長、レーザ発光のトリガからレーザ発光までの遅延時間、パルスレーザ光の発光周波数、レーザ光源におけるレーザ結晶の励起エネルギー、レーザ光源における経時変化に関する条件、レーザ光源におけるレーザ結晶温度のうちの少なくとも１つを含むことができる。

関数は、事前に実測されたパルス光の時間波形のフィッティング曲線、ガウス関数、ローレンツ関数、若しくは擬フォークト関数、又は、これら関数を組み合わせたものとすることができる。

本発明の光音響信号処理装置が、被検体に対して出射されたパルス光のパルス時間幅を測定する光パルス幅測定手段を更に備え、光波形微分信号生成手段が、光パルス幅測定手段で測定されたパルス時間幅を関数に適用することとしてもよい。

上記の場合、被検体に出射されるパルス光の光源の光出射端から被検体に対してパルス光が出射されるまでの間の光路中に、被検体に出射されるパルス光の一部を分岐する分岐部を設けておき、光パルス幅測定手段が、光路から分岐されたパルス光の光強度に基づいて被検体に対して出射されたパルス光のパルス時間幅を測定するようにしてもよい。

本発明の光音響信号処理装置では、光波形微分信号逆畳込み手段が、サンプリングされた光音響信号をフーリエ変換する第１のフーリエ変換手段と、光波形微分信号を所定のサンプリングレートでサンプリングした信号をフーリエ変換する第２のフーリエ変換手段と、フーリエ変換された光波形微分信号の逆数を逆フィルタとして求める逆フィルタ演算手段と、フーリエ変換された光音響信号に逆フィルタを適用するフィルタ適用手段と、逆フィルタが適用された光音響信号をフーリエ逆変換するフーリエ逆変換手段とを有する構成を採用することができる。

本発明では、光音響信号が第１のサンプリングレートでサンプリングされ、光波形微分信号が第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、第１のサンプリングレートでサンプリングされた光音響信号を、第２のサンプリングレートでリサンプルするリサンプル手段を更に備え、第１のフーリエ変換手段が、リサンプル手段でリサンプルされた光音響信号をフーリエ変換する構成を採用できる。

上記に代えて、光音響信号が第１のサンプリングレートでサンプリングされ、光波形微分信号が第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、かつ、第１のフーリエ変換手段が第１のデータ点数でフーリエ変換を行い、第２のフーリエ変換手段が、第１のデータ点数よりも多い第２のデータ点数でフーリエ変換を行うものであり、フーリエ変換された光音響信号に対して、第１のデータ点数と第２のデータ点数との差の分だけ中央に０を付加するゼロパディングを行うゼロパディング手段を更に備え、フィルタ適用手段が、ゼロパディング手段でゼロパディングが施された信号に対して逆フィルタを適用する構成としてもよい。

また、光音響信号が第１のサンプリングレートでサンプリングされ、光波形微分信号が第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、かつ、第１のフーリエ変換手段が第１のデータ点数でフーリエ変換を行い、第２のフーリエ変換手段が、第１のデータ点数よりも多い第２のデータ点数でフーリエ変換を行うものであり、フーリエ変換された光波形微分信号から、第１のデータ点数と第２のデータ点数の差の分だけ高周波成分サンプル点を除去する高周波成分サンプル点除去手段を更に備え、逆フィルタ演算手段が、フーリエ変換された光波形微分信号から高周波成分サンプル点を除去した信号の逆数を逆フィルタとして求めるこことしてもよい。

光音響信号のサンプリング間隔は、被検体に出射された光のパルス時間幅よりも長くしてもよい。

本発明の光音響信号処理装置は、複数の検出素子で検出され、サンプリング手段でサンプリングされた光音響信号に基づいて光音響信号を再構成する光音響信号再構成手段を更に備え、光波形微分信号逆畳込み手段が、光音響信号再構成手段で再構成された光音響信号から光波形微分信号をデコンボリューションする構成を採用することができる。

光波形微分信号逆畳込み手段は、光音響信号から、光波形微分信号に装置応答関数を掛け合わせたものをデコンボリューションしてもよい。

本発明では、被検体に対する出射光が複数の波長の光を含み、サンプリング手段が各波長の光に対応した光音響信号をサンプリングし、光波形微分信号逆畳込み手段が各波長の光に対応した光音響信号から光波形微分信号をデコンボリューションした信号を生成するものであってもよく、その場合、各波長の光に対応した、光音響信号から光波形微分信号をデコンボリューションした信号を演算処理する２波長データ演算手段を更に備えた構成とすることができる。

光音響信号から光波形微分信号をデコンボリューションした信号に基づいて光音響画像を生成する光音響画像生成手段を更に備える構成としてもよい。

サンプリング手段が、更に被検体に対して送信された音響波に対する反射音響波をサンプリングするものであり、サンプリングされた反射音響波に基づいて反射音響波画像を生成する反射音響波画像生成手段と、光音響画像と反射音響波画像とを合成する画像合成手段とを更に備える構成を採用してもよい。

画像合成手段は、光音響画像と反射音響波画像とを重畳することで画像合成を行ってもよい。

本発明は、また、被検体に対するパルス光の出射後に、そのパルス光に起因して被検体内で発生した光音響信号を検出するステップと、パルス光の光強度の時間波形又はそれを微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、被検体に出射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成するステップと、検出された光音響信号から、生成された光波形微分信号をデコンボリューションするステップとを有する光音響信号処理方法を提供する。

本発明の光音響信号処理装置及び方法では、パルス光の光強度の時間波形、又はその時間波形を微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、被検体に照射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成する。パルス時間幅を独立変数とする関数を用いて光波形微分信号を生成することで、被検体に照射されるパルス光の光強度の時間波形が常に一定ではなく、パルス時間幅に変動が生じる場合でも、パルス時間幅に応じた光波形微分信号を生成することができる。そのように生成された光波形微分信号をデコンボリューションすることで、光音響信号の検出信号から吸収分布を正しく求めることができる。

光パルス微分波形デコンボリューションの基本アルゴリズムを示すブロック図。 再構成後の光音響信号を示す波形図。 ＦＦＴ後の光音響信号ＦＦＴを示す波形図。 光パルス微分波形（ｈ）を示す波形図。 ＦＦＴ後の光パルス微分波形ＦＦＴ（ｆｆｔ＿ｈ）を示す波形図。 光パルス微分波形ＦＦＴフィルタを示す波形図。 デコンボリューション後のＦＦＴ波形を示す波形図。 逆変換された光音響信号を示す波形図。 再構成後の光音響信号に基づいて生成した光音響画像を示す図。 デコンボリューション後の光音響信号に基づいて生成した光音響画像を示す図。 本発明の第１実施形態の光音響信号処理装置を含む光音響画像生成装置を示すブロック図。 光波形微分信号生成手段を示すブロック図。 光波形微分信号逆畳込み手段を示すブロック図。 光音響画像生成の動作手順を示すフローチャート。 本発明の第２実施形態における光波形微分信号逆畳込み手段を示すブロック図。 サンプリングレート４００ＭＨｚ相当の光波形微分信号を示す波形図。 サンプリングレート４０ＭＨｚ相当の光波形微分信号を示す波形図。 本発明の第３実施形態における光波形微分信号逆畳込み手段を示すブロック図。 光音響信号（周波数領域）を示すグラフ。 ゼロパディング後の光音響信号を示すグラフ。 本発明の第４実施形態における光波形微分信号逆畳込み手段を示すブロック図。 光波形微分信号（周波数領域）を示すグラフ。 高周波成分サンプル点が除去された光波形微分信号を示すグラフ。 本発明の第５実施形態の光音響画像生成装置を示すブロック図。 本発明の第６実施形態の光音響画像生成装置を示すブロック図。 本発明の変形例の光音響画像生成装置を示すブロック図。

本発明の実施の形態の説明に先立って、本発明の概要について説明する。光吸収体であるミクロ吸収粒子を考え、このミクロ吸収粒子がパルスレーザ光を吸収して圧力波（光音響圧力波）が生じることを考える。時刻をｔとして、位置ｒにあるあるミクロ吸収粒子から発生する光音響圧力波を、位置Ｒで観測した場合の圧力波形ｐ_{ｍｉｃｒｏ}（Ｒ，ｔ）は、[Phys. Rev. Lett. 86(2001)3550.]より、以下の球面波となる。
ここで、Ｉ（ｔ）は励起光の光強度の時間波形であり、係数ｋは、粒子が光を吸収して音響波を出力する際の変換係数であり、ｖ_ｓは被検体の音速である。また、位置ｒ、Ｒは、空間上の位置を示すベクトルである。ミクロ吸収粒子から発生する圧力は、上記式に示すように、光パルス微分波形に比例した球面波となる。

実際にイメージングする対象から得られる圧力波形は、よりマクロな吸収体のサイズを有しているため、上記のミクロ吸収波形を重ね合わせた波形になると考える（重ね合わせの原理）。ここで、マクロな光音響波を発する粒子の吸収分布をＡ（ｒ−Ｒ）とし、そのマクロな吸収体からの圧力の観測波形をｐ_{ｍａｃｒｏ}（Ｒ，ｔ）とする。観測位置Ｒでは、各時刻において、観測位置Ｒから半径ｖ_ｓｔに位置する吸収粒子からの光音響波が観測されることになるため、観測波形ｐ_{ｍａｃｒｏ}（Ｒ，ｔ）は、以下の圧力波形の式で示される。
上記式（１）からわかるように、観測波形は、光パルス微分のコンボリューション型を示す。そこで、本発明では、デコンボリューション処理において、光パルス微分項を考慮してデコンボリューションを行うこととする。

更に、本発明では、従来の超音波システムにおいても用いられている圧力分布を求める再構成（ＦＴＡ法、ＤｎＳ法、ＢＰ法など）を適用後に、再構成後の画像がｔ≠０の圧力分布、すなわち圧力波の伝播過程の圧力分布であることを認識の上で、これを吸収分布に変換することを考えた。圧力分布再構成の基本的な考え方としては、検出位置Ｒ＝（ｘ，ｙ，０）の再構成後の圧力分布ｐ_ｒｅｃ（Ｒ，ｔ）は、各時刻ｔにおける、Ｒの検出軸（ｒ−Ｒ）上の｜ｒ−Ｒ｜位置に存在する吸収体から発生する球面波を、周囲の圧電素子の信号も含めて足し合わせてその位置における圧力強度を計算して得られる。従って、検出軸（ｒ−Ｒ）に存在するミクロの吸収体から発生し、伝播する光音響波を重ねあわせたｐ_ｒｅｃ（Ｒ，ｔ）は以下のように表記できる。
このように吸収分布を１次元で考えて良くなることで、上記式のような圧力表記が可能となる。上記式（２）は、検出軸（ｒ−Ｒ）をｚ軸、検出素子からの距離｜ｒ−Ｒ｜をｚとすると、下記のように表すことができる。
更に、式（３）から積分とは関係ないｘ，ｙは表記を省略し、ｚ軸を時間で表記すると、上記式は、下記式のように表すことができる。
このように、（ｘ，ｙ，０）に位置する検出素子における１軸（時間軸又はｚ’軸）のコンボリューション表記が可能となる。

上記式（４）の両辺をフーリエ変換し、周波数軸において、圧力分布のフーリエ係数を光パルスの時間微分のフーリエ係数で割ることで、光パルス微分をデコンボリューションすることができる。
デコンボリューション後、得られた式を、フーリエ逆変換することで、Ａ（ｘ，ｙ，ｖ_ｓｔ）を求め、吸収分布を画像化することができる。ここで求めたＡ（ｘ，ｙ，ｖ_ｓｔ）には、検出素子受信角度依存性Ｄ（ｘ，ｙ，ｚ）や、プローブ帯域の固有振動が重畳されている可能性がある。例えば、装置関数Ｄ（ｘ，ｙ，ｚ）を事前に求めておいた上でその逆数をＡ（ｘ，ｙ，ｖ_ｓｔ）に掛けることで、検出素子受信角度依存性の影響を除去できる。また、帯域の固有振動に関してはヒルベルト変換や直交検波処理により強度画像化すれば、その影響を除去することができる。更に、検体に入射する光空間分布Ｌ（ｘ，ｙ，ｚ）を別に観測やシミュレーションにより求め、吸収係数に比例する画素値μ（ｘ，ｙ，ｚ）を、μ（ｘ，ｙ，ｚ）＝Ａ（ｘ，ｙ，ｖ_ｓｔ）／Ｌ（ｘ，ｙ，ｚ）により求めてもよい。この場合、より生体組織と密接な関係のある物理量である吸収係数の分布画像を得ることができる。

図１は、光パルス微分波形デコンボリューションの基本アルゴリズムを示す。再構成後の光音響信号を入力し、再構成後の光音響信号をＦＦＴ（Fast Fourier Transform）によりフーリエ変換する（ステップＳ１０１）。図２Ａに再構成後の光音響信号を示し、図２ＢにＦＦＴ後の光音響信号ＦＦＴを示す。フーリエ変換することで、図２Ａに示す時間領域の信号が、図２Ｂに示すような周波数領域の信号に変換される。なお、図２Ｂでは、光音響信号ＦＦＴの絶対値を示しているが、実際の処理では複素数のまま処理される。

光パルス微分波形ｈをＦＦＴによりフーリエ変換する（ステップＳ１０２）。図２Ｃに光パルス微分波形（ｈ）を示し、図２ＤにＦＦＴ後の光パルス微分波形ＦＦＴ（ｆｆｔ＿ｈ）を示す。フーリエ変換することで、図２Ｃに示す時間領域の信号（波形）が、図２Ｄに示す周波数領域の信号に変換される。なお、図２Ｃにおける黒丸は、光パルス微分波形におけるサンプリング点を表している。また、図２Ｄでは、光パルス微分波形ＦＦＴの絶対値を示しているが、実際の処理では複素数のまま処理される。

ステップＳ１０２で得られたＦＦＴ後の光パルス微分波形ＦＦＴ（ｆｆｔ＿ｈ）の逆数を、光パルス微分波形ＦＦＴフィルタ（逆フィルタ）として求める（ステップＳ１０３）。光パルス微分波形ＦＦＴフィルタは、具体的にはconj(fft_h)/abs(fft_h)²で求めることができる。ここで、conj(fft_h)はfft_hの共役複素数、abs(fft_h)はfft_hの絶対値を表す。図２Ｅに、光パルス微分波形ＦＦＴフィルタを示す。図２Ｄに示す光パルス微分波形ＦＦＴの逆数を求めることで、図２Ｅに示すような光パルス微分波形ＦＦＴフィルタを得ることができる。

上記のようにして求めた光パルス微分ＦＦＴフィルタと、再構成後の光音響信号ＦＦＴとを要素ごとに乗算し、光音響信号ＦＦＴから光パルス微分波形をデコンボリューションする（ステップＳ１０４）。図２Ｆに、デコンボリューション後のＦＦＴ波形を示す。図２Ｂに示す光音響信号ＦＦＴと図２Ｅに示す光パルス微分波形ＦＦＴフィルタとの乗算を行うことで、図２Ｆに示すＦＦＴ波形が得られる。

ステップＳ１０４で光パルス微分波形をデコンボリューションしたＦＦＴ波形を、逆ＦＦＴによりフーリエ逆変換し、周波数領域の信号を時間領域の信号に戻す（ステップＳ１０５）。図２Ｇは、逆変換された光音響信号を示す。図２Ｆに示すＦＦＴ波形（周波数領域の信号）を逆ＦＦＴすることで、図２Ｇに示すデコンボリューション後の光音響信号（時間領域の信号）が得られる。このデコンボリューション後の光音響信号は、光吸収分布に光パルス微分波形（図２Ｃ）がコンボリューションされた再構成後の光音響信号（図２Ａ）から、光パルス微分波形をデコンボリューションした吸収分布に相当する。

図３Ａに、再構成後の光音響信号（図２Ａ）に基づいて生成した光音響画像を示し、図３Ｂに、デコンボリューション後の光音響信号（図２Ｇ）に基づいて生成した光音響画像を示す。図３Ａに示す、再構成後の光音響信号に基づいて生成した光音響画像は、実質的に圧力分布を画像化したものであり、１本の血管が二重に表示されるなど、画像判定上、血管の位置が確認しづらい。これに対し、図３Ｂに示すデコンボリューション後の光音響信号に基づいて生成した光音響画像は、光パルス微分波形をデコンボリューションしていることで吸収体の分布を画像化できており、血管の位置を確認しやすくなっている。

以下、図面を参照し、本発明の実施の形態を詳細に説明する。図４は、本発明の第１実施形態の光音響信号処理装置を含む光音響画像生成装置を示す。光音響画像生成装置（光音響画像診断装置）１０は、超音波探触子（プローブ）１１、超音波ユニット１２、及び光源ユニット（レーザユニット）１３を備える。

レーザユニット１３は、被検体に照射するレーザ光を生成する。レーザ光の波長は、観察対象物に応じて適宜設定すればよい。レーザユニット１３が出射するレーザ光は、例えば光ファイバなどの導光手段を用いてプローブ１１まで導光され、プローブ１１から被検体に照射される。プローブ１１は、レーザユニット１３から出射した光が被検体に照射された後に、被検体内の光吸収体がレーザ光を吸収することで生じた超音波（光音響信号）を検出する。プローブ１１は、例えば一次元配列された複数の超音波振動子を有する。

超音波ユニット１２は、光音響信号処理装置に相当する。超音波ユニット１２は、受信回路２１、ＡＤ変換手段２２、光音響画像再構成手段２４、光波形微分信号生成手段２５、光波形微分信号逆畳込み手段２６、検波・対数変換手段２７、光音響画像構築手段２８、トリガ制御回路２９、及び制御手段３０を有する。受信回路２１は、プローブ１１で検出された光音響信号を受信する。ＡＤ変換手段２２は、サンプリング手段であり、受信回路２１が受信した光音響信号をサンプリングしてデジタル信号に変換する。ＡＤ変換手段２２は、例えば、外部から入力する所定周波数のＡＤクロック信号に基づいて、所定のサンプリング周期で光音響信号をサンプリングする。受信メモリ２３は、ＡＤ変換手段２２でサンプリングされた光音響信号を記憶する。

光音響画像再構成手段２４は、受信メモリ２３から光音響信号を読み出し、プローブ１１の複数の超音波振動子で検出された光音響信号に基づいて、光音響画像の各ラインのデータを生成する。光音響画像再構成手段２４は、例えばプローブ１１の６４個の超音波振動子からのデータを、超音波振動子の位置に応じた遅延時間で加算し、１ライン分のデータを生成する（遅延加算法）。光音響画像再構成手段２４は、遅延加算法に代えて、ＢＰ法（Back Projection）により再構成を行ってもよい。あるいは光音響画像再構成手段２４は、ハフ変換法又はフーリエ変換法を用いて再構成を行ってもよい。

光波形微分信号生成手段２５には、パルス光の光強度の時間波形を表す関数が実装されている。光波形微分信号生成手段２５は、その関数を使用して、被検体に照射されたパルス光（典型的にはパルスレーザ光）の光強度の時間波形の微分波形である光波形微分信号を生成する。パルスレーザ光の光強度の時間波形を表す関数は、パルスレーザ光のパルス時間幅を独立変数とする。光波形微分信号生成手段２５は、被検体に照射されたパルスレーザ光のパルス時間幅を関数に適用することで、パルスレーザ光の光強度の時間波形を求め、その波形を微分した微分波形を光波形微分信号として求める。

光波形微分信号生成手段２５は、パルスレーザ光の光強度の時間波形を表す関数に代えて、パルスレーザ光の光強度の時間波形を微分した光微分波形を表す関数を実装していてもよい。その場合、光波形微分信号生成手段２５は、被検体に照射されたパルスレーザ光のパルス時間幅を関数に適用することで、光波形微分信号を求めればよい。

光波形微分信号生成手段２５は、生成した光波形微分信号を光波形微分信号逆畳込み手段２６に入力する。光波形微分信号逆畳込み手段２６には、再構成された光音響信号から、光波形微分信号をデコンボリューションした信号を生成する。光波形微分信号をデコンボリューションすることで、ｔ≠０に再構成した圧力分布から、ｔ＝０に再構成した圧力分布、すなわち吸収分布を求めることができる。光波形微分信号逆畳込み手段２６は、再構成前の光音響信号に対してデコンボリューションを行ってもよい。

検波・対数変換手段２７は、デコンボリューション後の各ラインのデータの包絡線を求め、求めた包絡線を対数変換する。包絡線を求める検波手段としては、ヒルベルト変換や直交検波など従来から用いられている手法を用いることができる。これにより、超音波振動子の固有振動による帯域の影響が除去できる。光音響画像構築手段２８は、対数変換が施された各ラインのデータに基づいて、光音響画像を生成する。光音響画像構築手段２８は、例えば光音響信号（ピーク部分）の時間軸方向の位置を光音響層画像における深さ方向の位置に変換して光音響画像を生成する。

制御手段３０は、超音波ユニット１２内の各部を制御する。トリガ制御回路２９は、光音響画像生成に際して、レーザユニット１３にフラッシュランプトリガ信号を送る。また、フラッシュランプトリガ信号の出力後に、Ｑスイッチトリガ信号を送る。レーザユニット１３は、フラッシュランプ３１とＱスイッチ３２とを含む。レーザユニット１３は、フラッシュランプトリガ信号を受けてフラッシュランプ３１を点灯し、レーザ励起を開始する。レーザユニット１３は、Ｑスイッチトリガ信号が入力されるとＱスイッチをＯＮにし、パルスレーザ光を出射する。トリガ制御回路２９は、被検体に対するレーザ光照射と同期してＡＤ変換手段２２にサンプリングトリガ信号を送り、ＡＤ変換手段２２における光音響信号のサンプリング開始タイミングを制御する。

なお、光波形微分信号逆畳込み手段２６の後段に補正手段を設け、その補正手段が、光波形微分信号がデコンボリューションされた信号から、プローブ１１における超音波振動子の受信角度依存特性の影響を除去する構成としてもよい。補正手段が、受信角度依存特性に加えて、又はこれらに代えて、光波形微分信号がデコンボリューションされた信号から被検体における光の入射光分布の影響を除去することとしてもよい。

図５は、光波形微分信号生成手段２５を示す。光波形微分信号生成手段２５は、光パルス幅推定部２５１と、関数適用部２５２と、光パルス幅情報記憶部２５３とを有する。光パルス幅情報記憶部２５３は、パルスレーザ光に関連した発光条件（レーザ発光条件）と、その発光条件下で被検体に照射されるパルスレーザ光のパルス時間幅（光パルス幅）との対応関係を示す光パルス幅情報を記憶する。レーザ発光条件は、例えばパルスレーザ光の波長、レーザ光源におけるレーザ発光のトリガからレーザ発光までの遅延時間、パルスレーザ光の発光周波数、レーザ光源におけるレーザ結晶の励起エネルギー、レーザ光源における経時変化に関する条件、レーザ光源におけるレーザ結晶温度のうちの少なくとも１つを含む。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、光パルス幅情報と被検体に照射されたパルスレーザ光の発光条件とに基づいて、被検体に照射されたパルスレーザ光の光パルス幅を推定する。

関数適用部２５２は、推定された光パルス幅を、パルスレーザ光の光強度の時間波形、又はその微分波形を表す関数に適用し、光波形微分信号を生成する。パルスレーザ光の光強度の時間波形、又はその微分波形を表す関数には、事前に実測されたパルスレーザ光の時間波形又はその微分波形のフィッティング曲線、ガウス関数、ローレンツ関数、又は擬フォークト関数を用いることができる。あるいは、それら関数を組み合わせたものを用いてもよい。例えば時刻ｔ＝０からある時刻まではガウス関数を用い、ある時刻以降はフィッティング曲線を用いることも可能である。

例えば光パルス幅情報記憶部２５３は、パルスレーザ光の波長と、光パルス幅とを対応付けて記憶する。光パルス幅情報記憶部２５３は、例えば、波長λ_１と光パルス幅Ｔ_ｐ１、波長λ_２と光パルス幅Ｔ_ｐ２、波長λ_３と光パルス幅Ｔ_ｐ３を、それぞれ対応付けたテーブル形式で光パルス幅情報を記憶する。光パルス幅推定部２５１は、例えば制御手段３０（図４）から、被検体に対して照射されたパルスレーザ光の波長に関する情報を入手する。光パルス幅推定部２５１は、例えば被検体に対して照射されたパルスレーザ光の波長がλ_１のときは、光パルス幅をＴ_ｐ１と推定する。この場合、レーザユニット１３（図４）が出射するパルスレーザ光の光パルス幅が波長に応じて変化するときでも、各波長のパルスレーザ光の光パルス幅を正しく推定できる。

光パルス幅情報記憶部２５３は、レーザユニット１３に対するレーザ発光のトリガから実際のレーザ発光までの遅延時間を、光パルス幅と対応付けて記憶してもよい。光パルス幅情報記憶部２５３は、例えば、遅延時間と、各遅延時間の場合の光パルス幅とを対応付けたテーブル形式で光パルス幅情報を記憶する。光パルス幅推定部２５１は、被検体に対してパルスレーザ光が照射されたときに、レーザ発光のトリガから実際のレーザ発光までの遅延時間を入手する。レーザ発光のトリガから実際のレーザ発光までの遅延時間は、例えばレーザユニット１３（図４）に対してＱスイッチトリガ信号が出力されたタイミングから、プローブ１１の超音波検出素子に対してパルスレーザ光が入射し、焦電効果による信号が検出されたタイミングまでの間の時間差で求めることができる。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、入手した遅延時間に対応する光パルス幅を、被検体に照射されたパルスレーザ光のパルス幅として推定する。この場合、遅延時間の差に起因してパルスレーザ光の光パルス幅が変化するようなときでも、被検体に照射されたパルスレーザ光の光パルス幅を正しく推定することができる。

光パルス幅情報記憶部２５３は、パルスレーザ光の発光周波数を、光パルス幅と対応付けて記憶してもよい。光パルス幅情報記憶部２５３は、例えば、発光周波数と、各発光周波数遅延時間の場合の光パルス幅とを対応付けたテーブル形式で光パルス幅情報を記憶する。ここで、パルスレーザ光の発光周波数は、レーザユニット１３から繰り返し定期的にパルスレーザ光を出射する際のレーザ出射の頻度を表す。光パルス幅推定部２５１は、例えば制御手段３０から、パルスレーザ光の発光周波数を取得する。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、取得した発光周波数に対応する光パルス幅を、被検体に照射されたパルスレーザ光のパルス幅として推定する。例えば、レーザユニット１３がパルスレーザ光を出射する頻度が高いときと、頻度が低いときとで、パルスレーザ光の光パルス幅が変化することが考えられる。パルスレーザ光の発光周波数と光パルス幅とを対応付けた光パルス幅情報を用いることで、パルスレーザ光の発光周波数に起因してパルスレーザ光の光パルス幅が変化するようなときでも、被検体に照射されたパルスレーザ光の光パルス幅を正しく推定できる。

光パルス幅情報記憶部２５３は、レーザユニット１３におけるレーザ結晶の励起エネルギーに関連した情報を、光パルス幅と対応付けて記憶してもよい。光パルス幅情報記憶部２５３は、例えば、レーザ結晶の励起エネルギーに関連した情報と、各励起エネルギーの場合の光パルス幅とを対応付けたテーブル形式で光パルス幅情報を記憶する。励起エネルギーに関連した情報としては、レーザユニット１３のフラッシュランプ３１の印加電圧が考えられる。光パルス幅推定部２５１は、被検体に対してパルスレーザ光が照射されたときに、例えばレーザユニット１３から、フラッシュランプ３１に印加された電圧に関する情報を、励起エネルギーに関連した情報として入手する。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、入手した印加電圧（励起エネルギー）に対応する光パルス幅を、被検体に照射されたパルスレーザ光のパルス幅として推定する。レーザユニット１３におけるレーザ結晶の励起エネルギーは、フラッシュランプの印加電圧に応じて変動するものと考えられ、励起エネルギーが変化すると、パルスレーザ光の光パルス幅に変動が生じると考えられる。レーザ結晶の励起エネルギーと光パルス幅とを対応付けた光パルス幅情報を用いることで、励起エネルギーの変化に起因してパルスレーザ光の光パルス幅が変化するようなときでも、被検体に照射されたパルスレーザ光の光パルス幅を正しく推定できる。

光パルス幅情報記憶部２５３は、レーザユニット１３における経時変化に関する条件に関連した情報と、光パルス幅とを対応付けて記憶してもよい。経時変化に関する条件に関連した情報としては、例えばレーザユニット１３におけるフラッシュランプ３１の累積発光回数、累積発光時間などが考えられる。光パルス幅情報記憶部２５３は、例えば、フラッシュランプの累積発光回数／発光時間と、各累積発光回数／発光時間の場合の光パルス幅とを対応付けたテーブル形式で光パルス幅情報を記憶する。光パルス幅推定部２５１は、被検体に対してパルスレーザ光が照射されたときに、例えばレーザユニット１３から、フラッシュランプ３１の累積発光回数又は発光時間を入手する。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、入手した累積発光回数／発光時間に対応する光パルス幅を、被検体に照射されたパルスレーザ光のパルス幅として推定する。フラッシュランプ３１の性能は、発光時間／発光回数が増えるに従って低下していくものと考えられ、それに伴って、レーザ結晶の励起エネルギーも低下していくと考えられる。そして、レーザ結晶の励起エネルギーの低下に伴い、パルスレーザ光の光パルス幅が変化すると考えられる。レーザユニット１３における経時変化に関する条件と光パルス幅とを対応付けた光パルス幅情報を用いることで、経時変化に起因してパルスレーザ光の光パルス幅が変化するようなときでも、被検体に照射されたパルスレーザ光の光パルス幅を正しく推定できる。

光パルス幅情報記憶部２５３は、レーザユニット１３におけるレーザ結晶温度に関連した情報と、光パルス幅とを対応付けて記憶してもよい。レーザ結晶温度に関連した情報としては、例えばレーザユニット１３における固体レーザ結晶の冷却水の温度などが考えられる。光パルス幅情報記憶部２５３は、例えば、固体レーザ結晶の冷却水温度と、各冷却水温度の場合の光パルス幅とを対応付けたテーブル形式で光パルス幅情報を記憶する。レーザユニット１３は、固体レーザ結晶の冷却水の温度をモニターしている。光パルス幅推定部２５１は、被検体に対してパルスレーザ光が照射されたときに、例えばレーザユニット１３から、固体レーザ結晶の冷却水温度を入手する。光パルス幅推定部２５１は、光パルス幅情報記憶部２５３を参照し、入手した冷却水温度に対応する光パルス幅を、被検体に照射されたパルスレーザ光のパルス幅として推定する。この場合、レーザユニット１３が出射するパルスレーザ光の光パルス幅がレーザ結晶温度に応じて変化するときでも、各波長のパルスレーザ光の光パルス幅を正しく推定できる。

なお、上記したレーザ発光条件は、組み合わせて使用することもできる。例えば、パルスレーザ光の波長と発光周波数との２つを組み合わせ、それら２つの組み合わせから光パルス幅を推定してもよい。その場合、光パルス幅情報記憶部２５３には、例えば波長及び発光周波数と、光パルス幅とを対応付けた二次元テーブルを記憶しておけばよい。

また、光パルス幅の推定に際して、被検体に照射されたパルスレーザ光のレーザ発光条件に合致するレーザ発光条件が光パルス幅情報（テーブル）に存在しないときは、テーブルに存在するパルス時間幅に基づいて、補間により、被検体に照射されたパルスレーザ光のパルス時間幅を推定すればよい。例えば、被検体に対して照射されたパルスレーザ光の波長と等しい波長に対するパルス時間幅が光パルス幅情報記憶部２５３に記憶されていないときは、テーブルに存在する前後の波長に対するパルス時間幅に基づいて、被検体に照射されたパルスレーザ光の波長に対するパルス時間幅を補間推定してもよい。具体的に、被検体に対して照射されたパルスレーザ光の波長がλ_１とλ_２との間の波長である場合には、光パルス幅Ｔｐ_１とＴｐ_２との間の値を光パルス幅として推定すればよい。

図６に、光波形微分信号逆畳込み手段２６を示す。光波形微分信号逆畳込み手段２６は、フーリエ変換手段４１、４２と、逆フィルタ演算手段４３と、フィルタ適用手段４４と、フーリエ逆変換手段４５とを有する。フーリエ変換手段（第１のフーリエ変換手段）４１は、離散フーリエ変換により、再構成された光音響信号を時間領域の信号から周波数領域の信号へと変換する。フーリエ変換手段（第２のフーリエ変換手段）４２は、離散フーリエ変換により、光波形微分信号を時間領域の信号から周波数領域の信号へと変換する。フーリエ変換のアルゴリズムにはＦＦＴを用いることができる。

本実施形態においては、光音響信号のサンプリングレートと光波形微分信号のサンプリングレートとは等しいものとする。例えば光音響信号は、Ｆｓ＝４０ＭＨｚのサンプリングクロックに同期してサンプリングされている。光波形微分信号は、Ｆｓ＿ｈ＝４０ＭＨｚのサンプリングレートでサンプリングされた信号に相当する。フーリエ変換手段４１は、４０ＭＨｚでサンプリングされた光音響信号を、例えば１０２４点のフーリエ変換でフーリエ変換する。また、フーリエ変換手段４２は、４０ＭＨｚでサンプリングされた信号に相当する光波形微分信号を１０２４点のフーリエ変換でフーリエ変換する。

逆フィルタ演算手段４３は、フーリエ変換された光波形微分信号の逆数を逆フィルタとして求める。例えば逆フィルタ演算手段４３は、光波形微分信号ｈをフーリエ変換した信号をｆｆｔ＿ｈとしたとき、conj(fft_h)/abs(fft_h)²を逆フィルタとして求める。フィルタ適用手段４４は、フーリエ変換手段４１でフーリエ変換された光音響信号に、逆フィルタ演算手段４３で求められた逆フィルタを適用する。フィルタ適用手段４４は、例えば、要素ごとに、光音響信号のフーリエ係数と逆フィルタのフーリエ係数とを乗算する。逆フィルタが適用されることで、周波数領域の信号において、光波形微分信号がデコンボリューションされる。フーリエ逆変換手段４５は、フーリエ逆変換により、逆フィルタが適用された光音響信号を、周波数領域の信号から時間領域の信号へと変換する。フーリエ逆変換により、時間領域の吸収分布信号が得られる。

図７は、動作手順を示す。トリガ制御回路２９は、レーザユニット１３に対してフラッシュランプトリガ信号を出力する。レーザユニット１３は、フラッシュランプトリガ信号を受けてフラッシュランプ３１を点灯する。トリガ制御回路２９は、所定のタイミングでＱスイッチトリガ信号を出力する。レーザユニット１３は、Ｑスイッチトリガ信号が入力されると、Ｑスイッチ３２をＯＮにし、パルスレーザ光を出射する。出射したパルスレーザ光は、例えばプローブ１１まで導光され、プローブ１１から被検体に照射される（ステップＳ１）。

プローブ１１は、パルスレーザ光の照射後、パルスレーザ光の照射により被検体内で発生した光音響信号を検出する（ステップＳ２）。超音波ユニット１２の受信回路２１は、プローブ１１で検出された光音響信号を受信する。トリガ制御回路２９は、被検体に対する光照射のタイミングに合わせてＡＤ変換手段２２にサンプリングトリガ信号を送る。ＡＤ変換手段２２は、サンプリングトリガ信号を受けて光音響信号のサンプリングを開始し、光音響信号のサンプリングデータを受信メモリ２３に格納する。光音響画像再構成手段２４は、受信メモリ２３から光音響信号のサンプリングデータを読み出し、読み出した光音響信号のサンプリングデータに基づいて、光音響信号を再構成する（ステップＳ３）。

光波形微分信号生成手段２５は、例えばレーザ発光条件に基づいて光パルス幅を推定し、推定した光パルス幅を、パルスレーザ光の光強度の時間波形又はその微分波形を表す関数に適用することで、光波形微分信号を生成する（ステップＳ４）。光波形微分信号逆畳込み手段２６は、ステップＳ３で再構成された光音響信号から、ステップＳ４で生成された光波形微分信号をデコンボリューションする（ステップＳ５）。このデコンボリューションにより、吸収分布を示す光音響信号が得られる。

検波・対数変換手段２７は、デコンボリューションが行われた光音響信号の包絡線を求め、求めた包絡線を対数変換する。光音響画像構築手段２８は、対数変換が施された各ラインのデータに基づいて、光音響画像を生成する（ステップＳ６）。この光音響信号は、吸収分布を画像化した吸収分布画像である。画像表示手段１４は、表示画面上に、吸収分布画像である光音響画像を表示する（ステップＳ７）。

本実施形態では、ひとまず、光音響画像再構成手段２４にて、通常の再構成法により発光時刻（ｔ＝０）の圧力分布として光音響信号（光音響画像）を再構成する。次に、光の発光時間は実際には有限の長さを持つことから、再構成時にｔ＝０としていた時刻を有限の長さの時間と考え、光波形微分信号逆畳込み手段２６にて、再構成後の光音響画像から光波形微分信号をデコンボリューションする。光波形微分信号をデコンボリューションすることで、吸収分布を得ることができ、吸収分布画像を生成することができる。このような手法を採用することで、実用的な光パルス幅と実用的な超音波システム、或いは実際の生体を観測した場合でも、吸収分布を画像化することができる。これは、現状システムの検出器の帯域やＡＤサンプリングを使用できる利点がある。また、本実施形態においては光音響画像の再構成で圧力分布を一度出しているため、既存の超音波アルゴリズム、装置との親和性が高い。

本実施形態では、光波形微分信号生成手段２５に、光パルスの光強度の時間波形又はその微分波形を表し、光パルス幅を独立変数とする関数を実装しており、その関数を用いて被検体に照射されたパルスレーザ光の光強度の時間波形を微分した光波形微分信号を生成している。関数を用いて光波形微分信号を生成することで、実際に被検体に対して照射されたパルスレーザ光の時間波形の微分波形に近い光波形微分信号を生成することができる。そのように生成された光波形微分信号をデコンボリューションに用いることで、固定的な光波形微分信号をデコンボリューションする場合に比して、より正確に吸収分布を得ることができる。

また、本実施形態では、レーザ発光条件と光パルス幅とを対応付けた光パルス幅情報を用いて光パルス幅を推定している。光パルス幅情報を用い、レーザ発光条件に基づいて被検体に対して照射されたパルスレーザ光の光パルス幅を推定することで、例えばレーザ波長、レーザ発光周波数、発光のトリガから実際のレーザ発光までの遅延時間などのレーザ発光条件の変化に起因して変化するパルスレーザ光の光パルス幅を正しく推定することができる。そのように推定された光パルス幅を関数に適用して光波形微分信号を生成することで、デコンボリューションに用いられる光波形微分信号を、被検体に対して実際に照射されたパルスレーザ光の時間波形の微分波形に近づけることができる。

次いで、本発明の第２実施形態を説明する。第１実施形態では、光音響信号のサンプリングレートと光波形微分信号のサンプリングレートとが一致しており、双方の信号を同じデータ点数でフーリエ変換した。本実施形態では、光音響信号を低速サンプリングする一方で、光パルス波形微分信号を高速サンプリングの信号とする。つまり、光波形微分信号を、光音響信号のサンプリングレートよりも高いサンプリングレートでサンプリングされた信号に相当する信号とする。例えば光音響信号のサンプリング間隔（サンプリングレートの逆数）は、被検体に照射される光のパルス時間幅よりも長く設定される。フーリエ変換に際しては、低サンプリングレートの光音響信号を、光波形微分信号のサンプリングレートと同じサンプリングレートでリサンプル（アップサンプル）した上で、フーリエ変換を行う。その他の点は、第１実施形態と同様でよい。

図８は、本実施形態における光波形微分信号逆畳込み手段２６ａを示す。本実施形態における光波形微分信号逆畳込み手段２６ａは、図６に示す第１実施形態における光波形微分信号逆畳込み手段２６の構成に加えて、リサンプル手段４６及び４７を有する。リサンプル手段４６は、アップサンプル手段であり、低いサンプリングレートでサンプリングされた光音響信号のサンプリングデータを、光波形微分信号のサンプリングレートと同じサンプリングレートでアップサンプルする。リサンプル手段４６は、例えば、低サンプリングレートでサンプリングされた光音響信号のサンプル点間にゼロを付加し、アップサンプル前のナイキスト周波数でカットするローパスフィルタをかけることでアップサンプルを行う。

例えば、ＡＤ変換手段２２（図４）における光音響信号のサンプリングレート（第１のサンプリングレート）が４０ＭＨｚであり、光波形微分信号のサンプリングレート（第２のサンプリングレート）が４００ＭＨｚであったとする。この場合、リサンプル手段４６は、４０ＭＨｚの光音響信号を４００ＭＨｚの信号にアップサンプルする。フーリエ変換手段４１は、リサンプル手段４６でアップサンプルされた光音響信号をフーリエ変換する。光音響信号をフーリエ変換するフーリエ変換手段４１と、光波形微分信号をフーリエ変換するフーリエ変換手段４２とは、同じデータ点数でフーリエ変換を行う。例えばフーリエ変換手段４１は光音響信号を８１９２点の周波数領域の信号に変換し、フーリエ変換手段４２は光波形微分信号を８１９２点の周波数領域の信号に変換する。

フィルタ適用手段４４は、アップサンプルされた光音響信号をフーリエ変換した信号に対して逆フィルタを適用する。フーリエ逆変換手段４５は、逆フィルタが適用された信号を、周波数領域の信号から時間領域の信号（吸収分布）へと変換する。時間領域の信号に戻された吸収分布信号は、例えば４００ＭＨｚにアップサンプルされた状態の信号となっている。リサンプル手段４７は、吸収分布信号が、光音響信号の元のサンプルリングレートでサンプリングされた信号となるように、吸収信号をダウンサンプルする。リサンプル手段４７は、例えば４００ＭＨｚの吸収信号を４０ＭＨｚの吸収信号にダウンサンプルする。ダウンサンプリングは、例えばダウンサンプル後のナイキスト周波数でカットするローパスフィルタをかけた後に、サンプル点を間引くことで行う。

図９Ａに、サンプリングレート４００ＭＨｚに相当する光波形微分信号を示し、図９Ｂに、サンプリングレート４０ＭＨｚに相当する光波形微分信号を示す。サンプリングレート４００ＭＨｚでは、図９Ａに示すように、パルスレーザ光の光強度の時間波形を微分した波形を正確に再現できる。一方、光波形微分信号のサンプリングレートを光音響信号のサンプリングレートに合わせ、４０ＭＨｚ相当の信号とすると、図９Ｂに示すように、光微分波形を正確に再現できなくなる。

フィルタ適用手段４４にて光音響信号をフーリエ変換した信号に逆フィルタを適用する際には、双方のデータ点数が揃っている必要がある。光音響信号のサンプリングレートに合わせて光波形微分信号のサンプリングレートを設定すると、図９Ｂに示したように、波形変化に対してサンプリング周波数が低すぎ、光微分波形が正確に再現できない。このような光波形微分信号から求めた逆フィルタを適用した場合、光パルス微分項を正確にデコンボリューションできずに、吸収分布を正しく求められないこともある。

一方、光微分波形を正確に再現するために光波形微分信号を４００ＭＨｚ相当の信号とし、光音響信号のサンプリングレートを４００ＭＨｚに合わせるとした場合は、光パルス微分項を正確にデコンボリューションでき、吸収分布を正しく求めることができる。しかしながら、その場合、ＡＤ変換手段２２には高速なＡＤ変換器が要求され、また、サンプリングデータの総数が増えることから、受信メモリ２３（図４）に要求されるメモリ容量が増大する。更に、光音響画像再構成手段２４で取り扱うデータが増えるため、再構成に要する時間も長くなる。

本実施形態では、リサンプル手段４６で、事後的に光音響信号のサンプリングデータをリサンプルする。本実施形態では、検出後の光音響信号を信号処理でアップサンプルしているため、光音響の検出から再構成までは低速サンプリングしつつも、光パルス微分項を正確にデコンボリューションすることができる。本実施形態では、ＡＤ変換手段２２に高速なＡＤ変換器は不要であり、受信メモリ２３に必要なメモリ容量も増大しない。また、光音響信号の再構成に要する時間も増大せず、光音響信号の検出時に高いサンプリングレートでサンプリングする場合に比して、処理時間を短縮することができる。

続いて、本発明の第３実施形態を説明する。本実施形態においても、第２実施形態と同様に、光波形微分信号のサンプリングレートを光音響信号のサンプリングレートよりも高く設定する。第２実施形態では、低サンプリングレートでサンプリングされた光音響信号をアップサンプルし、双方の信号を同じデータ点数でフーリエ変換した。本実施形態では、光波形微分信号のフーリエ変換を、光音響信号のフーリエ変換のデータ点数よりも多いデータ点数で行い、フーリエ変換された光音響信号に対して、データ点数の差の分だけ中央（高周波成分領域）にゼロ点を付加する。その他の点は、第１実施形態と同様でよい。

図１０は、本実施形態における光波形微分信号逆畳込み手段２６ｂを示す。本実施形態における光波形微分信号逆畳込み手段２６ｂは、図６に示す第１実施形態における光波形微分信号逆畳込み手段２６の構成に加えて、ゼロパディング手段４８とゼロ点除去手段４９とを有する。例えば、光音響信号のサンプリングレート（第１のサンプリングレート）は４０ＭＨｚであり、光波形微分信号のサンプリングレート（第２のサンプリングレート）は３２０ＭＨｚであるとする。フーリエ変換手段４１は、例えば４０ＭＨｚの光音響信号を１０２４点（第１のデータ点数）の周波数領域の信号に変換し、フーリエ変換手段４２は、３２０ＭＨｚの光波形微分信号を８１９２点（第２のデータ点数）の周波数領域の信号に変換する。第２のデータ点数は、第１のデータ点数に、第２のサンプリングレートと第１のサンプリングレートとの比を乗じたデータ点数と等しいか、又はそれよりも多い。

ゼロパディング手段４８は、フーリエ変換手段４１から周波数領域の信号に変換された光音響信号を入力する。ゼロパディング手段４８は、フーリエ変換された光音響信号に対して、フーリエ変換後の光音響信号と光波形微分信号のデータ点数の差の分だけ中央にゼロ点（信号値ゼロの点）を付加する。ゼロパディング手段４８は、例えばデータ点数１０２４点の光音響信号（周波数領域）を、周波数帯域の中心周波数で２つに分割し、分割した２つの周波数領域の間にデータ点数の差の分だけゼロ点を付加し、光波形微分信号（周波数領域）のデータ点数と同じデータ点数８１９２点の光音響信号を生成する。ゼロ点の付加は、周波数領域におけるアップサンプリングに相当する。

フィルタ適用手段４４は、ゼロパディング手段４８でゼロパディングが施された信号に対して逆フィルタを適用する。ゼロ点除去手段４９は、逆フィルタが適用された信号からゼロパディング手段４８で“０”が付加された周波数帯域を除去する。例えばゼロパディング手段４８にてデータ点数１０２４点の光音響信号（周波数領域）がデータ点数８１９２点の信号に変換されていたとき、ゼロ点除去手段４９は、フィルタ適用後の信号（データ点数８１９２点）をデータ点数１０２４点の信号に戻す。ゼロ点の除去は、周波数領域におけるダウンサンプリングに相当する。フーリエ逆変換手段４５は、データ点数１０２４点に戻された信号を、周波数領域の信号から時間領域の信号へと変換する。

図１１Ａに、フーリエ変換された光音響信号を示し、図１１Ｂに、ゼロパディング後の光音響信号を示す。例えば、ＡＤ変換手段２２（図４）における光音響信号のサンプリングレートが４０ＭＨｚであるとき、その光音響信号をフーリエ変換した信号は、図１１Ａに示すように、０ＭＨｚから４０ＭＨｚまでの周波数帯域の信号となる。この信号を、ナイキスト周波数（サンプリング周波数の１／２）である２０ＭＨｚを境に２つの領域Ａ、Ｂに２分割する。ゼロパディング手段４８は、図１１Ｂに示すように、２つの領域の間にゼロ点を８１９２−１０２４＝７１６８個挿入する。ゼロ点が付加された結果、領域Ｂの信号は、３００ＭＨｚから３２０ＭＨｚの周波数領域に対応した信号となる。

本実施形態では、低サンプリングレートでサンプリングされた光音響信号を周波数領域の信号に変換し、変換された周波数領域の信号の高周波成分の領域のゼロ点を付加する。本実施形態と第２実施形態との相違点は、第２実施形態では、光音響信号をアップサンプルするのに対し、本実施形態では、光音響信号を周波数領域でアップサンプルする点である。時間領域に代え、周波数領域において、双方の信号の帯域差を埋めるようにリサンプル（アップサンプル）を行う場合も、第２実施形態と同様に、光音響の検出から再構成までは低速サンプリングしつつも、光パルス微分項を正確にデコンボリューションすることができる。

引き続いて、本発明の第４実施形態を説明する。本実施形態においても、第２及び第３実施形態と同様に、光波形微分信号のサンプリングレートを光音響信号のサンプリングレートよりも高く設定する。本実施形態では、光波形微分信号のフーリエ変換を、光音響信号のフーリエ変換のデータ点数よりも多いデータ点数で行い、フーリエ変換された光波形微分信号から高周波成分サンプル点を除去し、その逆数を逆フィルタとして求める。その他の点は、第１実施形態と同様でよい。

図１２は、本実施形態における光波形微分信号逆畳込み手段２６ｃを示す。本実施形態における光波形微分信号逆畳込み手段２６ｃは、図６に示す第１実施形態における光波形微分信号逆畳込み手段２６の構成に加えて、高周波成分サンプル点除去手段５０を有する。例えば、光音響信号のサンプリングレート（第１のサンプリングレート）は４０ＭＨｚであり、光波形微分信号のサンプリングレート（第２のサンプリングレート）は３２０ＭＨｚであるとする。フーリエ変換手段４１は、例えば４０ＭＨｚの光音響信号を１０２４点（第１のデータ点数）の周波数領域の信号に変換し、フーリエ変換手段４２は、３２０ＭＨｚ相当の光波形微分信号を８１９２点（第２のデータ点数）の周波数領域の信号に変換する。第２のデータ点数は、第１のデータ点数に、第２のサンプリングレートと第１のサンプリングレートとの比を乗じたデータ点数と等しいか、又はそれよりも多い。

高周波成分サンプル点除去手段５０は、フーリエ変換手段４２から周波数領域の信号に変換された光波形微分信号を入力する。高周波成分サンプル点除去手段５０は、フーリエ変換された光波形微分信号から、フーリエ変換後の光音響信号と光波形微分信号のデータ点数の差の分だけ高周波成分サンプル点を除去する。高周波成分サンプル点除去手段５０は、例えばデータ点数８１９２点の光波形微分信号（周波数領域）から高周波成分に相当する中央のデータ点を削除し、光音響信号（周波数領域）のデータ点数と同じデータ点数１０２４点の光波形微分信号を生成する。高周波成分サンプル点の除去は、周波数領域における光波形微分信号のダウンサンプリングに相当する。

図１３Ａに、フーリエ変換された光波形微分信号を示し、図１３Ｂに、高周波成分サンプル点が除去された光波形微分信号を示す。例えば、光波形微分信号のサンプリングレートが３２０ＭＨｚであるとき、その光波形微分信号をフーリエ変換した信号（データ点数８１９２点）は、図１３Ａに示すように、０ＭＨｚから３２０ＭＨｚまでの周波数帯域の信号となる。この信号を、１番目のデータ点から５１２番目までの領域（領域Ａ）、５１３番目のデータ点から７６８０番目のデータ点までの領域（領域Ｂ）、及び、７６８１番目のデータ点から８１９２番目のデータ点までの領域（領域Ｃ）の３つの領域に分け、領域Ｂのデータ点を除去する。図１３Ｂに示すように、領域Ａと領域Ｃとをつなげることで、０ＭＨｚから４０ＭＨｚまでの周波数帯域に対応したデータ点数１０２４点の光波形微分信号（周波数領域）が得られる。

逆フィルタ演算手段４３は、高周波成分サンプル点が除去された光波形微分信号（周波数領域）の逆数を逆フィルタとして求める。逆フィルタ演算手段４３は、例えばデータ点が８１９２点から１０２４点に削減された光波形微分信号の逆数を逆フィルタとして求める。フィルタ適用手段４４は、例えばデータ点数１０２４点の光音響信号（周波数領域）と逆フィルタとを要素ごとに乗算する。フーリエ逆変換手段４５は、逆フィルタが適用された信号を、周波数領域の信号から時間領域の信号へと変換する。

ここで、第３実施形態では、フィルタ適用手段４４は、図１１Ｂに示す高周波成分の領域にゼロ点が付加された光音響信号（周波数領域）と、図１３Ａに示す光波形微分信号（周波数領域）の逆数とを乗算する。光音響信号の高周波成分領域の値は“０”であるため、光波形微分信号の高周波成分（図１３Ａの領域Ｂ）は、逆フィルタ適用後の光音響信号に影響を与えない。従って、本実施形態のように、光波形微分信号の周波数領域の信号から高周波成分サンプル点を除去し、高周波成分を除去した光波形微分信号から逆フィルタを求め、求めた逆フィルタを光音響信号（周波数領域）に適用しても、得られる結果は第３実施形態と同じ結果となる。つまり、本実施形態においても、第３実施形態と同様な効果が得られる。

続いて、本発明の第５実施形態を説明する。図１４は、本発明の第５実施形態の光音響画像生成装置を示す。本実施形態の光音響画像生成装置１０ａにおける超音波ユニット１２ａは、図４に示す第１実施形態の光音響画像生成装置１０における超音波ユニット１２の構成に加えて、送信制御回路３３、データ分離手段３４、超音波画像再構成手段３５、検波・対数変換手段３６、超音波画像構築手段３７、及び画像合成手段３８を備える。本実施形態の光音響画像生成装置１０ａは、光音響画像に加えて、超音波画像の生成を行う点で、第１実施形態の光音響画像生成装置１０と相違する。なお、本実施形態においては音響波として超音波を用いるが、被検対象や測定条件等に応じて適切な周波数を選択することにより、可聴周波数の音響波であっても良い。また、本実施形態を第２から第４実施形態の何れかと組み合わせ、それら実施形態において超音波画像の生成を行うこととしてもよい。

本実施形態では、プローブ１１は、光音響信号の検出に加えて、被検体に対する音響波（超音波）の出力（送信）、及び送信した超音波に対する被検体からの反射超音波の検出（受信）を行う。トリガ制御回路２９は、超音波画像（反射音響波画像）の生成時は、送信制御回路３３に超音波送信を指示する旨の超音波送信トリガ信号を送る。送信制御回路３３は、トリガ信号を受けると、プローブ１１から超音波を送信させる。プローブ１１は、超音波の送信後、被検体からの反射超音波を検出する。超音波の送受信は分離してもよい。例えばプローブ１１とは異なる位置から超音波の送信を行い、その送信された超音波に対する反射超音波をプローブ１１で受信してもよい。

プローブ１１が検出した反射超音波は、受信回路２１を介してＡＤ変換手段２２に入力される。トリガ制御回路２９は、超音波送信のタイミングに合わせてＡＤ変換手段２２にサンプリグトリガ信号を送り、反射超音波のサンプリングを開始させる。ここで、反射超音波はプローブ１１と超音波反射位置との間を往復するのに対し、光音響信号はその発生位置からプローブ１１までの片道である。反射超音波の検出には、同じ深さ位置で生じた光音響信号の検出に比して２倍の時間がかかるため、ＡＤ変換手段２２のサンプリングクロックは、光音響信号サンプリング時の半分、例えば２０ＭＨｚとしてもよい。ＡＤ変換手段２２は、反射超音波のサンプリングデータを受信メモリ２３に格納する。光音響信号の検出（サンプリング）と、反射超音波の検出（サンプリング）とは、どちらを先に行ってもよい。

データ分離手段３４は、受信メモリ２３に格納された光音響信号のサンプリングデータと反射超音波のサンプリングデータとを分離する。データ分離手段３４は、分離した光音響信号のサンプリングデータを光音響画像再構成手段２４に入力する。光波形微分信号の生成及びデコンボリューションを含む光音響画像（吸収分布画像）の生成は、第１実施形態と同様である。データ分離手段３４は、分離した反射超音波のサンプリングデータを、超音波画像再構成手段３５に入力する。

超音波画像再構成手段３５は、プローブ１１の複数の超音波振動子で検出された反射超音波（そのサンプリングデータ）に基づいて、超音波画像の各ラインのデータを生成する。各ラインのデータの生成には、光音響画像再構成手段２４における各ラインのデータの生成と同様に、遅延加算法などを用いることができる。検波・対数変換手段３６は、超音波画像再構成手段３５が出力する各ラインのデータの包絡線を求め、求めた包絡線を対数変換する。

超音波画像構築手段３７は、対数変換が施された各ラインのデータに基づいて、超音波画像を生成する。超音波画像再構成手段３５、検波・対数変換手段３６、及び超音波画像構築手段３７は、反射超音波に基づいて超音波画像を生成する超音波画像生成手段（反射音響波画像生成手段）を構成する。画像合成手段３８は、光音響画像と超音波画像とを合成する。画像合成手段３８は、例えば光音響画像と超音波画像とを重畳することで画像合成を行う。合成された画像は、画像表示手段１４に表示される。画像合成を行わずに、画像表示手段１４に、光音響画像と超音波画像とを並べて表示し、或いは光音響画像と超音波画像とを切り替えてすることも可能である。

本実施形態では、光音響画像生成装置は、光音響画像に加えて超音波画像を生成する。超音波画像を参照することで、光音響画像では画像化することができない部分を観察することができる。被検体に照射されたパルスレーザ光の微分波形を表す光波形微分信号を生成し、観測波形から光波形微分信号をデコンボリューションすることで吸収分布を画像化できる点は、第１実施形態と同様である。また、超音波画像の生成と光音響画像の生成とで、画像再構成や検波・対数変換などのアルゴリズムの大部分を共通化でき、ＦＰＧＡ回路構成やソフトの簡略化が可能であるなどの実用上のメリットを有する。

引き続いて、本発明の第６実施形態を説明する。図１５は、本発明の第６実施形態の光音響画像生成装置を示す。本実施形態は、被検体に対して複数の波長の光を照射する点で、第１実施形態と相違する。本実施形態の光音響画像生成装置１０ｂにおける超音波ユニット１２ｂは、図４に示す第１実施形態の光音響画像生成装置１０における超音波ユニット１２の構成に加えて、複数波長の光に対する光音響信号（光音響画像）の演算を行う２波長データ演算手段４０を備える。なお、本実施形態を第２から第５実施形態の何れかと組み合わせ、それら実施形態において複数波長の光を照射し、複数波長に対する光音響信号（光音響画像）の演算することとしてもよい。

本実施形態では、レーザユニット１３は、複数の波長の光を切り替えて出射可能に構成されている。レーザユニット１３は、例えば波長７５０ｎｍのパルスレーザ光と波長８００ｎｍのパルスレーザ光を切り替えて出射する。プローブ１１は、各波長のパルスレーザ光の出射後に被検体からの光音響信号を検出し、受信メモリ２３には、各波長に対応した光音響信号のサンプリングデータが格納される。格納された各波長に対応した光音響信号は、それぞれ光音響画像再構成手段で再構成される。

光波形微分信号逆畳込み手段２６は、光音響画像再構成手段２４による再構成後、各波長に対応した光音響信号（光音響画像）から、光波形微分信号生成手段２５が生成した各波長に対応した光波形微分信号をそれぞれデコンボリューションする。各波長に対応した光波形微分信号がデコンボリューションされた光音響信号は、２波長データ演算手段４０にて処理される。

ここで、生体組織の多くは光吸収特性が光の波長に応じて変わり、また一般に、その光吸収特性も組織ごとに特有のものとなっている。例えば、ヒトの動脈に多く含まれる酸素化ヘモグロビン（酸素と結合したヘモグロビン:oxy-Hb）の波長７５０ｎｍにおける分子吸収係数は波長８００ｎｍにおけるそれも低い。また、静脈に多く含まれる脱酸素化ヘモグロビン（酸素と結合していないヘモグロビンdeoxy-Hb)の波長７５０ｎｍにおける分子吸収係数は波長８００ｎｍにおけるそれよりも高い。この性質を利用し、波長８００ｎｍで得られた光音響信号に対して、波長７５０ｎｍで得られた光音響信号が相対的に大きいのか小さいのかを調べることで、動脈からの光音響信号と静脈からの光音響信号とを判別することができる。

２波長データ演算手段４０は、例えば複数波長に対応した光音響信号間の相対的な大小関係を比較する。具体的には、２波長データ演算手段４０は、波長７５０ｎｍの光が照射されたときに検出された光音響信号と、波長８００ｎｍの光が照射されたときに検出された光音響信号とを比較し、どちらがどれだけ大きいかを調べる。画像表示に際しては、波長７５０ｎｍの光が照射されたときに検出された光音響信号が大きければ静脈からの光音響信号と判断できるため、その部分を青色で表示するとよい。また、波長８００ｎｍの光が照射されたときに検出された光音響信号が大きければ動脈からの光音響信号と判断できるため、その部分を赤色で表示するとよい。

本実施形態では、２波長データ演算手段４０は、光波形微分信号のデコンボリューション後に、複数波長に対応した光音響信号の演算を行う。複数の波長の光を被検体に照射する場合、例えば１つ目の波長の光を被検体に照射した光音響信号を検出した後、２つ目の光を被検体に照射して光音響信号を検出するとき、体動などの影響で、波長間で位置ずれが生じることがある。複数波長間での光音響信号を比較する際には、同じ場所から発生した光音響信号同士を比較することが好ましい。光波形微分信号をデコンボリューションしない場合は、図３Ａに示したように１本の血管が二重に表示されるなど、画像判定上、血管の位置が確認しづらく、位置ずれ補正がしにくい。光波形微分信号をデコンボリューションすることで、図３Ｂに示したように光吸収分布を画像化でき、血管の位置が確認しやすくなり、位置ずれ補正も容易となる。

また、図３Ａに示したように１本の血管が二重に表示される場合は血管の内部に相当する部分には信号が存在せず（信号レベルが所定レベルよりも低い）、複数波長の画像間で位置ずれが生じていると、信号が存在する部分すなわち血管の重複部分が少なくなる。この場合、複数波長間で光音響信号の相対的な大小関係を適切に比較することが難しくなる。これに対し、図３Ｂに示したように光吸収分布を画像化した場合は、血管内の多くの部分が信号が存在する部分となり、多少位置ずれが生じたとしても、複数波長の画像間で多くの部分が重複する。従って、本実施形態において、双方の画像を位置合わせずに比較する場合でも、位置ずれの影響を軽減することができる。

なお、上記各実施形態では、光音響信号及び光波形微分信号を周波数領域の信号に変換し、周波数領域でデコンボリューション後に時間領域の信号に戻しているが、これには限定されない。光波形微分信号のデコンボリューションを時間領域で行うことも可能である。また、光波形微分信号逆畳込み手段２６は、デコンボリューション時に、光音響信号に対して何らかのフィルタをかける処理を行ってもよい。例えば光波形微分信号逆畳込み手段２６が、デコンボリューション時に、ノイズ増幅周波数帯をフィルタリングするようにしてもよい。あるいは、光波形微分信号に装置応答関数（圧電素子や電気系の周波数フィルタ）を掛け合わせたものをデコンボリューションしてもよい。

上記各実施形態では、光音響信号から光波形微分信号をデコンボリューションした後に光音響画像（吸収分布画像）を生成することとしているが、これに加えて、又はこれに代えて、光波形微分信号をデコンボリューションせずに光音響画像（圧力分布画像）を生成してもよい。例えば、ユーザが、スイッチや表示モニタ上で操作を行うことで、デコンボリューション処理の有無を選択できるようにしておき、ユーザがデコンボリューション処理の実施を選択したときには光波形微分信号のデコンボリューションを行った上で光音響画像を生成し、ユーザがデコンボリューション処理の不実施を選択したときは光波形微分信号のデコンボリューションを行わずに光音響画像を生成してもよい。例えば、光波形微分信号のデコンボリューションを行ったときは、光音響信号を赤・黒の色に対応付けて表示し、デコンボリューションなしのときは、光音響信号を青・黒の色に対応付けて表示してもよい。

また、デコンボリューションなしの場合の光音響画像を生成し、コンピュータがその光音響画像を解析することで、血管部分が２本に分かれているか否かを判定し、血管が２本に分かれていると判定されたときに、その血管部分のみを対象に光波形微分信号のデコンボリューション処理を行うようにしてもよい。その際、デコンボリューション処理を実施した血管部分の表示色を、他の未処理の血管部分の表示色とは異なる色とし、デコンボリューション処理が行われた血管と、他の未処理の血管とが容易に判別可能になるようにしてもよい。

また、上記各実施形態では、レーザ発光条件に基づいて被検体に照射されたパルスレーザ光の光パルス幅を推定したが、これには限定されない。例えば、光パルス幅を測定し、実測された光パルス幅を関数に適用することで光波形微分信号を生成してもよい。図１６は、変形例の光音響画像生成装置を示す。変形例の光音響画像生成装置１０ｃにおける超音波ユニット１２ｃは、図４に示す超音波ユニット１２の構成に加えて光パルス幅測定手段３９を更に備える。光パルス幅測定手段３９は、被検体に対して照射されたパルスレーザ光のパルス時間幅を測定する。光波形微分信号生成手段２５は、推定された光パルス幅に代えて、光パルス幅測定手段３９によって実測された光パルス幅を関数に適用し、光波形微分信号を生成する。このように、実測された光パルス幅と関数とにより光パルス波形微分信号を生成するようにしてもよい。

光パルス幅測定手段３９は、例えば、光を検出するフォトディテクタと、フォトディテクタが光を検出する時間の長さを測定する演算部とを含む。例えば、レーザユニット１３（図４）の光出射端から、被検体に対してパルスレーザ光が照射されるまでの間の光路（例えば光ファイバ）には、被検体に照射されるべきパルスレーザ光の一部を分岐する分岐部が設けられている。フォトディテクタは、例えば分岐部で分岐されたパルスレーザ光の一部が所定の光強度以上のとき、その旨を示す電気信号を出力する。演算部は、フォトディテクタの出力信号に基づいて、フォトディテクタが検出する光の光強度がしきい値以上となっている時間を、光パルス幅として測定する。

上記の分岐部には、例えば透過率が高いミラー、具体的には透過率９５％以上のミラーを用いることができる。そのようなミラーを、例えば被検体方向に向かうレーザ本流に対して４５°の角度で配置し、ミラーの反射成分を光パルス幅測定手段３９のフォトディテクタの方向に分岐する。分岐部に、ミラーに代えて透明ガラスを用い、透明ガラスの反射成分を光パルス幅測定手段３９のフォトディテクタの方向に分岐するようにしてもよい。また、別例として、例えば多数の光ファイバを束ねたバンドルファイバなどを用いてレーザ光を導光しているような場合であれば、複数本の光ファイバのうちのいくつか（例えば一本）を用いて、光パルス幅測定手段３９のフォトディテクタにパルスレーザ光の一部を導光してもよい。

以上、本発明をその好適な実施形態に基づいて説明したが、本発明の光音響信号処理装置及び方法は、上記実施形態にのみ限定されるものではなく、上記実施形態の構成から種々の修正及び変更を施したものも、本発明の範囲に含まれる。

１０：光音響画像生成装置（光音響画像診断装置）
１１：プローブ
１２：超音波ユニット
１３：レーザユニット
１４：画像表示手段
２１：受信回路
２２：ＡＤ変換手段
２３：受信メモリ
２４：光音響画像再構成手段
２５：光波形微分信号生成手段
２６：光波形微分信号逆畳込み手段
２７：検波・対数変換手段
２８：光音響画像構築手段
２９：トリガ制御回路
３０：制御手段
３１：フラッシュランプ
３２：Ｑスイッチ
３３：送信制御回路
３４：データ分離手段
３５：超音波画像再構成手段
３６：検波・対数変換手段
３７：超音波画像構築手段
３８：画像合成手段
３９：光パルス幅測定手段
４０：２波長データ演算手段
４１、４２：フーリエ変換手段
４３：逆フィルタ演算手段
４４：フィルタ適用手段
４５：フーリエ逆変換手段
４６、４７：リサンプル手段
４８：ゼロパディング手段
４９：ゼロ点除去手段
５０：高周波成分サンプル点除去手段
２５１：光パルス幅推定部
２５２：関数適用部
２５３：光パルス幅情報記憶部

Claims (21)

1. 被被検体に対するパルス光の出射後に、該パルス光より被検体内で生じた光音響信号をサンプリングするサンプリング手段と、
   前記パルス光の光強度の時間波形又は該時間波形を微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、前記被検体に対して出射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成する光波形微分信号生成手段と、
   前記サンプリングされた光音響信号から前記光波形微分信号生成手段で生成された光波形微分信号をデコンボリューションした信号を生成する光波形微分信号逆畳込み手段とを備えたことを特徴とする光音響信号処理装置。
2. 前記光波形微分信号生成手段が、
   パルス光に関連した発光条件と、該発光条件下で被検体に出射されるパルス光のパルス時間幅との対応関係を示す光パルス幅情報を参照し、該光パルス幅情報と前記被検体に出射されたパルス光の発光条件とに基づいて、前記被検体に出射されたパルス光のパルス時間幅を推定する光パルス幅推定部と、
   前記推定されたパルス時間幅を前記関数に適用する関数適用部とを含むことを特徴とする請求項１に記載の光音響信号処理装置。
3. 前記光パルス幅情報が、前記発光条件と前記パルス時間幅とを対応付けたテーブルを含むことを特徴とする請求項２に記載の光音響信号処理装置。
4. 前記光パルス幅推定部が、前記被検体に出射されたパルス光の発光条件に合致する発光条件が前記テーブルに存在しないときは、前記テーブルに存在する前記パルス時間幅に基づいて補間により前記被検体に出射されたパルス光のパルス時間幅を推定するものであることを特徴とする請求項３に記載の光音響信号処理装置。
5. 前記パルス光がパルスレーザ光であることを特徴とする請求項２から４何れかに記載の光音響信号処理装置。
6. 前記発光条件が、前記パルスレーザ光の波長、レーザ発光のトリガからレーザ発光までの遅延時間、前記パルスレーザ光の発光周波数、レーザ光源におけるレーザ結晶の励起エネルギー、レーザ光源における経時変化に関する条件、レーザ光源におけるレーザ結晶温度のうちの少なくとも１つを含むことを特徴とする請求項５に記載の光音響画像信号処理装置。
7. 前記関数が、事前に実測されたパルス光の時間波形のフィッティング曲線、ガウス関数、ローレンツ関数、若しくは擬フォークト関数、又は、これら関数を組み合わせたものであることを特徴とする請求項１から６何れかに記載の光音響信号処理装置。
8. 前記被検体に対して出射されたパルス光のパルス時間幅を測定する光パルス幅測定手段を更に備え、前記光波形微分信号生成手段が、前記光パルス幅測定手段で測定されたパルス時間幅を前記関数に適用することを特徴とする請求項１に記載の光音響信号処理装置。
9. 被検体に出射されるパルス光の光源の光出射端から前記被検体に対してパルス光が出射されるまでの間の光路中に、前記被検体に出射されるパルス光の一部を分岐する分岐部が設けられており、前記光パルス幅測定手段が、前記光路から分岐されたパルス光の光強度に基づいて前記被検体に対して出射されたパルス光のパルス時間幅を測定するものであることを特徴とする請求項８に記載の光音響信号処理装置。
10. 前記光波形微分信号逆畳込み手段が、
    前記サンプリングされた光音響信号をフーリエ変換する第１のフーリエ変換手段と、
    前記光波形微分信号を所定のサンプリングレートでサンプリングした信号をフーリエ変換する第２のフーリエ変換手段と、
    前記フーリエ変換された光波形微分信号の逆数を逆フィルタとして求める逆フィルタ演算手段と、
    前記フーリエ変換された光音響信号に前記逆フィルタを適用するフィルタ適用手段と、
    前記逆フィルタが適用された光音響信号をフーリエ逆変換するフーリエ逆変換手段と、
    を有するものであることを特徴とする請求項１から９何れかに記載の光音響信号処理装置。
11. 前記光音響信号が第１のサンプリングレートでサンプリングされ、前記光波形微分信号が前記第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、
    前記第１のサンプリングレートでサンプリングされた光音響信号を、前記第２のサンプリングレートでリサンプルするリサンプル手段を更に備え、
    前記第１のフーリエ変換手段が、前記リサンプル手段でリサンプルされた光音響信号をフーリエ変換するものであることを特徴とする請求項１０に記載の光音響信号処理装置。
12. 前記光音響信号が第１のサンプリングレートでサンプリングされ、前記光波形微分信号が前記第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、かつ、前記第１のフーリエ変換手段が第１のデータ点数でフーリエ変換を行い、前記第２のフーリエ変換手段が、前記第１のデータ点数よりも多い第２のデータ点数でフーリエ変換を行うものであり、
    前記フーリエ変換された光音響信号に対して、前記第１のデータ点数と前記第２のデータ点数との差の分だけ中央に０を付加するゼロパディングを行うゼロパディング手段を更に備え、
    前記フィルタ適用手段が、前記ゼロパディング手段でゼロパディングが施された信号に対して前記逆フィルタを適用するものであることを特徴とする請求項１０に記載の光音響信号処理装置。
13. 前記光音響信号が第１のサンプリングレートでサンプリングされ、前記光波形微分信号が前記第１のサンプリングレートよりも高い第２のサンプリングレートでサンプリングされた信号に相当し、かつ、前記第１のフーリエ変換手段が第１のデータ点数でフーリエ変換を行い、前記第２のフーリエ変換手段が、前記第１のデータ点数よりも多い第２のデータ点数でフーリエ変換を行うものであり、
    前記フーリエ変換された光波形微分信号から、前記第１のデータ点数と前記第２のデータ点数の差の分だけ高周波成分サンプル点を除去する高周波成分サンプル点除去手段を更に備え、
    前記逆フィルタ演算手段が、前記フーリエ変換された光波形微分信号から高周波成分サンプル点を除去した信号の逆数を逆フィルタとして求めるものであることを特徴とする請求項１０に記載の光音響信号処理装置。
14. 前記光音響信号のサンプリング間隔が、前記被検体に出射された光のパルス時間幅よりも長いことを特徴とする請求項１１から１３何れかに記載の光音響信号処理装置。
15. 複数の検出素子で検出され、前記サンプリング手段でサンプリングされた光音響信号に基づいて前記光音響信号を再構成する光音響信号再構成手段を更に備え、
    前記光波形微分信号逆畳込み手段が、前記光音響信号再構成手段で再構成された光音響信号から前記光波形微分信号をデコンボリューションするものであることを特徴とする請求項１から１４の何れかに記載の光音響信号処理装置。
16. 前記光波形微分信号逆畳込み手段が、前記光音響信号から、前記光波形微分信号に装置応答関数を掛け合わせたものをデコンボリューションすることを特徴とする請求項１から１５何れかに記載の光音響信号処理装置。
17. 被検体に対する出射光が複数の波長の光を含み、前記サンプリング手段が各波長の光に対応した光音響信号をサンプリングし、前記光波形微分信号逆畳込み手段が各波長の光に対応した光音響信号から光波形微分信号をデコンボリューションした信号を生成するものであり、
    各波長の光に対応した、光音響信号から光波形微分信号をデコンボリューションした信号を演算処理する２波長データ演算手段を更に備えたことを特徴とする請求項１から１６何れかに記載の光音響信号処理装置。
18. 前記光音響信号から前記光波形微分信号をデコンボリューションした信号に基づいて光音響画像を生成する光音響画像生成手段を更に備えること特徴とする請求項１から１７何れかに記載の光音響信号処理装置。
19. 前記サンプリング手段が、更に前記被検体に対して送信された音響波に対する反射音響波を更にサンプリングするものであり、
    前記サンプリング反射音響波に基づいて反射音響波画像を生成する反射音響波画像生成手段と、
    前記光音響画像と前記反射音響波画像とを合成する画像合成手段とを更に備えることを特徴とする請求項１８に記載の光音響信号処理装置。
20. 前記画像合成手段が、前記光音響画像と反射音響波画像とを重畳することで画像合成を行うものであることを特徴とする請求項１９に記載の光音響信号処理装置。
21. 被検体へのパルス光の出射後に、該パルス光に起因して被検体内で発生した光音響信号を検出するステップと、
    前記パルス光の光強度の時間波形又は該時間波形を微分した微分波形を表し、パルス光のパルス時間幅を独立変数とする関数を使用して、前記被検体に出射されたパルス光の光強度の時間波形の微分波形である光波形微分信号を生成するステップと、
    前記検出された光音響信号から、前記生成された光波形微分信号をデコンボリューションするステップとを有する光音響信号処理方法。