紫外域多纵模高光谱分辨率激光雷达探测气溶胶的技术实现和系统仿真∗
2018-03-26高飞南恒帅黄波汪丽李仕春王玉峰刘晶晶闫庆宋跃辉华灯鑫
高飞 南恒帅 黄波 汪丽 李仕春 王玉峰 刘晶晶闫庆 宋跃辉 华灯鑫
(西安理工大学机械与精密仪器工程学院,西安 710048)
1 引 言
气溶胶是大气污染物的重要组成部分,成分十分复杂,物理化学特性差异显著[1,2],是全球气候变化研究中最不确定的影响因子之一[3].气溶胶随气流漂浮或沉降,在大气辐射的吸收和散射、云雾降水的形成、大气污染以及大气光电现象的产生过程中都起着非常重要的作用[4].测量气溶胶的物理化学特性及其演化规律,对于揭示气溶胶时空分布与近地表环境的作用关系,分析区域环境和气候变化具有非常重要的作用.
高光谱分辨率激光雷达是激光雷达的一种重要形式,在大气污染及环境监测中扮演着重要的角色,其可精确提供表征气溶胶特性的物理参量––气溶胶消光和后向散射系数[5].其探测原理主要是利用超窄带滤波器分离单纵模(单频)脉冲激光的气溶胶米散射和大气分子瑞利散射光谱信号,从而实现在无需额外假设的条件下反演气溶胶消光、后向散射和激光雷达比等光学参量[6].单纵模高光谱分辨率激光雷达技术实现的关键主要有以下两点.1)高精度的激光锁频技术,以提供单纵模脉冲激光.Zhao等[7]设计的锁频系统实现了激光80 MHz频率稳定精度.Cheng等[8]利用多谐波外差技术实现了宽视场角迈克耳孙干涉仪的锁频技术.2)超窄带滤波器的设计与应用,以实现对米散射或者瑞利散射信号高达三个数量级的抑制率.目前主要有两种形式的滤波器成功应用在高光谱分辨率激光雷达上:原子/分子吸收滤波器和法布里-珀罗标准具滤波器.钡原子吸收滤波器和碘分子吸收滤波器是典型的原子/分子吸收滤波器,其中钡原子吸收滤波器利用特定温度下原子吸收谱线与激光大气回波的米散射光谱重合的原理,分离气溶胶米散射和大气分子瑞利散射信号[9].虽然钡原子吸收滤波器可以提供气溶胶米散射信号的高抑制比,但是激光器必须选择染料激光器,限制了其应用范围.碘分子吸收滤波器在普通室温条件下就可提供精细的超窄带吸收光谱,实现对米散射信号高达三个数量级的抑制率[10,11].以碘分子吸收滤波器为光学分光器的高光谱分辨率激光雷达可选择Nd:YAG激光器的二倍频输出532 nm作为激励光源,因此极大地推动了高光谱分辨率激光雷达在气溶胶探测中的应用[12,13].随着光学加工工艺的提高,法布里-珀罗标准具通过提供高精细的窄带干涉条纹谱线,成功应用于高光谱分辨率激光雷达中[14].通过调节法布里-珀罗标准具窄带反射条纹谱线的中心波长,使之与激光发射波长重合,从而利用窄带干涉条纹抑制大气分子瑞利散射信号而透过气溶胶米散射信号,实现气溶胶米散射和大气分子瑞利散射信号的分离[15,16].Hua等[17]在温度测量高光谱分辨率激光雷达中设计了多通道的法布里-珀罗标准具分离光谱信号,而气溶胶作为其副产品完全可以实现精确测量.
考虑到单纵模高光谱分辨率激光雷达系统建设的复杂性,日本国立环境研究所的Jin等[18]和阿根廷División Lidar的Ristori等[19]提出了多纵模高光谱分辨率激光雷达的思想.该思想摒弃了高光谱分辨率激光雷达必须选择单纵模脉冲激光的思路,利用可调谐的干涉仪来分离多纵模激光的气溶胶米散射和大气分子瑞利散射信号.浙江大学的刘东团队借助这种思想,利用设计的宽视场迈克耳孙干涉仪开展近红外波段多纵模高光谱分辨率激光雷达的气溶胶探测技术和理论研究[20].本文从激光器的多纵模模式出发,分析紫外域激光多纵模模式在大气中传输的气溶胶米散射光谱和瑞利散射光谱,探讨紫外域多纵模高光谱分辨率激光雷达建设的可能性,并利用大气模型对所设计的激光雷达系统进行系统仿真.
2 激光器的多纵模模式及其在大气传输中的散射光谱
Nd:YAG固体激光器可提供基频1064 nm红外光及其二倍频532 nm绿光、三倍频355 nm紫外光,已经成为激光雷达系统激励光源的的主流选择.考虑到对太阳背景光的抑制以及大气分子散射强度与波长的四次方成反比的关系,本文选择Nd:YAG脉冲激光器的三倍频输出355 nm作为研究对象,分析其多纵模模式及其在大气传输中的光谱信号.图1给出了为太阳背景光在不同波长下的强度分布,可以看出,紫外域波长355 nm相比于1064 nm和532 nm,太阳背景光的强度要弱些,使得以其作为激励波长的激光雷达更容易实现全天时探测.
图1 太阳背景光在不同波长下的强度分布Fig.1.Solar spectral irradiance at different wavelengths.
Nd:YAG固体激光器谐振腔的谐振频率νq为
式中,c表示光速,λ表示激光波长,n表示激光谐振腔的折射率,L表示谐振腔长度.激光谐振腔具有选频的作用,即从频带很宽的光波中,选出满足谐振频率的光波,相对应的模式称为纵模.相邻两个纵模的频率之差Δνd称为纵模间隔,为
则激光器辐射线宽Δν内的纵模个数N为
以美国Continuum公司生产的Surelite系列的Nd:YAG固体激光器为例,在无种子注入的情况下,其辐射线宽Δν为1cm−1(30 GHz),当激光谐振腔的有效长度为40cm时,纵模间隔为375 MHz,辐射线宽内的纵模个数为80.该系列激光器的纵模频谱分布图如图2所示.
激光束在大气中传输时,会与大气中的物质(大气分子、气溶胶粒子)发生相互作用产生大气散射回波光谱信号,包括由气溶胶粒子产生的米散射、大气分子产生的瑞利散射以及分子振动和转动产生的拉曼散射和荧光等.
图2 典型Nd:YAG固体激光器的多纵模模式及其频谱分布Fig.2.Multi-longitudinal modes and its spectral distribution of typical Nd:YAG laser.
气溶胶米散射与气溶胶的尺度谱、光学折射率以及激光波长等密切相关,散射机制较为复杂,其光谱宽度本质上是气溶胶粒子的布朗运动所造成的多普勒展宽,其半高全宽与激光器的线宽有关,一般为百兆赫兹.大气分子瑞利散射可看作是米散射在散射粒子的直径远远小于入射激光波长情况下的近似,其散射光谱RM(ν)是由大气分子热运动造成的多普勒展宽,其可简化为高斯模型表示:
式中,ν为频率,ν0为激光发射频率,m为分子质量,k为玻尔兹曼常数.大气分子瑞利散射强度分布I(ν)可近似为高斯分布:
式中,δν=7.16×10−7ν0(T/M)1/2为瑞利散射光谱宽度,T为大气温度,M为大气分子的平均摩尔质量(28.96 g/mol),I0表示ν0处的归一化强度.在T为300 K,ν0为355 nm时,光谱宽度约为4 GHz.激光大气后向散射信号的光谱分布形式,就是气溶胶后向散射信号叠加在大气分子后向散射的多普勒展宽信号上(图3).
而多纵模的脉冲激光在大气中传输时,相当于多个有序单纵模所产生的气溶胶米散射和大气分子瑞利散射信号的叠加,也就是说,多纵模激光的气溶胶米散射和大气分子瑞利散射光谱为多纵模脉冲激光的频谱与单纵模脉冲激光的气溶胶米散射和瑞利散射光谱的卷积.典型Nd:YAG多纵模脉冲激光器紫外光输出的气溶胶米散射和大气分子瑞利散射光谱如图3所示.
图3 典型Nd:YAG多纵模激光器紫外光输出的气溶胶米散射和大气分子瑞利散射光谱,右上角图为单纵模激光的米散射和瑞利散射光谱Fig.3.Spectra of aerosol Mie scattering and molecular Rayleigh scattering excited by typical Nd:YAG UV multi-mode lasers.The up-right layer displays the spectral of Mie and Rayleigh scattering by a single frequency laser.
3 紫外域多纵模高光谱分辨率激光雷达系统
为了实现对紫外域多纵模激光的气溶胶米散射信号和大气分子瑞利散射信号的接收和分离,设计的紫外域多纵模高光谱分辨率激光雷达系统结构图如图4所示.激光发射子系统选用高功率Nd:YAG脉冲激光器的三倍频输出355 nm作为激励波长,其辐射线宽为30 GHz.在谐振腔腔长为40cm的情况下,在辐射线宽内有80个纵模输出.多纵模激光经准直扩束后进入大气,所产生的后向散射光被卡塞格林望远镜接收,经多模光纤传输进入光学分光子系统.同时,发射激光在射向大气前,部分光被光束分离镜所截取并耦合在多模光纤以传输至接收望远镜中,这部分光束作为参考信号用来评价激光的光谱特性,并和大气回波信号一起被高速数据采集卡(采样频率>100 MHz)进行信号采集.
光学分光系统由窄带干涉滤光片和可调谐马赫-曾德尔干涉仪组成.窄带干涉滤光片主要用于抑制太阳背景光,以实现气溶胶时空变化的全天时观测.其中心波长为354.7 nm,带宽应综合考虑脉冲激光器的辐射线宽和太阳背景光的影响,系统选择滤光片的带宽为1 nm.可调谐的马赫-曾德尔干涉仪用于分离多纵模激光的气溶胶米散射光谱和大气分子瑞利散射光谱,其主要由一个安置在压电陶瓷上的屋脊反射镜和两个光束分离镜构成.多模光纤的输出经准直透镜后变成平行光,经窄带干涉滤光片IF后被光束分离镜BS1分成能量相等的两束光.其中反射光束射向可移动的屋脊反射镜并经折返后,与光束分离镜BS1的透射光束经光束分离镜BS2后相遇而发生干涉现象.光束分离镜BS2的输出为两束相位相反的叠加光束,并通过光电倍增管(PMT)进行光电转换.
图4 紫外域多纵模高光谱分辨率激光雷达系统结构图(IF,窄带干涉滤光片;BS,光束分离镜;PMT,光电倍增管)Fig.4.Schematic of ultraviolet multi-mode high-spectral-resolution lidar.IF,interference filter;BS,beam splitter;PMT,photomultiplier tube.
马赫-曾德尔干涉仪的透过率函数TMZ可表示为[21]
式中,ΔνMZ(=c/(2n×OPD))是与光程差相关的光谱宽度.为了实现对多纵模激光的气溶胶米散射和大气分子瑞利散射的精细光谱分离,需要设置马赫-曾德尔干涉仪的光程差,使其自由光谱范围与激光纵模间隔相一致.马赫-曾德尔干涉仪光程差的调节需要控制屋脊反射镜上的压电陶瓷,进而移动屋脊镜的位置.光程差的改变会导致相长干涉与相消干涉的变化.当屋脊反射镜的移动距离为Δd时,光程差的改变量为2Δd,因此在紫外域激光355 nm的测量条件下,每移动177.5 nm相长干涉与相消干涉就会发生一次周期变化.若控制屋脊反射镜的移动距离在纳米级运动,即可实现最佳光程差的调控.
对于具有40cm长度谐振腔的高能量Nd:YAG脉冲激光器来说,其在1cm−1(30 GHz)的辐射线宽内具有80个纵模.由于马赫-曾德尔干涉仪的输出是互补型的干涉圆环或者干涉条纹,因此当一路输出为相长干涉时,另外一路输出为相消干涉.通过调谐马赫-曾德尔干涉仪,可以确定分离米散射和瑞利散射的最佳光程差.在调谐过程中,通过包络分析参考信号(发射激光)来确定最优光程差.当参考信号的一路输出为最小,另一路输出为最大时,则每一个激光纵模的中心波长就被锁定在最优光程差处.在此激光雷达中,光电检测系统所测量的数据(两路PMT输出)分别为激光雷达回波中的多纵模脉冲激光的气溶胶米散射和大气分子瑞利散射信号与马赫-曾德尔干涉仪透过率函数的卷积.因此在测量激光大气回波信号时,其中一路输出就表现为气溶胶米散射信号的透过,称为米通道;而另外一路输出则表现为对气溶胶米散射信号的抑制,称为瑞利通道.图5为基于马赫-曾德尔干涉仪的多纵模激光米散射和瑞利散射信号的分离提取原理,其中图5(a)为瑞利通道透过率函数,图5(b)为米通道透过率函数.由于马赫-曾德尔干涉仪双通道输出的互补性,因此米通道和瑞利通道信号之和为总的激光雷达回波信号.
图5 基于马赫-曾德尔干涉仪的多纵模激光米散射和瑞利散射信号的分离提取技术(a)瑞利通道透过率函数;(b)米通道透过率函数Fig.5.Separation and extraction of multi-mode laser Mie scattering and Rayleigh scattering signals based on Mach-Zehnder interferometer:(a)Transmittance function of Rayleigh channel;(b)transmittance function of Mie channel.
4 气溶胶光学参量的数据反演
弹性散射激光雷达方程在分别考虑到气溶胶和大气分子在大气后向散射β(r)和α(r)消光中的贡献时,可表示为
式中,P(r)表示激光雷达的瞬时接收功率,C为激光雷达系统常数,P0为激光峰值发射功率,βa(r)和βm(r)分别表示气溶胶和大气分子的后向散射系数,αa(r)和αm(r)分别表示气溶胶和大气分子的消光系数.考虑到马赫-曾德尔干涉仪双通道输出的互补性,米通道和瑞利通道的激光雷达方程可以分别表示为
式中,Taa表示米通道的米散射信号透过率,Tam表示米通道的瑞利散射信号透过率,Tma表示瑞利通道的米散射信号透过率,Tmm表示瑞利通道的瑞利散射信号透过率.由于两个通道是互补输出的,根据能量守恒定律,在不考虑马赫-曾德尔干涉仪光学元件对光的吸收以及光能量损耗的情况下,
同时,考虑到大气分子瑞利散射谱宽远超过纵模间隔,因此多纵模激光脉冲的瑞利散射谱近似为一条30 GHz的方波谱(图3),因此在不考虑马赫-曾德尔干涉仪光学元器件的吸收以及能量损耗的情况下,则米散射通道的瑞利散射透过率(Tam)与瑞利散射通道的瑞利散射透过率(Tmm)相等,即Tam=Tmm=1/2.
若(9)与(10)式相加,则可得到激光雷达总的大气回波信号,即
若(9)与(10)式相减,则方程可表示为
即激光雷达方程中后向散射系数项中仅留有气溶胶后向散射项,这与拉曼激光雷达方程的形式非常相似,而拉曼激光雷达的后向散射中仅有大气分子的贡献.气溶胶后向散射与气溶胶粒子的粒子谱、折射率和形状特性有关,因此还是很难求解得到气溶胶消光系数.
若(9)与(10)式相比,则方程可表示为
在气溶胶光学参量研究中,气溶胶后向散射比Ra定义为
因此,(14)式可以改写为
米散射通道的米散射信号透过率Taa或者瑞利散射通道的米散射信号透过率Tma可通过理论计算以及实验的方法测定得到,因此通过马赫-曾德尔干涉仪互补双通道信号的比值可以准确反演气溶胶后向散射比,进而利用(13)和(15)式反演得到气溶胶后向散射系数和气溶胶消光系数.
5 系统仿真
为了验证紫外域多纵模高光谱分辨率激光雷达的探测性能,本文将从气溶胶光学参量数据反演的逆向出发,在假设气溶胶后向散射比的情况下,对多纵模高光谱分辨率激光雷达系统进行仿真,分别计算白天条件下米散射通道和瑞利散射通道的信号强度以及信噪比.激光雷达系统仿真参数如表1所列.系统仿真采用激光雷达实测数据的大气状态模型,如图6所示.其中图6(a)分别表示大气分子和气溶胶的后向散射系数,图6(b)表示气溶胶后向散射比,在5–6 km的高度上,存在一个明显的气溶胶/云层,使得气溶胶的后向散射系数和后向散射比有一个明显的跃变过程.
表1 紫外域多纵模高光谱分辨率激光雷达系统参数Table 1.Speci fications of the UV multi-mode high-spectral-resolution lidar system.
图6 系统仿真所用大气模型(a)气溶胶和大气分子后向散射系数;(b)气溶胶后向散射比Fig.6.The atmospheric model for system simulation study:(a)Aerosol and molecular backscatter coeff-cients;(b)aerosol backscatter ratio.
太阳背景光是全天时激光雷达白天探测的绝对干扰,由于望远镜视场角的存在,激光雷达所接收到的太阳背景光Pb计算公式为
式中,Sb为太阳辐射在地表的能量密度,此处取355 nm处的最大太阳辐射能量密度(Sb=0.3×109W/(m2·sr·nm));Ar为望远镜的有效接收面积;Δλ为窄带干涉滤光片的带宽;θ为望远镜的接收视场角.
根据激光雷达方程以及表1所列的激光雷达系统参数,分别计算得到紫外域多纵模高光谱分辨率激光雷达的米散射通道、瑞利散射通道和太阳背景光的信号强度如图7所示.在气溶胶/云层存在的地方,激光雷达回波信号的强度有明显的跃变.
在综合考虑太阳背景光、PMT本身暗电流等噪音的情况下,分别计算多纵模高光谱分辨率激光雷达的米散射通道和瑞利散射通道的信噪比为
图7 米散射通道和瑞利散射通道以及太阳背景光的回波信号强度分布Fig.7.The intensity distribution of return signals of Mie channel and Rayleigh channel as well as solar background light.
式中Ni米散射通道或者瑞利散射通道所接收到的光子数,n为采样累积次数,Nb表示太阳背景光Pb在PMT上所产生的光电子数,Nd表示PMT本身的暗电流所产生的暗计数.
紫外域多纵模高光谱分辨率激光雷达的米散射和瑞利散射通道的信噪比计算结果如图8所示.激光雷达的米散射通道和瑞利散射通道的信噪比在10 km的探测范围内均大于100,说明系统即使在模拟探测模式下,也可以实现大气边界层乃至整个对流层的气溶胶探测.在气溶胶/云层存在的高度上,回波信号强度增加,信噪比提升.
图8 米散射通道和瑞利散射通道的信噪比随高度的变化Fig.8.Signal-to-noise ratios of Mie channel and Rayleigh channel versus height.
6 讨 论
从原理上来讲,多纵模高光谱分辨率激光雷达利用可调谐的高光谱马赫-曾德尔干涉仪及其双通道输出的互补性,可精细分离提取多纵模激光的气溶胶米散射信号和大气分子瑞利散射信号,实现气溶胶光学特性的精细探测.然而,多纵模高光谱分辨率激光雷达的实现必须要仔细考虑以下两个因素.1)多纵模脉冲激光纵模模式输出的稳定性.只有干涉仪的自由光谱范围与多纵模脉冲激光的纵模模式相互匹配,才能实现多纵模激光的气溶胶米散射和大气分子瑞利散射信号的分离提取.而激光的多纵模模式受激光腔长、工作物质的增益线宽和激励水平等因素所决定,只有处于增益线宽内的纵模频率才能真正起振,形成多纵模振荡.因此,多纵模高光谱分辨率激光雷达的系统实现需要仔细研究多纵模激光的特性,尤其是脉冲激光的多纵模稳定特性.2)可调谐马赫-曾德尔干涉仪是分离提取多纵模激光的米散射和瑞利散射信号的关键.马赫-曾德尔干涉仪的输出稳定性受光束发散角、环境温度、空气湍流及机械振动等因素影响,选择合适的光程差,并合理分析这些因素对马赫-曾德尔干涉仪透过率曲线的影响,进而分析对气溶胶测量精度的影响,提出精确控制这些因素的方法也是多纵模高光谱分辨率激光雷达系统实现的重要保障.同时在分析干涉仪透过率影响的基础上,分析气溶胶测量结果的不确定度,为优化激光雷达系统性能提供可循依据.今后将注重这两方面的研究工作.
7 结 论
多纵模高光谱分辨率激光雷达是高光谱分辨率激光雷达的一个新概念,其摒弃了传统的高光谱分辨率激光雷达必须采用单纵模脉冲激光器的思想,采用普通的高功率Nd:YAG固体激光器作为激励光源,利用高光谱的马赫-曾德尔干涉仪分离提取多纵模激光的气溶胶米散射和大气分子瑞利散射光谱,进而实现气溶胶光学参量的精细探测.
本文以典型Nd:YAG固体脉冲激光器为例,精确计算Nd:YAG激光器的多纵模模式,分析其在大气中传输的气溶胶米散射和大气分子瑞利散射光谱信号.当激光器谐振腔的有效长度为40cm时,纵模间隔为375 MHz,在1cm−1(30 GHz)的辐射线宽内纵模个数为80个.多纵模激光的气溶胶米散射光谱和大气分子瑞利散射光谱是多纵模脉冲激光的频谱与单纵模脉冲激光的气溶胶米散射和瑞利散射光谱的卷积.
如何分离提取多纵模激光在大气中传输时产生的气溶胶米散射和瑞利散射光谱是多纵模高光谱分辨率激光雷达实现的关键.本文采用窄带干涉滤光片滤除太阳背景光,设计可调谐的马赫-曾德尔干涉仪分离提取多纵模激光的大气回波散射信号,并利用马赫-曾德尔干涉仪双通道输出互补性的原理,推导了气溶胶光学参量的反演过程.系统仿真结果表明,本文所设计的紫外域多纵模高光谱分辨率激光雷达,能够实现全天时条件下10 km高度范围内的气溶胶光学特性精细探测.
[1]Mao J T,Zhang J H,Wang M H 2002Acta Meteorolog.Sin.60 625(in Chinese)[毛节泰,张军华,王美华2002气象学报60 625]
[2]Cheng C L,Wang G H,Zhou B H,Meng J J,Liu J J,Cao J J,Xiao S 2013 Atmos.Environ.81 443
[3]Andreae M O,Rosenfeld D 2008 Earth-Sci.Rev.89 13
[4]Yu H B,Liu S C,Dickinson R E 2002 J.Geophys.Res.107 4142
[5]James D K 1981 Appl.Opt.20 211
[6]Fiocco G,Dewolf J B 1968 J.Atmos.Sci.25 488
[7]Zhao M,Xie C B,Zhong Z Q,Wang B X,Wang Z Z,Dai P D,Shang Z,Tan M,Liu D,Wang Y J 2015 J.Opt.Soc.Korea 19 695
[8]Cheng Z T,Liu D,Zhou Y D,Yang Y Y,Luo J,Zhang Y P,Shen Y B,Liu C,Bai J,Wang K W,Su L,Yang L M 2016 Opt.Lett.41 3916
[9]Shimizu H,Lee S A,She C Y 1983 Appl.Opt.22 1373
[10]Hair J W,Caldwell L M,Krueger D A,She C Y 2001 Appl.Opt.40 5280
[11]Liu J T,Chen W B,Song X Q 2010 Acta Opt.Sin.30 1548(in Chinese)[刘金涛,陈卫标,宋小全2010光学学报30 1548]
[12]Su W Y,Schuster G L,Loeb N G,Rogers R R,Ferrare R A,Hostetler C A,Hair J W,Obland M D 2008 J.Geophy.Res.Atmos.113 202
[13]Rogers R R,Hostetler C A,Hair J W,Ferrare R A,Liu Z,Obland M D,Harper D B,Cook A L,Powell K A,Vaughan M A,Winker D M 2011 Atmos.Chem.Phys.11 1295
[14]Hoffman D S,Repasky K S,Reagan J A,Carlsten J L 2012 Appl.Opt.51 6233
[15]Imaki M,Takegoshi Y,Kobayashi T 2005 Jpn.J.Appl.Phys.44 3063
[16]Di H G,Zhang Z F,Hua H B,Zhang J Q,Hua D X,Wang Y F,He T Y 2017 Opt.Express 25 5068
[17]Hua D X,Uchida M,Kobayashi T 2004 Opt.Lett.29 1063
[18]Jin Y,Sugimoto N,Nishizawa T,Ristori P,Otero L 2016 Proceeding of the 27th International Laser Radar Conference New York City,USA,July 5–10,2015 p02006
[19]Ristori P,Otero L,Jin Y,Sugimoto N,Nishizawa T,Quel E 2016 Proceeding of the 27th International Laser Radar Conference New York City,USA,July 5–10,2015 p06005
[20]Cheng Z T,Liu D,Zhang Y P,Liu C,Bai J,Wang D,Wang N C,Zhou Y D,Luo J,Yang Y Y,Shen Y B,Su L,Yang L M 2017 Opt.Express 25 979
[21]Liu Z Y,Kobayashi T 1996 Opt.Rev.3 47