基于GIS空间关系的温度序列复原方法

沈军^{1, 2}，杨敏华¹，钟荣华²

(1. 中南大学 地球科学与信息物理学院，湖南 长沙，410083；

2. 湖南益阳市气象局，湖南 益阳，413000)

摘要：在分析测量温度序列特征的基础上，结合GIS空间关系模型和空气热传导方程解的特点，建立温度序列复原的目标函数。利用多尺度小波分析方法和基于时-频变换稀疏理论提出基于GIS空间关系的气象温度序列复原方法。研究结果表明：在温度序列模型未知且各温度序列不完全独立的复杂条件下，该方法对于多源温度序列复原具有比较理想的效果。

关键词：盲信号复原；多尺度分析；稀疏矩阵；质量控制；气象观测

中图分类号：P203            文献标志码：A         文章编号：1672-7207(2012)03-1047-10

Recovering of temperature series data based on GIS

SHEN Jun^{1, 2}, YANG Min-hua¹, ZHONG Rong-hua²

(1. School of Geosciences and Info-Physics, Central South University, Changsha 410086, China;

2. Meteorological Bureau of Yiyang, Yiyang 413000, China)

Abstract: A new method for the recovering of temperature series from some adjacent automatic weather stations (AWS) was proposed with a GIS platform. By combining the short Fourier transformation with a thermodynamics equation and multi-scale analysis method, a goal function of the blind source separation was derived to get the optimal result from a blind source seperation(BSS) module. The results identify the usefulness of the propose method.

Key words: blind source separation; multi-scale analysis; sparsity matrix; quality control; meteorological observation

温度序列监控广泛应用于设备监测、医疗监测、植被培育、环境监测和气象观测等领域。基于无线传感器网络技术的气象自动站以观测连续性、实时性以及规范性等优点被广泛应用于气象观测中并在灾害监测和临近预报等行业中起着不可替代的作用^[1]。地球表面热辐射能力不一致使得某一自动观测站记录的温度为周边不同热辐射地表相互作用的结果。因此，从混合的观测温度序列中复原出原温度序列在气候研究和环境监测等领域中均有着重要意义^[2]。Ilin等^[3]利用多源降噪分离法(Denoise source separation)对降水数据、温度数据以及气压的观测混合数据进行分析，得到了厄尔尼诺现象动力因子和该现象的运动周期，并基于频谱分析技术和多源降噪分离法对20世纪五六十年代的逐日观测资料进行分析得到厄尔尼诺-南方涛动现象的发生周期^[4]。尽管气象自动观测站是无线传感器技术与气象观测技术相结合的新观测设备，相关数据处理研究研究成果暂未形成，但与其相关的盲信号复原问题在其他领域的应用需求不断被提出并受到了国内外学者的极大关注。Bell等^[5]利用递归神经网络和Hebb学习规则实现源信号分离。此后各种基于主成分分析方法(ICA)的多源信号复原方法大量被提出，如：Bell等^[5]提出最大化信息法^[5]；Aapo等^[6-7]提出基于独立信号的最大似然法；Mendel等^[8-9]提出的用于高斯噪声的高阶统计法等。该类方法均假定待估计信号的分布函数已知且累计误差服从期望为0的高斯分布。近来，Bobin等^[10-11]在可分多源信号和稀疏性的约束条件下提出了基于稀疏概念的信号复原方法并将该方法与形态学方法相结合在图像分割中取得了很好的纹理分离效果。在不确定的多源盲信号复原研究中，Abrard等^[12-13]分别在时间域和频率域提出了多种基于短时傅里叶变换的不确定信号复原方法如LI-TIFCORR法和LI-TICOM方法。该类方法在信号源个数与观测点不相等时能取得较好复原效果，但该算法仍要求各信号源在时间-空间域中有较好的可区分性。然而，在气象观测领域，受大气热传导方程等基本规律的影响，各站点的温度序列在时间域和空间域中均有一定相关性，具有混沌特性的温度序列难以用准确的数学模型进行表述或识别。此外，受大气扰动和其他气象条件如风速风向等变化的干扰，温度变化的差值序列也不一定服从高斯分布。传统复原方法仅仅依靠温度源的频率差异进行分离因而难取得较好复原效果。在此，本文作者在分析自动站测量温度序列特征的基础上，利用多尺度分析技术并结合GIS空间关系模型和热传导方程解的特点，建立温度序列复原的目标函数，基于目标函数最大化原理提出多源温度序列复原方法。该方法克服了传统方法的不足，将温度序列的热动力模型以及频谱特征相结合，因而表现出较好的复原特点，能有效应用于有相关性的多源信号的复原应用中。

1  基于GIS模型的温度序列BSS建模

1.1  温度序列特征分析与多源温度序列复原建模

图1(a)所示为基于地理信息系统(GIS)的地面气象观测站网络在某一时刻的局部区域温度分布，图1(b)所示为与之对应的观测站24 h同步温度序列。该温度观测频次为288个/(5 min)，观测精度为0.1 ℃。从图1(b) 可以看出：相邻观测站点，尽管总体温度变化趋一致，但在不同的地表条件下，相关的观测站温度变化率有很大不同。另外，在相同纬度、相同时间段内，局部区域温度变化受地理环境如城市热岛效应和河流、森林植被率以及太阳高度角等的影响，尤其受局部地形和植被指数以及地理环境等综合反辐射率的影响比较大，针对该特点，本文提出基于GIS空间关系的温度模型(如图2所示)。其中热源为不同反射率地表对外热辐射中心，该热源经过混合矩阵A作用   后，进入传感器，然后在加性噪声作用下得到观测序列v_i。

在一定空间范围内，利用传感器空间位置关系和热传导过程建立如下模型：用v_j表示与第i个热源对应的观测值，a_i,j表示第i 个热源对第j个传感器产生的影响系数，u_i为第i个热源对测量站的贡献值。以上模型可以用如下方程表示：

       (1)

其中：表示与第i个传感器对应的噪声信号，满足，。将式(1)表示成矩阵的形式

图1  基于GIS的观测网络以及相关温度序列

Fig.1  A GIS based observation net and correspond temperature series

图2  基于GIS的气象自动站温度复原模型

Fig.2  A model of blind source of temperature series based on GIS

有：。其中：

；；

。

根据站点间的空间位置关系，同时降低问题的复杂度，不妨假定矩阵为对称矩阵，即。在以上模型中，向量U和V是随时间变化的变量。但是，在一定的气象条件下(如无大风等强对流天气)，权矩阵A可近似为常系数矩阵。基于稳定天气的假设条件下并假定无噪声和无奇异值干扰，建立基于气象温度场的温度序列复原模型，即：

V=AU                   (2)

已知观测向量序列V，在矩阵A未知的情况下，求U的问题属于盲信号分离(blind source separation)范畴。图1(b) 描述了相邻空间域中不同地面反射率的5个观测站的温度序列。从图1(b)可以看出温度序列与传统信号有如下几点不同：

(1) 空气下垫面气温变化受天气系统变化以及周边环境的影响，局部区域温度及其相互影响只能在有限的时间内保持相对固定的关系，即矩阵A在有限的时间内为常系数矩阵；

(2) 不同于心电图(ECG)信号和声音信号，温度序列不是以某一频率为中心的周期性信号。

(3) 温度变化序列表现为在一定规律下的震荡性变化过程。

(4) 温度的差值序列不是以期望为0的周期性时间序列。

温度序列的以上特点决定了该类数据复原方法不同于传统的语音等信号复原方法，有着特殊性和不确定性。为了提高信号复原的效果，本文基于温度序列的物理特征方程构造1个特定的目标函数。

1.2  建立复原目标函数

为了在转换矩阵A未知的条件下，根据观察序列尽可能的恢复出原始的温度变化，确保温度序列在一个合理的范围内，假定如下：

(1) 在进行温度序列复原的过程中，混合矩阵A是1个常数矩阵，且其阶数等于信号源个数；

(2) 在时间段i内，各温度源相互独立，即温度序列为相互独立的信号源且相互不相关，即；

(3) 在任何子时间段内，各温度序列不存在对应比例关系，即不存在2个同类的温度序列；

假设由条件(1)可以得到混合矩阵是一满秩矩阵，即存在A的逆矩阵A^-1，满足：U=A^-1V。条件(2)确保了本问题被约束在线性复原范畴内；条件(3)确保了式(1)有唯一解。根据以上模型和假定条件，基于观测数据序列求出矩阵A的估计值并对其求逆可以得到本文问题的解。因此，温度序列复原的本质就是根据观测矩阵V得到最优的估计混合矩阵A，并利用得到的矩阵求出与观测序列对应的原序列U。为了得到与式(1)相同系数的解，引入相关数描述2个信号之间的互信息量。令信号

定义它们的比值函数为：

  (3)

设表示第j个信号和第i个信号在时间t的比值, 方差为。若该序列有k个子空间，第t个空间的方程为，令，当即随机变量的方差为0时，在某时间区域T内，为常数，不妨令 ，由假定条件(3)可知，此时间域内求得该混合矩阵的1个元素等价于求观测矩阵标准化后原系数的有效估计值。为此，定义1个目标函数：

         (4)

另外，从文献[14]可知大气热传递方式主要有长波辐射、空气对流、热传导等形式。由于气象自动站以及温度传感器建立在离地面高度1.5 m的百叶箱  内^[15]，在静风条件下，地面对温度传感器的效应主要表现以空气对流和热传导形式。不同距离的2个热源(u₁, u₂)对同一传感器的贡献率可以用热对流(热平衡)方程进行分析：

令热源u₁和u₂到传感器的距离分别为r₁和r₂    (r₁＜r )，并假定有如下偏微分方程成立：

满足边界条件

       (5)

在频率域中，该方程解的结构形式可表示为：。因此，

其中：c₁和c₂为系数；和为与热源到观测站距离差(r₂-r₁)相关的参数。令为时间域内相应的

解，其傅里叶算子为：。令，其中，为热源估计值。总目标函数最小化表示为：

            (6)

为了能实现时间序列目标函数的最小化，可以将时间序列分为若干个活动窗口，对每个活动窗口区间内的观测值求序列之间的系数最小估计值。根据以上假定和分析可知：目标函数成立的必要条件是多源温度源序列在某个区域内可以分辨，即要找到某个时间域，使得在这个时间域内，多源温度序列的是呈稀疏分布且相互可分的。然而，对于n组有限长度的随机信号序列，要找到n个活动窗口使得每一窗口仅对应1个温度观测序列的概率极低^[16]。

2  基于小波的STFT温度序列复原

根据气象知识，地面某观测点空气温度变化主要受太阳辐射的影响而呈周期性变化。一方面，在无强对流天气过程的情况下，周围空气温度也有着同样的时间周期性。在同等天气条件下，温度序列是以24 h为周期的低频函数；另一方面，由于观测站的空气温度记录值是多个热源相互作用的结果，是一个动态变化的热交换过程，因而表现出震荡波形(见图1(b))，且其震荡频率和幅值与观测场周边环境有关。若用x(t)表示某一温度序列对应的函数，不失一般性， 采用汉宁(Hamming widow)对x(t)进行加窗傅里叶变换。相应的短时傅里叶变换(STFT)为：

        (7)

其中：h(t)为加窗函数。由于温度是一个连续场，对连续变化的温度序列，当采样间隔较短时(1~5 min)，对于多数观测样本，有，即，所以，有

             (8)

式(8)表明：基于温度序列的加窗傅里叶变换得到的频率谱是集中分布的。

图3所示为对图1(b)中5个同步温度序列进行STFT变换得到频率域分布和时-频域分布，其中，图3(a)所示为对图1(b)所示温度序列进行短时短时傅里叶变换后的频谱分布状况。从图3(a)可见：在时间域内，大多数点都积聚在原点某一小邻域内；尽管在时间域中多源温度序列有不同的变化过程，但在频率域中，它们拥有相同的相空间分布。图3(b)描述了与图3(a)对应的功率谱分布。从图3(b)可见：温度序列的傅立叶功率谱也是高度集中。

图3  运用STFT对源数据进行处理的时间域和频率域结果

Fig.3  Output of Figure 1(b) with STFT in Temporary and frequency

为了提高频率域中温度序列复原的效率和观测序列之间的独立性，对以上观测序列引入多尺度分析技术。在多尺度分析过程中，低频信号对应于一段时间内信号的变化趋势，高频部分对应于局部时域内信号的变化过程。为了将不同温度序列的变化情况进行对比分析，取温度序列的高频部分进行分析。图4所示为对图1(b)所列出的5个相邻站点观测序列进行小波变换后的高频部分。从图4可以看出：图1(b)所列出的温度序列尽管有相同的变化过程，但在高频变化部分，不同温度序列反映了不同的地表辐射过程。图4中，星号虚线对应于图1(b)中的热源5。该序列所在的观测站位于植被较少且比热较小的岩石区域，因而表现为幅度变化较大的温度变化过程。对这些高频序列进行STFT变换，得到该类序列在时间域和频率域中的分布，如图5所示。

图5所示为观测序列经过小波分解后进行短时傅里叶变换结果，其中：图5(a)所示为STFT变换后复数域分布情况，可以看出在同一复数域空间中，5个温度序列相互依赖性大大降低，点列分布均匀。相对于图3(b)，图5(b)中高能量频率区(红色条纹)也均匀分布于时-频空间，因此，相对于图3，对于可分性较差的多源低频信号，其高频部分在频率域中有比较好的分辨性。本文基于时-频变换稀疏理论提出了基于多尺度的温度序列复原方法，将观测信号分成高频部分和低频部分：在低频部分，对奇异值进行检测并对明显奇异值进行检测并修正；在高频部分，用盲信号复原(HFBSS)方法进行多源信号复原。基本流程如图6所示。

图4  对图1(b) 进行sym2小波分解得到的高频温度序列

Fig.4  Resolution of sym set based on wavelet with high frequency

图5  基于STFT的高频温度序列时域分布和频率分布

Fig.5  Result of higher part data series with STFT

图6  系统处理流程图

Fig.6  Process of proposed method

此外，由于本文假定观测序列不受奇异值和噪声干扰，因此，在本文实验部分省去了低频的奇异值检测过程而直接对高频数据进行复原处理。

2.1  基于小波分解的温度序列预处理

根据小波变换的“线性不变”的特点可知：基于高频信号复原得到的混合矩阵与源信号复原的混合矩阵等价，即基于高频信号估计的混合矩阵与源信号混合矩阵等价。不妨令Y^h和X^h分别代表向量Y和X小波分解后对应的高频部分，本复原问题就变成了基于观测向量Y^h，求与之相对应的源信号向量X^h的过程。该过程的直接目的对混合矩阵A的合理估计，并使得其逆矩阵A^#满足X^h=A^#*Y^h。

2.2  基于短时傅里叶变换的高频温度序列复原

经过小波分解得到的高频温度序列虽然类似于声音等周期信号，但在时间域内对混合矩阵估计所需“稀疏性”等价于在空间域中要求各输入信号在对应的频率内有很强的“稀疏性”和“局部性”。然而，在实际问题中，多源信号在频率域中的任何频率段均能保持互不相关且呈聚集性(好的分辨特性)的概率很小^[16]；另一方面，根据加窗傅里叶变换的线性不变特性很容易找到某个频率，使得在该频率内某一信号相对于其他信号表现更加明显，即在该频率域中的某个时间段内各信号源是可区分的。为了能找到这个频率值并得到实用的稀疏模型，将原时域内的信号进行短时傅里叶变换(Short-time frequency transfer)处理，以提高信号在有限时间段内相互独立的概率。类似于式(3)定义时间-频率域中的目标函数：

         (9)

其中：表示在时间-频率域中信号第i个信号在时间t和频率w所对应的值；为由函数f(x)所决定的新的变量方差。由于和均大于0，当目标函数达到最小为0时，必有,且；当时，可估计出混合系数。

3  试验研究

为了定量描述估计矩阵与目标矩阵的差异，根据多源信号间的相对稀疏特征计算真实序列与估计序列的相似度(S_IRI) ，根据目标矩阵与估计矩阵差值定义距离指数。

3.1  信号干扰比

信号干扰比(Single interference ration, R_SI)是1个表征主分量与混合分量的比值，即对于给定的传感器i，起主要作用的热源与其他热源贡献率的比值：

其中：分别表示周围热源j 对热源i附近的传感器的贡献因子。该统计量表明在一个长度为t的序列中，主要热源与非主要热源的平均值比值。

根据该定义，分别对源信号和估计信号建立信号干扰比。源信号干扰比(S_{IR_in})和估计信号干扰比(S_{IR_out})为：

     (10)

    (11)

其中：和分别表示实际混合矩阵和估计混合矩阵的第i行j列元素；为了评估复原算法的优劣，可

以用来表示估计算法与源信号本

身的一致性。该定义表明：S_IRI绝对值越小，表明目标矩阵与估计矩阵越一致，当两者完全一致时，S_IRI等于0；相反，2个序列相差越大。

3.2  矩阵距离

在信号复原仿真实验中，混合矩阵的估计误差决定了信号估计的有效度。为此，定义了2个矩阵之间的距离为：

其中：和分别表示真实矩阵与估计矩阵对应的元素。

3.3  实验过程

为了证明本方法的有效性，利用该方法进行仿真实验和真实数据实验。在第1个实验中， 选用经过检验的气象数据进行仿真实验。对4个不同地区的观测站得到的温度序列进行模拟。在本次实验中，实验数据序列长度8 732。对该序列经过小波分频后，得到长度为4 456的高频温度序列(如图7中(a)，(d)和(g)所示)。将这4个高频温度序列按给定混合矩阵P进行混合得到模拟的“观测”数据序列。矩阵P满足热传导基本规律特征并设定如下：

图7  信号复原的结果与源信号的比较

Fig.7  Comparison of source data, mixed data and recovered data

混合后的序列对应于图7中(a1)，(b1)和(c1) 3个运算序列。为了简化计算，在本次实验中，将式(5)的定为常数，然后，对各序列分别进行短时傅里叶变换。实验结果表明：变换重叠率因子选为0.25~0.55，窗口内随机样本数选择为50~60时有比较稳定的结果。根据采样数据特点，在13个不同频率空间内对观测序列进行运算。最后求得估计混合矩阵为为：

对比P和可以看出，估计矩阵与原矩阵有相同的矩阵结构，且矩阵的距离D_M=0.221 4。

图7中有3个数据序列(a)，(d)和(g)经过矩阵P混合后生成了(b)，(e)和(h)序列，(c)，(f)和(i)为该组序列的复原结果。从图7可以看出：源信号经过混合后，在整体上数据特征基本不变，但在其他数据点处尤其是“零点”值附近的许多数据发生了变化：如在第一行组图中，许多原来为0的观测值经过变换后均变成了非零值(如从500~1 000和2 700~3 000的观测区域内，相应的零点均发生了改变)，但经过复原处理后，许多点基本恢复到0或更小的领域内。此外，源序列经过其他序列干扰后，某些特征点发生了变化，如在第492个时间历元处以及在2 000历元处，经过混合后，相应的特征数据消失，经过复原后，该特征点均在图(c)中得到恢复；在第2行组图中的原数据序列300~310个观测区间中，数据序列a₀，该值为-0.5~0.3，源数据b₀对应的值为-0.2~0.2，数据序列c₀对应的值为-0.2~0.6。经过矩阵P的混合后,在同样区间中被变成了-0.5~+0.5(如图7(e)所示)，使得观测值与实际值有明显区别，在与之对应的复原效果b₂中(如图7(f)所示)，得到的值为-0.2~+0.2。经进一步分析可知：在观测区间500~2 000的范围内，混合后的序列均发生了很大变化，主要表现在：“非零点“数据明显增多；某些特征点在混合后产生了新的变化；如在3 500附近的特征点经过混合运算后，大小和形状均发生了变化。以上观测点经过复原处理后与原数据基本一致。

在第3行组图中，观测数据(混合后数据)均出现了与第1行组图相同的情况，其中最明显的部分是在观测历元约为3 800处，特征值被明显拉长。经过复原运算后，该值在图7(i)中基本得到恢复。

从本实验运算效果来看，估计矩阵和原混合矩阵有很大相似性，且估计信号与原信号严格一致。

为了说明本方法复原的优势，在同样实验数据条件下，用LI-TIFCORR方法、对称LI-TIFCORR方    法^[13]和基于ICA全局最优法^[17]进行复原，得到输出信号干扰比S_{RI_out}、混合矩阵相似度S_IRI和矩阵距离D_M，如表1所示。

表1  运用多种方法对比效果

Table 1  Different result of many ways

从表1可见：在输出矩阵干扰(S_IRI)比以及矩阵对角相似度(S_IRI)中，基于ICA全局最优方法获得的估计矩阵信号干扰比为0.8 dB和5.7 dB。该2项参数表明利用ICA方法得到的矩阵有“对角占优”特点。然而，针对特定物理模型的信号复原过程中，本文提出的方法使得估计矩阵与实际矩阵距离最小。

表1所示结果表明：运用基于ICA的全局最优方法有比较好的主成分解，即该方法在源序列彼此独立的情况下，有比较好的复原效果。然而，该方法没考虑到实际的物理模型，忽略了多源信号相互干扰的情况，因而，在观测信号有一定相关性时失去了其复原的优势。由于本文考虑了大气热传导方程的特点，对目标函数提出了相应的约束，因而，在估计矩阵与真实矩阵一致性方面，该方法有更好的效果。

在第2个实验中，对图1(b)列出的温度序列进行复原处理。在本次试验中，采用sym2小波将各信号进行分解，对高频部分进行STFT变换处理。在STFT变化处理中，采样率上限为4 000，活动窗口长度为512，重叠区域选择35%；间隔样本选择512，得到时间-频率滤波的分布图和相对应的光谱分布图，见图5(a)和5(b)。图8所示为该序列的复原结果。各序列长度为146。图中实线表示对图1(b)进行小波变换后的高频部分。虚线为基于本方法的信号复原效果。从图8可以看出：对温度序列1和4，所得结果相对原序列较稳定；对温度序列2，3和5平稳序列所得结果有更大的方差。为了具体说明图8的复原结果，对复原前后的数据序列进行相关性分析，结果分别见表2和表3。

图8  对某一局部相邻5个自动站温度序列复原效果

Fig.8  Results of temperature series based on 5 adjacent AWS

表2  原数据之间的相关性

Table 2  Correlation between source series

表3  目标数据之间的相关性

Table 3  Correlation between the recovered data set

由表2和表3可知：2个观测序列的相关性与它们之间的距离成反比。其中在观测序列中，相关性最大的2个观测序列为温度序列3和温度序列4，相关系数为0.692 3；相关性最小的观测序列为温度序列1和温度序列5，其相关系数为0.045 3。经过运用本方法进行复原后，多源信号之间相关性有所降低，最大相关性为0.039 7，最小相关性为0.001 9。表2说明：运用本方法处理后的观测序列之间的相关性明显减少，信号之间的独立性增强。

4  结论

(1) 对气象自动观测站测量的温度序列进行分析，得出相邻观测站点温度序列之间的相互关系并建立相应的温度转换模型。

(2) 根据不同地表环境热辐射各异的特性，利用多尺度分析方法对各温度序列进行分解，得到了相关性更小的高频温度序列。

(3) 将热传导方程解的结构特征与GIS站点空间距离相结合建立混合温度序列复原的目标函数，直观描述了不同物理辐射率对同一观测站点的影响，因而在复原过程中约束了复原解的空间。

(4) 本方法对于有一定相关性的多源混合序列有比较好的复原效果。

然而，本方法是基于时间-频率域对输入信号进行统计分析，该计算结果受样本在时间-频率域中分布特征的影响，因此，频率区子空间的划分、样本点重叠率等参数对不同的温度序列有不同的影响，相关参数的确定方法有待进一步研究。

参考文献：

[1] Cardinali C. Monitoring the observation impact on the short-range forecast[J]. Quarterly Journal of the Royal Meteorological Society, 2009, 135(638): 239-250.

[2] Watkins R. Increased temperature and intensification of the urban heat island: implications for human comfort and urban design[J]. Built Environment, 2007, 33(1): 85-96.

[3] Ilin A, Valpola H, Oja E. Semiblind source separation of climate data detects El Nino as the component with the highest interannual variability[C]//Proceedings of International Joint Conference on Neural Networks (IJCNN’2005). Montreal, Canada: IEEE, 2005: 1722-1727.

[4] Ilin A, Valpola H, Oja E. Exploratory analysis of climate data using source separation methods[J]. Neural Networks, 2006, 19(2): 155-167.

[5] Bell A J, Sejnowski T. An information-maximization approach to blind separation and blind de-convolution[J]. Neural Computation, 1995, 7(6): 1129-1159.

[6] Hyv?rinen A, Oja E. Independent component analysis: Algorithms and applications[J]. Neural Networks, 2000, 3(4/5): 411-430.

[7] Lee T W, Girolami M, Sejnowski T J. Independent component analysis using an extended infomation algorithm for mixed sub-Gaussian and super-Gaussian sources[J]. Neural Computation, 1999, 11(2): 417-441.

[8] Mendel J M. Tutorial on higher-order statistics (spectra) in signal processing and system theory: Theoretical results and some applications[J]. Proceedings of the IEEE, 1991, 79(3): 278-305.

[9] Nikias C L, Mendel J M. Signal processing with higher-order spectra[J]. Signal Processing Magazine, IEEE, 1993, 10(3): 10-37.

[10] Bobin J, Strarck J L, Fadili J M, et al. Morphological component analysis: An adaptive thresholding strategy[J]. Image Processing, IEEE Transactions on, 2007, 16(11): 2675-2681.

[11] Bobin J, Moudden Y, Fadili J M. Morphological diversity and sparsity for multichannel data restoration[J]. Journal of Mathematical Imaging and Vision, 2009, 33(2): 149-168.

[12] Abrard F, Deville Y. A time-frequency blind signal separation method applicable to underdetermined mixtures of dependent sources[J]. Signal Processing, 2005, 85(7): 1389-1403.

[13] Deville Y, Puigt M. Temporal and time-frequency correlation- based blind source separation methods. Part I: Determined and underdetermined linear instantaneous mixtures[J]. Signal Processing, 2007, 87(3): 374-407.

[14] 沈钟平, 张华. 影响地面太阳辐射及其谱分布的因子分析[J]. 太阳能学报, 2009, 30(10): 1209-1215.
SHEN Zhong-ping, ZHANG Hua. Analysis on the factors affecting surface solar radiation and its special distribution[J]. Act Energe Soloaris Sinia, 2009, 30(9): 1209-1215.

[15] 中国气象局. 地面气象观测规范[M]. 北京: 气象出版社, 2003: 5-8.
China Meteorologial Adimistrator. Ground observation practice[M]. Beijing: Meteorological Press, 2003: 5-8.

[16] 白树忠, 刘琚, 孙国霞. 基于最小均方误差和稀疏特征的欠定盲源分离[J]. 山东大学学报, 2008, 38(4): 43-46.
BAI Shu-zhong, LIU Ju, SUN Guo-xia. An algorithm for under-determined blind source separation based on the least-mean-square error and sparse features[J]. Journal of Shandong University, 2008, 38(4): 43-46.

[17] Cheung Y M, Liu H L. A new approach to blind source separation with global optimal property[C]//Proceedings of the IASTED International Conference of Neural Networks and Computational Intelligence. Grindelwald, Switzerland: IASTED, 2004: 137-141.

(编辑 陈灿华)

收稿日期：2011-06-25；修回日期：2011-09-02

基金项目：国家自然科学基金资助项目(30570279)；国家自然科学基金青年科学基金资助项目(40901171)；湖南省气象局科研项目(201201)

通信作者：沈军(1972-)，男，湖南安乡人，博士研究生，从事测量数据处理、遥感影像处理研究；电话：0737-4224210；E-mail: csuva@yahoo.cn