基于R语言的风廓线雷达产品数据简单处理文档格式.docx
《基于R语言的风廓线雷达产品数据简单处理文档格式.docx》由会员分享,可在线阅读,更多相关《基于R语言的风廓线雷达产品数据简单处理文档格式.docx(5页珍藏版)》请在冰豆网上搜索。
虽然有着处理速度快,计算结果可靠等优点。
但是其使用并不便:
需要掌握编程语法和算法,不同计算机平台可能需要不同的编译系统等。
所以本文针对风廓线雷达资料,使用R语言来对其进行简单处理。
R语言的全称:
为统计计算和图形展示而设计的一种编程语言和统计环境。
目前有一个R核心开发团队对其进行定期维护和更新[5]。
其主要的特点在于:
(1)完全免费:
可以在其官方网站(http:
//cran.r-project.org)下载到完整的安装包并免费使用;
(2)开源软件:
R语言源代码完全公开,任何人均能提供各种数据处理模块,下载相应模块后能迅速地完成数据处理等工作;
(3)图形显示:
可以利用R语言处理数据后,直接获取各种统计分析图形;
(4)多平台使用:
R语言可以在Windows、mac、Unix这些操作系统上安装,不需要重复编译,其脚本在不同操作系统间可以任意使用。
针对风廓线雷达资料,我们可以利用R语言对其进行简要的数据处理,统计分析,然后利用图形显示系统得到各种分析图片。
本文将在第二部分介绍R语言读取和储存风廓线雷达产品数据资料;
文章第三部分对资料进行简要处理,最后会利用R语言的图形显示功能作图。
1雷达产品数据的读取
测站观测的风廓线雷达资料有两种数据[6]:
原始数据和产品数据。
其中原始数据为二进制格式,主要包含功率谱数据文件,瞬时径向谱数据文件。
产品数据文件为文本格式,包括实时的采样高度上的、半小时平均的采样高度上的、一小时平均的采样高度上的产品数据文件。
本文所处理资料为实时的采样高度上的产品数据文件。
R语言也能读取处理二进制格式文件,对雷达原始数据资料的处理会在以后的工作中进行。
1.1单个产品数据的读取
R语言有许多函数能够直接从文本文件中读取数据,比较常用的有:
read.table(),read.csv(),read.fwf()。
其中和Fortran比较接近的是read.fwf(),可以指定读取数据的长度和格式。
于Fortran不同的是,R语言在读取数据的时候,不用先给定数据类型,程序会直接读取数据,并存储到一个数据框(data.frame)里。
例如针对本文要处理的雷达产品数据,可以直接使用read.table()命令读取产品数据:
raw_data&
lt;
-read.table(fname,fill=TRUE)
其中raw_data为一个数据框(data.frame),用于储存我们需要处理的雷达产品数据,fname为要读取的雷达数据文件名字;
“header=FALSE”表示该实体数据中没有数据说明头文件;
由于本文要处理的雷达产品数据不是规则的表格形式,所以需要使用参数:
“fill=TRUE”,来自动填满不是数据表格的部分。
产品数据文件的前三行为测站基本参数,最后一行为结束行,中间部分为实际数据,包含:
采样高度,水平风向,水平风速,垂直风速,水平方向可信度,垂直方向可信度,Cn2。
为了便于资料处理,将实体数据单独储存到一个名为r_data的新数据框中:
r_data&
-raw_data[4:
(length(raw_data[,1])-1),]
其中使用“length(raw_data[,1])”函数判断raw_datpr_data&
-lapply(r_data,as.numeric)
-as.data.frame(r_data)
先使用函数“lapply”和“as.character”把每列元素转换为字符型,再使用as.numeric转换为数值型,最后再用as.data.frame把r_data转换成一个数据框。
从这些处理过程中不难看出,R语言能十分方便的利用“lapply”函数实现整列(整行)的数据处理,相应的apply类函数还有许多,从而省去了程序中的循环语句编写;
而且R的数据类型使用十分灵活,能够方便的将其转换为不同的数据类型,类似的函数还有as.matrix,as.logical等。
R在使用as.numeric时,会自动将无法转换成数值的字符,转换成R的缺测值(NA),本数据中的“//////”在经过转换后全部转换为了“NA”。
“NA”在R语言中可以参与计算,也可以使用一个简单的函数将数据中的“NA”去除。
例如针对本文中的数据:
m_data&
-r_data[complete.cases(r_data),]
其中complete.case函数能够获取数据中不含缺测的所有列,进而赋值后的数据框m_data中剔除了r_data的缺测值。
1.2批量数据文件的读取
通过1.1中几个简单的语句即可完成单个雷达产品数据读取和初始化,而在业务运行中,雷达产品数据肯定是大批量的生成,也需要程序脚本具有批处理功能。
我们将1.1中单个文件的雷达产品数据处理过程,整合成一个R语言函数readradar:
readradar&
-function(dir,fname){……}
readradar只需要给定路径和文件名,即可读取雷达产品数据资料并去除缺测值。
函数最后返还一个数据框,其中包含该文件中所有非缺测产品数据。
给定所有需要处理的路径和文件名后,即可完成资料的批处理。
而R语言能够十分便利的获取文件名,因为其能通过脚本语言进入当前运行的操作系统。
简单的说,R语言能够在程序内部完成操作系统的文件、文件夹处理、安装包的安装等。
例如可以直接使用file.create(create.dir)函数直接在R语言中生成新的文件(文件夹)等。
本文只需要R语言读取要处理的文件名:
fname&
-list.files(path=“./Qingdao/”
fname中存储了所有“./Qingdao”下的雷达数据文件名。
通过循环即可完成雷达数据的批量读取:
for(iin1:
length(fname)){
-readradar(dir,fname[i])
}
别的批量处理过程(资料简要处理,雷达风廓线图等)和读取类似,只需要加入循环即可。
2雷达资料的简要处理和画图
第二节中给出了R语言对雷达产品数据文件的批读取。
R语言最为实用的优点在于其计算和画图功能。
2.1雷达风廓线的简要计算
为了便于处理资料。
我们利用names函数将m_data数据框的每一列分别命名:
names(m_data)&
-c(“hgt”,”h_dir”,”h_speed”,”w_speed”,”h_rel”,”w_rel”,”cn2”)
每一列表示对应的资料的采样高度,水平风向,水平风速,垂直风速,水平方向可信度,垂直方向可信度,Cn2。
水平风的u分量计算如下:
m_data$u&
-m_data$h_speed*cos((270-m_data$h_dir)*pi/180)
其中pi为R语言中自带的圆周率的取值3.141593。
而计算出来的u风速可以直接存储到数据框m_data的新的一列中。
同理可以计算出水平风速的v分量:
m_data$v&
-m_data$h_speed*sin((270-m_data$h_dir)*pi/180)
上述计算过程可以看出,R语言中数据框的计算处理十分方便,不需要使用循环语句计算不同高度上的u/v分量,而且计算结果可以直接作为数据框新的一列存到原数据框中。
而R语言除了基本的计算外,还能使用内部函数做各种复杂数学计算,例如矩阵的求逆、线性回归分析、抽样分布、显著性试验等。
如果R语言自带的函数已经不能解决当前数学问题,还能上网搜索下载对应函数包。
2.2雷达风廓线图
R语言主要有三种绘图函数:
高级、低级和交互式。
通过调用高级绘图函数,能在R语言中直接绘制各种统计图;
低级绘图函数能够对现有的图进行修改;
交互式绘图能够让用户直接利用鼠标修改图形。
大量的内置函数让绘图变得十分的简易。
例如本文需要的u风场随高度变化图可以通过以下函数实现:
plot(m_data$u,m_data$hgt)
在plot函数中增加各种参数,能够优化所绘图形。
例如本文中使用如下命令获取风廓线图:
plot(m_data$u,m_data$hgt,type=“b”,
main=title,xlim=c(-20,20),ylim=c(0,5000),
xlab=“风速(m/s)”,ylab=“采样高度(m)”,pch=16,col=2)
其中tpye=“b”表示曲线为点画线;
main为主标题;
xlim和ylim表示横纵坐标取值范围;
xlab、ylab表示横、纵坐标标题;
pch=16表示点为实心圆圈,col=12表示颜色为红色。
通过使用低级绘图命令如points、lines等,可以在刚画的u风场廓线后增加v风场廓线:
points(m_data$v,m_data$hgt)
lines(m_data$v,m_data$hgt)
图例则可以用lengend函数添加:
legend(“topleft”,pch=c(16,17),col=c(2,19),lty=c(1,1),legend=c(“u”,”v”),bty=“n”)
R语言还能使用par函数对图片进行设置,例如本文中使用如下函数将2014年4月26日每隔3小时的风廓线显示到同一图片中:
par(mfrow=c(2,4),mar=c(5,4,4,2)+0.1,oma=c(0.1,0.1,0.4,0.1))
其中mfrow=c(2,4)表示图片分割为两行四列,参数mar和oma指定了图片的间距。
将上述命令和第二节的数据读取函数整合到一个R脚本中,运行脚本即可得到如图1风廓线图。
R语言的绘图函数很多,除了文中绘制风廓线图外还能绘制:
直方图、散点图、饼图等。
在安装绘图包后还能绘制3D图等。
其绘图功能,能满足大部分气象资料统计分析的出图需求。
3结语和讨论
免费的开源软件R语言,能够用于数据分析和图形显示。
文中使用其对风廓线雷达资料进行了批量读取、计算和绘图。
使用read.table函数简洁的实现资料读取。
R语言的操作系统函数,能便利的实现资料的批处理。
R语言的各种绘图函数,能迅速地绘制较为美观的风廓线图。
R语言的统计计算功能十分强大,本文只是简单的使用基础计算,在随后的工作中,可以使用R语言实现,雷达数据的质量控制,统计分析等。
本文只分析了雷达产品数据,R语言也能处理二进制数据资料,以后的工作中可以利用R语言分析其他类型的气象资料。
参考文献
[1]Green,J.L.,AtmosphericmeasurementsbyVHFpulsedDopplerradar.IEEETrans.OnGeoscienceelectronics.GE-17(4):
262-280.1979
[2]何平.相控阵风廓线雷达[M].气象出版社,2006:
104—122.
[3]王欣,卞林根,彭浩,李剑东.风廓线仪系统探测试验与应用[J].应用气象学报,2005,16(5):
693-698.
[4]胡明宝.风廓线雷达数据处理和应用研究[D].南京信息工程大学博士论文,2012.
[5]JohnM.Chambers.FacetsofR.TheRJournal,1
(1):
5-8,2009.
[6]关于进行风廓线雷达数据传输的通知.气测函[2011]61号.内部资料