[Soiltexture-commits] r95 - / pkg/soiltexture pkg/soiltexture/R pkg/soiltexture/inst pkg/soiltexture/man pkg/soiltexture/tests pkg/soiltexture/vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Jan 9 15:40:24 CET 2014
Author: jmoeys
Date: 2014-01-09 15:40:24 +0100 (Thu, 09 Jan 2014)
New Revision: 95
Removed:
1_RCMDBUILD_novignette.bat
1_RCMDBUILD_novignette.sh
pkg/soiltexture/R/text.transf.R
Modified:
2_RCMDcheck_noexamples.bat
3_RCMDINSTALL_build.BAT
pkg/soiltexture/DESCRIPTION
pkg/soiltexture/NEWS
pkg/soiltexture/R/onAttach.R
pkg/soiltexture/R/soiltexture.R
pkg/soiltexture/inst/polish_language_ANSI.r
pkg/soiltexture/man/TT.classes.tbl.Rd
pkg/soiltexture/man/TT.plot.Rd
pkg/soiltexture/man/TT.points.in.classes.Rd
pkg/soiltexture/man/TT.vertices.plot.Rd
pkg/soiltexture/man/TT.vertices.tbl.Rd
pkg/soiltexture/man/soiltexture-package.Rd
pkg/soiltexture/tests/TT.plot.R
pkg/soiltexture/tests/TT.points.in.classes.R
pkg/soiltexture/vignettes/soiltexture_vignette.Rnw
soiltexture_compile.R
Log:
FAO texture triangle renamed HYPRES texture triangle (name was erroneous)
Deleted: 1_RCMDBUILD_novignette.bat
===================================================================
--- 1_RCMDBUILD_novignette.bat 2014-01-09 10:12:27 UTC (rev 94)
+++ 1_RCMDBUILD_novignette.bat 2014-01-09 14:40:24 UTC (rev 95)
@@ -1,5 +0,0 @@
-c:
-cd \
-cd "C:\_R_PACKAGES\soiltexture\pkg"
-R CMD build --no-vignettes soiltexture
-pause
Deleted: 1_RCMDBUILD_novignette.sh
===================================================================
--- 1_RCMDBUILD_novignette.sh 2014-01-09 10:12:27 UTC (rev 94)
+++ 1_RCMDBUILD_novignette.sh 2014-01-09 14:40:24 UTC (rev 95)
@@ -1,5 +0,0 @@
-cd "/home/jules/Documents/_WORK/R_PACKAGES/soiltexture/pkg"
-R CMD build --no-vignettes soiltexture
-echo "press any key to continue"; read line
-
-
Modified: 2_RCMDcheck_noexamples.bat
===================================================================
--- 2_RCMDcheck_noexamples.bat 2014-01-09 10:12:27 UTC (rev 94)
+++ 2_RCMDcheck_noexamples.bat 2014-01-09 14:40:24 UTC (rev 95)
@@ -2,8 +2,8 @@
cd /D "%rPackagesDir%\%pkgname%\pkg"
-R CMD check --no-examples %pkgname%
+R CMD check --no-examples --as-cran %pkgname%_1.2.12.tar.gz
- at REM --as-cran _0.2.11.tar.gz
+ at REM --as-cran
pause
Modified: 3_RCMDINSTALL_build.BAT
===================================================================
--- 3_RCMDINSTALL_build.BAT 2014-01-09 10:12:27 UTC (rev 94)
+++ 3_RCMDINSTALL_build.BAT 2014-01-09 14:40:24 UTC (rev 95)
@@ -2,6 +2,6 @@
cd /D "%rPackagesDir%\%pkgname%\pkg"
-R CMD INSTALL --build --compact-docs --byte-compile %pkgname%
+R CMD INSTALL --build --compact-docs --byte-compile %pkgname%_1.2.12.tar.gz
pause
Modified: pkg/soiltexture/DESCRIPTION
===================================================================
--- pkg/soiltexture/DESCRIPTION 2014-01-09 10:12:27 UTC (rev 94)
+++ pkg/soiltexture/DESCRIPTION 2014-01-09 14:40:24 UTC (rev 95)
@@ -1,6 +1,6 @@
Package: soiltexture
-Version: 1.2.11
-Date: 2014-01-08
+Version: 1.2.12
+Date: 2014-01-09
Title: Functions for soil texture plot, classification and transformation
Author: Julien MOEYS <jules_m78-soiltexture at yahoo.fr>, contributions from Wei Shangguan.
Maintainer: Julien MOEYS <jules_m78-soiltexture at yahoo.fr>
Modified: pkg/soiltexture/NEWS
===================================================================
--- pkg/soiltexture/NEWS 2014-01-09 10:12:27 UTC (rev 94)
+++ pkg/soiltexture/NEWS 2014-01-09 14:40:24 UTC (rev 95)
@@ -1,16 +1,31 @@
The Soil Texture Wizard
change log. From 2009/10/09. Most recent changes first / on top.
Julien MOEYS
------------------------------------------------------------------------------------------
+-----------------------------------------------------------------
+VERSION 1.2.12
+
+ 2014/01/09 The triangle named FAO50.TT was in fact erroneously
+ named, as it is the HYPRES / European Soil Map
+ texture triangle. The option class.sys = "FAO50.TT"
+ is now giving a warning, and the new option is
+ class.sys = "HYPRES.TT". See the vignette for
+ details.
+
VERSION 1.2.11
2014/01/08 Compiled under R 3.0.2
- The functions: TT.check.ps.lim.Xm() and TT.text.transf.Xm()
- have been moved to the package
+ The functions: TT.check.ps.lim.Xm() and
+ TT.text.transf.Xm() have been moved to the package
soiltexturetransformation, available on r-forge
(https://r-forge.r-project.org/R/?group_id=740)
+
+ The polish triangle is now loaded at startup, but
+ loading occur within a try() statement, to (try to)
+ prevent possible encoding problems on some platforms.
+ A message is issued if the triangle was not loaded
+ properly (but no error or warning).
VERSION 1.2.10
Modified: pkg/soiltexture/R/onAttach.R
===================================================================
--- pkg/soiltexture/R/onAttach.R 2014-01-09 10:12:27 UTC (rev 94)
+++ pkg/soiltexture/R/onAttach.R 2014-01-09 14:40:24 UTC (rev 95)
@@ -20,6 +20,9 @@
package = pkgname ) ) )
}
+
+ # Extend language parameter if the polish triangle was loaded
+ # successfully
if( !("try-error" %in% class( tryRes )) ){
lang.par <- TT.get( "lang.par" )
@@ -27,12 +30,23 @@
lang.par,
lang.par2
)
+
+ TT.set( "lang.par" = lang.par )
+ }else{
+ packageStartupMessage( "soiltexture: The polish triangle could not be loaded" )
}
- TT.set( "lang.par" = lang.par )
- packageStartupMessage(
- paste("'", pkgname, "' loaded.", sep = "" )
- ) #
-} #
+ # Welcome message
+ if( interactive() ){
+ msg <- sprintf(
+ "%s %s For help type: help(pack='%s')",
+ pkgname,
+ as.character( packageVersion( pkgname ) ),
+ pkgname )
+
+ packageStartupMessage( msg )
+ }
+
+}
Modified: pkg/soiltexture/R/soiltexture.R
===================================================================
--- pkg/soiltexture/R/soiltexture.R 2014-01-09 10:12:27 UTC (rev 94)
+++ pkg/soiltexture/R/soiltexture.R 2014-01-09 14:40:24 UTC (rev 95)
@@ -71,7 +71,7 @@
# | - 7 pre-parameterized national Soil Texture Triangles/Classes |
# | systems, including classes polygons, classes and axis names, |
# | and geometric settings of the triangle. |
-# | (FAO, USDA, French-Aisne, French-GEPPA, German, UK and |
+# | (EU, USDA, French-Aisne, French-GEPPA, German, UK and |
# | Australian triangles) |
# | - Capable of projecting any Soil Texture Triangles/Classes systems|
# | in another family of ternary plot (custom angles and clock), |
@@ -337,7 +337,7 @@
at = seq(from=.1,to=.9,by=.1),
#
# CLASSES SYSTEM (polygons / Texture Triangle) used by default:
- class.sys = "FAO50.TT",
+ class.sys = "HYPRES.TT",
class.lab.col = NULL, # Color of classes names (abreviation)
class.p.bg.col = FALSE, # Fill classes polygon with color gradient ???
class.p.bg.hue = 0.04, # Hue (unique) of the classes polygon color gradient
@@ -432,14 +432,14 @@
trsf.add.opt1 = NA, # Additionnal option 1
trsf.add.opt2 = NA, # Additionnal option 2
#
- # +---------------------------------+
- # | HYPRES TEXTURE TRIANGLE |
- # +---------------------------------+
+ # +-------------------------------------------------+
+ # | HYPRES TEXTURE TRIANGLE -- ORIGINAL, WRONG NAME |
+ # +-------------------------------------------------+
#
FAO50.TT = list( # FAO TRIANGLE PARAMETERS :
- #
- main = "FAO/HYPRES",
- #
+
+ main = "HYPRES (please use 'HYPRES.TT' instead)",
+
# The list below specify the CSS coordinates of the different POINTS
# that are used to draw soil texture classes. One points can be
# used by several classes :
@@ -502,7 +502,82 @@
# MEDIUM-COARSE
#
), #
- #
+
+ # +----------------------------------+
+ # | HYPRES TEXTURE TRIANGLE, RENAMED |
+ # +----------------------------------+
+
+ HYPRES.TT = list( # EU SOIL MAP TRIANGLE PARAMETERS :
+
+ main = "HYPRES / European Soil Map",
+
+ # The list below specify the CSS coordinates of the different POINTS
+ # that are used to draw soil texture classes. One points can be
+ # used by several classes :
+ # =-P01- P02 P03 P04 P05 P06 P07 P08 -P09- P10 P11 -P12-
+ "tt.points" = data.frame(
+ "CLAY" = c( 1.000, 0.600, 0.600, 0.350, 0.350, 0.350, 0.180, 0.180, 0.000, 0.000, 0.000, 0.000),
+ "SILT" = c( 0.000, 0.000, 0.400, 0.000, 0.500, 0.650, 0.000, 0.170, 0.000, 0.350, 0.850, 1.000),
+ "SAND" = c( 0.000, 0.400, 0.000, 0.650, 0.150, 0.000, 0.820, 0.650, 1.000, 0.650, 0.150, 0.000)
+ ), #
+
+ # Abreviations; Names of the texture cl; Points marking the class limits (points specified above)
+ "tt.polygons" = list(
+ "VF" = list( "name" = "Very fine", "points" = c(02,01,03) ),
+ "F" = list( "name" = "Fine", "points" = c(04,02,03,06) ),
+ "M" = list( "name" = "Medium", "points" = c(07,04,05,11,10,08) ),
+ "MF" = list( "name" = "Medium fine", "points" = c(11,05,06,12) ),
+ "C" = list( "name" = "Coarse", "points" = c(09,07,08,10) )
+ ), #
+
+ # Traingle specific parameters for triangle geometry / appearance
+ # See general parameters above for detailed description of them
+ blr.clock = rep(T,3),
+ tlr.an = c(60,60,60),
+ #
+ blr.tx = c("SAND","CLAY","SILT"),
+ #
+ base.css.ps.lim = c(0,2,50,2000),
+ tri.css.ps.lim = c(0,2,50,2000),
+ #
+ unit.ps = quote(bold(mu) * bold('m')),
+ unit.tx = quote(bold('%')),
+ #
+ text.sum = 100
+ #
+ # In fact it is the FAO soil texture classes: Info from SysCan
+ # http://sis.agr.gc.ca/cansis/nsdb/lpdb/faotext.html
+ # FAO Soil Texture
+ # Texture is the relative proportion of sand, silt and clay of the dominant
+ # soil for each soil map polygon. Texture classes are:
+ #
+ # Coarse texture: sands, loamy sand and sandy loams with less than 18 % clay,
+ # and more than 65 % sand.
+ #
+ # Medium texture: sandy loams, loams, sandy clay loams, silt loams with less
+ # than 35 % clay and less than 65 % sand; the sand fractions may be as high as 82 % if a minimum of 18 % clay is present.
+ #
+ # Fine texture: clays, silty clays, sandy clays, clay loams and silty clay loams
+ # with more than 35 % clay.
+ #
+ # Where two or three texture names appear, this means that all named textures
+ # are present in the map unit.
+ #
+ # Texture Codeset
+ # COARSE
+ # FINE
+ # FINE-COARSE
+ # FINE-MED-CRS
+ # FINE-MEDIUM
+ # MEDIUM
+ # MEDIUM-COARSE
+ #
+ ), #
+
+ # +-----------------+
+ # | OTHER TRIANGLES |
+ # +-----------------+
+
USDA.TT = list( # USDA Triangle parameters
#
main = "USDA",
@@ -2060,7 +2135,11 @@
base.css.ps.lim = NULL
){ #
if( is.null(class.sys) ){ class.sys <- TT.get("class.sys") }
- #
+
+ if( class.sys == "FAO50.TT" ){
+ warning( "class.sys = 'FAO50.TT' must be replaced by class.sys = 'HYPRES.TT'. See the package vignette." )
+ }
+
geo.par <- c("blr.clock","tlr.an","blr.tx","text.sum","base.css.ps.lim")
#
p.env <- environment()
@@ -2071,11 +2150,11 @@
is.null( get(x=X,envir=p.env) )
} #
) ) #
- #
+
# +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~+
# Attributes either classes-system or default values
# to triangle geometry parameters (and others)
- #
+
if( any( null.geo.par ) )
{ #
geo.par <- geo.par[ null.geo.par ]
@@ -4329,7 +4408,7 @@
### clay silt sand coordinates of the triangle classes vertices.
### See also TT.vertices.plot().
- class.sys = "FAO50.TT",
+ class.sys = "HYPRES.TT",
collapse = NULL
){ #
TT.data <- TT.get( class.sys )
@@ -4368,7 +4447,7 @@
### TT.classes.tbl() to know the vertices that bounds each texture
### class. See also TT.vertices.plot().
- class.sys = "FAO50.TT"
+ class.sys = "HYPRES.TT"
){ #
TT.data <- TT.get( class.sys )
#
@@ -4398,7 +4477,7 @@
### for a non graphic, tabular equivalent of the plot.
geo,
- class.sys = "FAO50.TT",
+ class.sys = "HYPRES.TT",
fg = NULL,
col = NULL,
cex = NULL,
@@ -5157,9 +5236,8 @@
### geometry and particle class size system of the plot, unless
### higher level options are chosen (see the function definition).
### Possible values are "none" (no classification plotted), "USDA.TT"
-### (USDA texture triangle), "FAO50.TT" (FAO texture triangle with a 50
-### microns silt-sand limit. DEFAULT VALUE), "FR.AISNE.TT" (French
-### texture triangle
+### (USDA texture triangle), "HYPRES.TT" (texture triangle of the
+### European Sil Map), "FR.AISNE.TT" (French texture triangle
### of the Aisne region soil survey), "FR.GEPPA.TT" (French GEPPA
### texture triangle), "DE.BK94.TT" (German texture triangle),
### "UK.SSEW.TT" (Soil Survey of England and Wales), "AU.TT"
@@ -5700,9 +5778,9 @@
### Single text string. Text code of the texture classification
### system to be used for the classification of 'tri.data'.
### Possible values are "none" (no classification plotted), "USDA.TT"
-### (USDA texture triangle), "FAO50.TT" (FAO texture triangle with a 50
-### microns silt-sand limit. DEFAULT VALUE), "FR.AISNE.TT" (French
-### texture triangle of the Aisne region soil survey), "FR.GEPPA.TT" (French GEPPA
+### (USDA texture triangle), "HYPRES.TT" (texture triangle of the
+### European Soil Map), "FR.AISNE.TT" (French texture triangle of
+### the Aisne region soil survey), "FR.GEPPA.TT" (French GEPPA
### texture triangle), "DE.BK94.TT" (German texture triangle),
### "UK.SSEW.TT" (Soil Survey of England and Wales), "AU.TT"
### (Australian texture triangle), "BE.TT" (Belgium texture triangle),
@@ -7081,7 +7159,7 @@
# TT.contour( x = kde.res, geo = geo, add = FALSE, lwd = 2 )
# TT.plot(
-# class.sys = "FAO50.TT",
+# class.sys = "HYPRES.TT",
# geo = geo,
# grid.show = FALSE,
# add = TRUE
@@ -7317,7 +7395,7 @@
# TT.contour( x = kde.res, geo = geo, add = TRUE, lwd = 2 )
# TT.plot(
-# class.sys = "FAO50.TT",
+# class.sys = "HYPRES.TT",
# geo = geo,
# grid.show = FALSE,
# add = TRUE
@@ -7341,7 +7419,7 @@
# TT.contour( x = iwd.res, geo = geo, add = TRUE, lwd = 2 )
# TT.plot(
-# class.sys = "FAO50.TT",
+# class.sys = "HYPRES.TT",
# geo = geo,
# grid.show = FALSE,
# add = TRUE
Deleted: pkg/soiltexture/R/text.transf.R
===================================================================
--- pkg/soiltexture/R/text.transf.R 2014-01-09 10:12:27 UTC (rev 94)
+++ pkg/soiltexture/R/text.transf.R 2014-01-09 14:40:24 UTC (rev 95)
@@ -1,559 +0,0 @@
-# source( "C:/_RTOOLS/SWEAVE_WORK/SOIL_TEXTURES/rforge/pkg/soiltexture/R/text.transf.r" )
-
-
-
-TT.check.ps.lim.Xm <- function(# Internal. Check the consistency between 'base.ps.lim' and 'dat.ps.lim'.
-### Check the consistency between 'base.ps.lim' and 'dat.ps.lim'.
-### 4 tests performed.
-
- base.ps.lim,
-
- dat.ps.lim,
-
- ps.lim.length=c(4,4)
-### vector of 2 integers. Number of particle size classes + 1. c(base,dat)
-
-##author<<Wei Shangguan
-
-){ #
- # if( length( base.ps.lim ) != length( dat.ps.lim ) )
- # { #
- # stop( paste(
- # sep="",
- # "The length of the 'base' particle size classes limits must be equal to\n",
- # "the length of the 'dat' particle size classes limits.\n",
- # "Either check the 'base' particle size classes limits vector,\n",
- # "or check number of column in tri.data.\n"
- # ) ) #
- # } #
- #
- if( length( base.ps.lim ) != ps.lim.length[1] )
- { #
- stop( paste(
- sep="",
- "The length of the 'base' particle size classes limits must be equal to\n",
- ps.lim.length[1], " (number of particle size classes+1; from ps min to ps.max)\n",
- "Actual value: ", length( base.ps.lim ), ".\n",
- "Either check the 'base' particle size classes limits,\n",
- "or check number of column in tri.data.\n"
- ) ) #
- } #
- #
- if( length( dat.ps.lim ) != ps.lim.length[2] )
- { #
- stop( paste(
- sep="",
- "The length of the 'dat' particle size classes limits must be equal to\n",
- ps.lim.length[2], " (number of particle size classes +1; from ps min to ps.max)\n",
- "Actual value: ", length( dat.ps.lim ), ".\n",
- "Either check the 'dat' particle size classes limits,\n",
- "or check number of column in tri.data.\n"
- ) ) #
- } #
- #
- if( base.ps.lim[1] != dat.ps.lim[1] )
- { #
- stop( paste(
- sep="",
- "The first value of the 'dat' particle size classes limits must be equal to\n",
- "the first value of the 'base' particle size classes limits.\n",
- "Actual value, base: ", base.ps.lim[1], ", dat: ", dat.ps.lim[1]
- ) ) #
- } #
- #
- if( base.ps.lim[ps.lim.length[1]] != dat.ps.lim[ps.lim.length[2]] )
- { #
- stop( paste(
- sep="",
- "The last value of the 'dat' particle size classes limits must be equal to\n",
- "the last value of the 'base' particle size classes limits.\n",
- "Actual value, base: ", base.ps.lim[ps.lim.length[1]], ", dat: ", dat.ps.lim[ps.lim.length[2]]
- ) ) #
- } #
- #
-# if( base.ps.lim[1] == 0 )
-# { #
-# if( base.ps.lim[2] < dat.ps.lim[2] )
-# stop( paste(
-# sep="",
-# "When the 1st value of 'dat' and 'base' particle size classes limits is 0\n",
-# "The 2nd value of the 'base' particle size classes limits must higher or equal to\n",
-# "the 2nd value of the 'dat' particle size classes limits.\n"
-# ) ) #
-# } #
-} #
-
-
-
-
-
-
-TT.text.transf.Xm <- function(# Transformations of a soil texture data table between 2 particle size systems (X classes), various methods.
-### using various Particle Size Distribution (PSD) models including Anderson (AD), Fredlund4P (F4P), Fredlund3P (F3P),
-### modified logistic growth (ML), Offset-Nonrenormalized Lognormal (ONL), Offset-Renormalized Lognormal (ORL),
-### Skaggs (S), van Genuchten type(VG),van Genuchten modified(VGM), Weibull (W), Logarithm(L),
-### Logistic growth (LG), Simple Lognormal (SL),Shiozawa and Compbell (SC).
-### The performance of PSD models is influenced by many aspects like soil texture class,
-### number and position (or closeness) of observation points, clay content etc.
-### The latter four PSD models perform worse than the former ten.
-### The AD, F4P, S, and W model is recommended for most of texture classes.
-### And it will be even better to compare several different PSD models and using the results of the model
-### with the minimum residual sum of squares (or other measures).
-### Sometimes, the fitting will failed for the iteration is not converged or some errors happened.
-### Transformation of a soil texture data table
-### ('tri.data') from one
-### particle size system ('dat.ps.lim') into another
-### ('base.ps.lim'). No limit in the number of texture classes
-### in the input and output texture tables. See TT.text.transf
-### for transformation involving only 3 particle classes. 'tri.data'
-### can only contain texture data.
-### This function requires the "drc" package to be installed.
-##author<<Wei Shangguan
-
-tri.data,
-
-base.ps.lim,
-
-dat.ps.lim,
-
-text.sum= NULL,
-
-text.tol= NULL,
-
-tri.sum.tst= NULL,
-
-tri.pos.tst= NULL,
-
-psdmodel= "AD",
-
-omethod= "all",
-### see optim for available methods, the default "all" is to run all methods and
-### choose the best results with minimum residual sum of squares (RSS).
-
-tri.sum.norm=FALSE
-###Weather the sum of is
-
-){#
- TT.auto.set( set.par = FALSE )
- #
- TT.data.test.X(
- tri.data = tri.data,
- text.sum = text.sum,
- text.tol = text.tol,
- tri.sum.tst = tri.sum.tst,
- tri.pos.tst = tri.pos.tst
- ) #
- #
- tri.data <- t( apply(
- X = tri.data,
- MARGIN = 1,
- FUN = function(X){
- cumsum(X)/100
- } #
- ) ) #
- #
- ps.end <- dim( tri.data )[2] + 1
- #
- TT.check.ps.lim.Xm(
- base.ps.lim = base.ps.lim,
- dat.ps.lim = dat.ps.lim,
- ps.lim.length = c(length(base.ps.lim),ps.end)
- ) #
- #
- if( base.ps.lim[1] != 0 )
- { #
- tri.data <- cbind(
- "ZERO" = rep(0,dim(tri.data)[1]),
- tri.data
- ) #
- #
- ps.start <- 1
- }else{
- ps.start <- 2
- } #
- #Particle size diameter in micro-meters (will be converted in milli-meters)
- base.ps.lim <- base.ps.lim/1000
- dat.ps.lim <- dat.ps.lim/1000
-# base.ps.lim2 <- TT.dia2phi(base.ps.lim)
-# dat.ps.lim2 <- TT.dia2phi(dat.ps.lim)
- #
-# old.col.nm <- colnames( tri.data )
-
- # # Added 2010/06/13 Julien Moeys # Removed on 2012/03/07 by Julien Moeys
- # if( !"drc" %in% as.character( installed.packages()[,1] ) )
- # { #
- # stop( "The function 'TT.text.transf.Xm' needs the package 'drc'\n Please install it ( install.packages(\"drc\") )" )
- # }else{
- # require( "drc" )
- # } #
- #
- require( "drc" )
- #
- fitpsd <- function(
- y,
- xin,
- xout,
- psdmodel,
- method)
- {
- require( "drc" ) # Added 2010/08/11 by JM
- #
- #default max and min of initial parameters
- maxspa1 <- 1
- minspa1 <- 0.1
- maxspa2 <- 1
- minspa2 <- 0.1
- #erf error function for ONL and ORL model
- erf <- function(x) 2 * pnorm(x * sqrt(2)) - 1
- #r0 for S model
- r0 <- xin[1]
- #dmax, dmin for W model
- dmax <- xin[length(xin)]
- dmin <- r0
- # Particle Size Distribution (PSD) models
- AD <- function(dose, parm)
- {parm[, 1] + parm[, 2]*atan(parm[, 3]*log10(dose/parm[, 4]))}
- F4P <- function(dose, parm)
- {(1-(log(1+parm[,1]/dose)/log(1+parm[,1]/0.0001))^7)/(log(exp(1)+(parm[,2]/dose)^parm[,3]))^parm[,4]}
- F3P <- function(dose, parm)
- {(1-(log(1+0.001/dose)/log(1+0.001/0.0001))^7)/(log(exp(1)+(parm[,1]/dose)^parm[,2]))^parm[,3]}
- ML <- function(dose, parm)
- {1/(1+parm[,1]*exp(-parm[,2]*dose^parm[,3]))}
- ONL <- function(dose, parm)
- {
- t <- (-1)^(log(dose) >= parm[,1]+1)
- (1+t*erf((log(dose)+parm[,1])/parm[,2]*2^0.5))/2+parm[,3]
- }
- ORL <- function(dose, parm)
- {
- t <- (-1)^(log(dose) >= parm[,1]+1)
- (1-parm[,3])*(1+t*erf((log(dose)+parm[,1])/parm[,2]*2^0.5))/2+parm[,3]
- }
- # no ability to predict the content below r0
- S <- function(dose, parm)
- {1/(1+(1/y[1]-1)*exp(-parm[,1]*((dose-r0)/r0)^parm[,2]))}
- VG <- function(dose, parm)
- {(1+(parm[,1]/dose)^parm[,2])^(1/parm[,2]-1)}
- # old form is right
- VGM <- function(dose, parm)
- {y[1]+(1-y[1])*(1+(parm[,1]*dose)^(-parm[,2]))^(-1/parm[,2]-1)}
- #This is the wrong form of VGM
-# VGM <- function(dose, parm)
-# {y[1]+(1-y[1])*(1+(parm[,1]*dose)^(-parm[,2]))^(1/parm[,2]-1)}
- # no ability to predict the content below dmin
- W <- function(dose, parm)
- {parm[,3]+(1-parm[,3])*(1-exp(-parm[,1]*((dose-dmin)/(dmax-dmin))^parm[,2]))}
- L <- function(dose, parm)
- {parm[,1]*log(dose)+parm[,2]}
- LG <- function(dose, parm)
- {1/(1+parm[,1]*exp(-parm[,2]*dose))}
- SC <- function(dose, parm)
- {
- t <-(-1)^(log(dose) >= parm[,1]+1)
- (1-parm[,3])*(1+t*erf((log(dose)+parm[,1])/parm[,2]*2^0.5))/2+parm[,3]*(1+t*erf((log(dose)+1.96)/1*2^0.5))/2
- }
- SL <- function(dose, parm)
- {
- t <-(-1)^(log(dose) >= parm[,1]+1)
- (1+t*erf((log(dose)+parm[,1])/parm[,2]*2^0.5))/2
- }
-
- if( psdmodel == "AD" )
- {
- logi <- AD
- pn <- 4
- pname <- c("f0", "b", "c", "r0")
- }
- else if ( psdmodel == "F4P" )
- {
- logi <- F4P
- pn <- 4
- pname <- c("df","a","n","m")
- }
- else if ( psdmodel == "F3P" )
- {
- logi <- F3P
- pn <- 3
- pname <- c("a","n","m")
- }
- else if ( psdmodel == "ML" )
- {
- logi <- ML
- pn <- 3
- pname <- c("a","b","c")
- }
- else if ( psdmodel == "ONL" )
- {
- logi <- ONL
- pn <- 3
- pname <- c("u","o","c")
- }
- else if ( psdmodel == "ORL" )
- {
- logi <- ORL
- pn <- 3
- pname <- c("u","o","e")
- }
- else if ( psdmodel == "S" )
- {
- logi <- S
- pn <- 2
- pname <- c("u","c")
- #S model can not deal with first texture data with zero value
- if(y[1] == 0) y[1] <- 0.0001
- }
- else if ( psdmodel == "VG" )
- {
- logi <- VG
- pn <- 2
- pname <- c("dg","n")
- maxspa2 <- 2
- minspa2 <- 1
- }
- else if ( psdmodel == "VGM" )
- {
- logi <- VGM
- pn <- 2
- pname <- c("dg","n")
- maxspa1 <- 200
- minspa1 <- 4
- maxspa2 <- 2
- minspa2 <- 0.5
- }
- else if ( psdmodel == "W" )
- {
- logi <- W
- pn <- 3
- pname <- c("a","b","c")
- }
- else if ( psdmodel == "L" )
- {
- logi <- L
- pn <- 2
- pname <- c("a","b")
- }
- else if ( psdmodel == "LG" )
- {
- logi <- LG
- pn <- 2
- pname <- c("a","b")
- }
- else if ( psdmodel == "SC" )
- {
- logi <- SC
- pn <- 3
- pname <- c("u","o","e")
- }
- else if ( psdmodel == "SL" )
- {
- logi <- SL
- pn <- 2
- pname <- c("u","o")
- }
- #default lower and upper limit for drc::drm, these values should not set
- #at the beginning of the function for pn is set later
- lowerl <- rep(10e-9,times=pn)
- upperl <- rep(10e+5,times=pn)
- #Initailize spa for drc::drm
- spa <- c(1,1,1,1)
- #methods used in optim() of drc::drm
- meth <- c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN")
-
- mdev <- 100
- for( i in 1:5 ) # The nonlinear optimization runs were carried out using at least
- # five random initial parameter estimates for all soils.
- #When the final solution for each soil converged to different parameter values,
- #the parameter values with the best fitting statistics (RSS) were kept.
- {
- if( method == "all" )# using all optim methods
- {
- for( j in 1:5 )
- {
- countR <- 0 # Added by Julien Moeys on 2011/11/01
- repeat
- {
- countR <- countR + 1 # Added by Julien Moeys on 2011/11/01
-
- spa[1:pn-1] <- runif(n = pn-1,max = maxspa1,min = minspa1)
- spa[pn] <- runif(n = 1,max = maxspa2,min = minspa2)
- tt<- try( drc:::drm(y ~ xin, fct = list( logi, NULL,pname ), # JM:2010/08/11 changed drc::drm to drm alone
- start = spa[1:pn],
- #roust = "median",
- lowerl = lowerl,
- upperl = upperl,
- control = drc:::drmc(constr = TRUE,maxIt = 500, # JM:2010/08/11 changed drc::drmc to drmc alone
- noMessage = TRUE,
- method = meth[j],
- # method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN"),
- #trace = TRUE
- )
- )
- , silent = TRUE
- )
- if( !inherits(tt, "try-error") )
- {
- dev <- sum(residuals(tt)^2)
- if( mdev > dev )
- {
- mdev <- dev
- ttbest <- tt
- }
- break
- } #
-
- if( countR >= 100 ){ break } # Added by Julien Moeys on 2011/11/01
- }
- }
- }
- else
- {
- countR <- 0 # Added by Julien Moeys on 2011/11/01
- repeat
- {
- countR <- countR + 1 # Added by Julien Moeys on 2011/11/01
-
- spa[1:pn-1] <- runif(n=pn-1,max = maxspa1,min = minspa1)
- spa[pn] <- runif(n=1,max = maxspa2,min = minspa2)
- tt <- try( drc:::drm(y ~ xin, fct = list(logi, NULL, pname), # JM:2010/08/11 changed drc::drm to drm alone
- start = spa[1:pn],
- #roust = "median",
- lowerl = lowerl,
- upperl = upperl,
- control = drc:::drmc(constr = TRUE,maxIt = 500,# JM:2010/08/11 changed drc::drmc to drmc alone
- noMessage = TRUE,
- method = method,
- # method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN"),
- #trace = TRUE
- )
- )
- , silent = TRUE
- )
- if( !inherits(tt, "try-error") )
- {
- dev<-sum(residuals(tt)^2)
- if( mdev>dev )
- {
- mdev <- dev
- ttbest <- tt
- }
- break
- }
-
- if( countR >= 100 ){ break } # Added by Julien Moeys on 2011/11/01
- }
- }
- #when the residual sum of error (deviance) is very small, the iteration is stopped to save time
- if(mdev < 0.0001) break
- }
- #predict() has some bug for PSD model to predict the target values
- if( psdmodel == "AD" )
- {
- pre <- coef(ttbest)[1] + coef(ttbest)[2]*atan(coef(ttbest)[3]*log10(xout/coef(ttbest)[4]))
- }
- else if( psdmodel == "F4P" )
- {
- pre <- (1-(log(1+coef(ttbest)[1]/xout)/log(1+coef(ttbest)[1]/0.0001))^7)/(log(exp(1)+(coef(ttbest)[2]/xout)^coef(ttbest)[3]))^coef(ttbest)[4]
- }
- else if( psdmodel == "F3P" )
- {
- pre <- (1-(log(1+0.001/xout)/log(1+0.001/0.0001))^7)/(log(exp(1)+(coef(ttbest)[1]/xout)^coef(ttbest)[2]))^coef(ttbest)[3]
- }
- else if( psdmodel == "ML" )
- {
- pre <- 1/(1+coef(ttbest)[1]*exp(-coef(ttbest)[2]*xout^(coef(ttbest)[3])))
- }
- else if( psdmodel == "ONL" )
- {
- t <- (-1)^(log(xout) >= coef(ttbest)[1]+1)
- pre <- (1+t*erf((log(xout)+coef(ttbest)[1])/coef(ttbest)[2]*2^0.5))/2+(coef(ttbest)[3])
- }
- else if( psdmodel == "ORL" )
- {
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/soiltexture -r 95
More information about the Soiltexture-commits
mailing list