专利名称：计算关于局部血液流量循环指数的方法和装置的制作方法
技术领域：
本发明涉及一种关于脑组织中局部血液流量循环指数的计算方法和装置。
背景技术：
在X射线CT检查中，从简单CT图像中构造信息，和通过对比CT，在动态扫描中可以得到疾病位置周围的血液流量循环信息作为可视信息。近几年，通过多切片的高速扫描变为可能和认为对比CT动态扫描的使用范围逐渐扩大。
一般地说，存在一种称为CBP学习的方法用来计算关于脑组织中毛细管的血液流量循环指数。所述CBP学习包含得到指数例如CBP、CBV、WTT和Err，其量化表示组织中局部血液流量循环，也就是，局部组织中流经毛细管的血液循环；和输出这些指数图。
CBP表示脑组织毛细管中每单位体积和时间的血液流速[ml/100ml/min]；CBV表示脑组织中每单位体积的血液量[ml/100ml]；MTT表示毛细管的血液平均通行时间[秒]；和Err表示调制转换函数近似值中剩余误差的和或剩余误差的平方和的平方根。
CBP、CBV、MTT指数量化地表示了脑组织中毛细管的血液流量循环，与脑缺血性中风发展后的通行时间信息一起作为有用信息，用于区分缺血脑血管失调的病体，判断毛细管出现/没有扩张，或估计血液流速。例如，通常，在缺血脑血管失调中，提供的动脉血压下降，和血管内的流速降低。作为结果，甚至当CBV是常数时，MTT扩大和CBP降低。而且，在脑梗塞过急性阶段中，通过血压降低来补偿血液流速的降低，存在一种用于扩张毛细管和增加血液流速CBP的自动调整。因此，由于MTT扩大，甚至随着CBP下降、CBV增加，所述信息暗示毛细管的血管再形成的可能性。
在CBP学习中，不具有脑血管渗透性的造影剂，例如碘化造影剂用作示踪剂。经由肘脉注射碘化造影剂，例如通过注射器。通过注射器注射进静脉的碘化造影剂经由心脏和肺流进脑动脉。然后，造影剂经过脑组织中的毛细管从脑动脉中流出到脑静脉。碘化造影剂流经正常脑组织中的毛细管而没有任何血管外的泄漏。图1简要表示了这种状态。
通过动态CT扫描造影剂的通行状态，脑动脉上像素的时间-密度曲线Ca(t)、脑组织(毛细管)上像素的时间-密度曲线Ci(t)、和脑静脉上像素的时间-密度曲线Csss(t)可从连续图像中测量出。
这里，在CBP学习中，在脑动脉的时间-密度曲线Ca(t)和脑组织的时间-密度曲线Ci(t)之间建立的理想关系用作一种分析模型。假定如果就在脑组织之前经由血管注射造影剂，脑组织的单位体积(一个像素)里的时间-密度曲线垂直上升，停在一常数值上，并随后倾斜下降。这可以用矩形函数来近似(box-调制转换函数方法box-MTF方法)。
也就是，脑动脉的时间-密度曲线Ca(t)用作输入函数，用于脑组织的时间-密度曲线Ci(t)用作输出函数，和通过矩形函数来近似输入和输出函数间的调制转换函数。调制转换函数表示示踪剂经过毛细管的通行过程CBP学习具有下列问题。
由于计算各个指数CBP、CBV、MTT和Err用于每个像素(x，y，z)，可以构造使用所述值作为像素值的图像，和该图像称作为图。例如，当得到R类型的指数时，可以构造R个图。用这种方式准备的R个图可以看作为一幅图(向量值图)，其中每个像素具有一个向量值。也就是，所述图可以表示如下。
Vk(x，y，z)＝<Pk，1(x，y，z)，Pk，2(x，y，z)，…，Pk，R(x，y，z)>
例如，可以用这种方式来构造CBP学习，就是假定典型地R＝4，Pk，1(x，y，z)表示CBP值，Pk，2(x，y，z)表示CBV值，Pk，3(x，y，z)表示MTT值，Pk，4(x，y，z)表示剩余错误误差Err值。
这个向量值图Vk准备用于所述参考脑动脉的每个时间-密度曲线Ca(t)k。例如，假定从左和右半球的中间、前和后脑动脉中得到时间-密度曲线。在这种情况下，K＝6。而且，假定从影响区域的外围中动脉的几个部分中得到脑动脉的时间-密度曲线，K＝约10-15。
当脑动脉的时间-密度曲线Ca(t)k的数量K大(K＝1，2，…，K)时，作为结果得到的向量值图的数量Vk(K＝1，2，…，K)是大的，因此这是便于观察的。也就是，当所述图看作为正常灰度级图像或彩色级图像时，一个图由R个图像构成，存在K个图，和因此总共K×R个图像必须进行比较。而且，由动脉滋养的部分和所述动脉不必显出，和必须使用解剖知识来判断每个区域的要观察的图Vk(k＝1，2，…k)。尤其是在大脑血管失调发展的情况下，例如大脑梗塞，基于组织的动脉判断与解剖知识不一致，并且异常依赖关系频繁可见。这些缺陷产生的问题是难以解释向量值图的射线照相。
而且，在由多切片或量CT扫描的动态CT图像中，进一步观察到大量的动脉。这是因为在多个切片中可以观察到相同的动脉。如果脑动脉的时间-密度曲线准备用于这些动脉的所有X射线断层摄影图像，曲线的数量变得非常庞大。
而且CBP学习也有下列问题。由于经由肘脉进行丸药注射，对于用CT观察的对比增加效果来说，血液CT数上升到最大几百个HU(当不执行对比成像时，几十个HU)。然而，为了有效地分析脑血液流量，对比增强效果必须仅用几个百分点或更少误差进行测量。也就是，甚至当对比增强效果(CT数上升)约为20-40HU，必须检测对比增加效果。
单位体积的脑组织中毛细管的体积比最大约为3-4％。因此，当血液CT数上升20-40HU时，脑组织的平均CT数只上升约为0.5-1.5HU。
在CT图像中，噪声的标准偏差(sd)与X射线辐射剂量的平方根成反比，和在典型的辐射条件下，sd例如约为5-10HU。因此，为了检测0.5HU的对比增强效果，X射线辐射剂量必须增加约10-100倍，这意味着患者的曝光剂量相当大。而且，由于在动态CT中同一位置扫描几十次，扫描位置中的皮肤曝光量是正常曝光量的几百到几千倍，和考虑到辐射问题例如炎症、秃头症、坏死和致癌，这是不现实的。
而在动态CT中，与正常扫描相比，X射线辐射剂量必须减少。通常，每次扫描的X射线辐射剂量减少到例如约为正常剂量的1/2-1/10。因此，与正常一次CT扫描相比，只产生约几倍到20倍的X射线曝光量，不会出现任何辐射问题。然而，在CT图像中，其中X射线辐射剂量减少，sd例如约为15-20HU，几乎检测不到约0.5-1.5HU的对比增强效果。
因此，抑制图像的噪声成分是CBP学习中的重要问题之一。为此，1)增加切片厚度，2)平均相邻像素，和3)对图像进行滤波处理是通常采取的措施。然而，这些具有下列问题。
为了“增加切片厚度”，在扫描期间，切片厚度设定为大，或将连续薄切片的图像数量进行平均和产生厚切片图像。由于每个像素的X射线辐射剂量与切片厚度成比例地增加，因此图象噪音的sd与切片厚度成比例地减小。然而，当切片厚度增加时，产生部分体积效果。也就是，一个像素不表示相同脑组织，像素表示许多组织(白色物质、灰色物质、血管、脑沟、脑室等等)的平均CT数的概率是错误的，和作为分析结果得到的脑血液流速值变得不准确。
特别地，不可能正常分析包括血管影响的像素。因此，由于增加的切片厚度，只能得到包括许多像素的非常低质量结果，其不能用于分析。
通过对相邻的像素进行平均，空间分辨率受到一定程度的损失。例如得到正方形区域(包括n×n像素)的平均值作为整个正方形的平均CT数，正方形一边包括n个像素，该正方形看作为像素，和布置这些正方形来构成一个“像素捆绑图像”。例如，假定原始图像在一边包括512个像素(包括512×512个图像)，和n＝2，“像素捆绑图像”在一边由(512/2)个像素构成(该图像包括256×256个像素)。根据这种方法，与n成反比减少噪声。而且，作为分析目标的像素数量通过1/(n×n)的因数来增加，和因此具有计算量减少的优点。
然而，当n增加时，空间分辨率降低，并因此出现部分体积效果。也就是，一个像素不表示相同的脑组织，像素表示许多组织(白色物质、灰色物质、血管、脑沟、脑室等等)的平均CT数的概率变得不准确，和作为分析结果得到的脑血液流速值变得不准确。特别地，不可能正常分析包括血管影响的像素。因此，由于增加的n，空间分辨率低，和因此，只能得到包括许多像素的不准确和非常低质量的结果，其不能用于分析。因此，在实际使用中，n＝约2-4是极限，只使用这个措施，不能得到充分的噪声抑制效果。
而且，当图像被平滑时，也就是当使用对每个CT图像进行二维空间滤波和滤波所述图像的方法，得到充分的噪声抑制效果时，但是空间分辨率明显削弱。特别地，存在厚血管(动脉/静脉)的区域附近中的像素受到对比增加效果的影响，其产生在厚血管中，这些像素的时间-密度曲线是不准确的。因此，只能略微地执行滤波。这里，在只执行轻微滤波中，重要的是显著地减少图像滤波器的尺寸，例如，约为3×3。当使用3×3平滑滤波器得到最大图像噪声抑制效果时，上限是减少噪声sd到1/3，和不可能进一步抑制噪声。因此，不能得到充分的噪声抑制效果。
另一方面，当沿着时间轴滤波时，也就是，使用把为每个像素得到的时间-密度曲线看作为一个曲线和用一维滤波器滤波所述曲线的方法，和得到充分的噪声抑制效果，时间分辨率明显削弱。在CBP学习中，通过在短采样周期内初次执行动态CT来得到高时间分辨率，和精确地测量时间-密度曲线的微小和快速变化(尤其是产生于生理结构的滤波效果程度)，和随时间滤波是不适当的。
发明内容
本发明的一个目的是提高通过CBP学习得到的图的诊断效率。
本发明的另一个目的是抑制噪声而不降低空间和时间分辨率，和由此提高CBP学习的分析精度。
根据本发明的第一方面，提供一种方法从涉及注射有造影剂的受试者特定区域的许多连续图像中，准备关于特定区域中动脉的第一时间-密度曲线和关于特定区域中组织的第二时间-密度曲线；通过曲线拟合来计算调制转换函数，其表示组织中关于动脉的局部血液流量循环，由此最小化第二时间-密度曲线关于调制转换函数与第一时间-密度曲线卷积的剩余误差；从调制转换函数中计算关于各个动脉的局部血液流量循环的指数；准备用于所述动脉的指数图；和根据第一时间-密度曲线的剩余误差，把所述指数图合成为一个图。
根据本发明的第二方面，提供一种方法从涉及注射有造影剂的受试者特定区域的许多连续图像中，准备关于特定区域中动脉的第一时间-密度曲线和关于特定区域中组织的第二时间-密度曲线；选择一个与各个第二时间-密度曲线最相称的第一时间-密度曲线，来指定一根动脉，其具有最高可能性就是所述组织中每个局部血液流量循环是依赖的；和为每个组织准备一张图，来根据指定的一根动脉来区分动脉的非独立区。
根据本发明的第三方面，提供一种方法根据像素间的类似性使用权来对许多图像进行滤波，来减少来自涉及注射有造影剂的受试者特定区域的多个连续图像的噪声；从具有减少噪声的许多连续图像中，准备关于特定区域中动脉的第一时间-密度曲线和关于特定区域中组织的第二时间-密度曲线；通过曲线拟合来计算调制转换函数，其表示组织中关于各个动脉的局部血液流量循环，由此最小化第二时间-密度曲线关于调制转换函数与第一时间-密度曲线卷积的剩余误差；和从调制转换函数中计算关于各个动脉的局部血液流量循环的指数。
本发明的其它目的和优点将在其后的说明书中阐述，和在说明书中部分将是显而易见的，或可以通过本发明的实践来获知。依靠下文中提出的手段和结合可以实现和得到本发明的目的和优点。
与说明书组合并构成说明书的一部分的这些附图，图解了本发明目前优选实施例，和与上文给出的概括说明和下文给出的优选实施例的详细说明一起，用来解释本发明的原理。
图1是CBP学习的原理说明图；图2是表示根据本发明实施例关于脑组织中毛细管的血液流量循环的指数计算装置构造的框图；图3A、3B、3C是由本发明的相干滤波器执行的图像处理的说明图；图4A、4B是表示由本发明中相干滤波器执行的噪声抑制流程的流程图；图5是本实施例中整个指数计算处理前半段的流程图；图6是本实施例中整个指数计算处理后半段的流程图；图7是表示图5步骤S3划分线的一个实例的图；图8是表示在Voronoy图上每根脑动脉的影响范围的图，用于减少图6步骤S12的处理步骤的数量；图9A、9B、9C表示图5步骤S4的AT、PT和TT图；图10是表示图5步骤S4的AT、PT、TT的曲线图；图11是表示脑动脉ROI的曲线图，在图5步骤S6中的切片中相同；图12是表示在图5步骤S7中设置的上纵向窦性静脉ROI的图；图13是关于图5步骤S10的脑动脉的时间-密度曲线修正的补充图；图14A、14B表示在图5步骤S10、S11中准备的脑动脉的时间-密度曲线Ca(t)和用于脑组织的时间-密度曲线Ci(t)的一个实施例；图15是图6步骤S12的box-MTF方法的原理说明图；图16是图6步骤S12的box-MTF处理的说明图；图17表示设定图6步骤S14的每个指数的屏幕输出范围的一个实例；图18A至18D是表示在图6步骤S16中准备的每个CBP、CBV、MTT、Err图的一个实例的图；
图19表示在图6步骤S16中准备用于每根脑动脉的一列CBP、CBV、MTT和Err图；图20是图6步骤S17的图谱合成方法的说明图；图21A至21D是表示在图6步骤S19中计算的平均值显示实例的图；图22表示图(控制图)的产生方法，其中在图6步骤S17中具有高优势可能性的脑动脉(ACA、MCA、PCA)由每个像素(局部组织)的标签来识别；图23是表示由图22的产生方法产生的控制图的实例图；和图24A至24C表示使用控制图过滤的CBP图。
具体实施例方式
下文中将参考附图描述本发明的优选实施例。
本实施例的特征在于一种方法把许多指数图合成一个图，指数图通过CBP学习而产生，由此提高了CBP学习的诊断效率；并进一步使用相干滤波器来减少噪声和禁止空间与时间分辨率的降低以便提高指数的精确度。
而且，本实施例涉及一种方法和装置，其中从多个图像中计算表示局部血液流量循环的指数，该图像是关于物体的特殊区域并和时间有关，和产生多个图像作为目标的形式不局限于一个专门装置。该装置的实例包括X射线计算机断层照相术装置(X射线CT装置)、单光子发射照相术装置(SPECT)、正电子发射照相术装置(PET)、和磁共振照相装置(MRI)。这里，作为实例将描述X射线CT装置。
(装置构造)图2表示根据本实施例X射线CT装置的构连X射线CT装置由一构台部件10和计算机装置20构成。构台部件10包括一X射线管101、高压发生装置101a、X射线检测器102和数据获得系统(DAS)103。X射线管101和X射线检测器102经由物体P相互相对布置在旋转环(未示)上，旋转环连续高速旋转。
计算机装置20包括图像处理装置30、图像显示部件107和输入部件109。图像处理装置30包括一控制部件108作为中央单元，和还包括一预处理部件104用于通过修正处理，把从数据获得系统103中输出的原始数据转换成投影数据；一内存部件105用于保存投影数据；一图像重新构成部件106用于从投影数据中重新组成CT图像数据；一存储装置10M用于保存CT图像数据；一相干滤波处理部件110用于执行关于CT图像数据的相干滤波处理；和一CBP学习处理部件120用于使用进行了相干滤波处理的CT图像数据来执行CBP学习处理。
相干滤波处理部件110包括一离散值估计部件111、加权函数计算部件112和像素值计算部件(相干滤波部件)113。在下文相干滤波处理的详细说明中将描述这些离散值估计部件111、加权函数计算部件112和像素值计算部件113的功能。
CBP学习处理部件120包括一ROI设定支持部件121、时间-密度曲线预备部件122、脑动脉时间-密度曲线修正部件123、MTF处理部件124、指数计算部件125、图预备部件126和图合成部件127。
ROI设定支持部件121准备和提供信息(用于脑动脉ROI的AT、PT和TT图)用于支持在CT图像上设定关于脑动脉和静脉的ROI兴趣区域。
而且，脑动脉ROI分散设定在左脑和右脑区域，例如，关于前脑动脉(ACA)、中间脑动脉(MCA)、和后脑动脉(PCA)作为目标。因此，在该实例中，每侧三个，也就是，总共设定六个脑动脉ROI。而且，为了修正脑动脉的时间-密度曲线Ca(t)，使用另一个时间-密度曲线Csss(t)。该时间-密度曲线Csss(t)相对于设定在足够厚血管上的ROI来准备，在该血管中存在不包括部分体积的像素。设定Csss(t)的ROI，例如，在脑动脉中的最厚的上纵向窦性静脉。
时间-密度曲线预备部件122从保存在存储装置10M中的动态CT图像数据(多个时间连续图像数据)中准备关于脑动脉、大脑静脉和脑组织(毛细管)的时间-密度曲线。而且，分开准备脑动脉的时间-密度曲线Ca(t)，例如，用于六套脑动脉ROI中的每套。相对于设定在上纵向窦性静脉中的大脑静脉ROI，准备时间-密度曲线Csss(t)。而且，为脑组织上的所有像素中的每个准备脑组织的时间-密度曲线Ci(t)作为目标。
为了消除噪声影响和部分体积影响，脑动脉时间-密度曲线修正部件123根据上纵向窦性静脉的时间-密度曲线Csss(t)来校正脑动脉的时间-密度曲线Ca(t)。后面将描述该修正方法。MTF处理部件124根据通过box-MTF方法修正的脑动脉的时间-密度曲线Ca(t)和脑组织的时间-密度曲线Ci(t)来计算调制传递功能MTF，用于脑组织区域中所有像素中的每个作为目标。
指数计算部件125从调制传递功能MTF中计算表示血液流量循环的指数(CBP、CBV、MTF、Err)，用于脑组织区域中所有像素中的每个作为目标。图预备部件126产生计算的指数的图用于每根脑动脉(ACA、MCA、PCA)。关于每个切片，通过指数的类型(＝4)X脑动脉的数量(3，包括ACA、MCA、PCA)＝12种类型，来准备图，在多切片中，准备切片数量乘以图的类型。布置图合成部件127通过合成处理来减少图的巨大数量，并提高诊断效率。
在下文中将依次描述相干滤波处理和CBP学习处理。
以下简要描述相干滤波器的原理。例如，3×3相邻的局部像素被加权平均，和这个加权平均值用作局部中心像素值。根据中心像素和外围像素之间的类似性来改变外围像素的权。这里，类似性是一种表示可能性程度的指数，即在相同脑动脉的控制下，在解剖学上像素接近组织，具体地是脑组织(毛细管)。把高加权给具有高类似性的像素，和相反把接近零的低加权给具有低类似性的像素，由此抑制了噪声，空间分辨率的降低得到抑制。而且，在下文中，类似性可以适当地用保真度或风险率来代替。
在本实施例中，在、本课题的大脑，其中造影剂不具有大脑血管通透性，例如，注射(静脉注射)的碘化造影剂用作扫描目标，和多个连续获得的CT图像(动态CT图像)用来通过比较各个像素的时间-密度曲线来计算类似性。因此，类似性的确定依赖于和由采样频率、每单位时间里的图像数量和采样数量也就是所有图像的数量来确定。因此，有效的是把扫描间隔减少到例如0.5秒。
(相干滤波器)(相干滤波器的概述)(像素值V(X))通常，经由扫描设备例如照相机和CT扫描仪获得的数字图像由多个像素构成(或图像可以被认为是一组像素)。在下面说明中，像素的位置表示为向量X(也就是，一个坐标值的向量)，且像素X的值(也就是，表示暗纹、CT数HU的数字值)被表示为K维向量。对于二维图像，像素X是表示代表图像位置的坐标值(x，y)的二维向量。关于某一像素X定义的“像素值v(x)”表示如下v(x)＝(v1(x)，v2(x)， …，vk(x))…(1)其中在方程式(1)右侧中的v1(x)，v2(x)，…，vk(x)在下文中将称为关于像素X的“标量值”。
例如，当图像是“彩色图像”时，每个像素具有三个主要颜色(红、绿、蓝)的亮度，和因此每个像素的像素值v(x)可以看作为K＝3维向量。也就是，方程式(1)右侧中的下标与例如“红”、“绿”、“蓝”相关。而且，例如，当图像是由K静态图像构成的动态图像时，和第n个图像的像素具有标量值vn(x)，由相同公共点(相同坐标)的像素X拥有的像素值(标量值)布置在K静态图像上来构成K维向量值vn(x)＝(v1(x)，v2(x)，…，vk(x))，和这个值是像素值作为下文中描述的向量值。
(类似性(保真度或风险率)和加权)一组N(x)适当的外围像素被认为是关于像素X的(该组N(x)包括像素X)。因此，认为外围像素Y的权w(p(x，y))作为关于中心像素X的N(x)元素。这个权w(p(x，y))具有下列特性。
(类似性p(x，y))首先，将描述影响权w(p(x，y))值的函数p(x，y)的意义。该p(x，y)的意义是用于量化本实施例中提到的“类似性”。通常，这个值表示一个代表中心像素X和外围像素YεN(x)之间某种意义上的类似性程度的具体数字值。(例如，两个像素x和y的像素值v(x)和v(y)之间的公认的统计差程度)。
更具体地，例如，对于表示大值的p(x，y)来说，像素x和y具有“统计上非显著性差异(也就是高类似性)”和相互类似的可能性认为是高的。对于表示小值的p(x，y)来说，认为像素x和y“具有统计上显著性差异(也就是低类似性)”。
而且，像素值v(x)和v(y)(或标量值v1(x)，…，vk(x)和v1(y)，…，vk(y))必然包括噪声。例如，假定图像由CCD图像获得装置获取，装置中的暗电流产生的噪声和从外界入射的光量的不规则波动存在于构成图像的各个像素。
对于所有像素来说，噪声通常具有各种值。因此，即使像素x和y反映相同的物体(在外界)，某些情况下在实际观察的图像上的像素具有不相同的值。相反，假定在反映同一物体的像素x和y中，环境远离各个噪声，这些显示(＝公认为这样)在作为同一物体表示的图像上，或二者实质上具有相同(非常接近)像素值。
然后，考虑到上述噪声特性，和使用在统计检查方法中公知的“虚假设”观念。然后，这个像素p(x，y)具体表示如下。也就是，虚假设H“当消除各个噪声时，像素x和y具有相同像素值”，换句话说“v(x)＝v(y)，这里消除了产生于两个像素噪声的差”(也就是，当确立了这个建议时，假定“像素x和y之间的类似性假定是高的(保真度大)。然后，当否决虚假设H时，函数p(x，y)可以构成为风险率(或显著水平)(在这种情况下，p(x，y)定义为一个函数，使得值范围是
(p(x，y)ε
))。
因此，当风险率p(x，y)大时，也就是，具有高风险率的拒绝是错误的，虚假设H满足的可能性高。然而，相反，当风险率小时，也就是，具有低风险率的拒绝是错误的，虚假设H不能满足的可能性高(应当注意，由于统计确认中公知的问题，甚至具有没有“拒绝”的虚假设H，这并不意味着虚假设是“真实的”。在这种情况下，这只意味着由虚假设H表示的建议不能被否决)。
(权w(p(x，y)))由于权w(p(x，y))表示方法是显而易见的，如上所述权w(p(x，y))是风险率p(x，y)的函数。通常，权是保真度的函数。假定保真度是p(x，y)，函数可以构成为w(ρ(x，y))。而且，通常来讲，作用于风险率p(x，y)的权函数W具有实现“拒绝”的作用，为了得到权w(p(x，y))，关于x和y的组合得到风险率p(x，y)。相反，当风险率p(x，y)大时，权函数W，也就是权w(p(x，y))表示大的正值。在相反情况下，调整函数具有小正值(或“0”)。将在下文中描述权函数W的具体形式。也就是，当像素x和y满足由假设H表示的命题时，权w(p(x，y))表示大值。在相反的情况下，权具有小值。作为一个实例，可以显著地构成W表示两个值；“0”；和不同于“0”的常数值。
进一步地，将概述上述的假设H、风险率p(x，y)和权w(p(x，y))间的关系。当虚假设H具有正确的高可能性，类似性P也增加，和提高了给像素的权W。另一方面，虚假设H是正确的可能性低，因而类似性P是低的，和降低了给像素的权W。当根据类似性以这种方式改变加权平均值的贡献(权)时，抑制分辨率的降低和有效抑制噪声是可能的。而且，通常，权函数w(t)可以是“由tε
定义的正/负单调增函数”，和w(t)可以满足至少上述特性。
(相干滤波处理)通过上文说明，得到“相干滤波器”如下。也就是，首先关于构成图像的某一像素X来计算组N(x)元素的所有像素的上述权w(p(x，y))。因此，在下面方程式(2)中，使用这些权w(p(x， y))来计算构成像素X的新标量值V′k(x)。v′K(x)=Σy∈N(x)vk(y)w(p(x,y))Σy∈N(x)w(p(x,y))---(2)]]>
其中k＝1，2，…，k。因此由该公式得到的V′k(x)用来构成像素x的修改的像素值(新像素值)V′(X)如下。
V′(X)＝(V′1(x)，V′2(x)，…，V′k(x)) …(3)这里，用于转换像素值V(y)＝(V1(y)，V2(y)，…，Vk(y))(包括y＝x的情况)为V′(X)＝(V′1(x)，V′2(x)，…， V′k(x))的滤波器是“相干滤波器”的形式。该公式显而易见，这表示构成像素值的标量值Vk(y)的加权平均值。
该处理产生下列结果。也就是，这表示由加权平均值V′k(x)构成的向量，V′k(x)具有升高的像素y贡献，像素y似乎确实具有与像素x的像素值V′(X)相同的像素值，像素x的噪声消除(＝其具有满足假设H命题的高可能性)。而且，如果存在足够数量的像素y，像素值V′(X)具有不降低分辨率的值，而不偏离由像素x实质处理的真实值，并通过上述平均作用来抑制噪声。
进一步地，甚至当风险率p(x，y)小时，因此“拒绝”虚假设H，和减少权w(p(x，y))，正如从上面说明中所看到的，假设不是必须完全“拒绝”。这依赖于下文描述的权函数w的具体形式，但是甚至当风险率p(x，y)接近“0”时(＝0％)，可以设定w(p(x，y))≠0。然而，与p(x，y)接近“1”相比，权具有更小的正值。应当注意p(x，y)＝1相应于下文描述的v(x)＝v(y)。
也就是，可以承认小贡献，而不是完全拒绝该假设。应当注意，在这种情况下，w(p(x，y))＝0具有与完全拒绝相同的意义。
该处理通常进行如下。也就是，当存在构成某一图像的多个像素时，具有某一任意像素y(在上面说明中设定设定y∈N(x))的该像素x的保真度p(x，y)得到量化。此时量化上述h。当保真度大时，在使用像素值v(y)的加权平均处理中承认像素Y的大贡献。当保真度小时，只承认小贡献。因此，可以是在图像处理方法中有效抑制像素X的噪声。当像素x和y是所谓的“相互类似”时，允许像素y进一步进行平均处理。换句话说，当像素是“非相互类似”时，全部或几乎忽略像素y(权设定为零或近似值)。
当整个图像进行该处理时，图像的模糊，也就是空间分辨率的降低几乎不出现，和可以实现显著高的噪声抑制效果。而且，该使用不局限于噪声抑制。例如，甚至在图案识别领域中，以更好的具体形式来设定权函数或相干滤波。然后，可以实现更好的效果。
这里，上述“动态CT”扫描是一种扫描系统，其中X射线管101和X射线检测器102重复扫描物体P的相同部分(经常通过重复扫描中的连续旋转和连续旋转CT装置来进行重复扫描)来连续获得投影数据，根据投影数据连续执行还原处理，和获得序列图像。在这种情况下，图像显示部件107的图像显示受到控制，例如，通过计数器(未图示)，由此在扫描开始或终点之后，执行恒定时间的显示以用于收集投影数据作为原始图像。
因此，以这种方式获得/显示的图像是所谓的动态图像，其包括类似于电影中的多个序列静态图像。应当注意这种扫描系统典型地用于把造影剂注入进物体P，观察/分析随时间消逝的变化，和分析受影响区域的病理状态，例如组织中的狭窄和闭塞。而且，甚至在使用造影剂之前和之后，仅执行两次相同部分的CT扫描的系统认为是广泛的动态CT扫描。
而且，在传统技术中，在上述“动态CT”扫描期间，例如，当物体P在K次扫描操作期间具有某些变化时(例如，通常考虑到的是造影剂的浓度变化或呼吸操作)，为了抑制图像噪声而不损失空间分辨率，在时间方向上必须执行滤波。因此，损失时间分辨率的问题不能避免。
然而，由动态CT扫描获得的图像是如上所述的动态图像，并被扫描用于详细观察随时间变化的目的。因此，时间分辨率的损失实质上不是优选的情形。
由于使用相干滤波器，可以执行下面的动态相干滤波处理，其中不损失时间分辨率和可以从所有K静态图像(多个图像)中抑制噪声。
首先，对于像素X来说，其关于K静态图像定义为得到的如上所述的动态图像的，下式可以构成为如上所述的像素值v(x)。
v(X)＝(v1(x)，v2(x)，…，vk(x)) …(1)这里，在右侧项中的下标1，2，…，k是K静态图像的序号。
接下来，权函数w1的具体形式在这种情形下例如，用下面方程式(4)表示。wl(p(x,y))=exp[-{Σk=1k{vk(x)-vk(y)}2K(2σK)2}C]---(4)]]>其中y∈N(x)，和这个组N(x)可以关于像素x(＝也可以通过任何标准来设定)任意设定。然而，事实上，通过像素x和远离像素x的像素y满足假设“v(x)＝v(y)，这里由两个像素噪声产生的差得到消除”的可能性通常可能是低的。因此，在假设基础上限定组N(x)，假设组N(x)是一组布置在像素x附近的像素，具有例如计算速度的提高的实际意义。
因此，这里，作为一个实例，假定组N(x)是一组包括在以像素x为中心的矩形外围区域里的像素。更具体地，例如，对于组N(x)来说，当128×128个像素整体构成一个显著静态图像时，定义用于以像素x为中心的3×3像素区域。而且，当512×512个像素构成静态图像时，可以定义用于以像素x为中心的13×13像素区域。
而且，上面方程式(4)中的σk是在假设每个K图像具有恒定的普通程度的噪声的基础上，所估计的噪声的标准偏差，假设另一方面，C是一个可调节参数，其确定分配给上面方程式(4)的权w1(p(x，y))的作用。在下文中依次描述σk和C。
首先，将描述方程式(4)中的σk(下文中描述为散布的σk2)。如上所述，这个σk2是由K静态图像上每个像素的标量值所拥有的散布噪声组份。而且，上面方程式(4)中散布的σk2在假设基础上得到估计，假设K图像上每个像素的标量值包括具有恒值散布σk2的噪声。通常，这个假设在下列背景中具有充分的合理性。
当受试者P的尺寸和X射线管101、X射线检测器102和重构部件106的结构不变，和X射线的辐射能设定为常数，CT图像的噪声由X射线辐射剂量来确定，也就是，X射线管101的管电流和发光时间的乘积成比例关系(所谓的电流时间乘积(mA·s))。
另一方面，CT图像的噪声是附加的和实质上遵循高斯分布也是公知的。也就是，假定实际值(该值具有消除的噪声贡献)是Vn0(x)，关于构成某一像素x的像素值v(x)的任意标量值Vn(x)(n＝1，2，…，k)，Vn(x)-Vn0(x)的差值实质上遵循高斯分布，具有0和σk2的中间值。应当注意X射线剂量或电流时间乘积mAs实质上与噪声的散布σk2成反比。
而且，散布σk2依赖于像素x的位置(例如，如上所述每个坐标值x＝(x，y))。在正常的X射线CT扫描仪100中，用于调节X射线辐射剂量的物理X射线滤波器(例如，称为“楔”或“X射线滤波器”，由铜箔或固体金属片形成)布置在X射线管101和X射线检测器102之间。因此可以忽略散布。楔具有以下功能，物体P用实质上与水的密度相同的材料制造，和递归一部分X射线辐射剂量，由此在任何X射线检测器102中检测到相同程度的X射线辐射剂量。因此，根据该楔，作为结果，产生的效果是噪声的散布σk2实质上具有恒值，不管像素x的位置。通常，布置该楔是为了有效地使用X射线检测器102的动态范围的基本目的。
如上所述，它是适当的就是在通过动态CT扫描获得的K静态图像中估计散布σk2实质上关于K静态图像上所有像素是恒定的。当然，很容易想到扩展本实施例用于每个像素的散布差。
进一步地，为了具体计算上面方程式(4)，分配给散布σk2的数字值是一个问题。这是因为通常可以假定噪声分布的形式(上面是高斯分布)，但是在许多情况下散布σk2的具体值是不清楚的。
而且，通常，每次扫描可以改变辐射剂量(X射线管电流X辐射时间(mAs))。当用这种方式进行扫描时，在这种情况下，散布σk2与每个图像不一致。
另外，假定在第k个图像(k＝1，2，…，k)中每个像素标量值的噪声散布是σk2和在第k个图像的扫描中使用的辐射剂量是Rk，σk2与Rk成比例。因此，当关于至少一个k＝k0指定σk02时，关于另一个k也可精确地估计σk2如下。σk2=σk02RkRk0---(5)]]>在本实施例中，通过下面方法，可以关于至少一个k估计σk2的具体数字值。
该方法是有效的，就是使用N个图像(1＜N≤k)，假定该图象在K次扫描操作中物体P几乎不变化，和通过实际测量得到关于散布σk2的期望值E[σk2]。为了简化下面的说明，在N个图像中的辐射剂量是相同的。因此，假定σk2关于k＝1，2，…，N是常数(写成σ2)。期望包括在各自的标量值v1(xf)，v2(xf)，…，vk(xf)中的噪声遵循高斯分布，具有如上所述的0和σ2散布的平均值，标量值构成这些N个图像中某一像素的像素值v(xf)。使用下面方程式(6)得到平均值v*(xf)=1NΣK=1Nvk(xf)---(6)]]>然后，关于实际散布σ2可以得到期望值E[σ2]如下。E[σ2]=1N-1ΣK=1N{vk(xf)-v*(xf)}2---(7)]]>而且，可以认为如上所述散布的期望值E[σ2]关于所有K静态图像上的所有像素x是适当的。期望值用于替代实际散布σ2可能性通过不小于常数程度得到保证。因此，在上面方程式(4)的实际计算中，这个E[σ2]可以指定给方程式(4)的σ2。
进一步说，通过实际测量可以得到更具体的E[σ2]，例如，根据K静态图像中的第一和第二静态图像，这符合在上面方程式(6)和(7)中的N＝1的假设。而且，对于在上面方程式(6)和(7)的实际计算中使用的像素xf，来说，例如，只可以设计选择适当像素xf，其排斥一部分，在该部分中扫描了空气和骨头。当选择多个像素时，平均所有得到的E[σ2]值。而且，通常，通过物体P的移动也可以设计更好地抑制影响。
甚至由于这些N个图像扫描中的易变的辐射剂量，很容易推断σk2对Rk的比例用来准确地估计σk2。
随后，将描述上面方程式(4)中的参数C。首先，在方程式(4)中，在上述平常方式中的风险率p(x，y)的概念包括如下。也就是，在方程式(4)的右侧分子中的根号里面与假定遵循所谓的X平方分布的X2值一致，并除以(2σ)2。整个圆括号放在e的肩部，和该值是风险率p1(x，y)。也就是，下面方程式的结果。p1(x,y)=Aexp{-ΣK=1K{vk(x)-vk(y)}2K(2σk)2}---(8)]]>而且，上面方程式(4)转换为关于p1(x，y)并表示在方程式(8)中的下列方程式。w1(p(x,y))={p1(x,y)}C2---(9)]]>其中A是标准化的常数，由此p1表示一个值(0至1)。
最后，在方程式(4)中，上面以正常模式描述的风险率p(x，y)没有肯定地显示，但是如上所述，权w1(p(x，y))的模式的确可以看作为风险率(＝p1(x，y))的函数(方程式(9))，也就是，模式是“保真度函数”。如上所述，风险率和保真度是这样一种关系就是随着它们中的一个增加，另一个也增加。
而且，如从上面方程式(9)中所看到的，参数C有效地确定权w1(p(x，y))对风险率p1(x，y)的灵敏反作用程度。也就是，当C增加时，p1(x，y)轻微增加，和w1(p(x，y))接近0。进一步说，当C减少时，可以抑制这种灵敏反作用。而且，C可以具体设定约为1-10，和优选地C设定为3。
在本实施例中，可判断中心像素x和外围像素y之间的类似性。换句话说，通过一种所谓的上文中显而易见的X平方检查方法(统计检查方法)，根据风险率p1(x，y)来确定关于两个像素x和y的上述虚假设H的否定。
而且，正如从上面方程式(4)中所看到的，在本发明中，不是必须执行计算关于各个x，y组合的风险率p(x，y)和随后得到权w(p(x，y))的程序。没有具体得到风险率p(x，y)时，在构造中可以直接计算综合功能(wop)。
如上所述，当估计散布σ2时(例如方程式(7)的E[σ2])，和适当地确定参数C(例如C＝3)，可以使用方程式(4)得到关于所有像素y的具体权w1(p(x，y))，像素y包括在关于某一像素x定义的组N(x)中(例如上述以像素x为中心的3×3像素区域)。然后，在方程式(2)中代替w(p(x，y))，可以使用这个w1(p(x，y))来执行相干滤波器的具体数字值计算。而且，作为结果，没有削弱时间分辨率以及空间分辨率，像素值具有强抑制噪声V′(X)＝(V′1(x)，V′2(x)，…，V′k(x))(＝方程式(3))，也就是，可以得到具有强抑制噪声的K静态或动态图像。
为了便于理解，在图3A-3C中示意性地表示这个图像处理。也就是，首先在图3A中，在1，2，…，k静态图像中，关于某一像素x，假设以像素x为中心的3×3像素的矩形区域N3×3(x)。假设在矩形区域N3×3(x)的左角上的像素是y1，像素y1具有像素值v(y1)。
而且，通过构成像素值v(y1)的标量值V1(y1)，V2(y2)，…，Vk(y1)和像素值v(x)中的标量值V1(x)，V2(x)，…，Vk(x)，根据上面方程式(4)(图3B)来计算权w1(p(x，y1))。而且，关于矩形区域3X3(x)的剩余像素y2，…y8，最后如图3B所示，可以类似地得到w1(p(x，y1))，…，w1(p(x，y8))和w1(p(x，x))。在这种情况下，从方程式(8)中得出，风险率p(x，y)是“1”，和因此从方程式(9)中得出，权w1(p(x，y)也是“1”(＝应用最大权)。
随后，用这种方式得到的权w1(p(x，y1))，…，w1(p(x，y8))和w1(p(x，x))与第k个图像中相应像素的标量值Vk(y1)，Vk(y2)，…，Vk(y8)，Vk(x)相乘得到总和(与上面方程式(2)的分子相对应)。该和除以关于矩形区域N3×3(x)的权w1的和(类似地对应于上面方程式(2)的分母)。然后，关于第k个图像中的像素x，可以得到具有抑制噪声的标量值V′k(x)(图3C)。而且，关于所有图像k＝1，2，…，k，使用相同的权w1(p(x，y1))，…，w1(p(x，y8))和w1(p(x，x))得到具有抑制噪声的标量值V′k(x)，和由此得到像素x中具有抑制噪声的像素值V′k(x)＝(V′1(x)，V′2(x)，…，V′k(x))。当关于所有像素x重复上面的计算时，得到具有抑制噪声的K图像。
在由用相干滤波计算的像素值V′(x)构成的图像中，在原始图像中可见的随机噪声得到有效抑制。
而且，可以执行上述步骤，例如，参考图4A、4B中所示的流程。而且，为了在实际的X射线CT扫描仪100上实现各个步骤的计算和图像显示，例如，如图2所示，可以布置包括分散值评价部件111、权计算部件112和像素值计算部件113的图像处理部件110来执行这些步骤。
在这种情况下，配置权计算部件112直接从如上所述的像素值v(x)和v(y)中得到权w1(p(x，y))。因此，计算部件112是一种直接得到权而不具体得到风险率p1(x，y)的装置。此外，代替上述构造，一种构造可以由风险率计算部件(保真度量化部件)来执行两步骤的程序，用于具体得到风险率p1(x，y)的值，和权计算部件用于根据风险率计算部件的输出得到权w1(p(x，y))。在任何情况下，权计算部件112使用散布σ2、v(x)、v(y)来计算权w1(p(x，y))，散布σ2由分散值评价部件111来估计。
而且，像素值计算部件113使用像素值v(x)和v(y)和权w1(p(x，y))来计算像素值v′(x)，权w1(p(x，y))的数字值由权计算部件112来计算。也就是，计算部件113实际上执行抑制原始图像噪声的步骤，也就是，应用相干滤波器(下文中称为“应用相干滤波器”)。
当相干滤波器被应用到动态图像中时，包括上述动态相干滤波处理中的K静态图像，图像处理部件110中的处理可以包括首先还原所有静态图像；随后把图像保存在存储装置10M中；和把相干滤波器应用到所述图像作为后处理。然而，本发明不局限于这种模式。在上述的连续扫描、连续投影数据收集、连续还原和连续显示的流程中，可以实时执行应用相干滤波器的步骤(下文中称为“实时相干滤波处理”)。
在实时相干滤波处理的一个优选实施例中，每当扫描和还原新图像时，执行下面处理。在第一所获图像(图像号为1)到最后的图像(图像号为M)中，公共相同点(相同坐标)的像素x的像素值(标量值)排列在具有图像号M，M-1，…，M-K+1的静态图像上，来构造K维向量值v(x)＝(Vm(x)Vm-1(x)，…，Vm-k+1(x))。照这样，能够以与上述“动态相干滤波处理”相同的方式应用相干滤波器。实际上，像素值计算部件113根据最近图像(图像号为M)只计算标量值v′m(x)，代替计算所有元素的像素值v′(x)。作为结果，由于提高了计算速度，可以实时显示具有抑制噪声的最近图像。
作为“实时相干滤波处理”的另一优选实施例，其构造可以包括以和上述相同的方式应用相干滤波器得到某一时间上的V1′(x)，…，Vk′(x)，在该时间首先得到K图像；随后通过v(x)＝(Vm(x)，vm-1(x)，…，vm-k+1(x))，使用具有图像号M，M-1，…，M-K+1的K静态图像来构造K维向量值；和关于该值应用上述实时相干滤波处理。另外，在实时相干滤波处理期间，像素值向量维数K可以通过手动或自动设定按照需要而改变，并且这种构造是方便的。
以这种方式，通过相干滤波器，使用CT图像来执行CBP学习、组织中局部血液流量循环，也就是，流经局部组织中毛细管的血液流量循环得到量化分析，CT图像中仅仅有效抑制了噪声而不恶化空间或时间分辨率，和得到表示局部血液流量循环的指数(CBP、CBV、MTT、Err)，由此可以期望精确度和可靠性的增加。
如上所述，关于所述图像执行CBP学习处理，在该图像中抑制了分辨率的降低和消除了噪声。
(CBP学习)如上所述，CBP学习包括获得CBP、CBV、MTT和Err指数，这些指数量化地表示了脑组织中“流经毛细管的血液流量”的循环；以及，输出表示这些指数空间分布的图CBP脑组织毛细管中每单位体积和时间的血液流速[ml/100ml/min]；CBV脑组织中每单位体积的血液量[ml/100ml]；MTT毛细管的血液平均通行时间[秒]；和Err实际测量值与分析模型的偏差剩余错误误差指数。
而且，当剩余错误误差Err低时，意味着依赖参考脑动脉的高的可能性。相反，当剩余错误误差Err高时，意味着依赖参考脑动脉的低的可能性。
在CBP学习中，采用不具有大脑血管渗透性的造影剂，例如碘化造影剂作为示踪剂。经由肘脉通过注射器快速注入的碘化造影剂经由心脏和肺流进脑动脉。而且，造影剂经由脑组织中的毛细管从脑动脉中流出到肘脉。在这种情况下，造影剂不具有大脑血管渗透性，例如碘化造影剂流经正常脑组织中的毛细管而没有任何血管外的泄漏。
通过动态CT，连续地扫描造影剂的流经状态，和从连续图像中测量脑动脉上像素的时间-密度曲线Ca(t)、包括毛细管的脑组织上像素的时间-密度曲线Ci(t)、和大脑静脉上像素的时间-密度曲线Csss(t)。
在CBP学习中，在脑组织附近的大脑血管血液浓度的时间曲线Ca(t)和毛细管血液浓度的时间曲线Ci(t)之间建立的关于造影剂血液浓度的理想关系用作分析模型。也就是，假定如果就在脑组织之前经由血管注射造影剂，用于包括毛细管的脑组织的单位体积(一个像素)里的时间-密度曲线垂直上升，停在一常数值上，和稍微倾斜下落。这可以用矩形函数来近似(框-调制转换函数方法box-MTF方法)。
脑动脉血液的时间-密度曲线Ca(t)被用作输入函数，用于脑组织的时间-密度曲线Ci(t)被用作输出函数，且可以得到作为调制转换函数的流经毛细管的过程特征，该函数用矩形函数来表示。
(具体过程)图5和6表示根据本实施例CBP学习的典型程序。首先，执行组织例如肘脉的大丸药注射(造影剂只给予一次)，和在注射之后或之前立即执行动态CT(重复扫描相同位置)。
作为一种最典型的方法，当执行肘脉的大丸药注射时，扫描被重复，例如，以0.5-2秒的间隔持续约20-40秒。通过动态CT得到的N个CT图像中的j个图像的每个像素(x，y)的CT号假定为(x，y，j)。这只是像素(x，y)中的时间-密度曲线(光滑曲线)f(t，x，y)的采样值。
首先，在步骤S1中处理各个CT图像时，外观上可明显判断出对应于脑组织以外的组织的像素被从分析目标中排除。也就是，表示推断为脑组织CT数(例如10-60HU的CT数)的范围之外的值的像素，是与空气、骨头或脂肪相应的像素，它们与固定量大脑血液流量无关，并可以忽略。该分析范围设定为10-60HU作为默认值，但可以经由输入部件109任意设定。
而且，由于在步骤S2中的处理，差异增强效果被初始化。为了在每个像素中(CT数上升)得到差异增加效果，在造影剂到达对应于每个像素(x，y)的组织之前的图像(通常获得多个图像)用序号1，2，…，k来表示和得到时间平均值如下。(Σj=1Kv(x,y,j))/K---(10)]]>该值设定为b(x，y)。而且，每个图像j＝k+1，k+2，…，N的像素值v(x，y，j)计算如下。
q(x，y，j)＝v(x，y，j)-b(x，y)当j＜K时，可以设定下式。
q(x，y，j)＝0为了简化处理，关于任何像素可以使用相同K。这样得到的q(x，y，j)可以认为只是时间-密度曲线q(t，x，y)的采样值，其作为t＝t1，t2…，tN中的光滑连续曲线。使用这个q(t，x，y)来执行大脑血液流量的量化分析。
在量化分析中，首先在CT图像上可以将右半球和左半球分开。在上述CBP学习中，得到毛细管血液流量循环的状态，作为脑组织的时间-密度曲线Ci(t)关于脑动脉的时间-密度曲线Ca(t)的调制转换函数MTF。因此，除非作为分析目标的脑组织依赖于参考曲线Ca(t)的脑动脉，计算是无用的。单个脑动脉时间-密度曲线Ca(t)至少用于左脑和右脑来单独分析目标。也就是，左脑的脑动脉时间-密度曲线Ca(t)只在分析同一左脑的脑组织中使用。类似地，右脑的脑动脉时间-密度曲线Ca(t)只在分析同一右脑的脑组织中使用。这有效地减少了无效计算。
为了把大脑分成左半球和右半球，如图7中所示，划分线被添加和显示在屏幕中的CT图像上(S3)。可以首先构造划分线来显示在图像中间。操作员参考图像，移动划分线，移动构成划分线的多个点，任意弯曲该线，和由此分成左和右区。
当大脑分成左和右半球且分析范围被限定在各自区域中时，可以减少分析处理步骤的数量。也就是，左脑的脑动脉的时间-密度曲线Ca(t)(左ACA、MCA、PCA)只用来分析左脑区(调制转换函数最优化处理)。类似地，右脑的脑动脉的时间-密度曲线Ca(t)(右ACA、MCA、PCA)只用来分析同一右脑的脑组织。为了减少分析处理步骤的数量，左和右半球还分成区域，和分析处理可以限定在更窄区域。
为了划分区域，可以使用几何结构方法例如Voronoy方法。众所周知，Voronoy方法是一种频繁用于工厂最佳布置领域中的技术，例如医院、商场和消防队，和其特征在于根据布置在平面上的大量点(对应于商店、母线)距离把平面分成多个影响区域。
如图8中所示，在本实施例中，把Voronoy方法单独应用到左和右半球。左ACA、MCA、PCA用作三根母线把左脑区分成左ACA、MCA、PCA影响区域。Voronoy点设定在圆的中心，该圆经过三根与左ACA、MCA、PCA相应的母线。以Voronoy点为中心，左ACA和MCA的两根母线的垂直平分线、左MCA和PCA的两根母线的垂直平分线、和左ACA和PCA的两根母线的垂直平分线相互连接在一起。通过三根垂直平分线把左脑区分成三个影响区域。类似地，右ACA、MCA、PCA用作三根母线把右脑区分成右ACA、MCA、PCA影响区域。
脑组织的时间-密度曲线Ci(t)对左ACA的时间-密度曲线Ca(t)的调制转换函数MTF限定到左ACA的影响区域，并得到用于每个像素。类似地，调制转换函数MTF限定到关于左MCA和PCA、以及右ACA、MCA和PCA的影响区域，并得到用于每个像素。
当把左和右半球分成多个影响区域和分析范围限定到每个影响区域时，可以进一步减少分析处理步骤的数量。
随后，脑动脉ROI设定在CT图像上的脑动脉上。为了增加设定的精确性和帮助所述设定，ROI设定支持部件121准备一张支持图，和该图从或在CT图像(S4)中单独显示。如图9A-9C所示，支持图的实例包括出现时间图(AT)、最高时间图(PT)和通行时间图(TT)。关于各个像素，如图10所示，时间AT从反差成像之前的任意时间TO(例如数据收集启动时间)直到造影剂浓度达到最高的几个百分点(例如1％)，时间(最高时间)PT从时间TO直到造影剂浓度达到最高，或TT表示造影剂的移动时间，例如具有半值宽度，计算它们并产生和显示为图。缺省地，包括这些AT、PT和TT的所有类型得到产生和显示，但是操作员可以任意选择一种类型或两种类型。
与脑动脉中其它组织相比，这些数值趋于出现大值。因此，由彩色显示器经过彩色查找表，设定查找表只显示具有在该值中心的值的像素，很容易识别脑动脉位置，和可以准确设定脑动脉ROI(S5)。典型地，在每个左和右半球中，脑动脉ROI设定为前脑动脉(ACA)、中间脑动脉(MCA)和后脑动脉(PCA)三个位置。
而且，由于多切片扫描，例如，当相互相邻布置的四个切片是分析目标时，如图11所示，用于各个切片的各个脑动脉ROT的设定不仅在操作负担中是巨大的，而且在执行分析中是不必要的操作。因此，在某一任意切片中设定的脑动脉ROI也用于其它切片。而且，后面描述的相干递归方法也可以用于准备脑动脉的时间-密度曲线Ca(t)，其可以共有地用于所有切片。
随后，时间-密度曲线准备部件122从连续图像数据中通过关于每组脑动脉ROI的动态CT来准备时间-密度曲线Ca(t)(S6)。
这里，存在许多脑动脉，与像素大小进行比较，它们非常薄。另外，通常，脑动脉与CT切片不垂直相交。因此，图像上的每个像素不能准确地表示动脉血液的CT数，由脑动脉和其它组织的混合存在构成一个像素，和因为大部分情况下的部分体积影响，因此只表示较低的对比增强效果。而且，在这些包括任意部分体积的像素中，图像噪声是巨大的。特别是在部分产生大脑梗塞的动脉中由于对比增强效果是显著地小，噪声的影响是非常巨大。通过上述的相干滤波器来抑制图像噪声，但是还保持部分体积效果的影响。
对于脑动脉的时间-密度曲线来说，代替在单个切片图像中测量曲线，使用包括动脉的立体中的像素来应用后面描述的相干递归方法，由此可以解决问题。因此，代替上述相干滤波方法，在该阶段可以应用相干递归方法。
而且，根据这种方法，得到与每个动脉相应的仅仅一个切片图像的时间-密度曲线，和因此可以用于分析扫描容量里的所有切片中的任意部分。由此，关于特殊动脉，选择切片，通过该切片非常清楚地得到脑动脉时间-密度曲线，脑动脉时间-密度曲线可以应用到所有切片，和可以减少脑动脉时间-密度曲线的数量。
(相干递归方法)当准备时间-密度曲线时，重要的是消除部分体积效果和任意噪声的影响。首先，“时间-密度曲线”是一种表示动态CT图像特殊部分中CT数(密度值)随时间流失的变化的曲线。特别地，在上述的医学图像诊断装置中，为了检查局部血液流量循环和人体组织中新陈代谢功能的细节，人体特殊组织中的造影剂浓度随时间的变化已经测量为时间-密度曲线。而且，在天文观察中，为了分析特殊天体的发光变化，使用时间-密度曲线。将更正式和清楚地描述该曲线。也就是，时间-密度曲线表示为一对序列{<tk，dk>(k＝1，2，…，K)}，这里某一部分在时间tk中的密度值是dk。而且，在时间-密度曲线的许多使用中，不需要要求dk的绝对值，和只得到增量(dk-d1)是相当足够的，在增量中第一图像1用作参考。进一步说，在许多使用中，只得到数据A(dk-d1)是充分的(这里A是一个未知常数)，该数据与(dk-d1)成简单比例。因此，在这种情况下，这对序列{<tk，A(dk-d1)>(k＝1，2，…，K)}是所需的时间-密度曲线。
为了得到时间-密度曲线，原理上，包括在部分中的像素x的标量值vk(x)可以用来构成这对序列{<tk，vk(x)>}或{<tk，A(vk(x)-v1(x))>}，在该部分中时间-密度曲线在构成动态CT图像的各个图像(k＝1，2，…，K)中进行测量。
然而，在实际使用中，由医学图像诊断装置扫描的动态CT图像包括随机噪声，和因此存在的问题是实质上测量的时间-密度曲线不能准确得到。
此外，在实际使用中，在这些动态CT图像中，产生一种所谓的“部分体积效果”。部分体积效果是一种现象，其中物体中的微量材料的图像由图像上少量的像素来表示，但是少量像素也受到物体中邻近材料的图像的影响，和少量像素的像素值(本质上与要测量的密度值的波动成比例)只表示为一个相对小的波动。换句话说，少量像素的像素值只包括一个微小信号。因此，当产生部分体积效果时，对于任何像素x，这对序列{<tk，vk(x)>(k＝1，2，…，K)}具有非常低的信号电平，像素受到本质上没有测量的组织中密度值变化的影响，和还存在随机噪声。因此存在的问题是实质要测量的时间-密度曲线{<tk，dk>}不能准确得到。
为了解决该问题，迄今已经使用时间或空间滤波来抑制随机噪声。然而，由于时间平均，削弱了时间分辨率。而且，由于空间平均，存在的问题是实质要测量部分以外的部分密度随时间的变化混合到测量值中。为了解决该问题和得到更准确的时间-密度曲线，使用相干滤波器。
首先，将描述在本实施例的相干滤波器中使用的虚假设。假设要测量部分中的实际时间-密度曲线是{<tk，dk>(k＝1，2，…，K)}，为了测量主要转换{<tk，A(dk-d1)>(k＝1，2，…，K)}(这里A是一个未知系数)，设定一组R实质与要测量部分相应的像素，关于作为组R元素的任意像素xεR，在条件Q下即“像素x准确反映上述实际时间-密度曲线和几乎不受其它部分密度随时间变化的影响”，关于像素值v(x)＝(v1(x)，v2(x)，…，vk(x))作为向量值，考虑到部分体积效果和任意噪声的影响，可以假设建立下式。
Vk(x)＝p(x)dk+q(x)γk(x) …(11)
(k＝1，2，…，K)这里，p(x)和q(x)是未知系数，其随每个像素x不同，但是不随图像数k而变化(也就是扫描时间tk)，和模拟部分体积效果。而且，γk(x)模拟随机噪声，该值随每个像素x和图像号k而不同，但是期望值是0，和统计分布不依赖于像素x或图像号k。
根据上述假设，可以证明下式。如果关于两个任意像素x，y建立命题“两个像素x，y满足上述条件Q”作为组R的元素，建立下列关系。
Vk(x)＝a1vk(y)+a2+ξk(k＝1，2，…，K)…(12)这里，a1和a2是未知系数，其随每个像素组x，y而不同但是不随图像数k而变化(也就是扫描时间tk)。而且，ξk是随机噪声，该值随每个像素x，y和图像号k而不同，但是期望值是0。
得到方程式(12)如下。
Vk(y)＝p(y)dk+q(y)γk(y) …(13)也就是，当对通过给x赋值y得到的上述方程式进行修改时，得到下列方程式。dk=vk(y)-q(y)-γk(y)p(y)---(14)]]>当其赋值给方程式(20)时，得到下式。vk(x)=p(x)p(y)vk(y)+(q(x)-p(x)p(y)q(y))+(γk(x)-p(x)p(y)γk(y))---(15)]]>因此，使用下式，得到方程式(12)。a1=p(x)p(y)]]>a2=q(x)-p(x)p(y)---(16)]]>ξk=γk(x)-p(x)p(y)γk(y)]]>这里方程式(16)的a1和a2是表示部分体积效果的参数，和方程式(16)的ξk表示随机噪声。
上面已经表示，命题“像素x，y都满足条件Q”相当于虚假设H0′“Vk(x)＝a1vk(y)+a2+ξk(k＝1，2，…，K)”。
接下来将描述把虚假设H0 ′“Vk(x)＝a1vk(y)+a2+ξk(k＝1，2，…，K)”转换为一种实质等价命题形式，其实际上可以得到检验。为了重新严格地数学表示虚假设，虚假设H0′表示“由于存在某些常数a1和a2，ξk＝Vk(x)-a1vk(y)-a2(k＝1，2，…，K)遵循正常分布，其具有0和散布σ2h(a1)的中间值”。这里，系数h(a1)如下。
h(a1)＝1+a12…(17)从表示a1和ξk定义的方程式(16)和关于自由变量的散布通常特性中迅速得到方程式(17)。而且，可以简单地估计上述散布σ2的值充分用于实际使用中。
如上所述，如果可以确定常数a1和a2，检验上述虚假设H0′是可能的。而且，实际上，得到这些常数的最佳估计值a1-和a2-是足够的。
同样可以在计算常数a1和a2的最佳估计值中使用公知的拟合方法。然后，作为拟合方法的典型具体方法，将描述一种在线性化最小平方方法中使用的要点。为了把线性化最小平方方法应用到本实施例，假设虚假设ξk的平方和简化为s(a)，定义如下。S(a)ΣK=1K{vk(x)-a1vk(y)-a2}2---(18)]]>s(a)的值依赖于常数向量a＝(a1，a2)，也就是，常数a1和a2的值。使用常数向量a的计算，其中s(a)表示最小值，在无偏估计意义中，可以得到关于常数a1和a2的最优估计值a1-和a2-。而且，作为线性化最小平方方法的具体计算方法，可以使用各种已知的方法。另外，这些公知计算方法都是非常简单，和需要的计算时间非常少。常数a1和a2的最优估计值a1-和a2-用这种方式来计算。作为结果，由下式定义的剩余错误误差可以得到具体计算。
rk-(x，y)＝vk(x)-a1-vk(y)-a2-…(19)因此，可以使用该剩余错误差rk-来改写上述虚假设H0′作为实质等价的虚假设H0″“rk-(x，y)(k＝1，2，…，K)遵循具有0和散布(1+(a1-)2)σ2的中间值的正常分布”。这就是具体的命题，通过该命题可以实际执行检验计算。
而且，引入向量的下列表示a-＝(a1-，a2-)ζ＝(r1～，…，rk～)…(20)其中向量a和ζ依赖于像素组x，y。
f(a～，v(y))＝a1～，v(y)+a2-…(21)而且，由上式定义的向量值函数f用来改写虚假设H0′作为虚假设H0″“v(x)＝f(a～，v(y))+ζ(然而，ζ遵循具有0和散布(1+(a1-)2)σ2的中间值的正常分布”，和这是与上述虚假设HO的形式相同。也就是，显而易见本实施例是相干滤波器的一个改进实例。这里，f(a～，v(y))意味着使用像素x的像素值v(x)，表示部分体积效果的参数a可以调整到最优化和转换像素y的像素值v(y)以便具有最高保真度。
接下来将描述通过本实施例中的相干滤波器使用虚假设H0″得到时间-密度曲线的方法。对于通常与要测量的部分对应的组R来说，关于包括在组R中的某一像素x∈R，关于作为组R的元素的所有像素y∈R执行下列计算。也就是，使用上述方法来实际计算剩余错误误差rk-(x，y)(k＝1，2，…，K)。随后，具体计算风险率p(x，y)和权w(p(x，y))用于否决上面的虚假设H0″“rk-(x，y)(k＝1，2，…，K)遵循具有0和散布(1+(a1～)2)σ2的中间值的正常分布”。而且由下列方程式(22)和像素x中的时间-密度曲线{<tk，v′k(x)-v′1(x)>(k＝1，2，…，K)}来计算加权平均v′k(x)。v′(x)=Σy∈Rw(p(x,y)f(a~,-v(y))Σy∈Rw(p(x,y))---(22)]]>用这种方式得到的时间-密度曲线表示一个测量值，其近似于像素x{<tk，dk>}中的实际时间-密度曲线的初始转换{tk，A(dk-d1)>}(其中A是未知系数)。另外，通过加权平均的效果来抑制任意噪声。而且，关于另一像素y的像素值向量，正如从方程式中显而易见，使用所述曲线，在该曲线中部分体积效果的影响得到修正。而且，本发明具有相干滤波器的公共特征“不使用时间平均值，和根据具有像素x的保真度使用权来计算空间平均值”作为特性。因此，根据本实施例，可以得到时间-密度曲线，其中部分体积效果的影响得到抑制而不削弱时间分辨率和随机噪声得到抑制。而且，用这种方式得到时间-密度曲线的方法专门称为“相干递归方法”。
接下来将具体描述在通过医疗X射线CT中的动态CT扫描得到的动态CT图像中，时间-密度曲线的临床使用的一个实例。在该应用实例中，当造影剂快速注射进组织中时，执行扫描例如动态CT，测量存在于人体组织中的动脉图像的密度变化作为时间-密度曲线，和组织中的血液流量循环得到诊断。
在该应用实例中，在许多情况下，由于人体组织中的动脉通常非常薄，出现在CT摄影图像上的动脉图像产生部分体积效果。而且，不用说，该图像包含随机噪声。因此，在传统方法中得到关于动脉的充分准确的时间-密度曲线是困难的。当强迫执行测量时，仅仅得到测量值<tk，vk(x)-v1(x)>。该值近似于一定程度的像素值<tk，A(Dk-D1)>(这里Dk表示一组像素的像素值(作为标量值)，该组像素与时间tk k＝1，2，…，K中动脉图像相对应)作为关于动脉的实际时间-密度曲线<tk，Dk>的主转换。这个测量值包含随机噪声。而且，因为部分体积效果的影响，系数A保持未知。
然后，可能得到测量值<tk，(v′k(x)-v′1(x))>(k＝1，2，…，K)。其充分近似<tk，A(Dk-D1)>。另一方面，显著厚的血管存在于这些血管中，在相同的摄影图像上可以看到这些血管。因此，关于血管，在传统的方法中，可以得到时间-密度曲线的足够好的近似值<tk，(Jk-J1)>(k＝1，2，…，K)。这里Jk表示一组像素的像素值，其与时间tk中的血管图像相对应。
而且，在关于血液循环的时间-密度曲线中，公知的是建立命题S“如果时间t1中血液中的造影剂浓度、关于任何血管d的任何时间-密度曲线<tk，(dk-d1)>在曲线(AUC)下的区域里相匹配”作为特性。这里，曲线下的区域意味着时间-密度曲线<tk，(dk-d1)>的关于时间t的积分。
因此，可以近似计算关于某一组织的时间-密度曲线<tk，(dk-d1)>的曲线(AUC)下的面积，例如，通过下列方程式。AUC(d)≈ΣK=1K-1{(dk+1-d1)+(dk-d1)}(tk+1-tk)2---(23)]]>因此，可以使用方程式(22)计算关于时间-密度曲线{tk，(Jk-J1)>}的曲线AUC(J)下的面积，在传统方法中关于血管得到时间-密度曲线{tk，(Jk-J1)>}。J可以赋值给d。而且，如果关于动脉的时间-密度曲线{tk，(Dk-D1)>}是公知的，使用方程式(18)类似地计算曲线AUC(D)下的面积。另外，根据命题S，必须建立下式。
AUC(D)≌AUC(J) …(24)然而，实际上，由于时间-密度曲线{tk，(Dk-D1)>}是未知的，不能计算AUC(D)。
另一方面，在根据本实施例的方法中得到的时间-密度曲线<tk，(v′k(x)-v′1(x))>近似于<t k，A(Dk-D1)>，和后者包含未知系数A。因此，AUC(V′)下的面积必须是A与面积AUC(D)的乘积，使用方程式(23)从{<tk，(v′k(x)-v′1(x)>)}中可以具体计算AUC(V′)下的面积。也就是，下面结果。
AUC(V′)≌AUC(D)…(25)也就是，从方程式(24)和(25)中，建立下列关系。
A≌AUC(V′)/AUC(J)…(26)使用方程式(23)可以具体计算方程式(26)的右侧，和因此可以具体确定未知系数A的值然后，当系数A的这个值用来构成时间-密度曲线<tk，(v′k(x)-v′1(x)>)/A>，这只是动脉的时间-密度曲线<tk，(Dk-D1)>的近似。使用曲线下的面积构成时间-密度曲线的方法称作为“AUC方法”，该时间-密度曲线具有未知比例系数A的确定值。
因此，AUC方法还与相干递归方法相结合，其用于由动态CT扫描得到的动态CT图像中的时间-密度曲线的临床使用中，由此，消除了部分体积效果和随机噪声的影响，和甚至关于薄动脉的时间-密度曲线得到包括未知比例系数A的测量值。在传统方法中曲线测量已经是困难的或不可能。
而且，当然，AUC方法也可以应用到关于通过传统方法单独测量的动脉的时间-密度曲线<tk，(v′k(x)-v′1(x)>。不能消除部分体积效果和随机噪声的影响，但是可以构成时间-密度曲线，其具有未知比例系数A的确定值。
(使用上纵向窦性静脉的时间密度曲线校正脑动脉的时间密度曲线(抑制部分体积效果的影响))为了抑制部分体积的影响，代替或排除相干递归，可以使用上纵向窦性静脉的时间密度曲线Csss(t)校正脑动脉的时间密度曲线Ca(t)。
首先，在图5的步骤S7中，如图12所示，设定稍微大的上纵向窦性静脉ROI以便包围CT图像上的上纵向窦性静脉。与所述动脉相比，上纵向窦性静脉是大的和其位置相对固定，容易设定上纵向窦性静脉ROI。该稍微大的上纵向窦性静脉ROI包括多个像素。
随后，减少或处理所述上纵向窦性静脉ROI以便上纵向窦性静脉ROI的所有像素包括在整个区域上的上纵向窦性静脉中(S8)。例如，减少处理包括首先关于上纵向窦性静脉ROI中的每个像素执行极限值处理(二进制化)；和在ROI中准备二进制图(“0”，“1”)。所述极限值设定为一个值，该值把上纵向窦性静脉的图像从外围组织和骨头的图像中分离。“1”表示上纵向窦性静脉图像上的像素，和“0”表示外围组织和血液图像上的像素。根据附近中的四个或八个像素的值来更换二进制图的每个像素(中心像素)。只有当中心像素是“1”和所有附近中的四个或八个像素的值是“1”时，中心像素值保持在“1”。也就是，当然，当中心像素是“0”时，甚至当中心像素是“1”时，和所有附近中的四个或八个像素中任何一个也表示“0” 时，中心像素值用“0”取代。因此，与上纵向窦性静脉图像的外形相比，通过至少一个像素来减少上纵向窦性静脉ROI。由此，受到减少处理的上纵向窦性静脉ROI中所有像素是上纵向窦性静脉图像上的像素的条件，可以实现并具有很高程度确定性。
而且，代替这种技术，可以使用时间-密度曲线的曲线AUC下的面积来修正所述上纵向窦性静脉ROI。在这种情况下，使用稍微大的ROI作为搜索范围，和在该范围内关于每个像素计算时间-密度曲线的曲线AUC下的面积。通过对比增加效果，与外围像素值相比，上纵向窦性静脉图像上的像素的曲线AUC下的面积显然表示一个高值。因此，当关于曲线AUC下的面积执行极限值处理时，从ROI中能够选择的只有上纵向窦性静脉图像上的像素。
关于在很大程度肯定在上纵向窦性静脉图像上的多个像素，在这种方式中使用任何一种或两种技术通过AND条件获得所述图像，平均每个像素的时间-密度曲线，和准备上纵向窦性静脉的时间-密度曲线Csss(t)(S9)。这里，由于碘化造影剂不流经血脑障碍，原理上，碘浓度随脑动脉和静脉不改变。也就是，上纵向窦性静脉的时间-密度曲线Csss(t)的曲线AUC下的面积几乎等于S6中准备的脑动脉的时间-密度曲线Ca(t)的曲线AUC下的面积。因此，如图13所示，通过把S6中准备的脑动脉的时间-密度曲线Ca(t)的每个时间值与AUC(sss/AUCa)相乘来修正时间-密度曲线Ca(t)，由此S6中准备的脑动脉的时间-密度曲线Ca(t)的曲线AUCa下的面积几乎等于上纵向窦性静脉的时间-密度曲线Csss(t)的曲线AUCsss下的面积(S10)。
随后，图14A中所示的脑动脉的时间-密度曲线Ca(t)用来量化脑组织(毛细管)的血液流量循环，其中如上所述噪声和部分体积效果得到抑制。为此，首先关于脑组织上的每个像素准备图14B中表示的时间-密度曲线Ci(t)(S11)。
随后，如S12所示，单个脑动脉时间-密度曲线Ca(t)用于左区和右区，和脑动脉时间-密度曲线Ca(t)用作每个像素的输入函数和脑组织的时间-密度曲线Ci(t)用作每个像素的输出函数，得到示踪剂经过毛细管的经过过程特征，作为调制转换函数MTF。也就是，关于左区脑组织的时间-密度曲线Ci(t)，使用相同区的脑动脉时间-密度曲线Ca(t)。而且，关于右区脑组织的时间-密度曲线Ci(t)，使用同一右区的脑动脉时间-密度曲线Ca(t)来得到调制转换函数MTF。而且，由于如上所述为每个ACA、MCA、PCA准备脑动脉时间-密度曲线Ca(t)，每个Ca(t)重复计算调制转换函数MTF。
这里，如图15所示，在脑动脉时间-密度曲线Ca(t)和毛细管的时间-密度曲线Ci(t)之间建立的理想关系用作一种分析模型，其应用box-MTF方法。
图16表示box-MTF方法的原理。该方法包括评估脑动脉时间-密度曲线Ca(t)与用矩形函数表示的调制转换函数box-MTF的卷积和S11中准备的实际测量Ci(t)之间的剩余错误误差；和校正调制转换函数box-MTF以便减少剩余错误误差的平方和。重复该步骤程序来最小化剩余错误误差。
如图14所示，根据调制转换函数box-MTF来计算CBP、CBV、MTT，该函数最小化了剩余错误误差(S13)，和在S12中最小化的剩余错误误差平方和输出为Err。严格地讲，用下式执行修正。
CBP＝CBPCBV＝(1-Ht)/(1-b*Ht)*CBV′MTT＝(1-Ht)/(1-b*Ht)*MTT′这里，Ht是主要血管的血球容积计值，和b*Ht是外围血管的血球容积计值(通常b约为0.7)。
而且，由yi(x)-f(ti，x)给出剩余错误误差。Yi(x)表示时间ti中的框细胞x的标量值，和与脑组织的时间-密度曲线一致。f(ti，x)表示一模型的时间ti中的标量值，该模型适合于框细胞x的向量像素值，和符合调制转换函数与脑动脉时间-密度曲线的卷积。在近似调制转换函数中Err是剩余错误误差平方和的平方根，和计算，例如用下列方程式表示。Err=Σi=1N(yi(x)-f(ti,x))2(w(ti))2S]]>这里，S是常数，例如S＝N-p。p表示自由度，也就是包括在近似模型f中参数的数量。w(ti)是权系数，其确定时间ti-Err中剩余错误误差的贡献程度。例如，w(ti)不必要依赖于i，和可以具有固定值例如1。作为选择，w(ti)＝αe-ti2/β是合意的，构造其使得权w随着|ti|的增加适度地减少。随着时间的流失，剩余错误误差不必要作用到Err。
关于从调制转换函数中计算的CBP、CBV、MTT、Err，在这种方式中通过box-MTF方法得到调制转换函数，如图17中所示，设定各个输出范围(适当范围)(S14)。关于具有相应输出范围值的像素，保持该值。对于具有输出范围外值的像素来说，替换该值，例如，用与显示器上黑色相应的值(S15)。
随后，如图1SA-18D所示，从经过输出最优化的CBP、CBV、MTT、Err中产生各自图(S16)。关于前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA，单独计算用于各个切片的CBP、CBV、MTT、Err指数。因此，如图19所示，甚至使用一个切片，得到4×3＝12个图。当设定的脑动脉数量增加时，图增加到动脉增加数量的4倍数量。通常估计这样许多图是不现实的。然后，为了减少图的数量，合成这些图(S17)。
如图20所示，合成方法包括根据前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的剩余错误误差Err来合成前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的CBP图。例如，当从前脑动脉ACA的时间-密度曲线Ca(t)和用于控制下的脑组织的时间-密度曲线Ci(t)中得到调制转换函数MTF时，剩余错误误差Err相对小。相反，当从用于非控制下的脑组织的时间-密度曲线Ci(t)中得到调制转换函数MTF时，剩余错误误差Err相对大。也就是，剩余错误误差Err表示每根脑动脉的控制可能性。
因此，对于每个像素来说，从前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的CBP值中，选择与最低剩余错误误差Err值相应的CBP值作为像素值。用这种方式合成的图由脑组织的CBP值构成，其具有在前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的控制下的高可能性。这也可以应用到其它指数CBV、MTT的图合成。
这里，下文中将详细描述图合成。从存在于动脉相应位置的像素中得到的时间-密度曲线反映动脉血液中造影剂浓度，和应用上述相干递归方法，由此可以得到修正的动脉血液造影剂浓度的时间-密度曲线。可以准备这种脑动脉时间-密度曲线用于每根动脉，和依赖血液循环状态而变化。特别地，在大脑血管失调的情况下，有时差别是显著的。例如，在K个位置动脉中得到的脑动脉时间-密度曲线表示为Ak(t)(k＝1，2，…，K)。
当某一组织的时间-密度曲线与滋养(控制)组织的动脉的脑动脉的时间-密度曲线相比时，可能得到指数例如CBP，其反映组织中的微循环(毛细管系统的结构、功能)。由于计算这些指数用于每个部分(x，y，z)，可以用作为像素值的值来构造该图，和该图是指数图。例如，当得到R类型指数(如上所述典型的是四种类型CBP、CBV、MTT、Err)时，可以构造R图。用这种方式准备的R图可以看作为一个图(向量值图)，其中每个像素具有一向量值。也就是，下列结果。
Vk(x，y，z)＝<Pk，1(x，y，z)，Pk，2(x，y，z)，…，Pk，R(x，y，z)，>例如，在CBP学习中，如上所述典型地设定为R＝4，和可以构造Pk，1(x，y，z)来表示CBP的值、构造Pk，2(x，y，z)来表示CBV的值、构造Pk，3(x，y，z)来表示MTT的值、和构造Pk，4(x，y，z)来表示剩余错误误差Err的值。
在部分(x，y，z)中，从开始，作为分析目标，将显然与内部器官不相应的部分从分析目标中排除，和表示从分析目标中排除的特殊值可以分配给Pk，r(x，y)(步骤S14、S15)。使用绝对值大的负值作为该值是方便的。作为选择，作为要增加到向量Vk(x，y，z)的另一元素，也可以准备图Pk，R+1(x，y，z)＝(0当(x，y，z)从分析目标中排除时，否则1)。该图称为“掩模”。
为涉及的脑动脉的每个时间-密度曲线Ak准备这个向量值图Vk。例如，如果从左和右中间、前、后脑动脉中得到的脑动脉时间-密度曲线时，K＝6。而且，当从影响区域的外围中的几根动脉中得到的脑动脉时间-密度曲线时，K＝约10-15。
因此，从右半球动脉中得到的脑动脉时间-密度曲线必须只用在分析属于右半球的部分(x，y，z)中，然而从左半球动脉中得到的脑动脉时间-密度曲线必须只用在分析属于左半球的部分(x，y，z)中。然后，操作员指定右和左半球的边界(中线)为直线、曲线、弯线、平面或曲面，由此在构造中更好地准备图用于每个半球。然而，脑动脉时间-密度曲线的数量是K＝约3-10，脑动脉还存在于每个半球中。
当脑动脉时间-密度曲线Ak(k＝1，2，…K)的数量K大时，作为结果得到的向量值图Vk(k＝1，2，…K)的数量大，并且这不便于观察。也就是，对于正常的灰度级图像或彩色级图像的观察来说，一个图由R个图形构成，存在K个图，和因此总共K×R个图像必须相互比较。而且，由动脉滋养的部分和相应的动脉是不必要显出的，和必须使用解剖知识来判断每部分要观察的图Vk(k＝1，2，…K)。尤其是在大脑血管失调发展的情况下，例如大脑梗塞，控制组织的动脉判断与解剖知识不一致，和异常控制频繁可见。依据这些问题，存在的问题是难以解释向量值图的射线照相。
为解决所述问题，合成这些图。也就是，剩余错误误差图用来把K个向量值图Vk(k＝1，2，…K)放在一起组成一个向量值图V。例如，当Pk，R(x，y，z)是剩余错误误差图时，V(x，y，z)＝Vk(x，y，z)，其中K是在k＝1，2，…k中的使Pk，R(x，y，z)|最小化的值。
而且，可以增加图P0(x，y，z)＝(K，在k＝1，2，…k中使得|Pk，R(x，y，z)|最小化)来表示在k＝1，2，…k中每个部分中使用的K。
根据这种方法，也就是，在观察普通的灰度级图像或彩色级图像期间，可以观察到R或R+1个图像。
根据该方法，存在可能性就是错误使用Ak的计算结果用在要实质计算的部分(x，y，z)中，其使用脑动脉时间-密度曲线Aj。然而，要犯这种错误，从V(x，y，z)的定义中看出，关系|Pk，R(x，y，z)|<|Pj，R(x，y，z)|是必要的。只有当Aj显著类似于Ak时，才会建立这种关系。因此，认为在部分(x，y，z)中Vk(x，y，z)最初类似于Vj(x，y，z)，和在解释Vk(x，y，z)中犯这种错误几乎没有可能。
当实际应用这种方法时，和只有当Aj显著类似于Ak时，在表面上几乎均匀的组织中，按照每个部分，P0(x，y，z)＝K或P0(x，y，z)＝j。在这种情况下，Pk，r(x，y，z)≡Pj，r(x，y，z)(r＝1，2，…R)。可以看到甚至使用任一值几乎没有差别。
相反，关于受到特殊动脉控制的组织，当相应脑动脉的时间-密度曲线Ak不类似于另一曲线时，通过使用本方法，在所述组织的部分(x，y，z)中几乎安全地和自动地选择Vk(x，y，z)。因此，当看到P0(x，y，z)时，可以观察到组织中的控制和控制动脉而无需任何解剖知识。
这里，返回到图6。如图21A-21D中所示，关于如上所述的合成图，或每个脑动脉中的单个CBP、CBV、MTT或Err图，设定包括许多像素的影响ROI区域(S18)，计算ROI中的像素值(CBP、CBV、MTT、Err值)的平均值(CBP、CBV、MTT、Err平均值)(S19)，和有时该平均值用作诊断材料。在这个平均中，设定适当范围用于步骤S14中的CBP、CBV、MTT、Err，保持该范围内的值，和该范围外的值用最小值来取代，例如与黑色表示相应的值。因此，平均包括替代值的值，该平均值自然包括一个误差。因此，对于平均处理来说，需要选择只在适当范围内的值，或排除替代值，和执行平均处理。
而且，为了设定影响ROI的区域用于平均，当在CBP、CBV、MTT、Err图的任何图上设定影响ROI的区域时，ROI可以共有地用于其它图。简化ROI设定操作，和也可能计算关于相同ROI(CBP、CBV、MTT、Err平均值)的平均值。
这里，如上所述，某一组织的时间-密度曲线关于某一脑动脉的时间-密度曲线的最小剩余错误误差Err，表示组织中受脑动脉控制的程度，也就是，一定程度的血液流量通过脑动脉供应给组织。换句话说，所述误差表示属于脑动脉组织中的程度，也就是，从脑动脉供应给脑组织的一定程度的血液流量。与小剩余错误误差Err相应的脑动脉表示关于所述像素的脑组织的高控制可能性，和与大剩余错误误差Err相应的脑动脉表示关于所述像素的脑组织的低控制可能性。因此，可能产生一个图，其使用每个像素的标签从剩余错误误差Err中识别具有高控制可能性的脑动脉，也就是，一个控制图，其识别前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的具有高控制可能性的区域。
如图22所示，关于脑子左和右区，前脑动脉ACA、中间脑动脉MCA和后脑动脉PCA的剩余错误差Err相互比较用于每个像素。指出表示最小剩余错误误差Err的脑动脉(ACA、MCA或PCA)具有控制所述像素的脑组织的最高可能性。对于每个像素来说，指定具有最高控制可能性的脑动脉，也就是剩余错误误差Err的最小值。相应于指定的脑动脉的标签给每个像素。
图23表示产生的控制图实例。在这个控制图中，用颜色和底纹来区分和显示标签。而且，当用任意标签过滤指数图时，如图24A、24B、24C中所示，可能从指数图中提取具有高控制可能性的区域用于每根脑动脉(ACA、MCA、或PCA)。
如上所述，根据本实施例，通过相干滤波或相干递归，抑制了空间和时间分辨率的降低，和抑制了噪声，由此可以增加CBP学习的分析精度。
对于本领域技术人员来说，其它优点和改良将容易发生。因此，本发明的方面不局限于本文中所表示和描述的详细细节和典型实施例。因此，可以进行各种修改而不超出由附属的权利要求和它们的等价物所限定的总的发明精神或范围。
权利要求
1.一种指数计算方法，包括从有关注射有造影剂的受试者的一个特定区域的多个连续图像，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；通过曲线拟合来计算表示组织中关于动脉的局部血液流循环的调制转换函数，从而使第二时间-密度曲线的剩余误差相对于所述调制转换函数与第一时间-密度曲线的卷积得到最小化；从所述调制转换函数计算相对于各个动脉的局部血液流量循环的指数；制备所述动脉的指数图；和根据第一时间-密度曲线的所述剩余误差，把所述指数图合成为一个图。
2.根据权利要求1所述的指数计算方法，其中特定区域是脑部。
3.根据权利要求2所述的指数计算方法，还包含计算脑组织的每单位体积和时间的血液流速、脑组织中每单位体积的血液量和脑组织的血液流量平均通行时间中至少一个作为所述指数。
4.根据权利要求2所述的指数计算方法，还包含分别计算脑部左和右区的调制转换函数。
5.根据权利要求2所述的指数计算方法，还包含分别计算动脉影响区域的调制转换函数，所述区域通过用几何技术划分脑部区域而得到。
6.根据权利要求5所述的指数计算方法，其中几何技术是Voronoy技术。
7.根据权利要求2所述的指数计算方法，还包含根据从所述图像中准备的关于静脉的时间-密度曲线，修正关于所述动脉的第一时间-密度曲线。
8.根据权利要求7所述的指数计算方法，还包含修正第一时间-密度曲线，以使所述第一时间-密度曲线下的面积等于静脉的时间-密度曲线下的面积。
9.根据权利要求7所述的指数计算方法，还包含从关于静脉区域中的多个像素的平均密度中，准备上纵向窦性静脉的时间-密度曲线。
10.根据权利要求1所述的指数计算方法，还包含计算构成所述图像的多个像素的每个时间-密度曲线的特征值；和为了支持动脉设定显示特征值图。
11.根据权利要求10所述的指数计算方法，其中所述特征值是出现时间、最高峰时间或通行时间。
12.根据权利要求1所述的指数计算方法，还包含当所述指数值超出预定范围时，用预定值来替换所述指数值，该预定值与近似黑色的颜色或射线照相解释不使用的颜色相对应。
13.一种方法包含从有关注射有造影剂的受试者的一个特定区域的许多连续图像，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；选择所述第一时间-密度曲线中与相应的第二时间-密度曲线最拟合的一条第一时间密度曲线，以指定一根动脉，该动脉是所述组织中的每根动脉的局部血液流量循环所最可能依赖的；根据为每一组织指定的一条动脉，制备用于区分该动脉的一个依赖区的一个图。
14.根据权利要求13所述的方法，还包含用颜色区分和显示所述图的动脉非独立区。
15.根据权利要求13所述的方法，还包含根据所述图合成关于指数的指数图，所述指数表示相应动脉的局部血液流量循环。
16.一种指数计算方法，包含利用根据像素间的类似性的一个权来对多个图像进行滤波，以减少来自涉及注射有造影剂的受试者的一个特定区域的该多个连续图像的噪声；从带有被减少的噪声的该多个连续图像中，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；通过曲线拟合来计算表示组织中关于各个动脉的局部血液流量循环的调制转换函数，从而使第二时间-密度曲线的剩余误差相对于所述调制转换函数与第一时间-密度曲线的卷积得到最小化；以及从所述调制转换函数计算关于各个动脉的局部血液流量循环的指数，
17.根据权利要求16所述的方法，还包含根据所述像素间的时间-密度曲线的差别来确定所述像素间的类似性。
18.根据权利要求16所述的方法，还包含根据所述像素间的时间-密度曲线的差别和所述噪声的散布期望值来确定所述像素间的类似性。
19.根据权利要求16所述的方法，还包含用一个外围像素通过加权来更换所述像素的每个像素值，所述像素构成具有加权平均的图像。
20.根据权利要求16所述的方法，其中所述调制转换函数是一个矩形函数。
21.一种指数计算设备，包括时间-密度曲线准备部件，用于从有关注射有造影剂的受试者的一个特定区域的多个连续图像，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；调制转换函数计算部件，用于通过曲线拟合来计算表示组织中关于动脉的局部血液流量循环的调制转换函数，从而使第二时间-密度曲线的剩余误差相对于所述调制转换函数与第一时间-密度曲线的卷积得到最小化；指数计算部件，用于从所述调制转换函数中计算关于各个动脉的局部血液流量循环的指数；图制备部件，用于制备所述动脉的指数图；和图合成部件，用于根据第一时间-密度曲线的剩余误差，把所述指数图合成为一个图。
22.一种设备，包括时间-密度曲线准备部件，用于从有关注射有造影剂的受试者的一个特定区域的多个连续图像，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；用于选择所述第一时间-密度曲线中与相应的第二时间-密度曲线最拟合的一条第一时间-密度曲线以指定一根动脉的装置，该被指定的动脉具有一个最高可能性，该动脉是所述组织中的每根动脉的局部血液流量循环所最可能依赖的；和图准备部件，它根据为每个组织指定的一根动脉来制备用于区分该动脉的一个依赖区的一个图。
23.一种指数计算设备，包含滤波器部件，它利用根据像素间的类似性的一个权来对多个图像进行滤波，以减少来自涉及注射有造影剂的受试者的一个特定区域的该多个连续图像的噪声；时间-密度曲线准备部件，其从带有被减少噪声的该多个连续图像中，制备关于该特定区域中的动脉的第一时间-密度曲线和关于该特定区域中的组织的第二时间-密度曲线；调制转换函数计算部件，用于通过曲线拟合来计算表示组织中关于各个动脉的局部血液流量循环的调制转换函数，从而使第二时间-密度曲线的剩余误差相对于所述调制转换函数与第一时间-密度曲线的卷积得到最小化；和一指数计算部件，用于从所述调制转换函数计算关于各个动脉的局部血液流量循环的指数。
全文摘要
从涉及注射有造影剂的受试者特定区域的许多连续图像中，准备关于特定区域中动脉的第一时间－密度曲线和关于特定区域中组织的第二时间－密度曲线。通过曲线拟合来计算调制转换函数，其表示组织中关于动脉的局部血液流量循环，由此最小化第二时间－密度曲线关于调制转换函数与第一时间－密度曲线卷积的剩余误差。从调制转换函数中计算关于各个动脉的局部血液流量循环的指数。准备所述动脉的指数图，和根据第一时间－密度曲线的剩余误差，把所述指数图合成为一个图。
文档编号A61B6/00GK1413559SQ02154738
公开日2003年4月30日 申请日期2002年10月16日 优先权日2001年10月16日
发明者南部恭二郎, 池田佳弘 申请人:株式会社东芝