Tiller straightening using adelwheatdyn

README.md

@@ -30,7 +30,7 @@ gramineaous plant.
 ## Content
 
 The package hosts generic data structure and simulation tools for
-gramineaous plants(Fournier & Pradal, unpublished), the Adel-Maize
+gramineous plants(Fournier & Pradal, unpublished), the Adel-Maize
 (Fournier & Andrieu, 1998), Adel-Wheat (Fournier et al. 2003) models,
 together with the wheat parameterization model of Abichou et al. (2013)
 and the plastic leaf model of Fournier & Pradal (2012)
@@ -48,6 +48,6 @@ mamba env create -n adel -c openalea3 -c conda-forge openalea.adel
 
 ```bash
 git clone 'https://github.com/openalea/adel.git'
-cd caribu
+cd adel
 mamba env create -n adel_dev -f ./conda/environment.yml
 ```

src/openalea/adel/Adel.R

@@ -393,6 +393,19 @@ kinLvis <- function(kinlist,pars=NULL) {
 # stem inclination is beared by visible sheaths (option 0) or by nodes (option 1)
 # leaf inclination indicates whether leaf base angle is dynamic or not
 #
+#compute incT
+getIncT <- function(axename, HS, incBase, start_incT=1, incT_rate=30) {
+  if (axename == "MS") {
+    incT <- incBase
+  } else {
+    if (HS > start_incT) {
+      incT <- max(3,min(incBase, incT_rate * (HS - start_incT)))
+    } else {
+      incT <- 3
+    }
+  }
+  incT
+}
 #
 #returns stack of visible elements of the stem of an axe
 stemElements <- function(desc) {
@@ -416,29 +429,14 @@ stemElements <- function(desc) {
 #
 # Compute inclinations of stem elements
 #
-axe_inclination <- function(dat, HS, ht, axename, incBase, dredT, start_incT=1, incT_rate=30,epsillon=1e-6) {
+axe_inclination <- function(dat, incT, dredT, epsillon=1e-6) {
   nbphy <- nrow(dat)#inclus ear,ped et awn
       # Calcul des inclinaisons de tiges
       # 1er phyto = entrenoeud a incT
   Einc <- rep(0,nbphy)
   Ginc <- rep(0,nbphy)
-  if (axename == "MS") {
-    incT <- incBase
-  } else {
-    if (HS > start_incT) {
-      incT <- max(3,min(incBase, incT_rate * (HS - start_incT)))
-    } else {
-      incT <- 3
-    }
-  }
   Einc[1] <- incT
-  if (axename != "MS" && incT <= 3) { #Do not represent basal part of first metamer for non inclining tillers
-    dat$Lv[1] =  min(dat$Ll[1],max(0, dat$Ll[1] + dat$Gl[1] + dat$El[1] -  ht))
-    dat$Lr[1] = min(dat$Lv[1],dat$Lr[1])
-    dat$Gv[1] =  min(dat$Gl[1],max(0, dat$Gl[1] + dat$El[1] -  ht))
-    dat$Ev[1] =  min(dat$El[1],max(0, dat$El[1] -  ht))
-  }
-                                        # redressement (if any)
+      # redressement (if any)
   if (dredT > 0 & sum(dat$Ev+dat$Gv) > epsillon) {
           #distance inserton talle -> extremite stemElements
     stem <- stemElements(dat)
@@ -519,7 +517,8 @@ getdesc <- function(kinlist,plantlist,pars=list("senescence_leaf_shrink" = 0.5,"
           axilrank <- ms_pos(axename)
           ht <- kin$MS[t,axilrank+1,"ht"] # length of the tube the axe emerge from
         }
-        dat <- axe_inclination(dat, HS_axe, ht, axename, dataxe$incT, dataxe$dredT, start_incT, incT_rate)
+        incT <- getIncT(axename, HS_axe, dataxe$incT, start_incT, incT_rate)
+        dat <- axe_inclination(dat, incT, dataxe$dredT)
 
         #azimuts : Attention new 21 fev 2011 : azimuts en relatif / phytomere precedent !
         Laz <- datp$Azim
@@ -639,4 +638,5 @@ getAxeT <- function(plants) do.call('rbind', mapply(function(idpl,pl) {df=pl$axe
 getPhenT <- function(plants, axe='MS') do.call('rbind', mapply(function(idpl,pl) {df=pl$pheno[[axe]];df$plant=idpl;df},seq(plants),plants,SIMPLIFY=FALSE))
 #
 getPhytoT <- function(plants, axe='MS') do.call('rbind', mapply(function(idpl,pl) {df=data.frame(pl$phytoT[,,axe]);df$plant=idpl;df$axe=axe;df$n=seq(nrow(df));df},seq(plants),plants,SIMPLIFY=FALSE))
-
+#
+AxeDfList <- function(axeArray) setNames(lapply(seq_len(dim(axeArray)[3]), function(i) {as.data.frame(axeArray[, , i])}),dimnames(axeArray)[[3]])

src/openalea/adel/AdelR.py

@@ -58,6 +58,7 @@ def get_rcode(file_name):
 RgetAxeT = robj.globalEnv["getAxeT"]
 RgetPhenT = robj.globalEnv["getPhenT"]
 RgetPhytoT = robj.globalEnv["getPhytoT"]
+RAxeDfList = robj.globalEnv["AxeDfList"]
 # RgetLeafT = robj.globalEnv['getLeafT']
 
 
@@ -116,6 +117,32 @@ def dataframeAsdict(df):
             )  # r delegator is replaced by rx/rx2 in new rpy2
     return d
 
+def axeDfdict(axe_array):
+    """convert a [,,axe] Rarray to python axe list of dict"""
+    if r["is.null"](axe_array)[0]:
+        return None
+    return RlistAsDict(RAxeDfList(axe_array))
+
+
+def setAdelpars(Rpars):
+    """convert output of setAdel R object to python readable object"""
+    Rplants = RlistAsDict(Rpars)
+    pyplants = {}
+    for p in Rplants:
+        rout = RlistAsDict(Rplants[p])
+        pyout = {}
+        pyout['refp'] = rout['refp'][0]
+        pyout['axeT'] = pandas.DataFrame(dataframeAsdict(rout['axeT']))
+        _phytoT = axeDfdict(rout['phytoT'])
+        pyout['phytoT'] = {axe: pandas.DataFrame(dataframeAsdict(df)) for axe, df in _phytoT.items()}
+        _pheno = RlistAsDict(rout['pheno'])
+        pyout['pheno'] = {axe: pandas.DataFrame(dataframeAsdict(df)) for axe, df in _pheno.items()}
+        _pedT = RlistAsDict(rout['pedT'])
+        pyout['pedT'] = {axe: pandas.DataFrame(dataframeAsdict(df)) for axe, df in _pedT.items()}
+        if 'ssisenT' in rout:
+            pyout['ssisenT'] = pandas.DataFrame(dataframeAsdict(rout['ssisenT']))
+        pyplants[p] = pyout
+    return pyplants
 
 def _is_iterable(x):
     try:

src/openalea/adel/adel.py

@@ -22,7 +22,7 @@
     plot_statistics,
     midrib_statistics,
 )
-from openalea.adel.newmtg import exposed_areas, exposed_areas2canS, duplicate, mtg_factory
+from openalea.adel.newmtg import exposed_areas, exposed_areas2canS, stem_elements, stem_bases, duplicate, mtg_factory
 
 
 def flat_list(nested_list):
@@ -342,6 +342,25 @@ def get_exposed_areas(self, g, convert=False, TT=None):
         areas["TT"] = TT
         return areas
 
+    def get_stem_bases(self, g):
+        return stem_bases(g)
+
+    def get_stem_elements(self, g):
+        bases = self.get_stem_bases(g)
+        desc = self.get_exposed_areas(g)
+        desc = pandas.concat([desc, bases]).drop_duplicates().sort_values('vid')
+
+        result = []
+
+        for (pid, plant, axe), axe_desc in desc.groupby(["refplant_id", "plant", "axe"]):
+            stem = stem_elements(axe_desc, axe=axe)
+            stem["axe"] = axe
+            stem["plant"] = plant
+            stem["refplant_id"] = str(int(pid))
+            result.append(stem)
+
+        return pandas.concat(result, ignore_index=True)
+
     def axis_statistics(self, g):
         meta = self.meta_informations(g)
         df_lai = self.get_exposed_areas(g, convert=True)

src/openalea/adel/adelwheat_dynamic.py

@@ -106,25 +106,29 @@ def build_mtg(self, parameters, stand, **kwds):
         return g
 
     def update_geometry(
-        self, g, SI_units=False, properties_to_convert={"lengths": [], "areas": []}
+        self, g, SI_units=None, properties_to_convert={"lengths": [], "areas": []}
     ):
         """Update MTG geometry.
 
         :Parameters:
             - `g` (:class:`openalea.mtg.mtg.MTG`) - The MTG to update the geometry.
-            - `SI_units` (:class:`bool`) - A boolean indicating whether the MTG properties are expressed in SI units.
+            - `SI_units` (:class:`bool`) - deprecated : use scene _unit to specify unit for scene and properties.
             - `properties_to_convert` (:class:`dict` of :class:`pandas.DataFrame`) - A dictionnary with the list of length properties area properties to be converted.
         :Returns:
             MTG with updated geometry
         :Returns Type:
             :class:`openalea.mtg.mtg.MTG`
         """
 
-        if SI_units:
-            self.convert_to_ADEL_units(g, properties_to_convert)
+        if SI_units is not None:
+            raise ValueError('SI_units is deprecated, use scene _unit to specify unit for both scene and mtg properties, use dev_T unit to provide adel inputs in another unit')
+            #self.convert_to_ADEL_units(g, properties_to_convert)
 
         # update elements
         g = update_organ_elements(g, self.leaves, self.split, self.phyllochron())
+        inclinations = self.axe_inclination(g)
+        newinc = inclinations.set_index("vid")["inclination"].to_dict()
+        g.property('inclination').update(newinc)
         g = mtg_interpreter(g, self.leaves, min_length=self.min_length, face_up=self.face_up, classic=self.classic)
         pos = g.property("position")
         az = g.property("azimuth")