CN113177306A - 一种基于b样条函数的欧拉解集概率密度成像方法 - Google Patents

Abstract

本发明提供一种基于B样条函数的欧拉解集概率密度成像方法，首先确定待测区域范围，测量重力矢量数据及重力梯度张量数据；构建三维重力梯度张量欧拉反褶积方程，确定滑动窗口大小，利用滑动窗口内的数据构建线性方程组；获取线性方程组的欧拉解，输出欧拉解集；将欧拉解解集拆分成不同维度的组合，并估算估计区间的上/下界及估计网格带宽，以及确定估计网格的大小；计算并存储B样条函数，并对B样条函数和对应的系数进行卷积，获取欧拉解概率密度分布结果，根据欧拉解概率密度分布结果确定不同维度数据的异常源空间位置及类型。本发明减少了欧拉解缪解处理过程，优化了确定欧拉解优解的机制，实现了欧拉解簇的分离，进而有效地标识地下异常源。

Description

一种基于B样条函数的欧拉解集概率密度成像方法

技术领域

本发明涉及重力勘探技术领域，特别是涉及一种基于B样条函数的欧拉解集概率密度成像方法。

背景技术

欧拉反褶积是一种半自动/自动估算场源位置的位场数据解释方法，尤其适合于大面积位场数据的分析和解释，源于其使用参数少，且无需场源物性的先验信息即可有效地圈定地质体的范围，因而具有较强的适应性和灵活性。

在众多欧拉解中如何剔除缪解和确定优解一直是困扰人们的一个难题，如欧拉解在一处大量聚集，而在它处呈零星分布，传统欧拉反褶积图示方法上两者几乎没有差别。为此，许多学者将聚类分析引入至欧拉解的后处理，利用欧拉解的聚集特性来标识地质体的空间信息。如在算法Stavrev的基础上采用非聚类和聚类的方式对欧拉解进行统计分析，但其过度依赖于人工确定地下异常源的个数。

发明内容

鉴于以上所述现有技术的缺点，本发明的目的在于提供一种基于B样条函数的欧拉解集概率密度成像方法，采用基于规范化B样条函数的非参数概率密度估计方法对欧拉解集进行处理，在无先验信息的情况下，以各个欧拉解间的相似性及其集聚程度，作为计算欧拉解概率密度值的依据，进而区分由欧拉解组成的多个簇，以实现利用欧拉解标识多个相邻异常源的目的，算法验证及模型试验表明了B样条密度估计的可靠性，能够解决现有技术中存在的技术问题。

为实现上述目的及其他相关目的，本发明提供一种基于B样条函数的欧拉解集概率密度成像方法，包括有：

确定待测区域范围；测量重力矢量数据及重力梯度张量数据；或测量重力数据，通过离散余弦变换或傅里叶变换将重力数据换算到重力矢量数据及重力梯度张量数据；

根据所述重力矢量数据或所述重力梯度张量数据构建三维重力梯度张量欧拉反褶积方程；

确定滑动窗口大小，并利用滑动窗口内的数据构建线性方程组；

利用奇异值分解所述线性方程组，获取所述线性方程组的欧拉解，并输出对应的欧拉解集；其中，所述欧拉解包括异常源空间位置、异常源类型和奇异值分解误差；

将欧拉解解集拆分成不同维度的组合，并根据所述欧拉解解集组合估算估计区间的上/下界及估计网格带宽，以及确定所述估计网格的大小；

计算并存储B样条函数，并将欧拉解数据作为独立同分布采样投影至所述估计网格，计算预设权重以及所述B样条函数的系数；

对所述B样条函数以及对应的B样条函数的系数进行卷积，获取欧拉解概率密度分布结果，并根据所述欧拉解概率密度分布结果确定不同维度数据的异常源空间位置及类型。

可选地，所述根据所述重力矢量数据或所述重力梯度张量数据构建三维重力梯度张量欧拉反褶积方程，包括：

式中，B_α为异常背景场，g_α为重力梯度，g_αβ为重力梯度张量分量，{α,β∈x,y,z}；(x₀,y₀,z₀)为观测点位置坐标；(x,y,z)为待求解位场异常源位置；构造指数N是位场异常随场源深度衰减变化陡缓的量度。

可选地，所述确定滑动窗口大小，并利用滑动窗口内的数据构建线性方程组，以及利用奇异值分解所述线性方程组，获取所述线性方程组的欧拉解，并输出对应的欧拉解集，包括：

确定滑动窗口大小为Wx×Wy；

将滑窗中的Wx×Wy个测点数据代入(1)～(3)构建线性方程组，如下：

Am＝b (4)；

式中，

n_w＝w_x×w_y；

m＝[m¹,...,mⁱ,...]^T为欧拉反褶积解集；

第i个欧拉反褶积解可以表示为mⁱ＝[xⁱ,yⁱ,zⁱ,Nⁱ]^T；

利用奇异值分解求取异常源空间位置[x，y，z]及构造指数N；

逐一移动滑动窗口，构建线性方程组，并获得相应滑动窗口的欧拉解mⁱ，直至滑动窗口覆盖整个测区，构建欧拉解集{mⁱ}。

可选地，将样本Xⁱ投影至大小为λ＝[λ₁,...,λ_j,...,λ_d]^T且带宽为h＝[h₁,...,h_j,...,h_d]^T的估计网格χ，则其各维度方向的估计点向量分别为χ₁,...,χ_j,...,χ_d；

利用c次均匀B样条函数

表达d维B样条概率密度估计

则有：

由B样条密度估计的定义有：

其中，第j维的带宽h_j可以写为：

式中，

为第j维估计区间的上界，

为第j维估计区间的下界，

为第j维数据的最大值；

为第j维数据的最小值。

可选地，设定第一个网格点

最后一个网格点

λ_j为第j维的网格数；则有：

其中，δ_j为第j维数据的标准差，γ_j是一个与维度d和物理内存大小相关的经验系数，n为样本数，

为向下取整运算符，μ_j为第j维数据X_j算数平均值。

可选地，对公式(5)采用基于规范化的三次B样条基函数，则有：

式中，

为三次均匀规范B样条基函数；

系数

可以表示为：

可选地，对于各维度方向的带宽相同的情况，式(10)可以改写为：

可选地，在估计网格上，定义在j维度方向存在边

并将其简写为

若Xⁱ位于由边

构成的超矩形中，则利用

在

上的投影表示u_j，即：

其中，

为向下取整运算符，

为第j维数据X_j的第i个值。

如上所述，本发明提供一种基于B样条函数的欧拉解集概率密度成像方法，具有以下有益效果：

1、在处理如欧拉解在一处大量聚集，而在它处呈零星分布两种情况，传统欧拉反褶积图示方法于两者几乎没有差别，而本发明可通过欧拉解概率密度差异，剔除欧拉解缪解、确定欧拉解优解、区分欧拉解簇，以标识异常源；

2、当观测数据样本过大或估计网格过大时，传统B样条概率密度估计方法存在内存消耗过大，遍历估计网格耗时过多的问题，本发明可不构建大小为λ＝[λ₁,...,λ_j,...,λ_d]^T的估计网格χ；避免遍历估计网格存在多次重复计算B样条基函数，计算并存储B³(x₁),...,B³(x_j),...,B³(x_d)，通过计算

的卷积；

3、在某一维度方向的估计网格χ_j采用等间距的估计点，可以将观测数据样本{X_j}整体投影至χ_j，从而易于构建矢量化(于Matlab而言)/高效并行的算法(于C++等编程语言而言)，继而快速获得观测数据样本于估计网格的权重；

4、采用一类基于规范化B样条函数的非参数概率密度估计方法对欧拉解集进行处理，在无先验信息的情况下，以各个欧拉解间的相似性及其集聚程度，作为计算欧拉解概率密度值的依据，进而区分由欧拉解组成的多个簇。

附图说明

图1为一实施例提供的基于B样条函数的欧拉解概率密度成像的方法流程示意图；

图2为一实施例提供的立方体重力梯度张量正演图；

图3为另一实施例提供的立方体重力梯度张量正演图；

图4为一实施例提供的立方体全张量欧拉反褶积散点图；

图5为另一实施例提供的立方体全张量欧拉反褶积散点图；

图6为又一实施例提供的立方体全张量欧拉反褶积散点图；

图7为一实施例提供的一维B样条密度估计的概率密度曲线图；

图8为另一实施例提供的一维B样条密度估计的概率密度曲线图；

图9为一实施例提供的二维欧拉解集{y,x}的二维B样条密度估计的概率密度分布图；

图10为一实施例提供的二维欧拉解集{z,x}的二维B样条密度估计的概率密度分布图；

图11为一实施例提供的二维欧拉解集{y,z}的二维B样条密度估计的概率密度分布图；

图12为一实施例提供的二维欧拉解集{N,x}的二维B样条密度估计的概率密度分布图；

图13为一实施例提供的二维欧拉解集{N,y}的二维B样条密度估计的概率密度分布图；

图14为一实施例提供的二维欧拉解集{N,z}的二维B样条密度估计的概率密度分布图；

图15为一实施例提供的立方体三维欧拉解集{x,y,z}的三维B样条密度估计的概率密度等值面图；

图16为一实施例提供的立方体三维欧拉解集{x,y,N}的三维B样条密度估计的概率密度等值面图；

图17为一实施例提供的立方体三维欧拉解集{x,z,N}的三维B样条密度估计的概率密度等值面图；

图18为一实施例提供的立方体三维欧拉解集{y,z,N}的三维B样条密度估计的概率密度等值面图。

具体实施方式

以下通过特定的具体实例说明本发明的实施方式，本领域技术人员可由本说明书所揭露的内容轻易地了解本发明的其他优点与功效。本发明还可以通过另外不同的具体实施方式加以实施或应用，本说明书中的各项细节也可以基于不同观点与应用，在没有背离本发明的精神下进行各种修饰或改变。需说明的是，在不冲突的情况下，以下实施例及实施例中的特征可以相互组合。

需要说明的是，以下实施例中所提供的图示仅以示意方式说明本发明的基本构想，遂图式中仅显示与本发明中有关的组件而非按照实际实施时的组件数目、形状及尺寸绘制，其实际实施时各组件的型态、数量及比例可为一种随意的改变，且其组件布局型态也可能更为复杂。

请参阅图1至图18所示，本发明提供一种基于B样条函数的欧拉解概率密度成像方法，包括以下步骤：

确定待测区域范围；测量重力矢量数据及重力梯度张量数据；或测量重力数据，通过离散余弦变换或傅里叶变换将重力数据换算到重力矢量数据及重力梯度张量数据；

根据所述重力矢量数据或所述重力梯度张量数据构建三维重力梯度张量欧拉反褶积方程；

确定滑动窗口大小，并利用滑动窗口内的数据构建线性方程组；

利用奇异值分解所述线性方程组，获取所述线性方程组的欧拉解，并输出对应的欧拉解集；其中，所述欧拉解包括异常源空间位置、异常源类型和奇异值分解误差；

将欧拉解解集拆分成不同维度的组合，并根据所述欧拉解解集估算估计区间的上/下界及估计网格带宽，以及确定所述估计网格的大小；

计算并存储B样条函数，并将欧拉解数据作为独立同分布采样投影至所述估计网格，计算预设权重以及所述B样条函数的系数；

对所述B样条函数以及对应的B样条函数的系数进行卷积，获取欧拉解概率密度分布结果，并根据所述欧拉解概率密度分布结果确定不同维度数据的异常源空间位置及类型。

本方法可通过欧拉解概率密度差异，剔除欧拉解缪解、确定欧拉解优解、区分欧拉解簇，以标识异常源；本方法可不构建大小为λ＝[λ₁,...,λ_j,...,λ_d]^T的估计网格χ，避免遍历估计网格存在多次重复计算B样条基函数，计算并存储B³(x₁),...,B³(x_j),...,B³(x_d)，通过计算

的卷积；本方法在某一维度方向的估计网格χ_j采用等间距的估计点，可以将观测数据样本{X_j}整体投影至χ_j，从而易于构建矢量化(于Matlab而言)/高效并行的算法(于C++等编程语言而言)，继而快速获得观测数据样本于估计网格的权重；本方法采用一类基于规范化B样条函数的非参数概率密度估计方法对欧拉解集进行处理，在无先验信息的情况下，以各个欧拉解间的相似性及其集聚程度，作为计算欧拉解概率密度值的依据，进而区分由欧拉解组成的多个簇。

根据上述记载，本发明提供一种基于规范化B样条函数重力/重力张量欧拉解概率密度成像方法，包括以下步骤：

1)在传统的欧拉反褶积方法中，三维欧拉方程可表示为：

式中，T＝T(x,y,z)为位场函数；(x₀,y₀,z₀)为观测点位置坐标；(x,y,z)为待求异常源的空间位置；N为异常幅值随距离变化的衰减率，特定的地质体对应特定的衰减率，即构造指数。重磁场的一阶导数、二阶导数及延拓后的异常同样满足欧拉方程。

为消除区域背景场的影响，通常引入一个代表区域背景异常的参数B，(1-1)式可重写为：

其中，ε为失真率，当且仅当N＝0是存在。

将重力梯度张量三个方向的分量代入式(2-1)即得到三维重力梯度张量欧拉反褶积方程：

式中，α,β＝x,y,z，B_α为异常背景场，g_α为重力梯度，g_αβ为重力梯度张量分量。

2)利用一个大小为w＝w_x＝w_y的正方形滑动窗口逐步覆盖整个测区，每个滑动窗口可以由方程(3-1)-(5-1)组成的如下矩阵形式方程组：

Am＝b (6-1)；

式中，

n_w＝w_x×w_y，

m＝[m¹,...,mⁱ,...]^T为欧拉反褶积解集，第i个欧拉反褶积解可以表示为mⁱ＝[xⁱ,yⁱ,zⁱ,Nⁱ]^T，即为待求异常源空间位置及构造指数。

3)针对传统欧拉反褶积优化策略仅利用欧拉反褶积的构造指数，难于表征欧拉反褶积解聚集程度，亦难于判定单个欧拉解的可信度，因而难于区分异常及各异常之间存在的伪迹，为此引入B样条密度估计算法对欧拉反褶积解集的数据组合{x}，{x,y}，{x,y,z}，{x,y,z,N}等分别进行估计。

设X¹,...,Xⁱ,...,Xⁿ是d维概率密度f函数的独立同分布随机样本，第i个样本可以表示为

特别地，在本文中可以是欧拉解xⁱ,yⁱ,zⁱ,Nⁱ的不同组合，如xⁱ，[xⁱ,yⁱ]^T，[xⁱ,yⁱ,zⁱ]^T，[xⁱ,yⁱ,zⁱ,Nⁱ]^T等，1#d 4。将样本Xⁱ投影至大小为λ＝[λ₁,...,λ_j,...,λ_d]^T且带宽为h＝[h₁,...,h_j,...,h_d]^T的估计网格x，其各维度方向的估计点向量分别为x₁,...,x_j,...,x_d，则可以利用c次均匀B样条函数

表达d维B样条概率密度估计

即：

由B样条密度估计的定义，有：

其中，第j维的带宽h_j可以写为：

式中，

和

为第j维估计区间的上下界，

分别为第j维数据的最大值和最小值。一般第一个网格点

最后一个网格点

λ_j为第j维的网格数，由公式(10-1)计算：

其中δ_j为第j维数据的标准差，γ_j是一个与维度d和物理内存大小相关的经验系数，n为样本数，

为向下取整运算符。

实际应用中一般采用次数较低的B样条基函数，采用基于规范化的三次(四阶)B样条基函数，于是公式(7-1)可写为：

式中，

为三次均匀规范B样条基函数。

对于各维度方向的带宽相同的情况，式(12-1)可以改写为：

4)就对估计网格上每一点进行遍历计算式(13-1)而言，存在计算量过大且内存消耗过多的问题，尤其当估计网格尺寸过大时，这种瓶颈现象尤为严重。为此，对于各维度方向的带宽相同的情况，避免内存消耗过大，不构建大小为λ＝[λ₁,...,λ_j,...,λ_d]^T的估计网格x；避免遍历估计网格存在多次重复计算b样条基函数，计算并存储B³(x₁),...,B³(x_j),...,B³(x_d)，通过计算

的卷积，快速实现式(13-1)的概率密度值获取。

5)式(7-1)中系数

可以表示为：

当观测数据点较大时，难于估算

为计算式(12)中权重u，定义在j维度方向存在边

并将其简写为

假定Xⁱ恰位于由边

构成的超矩形中，则可以利用

在

上的投影表示u_j，即：

其中，

为向下取整运算符。当观测数据点较大时，可以利用式(15-1)将观测样本某一维度信息{X_j}直接投影至估计网格x_j，从而可以于Matlab计算环境实现矢量化或于C/C++等计算环境实现并行计算权重u。

根据上述记载，在一示例性实施例中，本发明结合具体实施例进行说明，按照图1所示的流程示意图建立理论模型，在本申请实施例中，模型的立方体大小1000m×1000m×1000m，质心为(-1000，-2000，1500)，剩余密度均为0.36g/cm³，观测高度为50m，测区范围为x：-10000m～10000m，y：-10000m～10000m；测网大小为100m×100m。由理论模型产生的地面重力数据共有200×200＝40000个采样点。

根据上述记载，本发明针对这个模型实例进行详细说明，如下：

a)确定三维欧拉方程，通过测量或位场转换获得重力矢量/重力梯度张量值，如图2和图3所示。

b)确定滑动窗口大小w＝w_x＝w_y，利用该正方形滑动窗口逐步覆盖整个测区，每个滑动窗口可以由方程(3-1)-(5-1)组成的如下矩阵形式方程组：

Am＝b；

式中，

n_w＝w_x×w_y，

m＝[m¹,...,mⁱ,...]^T为欧拉反褶积解集，第i个欧拉反褶积解可以表示为mⁱ＝[xⁱ,yⁱ,zⁱ,Nⁱ]^T，即为待求异常源空间位置及构造指数，如图4至图6所示。

c)将不同欧拉反褶积解集的数据组合{x}，{x,y}，{x,y,z}，{x,y,z,N}等，分别对不同数据组合进行估计，确定估计带宽h＝[h₁,...,h_j,...,h_d]^T，估计网格大小λ＝[λ₁,...,λ_j,...,λ_d]^T，d为数据组合的维度。

d)为此，对于各维度方向的带宽相同的情况，避免内存消耗过大，不构建大小为λ＝[λ₁,...,λ_j,...,λ_d]^T的估计网格χ；避免遍历估计网格存在多次重复计算B样条基函数，计算并存储B³(χ₁),...,B³(χ_j),...,B³(χ_d)；

e)当观测数据点较大时，难于估算

利用式(15-1)将观测样本某一维度信息{X_j}直接投影至估计网格χ_j，从而可以于Matlab计算环境实现矢量化或于C/C++等计算环境实现并行计算权重u。

f)通过计算

的卷积，快速实现式(13-1)的概率密度值获取。

g)对式(15-1)的计算获得概率密度结果，以概率密度曲线图(即图7和图8)、概率密度分布图(及图9至图14)或概率密度等值面图(即图15至图18)显示，并对这些结果进行解释，剔除欧拉解集中的缪解，分离欧拉解所构成的簇，以标识相邻异常源。

综上所述，本发明提供一种基于B样条函数的欧拉解集概率密度成像方法，首先通过测量或位场变换获得重力矢量及重力梯度张量数据，构建三维重力/重力梯度张量欧拉反褶积方程，逐一移动滑动窗口，获得一系列含异常源空间位置[x,y,z,N]欧拉解集；针对传统欧拉反褶积优化策略仅利用欧拉反褶积的构造指数，抑或不同大小“圈”难于表征欧拉反褶积解的聚集程度，亦难于判定单个欧拉解的可信度，因而难于区分异常及各异常之间存在的伪迹。为此，本发明以各个欧拉解间的相似性及其集聚程度，作为欧拉解分类及分离异常的标准，减少了欧拉解缪解处理过程，优化了确定欧拉解优解的机制，实现了欧拉解簇的分离，进而有效地标识地下异常源。本发明可通过欧拉解概率密度差异，剔除欧拉解缪解、确定欧拉解优解、区分欧拉解簇，以标识异常源；本发明可不构建大小为λ＝[λ₁,...,λ_j,...,λ_d]^T的估计网格χ，避免遍历估计网格存在多次重复计算B样条基函数，计算并存储B³(x₁),...,B³(x_j),...,B³(x_d)，通过计算

的卷积；本发明在某一维度方向的估计网格χ_j采用等间距的估计点，可以将观测数据样本{X_j}整体投影至χ_j，从而易于构建矢量化(于Matlab而言)/高效并行的算法(于C++等编程语言而言)，从而快速获得观测数据样本于估计网格的权重；本发明采用一类基于规范化B样条函数的非参数概率密度估计方法对欧拉解集进行处理，在无先验信息的情况下，以各个欧拉解间的相似性及其集聚程度，作为计算欧拉解概率密度值的依据，进而区分由欧拉解组成的多个簇，以实现利用欧拉解标识多个相邻异常源的目的，算法验证及模型试验表明了B样条密度估计的可靠性。

上述实施例仅例示性说明本发明的原理及其功效，而非用于限制本发明。任何熟悉此技术的人士皆可在不违背本发明的精神及范畴下，对上述实施例进行修饰或改变。因此，举凡所属技术领域中具有通常知识者在未脱离本发明所揭示的精神与技术思想下所完成的一切等效修饰或改变，仍应由本发明的权利要求所涵盖。

Claims (8)

1.一种基于B样条函数的欧拉解集概率密度成像方法，其特征在于，包括有：

确定待测区域范围；测量重力矢量数据及重力梯度张量数据；或测量重力数据，通过离散余弦变换或傅里叶变换将重力数据换算到重力矢量数据及重力梯度张量数据；

根据所述重力矢量数据或所述重力梯度张量数据构建三维重力梯度张量欧拉反褶积方程；

确定滑动窗口大小，并利用滑动窗口内的数据构建线性方程组；

利用奇异值分解所述线性方程组，获取所述线性方程组的欧拉解，并输出对应的欧拉解集；其中，所述欧拉解包括异常源空间位置、异常源类型和奇异值分解误差；

将欧拉解解集拆分成不同维度的组合，并根据所述欧拉解解集组合估算估计区间的上/下界及估计网格带宽，以及确定所述估计网格的大小；

计算并存储B样条函数，并将欧拉解数据组合作为独立同分布采样投影至所述估计网格，计算预设权重以及所述B样条函数的系数；

对所述B样条函数以及对应的B样条函数的系数进行卷积，获取欧拉解概率密度分布结果，并根据所述欧拉解概率密度分布结果确定不同维度数据的异常源空间位置及类型。

2.根据权利要求1所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，所述根据所述重力矢量数据或所述重力梯度张量数据构建三维重力梯度张量欧拉反褶积方程，包括：

式中，B_α为异常背景场，g_α为重力梯度，g_αβ为重力梯度张量分量，{α,β∈x,y,z}；(x₀,y₀,z₀)为观测点位置坐标；(x,y,z)为待求解位场异常源位置；构造指数N是位场异常随场源深度衰减变化陡缓的量度。

3.根据权利要求2所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，所述确定滑动窗口大小，并利用滑动窗口内的数据构建线性方程组，以及利用奇异值分解所述线性方程组，获取所述线性方程组的欧拉解，并输出对应的欧拉解集，包括：

确定滑动窗口大小为Wx×Wy；

将滑窗中的Wx×Wy个测点数据代入(1)～(3)构建线性方程组，如下：

Am＝b (4)；

式中，

n_w＝w_x×w_y；

m＝[m¹,...,mⁱ,...]^T为欧拉反褶积解集；

第i个欧拉反褶积解可以表示为mⁱ＝[xⁱ,yⁱ,zⁱ,Nⁱ]^T；

利用奇异值分解求取异常源空间位置[x,y,z]及构造指数N；

逐一移动滑动窗口，构建线性方程组，并获得相应滑动窗口的欧拉解mⁱ，直至滑动窗口覆盖整个测区，构建欧拉解集{mⁱ}。

4.根据权利要求3所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，将样本Xⁱ投影至大小为λ＝[λ₁,...,λ_j,...,λ_d]^T且带宽为h＝[h₁,...,h_j,...,h_d]^T的估计网格χ，则其各维度方向的估计点向量分别为χ₁,...,χ_j,...,χ_d；

利用c次均匀B样条函数

表达d维B样条概率密度估计

则有：

由B样条密度估计的定义有：

其中，第j维的带宽h_j可以写为：

式中，

为第j维估计区间的上界，

为第j维估计区间的下界，

为第j维数据的最大值；

为第j维数据的最小值。

5.根据权利要求4所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，设定第一个网格点

最后一个网格点

λ_j为第j维的网格数；则有：

其中，δ_j为第j维数据的标准差，γ_j是一个与维度d和物理内存大小相关的经验系数，n为样本数，

为向下取整运算符，μ_j为第j维数据X_j算数平均值。

6.根据权利要求4或5所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，对公式(5)采用基于规范化的三次B样条基函数，则有：

式中，

为三次均匀规范B样条基函数；

系数

可以表示为：

7.根据权利要求6所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，对于各维度方向的带宽相同的情况，式(10)可以改写为：

8.根据权利要求6所述的基于B样条函数的欧拉解集概率密度成像方法，其特征在于，在估计网格上，定义在j维度方向存在边

并将其简写为

若Xⁱ位于由边

构成的超矩形中，则利用

在

上的投影表示u_j，即：

其中，

为向下取整运算符，

为第j维数据X_j的第i个值。

Images