[2] | 1 | #!/usr/bin/env python |
---|
| 2 | |
---|
| 3 | import os, sys, tempfile |
---|
| 4 | |
---|
| 5 | assert sys.version_info[:2] >= (2.4) |
---|
| 6 | |
---|
| 7 | def stop_err( msg ): |
---|
| 8 | sys.stderr.write( "%s\n" % msg ) |
---|
| 9 | sys.exit() |
---|
| 10 | |
---|
| 11 | def check_nib_file( dbkey, GALAXY_DATA_INDEX_DIR ): |
---|
| 12 | nib_file = "%s/alignseq.loc" % GALAXY_DATA_INDEX_DIR |
---|
| 13 | nib_path = '' |
---|
| 14 | nibs = {} |
---|
| 15 | for i, line in enumerate( file( nib_file ) ): |
---|
| 16 | line = line.rstrip( '\r\n' ) |
---|
| 17 | if line and not line.startswith( "#" ): |
---|
| 18 | fields = line.split( '\t' ) |
---|
| 19 | if len( fields ) < 3: |
---|
| 20 | continue |
---|
| 21 | if fields[0] == 'seq': |
---|
| 22 | nibs[( fields[1] )] = fields[2] |
---|
| 23 | if nibs.has_key( dbkey ): |
---|
| 24 | nib_path = nibs[( dbkey )] |
---|
| 25 | return nib_path |
---|
| 26 | |
---|
| 27 | def check_twobit_file( dbkey, GALAXY_DATA_INDEX_DIR ): |
---|
| 28 | twobit_file = "%s/twobit.loc" % GALAXY_DATA_INDEX_DIR |
---|
| 29 | twobit_path = '' |
---|
| 30 | twobits = {} |
---|
| 31 | for i, line in enumerate( file( twobit_file ) ): |
---|
| 32 | line = line.rstrip( '\r\n' ) |
---|
| 33 | if line and not line.startswith( "#" ): |
---|
| 34 | fields = line.split( '\t' ) |
---|
| 35 | if len( fields ) < 2: |
---|
| 36 | continue |
---|
| 37 | twobits[( fields[0] )] = fields[1] |
---|
| 38 | if twobits.has_key( dbkey ): |
---|
| 39 | twobit_path = twobits[( dbkey )] |
---|
| 40 | return twobit_path |
---|
| 41 | |
---|
| 42 | def __main__(): |
---|
| 43 | # I/O |
---|
| 44 | source_format = sys.argv[1] # 0: dbkey; 1: upload file |
---|
| 45 | target_file = sys.argv[2] |
---|
| 46 | query_file = sys.argv[3] |
---|
| 47 | output_file = sys.argv[4] |
---|
| 48 | min_iden = sys.argv[5] |
---|
| 49 | tile_size = sys.argv[6] |
---|
| 50 | one_off = sys.argv[7] |
---|
| 51 | |
---|
| 52 | try: |
---|
| 53 | float(min_iden) |
---|
| 54 | except: |
---|
| 55 | stop_err('Invalid value for minimal identity.') |
---|
| 56 | |
---|
| 57 | try: |
---|
| 58 | test = int(tile_size) |
---|
| 59 | assert test >= 6 and test <= 18 |
---|
| 60 | except: |
---|
| 61 | stop_err('Invalid value for tile size. DNA word size must be between 6 and 18.') |
---|
| 62 | |
---|
| 63 | try: |
---|
| 64 | test = int(one_off) |
---|
| 65 | assert test >= 0 and test <= int(tile_size) |
---|
| 66 | except: |
---|
| 67 | stop_err('Invalid value for mismatch numbers in the word') |
---|
| 68 | |
---|
| 69 | GALAXY_DATA_INDEX_DIR = sys.argv[8] |
---|
| 70 | |
---|
| 71 | all_files = [] |
---|
| 72 | if source_format == '0': |
---|
| 73 | # check target genome |
---|
| 74 | dbkey = target_file |
---|
| 75 | nib_path = check_nib_file( dbkey, GALAXY_DATA_INDEX_DIR ) |
---|
| 76 | twobit_path = check_twobit_file( dbkey, GALAXY_DATA_INDEX_DIR ) |
---|
| 77 | if not os.path.exists( nib_path ) and not os.path.exists( twobit_path ): |
---|
| 78 | stop_err("No sequences are available for %s, request them by reporting this error." % dbkey) |
---|
| 79 | |
---|
| 80 | # check the query file, see whether all of them are legitimate sequence |
---|
| 81 | if nib_path and os.path.isdir( nib_path ): |
---|
| 82 | compress_files = os.listdir(nib_path) |
---|
| 83 | target_path = nib_path |
---|
| 84 | elif twobit_path: |
---|
| 85 | compress_files = [twobit_path] |
---|
| 86 | target_path = "" |
---|
| 87 | else: |
---|
| 88 | stop_err("Requested genome build has no available sequence.") |
---|
| 89 | |
---|
| 90 | for file in compress_files: |
---|
| 91 | file = "%s/%s" % ( target_path, file ) |
---|
| 92 | file = os.path.normpath(file) |
---|
| 93 | all_files.append(file) |
---|
| 94 | else: |
---|
| 95 | all_files = [target_file] |
---|
| 96 | |
---|
| 97 | for detail_file_path in all_files: |
---|
| 98 | output_tempfile = tempfile.NamedTemporaryFile().name |
---|
| 99 | command = "blat %s %s %s -oneOff=%s -tileSize=%s -minIdentity=%s -mask=lower -noHead -out=pslx 2>&1" % ( detail_file_path, query_file, output_tempfile, one_off, tile_size, min_iden ) |
---|
| 100 | os.system( command ) |
---|
| 101 | os.system( 'cat %s >> %s' % ( output_tempfile, output_file ) ) |
---|
| 102 | os.remove( output_tempfile ) |
---|
| 103 | |
---|
| 104 | if __name__ == '__main__': __main__() |
---|