CN106600695B - 一种三维体重构方法和系统 - Google Patents

一种三维体重构方法和系统 Download PDF

Info

Publication number
CN106600695B
CN106600695B CN201611247431.3A CN201611247431A CN106600695B CN 106600695 B CN106600695 B CN 106600695B CN 201611247431 A CN201611247431 A CN 201611247431A CN 106600695 B CN106600695 B CN 106600695B
Authority
CN
China
Prior art keywords
probe
roi
calculating
axis coordinate
axis
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
Application number
CN201611247431.3A
Other languages
English (en)
Other versions
CN106600695A (zh
Inventor
唐艳红
艾金钦
孙慧
刘旭江
王雅儒
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sonoscape Medical Corp
Original Assignee
Sonoscape Medical Corp
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Sonoscape Medical Corp filed Critical Sonoscape Medical Corp
Priority to CN201611247431.3A priority Critical patent/CN106600695B/zh
Publication of CN106600695A publication Critical patent/CN106600695A/zh
Application granted granted Critical
Publication of CN106600695B publication Critical patent/CN106600695B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/10Constructive solid geometry [CSG] using solid primitives, e.g. cylinders, cubes
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20092Interactive image processing based on input by user
    • G06T2207/20104Interactive definition of region of interest [ROI]

Landscapes

  • Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Engineering & Computer Science (AREA)
  • Computer Graphics (AREA)
  • Software Systems (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Ultra Sonic Daignosis Equipment (AREA)
  • Dental Tools And Instruments Or Auxiliary Dental Instruments (AREA)

Abstract

本申请提供一种三维体重构方法和系统,方法包括:建立XYZ直角坐标系,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,原点O为探头晶体圆心;设三维体数据中任意一点P在z=0平面上的投影点为S,过S点与转动轴L垂直的直线为L”,L与L”的交点为T;延长直线TS到点W,点W满足条件:|TW|=|TP|,点P由点W旋转phi角度得到,phi为探头的摆动角度;依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构;无需构建复杂的体重构查询表,因此,降低了用户工作量,同时也无需花费大量的时间由查询表中进行数据查询并取值,从而提高了重构速率。

Description

一种三维体重构方法和系统
技术领域
本发明涉及图像处理技术领域,具体涉及一种三维体重构方法和系统。
背景技术
三维重构是指对三维物体构建适合计算机表示和处理的数学模型,是在计算机环境下对其进行处理、操作和分析其性质的基础,也是在计算机中建立表达客观世界的虚拟现实的关键技术。
迈瑞公司申请的美国专利US7,988,632 B2中,提到了一种三维体重构方法:事先构建重构表,实时计算时,根据当前采样点的直角坐标查表得到其极坐标,插值计算采样点的灰度值。其中重构表保存从直角坐标系(x,y,z)到极坐标(R,theta,phi)的变换关系为:
Figure BDA0001197400450000011
phi=tan-1(z/y)
Figure BDA0001197400450000012
theta=tan-1(x/(pr-fr+R′))
申请人通过分析发现,事先构建的体重构查询表数据量很大,查询表构建过程耗时较长,且在实时重构时需花费较多时间从查询表中取值。
发明内容
有鉴于此,本发明实施例提供一种三维体重构方法和系统,以降低三维体重构过程中用户的工作量、以及提交重构速度。
为实现上述目的,本发明实施例提供如下技术方案:
一种三维体重构方法,包括:
扫描用户选择的感兴趣区域ROI;
建立XYZ直角坐标系,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
设三维体数据中任意一点P(x,y,z)在z=0平面上的投影点为S(x,y,0),过S点与转动轴L垂直的直线为L”,转动轴L与直线L”的交点为T,所述转动轴L为探头的摆动轴;
延长直线TS到点W(x1,y1,0),所述点W满足条件:|TW|=|TP|,|TW|指的是点T到点W之间的距离,|TP|是点T到点P之间的距离,且点P由点W旋转phi角度得到,所述phi为探头在马达带动下的摆动角度;
依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构,其中,R为扫描线方向上的半径,theta为扫描线所在的角度位置;
phi=tan-1(-z/|TS|) (公式1)
Figure BDA0001197400450000021
Figure BDA0001197400450000022
其中,所述|TS|表示点S到转动轴L的距离;其中所述Rprobe为探头晶体半径。
优选的,上述三维体重构方法中,依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式之前,还可以包括:
依据公式4计算得到|TS|;
|TS|=|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||(公式4)
其中,Rswing为探头摆动半径、Δθ为Y轴和直角坐标系X’Y’Z中的Y’轴的夹角,X’轴与探头的支柱方向一致;
采用公式4,对公式1进行变换,变换后的公式1为:phi=tan-1(-z/|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||)。
优选的,上述三维体重构方法中,延长直线TS到点W(x1,y1,0)之前,包括:
依据公式5计算得到点W的坐标;
Figure BDA0001197400450000031
其中,所述Δx为,Δy为,
Figure BDA0001197400450000032
优选的,上述三维体重构方法中,所述依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,包括:
采用计算着色器依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构。
一种三维体重构方法,包括:
扫描用户选择的感兴趣区域ROI;
建立XYZ坐标系,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
计算XY平面内的第一关键参数,所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam;
计算所述感兴趣区域ROI的四个顶点的坐标;
计算所述感兴趣区域ROI偏移后的四个顶点的坐标;
计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数;
计算ZY平面内的第二关键参数,所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴、探头摆动的轴向长度、探头摆动的横向长度、探头摆动时与Y轴的偏移量;
计算XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
依据所述第一关键参数、第二关键参数、感兴趣区域POI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构。
优选的,上述三维体重构方法中,所述计算XY平面内的第一关键参数,包括:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6)
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/2.0)) (公式7)
优选的,上述三维体重构方法中,所述计算所述感兴趣区域ROI的四个顶点的坐标,包括:
依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点的坐标Pta、Ptb、Ptc和Ptd;
Figure BDA0001197400450000041
Figure BDA0001197400450000042
Figure BDA0001197400450000043
Figure BDA0001197400450000044
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐,标Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
优选的,上述三维体重构方法中,所述计算所述感兴趣区域ROI偏移后的四个顶点的坐标,包括:
将计算所述感兴趣区域ROI的四个顶点的坐标绕原点旋转矩阵
Figure BDA0001197400450000051
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
优选的,上述三维体重构方法中,计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数,包括,
求取ROI感兴趣区域偏移后的四个顶点的坐标最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为外接矩形ABCD;如果是,比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的值记为min_y,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|。
优选的,上述三维体重构方法中,计算ZY平面内的第二关键参数,包括,
依据公式9.1计算得到探头摆动半径MotorCurve;
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemoror/2.0) (公式9.7)。
优选的,上述三维体重构方法中,进行三维重构,包括:
遍历每一张切面,利用插值系数表TableXY执行双线性插值;
对XY平面内每一张切片执行旋转操作,获得旋转体;
将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
利用ZY平面内的第二关键参数以及插值系数表TableZY执行双线性插值。
一种三维体重构系统,包括:
扫面单元,用于扫描用户选择的感兴趣区域ROI;
坐标系建立单元,用于建立XYZ直角坐标系和直角坐标系X’Y’Z,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心,X’轴与探头的支柱方向一致;
重构单元,用于:设三维体数据中任意一点P(x,y,z)在z=0平面上的投影点为S(x,y,0),过S点与转动轴L垂直的直线为L”,转动轴L与直线L”的交点为T,所述转动轴L为探头的摆动轴;延长直线TS到点W(x1,y1,0),所述点W满足条件:|TW|=|TP|,|TW|指的是点T到点W之间的距离,|TP|是点T到点P之间的距离,且点P由点W旋转phi角度得到,所述phi为探头在马达带动下的摆动角度;依据公式一、公式二和公式三得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构,其中,R为扫描线方向上的半径,theta为扫描线所在的角度位置;
phi=tan-1(-z/|TS|) (公式一)
Figure BDA0001197400450000071
Figure BDA0001197400450000072
其中,所述|TS|表示点S到转动轴L的距离;其中所述Rprobe为探头晶体半径。
优选的,上述三维体重构方法中,所述重构单元,还包括:
依据公式4计算得到|TS|;
Figure BDA0001197400450000073
其中,A为、B为、C为、Rswing为探头摆动半径、Δθ为Y轴和直角坐标系X’Y’Z中的Y’轴的夹角,X’轴与探头的支柱方向一致;
采用公式4,对公式1进行变换,变换后的公式1为:phi=tan-1(-z/|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||)。
优选的,上述三维体重构方法中,所述重构单元,还包括:
第二计算单元,用于依据公式5计算得到点W的坐标;
Figure BDA0001197400450000074
其中,所述Δx为,Δy为,
Figure BDA0001197400450000075
优选的,上述三维体重构方法中,所述重构单元,包括:
子重构单元,用于采用计算着色器依据公式一、公式二和公式三得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构。
一种三维体重构系统,包括:
扫描单元,用于扫描用户选择的感兴趣区域ROI;
坐标系建立单元,用于建立XYZ坐标系,所述感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
第一参数计算单元,用于计算XY平面内的第一关键参数,所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam;
第一ROI顶点计算单元,用于计算所述感兴趣区域ROI的四个顶点的坐标;
第二ROI顶点计算单元,用于计算所述感兴趣区域ROI偏移后的四个顶点的坐标;
外接矩形尺寸计算单元,用于计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数;
第二参数计算单元,用于计算ZY平面内的第二关键参数,所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴、探头摆动的轴向长度、探头摆动的横向长度、探头摆动时与Y轴的偏移量;
二维表建立单元,用于计算XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
重构单元,用于依据所述第一关键参数、第二关键参数、感兴趣区域POI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构。
优选的,上述三维体重构系统中,所述第一参数计算单元,具体用于:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6)
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/2.0)) (公式7)
优选的,上述三维体重构系统中,所述第一ROI顶点计算单元,具体用于:
依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点的坐标Pta、Ptb、Ptc和Ptd;
Figure BDA0001197400450000091
Figure BDA0001197400450000092
Figure BDA0001197400450000093
Figure BDA0001197400450000094
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐,标Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
优选的,上述三维体重构系统中,所述第二ROI顶点计算单元,具体用于:
将计算所述感兴趣区域ROI的四个顶点的坐标绕原点旋转矩阵
Figure BDA0001197400450000095
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
优选的,上述三维体重构系统中,所述外接矩形尺寸计算单元,具体用于:
求取ROI感兴趣区域偏移后的四个顶点的坐标最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为接矩形ABCD;如果是,比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的坐标记为min_y,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|。
优选的,上述三维体重构系统中,所述第二参数计算单元,具体用于:
依据公式9.1计算得到探头摆动半径MotorCurve;
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemotor/2.0) (公式9.7)。
优选的,上述三维体重构系统中,所述第二参数计算单元,具体用于:
遍历每一张切面,利用插值系数表TableXY执行双线性插值;
对XY平面内每一张切片执行旋转操作,获得旋转体;
将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
利用ZY平面内的第二关键参数以及插值系数表TableZY执行双线性插值。
基于上述技术方案,本发明实施例提供的三维体重构方法和系统,通过建立XYZ直角坐标系,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;设三维体数据中任意一点P(x,y,z)在z=0平面上的投影点为S(x,y,0),过S点与转动轴L垂直的直线为L”,转动轴L与直线L”的交点为T,所述转动轴L为探头的摆动轴;延长直线TS到点W(x1,y1,0),所述点W满足条件:|TW|=|TP|,|TW|指的是点T到点W之间的距离,|TP|是点T到点P之间的距离,且点P由点W旋转phi角度得到,所述phi为探头在马达带动下的摆动角度;依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构,无需构建体重构查询表,因此,降低了用户工作量,同时也无需花费大量的时间由查询表中进行数据查询并取值,从而提高了重构速率。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据提供的附图获得其他的附图。
图1为本申请实施例公开的一种三维体重构方法的流程示意图;
图2为本申请实施例公开的一种三维体重构方法的原理示意图;
图3为本申请实施例公开的一种三维体重构系统的结构示意图;
图4为本申请另一实施例公开的一种三维体重构方法的流程示意图;
图5为本申请实施例公开的一种三维体重构方法中XY平面内参数计算方法原理示意图;
图6为本申请实施例公开的一种三维体重构方法中二维查找表构建原理图;
图7为本申请实施例公开的一种三维体重构方法中ROI与坐标Y轴相交和不相交时的示意图;
图8为本申请实施例公开的一种三维体重构方法中ZY平面内摆动示意图;
图9a为本申请实施例公开的一种三维体重构方法中XY平面内扫描变换示意图;
图9b为本申请实施例公开的一种三维体重构方法中体旋转示意图;
图9c为本申请实施例公开的一种三维体重构方法中体数据存储顺序变换示意图;
图9d为本申请实施例公开的一种三维体重构方法中ZY平面内扫描变换示意图;
图10为本申请另一实施例公开的一种三维体重构系统的结构示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
申请人经研究还发现,US7,988,632 B2采用的技术方案(从直角坐标到极坐标变换关系)只适用于用户选择的感兴趣区域RIO中心线与容积探头晶体中心线重合的情况,实时三维成像中通常会采用凸阵容积探头获取三维体积数据,由于对象条件的限制,用户选择的感兴趣区域ROI的中心线常常与探头扫描中心线偏离,例如,实际扫描的六面体左右两面应为圆锥面,而采用上述变换关系得到的重构体左右平面为斜平面,与实际扫查对象形状不符,产生变形,使得被观察对象测量大小存在误差,容易造成医生误判。
针对于上述问题,本申请公开了一种三维体重构方法和系统,参见图1,方法可以包括:
步骤S101:扫描用户选择的感兴趣区域ROI;
当应用于本方法的系统处于在3D/4D预激活状态时,起始扫查切面内用户选择的感兴趣区域ROI(abcd)。
步骤S102:构建直角坐标系;
在本步骤中,依据所述感兴趣区域ROI(abcd)建立XYZ直角坐标系和X’Y’Z直角坐标系,在z=0的平面内直角坐标系为XOY,且O点为探头晶体圆心;所述X’Y’Z直角坐标系中的X’轴与探头的支柱方向一致,所述感兴趣区域ROI与所述XYZ直角坐标系和X’Y’Z直角坐标系的位置关系如图2所示,其中,所述感兴趣区域ROI的中心线与Y轴重合,与X’Y’Z直角坐标系中的ROI固定中心线位置(Y′轴)不重合,图2中的m为探头摆动轴圆心;
步骤S103:在三维体数据中任意一点P(x,y,z),对点P在z=0平面上进行投影,得到投影点S(x,y,0),设过S点与转动轴L垂直的直线为L”,转动轴L与直线L”的交点为T,所述转动轴L为探头的摆动轴;
步骤S104:延长直线TS到点W(x1,y1,0);
其中,所述点W满足条件:|TW|=|TP|,|TW|指的是点T到点W之间的距离,|TP|指的是点T到点P之间的距离,且点P由点W旋转phi角度得到,所述phi为探头在马达带动下的摆动角度;
步骤S105:依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构,其中,R为扫描线方向上的半径,theta为扫描线所在的角度位置;
phi=tan-1(-z/|TS|) (公式1)
Figure BDA0001197400450000131
Figure BDA0001197400450000132
其中,所述|TS|表示点S到转动轴L的距离,可由点到直线的距离公式进行计算;其中所述Rprobe为探头晶体半径。
在本申请上述实施例公开的方法进行三维体重构时,无需构建体重构查询表,因此,降低了用户工作量,同时也无需花费大量的时间由查询表中进行数据查询并取值,从而提高了重构速率。
在本申请另外一实施例公开的技术方案中,在依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式之前,为了减小系统工作量,本申请还可对所述公式1进行转换,即采用转换后的公式1公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,因此,其具体形式为:步骤S105之前还可以包括:
依据公式4计算得到|TS|;
|TS|=|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||(公式4)
其中,Rswing为探头摆动半径、Δθ为XYZ直角坐标系中Y轴和直角坐标系X’Y’Z中的Y’轴的夹角;
将公式4带入公式1,对公式1进行变换,变换后的公式1为:phi=tan-1(-z/|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||) (公式1.1)
在本申请一实施例公开的计算方案中,点W的坐标位置的确立过程可以包括:
依据公式5计算得到点W的坐标;
Figure BDA0001197400450000141
其中,所述
Figure BDA0001197400450000142
z为直角坐标(x,y,z)中的坐标值z。
在本申请另一实施例公开的技术方案中,上述方法可以应用于具有独立图像处理器(GPU)的系统中,可利用GPU能够海量数据运算速度的优势,采用计算着色器依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构。
针对于上述三维体重构方法,本申请公开了一种三维体重构系统,上述方法与系统可以相互借鉴,参见图3,该系统可以包括:
扫面单元100,用于扫描用户选择的感兴趣区域ROI;
坐标系建立单元200,用于依据所述感兴趣区域RIO构建XYZ直角坐标系和X’Y’Z直角坐标系;
其中,在z=0的平面内直角坐标系为XOY,且O点为探头晶体圆心;所述X’Y’Z直角坐标系中的X’轴与探头的支柱方向一致,所述感兴趣区域ROI与所述XYZ直角坐标系和X’Y’Z直角坐标系的位置关系如图2所示,所述感兴趣区域ROI的中心线与Y轴重合,与X’Y’Z直角坐标系中的ROI固定中心线位置(Y′轴)不重合,图2中的m为探头摆动轴圆心;
重构单元300,用于:在三维体数据中任意取一点P(x,y,z),对点P在z=0平面上进行投影,得到投影点S(x,y,0),设过S点与转动轴L垂直的直线为L”,转动轴L与直线L”的交点为T,所述转动轴L为探头的摆动轴;延长直线TS到点W(x1,y1,0),所述点W满足条件:|TW|=|TP|,|TW|指的是点T到点W之间的距离,|TP|指的是点T到点P之间的距离,且点P由点W旋转phi角度得到,所述phi为探头在马达带动下的摆动角度;依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构,其中,R为扫描线方向上的半径,theta为扫描线所在的角度位置;
phi=tan-1(-z/|TS|) (公式1)
Figure BDA0001197400450000151
Figure BDA0001197400450000152
其中,所述|TS|表示点S到转动轴L的距离;其中所述Rprobe为探头晶体半径。
与上述方法相对应,所述重构单元300,还包括公式1变电换单元,用于依据公式4计算得到|TS|;
|TS|=|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing|| (公式4)
其中,Rswing为探头摆动半径、Δθ为Y轴和直角坐标系X’Y’Z中的Y’轴的夹角,X’轴与探头的支柱方向一致;
将公式4带入公式,对公式1进行变换,变换后的公式1为:phi=tan-1(-z/|sinΔθ*x-cosΔθ*y-|Rprobe-Rswing||) (公式1.1);
此时,依据公式1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换等效为依据公式1.1、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换。
与上述方法相对应,所述重构单元300,还包括:
第二计算单元,用于依据公式5计算得到点W的坐标;
Figure BDA0001197400450000161
其中,
Figure BDA0001197400450000162
与上述方法相对应,当所述三维体重构系统应用于具有独立图像处理器(GPU)的系统中时,所述重构单元300,可以包括:
子重构单元,用于采用计算着色器依据公式1(或1.1)、公式2和公式3得到(R,theta,phi)与直角坐标(x,y,z)之间的转换公式,并进行三维体重构。
针对于,背景技术中提出的上述问题,本申请还公开了一种三维体重构方法,该方法可以应用于通过CPU实现三维重构计算的系统中,参见图4,该三维体重构方法包括:
步骤S401:扫描用户选择的感兴趣区域ROI;
步骤S402:构建XYZ坐标系;
依据所述感兴趣区域ROI建立XYZ坐标系,所述感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
本步骤中,XY平面内坐标系的建立如图5所示,用户感兴趣区域ROI的中心线为Y轴,感兴趣区域ROI的中心线过探头晶体圆心O,X轴垂直于Y轴;
步骤S403:计算XY平面内的第一关键参数;
所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam,所述第一关键参数和感兴趣区域ROI用于XY平面内插值;
步骤S404:计算所述感兴趣区域ROI的四个顶点(abcd)的坐标;
步骤S405:计算所述感兴趣区域ROI偏移后的四个顶点(abcd)的坐标;
步骤S406:计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数;
由于重构体数据为规则体数据,因此在本方法中,需要通过的便宜后的感兴趣区域的ROI的外接矩形限定XY平面内进行插值后的平面的大小;
步骤S407:计算ZY平面内的第二关键参数;
所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,探头摆动的轴向长度、探头摆动的横向长度、探头摆动时与Y轴的偏移量以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴,所述第二关键参数用于ZY平面内差值;
步骤S408:构建XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
具体的,构建查找表TableXY和TableZY的构建过程如图6所示,在计算两个平面(XY平面和ZY平面)的二维查找表TableXY和TableZY时,可采用双线性计算插值表进行计算,该查找表TableZY内存储源数据点坐标以及感兴趣区域RIO的四个顶点的权重值;
两个查找表大小分别为TableXYSize和TableZYSize;
在本申请公开的技术方案中,TableXYSize=20*dstRow*dstCol、TableZYSize=20*dstRow*dstSlice;其中,dstRow、dstCol分别代表在XY平面插值后,XY平面的高度和宽度。dstRow、dstSlice分别代表在ZY平面插值后,ZY平面的高度和宽度。
步骤S409:依据所述第一关键参数、第二关键参数、感兴趣区域ROI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构;
在本步骤中依据用户感兴趣区域ROI和计算得到的第一关键参数,进行XY平面内插值,通过外接矩形的尺寸限定XY平面内插值后的平面大小,通过第二关键参数进行ZY差值从而实现三维重构。
在本申请上述实施例公开的方案中,无需构建体重构查询表,因此,降低了用户工作量,同时也无需花费大量的时间由查询表中进行数据查询并取值,从而提高了计算速率。
在本申请一实施例公开的方案中,还对各个第一关键参数的是如何得到的进行了说明,即所述计算XY平面内的第一关键参数,可以包括:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6);
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/22.0)) (公式7);
Angleprobe/2.0表示ROI张角的一半。
在本申请一实施例公开的方案中,还对感兴趣区域ROI的四个顶点的坐标的是如何得到的进行了说明,所述计算所述感兴趣区域ROI的四个顶点的坐标,可以包括:
感兴趣区域ROI的存储坐标与数学计算笛卡尔坐标之间的关系如图5所示,依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点(abcd)的坐标Pta、Ptb、Ptc和Ptd;
Figure BDA0001197400450000191
Figure BDA0001197400450000192
Figure BDA0001197400450000193
Figure BDA0001197400450000194
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐,标Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
针对于上述实施例公开的技术方案,考虑ROI与Y轴有偏转,需要计算顶点(a,b,c,d)偏转后的直角坐标。具体的,可以利用二维平面内的四个点绕坐标原点旋转矩阵
Figure BDA0001197400450000195
计算获得感兴趣区域ROI的四个顶点在坐标系中的实际位置:Pta1,Ptb1,Ptc1,Ptd1,其旋转后的Pta1、Ptb1、Ptc1和Ptd1如图7所示。即,所述计算所述感兴趣区域ROI偏移后的四个顶点的坐标,可以包括:
将计算所述感兴趣区域ROI的四个顶点的坐标(Pta、Ptb、Ptc和Ptd)绕原点旋转矩阵
Figure BDA0001197400450000201
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
在本申请上述实施例公开的计算方案中,所述计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数,可以包括,
计算ROI感兴趣区域偏移后的四个顶点的坐标中的最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,通过所述最大X轴坐标和最小X轴坐标,参见图8,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为外接矩形ABCD;如果相交,则比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的值记为min_y,此时,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|,所述|AC|=|min_y-max_y|=BoxAxisMm,max_y为Y轴最大坐标,当ROI与Y轴不相交时,min_y为Y轴最小坐标;
所述|AB|=|min_x-max_x|=BoxLateralMm,max_x为Y轴最大坐标,min_x为X轴最小坐标;所述外接矩形ABCD靠近X轴的边与X轴之间的距离为D_ABX=|max_y|。
在本申请上述实施例公开的技术方案中,由于机械设计尺寸限制,容积探头摆动半径原点和探头扫查半径原点并不重合,如图8所示,O点为探头晶体圆心,O’为探头摆动圆心。在摆动ZY平面内外接矩形ABCD绕轴L转动。所述计算ZY平面内的第二关键参数可以包括:
依据公式9.1计算得到探头摆动半径MotorCurve(探头摆动圆心O’到外接矩形ABCD靠近X轴的边的距离);
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemotor/2.0) (公式9.7)。
在本申请上述实施例公开的技术方案中进行三维重构的具体过程可以包括:
参见图9a,遍历每一张切面,利用插值系数表TableXY执行双线性插值;
其中,图9a中所示(inpSlice,inpRow,inpCols)分别表示输入体数据的点数、线数和切面数量;(intCol1,intRow1,intSlice1)分别表示执行XY平面内二维执行双线性插值后的体数据点数、线数和切面数量;
参见图9b,对XY平面内每一张切片执行旋转操作,获得旋转体;
其中,图9b中(rotCol,rotRow,rotSlice)分别表示对每一张切面进行旋转后的点数、线数和切面数量;
参见图9c,将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
其中,图9c中(intCol1,intRow1,intSlice1)分别表示变换数据存储顺序后的点数、线数和切面数量;
参见图9d,利用ZY平面内的第二关键参数以及插值系数表TableZY执行双线性插值,即可实现三维体重构;
其中,图9d中(intCol3,intRow3,intSlice3)分别表示执行第二次二维插值后的体数据点数、线数和切面数量。
与上述方法相对应,本申请还公开了一种三维体重构系统A,参见图10,系统A可以包括:
扫描单元110,用于扫描用户选择的感兴趣区域ROI;
坐标系建立单元120,用于依据所述感兴趣区域ROI建立XYZ坐标系,所述感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;用户感兴趣区域ROI的中心线为Y轴,感兴趣区域ROI的中心线过探头晶体圆心O,X轴垂直于Y轴;
第一参数计算单元130,用于计算XY平面内的第一关键参数,所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam;
第一ROI顶点计算单元140,用于计算所述感兴趣区域ROI的四个顶点的坐标;
第二ROI顶点计算单元150,用于计算所述感兴趣区域ROI偏移后的四个顶点的坐标;
外接矩形尺寸计算单元160,用于计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数;
第二参数计算单元170,用于计算ZY平面内的第二关键参数,所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴、探头摆动的轴向长度、探头摆动的横向长度、探头摆动时与Y轴的偏移量;
二维表建立单元180,用于构建XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
重构单元190,用于依据所述第一关键参数、第二关键参数、感兴趣区域POI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构。
与上述方法相对应,所述第一参数计算单元130,具体用于:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6)
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/2.0)) (公式7)。
与上述方法相对应,所述第一ROI顶点计算单元140,具体用于:
依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点的坐标Pta、Ptb、Ptc和Ptd;
Figure BDA0001197400450000231
Figure BDA0001197400450000232
Figure BDA0001197400450000233
Figure BDA0001197400450000234
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐,标Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
与上述方法相对应,所述第二ROI顶点计算单元150,具体用于:
将计算所述感兴趣区域ROI的四个顶点的坐标绕原点旋转矩阵
Figure BDA0001197400450000241
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
与上述方法相对应,所述外接矩形尺寸计算单元160,具体用于:
计算ROI感兴趣区域偏移后的四个顶点的坐标中的最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,通过所述最大X轴坐标和最小X轴坐标,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为外接矩形ABCD;如果相交,则比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的值记为min_y,此时,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|,所述|AC|=|min_y-max_y|=BoxAxisMm,max_y为Y轴最大坐标,当ROI与Y轴不相交时,min_y为Y轴最小坐标;
所述|AB|=|min_x-max_x|=BoxLateralMm,max_x为Y轴最大坐标,min_x为X轴最小坐标;所述外接矩形ABCD靠近X轴的边与X轴之间的距离为D_ABX=|max_y|。
与上述方法相对应,所述第二参数计算单元170,具体用于:
依据公式9.1计算得到探头摆动半径MotorCurve(探头摆动圆心O’到外接矩形ABCD靠近X轴的边的距离);
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemotor/2.0) (公式9.7)
与上述方法相对应,所述重构单元190,具体用于:
遍历每一张切面,利用插值系数表TableXY执行双线性插值;
对XY平面内每一张切片执行旋转操作,获得旋转体;
将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
利用ZY平面内的第二关键参数以及插值系数表TableZY执行双线性插值,以实现三维体重构。
为了描述的方便,描述以上系统时以功能分为各种模块分别描述。当然,在实施本申请时可以把各模块的功能在同一个或多个软件和/或硬件中实现。
本说明书中的各个实施例均采用递进的方式描述,各个实施例之间相同相似的部分互相参见即可,每个实施例重点说明的都是与其他实施例的不同之处。尤其,对于系统或系统实施例而言,由于其基本相似于方法实施例,所以描述得比较简单,相关之处参见方法实施例的部分说明即可。以上所描述的系统及系统实施例仅仅是示意性的,其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的,作为单元显示的部件可以是或者也可以不是物理单元,即可以位于一个地方,或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性劳动的情况下,即可以理解并实施。
专业人员还可以进一步意识到,结合本文中所公开的实施例描述的各示例的单元及算法步骤,能够以电子硬件、计算机软件或者二者的结合来实现,为了清楚地说明硬件和软件的可互换性,在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行,取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能,但是这种实现不应认为超出本发明的范围。
结合本文中所公开的实施例描述的方法或算法的步骤可以直接用硬件、处理器执行的软件模块,或者二者的结合来实施。软件模块可以置于随机存储器(RAM)、内存、只读存储器(ROM)、电可编程ROM、电可擦除可编程ROM、寄存器、硬盘、可移动磁盘、CD-ROM、或技术领域内所公知的任意其它形式的存储介质中。
还需要说明的是,在本文中,诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下,在其它实施例中实现。因此,本发明将不会被限制于本文所示的这些实施例,而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (12)

1.一种三维体重构方法,其特征在于,包括:
扫描用户选择的感兴趣区域ROI;
建立XYZ坐标系,感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
计算XY平面内的第一关键参数,所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam,所述第一关键参数和感兴趣区域ROI用于XY平面内插值;
计算所述感兴趣区域ROI的四个顶点的坐标;
计算所述感兴趣区域ROI偏移后的四个顶点的坐标;
计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数,所述外接矩形ABCD的尺寸参数用于限定XY平面内插值后的平面大小;
计算ZY平面内的第二关键参数,所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴,探头摆动的轴向长度,探头摆动的横向长度,探头摆动时与Y轴的偏移量,所述第二关键参数用于ZY平面内插值;
计算XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
依据所述第一关键参数、第二关键参数、感兴趣区域ROI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构;
所述三维重构的具体过程包括:
遍历每一张切面,利用二维查找表TableXY执行双线性插值;
对XY平面内每一张切片执行旋转操作,获得旋转体;
将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
利用ZY平面内的第二关键参数以及二维查找表TableZY执行双线性插值,实现三维体重构。
2.根据权利要求1所述的三维体重构方法,其特征在于,所述计算XY平面内的第一关键参数,包括:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6)
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/2.0)) (公式7)。
3.根据权利要求2所述的三维体重构方法,其特征在于,所述计算所述感兴趣区域ROI的四个顶点的坐标,包括:
依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点的坐标Pta、Ptb、Ptc和Ptd;
Figure FDA0002221292120000021
Figure FDA0002221292120000022
Figure FDA0002221292120000023
Figure FDA0002221292120000024
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐标,Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
4.根据权利要求3所述的三维体重构方法,其特征在于,所述计算所述感兴趣区域ROI偏移后的四个顶点的坐标,包括:
将计算所述感兴趣区域ROI的四个顶点的坐标绕原点旋转矩阵
Figure FDA0002221292120000031
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
5.根据权利要求4所述的三维体重构方法,其特征在于,计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数,包括,
求取ROI感兴趣区域偏移后的四个顶点的坐标最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为外接矩形ABCD;如果是,比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的值记为min_y,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|。
6.根据权利要求5所述的三维体重构方法,其特征在于,计算ZY平面内的第二关键参数,包括,
依据公式9.1计算得到探头摆动半径MotorCurve;
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemotor/2.0) (公式9.7)。
7.一种三维体重构系统,其特征在于,包括:
扫描单元,用于扫描用户选择的感兴趣区域ROI;
坐标系建立单元,用于建立XYZ坐标系,所述感兴趣区域的中心线与XYZ坐标系中的Y轴重合,在z=0平面内直角坐标系为XOY,且O点为探头晶体圆心;
第一参数计算单元,用于计算XY平面内的第一关键参数,所述第一关键参数包括:所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve和外接矩形框在Y方向上的长度ProbeBeam,所述第一关键参数和感兴趣区域ROI用于XY平面内插值;
第一ROI顶点计算单元,用于计算所述感兴趣区域ROI的四个顶点的坐标;
第二ROI顶点计算单元,用于计算所述感兴趣区域ROI偏移后的四个顶点的坐标;
外接矩形尺寸计算单元,用于计算所述感兴趣区域ROI的外接矩形ABCD的尺寸参数,所述外接矩形ABCD的尺寸参数用于限定XY平面内插值后的平面大小;
第二参数计算单元,用于计算ZY平面内的第二关键参数,所述第二关键参数包括:探头摆动半径,探头摆动轴的有效长度,探头摆动轴总长度,以及ZY平面内摆动外接矩形与转动轴L的偏移量,所述转动轴L为探头的摆动轴,探头摆动的轴向长度,探头摆动的横向长度,探头摆动时与Y轴的偏移量,所述第二关键参数用于ZY平面内插值;
二维表建立单元,用于计算XY平面的二维查找表TableXY和ZY平面内的二维查找表TableZY;
重构单元,用于依据所述第一关键参数、第二关键参数、感兴趣区域ROI的四个顶点坐标、感兴趣区域偏移后的四个顶点坐标以及感兴趣区域ROI的外接矩形ABCD的尺寸参数进行三维重构;
所述三维重构的具体过程包括:
遍历每一张切面,利用二维查找表TableXY执行双线性插值;
对XY平面内每一张切片执行旋转操作,获得旋转体;
将旋转体数据存储顺序由存储顺序XYZ变为ZYX;
利用ZY平面内的第二关键参数以及二维查找表TableZY执行双线性插值,实现三维体重构。
8.根据权利要求7所述的三维体重构系统,其特征在于,所述第一参数计算单元,具体用于:
依据公式(6)计算所述感兴趣区域ROI外接矩形框与X轴的最小偏移量ProbeCurve;
ProbeCurve=Rprobe*cos(Angleprobe/2.0) (公式6)
其中,所述Angleprobe为感兴趣区域ROI的扇形角,Rprobe为探头晶体半径;
依据公式(7)计算得到所述外接矩形框在Y方向上的长度ProbeBeam;
ProbeBeam=|ef|+Rprobe*(1.0-cos(Angleprobe/2.0)) (公式7)。
9.根据权利要求8所述的三维体重构系统,其特征在于,所述第一ROI顶点计算单元,具体用于:
依据公式8.1、8.2、8.3和8.4计算得到所述感兴趣区域ROI的四个顶点的坐标Pta、Ptb、Ptc和Ptd;
Figure FDA0002221292120000061
Figure FDA0002221292120000062
Figure FDA0002221292120000063
Figure FDA0002221292120000064
其中,Pta.x为点Pta的X轴坐标,Pta.y为点Pta的Y轴坐,标Ptb.x为点Ptb的X轴坐标,Ptb.y为点Ptb的Y轴坐标,Ptc.x为点Ptc的X轴坐标,Ptc.y为点Ptc的Y轴坐标,Ptd.x为点Ptd的X轴坐标,Ptd.y为点Ptd的Y轴坐标。
10.根据权利要求9所述的三维体重构系统,其特征在于,所述第二ROI顶点计算单元,具体用于:
将计算所述感兴趣区域ROI的四个顶点的坐标绕原点旋转矩阵
Figure FDA0002221292120000065
得到所述感兴趣区域ROI偏移后的四个顶点的坐标。
11.根据权利要求10所述的三维体重构系统,其特征在于,所述外接矩形尺寸计算单元,具体用于:
求取ROI感兴趣区域偏移后的四个顶点的坐标最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值,判断所述感兴趣区域ROI与Y轴是否相交,如果否,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值和最小Y轴坐标值围成的矩形区域即为外接矩形ABCD;如果是,比较最小Y轴坐标值和-(ProbeCurve+ProveBeam)的大小,将两者之间较小的坐标记为min_y,所述最大X轴坐标值、最小X轴坐标值、最大Y轴坐标值以及min_y所围成的区域即为外接矩形ABCD;
计算得到所述外接矩形ABCD的高度|AC|和宽度|AB|。
12.根据权利要求11所述的三维体重构系统,其特征在于,所述第二参数计算单元,具体用于:
依据公式9.1计算得到探头摆动半径MotorCurve;
MotorCurve=D_ABX-(Rprobe-Rswing) (公式9.1)
其中,Rprobe为探头晶体半径,Rswing为探头摆动半径,D_ABX外接矩形ABCD的最大Y轴坐标;
依据公式9.2计算得到探头摆动轴的有效长度MotorBeam;
MotorBeam=|AC| (公式9.2)
依据公式9.3计算得到探头摆动轴总长度MotorTRadius;
MotorTRadius=MotorCurve+MotorBeam; (公式9.3)
依据公式9.4计算得到ZY平面内摆动外接矩形与L轴的偏移量为MotorAxialOffset;
MotorAxialOffset=MotorCurve*cos(Anglemotor/2.0) (公式9.4)
依据公式9.5计算得到探头摆动的轴向长度MotorAxial;
MotorAxial=(MotorCurve-MotorAxialOffset)+MotorBeam (公式9.5)
依据公式9.6计算得到探头摆动的横向长度MotorLateral;
MotorLateral=2.0*MotorTRadius*sin(Anglemotor/2.0) (公式9.6)
依据公式9.7计算得到探头与Y轴的偏移量MotorLateralOffset:
MotorLateralOffset=MotorTRadius*sin(Anglemotor/2.0) (公式9.7)。
CN201611247431.3A 2016-12-29 2016-12-29 一种三维体重构方法和系统 Active CN106600695B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201611247431.3A CN106600695B (zh) 2016-12-29 2016-12-29 一种三维体重构方法和系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201611247431.3A CN106600695B (zh) 2016-12-29 2016-12-29 一种三维体重构方法和系统

Publications (2)

Publication Number Publication Date
CN106600695A CN106600695A (zh) 2017-04-26
CN106600695B true CN106600695B (zh) 2020-04-10

Family

ID=58605120

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201611247431.3A Active CN106600695B (zh) 2016-12-29 2016-12-29 一种三维体重构方法和系统

Country Status (1)

Country Link
CN (1) CN106600695B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111053574B (zh) * 2019-12-05 2022-10-04 深圳开立生物医疗科技股份有限公司 一种三维表面重建方法、装置及超声设备和存储介质

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1891159A (zh) * 2005-07-01 2007-01-10 深圳迈瑞生物医疗电子股份有限公司 用于三维超声成像的快速扫描变换方法
CN101292883A (zh) * 2007-04-23 2008-10-29 深圳迈瑞生物医疗电子股份有限公司 超声三维快速成像方法及其装置
CN102106741A (zh) * 2009-12-25 2011-06-29 东软飞利浦医疗设备系统有限责任公司 一种二维超声图像的三维重建方法
CN104094319A (zh) * 2012-01-19 2014-10-08 株式会社东芝 图像处理设备、立体图像显示设备和图像处理方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR100751852B1 (ko) * 2003-12-31 2007-08-27 주식회사 메디슨 대상체의 3차원 초음파 데이터를 이용하여 그 단면을디스플레이하는 장치 및 방법

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1891159A (zh) * 2005-07-01 2007-01-10 深圳迈瑞生物医疗电子股份有限公司 用于三维超声成像的快速扫描变换方法
CN101292883A (zh) * 2007-04-23 2008-10-29 深圳迈瑞生物医疗电子股份有限公司 超声三维快速成像方法及其装置
CN102106741A (zh) * 2009-12-25 2011-06-29 东软飞利浦医疗设备系统有限责任公司 一种二维超声图像的三维重建方法
CN104094319A (zh) * 2012-01-19 2014-10-08 株式会社东芝 图像处理设备、立体图像显示设备和图像处理方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Three-dimensional ultrasound imaging;Aaron Fenster等;《PHYSICS IN MEDICINE AND BIOLOGY》;20011231;全文 *
一种扇形扫描三维超声成像系统三维重建方法;杨金耀等;《中国医疗器械信息》;20121231;全文 *

Also Published As

Publication number Publication date
CN106600695A (zh) 2017-04-26

Similar Documents

Publication Publication Date Title
Åkerblom et al. Analysis of geometric primitives in quantitative structure models of tree stems
CN101441780B (zh) 三维网格模型的剖切方法
Choi et al. Continuous collision detection for ellipsoids
US9916684B2 (en) Generating portable three-dimensional print-preview renderings of three-dimensional models
CN109613540A (zh) 一种基于WebGL的多普勒天气雷达三维可视化方法
US20100271369A1 (en) Curved surface rendering system and method
Surazhsky et al. Intrinsic morphing of compatible triangulations
CN106600695B (zh) 一种三维体重构方法和系统
Kim et al. Efficient offset trimming for planar rational curves using biarc trees
CN107053677B (zh) 一种面向3d打印的简单多边形多层旋转体生成方法及系统
Ajith Kumar et al. Near perfect path generation of corners chamfered rectangle and single synthesis cam-link mechanism to generate special-slot path
CN113643423B (zh) 一种具有轴对称特性的Waterbomb衍生折纸结构造型方法及折纸结构
JP5044144B2 (ja) 医用画像形成装置
JP3928016B2 (ja) 最大角度法を用いた三角形メッシュ生成方法及びプログラム
You et al. Boundary constrained swept surfaces for modelling and animation
EP2372659B1 (en) Method and system for rapidly vectorizing image by gradient meshes based on parameterization
Qiu et al. C 1 triangular Coons surface construction and image interpolation based on new side-side and side-vertex interpolation operators
US20040233193A1 (en) Method for visualising a spatially resolved data set using an illumination model
Au A simple algorithm for medial axis transform computation
Farouki et al. Computational methods for rapid prototyping of analytic solid models
Andersen et al. Height and tilt geometric texture
Kjenstad Construction and computation of geometries on the ellipsoid
Rypl Sweeping of unstructured meshes over generalized extruded volumes
EP3806028B1 (en) Ultrasound image display method and apparatus, and storage medium and electronic apparatus
Zhuo et al. Fleshing: Spine‐driven Bending with Local Volume Preservation

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