CN102510282B - 一种时间分辨单光子计数二维成像系统及方法 - Google Patents
一种时间分辨单光子计数二维成像系统及方法 Download PDFInfo
- Publication number
- CN102510282B CN102510282B CN201110328462.2A CN201110328462A CN102510282B CN 102510282 B CN102510282 B CN 102510282B CN 201110328462 A CN201110328462 A CN 201110328462A CN 102510282 B CN102510282 B CN 102510282B
- Authority
- CN
- China
- Prior art keywords
- time
- single photon
- dmd
- optical
- utmost point
- 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.)
- Expired - Fee Related
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 75
- 238000000034 method Methods 0.000 title claims abstract description 43
- 238000005259 measurement Methods 0.000 claims abstract description 45
- 238000005070 sampling Methods 0.000 claims abstract description 32
- 238000001514 detection method Methods 0.000 claims abstract description 15
- 230000001960 triggered effect Effects 0.000 claims abstract description 7
- 230000003287 optical effect Effects 0.000 claims description 76
- 238000004422 calculation algorithm Methods 0.000 claims description 72
- 239000011159 matrix material Substances 0.000 claims description 61
- 238000005457 optimization Methods 0.000 claims description 49
- 230000008859 change Effects 0.000 claims description 36
- 238000005516 engineering process Methods 0.000 claims description 28
- 238000012634 optical imaging Methods 0.000 claims description 21
- 230000000694 effects Effects 0.000 claims description 16
- 230000006835 compression Effects 0.000 claims description 15
- 238000007906 compression Methods 0.000 claims description 15
- 238000010168 coupling process Methods 0.000 claims description 15
- 239000000835 fiber Substances 0.000 claims description 10
- 238000012360 testing method Methods 0.000 claims description 7
- 230000036962 time dependent Effects 0.000 claims description 6
- 230000009466 transformation Effects 0.000 claims description 6
- 230000008878 coupling Effects 0.000 claims description 5
- 238000005859 coupling reaction Methods 0.000 claims description 5
- 238000013461 design Methods 0.000 claims description 5
- 238000001914 filtration Methods 0.000 claims description 5
- 230000004313 glare Effects 0.000 claims description 5
- 238000005375 photometry Methods 0.000 claims description 4
- 229920006395 saturated elastomer Polymers 0.000 claims description 4
- 230000001360 synchronised effect Effects 0.000 claims description 4
- 238000013459 approach Methods 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 230000001186 cumulative effect Effects 0.000 claims description 3
- 239000013307 optical fiber Substances 0.000 claims description 3
- 230000004913 activation Effects 0.000 claims description 2
- 238000003860 storage Methods 0.000 claims description 2
- 230000008569 process Effects 0.000 abstract description 14
- 230000003044 adaptive effect Effects 0.000 abstract description 3
- 230000008901 benefit Effects 0.000 abstract description 3
- 238000000701 chemical imaging Methods 0.000 abstract description 2
- 238000004891 communication Methods 0.000 abstract description 2
- 238000002059 diagnostic imaging Methods 0.000 abstract 1
- 230000006872 improvement Effects 0.000 description 7
- 238000012545 processing Methods 0.000 description 7
- 239000000523 sample Substances 0.000 description 7
- 230000035945 sensitivity Effects 0.000 description 6
- 238000004088 simulation Methods 0.000 description 6
- 210000004027 cell Anatomy 0.000 description 5
- 238000002474 experimental method Methods 0.000 description 5
- 230000006870 function Effects 0.000 description 5
- 206010052128 Glare Diseases 0.000 description 4
- 230000015572 biosynthetic process Effects 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 238000001228 spectrum Methods 0.000 description 3
- 230000001052 transient effect Effects 0.000 description 3
- 238000000018 DNA microarray Methods 0.000 description 2
- 230000001174 ascending effect Effects 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 230000001066 destructive effect Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 238000009826 distribution Methods 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 230000036039 immunity Effects 0.000 description 2
- 241000894007 species Species 0.000 description 2
- 241000282552 Chlorocebus aethiops Species 0.000 description 1
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 238000010009 beating Methods 0.000 description 1
- 229910002056 binary alloy Inorganic materials 0.000 description 1
- 230000029918 bioluminescence Effects 0.000 description 1
- 238000005415 bioluminescence Methods 0.000 description 1
- 238000007697 cis-trans-isomerization reaction Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007123 defense Effects 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000000368 destabilizing effect Effects 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000005281 excited state Effects 0.000 description 1
- 210000003292 kidney cell Anatomy 0.000 description 1
- 239000007791 liquid phase Substances 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000003825 pressing Methods 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000004064 recycling Methods 0.000 description 1
- 230000011514 reflex Effects 0.000 description 1
- 230000002787 reinforcement Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 238000007614 solvation Methods 0.000 description 1
- 238000004611 spectroscopical analysis Methods 0.000 description 1
- 230000001502 supplementing effect Effects 0.000 description 1
- 230000002459 sustained effect Effects 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Landscapes
- Investigating, Analyzing Materials By Fluorescence Or Luminescence (AREA)
Abstract
本发明提供一种时间分辨单光子计数成像系统及方法,属于极弱光探测技术领域,通过触发器2触发开始采样,每隔t时间间隔集中采样,在该间隔有进光就开始测量计数,从而实现对极弱光对象的时间分辨,生成时间序列图像。成像主要基于压缩传感(Compressive Sensing,简称CS)理论,利用数字微镜器件(DigitalMicromirror Device,简称DMD5)对可压缩的二维图像进行线性随机投影,经光学调制后用单光子计数器同步探测,仅需少量采样便可实现高分辨率的极弱光图像重建,其测量过程是线性的、非自适应的,而重建过程是非线性的,有更好的通用性、健壮性、可扩展性、叠加性和计算非对称性,可广泛应用于生命科学、医疗成像、数据采集、通信、天文、军事、超光谱成像、量子测量等领域。
Description
技术领域
本发明涉及极弱光探测的技术领域,特别涉及一种时间分辨单光子计数二维成像系统及方法。
背景技术
本发明是基于前人工作进行的改进和创新。在该领域,本研究所已有一份专利《一种单光子计数成像系统及其方法》(申请号或专利号:201110103559.3,申请人或专利权人:中国科学院空间科学与应用研究中心),该专利是本研究室前人所做工作,目的是做单光子计数成像,其特征在于,该单光子计数成像系统采用压缩传感理论和DLP技术,并以单光子计数器为探测元件,实现了单光子级别的极弱光对象的二维成像,但还存在一定的技术缺陷,如该专利所提供的系统没有设置触发器,仅能对静态物体成像,缺乏时间分辨,且所用算法较慢,没有考虑系统噪声对图像重建质量的影响,图像重建时间较长、效果较差,观察对象描述不明确,没有分类讨论,透射物体和反射物体的成像没有区分处理的方案,光衰减器和滤光片的设置条件没有明确说明,系统结构图尚还存在缺陷,没有考虑计数系统与DMD的同步问题,尚无DMD的驱动控制,系统装置尚不完善(有一些技术上的漏洞),实验条件不成熟,计数成像技术尚处在探索阶段。现基于此,特提出基于压缩传感理论的时间分辨单光子计数成像系统,以解决上述一系列的缺陷。此外,该系统与本所申请的另一项专利《一种时间分辨极弱光多光谱成像系统及方法》的区别在于,该系统采用了光衰减器和单光子计数器,成本较低,所成像为灰度图,适用于一些只需知道被测物体大致轮廓而不做颜色分析的场合,针对性更强一些,如在天文、军事等领域都会大展身手。
所谓时间分辨就是能观察物理和化学的瞬态过程并能分辨其时间,在液相中,很多物理和化学过程,如分子的顺-反异构和定向弛豫、电荷和质子的转移、激发态分子碰撞预解离、能量传递和荧光寿命以及电子在水中溶剂化等,仅需10-8秒就能完成。只有皮秒激光脉冲才有可能及时地观察这些极快的过程。在本发明中,就想在极端的时间间隔内进行单光子计数探测,实现高速采样。
当光强衰减到一定程度,达到单光子水平,就变成了离散脉冲信号。单光子是一种极微弱光,被认为是光不可分割的最小能量单位,是可以探测的极限。单光子探测技术应用于生物自发光、医疗诊断、非破坏性物质分析、天文观测、空间科学、高速现象检测、高分辨光谱测量、量子光学等领域,并在其中扮演着重要角色。研究极弱光成像探测技术有着很好的发展前景。
单光子计数方法利用弱光照射下光子探测器输出电信号自然离散的特点,采用脉冲甄别技术和数字计数技术把极其微弱的信号识别并提取出来,这种方法受不稳定因素的影响较小,可消除探测器大部分热噪声的影响,大大提高了测量结果的信噪比,并且可以输出数字信号,适合与微机接口连接进行数字数据处理。
而单光子计数成像就是一种极弱光探测技术,通常它通过记录成像位置的光子计数以及探测到光子的概率,在数据处理端进行累计和融合获得一幅图像,其核心是面元探测器,面元探测器规模(阵列大小)、灵敏度范围、以及响应波段直接影响能否获得单光子水平的图像获取质量。但是,用于单光子探测水平的面元探测器不但价格昂贵,只能够在少数波段可以实现,且面元探测器灵敏度低,存在技术不甚成熟与极弱光二维成像的强烈需求之间的矛盾。而点探测器无论在探测灵敏度、波长范围具有更宽的选择范围,成本优势明显,利用点探测器实现单光子计数成像成为未来单光子水平成像的重要发展趋势。
压缩传感理论(CS理论)是由E.J.Candès等人提出的,它打破了传统的线性采样模式,表明可压缩信号的少量线性随机投影中包含足够的信息来重建原信号。
CS理论包括两部分:压缩采样和稀疏重建。
1)假设x∈Rn是被测数据,y∈Rk是观测数据,Φ∈Rk×n是随机投影矩阵(k<<n),e∈Rk是测量噪声,那么,压缩采样的过程可以描述为(1)式:
y=Φx+e (1)
其中Φ满足RIP准则;如果x是变换域稀疏的,即θ=Ψx,Ψ是稀疏变换矩阵,那么(1)式变化为(2)式:
y=ΦΨ-1θ+e (2)
Φ与Ψ越不相关,采样所需的测量数k越小,计算量越小,且ΦΨ-1需满足RIP准则,所以在本发明中,(2)式中Ψ为小波变换矩阵,Φ为伪随机二值矩阵;
2)假设测量数为k,二维图像的像素总数为n,则(1)式中的测量矩阵则为Φ={Φ1,…,Φi,…,Φk},Φi是Φ的第i行,把a×b像素的二维图像的列首尾相连,化成n×1(其中n=a×b)的一维列向量,对应(1)式中的x,其中的每一个元素代表相应位置处的光子密度。主流的DMD由1024×768的阵列构成,它的列首尾相连,化成1×n的一维行向量,对应测量矩阵Φ中的一行,其中的每个元素代表相应位置处光子透射到聚焦系统的概率,而测量矩阵Φ共计k维,即k行n列的矩阵;
3)假设测量周期为t,即每段的时间间隔,在这段时间间隔内,驱动控制模块确保微镜阵列的随机翻转,单光子计数器探测到的光子数为N,相当于光子密度图像与DMD上的随机测量阵列的内积值,对应于(1)式中观察向量y的一个元素,式中,Φi,j、xj分别是Φi和x的第j个元素。根据测量矩阵,由驱动控制模块控制每次测量时DMD微镜阵列的排布,重复k次测量,就可以得到该t时间间隔内整组观测数据y(y为k×1的一维列向量);
4)所述的稀疏重建是在已知观测数据y和测量矩阵Φ的条件下求解(1)式中的x,一般用最优化算法求解,可描述为(3)式:
如果x是变换域稀疏的,对应于(2)式的重建问题可以描述为(4)式:
(3)式和(4)式中,第一项是最小二乘约束,记为f(x);第二项是对x稀疏度的一种约束;两项之和是目标函数,记为
上述“先采样压缩,后重建”的思想使得将二维信号转换为随时间分布的一维信号,使得时间分辨单光子计数成像成为可能。
DLP(Digital Light Processing)数字光处理技术是美国德州仪器公司(TI)提出的一项技术,即先将影像信号经过数字处理,然后再把光投影出来。其核心是DLP芯片——数字微镜器件(Digital Micro-mirror Device,简称DMD),这目前是世界上最精密的光开关。DMD是一种极小的反射镜,它是包含有成千上万个安装在铰链上的微镜的矩阵(主流的DMD由1024×768的阵列构成),每个微镜的大小小于人的头发丝的五分之一,每一个镜片可以通断一个像素的光,这些微镜皆悬浮着,并均可以静电方式向两侧倾斜10-12°左右(这里取+12°和-12°),把这两种状态记为1和0,分别对应“开”和“关”,当镜片不工作时,它们处于0°的“停泊”状态。对每一个镜片下的存储单元都以二进制平面信号进行电子化寻址。决定每个镜片倾斜在哪个方向上为多长时间的技术被称为脉冲宽度调制(PWM)。镜片可以在一秒内开关1000多次,这一相当快的速度允许数字灰度和颜色再现。
发明内容
本发明的目的在于,为解决能时间分辨地观察物理化学生物瞬态过程的强烈需求、目前面元探测器灵敏度低、技术不甚成熟与极弱光二维成像的强烈需求之间的矛盾,从而提供一种时间分辨单光子计数二维成像系统及方法。
为实现上述目的,本发明提供一种时间分辨单光子计数二维成像系统,该系统主要基于压缩传感理论,用于对随时间动态变化的物体成像,输出按时间序列排列的连续灰度视频图像帧,所述系统包含:触发器、光学成像系统、DMD、光学聚焦收集系统、光衰减器、单光子计数器、驱动控制模块和最优化算法模块;
所述触发器由位于其前端的极弱光触发源触发,该触发器输出端与驱动控制模块的输入端相连,当触发器被触发时驱动控制模块将输出驱动控制信号触发与其输出端相连的DMD和单光子计数器开始工作,则该DMD开始翻转,单光子计数器同时开始计数,光衰减器用于对光强度进行衰减,光学聚焦收集系统用于对光线的聚焦收集,单光子计数器是对光进行单光子计数;
所述单光子计数器的输出端和驱动控制模块的一个输出端,均与所述最优化算法模块的输入端相连,作为最优化算法模块的输入,则该最优化算法模块的作用即根据输入的单光子计数器的测量值与从驱动控制模块里导出的选定区域的随机测量矩阵来重建出稀疏信号,反演出光子密度图像,经过M个时间间隔,就能重建出按时间序列排列的可时间分辨的一系列二维灰度图像视频帧;
所述DMD为数字微镜器件,包括微镜阵列和相关集成电路部分(属该器件内部结构),都封装在一块板子上;
其中,所述驱动控制模块基于DLP数字光处理技术,驱动控制模块在选定DMD工作区域后,便下载伪随机测量矩阵来驱动控制所述DMD微镜阵列的翻转;所述DMD微镜阵列在翻转的同时会向所述单光子计数器发送同步信号,保证DMD微镜阵列与单光子计数器之间同步,即DMD微镜阵列每翻转一次,单光子计数器累积计数在该次翻转的时间间隔内的光单子数,DMD微镜阵列翻转完成后,单光子计数器清零重新开始累积计数,所有计数和该选定区域的随机测量矩阵均传送至最优化算法模块中。
可选的,所述极弱光触发源为:极弱光光源或自发光生物体。
可选的,当所述极弱光触发源采用极弱光光源时,其投射方式主要有两种:1)极弱光光源直接透射物体对象上,光源与光学成像系统光轴在同一条直线上,物体对象可以是半透明的或者镂空的,通过透射光成像;2)极弱光光源斜着打在物体上,光源与光学成像系统光轴不在同一条直线上,通过物体对象反射极弱光来成像;
当所述极弱光光源比较强时采用滤光片滤除该极弱光的杂光,若极弱光光源的光强极其微弱,且其波长在本发明探测所要求的波长范围内,则不需要再设置滤光片。
所述光学成像系统和光学聚焦收集系统均采用光学透镜组,分别负责光学成像和光学聚焦,极弱光通过光学成像系统后,可在DMD上成等大或缩小或放大的像,按实际需求进行成像调整;
其中,所述光学聚焦收集系统采用光纤耦合技术,即将经由所述分光光度计分光后的光束耦合到光纤中,利用光纤耦合技术将分光分别收集到对应的单光子计数器上。
可选的,所述的光衰减器是由多块衰减片组成,放置在光学聚焦收集系统和单光子计数器之间的光路上或光源处,用于将光衰减到单光子计数器的工作范围;该光衰减器的设计是为了防止被测光子密度过大和单光子计数器的门控时间过长引起的饱和,若被测光已经达到单光子级别,则不需要再设置光衰减器了;若将光衰减片放置在光源处,这可将光源变成单光子源,这种方式噪声会略大。
基于上述系统本发明还提供一种时间分辨单光子计数二维成像方法,该方法对随时间动态变化的物体成像,输出按时间序列排列的连续灰度视频图像帧,所述方法包含如下步骤:
步骤1,用于采用极弱光触发源触发触发器进行触发同步的步骤,该步骤实现了时间分辨的效果;
步骤2,用于将被测信号由高维向低维进行映射的步骤,该步骤采用压缩和采样对被测信号进行由高维向低维进行映射;
步骤3,用于对待输出的视频帧进行稀疏重建最优化,输出按时间序列排列的连续灰度视频图像帧的步骤。
所述步骤1进一步包含如下子步骤:
用于触发的步骤,所述极弱光光源或发光生物体触发触发器,该触发器进而触发激活整个成像系统的工作,所述触发器触发之后每隔t时间间隔进行一组单光子探测,最优化算法模块自动识别该t时间间隔内的极弱光对象二维图像的输入,分别累积测量该t时间间隔内的单光子数;
用于同步步骤,DMD在翻转的同时向单光子计数器发送同步信号,即DMD每翻转一次,单光子计数器累积计数在该次翻转的时间间隔内的单光子数,DMD翻转完成后,单光子计数器清零重新开始累积计数,所有的计数被传至最优化算法模块中;
所述的时间分辨,可以针对连续变化观测物体和周期变化观测物体,对于前者而言,采用上述每t时间间隔采样k次并持续测量M×t时间的方法即可实现,对于后者而言,通常采用的观测物体的变化周期极短,假设周期为T,将该时间周期等分为d个时间间隔,记做t1,t2,t3,…,td,在该周期T内保持对应随机矩阵不变,到下一个周期随机矩阵才发生改变,在每个小的时间间隔t1,t2,t3,…,td内分别计数,即对落在ti时间间隔内的单光子进行累积计数,依靠触发器2的精确触发,保证计数和随机矩阵的严格对应关系,测量k个周期,即对每个微小时间间隔测量了k次,分别对这d个时间间隔做最优化重建,便可反演出一个时间周期序列内的灰度观测对象的变化情况。
所述步骤2进一步包含如下子步骤:
用于压缩的步骤,DMD微镜阵列将可压缩的随时间动态变化的灰度物体的图像帧数据随机反射,当DMD中单个微镜+12°翻转时反射的光被单光子计数器接收,当DMD中单个微镜-12°翻转时反射光不能被单光子探测器接收,从而完成被对随时间动态变化的灰度物体的被测信号的压缩,同时确保DMD微镜阵列明暗阵列的最大随机性,进而控制极弱光被反射至光学聚焦收集系统的概率是随机的;
用于分光采样的步骤,经光学聚焦系统聚焦后的极弱光进入光衰减器衰减后再进入单光子计数器,单光子计数器再对极弱光进行探测采样。
所述步骤3进一步包含如下子步骤:
单光子计数器每t时间间隔内进行计数,将该计数值作为单光子技术器的测量值输入最优化算法模块;
最优化算法模块根据上步测量值、DMD微镜阵列上的测量矩阵和加载在二维图像上的稀疏矩阵,分别通过最优化算法重建出光子密度图像,反演出该时间间隔t内随时间动态变化的物体的灰度图;
重复执行上述两个步骤共M次,得到M×t时间段的M幅二维灰度图像帧序列,输出视频帧。
所述最优化算法为运用小波进行稀疏变换的可分离逼近的稀疏重建算法SpaRSA-DWT,该算法是基于压缩传感理论编写的,在现有的IST算法(该算法为迭代算法)基础上进行改进,利用小波变换对单光子计数测量值进行稀疏化,并调整原有IST算法每次迭代的步长系数αt,使αtI逼近f(x)在xt处的Hessian矩阵,通过自适应的步长系数提高了原有算法的迭代速度,并自适应地修改阈值,依靠反复迭代运算求解出对应的稀疏信号,最后反演出二维图像,经过M×t时间后,共读入M幅极弱光对象动态二维灰度图像,对应输出M幅随时间变化的二维灰度图像视频帧,以观察原对象变化情况,并计算出相关系数和峰值信噪比,用于比较原图和重建图稀疏系数的吻合程度。本发明比较了IST、TV、OMP、MP、StOMP、CoSaMP、LBI、SP、Bayesian、l1_ls、smp等算法,发现本算法的适用性更好、抗噪性更好、重建所用时间更短、重建图像的清晰度和对比度更高;
上述最优化算法在最优化算法模块里执行,输入为计数测量值、驱动控制模块里导出的随机测量矩阵,输出为M幅随时间变化的灰度二维重建图像视频帧,核心重建问题可描述为:
与现有技术相比本发明的优点在于:本发明综合触发器触发技术、压缩传感(Compressive Sensing,简称CS)理论、DLP数字光处理(Digital Light Processing,简称DLP)技术、光纤耦合技术和单光子计数器探测技术后提出的方案可以解决利用点探测器实现高探测灵敏度的时间分辨成像的问题,其灵敏度可以达到单光子水平,分辨率与DMD直接相关,而DMD目前可以达到很高的分辨率,有的已经可达2048×1152的分辨率。本发明解决了目前该领域中焦平面传感器灵敏度低、阵列规模小、探测波长范围相对单一与时间分辨极弱光二维成像的强烈需求之间的矛盾,改变了以往单一的二维成像,创新地加入第三维时间轴,变成按时间序列排列的连续二维图像帧输出,由于时间间隔很短,可制作成视频输出,实现了时间分辨,满足了目前对物理、化学、生物瞬态过程的时间分辨观察的强烈需求。本发明添置了触发器,能触发DMD和计数器同时开始运作,并可实现对连续变化的物体成像,具有时间分辨率。本发明改进了原有算法(IST、TV、OMP、MP、StOMP、CoSaMP、LBI、SP、Bayesian、l1_ls、smp等算法),使其通用性、鲁棒性更强,运算更快,对比度更高,抗噪性能更佳。本发明新置极弱光光源,并对不同观测对象都有明确的观测方案,更具有针对性。本发明中,光衰减器和滤光片的灵活设置也更能贴近实际需求。此外,本发明中的光学成像系统和光学聚焦收集系统,皆为光学透镜组,不再局限于原有技术的简单的单个透镜成像和聚焦,适用范围更广,所成像面积更大,可在DMD上成等大或缩小或放大的像,并创新地采用光纤耦合技术收集光线,使得成像质量更高。本发明还实现了对DMD的驱动控制,填补了原有技术的空白,解决了计数系统与DMD之间的同步问题。本发明通过一系列的实验,也进一步完善了计数成像技术。本发明可广泛应用于生物自发光检测、医疗成像、数据采集、遥感通信、非破坏性物质分析、天文观测、国防军事、超光谱测量、量子测量等领域。
附图说明
图1是本发明基于压缩传感理论的时间分辨单光子计数成像系统的结构示意图;
图2是当观测对象为自发光的物体时的时间分辨单光子计数成像系统的示意图;
图3(a)是对本发明的软件模拟实验结果;
图3(b)是对本发明的软件模拟实验结果;
图3(c)是对本发明的软件模拟实验结果;
图3(d)是对本发明的软件模拟实验结果;
图4(a)是对本发明的软件模拟实验结果;
图4(b)是对本发明的软件模拟实验结果;
图4(c)是对本发明的软件模拟实验结果;
图4(d)是对本发明的软件模拟实验结果;
图5(a)是对本发明的软件模拟实验结果;
图5(b)是对本发明的软件模拟实验结果;
图5(c)是对本发明的软件模拟实验结果;
图5(d)是对本发明的软件模拟实验结果;
图6(a)是对本发明的软件模拟实验结果;
图6(b)是对本发明的软件模拟实验结果;
图6(c)是对本发明的软件模拟实验结果;
图6(d)是对本发明的软件模拟实验结果;
图7(a)是对本发明的软件模拟实验结果;
图7(b)是对本发明的软件模拟实验结果;
图7(c)是对本发明的软件模拟实验结果;
图7(d)是对本发明的软件模拟实验结果。
附图标识
1、极弱光光源 2、触发器 3、滤光片
4、光学成像系统 5、DMD微镜阵列 6、光学聚焦收集系统
7、光衰减器 8、单光子计数器 9、驱动控制模块
10、最优化算法模块
具体实施方式
以下结合附图对本发明作进一步的详细说明。
本发明技术方案以压缩传感(Compressive Sensing,CS)理论为基础,在现有技术的基础上,创新地采用触发器进行触发,并控制时间间隔,利用DLP技术将图像信号随机投影,转化为随机的光强信号,再利用单光子计数器作为探测元件,探测到的光子数作为测量值,最后最优化算法模块对其进行重建,由于时间间隔很短,重建时间稍长,采用先集中采样后批量重建计算的方式,实现了用点探测器对对象进行时间分辨的极弱光二维成像。
为实现上述目的,本发明构建了一种新的时间分辨单光子计数成像系统和方法。
这种时间分辨单光子计数成像系统的特征在于,该系统基于触发器触发技术、压缩传感理论、DLP数字光处理技术、光纤耦合技术和单光子计数器探测技术,创新地由触发器2触发之后才开始探测计数,每间隔t时间间隔(按实际需求设定,在本发明中,创新之处就在于能在极端的时间间隔内进行单光子计数探测,比如纳秒级或皮秒级时间间隔,即实现高速采样,并具备时间分辨性)采样一组数据(一组数据对应一幅图像),对其二维图像进行投影,运用最优化算法进行重建,实现单光子级别的极弱光二维成像的时间分辨,以便观测对象的动态变化和后续的研究分析。
所述的时间分辨单光子计数成像系统主要由极弱光光源1、触发器2、滤光片3、光学成像系统4、DMD5、光学聚焦收集系统6、光衰减器7、单光子计数器8、驱动控制模块9和最优化算法模块10构成。
其中,新添加的极弱光光源1十分必要,它能将极弱光打在物体上,主要有两种方式:1)极弱光光源1直接透射物体对象上,光源与光学成像系统4光轴在同一条直线上,物体对象可以是半透明的或者镂空的,这样通过光子密度变化便能成像;2)极弱光光源1斜着打在物体上,光源与光学成像系统4光轴不在同一条直线上,通过物体对象反射极弱光来成像。
极弱光光源1采用普通极弱光光源或生物荧光或星光等,生物芯片是其中一种典型的极弱光源,目前主要通过荧光标记的方法使其便于观察,实际上生物都有自发光的特性,并且自发光光谱包含很多重要的信息,采用时间分辨单光子计数成像技术就可以直接观测,因而当观测对象为自发光的生物时,就可以去掉时间分辨单光子计数成像系统中的极弱光光源1,将触发器2和滤光片3平移到自发光的连续随时间变化的生物体右侧的光路上,进行直接观测,这也是对现有技术的补充。
新添置的触发器2与驱动控制模块9相连,触发器2当有极弱光射进来便开始触发,驱动控制模块9便发送驱动控制信号来使DMD5运作,单光子计数器8开始探测计数。简言之,触发器2的作用即告知驱动控制模块9触发信号,DMD5和计数器开始同时工作。
本发明又一改进之处在于,滤光片3的作用是滤除极弱光的杂光,使进入后续系统的极弱光在探测所需的波长范围内,一般在光强比较强的时候使用。若极弱光光源1的光强极其微弱,且其波长在本发明探测所要求的波长范围内,则不需要再设置滤光片3。
在现有技术上,本发明另一改进创新是添加了驱动控制模块9,这是基于DLP数字光处理技术,针对DMD5数字微镜器件的,即通过驱动控制模块9里的驱动程序编写伪随机测量矩阵,驱动其微镜的翻转,通过设置延时可调节灰度值。DMD5可以翻转+12°和-12°(有些型号的DMD5可以翻转+10°和-10°),在本系统中,设置+12°为能接收到的反射角度,-12°翻转能进入最后的单光子计数器8的反射光微乎其微,可忽略不计,因而光学聚焦收集系统6的主光轴与光学成像系统的主光轴所成夹角为24°。本发明可以自动生成601帧(帧数可设置)随机数文件,选定DMD5工作区域后,驱动控制模块9里的驱动控制程序下载并生成该区域的随机数文件,控制DMD5的随机翻转,DMD5在翻转的同时会向计数器发送同步信号,这保证了DMD5与单光子计数器8之间同步,即DMD5每翻转一次,单光子计数器8累积计数在该次翻转的时间间隔内的光单子数,DMD5翻转完成后,计数器清零重新开始累积计数,所有计数都会通过数据线传到最优化算法模块10上,存在txt文档中。
本发明又一改进创新之处在于,所述的光学成像系统4和光学聚焦收集系统6,皆为光学透镜组,分别负责光学成像和光学聚焦,而不再局限于原有技术的简单地用单个透镜进行成像和聚焦,因而本发明的适用范围更广,在DMD5上成像面积更大,最大可至768×1024像素,不同于现有的技术,在所述的光学成像系统4中,极弱光通过光学成像系统4后,将在DMD5上成等大或缩小或放大的像,并可按实际需求进行成像设置。而后续的光路收集相比现有技术更为复杂,涉及到光纤耦合技术,即将聚焦后的光束耦合到光纤里。创新点在于使用光纤耦合技术将聚焦光收集到单光子计数器8上,耦合的好坏直接影响成像质量,因而调节耦合也变得更加困难。
以单光子计数器8(即计数型单光子探测器)作为探测元件,在每t时间内对光子进行计数,传输到最优化算法模块10上,将这计数作为测量值,然后由压缩感知理论中的最优化算法重建出稀疏信号,反演出光子密度图像,重建出按时间序列排列的可时间分辨的一系列二维图像帧,以达到实时观测对象变化的目的。
对原有技术还有一个改进在于,所采用的光衰减器7是由多块衰减片组成,放置在光学聚焦收集系统6和单光子计数器8之间的光路上或光源处,用于将光衰减到单光子计数器的工作范围。该光衰减器7的设计是为了防止被测光子密度过大和单光子计数器8的门控时间过长引起的饱和。若被测光已经达到单光子级别,则不需要再设置光衰减器7了。若将光衰减片放置在光源处,这可将光源变成单光子源,这种方式噪声会略大。
为实现上述的另一发明目的,本发明还提供了一种新的时间分辨单光子计数成像方法,该方法采用了触发器触发技术、压缩传感理论、DLP数字光处理技术、光纤耦合技术和单光子计数器探测技术。利用最优化算法重建出连续变化的图像帧,实现了单光子级别的时间分辨极弱光二维成像。该方法可以远低于奈奎斯特采样频率进行采样,将压缩和采样合并在一起,采样次数大大减少,将后续的重建工作留给最优化算法,这种非对称的结构保证了时间分辨单光子计数成像的健壮性、可扩展性、叠加性,其测量过程是线性的、非自适应的,而重建过程是非线性的。
其具体步骤包括:
步骤1触发同步:通过触发器2触发来激活整个系统的工作,一旦触发后,系统每隔t时间周期进行一组单光子计数,最优化算法模块10自动识别该t时间间隔内的极弱光二维图像的输入,累积测量该t时间间隔的单光子数,这是至关重要的一个环节,关系到整个系统的同步问题;而DMD5在翻转的同时会向计数器发送同步信号,这保证了DMD5与单光子计数器8之间同步,即DMD5每翻转一次,单光子计数器8累积计数在该次翻转的时间间隔内的光单子数,DMD5翻转完成后,计数器清零重新开始累积计数,所有计数都会通过数据线传到最优化算法模块10上。所述的时间分辨,可以针对连续变化观测物体和周期变化观测物体,对于前者而言,采用上述每t时间间隔采样k次并持续测量M×t时间的方法即可实现,对于后者而言,通常采用的观测物体的变化周期极短,假设周期为T,将该时间周期等分为d个时间间隔,记做t1,t2,t3,…,td,在该周期T内保持对应随机矩阵不变,到下一个周期随机矩阵才发生改变,在每个小的时间间隔t1,t2,t3,…,td内分别计数,即对落在ti时间间隔内的单光子进行累积计数,依靠触发器2的精确触发,保证计数和随机矩阵的严格对应关系,测量k个周期,即对每个微小时间间隔测量了k次,分别对这d个时间间隔做最优化重建,便可反演出一个时间周期序列内的灰度观测对象的变化情况;
步骤2压缩采样:DMD微镜阵列将可压缩的随时间动态变化的物体的图像帧数据随机反射,当DMD中单个微镜+12°翻转时反射的光被单光子计数器接收,当DMD中单个微镜-12°翻转时反射光不能被单光子探测器接收,从而完成被对随时间动态变化的物体的被测信号的压缩,同时确保DMD微镜阵列明暗阵列的最大随机性,进而控制极弱光被反射至光学聚焦收集系统的概率是随机的;
步骤3稀疏重建:单光子计数器8在每t时间间隔内对光子进行计数,该计数作为测量值,最优化算法模块10根据这组测量值和DMD5上的测量矩阵、加载在原二维图像上的稀疏矩阵,通过最优化算法重建光子密度图像,反演出二维图像,M×t时间后,也就重建出了M幅二维图像,这样也就实现了时间分辨的效果。
作为上述技术方案的一种改进,所述的压缩采样和稀疏重建方法包括如下具体步骤:
压缩采样,是被测信号由高维向低维映射的过程,而稀疏重建是最优化的问题。
1)假设x∈Rn是被测数据,y∈Rk是观测数据,Φ∈Rk×n是随机投影矩阵(k<<n),e∈Rk是测量噪声,那么,压缩采样的过程可以描述为(1)式:
y=Φx+e (1)
其中Φ满足RIP准则;如果x是变换域稀疏的,即θ=Ψx,Ψ是稀疏变换矩阵,那么(1)式就变化为(2)式:
y=ΦΨ-1θ+e (2)
在本发明方案中,(2)式中Ψ为小波变换矩阵,Φ为伪随机二值矩阵,ΦΨ-1满足RIP准则;
2)假设测量数为k,二维图像的像素总数为n,则(1)式中的测量矩阵则为Φ={Φ1,…,Φi,…,Φk},Φi是Φ的第i行,把a×b像素的二维图像的列首尾相连,化成n×1(其中n=a×b)的一维列向量,对应(1)式中的x,其中的每一个元素代表相应位置处的光子密度。主流的DMD5由1024×768的阵列构成,它的列首尾相连,化成1×n的一维行向量,对应测量矩阵Φ中的一行,其中的每个元素代表相应位置处光子透射到聚焦系统的概率,而测量矩阵Φ共计k维,即k行n列的矩阵;
3)假设测量周期为t,即每段的时间间隔,在这段时间间隔内,驱动控制模块9确保微镜阵列的随机翻转,单光子计数器8探测到的光子数为N,则N相当于光子密度图像与DMD5上的随机测量阵列的内积值,对应于(1)式中观察向量y的一个元素,式中,Φi,j、xj分别是Φi和x的第j个元素。根据测量矩阵,由驱动控制模块9控制每次测量时DMD5微镜的排列,重复k次测量,就可以得到该t时间间隔内整组观测数据y(y为k×1的一维列向量);
4)所述的稀疏重建是在已知观测数据y和测量矩阵Φ的条件下求解(1)式中的x,一般用最优化算法求解,可描述为(3)式:
如果x是变换域稀疏的,对应于(2)式的重建问题可以描述为(4)式:
(3)式和(4)式中,第一项是最小二乘约束,记为f(x);第二项是对x稀疏度的一种约束;两项之和是目标函数,记为
作为原有技术方案的进一步的改进,本发明做了软件模拟,创新之处在于,所述最优化算法为运用小波进行稀疏变换的可分离逼近的稀疏重建算法SpaRSA-DWT,该算法是基于压缩传感理论编写的,在现有的IST算法(该算法为迭代算法)基础上进行改进,利用小波变换对单光子计数测量值进行稀疏化,并调整原有IST算法每次迭代的步长系数αt,使αtI逼近f(x)在xt处的Hessian矩阵,通过自适应的步长系数提高了原有算法的迭代速度,并自适应地修改阈值,依靠反复迭代运算求解出对应的稀疏信号,最后反演出二维图像,经过M×t时间后,共读入M幅极弱光对象动态二维灰度图像,对应输出M幅随时间变化的二维灰度图像视频帧,以观察原对象变化情况,并计算出相关系数和峰值信噪比,用于比较原图和重建图稀疏系数的吻合程度。本发明比较了IST、TV、OMP、MP、StOMP、CoSaMP、LBI、SP、Bayesian、l1_ls、smp等算法,发现本算法的适用性更好、抗噪性更好、重建所用时间更短、重建图像的清晰度和对比度更高;为了说明本发明设计方案的可行性,这里计算重建图像与原图像的相似程度,用相关系数来评价:
并且计算重建图像的信噪比,用PSNR(Peak Signal Noise Ratio)来评价:
其中,Xi,j和分别表示原图像和重建图像第i行第j列的像素点光子密度值,n2表示图像大小(这里假设图像像素是n×n的,这里的n与上述第2)中的n属不同概念,只是个记号),PSNR值越大,就代表失真越少。
根据上述时间分辨单光子计数成像方法,本发明具体实施方式如下:
如图1所示,极弱光光源1打出极弱光触发触发器2,使得后续的装置开始正常工作,然后每隔t时间间隔(该时间间隔可以极其短,期望达到皮秒甚至纳秒级别,以便实现时间分辨)进行一次集中采样,每次集中采样k次。每次时间间隔内,该极弱光经由滤光片3滤除杂光之后打在连续变化的对象上(直接透射则后续是透射光,斜射则后续是反射光),接着经光学成像系统4在DMD5上成像,利用DLP数字光处理技术调制DMD5随机明暗矩阵,来控制极弱光光子被反射到后续光学聚焦收集系统6的概率,使其尽可能随机,每t时间间隔调制k次(即每幅图像的测量数),这些经+12°反射的极弱光光子经由光学聚焦收集系统6聚焦到一点,设置在单光子计数器8之前的光衰减器7的作用就是当光过强时将光衰减到单光子计数器的工作范围,然后由单光子计数器8在一定时间内对光子进行计数,该计数作为测量值。最后由最优化算法模块10根据测量值、DMD5上k行n列的测量矩阵以及稀疏变换矩阵经过最优化算法重建出光子密度图像。最后形成时间分辨的重建二维图像帧序列。其中,光衰减片放置在光学聚焦收集系统6和单光子计数器8之间的光路上或光源处,该设计是为了防止被测光子密度过大和单光子计数器8的门控时间过长引起的饱和,可将光衰减到单光子计数器的工作范围,一般组合多块衰减片对光进行衰减。若被测光已经达到单光子级别,则不需要再设置光衰减片了。若将光衰减片放置在光源处,这可将光源变成单光子源,这种方式噪声会略大。
需要说明的是,当观测对象为自发光的生物时,就可以去掉时间分辨单光子计数成像系统中的极弱光光源1,将触发器2和滤光片3一起平移到自发光的连续随时间变化的生物体右侧的光路上,利用时间分辨单光子计数成像方法直接观测,如图2所示。
其中,滤光片3的作用是滤除极弱光的杂光,使进入后续系统的极弱光在探测所需的波长范围内,一般在光强比较强的时候使用。若极弱光光源1的光强极其微弱,且其波长在本发明探测所要求的波长范围内,则不需要再设置滤光片3。
所述的时间分辨,可以针对连续变化观测物体和周期变化观测物体。对于前者而言,采用上述每t时间间隔采样k次并持续测量M×t时间的方法即可实现。对于后者而言,通常采用的观测物体的变化周期极短,假设周期为T,将该时间周期等分为d个时间间隔,记做t1,t2,t3,…,td,在该周期T内保持对应随机矩阵不变,到下一个周期随机矩阵才发生改变,在每个小的时间间隔t1,t2,t3,…,td内分别计数,即对落在ti时间间隔内的单光子进行累积计数,依靠触发器2的精确触发,保证计数和随机矩阵的严格对应关系。测量k个周期,即对每个微小时间间隔测量了k次,分别对这d个时间间隔做最优化重建,便可反演出一个时间周期序列内的观测对象变化情况。
另外需要说明的是DMD5的反射机制,正如图1图2中所示,当入射光线与DMD5单个微镜法线成24°时,反射光线也与法线成24°,但当微镜翻转+12°时,图中DMD5微镜的法线也顺时针翻转+12°,那根据反射定律反射光线需顺时针翻转+24°,即与初始位置时的法线在同一直线上,由于本发明的时间分辨单光子计数成像系统的DMD5前后光学成像系统4和光学聚焦收集系统6的光轴之间设定的角度正是24°,那么微镜翻转+12°时,光线能正常反射到后续光学聚焦收集系统6。同理,当微镜翻转-12°时,这时的反射光线与初始位置时的法线成-48°,几乎不进入后续的光学聚焦收集系统6,因而微镜翻转-12°时的反射光可以忽略不计。这里取顺时针翻转为正,逆时针为负。
假设测量数为k,二维图像的像素总个数为n,则(1)式中的测量矩阵则为Φ={Φ1,…,Φi,…,Φk},Φi是Φ的第i行,把a×b像素的二维图像的列首尾相连,化成n×1(其中n=a×b)的一维列向量,对应(1)式中的x,其中的每一个元素代表相应位置处的光子密度。令DMD5微镜阵列的列首尾相连,化成1×n的一维行向量,对应测量矩阵Φ中的一行,其中的每个元素代表相应位置处光子透射到聚焦系统的概率,而测量矩阵Φ共计k维,即k行n列的矩阵;假设测量周期为t,即每段的时间间隔,在这段时间间隔内,驱动控制模块9确保微镜阵列的随机翻转,单光子计数器8探测到的光子数为N,则N就相当于光子密度图像与DMD5上的随机测量阵列的内积值,对应于(1)式中观察向量y的一个元素,式中,Φi,j、xj分别是Φi和x的第j个元素。根据测量矩阵,由驱动控制模块9控制每次测量时DMD5微镜的排列,重复k次测量,就可以得到该t时间间隔内整组观测数据y(y为k×1的一维列向量),在物理上实现(1)式的过程。
根据光子学的知识,在一个元面积dA内,任意时刻在r点观察到一个光子的概率p(r)dA正比于该处光强。
现在用软件模拟光子密度图像,采用的原始随时间连续变化的观测对象是根据Sara A Jones,Sang-Hee Shim,Jiang He & Xiaowei Zhuang所著的《Fast,three-dimensional super-resolution imaging of live Cells》一文中所提供的nmeth.1605-S2.mov截取的,截取时采用先截视频帧然后matlab精确定位截取某个观测区域块的方式来获取所需要的模拟实验对象数据。该视频是在STORM成像条件下获取的一个细胞的DIC图像序列,而这个细胞是BS-C-1非洲绿猴肾细胞株,视频记录了该细胞在657nm激光15kW/cm2的持续辐射下20分钟内的连续变化。实际上许多生物体都有自发光的特性,并且自发光光谱包含很多重要的信息。采用时间分辨光子计数成像技术就可以直接观测。生物芯片就是个典型的例子,目前主要通过荧光标记的方法使其便于观察。因而上述获取的图像序列信息是可靠可用的,符合本发明观测对象的特点,特借此作为模拟用实验数据。实验调试软件环境是matlab,所截取的块分辨率大小均为109×91,先将其分辨率调整到为64×64的灰度图像,由于是用unit8存储,因而灰度为256级,最高灰度级对应光子数为4.0×102s-1。在假设不知道原图像序列的情况下,采用伪随机二值矩阵进行压缩采样,SpaRSA稀疏重建算法进行连续二维图像重建,经过反复调试,发现每个时间间隔内的采样次数越多其重建效果越好,这里全部取k为3000次(大大小于传统奈奎斯特采样频率),而且采样次数越多所需的迭代次数越少,但k越大测量矩阵的维数也越大,计算量也越大,因而综合以上因素不断调试后将迭代次数设置为50,先一次性读入12组连续的图像,每隔0.01秒读入一幅图,然后批量重建,得到图3(a)、3(b)、3(c)、3(d)、4(a)、4(b)、4(c)、4(d)、5(a)、5(b)、5(c)、5(d)、6(a)、6(b)、6(c)、6(d)、7(a)、7(b)、7(c)、7(d)、所示的结果。图3(a)(b)(c)(d)、图4(a)(b)(c)(d)、图5(a)(b)(c)(d)、图6(a)(b)(c)(d)和图7(a)(b)(c)(d)是针对各t时间间隔内的模拟,其中,图3(a)、图4(a)、图5(a)、图6(a)和图7(a)是原始光子密度图像;图3(b)、图4(b)、图5(b)、图6(b)和图7(b)是SpaRSA算法的重建图像;图3(c)、图4(c)、图5(c)、图6(c)和图7(c)是SpaRSA算法的残差图像,图3(d)、图4(d)、图5(d)、图6(d)和图7(d)是psnr峰值信噪比上升曲线。psnr峰值信噪比上升曲线的横坐标为SpaRSA算法迭代次数,纵坐标为psnr的dB数。
图3至图7所示的图3(b)、图4(b)、图5(b)、图6(b)和图7(b)重建图像分别与图3(a)、图4(a)、图5(a)、图6(a)和图7(a)原始图像的相关系数Cov和psnr峰值信噪比如表1所示:
表1
Cov | Psnr | |
图3 | 0.9959 | 27.9733 |
图4 | 0.9961 | 28.124 |
图5 | 0.9956 | 27.6915 |
图6 | 0.9954 | 27.3899 |
图7 | 0.9961 | 28.0811 |
其中,psnr单位都是dB。PSNR值越大,就代表失真越少。从相关系数Cov和psnr峰值可以看出,图像重建效果很好,生成了与原图像相对应的可时间分辨的二维重建图像帧序列,完全满足对原随时间连续变化对象的观测需求,符合本发明的意图。实验结果表明,本发明硬件可以实现压缩采样,软件算法可以实现图像重建。
另外需要说明的是,本发明另一改进之处在于,实验模拟的噪声环境是用matlab模拟的高斯白噪声,所谓的高斯白噪声中的高斯是指概率分布是正态函数,它的幅度分布服从高斯分布,而它的功率密度又是均匀分布的,这里的高斯白噪声服从N(0,1)。实验中为了更好地模拟真实环境中的噪声,特意做了加噪测试,测试中用的函数是noise=randn(k,1),y=R*x+noise,其中k为测量次数,现在对noise加倍,实验结果如表2所示:
表2
Cov | Psnr | |
y=R*x+noise | 0.9988 | 30.7256 |
y=R*x+10*noise | 0.997 | 26.6837 |
y=R*x+50*noise | 0.9643 | 16.3655 |
y=R*x+70*noise | 0.9311 | 13.6902 |
y=R*x+100*noise | 0.8836 | 11.4918 |
y=R*x+150*noise | 0.7787 | 8.6116 |
y=R*x+200*noise | 0.6863 | 6.9938 |
y=R*x+250*noise | 0.6362 | 6.1363 |
y=R*x+300*noise | 0.5837 | 5.4176 |
其中,psnr单位都是dB,noise所乘上的倍数越大,噪声影响越大,重建效果越差。由于本发明实际输入信号直接是测得的单光子计数值,而没有原图,因而无法用由重建图像计算所得的信噪比来作为评价标准,为了更真实的模拟,在这里采用输入信号的信噪比作为评价标准,又作了一次加噪测试,结果如表3所示:
表3
其中,信噪比单位都是dB,可见加性噪声对重建图像效果有直接影响,但若信号加强,噪声的影响便会减弱。此外,也测试了在x(即原图对象)处加随机噪声,再乘上随机测量矩阵,会削弱一部分噪声影响,同样噪声很大时,重建效果也会变差。因而提高信噪比,对成像质量有很大帮助。
最后需要说明的是,具体实施方式中所述的实验用图仅用来说明本发明的技术方案软件算法的可行性而非局限于此例,算法已经经过大量实验数据验证,是真实可靠的,搭配硬件的功能便可实现本发明的技术方案。尽管参照实施例对本发明进行了详细说明,本领域的普通技术人员应当理解,对本发明的技术方案进行修改或者等同替换,都不脱离本发明技术方案的精神和范围,其均应涵盖在本发明的权利要求范围当中。
Claims (9)
1.一种时间分辨单光子计数二维成像系统,该系统主要基于压缩传感理论,用于对随时间动态变化的物体成像,输出按时间序列排列的连续灰度视频图像帧,所述系统包含:触发器、光学成像系统、DMD、光学聚焦收集系统、光衰减器、单光子计数器、驱动控制模块和最优化算法模块;
所述触发器由位于其前端的极弱光触发源触发,该触发器输出端与驱动控制模块的输入端相连,当触发器被触发时驱动控制模块将输出驱动控制信号触发与驱动控制模块输出端相连的DMD和单光子计数器开始工作,则该DMD开始翻转,单光子计数器同时开始计数,光衰减器用于对光强度进行衰减,光学聚焦收集系统用于对光线的聚焦收集,单光子计数器是对光进行单光子计数;
所述单光子计数器的输出端与所述最优化算法模块的一个输入端相连,将单光子计数器的计数值作为最优化算法模块的一个输入参数,所述最优化算法模块的另一输入端与所述驱动控制模块的一输出端相连,用于接收该驱动控制模块存储的选定工作区域的随机测量矩阵作为其另一输入参数;所述最优化算法模块根据输入的单光子计数器的测量值与所述随机测量矩阵重建出稀疏信号,并反演出光子密度图像,经过M个t时间间隔,就能重建出按时间序列排列的可时间分辨的一系列连续灰度视频图像帧;
所述DMD为数字微镜器件;
其中,所述驱动控制模块的另一输出端与所述DMD的输入端相连,用于驱动控制DMD的翻转;该驱动控制模块基于DLP数字光处理技术,该驱动控制模块在DMD选定工作区域后,下载随机测量矩阵驱动控制所述DMD的翻转;所述DMD在翻转的同时会向所述单光子计数器发送同步信号,保证DMD与单光子计数器之间同步,即DMD每翻转一次,单光子计数器累积计数在该次翻转的时间间隔内的单光子数,DMD翻转完成后,单光子计数器清零重新开始累积计数,所有计数和该选定工作区域的随机测量矩阵均传送至最优化算法模块中。
2.根据权利要求1所述的时间分辨单光子计数二维成像系统,其特征在于,所述极弱光触发源为:极弱光光源或自发光生物体。
3.根据权利要求2所述的时间分辨单光子计数二维成像系统,其特征在于,当所述极弱光触发源采用极弱光光源时,其投射方式有两种:1)极弱光光源直接透射至物体上,光源与光学成像系统光轴在同一条直线上,物体对象是半透明的或者镂空的,通过透射光成像;2)极弱光光源斜着打在物体上,光源与光学成像系统光轴不在同一条直线上,通过物体对象反射极弱光来成像;
当所述极弱光光源比较强时采用滤光片滤除该极弱光的杂光,若极弱光光源的光强极其微弱,且其波长在探测所要求的波长范围内,则不需要再设置滤光片。
4.根据权利要求1所述的时间分辨单光子计数二维成像系统,其特征在于,所述光学成像系统和光学聚焦收集系统均采用光学透镜组,分别负责光学成像和光学聚焦,极弱光通过光学成像系统后,可在DMD上成等大或缩小或放大的像,按实际需求进行成像调整;
其中,所述光学聚焦收集系统采用光纤耦合技术,即将经由分光光度计分光后的光束耦合到光纤中,利用光纤耦合技术将分光分别收集到对应的单光子计数器上。
5.根据权利要求1所述的时间分辨单光子计数二维成像系统,其特征在于,所述的光衰减器是由多块衰减片组成,放置在光学聚焦收集系统和单光子计数器之间的光路上或光源处,用于将光衰减到单光子计数器的工作范围;该光衰减器的设计是为了防止被测光子密度过大和单光子计数器的门控时间过长引起的饱和,若被测光已经达到单光子级别,则不需要再设置光衰减器了;若将光衰减片放置在光源处,这可将光源变成单光子源。
6.一种时间分辨单光子计数二维成像方法,该方法基于权利要求1-5任意一条权利要求记载的系统实现对随时间动态变化的物体成像,进而输出按时间序列排列的连续灰度视频图像帧,所述方法包含如下步骤:
步骤1,用于采用极弱光触发源触发触发器进行触发同步的步骤,该步骤实现了时间分辨的效果;
步骤2,用于将被测信号由高维向低维进行映射的步骤,该步骤采用压缩和采样对被测信号进行由高维向低维进行映射;
步骤3,用于对待输出的视频帧进行稀疏重建最优化,输出按时间序列排列的连续灰度视频图像帧的步骤;
所述步骤1进一步包含如下子步骤:
用于触发的步骤,极弱光光源或自发光生物体触发触发器,该触发器进而触发激活整个成像系统的工作,所述触发器触发之后每隔t时间间隔进行一组单光子探测,最优化算法模块自动识别该t时间间隔内的极弱光二维图像的输入,分别累积测量该t时间间隔内的单光子数;
用于同步的步骤,DMD在翻转的同时向单光子计数器发送同步信号,即DMD每翻转一次,单光子计数器累积计数在该次翻转的时间间隔内的单光子数,DMD翻转完成后,单光子计数器清零重新开始累积计数,所有的计数被传至最优化算法模块中;
其中,所述的时间分辨能够针对连续变化观测物体和周期变化观测物体,所述连续变化观测物体的时间分辨采用每t时间间隔采样k次并持续测量M×t时间的方法能实现;所述周期变化观测物体的时间分辨采用的观测物体的变化周期极短,假设周期为T,将时间周期等分为d个时间间隔,记做t1,t2,t3,…,td,在该周期T内保持对应随机测量矩阵不变,到下一个周期随机测量矩阵改变,在每个小的时间间隔t1,t2,t3,…,td内分别计数,即对落在ti时间间隔内的单光子进行累积计数,依靠触发器的触发,保证计数和随机测量矩阵的严格对应关系,测量k个周期,即对每个微小时间间隔测量了k次,分别对这d个时间间隔做最优化重建,便可反演出一个时间周期序列内的灰度观测对象的变化情况;所述M表示测量次数。
7.根据权利要求6所述的时间分辨单光子计数二维成像方法,其特征在于,所述步骤2进一步包含如下子步骤:
用于压缩的步骤,DMD将可压缩的随时间动态变化的灰度物体的图像帧数据随机反射,当DMD中单个微镜+12°翻转时反射的光被单光子计数器接收,当DMD中单个微镜-12°翻转时反射光不能被单光子计数器接收,从而完成对随时间动态变化的灰度物体的被测信号的压缩,同时确保DMD明暗阵列的最大随机性,进而控制极弱光被反射至光学聚焦收集系统的概率是随机的;
用于分光采样的步骤,经光学聚焦收集系统聚焦后的极弱光进入光衰减器衰减后再进入单光子计数器,单光子计数器再对极弱光进行探测采样。
8.根据权利要求6所述的时间分辨单光子计数二维成像方法,其特征在于,所述步骤3进一步包含如下子步骤:
单光子计数器每t时间间隔内进行计数,将该计数值作为单光子计数器的测量值输入最优化算法模块;
最优化算法模块根据上步测量值、DMD上的随机测量矩阵和加载在二维图像上的稀疏矩阵,分别通过最优化算法重建出光子密度图像,反演出时间间隔t内随时间动态变化的物体的灰度图;
重复执行上述两个子步骤共M次,得到M×t时间段的M幅连续灰度视频图像帧序列,输出视频帧。
9.根据权利要求6所述的时间分辨单光子计数二维成像方法,其特征在于,所述最优化算法为运用小波进行稀疏变换的可分离逼近的稀疏重建算法SpaRSA-DWT,该算法基于压缩传感理论,利用小波变换对单光子计数测量值进行稀疏化,并调整IST算法每次迭代的步长系数αt,使αtI逼近f(x)在xt处的Hessian矩阵,并自适应地修改阈值,依靠反复迭代运算求解出对应的稀疏信号,最后反演出二维灰度图像;
经过M×t时间后,共读入M幅极弱光对象二维变化图像,对应输出M幅随时间变化的连续灰度视频图像帧;
其中,所述M表示测量次数。
Priority Applications (6)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201110328462.2A CN102510282B (zh) | 2011-10-25 | 2011-10-25 | 一种时间分辨单光子计数二维成像系统及方法 |
EP12843969.2A EP2755327A4 (en) | 2011-10-25 | 2012-05-14 | MULTI-DIMENSIONAL MULTI-DIMENSIONAL IMAGING SYSTEM AND METHOD WITH TEMPORARY RESOLUTION OR ULTRA-LIGHT LIGHT |
US14/351,028 US9448162B2 (en) | 2011-10-25 | 2012-05-14 | Time-resolved single-photon or ultra-weak light multi-dimensional imaging spectrum system and method |
CN201280048647.0A CN104054266B (zh) | 2011-10-25 | 2012-05-14 | 一种时间分辨单光子或极弱光多维成像光谱系统及方法 |
JP2014537455A JP6002232B2 (ja) | 2011-10-25 | 2012-05-14 | 時間分解単一光子計数イメージングスペクトルシステム |
PCT/CN2012/075444 WO2013060134A1 (zh) | 2011-10-25 | 2012-05-14 | 一种时间分辨单光子或极弱光多维成像光谱系统及方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201110328462.2A CN102510282B (zh) | 2011-10-25 | 2011-10-25 | 一种时间分辨单光子计数二维成像系统及方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN102510282A CN102510282A (zh) | 2012-06-20 |
CN102510282B true CN102510282B (zh) | 2014-07-09 |
Family
ID=46222340
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201110328462.2A Expired - Fee Related CN102510282B (zh) | 2011-10-25 | 2011-10-25 | 一种时间分辨单光子计数二维成像系统及方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN102510282B (zh) |
Families Citing this family (30)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104054266B (zh) * | 2011-10-25 | 2016-11-23 | 中国科学院空间科学与应用研究中心 | 一种时间分辨单光子或极弱光多维成像光谱系统及方法 |
CN102769460B (zh) * | 2012-07-27 | 2015-05-06 | 中国科学院空间科学与应用研究中心 | 一种互补测量的时间分辨单光子计数成像系统及方法 |
CN103115681B (zh) * | 2013-01-24 | 2015-02-04 | 中国科学院空间科学与应用研究中心 | 一种超灵敏时间分辨光谱仪及其时间分辨方法 |
CN103115680B (zh) * | 2013-01-24 | 2014-11-12 | 中国科学院空间科学与应用研究中心 | 一种超灵敏光谱仪以及光谱检测方法 |
ES2696748T3 (es) * | 2013-12-09 | 2019-01-17 | Hamamatsu Photonics Kk | Elemento de conteo de fotones bidimensional |
CN104019898B (zh) * | 2014-05-28 | 2017-01-11 | 中国科学院空间科学与应用研究中心 | 一种超灵敏光谱成像天文望远镜及天文光谱成像方法 |
CN104702830A (zh) * | 2015-03-13 | 2015-06-10 | 四川大学 | 一种逐点扫描数字微镜阵列相机 |
CN104796674A (zh) * | 2015-04-17 | 2015-07-22 | 南京理工大学 | 基于压缩感知的彩色成像装置及方法 |
CN104914446B (zh) * | 2015-06-19 | 2017-06-27 | 南京理工大学 | 基于光子计数的三维距离图像时域实时去噪方法 |
DE102015112628A1 (de) * | 2015-07-31 | 2017-02-02 | Carl Zeiss Microscopy Gmbh | Verfahren zur Erzeugung eines digitalen Fluoreszenzbildes |
CN105606228B (zh) * | 2016-02-04 | 2018-07-17 | 北京理工大学 | 基于编码变换的双波长温度场成像设备、系统及方法 |
CN105675146B (zh) * | 2016-02-04 | 2018-07-17 | 北京理工大学 | 基于压缩感知的双波长三维温度场成像设备、系统及方法 |
CN105737992B (zh) * | 2016-02-04 | 2018-07-17 | 北京理工大学 | 基于压缩感知的双波长温度场成像设备、系统及方法 |
CN105915869A (zh) * | 2016-04-22 | 2016-08-31 | 南京理工大学 | 一种彩色自适应压缩计算鬼成像系统及方法 |
CN105915868A (zh) * | 2016-04-22 | 2016-08-31 | 南京理工大学 | 基于扩展小波树的彩色成像系统及方法 |
CN106770131A (zh) * | 2017-01-17 | 2017-05-31 | 清华大学 | 三维超光谱显微成像系统及成像方法 |
CN107121830A (zh) * | 2017-06-21 | 2017-09-01 | 山西大学 | 一种使用白光实现彩色显示的装置及方法 |
CN107729990B (zh) | 2017-07-20 | 2021-06-08 | 上海寒武纪信息科技有限公司 | 支持离散数据表示的用于执行正向运算的装置及方法 |
EP3467598B1 (en) * | 2017-10-04 | 2021-09-29 | TTTech Computertechnik AG | Method and apparatus for the determination of the slot-duration in a time-triggered control system |
CN108469673B (zh) * | 2018-01-16 | 2019-10-29 | 南昌大学 | 纠缠光子对时间和位置同步符合的量子成像装置及方法 |
CN108549275B (zh) * | 2018-03-01 | 2021-01-19 | 南昌大学 | 一种单光子压缩成像的控制装置及控制方法 |
CN109029727B (zh) * | 2018-06-11 | 2020-06-02 | 中国科学院紫金山天文台 | 基于编码孔径的高灵敏度太赫兹超导光谱成像系统及系统成像方法 |
CN109213037B (zh) * | 2018-08-28 | 2021-03-26 | 南昌大学 | 采样时间自适应单光子压缩成像控制方法及控制装置 |
CN109141654B (zh) * | 2018-09-21 | 2020-05-05 | 上海交通大学 | 基于单光子成像器件的光子级空间映射关联性测量方法 |
CN109361833B (zh) * | 2018-10-08 | 2020-08-11 | 南昌大学 | 一种单光子压缩视频传输装置的传输方法 |
CN109613556B (zh) * | 2018-11-26 | 2021-05-18 | 武汉大学 | 基于稀疏表征的光子计数激光三维探测成像方法 |
CN112484857B (zh) * | 2020-11-04 | 2023-04-07 | 西北工业大学宁波研究院 | 一种基于dmd的光谱成像系统及方法 |
CN113885023B (zh) * | 2021-08-18 | 2024-04-30 | 西安电子科技大学 | 基于lcim算法的三维模型微波光子回波成像方法 |
CN113674248B (zh) * | 2021-08-23 | 2022-08-12 | 广州市番禺区中心医院(广州市番禺区人民医院、广州市番禺区心血管疾病研究所) | 磁共振酰胺质子转移成像磁化转移率检测方法及相关设备 |
CN115479928B (zh) * | 2022-09-07 | 2024-06-07 | 山东大学 | 基于labview的低强度光信号的同步探测控制系统及控制方法 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060054778A1 (en) * | 2003-02-07 | 2006-03-16 | Klaus Suhling | Photon arrival time detection |
CN101387548A (zh) * | 2007-09-11 | 2009-03-18 | 中国科学院西安光学精密机械研究所 | 单光子计数成像仪 |
CN101881658A (zh) * | 2009-05-07 | 2010-11-10 | 中国科学院西安光学精密机械研究所 | 一种高分辨率位敏阳极探测器及其阳极解码方法 |
-
2011
- 2011-10-25 CN CN201110328462.2A patent/CN102510282B/zh not_active Expired - Fee Related
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060054778A1 (en) * | 2003-02-07 | 2006-03-16 | Klaus Suhling | Photon arrival time detection |
CN101387548A (zh) * | 2007-09-11 | 2009-03-18 | 中国科学院西安光学精密机械研究所 | 单光子计数成像仪 |
CN101881658A (zh) * | 2009-05-07 | 2010-11-10 | 中国科学院西安光学精密机械研究所 | 一种高分辨率位敏阳极探测器及其阳极解码方法 |
Non-Patent Citations (2)
Title |
---|
Single-Pixel Imaging Via Compressive Sampling;Marco F.Duarte等;《IEEE SIGNAL PROCESSING MAGZINE》;20080331;第84页左栏第1段至左栏第2段 * |
基于压缩传感的光子计数成像系统研究;杜克铭等;《二〇一〇国防空天信息技术前沿论坛论文集》;20101023;第2页左栏第1段至第4页右栏第1段 * |
Also Published As
Publication number | Publication date |
---|---|
CN102510282A (zh) | 2012-06-20 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN102510282B (zh) | 一种时间分辨单光子计数二维成像系统及方法 | |
CN102393248B (zh) | 一种时间分辨极弱光多光谱成像系统及方法 | |
Faccio et al. | Non-line-of-sight imaging | |
Tachella et al. | Real-time 3D reconstruction from single-photon lidar data using plug-and-play point cloud denoisers | |
US9727959B2 (en) | System and processor implemented method for improved image quality and generating an image of a target illuminated by quantum particles | |
US9131128B2 (en) | System and processor implemented method for improved image quality and generating an image of a target illuminated by quantum particles | |
Chen et al. | Steady-state non-line-of-sight imaging | |
Peng et al. | Photon-efficient 3d imaging with a non-local neural network | |
Gupta et al. | Reconstruction of hidden 3D shapes using diffuse reflections | |
Radwell et al. | Deep learning optimized single-pixel LiDAR | |
US9064315B2 (en) | System and processor implemented method for improved image quality and enhancement | |
CN102769460B (zh) | 一种互补测量的时间分辨单光子计数成像系统及方法 | |
CN102768069B (zh) | 一种互补测量的单光子光谱计数成像系统及方法 | |
CN102901564B (zh) | 一种互补测量的时间分辨单光子光谱计数成像系统及方法 | |
Edgar et al. | Real-time 3D video utilizing a compressed sensing time-of-flight single-pixel camera | |
CN109613556B (zh) | 基于稀疏表征的光子计数激光三维探测成像方法 | |
Edgar et al. | Real-time computational photon-counting LiDAR | |
US20230125131A1 (en) | Ultrafast light field tomography | |
CN102564614A (zh) | 激光光斑动态测量方法及测量仪 | |
CN102768070B (zh) | 一种互补测量的单光子计数成像系统及方法 | |
CN105737981B (zh) | 一种双向光谱光传输采集探测方法 | |
Chen et al. | A Compact Upconversion Single-Photon Imager for Full-Range and Accurate 3-D Imaging | |
CN113099209A (zh) | 基于光电倍增管阵列的非视域成像装置和方法 | |
Malik et al. | Flying With Photons: Rendering Novel Views of Propagating Light | |
CN105380638A (zh) | 一种用于激光散斑血流速度的定量成像装置及其方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
CP01 | Change in the name or title of a patent holder |
Address after: 100190 No. two south of Zhongguancun, Haidian District, Beijing 1 Patentee after: NATIONAL SPACE SCIENCE CENTER, CAS Address before: 100190 No. two south of Zhongguancun, Haidian District, Beijing 1 Patentee before: NATIONAL SPACE SCIENCE CENTER, CHINESE ACADEMY OF SCIENCES |
|
CP01 | Change in the name or title of a patent holder | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20140709 |
|
CF01 | Termination of patent right due to non-payment of annual fee |