isoton <-
function(x, wt = -1)
{
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# Performs a monotone smoothing on vector x, optionally using weights
# in parallel vector wt
#----------------------------------------------------------------------
        
        if(!is.numeric(x) || length(x) < 2)
                stop("x must be numeric  of length > 1")
        if(missing(wt)) {
                z <- .Fortran("monosm",
                        ans = as.double(x),
                        as.double(1),
                        as.integer(length(x)),
                        as.logical(FALSE),
                        integer(length(x)))
        }
        else {
                if(!is.numeric(wt) || length(wt) != length(x))
                        stop("wt must be numeric length same as x")
                z <- .Fortran("monosm",
                        ans = as.double(x),
                        as.double(wt),
                        as.integer(length(x)),
                        as.logical(TRUE),
                        integer(length(x)))
        }
        return(z$ans)
}


stepfun <-
function(datax, datay, type = "left")
{
# augment a set of points so that it plots as a left-continuous function
# allow both (x,y) and (structure with $x $y)  input
        if(missing(datay)) x <- datax$x else x <- datax
        if(missing(datay))
                y <- datax$y
        else y <- datay
        n <- length(x)
        type <- charmatch(type, c("left", "right"))
        if(is.na(type))
                stop("The type must be 'left' or 'right' continuous")
        if(any(diff(x) < 0))
                stop("The x vector must be sorted")
        if(type == 2) {
                x <- rev(x)
                y <- rev(y)
        }
        if(n > 2) {
# remove redundant points
                dupy <- c(T, diff(y[ - n]) != 0, T)
                dupx <- c(T, diff(x[ - n]) != 0, T)
                x <- x[dupx & dupy]
                y <- y[dupx & dupy]
                n <- length(x)
        }
#create the step function
        xrep <- rep(x[2:n], rep(2, n - 1))
        yrep <- rep(y[1:(n - 1)], rep(2, n - 1))
        if(type == 1)
                list(x = c(x[1], xrep), y = c(yrep, y[n]))
        else list(x = c(rev(xrep), x[1]), y = c(y[n], rev(yrep)))
}


dlwb <- 
function(x, y, rw = c(rep(1, length(y))), f = 1/3)
{
        nk <- order(x)
        y <- y[nk]
        x <- x[nk]
        n <- length(y)
        fit <- .Fortran("dlwb",
                x = as.double(x),
                y = as.double(y),
                as.integer(n),
                as.double(f),
                 ys = double(n),
                rw = as.double(rw),
                work = double(n))
        list(x = fit$x, ys = fit$ys)
}