R/gb2fasta.R

In seqinr: Biological Sequences Retrieval and Analysis

########################################################################
#
#                                   gb2fasta
#
#                   single entry genbank to fasta conversion
#
########################################################################
gb2fasta <- function(source.file, destination.file)
{
    input <- readLines(source.file)
    head <- input[1]
    head <- unlist(strsplit(head, split=" "))
    head <- head[nchar(head) > 0]
    seqname <- head[2]
    seqsize <- as.integer(head[3])
    outheader <- sprintf(">%s %d bp", seqname, seqsize)
    
    confile <- file(destination.file, open="w")
    writeLines(outheader, confile)
    #
    # Look for sequence position:
    #
    debut <- which(substring(input,1,6) == "ORIGIN") + 1
    if( length( debut ) > 1 )
        stop("Multiple entries not yet implemented !")
    fin <- which(substring(input,1,2) == "//") - 1
    if( length( fin ) > 1 )
        stop("Multiple entries not yet implemented !")
    
    input <- input[debut:fin]
    input <- sapply(input, function(x) {
        return(paste(substr(x,11,20),substr(x,22,31),substr(x,33,42),substr(x,44,53),
                     substr(x,55,64),substr(x,66,75),sep="",collapse="")) } )
    names(input)<-NULL
    writeLines(input, confile )
    close(confile)
}