CN102254090B - 采用核密度估计的闪电密度分布及ng值估计方法 - Google Patents

Abstract

本发明公布了一种采用核密度估计的闪电密度分布及NG值估计方法。本发明本方法实现：①结合地理信息系统（GIS）技术，构建和维护雷电空间数据库，并建立雷电实时定位监测系统来获取雷电发生的信息源，以空间要素形式实时地呈现闪电的定位分布；②利用核密度估计方法，求得区域闪电密度分布，再将用户选择的图形区域作为条件层作切割，形成区域密度分布图，并且能够消除边缘效应。③计算核密度估计的区域总面积（公里）以及年总闪电次数，求得对应区域的NG值，形成闪电密度图和数值结合的气象产品，最后通过共享网络来发布和共享雷电资料分析产品，通过以上闪电密度统计方法实现雷电业务服务的根本目标。

Description

采用核密度估计的闪电密度分布及NG值估计方法

技术领域

本发明属于GIS与气象信息技术交叉领域，地理信息系统(GIS)强大之处是能将空间地理分布状况及所具有的属性进行数字存储，建立数据管理系统，同时开发各种分析和处理功能，以便快速获取信息，并将处理结果以地图、图元及数据的形式表示出来。因此，处理海量且庞杂的实时雷电数据信息，需要通过准确快捷的方法求得区域闪电密度分布及闪电NG值，为雷电防护的分析和决策提供了科学手段，大大提高了以多种形式对雷电数据资料的统计效率。

背景技术

自然界雷电的发生具有很强的随机性、瞬时性和危险性，它与大气强对流系统密切相关。雷电监测网络主要由地面大气电场监测网、闪电定位监测网、地面自动观测站组成，这个网络包含3-5个电场仪、覆盖全省的9个闪电定位观测点和1-2个自动气象观测站。在当天雷电事件发生后，每个仪器都能即时地监测到的每次雷电的原始数据信息(包括序号、时间、经纬度、强度等信息)，并且以文本信息即时地更新传到服务器中，有时每天的雷电次数少则几条，而梅雨天的次数可能多达数万条。因此，当雷电的原始数据资料以文本记录形式存放时，空间庞大且信息冗杂，而且其灾害分布形势不具有直观性和鲜明性。

地理信息系统(GIS)是20世纪60年代中期兴起的一门交叉边缘学科。GIS利用计算机建立地理数据库，将空间地理分布状况及所具有的属性进行数字存储，建立数据管理系统，同时开发各种分析和处理功能，以便快速获取信息，并将处理结果以地图、图元及数据的形式表示出来。因此，结合雷电监测网络的实时雷电数据信息，并通过GIS的强大的空间查询分析和图形可视化功能为雷电防护的分析和决策提供了科学手段，一方面能够利用GIS的地图及图层加载功能，实现雷电精确定位及资料回放显示效果；另一方面能够利用GIS的空间查询分析及专题图制作功能，实现雷电资料查询显示以及闪电密度分析功能，大大提高了以多种形式对雷电数据资料的统计效率。

参考文献

[1]陈家宏，樊灵孟.雷电定位系统和雷电参数[C].The 3rd Asian Light-ningProtection Forum.Taisei Corp，Tokyo，Japan，2005.

[2]冯民学，焦雪等.江苏省雷电分布特征分析[J].气象科学，2009，29(2)：246-250

[3]刘岩，王振会，张慧良等.2007年浙江地区梅雨期闪电特征的分析[J].气象科学，2009，29(2)：225-228.

[4]刘艳，马劲松.核密度估计法在西藏人口空间分布研究中的应用，西藏科技，2007(4)，1498-1504

[5]杨宜东，孙志辉，张净.基于核密度估计的分布数据流离群点检测.计算机研究与发展，2005，42(9)

[6]国网武汉高压研究院.福建省雷电参数研究[R].武汉：国网武汉高压研究院，2006.

[7]Krider E，Pifer A E，Vance D L.Lightuing direction finding system for forest firedetection.

Ball Amer Metor SoC，1980，61：980-986.

[8]B.W Silverman，Density Estimation for Statistics and Data Analysis，Momographson Statistics and Applied Probability，1986.

发明内容

本发明目的是针对现有技术存在的缺陷提供一种，提供采用核密度估计的闪电密度分布及NG值统计方法，该方法有利于在雷电业务平台中进行闪电分布密度图及灾害指数NG的计算。

本发明为实现上述目的，采用如下技术方案：

本发明采用核密度估计的闪电密度分布及NG值估计方法，包括以下步骤：

(1)读取闪电资料数据库中的时间、位置信息、属性信息，形成点空间元素集合；

(2)确定平滑带宽，计算闪电密度计算尺度及密度函数，求出闪电密度集体；

(3)将闪电点矢量图层按照密度值进行渲染，按照条件区域进行适量图层切割；

(4)将矢量图层栅格化，并且按照上述条件区域进行栅格图层切割；

(5)计算用户选择的区域面积及区域闪电次数，求出闪电NG值，最终密度图与数字综合展示。

优选地，步骤(1)的实现方法是：输入时间跨度条件后，读取闪电资料数据库中相应闪电的信息，形成对应的二维空间链表；再利用ArcGIS Engine 9.2组件读取二维空间链表，将数据信息转成GIS地图元素，最终形成闪电点图层。

优选地，步骤(2)的实现方法是：闪电集合S是闪电S_i空间R上的离散点集，S＝{s₁，s₂，...，s_n}，给定带宽bandwidth或平滑量τ，计算点密度

${\mathrm{\λ}}_{\mathrm{\τ}}\left(s\right)=\frac{1}{{\mathrm{\δ}}_{\mathrm{\τ}}\left(s\right)}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\frac{1}{{\mathrm{\τ}}^2}\×k\left(\frac{s-s_i}{\mathrm{\τ}}\right)$

式中，带宽参数τ＞0，是一个以s为中心的圆的半径；δ_τ(s)代表了密度计算的尺度，k是一个概率密度函数即核心kernel密度函数；并将点密度添加到相应闪电元素的字段density中，i＝1，2，...n。

优选地，步骤(3)的实现方法是：根据闪电元素字段density，选择色带、区间值参数，利用ArcGIS Engine 9.2组件进行渲染，形成矢量面图层。

优选地，步骤(4)的实现方法是：利用ArcGIS Engine 9.2组件，设定分析环境参数，包括地理分析环境参数，通过执行IGeoAnalysisEnvironment接口实现；栅格分析环境参数，通过执行IRasterAnalysisEnvironment接口实现对矢量图层的栅格化。

优选地，步骤(5)的实现方法是：将所有闪电与条件区域作交集运算，求得该区域的闪电总次数，并计算条件区域的总面积，即可求出闪电NG值：

NG＝M/Y/P

其中，M为闪电总次数，Y为计算年限，P为区域面积。

本发明可以对复杂且庞大的雷电灾害数据进行时域和空间分布特征分析，利用核密度估计法求得区域内闪电密度分布以及年平均地面落雷密度(NG)值，并借助强大的GIS空间数据表现形式，了解区域内闪电的分布特征对确定闪电密度有重要意义，以便气象预报人员能够出色地完成闪电观测及分析工作。通常雷电地基监测网络由地面大气电场监测网、闪电定位监测网、地面自动观测站组成，但是粗糙的原始数据依然无法及时、便捷地提供给雷电防护人员观测分析。

附图说明

图1：本发明原理图。

具体实施方式

如图1所示，本发明采用核密度估计的闪电密度分布及NG值统计方法，包括以下步骤：

(1)读取闪电资料数据库中的时间、位置信息、属性等信息，形成点空间元素集合；

(2)确定平滑带宽，计算闪电密度计算尺度及密度函数，求出闪电密度集体；

(3)将闪电点矢量图层按照密度值进行渲染，按照条件区域进行适量图层切割；

(4)将矢量图层栅格化，并且按照上述条件区域进行栅格图层切割；

(5)计算用户选择的区域面积及区域闪电次数，求出闪电NG值，最终密度图与数字综合展示。

步骤(1)的实现方法是：闪电资料数据中存放的是有规则的数据，包括行号、时间、经纬度、强度、陡度、省、市、县等信息，当然也包括冗余和噪声数据，需要消除冗余和噪声数据形成规整的记录型数据。

表1闪电资料数据库字段

输入时间跨度条件后，读取闪电资料数据库中相应闪电的信息，形成对应的二维空间链表；再利用ArcGIS Engine 9.2组件读取二维空间链表，将数据信息转成GIS地图元素IPointFeature，最终形成闪电内存点图层IPointFeatureLayer。

步骤(2)的实现方法是：闪电集合S是闪电Si空间R上的离散点集，S＝{s₁，s₂，...，s_n}，给定带宽(bandwidth)或平滑量τ，在二维空间内计算点密度

${\mathrm{\λ}}_{\mathrm{\τ}}(x,y)=\frac{1}{{\mathrm{\δ}}_{\mathrm{\τ}}(x,y)}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\frac{1}{{\mathrm{\τ}}^2}\×k(\frac{x-x_i}{\mathrm{\τ}},\frac{y-y_i}{\mathrm{\τ}})$

式中，带宽参数τ＞0，实际上就是一个以s为中心的圆的半径，δ_τ(s)代表了密度计算的尺度，如下式：

${\mathrm{\δ}}_{\mathrm{\τ}}(x,y)=\∫\underset{R}{\∫}\frac{1}{{\mathrm{\τ}}^2}\×k(\frac{x-x_u}{\mathrm{\τ}},\frac{y-y_u}{\mathrm{\τ}}){\mathrm{dx}}_u{\mathrm{dy}}_u$

k是一个概率密度函数，即核心(kernel)密度函数，如下式：

$k(x,y)=\left\{\begin{array}{cc}\frac{3}{\mathrm{\π}}(1-{(x^2+y^2)}^2)& u^{T}u\≤1,\\ 0& u^{T}u>1.\end{array}\right.$

并将点密度添加到相应闪电元素的字段density中；

当满足条件(x-x_i)²+(y-y_i)²≤τ²的点集P(x，y)全部位于区域R上时，函数δ_τ(x，y)＝1，否则δ_τ(x，y)＜1，起到消除边缘效应的作用。

在核密度估计中，带宽τ的确定或选择对于计算的结果影响很大。一般而言，随着τ的增加，空间上点密度的变化更为光滑；当τ减小时，点密度变化突兀不平。

步骤(3)的实现方法是：根据闪电元素字段density，选择色带、区间值等参数，利用ArcGIS Engine 9.2组件对FeatureLayer进行渲染，形成矢量面图层；

步骤(4)的实现方法是：利用ArcGIS Engine 9.2组件IConversionOp接口下的RasterDataToLineFeatureData、RasterDataToPointFeatureData、RasterDataToPolygonFeatureData方法可以实现从栅格图层中提取出点状、线状、面状矢量要素，在AE中，可以实现提取矢量要素的栅格对象有Raster，RasterDataset，RasterBand和RasterDescriptor；ToRasterDataset方法可以实现将点、线、面、矢量要素栅格化，在AE中，可以进行栅格化的矢量要素对象有featureclass，featureclassdescriptor，featurelayer；ToRasterDataset方法可以实现矢量要素之间的互相转换。

设定接口方法中需要的环境参数，包括地理分析环境参数，通过执行IGeoAnalysisEnvironment接口实现；栅格分析环境参数，通过执行IRasterAnalysisEnvironment接口实现对闪电矢量图层的栅格化。

步骤(5)的实现方法是：将所有闪电与条件区域作交集运算，求得该区域的闪电总次数，并计算条件区域的总面积(km2)，即可求出闪电NG值，公式如下：

NG＝M/Y/P

其中，M为闪电总次数，Y为计算年限，P为区域面积。从而，得到区域闪电密度分布图和精确的闪电年平均每平方公里闪电次数NG值。

以上是本发明的较佳实施例，凡依本发明技术方案所作的改变，所产生的功能作用未超出本发明技术方案时，均发球本发明的保护范围。

Claims (5)

1.一种采用核密度估计的闪电密度分布及NG值估计方法，其特性在于，包括以下步骤：

（1）读取闪电资料数据库中的时间、位置信息、属性信息，形成点空间元素集合，最终形成闪电点图层；

（2）确定平滑带宽，计算闪电密度计算尺度及密度函数，求出闪电密度集体；

实现方法是：闪电集合S是闪电s_i在空间R上的离散点集，S={s₁,s₂,...,s_n}，给定带宽或平滑量τ，在二维空间内计算点密度

${\mathrm{\λ}}_{\mathrm{\τ}}(x,y)=\frac{1}{{\mathrm{\δ}}_{\mathrm{\τ}}(x,y)}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\frac{1}{{\mathrm{\τ}}^2}\×k(\frac{x-x_i}{\mathrm{\τ}},\frac{y-y_i}{\mathrm{\τ}})$

式中，参数τ＞0，实际上就是一个以s_i为中心的圆的半径，δ_τ(s)代表了密度计算的尺度，如下式：

${\mathrm{\δ}}_{\mathrm{\τ}}(x,y)=\underset{R}{\∫\∫}\frac{1}{{\mathrm{\τ}}^2}\×k(\frac{x-x_i}{\mathrm{\τ}},\frac{y-y_i}{\mathrm{\τ}}){\mathrm{dx}}_i{\mathrm{dy}}_i$

k是一个概率密度函数，即核心密度函数，如下式：

$k(x,y)=\left\{\begin{array}{cc}\frac{3}{\mathrm{\π}}(1-{(x^2+y^2)}^2)& u^{T}u\≤1,\\ 0& u^{T}u>1.\end{array}\right.$

并将点密度添加到相应闪电元素的字段density中，i=1，2，…n；

当满足条件(x-x_i)²+(y-y_i)²≤τ²的点集P(x,y)全部位于区域R上时，函数δ_τ(x,y)=1，否则δ_τ(x,y)<1，起到消除边缘效应的作用；

（3）将闪电点图层按照密度值进行渲染，按照条件区域进行适量图层切割，形成矢量面图层；

（4）将矢量面图层栅格化，并且按照上述条件区域进行栅格图层切割；

（5）计算用户选择的区域面积及区域闪电次数，求出闪电NG值，最终密度图与数字综合展示。

2.根据权利要求1所述的采用核密度估计的闪电密度分布及NG值估计方法,其特征在于，步骤(1)的实现方法是：输入时间跨度条件后,读取闪电资料数据库中相应闪电的信息,形成对应的二维空间链表；再利用ArcGIS Engine9.2组件读取二维空间链表，将数据信息转成GIS地图元素，最终形成闪电点图层。

3.根据权利要求1所述的采用核密度估计的闪电密度分布及NG值估计方法,其特征在于，步骤(3)的实现方法是：根据闪电元素字段density，选择色带、区间值参数，利用ArcGIS Engine9.2组件进行渲染，形成矢量面图层。

4.根据权利要求1所述的采用核密度估计的闪电密度分布及NG值估计方法,其特征在于，步骤(4)的实现方法是：利用ArcGIS Engine9.2组件，设定分析环境参数，包括地理分析环境参数，通过执行IGeoAnalysisEnvironment接口实现；栅格分析环境参数，通过执行IRasterAnalysisEnvironment接口实现对矢量面图层的栅格化。

5.根据权利要求1所述的采用核密度估计的闪电密度分布及NG值估计方法,其特征在于，步骤(5)的实现方法是：将所有闪电与条件区域作交集运算，求得该区域的闪电总次数，并计算条件区域的总面积，即可求出闪电NG值：

NG=M/Y/P

其中，M为闪电总次数，Y为计算年限，P为区域面积。