[Distr-commits] r729 - branches/distr-2.4/pkg/distr/R branches/distr-2.4/pkg/distr/man branches/distr-2.4/pkg/distrEx/src pkg/distr pkg/distr/R pkg/distr/man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Thu Sep 1 03:35:32 CEST 2011


Author: ruckdeschel
Date: 2011-09-01 03:35:31 +0200 (Thu, 01 Sep 2011)
New Revision: 729

Removed:
   branches/distr-2.4/pkg/distrEx/src/distrEx.dll
Modified:
   branches/distr-2.4/pkg/distr/R/AllClasses.R
   branches/distr-2.4/pkg/distr/R/AllInitialize.R
   branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R
   branches/distr-2.4/pkg/distr/R/LatticeDistribution.R
   branches/distr-2.4/pkg/distr/R/internalUtils.R
   branches/distr-2.4/pkg/distr/man/internals.Rd
   pkg/distr/DESCRIPTION
   pkg/distr/R/AllClasses.R
   pkg/distr/R/AllInitialize.R
   pkg/distr/R/DiscreteDistribution.R
   pkg/distr/R/LatticeDistribution.R
   pkg/distr/R/internalUtils.R
   pkg/distr/man/internals.Rd
Log:
[pkg distr]
+ revised initialize method convolution method and generating function for LatticeDistribution -> new routine to determine the smallest common grid for convolution
+ commented out require command in .onAttach
+ corrected error in definition of Beta() --- d slot was wrong
+ corrected type (inifinite...)


Modified: branches/distr-2.4/pkg/distr/R/AllClasses.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/AllClasses.R	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/AllClasses.R	2011-09-01 01:35:31 UTC (rev 729)
@@ -1,8 +1,9 @@
 ############# preparations ################
 # (.onload, .onattach ...)
-.onLoad <- function(lib, pkg) { # extended 03-28-06: P.R. 
-    require("methods", character = TRUE, quietly = TRUE)
-}
+### left out P.R. 01-09-11
+##.onLoad <- function(lib, pkg) { # extended 03-28-06: P.R. 
+##    require("methods", character = TRUE, quietly = TRUE)
+##}
 
 
 
@@ -722,7 +723,7 @@
                                      },
                       d = function(x, log = FALSE){
                               dbeta(x,  shape1 = 1, shape2 = 1, ncp = 0,
-                                    lower.tail = lower.tail, log.p = log.p)
+                                    log = log)
                                           },
                       p = function(q, lower.tail = TRUE, log.p = FALSE ){
                               pbeta(q,  shape1 = 1, shape2 = 1, ncp = 0,

Modified: branches/distr-2.4/pkg/distr/R/AllInitialize.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/AllInitialize.R	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/AllInitialize.R	2011-09-01 01:35:31 UTC (rev 729)
@@ -279,13 +279,14 @@
             
              OS  <- D at support 
 
-             if(is.null(lattice))  
-               {  if(! .is.vector.lattice(OS))
-                      stop("Support as given/generated is not a lattice.")
-                  .Object at lattice <- .make.lattice.es.vector(OS)
-             }else{
-                  .Object at lattice <- lattice 
-             }
+             #if(is.null(lattice))  
+             #  {  if(! .is.vector.lattice(OS))
+             #         stop("Support as given/generated is not a lattice.")
+             #     .Object at lattice <- .make.lattice.es.vector(OS)
+             #}else{
+                  .Object at lattice <- if(is.null(lattice )) 
+                          new("Lattice") else lattice
+             #}
 
 
             .Object at support <- OS

Modified: branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R	2011-09-01 01:35:31 UTC (rev 729)
@@ -15,7 +15,7 @@
     if(!is.numeric(supp))
         stop("'supp' is no numeric vector")
     if(any(!is.finite(supp)))   # admit +/- Inf?
-        stop("inifinite or missing values in supp")
+        stop("infinite or missing values in supp")
     len <- length(supp)
     if(missing(prob)){
         prob <- rep(1/len, len)
@@ -23,7 +23,7 @@
         if(len != length(prob))
             stop("'supp' and 'prob' must have equal lengths")
         if(any(!is.finite(prob)))
-            stop("inifinite or missing values in prob")
+            stop("infinite or missing values in prob")
         if(!identical(all.equal(sum(prob), 1,
                           tolerance = 2*getdistrOption("TruncQuantile")), TRUE))
             stop("sum of 'prob' has to be (approximately) 1")
@@ -226,78 +226,15 @@
             e1.L <- as(e1, "LatticeDistribution")
             e2.L <- as(e2, "LatticeDistribution")
             if(is(e1.L, "LatticeDistribution") & is(e2.L, "LatticeDistribution"))
-                {w1 <- width(lattice(e1.L))
-                 w2 <- width(lattice(e2.L))
-                 W <- sort(abs(c(w1,w2)))
-                 if (abs(abs(w1)-abs(w2))<getdistrOption("DistrResolution") ||
-                     W[2] %% W[1] < getdistrOption("DistrResolution") )
-                     return(e1.L + e2.L)
-                } 
-            convolutedsupport <- rep(support(e1), each = length(support(e2))) +
-                                 support(e2)
+                  {w1 <- width(lattice(e1.L))
+                   w2 <- width(lattice(e2.L))
+                   W <- sort(abs(c(w1,w2)))
+                   if (abs(abs(w1)-abs(w2))<getdistrOption("DistrResolution") ||
+                       W[2] %% W[1] < getdistrOption("DistrResolution") )
+                       return(e1.L + e2.L)
+                  } 
+            .convDiscrDiscr(e1,e2)})
 
-            gridvalues1 <- d(e1)(support(e1)); gridvalues2 <- d(e2)(support(e2))
-            convolutedvalues <- rep(gridvalues1, each = length(support(e2))) *
-                                gridvalues2
-            rm(gridvalues1,gridvalues2)
-
-            tmptable <- data.frame(x = convolutedsupport, dx = convolutedvalues)
-            rm(convolutedsupport,convolutedvalues)
-            tmp <- tapply(tmptable$dx, tmptable$x, sum)
-            rm(tmptable)
-
-            supp.u <- as.numeric(names(tmp))
-            prob.u <- as.numeric(tmp)
-
-            o <- order(supp.u)
-            supp <- supp.u[o]
-            prob <- prob.u[o]
-
-            #supp.u <- unique(supp)
-
-            len <- length(supp)
-
-            if(len > 1){
-               if (min(diff(supp))< getdistrOption("DistrResolution")){
-                   if (getdistrOption("DistrCollapse")){
-                       erg <- .DistrCollapse(supp, prob, 
-                                   getdistrOption("DistrResolution"))
-                       if ( len > length(erg$prob) && 
-                                getdistrOption("DistrCollapse.Unique.Warn") )
-                            warning("collapsing to unique support values")         
-                       prob <- erg$prob
-                       supp <- erg$supp
-                   }else
-                    stop("grid too narrow --> change DistrResolution")
-                }
-            }
-
-            rm(tmp, len)
-
-            .withSim <- e1 at .withSim || e2 at .withSim
-
-            rfun <- function(n) {}
-            body(rfun) <- substitute({ f(n) + g(n) },
-                                         list(f = e1 at r, g = e2 at r))
-
-            dfun <- .makeDNew(supp, prob, Cont = FALSE)
-            pfun <- .makePNew(supp, prob, .withSim, Cont = FALSE)
-            qfun <- .makeQNew(supp, cumsum(prob), rev(cumsum(rev(prob))),
-                      .withSim, min(supp), max(supp), Cont = FALSE)
-
-            object <- new("DiscreteDistribution", r = rfun, d = dfun, p = pfun,
-                           q = qfun, support = supp,
-                           .withSim = .withSim, .withArith = TRUE)
-            rm(rfun, dfun, qfun, pfun)
-
-            if(is(e1 at Symmetry,"SphericalSymmetry")&& 
-               is(e2 at Symmetry,"SphericalSymmetry"))
-               object at Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
-                                                     SymmCenter(e2 at Symmetry))   
-
-            object
-          })
-
 setMethod("+", c("Dirac","DiscreteDistribution"),
       function(e1,e2){e2+location(e1)})
 
@@ -569,3 +506,6 @@
                             .lowerExact = .lowerExact(object), 
                             .logExact = .logExact(object)))}
                   )
+
+
+          
\ No newline at end of file

Modified: branches/distr-2.4/pkg/distr/R/LatticeDistribution.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/LatticeDistribution.R	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/LatticeDistribution.R	2011-09-01 01:35:31 UTC (rev 729)
@@ -12,12 +12,12 @@
         {  D <- DiscreteDistribution
            if (is(lattice, "Lattice")) 
              { ### check consistency with support of DiscreteDistribution} 
-              if( !.is.consistent(lattice, support(D), eq.space = FALSE)){         
-                 if (check)
+              if (check){
+                 if( !.is.consistent(lattice, support(D), eq.space = FALSE))         
                      stop(paste("Argument 'lattice' is inconsistent to",
                             " the support of argument 'DiscreteDistribution'." , 
                             sep = ""))
-                 else return(D)}           
+              }           
               return(new("AffLinLatticeDistribution", r = D at r, d = D at d, 
                           q = D at q, p = D at p, support = D at support, 
                           a = D at a, b = D at b, X0 = D at X0,
@@ -25,11 +25,11 @@
                           .withSim = .withSim, img = D at img,
                           param = D at param, Symmetry = Symmetry))
               }else{
-               if( !.is.vector.lattice(support(D))){ 
-                   if (check)
+               if (check){
+                   if( !.is.vector.lattice(support(D)))
                        stop(paste("Support of argument 'DiscreteDistribution' ",
                               "is not a lattice.", sep = ""))
-                   else return(D)}           
+               }           
                return(new("AffLinLatticeDistribution", r = D at r, d = D at d, 
                           q = D at q, p = D at p, support = D at support, 
                           lattice = .make.lattice.es.vector(D at support), 
@@ -44,23 +44,23 @@
         {  D <- DiscreteDistribution
            if (is(lattice, "Lattice")) 
              { ### check consistency with support of DiscreteDistribution} 
-              if( !.is.consistent(lattice, support(D), eq.space = FALSE)){         
-                 if (check)
+              if (check){
+                  if( !.is.consistent(lattice, support(D), eq.space = FALSE))         
                      stop(paste("Argument 'lattice' is inconsistent to the",
                             " support of argument 'DiscreteDistribution'." , 
                             sep = ""))
-                 else return(D)}           
+              }           
               return(new("LatticeDistribution", r = D at r, d = D at d, 
                           q = D at q, p = D at p, support = D at support, 
                           lattice = lattice, .withArith = .withArith, 
                           .withSim = .withSim, img = D at img,
                           param = D at param, Symmetry = Symmetry))
               }else{
-               if( !.is.vector.lattice(support(D))){ 
-                 if (check)
+               if (check){
+                   if( !.is.vector.lattice(support(D)))
                      stop(paste("Support of argument 'DiscreteDistribution' is",
                               "not a lattice.", sep = " "))
-                 else return(D)}           
+               }           
  
                return(new("LatticeDistribution", r = D at r, d = D at d, 
                           q = D at q, p = D at p, support = D at support, 
@@ -76,10 +76,10 @@
                                    .withArith = .withArith, 
                                    .withSim = .withSim, Symmetry = Symmetry )
         
-        if( !.is.consistent(lattice, supp, eq.space = FALSE)){         
-            if (check)     
+        if (check){
+            if( !.is.consistent(lattice, supp, eq.space = FALSE))         
                 stop("Argument 'lattice' is inconsistent to argument 'supp'.")
-            else return(D)}
+        }
         
         return(new("LatticeDistribution", r = r(D), d = d(D), 
                     q = q(D), p = p(D), support = supp, 
@@ -103,9 +103,10 @@
                           lattice = lattice, .withArith = .withArith, 
                           .withSim = .withSim, Symmetry = Symmetry))
                   }else{ 
-                   if (check)
+                   #if (check)
                        stop("Lengths of lattice and probabilities differ.")
-                   else return(D)}    
+                   #else return(D)
+                   }    
               }else {if (is.null(prob))
                         stop(paste("Insufficient information given to ",
                                    "determine distribution.", sep = ""))
@@ -126,18 +127,16 @@
             {if (is.null(prob)) prob <- supp*0+1/length(supp)
              D <- DiscreteDistribution(supp, prob, .withArith = .withArith, 
                                        .withSim = .withSim, Symmetry = Symmetry )
-             if (!.is.vector.lattice (supp)){
-                 if (check)
+             if (check){
+                 if (!.is.vector.lattice (supp))
                      stop("Argument 'supp' given is not a lattice.")
-                 else return (D)    
-             }else{ 
-                  return(new("LatticeDistribution", r = D at r, d = D at d, 
+             }    
+             return(new("LatticeDistribution", r = D at r, d = D at d, 
                              q = D at q, p = D at p, support = D at support, 
                              lattice = .make.lattice.es.vector(D at support), 
                              .withArith = D at .withArith, 
                              .withSim = D at .withSim, img = D at img,
                              param = D at param, Symmetry = Symmetry))                           
-                 }
             }else 
              stop("Insufficient information given to determine distribution.")
 }
@@ -183,79 +182,100 @@
 
 setMethod("+", c("LatticeDistribution", "LatticeDistribution"),
 function(e1,e2){
-            ### Step 1
+            if(length(support(e1))==1) return(e2+support(e1))
+            if(length(support(e2))==1) return(e1+support(e2))
 
-#            e1 <- as(e1, "LatticeDistribution")
-#            e2 <- as(e2, "LatticeDistribution")
-#                  ### casting necessary due to setIs
+### Lattice calculations:
 
-            ### Lattice Calculations:
-            w1 <- width(lattice(e1))
-            w2 <- width(lattice(e2))
             sup1 <- support(e1)
             sup2 <- support(e2)
-            maxl <- length(sup1)*length(sup2)
-                              ### length of product grid
-            csup <- unique(sort(c(sup1,sup2)))
-                              ### grid width of convolution grid
+            # left and right endpoint of convolution support
+            su12.l <- sup1[1]+sup2[1]
+            su12.r <- (rev(sup1))[1]+(rev(sup2))[1]
 
-            w <- min(diff(csup))
-            commonsup <- unique(sort(c(outer(sup1,sup2,"+"))))
-                              ### grid width of convolution grid
-            dcs <- abs(diff(commonsup))
-            mw <- min(dcs[dcs>getdistrOption("DistrResolution")])
-            if (abs(abs(w1)-abs(w2)) > getdistrOption("DistrResolution")){
-                   W <- sort(abs(c(w1,w2)))
-                   if (W[2] %% W[1] > getdistrOption("DistrResolution")){
+            l1 <- length(sup1)
+            l2 <- length(sup2)
 
-                         ## check whether arrangement on common grid really
-                         ## saves something
+            lat1 <- lattice(e1)
+            lat2 <- lattice(e2)
+            L1 <- Length(lat1)
+            L2 <- Length(lat2)
+            w1 <- width(lat1)
+            w2 <- width(lat2)
 
-                         prob1 <- d(e1)(sup1)
-                         prob2 <- d(e2)(sup2)
-                              ###  convolutional grid
-                         comsup <- seq(min(commonsup),max(commonsup), by = mw)
 
-                         fct <- function(sup0, prob0, bw){
-                              ### expand original grid,prob onto new width:
-                                    sup00 <- seq(min(sup0), max(sup0), by = mw)
-                                    prb0 <- 0 * sup00
-                                    ind0 <- sup00 %in% sup0
-                                    prb0[ind0] <- prob0
-                                    return(LatticeDistribution(supp = sup00,
-                                                               prob = prb0))
-                                    }
-                        if(length(comsup) < maxl)
-                           return( fct(sup1,prob1,mw)  + fct(sup2,prob2,mw))
-                        else
-                           return(as(e1, "DiscreteDistribution") +
-                                  as(e2, "DiscreteDistribution"))
-                   }
-                   else
-                       w <- mw #generate common lattice / support
+            ### take care if lattice is infinite
+            L.inf <- !(is.finite(L1)&&is.finite(L2))
+            if(L.inf){
+               if(is.finite(L2)){
+                  if(w1>0)
+                     L.lr <- +1
+                  else
+                     L.lr <- -1   
+               }else{    
+                  if(is.finite(L1)){
+                     if(w2>0)
+                        L.lr <- +1
+                     else
+                        L.lr <- -1   
+                  }else{
+                     if(w1*w2>0) L.lr <- if(w1>0) +1 else -1
+                     if(w1*w2<0) L.lr <- if(abs(su12.l)<abs(su12.r)) +1 else -1
                   }
+               }
+            }   
 
-            newlat <- NULL
-            ### if any lattice is infinite: see if we can keep this in mind:
-            if( ! (is.finite(Length(lattice(e1))) &&
-                   is.finite(Length(lattice(e2)))   )  &&
-                width(lattice(e1)) * width(lattice(e2)) > 0 )
-                {p1 <- pivot(lattice(e1))
-                 p2 <- pivot(lattice(e2))
-                 p <- p1 + p2
-                 w0 <- if (width(lattice(e1))>0) w else -w
-                 newlat <- Lattice(pivot = p, width = w0, Length = Inf)
-                }
-            ### end Lattice Calculations
 
+            e0 <- NULL
+            tol0 <- .distroptions$DistrResolution/1000
+            
+            ## treat case separately when Discr + Discr is "faster"
+            if(l1*l2 < 100){ 
+                     d0 <- .convDiscrDiscr(e1,e2)
+                     sup0 <- support(d0)
+                     md <- min(diff(sup0))
+                     sup00 <- seq(from=min(sup0),to=max(sup0),by=md)
+                     sup0s <- intersect(sup00,sup0)
+                     sup01 <- .inWithTol(sup00, sup0s, tol=tol0)
+                     sup10 <- .inWithTol(sup0, sup0s, tol=tol0)
+                     if(!all(sup10)) return(d0)
+                     pr0 <- sup00*0
+                     pr0[sup01] <- (prob(d0))[sup10]
+                     pr0 <- pr0/sum(pr0)                              
+                     lat <- Lattice(pivot = sup00[1], width = md, 
+                                    Length = length(sup00))
+                     e0 <- LatticeDistribution(supp = sup00, prob = pr0, 
+                                               lattice = lat, check = FALSE)           
+                     if(L.inf){
+                         e0 at lattice <- if(L.lr>0){ 
+                            Lattice(pivot = su12.l, width = wa, Length = Inf) 
+                                }else{ 
+                            Lattice(pivot = su12.r, width = -wa, Length = Inf)}
+                     }       
+               }
+ 
+            ## step 1 common width
+            wa <- .getCommonWidth(abs(w1),abs(w2),
+                      tol=tol0)
+
+
+            ## treat case separately when no common support, i.e. when 
+            ## w1/w2 is not "rational" enough
+            
+            if(is.null(wa))  return(.convDiscrDiscr(e1,e2))
+            
+            
+            w0 <- ifelse(w1<0,-1,1) * wa
+            pi1 <- pivot(lat1)
+            pi2 <- pivot(lat2)                        
             ### Step 2
-            supp0 <- seq(by = abs(w),
-                         from = min(support(e1), support(e2)),
-                         to   = max(support(e1), support(e2)))
+            supp0 <- seq(by = wa, from = min(sup1-pi1, sup2-pi2), to = max(sup1-pi1, sup2-pi2))
+            s1 <- .inWithTol(supp0,sup1-pi1,tol0)
+            s2 <- .inWithTol(supp0,sup2-pi2,tol0)
+            d1 <- d2 <- 0*supp0
+            d1[s1] <- prob(as(e1,"DiscreteDistribution"))
+            d2[s2] <- prob(as(e2,"DiscreteDistribution"))
 
-            d1 <- d(e1)(supp0)
-            d2 <- d(e2)(supp0)
-
             L <- length(supp0)
             Ln <- 2^(ceiling(log(L)/log(2))+1)
 
@@ -270,13 +290,14 @@
 
             ## convolution theorem for DFTs
             newd <- (Re(fft(ftde1*ftde2, inverse = TRUE)) / Ln)[1:(2*L+1)]
-            newd <- (newd >= .Machine$double.eps)*newd
+            newd <- (newd >= .Machine$double.eps^1.5)*newd
 
 
             ## reduction to relevant support
-            supp1 <- seq(by = abs(w),
-                         from = 2 * min(support(e1), support(e2)),
-                         to   = 2 * max(support(e1), support(e2)))
+            supp1 <- seq(by = wa,
+                         from = 2 * min(sup1-pi1, sup2-pi2),
+                         to   = 2 * max(sup1-pi1, sup2-pi2))+pi1+pi2
+
             L1 <- length(supp1)
             newd <- newd[1:L1]
 
@@ -287,9 +308,6 @@
                            )
                 rsum.l <- 1 + sum( cumsum(newd) <
                                    getdistrOption("TruncQuantile")/2)
-                newd <- newd[rsum.l:rsum.u]
-                newd <- newd/sum(newd)
-                supp1 <- supp1[rsum.l:rsum.u]
             }else{
                 rsum.u <- min( sum( rev(cumsum(rev(newd))) >=
                                     .Machine$double.eps),
@@ -297,12 +315,15 @@
                            )
                 rsum.l <- 1 + sum( cumsum(newd) < .Machine$double.eps)
 
-                newd <- newd[rsum.l:rsum.u]
-                newd <- newd/sum(newd)
-                supp1 <- supp1[rsum.l:rsum.u]
             }
-            supp2 <- supp1[newd > getdistrOption("TruncQuantile")]
-            newd2 <- newd[newd  > getdistrOption("TruncQuantile")]
+            wi1 <- rsum.l:rsum.u
+            newd <- newd[wi1]
+            newd <- newd/sum(newd)
+            supp1 <- supp1[wi1]
+
+            wi2 <- newd > getdistrOption("TruncQuantile")
+            supp2 <- supp1[wi2]
+            newd2 <- newd[wi2]
             newd2 <- newd2/sum(newd2)
 
             Symmetry <- NoSymmetry()
@@ -311,13 +332,24 @@
                Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
                                               SymmCenter(e2 at Symmetry))   
 
-            
-            if( length(supp1) >= 2 * length(supp2))
+            if( length(supp1) >= 2 * length(supp2)){
                return(DiscreteDistribution(supp = supp2, prob = newd2,
                                            .withArith = TRUE, Symmetry = Symmetry))
-            else
-               return(LatticeDistribution(supp = supp1, prob = newd,
-                                          .withArith = TRUE, Symmetry = Symmetry))
+            }else{
+               lat <- Lattice(pivot=supp1[1],width=wa, Length=length(supp1))
+
+               e0 <- LatticeDistribution(supp = supp1, prob = newd,
+                                         lattice = lat,
+                                         .withArith = TRUE, Symmetry = Symmetry,
+                                         check = FALSE)
+               if(L.inf){
+                  e0 at lattice <- if(L.lr>0){ 
+                            Lattice(pivot = su12.l, width = wa, Length = Inf) 
+                         }else{ 
+                            Lattice(pivot = su12.r, width = -wa, Length = Inf)}
+               }
+               return(e0)
+            }
           })
 
 ## extra methods
@@ -336,7 +368,8 @@
                  Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+e2)   
               
               LatticeDistribution(lattice = L, 
-                     DiscreteDistribution = Distr, Symmetry = Symmetry)                                        
+                     DiscreteDistribution = Distr, Symmetry = Symmetry, 
+                     check = FALSE)                                        
               })       
 
 setMethod("*", c("LatticeDistribution", "numeric"),
@@ -356,7 +389,8 @@
                      Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry) * e2)   
               
                   return(LatticeDistribution(lattice = L, 
-                          DiscreteDistribution = Distr, Symmetry = Symmetry))
+                          DiscreteDistribution = Distr, Symmetry = Symmetry, 
+                          check = FALSE))
                 }
              }
           )              
@@ -371,7 +405,8 @@
               LatticeDistribution(lattice = L, 
                      DiscreteDistribution = 
                         as(e1, "AffLinDiscreteDistribution") + e2,
-                        Symmetry = Symmetry)                     
+                        Symmetry = Symmetry, 
+                     check = FALSE)                     
               })       
 
 setMethod("*", c("AffLinLatticeDistribution", "numeric"),
@@ -388,7 +423,8 @@
                   return(LatticeDistribution(lattice = L, 
                           DiscreteDistribution = 
                              as(e1, "AffLinDiscreteDistribution") * 
-                             e2, Symmetry = Symmetry))
+                             e2, Symmetry = Symmetry, 
+                             check = FALSE))
                 }
              }
           )              

Modified: branches/distr-2.4/pkg/distr/R/internalUtils.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/internalUtils.R	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/internalUtils.R	2011-09-01 01:35:31 UTC (rev 729)
@@ -48,6 +48,114 @@
 }
 
 #------------------------------------------------------------------------------
+### internal help function to determin common lattice width
+#------------------------------------------------------------------------------
+
+.EuclidAlgo <- function(n1,n2){
+   r<- 2
+   m <- round(max(n1,n2))
+   n <- round(min(n1,n2))
+   if (n==1) return(1)
+   while (r>1){
+      r <- m %% n
+      m <- n
+      if(r==1) return(1)
+      if(r==0) break
+      n <- r
+      }
+   return(n)
+}
+.getCommonWidth <- function(x1,x2, tol=.Machine$double.eps){
+  rI <- function(x) .isInteger(x, tol = x * tol)
+  if(rI(x1)&&rI(x2)) return(.EuclidAlgo(x1,x2))
+  if(rI(10000*x1) && rI(10000*x2)) return(.EuclidAlgo(10000*x1,10000*x2)/10000)
+  n <- max(x1,x2); m <- min(x1,x2)
+  if(rI(n/m)) return(m)
+  vc <- 1:10000
+  vecT <- sapply(vc,function(x) rI(x*(n/m)))
+  vec <- m/vc
+#  print(min(vc[vecT]))
+  if(any(vecT)) return((vec/min(vc[vecT]))[1])
+  return(NULL)  
+}
+
+#------------------------------------------------------------------------------
+### %in% for numerics with tolerance
+#------------------------------------------------------------------------------
+.inWithTol <- function(x,y,tol=.Machine$double.eps){
+   sapply(x, function(u) any(abs(u-y)<tol))
+}                                                                         
+
+#------------------------------------------------------------------------------
+### brute force convolution
+#------------------------------------------------------------------------------
+.convDiscrDiscr <- function(e1,e2){
+            convolutedsupport <- rep(support(e1), each = length(support(e2))) +
+                                 support(e2)
+
+            gridvalues1 <- d(e1)(support(e1)); gridvalues2 <- d(e2)(support(e2))
+            convolutedvalues <- rep(gridvalues1, each = length(support(e2))) *
+                                gridvalues2
+            rm(gridvalues1,gridvalues2)
+
+            tmptable <- data.frame(x = convolutedsupport, dx = convolutedvalues)
+            rm(convolutedsupport,convolutedvalues)
+            tmp <- tapply(tmptable$dx, tmptable$x, sum)
+            rm(tmptable)
+
+            supp.u <- as.numeric(names(tmp))
+            prob.u <- as.numeric(tmp)
+
+            o <- order(supp.u)
+            supp <- supp.u[o]
+            prob <- prob.u[o]
+
+            #supp.u <- unique(supp)
+
+            len <- length(supp)
+
+            if(len > 1){
+               if (min(diff(supp))< getdistrOption("DistrResolution")){
+                   if (getdistrOption("DistrCollapse")){
+                       erg <- .DistrCollapse(supp, prob, 
+                                   getdistrOption("DistrResolution"))
+                       if ( len > length(erg$prob) && 
+                                getdistrOption("DistrCollapse.Unique.Warn") )
+                            warning("collapsing to unique support values")         
+                       prob <- erg$prob
+                       supp <- erg$supp
+                   }else
+                    stop("grid too narrow --> change DistrResolution")
+                }
+            }
+
+            rm(tmp, len)
+
+            .withSim <- e1 at .withSim || e2 at .withSim
+
+            rfun <- function(n) {}
+            body(rfun) <- substitute({ f(n) + g(n) },
+                                         list(f = e1 at r, g = e2 at r))
+
+            dfun <- .makeDNew(supp, prob, Cont = FALSE)
+            pfun <- .makePNew(supp, prob, .withSim, Cont = FALSE)
+            qfun <- .makeQNew(supp, cumsum(prob), rev(cumsum(rev(prob))),
+                      .withSim, min(supp), max(supp), Cont = FALSE)
+
+            object <- new("DiscreteDistribution", r = rfun, d = dfun, p = pfun,
+                           q = qfun, support = supp,
+                           .withSim = .withSim, .withArith = TRUE)
+            rm(rfun, dfun, qfun, pfun)
+
+            if(is(e1 at Symmetry,"SphericalSymmetry")&& 
+               is(e2 at Symmetry,"SphericalSymmetry"))
+               object at Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
+                                                     SymmCenter(e2 at Symmetry))   
+
+            object
+          }
+
+#------------------------------------------------------------------------------
 ### .fm, .fM, .fM2 functions
 #------------------------------------------------------------------------------
 

Modified: branches/distr-2.4/pkg/distr/man/internals.Rd
===================================================================
--- branches/distr-2.4/pkg/distr/man/internals.Rd	2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/man/internals.Rd	2011-09-01 01:35:31 UTC (rev 729)
@@ -17,6 +17,8 @@
 \alias{.make.lattice.es.vector}
 \alias{.presubs}
 \alias{.inArgs}
+\alias{.EuclidAlgo}
+\alias{.getCommonWidth}
 \alias{.isEqual}
 \alias{.isIn}
 \alias{.setEqual}
@@ -55,6 +57,8 @@
 \alias{.trunc.low}
 \alias{.modifyqgaps}
 \alias{.DistrCollapse}
+\alias{.convDiscrDiscr}
+\alias{.inWithTol}
 \alias{devNew}
 
 \title{Internal functions of package distr}
@@ -125,6 +129,10 @@
 .trunc.low(object, lower)
 .modifyqgaps(pfun, qfun, gaps, leftright = "left")
 .DistrCollapse(support, prob, eps = getdistrOption("DistrResolution"))
+.EuclidAlgo(n1,n2)
+.getCommonWidth(x1,x2, tol=.Machine$double.eps)
+.convDiscrDiscr(e1,e2)
+.inWithTol(x,y,tol=.Machine$double.eps)
 devNew(...)
 }
 
@@ -228,6 +236,10 @@
                   the left continuous version; for slot \code{p}: if partially 
                   matched to \code{"left"} the left continuous version, else
                   the right continuous version;}
+ \item{n1}{integer argument for \code{.EuclidAlgo}}
+ \item{n2}{integer argument for \code{.EuclidAlgo}}
+ \item{x1}{width argument for \code{.getCommonWidth}}
+ \item{x2}{width argument for \code{.getCommonWidth}}
  \item{...}{arguments passed through to other functions}
 }
 
@@ -366,7 +378,14 @@
      \code{gaps} is not \code{NULL}. If argument \code{leftright} does not
      partially match \code{"right"} (default) returns the left continuous
      version of the quantile function, else the right continuous one.
-
+ 
+\code{.EuclidAlgo} computes the greatest common divisor of two integers by 
+      means of the Euclidean algorithm.
+\code{.getCommonWidth} for two lattices with widths \code{x1} and \code{x2}
+      computes the smallest common lattice width for convolution.
+\code{.convDiscrDiscr} computes the convolution of two discrete distributions by
+    brute force.
+\code{.inWithTol} works like \code{\%in\%} but with a given tolerance.
 \code{devNew} opens a new device. This function is for back compatibility
 with R versions < 2.8.0.
 }
@@ -433,6 +452,11 @@
        Criterium for collapsing: a distance smaller than argument
        \code{eps}.
        }
+\item{.EuclidAlgo}{returns the greatest common divisor (an integer).}
+\item{.getCommonWidth}{returns the  smallest common lattice width (a numeric).}
+\item{.convDiscrDiscr}{returns the convolution of two discrete distributions.}
+\item{.inWithTol}{returns a logical vector of same lenght as \code{x} for the
+    matches (up to tolerance) with vector \code{y}.}
 \item{devNew}{returns the return value of the device opened, 
 usually invisible \code{NULL}.}
 }

Deleted: branches/distr-2.4/pkg/distrEx/src/distrEx.dll
===================================================================
(Binary files differ)

Modified: pkg/distr/DESCRIPTION
===================================================================
--- pkg/distr/DESCRIPTION	2011-07-08 13:13:36 UTC (rev 728)
+++ pkg/distr/DESCRIPTION	2011-09-01 01:35:31 UTC (rev 729)
@@ -1,6 +1,6 @@
 Package: distr
-Version: 2.3.2
-Date: 2011-03-22
+Version: 2.3.3
+Date: 2011-08-31
 Title: Object oriented implementation of distributions
 Description: S4 Classes and Methods for distributions
 Author: Florian Camphausen, Matthias Kohl, Peter Ruckdeschel, Thomas Stabla

Modified: pkg/distr/R/AllClasses.R
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/distr -r 729


More information about the Distr-commits mailing list