source: trunk/exports/G2cif.py @ 1012

Last change on this file since 1012 was 1012, checked in by vondreele, 8 years ago

add svn header to G2cif.py and GSASIIpath to imports
add DEBUG to skip prnting of profile & reflections
add profile parm cif output
work on MC/SA scrolling

  • Property svn:eol-style set to native
File size: 56.3 KB
Line 
1#!/usr/bin/env python
2# -*- coding: utf-8 -*-
3#G2cif
4########### SVN repository information ###################
5# $Date: 2013-07-22 20:57:37 -0500 (Mon, 22 Jul 2013) $
6# $Author: toby $
7# $Revision: 1006 $
8# $URL: https://subversion.xray.aps.anl.gov/pyGSAS/trunk/exports/G2cif.py $
9# $Id: G2cif.py 1006 2013-07-23 01:57:37Z toby $
10########### SVN repository information ###################
11'''Development code to export a GSAS-II project as a CIF
12The heavy lifting is done in method export
13'''
14
15# TODO: need a mechanism for editing of instrument names, bond pub flags, templates,...
16
17import datetime as dt
18import os.path
19import numpy as np
20import GSASIIpath
21GSASIIpath.SetVersionNumber("$Revision: 1006 $")
22import GSASIIIO as G2IO
23#reload(G2IO)
24import GSASIIgrid as G2gd
25import GSASIIstrIO as G2stIO
26#reload(G2stIO)
27#import GSASIImapvars as G2mv
28import GSASIImath as G2mth
29reload(G2mth)
30import GSASIIlattice as G2lat
31import GSASIIspc as G2spg
32#reload(G2spg)
33
34DEBUG = True    #True to skip printing of reflection/powder profile lists
35
36def getCallerDocString(): # for development
37    "Return the calling function's doc string"
38    import inspect as ins
39    for item in ins.stack()[1][0].f_code.co_consts:
40        if type(item) is str:
41            return item
42    else:
43        return '?'
44
45class ExportCIF(G2IO.ExportBaseclass):
46    def __init__(self,G2frame):
47        super(self.__class__,self).__init__( # fancy way to say <parentclass>.__init__
48            G2frame=G2frame,
49            formatName = 'full CIF',
50            longFormatName = 'Export project as complete CIF'
51            )
52        self.author = ''
53
54    def export(self,mode='full'):
55        '''Export a CIF
56
57        :param str mode: "full" (default) to create a complete CIF of project,
58          "simple" for a simple CIF with only coordinates
59        '''
60   
61        def WriteCIFitem(name,value=''):
62            if value:
63                if "\n" in value or len(value)> 70:
64                    if name.strip(): print name
65                    print '; '+value
66                    print '; '
67                elif " " in value:
68                    if len(name)+len(value) > 65:
69                        print name,'\n   ','"' + str(value) + '"'
70                    else:
71                        print name,'  ','"' + str(value) + '"'
72                else:
73                    if len(name)+len(value) > 65:
74                        print name,'\n   ',value
75                    else:
76                        print name,'  ',value
77            else:
78                print name
79
80        def WriteAudit():
81            WriteCIFitem('_audit_creation_method',
82                         'created in GSAS-II')
83            WriteCIFitem('_audit_creation_date',self.CIFdate)
84            if self.author:
85                WriteCIFitem('_audit_author_name',self.author)
86            WriteCIFitem('_audit_update_record',
87                         self.CIFdate+'  Initial software-generated CIF')
88
89        def WriteOverall():
90            '''Write out overall refinement information.
91
92            More could be done here, but this is a good start.
93            '''
94            WriteCIFitem('_pd_proc_info_datetime', self.CIFdate)
95            WriteCIFitem('_pd_calc_method', 'Rietveld Refinement')
96            #WriteCIFitem('_refine_ls_shift/su_max',DAT1)
97            #WriteCIFitem('_refine_ls_shift/su_mean',DAT2)
98            WriteCIFitem('_computing_structure_refinement','GSAS-II')
99            try:
100                vars = str(len(self.OverallParms['Covariance']['varyList']))
101            except:
102                vars = '?'
103            WriteCIFitem('_refine_ls_number_parameters',vars)
104            try:
105                GOF = G2mth.ValEsd(self.OverallParms['Covariance']['Rvals']['GOF'],-0.009)
106            except:
107                GOF = '?'
108            WriteCIFitem('_refine_ls_goodness_of_fit_all',GOF)
109
110            # get restraint info
111            # restraintDict = self.OverallParms.get('Restraints',{})
112            # for i in  self.OverallParms['Constraints']:
113            #     print i
114            #     for j in self.OverallParms['Constraints'][i]:
115            #         print j
116            #WriteCIFitem('_refine_ls_number_restraints',TEXT)
117           
118            # other things to consider reporting
119            # _refine_ls_number_reflns
120            # _refine_ls_goodness_of_fit_obs
121            # _refine_ls_R_factor_all
122            # _refine_ls_R_factor_obs
123            # _refine_ls_wR_factor_all
124            # _refine_ls_wR_factor_obs
125            # _refine_ls_restrained_S_all
126            # _refine_ls_restrained_S_obs
127
128            # include an overall profile r-factor, if there is more than one powder histogram
129            if len(self.powderDict) > 1:
130                WriteCIFitem('\n# OVERALL POWDER R-FACTOR')
131                try:
132                    R = str(self.OverallParms['Covariance']['Rvals']['Rwp'])
133                except:
134                    R = '?'
135                WriteCIFitem('_pd_proc_ls_prof_wR_factor',R)
136                #WriteCIFitem('_pd_proc_ls_prof_R_factor',TEXT(11:20)) # who cares!
137            WriteCIFitem('_refine_ls_matrix_type','full')
138            #WriteCIFitem('_refine_ls_matrix_type','userblocks')
139
140        def WritePubTemplate():
141            '''TODO: insert the publication template ``template_publ.cif`` or some modified
142            version for this project. Store this in the GPX file?
143            '''
144            print getCallerDocString()
145
146        def WritePhaseTemplate():
147            '''TODO: insert the phase template ``template_phase.cif`` or some modified
148            version for this project
149            '''
150            print getCallerDocString()
151
152        def WritePowderTemplate():
153            '''TODO: insert the phase template ``template_instrument.cif`` or some modified
154            version for this project
155            '''
156            print getCallerDocString()
157
158        def WriteSnglXtalTemplate():
159            '''TODO: insert the single-crystal histogram template
160            for this project
161            '''
162            print getCallerDocString()
163
164        def FormatSH(phasenam):
165            'Format a full spherical harmonics texture description as a string'
166            phasedict = self.Phases[phasenam] # pointer to current phase info           
167            pfx = str(phasedict['pId'])+'::'
168            s = ""
169            textureData = phasedict['General']['SH Texture']   
170            if textureData.get('Order'):
171                s += "Spherical Harmonics correction. Order = "+str(textureData['Order'])
172                s += " Model: " + str(textureData['Model']) + "\n    Orientation angles: "
173                for name in ['omega','chi','phi']:
174                    aname = pfx+'SH '+name
175                    s += name + " = "
176                    sig = self.sigDict.get(aname,-0.09)
177                    s += G2mth.ValEsd(self.parmDict[aname],sig)
178                    s += "; "
179                s += "\n"
180                s1 = "    Coefficients:  "
181                for name in textureData['SH Coeff'][1]:
182                    aname = pfx+name
183                    if len(s1) > 60:
184                        s += s1 + "\n"
185                        s1 = "    "
186                    s1 += aname + ' = '
187                    sig = self.sigDict.get(aname,-0.0009)
188                    s1 += G2mth.ValEsd(self.parmDict[aname],sig)
189                    s1 += "; "
190                s += s1
191            return s
192
193        def FormatHAPpo(phasenam):
194            '''return the March-Dollase/SH correction for every
195            histogram in the current phase formatted into a
196            character string
197            '''
198            phasedict = self.Phases[phasenam] # pointer to current phase info           
199            s = ''
200            for histogram in sorted(phasedict['Histograms']):
201                if histogram.startswith("HKLF"): continue # powder only
202                Histogram = self.Histograms.get(histogram)
203                if not Histogram: continue
204                hapData = phasedict['Histograms'][histogram]
205                if hapData['Pref.Ori.'][0] == 'MD':
206                    aname = str(phasedict['pId'])+':'+str(Histogram['hId'])+':MD'
207                    if self.parmDict.get(aname,1.0) != 1.0: continue
208                    sig = self.sigDict.get(aname,-0.009)
209                    if s != "": s += '\n'
210                    s += 'March-Dollase correction'
211                    if len(self.powderDict) > 1:
212                        s += ', histogram '+str(Histogram['hId']+1)
213                    s += ' coef. = ' + G2mth.ValEsd(self.parmDict[aname],sig)
214                    s += ' axis = ' + str(hapData['Pref.Ori.'][3])
215                else: # must be SH
216                    if s != "": s += '\n'
217                    s += 'Simple spherical harmonic correction'
218                    if len(self.powderDict) > 1:
219                        s += ', histogram '+str(Histogram['hId']+1)
220                    s += ' Order = '+str(hapData['Pref.Ori.'][4])+'\n'
221                    s1 = "    Coefficients:  "
222                    for item in hapData['Pref.Ori.'][5]:
223                        print item
224                        aname = str(phasedict['pId'])+':'+str(Histogram['hId'])+':'+item
225                        print aname
226                        if len(s1) > 60:
227                            s += s1 + "\n"
228                            s1 = "    "
229                        s1 += aname + ' = '
230                        sig = self.sigDict.get(aname,-0.0009)
231                        s1 += G2mth.ValEsd(self.parmDict[aname],sig)
232                        s1 += "; "
233                    s += s1
234            return s
235        def FormatBackground(bkg):
236            '''Display the Background information as a descriptive text string.
237           
238            TODO: this needs to be expanded to show the diffuse peak and
239            Debye term information as well.
240
241            :returns: the text description (str)
242            '''
243            fxn, bkgdict = bkg
244            terms = fxn[2]
245            txt = 'Background function: "'+fxn[0]+'" function with '+str(terms)+' terms:\n'
246            l = "   "
247            for v in fxn[3:]:
248                if len(l) > 60:
249                    txt += l + '\n'
250                    l = '   '
251                l += G2mth.ValEsd(v,-.009)+', '
252            txt += l
253            return txt
254
255        def FormatInstProfile(instparmdict):
256            '''Format the instrumental profile parameters with a
257            string description. Will only be called on PWDR histograms
258            '''
259            s = ''
260            inst = instparmdict[0]
261            if 'C' in inst['Type'][0]:
262                s = 'Finger-Cox-Jephcoat function parameters U, V, W, X, Y, SH/L:\n   '
263                for item in ['U','V','W','X','Y','SH/L']:
264                    s += G2mth.ValEsd(inst[item][1],-.009)+', '                   
265            elif 'T' in inst['Type'][0]:    #to be tested after TOF Rietveld done
266                s = 'Von Dreele-Jorgenson-Windsor function parameters\n'+ \
267                    '   alpha, beta-0, beta-1, beta-q, sig-0, sig-1, sig-q, X, Y:\n   '
268                for item in ['alpha','bet-0','bet-1','bet-q','sig-0','sig-1','sig-q','X','Y']:
269                    s += G2mth.ValEsd(inst[item][1],-.009)+', '
270            return s
271
272        def FormatPhaseProfile(phasenam):
273            '''Format the phase-related profile parameters (size/strain)
274            with a string description.
275            return an empty string or None there are no
276            powder histograms for this phase.
277            '''
278            s = ''
279            phasedict = self.Phases[phasenam] # pointer to current phase info           
280            for histogram in sorted(phasedict['Histograms']):
281                if histogram.startswith("HKLF"): continue # powder only
282                Histogram = self.Histograms.get(histogram)
283                if not Histogram: continue
284                hapData = phasedict['Histograms'][histogram]
285            return 'TODO: Phase profile goes here'
286       
287        def FmtAtomType(sym):
288            'Reformat a GSAS-II atom type symbol to match CIF rules'
289            sym = sym.replace('_','') # underscores are not allowed: no isotope designation?
290            # in CIF, oxidation state sign symbols come after, not before
291            if '+' in sym:
292                sym = sym.replace('+','') + '+'
293            elif '-' in sym:
294                sym = sym.replace('-','') + '-'
295            return sym
296           
297        def PutInCol(val,wid):
298            '''Pad a value to >=wid+1 columns by adding spaces at the end. Always
299            adds at least one space
300            '''
301            val = str(val).replace(' ','')
302            if not val: val = '?'
303            fmt = '{:' + str(wid) + '} '
304            return fmt.format(val)
305
306        def MakeUniqueLabel(lbl,labellist):
307            'Make sure that every atom label is unique'
308            lbl = lbl.strip()
309            if not lbl: # deal with a blank label
310                lbl = 'A_1'
311            if lbl not in labellist:
312                labellist.append(lbl)
313                return lbl
314            i = 1
315            prefix = lbl
316            if '_' in lbl:
317                prefix = lbl[:lbl.rfind('_')]
318                suffix = lbl[lbl.rfind('_')+1:]
319                try:
320                    i = int(suffix)+1
321                except:
322                    pass
323            while prefix+'_'+str(i) in labellist:
324                i += 1
325            else:
326                lbl = prefix+'_'+str(i)
327                labellist.append(lbl)
328
329        def WriteAtomsNuclear(phasenam):
330            'Write atom positions to CIF'
331            phasedict = self.Phases[phasenam] # pointer to current phase info
332            General = phasedict['General']
333            cx,ct,cs,cia = General['AtomPtrs']
334            Atoms = phasedict['Atoms']
335            cfrac = cx+3
336            fpfx = str(phasedict['pId'])+'::Afrac:'       
337            for i,at in enumerate(Atoms):
338                fval = self.parmDict.get(fpfx+str(i),at[cfrac])
339                if fval != 0.0:
340                    break
341            else:
342                WriteCIFitem('\n# PHASE HAS NO ATOMS!')
343                return
344               
345            WriteCIFitem('\n# ATOMIC COORDINATES AND DISPLACEMENT PARAMETERS')
346            WriteCIFitem('loop_ '+
347                         '\n\t_atom_site_label'+
348                         '\n\t_atom_site_type_symbol'+
349                         '\n\t_atom_site_fract_x'+
350                         '\n\t_atom_site_fract_y'+
351                         '\n\t_atom_site_fract_z'+
352                         '\n\t_atom_site_occupancy'+
353                         '\n\t_atom_site_adp_type'+
354                         '\n\t_atom_site_U_iso_or_equiv'+
355                         '\n\t_atom_site_symmetry_multiplicity')
356
357            varnames = {cx:'Ax',cx+1:'Ay',cx+2:'Az',cx+3:'Afrac',
358                        cia+1:'AUiso',cia+2:'AU11',cia+3:'AU22',cia+4:'AU33',
359                        cia+5:'AU12',cia+6:'AU13',cia+7:'AU23'}
360            self.labellist = []
361           
362            pfx = str(phasedict['pId'])+'::'
363            # loop over all atoms
364            naniso = 0
365            for i,at in enumerate(Atoms):
366                s = PutInCol(MakeUniqueLabel(at[ct-1],self.labellist),6) # label
367                fval = self.parmDict.get(fpfx+str(i),at[cfrac])
368                if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
369                s += PutInCol(FmtAtomType(at[ct]),4) # type
370                if at[cia] == 'I':
371                    adp = 'Uiso '
372                else:
373                    adp = 'Uani '
374                    naniso += 1
375                    # compute Uequiv crudely
376                    # correct: Defined as "1/3 trace of diagonalized U matrix".
377                    # SEE cell2GS & Uij2Ueqv to GSASIIlattice. Former is needed to make the GS matrix used by the latter.
378                    t = 0.0
379                    for j in (2,3,4):
380                        var = pfx+varnames[cia+j]+":"+str(i)
381                        t += self.parmDict.get(var,at[cia+j])
382                for j in (cx,cx+1,cx+2,cx+3,cia+1):
383                    if j in (cx,cx+1,cx+2):
384                        dig = 11
385                        sigdig = -0.00009
386                    else:
387                        dig = 10
388                        sigdig = -0.009
389                    var = pfx+varnames[j]+":"+str(i)
390                    dvar = pfx+"d"+varnames[j]+":"+str(i)
391                    if dvar not in self.sigDict:
392                        dvar = var
393                    if j == cia+1 and adp == 'Uani ':
394                        val = t/3.
395                        sig = sigdig
396                    else:
397                        #print var,(var in self.parmDict),(var in self.sigDict)
398                        val = self.parmDict.get(var,at[j])
399                        sig = self.sigDict.get(dvar,sigdig)
400                    s += PutInCol(G2mth.ValEsd(val,sig),dig)
401                s += adp
402                s += PutInCol(at[cs+1],3)
403                WriteCIFitem(s)
404            if naniso == 0: return
405            # now loop over aniso atoms
406            WriteCIFitem('\nloop_' + '\n\t_atom_site_aniso_label' + 
407                         '\n\t_atom_site_aniso_U_11' + '\n\t_atom_site_aniso_U_12' +
408                         '\n\t_atom_site_aniso_U_13' + '\n\t_atom_site_aniso_U_22' +
409                         '\n\t_atom_site_aniso_U_23' + '\n\t_atom_site_aniso_U_33')
410            for i,at in enumerate(Atoms):
411                fval = self.parmDict.get(fpfx+str(i),at[cfrac])
412                if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
413                if at[cia] == 'I': continue
414                s = PutInCol(self.labellist[i],6) # label
415                for j in (2,3,4,5,6,7):
416                    sigdig = -0.0009
417                    var = pfx+varnames[cia+j]+":"+str(i)
418                    val = self.parmDict.get(var,at[cia+j])
419                    sig = self.sigDict.get(var,sigdig)
420                    s += PutInCol(G2mth.ValEsd(val,sig),11)
421                WriteCIFitem(s)
422
423        def HillSortElements(elmlist):
424            '''Sort elements in "Hill" order: C, H, others, (where others
425            are alphabetical).
426
427            :params list elmlist: a list of element strings
428
429            :returns: a sorted list of element strings
430            '''
431            newlist = []
432            oldlist = elmlist[:]
433            for elm in ('C','H'):
434                if elm in elmlist:
435                    newlist.append(elm)
436                    oldlist.pop(oldlist.index(elm))
437            return newlist+sorted(oldlist)
438
439        def WriteComposition(phasenam):
440            '''determine the composition for the unit cell, crudely determine Z and
441            then compute the composition in formula units
442            '''
443            phasedict = self.Phases[phasenam] # pointer to current phase info
444            General = phasedict['General']
445            Z = General.get('cellZ',0.0)
446            cx,ct,cs,cia = General['AtomPtrs']
447            Atoms = phasedict['Atoms']
448            fpfx = str(phasedict['pId'])+'::Afrac:'       
449            cfrac = cx+3
450            cmult = cs+1
451            compDict = {} # combines H,D & T
452            sitemultlist = []
453            massDict = dict(zip(General['AtomTypes'],General['AtomMass']))
454            cellmass = 0
455            for i,at in enumerate(Atoms):
456                atype = at[ct].strip()
457                if atype.find('-') != -1: atype = atype.split('-')[0]
458                if atype.find('+') != -1: atype = atype.split('+')[0]
459                atype = atype[0].upper()+atype[1:2].lower() # force case conversion
460                if atype == "D" or atype == "D": atype = "H"
461                fvar = fpfx+str(i)
462                fval = self.parmDict.get(fvar,at[cfrac])
463                mult = at[cmult]
464                if not massDict.get(at[ct]):
465                    print 'No mass found for atom type '+at[ct]
466                    print 'Will not compute cell contents for phase '+phasenam
467                    return
468                cellmass += massDict[at[ct]]*mult*fval
469                compDict[atype] = compDict.get(atype,0.0) + mult*fval
470                if fval == 1: sitemultlist.append(mult)
471            if len(compDict.keys()) == 0: return # no elements!
472            if Z < 1: # Z has not been computed or set by user
473                Z = 1
474                for i in range(2,min(sitemultlist)+1):
475                    for m in sitemultlist:
476                        if m % i != 0:
477                            break
478                        else:
479                            Z = i
480                General['cellZ'] = Z # save it
481
482            # when scattering factors are included in the CIF, this needs to be
483            # added to the loop here but only in the one-block case.
484            # For multiblock CIFs, scattering factors go in the histogram
485            # blocks  (for all atoms in all appropriate phases)
486
487            #if oneblock: # add scattering factors for current phase here
488            WriteCIFitem('\nloop_  _atom_type_symbol _atom_type_number_in_cell')
489            formula = ''
490            reload(G2mth)
491            for elem in HillSortElements(compDict.keys()):
492                WriteCIFitem('  ' + PutInCol(elem,4) +
493                             G2mth.ValEsd(compDict[elem],-0.009,True))
494                if formula: formula += " "
495                formula += elem
496                if compDict[elem] == Z: continue
497                formula += G2mth.ValEsd(compDict[elem]/Z,-0.009,True)
498            WriteCIFitem( '\n# Note that Z affects _cell_formula_sum and _weight')
499            WriteCIFitem( '_cell_formula_units_Z',str(Z))
500            WriteCIFitem( '_chemical_formula_sum',formula)
501            WriteCIFitem( '_chemical_formula_weight',
502                          G2mth.ValEsd(cellmass/Z,-0.09,True))
503
504        def WriteDistances(phasenam,SymOpList,offsetList,symOpList,G2oprList):
505            '''Report bond distances and angles for the CIF
506
507            Note that _geom_*_symmetry_* fields are values of form
508            n_klm where n is the symmetry operation in SymOpList (counted
509            starting with 1) and (k-5, l-5, m-5) are translations to add
510            to (x,y,z). See
511            http://www.iucr.org/__data/iucr/cifdic_html/1/cif_core.dic/Igeom_angle_site_symmetry_.html
512
513            TODO: need a method to select publication flags for distances/angles
514            '''
515            phasedict = self.Phases[phasenam] # pointer to current phase info           
516            Atoms = phasedict['Atoms']
517            cx,ct,cs,cia = phasedict['General']['AtomPtrs']
518            fpfx = str(phasedict['pId'])+'::Afrac:'       
519            cfrac = cx+3
520            # loop over interatomic distances for this phase
521            WriteCIFitem('\n# MOLECULAR GEOMETRY')
522            WriteCIFitem('loop_' + 
523                         '\n\t_geom_bond_atom_site_label_1' +
524                         '\n\t_geom_bond_atom_site_label_2' + 
525                         '\n\t_geom_bond_distance' + 
526                         '\n\t_geom_bond_site_symmetry_1' + 
527                         '\n\t_geom_bond_site_symmetry_2' + 
528                         '\n\t_geom_bond_publ_flag')
529
530            # Note that labels should be read from self.labellist to correspond
531            # to the values reported in the atoms table and zero occupancy atoms
532            # should not be included
533            fpfx = str(phasedict['pId'])+'::Afrac:'       
534            for i,at in enumerate(Atoms):
535                if self.parmDict.get(fpfx+str(i),at[cfrac]) == 0.0: continue
536                lbl = self.labellist[i]
537
538
539            # loop over interatomic angles for this phase
540            WriteCIFitem('loop_' + 
541                         '\n\t_geom_angle_atom_site_label_1' + 
542                         '\n\t_geom_angle_atom_site_label_2' + 
543                         '\n\t_geom_angle_atom_site_label_3' + 
544                         '\n\t_geom_angle' + 
545                         '\n\t_geom_angle_site_symmetry_1' +
546                         '\n\t_geom_angle_site_symmetry_2' + 
547                         '\n\t_geom_angle_site_symmetry_3' + 
548                         '\n\t_geom_angle_publ_flag')
549
550
551        def WritePhaseInfo(phasenam):
552            WriteCIFitem('\n# phase info for '+str(phasenam) + ' follows')
553            phasedict = self.Phases[phasenam] # pointer to current phase info           
554            WriteCIFitem('_pd_phase_name', phasenam)
555            pfx = str(phasedict['pId'])+'::'
556            A,sigA = G2stIO.cellFill(pfx,phasedict['General']['SGData'],self.parmDict,self.sigDict)
557            cellSig = G2stIO.getCellEsd(pfx,
558                                       phasedict['General']['SGData'],A,
559                                       self.OverallParms['Covariance'])  # returns 7 vals, includes sigVol
560            cellList = G2lat.A2cell(A) + (G2lat.calc_V(A),)
561            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
562            names = ['length_a','length_b','length_c',
563                     'angle_alpha','angle_beta ','angle_gamma',
564                     'volume']
565            prevsig = 0
566            for lbl,defsig,val,sig in zip(names,defsigL,cellList,cellSig):
567                if sig:
568                    txt = G2mth.ValEsd(val,sig)
569                    prevsig = -sig # use this as the significance for next value
570                else:
571                    txt = G2mth.ValEsd(val,min(defsig,prevsig),True)
572                WriteCIFitem('_cell_'+lbl,txt)
573                   
574            WriteCIFitem('_symmetry_cell_setting',
575                         phasedict['General']['SGData']['SGSys'])
576
577            spacegroup = phasedict['General']['SGData']['SpGrp'].strip()
578            # regularize capitalization and remove trailing H/R
579            spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
580            WriteCIFitem('_symmetry_space_group_name_H-M',spacegroup)
581
582            # generate symmetry operations including centering and center of symmetry
583            SymOpList,offsetList,symOpList,G2oprList = G2spg.AllOps(
584                phasedict['General']['SGData'])
585            WriteCIFitem('loop_ _space_group_symop_id _space_group_symop_operation_xyz')
586            for i,op in enumerate(SymOpList,start=1):
587                WriteCIFitem('   {:3d}  {:}'.format(i,op.lower()))
588
589            # loop over histogram(s) used in this phase
590            if not oneblock and not self.quickmode:
591                # report pointers to the histograms used in this phase
592                histlist = []
593                for hist in self.Phases[phasenam]['Histograms']:
594                    if self.Phases[phasenam]['Histograms'][hist]['Use']:
595                        if phasebyhistDict.get(hist):
596                            phasebyhistDict[hist].append(phasenam)
597                        else:
598                            phasebyhistDict[hist] = [phasenam,]
599                        blockid = datablockidDict.get(hist)
600                        if not blockid:
601                            print "Internal error: no block for data. Phase "+str(
602                                phasenam)+" histogram "+str(hist)
603                            histlist = []
604                            break
605                        histlist.append(blockid)
606
607                if len(histlist) == 0:
608                    WriteCIFitem('# Note: phase has no associated data')
609                else:
610                    WriteCIFitem('loop_  _pd_block_diffractogram_id')
611
612            # report atom params
613            if phasedict['General']['Type'] == 'nuclear':        #this needs macromolecular variant, etc!
614                WriteAtomsNuclear(phasenam)
615            else:
616                raise Exception,"no export for mm coordinates implemented"
617            # report cell contents
618            WriteComposition(phasenam)
619            if not self.quickmode:      # report distances and angles
620                WriteDistances(phasenam,SymOpList,offsetList,symOpList,G2oprList)
621
622        def WritePowderData(histlbl):
623            text = '?'
624            histblk = self.Histograms[histlbl]
625            inst = histblk['Instrument Parameters'][0]
626            hId = histblk['hId']
627            pfx = ':' + str(hId) + ':'
628            print 'TODO: powder here data for',histblk["Sample Parameters"]['InstrName']
629            # see wrpowdhist.for & wrreflist.for
630           
631            if 'Lam1' in inst:
632                ratio = self.parmDict.get('I(L2)/I(L1)',inst['I(L2)/I(L1)'][1])
633                sratio = self.sigDict.get('I(L2)/I(L1)',-0.0009)
634                lam1 = self.parmDict.get('Lam1',inst['Lam1'][1])
635                slam1 = self.sigDict.get('Lam1',-0.00009)
636                lam2 = self.parmDict.get('Lam2',inst['Lam2'][1])
637                slam2 = self.sigDict.get('Lam2',-0.00009)
638                # always assume Ka1 & Ka2 if two wavelengths are present
639                WriteCIFitem('loop_' + 
640                             '\n\t_diffrn_radiation_wavelength' +
641                             '\n\t_diffrn_radiation_wavelength_wt' + 
642                             '\n\t_diffrn_radiation_type' + 
643                             '\n\t_diffrn_radiation_wavelength_id')
644                WriteCIFitem('  ' + PutInCol(G2mth.ValEsd(lam1,slam1),15)+
645                             PutInCol('1.0',15) + 
646                             PutInCol('K\\a~1~',10) + 
647                             PutInCol('1',5))
648                WriteCIFitem('  ' + PutInCol(G2mth.ValEsd(lam2,slam2),15)+
649                             PutInCol(G2mth.ValEsd(ratio,sratio),15)+
650                             PutInCol('K\\a~2~',10) + 
651                             PutInCol('2',5))               
652            else:
653                lam1 = self.parmDict.get('Lam',inst['Lam'][1])
654                slam1 = self.sigDict.get('Lam',-0.00009)
655                WriteCIFitem('_diffrn_radiation_wavelength',G2mth.ValEsd(lam1,slam1))
656
657
658            if not oneblock:
659                if not phasebyhistDict.get(histlbl):
660                    WriteCIFitem('\n# No phases associated with this data set')
661                else:
662                    WriteCIFitem('\n# PHASE TABLE')
663                    WriteCIFitem('loop_' +
664                                 '\n\t_pd_phase_id' + 
665                                 '\n\t_pd_phase_block_id' + 
666                                 '\n\t_pd_phase_mass_%')
667                    wtFrSum = 0.
668                    for phasenam in phasebyhistDict.get(histlbl):
669                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
670                        General = self.Phases[phasenam]['General']
671                        wtFrSum += hapData['Scale'][0]*General['Mass']
672
673                    for phasenam in phasebyhistDict.get(histlbl):
674                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
675                        General = self.Phases[phasenam]['General']
676                        wtFr = hapData['Scale'][0]*General['Mass']/wtFrSum
677                        pfx = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':'
678                        if pfx+'Scale' in self.sigDict:
679                            sig = self.sigDict[pfx+'Scale']*wtFr/hapData['Scale'][0]
680                        else:
681                            sig = -0.0001
682                        WriteCIFitem(
683                            '  '+
684                            str(self.Phases[phasenam]['pId']) +
685                            '  '+datablockidDict[phasenam]+
686                            '  '+G2mth.ValEsd(wtFr,sig)
687                            )
688
689            # TODO: this will need help from Bob
690            # WriteCIFitem('_pd_proc_ls_prof_R_factor','?')
691            # WriteCIFitem('_pd_proc_ls_prof_wR_factor','?')
692            # WriteCIFitem('_pd_proc_ls_prof_wR_expected','?')
693            # WriteCIFitem('_refine_ls_R_Fsqd_factor','?')
694
695            if histblk['Instrument Parameters'][0]['Type'][1][1] == 'X':
696                WriteCIFitem('_diffrn_radiation_probe','x-ray')
697                pola = histblk['Instrument Parameters'][0].get('Polariz.')
698                if pola:
699                    pfx = ':' + str(hId) + ':'
700                    sig = self.sigDict.get(pfx+'Polariz.',-0.0009)
701                    txt = G2mth.ValEsd(pola[1],sig)
702                    WriteCIFitem('_diffrn_radiation_polarisn_ratio',txt)
703            elif histblk['Instrument Parameters'][0]['Type'][1][1] == 'N':
704                WriteCIFitem('_diffrn_radiation_probe','neutron')
705            # TOF (note that this may not be defined)
706            #if histblk['Instrument Parameters'][0]['Type'][1][2] == 'T':
707            #    WriteCIFitem('_pd_meas_2theta_fixed',text)
708           
709
710            # TODO: this will need help from Bob
711            #if not oneblock:
712            #WriteCIFitem('\n# SCATTERING FACTOR INFO')
713            #WriteCIFitem('loop_  _atom_type_symbol')
714            #if histblk['Instrument Parameters'][0]['Type'][1][1] == 'X':
715            #    WriteCIFitem('      _atom_type_scat_dispersion_real')
716            #    WriteCIFitem('      _atom_type_scat_dispersion_imag')
717            #    for lbl in ('a1','a2','a3', 'a4', 'b1', 'b2', 'b3', 'b4', 'c'):
718            #        WriteCIFitem('      _atom_type_scat_Cromer_Mann_'+lbl)
719            #elif histblk['Instrument Parameters'][0]['Type'][1][1] == 'N':
720            #    WriteCIFitem('      _atom_type_scat_length_neutron')
721            #WriteCIFitem('      _atom_type_scat_source')
722
723            WriteCIFitem('_pd_proc_ls_background_function',FormatBackground(histblk['Background']))
724
725            #WriteCIFitem('_exptl_absorpt_process_details','?')
726            #WriteCIFitem('_exptl_absorpt_correction_T_min','?')
727            #WriteCIFitem('_exptl_absorpt_correction_T_max','?')
728            #C extinction
729            #WRITE(IUCIF,'(A)') '# Extinction correction'
730            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_min',TEXT(1:10))
731            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_max',TEXT(11:20))
732
733            if not oneblock:                 # instrumental profile terms go here
734                WriteCIFitem('_pd_proc_ls_profile_function', 
735                    FormatInstProfile(histblk["Instrument Parameters"]))
736
737            #refprx = '_refln.' # mm
738            refprx = '_refln_' # normal
739            WriteCIFitem('\n# STRUCTURE FACTOR TABLE')           
740            # compute maximum intensity reflection
741            Imax = 0
742            for phasenam in histblk['Reflection Lists']:
743                scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
744                Icorr = np.array([refl[13] for refl in histblk['Reflection Lists'][phasenam]])
745                FO2 = np.array([refl[8] for refl in histblk['Reflection Lists'][phasenam]])
746                I100 = scale*FO2*Icorr
747                Imax = max(Imax,max(I100))
748
749            WriteCIFitem('loop_')
750            if len(histblk['Reflection Lists'].keys()) > 1:
751                WriteCIFitem('\t_pd_refln_phase_id')
752            WriteCIFitem('\t' + refprx + 'index_h' + 
753                         '\n\t' + refprx + 'index_k' + 
754                         '\n\t' + refprx + 'index_l' + 
755#                         '\n\t_pd_refln_wavelength_id' +
756#                         '\n\t' + refprx + 'status' +
757                         '\n\t' + refprx + 'F_squared_meas' + 
758#                         '\n\t' + refprx + 'F_squared_sigma' +
759                         '\n\t' + refprx + 'F_squared_calc' + 
760                         '\n\t' + refprx + 'phase_calc' + 
761                         '\n\t_pd_refln_d_spacing')
762            if Imax > 0:
763                WriteCIFitem('\t_gsas_i100_meas')
764
765            refcount = 0
766            hklmin = None
767            hklmax = None
768            dmax = None
769            dmin = None
770            for phasenam in histblk['Reflection Lists']:
771                scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
772                phaseid = self.Phases[phasenam]['pId']
773                refcount += len(histblk['Reflection Lists'][phasenam])
774                for ref in histblk['Reflection Lists'][phasenam]:
775                    if DEBUG:
776                        print 'DEBUG: skip reflection list'
777                        break
778                    if hklmin is None:
779                        hklmin = ref[0:3]
780                        hklmax = ref[0:3]
781                        dmax = dmin = ref[4]
782                    if len(histblk['Reflection Lists'].keys()) > 1:
783                        s = PutInCol(phaseid,2)
784                    else:
785                        s = ""
786                    for i,hkl in enumerate(ref[0:3]):
787                        hklmax[i] = max(hkl,hklmax[i])
788                        hklmin[i] = min(hkl,hklmin[i])
789                        s += PutInCol(int(hkl),4)
790                    for I in ref[8:10]:
791                        s += PutInCol(G2mth.ValEsd(I,-0.0009),10)
792                    s += PutInCol(G2mth.ValEsd(ref[10],-0.9),7)
793                    dmax = max(dmax,ref[4])
794                    dmin = min(dmin,ref[4])
795                    s += PutInCol(G2mth.ValEsd(ref[4],-0.009),8)
796                    if Imax > 0:
797                        I100 = 100.*scale*ref[8]*ref[13]/Imax
798                        s += PutInCol(G2mth.ValEsd(I100,-0.09),6)
799                    WriteCIFitem("  "+s)
800
801            WriteCIFitem('_reflns_number_total', str(refcount))
802            if hklmin is not None and len(histblk['Reflection Lists']) == 1: # hkl range has no meaning with multiple phases
803                WriteCIFitem('_reflns_limit_h_min', str(int(hklmin[0])))
804                WriteCIFitem('_reflns_limit_h_max', str(int(hklmax[0])))
805                WriteCIFitem('_reflns_limit_k_min', str(int(hklmin[1])))
806                WriteCIFitem('_reflns_limit_k_max', str(int(hklmax[1])))
807                WriteCIFitem('_reflns_limit_l_min', str(int(hklmin[2])))
808                WriteCIFitem('_reflns_limit_l_max', str(int(hklmax[2])))
809            if hklmin is not None:
810                WriteCIFitem('_reflns_d_resolution_high', G2mth.ValEsd(dmin,-0.009))
811                WriteCIFitem('_reflns_d_resolution_low', G2mth.ValEsd(dmax,-0.0009))
812
813            WriteCIFitem('\n# POWDER DATA TABLE')
814            # is data fixed step? If the step varies by <0.01% treat as fixed step
815            steps = histblk['Data'][0][1:] - histblk['Data'][0][:-1]
816            if abs(max(steps)-min(steps)) > abs(max(steps))/10000.:
817                fixedstep = False
818            else:
819                fixedstep = True
820
821            if fixedstep: # and not TOF
822                WriteCIFitem('_pd_meas_2theta_range_min', G2mth.ValEsd(histblk['Data'][0][0],-0.00009))
823                WriteCIFitem('_pd_meas_2theta_range_max', G2mth.ValEsd(histblk['Data'][0][-1],-0.00009))
824                WriteCIFitem('_pd_meas_2theta_range_inc', G2mth.ValEsd(steps.sum()/len(steps),-0.00009))
825                # zero correct, if defined
826                zero = None
827                zerolst = histblk['Instrument Parameters'][0].get('Zero')
828                if zerolst: zero = zerolst[1]
829                zero = self.parmDict.get('Zero',zero)
830                # TODO: Bob is zero added or subtracted?
831                if zero:
832                    WriteCIFitem('_pd_proc_2theta_range_min', G2mth.ValEsd(histblk['Data'][0][0]-zero,-0.00009))
833                    WriteCIFitem('_pd_proc_2theta_range_max', G2mth.ValEsd(histblk['Data'][0][-1]-zero,-0.00009))
834                    WriteCIFitem('_pd_proc_2theta_range_inc', G2mth.ValEsd(steps.sum()/len(steps),-0.00009))
835               
836            if zero:
837                WriteCIFitem('_pd_proc_number_of_points', str(len(histblk['Data'][0])))
838            else:
839                WriteCIFitem('_pd_meas_number_of_points', str(len(histblk['Data'][0])))
840            WriteCIFitem('\nloop_')
841            #            WriteCIFitem('\t_pd_proc_d_spacing') # need easy way to get this
842            if not fixedstep:
843                if zero:
844                    WriteCIFitem('\t_pd_proc_2theta_corrected')
845                else:
846                    WriteCIFitem('\t_pd_meas_2theta_scan')
847            # at least for now, always report weights. TODO: Should they be multiplied by weights?
848            #if countsdata:
849            #    WriteCIFitem('\t_pd_meas_counts_total')
850            #else:
851            WriteCIFitem('\t_pd_meas_intensity_total')
852            WriteCIFitem('\t_pd_calc_intensity_total')
853            WriteCIFitem('\t_pd_proc_intensity_bkg_calc')
854            WriteCIFitem('\t_pd_proc_ls_weight')
855            # TODO: are data excluded?
856            for x,yobs,yw,ycalc,ybkg in zip(histblk['Data'][0],
857                                            histblk['Data'][1],
858                                            histblk['Data'][2],
859                                            histblk['Data'][3],
860                                            histblk['Data'][4]):
861                if DEBUG:
862                    print 'DEBUG: skip reflection list'
863                    break
864                if fixedstep:
865                    s = ""
866                else:
867                    s = PutInCol(G2mth.ValEsd(x-zero,-0.00009),10)
868                s += PutInCol(G2mth.ValEsd(yobs,-0.00009),14)
869                s += PutInCol(G2mth.ValEsd(ycalc,-0.00009),14)
870                s += PutInCol(G2mth.ValEsd(ybkg,-0.00009),14)
871                s += PutInCol(G2mth.ValEsd(yw,-0.000009),14)
872                WriteCIFitem("  "+s)
873                         
874        def WriteSingleXtalData(histlbl):
875            histblk = self.Histograms[histlbl]
876            print 'TODO: single xtal here data for',histblk["Instrument Parameters"][0]['InstrName']
877
878            #refprx = '_refln.' # mm
879            refprx = '_refln_' # normal
880
881            print histblk.keys()
882           
883            WriteCIFitem('\n# STRUCTURE FACTOR TABLE')           
884            WriteCIFitem('loop_' + 
885                         '\n\t' + refprx + 'index_h' + 
886                         '\n\t' + refprx + 'index_k' + 
887                         '\n\t' + refprx + 'index_l' +
888                         '\n\t' + refprx + 'F_squared_meas' + 
889                         '\n\t' + refprx + 'F_squared_sigma' + 
890                         '\n\t' + refprx + 'F_squared_calc' + 
891                         '\n\t' + refprx + 'phase_calc' +
892                         '\n\t_pd_refln_d_spacing'
893                         )
894
895            hklmin = None
896            hklmax = None
897            dmax = None
898            dmin = None
899            refcount = len(histblk['Data'])
900            for ref in histblk['Data']:
901                s = "  "
902                if hklmin is None:
903                    hklmin = ref[0:3]
904                    hklmax = ref[0:3]
905                    dmax = dmin = ref[4]
906                for i,hkl in enumerate(ref[0:3]):
907                    hklmax[i] = max(hkl,hklmax[i])
908                    hklmin[i] = min(hkl,hklmin[i])
909                    s += PutInCol(int(hkl),4)
910                sig = ref[6] * ref[8] / ref[5]
911                s += PutInCol(G2mth.ValEsd(ref[8],-abs(sig/10)),12)
912                s += PutInCol(G2mth.ValEsd(sig,-abs(sig)/10.),10)
913                s += PutInCol(G2mth.ValEsd(ref[9],-abs(sig/10)),12)
914                s += PutInCol(G2mth.ValEsd(ref[10],-0.9),7)
915                dmax = max(dmax,ref[4])
916                dmin = min(dmin,ref[4])
917                s += PutInCol(G2mth.ValEsd(ref[4],-0.009),8)
918                WriteCIFitem(s)
919            WriteCIFitem('_reflns_number_total', str(refcount))
920            if hklmin is not None:
921                WriteCIFitem('_reflns_limit_h_min', str(int(hklmin[0])))
922                WriteCIFitem('_reflns_limit_h_max', str(int(hklmax[0])))
923                WriteCIFitem('_reflns_limit_k_min', str(int(hklmin[1])))
924                WriteCIFitem('_reflns_limit_k_max', str(int(hklmax[1])))
925                WriteCIFitem('_reflns_limit_l_min', str(int(hklmin[2])))
926                WriteCIFitem('_reflns_limit_l_max', str(int(hklmax[2])))
927                WriteCIFitem('_reflns_d_resolution_high', G2mth.ValEsd(dmin,-0.009))
928                WriteCIFitem('_reflns_d_resolution_low', G2mth.ValEsd(dmax,-0.0009))
929
930        #============================================================
931        # the export process starts here
932        # also load all of the tree into a set of dicts
933        self.loadTree()
934        #self.dumpTree()
935        # create a dict with refined values and their uncertainties
936        self.loadParmDict()
937        #
938
939        # get restraint info
940        #restraintDict = self.OverallParms.get('Restraints',{})
941        #for i in  self.OverallParms['Constraints']:
942        #    print i
943        #    for j in self.OverallParms['Constraints'][i]:
944        #        print j
945        #return
946
947        self.CIFdate = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%dT%H:%M")
948        # count phases, powder and single crystal histograms
949        self.nphase = len(self.Phases)
950        self.powderDict = {}
951        self.xtalDict = {}
952        for hist in self.Histograms:
953            i = self.Histograms[hist]['hId']
954            if hist.startswith("PWDR"): 
955                self.powderDict[i] = hist
956            elif hist.startswith("HKLF"): 
957                self.xtalDict[i] = hist
958        # is there anything to export?
959        if self.nphase + len(self.powderDict) + len(self.xtalDict) == 0: 
960            self.G2frame.ErrorDialog(
961                'Empty project',
962                'No data or phases to include in CIF')
963            return
964        # is there a file name defined?
965        self.CIFname = os.path.splitext(
966            os.path.split(self.G2frame.GSASprojectfile)[1]
967            )[0]
968        self.CIFname = self.CIFname.replace(' ','')
969        if not self.CIFname:
970            self.G2frame.ErrorDialog(
971                'No GPX name',
972                'Please save the project to provide a name')
973            return
974        # test for quick CIF mode or no data
975        self.quickmode = False
976        phasenam = phasenum = None # include all phases
977        if mode != "full" or len(self.powderDict) + len(self.xtalDict) == 0:
978            self.quickmode = True
979            oneblock = True
980            if self.nphase == 0:
981                self.G2frame.ErrorDialog(
982                    'No phase present',
983                    'Cannot create a coordinates CIF with no phases')
984                return
985            elif self.nphase > 1: # quick mode: choose one phase
986                choices = sorted(self.Phases.keys())
987                phasenum = G2gd.ItemSelector(choices,self.G2frame)
988                if phasenum is None: return
989                phasenam = choices[phasenum]
990        # will this require a multiblock CIF?
991        elif self.nphase > 1:
992            oneblock = False
993        elif len(self.powderDict) + len(self.xtalDict) > 1:
994            oneblock = False
995        else: # one phase, one dataset, Full CIF
996            oneblock = True
997
998        # make sure needed infomation is present
999        # get CIF author name -- required for full CIFs
1000        try:
1001            self.author = self.OverallParms['Controls'].get("Author",'').strip()
1002        except KeyError:
1003            pass
1004        while not (self.author or self.quickmode):
1005            dlg = G2gd.SingleStringDialog(self.G2frame,'Get CIF Author','Provide CIF Author name (Last, First)')
1006            if not dlg.Show(): return # cancel was pressed
1007            self.author = dlg.GetValue()
1008            dlg.Destroy()
1009        try:
1010            self.OverallParms['Controls']["Author"] = self.author # save for future
1011        except KeyError:
1012            pass
1013        self.shortauthorname = self.author.replace(',','').replace(' ','')[:20]
1014
1015        # check the instrument name for every histogram
1016        if not self.quickmode:
1017            dictlist = []
1018            keylist = []
1019            lbllist = []
1020            invalid = 0
1021            key3 = 'InstrName'
1022            for hist in self.Histograms:
1023                if hist.startswith("PWDR"): 
1024                    key2 = "Sample Parameters"
1025                    d = self.Histograms[hist][key2]
1026                elif hist.startswith("HKLF"): 
1027                    key2 = "Instrument Parameters"
1028                    d = self.Histograms[hist][key2][0]
1029                   
1030                lbllist.append(hist)
1031                dictlist.append(d)
1032                keylist.append(key3)
1033                instrname = d.get(key3)
1034                if instrname is None:
1035                    d[key3] = ''
1036                    invalid += 1
1037                elif instrname.strip() == '':
1038                    invalid += 1
1039            if invalid:
1040                msg = ""
1041                if invalid > 3: msg = (
1042                    "\n\nNote: it may be faster to set the name for\n"
1043                    "one histogram for each instrument and use the\n"
1044                    "File/Copy option to duplicate the name"
1045                    )
1046                if not G2gd.CallScrolledMultiEditor(
1047                    self.G2frame,dictlist,keylist,
1048                    prelbl=range(1,len(dictlist)+1),
1049                    postlbl=lbllist,
1050                    title='Instrument names',
1051                    header="Edit instrument names. Note that a non-blank\nname is required for all histograms"+msg,
1052                    ): return
1053
1054        if oneblock and not self.quickmode:
1055            # select a dataset to use (there should only be one set in one block,
1056            # but take whatever comes 1st
1057            for hist in self.Histograms:
1058                histblk = self.Histograms[hist]
1059                if hist.startswith("PWDR"): 
1060                    instnam = histblk["Sample Parameters"]['InstrName']
1061                    break # ignore all but 1st data histogram
1062                elif hist.startswith("HKLF"): 
1063                    instnam = histblk["Instrument Parameters"][0]['InstrName']
1064                    break # ignore all but 1st data histogram
1065        #======================================================================
1066        # Start writing the CIF - single block
1067        #======================================================================
1068        if oneblock:
1069            WriteCIFitem('data_'+self.CIFname)
1070            if phasenam is None: # if not already selected, select the first phase (should be one)
1071                phasenam = self.Phases.keys()[0]
1072            #print 'phasenam',phasenam
1073            phaseblk = self.Phases[phasenam] # pointer to current phase info
1074            if not self.quickmode:
1075                instnam = instnam.replace(' ','')
1076                WriteCIFitem('_pd_block_id',
1077                             str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
1078                             str(self.shortauthorname) + "|" + instnam)
1079                WriteAudit()
1080                WritePubTemplate()
1081                WriteOverall()
1082                WritePhaseTemplate()
1083            # report the phase info
1084            WritePhaseInfo(phasenam)
1085            if hist.startswith("PWDR") and not self.quickmode:
1086                # preferred orientation
1087                SH = FormatSH(phasenam)
1088                MD = FormatHAPpo(phasenam)
1089                if SH and MD:
1090                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', SH + '\n' + MD)
1091                elif SH or MD:
1092                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', SH + MD)
1093                else:
1094                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', 'none')
1095                    # report profile, since one-block: include both histogram and phase info
1096                WriteCIFitem('_pd_proc_ls_profile_function',
1097                             FormatInstProfile(histblk["Instrument Parameters"])
1098                             + '\n' +
1099                             FormatPhaseProfile(phasenam))
1100                WritePowderTemplate()
1101                WritePowderData(hist)
1102            elif hist.startswith("HKLF") and not self.quickmode:
1103                WriteSnglXtalTemplate()
1104                WriteSingleXtalData(hist)
1105        else:
1106        #======================================================================
1107        # Start writing the CIF - multiblock
1108        #======================================================================
1109            # publication info
1110            WriteCIFitem('\ndata_'+self.CIFname+'_publ')
1111            WriteAudit()
1112            WriteCIFitem('_pd_block_id',
1113                         str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
1114                         str(self.shortauthorname) + "|Overall")
1115            WritePubTemplate()
1116            # overall info
1117            WriteCIFitem('data_'+str(self.CIFname)+'_overall')
1118            WriteOverall()
1119            #============================================================
1120            WriteCIFitem('# POINTERS TO PHASE AND HISTOGRAM BLOCKS')
1121            datablockidDict = {} # save block names here -- N.B. check for conflicts between phase & hist names (unlikely!)
1122            # loop over phase blocks
1123            if self.nphase > 1:
1124                loopprefix = ''
1125                WriteCIFitem('loop_   _pd_phase_block_id')
1126            else:
1127                loopprefix = '_pd_phase_block_id'
1128           
1129            for phasenam in sorted(self.Phases.keys()):
1130                i = self.Phases[phasenam]['pId']
1131                datablockidDict[phasenam] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
1132                             'phase_'+ str(i) + '|' + str(self.shortauthorname))
1133                WriteCIFitem(loopprefix,datablockidDict[phasenam])
1134            # loop over data blocks
1135            if len(self.powderDict) + len(self.xtalDict) > 1:
1136                loopprefix = ''
1137                WriteCIFitem('loop_   _pd_block_diffractogram_id')
1138            else:
1139                loopprefix = '_pd_block_diffractogram_id'
1140            for i in sorted(self.powderDict.keys()):
1141                hist = self.powderDict[i]
1142                histblk = self.Histograms[hist]
1143                instnam = histblk["Sample Parameters"]['InstrName']
1144                instnam = instnam.replace(' ','')
1145                i = histblk['hId']
1146                datablockidDict[hist] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
1147                                         str(self.shortauthorname) + "|" +
1148                                         instnam + "_hist_"+str(i))
1149                WriteCIFitem(loopprefix,datablockidDict[hist])
1150            for i in sorted(self.xtalDict.keys()):
1151                hist = self.xtalDict[i]
1152                histblk = self.Histograms[hist]
1153                instnam = histblk["Instrument Parameters"][0]['InstrName']
1154                instnam = instnam.replace(' ','')
1155                i = histblk['hId']
1156                datablockidDict[hist] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
1157                                         str(self.shortauthorname) + "|" +
1158                                         instnam + "_hist_"+str(i))
1159                WriteCIFitem(loopprefix,datablockidDict[hist])
1160            #============================================================
1161            # loop over phases, exporting them
1162            phasebyhistDict = {} # create a cross-reference to phases by histogram
1163            for j,phasenam in enumerate(sorted(self.Phases.keys())):
1164                i = self.Phases[phasenam]['pId']
1165                WriteCIFitem('\ndata_'+self.CIFname+"_phase_"+str(i))
1166                print "debug, processing ",phasenam
1167                WriteCIFitem('# Information for phase '+str(i))
1168                WriteCIFitem('_pd_block_id',datablockidDict[phasenam])
1169                # report the phase
1170                WritePhaseTemplate()
1171                WritePhaseInfo(phasenam)
1172                # preferred orientation
1173                SH = FormatSH(phasenam)
1174                MD = FormatHAPpo(phasenam)
1175                if SH and MD:
1176                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', SH + '\n' + MD)
1177                elif SH or MD:
1178                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', SH + MD)
1179                else:
1180                    WriteCIFitem('_pd_proc_ls_pref_orient_corr', 'none')
1181                # report sample profile terms
1182                PP = FormatPhaseProfile(phasenam)
1183                if PP:
1184                    WriteCIFitem('_pd_proc_ls_profile_function',PP)
1185                   
1186            #============================================================
1187            # loop over histograms, exporting them
1188            for i in sorted(self.powderDict.keys()):
1189                hist = self.powderDict[i]
1190                histblk = self.Histograms[hist]
1191                if hist.startswith("PWDR"): 
1192                    WriteCIFitem('\ndata_'+self.CIFname+"_pwd_"+str(i))
1193                    #instnam = histblk["Sample Parameters"]['InstrName']
1194                    # report instrumental profile terms
1195                    WriteCIFitem('_pd_proc_ls_profile_function',
1196                        FormatInstProfile(histblk["Instrument Parameters"]))
1197                    WriteCIFitem('# Information for histogram '+str(i)+': '+hist)
1198                    WriteCIFitem('_pd_block_id',datablockidDict[hist])
1199                    WritePowderTemplate()
1200                    WritePowderData(hist)
1201            for i in sorted(self.xtalDict.keys()):
1202                hist = self.xtalDict[i]
1203                histblk = self.Histograms[hist]
1204                if hist.startswith("HKLF"): 
1205                    WriteCIFitem('\ndata_'+self.CIFname+"_sx_"+str(i))
1206                    #instnam = histblk["Instrument Parameters"][0]['InstrName']
1207                    WriteCIFitem('# Information for histogram '+str(i)+': '+hist)
1208                    WriteCIFitem('_pd_block_id',datablockidDict[hist])
1209                    WriteSnglXtalTemplate()
1210                    WriteSingleXtalData(hist)
1211
1212        WriteCIFitem('#--' + 15*'eof--' + '#')
1213
Note: See TracBrowser for help on using the repository browser.