(58)【調査した分野】(Int.Cl.,DB名)
【発明を実施するための形態】
【0024】
以下、添付図面を参照して、本発明を実施するための形態(以下、「実施の形態」という)を説明する。
【0025】
図1は、本発明の一実施の形態に係る超音波観測装置の構成を示すブロック図である。同図に示す超音波観測装置1は、超音波を用いて観測対象を観測するための装置である。
【0026】
超音波観測装置1は、観測対象へ超音波パルスを出力するとともに、観測対象によって反射された超音波エコーを受信する超音波探触子2と、超音波探触子2との間で電気信号の送受信を行う送受信部3と、超音波エコーを電気信号に変換した電気的なエコー信号に対して所定の演算を施す演算部4と、電気的なエコー信号に対応する画像データの生成を行う画像処理部5と、キーボード、マウス、タッチパネル等のユーザインタフェースを用いて実現され、各種情報の入力を受け付ける入力部6と、液晶または有機EL(Electro Luminescence)等からなる表示パネルを用いて実現され、画像処理部5が生成した画像を含む各種情報を表示する表示部7と、超音波観測に必要な各種情報を記憶する記憶部8と、超音波観測装置1の動作制御を行う制御部9と、を備える。
【0027】
超音波観測装置1は、超音波振動子21が設けられる超音波探触子2と、超音波探触子2が着脱可能に接続され、超音波探触子2以外の上記部分が設けられる処理装置(プロセッサ)とによって構成される。ここで、観測対象が生体組織である場合、超音波探触子2は、生体の体表から超音波を照射する体外式探触子の形態、消化管、胆膵管、血管等の管腔内に挿入する長軸の挿入部を備えたミニチュア超音波プローブの形態、管腔内超音波プローブに光学系をさらに備えた超音波内視鏡の形態、のいずれの形態であってもよい。このうち、超音波内視鏡の形態をとった場合には、管腔内超音波プローブの挿入部の先端側に超音波振動子21が設けられ、管腔内超音波プローブは基端側で処理装置と着脱可能に接続する。
【0028】
超音波振動子21は、送受信部3から受信した電気的なパルス信号を超音波パルス(音響パルス)に変換するとともに、外部の観測対象で反射された超音波エコーを電気的なエコー信号に変換する。超音波探触子2は、超音波振動子21をメカ的に走査させるものであってもよいし、超音波振動子21として複数の素子をアレイ状に設け、送受信にかかわる素子を電子的に切り替えたり、各素子の送受信に遅延をかけたりすることで、電子的に走査させるものであってもよい。本実施の形態では、超音波探触子2として、互いに異なる複数種類のいずれかの超音波探触子2を選択して使用することが可能である。
【0029】
送受信部3は、超音波探触子2と電気的に接続され、電気的なパルス信号を超音波探触子2へ送信するとともに、超音波探触子2から電気的な受信信号であるエコー信号を受信する。具体的には、送受信部3は、予め設定された波形および送信タイミングに基づいて電気的なパルス信号を生成し、この生成したパルス信号を超音波探触子2へ送信する。
【0030】
送受信部3は、エコー信号を増幅する信号増幅部31を有する。信号増幅部31は、受信深度が大きいエコー信号ほど高い増幅率で増幅するSTC(Sensitivity Time Control)補正を行う。
図2は、信号増幅部31が行うSTC補正処理における受信深度と増幅率との関係を示す図である。
図2に示す受信深度zは、超音波の受信開始時点からの経過時間に基づいて算出される量である。
図2に示すように、増幅率β(dB)は、受信深度zが閾値z
thより小さい場合、受信深度zの増加に伴ってβ
0からβ
th(>β
0)へ線型に増加する。また、増幅率β(dB)は、受信深度zが閾値z
th以上である場合、一定値β
thをとる。閾値z
thの値は、観測対象から受信する超音波信号がほとんど減衰してしまい、ノイズが支配的になるような値である。より一般に、増幅率βは、受信深度zが閾値z
thより小さい場合、受信深度zの増加に伴って単調増加すればよい。
【0031】
送受信部3は、信号増幅部31によって増幅されたエコー信号に対してフィルタリング等の処理を施した後、A/D変換することによって時間ドメインのデジタル高周波(RF:Radio Frequency)信号を生成して出力する。なお、超音波探触子2が複数の素子をアレイ状に設けた超音波振動子21を電子的に走査させるものである場合、送受信部3は、複数の素子に対応したビーム合成用の多チャンネル回路を有する。
【0032】
演算部4は、送受信部3が生成したデジタルRF信号に対して受信深度によらず増幅率βを一定とするよう増幅補正を行う増幅補正部41と、増幅補正を行ったデジタルRF信号に高速フーリエ変換(FFT:Fast Fourier Transfom)を施して周波数解析を行うことにより周波数スペクトルを算出する周波数解析部42と、周波数スペクトルの特徴量を算出する特徴量算出部43と、を有する。演算部4は、CPU(Central Proccesing Unit)や各種演算回路等を用いて実現される。
【0033】
図3は、増幅補正部41が行う増幅補正処理における受信深度と増幅率との関係を示す図である。
図3に示すように、増幅補正部41が行う増幅処理における増幅率β(dB)は、受信深度zがゼロのとき最大値β
th−β
0をとり、受信深度zがゼロから閾値z
thに達するまで線型に減少し、受信深度zが閾値z
th以上のときゼロである。このように定められる増幅率によって増幅補正部41がデジタルRF信号を増幅補正することにより、信号増幅部31におけるSTC補正の影響を相殺し、一定の増幅率β
thの信号を出力することができる。なお、増幅補正部41が行う受信深度zと増幅率βの関係は、信号増幅部31における受信深度と増幅率の関係に応じて異なることは勿論である。
【0034】
このような増幅補正を行う理由を説明する。STC補正は、アナログ信号波形の振幅を全周波数帯域にわたって均一に、かつ、深度に対しては単調増加する増幅率で増幅させることで、アナログ信号波形の振幅から減衰の影響を排除する補正処理である。このため、エコー信号の振幅を輝度に変換して表示するBモード画像を生成する場合、かつ、一様な組織を走査した場合には、STC補正を行うことによって深度によらず輝度値が一定になる。すなわち、Bモード画像の輝度値から減衰の影響を排除する効果を得ることができる。
【0035】
一方、本実施の形態のように超音波の周波数スペクトルを算出して解析した結果を利用する場合、STC補正でも超音波の伝播に伴う減衰の影響を正確に排除できるわけではない。なぜなら、一般に減衰量は周波数によって異なるが(後述する式(1)を参照)、STC補正の増幅率は距離だけに応じて変化し、周波数依存性がないためである。
【0036】
上述した問題、すなわち、超音波の周波数スペクトルを算出して解析した結果を利用する場合、STC補正でも超音波の伝播に伴う減衰の影響を正確に排除できるわけではない、という問題を解決するには、Bモード画像を生成する際にSTC補正を施した受信信号を出力する一方、周波数スペクトルに基づいた画像を生成する際に、Bモード画像を生成するための送信とは異なる新たな送信を行い、STC補正を施していない受信信号を出力することが考えられる。ところがこの場合には、受信信号に基づいて生成される画像データのフレームレートが低下してしまうという問題がある。
【0037】
そこで、本実施の形態では、生成される画像データのフレームレートを維持しつつ、Bモード画像用にSTC補正を施した信号に対してSTC補正の影響を排除するために、増幅補正部41によって増幅率の補正を行う。
【0038】
周波数解析部42は、エコー信号に基づくデジタルRF信号を増幅補正した信号の各音線(ラインデータ)を、所定の時間間隔でサンプリングした振幅データ群を高速フーリエ変換することによって音線上の複数の箇所(データ位置)における周波数スペクトルを算出する。
【0039】
図4は、超音波信号の1つの音線におけるデータ配列を模式的に示す図である。同図に示す音線SR
kにおいて、白または黒の長方形は、1つのデータを意味している。音線SR
kは、送受信部3が行うA/D変換におけるサンプリング周波数(例えば50MHz)に対応した時間間隔で離散化されている。
図4では、番号kの音線SR
kの1番目のデータ位置を受信深度zの方向の初期値Z
(k)0として設定した場合を示しているが、初期値の位置は任意に設定することができる。周波数解析部42による算出結果は複素数で得られ、記憶部8に格納される。
【0040】
図4に示すデータ群F
j(j=1、2、・・・、K)は、高速フーリエ変換の対象となる振幅データ群である。一般に、高速フーリエ変換を行うためには、振幅データ群が2のべき乗のデータ数を有している必要がある。この意味で、振幅データ群F
j(j=2、・・・、K−1)はデータ数が16(=2
4)で正常なデータ群である一方、振幅データ群F
1、F
Kは、それぞれデータ数が9、12であるため異常なデータ群である。異常なデータ群に対して高速フーリエ変換を行う際には、不足分だけゼロデータを挿入することにより、正常な振幅データ群を生成する処理を行う。この点については、周波数解析部42の処理を説明する際に詳述する(
図9を参照)。
【0041】
図5は、周波数解析部42が算出する周波数スペクトルの例を示す図である。ここでいう「周波数スペクトル」とは、振幅データ群を高速フーリエ変換(FFT演算)することによって得られた「ある受信深度zにおける強度の周波数分布」を意味する。また、ここでいう「強度」とは、例えばエコー信号の電圧、エコー信号の電力、超音波エコーの音圧、超音波エコーの音響エネルギー等のパラメータ、これらパラメータの振幅や時間積分値やその組み合わせのいずれかを指す。
【0042】
図5では、横軸に周波数fを取っている。また、
図5では、縦軸に、強度I
0を基準強度I
c(定数)で除した量の常用対数(デシベル表現)I=10log
10(I
0/I
c)を取っている。
図5において、受信深度zは一定である。
図5に示す直線L
10については後述する。なお、本実施の形態において、曲線および直線は、離散的な点の集合からなる。
【0043】
図5に示す周波数スペクトルC
1において、以後の演算に使用する周波数帯域の下限周波数f
Lおよび上限周波数f
Hは、超音波振動子21の周波数帯域、送受信部3が送信するパルス信号の周波数帯域などをもとに決定されるパラメータであり、例えばf
L=3MHz、f
H=10MHzである。以下、
図5において、下限周波数f
Lおよび上限周波数f
Hによって定まる周波数帯域を「周波数帯域F」という。
【0044】
一般に、周波数スペクトルは、観測対象が生体組織である場合、超音波が走査された生体組織の性状(属性)によって異なる傾向を示す。これは、周波数スペクトルが、超音波を散乱する散乱体の大きさ、数密度、音響インピーダンス等と相関を有しているためである。ここでいう「生体組織の性状」とは、例えば悪性腫瘍(癌)、良性腫瘍、内分泌腫瘍、粘液性腫瘍、正常組織、脈管などのことである。
【0045】
特徴量算出部43は、複数の周波数スペクトルの特徴量をそれぞれ算出し、超音波が観測対象を伝播する際の互いに異なる減衰特性を与える複数の減衰率候補値の各々において、各周波数スペクトルの特徴量(以下、補正前特徴量という)に対して超音波の減衰の影響を排除する減衰補正を行うことによって各周波数スペクトルの補正特徴量を算出し、該補正特徴量を用いて複数の減衰率候補値の中から観測対象に最適な減衰率を設定する。特徴量算出部43は、周波数スペクトルを直線で近似することによって周波数スペクトルの補正前特徴量を算出する近似部431と、近似部431が算出した補正前特徴量に対して複数の減衰率候補値の各々に基づいた減衰補正を行うことによって補正特徴量を算出する減衰補正部432と、減衰補正部432がすべての周波数スペクトルに対して算出した補正特徴量の統計的なばらつきに基づいて複数の減衰率候補値の中から最適な減衰率を設定する最適減衰率設定部433と、を有する。
【0046】
近似部431は、所定周波数帯域における周波数スペクトルの回帰分析を行って周波数スペクトルを一次式(回帰直線)で近似することにより、この近似した一次式を特徴付ける補正前特徴量を算出する。例えば、
図5に示す周波数スペクトルC
1の場合、近似部431は、周波数帯域Fで回帰分析を行い周波数スペクトルC
1を一次式で近似することによって回帰直線L
10を得る。換言すると、近似部431は、回帰直線L
10の傾きa
0、切片b
0、および周波数帯域Fの中心周波数f
M=(f
L+f
H)/2の回帰直線上の値であるミッドバンドフィット(Mid-band fit)c
0=a
0f
M+b
0を補正前特徴量として算出する。
【0047】
3つの補正前特徴量のうち、傾きa
0は、超音波の散乱体の大きさと相関を有し、一般に散乱体が大きいほど傾きが小さな値を有すると考えられる。また、切片b
0は、散乱体の大きさ、音響インピーダンスの差、散乱体の数密度(濃度)等と相関を有している。具体的には、切片b
0は、散乱体が大きいほど大きな値を有し、音響インピーダンスの差が大きいほど大きな値を有し、散乱体の数密度が大きいほど大きな値を有すると考えられる。ミッドバンドフィットc
0は、傾きa
0と切片b
0から導出される間接的なパラメータであり、有効な周波数帯域内の中心におけるスペクトルの強度を与える。このため、ミッドバンドフィットc
0は、散乱体の大きさ、音響インピーダンスの差、散乱体の数密度に加えて、Bモード画像の輝度とある程度の相関を有していると考えられる。なお、特徴量算出部43は、回帰分析によって二次以上の多項式で周波数スペクトルを近似するようにしてもよい。
【0048】
減衰補正部432が行う補正について説明する。一般に、超音波の減衰量A(f,z)は、超音波が受信深度0と受信深度zとの間を往復する間に生じる減衰であり、往復する前後の強度変化(デシベル表現での差)として定義される。減衰量A(f,z)は、一様な組織内では周波数に比例することが経験的に知られており、以下の式(1)で表現される。
A(f,z)=2αzf ・・・(1)
ここで、比例定数αは減衰率と呼ばれる量である。また、zは超音波の受信深度であり、fは周波数である。減衰率αの具体的な値は、観測対象が生体である場合、生体の部位に応じて定まる。減衰率αの単位は、例えばdB/cm/MHzである。本実施の形態において、減衰補正部432は、観測対象に最も適合する減衰率(最適な減衰率)を設定するために、複数の減衰率候補値に対してそれぞれ減衰補正を行う。複数の減衰率候補値の詳細については、
図8および
図10を参照して後述する。
【0049】
減衰補正部432は、近似部431が抽出した補正前特徴量(傾きa
0、切片b
0、ミッドバンドフィットc
0)に対し、以下に示す式(2)〜(4)にしたがって減衰補正を行うことにより、補正特徴量a、b、cを算出する。
a=a
0+2αz ・・・(2)
b=b
0 ・・・(3)
c=c
0+A(f
M,z)=c
0+2αzf
M(=af
M+b) ・・・(4)
式(2)、(4)からも明らかなように、減衰補正部432は、超音波の受信深度zが大きいほど、補正量が大きい補正を行う。また、式(3)によれば、切片に関する補正は恒等変換である。これは、切片が周波数0(Hz)に対応する周波数成分であって減衰の影響を受けないためである。
【0050】
図6は、減衰補正部432が補正した補正特徴量a、b、cをパラメータとして有する直線を示す図である。直線L
1の式は、
I=af+b=(a
0+2αz)f+b
0 ・・・(5)
で表される。この式(5)からも明らかなように、直線L
1は、減衰補正前の直線L
10と比較して、傾きが大きく(a>a
0)、かつ切片が同じ(b=b
0)である。
【0051】
最適減衰率設定部433は、減衰補正部432がすべての周波数スペクトルに対して減衰率候補値ごとに算出した補正特徴量の統計的なばらつきが最小である減衰率候補値を最適な減衰率として設定する。本実施の形態では、統計的なばらつきを示す量として分散を適用する。この場合、最適減衰率設定部433は、分散が最小となる減衰率候補値を最適な減衰率として設定する。上述した3つの補正特徴量a、b、cのうち独立なのは2つである。加えて、補正特徴量bは減衰率に依存しない。したがって、補正特徴量a、cに対して最適な減衰率を設定する場合、最適減衰率設定部433は、補正特徴量aおよびcのいずれか一方の分散を算出すればよい。
【0052】
ただし、最適減衰率設定部433が最適な減衰率を設定する際に用いる補正特徴量は、特徴量画像データ生成部52が特徴量画像データを生成する際に用いる補正特徴量と同じ種類であることが好ましい。すなわち、特徴量画像データ生成部52が補正特徴量として傾きを用いて特徴量画像データを生成する場合は補正特徴量aの分散を適用し、特徴量画像データ生成部52が補正特徴量としてミッドバンドフィットを用いて特徴量画像データを生成する場合は補正特徴量cの分散を適用するのがより好ましい。これは、減衰量A(f,z)を与える式(1)があくまで理想的なものに過ぎず、現実には以下の式(6)の方が適切であることによる。
A(f,z)=2αzf+2α
1z ・・・(6)
式(6)の右辺第2項のα
1は、超音波の受信深度zに比例して信号強度が変化する大きさを表す係数であり、観測対象の組織が不均一であることや、ビーム合成時のチャンネル数の変更などに起因して発生する信号強度の変化を表す係数である。式(6)の右辺第2項が存在するため、補正特徴量cを用いて特徴量画像を生成する場合は、補正特徴量cの分散を適用した方が正確に減衰を補正することができる(式(4)を参照)。一方、周波数fに比例する係数である補正特徴量aを用いて特徴量画像を生成する場合は、補正特徴量aの分散を適用した方が、右辺第2項の影響を排除して正確に減衰を補正することができる。例えば、減衰率αの単位がdB/cm/MHzである場合、係数α
1の単位はdB/cmである。
【0053】
ここで、統計的なばらつきに基づいて最適な減衰率を設定することができる理由を説明する。観測対象に最適な減衰率を適用した場合、観測対象と超音波振動子21との距離に関わらず、特徴量は観測対象に固有の値へ収束し、統計的なばらつきが小さくなると考えられる。その一方で、観測対象に適合しない減衰率候補値を最適な減衰率とした場合、減衰補正が過剰であるかまたは不足するため、超音波振動子21との距離に応じて特徴量にずれが生じ、特徴量の統計的なばらつきが大きくなると考えられる。したがって、統計的なばらつきが最も小さい減衰率候補値が、観察対象にとって最適な減衰率であるということができる。
【0054】
図7は、同じ観測対象に対して2つの異なる減衰率候補値に基づいてそれぞれ減衰補正された補正特徴量の分布例を模式的に示す図である。
図7では、横軸を補正特徴量とし、縦軸を頻度としている。
図7に示す2つの分布曲線N
1、N
2は、頻度の総和が同じである。
図7に示す場合、分布曲線N
1は、分布曲線N
2と比較して特徴量の統計的なばらつきが小さく(分散が小さく)、山が急峻な形状をなす。したがって、最適減衰率設定部433は、この2つの分布曲線N
1、N
2に対応する2つの減衰率候補値から最適な減衰率を設定する場合、分布曲線N
1に対応する減衰率候補値を最適な減衰率として設定する。
【0055】
画像処理部5は、エコー信号の振幅を輝度に変換して表示する超音波画像であるBモード画像データを生成するBモード画像データ生成部51と、最適減衰率設定部433が設定した最適な減衰率に基づく特徴量を視覚情報と関連づけてBモード画像とともに表示する特徴量画像データを生成する特徴量画像データ生成部52と、を有する。
【0056】
Bモード画像データ生成部51は、デジタル信号に対してバンドパスフィルタ、対数変換、ゲイン処理、コントラスト処理等の公知の技術を用いた信号処理を行うとともに、表示部7における画像の表示レンジに応じて定まるデータステップ幅に応じたデータの間引き等を行うことによってBモード画像データを生成する。Bモード画像は、色空間としてRGB表色系を採用した場合の変数であるR(赤)、G(緑)、B(青)の値を一致させたグレースケール画像である。
【0057】
特徴量画像データ生成部52は、特徴量算出部43が算出した特徴量に関連する視覚情報をBモード画像データにおける画像の各画素に対して重畳することによって特徴量画像データを生成する。特徴量画像データ生成部52は、例えば
図4に示す1つの振幅データ群F
j(j=1、2、・・・、K)のデータ量に対応する画素領域に対し、その振幅データ群F
jから算出される周波数スペクトルの特徴量に対応する視覚情報を割り当てる。特徴量画像データ生成部52は、例えば上述した傾き、切片、ミッドバンドフィットのいずれか一つに視覚情報としての色相を対応付けることによって特徴量画像を生成する。なお、特徴量画像データ生成部52が、傾き、切片、ミッドバンドフィットから選択される2つの特徴量の一方に色相を対応付けるとともに、他方に明暗を対応付けることによって特徴量画像データを生成するようにしてもよい。特徴量に関連する視覚情報としては、例えば色相、彩度、明度、輝度値、R(赤)、G(緑)、B(青)などの所定の表色系を構成する色空間の変数を挙げることができる。
【0058】
記憶部8は、減衰補正部432が減衰率候補値に応じて周波数スペクトルごとに算出した複数の特徴量、および該複数の特徴量の統計的なばらつきを与える分散を減衰率候補値と対応づけて記憶する特徴量情報記憶部81を有する。
【0059】
記憶部8は、上記以外にも、例えば増幅処理に必要な情報(
図2に示す増幅率と受信深度との関係)、増幅補正処理に必要な情報(
図3に示す増幅率と受信深度との関係)、減衰補正処理に必要な情報(式(1)参照)、周波数解析処理に必要な窓関数(Hamming、Hanning、Blackman等)の情報等を記憶する。
【0060】
また、記憶部8は、超音波観測装置1の作動方法を実行するための作動プログラムを含む各種プログラムを記憶する。作動プログラムは、ハードディスク、フラッシュメモリ、CD−ROM、DVD−ROM、フレキシブルディスク等のコンピュータ読み取り可能な記録媒体に記録して広く流通させることも可能である。なお、上述した各種プログラムは、通信ネットワークを介してダウンロードすることによって取得することも可能である。ここでいう通信ネットワークは、例えば既存の公衆回線網、LAN(Local Area Network)、WAN(Wide Area Network)などによって実現されるものであり、有線、無線を問わない。
【0061】
以上の構成を有する記憶部8は、各種プログラム等が予めインストールされたROM(Read Only Memory)、および各処理の演算パラメータやデータ等を記憶するRAM(Random Access Memory)等を用いて実現される。
【0062】
制御部9は、演算および制御機能を有するCPU(Central Proccesing Unit)や各種演算回路等を用いて実現される。制御部9は、記憶部8が記憶、格納する情報を記憶部8から読み出し、超音波観測装置1の作動方法に関連した各種演算処理を実行することによって超音波観測装置1を統括して制御する。なお、制御部9と演算部4を、共通のCPU等を用いて構成することも可能である。
【0063】
図8は、以上の構成を有する超音波観測装置1が行う処理の概要を示すフローチャートである。超音波観測装置1は、まず超音波探触子2によって新規の観測対象の測定を行う(ステップS1)。具体的には、超音波探触子2の超音波振動子21は、電気的なパルス信号を超音波パルスへ変換し、観測対象へ順次送信する。超音波パルスは観測対象によってそれぞれ反射され、超音波エコーが生じる。超音波振動子21は、超音波エコーを電気的なエコー信号に変換する。この際、パルス信号の周波数帯域は、超音波振動子21におけるパルス信号の超音波パルスへの電気音響変換の線形応答周波数帯域をほぼカバーする広帯域にするとよい。それにより、後述する周波数スペクトルの近似処理において、精度のよい近似を行うことが可能となる。
【0064】
超音波探触子2からエコー信号を受信した信号増幅部31は、そのエコー信号の増幅を行う(ステップS2)。ここで、信号増幅部31は、例えば
図2に示す増幅率と受信深度との関係に基づいてエコー信号の増幅(STC補正)を行う。この際、信号増幅部31におけるエコー信号の各種処理周波数帯域は、超音波振動子21による超音波エコーのエコー信号への音響電気変換の線型応答周波数帯域をほぼカバーする広帯域にするとよい。これも、後述する周波数スペクトルの近似処理において精度のよい近似を行うことを可能とするためである。
【0065】
続いて、Bモード画像データ生成部51は、信号増幅部31が増幅したエコー信号を用いてBモード画像データを生成する(ステップS3)。その後、制御部9は、生成されたBモード画像データに対応するBモード画像を表示部7に表示させる(ステップS4)。
【0066】
増幅補正部41は、送受信部3から出力された信号に対して受信深度によらず増幅率が一定となる増幅補正を行う(ステップS5)。ここで、増幅補正部41は、例えば
図3に示す増幅率と受信深度との関係に基づいて増幅補正を行う。
【0067】
この後、周波数解析部42は、FFT演算による周波数解析を行うことによって全ての振幅データ群に対する周波数スペクトルを算出する(ステップS6)。
図9は、ステップS6において周波数解析部42が実行する処理の概要を示すフローチャートである。以下、
図9に示すフローチャートを参照して、周波数解析処理を詳細に説明する。
【0068】
まず、周波数解析部42は、解析対象の音線を識別するカウンタkをk
0とする(ステップS21)。
【0069】
続いて、周波数解析部42は、FFT演算用に取得する一連のデータ群(振幅データ群)を代表するデータ位置(受信深度に相当)Z
(k)の初期値Z
(k)0を設定する(ステップS22)。例えば、
図4では、上述したように、音線SR
kの1番目のデータ位置を初期値Z
(k)0として設定した場合を示している。
【0070】
その後、周波数解析部42は、データ位置Z
(k)が属する振幅データ群を取得し(ステップS23)、取得した振幅データ群に対し、記憶部8が記憶する窓関数を作用させる(ステップS24)。このように振幅データ群に対して窓関数を作用させることにより、振幅データ群が境界で不連続になることを回避し、アーチファクトが発生するのを防止することができる。
【0071】
続いて、周波数解析部42は、データ位置Z
(k)の振幅データ群が正常なデータ群であるか否かを判定する(ステップS25)。
図4を参照した際に説明したように、振幅データ群は、2のべき乗のデータ数を有している必要がある。以下、正常な振幅データ群のデータ数を2
n(nは正の整数)とする。本実施の形態では、データ位置Z
(k)が、できるだけZ
(k)が属する振幅データ群の中心になるよう設定される。具体的には、振幅データ群のデータ数は2
nであるので、Z
(k)はその振幅データ群の中心に近い2
n/2(=2
n-1)番目の位置に設定される。この場合、振幅データ群が正常であるとは、データ位置Z
(k)の前方に2
n-1−1(=Nとする)個のデータがあり、データ位置Z
(k)の後方に2
n-1(=Mとする)個のデータがあることを意味する。
図4に示す場合、振幅データ群F
2、F
3はともに正常である。なお、
図4ではn=4(N=7,M=8)の場合を例示している。
【0072】
ステップS25における判定の結果、データ位置Z
(k)の振幅データ群が正常である場合(ステップS25:Yes)、周波数解析部42は、後述するステップS27へ移行する。
【0073】
ステップS25における判定の結果、データ位置Z
(k)の振幅データ群が正常でない場合(ステップS25:No)、周波数解析部42は、不足分だけゼロデータを挿入することによって正常な振幅データ群を生成する(ステップS26)。ステップS25において正常でないと判定された振幅データ群(例えば
図4の振幅データ群F
1、F
K)は、ゼロデータを追加する前に窓関数が作用されている。このため、振幅データ群にゼロデータを挿入してもデータの不連続は生じない。ステップS26の後、周波数解析部42は、後述するステップS27へ移行する。
【0074】
ステップS27において、周波数解析部42は、振幅データ群を用いてFFT演算を行うことにより、振幅の周波数分布である周波数スペクトルを得る(ステップS27)。
図5に示す周波数スペクトルC
1は、ステップS27の結果として得られる周波数スペクトルの一例である。
【0075】
続いて、周波数解析部42は、データ位置Z
(k)をステップ幅Dで変化させる(ステップS28)。ステップ幅Dは、記憶部8が予め記憶しているものとする。
図4では、D=15の場合を例示している。ステップ幅Dは、Bモード画像データ生成部51がBモード画像データを生成する際に利用するデータステップ幅と一致させることが望ましいが、周波数解析部42における演算量を削減したい場合には、ステップ幅Dとしてデータステップ幅より大きい値を設定してもよい。
【0076】
その後、周波数解析部42は、データ位置Z
(k)が音線SR
kにおける最大値Z
(k)maxより大きいか否かを判定する(ステップS29)。データ位置Z
(k)が最大値Z
(k)maxより大きい場合(ステップS29:Yes)、周波数解析部42はカウンタkを1増加させる(ステップS30)。これは、処理をとなりの音線へ移すことを意味する。一方、データ位置Z
(k)が最大値Z
(k)max以下である場合(ステップS29:No)、周波数解析部42はステップS23へ戻る。このようにして、周波数解析部42は、音線SR
kに対して、[(Z
(k)max−Z
(k)0+1)/D+1]個の振幅データ群に対するFFT演算を行う。ここで、[X]は、Xを超えない最大の整数を表す。
【0077】
ステップS30の後、周波数解析部42は、カウンタkが最大値k
maxより大きいか否かを判定する(ステップS31)。カウンタkがk
maxより大きい場合(ステップS31:Yes)、周波数解析部42は一連のFFT処理を終了する。一方、カウンタkがk
max以下である場合(ステップS31:No)、周波数解析部42はステップS22に戻る。
【0078】
このようにして、周波数解析部42は、解析対象領域内の(k
max−k
0+1)本の音線の各々について複数回のFFT演算を行う。
【0079】
なお、以上の説明では、周波数解析部42が超音波信号を受信したすべての領域に対して周波数解析処理を行うものとしたが、入力部6が特定の深度幅および音線幅で区切られる関心領域の設定入力を受け付け可能な構成とし、設定された関心領域内においてのみ周波数解析処理を行うようにすることも可能である。
【0080】
以上説明したステップS6の周波数解析処理に続いて、特徴量算出部43は、複数の周波数スペクトルの補正前特徴量をそれぞれ算出し、超音波が観測対象を伝播する際の互いに異なる減衰特性を与える複数の減衰率候補値の各々において、各周波数スペクトルの補正前特徴量に対して超音波の減衰の影響を排除する減衰補正を行うことによって各周波数スペクトルの補正特徴量を算出し、該補正特徴量を用いて複数の減衰率候補値の中から観測対象に最適な減衰率を設定する(ステップS7〜S13)。以下、ステップS7〜S13の処理を詳細に説明する。
【0081】
ステップS7において、近似部431は、周波数解析部42が算出した複数の周波数スペクトルをそれぞれ回帰分析することにより、各周波数スペクトルに対応する補正前特徴量を算出する(ステップS7)。具体的には、近似部431は、各周波数スペクトルを回帰分析することによって一次式で近似し、補正前特徴量として傾きa
0、切片b
0、ミッドバンドフィットc
0を算出する。例えば、
図5に示す直線L
10は、近似部431が周波数帯域Fの周波数スペクトルC
1に対し回帰分析によって近似した回帰直線である。
【0082】
この後、最適減衰率設定部433は、後述する減衰補正を行う際に適用する減衰率候補値αの値を所定の初期値α
0に設定する(ステップS8)。この初期値α
0の値は、予め記憶部8が記憶しておき、最適減衰率設定部433が記憶部8を参照するようにすればよい。
【0083】
続いて、減衰補正部432は、近似部431が各周波数スペクトルに対して近似した補正前特徴量に対し、減衰率候補値をαとして減衰補正を行うことにより、補正特徴量を算出し、減衰率候補値αとともに特徴量情報記憶部81に格納する(ステップS9)。
図6に示す直線L
1は、減衰補正部432が減衰補正処理を行うことによって得られる直線の例である。
【0084】
ステップS9において、減衰補正部432は、上述した式(2)、(4)における受信深度zに、超音波信号の音線のデータ配列を用いて得られるデータ位置Z=(f
sp/2v
s)Dnを代入することによって算出する。ここで、f
spはデータのサンプリング周波数、v
sは音速、Dはデータステップ幅、nは処理対象の振幅データ群のデータ位置までの音線の1番目のデータからのデータステップ数である。例えば、データのサンプリング周波数f
spを50MHzとし、音速v
sを1530m/secとし、
図4に示すデータ配列を採用してステップ幅Dを15とすると、z=0.2295n(mm)となる。
【0085】
最適減衰率設定部433は、減衰補正部432が各周波数スペクトルに対して減衰補正することによって得られた複数の補正特徴量のうち代表となる補正特徴量の分散を算出し、減衰率候補値αと対応づけて特徴量情報記憶部81へ格納する(ステップS10)。補正特徴量が傾きa、ミッドバンドフィットcである場合、上述したように、最適減衰率設定部433は、補正特徴量aおよびcのいずれか一方の分散を算出する。このステップS10において、特徴量画像データ生成部52が、傾きを用いて特徴量画像データを生成する場合は補正特徴量aの分散を適用し、ミッドバンドフィットを用いて特徴量画像データを生成する場合は補正特徴量cの分散を適用するのが好ましい。
【0086】
この後、最適減衰率設定部433は、減衰率候補値αの値をΔαだけ増加させ(ステップS11)、増加後の減衰率候補値αと所定の最大値α
maxとの大小を比較する(ステップS12)。ステップS12における比較の結果、減衰率候補値αが最大値α
maxより大きい場合(ステップS12:Yes)、超音波観測装置1はステップS13へ移行する。一方、ステップS12における比較の結果、減衰率候補値αが最大値α
max以下である場合(ステップS12:No)、超音波観測装置1はステップS9へ戻る。
【0087】
ステップS13において、最適減衰率設定部433は、特徴量情報記憶部81が記憶する減衰率候補値ごとの分散を参照し、分散が最小である減衰率候補値を最適な減衰率として設定する(ステップS13)。
【0088】
図10は、最適減衰率設定部433が行う処理の概要を示す図である。α
0=0(dB/cm/MHz)、α
max=1.0(dB/cm/MHz)、Δα=0.2(dB/cm/MHz)とした場合の減衰率候補値αと分散S(α)との関係の例を示す図である。
図10に示す場合、減衰率候補値αが0.2(dB/cm/MHz)のときに分散が最小値S(α)
minをとる。したがって、
図10に示す場合、最適減衰率設定部433は、α=0.2(dB/cm/MHz)を最適な減衰率として設定する。
【0089】
特徴量画像データ生成部52は、Bモード画像データ生成部51が生成したBモード画像データにおける各画素に対し、ステップS13で設定された最適な減衰率に基づく補正特徴量に関連づけた視覚情報(例えば色相)を重畳するとともに、最適な減衰率の情報を加えることによって特徴量画像データを生成する(ステップS14)。
【0090】
この後、表示部7は、制御部9の制御のもと、特徴量画像データ生成部52が生成した特徴量画像データに対応する特徴量画像を表示する(ステップS15)。
図11は、表示部7における特徴量画像の表示例を模式的に示す図である。同図に示す特徴量画像101は、Bモード画像に特徴量に関する視覚情報が重畳された画像を表示する重畳画像表示部102と、観測対象の識別情報および最適な減衰率として設定された減衰率候補値の情報を表示する情報表示部103とを有する。なお、情報表示部103に、特徴量の情報、近似式の情報、ゲインやコントラスト等の画像情報等をさらに表示するようにしてもよい。また、特徴量画像に対応するBモード画像を特徴量画像と並べて表示してもよい。また、減衰率候補値の情報を表示するか否かの指示信号を入力部6が受け付け可能な構成としてもよい。
【0091】
以上説明してきた一連の処理(ステップS1〜S15)において、ステップS4の処理とステップS5〜S13の処理とを並行して行うようにしてもよい。
【0092】
以上説明した本発明の一実施の形態によれば、超音波が観測対象を伝播する際の互いに異なる減衰特性を与える複数の減衰率候補値の中から観測対象に最適な減衰率を設定し、該最適な減衰率を用いて減衰補正を行うことによって複数の周波数スペクトルの各々の特徴量を算出するため、観測対象に適合した超音波の減衰特性を簡易な計算によって求めることができるとともに、その減衰特性を利用した観測を行うことができる。
【0093】
また、本実施の形態によれば、各周波数スペクトルを減衰補正した補正特徴量の統計的なばらつきに基づいて最適な減衰率を設定するため、複数の減衰モデルとフィッティングを行う従来技術と比較して、計算量を削減することができる。
【0094】
また、本実施の形態によれば、観測対象に適合する減衰率が未知の場合であっても、最適な減衰率を設定することができる。
【0095】
(実施の形態の変形例1)
図12は、本実施の形態の変形例1に係る超音波観測装置の最適減衰率設定部が行う処理の概要を示す図である。
図12では、α
0=0(dB/cm/MHz)、α
max=1.0(dB/cm/MHz)、Δα=0.2(dB/cm/MHz)とした場合の減衰率候補値αと分散S(α)との関係の例を示しており、減衰率候補値α=0、0.2、0.4、0.6、0.8、1.0(いずれもdB/cm/MHz)における分散S(α)の値は、
図10とそれぞれ同じである。本変形例1では、最適減衰率設定部433が最適な減衰率を設定する前に、近似部431が回帰分析を行うことによって減衰率候補値αにおける分散S(α)の値を補間する曲線Rを算出する。その後、最適減衰率設定部433は、この曲線Rに対し、0(dB/cm/MHz)≦α≦1.0(dB/cm/MHz)における最小値S(α)’
minを算出し、そのときの減衰率候補値の値α’を最適な減衰率として設定する。
図12に示す場合、最適な減衰率α’は、0(dB/cm/MHz)と0.2(dB/cm/MHz)の間の値となる。
【0096】
(実施の形態の変形例2)
次に、本発明の実施の形態の変形例2について説明する。本変形例では、最適減衰率設定部433が、特徴量画像として表示する際のダイナミックレンジよりも広いダイナミックレンジで最適な減衰率を設定する。
【0097】
具体的には、特徴量画像データ生成部52が生成する画像の表示ダイナミックレンジを70dBとすると、特徴量算出部43は、このダイナミックレンジ(70dB)よりも大きなダイナミックレンジ(例えば100dB)で減衰演算処理を行う。例えば、特徴量画像データ生成部52が、8bitの固定小数点方式を用いるのに対し、特徴量算出部43は、32bitの浮動小数点方式を用いて特徴量の算出から最適な減衰率の設定までを含む減衰演算処理を行う。
【0098】
本変形例2によれば、固定小数点方式を用いた減衰演算処理と比して、演算精度を向上させることができる。補正前特徴量の演算から分散に基づく二次曲線の生成を一層高精度に行うことで、最適な減衰率を高精度に算出することができる。
【0099】
ここまで、本発明を実施するための形態を説明してきたが、本発明は上述した実施の形態によってのみ限定されるべきものではない。例えば、最適減衰率設定部433は、超音波画像の全てのフレームで最適な減衰率に相当する最適減衰率相当値をそれぞれ算出し、最新のフレームにおける最適減衰率相当値を含む所定数の最適減衰率相当値の平均値、中央値または最頻値を最適な減衰率として設定してもよい。この場合には、各フレームで最適な減衰率を設定する場合と比較して、最適な減衰率の変化が少なくなってその値を安定させることができる。
【0100】
また、最適減衰率設定部433は、超音波画像の所定のフレーム間隔で最適な減衰率を設定するようにしてもよい。これにより、計算量を大幅に削減することができる。この場合には、次に最適な減衰率を設定するまでの間、最後に設定した最適な減衰率の値を使用すればよい。
【0101】
また、統計的なばらつきを算出する対象領域を音線ごととしてもよいし、受信深度が所定値以上の領域としてもよい。これらの領域の設定を入力部6が受け付け可能な構成としてもよい。
【0102】
また、最適減衰率設定部433が、設定された関心領域内とその関心領域外とで個別に最適な減衰率を設定するようにしてもよい。
【0103】
また、入力部6が減衰率候補値の初期値α
0の設定変更の入力を受け付け可能な構成としてもよい。
【0104】
また、統計的なばらつきを与える量として、例えば標準偏差、母集団における特徴量の最大値と最小値の差、特徴量の分布の半値幅のいずれかを適用することも可能である。なお、統計的なばらつきを与える量として分散の逆数を適用する場合も考えられるが、この場合には、その値が最大となる減衰率候補値が最適な減衰率となることはいうまでもない。
【0105】
また、最適減衰率設定部433が、複数種類の補正特徴量の統計的なばらつきをそれぞれ算出し、統計的なばらつきが最小である場合の減衰率候補値を最適な減衰率として設定することも可能である。
【0106】
また、減衰補正部432が複数の減衰率候補値を用いて周波数スペクトルを減衰補正した後、近似部431が減衰補正後の各周波数スペクトルに対して回帰分析を行うことによって補正特徴量を算出するようにしてもよい。
【0107】
このように、本発明は、特許請求の範囲に記載した技術的思想を逸脱しない範囲内において、様々な実施の形態を含みうるものである。