水文预报课程设计.docx

水文预报课程设计澴河流域洪水预报方案编制姓名：韩夕林学号：2014*组号：15目录一、设计目的1二、课程设计任务2任务一、A 站降雨径流相关图制作与预报21流域自然地理特征22水文气象特征33预报站基本情况34水利工程概况35预报方案编制说明45.1流域特征提取45.2水文数据整理55.2.1雨量站权重55.2.2原始流量资料格式处理55.2.3历年最枯流量平均值65.2.4流域平均蓄水容量WM65.2.5流域月蒸发能力Em75.2.6消退系数K75.2.7退水曲线85.2.8蓄泄关系曲线105.3产流分析115.3.1产流模式115.3.2洪峰流量频率分析125.3.3产流量计算145.4水源划分185.4.1两场洪水分析Fc185.4.2Fc分析215.5汇流分析225.5.1单位线225.5.2线性水库245.5.3水位流量关系图246方案的检验、评定256.1产流方案检验评定256.2汇流方案检验评定277预报方案的应用和注意事项288方案作业预报应用28任务二、A 站与 H 站相应水位流量法预报291预报站基本情况292水利工程概况293预报方案编制说明303.1基本资料引用303.2采用方法及有关处理技术304方案的检验、评定305方案作业预报应用30三、课程设计的总结与心得32附录33一、设计目的 流域洪水预报方案的用途：洪水预报方案是现代实时洪水预报调度系统的核心部分，是提高预报精度和增长预见期的关键技术。对水资源可持续利用：流域水文模型是水资源评价、开发、利用和管理的理论基础。对水环境和生态系统保护：流域水文模型是构建面污染模型和生态评价模型的主要平台。本次课程设计的目的是通过一个具体的降雨径流预报方案的制作，使学生了解生产单位对预报任务的要求。 洪水预报方案是现代实时洪水预报调度系统的核心部分，是提高预报精度和增长预见期的关键技术。对水资源可持续利用而言，流域水文模型是水资源评价、开发、利用和管理的理论基础。对水环境和生态系统保护而言，流域水文模型是构建面污染模型和生态评价模型的主要平台。流域水文模型亦为分析研究气候变化和人类活动对洪水、水资源和水环境影响的有效工具。通过课程设计，要求掌握如下内容： 流域综合退水曲线、地下水退水曲线的制作以及次洪分割方法； 掌握次洪径流深及地面、地下流量分割方法； 掌握降雨径流相关图（API 模型）编制的完整过程； 流域经验单位线的推求方法； 洪水预报方案精度评定方法； 利用预报方案进行实时洪水预报方法； 相应水位流量预报方案的制作方法； 利用相应水位流量法进行河道洪水预报。 水文预报课程设计对应工程教育认证要求：毕业要求 4. 研究：能够基于科学原理并采用科学方法对复杂工程问题进行研究，包括方案设计、分析与解释数据、并通过信息综合得到合理有效的结论，具有归纳、整理、分析研究结果，撰写报告能力。 毕业要求5. 使用现代工具：能够针对复杂工程问题，掌握计算机、文献检索、科技方法等现代工程工具和信息技术工具，实现对复杂工程问题的预测与模拟，并能够理解其局限性。 支撑毕业要求指标点：能够独立对实验数据进行有效分析、整理及撰写实验报告、研究报告或论文。 研究报告结论合理有效，能解决某个具体复杂工程任务。 针对具体复杂工程问题，能选择、使用与开发恰当的工具。 能对复杂工程问题进行预测与模拟。 二、课程设计任务 任务一、A 站降雨径流相关图制作与预报1 流域自然地理特征澴河，古称澴水，发源于河南省罗山县与湖北省大悟县接壤的灵山。澴河流域界于东经1131211433，北纬30263129之间，主要在湖北省孝感市境内，干流全长150.8公里，流域面积3618平方公里。距河口57km处有A水文站，控制面积为2591km，占流域面积的71.6% ，干流长130.5公里，河床平均坡度0.97。A站控制流域水系图及A站站网一览表见下。图 1 A站控制流域水系图图 2 流域水系地图表 1 A站站网一览表A站站网一览表站名站别集水面积（km）坐标设立日期资料系列长（年）东经北纬A水文站2571135800311500195652B水库站27211413003150001972.555C雨量站2751141700313900195955D水文站4591140700313200195855E水文站4171140011313712196050F水文站43311350303136421931.155G雨量站47811356123127431959552 水文气象特征流域属于副热带大陆性季风气候，气候温和，四季分明。一月份是本流域最冷的月份，多年最低气温南部在零下68，北部在零下79。七月份一般是本流域最热的月份，多年最高气温南部在37左右，北部在37以上。流域内雨量充沛，多年平均雨量在1150mm左右。各年雨量的变化比较大，一般在均值上下200mm以内，最多的1969年 （13001600mm）比最少的1966年（650800mm）多一倍。全流域各地雨量多集中在59月，最大一日雨量：1968年7月13日为355.2mm。历史实测最大洪水位1996年7月17日，洪峰流量6060m/s，洪峰水位39.95m。3 预报站基本情况A水文站于1956年设立，属国家驻测重要站。经常性监测项目有降水、蒸发、水位、流量、水质、泥沙、地下水位、土壤墒情等。该站位于东经11358，北纬3115，引用的绝对基准面为黄海基面，采用基面为冻结基面。A站设防水位36.50m，警戒水位37.50m，保证水位39.95m。水位级别划分：38.5m以上为高水，34.00m38.50m为中水，34.00m以下为低水。预报标准：当流量大于1000m/s时，需要作业预报。A站系河流中下游控制站，集水面积为257km。4 水利工程概况该流域内有中型水库6座，总库容2.1789108m，控制面积225.8km，占流域面积的8.7%；有小型水库306座，控制面积424km，占流域面积的16.4%。水库总控制面积占流域面积的25.1%。A站流域内水库基本情况一览表见下。表 2 A站流域内水库基本情况一览表A站流域内水库基本情况一览表序号水库名称水库类型集雨面积（km）所在流域面积（km）总库容（万m）1R1中型70259158402R2中3中型50259177644R4中型32.2259121625R5中型37259132826R6中型19.6259111805 预报方案编制说明 5.1 流域特征提取本文使用ArcGIS软件对流域特征进行提取。ArcGIS是由ESRI出品的一个地理信息系统系列软件的总称。ArcGIS产品线为用户提供一个可伸缩的，全面的GIS平台。笔者用ArcGIS对流域特征进行提取的方法流程如下。数据结果见图 1 A站控制流域水系图和表 2 A站流域内水库基本情况一览表。图 3 流域特征提取方法流程图5.2 水文数据整理5.2.1 雨量站权重本文降雨资料为7个站点的逐日降雨数据和部分降雨摘录表。逐日降雨数据提供了各个雨量站点每天总的降雨量，降雨部分摘录表提供了时段降雨量，在汛期雨量较大时时段比较短，在枯水期雨量较少时时段比较长。预报所需的降雨应为A站控制流域内的总的降雨数据，本文采用泰森多边形法计算流域平均雨量。笔者通过ArcGIS制作泰森多边形并计算各站点集水面积方法流程如下。各站点集水面积数据见表 3 雨量站面积权重表。图 4 集水面积计算方法流程图再通过泰森多边形面积对雨量站雨量数据进行面积加权平均，得到各雨量站雨量权重情况见下表。表 3 雨量站面积权重表雨量站面积权重表站名ABCDEFG集水面积（km）257272275459417433478面积权重0.100.110.110.180.160.170.185.2.2 原始流量资料格式处理流量资料从A站水文要素摘录表中可以得到。需要注意的是，年、月、日、时时间数据列需对格式进行转化，使其在一列单元格中体现。方法流程如下。图 5 日期格式转换方法过程5.2.3 历年最枯流量平均值统计历年最枯流量数据，统计表格如下。经分析计算历年最枯流量的平均值为2.22m/s。表 4 历年最枯流量数据统计表历年最枯流量数据年份1969197019711972197319741975流量（m/s）4.781.233.71.724.890.038年份1976197719781979198019811982流量（m/s）1.32.50.300.0271.740.803年份1983198419851986198719881989流量（m/s）0.412.952.61.4510.20.72.245.2.4 流域平均蓄水容量WMWM是流域平均蓄水容量。本文从所有的历史实测资料中选取久旱无雨之后突然一场暴雨引起的洪水来计算WM。通常逐日降雨量小于10mm的天数超过10天即可认为是久旱。本文采用的是1976年7月14日为起涨日的洪水，命名为“760714洪水”。13日和14日各站降雨量如下。表 5 760714洪水降水数据760714洪水降雨数据（mm）时间ABCDEFG流域平均13日153.842.582.6104.185.245.49886.3014日103.415.160.2102.316.610.79157.60从A站水文要素摘录表可以知道洪水起涨时间大约是1976年7月14日凌晨2点。则5.2.5 流域月蒸发能力Em蒸发资料从A站逐日水面蒸发数据可以获取。原数据是历史实测蒸发数据，本文需要推求流域蒸散发能力。寻找每年各月蒸发量的最大值近似作为本月的蒸发能力，再取21年的平均值，即可得到该流域各月份的蒸发能力。提取及计算MATLAB程序如下：clear;clcdata=xlsread(A站逐日水面蒸发量表,A2:D7423);year=(data(:,1)-1969).*12;way=year+data(:,2);data1=data(1:321,:);%1969年无12月数据，分开算data(1:321,:)=;way1=way(1:321);way(1:321)=;num=zeros(252,1);%21年共252个月for i=1:11 num(i)=max(data1(find(way1=i),4);endfor i=13:252 num(i)=max(data(find(way=i),4);end %提取各月的最大蒸发量作为该月蒸发能力k=reshape(num,12,21);ave=mean(k,2);ave(12)=mean(k(12,2:21),2);%21年的平均蒸发能力蒸散发能力数据结果见下表。表 6 多年平均蒸散发能力Em多年平均蒸散发能力Em月份JanFebMarAprMayJunJulAugSepOctNovDec蒸发能力（mm）3.173.775.247.498.478.849.178.98.044.953.762.675.2.6 消退系数KK是消退系数。在考虑到各月蒸发条件和蒸发能力不同时，消退系数K按月变化，其计算公式为 （1）其中 为，是土壤最大缺水量，其值近似为土壤蓄水容量,从而计算得到消退系数各月数值如下。表 7 K值表K值表月份JanFebMarAprMayJunJulAugSepOctNovDecK值0.970.960.950.930.920.910.910.910.920.950.960.975.2.7 退水曲线本文从1969年至1989年21年的时间里选取了16场退水数据，时间间隔为24小时。选取退水数据如下。表 8 16场退水流量数据16场退水流量数据（m/s）时012345678910134.0333.9833.9233.8833.8433.833.7733.7733.7433.7233.672867385155814832.124.319.916.815.30361121910967.750.838.529.422.117.8140443.931.62419.215.313.912.610.810.29.78.61533.630.124.92320.31714.912.210.810.10635.3134.3133.9533.7733.6633.533.4433.3933.3533.3233.3744812861.334.322.715.611.38.997.396.470848118710563.844.235.930.326.12219.116.4972595445403428262523010308122806555474237353233.291130414569.245.832.82623.922.519.81613.11227.520.211.97.115.224.443.22.922.391.381.221355.651.235.726.922.417.314.613.412.111.201415310789.170.952.54236.433.1933.1333.0601523010562.844.933.828.125.723.320.218.820.21613180.845.130.921.31612.3118.467.45.82由于数据系列不够长，因此对上述数据进行插值处理，使时间间隔为6h，插值程序如下。data=xlsread(16场洪水退水选取数据,B1:M17);data=data;t=data(:,1);q=data(:,2:17);tx=t(1):0.25:t(end);qx=zeros(45,16);for i=1:16 qx(:,i)=interp1(t,q(:,i),tx);endxlswrite(退水插值数据,qx);再根据这已知的16场插值洪水数据在Excel中绘制退水曲线，见图 6 16场退水图。通过分析退水曲线，可定退水规律符合指数函数的形式。图 6 16场退水图绘制相邻时段流量关系图见图 7 相邻时段流量关系图。从图可以看出，下部曲线大致重合，且接近直线，反映了地下退水特性，因此可以用外包线对其进行概括。上部分散的曲线簇用平均曲线进行代替，反映平均退水情况。图 7 相邻时段流量关系图其中，黑线是外包线，其斜率为1.0570，则。因为，则令，计算、，根据计算结果绘制流域的退水曲线见下图。图 8 流域标准退水曲线5.2.8 蓄泄关系曲线根据流域退水曲线建立退水段流量与相应的退水径流深之间的相关关系，点绘成如下图所示的关系曲线。图 9 蓄泄关系曲线此为蓄泄关系曲线，可以用来计算次洪径流深。笔者选取1983数据年内两场洪水进行计算示例。本文选取“830721洪水”和“830915洪水”进行计算。对两场洪水流量数据插值后点绘流量过程线见图 10 两场洪水流量过程线。根据公式（2）分别计算两场洪水径流深。其中为本次洪水退水流量所相应的退水径流深；为前次洪水退水流量所相应的退水径流深，都由图 9 蓄泄关系曲线查得。计算得到次洪径流深分别是图 10 两场洪水流量过程线5.3 产流分析 5.3.1 产流模式根据每年各站逐日降雨量和面积权重得到从1969年到1989年共21年在该流域的平均降雨量计算表如下。表 9 多年平均降雨量计算表多年平均降雨量计算表雨量站点累计雨量（mm）平均雨量（mm）权重流域平均雨量（mm）A24059.00 1145.67 0.10 1056.55B24371.60 1160.55 0.11 C23261.20 1107.68 0.11 D23922.00 1139.14 0.18 E21016.50 1000.79 0.16 F19703.30 938.25 0.17 G20913.30 995.87 0.18 可以看到，多年平均降雨量为1056.55mm，超过了1000mm。从所给流量资料里选取两场洪水，绘制洪水过程线见图 11 澴河流域洪水过程代表。可以看到该流域的洪水过程线不具有对称性。再结合该流所处地带，澴水流域属于我国南方地区，土壤颗粒较大、结构疏松、植被较高、地下水也比较丰富，因此属于蓄满产流区。图 11 澴河流域洪水过程代表该流域属于湿润地区,雨量充沛,植被良好,具有蓄满产流的特点,但水利工程和人类活动影响较大,产流方案采用以前期影响雨量 Pa 为参数的降雨径流经验相关法建立。 5.3.2 洪峰流量频率分析本文采用年最大值法对所给21年的流量数据进行样本选择，即寻找各年最大流量作为洪峰流量设计的样本。利用MATLAB对原始资料进行提取，得到各年最大流量数据如下。表 10 各年最大流量各年最大流量（m/s）年份1969197019711972197319741975流量3580276094451914309282630年份1976197719781979198019811982流量2050120020025934702321090年份1983198419851986198719881989流量2940110012908443750709879MATLAB提取程序如下：clear;clcnum=zeros(21,1);for temp=1969:1989t1=xlsread(A站洪水水文要素摘录表,num2str(temp),E6:E55);t2=xlsread(A站洪水水文要素摘录表,num2str(temp),J6:J55);t3=xlsread(A站洪水水文要素摘录表,num2str(temp),O6:O55);t4=xlsread(A站洪水水文要素摘录表,num2str(temp),E61:E110);t5=xlsread(A站洪水水文要素摘录表,num2str(temp),J61:J110);t6=xlsread(A站洪水水文要素摘录表,num2str(temp),O61:O110);t7=xlsread(A站洪水水文要素摘录表,num2str(temp),E116:E165);t8=xlsread(A站洪水水文要素摘录表,num2str(temp),J116:J165);t9=xlsread(A站洪水水文要素摘录表,num2str(temp),O116:O165);t=t1;t2;t3;t4;t5;t6;t7;num(temp-1968)=max(t);endv=reshape(num,7,3);v=v;xlswrite(年最大洪峰流量,v);本文不考虑特大洪水的影响，将流量数据降序排列，并按照统一处理法计算频率，计算公式如下（3）其中，N是样本容量，在这里是21，n是序列数。计算频率如下表。表 11 计算频率表计算频率表序号1234567流量3750358034702940276026302050频率0.050.090.140.180.230.270.32序号891011121314流量14301290120011001090944928频率0.360.410.450.500.550.590.64序号15161718192021流量879844709519259232200频率0.680.730.770.820.860.910.95本文采用武汉大学水文频率分布曲线适线软件（徐磊改版）-连续系列进行适线。适线结果图如下。获取两年一遇洪峰设计值为1318.02m/s。图 12 洪峰流量频率曲线笔者的数据年是1983年，从1983年的A站水文摘录表中选择洪峰为1370m/s的洪水进行分析，其起涨日为1983年9月15日，故命名为“830915洪水”。对洪水流量进行插值，间隔为6个小时，插值主要程序如下：clc;clear;a=textread(插值数据.txt);b=interp1(a(:,1),a(:,2),a(1,1):0.25:a(end,1);绘制流量过程线如下：图 13 830915洪水流量过程线5.3.3 产流量计算前期雨量指数Pa的计算公式为（4）其中， 是第t+1日的前期雨量指数，k是月消退系数， 是第t日的降雨量， 是第t日的前期雨量指数。本文选取1983年9月15日前的75天作为计算9月15日前期雨量指数的预热期，即7月2日至9月14日。先取各站逐日降雨量计算流域平均降雨量，假设7月2日的前期雨量指数，再按照公式（4）逐个计算各日的前期雨量指数。Pa值计算趋势如图 14 Pa值计算趋势图。可以看到，取不同的初始Pa值对计算9月15日的前期雨量指数并没有影响，都是32.1mm，这就是预热期的作用。图 14 Pa值计算趋势图蓄满产流模型：（5）（6）（7）（8）其中，b是常数，反映流域包气带蓄水容量分布的不均匀性。W是初始土壤含水量，在这里数值上可以理解成前期雨量。WM是流域平均蓄水容量。当b=0.4时，假设初始土壤含水量W=20mm，那么如果降雨量PE=20mm，由于a+PE=40.63FC)./PE(find(PEFC)+sum(r(PE10 q1=q2; q1(1)=;q1(end)=; Qtemp=Q-q1.*p(1)./10;0;0+0;0;q1.*p(3)./10; q2=Qtemp.*10./p(2); t=sum(abs(q2-0;q1;0); i=i+1; if i1000 break endend得到q1，再通过q1得到的出流过程与实际出流过程进行比较调整，调整对比结果图如下。绘制830915洪水单位线如下图。图 24 830915洪水单位线5.5.2 线性水库地下径流汇流过程采用线性水库计算，根据前文得到（11）其中，。时段为6h。那么 （12）则得到蓄泄常数 5.5.3 水位流量关系图经对比分析，A站洪峰水位流量关系逐年变化，即右移，故需要及时修改。提取历年最高水位、最大洪峰见表 18 历年最高水位、最大洪峰数据统计表。率定洪峰水位关系线见图 25 水位流量关系图。表 18 历年最高水位、最大洪峰数据统计表历年最高水位、最大洪峰数据统计表年份最高水位（m）最大流量（m/s）196939.23580197038.142760197136.23944197235.54519197337.091430197436.099281975992630197695205019779912001978992001979992591980833470198134.3232198235.75725198338.672940198436.391100198536.71290198636.73844198737.912270198835.66709198935.96879经过对数据分析，除掉异常点，点绘水位流量关系如下：图 25 水位流量关系图6 方案的检验、评定6.1 产流方案检验评定方案的评定、检验标准按水文情报预报规范（GB/T 22482-2008）中许可误差标准进行。该站为河道站，预报项目主要为洪峰流量及峰现时间。表 19 降雨径流相关图预报检验方案检验评定表序号降雨量P(mm)前期影响雨量Pa(mm)实测径流深R(mm)计算径流深R(mm)误差R-R(mm)许可误差(mm)评定1113.597.511011110.01合格282.551.447.133.9-13.2-0.28合格392.969.451.462.310.90.21不合格4147.121.570.668.6-2-0.03合格5113.470.589.783.9-5.8-0.06合格652.782.23534.9-0.10.00合格750.675.229.325.8-3.5-0.12合格886.835.330.422.1-8.3-0.27合格983.787.975.371.6-3.7-0.05合格10144.837.2788240.05合格11117.152.964.8705.20.08合格1288.349.43537.72.70.08合格1362.577.43339.96.90.21不合格1416161.9110122.912.90.12合格1574.899.768.674.55.90.09合格16208.461.9173.9170.3-3.6-0.02合格17105.244.854.450-4.4-0.08合格1892.478.467.270.83.60.05合格19121.864.865.386.621.30.33不合格20107.67477.781.63.90.05合格21130.447.96978.39.30.13合格2290.346.532.736.84.10.13合格2388.175.767.763.8-3.9-0.06合格2490.744.23534.9-0.10.00合格25104.840.44545.20.20.00合格2682.48068.662.4-6.2-0.09合格27110.350.738.76122.30.58不合格28182.146.1126.2128.220.02合格29128.8100120128.88.80.07合格30151.494.7150146.1-3.9-0.03合格31124.526.543.9517.10.16合格32122.13747.859.111.30.24不合格33201.627.9121.4129.58.10.07合格3474.865.837.640.630.08合格35101.685.584.787.12.40.03合格36176.759.1117.3135.818.50.16合格37106.159.376.165.4-10.7-0.14合格3885.995.364.881.216.40.25不合格3910379.668.982.613.70.20合格40115.4838598.413.40.16合格41195.652.6137.9148.210.30.07合格合格率85.4%6.2 汇流方案检验评定选760713洪水进行检验评定，地面净雨数据如下：表 20 760713洪水地面净雨地面净雨（mm）3.4721.7318.722.7616.414.94通过830915洪水推求的单位线计算出流过程结果见图 26 汇流预报出流。图 26 汇流预报出流7 预报方案的应用和注意事项在推求地面径流过程时，单位线为多条，当各时段地面净雨强度分布不均匀时，需要考查暴雨中心在流域的位置，即哪一个站点的降雨强度最大，位置分为暴雨中心在上游、暴雨中心在中游、暴雨中心在下游。分别采用对应相应暴雨中心位置的单位线进行推流。人类活动影响，本方案只考虑中型水库，对水库人为的蓄量应该进行还原。近几年来河道采砂严重，河床有所下降，如何预报洪峰水位。8 方案作业预报应用进入作业预报工作模式，设定相应降雨开始时间、初始Pa、预报根据时间等参数即可进行相应的作业预报，分析计算出该场洪水的洪峰流量。洪峰水位预报：通过建立的洪峰水位、洪峰流量相关图，用预报的洪峰流量查线推求相应的洪峰水位。表 21 A站控制流域雨量权重A站流域泰森多边形雨量权重雨量站ABCDEFG权重0.100.110.110.180.160.170.18预报检验表见附录页表格。任务二、A 站与 H 站相应水位流量法预报1 预报站基本情况 H 站位于东经 11353，北纬 3056，是河下游控制站，控制面积5157km2。与上游 A 站相距 56.8km，区间面积2566km2，占该站以上集水面积的49.8%。该站于1956 年 6 月设立为水文站，1957 年改为水位站至今。引用的绝对基面为黄海基面，本站所用为冻结基面，冻结基面以上水位减 2.125m 等于黄海基面以上水位。该站设防水位 26.50m，警戒水位 27.50m，保证水位 31.00m。历年实测最高水位31.67m。 H站控制流域水系图如下，可以看到，水系地图南端比较窄，支流很少，适合用河道洪水演算进行推流。图 27 H站控制流域水系2 水利工程概况 A 站至 H 站区间流域内有大、中型水库 3 座，总集水面积 135.9km2，总库容1.5123108m3，设计最大总泄洪流量 1392m3/s。沿岸建排灌闸 17 座，最大排灌流量总计 211m3/s。建泵站 5 处，总装机容量 1960KW。 表 22 A 站至 H 站区间流域内水库基本情况一览表A 站至 H 站区间流域内水库基本情况一览表序号水库名称水库类型集雨面积（km2）所在流域面积（km2）工程影响情况(%)总库容