技术领域
[0001] 本
发明涉及
计算机图形学技术领域，特别是涉及基于VR及手势交互技术的全球气候矢量场
数据可视化方法。
背景技术
[0002] 全球气候矢量数据可视化方法为传统气象矢量数据的显示提供了另一种直观清晰的表达方式，为相关人员提供了对气候数据进行分析的便利工具。
[0003] 气候数据一般有
风、洋流等矢量场数据，针对三维矢量场数据，一般提供箭矢图、流图、线积分卷积等传统矢量场可视化方法，常用的还有图标法、几何法以及纹理法等。对于预报员而言，通过这些复杂的图文来获得对整个大气情况的理解是很难的。三维大气数据可视化方法，尤其是矢量数据的可视化，有设计粒子
对流方法用动画表达矢量场数据，Mei HH,Chen HD,Zhao X,Liu HN,Zhu B,Chen W.Visualization system of 3D global scale meteorological data.Ruan Jian Xue Bao/Journal of Software,2016,27(5):1140-1150(inChinese)。本文在粒子对流方法模拟矢量风数据的同时，并且结合当前火热的VR，通过使用oculus rift DK2设备在将气象数据可视化在沉浸式
虚拟环境中，并通过预先设定leap motion中手势动作的定义，增强了交互性。
[0004] Leap Motion是面向PC以及Mac的体感
控制器制造公司Leap于2013年2月27日发布的体感控制器。支持Windows 7、Windows 8以及Mac OS X 10.7及10.8，Leap Motion控制器可追踪全部10只
手指，
精度高达1/100毫米，以超过每秒200
帧的速度追踪您的手部移动。
[0005] 近年来，随着
虚拟现实(VR)被广泛接受，其描述为用户被三维计算机生成的表示包围的一种体验，用户能够在其中移动，从不同的
角度查看它，并有能
力重塑它。除了直观的交互，虚拟环境也可以允许用户的更大程度的移动。通过使用VR，数据可以变得更容易理解，记住和参考。用户操作或观察，同时沉浸，可以形成用户的叙事，从而增强召回的能力数据。用户命令和输入可以立即
修改虚拟世界并描绘每个变化的影响，从而创造一个动态模型。实时交互是一个VR的关键特征，例如气象学家可以使用它来形成对气象数据的动态分析，比单独使用图表更快。VR对于进一步提高可视化数据的可读性有显著效果，大量可视化研究中采用VR，Irene Katsouri,Aimilia Tzanavari,Kyriakos Herakleous,and Charalambos Poullis.2015.Visualizing and assessing hypotheses for marine archaeology in a VR CAVE environment.ACM J.Comput.Cult.Herit.8,2,Article10(March2015),18 pages开发了一个使用VR CAVE的沉浸式3D可视化应用程序，旨在使研究人员能够挖掘这个古代沉船所提供的丰富的信息,也有对考古数据在虚拟现实环境中可视化的研究，但在虚拟现实环境中对三维大气数据可视化的研究鲜有。
发明内容
[0006] 针对
现有技术中存在的上述问题，本发明的目的在于提供基于VR及手势交互技术的全球气候矢量场数据可视化方法。
[0007] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于该方法包括以下步骤：
[0008] 步骤1)搭建3D环境，设计一个场景并根据球面体绘制方法构造一个球形代理几何；
[0009] 步骤2)获取由Natural Earth提供的地图数据，并对地图数据进行预处理，再绘制到步骤1)设计场景中的球体上；
[0010] 步骤3)获取由美国国家气象局操作的全球预报系统生成的天气数据，并对天气数据进行预处理，预测每天生成四次，可从NOMADS下载，这些文件是GRIB2格式，包含超过300条记录；
[0011] 步骤4)对步骤3)中的数据进行风向快速参考，用于转换向量风分量、风速和风向；
[0012] 步骤5)对经步骤4)处理后获得转换后的向量风分量再进行地转风近似求解；
[0013] 步骤6)将步骤5)处理后的数据投影到步骤2)的球体上；
[0014] 步骤7)有限差分近似用于在步骤6)插值过程期间估计该失真，确保风粒子路径被正确
渲染，得到可视化结果；
[0015] 步骤8)交互设计，根据需求采用手势交互技术定义各种操作手势，实现全球气候矢量场数据可视化，操作手势包括
水平竖直移动、旋转、放大及缩小。
[0016] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤1)的球面体绘制方法如下：球面体绘制方法基于传统的光线投射
算法，采用球面坐标表达构造一个球形代理几何，并以球面
坐标系的三个空间维度坐标作为纹理坐标，以归一化后的三维体数据作为纹理坐标，以所述球形代理几何作为载体构造三维体数据的三维体纹理，而在绘制过程前，又将球面坐标转换回直角坐标表达，通过求解二次方程计算投射光线与代理几何的交点，具体为：
[0017] 设当前相机视点
位置为(x0,y0,z0),光线方位为单位向量(xd,yd,zd)，则得到其参数化方程：
[0018] x＝x0+xd×t
[0019] y＝y0+yd×t
[0020] z＝z0+zd×t
[0021] 由于球形代理几何中心点为(0,0,0)，可以由方程 来表达，将投射光线参数化方程带入，得到：A×t2+B×t2+C＝0，其中，
[0022] 若二次方程无解，则相机视点在球体外且投射光线与代理几何不相交；若有解，则根据解得情况计算不同的光线积分的初始点和光线积分的终止点，获得当前光线位置转化为球面坐标。
[0023] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤2)中的地图数据预处理方法如下：将所获取的地图数据为GeoJSON格式转换为D3.js适用的TopoJSON格式，即首先在processing中安装GDAL和TopoJSON，然后将GeoJSON数据文件中边界线记录次数改为一次，且将浮点数据转换为整数形式。
[0024] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤4)中的风向快速参考方法如下：由于所得到的气象风数据以角度表示，为方便计算，将角度转
化成度，通过乘以DperR(180/π＝57.29578)将角度转换为度，并通过乘以RperD(π/180＝0.01745329)将角度转换为弧度：
[0025] 地理风坐标系：Ugeo，Vgeo；
[0026] +Ugeo分量代表风向东，即西风；-Ugeo表示风向西，即东风；
[0027] +Vgeo是风向北，即南风；
[0028] 地理风向：Dirgeo是相对于真正北方的方向，其中0＝北，90＝东，180＝南，270＝西，Dirgeo＝atan2(-Ugeo，-Vgeo)*DperR＝270-(atan2(Vgeo，Ugeo)*DperR)；
[0029] 水平风速Spd：Spd＝sqrt(Ugeo*Ugeo+Vgeo*Ugeo)；要从Dirgeo和Spd转换为Ugeo和Vgeo：Ugeo＝-Spd*sin(Dirgeo*RperD)；Vgeo＝-Spd*cos(Dirgeo*RperD)。
[0030] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤5)中的地转风近似求解如下：步骤4)的风向快速参考转换了向量风分量，近似被分解成其两个水平分量U和V，U分量表示风的东西分量，V分量表示南北分量，风的U分量和V分量的两个公式如下：
[0031] U＝((Y1-Y0)×R×Tavg×ln(P0/P1))/(2×omega×sin(lat)×(D^(2)))，
[0032] V＝((X1-X0)×R×Tavg×ln(P1/P0))/(2*omega×sin(lat)×(D^(2)))
[0033] 其中：(Y1-Y0)为当前位置和另一个参考点之间的距离的南北分量的差值，以米为单位；
[0034] (X1-X0)为当前位置和另一个参考点之间的东西向分量的差值，以米为单位；
[0035] R为干燥空气的平均气体常数，为287焦
耳/千克*开氏度；
[0036] Tavg为当前位置和另一个参考点之间的平均
温度，单位为开氏温标；
[0038] P1为另一个参考点处的大气压力；
[0039] omega为地球旋转的
角速度＝7.292×10^(-5)倒数秒；
[0040] Sin(lat)为当前所在位置的纬度的sin函数；
[0041] D为当前位置与另一个参考点之间的距离，以米为单位。
[0042] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤6)中的矢量风数据投影方法如下：采用双线性内插先在x轴上执行线性插值，然后再在y轴上执行线性插值，尽管每个步骤在
采样值和位置中是线性的，但是内插作为整体不是线性的，而是在采样位置的二次方，为了获取函数f在点(x，y)上的值，已知四个点Q11＝(x1，y1)，Q12＝(x1，y2)，Q21＝(x2，y1)和Q22＝(x2，y2)的f的值，我们首先在x方向进行线性插值，产生：
[0043]
[0044]
[0045] 再通过在y方向上内插来获得期望的估计：
[0046]
[0047] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤7)中的有限差分近似如下：首先投影坐标(x，y)相对于球面坐标(lat，lon)的所有四个一阶导数＝(phi，lambda)：dx/d(phi)，dx/d(lambda)，dy/d(phi)，dy/d(lambda),在距离h＝10-5.2弧度处使用一阶中心有限差分估计,要在点(phi，lambda)获得dx/d(phi)，投影点[0048] (phi-h/2,lambda)-->(x0,y0)
[0049] (phi+h/2,lambda)-->(x1,y1)
[0050] 使用估计值
[0051] dx/d(phi)＝(x1-x0)/h
[0052] dy/d(phi)＝(y1-y0)/h
[0053] 同样，投影点
[0054] (phi,lambda-h/2)-->(x2,y2)
[0055] (phi,lambda+h/2)-->(x3,y3)
[0056] 并使用估计
[0057] dx/d(lambda)＝(x3-x2)/h
[0058] dy/d(lambda)＝(y3-y2)/h，
[0059] 从上述衍
生物连同斯奈德公式，可以获得在(phi，lambda)的Tissot Indicatrix的轴的长度及其方向；根据地图的大小来确定比例因子，找到地图上的典型TI的大小，并缩放，使得那些TI将大约为地图的6％；以相同的量重新缩放所有的TI，因此它们可以被比较，并且将每个重新缩放在其自己的中心附近(其通过第五投影(phi，lambda)→(x，y)。
[0060] 所述的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其特征在于步骤8)中的手势交互技术如下：
[0061] 1)leap motion首先从双目摄像头采集操作者手势动作的左右视觉图像；
[0062] 2)经过立体标定后获取经过校准的立体图像对后，进行立体匹配，获得
视差图像；
[0063] 3)再利用摄像机的内参数及外参数进行三角计算获取
深度图像，然后对左或右视觉图像使用手势分割算法处理，分割出人手所在的初始位置信息；
[0064] 4)预先设定
阈值，当人手左右水平移动超过阈值，即对地球进行水平移动；当人手竖直移动超过阈值，即对地球进行竖直移动；人手五指握拳，实现缩小操作；人手由拳展开，实现放大操作；手掌旋转超过阈值角度，即对地球进行旋转。
[0065] 通过采用上述技术，与现有技术相比，本发明的有益效果如下:
[0066] 本发明针对现有的气象数据可视化的局限性，依据设计方案，通过采用leap motion增强了沉浸式环境下的交互性，实现了手势操作的相关交互设定，结合VR技术能够利用计算机生成较为真实的虚拟环境，借助于processing的可视化三维技术，更加逼真地模拟了全球
气候变化，突破了传统的二维显示的局限性，增加了模型和数据的
感知性和真实性，获得更高层次的信息细节，使用者不仅能够真实的感觉到研究的对象近在眼前，而且还能通过手势操作研究对象，从而构建了一个直观的空间数据可视化分析环境，它针对三维矢量场数据，设计了粒子对流方法用动画表达矢量场数据，最后实验可视化效果明显。本发明的方法更加逼真地模拟了全球气候变化，突破了传统的二维显示的局限性，增加了模型和数据的感知性和真实性，获得更高层次的信息细节，使用者不仅能够真实的感觉到研究的对象近在眼前，而且还能通过手势操作研究对象，从而构建了一个直观的空间数据可视化分析环境。
附图说明
[0068] 图2是本发明可视化案例示意图；
[0069] 图3是本发明手势定义图。
具体实施方式
[0070] 下面结合
说明书附图对本发明的基于VR及手势交互技术的全球气候矢量场数据可视化方法进行详细说明。
[0071] 如图1所示，本发明的基于VR及手势交互技术的全球气候矢量场数据可视化方法，其矢量场数据包括矢量风、洋流等，由于洋流的可视化方法与矢量风近似，下面仅以对矢量风可视化实施过程进行说明。它包括在processing中导入Oculus Rift库，然后搭建一个3D环境，设计包含一个地球和
宇宙背景的场景，将从Natural Earth获取的地图数据添加到所设计的地球上，然后从美国国家气象局操作的全球预报系统(GFS)上下载所需的矢量风数据，对获得的数据文件进行解码，随后对数据进行风向快速参考，用于转换向量风分量、速度和方向，进行地转风近似求解，网格点采用双线性插值，将矢量风投影到所设计的地球上，有限差分近似用于在插值过程期间估计该失真，确保风粒子路径被正确渲染，最后在processing中导入leapmotion的库，根据需要定义各种手势实现交互八个步骤。具体如下：
[0072] 步骤1)导入用于processing的Oculus Rift库：OculusRift.pde，该库目前是以PDE格式发布的。搭建3D环境，设计了一个场景并根据球面体绘制方法构造一个球形代理几何，所述球面体绘制方法如下：球面体绘制方法基于传统的光线投射算法，采用球面坐标表达构造一个球形代理几何，并以球面坐标系的三个空间维度坐标作为纹理坐标，以归一化后的三维体数据作为纹理坐标，以所述球形代理几何作为载体构造三维体数据的三维体纹理，而在绘制过程前，又将球面坐标转换回直角坐标表达，通过求解二次方程计算投射光线与代理几何的交点；
[0073] 步骤2)获取地图数据，并对地图数据进行预处理，再绘制到步骤1)设计场景中的球体上，对步骤1)设计场景中的球体，地图数据由Natural Earth提供，但必须转换为TopoJSON格式，同时需安装GDAL和TopoJSON；所述的地图数据预处理方法如下：将所获取的地图数据为GeoJSON格式转换为D3.js适用的TopoJSON格式，即首先在processing中安装GDAL和TopoJSON，然后将GeoJSON数据文件中边界线记录次数改为一次，且将浮点数据转换为整数形式；
[0074] 步骤3)获取可视化的矢量风数据；
[0075] 天气数据由美国国家气象局操作的全球预报系统(GFS)生成。预测每天生成四次，可从NOMADS下载。这些文件是GRIB2格式，包含超过300条记录。选择这些记录中的几个来可视化特定等压线上的
风力数据，对天气数据预处理，将步骤(3)中气象数据GRIB2文件的数据解码为JSON，需使用netCDF-Java GRIB
解码器；
[0076] 步骤4)对步骤3)中的数据进行风向快速参考，用于转换向量风分量，速度和方向，所述风向快速参考方法如下：由于所得到的气象风数据以角度表示，为方便计算，将角度转化成度，通过乘以DperR(180/π＝57.29578)将角度转换为度，并通过乘以RperD(π/180＝0.01745329)将角度转换为弧度：
[0077] 地理风坐标系：Ugeo，Vgeo；
[0078] +Ugeo分量代表风向东，即西风；-Ugeo表示风向西，即东风；
[0079] +Vgeo是风向北，即南风；
[0080] 地理风向：Dirgeo是相对于真正北方的方向，其中0＝北，90＝东，180＝南，270＝西，Dirgeo＝atan2(-Ugeo，-Vgeo)*DperR＝270-(atan2(Vgeo，Ugeo)*DperR)；
[0081] 水平风速Spd：Spd＝sqrt(Ugeo*Ugeo+Vgeo*Ugeo)；要从Dirgeo和Spd转换为Ugeo和Vgeo：Ugeo＝-Spd*sin(Dirgeo*RperD)；Vgeo＝-Spd*cos(Dirgeo*RperD)；
[0082] 步骤5)对经步骤4)处理后获得转换后的向量风分量再进行地转风近似求解，所述地转风近似求解如下：步骤4)的风向快速参考转换了向量风分量，近似被分解成其两个水平分量U和V，U分量表示风的东西分量，V分量表示南北分量，风的U分量和V分量的两个公式如下：
[0083] U＝((Y1-Y0)×R×Tavg×ln(P0/P1))/(2×omega×sin(lat)×(D^(2)))，
[0084] V＝((X1-X0)×R×Tavg×ln(P1/P0))/(2*omega×sin(lat)×(D^(2)))
[0085] 其中：(Y1-Y0)为当前位置和另一个参考点之间的距离的南北分量的差值，以米为单位；
[0086] (X1-X0)为当前位置和另一个参考点之间的东西向分量的差值，以米为单位；
[0087] R为干燥空气的平均气体常数，为287焦耳/千克*开氏度；
[0088] Tavg为当前位置和另一个参考点之间的平均温度，单位为开氏温标；
[0089] P0为当前所在位置的大气压力；
[0090] P1为另一个参考点处的大气压力；
[0091] omega为地球旋转的角速度＝7.292×10^(-5)倒数秒；
[0092] Sin(lat)为当前所在位置的纬度的sin函数；
[0093] D为当前位置与另一个参考点之间的距离，以米为单位；
[0094] 步骤6)将步骤5)处理后的数据投影到步骤2)的球体上，该的矢量风数据投影方法如下：采用双线性内插先在x轴上执行线性插值，然后再在y轴上执行线性插值，尽管每个步骤在采样值和位置中是线性的，但是内插作为整体不是线性的，而是在采样位置的二次方，为了获取函数f在点(x，y)上的值，已知四个点Q11＝(x1，y1)，Q12＝(x1，y2)，Q21＝(x2，y1)和Q22＝(x2，y2)的f的值，我们首先在x方向进行线性插值，产生：
[0095]
[0096]
[0097] 再通过在y方向上内插来获得期望的估计：
[0098]
[0099] 步骤7)有限差分近似用于在步骤(7)插值过程期间估计该失真，确保风粒子路径被正确渲染，得到可视化结果；
[0100] 需要投影坐标(x，y)相对于球面坐标(lat，lon)的所有四个一阶导数＝(phi，lambda)：dx/d(phi),dx/d(lambda)；dy/d(phi),dy/d(lambda)。关于Tissot Indicatrix(天梭指示线)的其它信息是根据这些(使用一些算术和三角函数：余弦，主反正弦和主反正切)计算的。计算需要描述地球的形状，为了获得最大精度，使用具有半长轴a和偏心率e的椭圆体基准；
[0101] 步骤7)的有限差分近似风粒子路径方法如下：
[0102] 步骤7-1)使用一阶中心有限差分估计在距离h＝10^(-5.2)弧度；
[0103] 步骤7-2)在点(phi，lambda)获得dx/d(phi)，让投影点
[0104] (phi-h/2,lambda)-->(x0,y0)
[0105] (phi+h/2,lambda)-->(x1,y1)
[0106] 使用估计值
[0107] dx/d(phi)＝(x1-x0)/h
[0108] dy/d(phi)＝(y1-y0)/h
[0109] 同样，投影点
[0110] (phi,lambda-h/2)-->(x2,y2)
[0111] (phi,lambda+h/2)-->(x3,y3)
[0112] 并使用估计
[0113] dx/d(lambda)＝(x3-x2)/h
[0114] dy/d(lambda)＝(y3-y2)/h
[0115] 步骤7-3)从上述衍生物连同斯奈德公式，可以获得在(phi，lambda)的Tissot Indicatrix的轴的长度及其方向；
[0116] 步骤7-4)根据地图的大小来确定比例因子，找到地图上的典型TI的大小，并缩放，使得那些TI将大约为地图的6％；
[0117] 步骤7-5)以相同的量重新缩放所有的TI，因此它们可以被比较，并且将每个重新缩放在其自己的中心附近(其通过第五投影(phi，lambda)→(x，y))；
[0118] 步骤8，交互设计，在processing中导入leapmotion的库，根据需求采用手势交互技术定义各种操作手势，实现全球气候矢量场数据可视化，操作手势包括水平竖直移动、旋转、放大及缩小；
[0119] 1)经过立体标定后获取经过校准的立体图像对后，进行立体匹配，获得视差图像；
[0120] 2)再利用摄像机的内参数及外参数进行三角计算获取深度图像，然后对左或右视觉图像使用手势分割算法处理，分割出人手所在的初始位置信息；
[0121] 3)预先设定阈值，当人手左右水平移动超过阈值，即对地球进行水平移动；当人手竖直移动超过阈值，即对地球进行竖直移动；人手五指握拳，实现缩小操作；人手由拳展开，实现放大操作；手掌旋转超过阈值角度，即对地球进行旋转。
[0122] leap motion手势设置如下：手掌左右水平移动控制地球水平移动；手掌上下竖直移动控制地球竖直移动；手掌五指握拳控制地球收缩；手掌五指展开控制地球放大；手掌转动控制地球旋转。
