source: trunk/exports/G2export_PDB.py @ 2152

Last change on this file since 2152 was 1299, checked in by toby, 11 years ago

Set conf flags consistently for .py files
  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 12.6 KB
Line 
1#!/usr/bin/env python
2# -*- coding: utf-8 -*-
3########### SVN repository information ###################
4# $Date: 2014-04-23 22:23:37 +0000 (Wed, 23 Apr 2014) $
5# $Author: vondreele $
6# $Revision: 1299 $
7# $URL: trunk/exports/G2export_PDB.py $
8# $Id: G2export_PDB.py 1299 2014-04-23 22:23:37Z vondreele $
9########### SVN repository information ###################
10'''
11*Module G2export_PDB: Macromolecular export*
12--------------------------------------------
13Code to export a phase into the venerated/obsolete (pick one)
14ASCII PDB format. Also defines exporter :class:`ExportPhaseCartXYZ`
15which writes atom positions in orthogonal coordinates for a phase.
16
17'''
18import numpy as np
19import os.path
20import GSASIIpath
21GSASIIpath.SetVersionNumber("$Revision: 1299 $")
22import GSASIIIO as G2IO
23#import GSASIIgrid as G2gd
24import GSASIIstrIO as G2stIO
25#import GSASIImath as G2mth
26import GSASIIlattice as G2lat
27import GSASIIspc as G2spc
28#import GSASIIphsGUI as G2pg
29#import GSASIIstrMain as G2stMn
30
31class ExportPhasePDB(G2IO.ExportBaseclass):
32    '''Used to create a PDB file for a phase
33
34    :param wx.Frame G2frame: reference to main GSAS-II frame
35    '''
36    def __init__(self,G2frame):
37        super(self.__class__,self).__init__( # fancy way to say <parentclass>.__init__
38            G2frame=G2frame,
39            formatName = 'PDB',
40            extension='.PDB',
41            longFormatName = 'Export phase as .PDB file'
42            )
43        self.exporttype = ['phase']
44        self.multiple = True
45
46    def Exporter(self,event=None):
47        '''Export as a PDB file
48        '''
49       
50        def PDBheader():
51            self.Write("HEADER phase "+str(phasenam)+" from "+str(self.G2frame.GSASprojectfile))
52            self.Write("TITLE")
53            self.Write("COMPND")
54            self.Write("SOURCE")
55            self.Write("KEYWDS")
56            self.Write("EXPDTA    X-RAY POWDER DIFFRACTION")
57            self.Write("REVDAT")
58            self.Write("JRNL")
59            self.Write("REMARK   1")
60            self.Write("REMARK   2")                                                                     
61            self.Write("REMARK   2 RESOLUTION. 2.66 ANGSTROMS.")                                         
62            self.Write("REMARK   2 POWDER DIFFRACTION MINIMUM D-SPACING.")
63           
64        def PDBremark250():                               
65            self.Write('REMARK 250')                                                                     
66            self.Write('REMARK 250 REFINEMENT.')                                                         
67            self.Write('REMARK 250   PROGRAM     : GSAS-II')                                                 
68            self.Write('REMARK 250   AUTHORS     : TOBY & VON DREELE')
69            self.Write('REMARK 250   REFRENCE    : J. APPL. CRYST. 46, 544-549(2013)')                             
70            self.Write('REMARK 250')
71            self.Write('REMARK 250  DATA USED IN REFINEMENT')                                             
72            self.Write('REMARK 250   RESOLUTION RANGE HIGH (ANGSTROMS) :  x.xx')                         
73            self.Write('REMARK 250   RESOLUTION RANGE LOW  (ANGSTROMS) : xx.xx')                         
74            self.Write('REMARK 250   POWDER DIFFRACTION DATA.')                                           
75            self.Write('REMARK 250')
76            self.Write('REMARK 250  FIT TO DATA USED IN REFINEMENT')                                     
77            self.Write('REMARK 250   NUMBER OF POWDER PATTERNS         :     x')                         
78            self.Write('REMARK 250   PROFILE R VALUES              (%) :  x.xx')                         
79            self.Write('REMARK 250   WEIGHTED PROFILE R VALUES     (%) :  x.xx')                         
80            self.Write('REMARK 250   F**2 R VALUES                 (%) : xx.xx')                         
81            self.Write('REMARK 250   NUMBERS OF POWDER PATTERN POINTS  :  xxxx')                         
82            self.Write('REMARK 250   NUMBERS OF REFLECTIONS            :  xxxx')                         
83            self.Write('REMARK 250   TOTAL NUMBER OF POWDER POINTS     :  xxxx')                         
84            self.Write('REMARK 250')
85            self.Write('REMARK 250  NUMBER OF NON-HYDROGEN ATOMS USED IN REFINEMENT.')                   
86            self.Write('REMARK 250   PROTEIN ATOMS       :      xxxx')                                             
87            self.Write('REMARK 250   NUCLEIC ACID ATOMS  :      xxxx')                                             
88            self.Write('REMARK 250   HETEROGEN ATOMS     :      xxxx')                                             
89            self.Write('REMARK 250   SOLVENT ATOMS       :      xxxx')                                             
90            self.Write('REMARK 250')
91            self.Write('REMARK 250  MODEL REFINEMENT.')                                                   
92            self.Write('REMARK 250   NUMBER OF LEAST-SQUARES PARAMETERS :  xxxx')                         
93            self.Write('REMARK 250   NUMBER OF RESTRAINTS               :  xxxx')                         
94            self.Write('REMARK 250')
95            self.Write('REMARK 250  RMS DEVIATIONS FROM RESTRAINT TARGET VALUES. NUMBER.')               
96            self.Write('REMARK 250   BOND ANGLES                      (DEG) : x.xx   xxx')               
97#            self.Write('REMARK 250   ANTI-BUMPING DISTANCE RESTRAINTS   (A) :x.xxx   xxx')               
98#            self.Write('REMARK 250   HYDROGEN BOND DISTANCE RESTRAINTS  (A) :x.xxx   xxx')               
99            self.Write('REMARK 250   INTERATOMIC DISTANCES              (A) :x.xxx   xxx')               
100            self.Write('REMARK 250   DISTANCES FROM RESTRAINT PLANES    (A) :x.xxx   xxx')               
101            self.Write('REMARK 250   TORSION PSEUDOPOTENTIAL RESTRAINTS (E) : x.xx   xxx')               
102            self.Write('REMARK 250   TORSION ANGLE RESTRAINTS           (E) : x.xx   xxx')               
103            self.Write('REMARK 250')
104            self.Write('REMARK 200')                                                                     
105            self.Write('DBREF')
106           
107        def PDBseqres(seqList):
108            chains = seqList.keys()
109            chains.sort()
110            nSeq = 0
111            for chain in chains:
112                nres = len(seqList[chain])
113                nrec = (nres-1)/13+1
114                iB = 0
115                for irec in range(nrec):
116                    iF = min(iB+13,nres)
117                    text = 'SEQRES {:3d}{:2s}{:5d}  '+(iF-iB)*'{:4s}'
118                    self.Write(text.format(irec+1,chain,nres,*seqList[chain][iB:iF]))
119                    nSeq += 1
120                    iB += 13
121            return nSeq
122           
123        # the export process starts here
124        self.InitExport(event)
125        # load all of the tree into a set of dicts
126        self.loadTree()
127        # create a dict with refined values and their uncertainties
128        self.loadParmDict()
129        if self.ExportSelect():    # set export parameters; ask for file name
130            return
131        filename = self.filename
132        for phasenam in self.phasenam:
133            phasedict = self.Phases[phasenam] # pointer to current phase info
134            General = phasedict['General']
135            if General['Type'] != 'macromolecular':
136                print 'phase '+str(phasenam)+' not macromolecular, skipping'
137                continue
138            i = self.Phases[phasenam]['pId']
139            if len(self.phasenam) > 1: # if more than one filename is included, add a phase #
140                self.filename = os.path.splitext(filename)[1] + "_" + str(i) + self.extension
141            fp = self.OpenFile()
142            Atoms = phasedict['Atoms']
143            cx,ct,cs,cia = General['AtomPtrs']
144            seqList = {}
145            AA3letter = ['ALA','ARG','ASN','ASP','CYS','GLN','GLU','GLY','HIS','ILE',
146                'LEU','LYS','MET','PHE','PRO','SER','THR','TRP','TYR','VAL','MSE']
147            seq = 0
148            for atom in Atoms:
149                if atom[ct-3] in AA3letter and int(atom[ct-4]) != seq:
150                    if atom[ct-2] not in seqList:
151                        seqList[atom[ct-2]] = []
152                    seqList[atom[ct-2]].append(atom[ct-3])
153                    seq = int(atom[ct-4])
154            PDBheader()
155            PDBremark250()
156            nSeq = PDBseqres(seqList)
157           
158            # get cell parameters
159            Cell = General['Cell'][1:7]
160            line = "CRYST1 {:8.3f} {:8.3f} {:8.3f} {:6.2f} {:6.2f} {:6.2f} ".format(*Cell)
161            line += General['SGData']['SpGrp'].ljust(13)
162            line += '%2d'%(len(General['SGData']['SGOps'])*len(General['SGData']['SGCen']))
163            self.Write(line)
164            self.Write('ORIGX1      1.000000  0.000000  0.000000        0.00000')
165            self.Write('ORIGX2      0.000000  1.000000  0.000000        0.00000')
166            self.Write('ORIGX3      0.000000  0.000000  1.000000        0.00000')
167            A,B = G2lat.cell2AB(Cell)
168            self.Write('SCALE1     {:9.6f} {:9.6f} {:9.6f}        0.00000'.format(*B[0]))
169            self.Write('SCALE2     {:9.6f} {:9.6f} {:9.6f}        0.00000'.format(*B[1]))
170            self.Write('SCALE3     {:9.6f} {:9.6f} {:9.6f}        0.00000'.format(*B[2]))
171            iatom = 1
172            nHet = 0
173            nTer = 0
174            fmt = '{:6s}{:5d}  {:4s}{:3s} {:1s}{:4s}    '+3*'{:8.3f}'+2*'{:6.2f}'+'{:s}'
175            for atom in Atoms:
176                if atom[cia] == 'I':    #need to deal with aniso thermals for proteins = "ANISOU" records
177                    Biso = atom[cia+1]*8.*np.pi**2
178                xyz = np.inner(A,np.array(atom[cx:cx+3]))
179                if atom[ct-3] in AA3letter:
180                    self.Write(fmt.format('ATOM  ',iatom,atom[ct-1],atom[ct-3].strip(),    \
181                        atom[ct-2].strip(),atom[ct-4].rjust(4),xyz[0],xyz[1],xyz[2],atom[cx+3], \
182                        Biso,atom[ct].rjust(12)))
183                    if atom[ct-1] == 'OXT':
184                        iatom += 1
185                        self.Write('{:6s}{:5d}  {:4s}{:3s}'.format('TER   ',iatom,atom[ct-1],atom[ct-3].strip()))
186                        nTer += 1
187                else:
188                    nHet += 1
189                    self.Write(fmt.format('HETATM',iatom,atom[ct-1],atom[ct-3].strip(),    \
190                        atom[ct-2].strip(),atom[ct-4].rjust(4),xyz[0],xyz[1],xyz[2],atom[cx+3], \
191                        Biso,atom[ct].rjust(12)))
192                iatom += 1
193           
194            vals = [3,0,nHet,0,0,0,6,len(Atoms),nTer,0,nSeq]
195            fmt = 'MASTER'+11*'{:5d}'
196            self.Write(fmt.format(*vals))
197            self.Write('END')
198            self.CloseFile()
199            print('Phase '+str(phasenam)+' written to PDB file '+str(self.fullpath))
200
201class ExportPhaseCartXYZ(G2IO.ExportBaseclass):
202    '''Used to create a Cartesian XYZ file for a phase
203
204    :param wx.Frame G2frame: reference to main GSAS-II frame
205    '''
206    def __init__(self,G2frame):
207        super(self.__class__,self).__init__( # fancy way to say <parentclass>.__init__
208            G2frame=G2frame,
209            formatName = 'Cartesian XYZ',
210            extension='.XYZ',
211            longFormatName = 'Export phase with Cartesian coordinates as .XYZ file'
212            )
213        self.exporttype = ['phase']
214        self.multiple = True
215
216    def Exporter(self,event=None):
217        '''Export as a XYZ file
218        '''
219        # the export process starts here
220        self.InitExport(event)
221        # load all of the tree into a set of dicts
222        self.loadTree()
223        # create a dict with refined values and their uncertainties
224        self.loadParmDict()
225        if self.ExportSelect():    # set export parameters; ask for file name
226            return
227        filename = self.filename
228        for phasenam in self.phasenam:
229            phasedict = self.Phases[phasenam] # pointer to current phase info
230            General = phasedict['General']
231            i = self.Phases[phasenam]['pId']
232            Atoms = phasedict['Atoms']
233            if not len(Atoms):
234                print('**** ERROR - Phase '+str(phasenam)+' has no atoms! ****')
235                continue
236            if len(self.phasenam) > 1: # if more than one filename is included, add a phase #
237                self.filename = os.path.splitext(filename)[1] + "_" + str(i) + self.extension
238            fp = self.OpenFile()
239            cx,ct,cs,cia = General['AtomPtrs']
240            Cell = General['Cell'][1:7]
241            A,B = G2lat.cell2AB(Cell)
242            fmt = '{:4s}'+3*'{:12.4f}'
243            self.Write('{:6d}'.format(len(Atoms)))
244            self.Write(' ')
245            for atom in Atoms:
246                xyz = np.inner(A,np.array(atom[cx:cx+3]))
247                self.Write(fmt.format(atom[ct],*xyz))
248            self.CloseFile()
249            print('Phase '+str(phasenam)+' written to XYZ file '+str(self.fullpath))
250   
Note: See TracBrowser for help on using the repository browser.