一、匯流原理

二、操作步驟（流向—產(chǎn)流比例—流程文件—水文站點位置—全局參數(shù)文件—匯流模型運行）

1.產(chǎn)流比例文件（fraction.txt）

文件內(nèi)容：里面是邊界網(wǎng)格被包含到實際的流域邊界內(nèi)的部分面積占該網(wǎng)格面積的比例，可以用ArcGIS制作。

（1）柵格轉(zhuǎn)點+點轉(zhuǎn)柵格：在ArcGIS中，導(dǎo)入研究區(qū)柵格數(shù)據(jù)"hhy_grid"，進行"柵格轉(zhuǎn)點"+"點轉(zhuǎn)柵格(pointed字段)"處理，得到分網(wǎng)格的柵格數(shù)據(jù)。

（2）柵格轉(zhuǎn)面：取消勾選"簡化面"，輸出"HHY_rout_1.shp"

（3）柵格轉(zhuǎn)面：導(dǎo)入"dem_one"數(shù)據(jù)，進行"柵格轉(zhuǎn)面"處理，輸出"HHY_rout_2.shp"

以上兩個步驟輸出數(shù)據(jù)的屬性表里應(yīng)有"面積"字段，可我輸出的沒有，重新運行，不修改輸出路徑，按GIS默認路徑輸出就沒有這個問題，這是什么bug???

（4）相交：將"HHY_rout_1.shp"和"HHY_rout_2.shp"進行相交處理

（5）導(dǎo)出屬性表：所需數(shù)據(jù)（紅框內(nèi)，數(shù)據(jù)表列數(shù)與教程中稍有出入，不影響結(jié)果）

（6）在RStudio中，打開代碼"1_Fraction"，修改輸入（上步屬性表）輸出路徑后運行

（7）柵格轉(zhuǎn)ASCII：在ArcGIS中，將研究區(qū)柵格數(shù)據(jù)"hhy_grid"進行柵格轉(zhuǎn)ASCII處理

（8）在RStudio中，打開代碼"2_precision"，修改輸入輸出路徑、行列數(shù)后運行

Error1：

解決：創(chuàng)建ckgssj副本，命名為ckgssj1和ckgssj2，將ckgssj1表頭刪除，并將R中輸入文件改為ckgssj1

Error2：同上

解決：修改行列數(shù)（寫反了）

產(chǎn)流比例文件

2.流向文件（flowdir.txt）

文件內(nèi)容：里面存放的是每個網(wǎng)格內(nèi)水流的方向，用數(shù)宇1-8表示八個方位，這里一定不能用ArcGIS直接做流向，因為會丟失很多細節(jié)（重采樣丟細節(jié)），而且容易造成封閉流向，所采用的方法是利用每個網(wǎng)格的最大匯流累積量代替高程，跑出來也基本上是對的。

（1）柵格轉(zhuǎn)點+點轉(zhuǎn)柵格：在ArcGIS中，導(dǎo)入研究區(qū)柵格數(shù)據(jù)"hhy_grid"，進行"柵格轉(zhuǎn)點"+"點轉(zhuǎn)柵格(pointed字段)"處理，得到分網(wǎng)格的柵格數(shù)據(jù)

?（2）重采樣：分辨率大小和dem的相同

（3）按掩膜提取

（4）柵格轉(zhuǎn)點：對第二步"重采樣"和第三步"按掩膜提取輸出"的數(shù)據(jù)進行柵格轉(zhuǎn)點處理（耗時較久，20min左右）

（5）添加連接：將第四步輸出的數(shù)據(jù)用"pointid"字段進行"連接"

（6）導(dǎo)出屬性表（10min）

（7）在RStudio中，打開代碼"3_flowdir"，修改輸入（上步屬性表）輸出路徑后運行

（8）打開代碼"4_precision"，修改輸入輸出路徑、行列數(shù)后運行

（9）打開代碼"5_flowdir2"，修改輸入輸出路徑后運行

（10）將"F_div"中的表頭拷貝至"hhy_F_dir1"

最終流向文件展示：

3.流程文件的制作（xmask）

（1）柵格轉(zhuǎn)點+點轉(zhuǎn)柵格+柵格轉(zhuǎn)面：在ArcGIS中，導(dǎo)入研究區(qū)柵格數(shù)據(jù)"hhy_grid"，進行"柵格轉(zhuǎn)點"+"點轉(zhuǎn)柵格(pointed字段)"+"柵格轉(zhuǎn)面"處理，分別輸出xmask_1/2/3。

（2）在RStudio中，打開代碼"7_Xmask"，先完成①②，后根據(jù)此修改代碼輸入輸出路徑等

①測量x_mask3中網(wǎng)格橫、豎、對角長度（流向）

②拷貝flowdir輸出結(jié)果"hhy_F_dir1"至x_mask文件夾，并刪除其表頭

③在"7_Xmask"中修改代碼輸入輸出路徑等，運行代碼輸出結(jié)果，給結(jié)果文件添加表頭

4.確定水文站在模擬柵格中的位置

計算公式如下：根據(jù)已知水文點經(jīng)緯度、xllcorner、yllcorner、cellsize，求RN(行數(shù))和CN(列數(shù))

其中xllcorner、yllcorner、cellsize（網(wǎng)格分辨率）值如下：

計算得出：

據(jù)此，制備水文站位置文件"F_F_stnloc.txt"

5.全局參數(shù)文件的制作

（1）將官網(wǎng)例子"stehekin\params\rout"中的全局參數(shù)文件"rout_input.STEHE"+"UH.all"和第1-4步制備的流向、產(chǎn)流比例、流程、水文站位置文件拷貝至"F:\VIC_modle\result\F_F\rout\param"中，將全局參數(shù)文件格式改為"txt"格式；將官網(wǎng)例子"route_1.0"文件夾拷貝至"F:\VIC_modle\result\F_F\rout"中

（2）修改全局參數(shù)文件輸入輸出路徑

（3）將"rout_1.0"中的"Makefile"文件中第25行"FC=g77"修改為"FC=gfortran"

6.匯流模型的運行

（1）打開"Cygwin64 Terminal"軟件，進行編譯、運行，輸出結(jié)果為"rout.exe"，將其拷貝至Cygwin安裝軟件目錄bin里

（2）繼續(xù)在"Cygwin64 Terminal"軟件中，進行編譯（全局參數(shù)）

Error1：

解決方法：參考https://blog.csdn.net/weixin_54697010/article/details/129163963?

修改：

Error2：

可能原因：產(chǎn)流（ArcGIS導(dǎo)出）、匯流（R）坐標(biāo)不同

解決步驟：

①重采樣：在ArcGIS中打開網(wǎng)格柵格數(shù)據(jù)，進行重采樣

②柵格轉(zhuǎn)ASCII

③修改匯流文件的坐標(biāo)：用輸出的ASCII文件的表頭，替換以下三個文件的表頭

④柵格轉(zhuǎn)點：給轉(zhuǎn)點數(shù)據(jù)的屬性表添加兩個字段"x""y"，并通過計算幾何賦值，最后導(dǎo)出屬性表

⑤重新制備forcing文件：打開制備氣象驅(qū)動數(shù)據(jù)的最后一個代碼"Conbination_PMNW"，修改經(jīng)緯度參考坐標(biāo)文件的路徑，刪除原有輸出文件的數(shù)據(jù)，重新運行代碼，將輸出結(jié)果拷貝至"F:\VIC_modle\result\F_F\vic\forcing"

⑥修改土壤參數(shù)文件的坐標(biāo)：打開"8modifySumSoilParam"代碼，修改輸入輸出路徑后運行

⑦重新制備產(chǎn)流文件：打開"cmd"命令提示符，運行vic模型

運行成功（比較耗時-1h）

⑧重新制備匯流文件：打開"cmd"命令提示符，運行rout模型

依然報錯！??！

不要慌！不要急！仔細檢查！繼續(xù)修改！問題就是一個個慢慢解決的??！

Error3：水文站點位置問題

解決方法：經(jīng)反復(fù)檢查，在第四步"4.確定水文站在模擬柵格中的位置"中計算得出的數(shù)據(jù)沒有錯誤，可能是數(shù)據(jù)稍有偏差，因此我嘗試小幅度修改行列值，在修改為"52，40"時成功了

重新運行rout模型，水文站點位置問題成功解決?。?/strong>

生成的過程線文件（TNH? .uh_s）1686行，此前只有1行

為方便后續(xù)調(diào)參，可在"F_F_stnloc"中，將運行出的正確的過程線文件地址輸入（替換"NONE"）

Error4（與Error2一樣）：

自動識別出的經(jīng)緯度坐標(biāo)與徑流輸出經(jīng)緯度不匹配（參數(shù)文件問題）制作的文件存在偏移

我通過對比報錯經(jīng)緯度和fluxes_原數(shù)據(jù)，發(fā)現(xiàn)：

①經(jīng)度：0.0049——0.0050

②緯度：0.0016 —— 0.0017

解決方法：修改參數(shù)（流向、產(chǎn)流比例、流程文件表頭的xllcorner和yllcorner）

（要確定此前用的同一個網(wǎng)格數(shù)據(jù)，否則調(diào)參無效，在Error2中已重做修改）

第一次修改：

失敗原因：緯度小了0.0001（第一次改動幅度太大了，把小數(shù)點后面刪了很多，吸取教訓(xùn)）

第二次修改：

失敗原因：經(jīng)度大了0.0001

第三次修改：

失敗原因：經(jīng)度小了0.0001

（本來以為這次成功了，前面非常絲滑，結(jié)果從第805個網(wǎng)格開始報錯/(ㄒoㄒ)/~~）

第四次修改：

失敗原因：和第三次一樣，經(jīng)度依然小了0.0001

第五次修改：

這次離成功很近了?。。≈皇Ｏ?5個not found ?。?！

第六次修改......第n次修改......放棄......

我決定從根本上解決問題，因為每次都是那25行經(jīng)度為_100.0417的問題。經(jīng)檢查，是經(jīng)緯度原數(shù)據(jù)中"100.041650164299995"的問題，將其全部更改為"100.041640164299995"，即保證輸出為_100.0416，然后重新制備產(chǎn)流文件（Error2中的第⑤⑥⑦步）

我想這次應(yīng)該是沒問題了，再次利用第五次修改中的坐標(biāo)數(shù)據(jù)，進行匯流

沒有報錯?。。R流成功?。?！

7月中旬至今，經(jīng)過間斷性學(xué)習(xí)，歷經(jīng)重重報錯，終于跑通了整個VIC模型，輕舟已過萬重山，當(dāng)然，后續(xù)還需要繼續(xù)學(xué)習(xí)參數(shù)率定，以及需要從頭來過再跑一遍VIC模型......因為

記得剛開始學(xué)習(xí)時，我有一個遺留問題，我搜索"VIC模型多站點"等關(guān)鍵詞并沒有得到答案，于是我天真的以為，VIC模型想輸出多個站點數(shù)據(jù)，每個站點都要重新跑一遍，于是我只跑了第一個站點，直到完成匯流，我才在某篇文章中發(fā)現(xiàn)，匯流時可以輸入多個站點數(shù)據(jù)......這意味著，這一個半月只是一次練習(xí)，在接下來的一周或者一個月，我得全部從頭開始再跑一遍數(shù)據(jù)輸出所有我需要的站點數(shù)據(jù)......哭了/(ㄒoㄒ)/~~......

參考文章：https://blog.csdn.net/qq_32611933/article/details/50847191?

我為什么沒有早點看到??。。?？？？

（今天才知道原來bilibili文章圖片不能超過100張，刪減了一下才得以提交）

標(biāo)簽：