本申请涉及科学数据可视化技术领域，尤其涉及一种基于大规模流线的矢量场实时仿真方法及系统。
背景技术
我们从维度上可以分为二维、三维以及多维(时间)的矢量场数据可视化，主要的应用领域包括在气象数据和计算机中的物理流体动力学速度场的可视化。
基于矢量场中每个采样点的数据是一个向量，其主要目标为展示其导向趋势信息、场中的模式、识别关键特征区域。而且二维或者三维流场可视化记录了水流、空气等流动过程中的方向信息，是应用最广泛、研究最深入的矢量场，因此流场的可视化是矢量场可视化的重要组成部分。这其中的实现方法包括速度场中流线、气体中的时线、液体中的脉线、固体中的迹线的静态表现方式，还有包括粒子运动、实时可视化数据真实效果的动态表现方式。这些可视化技术在实际运用过程中，能够在物理、气象、医学等诸多领域的应用中都很有成效。
目前，现有技术中通常采用实体点图标绘制技术和粒子对流方法。其中，实体点图标绘制技术是对二维矢量场中的图标法的沿用，改用三维箭头图标代替二维图标。其主要内容即通过箭头表达矢量的方向、线段长度表达矢量大小来表示矢量场信息，简单、直观、快速的实现方式，但是此方法实现不具动态性且表达方式太过单一，无法实现其矢量场的导向性并且杂乱无章，无法观测，在对其内部数据的提取过程中，进行特征数据块可视化，无法对整个矢量场有总体的趋势信息进行可视化。粒子对流方法是通过粒子在矢量场中的流动，以动画的方式表达某一时刻的矢量场数据，直观动态的展示矢量场特征。其通过在矢量场中播撒种子粒子，在每一帧绘制时，粒子沿所在位置的矢量方向移动，移动速度和矢量场成正比，粒子连续移动形成轨迹曲线，轨迹随时间逐渐淡化，形成渐变曲线，并使用gpu加速降低视角切换时重新计算矢量场分布所需时间，但是此方法在矢量场可视化中，其基于常规的cpu运算限制了其计算能力，不能够大规模提升可视化粒子数量并且其尾迹不能够大量的绘制。
因此，如何有效的进行矢量场实时仿真是一项亟待解决的问题。
技术实现要素：
有鉴于此，本申请提供了一种基于大规模流线的矢量场实时仿真方法，能够通过对原始数据的类矢量场处理和粒子运动并行化的设计，通过大规模的粒子流线来进行矢量场实时仿真，解决了现有技术中矢量场趋势表达不明显、无法大规模的粒子实现流场的问题。
本申请提供了一种基于大规模流线的矢量场实时仿真方法，所述方法包括：
将原始矢量场数据处理为高维数组；
对所述高维数组的矢量场进行粒子数量补充；
更新所述高维数组的矢量场中粒子的状态。
优选地，所述对所述高维数组的矢量场进行粒子数量补充包括：
计算所述高维数组的矢量场边界值；
基于所述高维数组的矢量场边界值确定出最大粒子数，并在第一帧生成最大粒子数个粒子；
获取当前场景粒子数；
当所述当前场景粒子数小于所述最大粒子数时，生成所述最大粒子数与所述当前场景粒子数的差值个粒子。
优选地，所述计算所述高维数组的矢量场边界值包括：
计算所述高维数组的矢量场的横、纵、垂直三个方向的格点数和格距；
基于所述格点数和格距，以空间坐标原点为基准，计算出所述高维数组的矢量场边界值。
优选地，所述更新所述高维数组的矢量场中粒子的状态包括：
分配分别用于缓存矢量场数据和粒子状态信息的两块缓存区；
分配线程，以使每个粒子的下一时刻状态计算能够分配一个线程；
基于当前场景中粒子的坐标信息，在矢量场中对粒子进行快速定位，得到粒子当前所处的网格的顶点矢量数据；
基于所述粒子当前所处的网格的顶点矢量数据进行矢量插值，得到粒子的当前速度和位置；
基于所述粒子的当前速度和位置，更新缓存所述粒子状态信息的缓存区中的粒子速度和位置。
优选地，所述更新所述高维数组的矢量场中粒子的状态还包括：
基于所述粒子的当前速度的模与所述矢量场的所有矢量模的最大值的比例计算出当前粒子的颜色值，并根据所述颜色值和配色表的对应关系更新当前粒子的颜色。
一种基于大规模流线的矢量场实时仿真系统，包括：
处理模块，用于将原始矢量场数据处理为高维数组；
补充模块，用于对所述高维数组的矢量场进行粒子数量补充；
更新模块，用于更新所述高维数组的矢量场中粒子的状态。
优选地，所述补充模块在对所述高维数组的矢量场进行粒子数量补充时具体用于：
计算所述高维数组的矢量场边界值；
基于所述高维数组的矢量场边界值确定出最大粒子数，并在第一帧生成最大粒子数个粒子；
获取当前场景粒子数；
当所述当前场景粒子数小于所述最大粒子数时，生成所述最大粒子数与所述当前场景粒子数的差值个粒子。
优选地，所述补充模块在计算所述高维数组的矢量场边界值时具体用于：
计算所述高维数组的矢量场的横、纵、垂直三个方向的格点数和格距；
基于所述格点数和格距，以空间坐标原点为基准，计算出所述高维数组的矢量场边界值。
优选地，所述更新模块在更新所述高维数组的矢量场中粒子的状态时具体用于：
分配分别用于缓存矢量场数据和粒子状态信息的两块缓存区；
分配线程，以使每个粒子的下一时刻状态计算能够分配一个线程；
基于当前场景中粒子的坐标信息，在矢量场中对粒子进行快速定位，得到粒子当前所处的网格的顶点矢量数据；
基于所述粒子当前所处的网格的顶点矢量数据进行矢量插值，得到粒子的当前速度和位置；
基于所述粒子的当前速度和位置，更新缓存所述粒子状态信息的缓存区中的粒子速度和位置。
优选地，所述更新模块在更新所述高维数组的矢量场中粒子的状态时还用于：
基于所述粒子的当前速度的模与所述矢量场的所有矢量模的最大值的比例计算出当前粒子的颜色值，并根据所述颜色值和配色表的对应关系更新当前粒子的颜色。
综上所述，本申请公开了一种基于大规模流线的矢量场实时仿真方法，当需要对矢量场进行实时仿真时，首先将原始矢量场数据处理为高维数组，然后对高维数组的矢量场进行粒子数量补充，更新高维数组的矢量场中粒子的状态。本申请能够通过对原始数据的类矢量场处理和粒子运动并行化的设计，通过大规模的粒子流线来进行矢量场实时仿真，解决了现有技术中矢量场趋势表达不明显、无法大规模的粒子实现流场的问题。
附图说明
为了更清楚地说明本申请实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本申请公开的一种基于大规模流线的矢量场实时仿真方法实施例1的流程图；
图2为本申请公开的粒子矢量场高维数组索引查找示意图；
图3为本申请公开的粒子所处位置的矢量场立方体网格示意图；
图4为本申请公开的一种基于大规模流线的矢量场实时仿真系统实施例1的结构示意图。。
具体实施方式
下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。
如图1所示，为本申请公开的一种基于大规模流线的矢量场实时仿真方法流程图，方法可以包括以下步骤：
s101、将原始矢量场数据处理为高维数组；
当需要对矢量场实时仿真时，首先对原始矢量场数据进行预处理，将原始矢量场数据处理为高维数组。即，如果原始矢量场为三维矢量场，则将三维矢量场处理为四维数组，第一维度为横向网格顶点索引u，第二维为纵向网格顶点索引v，第三维为垂直方向网格顶点索w，第四维为数据且长度为3，存储当前对应顶点的矢量值(x,y,z)。需要说明的是，原始矢量场可以为二维矢量场、三维矢量场等，二维矢量场数据的处理方法与上述举例的三维矢量场数据的处理原理相同，在此不再赘述。
s102、对高维数组的矢量场进行粒子数量补充；
当将原始矢量场数据处理为高维数组后，进一步对高维数组的矢量场进行粒子数量补充，及时补充场景中生命周期结束的粒子。
s103、更新高维数组的矢量场中粒子的状态。
然后，对高维数组的矢量场中粒子的状态进行更新，其中，粒子的状态包括粒子的速度、尾迹颜色、生命周期等，从而实现基于大规模流线的矢量场实时仿真。
综上所述，在上述实施例中，当需要对矢量场进行实时仿真时，首先将原始矢量场数据处理为高维数组，然后对高维数组的矢量场进行粒子数量补充，更新高维数组的矢量场中粒子的状态。能够通过对原始数据的类矢量场处理和粒子运动并行化的设计，通过大规模的粒子流线来进行矢量场实时仿真，解决了现有技术中矢量场趋势表达不明显、无法大规模的粒子实现流场的问题。
具体的，在上述实施例中，步骤s102对高维数组的矢量场进行粒子数量补充的其中一种实现方式可以是：计算高维数组的矢量场的横、纵、垂直三个方向的格点数和格距，根据格点数和格距，以空间坐标原点(0，0，0)为基准，计算出高维数组的矢量场边界值，基于高维数组的矢量场边界值确定出最大粒子数maxparticle，并在第一帧生成最大粒子数个粒子，并获取当前场景中存在的粒子数countparticle，当当前场景粒子数countparticle小于最大粒子数maxparticle时，生成的粒子个数为最大粒子数maxparticle与当前场景粒子数countparticle的差值，然后将生成的粒子及时补充场景中生命周期结束的粒子。
具体的，在上述实施例中，步骤s103更新高维数组的矢量场中粒子的状态的其中一种实现方式可以是：根据当前场景中粒子的坐标信息，进行粒子在矢量场中的快速定位，进而得到当前所处的网格的顶点矢量数据，然后进行矢量插值，得到当前粒子的速度和位置，计算出标量值与颜色表进行映射，继而更新粒子的速度、尾迹颜色、生命周期等。具体可以采用以下步骤：
步骤a、分配两块缓存区，矢量场数据通过缓存区与gpu(graphicsprocessingunit，图形处理器)之间进行数据交互。共两部分，一部分是缓存矢量场数据的vectorfieldbuffer缓存区，一部分是缓存粒子状态信息的particlestatebuffer缓存区。
步骤b、分配线程，主要根据gpu并行计算的线程分配机制，由当前最大粒子数maxparticle决定，通过gpu产生与当前场景最大粒子数个数相同的粒子，利用gpu的线程索引与缓存区中对应的粒子相互绑定，保证每个粒子的下一时刻状态计算能够分配一个线程。
步骤c、粒子快速定位，根据线程索引在步骤a中的粒子状态信息缓存particlestatebuffer中取出其空间坐标和矢量场在x、y、z方向上的格点位置，如图2所示计算出当前粒子在矢量场中最小顶点索引minindexu、minindexv、minindexw，然后以此minindexu+offset、minindexv+offset、minindexw+offset，其中offset∈[0,1]，计算出矢量场立方体网格的另外七个顶点索引。
步骤d、粒子速度插值，依照步骤c中的索引值minindexu、minindexv、minindexw和另外七个顶点的索引值在步骤a中矢量场数据的vectorfieldbuffer缓存区中取出粒子所处单位立方体网格的八个顶点矢量值v1–v8如图3所示，通过当前粒子的空间坐标p和八个顶点的空间坐标p1–p8的相互关系进行矢量三线性插值如公式(1)、(2)、(3)所示，计算得到当前粒子的速度v，继而计算出位置。
x方向：
z方向：
y方向：
其中，px、py、pz粒子当前位置的三维坐标，pix、piy、piz是第i个格点位置的坐标，vi是第i个格点的速度，vi,i+1是两个格点的速度插值结果，vi,i+1,i+2,i+3是四个格点的速度插值结果，v是八个格点的速度插值最终结果gridintervalx、gridintervaly、gridintervalz分别是x、y、z方向的格距。
步骤e、更新缓存粒子状态信息的particlestatebuffer缓存区中的速度和位置，在gpu计算的结果的缓存中取出其速度和位置信息，并根据取出的速度和位置信息更新粒子的速度和颜色等状态信息，完成粒子状态更新。
综上所述，对于其大规模矢量数据的仿真能够清晰的描述其矢量场的导向趋势信息、识别关键特征区域，并且通过此种大规模流线的形式还能够描述其矢量场的动态信息。在此实现方式中采用对单个粒子的并行化处理，为其分配单个线程进行计算，有效提高了其实时仿真帧率。
为了更加清楚的说明本申请公开的基于大规模流线的矢量场实时仿真方法，下面以具体的两个实例进行说明。
实例1、二维矢量场全球洋流实时仿真实例
采用全球洋流数据，即144*73网格二维矢量数据。首先进行数据处理，原始数据为netcdf数据，利用动态链接库文件解析，存储在[144,73,2]数组形式。
基于粒子系统构造矢量场。以空间坐标(0，0，0)为矢量场左上角原点，以矢量场格点数和格距计算出矢量场空间坐标边界值，此实例格距为2、横向格点数为144、纵向格点数为73，以此计算出边界顶点值为(288，0，0)、(0，0，146)、(288，0，146)。在此边界值范围内生成6*5*6*5*6*5(27000)个粒子，如果当前场景粒子数countparticle少于27000，则在边界生成(27000–countparticle)个粒子，以此保证矢量场中充足的粒子数。
更新粒子状态。分配缓冲区进行cpu与gpu的数据交互，一个是矢量场数据vectorfieldbuffer，另一个是粒子状态particlestatebuffer。
并行化设计。gpu分配6*6*6个线程组，每个线程组分配5*5*5个线程，此时我们根据线程索引从particlestatebuffer取出粒子状态信息，为每个粒子分配一个线程。
粒子快速定位。当前粒子所处位置(x，y，z)，依据格距(int)x/2、(int)(z/2)计算出在数组中的一维索引xindex、二维索引zindex,则网格的四个顶点的索引是一维索引xindex+offset、二维索引zindex+offset，offset∈[0，1]，得到粒子在矢量场中所处网格。
粒子速度插值。根据四个顶点的索引，从vectorfieldbuffer中取出矢量分量值，得到其矢量为(vectorfieldbuffer[xindex,zindex,0],0,vectorfieldbuffer[xindex,zindex,1])。依次得出四个顶点矢量v1、v2、v3、v4。采用双线性插值处理四个矢量得到当前粒子的矢量v，并以此计算出位置。
更新particlestatebuffer中的速度和位置，在gpu中取出其状态信息更新粒子的速度和颜色等状态信息，完成粒子状态更新。
实例2、三维矢量场飓风场实时仿真实例
采用经度(-100，100)，纬度(-100，100)，气压(650.0，1100.0)中19个等压面的三维矢量场数据，即200*200*19立方体网格三维矢量数据。首先进行数据处理，原始数据为netcdf数据，利用动态链接库文件解析，存储在[200,200,19,3]数组形式。
基于粒子系统构造矢量场。以空间坐标(0，0，0)为矢量场左上角原点，以矢量场格点数和格距计算出矢量场空间坐标边界值，此实例x，z格距为2、y格距为10，横向格点数为200、纵向格点数为200、垂直格点数19，以此计算出边界顶点值为(400，0，0)、(0，0，400)、(0，190，0)三个以及其他5个顶点.在此边界值范围内生成6*5*6*5*6*5(27000)个粒子，如果当前场景粒子数countparticle少于27000，则在边界内生成(27000–countparticle)个粒子，以此保证矢量场中充足的粒子数。
更新粒子状态。分配缓冲区进行cpu与gpu的数据交互，一个是矢量场数据vectorfieldbuffer，另一个是粒子状态particlestatebuffer。
并行化设计。gpu分配6*6*6个线程组，每个线程组分配5*5*5个线程，此时我们根据线程索引从particlestatebuffer取出粒子状态信息，为每个粒子分配一个线程。
粒子快速定位。当前粒子所处位置(x，y，z)，依据格距(int)x/2、(int)(z/2)、(int)(y/10)计算出在数组中的一维索引xindex、二维索引zindex,三维索引yindex，则立方体网格的八个顶点的索引是一维索引xindex+offset、二维索引zindex+offset、三维索引yindex+offset，offset∈[0，1]，得到粒子位置在矢量场中所处的立方体网格。
粒子速度插值。根据四个顶点的索引，从vectorfieldbuffer中取出矢量分量值，得到：(vectorfieldbuffer[xindex,zindex,yindex,0]，vectorfieldbuffer[xindex,zindex,yindex,1],vectorfieldbuffer[xindex,zindex,yindex,2])，依次得出8个顶点矢量v1～v8。采用三线性插值处理8个矢量得到当前粒子的矢量v，并以此计算出位置。
更新particlestatebuffer中的速度和位置，在gpu中取出其状态信息更新粒子的速度和颜色等状态信息，完成粒子状态更新。
如图4所示，为本申请公开的一种基于大规模流线的矢量场实时仿真系统的结构示意图，系统可以包括：
处理模块401，用于将原始矢量场数据处理为高维数组；
当需要对矢量场实时仿真时，首先对原始矢量场数据进行预处理，将原始矢量场数据处理为高维数组。即，如果原始矢量场为三维矢量场，则将三维矢量场处理为四维数组，第一维度为横向网格顶点索引u，第二维为纵向网格顶点索引v，第三维为垂直方向网格顶点索w，第四维为数据且长度为3，存储当前对应顶点的矢量值(x,y,z)。需要说明的是，原始矢量场可以为二维矢量场、三维矢量场等，二维矢量场数据的处理方法与上述举例的三维矢量场数据的处理原理相同，在此不再赘述。
补充模块402，用于对高维数组的矢量场进行粒子数量补充；
当将原始矢量场数据处理为高维数组后，进一步对高维数组的矢量场进行粒子数量补充，及时补充场景中生命周期结束的粒子。
更新模块403，用于更新高维数组的矢量场中粒子的状态。
然后，对高维数组的矢量场中粒子的状态进行更新，其中，粒子的状态包括粒子的速度、尾迹颜色、生命周期等，从而实现基于大规模流线的矢量场实时仿真。
综上所述，在上述实施例中，当需要对矢量场进行实时仿真时，首先将原始矢量场数据处理为高维数组，然后对高维数组的矢量场进行粒子数量补充，更新高维数组的矢量场中粒子的状态。本申请能够通过对原始数据的类矢量场处理和粒子运动并行化的设计，通过大规模的粒子流线来进行矢量场实时仿真，解决了现有技术中矢量场趋势表达不明显、无法大规模的粒子实现流场的问题。
具体的，在上述实施例中，补充模块402对高维数组的矢量场进行粒子数量补充的其中一种实现方式可以是：计算高维数组的矢量场的横、纵、垂直三个方向的格点数和格距，根据格点数和格距，以空间坐标原点(0，0，0)为基准，计算出高维数组的矢量场边界值，基于高维数组的矢量场边界值确定出最大粒子数maxparticle，并在第一帧生成最大粒子数个粒子，并获取当前场景中存在的粒子数countparticle，当当前场景粒子数countparticle小于最大粒子数maxparticle时，生成的粒子个数为最大粒子数maxparticle与当前场景粒子数countparticle的差值，然后将生成的粒子及时补充场景中生命周期结束的粒子。
具体的，在上述实施例中，步骤更新模块403更新高维数组的矢量场中粒子的状态的其中一种实现方式可以是：根据当前场景中粒子的坐标信息，进行粒子在矢量场中的快速定位，进而得到当前所处的网格的顶点矢量数据，然后进行矢量插值，得到当前粒子的速度和位置，计算出标量值与颜色表进行映射，继而更新粒子的速度、尾迹颜色、生命周期等。具体可以采用以下步骤：
步骤a、分配两块缓存区，矢量场数据通过缓存区与gpu(graphicsprocessingunit，图形处理器)之间进行数据交互。共两部分，一部分是缓存矢量场数据的vectorfieldbuffer缓存区，一部分是缓存粒子状态信息的particlestatebuffer缓存区。
步骤b、分配线程，主要根据gpu并行计算的线程分配机制，由当前最大粒子数maxparticle决定，通过gpu产生与当前场景最大粒子数个数相同的粒子，利用gpu的线程索引与缓存区中对应的粒子相互绑定，保证每个粒子的下一时刻状态计算能够分配一个线程。
步骤c、粒子快速定位，根据线程索引在步骤a中的粒子状态信息缓存particlestatebuffer中取出其空间坐标和矢量场在x、y、z方向上的格点位置，如图2所示计算出当前粒子在矢量场中最小顶点索引minindexu、minindexv、minindexw，然后以此minindexu+offset、minindexv+offset、minindexw+offset，其中offset∈[0,1]，计算出矢量场立方体网格的另外七个顶点索引。
步骤d、粒子速度插值，依照步骤c中的索引值minindexu、minindexv、minindexw和另外七个顶点的索引值在步骤a中矢量场数据的vectorfieldbuffer缓存区中取出粒子所处单位立方体网格的八个顶点矢量值v1–v8如图3所示，通过当前粒子的空间坐标p和八个顶点的空间坐标p1–p8的相互关系进行矢量三线性插值如公式(1)、(2)、(3)所示，计算得到当前粒子的速度v，继而计算出位置。
x方向：
z方向：
y方向：
其中，px、py、pz粒子当前位置的三维坐标，pix、piy、piz是第i个格点位置的坐标，vi是第i个格点的速度，vi,i+1是两个格点的速度插值结果，vi,i+1,i+2,i+3是四个格点的速度插值结果，v是八个格点的速度插值最终结果，gridintervalx、gridintervaly、gridintervalz分别是x、y、z方向的格距。
步骤e、更新缓存粒子状态信息的particlestatebuffer缓存区中的速度和位置，在gpu中取出其状态信息更新粒子的速度和颜色等状态信息，完成粒子状态更新。
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。
专业人员还可以进一步意识到，结合本文中所公开的实施例描述的各示例的单元及算法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，为了清楚地说明硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本申请的范围。
结合本文中所公开的实施例描述的方法或算法的步骤可以直接用硬件、处理器执行的软件模块，或者二者的结合来实施。软件模块可以置于随机存储器(ram)、内存、只读存储器(rom)、电可编程rom、电可擦除可编程rom、寄存器、硬盘、可移动磁盘、cd-rom、或技术领域内所公知的任意其它形式的存储介质中。
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本申请。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本申请的精神或范围的情况下，在其它实施例中实现。因此，本申请将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。