CN108364326A - 一种ct成像方法 - Google Patents
一种ct成像方法 Download PDFInfo
- Publication number
- CN108364326A CN108364326A CN201810128432.9A CN201810128432A CN108364326A CN 108364326 A CN108364326 A CN 108364326A CN 201810128432 A CN201810128432 A CN 201810128432A CN 108364326 A CN108364326 A CN 108364326A
- Authority
- CN
- China
- Prior art keywords
- point
- parameter
- straight line
- quasi
- solution
- 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
- 238000013170 computed tomography imaging Methods 0.000 title claims abstract description 34
- 238000000034 method Methods 0.000 title claims abstract description 31
- 238000003384 imaging method Methods 0.000 claims abstract description 44
- 238000004364 calculation method Methods 0.000 claims abstract description 37
- 230000000875 corresponding effect Effects 0.000 claims description 10
- 238000010521 absorption reaction Methods 0.000 claims description 9
- 238000006243 chemical reaction Methods 0.000 claims description 9
- 230000005684 electric field Effects 0.000 claims description 6
- 238000005259 measurement Methods 0.000 claims description 5
- 238000012512 characterization method Methods 0.000 claims description 3
- 230000002596 correlated effect Effects 0.000 claims description 3
- 230000000704 physical effect Effects 0.000 claims description 3
- 238000012805 post-processing Methods 0.000 claims description 2
- 238000002247 constant time method Methods 0.000 abstract description 4
- 238000010276 construction Methods 0.000 abstract description 4
- 230000001186 cumulative effect Effects 0.000 abstract description 3
- 238000004422 calculation algorithm Methods 0.000 description 34
- 238000010586 diagram Methods 0.000 description 7
- 238000005516 engineering process Methods 0.000 description 7
- 238000005553 drilling Methods 0.000 description 6
- 238000010291 electrical method Methods 0.000 description 4
- 241001212149 Cathetus Species 0.000 description 3
- 238000001514 detection method Methods 0.000 description 3
- 239000003814 drug Substances 0.000 description 3
- 239000011159 matrix material Substances 0.000 description 3
- 238000009659 non-destructive testing Methods 0.000 description 3
- 241000208340 Araliaceae Species 0.000 description 2
- 241001269238 Data Species 0.000 description 2
- 235000005035 Panax pseudoginseng ssp. pseudoginseng Nutrition 0.000 description 2
- 235000003140 Panax quinquefolius Nutrition 0.000 description 2
- 238000002591 computed tomography Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 235000008434 ginseng Nutrition 0.000 description 2
- 241001270131 Agaricus moelleri Species 0.000 description 1
- 102000011990 Sirtuin Human genes 0.000 description 1
- 108050002485 Sirtuin Proteins 0.000 description 1
- 230000004308 accommodation Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000005611 electricity Effects 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000012067 mathematical method Methods 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000000149 penetrating effect Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000003325 tomography Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
- G06T11/005—Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2211/00—Image generation
Abstract
一种CT成像方法,包括:建立覆盖被探测区域、源点和观测点的坐标系;获取源点和观测点坐标设为参数1,计算源点‑观测点直线方程设为参数2,整理每对源点‑观测点直线上的测量计算结果设为参数3;按参数1、参数2、参数3建立全部参数数据库;确定解算区域,将解算区域离散成一个个拟解算像点,获取每个像点的坐标参数4;分别取每一个拟解算像点,设定录取规则录取直线,计算每个像点到录取直线的距离Ri同时记录录取直线总数M;按距离反比原则,以Ri为关联的权重参数作用在该像点对应的参数3上进行累加后用M值平均,作为该像点的参数值5;以参数值5进行全域成像。本发明无需用传统CT方法建设和求解方程组,无需迭代计算,应用领域广泛。
Description
技术领域
本发明属于计算机成像技术领域,尤其涉及一种CT成像方法,可广泛应用于电磁波CT、超声波CT、传导电场CT、地震波CT等技术采集数据的成像,适用于医学、工业、安检、无损检测和地球物理勘探等众多领域。
背景技术
CT(ComputedTomography,简称CT)是“计算机断层摄影术”的简称,由HounsfieldG.N.1969年设计成功,1972年问世。CT是用X线束从多个方向对人体检查部位具有一定厚度的层面进行扫描,由探测器而不用胶片接收透过该层面的X线,转变为可见光后,由光电转换器转变为电信号,再经模拟/数字转换器转为数字,输入计算机处理。图像处理时将选定层面分成若干个体积相同的立方体,称之为体素(voxel)。扫描所得数据经计算而获得每个体素的X线衰减系数或称吸收系数,再排列成矩阵,即构成数字矩阵。数字矩阵中的每个数字经数字/模拟转换器转为由黑到白不等灰度的小方块,称之为像素(pixel),并按原有矩阵顺序排列,即构成CT图像。所以,CT图像是由一定数目像素组成的灰阶图像或RGB图像,是数字图像,是重建的断层图像。每个体素X线吸收系数可通过不同的数学方法算出。
由于现实的CT探测无法直接得到充分完备的投影数据,即观测得到的投影数据的数量及结构无法达到解析成像技术要求的理想条件。在CT成像过程中,为得到相对清晰的图像,一般采用解析重建算法和迭代重建算法。
解析重建算法对于二维平行束、扇束CT成像,常用的基于中心切片定理(归功于Bracewell)的滤波反投影重建算法;对于三维锥束CT重建,目前工程上应用最多是近似重建算法,其中典型的是FDK近似重建算法。迭代重建算法又分为代数重建算法和统计重建算法,代数重建中典型的算法有ART算法、SIRT算法、SART算法、MART算法等;统计迭代重建中,典型的算法有OSEM算法、MLEM算法、MAP算法等。实际工程应用中前者是主流,后者尽管仿真结果和试验结果是有效的,但需要大量的迭代计算,花费时间较长,无法满足实时性的要求。目前,CT成像算法的精度和快速性仍有待提升,成像算法的改进有助于提高CT成像装备的效能。
CT技术在源和接收传感器变通后又诞生了在工程无损检测中广泛应用的超声波CT、电磁波CT以及在地球物理勘探中广泛应用的弹性波CT(又称地震波CT)和传导类电法CT。这些CT方法总体上沿用X射线CT成像的基本原理,需要网格剖分拟成像的目标,并建立关于网格点(或体素)吸收系数、慢度(速度的导数)、电阻率等参数的大型方程组,通过最优化算法来实现每一个网格点相关参数的求解。运算量巨大,计算结果带来许多噪声干扰,需要后续算法对图像进行精细的处理。
目前,CT成像硬件发展较快,但算法研究进展较慢,也制约了相关领域的应用。以下是对CT成像算法改进的相关文献报道:
对比文件1:CN106846427A公开一种基于重加权各向异性全变分的有限角度CT重建方法,包括根据CT成像原理建立CT图像重建方程式,初始化各项参数,迭代求解数据保真项和重加权各向异性全变分极小化项,直至满足迭代停止条件。该对比文件提供的CT重建方法虽能提高重建效果,但是仍然需要通过迭代计算实现,还是需要花费较长的时间。
对比文件2:CN106037781A公开一种CT成像方法,该技术方案不对图像进行预重建,算法中也不涉及迭代步骤,只采用多个视角的投影图来定位被检测物体,实现对投影数据进行修复。在此基础上,对第一次修复后的数据进一步依据PET投影数据进行比对修复,而后再结合滤波手段和负对数运算校正实现了对CT图像的全面修正。该对比文件提供的CT成像方法,虽然算法中不涉及迭代步骤,能缩短计算时间,但是重建效果不佳。
本发明拟提供一种能广泛应用于包括电磁波CT、超声波CT、地震波CT、传导类电法CT在内的各种需要计算机成像领域的,且不会局限在断层扫描,特别适用于三维结构体的容积CT成像的方法。
发明内容
本发明要解决的技术问题是提供一种CT成像方法,该方法无需用传统CT方法建设和求解方程组,可广泛应用于电磁波CT、超声波CT、传导电场CT、地震波CT等技术采集数据的成像,适用于医学、工业、安检、无损检测和地球物理勘探等众多领域。
为了解决上述技术问题,本发明采用以下技术方案:
提供一种CT成像方法,具体包括如下步骤:
S1、建立一个覆盖被探测区域、源点和观测点的坐标系,被探测区域为小于或等于所有源点-观测点线段经过的面或容积;
S2、获取源点和观测点坐标,设为参数1;计算源点-观测点的直线方程,设为参数2;整理每对源点-观测点直线上的测量计算结果直接设为参数3,或者将测量计算结果进行转换后设为参数3;
S3、按参数1、参数2、参数3建立全部参数的数据库;
S4、确定拟解算区域,选定整个被探测区域作为拟解算区域或者在被探测区域内选定合适的某块区域作为拟解算区域,将拟解算区域离散成一个个独立的拟解算像点,并获取每个拟解算像点的坐标,设为参数4;
S5、分别取每一个拟解算像点,按照设定的录取规则录取源点-观测点直线,记为录取直线,并用点到直线的距离公式计算拟解算像点到各条录取直线的距离,设为Ri,同时记录属于该拟解算像点的所有录取直线的总数M,i=1~M;
S6、按照距离反比原则,即:距离越近贡献越大,距离越远贡献越小,以对应录取直线的Ri为关联的权重参数作用在该拟解算像点对应录取直线的参数3上进行累加后用M值平均,作为该拟解算像点的成像参数值5;
S7、串行或并行计算所有拟解算像点的成像参数值5,以这些参数值5进行全域成像。
进一步地,
还包括步骤S8、对步骤S7中以参数值5进行全域成像的结果根据专业领域或者用途进行后处理,获取更好表征被探测介质物性解释的图像结果。
进一步地,步骤S1中:
所述被探测区域为二维的剖面、断面或截面;或者为三维立体结构。
进一步地,步骤S1中:
当被探测区域为二维剖面、断面、截面时,所述坐标系为二维笛卡儿坐标系;步骤S2中对应的源点-观测点的直线方程为平面直线方程。
进一步地,
当被探测区域为三维立体结构时,所述坐标系为三维笛卡儿坐标系;步骤S2中对应的源点-观测点的直线方程为空间直线方程。
进一步地,步骤S1中:
所述坐标系为球坐标、柱坐标或其它坐标系,或为方便运算而特殊定义的坐标系。
进一步地,
当步骤S1中经过源点和观测点的直线方程也可以是曲线方程或者曲面方程,对应地,步骤S5中拟解算像点到直线的距离转化为点到曲线或者曲面的距离,拟解算像点的成像参数5的计算原则不变。
进一步地,步骤S2中,
对参数3进行转换包括:
对于电磁波CT,是将参数3进行专业换算后转换为衰减系数或吸收系数;
对于超声波CT或地震波CT,是将参数3进行专业换算后转换为波速或慢度;
对于传导类电场CT或其它低频电磁波场CT,是将参数3进行专业换算后转换为视电阻率或电导率。
进一步地,
步骤S5中,直线的录取规则,具体是:
设定以拟解算像点为中心的圆域或球域半径r,用点到直线的距离公式计算拟解算像点到各条源点-观测点直线的距离,如果距离小于或等于r,则录取该直线,并记为录取直线;
或者,设定以拟解算像点与源点连线为对称中心线的角域范围,录取位于角域范围内的源点-观测点直线,并记为录取直线。
步骤S5中直线的录取原则设定的原理在于:源点辐射到观测点的物理场,随着观测点的远离而存在扩散衰减和介质吸收衰竭,衰减速度很快。基于这个原理,可以认为远离拟解算像点的直线对成像参数5的影响不大,去掉一些影响不大的直线,同时参与计算的直线少,如此也可以节省计算资源。
进一步地,步骤S5中,
具体操作过程中,通过改变录取方式而放宽或减小直线录取率,使适当的结果参数3通过各自权重参与成像。
例如:通过加大或减小设定的圆域或球域半径r来提高或者降低直线的录取率,或者通过扩大或缩小设定的角域范围来提高或者降低直线的录取率。
进一步地,步骤S6中:
以距离Ri为关联的权重参数为:距离Ri的倒数,或者为距离Ri的其它正相关或负相关计算结果。
进一步地,
所述方法还包括:根据需要对步骤S6中的距离Ri的其它正相关或负相关计算结果进行归一化处理。
进一步地,
当观测系统需要源点与观测点互换时,直线方程保持不变,将互换后的测量数据作为该直线方程的另一个测量计算结果,两个测量结果取平均值作为新的参数3。
本发明的有益效果:
1、本发明的成像方法并不局限在断层扫描,而且尤其适用于三维结构体的容积CT成像。本方法为数据庞大的各种CT探测手段提供一种全新的思路,算法简单清晰,便于在计算机上实现并行计算;无需建立欠定的巨大线性方程组,能避免最优化算法的不利因素。包括X射线CT在内的各种成像问题,当观测装置几何形状、源点或观测点运动参数一定时,相关几何参数(源点及接收点坐标、拟解算区域内拟解算像点坐标乃至每个像点到各直线的距离等)都在数据库中事先清晰固定,从而可以大大节省计算时间。
2、本发明提供的CT成像算法无需用传统CT方法建设和求解方程组,不需要进行迭代计算,在医学CT、工业CT、安防CT以及地球物理勘探CT领域都有广阔的应用前景,既可以作为独立的成像算法对观测数据进行成像,也可以为其它更优秀的反演成像算法提供合理的初始模型。
3、本发明提供的泛CT成像算法,将源点和观测点坐标,设为参数1,计算源-观测点的直线方程,设为参数2,整理每对源-观测点的测量结果,设为参数3,并按参数1、参数2、参数3建立全部参数数据库,方便后期索引计算。如果观测系统相对固定,这些参数也基本固定,有利于节省计算资源。
4、本发明提供的成像方法,能广泛适用于电磁波CT、超声波CT、传导电场CT和地震波CT等领域的成像,行业应用前景广阔。
5、本发明提供的成像方法中,所述的权重表达方式或者测量计算结果(参数3)转换后的累加方式变化时,成像结果会有不同的针对性效果,在不同行业或不同探测对象的识别中可以分别使用。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1为本发明CT成像方法的流程示意图;
图2为本发明实施例1第四代人体X射线CT扫描方法及CT成像原理示意图;
图3为本发明实施例2和实施例3的CT成像原理示意图;
图4为本发明实施例4中地震CT成像原理示意图。
附图标记说明:
图2中:1-源点(X射线管);2-观测点(传感器接收单元);3-被探测待成像区域;4-源点-观测点直线;5-拟解算像点;图中的箭头代表源点的旋进方向;R1-R3为拟解算像点到源点-观测点直线的距离;
图3中:1’-源点;2’-观测点;3’-被探测区域;4’-源点-观测点直线;5’拟解算像点;S01-SN为源点编号;G01-GN为观测点编号;R1-R6为拟解算像点到源点-观测点直线的距离;
图4中:1”-源点(震源点);2”-观测点(地震接收点);3”-被探测区域;4”-源点-观测点直线;P为拟解算像点。
具体实施方式
下面结合实施例及附图对发明进一步说明,但不用来限制本发明的范围。
本发明提供一种CT成像方法,可广泛适用于电磁波CT、超声波CT、传导电场CT和地震波CT等成像领域。
应用实施例1人体医学X射线CT成像领域
如图2所示的第四代医学X射线扫描设备成像算法示意图,并参考图1。X射线管为移动源点1,环形阵列传感器接收单元为相对固定的一个个观测点2,构成一个圆环。X射线管1旋转一周(如图中箭头方向旋转)(多个源点)利用对应的部分阵列传感器完成对应源点个数的角域射线测量,测量计算参数为对应射线上的吸收系数,本实施例的成像方法如下:
S1、在圆环平面内以圆环中心为原点,建立覆盖被探测区域3、源点1和观测点2的笛卡儿坐标系;
S2、获取源点1和观测点2坐标,设为参数1,计算源点-观测点的直线方程,设为参数2,整理每对源点-观测点直线4上的测量计算结果,并将测量计算结果经专业换算后转换得到吸收系数,设为参数3;
S3、按参数1、参数2、参数3建立全部数据数据库;
S4、确定拟解算区域(被探测区域3),将拟解算区域离散成一个个独立的拟解算像点,点距为1mm,获取每个像点的坐标,设为参数4;
S5、分别对步骤S4中的每一个拟解算像点,设定每个X射线源点为角顶点,将角顶点与拟解算像点连线,以连线为中央线,两侧角域各15度范围内的数据库中直线为录取直线,并记录录取直线的总数M;并用点到直线的距离公式计算每个拟解算像点到录取直线的距离Ri,i=1~M;
S6、按照距离平方反比原则将M条录取直线上的测量计算结果Si,以1/Ri 2为权重参数采用加权求和平均的方式赋值到对应拟解算像点的吸收系数,即赋值作为该像点的成像用参数5;
S7、以这些参数值5进行全域吸收系数成像。
应用实施例2超声波CT或地震波CT成像领域
如图3所示的跨孔超声波CT或地震波CT无损检测成像模型示意图,并参考图1计算流程图。跨孔距离为24m,孔深200m。左侧观测孔中布设有众多声源点1’,声源点距离D=2m,接收传感器分布在右侧孔中作为一个个观测点2’,观测点点距为2m,整体构成一个覆盖被探测区域3’的矩形断面。本实施例的成像方法如下:
S1、在矩形断面上以左下角为原点0,建立覆盖被探测区域3’、源点1’和观测点2’的平面笛卡儿坐标系;
S2、获取源点1’和观测点2’坐标,设为参数1,计算源点-观测点的直线方程,设为参数2,整理每对源点-观测点直线4’上的测量计算结果,经专业换算后转换得到慢度,设为参数3,即平均波速的倒数;
S3、按参数1、参数2、参数3建立全部数据数据库;
S4、确定拟解算区域(被探测区域3’),将拟解算区域按像点距0.5m均匀离散成拟解算像点5’,共得到获49*401=19649个像点,获取每个像点的坐标参数4;
S5、分别取每一个拟解算像点5’,以拟解算像点5’为中心确定半径r=5m,若直线到拟解算像点5’的距离小于等于r,则录取该直线,保留该距离参数Ri,并提取该直线的参数3;录取总数M,i=1~M。
S6、按照距离反比原则将M条录取直线上的测量结果Ti(参数3)采用加权求和的方式赋值到对应拟解算像点的慢度,即赋值作为该像点的成像用参数5;
S7、以这些参数值5分别取倒数获得新的数据参数1/T,利用该参数进行全域速度成像。
应用实施例3对于传导类电法(地球物理勘探)电阻率CT成像领域
如图3所示的跨孔传导类电法CT地球物理勘探成像模型示意图,并参考图1所示成像计算流程图。跨孔距离为12m,孔深100m。左侧观测钻孔中布设有众多供电源点(源点1’),供电点点距为1m,接收电极(或传感器)分布在右侧钻孔中作为一个个观测点2’,观测点2’点距为1m;整体构成一个覆盖被探测区域3’的矩形断面。本实施例的成像方法如下:
S1、在矩形断面上以左下角为原点0,建立覆盖被探测区域3’、源点1’和观测点2’的笛卡儿坐标系;
S2、获取源点1’和观测点2’坐标,设为参数1,计算源点-观测点的直线方程,设为参数2,整理每对源点-观测点直线4’上的测量计算结果,经专业换算后得到视电阻率设为参数3;
S3、按参数1、参数2、参数3建立全部数据数据库;
S4、确定拟解算区域(被探测区域3’),将被解算矩形区域按像点点距0.25m进行均匀离散,共得到获49*401=19649个拟解算像点5’,获取每个像点的坐标参数4;
S5、分别取每一个拟解算像点5’,以拟解算像点5’纵坐标为参考值,数据库中各直线上源点1’和观测点2’纵坐标分别减去参考值之和的绝对值小于30米的直线为录取直线,并记录录取直线的总数M;用点到直线的距离公式计算每个拟解算像点5’到录取直线的距离Ri,i=1~M;
S6、按照距离反比原则将M条录取直线上的测量计算结果PSi(参数3)采用求和平均的累加方式赋值到对应像点,即赋值作为该像点的成像用电阻率参数5,这里Rmax为归一化参数,取值为100m,即钻孔深度。
S7、利用该参数5进行全域电阻率成像。
上述实施例1-3提供的是不同应用领域二维面的CT成像方法实施例,本发明算法还可以推广应用到三维空间结构体的CT成像,届时所述的坐标系为三维笛卡儿坐标系,所述的源点、观测点坐标都是三维坐标,所述的源点-观测点的直线为空间直线,拟成像点到直线的距离为空间中一点到空间中直线的距离,但是直线的录取原则以及相应的加权平均算法原则不会发生本质的改变。以下以实施例4为例做出具体说明:
实施例4:探测地下三维地质结构的地震波CT成像方法
如图4所示,水平地面上已完成四个深度均为30m正方形分布的钻孔ZK01、ZK02、ZK03、ZK04,正方形边长为15米。任选一对钻孔开展地震波CT探测作业,共有6对钻孔。一个做发射孔,发射源点(源点1”)点距1m,共31个点;一个做接收孔,接收点(观测点2”)点距1m,也是31个点。发射源点(源点1”)和接收点(观测点2”)同步或交错逐点移动,用地震仪进行地震波旅行时差记录,换算成该各条地震波旅行线段上的速度,用速度取倒数计算出指定连线上的慢度,作为参数3;6对钻孔测量完成后,得到6*31*31=5766个计算参数3(慢度),并分别属于各自源点-观测点直线(连线)4”。建立三维笛卡尔坐标系,原点0定为ZK01孔的底部,X方向指向ZK02,Y方向指向ZK04。于是,每个震源点(源点1”)和每个观测点2”的三维坐标可以求出作为参数1,每对源点-观测点的空间直线方程可以求出作为参数2;将参数1、参数2和参数3纳入全部参数数据库。
假设拟成像区域(被探测区域3”)为10m*10m*20m立方体,八个顶点的坐标分别为:(2.5,2.5,5)、(12.5,2.5,5)、(12.5,12.5,5)、(2.5,12.5,5)、(2.5,2.5,25)、(12.5,2.5,25)、(12.5,12.5,25)、(2.5,12.5,25),将该立方体均匀剖分成0.25m*0.25m*0.25m的拟解算像点,取任意一个拟解算像点中心的坐标P(x,y,z),按空间点到空间直线的距离公式可以计算出P点到所有观测连线的距离,设点Pi到观测连线的距离为Ri(i取值范围为1-5766),将所有P点到观测连线的距离中的极大值设为Rmax,则P点上的成像参数5为利用并行计算方式计算出所有拟解算像点的成像参数5;利用该参数5进行立方体介质的速度成像。
本发明还可以对上述实施例中获得的全域成像结果,根据专业领域或者用途进行后处理,以获取更好表征被探测介质物性解释的图像结果。
以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解,本发明不受上述实施例的限制,上述实施例和说明书中描述的只是说明本发明的原理,在不脱离本发明精神和范围的前提下本发明还会有各种变化和改进,这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等同物界定。
Claims (10)
1.一种CT成像方法,其特征在于,具体包括如下步骤:
S1、建立一个覆盖被探测区域、源点和观测点的坐标系,被探测区域为小于或等于所有源点-观测点线段经过的面或容积;
S2、获取源点和观测点坐标,设为参数1;计算源点-观测点的直线方程,设为参数2;整理每对源点-观测点直线上的测量计算结果直接设为参数3,或者将测量计算结果进行转换后设为参数3;
S3、按参数1、参数2、参数3建立全部参数的数据库;
S4、确定拟解算区域,将拟解算区域离散成一个个独立的拟解算像点,并获取每个拟解算像点的坐标,设为参数4;
S5、分别取每一个拟解算像点,设定录取规则录取源点-观测点直线,记为录取直线,并用点到直线的距离公式计算拟解算像点到各条录取直线的距离,设为Ri,同时记录属于该拟解算像点的所有录取直线的总数M,i=1~M;
S6、按照距离反比原则,以对应录取直线的Ri为关联的权重参数作用在该拟解算像点对应录取直线的参数3上进行累加后用M值平均,作为该拟解算像点的成像参数值5;
S7、串行或并行计算所有拟解算像点的成像参数值5,以这些参数值5进行全域成像。
2.根据权利要求1所述的一种CT成像方法,其特征在于,步骤S1中:
所述被探测区域为二维的剖面、断面或截面,对应地,所述的坐标系为二维笛卡儿坐标系,源点-观测点的直线方程为平面直线方程;
或者所述被探测区域为三维立体结构,对应地,所述坐标系为三维笛卡儿坐标系,源点-观测点的直线方程为空间直线方程。
3.根据权利要求1所述的一种CT成像方法,其特征在于,
步骤S1中经过源点和观测点的直线方程为曲线方程或者曲面方程,对应地,步骤S5中拟解算像点到直线的距离转化为点到曲线或者曲面的距离,拟解算像点的成像参数5的计算原则不变。
4.根据权利要求1-3任意一项所述的一种CT成像方法,其特征在于,
步骤S2中,对参数3进行转换包括:
对于电磁波CT,是将参数3进行专业换算后转换为衰减系数或吸收系数;
对于超声波CT或地震波CT,是将参数3进行专业换算后转换为波速或慢度;
对于传导类电场CT或其它低频电磁波场CT,是将参数3进行专业换算后转换为视电阻率或电导率。
5.根据权利要求1-3任意一项所述的一种CT成像方法,其特征在于,
步骤S5中,直线的录取规则,具体是:
设定以拟解算像点为中心的圆域或球域半径r,用点到直线的距离公式计算拟解算像点到各条源点-观测点直线的距离,如果距离小于或等于r,则录取该直线,并记为录取直线;
或者,设定以拟解算像点与源点连线为对称中心线的角域范围,录取位于角域范围内的源点-观测点直线,并记为录取直线。
6.根据权利要求5所述的一种CT成像方法,其特征在于,
通过加大设定的半径r来提高直线的录取率或者减小设定的半径r来降低直线的录取率;或者,通过扩大设定的角域范围来提高直线的录取率或缩小设定的角域范围来降低直线的录取率,使适当的结果参数3通过各自权重参与成像。
7.根据权利要求1-3任意一项所述的一种CT成像方法,其特征在于,
步骤S6中,以距离Ri为关联的权重参数为:距离Ri的倒数,或者为距离Ri的其它正相关或负相关计算结果。
8.根据权利要求7所述的一种CT成像方法,其特征在于,
根据需要进一步对步骤S6中的距离Ri的其它正相关或负相关计算结果进行归一化处理。
9.根据权利要求1-3任意一项所述的一种CT成像方法,其特征在于,
还包括步骤S8、对步骤S7中以参数值5进行全域成像的结果根据专业领域或者用途进行后处理,获取更好表征被探测介质物性解释的图像结果。
10.根据权利要求1-3任意一项所述的一种CT成像方法,其特征在于,
当观测系统需要源点与观测点互换时,直线方程保持不变,将互换后的测量数据作为该直线方程的另一个测量计算结果,两个测量计算结果取平均值作为新的参数3。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810128432.9A CN108364326B (zh) | 2018-02-08 | 2018-02-08 | 一种ct成像方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810128432.9A CN108364326B (zh) | 2018-02-08 | 2018-02-08 | 一种ct成像方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN108364326A true CN108364326A (zh) | 2018-08-03 |
CN108364326B CN108364326B (zh) | 2021-04-02 |
Family
ID=63004929
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810128432.9A Expired - Fee Related CN108364326B (zh) | 2018-02-08 | 2018-02-08 | 一种ct成像方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN108364326B (zh) |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110133644A (zh) * | 2019-06-03 | 2019-08-16 | 中南大学 | 基于插值尺度函数法的探地雷达三维正演方法 |
CN113960695A (zh) * | 2021-09-18 | 2022-01-21 | 山东大学 | 一种复杂城市环境富水岩溶精细探查方法 |
CN116468859A (zh) * | 2023-06-19 | 2023-07-21 | 之江实验室 | 适用于声速非均匀分布的三维光声断层成像方法、装置 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2002067201A1 (en) * | 2001-02-15 | 2002-08-29 | The Regents Of The University Of Michigan | Statistically reconstructing an x-ray computed tomography image with beam hardening corrector |
CN1722178A (zh) * | 2004-06-30 | 2006-01-18 | 通用电气公司 | 图像的三维重建方法和系统 |
CN1742296A (zh) * | 2003-02-05 | 2006-03-01 | 西门子公司 | 锥形光束计算机断层成像的立体像素驱动的螺旋重建方法 |
CN1751323A (zh) * | 2003-02-19 | 2006-03-22 | 皇家飞利浦电子股份有限公司 | 利用螺旋相对运动和锥形光束的计算机化x线断层摄影方法 |
CN101082991A (zh) * | 2006-05-31 | 2007-12-05 | 西门子公司 | 利用投影重建对象的方法以及执行该方法的装置 |
-
2018
- 2018-02-08 CN CN201810128432.9A patent/CN108364326B/zh not_active Expired - Fee Related
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2002067201A1 (en) * | 2001-02-15 | 2002-08-29 | The Regents Of The University Of Michigan | Statistically reconstructing an x-ray computed tomography image with beam hardening corrector |
CN1742296A (zh) * | 2003-02-05 | 2006-03-01 | 西门子公司 | 锥形光束计算机断层成像的立体像素驱动的螺旋重建方法 |
CN1751323A (zh) * | 2003-02-19 | 2006-03-22 | 皇家飞利浦电子股份有限公司 | 利用螺旋相对运动和锥形光束的计算机化x线断层摄影方法 |
CN1722178A (zh) * | 2004-06-30 | 2006-01-18 | 通用电气公司 | 图像的三维重建方法和系统 |
CN101082991A (zh) * | 2006-05-31 | 2007-12-05 | 西门子公司 | 利用投影重建对象的方法以及执行该方法的装置 |
Non-Patent Citations (1)
Title |
---|
刘征宇 等: "基于距离加权约束算法的围岩三维电阻率CT反演成像研究", 《岩土工程学报》 * |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110133644A (zh) * | 2019-06-03 | 2019-08-16 | 中南大学 | 基于插值尺度函数法的探地雷达三维正演方法 |
CN110133644B (zh) * | 2019-06-03 | 2023-03-31 | 中南大学 | 基于插值尺度函数法的探地雷达三维正演方法 |
CN113960695A (zh) * | 2021-09-18 | 2022-01-21 | 山东大学 | 一种复杂城市环境富水岩溶精细探查方法 |
CN116468859A (zh) * | 2023-06-19 | 2023-07-21 | 之江实验室 | 适用于声速非均匀分布的三维光声断层成像方法、装置 |
CN116468859B (zh) * | 2023-06-19 | 2023-09-08 | 之江实验室 | 适用于声速非均匀分布的三维光声断层成像方法、装置 |
Also Published As
Publication number | Publication date |
---|---|
CN108364326B (zh) | 2021-04-02 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Li et al. | Refraction corrected transmission ultrasound computed tomography for application in breast imaging | |
CN104809750B (zh) | 一种直线扫描ct系统及图像重建方法 | |
CN108577876A (zh) | 一种多边形静止ct及其工作方法 | |
CN108364326A (zh) | 一种ct成像方法 | |
WO2020206755A1 (zh) | 一种基于射线理论的超声ct图像重建方法及系统 | |
Fadakar Alghalandis et al. | The RANSAC method for generating fracture networks from micro-seismic event data | |
CN107328798A (zh) | 一种新型icl系统及实现方法 | |
CN106251381A (zh) | 图像重建方法 | |
CN107389793A (zh) | 一种铝合金多孔结构孔边腐蚀损伤监测方法 | |
CN101673413A (zh) | 基于表面数据的植物枝体三维形态建模系统及方法 | |
CN107392977A (zh) | 单视图切伦科夫发光断层成像重建方法 | |
CN102749053A (zh) | 基于三维可视化和蒙特卡罗方法的体积测量方法 | |
Lyubimov et al. | Application of the photon average trajectories method to real-time reconstruction of tissue inhomogeneities in diffuse optical tomography of strongly scattering media | |
CN103393410B (zh) | 一种基于交替迭代运算的荧光分子断层成像重建方法 | |
CN104268938A (zh) | 一种三维温度场的重建方法 | |
Zhang et al. | Fast projection algorithm for LIM-based simultaneous algebraic reconstruction technique and its parallel implementation on GPU | |
CN103908278A (zh) | 图像重建方法及装置、医疗成像系统 | |
CN107942326A (zh) | 一种具有高通用性的二维主动毫米波成像方法 | |
CN103620393A (zh) | 成像装置 | |
CN109598769A (zh) | 总变差正则化约束的超声成像同步代数迭代重建方法 | |
CN101639453B (zh) | 一种针对多相流动系统中ct射线硬化的校正方法 | |
Fan et al. | Sorted L1 regularization method for damage detection based on electrical impedance tomography | |
CN109102552A (zh) | 一种非均匀形状约束的像素值域滤波超声成像重建方法 | |
Cai et al. | Three-dimensional combustion diagnostics based on computed tomography of chemiluminescence | |
CN114741951A (zh) | 一种基于卷积神经网络的介质目标电磁探测方法 |
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 | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20210402 |
|
CF01 | Termination of patent right due to non-payment of annual fee |