PCSE机理剖析（一）

snackpdd

已于 2024-08-29 15:43:49 修改

阅读量606

点赞数 10

分类专栏：笔记文章标签： python

于 2024-08-29 15:38:59 首次发布

本文链接：https://blog.csdn.net/snackpdd/article/details/141643465

版权

笔记专栏收录该内容

1 篇文章 0 订阅

订阅专栏

PCSE机理剖析（一）

生育期算法

1.春化作用
春化作用指低温条件下的一种适应性发育过程，特指冬小麦等作物。其需要在一个最优温度范围内度过一定天数才能完成春化需求，在满足春化需求之前，生长发育将会延迟。

春化速率（VERNR）：由温度响应函数VERNRTB定义。在最优温度范围内，每天会增加一天到春化状态（VERN）。
发育速率减缩：春化的影响通过基础春化需求（VERNBASE）和饱和春化需求（VERNSAT）来计算。发育的减缩因子（VERNFAC）在VERNBASE和VERNSAT之间线性缩放。
关键发育阶段（VERNDVS）：使用一个关键的发育阶段来停止春化的影响，这个阶段被称为VERNDVS。这样做是为了提高模型的稳定性，以避免由于VERNSAT设定过高而永远达不到开花期（Anthesis）。如果发生这种情况，模型会在日志文件中写入警告。
状态变量定义：

    # Helper variable to indicate that DVS > VERNDVS
    _force_vernalisation = Bool(False)

    class Parameters(ParamTemplate):
        VERNSAT = Float(-99.)     # Saturated vernalisation requirements
        VERNBASE = Float(-99.)     # Base vernalisation requirements
        VERNRTB = AfgenTrait()    # Vernalisation temperature response
        VERNDVS = Float(-99.)     # Critical DVS for vernalisation fulfillment

    class RateVariables(RatesTemplate):
        VERNR = Float(-99.)        # Rate of vernalisation
        VERNFAC = Float(-99.)      # Red. factor for phenol. devel.

    class StateVariables(StatesTemplate):
        VERN = Float(-99.)              # Vernalisation state
        DOV = Instance(datetime.date)  # Day when vernalisation
                                            # requirements are fulfilled
        ISVERNALISED =  Bool()              # True when VERNSAT is reached and
                                            # Forced when DVS > VERNDVS

模拟参数

VERNSAT——饱和春化需求；VERNBASE——基础春化需求；VERNRTB——日平均温度函数下的春化速率；VERNDVS——达到此发育阶段后春化效果停止
状态变量

VERN——春化状态；DOV——满足春化需求日期；ISVERNALISED——标志位，表示已达到春化需求
速率变量

VERNR——春化速率；VERNFAC——由于春化效应对发育速率的减缩因子
外部依赖

DVS——发育阶段，仅用于判断是否达到春化的临界发育阶段（VERNDVS）
春化模块初始化

    # Helper variable to indicate that DVS > VERNDVS
    _force_vernalisation = Bool(False)

    class Parameters(ParamTemplate):
        VERNSAT = Float(-99.)     # Saturated vernalisation requirements
        VERNBASE = Float(-99.)     # Base vernalisation requirements
        VERNRTB = AfgenTrait()    # Vernalisation temperature response
        VERNDVS = Float(-99.)     # Critical DVS for vernalisation fulfillment

    class RateVariables(RatesTemplate):
        VERNR = Float(-99.)        # Rate of vernalisation
        VERNFAC = Float(-99.)      # Red. factor for phenol. devel.

    class StateVariables(StatesTemplate):
        VERN = Float(-99.)              # Vernalisation state
        DOV = Instance(datetime.date)  # Day when vernalisation
                                            # requirements are fulfilled
        ISVERNALISED =  Bool()              # True when VERNSAT is reached and
                                            # Forced when DVS > VERNDVS

参数实例化：使用从parvalues提供的参数值来实例化Parameters类。这将设定春化需求的基础和饱和值、温度响应函数和临界发育阶段等。
速率变量初始化：初始化RateVariables类的实例，这些变量将用来计算春化速率（VERNR）和发育速率减缩因子（VERNFAC）。这里指定了VERNFAC作为发布变量，意味着它可以被模型的其他部分访问。
状态变量初始化：设置StateVariables类的初始状态。其中，春化状态（VERN）初始化为0（表示没有春化天数累积），春化减缩因子（VERNFAC）也设为0（无减缩效果），春化完成日期（DOV）为None（未完成），以及春化是否完成的标志（ISVERNALISED）设为False。

计算VERNR和VERNFAC

    def calc_rates(self, day, drv):
        rates = self.rates
        states = self.states
        params = self.params

        DVS = self.kiosk["DVS"]
        if not states.ISVERNALISED:
            if DVS < params.VERNDVS:
                rates.VERNR = params.VERNRTB(drv.TEMP)
                r = (states.VERN - params.VERNBASE)/(params.VERNSAT-params.VERNBASE)
                rates.VERNFAC = limit(0., 1., r)
            else:
                rates.VERNR = 0.
                rates.VERNFAC = 1.0
                self._force_vernalisation = True
        else:
            rates.VERNR = 0.
            rates.VERNFAC = 1.0

春化状态判断：如果作物尚未完成春化（ISVERNALISED为False）：
(1)判断发育阶段：如果当前发育阶段（DVS）小于临界发育阶段（VERNDVS），则继续春化过程。
计算春化速率（VERNR）：根据当日平均温度（drv.TEMP）和温度响应函数（VERNRTB）计算春化速率。
计算减缩因子（VERNFAC）：计算发育减缩因子，基于春化状态（VERN）相对于基础春化需求（VERNBASE）和饱和春化需求（VERNSAT）的比例。该因子在0和1之间线性变化。
(2)达到临界发育阶段：如果DVS达到或超过VERNDVS，停止春化速率的增加，将春化速率（VERNR）设为0，减缩因子（VERNFAC）设为1.0，并设置强制春化标志（_force_vernalisation）为True。
春化状态判断：如果春化已经完成（即ISVERNALISED为True）：
(1)春化速率（VERNR）和减缩因子（VERNFAC）都保持不变，分别为0和1.0

更新春化状态和相关变量

    def integrate(self, day, delt=1.0):
        states = self.states
        rates = self.rates
        params = self.params
        
        states.VERN += rates.VERNR
        
        if states.VERN >= params.VERNSAT:  # Vernalisation requirements reached
            states.ISVERNALISED = True
            if states.DOV is None:
                states.DOV = day
                msg = "Vernalization requirements reached at day %s."
                self.logger.info(msg % day)

        elif self._force_vernalisation:  # Critical DVS for vernalisation reached
            # Force vernalisation, but do not set DOV
            states.ISVERNALISED = True

            # Write log message to warn about forced vernalisation
            msg = ("Critical DVS for vernalization (VERNDVS) reached " +
                   "at day %s, " +
                   "but vernalization requirements not yet fulfilled. " +
                   "Forcing vernalization now (VERN=%f).")
            self.logger.warning(msg % (day, states.VERN))

        else:  # Reduction factor for phenologic development
            states.ISVERNALISED = False

春化需求达成：
(1)如果春化状态(VERN)达到或超过了饱和春化需求(VERNSAT)：
设置春化完成标志(ISVERNALISED)为True。
如果春化完成日期(DOV)尚未设置，记录该日期并记录相应日志。
强制春化：
(1)如果达到临界发育阶段(VERNDVS)而春化需求尚未满足，将执行强制春化：
设置春化完成标志(ISVERNALISED)为True，但不设置春化完成日期(DOV)。
输出警告日志，说明已达到临界发育阶段但春化需求未满足，现在强制执行春化。
未达到春化需求：
(1)如果未达到上述条件，继续正常的春化过程，不改变春化完成标志(ISVERNALISED)。

2.生理发育
详细说明作物如何在环境因素如温度、日长和春化作用的影响下经历不同的生长发育阶段。
0-出苗；1-开花；2-成熟。此定义特别适用于谷类作物，对某些作物（土豆）而言，第一阶段实际上表示块茎开始形成，而非开花。
温度是主导生理发育的关键因素，在开花前，日长和春化效应会对作物生长发育产生影响。开花后，发育速率仅受温度控制。

模拟参数

TSUMEM——从播种到出苗所需的温度总和；TBASEM——出苗的基础温度；TEFFMX——出苗的最大有效温度；TSUM1——从出苗到开花所需的温度总和；TSUM2——从开花到成熟所需的温度总和；IDSL——生理发育模式的开关，用于选择只考虑温度（IDSL=0）、包括日长（IDSL=1）或同时包括日长和春化（IDSL>=2）的情况；DLO——最佳日长，影响生理发育；DLC——关键日长，影响生理发育；出苗时的初始发育阶段，通常为零，但对于一些移栽作物（如水稻）可能更高；DVSEND——最终发育阶段；DTSMTB——根据日平均温度计算的温度总和的日增量
状态变量

DVS——发育阶段；TSUM：温度总和；TSUME——出苗所需的温度总和；DOS——播种日；DOE——出苗日；DOA——开花日；DOM——成熟日；DOH——收获日；STAGE——当前生理阶段，可能的值包括：emerging（出苗期）、vegetative（营养生长期）、reproductive（生殖生长期）、mature（成熟期）
速率变量

DTSUME——出苗阶段的温度总和日增量；DTSUM——开花或成熟阶段的温度总和日增量；DVR——发育速率
生理发育初始化

    def initialize(self, day, kiosk, parvalues):
        """
        :param day: start date of the simulation
        :param kiosk: variable kiosk of this PCSE  instance
        :param parvalues: `ParameterProvider` object providing parameters as
                key/value pairs
        """

        self.params = self.Parameters(parvalues)
        self.rates = self.RateVariables(kiosk)
        self.kiosk = kiosk

        self._connect_signal(self._on_CROP_FINISH, signal=signals.crop_finish)

        # Define initial states
        DVS, DOS, DOE, STAGE = self._get_initial_stage(day)
        self.states = self.StateVariables(kiosk, publish="DVS",
                                          TSUM=0., TSUME=0., DVS=DVS,
                                          DOS=DOS, DOE=DOE, DOA=None, DOM=None,
                                          DOH=None, STAGE=STAGE)

        # initialize vernalisation for IDSL=2
        if self.params.IDSL >= 2:
            self.vernalisation = Vernalisation(day, kiosk, parvalues)

状态变量初始化：
调用 _get_initial_stage 方法确定模拟开始时作物的初始发育阶段，该方法根据开始日和作物的起始类型（种植或出苗）来设置。
初始化 StateVariables 类的实例，设置初始的发育阶段（DVS）、种植日（DOS）、出苗日（DOE）等，并将初始温度总和设置为0。

春化模块初始化：
如果 IDSL 参数设置为2或更高（表示模型中包含春化效应），则初始化 Vernalisation 类的实例，为接下来的模拟提供春化处理能力。

确定初始发育阶段（_get_initial_stage）

    def _get_initial_stage(self, day):
        """"""
        p = self.params

        # Define initial stage type (emergence/sowing) and fill the
        # respective day of sowing/emergence (DOS/DOE)
        if p.CROP_START_TYPE == "emergence":
            STAGE = "vegetative"
            DOE = day
            DOS = None
            DVS = p.DVSI

            # send signal to indicate crop emergence
            self._send_signal(signals.crop_emerged)

        elif p.CROP_START_TYPE == "sowing":
            STAGE = "emerging"
            DOS = day
            DOE = None
            DVS = -0.1

        else:
            msg = "Unknown start type: %s" % p.CROP_START_TYPE
            raise exc.PCSEError(msg)
            
        return DVS, DOS, DOE, STAGE

发育阶段确定：
出苗（emergence）：如果作物的起始类型是出苗，那么模拟从作物已经出苗的状态开始。设定发育阶段为“vegetative”（营养生长期），记录出苗日期为模拟的开始日期，而种植日期（DOS）则不适用。此外，发育阶段的指标（DVS）设置为参数 DVSI（通常为0，除非是特殊作物如移栽作物）。
播种（sowing）：如果起始类型是播种，则模拟从播种阶段开始。设定发育阶段为“emerging”（出苗期），种植日期记录为模拟开始的日期，出苗日期（DOE）则未知，DVS 初始值设为 -0.1，表示作物尚未出苗。

状态发送：
当作物被确定为已经出苗时，会发送一个信号（crop_emerged），通知模型其他部分作物已经进入营养生长期。

发育速率计算

    def calc_rates(self, day, drv):
        """Calculates the rates for phenological development
        """
        p = self.params
        r = self.rates
        s = self.states

        # Day length sensitivity
        DVRED = 1.
        if p.IDSL >= 1:
            DAYLP = daylength(day, drv.LAT)
            DVRED = limit(0., 1., (DAYLP - p.DLC)/(p.DLO - p.DLC))

        # Vernalisation
        VERNFAC = 1.
        if p.IDSL >= 2:
            if s.STAGE == 'vegetative':
                self.vernalisation.calc_rates(day, drv)
                VERNFAC = self.kiosk["VERNFAC"]

        # Development rates
        if s.STAGE == "emerging":
            r.DTSUME = limit(0., (p.TEFFMX - p.TBASEM), (drv.TEMP - p.TBASEM))
            r.DTSUM = 0.
            r.DVR = 0.1 * r.DTSUME/p.TSUMEM

        elif s.STAGE == 'vegetative':
            r.DTSUME = 0.
            r.DTSUM = p.DTSMTB(drv.TEMP) * VERNFAC * DVRED
            r.DVR = r.DTSUM/p.TSUM1

        elif s.STAGE == 'reproductive':
            r.DTSUME = 0.
            r.DTSUM = p.DTSMTB(drv.TEMP)
            r.DVR = r.DTSUM/p.TSUM2
        elif s.STAGE == 'mature':
            r.DTSUME = 0.
            r.DTSUM = 0.
            r.DVR = 0.
        else:  # Problem: no stage defined
            msg = "Unrecognized STAGE defined in phenology submodule: %s"
            raise exc.PCSEError(msg, self.states.STAGE)
        
        msg = "Finished rate calculation for %s"
        self.logger.debug(msg % day)

环境参数的获取和处理：
（1）光周期影响因子

DVRED=limit(0,1, （DAYLP-DLC/DLO-DLC ）)

（2）春化影响
如果参数 IDSL 表明模型包含春化影响（值>=2），并且作物当前处于营养生长阶段，则调用 vernalisation.calc_rates 方法计算春化的速率和影响因子（VERNFAC）。VERNFAC 通过春化模块计算得到，用于调整由温度累积引起的发育速率
（3）温度总和日增量计算（日有效积温）

出苗阶段：DTSUME=limit(0,TEFFMX−TBASEM,drv.TEMP−TBASEM)
营养生长和生殖生长阶段：DTSUM=DTSMTB(drv.TEMP)×VERNFAC×DVRED

（4）发育速率（DVR）
发育速率是基于温度总和的日增量（日有效积温）计算，与不同生长阶段所需的温度总和有关。

出苗阶段（Emerging）： DVR=0.1×(DTSUME/TSUMEM)
营养生长阶段（Vegetative）： DVR=DTSUM/TSUM1
生殖生长阶段（Reproductive）： DVR=DTSUM/TSUM2
成熟阶段（Mature）： DVR=0

发育状态更新

    def integrate(self, day, delt=1.0):
        """Updates the state variable and checks for phenologic stages
        """

        p = self.params
        r = self.rates
        s = self.states

        # Integrate vernalisation module
        if p.IDSL >= 2:
            if s.STAGE == 'vegetative':
                self.vernalisation.integrate(day, delt)
            else:
                self.vernalisation.touch()

        # Integrate phenologic states
        s.TSUME += r.DTSUME
        s.DVS += r.DVR
        s.TSUM += r.DTSUM

        # Check if a new stage is reached
        if s.STAGE == "emerging":
            if s.DVS >= 0.0:
                self._next_stage(day)
                s.DVS = 0.
        elif s.STAGE == 'vegetative':
            if s.DVS >= 1.0:
                self._next_stage(day)
                s.DVS = 1.0
        elif s.STAGE == 'reproductive':
            if s.DVS >= p.DVSEND:
                self._next_stage(day)
                s.DVS = p.DVSEND
        elif s.STAGE == 'mature':
            pass
        else: # Problem no stage defined
            msg = "No STAGE defined in phenology submodule"
            raise exc.PCSEError(msg)
            
        msg = "Finished state integration for %s"
        self.logger.debug(msg % day)

集成春化模块：
如果春化被包括在模型中（p.IDSL >= 2），并且作物处于营养生长阶段（vegetative），则调用 vernalisation.integrate 方法来更新春化状态。
对于非营养生长阶段，调用 vernalisation.touch 方法来更新内部状态，但不改变春化状态。

集成生理状态：
更新出苗所需的温度总和（s.TSUME）。
更新当前发育阶段指标（s.DVS）。
更新整个生长周期的温度总和（s.TSUM）。

检查新的发育阶段：
出苗阶段（emerging）：如果发育阶段指标（DVS）达到或超过0，调用 _next_stage 方法转入下一个阶段，并重置 DVS。
营养生长阶段（vegetative）：如果 DVS 达到或超过1.0，同样调用 _next_stage 转入生殖生长阶段。
生殖生长阶段（reproductive）：如果 DVS 达到或超过最终阶段指标（DVSEND），调用 _next_stage 进入成熟阶段。
成熟阶段（mature）：不进行任何更新，因为作物已完成生长周期。

发育阶段切换（_next_stage）：

    def _next_stage(self, day):
        """Moves states.STAGE to the next phenological stage"""
        s = self.states
        p = self.params

        current_STAGE = s.STAGE
        if s.STAGE == "emerging":
            s.STAGE = "vegetative"
            s.DOE = day
            # send signal to indicate crop emergence
            self._send_signal(signals.crop_emerged)
            
        elif s.STAGE == "vegetative":
            s.STAGE = "reproductive"
            s.DOA = day
                        
        elif s.STAGE == "reproductive":
            s.STAGE = "mature"
            s.DOM = day
            if p.CROP_END_TYPE in ["maturity","earliest"]:
                self._send_signal(signal=signals.crop_finish,
                                  day=day, finish_type="maturity",
                                  crop_delete=True)
        elif s.STAGE == "mature":
            msg = "Cannot move to next phenology stage: maturity already reached!"
            raise exc.PCSEError(msg)

        else: # Problem no stage defined
            msg = "No STAGE defined in phenology submodule."
            raise exc.PCSEError(msg)
        
        msg = "Changed phenological stage '%s' to '%s' on %s"
        self.logger.info(msg % (current_STAGE, s.STAGE, day))

3.总结
作物的营养生长和生殖生长两个阶段的划分参数为有效积温 TSUM1 和TSUM2。物候期的生长依据日发育速率计算（Daily Development Rate），日发育速率受到春化、日照长度和日有效温度的影响。
生长发育速率——营养生长和生殖生长阶段：t 时刻发育速率 = t时刻春化影响因子 * t时刻光周期影响因子 * 每日有效温度 * 完成阶段 i 需要的积温（DTSUM=DTSMTB(drv.TEMP)×VERNFAC×DVRED；DVR=DTSUM/TSUM1； DVR=DTSUM/TSUM2）
生长发育进程——DVS = t时刻DVS + t 时刻发育速率 * 时间步长