CN101190132B - 超声成像的预处理方法与装置 - Google Patents
超声成像的预处理方法与装置 Download PDFInfo
- Publication number
- CN101190132B CN101190132B CN2006101469892A CN200610146989A CN101190132B CN 101190132 B CN101190132 B CN 101190132B CN 2006101469892 A CN2006101469892 A CN 2006101469892A CN 200610146989 A CN200610146989 A CN 200610146989A CN 101190132 B CN101190132 B CN 101190132B
- Authority
- CN
- China
- Prior art keywords
- voxel
- line
- volume data
- voxels
- ultra sonic
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 49
- 238000002203 pretreatment Methods 0.000 title abstract description 4
- 238000000034 method Methods 0.000 claims abstract description 59
- 239000013598 vector Substances 0.000 claims abstract description 36
- 238000004364 calculation method Methods 0.000 claims abstract description 11
- 238000006243 chemical reaction Methods 0.000 claims abstract description 7
- 238000002604 ultrasonography Methods 0.000 claims description 23
- 230000003321 amplification Effects 0.000 claims description 13
- 238000003199 nucleic acid amplification method Methods 0.000 claims description 13
- 238000012545 processing Methods 0.000 claims description 9
- 238000005070 sampling Methods 0.000 claims description 9
- 238000011946 reduction process Methods 0.000 claims description 4
- 238000009499 grossing Methods 0.000 abstract description 10
- 238000010276 construction Methods 0.000 abstract 1
- 230000000694 effects Effects 0.000 description 10
- 238000001914 filtration Methods 0.000 description 6
- 239000000523 sample Substances 0.000 description 5
- 238000013507 mapping Methods 0.000 description 3
- 239000002245 particle Substances 0.000 description 3
- 238000013461 design Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 230000000007 visual effect Effects 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 1
- 238000002059 diagnostic imaging Methods 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 239000008187 granular material Substances 0.000 description 1
- 230000005865 ionizing radiation Effects 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000000750 progressive effect Effects 0.000 description 1
- 238000013139 quantization Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000012285 ultrasound imaging Methods 0.000 description 1
- 238000012800 visualization Methods 0.000 description 1
- 238000007794 visualization technique Methods 0.000 description 1
Images
Classifications
-
- G06T5/70—
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T15/00—3D [Three Dimensional] image rendering
- G06T15/08—Volume rendering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/20—Image enhancement or restoration by the use of local operators
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S15/00—Systems using the reflection or reradiation of acoustic waves, e.g. sonar systems
- G01S15/88—Sonar systems specially adapted for specific applications
- G01S15/89—Sonar systems specially adapted for specific applications for mapping or imaging
- G01S15/8906—Short-range imaging systems; Acoustic microscope systems using pulse-echo techniques
- G01S15/8993—Three dimensional imaging systems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2200/00—Indexing scheme for image data processing or generation, in general
- G06T2200/04—Indexing scheme for image data processing or generation, in general involving 3D image data
-
- 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/20172—Image enhancement details
- G06T2207/20192—Edge enhancement; Edge preservation
Abstract
本发明公开了一种超声成像的预处理方法和模块,所述方法包括:计算步骤,用于构建一个多值矢量场;以及平滑步骤,用于对整个体数据进行平滑处理。还可包括判断步骤,用于根据体素的梯度判断它是否位于物体的表面附近,然后只对表面附近的体素进行平滑处理,对其他区域的体素只进行均值滤波处理;以及缩、放步骤,用于对扫描转换前的体数据进行缩小处理,并用于将平滑处理之后的体数据放大到原始大小。所述模块包括计算单元和平滑单元,还可包括判断单元以及缩小和放大单元。按照本发明的超声成像的预处理方法,通过计算沿表面分布的多个结点的均值,可以有效地滤除斑点噪声,达到平滑的目的。因此,这种算法可以获得平滑和细节的两全。
Description
技术领域
本发明涉及一种图像信号的处理方法与装置,特别是涉及一种超声成像的预处理方法与可以实现该方法模块。
背景技术
传统的医疗影像设备只能提供人体内部的二维图像。医生们只能凭经验由多幅二维图像去估计病灶的大小及形状,以此“构思”病灶与其周围组织的三维几何关系,这给治疗带来了困难。而三维可视化技术可以由一系列二维图像重构出三维形体,并在终端显示出来。因此不仅能得到有关成像物体直观、形象的整体概念,而且还可保存许多重要的三维信息。由于超声成像较CT、MRI具有无创、无电离辐射以及操作灵活等明显优势,因此超声三维成像势必会在医学临床上得到广泛的应用,开展超声领域中的三维可视化研究显得十分必要。
目前有两种获取三维超声体数据的方法:一种是利用现有的二维超声诊断设备,结合某种定位机械,获取一系列空间位置已知的二维组织超声图像,进而以离线方式获取三维体数据;另一种是利用二维面阵探头发射金字塔形体积超声束,从而获得实时的三维体数据。
三维体数据的可视化已有比较成熟的方法,但是超声图像特有的斑点噪声成为影响成像效果、进而制约其实际应用的重要因素。三维超声体数据是由一系列超声图像组成的,而在超声图像形成的过程中,距离小于波长的反射点的反射回声相互干扰,由此便形成了斑点噪声。它的存在会掩盖图像特征,干扰观察和分析。若不对其进行处理,而直接以这样的体数据进行三维成像,其视觉效果是无法令人接受的。预处理的目的就是对超声体数据进行滤波或平滑,从而抑制噪声,以便在此基础上进行可视化,获得良好的成像质量,从而发挥三维超声成像的优势。
专利文献US005409007A公开了一种用于抑制超声成像中的斑点伪像的滤波器(Filter to Reduce Speckle Artifact in Ultrasound Imaging)。该专利提到的预处理方法简述如下:对超声图像(扫描转换前,即DSC前)的每一个像素,取以之为中心的5点菱形模板,计算这5个像素的灰度中值,以此取代中心像素的原灰度值,如图1和图2所示。
上述专利技术的缺点是不能直接用于三维,即使将其扩展到三维并对模板尺寸进行调整,要么噪声仍很突出,要么成像结果变得模糊,始终无法兼顾平滑和细节。因此,不能满足实际应用需要。
发明内容
本发明的目的就是要克服现有技术中的这些缺点,提供一种在不损失细节的前提下进行超声体数据平滑的方法。
按照本发明的第一方面,提供一种超声成像预处理方法,包括:计算步骤,用于计算在以每个体素为中心的半径为R的邻域内、以当前体素为中心的线段的均匀程度,以便从中选取N条最均匀的方向线作为2N个切线矢量,并且从中选取一条最不均匀的方向线作为法线;当对所有体素计算完毕后,即可得到一个多值矢量场;以及平滑步骤,对每个体素用其法线所代表的最大不均匀度作为梯度,先对整个梯度数据以半径为r的邻域进行一次均值滤波;若当前体素均值滤波后的梯度小于预定阈值,则计算其邻域的灰度均值并赋给所述当前体素;否则,沿着多值矢量场在所述当前体素位置的2N个方向上分别行进,其中每一步后又产生2N个分支,由此生成一个树形结构,计算树上的所有体素的灰度均值,赋给所述当前体素;所有体素计算完毕后,即完成了整个体数据的平滑。其中2N个切线矢量为N条方向线的两端的指向。其中R为1,r可取1或2,对于三维超声图像N为2,对于二维超声图像N为1,所述预定阈值可取分割背景与物体的经验值。
按照本发明第一方面的超声成像预处理方法,其中所述方向线的不均匀度V按照下式来计算:
式中G0是当前体素的灰度值,Gj是方向线上第j个体素的灰度值,其中j取负数和正数分别代表方向线的两端的指向。
按照本发明第一方面的超声成像预处理方法,可选地是还包括判断步骤,用于根据体素的梯度判断它是否位于物体的表面附近,然后只对表面附近的体素通过所述平滑步骤进行处理,对其他区域的体素只进行均值滤波处理。
按照本发明第一方面的超声成像预处理方法,可选地是还包括:缩小步骤,用于对扫描转换前的体数据进行缩小处理,然后再通过所述计算步骤和所述平滑步骤进行处理;以及放大步骤,用于将平滑处理之后的体数据放大到原始大小。优选地是,所述缩小步骤为:对于三维超声图像,可对体数据做X·Y·Z的均值滤波,再做X·Y·Z的采样,其中X、Y、Z为三个方向的缩小倍数;对于二维超声图像,可对体数据做X·Y的均值滤波,再做X·Y的采样,其中X和Y为两个方向的缩小倍数。优选地是,所述放大步骤为:对于三维超声图像,将体数据中的每个体素直接复制成X·Y·Z个,再做X·Y·Z的均值滤波,其中X、Y、Z为三个方向的放大倍数;对于二维超声图像,将体数据中的每个体素直接复制成X·Y个,再做X·Y的均值滤波,其中X、Y为两个方向的放大倍数。其中缩小和放大时做均值滤波的窗口大小可以取不同于X·Y·Z或X·Y的值。
按照本发明第一方面的超声成像预处理方法,其中对于三维超声图像,所述邻域为以每个体素为中心、以R为半径的立方体,其中R为立方体中心到表面的距离;所述方向线为以当前体素为中心并经过半径为R的立方体邻域表面的两个体素的线段,可取13条方向线或9条方向线。对于二维超声图像,所述邻域为以每个体素为中心、以R为半径的正方形,其中R为正方形中心到其边的距离;所述方向线为以当前体素为中心并经过半径为R的正方形邻边的两个体素的线段,共有4条方向线。
按照本发明的第二方面,提供一种超声成像预处理模块,包括:计算单元,用于计算在以每个体素为中心的半径为R的邻域内、以当前体素为中心的线段的均匀程度,以便从中选取N条最均匀的方向线作为2N个切线矢量,并且从中选取一条最不均匀的方向线作为法线;当对所有体素计算完毕后,即可得到一个多值矢量场;以及平滑单元,对每个体素用其法线所代表的最大不均匀度作为梯度,先对整个梯度数据以半径为r的邻域进行一次均值滤波;若当前体素均值滤波后的梯度小于预定阈值,则计算其邻域的灰度均值并赋给所述当前体素;否则,沿着多值矢量场在所述当前体素位置的2N个方向上分别行进,其中每一步后又产生2N个分支,由此生成一个树形结构,计算树上的所有体素的灰度均值,赋给所述当前体素;所有体素计算完毕后,即完成了整个体数据的平滑。其中2N个切线矢量为N条方向线的两端的指向。其中R为1,r可取1或2,对于三维超声图像N为2,对于二维超声图像N为1,所述预定阈值可取分割背景与物体的经验值。
按照本发明第二方面的超声成像预处理方法,其中所述方向线的不均匀度V按照下式来计算:
式中G0是当前体素的灰度值,Gj是方向线上第j个体素的灰度值,其中j取负数和正数分别代表方向线的两端的指向。
按照本发明第二方面的超声成像预处理模块,可选地是还包括判断单元,用于根据体素的梯度判断它是否位于物体的表面附近,然后只对表面附近的体素通过所述平滑单元进行处理,对其他区域的体素只进行均值滤波处理。
按照本发明第二方面的超声成像预处理模块,可选地是还包括:缩小单元,用于对扫描转换前的体数据进行缩小处理,然后再由所述计算单元和所述平滑单元对缩小后的体数据进行处理;以及放大单元,用于将平滑处理之后的体数据放大到原始大小。优选地是,所述缩小单元执行:对于三维超声图像,可对体数据做X·Y·Z的均值滤波,再做X·Y·Z的采样,其中X、Y、Z为三个方向的缩小倍数;对于二维超声图像,可对体数据做X·Y的均值滤波,再做X·Y的采样,其中X和Y为两个方向的缩小倍数。优选地是,所述放大单元执行:对于三维超声图像,将体数据中的每个体素直接复制成X·Y·Z个,再做X·Y·Z的均值滤波,其中X、Y、Z为三个方向的放大倍数;对于二维超声图像,将体数据中的每个体素直接复制成X·Y个,再做X·Y的均值滤波,其中X、Y为两个方向的放大倍数。其中缩小和放大时做均值滤波的窗口大小可以取不同于X·Y·Z或X·Y的值。
按照本发明第二方面的超声成像预处理方法,其中对于三维超声图像,所述邻域为以每个体素为中心、以R为半径的立方体,其中R为立方体中心到表面的距离;所述方向线为以当前体素为中心并经过半径为R的立方体邻域表面的两个体素的线段,可取13条方向线或9条方向线。对于二维超声图像,所述邻域为以每个体素为中心、以R为半径的正方形,其中R为正方形中心到其边的距离;所述方向线为以当前体素为中心并经过半径为R的正方形邻边的两个体素的线段,共有4条方向线。
采用本发明技术方案的超声成像预处理方法,在体数据内的物体表面处,切线树的树枝将尽可能沿表面的方向分布,而不会跨越表面或者脱离表面,因此可以有效地保护边缘。同时,通过计算沿表面分布的多个结点的均值,可以有效地滤除斑点噪声,达到平滑的目的。所以,这种算法可以获得平滑和细节的两全。
附图说明
图1是说明由射束成形器在矢量扫描数据的采集与显示中使用的格式的图(现有技术);
图2是说明在极坐标域对矢量扫描数据执行中值滤波的概念的图(现有技术);
图3是按照本发明对三维超声图像预处理时所取的13条方向线的示意图,其中R=1,白点是当前体素,黑点是方向线的端点;
图4是按照本发明对三维超声图像预处理时的N条切线和1条法线的示例,其中R=1、N=2;
图5是按照本发明的切线树的示例,其中2N=4、S=3,并且用不同的线段长度表示树的层次关系;
图6是按照本发明对二维超声图像预处理时所取的4条方向线的示意图;
图7是按照本发明的超声图像预处理模块的结构框图;以及
图8是本发明在实际应用中的示意图,描述了从探头扫描到显示三维图像的总流程。
具体实施方式
1.超声图像预处理方法
下面以三维超声图像为例,来具体说明本发明的超声成像预处理方法。
概括地讲,按照本发明的超声成像预处理方法,首先对体数据在扫描转换(即DSC)前进行缩小操作,然后执行算法的核心部分,最后放大到原来的尺寸。在算法的核心部分中,首先对于每个体素,在以其为中心的“半径”为R的立方体邻域(这里半径的定义为立方体中心到表面的距离,那么半径为R意味着边长为2R+1)内,计算以当前体素为中心的13条线段的均匀程度,并从中选取N条最均匀的线段,由此产生2N个最佳方向;所有体素计算完毕后,即得到一个多值矢量场。然后对于每个体素,若它具有较小的梯度幅值,则计算它的立方体邻域的灰度均值,赋给当前体素;否则,将沿着多值矢量场在该位置的2N个方向分别行进S步,在每一步后又将产生2N个分支,由此生成一个树形结构。计算树上的所有体素的灰度均值,赋给当前体素。所有体素计算完毕后,即完成了整个体数据的平滑。
为了提高算法的速度和稳定性,将在算法核心部分的前后对体数据采取缩小和放大处理。一般来说,三个方向的缩放倍数X、Y、Z均为整数。放缩倍数XYZ与扫描的空间状态有关,比如,若X方向对应的实际位置排列较密,则可取较大的X值;这也与算法对速度的要求有关,要求越高则可以采用越大的XYZ值。缩小时,可对体数据做X·Y·Z的均值滤波,再做X·Y·Z的采样来实现;这种方法无论速度还是效果都优于后向映射结合三线性插值的方法。放大时,可将体数据中的每个体素直接复制成X·Y·Z个,再做X·Y·Z的均值滤波来实现;这种方法的效果与后向映射结合三线性插值的方法相同,而速度优于该方法。
事实上,缩小和放大时做均值滤波的窗口大小不一定必须取X·Y·Z,也可以根据经验、所处理体数据的特性(如数据的尺寸、扫描的深度、目标的亮度、噪声的特性等)或实际需要而调整为不同的常数或者自适应的数值,例如取其他的X’·Y’·Z’和X”·Y”·Z”值,这样可以灵活地实现不同的效果。一般应有X’>=X,Y’>=Y,Z’>=Z,X”>=X,Y”>=Y,Z”>=Z。这6个不等式是独立的,等号成立时,无特殊效果,大于号成立时,有一个额外的平滑效果。这也与数据中噪声颗粒的平均大小有关,颗粒普遍较大时可取大于号成立。例如,一种可能的组合是:X=X’=X”=2,Y=Y’=Y”=2,Z=1,Z’=Z”=2。另外,缩、放方法还可采用普通的后向映射结合三线性插值的方法。
在算法的核心部分中,首先需要构建一个多值矢量场,方法如下:
如图3所示,对于当前体素,可作以其为中心的13条线段,每条线段都经过半径为1的立方体邻域表面的两个体素,并延此方向依次穿过半径从1到R的立方体邻域表面。这样会有13条线段,称这13条线段为“方向线”(如果不考虑立方体的4个对角方向,总的方向线数可由13条变为9条)。连同当前体素在内,每条方向线都包括2R+1个体素。这里R=1,白点是当前体素,黑点是方向线的端点,见图3。
对于当前体素,需刻划每条方向线的均匀程度,以便比较。这里可以采用传统的方差计算公式或梯度公式,然而本发明优选地采用如下公式进行计算:
式中G0是当前体素的灰度值,Gj是方向线上第j个体素的灰度值,其中j取负数和正数分别代表方向线的两端的指向。V表示当前方向线的“不均匀度”,其意义与方差类似。式(1)与方差计算相比,更突出了当前体素的作用。对于在本发明中的应用,它的准确程度并不亚于普通的方差计算公式,而且计算量要小得多。
随后,在当前体素的13条方向线中取出N个最均匀的和1个最不均匀的方向线,即V值最小的N条方向线和V值最大的1条方向线。可以近似认为,它们构成了体数据在当前体素位置的N条切线和1条法线。如图4所示,这里R=1,N=2,以P线表示切线,Q线表示法线。理论上,法线应与切线垂直,但由于体数据的量化效应以及不均匀度公式的设计,它们的方向往往并不严格垂直,这点不必深究。
现在暂时只考虑切线,后面再介绍法线的应用。由于每条方向线的两端都指向两个方向,所以实际获得了2N个“切线矢量”。若只考虑矢量的方向而不考虑长度,则每个切线矢量可以用它上面的最接近当前体素的那个体素(即该切线矢量与半径为1的立方体邻域表面的交点)来代表,称那个体素为当前体素的当前切线矢量的“定位点”。
对于当前体素靠近体数据边界的情况,为保证方向线上的体素完全落在体数据范围内,可以采取减小方向线跨度等方法。当所有体素的所有切线矢量计算完毕后,即得到一个多值矢量场,这时开始执行平滑操作,方法是:对当前体素,根据多值矢量场在那里的取值,可沿着它的2N个切线方向分别行进1步,到达2N个定位点。随后,对每个定位点,又可以根据多值矢量场在那里的取值,行进到它的2N个定位点。以此类推,将形成一个2N个叉的树形结构。
定义以下术语:称这个树形结构为“切线树”,记之为T;称当前处理的体素O为T的“根结点”;T经过的其他体素为T的“分支结点”。若从O开始经过d步行进可以到达T的某个结点C,称C的“层次”为d;称从O到C依次经过的d+1个结点为C的“路径”。若规定路径的最大层次为S,则称S为T的“深度”,即有0≤d≤S。显然,层次为0的唯一结点即为O,层次为d的结点共有(2N)d个。因此,T的总结点数为:
图5给出了一个切线树的例子,这里2N=4,S=3。例中N和S的取值正是实验给出的经验值,事实上,若N和S越小,算法速度越快,但不够平滑;反之则越慢而越平滑,甚至因过度平滑而损失细节。一般N可固定为2,而S的取值跟算法对速度的要求,以及数据中噪声颗粒的平均大小有关,越大则需要越大的S。图中,从根结点O开始,行进到它的4个定位点,从而扩展出O1、O2、O3、O4等4个层次为1的结点;再从O1开始,行进到它的4个定位点,扩展出O11、O12、O13、O14等4个层次为2的结点;又从O11开始,行进到它的4个定位点,扩展出O111、O112、O113、O114等4个层次为3的结点;由于O111等的层次已经达到设定的深度S,因此不再对它们进行扩展。以此类推,生成整个树T,共85个结点。这里指出,图中不同的线段长度只是便于表示树的层次关系,并不代表每次行进的距离。而且,在沿树前进的过程中,可能会折返到已经经过的位置。但是,这种情况不必特殊处理,即使数次经过同一个位置,仍在逻辑上将它们当作树上不同的结点看待。
按照这样的方法,构建树T并求其所有结点的灰度平均值,赋给结点O,即完成了当前体素的平滑。平滑后的灰度为:
式(3)中的M是树T的总结点数,一般可按照式(2)计算。对于当前体素靠近体数据边界的情况,树T将不能充分扩展,其结点数M将小于式(2)计算的结果,但式(3)仍然成立。因此,可计算实际扩展出的所有结点的灰度平均值,赋给结点O。
如此遍历整个体数据,则完成了体数据的平滑。称这种算法为“切线树平滑算法”。它可以很好地保护边缘,在不损失细节的前提下滤除噪声。但我们注意到,切线树平滑算法虽然具有理想的平滑效果,但复杂度稍高。若对所有体素都执行该算法,耗时将稍长。因此对全局方案进行优化,方法如下:
首先,根据当前体素的梯度幅值(以下简称梯度)判断它是否位于物体的表面附近。那么可以只在表面附近(即梯度较大的区域)进行切线树操作,而在其他区域(物体内部或背景处,即梯度较小的区域)只简单地进行均值滤波操作。因为均值滤波的速度是很快的,而位于物体表面附近的体素在体数据中所占的比例通常很小,因此这样可以显著地提高算法的速度。又由于只有物体表面附近的体素才具有较多的特征,也只有它们才对成像结果起主导作用,因此平滑的效果几乎不受影响。当然,小梯度情况下仍可执行切线树算法,此时可以减小当前切线树的深度S,也可以保持当前切线树的深度S不变。
考虑到前面构建多值矢量场时已经获取了每个体素的法线,而且法线的V值代表当前体素所有方向线的最大不均匀度,那么可以认为它近似地代表当前体素的梯度。对于本发明,这种方法判断表面的准确程度并不亚于普通的梯度计算公式,而且可以利用前面已经得到的结果,从而大大节省了计算量。
为了避免噪声和尖锐的边缘对梯度计算的影响,在求出体数据内所有体素的梯度之后,先对整个梯度数据以半径为r的立方体邻域进行一次均值滤波。这时,可以经验地设定一个阈值H。若当前体素的梯度(滤波后)小于H,则以体数据在当前体素的半径为R’的立方体邻域的灰度均值作为当前体素平滑的结果。根据经验,r可取1或2(r的值跟数据中噪声颗粒的平均大小有关,越大则需越大的r,但一般不超过2),R’可取方向线的跨度R(即邻域的半径R),H可取分割背景与物体的经验阈值。若当前体素的梯度(滤波后)不小于H,则执行前面介绍的切线树平滑算法。
对体数据内的所有体素进行上述的判断和分别处理,则完成了平滑操作。这时,完成了体数据的整个预处理。
当处理对象是二维超声图像时,所取的邻域变为正方形,备选的13条方向线将变为4条,而最佳方向线条数N应变为1,其他变更就是简单的降维,如图6所示。
2.超声图像预处理模块
图7为按照本发明的超声图像预处理模块,包括依次连接的缩小单元700、计算单元702、判断单元704、平滑单元706和放大单元708,以及一个均值滤波处理的分支。对原始的体数据,经过缩小单元702,得到缩小的体数据之后,由计算单元704计算以每个体素为中心的半径为R的邻域内、以当前体素为中心的线段的均匀程度,以从中选取N条最均匀的方向线作为2N个切线矢量,并且从中选取一条最不均匀的方向线作为法线;当对所有体素计算完毕后,即可得到一个多值矢量场。再经由判断单元700判断体素是否位于物体的表面附近,然后只对表面附近的体素经过平滑单元706进行平滑处理,而对其他区域的体素只进行均值滤波处理。所有体素处理完毕后,即完成了缩小的体数据的平滑。最后,放大单元708将该结果放大到原始大小,即得到平滑的体数据。
图8是本发明在实际应用中的示意图,描述了从探头扫描到显示三维图像的总流程。在前端,一组通过延迟聚焦的脉冲通过发射电路发送到探头,以发射出超声波;在经过一段延时后,探头接收到人体反射回的超声波,通过压电效应转化为电信号,经过前端放大,A/D转换模拟信号为数字信号,从而获取原始的体数据,并保存到内存中。在后端,原始体数据经过预处理模块,得到平滑后的体数据;再经过三维成像模块,就可以通过显示设备进行显示。
Claims (20)
1.一种超声成像预处理方法,其特征在于,包括以下步骤:
计算步骤,用于计算在以每个体素为中心的半径为R的邻域内、以当前体素为中心的线段的均匀程度,以便从中选取N条最均匀的方向线作为2N个切线矢量,并且从中选取一条最不均匀的方向线作为法线;当对所有体素计算完毕后,即可得到一个多值矢量场;以及
平滑步骤,对每个体素用其法线所代表的最大不均匀度作为梯度,先对整个梯度数据以半径为r的邻域进行一次均值滤波;若当前体素均值滤波后的梯度小于预定阈值,则计算其邻域的灰度均值并赋给所述当前体素;否则,沿着多值矢量场在所述当前体素位置的2N个切线矢量方向上分别行进,其中每一步又产生2N个分支,由此生成一个树形结构,计算树上的所有体素的灰度均值,赋给所述当前体素;所有体素计算完毕后,即完成了整个体数据的平滑。
2.如权利要求1所述的超声成像预处理方法,其特征在于,还包括:
判断步骤,用于根据体素的梯度判断它是否位于物体的表面附近,然后只对表面附近的体素通过所述平滑步骤进行处理,对其他区域的体素只进行均值滤波处理。
3.如权利要求1所述的超声成像预处理方法,其特征在于,还包括:
缩小步骤,用于对扫描转换前的体数据进行缩小处理,然后再通过所述计算步骤和所述平滑步骤进行处理;以及
放大步骤,用于将通过平滑步骤进行处理之后的体数据放大到原始大小。
4.如权利要求3所述的超声成像预处理方法,其特征在于:
其中所述缩小步骤为:对于三维超声图像,可对体数据做X·Y·Z的均值滤波,再做X·Y·Z的采样,其中X、Y、Z为三个方向的缩小倍数;对于二维超声图像,可对体数据做X·Y的均值滤波,再做X·Y的采样,其中X和Y为两个方向的缩小倍数;
其中所述放大步骤为:对于三维超声图像,将体数据中的每个体素直接复制成X·Y·Z个,再做X·Y·Z的均值滤波,其中X、Y、Z为三个方向的放大倍数;对于二维超声图像,将体数据中的每个体素直接复制成X·Y个,再做X·Y的均值滤波,其中X、Y为两个方向的放大倍数。
5.如权利要求4所述的超声成像预处理方法,其特征在于:
其中缩小和放大时做均值滤波的窗口大小可以取不同于X·Y·Z或X·Y的值。
6.如权利要求1至5中的任何一项所述的超声成像预处理方法,其特征在于:
其中2N个切线矢量为N条方向线的两端的指向。
7.如权利要求1至5中的任何一项所述的超声成像预处理方法,其特征在于:
R可取1;r可取1或2;对于三维超声图像N为2,对于二维超声图像N为1;所述预定阈值可取分割背景与物体的经验值。
8.如权利要求1至5中的任何一项所述的超声成像预处理方法,其特征在于,所述方向线的不均匀度V按照下式来计算:
式中G0是当前体素的灰度值,Gj是方向线上第j个元素的灰度值,其中j取负数和正数分别代表方向线的两端的指向。
9.如权利要求1至5中的任何一项所述的超声成像预处理方法,其特征在于:
对于三维超声图像,所述半径为R的邻域为以每个体素为中心、以R为半径的立方体,其中R为立方体中心到表面的距离;对于二维超声图像,所述半径为R的邻域为以每个体素为中心、以R为半径的正方形,其中R为正方形中心到其边的距离。
10.如权利要求9所述的超声成像预处理方法,其特征在于:
对于三维超声图像,所述方向线为以当前体素为中心并经过半径为R的立方体邻域表面的两个体素的线段,可取13条方向线或9条方向线;
对于二维超声图像,所述方向线为以当前体素为中心并经过半径为R的正方形邻边的两个体素的线段,共有4条方向线。
11.一种超声成像预处理模块,其特征在于,包括以下单元:
计算单元,用于计算在以每个体素为中心的半径为R的邻域内、以当前体素为中心的线段的均匀程度,以便从中选取N条最均匀的方向线作为2N个切线矢量,并且从中选取一条最不均匀的方向线作为法线;当对所有体素计算完毕后,即可得到一个多值矢量场;以及
平滑单元,对每个体素用其法线所代表的最大不均匀度作为梯度,先对整个梯度数据以半径为r的邻域进行一次均值滤波;若当前体素均值滤波后的梯度小于预定阈值,则计算其邻域的灰度均值并赋给所述当前体素;否则,沿着多值矢量场在所述当前体素位置的2N个切线矢量方向上分别行进,其中每一步又产生2N个分支,由此生成一个树形结构,计算树上的所有体素的灰度均值,赋给所述当前体素;所有体素计算完毕后,即完成了整个体数据的平滑。
12.如权利要求11所述的超声成像预处理模块,其特征在于,还包括:
判断单元,用于根据体素的梯度判断它是否位于物体的表面附近,然后只对表面附近的体素通过所述平滑单元进行处理,对其他区域的体素只进行均值滤波处理。
13.如权利要求11所述的超声成像预处理模块,其特征在于,还包括:
缩小单元,用于对扫描转换前的体数据进行缩小处理,然后再通过所述计算单元和所述平滑单元对缩小后的体数据进行处理;以及
放大单元,用于将通过平滑单元处理之后的体数据放大到原始大小。
14.如权利要求13所述的超声成像预处理模块,其特征在于:
其中所述缩小单元执行:对于三维超声图像,对体数据做X·Y·Z的均值滤波,再做X·Y·Z的采样,其中X、Y、Z为三个方向的缩小倍数;对于二维超声图像,对体数据做X·Y的均值滤波,再做X·Y的采样,其中X和Y为两个方向的缩小倍数;
其中所述放大单元执行:对于三维超声图像,将体数据中的每个体素直接复制成X·Y·Z个,再做X·Y·Z的均值滤波,其中X、Y、Z为三个方向的放大倍数;对于二维超声图像,将体数据中的每个体素直接复制成X·Y个,再做X·Y的均值滤波,其中X、Y为两个方向的放大倍数。
15.如权利要求14所述的超声成像预处理模块,其特征在于:
其中缩小和放大时做均值滤波的窗口大小可以取不同于X·Y·Z或X·Y的值。
16.如权利要求11至15中的任何一项所述的超声成像预处理模块,其特征在于:
其中2N个切线矢量为N条方向线的两端的指向。
17.如权利要求11至15中的任何一项所述的超声成像预处理模块,其特征在于:
R可取1;r可取1或2;对于三维超声图像N为2,对于二维超声图像N为1;所述预定阈值可取分割背景与物体的经验值。
18.如权利要求11至15中的任何一项所述的超声成像预处理模块,其特征在于,所述方向线的不均匀度V按照下式来计算:
式中G0是当前体素的灰度值,Gj是方向线上第j个元素的灰度值,其中j取负数和正数分别代表方向线的两端的指向。
19.如权利要求11至15中的任何一项所述的超声成像预处理模块,其特征在于:
对于三维超声图像,所述半径为R的邻域为以每个体素为中心、以R为半径的立方体,其中R为立方体中心到表面的距离;对于二维超声图像,所述半径为R的邻域为以每个体素为中心、以R为半径的正方形,其中R为正方形中心到其边的距离。
20.如权利要求19所述的超声成像预处理模块,其特征在于:
对于三维超声图像,所述方向线为以当前体素为中心并经过半径为R的立方体邻域表面的两个体素的线段,可取13条方向线或9条方向线;
对于二维超声图像,所述方向线为以当前体素为中心并经过半径为R的正方形邻边的两个体素的线段,共有4条方向线。
Priority Applications (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN2006101469892A CN101190132B (zh) | 2006-11-28 | 2006-11-28 | 超声成像的预处理方法与装置 |
US11/869,558 US8059914B2 (en) | 2006-11-28 | 2007-10-09 | Method and apparatus for preprocessing ultrasound imaging |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN2006101469892A CN101190132B (zh) | 2006-11-28 | 2006-11-28 | 超声成像的预处理方法与装置 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN101190132A CN101190132A (zh) | 2008-06-04 |
CN101190132B true CN101190132B (zh) | 2010-12-08 |
Family
ID=39463791
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN2006101469892A Active CN101190132B (zh) | 2006-11-28 | 2006-11-28 | 超声成像的预处理方法与装置 |
Country Status (2)
Country | Link |
---|---|
US (1) | US8059914B2 (zh) |
CN (1) | CN101190132B (zh) |
Families Citing this family (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US8553959B2 (en) * | 2008-03-21 | 2013-10-08 | General Electric Company | Method and apparatus for correcting multi-modality imaging data |
JP5632680B2 (ja) | 2010-08-25 | 2014-11-26 | 日立アロカメディカル株式会社 | 超音波画像処理装置 |
KR101881924B1 (ko) * | 2011-12-29 | 2018-08-27 | 삼성전자주식회사 | 초음파 영상 처리 방법 및 장치 |
ES2706749T3 (es) | 2013-11-04 | 2019-04-01 | Cyrill Gyger | Método para procesar datos de imagen que representan un volumen tridimensional |
CN105427261A (zh) * | 2015-11-27 | 2016-03-23 | 努比亚技术有限公司 | 一种去除图像彩色噪声的方法、装置及移动终端 |
CN106846445B (zh) * | 2016-12-08 | 2019-08-27 | 华中科技大学 | 一种基于cpu的三维超声图像体绘制方法 |
US11308609B2 (en) * | 2019-12-04 | 2022-04-19 | GE Precision Healthcare LLC | System and methods for sequential scan parameter selection |
CN112826483B (zh) * | 2021-01-08 | 2022-03-08 | 中国科学院自动化研究所 | 基于指尖视频的心率检测方法、系统及装置 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2000041982A (ja) * | 1998-07-29 | 2000-02-15 | Shimadzu Corp | 超音波診断装置 |
CN1867295A (zh) * | 2003-10-17 | 2006-11-22 | 松下电器产业株式会社 | 超声波多普勒血流测定装置 |
Family Cites Families (18)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
DE9300676U1 (zh) * | 1993-01-20 | 1993-03-11 | Eckhorn, Reinhard, Prof. Dr.-Ing., 3575 Kirchhain, De | |
WO1996021196A1 (en) * | 1994-12-30 | 1996-07-11 | Acuson Corporation | Adaptive temporal filtering to enhance fluid flow or tissue motion imaging |
US5988862A (en) * | 1996-04-24 | 1999-11-23 | Cyra Technologies, Inc. | Integrated system for quickly and accurately imaging and modeling three dimensional objects |
US5827189A (en) * | 1996-12-30 | 1998-10-27 | General Electric Company | Method and apparatus for preventing axial spatial aliasing in ultrasound imager having complex signal detector |
US6289135B1 (en) * | 1997-03-20 | 2001-09-11 | Inria Institut National De Recherche En Informatique Et En Antomatique | Electronic image processing device for the detection of motions |
US5860928A (en) * | 1997-08-07 | 1999-01-19 | Siemens Medical Systems, Inc. | Method and system for selectively smoothing color flow images in an ultrasound system |
US6847737B1 (en) * | 1998-03-13 | 2005-01-25 | University Of Houston System | Methods for performing DAF data filtering and padding |
US6001063A (en) * | 1998-06-23 | 1999-12-14 | Acuson Corporation | Ultrasonic imaging method and apparatus for providing doppler energy correction |
US6690816B2 (en) * | 2000-04-07 | 2004-02-10 | The University Of North Carolina At Chapel Hill | Systems and methods for tubular object processing |
WO2003021532A2 (en) * | 2001-09-06 | 2003-03-13 | Koninklijke Philips Electronics N.V. | Method and apparatus for segmentation of an object |
US6685641B2 (en) * | 2002-02-01 | 2004-02-03 | Siemens Medical Solutions Usa, Inc. | Plane wave scanning reception and receiver |
CN1190753C (zh) | 2003-02-28 | 2005-02-23 | 清华大学 | 基于△-σ变换的超声动态接收变迹方法 |
US6790759B1 (en) | 2003-07-31 | 2004-09-14 | Freescale Semiconductor, Inc. | Semiconductor device with strain relieving bump design |
CN1929781A (zh) * | 2003-08-21 | 2007-03-14 | 依斯克姆公司 | 用于脉管斑块检测和分析的自动化方法和系统 |
US7778493B2 (en) * | 2003-10-09 | 2010-08-17 | The Henry M. Jackson Foundation For The Advancement Of Military Medicine Inc. | Pixelation reconstruction for image resolution and image data transmission |
US8014625B2 (en) * | 2004-11-10 | 2011-09-06 | Agfa Healthcare | Method of performing measurements on digital images |
CN101112320B (zh) * | 2006-07-28 | 2010-07-21 | 深圳迈瑞生物医疗电子股份有限公司 | 波束合成的接收聚焦参数的实时计算方法及其装置 |
CN100595789C (zh) | 2006-09-21 | 2010-03-24 | 复旦大学 | 一种基于二维模糊聚类的超声图像分割方法 |
-
2006
- 2006-11-28 CN CN2006101469892A patent/CN101190132B/zh active Active
-
2007
- 2007-10-09 US US11/869,558 patent/US8059914B2/en active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2000041982A (ja) * | 1998-07-29 | 2000-02-15 | Shimadzu Corp | 超音波診断装置 |
CN1867295A (zh) * | 2003-10-17 | 2006-11-22 | 松下电器产业株式会社 | 超声波多普勒血流测定装置 |
Also Published As
Publication number | Publication date |
---|---|
US8059914B2 (en) | 2011-11-15 |
US20080123992A1 (en) | 2008-05-29 |
CN101190132A (zh) | 2008-06-04 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN101190132B (zh) | 超声成像的预处理方法与装置 | |
DE69732511T2 (de) | Verarbeitungsverfahren für Signale von Objekten mit sich bewegenden Teilen und Echographie-Vorrichtung dafür | |
Treece et al. | Surface interpolation from sparse cross sections using region correspondence | |
CN102106741A (zh) | 一种二维超声图像的三维重建方法 | |
DE69809538D1 (de) | Verfahren und Vorrichtung zur Bild-Aufnahme mit Ultraschall | |
CN105411624B (zh) | 一种超声三维流体成像与测速方法 | |
CN105982696A (zh) | 实时宽景超声成像装置及方法 | |
CN103156637A (zh) | 超声体积图像数据处理方法和设备 | |
Srinivasa et al. | Detection of edges from projections | |
CN110060337A (zh) | 颈动脉超声扫查三维重建方法及系统 | |
JP3944059B2 (ja) | 超音波診断装置 | |
Zhang et al. | Direct surface extraction from 3D freehand ultrasound images | |
CN101996415B (zh) | 眼球的三维建模方法 | |
Hafizah et al. | 3D ultrasound image reconstruction based on VTK | |
Singh et al. | Synthetic models of ultrasound image formation for speckle noise simulation and analysis | |
Lefebvre et al. | Automatic three-dimensional reconstruction and characterization of articular cartilage from high-resolution ultrasound acquisitions | |
CN106023293A (zh) | 一种基于c扫描超声图像的三维重建方法 | |
Liu et al. | A unified method based on wavelet transform and CV model for crack segmentation of 3D industrial CT images | |
CN110969694B (zh) | 无约束式扫描和基于体素的三维实时骨骼成像方法 | |
JPH08299341A (ja) | 超音波体積演算装置 | |
Thune et al. | Edge detection in noisy data using finite mixture distribution analysis | |
Sun et al. | X-FMAS Beamforming Method for 3D Imaging Using Row-Column Addressed Array | |
Chen et al. | Three-Dimensional Reconstruction of Ultrasound Images Based on Area Iteration | |
Ohashi et al. | Boundary estimation method for ultrasonic 3D imaging | |
de Jong et al. | Multiple photoacoustic sources localization using deep learning |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
EE01 | Entry into force of recordation of patent licensing contract | ||
EE01 | Entry into force of recordation of patent licensing contract |
Application publication date: 20080604 Assignee: Shenzhen Mindray Animal Medical Technology Co.,Ltd. Assignor: SHENZHEN MINDRAY BIO-MEDICAL ELECTRONICS Co.,Ltd. Contract record no.: X2022440020009 Denomination of invention: Ultrasound imaging preprocessing method and device Granted publication date: 20101208 License type: Common License Record date: 20220804 |