CN110428370A - 一种利用偏心旋转提高锥形束spect成像分辨率的方法 - Google Patents

一种利用偏心旋转提高锥形束spect成像分辨率的方法 Download PDF

Info

Publication number
CN110428370A
CN110428370A CN201910588491.9A CN201910588491A CN110428370A CN 110428370 A CN110428370 A CN 110428370A CN 201910588491 A CN201910588491 A CN 201910588491A CN 110428370 A CN110428370 A CN 110428370A
Authority
CN
China
Prior art keywords
resolution
detector
image
under
projected image
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
Application number
CN201910588491.9A
Other languages
English (en)
Other versions
CN110428370B (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.)
Beijing University of Technology
Beijing Institute of Technology BIT
Original Assignee
Beijing University of Technology
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 Beijing University of Technology filed Critical Beijing University of Technology
Priority to CN201910588491.9A priority Critical patent/CN110428370B/zh
Publication of CN110428370A publication Critical patent/CN110428370A/zh
Application granted granted Critical
Publication of CN110428370B publication Critical patent/CN110428370B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/02Arrangements for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03Computed tomography [CT]
    • A61B6/037Emission tomography
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/52Devices using data or image processing specially adapted for radiation diagnosis
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00Geometric image transformations in the plane of the image
    • G06T3/60Rotation of whole images or parts thereof
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/90Dynamic range modification of images or parts thereof
    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16HHEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H30/00ICT specially adapted for the handling or processing of medical images
    • G16H30/40ICT specially adapted for the handling or processing of medical images for processing medical images, e.g. editing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10108Single photon emission computed tomography [SPECT]

Landscapes

  • Engineering & Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Medical Informatics (AREA)
  • Physics & Mathematics (AREA)
  • Public Health (AREA)
  • General Health & Medical Sciences (AREA)
  • Radiology & Medical Imaging (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Pathology (AREA)
  • Animal Behavior & Ethology (AREA)
  • Optics & Photonics (AREA)
  • Biophysics (AREA)
  • Veterinary Medicine (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Epidemiology (AREA)
  • Primary Health Care (AREA)
  • Nuclear Medicine (AREA)
  • Image Processing (AREA)

Abstract

本发明涉及一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,包括如下步骤:步骤一,获得N个旋转角度下的低分辨率投影图像,在每个旋转角度下获得若干幅低分辨率投影图像,所述低分辨率投影图像的采集要保证检测器焦点始终为同一点,所述低分辨率投影图像与该旋转角度下高分辨率投影图像之间存在几何关系;步骤二,对单个旋转角度下的若干低分辨率投影图像进行处理,得到该旋转角度下对应高分辨率投影图像;步骤三,重复步骤二操作,得到N个旋转角度下的高分辨率投影图像,最终得到真实高分辨率图像。本发明在不改变检测器硬件的前提下,能够提高锥形束SPECT图像的分辨率。

Description

一种利用偏心旋转提高锥形束SPECT成像分辨率的方法
技术领域
本发明涉及一种SPECT成像方式,特别涉及一种利用偏心旋转提高锥形束 SPECT成像分辨率的方法,属于医学图像重建领域。
背景技术
单光子发射计算机断层成像(single photon emission computed tomography,SPECT)是核医学成像中的一项重要技术,放射性示踪剂通过辐射γ光子,穿过人体后被检测器检测,最终重建得到断层图像,在肿瘤诊断、心血管疾病诊断和脑神经科学研究等领域得到广泛应用。与其他医学成像方式(如X射线计算机断层扫描(CT)和磁共振成像(MRI))相比,SPECT具有分辨率较低的局限性。在SPECT中,限制图像空间分辨率和灵敏度的关键因素为检测器的准直孔径的大小,准直器孔径过大或过小都会影响成像质量的好坏。
超分辨率重建是利用同一场景的一幅或多幅低分辨率图像的附加信息重建出一幅或多幅高分辨率图像的技术,作为一种有效的提高图像分辨率的方法,广泛应用于遥感卫星,医学图像等领域。近年来超分辨率重建被应用到CT、PET、 SPECT、磁共振成像领域,以提高图像分辨率,这些方法主要应用在图像域或平行及扇形检测器投影重建领域。
锥形投影SPECT***与平行束、扇形束SPECT相比,具有高分辨率的特点,有利于对小目标进行检测及重建,利用检测器的旋转移动等操作,将锥形束投影重建和超分辨率重建相结合,将会极大提高SPECT的成像分辨率。
发明内容
本发明提供一种利用偏心旋转提高锥形束SPECT成像分辨率的方法。该方法在不改变检测器中准直器孔径大小、且保证每个旋转角度下检测器焦点始终为同一点的前提下,利用低分辨率检测器的偏心旋转,提高SPECT图像的分辨率。
本发明在每个角度下采集的低分辨率投影图像个数为M,与实际低分辨率投影图像大小,目标高分辨率投影图像大小有关;若低分辨率检测器大小为Sl×Sl个像素,目标高分辨率投影图像大小为Sh×Sh个像素,R为高分辨率图像大小与低分辨率图像大小像素之比,则M、N、R、Sh、Sl可由公式(1)、(2)、(3)给出。
N≥Sh (1)
M≥R×R (3)
本发明是通过下述技术方案实现的。
本发明公开的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,包括如下步骤:
步骤一,获得N个旋转角度下的低分辨率投影图像,在每个旋转角度下获得若干幅低分辨率投影图像,所述每个旋转角度下低分辨率投影图像的采集要保证检测器焦点始终为同一点,所述低分辨率图像与该旋转角度下高分辨率图像之间存在几何关系;
低分辨率投影图像采集,通过如下方法实现;
低分辨率检测器绕物体中心进行旋转,每次旋转角度,N为检测器绕物体旋转的次数;在每个旋转角度下,保证检测器焦点始终为同一点的前提下,对检测器进行M次偏心旋转;
采集得到N次旋转的每个位置处的M个低分辨率投影图像;记录每次偏心旋转的距离与角度;
记未偏心旋转前的初始检测器平面为基准平面,延长检测器焦点与偏心旋转后的检测器中心的连线,使连线与基准平面相交于一点,该点为交点;
所述检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的M 个交点,与基准平面中心点的距离为高分辨率像素大小的任意倍数,M个距离可以相等,也可以不相等,但要记录每次偏心旋转的距离;
所述检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的M 个交点,可以均匀分布,也可以不均匀分布,但要记录每个交点的位置坐标;
所述单个旋转角度下,单次偏心旋转得到的低分辨率投影图像与高分辨率投影图像之间存在几何关系:
设检测器焦点与第m次偏心旋转后的检测器中心连线的延长线与基准平面的交点和基准平面中心的距离为d(m),焦点到检测器中心的距离(即焦距)为 sid,高分辨率投影图像在基准平面内,以高分辨率投影图像中心为高分辨率投影图像坐标原点,以偏心旋转后低分辨率检测器中心为低分辨率投影图像坐标原点,对于低分辨率检测器中的任意像素点,取像素边界的四个点坐标为(i,j), (i,j+1),(i+1,j),(i+1,j+1),可得到其在高分辨率投影图像中对应区域块S*的边界点坐标,以坐标(i,j)为例;
φ(j)=θ+δ(j) (6)
g(j)=sid×tan(φ(j)) (7)
其中,θ为基准平面焦距与偏心旋转后检测器焦距之间的夹角,δ(j)为低分辨率投影图像中坐标j与偏心旋转后焦距之间的夹角,φ(j)为低分辨率投影图像中坐标j与基准平面焦距之间的夹角,g(j)为低分辨率投影图像中坐标j映射到高分辨率投影图像中的坐标,g(i)为低分辨率投影图像中坐标i映射到高分辨率投影图像中的坐标;
对于高分辨率投影图像中的任意像素点(a,b),判断其是否属于(g(i),g(j)), (g(i),g(j+1)),(g(i+1),g(j)),(g(i+1),g(j+1))组成的区域块S*内,若满足则对应 w(a,b)=1,这里w(a,b)表示S*内各高分辨率像素点对低分辨率像素点(a,b)的系数,由此得到低分辨率投影图像中任一像素点与高分辨率投影图像各像素点之间的映射关系;
步骤二、对单个旋转角度下的M个低分辨率投影图像进行处理,得到该旋转角度下对应高分辨率投影图像;
所述对低分辨率图像进行处理的方法包括如下步骤:
记PH0为初始假设高分辨率投影图像,即第1个旋转角度下的初始高分辨率投影图像,高分辨率投影图像设在基准平面内;PHnm为第n个旋转角度下,在检测器第m次偏心旋转位置对应的初始高分辨率投影图像,其中:n的取值范围为1-N,m的取值范围为1-M;PHnm1为第n个旋转角度下,该旋转角度下的初始高分辨率投影图像PHnm将其坐标系绕图像中心旋转γ角度得到的高分辨率投影图像,所述γ角度是第m个交点与基准平面中心点的连线和第一个交点与基准平面中心点的连线形成的夹角;
步骤2.1:在第n个旋转角度下,PHnm1可由式(9)求得:
其中:(x,y)为PHnm中像素点的坐标,(x1,y1)为(x,y)旋转坐标系后对应像素点位置;并由双线性插值公式求得PHnm1插值后每个像素点的坐标;
其中,(x0,y0)为(x1,y1)向下取整的第一个整数值坐标,(x2,y2)为(x1,y1)插值后PHnm1的像素点坐标值;
步骤2.2:记为第n个旋转角度下第m个低分辨率检测器测量得到的真实投影图像,PLnm为第n个旋转角度下对PHnm1进行按几何关系降采样得到的低分辨率投影图像,即将每个低分辨率投影图像像素点对应的高分辨率区域块S*的所有像素点按照公式(11)叠加:
其中,(x3,y3)为PLnm中各像素点的坐标,其中低分辨率图像中每个点对应高分辨率图像中的区域块S*,w(x2,y2)为S*内各高分辨率像素点对低分辨率像素点 (x3,y3)的系数;
将式(11)写为完整的向量形式即:
pLnm=WpHnm1 (12)
其中,pLnm表示低分辨率投影图像的列向量表示,长度为(Sl×Sl)×1,pHnm1为高分辨率投影图像的列向量表示,长度为(Sh×Sh)×1,W表示低分辨率投影图像中像素点与高分辨率投影图像中像素点的关系矩阵,其大小为(Sl×Sl)×(Sh×Sh);
步骤2.3:将步骤2.2得到的第n个角度下对第m个低分辨率投影图像PLnm与步骤一直接采集到的第n个角度下低分辨率图像按照公式(13)进行比较,求得真实低分辨率投影与目前低分辨率投影的差值ΔPLnm
步骤2.4:将经步骤2.3得到的对应像素点的差值ΔPLnm写为向量形式,即ΔpLnm,利用关系矩阵W,乘以每次调节的步长作为调节权重,得到对高分辨率图像pHnm1中每个点的调节权重ΔpHnm1
ΔpHnm1=WTΔpLnm×step (14)
其中:step为每次的调节步长,调节步长在0到1之间取值;
步骤2.5:由步骤2.4得到的调节权重ΔpHnm1还原为矩阵形式ΔPHnm1,并将其坐标系绕图像中心反向旋转γ度,由公式(15):
其中(x4,y4)为ΔPHnm1中像素点的坐标,(x5,y5)为(x4,y4)反向旋转坐标系后对应像素点的位置,并根据二维线性插值,得到PHnm的调节权重ΔPHnm
其中,(x50,y50)为(x5,y5)向下取整的第一个整数值坐标,(x6,y6)为(x5,y5)插值后ΔPHnm的像素点坐标值;
步骤2.6:将调节权重ΔPHnm叠加到在检测器第m次偏心旋转位置的初始高分辨率投影数据PHnm上,得到第(m+1)次偏心旋转位置的初始高分辨率投影数据 PHn(m+1),由公式(17)更新PHn(m+1)
PHn(m+1)=PHnm+ΔPHnm (17)
步骤2.7;此时得到第n个旋转角度下,第m次偏心旋转位置的高分辨率投影图像,即第m+1次偏心旋转位置的初始高分辨率投影图像;重复步骤2、1-2、 6,即能够得到第n个旋转角度下,经过M次处理后的高分辨率图像PHnM
此时完成全部M次调节;
再重复步骤2、1—2、6,继续从第n个旋转角度下的第1个位置开始调节图像,直到达到要求的迭代结束条件,得到满足投影分辨率要求的高分辨率投影图像,结束迭代过程,即得到检测器绕物体旋转时单个旋转角度下的高分辨率投影图像;
所述的迭代结束条件为:计算前后两次迭代过程高分辨率投影图像的差值的二范数,并与结束迭代阈值条件进行比较,从而确定是否结束迭代;所述的迭代阈值条件根据需要获得的高分辨率投影图像精度而定;
步骤三、重复步骤二操作,得到N个角度下的高分辨率投影图像,最终得到真实高分辨率图像;
有益效果:
1、本发明公开的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,由于每个旋转角度下低分辨率投影图像的采集保证检测器焦点始终为同一点,采集的多幅低分辨率投影图像都对同一目标进行成像,由此能够根据步骤二、三的处理方法从每个旋转角度下多幅低分辨率投影图像来获取该角度下的高分辨率投影图像,利用重建方法即可获得相应的高分辨率重建图,提高SPECT图像的分辨率;
2、本发明公开的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,利用检测器偏心旋转的方式,可以得到更多包含不同数据的低分辨率投影图像,尤其是得到更多图像中心点区域的数据,可以更好地提高整体成像质量。
附图说明
图1为锥形SPECT检测器模型示意图;
图2为锥形检测器绕物体中心旋转采集示意图,实线为检测器初始位置,虚线为检测器绕物体旋转一定角度到达的位置,检测器会绕着图中的虚线圆旋转一周;
图3为在单个旋转角度下,低分辨率检测器偏心旋转示意图,以圆形偏心旋转为例;
图4为在单个旋转角度下,检测器焦点与偏心旋转后的检测器中心连线的延长线与高分辨率投影图像(基准平面)的M个交点示意图,以圆形偏心旋转为例;
图5为在单个旋转角度下,高分辨率投影图像的坐标系绕图像中心旋转γ角度,便于坐标计算;
图6为在单个旋转角度下,低分辨率检测器偏心旋转示意图,以及低分辨率投影图像中任意像素点与高分辨率投影图像中像素点区域块的几何关系,以圆形偏心旋转为例;
图7为在单个旋转角度下,由多幅低分辨率投影图像获得一幅高分辨率投影图像流程图;
图8为由低分辨率投影图像获得高分辨率投影图像时单次迭代调整过程示意图,(a)为假设高分辨率投影图像,(b)为高分辨率投影图像变换坐标系后投影图像,(c)为由高分辨率投影图像(b)经过关系矩阵降采样后得到的假设低分辨率投影图像,(d)为采集到的实际低分辨率投影图像,(e)为由(c)和 (d)相减的差值得到的图像,将得到的(e)作为调节依据,对高分辨率投影图像(b)进行调节,得到高分辨率投影图像(f),(g)为将(f)反向变换坐标系,将图像移回初始位置,最终得到调节后的高分辨率投影图像(h),然后获得的(h)进入下一次调节过程;
图9为不同层原始图像、模拟高分辨率投影经FDK重建图像、由低分辨率投影采用本专利方法获得高分辨率投影经FDK重建图像、低分辨率投影经FDK 重建图像;其中检测器中心位于物体44层,检测器焦点到检测器中心的距离为 380个像素,焦点到物体中心的距离为310个像素;,数据均以圆形偏心旋转采集为例。
具体实施方式
下面结合附图和具体实施例对本发明方法的实施方式做详细说明。
实施例1:
检测器绕物体旋转的同时,在每个旋转角度下,检测器需进行偏心旋转操作,得到若干低分辨率投影图像。本实例以低分辨率检测器分辨率为32×32个像素,目标获得的高分辨率像素为128×128,N=128,R=4,M=16为例进行说明,但是并不说明该发明仅限于这些条件,本发明中所涉及的范围仅由权利说明书限定。
本实施例公开的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,其具体实施步骤如下:
步骤一,获得128个旋转角度下的低分辨率投影图像,在每个旋转角度下获得16幅低分辨率投影图像,所述每个旋转角度下低分辨率投影图像的采集要保证检测器焦点始终为同一点,所述低分辨率图像与该旋转角度下高分辨率图像之间存在几何关系;
低分辨率投影图像采集采用圆形偏心旋转的方式,通过如下方法实现;
低分辨率检测器绕物体中心进行旋转,每次旋转角度,128为检测器绕物体旋转的次数;图1为锥形束SPECT检测器模型示意图,锥形束SPECT检测器为一平板,且其所有准直器孔对应一个相同的焦点。被检测物***于检测器及检测器焦点中间,从物体发出的伽马射线可以被检测器接收。图2为检测器绕物体旋转进行光子采集示意图,实线为检测器初始位置,虚线为检测器绕物体旋转一定角度到达的位置,检测器会绕着图中的虚线圆旋转一周;
在每个旋转角度下,保证检测器焦点始终为同一点的前提下,对检测器进行16次偏心旋转;
记未偏心旋转前的初始检测器平面为基准平面,延长检测器焦点与偏心旋转后的检测器中心的连线,使连线与基准平面相交于一点,该点为交点,每个旋转角度下共16个交点,交点与基准平面中心点的距离不变,即构成圆形偏心,如图3所示;
采集得到128次旋转角度下的低分辨率投影图像,每个旋转角度下有16幅低分辨率投影图像;记录每次偏心转的距离与角度;
所述检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的16 个交点,与基准平面中心点的距离为一个高分辨率像素大小;
所述检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的16 个交点为均匀分布,如图3所示,图中的小黑点为16个交点;
所述单个旋转角度下,单次偏心旋转得到的低分辨率投影图像与高分辨率投影图像之间存在几何关系如下:
设检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的交点和基准平面中心的距离为d,d=1px,焦点到检测器中心的距离(即焦距)为sid, sid=380px,高分辨率投影图像在基准平面内,以高分辨率投影图像中心为高分辨率投影图像坐标原点,以偏心旋转后低分辨率检测器中心为低分辨率投影图像坐标原点,对于低分辨率检测器中的任意像素点,取像素边界的四个点坐标为(i,j),(i,j+1),(i+1,j),(i+1,j+1),可得到其在高分辨率投影图像中对应区域块S*的边界点坐标;
对于高分辨率投影图像中的任意像素点(a,b),判断其是否属于(g(i),g(j)), (g(i),g(j+1)),(g(i+1),g(j)),(g(i+1),g(j+1))组成的区域块S*内,若满足则对应 w(a,b)=1,这里w(a,b)表示S*内各高分辨率像素点对低分辨率像素点(a,b)的系数,由此得到低分辨率投影图像中任一像素点与高分辨率投影图像各像素点之间的映射关系;
步骤二、对单个旋转角度下的16幅低分辨率投影图像进行处理,得到该旋转角度下对应高分辨率投影图像;
所述对低分辨率图像进行处理的方法包括如下步骤:
记PH0为初始假设高分辨率投影图像,即第1个旋转角度下的初始高分辨率投影图像,高分辨率投影图像设在基准平面内;PHnm为第n个旋转角度下,在检测器第m次偏心旋转位置对应的初始高分辨率投影图像,其中:n的取值范围为1-128,m的取值范围为1-16;PHnm1为第n个旋转角度下,该旋转角度下的初始高分辨率投影图像PHnm将其坐标系绕图像中心旋转γ角度得到的高分辨率投影图像,所述γ角度是第m个交点与基准平面中心点的连线和第一个交点与基准平面中心点的连线形成的夹角,
步骤2.1:在第n个旋转角度下,为便于坐标计算,将PHnm坐标系绕其图像中心旋转γ度,并根据二维线性插值得到PHnm1,如图5所示;
步骤2.2:记为第n个旋转角度下第m个低分辨率检测器测量得到的真实投影图像,PLnm为第n个旋转角度下对PHnm1进行按几何关系降采样得到的低分辨率投影图像,并记录关系矩阵W,低分辨率投影图像中像素点与高分辨率投影图像像素点区域块的对应关系如图6所示;
步骤2.3:将步骤2.2得到的第n个角度下对第m个低分辨率投影图像PLnm与步骤一直接采集到的第n个角度下低分辨率图像按照公式(13)进行比较,求得真实低分辨率投影与目前低分辨率投影的差值ΔPLnm
步骤2.4:将经步骤2.3得到的对应像素点的差值ΔPLnm写为向量形式,即ΔpLnm,利用关系矩阵W,乘以每次调节的步长作为调节权重,得到对高分辨率图像pHnm1中每个点的调节权重ΔpHnm1
步骤2.5:由步骤2.4得到的调节权重ΔpHnm1还原为矩阵形式ΔPHnm1,并将其坐标系绕图像中心反向旋转γ度,并根据二维线性插值,得到PHnm的调节权重ΔPHnm
步骤2.6:将调节权重ΔPHnm叠加到在检测器第m次偏心旋转位置的初始高分辨率投影数据PHnm上,得到第(m+1)次偏心旋转位置的初始高分辨率投影数据 PHn(m+1),更新PHn(m+1)
步骤2.7;此时得到第n个旋转角度下,第m次偏心旋转位置的高分辨率投影图像,即第m+1次偏心旋转位置的初始高分辨率投影图像;重复步骤2、1-2、 6,即能够得到第n个旋转角度下,经过16次处理后的高分辨率图像PHnM
此时完成全部16次调节;
再重复步骤2、1—2、6,继续从第n个旋转角度下的第1个位置开始调节图像,直到达到要求的迭代结束条件,得到满足投影分辨率要求的高分辨率投影图像,结束迭代过程,即得到检测器绕物体旋转时单个旋转角度下的高分辨率投影数据;
图7为由16幅32×32像素的低分辨率投影图像获得1幅128×128像素高分辨率投影图像流程图,对应过程如图8所示;
流程图中初始高分辨率投影图像获得说明:
初始高分辨率投影图像获得选用下述两种方法之一:(1),对高分辨率投影图像随机赋值,此种方法在最终迭代后会有大量的随机噪声;(2),由采集到的低分辨率图像插值成高分辨率投影图像大小,并将此图像作为原始高分辨率投影图像,此种方法在最终迭代后获得高分辨率投影图像较好,且收敛速度较快;
迭代结束条件:迭代结束条件可有两种方法,这两种方法实际上是一样的, (1),通过计算图8中(a)与(h)的差值的二范数并与结束迭代阈值条件进行比较,从而确定是否结束迭代;(2),计算图8中(e)矩阵的二范数并与结束迭代阈值条件进行比较,从而确定是否结束迭代;
除圆形偏心旋转外,其他形式的偏心旋转,如椭圆偏心旋转,同样适用;
步骤三、重复步骤二操作,得到128个角度下的高分辨率投影图像,最终得到真实高分辨率图像;
根据步骤二、三得到的高分辨率投影根据FDK重建算法重建出高分辨率 SPECT图像。图9为利用FDK重建算法对投影图像重建所得结果。对比低分辨率投影直接重建结果,改善明显,与直接高分辨率投影直接重建结果基本接近,专利方法可有效提高分辨率。
该方法不仅可以用于SPECT超分重建,但该方法可以应用到其他超分重建算法中,如PET超分重建、CT超分重建等。
以上结合具体实施例对本发明的技术方案和具体实施方式作了说明,但这些说明不能被理解为限制了本发明的范围,这些仅是举例说明,可以对这些实施方式做出多种变更或修改,而不背离本发明的原理和实质。本发明的保护范围由随附的权利要求书限定,任何在本发明权利要求基础上的改动都是本发明的保护范围。

Claims (4)

1.一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,其特征在于:包括如下步骤:
步骤一,获得N个旋转角度下的低分辨率投影图像,在每个旋转角度下获得若干幅低分辨率投影图像,所述每个旋转角度下低分辨率投影图像的采集要保证检测器焦点始终为同一点,所述低分辨率图像与该旋转角度下高分辨率图像之间存在几何关系;
低分辨率投影图像采集,通过如下方法实现;
低分辨率检测器绕物体中心进行旋转,每次旋转角度,N为检测器绕物体旋转的次数;在每个旋转角度下,保证检测器焦点始终为同一点的前提下,对检测器进行M次偏心旋转;
采集得到N次旋转的每个位置处的M个低分辨率投影图像;记录每次偏心旋转的距离与角度;
步骤二、对单个旋转角度下的M个低分辨率投影图像进行处理,得到该旋转角度下对应高分辨率投影图像;
所述对低分辨率图像进行处理的方法包括如下步骤:
记未偏心旋转前的初始检测器平面为基准平面,延长检测器焦点与偏心旋转后的检测器中心的连线,使连线与基准平面相交于一点,该点为交点,每个旋转角度下共M个交点;
记PH0为初始假设高分辨率投影图像,即第1个旋转角度下的初始高分辨率投影图像,高分辨率投影图像设在基准平面内;PHnm为第n个旋转角度下,在检测器第m次偏心旋转位置对应的初始高分辨率投影图像,其中:n的取值范围为1-N,m的取值范围为1-M;PHnm1为第n个旋转角度下,该旋转角度下的初始高分辨率投影图像PHnm将其坐标系绕图像中心旋转γ角度后得到的高分辨率投影图像,所述γ角度是第m个交点与基准平面中心点的连线和第一个交点与基准平面中心点的连线形成的夹角;
步骤2.1:在第n个旋转角度下,PHnm1由式(9)求得:
其中:(x,y)为PHnm中像素点的坐标,(x1,y1)为(x,y)旋转坐标系后对应像素点位置;并由双线性插值公式求得PHnm1插值后每个像素点的坐标;
其中,(x0,y0)为(x1,y1)向下取整的第一个整数值坐标,(x2,y2)为(x1,y1)插值后PHnm1的像素点坐标值;
步骤2.2:记为第n个旋转角度下第m个低分辨率检测器测量得到的真实投影图像,PLnm为第n个旋转角度下对PHnm1进行按几何关系降采样得到的低分辨率投影图像,即将每个低分辨率投影图像像素点对应的高分辨率区域块S*的所有像素点按照公式(11)叠加:
其中,(x3,y3)为PLnm中各像素点的坐标,其中低分辨率图像中每个点对应高分辨率图像中的区域块S*,w(x2,y2)为S*内各高分辨率像素点对低分辨率像素点(x3,y3)的系数;
将式(11)写为完整的向量形式即:
pLnm=WpHnm1 (12)
其中,pLnm表示低分辨率投影图像的列向量表示,长度为(Sl×Sl)×1,pHnm1为高分辨率投影图像的列向量表示,长度为(Sh×Sh)×1,W表示低分辨率投影图像中像素点与高分辨率投影图像中像素点的关系矩阵,其大小为(Sl×Sl)×(Sh×Sh);
低分辨率检测器大小为Sl×Sl个像素,目标高分辨率投影图像大小为Sh×Sh个像素;
步骤2.3:将步骤2.2得到的第n个角度下对第m个低分辨率投影图像PLnm与步骤一直接采集到的第n个角度下低分辨率图像按照公式(13)进行比较,求得真实低分辨率投影与目前低分辨率投影的差值ΔPLnm
步骤2.4:将经步骤2.3得到的对应像素点的差值ΔPLnm写为向量形式,即ΔpLnm,利用关系矩阵W,乘以每次调节的步长作为调节权重,得到对高分辨率图像pHnm1中每个点的调节权重ΔpHnm1
ΔpHnm1=WTΔpLnm×step (14)
其中:step为每次的调节步长,调节步长在0到1之间取值;
步骤2.5:由步骤2.4得到的调节权重ΔpHnm1还原为矩阵形式ΔPHnm1,并将其坐标系绕图像中心反向旋转γ度,由公式(15):
其中(x4,y4)为ΔPHnm1中像素点的坐标,(x5,y5)为(x4,y4)反向旋转坐标系后对应像素点的位置,并根据二维线性插值,得到PHnm的调节权重ΔPHnm
其中,(x50,y50)为(x5,y5)向下取整的第一个整数值坐标,(x6,y6)为(x5,y5)插值后ΔPHnm的像素点坐标值;
步骤2.6:将调节权重ΔPHnm叠加到在检测器第m次偏心旋转位置的初始高分辨率投影数据PHnm上,得到第(m+1)次偏心旋转位置的初始高分辨率投影数据PHn(m+1),由公式(17)更新PHn(m+1)
PHn(m+1)=PHnm+ΔPHnm (17)
步骤2.7;此时得到第n个旋转角度下,第m次偏心旋转位置的高分辨率投影图像,即第m+1次偏心旋转位置的初始高分辨率投影图像;重复步骤2、1-2、6,即能够得到第n个旋转角度下,经过M次处理后的高分辨率图像PHnM
此时完成全部M次调节;
再重复步骤2、1—2、6,继续从第n个旋转角度下的第1个位置开始调节图像,直到达到要求的迭代结束条件,得到满足投影分辨率要求的高分辨率投影图像,结束迭代过程,即得到检测器绕物体旋转时单个旋转角度下的高分辨率投影图像;
所述的迭代结束条件为:计算前后两次迭代过程高分辨率投影图像的差值的二范数,并与结束迭代阈值条件进行比较,从而确定是否结束迭代;所述的迭代阈值条件根据需要获得的高分辨率投影图像精度而定;
步骤三、重复步骤二操作,得到N个角度下的高分辨率投影图像,最终得到真实高分辨率图像。
2.如权利要求1所述的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,其特征在于:检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的M个交点,与基准平面中心点的距离为高分辨率像素大小的任意倍数。
3.如权利要求1所述的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,其特征在于:检测器焦点与偏心旋转后的检测器中心连线的延长线与基准平面的M个交点,可以均匀分布,也可以不均匀分布,但要记录每个交点的位置坐标。
4.如权利要求1所述的一种利用偏心旋转提高锥形束SPECT成像分辨率的方法,其特征在于:单个旋转角度下,单次偏心旋转得到的低分辨率投影图像与高分辨率投影图像之间存在几何关系:
设检测器焦点与第m次偏心旋转后的检测器中心连线的延长线与基准平面的交点和基准平面中心的距离为d(m),检测器焦点到检测器中心的距离(即焦距)为sid,高分辨率投影图像在基准平面内,以高分辨率投影图像中心为高分辨率投影图像坐标原点,以偏心旋转后低分辨率检测器中心为低分辨率投影图像坐标原点,对于低分辨率检测器中的任意像素点,取像素边界的四个点坐标为(i,j),(i,j+1),(i+1,j),(i+1,j+1),得到其在高分辨率投影图像中对应区域块S*的边界点坐标,以坐标(i,j)为例;
φ(j)=θ+δ(j) (6)
g(j)=sid×tan(φ(j)) (7)
其中,θ为基准平面焦距与偏心旋转后检测器焦距之间的夹角,δ(j)为低分辨率投影图像中坐标j与偏心旋转后检测器焦距之间的夹角,φ(j)为低分辨率投影图像中坐标j与基准平面焦距之间的夹角,g(j)为低分辨率投影图像中坐标j映射到高分辨率投影图像中的坐标,g(i)为低分辨率投影图像中坐标i映射到高分辨率投影图像中的坐标;
对于高分辨率投影图像中的任意像素点(a,b),判断其是否属于(g(i),g(j)),(g(i),g(j+1)),(g(i+1),g(j)),(g(i+1),g(j+1))组成的区域块S*内,若满足则对应w(a,b)=1,这里w(a,b)表示S*内各高分辨率像素点对低分辨率像素点(a,b)的系数,由此得到低分辨率投影图像中任一像素点与高分辨率投影图像各像素点之间的映射关系。
CN201910588491.9A 2019-07-01 2019-07-01 一种利用偏心旋转提高锥形束spect成像分辨率的方法 Active CN110428370B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910588491.9A CN110428370B (zh) 2019-07-01 2019-07-01 一种利用偏心旋转提高锥形束spect成像分辨率的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910588491.9A CN110428370B (zh) 2019-07-01 2019-07-01 一种利用偏心旋转提高锥形束spect成像分辨率的方法

Publications (2)

Publication Number Publication Date
CN110428370A true CN110428370A (zh) 2019-11-08
CN110428370B CN110428370B (zh) 2021-11-23

Family

ID=68409988

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910588491.9A Active CN110428370B (zh) 2019-07-01 2019-07-01 一种利用偏心旋转提高锥形束spect成像分辨率的方法

Country Status (1)

Country Link
CN (1) CN110428370B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113100803A (zh) * 2021-04-20 2021-07-13 西门子数字医疗科技(上海)有限公司 用于显示静脉血栓的方法、装置、计算机设备和介质
CN116485789A (zh) * 2023-06-16 2023-07-25 新创碳谷集团有限公司 一种碳纤维劈丝缺陷检测方法、设备及存储介质
CN116483025A (zh) * 2023-04-23 2023-07-25 赛诺威盛科技(北京)股份有限公司 飞焦点模式下的数据采集***、方法、电子设备及介质

Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101185577A (zh) * 2006-11-24 2008-05-28 通用电气公司 利用多点发射源进行ct成像的方法和***
WO2010133983A1 (en) * 2009-05-18 2010-11-25 Koninklijke Philips Electronics, N.V. Interpolation free fan-to-parallel beam re-binning
CN102743182A (zh) * 2012-01-12 2012-10-24 北京理工大学 一种提高扇形束spect成像分辨率的方法
CN102842141A (zh) * 2012-07-03 2012-12-26 东南大学 一种旋转x射线造影图像迭代重建方法
CN102918565A (zh) * 2010-05-27 2013-02-06 皇家飞利浦电子股份有限公司 用于采用偏心平板检测器的锥形射束计算机断层摄影成像的改进重建
US20150043795A1 (en) * 2013-08-07 2015-02-12 Toshiba Medical Systems Corporation Image domain pansharpening method and system for spectral ct with large pixel energy discriminating detectors
CN106373087A (zh) * 2016-08-23 2017-02-01 大连理工大学 一种改进初始估计的图像超分辨率重建方法
CN107157505A (zh) * 2017-06-09 2017-09-15 北京理工大学 一种提高锥形束spect成像分辨率的方法
CN108283503A (zh) * 2018-02-12 2018-07-17 沈阳晟诺科技有限公司 一种ct机、扫描方法及图像重建方法
CN109187591A (zh) * 2018-06-04 2019-01-11 东南大学 一种x射线超分辨成像方法及其应用
CN109171792A (zh) * 2018-09-29 2019-01-11 江苏影医疗设备有限公司 成像方法及使用该成像方法的ct成像***
CN109325931A (zh) * 2018-08-22 2019-02-12 中北大学 基于生成对抗网络和超分辨率网络的多模态图像融合方法

Patent Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101185577A (zh) * 2006-11-24 2008-05-28 通用电气公司 利用多点发射源进行ct成像的方法和***
WO2010133983A1 (en) * 2009-05-18 2010-11-25 Koninklijke Philips Electronics, N.V. Interpolation free fan-to-parallel beam re-binning
CN102918565A (zh) * 2010-05-27 2013-02-06 皇家飞利浦电子股份有限公司 用于采用偏心平板检测器的锥形射束计算机断层摄影成像的改进重建
CN102743182A (zh) * 2012-01-12 2012-10-24 北京理工大学 一种提高扇形束spect成像分辨率的方法
CN102842141A (zh) * 2012-07-03 2012-12-26 东南大学 一种旋转x射线造影图像迭代重建方法
US20150043795A1 (en) * 2013-08-07 2015-02-12 Toshiba Medical Systems Corporation Image domain pansharpening method and system for spectral ct with large pixel energy discriminating detectors
CN106373087A (zh) * 2016-08-23 2017-02-01 大连理工大学 一种改进初始估计的图像超分辨率重建方法
CN107157505A (zh) * 2017-06-09 2017-09-15 北京理工大学 一种提高锥形束spect成像分辨率的方法
CN108283503A (zh) * 2018-02-12 2018-07-17 沈阳晟诺科技有限公司 一种ct机、扫描方法及图像重建方法
CN109187591A (zh) * 2018-06-04 2019-01-11 东南大学 一种x射线超分辨成像方法及其应用
CN109325931A (zh) * 2018-08-22 2019-02-12 中北大学 基于生成对抗网络和超分辨率网络的多模态图像融合方法
CN109171792A (zh) * 2018-09-29 2019-01-11 江苏影医疗设备有限公司 成像方法及使用该成像方法的ct成像***

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
ZIYE YAN ET AL.: "Super resolution SPECT reconstruction with non-uniform attenuation", 《COMPUTERS IN BIOLOGY AND MEDICINE》 *
曾凯等: "图像超分辨率重建的研究进展", 《计算机工程与应用》 *
陈云斌等: "非均匀衰减锥形投影SPECT局部重建", 《仪器仪表学报》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113100803A (zh) * 2021-04-20 2021-07-13 西门子数字医疗科技(上海)有限公司 用于显示静脉血栓的方法、装置、计算机设备和介质
CN116483025A (zh) * 2023-04-23 2023-07-25 赛诺威盛科技(北京)股份有限公司 飞焦点模式下的数据采集***、方法、电子设备及介质
CN116483025B (zh) * 2023-04-23 2024-03-22 赛诺威盛科技(北京)股份有限公司 飞焦点模式下的数据采集***、方法、电子设备及介质
CN116485789A (zh) * 2023-06-16 2023-07-25 新创碳谷集团有限公司 一种碳纤维劈丝缺陷检测方法、设备及存储介质
CN116485789B (zh) * 2023-06-16 2023-08-25 新创碳谷集团有限公司 一种碳纤维劈丝缺陷检测方法、设备及存储介质

Also Published As

Publication number Publication date
CN110428370B (zh) 2021-11-23

Similar Documents

Publication Publication Date Title
CN104107065B (zh) 3d图像集在不同空间之间的最佳变换
Kennedy et al. Super-resolution in PET imaging
JP5734664B2 (ja) 希薄化制約補正を用いた画像復元法
US7711087B2 (en) Patient setup using tomosynthesis techniques
JP5123191B2 (ja) 高度に限定された画像再構成法を使用する拡散テンソル・イメージング
CN107111867B (zh) 多模态成像***及方法
CN106725565B (zh) 一种稀疏投影下的锥束xct成像质量评估方法
Li et al. Radiation dose reduction in four‐dimensional computed tomography
EP1800264B1 (en) Image reconstruction with voxel dependent interpolation
CN102525527B (zh) Ct成像的投影数据加权方法
CN109187591B (zh) 一种x射线超分辨成像方法及其应用
CN105190696A (zh) 用于同时的图像伪像减少和断层扫描重建的***和方法
CN107157505B (zh) 一种提高锥形束spect成像分辨率的方法
CN100581471C (zh) 用于检查周期性运动的对象的ct方法
CN103608839A (zh) 对比度相关分辨率图像
US20170215818A1 (en) High-resolution computed tomography or c-arm imaging
CN107249465B (zh) 用于稀疏角度采样的断层摄影成像装置和方法
CN109146987B (zh) 一种基于gpu的快速锥束计算机断层成像重建方法
CN110428370A (zh) 一种利用偏心旋转提高锥形束spect成像分辨率的方法
Park et al. A fully GPU-based ray-driven backprojector via a ray-culling scheme with voxel-level parallelization for cone-beam CT reconstruction
US7058156B2 (en) Imaging method for a multi-slice spiral CT scan with 3D reconstruction, and a computed tomography unit for carrying out this method
Zhong et al. A dual‐view digital tomosynthesis imaging technique for improved chest imaging
CN110264536A (zh) 一种在平行束超分重建中计算高低分辨率投影关系的方法
Cheng et al. Super-resolution acquisition and reconstruction for cone-beam SPECT with low-resolution detector
Hoppe Accurate cone-beam image reconstruction in C-Arm computed tomography

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