1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/pair_cores/rasmol_output.py	Fri Apr 29 16:21:23 2011 +0400
     1.3 @@ -0,0 +1,85 @@
     1.4 +#!/usr/bin/python
     1.5 +
     1.6 +import os
     1.7 +import sys
     1.8 +from gzip import GzipFile
     1.9 +from copy import copy
    1.10 +
    1.11 +from protein_pdb import Alignment, Block, Monomer, Sequence
    1.12 +from allpy import structure
    1.13 +
    1.14 +def load(pdb_code):
    1.15 +    path = '/home/students/y08/bnagaev/pdb/%s/pdb%s.ent.gz' % (pdb_code[1:3], pdb_code)
    1.16 +    try:
    1.17 +        return GzipFile(path)
    1.18 +    except:
    1.19 +        return structure.download_pdb(pdb_code)
    1.20 +
    1.21 +def monomer2str(sequence, m):
    1.22 +    code = structure.pdb_id_parse(sequence.name)['code']
    1.23 +    model = structure.pdb_id_parse(sequence.name)['model']
    1.24 +    chain = structure.pdb_id_parse(sequence.name)['chain']
    1.25 +    resi = m.pdb_residue.get_id()[1]
    1.26 +    return '/%s/%i/%s/%i' % (code, model, chain, resi)
    1.27 +
    1.28 +def columns_for_align(columns, sequence, sequence2):
    1.29 +    cc = [c for c in columns if sequence in c and sequence2 in c and hasattr(c[sequence], 'pdb_residue') and hasattr(c[sequence2], 'pdb_residue')]
    1.30 +    str = ' or '.join(monomer2str(sequence, column[sequence]) for column in cc)
    1.31 +    return '((%s) and name ca)' % str
    1.32 +
    1.33 +aln_file = sys.argv[1]
    1.34 +blocks_file = sys.argv[2]
    1.35 +out_pdb = sys.argv[3]
    1.36 +
    1.37 +alignment = Alignment().append_file(open(aln_file))
    1.38 +blocks = alignment.blocks_from_file(open(blocks_file))
    1.39 +
    1.40 +for sequence in alignment.sequences:
    1.41 +    sequence.auto_pdb(pdb_getter=load)
    1.42 +
    1.43 +columns = set()
    1.44 +for b in blocks:
    1.45 +    columns |= set(b.columns)
    1.46 +columns = list(columns)
    1.47 +
    1.48 +_sequences = copy(alignment.sequences)
    1.49 +base_seq = _sequences.pop()
    1.50 +
    1.51 +alignment.superimpose(columns, extra_columns=True)
    1.52 +map = alignment.save_pdb(open(out_pdb, 'w'))
    1.53 +
    1.54 +for sequence, chain_model in map.items():
    1.55 +    print 'define str_%s *:%s' % (sequence.name, chain_model[0])
    1.56 +
    1.57 +print 'select all'
    1.58 +print 'color chain'
    1.59 +print 'restrict none'
    1.60 +print 'background [230, 240, 250]'
    1.61 +print 'select none'
    1.62 +
    1.63 +n = 0
    1.64 +
    1.65 +column2pos = alignment.column2pos()
    1.66 +
    1.67 +for block in blocks:
    1.68 +    selection = []
    1.69 +    for b in block.continuous_blocks():
    1.70 +        n += 1
    1.71 +        selections = []
    1.72 +        for sequence in b.sequences:
    1.73 +            try:
    1.74 +                chain = map[sequence][0]
    1.75 +                resi_from = b.columns[0][sequence].pdb_residue.get_id()[1]
    1.76 +                resi_to = b.columns[-1][sequence].pdb_residue.get_id()[1]
    1.77 +                selections.append('%i-%i:%s' % (resi_from, resi_to, chain))
    1.78 +            except:
    1.79 +                pass
    1.80 +        print ''
    1.81 +        print 'backbone'
    1.82 +        print 'define block_%i %s' % (n, ','.join(selections))
    1.83 +        print 'select block_%i' % n
    1.84 +        print 'echo block_%i %s %i-%i' % (n, ', '.join([b.sequences[0].name, b.sequences[1].name]),
    1.85 +            column2pos[b.columns[0]], column2pos[b.columns[-1]])
    1.86 +        print 'backbone 100'
    1.87 +        print 'pause'
    1.88 +