雨量預(yù)報(bào)分析的評(píng)價(jià)模型 數(shù)學(xué)建模

1、雨量預(yù)報(bào)分析的評(píng)價(jià)模型一、摘要我們將FORECAST文件夾中的數(shù)據(jù)按日期先后順序?qū)隡atlab，建立53×47×164的三維矩陣rain1和rain2；把MEASURING文件夾中的數(shù)據(jù)以同樣方法導(dǎo)入91×7×41的三維矩陣temp中，然后建立循環(huán)將temp矩陣中每一層的后4列提取，另存入一個(gè)91×164的rain3矩陣；在命令窗中直接導(dǎo)入預(yù)測點(diǎn)的經(jīng)度和緯度存入矩陣lon和lat中，導(dǎo)入實(shí)測點(diǎn)的經(jīng)度和緯度存入矩陣lon1和lat1中，并對其作圖，得到實(shí)測點(diǎn)和預(yù)測點(diǎn)的經(jīng)緯度圖。整理得到91個(gè)觀測點(diǎn)41天的預(yù)測值和測量值對應(yīng)的兩個(gè)91×

2、164矩陣，根據(jù)氣象部門將降雨的等級(jí)分為6個(gè)等級(jí)的分法，把矩陣中相應(yīng)的降雨量值轉(zhuǎn)化為其所對應(yīng)等級(jí)值，其中，預(yù)測中的零全部記為0，得到兩個(gè)預(yù)報(bào)等級(jí)矩陣。針對問題（1），利用插值基點(diǎn)為散亂節(jié)點(diǎn)的插值函數(shù) 1在Matlab中進(jìn)行三次樣條插值處理，將91個(gè)觀測站點(diǎn)41天164個(gè)時(shí)段的雨量情況進(jìn)行預(yù)測。利用殘差平方和以及平均誤差來作為評(píng)價(jià)的標(biāo)準(zhǔn)。殘差平方和與平均誤差值較小的一種預(yù)測方法作為較好的預(yù)報(bào)方法。殘差平方和以及平均誤差數(shù)值越小，表明預(yù)報(bào)越準(zhǔn)確度越高。預(yù)測方法一的殘差平方和為174290.00，平均誤差為0.4553。預(yù)測方法二的殘差平方和為195580.00，平均誤差為0.4753。雨量預(yù)報(bào)方法

3、一的準(zhǔn)確性更高一些。針對問題（2），兩個(gè)預(yù)報(bào)等級(jí)矩陣，繼續(xù)利用殘差平方和以及平均誤差來作為評(píng)價(jià)的標(biāo)準(zhǔn)。殘差平方和以及平均誤差數(shù)值越小，表明預(yù)報(bào)越準(zhǔn)確度越高，相應(yīng)公眾感受就越好。預(yù)測方法一的殘差平方和為2774，平均誤差為0.1730。預(yù)測方法二的殘差平方和為2806，平均誤差為0.1745。雨量預(yù)報(bào)方法一的準(zhǔn)確性更高一些。由于殘差平方和與平均誤差難以反映真實(shí)匯報(bào)的準(zhǔn)確度，我們將模型改進(jìn)優(yōu)化。把矩陣中相應(yīng)的降雨量值轉(zhuǎn)化為其所對應(yīng)等級(jí)值，得到兩個(gè)預(yù)報(bào)等級(jí)矩陣，將兩個(gè)預(yù)報(bào)等級(jí)矩陣與實(shí)測等級(jí)矩陣做差值運(yùn)算，得到兩個(gè)等級(jí)差矩陣，對等級(jí)差作絕對值處理，進(jìn)行等級(jí)差統(tǒng)計(jì)。我們利用預(yù)測準(zhǔn)確度檢驗(yàn)法對兩種預(yù)報(bào)進(jìn)行

4、評(píng)價(jià)。預(yù)測準(zhǔn)確度（）等于預(yù)報(bào)正確次數(shù)（）（即運(yùn)算之差為0的情況）和預(yù)測次數(shù)（）之比，即。準(zhǔn)確度越高，表明預(yù)報(bào)準(zhǔn)確度越高，相應(yīng)公眾感受就越好。預(yù)報(bào)1的預(yù)報(bào)準(zhǔn)確度為83.26%高于預(yù)報(bào)2的準(zhǔn)確度83.11%，公眾更易接受第一種預(yù)報(bào)方法。關(guān)鍵字：散亂節(jié)點(diǎn)插值 殘差平方和 平均誤差 預(yù)報(bào)等級(jí)矩陣 預(yù)測準(zhǔn)確度二、問題重述雨量預(yù)報(bào)對農(nóng)業(yè)生產(chǎn)和城市工作和生活有重要作用，但準(zhǔn)確、及時(shí)地對雨量作出預(yù)報(bào)是一個(gè)十分困難的問題，廣受世界各國關(guān)注。我國某地氣象臺(tái)和氣象研究所正在研究6小時(shí)雨量預(yù)報(bào)方法，即每天晚上20點(diǎn)預(yù)報(bào)從21點(diǎn)開始的4個(gè)時(shí)段（21點(diǎn)至次日3點(diǎn)，次日3點(diǎn)至9點(diǎn)，9點(diǎn)至15點(diǎn)，15點(diǎn)至21點(diǎn)）在某些位置的

5、雨量，這些位置位于東經(jīng)120度、北緯32度附近的53×47的等距網(wǎng)格點(diǎn)上。同時(shí)設(shè)立91個(gè)觀測站點(diǎn)實(shí)測這些時(shí)段的實(shí)際雨量，由于各種條件的限制，站點(diǎn)的設(shè)置是不均勻的。氣象部門希望建立一種科學(xué)評(píng)價(jià)預(yù)報(bào)方法好壞的數(shù)學(xué)模型與方法。氣象部門提供了41天的用兩種不同方法的預(yù)報(bào)數(shù)據(jù)和相應(yīng)的實(shí)測數(shù)據(jù)。預(yù)報(bào)數(shù)據(jù)在文件夾FORECAST中，實(shí)測數(shù)據(jù)在文件夾MEASURING中，其中的文件都可以用Windows系統(tǒng)的“寫字板”程序打開閱讀。FORECAST中的文件lon.dat和lat.dat分別包含網(wǎng)格點(diǎn)的經(jīng)緯度，其余文件名為<f日期i>_dis1和<f日期i>_dis2，例如f6

6、181_dis1中包含2002年6月18日晚上20點(diǎn)采用第一種方法預(yù)報(bào)的第一時(shí)段數(shù)據(jù)（其2491個(gè)數(shù)據(jù)為該時(shí)段各網(wǎng)格點(diǎn)的雨量），而f6183_dis2中包含2002年6月18日晚上20點(diǎn)采用第二種方法預(yù)報(bào)的第三時(shí)段數(shù)據(jù)。MEASURING中包含了41個(gè)名為<日期>.SIX的文件，如020618.SIX表示2002年6月18日晚上21點(diǎn)開始的連續(xù)4個(gè)時(shí)段各站點(diǎn)的實(shí)測數(shù)據(jù)（雨量），這些文件的數(shù)據(jù)格式是：站號(hào) 緯度 經(jīng)度 第1段 第2段 第3段 第4段 58138 32.9833 118.5167 0.0000 0.2000 10.1000 3.100058139 33.3000 118

7、.8500 0.0000 0.0000 4.6000 7.400058141 33.6667 119.2667 0.0000 0.0000 1.1000 1.400058143 33.8000 119.8000 0.0000 0.0000 0.0000 1.800058146 33.4833 119.8167 0.0000 0.0000 1.5000 1.9000雨量用毫米做單位，小于0.1毫米視為無雨。(1) 請建立數(shù)學(xué)模型來評(píng)價(jià)兩種6小時(shí)雨量預(yù)報(bào)方法的準(zhǔn)確性；(2) 氣象部門將6小時(shí)降雨量分為6等：0.12.5毫米為小雨，2.66毫米為中雨，6.112毫米為大雨，12.125毫米為暴雨，2

8、5.160毫米為大暴雨，大于60.1毫米為特大暴雨。若按此分級(jí)向公眾預(yù)報(bào)，如何在評(píng)價(jià)方法中考慮公眾的感受？三、名詞和符號(hào)說明表示預(yù)測點(diǎn)的經(jīng)度值表示預(yù)測點(diǎn)的緯度值表示預(yù)測點(diǎn)已知的預(yù)測值表示觀察站點(diǎn)的經(jīng)度值表示觀察站點(diǎn)的緯度值表示三次樣條插值的參數(shù)選項(xiàng)觀察站點(diǎn)的預(yù)測值表示兩種預(yù)測值矩陣中第i個(gè)元素的測量值表示實(shí)測矩陣中第i個(gè)元素的測量值表示殘差平方和表示均值誤差 表示預(yù)測準(zhǔn)確度預(yù)報(bào)正確次數(shù)預(yù)測次數(shù)四、模型假設(shè)：假設(shè)題目中全部數(shù)據(jù)真實(shí)可靠，忽略誤差；：假設(shè)觀測站所在位置的經(jīng)緯度準(zhǔn)確無誤；：假設(shè)天氣預(yù)報(bào)針對的位置在所給網(wǎng)格點(diǎn)附近；：假設(shè)雨量在各網(wǎng)點(diǎn)之間的變動(dòng)是連續(xù)的；五、問題分析針對問題1，我們將兩種

9、預(yù)測方法的所有預(yù)測值構(gòu)造成兩個(gè)以有序時(shí)間段對應(yīng)的預(yù)測值為列,以網(wǎng)格點(diǎn)的個(gè)數(shù)為行的2491×164矩陣，對于91觀測站點(diǎn)41天的實(shí)測值做同樣的處理，構(gòu)造成91×164的矩陣。這樣，繁瑣的數(shù)據(jù)經(jīng)過預(yù)處理后就整理成了三個(gè)矩陣。由于觀測站點(diǎn)相應(yīng)位置沒有兩種預(yù)測方法對應(yīng)的預(yù)測值，無法直接進(jìn)行評(píng)價(jià)，我們采用了三次樣條插值的方法進(jìn)行插值預(yù)處理，到了91個(gè)觀測站點(diǎn)兩種預(yù)測方法的相應(yīng)時(shí)刻的預(yù)測值，然后將兩種預(yù)測方法雨量預(yù)測值與雨量實(shí)測值進(jìn)行比較，從而判斷出兩種預(yù)測方法的準(zhǔn)確性。針對問題2，我們根據(jù)要求的雨量分級(jí)方式來考慮觀眾的感受。我們將問題1中91個(gè)觀測站點(diǎn)預(yù)測處理后雨量預(yù)測值構(gòu)成的兩個(gè)9

10、1×164矩陣和實(shí)際雨量觀測值構(gòu)成的91×164這個(gè)三個(gè)矩陣分別采用雨量等級(jí)記法構(gòu)造出三個(gè)新的矩陣，然后分別把兩個(gè)預(yù)測值構(gòu)成的降雨量等級(jí)矩陣和觀測值構(gòu)成的等級(jí)矩陣對應(yīng)元素相減并取絕對值，并進(jìn)行等級(jí)統(tǒng)計(jì)，再利用預(yù)測準(zhǔn)確度檢驗(yàn)法進(jìn)行判斷，準(zhǔn)確度越高說明我們預(yù)報(bào)的誤差越小，表明預(yù)測方法更準(zhǔn)確。六、模型建立1、數(shù)據(jù)預(yù)處理（1）針對問題1根據(jù)上面的分析，我們先對數(shù)據(jù)進(jìn)行預(yù)處理。處理方法為：把FORECAST文件夾中的第一種和第二種預(yù)測方式得到的數(shù)據(jù)分開兩個(gè)文件夾，分別以記事本格式按照日期的先后順序有序的導(dǎo)入Matlab的workspace工作空間中，然后建立m文件編輯公式將兩部分的數(shù)

11、據(jù)導(dǎo)入53×47×164的三維矩陣rain1和rain2中；把MEASURING文件夾中的數(shù)據(jù)以同樣方法導(dǎo)入91×7×41的三維矩陣temp中，然后建立循環(huán)將temp矩陣中每一層的后4列提取，另存入一個(gè)91×164的rain3矩陣；在命令窗中直接導(dǎo)入預(yù)測點(diǎn)的經(jīng)度和緯度存入矩陣lon和lat中，導(dǎo)入實(shí)測點(diǎn)的經(jīng)度和緯度存入矩陣lon1和lat1中，并對其作圖，如圖5-1。實(shí)現(xiàn)的matlab語句已呈現(xiàn)在附錄2.1中。圖5-1 預(yù)測點(diǎn)（彩色實(shí)線）與實(shí)測點(diǎn)（藍(lán)色孤點(diǎn)）由于三維矩陣無法用表格的形式呈現(xiàn)，我們分別截取了rain1和rain2矩陣的第一層呈現(xiàn)在

12、下表5-1和5-2中，rain3是二維矩陣，將其數(shù)據(jù)呈現(xiàn)在表5-3中：表5-1 預(yù)測方法一在6月18日第一個(gè)時(shí)間段的預(yù)測值構(gòu)造成的矩陣預(yù)測值B1B2B3.B47A10.03100.03210.0313.0.0081A20.03020.03100.0315.0.0078A30.03400.03410.0337.0.0079.A520.00880.00890.0084.0.0061A530.00890.00860.0086.0.0062表5-2 預(yù)測方法二在6月18日第一個(gè)時(shí)間段的預(yù)測值構(gòu)造成的矩陣預(yù)測值B1B2B3.B47A10.03190.02930.0322.0.0077A20.02880.

13、03110.0303.0.0085A30.03250.03230.0313.0.0082.A520.00950.00890.0079.0.0056A530.00850.00940.0093.0.0057表5-3 91個(gè)站點(diǎn)在41天共164個(gè)時(shí)段的雨量實(shí)測值預(yù)測值D1D2D3.D164P1000.0P2000.0P3000.0.P90004.1.0P91000.2.0（2）針對問題2由問題1我們可以整理得到91個(gè)觀測點(diǎn)41天的預(yù)測值和測量值對應(yīng)的兩個(gè)91×164矩陣，再根據(jù)問題2中氣象部門將降雨的等級(jí)分為6個(gè)等級(jí)的分法，把矩陣中相應(yīng)的降雨量值轉(zhuǎn)化為其所對應(yīng)等級(jí)值，其中，預(yù)測中的零全部記

14、為0，得到兩個(gè)預(yù)報(bào)等級(jí)矩陣，如下表5-4和5-5：表5-4 預(yù)測方法一在6月18日第一個(gè)時(shí)間段的預(yù)報(bào)等級(jí)構(gòu)造成的矩陣預(yù)報(bào)等級(jí)B1B2B3.B47A1000.0A2000.0A3000.0.A52000.0A53000.0表5-5 預(yù)測方法二在6月18日第一個(gè)時(shí)間段的預(yù)報(bào)等級(jí)構(gòu)造成的矩陣預(yù)報(bào)等級(jí)B1B2B3.B47A1000.0A2000.0A3000.0.A52000.0A53000.02、模型建立與求解（1）針對問題1由于91個(gè)觀測站點(diǎn)沒有相應(yīng)的預(yù)測值，因此不能夠直接對實(shí)測值進(jìn)行評(píng)價(jià)，屬于離散的散亂節(jié)點(diǎn)，我們利用插值基點(diǎn)為散亂節(jié)點(diǎn)的插值函數(shù)1在Matlab中進(jìn)行三次樣條插值處理，插值函數(shù)為：

15、 1 其中表示預(yù)測點(diǎn)的經(jīng)度值，表示預(yù)測點(diǎn)的緯度值，表示預(yù)測點(diǎn)的已知的預(yù)測值，表示觀察站點(diǎn)的緯度值，表示觀察站點(diǎn)的經(jīng)度值，表示三次樣條插值的參數(shù)選項(xiàng)，觀察站點(diǎn)的預(yù)測值。將91個(gè)觀測站點(diǎn)41天164個(gè)時(shí)段的雨量情況進(jìn)行預(yù)測之后，我們可以建立模型來評(píng)價(jià)這兩種預(yù)測方法。這里我們利用殘差平方和2以及平均誤差3來作為評(píng)價(jià)的標(biāo)準(zhǔn)： 2 3最后，我們根據(jù)和的值進(jìn)行評(píng)價(jià)，取值越大，表明預(yù)報(bào)的值準(zhǔn)確性越低。因此，殘差平方和與平均誤差值較小的一種預(yù)測方法作為較好的預(yù)報(bào)方法。利用公式（1），我們在Matlab中應(yīng)用編程求解，程序代碼見附錄2.2。求解之后得到91個(gè)觀測點(diǎn)41天164個(gè)時(shí)段的預(yù)測值，整理成91×

16、;164矩陣，然后把預(yù)測矩陣和實(shí)測矩陣對應(yīng)元素值相減取平方作殘差平方和，再作平均誤差，最后結(jié)果如下表5-6：表5-6 雨量預(yù)測計(jì)算結(jié)果數(shù)據(jù)表結(jié)果殘差平方和平均誤差預(yù)測方法一174290.000.4553預(yù)測方法二195580.000.4753由此可知，雨量預(yù)報(bào)方法一的準(zhǔn)確性更高一些。（2）針對問題2由問題1我們可以整理得到91個(gè)觀測點(diǎn)41天的預(yù)測值和測量值對應(yīng)的兩個(gè)91×164矩陣，再根據(jù)問題2中氣象部門將降雨的等級(jí)分為6個(gè)等級(jí)的分法，把矩陣中相應(yīng)的降雨量值轉(zhuǎn)化為其所對應(yīng)等級(jí)值，其中，預(yù)測中的零全部記為0，得到兩個(gè)預(yù)報(bào)等級(jí)矩陣，繼續(xù)利用殘差平方和以及平均誤差來作為評(píng)價(jià)的標(biāo)準(zhǔn)。殘差平方

17、和以及平均誤差數(shù)值越小，表明預(yù)報(bào)越準(zhǔn)確度越高，相應(yīng)公眾感受就越好。表5-7 雨量等級(jí)計(jì)算結(jié)果數(shù)據(jù)表結(jié)果殘差平方和平均誤差預(yù)測方法一27740.1730預(yù)測方法二28060.1745由結(jié)果可知，預(yù)測方法一的準(zhǔn)確度更高。七、模型優(yōu)化針對問題2，由于殘差平方和與平均誤差難以反映真實(shí)匯報(bào)的準(zhǔn)確度，我們將模型改進(jìn)優(yōu)化。由問題1我們可以整理得到91個(gè)觀測點(diǎn)41天的預(yù)測值和測量值對應(yīng)的兩個(gè)91×164矩陣，再根據(jù)問題2中氣象部門將降雨的等級(jí)分為6個(gè)等級(jí)的分法，把矩陣中相應(yīng)的降雨量值轉(zhuǎn)化為其所對應(yīng)等級(jí)值，得到兩個(gè)預(yù)報(bào)等級(jí)矩陣，將兩個(gè)預(yù)報(bào)等級(jí)矩陣與實(shí)測等級(jí)矩陣做差值運(yùn)算，得到兩個(gè)等級(jí)差矩陣，對等級(jí)差作

18、絕對值處理后，我們就可以從中進(jìn)行等級(jí)差統(tǒng)計(jì)。我們利用預(yù)測準(zhǔn)確度檢驗(yàn)法對兩種預(yù)報(bào)進(jìn)行評(píng)價(jià)。預(yù)測準(zhǔn)確度（）等于預(yù)報(bào)正確次數(shù)（）（即運(yùn)算之差為0的情況）和預(yù)測次數(shù)（）之比，即 4準(zhǔn)確度越高，表明預(yù)報(bào)越準(zhǔn)確度越高，相應(yīng)公眾感受就越好。雨量分級(jí)與統(tǒng)計(jì)程序見附錄4，數(shù)據(jù)統(tǒng)計(jì)結(jié)果見下表7-1：表7-1 兩種預(yù)報(bào)等級(jí)與實(shí)測等級(jí)的級(jí)差與比例表雨量等級(jí)差預(yù)報(bào)1等級(jí)差個(gè)數(shù)預(yù)報(bào)1等級(jí)差個(gè)數(shù)比例預(yù)報(bào)2等級(jí)差個(gè)數(shù)預(yù)報(bào)2等級(jí)差個(gè)數(shù)比例0124250.8326124040.8311121220.142221280.14262550.0037520.00353100.0007100.0007410.000120.0001500

19、0060000由表7-1可見預(yù)報(bào)1的預(yù)報(bào)準(zhǔn)確度為83.26%高于預(yù)報(bào)2的準(zhǔn)確度83.11%。預(yù)報(bào)1與實(shí)測雨量的符合度高，那么兩種預(yù)報(bào)方法中，公眾更易接受第一種預(yù)報(bào)方法。八、模型的評(píng)價(jià)與推廣優(yōu)點(diǎn)：本文思路自然流暢，較多地使用了原始數(shù)據(jù)，對91個(gè)站點(diǎn)的預(yù)測值比較準(zhǔn)確。缺點(diǎn)：認(rèn)為降雨量在不同的地區(qū)之間的變化有連續(xù)性，過于理想化。推廣：本文對降雨量的預(yù)測方法也可以運(yùn)用在工農(nóng)業(yè)生產(chǎn)中其他方面，比如一段時(shí)期某地區(qū)的氣溫變化、全國糧食產(chǎn)量預(yù)測等。本文應(yīng)用的散亂節(jié)點(diǎn)插值和殘差分析法可以應(yīng)用在很多方面。九、附錄附錄1：數(shù)據(jù)初始化>> File_FORECAST1=dir('預(yù)測方法一'

20、;)File_FORECAST1 =166x1 struct array with fields: name date bytes isdir>> File_FORECAST2=dir('預(yù)測方法二')File_FORECAST2 = 166x1 struct array with fields: name date bytesisdir>> File_MEASURING=dir('MEASURING')File_MEASURING = 43x1 struct array with fields: name date bytes isdi

21、r【mypro1.m】rain1=zeros(53,47,164);for n=1:164 filename=File_FORECAST1(n+2).name; rain1(:,:,n)=importdata(filename);end【mypro2.m】rain2=zeros(53,47,164);for n=1:164 filename=File_FORECAST2(n+2).name; rain2(:,:,n)=importdata(filename);end【mypro3.m】temp=zeros(91,7,41);for n=1:41 filename=File_MEASURING(

22、n+2).name; temp(:,:,n)=importdata(filename);endrain3=zeros(91,164);i=1;for j=1:41 for k=4:7 rain3(:,i)=temp(:,k,j); i=i+1; endend附錄2：問題1模型的建立與求解>> calc1pjwc = 0.4553ccpfh = 1.7429e+005>> calc2pjwc = 0.4753ccpfh = 1.9558e+005【calc1.m】lat;lon;lat1;lon1;ccpfh=0;pjwc=0;for i=1:164 wea=rain1(

23、:,:,i); wear=rain3(:,i); weap=griddata(lat,lon,wea,lat1,lon1); chazhi=wear-weap; pingfang=chazhi.2; juedui=abs(chazhi); ccpfh=ccpfh+sum(pingfang,1); pjwc=pjwc+sum(juedui);endpjwc=pjwc/(91*164);pjwcccpfh【calc2.m】lat;lon;lat1;lon1;ccpfh=0;pjwc=0;for i=1:164 wea=rain2(:,:,i); wear=rain3(:,i); weap=grid

24、data(lat,lon,wea,lat1,lon1); chazhi=wear-weap; pingfang=chazhi.2; juedui=abs(chazhi); ccpfh=ccpfh+sum(pingfang,1); pjwc=pjwc+sum(juedui);endpjwc=pjwc/(91*164);pjwcccpfh附錄3：問題2模型的建立與求解>> calc3pjwc1 = 0.1730ccpfh1 = 2774pjwc2 = 0.1745ccpfh2 = 2806【calc3.m】rank1=zeros(91,164);rank2=zeros(91,164);

25、rank=zeros(91,164);ccpfh1=0;pjwc1=0;ccpfh2=0;pjwc2=0;for i=1:164 wea1=rain1(:,:,i); wea2=rain2(:,:,i); wear=rain3(:,i); weap1(:,i)=griddata(lat,lon,wea1,lat1,lon1,'cubic'); weap2(:,i)=griddata(lat,lon,wea2,lat1,lon1,'cubic');endfor i=1:164 for j=1:91 if weap1(j,i)<0.1 rank1(j,i)=0

26、; elseif weap1(j,i)<=2.5 rank1(j,i)=1; elseif weap1(j,i)<=6.0 rank1(j,i)=2; elseif weap1(j,i)<=12.0 rank1(j,i)=3; elseif weap1(j,i)<=25.0 rank1(j,i)=4; elseif weap1(j,i)<=60.0 rank1(j,i)=5; else rank1(j,i)=6; end endendfor i=1:164 for j=1:91 if weap2(j,i)<0.1 rank2(j,i)=0; elseif w

27、eap2(j,i)<=2.5 rank2(j,i)=1; elseif weap2(j,i)<=6.0 rank2(j,i)=2; elseif weap2(j,i)<=12.0 rank2(j,i)=3; elseif weap2(j,i)<=25.0 rank2(j,i)=4; elseif weap2(j,i)<=60.0 rank2(j,i)=5; else rank2(j,i)=6; end endendfor i=1:164 for j=1:91 if rain3(j,i)<0.1 rank(j,i)=0; elseif rain3(j,i)&l

28、t;=2.5 rank(j,i)=1; elseif rain3(j,i)<=6.0 rank(j,i)=2; elseif rain3(j,i)<=12.0 rank(j,i)=3; elseif rain3(j,i)<=25.0 rank(j,i)=4; elseif rain3(j,i)<=60.0 rank(j,i)=5; else rank(j,i)=6; end endendchazhi1=rank1-rank;chazhi2=rank2-rank;pingfang1=chazhi1.2;pingfang2=chazhi2.2;juedui1=abs(chazhi1);juedui2=abs(chazhi2);ccpfh1=ccpfh1+sum(sum(pingfang1,1),2);ccpfh2=ccpfh2+sum(sum(pingfang2,1),2);pjwc1=pjwc1+sum(sum(juedui1),2);pjwc2=pjwc2+sum(sum(juedui2),2);pjwc1=pjwc1/(91*164);pjwc2=pjw