一种基于潮汐和径流的河口咸潮预报方法

摘要

本发明公开了一种基于潮汐和径流的河口咸潮预报方法，具体包括以下步骤：步骤A00：选取参数；步骤B00：量化参数；步骤C00：预估盐度和潮差统计模型；步骤D00：初次预估盐度和径流量统计模型；步骤E00：初次预估盐度统计模型；步骤F00：修正盐度统计模型；步骤G00：得到盐度预测预报模型。本基于潮汐和径流的河口咸潮预报方法建立了基于含氯度与径流、潮差关系的咸潮上溯预测模型，可提前1-7天预测，实现了中短期咸潮预测预报；预测模型只有一个主要未知参数径流，把咸潮预测复杂问题简单化，具有简便实用的突出特点。

说明书

技术领域

本发明涉及一种预报方法，尤其涉及的是一种基于潮汐和径流的河口咸潮预报方法。

背景技术

在建模过程中，因对咸潮入侵的动力过程缺少深入了解，现有的河口咸潮入侵统计模型建立的盐度与潮差和径流量关系不够理想，预测结果与实测资料比较误差较大，难于在涉水行业内得到应用。在数值预报模式中，给出的数值和边界条件复杂，计算量大，难于快速高效地预报河口取水口的盐度。

因此，现有技术还有待于改进和发展。

发明内容

本发明的目的在于提供一种基于潮汐和径流的河口咸潮预报方法，旨在解决现有的咸潮预报预见期短，预测咸潮边界条件复杂、参数多和预测模型计算复杂的问题，紧密结合河口咸潮入侵动力过程，建立河口咸潮入侵统计预报模型，根据预测预报结果采取避咸取淡措施，保障工农业生产和饮用水供水安全。

本发明的技术方案如下：

一种基于潮汐和径流的河口咸潮预报方法，其中，具体包括以下步骤：

步骤A00：选取参数：选取长时段的预测点盐度日特征值、口门区某站点的潮汐日特征值和上游来水量；

步骤B00：量化参数：选择日最大潮差值TR为潮差，径流量Q等于上游日总来水量；

步骤C00：预估盐度和潮差统计模型：通过分析盐度和潮差的散点分布形式，建立盐度和潮差拟合关系式S0（TR）；

步骤D00：初次预估盐度和径流量统计模型：将盐度差DS进行潮形归一化, 得到归一化盐度差SDS，对归一化盐度差SDS和径流量Q进行拟合，得到径流量盐度差拟合函数SDS_M，对径流量盐度差拟合函数SDS_M去归一化, 得到盐度差DS与径流量Q拟合关系式S0（Q）；

步骤E00：初次预估盐度统计模型：根据盐度和潮差拟合关系式S0（TR）与盐度差和径流量拟合关系式S0（Q），初次得出预估盐度统计模型；

步骤F00：修正盐度统计模型：消除径流量对实测盐度值S的影响，得到滤去径流量影响的盐度值S1，将S1与潮差作分析，得到修正的盐度潮差拟合关系式S1（TR）；将实测盐度值S消去修正的盐度潮差拟合关系式S1（TR）, 得到与径流量Q相关的盐度差，并作归一化处理得到归一化盐度差SDS1，将归一化盐度差SDS1和径流量Q进行拟合，得到径流量盐度差拟合函数SDS1_M，对径流量盐度差拟合函数SDS1_M去归一化, 得到修正盐度差与径流量拟合关系S1（Q）；

步骤G00：得到盐度预测预报模型：根据修正后盐度和潮差的关系式S1（TR）与修正后盐度差和径流量的关系式S1（Q），得到盐度与潮差和径流量Q的关系式S（Y）；

所述日特征值为日平均值、日最大值和日最小值，所以盐度日特征值是指盐度日平均值、日最大值和日最小值；潮汐特征值是指潮汐日平均值、日最大值和日最小值。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述步骤C00中，通过观察盐度和潮差的散点分布形式，建立盐度和潮差关系式S0（TR）：                                                （1）,（1）式中S0为盐度与潮差预估值，TR为潮差，、为拟合系数。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述步骤D00中，所述归一化盐度差的表达式SDS为：（2），（3），其中S为盐度实测值，DS为盐度差，SDS为归一化盐度差。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述径流量盐度差拟合函数SDS_M: （4），（4）中Q为径流量，SDS_M为归一化盐度差与径流量的拟合值，、、为拟合系数；所述盐度差与径流量拟合关系式S0(Q)：（5）。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述步骤F00中，滤去径流量影响的盐度S1的表达式为： ，其中S为实测盐度值，S(Q)为盐度差与径流量拟合关系式；修正后的盐度潮差拟合关系式为：（6），（6）式中为盐度与潮差修正值, 、为拟合系数。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述归一化盐度差SDS1的表达式为：（7）。

所述的咸潮预测预报方法，其中，将归一化盐度差SDS1与径流量Q采用反比例函数和最小二乘法拟合, 得到修正归一化盐度差SDS1_M与径流量Q的拟合关系式： （8），其中、、C0为拟合系数；修正后盐度差与径流量拟合关系:（9）。

所述的基于潮汐和径流的河口咸潮预报方法，其中，所述步骤G00中，所述盐度与潮差和径流量Q的关系式S（Y）：，，S(Y)=S1(TR)+ S1(TR)×S1(Q)= A*exp(B*TR)+ A*exp(B*TR)*[ C/(Q-D)-E] （10）；其中，S(Y)为预测盐度，TR为潮差，Q为流量，A、B、C、D、E均为拟合系数。

本发明的有益效果：本发明通过提供一种基于潮汐和径流的河口咸潮预报方法，本咸潮预测预报方法建立了基于含氯度与径流、潮差关系的咸潮上溯预测模型，可提前1-7天预测，实现了中短期咸潮预测预报；预测模型只有一个主要未知参数径流，把咸潮预测复杂问题简单化，具有简便实用的突出特点；特别适合使用单位各级管理、技术人员包括一线生产人员使用；为口门受咸潮影响地区避咸取淡工作提供了技术支撑。

附图说明

图1是本发明中基于潮汐和径流的河口咸潮预报方法的步骤流程图。

图2是本发明中原始资料处理后盐度、潮差和径流量的逐日变化图。

图3是本发明中潮差和盐度相位差相关系数图。

图4是本发明中除去径流量影响的日平均盐度与潮差分布和拟合曲线图。

图5是本发明中，修正归一化盐度差与径流量分布和拟合曲线图。

图6是本发明中2009年10月-12月平岗日平均含氯度预测预报结果与实测日平均含氯度结果对比图。

图7是本发明中2009年10月-12月平岗日最大含氯度预测预报结果与实测日最大含氯度结果对比图。

图8是本发明中2011年1月-2月平岗日平均含氯度预测预报结果与实测日平均含氯度结果对比图。

图9是本发明中2011年1月-2月平岗日最大含氯度预测预报结果与实测日最大含氯度结果对比图。

图10是本发明中2011年1月-2月广昌日平均含氯度预测预报结果与实测日平均含氯度结果对比图。

图11是本发明中2011年1月-2月广昌日最大含氯度预测预报结果与实测日最大含氯度结果对比图。

具体实施方式

为使本发明的目的、技术方案及优点更加清楚、明确，以下参照附图并举实施例对本发明进一步详细说明。

如图1所示，所述基于潮汐和径流的河口咸潮预报方法具体包括以下步骤：

步骤A00：选取参数：选取长时段的预测点盐度日特征值、口门区某站点的潮汐日特征值和上游来水量；

步骤B00：量化参数：选择日最大潮差值TR为潮差，径流量Q等于上游日总来水量；

步骤C00：预估盐度和潮差统计模型：通过分析盐度和潮差的散点分布形式，建立盐度和潮差拟合关系式S0（TR）；

步骤D00：初次预估盐度和径流量统计模型：将盐度差DS进行潮形归一化, 得到归一化盐度差SDS，对归一化盐度差SDS和径流量Q进行拟合，得到径流量盐度差拟合函数SDS_M，对径流量盐度差拟合函数SDS_M去归一化, 得到盐度差DS与径流量Q拟合关系式S0（Q）；

步骤E00：初次预估盐度统计模型：根据盐度和潮差拟合关系式S0（TR）与盐度差和径流量拟合关系式S0（Q），初次得出预估盐度统计模型；

步骤F00：修正盐度统计模型：消除径流量对实测盐度值S的影响，得到滤去径流量影响的盐度值S1，将S1与潮差作分析，得到修正的盐度潮差拟合关系式S1（TR）；将实测盐度值S消去修正的盐度潮差拟合关系式S1（TR）, 得到与径流量Q相关的盐度差，并作归一化处理得到归一化盐度差SDS1，将归一化盐度差SDS1和径流量Q进行拟合，得到径流量盐度差拟合函数SDS1_M，对径流量盐度差拟合函数SDS1_M去归一化, 得到修正盐度差与径流量拟合关系S1（Q）；

步骤G00：得到盐度预测预报模型：根据修正后盐度和潮差的关系式S1（TR）与修正后盐度差和径流量的关系式S1（Q），得到盐度与潮差和径流量Q的关系式S（Y）；

所述日特征值为日平均值、日最大值和日最小值，所以盐度日特征值是指盐度日平均值、日最大值和日最小值；潮汐特征值是指潮汐日平均值、日最大值和日最小值。

如图2-5所示，所述步骤C00中，因河口盐水入侵在半月时间尺度内主要受潮汐影响, 潮汐的强弱可以用潮差描述，故需先建立盐度和潮差的拟合关系。通过观察盐度和潮差的散点分布形式，可看到盐度与潮差总体上具有指数关系，因此可建立如下形式的盐度和潮差关系式：（1）（（1）式中S0为盐度与潮差预估值，TR为潮差，、为拟合系数）。

所述步骤D00中，观测盐度减去盐度与潮差拟合值得到一个盐度差, 这个差值可认为主要是径流量变化引起的, 通过拟合方法可找出盐度差与径流量的关系式。若对此盐度差与径流量关系式直接进行拟合建模，会发现拟合效果不理想，分析其原因，发现此盐度差还是受到潮形的影响, 相同流量下大潮期间由于盐度大，对应的盐度差大; 小潮期间对应的盐度值小，盐度差小。为了得到更为真实的径流量对盐度的影响，需将盐度差进行潮形归一化, 得到归一化盐度差：（2），（3）（其中S为盐度实测值，DS为盐度差，SDS为归一化盐度差）。

通过分析归一化盐度差与径流量之间的关系，发现二者具有反比例函数关系。对归一化盐度差和径流量进行拟合，得到流量盐度差拟合函数: 

（4）（（4）中Q为径流量，SDS_M为归一化盐度差与径流量的拟合值，、、为拟合系数）。

将归一化盐度差与径流量的拟合值去归一化, 得到盐度差与径流量拟合关系：（5）。

所述步骤F00中，在步骤C00中预估拟合盐度与潮差关系时, 实测盐度资料中包含了径流的贡献。现在已经知道了径流量对盐度的贡献, 那么将实测盐度S减去径流量对盐度的贡献, 得到滤去径流量影响的盐度S1, ，将S1与潮差画在图上, 可看到滤去径流量影响的盐度与潮差同样具有指数关系，建立两者间的指数关系, 用最小二乘法拟合,得到修正的盐度潮差拟合关系式：（6）（（6）式中为盐度与潮差修正值, 、为拟合系数）。

与步骤D00预估拟合盐度与径流量关系相同方法, 将实测盐度S减去由（6）式确定的, 得到与径流量相关的盐度差, 并作归一化处理：（7）；

将归一化盐度差SDS1与径流量采用反比例函数, 最小二乘法拟合, 得到修正归一化盐度差SDS1_M与径流量Q的拟合关系式： （8）（公式中、、C0为拟合系数）。

将修正归一化盐度差与径流量的拟合关系式去归一化, 得到修正盐度差与径流量拟合关系:（9）。

所述步骤G00中，经过公式整理后，得到以下最终关系式：，，S(Y)=S1(TR)+ S1(TR)×S1(Q)= A*exp(B*TR)+ A*exp(B*TR)*[ C/(Q-D)-E] （10），（其中，S(Y)为预测盐度 ，TR为潮差，Q为流量，A、B、C、D、E均为拟合系数）。

根据上述步骤A00——步骤G00，可分别拟合出不同站点盐度日特征值的预报模型。盐度差与径流量的拟合关系（式（9））根据不同站点具体的不同关系表达式。

本基于潮汐和径流的河口咸潮预报方法基于河口咸潮上溯主要决定于径流量和潮汐的动力机制，分离短时间尺度潮差变化和长时间尺度径流量变化对咸潮上溯不同影响，建立口门站盐度日特征值（日平均、日最大、日最小值）与潮差和前15天平均径流量统计模型。盐度和潮差资料分布表明两者存在指数关系，采用最小二乘法先拟合盐度和潮差的预估关系。在找出盐度与潮差关系后, 观测盐度减去盐度与潮差拟合值, 得到盐度差值, 归一化后通过拟合方法可得出盐度差与径流量的预估关系。盐度与潮差的拟合关系加上盐度差与径流量的拟合关系, 可得出预估盐度统计模型。在预估拟合盐度与潮差关系时, 实测盐度资料中包含了径流的贡献，减除这部分贡献可得出修正盐度统计模型。经误差分析，修正盐度统计模型比预估盐度统计模型具有更高的精度。

根据上述基于潮汐和径流的河口咸潮预报方法，提供以下应用实例：

对本发明在珠江流域磨刀门水道进行了预测预测预报应用，并提供预测预报结果与实测结果对比分析（为适应日常对咸潮程度表征习惯，以下盐度改为日常用语含氯度进行表述）。

实测含氯度选用磨刀门水道的平岗泵站和广昌泵站的2009年10月1日-12月31日和2011年1月1日-2月28日实测逐时含氯度值，并计算含氯度日特征值：日均值、日最大值和日最小值。

预测预报使用的潮差选用灯笼山水文站潮汐预报表的潮位值，计算日最大潮差TR（为达到咸潮预测预报的预见期，此处没有选用测量实际潮位得到的实测潮差）。

预测预报使用的径流量为上游马口水文站前15日径流量平均值（马口站每日有实测径流量；由于马口径流量可由上游广西梧州水文站及上游水库放水控制，枯季流量波动不大，可提前1-7天预测马口径流量，并且使用前15天径流量均值，径流量均值较为稳定。）

利用上述基于潮汐和径流的河口咸潮预报方法应用实例验证平岗泵站和广昌泵站两站预测预报模型，得出日平均含氯度与径流量、潮差模型表达式，详见表1-1所示，

表1-1平岗、广昌日平均含氯度与径流量、潮差模型表达式

经过整理，最后得出最终表达式：

Y=(Y3+Y3×Y4)*1000+M（11），其中，（11）式和表1-1中，Y3（单位是g/L）为预估(初次)/修正潮差含氯度模型，TR为潮差；Y4（单位是g/L）为预估(初次)/修正径流量含氯度模型，Q为流量，M为常数。

式（11）中Y为咸度（含氯度）预测值，单位是mg/L；表1-1中Y0为初次模型含氯度预测值，Y1为修正模型含氯度预测值。

将上述咸潮预测预报结果与实测结果绘图对比，详见图6-11所示，可见咸潮模型能较好预测站点的含氯度变化过程。

本基于潮汐和径流的河口咸潮预报方法操作简单，实用性强：1、建立了基于盐度（含氯度）与径流、潮差关系的咸潮上溯预测模型，可提前1-7天预测，实现了中短期咸潮预测预报；2、预测预报所使用的潮差可用提前一年已经预测预报的潮汐表计算而得，可较易获得潮差资料。3、预测模型只有一个主要未知参数径流，把咸潮预测复杂问题简单化，具有简便实用的突出特点。

应当理解的是，本发明的应用不限于上述的举例，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，所有这些改进和变换都应属于本发明所附权利要求的保护范围。