CN117078791B - 一种ct环形伪影校正方法、装置、电子设备及存储介质 - Google Patents
一种ct环形伪影校正方法、装置、电子设备及存储介质 Download PDFInfo
- Publication number
- CN117078791B CN117078791B CN202311329514.7A CN202311329514A CN117078791B CN 117078791 B CN117078791 B CN 117078791B CN 202311329514 A CN202311329514 A CN 202311329514A CN 117078791 B CN117078791 B CN 117078791B
- Authority
- CN
- China
- Prior art keywords
- original
- gradient
- image
- polar coordinate
- stripe
- 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
- 238000000034 method Methods 0.000 title claims abstract description 77
- 238000012937 correction Methods 0.000 title claims abstract description 41
- 230000009466 transformation Effects 0.000 claims abstract description 57
- 238000009499 grossing Methods 0.000 claims abstract description 48
- 238000012545 processing Methods 0.000 claims description 26
- 238000001914 filtration Methods 0.000 claims description 16
- 230000003044 adaptive effect Effects 0.000 claims description 11
- 238000000354 decomposition reaction Methods 0.000 claims description 8
- 230000006870 function Effects 0.000 claims description 4
- 230000001131 transforming effect Effects 0.000 claims description 3
- 230000000694 effects Effects 0.000 abstract description 19
- 239000002699 waste material Substances 0.000 abstract description 9
- 230000000007 visual effect Effects 0.000 abstract description 6
- 230000008569 process Effects 0.000 description 22
- 238000007781 pre-processing Methods 0.000 description 9
- 238000011156 evaluation Methods 0.000 description 7
- 238000004364 calculation method Methods 0.000 description 5
- 238000006243 chemical reaction Methods 0.000 description 5
- 238000012805 post-processing Methods 0.000 description 5
- 238000004422 calculation algorithm Methods 0.000 description 3
- 238000001514 detection method Methods 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 238000002591 computed tomography Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000013139 quantization Methods 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 230000000717 retained 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/008—Specific post-processing after tomographic reconstruction, e.g. voxelisation, metal artifact correction
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
- Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
- Y02P90/30—Computing systems specially adapted for manufacturing
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本申请提供了一种CT环形伪影校正方法、装置、电子设备及存储介质,对采用笛卡尔坐标系的原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像;确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到平滑垂直条带梯度,其中,平滑垂直条带梯度包括的伪影信息小于原始垂直条带梯度包括的伪影信息;根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到校正极坐标图像;对校正极坐标图像进行笛卡尔坐标变换,得到采用笛卡尔坐标系的校正图像。通过上述方法,去除效果较为直观,并且无需对投影数据进行额外的存储和处理,避免了资源浪费。
Description
技术领域
本申请涉及图像处理领域,具体涉及一种CT环形伪影校正方法、装置、电子设备及存储介质。
背景技术
由于电子计算机断层扫描(Computed Tomography,CT)系统在实际应用中,很难保证每个探测通道都保持一样的精度,并且CT系统在使用后探测单元都会受到不同程度的损耗,故而CT系统在实际应用中很难实现各探测器精度的一致,即各探测器会存在不同程度的误差,这些不同程度的误差会使得后续重建的CT图像中带有环形伪影。
在工业领域,CT图像中的环形伪影一方面会使得CT图像的原始信息被覆盖,无法正确表达CT图像所对应的真实信息,另一方面会对图像后续处理造成干扰,例如会对图像分割、自动识别和量化分析等图像后续处理造成干扰。
相关技术中主要通过投影的正弦图的前处理方式对环形伪影进行校正,具体的,投影的正弦图的前处理方式是在收集到CT系统的投影数据后,在图像重建之前对投影正弦图进行处理,从而避免在图像重建后出现环形伪影。
虽然通过投影的正弦图的前处理方式能够对环形伪影进行一定程度的校正,但是由于投影的正弦图的前处理方式是在图像重建之前实施的,因此全部的投影数据都需要以正弦图的形式进行额外的存储和处理,从而导致消耗更多内存,需要通过图像重建才能验证去除效果,即去除效果不够直观。
发明内容
本申请实施例提供了一种CT环形伪影校正方法、装置、电子设备及存储介质,针对包括环形伪影的原始图像,能够在进行极坐标变换将环形伪影变换为条带伪影之后,通过小波傅里叶方式对与条带伪影垂直的方向上的梯度进行快速平滑来对环形伪影进行校正,去除效果较为直观,并且在此过程中,无需对投影数据进行额外的存储和处理,避免了环形伪影的校正过程中的资源浪费。
有鉴于此,本申请实施例第一方面提供一种CT环形伪影校正方法,所述方法包括:
对原始图像进行极坐标变换,得到所述原始图像对应的原始极坐标图像;其中所述原始图像采用笛卡尔坐标系且包括环形伪影,所述原始极坐标图像采用极坐标系且包括与所述环形伪影对应的条带伪影;
确定所述原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;其中所述原始垂直条带梯度是指所述原始极坐标图像在与所述条带伪影垂直的方向上的梯度,所述原始平行条带梯度是指所述原始极坐标图像在与所述条带伪影平行的方向上的梯度;
通过小波傅里叶方式对所述原始垂直条带梯度进行平滑处理,得到所述原始垂直条带梯度对应的平滑垂直条带梯度;其中所述平滑垂直条带梯度包括的伪影信息小于所述原始垂直条带梯度包括的伪影信息;
根据所述平滑垂直条带梯度和所述原始平行条带梯度,去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像;
对所述校正极坐标图像进行笛卡尔坐标变换,得到所述原始图像对应的校正图像;所述校正图像采用笛卡尔坐标系。
本申请实施例第二方面提供一种CT环形伪影校正装置,所述装置包括:
第一变换单元,用于对原始图像进行极坐标变换,得到所述原始图像对应的原始极坐标图像;其中所述原始图像采用笛卡尔坐标系且包括环形伪影,所述原始极坐标图像采用极坐标系且包括与所述环形伪影对应的条带伪影;
确定单元,用于确定所述原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;其中所述原始垂直条带梯度是指所述原始极坐标图像在与所述条带伪影垂直的方向上的梯度,所述原始平行条带梯度是指所述原始极坐标图像在与所述条带伪影平行的方向上的梯度;
平滑单元,用于通过小波傅里叶方式对所述原始垂直条带梯度进行平滑处理,得到所述原始垂直条带梯度对应的平滑垂直条带梯度;其中所述平滑垂直条带梯度包括的伪影信息小于所述原始垂直条带梯度包括的伪影信息;
去除单元,用于根据所述平滑垂直条带梯度和所述原始平行条带梯度,去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像;
第二变换单元,用于对所述校正极坐标图像进行笛卡尔坐标变换,得到所述原始图像对应的校正图像;所述校正图像采用笛卡尔坐标系。
本申请实施例第三方面提供一种电子设备,包括:
存储器,用于存储可执行指令;
处理器,用于执行所述存储器中存储的可执行指令时,实现本申请实施例提供的CT环形伪影校正方法。
本申请实施例第四方面提供一种计算机可读介质,存储有可执行指令,用于被处理器执行时,实现本申请实施例提供的CT环形伪影校正方法。
综上所述,本申请实施例提供了一种CT环形伪影校正方法、装置、电子设备及存储介质,针对包括环形伪影的原始图像,对采用笛卡尔坐标系的原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像,原始极坐标图像采用极坐标系,通过极坐标变换,能够将原始图像中的环形伪影转换为原始极坐标图像中的条带伪影;确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度,其中,原始垂直条带梯度是指原始极坐标图像在与条带伪影垂直的方向上的梯度,原始平行条带梯度是指原始极坐标图像在与条带伪影平行的方向上的梯度;由于条带伪影对与条带伪影垂直的方向上的梯度有明显影响,而对与条带伪影平行的方向上的梯度基本没有影响,故可以通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度,其中,平滑垂直条带梯度包括的伪影信息小于原始垂直条带梯度包括的伪影信息,通过小波傅里叶方式能够快速对原始垂直条带梯度进行平滑,避免平滑时多次迭代造成的资源浪费,并为了减少计算量,对原始平行条带梯度不进行平滑处理;由于平滑垂直条带梯度中的伪影信息得到了有效去除,故可以根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;对校正极坐标图像进行笛卡尔坐标变换,得到原始图像对应的采用笛卡尔坐标系的校正图像。通过上述方法,针对包括环形伪影的原始图像,能够在进行极坐标变换将环形伪影变换为条带伪影之后,通过小波傅里叶方式对与条带伪影垂直的方向上的梯度进行快速平滑来对环形伪影进行校正,去除效果较为直观,并且在此过程中,无需对投影数据进行额外的存储和处理,避免了环形伪影的校正过程中的资源浪费。
附图说明
为了更清楚地说明本公开实施例的技术方案,下面将对实施例中所需要使用的附图作简单地介绍。此处的附图并入说明书中并构成本说明书的一部分,这些附图示出了符合本公开的实施例,并与说明书一起用于说明本公开的技术方案。应当理解,附图仅示出了本公开的某些实施例,因此不应被看作是对保护范围的限制,对于本领域普通技术人员而言,在不付出创造性劳动的前提下,还可以根据这些附图获得其他相关的附图。而且在整个附图中,用相同的标号表示相同的部件。在附图中:
图1为本申请实施例提供到的一种CT环形伪影校正方法的方法流程图;
图2为本申请实施例提供的一种对水平梯度进行平滑处理的流程图;
图3为本申请实施例提供的一种仿真数据的环形伪影校正结果对比图;
图4为本申请实施例提供的一种实测数据的环形伪影校正结果对比图;
图5为本申请实施例提供的一种CT环形伪影校正方法的具体流程图;
图6为本申请实施例提供的一种CT环形伪影校正装置的装置示意图。
具体实施方式
为了使本技术领域的人员更好地理解本申请方案,下面将结合本申请实施例中的附图,对本申请实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅是本申请一部分实施例,而不是全部的实施例。基于本申请中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本申请保护的范围。
本申请的说明书和权利要求书及上述附图中的术语“第一”、“第二”、“第三”、“第四”等(如果存在)是用于区别类似的对象,而不必用于描述特定的顺序或先后次序。应该理解这样使用的数据在适当情况下可以互换,以便这里描述的本申请的实施例能够以除了在这里图示或描述的那些以外的顺序实施。此外,术语“包括”和“具有”以及他们的任何变形,意图在于覆盖不排他的包含,例如,包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元,而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。
由于在工业CT系统中,对目标对象进行测量和识别都依赖于CT图像的质量,只有正确和清晰的CT图像才能让检测者做出正确的测量与识别。但是多种的伪影影响了图像质量,给图像后续处理带来了干扰。
环形伪影作为常见伪影的一种,已成为影响CT图像质量的一个重要因素,环形伪影的产生原因主要是由于CT系统中各探测器存在的不同程度的误差。
对此,CT图像的环形伪影校正可以提高对检测结果的判断准确率,即通过软件校正方法可以弥补硬件校正后的不足,更好的提升CT图像的质量。
相关技术中主要通过投影的正弦图的前处理方式对环形伪影进行校正,具体的,投影的正弦图的前处理方式是在收集到CT系统的投影数据后,在图像重建之前对投影正弦图进行处理,从而避免在图像重建后出现环形伪影。
此外,还可以通过重建图像的后处理方式对环形伪影进行校正,重建图像的后处理方式是指直接针对重建的CT图像进行处理,具体的,通过极坐标变换将环形伪影转换为直线特征再校正,环形伪影在投影正弦图或极坐标图上表现为多条平行直线或多条垂直平行直线,且偏向于高频信号,故使用简单的低通滤波器过滤掉图像的高频部分可以实现环形伪影去除。重建图像的后处理方式还包括小波傅里叶方式(Wavelet Fourier,WF)、中值均值结合滤波方式(Mean and Median Filter,MMF)和单向变分模型(UnidirectionalTotal Variation-Stokes,UTV)。小波傅里叶方式是指小波分解和傅里叶变换相结合的方式,小波傅里叶方式先是将图像小波分解为低频分量和高频分量,然后再对垂直高频带信息进行傅里叶变换和高斯低通滤波去除伪影分量,最后通过小波逆变换得到校正后的图像。中值均值结合滤波方式可以对伪影进行校正,对伪影垂直方向进行中值滤波,然后对滤波前后差分图像进行均值滤波获得伪影图像,最后减去该伪影图像得到校正后图像。单向变分模型可以对伪影进行校正,在极坐标图像下通过单向变分模型只对与环形伪影垂直的方向的梯度进行平滑,从而得到光滑的法向量,再用图像恢复模型得到去噪后的校正图像。
虽然通过投影的正弦图的前处理方式和重建图像的后处理方式能够对环形伪影进行一定程度的校正,但是对于投影的正弦图的前处理方式,由于投影的正弦图的前处理方式是在图像重建之前实施的,因此全部的投影数据都需要以正弦图的形式进行额外的存储和处理,从而导致消耗更多内存,需要通过图像重建才能验证去除效果,即去除效果不够直观。而对于重建图像的后处理方式,低通滤波器虽然可以消除环形伪影,但是也会在一定程度上影响图像中的其他高频信号,即通过简单的低通滤波器得到的CT图像质量不高;小波傅里叶方式虽然拥有较快的环形伪影校正效率,但是通过小波傅里叶对图像进行校正后得到的图像中会生成一些新的伪影,例如在灰度变化较为明显的区域附近会很容易生成新的伪影;中值均值结合滤波方式则涉及到滤波器窗口调节,滤波器窗口大小需要从去除伪影和保留信息两个方面进行考虑,如果图像中的伪影强度存在一定的差异,则强度更大的伪影将会被部分保留,此时,若要尽可能去除伪影则会导致图像质量明显下降;单向变分模型可以很好的去除环形伪影,其过程涉及图像切向量平滑和图像恢复,由于图像切向量平滑和图像恢复都使用了较多的迭代计算,故单向变分模型在处理高分辨率图像时速度很慢,效率较低。
鉴于此,本申请实施例提供了一种CT环形伪影校正方法、装置、电子设备及存储介质,针对包括环形伪影的原始图像,能够在进行极坐标变换将环形伪影变换为条带伪影之后,通过小波傅里叶方式对与条带伪影垂直的方向上的梯度进行快速平滑来对环形伪影进行校正,去除效果较为直观,并且在此过程中,无需对投影数据进行额外的存储和处理,避免了环形伪影的校正过程中的资源浪费。
本申请实施例所提供的CT环形伪影校正方法可以通过计算机设备实施,该计算机设备可以是终端设备或服务器,其中,服务器可以是独立的物理服务器,也可以是多个物理服务器构成的服务器集群或者分布式系统,还可以是提供云计算服务的云服务器。终端设备包括但不限于手机、电脑、智能语音交互设备、智能家电、车载终端、飞行器等。终端设备以及服务器可以通过有线或无线通信方式进行直接或间接地连接,本申请在此不做限制。
下面通过方法实施例来对本申请提供的一种CT环形伪影校正方法进行说明,如图1所示,图1为本申请实施例提供到的一种CT环形伪影校正方法的方法流程图,前述的计算机设备可以为服务器,该方法包括:
S101、对原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像。
原始图像是指包括环形伪影的图像,例如,若通过CT系统获取的CT图像中包括环形伪影,该CT图像可以作为原始图像。
原始图像采用笛卡尔坐标系,笛卡尔坐标系是一种正交坐标系,也称为直角坐标系,原始图像中的环形伪影表现为一系列圆心与原始图像的图像中心重合的同心圆环,即环形伪影会对图像的各个方向都造成影响,难以被检测和提取。
对此,服务器可以对原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像,原始极坐标图像采用极坐标系,原始图像中的环形伪影在原始极坐标图像中会变换为条带伪影,极坐标图像中的条带伪影表现为相互平行的条纹,这种相互平行的条带伪影比环形伪影更易检测和提取。
在实际应用中,在得到包括环形伪影的CT图像时,可以将该CT图像作为原始图像,为了便于后续的处理,可以将原始图像归一化至0-1,针对笛卡尔坐标下的原始图像F(x,y),可以通过下述公式进行极坐标变换,得到对应的原始极坐标图像:
其中,x为笛卡尔坐标系下的横坐标,y为笛卡尔坐标系下的纵坐标,
ρ为极坐标系下的半径坐标,表示与极点的距离,θ为极坐标系下的角坐标,表示按逆时针方向距离0°射线(有时也称作极轴)的角度。
需要说明的是,由于环形伪影表现为一系列圆心与原始图像的图像中心重合的同心圆环,故可以以原始图像的坐标(x,y)为中心,进行极坐标变换,变换的最大半径距离r为原始图像的坐标(x,y)到原始图像四角的距离最大值,相应的公式如下:
其中,r表示最大半径距离,(x,y)表示原始图像的坐标,M×N表示原始图像的尺寸。
在得到变换的最大半径距离之后,原始极坐标图像对应的半径坐标的取值需要不大于最大半径距离,即半径坐标的取值范围为。
在一种可能的实现方式中,S101对原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像,包括:
对原始图像进行极坐标变换,得到初始极坐标图像;
通过双线性差值方式对初始极坐标图像进行像素补偿,得到原始极坐标图像。
由于在对原始图像进行极坐标变换的实际过程中,通过极坐标变换得到的初始极坐标图像将损失较多的分辨率,故在本实施例中,服务器可以采用合适的插值算法,由目标像素点周围已知的像素点来计算目标像素点所需的像素值,从而减少像素损失,即服务器可以采用双线性插值方式对初始极坐标图像进行像素补偿,得到所需的原始极坐标图像。
具体的,假设已知像素点(X1,Y1)、(X1,Y2)、(X2,Y1)、(X2,Y2)的像素值分别为A11、A12、A21、A22,需要(X,Y)计算处的像素值P,服务器可以先在x方向上进行线性插值计算出像素点(X,Y1)和(X,Y2)处的像素值Q1、Q2,对应的公式如下所示:
其中,Q1、Q2分别表示像素点(X,Y1)和(X,Y2)处的像素值,A11、A12、A21、A22分别表示像素点(X1,Y1)、(X1,Y2)、(X2,Y1)、(X2,Y2)的像素值。
然后服务器可以在y方向上进行线性插值,得到像素点(X,Y)处的像素值P,对应的公式如下所示:
其中,P表示像素点(X,Y)处的像素值,Q1、Q2分别表示像素点(X,Y1)和(X,Y2)处的像素值。
S102、确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度。
原始垂直条带梯度是指原始极坐标图像在与条带伪影垂直的方向上的梯度,例如,当条带伪影为相互平行的垂直条纹时,原始垂直条带梯度是指原始极坐标图像在水平方向上的梯度。
原始平行条带梯度是指原始极坐标图像在与条带伪影平行的方向上的梯度,例如,当条带伪影为相互平行的垂直条纹时,原始平行条带梯度是指原始极坐标图像在垂直方向上的梯度。
服务器可以确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度,在本实施例的实际应用中,服务器可以通过图像有限差方式来确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度,例如当条带伪影为相互平行的竖直条纹时,服务器可以通过图像有限差方式来确定原始极坐标图像在水平方向上的水平梯度和在垂直方向上的垂直梯度,假定原始极坐标图像u的大小是M×N,可以通过下述公式确定原始极坐标图像中某个像素u(i, j)的水平前向差分、水平后向差分、垂直前向差分、垂直后向差分,水平前向差分和水平后向差分即是水平梯度,垂直前向差分和垂直后向差分即是垂直梯度:
其中,表示水平前向差分,/>表示水平后向差分,/>表示垂直前向差分,/>表示垂直后向差分,i表示M行中的第i行,j表示N列中的第j列,/>表示图像中第α行和第β列对应的像素点的像素值。
S103、通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度。
由于条带伪影在与条带伪影平行的方向上的梯度基本为零,而条带伪影在与条带伪影垂直的方向上的梯度相比于与条带伪影平行的方向来说较大,即条带伪影对与条带伪影垂直的方向上的梯度有明显影响,而对与条带伪影平行的方向上的梯度基本没有影响。
对此,可以只对与条带伪影垂直的方向上的梯度进行平滑来去除条带伪影,即服务器可以通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度,由于对原始垂直条带梯度进行了平滑处理,故得到的平滑垂直条带梯度包括的伪影信息会小于原始垂直条带梯度包括的伪影信息。
需要说明的是,在本实施例中,为了能够提高校正效率,服务器采用了小波傅里叶方式来对垂直条带梯度进行快速平滑处理,与单向变分模型需要通过多次迭代计算才能对切向量进行平滑相比,小波傅里叶方式只需要通过一次处理即可对垂直条带梯度进行平滑,即通过加快垂直条带梯度的平滑速度提高了算法的校正效率。
同时,在本实施例中,并未采用小波傅里叶方式直接对原始极坐标图像进行校正处理,而是采用小波傅里叶方式对原始极坐标图像对应的原始垂直条带梯度进行处理,从而避免小波傅里叶方式直接对图像进行处理会在图像中生成新的伪影的问题,进而保证后续步骤中得到的校正图像的图像质量。
在一种可能的实现方式中,S103中通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度,包括:
S11、对原始垂直条带梯度进行小波分解,得到原始垂直条带梯度对应的原始低频分量、原始对角高频分量、原始垂直高频分量和原始平行高频分量;
S12、对原始平行高频分量进行傅里叶变换并中心化,得到原始平行高频分量对应的原始频域信息;
S13、对原始频域信息进行低通滤波,得到原始频域信息对应的平滑频域信息;平滑频域信息包括的高频信息小于原始频域信息包括的高频信息;
S14、对平滑频域信息进行傅里叶逆变换,得到原始平行高频分量对应的平滑平行高频分量;
S15、基于平滑平行高频分量、原始低频分量、原始对角高频分量和原始垂直高频分量,通过小波重构得到平滑垂直条带梯度。
具体的,服务器可以对原始垂直条带梯度进行小波分解,得到原始垂直条带梯度对应的原始低频分量、原始对角高频分量、原始垂直高频分量和原始平行高频分量,其中,原始低频分量是指原始垂直条带梯度对应的频率较低的分量,原始对角高频分量是指原始垂直条带梯度对应的在对角线方向的频率较高的分量,原始垂直高频分量是指原始垂直条带梯度对应的在与条带伪影垂直的方向上的分量,原始平行高频分量是指原始垂直条带梯度对应的与条带伪影平行的方向上的分量,例如,如图2所示,当条带伪影为相互平行的垂直条纹、原始垂直条带梯度为水平梯度时,对水平梯度进行小波分解,可以得到低频带对应的原始低频分量、对角高频带对应的原始对角高频分量、水平高频带对应的原始垂直高频分量和垂直高频带对应的原始平行高频分量。
由于条带伪影表现为相互平行的条纹,故对原始垂直条带梯度进行小波分解之后,伪影信息主要在原始平行高频分量之中,例如,如图2所示,当条带伪影为相互平行的垂直条纹、原始垂直条带梯度为水平梯度时,伪影信息主要在垂直高频带对应的原始平行高频分量之中,故服务器可以对原始平行高频分量进行傅里叶变换并中心化,将空间域中的原始平行高频分量转换到频率域中,得到原始平行高频分量对应的原始频域信息,原始频域信息可以是频谱图的形式。
为了去除原始平行高频分量所携带的伪影信息,服务器可以对原始频域信息进行低通滤波,得到原始频域信息对应的平滑频域信息,低通滤波后的平滑频域信息包括的高频信息小于原始频域信息包括的高频信息。
在得到平滑频域信息之后,服务器可以对平滑频域信息进行傅里叶逆变换,将频率域中的平滑频域信息转换到空间域中,得到原始平行高频分量对应的平滑平行高频分量。
在得到平滑平行高频分量之后,服务器可以基于平滑平行高频分量、原始低频分量、原始对角高频分量和原始垂直高频分量,通过小波重构得到平滑垂直条带梯度,如图2所示,当条带伪影为相互平行的垂直条纹、原始垂直条带梯度为水平梯度时,可以通过小波重构得到平滑的水平梯度图。由于伪影信息主要在原始平行高频分量之中,通过对原始平行高频分量进行低通滤波能够去除相应的伪影信息,即平滑平行高频分量中的伪影信息已被有效去除,故平滑垂直条带梯度中的伪影信息也被有效去除。
在一种可能的实现方式中,当原始平行高频分量对应于竖直方向时,S12中对原始平行高频分量进行傅里叶变换并中心化,得到原始平行高频分量对应的原始频域信息,包括:
对原始平行高频分量按列进行傅里叶变换并中心化,得到频域信息。
当原始平行高频分量对应于竖直方向时,为了能够对原始平行高频分量进行准确处理,避免产生新的伪影信息,如图2所示,服务器可以对原始平行高频分量按列进行傅里叶变换并中心化,得到对应的频域信息。
同理,当原始平行高频分量对应于水平方向时,服务器可以对原始平行高频分量按行进行傅里叶变换并中心化,得到对应的频域信息。
在一种可能的实现方式中,S13中对原始频域信息进行低通滤波,得到原始频域信息对应的平滑频域信息,包括:
基于高斯函数,得到原始频域信息对应的衰减系数;
根据衰减系数,对原始频域信息进行系数衰减处理,得到平滑频域信息。
具体的,服务器可以基于高斯函数,得到原始频域信息对应的衰减系数,相应的公式如下所示:
其中,coef表示衰减系数,σ表示常量系数,h表示原始平行高频分量的行数,i表示原始平行高频分量的第i行,j表示原始平行高频分量的第j列。
在得到衰减系数之后,如图2所示,服务器可以根据衰减系数,对原始频域信息进行系数衰减处理,得到平滑频域信息,相应的公式如下所示:
其中,FCprocessed表示平滑频域信息,FCv表示原始频域信息,coef表示衰减系数。
S104、根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像。
由于条带伪影对与条带伪影平行的方向上的梯度基本没有影响,故在S103得到平滑垂直条带梯度之后,服务器可以根据平滑垂直条带梯度和原始平行条带梯度,得到与平滑垂直条带梯度和原始平行条带梯度相匹配的校正极坐标图像,也就是去除了条带伪影的校正极坐标图像。
需要说明的是,由于平滑垂直条带梯度保留了原始极坐标图像对应的图像信息,故根据平滑垂直条带梯度和原始平行条带梯度得到的校正极坐标图像能够在去除条带伪影的同时保留原始极坐标图像对应的图像信息。
在一种可能的实现方式中,S104中根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像,包括:
根据平滑垂直条带梯度和原始平行条带梯度,得到平滑法向量;
以原始极坐标图像作为迭代约束,将平滑法向量作为迭代目标,通过K次迭代去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;针对K次迭代中的第k+1次迭代,第k+1次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异小于第k次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异,其中,K为整数,K>1,k为整数,K>k≥1。
对于任意一幅灰度图f,可以通过下述公式表示切向量和法向量:
其中,t表示切向量,n表示法向量,为水平方向上的梯度算子,/>为垂直方向上的梯度算子,/>表示图像在水平方向上的梯度,/>表示图像的在垂直方向上的梯度。
为了能够得到与平滑垂直条带梯度和原始平行条带梯度相匹配的校正极坐标图像,服务器可以先根据平滑垂直条带梯度和原始平行条带梯度,得到对应的平滑法向量,例如,当条带伪影为相互平行的垂直条纹时,原始平行条带梯度是垂直梯度,原始垂直条带梯度是水平梯度,平滑垂直条带梯度是平滑水平梯度,根据垂直梯度和平滑水平梯度可以得到对应的平滑法向量。
在得到平滑法向量之后,服务器可以以原始极坐标图像作为迭代约束,将平滑法向量作为迭代目标,通过K次迭代去除原始极坐标图像包括的条带伪影,其中,K为整数,K>1,得到原始极坐标图像对应的校正极坐标图像。
需要说明的是,以原始极坐标图像作为迭代约束,是为了在一次次迭代过程中尽可能的保留原始极坐标图像对应的图像信息,将平滑法向量作为迭代目标,是为了在一次次迭代过程中去除伪影信息,以K次迭代中的第k+1次迭代为例进行说明,其中,k为整数,K>k≥1,中间极坐标图像是指每一次迭代对应的极坐标图像,第k+1次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异会小于第k次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异,说明第k+1次迭代对应的中间极坐标图像与第k次迭代对应的中间极坐标图像相比,与平滑法向量匹配程度更高,从而说明可以通过一次次迭代来得到与平滑法向量越来越匹配的极坐标图像,以便通过K次迭代能够得到与平滑法向量相匹配的校正极坐标图像。
在一种可能的实现方式中,以原始极坐标图像作为迭代约束,将平滑法向量作为迭代目标,通过K次迭代去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像,包括:
针对K次迭代中的第k+1次迭代,基于第k次迭代对应的中间极坐标图像的法向量与平滑法向量的差异,确定第k次迭代对应的适应平滑法向量项;
基于第k次迭代对应的中间极坐标图像与原始极坐标图像在平行于所述条带伪影的方向上的梯度差异,确定第k次迭代对应的平行条带梯度约束项;
基于第k次迭代对应的中间极坐标图像与原始极坐标图像的图像差异,确定第k次迭代对应的图像约束项;
基于第k次迭代对应的中间极坐标图像在垂直于条带伪影的方向上的梯度,确定第k次迭代对应的垂直条带梯度惩罚项;
根据第k次迭代对应的适应平滑法向量项、第k次迭代对应的平行条带梯度约束项、第k次迭代对应的图像约束项和第k次迭代对应的垂直条带梯度惩罚项,以第k次迭代对应的中间极坐标图像为基础,得到第k+1次迭代对应的中间极坐标图像。
具体的,以K次迭代中第k次迭代为例对迭代过程进行说明,服务器可以基于第k次迭代对应的中间极坐标图像的法向量与平滑法向量的差异,确定第k次迭代对应的适应平滑法向量项,适应平滑法向量项用于通过一次次迭代来得到与平滑法向量越来越匹配的极坐标图像,即将平滑法向量作为迭代目标。
此外,在本实施例中,为了更好的去除条带伪影,在确定适应平滑法向量项的基础上,由于条带伪影对在垂直于条带伪影的方向上的梯度影响较大,故服务器可以基于第k次迭代对应的中间极坐标图像在垂直于条带伪影的方向上的梯度,确定第k次迭代对应的垂直条带梯度惩罚项。
为了能够将原始极坐标图像作为迭代约束,使得得到的校正极坐标图像尽可能的保留原始极坐标图像对应的图像信息,由于条带伪影对在平行于条带伪影的方向上的梯度基本上无影响,故服务器可以基于第k次迭代对应的中间极坐标图像与原始极坐标图像在平行于条带伪影的方向上的梯度差异,确定第k次迭代对应的平行条带梯度约束项,从而通过平行条带梯度约束项让每一次迭代对应的中间极坐标与原始极坐标图像在平行于条带伪影的方向上的梯度尽可能的保持一致。
同时,服务器可以基于第k次迭代对应的中间极坐标图像与原始极坐标图像的图像差异,确定第k次迭代对应的图像约束项,从而通过图像约束项让每一次迭代对应的中间极坐标与原始极坐标图像尽可能的保持一致。
根据第k次迭代对应的适应平滑法向量项、第k次迭代对应的平行条带梯度约束项、第k次迭代对应的图像约束项和第k次迭代对应的垂直条带梯度惩罚项,服务器可以以第k次迭代对应的中间极坐标图像为基础,得到第k+1次迭代对应的中间极坐标图像。
当条带伪影为相互平行的垂直条纹、原始平行条带梯度是垂直梯度、原始垂直条带梯度是水平梯度、平滑垂直条带梯度是平滑水平梯度时,相应的公式如下所示:
其中,S*k+1表示第k+1次迭代对应的中间极坐标图像,S*k表示第k次迭代对应的中间极坐标图像,为学习长度,λ1、λ2、λ3为正则化参数,ε为很小的正数避免除数为0情况,I*为原始极坐标图,/>表示水平前向差分、/>表示水平后向差分、/>和/>都可以表示水平梯度,/>表示垂直前向差分、/>表示垂直后向差分,/>和/>都可以表示垂直梯度,n*表示平滑法向量,/>。
在上述公式中,是对应的适应平滑法向量项,是对应的垂直条带梯度惩罚项,/>是对应的平行条带梯度约束项,/>是对应的图像约束项。
S105、对校正极坐标图像进行笛卡尔坐标变换,得到原始图像对应的校正图像。
在S104到校正极坐标图像之后,需要将在极坐标下展示的校正极坐标图像变换为在笛卡尔坐标系下显示,此时,服务器可以对校正极坐标图像进行笛卡尔坐标变换,得到在笛卡尔坐标系下的校正图像,由于校正极坐标图像中去除了条带伪影,故校正图像中也被去除了环形伪影。
需要说明的是,由于校正极坐标图像能够在去除条带伪影的同时保留原始极坐标图像对应的图像信息,故校正图像也能够在去除环形伪影的同时保留原始图像对应的图像信息。
在实际应用中,可以将S101中极坐标变换的转换中心和最大半径距离作为S105的笛卡尔坐标变换的转换中心和最大半径距离,进行笛卡尔坐标变换,相应的公式如下:
其中,x为笛卡尔坐标系下的横坐标,y为笛卡尔坐标系下的纵坐标,
ρ为极坐标系下的半径坐标,θ为极坐标系下的角坐标。
由此可见,本申请实施例提供了一种CT环形伪影校正方法,针对包括环形伪影的原始图像,对采用笛卡尔坐标系的原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像,原始极坐标图像采用极坐标系,通过极坐标变换,能够将原始图像中的环形伪影转换为原始极坐标图像中的条带伪影;确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度,其中,原始垂直条带梯度是指原始极坐标图像在与条带伪影垂直的方向上的梯度,原始平行条带梯度是指原始极坐标图像在与条带伪影平行的方向上的梯度;由于条带伪影对与条带伪影垂直的方向上的梯度有明显影响,而对与条带伪影平行的方向上的梯度基本没有影响,故可以通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度,其中,平滑垂直条带梯度包括的伪影信息小于原始垂直条带梯度包括的伪影信息,通过小波傅里叶方式能够快速对原始垂直条带梯度进行平滑,避免平滑时多次迭代造成的资源浪费,并为了减少计算量,对原始平行条带梯度不进行平滑处理;由于平滑垂直条带梯度中的伪影信息得到了有效去除,故可以根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;对校正极坐标图像进行笛卡尔坐标变换,得到原始图像对应的采用笛卡尔坐标系的校正图像。通过上述方法,针对包括环形伪影的原始图像,能够在进行极坐标变换将环形伪影变换为条带伪影之后,通过小波傅里叶方式对与条带伪影垂直的方向上的梯度进行快速平滑来对环形伪影进行校正,去除效果较为直观,并且在此过程中,无需对投影数据进行额外的存储和处理,避免了环形伪影的校正过程中的资源浪费。
此外,通过小波傅里叶方式对原始垂直条带梯度进行平滑处理可以加快原始垂直条带梯度的平滑速度,提高校正效率,同时,并未通过小波傅里叶方式对原始极坐标图像直接进行处理,而是只对原始极坐标图像对应的原始垂直条带梯度进行处理,从而能够在保证去除效果的同时避免出现新的伪影。
下面通过评价指标来对本申请提供的环形伪影校正方法的进行评价,具体的,在本申请中,评价指标包括均方误差(Mean Squared Error,MSE)、峰值信噪比(Peak signal-to-noise ratio,PSNR)、结构相似度(Structural Similarity Index,SSIM)。
(1)均方误差
均方误差是指原始图像和校正图像之间的每个像素值的差,是反映估计量与被估计量之间差异程度的一种度量,均方误差可以通过下述公式确定:
其中,MSE表示均方误差、I(x,y)表示原始图像,S(x,y)表示校正图像,M×N表示原始图像的尺寸。
(2)峰值信噪比
峰值信噪比是通过计算校正图像与原始图像的像素值误差来评估处理结果与原始图像的差距,峰值信噪比可以通过下述公式确定:
其中,PSNR表示峰值信噪比,MAXI表示原始图像的像素值范围,MSE表示均方误差。
(3)结构相似度
结构相似度可以用来度量两个给定图像之间的相似性,从图像组成的角度将结构信息定义为独立于亮度、对比度的,反映场景中物体结构的属性,并将失真建模为亮度、对比度和结构三个不同因素的组合。用均值作为亮度的估计,标准差作为对比度的估计,协方差作为结构相似程度的度量,结构相似度SSIM可以通过下述公式确定:
其中,I表示原始图像,S表示校正图像,表示原始图像的平均值,/>表示校正图像的平均值,/>表示原始图像的方差,/>表示校正图像的方差,/>表示原始图像和校正图像的协方差,c1、c2是较小的常数。
在确定了评价指标之后,服务器可以将模型仿真图像和实测工件图像分别作为评价数据,将小波傅里叶方式和单向变分模型作为对比方法来对本申请的环形伪影校正方法进行评价,得到的对比结果如表1、图3和图4所示。
由表1可以得出,本申请提供的方法与小波傅里叶方式相比,在均方误差、峰值信噪比和结构相似度这三个评价指标上均有了明显提高,说明与小波傅里叶方式相比,本申请提供的方法能够对环形伪影起到更好的去除效果;本申请提供的方法与单向变分模型相比,虽然在均方误差、峰值信噪比和结构相似度这三个评价指标上相差不大,但是在运行时间上有了明显提高,说明与单向变分模型相比,本申请能够在保证环形伪影的去除效果的同时,提高校正效率。
表1 环形伪影校正评价指标对比结果
如图3所示,图3为本申请实施例提供的一种仿真数据的环形伪影校正结果对比图,图3中的(a)表示仿真得到的无伪影原始图像,图3中的(b)表示仿真得到的有伪影原始图像,图3中的(c)表示通过小波傅里叶方式得到的校正图像,图3中的(d)表示通过单向变分模型得到的校正图像,图3中的(e)表示通过本申请的环形伪影校正方法得到的校正图像,由图3可知,通过小波傅里叶方式得到的校正图像中出现了新的伪影信息,通过单向变分模型得到的校正图像中环形伪影的去除效果较差,而通过本申请的环形伪影校正方法得到的校正图像中对环形伪影的去除效果较好而且没有出现新的伪影。
如图4所示,图4为本申请实施例提供的一种实测数据的环形伪影校正结果对比图,图4中提供了对四组实测数据进行环形伪影校正的处理结果,以第一组为例,图4中的(a1)表示的是第一组有伪影的实测原始图像,图4中的(b1)表示的是第一组通过小波傅里叶方式得到的校正图像,图4中的(c1)表示的是第一组通过单向变分模型得到的校正图像,图4中的(d1)表示的是第一组通过本申请的环形伪影校正方法得到的校正图像,由图4可知,通过小波傅里叶方式得到的校正图像对环形伪影的去除效果最差,通过本申请的环形伪影校正方法得到的校正图像中对环形伪影的去除效果最好。
为了便于理解,下面基于图5来对本申请提供的一种CT环形伪影校正方法的具体流程进行说明:
输入含有环形伪影的CT图像,即输入原始图像,原始图像是灰度图,为了便于后续计算,将原始图像归一化至0-1。
由环形伪影中心将原始图像变换成极坐标图,得到原始极坐标图像。
确定原始极坐标图像的水平梯度和垂直梯度,在图5对应的原始极坐标图像中条带伪影为相互平行的垂直条纹,垂直梯度是原始平行条带梯度,水平梯度是原始垂直条带梯度。
通过小波傅里叶方式对水平梯度进行平滑处理,得到平滑水平梯度,并根据平滑水平梯度和垂直梯度得到光滑的法向量。
基于光滑的法向量,利用变分模型去除原始极坐标图像中的伪影信息,得到校正极坐标图像。
将校正极坐标图像通过笛卡尔坐标变换,得到对应的校正图像,校正图像即是去除伪影的CT图像。
下面通过装置实施例来对本申请提供的CT环形伪影校正装置进行说明,如图6所示,CT环形伪影校正装置600包括:
第一变换单元601,用于对原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像;其中原始图像采用笛卡尔坐标系且包括环形伪影,其中原始极坐标图像采用极坐标系且包括与环形伪影对应的条带伪影;
确定单元602,用于确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;其中原始垂直条带梯度是指原始极坐标图像在与条带伪影垂直的方向上的梯度,原始平行条带梯度是指原始极坐标图像在与条带伪影平行的方向上的梯度;
平滑单元603,用于通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度;其中平滑垂直条带梯度包括的伪影信息小于原始垂直条带梯度包括的伪影信息;
去除单元604,用于根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;
第二变换单元605,用于对校正极坐标图像进行笛卡尔坐标变换,得到原始图像对应的校正图像;校正图像采用笛卡尔坐标系。
在一种可能的实现方式中,平滑单元603,用于:
对原始垂直条带梯度进行小波分解,得到原始垂直条带梯度对应的原始低频分量、原始对角高频分量、原始垂直高频分量和原始平行高频分量;
对原始平行高频分量进行傅里叶变换并中心化,得到原始平行高频分量对应的原始频域信息;
对原始频域信息进行低通滤波,得到原始频域信息对应的平滑频域信息;平滑频域信息包括的高频信息小于原始频域信息包括的高频信息;
对平滑频域信息进行傅里叶逆变换,得到原始平行高频分量对应的平滑平行高频分量;
基于平滑平行高频分量、原始低频分量、原始对角高频分量和原始垂直高频分量,通过小波重构得到平滑垂直条带梯度。
在一种可能的实现方式中,平滑单元603,用于:
当原始平行高频分量对应于竖直方向时,对原始平行高频分量按列进行傅里叶变换并中心化,得到频域信息。
在一种可能的实现方式中,平滑单元603,用于:
基于高斯函数,得到原始频域信息对应的衰减系数;
根据衰减系数,对原始频域信息进行系数衰减处理,得到平滑频域信息。
在一种可能的实现方式中,去除单元604,用于:
根据平滑垂直条带梯度和原始平行条带梯度,得到平滑法向量;
以原始极坐标图像作为迭代约束,将平滑法向量作为迭代目标,通过K次迭代去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;针对K次迭代中的第k+1次迭代,第k+1次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异小于第k次迭代对应的中间极坐标图像的法向量与平滑法向量之间的差异,其中,K为整数,K>1,k为整数,K>k≥1。
在一种可能的实现方式中,去除单元604,用于:
针对K次迭代中的第k+1次迭代,基于第k次迭代对应的中间极坐标图像的法向量与平滑法向量的差异,确定第k次迭代对应的适应平滑法向量项;
基于第k次迭代对应的中间极坐标图像在垂直于条带伪影的方向上的梯度,确定第k次迭代对应的垂直条带梯度惩罚项;
基于第k次迭代对应的中间极坐标图像与原始极坐标图像在平行于条带伪影的方向上的梯度差异,确定第k次迭代对应的平行条带梯度约束项;
基于第k次迭代对应的中间极坐标图像与原始极坐标图像的图像差异,确定第k次迭代对应的图像约束项;
根据第k次迭代对应的适应平滑法向量项、第k次迭代对应的垂直条带梯度惩罚项、第k次迭代对应的平行条带梯度约束项和第k次迭代对应的图像约束项,以第k次迭代对应的中间极坐标图像为基础,得到第k+1次迭代对应的中间极坐标图像。
在一种可能的实现方式中,第一变换单元601,用于
对原始图像进行极坐标变换,得到初始极坐标图像;
通过双线性差值方式对初始极坐标图像进行像素补偿,得到原始极坐标图像。
需要说明的是,本申请上述实施例提供的各个模块的具体工作过程可相应地参考上述方法实施例中的相应的实施方式,此处不再赘述。
由此可见,针对包括环形伪影的原始图像,对采用笛卡尔坐标系的原始图像进行极坐标变换,得到原始图像对应的原始极坐标图像,原始极坐标图像采用极坐标系,通过极坐标变换,能够将原始图像中的环形伪影转换为原始极坐标图像中的条带伪影;确定原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度,其中,原始垂直条带梯度是指原始极坐标图像在与条带伪影垂直的方向上的梯度,原始平行条带梯度是指原始极坐标图像在与条带伪影平行的方向上的梯度;由于条带伪影对与条带伪影垂直的方向上的梯度有明显影响,而对与条带伪影平行的方向上的梯度基本没有影响,故可以通过小波傅里叶方式对原始垂直条带梯度进行平滑处理,得到原始垂直条带梯度对应的平滑垂直条带梯度,其中,平滑垂直条带梯度包括的伪影信息小于原始垂直条带梯度包括的伪影信息,通过小波傅里叶方式能够快速对原始垂直条带梯度进行平滑,避免平滑时多次迭代造成的资源浪费,并为了减少计算量,对原始平行条带梯度不进行平滑处理;由于平滑垂直条带梯度中的伪影信息得到了有效去除,故可以根据平滑垂直条带梯度和原始平行条带梯度,去除原始极坐标图像包括的条带伪影,得到原始极坐标图像对应的校正极坐标图像;对校正极坐标图像进行笛卡尔坐标变换,得到原始图像对应的采用笛卡尔坐标系的校正图像。通过上述装置,针对包括环形伪影的原始图像,能够在进行极坐标变换将环形伪影变换为条带伪影之后,通过小波傅里叶方式对与条带伪影垂直的方向上的梯度进行快速平滑来对环形伪影进行校正,去除效果较为直观,并且在此过程中,无需对投影数据进行额外的存储和处理,避免了环形伪影的校正过程中的资源浪费。
本申请另一实施例提供了一种电子设备,包括:
存储器,用于存储可执行指令;
处理器,用于执行存储器中存储的可执行指令时,实现本申请实施例上述方法实施例中方法。
本申请另一实施例提供了一种计算机可读存储介质,存储有可执行指令,用于被处理器执行时,实现本申请实施例上述的方法实施例中方法。
专业人员还可以进一步意识到,结合本文中所公开的实施例描述的各示例的单元及算法步骤,能够以电子硬件、计算机软件或者二者的结合来实现,为了清楚地说明硬件和软件的可互换性,在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行,取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能,但是这种实现不应认为超出本申请的范围。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本申请。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本文中所定义的一般原理可以在不脱离本申请的精神或范围的情况下,在其它实施例中实现。因此,本申请将不会被限制于本文所示的这些实施例,而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。
Claims (9)
1.一种CT环形伪影校正方法,其特征在于,所述方法包括:
对原始图像进行极坐标变换,得到所述原始图像对应的原始极坐标图像;其中所述原始图像采用笛卡尔坐标系且包括环形伪影,所述原始极坐标图像采用极坐标系且包括与所述环形伪影对应的条带伪影;
确定所述原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;其中所述原始垂直条带梯度是指所述原始极坐标图像在与所述条带伪影垂直的方向上的梯度,所述原始平行条带梯度是指所述原始极坐标图像在与所述条带伪影平行的方向上的梯度;
通过小波傅里叶方式对所述原始垂直条带梯度进行平滑处理,得到所述原始垂直条带梯度对应的平滑垂直条带梯度;其中所述平滑垂直条带梯度包括的伪影信息小于所述原始垂直条带梯度包括的伪影信息;
所述通过小波傅里叶方式对所述原始垂直条带梯度进行平滑处理,得到所述原始垂直条带梯度对应的平滑垂直条带梯度,包括:
对所述原始垂直条带梯度进行小波分解,得到所述原始垂直条带梯度对应的原始低频分量、原始对角高频分量、原始垂直高频分量和原始平行高频分量;对所述原始平行高频分量进行傅里叶变换并中心化,得到所述原始平行高频分量对应的原始频域信息;对所述原始频域信息进行低通滤波,得到所述原始频域信息对应的平滑频域信息;所述平滑频域信息包括的高频信息小于所述原始频域信息包括的高频信息;对所述平滑频域信息进行傅里叶逆变换,得到所述原始平行高频分量对应的平滑平行高频分量;基于所述平滑平行高频分量、所述原始低频分量、所述原始对角高频分量和所述原始垂直高频分量,通过小波重构得到所述平滑垂直条带梯度;
根据所述平滑垂直条带梯度和所述原始平行条带梯度,去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像;
对所述校正极坐标图像进行笛卡尔坐标变换,得到所述原始图像对应的校正图像;所述校正图像采用笛卡尔坐标系。
2.根据权利要求1所述的方法,其特征在于,当所述原始平行高频分量对应于竖直方向时,所述对所述原始平行高频分量进行傅里叶变换并中心化,得到所述原始平行高频分量对应的原始频域信息,包括:
对所述原始平行高频分量按列进行傅里叶变换并中心化,得到所述频域信息。
3.根据权利要求1所述的方法,其特征在于,所述对所述原始频域信息进行低通滤波,得到所述原始频域信息对应的平滑频域信息,包括:
基于高斯函数,得到所述原始频域信息对应的衰减系数;
根据所述衰减系数,对所述原始频域信息进行系数衰减处理,得到所述平滑频域信息。
4.根据权利要求1所述的方法,其特征在于,所述根据所述平滑垂直条带梯度和所述原始平行条带梯度,去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像,包括:
根据所述平滑垂直条带梯度和所述原始平行条带梯度,得到平滑法向量;
以所述原始极坐标图像作为迭代约束,将所述平滑法向量作为迭代目标,通过K次迭代去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像;针对所述K次迭代中的第k+1次迭代,所述第k+1次迭代对应的中间极坐标图像的法向量与所述平滑法向量之间的差异小于所述第k次迭代对应的中间极坐标图像的法向量与所述平滑法向量之间的差异,其中,K为整数,K>1,k为整数,K>k≥1。
5.根据权利要求4所述的方法,其特征在于,所述以所述原始极坐标图像作为迭代约束,将所述平滑法向量作为迭代目标,通过K次迭代去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像,包括:
针对所述K次迭代中的第k+1次迭代,基于所述第k次迭代对应的中间极坐标图像的法向量与所述平滑法向量的差异,确定所述第k次迭代对应的适应平滑法向量项;
基于所述第k次迭代对应的中间极坐标图像在垂直于所述条带伪影的方向上的梯度,确定所述第k次迭代对应的垂直条带梯度惩罚项;
基于所述第k次迭代对应的中间极坐标图像与所述原始极坐标图像在平行于所述条带伪影的方向上的梯度差异,确定所述第k次迭代对应的平行条带梯度约束项;
基于所述第k次迭代对应的中间极坐标图像与所述原始极坐标图像的图像差异,确定所述第k次迭代对应的图像约束项;
根据所述第k次迭代对应的适应平滑法向量项、所述第k次迭代对应的垂直条带梯度惩罚项、所述第k次迭代对应的平行条带梯度约束项和所述第k次迭代对应的图像约束项,以所述第k次迭代对应的中间极坐标图像为基础,得到所述第k+1次迭代对应的中间极坐标图像。
6.根据权利要求1所述的方法,其特征在于,所述对原始图像进行极坐标变换,得到所述原始图像对应的原始极坐标图像,包括:
对所述原始图像进行极坐标变换,得到初始极坐标图像;
通过双线性差值方式对所述初始极坐标图像进行像素补偿,得到所述原始极坐标图像。
7.一种CT环形伪影校正装置,其特征在于,所述装置包括:
第一变换单元,用于对原始图像进行极坐标变换,得到所述原始图像对应的原始极坐标图像;其中所述原始图像采用笛卡尔坐标系且包括环形伪影,所述原始极坐标图像采用极坐标系且包括与所述环形伪影对应的条带伪影;
确定单元,用于确定所述原始极坐标图像对应的原始垂直条带梯度和原始平行条带梯度;其中所述原始垂直条带梯度是指所述原始极坐标图像在与所述条带伪影垂直的方向上的梯度,所述原始平行条带梯度是指所述原始极坐标图像在与所述条带伪影平行的方向上的梯度;
平滑单元,用于通过小波傅里叶方式对所述原始垂直条带梯度进行平滑处理,得到所述原始垂直条带梯度对应的平滑垂直条带梯度;其中所述平滑垂直条带梯度包括的伪影信息小于所述原始垂直条带梯度包括的伪影信息;
所述平滑单元,具体用于对所述原始垂直条带梯度进行小波分解,得到所述原始垂直条带梯度对应的原始低频分量、原始对角高频分量、原始垂直高频分量和原始平行高频分量;对所述原始平行高频分量进行傅里叶变换并中心化,得到所述原始平行高频分量对应的原始频域信息;对所述原始频域信息进行低通滤波,得到所述原始频域信息对应的平滑频域信息;所述平滑频域信息包括的高频信息小于所述原始频域信息包括的高频信息;对所述平滑频域信息进行傅里叶逆变换,得到所述原始平行高频分量对应的平滑平行高频分量;基于所述平滑平行高频分量、所述原始低频分量、所述原始对角高频分量和所述原始垂直高频分量,通过小波重构得到所述平滑垂直条带梯度;
去除单元,用于根据所述平滑垂直条带梯度和所述原始平行条带梯度,去除所述原始极坐标图像包括的条带伪影,得到所述原始极坐标图像对应的校正极坐标图像;
第二变换单元,用于对所述校正极坐标图像进行笛卡尔坐标变换,得到所述原始图像对应的校正图像;所述校正图像采用笛卡尔坐标系。
8.一种电子设备,其特征在于,包括:
存储器,用于存储可执行指令;
处理器,用于执行所述存储器中存储的可执行指令时,实现权利要求1至6任一项所述的CT环形伪影校正方法。
9.一种计算机可读存储介质,其特征在于,存储有可执行指令,用于被处理器执行时,实现权利要求1至6任一项所述的CT环形伪影校正方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311329514.7A CN117078791B (zh) | 2023-10-13 | 2023-10-13 | 一种ct环形伪影校正方法、装置、电子设备及存储介质 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311329514.7A CN117078791B (zh) | 2023-10-13 | 2023-10-13 | 一种ct环形伪影校正方法、装置、电子设备及存储介质 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN117078791A CN117078791A (zh) | 2023-11-17 |
CN117078791B true CN117078791B (zh) | 2024-01-12 |
Family
ID=88717482
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311329514.7A Active CN117078791B (zh) | 2023-10-13 | 2023-10-13 | 一种ct环形伪影校正方法、装置、电子设备及存储介质 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117078791B (zh) |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105321155A (zh) * | 2015-10-29 | 2016-02-10 | 北京理工大学 | 一种cbct图像环形伪影消除方法 |
WO2018103015A1 (zh) * | 2016-12-07 | 2018-06-14 | 深圳先进技术研究院 | 一种环形伪影修正的方法及装置 |
CN109636872A (zh) * | 2018-12-10 | 2019-04-16 | 合肥中科离子医学技术装备有限公司 | 基于滑动窗口差值和条带噪声检测的cbct环形伪影消除方法 |
CN111047659A (zh) * | 2019-11-08 | 2020-04-21 | 湖北科技学院 | 结合滤波方法的ct环形伪影的校正方法 |
CN115526945A (zh) * | 2022-09-27 | 2022-12-27 | 苏州一目万相科技有限公司 | 环状伪影校正算法、系统、电子设备、存储介质和芯片 |
CN116485925A (zh) * | 2023-03-24 | 2023-07-25 | 北京航空航天大学 | 一种ct图像环状伪影抑制方法、装置、设备和存储介质 |
-
2023
- 2023-10-13 CN CN202311329514.7A patent/CN117078791B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105321155A (zh) * | 2015-10-29 | 2016-02-10 | 北京理工大学 | 一种cbct图像环形伪影消除方法 |
WO2018103015A1 (zh) * | 2016-12-07 | 2018-06-14 | 深圳先进技术研究院 | 一种环形伪影修正的方法及装置 |
CN109636872A (zh) * | 2018-12-10 | 2019-04-16 | 合肥中科离子医学技术装备有限公司 | 基于滑动窗口差值和条带噪声检测的cbct环形伪影消除方法 |
CN111047659A (zh) * | 2019-11-08 | 2020-04-21 | 湖北科技学院 | 结合滤波方法的ct环形伪影的校正方法 |
CN115526945A (zh) * | 2022-09-27 | 2022-12-27 | 苏州一目万相科技有限公司 | 环状伪影校正算法、系统、电子设备、存储介质和芯片 |
CN116485925A (zh) * | 2023-03-24 | 2023-07-25 | 北京航空航天大学 | 一种ct图像环状伪影抑制方法、装置、设备和存储介质 |
Non-Patent Citations (2)
Title |
---|
Stripe and ring artifact removal with combined wavelet-Fourier filtering;Beat Munch 等;《Optics Express》;第17卷(第10期);第8567-8591页 * |
基于变分和低秩矩阵分解的CT图像环形伪影校正算法研究;吴焕堂;《中国优秀硕士学位论文全文数据库信息科技辑》(第7期);第15-31页 * |
Also Published As
Publication number | Publication date |
---|---|
CN117078791A (zh) | 2023-11-17 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Ongie et al. | Off-the-grid recovery of piecewise constant images from few Fourier samples | |
KR101598265B1 (ko) | X선 컴퓨터 단층 촬영 장치 및 화상 재구성 방법 | |
CN109840927B (zh) | 一种基于各向异性全变分的有限角度ct重建算法 | |
CN109949349B (zh) | 一种多模态三维图像的配准及融合显示方法 | |
KR20130038794A (ko) | 디지털 엑스레이 프레임들 시리즈에서의 잡음 감소의 방법 | |
CN110133556B (zh) | 一种磁共振图像处理方法、装置、设备及存储介质 | |
US11035919B2 (en) | Image reconstruction using a colored noise model with magnetic resonance compressed sensing | |
CN111553960B (zh) | 一种基于投影均值图像的环状伪影快速校正方法 | |
CN105723416B (zh) | 图像降噪方法 | |
CN107845120B (zh) | Pet图像重建方法、系统、终端和可读存储介质 | |
CN104700440B (zh) | 磁共振部分k空间图像重建方法 | |
CN112017130B (zh) | 基于自适应各向异性全变分正则化的图像复原方法 | |
Liu et al. | Rician noise and intensity nonuniformity correction (NNC) model for MRI data | |
Zhang et al. | DSPI filtering evaluation method based on Sobel operator and image entropy | |
Chen et al. | Fast linearized alternating direction minimization algorithm with adaptive parameter selection for multiplicative noise removal | |
CN111091517A (zh) | 一种残差加权成像方法和装置 | |
US8989462B2 (en) | Systems, methods and computer readable storage mediums storing instructions for applying multiscale bilateral filtering to magnetic resonance (RI) images | |
CN117078791B (zh) | 一种ct环形伪影校正方法、装置、电子设备及存储介质 | |
CN107886478B (zh) | 一种ct图像重建方法及系统、终端及可读存储介质 | |
CN113506212B (zh) | 一种改进的基于pocs的高光谱图像超分辨率重建方法 | |
Dhiyanesh et al. | Image inpainting and image denoising in wavelet domain using fast curve evolution algorithm | |
Zhang et al. | A global sparse gradient based coupled system for image denoising | |
CN110853113B (zh) | 基于bpf的tof-pet图像重建算法及重建系统 | |
Taubmann et al. | Sharp as a tack: Measuring and Comparing Edge Sharpness in Motion-Compensated Medical Image Reconstruction | |
Gao et al. | Positron emission tomography image reconstruction using feature extraction |
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 |