CN109614693A - 基于bp神经网络模型计算各月份天文辐射量的方法 - Google Patents

Abstract

本发明涉及一种基于BP神经网络模型计算各月份天文辐射量的方法，包括以下步骤：步骤S1：采集待计算地区DEM影像，对DEM影像进行预处理；步骤S2：提取DEM影像对应的输入矩阵；步骤S3：对输入矩阵进行归一化处理；步骤S4：构建并训练BP神经网络模型；步骤S5：将归一化后的输入矩阵输入训练后的BP神经网络模型，得到输出矩阵；步骤S6：对输出矩阵进行反归一化，得到待计算地区的月天文辐射数值矩阵。本发明计算月天文辐射数值方法简单，时间短，精度高。

Description

基于BP神经网络模型计算各月份天文辐射量的方法

技术领域

本发明涉及地理信息系统领域，具体涉及一种基于BP神经网络模型计算各月份天文辐射量的方法。

背景技术

太阳辐射是地球生命活动的能量源泉，也是大气中一切物理过程或现象的基本动力因子，其变化对全球地表能量流动、物质循环、生物生长以及人类活动等都具有决定性的影响。天文辐射是指地球表面不考虑大气影响所能接收到的太阳辐射，是地表实际入射太阳辐射的基础背景，也是辐射计算、太阳能资源评估和农业生产潜力估算等方面的重要天文参量。

目前天文辐射的主要计算方法如下：1)，使用数字高程模型(Digital ElevationModel，DEM)并根据纬度和日期(积日)信息构建数学模型计算天文辐射。但构建出的数学模型模型复杂，计算量非常大，效率非常低。2)使用BP(Back Propagation)神经网络模型估算其它气象站点的天文辐射。模型是以气象站点的观测值(如气温)作为输入值以相应站点的天文辐射观测值为输出值进行训练的。但所得到的BP神经网络模型仅可使用于数量极为有限气象站点，而在非气象站点的广大地区无法使用，极大地限制了对天文辐射的研究以及应用。

发明内容

有鉴于此，本发明的目的在于提供一种基于BP神经网络模型计算各月份天文辐射量的方法，解决了以往计算方法中模型复杂、计算量大的问题，以及难以计算任意地点天文辐射的问题。

为实现上述目的，本发明采用如下技术方案：

一种基于BP神经网络模型计算各月份天文辐射量的方法，包括以下步骤：

步骤S1：采集待计算地区DEM影像，并对DEM影像进行预处理；

步骤S2：根据预处理后的DEM影像，提取对应的输入矩阵；

步骤S3：对输入矩阵进行归一化处理；

步骤S4：构建并训练BP神经网络模型；

步骤S5：将归一化后的输入矩阵输入训练后的BP神经网络模型，得到BP神经网络模型的输出矩阵；

步骤S6：对输出矩阵进行反归一化，得到待计算地区的月天文辐射数值矩阵。

进一步的，所述待计算地区DEM影像为任意空间分辨率的DEM影像，为TIFF格式，所包含像素的数量为p1×p2，像素值表示该像素所对应位置的高程值。根据地理信息系统专业的表达习惯，以下将DEM影像的像素称为栅格。为了方便描述，将每个栅格进行编号(编号从1到p1×p2)，编号按照从上到下，后从左到右的顺序进行。

进一步的，所述步骤S1具体为：

步骤S11：采集待计算地区DEM影像；

步骤S12:将采集的DEM影像转换为数值矩阵；

步骤S12:填充数值矩阵中的空缺值,得到预处理后的DEM影像。

进一步的，所述步骤S2具体为：

步骤S21：提取待计算地区DEM影像每个栅格的高程、坡度、坡向因子以及每个栅格上下左右四个方向邻域栅格的高程、坡度和坡向因子；

步骤S22：将每个栅格的经度、纬度、高程、坡度、坡向、四个方向邻域栅格的高程、坡度、坡向以及所要计算的月份作为该栅格的18个因子，并将每个栅格对应的18个因子作为一个列向量；

步骤S23：将p1×p2个栅格对应的列向量依次按照栅格的编号组合为一个矩阵，该矩阵即为待计算地区DEM影像对应的输入矩阵Input_m×n

其中，m为输入矩阵Input的行数，n为Input的列数；i＝1，2，3，…，m，该矩阵中m＝18；j为DEM栅格的编号，j＝1，2，3，…，n，该矩阵中n＝p1×p2；

x_1，j和x_2，j为第j个栅格的经度和纬度，取该幅DEM影像中心栅格的经度和纬度；

x_3，j、x_4，j、x_5，j、x_6，j、x_7，j分别为第j个栅格和该栅格四个方向邻域栅格的高程；

x_8，j、x_9，j、x_10，j、x_11，j、x_12，j分别为第j个栅格和该栅格四个方向邻域栅格的坡度；

x_13，j、x_14，j、x_15，j、x_16，j、x_17，j分别为第j个栅格和该栅格四个方向邻域栅格的坡向；

x_18，j为所要计算的月份的编号，即1月份为1，2月份为2，...，12月份为12。

进一步的，所述步骤S3具体为：

对输入矩阵进行归一化，使每一行的元素值限定在[-1,1]区间范围内；

可使用如下归一化公式：

式中，i＝1，2，3，…，18；j＝1，2，3，…，n；x_i，j是矩阵Input第i行第j列元素的值；是矩阵Input第i行所有元素中的最大值，是矩阵Input第i行所有元素中的最小值，P_i，j是归一化后的输入矩阵P第i行第j列元素的值，取值范围是[-1，1]。

进一步的，所述步骤S4具体为：

步骤S41：选取训练样本，包括N幅DEM影像，编号为1-N，并在每个DEM影像随机抽取M个栅格做为采样栅格；

步骤S42：构建BP神经网络模型输入矩阵：

构建编号为1的DEM影像含有1月份编号的矩阵l₁ ¹，矩阵l₁ ¹的第18行元素均为1(对应1月份的编号)；

式中，l₁ ¹为18行M列的矩阵；l₁ ¹的下标表示DEM影像的编号，上标表示月份编号；d＝1，2，3，...，M；矩阵l₁ ¹中第d列的元素x_1d，x_2d，...，x_17,d依次为第d个采样栅格的经度值、纬度值、高程值、四个方向邻域栅格的高程值、坡度值、四个方向邻域栅格的坡度值、坡向值、四个方向邻域栅格的坡向值；

构建上述所有采样栅格含有2月份编号的矩阵l₁ ²，矩阵l₁ ²的第18行元素均为2(对应2月份的编号)；

上述所有采样栅格所对应其余月份的矩阵仅改变矩阵第18行的月份编号，得到l₁ ³，l₁ ⁴，...，l₁ ¹²；

将矩阵l₁ ¹，l₁ ²，...，l₁ ¹²按如下方式进行组合，得到编号为1的DEM影像中采样栅格对应的矩阵l₁，

构建全部N幅DEM影像中采样栅格的矩阵l₁，l₂，...，l_N，将其按如下方式进行组合，得到BP神经网络模型的输入矩阵P，

P＝[l₁，l₂，...，l_N]

式中，P为18行，M×N列的矩阵；

步骤S43：构建BP神经网络模型期望输出矩阵；

该期望输出矩阵的每一个元素值为该元素所对应栅格的月天文辐射值，是通过数学模型推导计算得到的；

构建编号为1的DEM影像中所有采样栅格1月份的天文辐射矩阵t₁ ¹，

式中，t₁ ¹为1行M列的矩阵，下标表示DEM影像的编号，上标表示月份；y₁，y₂，...，y_M分别为上述所有采样栅格1月份的天文辐射值；

构建上述所有采样栅格其余月份的天文辐射矩阵t₁ ²，...，t₁ ¹²，将其按如下方式进行组合，得到编号为1的DEM影像中所有采样栅格12个月份的天文辐射矩阵t₁，

构建全部N幅DEM影像所有采样栅格的天文辐射矩阵t₁，t₂，...，t_N，将其按如下方式进行组合，得到BP神经网络模型的输出矩阵T，

T＝[t₁，t₂，...，t_N]

式中，T为1行M×N列的矩阵；

步骤S44：构建由一个输入层、三个隐含层、一个输出层构成的BP神经网络模型框架，设置网络基本参数：如最大训练次数、训练目标最小误差、学习率、最大失败次数；

步骤S45：训练BP神经网络模型，得到训练后的BP神经网络模型。

进一步的，所述步骤S45具体为：

步骤S451：初始化各层权值矩阵和阈值向量；

步骤S452：将步骤S42中构成的输入矩阵和步骤S43构成的期望输出矩阵分别归一化，同时需保存期望输出矩阵归一化之前的最大值和最小值，记为y_max和y_min；

步骤S453：将归一化后的输入矩阵输入BP神经网络模型的输入层，计算出其输出矩阵，求出训练目标误差E，其公式为：

式中，n为训练样本(即输入矩阵的列数)的数量；τ_k为归一化后的期望输出矩阵第k个元素的值；o_k为BP神经网络模型的输出矩阵第k个元素的值。

步骤S454：判断误差E是否达到误差要求，若E达到误差要求，即训练结束，转向步骤S455；若没有达到要求，则反向修正BP神经网络模型的各层权值矩阵和阈值向量，再转向步骤S453；

步骤S455：保存各层权值矩阵和阈值向量，得到训练后的BP神经网络模型。

进一步的，所述步骤S6具体为：

步骤S61：根据步骤S3的归一化方法，对BP神经网络模型的输出矩阵进行反归一化；

步骤S62：将单行矩阵Output重新组合为p1行p2列的矩阵，组合时每p1个元素为一列，即第1到第p1个元素为第一列，第p1+1到第2*p1个元素为第二列，以此类推，组合成一个p1行p2列的矩阵。

本发明与现有技术相比具有以下有益效果：

本发明构建出一种BP神经网络模型，能够快速计算出任意地点的天文辐射，解决了以往计算方法中模型复杂、计算量大的问题，以及难以计算任意地点天文辐射的问题；从而为深入了解天文辐射的空间分布特征，为计算地表实际入射太阳辐射(直接辐射、散射辐射、反射辐射)奠定基础，具有重要的科研价值和实用价值。

附图说明

图1是本发明每个栅格与四个方向邻域栅格的方位关系图；

图2是本发明的操作流程图；

图3是本发明BP神经网络模型的结构图；

图4是本发明BP神经网络模型训练流程图；

图5是本发明实施例中北京市顺义区DEM影像图；

图6是本发明实施例中Output_3601×3601的图像。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

本实施例提供一种基于BP神经网络模型计算各月份天文辐射量的方法，包括以下步骤：

步骤S1：采集待计算地区DEM影像，并对DEM影像进行预处理；

步骤S2：根据预处理后的DEM影像，提取对应的输入矩阵；

步骤S3：对输入矩阵进行归一化处理；

步骤S4：构建并训练BP神经网络模型；

步骤S5：将归一化后的输入矩阵输入训练后的BP神经网络模型，得到BP神经网络模型的输出矩阵；

步骤S6：对输出矩阵进行反归一化，得到待计算地区的月天文辐射数值矩阵。

在本实施例中，所述待计算地区DEM影像为30m空间分辨率的DEM影像，为TIFF格式，每幅DEM影像的经度与纬度间隔均为1°，所包含像素的数量为3601×3601，像素值表示该像素所对应位置的高程值；以下根据地理信息系统专业的表达习惯，将DEM影像的像素称为栅格，一幅DEM影像共有12967201(3601×3601)个栅格。为了方便描述，将每个栅格进行编号(编号从1到12967201)，编号按照先从上到下，后从左到右的顺序进行(第1行第1列的栅格编号为1，第2行第1列的栅格编号为2，第3行第1列的栅格编号为3，…，第3601行3601列的栅格编号为12967201)。另外，计算各月天文辐射时需要DEM影像每个栅格的经纬度，由于每幅DEM影像的经度和纬度间隔仅为1°，为方便计算，各栅格的经纬度坐标均取DEM影像中心栅格的经纬度坐标。

举例说明：某幅DEM影像中心栅格的经纬度坐标为(2.0333E，0.7069N)(弧度制坐标)，计算中该幅DEM影像的每个栅格的经度和纬度均取2.0333和0.7069。

在本实施例中，所述步骤S1中对DEM影像进行预处理的具体步骤为：

步骤S11：将采集的DEM影像转换为数值矩阵，使用MATLAB软件，则调用imread函数进行读取即可；

步骤S12：采用邻近值填充法填充数值矩阵中的空缺值。

在本实施例中，参照图1，所述步骤S2具体为：

步骤S21：使用MATLAB软件构建数学模型，提取待计算地区DEM影像每个栅格的高程、坡度、坡向因子以及每个栅格上下左右四个方向邻域栅格的高程、坡度和坡向因子；每个栅格与四个方向邻域栅格的方位关系见图1；

图1中，I，J为DEM每个栅格的行列号，I，J＝1，2，3，...，3601；Ω_IJ为第I行J列栅格的高程值；Ω_I，J-1，Ω_I，J+1，Ω_I-1，J，Ω_I+1，J分别为第I行J列栅格四个方向邻域栅格的高程值；

步骤S22：将每个栅格的经度、纬度、高程、坡度、坡向、四个方向邻域栅格的高程、坡度、坡向以及所要计算的月份(1月份为1，2月份为2，...，12月份为12)作为该栅格的18个因子，并将每个栅格对应的18个因子作为一个列向量；

步骤S23：将12967201个栅格对应的列向量依次按照栅格的编号组合为一个矩阵(编号为1的栅格对应的列向量为矩阵的第1列，编号为2的栅格对应的列向量为矩阵的第2列，…，编号为12967201的栅格对应的列向量为矩阵的第12967201列)，该矩阵即为待计算地区DEM影像对应的输入矩阵Input_m×n

其中，m为输入矩阵Input的行数，n为输入矩阵Input的列数；i＝1，2，3，…，m，该矩阵中m＝18；j为DEM栅格的编号，j＝1，2，3，…，n，该矩阵中n＝12967201；

x_1，j和x_2，j为第j个栅格的经度和纬度，取该幅DEM影像中心栅格的经度和纬度；

x_3，j、x_4，j、x_5，j、x_6，j、x_7，j分别为第j个栅格和该栅格四个方向邻域栅格的高程；

x_8，j、x_9，j、x_10，j、x_11，j、x_12，j分别为第j个栅格和该栅格四个方向邻域栅格的坡度；

x_13，j、x_14，j、x_15，j、x_16，j、x_17，j分别为第j个栅格和该栅格四个方向邻域栅格的坡向；

x_18，j为所要计算的月份的编号，即1月份为1，2月份为2，...，12月份为12。每个地区的月天文辐射需逐月进行计算，即每个地区需要构建12个输入矩阵，这12个输入矩阵只需更改第18行的月份编号即可，然后依次输入BP神经网络模型，即可得到地区各月的天文辐射。

参照图3，在本实施例中，所述步骤S3具体为：

矩阵中元素的归一化公式为：

式中，i＝1，2，3，…，18；j＝1，2，3，…，12967201；x_i，j是矩阵Input第i行第j列元素的值；是矩阵Input第i行所有元素中的最大值，是矩阵Input第i行所有元素中的最小值，P_i，j是归一化后的输入矩阵P第i行第j列元素的值，取值范围是[-1，1]。

以上归一化过程可以使用MATLAB软件的mapminmax函数完成，函数的调用格式如下：

[P，Pn]＝mapminmax(Input，pmin，pmax)

式中Input为输入矩阵，P为输出矩阵，Pn为调用mapminmax函数时得到的参数，主要包含了归一化前矩阵每一行全部元素的最小值和最大值，以及输出矩阵P的最小值pmin和最大值pmax(pmin默认为-1，pmax默认为1)。

在本实施例中，参照图3，所述步骤S4具体为：

步骤S41：选取训练样本，选取200幅DEM影像，编号为1-200，并在每个DEM影像随机抽取10个栅格，此后称抽取的栅格为‘采样栅格’；

步骤S42：构建BP神经网络模型输入矩阵：

构建编号为1的DEM影像含有1月份编号的矩阵l₁ ¹，矩阵l₁ ¹的第18行元素均为1(对应1月份的编号)；

式中，l₁ ¹为18行10列的矩阵；l₁ ¹的下标表示DEM影像的编号，上标表示月份编号；d＝1，2，3，...，10；矩阵l₁ ¹中第d列的元素x_1d，x_2d，...，x_17,d依次为第d个采样栅格的经度值、纬度值、高程值、四个方向邻域栅格的高程值、坡度值、四个方向邻域栅格的坡度值、坡向值、四个方向邻域栅格的坡向值；

构建上述所有采样栅格含有2月份编号的矩阵l₁ ²，矩阵l₁ ²的第18行元素均为2(对应2月份的编号)；

上述所有采样栅格所对应其余月份的矩阵仅改变矩阵第18行的月份编号，得到l₁ ³，l₁ ⁴，...，l₁ ¹²；

将矩阵l₁ ¹，l₁ ²，...，l₁ ¹²按如下方式进行组合，得到编号为1的DEM影像中采样栅格对应的矩阵l₁，

构建全部200幅DEM影像中所有采样栅格的矩阵，得到l₁，l₂，...，l₂₀₀，将其按照如下方式进行组合，得到BP神经网络模型的输入矩阵P，

P＝[l₁，l₂，...，l₁₂]

式中，P为18行，24000列的矩阵；

步骤S43：构建BP神经网络模型期望输出矩阵；

该期望输出矩阵的每一个元素值为该元素所对应栅格的月天文辐射值，是通过数学模型推导计算得到的；

构建编号为1的DEM影像中所有采样栅格1月份的天文辐射矩阵t₁ ¹，

式中，t₁ ¹为1行10列的矩阵，下标表示DEM影像的编号，上标表示月份；y₁，y₂，...，y₁₀分别为编号为1的DEM影像中所有采样栅格1月份的天文辐射值；

构建上述所有采样栅格其余月份的天文辐射矩阵t₁ ²，...，t₁ ¹²，将其按如下方式进行组合，得到编号为1的DEM影像中所有采样栅格12个月份的天文辐射矩阵t₁，

构建全部200幅DEM影像所有采样栅格的天文辐射矩阵t₁，t₂，...，t₂₀₀，将其按如下方式进行组合，得到BP神经网络模型的输出矩阵T，

T＝[t₁，t₂，...，t₂₀₀]

式中，T为1行24000列的矩阵；

步骤S44：构建由一个输入层、三个隐含层、一个输出层构成的BP神经网络模型框架，设置网络基本参数：最大训练次数为10000；训练目标最小误差设置为100；学习率设置为0.01；最大失败次数设置为20；

步骤S45：训练BP神经网络模型，得到训练后的BP神经网络模型；

在本实施例中，参照图4，所述步骤S45具体为：

步骤S451：初始化各层权值矩阵和阈值向量；

步骤S452：将上述输入矩阵和期望输出矩阵分别归一化，同时需保存期望输出矩阵归一化之前的最大值和最小值，记为y_max和y_min。使用MATLAB软件，则分别调用归一化函数mapminmax，同时需保存期望输出矩阵归一化时所获得的参数，记为Tn，该参数中包含了期望输出矩阵归一化之前所有元素的最大值和最小值；

步骤S453：将归一化后的输入矩阵输入BP神经网络模型的输入层，计算出其输出矩阵，求出训练目标误差E，其公式为：

式中，n为训练样本(即输入矩阵的列数)的数量；τ_k为归一化后的期望输出矩阵第k个元素的值；o_k为BP神经网络模型的输出矩阵第k个元素的值。

步骤S454：判断误差E是否达到误差要求，在本实施例中，将最小误差e设置为100，若E<e，即训练结束，转向步骤S455；若E>＝e，则反向修正BP神经网络模型的各层权值矩阵和阈值向量，再转向步骤S453；

步骤S455：保存各层权值矩阵和阈值向量，得到训练后的BP神经网络模型。

在本实施例中，所述步骤S5具体为：参照图3，

图2中，P为归一化后的输入矩阵；W_l和B_l为BP神经网络模型第ik层隐含层的权值矩阵和阈值向量，r为输入层神经元数量即P的行数，该模型中r＝18；s₁，s₂，s₃为隐含层神经元数量，该模型中s₁＝s₂＝s₃＝15；s₄为输出层神经元数量，该模型中s₄＝1；q为P的列数即DEM中栅格的个数，该模型中q＝12967201；h1＝1，2，3，…，s₁；h2＝1，2，3，…，s₂；h3＝1，2，3，…，s₃；h4＝1；

F₁，F₂，F₃为隐含层的激活函数，是非线性函数：

式中l＝1，2，3；A_l为BP神经网络模型第l层隐含层的输出矩阵，A_l-1为BP神经网络模型第l层隐含层的输入矩阵，二者的取值范围是[-1，1]，当l＝1时，A₀＝P，P为BP神经网络模型输入层的输入矩阵；

F₄为输出层的激活函数，是线性函数，公式为：

Y＝F₄(W_l×X+B_l)＝W_l×X+B_l

式中，Y为BP神经网络模型输出层的输出矩阵，X为BP神经网络模型输出层的输入矩阵；该式中l＝4，则O＝W₄*A₃+B₄＝A₄

O和A₄为BP神经网络模型输出层的输出矩阵(1行12967021列)，A₃为BP神经网络模型输出层的输入矩阵(1行12967021列)。

在本实施例中，所述步骤S6具体为：

步骤S61：反归一化公式为：

式中j＝1，2，3，…，12967201；output_j为反归一化后的值；O_j为BP神经网络模型的输出值；y_max和y_min是BP神经网络训练过程中，期望输出矩阵中的最大值和最小值。使用MATLAB软件，则调用mapminmax函数完成反归一化，函数的调用格式如下：

Output＝mapminmax(‘reverse’，O，Tn)

式中，Output_1×12967201为反归一化后的矩阵；O_1×12967201为归一化前的矩阵；参数Tn为步骤S452中所获参数，它包含了训练过程中期望输出矩阵归一化之前的最大值和最小值；

步骤S62：将单行矩阵Output_1×12967201重新组合为3601×3601的矩阵，组合时每3601个元素作为一列(第1到第3601个元素作为第一列，第3602到第7202个元素作为第二列，…，第12963421到第12967021个元素作为第3601列)。

实施例1：

参照图5和图6，假设需计算北京市顺义区8月份的天文辐射。该区域的经纬度范围约为40°N～40°18＇N，116°28＇E～116°58＇E(见图4中边线所围成的区域)；

具体计算步骤如下：

1)根据该区域的经纬度范围，下载包含该区域的DEM影像，如图4所示。由于DEM影像的经纬度范围为40°N～41°N，116°E～117°E，跨度仅为1°，为方便计算，各栅格的经纬度坐标均取DEM影像中心栅格的经度和纬度，即2.0333和0.7069。将DEM影像转换为数值矩阵；随后需对DEM影像中存在高程空缺值的栅格进行填充(此处使用邻近法进行填充，即将空缺值处填入最邻近栅格的高程值)。

2)如步骤S2所述，提取输入矩阵Input_18×12967201；

3)如步骤S3所述，对Input_18×12967201进行归一化，得到P_18×12967201矩阵；

4)如步骤S5所述，将P_18×12967201输入如图3所示训练后的BP神经网络模型，计算出输出矩阵O_1×12967201。

5)如步骤S6所述，对输出矩阵O_1×12967201进行反归一化得到Output_1×12967201。将单行Output_1×12967201矩阵重新组合为3601×3601的矩阵，得到该地区的8月份天文辐射矩阵Output_3601×3601。本例中Output_3601×3601以图像方式给出，见图6。

以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本发明的涵盖范围。

Claims (8)

1.一种基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于，包括以下步骤：

步骤S1：采集待计算地区DEM影像，并对DEM影像进行预处理；

步骤S2：根据预处理后的DEM影像，提取对应的输入矩阵；

步骤S3：对输入矩阵进行归一化处理；

步骤S4：构建并训练BP神经网络模型；

步骤S5：将归一化后的输入矩阵输入训练后的BP神经网络模型，得到BP神经网络模型的输出矩阵；

步骤S6：对输出矩阵进行反归一化，得到待计算地区的月天文辐射数值矩阵。

2.根据权利要求1所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述待计算地区DEM影像为任意空间分辨率的DEM影像，为TIFF格式，所包含像素的数量为p1×p2，像素值表示该像素所对应位置的高程值。

3.根据权利要求1所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S1具体为：

步骤S11：采集待计算地区DEM影像；

步骤S12:将采集的DEM影像转换为数值矩阵；

步骤S12:填充数值矩阵中的空缺值,得到预处理后的DEM影像。

4.根据权利要求2所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S2具体为：

步骤S21：提取待计算地区DEM影像每个栅格的高程、坡度、坡向因子以及每个栅格上下左右四个方向邻域栅格的高程、坡度和坡向因子；

步骤S22：将每个栅格的经度、纬度、高程、坡度、坡向、四个方向邻域栅格的高程、坡度、坡向以及所要计算的月份作为该栅格的18个因子，并将每个栅格对应的18个因子作为一个列向量；

步骤S23：将p1×p2个栅格对应的列向量依次按照栅格的编号组合为一个矩阵，该矩阵即为待计算地区DEM影像对应的输入矩阵Input_m×n

其中，m为输入矩阵Input的行数，n为Input的列数；i＝1，2，3，…，m，该矩阵中m＝18；j为DEM栅格的编号，j＝1，2，3，…，n，该矩阵中n＝p1×p2；

x_1，j和x_2，j为第j个栅格的经度和纬度，取该幅DEM影像中心栅格的经度和纬度；

x_3，j、x_4，j、x_5，j、x_6，j、x_7，j分别为第j个栅格和该栅格四个方向邻域栅格的高程；

x_8，j、x_9，j、x_10，j、x_11，j、x_12，j分别为第j个栅格和该栅格四个方向邻域栅格的坡度；

x_13，j、x_14，j、x_15，j、x_16，j、x_17，j分别为第j个栅格和该栅格四个方向邻域栅格的坡向；

x_18，j为所要计算的月份的编号，即1月份为1，2月份为2，...，12月份为12。

5.根据权利要求1所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S3具体为：对输入矩阵进行归一化，使每一行的元素值限定在[-1,1]区间范围内。

6.根据权利要求1所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S4具体为：

步骤S41：选取训练样本，包括N幅DEM影像，编号为1-N，并在每个DEM影像中随机抽取M个栅格做为采样栅格；

步骤S42：构建BP神经网络模型输入矩阵：

构建编号为1的DEM影像中所有采样栅格含有1月份编号的矩阵l₁ ¹，矩阵l₁ ¹的第18行元素均为1；

式中，l₁ ¹为18行M列的矩阵；l₁ ¹的下标表示DEM影像的编号，上标表示月份编号；d＝1，2，3，...，M；矩阵l₁ ¹中第d列的元素x_1d，x_2d，...，x_17,d依次为第d个采样栅格的经度值、纬度值、高程值、四个方向邻域栅格的高程值、坡度值、四个方向邻域栅格的坡度值、坡向值、四个方向邻域栅格的坡向值；

构建上述所有采样栅格含有2月份编号的矩阵l₁ ²，矩阵l₁ ²的第18行元素均为2(对应2月份的编号)；

上述所有采样栅格所对应其余月份的矩阵仅改变矩阵第18行的月份编号，得到l₁ ³，l₁ ⁴，...，l₁ ¹²；

将矩阵l₁ ¹，l₁ ²，...，l₁ ¹²按如下方式进行组合，得到编号为1的DEM影像中所有采样栅格对应的矩阵l₁，

构建得到全部N幅DEM影像中所有采样栅格的矩阵l₁，l₂，...，l_N，将其按如下方式进行组合，得到BP神经网络模型的输入矩阵P，

P＝[l₁，l₂，...，l_N]

式中，P为18行，M×N列的矩阵；

步骤S43：构建BP神经网络模型期望输出矩阵；

构建编号为1的DEM影像中所有采样栅格1月份的天文辐射矩阵t₁ ¹，

式中，t₁ ¹为1行M列的矩阵，下标表示DEM影像的编号，上标表示月份编号；y₁，y₂，...，y_M分别为上述中所有采样栅格1月份的天文辐射值；

构建上述所有采样栅格其余月份的天文辐射矩阵t₁ ²，...，t₁ ¹²，将其按如下方式进行组合，得到编号为1的DEM影像中所有采样栅格12个月份的天文辐射矩阵t₁，

构建全部N幅DEM影像所有采样栅格的天文辐射矩阵t₁，t₂，...，t_N，按如下方式进行组合，得到BP神经网络模型的输出矩阵T，

T＝[t₁，t₂，...，t_N]

式中，T为1行M×N列的矩阵；

步骤S44：构建由一个输入层、三个隐含层、一个输出层构成的BP神经网络模型框架，设置网络基本参数：如最大训练次数、训练目标最小误差、学习率、最大失败次数；

步骤S45：训练BP神经网络模型，得到训练后的BP神经网络模型。

7.根据权利要求6所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S45具体为：

步骤S451：初始化各层权值矩阵和阈值向量；

步骤S452：将步骤S42中构成的输入矩阵和步骤S43构成的期望输出矩阵分别归一化，同时需保存期望输出矩阵归一化前的最大值和最小值，记为y_max和y_min；

步骤S453：将归一化后的输入矩阵输入BP神经网络模型的输入层，计算出其输出矩阵，求出训练目标误差E，其公式为：

式中，n为训练样本(即输入矩阵的列数)的数量；τ_k为归一化后的期望输出矩阵第k个元素的值；o_k为BP神经网络模型的输出矩阵第k个元素的值。

步骤S454：判断误差E是否达到误差要求，若E达到误差要求，即训练结束，转向步骤S455；若没有达到要求，则反向修正BP神经网络模型的各层权值矩阵和阈值向量，再转向步骤S453；

步骤S455：保存各层权值矩阵和阈值向量，得到训练后的BP神经网络模型。

8.根据权利要求1所述的基于BP神经网络模型计算各月份天文辐射量的方法，其特征在于：所述步骤S6具体为：

步骤S61：根据步骤S3的归一化方法，对BP神经网络模型的输出矩阵进行反归一化；

步骤S62：将单行矩阵Output重新组合为p1行p2列的矩阵，组合时每p1个元素为一列。