← Back to branch summary

~ubuntu-branches/ubuntu/quantal/genometools/quantal-backports

« back to all changes in this revision

Viewing changes to scripts/splitmultifasta.rb

  • Committer: Package Import Robot
  • Author(s): Sascha Steinbiss
  • Date: 2012-07-09 14:10:23 UTC
  • Revision ID: package-import@ubuntu.com-20120709141023-juuu4spm6chqsf9o
Tags: upstream-1.4.1
Import upstream version 1.4.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
scripts/splitmultifasta.rb
 
1
#!/usr/bin/env ruby
 
2
 
 
3
# Split a fasta file into different files
 
4
# Stefan Kurtz, October 10, 2009.
 
5
 
 
6
# first parameter is prefix of files to generate. 
 
7
# second parameter is number of files to generate, 0 means to generate
 
8
# one file per sequence
 
9
# third parameter is input file.
 
10
 
 
11
# let infile be the name of the input file in fasta format. Suppose 
 
12
# that it contain 1952 sequences. Then
 
13
 
 
14
# splitmultifasta.rb tmp 0 infile
 
15
 
 
16
# generates 1952 files named tmp-0000, tmp-0001, ..., tmp-1951.
 
17
# each containing one sequence
 
18
# To check that the split was correct, execute the following commands:
 
19
# cat tmp-* > ALL
 
20
# diff ALL infile
 
21
# when nothing is reported by diff, then everything is fine. Otherwise
 
22
# check if there are other (not generated files) that begin with
 
23
# tmp-
 
24
 
 
25
def countnumofsequences(inputfile)
 
26
  seqcount = 0
 
27
  File.open(inputfile).each_line do |line|
 
28
    if line.match(/^>/)
 
29
      seqcount+=1
 
30
    end
 
31
  end
 
32
  return seqcount
 
33
end
 
34
 
 
35
def openoutfile(filename)
 
36
begin
 
37
  outfp = File.new(filename,"w")
 
38
rescue => error
 
39
  STDERR.puts "#{$0}: cannot open \"#{filename}\": #{error}"
 
40
  exit 1
 
41
end
 
42
  return outfp
 
43
end
 
44
 
 
45
def log10func(n)
 
46
  return (Math.log(n.to_f)/Math.log(10.0)).to_i
 
47
end
 
48
 
 
49
def splitfiles(inputfile,splitprefix,numoffiles,numofsequences)
 
50
  seqcount = 0
 
51
  filenum = 0
 
52
  fh = nil
 
53
  maxseqnum = nil
 
54
  numwidth = nil
 
55
 
 
56
  if numoffiles == 0
 
57
    maxseqnum = 1
 
58
    numwidth = 1+log10func(numofsequences-1)
 
59
  else
 
60
    maxseqnum = numofsequences/numoffiles + numofsequences % numoffiles
 
61
    numwidth = 1+log10func(numoffiles-1)
 
62
  end
 
63
  File.open(inputfile).each_line do |line|
 
64
    if line.match(/^\s*$/)     # discard blank line
 
65
      next
 
66
    elsif line.match(/^\s*#/)  # discard comment line
 
67
      next
 
68
    elsif line.match(/^>/)
 
69
      if seqcount >= maxseqnum
 
70
        seqcount = 0
 
71
      end
 
72
      if seqcount == 0
 
73
        outfilename = sprintf("%s-%0*d",splitprefix,numwidth,filenum)
 
74
        fh = openoutfile(outfilename)
 
75
        filenum+=1
 
76
      end
 
77
      seqcount += 1
 
78
    end
 
79
    fh.print line
 
80
  end
 
81
end
 
82
 
 
83
if ARGV.length != 3
 
84
  STDERR.puts "Usage: #{$0} <splitprefix> <numoffiles> <fastafile>"
 
85
  exit 1
 
86
end
 
87
 
 
88
splitprefix = ARGV[0]
 
89
numoffiles = ARGV[1].to_i
 
90
inputfile = ARGV[2]
 
91
 
 
92
numofsequences = countnumofsequences(inputfile)
 
93
splitfiles(inputfile,splitprefix,numoffiles,numofsequences)