CN111681740A - 一种基于活体超声图像的呼吸分离式应变成像方法 - Google Patents
一种基于活体超声图像的呼吸分离式应变成像方法 Download PDFInfo
- Publication number
- CN111681740A CN111681740A CN202010731441.4A CN202010731441A CN111681740A CN 111681740 A CN111681740 A CN 111681740A CN 202010731441 A CN202010731441 A CN 202010731441A CN 111681740 A CN111681740 A CN 111681740A
- Authority
- CN
- China
- Prior art keywords
- image sequence
- state
- image
- strain
- gamma
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 32
- 230000029058 respiratory gaseous exchange Effects 0.000 title claims abstract description 22
- 238000000926 separation method Methods 0.000 title abstract description 10
- 238000000034 method Methods 0.000 claims abstract description 49
- 230000033001 locomotion Effects 0.000 claims abstract description 44
- 238000006073 displacement reaction Methods 0.000 claims abstract description 43
- 238000004364 calculation method Methods 0.000 claims abstract description 34
- 238000012545 processing Methods 0.000 claims abstract description 19
- 238000009826 distribution Methods 0.000 claims abstract description 17
- 238000012216 screening Methods 0.000 claims abstract description 9
- 238000010183 spectrum analysis Methods 0.000 claims abstract description 8
- 230000008569 process Effects 0.000 claims description 25
- 238000002604 ultrasonography Methods 0.000 claims description 21
- 230000003434 inspiratory effect Effects 0.000 claims description 12
- BTCSSZJGUNDROE-UHFFFAOYSA-N gamma-aminobutyric acid Chemical compound NCCCC(O)=O BTCSSZJGUNDROE-UHFFFAOYSA-N 0.000 claims description 9
- 238000010438 heat treatment Methods 0.000 claims description 7
- 238000001228 spectrum Methods 0.000 claims description 6
- 238000009499 grossing Methods 0.000 claims description 5
- 230000036391 respiratory frequency Effects 0.000 claims description 3
- 238000001727 in vivo Methods 0.000 claims description 2
- 230000000737 periodic effect Effects 0.000 claims description 2
- 230000007547 defect Effects 0.000 abstract description 2
- 230000035565 breathing frequency Effects 0.000 abstract 1
- 230000001121 heart beat frequency Effects 0.000 abstract 1
- 210000001519 tissue Anatomy 0.000 description 39
- 239000000523 sample Substances 0.000 description 12
- 238000002679 ablation Methods 0.000 description 7
- 230000001629 suppression Effects 0.000 description 7
- 230000008859 change Effects 0.000 description 6
- 238000010586 diagram Methods 0.000 description 6
- 230000000241 respiratory effect Effects 0.000 description 6
- 238000012544 monitoring process Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 4
- 230000009471 action Effects 0.000 description 3
- 230000001186 cumulative effect Effects 0.000 description 3
- 238000002091 elastography Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000001914 filtration Methods 0.000 description 3
- 230000007246 mechanism Effects 0.000 description 3
- 238000003672 processing method Methods 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 230000006835 compression Effects 0.000 description 2
- 238000007906 compression Methods 0.000 description 2
- 238000002591 computed tomography Methods 0.000 description 2
- 238000002059 diagnostic imaging Methods 0.000 description 2
- 239000011159 matrix material Substances 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 238000003786 synthesis reaction Methods 0.000 description 2
- 238000002560 therapeutic procedure Methods 0.000 description 2
- 238000012285 ultrasound imaging Methods 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 241000282887 Suidae Species 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 210000000577 adipose tissue Anatomy 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000005452 bending Methods 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000000747 cardiac effect Effects 0.000 description 1
- 239000000919 ceramic Substances 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 239000013078 crystal Substances 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 230000005284 excitation Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000005247 gettering Methods 0.000 description 1
- 238000000608 laser ablation Methods 0.000 description 1
- 238000002595 magnetic resonance imaging Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000035790 physiological processes and functions Effects 0.000 description 1
- 238000007674 radiofrequency ablation Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000003938 response to stress Effects 0.000 description 1
- 230000002269 spontaneous effect Effects 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H30/00—ICT specially adapted for the handling or processing of medical images
- G16H30/20—ICT specially adapted for the handling or processing of medical images for handling medical images, e.g. DICOM, HL7 or PACS
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10132—Ultrasound image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20092—Interactive image processing based on input by user
- G06T2207/20104—Interactive definition of region of interest [ROI]
Landscapes
- Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Radiology & Medical Imaging (AREA)
- General Health & Medical Sciences (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Epidemiology (AREA)
- Primary Health Care (AREA)
- Public Health (AREA)
- Quality & Reliability (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Ultra Sonic Daignosis Equipment (AREA)
Abstract
本发明公开了一种基于活体超声图像的呼吸分离式应变成像方法，属于医学图像处理领域。本发明的方法为采集数字超声图像序列后，通过二维互相关计算和频谱分析得到活体的呼吸、心跳频率，随后在互相关曲线中进行交替极值检索，以对图像序列进行“呼气”、“吸气”状态划分，并进行初步的图像筛选，再分别提取运动状态最匹配的图像序列及其对应的运动补偿量，并据此分别计算组织的空间位移和空间应变图像序列。最后，将两个状态的位移和应变图像序列分别合并。本发明的目的在于克服现有技术中，活体呼吸、心跳等生理运动会影响斑点追踪精度，导致空间应变图像存在较大伪影和误差的不足，本发明可以得到精确的组织内位移和应变分布图像。
Description
技术领域
本发明涉及医学图像处理领域，更具体地说，涉及一种基于活体超声图像的呼吸分离式应变成像方法。
背景技术
随着医学影像技术的发展和治疗方式微创、无创化的趋势，医学影像在无创和介入式治疗的监控中应用愈加广泛，典型如局部组织热消融中的温度成像、弹性成像等功能成像技术。常用的成像方式包括磁共振成像、计算机断层扫描成像(CT)和超声成像等。以B超为主的超声成像设备由于兼具低成本、便携和实时性强等优点，已广泛应用于临床。
在利用B超进行活体成像时，超声图像由于活体的呼吸、心跳等生理过程表现出复杂的动态变化。基于此类超声图像进行超声热应变成像、弹性成像时，其结果中会存在较大的伪影和误差。为减小生理运动的影响，现有的运动抑制方法主要有图像配准法和呼吸门控法。但是，生理运动时组织的弯曲和压缩使图像配准法的运动抑制效果十分有限。对于呼吸门控法，由于个体的呼吸幅度、周期等存在较大差异，其运动抑制效果不具有很好的效果和普适性。此外，该方法一般不考虑心跳运动的影响。
在现有技术中针对克服活体生理运动、计算局部组织位移和应变，提出了一些技术方案。例如发明创造名称为：基于热膨胀和门控算法测量生物组织温度变化的超声方法(申请日：2017年9月25日；申请号：201710876349)，该方案公开了一种基于热膨胀和门控算法测量生物组织温度变化的超声方法，其针对目前各种基于对体内组织加热的方式治疗疾病的方法中无法有效监测靶区温升的问题，建立了一种利用B超RF信号评估生物组织温度变化的方法。该方案对生物组织利用聚焦超声、射频、微波等方法进行局部加热，用B型超声对靶区进行成像并收集其RF信号，基于B超时序图像，选取目标帧，计算超声经过组织时的时间延迟图像并由此得到温度变化图像；根据加热区域外的图像计算自适应滤波器的系数，对得到的温度变化图像进行噪声抑制。该方案在温升18℃范围内误差不超过2℃，其将推动B超的温升监控技术在热疗中的应用，可显著提高热疗的安全性和有效性。但该方案不足之处在于仅能在“呼气”和“吸气”两种状态之间二选其一，即仅能输出一种呼吸状态时期内的结果；由于缺少运动状态、周期的判断，可能导致门控失效，导致位移和应变计算误差增大，甚至发生错误。
综上所述，如何对活体生理运动进行有效的抑制，获取精确的组织位移和应变分布图像，是现有技术亟需解决的问题。
发明内容
1.要解决的问题
本发明的目的在于克服现有技术中，由于活体生理运动的存在导致超声图像移动且难以被有效抑制，进一步导致位移和和应变分布图像的计算结果伴有伪影和较大误差的不足，本发明提供了一种基于活体超声图像的呼吸分离式应变成像方法，可以获得活体图像在不同呼吸状态下的热应变，并且可以对活体图像序列进行精确的生理运动识别和补偿，进一步以此计算得到精确的组织位移分布图像和组织应变分布图像。
2.技术方案
为了解决上述问题，本发明所采用的技术方案如下：
本发明的一种基于活体超声图像的呼吸分离式应变成像方法，包括：对活体的应变集中区域进行超声图像采集得到数字超声图像序列；对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γn(t)，根据γn(t)选取“吸气”状态的第一参考帧；再根据“吸气”状态的第一参考帧对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γm(t)，对γm(t)进行处理得到互相关系数—时间曲线γ'm(t)，并根据γ'm(t)选取“呼气”状态的第一参考帧；查找γ'm(t)的极值，并根据极值时间坐标划分“吸气”状态和“呼气”状态；其中，“吸气”状态和“呼气”状态分别包括若干个周期；针对每个周期分别设置阈值，根据阈值对图像序列进行筛选得到剩余图像序列，其中，“吸气”状态对应的剩余图像序列为Iinh,0；“呼气”状态对应的剩余图像序列为Iexh,0；对数字超声图像序列设置感兴趣区域，并对感兴趣区域进行带运动补偿的互相关检索得到处理图像序列，其中，“吸气”状态对应的处理图像序列为Iinh,1；“呼气”状态对应的处理图像序列为Iexh,1；对Iinh,1和Iexh,1分别进行计算得到各自对应的位移图像序列和应变图像序列，其中，Iinh,1对应位移图像序列Dinh和应变图像序列Sinh；Iexh,1对应位移图像序列Dexh和应变图像序列Sexh；将Iinh,1和Iexh,1各自对应的位移图像序列和应变图像序列分别进行合并得到最终的组织位移图像序列和组织应变图像序列。
更进一步地，对活体的应变区域进行图像采集的具体过程为：对活体的靶区加热或施加外力，使得靶区产生应变集中区域；对靶区的应变集中区域进行等时间间隔的连续成像，采集得到按时间排列的数字超声图像序列，相邻图像的采集时间间隔为T0。
更进一步地，对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γn(t)的具体过程中：数字超声图像序列共N帧，第a帧图像根据时间坐标记为aT0，a∈[1,2,…,N]；取数字超声图像序列的第1～M帧图像，其中，第1～M帧图像包括至少3个完整的呼吸过程；以第n帧为参考帧Bref，将所有N帧图像依次作为目标帧Bex，计算参考帧和目标帧的二维互相关系数γn，并获得M条互相关系数—时间曲线γn(t)；其中，γn(t)为γn随时间t变化的曲线，n∈[1,2,…,M]。
更进一步地，选取“呼气”状态和“吸气”状态的第一参考帧的具体过程为：
对γn(t)进行频谱分析得到频域的幅度谱Sn(f)，f为频率变量，记Sn(f)的最大值为An；若[A1,A2,…，AM]中的最大值为Am，m∈[1,2,…,M]，则将Sm(f)的峰值频率记为呼吸频率fres，并将第m帧选为“吸气”状态第一参考帧；
以第m帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γm；并获得互相关系数—时间曲线γm(t)；其中，γm(t)为γm随时间t变化的曲线；
对γm(t)进行平滑降噪处理得到互相关系数—时间曲线γ'm(t)；若γ'm(t)最小值所处的时间坐标t＝kT0，则将第k帧选为“呼气”状态的第一参考帧；其中，k为整数。
更进一步地，查找互相关系数—时间曲线中的极值的具体过程为：极值包括极小值和极大值，以t＝mT0为起点和第一个极大值位置，先向正方向以步长Δtc交替检索极小值和极大值；然后在负方向以步长Δtc交替检索极小值和极大值，直至互相关系数—时间曲线γ'm(t)被检索完毕；其中，0.5/fres<Δtc<1/fres。
更进一步地，利用傅里叶变换或线性调频z变换对γn(t)进行频谱分析得到频域的幅度谱Sn(f)。
更进一步地，根据极值位置划分“吸气”状态和“呼气”状态的具体过程为：将所有极小值位置标记为“吸气”状态的周期边界，任意两个相邻极小值之间为“吸气”状态；所有极大值位置标记为“呼气”状态的周期边界，任意两个相邻极大值之间为“呼气”状态。
更进一步地，针对每个周期设置阈值，并根据阈值对图像序列进行筛选得到剩余图像序列的具体过程为：在“吸气”状态的每个周期内，将γm(t)的所有值按从大到小排列，取周期内前C1％数据中的最小值为该周期的相关性阈值；在“呼气”状态的每个周期内，将γk(t)的所有值按从大到小排列，取周期内前C1％数据中的最小值为该周期的相关性阈值；C1的取值范围均为1～100；其中，γk为以第k帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γk，γk(t)为γk随时间t变化的曲线；在“吸气”状态各周期对应的时间内，将γm(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iinh,0；在“呼气”状态各周期对应的时间内，将γk(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iexh,0。
更进一步地，计算得到处理图像序列的具体过程为：
(1)对数字超声图像序列预设统一的感兴趣区域ROI，ROI的尺寸为I*J，左下角顶点位置为i0，j0，使得所有图像中靶区均处于ROI内，I、J、i0、j0的单位均为像素；
(2)在“吸气”状态中，将Iinh,0中的第1帧记为参考帧Rref；将Rref的ROI区域作为一个新图像，置入处理图像序列Iinh,1；
(3)将位于Rref之后[Δt1,Δτ]范围内的所有帧依次记为目标帧Rex；其中，Δt1为选帧范围下限，Δτ为选帧范围上限，根据应变量级确定Δτ，要求Δt2>Δτ>Δt1，Δt2为Δτ选择范围上限，且Δt1>0.1/fres，Δt2<1.5/fres；
(4)将每个Rex的ROI横向移动i1，纵向移动j1；其中，-I/2≤i1≤I/2，-J/2≤j1≤J/2，i1和j1为整数；对Rref和Rex的ROI区域计算二维互相关系数γR(i1,j1)；若γR(i1,j1)的最大值γR,max对应于i1＝i'1，j1＝j'1，则当前目标帧的运动补偿量为i'1和j'1；
(5)对所有目标帧Rex计算γR,max和[i'1,j'1]；选择γR,max的最大值对应的目标帧记为新参考帧，并将新参考帧经运动补偿后的ROI区域作为一个新图像，置于Iinh,1中已有帧之后；其中，经运动补偿后的ROI区域的横向范围为：i0+i'1～i0+i'1+I-1，纵向范围为j0+j'1～j0+j'1+J-1；
(6)重复步骤(3)～(5)，直至Iinh,0中的所有帧被处理完毕，获得处理图像序列Iinh,1；
依据步骤(2)～(6)，对“呼气”状态对应的图像序列处理得到处理图像序列为Iexh,1。
更进一步地，合并得到最终的组织位移图像序列和组织应变图像序列的具体过程为：
(i)将Dinh中所有图像的时间坐标记为集合tinh，Dexh中所有图像的时间坐标记为集合texh；将Dinh和Dexh分别进行时域插值，使插值后的图像序列D'inh和D'exh中所有图像时间坐标的集合均为tinh和texh的合集tall；
(ii)将D'inh和D'exh中时间坐标相同的帧进行求和平均计算，将计算结果以同一时间坐标置入组织位移图像序列Dall；
依据步骤(i)和步骤(ii)对应变图像序列Sinh和Sexh进行处理得到组织应变图像序列Sall；
以tall为时序的组织位移图像序列Dall和组织应变图像序列Sall即为最终获得的组织位移分布图像和组织应变分布图像。
3.有益效果
相比于现有技术，本发明的有益效果为：
本发明的一种基于活体超声图像的呼吸分离式应变成像方法，采用了呼吸分离及再合成的处理方法，可获得“呼气”和“吸气”两种状态下的热应变，从而使得数据输出的时间间隔更小、更均匀，应用的实时性更好；进一步在图像运动状态匹配及位移和应变的计算中引入运动补偿机制，其结果可具有更好的精度，从而可以得到精确的组织位移分布图像和组织应变分布图像。
附图说明
图1为基于活体超声图像的呼吸分离式应变成像方法的流程示意图；
图2为实施例2中对生物组织加热过程采集B超图像的示意图；
图3为实施例2中通过交替检索峰值得到的呼吸、心跳周期划分结果示意图；其中，图3a为呼吸周期划分结果示意图，图3b为心跳周期划分结果示意图；
图4a为实施例2中“吸气”状态下不同时间的B超灰度图像示意图；
图4b为实施例2中“呼气”状态下不同时间的B超灰度图像示意图；
图5a为实施例2中“吸气”状态下的二维累积应变分布示意图；
图5b为实施例2中“呼气”状态下的二维累积应变分布示意图；
图6为实施例2中“吸气”状态、“呼气”状态通过插值后平均得到的应变-时间最大值曲线的对比示意图。
附图标记：100、超声成像系统；200、成像探头；300、计算机；400、微波消融针；500、微波加热器。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整的描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例；而且，各个实施例之间不是相对独立的，根据需要可以相互组合，从而达到更优的效果。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
为进一步了解本发明的内容，结合附图和实施例对本发明作详细描述。
实施例1
结合图1所示，本发明的一种基于活体超声图像的呼吸分离式应变成像方法，包括以下步骤：
1)采集数字超声图像序列
对活体的应变集中区域进行超声图像采集得到数字超声图像序列；具体地，先对活体的靶区加热或施加外力作用，使靶区内形成集中的应变；需要说明的是，本发明适用多种应用场景需求，如组织局部温升监控、弹性成像、组织散射体运动追踪等。进一步值得说明的是，对生物组织施加外部作用时，作用的具体形式不限，可为微波消融、射频消融、超声消融、激光消融和红外线消融，也可为对组织内部或外部施加静态或准静态的机械激励，或是利用组织内部自身存在的压力产生应力响应等。
进一步地，对靶区的应变集中区域进行超声监控成像，采集得到按时间排列的数字超声图像序列。值得说明的是，监控成像时不限定活体呼吸和心跳的方式，进一步地，呼吸可为自主呼吸、无创或有创呼吸机支持的呼吸、人工支持的呼吸等；心跳可为自主心跳、心脏起搏器支持的心跳、体外医疗设备支持的心跳等。进一步地，本发明不限定超声图像数据的形式，具体可为原始超声信号(RF)图像、原始信号经正交解调后的图像、对前述两类信号改变抽样率的图像或进行希尔伯特变换后形成的复值图像、对前述信号进行对数压缩后的灰度图像，以及对前述图像进行其它压缩、滤波处理以改善传输速率、显示质量等之后形成的图像等。
本实施例采用B超探头结合B超仪主机对靶区进行监控成像，采集并输出按时间排列的数字超声图像序列，要求输出图像的帧率保持稳定，即相邻两帧图像的时间间隔T0不变；所有图像具有相同的尺寸I0(横向)×J0(纵向)，单位为像素。需要说明的是，本实施例中根据待成像的活体类型、部位和靶区深度，按B超成像的相关技术规范选用不同类型的B超探头，设置B超主机的具体成像方式。此外需要说明的是，本发明不限定使用的B超探头类型和阵元数量，具体根据其外形可为线阵式探头、凸阵式探头、机械扇扫式探头、相控阵式探头、面阵式超声探头等，根据其阵元制作的材料可为压电陶瓷型、压电晶体型、压电复合材料型、压电式微机电型(PMUT)和电容式微机电型(CMUT)等。
2)互相关系数-时间曲线的计算
对数字超声图像序列以不同帧为参考帧计算得到一系列互相关系数—时间曲线；先对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γn(t)，具体地，数字超声图像序列共N帧，第a帧图像根据时间坐标记为aT0，a∈[1,2,…,N]；取数字超声图像序列的第1～M帧图像，其中，第1～M帧图像包括至少3个完整的呼吸过程；以第n帧为参考帧Bref，将N帧图像依次作为目标帧Bex，计算参考帧和目标帧的二维互相关系数γn，获得M条互相关系数—时间曲线γn(t)；其中，γn(t)为γn随时间变化的曲线，n∈[1,2,…,M]。值得说明的是，本发明适用的互相关计算的方法包括但不限于归一化互相关、绝对差求和、归一化协方差、平方差求和、非归一化互相关、Hybrid-sign相关、Polarity-coincidence相关或Meyr-Spies方法。以归一化互相关为例，令参考区域像素点矩阵为X，目标区域像素点矩阵为Y，则其互相关系数(皮尔逊相关系数)γ计算如下：
其中，COV(X,Y)表示X,Y的协方差，σX,σY表示X,Y的标准差。
3)选取“呼气”状态和“吸气”状态的第一参考帧
进一步地，根据互相关系数—时间曲线γn(t)选取“吸气”状态的第一参考帧；对γn(t)进行频谱分析得到频域的幅度谱Sn(f)，f为频率变量；记Sn(f)的最大值为An，若[A1,A2,…，AM]中的最大值为Am，m∈[1,2,…,M]，则将Sm(f)的峰值频率记为呼吸频率fres，将第m帧选为“吸气”状态第一参考帧；值得说明的是，本发明利用傅里叶变换或线性调频z变换对γn(t)进行频谱分析得到频域幅度谱Sn(f)。
进一步地，再根据“吸气”状态的第一参考帧对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γm(t)，具体地，以第m帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γm；并获得互相关系数—时间曲线γm(t)；其中，γm(t)为γm随时间t变化的曲线；对互相关系数—时间曲线γm(t)进行处理以识别“呼气”状态参考帧；具体地，对互相关系数—时间曲线γm(t)进行平滑降噪处理得到γ'm(t)；若γ'm(t)最小值所处的时间坐标为kT0，k为整数，则将第k帧选为“呼气”状态的第一参考帧。值得说明的是，可以采用Savitzky-Golay滤波，移动平均平滑处理，小波分析平滑，频域滤波平滑等方法对互相关系数—时间曲线进行平滑处理。
4)查找互相关系数-时间曲线的极值，划分“吸气”和“呼气”状态
进一步地，查找互相关系数—时间曲线中的极值；具体地，令γ'm(t)上第一个极大值位置的时间坐标t＝mT0为起点，先向正方向以步长Δtc交替检索极小值和极大值；0.5/fres<Δtc<1/fres，即以大于半个呼吸周期的步长进行极值检索。第一次在t＝mT0～mT0+Δtc的范围内找到极小值位于t＝pT0；第二次在t＝pT0～pT0+Δtc找到极大值位于t＝qT0；第三次在t＝pT0～pT0+Δtc寻找极小值，以此类推；然后在负方向以步长Δtc交替检索极小值和极大值，直至整条曲线被检索完毕。
进一步地，根据极值位置划分“吸气”状态和“呼气”状态；具体地，将所有极小值位置标记为“吸气”状态的周期边界，任意两个相邻极小值之间为“吸气”状态；所有极大值位置标记为“呼气”状态的周期边界，任意两个相邻极大值之间为“呼气”状态。值得说明的是，本发明通过划分“吸气”和“呼气”状态，进一步可以获得“呼气”和“吸气”两种状态下的热应变，从而使得结合后的累积位移、应变结果输出的时间间隔更小、更均匀，且由于其适合并行结构编程，应用的实时性更好。
5)分状态设置各周期阈值，初步筛选图像
“吸气”状态和“呼气”状态分别包括若干个周期；针对每个周期分别设置阈值，根据阈值对图像序列进行筛选得到剩余图像序列，其中，“吸气”状态对应的剩余图像序列为Iinh,0；“呼气”状态对应的剩余图像序列为Iexh,0；具体过程如下：
S100、在“吸气”状态的每个周期内，将γm(t)的所有值按从大到小排列，取周期内前C1％数据中的最小值为该周期的相关性阈值；
S200、在“呼气”状态的每个周期内，将γk(t)的所有值按从大到小排列，取前C1％数据中的最小值为该周期的相关性阈值；C1的取值范围均为1～100；γk为以第k帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γk，γk(t)为γk随时间t变化的曲线；
S300、令γ'm(t)上距离t＝T0最近的一个极值点位于t＝rT0，且为极小值，则该点之前标记为一个不完整的“吸气”状态，当该状态时长Thead≥0.9/fres时，根据步骤S100选取阈值；若0.9/fres>Thead≥0.5/fres，依据步骤S100选取前C2％数据中的最小值为阈值；若Thead<0.5/fres，采用相邻“吸气”状态的阈值。若距离t＝T0最近的一个极值点为极大值，则该点之前标记为一个不完整的“呼气”状态。当Thead≥0.9/fres时，根据步骤S200选取阈值；若0.9/fres>Thead≥0.5/fres，根据步骤S200选取前C2％数据中的最小值为阈值；若Thead<0.5/fres，采用相邻“呼气”状态的阈值。C1和C2的取值范围均为1～100，且C2>C1；
进一步地，考虑曲线γ'm(t)上距离t＝NT0最近的一个极值点，将该点之后标记为一个不完整状态，根据步骤S300设置阈值。
在“吸气”状态各周期对应的时间内，将γm(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iinh,0；在“呼气”状态各周期对应的时间内，将γk(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iexh,0。
6)筛选出高运动匹配的呼气、吸气状态图像(搜索处理图像序列)
对数字超声图像序列设置感兴趣区域，并对感兴趣区域进行带运动补偿的互相关检索得到处理图像序列，其中，“吸气”状态对应的处理图像序列为Iinh,1；“呼气”状态对应的处理图像序列为Iexh,1；具体过程如下：
(1)对数字超声图像序列预设统一的感兴趣区域ROI，ROI的尺寸为I*J，左下角顶点位置为i0，j0，使的所有图像中靶区均处于ROI内，I、J、i0、j0的单位均为像素；
(2)在“吸气”状态中，将Iinh,0中的第1帧记为参考帧Rref；将Rref的ROI区域作为一个新图像，置入处理图像序列Iinh,1；
(3)根据实际应变计算需求设置选帧范围上限Δτ，步骤4)中C1和C2的取值应使Iinh,0中位于Rref之后[Δt1,Δτ]范围内至少存在一个图像帧；将位于Rref之后[Δt1,Δτ]范围内的所有帧依次记为目标帧Rex；其中，Δt2>Δτ>Δt1，Δt1为选帧范围下限，Δt2为Δτ选择范围上限，且Δt1>0.1/fres，Δt2<1.5/fres；
(4)将每个Rex的ROI横向移动i1，纵向移动j1。其中，-I/2≤i1≤I/2，-J/2≤j1≤J/2，i1和j1为整数；对Rref和Rex的ROI计算二维互相关系数γR(i1,j1)。若γR(i1,j1)的最大值γR,max对应于i1＝i'1，j1＝j'1，则当前目标帧的运动补偿量为i'1(横)和j'1(纵)；
(5)对所有目标帧Rex计算γR,max和[i'1,j'1]；选择γR,max的最大值对应的目标帧记为新参考帧Rref，并将其经运动补偿后的ROI区域作为一个新图像，置于Iinh,1中已有帧之后；其中，经运动补偿后的ROI区域的横向范围为：i0+i'1～i0+i'1+I-1，纵向范围为j0+j'1～j0+j'1+J-1；值得说明的是，进行运动补偿的二维互相关检索时，需要明确，目标帧Rex的ROI横向和纵向移位量不必遍历-I/2≤i1≤I/2和-J/2≤j1≤J/2；该移位范围可根据应用场景减小，但在两个方向的正、负最大移位量至少应不小于组织的实际最大移位量。此外也可额外先对所有图像进行空间插值处理，以提高运动补偿的精度。
(6)重复步骤(3)～(5)，直至Iinh,0中的所有帧被处理完毕，获得处理图像序列Iinh,1；
依据步骤(2)～(6)，对“呼气”状态对应的图像序列处理得到处理图像序列为Iexh,1。值得说明的是，在进行运动状态匹配、位移计算、应变计算时都引入了运动补偿机制，使得计算结果更精确。
7)计算位移和应变
对Iinh,1和Iexh,1分别进行计算处理得到各自对应的位移图像序列和应变图像序列，该计算处理方法采用现有技术，例如发明创造名称为：一种基于低采样率B超图像计算热应变分布的方法(申请号：CN201910112033.8，公开日期：2019年4月12日)公开了计算处理得到各自对应的位移图像序列和应变图像序列的方法。本发明中Iinh,1对应位移图像序列Dinh和应变图像序列Sinh；Iexh,1对应位移图像序列Dexh和应变图像序列Sexh；
8)呼气、吸气状态的位移和应变合并
将Iinh,1和Iexh,1各自对应的位移图像序列和应变图像序列分别进行合并得到最终的组织位移图像序列和组织应变图像序列。具体过程如下：
(i)将Dinh中所有图像的时间坐标记为集合tinh，Dexh中所有图像的时间坐标记为集合texh；将Dinh和Dexh分别进行时域插值，使插值后的图像序列D'inh和D'exh中所有图像时间坐标的集合均为tinh和texh的合集tall；
(ii)将D'inh和D'exh中时间坐标相同的帧进行求和平均计算，将计算结果以同一时间坐标置入组织位移图像序列Dall；
依据步骤(i)和步骤(ii)对应变图像序列Sinh和Sexh进行处理得到组织应变图像序列Sall；
以tall为时序的组织位移图像序列Dall和组织应变图像序列Sall即为最终获得的组织位移分布图像和组织应变分布图像。
值得说明的是，本发明在完成“呼气”和“吸气”两种状态的识别后，对两种状态下图像数据的处理是独立的；因此，本发明适用于并行计算的编程结构，不但可以有效地提高数据处理的速度，也有利于提升结果输出的实时性。此外，本发明方法本质上是对时序排列的活体图像进行处理，因此需要明确，该方法同样可用于基于B超以外的其它物理手段获得的具有类似特征的二维活体图像序列。
进一步需要说明的是，现有技术中基于呼吸门控的方法进行运动抑制后计算组织位移和应变分布图像，会缺少运动状态、周期的判断，可能导致在“呼气”和“吸气”状态间的“过渡区间”进行图像选取，导致位移和应变计算误差增大，甚至发生错误；此外，在呼吸周期未知的情形下，难以明确在什么时间范围内寻找运动状态匹配的图像，可能导致运动抑制效果不佳，位移、应变的计算结果中包含较大累积误差，本发明基于两种呼吸状态计算组织位移和应变分布图像，即可克服上述缺点。具体地，本发明采用了呼吸分离及再合成的处理方法，且在图像运动状态匹配及位移和应变的计算中引入运动补偿机制，其结果可具有更短、更均匀的输出时间间隔和更好的精度。
实施例2
本实施例采用实施例1的方法，使用微波消融针对活猪的脂肪组织加热，计算其热应变，该方法具体步骤如下：
步骤一、如图2所示，本实施例的超声成像系统100采用采样率为40MHz的B超成像仪，成像探头200为128阵元、10.5MHz中心频率的线阵探头对活猪脂肪组织成像，成像深度4cm，微波消融针400的功率发射区域处于成像平面上；B超成像仪以平均50帧每秒的帧率，输出每帧尺寸为128像素*2048像素的RF图像，每次实验采用微波加热器500从微波开始加热起采集长度为20秒的图像序列。通过计算机300控制超声仪和微波机同时开始工作，将超声图像数据传输到计算机进行处理。
步骤二、对采集的RF图像序列，通过二维互相关计算和频谱分析，得到呼吸周期为2.67秒，“吸气”状态参考帧为m＝192帧，“呼气”状态参考帧为k＝633帧。
步骤三、在“吸气”状态参考帧对应的互相关系数曲线γm(t)上，从第192帧开始，交替检索极大值、极小值，对20s全部的图像划分“呼气”、“吸气”状态。如图3所示，图3(a)为“吸气”状态的周期划分结果，图3(b)为“呼气”状态的周期划分结果。
步骤四、取C1＝10，C2＝15，根据“呼气”、“吸气”状态下各周期的互相关系数分布，设置阈值并进行初步图像筛选。例如，第二个“吸气”周期的阈值设为0.87，第四个“呼气”周期的阈值设为0.63；
步骤五、根据加热区域在B超图像中的相对位置，设置感兴趣区域(ROI)的图像坐标范围为横向[i0,i0+I]＝[33,96]，纵向[j0,j0+J]＝[769,1794]；以下限Δt1＝5、上限Δτ＝30的选帧范围，带运动补偿地从步骤三中得到的图像序列中检索最优运动匹配ROI，得到的部分帧ROI图像如图4a和4b所示；图4a和4b表明：不同时刻下，ROI内的散斑位置和轮廓在“呼气”和“吸气”状态内分别具有很好的相似性，因此上述处理过程很好地克服了生理运动的影响。
步骤六、基于步骤五得到的高运动匹配的图像，在“呼气”、“吸气”两种状态下，分别计算位移和应变图像序列，如图5a，图5b分别为“吸气”状态和“呼气”状态的累积热应变分布，可以看到符合预期的热应变轮廓，且仅有较少的应变伪影。将两种状态的热应变通过时域插值后平均，其应变累积最大点在合并前后的应变-时间曲线如图6所示。本发明不需要使用额外的数据采集和运动补偿设备，使用B超图像可以独立地进行运动抑制以及位移和应变的计算，可以有效地减少诊断、治疗的成本。
在上文中结合具体的示例性实施例详细描述了本发明。但是，应当理解，可在不脱离由所附权利要求限定的本发明的范围的情况下进行各种修改和变型。详细的描述和附图应仅被认为是说明性的，而不是限制性的，如果存在任何这样的修改和变型，那么它们都将落入在此描述的本发明的范围内。此外，背景技术旨在为了说明本技术的研发现状和意义，并不旨在限制本发明或本申请和本发明的应用领域。
Claims (10)
1.一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，包括：
对活体的应变集中区域进行超声图像采集得到数字超声图像序列；
对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γn(t)，根据γn(t)选取“吸气”状态的第一参考帧；再根据“吸气”状态的第一参考帧对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γm(t)，对γm(t)进行处理得到互相关系数—时间曲线γ'm(t)，并根据γ'm(t)选取“呼气”状态的第一参考帧；
查找γ'm(t)的极值，并根据极值时间坐标划分“吸气”状态和“呼气”状态；其中，“吸气”状态和“呼气”状态分别包括若干个周期；
针对每个周期分别设置阈值，根据阈值对图像序列进行筛选得到剩余图像序列，其中，“吸气”状态对应的剩余图像序列为Iinh,0；“呼气”状态对应的剩余图像序列为Iexh,0；
对数字超声图像序列设置感兴趣区域，并对感兴趣区域进行带运动补偿的互相关检索得到处理图像序列，其中，“吸气”状态对应的处理图像序列为Iinh,1；“呼气”状态对应的处理图像序列为Iexh,1；
对Iinh,1和Iexh,1分别进行计算得到各自对应的位移图像序列和应变图像序列，其中，Iinh,1对应位移图像序列Dinh和应变图像序列Sinh；Iexh,1对应位移图像序列Dexh和应变图像序列Sexh；
将Iinh,1和Iexh,1各自对应的位移图像序列和应变图像序列分别进行合并得到最终的组织位移图像序列和组织应变图像序列。
2.根据权利要求1所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，对活体的应变区域进行图像采集的具体过程为：
对活体的靶区加热或施加外力，使得靶区产生应变集中区域；
对靶区的应变集中区域进行等时间间隔的连续成像，采集得到按时间排列的数字超声图像序列，相邻图像的采集时间间隔为T0。
3.根据权利要求1所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，对数字超声图像序列进行互相关计算得到互相关系数—时间曲线γn(t)的具体过程中：
数字超声图像序列共N帧，第a帧图像根据时间坐标记为aT0，a∈[1,2,…,N]；取数字超声图像序列的第1～M帧图像，其中，第1～M帧图像包括至少3个完整的呼吸过程；以第n帧为参考帧Bref，将所有N帧图像依次作为目标帧Bex，计算参考帧和目标帧的二维互相关系数γn，并获得M条互相关系数—时间曲线γn(t)；其中，γn(t)为γn随时间t变化的曲线，n∈[1,2,…,M]。
4.根据权利要求3所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，选取“呼气”状态和“吸气”状态的第一参考帧的具体过程为：
对γn(t)进行频谱分析得到频域的幅度谱Sn(f)，f为频率变量，记Sn(f)的最大值为An；若[A1,A2,…，AM]中的最大值为Am，m∈[1,2,…,M]，则将Sm(f)的峰值频率记为呼吸频率fres，并将第m帧选为“吸气”状态第一参考帧；
以第m帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γm；并获得互相关系数—时间曲线γm(t)；其中，γm(t)为γm随时间t变化的曲线；
对γm(t)进行平滑降噪处理得到互相关系数—时间曲线γ'm(t)；若γ'm(t)最小值所处的时间坐标t＝kT0，则将第k帧选为“呼气”状态的第一参考帧；其中，k为整数。
5.根据权利要求4所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，查找互相关系数—时间曲线中的极值的具体过程为：
极值包括极小值和极大值，以t＝mT0为起点和第一个极大值位置，先向正方向以步长Δtc交替检索极小值和极大值；然后在负方向以步长Δtc交替检索极小值和极大值，直至互相关系数—时间曲线γ'm(t)被检索完毕；其中，0.5/fres<Δtc<1/fres。
6.根据权利要求4所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，利用傅里叶变换或线性调频z变换对γn(t)进行频谱分析得到频域的幅度谱Sn(f)。
7.根据权利要求5所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，根据极值位置划分“吸气”状态和“呼气”状态的具体过程为：
将所有极小值位置标记为“吸气”状态的周期边界，任意两个相邻极小值之间为“吸气”状态；所有极大值位置标记为“呼气”状态的周期边界，任意两个相邻极大值之间为“呼气”状态。
8.根据权利要求4所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，针对每个周期设置阈值，并根据阈值对图像序列进行筛选得到剩余图像序列的具体过程为：
在“吸气”状态的每个周期内，将γm(t)的所有值按从大到小排列，取周期内前C1％数据中的最小值为该周期的相关性阈值；
在“呼气”状态的每个周期内，将γk(t)的所有值按从大到小排列，取周期内前C1％数据中的最小值为该周期的相关性阈值；C1的取值范围均为1～100；其中，γk为以第k帧为参考帧，以全部N帧图像依次为目标帧，计算得到二维互相关系数γk，γk(t)为γk随时间t变化的曲线；
在“吸气”状态各周期对应的时间内，将γm(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iinh,0；在“呼气”状态各周期对应的时间内，将γk(t)中小于阈值的点对应的目标帧从图像序列中移除，剩余图像序列记为Iexh,0。
9.根据权利要求4所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，计算得到处理图像序列的具体过程为：
(1)对数字超声图像序列预设统一的感兴趣区域ROI，ROI的尺寸为I*J，左下角顶点位置为i0，j0，使得所有图像中靶区均处于ROI内，I、J、i0、j0的单位均为像素；
(2)在“吸气”状态中，将Iinh,0中的第1帧记为参考帧Rref；将Rref的ROI区域作为一个新图像，置入处理图像序列Iinh,1；
(3)将位于Rref之后[Δt1,Δτ]范围内的所有帧依次记为目标帧Rex；其中，Δt1为选帧范围下限，Δτ为选帧范围上限，根据应变量级确定Δτ，要求Δt2>Δτ>Δt1，Δt2为Δτ选择范围上限，且Δt1>0.1/fres，Δt2<1.5/fres；
(4)将每个Rex的ROI横向移动i1，纵向移动j1；其中，-I/2≤i1≤I/2，-J/2≤j1≤J/2，i1和j1为整数；对Rref和Rex的ROI区域计算二维互相关系数γR(i1,j1)；若γR(i1,j1)的最大值γR,max对应于i1＝i'1，j1＝j'1，则当前目标帧的运动补偿量为i'1和j'1；
(5)对所有目标帧Rex计算γR,max和[i'1,j'1]；选择γR,max的最大值对应的目标帧记为新参考帧，并将新参考帧经运动补偿后的ROI区域作为一个新图像，置于Iinh,1中已有帧之后；其中，经运动补偿后的ROI区域的横向范围为：i0+i'1～i0+i'1+I-1，纵向范围为j0+j'1～j0+j'1+J-1；
(6)重复步骤(3)～(5)，直至Iinh,0中的所有帧被处理完毕，获得处理图像序列Iinh,1；
依据步骤(2)～(6)，对“呼气”状态对应的图像序列处理得到处理图像序列为Iexh,1。
10.根据权利要求1～9任一项所述的一种基于活体超声图像的呼吸分离式应变成像方法，其特征在于，合并得到最终的组织位移图像序列和组织应变图像序列的具体过程为：
(i)将Dinh中所有图像的时间坐标记为集合tinh，Dexh中所有图像的时间坐标记为集合texh；将Dinh和Dexh分别进行时域插值，使插值后的图像序列D'inh和D'exh中所有图像时间坐标的集合均为tinh和texh的合集tall；
(ii)将D'inh和D'exh中时间坐标相同的帧进行求和平均计算，将计算结果以同一时间坐标置入组织位移图像序列Dall；
依据步骤(i)和步骤(ii)对应变图像序列Sinh和Sexh进行处理得到组织应变图像序列Sall；
以tall为时序的组织位移图像序列Dall和组织应变图像序列Sall即为最终获得的组织位移分布图像和组织应变分布图像。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010731441.4A CN111681740B (zh) | 2020-07-27 | 2020-07-27 | 一种基于活体超声图像的呼吸分离式应变成像方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010731441.4A CN111681740B (zh) | 2020-07-27 | 2020-07-27 | 一种基于活体超声图像的呼吸分离式应变成像方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111681740A true CN111681740A (zh) | 2020-09-18 |
| CN111681740B CN111681740B (zh) | 2023-07-18 |
Family
ID=72438479
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010731441.4A Active CN111681740B (zh) | 2020-07-27 | 2020-07-27 | 一种基于活体超声图像的呼吸分离式应变成像方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111681740B (zh) |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP1772825A1 (en) * | 2005-10-06 | 2007-04-11 | Esaote S.p.A. | Method for registering images of a sequence of images, particularly ultrasound diagnostic images |
| US20080095417A1 (en) * | 2006-10-23 | 2008-04-24 | Gianni Pedrizzetti | Method for registering images of a sequence of images, particularly ultrasound diagnostic images |
| CN105869144A (zh) * | 2016-03-21 | 2016-08-17 | 常州大学 | 一种基于深度图像数据的非接触式呼吸监测方法 |
| CN106056589A (zh) * | 2016-05-24 | 2016-10-26 | 西安交通大学 | 一种呼吸运动补偿的超声造影灌注参量成像方法 |
-
2020
- 2020-07-27 CN CN202010731441.4A patent/CN111681740B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP1772825A1 (en) * | 2005-10-06 | 2007-04-11 | Esaote S.p.A. | Method for registering images of a sequence of images, particularly ultrasound diagnostic images |
| US20080095417A1 (en) * | 2006-10-23 | 2008-04-24 | Gianni Pedrizzetti | Method for registering images of a sequence of images, particularly ultrasound diagnostic images |
| CN105869144A (zh) * | 2016-03-21 | 2016-08-17 | 常州大学 | 一种基于深度图像数据的非接触式呼吸监测方法 |
| CN106056589A (zh) * | 2016-05-24 | 2016-10-26 | 西安交通大学 | 一种呼吸运动补偿的超声造影灌注参量成像方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111681740B (zh) | 2023-07-18 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP5031758B2 (ja) | 医用画像診断装置 | |
| CN107569256B (zh) | 基于热膨胀和门控算法测量生物组织温度变化的超声方法 | |
| Maurice et al. | Noninvasive vascular elastography for carotid artery characterization on subjects without previous history of atherosclerosis | |
| JP4857418B2 (ja) | 検体のモーションアーチファクト補正の方法およびシステム | |
| WO2013053000A1 (en) | Heart imaging method | |
| Jiang et al. | A novel performance descriptor for ultrasonic strain imaging: A preliminary study | |
| Soepriatna et al. | Three-dimensional myocardial strain correlates with murine left ventricular remodelling severity post-infarction | |
| CN111225617A (zh) | 超声成像系统和方法 | |
| Liu et al. | Asthma pattern identification via continuous diaphragm motion monitoring | |
| Hussain et al. | Direct and gradient-based average strain estimation by using weighted nearest neighbor cross-correlation peaks | |
| Yue et al. | Speckle tracking in intracardiac echocardiography for the assessment of myocardial deformation | |
| Mukherjee et al. | Computing myocardial motion in 4-dimensional echocardiography | |
| CN111272305B (zh) | 一种基于非线性热膨胀评估温度变化的超声方法及系统 | |
| Kolen et al. | Characterization of cardiovascular liver motion for the eventual application of elasticity imaging to the liver in vivo | |
| US10346989B2 (en) | Method and system for calculating a displacement of an object of interest | |
| CN111681740B (zh) | 一种基于活体超声图像的呼吸分离式应变成像方法 | |
| CN110931130B (zh) | 一种基于b超信号评估呼吸和心动周期的方法 | |
| Gatta et al. | Robust image-based ivus pullbacks gating | |
| JP2022168851A (ja) | ３ｄモデル構築のための自動フレーム選択 | |
| CN102525568A (zh) | 一种弹性减影成像方法 | |
| CN114240815A (zh) | 一种基于活体超声图像的多线程应变成像方法及装置 | |
| Isguder et al. | Manifold learning for image-based gating of intravascular ultrasound (IVUS) pullback sequences | |
| US20240046464A1 (en) | Intramyocardial tissue displacement and motion measurement and strain analysis from mri cine images using dense deep learning | |
| Szilágyi et al. | Volumetric analysis of the heart using echocardiography | |
| Shan et al. | Dynamic estimation of myocardial deformation using ultrasound RF-data: A preliminary study |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |