.packageName <- "DNAcopy"
CNA <- function(genomdat, chrom, maploc, data.type=c("logratio","binary"),
                   sampleid=NULL)
  {
    if (!is.numeric(genomdat)) stop("genomdat must be numeric")
    if (is.factor(chrom)) chrom <- as.character(chrom)
    if (!is.numeric(maploc)) stop("maploc must be numeric")
    data.type <- match.arg(data.type)
    if (sum(is.na(chrom)|is.na(maploc))>0)
      warning("markers with missing chrom and/or maploc removed\n")
    sortindex <- order(chrom, maploc, na.last=NA)
    if (is.vector(genomdat)) genomdat <- as.matrix(genomdat)
    if (!missing(sampleid)) {
      if (length(sampleid) != ncol(genomdat)) {
        warning("length(sampleid) and ncol(genomdat) differ, names ignored\n")
        sampleid <- paste("Sample", 1:ncol(genomdat))
      } 
    } else {
        sampleid <- paste("Sample", 1:ncol(genomdat))
    }
    colnames(genomdat) <- sampleid
    zzz <- data.frame(chrom=I(chrom), maploc=maploc, genomdat)
    zzz <- zzz[sortindex,]

    if(!all(sapply(unique(chrom),function(ichrom,chrom,maploc){
      length(maploc[chrom==ichrom])-length(unique(maploc[chrom==ichrom]))
    },chrom,maploc)==0)) {
      warning("array has repeated maploc positions\n")
    } else {
      rownames(zzz) <- paste(zzz$chrom,zzz$maploc,sep="-")
    }      
    attr(zzz, "data.type") <- data.type
    class(zzz) <- c("CNA","data.frame")
    zzz
  }

subset.CNA <- function(x, chromlist=NULL, samplelist=NULL, ...)
  {
    if (!inherits(x, 'CNA')) stop("First arg must be of class CNA")
    chrom <- x$chrom
    uchrom <- unique(chrom)
    if (missing(chromlist)) chromlist <- uchrom
    if (length(setdiff(chromlist, uchrom)) > 0)
      stop("chromlist contains chromosomes not in the data")
    if (length(chromlist) > length(unique(chromlist)))
      warning("duplicate chromosomes in chromlist removed")
    sampleid <- colnames(x)[-(1:2)]
    if (missing(samplelist)) samplelist <- sampleid
    nsample <- length(sampleid)
    if (length(setdiff(samplelist, 1:nsample)) > 0 & length(setdiff(samplelist, sampleid)) > 0)
      stop("samplelist should be a list of valid sample numbers or names")
    if (!is.integer(samplelist)) samplelist <- match(samplelist, names(x)) - 2
    if (length(samplelist) > length(unique(samplelist)))
      warning("duplicate samples in samplelist removed")
    samplelist <- unique(samplelist)
    y <- x[chrom %in% chromlist,c(1:2,samplelist+2)]
    attr(y, "data.type") <- attr(x, "data.type")
    y
  }

smooth.CNA <- function(x, smooth.region=2, outlier.SD.scale=4,
                       smooth.SD.scale=2, trim=0.025)
  {
    if (!inherits(x, 'CNA')) stop("First arg must be of class CNA")
    nsample <- ncol(x)-2
    chrom <- x$chrom
    uchrom <- unique(chrom)
    if(attr(x, "data.type")=="binary") stop("Not smoothing binary data ")
    for (isamp in 1:nsample) {
      genomdat <- x[,isamp+2]
      ina <- which(!is.na(genomdat) & !(abs(genomdat)==Inf))
      trimmed.SD <- sqrt(trimmed.variance(genomdat[ina], trim))
      outlier.SD <- outlier.SD.scale*trimmed.SD
      smooth.SD <- smooth.SD.scale*trimmed.SD
      k <- smooth.region
      for (i in uchrom) {
        ina <- which(!is.na(genomdat) & !(abs(genomdat)==Inf) & chrom==i)
        n <- length(genomdat[ina])
        smoothed.data <- sapply(1:n, function(i, x, n, nbhd, oSD, sSD) {
          xi <- x[i]
          nbhd <- i+nbhd
          xnbhd <- x[nbhd[nbhd>0 & nbhd <=n]]
          if (xi > max(xnbhd) + oSD) xi <- median(c(xi,xnbhd)) + sSD
          if (xi < min(xnbhd) - oSD) xi <- median(c(xi,xnbhd)) - sSD
          xi
        }, genomdat[ina], n, c(-k:-1, 1:k), outlier.SD, smooth.SD)
        x[,isamp+2][ina] <- smoothed.data
      }
    }
    x
  }

print.CNA <- function(x, ...)
  {
    if (!inherits(x, 'CNA')) stop("First arg must be of class CNA")
    cat("Number of Samples", ncol(x)-2,
        "\nNumber of Probes ", nrow(x),
        "\nData Type        ", attr(x,"data.type"),"\n")
  }

plot.DNAcopy <- function (x, plot.type=c("whole", "plateau", "samplebychrom",
                               "chrombysample"), xmaploc=FALSE, altcol=TRUE,
                          sbyc.layout=NULL, cbys.nchrom=1, cbys.layout=NULL,
                          include.means=TRUE, zeroline=TRUE, pt.pch=NULL,
                          pt.cex=NULL, pt.cols=NULL, segcol=NULL, zlcol=NULL,
                          ylim=NULL, lwd=NULL, ...)
{
  if (!inherits(x, "DNAcopy")) 
    stop("First arg must be the result of segment")
  xdat <- x$data
  nsample <- ncol(xdat)-2
  if(missing(ylim)) {
    uylim <- max(abs(xdat[,-(1:2)]), na.rm=TRUE)
    ylim <- c(-uylim, uylim)
  }
  xres <- x$output
  if(dev.cur() <= 1) get(getOption("device"))()
  int.dev <- dev.interactive()
  plot.type <- match.arg(plot.type)
  op <- par(no.readonly = TRUE)
  parask <- par("ask")
  if (int.dev & !parask & nsample>1) par(ask = TRUE)
  sampleid <- colnames(xdat)[-(1:2)]
  chrom0 <- xdat$chrom
  uchrom <- unique(chrom0)
  nchrom <- length(uchrom)
  if (xmaploc) {
    maploc0 <- as.numeric(xdat$maploc)
    if(length(uchrom)>1 & max(maploc0[chrom0==uchrom[1]]) > min(maploc0[chrom0==uchrom[2]])) {
      plen <- max(maploc0[chrom0==uchrom[1]])
      for(i in 2:nchrom) {
        maploc0[chrom0==uchrom[i]] <- plen + maploc0[chrom0==uchrom[i]]
        plen <- max(maploc0[chrom0==uchrom[i]])
      }
    }
  }
  if (missing(pt.pch)) pt.pch <- "."
  if (missing(pt.cex)) {
    if (pt.pch==".") { pt.cex <- 3}
    else {pt.cex <- 1}
  }
  wcol0 <- rep(1, length(chrom0))
  if (altcol) {
    j <- 0
    for (i in uchrom) {
      j <- (j+1) %% 2
      wcol0[chrom0==i] <- 1+j
    }
  }
  if (missing(pt.cols)) pt.cols <- c("black","green")
  if (missing(segcol)) segcol <- "red"
  if (missing(zlcol)) zlcol <- "grey"
  if (missing(lwd)) lwd <- 3
  if (plot.type == "chrombysample") {
    cat("Setting multi-figure configuration\n")
    par(mar = c(0, 4, 0, 2), oma = c(4, 0, 4, 0), mgp = c(2, 0.7, 0))
    if (missing(cbys.layout)) {
      nrow <- ncol <- ceiling(sqrt(nsample))
      if (nrow*ncol - nsample > 0) {
        nrow <- nrow - 1
        ncol <- ncol + 1
      }
      if (nrow*ncol - nsample >= nrow) ncol <- ncol - 1
      cbys.layout <- c(nrow, ncol)
    }
    lmat0 <- lmat1 <- c(1:nsample, rep(-cbys.nchrom*nsample, prod(cbys.layout) - nsample))
    for(i in 1:(cbys.nchrom-1)) {
      lmat1 <- c(lmat1,lmat0+nsample*i)
    }
    lmat1[lmat1<0] <- 0
    lmat <- matrix(lmat1, nrow = cbys.layout[1], ncol = cbys.nchrom*cbys.layout[2], byrow = FALSE)
    layout(lmat)
  }
  if (plot.type == "samplebychrom") {
    cat("Setting multi-figure configuration\n")
    par(mar = c(4, 4, 4, 2), oma = c(0, 0, 2, 0), mgp = c(2, 0.7, 0))
    if (missing(sbyc.layout)) {
      nrow <- ncol <- ceiling(sqrt(nchrom))
      if (nrow*ncol - nchrom > 0) {
        nrow <- nrow - 1
        ncol <- ncol + 1
      }
      if (nrow*ncol - nchrom > ncol) ncol <- ncol - 1
      sbyc.layout <- c(nrow, ncol)
    }
    lmat <- matrix(c(1:nchrom, rep(0,prod(sbyc.layout)-nchrom)),
                   nrow = sbyc.layout[1], ncol = sbyc.layout[2], byrow=TRUE)
    layout(lmat)
  }
  if (plot.type == "chrombysample") {
    atchrom <- 0.5/cbys.nchrom
    for (ichrom in uchrom) {
      for (isamp in 1:nsample) {
        genomdat <- xdat[chrom0==ichrom, isamp+2]
        ina <- which(!is.na(genomdat) & !(abs(genomdat) == Inf))
        genomdat <- genomdat[ina]
        ii <- cumsum(c(0, xres$num.mark[xres$ID == sampleid[isamp] & xres$chrom==ichrom]))
        mm <- xres$seg.mean[xres$ID == sampleid[isamp] & xres$chrom==ichrom]
        kk <- length(ii)
        zz <- cbind(ii[-kk] + 1, ii[-1])
        plot(genomdat, pch = pt.pch, cex=pt.cex, xaxt="n", ylim = ylim, ylab = sampleid[isamp])
        if(zeroline) abline(h=0, col=zlcol, lwd=lwd)
        if (isamp%%cbys.layout[1] == 0) {
          axis(1, outer=TRUE)
          title(xlab="Index")
        }
        if (include.means) {
          for (i in 1:(kk - 1))
            {
              lines(zz[i, ], rep(mm[i], 2), col = segcol, lwd=lwd)
            }
        }
      }
      mtext(paste("Chromosome",ichrom), side = 3, line = 1, at = atchrom, outer=TRUE, font=2)
      atchrom <- atchrom + 1/cbys.nchrom
      atchrom <- atchrom - floor(atchrom)
    }
  } else {
    for (isamp in 1:nsample)
      {
        genomdat <- xdat[, isamp+2]
        ina <- which(!is.na(genomdat) & !(abs(genomdat) == Inf))
        genomdat <- genomdat[ina]
        wcol <- wcol0[ina]
        chrom <- chrom0[ina]
        if (xmaploc) maploc <- maploc0[ina]
        ii <- cumsum(c(0, xres$num.mark[xres$ID == sampleid[isamp]]))
        mm <- xres$seg.mean[xres$ID == sampleid[isamp]]
        kk <- length(ii)
        zz <- cbind(ii[-kk] + 1, ii[-1])
        if(missing(ylim)) ylim <- range(c(genomdat, -genomdat))
        if (plot.type=="whole")
          {
            if (xmaploc) {
              plot(maploc, genomdat, pch = pt.pch, cex=pt.cex, col=pt.cols[wcol], main = sampleid[isamp], ylab = "", ylim = ylim)
              if(zeroline) abline(h=0, col=zlcol, lwd=lwd)
            } else {
              plot(genomdat, pch = pt.pch, cex=pt.cex, col=pt.cols[wcol], main = sampleid[isamp], ylab = "", ylim = ylim)
              if(zeroline) abline(h=0, col=zlcol, lwd=lwd)
            }
            if (include.means) {
              for (i in 1:(kk - 1))
                {
                  if (xmaploc) { 
                    lines(maploc[zz[i, ]], rep(mm[i], 2), col = segcol, lwd=lwd)
                  } else {
                    lines(zz[i, ], rep(mm[i], 2), col = segcol, lwd=lwd)
                  }
                }
            }
          }
        if (plot.type=="samplebychrom")
          {
            cc <- xres$chrom[xres$ID == sampleid[isamp]]
            for (ichrom in uchrom)
              {
                plot(genomdat[chrom == ichrom], pch = pt.pch, cex=pt.cex, ylab = "", main = paste("Chromosome", ichrom), ylim = ylim)
                if(zeroline) abline(h=0, col=zlcol, lwd=lwd)
                if (include.means) {
                  jj <- which(cc==ichrom)
                  jj0 <- min(jj)
                  for (i in jj)
                    {
                      lines(1+zz[i, ]-zz[jj0,1], rep(mm[i], 2), col = segcol, lwd=lwd)
                    }
                }
              }
            mtext(sampleid[isamp], side = 3, line = 0, outer = TRUE, font=2)
          }
        if (plot.type=="plateau")
          {
            omm <- order(mm)
            ozz <- zz[omm,]
            ina <- unlist(apply(ozz, 1, function(ii) ii[1]:ii[2]))
            plot(genomdat[ina], pch = pt.pch, cex=pt.cex, main = sampleid[isamp], ylab = "", ylim = ylim)
            if(zeroline) abline(h=0, col=zlcol, lwd=lwd)
            if (include.means) {
              ii <- cumsum(c(0, xres$num.mark[xres$ID == sampleid[isamp]][omm]))
              smm <- mm[omm]
              zz <- cbind(ii[-kk] + 1, ii[-1])
              for (i in 1:(kk-1)) lines(zz[i, ], rep(smm[i], 2), col = segcol, lwd=lwd)
            }
          }
      }
  }
  on.exit( if (plot.type=="chrombysample" | plot.type=="samplebychrom") {
    par(op)
  } else { if(int.dev & !parask & nsample>1) par(ask=parask) })
}

print.DNAcopy <- function(x, ...)
  {
    if (!inherits(x, "DNAcopy")) stop("Object is not the result of segment")
    if (!is.null(cl<- x$call))
      {
        cat("Call:\n")
        dput(cl)
        cat("\n")
      }
    print(x$output)
  }

subset.DNAcopy <- function(x, chromlist=NULL, samplelist=NULL, ...)
  {
    if (!inherits(x, 'DNAcopy')) stop("First arg must be of class DNAcopy")
    zdat <- x$data
    zres <- x$output
    chrom <- zdat$chrom
    uchrom <- unique(chrom)
    if (missing(chromlist)) chromlist <- uchrom
    if (length(setdiff(chromlist, uchrom)) > 0)
      stop("chromlist contains chromosomes not in the data")
    if (length(chromlist) > length(unique(chromlist)))
      warning("duplicate chromosomes in chromlist removed")
    sampleid <- colnames(zdat)[-(1:2)]
    if (missing(samplelist)) samplelist <- sampleid
    nsample <- length(sampleid)
    if (length(setdiff(samplelist, 1:nsample)) > 0 & length(setdiff(samplelist, sampleid)) > 0)
      stop("samplelist should be a list of valid sample numbers or names")
    if (!is.integer(samplelist)) samplelist <- match(samplelist, names(zdat)) - 2
    if (length(samplelist) > length(unique(samplelist)))
      warning("duplicate samples in samplelist removed")
    samplelist <- unique(samplelist)
    jj <- unlist(sapply(sampleid[samplelist], function(i, id) {which(id==i)}, zres$ID ))
    zres <- zres[jj,]
    y <- list()
    y$data <- zdat[chrom %in% chromlist,c(1:2,samplelist+2)]
    attr(y$data, "data.type") <- attr(zdat, "data.type")
    y$output <- zres[zres$chrom %in% chromlist,]
    class(y) <- "DNAcopy"
    y
  }

# Chromosome.Lengths <- c(263, 255, 214, 203, 194, 183, 171, 155, 145, 144, 144, 143, 114, 109, 106, 98, 92, 85, 67, 72, 50, 56, 164, 59)
# names(Chromosome.Lengths) <- c(as.character(1:22),"X","Y")
changepoints <- function(genomdat, data.type="logratio", alpha=0.01, sbdry,
                         sbn, nperm=10000, p.method="hybrid", window.size=NULL,
                         overlap=0.25, kmax=25, nmin=200, trimmed.SD = NULL,
                         undo.splits="none", undo.prune=0.05, undo.SD=3,
                         verbose=1, ngrid=100, tol=1e-6)
  {
    n <- length(genomdat)
    if (missing(trimmed.SD)) trimmed.SD <- mad(diff(genomdat))/sqrt(2)
    if (is.null(window.size)) window.size <- n
    seg.end <- c(0,n)
    k <- length(seg.end)
    change.loc <- NULL
    while (k > 1)
      {
        current.n <- seg.end[k]-seg.end[k-1]
        if (verbose>=3) cat(".... current segment:",seg.end[k-1]+1,"-",seg.end[k],"\n")
        if(current.n >= 4)
          {
            current.genomdat <- genomdat[(seg.end[k-1]+1):seg.end[k]]
            wsize <- min(current.n, window.size)
            winnum <- ceiling((current.n - wsize)/((1-overlap)*wsize)) + 1
            winloc <- round(seq(0,current.n - wsize,length=winnum))
            hybrid <- FALSE
            delta <- 0
            if ((p.method=="hybrid") & (nmin < current.n)) {
              hybrid <- TRUE
              delta <- (kmax+1)/current.n
            }
            zzz <- .Fortran("fndcpt",
                            n=as.integer(current.n),
                            w=as.integer(wsize),
                            twon=as.integer(2*wsize),
                            wn=as.integer(winnum),
                            wloc=as.integer(winloc),
                            x=as.double(current.genomdat),
                            px=double(current.n),
                            sx=double(wsize),
                            tx=double(2*wsize),
                            nperm=as.integer(nperm),
                            cpval=as.double(alpha),
                            ncpt=integer(1),
                            icpt=integer(2),
                            ibin=as.logical(data.type=="binary"),
                            hybrid=as.logical(hybrid),
                            hk=as.integer(kmax),
                            delta=as.double(delta),
                            ngrid=as.integer(ngrid),
                            sbn=as.integer(sbn),
                            sbdry=as.integer(sbdry),
                            tol= as.double(tol),
                            PACKAGE="DNAcopy")
          }
        else
          {
            zzz <- list()
            zzz$ncpt <- 0
          }
        if(zzz$ncpt==0) change.loc <- c(change.loc,seg.end[k])
        seg.end <- switch(1+zzz$ncpt,seg.end[-k],
                          c(seg.end[1:(k-1)],seg.end[k-1]+zzz$icpt[1],seg.end[k]),
                          c(seg.end[1:(k-1)],seg.end[k-1]+zzz$icpt,seg.end[k]))
        k <- length(seg.end)
        if(verbose>=3) cat(".... segments to go:",seg.end,"\n")
      }
    seg.ends <- rev(change.loc)
    nseg <- length(seg.ends)
    lseg <- diff(c(0,seg.ends))
    if (nseg > 1) {
        if (undo.splits == "prune") {
            lseg <- changepoints.prune(genomdat, lseg, undo.prune)
        }
        if (undo.splits == "sdundo") {
            lseg <- changepoints.sdundo(genomdat, lseg, trimmed.SD, undo.SD)
        }
    }
    segmeans <- 0*lseg
    ll <- uu <- 0
    for(i in 1:length(lseg)) {
      uu <- uu + lseg[i]
      segmeans[i] <- mean(genomdat[(ll+1):uu])
      ll <- uu
    }
    list("lseg" = lseg, "segmeans" = segmeans)
  }

changepoints.prune <- function(genomdat, lseg, change.cutoff=0.05) {
  n <- length(genomdat)
  nseg <- length(lseg)
  ncpt <- nseg-1
  zzz <- .Fortran("prune",
                  as.integer(n),
                  as.double(genomdat),
                  as.integer(nseg),
                  as.integer(lseg),
                  as.double(change.cutoff),
                  double(nseg),
                  as.integer(ncpt),
                  loc=integer(ncpt),
                  integer(2*ncpt),
                  pncpt=integer(1), PACKAGE="DNAcopy")
  pruned.ncpt <- zzz$pncpt
  pruned.cpts <- cumsum(lseg)[zzz$loc[1:pruned.ncpt]]
  pruned.lseg <- diff(c(0,pruned.cpts,n))
  pruned.lseg
}

changepoints.sdundo <- function(genomdat, lseg, trimmed.SD, change.SD=3) {
  change.SD <- trimmed.SD*change.SD
  cpt.loc <- cumsum(lseg)
  sdundo <- TRUE
  while(sdundo) {
    k <- length(cpt.loc)
    if (k>1) {
      segments0 <- cbind(c(1,1+cpt.loc[-k]),cpt.loc)
      segmed <- apply(segments0, 1, function(i,x) {median(x[i[1]:i[2]])}, genomdat)
      adsegmed <- abs(diff(segmed))
      if (min(adsegmed) < change.SD) {
        i <- which(adsegmed == min(adsegmed))
        cpt.loc <- cpt.loc[-i]
      } else {
        sdundo <- FALSE
      }
    } else {
      sdundo <- FALSE
    }
  }
  lseg.sdundo <- diff(c(0,cpt.loc))
  lseg.sdundo
}

trimmed.variance <- function(genomdat, trim=0.025)
  {
    n <- length(genomdat)
    n.keep <- round((1-2*trim)*(n-1))
    inflfact(trim)*sum((sort(abs(diff(genomdat)))[1:n.keep])^2 / (2*n.keep))
  }

inflfact <- function(trim)
  {
    a <- qnorm(1-trim)
    x <- seq(-a,a,length=10001)
    x1 <- (x[-10001] + x[-1])/2
    1/(sum(x1^2*dnorm(x1)/(1-2*trim))*(2*a/10000))
  }
getbdry <- function(eta, nperm, max.ones, tol= 1e-2) {
  bdry <- rep(0, max.ones*(max.ones+1)/2)
  zz <- .Fortran("getbdry",
                 as.double(eta),
                 as.integer(max.ones),
                 as.integer(nperm),
                 as.integer(max.ones*(max.ones+1)/2),
                 bdry=as.integer(bdry),
                 etastr=double(max.ones),
                 as.double(tol),
                 PACKAGE="DNAcopy")
#  list("eta.star"=zz$etastr, "boundary"=zz$bdry)
  zz$bdry
}
segment <- function(x, alpha=0.01, nperm=10000, p.method = c("hybrid","perm"),
                    kmax=25, nmin=200, eta=0.05, sbdry=NULL, window.size=NULL,
                    overlap=0.25, trim = 0.025, undo.splits=c("none","prune",
                      "sdundo"), undo.prune=0.05, undo.SD=3, verbose=1)
  {
    if (!inherits(x, 'CNA')) stop("First arg must be a copy number array object")
    call <- match.call()
    if (missing(sbdry)) {
      if (nperm==10000 & alpha==0.01 & eta==0.05) {
        sbdry <- default.DNAcopy.bdry
      } else {
        max.ones <- floor(nperm*alpha) + 1
        sbdry <- getbdry(eta, nperm, max.ones)
      }
    }
    sbn <- length(sbdry)
    nsample <- ncol(x)-2
    sampleid <- colnames(x)[-(1:2)]
    uchrom <- unique(x$chrom)
    data.type <- attr(x, "data.type")
    p.method <- match.arg(p.method)
    if (p.method=="hybrid") window.size <- NULL
    undo.splits <- match.arg(undo.splits)
    segres <- list()
    segres$data <- x
    allsegs <- list()
    allsegs$ID <- NULL
    allsegs$chrom <- NULL
    allsegs$loc.start <- NULL
    allsegs$loc.end <- NULL
    allsegs$num.mark <- NULL
    allsegs$seg.mean <- NULL
    for (isamp in 1:nsample) {
      if (verbose>=1) cat(paste("Analyzing:", sampleid[isamp],"\n"))
      genomdati <- x[,isamp+2]
      ina <- which(!is.na(genomdati) & !(abs(genomdati)==Inf))
      genomdati <- genomdati[ina]
      trimmed.SD <- sqrt(trimmed.variance(genomdati, trim))
      chromi <- x$chrom[ina]
#      maploci <- x$maploc[ina]
      sample.lsegs <- NULL
      sample.segmeans <- NULL
      for (ic in uchrom) {
        if (verbose>=2) cat(paste("  current chromosome:", ic, "\n"))
        segci <- changepoints(genomdati[chromi==ic], data.type, alpha, sbdry,
                              sbn, nperm, p.method, window.size, overlap, kmax,
                              nmin, trimmed.SD, undo.splits, undo.prune,
                              undo.SD, verbose)
        sample.lsegs <- c(sample.lsegs, segci$lseg)
        sample.segmeans <- c(sample.segmeans, segci$segmeans)
      }
      sample.nseg <- length(sample.lsegs)
      sample.segs.start <- ina[cumsum(c(1,sample.lsegs[-sample.nseg]))]
      sample.segs.end <- ina[cumsum(sample.lsegs)]
      allsegs$ID <- c(allsegs$ID, rep(isamp,sample.nseg))
      allsegs$chrom <- c(allsegs$chrom, x$chrom[sample.segs.end])
      allsegs$loc.start <- c(allsegs$loc.start, x$maploc[sample.segs.start])
      allsegs$loc.end <- c(allsegs$loc.end, x$maploc[sample.segs.end])
      allsegs$num.mark <- c(allsegs$num.mark, sample.lsegs)
      allsegs$seg.mean <- c(allsegs$seg.mean, sample.segmeans)
    }
    allsegs$ID <- sampleid[allsegs$ID]
    allsegs$seg.mean <- round(allsegs$seg.mean, 4)
    allsegs <- as.data.frame(allsegs)
    allsegs$ID <- as.character(allsegs$ID)
    segres$output <- allsegs
    segres$call <- call    
    class(segres) <- "DNAcopy"
    segres
  }
.First.lib <- function(lib, pkg)
{
    library.dynam("DNAcopy", pkg, lib)
}
