MODE目标检验改进版
<p>[TOC]</p>
<pre><code class="language-python">import meteva.base as meb
import meteva.product as mpd
import numpy as np
import pandas as pd
import meteva.method as mem
import math
import copy
import datetime
import xarray as xr
import matplotlib.pyplot as plt
import datetime</code></pre>
<pre><code class="language-python">#读取观测和预报数据
grid1 = meb.grid([100,120,0.05],[24,40,0.05])
filename_ob = r'H:\test_data\input\mem\mode\ob\rain03\20072611.000.nc'
filename_fo_03 = r'H:\test_data\input\mem\mode\ec\rain03\20072608.003.nc'
filename_fo_27 = r'H:\test_data\input\mem\mode\ec\rain03\20072508.027.nc'
grd_ob = meb.read_griddata_from_nc(filename_ob,grid = grid1,time = &quot;2020072611&quot;,dtime = 0,data_name = &quot;OBS&quot;)
grd_fo_03 = meb.read_griddata_from_nc(filename_fo_03,grid=grid1,time=&quot;2020072608&quot;,dtime=3,data_name=&quot;ECMWF&quot;)
grd_fo_27 = meb.read_griddata_from_nc(filename_fo_27,grid=grid1,time=&quot;2020072508&quot;,dtime=27,data_name=&quot;ECMWF&quot;)</code></pre>
<h1>MODE方法的问题</h1>
<p>在将MODE方法应用于不同时效或不同模式的检验对比时存在一个问题：相同的观测场和不同的预报场匹配时，得到观测目标是不相同的，这意味着基于MODE方法检验得到的不同时效或不同模式的预报误差没有严格的可对比性。通过如下实例可以看出观测和不同预报匹配的结果的不同。</p>
<pre><code class="language-python">print(&quot;观测和003H时效预报的匹配效果&quot;)
look_03 = mem.mode.feature_finder(grd_ob,grd_fo_03,smooth = 5,threshold = 5,minsize = 5)
look_match_03 = mem.mode.centmatch(look_03)
look_merge_03 = mem.mode.merge_force(look_match_03)
mem.mode.plot_label(look_merge_03)
print(&quot;————————————————————————————————&quot;)
print(&quot;观测和027H时效预报的匹配效果&quot;)
look_27 = mem.mode.feature_finder(grd_ob,grd_fo_27,smooth = 5,threshold = 5,minsize = 5)
look_match_27 = mem.mode.centmatch(look_27)
look_merge_27 = mem.mode.merge_force(look_match_27)
mem.mode.plot_label(look_merge_27)</code></pre>
<pre><code>观测和003H时效预报的匹配效果</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=4559f9f981615cbaf62c03220a537dfb" alt="" /></p>
<pre><code>————————————————————————————————
观测和027H时效预报的匹配效果</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=ad2b2481d6c928b4348a06c9fee883d4" alt="" /></p>
<h1>OBV概述</h1>
<hr />
<p><font face="黑体" color=Black size=4>
<br> OBV方法是在MODE方法的基础上改进后的检验方法，改进的点主要包括2个<br>
<font face="黑体" color=Black size=3>
<br> 1. MODE方法识别的1个目标是一片联通的区域，在OBV方法中，会进一步根据临近度将多个联通区域进行合并<br>
<br> 2. MODE方法，在匹配合并的时候，观测目标会进一步合并，在OBV方法中，目标匹配合并过程中，观测目标始终保持不变<br>
<font face="黑体" color=Black size=4><br>OBV 虽然对MODE方法中的识别和匹配有所改进，但是整个方法的体系仍和MODE是相同的，因此在MetEva中OBV并不构成一个独立模块，而是在MODE方法中增加了下面几项函数功能<br></p>
<h1>目标的识别合并</h1>
<p><font face="黑体" color=Blue size=3><strong>feature_finder_and_merge(grd_ob0,smooth, threshold, minsize,compare = ">=",near_dis = 200,near_rate = 0.3)</strong></font> </p>
<p>根据输入的观测和预报网格数据， 识别其中的联通目标，并将这些联通目标依次标记上序号。识别步骤包括 </p>
<ol>
<li>选用一个半径为smooth的圆盘形卷积内核，将观测观测场和预报场进行卷积平滑； </li>
<li>设定一个阈值threshold，将平滑后观测和预报格点场中的数值小于threshold的格点置为0； </li>
<li>通过联通域提取算法识别出观测和预报场中的目标； </li>
<li>设置minsize，将预报或观测场中目标的面积（格点数）小于minsize的删除； </li>
<li>判断两两目标之间的临近度，取多种临近度的最大值，若它大于near_rate,则将其对应的两个目标合并</li>
<li>重复步骤5，直至所有两两目标之间临近度小于near_rate算法终止</li>
</ol>
<p>具体的合并算法请参考文献：<br />
刘凑华, 牛若芸. 基于目标的降水检验方法及应用[J]. 气象, 2013, 39(6):681-690</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: center;">说明</th>
<th style="text-align: left;">备注</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>grd_ob</font></strong></td>
<td style="text-align: center;">网格观测或者预报数据</td>
<td style="text-align: left;"><a href="https://www.showdoc.com.cn/meteva?page_id=3975600815874861">网格数据格式</a></td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>smooth</font></strong></td>
<td style="text-align: center;">平滑系数,采用圆盘形的卷积内核对网格数据进行平滑，smooth相当圆盘的半径</td>
<td style="text-align: left;">非负整数，平滑系数为0时表示不做平滑</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>threshold</font></strong></td>
<td style="text-align: center;">雨量阈值，平滑后低于阈值的部分会被置0</td>
<td style="text-align: left;">实数</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>minsize</font></strong></td>
<td style="text-align: center;">目标的最小面积,面积小于该阈值的目标将会被删除</td>
<td style="text-align: left;">非负整数，面积（格点数）小于minsize的目标将被删除</td>
</tr>
<tr>
<td style="text-align: left;"><strong>compare</strong></td>
<td style="text-align: center;">比较方法，可选项包括">=",">","<=","<"</td>
<td style="text-align: left;">当关注目标是要素场大于阈值时,采用默认参数，如识别暴雨目标时选择参数">=",若关注对象是要素值小于阈值时，则可选取"<="，如根据能见度场识别大雾目标</td>
</tr>
<tr>
<td style="text-align: left;"><strong>near_dis</strong></td>
<td style="text-align: center;">邻域尺度，单位是km</td>
<td style="text-align: left;">该参数越大，目标被合并可能性越高，设为0时，表示所有的目标都不合并；邻域尺度参数应该根据研究者关注的预报对象的尺度来确定通常对大尺度目标的检验，near_dis 设置在100-300km较合适，对于中尺度目标的检验，near_dis 设置为30-100km</td>
</tr>
<tr>
<td style="text-align: left;"><strong>near_rate</strong></td>
<td style="text-align: center;">邻域比例， 表示目标落在相邻目标邻域内的比例阈值</td>
<td style="text-align: left;">取值范围是（0,1], 该阈值越低，相邻目标合并的可能性越大</td>
</tr>
</tbody>
</table>
<p><strong>————————————————————————————————————————————————————————————————————————————————————————————————————————————</strong></p>
<p><font face="黑体" color=green size=4><strong>feature_finder_and_merge函数返回结果内容说明</strong></font> </p>
<table>
<thead>
<tr>
<th style="text-align: left;">一级关键词</th>
<th style="text-align: left;">二级关键词</th>
<th style="text-align: center;">说明</th>
<th style="text-align: left;">备注</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>grd</font></strong></td>
<td style="text-align: left;"></td>
<td style="text-align: center;">原始观测场</td>
<td style="text-align: left;"><a href="https://www.showdoc.com.cn/meteva?page_id=3975600815874861">网格数据格式</a></td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>grd_label</font></strong></td>
<td style="text-align: left;"></td>
<td style="text-align: center;">目标标记场</td>
<td style="text-align: left;"><a href="https://www.showdoc.com.cn/meteva?page_id=3975600815874861">网格数据格式</a></td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>grd_features</font></strong></td>
<td style="text-align: left;"></td>
<td style="text-align: center;">标记的目标以及属性</td>
<td style="text-align: left;">字典</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>label_count</font></strong></td>
<td style="text-align: left;">网格场中目标的个数</td>
<td style="text-align: center;">整数</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>area</font></strong></td>
<td style="text-align: left;">各目标的面积</td>
<td style="text-align: center;">实数列表，长度=label_count, 其中一个目标的面积 = 目标格点数 <em> 径向格距 </em> 纬向格距</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>1</font></strong></td>
<td style="text-align: left;">观测场中1号目标的坐标</td>
<td style="text-align: center;">numpy数组的位置索引，其中记录了1号目标每一个点的坐标，具体内容为np.where(grd_ob_labeled.values.squeeze()==1)返回结果</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>2</font></strong></td>
<td style="text-align: left;">观测场中2号目标的坐标</td>
<td style="text-align: center;">numpy数组的位置索引,tuple</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=green size=3>...</font></strong></td>
<td style="text-align: left;">...</td>
<td style="text-align: center;">...</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">
look_ff_ob = mem.mode.feature_finder_and_merge(grd_ob,smooth = 5,threshold = 5,minsize = 10,
near_dis = 200,near_rate = 0.3)
look_ff_fo_03 = mem.mode.feature_finder_and_merge(grd_fo_03,smooth = 5,threshold = 5,minsize = 10,
near_dis = 100,near_rate = 0.5)
look_ff_fo_27 = mem.mode.feature_finder_and_merge(grd_fo_27,smooth = 5,threshold = 5,minsize = 10,
near_dis = 100,near_rate = 0.5)
print(&quot;*** 目标识别和合并完成 ***&quot;)</code></pre>
<pre><code>*** 目标识别和合并完成 ***</code></pre>
<p>改进后的MODE方法，一个观测同时多个预报场建立可对比的匹配关系，因此，有必要在同一张图上，同时显示一个观测场和多个预报场，以及场中的目标识别结果</p>
<h1>识别结果绘制</h1>
<h2>绘制要素场</h2>
<p><font face="黑体" color=Blue size=3><strong>plot_value_list(look_list, cmap ="rain_24h",clevs = None,save_path,show = True)</strong></font> </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_list</font></strong></td>
<td style="text-align: left;">识别合并或匹配后所得观测或者预报的目标对象的列表</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>cmap</font></strong></td>
<td style="text-align: left;">网格预报数据绘图colorbar的cmap参数，字符或者matplotlib 的cm对象。 缺省值”rainbow”对应彩虹色方案。 当cmap为字符时支持两种情况，一种是matplotlib 能够识别的字符串，例如”rainbow”,”bwr”等等，具体可以参考https://matplotlib.org/tutorials/colors/colormaps.html， 另外则是meteva定义的一些和预报检验有关的一些配色方案，具体可参考meb.def_cmap_clevs函数说明</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>clevs</font></strong></td>
<td style="text-align: left;">网格数据绘图colorbar的等级参数， 取值为列表或一维nump数值</td>
</tr>
<tr>
<td style="text-align: left;"><strong>save_path</strong></td>
<td style="text-align: left;">图片保存路径，缺省时不输出图片，而是以默认绘图窗口形式展示</td>
</tr>
<tr>
<td style="text-align: left;"><strong>show</strong></td>
<td style="text-align: left;">是否在屏幕显示图片，如果save_path 和save_dir 为都None时，程序内部会自动将show设置True</td>
</tr>
<tr>
<td style="text-align: left;"><strong>dpi</strong></td>
<td style="text-align: left;">图片的分辨率，效果同matplotlib的dpi参数</td>
</tr>
<tr>
<td style="text-align: left;"><strong>sup_fontsize</strong></td>
<td style="text-align: left;">图片标题的字体大小，其它字体将根据标题字体大小自动设置，其中坐标轴字体大小 = sup_fontsize <em> 0.9, 坐标刻度的字体大小 = sup_fontsize </em> 0.8</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>result</font></strong></td>
<td style="text-align: left;">无</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">mem.mode.plot_value_list([look_ff_ob,look_ff_fo_03,look_ff_fo_27],sup_fontsize = 8,
cmap = &quot;rain_3h&quot;,clevs = np.arange(0,61,3))</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=d18b70cc669baafccede416bb7cfc3e9" alt="" /></p>
<h2>绘制目标编号</h2>
<p><font face="黑体" color=Blue size=3><strong>plot_label_list(look_list,ncol = None,save_path = None,show = False)</strong></font> </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_list</font></strong></td>
<td style="text-align: left;">识别合并或匹配后所得观测或者预报的目标对象的列表</td>
</tr>
<tr>
<td style="text-align: left;"><strong>ncol</strong></td>
<td style="text-align: left;">子图列数，缺省时自动确定</td>
</tr>
<tr>
<td style="text-align: left;"><strong>save_path</strong></td>
<td style="text-align: left;">图片保存路径，缺省时不输出图片，而是以默认绘图窗口形式展示</td>
</tr>
<tr>
<td style="text-align: left;"><strong>show</strong></td>
<td style="text-align: left;">是否在屏幕显示图片，如果save_path 和save_dir 为都None时，程序内部会自动将show设置True</td>
</tr>
<tr>
<td style="text-align: left;"><strong>dpi</strong></td>
<td style="text-align: left;">图片的分辨率，效果同matplotlib的dpi参数</td>
</tr>
<tr>
<td style="text-align: left;"><strong>sup_fontsize</strong></td>
<td style="text-align: left;">图片标题的字体大小，其它字体将根据标题字体大小自动设置，其中坐标轴字体大小 = sup_fontsize <em> 0.9, 坐标刻度的字体大小 = sup_fontsize </em> 0.8</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>result</font></strong></td>
<td style="text-align: left;">无</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">mem.mode.plot_label_list([look_ff_ob,look_ff_fo_03,look_ff_fo_27],sup_fontsize = 8,ncol = 3)</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=a7cc30da619758f325e48650056522fd" alt="" /></p>
<p>在上面的示例中，经过初步的识别后，观测目标中包含多个联通区。预报的目标的识别合并时，有意将合并的临近尺度设得更低，临近比例阈值设得更高一些，这样目标更加零散，这样有利于后续匹配步骤中，进一步按照观测的目标的分布，重新组合预报目标，提升匹配的效果。</p>
<h2>绘制要素和编号</h2>
<p><font face="黑体" color=Blue size=3><strong>plot_value_and_label_list(look_list, cmap ="rain_24h",clevs = None,save_path,show = True)</strong></font> </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_list</font></strong></td>
<td style="text-align: left;">识别合并或匹配后所得观测或者预报的目标对象的列表</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>cmap</font></strong></td>
<td style="text-align: left;">网格预报数据绘图colorbar的cmap参数，字符或者matplotlib 的cm对象。 缺省值”rainbow”对应彩虹色方案。 当cmap为字符时支持两种情况，一种是matplotlib 能够识别的字符串，例如”rainbow”,”bwr”等等，具体可以参考https://matplotlib.org/tutorials/colors/colormaps.html， 另外则是meteva定义的一些和预报检验有关的一些配色方案，具体可参考meb.def_cmap_clevs函数说明</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>clevs</font></strong></td>
<td style="text-align: left;">网格数据绘图colorbar的等级参数， 取值为列表或一维nump数值</td>
</tr>
<tr>
<td style="text-align: left;"><strong>save_path</strong></td>
<td style="text-align: left;">图片保存路径，缺省时不输出图片，而是以默认绘图窗口形式展示</td>
</tr>
<tr>
<td style="text-align: left;"><strong>show</strong></td>
<td style="text-align: left;">是否在屏幕显示图片，如果save_path 和save_dir 为都None时，程序内部会自动将show设置True</td>
</tr>
<tr>
<td style="text-align: left;"><strong>dpi</strong></td>
<td style="text-align: left;">图片的分辨率，效果同matplotlib的dpi参数</td>
</tr>
<tr>
<td style="text-align: left;"><strong>sup_fontsize</strong></td>
<td style="text-align: left;">图片标题的字体大小，其它字体将根据标题字体大小自动设置，其中坐标轴字体大小 = sup_fontsize <em> 0.9, 坐标刻度的字体大小 = sup_fontsize </em> 0.8</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>result</font></strong></td>
<td style="text-align: left;">无</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">mem.mode.plot_value_and_label_list([look_ff_ob,look_ff_fo_03,look_ff_fo_27],sup_fontsize = 8,
cmap = &quot;rain_3h&quot;,clevs = np.arange(50))</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=aa0703d102e94914bb7fd5e5dff60740" alt="" /></p>
<h1>预报向观测匹配</h1>
<p><font face="黑体" color=Blue size=3><strong>unimatch(look_ob,look_fo,cover_dis,cover_rate)</strong></font><br />
将预报场中的目标向观测的目标匹配，匹配上的预报目标会和观测目标的编号相同，匹配过程中，预报目标会根据位置分布情况进一步合并，以提升和观测目标的重合度。<br />
具体的匹配算法请参考文献 ：<br />
刘凑华, 牛若芸. 基于目标的降水检验方法及应用[J]. 气象, 2013, 39(6):681-690</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: center;">说明</th>
<th style="text-align: left;">备注</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_ob</font></strong></td>
<td style="text-align: center;">feature_finder_and_merge函数结果识别出的观测目标</td>
<td style="text-align: left;">数据为dict格式</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_fo</font></strong></td>
<td style="text-align: center;">feature_finder_and_merge函数结果识别出的预报目标</td>
<td style="text-align: left;">数据为dict格式</td>
</tr>
<tr>
<td style="text-align: left;"><strong>cover_dis</strong></td>
<td style="text-align: center;">邻域尺度，单位是km</td>
<td style="text-align: left;">观测预报在位置上通常不是完全重合的，该参数是计算观测预报目标重合度时允许的距离偏差</td>
</tr>
<tr>
<td style="text-align: left;"><strong>cover_rate</strong></td>
<td style="text-align: center;">匹配度（即在位置无法范围内的重合度）阈值</td>
<td style="text-align: left;">观测和预报目标匹配度大于该阈值的记为成功匹配</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>return</font></strong></td>
<td style="text-align: center;">和feature_finder_and_merge识别的结果形式一致的预报目标变量</td>
<td style="text-align: left;">数据为dict格式</td>
</tr>
</tbody>
</table>
<pre><code class="language-python">matched_fo_03 = mem.mode.unimatch(look_ff_ob,look_ff_fo_03,cover_dis=150,cover_rate=0.5)
matched_fo_27 = mem.mode.unimatch(look_ff_ob,look_ff_fo_27,cover_dis=150,cover_rate=0.5)
mem.mode.plot_label_list([look_ff_ob,matched_fo_03,matched_fo_27],sup_fontsize =8,ncol = 3)</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=b12018d90bce8523e019ceb23721a22d" alt="" /></p>
<p>在上面的示例中，预报目标识别的结果非常零散，但经过和观测的匹配过程，预报目标重新聚合再一起，从不同时效预报和观测的匹配结果看出，它们和观测匹配关系是合理的。</p>
<h1>匹配结果合并</h1>
<p><font face="黑体" color=Blue size=5><strong>unimerge(look_ob,look_fo)</strong></font><br />
经过unimatch匹配后，观测的目标和预报的目标变量是分开的。为便于开展后续检验，采用unimerge函数将观测和预报目标集成为数据结构和mem.mode.merge_force返回结果相同的一个变量当中。</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: center;">说明</th>
<th style="text-align: left;">备注</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_ob</font></strong></td>
<td style="text-align: center;">feature_finder_and_merge 识别合并后的观测目标</td>
<td style="text-align: left;">数据为dict格式</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>look_fo</font></strong></td>
<td style="text-align: center;">unimatch 匹配后的预报目标</td>
<td style="text-align: left;">数据为dict格式</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>result</font></strong></td>
<td style="text-align: center;">匹配后的结果</td>
<td style="text-align: left;">数据结构同mem.mode.merge_force的结果相同</td>
</tr>
</tbody>
</table>
<pre><code class="language-python">look_merge_27 = mem.mode.unimerge(look_ff_ob,matched_fo_27)</code></pre>
<h1>检验综合函数</h1>
<p><font face="黑体" color=Blue size=3><strong>operate(grd_ob, grd_fo, smooth, threshold,minsize,compare = ">=",match_method = centmatch,near_dis =100,near_rate =0.3,cover_dis = 100,cover_rate =0.5,summary =True,save_dir = None,cmap = "rain_24h",clevs = None,show = False）</strong></font><br />
根据输入的观测和预报网格数据， 完成目标识别、匹配、合并和检验的完整步骤。 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: center;">说明</th>
<th style="text-align: left;">备注</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>grd_ob</font></strong></td>
<td style="text-align: center;">网格观测数据</td>
<td style="text-align: left;"><a href="https://www.showdoc.com.cn/meteva?page_id=3975600815874861">网格数据格式</a></td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>grd_fo</font></strong></td>
<td style="text-align: center;">网格预报数据</td>
<td style="text-align: left;"><a href="https://www.showdoc.com.cn/meteva?page_id=3975600815874861">网格数据格式</a></td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>smooth</font></strong></td>
<td style="text-align: center;">平滑系数,采用圆盘形的卷积内核对网格数据进行平滑，smooth相当圆盘的半径</td>
<td style="text-align: left;">可以为单个整数或列表形式，当其为单个整数时，观测和预报采用相同的参数，否则相应采用不同的参数，如smooth = [5, 4]。平滑系数为0时不做平滑</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>threshold</font></strong></td>
<td style="text-align: center;">雨量阈值，平滑后低于阈值的部分会被置0</td>
<td style="text-align: left;">可以为单个实数或列表形式，当其为单个实数时，观测和预报采用相同的参数，否则相应采用不同的参数，，如threshold = [9,10]</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>minsize</font></strong></td>
<td style="text-align: center;">目标的最小面积,面积小于该阈值的目标将会被删除</td>
<td style="text-align: left;">可以为单个整数或列表形式，当其为单个整数时，观测和预报采用相同的参数，否则相应采用不同的参数，<br>如minsize = [5, 4]，此列表中的数值表示的是格点数。</td>
</tr>
<tr>
<td style="text-align: left;"><strong>compare</strong></td>
<td style="text-align: center;">比较方法，可选项包括">=",">","<=","<"</td>
<td style="text-align: left;">当关注目标是要素场大于阈值时,采用默认参数，如识别暴雨目标时选择参数">=",若关注对象是要素值小于阈值时，则可选取"<="，如根据能见度场识别大雾目标</td>
</tr>
<tr>
<td style="text-align: left;"><strong>match_method</strong></td>
<td style="text-align: center;">目标匹配采用的方法</td>
<td style="text-align: left;">可选项包括 mem.mode.unimatch,<br>mem.mode.centmatch,<br>mem.mode.deltamm,<br>mem.mode.minboundmatch</td>
</tr>
<tr>
<td style="text-align: left;"><strong>near_dis</strong></td>
<td style="text-align: center;">观测和预报目标识别合并函数feature_finder_and_merge使用到的邻域尺度参数</td>
<td style="text-align: left;">可以为单个整数或列表形式，当其为单个整数时，观测和预报采用相同的参数，否则相应采用不同的参数，<br>如near_dis = [200,200]，单位km</td>
</tr>
<tr>
<td style="text-align: left;"><strong>near_rate</strong></td>
<td style="text-align: center;">观测和预报目标识别合并函数feature_finder_and_merge使用到的邻域比例参数</td>
<td style="text-align: left;">可以为单个整数或列表形式，当其为单个整数时，观测和预报采用相同的参数，否则相应采用不同的参数，<br>如near_rate = [0.3, 0.3]</td>
</tr>
<tr>
<td style="text-align: left;"><strong>cover_dis</strong></td>
<td style="text-align: center;">观测和预报目标匹配函数unimatch使用到的邻域尺度参数</td>
<td style="text-align: left;">实时，单位km</td>
</tr>
<tr>
<td style="text-align: left;"><strong>cover_dis</strong></td>
<td style="text-align: center;">观测和预报目标匹配函数unimatch使用到的覆盖比例参数</td>
<td style="text-align: left;">实时，取值范围是0-1</td>
</tr>
<tr>
<td style="text-align: left;"><strong>summary</strong></td>
<td style="text-align: center;">当该参数为True时返回相对重要的目标属性，否则返回完整的目标属性</td>
<td style="text-align: left;">布尔型</td>
</tr>
<tr>
<td style="text-align: left;"><strong>save_dir</strong></td>
<td style="text-align: center;">计算结果的输出目录,输出的结果包括4个部分,第一个部分是以netcdf格式输出预报和观测编号，第二部分是以图片形式输出输出原始场、平滑场和编号场，第三部分是以json格式输出检验结果概要，第四部分以图片形式输出检验结果概要</td>
<td style="text-align: left;">字符串</td>
</tr>
<tr>
<td style="text-align: left;"><strong>show</strong></td>
<td style="text-align: center;">是否在屏幕显示图片，如果save_dir 为都None时，程序内部会自动将show设置True</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>return</font></strong></td>
<td style="text-align: center;">成功匹配的目标的属性及检验结果</td>
<td style="text-align: left;">字典，内容如feature_merged_analyzer函数的返回结果</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">result = mem.mode.operate(grd_ob,grd_fo_27,smooth=5,threshold=5,minsize=5,match_method=mem.mode.unimatch,near_dis=[200,100],
near_rate=[0.3,0.5],cover_dis=150,cover_rate=0.5,
cmap = &quot;rain_3h&quot;,clevs = np.arange(0,61,3),show = True,save_dir = r&quot;H:\test_data\output\method\space\obv&quot;)</code></pre>
<pre><code>目标编号场已经输出至 H:\test_data\output\method\space\obv/label/2020072508.027.nc
图片已保存至H:\test_data\output\method\space\obv/png/2020072508.027.png</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=b9ef0ecc5f829f4de5045ad7f56246b1" alt="" /></p>
<pre><code>目标属性检验结果已经输出至 H:\test_data\output\method\space\obv/json/2020072508.027.txt
图片已保存至H:\test_data\output\method\space\obv/table/2020072508.027.png</code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=8f42c9b8d1d809ce50756a54127631a0" alt="" /></p>
<h1>长序列检验应用</h1>
<p>以上的程序都只是针对单个预报观测场进行检验，很多情况下，只有进行长时间序列的预报偏差进行统计分析才能得出有价值的结论。为此在以下提供一套将改进后的MODE用于长序列检验的示例程序和一些配套工具</p>
<pre><code class="language-python">#以下是批量读取观测和预报场数据并完成目标识别匹配和结果输出
dir_ob = r'H:\test_data\input\mem\mode\ob\rain03\YYMMDDHH.000.nc'
dir_fo = r'H:\test_data\input\mem\mode\ec\rain03\YYMMDDHH.TTT.nc'
save_dir = r&quot;H:\test_data\output\method\space\obv&quot;
time_s = datetime.datetime(2020,7,1,8,0)
time_e = datetime.datetime(2020,8,1,8,0)
time_fo = time_s
feature_summary_list = []
while time_fo &lt; time_e:
for dh in range(3,73,3):
time_ob = time_fo + datetime.timedelta(hours = dh)
path_ob = meb.get_path(dir_ob,time_ob)
grd_ob1 = meb.read_griddata_from_nc(path_ob,time = time_ob,dtime = 0,level =0,data_name = &quot;ob&quot;)
path_fo = meb.get_path(dir_fo,time_fo,dh)
grd_fo1 = meb.read_griddata_from_nc(path_fo,time = time_fo,dtime = dh,level =0,data_name = &quot;ECMWF&quot;)
if grd_ob1 is not None and grd_fo1 is not None:
feature_summary = mem.mode.operate(grd_ob1,grd_fo1,smooth = 5,threshold=5,minsize=5,
cmap = &quot;rain_3h&quot;,clevs = np.arange(0,61,3),
match_method=mem.mode.unimatch,near_dis=[200,100],
near_rate=[0.3,0.5],cover_dis=150,cover_rate=0.5,save_dir=save_dir)
feature_summary_list.append(feature_summary)
time_fo = time_fo + datetime.timedelta(hours = 12)</code></pre>
<pre><code>目标编号场已经输出至 H:\test_data\output\method\space\obv/label/2020070108.003.nc
图片已保存至H:\test_data\output\method\space\obv/png/2020070108.003.png
目标属性检验结果已经输出至 H:\test_data\output\method\space\obv/json/2020070108.003.txt
图片已保存至H:\test_data\output\method\space\obv/table/2020070108.003.png
目标编号场已经输出至 H:\test_data\output\method\space\obv/label/2020070108.006.nc
图片已保存至H:\test_data\output\method\space\obv/png/2020070108.006.png
目标属性检验结果已经输出至 H:\test_data\output\method\space\obv/json/2020070108.006.txt
图片已保存至H:\test_data\output\method\space\obv/table/2020070108.006.png
。。。</code></pre>
<h2>从文件中导入</h2>
<p><font face="黑体" color=Blue size=3><strong>load_feature_summary(filename）</strong></font><br />
从json格式文件中导入单个起报时刻预报时效的预报观测目标检验概要信息。 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>filename</font></strong></td>
<td style="text-align: left;">字符串形式文件路径，文件内容是json格式存放的检验概要信息</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>return</font></strong></td>
<td style="text-align: left;">检验概要信息，内容如feature_merged_analyzer函数在summary参数取True情况下的返回结果</td>
</tr>
</tbody>
</table>
<h2>从文件夹中批量导入</h2>
<p><font face="黑体" color=Blue size=3><strong>load_feature_summary_list(dir）</strong></font><br />
从json格式文件中导入单个起报时刻预报时效的预报观测目标检验概要信息。 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>dir</font></strong></td>
<td style="text-align: left;">字符串形式文件夹路径，文件夹里存放的是全部json格式文件，存放内容是检验概要信息</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>return</font></strong></td>
<td style="text-align: left;">检验概要信息列表，元素内容如feature_merged_analyzer函数在summary参数取True情况下的返回结果</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">feature_summary_list = mem.mode.load_feature_summary_list(r&quot;H:\test_data\output\method\space\obv\json&quot;)
print(&quot;*****批量导入feature_summary_list成功******&quot;) </code></pre>
<pre><code>*****批量导入feature_summary_list成功******</code></pre>
<h2>检验概要转DataFrame</h2>
<p><font face="黑体" color=Blue size=3><strong>features_list_to_df(feature_summary_list）</strong></font><br />
将feature_merged_analyzer函数返回的预报观测目标检验概要（字典）数据的列表转换成pandas的DataFrame格式，以方便进行统计分析。 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>feature_summary_list</font></strong></td>
<td style="text-align: left;">列表，其元素为feature_merged_analyzer函数返回结果的列表，其中summary需选择True</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=Blue size=5>return</font></strong></td>
<td style="text-align: left;">元组，包含4个元素，都为DataFrame格式，<br>其中元素1的列名称依次为'level', 'time', dtime', 'id', 'lon', 'lat', 'Hits', 'Misses', 'False alarms', 'Correct negatives';其中"id","lon","lat"列为缺省值999999.<br> 元素2的列名称依次为'level', 'time', 'dtime', 'id', 'lon', 'lat', 'ob_lenghts_MajorAxis','ob_lenghts_MinorAxis', 'ob_window_x0', 'ob_window_y0', 'ob_window_x1','ob_window_y1', 'ob_centroid_x', 'ob_centroid_y', 'ob_area','ob_intensity','ob_angle', 'fo_lenghts_MajorAxis', 'fo_lenghts_MinorAxis', 'fo_window_x0', 'fo_window_y0', 'fo_window_x1', 'fo_window_y1','fo_centroid_x', 'fo_centroid_y', 'fo_area', 'fo_intensity','fo_angle', 'cent_dist', 'angle_diff', 'area_ratio', 'int_area', 'bearing','bdelta', 'haus', 'medMiss', 'medFalseAlarm', 'msdMiss', 'msdFalseAlarm', 'ph', 'fom', 'minsep' ,其中"id","lon","lat"列分别为观测目标的编号，质心经度和质心纬度</td>
</tr>
</tbody>
</table>
<p><strong>调用示例</strong></p>
<pre><code class="language-python">df_hmfc,df_feature,_,_ = mem.mode.features_list_to_df(feature_summary_list)
print(&quot;df_hmfc:&quot;)
print(df_hmfc)
print(&quot;--------------&quot;)
print(df_feature.columns)</code></pre>
<pre><code>df_hmfc:
level time dtime id lon lat Hits Misses \
0 0 2020-07-01 08:00:00 3 999999 999999 999999 4 1
1 0 2020-07-01 08:00:00 6 999999 999999 999999 4 3
2 0 2020-07-01 08:00:00 9 999999 999999 999999 3 9
3 0 2020-07-01 08:00:00 12 999999 999999 999999 2 9
4 0 2020-07-01 08:00:00 15 999999 999999 999999 2 8
... ... ... ... ... ... ... ... ...
1483 0 2020-07-31 20:00:00 60 999999 999999 999999 7 6
1484 0 2020-07-31 20:00:00 63 999999 999999 999999 5 5
1485 0 2020-07-31 20:00:00 66 999999 999999 999999 2 2
1486 0 2020-07-31 20:00:00 69 999999 999999 999999 5 7
1487 0 2020-07-31 20:00:00 72 999999 999999 999999 8 9
False alarms Correct negatives
0 2 79
1 2 85
2 4 113
3 4 209
4 3 180
... ... ...
1483 4 120
1484 4 92
1485 10 184
1486 9 151
1487 7 241
[1488 rows x 10 columns]
--------------
Index(['level', 'time', 'dtime', 'id', 'lon', 'lat', 'ob_lenghts_MajorAxis',
'ob_lenghts_MinorAxis', 'ob_window_x0', 'ob_window_y0', 'ob_window_x1',
'ob_window_y1', 'ob_centroid_x', 'ob_centroid_y', 'ob_area',
'ob_intensity','ob_angle', 'fo_lenghts_MajorAxis', 'fo_lenghts_MinorAxis',
'fo_window_x0', 'fo_window_y0', 'fo_window_x1', 'fo_window_y1',
'fo_centroid_x', 'fo_centroid_y', 'fo_area', 'fo_intensity','fo_angle',
'ob_intensity_0', 'fo_intensity_0', 'ob_intensity_5', 'fo_intensity_5',
'ob_intensity_10', 'fo_intensity_10', 'ob_intensity_25',
'fo_intensity_25', 'ob_intensity_50', 'fo_intensity_50',
'ob_intensity_75', 'fo_intensity_75', 'ob_intensity_90',
'fo_intensity_90', 'ob_intensity_95', 'fo_intensity_95',
'ob_intensity_100', 'fo_intensity_100', 'cent_dist', 'angle_diff',
'area_ratio', 'int_area', 'bearing', 'bdelta', 'haus', 'medMiss',
'medFalseAlarm', 'msdMiss', 'msdFalseAlarm', 'ph', 'fom', 'minsep'],
dtype='object')</code></pre>
<pre><code class="language-python">#提取其中的中位数强度长度信息
intensity_50 = meb.sele_by_para(df_feature,member = [&quot;ob_intensity_50&quot;,&quot;fo_intensity_50&quot;])
print(intensity_50.iloc[0:10,:]) #打印前10行</code></pre>
<pre><code> level time dtime id lon lat \
0 0 2020-07-01 08:00:00 3 1 132.238746 29.386656
1 0 2020-07-01 08:00:00 3 2 107.618715 29.553073
2 0 2020-07-01 08:00:00 3 3 102.811404 24.855263
3 0 2020-07-01 08:00:00 3 4 117.340426 27.781915
4 0 2020-07-01 08:00:00 6 1 131.916587 28.759533
5 0 2020-07-01 08:00:00 6 3 107.002907 29.837209
6 0 2020-07-01 08:00:00 6 4 84.888298 15.404255
7 0 2020-07-01 08:00:00 6 5 97.244048 24.589286
8 0 2020-07-01 08:00:00 9 1 131.577514 28.161313
9 0 2020-07-01 08:00:00 9 12 96.571429 25.464286
ob_intensity_50 fo_intensity_50
0 13.300 8.499
1 4.013 5.121
2 2.766 7.545
3 5.951 7.305
4 13.475 7.137
5 5.806 4.021
6 6.799 5.503
7 7.366 4.223
8 10.625 6.466
9 12.673 3.044 </code></pre>
<h3>检验分析示例</h3>
<pre><code class="language-python">#应用meteva.product层的函数进行批量检验分析
result = mpd.score(intensity_50,mem.ob_fo_mean,g = &quot;ob_time&quot;,plot = &quot;line&quot;,ylabel = &quot;目标降水量中位数强度&quot;,
sup_fontsize = 16,grid = True) </code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=b686579d396af1471f465a2c74f9eed5" alt="" /></p>
<pre><code class="language-python">major_axis = meb.sele_by_para(df_feature,member = [&quot;ob_lenghts_MajorAxis&quot;,&quot;fo_lenghts_MajorAxis&quot;])
#提取部分地区目标进行检验，按观测目标的质心经纬度位置进行判别
result = mpd.score(major_axis,mem.ob_fo_mean,s ={&quot;lon&quot;:[100,125],&quot;lat&quot;:[20,40]},g = &quot;dtime&quot;,plot = &quot;line&quot;,
sup_fontsize = 16,grid = True,xlabel =&quot;时效（单位：小时）&quot;,ylabel = &quot;主轴长度（单位：100km）&quot;,
title =&quot;2020年7月中东部地区（100-125°E,20-40°N）\nECMWF模式逐3小时预报的降水目标主轴的平均长度同观测的对比&quot;) </code></pre>
<p><img src="https://www.showdoc.com.cn/server/api/attachment/visitFile?sign=6f30affcea1dfcb261c11dbb74820de7" alt="" /></p>