超聲多普勒血流成像的處理方法及處理系統(tǒng)與流程
本發(fā)明屬于醫(yī)療超聲技術(shù)領(lǐng)域，主要涉及一種超聲多普勒血流成像的處理方法及處理系統(tǒng)。
背景技術(shù)：
彩色超聲診斷儀(b超機(jī))的彩色血流成像，以其獨有的實時動態(tài)特性，成為現(xiàn)代醫(yī)學(xué)不可或缺的輔助診斷的手段之一，在臨床診斷中成為某些病癥的判斷標(biāo)準(zhǔn)。
當(dāng)前彩色血流成像普遍采用自相關(guān)技術(shù)；根據(jù)多普勒效應(yīng)，散射子運動導(dǎo)致的頻移大小正比于發(fā)射信號頻率和散射子運動速度的乘積。
傳統(tǒng)的一維自相關(guān)技術(shù)，通過在回波信號慢時間方向上的自相關(guān)來估算散射子在發(fā)射信號頻率上導(dǎo)致的相位差，從而估算出散射子的運動速度；該估計算法是建立在窄帶信號的模型基礎(chǔ)上，也就要求發(fā)射波形較長來保證符合窄帶模型，其速度估計性能，尤其是對快速運動的散射子的速度估計能力，隨發(fā)射信號的帶寬的增加而顯著惡化；同時，當(dāng)發(fā)射波形較長時，血流散射子的成像分辨率會下降。
另外，由于超聲在組織內(nèi)傳播時，存在的非線性頻率衰減，隨著聲速在組織內(nèi)的傳播距離增加導(dǎo)致的頻率成分不同衰減，傳統(tǒng)一維自相關(guān)血流成像技術(shù)因為是基于固定的發(fā)射信號中心頻率和固定帶寬而進(jìn)行計算，所以其估算速度的性能也隨傳播距離的增加而下降。
如上所述，對于傳統(tǒng)的一維自相關(guān)血流成像技術(shù)，由于技術(shù)理論基礎(chǔ)是基于窄帶信號模型而得出，這就需要發(fā)射信號脈寬較長，從而影響了圖像的分辨率。另外，超聲在組織內(nèi)存在非線性頻率衰減現(xiàn)象，隨著超聲傳播距離的增加，其信號的頻帶帶寬和中心頻率都會因此而偏移，傳統(tǒng)的自相關(guān)技術(shù)，只能在基于固定的發(fā)射信號中心頻率進(jìn)行速度估算，所以隨著深度增加時，其估算的速度誤差也會增大。
相應(yīng)的，當(dāng)前公開的互相關(guān)血流成像技術(shù)和蝶型搜索成像技術(shù)基于寬帶信號模型，解決了上述一維自相關(guān)技術(shù)與分辨率的矛盾，但也存在計算復(fù)雜度高，不易實時實現(xiàn)，和不能有效解決頻率衰減導(dǎo)致的速度偏差問題。
技術(shù)實現(xiàn)要素：
本發(fā)明的目的在于提供一種超聲多普勒血流成像的處理方法及處理系統(tǒng)。
為了實現(xiàn)上述發(fā)明目的之一，本發(fā)明一實施方式的超聲多普勒血流成像的處理方法，所述方法包括以下步驟：
s1、分別獲取每個掃查包下、每條掃查線上、每個采樣點對應(yīng)的回波rf信號；分別對獲取的回波rf信號做希爾伯特變換以形成正交iq復(fù)信號；
s2、在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，以獲得對應(yīng)的濾波變換信號；
s3、在快時間方向，分別對獲得的濾波變換信號做短時傅里葉變換，使其結(jié)果形成頻域矩陣；
s4、獲取當(dāng)前所述頻域矩陣中，每個快時間方向上的頻率分量；
s5、根據(jù)一維自相關(guān)算法，估計當(dāng)前頻域矩陣中每列數(shù)據(jù)在其對應(yīng)頻率分量上的平均速度和能量；
s6、對當(dāng)前頻域矩陣中對應(yīng)頻率分量的平均速度和能量進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述步驟s3具體包括：
p1、在同一深度下，對應(yīng)每個掃查包、在每條掃查線的快時間方向上分別獲取l個采樣點對應(yīng)的濾波變換信號，l小于等于n，n表示每條掃查線對應(yīng)的原始采樣點數(shù)量；
p2、對每條掃查線對應(yīng)的l個濾波變換信號序列進(jìn)行加窗補零后，做l_fft階短時傅里葉變換，使其結(jié)果形成對應(yīng)一個掃查包的m*l_fft頻域矩陣，m表示每個掃查包設(shè)置的掃查線數(shù)量；l_fft表示短時傅里葉變換的階數(shù)。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述步驟s4具體包括：
m1、獲取當(dāng)前深度下，m*l_fft頻域矩陣對應(yīng)的中心頻率位置；
m2、以獲取的中心頻率位置為中心，獲取當(dāng)前深度下，感興趣帶寬的截止頻率位置；
m3、根據(jù)獲取的感興趣帶寬的截止頻率位置獲取其對應(yīng)的m*k頻域矩陣，k表示感興趣帶寬的截止頻率所對應(yīng)的列的總和；
m4、獲取m*k頻域矩陣中，每列數(shù)據(jù)對應(yīng)的頻率分量。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述步驟s5包括：根據(jù)一維自相關(guān)算法，估計m*k頻域矩陣中每列數(shù)據(jù)在頻率分量f(k)上的平均速度和能量，并將其集合標(biāo)記為第一平均速度序列，和第一能量序列；
所述步驟s5和步驟s6之間，所述方法還包括：
n1、根據(jù)預(yù)設(shè)的高速度閾值，低速度閾值，噪聲能量閾值和組織能量閾值對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；
若相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值大于所述高速度閾值，對應(yīng)能量序列中的能量值小于噪聲能量閾值；
或相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值小于所述低速度閾值，對應(yīng)能量序列中的能量值大于組織能量閾值；
則剔除相同頻率分量對應(yīng)第一平均速度序列和第一能量序列對應(yīng)的中的平均速度值和能量值，形成第二平均速度序列和第二能量序列。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述步驟n1后，所述方法還包括：
n2、獲取第二平均速度序列對應(yīng)的方差值及速度平均值，并根據(jù)其對第二平均速度序列和第二能量序列進(jìn)行第二次預(yù)處理；
若第二平均速度序列中的任一平均速度值與第二平均速度序列的速度平均值差值的平方大于第二平均速度的方差值與第一預(yù)設(shè)常數(shù)值的乘積，則剔除相同頻率分量對應(yīng)第二平均速度序列和第二能量序列中的平均速度值和能量值，形成第三平均速度序列和第三能量序列。
為了實現(xiàn)上述發(fā)明目的之一，本發(fā)明一實施方式提供一種超聲多普勒血流成像的處理系統(tǒng)，所述系統(tǒng)包括：信號獲取模塊，用于分別獲取每個掃查包下、每條掃查線上、每個采樣點對應(yīng)的回波rf信號；分別對獲取的回波rf信號做希爾伯特變換以形成正交iq復(fù)信號；
濾波處理模塊，在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，以獲得對應(yīng)的濾波變換信號；
傅里葉處理模塊，在快時間方向，分別對獲得的濾波變換信號做短時傅里葉變換，使其結(jié)果形成頻域矩陣；
處理輸出模塊，用于獲取當(dāng)前所述頻域矩陣中，每個快時間方向上的頻率分量；
根據(jù)一維自相關(guān)算法，估計當(dāng)前頻域矩陣中每列數(shù)據(jù)在其對應(yīng)頻率分量上的平均速度和能量；
對當(dāng)前頻域矩陣中對應(yīng)頻率分量的平均速度和能量進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述傅里葉處理模塊具體用于：
在同一深度下，對應(yīng)每個掃查包、在每條掃查線的快時間方向上分別獲取l個采樣點對應(yīng)的濾波變換信號，l小于等于n，n表示每條掃查線對應(yīng)的原始采樣點數(shù)量；
對每條掃查線對應(yīng)的l個濾波變換信號序列進(jìn)行加窗補零后，做l_fft階短時傅里葉變換，使其結(jié)果形成對應(yīng)一個掃查包的m*l_fft頻域矩陣，m表示每個掃查包設(shè)置的掃查線數(shù)量；l_fft表示短時傅里葉變換的階數(shù)。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述處理輸出模塊具體用于：
獲取當(dāng)前深度下，m*l_fft頻域矩陣對應(yīng)的中心頻率位置；
以獲取的中心頻率位置為中心，獲取當(dāng)前深度下，感興趣帶寬的截止頻率位置；
根據(jù)獲取的感興趣帶寬的截止頻率位置獲取其對應(yīng)的m*k頻域矩陣，k表示感興趣帶寬的截止頻率所對應(yīng)的列的總和；
獲取m*k頻域矩陣中，每列數(shù)據(jù)對應(yīng)的頻率分量。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述處理輸出模塊還用于：
根據(jù)一維自相關(guān)算法，估計m*k頻域矩陣中每列數(shù)據(jù)在頻率分量f(k)上的平均速度和能量，并將其集合標(biāo)記為第一平均速度序列，和第一能量序列；
根據(jù)預(yù)設(shè)的高速度閾值，低速度閾值，噪聲能量閾值和組織能量閾值對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；
若相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值大于所述高速度閾值，對應(yīng)能量序列中的能量值小于噪聲能量閾值；
或相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值小于所述低速度閾值，對應(yīng)能量序列中的能量值大于組織能量閾值；
則剔除相同頻率分量對應(yīng)第一平均速度序列和第一能量序列對應(yīng)的中的平均速度值和能量值，形成第二平均速度序列和第二能量序列。
作為本發(fā)明一實施方式的進(jìn)一步改進(jìn)，所述處理輸出模塊還用于：
獲取第二平均速度序列對應(yīng)的方差值及速度平均值，并根據(jù)其對第二平均速度序列和第二能量序列進(jìn)行第二次預(yù)處理；
若第二平均速度序列中的任一平均速度值與第二平均速度序列的速度平均值差值的平方大于第二平均速度的方差值與第一預(yù)設(shè)常數(shù)值的乘積，則剔除相同頻率分量對應(yīng)第二平均速度序列和第二能量序列中的平均速度值和能量值，形成第三平均速度序列和第三能量序列。
與現(xiàn)有技術(shù)相比，本發(fā)明的超聲多普勒血流成像的處理方法及處理系統(tǒng)，考慮了組織衰減對超聲信號的影響，基于寬帶發(fā)射，提供更好的分辨率，同時，速度估計更穩(wěn)定準(zhǔn)確；提高了超聲成像設(shè)備的方便性和使用效率，提升了超聲圖像的質(zhì)量。
附圖說明
圖1是傳統(tǒng)的成像系統(tǒng)的整體模塊示意圖；
圖2是傳統(tǒng)的超聲多普勒血流成像的處理系統(tǒng)的模塊示意圖；
圖3是本發(fā)明一實施方式中超聲多普勒血流成像的處理方法的流程示意圖；
圖4是本發(fā)明一實施方式中超聲多普勒血流成像的處理系統(tǒng)的模塊示意圖；
圖5a、圖5b是本發(fā)明一具體示例中分別采用傳統(tǒng)的超聲多普勒血流成像的處理方法與本發(fā)明的超聲多普勒血流成像的處理方法獲得的血流圖像的對比示意圖；
圖5c、圖5d是本發(fā)明一具體示例中分別采用傳統(tǒng)的超聲多普勒血流成像的處理方法與本發(fā)明的超聲多普勒血流成像的處理方法獲得的速度譜圖的對比示意圖。
具體實施方式
以下將結(jié)合附圖所示的各實施方式對本發(fā)明進(jìn)行詳細(xì)描述。但這些實施方式并不限制本發(fā)明，本領(lǐng)域的普通技術(shù)人員根據(jù)這些實施方式所做出的結(jié)構(gòu)、方法、或功能上的變換均包含在本發(fā)明的保護(hù)范圍內(nèi)。
需要說明的是，本發(fā)明主要應(yīng)用于超聲設(shè)備，相應(yīng)的，所述待測物可為待測組織，在此不做詳細(xì)贅述。
結(jié)合圖1所示，多普勒成像系統(tǒng)的模塊示意圖；脈沖多普勒成像過程中；通過探頭向組織中發(fā)射脈沖信號，所述脈沖信號經(jīng)組織中反射形成超聲信號經(jīng)由探頭換能器的不同基元轉(zhuǎn)變?yōu)殡娔M信號，通過前放模塊放大，再由a/d數(shù)模轉(zhuǎn)換模塊轉(zhuǎn)換為數(shù)字信號；各個不同基元的數(shù)字信號經(jīng)過波束合成模塊，合成為射頻信號；射頻信號經(jīng)過固定頻率的正交解調(diào)后，將正交解調(diào)結(jié)果i/q信號送入相應(yīng)的處理模塊。
結(jié)合圖2所示，傳統(tǒng)的多普勒血流(cf)模式成像需經(jīng)過如下過程：在正交解調(diào)模塊中，回波射頻rf信號先與發(fā)射中心頻率為f0的正交信號相乘，再通過低通的基帶濾波器，獲得基帶iq信號，送入后面的壁濾波模塊，壁濾波一般為高通濾波器，主要用來慮除信號中低速運動的組織信號。壁濾波后的信號，送入后面的速度能量估計模塊計算速度及能量；再通過動態(tài)范圍壓縮模塊進(jìn)行對數(shù)壓縮，最后進(jìn)行血流顯示。
結(jié)合圖3所示，圖3為本發(fā)明一實施方式中超聲多普勒血流成像的處理方法的流程圖，所述方法包括：s1、分別獲取每個掃查包下、每條掃查線上、每個采樣點對應(yīng)的回波rf信號；分別對獲取的回波rf信號做希爾伯特變換以形成正交iq復(fù)信號；
本發(fā)明一具體示例中，為了方便描述，將每個掃查包下、每條掃查線上、每個采樣點分別對應(yīng)的正交iq復(fù)信號分別以iq(x,y)表示，所述(x,y)表示采樣點的坐標(biāo)。
需要理解的是，為了方便描述，以下出現(xiàn)的示例中僅以一個掃查包為例做具體介紹。相應(yīng)的，對于一個掃查包中的第m條掃查線，m取值為[1,2,3……m]，m表示一個掃查包內(nèi)掃查線的總數(shù)量，其對應(yīng)的正交iq復(fù)信號序列可表示為：{iq(m,1),iq(m,2),iq(m,3),…,iq(m,n)}，其中，n表示對應(yīng)任一條掃查線，其在快時間方向上的采樣點總數(shù)量，需要說明的是，在上述示例的實際應(yīng)用中，所述掃查包的數(shù)量以及每個掃查包中掃查線的數(shù)量均可以根據(jù)需要具體調(diào)節(jié)，在此不做詳細(xì)贅述。
進(jìn)一步的，本發(fā)明一實施方式中，所述方法還包括：s2、在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，以獲得對應(yīng)的濾波變換信號。
所述濾波處理通常為采用高通濾波器對數(shù)據(jù)進(jìn)行處理，其用于濾除正交iq復(fù)信號中低速運動的組織信號。接續(xù)上述示例，本發(fā)明一具體實施方式中，一個掃查包下的，m條掃查線對應(yīng)的正交iq復(fù)信號序列依次為：
{iq(1,1),iq(1,2),iq(1,3),…,iq(1,n)}
{iq(2,1),iq(2,2),iq(2,3),…,iq(2,n)}
{iq(3,1),iq(3,2),iq(3,3),…,iq(3,n)}
……
{iq(m,1),iq(m,2),iq(m,3),…,iq(m,n)}；
進(jìn)一步的，在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，獲得每條掃查線上的每個掃查點分別對應(yīng)的濾波變換信號序列依次為：
{w(1,1),w(1,2),w(1,3),…,w(1,n)}
{w(2,1),w(2,2),w(2,3),…,w(2,n)}
{w(3,1),w(3,2),w(3,3),…,w(3,n)}
……
{w(m,1),w(m,2),w(m,3),…,w(m,n)}，即獲取當(dāng)前數(shù)據(jù)包中的m*n數(shù)據(jù)矩陣對應(yīng)的濾波變換信號。
進(jìn)一步的，所述方法還包括：s3、在快時間方向，分別對獲得的濾波變換信號做短時傅里葉變換，使其結(jié)果形成頻域矩陣。
本發(fā)明一優(yōu)選實施方式中，所述步驟s3具體包括：p1、在同一深度下，對應(yīng)每個掃查包、在每條掃查線的快時間方向上分別獲取l個采樣點對應(yīng)的濾波變換信號，l小于等于n，n表示每條掃查線對應(yīng)的原始采樣點數(shù)量，l為偶數(shù)。
本發(fā)明優(yōu)選實施方式中，所述l＝k*fc/f0，其中，k為常數(shù)，其可以根據(jù)需要自行設(shè)置，一般其取值可為1、2或3；fc表示系統(tǒng)回波信號的采樣頻率，f0表示發(fā)射信號的中心頻率。
接續(xù)上述示例，本發(fā)明具體實施方式中，當(dāng)前深度以d表示，則對應(yīng)獲得每條掃查線上的l個掃查點分別對應(yīng)的濾波變換信號序列依次為：
{w(1,d-l/2+1),w(1,d-l/2+2),w(1,d-l/2+3),…,w(1,d-l/2)}
{w(2,d-l/2+1),w(2,d-l/2+2),w(2,d-l/2+3),…,w(2,d-l/2)}
{w(3,d-l/2+1),w(3,d-l/2+2),w(3,d-l/2+3),…,w(3,d-l/2)}
……
{w(m,d-l/2+1),w(m,d-l/2+2),w(m,d-l/2+3),…,w(m,d-l/2)}，即獲取當(dāng)前數(shù)據(jù)包中的m*l數(shù)據(jù)矩陣對應(yīng)的濾波變換信號。
p2、對每條掃查線對應(yīng)的l個濾波變換信號序列進(jìn)行加窗補零后，做短時傅里葉變換，使其結(jié)果形成對應(yīng)一個掃查包的m*l_fft頻域矩陣，m表示每個掃查包設(shè)置的掃查線數(shù)量；l_fft表示短時傅里葉變換的階數(shù)。通常情況下，短時傅里葉變換的階數(shù)為2的冪，其冪的值可根據(jù)需要自行設(shè)定，本發(fā)明優(yōu)選實施方式中，其冪取值范圍為[5,8]。
接續(xù)上述示例，本發(fā)明具體實施方式中，做l_fft階短時傅里葉變換，則獲得每個掃查包對應(yīng)的m*l_fft頻域矩陣為：
{s_fft(1,1),s_fft(1,2),s_fft(1,3),…,s_fft(1,l_fft)}
{s_fft(2,1),s_fft(2,2),s_fft(2,3),…,s_fft(2,l_fft)}
{s_fft(3,1),s_fft(3,2),s_fft(3,3),…,s_fft(3,l_fft)}
……
{s_fft(m,1),s_fft(m,2),s_fft(m,3),…,s_fft(m,l_fft)}。
進(jìn)一步的，本發(fā)明一實施方式中，所述方法還包括：s4、獲取當(dāng)前所述頻域矩陣中，每個快時間方向上的頻率分量。
本發(fā)明優(yōu)選實施方式中，所述步驟s4具體包括：m1、獲取當(dāng)前深度下，m*l_fft頻域矩陣對應(yīng)的中心頻率位置；本發(fā)明具體實施方式中，獲取m*l_fft頻域矩陣靠近中間一行的頻域數(shù)據(jù)，并獲取其對應(yīng)的幅值；本發(fā)明具體實施方式中，m為偶數(shù)，靠近中間一行的數(shù)據(jù)為第m/2行或m/2+1行。
繼續(xù)上述示例，以取第m/2行數(shù)據(jù)為例做具體介紹。相應(yīng)的，m*l_fft頻域矩陣中第m/2行的頻域數(shù)據(jù)表示為：
{s_fft(m/2,1),s_fft(m/2,2),s_fft(m/2,3),…,s_fft(m/2,l_fft)}，
其對應(yīng)的幅值序列為：
{p_fft(m/2,1),p_fft(m/2,2),p_fft(m/2,3),…,p_fft(m/2,l_fft)}；
進(jìn)一步的，以當(dāng)前的幅值序列為搜索基礎(chǔ)，在預(yù)設(shè)幅值范圍內(nèi)搜索幅值的最大值，將其確認(rèn)為m*l_fft頻域矩陣對應(yīng)的中心頻率位置，為了方便描述，將中心頻率位置對應(yīng)的列以n_fd表示。所述預(yù)設(shè)幅值范圍為[n_fl，n_f0]
n_fl＝(f0–b/2)*l_fft/fc；
n_f0＝f0*l_fft/fc；其中，b為預(yù)設(shè)的信號帶寬，f0、fc與上述示例表示相同。
m2、以獲取的中心頻率位置為中心，獲取當(dāng)前深度下，感興趣帶寬的截止頻率位置。
本發(fā)明具體實施方式中，以當(dāng)前的幅值序列為搜索基礎(chǔ)，以n_fd為搜索中心，向其兩側(cè)依次進(jìn)行搜索，至當(dāng)前幅值與中心頻率對應(yīng)的幅值的比值符合設(shè)置的要求或數(shù)據(jù)為止，分別獲得起止幅值，并將起止幅值分別對應(yīng)的列以n_frxl和n_frxh表示。
m3、根據(jù)獲取的感興趣帶寬的截止頻率位置獲取其對應(yīng)的m*k頻域矩陣，k表示感興趣帶寬的截止頻率所對應(yīng)的列的總和。
本發(fā)明具體實施方式中，以m*l_fft頻域矩陣為數(shù)據(jù)基礎(chǔ)，獲取從第n_frxl列起始至n_frxh列為止的m*k頻域矩陣；其中，k表示感興趣帶寬的截止頻率所對應(yīng)列的總和；其取值k＝n_frxl,n_frxl+1,…,n_frxh。
本發(fā)明具體實施方式中，為了方便描述，將所述m*k頻域矩陣中的任意一個數(shù)據(jù)以{s_fft-1(m,k)}表示，k＝[n_frxl,n_frxl+1,…,n_frxh]＝[1,2,3……k],相應(yīng)的，m*k頻域矩陣可表示為：
{s_fft-1(1,1),s_fft-1(1,2),s_fft-1(1,3),…,s_fft-1(1,k)}
{s_fft-1(2,1),s_fft-1(2,2),s_fft-1(2,),…,s_fft-1(2,k)}
{s_fft-1(3,1),s_fft-1(3,2),s_fft-1(3,3),…,s_fft-1(3,k)}
……
{s_fft-1(m,1),s_fft-1(m,2),s_fft-1(m,3),…,s_fft-1(m,k)}，
上述m*k頻域矩陣中的數(shù)據(jù)值與上述m*l_fft頻域矩陣中的第n_frxl列至第n_frxh列的數(shù)據(jù)值相同，例如：m*k頻域矩陣中的s_fft-1(1,1)，與上述m*l_fft頻域矩陣中s_fft(1,n_frxl)數(shù)值相同，上述m*k頻域矩陣中s_fft-1(m,k)，與上述m*l_fft頻域矩陣中s_fft(m,n_frxh)數(shù)值相同，在此不做一一贅述。
進(jìn)一步的，所述方法還包括：m4、獲取m*k頻域矩陣中，每列數(shù)據(jù)對應(yīng)的頻率分量。
本發(fā)明具體實施方式中，通過s_fft-1(m,k)對應(yīng)的具體數(shù)值為復(fù)數(shù)，可獲得m*k頻域矩陣中每個數(shù)據(jù)對應(yīng)的實部和虛部值。如下所示，
s_fft-1(m,k)＝rm(k)+i*im(k)，其中，rm(k)表示實部，im(k)表示虛部，i為虛部符號。
相應(yīng)的，m*k頻域矩陣中，第k列數(shù)據(jù)對應(yīng)的頻率分量可表示為：f(k)＝k*fc/l_fft。
進(jìn)一步的，本發(fā)明一實施方式中，所述方法還包括：s5、根據(jù)一維自相關(guān)算法，估計當(dāng)前頻域矩陣中每列數(shù)據(jù)在其對應(yīng)頻率分量上的平均速度和能量。
本發(fā)明優(yōu)選實施方式中，根據(jù)一維自相關(guān)算法，估計m*k頻域矩陣中每列數(shù)據(jù)在頻率分量f(k)上的平均速度和能量，并將其集合標(biāo)記為第一平均速度序列，和第一能量序列；相應(yīng)的，以公式表示如下：
如此，獲得m*k頻域矩陣中、k列數(shù)據(jù)分別在其對應(yīng)的頻率分量上的第一平均速度序列，和第一能量序列，其可以表示為：
{vf(1),vf(2),vf(3),……vf(k)}
{pf(1),pf(2),pf(3),……pf(k)}，
由于k＝[n_frxl,n_frxl+1,…,n_frxh]＝[1,2,3……k]，對于m*l_fft頻域矩陣，上式可表示為：
{vf(n_frxl),vf(n_frxl+1),vf(n_frxl+2),……vf(n_frxl)}
{pf(n_frxl),pf(n_frxl+1),pf(n_frxl+2),……pf(n_frxl)}。
進(jìn)一步的，本發(fā)明一實施方式中，所述方法還包括：s6、對當(dāng)前頻域矩陣中對應(yīng)頻率分量上的平均速度和能量進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
本發(fā)明一優(yōu)選實施方式中，為了使最終顯示的圖像更加清晰，所述步驟s6之前，所述方法還包括：n1、根據(jù)預(yù)設(shè)的高速度閾值，低速度閾值，噪聲能量閾值和組織能量閾值對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；若相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值大于所述高速度閾值，對應(yīng)能量序列中的能量值小于噪聲能量閾值；或相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值小于所述低速度閾值，對應(yīng)能量序列中的能量值大于組織能量閾值；則剔除相同頻率分量對應(yīng)第一平均速度序列和第一能量序列中的平均速度值和能量值，形成第二平均速度序列和第二能量序列。
本發(fā)明一具體實施方式中，接續(xù)上述示例，對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；
若vf(k)>tv_h且pf(k)<tp_l，剔除vf(k)和pf(k)；
若vf(k)<tv_l且pf(k)>tp_h，剔除vf(k)和pf(k)，
其中，tv_h表示高速度閾值，tv_l表示低速度閾值，tp_l表示噪聲能量閾值，tp_h表示組織能量閾值。
相應(yīng)的，形成的第二平均速度序列和第二能量序列，表示為：
{vf-1(1),vf-1(2),vf-1(3),……vf-1(p)}
{pf-1(1),pf-1(2),pf-1(3),……pf-1(p)}，其中，p小于等于k，第二平均速度序列中的任一平均速度以pf-1(p)表示，p＝1,2,3……p。
優(yōu)選的，本發(fā)明一實施方式中，所述方法還包括：n2，獲取第二平均速度序列對應(yīng)的方差值及速度平均值，并根據(jù)其對第二平均速度序列和第二能量序列進(jìn)行第二次預(yù)處理；若第二平均速度序列中的任一平均速度值與第二平均速度序列的速度平均值差值的平方大于第二平均速度的方差值與第一預(yù)設(shè)常數(shù)值的乘積，則剔除相同頻率分量對應(yīng)第二平均速度序列和第二能量序列中的平均速度值和能量值，形成第三平均速度序列和第三能量序列。
若(vf-1(p)-v_ave)2＞d_ave*k_d，剔除vf-1(p)和pf-1(p)；其中，vf-1(p)表示第二平均速度序列中的任一平均速度值，pf-1(p)表示第二能量序列中的任一能量值，v_ave表示第二平均速度序列的速度平均值，d_ave表示第二平均速度序列的方差值，k_d為預(yù)設(shè)常識值。
相應(yīng)的，形成的第三平均速度序列，和第三能量序列，表示為：
{vf-2(1),vf-2(2),vf-2(3),……vf-2(q)}
{pf-2(1),pf-2(2),pf-2(3),……pf-2(q)}，其中，q小于等于p，第三平均速度序列中的任一平均速度以pf-2(q)表示，q＝1,2,3……q。
進(jìn)一步的，對第三平均速度序列和第三能量序列進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
本實施方式中，用于最終的頻譜顯示的最終能量和最終速度可分別表示為：
進(jìn)一步的，將獲得的最終能量和最終速度進(jìn)行動態(tài)壓縮，以進(jìn)行輸出顯示，在此不做詳細(xì)贅述。
結(jié)合圖5a-圖5d所示，圖5a、圖5c為采用傳統(tǒng)的超聲多普勒血流成像的處理方法分別獲得的血流圖像和速度譜圖；圖5b、圖5d為采用本發(fā)明的超聲多普勒血流成像的處理方法分別獲得的血流圖像和速度譜圖；
通過比對圖5a和圖5b可知：對組織中兩支距離較近且細(xì)的血管進(jìn)行血流成像過程中，采用本發(fā)明獲得的血流圖像的分辨率更高，且能明顯區(qū)分出組織中兩支血管。
通過比對圖5c和圖5d可知：其中橫軸代表深度方向，縱軸代表速度大小，當(dāng)發(fā)射信號為較短時，對于同樣的2個散射子，采用傳統(tǒng)的超聲多普勒血流成像的處理方法獲得的速度估計不穩(wěn)定，在深度方向分布較寬且誤差較大，而采用本發(fā)明的超聲多普勒血流成像的處理方法獲得的速度估計更準(zhǔn)確且分辨率更高。
結(jié)合圖4所示，本發(fā)明一實施方式中提供的超聲多普勒血流成像的處理系統(tǒng)，所述系統(tǒng)包括：信號獲取模塊100、濾波處理模塊200、傅里葉處理模塊300、處理輸出模塊400。
信號獲取模塊100用于分別獲取每個掃查包下、每條掃查線上、每個采樣點對應(yīng)的回波rf信號；分別對獲取的回波rf信號做希爾伯特變換以形成正交iq復(fù)信號；
本發(fā)明一具體示例中，為了方便描述，將每個掃查包下、每條掃查線上、每個采樣點分別對應(yīng)的正交iq復(fù)信號分別以iq(x,y)表示，所述(x,y)表示采樣點的坐標(biāo)。
需要理解的是，為了方便描述，以下出現(xiàn)的示例中僅以一個掃查包為例做具體介紹。相應(yīng)的，對于一個掃查包中的第m條掃查線，m取值為[1,2,3……m]，m表示一個掃查包內(nèi)掃查線的總數(shù)量，其對應(yīng)的正交iq復(fù)信號序列可表示為：{iq(m,1),iq(m,2),iq(m,3),…,iq(m,n)}，其中，n表示對應(yīng)任一條掃查線，其在快時間方向上的采樣點總數(shù)量，需要說明的是，在上述示例的實際應(yīng)用中，所述掃查包的數(shù)量以及每個掃查包中掃查線的數(shù)量均可以根據(jù)需要具體調(diào)節(jié)，在此不做詳細(xì)贅述。
本發(fā)明一實施方式中，濾波處理模塊200用于在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，以獲得對應(yīng)的濾波變換信號。
所述濾波處理通常為采用高通濾波器對數(shù)據(jù)進(jìn)行處理，其用于濾除正交iq復(fù)信號中低速運動的組織信號。接續(xù)上述示例，本發(fā)明一具體實施方式中，一個掃查包下的，m條掃查線對應(yīng)的正交iq復(fù)信號序列依次為：
{iq(1,1),iq(1,2),iq(1,3),…,iq(1,n)}
{iq(2,1),iq(2,2),iq(2,3),…,iq(2,n)}
{iq(3,1),iq(3,2),iq(3,3),…,iq(3,n)}
……
{iq(m,1),iq(m,2),iq(m,3),…,iq(m,n)}；
進(jìn)一步的，濾波處理模塊200在慢時間方向上，分別對獲得的正交iq復(fù)信號做濾波處理，獲得每條掃查線上的每個掃查點分別對應(yīng)的濾波變換信號序列依次為：
{w(1,1),w(1,2),w(1,3),…,w(1,n)}
{w(2,1),w(2,2),w(2,3),…,w(2,n)}
{w(3,1),w(3,2),w(3,3),…,w(3,n)}
……
{w(m,1),w(m,2),w(m,3),…,w(m,n)}，即獲取當(dāng)前數(shù)據(jù)包中的m*n數(shù)據(jù)矩陣對應(yīng)的濾波變換信號。
本發(fā)明一實施方式中，傅里葉處理模塊300用于在快時間方向，分別對獲得的濾波變換信號做短時傅里葉變換，使其結(jié)果形成頻域矩陣。
本發(fā)明一優(yōu)選實施方式中，傅里葉處理模塊300具體用于：在同一深度下，對應(yīng)每個掃查包、在每條掃查線的快時間方向上分別獲取l個采樣點對應(yīng)的濾波變換信號，l小于等于n，n表示每條掃查線對應(yīng)的原始采樣點數(shù)量，l為偶數(shù)。
本發(fā)明優(yōu)選實施方式中，所述l＝k*fc/f0，其中，k為常數(shù)，其可以根據(jù)需要自行設(shè)置，一般其取值可為1、2或3；fc表示系統(tǒng)回波信號的采樣頻率，f0表示發(fā)射信號的中心頻率。
接續(xù)上述示例，本發(fā)明具體實施方式中，當(dāng)前深度以d表示，則對應(yīng)獲得每條掃查線上的l個掃查點分別對應(yīng)的濾波變換信號序列依次為：
{w(1,d-l/2+1),w(1,d-l/2+2),w(1,d-l/2+3),…,w(1,d-l/2)}
{w(2,d-l/2+1),w(2,d-l/2+2),w(2,d-l/2+3),…,w(2,d-l/2)}
{w(3,d-l/2+1),w(3,d-l/2+2),w(3,d-l/2+3),…,w(3,d-l/2)}
……
{w(m,d-l/2+1),w(m,d-l/2+2),w(m,d-l/2+3),…,w(m,d-l/2)}，即獲取當(dāng)前數(shù)據(jù)包中的m*l數(shù)據(jù)矩陣對應(yīng)的濾波變換信號。
傅里葉處理模塊300還用于：對每條掃查線對應(yīng)的l個濾波變換信號序列進(jìn)行加窗補零后，做短時傅里葉變換，使其結(jié)果形成對應(yīng)一個掃查包的m*l_fft頻域矩陣，m表示每個掃查包設(shè)置的掃查線數(shù)量；l_fft表示短時傅里葉變換的階數(shù)。通常情況下，短時傅里葉變換的階數(shù)為2的冪，其冪的值可根據(jù)需要自行設(shè)定，本發(fā)明優(yōu)選實施方式中，其冪取值范圍為[5,8]。
接續(xù)上述示例，本發(fā)明具體實施方式中，做l_fft階短時傅里葉變換，則獲得每個掃查包對應(yīng)的m*l_fft頻域矩陣為：
{s_fft(1,1),s_fft(1,2),s_fft(1,3),…,s_fft(1,l_fft)}
{s_fft(2,1),s_fft(2,2),s_fft(2,3),…,s_fft(2,l_fft)}
{s_fft(3,1),s_fft(3,2),s_fft(3,3),…,s_fft(3,l_fft)}
……
{s_fft(m,1),s_fft(m,2),s_fft(m,3),…,s_fft(m,l_fft)}。
本發(fā)明一實施方式中，處理輸出模塊400用于獲取當(dāng)前所述頻域矩陣中，每個快時間方向上的頻率分量。
本發(fā)明優(yōu)選實施方式中，處理輸出模塊400具體用于：獲取當(dāng)前深度下，m*l_fft頻域矩陣對應(yīng)的中心頻率位置；本發(fā)明具體實施方式中，處理輸出模塊400獲取m*l_fft頻域矩陣靠近中間一行的頻域數(shù)據(jù)，并獲取其對應(yīng)的幅值；本發(fā)明具體實施方式中，m為偶數(shù)，靠近中間一行的數(shù)據(jù)為第m/2行或m/2+1行。
繼續(xù)上述示例，以取第m/2行數(shù)據(jù)為例做具體介紹。相應(yīng)的，m*l_fft頻域矩陣中第m/2行的頻域數(shù)據(jù)表示為：
{s_fft(m/2,1),s_fft(m/2,2),s_fft(m/2,3),…,s_fft(m/2,l_fft)}，
其對應(yīng)的幅值序列為：
{p_fft(m/2,1),p_fft(m/2,2),p_fft(m/2,3),…,p_fft(m/2,l_fft)}；
進(jìn)一步的，處理輸出模塊400以當(dāng)前的幅值序列為搜索基礎(chǔ)，在預(yù)設(shè)幅值范圍內(nèi)搜索幅值的最大值，將其確認(rèn)為m*l_fft頻域矩陣對應(yīng)的中心頻率位置，為了方便描述，將中心頻率位置對應(yīng)的列以n_fd表示。所述預(yù)設(shè)幅值范圍為[n_fl，n_f0]
n_fl＝(f0–b/2)*l_fft/fc；
n_f0＝f0*l_fft/fc；其中，b為預(yù)設(shè)的信號帶寬，f0、fc與上述示例表示相同。
進(jìn)一步的，處理輸出模塊400還用于：以獲取的中心頻率位置為中心，獲取當(dāng)前深度下，感興趣帶寬的截止頻率位置。
本發(fā)明具體實施方式中，處理輸出模塊400以當(dāng)前的幅值序列為搜索基礎(chǔ)，以n_fd為搜索中心，向其兩側(cè)依次進(jìn)行搜索，至當(dāng)前幅值與中心頻率對應(yīng)的幅值的比值符合設(shè)置的要求或數(shù)據(jù)為止，分別獲得起止幅值，并將起止幅值分別對應(yīng)的列以n_frxl和n_frxh表示。
進(jìn)一步的，處理輸出模塊400還用于：根據(jù)獲取的感興趣帶寬的截止頻率位置獲取其對應(yīng)的m*k頻域矩陣，k表示感興趣帶寬的截止頻率所對應(yīng)的列的總和。
本發(fā)明具體實施方式中，處理輸出模塊400以m*l_fft頻域矩陣為數(shù)據(jù)基礎(chǔ)，獲取從第n_frxl列起始至n_frxh列為止的m*k頻域矩陣；其中，k表示感興趣帶寬的截止頻率所對應(yīng)列的總和；其取值k＝n_frxl,n_frxl+1,…,n_frxh。
本發(fā)明具體實施方式中，為了方便描述，將所述m*k頻域矩陣中的任意一個數(shù)據(jù)以{s_fft-1(m,k)}表示，k＝[n_frxl,n_frxl+1,…,n_frxh]＝[1,2,3……k],相應(yīng)的，m*k頻域矩陣可表示為：
{s_fft-1(1,1),s_fft-1(1,2),s_fft-1(1,3),…,s_fft-1(1,k)}
{s_fft-1(2,1),s_fft-1(2,2),s_fft-1(2,),…,s_fft-1(2,k)}
{s_fft-1(3,1),s_fft-1(3,2),s_fft-1(3,3),…,s_fft-1(3,k)}
……
{s_fft-1(m,1),s_fft-1(m,2),s_fft-1(m,3),…,s_fft-1(m,k)}，
上述m*k頻域矩陣中的數(shù)據(jù)值與上述m*l_fft頻域矩陣中的第n_frxl列至第n_frxh列的數(shù)據(jù)值相同，例如：m*k頻域矩陣中的s_fft-1(1,1)，與上述m*l_fft頻域矩陣中s_fft(1,n_frxl)數(shù)值相同，上述m*k頻域矩陣中s_fft-1(m,k)，與上述m*l_fft頻域矩陣中s_fft(m,n_frxh)數(shù)值相同，在此不做一一贅述。
進(jìn)一步的，處理輸出模塊400還用于：獲取m*k頻域矩陣中，每列數(shù)據(jù)對應(yīng)的頻率分量。
本發(fā)明具體實施方式中，處理輸出模塊400通過s_fft-1(m,k)對應(yīng)的具體數(shù)值為復(fù)數(shù)，可獲得m*k頻域矩陣中每個數(shù)據(jù)對應(yīng)的實部和虛部值。如下所示，
s_fft-1(m,k)＝rm(k)+i*im(k)，其中，rm(k)表示實部，im(k)表示虛部，i為虛部符號。
相應(yīng)的，m*k頻域矩陣中，第k列數(shù)據(jù)對應(yīng)的頻率分量可表示為：f(k)＝k*fc/l_fft。
本發(fā)明一實施方式中，處理輸出模塊400還用于：根據(jù)一維自相關(guān)算法，估計當(dāng)前頻域矩陣中每列數(shù)據(jù)在其對應(yīng)頻率分量上的平均速度和能量。
本發(fā)明優(yōu)選實施方式中，處理輸出模塊400根據(jù)一維自相關(guān)算法，估計m*k頻域矩陣中每列數(shù)據(jù)在頻率分量f(k)上的平均速度和能量，并將其集合標(biāo)記為第一平均速度序列，和第一能量序列；相應(yīng)的，以公式表示如下：
如此，獲得m*k頻域矩陣中、k列數(shù)據(jù)分別在其對應(yīng)的頻率分量上的第一平均速度序列，和第一能量序列，其可以表示為：
{vf(1),vf(2),vf(3),……vf(k)}
{pf(1),pf(2),pf(3),……pf(k)}，
由于k＝[n_frxl,n_frxl+1,…,n_frxh]＝[1,2,3……k]，對于m*l_fft頻域矩陣，上式可表示為：
{vf(n_frxl),vf(n_frxl+1),vf(n_frxl+2),……vf(n_frxl)}
{pf(n_frxl),pf(n_frxl+1),pf(n_frxl+2),……pf(n_frxl)}。
進(jìn)一步的，本發(fā)明一實施方式中，處理輸出模塊400還用于：對當(dāng)前頻域矩陣中對應(yīng)頻率分量上的平均速度和能量進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
本發(fā)明一優(yōu)選實施方式中，為了使最終顯示的圖像更加清晰，處理輸出模塊400還用于：根據(jù)預(yù)設(shè)的高速度閾值，低速度閾值，噪聲能量閾值和組織能量閾值對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；若相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值大于所述高速度閾值，對應(yīng)能量序列中的能量值小于噪聲能量閾值；或相同頻率分量對應(yīng)的第一平均速度序列中的平均速度值小于所述低速度閾值，對應(yīng)能量序列中的能量值大于組織能量閾值；則剔除相同頻率分量對應(yīng)第一平均速度序列和第一能量序列中的平均速度值和能量值，形成第二平均速度序列和第二能量序列。
本發(fā)明一具體實施方式中，接續(xù)上述示例，處理輸出模塊400對m*k頻域矩陣對應(yīng)的第一平均速度序列和第一能量序列進(jìn)行第一次預(yù)處理；
若vf(k)>tv_h且pf(k)<tp_l，剔除vf(k)和pf(k)；
若vf(k)<tv_l且pf(k)>tp_h，剔除vf(k)和pf(k)，
其中，tv_h表示高速度閾值，tv_l表示低速度閾值，tp_l表示噪聲能量閾值，tp_h表示組織能量閾值。
相應(yīng)的，形成的第二平均速度序列和第二能量序列，表示為：
{vf-1(1),vf-1(2),vf-1(3),……vf-1(p)}
{pf-1(1),pf-1(2),pf-1(3),……pf-1(p)}，其中，p小于等于k，第二平均速度序列中的任一平均速度以pf-1(p)表示，p＝1,2,3……p。
優(yōu)選的，本發(fā)明一實施方式中，處理輸出模塊400還用于：獲取第二平均速度序列對應(yīng)的方差值及速度平均值，并根據(jù)其對第二平均速度序列和第二能量序列進(jìn)行第二次預(yù)處理；若第二平均速度序列中的任一平均速度值與第二平均速度序列的速度平均值差值的平方大于第二平均速度的方差值與第一預(yù)設(shè)常數(shù)值的乘積，則剔除相同頻率分量對應(yīng)第二平均速度序列和第二能量序列中的平均速度值和能量值，形成第三平均速度序列和第三能量序列。
若(vf-1(p)-v_ave)2＞d_ave*k_d，剔除vf-1(p)和pf-1(p)；其中，vf-1(p)表示第二平均速度序列中的任一平均速度值，pf-1(p)表示第二能量序列中的任一能量值，v_ave表示第二平均速度序列的速度平均值，d_ave表示第二平均速度序列的方差值，k_d為預(yù)設(shè)常識值。
相應(yīng)的，形成的第三平均速度序列，和第三能量序列，表示為：
{vf-2(1),vf-2(2),vf-2(3),……vf-2(q)}
{pf-2(1),pf-2(2),pf-2(3),……pf-2(q)}，其中，q小于等于p，第三平均速度序列中的任一平均速度以pf-2(q)表示，q＝1,2,3……q。
進(jìn)一步的，處理輸出模塊400還用于：對第三平均速度序列和第三能量序列進(jìn)行線性擬合，獲取用于最終的頻譜顯示的最終能量和最終速度。
本實施方式中，用于最終的頻譜顯示的最終能量和最終速度可分別表示為：
進(jìn)一步的，本發(fā)明的超聲多普勒血流成像的處理系統(tǒng)還包括：動態(tài)壓縮模塊(未圖示)以及顯示模塊(未圖示)，所述動態(tài)壓縮模塊將獲得的最終能量和最終速度進(jìn)行動態(tài)壓縮，所述顯示模塊用于對經(jīng)過動態(tài)壓縮模塊處理的數(shù)據(jù)結(jié)果輸出顯示，在此不做詳細(xì)贅述。
所屬領(lǐng)域的技術(shù)人員可以清楚地了解到，為描述的方便和簡潔，上述描述的系統(tǒng)中模塊的具體工作過程，可以參考前述方法實施方式中的對應(yīng)過程，在此不再贅述。
綜上所述，本發(fā)明的超聲多普勒血流成像的處理方法及處理系統(tǒng)，考慮了組織衰減對超聲信號的影響，基于寬帶發(fā)射，提供更好的分辨率，同時，速度估計更穩(wěn)定準(zhǔn)確；提高了超聲成像設(shè)備的方便性和使用效率，提升了超聲圖像的質(zhì)量。
為了描述的方便，描述以上裝置時以功能分為各種模塊分別描述。當(dāng)然，在實施本申請時可以把各模塊的功能在同一個或多個軟件和/或硬件中實現(xiàn)。
應(yīng)當(dāng)理解，雖然本說明書按照實施方式加以描述，但并非每個實施方式僅包含一個獨立的技術(shù)方案，說明書的這種敘述方式僅僅是為清楚起見，本領(lǐng)域技術(shù)人員應(yīng)當(dāng)將說明書作為一個整體，各實施方式中的技術(shù)方案也可以經(jīng)適當(dāng)組合，形成本領(lǐng)域技術(shù)人員可以理解的其他實施方式。
上文所列出的一系列的詳細(xì)說明僅僅是針對本發(fā)明的可行性實施方式的具體說明，它們并非用以限制本發(fā)明的保護(hù)范圍，凡未脫離本發(fā)明技藝精神所作的等效實施方式或變更均應(yīng)包含在本發(fā)明的保護(hù)范圍之內(nèi)。
- 還沒有人留言評論。精彩留言會獲得點贊！