技术领域technical field
本发明属于农业遥感领域，具体涉及一种基于尺度转换和数据同化的农作物产量估测方法。The invention belongs to the field of agricultural remote sensing, and in particular relates to a crop yield estimation method based on scale conversion and data assimilation.
背景技术Background technique
利用数据同化技术把遥感反演参数融合到作物机理过程模型，是当前改进区域作物生长模拟精度和提高作物估产精度的重要途径。在田间尺度上，基于作物光合、呼吸、蒸腾、营养等机理过程的作物生长模型依靠其内在的物理过程和动力学机制，可以准确模拟作物对象以“天”为时间步长的连续演进和单点尺度上作物的生长发育状况及产量。而当作物模型扩展应用到区域尺度时，由于地表、近地表环境非均匀性，导致了作物模型中的初始条件、土壤参数、作物参数、气象强迫因子的空间分布的不确定性和资料获取的困难。卫星遥感具有空间连续和时间动态变化的优势，能够有效解决区域作物参数获取困难这一瓶颈。然而遥感对地观测由于受卫星时空分辨率等因素的制约，还不能真正揭示作物生长发育和产量形成的内在过程机理、个体生长发育状况及其与环境气象条件的关系，而这正是作物模型的优势所在。数据同化技术通过耦合遥感观测和作物模型，能够实现两者的优势互补，对提高区域作物产量估测具有巨大的潜力。Using data assimilation technology to integrate remote sensing inversion parameters into crop mechanism process models is an important way to improve the accuracy of regional crop growth simulation and crop yield estimation. On the field scale, crop growth models based on crop photosynthesis, respiration, transpiration, nutrition and other mechanistic processes rely on their inherent physical processes and dynamic mechanisms, which can accurately simulate the continuous evolution and single The growth and development status and yield of crops on the point scale. However, when the crop model is extended to the regional scale, due to the heterogeneity of the surface and near-surface environments, the uncertainty of the initial conditions, soil parameters, crop parameters, and meteorological forcing factors in the crop model and the uncertainty of data acquisition difficulty. Satellite remote sensing has the advantages of spatial continuity and temporal dynamic change, and can effectively solve the bottleneck of difficulty in obtaining regional crop parameters. However, due to the constraints of satellite space-time resolution and other factors, remote sensing earth observation cannot really reveal the internal process mechanism of crop growth and yield formation, the individual growth and development status and its relationship with environmental meteorological conditions, which is precisely the crop model. advantage. By coupling remote sensing observations and crop models, data assimilation technology can realize the complementary advantages of the two, and has great potential for improving regional crop yield estimation.
由于高时间分辨率遥感数据能扑捉作物生长发育变化以及物候信息，因此对农作物长势监测与估产具有非常重要的意义，然而，目前在轨运行的高时间分辨率遥感卫星数据往往空间分辨率都很低，例如AVHRR，MODIS，MERIS和SPOT vegetation空间分辨率为250m-1km尺度。这种低空间分辨率遥感数据，进一步增大了像元异质性，造成了遥感观测尺度与作物模型模拟尺度之间的不匹配，极大限制了数据同化模型的精度。另外，标准的MODIS LAI产品是面向全球尺度所有的植被类型开发，并不是针对农作物设计开发的。研究表明，MODIS LAI在冬小麦区域存在严重的低估现象。因此，直接同化标准MODIS LAI会导致同化后的LAI和产量不切实际的偏低。基于地统计（变异函数）尺度纠正的方法和像元分解降尺度的方法是处理像元异质性和尺度效应的两种手段。然而，这两种方法的一个共同的不足之处是要求一系列的作物生育期的中高空间分辨率数据，并且对异质性像元建模是一个非常复杂的问题，要求严格的近似和一些先验知识，这在实际应用中，不易获取。一种有潜力的方法是基于地面采样LAI、中高分辨率遥感数据和低分辨率遥感数据，构建二级尺度转换方法，生成尺度修正后的LAI。同化尺度修正后的LAI，将会较大提高区域冬小麦产量的估测。以克服直接同化MODIS LAI产品带来的误差是应用数据同化模型估产实践中亟待解决的关键问题。Because high-time-resolution remote sensing data can capture crop growth and development changes and phenological information, it is of great significance to crop growth monitoring and yield estimation. Very low, such as AVHRR, MODIS, MERIS and SPOT vegetation spatial resolution is 250m-1km scale. This low spatial resolution remote sensing data further increases the pixel heterogeneity, causing a mismatch between the remote sensing observation scale and the crop model simulation scale, which greatly limits the accuracy of the data assimilation model. In addition, the standard MODIS LAI product is developed for all vegetation types on a global scale, not for crops. The research shows that MODIS LAI is seriously underestimated in the winter wheat area. Thus, direct assimilation of standard MODIS LAI would result in unrealistically low assimilated LAI and yield. The method based on geostatistics (variogram) scale correction and the method of pixel decomposition and downscaling are two means to deal with pixel heterogeneity and scale effect. However, a common shortcoming of these two methods is that a series of crop growth period data with medium and high spatial resolution is required, and modeling heterogeneous pixels is a very complex problem, requiring strict approximation and some Prior knowledge, which is not easy to obtain in practical applications. A potential method is to construct a two-level scaling method based on ground-sampled LAI, medium-high resolution remote sensing data, and low-resolution remote sensing data to generate a scale-corrected LAI. The LAI after assimilation scale correction will greatly improve the estimation of regional winter wheat yield. To overcome the errors caused by direct assimilation of MODIS LAI products is a key problem to be solved in the practice of applying data assimilation models to estimate production.
发明内容Contents of the invention
为解决现有技术中存在的如下问题：“如何构建数据同化模型中的尺度转换方法，同化时间序列尺度修正后的LAI到作物生长模型，以减少遥感观测与作物模型之间尺度不匹配造成的数据同化模型误差”，本发明提供一种基于尺度转换和数据同化的农作物产量估测方法。In order to solve the following problems in the existing technology: "How to construct the scale conversion method in the data assimilation model, assimilate the time series scale corrected LAI to the crop growth model, so as to reduce the scale mismatch between the remote sensing observation and the crop model. Data assimilation model error", the present invention provides a crop yield estimation method based on scale conversion and data assimilation.
本发明提供一种基于尺度转换和数据同化的农作物产量估测方法，具体步骤如下：The present invention provides a crop yield estimation method based on scale conversion and data assimilation, the specific steps are as follows:
S1收集研究区内的气象参数、作物参数、土壤参数和管理参数并以之进行WOFOST作物模型的单点尺度标定；之后基于气象参数进行区域化，完成WOFOST作物模型的空间化；S1 Collect meteorological parameters, crop parameters, soil parameters and management parameters in the research area and use them to perform single-point scale calibration of WOFOST crop model; then perform regionalization based on meteorological parameters to complete the spatialization of WOFOST crop model;
S2基于多时相的遥感影像，结合地面调查获得的作物类型样本点，结合待测农作物的物候特征，构建分类规则，获取作物类型分布图和待测农作物纯度百分比图；S2 is based on multi-temporal remote sensing images, combined with the crop type sample points obtained from ground surveys, combined with the phenological characteristics of the crops to be tested, to construct classification rules, and to obtain the distribution map of crop types and the purity percentage map of the crops to be tested;
S3对待测农作物生育期时间序列MODIS LAI曲线进行滤波，以消除数据缺失和云污染的影响；基于地面实测样本点的实测LAI与多时相植被指数（TM VI）构建回归统计模型，获得区域30米的TM LAI；S3 Filter the time series MODIS LAI curve of the growth period of the crops to be measured to eliminate the influence of data loss and cloud pollution; build a regression statistical model based on the measured LAI of the sample points on the ground and the multi-temporal vegetation index (TM VI), and obtain an area of 30 meters TM LAI;
S4融合地面实测样本点LAI、多时相TM LAI和滤波后的MODIS LAI信息，构建二级尺度转换模型，生成时间序列尺度调整LAI；S4 integrates the ground-measured sample point LAI, multi-temporal TM LAI and filtered MODIS LAI information to construct a secondary scale conversion model and generate a time series scale-adjusted LAI;
S5基于地面样本点关键物候期的LAI，计算遥感观测和作物模型模拟的标准差，获得待测农作物关键物候期的遥感观测误差和作物模型误差；S5 calculates the standard deviation of remote sensing observations and crop model simulations based on the LAI of the key phenological periods of the ground sample points, and obtains the remote sensing observation errors and crop model errors of the key phenological periods of the crops to be measured;
S6以S4获得的尺度调整LAI和作物模型模拟LAI，并引入遥感观测误差和作物模型误差，构建四维变分代价函数，使用最优化算法最小化代价函数，通过多次迭代，获得优化作物模型参数；S6 uses the scale obtained in S4 to adjust the LAI and the crop model to simulate the LAI, and introduces the remote sensing observation error and the crop model error to construct a four-dimensional variational cost function, use the optimization algorithm to minimize the cost function, and obtain the optimized crop model parameters through multiple iterations ;
S7将S6获得的优化作物模型参数代入作物模型，选择待测农作物纯度大于60%的像元，逐像元单元运行作物模型模拟产量，然后汇总到县域行政单元，输出县域待测农作物单产，指导粮食生产。S7 substitutes the parameters of the optimized crop model obtained in S6 into the crop model, selects pixels with a purity greater than 60% of the crops to be tested, runs the crop model to simulate yields pixel by pixel, and then aggregates them to county administrative units to output the unit yield of the crops to be tested in the county. food production.
其中，S1所述气象参数为最高气温、最低气温、总辐射量、水汽压、风速、降水等参数。Wherein, the meteorological parameters mentioned in S1 are parameters such as maximum temperature, minimum temperature, total radiation, water vapor pressure, wind speed, and precipitation.
其中，S1所述作物参数为农作物分布参数、农作物的物候特征等。Wherein, the crop parameters described in S1 are crop distribution parameters, phenological characteristics of crops, and the like.
其中，S1所述基于气象参数进行区域化，完成作物模型的空间化是指结合遥感影像对收集的参数进行空间位置的匹配，对气象参数和作物模型所需的积温参数，采用Kriging插值算法，完成参数区域化标定。Among them, the regionalization based on the meteorological parameters mentioned in S1, and the completion of the spatialization of the crop model refers to the matching of the spatial position of the collected parameters in combination with the remote sensing image, and the Kriging interpolation algorithm is used for the meteorological parameters and the accumulated temperature parameters required by the crop model, Complete the parameter regionalization calibration.
其中，S2所述多时相的遥感影像为多时相的Landsat TM。Wherein, the multi-temporal remote sensing images mentioned in S2 are multi-temporal Landsat TM.
其中，S2所述获取作物类型分布图和待测农作物纯度百分比图为采用C5.0的决策树分类算法，结合农作物的物候特征和光谱特征，构建分类规则，获取作物类型分布图，采用1km格网套合，生成1km的待测农作物纯度百分比空间分布图，为了减少非待测农作物区域对数据同化结果的影响，选择待测农作物纯度大于60%的像元，作为后续县域尺度产量估测的像元。Among them, the acquisition of the crop type distribution map and the purity percentage map of the crops to be tested described in S2 is a decision tree classification algorithm using C5.0, combining the phenological characteristics and spectral characteristics of the crops, constructing classification rules, and obtaining the crop type distribution map, using a 1km grid The grid is combined to generate a 1km spatial distribution map of the purity percentage of the crops to be tested. In order to reduce the impact of the non-crops to be tested on the data assimilation results, the pixels with a purity of more than 60% of the crops to be tested are selected as the subsequent county-scale yield estimation. pixel.
其中，S3所述进行滤波，采用上包络线（Savitzky-Golay，S-G）滤波算法，其表达式见公式（1）：Among them, the filtering described in S3 adopts the upper envelope (Savitzky-Golay, S-G) filtering algorithm, and its expression is shown in formula (1):
式中：Yj+i表示原始LAI曲线上的一块窗口内的值，m为窗口的半径，N为卷积数目，窗口的宽度为2m+1；Yj*表示滤波后的LAI值；C表示第i个LAI值的滤波系数。In the formula: Yj+i represents the value in a window on the original LAI curve, m is the radius of the window, N is the number of convolutions, and the width of the window is 2m+1; Yj* represents the filtered LAI value; C Indicates the filter coefficient for the ith LAI value.
其中，S3所述进行滤波，在对S-G滤波改进的基础上，形成上包络线滤波方法，通过LAI时序变化趋势将LAI值分为真值和假值，用局部循环迭代的方式将S-G滤波值代替假值，与真值合成新的LAI平滑曲线，逐步拟合形成LAI时序数据的上包络线。Among them, the filtering described in S3 is based on the improvement of S-G filtering, and the upper envelope filtering method is formed. The LAI value is divided into true value and false value through the LAI time series change trend, and the S-G filtering is performed in a local loop iteration method. The value replaces the false value, synthesizes a new LAI smooth curve with the true value, and gradually fits to form the upper envelope of the LAI time series data.
其中，S4所述构建二级尺度转换模型由2个关键步骤构成：Among them, the construction of the secondary scale conversion model described in S4 consists of two key steps:
1）从S3生成区域的30米待测农作物TM LAI中选择3个关键生育期的TM LAI；1) Select 3 key growth period TM LAIs from the 30-meter crop TM LAIs in the S3 generation area;
2）利用已选的关键生育期TM LAI与对应时期的S-G MODIS LAI之间的比值系数，调整其他关键生育期S-G MODIS LAI，在待测农作物生育期的上升阶段和下降阶段，分别利用Logistic曲线进行拟合，2) Use the ratio coefficient between the selected key growth period TM LAI and the corresponding period S-G MODIS LAI to adjust the other key growth period S-G MODIS LAI, and use the Logistic curve respectively in the rising and falling stages of the growth period of the crops to be measured to fit,
Logistic曲线方程见公式（2）：Logistic curve equation see formula (2):
其中，t是MODIS LAI时间序列的索引，y(t)是t时间对应的LAI值，a和b是拟合参数，c+d是最大LAI值，d为LAI初始值，即LAI时间序列中的第一个值。Among them, t is the index of the MODIS LAI time series, y(t) is the LAI value corresponding to time t, a and b are the fitting parameters, c+d is the maximum LAI value, d is the initial value of LAI, that is, in the LAI time series the first value of .
其中，S5所述待测农作物关键物候期的遥感观测误差是基于实测LAI与尺度修正LAI之间的标准方差计算的；S5所述作物模型误差来源于优化的2个模型参数，TDWI和SPAN。Among them, the remote sensing observation error of the key phenological period of the crop to be measured in S5 is calculated based on the standard deviation between the measured LAI and the scale-corrected LAI; the crop model error in S5 comes from the optimized two model parameters, TDWI and SPAN.
其中，S6所述使用最优化算法最小化代价函数采用SCE_UA算法对四维变分代价函数进行最小化，四维变分代价函数按公式（3）计算：Among them, the optimization algorithm described in S6 minimizes the cost function using the SCE_UA algorithm to minimize the four-dimensional variational cost function, and the four-dimensional variational cost function is calculated according to formula (3):
其中，k表示代价函数中优化模型参数的个数，xk代表优化的WOFOST参数在区间范围的数值，xk0代表优化的WOFOST参数个数的经验值，B表示模型误差，以2个优化模型参数的误差表示，T表示矩阵的转置，N代表时间序列遥感观测数据的次数，yi表示遥感观测LAI，Hi(x)表示WOFOST模型模拟LAI，Q0代表遥感LAI的误差。Among them, k represents the number of optimized model parameters in the cost function, xk represents the value of the optimized WOFOST parameters in the interval range, xk0 represents the empirical value of the number of optimized WOFOST parameters, and B represents the model error. The parameter error representation, T represents the transposition of the matrix, N represents the number of time series remote sensing observation data, yi represents the remote sensing observation LAI, Hi (x) represents the WOFOST model simulation LAI, Q0 represents the error of the remote sensing LAI.
其中，S7具体步骤为：选择待测农作物纯度大于60%的像元通过不断迭代，重新初始化2个作物模型参数，使得作物模型输出的LAI和产量不断发生变化，然后在代价函数中对比MODIS LAI和作物模型输出LAI的最小差异下，当以下三个收敛条件满足其一即可结束同化。获得优化参数的数值：Among them, the specific steps of S7 are: select the pixel with a purity greater than 60% of the crop to be tested, and reinitialize the two crop model parameters through continuous iteration, so that the LAI and yield output by the crop model are constantly changing, and then compare the MODIS LAI in the cost function Under the minimum difference with crop model output LAI, the assimilation can be ended when one of the following three convergence conditions is met. Get the numerical value of the optimization parameter:
①连续5次循环后待优化参数值已收缩到指定的值域范围；① After 5 consecutive cycles, the value of the parameter to be optimized has shrunk to the specified value range;
②目标函数值在5次循环后无法提高0.0001%；②The objective function value cannot be increased by 0.0001% after 5 cycles;
③计算代价函数的次数超过10000次；③ Calculate the cost function more than 10,000 times;
将优化后的模型参数代入作物模型运行，获得产量，结合像元的待测农作物分布百分比，按县域行政边界汇总，输出区域上的单产结果。Substitute the optimized model parameters into the crop model to run to obtain the yield, combined with the distribution percentage of the crops to be tested in the pixel, summarize according to the administrative boundary of the county, and output the yield per unit area in the area.
其中，所述待测农作物优选为冬小麦。Wherein, the crop to be tested is preferably winter wheat.
本发明还提供所述基于尺度转换和数据同化的农作物产量估测方法在指导农作物生产中的应用。The invention also provides the application of the crop yield estimation method based on scale conversion and data assimilation in guiding crop production.
本发明与现有技术相比，有益效果为：Compared with the prior art, the present invention has beneficial effects as follows:
本发明克服了遥感观测像元与作物模型模拟单元之间尺度不匹配的问题，提高了数据同化模型的精度，适合于县域尺度的冬小麦产量估测。The invention overcomes the problem of scale mismatch between the remote sensing observation pixel and the crop model simulation unit, improves the accuracy of the data assimilation model, and is suitable for winter wheat yield estimation at the county level.
具体实施方式detailed description
下面结合实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。Below in conjunction with the examples, the specific implementation of the present invention will be further described in detail. The following examples are used to illustrate the present invention, but are not intended to limit the scope of the present invention.
实施例1Example 1
以河北保定和衡水地区的冬小麦产量估测为例进一步阐述本发明的技术方案。本实施例的基于尺度转换和数据同化的区域冬小麦产量估测方法的流程包括：The technical scheme of the present invention is further described by taking the estimation of winter wheat yield in Baoding and Hengshui, Hebei as examples. The flow of the regional winter wheat yield estimation method based on scale conversion and data assimilation in this embodiment includes:
S1收集研究区内的气象参数、作物参数、土壤参数和管理参数并以之进行WOFOST作物模型的单点尺度标定；之后基于气象参数进行区域化，完成作物模型的空间化。具体如下：S1 collects the meteorological parameters, crop parameters, soil parameters and management parameters in the study area and uses them to perform single-point scale calibration of the WOFOST crop model; then regionalizes based on the meteorological parameters to complete the spatialization of the crop model. details as follows:
选择河北保定和衡水地区为研究区域，该区为河北省重要的冬小麦产区。获取以下数据：选取6个气象站点的最高/最低气温、总辐射、水汽压、风速、降水模型所需的6个气象要素；研究区内农业气象试验站采集的土壤参数和作物参数等；2009年河北省分县冬小麦产量统计数据；MODIS15A2数据产品，统一所有数据的坐标投影为Albers Conical EqualArea，空间参考为WGS1984。The Baoding and Hengshui areas in Hebei Province were selected as the research area, which is an important winter wheat production area in Hebei Province. Obtain the following data: select the six meteorological elements required by the maximum/minimum temperature, total radiation, water vapor pressure, wind speed, and precipitation model of six meteorological stations; soil parameters and crop parameters collected by agricultural meteorological test stations in the study area; 2009 Statistics of winter wheat production in counties of Hebei Province in 2019; MODIS15A2 data product, the coordinate projection of all data is Albers Conical EqualArea, and the spatial reference is WGS1984.
对于WOFOST模型中不敏感的作物参数，或敏感程度较高但是变化范围较小的参数，通过查阅文献或者直接采用WOFOST模型的默认值来确定（如发育最适光长DLO，出苗最低温度下限TBASE，叶的维持呼吸最用RML等）；对于敏感程度高、变化范围大的作物参数，能够通过查阅本研究区的相关文献或者利用实测数据计算获得的（如积温TSUM1，TSUM2，叶同化物转化率CVL，比叶面积SLATB等）。也可以研究区的先验知识等确定参数可能的取值范围和数值（如叶片在35℃时的生命期SPAN，最大CO2同化速率AMAXTB，单叶片同化CO2的光能利用率EFFTB等）。For crop parameters that are not sensitive in the WOFOST model, or parameters that are highly sensitive but have a small range of variation, they can be determined by consulting the literature or directly using the default values of the WOFOST model (such as the optimal light length for development DLO, the lowest temperature limit for emergence TBASE , leaf maintenance respiration is most used RML, etc.); for crop parameters with high sensitivity and wide range of variation, they can be obtained by consulting relevant literature in this research area or using measured data (such as accumulated temperature TSUM1, TSUM2, leaf assimilate conversion rate CVL, specific leaf area SLATB, etc.). The prior knowledge of the study area can also be used to determine the possible value range and value of the parameters (such as the life span SPAN of the leaf at 35°C, the maximum CO2 assimilation rate AMAXTB, the light energy utilization rate of single leaf CO2 assimilation EFFTB, etc.).
利用研究区6个气象站点的观测数据，结合克里金空间插值方法对2008～2009年逐日最高气温、最低气温、总辐射量、水汽压、风速、降水量六个要素进行插值，并将其转换为WOFOST模型的标准格式，得到研究区每个格点相对应的气象数据。对于空间范围内容易获取的参数，例如播种期到出苗期的积温、出苗期到开花期的积温、开花期到成熟期的积温等，利用气象观测台站获得的日平均气温与冬小麦的生物学零度相减得到每日的积温，通过累加不同物候阶段的每日积温得到以上三个积温变量，同样采用与气象数据插值相同的克里金插值法，获得整个研究区内空间变化的积温参数。Using the observation data of 6 meteorological stations in the study area, combined with the Kriging spatial interpolation method, the six elements of daily maximum temperature, minimum temperature, total radiation, water vapor pressure, wind speed and precipitation from 2008 to 2009 were interpolated, and their Convert to the standard format of the WOFOST model to obtain the meteorological data corresponding to each grid point in the study area. For parameters that are easy to obtain in a spatial range, such as accumulated temperature from sowing to emergence, accumulated temperature from emergence to flowering, and accumulated temperature from flowering to maturity, the daily average temperature obtained from meteorological observation stations and the biological parameters of winter wheat The daily accumulated temperature was obtained by subtracting zero degrees, and the above three accumulated temperature variables were obtained by accumulating the daily accumulated temperature of different phenological stages. The same kriging interpolation method as the interpolation of meteorological data was also used to obtain the accumulated temperature parameters of the spatial variation in the entire study area.
S2利用研究区多时相的美国专题制图仪（Thematic Mapper，TM）影像，结合地面调查获得的作物类型样本点，采用C5.0的决策树分类算法，结合农作物的物候特征，构建分类规则，获取作物类型分布图，与1km格网套合，生成1km的冬小麦纯度百分比空间分布图，为了减少非冬小麦区域对数据同化结果的影响，选择冬小麦纯度大于60%的像元，作为后续县域尺度产量估测的基础。S2 uses the multi-temporal Thematic Mapper (TM) images of the study area, combined with the crop type sample points obtained from the ground survey, uses the C5.0 decision tree classification algorithm, and combines the phenological characteristics of the crops to construct classification rules to obtain The crop type distribution map is combined with the 1km grid to generate a 1km winter wheat purity percentage spatial distribution map. In order to reduce the impact of non-winter wheat areas on the data assimilation results, the pixels with winter wheat purity greater than 60% are selected as subsequent county-scale yield estimates. basis of measurement.
S3对实验区整个生育期内的MODIS LAI数据按时间序列合成，对每个栅格单元生成时间序列曲线。具体通过NASA提供的MRT投影转换工具，将河北地区的数据进行镶嵌、投影转换和格式转换，选取1到177天（儒略日），共23个时相影像进行叠加，生成基于栅格单元的时间序列曲线。时间序列MODIS LAI曲线进行滤波，以消除云污染造成的数据缺失。具体使用上包络线（Savitzky-Golay，S-G）滤波解决叶面积指数变化过程不连续的问题，其原理表达如下：S3 synthesized the MODIS LAI data in the whole growth period of the experimental area in time series, and generated time series curves for each grid unit. Specifically, through the MRT projection conversion tool provided by NASA, the data in the Hebei area were mosaiced, projected and format converted, and a total of 23 time-phase images were selected from 1 to 177 days (Julian day) for superimposition to generate a grid-based image. time series curve. The time-series MODIS LAI curves were filtered to remove missing data due to cloud pollution. Specifically, the upper envelope (Savitzky-Golay, S-G) filter is used to solve the problem of discontinuity in the change process of the leaf area index. The principle is expressed as follows:
式中：Yj+i表示原始LAI曲线上的一块窗口内的值，m为窗口的半径，N为卷积数目，窗口的宽度为2m+1；Yj*表示滤波后的LAI值；C表示第i个LAI值的滤波系数。In the formula: Yj+i represents the value in a window on the original LAI curve, m is the radius of the window, N is the number of convolutions, and the width of the window is 2m+1; Yj* represents the filtered LAI value; C Indicates the filter coefficient for the ith LAI value.
在对S-G滤波改进的基础上，形成上包络线滤波方法，通过LAI时序变化趋势将LAI值分为真值和假值，用局部循环迭代的方式将S-G滤波值代替假值，与真值合成新的LAI平滑曲线，逐步拟合形成LAI时序数据的上包络线。具体步骤如下：On the basis of improving the S-G filter, an upper envelope filter method is formed, and the LAI value is divided into true and false values through the LAI time series change trend, and the S-G filter value is replaced by the local loop iteration method. Synthesize a new LAI smooth curve, and gradually fit it to form the upper envelope of the LAI time series data. Specific steps are as follows:
（1）对有云污染的LAI值点进行线性插值。由于云覆盖或是大气效应的影响，会使一些像元值带有较大的噪声，对这些像元采用其周边可信像元值进行空间线性插值，生成初始的LAI数据。(1) Perform linear interpolation on the LAI value points with cloud pollution. Due to the influence of cloud coverage or atmospheric effects, some pixel values will have large noise, and the surrounding trusted pixel values are used to perform spatial linear interpolation on these pixels to generate initial LAI data.
（2）用S-G滤波进行LAI时序变化趋势的提取。根据假设，LAI时序数据应遵循植被物候动态变化的渐进特征，在LAI时序中突变被认为是由于云或大气产生的噪音造成的。因此，对研究区内每个冬小麦像元应用S-G滤波进行平滑，得到LAI时间序列的变化趋势数据LAItr。(2) Use S-G filter to extract the time series variation trend of LAI. According to the hypothesis, the LAI time series data should follow the gradual characteristics of vegetation phenology dynamic changes, and the sudden change in the LAI time series is considered to be caused by the noise generated by clouds or the atmosphere. Therefore, the S-G filter is applied to each winter wheat pixel in the study area for smoothing, and the change trend data LAItr of the LAI time series is obtained.
（3）确定LAI时间序列中每个值点的权重Wi，为后续的重构结果提供判断依据。Wi的计算公式为：(3) Determine the weight Wi of each value point in the LAI time series to provide a judgment basis for subsequent reconstruction results. The formula for calculating Wi is:
其中，dmax为di的最大值。in, dmax is the maximum value of di.
（4）生成新的LAI时间序列。由以上理论可知，含有噪声的LAI值要比时间变化趋势曲线上的LAItr值小，对这些点进行简单的替换生成新的LAI1时序曲线。(4) Generate a new LAI time series. It can be seen from the above theory that the LAI value containing noise is smaller than the LAItr value on the time trend curve, and these points are simply replaced to generate a new LAI1 time series curve.
（5）对LAI1进行S-G滤波，再将滤波后的数据与空间插值后的LAI值进行比较。同样，如果空间插值后的LAI值小于新生成的LAI j，用LAI j代替，生成新序列LAI j+1，迭代进行滤波重构，滤波效果由Fk的大小来评定：(5) Perform SG filtering on LAI1, and then compare the filtered data with the LAI value after spatial interpolation. Similarly, if the LAI value after spatial interpolation is smaller than the newly generated LAI j, replace it with LAI j to generate a new sequence LAI j+1, perform filtering reconstruction iteratively, and the filtering effect is evaluated by the size of Fk :
LAIik为第i个像元在第k次滤波重构后的LAI值，Wi为步骤（3）计算得到的权重，Fk为第k次滤波后的迭代终止判断参数。退出迭代的条件为Fk-1≥Fk≤Fk+1。LAIik is the LAI value of the i-th pixel after the k-th filtering reconstruction, Wi is the weight calculated in step (3), and Fk is the iteration termination judgment parameter after the k-th filtering. The condition for exiting the iteration is Fk-1 ≥Fk ≤Fk+1 .
获取了实验区3景无云的TM数据，时相分别为：2009-3-14、2009-5-17、2009-6-10，与对应的实测LAI分别构建回归统计模型。The TM data of Jing Wuyun in the experimental area 3 were obtained, and the time phases were: 2009-3-14, 2009-5-17, 2009-6-10, and the regression statistical models were constructed with the corresponding measured LAI.
2009年3月14日：由于3月份中旬冬小麦处于返青阶段，裸土对植被指数的影响较大，因此选择SAVI植被指数来建立LAI与TM VI之间的统计模型为：March 14, 2009: Since winter wheat is in the greening stage in mid-March, the bare soil has a greater impact on the vegetation index, so the SAVI vegetation index is selected to establish the statistical model between LAI and TM VI as follows:
决定系数R2=0.849。The coefficient of determination R2 =0.849.
2009年5月17日，建立实测LAI与TM NDVI之间的统计模型为：On May 17, 2009, the statistical model between the measured LAI and TM NDVI was established as follows:
决定系数R2=0.742。The coefficient of determination R2 =0.742.
2009年6月10日，建立实测LAI与TM NDVI之间的统计模型为:On June 10, 2009, the statistical model between the measured LAI and TM NDVI was established as:
决定系数R2=0.874。The coefficient of determination R2 =0.874.
基于上述的3个回归统计模型，获得研究区30米冬小麦的空间分布图。Based on the above three regression statistical models, the spatial distribution map of 30m winter wheat in the study area was obtained.
S4融合地面实测样本点LAI、多时相TM LAI和滤波后的MODISLAI信息，构建二级尺度转换模型，生成时间序列尺度调整LAI。S4 integrates the ground-measured sample point LAI, multi-temporal TM LAI and filtered MODISLAI information to construct a secondary scale conversion model and generate time series scale-adjusted LAI.
二级尺度转换模型包括2个步骤，首先，按S3步骤中生成的3个物候期的实测LAI和TM VI回归统计模型，生成区域30米冬小麦TMLAI；第二步，利用已有的关键生育期TM LAI与对应时期的S-G MODIS LAI之间的比值系数，调整其他关键生育期S-G MODIS LAI，在冬小麦生育期的上升阶段（返青到拔节期），获得4个物候期的调整LAI（孕穗期，拔节期，起身期，返青期）；在冬小麦的下降阶段（拔节到成熟期），获得4个物候期的调整LAI（孕穗期，抽穗期，开花期，成熟期），分别利用Logistic曲线进行拟合。The second-level scale conversion model includes two steps. First, according to the three phenological phases generated in step S3, the LAI and TM VI regression statistical models are used to generate the regional 30-meter winter wheat TMLAI; the second step is to use the existing key growth period The ratio coefficient between TM LAI and S-G MODIS LAI in the corresponding period was adjusted to adjust the S-G MODIS LAI in other key growth stages, and the adjusted LAI of the four phenological stages (booting stage, jointing stage, rising stage, greening stage); in the descending stage of winter wheat (jointing to maturity stage), the adjusted LAI of four phenological stages (booting stage, heading stage, flowering stage, maturity stage) were obtained, and the Logistic curves were used to simulate combine.
Logistic曲线方程如下：The Logistic curve equation is as follows:
其中，t是MODIS LAI时间序列的索引，y(t)是t时间对应的LAI值，a和b是拟合参数，c+d是最大LAI值，d为LAI初始值，即LAI时间序列中的第一个值。Among them, t is the index of the MODIS LAI time series, y(t) is the LAI value corresponding to time t, a and b are the fitting parameters, c+d is the maximum LAI value, d is the initial value of LAI, that is, in the LAI time series the first value of .
S5基于地面样本点关键物候期的LAI，计算遥感观测和作物模型模拟的标准差，获得冬小麦关键物候期的遥感观测误差和作物模型误差。S5 calculates the standard deviation of remote sensing observations and crop model simulations based on the LAI of key phenological periods of ground sample points, and obtains the errors of remote sensing observations and crop model errors of key phenological periods of winter wheat.
模型误差认为是来源于优化的2个参数（TDWI和SPAN）,经过分析计算获得尺度修正LAI在7个物候期的观测误差为0.34；模型误差设置为0.52。The model error is considered to come from the optimized two parameters (TDWI and SPAN). After analysis and calculation, the observation error of the scale-corrected LAI in the seven phenological periods is 0.34; the model error is set to 0.52.
S6以S4获得的尺度调整LAI和作物模型模拟LAI，并引入遥感观测误差和作物模型误差，构建四维变分代价函数，使用最优化算法最小化代价函数，通过多次迭代，获得优化作物模型参数。S6 uses the scale obtained in S4 to adjust the LAI and the crop model to simulate the LAI, and introduces the remote sensing observation error and the crop model error to construct a four-dimensional variational cost function, use the optimization algorithm to minimize the cost function, and obtain the optimized crop model parameters through multiple iterations .
在S1进行作物模型标定的基础上，在实验区域内，运行WOFOST作物模型模拟LAI，结合S4中生成的尺度修正LAI，引入S5中的遥感观测与模型模拟误差,构建四维变分代价函数，如下：On the basis of crop model calibration in S1, run the WOFOST crop model to simulate LAI in the experimental area, combine the scale correction LAI generated in S4, introduce the remote sensing observation and model simulation errors in S5, and construct a four-dimensional variational cost function, as follows :
其中，k表示代价函数中优化模型参数的个数，xk代表优化的WOFOST参数在区间范围的数值，xk0代表优化的WOFOST参数个数的经验值，B表示模型误差，以2个优化模型参数的误差表示，T表示矩阵的转置，N代表时间序列遥感观测数据的次数，yi表示遥感观测LAI，Hi(x)表示WOFOST模型模拟LAI，Q0代表遥感LAI的误差。Among them, k represents the number of optimized model parameters in the cost function, xk represents the value of the optimized WOFOST parameters in the interval range, xk0 represents the empirical value of the number of optimized WOFOST parameters, and B represents the model error. The parameter error representation, T represents the transposition of the matrix, N represents the number of time series remote sensing observation data, yi represents the remote sensing observation LAI, Hi (x) represents the WOFOST model simulation LAI, Q0 represents the error of the remote sensing LAI.
选择冬小麦纯度大于60%的像元，运用优化算法，不断迭代重新初始化2个WOFOST模型参数，使得模型输出的LAI和产量发生变化，然后在代价函数中对比MODIS LAI和WOFOST模型输出LAI的最小差异下，获得最优化参数。SCE_UA算法用来在初始参数空间中搜索全局最优解，使得代价函数快速收敛。当以下三个收敛条件满足其一即可结束同化，获得优化参数的数值。Select the pixel whose winter wheat purity is greater than 60%, use the optimization algorithm to iteratively re-initialize the two WOFOST model parameters, so that the LAI and yield output by the model change, and then compare the minimum difference between the MODIS LAI and the WOFOST model output LAI in the cost function Next, the optimal parameters are obtained. The SCE_UA algorithm is used to search for the global optimal solution in the initial parameter space, so that the cost function converges quickly. When one of the following three convergence conditions is satisfied, the assimilation can be ended and the value of the optimized parameter can be obtained.
（1）连续5次循环后待优化参数值已收缩到指定的值域范围；(1) After 5 consecutive cycles, the value of the parameter to be optimized has shrunk to the specified value range;
（2）目标函数值在5次循环后无法提高0.0001%；(2) The objective function value cannot be increased by 0.0001% after 5 cycles;
（3）计算代价函数的次数超过10000次。(3) Calculate the cost function more than 10,000 times.
S7将S6获得的优化作物模型参数代入作物模型，选择冬小麦纯度大于60%的像元，逐像元单元运行作物模型模拟产量，然后汇总到县域行政单元，输出县域冬小麦单产。将优化后的TDWI和SPAN参数代入WOFOST模型运行，获得产量。结合像元的冬小麦百分比，按县域行政边界汇总，输出区域上的冬小麦单产结果。S7 substitutes the optimized crop model parameters obtained in S6 into the crop model, selects pixels with a winter wheat purity greater than 60%, runs the crop model to simulate yields pixel by pixel, and then aggregates them to the county administrative units to output the county winter wheat yield. Substitute the optimized TDWI and SPAN parameters into the WOFOST model to run to obtain the yield. Combined with the winter wheat percentage of the cell, aggregated by county administrative boundaries, and output the winter wheat yield results for the region.
本发明通过利用S-G滤波算法对冬小麦生育期时间序列MODISLAI曲线进行滤波，消除了数据缺失和云污染的影响，生成了平滑的LAI时间序列曲线，为二级尺度转换模型提供了高精度的数据源。然后通过二级尺度转换模型生成了时间序列尺度调整LAI，有效融合了MODIS LAI时间序列的趋势信息及TM LAI信息，提高了LAI的精度，为同化过程提供了高精度的数据源。最后，通过同化调整后的LAI与作物模型模拟的LAI，提高了区域上冬小麦产量估测的精度。The invention filters the time series MODISLAI curve of winter wheat growth period by using the S-G filter algorithm, eliminates the influence of data loss and cloud pollution, generates a smooth LAI time series curve, and provides a high-precision data source for the secondary scale conversion model . Then, the time series scale adjustment LAI is generated through the secondary scale conversion model, which effectively integrates the trend information of MODIS LAI time series and TM LAI information, improves the accuracy of LAI, and provides a high-precision data source for the assimilation process. Finally, by assimilating the adjusted LAI with the LAI simulated by the crop model, the accuracy of winter wheat yield estimation in the region was improved.
粮食的生产者与消费者都需要及时准确地了解粮食产量信息，根据本方法，可以在冬小麦成熟期大面积地获得产量数据，为国家有关部门进行粮情判断、粮食调控等科学决策等提供重要的科学依据，并且可以作为粮食贸易的重要依据。Both grain producers and consumers need timely and accurate knowledge of grain production information. According to this method, production data can be obtained in a large area during the mature period of winter wheat, which provides important information for relevant state departments to make scientific decisions such as grain situation judgment and grain regulation. scientific basis, and can be used as an important basis for grain trade.
本发明的方法也可以用于其他农作物产量的估测。The method of the present invention can also be used to estimate the yield of other crops.
虽然，上文中已经用一般性说明及具体实施方案对本发明作了详尽的描述，但在本发明基础上，可以对之作一些修改或改进，这对本领域技术人员而言是显而易见的。因此，在不偏离本发明精神的基础上所做的这些修改或改进，均属于本发明要求保护的范围。Although the present invention has been described in detail with general descriptions and specific embodiments above, it is obvious to those skilled in the art that some modifications or improvements can be made on the basis of the present invention. Therefore, the modifications or improvements made on the basis of not departing from the spirit of the present invention all belong to the protection scope of the present invention.
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410155913.0ACN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410155913.0ACN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
| Publication Number | Publication Date |
|---|---|
| CN104134095A CN104134095A (en) | 2014-11-05 |
| CN104134095Btrue CN104134095B (en) | 2017-02-15 |
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201410155913.0AActiveCN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
| Country | Link |
|---|---|
| CN (1) | CN104134095B (en) |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106780323B (en)* | 2016-11-23 | 2020-03-17 | 湖北大学 | Agricultural condition acquisition and real-time updating method and system based on smart phone |
| CN107122739B (en)* | 2017-01-23 | 2020-07-17 | 东北农业大学 | Crop yield estimation model for reconstructing VI time series curve based on Extreme mathematical model |
| CN106951979A (en)* | 2017-02-20 | 2017-07-14 | 中国农业大学 | Crop Maturity Prediction Method Based on Remote Sensing, Crop Model and Meteorological Forecast |
| CN107016466B (en)* | 2017-04-11 | 2020-12-22 | 北京林业大学 | A method and system for optimizing the construction of landscape pattern |
| CN107274036B (en)* | 2017-07-25 | 2020-04-03 | 中国农业科学院农业信息研究所 | A method and system for predicting crop yield |
| CN108509836B (en)* | 2018-01-29 | 2021-10-08 | 中国农业大学 | Crop Yield Estimation Method Based on Data Assimilation of Dual Polarization Synthetic Aperture Radar and Crop Model |
| CN108304973A (en)* | 2018-02-11 | 2018-07-20 | 中国农业大学 | Area crops maturity period prediction technique based on accumulated temperature, radiation and soil moisture content |
| CN108133298B (en)* | 2018-03-08 | 2022-04-19 | 河南工业大学 | A Multiple Regression Model National Grain Consumption Forecast Method |
| CN108982369B (en)* | 2018-04-28 | 2020-09-01 | 中国农业大学 | Plot scale crop growth monitoring method integrating GF-1WFV and MODIS data |
| CN109033539A (en)* | 2018-07-02 | 2018-12-18 | 河北省科学院地理科学研究所 | The calculation method influenced based on plant growth mechanism model separation key factor pair crop phenology |
| CN108921351A (en)* | 2018-07-06 | 2018-11-30 | 北京兴农丰华科技有限公司 | Crop production forecast method based on trend yield and Meteorological Output |
| CN108984803B (en)* | 2018-10-22 | 2020-09-22 | 北京师范大学 | A method and system for spatializing crop yields |
| CN109800921B (en)* | 2019-01-30 | 2019-12-20 | 北京师范大学 | Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization |
| CN110008621B (en)* | 2019-04-15 | 2023-01-10 | 中国农业科学院农业资源与农业区划研究所 | Crop model remote sensing assimilation estimation method based on dual stream dependence set square root filtering assimilation algorithm |
| CN110378521B (en)* | 2019-06-28 | 2022-04-22 | 河南农业大学 | Construction and Application of Prediction Model of Winter Wheat Yield in Northeast Henan Province |
| CN110633841B (en)* | 2019-08-13 | 2022-04-01 | 中国农业大学 | Provincial range plot scale data assimilation yield prediction method based on set sampling |
| CN110766308B (en)* | 2019-10-17 | 2020-07-10 | 中国科学院地理科学与资源研究所 | A regional crop yield estimation method based on set assimilation strategy |
| CN110751094B (en)* | 2019-10-21 | 2020-08-28 | 北京师范大学 | Crop yield estimation method based on GEE comprehensive remote sensing image and deep learning method |
| CN113033262B (en)* | 2019-12-25 | 2022-08-05 | 中移(成都)信息通信科技有限公司 | Model training method and crop yield estimation method |
| US11935242B2 (en) | 2020-03-09 | 2024-03-19 | International Business Machines Corporation | Crop yield estimation |
| CN111915096B (en)* | 2020-08-14 | 2021-03-09 | 中国科学院地理科学与资源研究所 | Crop yield early-stage forecasting technology based on crop model, remote sensing data and climate forecasting information |
| CN112052988B (en)* | 2020-08-18 | 2024-03-22 | 中国农业大学 | Crop yield estimation method coupling multi-objective optimization and collection assimilation and application |
| CN115018105B (en)* | 2021-03-03 | 2025-04-04 | 中国农业大学 | A method and system for forecasting winter wheat meteorological yield |
| CN114254926A (en)* | 2021-12-22 | 2022-03-29 | 北京尚德智汇科技有限公司 | Crop yield remote sensing estimation method based on phenological index |
| CN114528282B (en)* | 2022-01-20 | 2024-07-19 | 重庆邮电大学 | Time-space Kriging improvement-based method for reconstructing MODIS NDVI time sequence in cloudy region |
| CN114529097B (en)* | 2022-02-26 | 2023-01-24 | 黑龙江八一农垦大学 | Multi-scale crop phenological period remote sensing dimensionality reduction prediction method |
| CN114692971B (en)* | 2022-03-30 | 2024-08-20 | 中国农业大学 | Crop yield prediction method and device based on yield difference |
| CN115753625B (en)* | 2022-11-02 | 2024-04-19 | 中国农业大学 | Regional crop yield estimation method, device, electronic device and storage medium |
| CN116432859B (en)* | 2023-05-08 | 2023-10-27 | 中山大学 | Crop yield statistical data downscaling method |
| CN117830860B (en)* | 2024-03-06 | 2024-12-06 | 江苏省基础地理信息中心 | A remote sensing automatic extraction method for winter wheat planting structure |
| CN118313160B (en)* | 2024-05-27 | 2024-11-12 | 中化现代农业有限公司 | A crop growth model optimization method and system based on phenotypic curve |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7702597B2 (en)* | 2004-04-20 | 2010-04-20 | George Mason Intellectual Properties, Inc. | Crop yield prediction using piecewise linear regression with a break point and weather and agricultural parameters |
| CN102323987A (en)* | 2011-09-09 | 2012-01-18 | 北京农业信息技术研究中心 | A method for assimilating crop leaf area index |
| CN102651096A (en)* | 2012-04-28 | 2012-08-29 | 中国农业大学 | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve |
| JP2012203875A (en)* | 2011-03-28 | 2012-10-22 | Tokyo Electric Power Co Inc:The | Yield estimation device and computer program |
| CN103345707A (en)* | 2013-06-04 | 2013-10-09 | 中国科学院遥感与数字地球研究所 | Crop maturation stage remote sensing prediction method based on multi-source remote sensing data |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7702597B2 (en)* | 2004-04-20 | 2010-04-20 | George Mason Intellectual Properties, Inc. | Crop yield prediction using piecewise linear regression with a break point and weather and agricultural parameters |
| JP2012203875A (en)* | 2011-03-28 | 2012-10-22 | Tokyo Electric Power Co Inc:The | Yield estimation device and computer program |
| CN102323987A (en)* | 2011-09-09 | 2012-01-18 | 北京农业信息技术研究中心 | A method for assimilating crop leaf area index |
| CN102651096A (en)* | 2012-04-28 | 2012-08-29 | 中国农业大学 | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve |
| CN103345707A (en)* | 2013-06-04 | 2013-10-09 | 中国科学院遥感与数字地球研究所 | Crop maturation stage remote sensing prediction method based on multi-source remote sensing data |
| Title |
|---|
| Estimating regional winter wheat yield by assimilation of time series of HJ-1 CCD NDVI into WOFOST-ACRM model with Ensemble Kalman Filter;Hongyuan Ma etal.;《Mathematical and Computer Modeling》;20130831;第58卷;第759-770页* |
| 基于遥感信息与作物模型集合卡尔曼滤波同化地区域冬小麦产量预测;黄健熙 等;《农业工程学报》;20120229;第28卷(第4期);第142-148页* |
| 条件植被温度指数的四维变分与集合卡尔曼同化方法;王维 等;《农业工程学报》;20111231;第27卷(第12期);第184-190页* |
| Publication number | Publication date |
|---|---|
| CN104134095A (en) | 2014-11-05 |
| Publication | Publication Date | Title |
|---|---|---|
| CN104134095B (en) | Crop yield estimation method based on scale transformation and data assimilation | |
| CN111104858B (en) | A large-scale crop phenology extraction method based on morphological model method | |
| CN102651096B (en) | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve | |
| Huang et al. | The improved winter wheat yield estimation by assimilating GLASS LAI into a crop growth model with the proposed Bayesian posterior-based ensemble Kalman filter | |
| CN110751094B (en) | Crop yield estimation method based on GEE comprehensive remote sensing image and deep learning method | |
| CN114062439B (en) | Method for jointly estimating salinity of soil profile by using time series remote sensing images | |
| CN108509836B (en) | Crop Yield Estimation Method Based on Data Assimilation of Dual Polarization Synthetic Aperture Radar and Crop Model | |
| CN109344865B (en) | A Data Fusion Method of Multiple Data Sources | |
| CN109800921B (en) | Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization | |
| Dan et al. | Delineating the rice crop activities in Northeast China through regional parametric synthesis using satellite remote sensing time-series data from 2000 to 2015 | |
| CN106951979A (en) | Crop Maturity Prediction Method Based on Remote Sensing, Crop Model and Meteorological Forecast | |
| CN105321120A (en) | Assimilation evapotranspiration and LAI (leaf area index) region soil moisture monitoring method | |
| CN108446999B (en) | Method for estimating ET (energy-efficiency) of different crops in irrigation area based on crown air temperature difference and remote sensing information | |
| CN110633841B (en) | Provincial range plot scale data assimilation yield prediction method based on set sampling | |
| Yan et al. | Assimilating multi-source remotely sensed data into a light use efficiency model for net primary productivity estimation | |
| CN113361191A (en) | Pixel scale winter wheat yield remote sensing estimation method based on multi-scenario simulation | |
| CN105160192A (en) | TRMM (Tropical Rainfall Measuring Mission) satellite rainfall data downscaling method based on M5-Local | |
| CN103955860A (en) | Regional crop yield estimation method based on ensemble Kalman filter assimilation | |
| CN102162850A (en) | Wheat yield remote sensing monitoring and forecasting method based on model | |
| CN105912836A (en) | Pure remote sensing data driven drainage basin water circulation simulation method | |
| CN107423850A (en) | Region corn maturity period Forecasting Methodology based on time series LAI curve integral areas | |
| CN111062526A (en) | Winter wheat yield per unit prediction method and system | |
| CN111611541A (en) | Method and system for deriving precipitation data in no-data area based on Copula function | |
| CN109033543A (en) | A kind of heterogeneous area's vegetation fraction estimation method, device and equipment of earth's surface | |
| Yu et al. | Improvement of medium-range forecasts using the analog-dynamical method |
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant |