|
1
|
1 # input a genome file and return a file genome.chrom.sizes to be associated with the custom build (or just have it as an output to be used later in the history.
|
|
|
2 # adapted from https://bioexpressblog.wordpress.com/2014/04/15/calculate-length-of-all-sequences-in-an-multi-fasta-file/
|
|
|
3 from sys import argv
|
|
|
4 # python calculating_chrom.sizes.py genome_input.fa output.chrom.sizes
|
|
6
|
5 fasta_source = str(argv[1])
|
|
|
6 genome = str(argv[2])
|
|
|
7 builtin = str(argv[3])
|
|
|
8 prefix = str(argv[4])
|
|
|
9 output = str(argv[5])
|
|
1
|
10 # genome = 'test-data/test.fasta'
|
|
|
11 # output = "test-data/test_chrom.sizes"
|
|
|
12
|
|
|
13 chromSizesoutput = open(output,"w")
|
|
|
14
|
|
|
15 records = []
|
|
|
16 record = False
|
|
|
17 for line in open(genome, 'r').readlines():
|
|
|
18 if line[0] == '>':
|
|
|
19 if record:
|
|
|
20 records.append(record)
|
|
|
21 record = [line.strip("\n").split(' ')[0][1:], 0]
|
|
|
22
|
|
|
23 else:
|
|
|
24 sequence = line.strip('\n')
|
|
|
25 record[1] += len(sequence)
|
|
3
|
26
|
|
|
27 if record not in records:
|
|
|
28 records.append(record)
|
|
|
29
|
|
1
|
30 for seq_record in records:
|
|
4
|
31 if prefix != 'none':
|
|
|
32 output_line = f"{prefix}{seq_record[0]}\t{seq_record[1]}\n"
|
|
|
33 else:
|
|
|
34 output_line = f"{seq_record[0]}\t{seq_record[1]}\n"
|
|
|
35
|
|
1
|
36 chromSizesoutput.write(output_line)
|
|
|
37
|
|
|
38 chromSizesoutput.close()
|
