news 2026/7/24 22:30:36

WRF模型性能评估

作者头像

张小明

前端开发工程师

1.2k 24
文章封面图
WRF模型性能评估

目录

一、数据下载（step1_data_download）

1. 筛选出需要的站点

（1）ArcGIS Pro 中打开 isd-history.csv

（2）按经纬度范围筛选（比如按 WRF/CMAQ 域）

2. 导出选中的站点

3. 下载ISD数据并合并

（1）把usaf和wban写入一个txt中

（2）下载数据

（3）文件整合

二、转换格式（step2_data_format）

三、转成metstat的输入格式（step3_data_to_metstat）

1. 准备两个输入文件：

2. 生成观测数据

四、运行metstat

五、结果整理

1. 步骤4的输出，放excel里看，用逗号分列

2. 用ncl整理一下，如daily结果的最后一列，是月均的结果，就可以直接拎出来

一、数据下载（step1_data_download）

1. 筛选出需要的站点

所用的数据来自NCDC：The Integrated Surface Database (ISD) is a global database that consists of hourly and synoptic surface observations compiled from numerous sources into a single common ASCII format and common data model.

数据网站：

https://www1.ncdc.noaa.gov/pub/data/noaa/、https://www.ncei.noaa.gov/pub/data/noaa/（一样）

wget ftp://ftp.ncdc.noaa.gov/pub/data/noaa/isd-history.txt
wget ftp://ftp.ncdc.noaa.gov/pub/data/noaa/isd-history.csv

下载isd-history.csv，记录了所有的站点

然后直接根据经纬度选择站点即可

也可以用arcmap看站点的分布（非必要）：打开ArcMAP，add data, 并根据经纬度显示，如下图

根据模拟域筛选出需要的站点，获取ISD station USAF IDs, 包括csv文件中的usaf和wban

（1）ArcGIS Pro 中打开 isd-history.csv

a. 把 CSV 加到工程里

  • 打开 ArcGIS Pro，进入你的项目（或新建一个）。
  • 打开 窗口- Catalog（目录） 窗格里：
    • 链接到文件夹，右键 → Add Folder Connection…，连接到 isd-history.csv 的文件夹。
    • 展开这个文件夹，你就能看到 isd-history.csv。

b. 把经纬度变成点（Add XY Data）

  • 确认 isd-history.csv 里面有经纬度列，一般是 LAT 和 LON（有时是 LATITUDE / LONGITUDE，看一下刚才打开的表头）。
  • 右击isd-history.csv，选择创建要素类 – 从XY表
  • 参数填写：
    • Input Table：选择 isd-history.csv。
    • X Field：选经度列，比如 LON。
    • Y Field：选纬度列，比如 LAT。
    • Z Field：可以空着。
    • Coordinate System：点右侧的地球图标 → 选 地理坐标系（Geographic Coordinate Systems）→ World → WGS 1984（GCS_WGS_1984）。

生成的文件：

（2）按经纬度范围筛选（比如按 WRF/CMAQ 域）

a. 窗口 - 搜索 – 按属性选择图层：Select By Attributes

  • 在弹出的窗口里：
    • 图层名称：XYisd-history（站点图层）
    • 选择类型：New selection

b. 写筛选语句（Expression）

  • "LAT" >= 10 AND "LAT" <= 60 AND "LON" >= 70 AND "LON" <= 150
  • c. 点 Apply / OK
  • 地图上就只高亮出这个框里的站点，如图

2. 导出选中的站点

1. 导出站点信息属性表

（1）在图层上右击打开属性表

（2）选择表选项 – 导出 所选站点的属性表

2. 导出站点图片

右击XYisd-histor图层，选择数据——导出数据——输出Export_Output.shp（如图）

3. 下载ISD数据并合并

（1）把usaf和wban写入一个txt中

比如命名为stn_ID.txt，需注意的是usaf是6位，wban是5位，如果不够需要补0在前面

（2）下载数据
cat downloadISD.csh #!/bin/csh -f set year=2020 mkdir -p ${year} set instr = "" foreach instr(`cat stn_ID.txt | tr -d ‘\r’`) # \r删除文件stn_ID.txt的回车字符 echo $instr set data = ftp://ftp.ncdc.noaa.gov/pub/data/noaa/${year}/${instr}-${year}.gz wget $data gzip -d ${instr}-${year}.gz mv ${instr}-${year} ./${year}/ end
（3）文件整合

如cat ./2020/* > isd2020

二、转换格式（step2_data_format）

从网址（https://www1.ncdc.noaa.gov/pub/data/noaa/）下载ishJava.java和ishJava.class

其实只用下载ishJava.java，class可以用该脚本生成。注意ishJava.java中有个地方要改：

注释掉539-542行，不要把WBAN的99999改成*****

然后生成class, 命令为：javac ishJava.java

建立超链接：ln -s路径/isd2020 isd2020 （ln -s ../step1_data_download/isd2013 ./）

运行：java -classpath . ishJava isd2020 isd2020.out

生成的文件格式在ish-abbreviated.txt里有介绍

Note

如果没有java, 下载jre, https://www.java.com/en/download/manual.jsp

如果没有javac, 下载jdk, https://www.oracle.com/java/technologies/downloads/

isd2020.out的格式为：

三、转成metstat的输入格式（step3_data_to_metstat）

1. 准备两个输入文件：

FILE_DATA_2013：记录的是第2步骤输出文件的完整路径

FILE_LAT_LON：记录的站点信息的文件的完整路径

FILE_LAT_LON.txt生成代码：

cat 1_writestnLatLonFile.ncl begin tmp = asciiread("/work/home/gzq22/apprepo/evaluate/WRF_evaluation/step1_data_download/isd-history.txt",(/29683/),"string") strs = tmp(22:29682) number = dimsizes(strs) tmp2 = asciiread("/work/home/gzq22/apprepo/evaluate/WRF_evaluation/step1_data_download/stn_ID_China_new.txt",-1,"string") USAFneed = str_get_cols(tmp2,0,5) WBANneed = str_get_cols(tmp2,7,11) numneed = dimsizes(USAFneed) print(""+numneed) USAF = str_get_cols(strs,0,5) WBAN = str_get_cols(strs,7,11) STATIONNAME = str_get_cols(strs,13,41) CTRY = str_get_cols(strs,43,44) State = str_get_cols(strs,48,49) LAT = str_get_cols(strs,57,63) LON = str_get_cols(strs,65,72) ELEV = str_get_cols(strs,74,80) ;722355 93999 VICKSBURG MUNICIPAL AIRPORT US MS KVKS +32.233 -090.933 +0031.4 20130319 20240422 ;------------ ------------------------------ -------------------------------------------------- ------------------------------ -------- --------- --------- header = (/"USAF-WBAN_ID STATION NAME COUNTRY STATE LATITUDE LONGITUDE ELEVATION", \ "------------ ------------------------------ -------------------------------------------------- ------------------------------ -------- --------- ---------"/) hlist = [/header/] fname = "FILE_LAT_LON.txt" system("rm -f "+fname) write_table(fname, "w", hlist, "%s") ;rm no lat lon LAT2 = str_sub_str(LAT," ","0") f_lat = stringtofloat(LAT2) do j=0,numneed-1 print(tostring(j)+": "+USAFneed(j)+" "+WBANneed(j)) index = ind(USAF.eq.USAFneed(j) .and. WBAN.eq.WBANneed(j)) ;print(index) if (.not. ismissing(index)) then i = index(0) write_table(fname, "a", [/USAF(i),WBAN(i),STATIONNAME(i),CTRY(i),State(i),LAT(i),LON(i),ELEV(i)/],"%s %s %s %s %s %s %s %s ") end if end do end

2. 生成观测数据

cat 2_batch_script #! /bin/csh -f #set year = (2013 2014 2015 2016 2017 2018 2019 2020 2021 2022 2023) set year = (2022 2023) set month = (01 02 03 04 05 06 07 08 09 10 11 12) set m_DIR = (JAN FEB MAR APR MAY JUN JUL AUG SEP OCT NOV DEC) set m_num = (1 2 3 4 5 6 7 8 9 10 11 12) foreach y (1 2 3 4 5 6 7 8 9 10 11) foreach m (1 2 3 4 5 6 7 8 9 10 11 12) echo "==============================" echo "Now running ${year[$y]}-${month[$m]}" echo "==============================" sed -e 's/MMMMMM/'${m_num[$m]}'/' \ -e 's/MONMON/'${month[$m]}'/g' \ -e 's/YYYYYY/'${year[$y]}'/g' dpp.s > dpp.s_${year[$y]}_${month[$m]} chmod +x dpp.s_${year[$y]}_${month[$m]} ./dpp.s_${year[$y]}_${month[$m]} FILE_DATA_${year[$y]} FILE_LAT_LON rm -f dpp.s_${year[$y]}_${month[$m]} end end

四、运行metstat （step4_run_metstat）

cat metstat.wrf #! /bin/csh -f set xw = (-1000. -1000.) set xe = ( 1000. 1000.) set ys = (-1000. -1000.) set yn = ( 1000. 1000.) set sub = (cn4_NCP cn9) set region = (CH) set domain = (d01 d02 d03) set year = (2024 2025) set month = (01 02 03 04 05 06 07 08 09 10 11 12) set s_day = (01 01 01 01 01 01 01 01 01 01 01 01) set e_day = (31 28 31 30 31 30 31 31 30 31 30 31) set e_day = (31 29 31 30 31 30 31 31 30 31 30 31) set wrf = /public/home/acez11b6ht/apprepo/wrf_wps/4.5-intelmpi2017/app/WRF/test foreach d (1) ## domain sub foreach y (2) ## year foreach m (7) ## month s_day e_day foreach i (1) ## region ./src_v3/metstat << IEOF Run note |WRF validation - $sub[$d] Hourly output file |./output/hrly_${year[$y]}_${month[$m]}_${domain[$d]}_$sub[$d]_$region[$i] Daily output file |./output/daly_${year[$y]}_${month[$m]}_${domain[$d]}_$sub[$d]_$region[$i] Input sfc obs |../step3_data_to_metstat/output/NCDC_${year[$y]}_${month[$m]}_${region[$i]} Tstart(yyyymmddhh) |${year[$y]}$month[$m]$s_day[$m]00 Tend |${year[$y]}$month[$m]$e_day[$m]23 Time Zone |0 Map Projection |LCP # sites to proc |-1 |$xw[$d] $xe[$d] |$ys[$d] $yn[$d] Met model output |WRF No met output files|$e_day[$m] Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-01_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-02_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-03_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-04_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-05_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-06_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-07_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-08_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-09_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-10_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-11_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-12_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-13_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-14_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-15_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-16_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-17_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-18_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-19_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-20_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-21_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-22_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-23_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-24_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-25_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-26_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-27_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-28_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-29_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-30_00:00:00 Met file name |$wrf/em_real_$year[$y]$month[$m]/WRFOUT/wrfout_${domain[$d]}_$year[$y]-$month[$m]-31_00:00:00 IEOF end end ## foreach regions end ## foreach m (1) month end ##foreach y (1)

五、结果整理 （step5_data_aggregate）

1. 步骤4的输出，放excel里看，用逗号分

2. 用ncl整理一下，如daily结果的最后一列，是月均的结果，就可以直接拎出来

cat aggregate.ncl function trim_left (s:string) local c, nvec, n, i0, i, c2 begin c = stringtochar(s) nvec = dimsizes(c) n = nvec(0) i0 = 0 do i = 0, n-1 if (c(i) .ne. " ") then i0 = i break end if end do c2 = c(i0:) return ( chartostring(c2) ) end function pad_right (s:string, width:integer) local c, nvec, n, t, pad, i begin c = stringtochar(s) nvec = dimsizes(c) n = nvec(0) if (n .ge. width) then return (s) end if t = s pad = width - n do i = 1, pad t = t + " " end do return (t) end function pad_left (s:string, width:integer) local c, nvec, n, t, pad, i begin c = stringtochar(s) nvec = dimsizes(c) n = nvec(0) if (n .ge. width) then return (s) end if t = "" pad = width - n do i = 1, pad t = t + " " end do t = t + s return (t) end ; ---------------- main program ---------------- begin region = (/"CH"/) year = (/"2018","2020","2021","2022","2023"/) dom = (/"d01"/) res = (/"cn4_NCP"/) ; 12 months month = (/"01","02","03","04","05","06","07","08","09","10","11","12"/) ; column index in daily files (your original nn) nn = (/ 35 , 32 , 35 , 34 , 35 , 34 , 35 , 35 , 34 , 35 , 34 , 35 /) ; only use Apr–Sep (4–9) m_start = 6 m_end = 8 nmon = m_end - m_start + 1 ; = 6 base_dir = "/work/home/gzq22/apprepo/evaluate/WRF_evaluation/step4_run_metstat/output/" out_dir = "/work/home/gzq22/apprepo/evaluate/WRF_evaluation/step5_data_aggregate/output/" ; column widths w_var = 10 w_stat = 18 w_unit = 8 w_month = 10 do yy = 0, dimsizes(year)-1 do ireg = 0, dimsizes(region)-1 ;========== 1. use April file to get row labels ========== mm0 = month(m_start-1) ; "04" file0 = base_dir + "daly_" + year(yy) + "_" + mm0 + "_" + dom + "_" + res + "_" + region(ireg) print(">>> using file for labels: " + file0) strs0 = asciiread(file0, -1, "string") title = strs0(0) ; line 1: WRF validation - ... lines = strs0(2:) ; from line 3: statistics nrow = dimsizes(lines) ; get Var / Stat / Unit (trim left spaces) var_en_raw = str_get_field(lines, 1, ",") stat_en_raw = str_get_field(lines, 2, ",") unit_raw = str_get_field(lines, 3, ",") var_en = new(nrow, "string") stat_en = new(nrow, "string") unit = new(nrow, "string") do j = 0, nrow-1 var_en(j) = trim_left(var_en_raw(j)) stat_en(j) = trim_left(stat_en_raw(j)) unit(j) = trim_left(unit_raw(j)) end do ;========== 2. read selected column from Apr–Sep daily files ========== vals = new((/nmon, nrow/), "string") do km = 0, nmon-1 m = m_start + km mm = month(m-1) filename = base_dir + "daly_" + year(yy) + "_" + mm + "_" + dom + "_" + res + "_" + region(ireg) print(">>> reading: " + filename) strs_m = asciiread(filename, -1, "string") lines_m = strs_m(2:) vals_raw = str_get_field(lines_m, nn(m-1), ",") do j = 0, nrow-1 vals(km,j) = trim_left(vals_raw(j)) end do end do ;========== 3. build aligned header ========== header = pad_right("Var", w_var) + " " + pad_right("Stat", w_stat) + " " + pad_right("Unit", w_unit) do km = 0, nmon-1 mm = month(m_start-1 + km) header = header + " " + pad_left(year(yy) + mm, w_month) end do ;========== 4. build each row ========== data_out = new(nrow, "string") do j = 0, nrow-1 line = pad_right(var_en(j), w_var) + " " + pad_right(stat_en(j), w_stat) + " " + pad_right(unit(j), w_unit) do km = 0, nmon-1 line = line + " " + pad_left(vals(km,j), w_month) end do data_out(j) = line end do ;========== 5. write aligned text file ========== outfname = out_dir + year(yy) + "_" + dom + "_" + res + region(ireg) + "_daly_summary.csv" system("rm -f " + outfname) out_all = new(nrow+2, "string") out_all(0) = title out_all(1) = header out_all(2:) = data_out asciiwrite(outfname, out_all) print(">>> wrote file: " + outfname) end do end do end
版权声明: 本文来自互联网用户投稿，该文观点仅代表作者本人，不代表本站立场。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如若内容造成侵权/违法违规/事实不符，请联系邮箱：809451989@qq.com进行投诉反馈，一经查实，立即删除！
网站建设 2026/7/24 22:29:04

HarmonyOS7 多媒体能力怎么放进一个页面：MediaPanelSwitcher 组合式案例解析

文章目录前言这个案例适合学什么完整代码先看页面骨架&#xff1a;一个状态切三块内容切换头部为什么用 ForEachWeb 面板讲的是“控制感”Canvas 面板讲的是“自定义表达”Image 面板讲的是“资源展示效率”这三个组件到底怎么选真实业务里怎么改造成可用页面这类组合页容易踩的…

作者头像 李华
网站建设 2026/7/24 22:29:02

AI大模型实战指南：低成本高效应用方案

1. 项目概述&#xff1a;AI大模型平民化应用指南去年帮朋友公司优化客服系统时&#xff0c;我用三个不同的大模型搭建了智能应答矩阵&#xff0c;成本比原方案降低67%的同时响应速度提升3倍。这个经历让我意识到&#xff0c;现在正是普通人驾驭顶尖AI的最佳时机——就像2007年智…

作者头像 李华
网站建设 2026/7/24 22:28:47

让你的桌面活起来：呆啵宠物DyberPet打造个性化AI陪伴体验

让你的桌面活起来&#xff1a;呆啵宠物DyberPet打造个性化AI陪伴体验 【免费下载链接】DyberPet Desktop Cyber Pet Framework based on PySide6 项目地址: https://gitcode.com/GitHub_Trending/dy/DyberPet 你是否曾梦想过让你喜爱的动漫角色、游戏伙伴或原创角色真正…

作者头像 李华
网站建设 2026/7/24 22:20:28

Claude Code 配置完全指南（一）：.claude 目录结构全景解析

Claude Code 配置完全指南&#xff08;一&#xff09;&#xff1a;.claude 目录结构全景解析系列第 1 篇 | 2026-07-22 配套仓库&#xff1a;C:\Users\zhang\.claude前言 Claude Code 的配置全部藏在 ~/.claude 目录下。这个目录不像 VS Code 的 settings.json 那样一目了然——…

作者头像 李华